LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
larg4::OpBoundaryProcessSimple Class Reference

Discrete process for reflection and diffusion at optical interfaces. More...

#include "OpBoundaryProcessSimple.hh"

Inheritance diagram for larg4::OpBoundaryProcessSimple:

Public Member Functions

 OpBoundaryProcessSimple (const G4String &processName="OpBoundary", G4ProcessType type=fOptical)
 
G4bool IsApplicable (const G4ParticleDefinition &aParticleType)
 
G4double GetMeanFreePath (const G4Track &, G4double, G4ForceCondition *condition)
 
G4VParticleChange * PostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
OpBoundaryProcessSimpleStatus GetStatus () const
 

Private Member Functions

G4bool G4BooleanRand (const G4double prob) const
 

Private Attributes

OpBoundaryProcessSimpleStatus fTheStatus
 
G4double fCarTolerance
 
int fVerbosity
 

Detailed Description

Discrete process for reflection and diffusion at optical interfaces.

See also
G4VDiscreteProcess

This class invokes a simplified model of optical reflections at boundaries between different materials. The relevant reflectivities are ultimately read from detinfo::LArProperties via larg4::MaterialPropertiesLoader.

The required parameters are total reflectance (detinfo::LArProperties::SurfaceReflectances()) and ratio of diffuse to specular reflectance (detinfo::LArProperties::SurfaceReflectanceDiffuseFractions()). Each photon crossing a boundary with a defined reflectance is randomly either reflected or absorbed and killed according to the supplied probability.

Every reflected photon with a defined diffuse reflection fraction is then randomly either diffusely or specularly reflected according to the supplied probability. All materials with no defined reflectance are assumed to be black and absorb all incident photons.

This physics process is loaded in larg4::OpticalPhysics physics constructor.

This class is based on the G4OpBoundaryProcess class in Geant4 and was adapted for LArSoft by Ben Jones, MIT, March 2010.

Definition at line 102 of file OpBoundaryProcessSimple.hh.

Constructor & Destructor Documentation

larg4::OpBoundaryProcessSimple::OpBoundaryProcessSimple ( const G4String processName = "OpBoundary",
G4ProcessType  type = fOptical 
)

Definition at line 67 of file OpBoundaryProcessSimple.cxx.

References fCarTolerance, fTheStatus, and larg4::Undefined.

68  : G4VDiscreteProcess(processName, type)
69  {
70  G4cout << GetProcessName() << " is created " << G4endl;
71 
72  SetProcessSubType(fOpBoundary);
73 
75 
76  fCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
77  }
OpBoundaryProcessSimpleStatus fTheStatus

Member Function Documentation

G4bool larg4::OpBoundaryProcessSimple::G4BooleanRand ( const G4double  prob) const
inlineprivate

Definition at line 132 of file OpBoundaryProcessSimple.hh.

Referenced by PostStepDoIt().

133  {
134  /* Returns a random boolean variable with the specified probability */
135  return (G4UniformRand() < prob);
136  }
G4double larg4::OpBoundaryProcessSimple::GetMeanFreePath ( const G4Track &  ,
G4double  ,
G4ForceCondition *  condition 
)

Definition at line 234 of file OpBoundaryProcessSimple.cxx.

237  {
238  *condition = Forced;
239 
240  return DBL_MAX;
241  }
OpBoundaryProcessSimpleStatus larg4::OpBoundaryProcessSimple::GetStatus ( ) const
inline

Definition at line 143 of file OpBoundaryProcessSimple.hh.

References fTheStatus.

144  {
145  return fTheStatus;
146  }
OpBoundaryProcessSimpleStatus fTheStatus
G4bool larg4::OpBoundaryProcessSimple::IsApplicable ( const G4ParticleDefinition &  aParticleType)
inline

Definition at line 138 of file OpBoundaryProcessSimple.hh.

139  {
140  return (&aParticleType == G4OpticalPhoton::OpticalPhoton());
141  }
G4VParticleChange * larg4::OpBoundaryProcessSimple::PostStepDoIt ( const G4Track &  aTrack,
const G4Step &  aStep 
)

Definition at line 82 of file OpBoundaryProcessSimple.cxx.

References fCarTolerance, fTheStatus, fVerbosity, G4BooleanRand(), larg4::NotAtBoundary, sim::LArG4Parameters::OpVerbosity(), larg4::SimpleAbsorbed, larg4::SimpleAbsorbedNoRefl, larg4::SimpleDiffuse, larg4::SimpleSpecular, larg4::StepTooSmall, and larg4::Undefined.

