Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // 26 // 27 // 27 // 28 // << 28 // 29 ////////////////////////////////////////////// 29 //////////////////////////////////////////////////////////////////////// 30 // Optical Photon Rayleigh Scattering Class Im 30 // Optical Photon Rayleigh Scattering Class Implementation 31 ////////////////////////////////////////////// 31 //////////////////////////////////////////////////////////////////////// 32 // 32 // 33 // File: G4OpRayleigh.cc 33 // File: G4OpRayleigh.cc 34 // Description: Discrete Process -- Rayleigh s 34 // Description: Discrete Process -- Rayleigh scattering of optical 35 // photons 35 // photons 36 // Version: 1.0 36 // Version: 1.0 37 // Created: 1996-05-31 37 // Created: 1996-05-31 38 // Author: Juliet Armstrong 38 // Author: Juliet Armstrong 39 // Updated: 2014-10-10 - This version cal << 39 // Updated: 2014-10-10 - This version calculates the Rayleigh scattering 40 // length for more materials than 40 // length for more materials than just Water (although the Water 41 // default is kept). To do this t 41 // default is kept). To do this the user would need to specify the 42 // ISOTHERMAL_COMPRESSIBILITY as 42 // ISOTHERMAL_COMPRESSIBILITY as a material property and 43 // optionally an RS_SCALE_LENGTH 43 // optionally an RS_SCALE_LENGTH (useful for testing). Code comes 44 // from Philip Graham (Queen Mary 44 // from Philip Graham (Queen Mary University of London). 45 // 2010-06-11 - Fix Bug 207; Than 45 // 2010-06-11 - Fix Bug 207; Thanks to Xin Qian 46 // (Kellogg Radiation Lab of Calt 46 // (Kellogg Radiation Lab of Caltech) 47 // 2005-07-28 - add G4ProcessType 47 // 2005-07-28 - add G4ProcessType to constructor 48 // 2001-10-18 by Peter Gumplinger 48 // 2001-10-18 by Peter Gumplinger 49 // eliminate unused variable warn 49 // eliminate unused variable warning on Linux (gcc-2.95.2) 50 // 2001-09-18 by mma 50 // 2001-09-18 by mma 51 // >numOfMaterials=G4Material::Ge << 51 // >numOfMaterials=G4Material::GetNumberOfMaterials() in BuildPhy 52 // 2001-01-30 by Peter Gumplinger 52 // 2001-01-30 by Peter Gumplinger 53 // > allow for positiv and negati 53 // > allow for positiv and negative CosTheta and force the 54 // > new momentum direction to be 54 // > new momentum direction to be in the same plane as the 55 // > new and old polarization vec 55 // > new and old polarization vectors 56 // 2001-01-29 by Peter Gumplinger 56 // 2001-01-29 by Peter Gumplinger 57 // > fix calculation of SinTheta 57 // > fix calculation of SinTheta (from CosTheta) 58 // 1997-04-09 by Peter Gumplinger 58 // 1997-04-09 by Peter Gumplinger 59 // > new physics/tracking scheme 59 // > new physics/tracking scheme >> 60 // mail: gum@triumf.ca 60 // 61 // 61 ////////////////////////////////////////////// 62 //////////////////////////////////////////////////////////////////////// 62 63 63 #include "G4OpRayleigh.hh" 64 #include "G4OpRayleigh.hh" >> 65 64 #include "G4ios.hh" 66 #include "G4ios.hh" 65 #include "G4PhysicalConstants.hh" 67 #include "G4PhysicalConstants.hh" 66 #include "G4SystemOfUnits.hh" 68 #include "G4SystemOfUnits.hh" 67 #include "G4OpticalParameters.hh" << 68 #include "G4OpProcessSubType.hh" 69 #include "G4OpProcessSubType.hh" 69 70 70 //....oooOO0OOooo........oooOO0OOooo........oo << 71 ///////////////////////// >> 72 // Class Implementation >> 73 ///////////////////////// >> 74 >> 75 ////////////// >> 76 // Operators >> 77 ////////////// >> 78 >> 79 // G4OpRayleigh::operator=(const G4OpRayleigh &right) >> 80 // { >> 81 // } >> 82 >> 83 ///////////////// >> 84 // Constructors >> 85 ///////////////// >> 86 71 G4OpRayleigh::G4OpRayleigh(const G4String& pro 87 G4OpRayleigh::G4OpRayleigh(const G4String& processName, G4ProcessType type) 72 : G4VDiscreteProcess(processName, type) << 88 : G4VDiscreteProcess(processName, type) 73 { 89 { 74 Initialise(); << 90 SetProcessSubType(fOpRayleigh); 75 SetProcessSubType(fOpRayleigh); << 76 thePhysicsTable = nullptr; << 77 91 78 if(verboseLevel > 0) << 92 thePhysicsTable = NULL; 79 { << 80 G4cout << GetProcessName() << " is created << 81 } << 82 } << 83 93 84 //....oooOO0OOooo........oooOO0OOooo........oo << 94 if (verboseLevel>0) { 85 G4OpRayleigh::~G4OpRayleigh() << 95 G4cout << GetProcessName() << " is created " << G4endl; 86 { << 96 } 87 // VI: inside this PhysicsTable all properti << 88 // it is not possible to destroy << 89 delete thePhysicsTable; << 90 } 97 } 91 98 92 //....oooOO0OOooo........oooOO0OOooo........oo << 99 // G4OpRayleigh::G4OpRayleigh(const G4OpRayleigh &right) 93 void G4OpRayleigh::PreparePhysicsTable(const G << 100 // { 94 { << 101 // } 95 Initialise(); << 102 96 } << 103 //////////////// >> 104 // Destructors >> 105 //////////////// 97 106 98 //....oooOO0OOooo........oooOO0OOooo........oo << 107 G4OpRayleigh::~G4OpRayleigh() 99 void G4OpRayleigh::Initialise() << 100 { 108 { 101 SetVerboseLevel(G4OpticalParameters::Instanc << 109 if (thePhysicsTable) { >> 110 thePhysicsTable->clearAndDestroy(); >> 111 delete thePhysicsTable; >> 112 } 102 } 113 } 103 114 104 //....oooOO0OOooo........oooOO0OOooo........oo << 115 //////////// 105 G4VParticleChange* G4OpRayleigh::PostStepDoIt( << 116 // Methods 106 << 117 //////////// 107 { << 108 aParticleChange.Initialize(aTrack); << 109 const G4DynamicParticle* aParticle = aTrack. << 110 << 111 if(verboseLevel > 1) << 112 { << 113 G4cout << "OpRayleigh: Scattering Photon!" << 114 << "Old Momentum Direction: " << aP << 115 << G4endl << "Old Polarization: " < << 116 << G4endl; << 117 } << 118 << 119 G4double cosTheta; << 120 G4ThreeVector oldMomDir, newMomDir; << 121 G4ThreeVector oldPol, newPol; << 122 G4double rand; << 123 G4double cost, sint, sinphi, cosphi; << 124 118 125 do << 119 // PostStepDoIt 126 { << 120 // ------------- 127 // Try to simulate the scattered photon mo << 121 // 128 // w.r.t. the initial photon momentum dire << 122 G4VParticleChange* 129 cost = G4UniformRand(); << 123 G4OpRayleigh::PostStepDoIt(const G4Track& aTrack, const G4Step& aStep) 130 sint = std::sqrt(1. - cost * cost); << 124 { 131 // consider for the angle 90-180 degrees << 125 aParticleChange.Initialize(aTrack); 132 if(G4UniformRand() < 0.5) << 133 cost = -cost; << 134 << 135 // simulate the phi angle << 136 rand = twopi * G4UniformRand(); << 137 sinphi = std::sin(rand); << 138 cosphi = std::cos(rand); << 139 << 140 // construct the new momentum direction << 141 newMomDir.set(sint * cosphi, sint * sinphi << 142 oldMomDir = aParticle->GetMomentumDirectio << 143 newMomDir.rotateUz(oldMomDir); << 144 << 145 // calculate the new polarization directio << 146 // The new polarization needs to be in the << 147 // momentum direction and the old polariza << 148 oldPol = aParticle->GetPolarization(); << 149 newPol = (oldPol - newMomDir.dot(oldPol) * << 150 << 151 // There is a corner case, where the new m << 152 // is the same as old polarization directi << 153 // random generate the azimuthal angle w.r << 154 if(newPol.mag() == 0.) << 155 { << 156 rand = G4UniformRand() * twopi; << 157 newPol.set(std::cos(rand), std::sin(rand << 158 newPol.rotateUz(newMomDir); << 159 } << 160 else << 161 { << 162 // There are two directions perpendicula << 163 if(G4UniformRand() < 0.5) << 164 newPol = -newPol; << 165 } << 166 << 167 // simulate according to the distribution << 168 cosTheta = newPol.dot(oldPol); << 169 // Loop checking, 13-Aug-2015, Peter Gumpl << 170 } while(std::pow(cosTheta, 2) < G4UniformRan << 171 126 172 aParticleChange.ProposePolarization(newPol); << 127 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle(); 173 aParticleChange.ProposeMomentumDirection(new << 174 128 175 if(verboseLevel > 1) << 129 if (verboseLevel>0) { 176 { << 130 G4cout << "Scattering Photon!" << G4endl; 177 G4cout << "New Polarization: " << newPol < << 131 G4cout << "Old Momentum Direction: " 178 << "Polarization Change: " << *(aPa << 132 << aParticle->GetMomentumDirection() << G4endl; 179 << G4endl << "New Momentum Directio << 133 G4cout << "Old Polarization: " 180 << "Momentum Change: " << *(aPartic << 134 << aParticle->GetPolarization() << G4endl; 181 << G4endl; << 135 } 182 } << 136 >> 137 G4double cosTheta; >> 138 G4ThreeVector OldMomentumDirection, NewMomentumDirection; >> 139 G4ThreeVector OldPolarization, NewPolarization; >> 140 >> 141 G4double rand, constant; >> 142 G4double CosTheta, SinTheta, SinPhi, CosPhi, unit_x, unit_y, unit_z; >> 143 >> 144 do { >> 145 // Try to simulate the scattered photon momentum direction >> 146 // w.r.t. the initial photon momentum direction >> 147 >> 148 CosTheta = G4UniformRand(); >> 149 SinTheta = std::sqrt(1.-CosTheta*CosTheta); >> 150 // consider for the angle 90-180 degrees >> 151 if (G4UniformRand() < 0.5) CosTheta = -CosTheta; >> 152 >> 153 // simulate the phi angle >> 154 rand = twopi*G4UniformRand(); >> 155 SinPhi = std::sin(rand); >> 156 CosPhi = std::cos(rand); >> 157 >> 158 // start constructing the new momentum direction >> 159 unit_x = SinTheta * CosPhi; >> 160 unit_y = SinTheta * SinPhi; >> 161 unit_z = CosTheta; >> 162 NewMomentumDirection.set (unit_x,unit_y,unit_z); >> 163 >> 164 // Rotate the new momentum direction into global reference system >> 165 OldMomentumDirection = aParticle->GetMomentumDirection(); >> 166 OldMomentumDirection = OldMomentumDirection.unit(); >> 167 NewMomentumDirection.rotateUz(OldMomentumDirection); >> 168 NewMomentumDirection = NewMomentumDirection.unit(); >> 169 >> 170 // calculate the new polarization direction >> 171 // The new polarization needs to be in the same plane as the new >> 172 // momentum direction and the old polarization direction >> 173 OldPolarization = aParticle->GetPolarization(); >> 174 constant = -NewMomentumDirection.dot(OldPolarization); >> 175 >> 176 NewPolarization = OldPolarization + constant*NewMomentumDirection; >> 177 NewPolarization = NewPolarization.unit(); >> 178 >> 179 // There is a corner case, where the Newmomentum direction >> 180 // is the same as oldpolariztion direction: >> 181 // random generate the azimuthal angle w.r.t. Newmomentum direction >> 182 if (NewPolarization.mag() == 0.) { >> 183 rand = G4UniformRand()*twopi; >> 184 NewPolarization.set(std::cos(rand),std::sin(rand),0.); >> 185 NewPolarization.rotateUz(NewMomentumDirection); >> 186 } else { >> 187 // There are two directions which are perpendicular >> 188 // to the new momentum direction >> 189 if (G4UniformRand() < 0.5) NewPolarization = -NewPolarization; >> 190 } >> 191 >> 192 // simulate according to the distribution cos^2(theta) >> 193 cosTheta = NewPolarization.dot(OldPolarization); >> 194 // Loop checking, 13-Aug-2015, Peter Gumplinger >> 195 } while (std::pow(cosTheta,2) < G4UniformRand()); >> 196 >> 197 aParticleChange.ProposePolarization(NewPolarization); >> 198 aParticleChange.ProposeMomentumDirection(NewMomentumDirection); >> 199 >> 200 if (verboseLevel>0) { >> 201 G4cout << "New Polarization: " >> 202 << NewPolarization << G4endl; >> 203 G4cout << "Polarization Change: " >> 204 << *(aParticleChange.GetPolarization()) << G4endl; >> 205 G4cout << "New Momentum Direction: " >> 206 << NewMomentumDirection << G4endl; >> 207 G4cout << "Momentum Change: " >> 208 << *(aParticleChange.GetMomentumDirection()) << G4endl; >> 209 } 183 210 184 return G4VDiscreteProcess::PostStepDoIt(aTra << 211 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep); 185 } 212 } 186 213 187 //....oooOO0OOooo........oooOO0OOooo........oo << 214 // BuildPhysicsTable for the Rayleigh Scattering process >> 215 // -------------------------------------------------------- 188 void G4OpRayleigh::BuildPhysicsTable(const G4P 216 void G4OpRayleigh::BuildPhysicsTable(const G4ParticleDefinition&) 189 { 217 { 190 if(thePhysicsTable) << 218 if (thePhysicsTable) { 191 { << 219 thePhysicsTable->clearAndDestroy(); 192 // thePhysicsTable->clearAndDestroy(); << 220 delete thePhysicsTable; 193 delete thePhysicsTable; << 221 thePhysicsTable = NULL; 194 thePhysicsTable = nullptr; << 195 } 222 } 196 223 197 const G4MaterialTable* theMaterialTable = G4 224 const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable(); 198 const size_t numOfMaterials = G4 << 225 const G4int numOfMaterials = G4Material::GetNumberOfMaterials(); 199 thePhysicsTable = ne << 200 226 201 for(size_t i = 0; i < numOfMaterials; ++i) << 227 thePhysicsTable = new G4PhysicsTable( numOfMaterials ); 202 { << 228 203 G4Material* material = (*the << 229 for( G4int iMaterial = 0; iMaterial < numOfMaterials; iMaterial++ ) 204 G4MaterialPropertiesTable* matProp = mater << 230 { 205 G4PhysicsFreeVector* rayleigh = nullptr; << 231 G4Material* material = (*theMaterialTable)[iMaterial]; 206 if(matProp) << 232 G4MaterialPropertiesTable* materialProperties = 207 { << 233 material->GetMaterialPropertiesTable(); 208 rayleigh = matProp->GetProperty(kRAYLEIG << 234 G4PhysicsOrderedFreeVector* rayleigh = NULL; 209 if(rayleigh == nullptr) << 235 if ( materialProperties != NULL ) { 210 rayleigh = CalculateRayleighMeanFreePa << 236 rayleigh = materialProperties->GetProperty( kRAYLEIGH ); 211 } << 237 if ( rayleigh == NULL ) rayleigh = 212 thePhysicsTable->insertAt(i, rayleigh); << 238 CalculateRayleighMeanFreePaths( material ); >> 239 } >> 240 thePhysicsTable->insertAt( iMaterial, rayleigh ); 213 } 241 } 214 } 242 } 215 243 216 //....oooOO0OOooo........oooOO0OOooo........oo << 244 // GetMeanFreePath() 217 G4double G4OpRayleigh::GetMeanFreePath(const G << 245 // ----------------- 218 G4Force << 246 // >> 247 G4double G4OpRayleigh::GetMeanFreePath(const G4Track& aTrack, >> 248 G4double , >> 249 G4ForceCondition* ) 219 { 250 { 220 auto rayleigh = static_cast<G4PhysicsFreeVec << 251 const G4DynamicParticle* particle = aTrack.GetDynamicParticle(); 221 (*thePhysicsTable)(aTrack.GetMaterial()- << 252 const G4double photonMomentum = particle->GetTotalMomentum(); 222 << 253 const G4Material* material = aTrack.GetMaterial(); >> 254 >> 255 G4PhysicsOrderedFreeVector* rayleigh = >> 256 static_cast<G4PhysicsOrderedFreeVector*> >> 257 ((*thePhysicsTable)(material->GetIndex())); >> 258 223 G4double rsLength = DBL_MAX; 259 G4double rsLength = DBL_MAX; 224 if(rayleigh) << 260 if( rayleigh != NULL ) rsLength = rayleigh->Value( photonMomentum ); 225 { << 226 rsLength = rayleigh->Value(aTrack.GetDynam << 227 idx_rslength); << 228 } << 229 return rsLength; 261 return rsLength; 230 } 262 } 231 263 232 //....oooOO0OOooo........oooOO0OOooo........oo << 264 // CalculateRayleighMeanFreePaths() 233 G4PhysicsFreeVector* G4OpRayleigh::CalculateRa << 265 // -------------------------------- 234 const G4Material* material) const << 266 // Private method to compute Rayleigh Scattering Lengths >> 267 G4PhysicsOrderedFreeVector* >> 268 G4OpRayleigh::CalculateRayleighMeanFreePaths( const G4Material* material ) const 235 { 269 { 236 G4MaterialPropertiesTable* MPT = material->G << 270 G4MaterialPropertiesTable* materialProperties = >> 271 material->GetMaterialPropertiesTable(); 237 272 238 // Retrieve the beta_T or isothermal compres 273 // Retrieve the beta_T or isothermal compressibility value. For backwards 239 // compatibility use a constant if the mater 274 // compatibility use a constant if the material is "Water". If the material 240 // doesn't have an ISOTHERMAL_COMPRESSIBILIT 275 // doesn't have an ISOTHERMAL_COMPRESSIBILITY constant then return 241 G4double betat; 276 G4double betat; 242 if(material->GetName() == "Water") << 277 if ( material->GetName() == "Water" ) 243 { << 278 betat = 7.658e-23*m3/MeV; 244 betat = 7.658e-23 * m3 / MeV; << 279 else if(materialProperties->ConstPropertyExists("ISOTHERMAL_COMPRESSIBILITY")) 245 } << 280 betat = materialProperties->GetConstProperty(kISOTHERMAL_COMPRESSIBILITY); 246 else if(MPT->ConstPropertyExists(kISOTHERMAL << 247 { << 248 betat = MPT->GetConstProperty(kISOTHERMAL_ << 249 } << 250 else 281 else 251 { << 282 return NULL; 252 return nullptr; << 253 } << 254 283 255 // If the material doesn't have a RINDEX pro 284 // If the material doesn't have a RINDEX property vector then return 256 G4MaterialPropertyVector* rIndex = MPT->GetP << 285 G4MaterialPropertyVector* rIndex = materialProperties->GetProperty(kRINDEX); 257 if(rIndex == nullptr) << 286 if ( rIndex == NULL ) return NULL; 258 return nullptr; << 259 287 260 // Retrieve the optional scale factor (scale << 288 // Retrieve the optional scale factor, (this just scales the scattering length 261 G4double scaleFactor = 1.0; 289 G4double scaleFactor = 1.0; 262 if(MPT->ConstPropertyExists(kRS_SCALE_FACTOR << 290 if( materialProperties->ConstPropertyExists( "RS_SCALE_FACTOR" ) ) 263 { << 291 scaleFactor= materialProperties->GetConstProperty(kRS_SCALE_FACTOR ); 264 scaleFactor = MPT->GetConstProperty(kRS_SC << 265 } << 266 292 267 // Retrieve the material temperature. For ba << 293 // Retrieve the material temperature. For backwards compatibility use a 268 // constant if the material is "Water" 294 // constant if the material is "Water" 269 G4double temperature; 295 G4double temperature; 270 if(material->GetName() == "Water") << 296 if( material->GetName() == "Water" ) 271 { << 297 temperature = 283.15*kelvin; // Temperature of water is 10 degrees celsius 272 temperature = << 273 283.15 * kelvin; // Temperature of wate << 274 } << 275 else 298 else 276 { << 277 temperature = material->GetTemperature(); 299 temperature = material->GetTemperature(); 278 } << 279 300 280 auto rayleighMFPs = new G4PhysicsFreeVector( << 301 G4PhysicsOrderedFreeVector* rayleighMeanFreePaths = >> 302 new G4PhysicsOrderedFreeVector(); 281 // This calculates the meanFreePath via the 303 // This calculates the meanFreePath via the Einstein-Smoluchowski formula 282 const G4double c1 = << 304 const G4double c1 = scaleFactor * betat * temperature * k_Boltzmann / 283 scaleFactor * betat * temperature * k_Bolt << 305 ( 6.0 * pi ); 284 306 285 for(size_t uRIndex = 0; uRIndex < rIndex->Ge << 307 for( size_t uRIndex = 0; uRIndex < rIndex->GetVectorLength(); uRIndex++ ) 286 { 308 { 287 const G4double energy = rIndex->Ene << 309 const G4double energy = rIndex->Energy( uRIndex ); 288 const G4double rIndexSquared = (*rIndex)[u << 310 const G4double rIndexSquared = (*rIndex)[uRIndex] * (*rIndex)[uRIndex]; 289 const G4double xlambda = h_Planck * << 311 const G4double xlambda = h_Planck * c_light / energy; 290 const G4double c2 = std::pow(tw << 312 const G4double c2 = std::pow(twopi/xlambda,4); 291 const G4double c3 = << 313 const G4double c3 = 292 std::pow(((rIndexSquared - 1.0) * (rInde << 314 std::pow(((rIndexSquared-1.0)*(rIndexSquared+2.0 )/3.0),2); 293 << 294 const G4double meanFreePath = 1.0 / (c1 * << 295 << 296 if(verboseLevel > 0) << 297 { << 298 G4cout << energy << "MeV\t" << meanFreeP << 299 } << 300 315 301 rayleighMFPs->InsertValues(energy, meanFre << 316 const G4double meanFreePath = 1.0 / ( c1 * c2 * c3 ); 302 } << 303 317 304 return rayleighMFPs; << 318 if( verboseLevel>0 ) 305 } << 319 G4cout << energy << "MeV\t" << meanFreePath << "mm" << G4endl; 306 320 307 //....oooOO0OOooo........oooOO0OOooo........oo << 321 rayleighMeanFreePaths->InsertValues( energy, meanFreePath ); 308 void G4OpRayleigh::SetVerboseLevel(G4int verbo << 322 } 309 { << 323 310 verboseLevel = verbose; << 324 return rayleighMeanFreePaths; 311 G4OpticalParameters::Instance()->SetRayleigh << 312 } 325 } 313 326