From 803b12abd75445e5749cc2779d38e6d7be9eb129 Mon Sep 17 00:00:00 2001
From: Markus Frank <markus.frank@cern.ch>
Date: Mon, 22 Sep 2014 18:34:38 +0000
Subject: [PATCH] Fix floating point exceptions and add generic input reader
template
---
DDG4/CMakeLists.txt | 1 +
DDG4/include/DDG4/Geant4InputAction.h | 119 ++++++++++++
DDG4/include/DDG4/Geant4ParticlePrint.h | 2 +-
DDG4/plugins/Geant4Factories.cpp | 4 +
DDG4/plugins/Geant4HepEventReader.cpp | 243 ++++++++++++++++++++++++
DDG4/src/Geant4Data.cpp | 1 -
DDG4/src/Geant4InputAction.cpp | 130 +++++++++++++
DDG4/src/Geant4ParticleHandler.cpp | 6 +-
DDG4/src/Geant4ParticlePrint.cpp | 17 +-
DDG4/src/Geant4PrimaryHandler.cpp | 3 +-
10 files changed, 514 insertions(+), 12 deletions(-)
create mode 100644 DDG4/include/DDG4/Geant4InputAction.h
create mode 100644 DDG4/plugins/Geant4HepEventReader.cpp
create mode 100644 DDG4/src/Geant4InputAction.cpp
diff --git a/DDG4/CMakeLists.txt b/DDG4/CMakeLists.txt
index 11b762163..78b7a760f 100644
--- a/DDG4/CMakeLists.txt
+++ b/DDG4/CMakeLists.txt
@@ -72,6 +72,7 @@ install(FILES
python/DD4hep.py
python/DDG4.py
python/DDG4Dict.C
+ python/SystemOfUnits.py
DESTINATION python)
install(TARGETS DD4hepG4 DD4hepG4Plugins DD4hepG4Legacy g4gdmlDisplay
diff --git a/DDG4/include/DDG4/Geant4InputAction.h b/DDG4/include/DDG4/Geant4InputAction.h
new file mode 100644
index 000000000..ca34c3b35
--- /dev/null
+++ b/DDG4/include/DDG4/Geant4InputAction.h
@@ -0,0 +1,119 @@
+// $Id: Geant4Converter.cpp 603 2013-06-13 21:15:14Z markus.frank $
+//====================================================================
+// AIDA Detector description implementation for LCD
+//--------------------------------------------------------------------
+//
+//====================================================================
+#ifndef DD4HEP_DDG4_GEANT4INPUTACTION_H
+#define DD4HEP_DDG4_GEANT4INPUTACTION_H
+
+// Framework include files
+#include "DDG4/Geant4Particle.h"
+#include "DDG4/Geant4GeneratorAction.h"
+
+// C/C++ include files
+#include <vector>
+
+// Forward declarations
+class G4Event;
+
+/*
+ * DD4hep namespace declaration
+ */
+namespace DD4hep {
+
+ /*
+ * Simulation namespace declaration
+ */
+ namespace Simulation {
+
+ /** @class Geant4EventReader Geant4EventReader.h DDG4/Geant4EventReader.h
+ *
+ * Base class to read input files containing simulation data.
+ *
+ * @author P.Kostka (main author)
+ * @author M.Frank (code reshuffeling into new DDG4 scheme)
+ * @version 1.0
+ */
+ class Geant4EventReader {
+
+ public:
+ typedef Geant4Particle Particle;
+ protected:
+ /// File name
+ std::string m_name;
+ /// Flag if direct event access is supported
+ bool m_directAccess;
+
+ public:
+ /// Initializing constructor
+ Geant4EventReader(const std::string& nam);
+ /// Default destructor
+ virtual ~Geant4EventReader();
+ /// File name
+ const std::string& name() const { return m_name; }
+ /// Flag if direct event access (by event sequence number) is supported (Default: false)
+ bool hasDirectAccess() const { return m_directAccess; }
+ /// Read an event and fill a vector of MCParticles.
+ virtual int readParticles(int event_number, std::vector<Particle*>& particles) = 0;
+ };
+
+ /** @class Geant4InputAction Geant4InputAction.h DDG4/Geant4InputAction.h
+ *
+ * Concrete implementation of the Geant4 generator action base class
+ * populating Geant4 primaries from Geant4 and HepStd files.
+ *
+ * @author P.Kostka (main author)
+ * @author M.Frank (code reshuffeling into new DDG4 scheme)
+ * @version 1.0
+ */
+ class Geant4InputAction : public Geant4GeneratorAction {
+
+ public:
+ typedef Geant4Particle Particle;
+
+ protected:
+ /// Property: input file
+ std::string m_input;
+ /// Property: SYNCEVT
+ int m_firstEvent;
+ /// Property; interaction mask
+ int m_mask;
+ /// Property: Momentum downscaler for debugging
+ double m_momScale;
+ /// Event reader object
+ Geant4EventReader* m_reader;
+
+ /// Read an event and return a LCCollectionVec of MCParticles.
+ int readParticles(int event_number, std::vector<Particle*>& particles);
+ /// helper to report Geant4 exceptions
+ std::string issue(int i) const;
+
+ public:
+ /// Standard constructor
+ Geant4InputAction(Geant4Context* context, const std::string& name);
+ /// Default destructor
+ virtual ~Geant4InputAction();
+ /// Callback to generate primary particles
+ virtual void operator()(G4Event* event);
+ };
+ } /* End namespace Simulation */
+} /* End namespace DD4hep */
+
+#include "DD4hep/Plugins.h"
+namespace {
+ /// Factory to create Geant4 physics constructions
+ template <typename P> class Factory<P, DD4hep::Simulation::Geant4EventReader*(std::string)> { public:
+ static void Func(void *ret, void*, const std::vector<void*>& a, void*)
+ { *(DD4hep::Simulation::Geant4EventReader**)ret = (DD4hep::Simulation::Geant4EventReader*)new P(*(std::string*)a[0]);}
+ };
+}
+/// Plugin defintion to create event reader objects
+#define DECLARE_GEANT4_EVENT_READER(name) \
+ PLUGINSVC_FACTORY_WITH_ID(name,std::string(#name),DD4hep::Simulation::Geant4EventReader*(std::string))
+
+/// Plugin defintion to create event reader objects
+#define DECLARE_GEANT4_EVENT_READER_NS(ns,name) typedef ns::name __##name##__; \
+ PLUGINSVC_FACTORY_WITH_ID(__##name##__,std::string(#name),DD4hep::Simulation::Geant4EventReader*(std::string))
+
+#endif /* DD4HEP_DDG4_GEANT4INPUTACTION_H */
diff --git a/DDG4/include/DDG4/Geant4ParticlePrint.h b/DDG4/include/DDG4/Geant4ParticlePrint.h
index 4429e408b..a9a44a697 100644
--- a/DDG4/include/DDG4/Geant4ParticlePrint.h
+++ b/DDG4/include/DDG4/Geant4ParticlePrint.h
@@ -53,7 +53,7 @@ namespace DD4hep {
/// Print tree of kept particles
void printParticleTree(const ParticleMap& particles) const;
/// Print particle table
- void makePrintout() const;
+ void makePrintout(int event_id) const;
public:
diff --git a/DDG4/plugins/Geant4Factories.cpp b/DDG4/plugins/Geant4Factories.cpp
index e4cfc7fa9..fa91f307e 100644
--- a/DDG4/plugins/Geant4Factories.cpp
+++ b/DDG4/plugins/Geant4Factories.cpp
@@ -91,6 +91,10 @@ DECLARE_GEANT4ACTION(Geant4InteractionMerger)
#include "DDG4/Geant4PrimaryHandler.h"
DECLARE_GEANT4ACTION(Geant4PrimaryHandler)
+//=============================
+#include "DDG4/Geant4InputAction.h"
+DECLARE_GEANT4ACTION(Geant4InputAction)
+
//=============================
#include "DDG4/Geant4TestActions.h"
namespace DD4hep { namespace Simulation {
diff --git a/DDG4/plugins/Geant4HepEventReader.cpp b/DDG4/plugins/Geant4HepEventReader.cpp
new file mode 100644
index 000000000..09918713b
--- /dev/null
+++ b/DDG4/plugins/Geant4HepEventReader.cpp
@@ -0,0 +1,243 @@
+// $Id: Geant4Converter.cpp 603 2013-06-13 21:15:14Z markus.frank $
+//====================================================================
+// AIDA Detector description implementation for LCD
+//--------------------------------------------------------------------
+//
+//====================================================================
+
+// Framework include files
+#include "DDG4/Geant4InputAction.h"
+
+// C/C++ include files
+#include <fstream>
+
+/*
+ * DD4hep namespace declaration
+ */
+namespace DD4hep {
+ /*
+ * lcio namespace declaration
+ */
+ namespace Simulation {
+
+ /** @class HepEventReader HepEventReader.h DDG4/HepEventReader.h
+ *
+ * Class to populate Geant4 primary particles and vertices from a
+ * file in StdHep format (ASCII)
+ *
+ * @author P.Kostka (main author)
+ * @author M.Frank (code reshuffeling into new DDG4 scheme)
+ * @version 1.0
+ */
+ struct Geant4HepEventReader : public Geant4EventReader {
+ protected:
+ std::ifstream m_input;
+ int m_format;
+
+ public:
+ /// Initializing constructor
+ Geant4HepEventReader(const std::string& nam);
+ /// Default destructor
+ virtual ~Geant4HepEventReader();
+ /// Read an event and fill a vector of MCParticles.
+ virtual int readParticles(int event_number, std::vector<Particle*>& particles);
+ };
+ } /* End namespace Simulation */
+} /* End namespace DD4hep */
+
+// $Id: Geant4Converter.cpp 603 2013-06-13 21:15:14Z markus.frank $
+//====================================================================
+// AIDA Detector description implementation for LCD
+//--------------------------------------------------------------------
+//
+//====================================================================
+// #include "DDG4/Geant4HepEventReader"
+
+#include "CLHEP/Units/SystemOfUnits.h"
+#include "CLHEP/Units/PhysicalConstants.h"
+
+#include <cerrno>
+
+using namespace CLHEP;
+using namespace DD4hep::Simulation;
+typedef DD4hep::ReferenceBitMask<int> PropertyMask;
+
+#define HEPEvt 1
+
+// Factory entry
+DECLARE_GEANT4_EVENT_READER(Geant4HepEventReader)
+
+/// Initializing constructor
+Geant4HepEventReader::Geant4HepEventReader(const std::string& nam)
+: Geant4EventReader(nam), m_input(), m_format(HEPEvt)
+{
+ // Need to set format from input specification!!!
+
+ // Now open the input file:
+}
+
+/// Default destructor
+Geant4HepEventReader::~Geant4HepEventReader() {
+ m_input.close();
+}
+
+/// Read an event and fill a vector of MCParticles.
+int Geant4HepEventReader::readParticles(int /* event_number */, std::vector<Particle*>& particles) {
+ //static const double c_light = 299.792;// mm/ns
+ //
+ // Read the event, check for errors
+ //
+ int NHEP; // number of entries
+ m_input >> NHEP;
+ if( m_input.eof() ) {
+ //
+ // End of File :: ??? Exception ???
+ // -> FG: EOF is not an exception as it happens for every file at the end !
+ return EIO;
+ }
+ //
+ // Loop over particles
+ int ISTHEP; // status code
+ int IDHEP; // PDG code
+ int JMOHEP1; // first mother
+ int JMOHEP2; // last mother
+ int JDAHEP1; // first daughter
+ int JDAHEP2; // last daughter
+ double PHEP1; // px in GeV/c
+ double PHEP2; // py in GeV/c
+ double PHEP3; // pz in GeV/c
+ double PHEP4; // energy in GeV
+ double PHEP5; // mass in GeV/c**2
+ double VHEP1; // x vertex position in mm
+ double VHEP2; // y vertex position in mm
+ double VHEP3; // z vertex position in mm
+ double VHEP4; // production time in mm/c
+
+ std::vector<int> daughter1;
+ std::vector<int> daughter2;
+
+ for( int IHEP=0; IHEP<NHEP; IHEP++ ) {
+ if ( m_format == HEPEvt)
+ m_input >> ISTHEP >> IDHEP >> JDAHEP1 >> JDAHEP2
+ >> PHEP1 >> PHEP2 >> PHEP3 >> PHEP5;
+ else
+ m_input >> ISTHEP >> IDHEP
+ >> JMOHEP1 >> JMOHEP2
+ >> JDAHEP1 >> JDAHEP2
+ >> PHEP1 >> PHEP2 >> PHEP3
+ >> PHEP4 >> PHEP5
+ >> VHEP1 >> VHEP2 >> VHEP3
+ >> VHEP4;
+
+ if(m_input.eof())
+ return EIO;
+ //
+ // Create a MCParticle and fill it from stdhep info
+ Particle* p = new Particle(IHEP);
+ PropertyMask status(p->status);
+ //
+ // PDGID
+ p->pdgID = IDHEP;
+ //
+ // Momentum vector
+ p->pex = p->psx = PHEP1*GeV;
+ p->pey = p->psy = PHEP2*GeV;
+ p->pez = p->psz = PHEP3*GeV;
+ //
+ // Mass
+ p->mass = PHEP5*GeV;
+ //
+ // Vertex
+ // (missing information in HEPEvt files)
+ p->vsx = 0.0;
+ p->vsy = 0.0;
+ p->vsz = 0.0;
+ p->vex = 0.0;
+ p->vey = 0.0;
+ p->vez = 0.0;
+ //
+ // Generator status
+ // Simulator status 0 until simulator acts on it
+ p->status = 0;
+ if ( ISTHEP == 0 ) status.set(G4PARTICLE_GEN_EMPTY);
+ if ( ISTHEP == 1 ) status.set(G4PARTICLE_GEN_STABLE);
+ if ( ISTHEP == 2 ) status.set(G4PARTICLE_GEN_DECAYED);
+ if ( ISTHEP == 3 ) status.set(G4PARTICLE_GEN_DOCUMENTATION);
+ //
+ // Creation time (note the units [1/c_light])
+ // (No information in HEPEvt files)
+ p->time = 0.0;
+ p->properTime = 0.0;
+ //
+ // Keep daughters information for later
+ daughter1.push_back(JDAHEP1);
+ daughter2.push_back(JDAHEP2);
+ //
+ // Add the particle to the collection vector
+ particles.push_back(p);
+ //
+ }// End loop over particles
+
+ //
+ // Now make a second loop over the particles, checking the daughter
+ // information. This is not always consistent with parent
+ // information, and this utility assumes all parents listed are
+ // parents and all daughters listed are daughters
+ //
+ for( int IHEP=0; IHEP<NHEP; IHEP++ ) {
+ struct ParticleHandler {
+ Particle* m_part;
+ ParticleHandler(Particle* p) : m_part(p) {}
+ void addParent(const Particle* p) {
+ m_part->parents.insert(p->id);
+ }
+ Particle* findParent(const Particle* p) {
+ return m_part->parents.find(p->id)==m_part->parents.end() ? 0 : m_part;
+ }
+ };
+
+ //
+ // Get the MCParticle
+ //
+ Particle* mcp = particles[IHEP];
+ //
+ // Get the daughter information, discarding extra information
+ // sometimes stored in daughter variables.
+ //
+ int fd = daughter1[IHEP] - 1;
+ int ld = daughter2[IHEP] - 1;
+
+ //
+ // As with the parents, look for range, 2 discreet or 1 discreet
+ // daughter.
+ if( (fd > -1) && (ld > -1) ) {
+ if(ld >= fd) {
+ for(int id=fd;id<ld+1;id++) {
+ //
+ // Get the daughter, and see if it already lists this particle as
+ // a parent. If not, add this particle as a parent
+ //
+ ParticleHandler part(particles[id]);
+ if ( !part.findParent(mcp) ) part.addParent(mcp);
+ }
+ }
+ else { // Same logic, discrete cases
+ ParticleHandler part_fd(particles[fd]);
+ if ( !part_fd.findParent(mcp) ) part_fd.addParent(mcp);
+
+ ParticleHandler part_ld(particles[ld]);
+ if ( !part_ld.findParent(mcp) ) part_ld.addParent(mcp);
+ }
+ }
+ else if(fd > -1) {
+ ParticleHandler part(particles[fd]);
+ if ( !part.findParent(mcp) ) part.addParent(mcp);
+ }
+ else if(ld > -1) {
+ ParticleHandler part(particles[ld]);
+ if ( !part.findParent(mcp) ) part.addParent(mcp);
+ }
+ } // End second loop over particles
+ return 0;
+}
+
diff --git a/DDG4/src/Geant4Data.cpp b/DDG4/src/Geant4Data.cpp
index a4a114bfd..c255e06dc 100644
--- a/DDG4/src/Geant4Data.cpp
+++ b/DDG4/src/Geant4Data.cpp
@@ -62,7 +62,6 @@ Geant4HitData::~Geant4HitData() {
/// Extract the MC contribution for a given hit from the step information
Geant4HitData::Contribution Geant4HitData::extractContribution(G4Step* step) {
Geant4StepHandler h(step);
- G4Track* trk = step->GetTrack();
double deposit =
(h.trackDef() == G4OpticalPhoton::OpticalPhotonDefinition()) ? h.trkEnergy() : h.totalEnergy();
const G4ThreeVector& pre = h.prePosG4();
diff --git a/DDG4/src/Geant4InputAction.cpp b/DDG4/src/Geant4InputAction.cpp
new file mode 100644
index 000000000..f833caa6c
--- /dev/null
+++ b/DDG4/src/Geant4InputAction.cpp
@@ -0,0 +1,130 @@
+// $Id: Geant4Converter.cpp 603 2013-06-13 21:15:14Z markus.frank $
+//====================================================================
+// AIDA Detector description implementation for LCD
+//--------------------------------------------------------------------
+//
+// @author P.Kostka (main author)
+// @author M.Frank (code reshuffeling into new DDG4 scheme)
+//
+//====================================================================
+
+// Framework include files
+#include "DD4hep/Printout.h"
+#include "DD4hep/Primitives.h"
+#include "DDG4/Geant4InputAction.h"
+#include "DDG4/Geant4Primary.h"
+#include "DDG4/Geant4Context.h"
+
+#include "G4Event.hh"
+
+using namespace std;
+using namespace DD4hep::Simulation;
+typedef DD4hep::ReferenceBitMask<int> PropertyMask;
+
+/// Initializing constructor
+Geant4EventReader::Geant4EventReader(const std::string& nam)
+ : m_name(nam), m_directAccess(false)
+{
+}
+
+/// Default destructor
+Geant4EventReader::~Geant4EventReader() {
+}
+
+/// Standard constructor
+Geant4InputAction::Geant4InputAction(Geant4Context* context, const string& name)
+ : Geant4GeneratorAction(context,name), m_reader(0)
+{
+ declareProperty("Input", m_input);
+ declareProperty("Sync", m_firstEvent=0);
+ declareProperty("Mask", m_mask = 0);
+ declareProperty("MomentumScale", m_momScale = 1.0);
+ m_needsControl = true;
+}
+
+/// Default destructor
+Geant4InputAction::~Geant4InputAction() {
+}
+
+/// helper to report Geant4 exceptions
+string Geant4InputAction::issue(int i) const {
+ stringstream str;
+ str << "Geant4InputAction[" << name() << "]: Event " << i << " ";
+ return str.str();
+}
+
+/// Read an event and return a LCCollection of MCParticles.
+int Geant4InputAction::readParticles(int evt_number, std::vector<Particle*>& particles) {
+ int evid = evt_number + m_firstEvent;
+ if ( 0 == m_reader ) {
+ if ( m_input.empty() ) {
+ throw runtime_error("InputAction: No inoput file declared!");
+ }
+ string err;
+ TypeName tn = TypeName::split(m_input,"|");
+ try {
+ m_reader = PluginService::Create<Geant4EventReader*>(tn.first,tn.second);
+ if ( 0 == m_reader ) {
+ abortRun(issue(evid)+"Error creating reader plugin.",
+ "Failed to create file reader of type %s. Cannot open dataset %s",
+ tn.first.c_str(),tn.second.c_str());
+ return 0;
+ }
+ }
+ catch(const exception& e) {
+ err = e.what();
+ }
+ if ( !err.empty() ) {
+ abortRun(issue(evid)+err,"Error when reading file %s",m_input.c_str());
+ return 0;
+ }
+ }
+ int status = m_reader->readParticles(evid,particles);
+ if ( 0 != status ) {
+ abortRun(issue(evid)+"Error when reading file - may be end of file.",
+ "Error when reading file %s",m_input.c_str());
+ }
+ return status;
+}
+
+/// Callback to generate primary particles
+void Geant4InputAction::operator()(G4Event* event) {
+ vector<Particle*> primaries;
+ Geant4Event& evt = context()->event();
+ Geant4PrimaryEvent* prim = evt.extension<Geant4PrimaryEvent>();
+ Geant4PrimaryInteraction* inter = new Geant4PrimaryInteraction();
+
+ int result = readParticles(event->GetEventID(),primaries);
+ if ( result != 0 ) { // handle I/O error, but how?
+ return;
+ }
+ prim->add(m_mask, inter);
+ // check if there is at least one particle
+ if ( primaries.empty() ) return;
+
+ print("+++ Particle interaction with %d generator particles ++++++++++++++++++++++++",
+ int(primaries.size()));
+ Geant4Vertex* vtx = new Geant4Vertex();
+ vtx->x = 0;
+ vtx->y = 0;
+ vtx->z = 0;
+ vtx->time = 0;
+ inter->vertices.insert(make_pair(m_mask,vtx));
+ // build collection of MCParticles
+ for(size_t i=0; i<primaries.size(); ++i ) {
+ Geant4ParticleHandle p(primaries[i]);
+ const double mom_scale = m_momScale;
+ PropertyMask status(p->status);
+ p->psx = mom_scale*p->psx;
+ p->psy = mom_scale*p->psy;
+ p->psz = mom_scale*p->psz;
+ if ( p->parents.size() == 0 ) {
+ if ( status.isSet(G4PARTICLE_GEN_EMPTY) || status.isSet(G4PARTICLE_GEN_DOCUMENTATION) )
+ vtx->in.insert(p->id); // Beam particles and primary quarks etc.
+ else
+ vtx->out.insert(p->id); // Stuff, to be given to Geant4 together with daughters
+ }
+ inter->particles.insert(make_pair(p->id,p));
+ p.dump3(outputLevel()-1,name(),"+->");
+ }
+}
diff --git a/DDG4/src/Geant4ParticleHandler.cpp b/DDG4/src/Geant4ParticleHandler.cpp
index 2b559ff09..67d3eb05a 100644
--- a/DDG4/src/Geant4ParticleHandler.cpp
+++ b/DDG4/src/Geant4ParticleHandler.cpp
@@ -427,12 +427,14 @@ void Geant4ParticleHandler::rebaseSimulatedTracks(int ) {
ParticleMap::const_iterator ipar;
while( (ipar=m_particleMap.find(g4_equiv)) == m_particleMap.end() ) {
TrackEquivalents::const_iterator iequiv = m_equivalentTracks.find(g4_equiv);
- if ( iequiv == iend )
+ if ( iequiv == iend ) {
break; // ERROR !! Will be handled by printout below because ipar==end()
+ }
g4_equiv = (*iequiv).second;
}
- if ( ipar != m_particleMap.end() )
+ if ( ipar != m_particleMap.end() ) {
equivalents[(*i).first] = (*ipar).second->id; // requires (1) !
+ }
else
error("+++ No Equivalent particle for track:%d last known is:%d",(*i).second,g4_equiv);
}
diff --git a/DDG4/src/Geant4ParticlePrint.cpp b/DDG4/src/Geant4ParticlePrint.cpp
index 3c7e9dc6b..af5b4e3ef 100644
--- a/DDG4/src/Geant4ParticlePrint.cpp
+++ b/DDG4/src/Geant4ParticlePrint.cpp
@@ -14,6 +14,8 @@
#include "DDG4/Geant4Data.h"
#include <cstring>
+#include "G4Event.hh"
+
using namespace std;
using namespace DD4hep;
using namespace DD4hep::Simulation;
@@ -37,10 +39,11 @@ Geant4ParticlePrint::~Geant4ParticlePrint() {
}
/// Print particle table
-void Geant4ParticlePrint::makePrintout() const {
+void Geant4ParticlePrint::makePrintout(int event_id) const {
Geant4ParticleMap* parts = context()->event().extension<Geant4ParticleMap>();
if ( parts ) {
const ParticleMap& particles = parts->particles();
+ print("+++ ******** MC Particle Printout for event ID:%d ********",event_id);
if ( (m_outputType&2) != 0 ) printParticleTree(particles); // Tree printout....
if ( (m_outputType&1) != 0 ) printParticles(particles); // Table printout....
return;
@@ -49,18 +52,18 @@ void Geant4ParticlePrint::makePrintout() const {
}
/// Generation action callback
-void Geant4ParticlePrint::operator()(G4Event* ) {
- if ( m_printGeneration ) makePrintout();
+void Geant4ParticlePrint::operator()(G4Event* e) {
+ if ( m_printGeneration ) makePrintout(e->GetEventID());
}
/// Pre-event action callback
-void Geant4ParticlePrint::begin(const G4Event* ) {
- if ( m_printBegin ) makePrintout();
+void Geant4ParticlePrint::begin(const G4Event* e) {
+ if ( m_printBegin ) makePrintout(e->GetEventID());
}
/// Post-event action callback
-void Geant4ParticlePrint::end(const G4Event* ) {
- if ( m_printEnd ) makePrintout();
+void Geant4ParticlePrint::end(const G4Event* e) {
+ if ( m_printEnd ) makePrintout(e->GetEventID());
}
void Geant4ParticlePrint::printParticle(const std::string& prefix, Geant4ParticleHandle p) const {
diff --git a/DDG4/src/Geant4PrimaryHandler.cpp b/DDG4/src/Geant4PrimaryHandler.cpp
index d60067efa..bdbd15065 100644
--- a/DDG4/src/Geant4PrimaryHandler.cpp
+++ b/DDG4/src/Geant4PrimaryHandler.cpp
@@ -63,7 +63,8 @@ Primaries getRelevant(set<int>& visited,
const Geant4Particle::Particles& dau = p->daughters;
int first_daughter = *(dau.begin());
Geant4ParticleHandle dp = pm[first_daughter];
- double me = p->mass / p.energy();
+ double en = p.energy();
+ double me = en > std::numeric_limits<double>::epsilon() ? p->mass / en : 0.0;
// fix by S.Morozov for real != 0
double proper_time = fabs(dp->time-p->time) * me;
double proper_time_Precision = pow(10.,-DBL_DIG)*me*fmax(fabs(p->time),fabs(dp->time));
--
GitLab