Functions
median_contours.C File Reference
#include "CAFAna/Analysis/Calcs.h"
#include "CAFAna/Analysis/Exposures.h"
#include "CAFAna/Fit/Fit.h"
#include "CAFAna/Analysis/Plots.h"
#include "3FlavorAna/Plotting/NuePlotStyle.h"
#include "CAFAna/Fit/FrequentistSurface.h"
#include "CAFAna/Analysis/Style.h"
#include "CAFAna/Experiment/Dmsq32Constraint.h"
#include "CAFAna/Experiment/MultiExperiment.h"
#include "CAFAna/Experiment/ReactorExperiment.h"
#include "CAFAna/Experiment/SingleSampleExperiment.h"
#include "CAFAna/FC/FCSurface.h"
#include "3FlavorAna/Prediction/PredictionSystJoint2018.h"
#include "CAFAna/Prediction/PredictionCombinePeriods.h"
#include "3FlavorAna/Systs/NumuSysts.h"
#include "CAFAna/Vars/FitVars.h"
#include "Utilities/rootlogon.C"
#include "OscLib/IOscCalc.h"
#include "3FlavorAna/Ana2020/joint_fit_2020_loader_tools.h"
#include "TCanvas.h"
#include "TMarker.h"
#include "TBox.h"
#include "TLatex.h"
#include "TColor.h"
#include "TGraph.h"
#include "TVectorD.h"
#include "TF1.h"
#include "TLegend.h"
#include "TLine.h"
#include "TStyle.h"
#include "TSystem.h"
#include "TGaxis.h"
#include <algorithm>
#include <vector>
#include <string>
#include <fstream>

Go to the source code of this file.

Functions

void DrawHieTag (int hie, bool th13=true, bool vert=false)
 
double sigma (TH1F *hist, double percentile)
 
bool is_file_exist (const char *fileName)
 
void median_contours (int universes=100, bool corrSysts=false, string folder="2020bf", TString options="joint_realData_both_onlyIH", bool dmsqSurf=false, bool th13Surf=false, bool zoomIn=true, bool fillContour=false)
 

Function Documentation

void DrawHieTag ( int  hie,
bool  th13 = true,
bool  vert = false 
)

Definition at line 458 of file median_contours.C.

References ana::kBlessedLabelFontSize, ana::kPrimColorIH, and ana::kPrimColorNH.

Referenced by median_contours().

458  {
459  if(hie>0){
460  TLatex * ltxNH = new TLatex(0.16,0.89,"#bf{NH}");
461  ltxNH->SetTextColor(kPrimColorNH);
462  ltxNH->SetNDC();
463  ltxNH->SetTextAlign(22);
464  ltxNH->SetTextSize(kBlessedLabelFontSize);
465  if(th13) ltxNH->Draw();
466  else if(!vert) ltxNH->DrawLatex(.85,.20,"#bf{NH}");
467  else ltxNH->DrawLatex(.85,.85,"#bf{NH}");
468  }
469  else{
470  TLatex * ltxIH = new TLatex (0.16,0.46,"#bf{IH}");
471  ltxIH->SetTextColor(kPrimColorIH);
472  ltxIH->SetNDC();
473  ltxIH->SetTextAlign(22);
474  ltxIH->SetTextSize(kBlessedLabelFontSize);
475  if(th13) ltxIH->Draw();
476  else if(!vert) ltxIH->DrawLatex(.85,.20,"#bf{IH}");
477  else ltxIH->DrawLatex(.85,.85,"#bf{IH}");
478  }
479 }
ana::kPrimColorIH
const Color_t kPrimColorIH
Definition: Style.h:64
ana::kPrimColorNH
const Color_t kPrimColorNH
Definition: Style.h:61
ana::kBlessedLabelFontSize
const Float_t kBlessedLabelFontSize
Definition: Style.h:90
th13
double th13
Definition: runWimpSim.h:98
bool is_file_exist ( const char *  fileName)

Definition at line 53 of file median_contours.C.

References om::cout, allTimeWatchdog::endl, and infile.

Referenced by median_contours().

54 {
55  std::ifstream infile(fileName);
56  std::cout<<fileName<<std::endl;
57  return infile.good();
58 }
plotROC.fileName
fileName
Definition: plotROC.py:78
infile
string infile
Definition: INukeNucleonCorr.cxx:60
om::cout
OStream cout
Definition: OStream.cxx:6
void median_contours ( int  universes = 100,
bool  corrSysts = false,
string  folder = "2020bf",
TString  options = "joint_realData_both_onlyIH",
bool  dmsqSurf = false,
bool  th13Surf = false,
bool  zoomIn = true,
bool  fillContour = false 
)

contour "minus sigma"

contour "minus sigma"

contour "plus sigma"

Definition at line 60 of file median_contours.C.

References ana::assert(), om::cout, gen_hdf5record::debug, DrawContourCanvas(), DrawHieTag(), APDHVSetting::dummy, allTimeWatchdog::endl, PandAna.Demos.tute_pid_validation::folder, MECModelEnuComparisons::g, ana::Gaussian2Sigma2D(), ana::Gaussian3Sigma2D(), ana::Gaussian68Percent2D(), MECModelEnuComparisons::i, infile, is_file_exist(), calib::j, ana::JoinGraphs(), ana::k2SigmaConfidenceColorIH, ana::k2SigmaConfidenceColorNH, ana::kBlessedLabelFontSize, kViolet, MECModelEnuComparisons::leg, ana::FrequentistSurface::LoadFrom(), m, median(), Munits::ms, fabricate::options, outdir, make_template_knob_config::outfilename, Preliminary(), nd_projection_maker::ps, save(), sigma(), PandAna.Demos.tute_pid_validation::specs, demo5::surf, tmp, art::to_string(), and ana::UniqueName().

