Functions | Variables
plot_predictions.C File Reference
#include "includes.h"
#include "numu_tools.h"
#include "3FlavorAna/Plotting/NuePlotStyle.h"

Go to the source code of this file.

Functions

void HornSettings (int quantile, double &pot, double &livetime, Sign::Sign_t &wrongs, double &maxy)
 
void plot_predictions (std::string anaName="fhc", std::string calculator="2018calc", bool nh=true, bool fake=true)
 

Variables

const double pot_fhc = kAna2018FHCPOT
 
const double pot_rhc = kAna2018RHCPOT
 
const double livetime_fhc = kAna2018FHCLivetime
 
const double livetime_rhc = kAna2018RHCLivetime
 
double pot = 0
 
double pot_one = pot_rhc
 
double pot_two = pot_fhc
 
double livetime = 0
 
double livetime_one = livetime_rhc
 
double livetime_two = livetime_fhc
 
double maxy = 0
 
double maxy_one = 1.999
 
double maxy_two = 3.999
 
Sign::Sign_t wrongs = Sign::kNu
 
Sign::Sign_t ws_one = Sign::kNu
 
Sign::Sign_t ws_two = Sign::kAntiNu
 
std::string extrapName = "extrap"
 
std::string horn_one = "rhc"
 
std::string horn_two = "fhc"
 
std::string period_one = "full"
 
std::string period_two = "full"
 
bool twobeams = true
 
int totquant = 8
 
bool DianaIsLazy = true
 
bool errorBand = true
 
Color_t color_band = kViolet-9
 
Color_t color_pred = kViolet+1
 
Color_t color_ws = kBlue-7
 
Color_t color_nc = kBlue-7
 
Color_t color_data = kBlack
 
Color_t color_cos = kGreen+1
 
Color_t color_bkg = kGray+1
 
Style_t line_ws = kSolid
 
Style_t line_nc = kSolid
 
Style_t line_data = kSolid
 
Style_t line_pred = kSolid
 
Style_t line_cos = kSolid
 
Style_t line_bkg = kSolid
 
Style_t marker_ws = kFullCircle
 
Style_t marker_nc = kFullSquare
 
Style_t marker_data = kFullCircle
 
Style_t marker_pred = kFullCircle
 
Style_t marker_cos = kFullTriangleUp
 
Style_t marker_bkg = kFullCircle
 

Function Documentation

void HornSettings ( int  quantile,
double &  pot,
double &  livetime,
Sign::Sign_t wrongs,
double &  maxy 
)

Definition at line 69 of file plot_predictions.C.

References livetime_one, livetime_two, maxy_one, maxy_two, pot_one, pot_two, ws_one, and ws_two.

69  {
70  if(quantile<4){
71  pot = pot_one;
73  wrongs = ws_one;
74  maxy = maxy_one;
75  }
76  else{
77  pot = pot_two;
79  wrongs = ws_two;
80  maxy = maxy_two;
81  }
82 
83 }// end HornSettings
double maxy
Sign::Sign_t ws_one
double pot_one
double maxy_one
double pot_two
double livetime_one
double maxy_two
double livetime
unsigned int quantile(std::vector< std::vector< double >> v, float E, float hadEfrac)
Definition: Toy_analyses.C:480
double livetime_two
Sign::Sign_t ws_two
double pot
Sign::Sign_t wrongs
void plot_predictions ( std::string  anaName = "fhc",
std::string  calculator = "2018calc",
bool  nh = true,
bool  fake = true 
)

Definition at line 87 of file plot_predictions.C.

References std::asin(), calc, canvas(), ana::CenterTitles(), color_band, color_bkg, color_cos, color_data, color_nc, color_pred, color_ws, osc::_NoOscillations< T >::Copy(), CornerLabel(), make_associated_cosmic_defs::cosmics, om::cout, ana::DefaultOscCalc(), DianaIsLazy, drawQuantLabel(), allTimeWatchdog::endl, HTMLTools::entry(), errorBand, extrapName, FakeData(), ana::FillWithDimColor(), ana::getAllAna2018Systs(), make_mec_shifts_plots::GetMaximum(), horn_one, horn_two, HornSettings(), MECModelEnuComparisons::i, ana::Flavors::kAll, ana::Flavors::kAllNuMu, ana::Sign::kAntiNu, ana::Sign::kBoth, ana::Current::kCC, ana::Current::kNC, ana::kNuMu, ana::kNumuAna2018, make_mec_shifts_plots::legend, line_bkg, line_cos, line_data, line_nc, line_pred, line_ws, livetime, livetime_fhc, livetime_one, livetime_two, ana::MakeHistCanvasReady_Quant(), marker_bkg, marker_cos, marker_data, marker_nc, marker_pred, marker_ws, maxy, maxy_one, runNovaSAM::outName, pad1, pad2, pad3, period_one, ana::PimpHist(), ana::PlotWithSystErrorBand_Quant(), pot, pot_fhc, pot_one, pot_two, preliminary(), osc::_IOscCalcAdjustable< T >::SetDmsq32(), ana::SystShifts::SetShift(), osc::_IOscCalcAdjustable< T >::SetTh23(), ana::SplitCanvasQuant(), std::sqrt(), string, systs, totquant, twobeams, wrongs, and ws_one.

