BPFEnergyEstimatorOnly_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 //
3 // \file BPFEnergyEstimatorOnly_module.cc
4 // \brief Module to produce the BPF Energy object.
5 // \version
6 // \author Michael Baird - mbaird42@fnal.gov
7 //
8 ////////////////////////////////////////////////////////////////////////
9 
10 // C/C++ includes
11 #include <cmath>
12 #include <iostream>
13 #include <vector>
14 
15 // ROOT includes
16 #include "TFile.h"
17 #include "TH2F.h"
18 #include "TMVA/Factory.h"
19 #include "TMVA/Tools.h"
20 #include "TMVA/Reader.h"
21 #include "TMVA/MethodKNN.h"
22 
23 // Framework includes
33 #include "fhiclcpp/ParameterSet.h"
35 
36 // NOvASoft includes
37 #include "RecoBase/CellHit.h"
38 #include "RecoBase/RecoHit.h"
39 #include "RecoBase/Cluster.h"
40 #include "RecoBase/Prong.h"
41 #include "RecoBase/Track.h"
42 #include "RecoBase/FitSum.h"
43 
45 #include "Utilities/AssociationUtil.h"
46 #include "Geometry/Geometry.h"
47 #include "Calibrator/Calibrator.h"
48 
53 
54 
55 
56 // BPFEnergyEstimatorOnly header
57 namespace bpfit {
59  public:
60  explicit BPFEnergyEstimatorOnly(fhicl::ParameterSet const& pset);
62 
63  void produce(art::Event& evt);
64  void reconfigure(const fhicl::ParameterSet& pset);
65  void beginRun(art::Run& run);
66  void resetVars();
67 
68  private:
70  std::string fProngLabel; ///< label for module that made the prongs
71  std::string fProngInstance; ///< instance label for the prongs to be used
72  std::string fTrackLabel; ///< label for module that made the tracks
73  std::string fBPFPIdLabel; ///< label for module that made the BPFPId object
74  std::string fLibPath; ///< location of all BPF library files listed below
75  std::string fEnXMLFile1; ///< name of the XML file for the first energy estimator
76  std::string fEnXMLFile2; ///< name of the XML file for the second energy estimator
77  std::string fEnXMLFile3; ///< name of the XML file for the third energy estimator
78  double fdxTOL; ///< lower limit for dx values to be used in dE/dx calculation
79 
80 
81  murem::TrackCleanUpAlg fTrackCleanUpAlg; ///< object to determine hadronic energy on the muon track
82 
83  bool fBookedKNNs; ///< Have we yet booked the kNNs?
84 
85  TMVA::Reader fEnergyEstm1; ///< TMVA reader object for the first energy estimator
86  TMVA::Reader fEnergyEstm2; ///< TMVA reader object for the second energy estimator
87  TMVA::Reader fEnergyEstm3; ///< TMVA reader object for the third energy estimator
88 
89 
90 
91  //
92  // Variables used with the kNNs
93  //
94 
95  // muon PID variables:
96  float fLength; ///< track length
97  float fChi2T; ///< total chi^2
98  float fdEdxLL; ///< dEdx LL value
99  float fHitRatio; ///< ratio of track hits to prong hits (used as a track/shower discriminator)
100 
101  // energy estimator variables:
102  // * target variables:
103  float fE1; ///< first estimator true energy
104  float fEsq1; ///< first estimator true energy squared (used to estiamte st. dev.)
105  float fE2; ///< second estimator true energy
106  float fEsq2; ///< second estimator true energy squared (used to estiamte st. dev.)
107  float fE3; ///< third estimator true energy
108  float fEsq3; ///< third estimator true energy squared (used to estiamte st. dev.)
109  // * input variables:
110  float fPmu; ///< reconstructed muon momentum from the track selected to be the muon
111  float fDirZmu; ///< muon track dirZ from the track selected to be the muon
112  float fN3Dprongs; ///< total number of 3D prongs
113  float fE3Dprongs; ///< sum of hit energies in 3D prongs (excluding the "muon" prong) - this does NOT account for dead material
114  float fEremain; ///< sum of remaining slice energy not in 3D prongs - this does NOT account for dead material
115  float fSumPE; ///< sum of total PE for all hits not on the muon track
116 
117  };
118 }
119 
120 
121 
122 namespace bpfit
123 {
124  //.......................................................................
126  fSlicerToken(consumes<std::vector<rb::Cluster>>(pset.get<std::string>("SlicerLabel"))),
127  fTrackCleanUpAlg(pset.get<fhicl::ParameterSet>("TrackCleanUpAlgPSet")),
128  fEnergyEstm1(),
129  fEnergyEstm2(),
130  fEnergyEstm3()
131  {
132  produces< std::vector<BPFEnergy> >();
133  produces< art::Assns<BPFEnergy, rb::Cluster> >();
134 
135  this->reconfigure(pset);
136 
137  fBookedKNNs = false;
138  }
139 
140  //......................................................................
142  {
143  //======================================================================
144  // Clean up any memory allocated by your module... ...if you dare...
145  //======================================================================
146  }
147 
148  //......................................................................
150  {
151  fProngLabel = pset.get<std::string> ("ProngLabel");
152  fProngInstance = pset.get<std::string> ("ProngInstance");
153  fTrackLabel = pset.get<std::string> ("TrackLabel");
154  fBPFPIdLabel = pset.get<std::string> ("BPFPIdLabel");
155  fLibPath = pset.get<std::string> ("LibPath");
156  fEnXMLFile1 = pset.get<std::string> ("EnXMLFile1");
157  fEnXMLFile2 = pset.get<std::string> ("EnXMLFile2");
158  fEnXMLFile3 = pset.get<std::string> ("EnXMLFile3");
159  fdxTOL = pset.get<double> ("dxTOL");
160  }
161 
162  //......................................................................
164  {
165 
166  // Interpret the env variable LibPath
168 
169  if(!fBookedKNNs) {
170 
171  //
172  // Set up the TMVA reader for the energy estimators.
173  //
174  fEnergyEstm1.AddVariable("Pmu",&fPmu);
175  fEnergyEstm1.AddVariable("DirZmu",&fDirZmu);
176  fEnergyEstm1.AddVariable("N3Dprongs",&fN3Dprongs);
177  fEnergyEstm1.AddVariable("E3Dprongs",&fE3Dprongs);
178  fEnergyEstm1.AddVariable("Eremain",&fEremain);
179  fEnergyEstm1.AddVariable("SumPE",&fSumPE);
180  fEnergyEstm1.BookMVA("KNN Method",libpath+fEnXMLFile1);
181 
182  fEnergyEstm2.AddVariable("Pmu",&fPmu);
183  fEnergyEstm2.AddVariable("DirZmu",&fDirZmu);
184  fEnergyEstm2.AddVariable("N3Dprongs",&fN3Dprongs);
185  fEnergyEstm2.AddVariable("E3Dprongs",&fE3Dprongs);
186  fEnergyEstm2.AddVariable("Eremain",&fEremain);
187  fEnergyEstm2.AddVariable("SumPE",&fSumPE);
188  fEnergyEstm2.BookMVA("KNN Method",libpath+fEnXMLFile2);
189 
190  fEnergyEstm3.AddVariable("Pmu",&fPmu);
191  fEnergyEstm3.AddVariable("DirZmu",&fDirZmu);
192  fEnergyEstm3.AddVariable("N3Dprongs",&fN3Dprongs);
193  fEnergyEstm3.AddVariable("E3Dprongs",&fE3Dprongs);
194  fEnergyEstm3.AddVariable("Eremain",&fEremain);
195  fEnergyEstm3.AddVariable("SumPE",&fSumPE);
196  fEnergyEstm3.BookMVA("KNN Method",libpath+fEnXMLFile3);
197 
198  fBookedKNNs = true;
199  }
200 
201  }
202 
203  //......................................................................
205  {
206  //
207  // PID variables:
208  //
209  fLength = -5.0;
210  fChi2T = 1.0e9;
211  fdEdxLL = -1.0e9;
212  fHitRatio = -5.0;
213 
214  //
215  // energy estimator variables:
216  //
217  fE1 = -5.0;
218  fEsq1 = -5.0;
219  fE2 = -5.0;
220  fEsq2 = -5.0;
221  fE3 = -5.0;
222  fEsq3 = -5.0;
223  fPmu = -5.0;
224  fDirZmu = -5.0;
225  fN3Dprongs = -5.0;
226  fE3Dprongs = -5.0;
227  fEremain = -5.0;
228  fSumPE = -5.0;
229  }
230 
231  //......................................................................
233  {
234 
235  // Get the slices and put them into a more convenient container.
237  evt.getByToken(fSlicerToken, slhandle);
239  for(unsigned int i=0; i < slhandle->size(); ++i) {
240  slices.push_back(art::Ptr<rb::Cluster>(slhandle,i));
241  }
242 
243  // Create the FindMany object for getting prongs.
245 
246  // Make the required service handles...
249 
250  // Make the containers for the BPFEnergy objects and their associations
251  std::unique_ptr< std::vector<BPFEnergy> >
252  energies(new std::vector<BPFEnergy>);
253  std::unique_ptr< art::Assns<BPFEnergy,rb::Cluster> >
254  EnSlAssns(new art::Assns<BPFEnergy,rb::Cluster>);
255 
256  // loop over all slices and get the fuzzyk 3D prongs
257  for(unsigned int islice = 0; islice < slhandle->size(); ++islice) {
258 
259  // As usual, skip the noise slice...
260  if((*slhandle)[islice].IsNoise()) continue;
261 
262  // NOTE: Currently, the fundamental object of analysis is the slice.
263  // If in the future we change things so that the fundamental object
264  // of analysis becomes the vertex, then we will have to get prongs
265  // associated with the vertex instead.
266 
267  // Reset all variables.
268  this->resetVars();
269 
270  // get the 3D prongs associated with this slice
271  std::vector< art::Ptr< rb::Prong > > prongs;
272  if(prongFMP.isValid()) prongs = prongFMP.at(islice);
273 
274  // create the FindMany object for getting tracks
275  art::FindManyP<rb::Track> trackFMP(prongs, evt, fTrackLabel);
276 
277 
278 
279  // Variables used to keep track with the best muon PID.
280  float bestPID = -5.0; // best muonPID value
281  int prongPick = -1; // Prong picked as having the longest track.
282 
283 
284 
285  //////////
286  //
287  // Step 1: Determine which prong is the most muon-like.
288  //
289  //////////
290 
291  // loop over all prongs and get the tracks
292  for(unsigned int iprong = 0; iprong < prongs.size(); ++iprong) {
293 
294  // get the tracks associated with this prong
295  std::vector< art::Ptr< rb::Track > > tracks;
296  if(trackFMP.isValid()) tracks = trackFMP.at(iprong);
297 
298  // create the FindMany object for getting FitSum objects
299  art::FindManyP<rb::FitSum> fitsumFMP(tracks, evt, fTrackLabel);
300 
301  art::FindManyP<bpfit::BPFPId> bpfpidFMP(tracks, evt, fBPFPIdLabel);
302  // Loop over all tracks...
303  for(unsigned int itrack = 0; itrack < tracks.size(); ++itrack) {
304 
305  // get the FitSum object associated with this track
306  std::vector< art::Ptr< rb::FitSum > > fitsums;
307  if(fitsumFMP.isValid()) fitsums = fitsumFMP.at(itrack);
308 
309  // get the BPFPId object associated with this track
310  std::vector< art::Ptr< bpfit::BPFPId > > bpfpids;
311  if(bpfpidFMP.isValid()) bpfpids = bpfpidFMP.at(itrack);
312 
313 
314  // There can be only one!!!
315  // (something went really wrong with BPF if there's not one and
316  // only one of these per track.)
317  if(fitsums.size() != 1) abort();
318 
319 
320  // (something went really wrong with BPF / the PID has not been produced
321  if(bpfpids.size() < 1){
322  continue;
323  }
324 
325  // For now, this energy estimator only uses information from the muon track fit.
326  if(fitsums[0]->PDG() != 13) continue;
327 
328  BPFPId bpfpid;
329  int bpfpid_index = -1;
330  for(size_t index=0;index<bpfpids.size();index++){
331  if(bpfpids.at(index)->Pdg()==13) bpfpid_index=index;
332  }
333  if(bpfpid_index==-1) continue;
334  if(bpfpidFMP.isValid()) bpfpid = *bpfpids[bpfpid_index];
335  fLength = bpfpid.GetLength();
336  fChi2T = bpfpid.GetChi2T();
337  fdEdxLL = bpfpid.GetdEdXLL();
338  fHitRatio = bpfpid.GetHitRatio();
339 
340 
341  // Set the variables for the track with the best muonPID
342  if(fLength > 0.0 && fdEdxLL > -1.0e9 && fChi2T < 1.0e9) {
343  // Get the PID value...
344  float pid = bpfpid.Value();
345 
346  if(pid > bestPID) {
347  bestPID = pid;
348  prongPick = (int)iprong;
349  fPmu = fitsums[0]->FourMom().P();
350  fDirZmu = tracks[itrack]->Dir().Z();
351  }
352  }
353 
354  } // end loop over tracks (itrack)
355 
356  } // end loop over prongs (iprong)
357 
358 
359 
360  //////////
361  //
362  // Step 2: Compute the values to be used for the energy estimator.
363  //
364  //////////
365 
366  //
367  // Add up the energy of the non-muon 3D prongs.
368  //
369 
370  float N3Dpr = 0.0;
371  float E3Dpr = 0.0;
372  float sumPE = 0.0;
373  rb::Cluster hits3Dprongs; // a container for all hits on 3D prongs
374 
375  for(unsigned int p = 0; p < prongs.size(); ++p) {
376 
377  // Skip 2D prongs.
378  if(prongs[p]->Is2D()) continue;
379  N3Dpr++;
380 
381  // Keep track of all hits in 3D prongs (to be used later when adding up
382  // the remaining slice energy.)
383  // NOTE: FuzzyK prongs can sometimes share hits, so we have to be careful
384  // not to double count.
385  for(unsigned int h = 0; h < prongs[p]->NCell(); ++h) {
386 
387  const rb::CellHit *pronghit = prongs[p]->Cell(h).get();
388 
389  // Check to see if this hit is already in the list.
390  bool onlist = false;
391  for(unsigned int hpr = 0; hpr < hits3Dprongs.NCell(); ++hpr) {
392  const rb::CellHit *pronghit3D = hits3Dprongs.Cell(hpr).get();
393  if((*pronghit) == (*pronghit3D)) {
394  onlist = true;
395  break;
396  }
397  } // end loop over hits in 3D prongs (hpr)
398 
399  // Add the hit if its not already on the list.
400  if(!onlist) hits3Dprongs.Add(prongs[p]->Cell(h));
401 
402  } // end loop over hits (h)
403 
404  // Skip this prong if it is the identified muon prong.
405  if((int)p == prongPick) continue;
406 
407  // Loop over all hits in this 3D, non-muon prong and add up their energies.
408  //
409  // NOTE: For the hit W position, I am currently using the prong W instead of
410  // the track W. Since I haven't determined which of the three BPF fits
411  // is the best one for this prong, this saves me from having to try
412  // to guess at which track to use. My assumption is that prong->W will
413  // be good enough since I am lumping all hits from 3D prongs together
414  // into one 3D prong energy number...
415  //
416  // WARNING: Since hits can be shared by fuzzyK prongs, there is the potential
417  // below to double (or triple) count energy on shared cell hits. For
418  // now, I am assuming that this is a small effect but something
419  // smarter could/should be put into place in the future.
420  for(unsigned int h = 0; h < prongs[p]->NCell(); ++h) {
421 
422  const rb::CellHit *chit = prongs[p]->Cell(h).get();
423  float W = prongs[p]->W(chit);
424  rb::RecoHit rhit(cal->MakeRecoHit(*chit,W));
425 
426  sumPE += chit->PE();
427  if(rhit.IsCalibrated()) E3Dpr += rhit.GeV();
428 
429  } // end loop over hits (h)
430 
431  } // end loop over prongs (p)
432 
433  // Set the prong variables.
434  if(bestPID >= 0.0) {
435  fN3Dprongs = N3Dpr;
436  fE3Dprongs = E3Dpr;
437  }
438 
439 
440 
441  //
442  // Add up the remaining energy in the slice.
443  //
444 
445  rb::Cluster hitsRemain; // a container for all hits NOT on 3D prongs
446 
447  // Loop over all hits in this slice.
448  for(unsigned int hsl = 0; hsl < (*slhandle)[islice].NCell(); ++hsl) {
449 
450  const rb::CellHit *slicehit = (*slhandle)[islice].Cell(hsl).get();
451 
452  // Check to see if this hit is in the list of hits on 3D prongs.
453  bool onlist = false;
454  for(unsigned int hpr = 0; hpr < hits3Dprongs.NCell(); ++hpr) {
455  const rb::CellHit *pronghit = hits3Dprongs.Cell(hpr).get();
456  if((*pronghit) == (*slicehit)) {
457  onlist = true;
458  break;
459  }
460  } // end loop over hits in 3D prongs (hpr)
461 
462  // If this hit is not in a 3D prong, add it to the list of remaining slice hits.
463  if(!onlist) {
464  hitsRemain.Add((*slhandle)[islice].Cell(hsl));
465  sumPE += slicehit->PE();
466  }
467 
468  } // end loop over hits in this slice (hsl)
469 
470  // Set the remaining slice energy variable.
471  if(bestPID >= 0.0) {
472  if(hitsRemain.NCell() > 0) fEremain = hitsRemain.TotalGeV();
473  else fEremain = 0.0;
474  }
475 
476  if(bestPID >= 0.0) {
477  fSumPE = sumPE;
478  }
479 
480 
481 
482  //
483  // Compute the hadronic contamination on the muon track and add that back in.
484  //
485  float muonTrackHadE = 0.0;
486 
487  // Only do this if we actually picked a track to be the muon.
488  if(prongPick >= 0) {
489 
490  std::vector< art::Ptr< rb::Track > > tracks;
491  if(trackFMP.isValid()) tracks = trackFMP.at(prongPick);
492  art::FindManyP<rb::FitSum> fitsumFMP(tracks, evt, fTrackLabel);
493 
494  // Figure out which BPF track was made with the muon assumption.
495  for(unsigned int itrack = 0; itrack < tracks.size(); ++itrack) {
496 
497  std::vector< art::Ptr< rb::FitSum > > fitsums;
498  if(fitsumFMP.isValid()) fitsums = fitsumFMP.at(itrack);
499  if(fitsums.size() != 1) abort();
500 
501  // Skip if not the muon track...
502  if(fitsums[0]->PDG() != 13) continue;
503 
504  // If it is the muon track, do the extra Ehad calculation
505  muonTrackHadE = fTrackCleanUpAlg.ExtraEOnTrackInGeV(*tracks[itrack],*slices[islice]);
506 
507  } // end loop over tracks (itrack)
508 
509  } // end if(prongPick >= 0)
510 
511  fEremain += muonTrackHadE;
512 
513 
514 
515  //
516  // Set the BPFEnergy object values and make the associations.
517  //
518  float Eres1 = -5.0;
519  float Eres2 = -5.0;
520  float Eres3 = -5.0;
521  if(bestPID >= 0.0) {
522  fE1 = (fEnergyEstm1.EvaluateRegression("KNN Method"))[0];
523  fEsq1 = (fEnergyEstm1.EvaluateRegression("KNN Method"))[1];
524  Eres1 = sqrt(fabs(fEsq1 - fE1*fE1));
525 
526  fE2 = (fEnergyEstm2.EvaluateRegression("KNN Method"))[0];
527  fEsq2 = (fEnergyEstm2.EvaluateRegression("KNN Method"))[1];
528  Eres2 = sqrt(fabs(fEsq2 - fE2*fE2));
529 
530  fE3 = (fEnergyEstm3.EvaluateRegression("KNN Method"))[0];
531  fEsq3 = (fEnergyEstm3.EvaluateRegression("KNN Method"))[1];
532  Eres3 = sqrt(fabs(fEsq3 - fE3*fE3));
533  }
534 
535  BPFEnergy energy(fE1,Eres1,fE2,Eres2,fE3,Eres3,
536  bestPID,fPmu,fDirZmu,
538  energies->push_back(energy);
539  util::CreateAssn(*this,evt,*(energies.get()),slices[islice],*(EnSlAssns.get()));
540 
541  } // end loop over slhandle (islice)
542 
543  evt.put(std::move(energies));
544  evt.put(std::move(EnSlAssns));
545 
546  }
547 
548 } // end namespace bpfit
549 ////////////////////////////////////////////////////////////////////////
550 
551 
552 
554 
555 namespace bpfit
556 {
558 }
float fPmu
reconstructed muon momentum from the track selected to be the muon
bool fBookedKNNs
Have we yet booked the kNNs?
static bool CreateAssn(art::EDProducer const &prod, art::Event &evt, std::vector< T > &a, art::Ptr< U > b, art::Assns< T, U > &assn, size_t indx=UINT_MAX, std::string const &instance=std::string())
Create a 1 to 1 association between a new product and one already in the event.
unsigned int NCell(geo::View_t view) const
Number of cells in view view.
Definition: Cluster.cxx:134
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
std::string fTrackLabel
label for module that made the tracks
TMVA::Reader fEnergyEstm3
TMVA reader object for the third energy estimator.
std::string fLibPath
location of all BPF library files listed below
float fEsq2
second estimator true energy squared (used to estiamte st. dev.)
const char * p
Definition: xmltok.h:285
T sqrt(T number)
Definition: d0nt_math.hpp:156
std::string EnvExpansion(const std::string &inString)
Function to expand environment variables.
Definition: EnvExpand.cxx:8
rb::RecoHit MakeRecoHit(rb::CellHit const &cellhit, double w)
A collection of associated CellHits.
Definition: Cluster.h:47
float fN3Dprongs
total number of 3D prongs
TMVA::Reader fEnergyEstm1
TMVA reader object for the first energy estimator.
void reconfigure(const fhicl::ParameterSet &pset)
DEFINE_ART_MODULE(TestTMapFile)
float fSumPE
sum of total PE for all hits not on the muon track
float fE3Dprongs
sum of hit energies in 3D prongs (excluding the "muon" prong) - this does NOT account for dead materi...
float GetdEdXLL() const
Definition: BPFPId.h:44
double Value() const
Definition: PID.h:22
Definition: Run.h:31
Calibrated quantities relying on position in the orthogonal view. To generate a rb::CellHit from a rb...
Definition: RecoHit.h:19
float fE2
second estimator true energy
virtual void Add(const art::Ptr< rb::CellHit > &cell, double weight=1)
Definition: Cluster.cxx:84
Calculates dE/dx in cells passed through by a track.
ProductID put(std::unique_ptr< PROD > &&product)
Definition: Event.h:102
float fDirZmu
muon track dirZ from the track selected to be the muon
std::string fProngInstance
instance label for the prongs to be used
unsigned short Cell() const
Definition: CellHit.h:40
float GetLength() const
Definition: BPFPId.h:42
murem::TrackCleanUpAlg fTrackCleanUpAlg
object to determine hadronic energy on the muon track
std::string fProngLabel
label for module that made the prongs
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:441
std::string fEnXMLFile3
name of the XML file for the third energy estimator
T get(std::string const &key) const
Definition: ParameterSet.h:231
BPFEnergyEstimatorOnly(fhicl::ParameterSet const &pset)
double energy
Definition: plottest35.C:25
float PE() const
Definition: CellHit.h:42
float GetHitRatio() const
Definition: BPFPId.h:45
TMVA::Reader fEnergyEstm2
TMVA reader object for the second energy estimator.
Perform a "2 point" Hough transform on a collection of hits.
Definition: run.py:1
float fEremain
sum of remaining slice energy not in 3D prongs - this does NOT account for dead material ...
const art::ProductToken< std::vector< rb::Cluster > > fSlicerToken
art::Ptr< rb::CellHit > Cell(geo::View_t view, unsigned int viewIdx) const
Get the ith cell from view view.
Definition: Cluster.cxx:145
float fE1
first estimator true energy
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
A rawdata::RawDigit with channel information decoded.
Definition: CellHit.h:27
std::string fEnXMLFile2
name of the XML file for the second energy estimator
T const * get() const
Definition: Ptr.h:321
std::string fBPFPIdLabel
label for module that made the BPFPId object
float fHitRatio
ratio of track hits to prong hits (used as a track/shower discriminator)
void geom(int which=0)
Definition: geom.C:163
float fEsq1
first estimator true energy squared (used to estiamte st. dev.)
double TotalGeV(EEnergyCalcScheme escheme=kRecomputeEnergy) const
Simple sum of the estimated GeV of all the hits.
Definition: Cluster.cxx:378
float fEsq3
third estimator true energy squared (used to estiamte st. dev.)
float fE3
third estimator true energy
std::string fEnXMLFile1
name of the XML file for the first energy estimator
bool getByToken(ProductToken< PROD > const &token, Handle< PROD > &result) const
Definition: DataViewImpl.h:387
double fdxTOL
lower limit for dx values to be used in dE/dx calculation
float GetChi2T() const
Definition: BPFPId.h:43
#define W(x)
float ExtraEOnTrackInGeV(const rb::Track &track, const rb::Cluster &slice, double *mipRangeHigh=NULL, double *mipRangeLow=NULL, double *mipValue=NULL, double *vertexRegionDeDxCutOff=NULL)
ProductToken< T > consumes(InputTag const &)
Encapsulate the geometry of one entire detector (near, far, ndos)