ParticleListAction.cxx
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 /// \brief Use Geant4's user "hooks" to maintain a list of particles generated by Geant4.
3 ///
4 /// \author seligman@nevis.columbia.edu, brebel@fnal.gov
5 ////////////////////////////////////////////////////////////////////////
6 // NOvA Includes
8 #include "nug4/G4Base/PrimaryParticleInformation.h"
9 #include "Simulation/Particle.h"
11 #include "Simulation/Simulation.h"
12 
13 
14 #include "nug4/G4Base/UserActionFactory.h"
15 USERACTIONREG3(g4n, ParticleListAction, g4n::ParticleListAction)
16 
17 // G4 includes
18 #include "Geant4/G4Event.hh"
19 #include "Geant4/G4Track.hh"
20 #include "Geant4/G4ThreeVector.hh"
21 #include "Geant4/G4ParticleDefinition.hh"
22 #include "Geant4/G4PrimaryParticle.hh"
23 #include "Geant4/G4DynamicParticle.hh"
24 #include "Geant4/G4VUserPrimaryParticleInformation.hh"
25 #include "Geant4/G4Step.hh"
26 #include "Geant4/G4StepPoint.hh"
27 #include "Geant4/G4VProcess.hh"
28 
29 // ROOT includes
30 #include <TLorentzVector.h>
31 
32 // C/C++ includes
33 #include <algorithm>
34 #include <string>
35 
36 // ART includes
39 
40 namespace g4n{
41 
42  // Initialize static members.
46  std::map<int, std::vector<float> > ParticleListAction::fTrueDepMap;
47 
48  //-------------------------------------------------------------
49  // Constructor.
51  bool manyParticles)
52  : fIsAborted (false)
53  , fEnergyCut (energythres*CLHEP::GeV)
54  , fManyParticles (manyParticles)
55  , fRecordFirstBrem (false)
56  , fRecordFirstPair (false)
57  , fStoreTrajOutsideBigBox(false)
58  , fisFirstBremRecorded (false)
59  , fisFirstPairRecorded (false)
60  {
61  // Create the particle list that we'll (re-)use during the course
62  // of the Geant4 simulation.
64  }
65 
66  //-------------------------------------------------------------
67  // Destructor.
69  {
70  // Delete anything that we created with "new'.
71  delete fParticleNav;
72  }
73 
74  //-------------------------------------------------------------
76  {
77  fEnergyCut = pset.get< double >("G4EnergyThreshold")*CLHEP::GeV;
78  fManyParticles = pset.get< bool >("ManyParticles");
79  fRecordFirstBrem = pset.get< bool >("RecordFirstBrem");
80  fRecordFirstPair = pset.get< bool >("RecordFirstPair");
81  fStoreTrajOutsideBigBox = pset.get< bool >("StoreTrajectoriesOutsideBigBox");
82  }
83 
84  //-------------------------------------------------------------
86  {
87  mf::LogInfo("ParticleListAction")
88  << "ParticleListAction::PrintConfig"
89  << "\n EnergyCut " << fEnergyCut
90  << "\n ManyParticles " << (fManyParticles?"true":"false")
91  << "\n RecordFirstBrem " << (fRecordFirstBrem?"true":"false")
92  << "\n RecordFirstPair " << (fRecordFirstPair?"true":"false")
93  << "\n StoreTrajOutsideBigBox " << (fStoreTrajOutsideBigBox?"true":"false");
94  }
95 
96  //-------------------------------------------------------------
98  {
99  // Clear any previous particle information.
100  fParticle = 0;
101  fParentIDMap.clear();
102  fParticleNav->clear();
103  fTrackIDToMCTruthIndex.clear();
104 
106 
107  //Karl::Clear the map for TrueDepMap
108  ParticleListAction::fTrueDepMap.clear();
109 
110  /// If it is configured to record the first of these
111  /// interactions, then set the variable to false.
112  /// Otherwise, treat it as if we already recorded it (even though we did not)
114  else fisFirstBremRecorded = true;
116  else fisFirstPairRecorded = true;
117  }
118 
119  //-------------------------------------------------------------
120  // figure out the ultimate parentage of the particle with track ID
121  // trackid
122  // assume that the current track id has already been added to
123  // the fParentIDMap
124  int ParticleListAction::GetParentage(int trackid) const
125  {
126  int parentid = sim::kNoParticleId;
127 
128  // search the fParentIDMap recursively until we have the parent id
129  // of the first EM particle that led to this one
130  std::map<int,int>::const_iterator itr = fParentIDMap.find(trackid);
131  while( itr != fParentIDMap.end() ){
132  LOG_DEBUG("ParticleListAction") << "parentage for " << trackid
133  << " " << (*itr).second;
134  // set the parentid to the current parent ID, when the loop ends
135  // this id will be the first EM particle
136  parentid = (*itr).second;
137  itr = fParentIDMap.find(parentid);
138  }
139  LOG_DEBUG("ParticleListAction") << "final parent ID " << parentid;
140 
141  return parentid;
142  }
143 
144  //-------------------------------------------------------------
145  // Create our initial sim::Particle object and add it to the sim::ParticleList.
147  {
148  // Karl::Reset the containment variables
150  fEnteredDet = fLeftDet = false;
151 
152  // get the track ID for this particle
153  const G4int trackID = track->GetTrackID() + fTrackIDOffset;
155  size_t mcTruthIndex = 0;
156 
157  // get the parent id from Geant for the current track:
158  G4int parentID = track->GetParentID() + fTrackIDOffset;
159 
160  // Is there an MCTruth object associated with this G4Track? We
161  // have to go up a "chain" of information to find out:
162  const G4ParticleDefinition* partdef = track->GetDefinition();
163  const G4int pdg = partdef->GetPDGEncoding();
164  const G4DynamicParticle* dp = track->GetDynamicParticle();
165  const G4PrimaryParticle* pp = dp->GetPrimaryParticle();
166 
167  LOG_DEBUG("ParticleListAction") << "preparing to track " << fCurrentTrackID
168  << " pdg " << pdg
169  << " with parent " << parentID;
170 
171  std::string process_name;
172  if ( pp != 0 ){
173  const G4VUserPrimaryParticleInformation* gppi = pp->GetUserInformation();
174  const g4b::PrimaryParticleInformation* ppi = dynamic_cast<const g4b::PrimaryParticleInformation*>(gppi);
175 
176  if ( ppi != 0 ){
177  // If we've made it this far, a PrimaryParticleInformation
178  // object exists and we are using a primary particle, set the
179  // process name accordingly.
180  process_name = "Primary";
181 
182  // primary particles should still have parentID = 0, even
183  // if there are multiple MCTruths for this event
184  parentID = 0;
185 
186  mcTruthIndex = ppi->MCTruthIndex();
187  }
188 
189  } // Is there a G4PrimaryParticle?
190  else{
191 
192  // figure out what process is making this track - skip it if it is
193  // one of pair production, compton scattering, photoelectric effect
194  // bremstrahlung, annihilation, any ionization - who wants to save
195  // a buttload of electrons that arent from a CC interaction?
196  process_name = track->GetCreatorProcess()->GetProcessName();
197  if( !fManyParticles && (process_name.find("conv") != std::string::npos
198  || process_name.find("LowEnConversion") != std::string::npos
199  ||(process_name.find("Pair") != std::string::npos && fisFirstPairRecorded)
200  || process_name.find("compt") != std::string::npos
201  || process_name.find("Compt") != std::string::npos
202  ||(process_name.find("Brem") != std::string::npos && fisFirstBremRecorded)
203  || process_name.find("phot") != std::string::npos
204  || process_name.find("Photo") != std::string::npos
205  || process_name.find("Ion") != std::string::npos
206  || process_name.find("annihil") != std::string::npos)){
207 
208  // find the ultimate parent of this particle that was not a secondary
209  // of the EM shower
210  // first add this track id and its parent to the fParentIDMap
212 
214 
215  LOG_DEBUG("ParticleListAction") << "current track ID " << fCurrentTrackID;
216 
217  // check that fCurrentTrackID is in the particle list - it is possible
218  // that this particle's parent is a particle that did not get tracked.
219  // An example is a parent that was made due to muMinusCaptureAtRest
220  // and the daughter was made by the phot process. The parent likely
221  // isn't saved in the particle list because it is below the energy cut
222  // which will put a bogus track id value into the sim::IDE object for
223  // the sim::SimChannel if we don't check it.
225 
227 
228  // set fParticle to 0 as we are not stepping this particle
229  // and adding trajectory points to it
230  fParticle = 0;
231 
232  LOG_DEBUG("ParticleListAction") << "killing TrackID: " << trackID << " bc EM daughter, "
233  << process_name << " " << pdg
234  << ", use track id " << fCurrentTrackID;
235 
236  return;
237  }// end if part of EM shower, but not primary particle
238 
239  // Check the energy of the particle. If it falls below the energy
240  // cut, don't add it to our list.
241  const G4double energy = track->GetKineticEnergy();
242  if ( energy < fEnergyCut ){
243  fParticle = 0;
244  LOG_DEBUG("ParticleListAction") << "killing TrackID: " << fCurrentTrackID << " energy/fEnergyCut";
245 
246  // do add the particle to the parent id map though
247  // and set the current track id to be it's ultimate parent
249 
251  return;
252  }
253 
254  // check to see if the parent particle has been stored in the particle navigator
255  // if not, then see if it is possible to walk up the fParentIDMap to find the
256  // ultimate parent of this particle. Use that ID as the parent ID for this
257  // particle
258  if( fParticleNav->find(parentID) == fParticleNav->end() ){
259  // do add the particle to the parent id map
260  // just in case it makes a daughter that we have to track as well
262  const int pid = this->GetParentage(parentID);
263 
264  // if we still can't find the parent in the particle navigator,
265  // we have to give up
266  if( fParticleNav->find(pid) == fParticleNav->end() ){
267  mf::LogWarning("ParticleListAction") << "can't find parent id: "
268  << parentID
269  << " in the particle navigator, or fParentIDMap."
270  << " Make " << parentID << " the mother ID for"
271  << " track ID " << fCurrentTrackID
272  << " in the hope that it will aid debugging.";
273  }
274  else
275  parentID = pid;
276  }
277 
278  // Attempt to find the MCTruth index corresponding to the
279  // current particle. If the fCurrentTrackID is not in the
280  // map try the parent ID, if that is not there, throw an
281  // exception
283  mcTruthIndex = fTrackIDToMCTruthIndex.at(fCurrentTrackID);
284  else if(fTrackIDToMCTruthIndex.count(parentID) > 0 )
285  mcTruthIndex = fTrackIDToMCTruthIndex.at(parentID);
286  else
287  throw cet::exception("ParticleListAction") << "Cannot find MCTruth index for track id "
288  << fCurrentTrackID << " or " << parentID;
289  }// end of if/else
290 
291 
292  // protect against potentially empty process names.
293  if ( process_name.empty() ) { process_name = "unknown";}
294  // if it's the brem photon, update the variable
295  else if( process_name.find("Brem") != std::string::npos){ fisFirstBremRecorded = true;}
296  else if( process_name.find("Pair") != std::string::npos){ fisFirstPairRecorded = true;}
297 
298 
299  // Create the sim::Particle object
300  fParticle = new sim::Particle(fCurrentTrackID, pdg, process_name,
301  parentID, dp->GetMass()/CLHEP::GeV);
302  //part.Print();
303 
304  // Polarization.
305  const G4ThreeVector& polarization = track->GetPolarization();
306  fParticle->SetPolarization( TVector3(polarization.x(),
307  polarization.y(),
308  polarization.z() ) );
309 
311 
313  LOG_WARNING("ParticleListAction") << "attempting to put " << fCurrentTrackID
314  << " into fTrackIDToMCTruthIndex map "
315  << " particle is\n" << *fParticle;
316 
317  fTrackIDToMCTruthIndex[fCurrentTrackID] = mcTruthIndex;
318 
319  return;
320  }
321 
322  //-------------------------------------------------------------
324  {
325  // Karl::Pushback new numbers for the containment variables
326  // --- I want loop through all of the parentage for this particle, and add the EDeps to them...
327  // --- Make an iterator to go through fParticleNav
329  while( ParItr != fParticleNav->end() ) {
330  sim::Particle TempPar = *((*ParItr).second);
331  int TempTrID = TempPar.TrackId();
332  // --- See if we already have an entry for this stored TrackID
333  if ( ParticleListAction::fTrueDepMap.count(TempTrID) == 0 ) {
334  std::vector<float> TempVec = { 0., 0., 0., 0. };
335  ParticleListAction::fTrueDepMap[TempTrID] = TempVec;
336  }
337  // --- Now check whether this particle entered the detector. If it did, change this parents vector!
338  if ( fEnteredDet ) {
339  // --- If looking at the actual particle, I want to set Entering & Leaving energy.
340  if ( TempTrID == track->GetTrackID() ) {
341  ParticleListAction::fTrueDepMap[ TempTrID ][0] = fThisEnStart;
342  ParticleListAction::fTrueDepMap[ TempTrID ][1] = fThisEnEscape;
343  }
344  // --- Increment the count for energy deposited in the detector
345  ParticleListAction::fTrueDepMap[ TempTrID ][2] += fThisEDep;
346  // --- Increment the count for escaping energy.
347  ParticleListAction::fTrueDepMap[ TempTrID ][3] += fThisEnEscape;
348  } else {
349  // If particle didn't enter the detector...Do I event care about them then?
350  }
351  // Finally, lets see if I have a parent to this particle!
352  ParItr = fParticleNav->find( TempPar.Mother() );
353  }
354  // Karl::Pushback new numbers for the containment variables
355 
356  if( !fParticle ) return;
357  // particle final position now saved in SteppingAction when
358  // filling the TrajectoryPoint, filled using the post step point.
359 
360  //RUN_ON_PARTICLELISTACTION_DEBUG(DEBUG_POSTTRACKING, track);
361 
362  const G4ThreeVector position = track->GetPosition();
363 
364  // If the final position wasn't inside the big box, then we wouldn't have
365  // written it as a trajectory point. But we are interested in where
366  // particles stop, so do it here.
368  !geom->isInsideDetectorBigBox(position.x() / CLHEP::cm,
369  position.y() / CLHEP::cm,
370  position.z() / CLHEP::cm)){
371  G4double time = track->GetGlobalTime();
372 
373  TLorentzVector fourPos(position.x() / CLHEP::cm,
374  position.y() / CLHEP::cm,
375  position.z() / CLHEP::cm,
376  time / CLHEP::ns);
377 
378  const G4ThreeVector momentum = track->GetMomentum();
379  const G4double energy = track->GetTotalEnergy();
380  TLorentzVector fourMom(momentum.x() / CLHEP::GeV,
381  momentum.y() / CLHEP::GeV,
382  momentum.z() / CLHEP::GeV,
383  energy / CLHEP::GeV);
384 
385  fParticle->AddTrajectoryPoint( fourPos, fourMom );
386  }
387  }
388 
389  //-------------------------------------------------------------
390  // With every step, add to the particle's trajectory - not
391  // done for nova - maybe we should?
393 
394  // Karl::Containment outputs
395  // --- Work out if PreStep point was in the detector volume.
396  bool PreStepInDet = geom->isInsideDetectorBigBox( step->GetPreStepPoint()->GetPosition().x()/CLHEP::cm,
397  step->GetPreStepPoint()->GetPosition().y()/CLHEP::cm,
398  step->GetPreStepPoint()->GetPosition().z()/CLHEP::cm
399  );
400  // --- Work out if the PostStep point was in the detector volume.
401  bool PosStepInDet = geom->isInsideDetectorBigBox( step->GetPostStepPoint()->GetPosition().x()/CLHEP::cm,
402  step->GetPostStepPoint()->GetPosition().y()/CLHEP::cm,
403  step->GetPostStepPoint()->GetPosition().z()/CLHEP::cm
404  );
405  // --- If this step took/kept the particle in the detector.
406  if ( PosStepInDet ) {
407  fThisEDep += step->GetTotalEnergyDeposit()/CLHEP::GeV;
408  // --- If the particle has entered the detector for the first time with this step.
409  if ( !fEnteredDet ) {
410  fEnteredDet = true;
411  // --- Take PreStepPoint, as this is the energy it had AS it entered.
412  fThisEnStart = step->GetPreStepPoint()->GetKineticEnergy()/CLHEP::GeV; // Do I want GetTotalEnergy() or GetKineticEnergy()??
413  }
414  // --- If the particle has left the detector before, and now come back into the detector I want to subtract this EDep
415  // --- from the amount of energy which I thought was not deposited in the detector.
416  if ( fLeftDet ) {
417  fThisEnEscape -= step->GetTotalEnergyDeposit()/CLHEP::GeV;
418  }
419  } else {
420  // --- If this is the first time the particle has left the detector with this step.
421  if ( PreStepInDet && !fLeftDet ) {
422  fLeftDet = true;
423  // --- Take PreStepPoint, as this is the energy it had AS it left.
424  fThisEnEscape = step->GetPreStepPoint()->GetKineticEnergy()/CLHEP::GeV; // Do I want GetTotalEnergy() or GetKineticEnergy??
425  }
426  }
427  // Karl::Containment outputs
428 
429  if ( !fParticle ) return;
430 
431  // For the most part, we just want to add the post-step
432  // information to the particle's trajectory. There's one
433  // exception: In PreTrackingAction, the correct time information
434  // is not available. So add the correct vertex information here.
435 
436  if ( fParticle->NumberTrajectoryPoints() == 0 ){
437  // Get the pre-step information from the G4Step.
438  const G4StepPoint* preStepPoint = step->GetPreStepPoint();
439 
440  const G4ThreeVector position = preStepPoint->GetPosition();
441  G4double time = preStepPoint->GetGlobalTime();
442 
443  // Remember that LArSoft uses cm, ns, GeV.
444  TLorentzVector fourPos(position.x() / CLHEP::cm,
445  position.y() / CLHEP::cm,
446  position.z() / CLHEP::cm,
447  time / CLHEP::ns);
448 
449  const G4ThreeVector momentum = preStepPoint->GetMomentum();
450  const G4double energy = preStepPoint->GetTotalEnergy();
451  TLorentzVector fourMom(momentum.x() / CLHEP::GeV,
452  momentum.y() / CLHEP::GeV,
453  momentum.z() / CLHEP::GeV,
454  energy / CLHEP::GeV);
455 
456  // Add the first point in the trajectory.
457  fParticle->AddTrajectoryPoint( fourPos, fourMom );
458  }
459 
460  // Get the post-step information from the G4Step.
461  const G4StepPoint* postStepPoint = step->GetPostStepPoint();
462 
463  const G4ThreeVector position = postStepPoint->GetPosition();
464  G4double time = postStepPoint->GetGlobalTime();
465 
466  // Remember that LArSoft uses cm, ns, GeV.
467  TLorentzVector fourPos(position.x() / CLHEP::cm,
468  position.y() / CLHEP::cm,
469  position.z() / CLHEP::cm,
470  time / CLHEP::ns );
471 
472  const G4ThreeVector momentum = postStepPoint->GetMomentum();
473  const G4double energy = postStepPoint->GetTotalEnergy();
474  TLorentzVector fourMom(momentum.x() / CLHEP::GeV,
475  momentum.y() / CLHEP::GeV,
476  momentum.z() / CLHEP::GeV,
477  energy / CLHEP::GeV );
478 
479  // Add another point in the trajectory.
481  geom->isInsideDetectorBigBox(position.x() / CLHEP::cm,
482  position.y() / CLHEP::cm,
483  position.z() / CLHEP::cm)){
484  fParticle->AddTrajectoryPoint( fourPos, fourMom );
485  }
486  }
487 
488  //-------------------------------------------------------------
489  std::vector<sim::Particle> ParticleListAction::GetList() const
490  {
491  // fill a vector<sim::Particle> using the sim::ParticleNavigator
492  std::vector<sim::Particle> plist;
493 
494  // make sure to set the fTrackIDOffset with the highest G4 track id + 100
495  // the 100 gives some cushion between lists
496  int highestID = 0;
497  for(sim::ParticleNavigator::iterator itr = fParticleNav->begin(); itr != fParticleNav->end(); ++itr){
498  plist.push_back(*((*itr).second));
499  if((*itr).first > highestID) highestID = (*itr).first;
500  }
501 
502  fTrackIDOffset = highestID + 100;
503 
504  LOG_DEBUG("ParticleListAction") << *fParticleNav << "\ntrack id offset is now " << fTrackIDOffset;
505 
506  return plist;
507  }
508 
509  //-------------------------------------------------------------
510  std::map<int, size_t> ParticleListAction::TrackIDToMCTruthIndexMap() const
511  {
512  return fTrackIDToMCTruthIndex;
513  }
514 
515  //-------------------------------------------------------------
516  /// Utility class for the EndOfEventAction method: update the
517  /// daughter relationships in the particle list.
518  class UpdateDaughterInformation : public std::unary_function<sim::ParticleNavigator::value_type, void>
519  {
520  public:
522  : particleNav(0)
523  {}
524  void SetParticleNav( sim::ParticleNavigator* p ) { particleNav = p; }
526  {
527  // We're looking at this Particle in the list.
528  sim::Particle* particle = particleNavEntry.second;
529 
530  // The parent ID of this particle.
531  int parentID = particle->Mother();
532 
533  // If the parentID <= 0, this is a primary particle.
534  if ( parentID <= 0 ) return;
535 
536  // If we get here, this particle is somebody's daughter. Add
537  // it to the list of daughter particles for that parent.
538 
539  // Get the parent particle from the list.
540  sim::ParticleNavigator::iterator parentEntry = particleNav->find( parentID );
541 
542  if ( parentEntry == particleNav->end() ){
543  // We have an "orphan": a particle whose parent isn't
544  // recorded in the particle list. This is not signficant;
545  // it's possible for a particle not to be saved in the list
546  // because it failed an energy cut, but for it to have a
547  // daughter that passed the cut (e.g., a nuclear decay).
548  return;
549  }
550 
551  // Add the current particle to the daughter list of the
552  // parent.
553  sim::Particle* parent = (*parentEntry).second;
554  parent->AddDaughter( particle->TrackId() );
555  }
556  private:
558  };
559 
560  //-------------------------------------------------------------
561  // There's one last thing to do: All the particles have their
562  // parent IDs set (in PostTrackingAction), but we haven't set the
563  // daughters yet. That's done in this method.
565  {
566  // If the event was aborted, we want to keep track of that fact
567  // by setting this flag
568  if (e->IsAborted()) fIsAborted = true;
569 
570  // Set up the utility class for the "for_each" algorithm. (We only
571  // need a separate set-up for the utility class because we need to
572  // give it the pointer to the particle list. We're using the STL
573  // "for_each" instead of the C++ "for loop" because it's supposed
574  // to be faster.
575  UpdateDaughterInformation updateDaughterInformation;
576  updateDaughterInformation.SetParticleNav( fParticleNav );
577 
578  // Update the daughter information for each particle in the list.
579  std::for_each(fParticleNav->begin(),
580  fParticleNav->end(),
581  updateDaughterInformation);
582 
583  LOG_DEBUG("ParticleListAction") << *fParticleNav;
584  }
585 
586 } // end namespace
TruthSlim – remove generated objects that don&#39;t contribute.
#define LOG_DEBUG(stream)
Definition: Messenger.h:149
static sim::Particle * fParticle
unsigned int NumberTrajectoryPoints() const
Definition: MCParticle.h:217
void operator()(sim::ParticleNavigator::value_type &particleNavEntry)
void AddDaughter(const int trackID)
Definition: MCParticle.h:264
float fThisEnStart
The energy which this particle has the first time in entered the active volume.
void AddTrajectoryPoint(TLorentzVector const &position, TLorentzVector const &momentum)
Definition: MCParticle.h:256
float fThisEDep
The energy which this particle deposited in the active volume.
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
void PreTrackingAction(const G4Track *)
int Mother() const
Definition: MCParticle.h:212
const char * p
Definition: xmltok.h:285
Definition: event.h:19
static constexpr double ns
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
static std::map< int, std::vector< float > > fTrueDepMap
list_type::iterator iterator
void SetParticleNav(sim::ParticleNavigator *p)
void Config(fhicl::ParameterSet const &pset)
int GetParentage(int trackid) const
int TrackId() const
Definition: MCParticle.h:209
bool fLeftDet
Boolean to state whether the particle has left the active volume before (after entering it)...
bool fIsAborted
Whether Geant4 was aborted while tracking a particle.
void PrintConfig(std::string const &opt)
static const int kNoParticleId
Definition: Simulation.h:14
static constexpr double cm
Definition: SystemOfUnits.h:99
std::map< int, size_t > fTrackIDToMCTruthIndex
map track ID to index of MCTruth in input list
ParticleListAction(double energythresh=0.0, bool manyParticles=false)
void SetPolarization(const TVector3 &p)
Definition: MCParticle.h:265
void PostTrackingAction(const G4Track *)
void EndOfEventAction(const G4Event *)
T get(std::string const &key) const
Definition: ParameterSet.h:231
double energy
Definition: plottest35.C:25
bool fRecordFirstPair
When true, saves the electron generated by the first Pair production.
static constexpr Double_t GeV
Definition: Munits.h:291
bool isInsideDetectorBigBox(const double x_cm, const double y_cm, const double z_cm) const
Is the particle inside the detector Big Box?
void BeginOfEventAction(const G4Event *)
G4double fEnergyCut
The minimum energy for a particle to be included in the list.
bool fRecordFirstBrem
When true, saves the gamma generated by the first Brem.
#define LOG_WARNING(category)
Int_t parentID
Definition: plot.C:85
std::vector< sim::Particle > GetList() const
sim::ParticleNavigator * fParticleNav
bool fManyParticles
When true, includes track ids from processes like compt and brem.
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
void SteppingAction(const G4Step *)
static int fTrackIDOffset
offset added to track ids when running over
bool fStoreTrajOutsideBigBox
Otherwise trajectory details outside the big box are thrown away to save space.
ParticleNavigator Add(const int &offset) const
std::map< int, size_t > TrackIDToMCTruthIndexMap() const
Int_t trackID
Definition: plot.C:84
std::map< int, int > fParentIDMap
key is current track ID, value is parent ID
Float_t e
Definition: plot.C:35
bool fEnteredDet
Boolean to state whether the particle has entered the active volume before.
static constexpr double GeV
iterator find(const key_type &key)
sim::ParticleNavigator * particleNav
list_type::value_type value_type
float fThisEnEscape
The energy which this particle had the last time it left the active volume.
art::ServiceHandle< geo::Geometry > geom
enum BeamMode string