14 #ifndef DD4HEP_GEANT4STEPHANDLER_H
15 #define DD4HEP_GEANT4STEPHANDLER_H
22 #include "G4StepPoint.hh"
23 #include "G4VTouchable.hh"
24 #include "G4VSensitiveDetector.hh"
25 #include "G4EmSaturation.hh"
26 #include "G4Version.hh"
32 namespace Simulation {
35 class Geant4StepHandler;
58 pre = s->GetPreStepPoint();
59 post = s->GetPostStepPoint();
60 track = s->GetTrack();
67 return track->GetDefinition();
70 return track->GetDefinition()->GetPDGEncoding();
73 static const char*
stepStatus(G4StepStatus status);
83 return step->GetTotalEnergyDeposit();
87 const G4ThreeVector& p =
pre->GetPosition();
88 return Position(p.x(), p.y(), p.z());
92 return pre->GetPosition();
96 const G4ThreeVector& p =
post->GetPosition();
97 return Position(p.x(), p.y(), p.z());
101 return post->GetPosition();
105 const G4ThreeVector& p1 =
pre->GetPosition();
106 const G4ThreeVector& p2 =
post->GetPosition();
107 G4ThreeVector r = 0.5*(p2+p1);
112 const G4ThreeVector& p1 =
pre->GetPosition();
113 const G4ThreeVector& p2 =
post->GetPosition();
114 G4ThreeVector r = (p2-p1);
123 const G4ThreeVector& p =
pre->GetMomentum();
124 return Momentum(p.x(), p.y(), p.z());
127 const G4ThreeVector& p =
post->GetMomentum();
128 return Momentum(p.x(), p.y(), p.z());
131 const G4ThreeVector& p =
track->GetMomentum();
132 return Momentum(p.x(), p.y(), p.z());
135 return step->GetTotalEnergyDeposit();
138 return step->GetStepLength();
141 return track->GetTrackID();
144 return track->GetParentID();
147 return track->GetGlobalTime();
150 return track->GetTotalEnergy();
153 return track->GetKineticEnergy();
156 return track->GetTrackStatus();
159 return track->GetTrackStatus() == fAlive;
162 return track->GetCreatorProcess();
165 return pre->GetTouchable();
168 return post->GetTouchable();
170 const char*
volName(
const G4StepPoint* p,
const char* undefined =
"")
const {
171 G4VPhysicalVolume*
v =
volume(p);
172 return v ? v->GetName().c_str() : undefined;
174 G4VPhysicalVolume*
volume(
const G4StepPoint* p)
const {
175 return p->GetTouchableHandle()->GetVolume();
177 G4VSolid*
solid(
const G4StepPoint* p)
const {
178 return p->GetTouchableHandle()->GetSolid();
180 G4VPhysicalVolume*
physvol(
const G4StepPoint* p)
const {
181 return p->GetPhysicalVolume();
183 G4LogicalVolume*
logvol(
const G4StepPoint* p)
const {
184 G4VPhysicalVolume* pv =
physvol(p);
185 return pv ? pv->GetLogicalVolume() : 0;
188 G4LogicalVolume* lv =
logvol(p);
189 return lv ? lv->GetSensitiveDetector() : 0;
191 const char*
sdName(
const G4StepPoint* p,
const char* undefined =
"")
const {
193 return s ? s->GetName().c_str() : undefined;
196 return lv ? (0 != lv->GetSensitiveDetector()) :
false;
199 return pv ?
isSensitive(pv->GetLogicalVolume()) :
false;
217 return step->IsFirstStepInVolume();
220 return step->IsLastStepInVolume();
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.
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
Helper class to ease the extraction of information from a G4Step object.
const char * postStepStatus() const
Returns the post-step status in form of a string.
const G4VTouchable * preTouchable() const
Position prePos() const
Returns the pre-step position.
double trkKineEnergy() const
G4VPhysicalVolume * postVolume() const
G4ParticleDefinition * trackDef() const
G4VSensitiveDetector * sd(const G4StepPoint *p) const
Simple container for a physics vector.
double stepLength() const
G4VPhysicalVolume * preVolume() const
ROOT::Math::XYZVector Position
Position direction() const
Direction calculated from the post- and pre-position ofthe step.
bool isSensitive(const G4VPhysicalVolume *pv) const
bool lastInVolume() 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.
G4LogicalVolume * logvol(const G4StepPoint *p) const
bool firstInVolume() 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
Geant4StepHandler & operator=(const Geant4StepHandler ©)=delete
Assignment operator inhibited. Should not be copied.
Geant4StepHandler(const G4Step *s)
Initializing constructor.