5 #include "CAFAna/Core/Binning.h" 6 #include "CAFAna/Core/Cut.h" 7 #include "CAFAna/Core/HistAxis.h" 16 #include "CAFAna/Cuts/NueCutsSecondAna.h" 53 TLegend*
leg =
new TLegend(x0, y0, x1, y1);
55 leg->SetTextSize(2/30.);
58 leg->SetHeader(
"#splitline{#bf{Neutrino Beam}}{#nu_{e} signal+WS}");
59 TH1*
dummy =
new TH1F(
"",
"", 1, 0, 1);
60 TH1* dummyFill =
new TH1F(
"",
"", 1, 0, 1);
61 leg->AddEntry(dummy->Clone(),
"",
"");
66 leg->AddEntry(dummy->Clone(),
"Nominal",
"l");
67 dummy->SetLineStyle(9);
69 leg->AddEntry(dummy->Clone(),
"Q^{2} reweighted",
"l");
70 dummy->SetLineStyle(8);
71 dummy->SetLineColor(kGray+2);
72 leg->AddEntry(dummy->Clone(),
"p_{T}/p reweighted",
"l");
73 dummy->SetLineStyle(3);
75 leg->AddEntry(dummy->Clone(),
"cos #theta reweighted",
"l");
82 TH1* coreperiTH1 = (TH1*)coreTH1->Clone(
"coreperi");
83 for(
int i = 18;
i < 28; ++
i){
84 coreperiTH1->SetBinContent(
i, periTH1->GetBinContent(
i));
93 gStyle->SetPadTickX(1);
94 gStyle->SetPadTickY(0);
96 TFile * file_no =
new TFile (filename_no.c_str(),
"READ");
97 if(file_no->IsZombie()){
std::cerr <<
"Can't find " << filename_no.c_str() <<
std::endl;
return;}
98 TFile * file_re_core =
new TFile (filename_re_core.c_str(),
"READ");
99 if(file_re_core->IsZombie()){
std::cerr <<
"Can't find " << filename_re_core.c_str() <<
std::endl;
return;}
100 TFile * file_re_peri =
new TFile (filename_re_peri.c_str(),
"READ");
101 if(file_re_peri->IsZombie()){
std::cerr <<
"Can't find " << filename_re_peri.c_str() <<
std::endl;
return;}
105 TFile*
file =
new TFile(
"FD_KinematicsCorrection_pTExtrap_FHC.root",
"RECREATE");
107 TDirectory* dpred_noextrap = file_no->GetDirectory(
"NOReweight/nueAxis_NoExtrap");
112 std::vector<const PredictionExtendToPeripheral* > pred_no_quant;
113 std::vector<const PredictionExtendToPeripheral* > pred_re_core_q2_quant;
114 std::vector<const PredictionExtendToPeripheral* > pred_re_core_PtP_quant;
115 std::vector<const PredictionExtendToPeripheral* > pred_re_core_Cos_quant;
116 std::vector<const PredictionExtendToPeripheral* > pred_re_peri_q2_quant;
117 std::vector<const PredictionExtendToPeripheral* > pred_re_peri_PtP_quant;
118 std::vector<const PredictionExtendToPeripheral* > pred_re_peri_Cos_quant;
120 const unsigned int nQuants = 3;
122 for (
unsigned int quant = 0; quant < nQuants; ++quant ){
124 TDirectory* dpred_no_quant = file_no->GetDirectory(Form(
"NOReweight/nueAxis_NueSignalExtrap_Quant%i",quant+1));
126 pred_no_quant.push_back(extend_nom);
129 TDirectory* dpred_re_core_q2_quant = file_re_core->GetDirectory(Form(
"CORE/nueAxis_NueSignalExtrap_Q2Weight_Quant%i", quant+1));
131 pred_re_core_q2_quant.push_back(extend_core_q2);
133 TDirectory* dpred_re_core_PtP_quant = file_re_core->GetDirectory(Form(
"CORE/nueAxis_NueSignalExtrap_PtPWeight_Quant%i", quant+1));
135 pred_re_core_PtP_quant.push_back(extend_core_PtP);
137 TDirectory* dpred_re_core_Cos_quant = file_re_core->GetDirectory(Form(
"CORE/nueAxis_NueSignalExtrap_CosWeight_Quant%i", quant+1));
139 pred_re_core_Cos_quant.push_back(extend_core_Cos);
142 TDirectory* dpred_re_peri_q2_quant = file_re_peri->GetDirectory(Form(
"PERI/nueAxis_NueSignalExtrap_Q2Weight_Quant%i", quant+1));
144 pred_re_peri_q2_quant.push_back(extend_peri_q2);
146 TDirectory* dpred_re_peri_PtP_quant = file_re_peri->GetDirectory(Form(
"PERI/nueAxis_NueSignalExtrap_PtPWeight_Quant%i", quant+1));
148 pred_re_peri_PtP_quant.push_back(extend_peri_PtP);
150 TDirectory* dpred_re_peri_Cos_quant = file_re_peri->GetDirectory(Form(
"PERI/nueAxis_NueSignalExtrap_CosWeight_Quant%i", quant+1));
152 pred_re_peri_Cos_quant.push_back(extend_peri_Cos);
158 Spectrum spec_nom = pred_no_quant.front()->PredictComponent(&mycalc,
164 Spectrum spec_core_q2 = pred_re_core_q2_quant.front()->PredictComponent(&mycalc,
168 spec_core_q2.
Clear();
170 Spectrum spec_core_PtP = pred_re_core_PtP_quant.front()->PredictComponent(&mycalc,
174 spec_core_PtP.
Clear();
176 Spectrum spec_core_Cos = pred_re_core_Cos_quant.front()->PredictComponent(&mycalc,
180 spec_core_Cos.
Clear();
182 Spectrum spec_peri_q2 = pred_re_peri_q2_quant.front()->PredictComponent(&mycalc,
186 spec_peri_q2.
Clear();
188 Spectrum spec_peri_PtP = pred_re_peri_PtP_quant.front()->PredictComponent(&mycalc,
192 spec_peri_PtP.
Clear();
194 Spectrum spec_peri_Cos = pred_re_peri_Cos_quant.front()->PredictComponent(&mycalc,
198 spec_peri_Cos.
Clear();
202 for (
unsigned int quant = 0; quant < nQuants; ++quant ){
203 spec_nom += pred_no_quant[quant]->PredictComponent(&mycalc,
208 spec_core_q2 += pred_re_core_q2_quant[quant]->PredictComponent(&mycalc,
213 spec_core_PtP += pred_re_core_PtP_quant[quant]->PredictComponent(&mycalc,
218 spec_core_Cos += pred_re_core_Cos_quant[quant]->PredictComponent(&mycalc,
223 spec_peri_q2 += pred_re_peri_q2_quant[quant]->PredictComponent(&mycalc,
228 spec_peri_PtP += pred_re_peri_PtP_quant[quant]->PredictComponent(&mycalc,
233 spec_peri_Cos += pred_re_peri_Cos_quant[quant]->PredictComponent(&mycalc,
242 const double kPOTFD = 13.5e20;
244 TH1* hAll_no = spec_nom.
ToTH1(kPOTFD);
246 TH1* hAll_core_q2 = spec_core_q2.
ToTH1(kPOTFD);
247 TH1* hAll_peri_q2 = spec_peri_q2.
ToTH1(kPOTFD);
250 TH1* hAll_core_PtP = spec_core_PtP.
ToTH1(kPOTFD);
251 TH1* hAll_peri_PtP = spec_peri_PtP.
ToTH1(kPOTFD);
254 TH1* hAll_core_Cos = spec_core_Cos.
ToTH1(kPOTFD);
255 TH1* hAll_peri_Cos = spec_peri_Cos.
ToTH1(kPOTFD);
259 hAll_re_q2->SetLineColor(
kOrange+1);
260 hAll_re_q2->SetLineStyle(9);
263 hAll_re_PtP->SetLineColor(kGray+2);
264 hAll_re_PtP->SetLineStyle(8);
265 hAll_re_PtP->SetLineWidth(1);
267 hAll_re_Cos->SetLineColor(kGray+2);
268 hAll_re_Cos->SetLineStyle(3);
269 hAll_re_Cos->SetLineWidth(1);
272 TCanvas *
c2 =
new TCanvas(
"c2",
"c2",500,400);
274 TPad *padAbs =
new TPad(
"padAbs",
"padAbs", 0, 0.3, 1, 1.0);
275 padAbs-> SetBottomMargin(0);
276 padAbs->SetFillStyle(0);
277 padAbs-> SetLeftMargin(0.15);
282 hAll_re_q2->GetXaxis()->CenterTitle();
283 hAll_re_q2->GetYaxis()->CenterTitle();
284 hAll_re_q2->GetXaxis()->SetDecimals();
285 hAll_re_q2->GetYaxis()->SetDecimals();
286 hAll_re_q2->GetYaxis()->SetTitleOffset(.8);
287 hAll_re_q2->GetYaxis()->SetTitleSize(.06);
288 hAll_re_q2->GetYaxis()->SetLabelSize(.06);
289 hAll_re_q2->GetXaxis()->SetLabelSize(.0);
290 hAll_re_q2->SetMaximum(1.1*hAll_re_q2->GetMaximum());
291 hAll_re_q2->GetYaxis()->SetTitle(
TString::Format(
"Events / %.02lf #times 10^{21} POT", kPOTFD/1E21));
293 hAll_no->Draw(
"hist");
294 hAll_re_PtP->Draw(
"same hist");
295 hAll_re_Cos->Draw(
"same hist");
296 hAll_re_q2->Draw(
"same hist");
299 hAll_re_q2->GetYaxis()->SetRangeUser(.001,20);
302 padAbs->RedrawAxis();
305 double legendxlow=0.55;
306 double legendxhigh=0.95;
308 Legend(legendxlow, .5, legendxhigh, .85);
310 TLatex* selTitle =
new TLatex(.9, .95,
"NOvA Simulation");
312 selTitle->SetTextColor(kGray+1);
313 selTitle->SetTextSize(2/30.);
314 selTitle->SetTextAlign(32);
318 TPad *padRatio =
new TPad(
"padRatio",
"padRatio", 0, 0.0, 1, 0.33);
319 padRatio -> SetTopMargin(0.11);
320 padRatio -> SetBottomMargin(0.33);
322 padRatio -> SetLeftMargin(0.15);
327 TH1* hRatio_q2 = (TH1*)hAll_re_q2 ->
Clone(
"hRatio_q2");
328 TH1* hRatio_PtP = (TH1*)hAll_re_PtP ->
Clone(
"hRatio_PtP");
329 TH1* hRatio_Cos = (TH1*)hAll_re_Cos ->
Clone(
"hRatio_Cos");
330 TLine* lOne =
new TLine(0,1.0,27.0,1.0);
333 hRatio_q2 ->
Divide(hAll_no);
334 hRatio_PtP ->
Divide(hAll_no);
335 hRatio_Cos ->
Divide(hAll_no);
337 hRatio_q2->GetXaxis()->CenterTitle();
338 hRatio_q2->GetYaxis()->CenterTitle();
339 hRatio_q2->GetYaxis()->SetLabelSize(.12);
340 hRatio_q2->GetXaxis()->SetLabelSize(.12);
341 hRatio_q2->GetXaxis()->SetTitleOffset(1.0);
342 hRatio_q2->GetXaxis()->SetTickSize(.1);
343 hRatio_q2->GetYaxis()->SetTitleSize(.12);
344 hRatio_q2->GetXaxis()->SetTitleSize(.12);
345 hRatio_q2->SetMarkerSize(.1);
346 hRatio_q2->GetXaxis()->SetDecimals();
347 hRatio_q2->GetYaxis()->SetDecimals();
348 hRatio_q2->GetYaxis()->SetTitleOffset(0.4);
349 hRatio_q2->GetYaxis()->SetTitle(
"#splitline{Ratio to}{nominal}");
351 hRatio_q2->GetYaxis()->SetNdivisions(502,kFALSE);
352 hRatio_q2 ->
Draw(
"hist");
353 hRatio_PtP ->
Draw(
"same hist");
354 hRatio_Cos ->
Draw(
"same hist");
355 lOne ->
Draw(
"same");
358 hRatio_q2->GetYaxis()->SetRangeUser(0.95,1.05);
360 c2->Print((
"plots/FDnue_pred_recoE_pTExtrap_FHC.png"));
361 c2->SaveAs(
"plots/FDnue_pred_recoE_pTExtrap_FHC.pdf",
"pdf");
363 hRatio_q2->
Write(
"trueQ2_FD_KinematicsCorrection_FHC");
364 hRatio_PtP->
Write(
"PtP_FD_KinematicsCorrection_FHC");
365 hRatio_Cos->
Write(
"Cos_FD_KinematicsCorrection_FHC");
369 double intAll_no = hAll_no->Integral();
370 double intAll_re_q2 = hAll_re_q2->Integral();
371 double intAll_re_PtP = hAll_re_PtP->Integral();
372 double intAll_re_Cos = hAll_re_Cos->Integral();
381 tfile.open(
"plots/FDextrap-TableOut-pTExtrap-FHC.tex",
ios::out);
383 tfile <<
"\\begin{tabular}{lll} \n";
384 tfile <<
"\\multicolumn{3}{c}{prop FD $\nu_e$ signal + WS} \\\\ \n";
385 tfile <<
"Sample & Signal+WS & /NOMINAL \\\\ \n";
386 tfile <<
"\\hline \n";
390 tfile <<
"NOMINAL" <<
" & " <<
round(100.0*intAll_no)/100.0 <<
" & \\\\ \n";
391 tfile <<
"trueQ2 REW." <<
" & " <<
round(100.0*intAll_re_q2)/100.0 <<
" & " <<
round(1000.0*((intAll_re_q2/intAll_no)-1.0))/10.0 <<
" \\% \\\\ \n";
392 tfile <<
"PtP REW." <<
" & " <<
round(100.0*intAll_re_PtP)/100.0 <<
" & " <<
round(1000.0*((intAll_re_PtP/intAll_no)-1.0))/10.0 <<
" \\% \\\\ \n";
393 tfile <<
"Cos REW." <<
" & " <<
round(100.0*intAll_re_Cos)/100.0 <<
" & " <<
round(1000.0*((intAll_re_Cos/intAll_no)-1.0))/10.0 <<
" \\% \\\\ \n";
394 tfile <<
"\\hline \n";
395 tfile <<
"\\end{tabular}";
void Nue2018ThreeBinAxis(THStack *axes, bool drawLabels, bool merged, bool coreOnly)
ratio_hxv Divide(hxv, goal_hxv)
Cuts and Vars for the 2020 FD DiF Study.
Float_t y1[n_points_granero]
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.
Float_t x1[n_points_granero]
const Color_t kTotalMCErrorBandColor
fvar< T > round(const fvar< T > &x)
ntuple SetFillStyle(1001)
Representation of a spectrum in any variable, with associated POT.
const Color_t kNueSignalColor
Charged-current interactions.
const Binning kNue2020Binning
Simple oscillation probability calculator that has no solar term or mass hierarchy or delta so it's s...
Both neutrinos and antineutrinos.
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
const Color_t kTotalMCColor