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
Geant4StepHandler.h
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 #ifndef DD4HEP_GEANT4STEPHANDLER_H
15 #define DD4HEP_GEANT4STEPHANDLER_H
16 
17 // Framework include files
18 #include "DDG4/Defs.h"
19 
20 // Geant4 include files
21 #include "G4Step.hh"
22 #include "G4StepPoint.hh"
23 #include "G4VTouchable.hh"
24 #include "G4VSensitiveDetector.hh"
25 #include "G4EmSaturation.hh"
26 #include "G4Version.hh"
27 
29 namespace DD4hep {
30 
32  namespace Simulation {
33 
34  // Forward declarations;
35  class Geant4StepHandler;
36 
38 
48  public:
49  const G4Step* step;
50  G4StepPoint* pre;
51  G4StepPoint* post;
52  G4Track* track;
55  Geant4StepHandler() = delete;
57  Geant4StepHandler(const G4Step* s) : step(s) {
58  pre = s->GetPreStepPoint();
59  post = s->GetPostStepPoint();
60  track = s->GetTrack();
61  applyBirksLaw = false;
62  }
65 
66  G4ParticleDefinition* trackDef() const {
67  return track->GetDefinition();
68  }
69  int trkPdgID() const {
70  return track->GetDefinition()->GetPDGEncoding();
71  }
73  static const char* stepStatus(G4StepStatus status);
75  const char* preStepStatus() const;
77  const char* postStepStatus() const;
79  double totalEnergy() const {
80  if(applyBirksLaw == true)
81  return this->birkAttenuation();
82  else
83  return step->GetTotalEnergyDeposit();
84  }
86  Position prePos() const {
87  const G4ThreeVector& p = pre->GetPosition();
88  return Position(p.x(), p.y(), p.z());
89  }
91  const G4ThreeVector& prePosG4() const {
92  return pre->GetPosition();
93  }
95  Position postPos() const {
96  const G4ThreeVector& p = post->GetPosition();
97  return Position(p.x(), p.y(), p.z());
98  }
100  const G4ThreeVector& postPosG4() const {
101  return post->GetPosition();
102  }
105  const G4ThreeVector& p1 = pre->GetPosition();
106  const G4ThreeVector& p2 = post->GetPosition();
107  G4ThreeVector r = 0.5*(p2+p1);
108  return Position(r.x(),r.y(),r.z());
109  }
111  Position direction() const {
112  const G4ThreeVector& p1 = pre->GetPosition();
113  const G4ThreeVector& p2 = post->GetPosition();
114  G4ThreeVector r = (p2-p1);
115  double len = r.r();
116  if ( len > 1e-15 ) {
117  r /= len;
118  return Position(r.x(),r.y(),r.z());
119  }
120  return Position();
121  }
122  Momentum preMom() const {
123  const G4ThreeVector& p = pre->GetMomentum();
124  return Momentum(p.x(), p.y(), p.z());
125  }
126  Momentum postMom() const {
127  const G4ThreeVector& p = post->GetMomentum();
128  return Momentum(p.x(), p.y(), p.z());
129  }
130  Momentum trkMom() const {
131  const G4ThreeVector& p = track->GetMomentum();
132  return Momentum(p.x(), p.y(), p.z());
133  }
134  double deposit() const {
135  return step->GetTotalEnergyDeposit();
136  }
137  double stepLength() const {
138  return step->GetStepLength();
139  }
140  int trkID() const {
141  return track->GetTrackID();
142  }
143  int parentID() const {
144  return track->GetParentID();
145  }
146  double trkTime() const {
147  return track->GetGlobalTime();
148  }
149  double trkEnergy() const {
150  return track->GetTotalEnergy();
151  }
152  double trkKineEnergy() const {
153  return track->GetKineticEnergy();
154  }
155  int trkStatus() const {
156  return track->GetTrackStatus();
157  }
158  bool trkAlive() const {
159  return track->GetTrackStatus() == fAlive;
160  }
161  const G4VProcess* trkProcess() const {
162  return track->GetCreatorProcess();
163  }
164  const G4VTouchable* preTouchable() const {
165  return pre->GetTouchable();
166  }
167  const G4VTouchable* postTouchable() const {
168  return post->GetTouchable();
169  }
170  const char* volName(const G4StepPoint* p, const char* undefined = "") const {
171  G4VPhysicalVolume* v = volume(p);
172  return v ? v->GetName().c_str() : undefined;
173  }
174  G4VPhysicalVolume* volume(const G4StepPoint* p) const {
175  return p->GetTouchableHandle()->GetVolume();
176  }
177  G4VSolid* solid(const G4StepPoint* p) const {
178  return p->GetTouchableHandle()->GetSolid();
179  }
180  G4VPhysicalVolume* physvol(const G4StepPoint* p) const {
181  return p->GetPhysicalVolume();
182  }
183  G4LogicalVolume* logvol(const G4StepPoint* p) const {
184  G4VPhysicalVolume* pv = physvol(p);
185  return pv ? pv->GetLogicalVolume() : 0;
186  }
187  G4VSensitiveDetector* sd(const G4StepPoint* p) const {
188  G4LogicalVolume* lv = logvol(p);
189  return lv ? lv->GetSensitiveDetector() : 0;
190  }
191  const char* sdName(const G4StepPoint* p, const char* undefined = "") const {
192  G4VSensitiveDetector* s = sd(p);
193  return s ? s->GetName().c_str() : undefined;
194  }
195  bool isSensitive(const G4LogicalVolume* lv) const {
196  return lv ? (0 != lv->GetSensitiveDetector()) : false;
197  }
198  bool isSensitive(const G4VPhysicalVolume* pv) const {
199  return pv ? isSensitive(pv->GetLogicalVolume()) : false;
200  }
201  bool isSensitive(const G4StepPoint* point) const {
202  return point ? isSensitive(volume(point)) : false;
203  }
204  G4VPhysicalVolume* preVolume() const {
205  return volume(pre);
206  }
208  return sd(pre);
209  }
210  G4VPhysicalVolume* postVolume() const {
211  return volume(post);
212  }
214  return sd(post);
215  }
216  bool firstInVolume() const {
217  return step->IsFirstStepInVolume();
218  }
219  bool lastInVolume() const {
220  return step->IsLastStepInVolume();
221  }
223 
224  Position localToGlobal(const Position& local) const;
226 
227  Position localToGlobal(const DDSegmentation::Vector3D& local) const;
229  Position localToGlobal(const G4ThreeVector& local) const;
231  Position localToGlobal(double x, double y, double z) const;
232 
234  Position globalToLocal(double x, double y, double z) const;
236  Position globalToLocal(const Position& global) const;
238  Position globalToLocal(const G4ThreeVector& global) const;
240  G4ThreeVector globalToLocalG4(double x, double y, double z) const;
242  G4ThreeVector globalToLocalG4(const G4ThreeVector& loc) const;
243 
245  double birkAttenuation() const;
247  void doApplyBirksLaw(void) {
248  applyBirksLaw = true;
249  }
250  };
251 
252  } // End namespace Simulation
253 } // End namespace DD4hep
254 
255 #endif // DD4HEP_GEANT4STEPHANDLER_H
Position localToGlobal(const Position &local) const
Coordinate transformation to global coordinates.
static const char * stepStatus(G4StepStatus status)
Returns the step status (argument) in form of a string.
Position globalToLocal(double x, double y, double z) const
Coordinate transformation to local coordinates.
G4ThreeVector globalToLocalG4(double x, double y, double z) const
Coordinate transformation to local coordinates.
bool isSensitive(const G4StepPoint *point) const
G4VSolid * solid(const G4StepPoint *p) const
G4VPhysicalVolume * volume(const G4StepPoint *p) const
Class of the Geant4 toolkit. See http://www-geant4.kek.jp/Reference.
Definition: Geant4Classes.h:48
G4VSensitiveDetector * postSD() const
Position postPos() const
Returns the post-step position.
const G4ThreeVector & postPosG4() const
Returns the post-step position as a G4ThreeVector.
bool isSensitive(const G4LogicalVolume *lv) const
const G4VProcess * trkProcess() const
G4VSensitiveDetector * preSD() const
void copy(Alignment from, Alignment to)
Copy alignment object from source object.
TGeoShape * s
Definition: Volumes.cpp:294
Helper class to ease the extraction of information from a G4Step object.
return e
Definition: Volumes.cpp:297
const char * postStepStatus() const
Returns the post-step status in form of a string.
Position Momentum
Definition: Defs.h:36
const G4VTouchable * preTouchable() const
Position prePos() const
Returns the pre-step position.
G4VPhysicalVolume * postVolume() const
G4ParticleDefinition * trackDef() const
G4VSensitiveDetector * sd(const G4StepPoint *p) const
Simple container for a physics vector.
Definition: Segmentation.h:41
G4VPhysicalVolume * preVolume() const
ROOT::Math::XYZVector Position
Definition: Objects.h:75
Position direction() const
Direction calculated from the post- and pre-position ofthe step.
bool isSensitive(const G4VPhysicalVolume *pv) const
const char * volName(const G4StepPoint *p, const char *undefined="") const
double totalEnergy() const
Returns total energy deposit.
void doApplyBirksLaw(void)
Set applyBirksLaw to ture.
View * v
Definition: MultiView.cpp:30
G4LogicalVolume * logvol(const G4StepPoint *p) const
const G4VTouchable * postTouchable() const
const G4ThreeVector & prePosG4() const
Returns the pre-step position as a G4ThreeVector.
Geant4StepHandler()=delete
Inhibit default constructor.
Position avgPosition() const
Average position vector of the step.
const char * sdName(const G4StepPoint *p, const char *undefined="") const
G4VPhysicalVolume * physvol(const G4StepPoint *p) const
const char * preStepStatus() const
Returns the pre-step status in form of a string.
double birkAttenuation() const
Apply BirksLaw.
ROOT::Math::XYZVector Position
Definition: Geant4Data.h:37
Geant4StepHandler & operator=(const Geant4StepHandler &copy)=delete
Assignment operator inhibited. Should not be copied.
Geant4StepHandler(const G4Step *s)
Initializing constructor.