Geant4  10.03
LXeSteppingAction.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: LXeSteppingAction.cc 73915 2013-09-17 07:32:26Z gcosmo $
27 //
30 //
31 //
32 #include "LXeSteppingAction.hh"
33 #include "LXeEventAction.hh"
34 #include "LXeTrackingAction.hh"
35 #include "LXeTrajectory.hh"
36 #include "LXePMTSD.hh"
39 #include "LXeSteppingMessenger.hh"
40 #include "LXeRecorderBase.hh"
41 
42 #include "G4SteppingManager.hh"
43 #include "G4SDManager.hh"
44 #include "G4EventManager.hh"
45 #include "G4ProcessManager.hh"
46 #include "G4Track.hh"
47 #include "G4Step.hh"
48 #include "G4Event.hh"
49 #include "G4StepPoint.hh"
50 #include "G4TrackStatus.hh"
51 #include "G4VPhysicalVolume.hh"
52 #include "G4ParticleDefinition.hh"
53 #include "G4ParticleTypes.hh"
54 
55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56 
58  : fRecorder(r),fOneStepPrimaries(false)
59 {
61 
63 }
64 
65 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
66 
68 
69 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
70 
72 
73  G4Track* theTrack = theStep->GetTrack();
74 
75  if ( theTrack->GetCurrentStepNumber() == 1 ) fExpectedNextStatus = Undefined;
76 
77  LXeUserTrackInformation* trackInformation
79  LXeUserEventInformation* eventInformation
82 
83  G4StepPoint* thePrePoint = theStep->GetPreStepPoint();
84  G4VPhysicalVolume* thePrePV = thePrePoint->GetPhysicalVolume();
85 
86  G4StepPoint* thePostPoint = theStep->GetPostStepPoint();
87  G4VPhysicalVolume* thePostPV = thePostPoint->GetPhysicalVolume();
88 
89  G4OpBoundaryProcessStatus boundaryStatus=Undefined;
90  static G4ThreadLocal G4OpBoundaryProcess* boundary=NULL;
91 
92  //find the boundary process only once
93  if(!boundary){
95  = theStep->GetTrack()->GetDefinition()->GetProcessManager();
96  G4int nprocesses = pm->GetProcessListLength();
97  G4ProcessVector* pv = pm->GetProcessList();
98  G4int i;
99  for( i=0;i<nprocesses;i++){
100  if((*pv)[i]->GetProcessName()=="OpBoundary"){
101  boundary = (G4OpBoundaryProcess*)(*pv)[i];
102  break;
103  }
104  }
105  }
106 
107  if(theTrack->GetParentID()==0){
108  //This is a primary track
109 
114 
115  //If we havent already found the conversion position and there were
116  //secondaries generated, then search for it
117  if(!eventInformation->IsConvPosSet() && tN2ndariesTot>0 ){
118  for(size_t lp1=(*fSecondary).size()-tN2ndariesTot;
119  lp1<(*fSecondary).size(); lp1++){
120  const G4VProcess* creator=(*fSecondary)[lp1]->GetCreatorProcess();
121  if(creator){
122  G4String creatorName=creator->GetProcessName();
123  if(creatorName=="phot"||creatorName=="compt"||creatorName=="conv"){
124  //since this is happening before the secondary is being tracked
125  //the Vertex position has not been set yet(set in initial step)
126  eventInformation->SetConvPos((*fSecondary)[lp1]->GetPosition());
127  }
128  }
129  }
130  }
131 
132  if(fOneStepPrimaries&&thePrePV->GetName()=="scintillator")
133  theTrack->SetTrackStatus(fStopAndKill);
134  }
135 
136  if(!thePostPV){//out of world
138  return;
139  }
140 
141  G4ParticleDefinition* particleType = theTrack->GetDefinition();
142  if(particleType==G4OpticalPhoton::OpticalPhotonDefinition()){
143  //Optical photon only
144 
145  if(thePrePV->GetName()=="Slab")
146  //force drawing of photons in WLS slab
147  trackInformation->SetForceDrawTrajectory(true);
148  else if(thePostPV->GetName()=="expHall")
149  //Kill photons entering expHall from something other than Slab
150  theTrack->SetTrackStatus(fStopAndKill);
151 
152  //Was the photon absorbed by the absorption process
153  if(thePostPoint->GetProcessDefinedStep()->GetProcessName()
154  =="OpAbsorption"){
155  eventInformation->IncAbsorption();
156  trackInformation->AddTrackStatusFlag(absorbed);
157  }
158 
159  boundaryStatus=boundary->GetStatus();
160 
161  //Check to see if the partcile was actually at a boundary
162  //Otherwise the boundary status may not be valid
163  //Prior to Geant4.6.0-p1 this would not have been enough to check
164  if(thePostPoint->GetStepStatus()==fGeomBoundary){
166  if(boundaryStatus!=StepTooSmall){
168  ed << "LXeSteppingAction::UserSteppingAction(): "
169  << "No reallocation step after reflection!"
170  << G4endl;
171  G4Exception("LXeSteppingAction::UserSteppingAction()", "LXeExpl01",
172  FatalException,ed,
173  "Something is wrong with the surface normal or geometry");
174  }
175  }
177  switch(boundaryStatus){
178  case Absorption:
179  trackInformation->AddTrackStatusFlag(boundaryAbsorbed);
180  eventInformation->IncBoundaryAbsorption();
181  break;
182  case Detection: //Note, this assumes that the volume causing detection
183  //is the photocathode because it is the only one with
184  //non-zero efficiency
185  {
186  //Triger sensitive detector manually since photon is
187  //absorbed but status was Detection
189  G4String sdName="/LXeDet/pmtSD";
190  LXePMTSD* pmtSD = (LXePMTSD*)SDman->FindSensitiveDetector(sdName);
191  if(pmtSD)pmtSD->ProcessHits_constStep(theStep,NULL);
192  trackInformation->AddTrackStatusFlag(hitPMT);
193  break;
194  }
195  case FresnelReflection:
198  case LobeReflection:
199  case SpikeReflection:
200  case BackScattering:
201  trackInformation->IncReflections();
203  break;
204  default:
205  break;
206  }
207  if(thePostPV->GetName()=="sphere")
208  trackInformation->AddTrackStatusFlag(hitSphere);
209  }
210  }
211 
212  if(fRecorder)fRecorder->RecordStep(theStep);
213 }
LXeSteppingMessenger * fSteppingMessenger
void SetTrackStatus(const G4TrackStatus aTrackStatus)
G4ParticleDefinition * GetDefinition() const
G4int GetParentID() const
virtual void UserSteppingAction(const G4Step *)
G4VUserEventInformation * GetUserInformation() const
Definition: G4Event.hh:199
Definition of the LXeRecorderBase class.
Definition of the LXeSteppingAction class.
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4StepStatus GetStepStatus() const
Definition of the LXeTrackingAction class.
Definition of the LXeUserEventInformation class.
G4OpBoundaryProcessStatus
void SetConvPos(const G4ThreeVector &p)
G4bool ProcessHits_constStep(const G4Step *, G4TouchableHistory *)
Definition: LXePMTSD.cc:96
#define G4ThreadLocal
Definition: tls.hh:89
int G4int
Definition: G4Types.hh:78
G4VUserTrackInformation * GetUserInformation() const
G4StepPoint * GetPreStepPoint() const
Definition of the LXePMTSD class.
G4int GetCurrentStepNumber() const
G4VPhysicalVolume * GetPhysicalVolume() const
const G4String & GetName() const
G4SteppingManager * fpSteppingManager
Definition of the LXeEventAction class.
Definition of the LXeTrajectory class.
Definition of the LXeUserTrackInformation class.
Definition: G4Step.hh:76
G4VSensitiveDetector * FindSensitiveDetector(G4String dName, G4bool warning=true)
Definition: G4SDManager.cc:128
G4int GetfN2ndariesAlongStepDoIt()
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4OpBoundaryProcessStatus fExpectedNextStatus
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
LXeSteppingAction(LXeRecorderBase *)
const G4VProcess * GetProcessDefinedStep() const
std::vector< G4Track * > G4TrackVector
G4ProcessManager * GetProcessManager() const
virtual void RecordStep(const G4Step *)
G4StepPoint * GetPostStepPoint() const
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40
static G4EventManager * GetEventManager()
#define G4endl
Definition: G4ios.hh:61
Definition of the LXeSteppingMessenger class.
static G4OpticalPhoton * OpticalPhotonDefinition()
const G4Event * GetConstCurrentEvent()
G4Track * GetTrack() const
G4int GetProcessListLength() const
virtual ~LXeSteppingAction()
G4TrackVector * GetfSecondary()
G4ProcessVector * GetProcessList() const
LXeRecorderBase * fRecorder