Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
WLSSteppingAction.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: WLSSteppingAction.cc 105203 2017-07-14 12:23:01Z gcosmo $
27 //
30 //
31 //
32 #include "G4Run.hh"
33 #include "G4Step.hh"
34 #include "G4Track.hh"
35 #include "G4StepPoint.hh"
36 #include "G4TrackStatus.hh"
37 #include "G4VPhysicalVolume.hh"
38 #include "G4ParticleDefinition.hh"
39 
40 #include "WLSSteppingAction.hh"
43 #include "WLSPhotonDetSD.hh"
44 
45 #include "G4ParticleTypes.hh"
46 
48 
49 #include "G4ProcessManager.hh"
50 #include "G4OpBoundaryProcess.hh"
51 
52 #include "G4RunManager.hh"
53 #include "G4SDManager.hh"
54 #include "G4UImanager.hh"
55 
56 #include "G4ThreeVector.hh"
57 #include "G4ios.hh"
58 #include "G4SystemOfUnits.hh"
59 #include <sstream>
60 
61 // Purpose: Save relevant information into User Track Information
62 
63 static const G4ThreeVector ZHat = G4ThreeVector(0.0,0.0,1.0);
64 
65 G4int WLSSteppingAction::fMaxRndmSave = 10000;
66 
67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
68 
70  : fDetector(detector)
71 {
72  fSteppingMessenger = new WLSSteppingActionMessenger(this);
73 
74  fCounterEnd = 0;
75  fCounterMid = 0;
76  fBounceLimit = 100000;
77 
78  fOpProcess = NULL;
79 
80  ResetCounters();
81 }
82 
83 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
84 
86 {
87  delete fSteppingMessenger;
88 }
89 
90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
91 
92 void WLSSteppingAction::SetBounceLimit(G4int i) {fBounceLimit = i;}
93 
94 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
95 
96 G4int WLSSteppingAction::GetNumberOfBounces() {return fCounterBounce;}
97 
98 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
99 
100 G4int WLSSteppingAction::GetNumberOfClad1Bounces() {return fCounterClad1Bounce;}
101 
102 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
103 
104 G4int WLSSteppingAction::GetNumberOfClad2Bounces() {return fCounterClad2Bounce;}
105 
106 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
107 
108 G4int WLSSteppingAction::GetNumberOfWLSBounces() {return fCounterWLSBounce;}
109 
110 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
111 
113  G4int temp = fCounterEnd; fCounterEnd = 0; return temp;
114 }
115 
116 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
117 
118 inline void WLSSteppingAction::SaveRandomStatus(G4String subDir)
119 // save the random status into a sub-directory
120 // Pre: subDir must be empty or ended with "/"
121 {
122 
123  // don't save if the maximum amount has been reached
124  if (WLSSteppingAction::fMaxRndmSave == 0) return;
125 
126  G4RunManager* theRunManager = G4RunManager::GetRunManager();
127  G4String randomNumberStatusDir = theRunManager->GetRandomNumberStoreDir();
128 
129  G4String fileIn = randomNumberStatusDir + "currentEvent.rndm";
130 
131  std::ostringstream os;
132 
133  os << "run" << theRunManager->GetCurrentRun()->GetRunID() << "evt"
134  << theRunManager->GetCurrentEvent()->GetEventID() << ".rndm" << '\0';
135 
136  G4String fileOut = randomNumberStatusDir + subDir + os.str();
137 
138  G4String copCmd = "/control/shell cp "+fileIn+" "+fileOut;
140 
141  WLSSteppingAction::fMaxRndmSave--;
142 }
143 
144 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
145 
147 {
148  G4Track* theTrack = theStep->GetTrack();
149  WLSUserTrackInformation* trackInformation
151 
152  G4StepPoint* thePrePoint = theStep->GetPreStepPoint();
153  G4StepPoint* thePostPoint = theStep->GetPostStepPoint();
154 
155  G4VPhysicalVolume* thePrePV = thePrePoint->GetPhysicalVolume();
156  G4VPhysicalVolume* thePostPV = thePostPoint->GetPhysicalVolume();
157 
158  G4String thePrePVname = " ";
159  G4String thePostPVname = " ";
160 
161  if (thePostPV) {
162  thePrePVname = thePrePV->GetName();
163  thePostPVname = thePostPV->GetName();
164  }
165 
166  //Recording data for start
167  if (theTrack->GetParentID()==0) {
168  //This is a primary track
169  if ( theTrack->GetCurrentStepNumber() == 1 ) {
170 // G4double x = theTrack->GetVertexPosition().x();
171 // G4double y = theTrack->GetVertexPosition().y();
172 // G4double z = theTrack->GetVertexPosition().z();
173 // G4double pz = theTrack->GetVertexMomentumDirection().z();
174 // G4double fInitTheta =
175 // theTrack->GetVertexMomentumDirection().angle(ZHat);
176  }
177  }
178 
179  // Retrieve the status of the photon
181 
182  G4ProcessManager* OpManager =
184 
185  if (OpManager) {
186  G4int MAXofPostStepLoops =
187  OpManager->GetPostStepProcessVector()->entries();
188  G4ProcessVector* fPostStepDoItVector =
190 
191  for ( G4int i=0; i<MAXofPostStepLoops; i++) {
192  G4VProcess* fCurrentProcess = (*fPostStepDoItVector)[i];
193  fOpProcess = dynamic_cast<G4OpBoundaryProcess*>(fCurrentProcess);
194  if (fOpProcess) { theStatus = fOpProcess->GetStatus(); break;}
195  }
196  }
197 
198  // Find the skewness of the ray at first change of boundary
199  if ( fInitGamma == -1 &&
200  (theStatus == TotalInternalReflection
201  || theStatus == FresnelReflection
202  || theStatus == FresnelRefraction)
203  && trackInformation->IsStatus(InsideOfFiber) ) {
204 
205  G4double px = theTrack->GetVertexMomentumDirection().x();
206  G4double py = theTrack->GetVertexMomentumDirection().y();
207  G4double x = theTrack->GetPosition().x();
208  G4double y = theTrack->GetPosition().y();
209 
210  fInitGamma = x * px + y * py;
211 
212  fInitGamma =
213  fInitGamma / std::sqrt(px*px + py*py) / std::sqrt(x*x + y*y);
214 
215  fInitGamma = std::acos(fInitGamma*rad);
216 
217  if ( fInitGamma / deg > 90.0) { fInitGamma = 180 * deg - fInitGamma;}
218  }
219  // Record Photons that missed the photon detector but escaped from readout
220  if ( !thePostPV && trackInformation->IsStatus(EscapedFromReadOut) ) {
221 // UpdateHistogramSuccess(thePostPoint,theTrack);
222  ResetCounters();
223 
224  return;
225  }
226 
227  // Assumed photons are originated at the fiber OR
228  // the fiber is the first material the photon hits
229  switch (theStatus) {
230 
231  // Exiting the fiber
232  case FresnelRefraction:
233  case SameMaterial:
234 
235  G4bool isFiber;
236  isFiber = thePostPVname == "WLSFiber"
237  || thePostPVname == "Clad1"
238  || thePostPVname == "Clad2";
239 
240  if ( isFiber ) {
241 
242  if (trackInformation->IsStatus(OutsideOfFiber))
243  trackInformation->AddStatusFlag(InsideOfFiber);
244 
245  // Set the Exit flag when the photon refracted out of the fiber
246  } else if (trackInformation->IsStatus(InsideOfFiber)) {
247 
248  // EscapedFromReadOut if the z position is the same as fiber's end
249  if (theTrack->GetPosition().z() == fDetector->GetWLSFiberEnd())
250  {
251  trackInformation->AddStatusFlag(EscapedFromReadOut);
252  fCounterEnd++;
253  }
254  else // Escaped from side
255  {
256  trackInformation->AddStatusFlag(EscapedFromSide);
257  trackInformation->SetExitPosition(theTrack->GetPosition());
258 
259 // UpdateHistogramEscape(thePostPoint,theTrack);
260 
261  fCounterMid++;
262  ResetCounters();
263  }
264 
265  trackInformation->AddStatusFlag(OutsideOfFiber);
266  trackInformation->SetExitPosition(theTrack->GetPosition());
267 
268  }
269 
270  return;
271 
272  // Internal Reflections
274 
275  // Kill the track if it's number of bounces exceeded the limit
276  if (fBounceLimit > 0 && fCounterBounce >= fBounceLimit)
277  {
278  theTrack->SetTrackStatus(fStopAndKill);
279  trackInformation->AddStatusFlag(murderee);
280  ResetCounters();
281  G4cout << "\n Bounce Limit Exceeded" << G4endl;
282  return;
283  }
284  break;
285 
286  case FresnelReflection:
287 
288  fCounterBounce++;
289 
290  if ( thePrePVname == "WLSFiber") fCounterWLSBounce++;
291 
292  else if ( thePrePVname == "Clad1") fCounterClad1Bounce++;
293 
294  else if ( thePrePVname == "Clad2") fCounterClad2Bounce++;
295 
296  // Determine if the photon has reflected off the read-out end
297  if (theTrack->GetPosition().z() == fDetector->GetWLSFiberEnd())
298  {
299  if (!trackInformation->IsStatus(ReflectedAtReadOut) &&
300  trackInformation->IsStatus(InsideOfFiber))
301  {
302  trackInformation->AddStatusFlag(ReflectedAtReadOut);
303 
304  if (fDetector->IsPerfectFiber() &&
305  theStatus == TotalInternalReflection)
306  {
307  theTrack->SetTrackStatus(fStopAndKill);
308  trackInformation->AddStatusFlag(murderee);
309 // UpdateHistogramReflect(thePostPoint,theTrack);
310  ResetCounters();
311  return;
312  }
313  }
314  }
315  return;
316 
317  // Reflection of the mirror
319  case LobeReflection:
320  case SpikeReflection:
321 
322  // Check if it hits the mirror
323  if ( thePostPVname == "Mirror" )
324  trackInformation->AddStatusFlag(ReflectedAtMirror);
325 
326  return;
327 
328  // Detected by a detector
329  case Detection:
330 
331  // Check if the photon hits the detector and process the hit if it does
332  if ( thePostPVname == "PhotonDet" ) {
333 
335  G4String SDname="WLS/PhotonDet";
336  WLSPhotonDetSD* mppcSD =
337  (WLSPhotonDetSD*)SDman->FindSensitiveDetector(SDname);
338 
339  if (mppcSD) mppcSD->ProcessHits_constStep(theStep,NULL);
340 
341  // Record Photons that escaped at the end
342 // if (trackInformation->IsStatus(EscapedFromReadOut))
343 // UpdateHistogramSuccess(thePostPoint,theTrack);
344 
345  // Stop Tracking when it hits the detector's surface
346  ResetCounters();
347  theTrack->SetTrackStatus(fStopAndKill);
348 
349  return;
350  }
351  break;
352 
353  default: break;
354 
355  }
356 
357  // Check for absorbed photons
358  if (theTrack->GetTrackStatus() != fAlive &&
359  trackInformation->IsStatus(InsideOfFiber))
360  {
361 // UpdateHistogramAbsorb(thePostPoint,theTrack);
362  ResetCounters();
363  return;
364  }
365 
366 }
void SetTrackStatus(const G4TrackStatus aTrackStatus)
G4int GetParentID() const
G4bool AddStatusFlag(TrackStatus s)
CLHEP::Hep3Vector G4ThreeVector
double x() const
const G4String & GetRandomNumberStoreDir() const
G4bool ProcessHits_constStep(const G4Step *, G4TouchableHistory *)
Definition of the WLSUserTrackInformation class.
void SetExitPosition(const G4ThreeVector &pos)
const G4ThreeVector & GetPosition() const
G4TrackStatus GetTrackStatus() const
G4bool IsStatus(TrackStatus s)
G4OpBoundaryProcessStatus
tuple x
Definition: test.py:50
static constexpr double rad
Definition: G4SIunits.hh:149
G4OpBoundaryProcessStatus GetStatus() const
int G4int
Definition: G4Types.hh:78
const G4Run * GetCurrentRun() const
Definition of the WLSDetectorConstruction class.
double z() const
static G4UImanager * GetUIpointer()
Definition: G4UImanager.cc:59
G4VUserTrackInformation * GetUserInformation() const
G4int GetEventID() const
Definition: G4Event.hh:151
G4StepPoint * GetPreStepPoint() const
virtual void UserSteppingAction(const G4Step *)
G4int entries() const
G4GLOB_DLL std::ostream G4cout
G4int GetCurrentStepNumber() const
G4VPhysicalVolume * GetPhysicalVolume() const
const G4String & GetName() const
bool G4bool
Definition: G4Types.hh:79
G4int GetRunID() const
Definition: G4Run.hh:76
virtual ~WLSSteppingAction()
Definition: G4Step.hh:76
G4VSensitiveDetector * FindSensitiveDetector(G4String dName, G4bool warning=true)
Definition: G4SDManager.cc:128
WLSSteppingAction(WLSDetectorConstruction *)
static G4RunManager * GetRunManager()
Definition: G4RunManager.cc:79
static G4OpticalPhoton * OpticalPhoton()
G4ProcessManager * GetProcessManager() const
G4StepPoint * GetPostStepPoint() const
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40
double y() const
Definition of the WLSPhotonDetSD class.
#define G4endl
Definition: G4ios.hh:61
static const G4ThreeVector ZHat
Definition of the WLSSteppingAction class.
void SetBounceLimit(G4int)
double G4double
Definition: G4Types.hh:76
static constexpr double deg
Definition: G4SIunits.hh:152
const G4Event * GetCurrentEvent() const
G4Track * GetTrack() const
const G4ThreeVector & GetVertexMomentumDirection() const
Definition of the WLSSteppingActionMessenger class.
G4int ApplyCommand(const char *aCommand)
Definition: G4UImanager.cc:447
G4ProcessVector * GetPostStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const