NuebarEnergyEstimatorHelper.h
Go to the documentation of this file.
2 #include "TH2D.h"
3 #include "TLegend.h"
4 
5 namespace ana
6 {
8  public:
9  InteractionSpectra(TFile* fin, bool _is_shifted,
10  bool _is_flat, std::string _fit_name)
11  : is_shifted(_is_shifted),
12  is_flat(_is_flat),
13  fit_name(_fit_name)
14  {
15  spec_name = "res_%s_";
17  if(is_shifted) spec_name += "_shifted";
18  spec_name += is_flat ? "_flat" : "_shape";
19 
20  spectra.push_back(Spectrum::LoadFrom(fin, TString::Format(spec_name.c_str(),
21  interaction_chns[0].name.c_str())).release());
22  POT = spectra[0]->POT();
23  hists.push_back(spectra.back()->ToTH1(POT));
24  for(int id = 1; id < kcNumIntChns; id++) {
25  int i = id;
26  if(i==1) i++; // hack
27  spectra.push_back(Spectrum::LoadFrom(fin, TString::Format(spec_name.c_str(),
28  interaction_chns[i].name.c_str())).release());
29  hists.push_back(spectra.back()->ToTH1(POT));
30  }
31  for(int id = 0; id < kcNumIntChns; id++) {
32  int i = id;
33  if(i==1) i++; //hack
34  hists[i]->SetLineColor(line_colors[i]);
35  }
36 
37  }
38  void MakePlot(std::string dir, TFile* fout)
39  {
40  gStyle->SetOptStat(0);
41  TDirectory* tmp = gDirectory;
42  TCanvas* c = new TCanvas();
43  c->cd();
44  std::string title = "Fractional Resolution by Interaction ";
45  title += is_shifted ? "(shift/" : "(noshift/";
46  title += is_flat ? "flat)" : "shape)";
47  hists[0]->SetTitle(title.c_str());
48  hists[0]->GetXaxis()->SetTitle("(Reco - True) / True #nu Energy");
49  hists[0]->Draw("hist");
50  for(int i = 1; i < kcNumIntChns; i++) {
51  // already cut on cc earlier. Don't plot it
52  if(interaction_chns[i].name == "cc") continue;
53  hists[i]->Draw("hist same");
54  }
55  hists[0]->Draw("hist same");
56 
57  TGraph* g = new TGraph();
58  g->SetPoint(0, 0, -10);
59  g->SetPoint(1, 0, 1e6);
60  g->SetLineColor(kRed);
61  g->SetLineStyle(2);
62  g->Draw("same");
63 
64  DrawLegend();
65  plot_name = "resolution_by_interaction_" + fit_name;
66  plot_name += is_shifted ? "_shift" : "_noshift";
67  plot_name += is_flat ? "_flat" : "_shape";
68  plot_name += ".pdf";
69  fout->cd();
70  c->SaveAs((dir + plot_name).c_str());
71  c->Write((dir + plot_name).c_str());
72 
73  delete c;
74  tmp->cd();
75  }
76 
77  void DrawLegend()
78  {
79  TLegend* leg = new TLegend();
80  for(int i = 0; i < kcNumIntChns; i++) {
81  // already cut on cc earlier. Don't plot it
82  if(interaction_chns[i].name == "cc") continue;
83  leg->AddEntry(hists[i], interaction_chns[i].name.c_str(), "l");
84  }
85  leg->Draw("same");
86  }
87 
88  std::vector<Spectrum*> spectra;
89  std::vector<TH1D*> hists;
90  std::vector<Int_t> line_colors = {kBlack, kGreen, kAzure,
91  kViolet, kCyan, kRed,
92  kOrange+8, kViolet+1};
93  double POT;
94  bool is_shifted;
95  bool is_flat;
99  };
100 
102  {
103  public:
104  ResolutionScan(TFile* fin, bool _is_shifted,
105  bool _is_flat, std::string _fit_name,
106  std::string _xvar)
107  : is_shifted(_is_shifted),
108  is_flat(_is_flat),
109  fit_name(_fit_name),
110  xvar(_xvar)
111  {
112  spec_name = "res_scan_" + xvar + fit_name;
113  spec_name += (is_shifted) ? "_shifted" : "";
114  spec_name += (is_flat) ? "_flat" : "_shape";
115  spec = Spectrum::LoadFrom(fin, spec_name).release();
116  POT = spec->POT();
117  hist = (TH2D*) spec->ToTH2(POT);
118  for(int i = 1; i < hist->GetNbinsX()+1; i++) {
119  cols.push_back(hist->ProjectionY(UniqueName().c_str(), i, i));
120  }
121  }
122 
123  void MakePlot(std::string dir, TFile* fout)
124  {
125  gStyle->SetPalette(kBlueYellow);
126  TDirectory* tmp = gDirectory;
127  std::string title = "Resolution Scan " + xvar + "(";
128  title += is_shifted ? "shifted/" : "noshift/";
129  title += is_flat ? "flat)" : "shape)";
130  hist->SetTitle(title.c_str());
131  hist->GetXaxis()->SetTitle(xvar.c_str());
132  hist->GetYaxis()->SetTitle("(Reco - True) / True #nu Energy");
133 
134  TGraph* g = new TGraph();
135  g->SetPoint(0, -1, 0);
136  g->SetPoint(1, 11, 0);
137  g->SetLineColor(kRed);
138  g->SetLineStyle(7);
139  g->SetLineWidth(3);
140 
141  TH1D* hMean = MeanScan();
142  hMean->SetLineColor(kCyan);
143 
144  TCanvas* c = new TCanvas();
145  c->cd();
146  hist->Draw("colz");
147  hMean->Draw("same el");
148  g->Draw("same");
149  fout->cd();
150  c->SaveAs((dir + spec_name + ".pdf").c_str());
151  c->Write((dir + spec_name).c_str());
152  delete c;
153  tmp->cd();
154  }
155 
156  TH1D* MeanScan()
157  {
158  TH1D* ret = new TH1D(UniqueName().c_str(),
159  "",
160  hist->GetNbinsX(),
161  hist->GetXaxis()->GetBinLowEdge(1),
162  hist->GetXaxis()->GetBinUpEdge(hist->GetNbinsX()));
163  for(int i = 1; i < ret->GetNbinsX()+1; i++) {
164  double mean = cols[i-1]->GetMean();
165  double rms = cols[i-1]->GetRMS();
166  int N = cols[i-1]->GetEntries();
167  if(mean != 0) {
168  ret->SetBinContent(i, mean);
169  ret->SetBinError(i, rms/sqrt(N));
170  }
171  }
172  return ret;
173  }
174 
175  TH1D* RMSScan()
176  {
177  TH1D* ret = new TH1D(UniqueName().c_str(),
178  "",
179  hist->GetNbinsX(),
180  hist->GetXaxis()->GetBinLowEdge(1),
181  hist->GetXaxis()->GetBinUpEdge(hist->GetNbinsX()));
182  for(int i = 1; i < ret->GetNbinsX()+1; i++) {
183  double rms = cols[i-1]->GetRMS();
184  ret->SetBinContent(i, rms);
185  ret->SetBinError(i, 0);
186  }
187  return ret;
188  }
189 
190  std::string GetXVar() {return xvar;}
191  bool IsShifted() {return is_shifted;}
193 
195  TH2D* hist;
196  double POT;
198  bool is_flat;
202  std::vector<TH1D*> cols;
203  };
204 }
::xsd::cxx::tree::id< char, ncname > id
Definition: Database.h:165
TString fin
Definition: Style.C:24
const XML_Char * name
Definition: expat.h:151
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
T sqrt(T number)
Definition: d0nt_math.hpp:156
const int kcNumIntChns
Definition: NueCCIncCuts.h:398
Float_t tmp
Definition: plot.C:36
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:33
static std::unique_ptr< Spectrum > LoadFrom(TDirectory *dir, const std::string &name)
Definition: Spectrum.cxx:607
const int cols[3]
const SelDef interaction_chns[kNumIntChns]
ResolutionScan(TFile *fin, bool _is_shifted, bool _is_flat, std::string _fit_name, std::string _xvar)
InteractionSpectra(TFile *fin, bool _is_shifted, bool _is_flat, std::string _fit_name)
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
void MakePlot(std::string dir, TFile *fout)
TDirectory * dir
Definition: macro.C:5
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
void MakePlot(std::string dir, TFile *fout)
std::vector< Spectrum * > spectra
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:30