91 {
92 
93 
94  bool preliminary = false;
95 
96  if(anaName=="rhc"){
97  totquant = 4;
98  twobeams = false;
99  }
100  // if fhc, change settings so the first 4 quantiles are fhc instead of rhc
101  // dont do anything if rhc as the first 4 quantiles are rhc already
102  if(anaName=="fhc"){
103  horn_one = "fhc";
104  period_one = "full";
105  pot_one = pot_fhc;
108  maxy_one = 3.999;
109  twobeams = false;
110  totquant = 4;
111  }
112 
113  // calculator
114  std::string calcName = "2018calc";
115  std::string hierarchy = "";
117  if(nh){
118  hierarchy = "nh";
119  calc->SetTh23(asin(sqrt(0.558)));
120  calc->SetDmsq32(0.00244);
121  }
122  if(!nh){
123  hierarchy = "ih";
124  calc->SetTh23(asin(sqrt(0.558)));
125  calc->SetDmsq32(-0.00244);
126  }
128  calcfake->SetTh23(asin(sqrt(0.4)));
129  calcfake->SetDmsq32(0.00244);
131  auto calcNoOsc = osc.Copy();
132 
133 
134  // print on screen /////////////////////////////////////////////////////////
135  std::cout << "\n\nmake table with predicted number of events\n" << std::endl;
136  std::cout << "ana name: " << anaName << std::endl;
137  std::cout << "calculator: " << calcName << ", hierarchy: " << hierarchy << std::endl;
138 
139 
140  //inputs
141  TString dout_name = "";
142  std::string ddata_name = "/nova/ana/nu_mu_ana/Ana2018/Data/";
143  std::string dcosm_name = "/nova/ana/nu_mu_ana/Ana2018/Cosmics/";
144  std::string dpred_name = "/nova/ana/nu_mu_ana/Ana2018/Predictions/";
145  std::string fpred_name = "";
146  std::string fpred_suffix = "";
147 
148 
149  // cosmics
150  std::string cosmics = "nocosmics";
151  std::vector<Spectrum> s_cosmics;
152  std::vector<TH1*> hCosmics;
153  std::cout << "\nloading cosmics" << std::endl;
154  std::cout << "" << horn_one << std::endl;
155  std::string fcosmics_name_one = dcosm_name + "cosmics_" + horn_one + "__numu2018.root";
156  TFile fcosmics_one(fcosmics_name_one.c_str());
157  hCosmics.push_back((TH1*)fcosmics_one.Get("cosmics_q0")); // remember, quant0 = all
158  hCosmics.push_back((TH1*)fcosmics_one.Get("cosmics_q1"));
159  hCosmics.push_back((TH1*)fcosmics_one.Get("cosmics_q2"));
160  hCosmics.push_back((TH1*)fcosmics_one.Get("cosmics_q3"));
161  hCosmics.push_back((TH1*)fcosmics_one.Get("cosmics_q4"));
162  auto *cosmics_quant1_one = (TH1*)fcosmics_one.Get("cosmics_q1");
163  auto *cosmics_quant2_one = (TH1*)fcosmics_one.Get("cosmics_q2");
164  auto *cosmics_quant3_one = (TH1*)fcosmics_one.Get("cosmics_q3");
165  auto *cosmics_quant4_one = (TH1*)fcosmics_one.Get("cosmics_q4");
166  s_cosmics.push_back(Spectrum (cosmics_quant1_one, pot_one, livetime_one));
167  s_cosmics.push_back(Spectrum (cosmics_quant2_one, pot_one, livetime_one));
168  s_cosmics.push_back(Spectrum (cosmics_quant3_one, pot_one, livetime_one));
169  s_cosmics.push_back(Spectrum (cosmics_quant4_one, pot_one, livetime_one));
170  if(twobeams){
171  std::string fcosmics_name_two = dcosm_name + "cosmics_" + horn_two + "__numu2018.root";
172  TFile fcosmics_two(fcosmics_name_two.c_str());
173  hCosmics.push_back((TH1*)fcosmics_two.Get("cosmics_q0")); // remember, quant0 = all
174  hCosmics.push_back((TH1*)fcosmics_two.Get("cosmics_q1"));
175  hCosmics.push_back((TH1*)fcosmics_two.Get("cosmics_q2"));
176  hCosmics.push_back((TH1*)fcosmics_two.Get("cosmics_q3"));
177  hCosmics.push_back((TH1*)fcosmics_two.Get("cosmics_q4"));
178  auto *cosmics_quant1_two = (TH1*)fcosmics_two.Get("cosmics_q1");
179  auto *cosmics_quant2_two = (TH1*)fcosmics_two.Get("cosmics_q2");
180  auto *cosmics_quant3_two = (TH1*)fcosmics_two.Get("cosmics_q3");
181  auto *cosmics_quant4_two = (TH1*)fcosmics_two.Get("cosmics_q4");
182  s_cosmics.push_back(Spectrum (cosmics_quant1_two, pot_two, livetime_two));
183  s_cosmics.push_back(Spectrum (cosmics_quant2_two, pot_two, livetime_two));
184  s_cosmics.push_back(Spectrum (cosmics_quant3_two, pot_two, livetime_two));
185  s_cosmics.push_back(Spectrum (cosmics_quant4_two, pot_two, livetime_two));
186  }
187  if(twobeams){
188  for(int quantId = 5; quantId < totquant+1; quantId++){
189  hCosmics[0]->Add(hCosmics[quantId]);
190  }
191  }
192 
193 
194 
195 
196  std::cout << "\nloading predictions" << std::endl;
197  const std::string dpred_name_fhc= dpred_name+"pred_numuconcat_fhc__numu2018.root";
198  const std::string dpred_name_rhc= dpred_name+"pred_numuconcat_rhc__numu2018.root";
199  std::string dpred_name_one = dpred_name_rhc;
200  std::string dpred_name_two = dpred_name_fhc;
201  if(anaName=="fhc") dpred_name_one = dpred_name_fhc;
202 
203  std::vector<PredictionInterp*> predictions;
205  int systSize = systs.size();
206 
207  //predictions here
208  std::cout << "\nloading predictions" << std::endl;
209  std::vector<PredictionSystJoint2018*> predictions;
210  ENu2018ExtrapType extrap = kNuMu;
211  for(int quant = 1; quant < totquant+1; quant++){
212  std::cout << "quantile " << quant << std::endl;
213  if(quant<5) predictions.push_back(new PredictionSystJoint2018(extrap, calc, horn_one, quant, systs, dpred_name_one));
214  else predictions.push_back(new PredictionSystJoint2018(extrap, calc, horn_two, quant-4, systs, dpred_name_two));
215  }
216 
217 
218  /*
219  //predictions old style
220  std::cout << "\nloading predictions" << std::endl;
221  for(int quant = 1; quant < totquant+1; quant++){
222  std::cout << "quantile " << quant << std::endl;
223  if(quant<5){
224  fpred_suffix = extrapName + "_stats__" + horn_one + "_" + calcName + "_" + hierarchy + "__" + period_one + "__numu2018.root";
225  fpred_name = horn_one + "/predictions/" + extrapName + "/" + Form("prediction_quant%d__", quant);
226  std::string input = dpred_name + fpred_name + fpred_suffix;
227  TFile* file = new TFile(input.c_str());
228  predictions.push_back(LoadFrom<PredictionInterp>(file, "prediction")).release();
229  file->Close();
230  }
231  else{
232  fpred_suffix = extrapName + "_stats__" + horn_two + "_" + calcName + "_" + hierarchy + "__" + period_two + "__numu2018.root";
233  fpred_name = horn_two + "/predictions/" + extrapName + "/" + Form("prediction_quant%d__", quant-4);
234  std::string input = dpred_name + fpred_name + fpred_suffix;
235  TFile* file = new TFile(input.c_str());
236  predictions.push_back(LoadFrom<PredictionInterp>(file, "prediction")).release();
237  file->Close();
238  }
239  }//loop quant
240  */
241 
242 
243  // real or fake fake data
244  std::cout << "\nData" << std::endl;
245  std::vector<Spectrum> s_predictions;
246  std::vector<Spectrum> s_fakedata;
247 
248  std::cout << "\nloading fake data" << std::endl;
249  for(int quant = 0; quant < totquant; quant++){
250  std::cout << "quantile " << quant+1 << std::endl;
251  HornSettings(quant, pot, livetime, wrongs, maxy);
252  s_predictions.push_back(predictions[quant]->Predict(calcfake));
253  s_fakedata.push_back(s_predictions[quant].FakeData(pot));
254  }
255  for(int quant = 0; quant < totquant; quant++) s_fakedata[quant] += s_cosmics[quant];
256 
257 
258  std::vector<TH1*> hDataUnscaled;
259  std::vector<TH1*> hData;
260  std::vector<TH1*> hNC;
261  std::vector<TH1*> hCC;
262  std::vector<TH1*> hNuMu;
263  std::vector<TH1*> hNoOsc;
264  std::vector<TH1*> hPred;
265  std::vector<TH1*> hBkg;
266  std::vector<TH1*> hAllBkg;
267  std::vector<TH1*> hWS;
268 
269 
271  hData.push_back(s_fakedata[0].ToTH1(pot));
272  hDataUnscaled.push_back(s_fakedata[0].ToTH1(pot));
273  hNC.push_back(predictions[0]->PredictComponent(calc, Flavors::kAll, Current::kNC, Sign::kBoth).ToTH1(pot));
274  hCC.push_back(predictions[0]->PredictComponent(calc, Flavors::kAll, Current::kCC, Sign::kBoth).ToTH1(pot));
275  hNuMu.push_back(predictions[0]->PredictComponent(calc, Flavors::kAllNuMu, Current::kCC, Sign::kBoth).ToTH1(pot));
276  hBkg.push_back(predictions[0]->Predict(calc).ToTH1(pot));
277  hAllBkg.push_back(predictions[0]->Predict(calc).ToTH1(pot));
278  hPred.push_back(predictions[0]->Predict(calc).ToTH1(pot));
279  hNoOsc.push_back(predictions[0]->Predict(calcNoOsc).ToTH1(pot));
280  hWS.push_back(predictions[0]->PredictComponent(calc, Flavors::kAllNuMu, Current::kCC, wrongs).ToTH1(pot));
281 
282  for(int quantId = 0; quantId < totquant; quantId++){
283  HornSettings(quantId, pot, livetime, wrongs, maxy);
284  hData.push_back(s_fakedata[quantId].ToTH1(pot));
285  hDataUnscaled.push_back(s_fakedata[quantId].ToTH1(pot));
286  hNC.push_back(predictions[quantId]->PredictComponent(calc, Flavors::kAll, Current::kNC, Sign::kBoth).ToTH1(pot));
287  hCC.push_back(predictions[quantId]->PredictComponent(calc, Flavors::kAll, Current::kCC, Sign::kBoth).ToTH1(pot));
288  hNuMu.push_back(predictions[quantId]->PredictComponent(calc, Flavors::kAllNuMu, Current::kCC, Sign::kBoth).ToTH1(pot));
289  hBkg.push_back(predictions[quantId]->Predict(calc).ToTH1(pot));
290  hAllBkg.push_back(predictions[quantId]->Predict(calc).ToTH1(pot));
291  hPred.push_back(predictions[quantId]->Predict(calc).ToTH1(pot));
292  hNoOsc.push_back(predictions[quantId]->Predict(calcNoOsc).ToTH1(pot));
293  hWS.push_back(predictions[quantId]->PredictComponent(calc, Flavors::kAllNuMu, Current::kCC, wrongs).ToTH1(pot));
294  if(quantId > 0){ // already added the first quantile
295  hDataUnscaled[0] ->Add(s_fakedata[quantId].ToTH1(pot));
296  hData[0] ->Add(s_fakedata[quantId].ToTH1(pot));
297  hNC[0] ->Add(predictions[quantId]->PredictComponent(calc, Flavors::kAll, Current::kNC, Sign::kBoth).ToTH1(pot));
298  hCC[0] ->Add(predictions[quantId]->PredictComponent(calc, Flavors::kAll, Current::kCC, Sign::kBoth).ToTH1(pot));
299  hNuMu[0] ->Add(predictions[quantId]->PredictComponent(calc, Flavors::kAllNuMu, Current::kCC, Sign::kBoth).ToTH1(pot));
300  hBkg[0] ->Add(predictions[quantId]->Predict(calc).ToTH1(pot));
301  hAllBkg[0]->Add(predictions[quantId]->Predict(calc).ToTH1(pot));
302  hPred[0] ->Add(predictions[quantId]->Predict(calc).ToTH1(pot));
303  hNoOsc[0] ->Add(predictions[quantId]->Predict(calcNoOsc).ToTH1(pot));
304  hWS[0] ->Add(predictions[quantId]->PredictComponent(calc, Flavors::kAllNuMu, Current::kCC, wrongs).ToTH1(pot));
305  }
306  }
307  // fake stack: draw bkg, cosmics, ws in this order
308  // add bkg and cosmics to ws
309  // add cosmics to bkg
310  for(int quantId = 0; quantId < totquant+1; quantId++){
311  hBkg[quantId]->Add(hNuMu[quantId], -1);
312  hAllBkg[quantId]->Add(hNuMu[quantId], -1);
313  hAllBkg[quantId]->Add(hCosmics[quantId]);
314  hWS[quantId]->Add(hAllBkg[quantId]);
315  hPred[quantId]->Add(hCosmics[quantId]);
316  }
317 
318 
319 
320  // scale everything, everything!
321  std::cout << "\nScaling histograms to 0.1 GeV" << std::endl;
322  for(int quantId = 0; quantId < totquant+1; quantId++){
323  hData[quantId] ->Scale(0.1, "width");
324  hCosmics[quantId]->Scale(0.1, "width");
325  hNC[quantId] ->Scale(0.1, "width");
326  hCC[quantId] ->Scale(0.1, "width");
327  hBkg[quantId] ->Scale(0.1, "width");
328  hAllBkg[quantId] ->Scale(0.1, "width");
329  hNuMu[quantId] ->Scale(0.1, "width");
330  hPred[quantId] ->Scale(0.1, "width");
331  hNoOsc[quantId] ->Scale(0.1, "width");
332  hWS[quantId] ->Scale(0.1, "width");
333  }
334 
335 
336 
337  std::vector<Spectrum> sPrediction;
338  for(int quantId = 0; quantId < totquant+1; quantId++){
339  if(quantId==0){
341  Spectrum prediction(hPred[quantId], pot, livetime);
342  sPrediction.push_back(prediction);
343  }
344  else{
345  HornSettings(quantId-1, pot, livetime, wrongs, maxy);
346  Spectrum prediction(hPred[quantId], pot, livetime);
347  sPrediction.push_back(prediction);
348  }
349  }
350 
351 
352  // save up and dn shifts for each quantile
353  // looks messy but really makes it less difficult
354  std::vector<Spectrum> ups[totquant+1], dns[totquant+1];
355  std::cout << "\nFilling ups and dns\n" << std::endl;
356 
357  if(errorBand){
358  for(int quantId = 0; quantId < totquant+1; quantId++){
359 
360  std::cout << "quantile " << quantId << std::endl;
361  // spectrums for syst band best fit
362  for(const ISyst* syst: systs){
363  SystShifts shifts;
364 
365  if(quantId==0){
367  shifts.SetShift(syst, +1);
368  TH1 * htempup = predictions[0]->PredictSyst(calc, shifts).ToTH1(pot);
369  for(int i = 1; i < totquant; i++){htempup->Add(predictions[i]->PredictSyst(calc, shifts).ToTH1(pot));}
370  htempup->Scale(0.1, "width");
371  htempup->Add(hCosmics[0]);
372  Spectrum specup(htempup, pot, livetime);
373  ups[quantId].push_back(specup);
374 
375  shifts.SetShift(syst, -1);
376  TH1 * htempdn = predictions[0]->PredictSyst(calc, shifts).ToTH1(pot);
377  for(int i = 1; i < totquant; i++){htempdn->Add(predictions[i]->PredictSyst(calc, shifts).ToTH1(pot));}
378  htempdn->Scale(0.1, "width");
379  htempdn->Add(hCosmics[0]);
380  Spectrum specdn(htempdn, pot, livetime);
381  dns[quantId].push_back(specdn);
382  }// all quantiles
383 
384  else{
385  HornSettings(quantId-1, pot, livetime, wrongs, maxy);
386  shifts.SetShift(syst, +1);
387  TH1* htempup = predictions[quantId-1]->PredictSyst(calc, shifts).ToTH1(pot);
388  htempup->Scale(0.1, "width");
389  htempup->Add(hCosmics[quantId]);
390  Spectrum specup(htempup, pot, livetime);
391  ups[quantId].push_back(specup);
392  shifts.SetShift(syst, -1);
393  TH1 * htempdn = predictions[quantId-1]->PredictSyst(calc, shifts).ToTH1(pot);
394  htempdn->Scale(0.1, "width");
395  htempdn->Add(hCosmics[quantId]);
396  Spectrum specdn(htempdn, pot, livetime);
397  dns[quantId].push_back(specdn);
398  }
399  } // shifts
400  }// for quantiles
401  }
402 
403 
404 
405  // pimping
406  for(int quantId = 0; quantId< 1; quantId++){
407  PimpHist(hData[quantId], line_data, color_data, 3, marker_data, color_data, 1);
408  PimpHist(hWS[quantId], line_ws, color_ws, 1, marker_ws, color_ws, 1);
409  PimpHist(hNC[quantId], line_nc, color_nc, 1, marker_nc, color_nc, 1);
410  PimpHist(hPred[quantId], line_pred, color_pred, 3, marker_pred, color_pred, 1);
411  PimpHist(hCosmics[quantId], line_cos, color_cos, 1, marker_cos, color_cos, 1);
412  PimpHist(hAllBkg[quantId], line_bkg, color_bkg, 1, marker_bkg, color_bkg, 1);
413  PimpHist(hBkg[quantId], line_bkg, color_bkg, 1, marker_bkg, color_bkg, 1);
414  FillWithDimColor(hWS[quantId], 0);
415  FillWithDimColor(hCosmics[quantId], 0);
416  FillWithDimColor(hAllBkg[quantId], 0);
417  }// for quantId
418  //root is not behaving for the quantiles
419  //set different line width and marker size
420  for(int quantId = 1; quantId< totquant+1; quantId++){
421  PimpHist(hData[quantId], line_data, color_data, 2, marker_data, color_data, 2);
422  PimpHist(hWS[quantId], line_ws, color_ws, 2, marker_ws, color_ws, 2);
423  PimpHist(hNC[quantId], line_nc, color_nc, 2, marker_nc, color_nc, 2);
424  PimpHist(hPred[quantId], line_pred, color_pred, 2, marker_pred, color_pred, 2);
425  PimpHist(hCosmics[quantId], line_cos, color_cos, 2, marker_cos, color_cos, 2);
426  PimpHist(hAllBkg[quantId], line_bkg, color_bkg, 2, marker_bkg, color_bkg, 2);
427  PimpHist(hBkg[quantId], line_bkg, color_bkg, 2, marker_bkg, color_bkg, 2);
428  FillWithDimColor(hWS[quantId], 0);
429  FillWithDimColor(hCosmics[quantId], 0);
430  FillWithDimColor(hAllBkg[quantId], 0);
431  }// for quantId
432 
433 
434  // for later use in the 4pad plots
435  TLatex* xTitle = new TLatex(0.5, 0.03, "Reconstructed Energy (GeV)");
436  xTitle->SetNDC();
437  xTitle->SetTextAlign(22);
438  TLatex* yTitle = new TLatex(0.02, 0.5, "Events / 0.1 GeV");
439  yTitle->SetNDC();
440  yTitle->SetTextAlign(22);
441  yTitle->SetTextAngle(90.);
442 
443 
444  const int cx = 600, cy = 500;
445  // all quantiles together
446  // equals quant 0
447  if(DianaIsLazy){
448  TString outName = dout_name + "prediction_spectra__" + anaName;
449  std::string idName = anaName + "_" + extrapName;
450  int quantId = 0;
451  TCanvas *canvas = new TCanvas("canvas","canvas", cx, cy);
452  canvas->cd();
453  ana::CenterTitles(hPred[0]);
454  hPred[0]->GetYaxis()->SetTitle("Events / 0.1 GeV");
455  hPred[0]->GetYaxis()->SetRangeUser(0.0, 1.3 * hPred[0]->GetMaximum());
457  if(errorBand) PlotWithSystErrorBand_Quant(quantId, sPrediction[quantId], ups[quantId], dns[quantId], pot, color_pred, color_band, hPred[quantId]->GetMaximum(), true);
458  else hPred[quantId]->Draw("hist same");
459  hWS[quantId]->Draw("hist same");
460  hAllBkg[quantId]->Draw("hist same");
461  hCosmics[quantId]->Draw("hist same");
462  hPred[quantId]->Draw("hist same");
463  //hAllBkg[quantId]->Draw("l p same");
464  //hCosmics[quantId]->Draw("hist p same");
465  TLegend *legend = new TLegend(0.58,0.60,0.78,0.85,"","NDC");
466  legend->AddEntry(hPred[quantId], "Prediction","l");
467  TLegendEntry *entry=legend->AddEntry("error","1-#sigma syst. range","bf");
468  entry->SetFillStyle(1001);
469  entry->SetFillColor(color_band);
470  entry->SetLineColor(color_band);
471  legend->AddEntry(hAllBkg[quantId], "Total bkgd.","bf");
472  legend->AddEntry(hCosmics[quantId], "Cosmic bkgd.","bf");
473  legend->AddEntry(hWS[quantId], "WS","bf");
474  //legend->AddEntry(hBkg[quantId], "Beam bkg.","lp");
475  //legend->AddEntry(hCosmics[quantId], "Cosmic bkg.","lp");
476  //legend->AddEntry(hWS[quantId], "WS","l");
477  legend->SetTextSize(0.05); //no border for legend
478  legend->SetBorderSize(0); //no border for legend
479  legend->SetFillStyle(0); //fill colour is white
480  legend->Draw();
481  CornerLabel(idName);
482  drawQuantLabel(0);
483  gPad->RedrawAxis();
484  gPad->Update();
485  canvas->cd();
486  if(errorBand){
487  xTitle->Draw();
488  yTitle->Draw();
489  }
490 
491  canvas->SaveAs(outName + ".eps");
492  canvas->SaveAs(outName + ".png");
493  canvas->SaveAs(outName + ".pdf");
494 
495  } // plot all quantiles
496 
497 
498  if(DianaIsLazy){
499 
500  TCanvas *canvas;
501  TPad *pad1, *pad2, *pad3, *pad4;
502  SplitCanvasQuant(canvas, pad1, pad2, pad3, pad4);
503  TString outName = dout_name + "prediction_spectra__" + anaName;
504  std::string idName = anaName + "_" + extrapName;
505  for(int quantId = 1; quantId < 5; quantId++){
506  if(quantId==1) pad1->cd();
507  if(quantId==2) pad2->cd();
508  if(quantId==3) pad3->cd();
509  if(quantId==4) pad4->cd();
510  HornSettings(quantId-1, pot, livetime, wrongs, maxy);
511  MakeHistCanvasReady_Quant(quantId, hPred[quantId], maxy);
512  if(errorBand) PlotWithSystErrorBand_Quant(quantId, sPrediction[quantId], ups[quantId], dns[quantId], pot, color_pred, color_band, maxy, true);
513  else hPred[quantId]->Draw("hist same");
514  hWS[quantId]->Draw("hist same");
515  hAllBkg[quantId]->Draw("hist same");
516  hCosmics[quantId]->Draw("hist same");
517  hPred[quantId]->Draw("hist same");
518  //hAllBkg[quantId]->Draw("l p same");
519  //hCosmics[quantId]->Draw("hist p same");
520  //grData[quantId] ->Draw("p0 same");
521  drawQuantLabel(quantId);
522  if(quantId==1){
523  pad1->cd();
524  TLegend *legend = new TLegend(0.09,0.70,0.27,0.88,"","NDC");
525  legend->AddEntry(hPred[quantId], "Prediction","l");
526  TLegendEntry *entry=legend->AddEntry("error","1-#sigma syst. range","bf");
527  entry->SetFillStyle(1001);
528  entry->SetFillColor(color_band);
529  entry->SetLineColor(color_band);
530  legend->AddEntry(hAllBkg[quantId], "Total bkgd.","bf");
531  legend->AddEntry(hCosmics[quantId], "Cosmic bkgd.","bf");
532  legend->AddEntry(hWS[quantId], "WS","bf");
533  //legend->AddEntry(hBkg[quantId], "Beam bkg.","lp");
534  //legend->AddEntry(hCosmics[quantId], "Cosmic bkg.","lp");
535  //legend->AddEntry(hWS[quantId], "WS","l");
536  legend->SetTextSize(0.025); //no border for legend
537  legend->SetBorderSize(0); //no border for legend
538  legend->SetFillStyle(0); //fill colour is white
539  legend->Draw();
540  CornerLabel(idName);
541  gPad->Update();
542  }
543  }// for quantiles
544 
545 
546  for(auto & p:{pad1,pad2,pad3,pad4}){
547  p->RedrawAxis();
548  }
549  canvas->Update();
550  canvas->cd();
551  xTitle->Draw();
552  yTitle->Draw();
553  canvas->SaveAs(outName + "_quant14.png");
554  canvas->SaveAs(outName + "_quant14.eps");
555  canvas->SaveAs(outName + "_quant14.pdf");
556 
557  } // plot quantiles 1 to 4
558 
559 
560 
561  if(DianaIsLazy && twobeams){
562  TCanvas *canvas;
563  TPad *pad1, *pad2, *pad3, *pad4;
564  SplitCanvasQuant(canvas, pad1, pad2, pad3, pad4);
565  TString outName = dout_name + "datapred_bestfit__" + anaName;
566  std::string idName = anaName + "_" + extrapName;
567  for(int quantId = 5; quantId <totquant+1; quantId++){ // revers prevents pad overlaping
568  if(quantId==5) pad1->cd();
569  if(quantId==6) pad2->cd();
570  if(quantId==7) pad3->cd();
571  if(quantId==8) pad4->cd();
572  HornSettings(quantId-1, pot, livetime, wrongs, maxy);
573  MakeHistCanvasReady_Quant(quantId, hPred[quantId], maxy);
574  if(errorBand) PlotWithSystErrorBand_Quant(quantId-4, sPrediction[quantId], ups[quantId], dns[quantId], pot, color_pred, color_band, maxy, true);
575  else hPred[quantId]->Draw("hist same");
576  hWS[quantId]->Draw("hist same");
577  hAllBkg[quantId]->Draw("hist same");
578  hCosmics[quantId]->Draw("hist same");
579  hPred[quantId]->Draw("hist same");
580  drawQuantLabel(quantId-4);
581  if(quantId==1){
582  pad1->cd();
583  TLegend *legend = new TLegend(0.09,0.70,0.27,0.88,"","NDC");
584  legend->AddEntry(hPred[quantId], "Prediction","l");
585  TLegendEntry *entry=legend->AddEntry("error","1-#sigma syst. range","bf");
586  entry->SetFillStyle(1001);
587  entry->SetFillColor(color_band);
588  entry->SetLineColor(color_band);
589  legend->AddEntry(hAllBkg[quantId], "Total bkgd.","bf");
590  legend->AddEntry(hCosmics[quantId], "Cosmic bkgd.","bf");
591  legend->AddEntry(hWS[quantId], "WS","bf");
592  //legend->AddEntry(hNC[quantId], "NC bkg.","lp");
593  //legend->AddEntry(hBkg[quantId], "Beam bkg.","lp");
594  //legend->AddEntry(hCosmics[quantId], "Cosmic bkg.","lp");
595  //legend->AddEntry(hWS[quantId], "WS","l");
596  legend->SetTextSize(0.025); //no border for legend
597  legend->SetBorderSize(0); //no border for legend
598  legend->SetFillStyle(0); //fill colour is white
599  legend->Draw();
600  CornerLabel(idName);
601  gPad->Update();
602  }
603  }// for quantiles
604 
605  for(auto & p:{pad1,pad2,pad3,pad4}){
606  p->RedrawAxis();
607  }
608  canvas->Update();
609  canvas->cd();
610  xTitle->Draw();
611  yTitle->Draw();
612  canvas->SaveAs(outName + "_quant58.png");
613  canvas->SaveAs(outName + "_quant58.eps");
614  canvas->SaveAs(outName + "_quant58.pdf");
615 
616  } // plot quantiles 1 to 4
617 
618 
619 
620 } // End of function
Style_t marker_pred
Color_t color_data
Color_t color_pred
Pass neutrinos through unchanged.
Definition: IOscCalc.h:45
const double pot_fhc
Antineutrinos-only.
Definition: IPrediction.h:50
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
void MakeHistCanvasReady_Quant(const int quant, TH1 *hist, double ymax)
Definition: Plots.cxx:1365
double maxy
int totquant
std::string extrapName
Loads shifted spectra from files.
Sign::Sign_t ws_one
Style_t line_data
Simple record of shifts applied to systematic parameters.
Definition: SystShifts.h:20
double pot_one
const double livetime_fhc
const char * p
Definition: xmltok.h:285
bool twobeams
T sqrt(T number)
Definition: d0nt_math.hpp:156
void CornerLabel(std::string &s)
Definition: numu_tools.h:145
Color_t color_bkg
void FakeData()
Definition: rootlogon.C:156
Color_t color_band
void HornSettings(int quantile, double &pot, double &livetime, Sign::Sign_t &wrongs, double &maxy)
Style_t line_cos
void CenterTitles(TH1 *histo)
Definition: Plots.cxx:1483
TCanvas * canvas(const char *nm, const char *ti, const char *a)
Definition: pass1_plots.C:36
double maxy_one
osc::OscCalcDumb calc
osc::IOscCalcAdjustable * DefaultOscCalc()
Create a new calculator with default assumptions for all parameters.
Definition: Calcs.cxx:49
Encapsulate code to systematically shift a caf::SRProxy.
Definition: ISyst.h:14
virtual void SetDmsq32(const T &dmsq32)=0
Style_t line_bkg
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
void PimpHist(TH1 *hist, Style_t linestyle, Color_t linecolor, int linewidth, Style_t markerstyle, Color_t markercolor, double markersize)
Pimp histogram once and for all.
Definition: Plots.cxx:1408
double pot_two
double livetime_one
void preliminary()
std::vector< const ISyst * > getAllAna2018Systs(const EAnaType2018 ana, const bool smallgenie, const BeamType2018 beam, bool isFit)
Charged-current interactions.
Definition: IPrediction.h:39
Style_t line_nc
double livetime
Color_t color_nc
virtual _IOscCalc< T > * Copy() const override
Definition: IOscCalc.h:49
Style_t marker_bkg
Style_t line_pred
std::string horn_two
Style_t marker_ws
Oscillation probability calculators.
Definition: Calcs.h:5
OStream cout
Definition: OStream.cxx:6
Style_t line_ws
Style_t marker_nc
Color_t color_ws
Style_t marker_data
double livetime_two
std::string period_one
bool errorBand
TGraphAsymmErrors * PlotWithSystErrorBand_Quant(const int quant, IPrediction *pred, const std::vector< const ISyst * > &systs, osc::IOscCalc *calc, double pot, int col, int errCol, float maxy, bool newaxis)
Definition: Plots.cxx:543
TPaveText * drawQuantLabel(int quantId=0)
Definition: numu_tools.h:4
virtual void SetTh23(const T &th23)=0
std::string horn_one
TPad * pad3
Definition: analysis.C:13
void SplitCanvasQuant(TCanvas *&canvas, TPad *&pad1, TPad *&pad2, TPad *&pad3, TPad *&pad4)
Definition: Plots.cxx:1439
double pot
bool DianaIsLazy
TPad * pad2
Definition: analysis.C:13
Color_t color_cos
All neutrinos, any flavor.
Definition: IPrediction.h:26
Style_t marker_cos
def entry(str)
Definition: HTMLTools.py:26
Sign::Sign_t wrongs
void SetShift(const ISyst *syst, double shift, bool force=false)
Definition: SystShifts.cxx:81
void FillWithDimColor(TH1 *h, bool usealpha, float dim)
TPad * pad1
Definition: analysis.C:13
T asin(T number)
Definition: d0nt_math.hpp:60
enum BeamMode string

