example_plot.C
Go to the documentation of this file.
1 // Plot the histograms made using example_macro.C
2 // Author: Cathal Sweeney - csweeney@fnal.gov
3 
4 
5 #include "CAFAna/Core/Spectrum.h"
6 
7 #include "TFile.h"
8 #include "TH1.h"
9 #include "TH2.h"
10 #include "TCanvas.h"
11 #include "TGraphAsymmErrors.h"
12 #include "TString.h"
13 #include "TColor.h"
14 #include "TGaxis.h"
15 #include "TROOT.h"
16 #include "TStyle.h"
17 #include "TLegend.h"
18 #include "TLegendEntry.h"
19 
20 // for plotting
21 #include "CAFAna/Analysis/Plots.h"
22 
23 namespace ana
24 {
25 
26  std::vector<std::pair<std::string,std::string>> systName_list = {
27  {"tot", "Total inelastic"}/*,
28  {"cex", "Charge exchange"},
29  {"dcex","Double charge exchange"},
30  {"qe", "Quasi-elastic"},
31  {"abs", "Absorption"},
32  {"prod","Pion production"}
33  */
34  };
35 
36  std::vector<std::pair<std::string,std::string>> cutName_list = {
37  {"all", "All"}
38 /* {"numucc", "#nu_{#mu} CC"},
39  {"numubarcc", "#bar{#nu}_{#mu} CC"},
40  {"nuecc", "#nu_{e} CC"},
41  {"nuebarcc", "#bar{#nu}_{e} CC"},
42  {"nc", "NC"}
43  */
44  };
45 
46 
47  std::vector<std::pair<std::string,std::string>> varName_list = {
48  {"Enu", "True E_{#nu}"},
49  {"Epi", "True E_{#pi}"}
50  };
51 
52 
53 std::pair<TH1*,TH1*> RatioPlots(TH1* hSig,
54  std::vector<TH1*> univs,
55  std::string myFullLabel)
56 {
57 
58  std::pair<TH1*, TH1*> ret;
59 
60  int down_idx = 0;
61  int up_idx = 1;
62 
63  TH1* hDown = univs[down_idx];
64  hDown->Divide(hSig);
65  hDown->GetXaxis()->SetTitle(myFullLabel.c_str());
66 
67  TH1* hUp = univs[up_idx];
68  hUp->Divide(hSig);
69  hUp->GetXaxis()->SetTitle(myFullLabel.c_str());
70 
71  ret = std::make_pair(hDown, hUp);
72 
73  return ret;
74 
75 }//end RatioPlots()
76 
77 
78 
79 
80 
81 }
82 
83 
84 using namespace ana;
85 
86 void example_plot()
87 {
88 
89 
90 
91 
92  TGaxis::SetMaxDigits(5);
93  std::string inName = "outFile.root";
94 
95 
96  TFile* inFile = TFile::Open(inName.c_str(), "read");
97 
98  for(uint iVar=0; iVar<varName_list.size(); ++iVar){
99 
100 
101  std::string myVarName = varName_list[iVar].first;
102  std::string myVarLabel = varName_list[iVar].second;
103 
104 
105  for(uint iCut=0; iCut<cutName_list.size(); ++iCut){
106 
107  std::string myCutName = cutName_list[iCut].first;
108  std::string myCutLabel = cutName_list[iCut].second;
109 
110 
111  std::string histTitle = myCutLabel;
112 
113  std::string base_dir = myVarName + "_" + myCutName ;
114 
115  // Load in nominal spectrum
116  Spectrum * sSig = Spectrum::LoadFrom( inFile->GetDirectory("foo") , (base_dir + "_nom").c_str() ).release();
117 
118 
119  double pot = sSig->POT();
120  TH1 * hSig = sSig->ToTH1(pot);
121  hSig->GetXaxis()->SetTitle(myVarName.c_str());
122 
123  for(uint iSyst=0; iSyst<systName_list.size(); ++iSyst){
124 
125 
126  std::string mySystName = systName_list[iSyst].first;
127  std::string mySystLabel = systName_list[iSyst].second;
128 
129  std::string myHistTitle = histTitle + " " + mySystLabel;
130 
131  std::string outName = "./g4rwgt_" + base_dir + "_" + mySystName;
132 
133  hSig->SetTitle((myHistTitle).c_str());
134 
135  std::vector<TH1*> hUnivs{
136  (Spectrum::LoadFrom( inFile->GetDirectory("foo") , (base_dir + "_" + mySystName + "_minus1sigma").c_str() ).release() )->ToTH1(pot),
137  (Spectrum::LoadFrom( inFile->GetDirectory("foo") , (base_dir + "_" + mySystName + "_plus1sigma").c_str() ).release() )->ToTH1(pot)
138  };
139 
140  TH1* hDown = hUnivs[0];
141  TH1* hUp = hUnivs[1];
142 
143  hDown->SetLineColor(kRed);
144  hUp->SetLineColor(kBlue);
145 
146  TH1* frac_down = (TH1*)hDown->Clone();
147  frac_down->Divide(hSig);
148 
149  TH1* frac_up = (TH1*)hUp->Clone();
150  frac_up->Divide(hSig);
151 
152  auto legend = new TLegend(0.6,0.5,0.95,0.9);
153  legend->AddEntry(hSig, "Nominal", "l");
154  legend->AddEntry(hDown, "Knob val = 0.8", "l");
155  legend->AddEntry(hUp, "Knob val = 1.2", "l");
156 
157  //.........................................
158  TString rCname = "rC";
159 
160  TCanvas* rC = new TCanvas(rCname, rCname);
161 
162  rC -> SetBottomMargin(0.);
163  double Spl = 0.3;
164  TPad* P1 = new TPad( "Temp_1", "", 0.0, Spl, 1.0, 1.0, 0 );
165  TPad* P2 = new TPad( "Temp_2", "", 0.0, 0.0, 1.0, Spl, 0 );
166  P2 -> SetRightMargin (.03);
167  P2 -> SetTopMargin (.00);
168  P2 -> SetBottomMargin(.3);
169  P2 -> SetLeftMargin (.13);
170  P2 -> Draw();
171  P1 -> SetRightMargin (.03);
172  P1 -> SetLeftMargin (.13);
173  P1 -> SetTopMargin (.1);
174  P1 -> SetBottomMargin(.00);
175  P1 -> Draw();
176  // Set some label sizes.
177  double Lb1 = 0.07;
178  double Lb2 = 0.13;
179  // --- First, draw the fracUncert so cd onto Pad2
180  P2 -> cd();
181 
182  // Set axis ranges etc.
183  frac_down->GetYaxis()->SetTitleSize( Lb2 );
184  frac_down->GetYaxis()->SetTitleOffset(0.4);
185  //frac_down->GetYaxis()->SetLabelOffset(0.05);
186  frac_down->GetYaxis()->SetLabelSize( Lb2 );
187  frac_down->GetXaxis()->SetTitleSize( Lb2 );
188  frac_down->GetXaxis()->SetLabelSize( Lb2 );
189  frac_down->GetYaxis()->SetRangeUser(0.9, 1.1);
190  frac_down->GetYaxis()->SetTitle("Ratio");
191 
192  frac_down->Draw("HIST");
193 
194  double x[2] = {frac_down->GetXaxis()->GetXmin(), frac_down->GetXaxis()->GetXmax()};
195  double y[2] = {1, 1};
196  TGraph *line = new TGraph(2, x, y);
197  line->SetLineColor(kGray);
198  line->SetLineWidth(3);
199  line->SetLineStyle(2);
200  line->Draw("l same");
201 
202 
203  frac_down->Draw("SAME HIST AXIS");
204  frac_up->Draw("SAME HIST");
205 
206  P1->cd();
207 
208  hSig->GetYaxis()->SetTitleSize( Lb1 );
209  hSig->GetYaxis()->SetLabelSize( Lb1 );
210  hSig->GetYaxis()->SetTitleOffset( 0.7 );
211  // Remove the x axis labels
212  hSig->GetXaxis()->SetLabelSize (0 );
213  hSig->GetXaxis()->SetLabelOffset(99);
214 
215  //.........................................
216 
217  hSig->Draw("hist ][");
218  hDown->Draw("SAME HIST");
219  hUp->Draw("SAME HIST");
220  legend->Draw();
221 
222  rC->SaveAs((outName + ".png").c_str());
223 
224  } // end for(iSyst)
225  } // end for(iCut)
226  }//end for(iVar)
227 
228  inFile->Close();
229 
230 }
231 
Draw
tree Draw("slc.nhit")
kRed
enum BeamMode kRed
ana
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
ana::Spectrum::ToTH1
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:148
make_syst_table_plots.line
line
Definition: make_syst_table_plots.py:268
make_pair
std::pair< Spectrum *, CheatDecomp * > make_pair(SpectrumLoaderBase &loader_data, SpectrumLoaderBase &loader_mc, HistAxis *axis, Cut *cut, const SystShifts &shift, const Var &wei)
Definition: DataMCLoad.C:336
ana::cutName_list
std::vector< std::pair< std::string, std::string > > cutName_list
Definition: example_plot.C:36
make_mec_shifts_plots.legend
legend
Definition: make_mec_shifts_plots.py:148
ana::Spectrum
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
submit_syst.y
y
Definition: submit_syst.py:249
inFile
ifstream inFile
Definition: AnaPlotMaker.h:34
submit_syst.x
x
Definition: submit_syst.py:249
ana::Spectrum::LoadFrom
static std::unique_ptr< Spectrum > LoadFrom(TDirectory *dir, const std::string &name)
Definition: Spectrum.cxx:535
pot
#define pot
Definition: template_basic.C:49
ana::varName_list
std::vector< std::pair< std::string, std::string > > varName_list
Definition: example_plot.C:47
ana::Spectrum::POT
double POT() const
Definition: Spectrum.h:227
example_plot
void example_plot()
Definition: example_plot.C:86
ana::RatioPlots
std::pair< TH1 *, TH1 * > RatioPlots(TH1 *hSig, std::vector< TH1 * > univs, std::string myFullLabel)
Definition: example_plot.C:53
kBlue
enum BeamMode kBlue
cd
c cd(1)
ana::systName_list
std::vector< std::pair< std::string, std::string > > systName_list
Definition: example_plot.C:26
string
enum BeamMode string
uint
unsigned int uint
Definition: DAQChannelMapBaseOLD.h:22
Generated on Sat Mar 6 2021 14:12:57 for NOvA by   doxygen 1.8.11