From db9d0e528be6874fc108b0c2647c198d231f3673 Mon Sep 17 00:00:00 2001 From: Markus Frank <markus.frank@cern.ch> Date: Wed, 28 May 2014 06:17:08 +0000 Subject: [PATCH] Please see release notes --- DDCore/src/IDDescriptor.cpp | 4 +- DDCore/src/Shapes.cpp | 6 +- DDG4/examples/CLICSidSimu.py | 15 +- DDG4/include/DDG4/Factories.h | 6 + DDG4/include/DDG4/Geant4Action.h | 12 ++ DDG4/include/DDG4/Geant4Converter.h | 5 + .../include/DDG4/Geant4DetectorConstruction.h | 16 +- DDG4/include/DDG4/Geant4Kernel.h | 18 ++ DDG4/include/DDG4/Geant4ParticleGun.h | 19 +- DDG4/include/DDG4/Geant4StepHandler.h | 25 +++ DDG4/plugins/Geant4EscapeCounter.cpp | 122 ++++++++++++ DDG4/plugins/Geant4SDActions.cpp | 46 +++-- DDG4/plugins/Geant4SensDet.cpp | 1 + DDG4/python/DDG4.py | 2 +- DDG4/src/Geant4Action.cpp | 28 +++ DDG4/src/Geant4Converter.cpp | 57 +++--- DDG4/src/Geant4DetectorConstruction.cpp | 22 ++- DDG4/src/Geant4Exec.cpp | 7 +- DDG4/src/Geant4Kernel.cpp | 31 ++- DDG4/src/Geant4OutputAction.cpp | 18 +- DDG4/src/Geant4ParticleGun.cpp | 32 ++- DDG4/src/Geant4TestActions.cpp | 40 ++-- DDG4/src/Geant4VolumeManager.cpp | 17 +- doc/release.notes | 10 + examples/CLICSiD/compact/compact.xml | 2 +- .../src/CylindricalEndcapCalorimeter_geo.cpp | 22 ++- examples/CLICSiD/src/ForwardDetector_geo.cpp | 30 ++- examples/CLICSiD/src/PolyconeSupport_geo.cpp | 17 +- .../ClientTests/compact/FCC_HcalBarrel.xml | 30 ++- .../ClientTests/src/FCC_HcalBarrel2_geo.cpp | 184 ++++++++++++++++++ .../ClientTests/src/FCC_HcalBarrel_geo.cpp | 59 +++--- .../ClientTests/src/ZylinderShell_geo.cpp | 51 +++++ 32 files changed, 779 insertions(+), 175 deletions(-) create mode 100644 DDG4/plugins/Geant4EscapeCounter.cpp create mode 100644 examples/ClientTests/src/FCC_HcalBarrel2_geo.cpp create mode 100644 examples/ClientTests/src/ZylinderShell_geo.cpp diff --git a/DDCore/src/IDDescriptor.cpp b/DDCore/src/IDDescriptor.cpp index 97677ed7b..c07163bc3 100644 --- a/DDCore/src/IDDescriptor.cpp +++ b/DDCore/src/IDDescriptor.cpp @@ -113,7 +113,9 @@ VolumeID IDDescriptor::encode(const std::vector<VolID>& ids) const { VolumeID id = 0; for (VolIds::const_iterator i = ids.begin(); i != ids.end(); ++i) { Field f = field((*i).first); - id |= f->value((*i).second << f->offset()) << f->offset(); + VolumeID vid = (*i).second; + vid = vid << f->offset(); + id |= f->value(vid) << f->offset(); } return id; } diff --git a/DDCore/src/Shapes.cpp b/DDCore/src/Shapes.cpp index 520526c42..16967a4ae 100644 --- a/DDCore/src/Shapes.cpp +++ b/DDCore/src/Shapes.cpp @@ -78,7 +78,7 @@ double Box::z() const { /// Constructor to be used when creating a new object Polycone::Polycone(double start, double delta) { - _assign(new TGeoPcon(start, delta, 0), "", "polycone", false); + _assign(new TGeoPcon(start*RAD_2_DEGREE, delta*RAD_2_DEGREE, 0), "", "polycone", false); } /// Constructor to be used when creating a new polycone object. Add at the same time all Z planes @@ -88,8 +88,8 @@ Polycone::Polycone(double start, double delta, const vector<double>& rmin, const if (rmin.size() < 2) { throw runtime_error("DD4hep: PolyCone::addZPlanes> Not enough Z planes. minimum is 2!"); } - params.push_back(start); - params.push_back(delta); + params.push_back(start*RAD_2_DEGREE); + params.push_back(delta*RAD_2_DEGREE); params.push_back(rmin.size()); for (size_t i = 0; i < rmin.size(); ++i) { params.push_back(z[i] ); diff --git a/DDG4/examples/CLICSidSimu.py b/DDG4/examples/CLICSidSimu.py index 1fc88a9f8..6551d470c 100644 --- a/DDG4/examples/CLICSidSimu.py +++ b/DDG4/examples/CLICSidSimu.py @@ -1,6 +1,7 @@ # # -import os, DDG4 +import os, time, DDG4 +from DDG4 import OutputLevel as Output from SystemOfUnits import * # # @@ -62,18 +63,18 @@ def run(): kernel.eventAction().add(evt1) kernel.eventAction().add(evt2) - + """ trk = DDG4.Action(kernel,"Geant4TrackPersistency/MonteCarloTruthHandler") mc = DDG4.Action(kernel,"Geant4MonteCarloRecordManager/MonteCarloRecordManager") kernel.registerGlobalAction(trk) kernel.registerGlobalAction(mc) trk.release() mc.release() - + """ # Configure I/O evt_root = DDG4.EventAction(kernel,'Geant4Output2ROOT/RootOutput') evt_root.Control = True - evt_root.Output = "simple.root" + evt_root.Output = "CLICSiD_"+time.strftime("%Y-%m-%d_%H-%M")+".root" evt_root.enableUI() kernel.eventAction().add(evt_root) @@ -84,10 +85,11 @@ def run(): # Setup particle gun gun = DDG4.GeneratorAction(kernel,"Geant4ParticleGun/Gun") - gun.energy = 0.5*GeV + gun.energy = 50*GeV gun.particle = 'e-' gun.multiplicity = 1 - gun.position = (0.15*mm,0.12*mm,0.1*cm) + gun.position = (0*mm,0*mm,0*cm) + gun.isotrop = True gun.enableUI() kernel.generatorAction().add(gun) """ @@ -136,6 +138,7 @@ def run(): seq = DDG4.SensitiveSequence(kernel,'Geant4SensDetActionSequence/SiTrackerEndcap') act = DDG4.SensitiveAction(kernel,'Geant4SimpleTrackerAction/SiTrackerEndcapHandler','SiTrackerEndcap') + #act.OutputLevel = Output.INFO seq.add(act) seq = DDG4.SensitiveSequence(kernel,'Geant4SensDetActionSequence/SiTrackerForward') diff --git a/DDG4/include/DDG4/Factories.h b/DDG4/include/DDG4/Factories.h index e17abb912..c7ab324a8 100644 --- a/DDG4/include/DDG4/Factories.h +++ b/DDG4/include/DDG4/Factories.h @@ -16,6 +16,10 @@ // Framework include files #include "DDG4/Defs.h" +#include "DD4hep/Printout.h" +#include "DD4hep/Primitives.h" + +// C/C++ include files #include <string> #include <map> @@ -150,6 +154,8 @@ namespace { template <typename P> class Factory<P, G4VUserPhysicsList*(DS::Geant4PhysicsListActionSequence*, int)> { public: static void Func(void *ret, void*, ARGS a, void*) { + DD4hep::printout(DD4hep::INFO,"PhysicsList","+++ Create physics list of type:%s", + DD4hep::typeName(typeid(P)).c_str()); *(G4VUserPhysicsList**) ret = (G4VUserPhysicsList*) new P((DS::Geant4PhysicsListActionSequence*) a[0], *(int*) a[1]); } }; diff --git a/DDG4/include/DDG4/Geant4Action.h b/DDG4/include/DDG4/Geant4Action.h index 4a18b3b75..e5b7fbe75 100644 --- a/DDG4/include/DDG4/Geant4Action.h +++ b/DDG4/include/DDG4/Geant4Action.h @@ -10,6 +10,7 @@ #define DD4HEP_DDG4_GEANT4ACTION_H // Framework include files +#include "DD4hep/Printout.h" #include "DDG4/Geant4Context.h" #include "DDG4/Geant4Callback.h" #include "DDG4/ComponentProperties.h" @@ -265,6 +266,12 @@ namespace DD4hep { PropertyManager& properties() { return m_properties; } + /// Access the output level + PrintLevel outputLevel() const { + return (PrintLevel)m_outputLevel; + } + /// Set the output level; returns previous value + PrintLevel setOutputLevel(PrintLevel new_level); /// Access to the UI messenger Geant4UIMessenger* control() const; /// Enable and install UI messenger @@ -287,6 +294,11 @@ namespace DD4hep { /// Install property control messenger if wanted void installPropertyMessenger(); + /// Support for messages with variable output level using output level + void print(const std::string& fmt, ...) const; + /// Support for messages with variable output level using output level + void print(const std::string& tag, const std::string& fmt, ...) const; + /// Support of debug messages. void debug(const std::string& fmt, ...) const; /// Support of info messages. diff --git a/DDG4/include/DDG4/Geant4Converter.h b/DDG4/include/DDG4/Geant4Converter.h index c23fb38ed..aff071739 100644 --- a/DDG4/include/DDG4/Geant4Converter.h +++ b/DDG4/include/DDG4/Geant4Converter.h @@ -10,6 +10,7 @@ #define DD4HEP_GEANT4CONVERTER_H // Framework include files +#include "DD4hep/Printout.h" #include "DDG4/Geant4Mapping.h" /* @@ -32,10 +33,14 @@ namespace DD4hep { class Geant4Converter : public Geometry::GeoHandler, public Geant4Mapping { public: bool m_checkOverlaps; + PrintLevel m_outputLevel; /// Initializing Constructor Geant4Converter(LCDD& lcdd); + /// Initializing Constructor + Geant4Converter(LCDD& lcdd, PrintLevel level); + /// Standard destructor virtual ~Geant4Converter(); diff --git a/DDG4/include/DDG4/Geant4DetectorConstruction.h b/DDG4/include/DDG4/Geant4DetectorConstruction.h index f4a7137ed..d67a6a9d8 100644 --- a/DDG4/include/DDG4/Geant4DetectorConstruction.h +++ b/DDG4/include/DDG4/Geant4DetectorConstruction.h @@ -8,6 +8,7 @@ #ifndef DD4HEP_GEANT4DETECTORCONSTRUCTION_H #define DD4HEP_GEANT4DETECTORCONSTRUCTION_H +#include "DD4hep/Printout.h" #include "G4VUserDetectorConstruction.hh" /* @@ -24,16 +25,25 @@ namespace DD4hep { * Simulation namespace declaration */ namespace Simulation { + class Geant4Kernel; class Geant4DetectorConstruction : public G4VUserDetectorConstruction { public: - + /// Initializing constructor for DDG4 + Geant4DetectorConstruction(Geant4Kernel& kernel); + /// Initializing constructor for other clients Geant4DetectorConstruction(Geometry::LCDD& lcdd); - virtual ~Geant4DetectorConstruction() { - } + /// Default destructor + virtual ~Geant4DetectorConstruction(); + /// Geometry construction callback G4VPhysicalVolume* Construct(); + private: + /// Printlevel used for the geometry conversion + PrintLevel m_outputLevel; + /// Reference to geometry object Geometry::LCDD& m_lcdd; + /// Reference to the world after construction G4VPhysicalVolume* m_world; }; } diff --git a/DDG4/include/DDG4/Geant4Kernel.h b/DDG4/include/DDG4/Geant4Kernel.h index 27cd4b740..ce56eeec2 100644 --- a/DDG4/include/DDG4/Geant4Kernel.h +++ b/DDG4/include/DDG4/Geant4Kernel.h @@ -10,6 +10,7 @@ #define DD4HEP_DDG4_GEANT4KERNEL_H // Framework include files +#include "DD4hep/Printout.h" #include "DDG4/Geant4Primitives.h" #include "DDG4/Geant4Action.h" @@ -114,6 +115,11 @@ namespace DD4hep { std::string m_uiName; /// Property: Number of events to be executed in batch mode long m_numEvent; + /// Property: Output level + int m_outputLevel; + /// Property: Client output levels + typedef std::map<std::string,int> ClientOutputLevels; + ClientOutputLevels m_clientLevels; /// Helper to register an action sequence template <typename C> bool registerSequence(C*& seq, const std::string& name); @@ -164,6 +170,8 @@ namespace DD4hep { PropertyManager& properties() { return m_properties; } + /// Print the property values + void printProperties() const; /// Access phase phases const Phases& phases() const { return m_phases; @@ -194,6 +202,16 @@ namespace DD4hep { bool hasProperty(const std::string& name) const; /// Access single property Property& property(const std::string& name); + /// Access the output level + PrintLevel outputLevel() const { + return (PrintLevel)m_outputLevel; + } + /// Set the global output level of the kernel object; returns previous value + PrintLevel setOutputLevel(PrintLevel new_level); + /// Fill cache with the global output level of a named object. Must be set before instantiation + void setOutputLevel(const std::string object, PrintLevel new_level); + /// Retrieve the global output level of a named object. + PrintLevel getOutputLevel(const std::string object) const; /// Register action by name to be retrieved when setting up and connecting action objects /** Note: registered actions MUST be unique. diff --git a/DDG4/include/DDG4/Geant4ParticleGun.h b/DDG4/include/DDG4/Geant4ParticleGun.h index 4a0659185..69bb5cd2a 100644 --- a/DDG4/include/DDG4/Geant4ParticleGun.h +++ b/DDG4/include/DDG4/Geant4ParticleGun.h @@ -16,6 +16,7 @@ // Forward declarations class G4ParticleDefinition; class G4ParticleGun; +class TRandom1; /* * DD4hep namespace declaration @@ -37,16 +38,22 @@ namespace DD4hep { protected: /// Position and shooting direction of the gun ROOT::Math::XYZVector m_position, m_direction; - /// Particle energy - double m_energy; - /// Desired multiplicity of the particles to be shot - int m_multiplicity; - /// Particle name - std::string m_particleName; /// Pointer to geant4 particle definition G4ParticleDefinition* m_particle; /// Pointer to the particle gun itself G4ParticleGun* m_gun; + /// Random number generator + TRandom1* m_rndm; + /// Particle energy + double m_energy; + /// Particle name + std::string m_particleName; + /// Desired multiplicity of the particles to be shot + int m_multiplicity; + /// Shot number in sequence + int m_shotNo; + /// Isotrope particles? + bool m_isotrop; public: /// Standard constructor Geant4ParticleGun(Geant4Context* context, const std::string& name); diff --git a/DDG4/include/DDG4/Geant4StepHandler.h b/DDG4/include/DDG4/Geant4StepHandler.h index f16c04973..fd774307a 100644 --- a/DDG4/include/DDG4/Geant4StepHandler.h +++ b/DDG4/include/DDG4/Geant4StepHandler.h @@ -54,6 +54,9 @@ namespace DD4hep { G4ParticleDefinition* trackDef() const { return track->GetDefinition(); } + int trkPdgID() const { + return track->GetDefinition()->GetPDGEncoding(); + } static const char* stepStatus(G4StepStatus status); const char* preStepStatus() const; const char* postStepStatus() const; @@ -73,6 +76,22 @@ namespace DD4hep { const G4ThreeVector& p = post->GetMomentum(); return Momentum(p.x(), p.y(), p.z()); } + Momentum trkMom() const { + const G4ThreeVector& p = track->GetMomentum(); + return Momentum(p.x(), p.y(), p.z()); + } + double deposit() const { + return step->GetTotalEnergyDeposit(); + } + int trkID() const { + return track->GetTrackID(); + } + double trkEnergy() const { + return track->GetTotalEnergy(); + } + double trkKineEnergy() const { + return track->GetKineticEnergy(); + } const G4VTouchable* preTouchable() const { return pre->GetTouchable(); } @@ -122,6 +141,12 @@ namespace DD4hep { G4VSensitiveDetector* postSD() const { return sd(post); } + bool firstInVolume() const { + return step->IsFirstStepInVolume(); + } + bool lastInVolume() const { + return step->IsLastStepInVolume(); + } }; } // End namespace Simulation diff --git a/DDG4/plugins/Geant4EscapeCounter.cpp b/DDG4/plugins/Geant4EscapeCounter.cpp new file mode 100644 index 000000000..9ec56ac4b --- /dev/null +++ b/DDG4/plugins/Geant4EscapeCounter.cpp @@ -0,0 +1,122 @@ +// $Id: Geant4Converter.cpp 603 2013-06-13 21:15:14Z markus.frank $ +//==================================================================== +// AIDA Detector description implementation for LCD +//-------------------------------------------------------------------- +// +// Author : M.Frank +// +//==================================================================== +#ifndef DD4HEP_DDG4_GEANT4ESCAPECOUNTER_H +#define DD4HEP_DDG4_GEANT4ESCAPECOUNTER_H + +#include "DD4hep/Detector.h" +#include "DDG4/Geant4SensDetAction.h" +#include "DDG4/Geant4SteppingAction.h" + +/* + * DD4hep namespace declaration + */ +namespace DD4hep { + + /* + * Simulation namespace declaration + */ + namespace Simulation { + + /** @class Geant4EscapeCounter Geant4EscapeCounter.h DDG4/Geant4EscapeCounter.h + * + * Measure escaping energy.... + * + * @author M.Frank + * @version 1.0 + */ + class Geant4EscapeCounter : /* virtual public Geant4SteppingAction, virtual */ public Geant4Sensitive { + /// Collection identifiers + size_t m_collectionID; + std::vector<std::string> m_detectorNames; + public: + /// Standard constructor + Geant4EscapeCounter(Geant4Context* context, const std::string& name, DetElement det, LCDD& lcdd); + /// Default destructor + virtual ~Geant4EscapeCounter(); + /// G4VSensitiveDetector interface: Method for generating hit(s) using the information of G4Step object. + virtual bool process(G4Step* step, G4TouchableHistory* history); + }; + + } // End namespace Simulation +} // End namespace DD4hep + +#endif /* DD4HEP_DDG4_GEANT4ESCAPECOUNTER_H */ + +// $Id: Geant4Converter.cpp 603 2013-06-13 21:15:14Z markus.frank $ +//==================================================================== +// AIDA Detector description implementation for LCD +//-------------------------------------------------------------------- +// +// Author : M.Frank +// +//==================================================================== +#include "DD4hep/Printout.h" +#include "DD4hep/InstanceCount.h" + +#include "DDG4/Geant4TouchableHandler.h" +#include "DDG4/Geant4TrackHandler.h" +#include "DDG4/Geant4StepHandler.h" +#include "DDG4/Geant4Mapping.h" +#include "DDG4/Geant4Data.h" + +#include "CLHEP/Units/SystemOfUnits.h" +#include "G4VProcess.hh" + +using namespace std; +using namespace CLHEP; +using namespace DD4hep; +using namespace DD4hep::Geometry; +using namespace DD4hep::Simulation; + +/// Standard constructor +Geant4EscapeCounter::Geant4EscapeCounter(Geant4Context* context, const string& nam, DetElement det, LCDD& lcdd) + // : Geant4SteppingAction(context, nam) + : Geant4Sensitive(context, nam, det, lcdd) +{ + string coll_name = name()+"Hits"; + m_needsControl = true; + declareProperty("Shells",m_detectorNames); + m_collectionID = defineCollection<SimpleTracker::Hit>(coll_name); + InstanceCount::increment(this); +} + +/// Default destructor +Geant4EscapeCounter::~Geant4EscapeCounter() { + InstanceCount::decrement(this); +} + +/// G4VSensitiveDetector interface: Method for generating hit(s) using the information of G4Step object. +bool Geant4EscapeCounter::process(G4Step* step, G4TouchableHistory* /* history */) { + typedef SimpleHit::Contribution HitContribution; + typedef vector<string> _V; + Geant4StepHandler h(step); + Geant4TrackHandler th(h.track); + Geant4TouchableHandler handler(step); + string path = handler.path(); + Position prePos = h.prePos(); + HitCollection* coll = collection(m_collectionID); + SimpleTracker::Hit* hit = new SimpleTracker::Hit(th.id(),th.pdgID(),h.deposit(),th.time()); + hit->cellID = volumeID(step); + hit->energyDeposit = th.energy(); + hit->position = prePos; + hit->momentum = h.trkMom(); + hit->length = 0; + coll->add(hit); + mark(h.track); + + print(name(),"+++ Track:%4d %8.2f MeV [%s] %s Geant4 path:%s", + h.trkID(),h.trkEnergy()/MeV,th.name().c_str(),th.creatorName().c_str(),path.c_str()); + // Kill track, so that it does no longer participate in the propagation + h.track->SetTrackStatus(fStopAndKill); + return true; +} + +#include "DDG4/Factories.h" +DECLARE_GEANT4SENSITIVE(Geant4EscapeCounter) + diff --git a/DDG4/plugins/Geant4SDActions.cpp b/DDG4/plugins/Geant4SDActions.cpp index 36ade55af..9cee4fb23 100644 --- a/DDG4/plugins/Geant4SDActions.cpp +++ b/DDG4/plugins/Geant4SDActions.cpp @@ -136,17 +136,13 @@ namespace DD4hep { double new_len = mean_length(h.preMom(),h.postMom())/hit_len; direction *= new_len/hit_len; } - printout(DEBUG,"SimpleTracker","%s> Add hit with deposit:%f Pos:%f %f %f", - c_name(),step->GetTotalEnergyDeposit(),position.X(),position.Y(),position.Z()); - Hit* hit = new Hit(h.track->GetTrackID(), - h.track->GetDefinition()->GetPDGEncoding(), - step->GetTotalEnergyDeposit(), - h.track->GetGlobalTime()); + print("SimpleTracker","%s> Add hit with deposit:%7.2f MeV Pos:%8.2f %8.2f %8.2f", + c_name(),step->GetTotalEnergyDeposit(),position.X(),position.Y(),position.Z()); + Hit* hit = new Hit(h.trkID(), h.trkPdgID(), h.deposit(), h.track->GetGlobalTime()); if ( hit ) { HitContribution contrib = Hit::extractContribution(step); - hit->cellID = volumeID( step ) ; - - hit->energyDeposit = contrib.deposit ; + hit->cellID = volumeID(step); + hit->energyDeposit = contrib.deposit; hit->position = position; hit->momentum = direction; hit->length = hit_len; @@ -156,11 +152,11 @@ namespace DD4hep { hit->cellID = volumeID( step ) ; throw runtime_error("Invalid CELL ID for hit!"); } - printout(INFO,"SimpleTracker","%s> Hit with deposit:%f Pos:%f %f %f ID=%016X", - c_name(),step->GetTotalEnergyDeposit(),position.X(),position.Y(),position.Z(), - (void*)hit->cellID); + print("SimpleTracker","%s> Hit with deposit:%f Pos:%f %f %f ID=%016X", + c_name(),step->GetTotalEnergyDeposit(),position.X(),position.Y(),position.Z(), + (void*)hit->cellID); Geant4TouchableHandler handler(step); - printout(INFO,"SimpleTracker","%s> Geant4 path:%s",c_name(),handler.path().c_str()); + print("SimpleTracker","%s> Geant4 path:%s",c_name(),handler.path().c_str()); return true; } throw runtime_error("new() failed: Cannot allocate hit object"); @@ -182,21 +178,23 @@ namespace DD4hep { Position pos = 0.5 * (h.prePos() + h.postPos()); HitContribution contrib = Hit::extractContribution(step); HitCollection* coll = collection(m_collectionID); - Hit* hit = coll->find<Hit>(PositionCompare<Hit>(pos)); + Hit* hit = 0;//coll->find<Hit>(PositionCompare<Hit>(pos)); if ( !hit ) { - hit = new Hit(pos) ; - hit->cellID = volumeID( step ) ; - coll->add(hit) ; - printout(DEBUG,"SimpleTracker","%s> CREATE hit with deposit:%f Pos:%f %f %f", - c_name(),contrib.deposit,pos.X(),pos.Y(),pos.Z()); + Geant4TouchableHandler handler(step); + //hit = new Hit(pos); + hit = new Hit(h.prePos()); + hit->cellID = volumeID(step); + coll->add(hit); + print("SimpleCalorimeter","%s> CREATE hit with deposit:%7.3f MeV Pos:%8.2f %8.2f %8.2f %s", + c_name(),contrib.deposit,pos.X(),pos.Y(),pos.Z(),handler.path().c_str()); if ( 0 == hit->cellID ) { - hit->cellID = volumeID( step ) ; + hit->cellID = volumeID(step); throw runtime_error("Invalid CELL ID for hit!"); } } else { - printout(DEBUG,"SimpleTracker","%s> UPDATE hit with deposit:%f Pos:%f %f %f", - c_name(),contrib.deposit,pos.X(),pos.Y(),pos.Z()); + print("SimpleCalorimeter","%s> UPDATE hit with deposit:%7.3f MeV Pos:%8.2f %8.2f %8.2f", + c_name(),contrib.deposit,pos.X(),pos.Y(),pos.Z()); } hit->truth.push_back(contrib); hit->energyDeposit += contrib.deposit; @@ -233,10 +231,10 @@ namespace DD4hep { Hit* hit = coll->find<Hit>(PositionCompare<Hit>(pos)); if ( !hit ) { hit = new Hit(pos); - hit->cellID = volumeID( step ) ; + hit->cellID = volumeID(step); coll->add(hit); if ( 0 == hit->cellID ) { - hit->cellID = volumeID( step ) ; + hit->cellID = volumeID(step); throw runtime_error("Invalid CELL ID for hit!"); } } diff --git a/DDG4/plugins/Geant4SensDet.cpp b/DDG4/plugins/Geant4SensDet.cpp index 66c29fb13..0c5878f38 100644 --- a/DDG4/plugins/Geant4SensDet.cpp +++ b/DDG4/plugins/Geant4SensDet.cpp @@ -54,6 +54,7 @@ namespace DD4hep { { Geant4Kernel& kernel = Geant4Kernel::access(lcdd); setContext(kernel.context()); + m_outputLevel = kernel.getOutputLevel(nam); _aquire(kernel.sensitiveAction(nam)); m_sequence->defineCollections(this); this->G4VSensitiveDetector::SetFilter(this); diff --git a/DDG4/python/DDG4.py b/DDG4/python/DDG4.py index f8e3e3581..e0e47ee77 100644 --- a/DDG4/python/DDG4.py +++ b/DDG4/python/DDG4.py @@ -58,7 +58,7 @@ def _registerGlobalFilter(self,filter): self.get().registerGlobalFilter(Interface.toAction(filter)) #--------------------------------------------------------------------------- def _getKernelProperty(self, name): - print '_getKernelProperty:',str(type(self)),name + #print '_getKernelProperty:',str(type(self)),name ret = Interface.getPropertyKernel(self.get(),name) if ret.status > 0: return ret.data diff --git a/DDG4/src/Geant4Action.cpp b/DDG4/src/Geant4Action.cpp index 3f377a06d..a8837abbf 100644 --- a/DDG4/src/Geant4Action.cpp +++ b/DDG4/src/Geant4Action.cpp @@ -39,6 +39,7 @@ TypeName TypeName::split(const string& type_name) { Geant4Action::Geant4Action(Geant4Context* context, const string& nam) : m_context(context), m_control(0), m_outputLevel(INFO), m_needsControl(false), m_name(nam), m_refCount(1) { InstanceCount::increment(this); + m_outputLevel = context ? context->kernel().getOutputLevel(nam) : (printLevel()-1); declareProperty("Name", m_name); declareProperty("name", m_name); declareProperty("OutputLevel", m_outputLevel); @@ -73,6 +74,13 @@ void Geant4Action::setContext(Geant4Context* context) { m_context = context; } +/// Set the output level; returns previous value +PrintLevel Geant4Action::setOutputLevel(PrintLevel new_level) { + int old = m_outputLevel; + m_outputLevel = new_level; + return (PrintLevel)old; +} + /// Set object properties Geant4Action& Geant4Action::setProperties(PropertyConfigurator& setup) { m_properties.set(m_name, setup); @@ -125,6 +133,26 @@ void Geant4Action::enableUI() { installMessengers(); } +/// Support for messages with variable output level using output level +void Geant4Action::print(const std::string& fmt, ...) const { + if ( outputLevel() >= printLevel() ) { + va_list args; + va_start(args, fmt); + DD4hep::printout((PrintLevel)outputLevel(),m_name, fmt, args); + va_end(args); + } +} + +/// Support for messages with variable output level using output level +void Geant4Action::print(const std::string& tag, const std::string& fmt, ...) const { + if ( outputLevel() >= printLevel() ) { + va_list args; + va_start(args, fmt); + DD4hep::printout((PrintLevel)outputLevel(),tag, fmt, args); + va_end(args); + } +} + /// Support of debug messages. void Geant4Action::debug(const string& fmt, ...) const { va_list args; diff --git a/DDG4/src/Geant4Converter.cpp b/DDG4/src/Geant4Converter.cpp index ae3be70e9..fa3af011d 100644 --- a/DDG4/src/Geant4Converter.cpp +++ b/DDG4/src/Geant4Converter.cpp @@ -249,6 +249,14 @@ namespace { Geant4Converter::Geant4Converter(LCDD& lcdd) : Geant4Mapping(lcdd), m_checkOverlaps(true) { this->Geant4Mapping::init(); + m_outputLevel = PrintLevel(printLevel() - 1); +} + +/// Initializing Constructor +Geant4Converter::Geant4Converter(LCDD& lcdd, PrintLevel level) + : Geant4Mapping(lcdd), m_checkOverlaps(true) { + this->Geant4Mapping::init(); + m_outputLevel = level; } /// Standard destructor @@ -277,7 +285,7 @@ void* Geant4Converter::handleElement(const string& name, const TGeoElement* elem } stringstream str; str << (*g4e); - printout(DEBUG, "Geant4Converter", "++ Created G4 %s", str.str().c_str()); + printout(m_outputLevel, "Geant4Converter", "++ Created G4 %s", str.str().c_str()); } data().g4Elements[element] = g4e; } @@ -331,7 +339,7 @@ void* Geant4Converter::handleMaterial(const string& name, const TGeoMedium* medi } stringstream str; str << (*mat); - printout(DEBUG, "Geant4Converter", "++ Created G4 %s", str.str().c_str()); + printout(m_outputLevel, "Geant4Converter", "++ Created G4 %s", str.str().c_str()); } data().g4Materials[medium] = mat; } @@ -511,7 +519,7 @@ void* Geant4Converter::handleVolume(const string& name, const TGeoVolume* volume } } - printout(DEBUG, "Geant4Converter", "++ Convert Volume %-32s: %p %s/%s assembly:%s sensitive:%s", n.c_str(), v, + printout(m_outputLevel, "Geant4Converter", "++ Convert Volume %-32s: %p %s/%s assembly:%s sensitive:%s", n.c_str(), v, s->IsA()->GetName(), v->IsA()->GetName(), yes_no(assembly), yes_no(det.isValid())); if (assembly) { @@ -526,11 +534,11 @@ void* Geant4Converter::handleVolume(const string& name, const TGeoVolume* volume throw runtime_error("G4Converter: No Geant4 material present for volume:" + n); } if (user_limits) { - printout(DEBUG, "Geant4Converter", "++ Volume + Apply LIMITS settings:%-24s to volume %s.", lim.name(), _v.name()); + printout(m_outputLevel, "Geant4Converter", "++ Volume + Apply LIMITS settings:%-24s to volume %s.", lim.name(), _v.name()); } G4LogicalVolume* vol = new G4LogicalVolume(solid, medium, n, 0, sd, user_limits); if (region) { - printout(DEBUG, "Geant4Converter", "++ Volume + Apply REGION settings: %s to volume %s.", reg.name(), _v.name()); + printout(m_outputLevel, "Geant4Converter", "++ Volume + Apply REGION settings: %s to volume %s.", reg.name(), _v.name()); vol->SetRegion(region); region->AddRootLogicalVolume(vol); } @@ -538,11 +546,11 @@ void* Geant4Converter::handleVolume(const string& name, const TGeoVolume* volume vol->SetVisAttributes(vis_attr); } if (sd) { - printout(DEBUG, "Geant4Converter", "++ Volume: + %s <> %s Solid:%s Mat:%s SD:%s", name.c_str(), vol->GetName().c_str(), + printout(m_outputLevel, "Geant4Converter", "++ Volume: + %s <> %s Solid:%s Mat:%s SD:%s", name.c_str(), vol->GetName().c_str(), solid->GetName().c_str(), medium->GetName().c_str(), sd->GetName().c_str()); } info.g4Volumes[v] = vol; - printout(DEBUG, "Geant4Converter", "++ Volume + %s converted: %p ---> G4: %p", n.c_str(), v, vol); + printout(m_outputLevel, "Geant4Converter", "++ Volume + %s converted: %p ---> G4: %p", n.c_str(), v, vol); } return 0; } @@ -588,7 +596,7 @@ void* Geant4Converter::handleAssembly(const std::string& name, const TGeoNode* n __FILE__, __LINE__, name.c_str(), d->GetName()); } g4->placeAssembly(d,(*assIt).second,transform); - printout(DEBUG, "Geant4Converter", "+++ Assembly: AddPlacedAssembly : dau:%s " + printout(m_outputLevel, "Geant4Converter", "+++ Assembly: AddPlacedAssembly : dau:%s " "to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f", dau_vol->GetName(), mot_vol->GetName(), transform.dx(), transform.dy(), transform.dz()); @@ -600,7 +608,7 @@ void* Geant4Converter::handleAssembly(const std::string& name, const TGeoNode* n __FILE__, __LINE__, name.c_str(), d->GetName()); } g4->placeVolume(d,(*volIt).second, transform); - printout(DEBUG, "Geant4Converter", "+++ Assembly: AddPlacedVolume : dau:%s " + printout(m_outputLevel, "Geant4Converter", "+++ Assembly: AddPlacedVolume : dau:%s " "to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f", dau_vol->GetName(), mot_vol->GetName(), transform.dx(), transform.dy(), transform.dz()); @@ -642,7 +650,7 @@ void* Geant4Converter::handlePlacement(const string& name, const TGeoNode* node) // -- placed volumes were already added before in "handleAssembly" // -- imprint cannot be made, because this requires a logical volume as a mother // - printout(DEBUG, "Geant4Converter", "+++ Assembly: **** : dau:%s " + printout(m_outputLevel, "Geant4Converter", "+++ Assembly: **** : dau:%s " "to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f", vol->GetName(), mot_vol->GetName(), transform.dx(), transform.dy(), transform.dz()); @@ -653,7 +661,7 @@ void* Geant4Converter::handlePlacement(const string& name, const TGeoNode* node) // Node is an assembly: // Imprint the assembly. The mother MUST already be transformed. // - printout(INFO, "Geant4Converter", "+++ Assembly: makeImprint: dau:%s in mother %s " + printout(m_outputLevel, "Geant4Converter", "+++ Assembly: makeImprint: dau:%s in mother %s " "Tr:x=%8.3f y=%8.3f z=%8.3f", node->GetName(), mot_vol->GetName(), transform.dx(), transform.dy(), transform.dz()); @@ -705,7 +713,7 @@ void* Geant4Converter::handleRegion(const TNamed* region, const set<const TGeoVo info->storeSecondaries = r.storeSecondaries(); g4->SetUserInformation(info); - printout(INFO, "Geant4Converter", "++ Converted region settings of:%s.", r.name()); + printout(m_outputLevel, "Geant4Converter", "++ Converted region settings of:%s.", r.name()); vector < string > &limits = r.limits(); for (vector<string>::const_iterator i = limits.begin(); i != limits.end(); ++i) { const string& nam = *i; @@ -849,7 +857,7 @@ void Geant4Converter::handleProperties(LCDD::Properties& prp) const { if (result != 1) { throw runtime_error("Failed to invoke the plugin " + tag + " of type " + type); } - printout(INFO, "Geant4Converter", "+++++ Executed Successfully Geant4 setup module *%s*.", type.c_str()); + printout(m_outputLevel, "Geant4Converter", "+++++ Executed Successfully Geant4 setup module *%s*.", type.c_str()); } } @@ -913,18 +921,18 @@ void* Geant4Converter::printPlacement(const string& name, const TGeoNode* node) stringstream str; str << "G4Cnv::placement: + " << name << " No:" << node->GetNumber() << " Vol:" << vol->GetName() << " Solid:" << sol->GetName(); - printout(DEBUG, "G4Placement", str.str().c_str()); + printout(m_outputLevel, "G4Placement", str.str().c_str()); str.str(""); str << " |" << " Loc: x=" << tr.x() << " y=" << tr.y() << " z=" << tr.z(); - printout(DEBUG, "G4Placement", str.str().c_str()); - printout(DEBUG, "G4Placement", printSolid(sol).c_str()); + printout(m_outputLevel, "G4Placement", str.str().c_str()); + printout(m_outputLevel, "G4Placement", printSolid(sol).c_str()); str.str(""); str << " |" << " Ndau:" << vol->GetNoDaughters() << " physvols." << " Mat:" << vol->GetMaterial()->GetName() << " Mother:" << (char*) (mot ? mot->GetName().c_str() : "---"); - printout(DEBUG, "G4Placement", str.str().c_str()); + printout(m_outputLevel, "G4Placement", str.str().c_str()); str.str(""); str << " |" << " SD:" << sd->GetName(); - printout(DEBUG, "G4Placement", str.str().c_str()); + printout(m_outputLevel, "G4Placement", str.str().c_str()); return g4; } @@ -965,17 +973,17 @@ Geant4Converter& Geant4Converter::create(DetElement top) { handle(this, geo.volumes, &Geant4Converter::collectVolume); handle(this, geo.solids, &Geant4Converter::handleSolid); - printout(INFO, "Geant4Converter", "++ Handled %ld solids.", geo.solids.size()); + printout(m_outputLevel, "Geant4Converter", "++ Handled %ld solids.", geo.solids.size()); handle(this, geo.vis, &Geant4Converter::handleVis); - printout(INFO, "Geant4Converter", "++ Handled %ld visualization attributes.", geo.vis.size()); + printout(m_outputLevel, "Geant4Converter", "++ Handled %ld visualization attributes.", geo.vis.size()); handleMap(this, geo.sensitives, &Geant4Converter::handleSensitive); - printout(INFO, "Geant4Converter", "++ Handled %ld sensitive detectors.", geo.sensitives.size()); + printout(m_outputLevel, "Geant4Converter", "++ Handled %ld sensitive detectors.", geo.sensitives.size()); handleMap(this, geo.limits, &Geant4Converter::handleLimitSet); - printout(INFO, "Geant4Converter", "++ Handled %ld limit sets.", geo.limits.size()); + printout(m_outputLevel, "Geant4Converter", "++ Handled %ld limit sets.", geo.limits.size()); handleMap(this, geo.regions, &Geant4Converter::handleRegion); - printout(INFO, "Geant4Converter", "++ Handled %ld regions.", geo.regions.size()); + printout(m_outputLevel, "Geant4Converter", "++ Handled %ld regions.", geo.regions.size()); handle(this, geo.volumes, &Geant4Converter::handleVolume); - printout(INFO, "Geant4Converter", "++ Handled %ld volumes.", geo.volumes.size()); + printout(m_outputLevel, "Geant4Converter", "++ Handled %ld volumes.", geo.volumes.size()); handleRMap(this, *m_data, &Geant4Converter::handleAssembly); // Now place all this stuff appropriately handleRMap(this, *m_data, &Geant4Converter::handlePlacement); @@ -987,5 +995,6 @@ Geant4Converter& Geant4Converter::create(DetElement top) { geo.setWorld(top.placement().ptr()); geo.valid = true; + printout(INFO, "Geant4Converter", "+++ Successfully converted geometry to Geant4."); return *this; } diff --git a/DDG4/src/Geant4DetectorConstruction.cpp b/DDG4/src/Geant4DetectorConstruction.cpp index 8be4b853e..4a7ae76ea 100644 --- a/DDG4/src/Geant4DetectorConstruction.cpp +++ b/DDG4/src/Geant4DetectorConstruction.cpp @@ -1,6 +1,15 @@ +// $Id$ +//==================================================================== +// AIDA Detector description implementation +//-------------------------------------------------------------------- +// +// Author : M.Frank +// +//==================================================================== #include "DDG4/Geant4DetectorConstruction.h" #include "DDG4/Geant4HierarchyDump.h" #include "DDG4/Geant4Converter.h" +#include "DDG4/Geant4Kernel.h" #include "DD4hep/LCDD.h" #include "TGeoManager.h" #include "G4PVPlacement.hh" @@ -15,14 +24,23 @@ using namespace DD4hep; using namespace DD4hep::Geometry; DD4hep::Simulation::Geant4DetectorConstruction::Geant4DetectorConstruction(Geometry::LCDD& lcdd) - : m_lcdd(lcdd), m_world(0) { + : m_outputLevel(PrintLevel(printLevel()-1)), m_lcdd(lcdd), m_world(0) { +} + +DD4hep::Simulation::Geant4DetectorConstruction::Geant4DetectorConstruction(DD4hep::Simulation::Geant4Kernel& kernel) + : m_lcdd(kernel.lcdd()), m_world(0) { + m_outputLevel = kernel.getOutputLevel("Geant4Converter"); +} + +/// Default destructor +DD4hep::Simulation::Geant4DetectorConstruction::~Geant4DetectorConstruction() { } G4VPhysicalVolume* DD4hep::Simulation::Geant4DetectorConstruction::Construct() { typedef Simulation::Geant4Converter Geant4Converter; Geant4Mapping& g4map = Geant4Mapping::instance(); DetElement world = m_lcdd.world(); - Geant4Converter conv(m_lcdd); + Geant4Converter conv(m_lcdd, m_outputLevel); Geant4GeometryInfo* info = conv.create(world).detach(); g4map.attach(info); m_world = info->world(); diff --git a/DDG4/src/Geant4Exec.cpp b/DDG4/src/Geant4Exec.cpp index 20a01b3c6..ecae73f3c 100644 --- a/DDG4/src/Geant4Exec.cpp +++ b/DDG4/src/Geant4Exec.cpp @@ -235,11 +235,12 @@ int Geant4Exec::configure(Geant4Kernel& kernel) { G4RunManager& runManager = kernel.runManager(); // Check if the geometry was loaded - if (lcdd.detectors().size() <= 2) { - cout << "Error, no geometry." << endl; + if (lcdd.detectors().size() <= 1) { + printout(INFO, "Geant4Exec", "+++ Only %d subdetectors present. " + "You sure you loaded the geometry properly?",int(lcdd.detectors().size())); } // Get the detector constructed - Geant4DetectorConstruction* detector = new Geant4DetectorConstruction(lcdd); + Geant4DetectorConstruction* detector = new Geant4DetectorConstruction(kernel); runManager.SetUserInitialization(detector); G4VUserPhysicsList* physics = 0; diff --git a/DDG4/src/Geant4Kernel.cpp b/DDG4/src/Geant4Kernel.cpp index 7febc78dc..5cf0bb546 100644 --- a/DDG4/src/Geant4Kernel.cpp +++ b/DDG4/src/Geant4Kernel.cpp @@ -24,7 +24,6 @@ #include "DDG4/Geant4StackingAction.h" #include "DDG4/Geant4GeneratorAction.h" #include "DDG4/Geant4SensDetAction.h" -#include "DDG4/Geant4DetectorConstruction.h" #include "DDG4/Geant4MonteCarloRecordManager.h" #include "DDG4/Geant4TrackPersistency.h" @@ -84,7 +83,9 @@ Geant4Kernel::Geant4Kernel(LCDD& lcdd) m_context = new Geant4Context(this); m_lcdd.addExtension < Geant4Kernel > (this); declareProperty("UI",m_uiName); + declareProperty("OutputLevel",m_outputLevel = DEBUG); declareProperty("NumEvents",m_numEvent = 10); + declareProperty("OutputLevels",m_clientLevels); m_controlName = "/ddg4/"; m_control = new G4UIdirectory(m_controlName.c_str()); m_control->SetGuidance("Control for all named Geant4 actions"); @@ -121,6 +122,15 @@ Geant4Kernel& Geant4Kernel::instance(LCDD& lcdd) { return obj; } +void Geant4Kernel::printProperties() const { + printout(ALWAYS,"Geant4Kernel","OutputLevel: %d",m_outputLevel); + printout(ALWAYS,"Geant4Kernel","UI: %s",m_uiName.c_str()); + printout(ALWAYS,"Geant4Kernel","NumEvents: %ld",m_numEvent); + for(ClientOutputLevels::const_iterator i=m_clientLevels.begin(); i!=m_clientLevels.end();++i) { + printout(ALWAYS,"Geant4Kernel","OutputLevel[%s]: %d",(*i).first.c_str(),(*i).second); + } +} + /// Check property for existence bool Geant4Kernel::hasProperty(const std::string& name) const { return m_properties.exists(name); @@ -141,6 +151,25 @@ Geant4Kernel& Geant4Kernel::access(LCDD& lcdd) { return *kernel; } +/// Fill cache with the global output level of a named object. Must be set before instantiation +void Geant4Kernel::setOutputLevel(const std::string object, PrintLevel new_level) { + m_clientLevels[object] = new_level; +} + +/// Retrieve the global output level of a named object. +DD4hep::PrintLevel Geant4Kernel::getOutputLevel(const std::string object) const { + ClientOutputLevels::const_iterator i=m_clientLevels.find(object); + if ( i != m_clientLevels.end() ) return (PrintLevel)(*i).second; + return DD4hep::PrintLevel(DD4hep::printLevel()-1); +} + +/// Set the output level; returns previous value +DD4hep::PrintLevel Geant4Kernel::setOutputLevel(PrintLevel new_level) { + int old = m_outputLevel; + m_outputLevel = new_level; + return (PrintLevel)old; +} + /// Access to the Geant4 run manager G4RunManager& Geant4Kernel::runManager() { if (m_runManager) diff --git a/DDG4/src/Geant4OutputAction.cpp b/DDG4/src/Geant4OutputAction.cpp index d0c61b839..fbba4ea29 100644 --- a/DDG4/src/Geant4OutputAction.cpp +++ b/DDG4/src/Geant4OutputAction.cpp @@ -8,6 +8,7 @@ //==================================================================== // Framework include files +#include "DD4hep/Printout.h" #include "DD4hep/InstanceCount.h" #include "DDG4/Geant4RunAction.h" #include "DDG4/Geant4OutputAction.h" @@ -42,13 +43,18 @@ void Geant4OutputAction::begin(const G4Event* /* event */) { void Geant4OutputAction::end(const G4Event* evt) { OutputContext < G4Event > ctxt(evt); G4HCofThisEvent* hce = evt->GetHCofThisEvent(); - int nCol = hce->GetNumberOfCollections(); - saveEvent(ctxt); - for (int i = 0; i < nCol; ++i) { - G4VHitsCollection* hc = hce->GetHC(i); - saveCollection(ctxt, hc); + if ( hce ) { + int nCol = hce->GetNumberOfCollections(); + saveEvent(ctxt); + for (int i = 0; i < nCol; ++i) { + G4VHitsCollection* hc = hce->GetHC(i); + saveCollection(ctxt, hc); + } + commit(ctxt); + return; } - commit(ctxt); + printout(WARNING,"Geant4OutputAction", + "+++ The value of G4HCofThisEvent is NULL. No collections saved!"); } /// Commit data at end of filling procedure diff --git a/DDG4/src/Geant4ParticleGun.cpp b/DDG4/src/Geant4ParticleGun.cpp index da647f25c..bfefea6fa 100644 --- a/DDG4/src/Geant4ParticleGun.cpp +++ b/DDG4/src/Geant4ParticleGun.cpp @@ -8,13 +8,17 @@ //==================================================================== // Framework include files +#include "DD4hep/Printout.h" #include "DD4hep/InstanceCount.h" #include "DDG4/Geant4ParticleGun.h" +#include "CLHEP/Units/SystemOfUnits.h" #include "G4ParticleGun.hh" #include "G4ParticleTable.hh" #include "G4ParticleDefinition.hh" +#include "TRandom1.h" + // C/C++ include files #include <stdexcept> #include <limits> @@ -26,7 +30,7 @@ using namespace DD4hep::Simulation; /// Standard constructor Geant4ParticleGun::Geant4ParticleGun(Geant4Context* context, const string& name) : Geant4GeneratorAction(context, name), m_position(0,0,0), m_direction(1,1,0.3), - m_particle(0), m_gun(0) + m_particle(0), m_gun(0), m_rndm(0), m_shotNo(0) { InstanceCount::increment(this); m_needsControl = true; @@ -35,12 +39,15 @@ Geant4ParticleGun::Geant4ParticleGun(Geant4Context* context, const string& name) declareProperty("multiplicity", m_multiplicity = 1); declareProperty("position", m_position); declareProperty("direction", m_direction); + declareProperty("isotrop", m_isotrop = false); } /// Default destructor Geant4ParticleGun::~Geant4ParticleGun() { if (m_gun) delete m_gun; + if ( m_rndm ) + delete m_rndm; InstanceCount::decrement(this); } @@ -53,16 +60,31 @@ void Geant4ParticleGun::operator()(G4Event* event) { G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable(); m_particle = particleTable->FindParticle(m_particleName); if (0 == m_particle) { - throw runtime_error("Bad particle type!"); + throw runtime_error("Bad particle type:"+m_particleName+"!"); } } - double r = m_direction.R(), eps = numeric_limits<float>::epsilon(); - if ( r > eps && std::fabs(r-1.0) > eps ) { - m_direction.SetXYZ(m_direction.X()/r, m_direction.Y()/r, m_direction.Z()/r); + if ( m_isotrop ) { + if ( 0 == m_rndm ) m_rndm = new TRandom1(); + double phi = 2*M_PI*m_rndm->Rndm(); + double theta = M_PI*m_rndm->Rndm(); + double x1 = std::sin(theta)*std::cos(phi); + double x2 = std::sin(theta)*std::sin(phi); + double x3 = std::cos(theta); + m_direction.SetXYZ(x1,x2,x3); + } + else { + double r = m_direction.R(), eps = numeric_limits<float>::epsilon(); + if ( r > eps ) { + m_direction.SetXYZ(m_direction.X()/r, m_direction.Y()/r, m_direction.Z()/r); + } } + print("Geant4ParticleGun","Shoot [%d] %.3f GeV %s pos:(%.3f %.3f %.3f)[mm] dir:(%6.3f %6.3f %6.3f)", + m_shotNo, m_energy/GeV, m_particleName.c_str(), m_position.X(), m_position.Y(), m_position.Z(), + m_direction.X(),m_direction.Y(), m_direction.Z()); m_gun->SetParticleDefinition(m_particle); m_gun->SetParticleEnergy(m_energy); m_gun->SetParticleMomentumDirection(G4ThreeVector(m_direction.X(), m_direction.Y(), m_direction.Z())); m_gun->SetParticlePosition(G4ThreeVector(m_position.X(), m_position.Y(), m_position.Z())); m_gun->GeneratePrimaryVertex(event); + ++m_shotNo; } diff --git a/DDG4/src/Geant4TestActions.cpp b/DDG4/src/Geant4TestActions.cpp index 4db4b11c7..c1853e4b5 100644 --- a/DDG4/src/Geant4TestActions.cpp +++ b/DDG4/src/Geant4TestActions.cpp @@ -54,21 +54,21 @@ Geant4TestRunAction::~Geant4TestRunAction() { } /// begin-of-run callback void Geant4TestRunAction::begin(const G4Run* run) { - printout(INFO, name(), "%s> calling begin(run_id=%d,num_event=%d)", m_type.c_str(), run->GetRunID(), - run->GetNumberOfEventToBeProcessed()); + print(name(), "%s> calling begin(run_id=%d,num_event=%d)", m_type.c_str(), run->GetRunID(), + run->GetNumberOfEventToBeProcessed()); } /// End-of-run callback void Geant4TestRunAction::end(const G4Run* run) { - printout(INFO, name(), "%s> calling end(run_id=%d, num_event=%d)", - m_type.c_str(), run->GetRunID(), run->GetNumberOfEvent()); + print(name(), "%s> calling end(run_id=%d, num_event=%d)", + m_type.c_str(), run->GetRunID(), run->GetNumberOfEvent()); } /// begin-of-event callback void Geant4TestRunAction::beginEvent(const G4Event* evt) { - printout(INFO, name(), "%s> calling beginEvent(event_id=%d)", m_type.c_str(), evt->GetEventID()); + print(name(), "%s> calling beginEvent(event_id=%d)", m_type.c_str(), evt->GetEventID()); } /// End-of-event callback void Geant4TestRunAction::endEvent(const G4Event* evt) { - printout(INFO, name(), "%s> calling endEvent(event_id=%d)", m_type.c_str(), evt->GetEventID()); + print(name(), "%s> calling endEvent(event_id=%d)", m_type.c_str(), evt->GetEventID()); } Geant4TestEventAction::Geant4TestEventAction(Geant4Context* c, const std::string& n) @@ -80,24 +80,24 @@ Geant4TestEventAction::~Geant4TestEventAction() { } /// begin-of-event callback void Geant4TestEventAction::begin(const G4Event* evt) { - printout(INFO, name(), "%s> calling begin(event_id=%d)", m_type.c_str(), evt->GetEventID()); + print(name(), "%s> calling begin(event_id=%d)", m_type.c_str(), evt->GetEventID()); } /// End-of-event callback void Geant4TestEventAction::end(const G4Event* evt) { - printout(INFO, name(), "%s> calling end(event_id=%d)", m_type.c_str(), evt->GetEventID()); + print(name(), "%s> calling end(event_id=%d)", m_type.c_str(), evt->GetEventID()); } /// begin-of-run callback void Geant4TestEventAction::beginRun(const G4Run* run) { - printout(INFO, name(), "%s> calling beginRun(run_id=%d,num_event=%d)", + print(name(), "%s> calling beginRun(run_id=%d,num_event=%d)", m_type.c_str(), run->GetRunID(), run->GetNumberOfEventToBeProcessed()); } /// End-of-run callback void Geant4TestEventAction::endRun(const G4Run* run) { - printout(INFO, name(), "%s> calling endRun(run_id=%d, num_event=%d)", + print(name(), "%s> calling endRun(run_id=%d, num_event=%d)", m_type.c_str(), run->GetRunID(), run->GetNumberOfEvent()); } @@ -111,14 +111,14 @@ Geant4TestTrackAction::~Geant4TestTrackAction() { } /// Begin-of-tracking callback void Geant4TestTrackAction::begin(const G4Track* trk) { - printout(INFO, name(), "%s> calling begin(track=%d, parent=%d, position=(%f,%f,%f))", + print(name(), "%s> calling begin(track=%d, parent=%d, position=(%f,%f,%f))", m_type.c_str(), trk->GetTrackID(), trk->GetParentID(), trk->GetPosition().x(), trk->GetPosition().y(), trk->GetPosition().z()); } /// End-of-tracking callback void Geant4TestTrackAction::end(const G4Track* trk) { - printout(INFO, name(), "%s> calling end(track=%d, parent=%d, position=(%f,%f,%f))", + print(name(), "%s> calling end(track=%d, parent=%d, position=(%f,%f,%f))", m_type.c_str(), trk->GetTrackID(), trk->GetParentID(), trk->GetPosition().x(), trk->GetPosition().y(), trk->GetPosition().z()); } @@ -132,14 +132,14 @@ Geant4TestStepAction::~Geant4TestStepAction() { } /// User stepping callback void Geant4TestStepAction::operator()(const G4Step*, G4SteppingManager*) { - printout(INFO, name(), "%s> calling operator()", m_type.c_str()); + print(name(), "%s> calling operator()", m_type.c_str()); } Geant4TestSensitive::Geant4TestSensitive(Geant4Context* c, const std::string& n, DetElement det, LCDD& lcdd) : Geant4Sensitive(c, n, det, lcdd), Geant4TestBase(this, "Geant4TestSensitive") { InstanceCount::increment(this); m_collectionID = defineCollection < TestHit > (n); - printout(INFO, name(), "%s> Collection ID is %d", m_type.c_str(), int(m_collectionID)); + print(name(), "%s> Collection ID is %d", m_type.c_str(), int(m_collectionID)); } Geant4TestSensitive::~Geant4TestSensitive() { InstanceCount::decrement(this); @@ -148,7 +148,7 @@ Geant4TestSensitive::~Geant4TestSensitive() { /// Begin-of-tracking callback void Geant4TestSensitive::begin(G4HCofThisEvent* hce) { Geant4HitCollection* c = collectionByID(m_collectionID); - printout(INFO, name(), "%s> calling begin(num_coll=%d, coll=%s)", + print(name(), "%s> calling begin(num_coll=%d, coll=%s)", m_type.c_str(), hce->GetNumberOfCollections(), c ? c->GetName().c_str() : "None"); } @@ -156,16 +156,16 @@ void Geant4TestSensitive::begin(G4HCofThisEvent* hce) { /// End-of-tracking callback void Geant4TestSensitive::end(G4HCofThisEvent* hce) { Geant4HitCollection* c = collection(m_collectionID); - printout(INFO, name(), "%s> calling end(num_coll=%d, coll=%s)", - m_type.c_str(), hce->GetNumberOfCollections(), - c ? c->GetName().c_str() : "None"); + print(name(), "%s> calling end(num_coll=%d, coll=%s)", + m_type.c_str(), hce->GetNumberOfCollections(), + c ? c->GetName().c_str() : "None"); } /// Method for generating hit(s) using the information of G4Step object. bool Geant4TestSensitive::process(G4Step* step, G4TouchableHistory*) { Geant4HitCollection* c = collection(m_collectionID); - printout(INFO, name(), "%s> calling process(track=%d, dE=%f, dT=%f len=%f, First,last in Vol=(%c,%c), coll=%s)", - m_type.c_str(), step->GetTrack()->GetTrackID(), + print(name(), "%s> calling process(track=%d, dE=%f, dT=%f len=%f, First,last in Vol=(%c,%c), coll=%s)", + m_type.c_str(), step->GetTrack()->GetTrackID(), step->GetTotalEnergyDeposit(), step->GetDeltaTime(), step->GetStepLength(), step->IsFirstStepInVolume() ? 'Y' : 'N', step->IsLastStepInVolume() ? 'Y' : 'N', diff --git a/DDG4/src/Geant4VolumeManager.cpp b/DDG4/src/Geant4VolumeManager.cpp index 03911c434..19536172c 100644 --- a/DDG4/src/Geant4VolumeManager.cpp +++ b/DDG4/src/Geant4VolumeManager.cpp @@ -114,7 +114,7 @@ namespace { PrintLevel print_res = DEBUG; printout(print_action,"Geant4VolumeManager","+++ Add path:%s vid:%016X", - DetectorTools::placementPath(nodes,true).c_str(),code); + DetectorTools::placementPath(nodes,false).c_str(),code); if (i == m_entries.end()) { path.reserve(nodes.size()); @@ -165,16 +165,17 @@ namespace { int(control.size()),DetectorTools::placementPath(control,true).c_str()); goto Err; } - printout(ERROR, "Geant4VolumeManager", "populate: Severe error: Duplicated Volume entry: %X %s", - code, " [THIS SHOULD NEVER HAPPEN]"); + printout(ERROR, "Geant4VolumeManager", "populate: Severe error: Duplicated Volume entry: %X" + " [THIS SHOULD NEVER HAPPEN]", code); Err: if ( i != m_entries.end() ) - printout(ERROR, "Geant4VolumeManager", " Geant4 path: %s",placementPath((*i).second).c_str()); - else if ( !path.empty() ) - printout(ERROR, "Geant4VolumeManager", " Geant4 path: %s",placementPath(path).c_str()); - printout(ERROR, "Geant4VolumeManager", " Offend.VolIDs: %s", - DetectorTools::toString(ids).c_str()); + printout(ERROR,"Geant4VolumeManager"," Known G4 path: %s",placementPath((*i).second).c_str()); + if ( !path.empty() ) + printout(ERROR,"Geant4VolumeManager"," New G4 path: %s",placementPath(path).c_str()); + if ( !nodes.empty() ) + printout(ERROR,"Geant4VolumeManager"," TGeo path: %s",DetectorTools::placementPath(nodes,false).c_str()); + printout(ERROR,"Geant4VolumeManager", " Offend.VolIDs: %s",DetectorTools::toString(ro.idSpec(),ids,code).c_str()); throw runtime_error("Failed to populate Geant4 volume manager!"); } }; diff --git a/doc/release.notes b/doc/release.notes index a0a7742e9..06ac1d091 100644 --- a/doc/release.notes +++ b/doc/release.notes @@ -1,6 +1,16 @@ DD4hep ---- Release Notes ================================= +2014/05/28 Markus Frank +----------------------- + - Several fixes and improvements in DDG4 + - Cache output levels in Geant4Kernel object + - Add shell handler to measure energy escape e.g. in calorimeters + - Fix bug in Polycone shapes concerning angles. + IMPORTANT: If you forked the repository examples/CLICSiD + your PolyconeSupport does not work anymore, because the angles + were given in degrees. They should be given in radians! + 2014/05/21 Markus Frank ----------------------- - Fix material creation from XML diff --git a/examples/CLICSiD/compact/compact.xml b/examples/CLICSiD/compact/compact.xml index c69ffdd24..b1485ee2e 100644 --- a/examples/CLICSiD/compact/compact.xml +++ b/examples/CLICSiD/compact/compact.xml @@ -1528,7 +1528,7 @@ </readout> <readout name="LumiCalHits"> <segmentation type="CartesianGridXY" grid_size_x="0.35*cm" grid_size_y="0.35*cm" /> - <id>system:8,layer:8,barrel:3,layer:8,slice:5,x:32:-16,y:-16</id> + <id>system:8,barrel:3,layer:8,slice:8,x:32:-16,y:-16</id> </readout> <readout name="BeamCalHits"> <segmentation type="CartesianGridXY" grid_size_x="0.35*cm" grid_size_y="0.35*cm" /> diff --git a/examples/CLICSiD/src/CylindricalEndcapCalorimeter_geo.cpp b/examples/CLICSiD/src/CylindricalEndcapCalorimeter_geo.cpp index 3aebcc9b7..3dd8d4e2d 100644 --- a/examples/CLICSiD/src/CylindricalEndcapCalorimeter_geo.cpp +++ b/examples/CLICSiD/src/CylindricalEndcapCalorimeter_geo.cpp @@ -35,7 +35,7 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { double layerWidth = 0; for(xml_coll_t l(x_layer,_U(slice)); l; ++l) layerWidth += xml_comp_t(l).thickness(); - for(int i=0, m=0, repeat=x_layer.repeat(); i<repeat; ++i, m=0) { + for(int i=0, m=0, repeat=x_layer.repeat(); i<repeat; ++i) { double zlayer = z; string layer_name = det_name + _toString(layer_num,"_layer%d"); Volume layer_vol(layer_name,Tube(rmin,rmax,layerWidth),air); @@ -69,15 +69,23 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { envelopeVol.setAttributes(lcdd,x_det.regionStr(),x_det.limitsStr(),x_det.visStr()); DetElement sdet(det_name,x_det.id()); + Assembly assembly(det_name+"_assembly"); Volume motherVol = lcdd.pickMotherVolume(sdet); - PlacedVolume phv = motherVol.placeVolume(envelopeVol,Position(0,0,zmin+totWidth/2)); - phv.addPhysVolID("system",sdet.id()) - .addPhysVolID("barrel",1); + PlacedVolume phv = motherVol.placeVolume(assembly); + phv.addPhysVolID("system",sdet.id()); sdet.setPlacement(phv); + + DetElement sdetA(sdet,det_name+(reflect ? "_A" : ""),x_det.id()); + phv = assembly.placeVolume(envelopeVol,Position(0,0,zmin+totWidth/2)); + phv.addPhysVolID("barrel",1); + sdetA.setPlacement(phv); + if ( reflect ) { - phv=motherVol.placeVolume(envelopeVol,Transform3D(RotationZ(M_PI),Position(0,0,-zmin-totWidth/2))); - phv.addPhysVolID("system",sdet.id()) - .addPhysVolID("barrel",2); + phv=assembly.placeVolume(envelopeVol,Transform3D(RotationZ(M_PI),Position(0,0,-zmin-totWidth/2))); + phv.addPhysVolID("barrel",2); + /// Create the detector element for the opposite side.... + DetElement sdetB(sdet,det_name+"_B",x_det.id()); + sdetB.setPlacement(phv); } return sdet; } diff --git a/examples/CLICSiD/src/ForwardDetector_geo.cpp b/examples/CLICSiD/src/ForwardDetector_geo.cpp index a76510afb..28394668c 100644 --- a/examples/CLICSiD/src/ForwardDetector_geo.cpp +++ b/examples/CLICSiD/src/ForwardDetector_geo.cpp @@ -159,17 +159,29 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { } sdet.setVisAttributes(lcdd, x_det.visStr(), envelopeVol); - PlacedVolume env_phv = motherVol.placeVolume(envelopeVol,Transform3D(RotationZ(M_PI),Position(0,0,zpos))); - env_phv.addPhysVolID("system", id); - env_phv.addPhysVolID("barrel", 1); - sdet.setPlacement(env_phv); // Reflect it. if ( reflect ) { - env_phv = motherVol.placeVolume(envelopeVol,Transform3D(RotationY(M_PI),Position(0,0,-zpos))); - env_phv.addPhysVolID("system", id); - env_phv.addPhysVolID("barrel", 2); - DetElement rdet(det_name+"_reflect",x_det.id()); - rdet.setPlacement(env_phv); + Assembly assembly(det_name+"_assembly"); + PlacedVolume pv = motherVol.placeVolume(assembly); + pv.addPhysVolID("system", id); + sdet.setPlacement(pv); + + DetElement sdetA(sdet,det_name+"_A",x_det.id()); + pv = assembly.placeVolume(envelopeVol,Transform3D(RotationZ(M_PI),Position(0,0,zpos))); + pv.addPhysVolID("barrel", 1); + sdetA.setPlacement(pv); + + DetElement sdetB = sdetA.clone(det_name+"_B",x_det.id()); + sdet.add(sdetB); + pv = assembly.placeVolume(envelopeVol,Transform3D(RotationY(M_PI),Position(0,0,-zpos))); + pv.addPhysVolID("barrel", 2); + sdetB.setPlacement(pv); + } + else { + PlacedVolume pv = motherVol.placeVolume(envelopeVol,Transform3D(RotationZ(M_PI),Position(0,0,zpos))); + pv.addPhysVolID("system", id); + pv.addPhysVolID("barrel", 1); + sdet.setPlacement(pv); } return sdet; } diff --git a/examples/CLICSiD/src/PolyconeSupport_geo.cpp b/examples/CLICSiD/src/PolyconeSupport_geo.cpp index f2634f0b7..b1aa834a6 100644 --- a/examples/CLICSiD/src/PolyconeSupport_geo.cpp +++ b/examples/CLICSiD/src/PolyconeSupport_geo.cpp @@ -13,10 +13,10 @@ using namespace DD4hep; using namespace DD4hep::Geometry; static Ref_t create_detector(LCDD& lcdd, xml_h e, Ref_t) { - xml_det_t x_det = e; - string name = x_det.nameStr(); - DetElement sdet (name,x_det.id()); - Material mat (lcdd.material(x_det.materialStr())); + xml_det_t x_det = e; + string name = x_det.nameStr(); + DetElement sdet (name,x_det.id()); + Material mat (lcdd.material(x_det.materialStr())); vector<double> rmin,rmax,z; int num = 0; @@ -29,10 +29,15 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, Ref_t) { if ( num < 2 ) { throw runtime_error("PolyCone["+name+"]> Not enough Z planes. minimum is 2!"); } - Polycone cone (0.,2.*M_PI*RAD_2_DEGREE,rmin,rmax,z); + Polycone cone (0,2*M_PI,rmin,rmax,z); Volume volume(name, cone, mat); volume.setVisAttributes(lcdd, x_det.visStr()); - sdet.setPlacement(lcdd.pickMotherVolume(sdet).placeVolume(volume)); + PlacedVolume pv = lcdd.pickMotherVolume(sdet).placeVolume(volume); + sdet.setPlacement(pv); + if ( x_det.hasAttr(_U(id)) ) { + int det_id = x_det.id(); + pv.addPhysVolID("system",det_id); + } return sdet; } diff --git a/examples/ClientTests/compact/FCC_HcalBarrel.xml b/examples/ClientTests/compact/FCC_HcalBarrel.xml index 53fa198fe..5d52bb835 100644 --- a/examples/ClientTests/compact/FCC_HcalBarrel.xml +++ b/examples/ClientTests/compact/FCC_HcalBarrel.xml @@ -21,6 +21,7 @@ <display> <vis name="Invisible" showDaughters="false" visible="false"/> <vis name="InvisibleWithChildren" showDaughters="true" visible="false"/> + <vis name="VisibleRed" r="1.0" g="0.0" b="0.0" showDaughters="true" visible="true"/> <vis name="VisibleBlue" r="0.0" g="0.0" b="1.0" showDaughters="false" visible="true"/> <vis name="VisibleGreen" alpha="1.0" r="0.0" g="1.0" b="0.0" drawingStyle="solid" lineStyle="solid" showDaughters="true" visible="true"/> <vis name="HcalStrip" r="0.0" g="0.0" b="1.0" showDaughters="true" visible="false"/> @@ -35,7 +36,7 @@ </limits> <detectors> - <detector id="1" name="HcalBarrel" type="FCC_HcalBarrel" readout="HcalBarrelHits" vis="InvisibleWithChildren"> + <detector id="1" name="HcalBarrel" type="FCC_HcalBarrel2" readout="HcalBarrelHits" vis="VisibleGreen" limits="cal_limits"> <comment>A barrel hadronic calorimeter inspired on the ATLAS Tile hadronic calorimeter</comment> <dimensions rmin="HcalBarrel_rmin" rmax="HcalBarrel_rmax" z="HcalBarrel_zmax" phiBins="64" /> <layer dr="10.00*cm" dz="3.00*cm" vis="HcalStrip"> @@ -51,14 +52,33 @@ <slice material="Polystyrene" dz="0.50*cm" sensitive="yes" limits="cal_limits" vis="HcalActive"/> </layer> </detector> - </detectors> + + <detector id="2" name="ContainmentShell" type="ZylinderShell" vis="VisibleRed" readout="ContainmentHits" > + <comment>Containment shell to measure calorimeter escapes</comment> + <material name="Air"/> + <module name="Barrel" id="0" vis="VisibleRed"> + <zplane rmin="HcalBarrel_rmax+20*cm" rmax="HcalBarrel_rmax+22*cm" z="-2*HcalBarrel_zmax"/> + <zplane rmin="HcalBarrel_rmax+20*cm" rmax="HcalBarrel_rmax+22*cm" z="2*HcalBarrel_zmax"/> + </module> + <module name="SideA" id="1" vis="VisibleRed"> + <zplane rmin="0" rmax="HcalBarrel_rmax+22*cm" z="2*HcalBarrel_zmax+10*cm"/> + <zplane rmin="0" rmax="HcalBarrel_rmax+22*cm" z="2*HcalBarrel_zmax+20*cm"/> + </module> + <module name="SideB" id="2" vis="VisibleRed"> + <zplane rmin="0" rmax="HcalBarrel_rmax+22*cm" z="-(2*HcalBarrel_zmax+10*cm)"/> + <zplane rmin="0" rmax="HcalBarrel_rmax+22*cm" z="-(2*HcalBarrel_zmax+20*cm)"/> + </module> + </detector> + </detectors> <readouts> <readout name="HcalBarrelHits"> - <segmentation type="CartesianGridXY" grid_size_x="3.0*cm" grid_size_y="3.0*cm" /> - <id>system:8,barrel:3,stave:8,module:6,layer:8,slice:5,x:46:-8,y:-8</id> + <id>system:8,barrel:3,module:7,layer:5,slice:10</id> + </readout> + <readout name="ContainmentHits"> + <id>system:8,barrel:3</id> </readout> - </readouts> + </readouts> <fields> <field name="GlobalSolenoid" type="solenoid" diff --git a/examples/ClientTests/src/FCC_HcalBarrel2_geo.cpp b/examples/ClientTests/src/FCC_HcalBarrel2_geo.cpp new file mode 100644 index 000000000..7b2436734 --- /dev/null +++ b/examples/ClientTests/src/FCC_HcalBarrel2_geo.cpp @@ -0,0 +1,184 @@ +/********************************* + * HcalBarrel_geo + * Implementing a detector + * + * Carlos.Solans@cern.ch + *********************************/ + +#include "DD4hep/DetFactoryHelper.h" +#include "XML/Layering.h" + +#include <vector> + +using namespace std; +using namespace DD4hep; +using namespace DD4hep::Geometry; + +static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { + + //XML detector object: DDCore/XML/XMLDetector.h + DD4hep::XML::DetElement x_det = e; + + //Create the DetElement for DD4hep + DetElement d_det(x_det.nameStr(),x_det.id()); + + //Pick the mothervolume + Volume det_vol = lcdd.pickMotherVolume(d_det); + + //XML dimension object: DDCore/XML/XMLDimension.h + DD4hep::XML::Dimension x_det_dim(x_det.dimensions()); + + //Tube: DDCore/DD4hep/Shapes.h + Tube calo_shape(x_det_dim.rmin(),x_det_dim.rmax(),x_det_dim.z()); + + //Create the detector mother volume + Volume calo_vol(x_det.nameStr()+"_envelope",calo_shape,lcdd.air()); + + //Set envelope volume attributes + calo_vol.setAttributes(lcdd,x_det.regionStr(),x_det.limitsStr(),x_det.visStr()); + + //Place inside the mother volume + PlacedVolume calo_plv = det_vol.placeVolume(calo_vol); + + calo_plv.addPhysVolID("system",x_det.id()); + calo_plv.addPhysVolID("barrel",0); + d_det.setPlacement(calo_plv); + + //Declare this sensitive detector as a calorimeter + sens.setType("calorimeter"); + + + int layer_num = 0; + float layer_pos_z = 0; + double tile_phi = 2*M_PI/x_det_dim.phiBins(); + float r = x_det_dim.rmin(); + + bool debug = false; + + std::map<std::string,int> vol_max; + + //Repeat layers until we reach the rmax + while(r<x_det_dim.rmax()){ + + //Loop over layers of type: XML Collection_t object: DDCore/XML/XMLElements.h + for(DD4hep::XML::Collection_t layerIt(x_det,_U(layer));layerIt; ++layerIt){ + + //Build a layer volume + DD4hep::XML::Component x_det_layer = layerIt; + + float dr = x_det_layer.dr(); + + string layer_name = x_det.nameStr()+_toString(layer_num,"_layer%d"); + + float x1 = r * tan(tile_phi/2.); + float x2 = (r+dr) * tan(tile_phi/2.); + float y1 = x_det_dim.z(); + float y2 = x_det_dim.z(); + float z = x_det_layer.dr(); + + if(debug){ + cout << " r:" << r + << " dr:" << dr + << " x1:" << x1 + << " x2:" << x2 + << " y1:" << y1 + << " y2:" << y2 + << " z:" << z + << endl; + } + + //Shape a Trapezoid (tile): DDCore/DD4hep/Shapes.h + Trapezoid layer_shape(x1,x2,y1,y2,z); + + //Create a volume with trapezoid shape + Volume layer_vol(layer_name, layer_shape, lcdd.air()); + layer_vol.setAttributes(lcdd,x_det.regionStr(),x_det.limitsStr(),x_det_layer.visStr()); + + //DetElement layer(layer_name,_toString(layer_num,"layer%d"),x_det.id()); + + //Fill the volume with tiles + + + int tile_number = 0; + vector<Volume> tiles; + + //Repeat slices until we reach the end of the calorimeter + for(xml_coll_t k(x_det_layer,_U(slice)); k; ++k) { + + DD4hep::XML::Component tile_xml = k; + string tile_name = layer_name + _toString(tile_number,"_slice%d"); + Material tile_material = lcdd.material(tile_xml.materialStr()); + float tile_thickness = tile_xml.dz(); + float tile_y1 = tile_thickness; + float tile_y2 = tile_thickness; + float tile_z = x_det_layer.dr(); + + //Shape a Trapezoid (tile): DDCore/DD4hep/Shapes.h + Trapezoid tile_shape(x1,x2,tile_y1,tile_y2,tile_z); + + //Create a volume with trapezoid shape + Volume tile_vol(tile_name,tile_shape,tile_material); + + if ( tile_xml.isSensitive() ) { + tile_vol.setSensitiveDetector(sens); + } + + //Set region, limitset, and visibility settings + tile_vol.setAttributes(lcdd,tile_xml.regionStr(),tile_xml.limitsStr(),tile_xml.visStr()); + + tiles.push_back(tile_vol); + tile_number++; + } + + //Place the same volumes inside the envelope + float tile_pos_z = -x_det_dim.z()/2.; + int slice_num = 0; + while(tile_pos_z<x_det_dim.z()/2.){ + tile_number=0; + for(xml_coll_t k(x_det_layer,_U(slice)); k; ++k) { + + DD4hep::XML::Component tile_xml = k; + float tile_thickness = tile_xml.dz(); + + //Place the tile inside the layer + PlacedVolume tile_plv = layer_vol.placeVolume(tiles.at(tile_number),Position(0,tile_pos_z,0)); + tile_plv.addPhysVolID("layer",layer_num); + tile_plv.addPhysVolID("slice",slice_num); + + vol_max["layer"] = std::max(vol_max["layer"],layer_num); + vol_max["slice"] = std::max(vol_max["slice"],slice_num); + + //Increment the z pos of the tile + tile_pos_z += tile_thickness; + tile_number++; + slice_num++; + } + } + + //Place the same layer around the beam axis phiBins times + double mod_x_off = r; + double mod_y_off = 0; + for(int i=0;i<x_det_dim.phiBins();i++){ + if(debug) cout << "Layer:" << i << " phi:" << tile_phi << " rotz:" << (tile_phi*i) << endl; + double layer_pos_x = mod_x_off * cos(tile_phi*i) - mod_y_off * sin(tile_phi*i); + double layer_pos_y = mod_x_off * sin(tile_phi*i) + mod_y_off * cos(tile_phi*i); + Transform3D tr(RotationZYX(M_PI*0.5,M_PI*0.5,0)*RotationZYX(0,tile_phi*i,0), + Translation3D(layer_pos_x,layer_pos_y,layer_pos_z)); + PlacedVolume pv = calo_vol.placeVolume(layer_vol,tr); + pv.addPhysVolID("module",i+1); + vol_max["module"] = std::max(vol_max["module"],i+1); + //DetElement sd = i==0 ? stave_det : stave_det.clone(_toString(i,"stave%d")); + } + + r += dr; + layer_num += 1; + } + } + for(map<string,int>::const_iterator i=vol_max.begin(); i!=vol_max.end();++i) + cout << "Volume ID:" << (*i).first << " max:" << (*i).second << endl; + //Place the calo inside the world + + return d_det; +} + +DECLARE_DETELEMENT(FCC_HcalBarrel2,create_detector) diff --git a/examples/ClientTests/src/FCC_HcalBarrel_geo.cpp b/examples/ClientTests/src/FCC_HcalBarrel_geo.cpp index 69ef731f7..c7cbfb3c5 100644 --- a/examples/ClientTests/src/FCC_HcalBarrel_geo.cpp +++ b/examples/ClientTests/src/FCC_HcalBarrel_geo.cpp @@ -24,19 +24,25 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { xml_dim_t x_det_dim(x_det.dimensions()); double inner_r = x_det_dim.rmin(); double outer_r = x_det_dim.rmax(); - // Tube: DDCore/DD4hep/Shapes.h - //Tube calo_shape(inner_r,outer_r,x_det_dim.z(),0.0,2*M_PI/x_det_dim.phiBins()); - // Create the detector mother volume - //Volume calo_vol(x_det.nameStr()+"_envelope",calo_shape,lcdd.air()); Assembly calo_vol(x_det.nameStr()+"_envelope"); + PlacedVolume pv; //Set envelope volume attributes calo_vol.setAttributes(lcdd,x_det.regionStr(),x_det.limitsStr(),x_det.visStr()); +#if 0 //Declare this sensitive detector as a calorimeter + Tube tub(inner_r,outer_r,x_det_dim.z()/2.0,0.0,2*M_PI); + //Volume tub_vol(x_det.nameStr()+"_tube",tub,lcdd.material("PyrexGlass")); + Volume tub_vol(x_det.nameStr()+"_tube",tub,lcdd.material("Iron")); + calo_vol.placeVolume(tub_vol); sens.setType("calorimeter"); + tub_vol.setSensitiveDetector(sens); + d_det.setAttributes(lcdd,tub_vol,x_det.regionStr(),x_det.limitsStr(),x_det.visStr()); +#endif +#if 1 int layer_num = 0; float layer_pos_z = 0; @@ -88,18 +94,17 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { //DetElement layer(layer_name,_toString(layer_num,"layer%d"),x_det.id()); //Fill the volume with tiles - - - int tile_number = 0; vector<Volume> tiles; //Assembly tile_seq(layer_name+"_seq"); Trapezoid tile_seq_shape(x1,x2,x_det_layer.dz(),x_det_layer.dz(),x_det_layer.dr()); Volume tile_seq(layer_name + "_seq",tile_seq_shape,lcdd.air()); - tile_seq.setVisAttributes(lcdd.visAttributes("VisibleGreen")); double total_thickness = 0; //Repeat slices until we reach the end of the calorimeter - for(xml_coll_t k(x_det_layer,_U(slice)); k; ++k) { + int slice_num = 0, tile_number = 0; + + tile_seq.setVisAttributes(lcdd.visAttributes("VisibleGreen")); + for(xml_coll_t k(x_det_layer,_U(slice)); k; ++k, ++slice_num) { xml_comp_t tile_xml = k; string tile_name = layer_name + _toString(tile_number,"_slice%d"); Material tile_material = lcdd.material(tile_xml.materialStr()); @@ -108,49 +113,36 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { float tile_y2 = tile_thickness; float tile_z = x_det_layer.dr(); - //Shape a Trapezoid (tile): DDCore/DD4hep/Shapes.h Trapezoid tile_shape(x1,x2,tile_y1,tile_y2,tile_z); - - //Create a volume with trapezoid shape Volume tile_vol(tile_name,tile_shape,tile_material); - tile_seq.placeVolume(tile_vol,Position(0,total_thickness,0)); + pv = tile_seq.placeVolume(tile_vol,Position(0,total_thickness,0)); + pv.addPhysVolID("slice",slice_num); total_thickness += tile_thickness; if ( tile_xml.isSensitive() ) { cout << "Set volume " << tile_name << " sensitive...." << endl; tile_vol.setSensitiveDetector(sens); } - //Set region, limitset, and visibility settings + // Set region, limitset, and visibility settings tile_vol.setAttributes(lcdd,tile_xml.regionStr(),tile_xml.limitsStr(),tile_xml.visStr()); tiles.push_back(tile_vol); tile_number++; } - //Place the same volumes inside the envelope + // Place the same volumes inside the envelope float tile_pos_z = -x_det_dim.z()/2.; + int tile_num = 0; while(tile_pos_z<x_det_dim.z()/2.){ - layer_vol.placeVolume(tile_seq,Position(0,tile_pos_z,0)); + pv = layer_vol.placeVolume(tile_seq,Position(0,tile_pos_z,0)); + pv.addPhysVolID("tile",tile_num); tile_pos_z += total_thickness; -#if 0 - for(xml_coll_t k(x_det_layer,_U(slice)); k; ++k) { - xml_comp_t tile_xml = k; - float tile_thickness = tile_xml.dz(); - - //Place the tile inside the layer - PlacedVolume tile_plv = layer_vol.placeVolume(tiles.at(tile_number),Position(0,tile_pos_z,0)); - //Increment the z pos of the tile - tile_pos_z += tile_thickness; - tile_number++; - } -#endif + tile_num++; } - //Place the same layer around the beam axis phiBins times + // Place the same layer around the beam axis phiBins times Transform3D tr(RotationZYX(M_PI*0.5,M_PI*0.5,0),Translation3D(r,0,layer_pos_z)); PlacedVolume pv = stave_vol.placeVolume(layer_vol,tr); - //pv.addPhysVolID("layer",i+1); - //DetElement sd = i==0 ? stave_det : stave_det.clone(_toString(i,"stave%d")); - + pv.addPhysVolID("layer",layer_num); r += dr; cout << "+++ R=" << r << endl; } @@ -168,9 +160,8 @@ static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sens) { pv.addPhysVolID("stave",i+1); } - cout << "Number of layers: " << layer_num << endl; - +#endif //Place the calo inside the world PlacedVolume calo_plv = lcdd.pickMotherVolume(d_det).placeVolume(calo_vol); calo_plv.addPhysVolID("system",x_det.id()); diff --git a/examples/ClientTests/src/ZylinderShell_geo.cpp b/examples/ClientTests/src/ZylinderShell_geo.cpp new file mode 100644 index 000000000..d1514460b --- /dev/null +++ b/examples/ClientTests/src/ZylinderShell_geo.cpp @@ -0,0 +1,51 @@ +// $Id: PolyconeSupport_geo.cpp 941 2013-12-12 18:47:03Z markus.frank@cern.ch $ +//==================================================================== +// AIDA Detector description implementation for LCD +//-------------------------------------------------------------------- +// +// Author : M.Frank +// +//==================================================================== +#include "DD4hep/DetFactoryHelper.h" + +using namespace std; +using namespace DD4hep; +using namespace DD4hep::Geometry; + +static Ref_t create_detector(LCDD& lcdd, xml_h e, SensitiveDetector sensitive) { + xml_det_t x_det = e; + string name = x_det.nameStr(); + DetElement sdet (name,x_det.id()); + Assembly assembly(name+"_assembly"); + Material mat (lcdd.material(x_det.materialStr())); + double gap = 200; + PlacedVolume pv; + + for(xml_coll_t m(e,_U(module)); m; ++m) { + xml_comp_t mod = m; + vector<double> rmin,rmax,z; + int num = 0; + for(xml_coll_t c(m,_U(zplane)); c; ++c, ++num) { + xml_comp_t dim(c); + rmin.push_back(dim.rmin()); + rmax.push_back(dim.rmax()); + z.push_back(dim.z()/2); + } + if ( num < 2 ) { + throw runtime_error("ZylinderShell["+name+"]> Not enough Z planes. minimum is 2!"); + } + Polycone cone (0.,2*M_PI,rmin,rmax,z); + Volume volume(name, cone, mat); + volume.setVisAttributes(lcdd, x_det.visStr()); + volume.setSensitiveDetector(sensitive); + pv = assembly.placeVolume(volume); + pv.addPhysVolID("barrel",mod.id()); + } + + pv = lcdd.pickMotherVolume(sdet).placeVolume(assembly); + pv.addPhysVolID("system",x_det.id()); + sdet.setPlacement(pv); + return sdet; +} + +DECLARE_DETELEMENT(ZylinderShell,create_detector); -- GitLab