Variable Documentation

Color_t color_band = kViolet-9

Definition at line 46 of file plot_predictions.C.

Referenced by plot_predictions().

Color_t color_bkg = kGray+1

Definition at line 52 of file plot_predictions.C.

Referenced by plot_predictions().

Color_t color_cos = kGreen+1

Definition at line 51 of file plot_predictions.C.

Referenced by plot_predictions().

Color_t color_data = kBlack

Definition at line 50 of file plot_predictions.C.

Referenced by plot_predictions().

Color_t color_nc = kBlue-7

Definition at line 49 of file plot_predictions.C.

Referenced by plot_predictions().

Color_t color_pred = kViolet+1

Definition at line 47 of file plot_predictions.C.

Referenced by plot_predictions().

Color_t color_ws = kBlue-7

Definition at line 48 of file plot_predictions.C.

Referenced by plot_predictions().

bool DianaIsLazy = true

Definition at line 41 of file plot_predictions.C.

bool errorBand = true

Definition at line 42 of file plot_predictions.C.

std::string extrapName = "extrap"

Definition at line 32 of file plot_predictions.C.

std::string horn_one = "rhc"

Definition at line 33 of file plot_predictions.C.

Referenced by plot_predictions().

std::string horn_two = "fhc"

Definition at line 34 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t line_bkg = kSolid

