nd_plot_bless.C
Go to the documentation of this file.
2 
4 #include "CAFAna/Core/Ratio.h"
6 #include "CAFAna/Core/Spectrum.h"
9 
10 #include "CAFAna/Vars/Vars.h"
11 #include "CAFAna/Vars/NueVars.h"
12 #include "CAFAna/Analysis/Style.h"
13 using namespace ana;
14 
16 
17 #include "TCanvas.h"
18 #include "TFile.h"
19 #include "TH1.h"
20 #include "TLatex.h"
21 #include "TLine.h"
22 #include "TLegend.h"
23 #include "TPad.h"
24 
25 //make a pretty legend
26 void Legend(double x0, double y0, double x1, double y1)
27 {
28  // Doesn't handle log-y well
29  // TLegend* leg = AutoPlaceLegend(x1-x0, y1-y0);
30  TLegend* leg = new TLegend(x0, y0, x1, y1);
31  leg->SetFillStyle(0);
32 
33  TH1* dummy = new TH1F("", "", 1, 0, 1);
34  TH1* dummyFill = new TH1F("", "", 1, 0, 1);
35  dummyFill->SetLineColor(kTotalMCColor);
36  dummyFill->SetFillColor(kTotalMCErrorBandColor);
37 
38  dummy->SetMarkerStyle(kFullCircle);
39 
40  leg->AddEntry(dummy->Clone(), "ND data", "lep");
41  dummy->SetLineColor(kGray+2);
42  dummy->SetMarkerColor(kGray+2);
43  // leg->AddEntry(dummy->Clone(), "ND data stagger", "lep");
44  dummy->SetLineColor(kTotalMCColor);
45  leg->AddEntry(dummy->Clone(), "Total MC", "l");
46  leg->AddEntry(dummyFill->Clone(), "Flux Uncert.", "fl");
47  dummy->SetLineColor(kNCBackgroundColor);
48  leg->AddEntry(dummy->Clone(), "NC", "l");
49  dummy->SetLineColor(kBeamNueBackgroundColor);
50  leg->AddEntry(dummy->Clone(), "Beam #nu_{e} CC", "l");
51  dummy->SetLineColor(kNumuBackgroundColor);
52  leg->AddEntry(dummy->Clone(), "#nu_{#mu} CC", "l");
53 
54  leg->Draw();
55 }
56 
57 //create a nova prelimary tag
59 {
60  TLatex* prelim = new TLatex(.935, .95, "NOvA Preliminary");
61  prelim->SetTextColor(kBlue);
62  prelim->SetNDC();
63  prelim->SetTextSize(2/30.);
64  prelim->SetTextAlign(32);
65  prelim->Draw();
66 }
67 
69 {
70  struct HistDef
71  {
73  };
74 
75  //which variables do you actually want to look at?
76  const int kNumVars =14;
77  const HistDef defs[kNumVars] = {
78  "nhits",
79  "inelast",
80  "vtxx",
81  "vtxy",
82  "vtxz",
83  "shw_len",
84  "shw_width",
85  "recoE",
86  "nshw",
87  "shw_e",
88  "had_e",
89  "ptp",
90  "cvne",
91  "cvne_zoom"
92  };
93 
94  //which cut tiers do you want to look at?
95  const int kNumSels = 2;
96  const std::string selNames[kNumSels] = {"presel","cvneloose"};
97 
99  IPrediction* preds[kNumSels][kNumVars];
102 
104 
105  //pick up spectra/prediction objects from the input root file
106  for(int selIdx = 0; selIdx < kNumSels; ++selIdx){
107  const char* selName = selNames[selIdx].c_str();
108  for(int varIdx = 0; varIdx < kNumVars; ++varIdx){
109  const char* varName = defs[varIdx].name.c_str();
110 
111  std::cout<<TString::Format("%s/spect_%s", selName, varName).Data()<<std::endl;
112  spects[selIdx][varIdx] = LoadFromFile<Spectrum>(fname, TString::Format("%s/spect_%s", selName, varName).Data()).release();
113 
114  preds[selIdx][varIdx] = LoadFromFile<IPrediction>(fname, TString::Format("%s/pred_%s", selName, varName).Data()).release();
115  ups[selIdx][varIdx] = LoadFromFile<IPrediction>(fname, TString::Format("%s/up_%s", selName, varName).Data()).release();
116  downs[selIdx][varIdx] = LoadFromFile<IPrediction>(fname, TString::Format("%s/down_%s", selName, varName).Data()).release();
117 
118  if(selIdx == 0 && varIdx == 0)
119  std::cout << "data, MC POT: "
120  << spects[selIdx][varIdx]->POT() << " "
121  << preds[selIdx][varIdx]->Predict(&noosc).POT() << std::endl;
122 
123  } // end for varIdx
124  } // end for selIdx
125 
126  //loop over the different selector/cut levels to plot each one
127  for(int selIdx = 0; selIdx < kNumSels; ++selIdx){
128  const std::string selName = selNames[selIdx];
129  for(int varIdx = 0; varIdx < kNumVars; ++varIdx){
130 
131  const std::string varName = defs[varIdx].name;
132 
133  //play with the binning choice for given variables:
134  int rebinFactor=1;
135 
136  if(varName=="recoE"){
137  rebinFactor=4;
138  }
139 
140  //create a canvas
141  TCanvas *c2 =new TCanvas("c2","c2",500,400);
142  c2->SetLeftMargin(0.15); //default is 10%
143  c2->SetRightMargin(0.05);
144 
145  gPad->SetFillStyle(0);
146 
147  double legendxlow=0.7;
148  double legendxhigh=0.95;
149 
150  //Pull down the variable we want
151  Spectrum* spect = spects[selIdx][varIdx];
152 
153  IPrediction* pred = preds[selIdx][varIdx];
154  IPrediction* up = ups[selIdx][varIdx];
155  IPrediction* down = downs[selIdx][varIdx];
156 
157  //Get components of the prediction and store them in histograms
158  TH1* hMC = pred->Predict(&noosc).ToTH1(spect->POT());
159  hMC->SetLineColor(kTotalMCColor);
160  hMC->Rebin(rebinFactor);
161 
162 
163  TH1* hPlus = up->Predict(&noosc).ToTH1(spect->POT());
164  hPlus->Rebin(rebinFactor);
165 
166  TH1* hMinus = down->Predict(&noosc).ToTH1(spect->POT());
167  hMinus->Rebin(rebinFactor);
168 
169  TH1* hNC = pred->PredictComponent(&noosc,
171  Current::kNC,
172  Sign::kBoth).ToTH1(spect->POT());
173  hNC->SetLineColor(kNCBackgroundColor);
174  hNC->Rebin(rebinFactor);
175 
176  TH1* hCC = pred->PredictComponent(&noosc,
178  Current::kCC,
179  Sign::kBoth).ToTH1(spect->POT());
180  hCC->SetLineColor(kNumuBackgroundColor);
181  hCC->Rebin(rebinFactor);
182 
183  TH1* hBeam = pred->PredictComponent(&noosc,
185  Current::kCC,
186  Sign::kBoth).ToTH1(spect->POT());
187  hBeam->SetLineColor(kBeamNueBackgroundColor);
188  hBeam->Rebin(rebinFactor);
189 
190  TH1* hData = spect->ToTH1(spect->POT());
191  hData->SetMarkerStyle(kFullCircle);
192  hData->Rebin(rebinFactor);
193 
194  hPlus ->SetFillColor(kTotalMCErrorBandColor);
195  hPlus ->SetLineColor(kTotalMCErrorBandColor);
196  hMinus->SetFillColor(10);
197  hMinus->SetLineColor(kTotalMCErrorBandColor);
198 
199  hPlus->GetXaxis()->CenterTitle();
200  hPlus->GetYaxis()->CenterTitle();
201  hPlus->GetXaxis()->SetDecimals();
202  hPlus->GetYaxis()->SetDecimals();
203  hPlus->GetYaxis()->SetTitleOffset(1.15);
204  hPlus->GetYaxis()->SetTitle(TString::Format("Events / %.02lf #times 10^{20} POT", spect->POT()/1E20));
205  if(varName=="had_e"){
206  hPlus->GetXaxis()->SetTitle("Hadronic calorimetric energy (GeV)");
207  }
208  if(varName=="cvne"){
209  hPlus->GetXaxis()->SetTitle("CVN #nu_{e} classifier");
210  }
211  if(varName=="cvne_zoom"){
212  hPlus->GetXaxis()->SetTitle("CVN #nu_{e} classifier");
213  }
214 
215  if(varName=="inelast"){
216  hPlus->GetXaxis()->SetTitle("Reconstructed inelasticity");
217  }
218 
219  if(varName=="recoE"){
220  hPlus->GetXaxis()->SetTitle("Reconstructed neutrino energy (GeV)");
221  }
222 
223  hPlus->GetXaxis()->SetNdivisions(405,kFALSE);
224 
225  //Draw the error band, MC, Data, MC components, and then Data one last time
226  if(varName=="vtxx"||varName=="vtxy"||varName=="vtxz"){
227  hPlus->Scale(1./1000.);
228  hMinus->Scale(1./1000.);
229  hMC->Scale(1./1000.);
230  hData->Scale(1./1000.);
231  hNC->Scale(1./1000.);
232  hCC->Scale(1./1000.);
233  hBeam->Scale(1./1000.);
234  hPlus->GetYaxis()->SetTitle(TString::Format("10^{3} Events / %.02lf #times 10^{20} POT", spect->POT()/1E20));
235  hPlus->SetMaximum(1.7*hPlus->GetMaximum());
236  }
237 
238  if(varName=="cvne_zoom"||varName=="shw_len"){
239  legendxlow=0.2;
240  legendxhigh=0.55;
241  }
242  hPlus ->SetMinimum(.001);
243 
244  if(varName=="shw_width"){
245  hPlus->GetXaxis()->SetRangeUser(0,15);
246  }
247 
248  if(varName=="had_e"){
249  hPlus->GetXaxis()->SetRangeUser(0,2);
250  }
251 
252  hPlus ->DrawCopy("hist");
253  hMinus->Draw("hist same ");
254 
255  hMC->Draw("hist same");
256  hData->Draw("ep same");
257 
258  hNC->Draw("hist same");
259  hCC->Draw("hist same");
260  hBeam->Draw("hist same");
261  hMC->Draw("hist same");
262  hData->Draw("ep same");
263 
264  //Fix up the axis
265  gPad->RedrawAxis();
266  gPad->Update();
267 
268  //Add the legend and preliminary
269  Legend(legendxlow, .6, legendxhigh, .85);
270  Preliminary();
271 
272  gPad->Print(("plots/"+selName+"_"+varName+"_bless.pdf").c_str());
273  gPad->Print(("plots/"+selName+"_"+varName+"_bless.eps").c_str());
274  gPad->Print(("plots/"+selName+"_"+varName+"_bless.png").c_str());
275  gPad->Print(("plots/"+selName+"_"+varName+"_bless.gif").c_str());
276  gPad->Print(("plots/"+selName+"_"+varName+"_bless.C").c_str());
277 
278  c2->SetLogy();
279  if(varName=="cvne"){
280  hData->SetMarkerSize(0.5);
281  }
282 
283  hPlus ->SetMinimum(10);
284  hPlus ->SetMaximum(10*hPlus ->GetMaximum());
285  hPlus ->DrawCopy("hist");
286  hMinus->Draw("hist same ");
287 
288  hMC->Draw("hist same");
289  hData->Draw("ep same");
290 
291  hNC->Draw("hist same");
292  hCC->Draw("hist same");
293  hBeam->Draw("hist same");
294  hMC->Draw("hist same");
295  hData->Draw("ep same");
296 
297  //Fix up the axis
298  gPad->RedrawAxis();
299  gPad->Update();
300 
301  //Add the legend and preliminary
302  Legend(legendxlow, .6, legendxhigh, .85);
303  Preliminary();
304 
305 
306  gPad->Print(("plots/"+selName+"_"+varName+"_bless_log.pdf").c_str());
307  gPad->Print(("plots/"+selName+"_"+varName+"_bless_log.eps").c_str());
308  gPad->Print(("plots/"+selName+"_"+varName+"_bless_log.png").c_str());
309  gPad->Print(("plots/"+selName+"_"+varName+"_bless_log.gif").c_str());
310  gPad->Print(("plots/"+selName+"_"+varName+"_bless_log.C").c_str());
311 
312  } // end for varIdx
313  } // end for selIdx
314 
315 }
Pass neutrinos through unchanged.
const int kNumVars
Definition: vars.h:14
Oscillation analysis framework, runs over CAF files outside of ART.
Float_t y1[n_points_granero]
Definition: compare.C:5
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
Float_t x1[n_points_granero]
Definition: compare.C:5
(&#39;beam &#39;)
Definition: IPrediction.h:15
const Color_t kTotalMCErrorBandColor
Definition: Style.h:18
const Color_t kNumuBackgroundColor
Definition: Style.h:31
virtual Spectrum Predict(osc::IOscCalculator *calc) const =0
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:37
c2
Definition: demo5.py:33
Charged-current interactions.
Definition: IPrediction.h:39
void nd_plot_bless(std::string fname)
Definition: nd_plot_bless.C:68
const Color_t kBeamNueBackgroundColor
Definition: Style.h:24
const HistDef defs[kNumVars]
Definition: vars.h:15
const int kNumSels
Definition: vars.h:43
osc::NoOscillations noosc
Definition: NueSystFuncs.h:44
TLatex * prelim
Definition: Xsec_final.C:133
double POT() const
Definition: Spectrum.h:263
std::string name
Definition: NuePlotLists.h:12
OStream cout
Definition: OStream.cxx:6
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
const std::map< std::string, RowDef > downs
void Preliminary()
Definition: nd_plot_bless.C:58
void Legend(double x0, double y0, double x1, double y1)
Definition: nd_plot_bless.C:26
Neutral-current interactions.
Definition: IPrediction.h:40
Both neutrinos and antineutrinos.
Definition: IPrediction.h:52
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:17
void spects(osc::IOscCalculatorAdjustable *calc, IPrediction *predFHC, IPrediction *predRHC)
Definition: sensitivity.C:100
All neutrinos, any flavor.
Definition: IPrediction.h:26
const std::string selNames[kNumSels]
Definition: vars.h:46
const Color_t kNCBackgroundColor
Definition: Style.h:22