Geant4  10.01
G4OpBoundaryProcess.hh
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 //
27 // $Id: G4OpBoundaryProcess.hh 84717 2014-10-20 07:39:47Z gcosmo $
28 //
29 //
31 // Optical Photon Boundary Process Class Definition
33 //
34 // File: G4OpBoundaryProcess.hh
35 // Description: Discrete Process -- reflection/refraction at
36 // optical interfaces
37 // Version: 1.1
38 // Created: 1997-06-18
39 // Modified: 2005-07-28 add G4ProcessType to constructor
40 // 1999-10-29 add method and class descriptors
41 // 1999-10-10 - Fill NewMomentum/NewPolarization in
42 // DoAbsorption. These members need to be
43 // filled since DoIt calls
44 // aParticleChange.SetMomentumChange etc.
45 // upon return (thanks to: Clark McGrew)
46 // 2006-11-04 - add capability of calculating the reflectivity
47 // off a metal surface by way of a complex index
48 // of refraction - Thanks to Sehwook Lee and John
49 // Hauptman (Dept. of Physics - Iowa State Univ.)
50 // 2009-11-10 - add capability of simulating surface reflections
51 // with Look-Up-Tables (LUT) containing measured
52 // optical reflectance for a variety of surface
53 // treatments - Thanks to Martin Janecek and
54 // William Moses (Lawrence Berkeley National Lab.)
55 // 2013-06-01 - add the capability of simulating the transmission
56 // of a dichronic filter
57 //
58 // Author: Peter Gumplinger
59 // adopted from work by Werner Keil - April 2/96
60 // mail: gum@triumf.ca
61 //
63 
64 #ifndef G4OpBoundaryProcess_h
65 #define G4OpBoundaryProcess_h 1
66 
68 // Includes
70 
71 #include "globals.hh"
72 #include "templates.hh"
73 #include "geomdefs.hh"
74 #include "Randomize.hh"
75 
76 #include "G4RandomTools.hh"
77 #include "G4RandomDirection.hh"
78 
79 #include "G4Step.hh"
80 #include "G4VDiscreteProcess.hh"
81 #include "G4DynamicParticle.hh"
82 #include "G4Material.hh"
84 #include "G4LogicalSkinSurface.hh"
85 #include "G4OpticalSurface.hh"
86 #include "G4OpticalPhoton.hh"
88 
89 // Class Description:
90 // Discrete Process -- reflection/refraction at optical interfaces.
91 // Class inherits publicly from G4VDiscreteProcess.
92 // Class Description - End:
93 
95 // Class Definition
97 
130 
132 {
133 
134 public:
135 
137  // Constructors and Destructor
139 
140  G4OpBoundaryProcess(const G4String& processName = "OpBoundary",
141  G4ProcessType type = fOptical);
143 
144 private:
145 
147 
149  // Operators
151 
153 
154 public:
155 
157  // Methods
159 
160  G4bool IsApplicable(const G4ParticleDefinition& aParticleType);
161  // Returns true -> 'is applicable' only for an optical photon.
162 
164  G4double ,
166  // Returns infinity; i. e. the process does not limit the step,
167  // but sets the 'Forced' condition for the DoIt to be invoked at
168  // every step. However, only at a boundary will any action be
169  // taken.
170 
171  G4VParticleChange* PostStepDoIt(const G4Track& aTrack,
172  const G4Step& aStep);
173  // This is the method implementing boundary processes.
174 
176  // Returns the current status.
177 
178 private:
179 
180  G4bool G4BooleanRand(const G4double prob) const;
181 
183  const G4ThreeVector& Normal) const;
184 
185  void DielectricMetal();
186  void DielectricDielectric();
187  void DielectricLUT();
188  void DielectricDichroic();
189 
190  void ChooseReflection();
191  void DoAbsorption();
192  void DoReflection();
193 
195  // Returns the incident angle of optical photon
196 
198  G4double E1_parl,
199  G4double incidentangle,
200  G4double RealRindex,
201  G4double ImaginaryRindex);
202  // Returns the Reflectivity on a metalic surface
203 
204  void CalculateReflectivity(void);
205 
206  void BoundaryProcessVerbose(void) const;
207 
208  // Invoke SD for post step point if the photon is 'detected'
209  G4bool InvokeSD(const G4Step* step);
210 
211 private:
212 
214 
217 
220 
223 
226 
228 
232 
235 
237 
239 
241 
243 
247 
249 
251 
253 
255 
256  size_t idx, idy;
258 };
259 
261 // Inline methods
263 
264 inline
266 {
267  /* Returns a random boolean variable with the specified probability */
268 
269  return (G4UniformRand() < prob);
270 }
271 
272 inline
274  aParticleType)
275 {
276  return ( &aParticleType == G4OpticalPhoton::OpticalPhoton() );
277 }
278 
279 inline
281 {
282  return theStatus;
283 }
284 
285 inline
287 {
288  G4double rand = G4UniformRand();
289  if ( rand >= 0.0 && rand < prob_ss ) {
292  }
293  else if ( rand >= prob_ss &&
294  rand <= prob_ss+prob_sl) {
296  }
297  else if ( rand > prob_ss+prob_sl &&
298  rand < prob_ss+prob_sl+prob_bs ) {
300  }
301  else {
303  }
304 }
305 
306 inline
308 {
310 
311  if ( G4BooleanRand(theEfficiency) ) {
312 
313  // EnergyDeposited =/= 0 means: photon has been detected
316  }
317  else {
319  }
320 
323 
324 // aParticleChange.ProposeEnergy(0.0);
326 }
327 
328 inline
330 {
331  if ( theStatus == LambertianReflection ) {
332 
335 
336  }
337  else if ( theFinish == ground ) {
338 
341  } else {
344  }
346  NewMomentum = OldMomentum - (2.*PdotN)*theFacetNormal;
347 
348  }
349  else {
350 
354  NewMomentum = OldMomentum - (2.*PdotN)*theFacetNormal;
355 
356  }
358  NewPolarization = -OldPolarization + (2.*EdotN)*theFacetNormal;
359 }
360 
361 #endif /* G4OpBoundaryProcess_h */
G4double condition(const G4ErrorSymMatrix &m)
G4double GetReflectivity(G4double E1_perp, G4double E1_parl, G4double incidentangle, G4double RealRindex, G4double ImaginaryRindex)
CLHEP::Hep3Vector G4ThreeVector
G4Physics2DVector * DichroicVector
G4ThreeVector G4LambertianRand(const G4ThreeVector &normal)
G4OpBoundaryProcess & operator=(const G4OpBoundaryProcess &right)
G4OpticalSurface * OpticalSurface
G4VParticleChange * PostStepDoIt(const G4Track &aTrack, const G4Step &aStep)
G4OpBoundaryProcessStatus
G4OpBoundaryProcess(const G4String &processName="OpBoundary", G4ProcessType type=fOptical)
G4OpBoundaryProcessStatus GetStatus() const
G4MaterialPropertyVector * PropertyPointer1
int G4int
Definition: G4Types.hh:78
G4OpticalSurfaceFinish theFinish
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4OpBoundaryProcessStatus theStatus
#define G4UniformRand()
Definition: Randomize.hh:95
G4bool G4BooleanRand(const G4double prob) const
bool G4bool
Definition: G4Types.hh:79
G4double GetMeanFreePath(const G4Track &, G4double, G4ForceCondition *condition)
void BoundaryProcessVerbose(void) const
Definition: G4Step.hh:76
G4OpticalSurfaceModel theModel
G4bool InvokeSD(const G4Step *step)
static G4OpticalPhoton * OpticalPhoton()
G4MaterialPropertyVector * PropertyPointer
G4OpticalSurfaceFinish
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289
G4bool IsApplicable(const G4ParticleDefinition &aParticleType)
G4OpticalSurfaceModel
G4MaterialPropertyVector * PropertyPointer2
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
G4ForceCondition
G4ThreeVector GetFacetNormal(const G4ThreeVector &Momentum, const G4ThreeVector &Normal) const
G4ProcessType