13 #include "cetlib_except/exception.h" 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" 27 #include <TLorentzVector.h> 58 bool storeTrajectories,
59 bool keepEMShowerDaughters,
60 bool keepMCParticleList,
61 bool storeDroppedMCParticles)
93 const std::map<int, int>* parentIDMap =
99 while (itr != parentIDMap->end()) {
100 MF_LOG_DEBUG(
"ParticleListAction") <<
"parentage for " << trackid <<
" " << (*itr).second;
104 parentid = (*itr).second;
105 itr = parentIDMap->find(parentid);
107 MF_LOG_DEBUG(
"ParticleListAction") <<
"final parent ID " << parentid;
117 G4ParticleDefinition* particleDefinition = track->GetDefinition();
118 G4int pdgCode = particleDefinition->GetPDGEncoding();
137 std::string process_name =
"unknown";
138 bool drop_shower_daughter =
false;
142 const G4DynamicParticle* dynamicParticle = track->GetDynamicParticle();
143 const G4PrimaryParticle* primaryParticle = dynamicParticle->GetPrimaryParticle();
145 if (primaryParticle != 0) {
146 const G4VUserPrimaryParticleInformation* gppi = primaryParticle->GetUserInformation();
155 process_name =
"primary";
170 process_name = track->GetCreatorProcess()->GetProcessName();
171 drop_shower_daughter =
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) {
204 G4double
energy = track->GetKineticEnergy();
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.";
247 double mass = dynamicParticle->GetMass() / CLHEP::GeV;
259 const G4ThreeVector& polarization = track->GetPolarization();
261 TVector3(polarization.x(), polarization.y(), polarization.z()));
265 if (drop_shower_daughter) {
284 aTrack->GetStep()->GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName();
324 globalTime = step->GetTrack()->GetGlobalTime();
326 velocity_step = step->GetStepLength() / step->GetDeltaTime();
327 if ((step->GetTrack()->GetDefinition()->GetPDGEncoding() == 0) &&
331 step->GetPostStepPoint()->SetGlobalTime(
globalTime - step->GetDeltaTime() +
343 const G4StepPoint* preStepPoint = step->GetPreStepPoint();
345 const G4ThreeVector position = preStepPoint->GetPosition();
346 G4double time = preStepPoint->GetGlobalTime();
349 TLorentzVector fourPos(position.x() / CLHEP::cm,
350 position.y() / CLHEP::cm,
351 position.z() / CLHEP::cm,
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);
372 G4String process = step->GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName();
373 G4bool ignoreProcess = process.contains(
"LArVoxel") || process.contains(
"OpDetReadout");
376 <<
": DEBUG - process='" << process <<
"'" 377 <<
" ignoreProcess=" << ignoreProcess <<
" fstoreTrajectories=" <<
fstoreTrajectories;
383 if (fstoreTrajectories && !ignoreProcess) {
389 const G4StepPoint* postStepPoint = step->GetPostStepPoint();
391 const G4ThreeVector position = postStepPoint->GetPosition();
392 G4double time = postStepPoint->GetGlobalTime();
395 TLorentzVector fourPos(position.x() / CLHEP::cm,
396 position.y() / CLHEP::cm,
397 position.z() / CLHEP::cm,
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);
422 int particleID = particleListEntry.first;
427 int parentID = particleList->GetMotherOf(particleID);
430 if (parentID <= 0)
return;
438 if (parentEntry == particleList->end()) {
446 if (!parentEntry->second)
return;
489 if ((*pn).first > highestID) highestID = (*pn).first;
495 if ((*pn).first > highestID) highestID = (*pn).first;
518 <<
"ParticleListAction::YieldList(): particle list not build by user request.\n";
525 if ((*pn).first > highestID) highestID = (*pn).first;
531 if ((*pn).first > highestID) highestID = (*pn).first;
546 <<
"ParticleListAction::YieldDroppedList(): dropped particle list not build by user " 554 if ((*pn).first > highestID) highestID = (*pn).first;
559 if ((*pn).first > highestID) highestID = (*pn).first;
568 TLorentzVector
const& mom,
569 std::string
const& process)
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)
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
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
int GetParentage(int trackid, bool useOrigTrackIDMap=false) const
list_type::iterator iterator
const sim::ParticleList * GetList() const
GeneratedParticleIndex_t truthInfoIndex() const
Returns the index of the particle in the generator truth record.
ParticleInfo_t fCurrentParticle
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.
constexpr GeneratedParticleIndex_t NoGeneratedParticleIndex
Constant representing the absence of generator truth information.
void clear()
Resets the information (does not release memory it does not own)
sim::ParticleList && YieldList()
Use Geant4's user "hooks" to maintain a list of particles generated by Geant4.
static int fTrackIDOffset
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
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).
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.
cet::exempt_ptr< simb::MCParticle > particle
Object representing particle.
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)
Particle list in DetSim contains Monte Carlo particle information.
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
std::size_t GeneratedParticleIndex_t
Type of particle index in the generator truth record (simb::MCTruth).
bool keep
if there was decision to keep
static int fCurrentTrackID