LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
PFParticleCosmicAna_module.cc
Go to the documentation of this file.
1 
9 
11 
12 #include "TTree.h"
13 
14 #include <string>
15 
16 //------------------------------------------------------------------------------------------------------------------------------------------
17 
18 namespace lar_pandora {
19 
24  public:
31 
35  virtual ~PFParticleCosmicAna();
36 
37  void beginJob();
38  void endJob();
39  void analyze(const art::Event& evt);
40  void reconfigure(fhicl::ParameterSet const& pset);
41 
42  private:
50  void FillRecoTree(const PFParticlesToHits& recoParticlesToHits,
51  const PFParticlesToTracks& recoParticlesToTracks,
52  const TracksToCosmicTags& recoTracksToCosmicTags);
53 
64  void FillTrueTree(const HitVector& hitVector,
65  const HitsToMCParticles& trueHitsToParticles,
66  const HitsToPFParticles& recoHitsToParticles,
67  const MCParticlesToMCTruth& particlesToTruth,
68  const PFParticlesToTracks& particlesToTracks,
69  const TracksToCosmicTags& tracksToCosmicTags);
70 
78  float GetCosmicScore(const art::Ptr<recob::PFParticle> particle,
79  const PFParticlesToTracks& recoParticlesToTracks,
80  const TracksToCosmicTags& recoTracksToCosmicTags) const;
81 
82  TTree* m_pRecoTree;
83  TTree* m_pTrueTree;
84 
85  int m_run;
86  int m_event;
87  int m_index;
88 
89  int m_self;
90  int m_pdgCode;
93  float m_cosmicScore;
94  int m_nTracks;
95  int m_nHits;
96 
97  float m_trackVtxX;
98  float m_trackVtxY;
99  float m_trackVtxZ;
100  float m_trackEndX;
101  float m_trackEndY;
102  float m_trackEndZ;
115 
118 
125 
132 
133  std::string m_hitfinderLabel;
134  std::string m_trackfitLabel;
135  std::string m_particleLabel;
136  std::string m_cosmicLabel;
137  std::string m_geantModuleLabel;
138 
141 
143  };
144 
146 
147 } // namespace lar_pandora
148 
149 //------------------------------------------------------------------------------------------------------------------------------------------
150 // implementation follows
151 
155 #include "art_root_io/TFileDirectory.h"
156 #include "art_root_io/TFileService.h"
159 #include "fhiclcpp/ParameterSet.h"
161 
167 
168 #include <iostream>
169 
170 namespace lar_pandora {
171 
173  {
174  this->reconfigure(pset);
175  }
176 
177  //------------------------------------------------------------------------------------------------------------------------------------------
178 
180 
181  //------------------------------------------------------------------------------------------------------------------------------------------
182 
184  {
185  m_cosmicLabel = pset.get<std::string>("CosmicTagModule", "cosmictagger");
186  m_particleLabel = pset.get<std::string>("PFParticleModule", "pandora");
187  m_trackfitLabel = pset.get<std::string>("TrackFitModule", "trackfit");
188  m_hitfinderLabel = pset.get<std::string>("HitFinderModule", "gaushit");
189  m_geantModuleLabel = pset.get<std::string>("GeantModule", "largeant");
190 
191  m_useDaughterPFParticles = pset.get<bool>("UseDaughterPFParticles", true);
192  m_useDaughterMCParticles = pset.get<bool>("UseDaughterMCParticles", true);
193 
194  m_cosmicContainmentCut = pset.get<double>("CosmicContainmentCut", 5.0);
195  }
196 
197  //------------------------------------------------------------------------------------------------------------------------------------------
198 
200  {
201  mf::LogDebug("LArPandora") << " *** PFParticleCosmicAna::beginJob() *** " << std::endl;
202 
203  //
205 
206  m_pRecoTree = tfs->make<TTree>("recoTree", "LAr Cosmic Reco Tree");
207  m_pRecoTree->Branch("run", &m_run, "run/I");
208  m_pRecoTree->Branch("event", &m_event, "event/I");
209  m_pRecoTree->Branch("index", &m_index, "index/I");
210  m_pRecoTree->Branch("self", &m_self, "self/I");
211  m_pRecoTree->Branch("pdgCode", &m_pdgCode, "pdgCode/I");
212  m_pRecoTree->Branch("isTrackLike", &m_isTrackLike, "isTrackLike/I");
213  m_pRecoTree->Branch("isPrimary", &m_isPrimary, "isPrimary/I");
214  m_pRecoTree->Branch("cosmicScore", &m_cosmicScore, "cosmicScore/F");
215  m_pRecoTree->Branch("trackVtxX", &m_trackVtxX, "trackVtxX/F");
216  m_pRecoTree->Branch("trackVtxY", &m_trackVtxY, "trackVtxY/F");
217  m_pRecoTree->Branch("trackVtxZ", &m_trackVtxZ, "trackVtxZ/F");
218  m_pRecoTree->Branch("trackEndX", &m_trackEndX, "trackEndX/F");
219  m_pRecoTree->Branch("trackEndY", &m_trackEndY, "trackEndY/F");
220  m_pRecoTree->Branch("trackEndZ", &m_trackEndZ, "trackEndZ/F");
221  m_pRecoTree->Branch("trackVtxDirX", &m_trackVtxDirX, "trackVtxDirX/F");
222  m_pRecoTree->Branch("trackVtxDirY", &m_trackVtxDirY, "trackVtxDirY/F");
223  m_pRecoTree->Branch("trackVtxDirZ", &m_trackVtxDirZ, "trackVtxDirZ/F");
224  m_pRecoTree->Branch("trackEndDirX", &m_trackEndDirX, "trackEndDirX/F");
225  m_pRecoTree->Branch("trackEndDirY", &m_trackEndDirY, "trackEndDirY/F");
226  m_pRecoTree->Branch("trackEndDirZ", &m_trackEndDirZ, "trackEndDirZ/F");
227  m_pRecoTree->Branch("trackLength", &m_trackLength, "trackLength/F");
228  m_pRecoTree->Branch("trackWidthX", &m_trackWidthX, "trackWidthX/F");
229  m_pRecoTree->Branch("trackWidthY", &m_trackWidthY, "trackWidthY/F");
230  m_pRecoTree->Branch("trackWidthZ", &m_trackWidthZ, "trackWidthZ/F");
231  m_pRecoTree->Branch("trackVtxDeltaYZ", &m_trackVtxDeltaYZ, "trackVtxDeltaYZ/F");
232  m_pRecoTree->Branch("trackEndDeltaYZ", &m_trackEndDeltaYZ, "trackEndDeltaYZ/F");
233  m_pRecoTree->Branch("trackVtxContained", &m_trackVtxContained, "trackVtxContained/I");
234  m_pRecoTree->Branch("trackEndContained", &m_trackEndContained, "trackEndContained/I");
235  m_pRecoTree->Branch("nTracks", &m_nTracks, "nTracks/I");
236  m_pRecoTree->Branch("nHits", &m_nHits, "nHits/I");
237 
238  m_pTrueTree = tfs->make<TTree>("trueTree", "LAr Cosmic True Tree");
239  m_pTrueTree->Branch("run", &m_run, "run/I");
240  m_pTrueTree->Branch("event", &m_event, "event/I");
241  m_pTrueTree->Branch("nHits", &m_nHits, "nHits/I");
242  m_pTrueTree->Branch("nNeutrinoHits", &m_nNeutrinoHits, "nNeutrinoHits/I");
243  m_pTrueTree->Branch(
244  "nNeutrinoHitsFullyTagged", &m_nNeutrinoHitsFullyTagged, "nNeutrinoHitsFullyTagged/I");
245  m_pTrueTree->Branch(
246  "nNeutrinoHitsSemiTagged", &m_nNeutrinoHitsSemiTagged, "nNeutrinoHitsSemiTagged/I");
247  m_pTrueTree->Branch(
248  "nNeutrinoHitsNotTagged", &m_nNeutrinoHitsNotTagged, "nNeutrinoHitsNotTagged/I");
249  m_pTrueTree->Branch("nNeutrinoHitsNotReconstructed",
251  "nNeutrinoHitsNotReconstructed/I");
252  m_pTrueTree->Branch(
253  "nNeutrinoHitsReconstructed", &m_nNeutrinoHitsReconstructed, "nNeutrinoHitsReconstructed/I");
254  m_pTrueTree->Branch("nCosmicHits", &m_nCosmicHits, "nCosmicHits/I");
255  m_pTrueTree->Branch(
256  "nCosmicHitsFullyTagged", &m_nCosmicHitsFullyTagged, "nCosmicHitsFullyTagged/I");
257  m_pTrueTree->Branch(
258  "nCosmicHitsSemiTagged", &m_nCosmicHitsSemiTagged, "nCosmicHitsSemiTagged/I");
259  m_pTrueTree->Branch("nCosmicHitsNotTagged", &m_nCosmicHitsNotTagged, "nCosmicHitsNotTagged/I");
260  m_pTrueTree->Branch("nCosmicHitsNotReconstructed",
262  "nCosmicHitsNotReconstructed/I");
263  m_pTrueTree->Branch(
264  "nCosmicHitsReconstructed", &m_nCosmicHitsReconstructed, "nCosmicHitsReconstructed/I");
265  }
266 
267  //------------------------------------------------------------------------------------------------------------------------------------------
268 
270 
271  //------------------------------------------------------------------------------------------------------------------------------------------
272 
274  {
275  std::cout << " *** PFParticleCosmicAna::analyze(...) *** " << std::endl;
276 
277  //
278  // Note: I've made this is MicroBooNE-only module
279  //
280 
281  m_run = evt.run();
282  m_event = evt.id().event();
283 
284  std::cout << " Run: " << m_run << std::endl;
285  std::cout << " Event: " << m_event << std::endl;
286 
287  // Collect True Particles
288  // ======================
289  HitVector hitVector;
290  MCTruthToMCParticles truthToParticles;
291  MCParticlesToMCTruth particlesToTruth;
292  MCParticlesToHits trueParticlesToHits;
293  HitsToMCParticles trueHitsToParticles;
294 
297  evt, m_geantModuleLabel, truthToParticles, particlesToTruth);
300  hitVector,
301  trueParticlesToHits,
302  trueHitsToParticles,
306 
307  // Collect Reco Particles
308  // ======================
309  PFParticleVector recoParticleVector;
310  PFParticlesToHits recoParticlesToHits;
311  HitsToPFParticles recoHitsToParticles;
312 
313  LArPandoraHelper::CollectPFParticles(evt, m_particleLabel, recoParticleVector);
317  recoParticlesToHits,
318  recoHitsToParticles,
321  LArPandoraHelper::kIgnoreDaughters));
322 
323  std::cout << " PFParticles: " << recoParticleVector.size() << std::endl;
324 
325  // Collect Reco Tracks
326  // ===================
327  TrackVector recoTrackVector;
328  PFParticlesToTracks recoParticlesToTracks;
329  LArPandoraHelper::CollectTracks(evt, m_trackfitLabel, recoTrackVector, recoParticlesToTracks);
330 
331  // Collect Cosmic Tags
332  // =====================
333  CosmicTagVector recoCosmicTagVector;
334  TracksToCosmicTags recoTracksToCosmicTags;
336  evt, m_cosmicLabel, recoCosmicTagVector, recoTracksToCosmicTags);
337 
338  // Analyse Reconstructed Particles
339  // ===============================
340  this->FillRecoTree(recoParticlesToHits, recoParticlesToTracks, recoTracksToCosmicTags);
341 
342  // Analyse True Hits
343  // =================
344  this->FillTrueTree(hitVector,
345  trueHitsToParticles,
346  recoHitsToParticles,
347  particlesToTruth,
348  recoParticlesToTracks,
349  recoTracksToCosmicTags);
350  }
351 
352  //------------------------------------------------------------------------------------------------------------------------------------------
353 
354  void PFParticleCosmicAna::FillRecoTree(const PFParticlesToHits& recoParticlesToHits,
355  const PFParticlesToTracks& recoParticlesToTracks,
356  const TracksToCosmicTags& recoTracksToCosmicTags)
357  {
358  // Set up Geometry Service
359  // =======================
361 
362  const double xmin(0.0);
363  const double xmax(2.0 * theGeometry->DetHalfWidth());
364  const double ymin(-theGeometry->DetHalfHeight());
365  const double ymax(+theGeometry->DetHalfHeight());
366  const double zmin(0.0);
367  const double zmax(theGeometry->DetLength());
368  const double xyzCut(m_cosmicContainmentCut);
369 
370  m_index = 0;
371 
372  m_self = 0;
373  m_pdgCode = 0;
374  m_isTrackLike = 0;
375  m_isPrimary = 0;
376  m_cosmicScore = 0.f;
377 
378  m_trackVtxX = 0.f;
379  m_trackVtxY = 0.f;
380  m_trackVtxZ = 0.f;
381  m_trackEndX = 0.f;
382  m_trackEndY = 0.f;
383  m_trackEndZ = 0.f;
384  m_trackVtxDirX = 0.f;
385  m_trackVtxDirY = 0.f;
386  m_trackVtxDirZ = 0.f;
387  m_trackEndDirX = 0.f;
388  m_trackEndDirY = 0.f;
389  m_trackEndDirZ = 0.f;
390  m_trackLength = 0.f;
391  m_trackWidthX = 0.f;
392  m_trackWidthY = 0.f;
393  m_trackWidthZ = 0.f;
394  m_trackVtxDeltaYZ = 0.f;
395  m_trackEndDeltaYZ = 0.f;
396 
399 
400  m_nTracks = 0;
401  m_nHits = 0;
402 
403  // Loop over Reco Particles
404  // ========================
405  for (PFParticlesToHits::const_iterator iter1 = recoParticlesToHits.begin(),
406  iterEnd1 = recoParticlesToHits.end();
407  iter1 != iterEnd1;
408  ++iter1) {
409  const art::Ptr<recob::PFParticle> recoParticle = iter1->first;
410 
411  const HitVector& hitVector = iter1->second;
412  if (hitVector.empty()) continue;
413 
414  PFParticlesToTracks::const_iterator iter2 = recoParticlesToTracks.find(recoParticle);
415  if (recoParticlesToTracks.end() == iter2) continue;
416 
417  const TrackVector& trackVector = iter2->second;
418  if (trackVector.empty()) continue;
419 
420  m_nHits = hitVector.size();
421  m_nTracks = trackVector.size();
422 
423  m_self = recoParticle->Self();
424  m_pdgCode = recoParticle->PdgCode();
425  m_isPrimary = recoParticle->IsPrimary();
427  m_cosmicScore =
428  this->GetCosmicScore(recoParticle, recoParticlesToTracks, recoTracksToCosmicTags);
429 
430  m_trackVtxX = 0.f;
431  m_trackVtxY = 0.f;
432  m_trackVtxZ = 0.f;
433  m_trackEndX = 0.f;
434  m_trackEndY = 0.f;
435  m_trackEndZ = 0.f;
436  m_trackVtxDirX = 0.f;
437  m_trackVtxDirY = 0.f;
438  m_trackVtxDirZ = 0.f;
439  m_trackEndDirX = 0.f;
440  m_trackEndDirY = 0.f;
441  m_trackEndDirZ = 0.f;
442  m_trackLength = 0.f;
443  m_trackWidthX = 0.f;
444  m_trackWidthY = 0.f;
445  m_trackWidthZ = 0.f;
446  m_trackVtxDeltaYZ = 0.f;
447  m_trackEndDeltaYZ = 0.f;
448 
451 
452  for (TrackVector::const_iterator iter3 = trackVector.begin(), iterEnd3 = trackVector.end();
453  iter3 != iterEnd3;
454  ++iter3) {
455  const art::Ptr<recob::Track> track = *iter3;
456  const float trackLength(track->Length());
457 
458  if (trackLength < m_trackLength) continue;
459 
460  m_trackLength = trackLength;
461 
462  const auto& trackVtxPosition = track->Vertex();
463  const auto& trackVtxDirection = track->VertexDirection();
464  const auto& trackEndPosition = track->End();
465  const auto& trackEndDirection = track->EndDirection();
466 
467  m_trackVtxX = trackVtxPosition.x();
468  m_trackVtxY = trackVtxPosition.y();
469  m_trackVtxZ = trackVtxPosition.z();
470  m_trackVtxDirX = trackVtxDirection.x();
471  m_trackVtxDirY = trackVtxDirection.y();
472  m_trackVtxDirZ = trackVtxDirection.z();
473  m_trackEndX = trackEndPosition.x();
474  m_trackEndY = trackEndPosition.y();
475  m_trackEndZ = trackEndPosition.z();
476  m_trackEndDirX = trackEndDirection.x();
477  m_trackEndDirY = trackEndDirection.y();
478  m_trackEndDirZ = trackEndDirection.z();
479 
480  m_trackWidthX = std::fabs(m_trackEndX - m_trackVtxX);
481  m_trackWidthY = std::fabs(m_trackEndY - m_trackVtxY);
482  m_trackWidthZ = std::fabs(m_trackEndZ - m_trackVtxZ);
483 
485  std::min((ymax - m_trackVtxY), std::min((m_trackVtxZ - zmin), (zmax - m_trackVtxZ)));
487  std::min((m_trackEndY - ymin), std::min((m_trackEndZ - zmin), (zmax - m_trackEndZ)));
488 
489  m_trackVtxContained = ((m_trackVtxX > xmin + xyzCut && m_trackVtxX < xmax - xyzCut) &&
490  (m_trackVtxY > ymin + xyzCut && m_trackVtxY < ymax - xyzCut) &&
491  (m_trackVtxZ > zmin + xyzCut && m_trackVtxZ < zmax - xyzCut));
492  m_trackEndContained = ((m_trackEndX > xmin + xyzCut && m_trackEndX < xmax - xyzCut) &&
493  (m_trackEndY > ymin + xyzCut && m_trackEndY < ymax - xyzCut) &&
494  (m_trackEndZ > zmin + xyzCut && m_trackEndZ < zmax - xyzCut));
495  }
496 
497  std::cout << " PFParticle: [" << m_index << "] nHits=" << m_nHits
498  << ", nTracks=" << m_nTracks << ", cosmicScore=" << m_cosmicScore << std::endl;
499 
500  m_pRecoTree->Fill();
501  ++m_index;
502  }
503  }
504 
505  //------------------------------------------------------------------------------------------------------------------------------------------
506 
508  const HitsToMCParticles& trueHitsToParticles,
509  const HitsToPFParticles& recoHitsToParticles,
510  const MCParticlesToMCTruth& particlesToTruth,
511  const PFParticlesToTracks& particlesToTracks,
512  const TracksToCosmicTags& tracksToCosmicTags)
513  {
514  m_nHits = 0;
515 
516  m_nNeutrinoHits = 0;
522 
523  m_nCosmicHits = 0;
529 
530  for (HitVector::const_iterator iter2 = hitVector.begin(), iterEnd2 = hitVector.end();
531  iter2 != iterEnd2;
532  ++iter2) {
533  const art::Ptr<recob::Hit> hit = *iter2;
534 
535  HitsToMCParticles::const_iterator iter3 = trueHitsToParticles.find(hit);
536  if (trueHitsToParticles.end() == iter3) continue;
537 
538  const art::Ptr<simb::MCParticle> trueParticle = iter3->second;
539 
540  MCParticlesToMCTruth::const_iterator iter4 = particlesToTruth.find(trueParticle);
541  if (particlesToTruth.end() == iter4)
542  throw cet::exception("LArPandora") << " PFParticleCosmicAna::analyze --- Found a true "
543  "particle without any ancestry information ";
544 
545  const art::Ptr<simb::MCTruth> truth = iter4->second;
546 
547  float cosmicScore(-0.2);
548 
549  HitsToPFParticles::const_iterator iter5 = recoHitsToParticles.find(hit);
550  if (recoHitsToParticles.end() != iter5) {
551  const art::Ptr<recob::PFParticle> particle = iter5->second;
552  cosmicScore = this->GetCosmicScore(particle, particlesToTracks, tracksToCosmicTags);
553  }
554 
555  ++m_nHits;
556 
557  if (truth->NeutrinoSet()) {
558  ++m_nNeutrinoHits;
559 
560  if (cosmicScore >= 0)
562  else
564 
565  if (cosmicScore > 0.51)
567  else if (cosmicScore > 0.39)
569  else
571  }
572  else {
573  ++m_nCosmicHits;
574 
575  if (cosmicScore >= 0)
577  else
579 
580  if (cosmicScore > 0.51)
582  else if (cosmicScore > 0.39)
584  else
586  }
587  }
588 
589  m_pTrueTree->Fill();
590  }
591 
592  //------------------------------------------------------------------------------------------------------------------------------------------
593 
595  const PFParticlesToTracks& recoParticlesToTracks,
596  const TracksToCosmicTags& recoTracksToCosmicTags) const
597  {
598  float cosmicScore(0.f);
599 
600  // Get cosmic tags associated with this particle
601  PFParticlesToTracks::const_iterator iter2 = recoParticlesToTracks.find(particle);
602  if (recoParticlesToTracks.end() != iter2) {
603  for (TrackVector::const_iterator iter3 = iter2->second.begin(),
604  iterEnd3 = iter2->second.end();
605  iter3 != iterEnd3;
606  ++iter3) {
607  const art::Ptr<recob::Track> track = *iter3;
608 
609  TracksToCosmicTags::const_iterator iter4 = recoTracksToCosmicTags.find(track);
610  if (recoTracksToCosmicTags.end() != iter4) {
611  for (CosmicTagVector::const_iterator iter5 = iter4->second.begin(),
612  iterEnd5 = iter4->second.end();
613  iter5 != iterEnd5;
614  ++iter5) {
615  const art::Ptr<anab::CosmicTag> cosmicTag = *iter5;
616  if (cosmicTag->CosmicScore() > cosmicScore) cosmicScore = cosmicTag->CosmicScore();
617  }
618  }
619  }
620  }
621 
622  return cosmicScore;
623  }
624 
625 } //namespace lar_pandora
static void BuildPFParticleHitMaps(const PFParticleVector &particleVector, const PFParticlesToSpacePoints &particlesToSpacePoints, const SpacePointsToHits &spacePointsToHits, PFParticlesToHits &particlesToHits, HitsToPFParticles &hitsToParticles, const DaughterMode daughterMode=kUseDaughters)
Build mapping between PFParticles and Hits using PFParticle/SpacePoint/Hit maps.
std::map< art::Ptr< recob::Hit >, art::Ptr< simb::MCParticle > > HitsToMCParticles
size_t Self() const
Returns the index of this particle.
Definition: PFParticle.h:88
void FillTrueTree(const HitVector &hitVector, const HitsToMCParticles &trueHitsToParticles, const HitsToPFParticles &recoHitsToParticles, const MCParticlesToMCTruth &particlesToTruth, const PFParticlesToTracks &particlesToTracks, const TracksToCosmicTags &tracksToCosmicTags)
Fill track-level variables using input maps between reconstructed objects.
std::map< art::Ptr< simb::MCParticle >, art::Ptr< simb::MCTruth > > MCParticlesToMCTruth
Length_t DetHalfWidth(TPCID const &tpcid=tpc_zero) const
Returns the half width of the active volume of the specified TPC.
Vector_t VertexDirection() const
Access to track direction at different points.
Definition: Track.h:166
int PdgCode() const
Return the type of particle as a PDG ID.
Definition: PFParticle.h:79
intermediate_table::const_iterator const_iterator
std::vector< art::Ptr< anab::CosmicTag > > CosmicTagVector
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.cc:6
Length_t DetLength(TPCID const &tpcid=tpc_zero) const
Returns the length of the active volume of the specified TPC.
std::map< art::Ptr< recob::PFParticle >, TrackVector > PFParticlesToTracks
float & CosmicScore()
Definition: CosmicTag.h:54
std::map< art::Ptr< recob::Track >, CosmicTagVector > TracksToCosmicTags
TFile f
Definition: plotHisto.C:6
void reconfigure(fhicl::ParameterSet const &pset)
double Length(size_t p=0) const
Access to various track properties.
Definition: Track.h:207
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:65
std::map< art::Ptr< simb::MCParticle >, HitVector > MCParticlesToHits
std::vector< art::Ptr< recob::PFParticle > > PFParticleVector
std::map< art::Ptr< simb::MCTruth >, MCParticleVector > MCTruthToMCParticles
T get(std::string const &key) const
Definition: ParameterSet.h:314
std::map< art::Ptr< recob::PFParticle >, HitVector > PFParticlesToHits
Point_t const & Vertex() const
Access to track position at different points.
Definition: Track.h:158
Provides recob::Track data product.
bool IsPrimary() const
Returns whether the particle is the root of the flow.
Definition: PFParticle.h:82
std::vector< art::Ptr< recob::Track > > TrackVector
static void CollectTracks(const art::Event &evt, const std::string &label, TrackVector &trackVector, PFParticlesToTracks &particlesToTracks)
Collect the reconstructed PFParticles and associated Tracks from the ART event record.
static void CollectPFParticles(const art::Event &evt, const std::string &label, PFParticleVector &particleVector)
Collect the reconstructed PFParticles from the ART event record.
void FillRecoTree(const PFParticlesToHits &recoParticlesToHits, const PFParticlesToTracks &recoParticlesToTracks, const TracksToCosmicTags &recoTracksToCosmicTags)
Fill event-level variables using input maps between reconstructed objects.
static void CollectHits(const art::Event &evt, const std::string &label, HitVector &hitVector)
Collect the reconstructed Hits from the ART event record.
Detector simulation of raw signals on wires.
std::vector< art::Ptr< recob::Hit > > HitVector
float GetCosmicScore(const art::Ptr< recob::PFParticle > particle, const PFParticlesToTracks &recoParticlesToTracks, const TracksToCosmicTags &recoTracksToCosmicTags) const
Get cosmic score for a PFParticle using track-level information.
std::map< art::Ptr< recob::Hit >, art::Ptr< recob::PFParticle > > HitsToPFParticles
MaybeLogger_< ELseverityLevel::ELsev_success, false > LogDebug
Vector_t EndDirection() const
Access to track direction at different points.
Definition: Track.h:167
static bool IsTrack(const art::Ptr< recob::PFParticle > particle)
Determine whether a particle has been reconstructed as track-like.
Definition: MVAAlg.h:12
EventNumber_t event() const
Definition: EventID.h:116
Point_t const & End() const
Access to track position at different points.
Definition: Track.h:159
Length_t DetHalfHeight(TPCID const &tpcid=tpc_zero) const
Returns the half height of the active volume of the specified TPC.
bool NeutrinoSet() const
Definition: MCTruth.h:78
TCEvent evt
Definition: DataStructs.cxx:8
RunNumber_t run() const
Definition: Event.cc:29
static void BuildMCParticleHitMaps(const art::Event &evt, const HitVector &hitVector, const SimChannelVector &simChannelVector, HitsToTrackIDEs &hitsToTrackIDEs)
Collect the links from reconstructed hits to their true energy deposits.
static void CollectMCParticles(const art::Event &evt, const std::string &label, MCParticleVector &particleVector)
Collect a vector of MCParticle objects from the ART event record.
Float_t track
Definition: plot.C:35
helper function for LArPandoraInterface producer module
PFParticleCosmicAna(fhicl::ParameterSet const &pset)
Constructor.
EventID id() const
Definition: Event.cc:23
art framework interface to geometry description
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
static void CollectCosmicTags(const art::Event &evt, const std::string &label, CosmicTagVector &cosmicTagVector, TracksToCosmicTags &tracksToCosmicTags)
Collect a vector of cosmic tags from the ART event record.