12 #include "CAFAna/Core/Binning.h" 18 #include "CAFAna/Cuts/NumuCuts.h" 23 #include "CAFAna/Vars/NumuVars.h" 46 TLatex*
prelim =
new TLatex(.9, .95,
"NO#nuA Preliminary");
47 prelim->SetTextColor(
kBlue);
49 prelim->SetTextSize(2/30.);
50 prelim->SetTextAlign(32);
77 if(
sr->
mc.
nnu == 0)
return false;
102 std::cout <<
"\nrun : ---- Running CAF NuMU ND Validation.\n";
111 TH1D* spillHist =
new TH1D(
"reco-spills",
";Detector;Spills", 3, 0, 3);
120 std::vector<Selection> selections_base;
122 if(isFD) selections_base.emplace_back(
kNumuFD,
"numuFD",
123 "Selected events pass the numu FD selection");
124 else selections_base.emplace_back(
kNumuND,
"numuND",
125 "Selected events pass the numu ND selection");
129 const Var kRecoMuE({
"energy.mutrkE.E",
"sel.remid.bestidx",
"energy.mutrkE"},
134 if(
sr->
energy.mutrkE.size() == 0)
return 0.
f;
156 std::vector<TangibleAxis> variables;
157 variables.emplace_back(kNumuCCShiftAxis,
"reco-numuEShift",
158 "Shifted CC energy estimator (GeV) ");
160 "CC energy estimator (GeV) ");
161 variables.emplace_back(axRecoMuE,
"reco-muE",
162 "Muon energy (GeV)");
164 "reco-slcNHit",
"Number of hits in slice. " );
165 variables.emplace_back(
HistAxis(
"Slice Calorimetric Energy (GeV)", kEnergyBinning,
kCaloE),
166 "reco-calE",
"Calorimetric energy of slice. " );
167 variables.emplace_back(
HistAxis(
"Hadronic Energy (GeV)", kHadronicEnergyBinning,
kHadE),
169 "Hadronic enery, i.e. numu energy estimate minus muon track energy. " );
170 variables.emplace_back(
HistAxis(
"Hadronic Energy (GeV)", kHadronicEnergyBinning, kHadEshift),
172 "Hadronic enery, i.e. numu energy estimate minus muon track energy. " );
173 variables.emplace_back(
HistAxis(
"Average Hadronic Energy Per Hit (GeV)",
176 "Average energy per hit in hadronic cluster, i.e. had_E/had_n_hit. " );
180 "Average energy per hit on primary track, i.e. trk_E/trk_n_hit. " );
182 "reco-hadNHit",
"Number of hits in hadronic cluster. ");
183 variables.emplace_back(
HistAxis(
"Calorimetric Hadronic Energy (GeV)",
185 "reco-hadcalE",
"Sum of calibrated energy deposits in hadronic cluster. " );
186 variables.emplace_back(
HistAxis(
"Slice Maximum Y (m)", kXYBins,
kMaxY),
187 "reco-maxy",
"Maximum Y position of slice. " );
189 "reco-nkal",
"Number of tracks in slice. " );
191 "reco-nhit",
"Number of hits in slice. " );
193 "reco-dirZ",
"Z-direction of muon track. " );
195 "reco-trkStartX",
"Track start x position. " );
197 "reco-trkStartY",
"Track start y position. " );
199 "reco-trkStartZ",
"Track start z position. " );
200 variables.emplace_back(
HistAxis(
"Track End X Position (m)", kXYBins,
kTrkEndX),
201 "reco-trkEndX",
"Track stop x position. " );
202 variables.emplace_back(
HistAxis(
"Track End Y Position (m)", kXYBins,
kTrkEndY),
203 "reco-trkEndY",
"Track stop y position. " );
204 variables.emplace_back(
HistAxis(
"Track End Z Position (m)", kZBins,
kTrkEndZ),
205 "reco-trkEndZ",
"Track stop z position. " );
207 "reco-sliceDuration",
"Slice duration. " );
209 "reco-trkNhits",
"Number of hits on primary track. ");
211 "reco-trkLength",
"Primary track length. " );
213 "reco-scatLL",
"ReMId scattering log log-likelihood for primary track. " );
215 "reco-dedxLL",
"ReMId dE/dx log log-likelihood for primary track. " );
216 variables.emplace_back(
HistAxis(
"Non-hadronic Plane Fraction",
218 "reco-nonHadPlaneFrac",
"ReMId Non-hadronic plane fraction. " );
220 "reco-remid",
"ReMId kNN score. " );
223 std::vector<UsefulHist>
hists;
224 hists.reserve(selections_base.size() * variables.size());
226 for(
const auto& sel_base:selections_base){
227 for(
const auto& variable:variables){
228 hists.emplace_back(sel_base, variable, loader);
233 std::cout <<
"\nrun : --- run loader.\n";
236 double pot = hists[0].fHist.POT();
238 std::cout <<
"\nrun : --- save output.\n";
240 TFile outputfile(kOutputFileName.c_str(),
"RECREATE");
242 for(
const auto&
hist:hists){
245 temp->Write(tempName.c_str());
250 TH1F *pothist =
new TH1F(
"meta-TotalPOT",
"TotalPOT", 1, 0, 1);
251 TH1F *evthist =
new TH1F(
"meta-TotalEvents",
"TotalEvents", 1, 0, 1);
252 pothist->SetBinContent(1, pot);
253 evthist->SetBinContent(1, spillHist->GetEntries());
254 pothist->Write(
"meta-TotalPOT");
255 evthist->Write(
"meta-TotalEvents");
257 std::cout <<
"\nrun : --- close output files.\n";
const Var kHadNHit([](const caf::SRProxy *sr){unsigned int nought=0;if(sr->trk.kalman.ntracks< 1) return nought;return sr->slc.nhit-sr->trk.kalman.tracks[0].nhit;})
_HistAxis< Var > HistAxis
Represent the binning of a Spectrum's x-axis.
Cuts and Vars for the 2020 FD DiF Study.
UsefulHist(Selection selA, TangibleAxis tanAxis, SpectrumLoader &loader, const Cut &bkg=kNoCut)
const Binning kRemidBinning
Binning for plotting remid attractively.
const Var kTrkStartY([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return sr->trk.kalman.tracks[0].start.Y()/100;})
Tangible< Cut > Selection
const Var kSliceDuration([](const caf::SRProxy *sr){return(sr->slc.endtime-sr->slc.starttime);})
void SetSpillCut(const SpillCut &cut)
const Var kDirZ([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-5.f;return sr->trk.kalman.tracks[0].dir.Z();})
const Var kTrkLength([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-5.f;return sr->trk.kalman.tracks[0].len/100;})
const Var kHadEPerNHit([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return 0.0f;int nHit=sr->slc.nhit-sr->trk.kalman.tracks[0].nhit;if(nHit<=0) return 0.0f;float hadE=sr->energy.numu.hadcalE;return hadE/nHit;})
Representation of a spectrum in any variable, with associated POT.
void Preliminary()
Put NOvA Preliminary on plots.
const Var kMaxY([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1000.0f;if(sr->vtx.elastic.fuzzyk.nshwlid==0) return-1000.0f;float maxyall=-999.0;for(unsigned int i=0;i< sr->vtx.elastic.fuzzyk.nshwlid;i++){maxyall=std::max(std::max(sr->vtx.elastic.fuzzyk.png[i].shwlid.start.Y(), sr->vtx.elastic.fuzzyk.png[i].shwlid.stop.Y()), maxyall);}return maxyall;})
const HistAxis kNumuCCAxis("Reconstructed Neutrino Energy (GeV)", kNumuEnergyBinning, kCCE)
virtual void AddSpillHistogram(TH1 *h, const SpillVar &var, const SpillCut &cut, const SpillVar &wei=kSpillUnweighted)
Uses include counting the total POT or spills in a run.
const Var kCaloE
Summed calorimetric energy of all hits in slice, uncorrected.
Tangible(const T &obj, const std::string &shortName, const std::string &blurb)
const Binning kHadronicEnergyBinning
const Var kTrkStartZ([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return sr->trk.kalman.tracks[0].start.Z()/100;})
SRRemid remid
Output from RecoMuonID (ReMId) package.
Tangible< HistAxis > TangibleAxis
const Binning kEnergyBinning
virtual void Go() override
Load all the registered spectra.
short nnu
Number of neutrinos in nu vector (0 or 1)
const Var kTrkStartX([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return sr->trk.kalman.tracks[0].start.X()/100;})
const Var kTrkEndZ([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return sr->trk.kalman.tracks[0].stop.Z()/100;})
The StandardRecord is the primary top-level object in the Common Analysis File trees.
const Var kTrkEndY([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return sr->trk.kalman.tracks[0].stop.Y()/100;})
const Var kTrkEndX([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return sr->trk.kalman.tracks[0].stop.X()/100;})
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
SRIDBranch sel
Selector (PID) branch.
assert(nhit_max >=nhit_nbins)
string shortName
THUMBNAIL BLOCK: We need to make a thumbnail for each.
SRTruthBranch mc
Truth branch for MC: energy, flavor, etc.
const SpillCut kStandardSpillCuts
Apply this unless you're doing something special.
void caf_numu_vars(std::string kInputFileName, std::string kOutputFileName, bool isFD)
Template for Var and SpillVar.
SREnergyBranch energy
Energy estimator branch.
static Binning Simple(int n, double lo, double hi, const std::vector< std::string > &labels={})
std::vector< SRNeutrino > nu
implemented as a vector to maintain mc.nu structure, i.e. not a pointer, but with 0 or 1 entries...
const Cut kNoCut
The simplest possible cut: pass everything, used as a default.
const Var kTrkNhits([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return 65535;return int(sr->trk.kalman.tracks[0].nhit);})