Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 96023 2016-03-08 08:25:42Z 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 
152  G4OpBoundaryProcess& operator=(const G4OpBoundaryProcess &right);
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  void SetInvokeSD(G4bool );
179  // Set flag for call to InvokeSD method.
180 
181 private:
182 
183  G4bool G4BooleanRand(const G4double prob) const;
184 
185  G4ThreeVector GetFacetNormal(const G4ThreeVector& Momentum,
186  const G4ThreeVector& Normal) const;
187 
188  void DielectricMetal();
189  void DielectricDielectric();
190  void DielectricLUT();
191  void DielectricDichroic();
192 
193  void ChooseReflection();
194  void DoAbsorption();
195  void DoReflection();
196 
197  G4double GetIncidentAngle();
198  // Returns the incident angle of optical photon
199 
200  G4double GetReflectivity(G4double E1_perp,
201  G4double E1_parl,
202  G4double incidentangle,
203  G4double RealRindex,
204  G4double ImaginaryRindex);
205  // Returns the Reflectivity on a metalic surface
206 
207  void CalculateReflectivity(void);
208 
209  void BoundaryProcessVerbose(void) const;
210 
211  // Invoke SD for post step point if the photon is 'detected'
212  G4bool InvokeSD(const G4Step* step);
213 
214 private:
215 
216  G4double thePhotonMomentum;
217 
218  G4ThreeVector OldMomentum;
219  G4ThreeVector OldPolarization;
220 
221  G4ThreeVector NewMomentum;
222  G4ThreeVector NewPolarization;
223 
224  G4ThreeVector theGlobalNormal;
225  G4ThreeVector theFacetNormal;
226 
227  G4Material* Material1;
228  G4Material* Material2;
229 
230  G4OpticalSurface* OpticalSurface;
231 
232  G4MaterialPropertyVector* PropertyPointer;
233  G4MaterialPropertyVector* PropertyPointer1;
234  G4MaterialPropertyVector* PropertyPointer2;
235 
236  G4double Rindex1;
237  G4double Rindex2;
238 
239  G4double cost1, cost2, sint1, sint2;
240 
241  G4OpBoundaryProcessStatus theStatus;
242 
243  G4OpticalSurfaceModel theModel;
244 
245  G4OpticalSurfaceFinish theFinish;
246 
247  G4double theReflectivity;
248  G4double theEfficiency;
249  G4double theTransmittance;
250 
251  G4double theSurfaceRoughness;
252 
253  G4double prob_sl, prob_ss, prob_bs;
254 
255  G4int iTE, iTM;
256 
257  G4double kCarTolerance;
258 
259  size_t idx, idy;
260  G4Physics2DVector* DichroicVector;
261 
262  G4bool fInvokeSD;
263 };
264 
266 // Inline methods
268 
269 inline
270 G4bool G4OpBoundaryProcess::G4BooleanRand(const G4double prob) const
271 {
272  /* Returns a random boolean variable with the specified probability */
273 
274  return (G4UniformRand() < prob);
275 }
276 
277 inline
279  aParticleType)
280 {
281  return ( &aParticleType == G4OpticalPhoton::OpticalPhoton() );
282 }
283 
284 inline
286 {
287  return theStatus;
288 }
289 
290 inline
292 {
293  fInvokeSD = flag;
294 }
295 
296 inline
297 void G4OpBoundaryProcess::ChooseReflection()
298 {
299  G4double rand = G4UniformRand();
300  if ( rand >= 0.0 && rand < prob_ss ) {
301  theStatus = SpikeReflection;
302  theFacetNormal = theGlobalNormal;
303  }
304  else if ( rand >= prob_ss &&
305  rand <= prob_ss+prob_sl) {
306  theStatus = LobeReflection;
307  }
308  else if ( rand > prob_ss+prob_sl &&
309  rand < prob_ss+prob_sl+prob_bs ) {
310  theStatus = BackScattering;
311  }
312  else {
313  theStatus = LambertianReflection;
314  }
315 }
316 
317 inline
318 void G4OpBoundaryProcess::DoAbsorption()
319 {
320  theStatus = Absorption;
321 
322  if ( G4BooleanRand(theEfficiency) ) {
323 
324  // EnergyDeposited =/= 0 means: photon has been detected
325  theStatus = Detection;
326  aParticleChange.ProposeLocalEnergyDeposit(thePhotonMomentum);
327  }
328  else {
330  }
331 
332  NewMomentum = OldMomentum;
333  NewPolarization = OldPolarization;
334 
335 // aParticleChange.ProposeEnergy(0.0);
337 }
338 
339 inline
340 void G4OpBoundaryProcess::DoReflection()
341 {
342  if ( theStatus == LambertianReflection ) {
343 
344  NewMomentum = G4LambertianRand(theGlobalNormal);
345  theFacetNormal = (NewMomentum - OldMomentum).unit();
346 
347  }
348  else if ( theFinish == ground ) {
349 
350  theStatus = LobeReflection;
351  if ( PropertyPointer1 && PropertyPointer2 ){
352  } else {
353  theFacetNormal =
354  GetFacetNormal(OldMomentum,theGlobalNormal);
355  }
356  G4double PdotN = OldMomentum * theFacetNormal;
357  NewMomentum = OldMomentum - (2.*PdotN)*theFacetNormal;
358 
359  }
360  else {
361 
362  theStatus = SpikeReflection;
363  theFacetNormal = theGlobalNormal;
364  G4double PdotN = OldMomentum * theFacetNormal;
365  NewMomentum = OldMomentum - (2.*PdotN)*theFacetNormal;
366 
367  }
368  G4double EdotN = OldPolarization * theFacetNormal;
369  NewPolarization = -OldPolarization + (2.*EdotN)*theFacetNormal;
370 }
371 
372 #endif /* G4OpBoundaryProcess_h */
G4double condition(const G4ErrorSymMatrix &m)
G4ThreeVector G4LambertianRand(const G4ThreeVector &normal)
G4VParticleChange * PostStepDoIt(const G4Track &aTrack, const G4Step &aStep)
G4OpBoundaryProcessStatus
G4OpBoundaryProcess(const G4String &processName="OpBoundary", G4ProcessType type=fOptical)
G4OpBoundaryProcessStatus GetStatus() const
int G4int
Definition: G4Types.hh:78
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
#define G4UniformRand()
Definition: Randomize.hh:97
bool G4bool
Definition: G4Types.hh:79
G4double GetMeanFreePath(const G4Track &, G4double, G4ForceCondition *condition)
Definition: G4Step.hh:76
static G4OpticalPhoton * OpticalPhoton()
G4OpticalSurfaceFinish
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289
G4bool IsApplicable(const G4ParticleDefinition &aParticleType)
G4OpticalSurfaceModel
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
G4ForceCondition
G4ProcessType