Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4IonsShenCrossSection Class Reference

#include <G4IonsShenCrossSection.hh>

Inheritance diagram for G4IonsShenCrossSection:
Collaboration diagram for G4IonsShenCrossSection:

Public Member Functions

 G4IonsShenCrossSection ()
 
virtual ~G4IonsShenCrossSection ()
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *aDP, G4int Z, const G4Material *)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *)
 
virtual G4double GetIsoCrossSection (const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
 
virtual void CrossSectionDescription (std::ostream &) const
 
- Public Member Functions inherited from G4VCrossSectionDataSet
 G4VCrossSectionDataSet (const G4String &nam="")
 
virtual ~G4VCrossSectionDataSet ()
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=0, const G4Material *mat=0)
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
virtual G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
virtual G4int GetVerboseLevel () const
 
virtual void SetVerboseLevel (G4int value)
 
G4double GetMinKinEnergy () const
 
void SetMinKinEnergy (G4double value)
 
G4double GetMaxKinEnergy () const
 
void SetMaxKinEnergy (G4double value)
 
const G4StringGetName () const
 

Additional Inherited Members

- Protected Member Functions inherited from G4VCrossSectionDataSet
void SetName (const G4String &)
 
- Protected Attributes inherited from G4VCrossSectionDataSet
G4int verboseLevel
 

Detailed Description

Definition at line 51 of file G4IonsShenCrossSection.hh.

Constructor & Destructor Documentation

G4IonsShenCrossSection::G4IonsShenCrossSection ( )

Definition at line 42 of file G4IonsShenCrossSection.cc.

43  : G4VCrossSectionDataSet("IonsShen"),
44  upperLimit( 10*GeV ),
45 // lowerLimit( 10*MeV ),
46  r0 ( 1.1 )
47 {}
G4VCrossSectionDataSet(const G4String &nam="")
static constexpr double GeV
Definition: G4SIunits.hh:217
G4IonsShenCrossSection::~G4IonsShenCrossSection ( )
virtual

Definition at line 49 of file G4IonsShenCrossSection.cc.

50 {}

Member Function Documentation

void G4IonsShenCrossSection::CrossSectionDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 53 of file G4IonsShenCrossSection.cc.

54 {
55  outFile << "G4IonsShenCrossSection calculates the total reaction cross\n"
56  << "section for nucleus-nucleus scattering using the Shen\n"
57  << "parameterization. It is valid for projectiles and targets of\n"
58  << "all Z, and projectile energies up to 1 TeV/n. Above 10 GeV/n"
59  << "the cross section is constant. Below 10 MeV/n zero cross\n"
60  << "is returned.\n";
61 }
G4double G4IonsShenCrossSection::GetElementCrossSection ( const G4DynamicParticle aParticle,
G4int  Z,
const G4Material  
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 70 of file G4IonsShenCrossSection.cc.

73 {
74  G4int A = G4lrint(G4NistManager::Instance()->GetAtomicMassAmu(Z));
75  return GetIsoCrossSection(aParticle, Z, A);
76 }
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
double A(double temperature)
int G4lrint(double ad)
Definition: templates.hh:163
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)

Here is the call graph for this function:

G4double G4IonsShenCrossSection::GetIsoCrossSection ( const G4DynamicParticle aParticle,
G4int  Z,
G4int  A,
const G4Isotope iso = 0,
const G4Element elm = 0,
const G4Material mat = 0 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 78 of file G4IonsShenCrossSection.cc.

84 {
85  G4double xsection = 0.0;
86 
87  G4int Ap = aParticle->GetDefinition()->GetBaryonNumber();
88  G4int Zp = G4lrint(aParticle->GetDefinition()->GetPDGCharge()/eplus);
89  G4double ke_per_N = aParticle->GetKineticEnergy() / Ap;
90  if ( ke_per_N > upperLimit ) { ke_per_N = upperLimit; }
91 
92  // Apply energy check, if less than lower limit then 0 value is returned
93  //if ( ke_per_N < lowerLimit ) { return xsection; }
94 
95  G4Pow* g4pow = G4Pow::GetInstance();
96 
97  G4double cubicrAt = g4pow->Z13(At);
98  G4double cubicrAp = g4pow->Z13(Ap);
99 
100  G4double Rt = 1.12 * cubicrAt - 0.94 * ( 1.0 / cubicrAt );
101  G4double Rp = 1.12 * cubicrAp - 0.94 * ( 1.0 / cubicrAp );
102 
103  G4double r = Rt + Rp + 3.2; // in fm
104  G4double b = 1.0; // in MeV/fm
105  G4double targ_mass = G4NucleiProperties::GetNuclearMass(At, Zt);
106 
107  G4double proj_mass = aParticle->GetMass();
108  G4double proj_momentum = aParticle->GetMomentum().mag();
109 
110  G4double Ecm = calEcmValue (proj_mass, targ_mass, proj_momentum);
111 
112  G4double B = 1.44 * Zt * Zp / r - b * Rt * Rp / ( Rt + Rp );
113  if(Ecm <= B) { return xsection; }
114 
115  G4double c = calCeValue ( ke_per_N / MeV );
116 
117  G4double R1 = r0 * (cubicrAt + cubicrAp + 1.85*cubicrAt*cubicrAp/(cubicrAt + cubicrAp) - c);
118 
119  G4double R2 = 1.0 * ( At - 2 * Zt ) * Zp / ( Ap * At );
120 
121 
122  G4double R3 = (0.176 / g4pow->A13(Ecm)) * cubicrAt * cubicrAp /(cubicrAt + cubicrAp);
123 
124  G4double R = R1 + R2 + R3;
125 
126  xsection = 10 * pi * R * R * ( 1 - B / Ecm );
127  xsection = xsection * millibarn; // mulitply xsection by millibarn
128 
129  return xsection;
130 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4double GetKineticEnergy() const
Definition: G4Pow.hh:56
G4ParticleDefinition * GetDefinition() const
double B(double temperature)
int G4int
Definition: G4Types.hh:78
tuple b
Definition: test.py:12
G4double Z13(G4int Z) const
Definition: G4Pow.hh:127
G4double GetMass() const
static constexpr double eplus
Definition: G4SIunits.hh:199
G4double A13(G4double A) const
Definition: G4Pow.hh:132
int G4lrint(double ad)
Definition: templates.hh:163
static constexpr double MeV
Definition: G4SIunits.hh:214
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13
G4double GetPDGCharge() const
double mag() const
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4ThreeVector GetMomentum() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4IonsShenCrossSection::IsElementApplicable ( const G4DynamicParticle aDP,
G4int  Z,
const G4Material  
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 63 of file G4IonsShenCrossSection.cc.

65 {
66  return (1 <= aDP->GetDefinition()->GetBaryonNumber());
67 }
G4ParticleDefinition * GetDefinition() const

Here is the call graph for this function:


The documentation for this class was generated from the following files: