G4TripathiLightCrossSection Class Reference

#include <G4TripathiLightCrossSection.hh>

Inheritance diagram for G4TripathiLightCrossSection:

G4VCrossSectionDataSet

Public Member Functions

 G4TripathiLightCrossSection ()
 ~G4TripathiLightCrossSection ()
virtual G4bool IsElementApplicable (const G4DynamicParticle *theProjectile, G4int Z, const G4Material *)
virtual G4double GetElementCrossSection (const G4DynamicParticle *theProjectile, G4int Z, const G4Material *mat=0)
void SetLowEnergyCheck (G4bool)

Detailed Description

Definition at line 83 of file G4TripathiLightCrossSection.hh.


Constructor & Destructor Documentation

G4TripathiLightCrossSection::G4TripathiLightCrossSection (  ) 

Definition at line 77 of file G4TripathiLightCrossSection.cc.

00078  : G4VCrossSectionDataSet("TripathiLightIons")
00079 {
00080   // Constructor only needs to instantiate the object which provides functions
00081   // to calculate the nuclear radius, and some other constants used to
00082   // calculate cross-sections.
00083 
00084   theWilsonRadius = new G4WilsonRadius();
00085   r_0             = 1.1  * fermi;
00086 
00087   // The following variable is set to true if
00088   // G4TripathiLightCrossSection::GetCrossSection is going to be called from
00089   // within G4TripathiLightCrossSection::GetCrossSection to check whether the
00090   // cross-section is behaviing anomalously in the low-energy region.
00091 
00092   lowEnergyCheck = false;
00093 }

G4TripathiLightCrossSection::~G4TripathiLightCrossSection (  ) 

Definition at line 96 of file G4TripathiLightCrossSection.cc.

00097 {
00098   //
00099   // Destructor just needs to delete the pointer to the G4WilsonRadius object.
00100   //
00101   delete theWilsonRadius;
00102 }


Member Function Documentation

G4double G4TripathiLightCrossSection::GetElementCrossSection ( const G4DynamicParticle theProjectile,
G4int  Z,
const G4Material mat = 0 
) [virtual]

Reimplemented from G4VCrossSectionDataSet.

Definition at line 125 of file G4TripathiLightCrossSection.cc.

References G4Pow::A13(), G4lrint(), G4DynamicParticle::Get4Momentum(), G4NistManager::GetAtomicMassAmu(), G4ParticleDefinition::GetBaryonNumber(), G4DynamicParticle::GetDefinition(), G4Pow::GetInstance(), G4DynamicParticle::GetKineticEnergy(), G4NucleiProperties::GetNuclearMass(), G4ParticleDefinition::GetPDGCharge(), G4WilsonRadius::GetWilsonRMSRadius(), G4NistManager::Instance(), G4INCL::Math::pi, G4DynamicParticle::SetKineticEnergy(), SetLowEnergyCheck(), and G4Pow::Z13().

