G4ProtonInelasticCrossSection Class Reference

#include <G4ProtonInelasticCrossSection.hh>

Inheritance diagram for G4ProtonInelasticCrossSection:

G4VCrossSectionDataSet

Public Member Functions

 G4ProtonInelasticCrossSection ()
 ~G4ProtonInelasticCrossSection ()
virtual G4bool IsElementApplicable (const G4DynamicParticle *aPart, G4int Z, const G4Material *)
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *)
G4double GetProtonCrossSection (G4double kineticEnergy, G4int Z)

Detailed Description

Definition at line 52 of file G4ProtonInelasticCrossSection.hh.


Constructor & Destructor Documentation

G4ProtonInelasticCrossSection::G4ProtonInelasticCrossSection (  ) 

Definition at line 41 of file G4ProtonInelasticCrossSection.cc.

References G4NistManager::Instance(), and G4Proton::Proton().

00042   : G4VCrossSectionDataSet("Axen-Wellisch"),thEnergy(19.8*CLHEP::GeV)
00043 {
00044   nist = G4NistManager::Instance();
00045   theProton = G4Proton::Proton();
00046 }

G4ProtonInelasticCrossSection::~G4ProtonInelasticCrossSection (  ) 

Definition at line 48 of file G4ProtonInelasticCrossSection.cc.

00049 {}


Member Function Documentation

G4double G4ProtonInelasticCrossSection::GetElementCrossSection ( const G4DynamicParticle ,
G4int  Z,
const G4Material  
) [virtual]

Reimplemented from G4VCrossSectionDataSet.

Definition at line 59 of file G4ProtonInelasticCrossSection.cc.

References G4DynamicParticle::GetKineticEnergy(), and GetProtonCrossSection().

00062 {
00063   return GetProtonCrossSection(aPart->GetKineticEnergy(), Z);
00064 }

G4double G4ProtonInelasticCrossSection::GetProtonCrossSection ( G4double  kineticEnergy,
G4int  Z 
)

Definition at line 66 of file G4ProtonInelasticCrossSection.cc.

References G4lrint(), G4NistManager::GetAtomicMassAmu(), and G4INCL::Math::pi.

Referenced by GetElementCrossSection(), and G4IonProtonCrossSection::GetElementCrossSection().

00068 {
00069   if(kineticEnergy <= 0.0) { return 0.0; }
00070  
00071   // constant cross section above ~20GeV
00072   if (kineticEnergy > thEnergy) { kineticEnergy = thEnergy; }
00073 
00074   G4double a = nist->GetAtomicMassAmu(Z);
00075   G4double a13 = std::pow(a,-0.3333333333);
00076   G4int nOfNeutrons = G4lrint(a) - Z;
00077   kineticEnergy /=GeV;
00078   G4double alog10E = std::log10(kineticEnergy);
00079 
00080   const G4double nuleonRadius=1.36E-15;
00081   const G4double fac=CLHEP::pi*nuleonRadius*nuleonRadius;
00082 
00083   G4double b0   = 2.247-0.915*(1 - a13);
00084   G4double fac1 = b0*(1 - a13);
00085   G4double fac2 = 1.;
00086   if(nOfNeutrons > 1) { fac2=std::log((G4double(nOfNeutrons))); }
00087   G4double crossSection = 1.0E31*fac*fac2*(1. + 1./a13 - fac1);
00088 
00089   // high energy correction
00090   crossSection *= (1 - 0.15*std::exp(-kineticEnergy))/(1.0 - 0.0007*a);
00091 
00092   // first try on low energies: rise
00093   G4double ff1= 0.70-0.002*a;  // slope of the drop at medium energies.
00094   G4double ff2= 1.00+1/a;  // start of the slope.
00095   G4double ff3= 0.8+18/a-0.002*a; // stephight
00096 
00097   G4double ff4= 1.0 - (1.0/(1+std::exp(-8*ff1*(alog10E + 1.37*ff2))));
00098 
00099   crossSection *= (1 + ff3*ff4);
00100 
00101   // low energy return to zero
00102 
00103   ff1=1. - 1./a - 0.001*a;       // slope of the rise
00104   ff2=1.17 - 2.7/a - 0.0014*a;   // start of the rise
00105  
00106   ff4=-8.*ff1*(alog10E + 2.0*ff2);
00107  
00108   crossSection *= millibarn/(1. + std::exp(ff4));
00109   return crossSection;
00110 }

G4bool G4ProtonInelasticCrossSection::IsElementApplicable ( const G4DynamicParticle aPart,
G4int  Z,
const G4Material  
) [virtual]

Reimplemented from G4VCrossSectionDataSet.

Definition at line 52 of file G4ProtonInelasticCrossSection.cc.

References G4DynamicParticle::GetDefinition().

00055 {
00056   return ((1 < Z) && (aPart->GetDefinition() == theProton));
00057 }


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