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
28 #include "art_root_io/TFileDirectory.h"
29 #include "art_root_io/TFileService.h"
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  EDProducer(pset),
127  fSlicerToken(consumes<std::vector<rb::Cluster>>(pset.get<std::string>("SlicerLabel"))),
128  fTrackCleanUpAlg(pset.get<fhicl::ParameterSet>("TrackCleanUpAlgPSet")),
129  fEnergyEstm1(),
130  fEnergyEstm2(),
131  fEnergyEstm3()
132  {
133  produces< std::vector<BPFEnergy> >();
134  produces< art::Assns<BPFEnergy, rb::Cluster> >();
135 
136  this->reconfigure(pset);
137 
138  fBookedKNNs = false;
139  }
140 
141  //......................................................................
143  {
144  //======================================================================
145  // Clean up any memory allocated by your module... ...if you dare...
146  //======================================================================
147  }
148 
149  //......................................................................
151  {
152  fProngLabel = pset.get<std::string> ("ProngLabel");
153  fProngInstance = pset.get<std::string> ("ProngInstance");
154  fTrackLabel = pset.get<std::string> ("TrackLabel");
155  fBPFPIdLabel = pset.get<std::string> ("BPFPIdLabel");
156  fLibPath = pset.get<std::string> ("LibPath");
157  fEnXMLFile1 = pset.get<std::string> ("EnXMLFile1");
158  fEnXMLFile2 = pset.get<std::string> ("EnXMLFile2");
159  fEnXMLFile3 = pset.get<std::string> ("EnXMLFile3");
160  fdxTOL = pset.get<double> ("dxTOL");
161  }
162 
163  //......................................................................
165  {
166 
167  // Interpret the env variable LibPath
169 
170  if(!fBookedKNNs) {
171 
172  //
173  // Set up the TMVA reader for the energy estimators.
174  //
175  fEnergyEstm1.AddVariable("Pmu",&fPmu);
176  fEnergyEstm1.AddVariable("DirZmu",&fDirZmu);
177  fEnergyEstm1.AddVariable("N3Dprongs",&fN3Dprongs);
178  fEnergyEstm1.AddVariable("E3Dprongs",&fE3Dprongs);
179  fEnergyEstm1.AddVariable("Eremain",&fEremain);
180  fEnergyEstm1.AddVariable("SumPE",&fSumPE);
181  fEnergyEstm1.BookMVA("KNN Method",libpath+fEnXMLFile1);
182 
183  fEnergyEstm2.AddVariable("Pmu",&fPmu);
184  fEnergyEstm2.AddVariable("DirZmu",&fDirZmu);
185  fEnergyEstm2.AddVariable("N3Dprongs",&fN3Dprongs);
186  fEnergyEstm2.AddVariable("E3Dprongs",&fE3Dprongs);
187  fEnergyEstm2.AddVariable("Eremain",&fEremain);
188  fEnergyEstm2.AddVariable("SumPE",&fSumPE);
189  fEnergyEstm2.BookMVA("KNN Method",libpath+fEnXMLFile2);
190 
191  fEnergyEstm3.AddVariable("Pmu",&fPmu);
192  fEnergyEstm3.AddVariable("DirZmu",&fDirZmu);
193  fEnergyEstm3.AddVariable("N3Dprongs",&fN3Dprongs);
194  fEnergyEstm3.AddVariable("E3Dprongs",&fE3Dprongs);
195  fEnergyEstm3.AddVariable("Eremain",&fEremain);
196  fEnergyEstm3.AddVariable("SumPE",&fSumPE);
197  fEnergyEstm3.BookMVA("KNN Method",libpath+fEnXMLFile3);
198 
199  fBookedKNNs = true;
200  }
201 
202  }
203 
204  //......................................................................
206  {
207  //
208  // PID variables:
209  //
210  fLength = -5.0;
211  fChi2T = 1.0e9;
212  fdEdxLL = -1.0e9;
213  fHitRatio = -5.0;
214 
215  //
216  // energy estimator variables:
217  //
218  fE1 = -5.0;
219  fEsq1 = -5.0;
220  fE2 = -5.0;
221  fEsq2 = -5.0;
222  fE3 = -5.0;
223  fEsq3 = -5.0;
224  fPmu = -5.0;
225  fDirZmu = -5.0;
226  fN3Dprongs = -5.0;
227  fE3Dprongs = -5.0;
228  fEremain = -5.0;
229  fSumPE = -5.0;
230  }
231 
232  //......................................................................
234  {
235 
236  // Get the slices and put them into a more convenient container.
238  evt.getByToken(fSlicerToken, slhandle);
240  for(unsigned int i=0; i < slhandle->size(); ++i) {
241  slices.push_back(art::Ptr<rb::Cluster>(slhandle,i));
242  }
243 
244  // Create the FindMany object for getting prongs.
246 
247  // Make the required service handles...
250 
251  // Make the containers for the BPFEnergy objects and their associations
252  std::unique_ptr< std::vector<BPFEnergy> >
253  energies(new std::vector<BPFEnergy>);
254  std::unique_ptr< art::Assns<BPFEnergy,rb::Cluster> >
255  EnSlAssns(new art::Assns<BPFEnergy,rb::Cluster>);
256 
257  // loop over all slices and get the fuzzyk 3D prongs
258  for(unsigned int islice = 0; islice < slhandle->size(); ++islice) {
259 
260  // As usual, skip the noise slice...
261  if((*slhandle)[islice].IsNoise()) continue;
262 
263  // NOTE: Currently, the fundamental object of analysis is the slice.
264  // If in the future we change things so that the fundamental object
265  // of analysis becomes the vertex, then we will have to get prongs
266  // associated with the vertex instead.
267 
268  // Reset all variables.
269  this->resetVars();
270 
271  // get the 3D prongs associated with this slice
272  std::vector< art::Ptr< rb::Prong > > prongs;
273  if(prongFMP.isValid()) prongs = prongFMP.at(islice);
274 
275  // create the FindMany object for getting tracks
276  art::FindManyP<rb::Track> trackFMP(prongs, evt, fTrackLabel);
277 
278 
279 
280  // Variables used to keep track with the best muon PID.
281  float bestPID = -5.0; // best muonPID value
282  int prongPick = -1; // Prong picked as having the longest track.
283 
284 
285 
286  //////////
287  //
288  // Step 1: Determine which prong is the most muon-like.
289  //
290  //////////
291 
292  // loop over all prongs and get the tracks
293  for(unsigned int iprong = 0; iprong < prongs.size(); ++iprong) {
294 
295  // get the tracks associated with this prong
296  std::vector< art::Ptr< rb::Track > > tracks;
297  if(trackFMP.isValid()) tracks = trackFMP.at(iprong);
298 
299  // create the FindMany object for getting FitSum objects
300  art::FindManyP<rb::FitSum> fitsumFMP(tracks, evt, fTrackLabel);
301 
302  art::FindManyP<bpfit::BPFPId> bpfpidFMP(tracks, evt, fBPFPIdLabel);
303  // Loop over all tracks...
304  for(unsigned int itrack = 0; itrack < tracks.size(); ++itrack) {
305 
306  // get the FitSum object associated with this track
307  std::vector< art::Ptr< rb::FitSum > > fitsums;
308  if(fitsumFMP.isValid()) fitsums = fitsumFMP.at(itrack);
309 
310  // get the BPFPId object associated with this track
311  std::vector< art::Ptr< bpfit::BPFPId > > bpfpids;
312  if(bpfpidFMP.isValid()) bpfpids = bpfpidFMP.at(itrack);
313 
314 
315  // There can be only one!!!
316  // (something went really wrong with BPF if there's not one and
317  // only one of these per track.)
318  if(fitsums.size() != 1) abort();
319 
320 
321  // (something went really wrong with BPF / the PID has not been produced
322  if(bpfpids.size() < 1){
323  continue;
324  }
325 
326  // For now, this energy estimator only uses information from the muon track fit.
327  if(fitsums[0]->PDG() != 13) continue;
328 
329  BPFPId bpfpid;
330  int bpfpid_index = -1;
331  for(size_t index=0;index<bpfpids.size();index++){
332  if(bpfpids.at(index)->Pdg()==13) bpfpid_index=index;
333  }
334  if(bpfpid_index==-1) continue;
335  if(bpfpidFMP.isValid()) bpfpid = *bpfpids[bpfpid_index];
336  fLength = bpfpid.GetLength();
337  fChi2T = bpfpid.GetChi2T();
338  fdEdxLL = bpfpid.GetdEdXLL();
339  fHitRatio = bpfpid.GetHitRatio();
340 
341 
342  // Set the variables for the track with the best muonPID
343  if(fLength > 0.0 && fdEdxLL > -1.0e9 && fChi2T < 1.0e9) {
344  // Get the PID value...
345  float pid = bpfpid.Value();
346 
347  if(pid > bestPID) {
348  bestPID = pid;
349  prongPick = (int)iprong;
350  fPmu = fitsums[0]->FourMom().P();
351  fDirZmu = tracks[itrack]->Dir().Z();
352  }
353  }
354 
355  } // end loop over tracks (itrack)
356 
357  } // end loop over prongs (iprong)
358 
359 
360 
361  //////////
362  //
363  // Step 2: Compute the values to be used for the energy estimator.
364  //
365  //////////
366 
367  //
368  // Add up the energy of the non-muon 3D prongs.
369  //
370 
371  float N3Dpr = 0.0;
372  float E3Dpr = 0.0;
373  float sumPE = 0.0;
374  rb::Cluster hits3Dprongs; // a container for all hits on 3D prongs
375 
376  for(unsigned int p = 0; p < prongs.size(); ++p) {
377 
378  // Skip 2D prongs.
379  if(prongs[p]->Is2D()) continue;
380  N3Dpr++;
381 
382  // Keep track of all hits in 3D prongs (to be used later when adding up
383  // the remaining slice energy.)
384  // NOTE: FuzzyK prongs can sometimes share hits, so we have to be careful
385  // not to double count.
386  for(unsigned int h = 0; h < prongs[p]->NCell(); ++h) {
387 
388  const rb::CellHit *pronghit = prongs[p]->Cell(h).get();
389 
390  // Check to see if this hit is already in the list.
391  bool onlist = false;
392  for(unsigned int hpr = 0; hpr < hits3Dprongs.NCell(); ++hpr) {
393  const rb::CellHit *pronghit3D = hits3Dprongs.Cell(hpr).get();
394  if((*pronghit) == (*pronghit3D)) {
395  onlist = true;
396  break;
397  }
398  } // end loop over hits in 3D prongs (hpr)
399 
400  // Add the hit if its not already on the list.
401  if(!onlist) hits3Dprongs.Add(prongs[p]->Cell(h));
402 
403  } // end loop over hits (h)
404 
405  // Skip this prong if it is the identified muon prong.
406  if((int)p == prongPick) continue;
407 
408  // Loop over all hits in this 3D, non-muon prong and add up their energies.
409  //
410  // NOTE: For the hit W position, I am currently using the prong W instead of
411  // the track W. Since I haven't determined which of the three BPF fits
412  // is the best one for this prong, this saves me from having to try
413  // to guess at which track to use. My assumption is that prong->W will
414  // be good enough since I am lumping all hits from 3D prongs together
415  // into one 3D prong energy number...
416  //
417  // WARNING: Since hits can be shared by fuzzyK prongs, there is the potential
418  // below to double (or triple) count energy on shared cell hits. For
419  // now, I am assuming that this is a small effect but something
420  // smarter could/should be put into place in the future.
421  for(unsigned int h = 0; h < prongs[p]->NCell(); ++h) {
422 
423  const rb::CellHit *chit = prongs[p]->Cell(h).get();
424  float W = prongs[p]->W(chit);
425  rb::RecoHit rhit(cal->MakeRecoHit(*chit,W));
426 
427  sumPE += chit->PE();
428  if(rhit.IsCalibrated()) E3Dpr += rhit.GeV();
429 
430  } // end loop over hits (h)
431 
432  } // end loop over prongs (p)
433 
434  // Set the prong variables.
435  if(bestPID >= 0.0) {
436  fN3Dprongs = N3Dpr;
437  fE3Dprongs = E3Dpr;
438  }
439 
440 
441 
442  //
443  // Add up the remaining energy in the slice.
444  //
445 
446  rb::Cluster hitsRemain; // a container for all hits NOT on 3D prongs
447 
448  // Loop over all hits in this slice.
449  for(unsigned int hsl = 0; hsl < (*slhandle)[islice].NCell(); ++hsl) {
450 
451  const rb::CellHit *slicehit = (*slhandle)[islice].Cell(hsl).get();
452 
453  // Check to see if this hit is in the list of hits on 3D prongs.
454  bool onlist = false;
455  for(unsigned int hpr = 0; hpr < hits3Dprongs.NCell(); ++hpr) {
456  const rb::CellHit *pronghit = hits3Dprongs.Cell(hpr).get();
457  if((*pronghit) == (*slicehit)) {
458  onlist = true;
459  break;
460  }
461  } // end loop over hits in 3D prongs (hpr)
462 
463  // If this hit is not in a 3D prong, add it to the list of remaining slice hits.
464  if(!onlist) {
465  hitsRemain.Add((*slhandle)[islice].Cell(hsl));
466  sumPE += slicehit->PE();
467  }
468 
469  } // end loop over hits in this slice (hsl)
470 
471  // Set the remaining slice energy variable.
472  if(bestPID >= 0.0) {
473  if(hitsRemain.NCell() > 0) fEremain = hitsRemain.TotalGeV();
474  else fEremain = 0.0;
475  }
476 
477  if(bestPID >= 0.0) {
478  fSumPE = sumPE;
479  }
480 
481 
482 
483  //
484  // Compute the hadronic contamination on the muon track and add that back in.
485  //
486  float muonTrackHadE = 0.0;
487 
488  // Only do this if we actually picked a track to be the muon.
489  if(prongPick >= 0) {
490 
491  std::vector< art::Ptr< rb::Track > > tracks;
492  if(trackFMP.isValid()) tracks = trackFMP.at(prongPick);
493  art::FindManyP<rb::FitSum> fitsumFMP(tracks, evt, fTrackLabel);
494 
495  // Figure out which BPF track was made with the muon assumption.
496  for(unsigned int itrack = 0; itrack < tracks.size(); ++itrack) {
497 
498  std::vector< art::Ptr< rb::FitSum > > fitsums;
499  if(fitsumFMP.isValid()) fitsums = fitsumFMP.at(itrack);
500  if(fitsums.size() != 1) abort();
501 
502  // Skip if not the muon track...
503  if(fitsums[0]->PDG() != 13) continue;
504 
505  // If it is the muon track, do the extra Ehad calculation
506  muonTrackHadE = fTrackCleanUpAlg.ExtraEOnTrackInGeV(*tracks[itrack],*slices[islice]);
507 
508  } // end loop over tracks (itrack)
509 
510  } // end if(prongPick >= 0)
511 
512  fEremain += muonTrackHadE;
513 
514 
515 
516  //
517  // Set the BPFEnergy object values and make the associations.
518  //
519  float Eres1 = -5.0;
520  float Eres2 = -5.0;
521  float Eres3 = -5.0;
522  if(bestPID >= 0.0) {
523  fE1 = (fEnergyEstm1.EvaluateRegression("KNN Method"))[0];
524  fEsq1 = (fEnergyEstm1.EvaluateRegression("KNN Method"))[1];
525  Eres1 = sqrt(fabs(fEsq1 - fE1*fE1));
526 
527  fE2 = (fEnergyEstm2.EvaluateRegression("KNN Method"))[0];
528  fEsq2 = (fEnergyEstm2.EvaluateRegression("KNN Method"))[1];
529  Eres2 = sqrt(fabs(fEsq2 - fE2*fE2));
530 
531  fE3 = (fEnergyEstm3.EvaluateRegression("KNN Method"))[0];
532  fEsq3 = (fEnergyEstm3.EvaluateRegression("KNN Method"))[1];
533  Eres3 = sqrt(fabs(fEsq3 - fE3*fE3));
534  }
535 
536  BPFEnergy energy(fE1,Eres1,fE2,Eres2,fE3,Eres3,
537  bestPID,fPmu,fDirZmu,
539  energies->push_back(energy);
540  util::CreateAssn(evt,*(energies.get()),slices[islice],*(EnSlAssns.get()));
541 
542  } // end loop over slhandle (islice)
543 
544  evt.put(std::move(energies));
545  evt.put(std::move(EnSlAssns));
546 
547  }
548 
549 } // end namespace bpfit
550 ////////////////////////////////////////////////////////////////////////
551 
552 
553 
555 
556 namespace bpfit
557 {
559 }
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
bool getByToken(ProductToken< PROD > const &, Handle< PROD > &result) const
Definition: DataViewImpl.h:462
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
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.h:20
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:21
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.
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:435
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)
int evt
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
A rawdata::RawDigit with channel information decoded.
Definition: CellHit.h:27
std::string fEnXMLFile2
name of the XML file for the second energy estimator
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.)
ProductToken< T > consumes(InputTag const &)
Definition: ModuleBase.h:55
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
double fdxTOL
lower limit for dx values to be used in dE/dx calculation
T const * get() const
Definition: Ptr.h:149
ProductID put(std::unique_ptr< PROD > &&edp, FullSemantic< Level::Run > const semantic)
Definition: DataViewImpl.h:730
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)
Encapsulate the geometry of one entire detector (near, far, ndos)
enum BeamMode string