ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4CEmcTestBeamSteppingAction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4CEmcTestBeamSteppingAction.cc
2 
4 
6 #include <g4main/PHG4Hit.h>
7 #include <g4main/PHG4Hitv1.h>
8 #include <g4main/PHG4Shower.h>
9 #include <g4main/PHG4SteppingAction.h> // for PHG4SteppingAction
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>
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 operator+, string, ope...
31 
32 class G4VPhysicalVolume;
33 class PHCompositeNode;
34 
35 using namespace std;
36 //____________________________________________________________________________..
38  PHG4SteppingAction(detector->GetName()),
39  detector_( detector ),
40  hits_(nullptr),
41  absorberhits_(nullptr),
42  hit(nullptr)
43 {}
44 
45 //____________________________________________________________________________..
47 {
48 
50  // get volume of the current step
52 
53  int whichactive = detector_->IsInCEmcTestBeam(volume);
54 
55  if (!whichactive)
56  {
57  return false;
58  }
59 
60  // collect energy and track length step by step
63 
64  const G4Track* aTrack = aStep->GetTrack();
65 
66  // if this block stops everything, just put all kinetic energy into edep
67  if (detector_->IsBlackHole())
68  {
69  edep = aTrack->GetKineticEnergy() / GeV;
70  G4Track* killtrack = const_cast<G4Track *> (aTrack);
71  killtrack->SetTrackStatus(fStopAndKill);
72  }
73 
74  int tower_id = touch->GetCopyNumber(2); // tower number
75  // make sure we are in a volume
76  if ( detector_->IsActive() )
77  {
78  bool geantino = false;
79  // the check for the pdg code speeds things up, I do not want to make
80  // an expensive string compare for every track when we know
81  // geantino or chargedgeantino has pid=0
82  if (aTrack->GetParticleDefinition()->GetPDGEncoding() == 0 &&
83  aTrack->GetParticleDefinition()->GetParticleName().find("geantino") != string::npos)
84  {
85  geantino = true;
86  }
87  G4StepPoint * prePoint = aStep->GetPreStepPoint();
88  G4StepPoint * postPoint = aStep->GetPostStepPoint();
89  // cout << "track id " << aTrack->GetTrackID() << endl;
90  // cout << "time prepoint: " << prePoint->GetGlobalTime() << endl;
91  // cout << "time postpoint: " << postPoint->GetGlobalTime() << endl;
92  switch (prePoint->GetStepStatus())
93  {
94  case fGeomBoundary:
95  case fUndefined:
96  hit = new PHG4Hitv1();
97  hit->set_layer((unsigned int)tower_id);
98  hit->set_scint_id(touch->GetCopyNumber(1)); // the copy number of the sandwich
99  //here we set the entrance values in cm
100  hit->set_x( 0, prePoint->GetPosition().x() / cm);
101  hit->set_y( 0, prePoint->GetPosition().y() / cm );
102  hit->set_z( 0, prePoint->GetPosition().z() / cm );
103  // time in ns
104  hit->set_t( 0, prePoint->GetGlobalTime() / nanosecond );
105  //set the track ID
106  {
107  hit->set_trkid(aTrack->GetTrackID());
108  if ( G4VUserTrackInformation* p = aTrack->GetUserInformation() )
109  {
110  if ( PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p) )
111  {
112  hit->set_trkid(pp->GetUserTrackId());
113  hit->set_shower_id(pp->GetShower()->get_id());
114  }
115  }
116  }
117 
118  //set the initial energy deposit
119  hit->set_edep(0);
120  hit->set_eion(0); // only implemented for v5 otherwise empty
121  PHG4HitContainer *hitcontainer;
122  // here we do things which are different between scintillator and absorber hits
123  if (whichactive > 0) // return of isinCEmcTestDetector, > 0 hit in scintillator, < 0 hit in absorber
124  {
125  hitcontainer = hits_;
126  }
127  else
128  {
129  hitcontainer = absorberhits_;
130  }
131  // here we set what is common for scintillator and absorber hits
132  hitcontainer->AddHit(tower_id, hit);
133  if ( G4VUserTrackInformation* p = aTrack->GetUserInformation() )
134  {
135  if ( PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p) )
136  {
137  pp->GetShower()->add_g4hit_id(hitcontainer->GetID(),hit->get_hit_id());
138  }
139  }
140  break;
141  default:
142  break;
143  }
144  // here we just update the exit values, it will be overwritten
145  // for every step until we leave the volume or the particle
146  // ceases to exist
147  hit->set_x( 1, postPoint->GetPosition().x() / cm );
148  hit->set_y( 1, postPoint->GetPosition().y() / cm );
149  hit->set_z( 1, postPoint->GetPosition().z() / cm );
150 
151  hit->set_t( 1, postPoint->GetGlobalTime() / nanosecond );
152  //sum up the energy to get total deposited
153  hit->set_edep(hit->get_edep() + edep);
154  hit->set_eion(hit->get_eion() + eion);
155  if (geantino)
156  {
157  hit->set_edep(-1); // only energy=0 g4hits get dropped, this way geantinos survive the g4hit compression
158  hit->set_eion(-1);
159  }
160  if (edep > 0)
161  {
162  if ( G4VUserTrackInformation* p = aTrack->GetUserInformation() )
163  {
164  if ( PHG4TrackUserInfoV1* pp = dynamic_cast<PHG4TrackUserInfoV1*>(p) )
165  {
166  pp->SetKeep(1); // we want to keep the track
167  }
168 
169 
170  }
171  }
172 
173  // hit->identify();
174  // return true to indicate the hit was used
175  return true;
176 
177  }
178  else
179  {
180  return false;
181  }
182 }
183 
184 //____________________________________________________________________________..
186 {
187 
188  string hitnodename;
189  string absorbernodename;
190  if (detector_->SuperDetector() != "NONE")
191  {
192  hitnodename = "G4HIT_" + detector_->SuperDetector();
193  absorbernodename = "G4HIT_ABSORBER_" + detector_->SuperDetector();
194  }
195  else
196  {
197  hitnodename = "G4HIT_" + detector_->GetName();
198  absorbernodename = "G4HIT_ABSORBER_" + detector_->GetName();
199  }
200 
201  //now look for the map and grab a pointer to it.
202  hits_ = findNode::getClass<PHG4HitContainer>( topNode , hitnodename.c_str() );
203  absorberhits_ = findNode::getClass<PHG4HitContainer>( topNode , absorbernodename.c_str() );
204 
205  // if we do not find the node it's messed up.
206  if ( ! hits_ )
207  {
208  std::cout << "PHG4CEmcTestBeamSteppingAction::SetTopNode - unable to find " << hitnodename << std::endl;
209  }
210  if ( ! absorberhits_)
211  {
212  if (Verbosity() > 0)
213  {
214  cout << "PHG4HcalSteppingAction::SetTopNode - unable to find " << absorbernodename << endl;
215  }
216  }
217 }