68 {
69 
70  bool debug = 1;
71  bool nueOnly = options.Contains("nueOnly");
72  bool numuOnly = options.Contains("numuOnly");
73  bool joint = options.Contains("joint");
74  assert (nueOnly || numuOnly || joint);
75 
76  bool FHCOnly = options.Contains("FHCOnly");
77  bool RHCOnly = options.Contains("RHCOnly");
78  bool both = options.Contains("both");
79  assert (FHCOnly || RHCOnly || both);
80 
81  bool fake2017 = options.Contains("fake2017");
82  bool realData = options.Contains("realData");
83 
84  bool onlyNH = options.Contains("onlyNH");
85  bool onlyIH = options.Contains("onlyIH");
86 
87  auto suffix = options;
88  if(dmsqSurf) suffix+= "_dmsq";
89  if(th13Surf) suffix+= "_th13";
90  if(corrSysts) suffix+="_systs";
91  else suffix+="_stats10";
92  string syst = (corrSysts?"syst":"stat");
93 
94  TString outdir="./";
95 
96  //////////////////////////////////////////////////////////////////////
97  // Make plots
98  //////////////////////////////////////////////////////////////////////
99  // load old analyses contours for comparisons
100  TFile * official =0;
101  string outdir2018 = "";
102  if(th13Surf) outdir2018 = "th13_delta/";
103  else if (dmsqSurf) outdir2018 = "th23_dmsq/";
104  else outdir2018 = "delta_th23/";
105  TString ana2018 ("/nova/ana/nu_e_ana/Ana2018/Results/RHC_and_FHC/contours/" + outdir2018 + "stat/hist_contours_2018_" + options + "combo_"+(dmsqSurf?"dmsq_":"")+"stats");
106  TFile * infile2018ana = new TFile (ana2018 + ".root", "read");
107  auto mins2018 =* (TVectorD*)infile2018ana->Get("overall_minchi");
108  double minchi23_2018ana = mins2018[0];
109 
110  TString ana2019 ("/nova/ana/nu_e_ana/Ana2019/Results/RHC_and_FHC/contours/" + outdir2018 + "stat/hist_contours_2019_" + options + "_"+(dmsqSurf?"dmsq_":"")+"stats");
111  TFile * infile2019ana = new TFile (ana2019 + ".root", "read");
112  auto mins2019 =* (TVectorD*)infile2019ana->Get("overall_minchi");
113  double minchi23_2019ana = mins2019[0];
114 
115  TString ana2020 ("/nova/ana/3flavor/Ana2020/Fits/RealData/root_contours/hist_contours_2020_" + options + "_"+(dmsqSurf?"dmsq_":"")+"stats");
116  TFile * infile2020ana = new TFile (ana2020 + ".root", "read");
117  auto mins2020 =* (TVectorD*)infile2020ana->Get("overall_minchi");
118  double minchi23_2020ana = mins2020[0];
119 
120  // load universes
121  std::vector <TString> surfNames;
122  if(!th13Surf && !dmsqSurf)surfNames = {"delta_th23"};
123  if (dmsqSurf) surfNames.push_back("th23_dm32");
124  if (th13Surf) surfNames.push_back("th13_delta");
125  for(TString surfName:surfNames){
126  for (int hie:{-1,+1}){
127  if((onlyNH && hie<0) || (onlyIH && hie>0)) continue;
128 
129  auto surf_2018ana = *FrequentistSurface::LoadFrom(infile2018ana, "surf_"+surfName+(hie>0?"_NH":"_IH"));
130  TH2 * surf_1Sigma_2018ana = Gaussian68Percent2D(surf_2018ana);
131  TH2 * surf_2Sigma_2018ana = Gaussian2Sigma2D(surf_2018ana);
132  TH2 * surf_3Sigma_2018ana = Gaussian3Sigma2D(surf_2018ana);
133  auto g1_2018ana = surf_2018ana.GetGraphs(surf_1Sigma_2018ana, minchi23_2018ana);
134  auto g2_2018ana = surf_2018ana.GetGraphs(surf_2Sigma_2018ana, minchi23_2018ana);
135  auto g3_2018ana = surf_2018ana.GetGraphs(surf_3Sigma_2018ana, minchi23_2018ana);
136 
137  auto surf_2019ana = *FrequentistSurface::LoadFrom(infile2019ana, "surf_"+surfName+(hie>0?"_NH":"_IH"));
138  TH2 * surf_1Sigma_2019ana = Gaussian68Percent2D(surf_2019ana);
139  TH2 * surf_2Sigma_2019ana = Gaussian2Sigma2D(surf_2019ana);
140  TH2 * surf_3Sigma_2019ana = Gaussian3Sigma2D(surf_2019ana);
141  auto g1_2019ana = surf_2019ana.GetGraphs(surf_1Sigma_2019ana, minchi23_2019ana);
142  auto g2_2019ana = surf_2019ana.GetGraphs(surf_2Sigma_2019ana, minchi23_2019ana);
143  auto g3_2019ana = surf_2019ana.GetGraphs(surf_3Sigma_2019ana, minchi23_2019ana);
144 
145  auto surf_2020ana = *FrequentistSurface::LoadFrom(infile2020ana, "surf_"+surfName+(hie>0?"_NH":"_IH"));
146  TH2 * surf_1Sigma_2020ana = Gaussian68Percent2D(surf_2020ana);
147  TH2 * surf_2Sigma_2020ana = Gaussian2Sigma2D(surf_2020ana);
148  TH2 * surf_3Sigma_2020ana = Gaussian3Sigma2D(surf_2020ana);
149  auto g1_2020ana = surf_2020ana.GetGraphs(surf_1Sigma_2020ana, minchi23_2020ana);
150  auto g2_2020ana = surf_2020ana.GetGraphs(surf_2Sigma_2020ana, minchi23_2020ana);
151  auto g3_2020ana = surf_2020ana.GetGraphs(surf_3Sigma_2020ana, minchi23_2020ana);
152 
153  TColor *c_coral2 = new TColor(102, 255/255.0, 169/255.0, 155/255.0);
154  TColor *c_coral = new TColor(100, 255/255.0, 228/255.0, 218/255.0);
155  TColor *c_wave = new TColor(101, 0/255.0, 167/255.0, 200/255.0);
156  Color_t kFillTestColor = ((hie>0)? 101:102);
157  Int_t kFillColorCont = TColor::GetColorTransparent(kFillTestColor, 0.3);
158  if (onlyIH) kFillColorCont = 100;
159  Color_t k2SigmaConfidenceColor = ((hie>0)? k2SigmaConfidenceColorNH:k2SigmaConfidenceColorIH);
160  Int_t kFillColor;
161  Int_t kFillColor1 = TColor::GetColorTransparent(k2SigmaConfidenceColor, 0.75);
162  Int_t kFillColor2 = TColor::GetColorTransparent(k2SigmaConfidenceColor, 0.38);
163  Int_t kFillColor3 = TColor::GetColorTransparent(k2SigmaConfidenceColor, 0.18);
164  Int_t kDarkColor = kGray+3;
165 
166  int thr = 10; // submit jobs with 10 in one
167  std::vector <TH2 *> hhist; // all contours
168  std::vector <double> chis; // all minima
169 
170  for(int j=0; j< universes; j++){
171  TString outfilename ("/nova/ana/users/lkolupae/mock/contours"+ folder+ "/" + syst + "_hist_contours_" + suffix + "_thr"+std::to_string(j) + ".root");
172  if(is_file_exist(outfilename)){
173  TFile * infile = new TFile (outfilename,"read");
174  for (int k=0; k<10; k++){
175  cout<<"univ "<<j<<endl;
176  auto mins =* (TVectorD*)infile->Get(("universe_"+std::to_string(k)+"/bf_"+std::to_string(k)).c_str());
177  double minchi23 = mins[0];
178  cout<<minchi23<<endl;
179  auto surf = *FrequentistSurface::LoadFrom(infile, "universe_"+std::to_string(k)+"/surf_"+surfName+"_univ"+std::to_string(k)+"_"+(hie>0?"NH":"IH"));
180  hhist.push_back(surf.ToTH2(minchi23));
181  chis.push_back(minchi23);
182  }
183  infile->Close();
184  }
185  }
186  TH2D * median;
187  TH2D * plsigma;
188  TH2D * minsigma;
189  // fill plus, minus and median histograms
190  if (surfName=="delta_th23"){
191  median = new TH2D("me", "", 30, 0-1/30.0, 2+1/30.0, 30, 0.3-0.2/30, 0.7+0.2/30);// plsigma, minsigma;
192  plsigma = new TH2D("p", "", 30, 0-1/30.0, 2+1/30.0, 30, 0.3-0.2/30, 0.7+0.2/30);
193  minsigma = new TH2D("mi", "", 30, 0-1/30.0, 2+1/30.0, 30, 0.3-0.2/30, 0.7+0.2/30);}
194  if (surfName=="th23_dm32" && onlyNH){
195  median = new TH2D("me", "", 20, 0.35-0.35/20, 0.7+0.35/20, 20, 2.1-0.6/20.0, 2.7+0.6/20.0);// plsigma, minsigma;
196  plsigma = new TH2D("p", "", 20, 0.35-0.35/20, 0.7+0.35/20, 20, 2.1-0.6/20.0, 2.7+0.6/20.0);
197  minsigma = new TH2D("mi", "", 20, 0.35-0.35/20, 0.7+0.35/20, 20, 2.1-0.6/20.0, 2.7+0.6/20.0);}
198  if (surfName=="th23_dm32" && onlyIH){
199  median = new TH2D("me", "", 20, 0.35-0.35/20, 0.7+0.35/20, 20, -2.7-0.5/20.0, -2.2+0.5/20.0);// plsigma, minsigma;
200  plsigma = new TH2D("p", "", 20, 0.35-0.35/20, 0.7+0.35/20, 20, -2.7-0.5/20.0, -2.2+0.5/20.0);
201  minsigma = new TH2D("mi", "", 20, 0.35-0.35/20, 0.7+0.35/20, 20, -2.7-0.5/20.0, -2.2+0.5/20.0);}
202  // find plus, minus and median values for each bin
203  for (int i = 0; i<=31; i++){
204  for (int k = 0; k<=31; k++){
205  // fill 1d histogram with all LL values for each universe
206  TH1F* m = new TH1F("tmp", "tmp", 500, 0, 200);
207  for(int j=0; j< hhist.size(); j++){
208  auto tmp = hhist[j]->GetBinContent(i, k);
209  m->Fill(tmp);
210  }
211  auto med = sigma(m, 0.5);
212  median->SetBinContent(i, k, med);
213  auto ps = sigma(m, 0.84);
214  plsigma->SetBinContent(i, k, ps);
215  auto ms = sigma(m, 0.16);
216  minsigma->SetBinContent(i, k, ms);
217  cout<<"m "<<med<<" ps "<<ps<<" min "<<ms<<endl;
218  }
219  }
220 
221  cout<<"test "<<median->GetMinimum()<<" "<<plsigma->GetMinimum()<<" "<<minsigma->GetMinimum()<<endl;
222  new TCanvas;
223  median->Draw("colz");
224  cout<<"minimum "<< median->GetMinimum()<<endl;
225  gPad->Print(("median_"+suffix + "_" + folder + ".pdf"));
226 
227  new TCanvas;
228  plsigma->Draw("colz");
229  cout<<"minimum "<< plsigma->GetMinimum()<<endl;
230  gPad->Print(("plsigma_"+suffix + "_" + folder + ".pdf"));
231 
232  new TCanvas;
233  minsigma->Draw("colz");
234  cout<<"minimum "<< minsigma->GetMinimum()<<endl;
235  gPad->Print(("minsigma_"+suffix + "_" + folder +".pdf"));
236  // for plus, minus and median surfaces find 1,2,3 sigma contours
237  for(auto &sigma:{1,2,3}){
238 
239  if(surfName.Contains("delta_th23")) {
240  if (zoomIn) DrawContourCanvas(surfName, 0., 2., (RHCOnly?0.15:0.275), (RHCOnly?0.8:0.725));
241  else DrawContourCanvas(surfName, 0., 2., 0., 1.);
242  }
243  if(surfName.Contains("th23_dm32")) {
244  if (hie>0) DrawContourCanvas(surfName, ((nueOnly||RHCOnly)?0.2:0.3), ((nueOnly||RHCOnly)?0.8:0.7), (RHCOnly?1.8:2.05), (RHCOnly?3.3:2.85));
245  else DrawContourCanvas(surfName, ((nueOnly||RHCOnly)?0.2:0.3), ((nueOnly||RHCOnly)?0.8:0.7), (RHCOnly?-3.5:-2.95), (RHCOnly?-2.0:-2.25));
246  }
247  if(surfName.Contains("th13_delta")) {
248  if (zoomIn) DrawContourCanvas(surfName, 0., 0.38, 0, 2);
249  else DrawContourCanvas(surfName, 0., 1.0, 0, 2);
250  }
251 
252  kFillColor = TColor::GetColorTransparent(k2SigmaConfidenceColor, 0.75);
253 
254  double value = 0;
255  if (sigma == 1) value = 2.3;
256  if (sigma == 2) value = 6.18;
257  if (sigma == 3) value = 11.83;
258 
259  TVirtualPad* save = gPad;
260  TCanvas tmp;
261  median->SetContour(1, &value);
262  median->Draw("CONT Z LIST");
263  tmp.Update();
264  tmp.Paint();
265  tmp.Print(("median_cont_" + folder + ".pdf").c_str());
266  std::vector<TGraph*> gr_m;
267  TCollection* specs = gROOT->GetListOfSpecials();
268  TIter nextSpec(specs);
269  while(TObject* spec = nextSpec()){
270  if(!spec->InheritsFrom(TObjArray::Class())) continue;
271  TObjArray* conts = (TObjArray*)spec;
272  if(conts->IsEmpty()) continue;
273  if(!conts->At(0)->InheritsFrom(TList::Class())) continue;
274  TList* cont = (TList*)conts->At(0);
275  TIter nextObj(cont);
276  while(TObject* obj = nextObj()){
277  if(!obj->InheritsFrom(TGraph::Class())) continue;
278  gr_m.push_back((TGraph*)obj->Clone(UniqueName().c_str()));
279  }
280  }
281 
282  TCanvas tmp2;
283  plsigma->SetContour(1, &value);
284  plsigma->Draw("CONT Z LIST");
285  tmp2.Update();
286  tmp2.Paint();
287  tmp2.Print(("plsigma_cont_" + folder + ".pdf").c_str());
288  std::vector<TGraph*> gr_p;
289  TCollection* specs_p = gROOT->GetListOfSpecials();
290  TIter nextSpec_p(specs_p);
291  while(TObject* spec_p = nextSpec_p()){
292  if(!spec_p->InheritsFrom(TObjArray::Class())) continue;
293  TObjArray* conts_p = (TObjArray*)spec_p;
294  if(conts_p->IsEmpty()) continue;
295  if(!conts_p->At(0)->InheritsFrom(TList::Class())) continue;
296  TList* cont_p = (TList*)conts_p->At(0);
297  TIter nextObj_p(cont_p);
298  while(TObject* obj_p = nextObj_p()){
299  if(!obj_p->InheritsFrom(TGraph::Class())) continue;
300  gr_p.push_back((TGraph*)obj_p->Clone(UniqueName().c_str()));
301  }
302  }
303 
304  TCanvas tmp3;
305  minsigma->SetContour(1, &value);
306  minsigma->Draw("CONT Z LIST");
307  tmp3.Update();
308  tmp3.Paint();
309  tmp3.Print(("minsigma_cont_" + folder + ".pdf").c_str());
310  std::vector<TGraph*> gr_min;
311  TCollection* specs_min = gROOT->GetListOfSpecials();
312  TIter nextSpec_min(specs_min);
313  while(TObject* spec_min = nextSpec_min()){
314  if(!spec_min->InheritsFrom(TObjArray::Class())) continue;
315  TObjArray* conts = (TObjArray*)spec_min;
316  if(conts->IsEmpty()) continue;
317  if(!conts->At(0)->InheritsFrom(TList::Class())) continue;
318  TList* cont = (TList*)conts->At(0);
319  TIter nextObj(cont);
320  while(TObject* obj = nextObj()){
321  if(!obj->InheritsFrom(TGraph::Class())) continue;
322  gr_min.push_back((TGraph*)obj->Clone(UniqueName().c_str()));
323  }
324  }
325  gPad=save;
326  cout<<"size of gr "<<gr_m.size()<<" p "<<gr_p.size()<<" min "<<gr_min.size()<<endl;
327 
328  if (folder == "2020max"){
329  if(surfName.Contains("delta_th23") && onlyNH && sigma > 2) { ///contour "minus sigma"
330  JoinGraphs(gr_min[0], gr_min[1], kFillColorCont)->Draw("f");
331  for(TGraph* g: gr_min){
332  g->SetLineWidth(3);
333  g->SetLineStyle(kSolid);
334  g->SetLineColor(kFillColorCont);
335  g->Draw("l");
336  }
337  }
338  else {
339  for(TGraph* g: gr_min){
340  g->SetLineWidth(3);
341  g->SetLineStyle(kSolid);
342  g->SetLineColor(kFillColorCont);
343  g->SetFillColor(kFillColorCont);
344  g->Draw("f");
345  g->Draw("l");
346  }
347  }
348  for(TGraph* g: gr_p){
349  g->SetLineWidth(3);
350  g->SetLineStyle(kSolid);
351  g->SetLineColor(kFillColorCont);
352  g->SetFillColor(kWhite);
353  g->Draw("f");
354  g->Draw("l");
355  }
356  }
357 
358  if (folder == "2020bf"){
359  if(surfName.Contains("delta_th23") && sigma >3 && onlyIH || surfName.Contains("delta_th23") && onlyNH) { ///contour "minus sigma"
360  JoinGraphs(gr_min[0], gr_min[1], kFillColorCont)->Draw("f");
361  for(TGraph* g: gr_min){
362  g->SetLineWidth(3);
363  g->SetLineStyle(kSolid);
364  g->SetLineColor(kFillColorCont);
365  g->Draw("l");
366  }
367  }
368  else {
369  for(TGraph* g: gr_min){
370  g->SetLineWidth(3);
371  g->SetLineStyle(kSolid);
372  g->SetLineColor(kFillColorCont);
373  g->SetFillColor(kFillColorCont);
374  g->Draw("f");
375  g->Draw("l");
376  }
377  }
378  if (surfName.Contains("delta_th23") && sigma >2 && onlyNH){ ///contour "plus sigma"
379  if(sigma == 3) JoinGraphs(gr_p[0], gr_p[1], kWhite)->Draw("f");
380  for(TGraph* g: gr_p){
381  g->SetLineWidth(3);
382  g->SetLineStyle(kSolid);
383  g->SetLineColor(kFillColorCont);
384  g->Draw("l");
385  }
386  }
387  else {
388  for(TGraph* g: gr_p){
389  g->SetLineWidth(3);
390  g->SetLineStyle(kSolid);
391  g->SetLineColor(kFillColorCont);
392  g->SetFillColor(kWhite);
393  g->Draw("f");
394  g->Draw("l");
395  }
396  }
397  }
398 
399  for(TGraph* g: gr_m){
400  g->SetLineWidth(3);
401  g->SetLineStyle(kSolid);
402  g->SetLineColor(kFillColor);
403  g->Draw("l");
404  }
405 
406  auto gr2018 = g1_2018ana;
407  auto gr2019 = g1_2019ana;
408  auto gr2020 = g1_2020ana;
409  if (sigma==1) cout<<"2018 size "<<gr2018.size()<< "2019 size "<<gr2019.size()<<" 2020 size "<<gr2020.size()<<endl;
410  if (sigma==2) {gr2018 = g2_2018ana; gr2019 = g2_2019ana; gr2020 = g2_2020ana; cout<<"2018 size "<<gr2018.size()<<" 2019 size "<<gr2019.size() << " 2020 size "<<gr2020.size() <<endl;}
411  if (sigma==3) {gr2018 = g3_2018ana; gr2019 = g3_2019ana; gr2020 = g3_2020ana; cout<<"2018 size "<<gr2018.size()<<" 2019 size "<<gr2019.size() << " 2020 size "<<gr2020.size() << endl; }
412 
413  for (auto & g:gr2018) { g->SetLineColor(kGray); g->Draw("l");}
414  for (auto & g:gr2019) { g->SetLineColor(kGray+2); g->Draw("l");}
415  for (auto & g:gr2020) { g->SetLineColor(kViolet-6); g->Draw("l");}
416 
417  TLegend * leg = new TLegend(0.16,0.17,0.8,0.28);
418  leg->SetTextSize(kBlessedLabelFontSize*0.75);
419  leg->SetFillStyle(0);
420  leg->SetMargin(1.3*leg->GetMargin());
421  leg->SetNColumns(3);
422  leg->SetMargin(1.4*leg->GetMargin());
423  TGraph * dummy = new TGraph;
424  dummy->SetLineWidth(3);
425  dummy->SetLineColor(kFillColor);
426  leg->AddEntry(dummy->Clone(),"median","l");
427  dummy->SetFillColor(kFillColorCont);
428  leg->AddEntry(dummy->Clone(),"68\% exp","f");
429  leg->SetHeader((std::to_string(sigma) + " sigma contour").c_str());
430 
431  TLegend * leg2 = new TLegend(0.16,0.8,0.9,0.85);
432  leg2->SetTextSize(kBlessedLabelFontSize*0.75);
433  leg2->SetNColumns(3);
434  dummy->SetLineColor(kGray);
435  leg2->AddEntry(dummy->Clone(),"2018 stat.fit","l");
436  dummy->SetLineColor(kGray+2);
437  leg2->AddEntry(dummy->Clone(),"2019 stat fit","l");
438  dummy->SetLineColor(kViolet-6);
439  leg2->AddEntry(dummy->Clone(),"2020 stat fit","l");
440  leg->Draw();
441  leg2->Draw();
442 
443  Preliminary();
444  DrawHieTag(hie,false);
445  gPad->RedrawAxis();
446 
447  gPad->Print(outdir + "contour_" + suffix + "_" + folder + "_" +
448  surfName + "_" + std::to_string(universes) + "_sigma" + std::to_string(sigma) + "_"+
449  (hie>0?"NH":"IH") +
450  (zoomIn?"_zoom":"") +
451  ".pdf");
452  }//end sigma
453  }//end hie
454  }//end surf
455 }//end all
Munits::ms
static constexpr Double_t ms
Definition: Munits.h:192
ana::Gaussian68Percent2D
TH2 * Gaussian68Percent2D(const FrequentistSurface &s)
Up-value surface for 68% confidence in 2D in gaussian approximation.
Definition: FrequentistSurface.cxx:485
gen_hdf5record.debug
bool debug
Definition: gen_hdf5record.py:290
ana::k2SigmaConfidenceColorNH
const Color_t k2SigmaConfidenceColorNH
Definition: Style.h:77
ana::Gaussian2Sigma2D
TH2 * Gaussian2Sigma2D(const FrequentistSurface &s)
Up-value surface for 2 sigma confidence in 2D in gaussian approximation.
Definition: FrequentistSurface.cxx:488
tmp
Float_t tmp
Definition: plot.C:36
ana::LoadFrom
std::unique_ptr< T > LoadFrom(TDirectory *dir, const std::string &label)
Definition: LoadFromFile.h:17
m
#define m
Definition: MakeGENIEReweightTable_module.cc:421
make_template_knob_config.outfilename
string outfilename
knobs that need extra care
Definition: make_template_knob_config.py:120
infile
string infile
Definition: INukeNucleonCorr.cxx:60
value
const XML_Char int const XML_Char * value
Definition: expat.h:331
DrawHieTag
void DrawHieTag(int hie, bool th13=true, bool vert=false)
Definition: median_contours.C:458
ana::Gaussian3Sigma2D
TH2 * Gaussian3Sigma2D(const FrequentistSurface &s)
Up-value surface for 3 sigma confidence in 2D in gaussian approximation.
Definition: FrequentistSurface.cxx:490
calib::j
const double j
Definition: BetheBloch.cxx:29
median
double median(TH1D *h)
Definition: absCal.cxx:524
save
void save(bool IsData, bool Bad)
Definition: cellShifts.C:453
Preliminary
void Preliminary()
Definition: draw_decomp_plots.C:55
sigma
double sigma(TH1F *hist, double percentile)
Definition: median_contours.C:46
DrawContourCanvas
void DrawContourCanvas(TString surfName, double minx=0, double maxx=2, double miny=0, double maxy=1)
Definition: joint_fit_2017_contours.C:51
ana::JoinGraphs
TGraph * JoinGraphs(TGraph *a, TGraph *b, int fillcol)
Join graphs and set a fill color. Useful for contours.
Definition: Plots.cxx:1531
om::cout
OStream cout
Definition: OStream.cxx:6
PandAna.Demos.tute_pid_validation.folder
folder
Definition: tute_pid_validation.py:180
kViolet
enum BeamMode kViolet
ana::assert
assert(nhit_max >=nhit_nbins)
PandAna.Demos.tute_pid_validation.specs
dictionary specs
Definition: tute_pid_validation.py:144
ana::kBlessedLabelFontSize
const Float_t kBlessedLabelFontSize
Definition: Style.h:90
nd_projection_maker.ps
list ps
Definition: nd_projection_maker.py:5
art::to_string
std::string to_string(ModuleType const mt)
Definition: ModuleType.h:34
ana::UniqueName
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:28
outdir
const std::string outdir
Definition: nus_error_barchart.C:10
ana::k2SigmaConfidenceColorIH
const Color_t k2SigmaConfidenceColorIH
Definition: Style.h:82
is_file_exist
bool is_file_exist(const char *fileName)
Definition: median_contours.C:53
demo5.surf
surf
Definition: demo5.py:35
double sigma ( TH1F *  hist,
double  percentile 
)

