G4WentzelOKandVIxSection Class Reference

#include <G4WentzelOKandVIxSection.hh>


Public Member Functions

 G4WentzelOKandVIxSection ()
virtual ~G4WentzelOKandVIxSection ()
void Initialise (const G4ParticleDefinition *, G4double CosThetaLim)
void SetupParticle (const G4ParticleDefinition *)
G4double SetupTarget (G4int Z, G4double cut=DBL_MAX)
G4double ComputeTransportCrossSectionPerAtom (G4double CosThetaMax)
G4ThreeVector SampleSingleScattering (G4double CosThetaMin, G4double CosThetaMax, G4double elecRatio=0.0)
G4double ComputeNuclearCrossSection (G4double CosThetaMin, G4double CosThetaMax)
G4double ComputeElectronCrossSection (G4double CosThetaMin, G4double CosThetaMax)
G4double SetupKinematic (G4double kinEnergy, const G4Material *mat)
void SetTargetMass (G4double value)
G4double GetMomentumSquare () const
G4double GetCosThetaNuc () const
G4double GetCosThetaElec () const

Detailed Description

Definition at line 71 of file G4WentzelOKandVIxSection.hh.

Constructor & Destructor Documentation

G4WentzelOKandVIxSection::G4WentzelOKandVIxSection (  ) 

Definition at line 66 of file G4WentzelOKandVIxSection.cc.

References DBL_MAX, G4Electron::Electron(), G4NistManager::GetA27(), G4Pow::GetInstance(), G4NistManager::Instance(), G4INCL::Math::pi, G4Positron::Positron(), G4Proton::Proton(), and G4Pow::Z13().

00066                                                    :
00067   numlimit(0.1),
00068   nwarnings(0),
00069   nwarnlimit(50),
00070   alpha2(fine_structure_const*fine_structure_const)
00071 {
00072   fNistManager = G4NistManager::Instance();
00073   fG4pow = G4Pow::GetInstance();
00074   theElectron = G4Electron::Electron();
00075   thePositron = G4Positron::Positron();
00076   theProton   = G4Proton::Proton();
00077   lowEnergyLimit = 1.0*eV;
00078   G4double p0 = electron_mass_c2*classic_electr_radius;
00079   coeff = twopi*p0*p0;
00080   particle = 0;
00081 
00082   // Thomas-Fermi screening radii
00083   // Formfactors from A.V. Butkevich et al., NIM A 488 (2002) 282
00084 
00085   if(0.0 == ScreenRSquare[0]) {
00086     G4double a0 = electron_mass_c2/0.88534; 
00087     G4double constn = 6.937e-6/(MeV*MeV);
00088 
00089     ScreenRSquare[0] = alpha2*a0*a0;
00090     ScreenRSquareElec[0] = ScreenRSquare[0];
00091     for(G4int j=1; j<100; ++j) {
00092       G4double x = a0*fG4pow->Z13(j);
00093       if(1 == j) { ScreenRSquare[j] = 0.5*alpha2*a0*a0; }
00094       else {
00095         ScreenRSquare[j] = 0.5*(1 + exp(-j*j*0.001))*alpha2*x*x;
00096         ScreenRSquareElec[j] = 0.5*alpha2*x*x;
00097       }
00098       x = fNistManager->GetA27(j);
00099       FormFactor[j] = constn*x*x;
00100     } 
00101   }
00102   currentMaterial = 0;
00103   elecXSRatio = factB = factD = formfactA = screenZ = 0.0;
00104   cosTetMaxElec = cosTetMaxNuc = invbeta2 = kinFactor = gam0pcmp = pcmp2 = 1.0;
00105 
00106   factB1= 0.5*CLHEP::pi*fine_structure_const;
00107 
00108   tkin = mom2 = momCM2 = factorA2 = mass = spin = chargeSquare = charge3 = 0.0;
00109   ecut = etag = DBL_MAX;
00110   targetZ = 0;
00111   cosThetaMax = 1.0;
00112   targetMass = proton_mass_c2;
00113   particle = 0;
00114 }

