draw_plots_util.h
Go to the documentation of this file.
6 #include "CAFAna/Core/Binning.h"
7 #include "CAFAna/Core/HistAxis.h"
16 #include "CAFAna/Vars/FitVars.h"
17 #include "CAFAna/Decomp/IDecomp.h"
19 #include "CAFAna/Extrap/IExtrap.h"
22 
23 #include "OscLib/OscCalcDumb.h"
24 #include "OscLib/IOscCalc.h"
25 
26 #include "TStyle.h"
27 #include "TFile.h"
28 #include "TH1.h"
29 #include "TCanvas.h"
30 #include "TGraph.h"
31 #include "TGraphAsymmErrors.h"
32 #include "TGaxis.h"
33 #include "TLatex.h"
34 #include "TPaveText.h"
35 #include "TLegend.h"
36 #include "TLine.h"
37 #include "TSystem.h"
38 
39 #include <iostream>
40 #include <fstream>
41 #include <iomanip>
42 
43 #include "tex_tables_util.h"
44 
45 static bool isFHC=true; // Global variable for FHC/RHC
46 
47 //static const double kFDPOT= ana::kAna2017POT;
48 const bool printratioplots=true;
49 
51 
52 bool usedumb = false;
53 bool mergePeripheral=false;
54 
55 osc::IOscCalc * GetCalculator (bool usingdumb){
56  double dCP_Pi = 1.48;
57  double ssth23 = 0.558;
58  double dmsq32 = 2.44e-3;
59 
61  if(usingdumb){
62  osc::OscCalcDumb dumb;
63  calc = dumb.Copy();
64  }
65  else{
66  auto calc2 = DefaultOscCalc();
67  kFitDeltaInPiUnits.SetValue(calc2, dCP_Pi);
68  kFitSinSqTheta23.SetValue(calc2, ssth23);
69  kFitDmSq32.SetValue(calc2, dmsq32);
70  calc = calc2->Copy();
71  }
72  return calc;
73 }
74 
75 
76 //////////////////////////////////////////////////////////////////////
77 // Histogram
78 //////////////////////////////////////////////////////////////////////
79 
80 TH1* PrettyRatio( TH1* shi, TH1* nom, int color, double alpha = 1.,TString titley = "Shift / Nom"){
81  auto htemp = (TH1*) nom->Clone(ana::UniqueName().c_str());
82  htemp->Divide(shi,nom);
83  htemp->SetMarkerColorAlpha(color,alpha);
84  htemp->SetLineColorAlpha(color,alpha);
85  htemp->SetMarkerStyle(34);
86  htemp->SetMarkerSize(1.4);
87  htemp->GetYaxis()->SetTitle(titley);
88  return htemp;
89 }
90 
91 //////////////////////////////////////////////////////////////////////
92 
93 double FindLimitY (std::vector<TH1*> histos, bool doMax = true)
94 {
95  double max = histos[0]->GetMaximum();
96  double min = histos[0]->GetMinimum();
97 
98  for(unsigned int i=0; i<histos.size(); i++){
99  if ( histos[i]->GetMaximum() > max ) max = histos[i]->GetMaximum();
100  if ( histos[i]->GetMinimum() < min ) min = histos[i]->GetMinimum();
101  }
102 
103  if ( doMax ) return max;
104  else return min;
105 }
106 
107 //////////////////////////////////////////////////////////////////////
108 void FormatErrorBand(TH1* hplus, TH1* hminus,
109  bool signal=false, bool fixbins=false)
110 {
111  int mccolor = signal ? kViolet -9 : kTotalMCErrorBandColor;
112  std::cout << "hplus: " << hplus->Integral() << std::endl;
113  std::cout << "hminus: " << hminus->Integral() << std::endl;
114  hplus ->SetLineColor(mccolor);
115  hminus->SetLineColor(mccolor);
116  hplus ->SetFillColor(mccolor);
117  hminus->SetFillColor(10);
118  hplus ->SetMarkerColor(kRed-7);
119  hminus->SetMarkerColor(kBlue-7);
120 
121  if (fixbins)
122  {
123  for(int i=1; i <= hplus->GetNbinsX(); i++)
124  {
125  double tplus = hplus->GetBinContent(i);
126  double tminus = hminus->GetBinContent(i);
127  hplus->SetBinContent(i,max(tplus, tminus));
128  hminus->SetBinContent(i,min(tplus,tminus));
129  }
130  }
131 }
132 
133 //////////////////////////////////////////////////////////////////////
134 TH1* HistNumuOpt (TH1* orig){
135  auto ret = (TH1*) orig->Clone(UniqueName().c_str());
136  ret->Scale(0.1,"width");
137  ret->GetYaxis()->SetTitle("Events / 0.1 GeV");
138  ret->SetBit(kCanDelete);
139  ret->SetDirectory(gROOT);
140  return ret;
141 }
142 
143 //////////////////////////////////////////////////////////////////////
144 // Adornments
145 //////////////////////////////////////////////////////////////////////
146 
147 void PrettyTag(TString pid, int color, double xndc, double yndc){
148  TLatex* tag = new TLatex(xndc, yndc, pid);
149  tag->SetTextColor(color);
150  tag->SetNDC();
151  tag->SetTextSize(0.04);
152  tag->SetTextFont(42);
153  tag->SetTextAlign(22);
154  tag->Draw();
155  tag->SetBit(kCanDelete);
156 }
157 
158 void PIDTag(TString pid) { PrettyTag (pid, kAzure+3, .83, .92); }
159 void HTag (TString hie) { PrettyTag (hie, (hie=="NH")?kAzure+3:kOrange+2, .18, .825);}
160 
161 //////////////////////////////////////////////////////////////////////
162 // Legends
163 //////////////////////////////////////////////////////////////////////
164 
165 const double kLegendDefaultX1 = 0.63;
166 const double kLegendDefaultX2 = 0.85;
167 const double kLegendDefaultY1 = 0.68;
168 const double kLegendDefaultY2 = 0.88;
169 
170 TLegend * CustomLegend(
171  std::vector<TH1*>h, std::vector<TString> titles, std::vector<TString> option,
172  double x1=kLegendDefaultX1, double y1=kLegendDefaultY1,
173  double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
174 {
175  TLegend* leg = new TLegend(x1, y1, x2, y2);
176  leg->SetFillStyle(0);
177  for (unsigned int v = 0; v < h.size(); ++v)
178  leg->AddEntry(h[v], titles[v], option[v]);
179  leg->SetBit(kCanDelete);
180  return leg;
181 }
182 
183 //////////////////////////////////////////////////////////////////////
184 
185 void FixLegend(TLegend * leg, TString opt="default")
186 {
187  leg->SetFillStyle(0);
188  if(opt=="default")
189  {
190  leg->SetX1(kLegendDefaultX1); leg->SetX2(kLegendDefaultX2);
191  leg->SetY1(kLegendDefaultY1); leg->SetY2(kLegendDefaultY2);
192  }
193  if(opt=="ExPID")
194  {
195  leg->SetX1(0.405); leg->SetX2(0.565);
196  leg->SetY1(0.63); leg->SetY2(0.82);
197  leg->SetBorderSize(0);
198  leg->SetTextSize(0.033);
199  }
200  if(opt=="ExPIDPads")
201  {
202  leg->SetX1(0.55); leg->SetX2(0.89);
203  leg->SetY1(0.75); leg->SetY2(0.89);
204  leg->SetTextSize(0.045);
205  }
206 }
207 //////////////////////////////////////////////////////////////////////
208 
209 TGraph * DummyGraph(int color=kBlack, int style =1, int linewidth=2, int markerstyle = kFullCircle)
210 {
211  auto gr = new TGraph();
212  gr->SetLineColor(color);
213  gr->SetLineStyle(style);
214  gr->SetLineWidth(linewidth);
215  gr->SetMarkerStyle(markerstyle);
216  gr->SetMarkerColor(color);
217  return gr;
218 }
219 //////////////////////////////////////////////////////////////////////
220 
222  bool sig=true, bool tot=true, bool data=false,
223  double x1=kLegendDefaultX1, double y1=kLegendDefaultY1,
224  double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
225 {
226  TLegend* leg = new TLegend(x1,y1,x2,y2);
227  leg->SetFillStyle(0);
228  if(sig) leg->AddEntry( DummyGraph(kNueSignalColor), "#nu_{e} Signal","l");
229  if(data) leg->AddEntry( DummyGraph(), "ND Data","epl");
230  if(tot) leg->AddEntry( DummyGraph(kTotalMCColor), "Total Bkg","l");
231  leg->AddEntry( DummyGraph(kBeamNueBackgroundColor), "Beam #nu_{e}","l");
232  leg->AddEntry( DummyGraph(kNCBackgroundColor), "NC","l");
233  leg->AddEntry( DummyGraph(kNumuBackgroundColor), "#nu_{#mu} CC","l");
234  leg->SetBit(kCanDelete);
235  return leg;
236 
237 }
238 //////////////////////////////////////////////////////////////////////
239 
240 TLegend * DefaultNumuLegend (bool sig=true, bool tot=true,bool data=false,
241  double x1=kLegendDefaultX1, double y1=kLegendDefaultY1,
242  double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
243 {
244  TLegend* leg = new TLegend(x1,y1,x2,y2);
245  leg->SetFillStyle(0);
246  if(sig) leg->AddEntry( DummyGraph(kNueSignalColor), "#nu_{#mu} Signal","l");
247  if(data) leg->AddEntry( DummyGraph(), "ND Data","epl");
248  if(tot) leg->AddEntry( DummyGraph(kTotalMCColor), "Total Bkg","l");
249  if(isNumuAna) leg->AddEntry( DummyGraph(kNumuBackgroundColor), "#nu_{#mu} app","l");
250  else leg->AddEntry( DummyGraph(kNumuBackgroundColor), "#nu_{#mu} CC" ,"l");
251  leg->AddEntry( DummyGraph(kNCBackgroundColor), "NC","l");
252  leg->AddEntry( DummyGraph(kBeamNueBackgroundColor), "#nu_{e} CC","l");
253  leg->AddEntry( DummyGraph(kNuTauBackgroundColor), "#nu_{#tau} CC","l");
254  leg->SetBit(kCanDelete);
255  return leg;
256 
257 }
258 //////////////////////////////////////////////////////////////////////
259 
260 void PrintCanvasAll (TCanvas * c, const TString outname, const TString outdir = "plots/",
261  const std::vector <TString> exts = {"pdf","png","root"} )
262 {
263  for(auto const & ext:exts){
264  gSystem -> MakeDirectory( outdir );
265  gSystem -> MakeDirectory( outdir+"rootfiles/" );
266  if(ext=="png") c->SaveAs(outdir+"png_"+outname + "." + ext);
267  else if(ext=="root") c->SaveAs(outdir+"rootfiles/"+outname + "." + ext);
268  else c->SaveAs(outdir+outname + "." + ext);
269  }
270  delete c;
271 }
272 
273 //////////////////////////////////////////////////////////////////////
274 // Nue 3 bins functions
275 //////////////////////////////////////////////////////////////////////
276 
277 // NB, loses all the style info, so call early
278 TH1* RebinToShortAxis(TH1* h)
279 {
280  TH1* ret = new TH1D(UniqueName().c_str(), h->GetTitle(),
281  18, 0, 18);
282  ret->GetXaxis()->SetTitle(h->GetXaxis()->GetTitle());
283  ret->GetYaxis()->SetTitle(h->GetYaxis()->GetTitle());
284 
285  std::vector<std::pair<int, int>> binMap = {{3, 2},
286  {4, 3},
287  {5, 4},
288  {6, 5}, // end of 1st spectrum
289  {13, 8},
290  {14, 9},
291  {15, 10},
292  {16, 11}, // end of 2nd
293  {23, 14},
294  {24, 15},
295  {25, 16},
296  {26, 17}};
297 
298  for(auto it: binMap) ret->SetBinContent(it.second, h->GetBinContent(it.first));
299 
300  return ret;
301 }
302 
303 //////////////////////////////////////////////////////////////////////
304 void NuePID3bin(TH1* h, TString pidtag, bool binlabels=true)
305 {
307  if (binlabels) Nue2018ThreeBinLabels();
308  //if (h) PID2DAxis(h);
309 }
310 //////////////////////////////////////////////////////////////////////
311 TCanvas * ExPIDPlot (std::vector<TH1*> topHistos,
312  std::vector<TString> topOption)
313 {
314  TCanvas *cpid = new TCanvas("cpid", "cpid", 410, 850);
315  float Ebins = 10;
316  int pidBins = 3;
317 
318  TPad* p[pidBins];
319  p[0] = new TPad("", "", 0, 0, 1, 1);
320  p[1] = new TPad("", "", 0, 0, 1, 1);
321  p[2] = new TPad("", "", 0, 0, 1, 1);
322  p[0]->SetTopMargin(0.1);
323  p[0]->SetBottomMargin(.634);
324  p[1]->SetTopMargin(.369);
325  p[1]->SetBottomMargin(.369);
326  p[2]->SetTopMargin(.634);
327  p[2]->SetBottomMargin(0.1);
328  p[0]->SetFillStyle(0);
329  p[1]->SetFillStyle(0);
330  p[2]->SetFillStyle(0);
331  p[0]->Draw();
332  p[1]->Draw();
333  p[2]->Draw();
334 
335  TH1* he[pidBins];
336  TH1* h[topHistos.size()][pidBins];
337  for(int k=pidBins-1; k>=0; --k){
338  he[k] = (TH1*) topHistos[0]->Clone();
339  he[k]->GetYaxis()->SetRangeUser(FindLimitY(topHistos,0) ,1.2*FindLimitY(topHistos,1));
340  he[k]->GetXaxis()->SetRangeUser(Ebins*k,Ebins*(k+1));
341  if ( k>0 ) he[k]->GetXaxis()->SetLabelSize(0);
342  p[pidBins-1-k]->cd();
343  he[k]->Draw(topOption[0]);
344  he[k]->GetYaxis()->SetTitleSize(24);
345  he[k]->GetYaxis()->SetTitleFont(43);
346  he[k]->GetYaxis()->SetTitleOffset(1.25);
347  he[k]->GetYaxis()->SetLabelFont(43);
348  he[k]->GetYaxis()->SetLabelSize(14);
349  he[k]->GetXaxis()->SetLabelOffset(-0.008);
350  he[k]->GetXaxis()->SetTitleOffset(0.5);
351  p[pidBins-1-k]->SetGridy();
352  }
353  for(unsigned int ii=0;ii<topHistos.size();ii++){
354  for(int j=pidBins-1; j>=0; --j){
355  h[ii][j] = (TH1*) topHistos[ii]->Clone();
356  h[ii][j]->GetXaxis()->SetRangeUser(Ebins*j,Ebins*(j+1));
357  p[pidBins-1-j]->cd();
358  h[ii][j]->Draw("same "+topOption[ii]);
359  }
360  }
361  for(int h=pidBins-1; h>=0; --h){
362  he[h]->SetStats(0);
363  p[pidBins-1-h]->cd();
364  he[h]->Draw("same "+topOption[0]);
365  }
366  p[0]->RedrawAxis(); //needed for white histograms e.g. beam
367  p[1]->RedrawAxis();
368  p[2]->RedrawAxis();
369  return cpid;
370 }
371 
372 //////////////////////////////////////////////////////////////////////
373 // Pretty Canvas
374 //////////////////////////////////////////////////////////////////////
375 TCanvas * RatioPlot (std::vector<TH1*> topHistos,
376  std::vector<TString> topOption,
377  std::vector<TH1*> bottomHistos,
378  std::vector<TString> bottomOption,
379  TString pidtag,
380  AxisType pidaxis = kNue1bin)
381 {
382  gStyle->SetTitleStyle(0);
383  TCanvas *c = new TCanvas("c", "canvas", 600, 700);
384  if(pidaxis == kNue3bin) c->SetCanvasSize(750,650);
385  //define top and bottom pads
386  TPad* p1 = new TPad("", "", 0, 0, 1, 1);
387  TPad* p2 = new TPad("", "", 0, 0, 1, 1);
388  p1->SetBottomMargin(.3);
389  p2->SetTopMargin(.7);
390 
391  for(auto p:{p1,p2}){
392  p->SetFillStyle(0);
393  p->Draw();
394  }
395 
396  // format histograms
397  auto h1 = (TH1*) topHistos[0]->Clone();
398  auto h3 = (TH1*) bottomHistos[0]->Clone();
399 
400  for(auto & h:{h1,h3}){
401  h->SetStats(0);
402 
403  h->GetYaxis()->SetTitleSize(26);
404  h->GetYaxis()->SetTitleFont(43);
405  h->GetYaxis()->SetTitleOffset(1.3);
406  h->GetYaxis()->SetLabelFont(43);
407 
408  h->GetXaxis()->SetTitleOffset(1.2);
409  h->GetXaxis()->SetTitleSize(28);
410  h->GetXaxis()->SetTitleFont(43);
411  h->GetXaxis()->SetLabelFont(43);
412  }
413 
414  h1->GetYaxis()->SetRangeUser(FindLimitY(topHistos,0), 1.2*FindLimitY(topHistos,1));
415  h1->GetYaxis()->SetLabelSize(18);
416 
417  h1->GetXaxis()->SetLabelSize(0);
418  h1->GetXaxis()->SetTitleSize(0);
419 
420  h3->SetTitle("");
421  h3->GetYaxis()->SetRangeUser(max(0.,2-FindLimitY(bottomHistos,1)),
422  1.01*FindLimitY(bottomHistos,1));
423  h3->GetYaxis()->SetLabelSize(15);
424  h3->GetYaxis()->SetDecimals();
425 
426  h3->GetXaxis()->SetTitle(h1->GetXaxis()->GetTitle());
427  h3->GetXaxis()->SetLabelSize(18);
428  h3->GetXaxis()->SetLabelOffset(0.005);
429 
430  p1->cd();
431 
432  h1->Draw(topOption[0]);
433 
434  for(unsigned int ii = 1; ii < topHistos.size(); ++ii) {
435  topHistos[ii]->Draw("same "+topOption[ii]);}
436  h1->Draw("same "+topOption[0]);
437 
438  p1->RedrawAxis(); //white histograms cover axes
439 
440  if(pidaxis==kNue3bin) NuePID3bin(0,pidtag);
441 
442  p2->cd();
443  h3->Draw(bottomOption[0]);
444  for(unsigned int ii = 1; ii < bottomHistos.size(); ++ii){
445  bottomHistos[ii]->Draw("same " + bottomOption[ii]);}
446  h3->Draw("same " + bottomOption[0]);
447 
448  auto lone = new TLine();
449  lone->SetLineStyle(3);
450  lone->SetLineColor(kGray+3);
451  p2->Update();
452  lone->DrawLine(p2->GetUxmin(),1.0,p2->GetUxmax(),1.0);
453  lone->DrawLine(p2->GetUxmin(),1.2,p2->GetUxmax(),1.2);
454  lone->DrawLine(p2->GetUxmin(),0.8,p2->GetUxmax(),0.8);
455  float p2Y = p2->GetUymax()-p2->GetUymin();
456  h3->GetYaxis()->SetRangeUser(p2->GetUymin()-(0.01*p2Y),p2->GetUymax()+(0.01*p2Y));
457  p2->RedrawAxis();
458  if(pidaxis==kNue3bin) NuePID3bin(h3,pidtag, false);
459 
460  return c;
461 }
462 //////////////////////////////////////////////////////////////////////
463 // Pretty plot styles: Data MC, MC components, error bands
464 //////////////////////////////////////////////////////////////////////
465 
466 void PlotDataMC(std::vector<TH1*>vnom, std::vector<bool>isdata,
467  TLegend * leg, TString pidtag="",TString htag="",
468  TString out_name="plot_FD",bool ratioplot=false,
469  bool ratioerrors=false, AxisType pidaxis = kNue1bin)
470 {
471  gStyle->SetPadLeftMargin(0.15);
472  gStyle->SetPadTopMargin(0.08);
473 
474  // Fill empty bins with 0 if in range and <0 if out of range
475  for ( int bin = 1; bin<=vnom[0]->GetNbinsX(); ++bin){
476  bool filledbin=false;
477  for ( unsigned int v = 0; v<vnom.size(); ++v){
478  if ( !isdata[v] )
479  if (vnom[v]->GetBinContent(bin) > 0 ) filledbin=true;
480  }
481  for ( uint v = 0; v<vnom.size(); ++v){
482  if (vnom[v]->GetBinContent(bin) == 0 ){
483  if ( !filledbin ) vnom[v]->SetBinContent(bin,-0.01);
484  else vnom[v]->SetBinContent(bin,0.001);
485  }
486  vnom[v]->GetYaxis()->SetTitleOffset(1.05);
487  }
488  }
489 
490  // Draw marker w/error or line when appropriate
491  int datahist=-1;
492  for ( unsigned int v = 0; v<vnom.size(); ++v)
493  if (vnom[v]->GetLineColor() == kBlack ) datahist=v;
494 
495  // Option to make a ratio plot given
496  if(!ratioplot){
497  auto c2 = new TCanvas ("c2","c2",715,500);
498  for ( unsigned int v = 0; v<vnom.size(); ++v){
499  TGraphAsymmErrors* gr = GraphWithPoissonErrors(vnom[v], 0, 1);
500  if ( isdata[v] ){
501  if (v == datahist ) gr->Draw("ep same");
502  else vnom[v]->Draw("p same");
503  }
504  else vnom[v]->Draw("hist same");
505  }
506  gPad->RedrawAxis();
507 
508  // Adjust for ExPID axis
509  if(pidaxis==kNue3bin)
510  {
511  FixLegend(leg);
514  }
515 
516  // pot label
517  PrettyTag(TString::Format("%g#times10^{20} POT equiv.", 8.85).Data(),
518  kBlack, .175, .74);
519  HTag(htag);
520  leg->Draw("same");
521 
522  gPad->SetFillStyle(0);
523  c2->SetFillStyle(0);
524  PrintCanvasAll (c2,out_name+"_datamc","plots/",{"pdf","eps","png","C","root"});
525  PrintLatexFigure(out_name+"_datamc");
526  //removed: prelim tag
527  }
528  else{
529  std::vector<TH1*> vratio;
530  for(unsigned int ii=1;ii<vnom.size();ii++){
531  auto htemp = PrettyRatio (vnom[datahist], vnom[ii], vnom[ii]->GetLineColor(), 1.,
532  "Data / MC");
533  htemp->GetYaxis()->SetRangeUser(0.5,1.5);
534  vratio.push_back(htemp);
535  }
536  std::vector<TString> topOption = {" ","hist"};
537  std::vector<TString> bottomOption (vratio.size(),ratioerrors ?"p":"hist p");
538  auto cratio2= RatioPlot(vnom, topOption,vratio,bottomOption,pidtag,pidaxis);
539  cratio2->cd();
540  leg->Draw("same");
541  HTag(pidtag);
542  PrintCanvasAll(cratio2,out_name+"_datamc_ratio");
543  PrintLatexFigure("/png_"+out_name+"_datamc_ratio");
544  }
545 }
546 
547 //////////////////////////////////////////////////////////////////////
548 
549 void PlotNDDataTotalMCComparison(TH1* hdata,std::vector<TH1*> htots,
550  TLegend * leg,TString pidtag="",
551  TString out_name="plot_nd",
552  bool ratioplot=false, AxisType pidaxis = kNue1bin)
553 {
554  if(pidaxis==kNue3bin)
555  FixLegend(leg);
556  std::vector<TString> h_opt;
557  for(auto const & h:htots) h_opt.push_back("hist");
558  h_opt.push_back("p");
559 
560  if(!ratioplot)
561  {
562  auto c1 = new TCanvas ("c1","c1");
563  htots[0]->SetMaximum(1.3*htots[0]->GetMaximum());
564  for(auto const & h:htots) h->Draw("hist same");
565  htots[0]->Draw("hist same"); //need to redraw for error bands
566  hdata->Draw("p same");
567  PIDTag(pidtag);
568  if(pidaxis==kNue3bin) NuePID3bin (htots[0],pidtag,true);
569  leg->Draw("same");
570  PrintCanvasAll(c1, out_name+"_datamc");
571  PrintLatexFigure("/png_"+out_name+"_datamc");
572  }
573  else
574  {
575  std::vector<TH1*> ratios;
576  std::vector<TString> r_opt;
577  for(auto const & h:htots){
578  auto ratio = PrettyRatio(hdata,h, h->GetLineColor(), 1., "Data / MC");
579  ratio->GetYaxis()->SetRangeUser(0.5,1.5);
580  ratios.push_back(ratio);
581  r_opt.push_back("e1");
582  }
583  htots.push_back(hdata);
584  auto cratio1 = RatioPlot(htots,h_opt,ratios,r_opt,pidtag,pidaxis);
585  cratio1->cd();
586  leg->Draw("same");
587  PIDTag(pidtag);
588 
589  PrintCanvasAll(cratio1,out_name+"_datamc_ratio");
590  PrintLatexFigure("/png_"+out_name+"_datamc_ratio");
591  }
592 }
593 //////////////////////////////////////////////////////////////////////
595  std::vector<TH1*>vnom, std::vector<TH1*>vshi,
596  TLegend * leg, TString pidtag="",
597  TString out_name="plot_nd"/* Why ND? */, bool ratioplot=false,
598  bool ratioerrors=false, AxisType pidaxis=kNue1bin)
599 {
600  if(pidaxis==kNue3bin) FixLegend(leg);
601 
602  if(!ratioplot)
603  {
604  auto c2 = new TCanvas (UniqueName().c_str(), "c2");
605  for(auto const & v:vnom) v->Draw("hist same");
606  for(auto const & v:vshi) v->Draw("hist same");
607  if(pidaxis==kNue3bin) NuePID3bin (vnom[0],pidtag,true);
608  PIDTag(pidtag);
609  leg->Draw("same");
610  PrintCanvasAll(c2, out_name+"_mccomp");
611  PrintLatexFigure("/png_"+out_name+"_mccomp");
612  c2->Close();
613  }
614  else
615  {
616  std::vector<TH1*> vratio;
617  for(unsigned int ii=0;ii<vshi.size();ii++){
618  auto htemp = PrettyRatio(vshi[ii], vnom[ii],
619  vnom[ii]->GetLineColor(), ii>0 ? 0.4:1.);
620  htemp->GetYaxis()->SetRangeUser(0.5, 1.5);
621  vratio.push_back (htemp);
622  }
623  vnom.insert(vnom.end(),vshi.begin(),vshi.end());
624  std::vector<TString> topOption(vnom.size(), "hist");
625  std::vector<TString> bottomOption(vratio.size(),
626  ratioerrors ? "p":"hist p");
627  auto cratio2= RatioPlot(vnom, topOption, vratio,
628  bottomOption, pidtag, pidaxis);
629  cratio2->cd();
630  leg->Draw("same");
631  PIDTag(pidtag);
632  PrintCanvasAll(cratio2, out_name+"_mccomp_ratio");
633  PrintLatexFigure("/png_"+out_name+"_mccomp_ratio");
634  }
635  //if(pidaxis==kNue3bin) {
636  // vnom.insert(vnom.end(),vshi.begin(),vshi.end());
637  // std::vector<TString> topOption (vnom.size(),"hist");
638  // auto cpads =ExPIDPlot(vnom,topOption);
639  // cpads->cd();
640  // FixLegend(leg,"default");
641  // leg->Draw("same");
642  // PIDTag(pidtag);
643  // PrintCanvasAll(cpads, out_name+"_mccomp_pads");
644  // PrintLatexFigure("/png_"+out_name+"_mccomp_pads");
645  //}
646 }
647 
648 //////////////////////////////////////////////////////////////////////
649 
650 void PlotMCComponentsErrorBand(std::vector<TH1*>vnom,std::vector<TH1*>vshi,
651  TLegend * leg, TString pidtag="",
652  TString out_name="plot_nd",bool ratioplot=false,
653  bool ratioerrors=false,AxisType pidaxis = kNue1bin)
654 {
655  auto hplu_plot = (TH1*) vshi[0]->Clone();
656  auto hmin_plot = (TH1*) vshi[1]->Clone();
657  FormatErrorBand(hplu_plot,hmin_plot,out_name.Contains("_sig"),true);
658  std::vector<TH1*>vshi_plot={hplu_plot,hmin_plot};
659 
660  if(pidaxis==kNue3bin)
661  FixLegend(leg);
662 
663  if(!ratioplot){
664  auto c2 = new TCanvas ("c2","c2");
665  vnom[0]->SetMaximum(1.2*vnom[0]->GetMaximum());
666  vnom[0]->Draw("hist");
667  for(auto v:vshi_plot) v->Draw("hist same");
668  for(auto v:vnom) v->Draw("hist same");
669  if(pidaxis==kNue3bin) NuePID3bin (vnom[0],pidtag,true);
670  PIDTag(pidtag);
671  leg->Draw("same");
672  PrintCanvasAll(c2, out_name+"_mccomp");
673  PrintLatexFigure("/png_"+out_name+"_mccomp");
674  c2->Close();
675  }
676  else{
677  std::vector<TH1*> vratio;
678  for(unsigned int ii=0;ii<vshi.size();ii++){
679  auto htemp = PrettyRatio(vshi[ii],vnom[0],vshi[ii]->GetMarkerColor());
680  vratio.push_back(htemp);
681  }
682  vratio[0]->GetYaxis()->SetRangeUser(0.5,1.5);
683  vnom.insert(vnom.begin()+1,vshi_plot.begin(),vshi_plot.end());
684  std::vector<TString> topOption (vnom.size(),"hist");
685  std::vector<TString> bottomOption (vratio.size(),ratioerrors ?"p":"hist p");
686  auto cratio2= RatioPlot(vnom, topOption,vratio,bottomOption,pidtag,pidaxis);
687  cratio2->cd();
688  leg->Draw("same");
689  PIDTag(pidtag);
690  PrintCanvasAll(cratio2, out_name+"_mccomp_ratio");
691  PrintLatexFigure("/png_"+out_name+"_mccomp_ratio");
692  }
693 }
694 
695 //////////////////////////////////////////////////////////////////////
696 
697 
698 
699 //////////////////////////////////////////////////////////////////////
700 // Define container for histograms + format
701 //////////////////////////////////////////////////////////////////////
702 
703 namespace ana
704 {
705 
706  struct DataMCComponents
707  {
708  TH1* data;
709  TH1* sig, *wrongSign, *bkg, *beam, *numu, *nc, *tau;
711  };
712 
713  void DefaultFormatNue (DataMCComponents & comp, int linestyle = 1){
714  std::vector <std::pair <TH1*, int>> pairs = {
715  {comp.sig, kNueSignalColor},
716  {comp.wrongSign, kNueSignalWSColor},
717  {comp.beam, kBeamNueBackgroundColor},
718  {comp.nc, kNCBackgroundColor},
719  {comp.numu, kNumuBackgroundColor},
720  {comp.bkg, kTotalMCColor},
721  {comp.nuetonumu, kNumuBackgroundColor},
722  {comp.numutonumu, kNumuBackgroundColor},
723  {comp.nuetonutau, kGray},
724  {comp.numutonutau, kGray}
725  };
726  for (auto const & mypair:pairs){
727  if(!mypair.first) continue;
728  mypair.first->SetLineColor(mypair.second);
729  mypair.first->SetLineStyle(linestyle);
730  ana::CenterTitles(mypair.first);
731  }
732  if(comp.data) comp.data->SetMarkerStyle(kFullCircle);
733  return;
734  }
735 //////////////////////////////////////////////////////////////////////
736 
738  {
739 
740  // Get correct POT for horn current
741  double ndPOT;
742 
743  if (isFHC)
744  ndPOT = kAna2018FHCPOT;
745  else
746  ndPOT = kAna2018RHCPOT;
747 
748  // Get components from decomp
749  Spectrum nueComp = decomp->NueComponent();
750  std::cout << "Loaded NueComponent..." << std::endl;
751  Spectrum antinueComp = decomp->AntiNueComponent();
752  std::cout << "Loaded AnitNueComponent..." << std::endl;
753  Spectrum numuComp = decomp->NumuComponent();
754  std::cout << "Loaded NumuComponent..." << std::endl;
755  Spectrum antinumuComp = decomp->AntiNumuComponent();
756  std::cout << "Loaded AntiNumuComponent..." << std::endl;
757  Spectrum ncTotComp = decomp->NCTotalComponent();
758  std::cout << "Loaded NCTotalComponent..." << std::endl;
759  Spectrum dataComp = decomp->Data_Component();
760  std::cout << "Loaded Data_Component..." << std::endl;
761 
762  auto beam = nueComp + antinueComp;
763  auto numu = numuComp + antinumuComp;
764 
766  ret.beam = beam.ToTH1(ndPOT);
767  std::cout << "beam component: " << ret.beam->Integral()
768  << ", in " << ret.beam->GetNbinsX() << " bins" << std::endl;
769  ret.numu = numu.ToTH1(ndPOT);
770  std::cout << "numu component: " << ret.numu->Integral()
771  << ", in " << ret.beam->GetNbinsX() << " bins" << std::endl;
772  ret.nc = ncTotComp.ToTH1(ndPOT);
773  std::cout << "nc component: " << ret.nc->Integral()
774  << ", in " << ret.beam->GetNbinsX() << " bins" << std::endl;
775  ret.data = dataComp.ToTH1(ndPOT);
776  std::cout << "data component: " << ret.data->Integral()
777  << ", in " << ret.beam->GetNbinsX() << " bins" << std::endl;
778  ret.bkg = reinterpret_cast<TH1*>(ret.beam->Clone());
779  ret.bkg->Add(ret.numu);
780  ret.bkg->Add(ret.nc);
781  std::cout << "bkg component: " << ret.bkg->Integral()
782  << ", in " << ret.beam->GetNbinsX() << " bins" << std::endl;
783  ret.sig = 0;
784  ret.tau = 0;
785  ret.nuetonumu = 0;
786  ret.nuetonutau = 0;
787  ret.numutonumu = 0;
788  ret.numutonutau = 0;
789 
790  std::cout << "Loaded ND components" << std::endl;
791  return ret;
792  }
793 
794 //////////////////////////////////////////////////////////////////////
795 
796  struct DataMCComponents GetFDMCComponents(
797  osc::IOscCalc * calc, IPrediction * pred_no, TString
798  output_name="nue", int linestyle=1, bool bkgdetails=false)
799  {
801  // Set some defaults for the components.
802  Flavors::Flavors_t MySig, MyBeam, MyNuMu, MyNC, MyTau;
803  MySig = MyBeam = MyNuMu = MyNC = MyTau = Flavors::kAll;
804  // Change the components depending on what I'm looking at.
805  if (!isNumuAna)
806  { // Nu e
807  MySig = Flavors::kNuMuToNuE;
808  MyBeam = Flavors::kNuEToNuE;
809  MyNuMu = Flavors::kAllNuMu;
810  MyNC = Flavors::kAll;
811  MyTau = Flavors::kAllNuTau;
812  }
813  else
814  { // Nu mu
815  MySig = Flavors::kNuMuToNuMu;
816  MyBeam = Flavors::kAllNuE;
817  MyNuMu = Flavors::kNuEToNuMu;
818  MyNC = Flavors::kAll;
819  MyTau = Flavors::kAllNuTau;
820  }
821 
822  // Specify sign depending on polarity
824  sign = wrongSign = Sign::kBoth;
825  if (!isNumuAna) {
826  if (isFHC)
827  {
828  sign = Sign::kNu;
829  wrongSign = Sign::kAntiNu;
830  }
831  else
832  {
833  sign = Sign::kAntiNu;
834  wrongSign = Sign::kNu;
835  }
836  }
837  const double pot = (isFHC) ? kAna2018FHCPOT:kAna2018RHCPOT;
838 
839  // Predict components
840  ret.sig = pred_no->PredictComponent(
841  calc, MySig, Current::kCC, sign).ToTH1(pot);
842  ret.wrongSign = pred_no->PredictComponent(
843  calc, MySig, Current::kCC, wrongSign).ToTH1(pot);
844  ret.beam = pred_no->PredictComponent(
845  calc, MyBeam, Current::kCC, Sign::kBoth).ToTH1(pot);
846  ret.numu = pred_no->PredictComponent(
847  calc, MyNuMu, Current::kCC, Sign::kBoth).ToTH1(pot);
848  ret.nc = pred_no->PredictComponent(
849  calc, MyNC, Current::kNC, Sign::kBoth).ToTH1(pot);
850  ret.tau = pred_no->PredictComponent(
851  calc, MyTau, Current::kCC, Sign::kBoth).ToTH1(pot);
852 
853  ret.bkg = (pred_no->Predict(calc) - pred_no->PredictComponent(
854  calc, MySig, Current::kCC, sign)).ToTH1(pot);
855 
856  if(bkgdetails)
857  {
858  ret.nuetonumu = pred_no->PredictComponent(
860  ret.numutonumu = pred_no->PredictComponent(
862  ret.nuetonutau = pred_no->PredictComponent(
864  ret.numutonutau = pred_no->PredictComponent(
866  }
867  else
868  {
869  ret.nuetonumu = 0; ret.numutonumu = 0;
870  ret.nuetonutau = 0; ret.numutonutau = 0;
871  }
872 
873  ret.data = 0;
874  DefaultFormatNue(ret, linestyle);
875 
876  return ret;
877  }
878 
879 //////////////////////////////////////////////////////////////////////
880 // Compare ND Data MC
881 //////////////////////////////////////////////////////////////////////
882 
884  TString plottitle, TString out_name, TString pidtag,
885  AxisType pidaxis, bool printtable=true)
886  {
887  h_sh.bkg->SetLineColor(kTotalMCErrorBandColor);
888  h_sh.bkg->SetLineStyle(1);
889 
890  h_no.bkg->SetMaximum(FindLimitY({h_no.bkg,h_no.data}));
891  h_no.bkg->SetTitle(plottitle);
892 
893  TGaxis::SetMaxDigits(3);
894  h_no.bkg->GetXaxis()->SetDecimals();
895 
896  auto leg1 = CustomLegend({h_no.data, h_no.bkg, h_sh.bkg},{"ND Data","Nominal","Shift"},{"epl","l","l"});
897 
898  auto leg2 = DefaultNueLegend(false,true,false);
899 
900  TH1 * mydata = h_no.data;
901  std::vector <TH1*> comp_plot = {h_no.bkg, h_sh.bkg};
902  std::vector <TH1*> vnom_plot = {h_no.bkg, h_no.beam,h_no.nc,h_no.numu};
903  std::vector <TH1*> vshi_plot = {h_sh.bkg, h_sh.beam,h_sh.nc,h_sh.numu};
904 
905  std::vector <TH1*> vnom = {h_no.bkg,h_no.beam,h_no.nc,h_no.numu};
906  std::vector <TH1*> vshi = {h_sh.bkg,h_sh.beam,h_sh.nc,h_sh.numu};
907  std::vector<TString> labels = {"Total MC","Beam $\\nu_{e}$","NC", "$\\nu_{\\mu}$ CC"};
908 
909  if(pidaxis == kNumuOpt){
910  leg2 = DefaultNumuLegend(false,true,false);
911  mydata = HistNumuOpt(h_no.data);
912  comp_plot = { HistNumuOpt(h_no.bkg), HistNumuOpt(h_sh.bkg) };
913  vnom_plot = { HistNumuOpt(h_no.bkg), HistNumuOpt(h_no.numu), HistNumuOpt(h_no.nc)};
914  vshi_plot = { HistNumuOpt(h_sh.bkg), HistNumuOpt(h_sh.numu), HistNumuOpt(h_sh.nc)};
915 
916  vnom = {h_no.bkg,h_no.numu,h_no.nc,h_no.beam};
917  vshi = {h_sh.bkg,h_sh.numu,h_sh.nc,h_sh.beam};
918  labels = {"Total MC","$\\nu_{\\mu}$ CC","NC", "$\\nu_{e}$ CC"};
919  }
920 
921  PlotNDDataTotalMCComparison(mydata, comp_plot,
922  leg1, pidtag, out_name, printratioplots, pidaxis);
923 
924  h_sh.bkg->SetLineColor(kTotalMCColor+1);
925  h_sh.bkg->SetLineStyle(7);
926 
927  leg2->AddEntry(DummyGraph(kGray+1,2),"Shifted","l");
928 
929  PlotMCComponentsComparison(vnom_plot, vshi_plot,
930  leg2, pidtag, out_name, printratioplots, false, pidaxis);
931 
932  if(printtable){
933 
934  ComparisonTableOne (vnom, vshi, labels, out_name);
935 
936  int nbins = 1;
937  if (pidaxis== kNue3bin) nbins = 3;
938  ComparisonTableOneNbins({h_no.bkg,h_no.beam,h_no.nc,h_no.numu},
939  {h_sh.bkg,h_sh.beam,h_sh.nc,h_sh.numu},
940  labels,
941  out_name,
942  nbins);
943  }
944 
945  }
946 //////////////////////////////////////////////////////////////////////
947 
948  void CompareNDDataTwoMC(
950  TString plottitle, TString out_name, TString pidtag,
951  AxisType pidaxis, bool printtable=true)
952  {
953  //FormatErrorBand(hplu.bkg, hmin.bkg);
954 
955  /*
956  std::cout << "DEBUG[CompareNDDataTwoMC]: hnom.bkg "
957  << hnom.bkg->Integral() << std::endl;
958  std::cout << "DEBUG[CompareNDDataTwoMC]: hplu.bkg "
959  << hplu.bkg->Integral() << std::endl;
960  std::cout << "DEBUG[CompareNDDataTwoMC]: hmin.bkg "
961  << hmin.bkg->Integral() << std::endl;
962 
963  std::cout << "DEBUG[CompareNDDataTwoMC]: hnom.data "
964  << hnom.data->Integral() << std::endl;
965  std::cout << "DEBUG[CompareNDDataTwoMC]: hplu.data "
966  << hplu.data->Integral() << std::endl;
967  std::cout << "DEBUG[CompareNDDataTwoMC]: hmin.data "
968  << hmin.data->Integral() << std::endl;*/
969 
970  hnom.bkg->SetMaximum(FindLimitY({hnom.bkg, hnom.data, hplu.bkg, hmin.bkg}));
971  hnom.bkg->SetTitle(plottitle);
972 
973  TGaxis::SetMaxDigits(3);
974  hnom.bkg->GetXaxis()->SetDecimals();
975 
976  auto leg1 = CustomLegend(
977  {hnom.data, hnom.bkg, hplu.bkg},
978  {"ND Data","Nominal","Syst."}, {"epl","l","f"});
979  auto leg2 = DefaultNueLegend(false, true, false);
980 
981  TH1 * mydata = hnom.data;
982  std::vector <TH1*> comp_plot = {hnom.bkg, hplu.bkg, hmin.bkg};
983  std::vector <TH1*> vnom_plot = {hnom.bkg, hnom.beam, hnom.nc, hnom.numu};
984  std::vector <TH1*> vshi_plot = {hplu.bkg, hmin.bkg};
985 
986  std::vector <TH1*> vnom = {hnom.bkg, hnom.beam, hnom.nc, hnom.numu};
987  std::vector <TH1*> vplu = {hplu.bkg, hplu.beam, hplu.nc, hplu.numu};
988  std::vector <TH1*> vmin = {hmin.bkg, hmin.beam, hmin.nc, hmin.numu};
989  std::vector<TString> labels = {
990  "Total MC", "Beam $\\nu_{e}$", "NC", "$\\nu_{\\mu}$ CC"};
991 
992  if(pidaxis == kNumuOpt)
993  {
994  leg2 = DefaultNumuLegend(false, true, false);
995  mydata = HistNumuOpt(hnom.data);
996  comp_plot = {HistNumuOpt(hnom.bkg),
997  HistNumuOpt(hplu.bkg), HistNumuOpt(hmin.bkg)};
998 
999  vnom_plot = {HistNumuOpt(hnom.bkg),
1000  HistNumuOpt(hnom.numu), HistNumuOpt(hnom.nc)};
1001  vshi_plot = {HistNumuOpt(hplu.bkg), HistNumuOpt(hmin.bkg)};
1002 
1003  vnom = {hnom.bkg,hnom.numu,hnom.nc,hnom.beam};
1004  vplu = {hplu.bkg,hplu.numu,hplu.nc,hplu.beam};
1005  vmin = {hmin.bkg,hmin.numu,hmin.nc,hmin.beam};
1006  labels = {"Total MC", "$\\nu_{\\mu}$ CC", "NC", "$\\nu_{e}$ CC"};
1007  }
1008 
1009  PlotNDDataTotalMCComparison(mydata, comp_plot,
1010  leg1, pidtag, out_name, printratioplots, pidaxis);
1011 
1012  leg2->AddEntry(hplu.bkg,"Syst.","f");
1013  PlotMCComponentsErrorBand(vnom_plot, vshi_plot,
1014  leg2, pidtag, out_name, printratioplots, false, pidaxis);
1015 
1016  if(printtable){
1017  ComparisonTable (vnom, vplu, vmin, labels,
1018  out_name);
1019  int nbins = 1;
1020  if (pidaxis== kNue3bin) nbins = 3;
1021  ComparisonTableNbins(vnom, vplu, vmin, labels,
1022  out_name, nbins);
1023  }
1024  }
1025 
1026 //////////////////////////////////////////////////////////////////////
1027 
1029  {
1034  };
1035 
1036  const IDecomp* GetDecomp(IPrediction* prediction,
1037  EModExtrapComps modExtrapComp)
1038  {
1039  // First cast IPrediction as a PredictionExtrap object and check it extist
1040  //std::cout << "Checking if IPrediction exists" << std::endl;
1041  if (!prediction)
1042  abort();
1043 
1044  PredictionExtrap* predExtrap;
1045 
1046  // If numu then it's simmple;
1047  if (isNumuAna) {
1048  predExtrap = dynamic_cast<PredictionExtrap*>(prediction);
1049  } else {
1050  // If nue then we have to do some wrangling to get the correct object
1051  //std::cout << DemangledTypeName(prediction) << std::endl;
1052  // This will be a PredictionExtendToPeripheral so cast as that
1053  PredictionExtendToPeripheral* predExtend = dynamic_cast<PredictionExtendToPeripheral*>(prediction);
1054  //std::cout << DemangledTypeName(predExtend) << std::endl;
1055  //std::cout << "Checking if PredictionExtendToPeripheral exists" << std::endl;
1056  if (!predExtend)
1057  abort();
1058 
1059  // Now we want to get the core part of this
1060  IPrediction* predCore = predExtend->GetCore();
1061  //std::cout << "Checking if IPrediction (core) exists" << std::endl;
1062  if (!predCore)
1063  abort();
1064 
1065  // Now we want to cast this as a PredictionExtrap
1066  //std::cout << DemangledTypeName(predCore) << std::endl;
1067  predExtrap = dynamic_cast<PredictionExtrap*>(predCore);
1068  //std::cout << "Checking if PredictionExtrap exists" << std::endl;
1069  //std::cout << DemangledTypeName(predExtrap) << std::endl;
1070  }
1071  if (!predExtrap)
1072  abort();
1073  // Get the Extrap object
1074  IExtrap* extrap = predExtrap->GetExtrap();
1075 
1076  // Cast as modular extrap
1077  ModularExtrap* modExtrap = dynamic_cast<ModularExtrap*>(extrap);
1078 
1079  // Get mod extrap components
1080  std::vector<ModularExtrapComponent*> modExtrapComps =
1081  modExtrap->GetModExtrapComponents();
1082 
1083  // Get decomp for EEextrap
1084  // Use vector index to select appropriate set of extrap components
1085  std::cout << std::endl;
1086  std::cout << "*************************************" << std::endl;
1087  std::cout << "Getting modular extrap component: "
1088  << modExtrapComp << std::endl;
1089  std::cout << "*************************************" << std::endl;
1090  std::cout << std::endl;
1091  const IDecomp* decomp = modExtrapComps[modExtrapComp]->GetDecomp();
1092 
1093  return decomp;
1094  }
1095 
1096 //////////////////////////////////////////////////////////////////////
1097 
1099  PredictionSystJoint2018* predictionSyst, const ISyst* syst,
1100  EModExtrapComps modExtrapComp, TString plottitle, TString out_name,
1101  TString tag, AxisType pidaxis, bool printtable=true)
1102  {
1103  // Get nominal from PredictionSystJoint2018
1104  //std::cout << "Getting Nominal prediction" << std::endl;
1105  auto pred_no = predictionSyst->GetNominalPrediction();
1106  //std::cout << "In ComparePredictionsFromVector -- pred_no is " << pred_no << std::endl;
1107  // Get Nominal decomp
1108  const IDecomp* nomDecomp = GetDecomp(pred_no, modExtrapComp);
1109  //std::cout << "Got nomDecomp" << std::endl;
1110  DataMCComponents nomNDComps = GetNDComponents(nomDecomp);
1111  //std::cout << "Got nomNDComps" << std::endl;
1112  // Get shifted predictions from PredictionSystJoint2018
1113  const PredictionInterp::ShiftedPreds& shiftedPred =
1114  predictionSyst->GetShiftedPrediction(syst);
1115 
1116  //std::cout << "ShiftedPreds name: " << shiftedPred.systName << std::endl;
1117 
1118  int numSh = shiftedPred.preds.size();
1119  //std::cout << "Number of shifts: " << numSh << std::endl;
1120  //for (auto shift:shiftedPred.preds) { std::cout << "Func shift: " << shift << std::endl; }
1121 
1122  std::vector<IPrediction*> pred_sh = shiftedPred.preds;
1123 
1124  // create vecor of components
1125  std::vector<DataMCComponents> comps;
1126  for (auto pred:pred_sh)
1127  {
1128  const IDecomp* shiftedDecomp = GetDecomp(pred, modExtrapComp);
1129  DataMCComponents shiftedComps = GetNDComponents(shiftedDecomp);
1130  comps.push_back(shiftedComps);
1131  }
1132 
1133  //std::cout << "comps has length: " << comps.size() << std::endl;
1134 
1135  if (numSh == 3)
1136  {
1137  // shift 0 = -1
1138  // shift 1 = Nominal
1139  // shift 2 = +1
1140  //std::cout << "Entering shift = 3" << std::endl;
1142  nomNDComps, comps[2], comps[0], plottitle,
1143  out_name+"_1sigma", "#pm1#sigma", pidaxis, printtable);
1144  }
1145  else if (numSh == 5)
1146  {
1147  // shift 1 = -1
1148  // shift 2 = Nominal
1149  // shift 3 = +1
1150  //std::cout << "Entering shift = 5" << std::endl;
1152  nomNDComps, comps[3], comps[1], plottitle,
1153  out_name+"_1sigma", "#pm1#sigma", pidaxis, printtable);
1154  }
1155  else if (numSh == 2)
1156  {
1157  // shift 0 = Nominal
1158  // shift 1 = +1
1159  //std::cout << "Entering shift = 2" << std::endl;
1161  nomNDComps, comps[1], plottitle, out_name+"_1sigma",
1162  "#plus1#sigma", pidaxis, printtable);
1163  }
1164  }
1165 
1166 //////////////////////////////////////////////////////////////////////
1167 // Compare FD Predictions
1168 //////////////////////////////////////////////////////////////////////
1169 
1170  void CompareOneShiftPred(IPrediction * pred_no, IPrediction * pred_sh,
1171  TString plottitle, TString out_name, TString pidtag,
1172  bool printtable=true, AxisType pidaxis = kNue1bin, bool PrintParams=false){
1173 
1174  auto calc = GetCalculator (usedumb);
1175  auto h_no = GetFDMCComponents (calc, pred_no, out_name);
1176  auto h_sh = GetFDMCComponents (calc, pred_sh, out_name, 2);
1177 
1178  h_no.sig->SetMaximum(FindLimitY({h_no.sig,h_sh.sig}));
1179  h_no.bkg->SetMaximum(FindLimitY({h_no.bkg,h_sh.bkg}));
1180  h_no.sig->SetTitle(plottitle);
1181  h_no.bkg->SetTitle(plottitle);
1182  TGaxis::SetMaxDigits(3);
1183 
1184  auto leg = DefaultNueLegend(false,true,false);
1185  if (isNumuAna) leg = DefaultNumuLegend(false,true,false);
1186  leg->AddEntry(DummyGraph(kGray+1,2),"Shifted","l");
1187 
1188  std::vector<TString> legnames = {"#nu_{e} Signal","Shifted"};
1189  if (isNumuAna) legnames = {"#nu_{#mu} Signal","Shifted"};
1190  auto leg2 = CustomLegend({h_no.sig, h_sh.sig},legnames,{"l","l"});
1191 
1192  if(pidaxis == kNumuOpt){ //fix the plot but don't mess up the tables
1193 
1194  PlotMCComponentsComparison({ HistNumuOpt(h_no.bkg), HistNumuOpt(h_no.beam), HistNumuOpt(h_no.nc), HistNumuOpt(h_no.numu) },
1195  { HistNumuOpt(h_sh.bkg), HistNumuOpt(h_sh.beam), HistNumuOpt(h_sh.nc), HistNumuOpt(h_sh.numu) },
1196  leg,pidtag,out_name+"_bkg",printratioplots,false,pidaxis);
1198  { HistNumuOpt(h_sh.sig) },
1199  leg2,pidtag,out_name+"_sig",printratioplots,false,pidaxis);
1200  }
1201  else{
1202  PlotMCComponentsComparison({ h_no.bkg, h_no.beam, h_no.nc, h_no.numu },
1203  { h_sh.bkg, h_sh.beam, h_sh.nc, h_sh.numu },
1204  leg,pidtag, out_name+"_bkg", printratioplots, false, pidaxis);
1205  PlotMCComponentsComparison({ h_no.sig },
1206  { h_sh.sig },
1207  leg2,pidtag,out_name+"_sig",printratioplots,false,pidaxis);
1208  }
1209 
1210  if(PrintParams) PrintOscilationParameters (calc, usedumb);
1211 
1212  if(printtable){
1213  std::vector<TString> labels = {"$\\nu_e$ signal","Tot beam bkg","Beam $\\nu_{e}$ CC","NC", "$\\nu_{\\mu}$ CC", "$\\nu_{\\tau}$ CC"};
1214 
1215  if (isNumuAna) {
1216  labels = {"$\\nu_{\\mu}$ signal","Tot beam bkg","$\\nu_{e}$ CC","NC", "$\\nu_{\\mu} $App", "$\\nu_{\\tau}$ CC"};
1217  }
1218  ComparisonTableOne ({h_no.sig, h_no.bkg, h_no.beam, h_no.nc, h_no.numu, h_no.tau},
1219  {h_sh.sig, h_sh.bkg, h_sh.beam, h_sh.nc, h_sh.numu, h_sh.tau},
1220  labels,
1221  out_name
1222  );
1223  int nbins = 1;
1224  if (pidaxis== kNue3bin) nbins = 3;
1225  ComparisonTableOneNbins({h_no.sig, h_no.bkg, h_no.beam,h_no.nc,h_no.numu,h_no.tau},
1226  {h_sh.sig, h_sh.bkg, h_sh.beam,h_sh.nc,h_sh.numu,h_sh.tau},
1227  labels,
1228  out_name,
1229  nbins);
1230  }//end if table
1231  }
1232 //////////////////////////////////////////////////////////////////////
1233 
1234  void CompareTwoShiftPred(
1235  IPrediction * pred_no, IPrediction * pred_pl, IPrediction * pred_mi,
1236  TString plottitle, TString out_name, TString pidtag,
1237  bool printtable=true, AxisType pidaxis = kNue1bin, bool PrintParams=false)
1238  {
1239 
1240  auto calc = GetCalculator(usedumb);
1241 
1242  auto hnom = GetFDMCComponents(calc, pred_no, out_name);
1243  auto hplu = GetFDMCComponents(calc, pred_pl, out_name);
1244  auto hmin = GetFDMCComponents(calc, pred_mi, out_name);
1245 
1246  FormatErrorBand(hplu.bkg, hmin.bkg, false);
1247  hnom.bkg->SetMaximum(FindLimitY({hnom.bkg, hplu.bkg, hmin.bkg}));
1248  hnom.bkg->SetTitle(plottitle);
1249 
1250  FormatErrorBand(hplu.sig, hmin.sig, true);
1251  hnom.sig->SetMaximum(FindLimitY({hnom.sig, hplu.sig, hmin.sig}));
1252  hnom.sig->SetTitle(plottitle);
1253 
1254  TGaxis::SetMaxDigits(3);
1255 
1256  TLegend * leg, *leg2;
1257  if( pidaxis == kNumuOpt || pidaxis == kNumuOther)
1258  {
1259  leg = DefaultNumuLegend(false, true, false);
1260  leg2 = CustomLegend({hnom.sig, hplu.sig},
1261  {"#nu_{#mu} Signal","Shifted"},
1262  {"l","f"});
1263  }
1264  else
1265  {
1266  leg = DefaultNueLegend(false, true, false);
1267  leg2 = CustomLegend({hnom.sig, hplu.sig},
1268  {"#nu_{e} Signal","Shifted"},
1269  {"l","f"});
1270  }
1271  leg->AddEntry(hplu.bkg,"Syst.","f");
1272 
1273  std::vector <TH1*> nom_bkg_plot = {hnom.bkg, hnom.beam, hnom.nc, hnom.numu, hnom.tau};
1274  std::vector <TH1*> shi_bkg_plot = {hplu.bkg, hmin.bkg};
1275  std::vector <TH1*> nom_sig_plot = {hnom.sig};
1276  std::vector <TH1*> shi_sig_plot = {hplu.sig, hmin.sig};
1277 
1278  std::vector <TH1*> vnom = {
1279  hnom.sig, hnom.bkg, hnom.beam, hnom.nc, hnom.numu, hnom.tau};
1280  std::vector <TH1*> vplu = {
1281  hplu.sig, hplu.bkg, hplu.beam, hplu.nc, hplu.numu, hplu.tau};
1282  std::vector <TH1*> vmin = {
1283  hmin.sig, hmin.bkg, hmin.beam, hmin.nc, hmin.numu, hmin.tau};
1284 
1285  std::vector <TString> labels = {
1286  "$\\nu_e$ signal", "Tot beam bkg", "Beam $\\nu_{e}$ CC",
1287  "NC", "$\\nu_{\\mu}$ CC", "$\\nu_{\\tau}$ CC"};
1288  if(pidaxis == kNumuOpt || pidaxis == kNumuOther)
1289  {
1290  labels = {"$\\nu_{\\mu}$ signal", "Tot beam bkg", "NC",
1291  "$\\nu_{\\mu}$ App", "$\\nu_{e}$ CC", "$\\nu_{\\tau}$ CC"};
1292  vnom = {hnom.sig, hnom.bkg, hnom.nc, hnom.numu, hnom.beam, hnom.tau};
1293  vplu = {hplu.sig, hplu.bkg, hplu.nc, hplu.numu, hnom.beam, hplu.tau};
1294  vmin = {hmin.sig, hmin.bkg, hmin.nc, hmin.numu, hnom.beam, hmin.tau};
1295  }
1296  if(pidaxis == kNumuOpt)
1297  {
1298  nom_bkg_plot = {HistNumuOpt(hnom.bkg), HistNumuOpt(hnom.beam), HistNumuOpt(hnom.nc), HistNumuOpt(hnom.numu), HistNumuOpt(hnom.tau)};
1299  shi_bkg_plot = {HistNumuOpt(hplu.bkg), HistNumuOpt(hmin.bkg)};
1300  nom_sig_plot = {HistNumuOpt(hnom.sig)};
1301  shi_sig_plot = {HistNumuOpt(hplu.sig), HistNumuOpt(hmin.sig)};
1302  }
1303 
1304  PlotMCComponentsErrorBand(nom_bkg_plot, shi_bkg_plot,
1305  leg, pidtag, out_name + "_bkg", printratioplots, false, pidaxis);
1306  PlotMCComponentsErrorBand(nom_sig_plot, shi_sig_plot,
1307  leg2, pidtag, out_name + "_sig", printratioplots, false, pidaxis);
1308 
1309  if(PrintParams) PrintOscilationParameters (calc, usedumb);
1310 
1311  if(printtable){
1312 
1313  ComparisonTable (vnom, vplu, vmin, labels, out_name);
1314 
1315  int nbins = 1;
1316  if (pidaxis== kNue3bin) nbins = 3;
1317  ComparisonTableNbins (vnom, vplu, vmin, labels, out_name, nbins);
1318  }
1319  }
1320 
1321 //////////////////////////////////////////////////////////////////////
1322 
1324  PredictionSystJoint2018* predictionSyst, const ISyst* syst,
1325  TString plottitle, TString out_name,
1326  bool printtable=true, AxisType pidaxis=kNue3bin)
1327  {
1328  // Get nominal from PredictionSystJoint2018
1329  auto pred_no = predictionSyst->GetNominalPrediction();
1330  //std::cout << "In ComparePredictionsFromVector -- pred_no is " << pred_no << std::endl;
1331 
1332  // Get shifted predictions from PredictionSystJoint2018
1333  const PredictionInterp::ShiftedPreds& shiftedPred =
1334  predictionSyst->GetShiftedPrediction(syst);
1335 
1336  //std::cout << "ShiftedPreds name: " << shiftedPred.systName << std::endl;
1337 
1338  int numSh = shiftedPred.preds.size();
1339  //std::cout << "Number of shifts: " << numSh << std::endl;
1340  //for (auto shift:shiftedPred.preds) { std::cout << "Func shift: " << shift << std::endl; }
1341 
1342  std::vector<IPrediction*> pred_sh = shiftedPred.preds;
1343 
1344  if (numSh == 3)
1345  {
1346  // Just get rid of this.
1347  //if (out_name.Contains("RelativeCalib"))
1348  //{
1349  // shift 0 = -2
1350  // shift 1 = Nominal
1351  // shift 2 = +2
1352  // FOR FD only --> need
1353  //CompareTwoShiftPred(
1354  // pred_sh[1], pred_sh[0], pred_sh[2], plottitle,
1355  // out_name+"_2sigma", "#pm2#sigma", printtable, pidaxis);
1356  //}
1357  //else
1358  {
1359  // shift 0 = -1
1360  // shift 1 = Nominal
1361  // shift 2 = +1
1362  //std::cout << "Passing; " << pred_no << ", " << pred_sh[0] << ", " << pred_sh[2] << std::endl;
1364  pred_no, pred_sh[0], pred_sh[2], plottitle,
1365  out_name+"_1sigma", "#pm1#sigma", printtable, pidaxis);
1366  }
1367  }
1368  else if (numSh == 5) {
1369  // shift 1 = -1
1370  // shift 2 = Nominal
1371  // shift 3 = +1
1372  //std::cout << "Passing; " << pred_no << ", " << pred_sh[1] << ", " << pred_sh[3] << std::endl;
1374  pred_no, pred_sh[1], pred_sh[3], plottitle,
1375  out_name+"_1sigma", "#pm1#sigma", printtable, pidaxis);
1376  } else if (numSh == 2) {
1377  // shift 0 = Nominal
1378  // shift 1 = +1
1379  //std::cout << "Passing; " << pred_no << ", " << pred_sh[1] << std::endl;
1381  pred_no, pred_sh[1], plottitle, out_name+"_1sigma",
1382  "#plus1#sigma", printtable, pidaxis);
1383  }
1384  }
1385 }
struct DataMCComponets GetFDMCComponents(osc::IOscCalc *calc, IPrediction *pred_no, TString output_name="nue", int linestyle=1, bool bkgdetails=false)
const double kLegendDefaultY2
enum BeamMode kOrange
const double kLegendDefaultY1
const int pidBins
enum BeamMode kRed
string output_name
Pickle.
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
void ComparisonTableOneNbins(const std::vector< TH1 * > &mcnom, const std::vector< TH1 * > &mcshift, const std::vector< TString > &labels, TString ltxcommand, int N=3)
const IDecomp * GetDecomp(IPrediction *prediction, EModExtrapComps modExtrapComp)
virtual Spectrum AntiNueComponent() const =0
double ssth23
set< int >::iterator it
osc::IOscCalc * GetCalculator(bool usingdumb)
Antineutrinos-only.
Definition: IPrediction.h:50
bool mergePeripheral
Float_t y1[n_points_granero]
Definition: compare.C:5
(&#39; appearance&#39;)
Definition: IPrediction.h:18
void CompareOneShiftPred(IPrediction *pred_no, IPrediction *pred_sh, TString plottitle, TString out_name, TString pidtag, bool printtable=true, AxisType pidaxis=kNue1bin, bool PrintParams=false)
TH1D * ToTH1(double exposure, Color_t col=kBlack, Style_t style=kSolid, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Histogram made from this Spectrum, scaled to some exposure.
Definition: Spectrum.cxx:148
virtual Spectrum NumuComponent() const =0
Loads shifted spectra from files.
const FitDmSq32 kFitDmSq32
Definition: FitVars.cxx:18
const double kLegendDefaultX2
Float_t x1[n_points_granero]
Definition: compare.C:5
void PlotMCComponentsComparison(std::vector< TH1 * >vnom, std::vector< TH1 * >vshi, TLegend *leg, TString pidtag="", TString out_name="plot_nd", bool ratioplot=false, bool ratioerrors=false, bool pidaxis=0)
TH1F * h3
Definition: berger.C:36
void SetValue(osc::IOscCalcAdjustable *osc, double val) const override
Definition: FitVars.cxx:177
(&#39;beam &#39;)
Definition: IPrediction.h:15
const char * p
Definition: xmltok.h:285
const Color_t kNueSignalWSColor
Definition: Style.h:20
TH1 * ratio(TH1 *h1, TH1 *h2)
General interface to oscillation calculators.
Definition: StanTypedefs.h:22
virtual Spectrum Data_Component() const
Definition: IDecomp.cxx:16
TLegend * CustomLegend(std::vector< TH1 * >h, std::vector< TString > title, std::vector< TString > option, double x1=kLegendDefaultX1, double y1=kLegendDefaultY1, double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
void ComparisonTableOne(const std::vector< TH1 * > &mcnom, const std::vector< TH1 * > &mcshift, const std::vector< TString > &labels, const TString ltxcommand="")
void SetValue(osc::IOscCalcAdjustable *osc, double val) const override
Forward to wrapped Var&#39;s SetValue()
TLegend * leg1
Definition: plot_hist.C:105
void CompareTwoShiftPred(IPrediction *pred_no, IPrediction *pred_pl, IPrediction *pred_mi, TString plottitle, TString out_name, TString pidtag, bool printtable=true, AxisType pidaxis=kNue1bin, bool PrintParams=false)
TH1F * hplus
Definition: Xsec_final.C:75
void CenterTitles(TH1 *histo)
Definition: Plots.cxx:1481
TGraphAsymmErrors * GraphWithPoissonErrors(const TH1 *h, bool noErrorsXaxis, bool drawEmptyBins)
Calculate statistical errors appropriate for small Poisson numbers.
Definition: Plots.cxx:910
const Color_t kTotalMCErrorBandColor
Definition: Style.h:17
osc::OscCalcDumb calc
osc::IOscCalcAdjustable * DefaultOscCalc()
Create a new calculator with default assumptions for all parameters.
Definition: Calcs.cxx:49
void PrettyTag(TString pid, int color, double xndc, double yndc)
Encapsulate code to systematically shift a caf::SRProxy.
Definition: ISyst.h:14
const Color_t kNumuBackgroundColor
Definition: Style.h:30
void PlotMCComponentsErrorBand(std::vector< TH1 * >vnom, std::vector< TH1 * >vshi, TLegend *leg, TString pidtag="", TString out_name="plot_nd", bool ratioplot=false, bool ratioerrors=false, bool pidaxis=0)
const ShiftedPreds & GetShiftedPrediction(const ISyst *syst)
_IOscCalc< double > IOscCalc
Definition: IOscCalc.h:41
double FindLimitY(std::vector< TH1 * > histos, bool doMax=true)
struct DataMCComponets GetNDComponents(TDirectory *d_no, double kNDPOT, int linestyle=1)
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
const double kAna2018RHCPOT
Definition: Exposures.h:208
const Color_t kNueSignalColor
Definition: Style.h:19
void PlotNDDataTotalMCComparison(TH1 *hdata, std::vector< TH1 * > htots, TLegend *leg, TString pidtag="", TString out_name="plot_nd", bool ratioplot=false, bool pidaxis=false)
TLegend * DefaultNueLegend(bool sig=true, bool tot=true, bool data=false, double x1=kLegendDefaultX1, double y1=kLegendDefaultY1, double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
void PlotMCComponentsErrorBand(std::vector< TH1 * >vnom, std::vector< TH1 * >vshi, TLegend *leg, TString pidtag="", TString out_name="plot_nd", bool ratioplot=false, bool ratioerrors=false, AxisType pidaxis=kNue1bin)
const XML_Char const XML_Char * data
Definition: expat.h:268
c2
Definition: demo5.py:33
static bool isFHC
const Color_t kNuTauBackgroundColor
Definition: Style.h:26
const int nbins
Definition: cellShifts.C:15
Charged-current interactions.
Definition: IPrediction.h:39
TLegend * DefaultNumuLegend(bool sig=true, bool tot=true, bool data=false, double x1=kLegendDefaultX1, double y1=kLegendDefaultY1, double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
void CompareNDDataMCFromVector(TDirectory *d_no, std::vector< TDirectory * > d_sh, TString plottitle, TString out_name, TString tag, bool printtable=true, AxisType pidaxis=kNue1bin)
const bool printratioplots
bool isNumuAna
void PIDTag(TString pid)
void Nue2018ThreeBinDivisions(bool coreOnly, const int color, const int style)
void ComparisonTable(const std::vector< TH1 * > &mcnom, const std::vector< TH1 * > &mcplus, const std::vector< TH1 * > &mcminus, std::vector< TString > labels, TString ltxcommand="")
void FormatErrorBand(TH1 *hplus, TH1 *hminus, bool signal=false, bool fixbins=false)
void SetValue(osc::IOscCalcAdjustable *osc, double val) const override
void NuePID3bin(TH1 *h, TString pidtag, bool binlabels=true)
const Color_t kBeamNueBackgroundColor
Definition: Style.h:24
double Ebins[11]
Definition: Xdiff_gwt.C:8
#define pot
std::vector< ModularExtrapComponent * > GetModExtrapComponents() const
Definition: ModularExtrap.h:60
TGraph * DummyGraph(int color=kBlack, int style=1, int linewidth=2, int markerstyle=kFullCircle)
Simple oscillation probability calculator that has no solar term or mass hierarchy or delta so it&#39;s s...
Definition: OscCalcDumb.h:16
const double kLegendDefaultX1
std::vector< IPrediction * > preds
const double j
Definition: BetheBloch.cxx:29
double he
Definition: runWimpSim.h:113
void HTag(TString hie)
std::vector< std::string > comps
void FormatErrorBand(TH1 *hplus, TH1 *hminus, bool signal=false, bool fixbins=false)
virtual Spectrum AntiNumuComponent() const =0
void ComparePredictionsFromVector(TDirectory *dpred_no, std::vector< TDirectory * > dpred_sh, TString plottitle, TString out_name, bool printtable=true, AxisType pidaxis=kNue1bin)
float bin[41]
Definition: plottest35.C:14
void CompareNDDataOneMC(TDirectory *d_no, TDirectory *d_sh, TString plottitle, TString out_name, TString pidtag, bool printtable=true, AxisType pidaxis=kNue1bin)
Neutrinos-only.
Definition: IPrediction.h:49
void PrintLatexFigure(TString name, TString plotdir="plots/")
OStream cout
Definition: OStream.cxx:6
(&#39; survival&#39;)
Definition: IPrediction.h:19
virtual _IOscCalc * Copy() const =0
void CompareNDDataTwoMC(TDirectory *d_no, TDirectory *d_pl, TDirectory *d_mi, TString plottitle, TString out_name, TString pidtag, bool printtable=true, AxisType pidaxis=kNue1bin)
static float min(const float a, const float b, const float c)
Definition: absgeo.cxx:45
TH1F * h1
TLegend * CustomLegend(std::vector< TH1 * >h, std::vector< bool >isdata, std::vector< TString >title, double x1=kLegendDefaultX1, double y1=kLegendDefaultY1, double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
TH1 * PrettyRatio(TH1 *shi, TH1 *nom, int color, double alpha=1., TString titley="Shift / Nom")
Interface to extrapolation procedures.
Definition: IExtrap.h:8
AxisType
double dmsq32
TCanvas * ExPIDPlot(std::vector< TH1 * > topHistos, std::vector< TString > topOption)
const IExtrap * GetExtrap() const
TH1F * hminus
Definition: Xsec_final.C:76
void PrintOscilationParameters(osc::IOscCalc *calc, bool usingdumb)
TH1 * HistNumuOpt(TH1 *orig)
enum BeamMode kViolet
void ComparisonTableNbins(const std::vector< TH1 * > &mcnom, const std::vector< TH1 * > &mcplus, const std::vector< TH1 * > &mcminus, std::vector< TString > labels, TString ltxcommand="", int N=3)
virtual Spectrum NCTotalComponent() const
Definition: IDecomp.h:18
TH1 * RebinToShortAxis(TH1 *h)
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
Standard interface to all decomposition techniques.
Definition: IDecomp.h:13
void Nue2018ThreeBinLabels(const double yNDC, const double textsize, const int color, const bool nd)
Neutral-current interactions.
Definition: IPrediction.h:40
void PlotNDDataTotalMCComparison(TH1 *hdata, std::vector< TH1 * > htots, TLegend *leg, TString pidtag="", TString out_name="plot_nd", bool ratioplot=false, AxisType pidaxis=kNue1bin)
Both neutrinos and antineutrinos.
Definition: IPrediction.h:52
c1
Definition: demo5.py:24
Standard interface to all prediction techniques.
Definition: IPrediction.h:57
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
const Color_t kTotalMCColor
Definition: Style.h:16
void DefaultFormatNue(DataMCComponets &comp, int linestyle=1)
const double kAna2018FHCPOT
Definition: Exposures.h:207
Take the output of an extrapolation and oscillate it as required.
All neutrinos, any flavor.
Definition: IPrediction.h:26
(&#39; appearance&#39;)
Definition: IPrediction.h:16
TLegend * DefaultNueLegend(bool sig=true, bool tot=true, bool data=false, double x1=kLegendDefaultX1, double y1=kLegendDefaultY1, double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
bool usedumb
Extrapolate each component using a separate ModularExtrapComponent.
Definition: ModularExtrap.h:23
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
const FitVarWithPrior fitDeltaInPiUnits_UniformPriordCP & kFitDeltaInPiUnits
const Color_t kNCBackgroundColor
Definition: Style.h:22
enum BeamMode kBlue
void PlotDataMC(std::vector< TH1 * >vnom, std::vector< bool >isdata, TLegend *leg, TString pidtag="", TString htag="", TString out_name="plot_FD", bool ratioplot=false, bool ratioerrors=false, AxisType pidaxis=kNue1bin)
virtual IOscCalc * Copy() const override
Definition: OscCalcDumb.h:19
void PlotMCComponentsComparison(std::vector< TH1 * >vnom, std::vector< TH1 * >vshi, TLegend *leg, TString pidtag="", TString out_name="plot_nd", bool ratioplot=false, bool ratioerrors=false, AxisType pidaxis=kNue1bin)
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:29
TCanvas * RatioPlot(std::vector< TH1 * > topHistos, std::vector< TString > topOption, std::vector< TH1 * > bottomHistos, std::vector< TString > bottomOption, TString pidtag, AxisType pidaxis=kNue1bin)
void FixLegend(TLegend *leg, TString opt="default")
void PrintCanvasAll(TCanvas *c, const TString outname, const TString outdir="plots/", const std::vector< TString > exts={"pdf","png","root"})
const std::string outdir
virtual Spectrum NueComponent() const =0
def sign(x)
Definition: canMan.py:197
unsigned int uint