#include <G4PairProductionRelModel.hh>
Inheritance diagram for G4PairProductionRelModel:
Definition at line 61 of file G4PairProductionRelModel.hh.
G4PairProductionRelModel::G4PairProductionRelModel | ( | const G4ParticleDefinition * | p = 0 , |
|
const G4String & | nam = "BetheHeitlerLPM" | |||
) |
Definition at line 86 of file G4PairProductionRelModel.cc.
References currentZ, G4Electron::Electron(), fCoulomb, Fel, Finel, fParticleChange, G4Gamma::Gamma(), gLPM, G4NistManager::Instance(), lnZ, nist, phiLPM, G4Positron::Positron(), theElectron, theGamma, thePositron, xiLPM, z13, and z23.
00088 : G4VEmModel(nam), 00089 fLPMconstant(fine_structure_const*electron_mass_c2*electron_mass_c2/(4.*pi*hbarc)*0.5), 00090 fLPMflag(true), 00091 lpmEnergy(0.), 00092 use_completescreening(false) 00093 { 00094 fParticleChange = 0; 00095 theGamma = G4Gamma::Gamma(); 00096 thePositron = G4Positron::Positron(); 00097 theElectron = G4Electron::Electron(); 00098 00099 nist = G4NistManager::Instance(); 00100 00101 currentZ = z13 = z23 = lnZ = Fel = Finel = fCoulomb = phiLPM = gLPM = xiLPM = 0; 00102 00103 }
G4PairProductionRelModel::~G4PairProductionRelModel | ( | ) | [virtual] |
G4PairProductionRelModel::G4PairProductionRelModel | ( | const G4PairProductionRelModel & | ) | [protected] |
Definition at line 222 of file G4PairProductionRelModel.cc.
References facFel, gLPM, lnZ, logTwo, lpmEnergy, phiLPM, G4INCL::Math::pi, preS1, G4InuclParticleNames::s0, sqr(), xiLPM, and z23.
Referenced by ComputeRelDXSectionPerAtom(), and SampleSecondaries().
00223 { 00224 // *** calculate lpm variable s & sprime *** 00225 // Klein eqs. (78) & (79) 00226 G4double sprime = sqrt(0.125*k*lpmEnergy/(eplusEnergy*(k-eplusEnergy))); 00227 00228 G4double s1 = preS1*z23; 00229 G4double logS1 = 2./3.*lnZ-2.*facFel; 00230 G4double logTS1 = logTwo+logS1; 00231 00232 xiLPM = 2.; 00233 00234 if (sprime>1) 00235 xiLPM = 1.; 00236 else if (sprime>sqrt(2.)*s1) { 00237 G4double h = log(sprime)/logTS1; 00238 xiLPM = 1+h-0.08*(1-h)*(1-sqr(1-h))/logTS1; 00239 } 00240 00241 G4double s0 = sprime/sqrt(xiLPM); 00242 // G4cout<<"k="<<k<<" y="<<eplusEnergy/k<<G4endl; 00243 // G4cout<<"s0="<<s0<<G4endl; 00244 00245 // *** calculate supression functions phi and G *** 00246 // Klein eqs. (77) 00247 G4double s2=s0*s0; 00248 G4double s3=s0*s2; 00249 G4double s4=s2*s2; 00250 00251 if (s0<0.1) { 00252 // high suppression limit 00253 phiLPM = 6.*s0 - 18.84955592153876*s2 + 39.47841760435743*s3 00254 - 57.69873135166053*s4; 00255 gLPM = 37.69911184307752*s2 - 236.8705056261446*s3 + 807.7822389*s4; 00256 } 00257 else if (s0<1.9516) { 00258 // intermediate suppression 00259 // using eq.77 approxim. valid s0<2. 00260 phiLPM = 1.-exp(-6.*s0*(1.+(3.-pi)*s0) 00261 +s3/(0.623+0.795*s0+0.658*s2)); 00262 if (s0<0.415827397755) { 00263 // using eq.77 approxim. valid 0.07<s<2 00264 G4double psiLPM = 1-exp(-4*s0-8*s2/(1+3.936*s0+4.97*s2-0.05*s3+7.50*s4)); 00265 gLPM = 3*psiLPM-2*phiLPM; 00266 } 00267 else { 00268 // using alternative parametrisiation 00269 G4double pre = -0.16072300849123999 + s0*3.7550300067531581 + s2*-1.7981383069010097 00270 + s3*0.67282686077812381 + s4*-0.1207722909879257; 00271 gLPM = tanh(pre); 00272 } 00273 } 00274 else { 00275 // low suppression limit valid s>2. 00276 phiLPM = 1. - 0.0119048/s4; 00277 gLPM = 1. - 0.0230655/s4; 00278 } 00279 00280 // *** make sure suppression is smaller than 1 *** 00281 // *** caused by Migdal approximation in xi *** 00282 if (xiLPM*phiLPM>1. || s0>0.57) xiLPM=1./phiLPM; 00283 }
G4double G4PairProductionRelModel::ComputeCrossSectionPerAtom | ( | const G4ParticleDefinition * | , | |
G4double | kinEnergy, | |||
G4double | Z, | |||
G4double | A = 0. , |
|||
G4double | cut = 0. , |
|||
G4double | emax = DBL_MAX | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 289 of file G4PairProductionRelModel.cc.
References ComputeXSectionPerAtom(), fCoulomb, Fel, Finel, and SetCurrentElement().
00292 { 00293 // static const G4double gammaEnergyLimit = 1.5*MeV; 00294 G4double crossSection = 0.0 ; 00295 if ( Z < 0.9 ) return crossSection; 00296 if ( gammaEnergy <= 2.0*electron_mass_c2 ) return crossSection; 00297 00298 SetCurrentElement(Z); 00299 // choose calculator according to parameters and switches 00300 // in the moment only one calculator: 00301 crossSection=ComputeXSectionPerAtom(gammaEnergy,Z); 00302 00303 G4double xi = Finel/(Fel - fCoulomb); // inelastic contribution 00304 crossSection*=4.*Z*(Z+xi)*fine_structure_const*classic_electr_radius*classic_electr_radius; 00305 00306 return crossSection; 00307 }
G4double G4PairProductionRelModel::ComputeDXSectionPerAtom | ( | G4double | eplusEnergy, | |
G4double | totalEnergy, | |||
G4double | Z | |||
) | [protected] |
Definition at line 165 of file G4PairProductionRelModel.cc.
References DeltaMin(), fCoulomb, Fel, lnZ, Phi1(), Phi2(), and use_completescreening.
Referenced by ComputeXSectionPerAtom().
00168 { 00169 // most simple case - complete screening: 00170 00171 // dsig/dE+ = 4 * alpha * Z**2 * r0**2 / k 00172 // * [ (y**2 + (1-y**2) + 2/3*y*(1-y) ) * ( log (183 * Z**-1/3) + 1/9 * y*(1-y) ] 00173 // y = E+/k 00174 G4double yp=eplusEnergy/totalEnergy; 00175 G4double ym=1.-yp; 00176 00177 G4double cross = 0.; 00178 if (use_completescreening) 00179 cross = (yp*yp + ym*ym + 2./3.*ym*yp)*(Fel - fCoulomb) + 1./9.*yp*ym; 00180 else { 00181 G4double delta = 0.25*DeltaMin(totalEnergy)/(yp*ym); 00182 cross = (yp*yp + ym*ym)*(0.25*Phi1(delta) - lnZ/3. - fCoulomb) 00183 + 2./3.*ym*yp*(0.25*Phi2(delta) - lnZ/3. - fCoulomb); 00184 } 00185 return cross/totalEnergy; 00186 00187 }
G4double G4PairProductionRelModel::ComputeRelDXSectionPerAtom | ( | G4double | eplusEnergy, | |
G4double | totalEnergy, | |||
G4double | Z | |||
) | [protected] |
Definition at line 191 of file G4PairProductionRelModel.cc.
References CalcLPMFunctions(), DeltaMin(), fCoulomb, Fel, gLPM, lnZ, Phi1(), Phi2(), phiLPM, use_completescreening, and xiLPM.
Referenced by ComputeXSectionPerAtom().
00194 { 00195 // most simple case - complete screening: 00196 00197 // dsig/dE+ = 4 * alpha * Z**2 * r0**2 / k 00198 // * [ (y**2 + (1-y**2) + 2/3*y*(1-y) ) * ( log (183 * Z**-1/3) + 1/9 * y*(1-y) ] 00199 // y = E+/k 00200 G4double yp=eplusEnergy/totalEnergy; 00201 G4double ym=1.-yp; 00202 00203 CalcLPMFunctions(totalEnergy,eplusEnergy); // gamma 00204 00205 G4double cross = 0.; 00206 if (use_completescreening) 00207 cross = xiLPM*(2./3.*phiLPM*(yp*yp + ym*ym) + gLPM)*(Fel - fCoulomb); 00208 else { 00209 G4double delta = 0.25*DeltaMin(totalEnergy)/(yp*ym); 00210 cross = (1./3.*gLPM + 2./3.*phiLPM)*(yp*yp + ym*ym) 00211 *(0.25*Phi1(delta) - lnZ/3. - fCoulomb) 00212 + 2./3.*gLPM*ym*yp*(0.25*Phi2(delta) - lnZ/3. - fCoulomb); 00213 cross *= xiLPM; 00214 } 00215 return cross/totalEnergy; 00216 00217 }
G4double G4PairProductionRelModel::ComputeXSectionPerAtom | ( | G4double | totalEnergy, | |
G4double | Z | |||
) | [protected] |
Definition at line 121 of file G4PairProductionRelModel.cc.
References ComputeDXSectionPerAtom(), ComputeRelDXSectionPerAtom(), DeltaMax(), DeltaMin(), fLPMflag, CLHEP::detail::n, wgi, and xgi.
Referenced by ComputeCrossSectionPerAtom().
00122 { 00123 G4double cross = 0.0; 00124 00125 // number of intervals and integration step 00126 G4double vcut = electron_mass_c2/totalEnergy ; 00127 00128 // limits by the screening variable 00129 G4double dmax = DeltaMax(); 00130 G4double dmin = min(DeltaMin(totalEnergy),dmax); 00131 G4double vcut1 = 0.5 - 0.5*sqrt(1. - dmin/dmax) ; 00132 vcut = max(vcut, vcut1); 00133 00134 00135 G4double vmax = 0.5; 00136 G4int n = 1; // needs optimisation 00137 00138 G4double delta = (vmax - vcut)*totalEnergy/G4double(n); 00139 00140 G4double e0 = vcut*totalEnergy; 00141 G4double xs; 00142 00143 // simple integration 00144 for(G4int l=0; l<n; l++,e0 += delta) { 00145 for(G4int i=0; i<8; i++) { 00146 00147 G4double eg = (e0 + xgi[i]*delta); 00148 if (fLPMflag && totalEnergy>100.*GeV) 00149 xs = ComputeRelDXSectionPerAtom(eg,totalEnergy,Z); 00150 else 00151 xs = ComputeDXSectionPerAtom(eg,totalEnergy,Z); 00152 cross += wgi[i]*xs; 00153 00154 } 00155 } 00156 00157 cross *= delta*2.; 00158 00159 return cross; 00160 }
G4double G4PairProductionRelModel::DeltaMax | ( | ) | const [inline, protected] |
Definition at line 274 of file G4PairProductionRelModel.hh.
Referenced by ComputeXSectionPerAtom().
00275 { 00276 // k > 50 MeV 00277 G4double FZ = 8.*(lnZ/3. + fCoulomb); 00278 return std::exp( (42.24-FZ)/8.368 ) + 0.952; 00279 }
Definition at line 281 of file G4PairProductionRelModel.hh.
References z13.
Referenced by ComputeDXSectionPerAtom(), ComputeRelDXSectionPerAtom(), and ComputeXSectionPerAtom().
00282 { 00283 return 4.*136./z13*(CLHEP::electron_mass_c2/k); 00284 }
void G4PairProductionRelModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 112 of file G4PairProductionRelModel.cc.
References fParticleChange, G4VEmModel::GetParticleChangeForGamma(), and G4VEmModel::InitialiseElementSelectors().
00114 { 00115 if(!fParticleChange) { fParticleChange = GetParticleChangeForGamma(); } 00116 InitialiseElementSelectors(p, cuts); 00117 }
G4double G4PairProductionRelModel::LPMconstant | ( | ) | const [inline] |
Definition at line 166 of file G4PairProductionRelModel.hh.
References fLPMconstant.
00167 { 00168 return fLPMconstant; 00169 }
G4bool G4PairProductionRelModel::LPMflag | ( | ) | const [inline] |
Definition at line 182 of file G4PairProductionRelModel.hh.
References fLPMflag.
00183 { 00184 return fLPMflag; 00185 }
G4PairProductionRelModel& G4PairProductionRelModel::operator= | ( | const G4PairProductionRelModel & | right | ) | [protected] |
Definition at line 213 of file G4PairProductionRelModel.hh.
Referenced by ComputeDXSectionPerAtom(), ComputeRelDXSectionPerAtom(), and SampleSecondaries().
00214 { 00215 G4double screenVal; 00216 00217 if (delta > 1.) 00218 screenVal = 21.12 - 4.184*std::log(delta+0.952); 00219 else 00220 screenVal = 20.868 - delta*(3.242 - 0.625*delta); 00221 00222 return screenVal; 00223 }
Definition at line 227 of file G4PairProductionRelModel.hh.
Referenced by ComputeDXSectionPerAtom(), ComputeRelDXSectionPerAtom(), and SampleSecondaries().
00228 { 00229 G4double screenVal; 00230 00231 if (delta > 1.) 00232 screenVal = 21.12 - 4.184*std::log(delta+0.952); 00233 else 00234 screenVal = 20.209 - delta*(1.930 + 0.086*delta); 00235 00236 return screenVal; 00237 }
void G4PairProductionRelModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 312 of file G4PairProductionRelModel.cc.
References CalcLPMFunctions(), F10, F20, fLPMflag, fParticleChange, fStopAndKill, G4UniformRand, G4Element::GetfCoulomb(), G4Element::GetIonisation(), G4DynamicParticle::GetKineticEnergy(), G4IonisParamElm::GetlogZ3(), G4MaterialCutsCouple::GetMaterial(), G4DynamicParticle::GetMomentumDirection(), G4IonisParamElm::GetZ3(), gLPM, Phi1(), Phi2(), phiLPM, G4VParticleChange::ProposeTrackStatus(), ScreenFunction1(), ScreenFunction2(), G4VEmModel::SelectRandomAtom(), G4ParticleChangeForGamma::SetProposedKineticEnergy(), theElectron, theGamma, thePositron, and xiLPM.
00324 : Effects due to the breakdown of the Born approximation at 00325 // low energy are ignored. 00326 // Note 2 : The differential cross section implicitly takes account of 00327 // pair creation in both nuclear and atomic electron fields. 00328 // However triplet prodution is not generated. 00329 { 00330 const G4Material* aMaterial = couple->GetMaterial(); 00331 00332 G4double GammaEnergy = aDynamicGamma->GetKineticEnergy(); 00333 G4ParticleMomentum GammaDirection = aDynamicGamma->GetMomentumDirection(); 00334 00335 G4double epsil ; 00336 G4double epsil0 = electron_mass_c2/GammaEnergy ; 00337 if(epsil0 > 1.0) { return; } 00338 00339 // do it fast if GammaEnergy < 2. MeV 00340 static const G4double Egsmall=2.*MeV; 00341 00342 // select randomly one element constituing the material 00343 const G4Element* anElement = SelectRandomAtom(aMaterial, theGamma, GammaEnergy); 00344 00345 if (GammaEnergy < Egsmall) { 00346 00347 epsil = epsil0 + (0.5-epsil0)*G4UniformRand(); 00348 00349 } else { 00350 // now comes the case with GammaEnergy >= 2. MeV 00351 00352 // Extract Coulomb factor for this Element 00353 G4double FZ = 8.*(anElement->GetIonisation()->GetlogZ3()); 00354 if (GammaEnergy > 50.*MeV) FZ += 8.*(anElement->GetfCoulomb()); 00355 00356 // limits of the screening variable 00357 G4double screenfac = 136.*epsil0/(anElement->GetIonisation()->GetZ3()); 00358 G4double screenmax = exp ((42.24 - FZ)/8.368) - 0.952 ; 00359 G4double screenmin = min(4.*screenfac,screenmax); 00360 00361 // limits of the energy sampling 00362 G4double epsil1 = 0.5 - 0.5*sqrt(1. - screenmin/screenmax) ; 00363 G4double epsilmin = max(epsil0,epsil1) , epsilrange = 0.5 - epsilmin; 00364 00365 // 00366 // sample the energy rate of the created electron (or positron) 00367 // 00368 //G4double epsil, screenvar, greject ; 00369 G4double screenvar, greject ; 00370 00371 G4double F10 = ScreenFunction1(screenmin) - FZ; 00372 G4double F20 = ScreenFunction2(screenmin) - FZ; 00373 G4double NormF1 = max(F10*epsilrange*epsilrange,0.); 00374 G4double NormF2 = max(1.5*F20,0.); 00375 00376 do { 00377 if ( NormF1/(NormF1+NormF2) > G4UniformRand() ) { 00378 epsil = 0.5 - epsilrange*pow(G4UniformRand(), 0.333333); 00379 screenvar = screenfac/(epsil*(1-epsil)); 00380 if (fLPMflag && GammaEnergy>100.*GeV) { 00381 CalcLPMFunctions(GammaEnergy,GammaEnergy*epsil); 00382 greject = xiLPM*((gLPM+2.*phiLPM)*Phi1(screenvar) - gLPM*Phi2(screenvar) - phiLPM*FZ)/F10; 00383 } 00384 else { 00385 greject = (ScreenFunction1(screenvar) - FZ)/F10; 00386 } 00387 00388 } else { 00389 epsil = epsilmin + epsilrange*G4UniformRand(); 00390 screenvar = screenfac/(epsil*(1-epsil)); 00391 if (fLPMflag && GammaEnergy>100.*GeV) { 00392 CalcLPMFunctions(GammaEnergy,GammaEnergy*epsil); 00393 greject = xiLPM*((0.5*gLPM+phiLPM)*Phi1(screenvar) + 0.5*gLPM*Phi2(screenvar) - 0.5*(gLPM+phiLPM)*FZ)/F20; 00394 } 00395 else { 00396 greject = (ScreenFunction2(screenvar) - FZ)/F20; 00397 } 00398 } 00399 00400 } while( greject < G4UniformRand() ); 00401 00402 } // end of epsil sampling 00403 00404 // 00405 // fixe charges randomly 00406 // 00407 00408 G4double ElectTotEnergy, PositTotEnergy; 00409 if (G4UniformRand() > 0.5) { 00410 00411 ElectTotEnergy = (1.-epsil)*GammaEnergy; 00412 PositTotEnergy = epsil*GammaEnergy; 00413 00414 } else { 00415 00416 PositTotEnergy = (1.-epsil)*GammaEnergy; 00417 ElectTotEnergy = epsil*GammaEnergy; 00418 } 00419 00420 // 00421 // scattered electron (positron) angles. ( Z - axis along the parent photon) 00422 // 00423 // universal distribution suggested by L. Urban 00424 // (Geant3 manual (1993) Phys211), 00425 // derived from Tsai distribution (Rev Mod Phys 49,421(1977)) 00426 00427 G4double u; 00428 const G4double a1 = 0.625 , a2 = 3.*a1 , d = 27. ; 00429 00430 if (9./(9.+d) >G4UniformRand()) u= - log(G4UniformRand()*G4UniformRand())/a1; 00431 else u= - log(G4UniformRand()*G4UniformRand())/a2; 00432 00433 G4double TetEl = u*electron_mass_c2/ElectTotEnergy; 00434 G4double TetPo = u*electron_mass_c2/PositTotEnergy; 00435 G4double Phi = twopi * G4UniformRand(); 00436 G4double dxEl= sin(TetEl)*cos(Phi),dyEl= sin(TetEl)*sin(Phi),dzEl=cos(TetEl); 00437 G4double dxPo=-sin(TetPo)*cos(Phi),dyPo=-sin(TetPo)*sin(Phi),dzPo=cos(TetPo); 00438 00439 // 00440 // kinematic of the created pair 00441 // 00442 // the electron and positron are assumed to have a symetric 00443 // angular distribution with respect to the Z axis along the parent photon. 00444 00445 G4double ElectKineEnergy = max(0.,ElectTotEnergy - electron_mass_c2); 00446 00447 G4ThreeVector ElectDirection (dxEl, dyEl, dzEl); 00448 ElectDirection.rotateUz(GammaDirection); 00449 00450 // create G4DynamicParticle object for the particle1 00451 G4DynamicParticle* aParticle1= new G4DynamicParticle( 00452 theElectron,ElectDirection,ElectKineEnergy); 00453 00454 // the e+ is always created (even with Ekine=0) for further annihilation. 00455 00456 G4double PositKineEnergy = max(0.,PositTotEnergy - electron_mass_c2); 00457 00458 G4ThreeVector PositDirection (dxPo, dyPo, dzPo); 00459 PositDirection.rotateUz(GammaDirection); 00460 00461 // create G4DynamicParticle object for the particle2 00462 G4DynamicParticle* aParticle2= new G4DynamicParticle( 00463 thePositron,PositDirection,PositKineEnergy); 00464 00465 // Fill output vector 00466 fvect->push_back(aParticle1); 00467 fvect->push_back(aParticle2); 00468 00469 // kill incident photon 00470 fParticleChange->SetProposedKineticEnergy(0.); 00471 fParticleChange->ProposeTrackStatus(fStopAndKill); 00472 }
Definition at line 239 of file G4PairProductionRelModel.hh.
Referenced by SampleSecondaries().
00243 { 00244 G4double screenVal; 00245 00246 if (ScreenVariable > 1.) 00247 screenVal = 42.24 - 8.368*std::log(ScreenVariable+0.952); 00248 else 00249 screenVal = 42.392 - ScreenVariable*(7.796 - 1.961*ScreenVariable); 00250 00251 return screenVal; 00252 }
Definition at line 256 of file G4PairProductionRelModel.hh.
Referenced by SampleSecondaries().
00260 { 00261 G4double screenVal; 00262 00263 if (ScreenVariable > 1.) 00264 screenVal = 42.24 - 8.368*std::log(ScreenVariable+0.952); 00265 else 00266 screenVal = 41.405 - ScreenVariable*(5.828 - 0.8945*ScreenVariable); 00267 00268 return screenVal; 00269 }
void G4PairProductionRelModel::SetCurrentElement | ( | G4double | ) | [inline] |
Definition at line 189 of file G4PairProductionRelModel.hh.
References currentZ, facFel, facFinel, fCoulomb, Fel, Fel_light, Finel, Finel_light, G4VEmModel::GetCurrentElement(), G4Element::GetfCoulomb(), G4NistManager::GetLOGZ(), G4NistManager::GetZ13(), lnZ, nist, z13, and z23.
Referenced by ComputeCrossSectionPerAtom().
00190 { 00191 if(Z != currentZ) { 00192 currentZ = Z; 00193 00194 G4int iz = G4int(Z); 00195 z13 = nist->GetZ13(iz); 00196 z23 = z13*z13; 00197 lnZ = nist->GetLOGZ(iz); 00198 00199 if (iz <= 4) { 00200 Fel = Fel_light[iz]; 00201 Finel = Finel_light[iz] ; 00202 } 00203 else { 00204 Fel = facFel - lnZ/3. ; 00205 Finel = facFinel - 2.*lnZ/3. ; 00206 } 00207 fCoulomb=GetCurrentElement()->GetfCoulomb(); 00208 } 00209 }
void G4PairProductionRelModel::SetLPMconstant | ( | G4double | val | ) | [inline] |
Definition at line 158 of file G4PairProductionRelModel.hh.
References fLPMconstant.
00159 { 00160 fLPMconstant = val; 00161 }
void G4PairProductionRelModel::SetLPMflag | ( | G4bool | ) | [inline] |
Definition at line 174 of file G4PairProductionRelModel.hh.
References fLPMflag.
00175 { 00176 fLPMflag = val; 00177 }
void G4PairProductionRelModel::SetupForMaterial | ( | const G4ParticleDefinition * | , | |
const G4Material * | , | |||
G4double | ||||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 477 of file G4PairProductionRelModel.cc.
References fLPMconstant, G4Material::GetRadlen(), and lpmEnergy.
00479 { 00480 lpmEnergy = mat->GetRadlen()*fLPMconstant; 00481 // G4cout<<" lpmEnergy="<<lpmEnergy<<G4endl; 00482 }
G4double G4PairProductionRelModel::currentZ [protected] |
Definition at line 135 of file G4PairProductionRelModel.hh.
Referenced by G4PairProductionRelModel(), and SetCurrentElement().
const G4double G4PairProductionRelModel::facFel = log(184.15) [static, protected] |
Definition at line 147 of file G4PairProductionRelModel.hh.
Referenced by CalcLPMFunctions(), and SetCurrentElement().
const G4double G4PairProductionRelModel::facFinel = log(1194.) [static, protected] |
G4double G4PairProductionRelModel::fCoulomb [protected] |
Definition at line 134 of file G4PairProductionRelModel.hh.
Referenced by ComputeCrossSectionPerAtom(), ComputeDXSectionPerAtom(), ComputeRelDXSectionPerAtom(), DeltaMax(), G4PairProductionRelModel(), and SetCurrentElement().
G4double G4PairProductionRelModel::Fel [protected] |
Definition at line 134 of file G4PairProductionRelModel.hh.
Referenced by ComputeCrossSectionPerAtom(), ComputeDXSectionPerAtom(), ComputeRelDXSectionPerAtom(), G4PairProductionRelModel(), and SetCurrentElement().
const G4double G4PairProductionRelModel::Fel_light = {0., 5.31 , 4.79 , 4.74 , 4.71} [static, protected] |
G4double G4PairProductionRelModel::Finel [protected] |
Definition at line 134 of file G4PairProductionRelModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4PairProductionRelModel(), and SetCurrentElement().
const G4double G4PairProductionRelModel::Finel_light = {0., 6.144 , 5.621 , 5.805 , 5.924} [static, protected] |
G4double G4PairProductionRelModel::fLPMconstant [protected] |
Definition at line 129 of file G4PairProductionRelModel.hh.
Referenced by LPMconstant(), SetLPMconstant(), and SetupForMaterial().
G4bool G4PairProductionRelModel::fLPMflag [protected] |
Definition at line 130 of file G4PairProductionRelModel.hh.
Referenced by ComputeXSectionPerAtom(), LPMflag(), SampleSecondaries(), and SetLPMflag().
Definition at line 127 of file G4PairProductionRelModel.hh.
Referenced by G4PairProductionRelModel(), Initialise(), and SampleSecondaries().
G4double G4PairProductionRelModel::gLPM [protected] |
Definition at line 139 of file G4PairProductionRelModel.hh.
Referenced by CalcLPMFunctions(), ComputeRelDXSectionPerAtom(), G4PairProductionRelModel(), and SampleSecondaries().
G4double G4PairProductionRelModel::lnZ [protected] |
Definition at line 133 of file G4PairProductionRelModel.hh.
Referenced by CalcLPMFunctions(), ComputeDXSectionPerAtom(), ComputeRelDXSectionPerAtom(), DeltaMax(), G4PairProductionRelModel(), and SetCurrentElement().
const G4double G4PairProductionRelModel::logTwo = log(2.) [static, protected] |
G4double G4PairProductionRelModel::lpmEnergy [protected] |
Definition at line 138 of file G4PairProductionRelModel.hh.
Referenced by CalcLPMFunctions(), and SetupForMaterial().
G4NistManager* G4PairProductionRelModel::nist [protected] |
Definition at line 122 of file G4PairProductionRelModel.hh.
Referenced by G4PairProductionRelModel(), and SetCurrentElement().
G4double G4PairProductionRelModel::phiLPM [protected] |
Definition at line 139 of file G4PairProductionRelModel.hh.
Referenced by CalcLPMFunctions(), ComputeRelDXSectionPerAtom(), G4PairProductionRelModel(), and SampleSecondaries().
const G4double G4PairProductionRelModel::preS1 = 1./(184.15*184.15) [static, protected] |
Definition at line 125 of file G4PairProductionRelModel.hh.
Referenced by G4PairProductionRelModel(), and SampleSecondaries().
G4ParticleDefinition* G4PairProductionRelModel::theGamma [protected] |
Definition at line 124 of file G4PairProductionRelModel.hh.
Referenced by G4PairProductionRelModel(), and SampleSecondaries().
Definition at line 126 of file G4PairProductionRelModel.hh.
Referenced by G4PairProductionRelModel(), and SampleSecondaries().
Definition at line 142 of file G4PairProductionRelModel.hh.
Referenced by ComputeDXSectionPerAtom(), and ComputeRelDXSectionPerAtom().
const G4double G4PairProductionRelModel::wgi [static, protected] |
Initial value:
{ 0.0506, 0.1112, 0.1569, 0.1813, 0.1813, 0.1569, 0.1112, 0.0506 }
Definition at line 144 of file G4PairProductionRelModel.hh.
Referenced by ComputeXSectionPerAtom().
const G4double G4PairProductionRelModel::xgi [static, protected] |
Initial value:
{ 0.0199, 0.1017, 0.2372, 0.4083, 0.5917, 0.7628, 0.8983, 0.9801 }
Definition at line 144 of file G4PairProductionRelModel.hh.
Referenced by ComputeXSectionPerAtom().
G4double G4PairProductionRelModel::xiLPM [protected] |
Definition at line 139 of file G4PairProductionRelModel.hh.
Referenced by CalcLPMFunctions(), ComputeRelDXSectionPerAtom(), G4PairProductionRelModel(), and SampleSecondaries().
G4double G4PairProductionRelModel::z13 [protected] |
Definition at line 133 of file G4PairProductionRelModel.hh.
Referenced by DeltaMin(), G4PairProductionRelModel(), and SetCurrentElement().
G4double G4PairProductionRelModel::z23 [protected] |
Definition at line 133 of file G4PairProductionRelModel.hh.
Referenced by CalcLPMFunctions(), G4PairProductionRelModel(), and SetCurrentElement().