VarVals.h
Go to the documentation of this file.
1 /*
2  * \file VarVals.h
3  * \brief Base container for the Vars that constitute an event.
4  *
5  * Created on: \date June 13, 2016
6  * Original author: \author J. Wolcott <jwolcott@fnal.gov>
7  *
8  */
9 
10 #ifndef COVARIANCEMATRIXFIT_CORE_VARVALS_H_
11 #define COVARIANCEMATRIXFIT_CORE_VARVALS_H_
12 
13 #include <string>
14 #include <cstdint>
15 
16 #include "cetlib_except/exception.h"
18 
22 
23 namespace cmf
24 {
25  struct TruthVars{
26 
28  : fTrueE (std::numeric_limits<float>::lowest())
29  , fTruePDG (std::numeric_limits<float>::lowest())
30  , fTrueCCNC (std::numeric_limits<float>::lowest())
31  , fTrueIntType (std::numeric_limits<float>::lowest())
32  , fTruePDGOrig (std::numeric_limits<float>::lowest())
33  , fTrueParentPDG (std::numeric_limits<float>::lowest())
34  , fTrueParentPT (std::numeric_limits<float>::lowest())
35  , fTrueParentPZ (std::numeric_limits<float>::lowest())
36  , fTrueParentDecay (std::numeric_limits<float>::lowest())
37  , fTrueParentTargetPDG(std::numeric_limits<float>::lowest())
38  , fTrueHitNuc (std::numeric_limits<float>::lowest())
39  , fTrueIntMode (std::numeric_limits<float>::lowest())
40  {}
41 
42  float fTrueE; ///< True nu energy
43  float fTruePDG; ///< true PDG
44  float fTrueCCNC; ///< true cc vs nc
45  float fTrueIntType; ///< true interaction type
46  float fTruePDGOrig; ///< true PDG of original neutrino
47  float fTrueParentPDG; ///< true PDG of neutrino parent
48  float fTrueParentPT; ///< true p_T of neutrino parent
49  float fTrueParentPZ; ///< true p_Z of neutrino parent
50  float fTrueParentDecay; ///< true parent decay mode
51  float fTrueParentTargetPDG; ///< true parent pdg code off the target
52  float fTrueHitNuc; ///< true hit nucleus
53  float fTrueIntMode; ///< true interaction mode
54  };
55 
56  class EventId{
57 
58  public:
60  : run (std::numeric_limits<int>::lowest())
61  , subrun(std::numeric_limits<int>::lowest())
62  , event (std::numeric_limits<int>::lowest())
63  , cycle (std::numeric_limits<int>::lowest())
64  , slice (std::numeric_limits<int>::lowest())
65  {}
66 
67  EventId(int r,
68  int s,
69  int e,
70  int c,
71  int slc)
72  : run (r)
73  , subrun(s)
74  , event (e)
75  , cycle (c)
76  , slice (slc)
77  {}
78 
79  int run;
80  int subrun;
81  int event;
82  int cycle;
83  int slice;
84 
85  friend std::ostream& operator << (std::ostream & o, cmf::EventId const& evid);
86 
87  bool operator<(EventId const& evid) const;
88  bool operator==(EventId const& evid) const;
89 
90  };
91 
92  struct DataVars{
93 
94  // initialize the Fake_Weight to be 1, it will be reset if Fake_Data are generated
96  : fNuE_CVN (std::numeric_limits<float>::lowest())
97  , fNuE_NumMichel (std::numeric_limits<float>::lowest())
98  , fNu_RecoE (std::numeric_limits<float>::lowest())
99  , fHad_RecoE (std::numeric_limits<float>::lowest())
100  , fLep_RecoE (std::numeric_limits<float>::lowest())
101  , fLep_RecoE_MCFrac(0.)
102  , fRecoQ2 (std::numeric_limits<float>::lowest())
103  , fFake_Weight (1.)
104  {}
105 
106  DataVars(float nuecvn,
107  float nuenummichel,
108  float nurece,
109  float hadrece,
110  float leprece,
111  float leprecemcfrac,
112  float recoQ2,
113  float fakewgt)
114  : fNuE_CVN (nuecvn)
115  , fNuE_NumMichel (nuenummichel)
116  , fNu_RecoE (nurece)
117  , fHad_RecoE (hadrece)
118  , fLep_RecoE (leprece)
119  , fLep_RecoE_MCFrac(leprecemcfrac)
120  , fRecoQ2 (recoQ2)
121  , fFake_Weight (fakewgt)
122  {}
123 
124  float fNuE_CVN; ///<
125  float fNuE_NumMichel; ///<
126  float fNu_RecoE; ///<
127  float fHad_RecoE; ///<
128  float fLep_RecoE; ///<
129  float fLep_RecoE_MCFrac; ///< fraction of leptonic energy in muon catcher
130  float fRecoQ2; ///< reconstructed Q^2
131  float fFake_Weight; ///< Weight for fake data events
132  };
133 
134  // set the default values to be unity, ie no weighting applied
135  struct WeightVars{
136 
138  : fXSecCVPPFX_Weight (1.)
139  , fRPACCQE_Weight (1.)
140  , fRPARES_Weight (1.)
141  , fDISvnCC1pi_Weight (1.)
142  , fEmpiricalMEC_Weight (1.)
143  , fEmpiricalMECtoGENIEQE_Weight (1.)
144  , fEmpiricalMECtoGENIERES_Weight(1.)
145  , fOtherDIS_Weight (1.)
146  , fXSecCV2020_Weight (1.)
147  , fPPFXFluxCV_Weight (1.)
148  {}
149 
150  float fXSecCVPPFX_Weight; ///< Deprecated, Was Tufts weight for SA
151  float fRPACCQE_Weight; ///< Deprecated
152  float fRPARES_Weight; ///< To be used for systematic evaluation ONLY
153  float fDISvnCC1pi_Weight; ///< Deprecated
154  float fEmpiricalMEC_Weight; ///< Deprecated
155  float fEmpiricalMECtoGENIEQE_Weight; ///< Deprecated
156  float fEmpiricalMECtoGENIERES_Weight; ///< Deprecated
157  float fOtherDIS_Weight; ///< Deprecated
158  float fXSecCV2020_Weight; ///< For 2020 Ana, store weights seperately
159  float fPPFXFluxCV_Weight; ///< For 2020 Ana, store weights seperately
160  };
161 
162  // struct to hold weights of systematic parameters with shifts related to
163  // different sigma shifts in the underlying parameter
164  struct SystVar {
165 
166  SystVar(uint8_t type)
167  : fType(type)
168  , fMinus2Sigma(1.)
169  , fMinus1Sigma(1.)
170  , fPlus1Sigma(1.)
171  , fPlus2Sigma(1.)
172  {}
173 
174  SystVar(uint8_t type,
175  float minus2,
176  float minus1,
177  float plus1,
178  float plus2)
179  : fType(type)
180  , fMinus2Sigma(minus2)
181  , fMinus1Sigma(minus1)
182  , fPlus1Sigma(plus1)
183  , fPlus2Sigma(plus2)
184  {}
185 
186  uint8_t Type() const { return fType; }
187  std::vector<float> SigmasAsVec() const;
188  float SigmaWeight(uint8_t const& sigmaLevel) const;
189 
190  bool operator==(cmf::SystVar const& other) const { return fType == other.Type(); }
191  bool operator <(cmf::SystVar const& other) const { return fType < other.Type(); }
192 
193  friend std::ostream& operator << (std::ostream & o, cmf::SystVar const& sv);
194 
195  uint8_t fType; ///< key for this parameter
196  float fMinus2Sigma; ///< weight associated with this change in the underlying parameter
197  float fMinus1Sigma; ///< weight associated with this change in the underlying parameter
198  float fPlus1Sigma; ///< weight associated with this change in the underlying parameter
199  float fPlus2Sigma; ///< weight associated with this change in the underlying parameter
200  };
201 
202  //----------------------------------------------------------------------------
203  inline std::vector<float> SystVar::SigmasAsVec() const
204  {
205  return std::move(std::vector<float>({fMinus2Sigma,
206  fMinus1Sigma,
207  fPlus1Sigma,
208  fPlus2Sigma}));
209  }
210 
211  //----------------------------------------------------------------------------
212  inline float SystVar::SigmaWeight(uint8_t const& sigmaLevel) const
213  {
214  if(sigmaLevel == cmf::kMinus2Sigma) return fMinus2Sigma;
215  else if(sigmaLevel == cmf::kMinus1Sigma) return fMinus1Sigma;
216  else if(sigmaLevel == cmf::kPlus1Sigma ) return fPlus1Sigma;
217  else if(sigmaLevel == cmf::kPlus2Sigma ) return fPlus2Sigma;
218 
219  // default to no change in weight
220  return 1.;
221  }
222 
223  typedef std::vector<cmf::SystVar> SystVarColl;
224 
225 
226  /// \brief Base container for the MC related Vars that constitute an event.
227  class MCVarVals
228  {
229  private:
233 
234  public:
236  {}
237 
239  cmf::WeightVars const& wv)
240  : fTruthVars (tv)
241  , fWeightVars(wv)
242  {}
243 
245  cmf::WeightVars const& wv,
246  cmf::SystVarColl const& gv)
247  : fTruthVars (tv)
248  , fWeightVars(wv)
249  , fSystematicShifts(gv)
250  {}
251 
252 
253  /// Access a Var by name.
254  float val_at(uint8_t const& varkey) const;
255 
256  std::map<float, float> SystVarWgts(uint8_t const& varkey) const;
257 
258  friend std::ostream& operator << (std::ostream & o, cmf::MCVarVals const& mcvv);
259 
260  // The following functions should only be used for making eventlist trees
261  cmf::TruthVars const& TruthVars() const { return fTruthVars; }
262  cmf::WeightVars const& WeightVars() const { return fWeightVars; }
263  cmf::SystVarColl const& SystematicShifts() const { return fSystematicShifts; }
264 
265  };
266 
267  /// \brief Base container for the MC related Vars that constitute an event.
269  {
270  private:
271 
272  float NuRecoEVal(cmf::MetaData const& md) const;
273 
275 
276  public:
277  DataVarVals () = default;
278 
280  : fDataVars(dv)
281  {}
282 
283  /// Access a Var by name.
284  /// Throws if the var doesn't exist
285  float val_at(uint8_t const& varkey,
286  cmf::MetaData const& md) const;
287 
288  void set_val_at(uint8_t const& varkey,
289  float const& val);
290 
291  friend std::ostream& operator << (std::ostream & o, cmf::DataVarVals const& dvv);
292 
293  // The following function should only be used for making eventlist trees
294  cmf::DataVars const& DataVars() const { return fDataVars; }
295 
296  };
297 
298  static float RecoEnergy(float leptonicEnergy,
299  float hadronicEnergy,
300  cmf::MetaData const& md);
301 
302  static float NuMuEnergy(float leptonicEnergy,
303  float hadronicEnergy);
304 
305  static float NuEEnergy(float leptonicEnergy,
306  float hadronicEnergy,
307  cmf::BeamType_t const& beamType);
308 
309  static float NCEnergy(float leptonicEnergy,
310  float hadronicEnergy,
311  cmf::BeamType_t const& beamType,
312  cmf::DetType_t const& detector);
313 
314 } /* namespace cmf */
315 
316 //----------------------------------------------------------------------------
317 inline float cmf::RecoEnergy(float leptonicEnergy,
318  float hadronicEnergy,
319  cmf::MetaData const& md)
320 {
321  if (md.IsNuMuSelected()) return cmf::NuMuEnergy(leptonicEnergy, hadronicEnergy );
322  else if(md.IsNCSelected() ) return cmf::NCEnergy (leptonicEnergy, hadronicEnergy, md.BeamType(), md.detector);
323  else if(md.IsNuESelected() ) return cmf::NuEEnergy (leptonicEnergy, hadronicEnergy, md.BeamType());
324  else
325  throw cet::exception("DataVarVals")
326  << "Unknown selection type "
327  << md.selectionType
328  << " requested to be used in calculating reconstructed energy";
329 
330  return std::numeric_limits<float>::lowest();
331 }
332 
333 //----------------------------------------------------------------------------
334 inline float cmf::NuMuEnergy(float leptonicEnergy,
335  float hadronicEnergy)
336 {
337  return leptonicEnergy + hadronicEnergy;
338 
339 }
340 
341 //----------------------------------------------------------------------------
342 inline float cmf::NuEEnergy(float leptonicEnergy,
343  float hadronicEnergy,
344  cmf::BeamType_t const& beamType)
345 {
346  // p0 and p3 are 0.0 for both FHC and RHC
347  // 2020
348  float p0 = 0.0;
349  float p1 = (beamType == kFHC) ? 1.01777 : 0.988258;
350  float p2 = (beamType == kFHC) ? 1.10868 : 1.20084;
351  float p3 = 0.0;
352  float p4 = (beamType == kFHC) ? 1.43541e-03 : 1.92904e-07;
353  float p5 = (beamType == kFHC) ? 1.09628e-01 : 1.20704e-07;
354  float fr = (beamType == kFHC) ? 1. / 1.0355622 : 1. / 1.0111873;
355 
356  return fr *(( p5 * hadronicEnergy * hadronicEnergy) +
357  ( p4 * leptonicEnergy * leptonicEnergy) +
358  ( p3 * leptonicEnergy * hadronicEnergy) +
359  ( p2 * hadronicEnergy) +
360  ( p1 * leptonicEnergy) +
361  p0 );
362 }
363 
364 //----------------------------------------------------------------------------
365 inline float cmf::NCEnergy(float leptonicEnergy,
366  float hadronicEnergy,
367  cmf::BeamType_t const& beamType ,
368  cmf::DetType_t const& detector )
369 {
370  // Note: leptonicEnergy must be the EM-like energy, as calculated by kEME_2020
371  // and hadronicEnergy must be the hadronic-like energy, as calculated by kHADE_2020
372  // If we don't have a tuned correction return the total calorimetric energy
373  if (beamType != kFHC ||
374  (detector != cmf::kNEARDET && detector != cmf::kFARDET))
375  return leptonicEnergy + hadronicEnergy;
376 
377  // Set up correction parameters
378  double pars[4][6] = { {1.049, 0.795, 0.8409, 0.17, 0.82, -1.00}, // ND FHC
379  {1.025, 0.797, 0.9162, 0.53, -0.26, -1.48}, // FD FHC
380  {1.000, 1.000, 1.0000, 0.00, 0.00, 0.00}, // ND RHC
381  {1.000, 1.000, 1.0000, 0.00, 0.00, 0.00} };// FD RHC
382 
383  // set detector/beam index
384  int detCur = (detector == cmf::kFARDET) + ((beamType == kRHC)<<1); // bitshift magic
385 
386  // return the energy estimator
387  return (leptonicEnergy / pars[detCur][0] + hadronicEnergy / pars[detCur][1])
388  / (pars[detCur][2] + pars[detCur][3] * std::pow(leptonicEnergy + hadronicEnergy + pars[detCur][4],2) * std::exp(pars[detCur][5] * (leptonicEnergy + hadronicEnergy)));
389 }
390 #endif /* COVARIANCEMATRIXFIT_CORE_VARVALS_H_ */
float fDISvnCC1pi_Weight
Deprecated.
Definition: VarVals.h:153
float fMinus2Sigma
weight associated with this change in the underlying parameter
Definition: VarVals.h:196
float fTruePDGOrig
true PDG of original neutrino
Definition: VarVals.h:46
static float NuMuEnergy(float leptonicEnergy, float hadronicEnergy)
Definition: VarVals.h:334
MCVarVals(cmf::TruthVars const &tv, cmf::WeightVars const &wv, cmf::SystVarColl const &gv)
Definition: VarVals.h:244
float fEmpiricalMECtoGENIEQE_Weight
Deprecated.
Definition: VarVals.h:155
static float NuEEnergy(float leptonicEnergy, float hadronicEnergy, cmf::BeamType_t const &beamType)
Definition: VarVals.h:342
bool operator==(cmf::SystVar const &other) const
Definition: VarVals.h:190
float fEmpiricalMEC_Weight
Deprecated.
Definition: VarVals.h:154
float fFake_Weight
Weight for fake data events.
Definition: VarVals.h:131
cmf::WeightVars fWeightVars
Definition: VarVals.h:231
float fRecoQ2
reconstructed Q^2
Definition: VarVals.h:130
EventId(int r, int s, int e, int c, int slc)
Definition: VarVals.h:67
cmf::DetType_t detector
Definition: Structs.h:114
float fPPFXFluxCV_Weight
For 2020 Ana, store weights seperately.
Definition: VarVals.h:159
float fTrueHitNuc
true hit nucleus
Definition: VarVals.h:52
float fTrueCCNC
true cc vs nc
Definition: VarVals.h:44
float fTrueParentDecay
true parent decay mode
Definition: VarVals.h:50
float fTruePDG
true PDG
Definition: VarVals.h:43
constexpr T pow(T x)
Definition: pow.h:75
float fNuE_CVN
Definition: VarVals.h:124
float fPlus2Sigma
weight associated with this change in the underlying parameter
Definition: VarVals.h:199
enum cmf::det_type DetType_t
float fTrueParentPT
true p_T of neutrino parent
Definition: VarVals.h:48
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
float fEmpiricalMECtoGENIERES_Weight
Deprecated.
Definition: VarVals.h:156
float fNu_RecoE
Definition: VarVals.h:126
uint8_t fType
key for this parameter
Definition: VarVals.h:195
cmf::DataVars const & DataVars() const
Definition: VarVals.h:294
float fRPARES_Weight
To be used for systematic evaluation ONLY.
Definition: VarVals.h:152
std::ostream & operator<<(std::ostream &o, cmf::Event const &event)
Definition: Event.cxx:27
cmf::SystVarColl fSystematicShifts
Definition: VarVals.h:232
float fTrueE
True nu energy.
Definition: VarVals.h:42
Base container for the MC related Vars that constitute an event.
Definition: VarVals.h:227
float fLep_RecoE_MCFrac
fraction of leptonic energy in muon catcher
Definition: VarVals.h:129
MCVarVals(cmf::TruthVars const &tv, cmf::WeightVars const &wv)
Definition: VarVals.h:238
const XML_Char * s
Definition: expat.h:262
float SigmaWeight(uint8_t const &sigmaLevel) const
Definition: VarVals.h:212
SystVar(uint8_t type, float minus2, float minus1, float plus1, float plus2)
Definition: VarVals.h:174
int subrun
Definition: VarVals.h:80
float fTrueParentPDG
true PDG of neutrino parent
Definition: VarVals.h:47
cmf::TruthVars fTruthVars
Definition: VarVals.h:230
float fRPACCQE_Weight
Deprecated.
Definition: VarVals.h:151
cmf::TruthVars const & TruthVars() const
Definition: VarVals.h:261
cmf::SelectionType_t selectionType
Definition: Structs.h:116
float fMinus1Sigma
weight associated with this change in the underlying parameter
Definition: VarVals.h:197
int cycle
Definition: VarVals.h:82
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:10
bool operator<(StanConfig::Verbosity a, StanConfig::Verbosity b)
Allow for comparing them, since kQuiet is definitely "less" verbose than kVerbose.
Definition: StanConfig.h:95
bool IsNuMuSelected() const
Definition: Structs.cxx:226
std::vector< float > SigmasAsVec() const
Definition: VarVals.h:203
float fXSecCV2020_Weight
For 2020 Ana, store weights seperately.
Definition: VarVals.h:158
bool operator==(FC &A, FC &B)
DataVars(float nuecvn, float nuenummichel, float nurece, float hadrece, float leprece, float leprecemcfrac, float recoQ2, float fakewgt)
Definition: VarVals.h:106
int slice
Definition: VarVals.h:83
Base container for the MC related Vars that constitute an event.
Definition: VarVals.h:268
float fXSecCVPPFX_Weight
Deprecated, Was Tufts weight for SA.
Definition: VarVals.h:150
float fTrueParentPZ
true p_Z of neutrino parent
Definition: VarVals.h:49
Definition: run.py:1
std::vector< cmf::SystVar > SystVarColl
Definition: VarVals.h:223
cmf::DataVars fDataVars
Definition: VarVals.h:274
Module to combine a set of results into a single file currently only does one data product type at a ...
Definition: Event.cxx:24
DataVarVals(cmf::DataVars const &dv)
Definition: VarVals.h:279
std::string pars("Th23Dmsq32")
static float RecoEnergy(float leptonicEnergy, float hadronicEnergy, cmf::MetaData const &md)
Definition: VarVals.h:317
cmf::WeightVars const & WeightVars() const
Definition: VarVals.h:262
float fNuE_NumMichel
Definition: VarVals.h:125
enum cmf::beam_type BeamType_t
TRandom3 r(0)
uint8_t Type() const
Definition: VarVals.h:186
cmf::BeamType_t BeamType() const
Definition: Structs.cxx:184
float fOtherDIS_Weight
Deprecated.
Definition: VarVals.h:157
SystVar(uint8_t type)
Definition: VarVals.h:166
float fHad_RecoE
Definition: VarVals.h:127
cmf::SystVarColl const & SystematicShifts() const
Definition: VarVals.h:263
Float_t e
Definition: plot.C:35
float fTrueIntType
true interaction type
Definition: VarVals.h:45
float fTrueIntMode
true interaction mode
Definition: VarVals.h:53
int event
Definition: VarVals.h:81
float fTrueParentTargetPDG
true parent pdg code off the target
Definition: VarVals.h:51
float fLep_RecoE
Definition: VarVals.h:128
bool IsNuESelected() const
Definition: Structs.cxx:232
bool IsNCSelected() const
Definition: Structs.cxx:238
float fPlus1Sigma
weight associated with this change in the underlying parameter
Definition: VarVals.h:198
static float NCEnergy(float leptonicEnergy, float hadronicEnergy, cmf::BeamType_t const &beamType, cmf::DetType_t const &detector)
Definition: VarVals.h:365