BPFEnergyEstimator_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 //
3 // \file BPFEnergyEstimator_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 
52 
53 
54 
55 // BPFEnergyEstimator header
56 namespace bpfit {
58  public:
59  explicit BPFEnergyEstimator(fhicl::ParameterSet const& pset);
61 
62  void produce(art::Event& evt);
63  void reconfigure(const fhicl::ParameterSet& pset);
64  void beginRun(art::Run& run);
65  void resetVars();
66 
67  private:
69  std::string fProngLabel; ///< label for module that made the prongs
70  std::string fProngInstance; ///< instance label for the prongs to be used
71  std::string fTrackLabel; ///< label for module that made the tracks
72  std::string fLibPath; ///< location of all BPF library files listed below
73  std::string fHistoFile; ///< name of dE/dx log-likelihood histogram file
74  std::string fPIDXMLFile; ///< name of the XML file for the muon PID
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  TH2F *fdEdxVSbg[2][2]; ///< histo of computed dE/dx for cell hits (vs. beta*gamma)
81  bpfit::dEdxCalculator fdEdxCalc; ///< helper class for computing dEdx values
82 
83  murem::TrackCleanUpAlg fTrackCleanUpAlg; ///< object to determine hadronic energy on the muon track
84 
85  bool fBookedKNNs; ///< Have we yet booked the kNNs?
86 
87  TMVA::Reader fMuonPID; ///< TMVA reader object for the muon PID
88  TMVA::Reader fEnergyEstm1; ///< TMVA reader object for the first energy estimator
89  TMVA::Reader fEnergyEstm2; ///< TMVA reader object for the second energy estimator
90  TMVA::Reader fEnergyEstm3; ///< TMVA reader object for the third energy estimator
91 
92 
93 
94  //
95  // Variables used with the kNNs
96  //
97 
98  // muon PID variables:
99  float fLength; ///< track length
100  float fChi2T; ///< total chi^2
101  float fdEdxLL; ///< dEdx LL value
102  float fHitRatio; ///< ratio of track hits to prong hits (used as a track/shower discriminator)
103 
104  // energy estimator variables:
105  // * target variables:
106  float fE1; ///< first estimator true energy
107  float fEsq1; ///< first estimator true energy squared (used to estiamte st. dev.)
108  float fE2; ///< second estimator true energy
109  float fEsq2; ///< second estimator true energy squared (used to estiamte st. dev.)
110  float fE3; ///< third estimator true energy
111  float fEsq3; ///< third estimator true energy squared (used to estiamte st. dev.)
112  // * input variables:
113  float fPmu; ///< reconstructed muon momentum from the track selected to be the muon
114  float fDirZmu; ///< muon track dirZ from the track selected to be the muon
115  float fN3Dprongs; ///< total number of 3D prongs
116  float fE3Dprongs; ///< sum of hit energies in 3D prongs (excluding the "muon" prong) - this does NOT account for dead material
117  float fEremain; ///< sum of remaining slice energy not in 3D prongs - this does NOT account for dead material
118  float fSumPE; ///< sum of total PE for all hits not on the muon track
119 
120  };
121 }
122 
123 
124 
125 namespace bpfit
126 {
127  //.......................................................................
129  fSlicerToken(consumes<std::vector<rb::Cluster>>(pset.get<std::string>("SlicerLabel"))),
130  fdEdxCalc(),
131  fTrackCleanUpAlg(pset.get<fhicl::ParameterSet>("TrackCleanUpAlgPSet")),
132  fMuonPID(),
133  fEnergyEstm1(),
134  fEnergyEstm2(),
135  fEnergyEstm3()
136  {
137  produces< std::vector<BPFEnergy> >();
138  produces< art::Assns<BPFEnergy, rb::Cluster> >();
139 
140  this->reconfigure(pset);
141 
142  fBookedKNNs = false;
143  }
144 
145  //......................................................................
147  {
148  //======================================================================
149  // Clean up any memory allocated by your module... ...if you dare...
150  //======================================================================
151  }
152 
153  //......................................................................
155  {
156  fProngLabel = pset.get<std::string> ("ProngLabel");
157  fProngInstance = pset.get<std::string> ("ProngInstance");
158  fTrackLabel = pset.get<std::string> ("TrackLabel");
159  fLibPath = pset.get<std::string> ("LibPath");
160  fHistoFile = pset.get<std::string> ("HistoFile");
161  fPIDXMLFile = pset.get<std::string> ("PIDXMLFile");
162  fEnXMLFile1 = pset.get<std::string> ("EnXMLFile1");
163  fEnXMLFile2 = pset.get<std::string> ("EnXMLFile2");
164  fEnXMLFile3 = pset.get<std::string> ("EnXMLFile3");
165  fdxTOL = pset.get<double> ("dxTOL");
166  }
167 
168  //......................................................................
170  {
171 
172  // Interpret the env variable LibPath
174 
175  //
176  // Read in the log-likelihood histos from the input file
177  //
178  const std::string filename = libpath+fHistoFile;
179  TFile infile(filename.c_str());
180 
181  fdEdxVSbg[0][0] = (TH2F*)infile.FindObjectAny("fdEdxVSbg_W0_X0");
182  fdEdxVSbg[0][1] = (TH2F*)infile.FindObjectAny("fdEdxVSbg_W0_X1");
183  fdEdxVSbg[1][0] = (TH2F*)infile.FindObjectAny("fdEdxVSbg_W1_X0");
184  fdEdxVSbg[1][1] = (TH2F*)infile.FindObjectAny("fdEdxVSbg_W1_X1");
185 
186  infile.Close();
187 
188 
189 
190  if(!fBookedKNNs) {
191  //
192  // Set up the TMVA reader for the muon PID.
193  //
194  fMuonPID.AddVariable("length",&fLength);
195  fMuonPID.AddVariable("chi2T",&fChi2T);
196  fMuonPID.AddVariable("dEdxLL",&fdEdxLL);
197  fMuonPID.AddVariable("hitRatio",&fHitRatio);
198  fMuonPID.BookMVA("KNN Method",libpath+fPIDXMLFile);
199 
200  //
201  // Set up the TMVA reader for the energy estimators.
202  //
203  fEnergyEstm1.AddVariable("Pmu",&fPmu);
204  fEnergyEstm1.AddVariable("DirZmu",&fDirZmu);
205  fEnergyEstm1.AddVariable("N3Dprongs",&fN3Dprongs);
206  fEnergyEstm1.AddVariable("E3Dprongs",&fE3Dprongs);
207  fEnergyEstm1.AddVariable("Eremain",&fEremain);
208  fEnergyEstm1.AddVariable("SumPE",&fSumPE);
209  fEnergyEstm1.BookMVA("KNN Method",libpath+fEnXMLFile1);
210 
211  fEnergyEstm2.AddVariable("Pmu",&fPmu);
212  fEnergyEstm2.AddVariable("DirZmu",&fDirZmu);
213  fEnergyEstm2.AddVariable("N3Dprongs",&fN3Dprongs);
214  fEnergyEstm2.AddVariable("E3Dprongs",&fE3Dprongs);
215  fEnergyEstm2.AddVariable("Eremain",&fEremain);
216  fEnergyEstm2.AddVariable("SumPE",&fSumPE);
217  fEnergyEstm2.BookMVA("KNN Method",libpath+fEnXMLFile2);
218 
219  fEnergyEstm3.AddVariable("Pmu",&fPmu);
220  fEnergyEstm3.AddVariable("DirZmu",&fDirZmu);
221  fEnergyEstm3.AddVariable("N3Dprongs",&fN3Dprongs);
222  fEnergyEstm3.AddVariable("E3Dprongs",&fE3Dprongs);
223  fEnergyEstm3.AddVariable("Eremain",&fEremain);
224  fEnergyEstm3.AddVariable("SumPE",&fSumPE);
225  fEnergyEstm3.BookMVA("KNN Method",libpath+fEnXMLFile3);
226 
227  fBookedKNNs = true;
228  }
229 
230  }
231 
232  //......................................................................
234  {
235  //
236  // PID variables:
237  //
238  fLength = -5.0;
239  fChi2T = 1.0e9;
240  fdEdxLL = -1.0e9;
241  fHitRatio = -5.0;
242 
243  //
244  // energy estimator variables:
245  //
246  fE1 = -5.0;
247  fEsq1 = -5.0;
248  fE2 = -5.0;
249  fEsq2 = -5.0;
250  fE3 = -5.0;
251  fEsq3 = -5.0;
252  fPmu = -5.0;
253  fDirZmu = -5.0;
254  fN3Dprongs = -5.0;
255  fE3Dprongs = -5.0;
256  fEremain = -5.0;
257  fSumPE = -5.0;
258  }
259 
260  //......................................................................
262  {
263 
264  // Get the slices and put them into a more convenient container.
266  evt.getByToken(fSlicerToken, slhandle);
268  for(unsigned int i=0; i < slhandle->size(); ++i) {
269  slices.push_back(art::Ptr<rb::Cluster>(slhandle,i));
270  }
271 
272  // Create the FindMany object for getting prongs.
274 
275  // Make the required service handles...
278 
279  // Make the containers for the BPFEnergy objects and their associations
280  std::unique_ptr< std::vector<BPFEnergy> >
281  energies(new std::vector<BPFEnergy>);
282  std::unique_ptr< art::Assns<BPFEnergy,rb::Cluster> >
283  EnSlAssns(new art::Assns<BPFEnergy,rb::Cluster>);
284 
285  // loop over all slices and get the fuzzyk 3D prongs
286  for(unsigned int islice = 0; islice < slhandle->size(); ++islice) {
287 
288  // As usual, skip the noise slice...
289  if((*slhandle)[islice].IsNoise()) continue;
290 
291  // NOTE: Currently, the fundamental object of analysis is the slice.
292  // If in the future we change things so that the fundamental object
293  // of analysis becomes the vertex, then we will have to get prongs
294  // associated with the vertex instead.
295 
296  // Reset all variables.
297  this->resetVars();
298 
299  // get the 3D prongs associated with this slice
300  std::vector< art::Ptr< rb::Prong > > prongs;
301  if(prongFMP.isValid()) prongs = prongFMP.at(islice);
302 
303  // create the FindMany object for getting tracks
304  art::FindManyP<rb::Track> trackFMP(prongs, evt, fTrackLabel);
305 
306 
307 
308  // Variables used to keep track with the best muon PID.
309  float bestPID = -5.0; // best muonPID value
310  int prongPick = -1; // Prong picked as having the longest track.
311 
312 
313 
314  //////////
315  //
316  // Step 1: Determine which prong is the most muon-like.
317  //
318  //////////
319 
320  // loop over all prongs and get the tracks
321  for(unsigned int iprong = 0; iprong < prongs.size(); ++iprong) {
322 
323  // get the tracks associated with this prong
324  std::vector< art::Ptr< rb::Track > > tracks;
325  if(trackFMP.isValid()) tracks = trackFMP.at(iprong);
326 
327  // create the FindMany object for getting FitSum objects
328  art::FindManyP<rb::FitSum> fitsumFMP(tracks, evt, fTrackLabel);
329 
330  // Loop over all tracks...
331  for(unsigned int itrack = 0; itrack < tracks.size(); ++itrack) {
332 
333  // get the FitSum object associated with this track
334  std::vector< art::Ptr< rb::FitSum > > fitsums;
335  if(fitsumFMP.isValid()) fitsums = fitsumFMP.at(itrack);
336 
337  // There can be only one!!!
338  // (something went really wrong with BPF if there's not one and
339  // only one of these per track.)
340  if(fitsums.size() != 1) abort();
341 
342  // For now, this energy estimator only uses information from the muon track fit.
343  if(fitsums[0]->PDG() != 13) continue;
344 
345  float LLtemp = 0.0;
346 
347  // fill in chi2 and NDOF values
348  float tempH = 0.0;
349  float tempA = 0.0;
350  float tempN = 0.0;
351  tempH = fitsums[0]->Chi2(0)+fitsums[0]->Chi2(1); // chi2 from the hits
352  tempA = fitsums[0]->Chi2(2)+fitsums[0]->Chi2(3); // chi2 from the angles
353  tempN = fitsums[0]->Ndof(0)+fitsums[0]->Ndof(1)-4; // NDOF (NOTE: currently 4 is NOT subtracted from NDOF when the fitsum object is made by BPF)
354 
355  if(tempN <= 0) continue;
356 
357  // Keep track of the relevant info...
358  fChi2T = tempH/tempN + tempA/tempN;
359  fLength = tracks[itrack]->TotalLength();
360 
361  //
362  // Calculate dE/dx for each cell the track passes through using
363  // the helper class.
364  //
365  fdEdxCalc.computeDEDX(tracks[itrack],fitsums[0]->PDG());
366 
367  std::vector<double> dE;
368  std::vector<double> dx;
369  std::vector<double> W;
370  std::vector<double> BG;
371  std::vector<double> s; // this is a throw away variable...
372 
373  fdEdxCalc.getResults(dE,dx,W,BG,s,fdxTOL);
374 
375  //
376  // Compute the de/dx likelihood values for all cells.
377  //
378  for(unsigned int a = 0; a < dE.size(); ++a) {
379 
380  unsigned int w = 0;
381  unsigned int x = 0;
382 
383  if(W[a] > 0.0) w = 1;
384  if(dx[a] > 8.0) x = 1;
385 
386  // For the NearDet, we will just the "close to the readout end" FarDet histos.
387  if(geom->DetId() == novadaq::cnv::kNEARDET) w = 1;
388 
389  double content = fdEdxVSbg[w][x]->GetBinContent(fdEdxVSbg[w][x]->FindBin(log(BG[a]),dE[a]/dx[a]));
390  // We can't take the log of zero, so put a minimum on it.
391  if(content < 1.0e-9) content = 1.0e-9;
392  LLtemp += log(content);
393 
394  }
395 
396  if(dE.size() > 0) {
397  // Keep track of the relevant info...
398  fdEdxLL = LLtemp/(float)dE.size();
399  // NOTE: For this ratio, the number of hits on the track is taken to be the
400  // number of hits for which a value of dE/dx was successfully calculated.
401  fHitRatio = (float)dE.size()/(float)prongs[iprong]->NCell();
402  }
403  // end computing dE/dx LL info...
404 
405 
406 
407  // Set the variables for the track with the best muonPID
408  if(fLength > 0.0 && fdEdxLL > -1.0e9 && fChi2T < 1.0e9) {
409  // Get the PID value...
410  float pid = fMuonPID.EvaluateMVA("KNN Method");
411  if(pid > bestPID) {
412  bestPID = pid;
413  prongPick = (int)iprong;
414  fPmu = fitsums[0]->FourMom().P();
415  fDirZmu = tracks[itrack]->Dir().Z();
416  }
417  }
418 
419  } // end loop over tracks (itrack)
420 
421  } // end loop over prongs (iprong)
422 
423 
424 
425  //////////
426  //
427  // Step 2: Compute the values to be used for the energy estimator.
428  //
429  //////////
430 
431  //
432  // Add up the energy of the non-muon 3D prongs.
433  //
434 
435  float N3Dpr = 0.0;
436  float E3Dpr = 0.0;
437  float sumPE = 0.0;
438  rb::Cluster hits3Dprongs; // a container for all hits on 3D prongs
439 
440  for(unsigned int p = 0; p < prongs.size(); ++p) {
441 
442  // Skip 2D prongs.
443  if(prongs[p]->Is2D()) continue;
444  N3Dpr++;
445 
446  // Keep track of all hits in 3D prongs (to be used later when adding up
447  // the remaining slice energy.)
448  // NOTE: FuzzyK prongs can sometimes share hits, so we have to be careful
449  // not to double count.
450  for(unsigned int h = 0; h < prongs[p]->NCell(); ++h) {
451 
452  const rb::CellHit *pronghit = prongs[p]->Cell(h).get();
453 
454  // Check to see if this hit is already in the list.
455  bool onlist = false;
456  for(unsigned int hpr = 0; hpr < hits3Dprongs.NCell(); ++hpr) {
457  const rb::CellHit *pronghit3D = hits3Dprongs.Cell(hpr).get();
458  if((*pronghit) == (*pronghit3D)) {
459  onlist = true;
460  break;
461  }
462  } // end loop over hits in 3D prongs (hpr)
463 
464  // Add the hit if its not already on the list.
465  if(!onlist) hits3Dprongs.Add(prongs[p]->Cell(h));
466 
467  } // end loop over hits (h)
468 
469  // Skip this prong if it is the identified muon prong.
470  if((int)p == prongPick) continue;
471 
472  // Loop over all hits in this 3D, non-muon prong and add up their energies.
473  //
474  // NOTE: For the hit W position, I am currently using the prong W instead of
475  // the track W. Since I haven't determined which of the three BPF fits
476  // is the best one for this prong, this saves me from having to try
477  // to guess at which track to use. My assumption is that prong->W will
478  // be good enough since I am lumping all hits from 3D prongs together
479  // into one 3D prong energy number...
480  //
481  // WARNING: Since hits can be shared by fuzzyK prongs, there is the potential
482  // below to double (or triple) count energy on shared cell hits. For
483  // now, I am assuming that this is a small effect but something
484  // smarter could/should be put into place in the future.
485  for(unsigned int h = 0; h < prongs[p]->NCell(); ++h) {
486 
487  const rb::CellHit *chit = prongs[p]->Cell(h).get();
488  float W = prongs[p]->W(chit);
489  rb::RecoHit rhit(cal->MakeRecoHit(*chit,W));
490 
491  sumPE += chit->PE();
492  if(rhit.IsCalibrated()) E3Dpr += rhit.GeV();
493 
494  } // end loop over hits (h)
495 
496  } // end loop over prongs (p)
497 
498  // Set the prong variables.
499  if(bestPID >= 0.0) {
500  fN3Dprongs = N3Dpr;
501  fE3Dprongs = E3Dpr;
502  }
503 
504 
505 
506  //
507  // Add up the remaining energy in the slice.
508  //
509 
510  rb::Cluster hitsRemain; // a container for all hits NOT on 3D prongs
511 
512  // Loop over all hits in this slice.
513  for(unsigned int hsl = 0; hsl < (*slhandle)[islice].NCell(); ++hsl) {
514 
515  const rb::CellHit *slicehit = (*slhandle)[islice].Cell(hsl).get();
516 
517  // Check to see if this hit is in the list of hits on 3D prongs.
518  bool onlist = false;
519  for(unsigned int hpr = 0; hpr < hits3Dprongs.NCell(); ++hpr) {
520  const rb::CellHit *pronghit = hits3Dprongs.Cell(hpr).get();
521  if((*pronghit) == (*slicehit)) {
522  onlist = true;
523  break;
524  }
525  } // end loop over hits in 3D prongs (hpr)
526 
527  // If this hit is not in a 3D prong, add it to the list of remaining slice hits.
528  if(!onlist) {
529  hitsRemain.Add((*slhandle)[islice].Cell(hsl));
530  sumPE += slicehit->PE();
531  }
532 
533  } // end loop over hits in this slice (hsl)
534 
535  // Set the remaining slice energy variable.
536  if(bestPID >= 0.0) {
537  if(hitsRemain.NCell() > 0) fEremain = hitsRemain.TotalGeV();
538  else fEremain = 0.0;
539  }
540 
541  if(bestPID >= 0.0) {
542  fSumPE = sumPE;
543  }
544 
545 
546 
547  //
548  // Compute the hadronic contamination on the muon track and add that back in.
549  //
550  float muonTrackHadE = 0.0;
551 
552  // Only do this if we actually picked a track to be the muon.
553  if(prongPick >= 0) {
554 
555  std::vector< art::Ptr< rb::Track > > tracks;
556  if(trackFMP.isValid()) tracks = trackFMP.at(prongPick);
557  art::FindManyP<rb::FitSum> fitsumFMP(tracks, evt, fTrackLabel);
558 
559  // Figure out which BPF track was made with the muon assumption.
560  for(unsigned int itrack = 0; itrack < tracks.size(); ++itrack) {
561 
562  std::vector< art::Ptr< rb::FitSum > > fitsums;
563  if(fitsumFMP.isValid()) fitsums = fitsumFMP.at(itrack);
564  if(fitsums.size() != 1) abort();
565 
566  // Skip if not the muon track...
567  if(fitsums[0]->PDG() != 13) continue;
568 
569  // If it is the muon track, do the extra Ehad calculation
570  muonTrackHadE = fTrackCleanUpAlg.ExtraEOnTrackInGeV(*tracks[itrack],*slices[islice]);
571 
572  } // end loop over tracks (itrack)
573 
574  } // end if(prongPick >= 0)
575 
576  fEremain += muonTrackHadE;
577 
578 
579 
580  //
581  // Set the BPFEnergy object values and make the associations.
582  //
583  float Eres1 = -5.0;
584  float Eres2 = -5.0;
585  float Eres3 = -5.0;
586  if(bestPID >= 0.0) {
587  fE1 = (fEnergyEstm1.EvaluateRegression("KNN Method"))[0];
588  fEsq1 = (fEnergyEstm1.EvaluateRegression("KNN Method"))[1];
589  Eres1 = sqrt(fabs(fEsq1 - fE1*fE1));
590 
591  fE2 = (fEnergyEstm2.EvaluateRegression("KNN Method"))[0];
592  fEsq2 = (fEnergyEstm2.EvaluateRegression("KNN Method"))[1];
593  Eres2 = sqrt(fabs(fEsq2 - fE2*fE2));
594 
595  fE3 = (fEnergyEstm3.EvaluateRegression("KNN Method"))[0];
596  fEsq3 = (fEnergyEstm3.EvaluateRegression("KNN Method"))[1];
597  Eres3 = sqrt(fabs(fEsq3 - fE3*fE3));
598  }
599 
600  BPFEnergy energy(fE1,Eres1,fE2,Eres2,fE3,Eres3,
601  bestPID,fPmu,fDirZmu,
603  energies->push_back(energy);
604  util::CreateAssn(*this,evt,*(energies.get()),slices[islice],*(EnSlAssns.get()));
605 
606  } // end loop over slhandle (islice)
607 
608  evt.put(std::move(energies));
609  evt.put(std::move(EnSlAssns));
610 
611  }
612 
613 } // end namespace bpfit
614 ////////////////////////////////////////////////////////////////////////
615 
616 
617 
619 
620 namespace bpfit
621 {
623 }
std::string fTrackLabel
label for module that made the tracks
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.
float fEsq2
second estimator true energy squared (used to estiamte st. dev.)
unsigned int NCell(geo::View_t view) const
Number of cells in view view.
Definition: Cluster.cxx:134
std::string fHistoFile
name of dE/dx log-likelihood histogram file
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
void reconfigure(const fhicl::ParameterSet &pset)
bpfit::dEdxCalculator fdEdxCalc
helper class for computing dEdx values
double fdxTOL
lower limit for dx values to be used in dE/dx calculation
TMVA::Reader fMuonPID
TMVA reader object for the muon PID.
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
std::string fLibPath
location of all BPF library files listed below
float fEremain
sum of remaining slice energy not in 3D prongs - this does NOT account for dead material ...
void computeDEDX(art::Ptr< rb::Track > &track, int pdg=13)
Compute dE/dx for all cells along this track and the fitsum object that went with it...
std::string fEnXMLFile3
name of the XML file for the third energy estimator
string filename
Definition: shutoffs.py:106
std::string fEnXMLFile1
name of the XML file for the first energy estimator
TMVA::Reader fEnergyEstm3
TMVA reader object for the third energy estimator.
float fPmu
reconstructed muon momentum from the track selected to be the muon
DEFINE_ART_MODULE(TestTMapFile)
float fE2
second estimator true energy
double dE
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 fEsq3
third estimator true energy squared (used to estiamte st. dev.)
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
unsigned short Cell() const
Definition: CellHit.h:40
TMVA::Reader fEnergyEstm1
TMVA reader object for the first energy estimator.
const XML_Char * s
Definition: expat.h:262
double dx[NP][NC]
float fE1
first estimator true energy
string infile
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:441
float fSumPE
sum of total PE for all hits not on the muon track
std::string fEnXMLFile2
name of the XML file for the second energy estimator
const double a
TH2F * fdEdxVSbg[2][2]
histo of computed dE/dx for cell hits (vs. beta*gamma)
T get(std::string const &key) const
Definition: ParameterSet.h:231
novadaq::cnv::DetId DetId() const
Prefer ds::DetectorService::DetId() instead.
Definition: GeometryBase.h:243
double energy
Definition: plottest35.C:25
float PE() const
Definition: CellHit.h:42
Near Detector in the NuMI cavern.
float fEsq1
first estimator true energy squared (used to estiamte st. dev.)
murem::TrackCleanUpAlg fTrackCleanUpAlg
object to determine hadronic energy on the muon track
Perform a "2 point" Hough transform on a collection of hits.
Definition: run.py:1
float fN3Dprongs
total number of 3D prongs
art::Ptr< rb::CellHit > Cell(geo::View_t view, unsigned int viewIdx) const
Get the ith cell from view view.
Definition: Cluster.cxx:145
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
A rawdata::RawDigit with channel information decoded.
Definition: CellHit.h:27
float fE3Dprongs
sum of hit energies in 3D prongs (excluding the "muon" prong) - this does NOT account for dead materi...
T const * get() const
Definition: Ptr.h:321
const art::ProductToken< std::vector< rb::Cluster > > fSlicerToken
void geom(int which=0)
Definition: geom.C:163
std::string fPIDXMLFile
name of the XML file for the muon PID
double TotalGeV(EEnergyCalcScheme escheme=kRecomputeEnergy) const
Simple sum of the estimated GeV of all the hits.
Definition: Cluster.cxx:378
bool fBookedKNNs
Have we yet booked the kNNs?
bool getByToken(ProductToken< PROD > const &token, Handle< PROD > &result) const
Definition: DataViewImpl.h:387
float fHitRatio
ratio of track hits to prong hits (used as a track/shower discriminator)
Float_t e
Definition: plot.C:35
BPFEnergyEstimator(fhicl::ParameterSet const &pset)
Float_t w
Definition: plot.C:20
std::string fProngInstance
instance label for the prongs to be used
#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)
float fE3
third estimator true energy
ProductToken< T > consumes(InputTag const &)
float fDirZmu
muon track dirZ from the track selected to be the muon
void getResults(std::vector< double > &dEtemp, std::vector< double > &dxtemp, std::vector< double > &Wtemp, std::vector< double > &bgtemp, std::vector< double > &stemp, double dxTol)
Return the results of the dE/dx calculation.
TMVA::Reader fEnergyEstm2
TMVA reader object for the second energy estimator.
Encapsulate the geometry of one entire detector (near, far, ndos)
std::string fProngLabel
label for module that made the prongs