DD4hep  1.30.0
Detector Description Toolkit for High Energy Physics
Geant4MaterialScanner.cpp
Go to the documentation of this file.
1 //==========================================================================
2 // AIDA Detector description implementation
3 //--------------------------------------------------------------------------
4 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
5 // All rights reserved.
6 //
7 // For the licensing terms see $DD4hepINSTALL/LICENSE.
8 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
9 //
10 // Author : M.Frank
11 //
12 //==========================================================================
13 
14 // Framework include files
15 #include <DD4hep/Objects.h>
16 #include <DDG4/Defs.h>
18 
19 // Forward declarations
20 class G4LogicalVolume;
21 
23 namespace dd4hep {
24 
26  namespace sim {
27 
29 
36  protected:
38  class StepInfo {
39  public:
43  const G4LogicalVolume* volume;
44 
46  StepInfo(const Position& pre, const Position& post, const G4LogicalVolume* volume);
48  StepInfo(const StepInfo& c);
50  ~StepInfo() {}
52  StepInfo& operator=(const StepInfo& c);
53  };
54  typedef std::vector<StepInfo*> Steps;
55 
56  double m_sumX0 = 0E0;
57  double m_sumLambda = 0E0;
58  double m_sumPath = 0E0;
60 
61  public:
63  Geant4MaterialScanner(Geant4Context* context, const std::string& name);
65  virtual ~Geant4MaterialScanner();
67  virtual void operator()(const G4Step* step, G4SteppingManager* mgr);
69  virtual void begin(const G4Track* track);
71  virtual void end(const G4Track* track);
73  void beginEvent(const G4Event* event);
74  };
75  }
76 }
77 
78 //====================================================================
79 // AIDA Detector description implementation
80 //--------------------------------------------------------------------
81 //
82 // Author : M.Frank
83 //
84 //====================================================================
85 
86 // Framework include files
87 #include <DD4hep/InstanceCount.h>
88 #include <DD4hep/Printout.h>
90 #include <DDG4/Geant4StepHandler.h>
91 #include <DDG4/Geant4EventAction.h>
93 #include <CLHEP/Units/SystemOfUnits.h>
94 #include <G4LogicalVolume.hh>
95 #include <G4Material.hh>
96 
97 using namespace dd4hep::sim;
98 
99 #include <DDG4/Factories.h>
101 
102 Geant4MaterialScanner::StepInfo::StepInfo(const Position& prePos, const Position& postPos, const G4LogicalVolume* vol)
104 : pre(prePos), post(postPos), volume(vol)
105 {
106 }
107 
110  : pre(c.pre), post(c.post), volume(c.volume)
111 {
112 }
113 
116  pre = c.pre;
117  post = c.post;
118  volume = c.volume;
119  return *this;
120 }
121 
124  : Geant4SteppingAction(ctxt,nam)
125 {
126  m_needsControl = true;
131 }
132 
136 }
137 
139 void Geant4MaterialScanner::operator()(const G4Step* step, G4SteppingManager*) {
140  Geant4StepHandler h(step);
141 #if 0
142  Geant4TouchableHandler pre_handler(step);
143  std::string prePath = pre_handler.path();
144  Geant4TouchableHandler post_handler(step);
145  std::string postPath = post_handler.path();
146 #endif
147  G4LogicalVolume* logVol = h.logvol(h.pre);
148  m_steps.emplace_back(new StepInfo(h.prePos(), h.postPos(), logVol));
149 }
150 
152 void Geant4MaterialScanner::beginEvent(const G4Event* /* event */) {
153  for_each(m_steps.begin(),m_steps.end(),detail::DestroyObject<StepInfo*>());
154  m_steps.clear();
155  m_sumX0 = 0;
156  m_sumLambda = 0;
157  m_sumPath = 0;
158 }
159 
161 void Geant4MaterialScanner::begin(const G4Track* track) {
162  printP2("Starting tracking action for track ID=%d",track->GetTrackID());
163  for_each(m_steps.begin(),m_steps.end(),detail::DestroyObject<StepInfo*>());
164  m_steps.clear();
165  m_sumX0 = 0;
166  m_sumLambda = 0;
167  m_sumPath = 0;
168 }
169 
171 void Geant4MaterialScanner::end(const G4Track* track) {
172  if ( !m_steps.empty() ) {
173  constexpr const char* line = " +--------------------------------------------------------------------------------------------------------------------------------------------------\n";
174  constexpr const char* fmt1 = " | %5d %-20s %3.0f %8.3f %8.4f %11.4f %11.4f %10.3f %8.2f %11.6f %11.6f (%7.2f,%7.2f,%7.2f)\n";
175  constexpr const char* fmt2 = " | %5d %-20s %3.0f %8.3f %8.4f %11.6g %11.6g %10.3f %8.2f %11.6f %11.6f (%7.2f,%7.2f,%7.2f)\n";
176  const Position& pre = m_steps[0]->pre;
177  const Position& post = m_steps[m_steps.size()-1]->post;
178 
179  ::printf("%s + Material scan between: x_0 = (%7.2f,%7.2f,%7.2f) [cm] and x_1 = (%7.2f,%7.2f,%7.2f) [cm] TrackID:%d: \n%s",
180  line,pre.X()/CLHEP::cm,pre.Y()/CLHEP::cm,pre.Z()/CLHEP::cm,post.X()/CLHEP::cm,post.Y()/CLHEP::cm,post.Z()/CLHEP::cm,track->GetTrackID(),line);
181  ::printf(" | \\ %-11s Atomic Radiation Interaction Path Integrated Integrated Material\n","Material");
182  ::printf(" | Num. \\ %-11s Number/Z Mass/A Density Length Length Thickness Length X0 Lambda Endpoint \n","Name");
183  ::printf(" | Layer \\ %-11s [g/mole] [g/cm3] [cm] [cm] [cm] [cm] [cm] [cm] ( cm, cm, cm)\n","");
184  ::printf("%s",line);
185  int count = 1;
186  for(Steps::const_iterator i=m_steps.begin(); i!=m_steps.end(); ++i, ++count) {
187  const G4LogicalVolume* logVol = (*i)->volume;
188  G4Material* material = logVol->GetMaterial();
189  const Position& prePos = (*i)->pre;
190  const Position& postPos = (*i)->post;
191  Position direction = postPos - prePos;
192  double length = direction.R()/CLHEP::cm;
193  double intLen = material->GetNuclearInterLength()/CLHEP::cm;
194  double radLen = material->GetRadlen()/CLHEP::cm;
195  double density = material->GetDensity()/(CLHEP::gram/CLHEP::cm3);
196  double nLambda = length / intLen;
197  double nx0 = length / radLen;
198  double Aeff = 0.0;
199  double Zeff = 0.0;
200  const char* fmt = radLen >= 1e5 ? fmt2 : fmt1;
201  const double* fractions = material->GetFractionVector();
202  for(size_t j=0; j<material->GetNumberOfElements(); ++j) {
203  Zeff += fractions[j]*(material->GetElement(j)->GetZ());
204  Aeff += fractions[j]*(material->GetElement(j)->GetA())/CLHEP::gram;
205  }
206  m_sumX0 += nx0;
207  m_sumLambda += nLambda;
208  m_sumPath += length;
209  ::printf(fmt,count,material->GetName().c_str(),
210  Zeff, Aeff, density, radLen, intLen, length,
212  postPos.X()/CLHEP::cm,postPos.Y()/CLHEP::cm,postPos.Z()/CLHEP::cm);
213  //cout << *m << endl;
214  }
215  for_each(m_steps.begin(),m_steps.end(),detail::DestroyObject<StepInfo*>());
216  m_steps.clear();
217  }
218 }
dd4hep::sim::Geant4Action::printP2
void printP2(const char *fmt,...) const
Support for messages with variable output level using output level+2.
Definition: Geant4Action.cpp:188
dd4hep::sim::Geant4MaterialScanner::StepInfo::~StepInfo
~StepInfo()
Default destructor.
Definition: Geant4MaterialScanner.cpp:50
dd4hep::sim::Geant4TouchableHandler
Helper class to ease the extraction of information from a G4Touchable object.
Definition: Geant4TouchableHandler.h:44
dd4hep::sim::Geant4Action::m_needsControl
bool m_needsControl
Default property: Flag to create control instance.
Definition: Geant4Action.h:123
Objects.h
dd4hep::sim::Geant4StepHandler::prePos
Position prePos() const
Returns the pre-step position.
Definition: Geant4StepHandler.h:84
dd4hep::sim::Geant4MaterialScanner::m_sumX0
double m_sumX0
Definition: Geant4MaterialScanner.cpp:56
dd4hep::sim::Geant4SteppingAction
Concrete implementation of the Geant4 stepping action sequence.
Definition: Geant4SteppingAction.h:40
dd4hep::sim::Geant4MaterialScanner::m_sumLambda
double m_sumLambda
Definition: Geant4MaterialScanner.cpp:57
dd4hep::sim::Geant4StepHandler::pre
G4StepPoint * pre
Definition: Geant4StepHandler.h:49
dd4hep::sim::Geant4MaterialScanner::Geant4MaterialScanner
Geant4MaterialScanner(Geant4Context *context, const std::string &name)
Standard constructor.
Definition: Geant4MaterialScanner.cpp:123
Geant4TrackingAction.h
Geant4EventAction.h
dd4hep::sim::Geant4MaterialScanner::StepInfo::volume
const G4LogicalVolume * volume
Reference to the logical volue.
Definition: Geant4MaterialScanner.cpp:43
dd4hep::InstanceCount::increment
static void increment(T *)
Increment count according to type information.
Definition: InstanceCount.h:98
dd4hep::sim::Geant4MaterialScanner::StepInfo
Structure to hold the information of one simulation step.
Definition: Geant4MaterialScanner.cpp:38
dd4hep::sim::Geant4TrackingActionSequence::callAtBegin
void callAtBegin(Q *p, void(T::*f)(const G4Track *), CallbackSequence::Location where=CallbackSequence::END)
Register Pre-track action callback.
Definition: Geant4TrackingAction.h:150
Geant4StepHandler.h
DECLARE_GEANT4ACTION
#define DECLARE_GEANT4ACTION(name)
Plugin defintion to create Geant4Action objects.
Definition: Factories.h:210
Geant4SteppingAction.h
dd4hep::sim::Geant4MaterialScanner
Class to perform directional material scans using Geantinos.
Definition: Geant4MaterialScanner.cpp:35
dd4hep::sim::Geant4StepHandler
Helper class to ease the extraction of information from a G4Step object.
Definition: Geant4StepHandler.h:46
dd4hep::sim::Geant4MaterialScanner::end
virtual void end(const G4Track *track)
End-of-tracking callback.
Definition: Geant4MaterialScanner.cpp:171
dd4hep::sim::Geant4Action::eventAction
Geant4EventActionSequence & eventAction() const
Access to the main event action sequence from the kernel object.
Definition: Geant4Action.cpp:283
dd4hep::sim::Geant4EventActionSequence::callAtBegin
void callAtBegin(Q *p, void(T::*f)(const G4Event *))
Register begin-of-event callback.
Definition: Geant4EventAction.h:152
dd4hep::sim::Geant4MaterialScanner::beginEvent
void beginEvent(const G4Event *event)
Registered callback on Begin-event.
Definition: Geant4MaterialScanner.cpp:152
dd4hep::InstanceCount::decrement
static void decrement(T *)
Decrement count according to type information.
Definition: InstanceCount.h:102
Geant4TouchableHandler.h
dd4hep::sim::Geant4MaterialScanner::m_sumPath
double m_sumPath
Definition: Geant4MaterialScanner.cpp:58
dd4hep::sim::Geant4Action::name
const std::string & name() const
Access name of the action.
Definition: Geant4Action.h:280
dd4hep::sim::Geant4MaterialScanner::operator()
virtual void operator()(const G4Step *step, G4SteppingManager *mgr)
User stepping callback.
Definition: Geant4MaterialScanner.cpp:139
dd4hep::sim::Geant4MaterialScanner::~Geant4MaterialScanner
virtual ~Geant4MaterialScanner()
Default destructor.
Definition: Geant4MaterialScanner.cpp:134
dd4hep::sim::Geant4StepHandler::postPos
Position postPos() const
Returns the post-step position.
Definition: Geant4StepHandler.h:93
dd4hep::sim::Geant4MaterialScanner::begin
virtual void begin(const G4Track *track)
Begin-of-tracking callback.
Definition: Geant4MaterialScanner.cpp:161
dd4hep::sim::Geant4MaterialScanner::Steps
std::vector< StepInfo * > Steps
Definition: Geant4MaterialScanner.cpp:54
dd4hep::sim::Geant4MaterialScanner::StepInfo::pre
Position pre
Pre-step and Post-step position.
Definition: Geant4MaterialScanner.cpp:41
dd4hep::Position
ROOT::Math::XYZVector Position
Definition: Objects.h:81
Factories.h
dd4hep::sim::Geant4MaterialScanner::StepInfo::operator=
StepInfo & operator=(const StepInfo &c)
Assignment operator.
Definition: Geant4MaterialScanner.cpp:115
dd4hep::sim
Namespace for the Geant4 based simulation part of the AIDA detector description toolkit.
Definition: Geant4Output2EDM4hep.cpp:49
dd4hep::sim::Geant4TouchableHandler::path
std::string path() const
Helper: Access the placement path of a Geant4 touchable object as a string.
Definition: Geant4TouchableHandler.cpp:58
dd4hep
Namespace for the AIDA detector description toolkit.
Definition: AlignmentsCalib.h:28
dd4hep::sim::Geant4MaterialScanner::m_steps
Steps m_steps
Definition: Geant4MaterialScanner.cpp:59
dd4hep::sim::Geant4MaterialScanner::StepInfo::post
Position post
Definition: Geant4MaterialScanner.cpp:41
InstanceCount.h
Defs.h
dd4hep::sim::Geant4MaterialScanner::StepInfo::StepInfo
StepInfo(const Position &pre, const Position &post, const G4LogicalVolume *volume)
Initializing constructor.
Definition: Geant4MaterialScanner.cpp:103
Printout.h
dd4hep::sim::Geant4Context
Generic context to extend user, run and event information.
Definition: Geant4Context.h:201
dd4hep::sim::Geant4StepHandler::logvol
G4LogicalVolume * logvol(const G4StepPoint *p) const
Definition: Geant4StepHandler.h:160
dd4hep::sim::Geant4Action::context
Geant4Context * context() const
Access the context.
Definition: Geant4Action.h:270
dd4hep::sim::Geant4Action::trackingAction
Geant4TrackingActionSequence & trackingAction() const
Access to the main tracking action sequence from the kernel object.
Definition: Geant4Action.cpp:293
dd4hep::sim::Geant4TrackingActionSequence::callAtEnd
void callAtEnd(Q *p, void(T::*f)(const G4Track *), CallbackSequence::Location where=CallbackSequence::END)
Register Post-track action callback.
Definition: Geant4TrackingAction.h:156