G4WentzelOKandVIxSection::~G4WentzelOKandVIxSection (  )  [virtual]

Definition at line 118 of file G4WentzelOKandVIxSection.cc.

00119 {}

Member Function Documentation

G4double G4WentzelOKandVIxSection::ComputeElectronCrossSection ( G4double  CosThetaMin,
G4double  CosThetaMax 
) [inline]

Definition at line 237 of file G4WentzelOKandVIxSection.hh.

Referenced by G4eCoulombScatteringModel::ComputeCrossSectionPerAtom().

00239 {
00240   G4double xsec = 0.0;
00241   G4double cost1 = std::max(cosTMin,cosTetMaxElec);
00242   G4double cost2 = std::max(cosTMax,cosTetMaxElec);
00243   if(cost1 > cost2) {
00244     xsec = kinFactor*(cost1 - cost2)/((1.0 - cost1 + screenZ)*(1.0 - cost2 + screenZ));
00245   }
00246   return xsec;
00247 }

G4double G4WentzelOKandVIxSection::ComputeNuclearCrossSection ( G4double  CosThetaMin,
G4double  CosThetaMax 
) [inline]

Definition at line 223 of file G4WentzelOKandVIxSection.hh.

Referenced by G4eCoulombScatteringModel::ComputeCrossSectionPerAtom().

00225 {
00226   G4double xsec = 0.0;
00227   if(cosTMax < cosTMin) {
00228     xsec = targetZ*kinFactor*(cosTMin - cosTMax)/
00229       ((1.0 - cosTMin + screenZ)*(1.0 - cosTMax + screenZ));
00230   }
00231   return xsec;
00232 }

G4double G4WentzelOKandVIxSection::ComputeTransportCrossSectionPerAtom ( G4double  CosThetaMax  ) 

Definition at line 202 of file G4WentzelOKandVIxSection.cc.

References G4cout, G4endl, and G4ParticleDefinition::GetParticleName().

Referenced by G4WentzelVIModel::ComputeCrossSectionPerAtom().

