32 #include "TGraphAsymmErrors.h" 41 TLegend*
Legend(
double x0,
double y0,
double x1,
double y1,
bool useData,
bool smartLeg=
false,
bool mc=
false)
43 if (!smartLeg) {leg =
new TLegend(x0, y0, x1, y1);}
46 leg->SetFillStyle(1001);
47 leg->SetFillColor(kWhite);
49 TH1*
dummy =
new TH1F(
"",
"", 1, 0, 1);
50 TH1* dummyFill =
new TH1F(
"",
"", 1, 0, 1);
53 dummy->SetMarkerStyle(kFullCircle);
54 if( useData) leg->AddEntry(dummy->Clone(),
"FD data",
"lep");
55 if(!useData) leg->AddEntry(dummy->Clone(),
"FD mock data",
"lep");
57 leg->AddEntry(dummy->Clone(),
"Best Fit prediction",
"l");
58 dummy->SetLineColor(kGray+1);
60 leg->AddEntry(dummy->Clone(),
"Total Background",
"bf");
61 dummy->SetLineColor(kAzure+1);
63 leg->AddEntry(dummy->Clone(),
"Cosmic Background",
"bf");}
68 leg->AddEntry(dummy->Clone(),
"signal #nu_{e} CC",
"bf");
69 dummy->SetLineColor(kPink+9);
71 leg->AddEntry(dummy->Clone(),
"beam #nu_{e} CC",
"bf");
72 dummy->SetLineColor(kAzure);
74 leg->AddEntry(dummy->Clone(),
"NC",
"bf");
75 TColor *color4 =
new TColor(232, 0.85, 0.75, 0.58);
76 dummy->SetLineColor(
kGreen+2);
78 leg->AddEntry(dummy->Clone(),
"#nu_{#mu} CC",
"bf");
79 dummy->SetLineColor(kMagenta-9);
81 leg->AddEntry(dummy->Clone(),
"#nu_{#tau} CC",
"bf");
82 dummy->SetLineColor(kMagenta+3);
84 leg->AddEntry(dummy->Clone(),
"#bar{#nu}_{e} CC",
"bf");
85 dummy->SetLineColor(kAzure+1);
87 leg->AddEntry(dummy->Clone(),
"cosmic",
"bf");
113 if(sideband) {end = 6;}
114 if(!sideband) {start = 6;}
128 calc->SetDmsq32(2.44
e-3);
129 for(
int selIdx = 0; selIdx <
kNumSels; ++selIdx){
130 const char* selName =
selNames[selIdx].c_str();
131 for(
int varIdx = 0; varIdx <
kNumVars; ++varIdx){
132 const char* varName =
defs[varIdx].
name.c_str();
133 spectMC[selIdx][varIdx] = LoadFromFile<PredictionNoExtrap>(fnameMC,
TString::Format(
"%s/%s", selName, varName).Data()).
release();
134 cosmic[selIdx][varIdx] = LoadFromFile<Spectrum>(fnameCo,
TString::Format(
"%s/cosm_%s", selName, varName).Data()).
release();
135 data[selIdx][varIdx] = LoadFromFile<Spectrum>(fnameDa,
TString::Format(
"%s/spect_%s", selName, varName).Data()).
release();
136 rock[selIdx][varIdx] = LoadFromFile<PredictionNoExtrap>(fnameRo,
TString::Format(
"%s/rock_%s", selName, varName).Data()).
release();
137 extrapCombo[selIdx][varIdx]=LoadFromFile<PredictionExtrap>(fnameMC,
TString::Format(
"%s/%s_extrap_combo", selName, varName).Data()).
release();
156 for(
int selIdx = start; selIdx <
end; ++selIdx){
158 for(
int varIdx = 0; varIdx <
kNumVars; ++varIdx){
159 const char* varName =
defs[varIdx].
name.c_str();
160 double nonsScale = 1./10.45;
161 double swapScale = 1./12.91;
163 auto data1 = data[selIdx][varIdx];
164 auto cosmic1 = cosmic[selIdx][varIdx];
165 auto rock1 = rock[selIdx][varIdx];
166 const double POT=data[selIdx][varIdx]->
POT();
167 double Livetime = data[selIdx][varIdx]->
Livetime();
171 auto noextrap = spectMC[selIdx][varIdx];
173 rock1 = rock[selIdx][21];
174 data1 = data[selIdx][21];
175 cosmic1 = cosmic[selIdx][21];
177 else {
auto spectMC2 = extrapCombo[selIdx][varIdx];
cout<<
"hi"<<
endl;}
186 hMC = spectMC1->Predict(
calc).ToTH1(POT);
188 bool useScales = !useExtrap &&
true;
190 std::cout <<
"Rescaling MC components to Extrap/No Extrap TODO";
191 hMC->Scale(57.5/57.1);
192 hNue->Scale(42.2/43.1);
193 hBeam->Scale(7.09/7.03);
194 hNC->Scale(6.62/5.52);
196 hTau->Scale(0.446/0.446);
197 hNue_bg->Scale(1.0/1.0);
199 hTotbkg = (TH1*)hNue_bg->Clone(
UniqueName().c_str());
205 hTotbkg->SetLineColor(kGray+1);
208 hRock = spectMC1->Predict(
calc).ToTH1(POT);
210 hRock->Add(hJustMC, -1);
214 hbkg = (TH1D*) hTotbkg->Clone(
UniqueName().c_str());
217 hCos->SetLineColor(kAzure+1);
221 std::cout <<
"\n\nWARNING useData=false; plotting mock data\n\n";
222 hData = spectMC1->Predict(
calc).ToTH1(POT,kBlack,kSolid);
237 if((varIdx==2 || varIdx==20 )&& selIdx==4){
238 for(
auto &
h:{hNue,
hbkg,hTotbkg,hCos,hMC}){
239 if(varIdx==2)
h->Rebin(10);
240 if(varIdx==20)
h->Rebin(15);
242 if(!mc) {
if(varIdx==2) hData->Rebin(10);
243 if(varIdx==20) hData->Rebin(15);}
246 if(!mc) {hData->SetMarkerStyle(kFullCircle);}
249 if (selIdx == 6 && varIdx == 2) {rebin = 2;}
250 if (varIdx==20) {rebin = 3;}
251 if (varIdx==19 || varIdx==6 || varIdx==7 || varIdx==8) {rebin = 4;}
252 if (varIdx==11 || varIdx==12 || varIdx==13 || varIdx==14 || varIdx==15 || varIdx==16 || varIdx==17) {rebin = 10;}
253 if (!mc) {hData->Rebin(rebin);}
254 for (
auto &
h:{hMC, hNue,
hbkg, hTotbkg, hNue_bg, hNC, hCC, hBeam, hTau, hRock, hCos}){
258 cout<<
"MC "<<hMC->Integral()<<
" signal "<<hNue->Integral()<<
" total beam "<<hbkg->Integral()<<
" cosmic "<<hCos->Integral()<<
"rock "<<hRock->Integral()<<
endl;
289 hMC->SetMaximum(1.1*hMC->GetMaximum());
290 if(!mc) {
int n = hData->GetNbinsX();
291 for (
int i=1;
i<=
n;
i++) {
if(hMC->GetBinContent(
i) == 0) hData->SetBinContent(
i,-1);}
293 if (varIdx==0) {mult = 1.13;}
294 if ( varIdx==7 || varIdx==8) {mult = 1.6;}
295 if (varIdx==6) {mult = 1.8;}
296 hMC->SetMaximum(mult*
max(hMC->GetMaximum(),(hData->GetMaximum()+
sqrt(hData->GetMaximum()))));}
297 hMC->GetXaxis()->CenterTitle();
298 if (varIdx==1) {hMC->GetXaxis()->SetRange(1,12);}
299 if (varIdx==5) {hMC->GetXaxis()->SetRange(1,20);}
300 if (selIdx==6 && varIdx==2) {hMC->GetXaxis()->SetRange(50.0/rebin+1,100.0/rebin);}
301 hMC->GetYaxis()->CenterTitle();
302 hMC->GetYaxis()->SetTitleOffset(1.15);
303 hMC->GetYaxis()->SetTitle(
"Events / 8.85 #times 10^{20} POT-equiv");
304 if(varIdx == 6){hMC->GetXaxis()->SetTitle(
"Vertex X (cm)");}
305 if(varIdx == 7){hMC->GetXaxis()->SetTitle(
"Vertex Y (cm)");}
306 if(varIdx == 8){hMC->GetXaxis()->SetTitle(
"Vertex Z (cm)");}
307 if(varIdx == 19){hMC->GetXaxis()->SetTitle(
"p_{T}/p");}
308 if(varIdx == 20){hMC->GetXaxis()->SetRange(1,80/rebin);}
312 TCanvas *
c =
new TCanvas(
"c",
"c",1000,800);
313 c->SetLeftMargin(0.15);
314 c->SetRightMargin(0.05);
316 gPad->SetFillStyle(0);
358 hTotbkg->SetLineColor(kAzure);
360 hTotbkg->Draw(
"hist same");
361 TH1* wonc = (TH1*)hTotbkg->Clone(
UniqueName().c_str());
362 wonc->SetLineColor(kPink+9);
365 wonc->Draw(
"hist same");
366 TH1* wobeam = (TH1*)wonc->Clone(
UniqueName().c_str());
367 wobeam->SetLineColor(
kGreen+2);
369 wobeam->Add(hBeam,-1);
370 wobeam->Draw(
"hist same");
371 TH1* wocc = (TH1*)wobeam->Clone(
UniqueName().c_str());
372 wocc->SetLineColor(kMagenta-9);
375 wocc->Draw(
"hist same");
376 TH1* wotau = (TH1*)wocc->Clone(
UniqueName().c_str());
377 wotau->SetLineColor(kMagenta+3);
380 wotau->Draw(
"hist same");
381 TH1* wown = (TH1*)wotau->Clone(
UniqueName().c_str());
382 wown->SetLineColor(kGray);
384 wown->Add(hNue_bg,-1);
385 wown->Draw(
"hist same");
386 TH1* worock = (TH1*)wown->Clone(
UniqueName().c_str());
387 worock->SetLineColor(kAzure+1);
390 worock->Draw(
"hist same");
394 hTotbkg->Draw(
"hist same");
395 hCos->Draw(
"hist same");
396 hMC->Draw(
"hist same");
398 gr->SetMarkerStyle(kFullCircle);
400 gr->Draw(
"ep same");}
404 double legendxlow=0.2;
405 double legendxhigh=0.5;
406 double legendylow=0.5;
407 double legendyhigh=0.7;
410 if (varIdx==0) { marker =
false; legendxlow=0.2; legendxhigh=0.55; legendylow=0.6; legendyhigh=0.8; }
411 if (selIdx == 5 && varIdx == 4 || selIdx == 4 && varIdx == 20 || selIdx == 4 && varIdx == 2 || selIdx == 0 && varIdx == 3 || selIdx == 0 && varIdx == 4 ||
412 selIdx == 2 && varIdx == 3) { marker =
false; legendxlow=0.17; legendxhigh=0.465; legendylow=0.6; legendyhigh=0.8; }
413 if (varIdx==20) { marker =
false; legendxlow=0.2; legendxhigh=0.5; legendylow=0.6; legendyhigh=0.8; }
414 if (varIdx==21) { marker =
false; legendxlow=0.2; legendxhigh=0.5; legendylow=0.5; legendyhigh=0.7; }
415 if (mc) {legendylow=0.5; legendyhigh=0.8;}
417 auto leg =
Legend(legendxlow, legendylow, legendxhigh, legendyhigh, useData, marker, mc);
418 leg->SetTextSize(hMC->GetXaxis()->GetLabelSize());
428 if (selName==
"Peripheral") {
auto *
tt =
new TText(gPad->GetUxmin(), gPad->GetUymax() * 1.05,
"Peripheral Sample");
430 if (selName==
"Core") {
auto *
tt =
new TText(gPad->GetUxmin(), gPad->GetUymax() * 1.05 ,
"Core Sample");
432 if (selName==
"Presel") {
auto *
tt =
new TText(gPad->GetUxmin(), gPad->GetUymax() * 1.05 ,
"Preselection cut tier");
435 if (selIdx == 5) {
auto *
tt =
new TText(gPad->GetUxmin(), gPad->GetUymax() * 1.02 ,
"High Energy Sideband (Core)");
tt->SetTextSize(0.03);
tt->Draw();}
436 if (selIdx == 4) {
auto *
tt =
new TText(gPad->GetUxmin(), gPad->GetUymax() * 1.02 ,
"Low PID Sideband (Core)");
tt->SetTextSize(0.03);
tt->Draw();}
437 if (selIdx == 0) {
auto *
tt =
new TText(gPad->GetUxmin(), gPad->GetUymax() * 1.01 ,
"High Energy, low PID Sideband (Peripheral)");
tt->SetTextSize(0.03);
tt->Draw(); }
438 if (selIdx == 2) {
auto *
tt =
new TText(gPad->GetUxmin(), gPad->GetUymax() * 1.02 ,
"Low PID Sideband (Peripheral)");
tt->SetTextSize(0.03);
tt->Draw();}
444 ofstream
caption((
"plots/sideband/" + selName +
"_" + varName +
".txt").c_str());
445 caption<<selName<<
" cut applied, sideband data/prediction study, with extrapolation from the ND and combo decomposition, 2017 Best Fit oscillation parameters.";
446 c->Print((
"plots/sideband/" + selName +
"_" + varName +
".pdf").c_str());}
448 ofstream
caption((
"plots/fdmc/fdmc_" + selName +
"_" + varName +
".txt").c_str());
449 caption<<selName<<
" cut applied, FDMC with extrapolation from the ND and combo decomposition, 2017 Best Fit oscillation parameters.";
450 c->Print((
"plots/fdmc/fdmc_" + selName +
"_" + varName +
".pdf").c_str());}
451 if(useExtrap & !mc & !sideband) {
452 ofstream
caption((
"plots/extrap/decomp_" + selName +
"_" + varName + (useData?
"_data":
"_fake") +
".txt").c_str());
453 caption<<selName<<
" cut applied, data/prediction comparisons, with extrapolation from the ND and combo decomposition, 2017 Best Fit oscillation parameters.";
454 c->Print((
"plots/extrap/decomp_" + selName +
"_" + varName + (useData?
"_data":
"_fake") +
".pdf").c_str());}
void Nue2017FourBinLabels(const double yNDC, const double textsize, const int color, const bool merged)
void Nue2017FourBinAxis(TH1 *axes, bool drawLabels, bool merged)
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.
virtual Spectrum Predict(osc::IOscCalc *calc) const =0
Float_t x1[n_points_granero]
TLegend * Legend(double x0, double y0, double x1, double y1, bool useData, bool smartLeg=false, bool mc=false)
TGraphAsymmErrors * GraphWithPoissonErrors(const TH1 *h, bool noErrorsXaxis, bool drawEmptyBins)
Calculate statistical errors appropriate for small Poisson numbers.
const Color_t kTotalMCErrorBandColor
osc::IOscCalcAdjustable * DefaultOscCalc()
Create a new calculator with default assumptions for all parameters.
Spectrum MockData(double pot, int seed=0) const
Mock data is FakeData with Poisson fluctuations applied.
Representation of a spectrum in any variable, with associated POT.
const XML_Char const XML_Char * data
Charged-current interactions.
Sum up livetimes from individual cosmic triggers.
const HistDef defs[kNumVars]
void fd_plot(std::string fnameMC, std::string fnameCo, std::string fnameDa, std::string fnameRo, bool useData=false)
std::vector< double > POT
Neutral-current interactions.
Both neutrinos and antineutrinos.
Standard interface to all prediction techniques.
TLegend * AutoPlaceLegend(double dx, double dy, double yPin)
Create a legend, maximizing distance from all histograms.
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
void Nue2017FourBinDivisions(const int color, const int style)
All neutrinos, any flavor.
const std::string selNames[kNumSels]
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
double Livetime() const
Seconds. For informational purposes only. No calculations use this.
std::string UniqueName()
Return a different string each time, for creating histograms.
void FillWithDimColor(TH1 *h, bool usealpha, float dim)
void rock(std::string suffix="full")