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

#include <G4ComponentGGNuclNuclXsc.hh>

Inheritance diagram for G4ComponentGGNuclNuclXsc:
Collaboration diagram for G4ComponentGGNuclNuclXsc:

Public Member Functions

 G4ComponentGGNuclNuclXsc ()
 
virtual ~G4ComponentGGNuclNuclXsc ()
 
virtual G4double GetTotalIsotopeCrossSection (const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4int A)
 
virtual G4double GetTotalElementCrossSection (const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
 
virtual G4double GetInelasticIsotopeCrossSection (const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4int A)
 
virtual G4double GetInelasticElementCrossSection (const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
 
virtual G4double GetElasticElementCrossSection (const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
 
virtual G4double GetElasticIsotopeCrossSection (const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4int A)
 
virtual G4double ComputeQuasiElasticRatio (const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4int A)
 
G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, const G4Material *)
 
G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *)
 
G4double GetZandACrossSection (const G4DynamicParticle *, G4int Z, G4int A)
 
G4double GetCoulombBarier (const G4DynamicParticle *, G4double Z, G4double A, G4double pR, G4double tR)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
virtual void CrossSectionDescription (std::ostream &) const
 
G4double GetRatioSD (const G4DynamicParticle *, G4double At, G4double Zt)
 
G4double GetRatioQE (const G4DynamicParticle *, G4double At, G4double Zt)
 
G4double GetHadronNucleonXsc (const G4DynamicParticle *, const G4Element *)
 
G4double GetHadronNucleonXsc (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double GetHadronNucleonXscPDG (const G4ParticleDefinition *, G4double sMand, const G4ParticleDefinition *)
 
G4double GetHadronNucleonXscNS (const G4ParticleDefinition *, G4double pTkin, const G4ParticleDefinition *)
 
G4double GetHNinelasticXscVU (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double CalculateEcmValue (const G4double, const G4double, const G4double)
 
G4double CalcMandelstamS (const G4double, const G4double, const G4double)
 
G4double GetElasticGlauberGribov (const G4DynamicParticle *, G4int Z, G4int A)
 
G4double GetInelasticGlauberGribov (const G4DynamicParticle *, G4int Z, G4int A)
 
G4double GetTotalGlauberGribovXsc ()
 
G4double GetElasticGlauberGribovXsc ()
 
G4double GetInelasticGlauberGribovXsc ()
 
G4double GetProductionGlauberGribovXsc ()
 
G4double GetDiffractionGlauberGribovXsc ()
 
G4double GetRadiusConst ()
 
G4double GetNucleusRadius (const G4DynamicParticle *, const G4Element *)
 
G4double GetNucleusRadius (G4double Zt, G4double At)
 
G4double GetNucleusRadiusGG (G4double At)
 
G4double GetNucleusRadiusDE (G4double Zt, G4double At)
 
G4double GetNucleusRadiusRMS (G4double Zt, G4double At)
 
void SetEnergyLowerLimit (G4double E)
 
- Public Member Functions inherited from G4VComponentCrossSection
 G4VComponentCrossSection (const G4String &nam="")
 
virtual ~G4VComponentCrossSection ()
 
G4double GetTotalElementCrossSection (const G4ParticleDefinition *, G4double kinEnergy, const G4Element *)
 
G4double GetInelasticElementCrossSection (const G4ParticleDefinition *, G4double kinEnergy, const G4Element *)
 
G4double GetElasticElementCrossSection (const G4ParticleDefinition *, G4double kinEnergy, const G4Element *)
 
virtual void Description () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
G4double GetMinKinEnergy () const
 
void SetMinKinEnergy (G4double value)
 
G4double GetMaxKinEnergy () const
 
void SetMaxKinEnergy (G4double value)
 
const G4StringGetName () const
 

Detailed Description

Definition at line 50 of file G4ComponentGGNuclNuclXsc.hh.

Constructor & Destructor Documentation

G4ComponentGGNuclNuclXsc::G4ComponentGGNuclNuclXsc ( )

Definition at line 40 of file G4ComponentGGNuclNuclXsc.cc.

41  : G4VComponentCrossSection("Glauber-Gribov nucleus nucleus"),
42 // fUpperLimit(100000*GeV),
43  fLowerLimit(0.1*MeV),
44  fRadiusConst(1.08*fermi), // 1.1, 1.3 ?
45  fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0), fProductionXsc(0.0),
46  fDiffractionXsc(0.0),
47  cacheDP(G4Proton::Proton(),G4ParticleMomentum(1.,0,0),0),
48  dProton(G4Proton::Proton(),G4ParticleMomentum(1.,0,0),0),
49  dNeutron(G4Neutron::Neutron(),G4ParticleMomentum(1.,0,0),0)
50 // , fHadronNucleonXsc(0.0)
51 {
52  theProton = G4Proton::Proton();
53  theNeutron = G4Neutron::Neutron();
54  hnXsc = new G4HadronNucleonXsc();
55 }
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
G4VComponentCrossSection(const G4String &nam="")
static constexpr double MeV
Definition: G4SIunits.hh:214
static constexpr double fermi
Definition: G4SIunits.hh:103
G4ThreeVector G4ParticleMomentum

Here is the call graph for this function:

G4ComponentGGNuclNuclXsc::~G4ComponentGGNuclNuclXsc ( )
virtual

Definition at line 58 of file G4ComponentGGNuclNuclXsc.cc.

59 {
60  delete hnXsc;
61 }

Member Function Documentation

virtual void G4ComponentGGNuclNuclXsc::BuildPhysicsTable ( const G4ParticleDefinition )
inlinevirtual

Reimplemented from G4VComponentCrossSection.

Definition at line 111 of file G4ComponentGGNuclNuclXsc.hh.

112  {}
G4double G4ComponentGGNuclNuclXsc::CalcMandelstamS ( const G4double  mp,
const G4double  mt,
const G4double  Plab 
)

Definition at line 857 of file G4ComponentGGNuclNuclXsc.cc.

860 {
861  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
862  G4double sMand = mp*mp + mt*mt + 2*Elab*mt ;
863 
864  return sMand;
865 }
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4ComponentGGNuclNuclXsc::CalculateEcmValue ( const G4double  mp,
const G4double  mt,
const G4double  Plab 
)

Definition at line 840 of file G4ComponentGGNuclNuclXsc.cc.

843 {
844  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
845  G4double Ecm = std::sqrt ( mp * mp + mt * mt + 2 * Elab * mt );
846  // G4double Pcm = Plab * mt / Ecm;
847  // G4double KEcm = std::sqrt ( Pcm * Pcm + mp * mp ) - mp;
848 
849  return Ecm ; // KEcm;
850 }
double G4double
Definition: G4Types.hh:76
G4double G4ComponentGGNuclNuclXsc::ComputeQuasiElasticRatio ( const G4ParticleDefinition aParticle,
G4double  kinEnergy,
G4int  Z,
G4int  A 
)
virtual

Reimplemented from G4VComponentCrossSection.

Definition at line 137 of file G4ComponentGGNuclNuclXsc.cc.

140 {
141  cacheDP.SetDefinition(aParticle);
142  cacheDP.SetKineticEnergy(kinEnergy);
143  fInelasticXsc = GetZandACrossSection(&cacheDP, Z, A);
144  G4double ratio = 0.;
145 
146  if(fInelasticXsc > 0.)
147  {
148  ratio = (fInelasticXsc - fProductionXsc)/fInelasticXsc;
149  if(ratio < 0.) ratio = 0.;
150  }
151  return ratio;
152 }
double A(double temperature)
void SetKineticEnergy(G4double aEnergy)
double G4double
Definition: G4Types.hh:76
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

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

Definition at line 157 of file G4ComponentGGNuclNuclXsc.cc.

158 {
159  outFile << "G4ComponentGGNuclNuclXsc calculates total, inelastic and\n"
160  << "elastic cross sections for nucleus-nucleus collisions using\n"
161  << "the Glauber model with Gribov corrections. It is valid for\n"
162  << "all incident energies above 100 keV./n";
163 }
virtual void G4ComponentGGNuclNuclXsc::DumpPhysicsTable ( const G4ParticleDefinition )
inlinevirtual

Reimplemented from G4VComponentCrossSection.

Definition at line 115 of file G4ComponentGGNuclNuclXsc.hh.

116  {G4cout << "G4NuclNuclCrossSection: uses Glauber-Gribov formula"<<G4endl;}
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
G4double G4ComponentGGNuclNuclXsc::GetCoulombBarier ( const G4DynamicParticle aParticle,
G4double  Z,
G4double  A,
G4double  pR,
G4double  tR 
)

Definition at line 294 of file G4ComponentGGNuclNuclXsc.cc.

296 {
297  G4double ratio;
298  G4double pZ = aParticle->GetDefinition()->GetPDGCharge();
299 
300  G4double pTkin = aParticle->GetKineticEnergy();
301  // G4double pPlab = aParticle->GetTotalMomentum();
302  G4double pM = aParticle->GetDefinition()->GetPDGMass();
303  // G4double tM = tZ*proton_mass_c2 + (tA-tZ)*neutron_mass_c2; // ~ 1% accuracy
305  G4double pElab = pTkin + pM;
306  G4double totEcm = std::sqrt(pM*pM + tM*tM + 2.*pElab*tM);
307  // G4double pPcm = pPlab*tM/totEcm;
308  // G4double pTcm = std::sqrt(pM*pM + pPcm*pPcm) - pM;
309  G4double totTcm = totEcm - pM -tM;
310 
312  bC /= pR + tR;
313  bC /= 2.; // 4., 2. parametrisation cof ??? vmg
314 
315  // G4cout<<"pTkin = "<<pTkin/GeV<<"; pPlab = "
316  // <<pPlab/GeV<<"; bC = "<<bC/GeV<<"; pTcm = "<<pTcm/GeV<<G4endl;
317 
318  if( totTcm <= bC ) ratio = 0.;
319  else ratio = 1. - bC/totTcm;
320 
321  // if(ratio < DBL_MIN) ratio = DBL_MIN;
322  if( ratio < 0.) ratio = 0.;
323 
324  // G4cout <<"ratio = "<<ratio<<G4endl;
325  return ratio;
326 }
G4double GetKineticEnergy() const
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
G4IonTable * GetIonTable() const
G4double GetIonMass(G4int Z, G4int A, G4int L=0, G4int lvl=0) const
Definition: G4IonTable.cc:1335
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetDiffractionGlauberGribovXsc ( )
inline

Definition at line 141 of file G4ComponentGGNuclNuclXsc.hh.

141 { return fDiffractionXsc; };
G4double G4ComponentGGNuclNuclXsc::GetElasticElementCrossSection ( const G4ParticleDefinition aParticle,
G4double  kinEnergy,
G4int  Z,
G4double  A 
)
virtual

Implements G4VComponentCrossSection.

Definition at line 113 of file G4ComponentGGNuclNuclXsc.cc.

116 {
117  cacheDP.SetDefinition(aParticle);
118  cacheDP.SetKineticEnergy(kinEnergy);
119  fInelasticXsc = GetZandACrossSection(&cacheDP, Z, G4int(A));
120  return fElasticXsc;
121 }
int G4int
Definition: G4Types.hh:78
double A(double temperature)
void SetKineticEnergy(G4double aEnergy)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetElasticGlauberGribov ( const G4DynamicParticle dp,
G4int  Z,
G4int  A 
)
inline

Definition at line 177 of file G4ComponentGGNuclNuclXsc.hh.

179 {
180  GetZandACrossSection(dp, Z, A);
181  return fElasticXsc;
182 }
double A(double temperature)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetElasticGlauberGribovXsc ( )
inline

Definition at line 138 of file G4ComponentGGNuclNuclXsc.hh.

138 { return fElasticXsc; };
G4double G4ComponentGGNuclNuclXsc::GetElasticIsotopeCrossSection ( const G4ParticleDefinition aParticle,
G4double  kinEnergy,
G4int  Z,
G4int  A 
)
virtual

Implements G4VComponentCrossSection.

Definition at line 125 of file G4ComponentGGNuclNuclXsc.cc.

128 {
129  cacheDP.SetDefinition(aParticle);
130  cacheDP.SetKineticEnergy(kinEnergy);
131  fInelasticXsc = GetZandACrossSection(&cacheDP, Z, A);
132  return fElasticXsc;
133 }
double A(double temperature)
void SetKineticEnergy(G4double aEnergy)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetElementCrossSection ( const G4DynamicParticle aParticle,
G4int  Z,
const G4Material  
)

Definition at line 187 of file G4ComponentGGNuclNuclXsc.cc.

189 {
190  G4int A = G4lrint(G4NistManager::Instance()->GetAtomicMassAmu(Z));
191  return GetZandACrossSection(aParticle, Z, A);
192 }
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
double A(double temperature)
int G4lrint(double ad)
Definition: templates.hh:163
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetHadronNucleonXsc ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 419 of file G4ComponentGGNuclNuclXsc.cc.

421 {
422  G4int At = G4lrint(anElement->GetN()); // number of nucleons
423  G4int Zt = G4lrint(anElement->GetZ()); // number of protons
424  return GetHadronNucleonXsc(aParticle, At, Zt);
425 }
G4double GetN() const
Definition: G4Element.hh:135
G4double GetZ() const
Definition: G4Element.hh:131
int G4int
Definition: G4Types.hh:78
G4double GetHadronNucleonXsc(const G4DynamicParticle *, const G4Element *)
int G4lrint(double ad)
Definition: templates.hh:163

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetHadronNucleonXsc ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 438 of file G4ComponentGGNuclNuclXsc.cc.

440 {
441  G4double xsection = 0.;
442 
444  GetIonTable()->GetIonMass(Zt, At);
445  targ_mass = 0.939*GeV; // ~mean neutron and proton ???
446 
447  G4double proj_mass = aParticle->GetMass();
448  G4double proj_momentum = aParticle->GetMomentum().mag();
449  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
450 
451  sMand /= GeV*GeV; // in GeV for parametrisation
452  proj_momentum /= GeV;
453  const G4ParticleDefinition* pParticle = aParticle->GetDefinition();
454 
455  if(pParticle == theNeutron) // as proton ???
456  {
457  xsection = G4double(At)*(21.70*G4Pow::GetInstance()->powA(sMand,0.0808) + 56.08*G4Pow::GetInstance()->powA(sMand,-0.4525));
458  }
459  else if(pParticle == theProton)
460  {
461  xsection = G4double(At)*(21.70*G4Pow::GetInstance()->powA(sMand,0.0808) + 56.08*G4Pow::GetInstance()->powA(sMand,-0.4525));
462  }
463 
464  xsection *= millibarn;
465  return xsection;
466 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4ParticleDefinition * GetDefinition() const
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
G4double GetMass() const
static G4ParticleTable * GetParticleTable()
static constexpr double GeV
Definition: G4SIunits.hh:217
double G4double
Definition: G4Types.hh:76
double mag() const
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4ThreeVector GetMomentum() const

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetHadronNucleonXscNS ( const G4ParticleDefinition pParticle,
G4double  pTkin,
const G4ParticleDefinition tParticle 
)

Definition at line 534 of file G4ComponentGGNuclNuclXsc.cc.

537 {
538  G4double xsection(0);
539  // G4double Delta; DHW 19 May 2011: variable set but not used
540  G4double A0, B0;
541  G4double hpXscv(0);
542  G4double hnXscv(0);
543 
544  G4double targ_mass = tParticle->GetPDGMass();
545  G4double proj_mass = pParticle->GetPDGMass();
546 
547  G4double proj_energy = proj_mass + pTkin;
548  G4double proj_momentum = std::sqrt(pTkin*(pTkin+2*proj_mass));
549 
550  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
551 
552  sMand /= GeV*GeV; // in GeV for parametrisation
553  proj_momentum /= GeV;
554  proj_energy /= GeV;
555  proj_mass /= GeV;
556 
557  // General PDG fit constants
558 
559  // G4double s0 = 5.38*5.38; // in Gev^2
560  // G4double eta1 = 0.458;
561  // G4double eta2 = 0.458;
562  // G4double B = 0.308;
563 
564  if( proj_momentum >= 373.)
565  {
566  return GetHadronNucleonXscPDG(pParticle,sMand,tParticle);
567  }
568  else if( proj_momentum >= 10. ) // high energy: pp = nn = np
569  // if( proj_momentum >= 2.)
570  {
571  // Delta = 1.; // DHW 19 May 2011: variable set but not used
572  // if (proj_energy < 40.) Delta = 0.916+0.0021*proj_energy;
573 
574  //AR-12Aug2016 if (proj_momentum >= 10.) {
575  B0 = 7.5;
576  A0 = 100. - B0*G4Log(3.0e7);
577 
578  xsection = A0 + B0*G4Log(proj_energy) - 11
579  + 103*G4Pow::GetInstance()->powA(2*0.93827*proj_energy + proj_mass*proj_mass+
580  0.93827*0.93827,-0.165); // mb
581  //AR-12Aug2016 }
582  }
583  else // low energy pp = nn != np
584  {
585  if(pParticle == tParticle) // pp or nn // nn to be pp
586  {
587  if( proj_momentum < 0.73 )
588  {
589  hnXscv = 23 + 50*( G4Pow::GetInstance()->powA( G4Log(0.73/proj_momentum), 3.5 ) );
590  }
591  else if( proj_momentum < 1.05 )
592  {
593  hnXscv = 23 + 40*(G4Log(proj_momentum/0.73))*
594  (G4Log(proj_momentum/0.73));
595  }
596  else // if( proj_momentum < 10. )
597  {
598  hnXscv = 39.0 +
599  75*(proj_momentum - 1.2)/(G4Pow::GetInstance()->powA(proj_momentum,3.0) + 0.15);
600  }
601  xsection = hnXscv;
602  }
603  else // pn to be np
604  {
605  if( proj_momentum < 0.8 )
606  {
607  hpXscv = 33+30*G4Pow::GetInstance()->powA(G4Log(proj_momentum/1.3),4.0);
608  }
609  else if( proj_momentum < 1.4 )
610  {
611  hpXscv = 33+30*G4Pow::GetInstance()->powA(G4Log(proj_momentum/0.95),2.0);
612  }
613  else // if( proj_momentum < 10. )
614  {
615  hpXscv = 33.3+
616  20.8*(G4Pow::GetInstance()->powA(proj_momentum,2.0)-1.35)/
617  (G4Pow::GetInstance()->powA(proj_momentum,2.50)+0.95);
618  }
619  xsection = hpXscv;
620  }
621  }
622  xsection *= millibarn; // parametrised in mb
623  return xsection;
624 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double GetPDGMass() const
static constexpr double GeV
Definition: G4SIunits.hh:217
double G4double
Definition: G4Types.hh:76
G4double GetHadronNucleonXscPDG(const G4ParticleDefinition *, G4double sMand, const G4ParticleDefinition *)
static constexpr double millibarn
Definition: G4SIunits.hh:106

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetHadronNucleonXscPDG ( const G4ParticleDefinition pParticle,
G4double  sMand,
const G4ParticleDefinition tParticle 
)

Definition at line 477 of file G4ComponentGGNuclNuclXsc.cc.

480 {
481  G4double xsection = 0.;
482  // G4bool pORn = (tParticle == theProton || nucleon == theNeutron );
483  G4bool proton = (tParticle == theProton);
484  G4bool neutron = (tParticle == theNeutron);
485 
486  // General PDG fit constants
487 
488  G4double s0 = 5.38*5.38; // in Gev^2
489  G4double eta1 = 0.458;
490  G4double eta2 = 0.458;
491  G4double B = 0.308;
492 
493  // const G4ParticleDefinition* pParticle = aParticle->GetDefinition();
494 
495  if(pParticle == theNeutron) // proton-neutron fit
496  {
497  if ( proton )
498  {
499  xsection = ( 35.80 + B*G4Pow::GetInstance()->powA(G4Log(sMand/s0),2.)
500  + 40.15*G4Pow::GetInstance()->powA(sMand,-eta1) - 30.*G4Pow::GetInstance()->powA(sMand,-eta2));
501  }
502  if ( neutron )
503  {
504  xsection = (35.45 + B*G4Pow::GetInstance()->powA(G4Log(sMand/s0),2.)
505  + 42.53*G4Pow::GetInstance()->powA(sMand,-eta1) - 33.34*G4Pow::GetInstance()->powA(sMand,-eta2)); // pp for nn
506  }
507  }
508  else if(pParticle == theProton)
509  {
510  if ( proton )
511  {
512  xsection = (35.45 + B*G4Pow::GetInstance()->powA(G4Log(sMand/s0),2.)
513  + 42.53*G4Pow::GetInstance()->powA(sMand,-eta1) - 33.34*G4Pow::GetInstance()->powA(sMand,-eta2));
514 
515  }
516  if ( neutron )
517  {
518  xsection = (35.80 + B*G4Pow::GetInstance()->powA(G4Log(sMand/s0),2.)
519  + 40.15*G4Pow::GetInstance()->powA(sMand,-eta1) - 30.*G4Pow::GetInstance()->powA(sMand,-eta2));
520  }
521  }
522  xsection *= millibarn; // parametrised in mb
523  return xsection;
524 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
double B(double temperature)
bool G4bool
Definition: G4Types.hh:79
G4double G4Log(G4double x)
Definition: G4Log.hh:230
double G4double
Definition: G4Types.hh:76
static constexpr double millibarn
Definition: G4SIunits.hh:106

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetHNinelasticXscVU ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 631 of file G4ComponentGGNuclNuclXsc.cc.

633 {
634  G4int PDGcode = aParticle->GetDefinition()->GetPDGEncoding();
635  G4int absPDGcode = std::abs(PDGcode);
636  G4double Elab = aParticle->GetTotalEnergy();
637  // (s - 2*0.88*GeV*GeV)/(2*0.939*GeV)/GeV;
638  G4double Plab = aParticle->GetMomentum().mag();
639  // std::sqrt(Elab * Elab - 0.88);
640 
641  Elab /= GeV;
642  Plab /= GeV;
643 
644  G4double LogPlab = G4Log( Plab );
645  G4double sqrLogPlab = LogPlab * LogPlab;
646 
647  //G4cout<<"Plab = "<<Plab<<G4endl;
648 
649  G4double NumberOfTargetProtons = Zt;
650  G4double NumberOfTargetNucleons = At;
651  G4double NumberOfTargetNeutrons = NumberOfTargetNucleons - NumberOfTargetProtons;
652 
653  if(NumberOfTargetNeutrons < 0.) NumberOfTargetNeutrons = 0.;
654 
655  G4double Xtotal = 0., Xelastic = 0., Xinelastic =0.;
656 
657  if( absPDGcode > 1000 ) //------Projectile is baryon --------
658  {
659  G4double XtotPP = 48.0 + 0. *G4Pow::GetInstance()->powA(Plab, 0. ) +
660  0.522*sqrLogPlab - 4.51*LogPlab;
661 
662  G4double XtotPN = 47.3 + 0. *G4Pow::GetInstance()->powA(Plab, 0. ) +
663  0.513*sqrLogPlab - 4.27*LogPlab;
664 
665  G4double XelPP = 11.9 + 26.9*G4Pow::GetInstance()->powA(Plab,-1.21) +
666  0.169*sqrLogPlab - 1.85*LogPlab;
667 
668  G4double XelPN = 11.9 + 26.9*G4Pow::GetInstance()->powA(Plab,-1.21) +
669  0.169*sqrLogPlab - 1.85*LogPlab;
670 
671  Xtotal = ( NumberOfTargetProtons * XtotPP +
672  NumberOfTargetNeutrons * XtotPN );
673 
674  Xelastic = ( NumberOfTargetProtons * XelPP +
675  NumberOfTargetNeutrons * XelPN );
676  }
677 
678  Xinelastic = Xtotal - Xelastic;
679  if(Xinelastic < 0.) Xinelastic = 0.;
680 
681  return Xinelastic*= millibarn;
682 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4double GetTotalEnergy() const
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
G4double G4Log(G4double x)
Definition: G4Log.hh:230
static constexpr double GeV
Definition: G4SIunits.hh:217
double G4double
Definition: G4Types.hh:76
double mag() const
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4ThreeVector GetMomentum() const

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetInelasticElementCrossSection ( const G4ParticleDefinition aParticle,
G4double  kinEnergy,
G4int  Z,
G4double  A 
)
virtual

Implements G4VComponentCrossSection.

Definition at line 101 of file G4ComponentGGNuclNuclXsc.cc.

104 {
105  cacheDP.SetDefinition(aParticle);
106  cacheDP.SetKineticEnergy(kinEnergy);
107  fInelasticXsc = GetZandACrossSection(&cacheDP, Z, G4int(A));
108  return fInelasticXsc;
109 }
int G4int
Definition: G4Types.hh:78
double A(double temperature)
void SetKineticEnergy(G4double aEnergy)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetInelasticGlauberGribov ( const G4DynamicParticle dp,
G4int  Z,
G4int  A 
)
inline

Definition at line 187 of file G4ComponentGGNuclNuclXsc.hh.

189 {
190  GetZandACrossSection(dp, Z, A);
191  return fInelasticXsc;
192 }
double A(double temperature)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetInelasticGlauberGribovXsc ( )
inline

Definition at line 139 of file G4ComponentGGNuclNuclXsc.hh.

139 { return fInelasticXsc; };
G4double G4ComponentGGNuclNuclXsc::GetInelasticIsotopeCrossSection ( const G4ParticleDefinition aParticle,
G4double  kinEnergy,
G4int  Z,
G4int  A 
)
virtual

Implements G4VComponentCrossSection.

Definition at line 89 of file G4ComponentGGNuclNuclXsc.cc.

92 {
93  cacheDP.SetDefinition(aParticle);
94  cacheDP.SetKineticEnergy(kinEnergy);
95  fInelasticXsc = GetZandACrossSection(&cacheDP, Z, A);
96  return fInelasticXsc;
97 }
double A(double temperature)
void SetKineticEnergy(G4double aEnergy)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetNucleusRadius ( const G4DynamicParticle ,
const G4Element anElement 
)

Definition at line 689 of file G4ComponentGGNuclNuclXsc.cc.

691 {
692  G4double At = anElement->GetN();
693  G4double oneThird = 1.0/3.0;
694  G4double cubicrAt = G4Pow::GetInstance()->powA (At, oneThird);
695 
696  G4double R; // = fRadiusConst*cubicrAt;
697  R = fRadiusConst*cubicrAt;
698 
699  G4double meanA = 21.;
700  G4double tauA1 = 40.;
701  G4double tauA2 = 10.;
702  G4double tauA3 = 5.;
703 
704  G4double a1 = 0.85;
705  G4double b1 = 1. - a1;
706 
707  G4double b2 = 0.3;
708  G4double b3 = 4.;
709 
710  if (At > 20.) // 20.
711  {
712  R *= ( a1 + b1*G4Exp( -(At - meanA)/tauA1) );
713  }
714  else if (At > 3.5)
715  {
716  R *= ( 1.0 + b2*( 1. - G4Exp( (At - meanA)/tauA2) ) );
717  }
718  else
719  {
720  R *= ( 1.0 + b3*( 1. - G4Exp( (At - meanA)/tauA3) ) );
721  }
722 
723  return R;
724 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4double GetN() const
Definition: G4Element.hh:135
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
double G4double
Definition: G4Types.hh:76
const G4double oneThird

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetNucleusRadius ( G4double  Zt,
G4double  At 
)

Definition at line 731 of file G4ComponentGGNuclNuclXsc.cc.

732 {
733  G4double R;
734  R = GetNucleusRadiusDE(Zt,At);
735  // R = GetNucleusRadiusRMS(Zt,At);
736 
737  return R;
738 }
double G4double
Definition: G4Types.hh:76
G4double GetNucleusRadiusDE(G4double Zt, G4double At)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetNucleusRadiusDE ( G4double  Zt,
G4double  At 
)

Definition at line 771 of file G4ComponentGGNuclNuclXsc.cc.

772 {
773  // algorithm from diffuse-elastic
774 
775  G4double R, r0, a11, a12, a13, a2, a3;
776 
777  a11 = 1.26; // 1.08, 1.16
778  a12 = 1.; // 1.08, 1.16
779  a13 = 1.12; // 1.08, 1.16
780  a2 = 1.1;
781  a3 = 1.;
782 
783  // Special rms radii for light nucleii
784 
785  if (A < 50.)
786  {
787  if (std::abs(A-1.) < 0.5) return 0.89*fermi; // p
788  else if(std::abs(A-2.) < 0.5) return 2.13*fermi; // d
789  else if(std::abs(Z-1.) < 0.5 && std::abs(A-3.) < 0.5) return 1.80*fermi; // t
790 
791  else if(std::abs(Z-2.) < 0.5 && std::abs(A-3.) < 0.5) return 1.96*fermi; // He3
792  else if(std::abs(Z-2.) < 0.5 && std::abs(A-4.) < 0.5) return 1.68*fermi; // He4
793 
794  else if(std::abs(Z-3.) < 0.5) return 2.40*fermi; // Li7
795  else if(std::abs(Z-4.) < 0.5) return 2.51*fermi; // Be9
796 
797  else if( 10. < A && A <= 16. ) r0 = a11*( 1 - G4Pow::GetInstance()->powA(A, -2./3.) )*fermi; // 1.08*fermi;
798  else if( 15. < A && A <= 20. ) r0 = a12*( 1 - G4Pow::GetInstance()->powA(A, -2./3.) )*fermi;
799  else if( 20. < A && A <= 30. ) r0 = a13*( 1 - G4Pow::GetInstance()->powA(A, -2./3.) )*fermi;
800  else r0 = a2*fermi;
801 
802  R = r0*G4Pow::GetInstance()->powA( A, 1./3. );
803  }
804  else
805  {
806  r0 = a3*fermi;
807 
808  R = r0*G4Pow::GetInstance()->powA(A, 0.27);
809  }
810  return R;
811 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
double A(double temperature)
double G4double
Definition: G4Types.hh:76
static constexpr double fermi
Definition: G4SIunits.hh:103

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetNucleusRadiusGG ( G4double  At)

Definition at line 743 of file G4ComponentGGNuclNuclXsc.cc.

744 {
745  G4double oneThird = 1.0/3.0;
746  G4double cubicrAt = G4Pow::GetInstance()->powA (At, oneThird);
747 
748  G4double R; // = fRadiusConst*cubicrAt;
749  R = fRadiusConst*cubicrAt;
750 
751  G4double meanA = 20.;
752  G4double tauA = 20.;
753 
754  if ( At > 20.) // 20.
755  {
756  R *= ( 0.8 + 0.2*G4Exp( -(At - meanA)/tauA) );
757  }
758  else
759  {
760  R *= ( 1.0 + 0.1*( 1. - G4Exp( (At - meanA)/tauA) ) );
761  }
762 
763  return R;
764 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
double G4double
Definition: G4Types.hh:76
const G4double oneThird

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetNucleusRadiusRMS ( G4double  Zt,
G4double  At 
)

Definition at line 819 of file G4ComponentGGNuclNuclXsc.cc.

820 {
821 
822  if (std::abs(A-1.) < 0.5) return 0.89*fermi; // p
823  else if(std::abs(A-2.) < 0.5) return 2.13*fermi; // d
824  else if(std::abs(Z-1.) < 0.5 && std::abs(A-3.) < 0.5) return 1.80*fermi; // t
825 
826  else if(std::abs(Z-2.) < 0.5 && std::abs(A-3.) < 0.5) return 1.96*fermi; // He3
827  else if(std::abs(Z-2.) < 0.5 && std::abs(A-4.) < 0.5) return 1.68*fermi; // He4
828 
829  else if(std::abs(Z-3.) < 0.5) return 2.40*fermi; // Li7
830  else if(std::abs(Z-4.) < 0.5) return 2.51*fermi; // Be9
831 
832  else return 1.24*G4Pow::GetInstance()->powA(A, 0.28 )*fermi; // A > 9
833 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
double A(double temperature)
static constexpr double fermi
Definition: G4SIunits.hh:103

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetProductionGlauberGribovXsc ( )
inline

Definition at line 140 of file G4ComponentGGNuclNuclXsc.hh.

140 { return fProductionXsc; };
G4double G4ComponentGGNuclNuclXsc::GetRadiusConst ( )
inline

Definition at line 142 of file G4ComponentGGNuclNuclXsc.hh.

142 { return fRadiusConst; };
G4double G4ComponentGGNuclNuclXsc::GetRatioQE ( const G4DynamicParticle aParticle,
G4double  At,
G4double  Zt 
)

Definition at line 374 of file G4ComponentGGNuclNuclXsc.cc.

375 {
376  G4double sigma, cofInelastic = 2.4, cofTotal = 2.0, nucleusSquare, ratio;
377 
378  G4double pZ = aParticle->GetDefinition()->GetPDGCharge();
379  G4double pA = aParticle->GetDefinition()->GetBaryonNumber();
380 
381  G4double pTkin = aParticle->GetKineticEnergy();
382  pTkin /= pA;
383 
384  G4double pN = pA - pZ;
385  if( pN < 0. ) pN = 0.;
386 
387  G4double tN = tA - tZ;
388  if( tN < 0. ) tN = 0.;
389 
390  G4double tR = GetNucleusRadius(tZ,tA);
391  G4double pR = GetNucleusRadius(pZ,pA);
392 
393  sigma = (pZ*tZ+pN*tN)*GetHadronNucleonXscNS(theProton, pTkin, theProton) +
394  (pZ*tN+pN*tZ)*GetHadronNucleonXscNS(theProton, pTkin, theNeutron);
395 
396  nucleusSquare = cofTotal*pi*( pR*pR + tR*tR ); // basically 2piRR
397  ratio = sigma/nucleusSquare;
398  fInelasticXsc = nucleusSquare*G4Log(1. + cofInelastic*ratio)/cofInelastic;
399 
400  // sigma = GetHNinelasticXsc(aParticle, tA, tZ);
401  ratio = sigma/nucleusSquare;
402  fProductionXsc = nucleusSquare*G4Log(1. + cofInelastic*ratio)/cofInelastic;
403 
404  if (fInelasticXsc > fProductionXsc) ratio = (fInelasticXsc-fProductionXsc)/fInelasticXsc;
405  else ratio = 0.;
406  if ( ratio < 0. ) ratio = 0.;
407 
408  return ratio;
409 }
G4double GetKineticEnergy() const
G4double GetHadronNucleonXscNS(const G4ParticleDefinition *, G4double pTkin, const G4ParticleDefinition *)
G4double GetNucleusRadius(const G4DynamicParticle *, const G4Element *)
G4ParticleDefinition * GetDefinition() const
G4double G4Log(G4double x)
Definition: G4Log.hh:230
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetRatioSD ( const G4DynamicParticle aParticle,
G4double  At,
G4double  Zt 
)

Definition at line 334 of file G4ComponentGGNuclNuclXsc.cc.

335 {
336  G4double sigma, cofInelastic = 2.4, cofTotal = 2.0, nucleusSquare, ratio;
337 
338  G4double pZ = aParticle->GetDefinition()->GetPDGCharge();
339  G4double pA = aParticle->GetDefinition()->GetBaryonNumber();
340 
341  G4double pTkin = aParticle->GetKineticEnergy();
342  pTkin /= pA;
343 
344  G4double pN = pA - pZ;
345  if( pN < 0. ) pN = 0.;
346 
347  G4double tN = tA - tZ;
348  if( tN < 0. ) tN = 0.;
349 
350  G4double tR = GetNucleusRadius(tZ,tA);
351  G4double pR = GetNucleusRadius(pZ,pA);
352 
353  sigma = (pZ*tZ+pN*tN)*GetHadronNucleonXscNS(theProton, pTkin, theProton) +
354  (pZ*tN+pN*tZ)*GetHadronNucleonXscNS(theProton, pTkin, theNeutron);
355 
356  nucleusSquare = cofTotal*pi*( pR*pR + tR*tR ); // basically 2piRR
357  ratio = sigma/nucleusSquare;
358  fInelasticXsc = nucleusSquare*G4Log(1. + cofInelastic*ratio)/cofInelastic;
359  G4double difratio = ratio/(1.+ratio);
360 
361  fDiffractionXsc = 0.5*nucleusSquare*( difratio - G4Log( 1. + difratio ) );
362 
363  if (fInelasticXsc > 0.) ratio = fDiffractionXsc/fInelasticXsc;
364  else ratio = 0.;
365 
366  return ratio;
367 }
G4double GetKineticEnergy() const
G4double GetHadronNucleonXscNS(const G4ParticleDefinition *, G4double pTkin, const G4ParticleDefinition *)
G4double GetNucleusRadius(const G4DynamicParticle *, const G4Element *)
G4ParticleDefinition * GetDefinition() const
G4double G4Log(G4double x)
Definition: G4Log.hh:230
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetTotalElementCrossSection ( const G4ParticleDefinition aParticle,
G4double  kinEnergy,
G4int  Z,
G4double  A 
)
virtual

Implements G4VComponentCrossSection.

Definition at line 77 of file G4ComponentGGNuclNuclXsc.cc.

80 {
81  cacheDP.SetDefinition(aParticle);
82  cacheDP.SetKineticEnergy(kinEnergy);
83  fInelasticXsc = GetZandACrossSection(&cacheDP, Z, G4int(A));
84  return fTotalXsc;
85 }
int G4int
Definition: G4Types.hh:78
double A(double temperature)
void SetKineticEnergy(G4double aEnergy)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetTotalGlauberGribovXsc ( )
inline

Definition at line 137 of file G4ComponentGGNuclNuclXsc.hh.

137 { return fTotalXsc; };
G4double G4ComponentGGNuclNuclXsc::GetTotalIsotopeCrossSection ( const G4ParticleDefinition aParticle,
G4double  kinEnergy,
G4int  Z,
G4int  A 
)
virtual

Implements G4VComponentCrossSection.

Definition at line 65 of file G4ComponentGGNuclNuclXsc.cc.

68 {
69  cacheDP.SetDefinition(aParticle);
70  cacheDP.SetKineticEnergy(kinEnergy);
71  fInelasticXsc = GetZandACrossSection(&cacheDP, Z, A);
72  return fTotalXsc;
73 }
double A(double temperature)
void SetKineticEnergy(G4double aEnergy)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetZandACrossSection(const G4DynamicParticle *, G4int Z, G4int A)

Here is the call graph for this function:

G4double G4ComponentGGNuclNuclXsc::GetZandACrossSection ( const G4DynamicParticle aParticle,
G4int  Z,
G4int  A 
)

Definition at line 204 of file G4ComponentGGNuclNuclXsc.cc.

206 {
207  G4double xsection;
208  G4double sigma;
209  G4double cofInelastic = 2.4;
210  G4double cofTotal = 2.0;
211  G4double nucleusSquare;
212  G4double cB;
213  G4double ratio;
214 
215  G4double pZ = aParticle->GetDefinition()->GetPDGCharge();
216  G4double pA = aParticle->GetDefinition()->GetBaryonNumber();
217 
218  G4double pTkin = aParticle->GetKineticEnergy();
219  pTkin /= pA;
220 
221  G4double pN = pA - pZ;
222  if( pN < 0. ) pN = 0.;
223 
224  G4double tN = tA - tZ;
225  if( tN < 0. ) tN = 0.;
226 
227  G4double tR = GetNucleusRadius( G4double(tZ),G4double(tA) );
228  G4double pR = GetNucleusRadius(pZ,pA);
229 
230  cB = GetCoulombBarier(aParticle, G4double(tZ), G4double(tA), pR, tR);
231 
232  if ( cB > 0. )
233  {
234  dProton.SetKineticEnergy(pTkin);
235  dNeutron.SetKineticEnergy(pTkin);
236 
237  sigma = (pZ*tZ+pN*tN)*hnXsc->GetHadronNucleonXscNS(&dProton, theProton);
238 
239  G4double ppInXsc = hnXsc->GetInelasticHadronNucleonXsc();
240 
241  sigma += (pZ*tN+pN*tZ)*hnXsc->GetHadronNucleonXscNS(&dNeutron, theProton);
242 
243  G4double npInXsc = hnXsc->GetInelasticHadronNucleonXsc();
244 
245  // G4cout<<"ppInXsc = "<<ppInXsc/millibarn<<"; npInXsc = "<<npInXsc/millibarn<<G4endl;
246  // G4cout<<"npTotXsc = "<<hnXsc->GetTotalHadronNucleonXsc()/millibarn<<"; npElXsc = "
247  // <<hnXsc->GetElasticHadronNucleonXsc()/millibarn<<G4endl;
248 
249  nucleusSquare = cofTotal*pi*( pR*pR + tR*tR ); // basically 2piRR
250 
251  ratio = sigma/nucleusSquare;
252  xsection = nucleusSquare*G4Log( 1. + ratio );
253  fTotalXsc = xsection;
254  fTotalXsc *= cB;
255 
256  fInelasticXsc = nucleusSquare*G4Log( 1. + cofInelastic*ratio )/cofInelastic;
257 
258  fInelasticXsc *= cB;
259  fElasticXsc = fTotalXsc - fInelasticXsc;
260 
261  // if (fElasticXsc < DBL_MIN) fElasticXsc = DBL_MIN;
262  /*
263  G4double difratio = ratio/(1.+ratio);
264 
265  fDiffractionXsc = 0.5*nucleusSquare*( difratio - G4Log( 1. + difratio ) );
266  */
267  // production to be checked !!! edit MK xsc
268 
269  //sigma = (pZ*tZ+pN*tN)*GetHadronNucleonXscMK(theProton, pTkin, theProton) +
270  // (pZ*tN+pN*tZ)*GetHadronNucleonXscMK(theProton, pTkin, theNeutron);
271 
272  sigma = (pZ*tZ+pN*tN)*ppInXsc + (pZ*tN+pN*tZ)*npInXsc;
273 
274  ratio = sigma/nucleusSquare;
275  fProductionXsc = nucleusSquare*G4Log( 1. + cofInelastic*ratio )/cofInelastic;
276 
277  if (fElasticXsc < 0.) fElasticXsc = 0.;
278  }
279  else
280  {
281  fInelasticXsc = 0.;
282  fTotalXsc = 0.;
283  fElasticXsc = 0.;
284  fProductionXsc = 0.;
285  }
286  return fInelasticXsc; // xsection;
287 }
G4double GetKineticEnergy() const
G4double GetNucleusRadius(const G4DynamicParticle *, const G4Element *)
G4ParticleDefinition * GetDefinition() const
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4ParticleDefinition *)
void SetKineticEnergy(G4double aEnergy)
G4double G4Log(G4double x)
Definition: G4Log.hh:230
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const
G4double GetCoulombBarier(const G4DynamicParticle *, G4double Z, G4double A, G4double pR, G4double tR)
G4double GetInelasticHadronNucleonXsc()

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 168 of file G4ComponentGGNuclNuclXsc.cc.

170 {
171  G4bool applicable = false;
172  G4double kineticEnergy = aDP->GetKineticEnergy();
173 
174  if (kineticEnergy >= fLowerLimit && Z > 1) applicable = true;
175  return applicable;
176 }
G4double GetKineticEnergy() const
bool G4bool
Definition: G4Types.hh:79
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

void G4ComponentGGNuclNuclXsc::SetEnergyLowerLimit ( G4double  E)
inline

Definition at line 151 of file G4ComponentGGNuclNuclXsc.hh.

151 {fLowerLimit=E;};

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