#include <G4SynchrotronRadiationInMat.hh>
Inheritance diagram for G4SynchrotronRadiationInMat:
Definition at line 68 of file G4SynchrotronRadiationInMat.hh.
G4SynchrotronRadiationInMat::G4SynchrotronRadiationInMat | ( | const G4String & | processName = "SynchrotronRadiation" , |
|
G4ProcessType | type = fElectromagnetic | |||
) |
Definition at line 123 of file G4SynchrotronRadiationInMat.cc.
References fSynchrotronRadiation, G4TransportationManager::GetPropagatorInField(), G4TransportationManager::GetTransportationManager(), and G4VProcess::SetProcessSubType().
00124 :G4VDiscreteProcess (processName, type), 00125 LowestKineticEnergy (10.*keV), 00126 HighestKineticEnergy (100.*TeV), 00127 TotBin(200), 00128 theGamma (G4Gamma::Gamma() ), 00129 theElectron ( G4Electron::Electron() ), 00130 thePositron ( G4Positron::Positron() ), 00131 GammaCutInKineticEnergy(0), 00132 ElectronCutInKineticEnergy(0), 00133 PositronCutInKineticEnergy(0), 00134 ParticleCutInKineticEnergy(0), 00135 fAlpha(0.0), fRootNumber(80), 00136 fVerboseLevel( verboseLevel ) 00137 { 00138 G4TransportationManager* transportMgr = G4TransportationManager::GetTransportationManager(); 00139 00140 fFieldPropagator = transportMgr->GetPropagatorInField(); 00141 SetProcessSubType(fSynchrotronRadiation); 00142 CutInRange = GammaCutInKineticEnergyNow = ElectronCutInKineticEnergyNow = 00143 PositronCutInKineticEnergyNow = ParticleCutInKineticEnergyNow = fKsi = 00144 fPsiGamma = fEta = fOrderAngleK = 0.0; 00145 }
G4SynchrotronRadiationInMat::~G4SynchrotronRadiationInMat | ( | ) | [virtual] |
Definition at line 636 of file G4SynchrotronRadiationInMat.cc.
References GetIntegrandForAngleK(), G4Integrator< T, F >::Laguerre(), and CLHEP::detail::n.
Referenced by GetAngleNumberAtGammaKsi().
00637 { 00638 fEta = eta; // should be > 0. ! 00639 G4int n; 00640 G4double result, a; 00641 00642 a = fAlpha; // always = 0. 00643 n = fRootNumber; // around default = 80 00644 00645 G4Integrator<G4SynchrotronRadiationInMat, G4double(G4SynchrotronRadiationInMat::*)(G4double)> integral; 00646 00647 result = integral.Laguerre(this, 00648 &G4SynchrotronRadiationInMat::GetIntegrandForAngleK, a, n); 00649 00650 return result; 00651 }
Definition at line 657 of file G4SynchrotronRadiationInMat.cc.
References GetAngleK().
00658 { 00659 G4double result, funK, funK2, gpsi2 = gpsi*gpsi; 00660 00661 fPsiGamma = gpsi; 00662 fEta = 0.5*fKsi*(1 + gpsi2)*std::sqrt(1 + gpsi2); 00663 00664 fOrderAngleK = 1./3.; 00665 funK = GetAngleK(fEta); 00666 funK2 = funK*funK; 00667 00668 result = gpsi2*funK2/(1 + gpsi2); 00669 00670 fOrderAngleK = 2./3.; 00671 funK = GetAngleK(fEta); 00672 funK2 = funK*funK; 00673 00674 result += funK2; 00675 result *= (1 + gpsi2)*fKsi; 00676 00677 return result; 00678 }
static G4double G4SynchrotronRadiationInMat::GetEnergyConst | ( | ) | [inline, static] |
Definition at line 599 of file G4SynchrotronRadiationInMat.cc.
References GetProbSpectrumSRforEnergy(), G4Integrator< T, F >::Laguerre(), CLHEP::detail::n, and G4INCL::Math::pi.
00600 { 00601 if (ksi <= 0.) return 1.0; 00602 fKsi = ksi; // should be > 0. ! 00603 G4int n; 00604 G4double result, a; 00605 00606 a = fAlpha; // always = 0. 00607 n = fRootNumber; // around default = 80 00608 00609 G4Integrator<G4SynchrotronRadiationInMat, G4double(G4SynchrotronRadiationInMat::*)(G4double)> integral; 00610 00611 result = integral.Laguerre(this, 00612 &G4SynchrotronRadiationInMat::GetProbSpectrumSRforEnergy, a, n); 00613 00614 result *= 9.*std::sqrt(3.)*ksi/8./pi; 00615 00616 return result; 00617 }
Definition at line 623 of file G4SynchrotronRadiationInMat.cc.
Referenced by GetAngleK().
00624 { 00625 G4double result, hypCos=std::cosh(t); 00626 00627 result = std::cosh(fOrderAngleK*t)*std::exp(t - fEta*hypCos); // fEta > 0. ! 00628 result /= hypCos; 00629 return result; 00630 }
Definition at line 560 of file G4SynchrotronRadiationInMat.cc.
References GetProbSpectrumSRforInt(), G4Integrator< T, F >::Laguerre(), CLHEP::detail::n, and G4INCL::Math::pi.
00561 { 00562 if (ksi <= 0.) return 1.0; 00563 fKsi = ksi; // should be > 0. ! 00564 G4int n; 00565 G4double result, a; 00566 00567 a = fAlpha; // always = 0. 00568 n = fRootNumber; // around default = 80 00569 00570 G4Integrator<G4SynchrotronRadiationInMat, G4double(G4SynchrotronRadiationInMat::*)(G4double)> integral; 00571 00572 result = integral.Laguerre(this, 00573 &G4SynchrotronRadiationInMat::GetProbSpectrumSRforInt, a, n); 00574 00575 result *= 3./5./pi; 00576 00577 return result; 00578 }
static G4double G4SynchrotronRadiationInMat::GetLambdaConst | ( | ) | [inline, static] |
G4double G4SynchrotronRadiationInMat::GetMeanFreePath | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4ForceCondition * | condition | |||
) | [virtual] |
Implements G4VDiscreteProcess.
Definition at line 174 of file G4SynchrotronRadiationInMat.cc.
References DBL_MAX, G4PropagatorInField::FindAndSetFieldManager(), G4cout, G4endl, G4DynamicParticle::GetDefinition(), G4FieldManager::GetDetectorField(), G4Track::GetDynamicParticle(), G4Field::GetFieldValue(), G4Track::GetGlobalTime(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGCharge(), G4Track::GetPosition(), G4DynamicParticle::GetTotalEnergy(), G4DynamicParticle::GetTotalMomentum(), G4Track::GetVolume(), and NotForced.
00177 { 00178 // gives the MeanFreePath in GEANT4 internal units 00179 G4double MeanFreePath; 00180 00181 const G4DynamicParticle* aDynamicParticle = trackData.GetDynamicParticle(); 00182 // G4Material* aMaterial = trackData.GetMaterial(); 00183 00184 //G4bool isOutRange ; 00185 00186 *condition = NotForced ; 00187 00188 G4double gamma = aDynamicParticle->GetTotalEnergy()/ 00189 aDynamicParticle->GetMass(); 00190 00191 G4double particleCharge = aDynamicParticle->GetDefinition()->GetPDGCharge(); 00192 00193 G4double KineticEnergy = aDynamicParticle->GetKineticEnergy(); 00194 00195 if ( KineticEnergy < LowestKineticEnergy || gamma < 1.0e3 ) MeanFreePath = DBL_MAX; 00196 else 00197 { 00198 00199 G4ThreeVector FieldValue; 00200 const G4Field* pField = 0; 00201 00202 G4FieldManager* fieldMgr=0; 00203 G4bool fieldExertsForce = false; 00204 00205 if( (particleCharge != 0.0) ) 00206 { 00207 fieldMgr = fFieldPropagator->FindAndSetFieldManager( trackData.GetVolume() ); 00208 00209 if ( fieldMgr != 0 ) 00210 { 00211 // If the field manager has no field, there is no field ! 00212 00213 fieldExertsForce = ( fieldMgr->GetDetectorField() != 0 ); 00214 } 00215 } 00216 if ( fieldExertsForce ) 00217 { 00218 pField = fieldMgr->GetDetectorField() ; 00219 G4ThreeVector globPosition = trackData.GetPosition(); 00220 00221 G4double globPosVec[4], FieldValueVec[6]; 00222 00223 globPosVec[0] = globPosition.x(); 00224 globPosVec[1] = globPosition.y(); 00225 globPosVec[2] = globPosition.z(); 00226 globPosVec[3] = trackData.GetGlobalTime(); 00227 00228 pField->GetFieldValue( globPosVec, FieldValueVec ); 00229 00230 FieldValue = G4ThreeVector( FieldValueVec[0], 00231 FieldValueVec[1], 00232 FieldValueVec[2] ); 00233 00234 00235 00236 G4ThreeVector unitMomentum = aDynamicParticle->GetMomentumDirection(); 00237 G4ThreeVector unitMcrossB = FieldValue.cross(unitMomentum) ; 00238 G4double perpB = unitMcrossB.mag() ; 00239 G4double beta = aDynamicParticle->GetTotalMomentum()/ 00240 (aDynamicParticle->GetTotalEnergy() ); 00241 00242 if( perpB > 0.0 ) MeanFreePath = fLambdaConst*beta/perpB; 00243 else MeanFreePath = DBL_MAX; 00244 } 00245 else MeanFreePath = DBL_MAX; 00246 } 00247 if(fVerboseLevel > 0) 00248 { 00249 G4cout<<"G4SynchrotronRadiationInMat::MeanFreePath = "<<MeanFreePath/m<<" m"<<G4endl; 00250 } 00251 return MeanFreePath; 00252 }
G4double G4SynchrotronRadiationInMat::GetPhotonEnergy | ( | const G4Track & | trackData, | |
const G4Step & | stepData | |||
) |
Definition at line 421 of file G4SynchrotronRadiationInMat.cc.
References G4PropagatorInField::FindAndSetFieldManager(), G4UniformRand, G4DynamicParticle::GetDefinition(), G4FieldManager::GetDetectorField(), G4Track::GetDynamicParticle(), G4Field::GetFieldValue(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGCharge(), G4Track::GetPosition(), G4DynamicParticle::GetTotalEnergy(), and G4Track::GetVolume().
00424 { 00425 G4int i ; 00426 G4double energyOfSR = -1.0 ; 00427 //G4Material* aMaterial=trackData.GetMaterial() ; 00428 00429 const G4DynamicParticle* aDynamicParticle=trackData.GetDynamicParticle(); 00430 00431 G4double gamma = aDynamicParticle->GetTotalEnergy()/ 00432 (aDynamicParticle->GetMass() ) ; 00433 00434 G4double particleCharge = aDynamicParticle->GetDefinition()->GetPDGCharge(); 00435 00436 G4ThreeVector FieldValue; 00437 const G4Field* pField = 0 ; 00438 00439 G4FieldManager* fieldMgr=0; 00440 G4bool fieldExertsForce = false; 00441 00442 if( (particleCharge != 0.0) ) 00443 { 00444 fieldMgr = fFieldPropagator->FindAndSetFieldManager( trackData.GetVolume() ); 00445 if ( fieldMgr != 0 ) 00446 { 00447 // If the field manager has no field, there is no field ! 00448 00449 fieldExertsForce = ( fieldMgr->GetDetectorField() != 0 ); 00450 } 00451 } 00452 if ( fieldExertsForce ) 00453 { 00454 pField = fieldMgr->GetDetectorField(); 00455 G4ThreeVector globPosition = trackData.GetPosition(); 00456 G4double globPosVec[3], FieldValueVec[3]; 00457 00458 globPosVec[0] = globPosition.x(); 00459 globPosVec[1] = globPosition.y(); 00460 globPosVec[2] = globPosition.z(); 00461 00462 pField->GetFieldValue( globPosVec, FieldValueVec ); 00463 FieldValue = G4ThreeVector( FieldValueVec[0], 00464 FieldValueVec[1], 00465 FieldValueVec[2] ); 00466 00467 G4ThreeVector unitMomentum = aDynamicParticle->GetMomentumDirection(); 00468 G4ThreeVector unitMcrossB = FieldValue.cross(unitMomentum) ; 00469 G4double perpB = unitMcrossB.mag(); 00470 if( perpB > 0.0 ) 00471 { 00472 // M-C of synchrotron photon energy 00473 00474 G4double random = G4UniformRand() ; 00475 for(i=0;i<200;i++) 00476 { 00477 if(random >= fIntegralProbabilityOfSR[i]) break ; 00478 } 00479 energyOfSR = 0.0001*i*i*fEnergyConst*gamma*gamma*perpB ; 00480 00481 // check against insufficient energy 00482 00483 if(energyOfSR <= 0.0) 00484 { 00485 return -1.0 ; 00486 } 00487 //G4double kineticEnergy = aDynamicParticle->GetKineticEnergy(); 00488 //G4ParticleMomentum 00489 //particleDirection = aDynamicParticle->GetMomentumDirection(); 00490 00491 // Gamma production cut in this material 00492 //G4double 00493 //gammaEnergyCut = (G4Gamma::GetCutsInEnergy())[aMaterial->GetIndex()]; 00494 00495 // SR photon has energy more than the current material cut 00496 // M-C of its direction 00497 00498 //G4double Teta = G4UniformRand()/gamma ; // Very roughly 00499 00500 //G4double Phi = twopi * G4UniformRand() ; 00501 } 00502 else 00503 { 00504 return -1.0 ; 00505 } 00506 } 00507 return energyOfSR ; 00508 }
Definition at line 584 of file G4SynchrotronRadiationInMat.cc.
Referenced by GetEnergyProbSR().
00585 { 00586 G4double result, hypCos=std::cosh(t); 00587 00588 result = std::cosh(5.*t/3.)*std::exp(t - fKsi*hypCos); // fKsi > 0. ! 00589 result /= hypCos; 00590 return result; 00591 }
Definition at line 544 of file G4SynchrotronRadiationInMat.cc.
Referenced by GetIntProbSR().
00545 { 00546 G4double result, hypCos2, hypCos=std::cosh(t); 00547 00548 hypCos2 = hypCos*hypCos; 00549 result = std::cosh(5.*t/3.)*std::exp(t-fKsi*hypCos); // fKsi > 0. ! 00550 result /= hypCos2; 00551 return result; 00552 }
Definition at line 514 of file G4SynchrotronRadiationInMat.cc.
References G4UniformRand, and position.
Referenced by PostStepDoIt().
00515 { 00516 G4int i, iMax; 00517 G4double energySR, random, position; 00518 00519 iMax = 200; 00520 random = G4UniformRand(); 00521 00522 for( i = 0; i < iMax; i++ ) 00523 { 00524 if( random >= fIntegralProbabilityOfSR[i] ) break; 00525 } 00526 if(i <= 0 ) position = G4UniformRand(); // 0. 00527 else if( i>= iMax) position = G4double(iMax); 00528 else position = i + G4UniformRand(); // -1 00529 // 00530 // it was in initial implementation: 00531 // energyOfSR = 0.0001*i*i*fEnergyConst*gamma*gamma*perpB ; 00532 00533 energySR = 0.0001*position*position*fEnergyConst*gamma*gamma*perpB; 00534 00535 if( energySR < 0. ) energySR = 0.; 00536 00537 return energySR; 00538 }
G4bool G4SynchrotronRadiationInMat::IsApplicable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 157 of file G4SynchrotronRadiationInMat.cc.
00158 { 00159 00160 return ( ( &particle == (const G4ParticleDefinition *)theElectron ) || 00161 ( &particle == (const G4ParticleDefinition *)thePositron ) ); 00162 00163 }
G4VParticleChange * G4SynchrotronRadiationInMat::PostStepDoIt | ( | const G4Track & | track, | |
const G4Step & | Step | |||
) | [virtual] |
Reimplemented from G4VDiscreteProcess.
Definition at line 259 of file G4SynchrotronRadiationInMat.cc.
References G4ParticleChange::AddSecondary(), G4VProcess::aParticleChange, fcos(), G4PropagatorInField::FindAndSetFieldManager(), fStopAndKill, fStopButAlive, G4cout, G4endl, G4UniformRand, G4Gamma::Gamma(), G4DynamicParticle::GetDefinition(), G4FieldManager::GetDetectorField(), G4Track::GetDynamicParticle(), G4Field::GetFieldValue(), G4Track::GetGlobalTime(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGCharge(), G4Track::GetPosition(), GetRandomEnergySR(), G4DynamicParticle::GetTotalEnergy(), G4Track::GetVolume(), G4ParticleChange::Initialize(), G4VDiscreteProcess::PostStepDoIt(), G4ParticleChange::ProposeEnergy(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChange::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), G4VParticleChange::SetNumberOfSecondaries(), and G4DynamicParticle::SetPolarization().
00262 { 00263 aParticleChange.Initialize(trackData); 00264 00265 const G4DynamicParticle* aDynamicParticle=trackData.GetDynamicParticle(); 00266 00267 G4double gamma = aDynamicParticle->GetTotalEnergy()/ 00268 (aDynamicParticle->GetMass() ); 00269 00270 if(gamma <= 1.0e3 ) 00271 { 00272 return G4VDiscreteProcess::PostStepDoIt(trackData,stepData); 00273 } 00274 G4double particleCharge = aDynamicParticle->GetDefinition()->GetPDGCharge(); 00275 00276 G4ThreeVector FieldValue; 00277 const G4Field* pField = 0 ; 00278 00279 G4FieldManager* fieldMgr=0; 00280 G4bool fieldExertsForce = false; 00281 00282 if( (particleCharge != 0.0) ) 00283 { 00284 fieldMgr = fFieldPropagator->FindAndSetFieldManager( trackData.GetVolume() ); 00285 if ( fieldMgr != 0 ) 00286 { 00287 // If the field manager has no field, there is no field ! 00288 00289 fieldExertsForce = ( fieldMgr->GetDetectorField() != 0 ); 00290 } 00291 } 00292 if ( fieldExertsForce ) 00293 { 00294 pField = fieldMgr->GetDetectorField() ; 00295 G4ThreeVector globPosition = trackData.GetPosition() ; 00296 G4double globPosVec[4], FieldValueVec[6] ; 00297 globPosVec[0] = globPosition.x() ; 00298 globPosVec[1] = globPosition.y() ; 00299 globPosVec[2] = globPosition.z() ; 00300 globPosVec[3] = trackData.GetGlobalTime(); 00301 00302 pField->GetFieldValue( globPosVec, FieldValueVec ) ; 00303 FieldValue = G4ThreeVector( FieldValueVec[0], 00304 FieldValueVec[1], 00305 FieldValueVec[2] ); 00306 00307 G4ThreeVector unitMomentum = aDynamicParticle->GetMomentumDirection(); 00308 G4ThreeVector unitMcrossB = FieldValue.cross(unitMomentum); 00309 G4double perpB = unitMcrossB.mag() ; 00310 if(perpB > 0.0) 00311 { 00312 // M-C of synchrotron photon energy 00313 00314 G4double energyOfSR = GetRandomEnergySR(gamma,perpB); 00315 00316 if(fVerboseLevel > 0) 00317 { 00318 G4cout<<"SR photon energy = "<<energyOfSR/keV<<" keV"<<G4endl; 00319 } 00320 // check against insufficient energy 00321 00322 if( energyOfSR <= 0.0 ) 00323 { 00324 return G4VDiscreteProcess::PostStepDoIt(trackData,stepData); 00325 } 00326 G4double kineticEnergy = aDynamicParticle->GetKineticEnergy(); 00327 G4ParticleMomentum 00328 particleDirection = aDynamicParticle->GetMomentumDirection(); 00329 00330 // M-C of its direction, simplified dipole busted approach 00331 00332 // G4double Teta = G4UniformRand()/gamma ; // Very roughly 00333 00334 G4double cosTheta, sinTheta, fcos, beta; 00335 00336 do 00337 { 00338 cosTheta = 1. - 2.*G4UniformRand(); 00339 fcos = (1 + cosTheta*cosTheta)*0.5; 00340 } 00341 while( fcos < G4UniformRand() ); 00342 00343 beta = std::sqrt(1. - 1./(gamma*gamma)); 00344 00345 cosTheta = (cosTheta + beta)/(1. + beta*cosTheta); 00346 00347 if( cosTheta > 1. ) cosTheta = 1.; 00348 if( cosTheta < -1. ) cosTheta = -1.; 00349 00350 sinTheta = std::sqrt(1. - cosTheta*cosTheta ); 00351 00352 G4double Phi = twopi * G4UniformRand() ; 00353 00354 G4double dirx = sinTheta*std::cos(Phi) , 00355 diry = sinTheta*std::sin(Phi) , 00356 dirz = cosTheta; 00357 00358 G4ThreeVector gammaDirection ( dirx, diry, dirz); 00359 gammaDirection.rotateUz(particleDirection); 00360 00361 // polarization of new gamma 00362 00363 // G4double sx = std::cos(Teta)*std::cos(Phi); 00364 // G4double sy = std::cos(Teta)*std::sin(Phi); 00365 // G4double sz = -std::sin(Teta); 00366 00367 G4ThreeVector gammaPolarization = FieldValue.cross(gammaDirection); 00368 gammaPolarization = gammaPolarization.unit(); 00369 00370 // (sx, sy, sz); 00371 // gammaPolarization.rotateUz(particleDirection); 00372 00373 // create G4DynamicParticle object for the SR photon 00374 00375 G4DynamicParticle* aGamma= new G4DynamicParticle ( G4Gamma::Gamma(), 00376 gammaDirection, 00377 energyOfSR ); 00378 aGamma->SetPolarization( gammaPolarization.x(), 00379 gammaPolarization.y(), 00380 gammaPolarization.z() ); 00381 00382 00383 aParticleChange.SetNumberOfSecondaries(1); 00384 aParticleChange.AddSecondary(aGamma); 00385 00386 // Update the incident particle 00387 00388 G4double newKinEnergy = kineticEnergy - energyOfSR ; 00389 00390 if (newKinEnergy > 0.) 00391 { 00392 aParticleChange.ProposeMomentumDirection( particleDirection ); 00393 aParticleChange.ProposeEnergy( newKinEnergy ); 00394 aParticleChange.ProposeLocalEnergyDeposit (0.); 00395 } 00396 else 00397 { 00398 aParticleChange.ProposeEnergy( 0. ); 00399 aParticleChange.ProposeLocalEnergyDeposit (0.); 00400 G4double charge = aDynamicParticle->GetDefinition()->GetPDGCharge(); 00401 if (charge<0.) 00402 { 00403 aParticleChange.ProposeTrackStatus(fStopAndKill) ; 00404 } 00405 else 00406 { 00407 aParticleChange.ProposeTrackStatus(fStopButAlive) ; 00408 } 00409 } 00410 } 00411 else 00412 { 00413 return G4VDiscreteProcess::PostStepDoIt(trackData,stepData); 00414 } 00415 } 00416 return G4VDiscreteProcess::PostStepDoIt(trackData,stepData); 00417 }
void G4SynchrotronRadiationInMat::SetEta | ( | G4double | eta | ) | [inline] |
void G4SynchrotronRadiationInMat::SetKsi | ( | G4double | ksi | ) | [inline] |
void G4SynchrotronRadiationInMat::SetOrderAngleK | ( | G4double | ord | ) | [inline] |
void G4SynchrotronRadiationInMat::SetPsiGamma | ( | G4double | psg | ) | [inline] |
void G4SynchrotronRadiationInMat::SetRootNumber | ( | G4int | rn | ) | [inline] |
void G4SynchrotronRadiationInMat::SetVerboseLevel | ( | G4int | v | ) | [inline] |