#include <G4TripathiLightCrossSection.hh>
Inheritance diagram for G4TripathiLightCrossSection:
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) |
Definition at line 83 of file G4TripathiLightCrossSection.hh.
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 }
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().