00127 {
00128   // Initialise the result.
00129   G4double result = 0.0;
00130 
00131   // Get details of the projectile and target (nucleon number, atomic number,
00132   // kinetic enery and energy/nucleon.
00133 
00134   G4double xAT= G4NistManager::Instance()->GetAtomicMassAmu(ZT);
00135   G4int    AT = G4lrint(xAT);
00136   G4double EA = theProjectile->GetKineticEnergy()/MeV;
00137   G4int    AP = theProjectile->GetDefinition()->GetBaryonNumber();
00138   G4double xAP= G4double(AP);
00139   G4double ZP = G4lrint(theProjectile->GetDefinition()->GetPDGCharge()/eplus);
00140   G4double E  = EA / xAP;
00141 
00142   G4Pow* g4pow = G4Pow::GetInstance();
00143   
00144   G4double AT13 = g4pow->Z13(AT);
00145   G4double AP13 = g4pow->Z13(AP);
00146 
00147   // Determine target mass and energy within the centre-of-mass frame.
00148 
00149   G4double mT = G4NucleiProperties::GetNuclearMass(AT, ZT);
00150   G4LorentzVector pT(0.0, 0.0, 0.0, mT);
00151   G4LorentzVector pP(theProjectile->Get4Momentum());
00152   pT += pP;
00153   G4double E_cm = (pT.mag()-mT-pP.m())/MeV;
00154 
00155   //G4cout << G4endl;
00156   //G4cout << "### EA= " << EA << " ZT= " << ZT << " AT= " << AT 
00157   //     << "  ZP= " << ZP << " AP= " << AP << " E_cm= " << E_cm 
00158   //     << " Elim= " << (0.8 + 0.04*ZT)*xAP << G4endl;
00159 
00160   if (E_cm <= 0.0) { return 0.; }
00161   if (E_cm <= (0.8 + 0.04*ZT)*xAP && !lowEnergyCheck) { return 0.; }
00162   
00163   G4double E_cm13 = g4pow->A13(E_cm);
00164 
00165   // Determine nuclear radii.  Note that the r_p and r_T are defined differently
00166   // from Wilson et al.
00167 
00168   G4double r_rms_p = theWilsonRadius->GetWilsonRMSRadius(xAP);
00169   G4double r_rms_t = theWilsonRadius->GetWilsonRMSRadius(xAT);
00170 
00171   G4double r_p = 1.29*r_rms_p;
00172   G4double r_t = 1.29*r_rms_t;
00173 
00174   G4double Radius = (r_p + r_t)/fermi + 1.2*(AT13 + AP13)/E_cm13;
00175 
00176   G4double B = 1.44 * ZP * ZT / Radius;
00177 
00178   // Now determine other parameters associated with the parametric
00179   // formula, depending upon the projectile and target.
00180 
00181   G4double T1 = 0.0;
00182   G4double D  = 0.0;
00183   G4double G  = 0.0;
00184 
00185   if ((AT==1 && ZT==1) || (AP==1 && ZP==1)) {
00186     T1 = 23.0;
00187     D  = 1.85 + 0.16/(1+std::exp((500.0-E)/200.0));
00188 
00189   } else if ((AT==1 && ZT==0) || (AP==1 && ZP==0)) {
00190     T1 = 18.0;
00191     D  = 1.85 + 0.16/(1+std::exp((500.0-E)/200.0));
00192 
00193   } else if ((AT==2 && ZT==1) || (AP==2 && ZP==1)) {
00194     T1 = 23.0;
00195     D  = 1.65 + 0.1/(1+std::exp((500.0-E)/200.0));
00196 
00197   } else if ((AT==3 && ZT==2) || (AP==3 && ZP==2)) {
00198     T1 = 40.0;
00199     D  = 1.55;
00200 
00201   } else if (AP==4 && ZP==2) {
00202     if      (AT==4 && ZT==2) {T1 = 40.0; G = 300.0;}
00203     else if (ZT==4)          {T1 = 25.0; G = 300.0;}
00204     else if (ZT==7)          {T1 = 40.0; G = 500.0;}
00205     else if (ZT==13)         {T1 = 25.0; G = 300.0;}
00206     else if (ZT==26)         {T1 = 40.0; G = 300.0;}
00207     else                     {T1 = 40.0; G = 75.0;}
00208     D = 2.77 - 8.0E-3*AT + 1.8E-5*AT*AT-0.8/(1.0+std::exp((250.0-E)/G));
00209   }
00210   else if (AT==4 && ZT==2) {
00211     if      (AP==4 && ZP==2) {T1 = 40.0; G = 300.0;}
00212     else if (ZP==4)          {T1 = 25.0; G = 300.0;}
00213     else if (ZP==7)          {T1 = 40.0; G = 500.0;}
00214     else if (ZP==13)         {T1 = 25.0; G = 300.0;}
00215     else if (ZP==26)         {T1 = 40.0; G = 300.0;}
00216     else                     {T1 = 40.0; G = 75.0;}
00217     D = 2.77 - 8.0E-3*AP + 1.8E-5*AP*AP-0.8/(1.0+std::exp((250.0-E)/G));
00218   }
00219 
00220   // C_E, S, deltaE, X1, S_L and X_m correspond directly with the original
00221   // formulae of Tripathi et al in his report.
00222   //G4cout << "E= " << E << " T1= " << T1 << "  AP= " << AP << " ZP= " << ZP 
00223   //     << " AT= " << AT << " ZT= " << ZT << G4endl;
00224   G4double C_E = D*(1.0-std::exp(-E/T1)) -
00225                  0.292*std::exp(-E/792.0)*std::cos(0.229*std::pow(E,0.453));
00226 
00227   G4double S = AP13*AT13/(AP13 + AT13);
00228 
00229   G4double deltaE = 0.0;
00230   G4double X1     = 0.0;
00231   if (AT >= AP)
00232   {
00233     deltaE = 1.85*S + 0.16*S/E_cm13 - C_E + 0.91*(AT-2*ZT)*ZP/(xAT*xAP);
00234     X1     = 2.83 - 3.1E-2*AT + 1.7E-4*AT*AT;
00235   }
00236   else
00237   {
00238     deltaE = 1.85*S + 0.16*S/E_cm13 - C_E + 0.91*(AP-2*ZP)*ZT/(xAT*xAP);
00239     X1     = 2.83 - 3.1E-2*AP + 1.7E-4*AP*AP;
00240   }
00241   G4double S_L = 1.2 + 1.6*(1.0-std::exp(-E/15.0));
00242   //JMQ 241110 bug fixed 
00243   G4double X_m = 1.0 - X1*std::exp(-E/(X1*S_L));
00244 
00245   //G4cout << "deltaE= " << deltaE << "  X1= " << X1 << " S_L= " << S_L << " X_m= " << X_m << G4endl;
00246 
00247   // R_c is also highly dependent upon the A and Z of the projectile and
00248   // target.
00249 
00250   G4double R_c = 1.0;
00251   if (AP==1 && ZP==1)
00252   {
00253     if      (AT==2 && ZT==1) R_c = 13.5;
00254     else if (AT==3 && ZT==2) R_c = 21.0;
00255     else if (AT==4 && ZT==2) R_c = 27.0;
00256     else if (ZT==3)          R_c = 2.2;
00257   }
00258   else if (AT==1 && ZT==1)
00259   {
00260     if      (AP==2 && ZP==1) R_c = 13.5;
00261     else if (AP==3 && ZP==2) R_c = 21.0;
00262     else if (AP==4 && ZP==2) R_c = 27.0;
00263     else if (ZP==3)          R_c = 2.2;
00264   }
00265   else if (AP==2 && ZP==1)
00266   {
00267     if       (AT==2 && ZT==1) R_c = 13.5;
00268     else if (AT==4 && ZT==2)  R_c = 13.5;
00269     else if (AT==12 && ZT==6) R_c = 6.0;
00270   }
00271   else if (AT==2 && ZT==1)
00272   {
00273     if       (AP==2 && ZP==1) R_c = 13.5;
00274     else if (AP==4 && ZP==2)  R_c = 13.5;
00275     else if (AP==12 && ZP==6) R_c = 6.0;
00276   }
00277   else if ((AP==4 && ZP==2 && (ZT==73 || ZT==79)) ||
00278            (AT==4 && ZT==2 && (ZP==73 || ZP==79))) R_c = 0.6;
00279 
00280   // Find the total cross-section.  Check that it's value is positive, and if
00281   // the energy is less that 10 MeV/nuc, find out if the cross-section is
00282   // increasing with decreasing energy.  If so this is a sign that the function
00283   // is behaving badly at low energies, and the cross-section should be
00284   // set to zero.
00285 
00286   G4double xr = r_0*(AT13 + AP13 + deltaE);
00287   result = pi * xr * xr * (1.0 - R_c*B/E_cm) * X_m;
00288   //G4cout << "       result= " << result << " E= " << E << "  check= "<< lowEnergyCheck << G4endl;
00289   if (result < 0.0) {
00290     result = 0.0;
00291 
00292   } else if (!lowEnergyCheck && E < 6.0) {
00293     G4double f  = 0.95;
00294     G4DynamicParticle slowerProjectile = *theProjectile;
00295     slowerProjectile.SetKineticEnergy(f * EA * MeV);
00296 
00297     G4bool savelowenergy = lowEnergyCheck;
00298     SetLowEnergyCheck(true);
00299     G4double resultp = GetElementCrossSection(&slowerProjectile, ZT);
00300     SetLowEnergyCheck(savelowenergy);
00301     //G4cout << "           newres= " << resultp << " f*EA= " << f*EA << G4endl;  
00302     if (resultp > result) { result = 0.0; }
00303   }
00304 
00305   return result;
00306 }

