plot_quantile_boundaries_twosamples_2020.C
Go to the documentation of this file.
1 //////////////////////////////////////////////////////////////////////////
2 // plot_quantbound_twosamples.C //
3 // //
4 // print the quantile boundaries from two samples in the same canvas //
5 // fhc = true prints fhc and rhc boundaries with fhc colz background //
6 // uses rhc colz if fhc = false //
7 // //
8 // colzbkg = false prints the two boundaries but no colz background //
9 //////////////////////////////////////////////////////////////////////////
10 
13 #include "Utilities/rootlogon.C"
14 
15 #include "TCanvas.h"
16 #include "TColor.h"
17 #include "TFile.h"
18 #include "TH2.h"
19 #include "TLegend.h"
20 #include "TStyle.h"
21 
22 #include <iostream>
23 
24 using namespace ana;
25 
26 void plot_quantile_boundaries_twosamples_2020(bool isFHC, bool colzbkg = true){
27 
30  std::string name_one, name_two;
31  Color_t color_one, color_two;
32  Style_t line_one = kSolid, line_two = kSolid;
33  int color_pal = kFuchsia;
34  //TColor::InvertPalette();
35 
36  if(isFHC){
37  horn_one = "fhc"; name_one = "FHC";
38  horn_two = "rhc"; name_two = "RHC";
39  color_one = kBlue-7; color_two = kMagenta+1;
40  } else {
41  horn_one = "rhc"; name_one = "RHC";
42  horn_two = "fhc"; name_two = "FHC";
43  color_one = kMagenta+1; color_two = kBlue-7;
44  }
45 
46 
47  const int nBins= 19;
48  const int nQuant= 4;
49 
50  std::string filename1 = (std::string)std::getenv("NUMUDATA_DIR")+"lib/ana2020/Quantiles/quantiles_"+horn_one+"_full_numu2020.root";
51  std::string filename2 = (std::string)std::getenv("NUMUDATA_DIR")+"lib/ana2020/Quantiles/quantiles_"+horn_two+"_full_numu2020.root";
52 
53  TFile* infile1= TFile::Open(pnfs2xrootd(filename1).c_str());
54  TFile* infile2= TFile::Open(pnfs2xrootd(filename2).c_str());
55  TH2* quanthist1 = (TH2*)infile1->Get( "FDSpec2D" );
56  TH2* quanthist2 = (TH2*)infile2->Get( "FDSpec2D" );
57  std::vector< std::vector<double> > QuantBounds1 = GetQuantileBins(quanthist1, kNumuCCOptimisedAxis, kHadEFracAxis, nQuant, 1);
58  std::vector< std::vector<double> > QuantBounds2 = GetQuantileBins(quanthist2, kNumuCCOptimisedAxis, kHadEFracAxis, nQuant, 1);
59 
60  std::vector< TH1* > hist_bound1;
61  std::vector< TH1* > hist_bound2;
62 
63  // fill histograms in different loops, otherwise script messes up
64  for(int thisQuant = 0; thisQuant < nQuant; thisQuant++){
65  TH1* h_temp1 = quanthist1->ProjectionX();
66  h_temp1-> Clear();
67  for(int bin = 1; bin <= nBins; bin++){
68  double bound1 = QuantBounds1[bin][thisQuant];
69  h_temp1-> SetBinContent(bin, bound1);
70  }
71  hist_bound1.push_back(h_temp1);
72  }
73  for(int thisQuant = 0; thisQuant < nQuant; thisQuant++){
74  TH1* h_temp2 = quanthist2->ProjectionX();
75  h_temp2-> Clear();
76  for(int bin = 1; bin <= nBins; bin++){
77  double bound2 = QuantBounds2[bin][thisQuant];
78  h_temp2-> SetBinContent(bin, bound2);
79  }
80  hist_bound2.push_back(h_temp2);
81  }
82 
83  gStyle->SetPalette(color_pal);
84  gStyle->SetPadRightMargin(0.11);
85  new TCanvas("canvas");
86  gPad->SetFillStyle(0);
87  TH2* axes = new TH2F("axes", "", 20, 0.0, 5.0, 10, 0., 1.);
88  std::string name_colz;
89  if(colzbkg){
90  name_colz = "_colzbkg";
91  quanthist1-> Draw("colz");
92  quanthist1->SetTitle("");
93  quanthist1->GetYaxis()->SetTitle("E_{had} / E_{#nu}");
94  quanthist1->GetXaxis()->CenterTitle();
95  quanthist1->GetYaxis()->CenterTitle();
96  }
97  else{
98  axes->SetTitle("");
99  axes->GetXaxis()->SetTitle("Reconstructed Neutrino Energy (GeV)");
100  axes->GetYaxis()->SetTitle("E_{had} / E_{#nu}");
101  axes->GetXaxis()->CenterTitle();
102  axes->GetYaxis()->CenterTitle();
103  axes->Draw();
104  }
105  for(int thisQuant = 1; thisQuant < nQuant; thisQuant++){
106  hist_bound1[thisQuant] -> SetLineColor(color_one);
107  hist_bound2[thisQuant] -> SetLineColor(color_two);
108  hist_bound1[thisQuant] -> SetLineStyle(line_one);
109  hist_bound2[thisQuant] -> SetLineStyle(line_two);
110  hist_bound1[thisQuant] -> Draw("hist same][");
111  hist_bound2[thisQuant] -> Draw("hist same][");
112  }
113  if (colzbkg) {
114  if (isFHC ) CornerLabel( "Neutrino beam" );
115  else CornerLabel( "Antineutrino beam" );
116  }
117  Simulation();
118  TLegend *legend = new TLegend(0.15,0.70,0.30,0.85);
119  legend->AddEntry(hist_bound1[1], Form("%s", name_one.c_str()),"l");
120  legend->AddEntry(hist_bound2[1], Form("%s", name_two.c_str()),"l");
121  legend->SetTextSize(0.05); //no border for legend
122  legend->SetBorderSize(0); //no border for legend
123  legend->SetFillStyle(0); //fill colour is transparent
124  legend->Draw();
125  gPad->Update();
126 
127 
128  TString outname = "quantiles_" + horn_one + "_" + horn_two + "_full_numu2020" + name_colz;
129  gPad-> Print(outname + ".pdf");
130  gPad-> Print(outname + ".png");
131 }
void Simulation()
Definition: tools.h:16
tree Draw("slc.nhit")
int nBins
Definition: plotROC.py:16
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
hmean SetLineStyle(2)
void CornerLabel(std::string &s)
Definition: numu_tools.h:145
void plot_quantile_boundaries_twosamples_2020(bool isFHC, bool colzbkg=true)
const HistAxis kHadEFracAxis("E_{had.} / E_{#nu}", Binning::Simple(200, 0, 1), kHadEFrac)
HistAxis that implements Hadronic Energy fraction binning used by L.Vinton to derive Hadronic Energy ...
Definition: HistAxes.h:30
std::vector< std::vector< double > > GetQuantileBins(TH2 *quantileHist, const HistAxis &independentAxis, const HistAxis &quantileAxis, const unsigned int &numQuantiles, const bool verbose)
Returns a 2D vector First index is the independentAxis bin number Second index is the high bin edge f...
std::string pnfs2xrootd(std::string loc, bool unauth)
Definition: Utilities.cxx:620
std::string getenv(std::string const &name)
hmean SetLineColor(4)
float bin[41]
Definition: plottest35.C:14
const HistAxis kNumuCCOptimisedAxis("Reconstructed Neutrino Energy (GeV)", kNumuCCEOptimisedBinning, kCCE)
HistAxis that implements optimised numuCCE from L. Vinton. See docdb 16332. This was close to &#39;custC&#39;...
Definition: HistAxes.h:24
static bool isFHC
cout<< "--"<< endl;for(Int_t iP=1;iP<=hyz->GetNbinsX();iP++){for(Int_t iC=1;iC<=hyz->GetNbinsY();iC++){if(hyv->GetBinContent(iP, iC)>-999){goal_hyv-> SetBinContent(iP, iC,-(dy[iP-1][iC-1]))
void Print(std::string prefix, std::string name, std::string suffix="")
Definition: nue_pid_effs.C:68
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
const int nQuant
Definition: varsandcuts.h:4
std::string horn_one
Definition: saveFDMCHists.C:33
std::string horn_two
Definition: saveFDMCHists.C:34