00203 {
00204   G4double xsec = 0.0;
00205   if(cosTMax >= 1.0) { return xsec; }
00206  
00207   G4double xSection = 0.0;
00208   G4double x = 0; 
00209   G4double y = 0;
00210   G4double x1= 0;
00211   G4double x2= 0;
00212   G4double xlog = 0.0;
00213 
00214   G4double costm = std::max(cosTMax,cosTetMaxElec); 
00215   G4double fb = screenZ*factB;
00216 
00217   // scattering off electrons
00218   if(costm < 1.0) {
00219     x = (1.0 - costm)/screenZ;
00220     if(x < numlimit) { 
00221       x2 = 0.5*x*x;
00222       y  = x2*(1.0 - 1.3333333*x + 3*x2);
00223       if(0.0 < factB) { y -= fb*x2*x*(0.6666667 - x); }
00224     } else { 
00225       x1= x/(1 + x);
00226       xlog = log(1.0 + x);  
00227       y = xlog - x1; 
00228       if(0.0 < factB) { y -= fb*(x + x1 - 2*xlog); }
00229     }
00230 
00231     if(y < 0.0) {
00232       ++nwarnings;
00233       if(nwarnings < nwarnlimit) {
00234         G4cout << "G4WentzelOKandVIxSection::ComputeTransportCrossSectionPerAtom scattering on e- <0"
00235                << G4endl;
00236         G4cout << "y= " << y 
00237                << " e(MeV)= " << tkin << " p(MeV/c)= " << sqrt(mom2) 
00238                << " Z= " << targetZ << "  " 
00239                << particle->GetParticleName() << G4endl;
00240         G4cout << " 1-costm= " << 1.0-costm << " screenZ= " << screenZ 
00241                << " x= " << x << G4endl;
00242       }
00243       y = 0.0;
00244     }
00245     xSection = y;
00246   }
00247   /* 
00248   G4cout << "G4WentzelVI:XS per A " << " Z= " << targetZ 
00249          << " e(MeV)= " << tkin/MeV << " XSel= " << xSection
00250          << " cut(MeV)= " << ecut/MeV  
00251          << " zmaxE= " << (1.0 - cosTetMaxElec)/screenZ 
00252          << " zmaxN= " << (1.0 - cosThetaMax)/screenZ 
00253          << " 1-costm= " << 1.0 - cosThetaMax << G4endl;
00254   */
00255   // scattering off nucleus
00256   if(cosTMax < 1.0) {
00257     x = (1.0 - cosTMax)/screenZ;
00258     if(x < numlimit) { 
00259       x2 = 0.5*x*x;
00260       y  = x2*(1.0 - 1.3333333*x + 3*x2); 
00261       if(0.0 < factB) { y -= fb*x2*x*(0.6666667 - x); }
00262     } else { 
00263       x1= x/(1 + x);
00264       xlog = log(1.0 + x);  
00265       y = xlog - x1; 
00266       if(0.0 < factB) { y -= fb*(x + x1 - 2*xlog); }
00267     }
00268 
00269     if(y < 0.0) {
00270       ++nwarnings;
00271       if(nwarnings < nwarnlimit) {
00272         G4cout << "G4WentzelOKandVIxSection::ComputeTransportCrossSectionPerAtom scattering on e- <0"
00273                << G4endl;
00274         G4cout << "y= " << y 
00275                << " e(MeV)= " << tkin << " Z= " << targetZ << "  " 
00276                << particle->GetParticleName() << G4endl;
00277         G4cout << " formfactA= " << formfactA << " screenZ= " << screenZ 
00278                << " x= " << " x1= " << x1 <<G4endl;
00279       }
00280       y = 0.0;
00281     }
00282     xSection += y*targetZ; 
00283   }
00284   xSection *= kinFactor;
00285   /*
00286   G4cout << "Z= " << targetZ << " XStot= " << xSection/barn 
00287          << " screenZ= " << screenZ << " formF= " << formfactA 
00288          << " for " << particle->GetParticleName() 
00289          << " m= " << mass << " 1/v= " << sqrt(invbeta2) << " p= " << sqrt(mom2)
00290          << " x= " << x 
00291          << G4endl;
00292   */
00293   return xSection; 
00294 }

G4double G4WentzelOKandVIxSection::GetCosThetaElec (  )  const [inline]

Definition at line 215 of file G4WentzelOKandVIxSection.hh.

00216 {
00217   return cosTetMaxElec;
00218 }

G4double G4WentzelOKandVIxSection::GetCosThetaNuc (  )  const [inline]

Definition at line 208 of file G4WentzelOKandVIxSection.hh.

00209 {
00210   return cosTetMaxNuc;
00211 }

G4double G4WentzelOKandVIxSection::GetMomentumSquare (  )  const [inline]

Definition at line 201 of file G4WentzelOKandVIxSection.hh.

Referenced by G4eCoulombScatteringModel::SampleSecondaries().

00202 {
00203   return mom2;
00204 }

void G4WentzelOKandVIxSection::Initialise ( const G4ParticleDefinition *  ,
G4double  CosThetaLim 
)

Definition at line 123 of file G4WentzelOKandVIxSection.cc.

References DBL_MAX, G4LossTableManager::FactorForAngleLimit(), G4LossTableManager::Instance(), and SetupParticle().

Referenced by G4WentzelVIModel::Initialise(), and G4eCoulombScatteringModel::Initialise().