G4bool G4TripathiLightCrossSection::IsElementApplicable ( const G4DynamicParticle theProjectile,
G4int  Z,
const G4Material  
) [virtual]

Reimplemented from G4VCrossSectionDataSet.

Definition at line 106 of file G4TripathiLightCrossSection.cc.

References G4lrint(), G4ParticleDefinition::GetBaryonNumber(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetKineticEnergy(), G4ParticleDefinition::GetPDGCharge(), and G4NistManager::Instance().

00108 {
00109   G4bool result = false;
00110   G4int AT = G4lrint(G4NistManager::Instance()->GetAtomicMassAmu(ZT));
00111   G4int ZP = G4lrint(theProjectile->GetDefinition()->GetPDGCharge()/eplus);
00112   G4int AP = theProjectile->GetDefinition()->GetBaryonNumber();
00113   if (theProjectile->GetKineticEnergy()/AP < 10.0*GeV &&
00114       ((AT==1 && ZT==1) || (AP==1 && ZP==1) ||
00115        (AT==1 && ZT==0) || (AP==1 && ZP==0) ||
00116        (AT==2 && ZT==1) || (AP==2 && ZP==1) ||
00117        (AT==3 && ZT==2) || (AP==3 && ZP==2) ||
00118        (AT==4 && ZT==2) || (AP==4 && ZP==2))) { result = true; }
00119   return result;
00120 }

void G4TripathiLightCrossSection::SetLowEnergyCheck ( G4bool   )  [inline]

Definition at line 107 of file G4TripathiLightCrossSection.hh.

Referenced by GetElementCrossSection().

00108 {
00109   lowEnergyCheck = aLowEnergyCheck;
00110 }


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