From 87ea56c102bbc1161ae58bff6ffe413df7fa5383 Mon Sep 17 00:00:00 2001
From: Markus Frank <markus.frank@cern.ch>
Date: Tue, 8 Mar 2016 14:00:07 +0000
Subject: [PATCH] Extend checks in the CLID random numbers example
---
DDG4/include/DDG4/Geant4Random.h | 6 +++
DDG4/python/DDG4.py | 5 +-
DDG4/python/DDG4Dict.C | 10 ++++
DDG4/src/Geant4Exec.cpp | 1 +
DDG4/src/Geant4IsotropeGenerator.cpp | 8 +--
DDG4/src/Geant4Random.cpp | 70 +++++++++++++++++---------
examples/CLICSiD/scripts/CLICRandom.py | 51 ++++++++++++++++---
7 files changed, 115 insertions(+), 36 deletions(-)
diff --git a/DDG4/include/DDG4/Geant4Random.h b/DDG4/include/DDG4/Geant4Random.h
index 6322ee7f5..1085153a4 100644
--- a/DDG4/include/DDG4/Geant4Random.h
+++ b/DDG4/include/DDG4/Geant4Random.h
@@ -133,6 +133,12 @@ namespace DD4hep {
/** Basic random generator functions */
+ /// Create flat distributed random numbers in the interval ]0,1] calling CLHEP
+ /** This is more or less a test function, since the result should be
+ * identical to calling rndm.
+ */
+ double rndm_clhep();
+
/// Create flat distributed random numbers in the interval ]0,1]
double rndm(int i=0);
/// Create a float array of flat distributed random numbers in the interval ]0,1]
diff --git a/DDG4/python/DDG4.py b/DDG4/python/DDG4.py
index 3ee3ac22e..15c9e50f8 100644
--- a/DDG4/python/DDG4.py
+++ b/DDG4/python/DDG4.py
@@ -37,10 +37,10 @@ except Exception,X:
print '| %-100s |'%('Try to compile AClick on the fly.',)
print '+--%-100s--+'%(100*'-',)
DD4hep = compileAClick(dictionary='DDG4Dict.C',g4=True)
+from ROOT import CLHEP as CLHEP
Core = DD4hep
Sim = DD4hep.Simulation
Simulation = DD4hep.Simulation
-
Kernel = Sim.KernelHandle
Interface = Sim.Geant4ActionCreation
LCDD = Geo.LCDD
@@ -225,6 +225,9 @@ _import_class('Sim','Geant4UserParticleHandler')
_import_class('Sim','Geant4UserInitialization')
_import_class('Sim','Geant4DetectorConstruction')
_import_class('Sim','Geant4GeneratorWrapper')
+_import_class('Sim','Geant4Random')
+_import_class('CLHEP','HepRandom')
+_import_class('CLHEP','HepRandomEngine')
#---------------------------------------------------------------------------
def _get(self, name):
diff --git a/DDG4/python/DDG4Dict.C b/DDG4/python/DDG4Dict.C
index dc3829f91..d55afefae 100644
--- a/DDG4/python/DDG4Dict.C
+++ b/DDG4/python/DDG4Dict.C
@@ -333,6 +333,16 @@ namespace {
#endif
+// CLHEP stuff
+#include "CLHEP/Random/Random.h"
+#if defined(__CINT__) || defined(__MAKECINT__) || defined(__CLING__) || defined(__ROOTCLING__)
+
+#pragma link C++ namespace CLHEP;
+#pragma link C++ class CLHEP::HepRandom;
+#pragma link C++ class CLHEP::HepRandomEngine;
+#endif
+
+
int Geant4Dict() {
return 0;
}
diff --git a/DDG4/src/Geant4Exec.cpp b/DDG4/src/Geant4Exec.cpp
index 253f90288..4eb99beaa 100644
--- a/DDG4/src/Geant4Exec.cpp
+++ b/DDG4/src/Geant4Exec.cpp
@@ -522,6 +522,7 @@ int Geant4Exec::configure(Geant4Kernel& kernel) {
/// Initialize the engine etc.
rndm->initialize();
}
+ Geant4Random::setMainInstance(rndm);
kernel.executePhase("configure",0);
// Construct the default run manager
diff --git a/DDG4/src/Geant4IsotropeGenerator.cpp b/DDG4/src/Geant4IsotropeGenerator.cpp
index b6e5bbbea..f8157a4bd 100644
--- a/DDG4/src/Geant4IsotropeGenerator.cpp
+++ b/DDG4/src/Geant4IsotropeGenerator.cpp
@@ -40,7 +40,7 @@ Geant4IsotropeGenerator::~Geant4IsotropeGenerator() {
/// Uniform particle distribution
void Geant4IsotropeGenerator::getParticleDirectionUniform(int, ROOT::Math::XYZVector& direction, double& momentum) const {
- Geant4Event& evt = context()->event();
+ Geant4Event& evt = context()->event();
Geant4Random& rnd = evt.random();
double phi = m_phiMin+(m_phiMax-m_phiMin)*rnd.rndm();
double theta = m_thetaMin+(m_thetaMax-m_thetaMin)*rnd.rndm();
@@ -54,7 +54,7 @@ void Geant4IsotropeGenerator::getParticleDirectionUniform(int, ROOT::Math::XYZVe
/// Particle distribution ~ cos(theta)
void Geant4IsotropeGenerator::getParticleDirectionCosTheta(int, ROOT::Math::XYZVector& direction, double& momentum) const {
- Geant4Event& evt = context()->event();
+ Geant4Event& evt = context()->event();
Geant4Random& rnd = evt.random();
double phi = m_phiMin+(m_phiMax-m_phiMin)*rnd.rndm();
double cos_theta = std::cos(m_thetaMin)+(std::cos(m_thetaMax)-std::cos(m_thetaMin))*rnd.rndm();
@@ -72,7 +72,7 @@ void Geant4IsotropeGenerator::getParticleDirectionEta(int, ROOT::Math::XYZVector
struct Distribution {
static double eta(double x) { return -1.0*std::log(std::tan(x/2.0)); }
};
- Geant4Event& evt = context()->event();
+ Geant4Event& evt = context()->event();
Geant4Random& rnd = evt.random();
// See https://en.wikipedia.org/wiki/Pseudorapidity
const double dmin = std::numeric_limits<double>::epsilon();
@@ -94,7 +94,7 @@ void Geant4IsotropeGenerator::getParticleDirectionFFbar(int, ROOT::Math::XYZVect
static double ffbar(double x) { double c = std::cos(x); return 1 + c*c; }
//static double integral(double x) { return 1.5*x + sin(2.*x)/4.0; }
};
- Geant4Event& evt = context()->event();
+ Geant4Event& evt = context()->event();
Geant4Random& rnd = evt.random();
double phi = m_phiMin+(m_phiMax-m_phiMin)*rnd.rndm();
diff --git a/DDG4/src/Geant4Random.cpp b/DDG4/src/Geant4Random.cpp
index ad663c8bf..c0fce221f 100644
--- a/DDG4/src/Geant4Random.cpp
+++ b/DDG4/src/Geant4Random.cpp
@@ -87,9 +87,16 @@ Geant4Random::Geant4Random(Geant4Context* ctxt, const std::string& nam)
/// Default destructor
Geant4Random::~Geant4Random() {
- if ( s_instance == this ) s_instance = 0;
- if ( !m_engineType.empty() ) deletePtr(m_engine);
+ // Only delete the engine if it is NOT the CLEP default one
+ // BUT: Just cannot delete the engine. Causes havoc with static destructors!
+ //CLHEP::HepRandomEngine* curr = CLHEP::HepRandom::getTheEngine();
+ //if ( !m_engineType.empty() && m_engine != curr ) deletePtr(m_engine);
+
+ // Set gRandom to the old value
if ( m_rootRandom == gRandom ) gRandom = m_rootOLD;
+ // Reset instance pointer
+ if ( s_instance == this ) s_instance = 0;
+ // Finally delete the TRandom instance wrapper
deletePtr(m_rootRandom);
InstanceCount::decrement(this);
}
@@ -103,26 +110,31 @@ Geant4Random* Geant4Random::instance(bool throw_exception) {
}
/// Make this random generator instance the one used by Geant4
- Geant4Random* Geant4Random::setMainInstance(Geant4Random* ptr) {
- if ( s_instance != ptr ) {
- if ( !ptr ) {
- throw runtime_error("Attempt to declare invalid Geant4Random instance.");
- }
- if ( !ptr->m_inited ) {
- throw runtime_error("Attempt to declare uninitialized Geant4Random instance.");
- }
- Geant4Random* old = s_instance;
- if ( ptr->m_engine != CLHEP::HepRandom::getTheEngine() ) {
- CLHEP::HepRandom::setTheEngine(ptr->m_engine);
- }
- if ( ptr->m_replace ) {
- ptr->m_rootOLD = gRandom;
- gRandom = ptr->m_rootRandom;
- }
- s_instance = ptr;
- return old;
- }
- return 0;
+Geant4Random* Geant4Random::setMainInstance(Geant4Random* ptr) {
+ if ( ptr && !ptr->m_inited ) {
+ ptr->initialize();
+ }
+ if ( s_instance != ptr ) {
+ if ( !ptr ) {
+ throw runtime_error("Attempt to declare invalid Geant4Random instance.");
+ }
+ if ( !ptr->m_inited ) {
+ throw runtime_error("Attempt to declare uninitialized Geant4Random instance.");
+ }
+ Geant4Random* old = s_instance;
+ CLHEP::HepRandomEngine* curr = CLHEP::HepRandom::getTheEngine();
+ if ( ptr->m_engine != curr ) {
+ ptr->printP2("Moving CLHEP random instance from %p to %p",curr,ptr->m_engine);
+ CLHEP::HepRandom::setTheEngine(ptr->m_engine);
+ }
+ if ( ptr->m_replace ) {
+ ptr->m_rootOLD = gRandom;
+ gRandom = ptr->m_rootRandom;
+ }
+ s_instance = ptr;
+ return old;
+ }
+ return 0;
}
#include "CLHEP/Random/DualRand.h"
@@ -216,11 +228,13 @@ void Geant4Random::showStatus() const {
if ( !m_file.empty() )
printP2(" Created from file: %s",m_file.c_str());
else if ( !m_engineType.empty() )
- printP2(" Special instance created of type:%s",m_engineType.c_str());
+ printP2(" Special instance created of type:%s @ 0x%p",
+ m_engineType.c_str(),m_engine);
else
- printP2(" Reused HepRandom instance @ 0x%p",m_engine);
+ printP2(" Reused HepRandom engine instance %s @ 0x%p",
+ m_engine ? m_engine->name().c_str() : "???", m_engine);
- if ( m_engine != CLHEP::HepRandom::getTheEngine() )
+ if ( m_engine == CLHEP::HepRandom::getTheEngine() )
printP2(" Instance is identical to Geant4's HepRandom instance.");
printP2(" Instance is %sidentical to ROOT's gRandom instance.",
@@ -232,6 +246,12 @@ void Geant4Random::showStatus() const {
m_engine->showStatus();
}
+/// Create flat distributed random numbers in the interval ]0,1]
+double Geant4Random::rndm_clhep() {
+ if ( !m_inited ) initialize();
+ return m_engine->flat();
+}
+
/// Create flat distributed random numbers in the interval ]0,1]
double Geant4Random::rndm(int i) {
if ( !m_inited ) initialize();
diff --git a/examples/CLICSiD/scripts/CLICRandom.py b/examples/CLICSiD/scripts/CLICRandom.py
index 338a532b2..247218dbc 100644
--- a/examples/CLICSiD/scripts/CLICRandom.py
+++ b/examples/CLICSiD/scripts/CLICRandom.py
@@ -5,16 +5,55 @@
@version 1.0
"""
+from ROOT import TRandom
+from ROOT import gRandom
+
if __name__ == "__main__":
- import CLICSid
+ import CLICSid, DDG4
sid = CLICSid.CLICSid()
+ print 'DEFAULT Engine:', DDG4.CLHEP.HepRandom.getTheEngine().name()
# <<-- See this function to know how it's done....
- rndm1=sid.setupRandom('R1',seed=987654321,type='MTwistEngine')
- rndm2=sid.setupRandom('R2',seed=1234321)
+ rndm1 = sid.setupRandom('R1',seed=987654321,type='RunluxEngine')
+ print 'R1:', rndm1.name, ' Default instance:', rndm1.instance().name()
+ print ' Engine: ', rndm1.engine().name()
+ print ' DEFAULT:', DDG4.CLHEP.HepRandom.getTheEngine().name()
+ rndm1.showStatus()
+
+ rndm2 = sid.setupRandom('R2',seed=1234321,type='MTwistEngine')
+ print 'R2:', rndm2.name, ' Default instance:', rndm2.instance().name()
+ print ' Engine:', rndm2.engine().name()
+ print ' DEFAULT:', DDG4.CLHEP.HepRandom.getTheEngine().name()
+ #rndm2.showStatus()
+
+ DDG4.Geant4Random.setMainInstance(rndm1.get())
+ rndm1.showStatus()
# Move main geant random instance from rndm1 to rndm2:
# See how gRandom and HepRandom instances move
- rndm2.setMainInstance(rndm2.get())
- rndm2.showStatus()
- rndm1.showStatus()
+ DDG4.Geant4Random.setMainInstance(rndm1.get())
+ print 'DEFAULT Engine:', DDG4.CLHEP.HepRandom.getTheEngine().name()
+ print 'DDG4 DEFAULT:', DDG4.Geant4Random.instance().engine().name()
+ rndm = DDG4.Geant4Random.instance()
+
+ rndm.setSeed(1234)
+ rndm.showStatus()
+ for i in xrange(10):
+ print rndm.name(), ' -- 0 gRandome.Rndm() -- Shoot random[',i,']= ',gRandom.Rndm()
+
+ rndm.setSeed(1234)
+ for i in xrange(10):
+ print rndm.name(), ' -- 1 gRandome.Rndm() -- Shoot random[',i,']= ',gRandom.Rndm()
+
+ rndm.setSeed(1234)
+ for i in xrange(10):
+ print rndm.name(), ' -- 2 Geant4Random(CLHEP) -- Shoot random[',i,']= ',rndm.rndm_clhep()
+
+ rndm.setSeed(1234)
+ for i in xrange(10):
+ print rndm.name(), ' -- 3 Geant4Random(CLHEP) -- Shoot random[',i,']= ',rndm.rndm_clhep()
+
+ rndm.setSeed(1234)
+ for i in xrange(10):
+ print rndm.name(), ' -- 4 HepRandomEngine(CLHEP) -- Shoot random[',i,']= ',rndm.engine().flat()
+
sid.test_run(have_geo=False)
--
GitLab