signal_count.C
Go to the documentation of this file.
1 #include "CAFAna/Core/Spectrum.h"
7 
9 
10 #include "TCanvas.h"
11 #include "TStyle.h"
12 #include "TLegend.h"
13 #include "TLegendEntry.h"
18 // temporary measure to use prod4 trained muonid
19 #include "NDAna/muonid/NDXSecMuonPID.h"
20 namespace nuebarccinc {
21  const ana::Cut kMuonIDProd4Cut = ana::kMuonIDProd4 < -0.2;
22 }
23 
25 {
26  const int ret = gSystem->Load("libNDAnamuonid");
27  if(ret != 0) exit(1);
28 }
29 
30 using namespace ana;
31 
32 void signal_count(bool make_plots=false,
33  std::string input_file_name="/nova/ana/users/ddoyle/NuebarResolutionStudy/signal_count.root",
34  std::string plot_dump="/nova/ana/users/ddoyle/NuebarResolutionStudy/signal_count_plots/")
35 {
36  if(!make_plots) {
38 
41 
42  const Var kCVWeight = kPPFXFluxCVWgt * kXSecCVWgt2020;
43 
44  const Cut kSelection =
47  // nuebarccinc::kRemid;
48 
49  Spectrum nue(lNDMC,
51  kSelection && nuebarccinc::kNueCC,
52  kNoShift,
53  kCVWeight);
54 
55 
56  Spectrum nuebar(lNDMC,
58  kSelection && nuebarccinc::kNuebarCC,
59  kNoShift,
60  kCVWeight);
61 
62  lNDMC.Go();
63 
64  TFile output("signal_count.root", "recreate");
65  nue.SaveTo(output.mkdir("nue"));
66  nuebar.SaveTo(output.mkdir("nuebar"));
67  output.Close();
68  }
69  else {
70  TFile * input = TFile::Open(input_file_name.c_str());
71  TCanvas * c = new TCanvas();
72 
73  Spectrum * nue1D = Spectrum::LoadFrom(input, "nue" )).release(;
74  Spectrum * nuebar1D = Spectrum::LoadFrom(input, "nuebar")).release(;
75  Spectrum * nue2D = new Spectrum(nue1D->ToTH1(kAna2020RHCPOT),
80  kAna2020RHCPOT, 0);
81  Spectrum * nuebar2D = new Spectrum(nuebar1D->ToTH1(kAna2020RHCPOT),
86  kAna2020RHCPOT, 0);
87 
88  // count plots
89  gStyle->SetPalette(kDeepSea);
90  gStyle->SetPaintTextFormat(".2f");
91  TH2 * nue = nue2D ->ToTH2(kAna2020RHCPOT);
92  nue->SetTitle("Signal #nu_{e}");
93  nue->SetMarkerColor(kWhite);
94  nue->GetXaxis()->SetRangeUser(0.5, 1);
95  nue->GetYaxis()->SetRangeUser(0, 6);
96  nue->Draw("colz text");
97  c->Print((plot_dump + "/nue.pdf").c_str());
98  TH2 * nuebar = nuebar2D->ToTH2(kAna2020RHCPOT);
99  nuebar->SetTitle("Signal #bar{#nu}_{e}");
100  nuebar->SetMarkerColor(kWhite);
101  nuebar->GetXaxis()->SetRangeUser(0.5, 1);
102  nuebar->GetYaxis()->SetRangeUser(0, 6);
103  nuebar->Draw("colz text");
104  c->Print((plot_dump + "/nuebar.pdf").c_str());
105 
106  TH2 * nunubar = (TH2*) nue->Clone("nunubar");
107  nunubar->Add(nuebar);
108  nunubar->SetTitle("Signal #bar{#nu}_{e} + #nu_{e}");
109  nunubar->SetMarkerColor(kWhite);
110  nunubar->GetXaxis()->SetRangeUser(0.5, 1);
111  nunubar->GetYaxis()->SetRangeUser(0, 6);
112  nunubar->Draw("colz text");
113  c->Print((plot_dump + "/nunubar.pdf").c_str());
114 
115 
116  // fraction plots
117  auto kNuebarColor = kGreen;
118  auto kNueColor = kMagenta;
119 
120  gStyle->SetPaintTextFormat("4.2g");
121  TH2 * nue_frac = (TH2*) nue ->Clone("nue_frac" );
122  TH2 * nuebar_frac = (TH2*) nuebar->Clone("nuebar_frac");
123  nue_frac->SetMarkerColor(kNueColor);
124  nuebar_frac->SetMarkerColor(kNuebarColor);
125 
126  TLegend * leg = new TLegend(0.2, 0.8, 0.3, 0.7, "", "NDC");
127 
128  auto entry = leg->AddEntry((TObject*)0, "#bar{#nu}_{e} Fraction", "");
129  entry->SetTextColor(kNuebarColor);
130  entry = leg->AddEntry((TObject*)0, "#nu_{e} Fraction" , "");
131  entry->SetTextColor(kNueColor);
132 
133  nue_frac->Divide(nunubar);
134  nuebar_frac->Divide(nunubar);
135 
136  nunubar ->Draw("colz");
137  // nuebar on top
138  nuebar_frac->SetBarOffset(0.2);
139  nue_frac ->SetBarOffset(-0.2);
140  nue_frac ->Draw("text same");
141  nuebar_frac->Draw("text same");
142  leg->Draw();
143  c->Print((plot_dump + "/nunubar_fracs.pdf").c_str());
144 
145 
146 
147  // uncertainty plots
148  TH2 * nuebar_uncert = (TH2*) nuebar->Clone("nuebar_uncert");
149  TH2 * nunubar_uncert = (TH2*) nunubar->Clone("nunubar_uncert");
150  auto NX = nuebar_uncert->GetNbinsX();
151  auto NY = nuebar_uncert->GetNbinsY();
152  for(auto x = 1; x <= NX; x++) {
153  for(auto y = 1; y <= NY; y++) {
154  auto nnuebar = nuebar->GetBinContent(x, y);
155  auto nnunubar = nunubar->GetBinContent(x, y);
156 
157  if(nnuebar)
158  nuebar_uncert->SetBinContent(x, y, 1/std::sqrt(nnuebar));
159  else
160  nuebar_uncert->SetBinContent(x, y, 0);
161 
162  if(nnunubar)
163  nunubar_uncert->SetBinContent(x, y, 1/std::sqrt(nnunubar));
164  else
165  nunubar_uncert->SetBinContent(x, y, 0);
166  }
167  }
168 
169  TH2 * uncert_diff = (TH2*) nunubar_uncert->Clone("uncert_diff");
170  uncert_diff->Add(nuebar_uncert, -1);
171 
172  nunubar_uncert->SetTitle("Signal #nu_{e} + #bar{#nu}_{e} Stat. Uncert.");
173  nuebar_uncert->SetTitle("Signal #bar{#nu}_{e} Stat. Uncert.");
174  uncert_diff->SetTitle("#Delta Stat. Uncert. (#bar{#nu} + #nu) - #bar{#nu}");
175 
176  nunubar_uncert->Draw("colz text");
177  c->Print((plot_dump + "/nunubar_uncert.pdf").c_str());
178 
179  nuebar_uncert->Draw("colz text");
180  c->Print((plot_dump + "/nuebar_uncert.pdf").c_str());
181 
182  uncert_diff->Draw("colz text");
183  c->Print((plot_dump + "/uncert_diff.pdf").c_str());
184 
185 
186  }
187 }
188 
ofstream output
TH2 * ToTH2(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Spectrum must be 2D to obtain TH2.
Definition: Spectrum.cxx:165
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
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
const ana::Cut kNuebarCC
const ana::HistAxis kTrueElectronEStandardAxis
T sqrt(T number)
Definition: d0nt_math.hpp:156
void load_libs_muonid()
Definition: signal_count.C:24
void signal_count(bool make_plots=false, std::string input_file_name="/nova/ana/users/ddoyle/NuebarResolutionStudy/signal_count.root", std::string plot_dump="/nova/ana/users/ddoyle/NuebarResolutionStudy/signal_count_plots/")
Definition: signal_count.C:32
void SetSpillCut(const SpillCut &cut)
const ana::HistAxis kTrueElectronEVsCosStandardAxis
std::unique_ptr< T > LoadFrom(TDirectory *dir, const std::string &label)
Definition: LoadFromFile.h:17
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
const ana::HistAxis kTrueElectronCosThetaStandardAxis
const std::string PROD5_MC_RHC_NOMINAL
virtual void Go() override
Load all the registered spectra.
void SaveTo(TDirectory *dir, const std::string &name) const
Definition: Spectrum.cxx:517
std::vector< float > Spectrum
Definition: Constants.h:570
const ana::Cut kPreselectionLoose
const double kAna2020RHCPOT
Definition: Exposures.h:235
const SystShifts kNoShift
Definition: SystShifts.cxx:21
const Cut kSelection
Definition: SINCpi0_Cuts.h:328
const std::vector< Binning > & GetBinnings() const
Definition: LabelsAndBins.h:69
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
exit(0)
const ana::Cut kMuonIDProd4Cut
Definition: datamc.C:28
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
const ana::Cut kNueCC
def entry(str)
Definition: HTMLTools.py:26
void make_plots(TFile *f, TH2 *hVsRun, TGraph *g, int run, std::string suffix, double xpos, TH2 *&cut)
const std::vector< std::string > & GetLabels() const
Definition: LabelsAndBins.h:68
const Var kXSecCVWgt2020
Definition: XsecTunes.h:106