NueVarsExtra.cxx
Go to the documentation of this file.
2 
7 
8 #include "CAFAna/Vars/Vars.h"
10 
11 #include "CAFAna/Vars/XsecTunes.h"
12 
14 
16 
17 #include "Utilities/func/MathUtil.h"
18 
19 #include <iostream>
20 
21 namespace ana
22 {
23 
24 #define ELASTICVAR(VAR) \
25  Var([](const caf::SRProxy* sr) \
26  { \
27  if(!sr->vtx.elastic.IsValid) return -1000.; \
28  return double(sr->vtx.elastic.vtx.VAR); \
29  })
30 
31 const Var kVtxX = ELASTICVAR(x);
32 const Var kVtxY = ELASTICVAR(y);
33 const Var kVtxZ = ELASTICVAR(z);
34 
35 const Var kEPerHit([](const caf::SRProxy* sr)
36  {
37  if(sr->slc.nhit>0) return 1000.0*(sr->slc.calE/sr->slc.nhit);
38  else return -5.;
39  });
40 
41 const Var kRecoEPerHit([](const caf::SRProxy* sr)
42  {
43  if(sr->slc.nhit>0) return 1000.0*(kNueEnergy2018(sr)/sr->slc.nhit);
44  else return -5.;
45  });
46 
47 const Var kemCalE([](const caf::SRProxy* sr)
48  {
49  if(!sr->vtx.elastic.IsValid) return -1.0;
50  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
51 
52  double CVNem_CalE = 0;
53  for(const caf::SRFuzzyKProngProxy& png: sr->vtx.elastic.fuzzyk.png){
54  double png_CalE = png.calE;
55 
56  double emPID = ( (double)png.cvnpart.photonid +
57  (double)png.cvnpart.pizeroid +
58  (double)png.cvnpart.electronid);
59  double haPID = ( (double)png.cvnpart.protonid +
60  (double)png.cvnpart.pionid +
61  (double)png.cvnpart.neutronid +
62  (double)png.cvnpart.otherid +
63  (double)png.cvnpart.muonid );
64 
65  if ( emPID <= 0 ) continue; // for no cvn scores
66  if ( emPID >= haPID ) CVNem_CalE += png_CalE;
67  } // png
68 
69  return CVNem_CalE * CalibrationBugCorrectionFactor(sr->hdr);
70  });
71 
72 const Var kEMNhit([](const caf::SRProxy* sr)
73  {
74  if(!sr->vtx.elastic.IsValid) return -1.0;
75  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
76 
77  double Nhit = 0;
78  for(const caf::SRFuzzyKProngProxy& png: sr->vtx.elastic.fuzzyk.png){
79  double hits = png.nhit;
80 
81  double emPID = ( (double)png.cvnpart.photonid +
82  (double)png.cvnpart.pizeroid +
83  (double)png.cvnpart.electronid);
84  double haPID = ( (double)png.cvnpart.protonid +
85  (double)png.cvnpart.pionid +
86  (double)png.cvnpart.neutronid +
87  (double)png.cvnpart.otherid +
88  (double)png.cvnpart.muonid );
89 
90  if ( emPID <= 0 ) continue; // for no cvn scores
91  if ( emPID >= haPID ) Nhit += hits;
92  } // png
93 
94  return Nhit;
95  });
96 
97 const Var kEMEPerHit([](const caf::SRProxy* sr)
98  {
99  if(kEMNhit(sr)>0) return 1000.0*(kemCalE(sr)/kEMNhit(sr));
100  else return -5.;
101  });
102 
103 const Var kHADEPerHit([](const caf::SRProxy* sr)
104  {
105  if((sr->slc.nhit - kEMNhit(sr))>0) return 1000.0*(sr->slc.calE*CalibrationBugCorrectionFactor(sr->hdr) - kemCalE(sr))/(sr->slc.nhit - kEMNhit(sr));
106  else return -5.;
107  });
108 
109 
110 #define SHWLIDVAR(VAR) \
111  Var( \
112  [](const caf::SRProxy* sr) \
113  { \
114  if(!sr->vtx.elastic.IsValid) return -1000.; \
115  if(sr->vtx.elastic.fuzzyk.nshwlid < 1) return -1000.; \
116  return double(sr->vtx.elastic.fuzzyk.png[0].shwlid.VAR); \
117  })
118 
122 
123 const Var kShwStopX = SHWLIDVAR(stop.x);
124 const Var kShwStopY = SHWLIDVAR(stop.y);
125 const Var kShwStopZ = SHWLIDVAR(stop.z);
126 
127 const Var kShwCalE = SHWLIDVAR(calE);
129 
130 const Var kShwWidth = SHWLIDVAR(width);
132 const Var kShwGap = SHWLIDVAR(gap);
133 
134 const Var kInelasticity([](const caf::SRProxy* sr)
135  {
136  if(!sr->vtx.elastic.IsValid) return 1.f;
137  if(sr->vtx.elastic.fuzzyk.nshwlid < 1) return 1.f;
138  return ((sr->slc.calE -
139  sr->vtx.elastic.fuzzyk.png[0].shwlid.calE)/(sr->slc.calE));
140  });
141 
142 const Var kHadCalE([](const caf::SRProxy* sr)
143  {
144  if(!sr->vtx.elastic.IsValid) return float(sr->slc.calE);
145  if(sr->vtx.elastic.fuzzyk.nshwlid < 1) return float(sr->slc.calE);
146  return ((sr->slc.calE -
147  sr->vtx.elastic.fuzzyk.png[0].shwlid.calE));
148  });
149 
150 const Var kemE2017([](const caf::SRProxy* sr)
151  {
152  if(!sr->vtx.elastic.IsValid) return -1.0;
153  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
154  return NueRecoE_2017FDFit(kCVNemE(sr), 0.0);
155  });
156 
157 const Var khadE2017([](const caf::SRProxy* sr)
158  {
159  if(!sr->vtx.elastic.IsValid) return -1.0;
160  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
161  return NueRecoE_2017FDFit(0.0, kCVNhadE(sr));
162  });
163 
164 const Var kemE2018([](const caf::SRProxy* sr)
165  {
166  if(!sr->vtx.elastic.IsValid) return -1.0;
167  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
168  if(kIsRHC(sr)) return NueRecoE_2018RHCFit(kRecoEME(sr), 0.0);
169  else return NueRecoE_2017FDFit(kRecoEME(sr), 0.0);
170  });
171 
172 const Var khadE2018([](const caf::SRProxy* sr)
173  {
174  if(!sr->vtx.elastic.IsValid) return -1.0;
175  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
176  if(kIsRHC(sr)) return NueRecoE_2018RHCFit(0.0, kRecoHADE(sr));
177  else return NueRecoE_2017FDFit(0.0, kRecoHADE(sr));
178  });
179 
180 const Var kNueHadEFrac([](const caf::SRProxy* sr)
181  {
182  if(!sr->vtx.elastic.IsValid) return -1.0;
183  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
184  if (kNueEnergy2018(sr)>0){
185  if(kIsRHC(sr)) return (NueRecoE_2018RHCFit(0.0, kRecoHADE(sr)))/(NueRecoE_2018RHCFit(kRecoEME(sr), kRecoHADE(sr)));
186  else return (NueRecoE_2017FDFit(0.0, kRecoHADE(sr)))/NueRecoE_2017FDFit(kRecoEME(sr), kRecoHADE(sr));
187  }
188  else return -5.;
189  });
190 
191 const Var kNueHadEFrac2([](const caf::SRProxy* sr)
192  {
193  if(!sr->vtx.elastic.IsValid) return -1.0;
194  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
195  if(kNueEnergy2018(sr) >0) return (kRecoHADE(sr)/kNueEnergy2018(sr));
196  else return -5.;
197  });
198 
199 const Var kNueHadEFrac2020([](const caf::SRProxy* sr)
200  {
201  if(!sr->vtx.elastic.IsValid) return -1.0;
202  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1.0;
203  if(kNueEnergy2020(sr) > 0.) return (kRecoHADE(sr)/kNueEnergy2020(sr));
204  else return -5.;
205  });
206 
207 
208 
209 const Var kNueHadCalEFrac([](const caf::SRProxy* sr)
210  {
211  return kHadCalE(sr)/sr->slc.calE;
212  });
213 
214 #define PRONGVAR(VAR) \
215  Var( \
216  [](const caf::SRProxy* sr) \
217  { \
218  if(!sr->vtx.elastic.IsValid) return -1000.; \
219  if(sr->vtx.elastic.fuzzyk.npng < 1) return -1000.; \
220  return double(sr->vtx.elastic.fuzzyk.png[0].cvnpart.VAR); \
221  })
222 
223 const Var kProtonCVN = PRONGVAR(protonid);
224 const Var kElectronCVN = PRONGVAR(electronid);
225 const Var kPhotonCVN = PRONGVAR(photonid);
226 const Var kPizeroCVN = PRONGVAR(pizeroid);
227 const Var kNeutronCVN = PRONGVAR(neutronid);
228 const Var kOtherCVN = PRONGVAR(otherid);
230 
231 #define ELECIDVAR(VAR) \
232  Var( \
233  [](const caf::SRProxy* sr) \
234  { \
235  if(!sr->vtx.elastic.IsValid) return -1000.; \
236  if(sr->vtx.elastic.fuzzyk.nshwlid < 1) return -1000.; \
237  return double(sr->vtx.elastic.fuzzyk.png[0].shwlid.lid.VAR); \
238  })
239 
240 const Var kCosTheta = ELECIDVAR(costheta);
241 const Var kPi0Mass = ELECIDVAR(pi0mass);
242 const Var kShwEFrac = ELECIDVAR(shwEFrac);
243 
244 const Var kMuLLL = ELECIDVAR(mulll);
245 const Var kMuLLT = ELECIDVAR(mullt);
246 const Var kELLL = ELECIDVAR(elll);
247 const Var kELLT = ELECIDVAR(ellt);
248 const Var kEGLLL = ELECIDVAR(eglll);
249 const Var kEGLLT = ELECIDVAR(egllt);
250 const Var kEPLLL = ELECIDVAR(eplll);
251 const Var kEPLLT = ELECIDVAR(epllt);
252 const Var kENLLL = ELECIDVAR(enlll);
253 const Var kENLLT = ELECIDVAR(enllt);
254 const Var kEMuLLL = ELECIDVAR(emulll);
255 const Var kEMuLLT = ELECIDVAR(emullt);
256 const Var kEPi0LLL = ELECIDVAR(epi0lll);
257 const Var kEPi0LLT = ELECIDVAR(epi0llt);
258 const Var kEPiLLL = ELECIDVAR(epilll);
259 const Var kEPiLLT = ELECIDVAR(epillt);
260 
261 const Var kFirstProngCosTheta( []( const caf::SRProxy* sr )
262 {
263  if ( !sr->vtx.elastic.IsValid ) return -1000.f;
264  if ( sr->vtx.elastic.fuzzyk.npng < 1 ) return -1000.f;
265  TVector3 prongDir = (TVector3)sr->vtx.elastic.fuzzyk.png[0].dir;
266  TVector3 beamDir = ana::NuMIBeamDirection( sr->hdr.det );
267  return float( prongDir.Dot( beamDir ) );
268 });
269 
270 const Var kFirstProngTheta( []( const caf::SRProxy* sr )
271 {
272  if ( !sr->vtx.elastic.IsValid ) return -1000.f;
273  if ( sr->vtx.elastic.fuzzyk.npng < 1 ) return -1000.f;
274  float ct = kFirstProngCosTheta( sr );
275  if ( fabs( ct ) > 1 ) return -1000.f;
276  return float( acos( ct ) * 180.0 / 3.1416 );
277 });
278 
279 const Var kCVNeQE([](const caf::SRProxy* sr){
280  float cvnOut=-1;
281  std::cout << "ERROR::kCVNeQE. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
282  abort();
283  //if(sr->sel.cvnProd3Train.noutput>1){
284  // cvnOut=sr->sel.cvnProd3Train.output[4];
285  //}
286  return cvnOut;});
287 
288 const Var kCVNeRES([](const caf::SRProxy* sr){
289  float cvnOut=-1;
290  std::cout << "ERROR::kCVNeRES. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
291  abort();
292  //if(sr->sel.cvnProd3Train.noutput>1){
293  // cvnOut=sr->sel.cvnProd3Train.output[5];
294  //}
295  return cvnOut;});
296 
297 const Var kCVNeDIS([](const caf::SRProxy* sr){
298  float cvnOut=-1;
299  std::cout << "ERROR::kCVNeDIS. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
300  abort();
301  //if(sr->sel.cvnProd3Train.noutput>1){
302  // cvnOut=sr->sel.cvnProd3Train.output[6];
303  //}
304  return cvnOut;});
305 
306 const Var kCVNeOTHER([](const caf::SRProxy* sr){
307  float cvnOut=-1;
308  std::cout << "ERROR::kCVNeOTHER. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
309  abort();
310  //if(sr->sel.cvnProd3Train.noutput>1){
311  // cvnOut=sr->sel.cvnProd3Train.output[7];
312  //}
313  return cvnOut;});
314 
315 
316 const Var kCVNmQE([](const caf::SRProxy* sr){
317  float cvnOut=-1;
318  std::cout << "ERROR::kCVNeQE. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
319  abort();
320  //if(sr->sel.cvnProd3Train.noutput>1){
321  // cvnOut=sr->sel.cvnProd3Train.output[0];
322  //}
323  return cvnOut;});
324 
325 const Var kCVNmRES([](const caf::SRProxy* sr){
326  float cvnOut=-1;
327  std::cout << "ERROR::kCVNmRES. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
328  abort();
329  //if(sr->sel.cvnProd3Train.noutput>1){
330  // cvnOut=sr->sel.cvnProd3Train.output[1];
331  //}
332  return cvnOut;});
333 
334 const Var kCVNmDIS([](const caf::SRProxy* sr){
335  float cvnOut=-1;
336  std::cout << "ERROR::kCVNmDIS. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
337  abort();
338  //if(sr->sel.cvnProd3Train.noutput>1){
339  // cvnOut=sr->sel.cvnProd3Train.output[2];
340  //}
341  return cvnOut;});
342 
343 const Var kCVNmOTHER([](const caf::SRProxy* sr){
344  float cvnOut=-1;
345  std::cout << "ERROR::kCVNmOTHER. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
346  abort();
347  //if(sr->sel.cvnProd3Train.noutput>1){
348  // cvnOut=sr->sel.cvnProd3Train.output[3];
349  //}
350  return cvnOut;});
351 
352 const Var kConfusion([](const caf::SRProxy* sr)
353  {
354  if(sr->sel.lem.pid < .47 && sr->sel.lid.ann < .63) return .5;
355  if(sr->sel.lem.pid < .47 && sr->sel.lid.ann >= .63) return 1.5;
356  if(sr->sel.lem.pid >= .47 && sr->sel.lid.ann < .63) return 2.5;
357  if(sr->sel.lem.pid >= .47 && sr->sel.lid.ann >= .63) return 3.5;
358  return -5.;
359  });
360 
361 const Var kConfusionCVNLEM([](const caf::SRProxy* sr)
362  {
363  std::cout << "ERROR::kConfusionCVNLEM. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
364  abort();
365  /*
366  if(sr->sel.cvnProd3Train.nueid < .75 && sr->sel.lem.pid < .47) return .5;
367  if(sr->sel.cvnProd3Train.nueid < .75 && sr->sel.lem.pid >= .47) return 1.5;
368  if(sr->sel.cvnProd3Train.nueid >= .75 && sr->sel.lem.pid < .47) return 2.5;
369  if(sr->sel.cvnProd3Train.nueid >= .75 && sr->sel.lem.pid >= .47) return 3.5;
370  */
371  return -5.;
372  });
373 
374 const Var kConfusionCVNLID([](const caf::SRProxy* sr)
375  {
376  std::cout << "ERROR::kConfusionCVNLID. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
377  abort();
378  /*
379  if(sr->sel.cvnProd3Train.nueid < .75 && sr->sel.lid.ann < .63) return .5;
380  if(sr->sel.cvnProd3Train.nueid < .75 && sr->sel.lid.ann >= .63) return 1.5;
381  if(sr->sel.cvnProd3Train.nueid >= .75 && sr->sel.lid.ann < .63) return 2.5;
382  if(sr->sel.cvnProd3Train.nueid >= .75 && sr->sel.lid.ann >= .63) return 3.5;
383  */
384  return -5.;
385  });
386 
387 const Var kCVNMode([](const caf::SRProxy* sr)
388  {
389  std::cout << "ERROR::kCVNMode. Looking for cvnProd3Train. Branch no longer exists." << std::endl;
390  abort();
391  if(sr->sel.cvn.nueid<0){
392  return -1.;
393  }
394  /*
395  double nueVals[4]={sr->sel.cvnProd3Train.output[4],
396  sr->sel.cvnProd3Train.output[5],
397  sr->sel.cvnProd3Train.output[6],
398  sr->sel.cvnProd3Train.output[7]};
399  double temp=0;
400  int arrayEntry=-1;
401  for(int i=0;i<4;i++)
402  {
403  if(nueVals[i]>temp){
404  temp=nueVals[i];
405  arrayEntry=i;
406  }
407  }
408  if(arrayEntry==-1) return -1.;
409  if(arrayEntry==0) return .5;
410  if(arrayEntry==1) return 1.5;
411  if(arrayEntry==2) return 2.5;
412  if(arrayEntry==3) return 3.5;
413  */
414  return -5.;
415  });
416 
417 const Var kPIDExp = SIMPLEVAR(sel.lem.pidexp);
418 const Var kMeany = SIMPLEVAR(sel.lem.meanyexp);
419 const Var kQFrac = SIMPLEVAR(sel.lem.meanqfracexp);
420 const Var kEDiff = SIMPLEVAR(sel.lem.energydiffexp);
421 const Var kEnrich = SIMPLEVAR(sel.lem.enrichfracexp);
422  //const Var kNElastic = SIMPLEVAR(vtx.elastic.IsValid);
423 
424 const Var kNElastic([](const caf::SRProxy* sr)
425  {
426  if(!sr->vtx.elastic.IsValid){
427  return 0.;
428  }
429  else{
430  return 1.;
431  }
432  });
433 
434 const Var kNShwLID ([](const caf::SRProxy* sr)
435  {
436  if(!sr->vtx.elastic.IsValid) return 0u;
437  return (unsigned int)sr->vtx.elastic.fuzzyk.nshwlid;
438  });
439 
440 const Var kPtpMassless([](const caf::SRProxy* sr)
441  {
442  // This happened in the MC
443  if(std::isnan(1.*sr->sel.nuecosrej.photptp) ||
444  std::isinf(1.*sr->sel.nuecosrej.photptp)) return -1.f;
445  return float(sr->sel.nuecosrej.photptp);
446  });
447 
448 
449 const Var kHitsPerPlane = SIMPLEVAR(sel.nuecosrej.hitsperplane);
450 
451 const Var kCosdang = SIMPLEVAR(sel.nuecosrej.cosdang);
452 
453 
455  [](const caf::SRProxy* sr)
456  {
457  if(!sr->vtx.elastic.IsValid) return -5.f;
458  if(sr->vtx.elastic.fuzzyk.nshwlid >= 2)
459  return (sr->vtx.elastic.fuzzyk.png[0].shwlid.dir).Dot(
460  sr->vtx.elastic.fuzzyk.png[1].shwlid.dir);
461  else return -5.f;
462  });
463 
464 //asymmetry between x and y hits of most energetic shower
465 const Var kShwHitXYAsymm(
466  [](const caf::SRProxy* sr)
467  {
468  if(!sr->vtx.elastic.IsValid) return -5.f;
469  if(sr->vtx.elastic.fuzzyk.nshwlid < 1) return -5.f;
470 
471  float xyhitdiff = sr->vtx.elastic.fuzzyk.png[0].shwlid.nhitx -
472  sr->vtx.elastic.fuzzyk.png[0].shwlid.nhity;
473 
474  float xyhitsum = (sr->vtx.elastic.fuzzyk.png[0].shwlid.nhitx +
475  sr->vtx.elastic.fuzzyk.png[0].shwlid.nhity);
476 
477  if(xyhitsum>0) return xyhitdiff/xyhitsum;
478  else return -5.f;
479  });
480 
481 //unsigned asymmetry of x vs y hits
483  [](const caf::SRProxy* sr)
484  {
485  return abs(kShwHitXYAsymm(sr));
486  });
487 
488 /// Fraction of hits in a 3D shower out of total slice
489 const Var kShwHitFrac(
490  [](const caf::SRProxy* sr)
491  {
492  if(sr->slc.nhit <1) return -5.f;
493  if(!sr->vtx.elastic.IsValid) return 0.f;
494  float nhitall = 0.0;
495  for(unsigned int j=0; j<sr->vtx.elastic.fuzzyk.nshwlid; j++){
496  nhitall += sr->vtx.elastic.fuzzyk.png[j].shwlid.nhit;
497  }
498  return (nhitall/sr->slc.nhit);
499  });
500 
501 const Var kNPlanesToFront = SIMPLEVAR(sel.contain.nplanestofront);
502 
503 ///Count number of prongs identified as protons by CVN prong.
504 ///Applying a CVN cut of 0.9, and protons within 20 cm of vertex
505 ///used in wrong-sign studies in docdb-25899
506 const Var kProtonCount(
507  [](const caf::SRProxy* sr)
508  {
509  int count = 0;
510  const caf::SRVector3DProxy& vtx = sr->vtx.elastic.vtx;
511  for (unsigned int i=0; i< sr->vtx.elastic.fuzzyk.npng; ++i){
512  if (sr->vtx.elastic.fuzzyk.png[i].cvnpart.pdgmax != 2212) continue;
513  const caf::SRVector3DProxy& st = sr->vtx.elastic.fuzzyk.png[i].start;
514  if (util::pythag(vtx.X()-st.X(),vtx.Y()-st.Y(),vtx.Z()-st.Z())>20)
515  if (sr->vtx.elastic.fuzzyk.png[i].cvnpart.maxval < 0.9) continue;
516  count++;
517  }
518  return count;
519  });
520 
521 ///return CVN proton score for most proton-like prong
522 const Var kProngCVNMaxProton(
523  [](const caf::SRProxy* sr)
524  {
525  double score = -5;
526  for (unsigned int i=0; i< sr->vtx.elastic.fuzzyk.npng; ++i){
527  if (sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid > score)
528  score = sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid;
529  }
530  return score;
531  });
532 
533 ///return distance from vertex for most proton-like prong
535  [](const caf::SRProxy* sr)
536  {
537  double score = -5;
538  double dist = -5;
539  const caf::SRVector3DProxy& vtx = sr->vtx.elastic.vtx;
540  for (unsigned int i=0; i< sr->vtx.elastic.fuzzyk.npng; ++i){
541  if (sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid > score){
542  score = sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid;
543  const caf::SRVector3DProxy& st = sr->vtx.elastic.fuzzyk.png[i].start;
544  dist = util::pythag(vtx.X()-st.X(),vtx.Y()-st.Y(),vtx.Z()-st.Z());
545  }
546  }
547  return dist;
548  });
549 
550 }//namespace
const XML_Char int len
Definition: expat.h:262
caf::Proxy< size_t > npng
Definition: SRProxy.h:2038
caf::Proxy< unsigned int > nshwlid
Definition: SRProxy.h:2040
caf::Proxy< float > pid
Definition: SRProxy.h:949
const Var kShwNHit
int isinf(const stan::math::var &a)
Definition: std_isinf.hpp:16
caf::Proxy< caf::SRFuzzyK > fuzzyk
Definition: SRProxy.h:2059
const Var kProtonCount([](const caf::SRProxy *sr){int count=0;const caf::SRVector3DProxy &vtx=sr->vtx.elastic.vtx;for(unsigned int i=0;i< sr->vtx.elastic.fuzzyk.npng;++i){if(sr->vtx.elastic.fuzzyk.png[i].cvnpart.pdgmax!=2212) continue;const caf::SRVector3DProxy &st=sr->vtx.elastic.fuzzyk.png[i].start;if(util::pythag(vtx.X()-st.X(), vtx.Y()-st.Y(), vtx.Z()-st.Z())>20) if(sr->vtx.elastic.fuzzyk.png[i].cvnpart.maxval< 0.9) continue;count++;}return count;})
Definition: NueVarsExtra.h:131
const Var kCVNemE([](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;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;}return CVNem_CalE *CalibrationBugCorrectionFactor(sr->hdr);})
Definition: NueEnergy2017.h:9
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Var kConfusionCVNLEM([](const caf::SRProxy *sr){std::cout<< "ERROR::kConfusionCVNLEM. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort();return-5.;})
Definition: NueVarsExtra.h:96
const Var kEMNhit([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;double Nhit=0;for(const caf::SRFuzzyKProngProxy &png:sr->vtx.elastic.fuzzyk.png){double hits=png.nhit;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 ) Nhit+=hits;}return Nhit;})
Definition: NueVarsExtra.h:18
const Var kHADEPerHit([](const caf::SRProxy *sr){if((sr->slc.nhit-kEMNhit(sr))>0) return 1000.0 *(sr->slc.calE *CalibrationBugCorrectionFactor(sr->hdr)-kemCalE(sr))/(sr->slc.nhit-kEMNhit(sr));else return-5.;})
Definition: NueVarsExtra.h:21
const Var kConfusionCVNLID([](const caf::SRProxy *sr){std::cout<< "ERROR::kConfusionCVNLID. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort();return-5.;})
Definition: NueVarsExtra.h:97
const Var kEPi0LLL
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
const Var kCVNmQE([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNeQE. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:90
const Var kProngCVNMaxProtonDist([](const caf::SRProxy *sr){double score=-5;double dist=-5;const caf::SRVector3DProxy &vtx=sr->vtx.elastic.vtx;for(unsigned int i=0;i< sr->vtx.elastic.fuzzyk.npng;++i){if(sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid > score){score=sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid;const caf::SRVector3DProxy &st=sr->vtx.elastic.fuzzyk.png[i].start;dist=util::pythag(vtx.X()-st.X(), vtx.Y()-st.Y(), vtx.Z()-st.Z());}}return dist;})
return distance from vertex for most proton-like prong
Definition: NueVarsExtra.h:137
const Var kFirstProngCosTheta([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1000.f;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1000.f;TVector3 prongDir=(TVector3) sr->vtx.elastic.fuzzyk.png[0].dir;TVector3 beamDir=ana::NuMIBeamDirection(sr->hdr.det);return float(prongDir.Dot(beamDir));})
Definition: NueVarsExtra.h:62
const Var kemCalE([](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;for(const caf::SRFuzzyKProngProxy &png:sr->vtx.elastic.fuzzyk.png){double png_CalE=png.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;}return CVNem_CalE *CalibrationBugCorrectionFactor(sr->hdr);})
Definition: NueVarsExtra.h:17
double CalibrationBugCorrectionFactor(const caf::SRHeaderProxy &hdr)
See docdb 23597.
Definition: Vars.cxx:141
double NueRecoE_2018RHCFit(double rawEM, double rawHA)
caf::Proxy< caf::SRHeader > hdr
Definition: SRProxy.h:2137
const Var kCosdang
const Var kCVNmDIS([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNmDIS. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:92
const Var kShwStopX
const Var kShwCalE
Proxy for caf::StandardRecord.
Definition: SRProxy.h:2126
const Var kEGLLL
const Var kCVNhadE([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;const double calE=sr->slc.calE *CalibrationBugCorrectionFactor(sr->hdr);return std::max(calE-kCVNemE(sr), 0.);})
Definition: NueEnergy2017.h:10
const Var kemE2017([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;return NueRecoE_2017FDFit(kCVNemE(sr), 0.0);})
Definition: NueVarsExtra.h:42
const Var kEGLLT
const Var kENLLT
const Var kShwStartY
nhit
Definition: demo1.py:25
const Var kNueHadEFrac2020([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;if(kNueEnergy2020(sr) > 0.) return(kRecoHADE(sr)/kNueEnergy2020(sr));else return-5.;})
Definition: NueVarsExtra.h:50
void abs(TH1 *hist)
const Var kNueHadEFrac2([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;if(kNueEnergy2018(sr) >0) return(kRecoHADE(sr)/kNueEnergy2018(sr));else return-5.;})
Definition: NueVarsExtra.h:49
const Var kNueHadEFrac([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;if(kNueEnergy2018(sr)>0){if(kIsRHC(sr)) return(NueRecoE_2018RHCFit(0.0, kRecoHADE(sr)))/(NueRecoE_2018RHCFit(kRecoEME(sr), kRecoHADE(sr)));else return(NueRecoE_2017FDFit(0.0, kRecoHADE(sr)))/NueRecoE_2017FDFit(kRecoEME(sr), kRecoHADE(sr));}else return-5.;})
Definition: NueVarsExtra.h:48
const Var kPizeroCVN
Proxy for caf::SRFuzzyKProng.
Definition: SRProxy.h:2003
const Var kCVNmRES([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNmRES. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:91
const Var kEPerHit([](const caf::SRProxy *sr){if(sr->slc.nhit >0) return 1000.0 *(sr->slc.calE/sr->slc.nhit);else return-5.;})
Definition: NueVarsExtra.h:14
const Var kVtxX
T acos(T number)
Definition: d0nt_math.hpp:54
A PID for muons.
Definition: FillPIDs.h:11
const Var kShwEFrac
const Var kHadCalE([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return float(sr->slc.calE);if(sr->vtx.elastic.fuzzyk.nshwlid< 1) return float(sr->slc.calE);return((sr->slc.calE- sr->vtx.elastic.fuzzyk.png[0].shwlid.calE));})
Definition: NueVarsExtra.h:40
_Var< caf::SRProxy > Var
Representation of a variable to be retrieved from a caf::StandardRecord object.
Definition: Var.h:74
int isnan(const stan::math::var &a)
Definition: std_isnan.hpp:18
const Var kMeany
const Var kShwWidth
double dist
Definition: runWimpSim.h:113
const Var kEnrich
const Var kPtpMassless([](const caf::SRProxy *sr){if(std::isnan(1.*sr->sel.nuecosrej.photptp)|| std::isinf(1.*sr->sel.nuecosrej.photptp)) return-1.f;return float(sr->sel.nuecosrej.photptp);})
Definition: NueVarsExtra.h:110
caf::Proxy< caf::SRElastic > elastic
Definition: SRProxy.h:2118
caf::Proxy< caf::SRNueCosRej > nuecosrej
Definition: SRProxy.h:1265
const Var kENLLL
const Var kHitsPerPlane
caf::Proxy< std::vector< caf::SRFuzzyKProng > > png
Definition: SRProxy.h:2043
const Var kEMuLLT
const Var kShwStopY
const Var kCVNeOTHER([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNeOTHER. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:89
caf::Proxy< unsigned int > nhit
Definition: SRProxy.h:1315
const Var kCVNmOTHER([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNmOTHER. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:93
const Var kMuLLT
const Var kShwStartX
const Var kCosAngleNextShower([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-5.f;if(sr->vtx.elastic.fuzzyk.nshwlid >=2) return(sr->vtx.elastic.fuzzyk.png[0].shwlid.dir).Dot(sr->vtx.elastic.fuzzyk.png[1].shwlid.dir);else return-5.f;})
Definition: NueVarsExtra.h:116
if(dump)
#define PRONGVAR(VAR)
Proxy for caf::SRVector3D.
Definition: SRProxy.h:78
const Var kNueEnergy2018([](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 kNueEnergyRHC(sr);else return kNueEnergyFHC(sr);})
Definition: NueEnergy2018.h:25
const Var kEPiLLL
caf::Proxy< caf::SRCVNResult > cvn
Definition: SRProxy.h:1253
caf::Proxy< float > photptp
Definition: SRProxy.h:1058
const Var kConfusion([](const caf::SRProxy *sr){if(sr->sel.lem.pid< .47 &&sr->sel.lid.ann< .63) return.5;if(sr->sel.lem.pid< .47 &&sr->sel.lid.ann >=.63) return 1.5;if(sr->sel.lem.pid >=.47 &&sr->sel.lid.ann< .63) return 2.5;if(sr->sel.lem.pid >=.47 &&sr->sel.lid.ann >=.63) return 3.5;return-5.;})
Definition: NueVarsExtra.h:95
caf::Proxy< float > nueid
Definition: SRProxy.h:906
double NueRecoE_2017FDFit(double rawEM, double rawHA)
#define ELECIDVAR(VAR)
const Var kEMuLLL
const Var kEPiLLT
const Var kEPi0LLT
const Var kRecoEME([](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(CVNem_CalE==0.0) CVNem_CalE=prim_png.shwlid.calE;return CVNem_CalE *CalibrationBugCorrectionFactor(sr->hdr);})
Definition: NueEnergy2018.h:17
const Var kShwGap
caf::StandardRecord * sr
const Var kMuLLL
const double j
Definition: BetheBloch.cxx:29
const Var kCVNeQE([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNeQE. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:86
const Var kShwStopZ
const Var kFirstProngTheta([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1000.f;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1000.f;float ct=kFirstProngCosTheta(sr);if(fabs(ct) > 1) return-1000.f;return float(acos(ct)*180.0/3.1416);})
Definition: NueVarsExtra.h:63
const Var kELLT
const Var kCVNMode([](const caf::SRProxy *sr){std::cout<< "ERROR::kCVNMode. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort();if(sr->sel.cvn.nueid< 0){return-1.;}return-5.;})
Definition: NueVarsExtra.h:99
const Var kEDiff
z
Definition: test.py:28
#define SHWLIDVAR(VAR)
const Var kELLL
caf::Proxy< caf::SRLem > lem
Definition: SRProxy.h:1260
const Var kRecoEPerHit([](const caf::SRProxy *sr){if(sr->slc.nhit >0) return 1000.0 *(kNueEnergy2018(sr)/sr->slc.nhit);else return-5.;})
Definition: NueVarsExtra.h:15
OStream cout
Definition: OStream.cxx:6
const Var kNueHadCalEFrac([](const caf::SRProxy *sr){return kHadCalE(sr)/sr->slc.calE;})
Definition: NueVarsExtra.h:52
const Var khadE2018([](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 NueRecoE_2018RHCFit(0.0, kRecoHADE(sr));else return NueRecoE_2017FDFit(0.0, kRecoHADE(sr));})
Definition: NueVarsExtra.h:46
caf::Proxy< caf::SRELid > lid
Definition: SRProxy.h:1261
const Var kPIDExp
caf::Proxy< float > ann
Definition: SRProxy.h:971
const Var kInelasticity([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return 1.f;if(sr->vtx.elastic.fuzzyk.nshwlid< 1) return 1.f;return((sr->slc.calE- sr->vtx.elastic.fuzzyk.png[0].shwlid.calE)/(sr->slc.calE));})
Definition: NueVarsExtra.h:38
const Cut kIsRHC([](const caf::SRProxy *sr){return sr->spill.isRHC;})
Definition: Vars.h:16
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
const Var kCVNeDIS([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNeDIS. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:88
const Var kPhotonCVN
const Var kMuonCVN
caf::Proxy< bool > IsValid
Definition: SRProxy.h:2058
const Var kProngCVNMaxProton([](const caf::SRProxy *sr){double score=-5;for(unsigned int i=0;i< sr->vtx.elastic.fuzzyk.npng;++i){if(sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid > score) score=sr->vtx.elastic.fuzzyk.png[i].cvnpart.protonid;}return score;})
return CVN proton score for most proton-like prong
Definition: NueVarsExtra.h:134
const Var kShwHitFrac([](const caf::SRProxy *sr){if(sr->slc.nhit< 1) return-5.f;if(!sr->vtx.elastic.IsValid) return 0.f;float nhitall=0.0;for(unsigned int j=0;j< sr->vtx.elastic.fuzzyk.nshwlid;j++){nhitall+=sr->vtx.elastic.fuzzyk.png[j].shwlid.nhit;}return(nhitall/sr->slc.nhit);})
Fraction of hits in a 3D shower out of total slice.
Definition: NueVarsExtra.h:125
const Var kVtxY
caf::Proxy< caf::SRSlice > slc
Definition: SRProxy.h:2142
const Var kEMEPerHit([](const caf::SRProxy *sr){if(kEMNhit(sr)>0) return 1000.0 *(kemCalE(sr)/kEMNhit(sr));else return-5.;})
Definition: NueVarsExtra.h:20
caf::Proxy< float > calE
Definition: SRProxy.h:1292
double pythag(double x, double y)
2D Euclidean distance
Definition: MathUtil.h:29
const Var kShwHitXYUnsignedAsymm([](const caf::SRProxy *sr){return abs(kShwHitXYAsymm(sr));})
Asymmetry of x vs y hits.
Definition: NueVarsExtra.h:122
const Var kPi0Mass
const Var kNPlanesToFront
#define SIMPLEVAR(CAFNAME)
For Vars where literally all you need is a single CAF variable.
Definition: Var.h:88
const Var kemE2018([](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 NueRecoE_2018RHCFit(kRecoEME(sr), 0.0);else return NueRecoE_2017FDFit(kRecoEME(sr), 0.0);})
Definition: NueVarsExtra.h:45
const Var kElectronCVN
const Var kShwHitXYAsymm([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-5.f;if(sr->vtx.elastic.fuzzyk.nshwlid< 1) return-5.f;float xyhitdiff=sr->vtx.elastic.fuzzyk.png[0].shwlid.nhitx-sr->vtx.elastic.fuzzyk.png[0].shwlid.nhity;float xyhitsum=(sr->vtx.elastic.fuzzyk.png[0].shwlid.nhitx+ sr->vtx.elastic.fuzzyk.png[0].shwlid.nhity);if(xyhitsum >0) return xyhitdiff/xyhitsum;else return-5.f;})
asymmetry between x and y hits of most energetic shower
Definition: NueVarsExtra.h:119
const Var kVtxZ
const Var kEPLLT
caf::Proxy< caf::SRVector3D > vtx
Definition: SRProxy.h:2073
#define ELASTICVAR(VAR)
TVector3 NuMIBeamDirection(caf::Det_t det)
Average direction of NuMI neutrinos in a given detector This function is a copy of geo::GeometryBase:...
Definition: Utilities.cxx:329
const Var kEPLLL
const Var kNShwLID([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return 0u;return(unsigned int) sr->vtx.elastic.fuzzyk.nshwlid;})
Definition: NueVarsExtra.h:108
caf::Proxy< caf::SRIDBranch > sel
Definition: SRProxy.h:2141
const Var kCosTheta
caf::Proxy< caf::SRVertexBranch > vtx
Definition: SRProxy.h:2146
const Var kNElastic([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid){return 0.;}else{return 1.;}})
Definition: NueVarsExtra.h:106
const Var kShwStartZ
const Var kCVNeRES([](const caf::SRProxy *sr){float cvnOut=-1;std::cout<< "ERROR::kCVNeRES. Looking for cvnProd3Train. Branch no longer exists."<< std::endl;abort(); return cvnOut;})
Definition: NueVarsExtra.h:87
const Var kOtherCVN
const Var kShwLen
const Var kNeutronCVN
const Var kProtonCVN
const Var khadE2017([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;return NueRecoE_2017FDFit(0.0, kCVNhadE(sr));})
Definition: NueVarsExtra.h:43
const Var kQFrac
caf::Proxy< caf::Det_t > det
Definition: SRProxy.h:232
const Var kRecoHADE([](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 *CalibrationBugCorrectionFactor(sr->hdr);return std::max(CVNha_CalE-kRecoEME(sr), 0.);})
Definition: NueEnergy2018.h:18