Definition at line 58 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t line_cos = kSolid

Definition at line 57 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t line_data = kSolid

Definition at line 55 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t line_nc = kSolid

Definition at line 54 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t line_pred = kSolid

Definition at line 56 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t line_ws = kSolid

Definition at line 53 of file plot_predictions.C.

Referenced by plot_predictions().

double livetime = 0

Definition at line 21 of file plot_predictions.C.

Referenced by plot_predictions().

const double livetime_fhc = kAna2018FHCLivetime

Definition at line 12 of file plot_predictions.C.

Referenced by plot_predictions().

double livetime_one = livetime_rhc

Definition at line 22 of file plot_predictions.C.

Referenced by HornSettings(), and plot_predictions().

const double livetime_rhc = kAna2018RHCLivetime

Definition at line 13 of file plot_predictions.C.

double livetime_two = livetime_fhc

Definition at line 23 of file plot_predictions.C.

Referenced by HornSettings(), and plot_predictions().

Style_t marker_bkg = kFullCircle

Definition at line 64 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t marker_cos = kFullTriangleUp

Definition at line 63 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t marker_data = kFullCircle

Definition at line 61 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t marker_nc = kFullSquare

Definition at line 60 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t marker_pred = kFullCircle

Definition at line 62 of file plot_predictions.C.

