DD4hep - The AIDA detector description toolkit for high energy physics experiments
DD4hep  Rev:Unversioneddirectory
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Geant4MaterialScanner.cpp
Go to the documentation of this file.
1 // $Id: $
2 //==========================================================================
3 // AIDA Detector description implementation for LCD
4 //--------------------------------------------------------------------------
5 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
6 // All rights reserved.
7 //
8 // For the licensing terms see $DD4hepINSTALL/LICENSE.
9 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
10 //
11 // Author : M.Frank
12 //
13 //==========================================================================
14 
15 // Framework include files
16 #include "DD4hep/Objects.h"
17 #include "DDG4/Defs.h"
19 
20 // Forward declarations
21 class G4LogicalVolume;
22 
24 namespace DD4hep {
25 
27  namespace Simulation {
28 
30 
37  protected:
39  class StepInfo {
40  public:
44  const G4LogicalVolume* volume;
45 
47  StepInfo(const Position& pre, const Position& post, const G4LogicalVolume* volume);
49  StepInfo(const StepInfo& c);
51  ~StepInfo() {}
53  StepInfo& operator=(const StepInfo& c);
54  };
55  typedef std::vector<StepInfo*> Steps;
56 
57  double m_sumX0;
58  double m_sumLambda;
59  double m_sumPath;
61 
62  public:
64  Geant4MaterialScanner(Geant4Context* context, const std::string& name);
66  virtual ~Geant4MaterialScanner();
68  virtual void operator()(const G4Step* step, G4SteppingManager* mgr);
70  virtual void begin(const G4Track* track);
72  virtual void end(const G4Track* track);
74  void beginEvent(const G4Event* event);
75  };
76  }
77 }
78 
79 // $Id: $
80 //====================================================================
81 // AIDA Detector description implementation for LCD
82 //--------------------------------------------------------------------
83 //
84 // Author : M.Frank
85 //
86 //====================================================================
87 
88 // Framework include files
89 #include "DD4hep/InstanceCount.h"
90 #include "DD4hep/Printout.h"
92 #include "DDG4/Geant4StepHandler.h"
93 #include "DDG4/Geant4EventAction.h"
95 #include "CLHEP/Units/SystemOfUnits.h"
96 #include "G4LogicalVolume.hh"
97 #include "G4Material.hh"
98 
99 using namespace std;
100 using namespace DD4hep::Simulation;
101 
102 #include "DDG4/Factories.h"
104 
105 Geant4MaterialScanner::StepInfo::StepInfo(const Position& prePos, const Position& postPos, const G4LogicalVolume* vol)
107 : pre(prePos), post(postPos), volume(vol)
108 {
109 }
110 
112 Geant4MaterialScanner::StepInfo::StepInfo(const StepInfo& c)
113 : pre(c.pre), post(c.post), volume(c.volume)
114 {
115 }
116 
119  pre = c.pre;
120  post = c.post;
121  volume = c.volume;
122  return *this;
123 }
124 
127 : Geant4SteppingAction(ctxt,nam)
128 {
129  m_needsControl = true;
134 }
135 
139 }
140 
142 void Geant4MaterialScanner::operator()(const G4Step* step, G4SteppingManager*) {
143  Geant4StepHandler h(step);
144 #if 0
145  Geant4TouchableHandler pre_handler(step);
146  string prePath = pre_handler.path();
147  Geant4TouchableHandler post_handler(step);
148  string postPath = post_handler.path();
149 #endif
150  G4LogicalVolume* logVol = h.logvol(h.pre);
151  m_steps.push_back(new StepInfo(h.prePos(), h.postPos(), logVol));
152 }
153 
155 void Geant4MaterialScanner::beginEvent(const G4Event* /* event */) {
156  for_each(m_steps.begin(),m_steps.end(),DestroyObject<StepInfo*>());
157  m_steps.clear();
158  m_sumX0 = 0;
159  m_sumLambda = 0;
160  m_sumPath = 0;
161 }
162 
164 void Geant4MaterialScanner::begin(const G4Track* track) {
165  printP2("Starting tracking action for track ID=%d",track->GetTrackID());
166  for_each(m_steps.begin(),m_steps.end(),DestroyObject<StepInfo*>());
167  m_steps.clear();
168  m_sumX0 = 0;
169  m_sumLambda = 0;
170  m_sumPath = 0;
171 }
172 
174 void Geant4MaterialScanner::end(const G4Track* track) {
175  using namespace CLHEP;
176  if ( !m_steps.empty() ) {
177  const char* line = " +--------------------------------------------------------------------------------------------------------------------------------------------------\n";
178  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";
179  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";
180  const Position& pre = m_steps[0]->pre;
181  const Position& post = m_steps[m_steps.size()-1]->post;
182 
183  ::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",
184  line,pre.X()/cm,pre.Y()/cm,pre.Z()/cm,post.X()/cm,post.Y()/cm,post.Z()/cm,track->GetTrackID(),line);
185  ::printf(" | \\ %-11s Atomic Radiation Interaction Path Integrated Integrated Material\n","Material");
186  ::printf(" | Num. \\ %-11s Number/Z Mass/A Density Length Length Thickness Length X0 Lambda Endpoint \n","Name");
187  ::printf(" | Layer \\ %-11s [g/mole] [g/cm3] [cm] [cm] [cm] [cm] [cm] [cm] ( cm, cm, cm)\n","");
188  ::printf("%s",line);
189  int count = 1;
190  for(Steps::const_iterator i=m_steps.begin(); i!=m_steps.end(); ++i, ++count) {
191  const G4LogicalVolume* logVol = (*i)->volume;
192  G4Material* material = logVol->GetMaterial();
193  const Position& prePos = (*i)->pre;
194  const Position& postPos = (*i)->post;
195  Position direction = postPos - prePos;
196  double length = direction.R()/cm;
197  double intLen = material->GetNuclearInterLength()/cm;
198  double radLen = material->GetRadlen()/cm;
199  double density = material->GetDensity()/(gram/cm3);
200  double nLambda = length / intLen;
201  double nx0 = length / radLen;
202  double Aeff = 0.0;
203  double Zeff = 0.0;
204  const char* fmt = radLen >= 1e5 ? fmt2 : fmt1;
205  const double* fractions = material->GetFractionVector();
206  for(size_t j=0; j<material->GetNumberOfElements(); ++j) {
207  Zeff += fractions[j]*(material->GetElement(j)->GetZ());
208  Aeff += fractions[j]*(material->GetElement(j)->GetA())/gram;
209  }
210  m_sumX0 += nx0;
211  m_sumLambda += nLambda;
212  m_sumPath += length;
213  ::printf(fmt,count,material->GetName().c_str(),
214  Zeff, Aeff, density, radLen, intLen, length,
216  postPos.X()/cm,postPos.Y()/cm,postPos.Z()/cm);
217  //cout << *m << endl;
218  }
219  for_each(m_steps.begin(),m_steps.end(),DestroyObject<StepInfo*>());
220  m_steps.clear();
221  }
222 }
void callAtEnd(Q *p, void(T::*f)(const G4Track *), CallbackSequence::Location where=CallbackSequence::END)
Register Post-track action callback.
Class to perform directional material scans using Geantinos.
void beginEvent(const G4Event *event)
Registered callback on Begin-event.
Geant4TrackingActionSequence & trackingAction() const
Access to the main tracking action sequence from the kernel object.
std::string path() const
Helper: Access the placement path of a Geant4 touchable object as a string.
#define DECLARE_GEANT4ACTION(name)
Plugin defintion to create Geant4Action objects.
Definition: Factories.h:166
const G4LogicalVolume * volume
Reference to the logical volue.
void callAtBegin(Q *p, void(T::*f)(const G4Track *), CallbackSequence::Location where=CallbackSequence::END)
Register Pre-track action callback.
static void decrement(T *)
Decrement count according to type information.
Position postPos() const
Returns the post-step position.
virtual ~Geant4MaterialScanner()
Default destructor.
Position pre
Pre-step and Post-step position.
Helper class to ease the extraction of information from a G4Step object.
virtual void begin(const G4Track *track)
Begin-of-tracking callback.
Geant4Context * context() const
Access the context.
Definition: Geant4Action.h:261
void printP2(const char *fmt,...) const
Support for messages with variable output level using output level+2.
const std::string & name() const
Access name of the action.
Definition: Geant4Action.h:271
void callAtBegin(Q *p, void(T::*f)(const G4Event *))
Register begin-of-event callback.
Concrete implementation of the Geant4 stepping action sequence.
Functor to delete objects from heap and reset the pointer.
Definition: Primitives.h:249
virtual void end(const G4Track *track)
End-of-tracking callback.
Position prePos() const
Returns the pre-step position.
virtual void operator()(const G4Step *step, G4SteppingManager *mgr)
User stepping callback.
static const double cm3
Definition: DD4hepUnits.h:75
Helper class to ease the extraction of information from a G4Touchable object.
ROOT::Math::XYZVector Position
Definition: Objects.h:75
static const double gram
Definition: DD4hepUnits.h:157
bool m_needsControl
Default property: Flag to create control instance.
Definition: Geant4Action.h:101
static const double cm
Definition: DD4hepUnits.h:73
static void increment(T *)
Increment count according to type information.
Definition: InstanceCount.h:98
Generic context to extend user, run and event information.
G4LogicalVolume * logvol(const G4StepPoint *p) const
Geant4EventActionSequence & eventAction() const
Access to the main event action sequence from the kernel object.
Geant4MaterialScanner(Geant4Context *context, const std::string &name)
Standard constructor.
Structure to hold the information of one simulation step.
StepInfo & operator=(const StepInfo &c)
Assignment operator.
StepInfo(const Position &pre, const Position &post, const G4LogicalVolume *volume)
Initializing constructor.