DD4hep - The AIDA detector description toolkit for high energy physics experiments
DD4hep  Rev:Unversioneddirectory
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Geant4ParticleHandler.cpp
Go to the documentation of this file.
1 // $Id: $
2 //==========================================================================
3 // AIDA Detector description implementation for LCD
4 //--------------------------------------------------------------------------
5 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
6 // All rights reserved.
7 //
8 // For the licensing terms see $DD4hepINSTALL/LICENSE.
9 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
10 //
11 // Author : M.Frank
12 //
13 //==========================================================================
14 
15 // Framework include files
16 #include "DD4hep/Primitives.h"
17 #include "DD4hep/InstanceCount.h"
18 #include "DDG4/Geant4StepHandler.h"
20 #include "DDG4/Geant4EventAction.h"
26 
27 // Geant4 include files
28 #include "G4Step.hh"
29 #include "G4Track.hh"
30 #include "G4Event.hh"
31 #include "G4TrackStatus.hh"
32 #include "G4PrimaryVertex.hh"
33 #include "G4PrimaryParticle.hh"
34 #include "G4TrackingManager.hh"
35 #include "G4ParticleDefinition.hh"
36 #include "CLHEP/Units/SystemOfUnits.h"
37 
38 // C/C++ include files
39 #include <set>
40 #include <stdexcept>
41 #include <algorithm>
42 
43 using CLHEP::MeV;
44 using namespace std;
45 using namespace DD4hep;
46 using namespace DD4hep::Simulation;
47 
49 
51 Geant4ParticleHandler::Geant4ParticleHandler(Geant4Context* ctxt, const string& nam)
52 : Geant4GeneratorAction(ctxt,nam), Geant4MonteCarloTruth(), m_userHandler(0), m_primaryMap(0)
53 {
55  //generatorAction().adopt(this);
62  declareProperty("PrintEndTracking", m_printEndTracking = false);
63  declareProperty("PrintStartTracking", m_printStartTracking = false);
64  declareProperty("KeepAllParticles", m_keepAll = false);
65  declareProperty("SaveProcesses", m_processNames);
66  declareProperty("MinimalKineticEnergy",m_kinEnergyCut = 100e0*CLHEP::MeV);
67  declareProperty("MinDistToParentVertex",m_minDistToParentVertex = 2.2e-14*CLHEP::mm);//default tolerance for g4ThreeVector isNear
68  m_needsControl = true;
69 }
70 
73 }
74 
77  clear();
80 }
81 
84  return *this;
85 }
86 
89  if ( action ) {
90  if ( !m_userHandler ) {
91  Geant4UserParticleHandler* h = dynamic_cast<Geant4UserParticleHandler*>(action);
92  if ( h ) {
93  m_userHandler = h;
95  return true;
96  }
97  except("Cannot add an invalid user particle handler object [Invalid-object-type].", c_name());
98  }
99  except("Cannot add an user particle handler object [Object-exists].", c_name());
100  }
101  except("Cannot add an invalid user particle handler object [NULL-object].", c_name());
102  return false;
103 }
104 
108  m_particleMap.clear();
109  m_equivalentTracks.clear();
110 }
111 
113 void Geant4ParticleHandler::mark(const G4Track* track, int reason) {
114  if ( track ) {
115  if ( reason != 0 ) {
116  PropertyMask(m_currTrack.reason).set(reason);
117  return;
118  }
119  }
120  except("Cannot mark the G4Track if the pointer is invalid!", c_name());
121 }
122 
124 void Geant4ParticleHandler::mark(const G4Step* step_value, int reason) {
125  if ( step_value ) {
126  mark(step_value->GetTrack(),reason);
127  return;
128  }
129  except("Cannot mark the G4Track if the step-pointer is invalid!", c_name());
130 }
131 
133 void Geant4ParticleHandler::mark(const G4Step* step_value) {
134  if ( step_value ) {
135  mark(step_value->GetTrack());
136  return;
137  }
138  except("Cannot mark the G4Track if the step-pointer is invalid!", c_name());
139 }
140 
142 void Geant4ParticleHandler::mark(const G4Track* track) {
146  // If yes, flag it, because it is a candidate for removal.
147  G4LogicalVolume* vol = track->GetVolume()->GetLogicalVolume();
148  G4VSensitiveDetector* g4 = vol->GetSensitiveDetector();
149  Geant4ActionSD* sd = dynamic_cast<Geant4ActionSD*>(g4);
150  string typ = sd ? sd->sensitiveType() : string();
151  if ( typ == "calorimeter" )
153  else if ( typ == "tracker" )
155  else // Assume by default "tracker"
157 
158  //Geant4ParticleHandle(&m_currTrack).dump4(outputLevel(),vol->GetName(),"hit created by particle");
159 }
160 
162 void Geant4ParticleHandler::operator()(G4Event* event) {
163  typedef Geant4MonteCarloTruth _MC;
164  debug("+++ Event:%d Add EVENT extension of type Geant4ParticleHandler.....",event->GetEventID());
165  context()->event().addExtension((_MC*)this, typeid(_MC), 0);
166  clear();
168  if ( m_userHandler ) {
169  m_userHandler->generate(event, this);
170  }
171 }
172 
174 void Geant4ParticleHandler::step(const G4Step* step_value, G4SteppingManager* mgr) {
175  typedef vector<const G4Track*> _Sec;
176  ++m_currTrack.steps;
178  //
179  // Tracks below the energy threshold are NOT stored.
180  // If these tracks produce hits or are selected due to another signature,
181  // this criterium will anyhow take precedence.
182  //
183  const _Sec* sec=step_value->GetSecondaryInCurrentStep();
184  if ( not sec->empty() ) {
186  }
187  }
189  if ( m_userHandler ) {
190  m_userHandler->step(step_value, mgr, m_currTrack);
191  }
192 }
193 
195 void Geant4ParticleHandler::begin(const G4Track* track) {
196  Geant4TrackHandler h(track);
197  double kine = h.kineticEnergy();
198  G4ThreeVector m = h.momentum();
199  const G4ThreeVector& v = h.vertex();
200  int reason = (kine > m_kinEnergyCut) ? G4PARTICLE_ABOVE_ENERGY_THRESHOLD : 0;
201  const G4PrimaryParticle* prim = h.primary();
202  Particle* prim_part = 0;
203 
204  if ( prim ) {
205  prim_part = m_primaryMap->get(prim);
206  if ( !prim_part ) {
207  except("+++ Tracking preaction: Primary particle without generator particle!");
208  }
210  m_particleMap[h.id()] = prim_part->addRef();
211  }
212 
213  if ( prim_part ) {
214  m_currTrack.id = prim_part->id;
215  m_currTrack.reason = prim_part->reason|reason;
216  m_currTrack.mask = prim_part->mask;
217  m_currTrack.status = prim_part->status;
218  m_currTrack.spin[0] = prim_part->spin[0];
219  m_currTrack.spin[1] = prim_part->spin[1];
220  m_currTrack.spin[2] = prim_part->spin[2];
221  m_currTrack.colorFlow[0] = prim_part->colorFlow[0];
222  m_currTrack.colorFlow[1] = prim_part->colorFlow[1];
223  m_currTrack.parents = prim_part->parents;
224  m_currTrack.daughters = prim_part->daughters;
225  m_currTrack.pdgID = prim_part->pdgID;
226  m_currTrack.mass = prim_part->mass;
227  }
228  else {
230  m_currTrack.reason = reason;
231  m_currTrack.mask = 0;
233  m_currTrack.spin[0] = 0;
234  m_currTrack.spin[1] = 0;
235  m_currTrack.spin[2] = 0;
236  m_currTrack.colorFlow[0] = 0;
237  m_currTrack.colorFlow[1] = 0;
238  m_currTrack.parents.clear();
239  m_currTrack.daughters.clear();
240  m_currTrack.pdgID = h.trackDef()->GetPDGEncoding();
241  m_currTrack.mass = h.trackDef()->GetPDGMass();
243  }
244  m_currTrack.steps = 0;
249  m_currTrack.vsx = v.x();
250  m_currTrack.vsy = v.y();
251  m_currTrack.vsz = v.z();
252  m_currTrack.vex = 0.0;
253  m_currTrack.vey = 0.0;
254  m_currTrack.vez = 0.0;
255  m_currTrack.psx = m.x();
256  m_currTrack.psy = m.y();
257  m_currTrack.psz = m.z();
258  m_currTrack.pex = 0.0;
259  m_currTrack.pey = 0.0;
260  m_currTrack.pez = 0.0;
261  // If the creator process of the track is in the list of process products to be kept, set the proper flag
262  if ( m_currTrack.process ) {
263  Processes::iterator i=find(m_processNames.begin(),m_processNames.end(),m_currTrack.process->GetProcessName());
264  if ( i != m_processNames.end() ) {
266  }
267  }
268  if ( m_keepAll ) {
270  }
271 
273  if ( m_userHandler ) {
275  }
276 }
277 
279 void Geant4ParticleHandler::end(const G4Track* track) {
280  Geant4TrackHandler h(track);
282  int g4_id = h.id();
283  int track_reason = m_currTrack.reason;
285  // Update vertex end point and final momentum
286  G4ThreeVector m = track->GetMomentum();
287  const G4ThreeVector& p = track->GetPosition();
288  ph->pex = m.x();
289  ph->pey = m.y();
290  ph->pez = m.z();
291  ph->vex = p.x();
292  ph->vey = p.y();
293  ph->vez = p.z();
294 
295 
296  // Set the simulator status bits
297  PropertyMask simStatus(m_currTrack.status);
298 
299  // check if the last step ended on the worldVolume boundary
300  const G4Step* theLastStep = track->GetStep();
301  G4StepPoint* theLastPostStepPoint = NULL;
302  if(theLastStep) theLastPostStepPoint = theLastStep->GetPostStepPoint();
303  if( theLastPostStepPoint &&
304  ( theLastPostStepPoint->GetStepStatus() == fWorldBoundary //particle left world volume
305  //|| theLastPostStepPoint->GetStepStatus() == fGeomBoundary
306  )
307  ) {
309  }
310 
311  if(track->GetKineticEnergy() <= 0.) {
312  simStatus.set(G4PARTICLE_SIM_STOPPED);
313  }
314 
316  if ( m_userHandler ) {
317  m_userHandler->end(track, m_currTrack);
318  }
319 
320  // These are candate tracks with a probability to be stored due to their properties:
321  // - primary particle
322  // - hits created
323  // - secondaries
324  // - above energy threshold
325  // - to be kept due to creator process
326  //
327  if ( !mask.isNull() ) {
328  m_equivalentTracks[g4_id] = g4_id;
329  ParticleMap::iterator ip = m_particleMap.find(g4_id);
330  if ( mask.isSet(G4PARTICLE_PRIMARY) ) {
331  ph.dump2(outputLevel()-1,name(),"Add Primary",h.id(),ip!=m_particleMap.end());
332  }
333  // Create a new MC particle from the current track information saved in the pre-tracking action
334  Particle* part = 0;
335  if ( ip==m_particleMap.end() ) part = m_particleMap[g4_id] = new Particle();
336  else part = (*ip).second;
337  part->get_data(m_currTrack);
338  }
339  else {
340  // These are tracks without any special properties.
341  //
342  // We will not store them on the record, but have to memorise the
343  // track identifier in order to restore the history for the created hits.
344  int pid = m_currTrack.g4Parent;
345  m_equivalentTracks[g4_id] = pid;
346  // Need to find the last stored particle and OR this particle's mask
347  // with the mask of the last stored particle
348  TrackEquivalents::const_iterator iequiv, iend = m_equivalentTracks.end();
349  ParticleMap::iterator ip;
350  for(ip=m_particleMap.find(pid); ip == m_particleMap.end(); ip=m_particleMap.find(pid)) {
351  if ((iequiv=m_equivalentTracks.find(pid)) == iend) break; // ERROR
352  pid = (*iequiv).second;
353  }
354  if ( ip != m_particleMap.end() )
355  (*ip).second->reason |= track_reason;
356  else
357  ph.dumpWithVertex(outputLevel()+3,name(),"FATAL: No real particle parent present");
358  }
359 }
360 
362 void Geant4ParticleHandler::beginEvent(const G4Event* event) {
364  info("+++ Event %d Begin event action. Access event related information.",event->GetEventID());
366  m_globalParticleID = interaction->nextPID();
367  m_particleMap.clear();
368  m_equivalentTracks.clear();
370  if ( m_userHandler ) {
371  m_userHandler->begin(event);
372  }
373 }
374 
376 void Geant4ParticleHandler::dumpMap(const char* tag) const {
377  for(ParticleMap::const_iterator iend=m_particleMap.end(), i=m_particleMap.begin(); i!=iend; ++i) {
378  Geant4ParticleHandle((*i).second).dump4(INFO,name(),tag);
379  }
380 }
381 
383 void Geant4ParticleHandler::endEvent(const G4Event* event) {
384  int count = 0;
385  int level = outputLevel();
386  do {
387  if ( level <= VERBOSE ) dumpMap("Particle");
388  debug("+++ Iteration:%d Tracks:%d Equivalents:%d",++count,m_particleMap.size(),m_equivalentTracks.size());
389  } while( recombineParents() > 0 );
390 
391  if ( level <= VERBOSE ) dumpMap("Recombined");
392  // Rebase the simulated tracks, so that they fit to the generator particles
394  if ( level <= VERBOSE ) dumpMap("Rebased");
395  // Consistency check....
398  if ( m_userHandler ) {
399  m_userHandler->end(event);
400  }
402 
403  // Now export the data to the final record.
406  m_primaryMap = 0;
407  clear();
408 }
409 
413  TrackEquivalents equivalents, orgParticles;
414  ParticleMap finalParticles;
415  ParticleMap::const_iterator ipar, iend, i;
416  int count;
417 
419  ParticleMap& pm = interaction->particles;
420 
421  // (1.0) Copy the pre-defined particle mapping for the simulated tracks
422  // It is assumed the mapping is ZERO based without holes.
423  for(count = 0, iend=pm.end(), i=pm.begin(); i!=iend; ++i) {
424  Particle* p = (*i).second;
425  orgParticles[p->id] = p->id;
426  finalParticles[p->id] = p;
427  if ( p->id > count ) count = p->id;
429  p->addRef();
430  }
431  }
432  // (1.1) Define the new particle mapping for the simulated tracks
433  for(++count, iend=m_particleMap.end(), i=m_particleMap.begin(); i!=iend; ++i) {
434  Particle* p = (*i).second;
436  //if ( orgParticles.find(p->id) == orgParticles.end() ) {
437  orgParticles[p->id] = count;
438  finalParticles[count] = p;
439  p->id = count;
440  ++count;
441  }
442  }
443  // (2) Re-evaluate the corresponding geant4 track equivalents using the new mapping
444  for(TrackEquivalents::iterator ie=m_equivalentTracks.begin(),ie_end=m_equivalentTracks.end(); ie!=ie_end; ++ie) {
445  int g4_equiv = (*ie).first;
446  while( (ipar=m_particleMap.find(g4_equiv)) == m_particleMap.end() ) {
447  TrackEquivalents::const_iterator iequiv = m_equivalentTracks.find(g4_equiv);
448  if ( iequiv == ie_end ) {
449  break; // ERROR !! Will be handled by printout below because ipar==end()
450  }
451  g4_equiv = (*iequiv).second;
452  }
453  TrackEquivalents::mapped_type equiv = (*ie).second;
454  if ( ipar != m_particleMap.end() ) {
455  equivalents[(*ie).first] = (*ipar).second->id; // requires (1) !
456  Geant4ParticleHandle p = (*ipar).second;
457  const G4ParticleDefinition* def = p.definition();
458  int pdg = int(fabs(def->GetPDGEncoding())+0.1);
459  if ( pdg != 0 && pdg<36 && !(pdg > 10 && pdg < 17) && pdg != 22 ) {
460  error("+++ ERROR: Geant4 particle for track:%d last known is:%d -- is gluon or quark!",equiv,g4_equiv);
461  }
462  pdg = int(fabs(p->pdgID)+0.1);
463  if ( pdg != 0 && pdg<36 && !(pdg > 10 && pdg < 17) && pdg != 22 ) {
464  error("+++ ERROR(2): Geant4 particle for track:%d last known is:%d -- is gluon or quark!",equiv,g4_equiv);
465  }
466  }
467  else {
468  error("+++ No Equivalent particle for track:%d last known is:%d",equiv,g4_equiv);
469  }
470  }
471 
472  // (3) Compute the particle's parents and daughters.
473  // Replace the original Geant4 track with the
474  // equivalent particle still present in the record.
475  for(iend=m_particleMap.end(), i=m_particleMap.begin(); i!=iend; ++i) {
476  Particle* p = (*i).second;
477  if ( p->g4Parent > 0 ) {
478  int equiv_id = equivalents[p->g4Parent];
479  if ( (ipar=finalParticles.find(equiv_id)) != finalParticles.end() ) {
480  Particle* q = (*ipar).second;
481  q->daughters.insert(p->id);
482  p->parents.insert(q->id);
483  }
484  else {
485  error("+++ Inconsistency in particle record: Geant4 parent %d "
486  "of particle %d (equiv:%d) not in record!",
487  p->g4Parent,p->id,equiv_id);
488  }
489  }
490  }
491  m_equivalentTracks = equivalents;
492  m_particleMap = finalParticles;
493 }
494 
497  PropertyMask mask(particle.reason);
498  bool secondaries = mask.isSet(G4PARTICLE_HAS_SECONDARIES);
499  bool tracker_track = mask.isSet(G4PARTICLE_CREATED_TRACKER_HIT);
500  bool calo_track = mask.isSet(G4PARTICLE_CREATED_CALORIMETER_HIT);
501  bool hits_produced = mask.isSet(G4PARTICLE_CREATED_HIT);
502  bool low_energy = !mask.isSet(G4PARTICLE_ABOVE_ENERGY_THRESHOLD);
503 
505  if ( mask.isNull() || (secondaries && low_energy && !hits_produced) ) {
506  return true;
507  }
509  else if ( !hits_produced && low_energy ) {
510  return true;
511  }
513  else if ( !tracker_track && calo_track && low_energy ) {
514  return true;
515  }
516  else {
517  //printout(INFO,name(),"+++ Track: %d should be kept for no obvious reason....",id);
518  }
519  return false;
520 }
521 
525  set<int> remove;
526 
528  for(ParticleMap::reverse_iterator i=m_particleMap.rbegin(); i!=m_particleMap.rend(); ++i) {
529  Particle* p = (*i).second;
530  PropertyMask mask(p->reason);
531  // Allow the user to force the particle handling either by
532  // or the reason mask with G4PARTICLE_KEEP_USER or
533  // to set the reason mask to NULL in order to drop it.
534  //
535  // If the mask entry is set to G4PARTICLE_FORCE_KILL
536  // or is set to NULL, the particle is ALWAYS removed
537  //
538  // Note: This may override all other decisions!
539  bool remove_me = m_userHandler ? m_userHandler->keepParticle(*p) : defaultKeepParticle(*p);
540 
541  // Now look at the property mask of the particle
542  if ( mask.isNull() || mask.isSet(G4PARTICLE_FORCE_KILL) ) {
543  remove_me = true;
544  }
545  else if ( mask.isSet(G4PARTICLE_KEEP_USER) ) {
548  continue;
549  }
550  else if ( mask.isSet(G4PARTICLE_PRIMARY) ) {
552  continue;
553  }
554  else if ( mask.isSet(G4PARTICLE_KEEP_ALWAYS) ) {
555  continue;
556  }
557  else if ( mask.isSet(G4PARTICLE_KEEP_PARENT) ) {
558  //continue;
559  }
560  else if ( mask.isSet(G4PARTICLE_KEEP_PROCESS) ) {
561  ParticleMap::iterator ip = m_particleMap.find(p->g4Parent);
562  if ( ip != m_particleMap.end() ) {
563  Particle* parent_part = (*ip).second;
564  PropertyMask parent_mask(parent_part->reason);
565  if ( parent_mask.isSet(G4PARTICLE_ABOVE_ENERGY_THRESHOLD) ) {
566  parent_mask.set(G4PARTICLE_KEEP_PARENT);
567  continue;
568  }
569  }
570  // Low energy stuff. Remove it. Reassign to parent.
571  //remove_me = true;
572  }
573 
575  if ( remove_me ) {
576  int g4_id = (*i).first;
577  ParticleMap::iterator ip = m_particleMap.find(p->g4Parent);
578  remove.insert(g4_id);
579  m_equivalentTracks[g4_id] = p->g4Parent;
580  if ( ip != m_particleMap.end() ) {
581  Particle* parent_part = (*ip).second;
582  PropertyMask(parent_part->reason).set(mask.value());
583  parent_part->steps += p->steps;
584  parent_part->secondaries += p->secondaries;
586  if ( m_userHandler ) {
587  m_userHandler->combine(*p, *parent_part);
588  }
589  }
590  }
591  }
592  for(set<int>::const_iterator r=remove.begin(); r!=remove.end();++r) {
593  ParticleMap::iterator ir = m_particleMap.find(*r);
594  if ( ir != m_particleMap.end() ) {
595  (*ir).second->release();
596  m_particleMap.erase(ir);
597  }
598  }
599  return int(remove.size());
600 }
601 
604  int num_errors = 0;
605 
607  for(ParticleMap::const_iterator j, i=m_particleMap.begin(); i!=m_particleMap.end(); ++i) {
608  Geant4ParticleHandle p((*i).second);
609  PropertyMask mask(p->reason);
610  PropertyMask status(p->status);
611  set<int>& daughters = p->daughters;
612  // For all particles, the set of daughters must be contained in the record.
613  for(set<int>::const_iterator id=daughters.begin(); id!=daughters.end(); ++id) {
614  int id_dau = *id;
615  if ( (j=m_particleMap.find(id_dau)) == m_particleMap.end() ) {
616  ++num_errors;
617  error("+++ Particle:%d Daughter %d is not in particle map!",p->id,id_dau);
618  }
619  }
620  // We assume that particles from the generator have consistent parents
621  // For all other particles except the primaries, the parent must be contained in the record.
622  if ( !mask.isSet(G4PARTICLE_PRIMARY) && !status.anySet(G4PARTICLE_GEN_GENERATOR) ) {
623  TrackEquivalents::const_iterator eq_it = m_equivalentTracks.find(p->g4Parent);
624  bool in_map = false, in_parent_list = false;
625  int parent_id = -1;
626  if ( eq_it != m_equivalentTracks.end() ) {
627  parent_id = (*eq_it).second;
628  in_map = (j=m_particleMap.find(parent_id)) != m_particleMap.end();
629  in_parent_list = p->parents.find(parent_id) != p->parents.end();
630  }
631  if ( !in_map || !in_parent_list ) {
632  char parent_list[1024];
633  parent_list[0] = 0;
634  ++num_errors;
635  p.dumpWithMomentum(ERROR,name(),"INCONSISTENCY");
636  for(set<int>::const_iterator ip=p->parents.begin(); ip!=p->parents.end();++ip)
637  ::snprintf(parent_list+strlen(parent_list),sizeof(parent_list)-strlen(parent_list),"%d ",*ip);
638  error("+++ Particle:%d Parent %d (G4id:%d) In record:%s In parent list:%s [%s]",
639  p->id,parent_id,p->g4Parent,yes_no(in_map),yes_no(in_parent_list),parent_list);
640  }
641  }
642  }
643 
644  if ( num_errors > 0 ) {
645  except("+++ Consistency check failed. Found %d problems.",num_errors);
646  }
647 }
648 
650 
652  ParticleMap::const_iterator iend, i;
653  for(iend=pm.end(), i=pm.begin(); i!=iend; ++i) {
654  Particle* p = (*i).second;
655 
656  if( p->parents.empty() ) {
657  continue;
658  }
659 
660  Geant4Particle *parent(pm[ *p->parents.begin() ]);
661  const double X( parent->vex - p->vsx );
662  const double Y( parent->vey - p->vsy );
663  const double Z( parent->vez - p->vsz );
664  if( sqrt(X*X + Y*Y + Z*Z) > m_minDistToParentVertex ){
666  }
667  }
668 
669 }
PrintLevel outputLevel() const
Access the output level.
Definition: Geant4Action.h:287
Geant4ParticleMap::ParticleMap ParticleMap
Geant4Event & event() const
Access the geant4 event – valid only between BeginEvent() and EndEvent()!
TrackEquivalents m_equivalentTracks
Map associating the G4Track identifiers with identifiers of existing MCParticles. ...
virtual ~Geant4ParticleHandler()
Default destructor.
Geant4Particle * addRef()
Increase reference count.
const G4VProcess * process
not persisten. ROOT cannot handle
static const double MeV
Definition: DD4hepUnits.h:146
int nextPID()
Access a new particle identifier within the interaction.
Processes m_processNames
Property: All the processes of which the decay products will be explicitly stored.
double m_kinEnergyCut
Property: Energy cut below which particles are not collected, but assigned to the parent...
Geant4TrackingActionSequence & trackingAction() const
Access to the main tracking action sequence from the kernel object.
Data structure to store the MC particle information.
void adopt(ParticleMap &pm, TrackEquivalents &equiv)
Adopt particle maps.
void dumpWithVertex(int level, const std::string &src, const char *tag) const
Output type 3:+++ "tag" ID: 0 e- status:00000014 type: 11 Vertex:(+0.00e+00,+0.00e+00,+0.00e+00) [mm] time: +0.00e+00 [ns] #Par: 0 #Dau: 4.
void set(const T &m)
Definition: Primitives.h:535
virtual void generate(G4Event *event, Geant4ParticleHandler *handler)
Event generation action callback.
ReferenceBitMask< int > PropertyMask
Class of the Geant4 toolkit. See http://www-geant4.kek.jp/Reference.
Definition: Geant4Classes.h:48
Default Interface class to handle monte carlo truth records.
bool m_printStartTracking
Property: Steer printout at tracking action begin.
int parent() const
Track's parent identifier.
bool isNull() const
Definition: Primitives.h:554
static void decrement(T *)
Decrement count according to type information.
const char * yes_no(bool value)
Helper function to print booleans in format YES/NO.
Definition: Printout.h:295
virtual void begin(const G4Event *event)
Pre-event action callback.
Interface class to access properties of the underlying Geant4 sensitive detector structure.
ParticleMap m_particleMap
Map with stored MC Particles.
void rebaseSimulatedTracks(int base)
Rebase the simulated tracks, so that they fit to the generator particles.
void callAtEnd(Q *p, void(T::*f)(const G4Event *))
Register end-of-event callback.
void dump4(int level, const std::string &src, const char *tag) const
double kineticEnergy() const
Track's kinetic energy.
Helper class to ease the extraction of information from a G4Track object.
void call(Q *p, void(T::*f)(const G4Step *, G4SteppingManager *))
Register stepping action callback. Types Q and T must be polymorph!
Class modelling a single interaction with multiple primary vertices and particles.
Definition: Geant4Primary.h:96
Geant4ParticleMap::TrackEquivalents TrackEquivalents
void dump2(int level, const std::string &src, const char *tag, int g4id, bool inrec) const
Output type 2:+++ "tag" 20 G4: 7 def:0xde4eaa8 [gamma , gamma] reason: 20 E:+3.304035e+01 in record:Y...
void dumpMap(const char *tag) const
Debugging: Dump Geant4 particle map.
virtual void beginEvent(const G4Event *event)
Pre-event action callback.
virtual const std::string & sensitiveType() const =0
Access to the sensitive type of the detector.
return e
Definition: Volumes.cpp:297
void checkConsistency() const
Check the record consistency.
Geant4ParticleHandler user extension action called by the particle handler.
void error(const char *fmt,...) const
Support of error messages.
Data structure to manipulate a bitmask held by reference and represented by an integer.
Definition: Primitives.h:526
const char * c_name() const
Access name of the action.
Definition: Geant4Action.h:275
Geant4Context * context() const
Access the context.
Definition: Geant4Action.h:261
Data structure to map primaries to particles.
Definition: Geant4Primary.h:65
const std::string & name() const
Access name of the action.
Definition: Geant4Action.h:271
void callAtBegin(Q *p, void(T::*f)(const G4Event *))
Register begin-of-event callback.
virtual void mark(const G4Track *track)
Mark a Geant4 track to be kept for later MC truth analysis. Default flag: CREATED_HIT.
void except(const char *fmt,...) const
Support of exceptions: Print fatal message and throw runtime_error.
void releasePtr(T &p)
Helper to delete objects from heap and reset the pointer. Saves many many lines of code...
Definition: Primitives.h:316
void setVertexEndpointBit()
set the endpointIsNotVertexOfParentFlag at the end of the event
Geant4PrimaryMap * m_primaryMap
Primary map.
virtual void combine(Particle &to_be_deleted, Particle &remaining_parent)
Callback when parent should be combined.
virtual void step(const G4Step *step, G4SteppingManager *mgr, Particle &particle)
User stepping callback.
Particles parents
The list of daughters of this MC particle.
void * addExtension(void *ptr, const std::type_info &info, destruct_t dtor)
Add an extension object to the detector element.
double m_minDistToParentVertex
Property: Minimal distance after which the vertexIsNotEndpointOfParent flag is set.
virtual void end(const G4Track *track)
Post-track action callback.
bool anySet(const T &m) const
Definition: Primitives.h:547
ParticleMap particles
The map of particles participating in this primary interaction.
Data structure to map particles produced during the generation and the simulation.
long addRef()
Increase reference count.
const G4ParticleDefinition * definition() const
Access the Geant4 particle definition object (expensive!)
Geant4Action to collect the MC particle information.
Geant4SteppingActionSequence & steppingAction() const
Access to the main stepping action sequence from the kernel object.
T * extension(bool alert=true)
Access to type safe extension object. Exception is thrown if the object is invalid.
static T::const_iterator find(const T &c, const string &s)
Definition: View.cpp:34
virtual void begin(const G4Track *track)
Pre-track action callback.
void callUpFront(Q *p, void(T::*f)(const G4Track *), CallbackSequence::Location where=CallbackSequence::END)
Register Pre-track action callback before anything else.
virtual void end(const G4Event *event)
Post-event action callback.
Particle m_currTrack
Local buffer about the 'current' G4Track.
Concrete implementation of the Geant4 generator action base class.
void releaseObjects(M &m)
Definition: Primitives.h:344
Geant4Particle * get(const G4PrimaryParticle *particle)
Access DDG4 particle by G4 primary particle.
bool m_keepAll
Property: Flag to keep all particles generated.
int recombineParents()
Recombine particles and associate the to parents with cleanup.
double globalTime() const
Track global time.
virtual bool keepParticle(Particle &particle)
Callback to be answered if the particle MUST be kept during recombination step.
void callAtFinal(Q *p, void(T::*f)(const G4Track *), CallbackSequence::Location where=CallbackSequence::END)
Register Post-track action callback.
bool m_needsControl
Default property: Flag to create control instance.
Definition: Geant4Action.h:101
G4ParticleDefinition * trackDef() const
Track's particle definition.
Data structure to access derived MC particle information.
Geant4ParticleHandler & operator=(const Geant4ParticleHandler &c)
No assignment operator.
static void increment(T *)
Increment count according to type information.
Definition: InstanceCount.h:98
virtual void step(const G4Step *step, G4SteppingManager *mgr)
User stepping callback.
static bool defaultKeepParticle(Particle &particle)
Default callback to be answered if the particle should be kept if NO user handler is installed...
Generic context to extend user, run and event information.
bool m_printEndTracking
Property: Steer printout at tracking action end.
View * v
Definition: MultiView.cpp:30
int m_globalParticleID
Global particle identifier. Obtained at the begin of the event.
const G4ThreeVector & vertex() const
Track's vertex position, where the track was created.
bool adopt(Geant4Action *action)
Adopt the user particle handler.
virtual void operator()(G4Event *event)
Event generation action callback.
Geant4UserParticleHandler * m_userHandler
User action pointer.
const G4VProcess * creatorProcess() const
Physical process of the track generation.
void info(const char *fmt,...) const
Support of info messages.
Geant4Particle & get_data(Geant4Particle &c)
Assignment operator.
Geant4EventActionSequence & eventAction() const
Access to the main event action sequence from the kernel object.
static const double mm
Definition: DD4hepUnits.h:69
virtual void endEvent(const G4Event *event)
Post-event action callback.
Geant4Action & declareProperty(const std::string &nam, T &val)
Declare property.
Definition: Geant4Action.h:358
bool isSet(const T &m) const
Definition: Primitives.h:544
Default base class for all Geant 4 actions and derivates thereof.
Definition: Geant4Action.h:91
int id() const
Track's identifier.
TGeoShape TGeoMedium * m
Definition: Volumes.cpp:294
void debug(const char *fmt,...) const
Support of debug messages.
const G4PrimaryParticle * primary() const
Access the primary particle of the track object (if present)
void dumpWithMomentum(int level, const std::string &src, const char *tag) const
Output type 3:+++ <tag> ID: 0 e- status:00000014 type: 11 Vertex:(+0.00e+00,+0.00e+00,+0.00e+00) [mm] time: +0.00e+00 [ns] #Par: 0 #Dau: 4.