get_eventlist2020.C
Go to the documentation of this file.
1 #pragma once
2 
4 #include "CAFAna/Core/Loaders.h"
5 #include "CAFAna/Core/Spectrum.h"
9 #include "CAFAna/Core/Cut.h"
10 
19 #include "NuXAna/Vars/NusVars.h"
20 #include "NuXAna/Cuts/NusCuts.h"
21 #include "NuXAna/Cuts/NusCuts20.h"
24 #include "CAFAna/Vars/XsecTunes.h"
25 #include "CAFAna/Cuts/Cuts.h"
26 #include "CAFAna/Cuts/SpillCuts.h"
27 
28 using namespace ana;
29 
31  std::string outlist="OUTLIST",
32  TString sel="SELECT",
33  std::string det="DETECTOR")
34 {
35 
36  std::cout
37  << "\nSAM definition: " << sam
38  << "\nOutput file : " << outlist
39  << "\nSelection : " << sel
40  << "\nDetector : " << det
41  << std::endl;
42 
43  // get quantile stuff
44  std::string infile_quant = "/cvmfs/nova.opensciencegrid.org/externals/numudata/v00.05/NULL/lib/ana2020/Quantiles/quantiles_fhc_full_numu2020.root";
45  TFile* infile = new TFile(infile_quant.c_str());
46  TH2* FDSpec2D = (TH2*)infile->FindObjectAny("FDSpec2D");
47  std::vector<Cut> HadEFracQuantCuts = QuantileCutsFromTH2(FDSpec2D, kNumuCCOptimisedAxis, kHadEFracAxis, 4);
48  infile->Close();
49  delete infile;
50 
51  const ana::Var kRunNumber = SIMPLEVAR(hdr.run);
52  const ana::Var kRun = SIMPLEVAR(hdr.run);
53  const ana::Var kSubRun = SIMPLEVAR(hdr.subrun);
54  const ana::Var kEvt = SIMPLEVAR(hdr.evt);
55  const ana::Var kCycle = SIMPLEVAR(hdr.cycle);
56  const ana::Var kSlice = SIMPLEVAR(hdr.subevt);
57  const ana::Var kSliceTime = SIMPLEVAR(hdr.unixtime);
58  const ana::Var kSliceMeanTime = SIMPLEVAR(hdr.subevtmeantime);
59  const ana::Var kSliceTsd = SIMPLEVAR(slc.tsd);
60  const ana::Var kAlwaysOnWgt = ana::kXSecCVWgt2020 * ana::kPPFXFluxCVWgt;
61 
62  const ana::Var kCCNC([](const caf::SRProxy* sr)
63  {
64  int nnu = sr->mc.nu.size();
65  if (nnu == 0) return -5;
66  else {
67  int ccnc = sr->mc.nu[0].iscc;
68  return ccnc;
69  }
70  });
71 
72  const ana::Var kRecoEnergy([sel](const caf::SRProxy* sr)
73  {
74  double en = -5;
75  if (sel == "numu")
76  en = ana::kCCE(sr);
77  if (sel == "nue")
78  en = ana::kNueEnergy2020(sr);
79  if (sel == "nc")
80  en = ana::kNus20Energy(sr);
81  return en;
82  });
83 
84  //
85  const auto kSpillCuts
87 
88  ana::Cut kCuts = ana::kNoCut;
89  if (sel == "numuq1")
90  kCuts = kCuts && HadEFracQuantCuts.at(0);
91  if (sel == "numuq2")
92  kCuts = kCuts && HadEFracQuantCuts.at(1);
93  if (sel == "numuq3")
94  kCuts = kCuts && HadEFracQuantCuts.at(2);
95  if (sel == "numuq4")
96  kCuts = kCuts && HadEFracQuantCuts.at(3);
97 
98  if (sel.Contains("numu")){
99  if (sam.find("genie") != std::string::npos){
100  if (det == "fd")
101  kCuts = kCuts && ana::kNumu2020FD_ML && ana::kInBeamSpill;
102  else
103  kCuts = kCuts && ana::kNumu2020ND;
104  }
105  else
107  }
108  else if (sel == "nue"){
109  if (sam.find("genie") != std::string::npos){
110  if (det == "fd"){
111  std::cout << __LINE__ << std::endl;
113  }
114  else
115  kCuts = ana::kNue2020ND;
116  }
117  else
119  }
120  else if (sel == "nc"){
121  if (sam.find("genie") != std::string::npos){
122  if (det == "fd")
124  else
125  kCuts = ana::kNus20NDCuts;
126  }
127  else
129  }
130  else throw std::logic_error("selection must be either numu, nue or nc");
131 
132  std::cout << sam
133  << " "
134  << outlist
135  << std::endl;
136 
137  std::vector<const Var*> intVars = {&kRun, &kSubRun, &kEvt, &kCycle, &kSlice};
138  std::vector<const Var*> floatVars;
139 
140  floatVars.push_back(&kEME_2020);
141  floatVars.push_back(&kHADE_2020);
142  floatVars.push_back(&kCaloE);
143  floatVars.push_back(&kRecoEnergy);
144  if (sam.find("genie") != std::string::npos){
145  floatVars.push_back(&kTrueE);
146  floatVars.push_back(&kTruePDG);
147  floatVars.push_back(&kCCNC);
148  floatVars.push_back(&kAlwaysOnWgt);
149  }
150 
151  MakeTextListFile(sam, {kCuts}, {outlist}, floatVars, intVars, &kSpillCuts);
152 }
153 
const Var kTruePDG
Definition: TruthVars.h:48
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
const Cut kCCNC([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;return(sr->mc.nu[0].iscc)?false:true;})
Proxy for caf::StandardRecord.
Definition: SRProxy.h:2126
const Var kNus20Energy([](const caf::SRProxy *sr){if(sr->hdr.det!=caf::kNEARDET &&sr->hdr.det!=caf::kFARDET) return(double) sr->slc.calE;double pars[4][6]={{1.049, 0.795, 0.8409, 0.17, 0.82,-1.00},{1.025, 0.797, 0.9162, 0.53,-0.26,-1.48},{1.000, 1.000, 1.0000, 0.00, 0.00, 0.00},{1.000, 1.000, 1.0000, 0.00, 0.00, 0.00}};int detCur=(sr->hdr.det==caf::kFARDET)+((sr->spill.isRHC)<< 1);double e_EM=ana::kEME_2020(sr);double e_Had=ana::kHADE_2020(sr);double e_Cal=sr->slc.calE;return(e_EM/pars[detCur][0]+e_Had/pars[detCur][1])/(pars[detCur][2]+pars[detCur][3]*std::pow(e_Cal+pars[detCur][4], 2)*std::exp(pars[detCur][5]*e_Cal));})
Definition: NusVars.h:64
caf::Proxy< std::vector< caf::SRNeutrino > > nu
Definition: SRProxy.h:618
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
const Cut kNumu2020ND
Definition: NumuCuts2020.h:57
const Var kSliceTime([](const caf::SRProxy *sr){return sr->slc.meantime/1000;})
Definition: NumuVars.h:34
void get_eventlist2020(std::string sam="SAMLIST", std::string outlist="OUTLIST", TString sel="SELECT", std::string det="DETECTOR")
const Cut kNue2020ND
Definition: NueCuts2020.h:178
const Cut kNus20NDCuts
Definition: NusCuts20.h:102
const Cut kInCosmicTimingWindow
Is the event far from the start and ends of the spill ? For FD cosmic selection.
Definition: TimingCuts.cxx:165
void MakeTextListFile(const std::string &wildcard, const std::vector< Cut > &cut, const std::vector< std::string > &output, const std::vector< const Var * > &floatVars, const std::vector< const Var * > &intVars, const SpillCut *spillCut)
Make a file listing all the events passing the specified cut.
Definition: EventList.cxx:193
const Var kTrueE([](const caf::SRProxy *sr){assert(sr->mc.nnu==1);return sr->mc.nu[0].E;})
Definition: Vars.cxx:85
const Var kEvt
Definition: Vars.cxx:23
string infile
const Var kCaloE
Summed calorimetric energy of all hits in slice, uncorrected.
Definition: Vars.cxx:52
const Var kCCE
Definition: NumuVars.h:21
caf::StandardRecord * sr
const Var kSubRun
const Var kCycle
Definition: Vars.cxx:22
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
const Var kHADE_2020([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;const double CVNha_CalE=sr->slc.calE;return std::max(CVNha_CalE-kEME_2020(sr), 0.);})
Definition: NueEnergy2020.h:15
OStream cout
Definition: OStream.cxx:6
caf::Proxy< caf::SRTruthBranch > mc
Definition: SRProxy.h:2138
const Var kNueEnergy2020([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;if(kIsRHC(sr)) return kNueEnergyRHC2020(sr);else return kNueEnergyFHC2020(sr);})
Nue energy with 3d prong info. only (old version of vars)
Definition: NueEnergy2020.h:38
#define SIMPLEVAR(CAFNAME)
For Vars where literally all you need is a single CAF variable.
Definition: Var.h:88
std::vector< Cut > QuantileCutsFromTH2(TH2 *quantileHist, const HistAxis &independentAxis, const HistAxis &quantileAxis, const unsigned int &numQuantiles, const bool verbose)
: Do the same as the QuantileCuts function but taking in the TH2 instead of making it...
const Cut kInBeamSpill([](const caf::SRProxy *sr){if(sr->spill.run > util::kLastBadTimingRun) return kInBeamSpill_main(sr);else return kInBeamSpill_main(sr)||kInBeamSpill_shifted(sr);}, [](const caf::SRSpillProxy *spill){if(spill->run > util::kLastBadTimingRun) return kInBeamSpill_main.Livetime(spill);else return kInBeamSpill_main.Livetime(spill)+kInBeamSpill_shifted.Livetime(spill);}, [](const caf::SRSpillProxy *spill) -> double{return spill->spillpot;})
Does the event fall inside the window we call the beam spill?
Definition: TimingCuts.h:8
const Var kSlice
const Cut kNus20FDCuts
Definition: NusCuts20.h:174
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
const Cut kNue2020FDAllSamples
Definition: NueCuts2020.h:84
const Var kRun
Definition: Vars.cxx:20
const Cut kNumu2020FD_ML
const Cut kNoCut
The simplest possible cut: pass everything, used as a default.
Definition: Cut.h:109
const Var kEME_2020([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;double CVNem_CalE=0.0;const caf::SRFuzzyKProngProxy &prim_png=sr->vtx.elastic.fuzzyk.png[0];for(const caf::SRFuzzyKProngProxy &png:sr->vtx.elastic.fuzzyk.png){double png_CalE=png.shwlid.calE;double emPID=((double) png.cvnpart.photonid+ (double) png.cvnpart.pizeroid+ (double) png.cvnpart.electronid);double haPID=((double) png.cvnpart.protonid+ (double) png.cvnpart.pionid+ (double) png.cvnpart.neutronid+ (double) png.cvnpart.otherid+ (double) png.cvnpart.muonid);if(emPID< 0) continue;if(emPID >=haPID ) CVNem_CalE+=png_CalE;else continue;}if(kLongestProng(sr) >=500) CVNem_CalE=prim_png.shwlid.calE;return CVNem_CalE;})
Definition: NueEnergy2020.h:14
const Var kXSecCVWgt2020
Definition: XsecTunes.h:106
enum BeamMode string