Referenced by plot_predictions().

Style_t marker_ws = kFullCircle

Definition at line 59 of file plot_predictions.C.

Referenced by plot_predictions().

double maxy = 0

Definition at line 24 of file plot_predictions.C.

double maxy_one = 1.999

Definition at line 25 of file plot_predictions.C.

double maxy_two = 3.999

Definition at line 26 of file plot_predictions.C.

std::string period_one = "full"

Definition at line 35 of file plot_predictions.C.

Referenced by plot_predictions().

std::string period_two = "full"

Definition at line 36 of file plot_predictions.C.

double pot = 0

Definition at line 18 of file plot_predictions.C.

Referenced by plot_datamcpred(), and plot_predictions().

const double pot_fhc = kAna2018FHCPOT

Definition at line 10 of file plot_predictions.C.

Referenced by plot_predictions().

double pot_one = pot_rhc

Definition at line 19 of file plot_predictions.C.

Referenced by HornSettings(), and plot_predictions().

const double pot_rhc = kAna2018RHCPOT

Definition at line 11 of file plot_predictions.C.

double pot_two = pot_fhc

Definition at line 20 of file plot_predictions.C.

Referenced by HornSettings(), and plot_predictions().

int totquant = 8

Definition at line 40 of file plot_predictions.C.

Referenced by plot_predictions().

bool twobeams = true

Definition at line 39 of file plot_predictions.C.

Referenced by plot_predictions().

Sign::Sign_t wrongs = Sign::kNu

Definition at line 27 of file plot_predictions.C.

Referenced by plot_predictions().

Sign::Sign_t ws_one = Sign::kNu

Definition at line 28 of file plot_predictions.C.

Referenced by HornSettings(), and plot_predictions().

Sign::Sign_t ws_two = Sign::kAntiNu

Definition at line 29 of file plot_predictions.C.

Referenced by HornSettings().