Changeset 1340 for trunk/source/processes/hadronic/cross_sections/src/G4BGGNucleonElasticXS.cc

Timestamp:

Nov 5, 2010, 3:45:55 PM (14 years ago)

Author:

garnier

Message:

update ti head

File:

• trunk/source/processes/hadronic/cross_sections/src/G4BGGNucleonElasticXS.cc (modified) (8 diffs)

trunk/source/processes/hadronic/cross_sections/src/G4BGGNucleonElasticXS.cc

@@ -24 +24 @@
 // ********************************************************************
 //
-// $Id: G4BGGNucleonElasticXS.cc,v 1.7 2009/11/19 19:40:45 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-04-beta-01 $
+// $Id: G4BGGNucleonElasticXS.cc,v 1.8 2010/10/12 06:02:41 dennis Exp $
+// GEANT4 tag $Name: hadr-cross-V09-03-12 $
 //
 // -------------------------------------------------------------------
@@ -51 +51 @@
 #include "G4NistManager.hh"
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 
 G4BGGNucleonElasticXS::G4BGGNucleonElasticXS(const G4ParticleDefinition*) 
@@ -75 +74 @@
 }
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-G4double G4BGGNucleonElasticXS::GetIsoZACrossSection(const G4DynamicParticle* dp, 
-                                                       G4double Z,
-                                                       G4double A, 
-                                                       G4double)
+
+G4double
+G4BGGNucleonElasticXS::GetZandACrossSection(const G4DynamicParticle* dp,
+                                            G4int Z, G4int A, G4double)
 {
   G4double cross = 0.0;
   G4double ekin = dp->GetKineticEnergy();
-  G4int iz = G4int(Z);
-  if(iz > 92) iz = 92;
+  if(Z > 92) Z = 92;
 
   if(ekin <= fLowEnergy) {
-    cross = theCoulombFac[iz];
+    cross = theCoulombFac[Z];
     if(isProton) { cross *= CoulombFactor(ekin, A); }
 
-  } else if(iz == 1) {
-    if( A < 1.5) {
+  } else if(Z == 1) {
+    if( A < 2) {
       fHadron->GetHadronNucleonXscNS(dp, G4Proton::Proton());
       cross = fHadron->GetElasticHadronNucleonXsc();
@@ -102 +98 @@
     }
   } else if(ekin > fGlauberEnergy) {
-    cross = theGlauberFac[iz]*fGlauber->GetElasticGlauberGribov(dp, Z, A);
+    cross = theGlauberFac[Z]*fGlauber->GetElasticGlauberGribov(dp, Z, A);
   } else {
     cross = fNucleon->GetElasticCrossSection(dp, Z, A);
@@ -159 +155 @@
   G4NistManager* nist = G4NistManager::Instance();
 
-  G4double A, csup, csdn;
+  G4double csup, csdn;
+  G4int A;
 
   if(verboseLevel > 0) G4cout << "### G4BGGNucleonElasticXS::Initialise for "
@@ -167 +164 @@
 
     G4double Z = G4double(iz);
-    A = nist->GetAtomicMassAmu(iz);
-
-    csup = fGlauber->GetElasticGlauberGribov(&dp, Z, A);
-    csdn = fNucleon->GetElasticCrossSection(&dp, Z, A);
+    A = G4lrint(nist->GetAtomicMassAmu(iz));
+
+    csup = fGlauber->GetElasticGlauberGribov(&dp, iz, A);
+    csdn = fNucleon->GetElasticCrossSection(&dp, iz, A);
 
     theGlauberFac[iz] = csdn/csup;
@@ -179 +176 @@
   fHadron->GetHadronNucleonXscNS(&dp, G4Proton::Proton());
   theCoulombFac[1] = fHadron->GetElasticHadronNucleonXsc();
-  if(isProton) { theCoulombFac[1] /= CoulombFactor(fLowEnergy, 1.0); }
+  if(isProton) { theCoulombFac[1] /= CoulombFactor(fLowEnergy, 1); }
 
   for(G4int iz=2; iz<93; iz++) {
 
     G4double Z = G4double(iz);
-    A = nist->GetAtomicMassAmu(iz);
-
-    theCoulombFac[iz] = fNucleon->GetElasticCrossSection(&dp, Z, A);
+    A = G4lrint(nist->GetAtomicMassAmu(iz));
+
+    theCoulombFac[iz] = fNucleon->GetElasticCrossSection(&dp, iz, A);
     if(isProton) { theCoulombFac[iz] /= CoulombFactor(fLowEnergy, A); }
     if(verboseLevel > 0) G4cout << "Z= " << Z <<  "  A= " << A 
@@ -195 +192 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 
-G4double G4BGGNucleonElasticXS::CoulombFactor(G4double kinEnergy, G4double A)
-{
-  G4double res= 0.0;
+G4double G4BGGNucleonElasticXS::CoulombFactor(G4double kinEnergy, G4int A)
+{
+  G4double res = 0.0;
   if(kinEnergy <= DBL_MIN) return res;
-  else if(A < 1.5) return kinEnergy*kinEnergy;
+  else if(A < 2) return kinEnergy*kinEnergy;
 
   G4double elog = std::log10(kinEnergy/GeV);
+  G4double aa = A;
 
   // from G4ProtonInelasticCrossSection
-  G4double f1 = 8.0  - 8.0/A - 0.008*A;
-  G4double f2 = 2.34 - 5.4/A - 0.0028*A;
+  G4double f1 = 8.0  - 8.0/aa - 0.008*aa;
+  G4double f2 = 2.34 - 5.4/aa - 0.0028*aa;
 
   res = 1.0/(1.0 + std::exp(-f1*(elog + f2)));
  
-  f1 = 5.6 - 0.016*A;
-  f2 = 1.37 + 1.37/A;
-  res *= ( 1.0 + (0.8 + 18./A - 0.002*A)/(1.0 + std::exp(f1*(elog + f2))));
+  f1 = 5.6 - 0.016*aa;
+  f2 = 1.37 + 1.37/aa;
+  res *= ( 1.0 + (0.8 + 18./aa - 0.002*aa)/(1.0 + std::exp(f1*(elog + f2))));
   return res;  
 }
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-
-
+