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

#include <G4BGGNucleonInelasticXS.hh>

Inheritance diagram for G4BGGNucleonInelasticXS:
Collaboration diagram for G4BGGNucleonInelasticXS:

Public Member Functions

 G4BGGNucleonInelasticXS (const G4ParticleDefinition *)
 
virtual ~G4BGGNucleonInelasticXS ()
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=0, const G4Material *mat=0)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4double GetIsoCrossSection (const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void CrossSectionDescription (std::ostream &) const
 
void SetLowestCrossSection (G4double val)
 
- Public Member Functions inherited from G4VCrossSectionDataSet
 G4VCrossSectionDataSet (const G4String &nam="")
 
virtual ~G4VCrossSectionDataSet ()
 
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 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 64 of file G4BGGNucleonInelasticXS.hh.

Constructor & Destructor Documentation

G4BGGNucleonInelasticXS::G4BGGNucleonInelasticXS ( const G4ParticleDefinition p)

Definition at line 65 of file G4BGGNucleonInelasticXS.cc.

66  : G4VCrossSectionDataSet("Barashenkov-Glauber")
67 {
68  verboseLevel = 0;
69  fGlauberEnergy = 91.*GeV;
70  fLowEnergy = 14.*MeV;
71  fHighEnergy = 5.*GeV;
72  fSAIDHighEnergyLimit = 1.3*GeV;
73  fLowestXSection = millibarn;
74  for (G4int i = 0; i < 93; ++i) {
75  theGlauberFac[i] = 0.0;
76  theCoulombFac[i] = 0.0;
77  theA[i] = 1;
78  }
79  fNucleon = 0;
80  fGlauber = 0;
81  fHadron = 0;
82  fSAID = 0;
83 
84  particle = p;
85  theProton= G4Proton::Proton();
86  isProton = false;
87  isInitialized = false;
88 }
G4VCrossSectionDataSet(const G4String &nam="")
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
static G4Proton * Proton()
Definition: G4Proton.cc:93
static constexpr double GeV
Definition: G4SIunits.hh:217
static constexpr double MeV
Definition: G4SIunits.hh:214
G4bool isInitialized()
static constexpr double millibarn
Definition: G4SIunits.hh:106

Here is the call graph for this function:

G4BGGNucleonInelasticXS::~G4BGGNucleonInelasticXS ( )
virtual

Definition at line 92 of file G4BGGNucleonInelasticXS.cc.

93 {
94  delete fSAID;
95  delete fHadron;
96  // The cross section registry will delete fNucleon
97  delete fGlauber;
98 }

Member Function Documentation

void G4BGGNucleonInelasticXS::BuildPhysicsTable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 201 of file G4BGGNucleonInelasticXS.cc.

202 {
203  if(&p == theProton || &p == G4Neutron::Neutron()) {
204  particle = &p;
205  } else {
206  G4cout << "### G4BGGNucleonInelasticXS WARNING: is not applicable to "
207  << p.GetParticleName()
208  << G4endl;
209  throw G4HadronicException(__FILE__, __LINE__,
210  "G4BGGNucleonElasticXS::BuildPhysicsTable is used for wrong particle");
211  return;
212  }
213 
214  if(isInitialized) { return; }
215  isInitialized = true;
216 
218  fGlauber = new G4ComponentGGHadronNucleusXsc();
219  fHadron = new G4HadronNucleonXsc();
220  fSAID = new G4ComponentSAIDTotalXS();
221 
222  fNucleon->BuildPhysicsTable(*particle);
223  fGlauber->BuildPhysicsTable(*particle);
224 
225  if(particle == theProton) {
226  isProton = true;
227  fSAIDHighEnergyLimit = 2*GeV;
228  fHighEnergy = 2*GeV;
229  }
230 
231  G4ThreeVector mom(0.0,0.0,1.0);
232  G4DynamicParticle dp(particle, mom, fGlauberEnergy);
233 
235  G4int A;
236 
237  G4double csup, csdn;
238 
239  if(verboseLevel > 0) {
240  G4cout << "### G4BGGNucleonInelasticXS::Initialise for "
241  << particle->GetParticleName() << G4endl;
242  }
243  for(G4int iz=2; iz<93; iz++) {
244 
245  A = G4lrint(nist->GetAtomicMassAmu(iz));
246  theA[iz] = A;
247 
248  csup = fGlauber->GetInelasticGlauberGribov(&dp, iz, A);
249  csdn = fNucleon->GetElementCrossSection(&dp, iz);
250 
251  theGlauberFac[iz] = csdn/csup;
252  if(verboseLevel > 0) {
253  G4cout << "Z= " << iz << " A= " << A
254  << " GlauberFactor= " << theGlauberFac[iz] << G4endl;
255  }
256  }
257  //const G4Material* mat = 0;
258 
259  dp.SetKineticEnergy(fSAIDHighEnergyLimit);
260  fHadron->GetHadronNucleonXscNS(&dp, theProton);
261  theCoulombFac[0] = fSAIDHighEnergyLimit*
262  (fSAID->GetInelasticIsotopeCrossSection(particle,fSAIDHighEnergyLimit,1,1)
263  /fHadron->GetInelasticHadronNucleonXsc() - 1);
264 
265  //G4cout << "Z=1 E(GeV)= " << fSAIDHighEnergyLimit/GeV
266  // << " xsNS(b)= " << fHadron->GetInelasticHadronNucleonXsc()/barn;
267  fHadron->GetHadronNucleonXscPDG(&dp, theProton);
268  //G4cout << " xsPDG(b)= " << fHadron->GetInelasticHadronNucleonXsc()/barn;
269  //G4cout << " xsSAID(b)= " << fSAID->GetInelasticIsotopeCrossSection(particle,fSAIDHighEnergyLimit,1,1)/barn << G4endl;
270 
271  dp.SetKineticEnergy(fHighEnergy);
272  fHadron->GetHadronNucleonXscPDG(&dp, theProton);
274 
275  //G4cout << "Z=1 E(GeV)= " << fHighEnergy/GeV
276  // << " xsPDG(b)= " << fHadron->GetInelasticHadronNucleonXsc()/barn;
277 
278  fHadron->GetHadronNucleonXscNS(&dp, theProton);
279  theCoulombFac[1] = fHighEnergy*((theCoulombFac[0]/fHighEnergy + 1)
280  *fHadron->GetInelasticHadronNucleonXsc()/x - 1);
281 
282  fHadron->GetHadronNucleonXscNS(&dp, theProton);
283  //G4cout <<" xsNS(b)= "<<fHadron->GetInelasticHadronNucleonXsc()/barn<<G4endl;
284 
285  if(verboseLevel > 0) {
286  G4cout << "Z=1 A=1" << " CoulombFactor[0]= " << theCoulombFac[0]
287  << " CoulombFactor[1]= " << theCoulombFac[1] << G4endl;
288  }
289  theCoulombFac[2] = 1.0;
290 
291  dp.SetKineticEnergy(fLowEnergy);
292  for(G4int iz=3; iz<93; iz++) {
293  theCoulombFac[iz] =
294  fNucleon->GetElementCrossSection(&dp, iz)/CoulombFactor(fLowEnergy, iz);
295 
296  if(verboseLevel > 0) {
297  G4cout << "Z= " << iz << " A= " << theA[iz]
298  << " CoulombFactor= " << theCoulombFac[iz] << G4endl;
299  }
300  }
301 }
G4VCrossSectionDataSet * GetCrossSectionDataSet(const G4String &name, G4bool warning=true)
virtual G4double GetInelasticIsotopeCrossSection(const G4ParticleDefinition *, G4double kinEnergy, G4int, G4int)
const char * p
Definition: xmltok.h:285
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
virtual G4double GetElementCrossSection(const G4DynamicParticle *aParticle, G4int Z, const G4Material *mat=0)
tuple x
Definition: test.py:50
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
const G4String & GetParticleName() const
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4ParticleDefinition *)
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
static G4CrossSectionDataSetRegistry * Instance()
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
int G4lrint(double ad)
Definition: templates.hh:163
G4double GetAtomicMassAmu(const G4String &symb) const
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4ParticleDefinition *)
#define G4endl
Definition: G4ios.hh:61
G4double GetInelasticGlauberGribov(const G4DynamicParticle *, G4int Z, G4int A)
double G4double
Definition: G4Types.hh:76
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
G4bool isInitialized()
G4double GetInelasticHadronNucleonXsc()

