Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4ConeSteppingAction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4ConeSteppingAction.cc
2 #include "PHG4ConeDetector.h"
3 
4 #include <g4main/PHG4Hit.h>
6 #include <g4main/PHG4Hitv1.h>
7 #include <g4main/PHG4Shower.h>
8 #include <g4main/PHG4SteppingAction.h> // for PHG4SteppingAction
10 
11 #include <phool/getClass.h>
12 
13 #include <Geant4/G4ParticleDefinition.hh> // for G4ParticleDefinition
14 #include <Geant4/G4Step.hh>
15 #include <Geant4/G4StepPoint.hh> // for G4StepPoint
16 #include <Geant4/G4StepStatus.hh> // for fGeomBoundary, fUndefined
17 #include <Geant4/G4String.hh> // for G4String
18 #include <Geant4/G4SystemOfUnits.hh> // for cm, nanosecond, GeV
19 #include <Geant4/G4ThreeVector.hh> // for G4ThreeVector
20 #include <Geant4/G4TouchableHandle.hh> // for G4TouchableHandle
21 #include <Geant4/G4Track.hh> // for G4Track
22 #include <Geant4/G4TrackStatus.hh> // for fStopAndKill
23 #include <Geant4/G4Types.hh> // for G4double
24 #include <Geant4/G4VTouchable.hh> // for G4VTouchable
25 #include <Geant4/G4VUserTrackInformation.hh> // for G4VUserTrackInformation
26 
27 #include <iostream>
28 #include <string> // for string, operator+, oper...
29 
30 class G4VPhysicalVolume;
31 class PHCompositeNode;
32 
33 using namespace std;
34 //____________________________________________________________________________..
36  : PHG4SteppingAction(detector->GetName())
37  , detector_(detector)
38  , hits_(nullptr)
39  , hit(nullptr)
40  , saveshower(nullptr)
41 {
42 }
43 
45 {
46  // if the last hit was a zero energie deposit hit, it is just reset
47  // and the memory is still allocated, so we need to delete it here
48  // if the last hit was saved, hit is a nullptr pointer which are
49  // legal to delete (it results in a no operation)
50  delete hit;
51 }
52 
53 //____________________________________________________________________________..
54 bool PHG4ConeSteppingAction::UserSteppingAction(const G4Step* aStep, bool)
55 {
56  // get volume of the current step
57  G4VPhysicalVolume* volume = aStep->GetPreStepPoint()->GetTouchableHandle()->GetVolume();
58 
59  // collect energy and track length step by step
60  G4double edep = aStep->GetTotalEnergyDeposit() / GeV;
61 
62  const G4Track* aTrack = aStep->GetTrack();
63 
64  int layer_id = detector_->get_Layer();
65  // make sure we are in a volume
66  if (detector_->IsInConeActive(volume))
67  {
68  bool geantino = false;
69  // the check for the pdg code speeds things up, I do not want to make
70  // an expensive string compare for every track when we know
71  // geantino or chargedgeantino has pid=0
72  if (aTrack->GetParticleDefinition()->GetPDGEncoding() == 0 &&
73  aTrack->GetParticleDefinition()->GetParticleName().find("geantino") != string::npos)
74  {
75  geantino = true;
76  }
77  G4StepPoint* prePoint = aStep->GetPreStepPoint();
78  G4StepPoint* postPoint = aStep->GetPostStepPoint();
79  // cout << "track id " << aTrack->GetTrackID() << endl;
80  // cout << "time prepoint: " << prePoint->GetGlobalTime() << endl;
81  // cout << "time postpoint: " << postPoint->GetGlobalTime() << endl;
82  switch (prePoint->GetStepStatus())
83  {
84  case fGeomBoundary:
85  case fUndefined:
86  if (!hit)
87  {
88  hit = new PHG4Hitv1();
89  }
90  //here we set the entrance values in cm
91  hit->set_x(0, prePoint->GetPosition().x() / cm);
92  hit->set_y(0, prePoint->GetPosition().y() / cm);
93  hit->set_z(0, prePoint->GetPosition().z() / cm);
94  // time in ns
95  hit->set_t(0, prePoint->GetGlobalTime() / nanosecond);
96  //set the track ID
97  hit->set_trkid(aTrack->GetTrackID());
98  //set the initial energy deposit
99  hit->set_edep(0);
100 
101  if (G4VUserTrackInformation* p = aTrack->GetUserInformation())
102  {
103  if (PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p))
104  {
105  hit->set_trkid(pp->GetUserTrackId());
106  hit->set_shower_id(pp->GetShower()->get_id());
107  saveshower = pp->GetShower();
108  }
109  }
110 
111  break;
112  default:
113  break;
114  }
115  // here we just update the exit values, it will be overwritten
116  // for every step until we leave the volume or the particle
117  // ceases to exist
118  hit->set_x(1, postPoint->GetPosition().x() / cm);
119  hit->set_y(1, postPoint->GetPosition().y() / cm);
120  hit->set_z(1, postPoint->GetPosition().z() / cm);
121 
122  hit->set_t(1, postPoint->GetGlobalTime() / nanosecond);
123  //sum up the energy to get total deposited
124  hit->set_edep(hit->get_edep() + edep);
125  if (geantino)
126  {
127  hit->set_edep(-1); // only energy=0 g4hits get dropped, this way geantinos survive the g4hit compression
128  }
129  if (edep > 0)
130  {
131  if (G4VUserTrackInformation* p = aTrack->GetUserInformation())
132  {
133  if (PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p))
134  {
135  pp->SetKeep(1); // we want to keep the track
136  }
137  }
138  }
139 
140  // if any of these conditions is true this is the last step in
141  // this volume and we need to save the hit
142  // postPoint->GetStepStatus() == fGeomBoundary: track leaves this volume
143  // postPoint->GetStepStatus() == fWorldBoundary: track leaves this world
144  // (not sure if this will ever be the case)
145  // aTrack->GetTrackStatus() == fStopAndKill: track ends
146  if (postPoint->GetStepStatus() == fGeomBoundary || postPoint->GetStepStatus() == fWorldBoundary || aTrack->GetTrackStatus() == fStopAndKill)
147  {
148  // save only hits with energy deposit (or -1 for geantino)
149  if (hit->get_edep())
150  {
151  hits_->AddHit(layer_id, hit);
152  if (saveshower)
153  {
155  }
156  // ownership has been transferred to container, set to null
157  // so we will create a new hit for the next track
158  hit = nullptr;
159  }
160  else
161  {
162  // if this hit has no energy deposit, just reset it for reuse
163  // this means we have to delete it in the dtor. If this was
164  // the last hit we processed the memory is still allocated
165  hit->Reset();
166  }
167  }
168  // return true to indicate the hit was used
169  return true;
170  }
171  else
172  {
173  return false;
174  }
175 }
176 
177 //____________________________________________________________________________..
179 {
180  string hitnodename;
181  if (detector_->SuperDetector() != "NONE")
182  {
183  hitnodename = "G4HIT_" + detector_->SuperDetector();
184  }
185  else
186  {
187  hitnodename = "G4HIT_" + detector_->GetName();
188  }
189 
190  //now look for the map and grab a pointer to it.
191  hits_ = findNode::getClass<PHG4HitContainer>(topNode, hitnodename);
192 
193  // if we do not find the node we need to scream.
194  if (!hits_)
195  {
196  std::cout << "PHG4ConeSteppingAction::SetTopNode - unable to find " << hitnodename << std::endl;
197  }
198 }