#include <G4DNABrownianTransportation.hh>
Inheritance diagram for G4DNABrownianTransportation:
Definition at line 49 of file G4DNABrownianTransportation.hh.
G4DNABrownianTransportation::G4DNABrownianTransportation | ( | const G4String & | aName = "DNABrownianTransportation" , |
|
G4int | verbosityLevel = 1 | |||
) |
Definition at line 83 of file G4DNABrownianTransportation.cc.
References G4NistManager::FindOrBuildMaterial(), fNistWater, fpWaterDensity, fUseMaximumTimeBeforeReachingBoundary, G4NistManager::Instance(), G4VProcess::SetProcessSubType(), and G4VProcess::verboseLevel.
00083 : 00084 G4ITTransportation(aName, verbosity), 00085 InitProcessState(fpBrownianState, fTransportationState) 00086 { 00087 //ctor 00088 SetProcessSubType(61); 00089 verboseLevel = 1; 00090 fUseMaximumTimeBeforeReachingBoundary = true; 00091 fNistWater = G4NistManager::Instance()->FindOrBuildMaterial("G4_WATER"); 00092 fpWaterDensity = 0; 00093 }
G4DNABrownianTransportation::~G4DNABrownianTransportation | ( | ) | [virtual] |
G4DNABrownianTransportation::G4DNABrownianTransportation | ( | const G4DNABrownianTransportation & | other | ) |
Definition at line 98 of file G4DNABrownianTransportation.cc.
References fNistWater, fpWaterDensity, fUseMaximumTimeBeforeReachingBoundary, G4VProcess::SetProcessSubType(), and G4VProcess::verboseLevel.
00098 : 00099 G4ITTransportation(right), 00100 InitProcessState(fpBrownianState, fTransportationState) 00101 { 00102 //copy ctor 00103 SetProcessSubType(61); 00104 verboseLevel = right.verboseLevel; 00105 fUseMaximumTimeBeforeReachingBoundary = right.fUseMaximumTimeBeforeReachingBoundary; 00106 fNistWater = right.fNistWater; 00107 fpWaterDensity = right.fpWaterDensity ; 00108 }
G4VParticleChange * G4DNABrownianTransportation::AlongStepDoIt | ( | const G4Track & | track, | |
const G4Step & | ||||
) | [virtual] |
Reimplemented from G4ITTransportation.
Definition at line 362 of file G4DNABrownianTransportation.cc.
References G4ITTransportation::AlongStepDoIt(), Diffusion(), and G4ITTransportation::fParticleChange.
00364 { 00365 G4ITTransportation::AlongStepDoIt(track,step); 00366 Diffusion(track); 00367 return &fParticleChange; 00368 }
G4double G4DNABrownianTransportation::AlongStepGetPhysicalInteractionLength | ( | const G4Track & | , | |
G4double | , | |||
G4double | , | |||
G4double & | , | |||
G4GPILSelection * | ||||
) | [virtual] |
Reimplemented from G4ITTransportation.
Definition at line 314 of file G4DNABrownianTransportation.cc.
References G4ITTransportation::AlongStepGetPhysicalInteractionLength(), fUseMaximumTimeBeforeReachingBoundary, G4UniformRand, G4Molecule::GetDiffusionCoefficient(), G4Track::GetGlobalTime(), GetMolecule(), and State.
00320 { 00321 G4double geometryStepLength = G4ITTransportation::AlongStepGetPhysicalInteractionLength(track, 00322 previousStepSize, 00323 currentMinimumStep, 00324 currentSafety, 00325 selection); 00326 00327 G4double diffusionCoefficient = GetMolecule(track)->GetDiffusionCoefficient(); 00328 // G4double diffusionCoefficient = GetITBrownianObject(track)->GetDiffusionCoefficient(track.GetMaterial()); 00329 00330 if(State(fGeometryLimitedStep)) 00331 { 00332 // 99 % of the space step distribution is lower than 00333 // d_99 = 8 * sqrt(D*t) 00334 // where t is the corresponding time step 00335 // so by inversion : 00336 if(fUseMaximumTimeBeforeReachingBoundary) 00337 { 00338 State(theInteractionTimeLeft) = (geometryStepLength*geometryStepLength)/(64 * diffusionCoefficient); 00339 } 00340 else // Will use a random time 00341 { 00342 State(theInteractionTimeLeft) = 1/(4*diffusionCoefficient) * pow(geometryStepLength/InvErfc(G4UniformRand()),2); 00343 } 00344 00345 State(fCandidateEndGlobalTime) = track.GetGlobalTime() + State(theInteractionTimeLeft); 00346 State(fPathLengthWasCorrected) = false; 00347 } 00348 else 00349 { 00350 geometryStepLength = 2*sqrt(diffusionCoefficient*State(theInteractionTimeLeft) ) *InvErf(G4UniformRand()); 00351 State(fPathLengthWasCorrected) = true; 00352 State(endpointDistance) = geometryStepLength; 00353 } 00354 00355 return geometryStepLength ; 00356 }
void G4DNABrownianTransportation::BuildPhysicsTable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4ITTransportation.
Definition at line 129 of file G4DNABrownianTransportation.cc.
References fpWaterDensity, G4cout, G4endl, G4DNAMolecularMaterial::GetDensityTableFor(), G4Material::GetMaterial(), G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), G4VProcess::GetProcessSubType(), G4DNAMolecularMaterial::Instance(), and G4VProcess::verboseLevel.
00130 { 00131 if(verboseLevel > 0) 00132 { 00133 G4cout << G4endl << GetProcessName() << ": for " 00134 << setw(24) << particle.GetParticleName() 00135 << "\tSubType= " << GetProcessSubType() << G4endl; 00136 } 00137 00138 // Initialize water density pointer 00139 fpWaterDensity = G4DNAMolecularMaterial::Instance()->GetDensityTableFor(G4Material::GetMaterial("G4_WATER")); 00140 }
void G4DNABrownianTransportation::ComputeStep | ( | const G4Track & | , | |
const G4Step & | , | |||
const | double, | |||
double & | ||||
) | [virtual] |
Reimplemented from G4ITTransportation.
Definition at line 142 of file G4DNABrownianTransportation.cc.
References DBL_MAX, FatalErrorInArgument, G4ITTransportation::fVerboseLevel, G4BestUnit, G4cout, G4endl, G4Exception(), G4Molecule::GetDiffusionCoefficient(), G4StepPoint::GetGlobalTime(), GetIT(), GetMolecule(), G4StepPoint::GetMomentumDirection(), G4Track::GetPosition(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4StepPoint::GetProcessDefinedStep(), G4Track::GetTrackID(), GREEN_ON_BLUE, G4VITProcess::ProposesTimeStep(), RESET, and State.
00146 { 00147 // G4cout << "G4ITBrownianTransportation::ComputeStep" << G4endl; 00148 00149 // If this method is called, this step 00150 // cannot be geometry limited. 00151 // In case the step is limited by the geometry, 00152 // this method should not be called. 00153 // The fTransportEndPosition calculated in 00154 // the method AlongStepIL should be taken 00155 // into account. 00156 // In order to do so, the flag IsLeadingStep 00157 // is on. Meaning : this track has the minimum 00158 // interaction length over all others. 00159 if(GetIT(track)->GetTrackingInfo()->IsLeadingStep()) 00160 { 00161 const G4VITProcess* ITProc = ((const G4VITProcess*) step.GetPostStepPoint()->GetProcessDefinedStep()); 00162 bool makeException = true; 00163 00164 if(ITProc && ITProc->ProposesTimeStep()) makeException = false; 00165 00166 if(makeException) 00167 { 00168 00169 G4ExceptionDescription exceptionDescription ; 00170 exceptionDescription << "ComputeStep is called while the track has the minimum interaction time"; 00171 exceptionDescription << " so it should not recompute a timeStep "; 00172 G4Exception("G4DNABrownianTransportation::ComputeStep","G4DNABrownianTransportation001", 00173 FatalErrorInArgument,exceptionDescription); 00174 } 00175 } 00176 00177 State(fGeometryLimitedStep) = false; 00178 // TODO : generalize this process to all kind of brownian objects 00179 // G4ITBrownianObject* ITBrown = GetITBrownianObject(track) ; 00180 // G4double diffCoeff = ITBrown->GetDiffusionCoefficient(track.GetMaterial()); 00181 G4Molecule* molecule = GetMolecule(track) ; 00182 G4double diffCoeff = molecule->GetDiffusionCoefficient(); 00183 00184 if(timeStep > 0) 00185 { 00186 spaceStep = DBL_MAX; 00187 00188 while(spaceStep > State(endpointDistance)) 00189 // Probably inefficient when the track is close close to boundaries 00190 // it goes with fUserMaximumTimeBeforeReachingBoundary == false 00191 // fUserMaximumTimeBeforeReachingBoundary == true, it should never loop 00192 { 00193 G4double x = G4RandGauss::shoot(0,sqrt(2*diffCoeff*timeStep)); 00194 G4double y = G4RandGauss::shoot(0,sqrt(2*diffCoeff*timeStep)); 00195 G4double z = G4RandGauss::shoot(0,sqrt(2*diffCoeff*timeStep)); 00196 00197 spaceStep = sqrt(x*x + y*y + z*z); 00198 } 00199 // State(fTransportEndPosition).set(x + track.GetPosition().x(), 00200 // y + track.GetPosition().y(), 00201 // z + track.GetPosition().z()); 00202 00203 State(fTransportEndPosition)= spaceStep*step.GetPostStepPoint()->GetMomentumDirection() + track.GetPosition(); 00204 } 00205 else 00206 { 00207 spaceStep = 0. ; 00208 State(fTransportEndPosition) = track.GetPosition() ; 00209 } 00210 00211 State(fCandidateEndGlobalTime) = step.GetPreStepPoint()->GetGlobalTime() + timeStep ; 00212 State(fEndGlobalTimeComputed) = true ; 00213 00214 #ifdef G4VERBOSE 00215 // DEBUG 00216 if(fVerboseLevel>1) 00217 { 00218 G4cout<< GREEN_ON_BLUE 00219 << "G4ITBrownianTransportation::ComputeStep() : " 00220 << " trackID : " << track.GetTrackID() 00221 << " : Molecule name: " << molecule-> GetName() 00222 << G4endl 00223 << "Diffusion length : " << G4BestUnit(spaceStep, "Length") 00224 << " within time step : " << G4BestUnit(timeStep,"Time") 00225 << RESET 00226 << G4endl<< G4endl; 00227 } 00228 #endif 00229 }
void G4DNABrownianTransportation::Diffusion | ( | const G4Track & | track | ) | [protected] |
Definition at line 253 of file G4DNABrownianTransportation.cc.
References G4ITTransportation::fParticleChange, fStopButAlive, G4ITTransportation::fVerboseLevel, G4BestUnit, G4cout, G4endl, G4UniformRand, G4Track::GetCurrentStepNumber(), G4Track::GetGlobalTime(), GetIT(), G4Track::GetLocalTime(), G4Track::GetMaterial(), GetMolecule(), G4Molecule::GetName(), G4IT::GetName(), G4Track::GetTrackID(), GREEN_ON_BLUE, G4INCL::Math::pi, G4ParticleChange::ProposeEnergy(), G4ParticleChange::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), RESET, G4ParticleChangeForTransport::SetMomentumChanged(), and State.
Referenced by AlongStepDoIt().
00255 { 00256 00257 #ifdef G4VERBOSE 00258 // DEBUG 00259 if (fVerboseLevel>1) 00260 { 00261 G4cout<< GREEN_ON_BLUE 00262 << setw(18)<< "G4DNABrownianTransportation::Diffusion :" 00263 << setw(8) << GetIT(track)->GetName() 00264 << "\t trackID:" << track.GetTrackID() <<"\t" 00265 << " Global Time = " << G4BestUnit(track.GetGlobalTime(),"Time") 00266 << RESET 00267 << G4endl<< G4endl; 00268 } 00269 #endif 00270 00271 G4Material* material = track.GetMaterial(); 00272 // if (material != fNistWater && material->GetBaseMaterial() != fNistWater) 00273 00274 G4double waterDensity = (*fpWaterDensity)[material->GetIndex()]; 00275 00276 if(waterDensity == 0.0) 00277 // if (material == nistwater || material->GetBaseMaterial() == nistwater) 00278 { 00279 G4cout << "A track is outside water material : trackID"<< track.GetTrackID() << " (" << GetMolecule(track)->GetName() <<")" << G4endl; 00280 G4cout << "Local Time : "<< (track.GetLocalTime()) /s<<G4endl; 00281 G4cout<< "Step Number :" << track.GetCurrentStepNumber() <<G4endl; 00282 00283 fParticleChange.ProposeEnergy(0.) ; 00284 fParticleChange.ProposeTrackStatus(fStopButAlive); 00285 00286 // Got pb with : 00287 // fParticleChange.ProposeTrackStatus(fStopAndKill); 00288 // It makes the tracks going straight without killing them 00289 00290 return ; // &fParticleChange is the final returned object 00291 } 00292 00293 G4double costheta = (2*G4UniformRand()-1); 00294 G4double theta = acos (costheta); 00295 G4double phi = 2*pi*G4UniformRand(); 00296 00297 G4double xMomentum = cos(phi)* sin(theta); 00298 G4double yMomentum = sin(theta)*sin(phi); 00299 G4double zMomentum = costheta; 00300 00301 fParticleChange.ProposeMomentumDirection(xMomentum, yMomentum, zMomentum); 00302 State(fMomentumChanged) = true; 00303 fParticleChange.SetMomentumChanged(true) ; 00304 00305 // G4cout << "BROWN : Propose new direction :" << G4ThreeVector(xMomentum, yMomentum, zMomentum) << G4endl; 00306 00307 // Alternative 00308 //fParticleChange.ProposeMomentumDirection(G4RandomDirection()); 00309 00310 return; // &fParticleChange is the final returned object 00311 }
G4DNABrownianTransportation & G4DNABrownianTransportation::operator= | ( | const G4DNABrownianTransportation & | other | ) |
Definition at line 110 of file G4DNABrownianTransportation.cc.
00111 { 00112 if (this == &rhs) return *this; // handle self assignment 00113 //assignment operator 00114 return *this; 00115 }
G4VParticleChange * G4DNABrownianTransportation::PostStepDoIt | ( | const G4Track & | track, | |
const G4Step & | ||||
) | [virtual] |
Reimplemented from G4ITTransportation.
Definition at line 231 of file G4DNABrownianTransportation.cc.
References G4ITTransportation::fParticleChange, G4ITTransportation::fVerboseLevel, G4BestUnit, G4cout, G4endl, G4Step::GetDeltaTime(), G4Track::GetGlobalTime(), GetMolecule(), G4Step::GetStepLength(), G4Track::GetTrackID(), GREEN_ON_BLUE, G4ITTransportation::PostStepDoIt(), and RESET.
00232 { 00233 G4ITTransportation::PostStepDoIt(track,step); 00234 00235 #ifdef G4VERBOSE 00236 // DEBUG 00237 if(fVerboseLevel>1) 00238 { 00239 G4cout<< GREEN_ON_BLUE 00240 << "G4ITBrownianTransportation::PostStepDoIt() :" 00241 << " trackID : " << track.GetTrackID() 00242 << " Molecule name: " << GetMolecule(track)-> GetName() << G4endl; 00243 G4cout<< "Diffusion length : "<< G4BestUnit(step.GetStepLength(),"Length") <<" within time step : " 00244 << G4BestUnit(step.GetDeltaTime(),"Time") << "\t" 00245 << " Current global time : " << G4BestUnit(track.GetGlobalTime(),"Time") 00246 << RESET 00247 << G4endl<< G4endl; 00248 } 00249 #endif 00250 return &fParticleChange ; 00251 }
void G4DNABrownianTransportation::StartTracking | ( | G4Track * | aTrack | ) | [virtual] |
Reimplemented from G4ITTransportation.
Definition at line 122 of file G4DNABrownianTransportation.cc.
References G4VITProcess::fpState, G4ITTransportation::SetInstantiateProcessState(), and G4ITTransportation::StartTracking().
00123 { 00124 fpState = new G4ITBrownianState(); 00125 SetInstantiateProcessState(false); 00126 G4ITTransportation::StartTracking(track); 00127 }
G4Material* G4DNABrownianTransportation::fNistWater [protected] |
Definition at line 92 of file G4DNABrownianTransportation.hh.
Referenced by G4DNABrownianTransportation().
G4ITBrownianState* const& G4DNABrownianTransportation::fpBrownianState [protected] |
Definition at line 89 of file G4DNABrownianTransportation.hh.
const std::vector<G4double>* G4DNABrownianTransportation::fpWaterDensity [protected] |
Definition at line 95 of file G4DNABrownianTransportation.hh.
Referenced by BuildPhysicsTable(), and G4DNABrownianTransportation().
Definition at line 91 of file G4DNABrownianTransportation.hh.
Referenced by AlongStepGetPhysicalInteractionLength(), and G4DNABrownianTransportation().