#include <G4OpRayleigh.hh>
Inheritance diagram for G4OpRayleigh:
Public Member Functions | |
G4OpRayleigh (const G4String &processName="OpRayleigh", G4ProcessType type=fOptical) | |
~G4OpRayleigh () | |
G4bool | IsApplicable (const G4ParticleDefinition &aParticleType) |
G4double | GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *) |
G4VParticleChange * | PostStepDoIt (const G4Track &aTrack, const G4Step &aStep) |
G4PhysicsTable * | GetPhysicsTable () const |
void | DumpPhysicsTable () const |
Protected Attributes | |
G4PhysicsTable * | thePhysicsTable |
Definition at line 76 of file G4OpRayleigh.hh.
G4OpRayleigh::G4OpRayleigh | ( | const G4String & | processName = "OpRayleigh" , |
|
G4ProcessType | type = fOptical | |||
) |
Definition at line 82 of file G4OpRayleigh.cc.
References fOpRayleigh, G4cout, G4endl, G4VProcess::GetProcessName(), G4VProcess::SetProcessSubType(), thePhysicsTable, and G4VProcess::verboseLevel.
00083 : G4VDiscreteProcess(processName, type) 00084 { 00085 SetProcessSubType(fOpRayleigh); 00086 00087 thePhysicsTable = 0; 00088 00089 DefaultWater = false; 00090 00091 if (verboseLevel>0) { 00092 G4cout << GetProcessName() << " is created " << G4endl; 00093 } 00094 00095 BuildThePhysicsTable(); 00096 }
G4OpRayleigh::~G4OpRayleigh | ( | ) |
Definition at line 106 of file G4OpRayleigh.cc.
References G4PhysicsTable::clearAndDestroy(), and thePhysicsTable.
00107 { 00108 if (thePhysicsTable!= 0) { 00109 thePhysicsTable->clearAndDestroy(); 00110 delete thePhysicsTable; 00111 } 00112 }
void G4OpRayleigh::DumpPhysicsTable | ( | ) | const [inline] |
Definition at line 163 of file G4OpRayleigh.hh.
References G4PhysicsOrderedFreeVector::DumpValues(), G4PhysicsTable::entries(), and thePhysicsTable.
00165 { 00166 G4int PhysicsTableSize = thePhysicsTable->entries(); 00167 G4PhysicsOrderedFreeVector *v; 00168 00169 for (G4int i = 0 ; i < PhysicsTableSize ; i++ ) 00170 { 00171 v = (G4PhysicsOrderedFreeVector*)(*thePhysicsTable)[i]; 00172 v->DumpValues(); 00173 } 00174 }
G4double G4OpRayleigh::GetMeanFreePath | ( | const G4Track & | aTrack, | |
G4double | , | |||
G4ForceCondition * | ||||
) | [virtual] |
Implements G4VDiscreteProcess.
Definition at line 262 of file G4OpRayleigh.cc.
References DBL_MAX, G4Track::GetDynamicParticle(), G4Material::GetIndex(), G4Track::GetMaterial(), G4Material::GetMaterialPropertiesTable(), G4Material::GetName(), G4MaterialPropertiesTable::GetProperty(), and G4DynamicParticle::GetTotalEnergy().
00265 { 00266 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle(); 00267 const G4Material* aMaterial = aTrack.GetMaterial(); 00268 00269 G4double thePhotonEnergy = aParticle->GetTotalEnergy(); 00270 00271 G4double AttenuationLength = DBL_MAX; 00272 00273 if (aMaterial->GetName() == "Water" && DefaultWater){ 00274 00275 G4bool isOutRange; 00276 00277 AttenuationLength = 00278 (*thePhysicsTable)(aMaterial->GetIndex())-> 00279 GetValue(thePhotonEnergy, isOutRange); 00280 } 00281 else { 00282 00283 G4MaterialPropertiesTable* aMaterialPropertyTable = 00284 aMaterial->GetMaterialPropertiesTable(); 00285 00286 if(aMaterialPropertyTable){ 00287 G4MaterialPropertyVector* AttenuationLengthVector = 00288 aMaterialPropertyTable->GetProperty("RAYLEIGH"); 00289 if(AttenuationLengthVector){ 00290 AttenuationLength = AttenuationLengthVector -> 00291 Value(thePhotonEnergy); 00292 } 00293 else{ 00294 // G4cout << "No Rayleigh scattering length specified" << G4endl; 00295 } 00296 } 00297 else{ 00298 // G4cout << "No Rayleigh scattering length specified" << G4endl; 00299 } 00300 } 00301 00302 return AttenuationLength; 00303 }
G4PhysicsTable * G4OpRayleigh::GetPhysicsTable | ( | ) | const [inline] |
Definition at line 176 of file G4OpRayleigh.hh.
References thePhysicsTable.
00177 { 00178 return thePhysicsTable; 00179 }
G4bool G4OpRayleigh::IsApplicable | ( | const G4ParticleDefinition & | aParticleType | ) | [inline, virtual] |
Reimplemented from G4VProcess.
Definition at line 157 of file G4OpRayleigh.hh.
References G4OpticalPhoton::OpticalPhoton().
00158 { 00159 return ( &aParticleType == G4OpticalPhoton::OpticalPhoton() ); 00160 }
G4VParticleChange * G4OpRayleigh::PostStepDoIt | ( | const G4Track & | aTrack, | |
const G4Step & | aStep | |||
) | [virtual] |
Reimplemented from G4VDiscreteProcess.
Definition at line 122 of file G4OpRayleigh.cc.
References G4VProcess::aParticleChange, G4cout, G4endl, G4UniformRand, G4Track::GetDynamicParticle(), G4ParticleChange::GetMomentumDirection(), G4DynamicParticle::GetMomentumDirection(), G4ParticleChange::GetPolarization(), G4DynamicParticle::GetPolarization(), G4ParticleChange::Initialize(), G4VDiscreteProcess::PostStepDoIt(), G4ParticleChange::ProposeMomentumDirection(), G4ParticleChange::ProposePolarization(), and G4VProcess::verboseLevel.
00123 { 00124 aParticleChange.Initialize(aTrack); 00125 00126 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle(); 00127 00128 if (verboseLevel>0) { 00129 G4cout << "Scattering Photon!" << G4endl; 00130 G4cout << "Old Momentum Direction: " 00131 << aParticle->GetMomentumDirection() << G4endl; 00132 G4cout << "Old Polarization: " 00133 << aParticle->GetPolarization() << G4endl; 00134 } 00135 00136 G4double cosTheta; 00137 G4ThreeVector OldMomentumDirection, NewMomentumDirection; 00138 G4ThreeVector OldPolarization, NewPolarization; 00139 00140 do { 00141 // Try to simulate the scattered photon momentum direction 00142 // w.r.t. the initial photon momentum direction 00143 00144 G4double CosTheta = G4UniformRand(); 00145 G4double SinTheta = std::sqrt(1.-CosTheta*CosTheta); 00146 // consider for the angle 90-180 degrees 00147 if (G4UniformRand() < 0.5) CosTheta = -CosTheta; 00148 00149 // simulate the phi angle 00150 G4double rand = twopi*G4UniformRand(); 00151 G4double SinPhi = std::sin(rand); 00152 G4double CosPhi = std::cos(rand); 00153 00154 // start constructing the new momentum direction 00155 G4double unit_x = SinTheta * CosPhi; 00156 G4double unit_y = SinTheta * SinPhi; 00157 G4double unit_z = CosTheta; 00158 NewMomentumDirection.set (unit_x,unit_y,unit_z); 00159 00160 // Rotate the new momentum direction into global reference system 00161 OldMomentumDirection = aParticle->GetMomentumDirection(); 00162 OldMomentumDirection = OldMomentumDirection.unit(); 00163 NewMomentumDirection.rotateUz(OldMomentumDirection); 00164 NewMomentumDirection = NewMomentumDirection.unit(); 00165 00166 // calculate the new polarization direction 00167 // The new polarization needs to be in the same plane as the new 00168 // momentum direction and the old polarization direction 00169 OldPolarization = aParticle->GetPolarization(); 00170 G4double constant = -1./NewMomentumDirection.dot(OldPolarization); 00171 00172 NewPolarization = NewMomentumDirection + constant*OldPolarization; 00173 NewPolarization = NewPolarization.unit(); 00174 00175 // There is a corner case, where the Newmomentum direction 00176 // is the same as oldpolariztion direction: 00177 // random generate the azimuthal angle w.r.t. Newmomentum direction 00178 if (NewPolarization.mag() == 0.) { 00179 rand = G4UniformRand()*twopi; 00180 NewPolarization.set(std::cos(rand),std::sin(rand),0.); 00181 NewPolarization.rotateUz(NewMomentumDirection); 00182 } else { 00183 // There are two directions which are perpendicular 00184 // to the new momentum direction 00185 if (G4UniformRand() < 0.5) NewPolarization = -NewPolarization; 00186 } 00187 00188 // simulate according to the distribution cos^2(theta) 00189 cosTheta = NewPolarization.dot(OldPolarization); 00190 } while (std::pow(cosTheta,2) < G4UniformRand()); 00191 00192 aParticleChange.ProposePolarization(NewPolarization); 00193 aParticleChange.ProposeMomentumDirection(NewMomentumDirection); 00194 00195 if (verboseLevel>0) { 00196 G4cout << "New Polarization: " 00197 << NewPolarization << G4endl; 00198 G4cout << "Polarization Change: " 00199 << *(aParticleChange.GetPolarization()) << G4endl; 00200 G4cout << "New Momentum Direction: " 00201 << NewMomentumDirection << G4endl; 00202 G4cout << "Momentum Change: " 00203 << *(aParticleChange.GetMomentumDirection()) << G4endl; 00204 } 00205 00206 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep); 00207 }
G4PhysicsTable* G4OpRayleigh::thePhysicsTable [protected] |
Definition at line 141 of file G4OpRayleigh.hh.
Referenced by DumpPhysicsTable(), G4OpRayleigh(), GetPhysicsTable(), and ~G4OpRayleigh().