Here is the call graph for this function:

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

Reimplemented from G4VCrossSectionDataSet.

Definition at line 346 of file G4BGGNucleonInelasticXS.cc.

347 {
348  outFile << "The Barashenkov-Glauber-Gribov cross section calculates inelastic\n"
349  << "scattering of protons and neutrons from nuclei using the\n"
350  << "Barashenkov parameterization below 91 GeV and the Glauber-Gribov\n"
351  << "parameterization above 91 GeV. It uses the G4HadronNucleonXsc\n"
352  << "cross section component for hydrogen targets, and the\n"
353  << "G4ComponentGGHadronNucleusXsc component for other targets.\n";
354 }
G4double G4BGGNucleonInelasticXS::GetElementCrossSection ( const G4DynamicParticle dp,
G4int  Z,
const G4Material mat = 0 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 121 of file G4BGGNucleonInelasticXS.cc.

123 {
124  // this method should be called only for Z > 1
125 
126  G4double cross = 0.0;
127  G4double ekin = dp->GetKineticEnergy();
128  G4int Z = ZZ;
129  if(1 == Z) {
130  cross = 1.0115*GetIsoCrossSection(dp,1,1);
131  } else if(2 == Z) {
132  if(ekin > fGlauberEnergy) {
133  cross = theGlauberFac[Z]*fGlauber->GetInelasticGlauberGribov(dp, Z, theA[Z]);
134  } else {
135  cross = fNucleon->GetElementCrossSection(dp, Z);
136  }
137 
138  } else {
139  if(Z > 92) { Z = 92; }
140 
141  if(ekin <= fLowEnergy) {
142  cross = theCoulombFac[Z]*CoulombFactor(ekin, Z);
143  } else if(ekin > fGlauberEnergy) {
144  cross = theGlauberFac[Z]*fGlauber->GetInelasticGlauberGribov(dp, Z, theA[Z]);
145  } else {
146  cross = fNucleon->GetElementCrossSection(dp, Z);
147  }
148  }
149 
150  if(verboseLevel > 1) {
151  G4cout << "G4BGGNucleonInelasticXS::GetCrossSection for "
152  << dp->GetDefinition()->GetParticleName()
153  << " Ekin(GeV)= " << dp->GetKineticEnergy()/CLHEP::GeV
154  << " in nucleus Z= " << Z << " A= " << theA[Z]
155  << " XS(b)= " << cross/barn
156  << G4endl;
157  }
158  //AR-18Dec2013 if(cross <= fLowestXSection) { cross = 0.0; }
159  return cross;
160 }
G4double GetKineticEnergy() const
virtual G4double GetElementCrossSection(const G4DynamicParticle *aParticle, G4int Z, const G4Material *mat=0)
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
static constexpr double GeV
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
#define G4endl
Definition: G4ios.hh:61
G4double GetInelasticGlauberGribov(const G4DynamicParticle *, G4int Z, G4int A)
double G4double
Definition: G4Types.hh:76
static constexpr double barn
Definition: G4SIunits.hh:105

Here is the call graph for this function:

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

Reimplemented from G4VCrossSectionDataSet.

Definition at line 165 of file G4BGGNucleonInelasticXS.cc.

170 {
171  // this method should be called only for Z = 1
172 
173  G4double cross = 0.0;
174  G4double ekin = dp->GetKineticEnergy();
175 
176  if(ekin <= fSAIDHighEnergyLimit) {
177  cross = fSAID->GetInelasticIsotopeCrossSection(particle, ekin, 1, 1);
178  } else if(ekin < fHighEnergy) {
179  fHadron->GetHadronNucleonXscNS(dp, theProton);
180  cross = (theCoulombFac[0]/ekin + 1)*fHadron->GetInelasticHadronNucleonXsc();
181  } else {
182  fHadron->GetHadronNucleonXscPDG(dp, theProton);
183  cross = (theCoulombFac[1]/ekin + 1)*fHadron->GetInelasticHadronNucleonXsc();
184  }
185  cross *= A;
186 
187  if(verboseLevel > 1) {
188  G4cout << "G4BGGNucleonInelasticXS::GetCrossSection for "
189  << dp->GetDefinition()->GetParticleName()
190  << " Ekin(GeV)= " << dp->GetKineticEnergy()/CLHEP::GeV
191  << " in nucleus Z= " << Z << " A= " << A
192  << " XS(b)= " << cross/barn
193  << G4endl;
194  }
195  //AR-18Dec2013 if(cross <= fLowestXSection) { cross = 0.0; }
196  return cross;
197 }
virtual G4double GetInelasticIsotopeCrossSection(const G4ParticleDefinition *, G4double kinEnergy, G4int, G4int)
G4double GetKineticEnergy() const
G4ParticleDefinition * GetDefinition() const
const G4String & GetParticleName() const
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4ParticleDefinition *)
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
static constexpr double GeV
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4ParticleDefinition *)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static constexpr double barn
Definition: G4SIunits.hh:105
G4double GetInelasticHadronNucleonXsc()

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4BGGNucleonInelasticXS::IsElementApplicable ( const G4DynamicParticle ,
G4int  Z,
const G4Material mat = 0 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Reimplemented in G4ParticleHPBGGNucleonInelasticXS.

Definition at line 102 of file G4BGGNucleonInelasticXS.cc.

104 {
105  return true;
106 }
G4bool G4BGGNucleonInelasticXS::IsIsoApplicable ( const G4DynamicParticle ,
G4int  Z,
G4int  A,
const G4Element elm = 0,
const G4Material mat = 0 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Reimplemented in G4ParticleHPBGGNucleonInelasticXS.

Definition at line 110 of file G4BGGNucleonInelasticXS.cc.

114 {
115  return (1 == Z);
116 }
void G4BGGNucleonInelasticXS::SetLowestCrossSection ( G4double  val)
inline

Definition at line 125 of file G4BGGNucleonInelasticXS.hh.

126 {
127  fLowestXSection = val;
128 }

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