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
Geant4ParticlePrint.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/Printout.h"
17 #include "DD4hep/Primitives.h"
18 #include "DD4hep/InstanceCount.h"
20 #include "DDG4/Geant4Data.h"
22 
23 // Geant4 include files
24 #include "G4Event.hh"
25 
26 // C/C++ include files
27 #include <cstring>
28 
29 using namespace std;
30 using namespace DD4hep;
31 using namespace DD4hep::Simulation;
32 
34 
36 Geant4ParticlePrint::Geant4ParticlePrint(Geant4Context* ctxt, const std::string& nam)
37  : Geant4EventAction(ctxt,nam)
38 {
39  declareProperty("OutputType",m_outputType=3);
40  declareProperty("PrintBegin",m_printBegin=false);
41  declareProperty("PrintEnd", m_printEnd=true);
42  declareProperty("PrintGen", m_printGeneration=true);
43  declareProperty("PrintHits", m_printHits=false);
45 }
46 
50 }
51 
53 void Geant4ParticlePrint::makePrintout(const G4Event* e) const {
55  if ( parts ) {
56  const ParticleMap& particles = parts->particles();
57  print("+++ ******** MC Particle Printout for event ID:%d ********",e->GetEventID());
58  if ( (m_outputType&2) != 0 ) printParticleTree(e,particles); // Tree printout....
59  if ( (m_outputType&1) != 0 ) printParticles(0,particles); // Table printout....
60  return;
61  }
62  except("No Geant4MonteCarloTruth object present in the event structure to access the particle information!", c_name());
63 }
64 
68 }
69 
71 void Geant4ParticlePrint::begin(const G4Event* e) {
72  if ( m_printBegin ) makePrintout(e);
73 }
74 
76 void Geant4ParticlePrint::end(const G4Event* e) {
77  if ( m_printEnd ) makePrintout(e);
78 }
79 
80 void Geant4ParticlePrint::printParticle(const std::string& prefix, const G4Event* e, Geant4ParticleHandle p) const {
81  char equiv[32];
82  PropertyMask mask(p->reason);
83  PropertyMask status(p->status);
84  string proc_name = p.processName();
85  string proc_type = p.processTypeName();
86  int parent_id = p->parents.empty() ? -1 : *(p->parents.begin());
87 
88  equiv[0] = 0;
89  if ( p->parents.end() == p->parents.find(p->g4Parent) ) {
90  ::snprintf(equiv,sizeof(equiv),"/%d",p->g4Parent);
91  }
92  print("+++ %s ID:%7d %12s %6d%-7s %7s %3s %5d %3s %+.3e %-4s %-7s %-3s %-3s %2d [%s%s%s] %c%c%c%c",
93  prefix.c_str(),
94  p->id,
95  p.particleName().c_str(),
96  parent_id,equiv,
99  int(p->daughters.size()),
101  p.energy(),
105  mask.isSet(G4PARTICLE_KEEP_PARENT) ? "YES" : "",
106  p.numParent(),
107  proc_name.c_str(),
108  p->process ? "/" : "",
109  proc_type.c_str(),
110  status.isSet(G4PARTICLE_GEN_EMPTY) ? 'E' : '.',
111  status.isSet(G4PARTICLE_GEN_STABLE) ? 'S' : '.',
112  status.isSet(G4PARTICLE_GEN_DECAYED) ? 'D' : '.',
113  status.isSet(G4PARTICLE_GEN_DOCUMENTATION) ? 'd' : '.'
114  );
115  if ( e && m_printHits ) {
117  G4HCofThisEvent* hc = e->GetHCofThisEvent();
118  for (int ihc=0, nhc=hc->GetNumberOfCollections(); ihc<nhc; ++ihc) {
119  G4VHitsCollection* c = hc->GetHC(ihc);
120  Geant4HitCollection* coll = dynamic_cast<Geant4HitCollection*>(c);
121  size_t nhits = coll->GetSize();
122  for(size_t i=0; i<nhits; ++i) {
123  Geant4HitData* h = coll->hit(i);
124  Geant4Tracker::Hit* trk_hit = dynamic_cast<Geant4Tracker::Hit*>(h);
125  if ( 0 != trk_hit ) {
126  Geant4HitData::Contribution& t = trk_hit->truth;
127  int trackID = t.trackID;
128  int trueID = truth->particleID(trackID);
129  if ( trueID == p->id ) {
130  print("+++ %20s %s[%d] (%+.2e,%+.2e,%+.2e)[mm]","",c->GetName().c_str(),i,
131  trk_hit->position.x(),trk_hit->position.y(),trk_hit->position.z());
132  }
133  }
134  Geant4Calorimeter::Hit* cal_hit = dynamic_cast<Geant4Calorimeter::Hit*>(h);
135  if ( 0 != cal_hit ) {
136  Geant4HitData::Contributions& contrib = cal_hit->truth;
137  for(Geant4HitData::Contributions::iterator j=contrib.begin(); j!=contrib.end(); ++j) {
139  int trackID = t.trackID;
140  int trueID = truth->particleID(trackID);
141  if ( trueID == p->id ) {
142  print("+++ %20s %s[%d] (%+.2e,%+.2e,%+.2e)[mm]","",c->GetName().c_str(),i,
143  cal_hit->position.x(),cal_hit->position.y(),cal_hit->position.z());
144  }
145  }
146  }
147  }
148  }
149  }
150 }
151 
153 void Geant4ParticlePrint::printParticles(const G4Event* e, const ParticleMap& particles) const {
154  int num_energy = 0;
155  int num_parent = 0;
156  int num_process = 0;
157  int num_primary = 0;
158  int num_secondaries = 0;
159  int num_calo_hits = 0;
160  int num_tracker_hits = 0;
161 
162  print("+++ MC Particles #Tracks:%7d ParticleType Parent/Geant4 "
163  "Primary Secondary Energy in [MeV] Calo Tracker Process/Par Details",
164  int(particles.size()));
165  for(ParticleMap::const_iterator i=particles.begin(); i!=particles.end(); ++i) {
166  Geant4ParticleHandle p = (*i).second;
167  PropertyMask mask(p->reason);
168  printParticle("MC Particle Track",e, p);
169  num_secondaries += int(p->daughters.size());
170  if ( mask.isSet(G4PARTICLE_ABOVE_ENERGY_THRESHOLD) ) ++num_energy;
171  if ( mask.isSet(G4PARTICLE_PRIMARY) ) ++num_primary;
172  if ( mask.isSet(G4PARTICLE_CREATED_TRACKER_HIT) ) ++num_tracker_hits;
173  if ( mask.isSet(G4PARTICLE_CREATED_CALORIMETER_HIT) ) ++num_calo_hits;
174  if ( mask.isSet(G4PARTICLE_KEEP_PARENT) ) ++num_parent;
175  else if ( mask.isSet(G4PARTICLE_KEEP_PROCESS) ) ++num_process;
176  }
177  print("+++ MC Particles #Tracks:%7d ParticleType Parent/Geant4 "
178  "Primary Secondary Energy Calo Tracker Process/Par",
179  int(particles.size()));
180  print("+++ MC Particle Summary: %7d %10d %7d %7d %9d %5d %6d",
181  num_primary, num_secondaries, num_energy,
182  num_calo_hits,num_tracker_hits,num_process,num_parent);
183 }
184 
185 void Geant4ParticlePrint::printParticleTree(const G4Event* e, const ParticleMap& particles, int level, Geant4ParticleHandle p) const {
186  char txt[32];
187  size_t len = sizeof(txt)-1;
188  // Ensure we do not overwrite the array
189  if ( level>int(len)-3 ) level=len-3;
190 
191  ::snprintf(txt,sizeof(txt),"%5d ",level);
192  ::memset(txt+6,' ',len-6);
193  txt[len-1] = 0;
194  txt[len-2] = '>';
195  txt[level+6]='+';
196  ::memset(txt+level+6+1,'-',len-level-3-6);
197 
198  printParticle(txt, e, p);
199  const set<int>& daughters = p->daughters;
200  // For all particles, the set of daughters must be contained in the record.
201  for(set<int>::const_iterator id=daughters.begin(); id!=daughters.end(); ++id) {
202  int id_dau = *id;
203  Geant4ParticleHandle dau = (*particles.find(id_dau)).second;
204  printParticleTree(e, particles, level+1, dau);
205  }
206 }
207 
209 void Geant4ParticlePrint::printParticleTree(const G4Event* e, const ParticleMap& particles) const {
210  print("+++ MC Particle Parent daughter relationships. [%d particles]",int(particles.size()));
211  print("+++ MC Particles %12s #Tracks:%7d %-12s Parent%-7s "
212  "Primary Secondary Energy %-8s Calo Tracker Process/Par Details",
213  "",int(particles.size()),"ParticleType","","in [MeV]");
214  for(ParticleMap::const_iterator i=particles.begin(); i!=particles.end(); ++i) {
215  Geant4ParticleHandle p = (*i).second;
216  PropertyMask mask(p->reason);
217  if ( mask.isSet(G4PARTICLE_PRIMARY) ) printParticleTree(e, particles, 0, p);
218  }
219 }
DDG4 tracker hit class used by the generic DDG4 tracker sensitive detector.
Definition: Geant4Data.h:228
void makePrintout(const G4Event *e) const
Print particle table.
void printParticle(const std::string &prefix, const G4Event *e, Geant4ParticleHandle p) const
Geant4Event & event() const
Access the geant4 event – valid only between BeginEvent() and EndEvent()!
const G4VProcess * process
not persisten. ROOT cannot handle
void printParticleTree(const G4Event *e, const ParticleMap &particles, int level, Geant4ParticleHandle p) const
Print tree of kept particles.
void printParticles(const G4Event *e, const ParticleMap &particles) const
Print record of kept particles.
std::vector< MonteCarloContrib > Contributions
Definition: Geant4Data.h:199
static void decrement(T *)
Decrement count according to type information.
DDG4 calorimeter hit class used by the generic DDG4 calorimeter sensitive detector.
Definition: Geant4Data.h:273
bool m_printGeneration
Property: Flag to indicate output type as part of the generator action.
const char * yes_no(bool value)
Helper function to print booleans in format YES/NO.
Definition: Printout.h:295
Geant4HitWrapper & hit(size_t which)
Access the hit wrapper.
int particleID(int track, bool throw_if_not_found=true) const
Access the equivalent track id (shortcut to the usage of TrackEquivalents)
Class of the Geant4 toolkit. See http://www-geant4.kek.jp/Reference.
Definition: Geant4Classes.h:40
const ParticleMap & particles() const
Access the particle map.
return e
Definition: Volumes.cpp:297
Concrete basic implementation of the Geant4 event action.
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
Contribution truth
Monte Carlo / Geant4 information.
Definition: Geant4Data.h:237
Contributions truth
Hit contributions by individual particles.
Definition: Geant4Data.h:278
void print(const char *fmt,...) const
Support for messages with variable output level using output level.
void except(const char *fmt,...) const
Support of exceptions: Print fatal message and throw runtime_error.
Geant4ParticleMap::ParticleMap ParticleMap
int m_outputType
Property: Flag to indicate output type: 1: TABLE, 2:TREE, 3:BOTH (default)
not persisten. ROOT cannot handle
Definition: Geant4Data.h:140
virtual void begin(const G4Event *event)
Pre-event action callback.
Particles parents
The list of daughters of this MC particle.
Position position
Hit position.
Definition: Geant4Data.h:231
Data structure to map particles produced during the generation and the simulation.
virtual void end(const G4Event *event)
Post-event action callback.
std::string particleName() const
Access to the Geant4 particle name.
T * extension(bool alert=true)
Access to type safe extension object. Exception is thrown if the object is invalid.
Generic hit container class using Geant4HitWrapper objects.
size_t numParent() const
Accessor to the number of particle parents.
virtual void operator()(G4Event *event)
Generation action callback.
Data structure to access derived MC particle information.
static void increment(T *)
Increment count according to type information.
Definition: InstanceCount.h:98
bool m_printBegin
Property: Flag to indicate output type at begin of event.
double energy() const
Scalar particle energy.
ReferenceBitMask< const int > PropertyMask
bool m_printHits
Property: Flag to indicate output of hit data in tree.
Generic context to extend user, run and event information.
static const double second
Definition: DD4hepUnits.h:112
bool m_printEnd
Property: Flag to indicate output type at end of event.
virtual size_t GetSize() const
Access the collection size.
std::string processName() const
Access to the creator process name.
virtual ~Geant4ParticlePrint()
Default destructor.
Geant4Action & declareProperty(const std::string &nam, T &val)
Declare property.
Definition: Geant4Action.h:358
bool isSet(const T &m) const
Definition: Primitives.h:544
Base class for geant4 hit structures used by the default DDG4 sensitive detector implementations.
Definition: Geant4Data.h:108
std::string processTypeName() const
Access to the creator process type name.