draw_plots_util.h
Go to the documentation of this file.
4 #include "CAFAna/Core/Binning.h"
5 #include "CAFAna/Core/Cut.h"
6 #include "CAFAna/Core/HistAxis.h"
8 #include "CAFAna/Core/Spectrum.h"
11 #include "CAFAna/Cuts/Cuts.h"
12 #include "CAFAna/Cuts/NumuCuts.h"
14 #include "CAFAna/Cuts/SpillCuts.h"
16 #include "CAFAna/Vars/NueVars.h"
17 #include "CAFAna/Vars/Vars.h"
19 
20 #include "TStyle.h"
21 #include "TFile.h"
22 #include "TH1.h"
23 #include "TCanvas.h"
24 #include "TGraph.h"
25 #include "TGraphAsymmErrors.h"
26 #include "TGaxis.h"
27 #include "TLatex.h"
28 #include "TLegend.h"
29 #include "TLine.h"
30 #include "TSystem.h"
31 #include <iostream>
32 #include <fstream>
33 #include <iomanip>
34 
35 namespace osc{
36  class OscCalculatorDumb;
37 };
38 
41 const bool printratioplots=true;
42 
43 
44 // Put a "NOvA Preliminary" tag in the corner
45 void Prelim()
46 {
47  TLatex* prelim = new TLatex(.9, .95, "NOvA Preliminary");
48  prelim->SetTextColor(kBlue);
49  prelim->SetNDC();
50  prelim->SetTextSize(2/36.);
51  prelim->SetTextAlign(32);
52  prelim->Draw();
53 }
54 
55 void CenterTitles(TH1* histo)
56 {
57  histo->GetXaxis()->CenterTitle();
58  histo->GetYaxis()->CenterTitle();
59  histo->GetZaxis()->CenterTitle();
60 }
61 
62 void PIDTag(TString pid)
63 {
64  TLatex* tag = new TLatex(.83, .92, pid);
65  tag->SetTextColor(kAzure+3);
66  tag->SetNDC();
67  tag->SetTextSize(0.038);
68  tag->SetTextFont(42);
69  tag->SetTextAlign(12);
70  tag->Draw();
71 }
72 
73 void HTag(TString pid)
74 {
75  TLatex* tag = new TLatex(.18, .825, pid);
76  if ( pid == "IH") tag->SetTextColor(kOrange+2);
77  if ( pid == "NH") tag->SetTextColor(kAzure+3);
78  tag->SetNDC();
79  tag->SetTextSize(0.043);
80  tag->SetTextFont(42);
81  tag->SetTextAlign(12);
82  tag->Draw();
83 }
84 
85 // NB, loses all the style info, so call early
86 TH1* RebinToShortAxis(TH1* h)
87 {
88  TH1* ret = new TH1D(UniqueName().c_str(), h->GetTitle(),
89  18, 0, 18);
90  ret->GetXaxis()->SetTitle(h->GetXaxis()->GetTitle());
91  ret->GetYaxis()->SetTitle(h->GetYaxis()->GetTitle());
92 
93  std::vector<std::pair<int, int>> binMap = {{3, 2},
94  {4, 3},
95  {5, 4},
96  {6, 5}, // end of 1st spectrum
97  {13, 8},
98  {14, 9},
99  {15, 10},
100  {16, 11}, // end of 2nd
101  {23, 14},
102  {24, 15},
103  {25, 16},
104  {26, 17}};
105 
106  for(auto it: binMap) ret->SetBinContent(it.second, h->GetBinContent(it.first));
107 
108  return ret;
109 }
110 
112  // Separators between the different PID ranges
113  // TGraph not TLine because TGraph is properly clipped by frame
114  TGraph* l1 = new TGraph;
115  l1->SetPoint(0, 10, 0);
116  l1->SetPoint(1, 10, 50000);
117  TGraph* l2 = new TGraph;
118  l2->SetPoint(0, 20, 0);
119  l2->SetPoint(1, 20, 50000);
120  for(TGraph* l: {l1, l2}){
121  l->SetLineWidth(1);
122  l->SetLineColor(kGray+3);
123  l->Draw("l same");
124  }
125 }
126 
128  // Put the labels near the top of the plot
129  gPad->Update();
130  const double y = gPad->GetUymax() * 0.955;
131 
132  const ana::Binning* bins = 0;
133  if(pid == "LID") bins = &ana::kLIDNLBinning;
134  if(pid == "LEM") bins = &ana::kLEMNLBinning;
135  if(pid == "CVN") bins = &ana::kCVNNLBinning;
136  assert(bins);
137 
138  for(int i = 0; i < 3; ++i){
139  TLatex* ltx = new TLatex(3+6*i, y, TString::Format("%g < %s < %g",
140  bins->Edges()[i],
141  pid.c_str(),
142  bins->Edges()[i+1]));
143  ltx->SetTextSize(0.032);
144  ltx->SetTextColor(kGray+3);
145  ltx->SetTextAlign(22);
146  ltx->Draw();
147  }
148 
149 }
150 
151 void PID2DAxis(TH1* axes, bool third = false)
152 {
153 
154  double xmax = axes->GetXaxis()->GetXmax();
155  double xmin = axes->GetXaxis()->GetXmin();
156  double xrange = xmax-xmin;
157 
158  double ymin = axes->GetMinimum();
159  TGaxis* gax1 = new TGaxis(0, ymin, 8, ymin, 0, 4, 4, "");
160  TGaxis* gax2 = new TGaxis(10, ymin, 18, ymin, 0, 4, 4, "");
161  TGaxis* gax3 = new TGaxis(20, ymin, 30, ymin, 0, 5, 5, "");
162 
163  float size = axes->GetXaxis()->GetLabelSize();
164  float offs = axes->GetXaxis()->GetLabelOffset();
165  int font = axes->GetXaxis()->GetLabelFont();
166  // And then hide its title
167  axes->GetXaxis()->SetLabelSize(0);
168  // and make all the tick marks the same size
169  axes->GetXaxis()->SetNdivisions(3000);
170  axes->GetXaxis()->SetTitle("Reconstructed Neutrino Energy (GeV)");
171  axes->GetYaxis()->SetTitle("Events / 6.05 x 10^{20} POT-equiv");
172  axes->GetYaxis()->SetTitleOffset(0.8);
173  axes->SetTitle("");
174  CenterTitles(axes);
175 
176  if ( third ){
177  gax1->SetLabelSize(size);
178  gax1->SetLabelOffset(offs);
179  gax1->SetLabelFont(font);
180  gax1->Draw();
181  }
182  else
183  for(TGaxis* ax: {gax1, gax2, gax3}){
184  // Copy settings from pre-existing axis
185  ax->SetLabelSize(size-0.001);
186  ax->SetLabelOffset(offs);
187  ax->SetLabelFont(font);
188  ax->Draw();
189  }
190 
191  TLine* xline = new TLine();
192  xline->SetLineWidth(2);
193  xline->DrawLine(xmin, 0, xmax, 0);
194  gPad->RedrawAxis();
195 }
196 //----------------------------------------------------------------------
198 {
199  // Separators between the different PID ranges
200  // TGraph not TLine because TGraph is properly clipped by frame
201  TGraph* l1 = new TGraph;
202  l1->SetPoint(0, 6, 0);
203  l1->SetPoint(1, 6, 50000);
204  TGraph* l2 = new TGraph;
205  l2->SetPoint(0, 12, 0);
206  l2->SetPoint(1, 12, 50000);
207  for(TGraph* l: {l1, l2}){
208  // l->SetLineStyle(3);
209  l->SetLineWidth(1);
210  l->SetLineColor(kGray+3);
211  l->Draw("l same");
212  }
213 
214  // Hardcoding maxy is obviously unsustainable
215  TGaxis* gax1 = new TGaxis(6, 0, 6, 18, 0, 18, 504, "U+-");
216  TGaxis* gax2 = new TGaxis(12, 0, 12, 18, 0, 18, 504, "U+-");
217 }
218 //----------------------------------------------------------------------
220 {
221  const std::string fPID = "CVN"; // HACK
222 
223  // Put the labels near the top of the plot
224  std::string pid = (std::string)fPID;
225  gPad->Update();
226  const double y = gPad->GetUymax() * .94;
227 
228  const ana::Binning* bins = 0;
229  if(fPID == "LID") bins = &ana::kLIDNLBinning;
230  if(fPID == "LEM") bins = &ana::kLEMNLBinning;
231  if(fPID == "CVN") bins = &ana::kCVNNLBinning;
232  assert(bins);
233 
234  for(int i = 0; i < 3; ++i){
235  TLatex* ltx = new TLatex(3+6*i, y, TString::Format("%g < %s < %g",
236  bins->Edges()[i],
237  pid.c_str(),
238  bins->Edges()[i+1]));
239  ltx->SetTextSize(.04);//0.035);
240  ltx->SetTextAlign(22);
241  ltx->Draw();
242  }
243 }
245 {
246  double ymin = hist->GetMinimum();
247  // Paint our own x-axis fragments
248  TGaxis* gax1 = new TGaxis(1, ymin, 5, ymin, 1, 3, 3, "");
249  TGaxis* gax2 = new TGaxis(7, ymin, 11, ymin, 1, 3, 3, "");
250  TGaxis* gax3 = new TGaxis(13, ymin, 17, ymin, 1, 3, 3, "");
251 
252  float tick = hist->GetXaxis()->GetTickLength();
253  float size = hist->GetXaxis()->GetLabelSize();
254  float offs = hist->GetXaxis()->GetLabelOffset();
255  int font = hist->GetXaxis()->GetLabelFont();
256  // And then hide its title
257  hist->GetXaxis()->SetLabelSize(0);
258  // and make all the tick marks the same size
259  hist->GetXaxis()->SetNdivisions(1800);
260  hist->GetXaxis()->SetTitle("Reconstructed neutrino energy (GeV)");
261  hist->GetYaxis()->SetTitle("Events / 0.5 GeV Bin");
262  hist->GetYaxis()->SetTitleOffset(0.8);
263  hist->GetXaxis()->CenterTitle();
264  hist->GetYaxis()->CenterTitle();
265  hist->SetTitle("");
266 
267  if (false){// fThird ){
268  gax1->SetTickSize(tick);
269  gax1->SetLabelSize(size);
270  gax1->SetLabelOffset(offs);
271  gax1->SetLabelFont(font);
272  gax1->Draw();
273  }
274  else{
275  for(TGaxis* ax: {gax1, gax2, gax3}){
276  // Copy settings from pre-existing axis
277  ax->SetTickSize(tick);
278  ax->SetLabelSize(size);
279  ax->SetLabelOffset(offs);
280  ax->SetLabelFont(font);
281  ax->Draw();
282  }
283  }
284 }
285 
286 double FindLimitY (std::vector<TH1*> histos, bool doMax = true)
287 {
288  double max = histos[0]->GetMaximum();
289  double min = histos[0]->GetMinimum();
290 
291  for(unsigned int i=0; i<histos.size(); i++){
292  if ( histos[i]->GetMaximum() > max ) max = histos[i]->GetMaximum();
293  if ( histos[i]->GetMinimum() < min ) min = histos[i]->GetMinimum();
294  }
295 
296  if ( doMax ) return max;
297  else return min;
298 
299 }
300 
301 namespace ana
302 {
303 
304  void ComparisonTable(std::vector<TH1*> mcnom, std::vector<TH1*> mcshift, std::vector<TString> labels){
305  unsigned int nentries = mcnom.size();
306  if(mcshift.size()!=nentries || labels.size()!=nentries){std::cerr << "incompatible sizes\n"; return;}
307  std::cout << "\\begin{tabular}{ l | r r r } \n"
308  << std::setw(20) <<" " << " & "
309  << std::setw(10) << "Nominal" << " & "
310  << std::setw(10) << "Shift" << " & "
311  << std::setw(12) << "\\% Diff."
312  <<"\t \\\\ \\hline \n";
313  for(unsigned int i=0;i<nentries;i++){
314  std::cout.unsetf(std::ios_base::floatfield);
315  std::cout << std::setprecision(3)
316  << std::setw(20) << labels[i] << " & "
317  << std::setw(10) << mcnom[i] ->Integral() << " & "
318  << std::setw(10) << mcshift[i]->Integral() << " & "
319  << std::setw(10) << std::setprecision(1) << std::fixed
320  << ( mcnom[i]->Integral()>0 ? (((mcshift[i]->Integral()/mcnom[i]->Integral()) -1.)*100.):0.)
321  << "\\% \t \\\\ \n";
322  }
323  std::cout << "\\end{tabular}\n";
324  }
325 
326  TString MakeLatexCommandName(TString mystring){
327  //There must be a better way
328  std::vector<TString> in = {"0","1","2","3","4","5",
329  "6","7","8","9","_"};
330  std::vector<TString> out = {"zero","one","two","three","four","five",
331  "six","seven","eight","nine",""};
332  for(unsigned int i=0;i<in.size();i++)
333  mystring.ReplaceAll(in[i],out[i]);
334  return mystring;
335 
336  }
337 
338  void ComparisonTable(std::vector<TH1*> mcnom, std::vector<TH1*> mcplus, std::vector<TH1*> mcminus,
339  std::vector<TString> labels,TString ltxcommand = ""){
340  unsigned int nentries = mcnom.size();
341  if(mcplus.size()!=nentries || mcminus.size()!=nentries || labels.size()!=nentries){std::cerr << "incompatible sizes\n"; return;}
342  if(ltxcommand != "") std::cout << "\\newcommand{\\" << MakeLatexCommandName(ltxcommand) << "}\n";
343  std::cout << "{ \\begin{tabular}{ l | r r r | r r } \n"
344  << std::setw(20) <<" " << " & "
345  << std::setw(10) << "Nominal" << " & "
346  << std::setw(10) << "Shift (+)" << " & "
347  << std::setw(10) << "Shift (-)" << " & "
348  << std::setw(12) << "\\% Diff. (+)" << " & "
349  << std::setw(12) << "\\% Diff. (-)"
350  <<"\t \\\\ \\hline \n";
351  for(unsigned int i=0;i<nentries;i++){
352  std::cout.unsetf(std::ios_base::floatfield);
353  std::cout << std::setprecision(3)
354  << std::setw(20) << labels[i] << " & "
355  << std::setw(10) << mcnom[i] ->Integral() << " & "
356  << std::setw(10) << mcplus[i] ->Integral() << " & "
357  << std::setw(10) << mcminus[i]->Integral() << " & "
358  << std::setw(10) << std::setprecision(1) << std::fixed
359  << ( mcnom[i]->Integral()>0 ? (((mcplus [i]->Integral()/mcnom[i]->Integral()) -1.)*100.):0.)
360  << "\\% & "
361  << std::setw(10)
362  << ( mcnom[i]->Integral()>0 ? (((mcminus[i]->Integral()/mcnom[i]->Integral()) -1.)*100.):0.)
363  << "\\% \t \\\\ \n";
364  }
365  std::cout << "\\end{tabular} } \n";
366  }
367 
368  ////////////////////////////////////////
369 
370  const double kLegendDefaultX1 = 0.63;
371  const double kLegendDefaultX2 = 0.85;
372  const double kLegendDefaultY1 = 0.68;
373  const double kLegendDefaultY2 = 0.88;
374 
375  void FixLegend(TLegend * leg, TString opt="default"){
376  leg->SetFillStyle(0);
377  if(opt=="default"){
378  leg->SetX1(kLegendDefaultX1); leg->SetX2(kLegendDefaultX2);
379  leg->SetY1(kLegendDefaultY1); leg->SetY2(kLegendDefaultY2);
380  }
381  if(opt=="ExPID"){
382  leg->SetX1(0.405); leg->SetX2(0.565);
383  leg->SetY1(0.59); leg->SetY2(0.82);
384  leg->SetBorderSize(0);
385  leg->SetTextSize(0.037);
386  }
387  if(opt=="ExPIDPads"){
388  leg->SetX1(0.55); leg->SetX2(0.89);
389  leg->SetY1(0.75); leg->SetY2(0.89);
390  leg->SetTextSize(0.045);
391  }
392  }
393  ////////////////////////////////////////
394  TLegend * DefaultNueLegend (bool sig=true, bool tot=true,bool data=false,
395  double x1=kLegendDefaultX1, double y1=kLegendDefaultY1,
396  double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
397  {
398  TLegend* leg = new TLegend(x1,y1,x2,y2);
399  leg->SetFillStyle(0);
400  TH1* hdata = new TH1F("", "", 1, 0, 1);
401  TH1* hsig = new TH1F("", "", 1, 0, 1);
402  TH1* hNC = new TH1F("", "", 1, 0, 1);
403  TH1* hCC = new TH1F("", "", 1, 0, 1);
404  TH1* hBeam = new TH1F("", "", 1, 0, 1);
405  TH1* htot = new TH1F("", "", 1, 0, 1);
406  hdata->SetMarkerStyle(kFullCircle);
407  hsig->SetLineColor(kNueSignalColor);
408  hNC->SetLineColor(kNCBackgroundColor);
409  hCC->SetLineColor(kNumuBackgroundColor);
410  hBeam->SetLineColor(kBeamNueBackgroundColor);
411  htot->SetLineColor(kTotalMCColor);
412  if(sig) leg->AddEntry(hsig,"#nu_{e} Signal","l");
413  if(data) leg->AddEntry(hdata,"ND Data","epl");
414  if(tot)leg->AddEntry(htot,"Total Bkg","l");
415  leg->AddEntry(hBeam,"Beam #nu_{e}","l");
416  leg->AddEntry(hNC,"NC","l");
417  leg->AddEntry(hCC,"#nu_{#mu} CC","l");
418  return leg;
419 
420  }
421  ////////////////////////////////////////
422  TLegend * CustomLegend (std::vector<TH1*>h, std::vector<bool>isdata,
423  std::vector<TString>title,
424  double x1=kLegendDefaultX1, double y1=kLegendDefaultY1,
425  double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
426  {
427  TLegend* leg = new TLegend(x1,y1,x2,y2);
428  leg->SetFillStyle(0);
429  for ( unsigned int v = 0; v<h.size(); ++v){
430  if (isdata[v]) TLegendEntry* le = leg->AddEntry(h[v],title[v],"epl");
431  else{
432  if (h[v]->GetLineColor() == kRed) TLegendEntry* le = leg->AddEntry(h[v],title[v],"l");
433  else TLegendEntry* le = leg->AddEntry(h[v],title[v],"f");
434  }
435  }
436  return leg;
437 
438  }
439  ////////////////////////////////////////
440  TCanvas * ExPIDPlot (std::vector<TH1*> topHistos,
441  std::vector<TString> topOption)
442  {
443  TCanvas *cpid = new TCanvas("cpid", "cpid", 410, 850);
444  float Ebins = 10;
445  int pidBins = 3;
446 
447  TPad* p[pidBins];
448  p[0] = new TPad("", "", 0, 0, 1, 1);
449  p[1] = new TPad("", "", 0, 0, 1, 1);
450  p[2] = new TPad("", "", 0, 0, 1, 1);
451  p[0]->SetTopMargin(0.1);
452  p[0]->SetBottomMargin(.634);
453  p[1]->SetTopMargin(.369);
454  p[1]->SetBottomMargin(.369);
455  p[2]->SetTopMargin(.634);
456  p[2]->SetBottomMargin(0.1);
457  p[0]->SetFillStyle(0);
458  p[1]->SetFillStyle(0);
459  p[2]->SetFillStyle(0);
460  p[0]->Draw();
461  p[1]->Draw();
462  p[2]->Draw();
463 
464  TH1* he[pidBins];
465  TH1* h[topHistos.size()][pidBins];
466  for(int k=pidBins-1; k>=0; --k){
467  he[k] = (TH1*) topHistos[0]->Clone();
468  he[k]->GetYaxis()->SetRangeUser(FindLimitY(topHistos,0) ,1.2*FindLimitY(topHistos,1));
469  he[k]->GetXaxis()->SetRangeUser(Ebins*k,Ebins*(k+1));
470  if ( k>0 ) he[k]->GetXaxis()->SetLabelSize(0);
471  p[pidBins-1-k]->cd();
472  he[k]->Draw(topOption[0]);
473  he[k]->GetYaxis()->SetTitleSize(26);
474  he[k]->GetYaxis()->SetTitleFont(43);
475  he[k]->GetYaxis()->SetTitleOffset(1.25);
476  he[k]->GetYaxis()->SetLabelFont(43);
477  he[k]->GetYaxis()->SetLabelSize(14);
478  he[k]->GetXaxis()->SetLabelOffset(-0.008);
479  he[k]->GetXaxis()->SetTitleOffset(0.5);
480  p[pidBins-1-k]->SetGridy();
481  }
482  for(unsigned int ii=0;ii<topHistos.size();ii++){
483  for(int j=pidBins-1; j>=0; --j){
484  h[ii][j] = (TH1*) topHistos[ii]->Clone();
485  h[ii][j]->GetXaxis()->SetRangeUser(Ebins*j,Ebins*(j+1));
486  p[pidBins-1-j]->cd();
487  h[ii][j]->Draw("same "+topOption[ii]);
488  }
489  }
490  for(int h=pidBins-1; h>=0; --h){
491  he[h]->SetStats(0);
492  p[pidBins-1-h]->cd();
493  he[h]->Draw("same "+topOption[0]);
494  }
495  p[0]->RedrawAxis(); //needed for white histograms e.g. beam
496  p[1]->RedrawAxis();
497  p[2]->RedrawAxis();
498  return cpid;
499  }
500  ////////////////////////////////////////
501  TCanvas * RatioPlot (std::vector<TH1*> topHistos,
502  std::vector<TString> topOption,
503  std::vector<TH1*> bottomHistos,
504  std::vector<TString> bottomOption,
505  TString pidtag,
506  bool pidaxis=false)
507  {
508  TCanvas *c = new TCanvas("c", "canvas", 600+(150*pidaxis), 700-(50*pidaxis));
509  auto h1 = (TH1*) topHistos[0]->Clone();
510  auto h3 = (TH1*) bottomHistos[0]->Clone();
511 
512  CenterTitles(h1);
513  CenterTitles(h3);
514 
515  TPad* p1 = new TPad("", "", 0, 0, 1, 1);
516  TPad* p2 = new TPad("", "", 0, 0, 1, 1);
517  p1->SetBottomMargin(.3);
518  p2->SetTopMargin(.7);
519  p1->SetFillStyle(0);
520  p2->SetFillStyle(0);
521  p1->Draw();
522  p2->Draw();
523 
524  p1->cd();
525  h1->SetStats(0);
526  h1->Draw(topOption[0]);
527  for(unsigned int ii=0;ii<topHistos.size();ii++){
528  topHistos[ii]->Draw("same "+topOption[ii]);}
529  h1->Draw("same "+topOption[0]);
530  p1->RedrawAxis(); //white histograms cover axes
531  h1->GetYaxis()->SetRangeUser(FindLimitY(topHistos,0) ,1.2*FindLimitY(topHistos,1));
532  h1->GetYaxis()->SetTitleSize(26);
533  h1->GetYaxis()->SetTitleFont(43);
534  h1->GetYaxis()->SetTitleOffset(1.2);
535  h1->GetYaxis()->SetLabelFont(43);
536  h1->GetYaxis()->SetLabelSize(18);
537  h1->GetXaxis()->SetLabelSize(0);
538  h1->GetXaxis()->SetTitleSize(0);
539  if (pidaxis){
540  PIDDivisions();
541  PIDBinLabels((std::string)pidtag);
542  }
543 
544  p2->cd();
545  h3->SetTitle("");
546  h3->GetYaxis()->SetRangeUser(max(0.,2-FindLimitY(bottomHistos,1)) ,1.01*FindLimitY(bottomHistos,1));
547  h3->GetYaxis()->SetTitleSize(26);
548  h3->GetYaxis()->SetTitleFont(43);
549  h3->GetYaxis()->SetTitleOffset(1.2);
550  h3->GetYaxis()->SetLabelFont(43);
551  h3->GetYaxis()->SetLabelSize(15);
552  h3->GetYaxis()->SetDecimals();
553 
554  h3->GetXaxis()->SetTitle(h1->GetXaxis()->GetTitle());
555  h3->GetXaxis()->SetTitleSize(26);
556  h3->GetXaxis()->SetTitleFont(43);
557  h3->GetXaxis()->SetLabelFont(43);
558  h3->GetXaxis()->SetLabelSize(18);
559  h3->GetXaxis()->SetLabelOffset(0.005);
560  h3->Draw(bottomOption[0]);
561  for(unsigned int ii=0;ii<bottomHistos.size();ii++){
562  bottomHistos[ii]->Draw("same "+bottomOption[ii]);}
563  h3->Draw("same "+bottomOption[0]);
564 
565  auto lone = new TLine();
566  lone->SetLineStyle(3);
567  p2->Update();
568  lone->DrawLine(p2->GetUxmin(),1.0,p2->GetUxmax(),1.0);
569  lone->DrawLine(p2->GetUxmin(),1.2,p2->GetUxmax(),1.2);
570  lone->DrawLine(p2->GetUxmin(),0.8,p2->GetUxmax(),0.8);
571  float p2Y = p2->GetUymax()-p2->GetUymin();
572  h3->GetYaxis()->SetRangeUser(p2->GetUymin()-(0.01*p2Y),p2->GetUymax()+(0.01*p2Y));
573 
574  if (pidaxis){
575  PIDDivisions();
576  PID2DAxis(h3);
577  }
578  return c;
579  }
580  ////////////////////////////////////////
581  void PlotDataMC(std::vector<TH1*>vnom, std::vector<bool>isdata,
582  TLegend * leg, TString pidtag="",TString htag="",
583  TString out_name="plot_FD",bool ratioplot=false,
584  bool ratioerrors=false, bool pidaxis = 0)
585  {
586  gStyle->SetPadLeftMargin(0.15);
587  gStyle->SetPadTopMargin(0.08);
588 
589  // Fill empty bins with 0 if in range and <0 if out of range
590  for ( unsigned int bin = 1; bin<=vnom[0]->GetNbinsX(); ++bin){
591  bool filledbin=false;
592  for ( unsigned int v = 0; v<vnom.size(); ++v){
593  if ( !isdata[v] )
594  if (vnom[v]->GetBinContent(bin) > 0 ) filledbin=true;
595  }
596  for ( unsigned int v = 0; v<vnom.size(); ++v){
597  if (vnom[v]->GetBinContent(bin) == 0 ){
598  if ( !filledbin ) vnom[v]->SetBinContent(bin,-0.01);
599  else vnom[v]->SetBinContent(bin,0.001);
600  }
601  vnom[v]->GetYaxis()->SetTitleOffset(1.05);
602  }
603  }
604 
605  // Draw marker w/error or line when appropriate
606  int datahist=-1;
607  for ( unsigned int v = 0; v<vnom.size(); ++v)
608  if (vnom[v]->GetLineColor() == kBlack ) datahist=v;
609 
610  // Option to make a ratio plot given
611  if(!ratioplot){
612  auto c2 = new TCanvas ("c2","c2",715,500);
613  for ( unsigned int v = 0; v<vnom.size(); ++v){
614  TGraphAsymmErrors* gr = GraphWithPoissonErrors(vnom[v], 0, 1);
615  if ( isdata[v] ){
616  if (v == datahist ) gr->Draw("ep same");
617  else vnom[v]->Draw("p same");
618  }
619  else vnom[v]->Draw("hist same");
620  }
621  gPad->RedrawAxis();
622 
623  // Adjust for ExPID axis
624  if (pidaxis){
625  FixLegend(leg,"ExPID");
626  NuePlotStyle style(pidtag,0);
627  style.PIDDivisions(true);
628  style.PID2DAxis(vnom[0],true);
629  style.PIDBinLabels(true);
630  }
631 
632  // pot label
633  TLatex* ltx = new TLatex(.175, 0.74,
634  TString::Format("%g#times10^{20} POT equiv.", 6.05).Data());
635  ltx->SetNDC();
636  ltx->SetTextAlign(11);
637  ltx->SetTextSize(0.0365);
638  ltx->Draw();
639  HTag(htag);
640  leg->Draw("same");
641 
642  // Save with & without prelim tag
643  gPad->SetFillStyle(0);
644  c2->SetFillStyle(0);
645  c2->SaveAs("plots/"+out_name+"_datamc.pdf");
646  c2->SaveAs("plots/"+out_name+"_datamc.eps");
647  c2->SaveAs("plots/"+out_name+"_datamc.C");
648  c2->SaveAs("plots/"+out_name+"_datamc.root");
649  Prelim();
650  c2->SaveAs("plots/"+out_name+"_pre_datamc.pdf");
651  c2->SaveAs("plots/"+out_name+"_pre_datamc.eps");
652  c2->SaveAs("plots/"+out_name+"_pre_datamc.C");
653  c2->SaveAs("plots/"+out_name+"_pre_datamc.root");
654  }
655  else{
656  std::vector<TH1*> vratio;
657  for(unsigned int ii=1;ii<vnom.size();ii++){
658  auto htemp = (TH1*) vnom[ii]->Clone(ana::UniqueName().c_str());
659  htemp->Divide(vnom[datahist],vnom[ii]);
660  htemp->SetMarkerColor(vnom[ii]->GetLineColor());
661  htemp->SetMarkerStyle(34);
662  htemp->SetMarkerSize(1.2);
663  htemp->GetYaxis()->SetTitle("Data / MC");
664  htemp->GetYaxis()->SetRangeUser(0.5,1.5);
665  vratio.push_back(htemp);
666  }
667  std::vector<TString> topOption = {" ","hist"};
668  std::vector<TString> bottomOption (vratio.size(),ratioerrors ?"p":"hist p");
669  auto cratio2= RatioPlot(vnom, topOption,vratio,bottomOption,pidtag,pidaxis);
670  cratio2->cd();
671  leg->Draw("same");
672  HTag(pidtag);
673  cratio2->SaveAs("plots/"+out_name+"_datamc_ratio.pdf");
674  cratio2->SaveAs("plots/png_"+out_name+"_datamc_ratio.png");
675  cratio2->SaveAs("plots/rootfiles/"+out_name+"_datamc_ratio.root");
676  }
677  }
678 
679  ////////////////////////////////////////
680  void PlotNDDataTotalMCComparison(TH1* hdata,std::vector<TH1*> htots,
681  TLegend * leg,TString pidtag="",
682  TString out_name="plot_nd",
683  bool ratioplot=false, bool pidaxis=false)
684  {
685  if(pidaxis) FixLegend(leg,"ExPID");
686  std::vector<TString> h_opt;
687  for(auto h:htots) h_opt.push_back("hist");
688  h_opt.push_back("p");
689 
690  if(!ratioplot){
691  auto c1 = new TCanvas ("c1","c1");
692  htots[0]->SetMaximum(1.3*htots[0]->GetMaximum());
693  for(auto h:htots) h->Draw("hist same");
694  htots[0]->Draw("hist same"); //need to redraw for error bands
695  hdata->Draw("p same");
696  PIDTag(pidtag);
697  if (pidaxis){
698  PIDDivisions();
699  PIDBinLabels((std::string)pidtag);
700  PID2DAxis(htots[0]);
701  }
702  leg->Draw("same");
703  c1->SaveAs("plots/"+out_name+"_datamc.pdf");
704  c1->SaveAs("plots/png_"+out_name+"_datamc.png");
705  c1->SaveAs("plots/rootfiles/"+out_name+"_datamc.root");
706  }
707  else{
708  std::vector<TH1*> ratios;
709  std::vector<TString> r_opt;
710  for(auto h:htots){
711  auto ratio = (TH1*) h->Clone();
712  ratio->Divide(hdata,h);
713  ratios.push_back(ratio);
714  r_opt.push_back("e1");
715  }
716  ratios[0]->GetYaxis()->SetTitle("Data / MC");
717  ratios[0]->GetYaxis()->SetRangeUser(0.5,1.5);
718  htots.push_back(hdata);
719  auto cratio1 =RatioPlot(htots,h_opt,ratios,r_opt,pidtag,pidaxis);
720  cratio1->cd();
721  leg->Draw("same");
722  PIDTag(pidtag);
723  cratio1->SaveAs("plots/"+out_name+"_datamc_ratio.pdf");
724  cratio1->SaveAs("plots/png_"+out_name+"_datamc_ratio.png");
725  cratio1->SaveAs("plots/rootfiles/"+out_name+"_datamc_ratio.root");
726  }
727  if(pidaxis){
728  htots.push_back(hdata);
729  auto cpads =ExPIDPlot(htots,h_opt);
730  cpads->cd();
731  FixLegend(leg,"ExPIDPads");
732  leg->Draw("same");
733  PIDTag(pidtag);
734  cpads->SaveAs("plots/"+out_name+"_datamc_pads.pdf");
735  cpads->SaveAs("plots/png_"+out_name+"_datamc_pads.png");
736  cpads->SaveAs("plots/rootfiles/"+out_name+"_datamc_pads.root");
737  }
738  }
739 
740  ////////////////////////////////////////
741  void PlotMCComponentsComparison(std::vector<TH1*>vnom,std::vector<TH1*>vshi,
742  TLegend * leg, TString pidtag="",
743  TString out_name="plot_nd",bool ratioplot=false,
744  bool ratioerrors=false, bool pidaxis = 0)
745  {
746 
747  if(pidaxis) FixLegend(leg,"ExPID");
748  if(!ratioplot){
749  auto c2 = new TCanvas ("c2","c2");
750  for(auto v:vnom) v->Draw("hist same");
751  for(auto v:vshi) v->Draw("hist same");
752  if (pidaxis){
753  PIDDivisions();
754  PIDBinLabels((std::string)pidtag);
755  PID2DAxis(vnom[0]);
756  }
757  PIDTag(pidtag);
758  leg->Draw("same");
759  c2->SaveAs("plots/"+out_name+"_mccomp.pdf");
760  c2->SaveAs("plots/png_"+out_name+"_mccomp.png");
761  c2->SaveAs("plots/rootfiles/"+out_name+"_mccomp.root");
762  c2->Close();
763  }
764  else{
765  std::vector<TH1*> vratio;
766  for(unsigned int ii=0;ii<vshi.size();ii++){
767  auto htemp = (TH1*) vnom[ii]->Clone(ana::UniqueName().c_str());
768  htemp->Divide(vshi[ii],vnom[ii]);
769  htemp->SetMarkerColor(vshi[ii]->GetLineColor());
770  htemp->SetMarkerStyle(34);
771  htemp->SetMarkerSize(1.2);
772  htemp->GetYaxis()->SetTitle("Shift / Nom");
773  htemp->GetYaxis()->SetRangeUser(0.5,1.5);
774  vratio.push_back(htemp);
775  }
776  vnom.insert(vnom.end(),vshi.begin(),vshi.end());
777  std::vector<TString> topOption (vnom.size(),"hist");
778  std::vector<TString> bottomOption (vratio.size(),ratioerrors ?"p":"hist p");
779  auto cratio2= RatioPlot(vnom, topOption,vratio,bottomOption,pidtag,pidaxis);
780  cratio2->cd();
781  leg->Draw("same");
782  PIDTag(pidtag);
783  cratio2->SaveAs("plots/"+out_name+"_mccomp_ratio.pdf");
784  cratio2->SaveAs("plots/png_"+out_name+"_mccomp_ratio.png");
785  cratio2->SaveAs("plots/rootfiles/"+out_name+"_mccomp_ratio.root");
786  }
787  if(pidaxis){
788  vnom.insert(vnom.end(),vshi.begin(),vshi.end());
789  std::vector<TString> topOption (vnom.size(),"hist");
790  auto cpads =ExPIDPlot(vnom,topOption);
791  cpads->cd();
792  FixLegend(leg,"ExPIDPads");
793  leg->Draw("same");
794  PIDTag(pidtag);
795  cpads->SaveAs("plots/"+out_name+"_mccomp_pads.pdf");
796  cpads->SaveAs("plots/png_"+out_name+"_mccomp_pads.png");
797  cpads->SaveAs("plots/rootfiles/"+out_name+"_mccomp_pads.root");
798 
799  }
800  }
801 
802  void CompareNDDataMC (TDirectory * d_no, TDirectory* d_sh,TString plottitle,TString out_name, TString pidtag, bool printtable=true, bool pidaxis = 0){
803  auto data = Spectrum::LoadFrom(d_no->GetDirectory("data_comp"));
804  auto kNDPOT = data->POT();
805  auto hdata = data->ToTH1(kNDPOT);
806  hdata ->SetMarkerStyle(kFullCircle);
807 
808  auto hnue_no = (*Spectrum::LoadFrom(d_no->GetDirectory("nue_comp")) +
809  *Spectrum::LoadFrom(d_no->GetDirectory("antinue_comp"))).ToTH1(kNDPOT, kBeamNueBackgroundColor, 1);
810  auto hnumu_no = (*Spectrum::LoadFrom(d_no->GetDirectory("numu_comp")) +
811  *Spectrum::LoadFrom(d_no->GetDirectory("antinumu_comp"))).ToTH1(kNDPOT, kNumuBackgroundColor, 1);
812  auto hnc_no = (*Spectrum::LoadFrom(d_no->GetDirectory("nc_tot_comp"))).ToTH1(kNDPOT, kNCBackgroundColor, 1);
813  auto htot_no = (TH1*)hnue_no->Clone();
814  htot_no->Add(hnumu_no);htot_no->Add(hnc_no);
815 
816  auto hnue_sh = (*Spectrum::LoadFrom(d_sh->GetDirectory("nue_comp")) +
817  *Spectrum::LoadFrom(d_sh->GetDirectory("antinue_comp"))).ToTH1(kNDPOT, kBeamNueBackgroundColor, 2);
818  auto hnumu_sh = (*Spectrum::LoadFrom(d_sh->GetDirectory("numu_comp")) +
819  *Spectrum::LoadFrom(d_sh->GetDirectory("antinumu_comp"))).ToTH1(kNDPOT, kNumuBackgroundColor, 2);
820  auto hnc_sh = (*Spectrum::LoadFrom(d_sh->GetDirectory("nc_tot_comp"))).ToTH1(kNDPOT, kNCBackgroundColor, 2);
821  auto htot_sh = (TH1*)hnue_sh->Clone();
822  htot_sh->Add(hnumu_sh);htot_sh->Add(hnc_sh);
823 
824  htot_no->SetLineColor(kTotalMCColor+1);
825  htot_sh->SetLineColor(kTotalMCColor+2);
826 
827  htot_no->SetMaximum(FindLimitY({htot_no,hdata}));
828  htot_no->SetTitle(plottitle);
829  hnumu_no->SetMaximum(FindLimitY({hnumu_no,hnc_no,hnue_no,hnumu_sh,hnc_sh,hnue_sh}));
830  hnumu_no->SetTitle(plottitle);
831 
832  CenterTitles(htot_no);
833  CenterTitles(hnue_no);
834  TGaxis::SetMaxDigits(3);
835  htot_no->GetXaxis()->SetDecimals();
836 
837  auto leg1=new TLegend(0.63,0.67,0.85,0.87);
838  leg1->SetFillStyle(0);
839  leg1->AddEntry(hdata, "ND Data","epl");
840  leg1->AddEntry(htot_no,"Nominal","l");
841  leg1->AddEntry(htot_sh,"Shift","l");
842 
843  PlotNDDataTotalMCComparison(hdata,{htot_no,htot_sh},leg1,pidtag,out_name,printratioplots,pidaxis);
844 
845  htot_sh->SetLineColor(kTotalMCColor+1);
846  std::vector<TH1*>vnom={hnumu_no,hnc_no,hnue_no};
847  std::vector<TH1*>vshi={hnumu_sh,hnc_sh,hnue_sh};
848 
849  auto leg2 = DefaultNueLegend(false,false,false);
850  TH1* htemp = new TH1F("", "", 1, 0, 1);
851  htemp->SetLineColor(kGray+1);
852  htemp->SetLineStyle(2);
853  leg2->AddEntry(htemp,"Shifted","l");
854 
855  PlotMCComponentsComparison(vnom,vshi,leg2,pidtag,out_name,printratioplots,false,pidaxis);
856 
857  if(printtable){
858  std::cout << std::string(60,'-') << std::endl;
859  std:: cout << "% - " << (plottitle.Remove(plottitle.First(";"))).Data() << " - ND - " << pidtag << " - " << std::endl;
860  ComparisonTable ({htot_no,hnue_no,hnc_no,hnumu_no},
861  {htot_sh,hnue_sh,hnc_sh,hnumu_sh},
862  {"Total MC","Beam $\\nu_{e}$","NC", "$\\nu_{\\mu}$ CC"}
863  );
864  std::cout << std::string(60,'-') << std::endl;
865  }
866 
867 
868  }
869 
870 
871  void ComparePredictions(IPrediction * pred_no, IPrediction * pred_sh, TString plottitle, TString out_name, TString pidtag, bool printtable=true, bool pidaxis = false){
872 
873  TH1* hNueSig_no = pred_no->PredictComponent(&calc, ana::Flavors::kNuMuToNuE,
875  TH1* hNueBeam_no = pred_no->PredictComponent(&calc, ana::Flavors::kNuEToNuE,
877  TH1* hNumu_no = pred_no->PredictComponent(&calc, ana::Flavors::kAllNuMu,
879  TH1* hNC_no = pred_no->PredictComponent(&calc, ana::Flavors::kAll,
881  TH1* hTau_no = pred_no->PredictComponent(&calc, ana::Flavors::kAllNuTau,
883 
884  TH1* hNueSig_sh = pred_sh->PredictComponent(&calc, ana::Flavors::kNuMuToNuE,
886  TH1* hNueBeam_sh = pred_sh->PredictComponent(&calc, ana::Flavors::kNuEToNuE,
888  TH1* hNumu_sh = pred_sh->PredictComponent(&calc, ana::Flavors::kAllNuMu,
890  TH1* hNC_sh = pred_sh->PredictComponent(&calc, ana::Flavors::kAll,
892  TH1* hTau_sh = pred_sh->PredictComponent(&calc, ana::Flavors::kAllNuTau,
894 
895  std::cout<<"Making plots"<<std::endl;
896  hNueBeam_no->SetLineColor(kBeamNueBackgroundColor);
897  hNumu_no->SetLineColor(kNumuBackgroundColor);
898  hNC_no->SetLineColor(kNCBackgroundColor);
899  hNueSig_no->SetLineColor(kNueSignalColor);
900 
901  hNueBeam_sh->SetLineColor(kBeamNueBackgroundColor);
902  hNumu_sh->SetLineColor(kNumuBackgroundColor);
903  hNC_sh->SetLineColor(kNCBackgroundColor);
904  hNueSig_sh->SetLineColor(kNueSignalColor);
905 
906  hNueBeam_sh->SetLineStyle(2);
907  hNumu_sh->SetLineStyle(2);
908  hNC_sh->SetLineStyle(2);
909  hNueSig_sh->SetLineStyle(2);
910 
911  hNueSig_no->SetMaximum(1.3*max(hNueSig_no->GetMaximum(),hNueSig_sh->GetMaximum()));
912  hNueSig_no->SetMaximum(FindLimitY({hNueSig_no,hNueBeam_no,hNC_no,hNumu_no,
913  hNueSig_sh,hNueBeam_sh,hNC_sh,hNumu_sh}));
914  hNueSig_no->SetTitle(plottitle);
915  CenterTitles(hNueSig_no);
916  TGaxis::SetMaxDigits(3);
917 
918  std::vector<TH1*>vnom={hNueSig_no,hNueBeam_no,hNC_no,hNumu_no};
919  std::vector<TH1*>vshi={hNueSig_sh,hNueBeam_sh,hNC_sh,hNumu_sh};
920  auto leg = DefaultNueLegend(true,false,false);
921  TH1* htemp = new TH1F("", "", 1, 0, 1);
922  htemp->SetLineColor(kGray+1);
923  htemp->SetLineStyle(2);
924  leg->AddEntry(htemp,"Shifted","l");
925 
926  PlotMCComponentsComparison(vnom,vshi,leg,pidtag,out_name,printratioplots,false,pidaxis);
927 
928  if(printtable){
929  auto hTotBkg_no = (TH1*) hNueBeam_no->Clone();
930  hTotBkg_no->Add(hNC_no); hTotBkg_no->Add(hNumu_no); hTotBkg_no->Add(hTau_no);
931 
932  auto hTotBkg_sh = (TH1*) hNueBeam_sh->Clone();
933  hTotBkg_sh->Add(hNC_sh); hTotBkg_sh->Add(hNumu_sh); hTotBkg_sh->Add(hTau_sh);
934 
935  std::cout << std::string(60,'-') << std::endl;
936  std::cout << "% - " << (plottitle.Remove(plottitle.First(";"))).Data() << " - " << pidtag << " - " << std::endl;
937  ComparisonTable ({hNueSig_no, hTotBkg_no, hNueBeam_no,hNC_no,hNumu_no,hTau_no},
938  {hNueSig_sh, hTotBkg_sh, hNueBeam_sh,hNC_sh,hNumu_sh,hTau_sh},
939  {"$\\nu_e$ signal","Total Beam Bkg","Beam $\\nu_{e}$CC","NC", "$\\nu_{\\mu}$CC", "$\\nu_{\\tau}$CC"}
940  );
941  std::cout << std::string(60,'-') << std::endl;
942  }
943 
944  }
945 
946  void FormatErrorBand(TH1* hplus, TH1* hminus, bool signal=false, bool fixbins=false){
947 
948  int mccolor = signal ? kViolet -9 : kTotalMCErrorBandColor;
949  hplus ->SetLineColor(mccolor);
950  hminus->SetLineColor(mccolor);
951  hplus ->SetFillColor(mccolor);
952  hminus->SetFillColor(10);
953  hplus ->SetMarkerColor(kRed-7); hminus->SetMarkerColor(kBlue-7);
954 
955  if (fixbins){
956  for(int i=1;i<=hplus->GetNbinsX();i++){
957  double tplus = hplus->GetBinContent(i);
958  double tminus = hminus->GetBinContent(i);
959  hplus->SetBinContent(i,max(tplus, tminus));
960  hminus->SetBinContent(i,min(tplus,tminus));
961  }
962  }
963 
964  }
965 
966  void PlotMCComponentsErrorBand(std::vector<TH1*>vnom,std::vector<TH1*>vshi,
967  TLegend * leg, TString pidtag="",
968  TString out_name="plot_nd",bool ratioplot=false,
969  bool ratioerrors=false, bool pidaxis = 0)
970  {
971  auto hplu_plot = (TH1*) vshi[0]->Clone();
972  auto hmin_plot = (TH1*) vshi[1]->Clone();
973  FormatErrorBand(hplu_plot,hmin_plot,out_name.Contains("_sig"),true);
974  std::vector<TH1*>vshi_plot={hplu_plot,hmin_plot};
975 
976  if(pidaxis) FixLegend(leg,"ExPID");
977  if(!ratioplot){
978  auto c2 = new TCanvas ("c2","c2");
979  vnom[0]->SetMaximum(1.2*vnom[0]->GetMaximum());
980  vnom[0]->Draw("hist");
981  for(auto v:vshi_plot) v->Draw("hist same");
982  for(auto v:vnom) v->Draw("hist same");
983  if (pidaxis){
984  PIDDivisions();
985  PIDBinLabels((std::string)pidtag);
986  PID2DAxis(vnom[0]);
987  }
988  PIDTag(pidtag);
989  leg->Draw("same");
990  c2->SaveAs("plots/"+out_name+"_mccomp.pdf");
991  c2->SaveAs("plots/png_"+out_name+"_mccomp.png");
992  c2->SaveAs("plots/rootfiles/"+out_name+"_mccomp.root");
993  c2->Close();
994  }
995  else{
996  std::vector<TH1*> vratio;
997  for(unsigned int ii=0;ii<vshi.size();ii++){
998  auto htemp = (TH1*) vshi[ii]->Clone(ana::UniqueName().c_str());
999  htemp->Divide(vshi[ii],vnom[0]);
1000  htemp->SetMarkerStyle(34);
1001  htemp->SetMarkerSize(1.2);
1002  htemp->GetYaxis()->SetTitle("Shift / Nom");
1003  vratio.push_back(htemp);
1004  }
1005  vratio[0]->GetYaxis()->SetRangeUser(0.5,1.5);
1006  vnom.insert(vnom.begin()+1,vshi_plot.begin(),vshi_plot.end());
1007  std::vector<TString> topOption (vnom.size(),"hist");
1008  std::vector<TString> bottomOption (vratio.size(),ratioerrors ?"p":"hist p");
1009  auto cratio2= RatioPlot(vnom, topOption,vratio,bottomOption,pidtag,pidaxis);
1010  cratio2->cd();
1011  leg->Draw("same");
1012  PIDTag(pidtag);
1013  cratio2->SaveAs("plots/"+out_name+"_mccomp_ratio.pdf");
1014  cratio2->SaveAs("plots/png_"+out_name+"_mccomp_ratio.png");
1015  cratio2->SaveAs("plots/rootfiles/"+out_name+"_mccomp_ratio.root");
1016  }
1017  }
1018 
1019 ////////////////////////////////////////
1021  TH1* sig, *bkg,*beam,*numu,*nc,*tau;
1022  TH1* nuetonumu,*numutonumu,*nuetonutau,*numutonutau;
1023  };
1024 
1025  struct NueMCComponents GetNDMCComponents(TDirectory * d_no, double kNDPOT,int linestyle =1){
1027  ret.beam = (*Spectrum::LoadFrom(d_no->GetDirectory("nue_comp")) +
1028  *Spectrum::LoadFrom(d_no->GetDirectory("antinue_comp"))).ToTH1(kNDPOT, kBeamNueBackgroundColor, linestyle);
1029  ret.numu = (*Spectrum::LoadFrom(d_no->GetDirectory("numu_comp")) +
1030  *Spectrum::LoadFrom(d_no->GetDirectory("antinumu_comp"))).ToTH1(kNDPOT, kNumuBackgroundColor, linestyle);
1031  ret.nc = (*Spectrum::LoadFrom(d_no->GetDirectory("nc_tot_comp"))).ToTH1(kNDPOT, kNCBackgroundColor, linestyle);
1032  ret.bkg = (TH1*)ret.beam->Clone();
1033  ret.bkg->Add(ret.numu);ret.bkg->Add(ret.nc);
1034  ret.sig=0; ret.tau=0;
1035  ret.nuetonumu=0; ret.nuetonutau=0;
1036  ret.numutonumu=0; ret.numutonutau=0;
1037  ret.bkg->SetLineColor(kTotalMCColor);
1038  return ret;
1039 
1040  }
1041 ////////////////////////////////////////
1042  void CompareNDDataMC (TDirectory * d_no, TDirectory* d_pl,TDirectory*d_mi,TString plottitle,TString out_name, TString pidtag, bool printtable=true, bool pidaxis = 0){
1043  auto data = Spectrum::LoadFrom(d_no->GetDirectory("data_comp"));
1044  auto kNDPOT = data->POT();
1045  auto hdata = data->ToTH1(kNDPOT);
1046  hdata ->SetMarkerStyle(kFullCircle);
1047 
1048  auto hnom = GetNDMCComponents(d_no,kNDPOT);
1049  auto hplu = GetNDMCComponents(d_pl,kNDPOT);
1050  auto hmin = GetNDMCComponents(d_mi,kNDPOT);
1051 
1052  FormatErrorBand(hplu.bkg,hmin.bkg);
1053 
1054  hnom.bkg->SetMaximum(FindLimitY({hnom.bkg,hdata,hplu.bkg,hmin.bkg}));
1055  hnom.bkg->SetTitle(plottitle);
1056 
1057  CenterTitles(hnom.bkg);
1058  TGaxis::SetMaxDigits(3);
1059  hnom.bkg->GetXaxis()->SetDecimals();
1060 
1061  auto leg1=new TLegend(kLegendDefaultX1,kLegendDefaultY1,kLegendDefaultX2,kLegendDefaultY2);
1062  leg1->SetFillStyle(0);
1063  leg1->AddEntry(hdata, "ND Data","epl");
1064  leg1->AddEntry(hnom.bkg,"Nominal","l");
1065  leg1->AddEntry(hplu.bkg,"Syst.","f");
1066 
1067  PlotNDDataTotalMCComparison(hdata,{hnom.bkg,hplu.bkg,hmin.bkg},
1068  leg1,pidtag,out_name,printratioplots,pidaxis);
1069 
1070 
1071  auto leg2 = DefaultNueLegend(false,true,false);
1072  leg2->AddEntry(hplu.bkg,"Syst.","f");
1073  PlotMCComponentsErrorBand({hnom.bkg,hnom.numu,hnom.nc,hnom.beam},{hplu.bkg,hmin.bkg},
1074  leg2,pidtag,out_name,printratioplots,false,pidaxis);
1075 
1076 
1077  if(printtable){
1078  std::cout << "% " << std::string(60,'-') << std::endl;
1079  std::cout << "% - " << (plottitle.Remove(plottitle.First(";"))).Data() << " - ND - " << pidtag << " - " << std::endl;
1080  ComparisonTable ({hnom.bkg,hnom.beam,hnom.nc,hnom.numu},
1081  {hplu.bkg,hplu.beam,hplu.nc,hplu.numu},
1082  {hmin.bkg,hmin.beam,hmin.nc,hmin.numu},
1083  {"Total MC","Beam $\\nu_{e}$","NC", "$\\nu_{\\mu}$ CC"},
1084  out_name
1085  );
1086  std::cout << "% " << std::string(60,'-') << std::endl;
1087  }
1088 
1089 
1090  }
1091 
1092  struct NueMCComponents GetFDMCComponents(IPrediction * pred_no,int linestyle=1){
1094  ret.sig = pred_no->PredictComponent(&calc, ana::Flavors::kNuMuToNuE,
1096  ret.beam = pred_no->PredictComponent(&calc, ana::Flavors::kNuEToNuE,
1098  ret.numu = pred_no->PredictComponent(&calc, ana::Flavors::kAllNuMu,
1100  ret.nc = pred_no->PredictComponent(&calc, ana::Flavors::kAll,
1102  ret.tau = pred_no->PredictComponent(&calc, ana::Flavors::kAllNuTau,
1104  ret.bkg = (TH1*) ret.beam->Clone();
1105  ret.bkg->Add(ret.nc); ret.bkg->Add(ret.numu); ret.bkg->Add(ret.tau);
1106  ret.bkg->SetLineColor(kTotalMCColor);
1107 
1108  ret.nuetonumu = pred_no->PredictComponent(&calc, ana::Flavors::kNuEToNuMu,
1110  ret.numutonumu = pred_no->PredictComponent(&calc, ana::Flavors::kNuMuToNuMu,
1112  ret.nuetonutau = pred_no->PredictComponent(&calc, ana::Flavors::kNuEToNuTau,
1114  ret.numutonutau = pred_no->PredictComponent(&calc, ana::Flavors::kNuMuToNuTau,
1116  return ret;
1117  }
1118 
1119  Double_t GetPercentage(IPrediction* pred_no, IPrediction* pred_pl, IPrediction* pred_mi, int compIdx){
1120 
1121  auto hnom = GetFDMCComponents(pred_no);
1122  auto hplu = GetFDMCComponents(pred_pl);
1123  auto hmin = GetFDMCComponents(pred_mi);
1124 
1125  Double_t percentShift[kNumComp];
1126  percentShift[signalevts]=100 * ( TMath::Abs( ( hnom.sig->Integral() - hplu.sig->Integral()) / hnom.sig->Integral() ) + TMath::Abs( (hnom.sig->Integral() - hmin.sig->Integral()) / hnom.sig->Integral() ) ) / 2.0;
1127  percentShift[background]=100 * ( TMath::Abs( ( hnom.bkg->Integral() - hplu.bkg->Integral()) / hnom.bkg->Integral() ) + TMath::Abs( (hnom.bkg->Integral() - hmin.bkg->Integral()) / hnom.bkg->Integral() ) ) / 2.0;
1128  percentShift[numucc]=100 * ( TMath::Abs( ( hnom.numu->Integral() - hplu.numu->Integral()) / hnom.numu->Integral() ) + TMath::Abs( (hnom.numu->Integral() - hmin.numu->Integral()) / hnom.numu->Integral() ) ) / 2.0;
1129  percentShift[beamnue]=100 * ( TMath::Abs( ( hnom.beam->Integral() - hplu.beam->Integral()) / hnom.beam->Integral() ) + TMath::Abs( (hnom.beam->Integral() - hmin.beam->Integral()) / hnom.beam->Integral() ) ) / 2.0;
1130  percentShift[neutralc]=100 * ( TMath::Abs( ( hnom.nc->Integral() - hplu.nc->Integral()) / hnom.nc->Integral() ) + TMath::Abs( (hnom.nc->Integral() - hmin.nc->Integral()) / hnom.nc->Integral() ) ) / 2.0;
1131  percentShift[nutau]=100 * ( TMath::Abs( ( hnom.tau->Integral() - hplu.tau->Integral()) / hnom.tau->Integral() ) + TMath::Abs( (hnom.tau->Integral() - hmin.tau->Integral()) / hnom.tau->Integral() ) ) / 2;
1132 
1133  return percentShift[compIdx];
1134  }
1135 
1136  TString Latexify (TString s){
1137  for (int i=0; i<s.Length(); i++){
1138  if (s[i]=='_'){
1139  s.Insert(i,"\\");
1140  i++;
1141  }
1142  }
1143  return s;
1144  }
1145 
1146  void PrintFullShiftTable(int nSyst, Double_t percentages[][kNumSels][kNumVars][kNumComp],int selIdx, int varIdx, TString pidtag, TString out_name, TString* shiftLabel, Bool_t isGenie, Bool_t* turnedOn){
1147  TString compLabel[kNumComp]={"Signal", "Tot. Background", "$\\nu_\\mu$ CC", "Beam $\\nu_e$", "N.C.", "$\\nu_\\tau$"};
1148 
1149  std::ofstream theFile;
1150  theFile.open(out_name.Data(),ios::out|ios::app);
1151  theFile << "\\newcommand{\\FullShiftTable" << pidtag << "}"<<std::endl;
1152 
1153  theFile << "{ \\begin{tabular}{ c | c | ";
1154  for (int i=1; i<kNumComp; ++i)
1155  theFile<<"c ";
1156  theFile<<"}\n";
1157 
1158  for (int i=0; i<kNumComp; ++i){
1159  theFile<<" & "<<compLabel[i];
1160  }
1161 
1162 
1163 
1164  theFile<<"\\\\\\hline\\hline\n";
1165  for (int iSyst=0; iSyst<nSyst; ++iSyst){
1166  if (!isGenie || turnedOn[iSyst]){
1167  theFile<<Latexify(shiftLabel[iSyst]);
1168  for (int i=0; i<kNumComp; ++i){
1169  theFile<<" & "<< std::fixed <<std::setprecision(2) <<percentages[iSyst][selIdx][varIdx][i] ;
1170  }
1171  theFile<<" \\\\\n";
1172  }
1173  }
1174  theFile << "\\hline \n";
1175  theFile << "\\end{tabular}} \n\n";
1176  theFile.close();
1177  }
1178 
1179 
1180 
1181  void ComparePredictions(IPrediction * pred_no, IPrediction * pred_pl, IPrediction * pred_mi,TString plottitle, TString out_name, TString pidtag, bool printtable=true, bool pidaxis = false){
1182 
1183  auto hnom = GetFDMCComponents(pred_no);
1184  auto hplu = GetFDMCComponents(pred_pl);
1185  auto hmin = GetFDMCComponents(pred_mi);
1186 
1187  TGaxis::SetMaxDigits(3);
1188 
1189  //FD Bkg Error band plot
1190  FormatErrorBand(hplu.bkg,hmin.bkg,false);
1191  hnom.bkg->SetMaximum(FindLimitY({hnom.bkg,hplu.bkg,hmin.bkg}));
1192  hnom.bkg->SetTitle(plottitle);
1193  CenterTitles(hnom.bkg);
1194  auto leg = DefaultNueLegend(false,true,false);
1195  leg->AddEntry(hplu.bkg,"Syst.","f");
1196 
1197  PlotMCComponentsErrorBand({hnom.bkg,hnom.beam,hnom.nc,hnom.numu},{hplu.bkg,hmin.bkg},
1198  leg,pidtag,out_name+"_bkg",printratioplots,false,pidaxis);
1199 
1200  //FD Sig Error band plot
1201  FormatErrorBand(hplu.sig,hmin.sig,true);
1202  hnom.sig->SetMaximum(FindLimitY({hnom.sig,hplu.sig,hmin.sig}));
1203  hnom.sig->SetTitle(plottitle);
1204  CenterTitles(hnom.sig);
1205  auto leg2 = new TLegend(kLegendDefaultX1,kLegendDefaultY1,kLegendDefaultX2,kLegendDefaultY2);
1206  leg2->AddEntry(hnom.sig,"#nu_{e} Signal","l");
1207  leg2->AddEntry(hplu.sig,"Syst.","f");
1208 
1209  PlotMCComponentsErrorBand({hnom.sig},{hplu.sig,hmin.sig},leg2,pidtag,out_name+"_sig",printratioplots,false,pidaxis);
1210 
1211  //Full comparison table
1212  if(printtable){
1213  std::cout << "% "<< std::string(60,'-') << std::endl;
1214  std::cout << "% - " << (plottitle.Remove(plottitle.First(";"))).Data() << " - " << pidtag << " - " << std::endl;
1215  ComparisonTable ({hnom.sig, hnom.bkg, hnom.beam, hnom.nc, hnom.numu, hnom.tau},
1216  {hplu.sig, hplu.bkg, hplu.beam, hplu.nc, hplu.numu, hplu.tau},
1217  {hmin.sig, hmin.bkg, hmin.beam, hmin.nc, hmin.numu, hmin.tau},
1218  {"$\\nu_e$ signal","Total Beam Bkg","Beam $\\nu_{e}$CC","NC", "$\\nu_{\\mu}$CC", "$\\nu_{\\tau}$CC"},
1219  out_name
1220  );
1221  std::cout << "% "<< std::string(60,'-') << std::endl;
1222  }
1223 
1224  }
1225 
1226  void MakeNueSystematicsFile(IPrediction * pred_no, IPrediction * pred_sh, TString pid, TString sigma, TFile* nue_syst_file){
1227  auto hnom = GetFDMCComponents(pred_no);
1228  auto hshi = GetFDMCComponents(pred_sh);
1229 
1230  std::vector<TString> channels = {"numutonumucc" ,"numutonuecc",
1231  "numutonutaucc","nuetonumucc",
1232  "nuetonuecc" ,"nuetonutaucc" ,
1233  "nc"};
1234  std::vector<TH1*> vshi = {hshi.numutonumu, hshi.sig,
1235  hshi.numutonutau, hshi.nuetonumu,
1236  hshi.beam, hshi.nuetonutau,
1237  hshi.nc};
1238  std::vector<TH1*> vnom = {hnom.numutonumu, hnom.sig,
1239  hnom.numutonutau, hnom.nuetonumu,
1240  hnom.beam, hnom.nuetonutau,
1241  hnom.nc};
1242  if(hshi.sig->GetNbinsX()!=30) {
1243  std::cerr << "\nWARNING\n"
1244  <<"MakeNueSystematicsFile doesn't know what to do with this binning\n"
1245  <<"Nothing will be saved\n";
1246  return;
1247  };
1248 
1249  const ana::Binning* bins = 0;
1250  if(pid == "LID") bins = &ana::kLIDNLBinning;
1251  if(pid == "LEM") bins = &ana::kLEMNLBinning;
1252  if(pid == "CVN") bins = &ana::kCVNNLBinning;
1253  assert(bins);
1254 
1255  Double_t pidEdges[]={bins->Edges()[0],bins->Edges()[1],bins->Edges()[2],bins->Edges()[3]};
1256  for(unsigned int i=0;i<vshi.size();i++){
1257  vshi[i]->Divide(vnom[i]);
1258  TString hname = "fdratio_"+pid+"_"+channels[i]+"_"+sigma;
1259  TString htitle = "FD "+channels[i]+";Reconstructed Energy (GeV);"+pid;
1260 
1261  TH2D* h_epid = new TH2D(hname,htitle,3,pidEdges,10,0,5);
1262  for(int xx=1;xx<=3;xx++){
1263  for(int yy=1;yy<=10;yy++){
1264  double ratio = vshi[i]->GetBinContent(yy+10*(xx-1));
1265  // double shift = ratio>0? ratio -1:0;
1266  if(ratio>0)h_epid->SetBinContent(xx,yy,ratio-1);
1267  }
1268  }
1269  nue_syst_file->cd();
1270  h_epid->Write();
1271  }
1272 
1273  }
1274 
1275  void PrintRawEventCounts(TDirectory * dpred, TString title){
1276 
1277  std::vector<int> ndcounts;
1278  std::vector<int> fdcounts;
1279  std::vector<TString> labels;
1280  //Data
1281  labels.push_back("Data - $\\nu_\\mu $");
1282  ndcounts.push_back(((TH1*)dpred->Get("extrap/MEextrap/Decomp/data_comp/hist"))->GetEntries());
1283  fdcounts.push_back(0);
1284  labels.push_back("Data - $\\nu_e $");
1285  ndcounts.push_back(((TH1*)dpred->Get("extrap/EEextrap/Decomp/data_comp/hist"))->GetEntries());
1286  fdcounts.push_back(0);
1287 
1288  //Signal
1289  labels.push_back("Signal - $\\nu_\\mu \\to \\nu_e $");
1290  labels.push_back("Signal - $\\bar{\\nu}_\\mu \\to \\bar{\\nu}_e $");
1291  for(TString xp:{"ME","MEAnti"}){
1292  ndcounts.push_back(((TH1*)dpred->Get("extrap/"+xp+"extrap/RecoToTrueND/hist"))->GetEntries());
1293  fdcounts.push_back(((TH1*)dpred->Get("extrap/"+xp+"extrap/TrueToRecoFD/hist"))->GetEntries());
1294  }
1295 
1296  //Major bkg
1297  labels.push_back("Bkg. - $\\nu_e \\to \\nu_e $");
1298  labels.push_back("Bkg. - $\\nu_\\mu \\to \\nu_\\mu $");
1299  labels.push_back("Bkg. - $NC \\to NC $");
1300  for(TString xp:{"EE","MM","NC"}){
1301  ndcounts.push_back(((TH1*)dpred->Get("extrap/"+xp+"extrap/RecoND/hist"))->GetEntries());
1302  fdcounts.push_back(((TH1*)dpred->Get("extrap/"+xp+"extrap/TrueToRecoFD/hist"))->GetEntries());
1303  }
1304 
1305  //Minor bkg
1306  labels.push_back("Bkg. - $\\nu_e \\to \\nu_\\tau $");
1307  labels.push_back("Bkg. - $\\nu_\\mu \\to \\nu_\\tau $");
1308  for(TString xp: {"ET","MT"}){
1309  ndcounts.push_back(0);
1310  fdcounts.push_back(((TH1*)dpred->Get("extrap/"+xp+"extrap/RecoFD/hist"))->GetEntries());
1311  }
1312  labels.push_back("Bkg. - Other");
1313  double fdcount_minor=0;
1314  for(TString xp: {"EM","EEAnti","EMAnti","ETAnti","MMAnti","MTAnti"}){
1315  fdcount_minor+=(((TH1*)dpred->Get("extrap/"+xp+"extrap/RecoFD/hist"))->GetEntries());
1316  }
1317  ndcounts.push_back(0);
1318  fdcounts.push_back(fdcount_minor);
1319 
1320  //Print
1321  std::cout << "% " << std::string(60,'-') << std::endl;
1322  std::cout << title << "\n";
1323  std::cout << "\\begin{tabular}{ l | r | r }\n"
1324  << std::setw(45) << " " << " & "
1325  << std::setw(10) << "ND" << " & "
1326  << std::setw(10) << "FD"
1327  << " \\\\ \\hline \n";
1328  for(unsigned int i=0;i<labels.size();i++){
1329  std::cout << std::setw(45) << labels[i] << " & "
1330  << std::setw(10)<< ndcounts[i] << " & "
1331  << std::setw(10)<< fdcounts[i] << " \\\\ ";
1332  if(i==1 || i==3 || i==6) std::cout << "\\hline \n";
1333  else std::cout << "\n";
1334  }
1335  std::cout << "\\end{tabular} \n";
1336  std::cout << "% " << std::string(60,'-') << std::endl;
1337  }
1338 
1339 
1340 
1341 }
1342 
1343 
1344 
T max(const caf::Proxy< T > &a, T b)
struct NueMCComponents GetNDMCComponents(TDirectory *d_no, double kNDPOT, int linestyle=1)
const int pidBins
const int kNumVars
Definition: vars.h:14
Represent the binning of a Spectrum&#39;s x-axis.
Definition: Binning.h:16
Double_t xx
Definition: macro.C:12
Oscillation analysis framework, runs over CAF files outside of ART.
TString Latexify(TString s)
double FindLimitY(std::vector< TH1 * > histos, bool doMax=true)
set< int >::iterator it
const double kLegendDefaultX2
Float_t y1[n_points_granero]
Definition: compare.C:5
(&#39; appearance&#39;)
Definition: IPrediction.h:18
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:553
virtual Spectrum PredictComponent(osc::IOscCalculator *calc, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const =0
TCanvas * ExPIDPlot(std::vector< TH1 * > topHistos, std::vector< TString > topOption)
const Binning kCVNNLBinning
Nonlinear binning for CVN.
Definition: Binning.h:121
Float_t x1[n_points_granero]
Definition: compare.C:5
TH1F * h3
Definition: berger.C:36
TLegend * DefaultNueLegend(bool sig=true, bool tot=true, bool data=false, double x1=kLegendDefaultX1, double y1=kLegendDefaultY1, double x2=kLegendDefaultX2, double y2=kLegendDefaultY2)
void FixLegend(TLegend *leg, TString opt="default")
(&#39;beam &#39;)
Definition: IPrediction.h:15
const char * p
Definition: xmltok.h:285
void PID2DShortAxis(TH1 *hist)
void PIDBinLabels(bool shortaxis)
void PlotNDDataTotalMCComparison(TH1 *hdata, std::vector< TH1 * > htots, TLegend *leg, TString pidtag="", TString out_name="plot_nd", bool ratioplot=false, AxisType pidaxis=kNue1bin)
OStream cerr
Definition: OStream.cxx:7
void CompareNDDataMC(TDirectory *d_no, TDirectory *d_pl, TDirectory *d_mi, TString plottitle, TString out_name, TString pidtag, bool printtable=true, bool pidaxis=0)
TLegend * leg1
Definition: plot_hist.C:105
TH1F * hplus
Definition: Xsec_final.C:75
TGraphAsymmErrors * GraphWithPoissonErrors(const TH1 *h, bool noErrorsXaxis, bool drawEmptyBins)
Calculate statistical errors appropriate for small Poisson numbers.
Definition: Plots.cxx:842
const Color_t kTotalMCErrorBandColor
Definition: Style.h:18
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)
void HTag(TString pid)
TH1F * hsig
Definition: Xdiff_gwt.C:59
cout<< t1-> GetEntries()
Definition: plottest35.C:29
const Color_t kNumuBackgroundColor
Definition: Style.h:31
const Binning bins
Definition: sa_fd_arrays.h:120
TH1 * RebinToShortAxis(TH1 *h)
const Color_t kNueSignalColor
Definition: Style.h:19
static const double kFDPOT
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)
TCanvas * RatioPlot(std::vector< TH1 * > topHistos, std::vector< TString > topOption, std::vector< TH1 * > bottomHistos, std::vector< TString > bottomOption, TString pidtag, AxisType pidaxis=kNue1bin)
const XML_Char const XML_Char * data
Definition: expat.h:268
void MakeNueSystematicsFile(IPrediction *pred_no, IPrediction *pred_sh, TString pid, TString sigma, TFile *nue_syst_file)
c2
Definition: demo5.py:33
const XML_Char * s
Definition: expat.h:262
struct DataMCComponets GetFDMCComponents(osc::IOscCalculator *calc, IPrediction *pred_no, TString output_name="nue", int linestyle=1, bool bkgdetails=false)
Charged-current interactions.
Definition: IPrediction.h:39
void CenterTitles(TH1 *histo)
std::map< ToFCounter, std::vector< unsigned int > > channels
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 PID2DAxis(TH1 *axes, bool third=false)
Double_t GetPercentage(IPrediction *pred_no, IPrediction *pred_pl, IPrediction *pred_mi, int compIdx)
Simple oscillation probability calculator that has no solar term or mass hierarchy or delta so it&#39;s s...
TString compLabel[kNumComp]
const Color_t kBeamNueBackgroundColor
Definition: Style.h:24
double Ebins[11]
Definition: Xdiff_gwt.C:8
void Prelim()
const double j
Definition: BetheBloch.cxx:29
double he
Definition: runWimpSim.h:113
TH2D * histo
void PrintFullShiftTable(int nSyst, Double_t percentages[][kNumSels][kNumVars][kNumComp], int selIdx, int varIdx, TString pidtag, TString out_name, TString *shiftLabel, Bool_t isGenie, Bool_t *turnedOn)
const bool printratioplots
const int kNumSels
Definition: vars.h:43
void FormatErrorBand(TH1 *hplus, TH1 *hminus, bool signal=false, bool fixbins=false)
float bin[41]
Definition: plottest35.C:14
TLatex * prelim
Definition: Xsec_final.C:133
Oscillation probability calculators.
Definition: Calcs.h:5
const std::vector< double > & Edges() const
Definition: Binning.h:30
OStream cout
Definition: OStream.cxx:6
(&#39; survival&#39;)
Definition: IPrediction.h:19
TString MakeLatexCommandName(TString mystring)
static float min(const float a, const float b, const float c)
Definition: absgeo.cxx:45
void PIDBinLabelsShortAxis()
const double kLegendDefaultY2
TH1F * h1
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
ifstream in
Definition: comparison.C:7
TChain * le
Definition: leana.C:12
const Binning kLIDNLBinning
Nonlinear binning for LID.
Definition: Binning.h:114
TH1F * hminus
Definition: Xsec_final.C:76
void PIDDivisionsShortAxis()
const Binning kLEMNLBinning
Nonlinear binning for LEM.
Definition: Binning.h:107
app
Definition: demo.py:189
void PIDBinLabels(std::string pid)
Neutral-current interactions.
Definition: IPrediction.h:40
const double kSecondAnaRunMatchedMCPOTScale
Definition: Exposures.h:96
void PIDTag(TString pid)
Both neutrinos and antineutrinos.
Definition: IPrediction.h:52
c1
Definition: demo5.py:24
const double kLegendDefaultY1
Standard interface to all prediction techniques.
Definition: IPrediction.h:57
Double_t ymin
Definition: plot.C:24
void ComparisonTable(std::vector< TH1 * > mcnom, std::vector< TH1 * > mcplus, std::vector< TH1 * > mcminus, std::vector< TString > labels, TString ltxcommand="")
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:17
std::unique_ptr< T > LoadFrom(TDirectory *dir)
Definition: LoadFromFile.h:18
All neutrinos, any flavor.
Definition: IPrediction.h:26
(&#39; appearance&#39;)
Definition: IPrediction.h:16
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)
const Color_t kNCBackgroundColor
Definition: Style.h:22
void PIDDivisions()
osc::OscCalculatorDumb calc
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:30
void ComparePredictions(IPrediction *pred_no, IPrediction *pred_pl, IPrediction *pred_mi, TString plottitle, TString out_name, TString pidtag, bool printtable=true, bool pidaxis=false)
const double kLegendDefaultX1
void PID2DAxis(TH1 *hist, bool shortaxis)
void PrintRawEventCounts(TDirectory *dpred, TString title)
h
Definition: demo3.py:41
void PIDDivisions(bool shortaxis)