84  {
85 
86  // Note - these have to be loaded here, since this object
87  // is constructed before LArG4, and hence before the
88  // LArG4 parameters are read from config
89 
91  aParticleChange.Initialize(aTrack);
92 
93  G4StepPoint* pPreStepPoint = aStep.GetPreStepPoint();
94  G4StepPoint* pPostStepPoint = aStep.GetPostStepPoint();
95 
96  if (pPostStepPoint->GetStepStatus() != fGeomBoundary) {
98  return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
99  }
100  if (aTrack.GetStepLength() <= fCarTolerance / 2) {
102  return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
103  }
104 
105  G4Material* Material1 = pPreStepPoint->GetMaterial();
106  G4Material* Material2 = pPostStepPoint->GetMaterial();
107 
108  if (Material1 != Material2) {
109 
111 
112  fVerbosity = lgp->OpVerbosity();
113 
114  G4ThreeVector NewMomentum;
115  G4ThreeVector NewPolarization;
116 
117  G4ThreeVector theGlobalNormal;
118  G4ThreeVector theFacetNormal;
119 
120  std::string Material1Name = Material1->GetName();
121  std::string Material2Name = Material2->GetName();
122 
123  const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
124  G4double thePhotonMomentum = aParticle->GetTotalMomentum();
125  G4ThreeVector OldMomentum = aParticle->GetMomentumDirection();
126  G4ThreeVector OldPolarization = aParticle->GetPolarization();
127 
128  if (fVerbosity > 9)
129  std::cout << "OpBoundaryProcessSimple Debug: Photon " << aTrack.GetTrackID() << " momentum "
130  << thePhotonMomentum << " Material1 " << Material1Name.c_str() << " Material 2 "
131  << Material2Name.c_str() << std::endl;
132 
133  G4ThreeVector theGlobalPoint = pPostStepPoint->GetPosition();
134 
135  G4Navigator* theNavigator =
136  G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking();
137 
138  G4ThreeVector theLocalPoint =
139  theNavigator->GetGlobalToLocalTransform().TransformPoint(theGlobalPoint);
140 
141  // Normal points back into volume
142  G4bool valid;
143  G4ThreeVector theLocalNormal = theNavigator->GetLocalExitNormal(&valid);
144 
145  if (valid) { theLocalNormal = -theLocalNormal; }
146  else {
147  G4cerr << " OpBoundaryProcessSimple/PostStepDoIt(): "
148  << " The Navigator reports that it returned an invalid normal" << G4endl;
149  }
150 
151  theGlobalNormal = theNavigator->GetLocalToGlobalTransform().TransformAxis(theLocalNormal);
152 
153  if (OldMomentum * theGlobalNormal > 0.0) {
154 #ifdef G4DEBUG_OPTICAL
155  G4cerr << " OpBoundaryProcessSimple/PostStepDoIt(): "
156  << " theGlobalNormal points the wrong direction " << G4endl;
157 #endif
158  theGlobalNormal = -theGlobalNormal;
159  }
160 
161  G4MaterialPropertiesTable* aMaterialPropertiesTable;
162  G4MaterialPropertyVector* Reflectance;
163  G4MaterialPropertyVector* DiffuseReflectanceFraction;
164 
165  aMaterialPropertiesTable = Material1->GetMaterialPropertiesTable();
166 
167  if (aMaterialPropertiesTable) {
168 
169  std::stringstream PropertyName;
170  PropertyName << "REFLECTANCE_" << Material2Name.c_str();
171  Reflectance = aMaterialPropertiesTable->GetProperty(PropertyName.str().c_str());
172 
173  PropertyName.str("");
174  PropertyName << "DIFFUSE_REFLECTANCE_FRACTION_" << Material2Name.c_str();
175  DiffuseReflectanceFraction =
176  aMaterialPropertiesTable->GetProperty(PropertyName.str().c_str());
177 
178  if (Reflectance) {
179  double theReflectance = Reflectance->Value(thePhotonMomentum);
180  if (G4BooleanRand(theReflectance)) {
181  if (DiffuseReflectanceFraction) {
182  double theDiffuseReflectanceFraction =
183  DiffuseReflectanceFraction->Value(thePhotonMomentum);
184  if (G4BooleanRand(theDiffuseReflectanceFraction)) { fTheStatus = SimpleDiffuse; }
185  else {
187  }
188  }
189  else {
191  }
192  }
193  else {
195  }
196  }
197  else {
199  }
200  }
201 
202  // Take action according to track status
203 
204  if (fTheStatus == SimpleDiffuse) {
205  NewMomentum = G4LambertianRand(theGlobalNormal);
206  theFacetNormal = (NewMomentum - OldMomentum).unit();
207  }
208 
209  else if (fTheStatus == SimpleSpecular) {
210  theFacetNormal = theGlobalNormal;
211  G4double PdotN = OldMomentum * theFacetNormal;
212  NewMomentum = OldMomentum - (2. * PdotN) * theFacetNormal;
213  }
214 
216  aParticleChange.ProposeTrackStatus(fStopAndKill);
217  }
218 
219  else if ((fTheStatus == NotAtBoundary) || (fTheStatus == StepTooSmall)) {
220  NewMomentum = OldMomentum;
221  }
222  else
223  aParticleChange.ProposeTrackStatus(fStopAndKill);
224 
225  NewMomentum = NewMomentum.unit();
226  aParticleChange.ProposeMomentumDirection(NewMomentum);
227  }
228  return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
229  }
OpBoundaryProcessSimpleStatus fTheStatus
int OpVerbosity() const
G4bool G4BooleanRand(const G4double prob) const

Member Data Documentation

G4double larg4::OpBoundaryProcessSimple::fCarTolerance
private

Definition at line 127 of file OpBoundaryProcessSimple.hh.

Referenced by OpBoundaryProcessSimple(), and PostStepDoIt().

OpBoundaryProcessSimpleStatus larg4::OpBoundaryProcessSimple::fTheStatus
private

Definition at line 126 of file OpBoundaryProcessSimple.hh.

Referenced by GetStatus(), OpBoundaryProcessSimple(), and PostStepDoIt().

int larg4::OpBoundaryProcessSimple::fVerbosity
private

Definition at line 129 of file OpBoundaryProcessSimple.hh.

Referenced by PostStepDoIt().


The documentation for this class was generated from the following files: