Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
LXeSteppingAction.cc
이 파일의 문서화 페이지로 가기
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 110138 2018-05-16 07:31:43Z gcosmo $
27 //
30 //
31 //
32 #include "LXeSteppingAction.hh"
33 #include "LXeEventAction.hh"
34 #include "LXeTrackingAction.hh"
35 #include "LXeTrajectory.hh"
36 #include "LXePMTSD.hh"
38 #include "LXeSteppingMessenger.hh"
39 
40 #include "G4SteppingManager.hh"
41 #include "G4SDManager.hh"
42 #include "G4EventManager.hh"
43 #include "G4ProcessManager.hh"
44 #include "G4Track.hh"
45 #include "G4Step.hh"
46 #include "G4Event.hh"
47 #include "G4StepPoint.hh"
48 #include "G4TrackStatus.hh"
49 #include "G4VPhysicalVolume.hh"
50 #include "G4ParticleDefinition.hh"
51 #include "G4ParticleTypes.hh"
52 
53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
54 
56  : fOneStepPrimaries(false),
57  fEventAction(ea)
58 {
60 
62 }
63 
64 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
65 
67 
68 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
69 
71 
72  G4Track* theTrack = theStep->GetTrack();
73 
74  if ( theTrack->GetCurrentStepNumber() == 1 ) fExpectedNextStatus = Undefined;
75 
76  LXeUserTrackInformation* trackInformation
78 
79  G4StepPoint* thePrePoint = theStep->GetPreStepPoint();
80  G4VPhysicalVolume* thePrePV = thePrePoint->GetPhysicalVolume();
81 
82  G4StepPoint* thePostPoint = theStep->GetPostStepPoint();
83  G4VPhysicalVolume* thePostPV = thePostPoint->GetPhysicalVolume();
84 
85  G4OpBoundaryProcessStatus boundaryStatus=Undefined;
86  static G4ThreadLocal G4OpBoundaryProcess* boundary = nullptr;
87 
88  //find the boundary process only once
89  if(!boundary){
91  = theStep->GetTrack()->GetDefinition()->GetProcessManager();
92  G4int nprocesses = pm->GetProcessListLength();
93  G4ProcessVector* pv = pm->GetProcessList();
94  G4int i;
95  for( i=0;i<nprocesses;i++){
96  if((*pv)[i]->GetProcessName()=="OpBoundary"){
97  boundary = (G4OpBoundaryProcess*)(*pv)[i];
98  break;
99  }
100  }
101  }
102 
103  if(theTrack->GetParentID()==0){
104  //This is a primary track
105 
110 
111  //If we havent already found the conversion position and there were
112  //secondaries generated, then search for it
113  if(!fEventAction->IsConvPosSet() && tN2ndariesTot>0 ){
114  for(size_t lp1=(*fSecondary).size()-tN2ndariesTot;
115  lp1<(*fSecondary).size(); lp1++){
116  const G4VProcess* creator=(*fSecondary)[lp1]->GetCreatorProcess();
117  if(creator){
118  G4String creatorName=creator->GetProcessName();
119  if(creatorName=="phot"||creatorName=="compt"||creatorName=="conv"){
120  //since this is happening before the secondary is being tracked
121  //the Vertex position has not been set yet(set in initial step)
122  fEventAction->SetConvPos((*fSecondary)[lp1]->GetPosition());
123  }
124  }
125  }
126  }
127 
128  if(fOneStepPrimaries&&thePrePV->GetName()=="scintillator")
129  theTrack->SetTrackStatus(fStopAndKill);
130  }
131 
132  if(!thePostPV){//out of world
134  return;
135  }
136 
137  G4ParticleDefinition* particleType = theTrack->GetDefinition();
138  if(particleType==G4OpticalPhoton::OpticalPhotonDefinition()){
139  //Optical photon only
140 
141  if(thePrePV->GetName()=="Slab")
142  //force drawing of photons in WLS slab
143  trackInformation->SetForceDrawTrajectory(true);
144  else if(thePostPV->GetName()=="expHall")
145  //Kill photons entering expHall from something other than Slab
146  theTrack->SetTrackStatus(fStopAndKill);
147 
148  //Was the photon absorbed by the absorption process
149  if(thePostPoint->GetProcessDefinedStep()->GetProcessName()
150  =="OpAbsorption"){
152  trackInformation->AddTrackStatusFlag(absorbed);
153  }
154 
155  boundaryStatus=boundary->GetStatus();
156 
157  //Check to see if the partcile was actually at a boundary
158  //Otherwise the boundary status may not be valid
159  //Prior to Geant4.6.0-p1 this would not have been enough to check
160  if(thePostPoint->GetStepStatus()==fGeomBoundary){
162  if(boundaryStatus!=StepTooSmall){
164  ed << "LXeSteppingAction::UserSteppingAction(): "
165  << "No reallocation step after reflection!"
166  << G4endl;
167  G4Exception("LXeSteppingAction::UserSteppingAction()", "LXeExpl01",
168  FatalException,ed,
169  "Something is wrong with the surface normal or geometry");
170  }
171  }
173  switch(boundaryStatus){
174  case Absorption:
175  trackInformation->AddTrackStatusFlag(boundaryAbsorbed);
177  break;
178  case Detection: //Note, this assumes that the volume causing detection
179  //is the photocathode because it is the only one with
180  //non-zero efficiency
181  {
182  //Triger sensitive detector manually since photon is
183  //absorbed but status was Detection
185  G4String sdName="/LXeDet/pmtSD";
186  LXePMTSD* pmtSD = (LXePMTSD*)SDman->FindSensitiveDetector(sdName);
187  if(pmtSD)pmtSD->ProcessHits_constStep(theStep, nullptr);
188  trackInformation->AddTrackStatusFlag(hitPMT);
189  break;
190  }
191  case FresnelReflection:
194  case LobeReflection:
195  case SpikeReflection:
196  case BackScattering:
197  trackInformation->IncReflections();
199  break;
200  default:
201  break;
202  }
203  if(thePostPV->GetName()=="sphere")
204  trackInformation->AddTrackStatusFlag(hitSphere);
205  }
206  }
207 }
G4OpBoundaryProcessStatus
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:45
Definition of the LXeTrackingAction class.
G4VSensitiveDetector * FindSensitiveDetector(G4String dName, G4bool warning=true)
Definition: G4SDManager.cc:128
G4int GetProcessListLength() const
G4StepPoint * GetPreStepPoint() const
#define G4endl
Definition: G4ios.hh:61
G4double IsConvPosSet()
LXeSteppingMessenger * fSteppingMessenger
G4VPhysicalVolume * GetPhysicalVolume() const
G4int GetCurrentStepNumber() const
G4ProcessVector * GetProcessList() const
#define G4ThreadLocal
Definition: tls.hh:69
G4ParticleDefinition * GetDefinition() const
Definition of the LXeEventAction class.
G4StepStatus GetStepStatus() const
LXeSteppingAction(LXeEventAction *)
G4bool ProcessHits_constStep(const G4Step *, G4TouchableHistory *)
Definition: LXePMTSD.cc:96
Definition of the LXeUserTrackInformation class.
G4Track * GetTrack() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:411
static G4OpticalPhoton * OpticalPhotonDefinition()
Definition of the LXeSteppingMessenger class.
Definition: G4Step.hh:76
G4StepPoint * GetPostStepPoint() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40
void SetConvPos(const G4ThreeVector &p)
int G4int
Definition: G4Types.hh:78
G4SteppingManager * fpSteppingManager
G4ProcessManager * GetProcessManager() const
std::vector< G4Track * > G4TrackVector
virtual ~LXeSteppingAction()
Definition of the LXePMTSD class.
virtual void UserSteppingAction(const G4Step *)
G4VUserTrackInformation * GetUserInformation() const
void SetTrackStatus(const G4TrackStatus aTrackStatus)
G4int GetfN2ndariesAlongStepDoIt()
const G4VProcess * GetProcessDefinedStep() const
Definition of the LXeSteppingAction class.
Definition of the LXeTrajectory class.
void IncBoundaryAbsorption()
const G4String & GetName() const
G4TrackVector * GetfSecondary()
G4int GetParentID() const
LXeEventAction * fEventAction
G4OpBoundaryProcessStatus fExpectedNextStatus