#include <G4eBremsstrahlungModel.hh>
Inheritance diagram for G4eBremsstrahlungModel:
Definition at line 64 of file G4eBremsstrahlungModel.hh.
G4eBremsstrahlungModel::G4eBremsstrahlungModel | ( | const G4ParticleDefinition * | p = 0 , |
|
const G4String & | nam = "eBrem" | |||
) |
Definition at line 85 of file G4eBremsstrahlungModel.cc.
References fParticleChange, G4Gamma::Gamma(), G4VEmModel::HighEnergyLimit(), G4VEmModel::LowEnergyLimit(), G4VEmModel::SetAngularDistribution(), and theGamma.
00087 : G4VEmModel(nam), 00088 particle(0), 00089 isElectron(true), 00090 probsup(1.0), 00091 MigdalConstant(classic_electr_radius*electron_Compton_length*electron_Compton_length*4.0*pi), 00092 LPMconstant(fine_structure_const*electron_mass_c2*electron_mass_c2/(4.*pi*hbarc)), 00093 isInitialised(false) 00094 { 00095 if(p) { SetParticle(p); } 00096 theGamma = G4Gamma::Gamma(); 00097 SetAngularDistribution(new G4ModifiedTsai()); 00098 highKinEnergy = HighEnergyLimit(); 00099 lowKinEnergy = LowEnergyLimit(); 00100 fParticleChange = 0; 00101 }
G4eBremsstrahlungModel::~G4eBremsstrahlungModel | ( | ) | [virtual] |
Definition at line 105 of file G4eBremsstrahlungModel.cc.
References CLHEP::detail::n.
00106 { 00107 size_t n = partialSumSigma.size(); 00108 if(n > 0) { 00109 for(size_t i=0; i<n; i++) { 00110 delete partialSumSigma[i]; 00111 } 00112 } 00113 }
G4double G4eBremsstrahlungModel::ComputeCrossSectionPerAtom | ( | const G4ParticleDefinition * | , | |
G4double | tkin, | |||
G4double | Z, | |||
G4double | , | |||
G4double | cut, | |||
G4double | maxE = DBL_MAX | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 484 of file G4eBremsstrahlungModel.cc.
References isElectron.
Referenced by CrossSectionPerVolume().
00493 { 00494 G4double cross = 0.0 ; 00495 if ( kineticEnergy < keV || kineticEnergy < cut) { return cross; } 00496 00497 static const G4double ksi=2.0, alfa=1.00; 00498 static const G4double csigh = 0.127, csiglow = 0.25, asiglow = 0.020*MeV ; 00499 static const G4double Tlim = 10.*MeV ; 00500 00501 static const G4double xlim = 1.2 ; 00502 static const G4int NZ = 8 ; 00503 static const G4int Nsig = 11 ; 00504 static const G4double ZZ[NZ] = 00505 {2.,4.,6.,14.,26.,50.,82.,92.} ; 00506 static const G4double coefsig[NZ][Nsig] = { 00507 // Z=2 00508 { 0.4638, 0.37748, 0.32249, -0.060362, -0.065004, 00509 -0.033457, -0.004583, 0.011954, 0.0030404, -0.0010077, 00510 -0.00028131}, 00511 00512 // Z=4 00513 { 0.50008, 0.33483, 0.34364, -0.086262, -0.055361, 00514 -0.028168, -0.0056172, 0.011129, 0.0027528, -0.00092265, 00515 -0.00024348}, 00516 00517 // Z=6 00518 { 0.51587, 0.31095, 0.34996, -0.11623, -0.056167, 00519 -0.0087154, 0.00053943, 0.0054092, 0.00077685, -0.00039635, 00520 -6.7818e-05}, 00521 00522 // Z=14 00523 { 0.55058, 0.25629, 0.35854, -0.080656, -0.054308, 00524 -0.049933, -0.00064246, 0.016597, 0.0021789, -0.001327, 00525 -0.00025983}, 00526 00527 // Z=26 00528 { 0.5791, 0.26152, 0.38953, -0.17104, -0.099172, 00529 0.024596, 0.023718, -0.0039205, -0.0036658, 0.00041749, 00530 0.00023408}, 00531 00532 // Z=50 00533 { 0.62085, 0.27045, 0.39073, -0.37916, -0.18878, 00534 0.23905, 0.095028, -0.068744, -0.023809, 0.0062408, 00535 0.0020407}, 00536 00537 // Z=82 00538 { 0.66053, 0.24513, 0.35404, -0.47275, -0.22837, 00539 0.35647, 0.13203, -0.1049, -0.034851, 0.0095046, 00540 0.0030535}, 00541 00542 // Z=92 00543 { 0.67143, 0.23079, 0.32256, -0.46248, -0.20013, 00544 0.3506, 0.11779, -0.1024, -0.032013, 0.0092279, 00545 0.0028592} 00546 00547 } ; 00548 00549 G4int iz = 0 ; 00550 G4double delz = 1.e6 ; 00551 for (G4int ii=0; ii<NZ; ii++) 00552 { 00553 G4double absdelz = std::abs(Z-ZZ[ii]); 00554 if(absdelz < delz) 00555 { 00556 iz = ii ; 00557 delz = absdelz; 00558 } 00559 } 00560 00561 G4double xx = log10(kineticEnergy/MeV); 00562 G4double fs = 1. ; 00563 00564 if (xx <= xlim) { 00565 00566 fs = coefsig[iz][Nsig-1] ; 00567 for (G4int j=Nsig-2; j>=0; j--) { 00568 00569 fs = fs*xx+coefsig[iz][j] ; 00570 } 00571 if(fs < 0.) { fs = 0.; } 00572 } 00573 00574 cross = Z*(Z+ksi)*(1.-csigh*exp(log(Z)/4.))*pow(log(kineticEnergy/cut),alfa); 00575 00576 if (kineticEnergy <= Tlim) 00577 cross *= exp(csiglow*log(Tlim/kineticEnergy)) 00578 *(1.+asiglow/(sqrt(Z)*kineticEnergy)); 00579 00580 if (!isElectron) 00581 cross *= PositronCorrFactorSigma(Z, kineticEnergy, cut); 00582 00583 cross *= fs/Avogadro ; 00584 00585 if (cross < 0.) { cross = 0.; } 00586 return cross; 00587 }
G4double G4eBremsstrahlungModel::ComputeDEDXPerVolume | ( | const G4Material * | , | |
const G4ParticleDefinition * | , | |||
G4double | kineticEnergy, | |||
G4double | cutEnergy | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 166 of file G4eBremsstrahlungModel.cc.
References G4Material::GetAtomicNumDensityVector(), G4Material::GetElectronDensity(), G4Material::GetElementVector(), G4Material::GetNumberOfElements(), isElectron, CLHEP::detail::n, G4InuclParticleNames::nn, and particle.
00171 { 00172 if(!particle) { SetParticle(p); } 00173 if(kineticEnergy < lowKinEnergy) { return 0.0; } 00174 00175 const G4double thigh = 100.*GeV; 00176 00177 G4double cut = std::min(cutEnergy, kineticEnergy); 00178 00179 G4double rate, loss; 00180 const G4double factorHigh = 36./(1450.*GeV); 00181 const G4double coef1 = -0.5; 00182 const G4double coef2 = 2./9.; 00183 00184 const G4ElementVector* theElementVector = material->GetElementVector(); 00185 const G4double* theAtomicNumDensityVector = material->GetAtomicNumDensityVector(); 00186 00187 G4double totalEnergy = kineticEnergy + electron_mass_c2; 00188 G4double dedx = 0.0; 00189 00190 // loop for elements in the material 00191 for (size_t i=0; i<material->GetNumberOfElements(); i++) { 00192 00193 G4double Z = (*theElementVector)[i]->GetZ(); 00194 G4double natom = theAtomicNumDensityVector[i]; 00195 00196 // loss for MinKinEnergy<KineticEnergy<=100 GeV 00197 if (kineticEnergy <= thigh) { 00198 00199 // x = log(totalEnergy/electron_mass_c2); 00200 loss = ComputeBremLoss(Z, kineticEnergy, cut) ; 00201 if (!isElectron) loss *= PositronCorrFactorLoss(Z, kineticEnergy, cut); 00202 00203 // extrapolation for KineticEnergy>100 GeV 00204 } else { 00205 00206 // G4double xhigh = log(thigh/electron_mass_c2); 00207 G4double cuthigh = thigh*0.5; 00208 00209 if (cut < thigh) { 00210 00211 loss = ComputeBremLoss(Z, thigh, cut) ; 00212 if (!isElectron) loss *= PositronCorrFactorLoss(Z, thigh, cut) ; 00213 rate = cut/totalEnergy; 00214 loss *= (1. + coef1*rate + coef2*rate*rate); 00215 rate = cut/thigh; 00216 loss /= (1.+coef1*rate+coef2*rate*rate); 00217 00218 } else { 00219 00220 loss = ComputeBremLoss(Z, thigh, cuthigh) ; 00221 if (!isElectron) loss *= PositronCorrFactorLoss(Z, thigh, cuthigh) ; 00222 rate = cut/totalEnergy; 00223 loss *= (1. + coef1*rate + coef2*rate*rate); 00224 loss *= cut*factorHigh; 00225 } 00226 } 00227 loss *= natom; 00228 00229 G4double kp2 = MigdalConstant*totalEnergy*totalEnergy 00230 * (material->GetElectronDensity()) ; 00231 00232 // now compute the correction due to the supression(s) 00233 G4double kmin = 1.*eV; 00234 G4double kmax = cut; 00235 00236 if (kmax > kmin) { 00237 00238 G4double floss = 0.; 00239 G4int nmax = 100; 00240 00241 G4double vmin=log(kmin); 00242 G4double vmax=log(kmax) ; 00243 G4int nn = (G4int)(nmax*(vmax-vmin)/(log(highKinEnergy)-vmin)) ; 00244 G4double u,fac,c,v,dv ; 00245 if(nn > 0) { 00246 00247 dv = (vmax-vmin)/nn ; 00248 v = vmin-dv ; 00249 00250 for(G4int n=0; n<=nn; n++) { 00251 00252 v += dv; 00253 u = exp(v); 00254 fac = u*SupressionFunction(material,kineticEnergy,u); 00255 fac *= probsup*(u*u/(u*u+kp2))+1.-probsup; 00256 if ((n==0)||(n==nn)) c=0.5; 00257 else c=1. ; 00258 fac *= c ; 00259 floss += fac ; 00260 } 00261 floss *=dv/(kmax-kmin); 00262 00263 } else { 00264 floss = 1.; 00265 } 00266 if(floss > 1.) floss = 1.; 00267 // correct the loss 00268 loss *= floss; 00269 } 00270 dedx += loss; 00271 } 00272 if(dedx < 0.) { dedx = 0.; } 00273 return dedx; 00274 }
G4double G4eBremsstrahlungModel::CrossSectionPerVolume | ( | const G4Material * | , | |
const G4ParticleDefinition * | , | |||
G4double | kineticEnergy, | |||
G4double | cutEnergy, | |||
G4double | maxEnergy | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 406 of file G4eBremsstrahlungModel.cc.
References ComputeCrossSectionPerAtom(), G4Material::GetAtomicNumDensityVector(), G4Material::GetElectronDensity(), G4Material::GetElementVector(), G4Material::GetNumberOfElements(), CLHEP::detail::n, G4InuclParticleNames::nn, and particle.
00412 { 00413 if(!particle) { SetParticle(p); } 00414 G4double cross = 0.0; 00415 G4double tmax = min(maxEnergy, kineticEnergy); 00416 G4double cut = min(cutEnergy, kineticEnergy); 00417 if(cut >= tmax) { return cross; } 00418 00419 const G4ElementVector* theElementVector = material->GetElementVector(); 00420 const G4double* theAtomNumDensityVector = material->GetAtomicNumDensityVector(); 00421 G4double dum=0.; 00422 00423 for (size_t i=0; i<material->GetNumberOfElements(); i++) { 00424 00425 cross += theAtomNumDensityVector[i] * ComputeCrossSectionPerAtom(p, 00426 kineticEnergy, (*theElementVector)[i]->GetZ(), dum, cut); 00427 if (tmax < kineticEnergy) { 00428 cross -= theAtomNumDensityVector[i] * ComputeCrossSectionPerAtom(p, 00429 kineticEnergy, (*theElementVector)[i]->GetZ(), dum, tmax); 00430 } 00431 } 00432 00433 // now compute the correction due to the supression(s) 00434 00435 G4double kmax = tmax; 00436 G4double kmin = cut; 00437 00438 G4double totalEnergy = kineticEnergy+electron_mass_c2 ; 00439 G4double kp2 = MigdalConstant*totalEnergy*totalEnergy 00440 *(material->GetElectronDensity()); 00441 00442 G4double fsig = 0.; 00443 G4int nmax = 100; 00444 G4double vmin=log(kmin); 00445 G4double vmax=log(kmax) ; 00446 G4int nn = (G4int)(nmax*(vmax-vmin)/(log(highKinEnergy)-vmin)); 00447 G4double u,fac,c,v,dv,y ; 00448 if(nn > 0) { 00449 00450 dv = (vmax-vmin)/nn ; 00451 v = vmin-dv ; 00452 for(G4int n=0; n<=nn; n++) { 00453 00454 v += dv; 00455 u = exp(v); 00456 fac = SupressionFunction(material, kineticEnergy, u); 00457 y = u/kmax; 00458 fac *= (4.-4.*y+3.*y*y)/3.; 00459 fac *= probsup*(u*u/(u*u+kp2))+1.-probsup; 00460 00461 if ((n==0)||(n==nn)) c=0.5; 00462 else c=1. ; 00463 00464 fac *= c; 00465 fsig += fac; 00466 } 00467 y = kmin/kmax ; 00468 fsig *=dv/(-4.*log(y)/3.-4.*(1.-y)/3.+0.5*(1.-y*y)); 00469 00470 } else { 00471 00472 fsig = 1.; 00473 } 00474 if (fsig > 1.) fsig = 1.; 00475 00476 // correct the cross section 00477 cross *= fsig; 00478 00479 return cross; 00480 }
void G4eBremsstrahlungModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Reimplemented in G4ePolarizedBremsstrahlungModel.
Definition at line 126 of file G4eBremsstrahlungModel.cc.
References DBL_MAX, fParticleChange, G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4VEmModel::GetParticleChangeForLoss(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4VEmModel::HighEnergyLimit(), G4VEmModel::LowEnergyLimit(), and G4InuclParticleNames::nn.
Referenced by G4ePolarizedBremsstrahlungModel::Initialise().
00128 { 00129 if(p) { SetParticle(p); } 00130 highKinEnergy = HighEnergyLimit(); 00131 lowKinEnergy = LowEnergyLimit(); 00132 const G4ProductionCutsTable* theCoupleTable= 00133 G4ProductionCutsTable::GetProductionCutsTable(); 00134 00135 if(theCoupleTable) { 00136 G4int numOfCouples = theCoupleTable->GetTableSize(); 00137 00138 G4int nn = partialSumSigma.size(); 00139 G4int nc = cuts.size(); 00140 if(nn > 0) { 00141 for (G4int ii=0; ii<nn; ii++){ 00142 G4DataVector* a=partialSumSigma[ii]; 00143 if ( a ) delete a; 00144 } 00145 partialSumSigma.clear(); 00146 } 00147 if(numOfCouples>0) { 00148 for (G4int i=0; i<numOfCouples; i++) { 00149 G4double cute = DBL_MAX; 00150 if(i < nc) cute = cuts[i]; 00151 const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(i); 00152 const G4Material* material = couple->GetMaterial(); 00153 G4DataVector* dv = ComputePartialSumSigma(material, 0.5*highKinEnergy, 00154 std::min(cute, 0.25*highKinEnergy)); 00155 partialSumSigma.push_back(dv); 00156 } 00157 } 00158 } 00159 if(isInitialised) return; 00160 fParticleChange = GetParticleChangeForLoss(); 00161 isInitialised = true; 00162 }
void G4eBremsstrahlungModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Reimplemented in G4ePolarizedBremsstrahlungModel.
Definition at line 643 of file G4eBremsstrahlungModel.cc.
References fParticleChange, fStopAndKill, G4lrint(), G4UniformRand, G4VEmModel::GetAngularDistribution(), G4Material::GetElectronDensity(), G4Element::GetIonisation(), G4DynamicParticle::GetKineticEnergy(), G4IonisParamElm::GetlogZ3(), G4MaterialCutsCouple::GetMaterial(), G4DynamicParticle::GetMomentumDirection(), G4Element::GetZ(), G4IonisParamElm::GetZ3(), G4IonisParamElm::GetZZ3(), G4VEmModel::LPMFlag(), G4VParticleChange::ProposeTrackStatus(), G4VEmAngularDistribution::SampleDirection(), G4VEmModel::SecondaryThreshold(), SelectRandomAtom(), G4ParticleChangeForLoss::SetProposedKineticEnergy(), G4ParticleChangeForLoss::SetProposedMomentumDirection(), and theGamma.
Referenced by G4ePolarizedBremsstrahlungModel::SampleSecondaries().
00650 : 00651 // cross-section values of Seltzer and Berger for electron energies 00652 // 1 keV - 10 GeV, 00653 // screened Bethe Heilter differential cross section above 10 GeV, 00654 // Migdal corrections in both case. 00655 // Seltzer & Berger: Nim B 12:95 (1985) 00656 // Nelson, Hirayama & Rogers: Technical report 265 SLAC (1985) 00657 // Migdal: Phys Rev 103:1811 (1956); Messel & Crawford: Pergamon Press (1970) 00658 // 00659 // A modified version of the random number techniques of Butcher&Messel is used 00660 // (Nuc Phys 20(1960),15). 00661 { 00662 G4double kineticEnergy = dp->GetKineticEnergy(); 00663 G4double tmax = min(maxEnergy, kineticEnergy); 00664 if(tmin >= tmax) { return; } 00665 00666 // 00667 // GEANT4 internal units. 00668 // 00669 static const G4double 00670 ah10 = 4.67733E+00, ah11 =-6.19012E-01, ah12 = 2.02225E-02, 00671 ah20 =-7.34101E+00, ah21 = 1.00462E+00, ah22 =-3.20985E-02, 00672 ah30 = 2.93119E+00, ah31 =-4.03761E-01, ah32 = 1.25153E-02; 00673 00674 static const G4double 00675 bh10 = 4.23071E+00, bh11 =-6.10995E-01, bh12 = 1.95531E-02, 00676 bh20 =-7.12527E+00, bh21 = 9.69160E-01, bh22 =-2.74255E-02, 00677 bh30 = 2.69925E+00, bh31 =-3.63283E-01, bh32 = 9.55316E-03; 00678 00679 static const G4double 00680 al00 =-2.05398E+00, al01 = 2.38815E-02, al02 = 5.25483E-04, 00681 al10 =-7.69748E-02, al11 =-6.91499E-02, al12 = 2.22453E-03, 00682 al20 = 4.06463E-02, al21 =-1.01281E-02, al22 = 3.40919E-04; 00683 00684 static const G4double 00685 bl00 = 1.04133E+00, bl01 =-9.43291E-03, bl02 =-4.54758E-04, 00686 bl10 = 1.19253E-01, bl11 = 4.07467E-02, bl12 =-1.30718E-03, 00687 bl20 =-1.59391E-02, bl21 = 7.27752E-03, bl22 =-1.94405E-04; 00688 00689 static const G4double tlow = 1.*MeV; 00690 00691 G4double gammaEnergy; 00692 G4bool LPMOK = false; 00693 const G4Material* material = couple->GetMaterial(); 00694 00695 // select randomly one element constituing the material 00696 const G4Element* anElement = SelectRandomAtom(couple); 00697 00698 // Extract Z factors for this Element 00699 G4double lnZ = 3.*(anElement->GetIonisation()->GetlogZ3()); 00700 G4double FZ = lnZ* (4.- 0.55*lnZ); 00701 G4double ZZ = anElement->GetIonisation()->GetZZ3(); 00702 00703 // limits of the energy sampling 00704 G4double totalEnergy = kineticEnergy + electron_mass_c2; 00705 00706 G4double xmin = tmin/kineticEnergy; 00707 G4double xmax = tmax/kineticEnergy; 00708 G4double kappa = 0.0; 00709 if(xmax >= 1.) { xmax = 1.; } 00710 else { kappa = log(xmax)/log(xmin); } 00711 G4double epsilmin = tmin/totalEnergy; 00712 G4double epsilmax = tmax/totalEnergy; 00713 00714 // Migdal factor 00715 G4double MigdalFactor = (material->GetElectronDensity())*MigdalConstant 00716 / (epsilmax*epsilmax); 00717 00718 G4double x, epsil, greject, migdal, grejmax, q; 00719 G4double U = log(kineticEnergy/electron_mass_c2); 00720 G4double U2 = U*U; 00721 00722 // precalculated parameters 00723 G4double ah, bh; 00724 G4double screenfac = 0.0; 00725 00726 if (kineticEnergy > tlow) { 00727 00728 G4double ah1 = ah10 + ZZ* (ah11 + ZZ* ah12); 00729 G4double ah2 = ah20 + ZZ* (ah21 + ZZ* ah22); 00730 G4double ah3 = ah30 + ZZ* (ah31 + ZZ* ah32); 00731 00732 G4double bh1 = bh10 + ZZ* (bh11 + ZZ* bh12); 00733 G4double bh2 = bh20 + ZZ* (bh21 + ZZ* bh22); 00734 G4double bh3 = bh30 + ZZ* (bh31 + ZZ* bh32); 00735 00736 ah = 1. + (ah1*U2 + ah2*U + ah3) / (U2*U); 00737 bh = 0.75 + (bh1*U2 + bh2*U + bh3) / (U2*U); 00738 00739 // limit of the screening variable 00740 screenfac = 00741 136.*electron_mass_c2/((anElement->GetIonisation()->GetZ3())*totalEnergy); 00742 G4double screenmin = screenfac*epsilmin/(1.-epsilmin); 00743 00744 // Compute the maximum of the rejection function 00745 G4double F1 = max(ScreenFunction1(screenmin) - FZ ,0.); 00746 G4double F2 = max(ScreenFunction2(screenmin) - FZ ,0.); 00747 grejmax = (F1 - epsilmin* (F1*ah - bh*epsilmin*F2))/(42.392 - FZ); 00748 00749 } else { 00750 00751 G4double al0 = al00 + ZZ* (al01 + ZZ* al02); 00752 G4double al1 = al10 + ZZ* (al11 + ZZ* al12); 00753 G4double al2 = al20 + ZZ* (al21 + ZZ* al22); 00754 00755 G4double bl0 = bl00 + ZZ* (bl01 + ZZ* bl02); 00756 G4double bl1 = bl10 + ZZ* (bl11 + ZZ* bl12); 00757 G4double bl2 = bl20 + ZZ* (bl21 + ZZ* bl22); 00758 00759 ah = al0 + al1*U + al2*U2; 00760 bh = bl0 + bl1*U + bl2*U2; 00761 00762 // Compute the maximum of the rejection function 00763 grejmax = max(1. + xmin* (ah + bh*xmin), 1.+ah+bh); 00764 G4double xm = -ah/(2.*bh); 00765 if ( xmin < xm && xm < xmax) grejmax = max(grejmax, 1.+ xm* (ah + bh*xm)); 00766 } 00767 00768 // 00769 // sample the energy rate of the emitted gamma for e- kin energy > 1 MeV 00770 // 00771 00772 do { 00773 if (kineticEnergy > tlow) { 00774 do { 00775 q = G4UniformRand(); 00776 x = pow(xmin, q + kappa*(1.0 - q)); 00777 epsil = x*kineticEnergy/totalEnergy; 00778 G4double screenvar = screenfac*epsil/(1.0-epsil); 00779 G4double F1 = max(ScreenFunction1(screenvar) - FZ ,0.); 00780 G4double F2 = max(ScreenFunction2(screenvar) - FZ ,0.); 00781 migdal = (1. + MigdalFactor)/(1. + MigdalFactor/(x*x)); 00782 greject = migdal*(F1 - epsil* (ah*F1 - bh*epsil*F2))/(42.392 - FZ); 00783 /* 00784 if ( greject > grejmax ) { 00785 G4cout << "### G4eBremsstrahlungModel Warning: Majoranta exceeded! " 00786 << greject << " > " << grejmax 00787 << " x= " << x 00788 << " e= " << kineticEnergy 00789 << G4endl; 00790 } 00791 */ 00792 } while( greject < G4UniformRand()*grejmax ); 00793 00794 } else { 00795 00796 do { 00797 q = G4UniformRand(); 00798 x = pow(xmin, q + kappa*(1.0 - q)); 00799 migdal = (1. + MigdalFactor)/(1. + MigdalFactor/(x*x)); 00800 greject = migdal*(1. + x* (ah + bh*x)); 00801 /* 00802 if ( greject > grejmax ) { 00803 G4cout << "### G4eBremsstrahlungModel Warning: Majoranta exceeded! " 00804 << greject << " > " << grejmax 00805 << " x= " << x 00806 << " e= " << kineticEnergy 00807 << G4endl; 00808 } 00809 */ 00810 } while( greject < G4UniformRand()*grejmax ); 00811 } 00812 gammaEnergy = x*kineticEnergy; 00813 00814 if (LPMFlag()) { 00815 // take into account the supression due to the LPM effect 00816 if (G4UniformRand() <= SupressionFunction(material,kineticEnergy, 00817 gammaEnergy)) 00818 LPMOK = true; 00819 } 00820 else LPMOK = true; 00821 00822 } while (!LPMOK); 00823 00824 // 00825 // angles of the emitted gamma. ( Z - axis along the parent particle) 00826 // use general interface 00827 // 00828 00829 G4ThreeVector gammaDirection = 00830 GetAngularDistribution()->SampleDirection(dp, totalEnergy-gammaEnergy, 00831 G4lrint(anElement->GetZ()), 00832 couple->GetMaterial()); 00833 00834 // create G4DynamicParticle object for the Gamma 00835 G4DynamicParticle* gamma = new G4DynamicParticle(theGamma,gammaDirection, 00836 gammaEnergy); 00837 vdp->push_back(gamma); 00838 00839 G4double totMomentum = sqrt(kineticEnergy*(totalEnergy + electron_mass_c2)); 00840 00841 G4ThreeVector direction = (totMomentum*dp->GetMomentumDirection() 00842 - gammaEnergy*gammaDirection).unit(); 00843 00844 // energy of primary 00845 G4double finalE = kineticEnergy - gammaEnergy; 00846 00847 // stop tracking and create new secondary instead of primary 00848 if(gammaEnergy > SecondaryThreshold()) { 00849 fParticleChange->ProposeTrackStatus(fStopAndKill); 00850 fParticleChange->SetProposedKineticEnergy(0.0); 00851 G4DynamicParticle* el = 00852 new G4DynamicParticle(const_cast<G4ParticleDefinition*>(particle), 00853 direction, finalE); 00854 vdp->push_back(el); 00855 00856 // continue tracking 00857 } else { 00858 fParticleChange->SetProposedMomentumDirection(direction); 00859 fParticleChange->SetProposedKineticEnergy(finalE); 00860 } 00861 }
const G4Element * G4eBremsstrahlungModel::SelectRandomAtom | ( | const G4MaterialCutsCouple * | couple | ) | [protected] |
Definition at line 865 of file G4eBremsstrahlungModel.cc.
References G4UniformRand, G4Material::GetElementVector(), G4MaterialCutsCouple::GetIndex(), G4MaterialCutsCouple::GetMaterial(), G4Material::GetNumberOfElements(), and G4VEmModel::SetCurrentElement().
Referenced by SampleSecondaries().
00867 { 00868 // select randomly 1 element within the material 00869 00870 const G4Material* material = couple->GetMaterial(); 00871 G4int nElements = material->GetNumberOfElements(); 00872 const G4ElementVector* theElementVector = material->GetElementVector(); 00873 00874 const G4Element* elm = 0; 00875 00876 if(1 < nElements) { 00877 00878 --nElements; 00879 G4DataVector* dv = partialSumSigma[couple->GetIndex()]; 00880 G4double rval = G4UniformRand()*((*dv)[nElements]); 00881 00882 elm = (*theElementVector)[nElements]; 00883 for (G4int i=0; i<nElements; ++i) { 00884 if (rval <= (*dv)[i]) { 00885 elm = (*theElementVector)[i]; 00886 break; 00887 } 00888 } 00889 } else { elm = (*theElementVector)[0]; } 00890 00891 SetCurrentElement(elm); 00892 return elm; 00893 }
Definition at line 131 of file G4eBremsstrahlungModel.hh.
Referenced by G4eBremsstrahlungModel(), Initialise(), G4ePolarizedBremsstrahlungModel::SampleSecondaries(), and SampleSecondaries().
G4bool G4eBremsstrahlungModel::isElectron [protected] |
Definition at line 133 of file G4eBremsstrahlungModel.hh.
Referenced by ComputeCrossSectionPerAtom(), and ComputeDEDXPerVolume().
const G4ParticleDefinition* G4eBremsstrahlungModel::particle [protected] |
Definition at line 129 of file G4eBremsstrahlungModel.hh.
Referenced by ComputeDEDXPerVolume(), and CrossSectionPerVolume().
G4ParticleDefinition* G4eBremsstrahlungModel::theGamma [protected] |
Definition at line 130 of file G4eBremsstrahlungModel.hh.
Referenced by G4eBremsstrahlungModel(), and SampleSecondaries().