LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
ParticleListAction.cxx
Go to the documentation of this file.
1 
9 #include "lardataobj/Simulation/sim.h" // sim::NoParticleId
12 
13 #include "cetlib_except/exception.h"
15 
16 #include "Geant4/G4DynamicParticle.hh"
17 #include "Geant4/G4ParticleDefinition.hh"
18 #include "Geant4/G4PrimaryParticle.hh"
19 #include "Geant4/G4Step.hh"
20 #include "Geant4/G4StepPoint.hh"
21 #include "Geant4/G4String.hh"
22 #include "Geant4/G4ThreeVector.hh"
23 #include "Geant4/G4Track.hh"
24 #include "Geant4/G4VProcess.hh"
25 #include "Geant4/G4VUserPrimaryParticleInformation.hh"
26 
27 #include <TLorentzVector.h>
28 
29 #include <algorithm>
30 #include <cassert>
31 
32 //const G4bool debug = false; // unused
33 
34 // Photon variables defined at each step, for use
35 // in temporary velocity bug fix. -wforeman
37 bool entra = true;
38 
39 namespace larg4 {
40 
41  // Initialize static members.
46 
47  //----------------------------------------------------------------------------
48  // Dropped particle test
49 
51  {
52  return !p || p->Trajectory().empty();
53  } // ParticleListAction::isDropped()
54 
55  //----------------------------------------------------------------------------
56  // Constructor.
58  bool storeTrajectories,
59  bool keepEMShowerDaughters,
60  bool keepMCParticleList, /* = true */
61  bool storeDroppedMCParticles)
62  : fenergyCut(energyCut * CLHEP::GeV)
63  , fparticleList(keepMCParticleList ? std::make_unique<sim::ParticleList>() : nullptr)
64  , fdroppedParticleList(storeDroppedMCParticles ? std::make_unique<sim::ParticleList>() :
65  nullptr)
66  , fstoreTrajectories(storeTrajectories)
67  , fKeepEMShowerDaughters(keepEMShowerDaughters)
68  {}
69 
70  //----------------------------------------------------------------------------
71  // Begin the event
73  {
74  // Clear any previous particle information.
76  if (fparticleList) fparticleList->clear();
78  fParentIDMap.clear();
82  fCurrentPdgCode = 0;
83  }
84 
85  //-------------------------------------------------------------
86  // figure out the ultimate parentage of the particle with track ID
87  // trackid
88  // assume that the current track id has already been added to
89  // the fParentIDMap
90  int ParticleListAction::GetParentage(int trackid, bool useOrigTrackIDMap) const
91  {
92  int parentid = sim::NoParticleId;
93  const std::map<int, int>* parentIDMap =
94  useOrigTrackIDMap ? &fParentIDMap_OrigTrackID : &fParentIDMap;
95 
96  // search the fParentIDMap recursively until we have the parent id
97  // of the first EM particle that led to this one
98  std::map<int, int>::const_iterator itr = parentIDMap->find(trackid);
99  while (itr != parentIDMap->end()) {
100  MF_LOG_DEBUG("ParticleListAction") << "parentage for " << trackid << " " << (*itr).second;
101 
102  // set the parentid to the current parent ID, when the loop ends
103  // this id will be the first EM particle
104  parentid = (*itr).second;
105  itr = parentIDMap->find(parentid);
106  }
107  MF_LOG_DEBUG("ParticleListAction") << "final parent ID " << parentid;
108 
109  return parentid;
110  }
111 
112  //----------------------------------------------------------------------------
113  // Create our initial simb::MCParticle object and add it to the sim::ParticleList.
115  {
116  // Particle type.
117  G4ParticleDefinition* particleDefinition = track->GetDefinition();
118  G4int pdgCode = particleDefinition->GetPDGEncoding();
119 
120  // Get Geant4's ID number for this track. This will be the same
121  // ID number that we'll use in the ParticleList.
122  // It is offset by the number of tracks accumulated from the previous Geant4
123  // runs (if any)
124  G4int trackID = track->GetTrackID() + fTrackIDOffset;
125  fCurrentTrackID = trackID;
126  fCurrentOrigTrackID = trackID;
127  fCurrentPdgCode = pdgCode;
128 
129  if (!fparticleList) {
130  // the rest is about adding a new particle to the list: skip
131  return; // note that fCurrentParticle is clear()'ed
132  }
133 
134  // And the particle's parent (same offset as above):
135  G4int parentID = track->GetParentID() + fTrackIDOffset;
136 
137  std::string process_name = "unknown";
138  bool drop_shower_daughter = false;
139 
140  // Is there an MCTruth object associated with this G4Track? We
141  // have to go up a "chain" of information to find out:
142  const G4DynamicParticle* dynamicParticle = track->GetDynamicParticle();
143  const G4PrimaryParticle* primaryParticle = dynamicParticle->GetPrimaryParticle();
145  if (primaryParticle != 0) {
146  const G4VUserPrimaryParticleInformation* gppi = primaryParticle->GetUserInformation();
148  dynamic_cast<const g4b::PrimaryParticleInformation*>(gppi);
149  if (ppi != 0) {
150  primaryIndex = ppi->MCParticleIndex();
151 
152  // If we've made it this far, a PrimaryParticleInformation
153  // object exists and we are using a primary particle, set the
154  // process name accordingly
155  process_name = "primary";
156 
157  // primary particles should have parentID = 0, even if there
158  // are multiple MCTruths for this event
159  parentID = 0;
160  } // end else no primary particle information
161  } // Is there a G4PrimaryParticle?
162  // If this is not a primary particle...
163  else {
164  // check if this particle was made in an EM shower, don't put it in the particle
165  // list as we don't care about secondaries, tertiaries, etc for these showers
166  // figure out what process is making this track - skip it if it is
167  // one of pair production, compton scattering, photoelectric effect
168  // bremstrahlung, annihilation, any ionization - who wants to save
169  // a buttload of electrons that arent from a CC interaction?
170  process_name = track->GetCreatorProcess()->GetProcessName();
171  drop_shower_daughter =
172  (!fKeepEMShowerDaughters && (process_name.find("conv") != std::string::npos ||
173  process_name.find("LowEnConversion") != std::string::npos ||
174  process_name.find("Pair") != std::string::npos ||
175  process_name.find("compt") != std::string::npos ||
176  process_name.find("Compt") != std::string::npos ||
177  process_name.find("Brem") != std::string::npos ||
178  process_name.find("phot") != std::string::npos ||
179  process_name.find("Photo") != std::string::npos ||
180  process_name.find("Ion") != std::string::npos ||
181  process_name.find("annihil") != std::string::npos));
182  if (drop_shower_daughter) {
183 
184  // figure out the ultimate parentage of this particle
185  // first add this track id and its parent to the fParentIDMap
186  fParentIDMap[trackID] = parentID;
187 
188  fCurrentTrackID = -1 * this->GetParentage(
189  trackID); // the real trackID remains stored in fCurrentOrigTrackID
190 
191  // check that fCurrentTrackID is in the particle list - it is possible
192  // that this particle's parent is a particle that did not get tracked.
193  // An example is a partent that was made due to muMinusCaptureAtRest
194  // and the daughter was made by the phot process. The parent likely
195  // isn't saved in the particle list because it is below the energy cut
196  // which will put a bogus track id value into the sim::IDE object for
197  // the sim::SimChannel if we don't check it.
199 
200  } // end if keeping EM shower daughters
201 
202  // Check the energy of the particle. If it falls below the energy
203  // cut, don't add it to our list.
204  G4double energy = track->GetKineticEnergy();
205  if (energy < fenergyCut) {
207 
208  // do add the particle to the parent id map though
209  // and set the current track id to be it's ultimate parent
210  fParentIDMap[trackID] = parentID;
211  fParentIDMap_OrigTrackID[trackID] = parentID;
212 
213  fCurrentTrackID = -1 * this->GetParentage(trackID);
214  fCurrentOrigTrackID = -1 * this->GetParentage(trackID, true);
215 
216  return;
217  }
218 
219  // check to see if the parent particle has been stored in the particle navigator
220  // if not, then see if it is possible to walk up the fParentIDMap to find the
221  // ultimate parent of this particle. Use that ID as the parent ID for this
222  // particle
223  if (!fparticleList->KnownParticle(parentID) &&
224  !(fdroppedParticleList && fdroppedParticleList->KnownParticle(parentID))) {
225  // do add the particle to the parent id map
226  // just in case it makes a daughter that we have to track as well
227  fParentIDMap[trackID] = parentID;
228  fParentIDMap_OrigTrackID[trackID] = parentID;
229  int pid = this->GetParentage(parentID);
230 
231  // if we still can't find the parent in the particle navigator,
232  // we have to give up
233  if (!fparticleList->KnownParticle(pid) &&
234  !(fdroppedParticleList && fdroppedParticleList->KnownParticle(parentID))) {
235  MF_LOG_WARNING("ParticleListAction")
236  << "can't find parent id: " << parentID << " in the particle list, or fParentIDMap."
237  << " Make " << parentID << " the mother ID for"
238  << " track ID " << fCurrentTrackID << " in the hope that it will aid debugging.";
239  }
240  else
241  parentID = pid;
242  }
243 
244  } // end if not a primary particle
245 
246  // This is probably the PDG mass, but just in case:
247  double mass = dynamicParticle->GetMass() / CLHEP::GeV;
248 
249  // Create the sim::Particle object.
252  new simb::MCParticle(trackID, pdgCode, process_name, parentID, mass);
253  fCurrentParticle.truthIndex = primaryIndex;
254 
255  // if we are not filtering, we have a decision already
256  if (!fFilter) fCurrentParticle.keep = true;
257 
258  // Polarization.
259  const G4ThreeVector& polarization = track->GetPolarization();
260  fCurrentParticle.particle->SetPolarization(
261  TVector3(polarization.x(), polarization.y(), polarization.z()));
262 
263  // if KeepEMShowerDaughters = False and we decided to drop this particle,
264  // record it before throwing it away.
265  if (drop_shower_daughter) {
266  fCurrentParticle.drop = true;
268  return;
269  }
270 
271  // Save the particle in the ParticleList.
273  }
274 
275  //----------------------------------------------------------------------------
276  void ParticleListAction::PostTrackingAction(const G4Track* aTrack)
277  {
278  if (!fCurrentParticle.hasParticle()) return;
279  assert(fparticleList);
280 
281  if (aTrack) {
282  fCurrentParticle.particle->SetWeight(aTrack->GetWeight());
283  G4String process =
284  aTrack->GetStep()->GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName();
285  fCurrentParticle.particle->SetEndProcess(process);
286  }
287 
288  // We would keep it if it wasn't for KeepEMShowerDaughters = false
289  // It goes into fDroppedParticleList in case LArG4 is configured
290  // to keep a minimal version of it.
291  if (fCurrentParticle.drop) {
295  return;
296  }
297 
298  // if we have found no reason to keep it, drop it!
299  // (we might still need parentage information though)
300  if (!fCurrentParticle.keep) {
301  fparticleList->Archive(fCurrentParticle.particle.get());
302  // after the particle is archived, it is deleted
304  return;
305  }
306 
307  // store truth record pointer, only if it is available
308  if (fCurrentParticle.isPrimary()) {
310  }
311 
312  return;
313  }
314 
315  //----------------------------------------------------------------------------
316  // With every step, add to the particle's trajectory.
317  void ParticleListAction::SteppingAction(const G4Step* step)
318  {
319 
320  if (!fCurrentParticle.hasParticle()) { return; }
321 
322  // Temporary fix for problem where DeltaTime on the first step
323  // of optical photon propagation is calculated incorrectly. -wforeman
324  globalTime = step->GetTrack()->GetGlobalTime();
325  velocity_G4 = step->GetTrack()->GetVelocity();
326  velocity_step = step->GetStepLength() / step->GetDeltaTime();
327  if ((step->GetTrack()->GetDefinition()->GetPDGEncoding() == 0) &&
328  fabs(velocity_G4 - velocity_step) > 0.0001) {
329  // Subtract the faulty step time from the global time,
330  // and add the correct step time based on G4 velocity.
331  step->GetPostStepPoint()->SetGlobalTime(globalTime - step->GetDeltaTime() +
332  step->GetStepLength() / velocity_G4);
333  }
334 
335  // For the most part, we just want to add the post-step
336  // information to the particle's trajectory. There's one
337  // exception: In PreTrackingAction, the correct time information
338  // is not available. So add the correct vertex information here.
339 
340  if (fCurrentParticle.particle->NumberTrajectoryPoints() == 0) {
341 
342  // Get the pre/along-step information from the G4Step.
343  const G4StepPoint* preStepPoint = step->GetPreStepPoint();
344 
345  const G4ThreeVector position = preStepPoint->GetPosition();
346  G4double time = preStepPoint->GetGlobalTime();
347 
348  // Remember that LArSoft uses cm, ns, GeV.
349  TLorentzVector fourPos(position.x() / CLHEP::cm,
350  position.y() / CLHEP::cm,
351  position.z() / CLHEP::cm,
352  time / CLHEP::ns);
353 
354  const G4ThreeVector momentum = preStepPoint->GetMomentum();
355  const G4double energy = preStepPoint->GetTotalEnergy();
356  TLorentzVector fourMom(momentum.x() / CLHEP::GeV,
357  momentum.y() / CLHEP::GeV,
358  momentum.z() / CLHEP::GeV,
359  energy / CLHEP::GeV);
360 
361  // Add the first point in the trajectory.
362  AddPointToCurrentParticle(fourPos, fourMom, "Start");
363 
364  } // end if this is the first step
365 
366  // At this point, the particle is being transported through the
367  // simulation. This method is being called for every voxel that
368  // the track passes through, but we don't want to update the
369  // trajectory information if we're just updating voxels. To check
370  // for this, look at the process name for the step, and compare it
371  // against the voxelization process name (set in PhysicsList.cxx).
372  G4String process = step->GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName();
373  G4bool ignoreProcess = process.contains("LArVoxel") || process.contains("OpDetReadout");
374 
375  MF_LOG_DEBUG("ParticleListAction::SteppingAction")
376  << ": DEBUG - process='" << process << "'"
377  << " ignoreProcess=" << ignoreProcess << " fstoreTrajectories=" << fstoreTrajectories;
378 
379  // We store the initial creation point of the particle
380  // and its final position (ie where it has no more energy, or at least < 1 eV) no matter
381  // what, but whether we store the rest of the trajectory depends
382  // on the process, and on a user switch
383  if (fstoreTrajectories && !ignoreProcess) {
384  // If the particle was marked to drop (EM shower daughters),
385  // then skip unless this step is last in volume.
386  //if (fCurrentParticle.drop && !step->IsLastStepInVolume()) return;
387 
388  // Get the post-step information from the G4Step.
389  const G4StepPoint* postStepPoint = step->GetPostStepPoint();
390 
391  const G4ThreeVector position = postStepPoint->GetPosition();
392  G4double time = postStepPoint->GetGlobalTime();
393 
394  // Remember that LArSoft uses cm, ns, GeV.
395  TLorentzVector fourPos(position.x() / CLHEP::cm,
396  position.y() / CLHEP::cm,
397  position.z() / CLHEP::cm,
398  time / CLHEP::ns);
399 
400  const G4ThreeVector momentum = postStepPoint->GetMomentum();
401  const G4double energy = postStepPoint->GetTotalEnergy();
402  TLorentzVector fourMom(momentum.x() / CLHEP::GeV,
403  momentum.y() / CLHEP::GeV,
404  momentum.z() / CLHEP::GeV,
405  energy / CLHEP::GeV);
406 
407  // Add another point in the trajectory.
408  AddPointToCurrentParticle(fourPos, fourMom, std::string(process));
409  }
410  }
411 
412  //----------------------------------------------------------------------------
416  public:
417  UpdateDaughterInformation() : particleList(0) {}
418  void SetParticleList(sim::ParticleList* p) { particleList = p; }
419  void operator()(sim::ParticleList::value_type& particleListEntry) const
420  {
421  // We're looking at this Particle in the list.
422  int particleID = particleListEntry.first;
423 
424  // The parent ID of this particle;
425  // we ask the particle list since the particle itself might have been lost
426  // ("archived"), but the particle list still holds the information we need
427  int parentID = particleList->GetMotherOf(particleID);
428 
429  // If the parentID <= 0, this is a primary particle.
430  if (parentID <= 0) return;
431 
432  // If we get here, this particle is somebody's daughter. Add
433  // it to the list of daughter particles for that parent.
434 
435  // Get the parent particle from the list.
436  sim::ParticleList::iterator parentEntry = particleList->find(parentID);
437 
438  if (parentEntry == particleList->end()) {
439  // We have an "orphan": a particle whose parent isn't
440  // recorded in the particle list. This is not signficant;
441  // it's possible for a particle not to be saved in the list
442  // because it failed an energy cut, but for it to have a
443  // daughter that passed the cut (e.g., a nuclear decay).
444  return;
445  }
446  if (!parentEntry->second) return; // particle archived, nothing to update
447 
448  // Add the current particle to the daughter list of the
449  // parent.
450  simb::MCParticle* parent = (*parentEntry).second;
451  parent->AddDaughter(particleID);
452  }
453 
454  private:
455  sim::ParticleList* particleList;
456  };
457 
458  //----------------------------------------------------------------------------
459  // There's one last thing to do: All the particles have their
460  // parent IDs set (in PostTrackingAction), but we haven't set the
461  // daughters yet. That's done in this method.
463  {
464  if (!fparticleList) return;
465 
466  // Set up the utility class for the "for_each" algorithm. (We only
467  // need a separate set-up for the utility class because we need to
468  // give it the pointer to the particle list. We're using the STL
469  // "for_each" instead of the C++ "for loop" because it's supposed
470  // to be faster.
471  UpdateDaughterInformation updateDaughterInformation;
472  updateDaughterInformation.SetParticleList(fparticleList.get());
473 
474  // Update the daughter information for each particle in the list.
475  std::for_each(fparticleList->begin(), fparticleList->end(), updateDaughterInformation);
476  }
477 
478  //----------------------------------------------------------------------------
479  // Returns the ParticleList accumulated during the current event.
481  {
482  if (!fparticleList) return nullptr;
483 
484  // check if the ParticleNavigator has entries, and if
485  // so grab the highest track id value from it to
486  // add to the fTrackIDOffset
487  int highestID = 0;
488  for (auto pn = fparticleList->begin(); pn != fparticleList->end(); pn++)
489  if ((*pn).first > highestID) highestID = (*pn).first;
490 
491  // If we have stored dropped particles,
492  // include them in the offset.
493  if (fdroppedParticleList) {
494  for (auto pn = fdroppedParticleList->begin(); pn != fdroppedParticleList->end(); pn++)
495  if ((*pn).first > highestID) highestID = (*pn).first;
496  }
497 
498  //Only change the fTrackIDOffset if there is in fact a particle to add to the event
499  if ((fparticleList->size()) != 0) { fTrackIDOffset = highestID + 1; }
500 
501  return fparticleList.get();
502  }
503 
504  //----------------------------------------------------------------------------
506  {
507  auto const iInfo = GetPrimaryTruthMap().find(trackId);
508  return (iInfo == GetPrimaryTruthMap().end()) ? simb::NoGeneratedParticleIndex : iInfo->second;
509  } // ParticleListAction::GetPrimaryTruthIndex()
510 
511  //----------------------------------------------------------------------------
512  // Yields the ParticleList accumulated during the current event.
514  {
515  if (!fparticleList) {
516  // check with `hasList()` before calling this method
517  throw cet::exception("ParticleListAction")
518  << "ParticleListAction::YieldList(): particle list not build by user request.\n";
519  }
520  // check if the ParticleNavigator has entries, and if
521  // so grab the highest track id value from it to
522  // add to the fTrackIDOffset
523  int highestID = 0;
524  for (auto pn = fparticleList->begin(); pn != fparticleList->end(); pn++)
525  if ((*pn).first > highestID) highestID = (*pn).first;
526 
527  // If we have stored dropped particles,
528  // include them in the offset.
529  if (fdroppedParticleList) {
530  for (auto pn = fdroppedParticleList->begin(); pn != fdroppedParticleList->end(); pn++)
531  if ((*pn).first > highestID) highestID = (*pn).first;
532  }
533 
534  //Only change the fTrackIDOffset if there is in fact a particle to add to the event
535  if ((fparticleList->size()) != 0) { fTrackIDOffset = highestID + 1; }
536 
537  return std::move(*fparticleList);
538  } // ParticleList&& ParticleListAction::YieldList()
539 
540  //----------------------------------------------------------------------------
541  // Yields the (dropped) ParticleList accumulated during the current event.
543  {
544  if (!fdroppedParticleList) {
545  throw cet::exception("ParticleListAction")
546  << "ParticleListAction::YieldDroppedList(): dropped particle list not build by user "
547  "request.\n";
548  }
549  // check if the ParticleNavigator has entries, and if
550  // so grab the highest track id value from it to
551  // add to the fTrackIDOffset
552  int highestID = 0;
553  for (auto pn = fparticleList->begin(); pn != fparticleList->end(); pn++)
554  if ((*pn).first > highestID) highestID = (*pn).first;
555 
556  // If we have stored dropped particles,
557  // include them in the offset.
558  for (auto pn = fdroppedParticleList->begin(); pn != fdroppedParticleList->end(); pn++)
559  if ((*pn).first > highestID) highestID = (*pn).first;
560 
561  //Only change the fTrackIDOffset if there is in fact a particle to add to the event
562  if ((fparticleList->size()) != 0) { fTrackIDOffset = highestID + 1; }
563 
564  return std::move(*fdroppedParticleList);
565  } // ParticleList&& ParticleListAction::YieldDroppedList()
566  //----------------------------------------------------------------------------
567  void ParticleListAction::AddPointToCurrentParticle(TLorentzVector const& pos,
568  TLorentzVector const& mom,
569  std::string const& process)
570  {
571 
572  // Add the first point in the trajectory.
573  fCurrentParticle.particle->AddTrajectoryPoint(pos, mom, process);
574 
575  // also see if we can decide to keep the particle
576  if (!fCurrentParticle.keep) fCurrentParticle.keep = fFilter->mustKeep(pos);
577 
578  } // ParticleListAction::AddPointToCurrentParticle()
579 
580  //----------------------------------------------------------------------------
581 
582 } // namespace LArG4
ParticleListAction(double energyCut, bool storeTrajectories=false, bool keepEMShowerDaughters=false, bool keepMCParticleList=true, bool storeDroppedMCParticles=false)
G4bool fstoreTrajectories
Whether to store particle trajectories with each particle.
void AddDaughter(const int trackID)
Definition: MCParticle.h:269
void AddPointToCurrentParticle(TLorentzVector const &pos, TLorentzVector const &mom, std::string const &process)
Adds a trajectory point to the current particle, and runs the filter.
static int fCurrentPdgCode
pdg code of current particle
virtual void BeginOfEventAction(const G4Event *)
G4UserEventAction interfaces.
std::map< int, GeneratedParticleIndex_t > fPrimaryTruthMap
Map: particle track ID -> index of primary information in MC truth.
std::unique_ptr< sim::ParticleList > fdroppedParticleList
const simb::MCTrajectory & Trajectory() const
Definition: MCParticle.h:254
std::map< int, int > fParentIDMap_OrigTrackID
key is current track ID, value is parent ID – for real G4 track ID tracking only ...
GeneratedParticleIndex_t truthIndex
Index of the particle in the original generator truth record.
list_type::value_type value_type
Definition: ParticleList.h:130
int GetParentage(int trackid, bool useOrigTrackIDMap=false) const
list_type::iterator iterator
Definition: ParticleList.h:131
Geant4 interface.
STL namespace.
const sim::ParticleList * GetList() const
intermediate_table::const_iterator const_iterator
double velocity_step
GeneratedParticleIndex_t truthInfoIndex() const
Returns the index of the particle in the generator truth record.
bool empty() const
Definition: MCTrajectory.h:167
std::map< int, int > fParentIDMap
key is current track ID, value is parent ID
virtual void PreTrackingAction(const G4Track *)
G4UserTrackingAction interfaces.
decltype(auto) constexpr end(T &&obj)
ADL-aware version of std::end.
Definition: StdUtils.h:77
constexpr GeneratedParticleIndex_t NoGeneratedParticleIndex
Constant representing the absence of generator truth information.
Definition: simb.h:34
void clear()
Resets the information (does not release memory it does not own)
sim::ParticleList && YieldList()
Use Geant4&#39;s user "hooks" to maintain a list of particles generated by Geant4.
static bool isDropped(simb::MCParticle const *p)
returns whether the specified particle has been marked as dropped
bool fKeepEMShowerDaughters
whether to keep EM shower secondaries, tertiaries, etc
virtual void PostTrackingAction(const G4Track *)
static const int NoParticleId
Definition: sim.h:21
bool entra
double energy
Definition: plottest35.C:25
void SetParticleList(sim::ParticleList *p)
sim::ParticleList && YieldDroppedList()
Yields the (dropped) ParticleList accumulated during the current event.
GeneratedParticleIndex_t GetPrimaryTruthIndex(int trackId) const
Returns the index of primary truth (sim::NoGeneratorIndex if none).
double velocity_G4
Monte Carlo Simulation.
std::unique_ptr< util::PositionInVolumeFilter > fFilter
filter for particles to be kept
std::unique_ptr< sim::ParticleList > fparticleList
bool isPrimary() const
Returns whether there is a particle.
double globalTime
cet::exempt_ptr< simb::MCParticle > particle
Object representing particle.
#define MF_LOG_DEBUG(id)
static int fCurrentOrigTrackID
except for EM shower particles where it always shows the original track ID
virtual void EndOfEventAction(const G4Event *)
bool hasParticle() const
Returns whether there is a particle.
#define MF_LOG_WARNING(category)
void operator()(sim::ParticleList::value_type &particleListEntry) const
Particle list in DetSim contains Monte Carlo particle information.
Float_t track
Definition: plot.C:35
virtual void SteppingAction(const G4Step *)
G4UserSteppingAction interface.
Tools and modules for checking out the basics of the Monte Carlo.
std::map< int, GeneratedParticleIndex_t > const & GetPrimaryTruthMap() const
GeneratedParticleIndex_t MCParticleIndex() const
Returns the index of the corresponding particle in truth record.
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
std::size_t GeneratedParticleIndex_t
Type of particle index in the generator truth record (simb::MCTruth).
Definition: simb.h:30
bool keep
if there was decision to keep