Definition at line 46 of file median_contours.C.

Referenced by ana::SystMakerLoaderShift::AddOnOff(), ana::SystMakerWeightShift::AddOnOff(), ana::SystMakerLoaderShift::AddSigma(), ana::SystMakerWeightShift::AddSigma(), ana::SystMakerSystShift::AddSigma(), ana::BeamIntensitySyst::BeamIntensitySyst(), ana::BeamNueScaleSyst::BeamNueScaleSyst(), rwgt::MakeGENIEReweightTable::beginJob(), ana::best_spline_find_mode(), ana::best_spline_mode_fcn(), ana::BolombMuESyst2019::BolombMuESyst2019(), CAFReweightExample(), ana::CalcChi2(), check_predinterp(), check_predinterp_numu(), ana::CompNormSyst::CompNormSyst(), jw::CompNormSyst::CompNormSyst(), ana::ContainmentSyst::ContainmentSyst(), ana::CorrMuEScaleSyst2019::CorrMuEScaleSyst2019(), ana::CorrMuEScaleSyst2020::CorrMuEScaleSyst2020(), d3sigma_calc(), d4sigma_calc(), DCMOffsetCalculator(), demo_nueNumuSysts(), demoStarPlots(), ana::DirectHadEScaleSyst2017::DirectHadEScaleSyst2017(), ana::DirectRelHadEScaleSyst2017::DirectRelHadEScaleSyst2017(), ana::DISHighWScaleSyst::DISHighWScaleSyst(), ana::SystMakerShift::DrawBkg(), ana::SystMakerShift::DrawShift(), ana::SystMakerShift::DrawSig(), rsim::ExcessNoiseMaker::DrawSmearedPE(), calib::DCMTimingOffset::endJob(), calib::DCMTimingOffset::endRun(), FCContour(), ana::FDHadEnergyScaleSyst::FDHadEnergyScaleSyst(), ana::FDMuEnergyScaleSyst::FDMuEnergyScaleSyst(), hough::Gaus(), ana::GenericSystComponentScale< T >::GetCut(), ana::SystMakerShift::GetLabel(), ana::getNue2018PredVar(), ana::SystMakerShift::GetShiftedSpectrum(), SystGroupHelper::GetShifts(), ana::SystMaker::GetShifts(), GetShifts2017(), ana::GetShifts2018(), ana::GetShifts2020(), rsim::ReadoutSim::GetSmearedPE(), ana::SystMaker::GetTotalShift(), ana::HadEnergyScaleSyst::HadEnergyScaleSyst(), HighlightSignPoints(), ana::KaonScaleSyst::KaonScaleSyst(), ana::ISyst::LatexName(), ana::LEMScaleSyst::LEMScaleSyst(), ana::LeptonAngleSyst::LeptonAngleSyst(), novarwgt::LinearInterp(), ana::ReactorExperiment::LoadFrom(), ana::SystMakerShift::LoadFrom(), ana::LoadMaps(), main(), make_nue_xsec_pred(), make_starplots(), make_time_profile(), ana::MakeLatexCommandName(), MakeName(), ana::MakeNueSystematicsFile(), MakeNus18Systs(), MakeNus2019CovMxSysts(), MakeSysts(), rsim::ExcessNoiseMaker::MakeTemplate(), ana::MECDoubleGaussEnhSyst::MECDoubleGaussEnhSyst(), ana::MECDoubleGaussEnhSystNux::MECDoubleGaussEnhSystNux(), median_contours(), ana::MichelTaggingSyst2018::MichelTaggingSyst2018(), ana::MichelTaggingSyst2020::MichelTaggingSyst2020(), ana::MuEnergyScaleSyst::MuEnergyScaleSyst(), ana::MuEScaleSyst2017::MuEScaleSyst2017(), ana::MyWrongSignScale::MyWrongSignScale(), ana::NCScaleSyst::NCScaleSyst(), ana::NCtest::NCtest(), ana::NeutronVisEScaleSyst2018::NeutronVisEScaleSyst2018(), ana::NormSyst::NormSyst(), ana::covmx::NormSyst::NormSyst(), ana::NOvARwgtSyst::NOvARwgtSyst(), ana::NueAcceptSystBkg2018FHC::NueAcceptSystBkg2018FHC(), ana::NueAcceptSystBkg2018RHC::NueAcceptSystBkg2018RHC(), ana::NueAcceptSystSignalKin2018FHC::NueAcceptSystSignalKin2018FHC(), ana::NueAcceptSystSignalKin2018RHC::NueAcceptSystSignalKin2018RHC(), ana::NueAcceptSystSignalKin2020FHC::NueAcceptSystSignalKin2020FHC(), ana::NueAcceptSystSignalKin2020RHC::NueAcceptSystSignalKin2020RHC(), ana::NueExtrapSystBkg2017::NueExtrapSystBkg2017(), ana::NueExtrapSystSignalKin2017::NueExtrapSystSignalKin2017(), ana::Nuetest::Nuetest(), ana::NumuCCScaleSyst::NumuCCScaleSyst(), ana::NumuEnergyScaleSyst::NumuEnergyScaleSyst(), ana::NumuGEANTNormSyst::NumuGEANTNormSyst(), ana::NumuGEANTScaleSyst::NumuGEANTScaleSyst(), ana::NumuNCScaleSyst::NumuNCScaleSyst(), ana::NumuNormSyst::NumuNormSyst(), ana::NumuRelEnergyScaleSyst::NumuRelEnergyScaleSyst(), ana::NumuRelNormSyst::NumuRelNormSyst(), ana::NumuTauContaminationSyst::NumuTauContaminationSyst(), ana::Numutest::Numutest(), Nus18SystsBeamTranspLoad(), Nus18SystsGENIELoad(), Nus18SystsNeutronLoad(), cmf::ParameterSetToLocation(), ana::PileupMuESyst2019::PileupMuESyst2019(), ana::PileupMuESyst2020::PileupMuESyst2020(), pion_energy_estimator(), predEventCountWithSystError(), ana::SingleSampleExperiment::PredHistIncCosmics(), ana::PredictionInterp::PredictionInterp(), ana::PredictionSyst3Flavor2020::PredictionSyst3Flavor2020(), ana::PredictionSystJoint2018::PredictionSystJoint2018(), ana::PredictionSystJointDemo::PredictionSystJointDemo(), ana::PredictionSystNue2017::PredictionSystNue2017(), ana::PredictionSystNumu2017::PredictionSystNumu2017(), ana::SystMaker::ProcessSyst(), calib::TimingCalibration::produce(), ana::ReinteractionSyst::ReinteractionSyst(), ana::RelMuEScaleSyst2017::RelMuEScaleSyst2017(), ana::RelNormSyst::RelNormSyst(), ana::RemIDScaleSyst::RemIDScaleSyst(), ana::RockMuonNormSyst::RockMuonNormSyst(), novaddt::UpMuAna::RootFitLLR(), runCheatDecomp(), ana::NuISyst::SampleName(), ana::NuISyst::SaveTo(), ana::SystMakerShiftSigma::SaveTo(), ana::SigmaDelta::SetFitVars(), earms::ElasticArms::SetHit(), ana::WrongSignContamination::Shift(), ana::WrongSignScale::Shift(), ana::NueExtrapSystBkg2017::Shift(), ana::NumuGEANTNormSyst::Shift(), ana::NueAcceptSystBkg2018FHC::Shift(), ana::NCtest::Shift(), ana::MuEScaleSyst2017::Shift(), ana::WrongSignEnergyBinRHC::Shift(), ana::HadEnergyScaleSyst::Shift(), ana::NormSyst::Shift(), ana::CorrMuEScaleSyst2020::Shift(), ana::NumuGEANTScaleSyst::Shift(), ana::CorrMuEScaleSyst2019::Shift(), ana::Nuetest::Shift(), ana::WrongSignEnergyBin::Shift(), ana::NumuSummedSmallGENIESyst::Shift(), ana::Numutest::Shift(), ana::RelMuEScaleSyst2017::Shift(), ana::NCScaleSyst::Shift(), ana::WrongSignScale100::Shift(), ana::MuEnergyScaleSyst::Shift(), ana::UnCorrFDMuEScaleSyst2019::Shift(), ana::UnCorrFDMuEScaleSyst2020::Shift(), ana::NueAcceptSystBkg2018RHC::Shift(), ana::BeamNueScaleSyst::Shift(), ana::NumuNCScaleSyst::Shift(), ana::NumuCCScaleSyst::Shift(), ana::DirectHadEScaleSyst2017::Shift(), ana::MyWrongSignScale::Shift(), ana::NumuTauContaminationSyst::Shift(), ana::UnCorrNDMuEScaleSyst2020::Shift(), ana::KaonScaleSyst::Shift(), ana::UnCorrNDMuEScaleSyst2019::Shift(), ana::NumuNormSyst::Shift(), ana::UnCorrMuCatMuESyst2020::Shift(), ana::NumuRelNormSyst::Shift(), ana::UnCorrMuCatMuESyst2019::Shift(), ana::LEMScaleSyst::Shift(), ana::RemIDScaleSyst::Shift(), ana::NumuEnergyScaleSyst::Shift(), ana::PileupMuESyst2020::Shift(), ana::PileupMuESyst2019::Shift(), ana::ContainmentSyst::Shift(), ana::DISHighWScaleSyst::Shift(), ana::CompNormSyst::Shift(), jw::CompNormSyst::Shift(), jw::MECInitStateNPFracShift::Shift(), cmf::ShifterAndWeighter::ShiftForParameter(), starPlot(), syst_test(), SystsBeamLoad(), SystsBeamTranspLoad17(), SystsGENIELoad(), SystsGENIELoad17(), demo::DemoSyst1::TruthShift(), ana::UnCorrFDMuEScaleSyst2019::UnCorrFDMuEScaleSyst2019(), ana::UnCorrFDMuEScaleSyst2020::UnCorrFDMuEScaleSyst2020(), ana::UnCorrMuCatMuESyst2019::UnCorrMuCatMuESyst2019(), ana::UnCorrMuCatMuESyst2020::UnCorrMuCatMuESyst2020(), ana::UnCorrNDMuEScaleSyst2019::UnCorrNDMuEScaleSyst2019(), ana::UnCorrNDMuEScaleSyst2020::UnCorrNDMuEScaleSyst2020(), ana::WrongSignContamination::WrongSignContamination(), ana::WrongSignEnergyBin::WrongSignEnergyBin(), ana::WrongSignEnergyBinRHC::WrongSignEnergyBinRHC(), ana::WrongSignScale::WrongSignScale(), ana::WrongSignScale100::WrongSignScale100(), xsec_extrap_plots(), and ana::SummedSyst::~SummedSyst().

46  {
47  double xq[1] = {percentile};
48  double yq[1];
49  hist->GetQuantiles(1, yq, xq);
50  return yq[0];
51 }
Generated on Wed Apr 14 2021 14:11:47 for NOvA by   doxygen 1.8.11