ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4EnvelopeSteppingAction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4EnvelopeSteppingAction.cc
2 #include "PHG4EnvelopeDetector.h"
3 
5 #include <g4main/PHG4Hit.h>
6 #include <g4main/PHG4Hitv1.h>
7 #include <g4main/PHG4Shower.h>
8 #include <g4main/PHG4SteppingAction.h> // for PHG4SteppingAction
9 
11 
12 #include <phool/getClass.h>
13 
14 #include <Geant4/G4ParticleDefinition.hh> // for G4ParticleDefinition
15 #include <Geant4/G4ReferenceCountedHandle.hh> // for G4ReferenceCountedHandle
16 #include <Geant4/G4Step.hh>
17 #include <Geant4/G4StepPoint.hh> // for G4StepPoint
18 #include <Geant4/G4StepStatus.hh> // for fGeomBoundary, fUndefined
19 #include <Geant4/G4String.hh> // for G4String
20 #include <Geant4/G4SystemOfUnits.hh> // for mm, m
21 #include <Geant4/G4ThreeVector.hh> // for G4ThreeVector
22 #include <Geant4/G4TouchableHandle.hh> // for G4TouchableHandle
23 #include <Geant4/G4Track.hh> // for G4Track
24 #include <Geant4/G4TrackStatus.hh> // for fStopAndKill
25 #include <Geant4/G4Types.hh> // for G4double
26 #include <Geant4/G4VTouchable.hh> // for G4VTouchable
27 #include <Geant4/G4VUserTrackInformation.hh> // for G4VUserTrackInformation
28 
29 #include <iostream>
30 #include <string> // for string, operator+, ope...
31 
32 class G4VPhysicalVolume;
33 class PHCompositeNode;
34 
35 using namespace std;
36 
37 //______________________________________________________________
39  PHG4SteppingAction(detector->GetName()),
40  detector_( detector ),
41  hits_(nullptr),
42  hit(nullptr)
43 {
44 }
45 
46 //____________________________________________________________________________..
48 {
51 
52  int whichactive = detector_->IsInEnvelope(volume);
53 
54  if ( !whichactive )
55  {
56  return false;
57  }
58 
59  int layer_id = detector_->get_Layer();
60  int tower_id = touch->GetCopyNumber();
61 
62  /* Get energy deposited by this step */
63  //G4double edep = aStep->GetTotalEnergyDeposit() / GeV;
64  G4double eion = (aStep->GetTotalEnergyDeposit() - aStep->GetNonIonizingEnergyDeposit()) / GeV;
65 
66  /* Get pointer to associated Geant4 track */
67  const G4Track* aTrack = aStep->GetTrack();
68 
69  //This detector is a black hole! Just put all kinetic energy into edep
70  G4double edep = aTrack->GetKineticEnergy() / GeV;
71  G4Track* killtrack = const_cast<G4Track *> (aTrack);
72  killtrack->SetTrackStatus(fStopAndKill);
73 
74  /* Make sure we are in a volume */
75  if ( detector_->IsActive() )
76  {
77  /* Check if particle is 'geantino' */
78  bool geantino = false;
79  if (aTrack->GetParticleDefinition()->GetPDGEncoding() == 0 && aTrack->GetParticleDefinition()->GetParticleName().find("geantino") != string::npos)
80  {
81  geantino = true;
82  }
83 
84  /* Get Geant4 pre- and post-step points */
85  G4StepPoint * prePoint = aStep->GetPreStepPoint();
86  G4StepPoint * postPoint = aStep->GetPostStepPoint();
87 
88  switch (prePoint->GetStepStatus())
89  {
90  case fGeomBoundary:
91  case fUndefined:
92  hit = new PHG4Hitv1();
93  // hit->set_layer(0);
94  hit->set_scint_id(tower_id);
95 
96  /* Set hit location (tower index) */
97  //hit->set_index_j(idx_j);
98  //hit->set_index_k(idx_k);
99  //hit->set_index_l(idx_l);
100  hit->set_index_j(0);
101  hit->set_index_k(0);
102  hit->set_index_l(0);
103 
104  /* Set hit location (space point) */
105  hit->set_x( 0, prePoint->GetPosition().x() / cm);
106  hit->set_y( 0, prePoint->GetPosition().y() / cm );
107  hit->set_z( 0, prePoint->GetPosition().z() / cm );
108 
109  /* Set momentum */
110  hit->set_x( 0, prePoint->GetMomentum().x() / GeV );
111  hit->set_y( 0, prePoint->GetMomentum().y() / GeV );
112  hit->set_z( 0, prePoint->GetMomentum().z() / GeV );
113 
114  /* Set hit time */
115  hit->set_t( 0, prePoint->GetGlobalTime() / nanosecond );
116 
117  /* set the track ID */
118  {
119  hit->set_trkid(aTrack->GetTrackID());
120  if ( G4VUserTrackInformation* p = aTrack->GetUserInformation() )
121  {
122  if ( PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p) )
123  {
124  hit->set_trkid(pp->GetUserTrackId());
125  hit->set_shower_id(pp->GetShower()->get_id());
126  }
127  }
128  }
129 
130  /* set intial energy deposit */
131  hit->set_edep( 0 );
132  hit->set_eion( 0 );
133 
134  hits_->AddHit(layer_id, hit);
135 
136  {
137  if ( G4VUserTrackInformation* p = aTrack->GetUserInformation() )
138  {
139  if ( PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p) )
140  {
141  pp->GetShower()->add_g4hit_id(hits_->GetID(),hit->get_hit_id());
142  }
143  }
144  }
145 
146  break;
147  default:
148  break;
149  }
150 
151  /* Update exit values- will be overwritten with every step until
152  * we leave the volume or the particle ceases to exist */
153  hit->set_x( 1, postPoint->GetPosition().x() / cm );
154  hit->set_y( 1, postPoint->GetPosition().y() / cm );
155  hit->set_z( 1, postPoint->GetPosition().z() / cm );
156 
157  hit->set_t( 1, postPoint->GetGlobalTime() / nanosecond );
158 
159  /* sum up the energy to get total deposited */
160  hit->set_edep(hit->get_edep() + edep);
161  hit->set_eion(hit->get_eion() + eion);
162 
163  if (geantino)
164  {
165  hit->set_edep(-1); // only energy=0 g4hits get dropped, this way geantinos survive the g4hit compression
166  hit->set_eion(-1);
167  }
168  if (edep > 0)
169  {
170  if ( G4VUserTrackInformation* p = aTrack->GetUserInformation() )
171  {
172  if ( PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p) )
173  {
174  pp->SetKeep(1); // we want to keep the track
175  }
176  }
177  }
178  return true;
179  }
180  else
181  {
182  return false;
183  }
184 }
185 
187 {
188  string hitnodename;
189 
190  if (detector_->SuperDetector() != "NONE")
191  {
192  hitnodename = "G4HIT_ENVELOPE_" + detector_->SuperDetector();
193  }
194  else
195  {
196  hitnodename = "G4HIT_ENVELOPE_" + detector_->GetName();
197  }
198 
199  hits_ = findNode::getClass<PHG4HitContainer>( topNode , hitnodename.c_str() );
200 
201  // if we do not find the node it's messed up.
202  if ( ! hits_ )
203  {
204  std::cout << "PHG4CrystalCalorimeterSteppingAction::SetTopNode - unable to find " << hitnodename << std::endl;
205  }
206 
207 
208 }