00125 {
00126   SetupParticle(p);
00127   tkin = mom2 = momCM2 = 0.0;
00128   ecut = etag = DBL_MAX;
00129   targetZ = 0;
00130   cosThetaMax = CosThetaLim;
00131   G4double a = G4LossTableManager::Instance()->FactorForAngleLimit()
00132     *CLHEP::hbarc/CLHEP::fermi;
00133   factorA2 = 0.5*a*a;
00134   currentMaterial = 0;
00135   
00136   //G4cout << "G4WentzelOKandVIxSection::Initialise  mass= " << mass
00137   //     << "  " << p->GetParticleName() << G4endl; 
00138   
00139 }

G4ThreeVector G4WentzelOKandVIxSection::SampleSingleScattering ( G4double  CosThetaMin,
G4double  CosThetaMax,
G4double  elecRatio = 0.0 
)

Definition at line 299 of file G4WentzelOKandVIxSection.cc.

References G4UniformRand.

Referenced by G4WentzelVIModel::SampleScattering(), and G4eCoulombScatteringModel::SampleSecondaries().

00302 {
00303   G4ThreeVector v(0.0,0.0,1.0);
00304  
00305   G4double formf = formfactA;
00306   G4double cost1 = cosTMin;
00307   G4double cost2 = cosTMax;
00308   if(elecRatio > 0.0) {
00309     if(G4UniformRand() <= elecRatio) {
00310       formf = 0.0;
00311       cost1 = std::max(cost1,cosTetMaxElec);
00312       cost2 = std::max(cost2,cosTetMaxElec);
00313     }
00314   }
00315   if(cost1 < cost2) { return v; }
00316 
00317   G4double w1 = 1. - cost1 + screenZ;
00318   G4double w2 = 1. - cost2 + screenZ;
00319   G4double z1 = w1*w2/(w1 + G4UniformRand()*(w2 - w1)) - screenZ;
00320 
00321   //G4double fm =  1.0 + formf*z1/(1.0 + (mass + tkin)*z1/targetMass);
00322   G4double fm =  1.0 + formf*z1;
00323   //G4double grej = (1. - z1*factB)/( (1.0 + z1*factD)*fm*fm );
00324   G4double grej = (1. - z1*factB + factB1*targetZ*sqrt(z1*factB)*(2 - z1))/( (1.0 + z1*factD)*fm*fm );
00325   // "false" scattering
00326   if( G4UniformRand() > grej ) { return v; }
00327     // } 
00328   G4double cost = 1.0 - z1;
00329 
00330   if(cost > 1.0)       { cost = 1.0; }
00331   else if(cost < -1.0) { cost =-1.0; }
00332   G4double sint = sqrt((1.0 - cost)*(1.0 + cost));
00333   //G4cout << "sint= " << sint << G4endl;
00334   G4double phi  = twopi*G4UniformRand();
00335   G4double vx1 = sint*cos(phi);
00336   G4double vy1 = sint*sin(phi);
00337 
00338   // only direction is changed
00339   v.set(vx1,vy1,cost);
00340   return v;
00341 }

void G4WentzelOKandVIxSection::SetTargetMass ( G4double  value  )  [inline]

Definition at line 193 of file G4WentzelOKandVIxSection.hh.

Referenced by G4eCoulombScatteringModel::SampleSecondaries(), and SetupTarget().

00194 {
00195   targetMass = value;
00196   factD = std::sqrt(mom2)/value;
00197 }

G4double G4WentzelOKandVIxSection::SetupKinematic ( G4double  kinEnergy,
const G4Material *  mat 
) [inline]

Definition at line 177 of file G4WentzelOKandVIxSection.hh.

References G4IonisParamMat::GetInvA23(), and G4Material::GetIonisation().

Referenced by G4WentzelVIModel::ComputeCrossSectionPerAtom(), G4eCoulombScatteringModel::ComputeCrossSectionPerAtom(), G4WentzelVIModel::ComputeGeomPathLength(), G4WentzelVIModel::ComputeTruePathLengthLimit(), and G4WentzelVIModel::ComputeTrueStepLength().

00178 {
00179   if(ekin != tkin || mat != currentMaterial) { 
00180     currentMaterial = mat;
00181     tkin  = ekin;
00182     mom2  = tkin*(tkin + 2.0*mass);
00183     invbeta2 = 1.0 +  mass*mass/mom2;
00184     factB = spin/invbeta2;
00185     cosTetMaxNuc = 
00186         std::max(cosThetaMax,1.-factorA2*mat->GetIonisation()->GetInvA23()/mom2);
00187   } 
00188   return cosTetMaxNuc;
00189 }

void G4WentzelOKandVIxSection::SetupParticle ( const G4ParticleDefinition *   ) 

Definition at line 143 of file G4WentzelOKandVIxSection.cc.

References G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGMass(), and G4ParticleDefinition::GetPDGSpin().

Referenced by Initialise(), and G4eCoulombScatteringModel::SetupParticle().

00144 {
00145   particle = p;
00146   mass = particle->GetPDGMass();
00147   spin = particle->GetPDGSpin();
00148   if(0.0 != spin) { spin = 0.5; }
00149   G4double q = std::fabs(particle->GetPDGCharge()/eplus);
00150   chargeSquare = q*q;
00151   charge3 = chargeSquare*q;
00152   tkin = 0.0;
00153   currentMaterial = 0;
00154   targetZ = 0;
00155 }

G4double G4WentzelOKandVIxSection::SetupTarget ( G4int  Z,
G4double  cut = DBL_MAX 
)

Definition at line 160 of file G4WentzelOKandVIxSection.cc.

References G4NistManager::GetAtomicMassAmu(), SetTargetMass(), and G4Pow::Z23().

Referenced by G4WentzelVIModel::ComputeCrossSectionPerAtom(), G4eCoulombScatteringModel::ComputeCrossSectionPerAtom(), and G4WentzelVIModel::SampleScattering().

00161 {
00162   G4double cosTetMaxNuc2 = cosTetMaxNuc;
00163   if(Z != targetZ || tkin != etag) {
00164     etag    = tkin; 
00165     targetZ = Z;
00166     if(targetZ > 99) { targetZ = 99; }
00167     SetTargetMass(fNistManager->GetAtomicMassAmu(targetZ)*CLHEP::amu_c2);
00168     //G4double tmass2 = targetMass*targetMass;
00169     //G4double etot = tkin + mass;
00170     //G4double invmass2 = mass*mass + tmass2 + 2*etot*targetMass;
00171     //momCM2 = mom2*tmass2/invmass2;
00172     //gam0pcmp = (etot + targetMass)*targetMass/invmass2;
00173     //pcmp2    = tmass2/invmass2;
00174 
00175     kinFactor = coeff*Z*chargeSquare*invbeta2/mom2;
00176 
00177     if(1 == Z) {
00178       screenZ = ScreenRSquare[targetZ]/mom2;
00179     } else if(mass > MeV) {
00180       screenZ = std::min(Z*1.13,1.13 +3.76*Z*Z*invbeta2*alpha2*chargeSquare)*
00181         ScreenRSquare[targetZ]/mom2;
00182     } else {
00183       G4double tau = tkin/mass;
00184       screenZ = std::min(Z*1.13,(1.13 +3.76*Z*Z
00185           *invbeta2*alpha2*std::sqrt(tau/(tau + fG4pow->Z23(targetZ)))))*
00186         ScreenRSquareElec[targetZ]/mom2;
00187     }
00188     if(targetZ == 1 && cosTetMaxNuc2 < 0.0 && particle == theProton) {
00189       cosTetMaxNuc2 = 0.0;
00190     }
00191     formfactA = FormFactor[targetZ]*mom2;
00192 
00193     cosTetMaxElec = 1.0;
00194     ComputeMaxElectronScattering(cut); 
00195   }
00196   return cosTetMaxNuc2;
00197 }

The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:59 2013 for Geant4 by  doxygen 1.4.7