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.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 "DDG4/Geant4StepHandler.h"
18 #include "DD4hep/DD4hepUnits.h"
19 
20 using namespace DD4hep::Simulation;
21 
23 const char* Geant4StepHandler::stepStatus(G4StepStatus status) {
24  switch (status) {
25  // Step reached the world boundary
26  case fWorldBoundary:
27  return "WorldBoundary";
28  // Step defined by a geometry boundary
29  case fGeomBoundary:
30  return "GeomBoundary";
31  // Step defined by a PreStepDoItVector
32  case fAtRestDoItProc:
33  return "AtRestDoItProc";
34  // Step defined by a AlongStepDoItVector
35  case fAlongStepDoItProc:
36  return "AlongStepDoItProc";
37  // Step defined by a PostStepDoItVector
38  case fPostStepDoItProc:
39  return "PostStepDoItProc";
40  // Step defined by the user Step limit in the logical volume
41  case fUserDefinedLimit:
42  return "UserDefinedLimit";
43  // Step defined by an exclusively forced PostStepDoIt process
44  case fExclusivelyForcedProc:
45  return "ExclusivelyForcedProc";
46  // Step not defined yet
47  case fUndefined:
48  default:
49  return "Undefined";
50  };
51 }
52 
54 const char* Geant4StepHandler::preStepStatus() const {
55  return stepStatus(pre ? pre->GetStepStatus() : fUndefined);
56 }
57 
59 const char* Geant4StepHandler::postStepStatus() const {
60  return stepStatus(post ? post->GetStepStatus() : fUndefined);
61 }
62 
65  return localToGlobal(G4ThreeVector(local.X / dd4hep::mm,local.Y / dd4hep::mm,local.Z / dd4hep::mm));
66 }
67 
70  return localToGlobal(G4ThreeVector(local.X(),local.Y(),local.Z()));
71 }
72 
74 Position Geant4StepHandler::localToGlobal(double x, double y, double z) const {
75  return localToGlobal(G4ThreeVector(x,y,z));
76 }
77 
79 Position Geant4StepHandler::localToGlobal(const G4ThreeVector& loc) const {
80  G4TouchableHandle t = step->GetPreStepPoint()->GetTouchableHandle();
81  G4ThreeVector p = t->GetHistory()->GetTopTransform().Inverse().TransformPoint(loc);
82  return Position(p.x(),p.y(),p.z());
83 }
84 
86 Position Geant4StepHandler::globalToLocal(double x, double y, double z) const {
87  G4ThreeVector p = globalToLocalG4(G4ThreeVector(x,y,z));
88  return Position(p.x(),p.y(),p.z());
89 }
90 
93  G4ThreeVector p = globalToLocalG4(G4ThreeVector(global.X(),global.Y(),global.Z()));
94  return Position(p.x(),p.y(),p.z());
95 }
96 
98 Position Geant4StepHandler::globalToLocal(const G4ThreeVector& global) const {
99  G4ThreeVector p = globalToLocalG4(global);
100  return Position(p.x(),p.y(),p.z());
101 }
102 
104 G4ThreeVector Geant4StepHandler::globalToLocalG4(double x, double y, double z) const {
105  return globalToLocalG4(G4ThreeVector(x,y,z));
106 }
107 
109 G4ThreeVector Geant4StepHandler::globalToLocalG4(const G4ThreeVector& global) const {
110  G4TouchableHandle t = step->GetPreStepPoint()->GetTouchableHandle();
111  return t->GetHistory()->GetTopTransform().TransformPoint(global);
112 }
113 
116 #if G4VERSION_NUMBER >= 1001
117  static G4EmSaturation s_emSaturation(0);
118 #else
119  static G4EmSaturation s_emSaturation();
120 #endif
121 
122  double energyDeposition = step->GetTotalEnergyDeposit();
123  double length = step->GetStepLength();
124  double niel = step->GetNonIonizingEnergyDeposit();
125  const G4Track* trk = step->GetTrack();
126  const G4ParticleDefinition* particle = trk->GetDefinition();
127  const G4MaterialCutsCouple* couple = trk->GetMaterialCutsCouple();
128  double engyVis = s_emSaturation.VisibleEnergyDeposition(particle,
129  couple,
130  length,
131  energyDeposition,
132  niel);
133  return engyVis;
134 }
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.
const char * postStepStatus() const
Returns the post-step status in form of a string.
Simple container for a physics vector.
Definition: Segmentation.h:41
ROOT::Math::XYZVector Position
Definition: Objects.h:75
static const double mm
Definition: DD4hepUnits.h:69
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