Compare_Spectra.C
Go to the documentation of this file.
4 #include "CAFAna/Core/Spectrum.h"
6 #include "CAFAna/Core/Loaders.h"
7 
8 #include "CAFAna/Vars/HistAxes.h"
11 
13 
14 #include "TFile.h"
15 
16 #include <vector>
17 #include <string>
18 
19 using namespace ana;
20 
22  std::cout << "Please provide;"
23  << "An input samweb definition. \n"
24  << "Whether you want to run all cuts or not. I will assume full cuts.\n"
25  << std::endl;
26 }
27 
28 void Compare_Spectra( std::string InputDef, bool IsCosmics = false, bool RunAllCuts=false, bool Debug=false )
29 {
30 
31  // What PID and Cos Rej Cut am I using?
32  Cut kMyPIDCut = kNumuPID2018;
33  Cut kMyCosRejCut = kNumuCosmicRej2018;
34 
35  // Add the detector specific ones.
36  Cut kDetCut = kNumuContainFD2017;
37  std::string sDet = "FD_";
38  if ( InputDef.find("nd") != std::string::npos) {
39  kDetCut = kNumuContainND2017;
40  sDet = "ND_";
41  }
42 
43  // What weights do I want to use?
44  Var MyWeight = kUnweighted;
45  if ( IsCosmics){
46  MyWeight = kTimingSidebandWeight;
47  std::cout << "Cosmics --> Sideband" << std::endl;
48  }
49 
50  // --- The timing cut.
51  Cut TimingCut = kInBeamSpill;
52  if (IsCosmics) {
53  TimingCut = kInTimingSideband;
54  std::cout << "Cosmics --> SidebandTime" << std::endl;
55  }
56 
57  // ----- Define my cuts....
58  std::vector<Cut> kDefinedCut;
59  unsigned int NumCuts = 1;
60  if (RunAllCuts) NumCuts=4;
61  std::string CutNames[NumCuts];
62  if (NumCuts == 4) { // --- NumCuts == 4
63  // --- Cut
64  kDefinedCut.push_back( TimingCut && kNumuQuality ); // Good quality NuMu events
65  kDefinedCut.push_back( TimingCut && kNumuQuality && kDetCut ); // That the event is contained.
66  kDefinedCut.push_back( TimingCut && kNumuQuality && kDetCut ); // Reject cosmics
67  kDefinedCut.push_back( TimingCut && kNumuQuality && kDetCut && kMyPIDCut ); // Reject Neutral current events
68  // --- CutNames
69  CutNames[0] = "_QualCuts";
70  CutNames[1] = "_DetCuts" ;
71  CutNames[2] = "_CosmCuts";
72  CutNames[3] = "_PIDCuts";
73  // --- CosmicRej
74  if ( InputDef.find("fd") != std::string::npos ) {
75  kDefinedCut[2] = kDefinedCut[2] && kMyCosRejCut;
76  kDefinedCut[3] = kDefinedCut[3] && kMyCosRejCut;
77  }
78  }
79  else if (NumCuts == 1) { // --- NumCuts == 1
80  // --- Cut
81  kDefinedCut.push_back( TimingCut && kNumuQuality && kDetCut && kMyPIDCut ); // Go straight to full cuts.
82  // --- CutNames
83  CutNames[0] = "_AllCuts";
84  // --- CosmicRej
85  if ( InputDef.find("fd") != std::string::npos ) {
86  kDefinedCut[0] = kDefinedCut[0] && kMyCosRejCut;
87  }
88 
89  //kDefinedCut[0] = kDefinedCut[0] && kNumu2017Veto;
90  }
91 
92  // --- Load the Quantile Cuts.
93  std::string sFHC = "fhc";
94  if (InputDef.find("rhc") != std::string::npos ) sFHC = "rhc";
95  // Figure out the file name should be...
96  std::string NUMUDATA_DIR = "/cvmfs/nova.opensciencegrid.org/externals/numudata/v00.00/NULL/";
97  std::string fdspecfile = NUMUDATA_DIR+"ana2018/Quantiles/quantiles__"+sFHC+"_full__numu2018.root";
98  //std::string fdspecfile = (std::string)std::getenv("NUMUDATA_DIR")+"ana2018/Quantiles/quantiles__"+sFHC+"_full__numu2018.root";
99 
100  // Load the file...
101  std::cerr << "fdspecfile is " << fdspecfile << std::endl;
102  TFile* inFile = TFile::Open( fdspecfile.c_str() );
103  std::cerr << "Opened file." << std::endl;
104  TH2 *FDSpec2D = (TH2*)inFile->FindObjectAny( "FDSpec2D" );
105  std::cerr << "Loaded object." << std::endl;
106  // Add get the cuts....
107  std::vector<Cut> HadEFracQuantCuts = QuantileCutsFromTH2(FDSpec2D, kNumuCCOptimisedAxis, kHadEFracAxis, 4);
108 
109  std::cerr << "Loaded the quartiles." << std::endl;
110 
111  // ----- Define my quantiles.
112  const unsigned int NumQuants = 5;
113  Cut kQuantCut[NumQuants] = { kNoCut , HadEFracQuantCuts[0], HadEFracQuantCuts[1], HadEFracQuantCuts[2], HadEFracQuantCuts[3] };
114  const std::string QuantNames[NumQuants] = { "_AllQuants", "_Quant1" , "_Quant2" , "_Quant3" , "_Quant4" };
115 
116  // --- Decide which kind of sample I am looking at....
117  // I need to add the period to the output name.
118  std::string sPer = "_full_";
119  if (InputDef.find("7d8b") != std::string::npos ) sPer = "_Epochs_7d_8b_";
120  /*
121  if (InputDef.find("period1") != std::string::npos ) sPer = "_Period1_";
122  else if (InputDef.find("period2") != std::string::npos ) sPer = "_Period2_";
123  else if (InputDef.find("period3") != std::string::npos ) sPer = "_Period3_";
124  else if (InputDef.find("period4") != std::string::npos ) sPer = "_Period4_";
125  else if (InputDef.find("period5") != std::string::npos ) sPer = "_Period5_";
126  else if (InputDef.find("period6") != std::string::npos ) sPer = "_Period6_";
127  else if (InputDef.find("period7") != std::string::npos ) sPer = "_Period7_";
128  else if (InputDef.find("period8") != std::string::npos ) sPer = "_Period8_";
129  */
130  // Need to flag whether NuMI data or cosmics.
131  std::string sSwap = "NuMI_";
132  if ( IsCosmics) sSwap = "Cosmics_";
133 
134  // --- What kind of files?
135  std::string sCAF = "CAF";
136  if (InputDef.find("sumdecaf") != std::string::npos ) sCAF = "concat";
137  else if (InputDef.find("restricted") != std::string::npos ) sCAF = "res-concat";
138 
139  // --- Make the name of my output file.
140  std::string OutputName = "CompPlots_" + sDet + sSwap + sFHC + sPer + sCAF + ".root";
141 
142  // ---- Write some output.
143  std::cout << "Before I declare anything, the following is set;\n"
144  << "\tInput definit is " << InputDef << "\n"
145  << "\tOutput file name " << OutputName << "\n"
146  << "\tNumber of Cuts " << NumCuts << "\n"
147  << "\tIs debug mode is " << Debug << "\n"
148  << std::endl;
149 
150 
151  // --- Declare all of my spectra, before setting them in a loop...
152  // *** Some key ones first.
153  Spectrum *sReconstEnergy[NumCuts][NumQuants]; // --- Reconstructed energy --- Vars.h ===> kCCE ---> sr->energy.numu.E
154  Spectrum *sSliceTimeFull[NumCuts][NumQuants]; // --- Mean time of slice Full --- NumuVars.h ===> kSliceTime ---> sr->slc.meantime/1000
155  Spectrum *sSliceTimeZoom[NumCuts][NumQuants]; // --- Mean time of slice Zoom --- NumuVars.h ===> kSliceTime ---> sr->slc.meantime/1000
156  Spectrum *sMuonEnergy [NumCuts][NumQuants]; // --- Muon energy in the slice --- Vars.h ===> kMuE --> energy.numu.recomuonE
157  Spectrum *sMuonEnPerHit [NumCuts][NumQuants]; // --- Muon energy per hit --- NumuVars.h ===> kTrkEPerNHit --> sr->trk.kalman.tracks[0].calE / sr->trk.kalman.tracks[0].nhit
158  Spectrum *sHadronEnergy [NumCuts][NumQuants]; // --- Hadronic energy --- NumuVars.h ===> kHadE --> sr->energy.numu.E - sr->energy.numu.recomuonE
159  Spectrum *sHadroEnPerHit[NumCuts][NumQuants]; // --- Hadronic energy per hit --- NumuVars.h ===> kHadEPerNHit --> function.
160  Spectrum *sHadFracEnergy[NumCuts][NumQuants]; // --- Hadronic energy fraction --- NumuVars.h ===> kHadEFrac --> kHadE / kCCE
161  Spectrum *sHitsPerSlice [NumCuts][NumQuants]; // --- Number of hits per slice --- Vars.h ===> kNHit --> slc.nhit
162  Spectrum *sRemIDScore [NumCuts][NumQuants]; // --- The RemID score --- Vars.h ===> kRemID --> sr->sel.remid.pid
163  Spectrum *sCVNCosmicScor[NumCuts][NumQuants]; // --- CVN Cosmic score --- Vars.h ===> sr->sel.cvn.output[390]
164  Spectrum *sCVNNuMuIDScor[NumCuts][NumQuants]; // --- CVN NuMuID score --- Vars.h ===> sr->sel.cvnProd3Train.numuid
165  Spectrum *sCVNNuMuID2017[NumCuts][NumQuants]; // --- CVN NuMuID score --- Vars.h ===> sr->sel.cvn2017.numuid
166  Spectrum *sNuMuContPID [NumCuts][NumQuants]; // --- Third analysis BDT --- NuMuVars.h ===> sr->sel.cosrej.numucontpid
167  Spectrum *sSANuMuContPID[NumCuts][NumQuants]; // --- Second analysis BDT --- Above.
168  // *** Some important Truth based ones.
169  Spectrum *sTrueNuEnergy [NumCuts][NumQuants]; // --- True neutrino energy --- Vars.h ===> kTrueE ---> sr->mc.nu[0].E
170  Spectrum *sReTrOverTrNuE[NumCuts][NumQuants]; // --- Reco - True / True Nu E --- Above.
171 
172  // *** Basic positional plots.
173  Spectrum *sTrkStartX[NumCuts][NumQuants]; // --- Start Pos of Track --- NumuVars.h ===> kTrkStartX --> trk.kalman.tracks[0].start.X()/100
174  Spectrum *sTrkStartY[NumCuts][NumQuants]; // --- Start Pos of Track --- NumuVars.h ===> kTrkStartY --> trk.kalman.tracks[0].start.Y()/100
175  Spectrum *sTrkStartZ[NumCuts][NumQuants]; // --- Start Pos of Track --- NumuVars.h ===> kTrkStartZ --> trk.kalman.tracks[0].start.Z()/100
176  Spectrum *sTrkEndX [NumCuts][NumQuants]; // --- End Pos of Track --- NumuVars.h ===> kTrkEndX --> trk.kalman.tracks[0].stop.X()/100
177  Spectrum *sTrkEndY [NumCuts][NumQuants]; // --- End Pos of Track --- NumuVars.h ===> kTrkEndY --> trk.kalman.tracks[0].stop.Y()/100
178  Spectrum *sTrkEndZ [NumCuts][NumQuants]; // --- End Pos of Track --- NumuVars.h ===> kTrkEndZ --> trk.kalman.tracks[0].stop.Z()/100
179  Spectrum *sTrkLenXY [NumCuts][NumQuants]; // --- Diff in Track Len X Y --- Above. ===> kTrkLenXY --> (tracks[0].start.X - tracks[0].stop.X) / 100
180  //*
181  // *** Track based ones...
182  Spectrum *sNumKalTracks [NumCuts][NumQuants]; // --- Number of Kal trk in slc --- NuMuVars.h ===> sr->trk.kalman.ntracks
183  Spectrum *sNumKalTrHits [NumCuts][NumQuants]; // --- Number of hits in Kal tr --- NuMuVars.h ===> sr->trk.kalman.tracks[0].nhit
184  Spectrum *sRatKalHitSlc [NumCuts][NumQuants]; // --- Ratio of Kal hits in slc --- Above.
185  Spectrum *sKalTrLength [NumCuts][NumQuants]; // --- Length of Kalman track --- NuMuVars.h ===> sr->trk.kalman.tracks[0].len / 100
186  Spectrum *sKalTrBeamAng [NumCuts][NumQuants]; // --- Angle Kal track & beam --- NuMuVars.h ===> kCosNumi --> sr->trk.kalman.tracks[0].dir.Dot(beamDirND)
187  Spectrum *sKalMostFwdCel[NumCuts][NumQuants]; // --- Most forward Kal tr hit --- Above.
188  Spectrum *sKalMostBakCel[NumCuts][NumQuants]; // --- Most backward Kal tr hit --- Above.
189  Spectrum *sKalTrVer_MaxY[NumCuts][NumQuants]; // --- Highest Y Kalman tr vert --- Above.
190  Spectrum *sKalTrVer_MaxZ[NumCuts][NumQuants]; // --- Highest Z Kalman tr vert --- Above.
191  Spectrum *sKalTrDir_Y [NumCuts][NumQuants]; // --- Y direction of Kal track --- NuMuVars.h ===> sr->trk.kalman.tracks[0].dir.Y()
192  Spectrum *sScattKalTrLen[NumCuts][NumQuants]; // --- Scattering over tr len --- Above.
193 
194  // *** Hit based ones...
195  Spectrum *sFirstHitCell [NumCuts][NumQuants]; // --- The first cell hit --- Above.
196  Spectrum *sLastHitCell [NumCuts][NumQuants]; // --- The last cell hit --- Above.
197  Spectrum *sMaxActivity_Y[NumCuts][NumQuants]; // --- Maximum Y pos of activty --- NuMuVars.h ===> sr->slc.boxmax.Y() / 100.
198  Spectrum *sMinActivity_Y[NumCuts][NumQuants]; // --- Minimum Y pos of activty --- NuMuVars.h ===> sr->slc.boxmin.Y() / 100.
199  Spectrum *sMaxActivity_Z[NumCuts][NumQuants]; // --- Maximum Z pos of activty --- NuMuVars.h ===> sr->slc.boxmax.Z() / 100.
200  Spectrum *sMinActivity_Z[NumCuts][NumQuants]; // --- Minimum Z pos of activty --- NuMuVars.h ===> sr->slc.boxmin.Z() / 100.
201  Spectrum *sMinCellToEdge[NumCuts][NumQuants]; // --- Minimum cells to edge --- Above.
202 
203  // *** RemID score / inputs
204  Spectrum *sRemIDScatLLH [NumCuts][NumQuants]; // --- RemID scattering angle --- NuMuVars.h ===> kReMIdScatLLH --> sel.remid.scatllh
205  Spectrum *sRemIDdEdxLLH [NumCuts][NumQuants]; // --- RemID dE/dx --- NuMuVars.h ===> kReMIdDEDxLLH --> sel.remid.dedxllh
206  Spectrum *sRemIDMeasFrac[NumCuts][NumQuants]; // --- RemID fraction of planes --- NuMuVars.h ===> kReMIdMeasFrac -> sel.remid.measfrac
207 
208  // *** Some calibration quantities
209  Spectrum *sCor_RecEn[NumCuts][NumQuants][4];
210  Spectrum *sCor_MuoEn[NumCuts][NumQuants][4];
211  Spectrum *sCor_HadEn[NumCuts][NumQuants][4];
212  Spectrum *sCor_HFrEn[NumCuts][NumQuants][4];
213  // Some special truth based ones.
214  Spectrum *sCor_RecTr[NumCuts][NumQuants][4];
215  Spectrum *sCor_MuoTr[NumCuts][NumQuants][4];
216  Spectrum *sCor_HadTr[NumCuts][NumQuants][4];
217  Spectrum *sCor_HFrTr[NumCuts][NumQuants][4];
218  Spectrum *sCor_ReTr [NumCuts][NumQuants][4];
219  //*/
220  const Cut kCorner[4] = { kW_TL, kW_TR, kW_BL, kW_BR };
221  const std::string CutName[4] = { "_TopL", "_TopR", "_BotL", "_BotR" };
222 
223  // --- Set my loader
224  SpectrumLoader loader( InputDef );
226 
227  // --- And now set what is going to go into each of my spectra.
228  for (unsigned int cc=0; cc<NumCuts; ++cc) {
229  // --- Loop through my quantiles.
230  for (unsigned int qq=0; qq<NumQuants; ++qq) {
231  // What is the cut this time?
232  const Cut kThisCut = kDefinedCut[cc] && kQuantCut[qq];
233 
234  // *** Some key ones first...
235  sReconstEnergy[cc][qq] = new Spectrum( "The reconstructed energy", OptEnBins , loader, kCCE , kThisCut, kNoShift, MyWeight );
236  sSliceTimeFull[cc][qq] = new Spectrum( "The time of slice - Full", TimeFBins , loader, kSliceTime , kThisCut, kNoShift, MyWeight );
237  sSliceTimeZoom[cc][qq] = new Spectrum( "The time of slice - Zoom", TimeZBins , loader, kSliceTime , kThisCut, kNoShift, MyWeight );
238  sMuonEnergy [cc][qq] = new Spectrum( "The muon energy " , LowEnBins , loader, kMuE , kThisCut, kNoShift, MyWeight );
239  sMuonEnPerHit [cc][qq] = new Spectrum( "The muon energy per hit" , EnPHitBins, loader, kTrkEPerNHit , kThisCut, kNoShift, MyWeight );
240  sHadronEnergy [cc][qq] = new Spectrum( "The hadronic energy" , LowEnBins , loader, kHadE , kThisCut, kNoShift, MyWeight );
241  sHadroEnPerHit[cc][qq] = new Spectrum( "The hadronic En per hit ", EnPHitBins, loader, kHadEPerNHit , kThisCut, kNoShift, MyWeight );
242  sHadFracEnergy[cc][qq] = new Spectrum( "The hadronic energy frac", RatioBins , loader, kHadEFrac , kThisCut, kNoShift, MyWeight );
243  sHitsPerSlice [cc][qq] = new Spectrum( "Number of hits per slice", NumHitBins, loader, kNHit , kThisCut, kNoShift, MyWeight );
244  sRemIDScore [cc][qq] = new Spectrum( "The RemID score" , RemIDBins , loader, kRemID , kThisCut, kNoShift, MyWeight );
245  sCVNCosmicScor[cc][qq] = new Spectrum( "CVN Cosmic score" , RatioBins , loader, kCVNcos , kThisCut, kNoShift, MyWeight );
246  sCVNNuMuIDScor[cc][qq] = new Spectrum( "CVN NuMu score" , RatioBins , loader, kCVNm , kThisCut, kNoShift, MyWeight );
247  sCVNNuMuID2017[cc][qq] = new Spectrum( "2017 CVN NuMu score" , RatioBins , loader, kCVNm2017 , kThisCut, kNoShift, MyWeight );
248  sNuMuContPID [cc][qq] = new Spectrum( "The NuMu containment PID", RatioBins , loader, kNumuContPID , kThisCut, kNoShift, MyWeight );
249  sSANuMuContPID[cc][qq] = new Spectrum( "SA NuMu containment PID" , RatioBins , loader, kSANumuContPID, kThisCut, kNoShift, MyWeight );
250  // *** Some important truth based ones.
251  sTrueNuEnergy [cc][qq] = new Spectrum( "The true neutrino energy", EnergyBins, loader, kTrueE , kThisCut, kNoShift, MyWeight );
252  sReTrOverTrNuE[cc][qq] = new Spectrum( "The R-T ov T Nu energy" , RangeBins , loader, kRTOvTNuE, kThisCut, kNoShift, MyWeight );
253 
254  // *** Positional plots
255  sTrkStartX[cc][qq] = new Spectrum( "Starting X Pos of Track", YRangeBins, loader, kTrkStartX, kThisCut, kNoShift, MyWeight );
256  sTrkStartY[cc][qq] = new Spectrum( "Starting Y Pos of Track", YRangeBins, loader, kTrkStartY, kThisCut, kNoShift, MyWeight );
257  sTrkStartZ[cc][qq] = new Spectrum( "Starting Z Pos of Track", ZRangeBins, loader, kTrkStartZ, kThisCut, kNoShift, MyWeight );
258  sTrkEndX[cc][qq] = new Spectrum( "Ending X Pos of Track" , YRangeBins, loader, kTrkEndX , kThisCut, kNoShift, MyWeight );
259  sTrkEndY[cc][qq] = new Spectrum( "Ending Y Pos of Track" , YRangeBins, loader, kTrkEndY , kThisCut, kNoShift, MyWeight );
260  sTrkEndZ[cc][qq] = new Spectrum( "Ending Z Pos of Track" , ZRangeBins, loader, kTrkEndZ , kThisCut, kNoShift, MyWeight );
261  sTrkLenXY[cc][qq] = new Spectrum( "Difference in XY Tr Len", ZRangeBins, loader, kTrkLenXY , kThisCut, kNoShift, MyWeight );
262  //*
263  // *** Track based ones...
264  sNumKalTracks [cc][qq] = new Spectrum( "Number of Kalman Tracks" , NumTrkBins, loader, kNKalman , kThisCut, kNoShift, MyWeight );
265  sNumKalTrHits [cc][qq] = new Spectrum( "Number of Kal Track hits", NumHitBins, loader, kTrkNhits , kThisCut, kNoShift, MyWeight );
266  sRatKalHitSlc [cc][qq] = new Spectrum( "Ratio of KalHit in Slice", RatioBins , loader, kRatOfKalHi, kThisCut, kNoShift, MyWeight );
267  sKalTrLength [cc][qq] = new Spectrum( "Length of Kalman Track" , LengthBins, loader, kTrkLength , kThisCut, kNoShift, MyWeight );
268  sKalTrBeamAng [cc][qq] = new Spectrum( "Angle of KalTr to Beam" , RatioBins , loader, kCosNumi , kThisCut, kNoShift, MyWeight );
269  sKalMostFwdCel[cc][qq] = new Spectrum( "Most forward KalTr Cell" , CellsBins , loader, kKalFwdCell, kThisCut, kNoShift, MyWeight );
270  sKalMostBakCel[cc][qq] = new Spectrum( "Most backward KalTr Cell", CellsBins , loader, kKalBakCell, kThisCut, kNoShift, MyWeight );
271  sKalTrVer_MaxY[cc][qq] = new Spectrum( "Most pos Y KalTr Vertex" , YRangeBins, loader, kMaxKalYPos, kThisCut, kNoShift, MyWeight );
272  sKalTrVer_MaxZ[cc][qq] = new Spectrum( "Most pos Z KalTr Vertex" , ZRangeBins, loader, kMaxKalZPos, kThisCut, kNoShift, MyWeight );
273  sKalTrDir_Y [cc][qq] = new Spectrum( "Y Direction of KalTr" , CompBins , loader, kDirY , kThisCut, kNoShift, MyWeight );
274  sScattKalTrLen[cc][qq] = new Spectrum( "Scatt over KalTr length" , TinyBins , loader, kScattTrLen, kThisCut, kNoShift, MyWeight );
275 
276  // *** Hit based ones...
277  sFirstHitCell [cc][qq] = new Spectrum( "The first cell hit in sl", CellsBins , loader, kFirstCell , kThisCut, kNoShift, MyWeight );
278  sLastHitCell [cc][qq] = new Spectrum( "The last cell hit in sl" , CellsBins , loader, kLastCell , kThisCut, kNoShift, MyWeight );
279  sMaxActivity_Y[cc][qq] = new Spectrum( "Largest Y loc with a hit", YRangeBins, loader, kSlcMaxY , kThisCut, kNoShift, MyWeight );
280  sMinActivity_Y[cc][qq] = new Spectrum( "Smalest Y loc with a hit", YRangeBins, loader, kSlcMinY , kThisCut, kNoShift, MyWeight );
281  sMaxActivity_Z[cc][qq] = new Spectrum( "Largest Z loc with a hit", ZRangeBins, loader, kSlcMaxZ , kThisCut, kNoShift, MyWeight );
282  sMinActivity_Z[cc][qq] = new Spectrum( "Largest Z loc with a hit", ZRangeBins, loader, kSlcMinZ , kThisCut, kNoShift, MyWeight );
283  sMinCellToEdge[cc][qq] = new Spectrum( "Min cells hit to edge" , CellsBins , loader, kMinCellEdg, kThisCut, kNoShift, MyWeight );
284 
285  // *** RemID score / inputs
286  sRemIDScatLLH [cc][qq] = new Spectrum( "RemID Scattering angle" , CompZBins , loader, kReMIdScatLLH , kThisCut, kNoShift, MyWeight );
287  sRemIDdEdxLLH [cc][qq] = new Spectrum( "RemID dE/dx" , CompBins , loader, kReMIdDEDxLLH , kThisCut, kNoShift, MyWeight );
288  sRemIDMeasFrac[cc][qq] = new Spectrum( "RemID fraction of planes", RatioBins , loader, kReMIdMeasFrac , kThisCut, kNoShift, MyWeight );
289 
290  // *** Calibration quantities
291  for (int cor=0; cor<4; ++cor) {
292  const Cut kCorCut = kThisCut && kCorner[cor];
293  sCor_RecEn[cc][qq][cor] = new Spectrum( "Reconstructed Energy with corner cut", EnergyBins, loader, kCCE , kCorCut, kNoShift, MyWeight);
294  sCor_MuoEn[cc][qq][cor] = new Spectrum( "Muon energy with corner cut" , EnergyBins, loader, kMuE , kCorCut, kNoShift, MyWeight);
295  sCor_HadEn[cc][qq][cor] = new Spectrum( "Hadronic energy with corner cut" , EnergyBins, loader, kHadE , kCorCut, kNoShift, MyWeight);
296  sCor_HFrEn[cc][qq][cor] = new Spectrum( "Hadronic energy frac with corner cut", RatioBins , loader, kHadEFrac, kCorCut, kNoShift, MyWeight);
297  // Additional truth based ones.
298  sCor_RecTr[cc][qq][cor] = new Spectrum( "True Energy with corner cut" , EnergyBins, loader, kTrueE , kCorCut, kNoShift, MyWeight);
299  sCor_MuoTr[cc][qq][cor] = new Spectrum( "True Mu En with corner cut" , EnergyBins, loader, kTrueMuonE, kCorCut, kNoShift, MyWeight);
300  sCor_HadTr[cc][qq][cor] = new Spectrum( "True Had En with corner cut" , EnergyBins, loader, kTrueHadE , kCorCut, kNoShift, MyWeight);
301  sCor_HFrTr[cc][qq][cor] = new Spectrum( "True HadEFrac with corner cut", RatioBins , loader, kTrueHFrE , kCorCut, kNoShift, MyWeight);
302  sCor_ReTr [cc][qq][cor] = new Spectrum( "Reco-True/True w. corner cut" , RangeBins , loader, kRTOvTNuE , kCorCut, kNoShift, MyWeight);
303  }
304  //*/
305  } // Quantile loop
306  } // Cut loop
307 
308  // --- Do it!
309  loader.Go();
310 
311  // --- Open my outfile so that I save histograms!
312  TFile *outFile = new TFile(OutputName.c_str(),"RECREATE");
313  std::cerr << "Made my output file...." << std::endl;
314 
315  // --- Finally, write them to disk!
316  for (unsigned int cc=0; cc<NumCuts; ++cc) {
317  // --- Loop through my quantiles.
318  for (unsigned int qq=0; qq<NumQuants; ++qq) {
319  std::string MyAppend = CutNames[cc] + QuantNames[qq];
320  // *** Some key ones first...
321  sReconstEnergy[cc][qq] -> SaveTo( outFile, TString("sReconstEnergy")+TString(MyAppend) ) ;
322  sSliceTimeFull[cc][qq] -> SaveTo( outFile, TString("sSliceTimeFull")+TString(MyAppend) ) ;
323  sSliceTimeZoom[cc][qq] -> SaveTo( outFile, TString("sSliceTimeZoom")+TString(MyAppend) ) ;
324  sMuonEnergy [cc][qq] -> SaveTo( outFile, TString("sMuonEnergy") +TString(MyAppend) ) ;
325  sMuonEnPerHit [cc][qq] -> SaveTo( outFile, TString("sMuonEnPerHit") +TString(MyAppend) ) ;
326  sHadronEnergy [cc][qq] -> SaveTo( outFile, TString("sHadronEnergy") +TString(MyAppend) ) ;
327  sHadroEnPerHit[cc][qq] -> SaveTo( outFile, TString("sHadroEnPerHit")+TString(MyAppend) ) ;
328  sHadFracEnergy[cc][qq] -> SaveTo( outFile, TString("sHadFracEnergy")+TString(MyAppend) ) ;
329  sHitsPerSlice [cc][qq] -> SaveTo( outFile, TString("sHitsPerSlice") +TString(MyAppend) ) ;
330  sRemIDScore [cc][qq] -> SaveTo( outFile, TString("sRemIDScore") +TString(MyAppend) ) ;
331  sCVNCosmicScor[cc][qq] -> SaveTo( outFile, TString("sCVNCosmicScor")+TString(MyAppend) ) ;
332  sCVNNuMuIDScor[cc][qq] -> SaveTo( outFile, TString("sCVNNuMuIDScor")+TString(MyAppend) ) ;
333  sCVNNuMuID2017[cc][qq] -> SaveTo( outFile, TString("sCVNNuMuID2017")+TString(MyAppend) ) ;
334  sNuMuContPID [cc][qq] -> SaveTo( outFile, TString("sNuMuContPID") +TString(MyAppend) ) ;
335  sSANuMuContPID[cc][qq] -> SaveTo( outFile, TString("sSANuMuContPID")+TString(MyAppend) ) ;
336  // *** Some important truth based ones.
337  sTrueNuEnergy [cc][qq] -> SaveTo( outFile, TString("sTrueNuEnergy") +TString(MyAppend) ) ;
338  sReTrOverTrNuE[cc][qq] -> SaveTo( outFile, TString("sReTrOverTrNuE")+TString(MyAppend) ) ;
339 
340  // *** Positional plots
341  sTrkStartX[cc][qq] -> SaveTo( outFile, TString("sTrkStartX")+TString(MyAppend) ) ;
342  sTrkStartY[cc][qq] -> SaveTo( outFile, TString("sTrkStartY")+TString(MyAppend) ) ;
343  sTrkStartZ[cc][qq] -> SaveTo( outFile, TString("sTrkStartZ")+TString(MyAppend) ) ;
344  sTrkEndX[cc][qq] -> SaveTo( outFile, TString("sTrkEndX") +TString(MyAppend) ) ;
345  sTrkEndY[cc][qq] -> SaveTo( outFile, TString("sTrkEndY") +TString(MyAppend) ) ;
346  sTrkEndZ[cc][qq] -> SaveTo( outFile, TString("sTrkEndZ") +TString(MyAppend) ) ;
347  sTrkLenXY[cc][qq] -> SaveTo( outFile, TString("sTrkLenXY") +TString(MyAppend) ) ;
348  //*
349  // *** Track based ones...
350  sNumKalTracks [cc][qq] -> SaveTo( outFile, TString("sNumKalTracks") +TString(MyAppend) ) ;
351  sNumKalTrHits [cc][qq] -> SaveTo( outFile, TString("sNumKalTrHits") +TString(MyAppend) ) ;
352  sRatKalHitSlc [cc][qq] -> SaveTo( outFile, TString("sRatKalHitSlc") +TString(MyAppend) ) ;
353  sKalTrLength [cc][qq] -> SaveTo( outFile, TString("sKalTrLength") +TString(MyAppend) ) ;
354  sKalTrBeamAng [cc][qq] -> SaveTo( outFile, TString("sKalTrBeamAng") +TString(MyAppend) ) ;
355  sKalMostFwdCel[cc][qq] -> SaveTo( outFile, TString("sKalMostFwdCel")+TString(MyAppend) ) ;
356  sKalMostBakCel[cc][qq] -> SaveTo( outFile, TString("sKalMostBakCel")+TString(MyAppend) ) ;
357  sKalTrVer_MaxY[cc][qq] -> SaveTo( outFile, TString("sKalTrVer_MaxY")+TString(MyAppend) ) ;
358  sKalTrVer_MaxZ[cc][qq] -> SaveTo( outFile, TString("sKalTrVer_MaxZ")+TString(MyAppend) ) ;
359  sKalTrDir_Y [cc][qq] -> SaveTo( outFile, TString("sKalTrDir_Y") +TString(MyAppend) ) ;
360  sScattKalTrLen[cc][qq] -> SaveTo( outFile, TString("sScattKalTrLen")+TString(MyAppend) ) ;
361 
362  // *** Hit based ones...
363  sFirstHitCell [cc][qq] -> SaveTo( outFile, TString("sFirstHitCell") +TString(MyAppend) ) ;
364  sLastHitCell [cc][qq] -> SaveTo( outFile, TString("sLastHitCell") +TString(MyAppend) ) ;
365  sMaxActivity_Y[cc][qq] -> SaveTo( outFile, TString("sMaxActivity_Y")+TString(MyAppend) ) ;
366  sMinActivity_Y[cc][qq] -> SaveTo( outFile, TString("sMinActivity_Y")+TString(MyAppend) ) ;
367  sMinCellToEdge[cc][qq] -> SaveTo( outFile, TString("sMinCellToEdge")+TString(MyAppend) ) ;
368  sMaxActivity_Z[cc][qq] -> SaveTo( outFile, TString("sMaxActivity_Z")+TString(MyAppend) ) ;
369  sMinActivity_Z[cc][qq] -> SaveTo( outFile, TString("sMinActivity_Z")+TString(MyAppend) ) ;
370 
371  // *** RemID score / inputs
372  sRemIDScatLLH [cc][qq] -> SaveTo( outFile, TString("sRemIDScatLLH") +TString(MyAppend) ) ;
373  sRemIDdEdxLLH [cc][qq] -> SaveTo( outFile, TString("sRemIDdEdxLLH") +TString(MyAppend) ) ;
374  sRemIDMeasFrac[cc][qq] -> SaveTo( outFile, TString("sRemIDMeasFrac")+TString(MyAppend) ) ;
375 
376  // *** Calibration quantities
377  for (int cor=0; cor<4; ++cor) {
378  sCor_RecEn[cc][qq][cor] -> SaveTo( outFile, TString("sCor_RecEn")+TString(MyAppend)+TString(CutName[cor]) ) ;
379  sCor_MuoEn[cc][qq][cor] -> SaveTo( outFile, TString("sCor_MuoEn")+TString(MyAppend)+TString(CutName[cor]) ) ;
380  sCor_HadEn[cc][qq][cor] -> SaveTo( outFile, TString("sCor_HadEn")+TString(MyAppend)+TString(CutName[cor]) ) ;
381  sCor_HFrEn[cc][qq][cor] -> SaveTo( outFile, TString("sCor_HFrEn")+TString(MyAppend)+TString(CutName[cor]) ) ;
382  // Additional truth based ones.
383  sCor_RecTr[cc][qq][cor] -> SaveTo( outFile, TString("sCor_RecTr")+TString(MyAppend)+TString(CutName[cor]) ) ;
384  sCor_MuoTr[cc][qq][cor] -> SaveTo( outFile, TString("sCor_MuoTr")+TString(MyAppend)+TString(CutName[cor]) ) ;
385  sCor_HadTr[cc][qq][cor] -> SaveTo( outFile, TString("sCor_HadTr")+TString(MyAppend)+TString(CutName[cor]) ) ;
386  sCor_HFrTr[cc][qq][cor] -> SaveTo( outFile, TString("sCor_HFrTr")+TString(MyAppend)+TString(CutName[cor]) ) ;
387  sCor_ReTr [cc][qq][cor] -> SaveTo( outFile, TString("sCor_ReTr" )+TString(MyAppend)+TString(CutName[cor]) ) ;
388  }
389  //*/
390  std::cerr << "\tWritten quantile set " << qq << " for cut set " << cc << " to disk!!!" << std::endl;
391  }
392  std::cerr << "\n\tWritten cut set " << cc << " to disk!!!" << std::endl;
393  }
394 
395 }
const Var kHadE
Definition: NumuVars.h:23
const Binning ZRangeBins
const Var kNKalman
Definition: NumuVars.cxx:540
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Var kMaxKalYPos([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return std::max(sr->trk.kalman.tracks[0].start.Y()/100, sr->trk.kalman.tracks[0].stop.Y()/100);})
const Var kReMIdScatLLH
Definition: NumuVars.cxx:555
void SaveTo(const osc::IOscCalc &x, TDirectory *dir, const std::string &name)
void Compare_Spectra(std::string InputDef, bool IsCosmics=false, bool RunAllCuts=false, bool Debug=false)
const Var kSlcMaxZ([](const caf::SRProxy *sr){return sr->slc.boxmax.Z()/100.;})
Definition: NumuVars.h:93
const Var kSlcMaxY([](const caf::SRProxy *sr){return sr->slc.boxmax.Y()/100.;})
Definition: NumuVars.h:92
const Var kNumuContPID
Definition: NumuVars.cxx:553
const Cut kNumuContainND2017([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return false;for(unsigned int i=0;i< sr->vtx.elastic.fuzzyk.nshwlid;++i){TVector3 start=sr->vtx.elastic.fuzzyk.png[i].shwlid.start;TVector3 stop=sr->vtx.elastic.fuzzyk.png[i].shwlid.stop;if(std::min(start.X(), stop.X())< -180.0) return false;if(std::max(start.X(), stop.X()) > 180.0) return false;if(std::min(start.Y(), stop.Y())< -180.0) return false;if(std::max(start.Y(), stop.Y()) > 180.0) return false;if(std::min(start.Z(), stop.Z())< 20.0) return false;if(std::max(start.Z(), stop.Z()) > 1525.0) return false;}if(sr->trk.kalman.ntracks< 1) return false;for(unsigned int i=0;i< sr->trk.kalman.ntracks;++i){if(i==sr->trk.kalman.idxremid) continue;else if(sr->trk.kalman.tracks[i].start.Z() > 1275||sr->trk.kalman.tracks[i].stop.Z() > 1275) return false;}return(sr->trk.kalman.ntracks > sr->trk.kalman.idxremid &&sr->slc.firstplane > 1 &&sr->slc.lastplane< 212 &&sr->trk.kalman.tracks[0].start.Z()< 1100 &&(sr->trk.kalman.tracks[0].stop.Z()< 1275 ||sr->sel.contain.kalyposattrans< 55) &&sr->sel.contain.kalfwdcellnd > 5 &&sr->sel.contain.kalbakcellnd > 10);})
Definition: NumuCuts2017.h:11
const Binning YRangeBins
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 Binning EnPHitBins
const Var kMinCellEdg([](const caf::SRProxy *sr){return std::min((sr->sel.contain.kalfwdcell+sr->sel.contain.kalbakcell), (sr->sel.contain.cosfwdcell+sr->sel.contain.cosbakcell));})
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;})
Definition: NumuVars.h:52
const Var kFirstCell([](const caf::SRProxy *sr){return sr->slc.firstcell;})
const Var kSliceTime([](const caf::SRProxy *sr){return sr->slc.meantime/1000;})
Definition: NumuVars.h:34
OStream cerr
Definition: OStream.cxx:7
const Binning LowEnBins
const Var kTrueHFrE([](const caf::SRProxy *sr){float HadE=kTrueHadE(sr);float TruE=kTrueE(sr);return(HadE/TruE);})
const Var kLastCell([](const caf::SRProxy *sr){return sr->slc.lastcell;})
const Binning TimeZBins
void SetSpillCut(const SpillCut &cut)
const Var kTrkEPerNHit([](const caf::SRProxy *sr){return(sr->trk.kalman.tracks[0].calE/sr->trk.kalman.tracks[0].nhit);})
const Var kReMIdMeasFrac
Definition: NumuVars.cxx:557
const Var kKalFwdCell([](const caf::SRProxy *sr){return sr->sel.contain.kalfwdcell;})
const Var kTrkLength([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-5.f;return sr->trk.kalman.tracks[0].len/100;})
Definition: NumuVars.h:65
const Var kScattTrLen([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return(float)((float) sr->sel.cosrej.scatt/(float) sr->trk.kalman.tracks[0].len);})
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;})
Definition: NumuVars.h:63
const Var kDirY([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-5.f;return sr->trk.kalman.tracks[0].dir.Y();})
Definition: NumuVars.h:38
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
const Var kSlcMinY([](const caf::SRProxy *sr){return sr->slc.boxmin.Y()/100.;})
Definition: NumuVars.h:88
const Cut kInTimingSideband([](const caf::SRProxy *sr){if(sr->spill.run > util::kLastBadTimingRun) return(kInTimingSideband_before(sr)|| kInTimingSideband_after(sr));else return(kInTimingSideband_before(sr)|| kInTimingSideband_afterA(sr)|| kInTimingSideband_afterB(sr));}, [](const caf::SRSpillProxy *spill){if(spill->run > util::kLastBadTimingRun) return(kInTimingSideband_before.Livetime(spill)+ kInTimingSideband_after.Livetime(spill));else return(kInTimingSideband_before.Livetime(spill)+ kInTimingSideband_afterA.Livetime(spill)+ kInTimingSideband_afterB.Livetime(spill));}, [](const caf::SRSpillProxy *spill){return 0;})
Definition: TimingCuts.h:12
const Cut kNumuCosmicRej2018([](const caf::SRProxy *sr){return(sr->sel.cosrej.anglekal > 0.5 && sr->sel.cosrej.numucontpid2019 > 0.53 && sr->slc.nhit< 400 && sr->sel.nuecosrej.pngptp< 0.9 );})
Definition: NumuCuts2018.h:19
const Binning RemIDBins
ifstream inFile
Definition: AnaPlotMaker.h:34
std::string sFHC
Definition: MakeCutFlow.C:35
const Var kTrueE([](const caf::SRProxy *sr){assert(sr->mc.nnu==1);return sr->mc.nu[0].E;})
Definition: Vars.cxx:85
const Var kRatOfKalHi([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return(float)((float) sr->trk.kalman.tracks[0].nhit/(float) sr->slc.nhit);})
const Binning EnergyBins
TFile * outFile
Definition: PlotXSec.C:135
const Cut kNumuContainFD2017
Definition: NumuCuts2017.h:21
const Var kRemID
PID
Definition: Vars.cxx:81
const Cut kW_TR([](const caf::SRProxy *sr){bool IsTrue(sr->trk.kalman.tracks[0].start.X() > XCent && sr->trk.kalman.tracks[0].start.Y() > YCent);return IsTrue;})
const Var kSANumuContPID([](const caf::SRProxy *sr){return sr->sel.cosrej.numuSAcontpid;})
const Var kHadEFrac
Definition: NumuVars.h:24
const Var kNHit
Definition: Vars.cxx:71
const Cut kW_BR([](const caf::SRProxy *sr){bool IsTrue(sr->trk.kalman.tracks[0].start.X() > XCent && sr->trk.kalman.tracks[0].start.Y()< YCent);return IsTrue;})
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;})
Definition: NumuVars.h:53
const Var kCCE
Definition: NumuVars.h:21
std::string sPer
Definition: MakeThePlots.C:111
virtual void Go() override
Load all the registered spectra.
const Binning CompBins
const std::string QuantNames[NumQuant]
Definition: DoThePlots.C:207
loader
Definition: demo0.py:10
const Binning CompZBins
std::vector< float > Spectrum
Definition: Constants.h:610
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 Cut kW_TL([](const caf::SRProxy *sr){bool IsTrue=(sr->trk.kalman.tracks[0].start.X()< XCent && sr->trk.kalman.tracks[0].start.Y() > YCent);return IsTrue;})
const Binning CellsBins
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;})
Definition: NumuVars.h:51
const SystShifts kNoShift
Definition: SystShifts.cxx:21
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;})
Definition: NumuVars.h:57
OStream cout
Definition: OStream.cxx:6
const Var kTrkLenXY([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;float XLen=(sr->trk.kalman.tracks[0].start.X()-sr->trk.kalman.tracks[0].stop.X())/100;float YLen=(sr->trk.kalman.tracks[0].start.Y()-sr->trk.kalman.tracks[0].stop.Y())/100;return XLen-YLen;})
const Var kKalBakCell([](const caf::SRProxy *sr){return sr->sel.contain.kalbakcell;})
const Cut kW_BL([](const caf::SRProxy *sr){bool IsTrue(sr->trk.kalman.tracks[0].start.X()< XCent && sr->trk.kalman.tracks[0].start.Y()< YCent);return IsTrue;})
const Var kRTOvTNuE
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;})
Definition: NumuVars.h:56
const Binning NumHitBins
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;})
Definition: NumuVars.h:55
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
const Binning LengthBins
void cc()
Definition: test_ana.C:28
const Var kTrueMuonE([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return 0.f;if(sr->mc.nu[0].prim.empty()) return 0.f;if(std::abs(sr->mc.nu[0].prim[0].pdg)!=13) return 0.f;return float(sr->mc.nu[0].prim[0].p.E);})
Definition: NumuVars.h:107
const Cut kNumuPID2018([](const caf::SRProxy *sr){std::cout<< "ERROR::kNumuPID2018, cutting on both cvnProd3Train and cvn2017."<< " Neither branch exists anymore. Returning False."<< std::endl;abort();return false;})
Definition: NumuCuts2018.h:22
const Var kTrueHadE
const Binning OptEnBins
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 Binning RangeBins
const Binning NumTrkBins
const Var kReMIdDEDxLLH
Definition: NumuVars.cxx:556
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
const Binning TimeFBins
const Var kMaxKalZPos([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return-10.0f;return std::max(sr->trk.kalman.tracks[0].start.Z()/100, sr->trk.kalman.tracks[0].stop.Z()/100);})
const Cut kNumuQuality
Definition: NumuCuts.h:18
const Var kMuE
Definition: NumuVars.h:22
void SystematicComp()
const Var kCVNm
PID
Definition: Vars.cxx:39
const Var kUnweighted
The simplest possible Var, always 1. Used as a default weight.
Definition: Var.h:96
const Cut kNoCut
The simplest possible cut: pass everything, used as a default.
Definition: Cut.h:109
const Var kCosNumi([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks > 0 &&sr->trk.kalman.idxremid!=999){if(sr->hdr.det==1){return sr->trk.kalman.tracks[0].dir.Dot(beamDirND);}if(sr->hdr.det==2){return sr->trk.kalman.tracks[0].dir.Dot(beamDirFD);}}return-5.f;})
Definition: NumuVars.h:43
const Var kSlcMinZ([](const caf::SRProxy *sr){return sr->slc.boxmin.Z()/100.;})
Definition: NumuVars.h:89
const Var kTrkNhits([](const caf::SRProxy *sr){if(sr->trk.kalman.ntracks< 1) return 65535;return int(sr->trk.kalman.tracks[0].nhit);})
Definition: NumuVars.h:59
const Binning TinyBins
const Binning RatioBins
std::vector< std::string > CutNames
Definition: MakeCutFlow.C:49
enum BeamMode string