Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Static Public Member Functions
G4GlauberGribovCrossSection Class Reference

#include <G4GlauberGribovCrossSection.hh>

Inheritance diagram for G4GlauberGribovCrossSection:
G4VCrossSectionDataSet

Public Member Functions

 G4GlauberGribovCrossSection ()
 
virtual ~G4GlauberGribovCrossSection ()
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *aDP, G4int Z, G4int A, const G4Element *elm=0, 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)
 
G4double GetRatioSD (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double GetRatioQE (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double GetHadronNucleonXsc (const G4DynamicParticle *, const G4Element *)
 
G4double GetHadronNucleonXsc (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double GetHadronNucleonXscPDG (const G4DynamicParticle *, const G4Element *)
 
G4double GetHadronNucleonXscPDG (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double GetHadronNucleonXscNS (const G4DynamicParticle *, const G4Element *)
 
G4double GetHadronNucleonXscNS (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double GetKaonNucleonXscVector (const G4DynamicParticle *, G4int At, G4int Zt)
 
G4double GetHNinelasticXsc (const G4DynamicParticle *, const G4Element *)
 
G4double GetHNinelasticXsc (const G4DynamicParticle *, G4int At, G4int Zt)
 
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 GetNucleusRadius (const G4DynamicParticle *, const G4Element *)
 
G4double GetNucleusRadius (G4int At)
 
virtual void CrossSectionDescription (std::ostream &) const
 
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 GetParticleBarCorTot (const G4ParticleDefinition *theParticle, G4int Z)
 
G4double GetParticleBarCorIn (const G4ParticleDefinition *theParticle, G4int Z)
 
void SetEnergyLowerLimit (G4double E)
 
- Public Member Functions inherited from G4VCrossSectionDataSet
 G4VCrossSectionDataSet (const G4String &nam="")
 
virtual ~G4VCrossSectionDataSet ()
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, 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 G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, 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
 

Static Public Member Functions

static const char * Default_Name ()
 

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 55 of file G4GlauberGribovCrossSection.hh.

Constructor & Destructor Documentation

G4GlauberGribovCrossSection::G4GlauberGribovCrossSection ( )

Definition at line 230 of file G4GlauberGribovCrossSection.cc.

References G4Alpha::Alpha(), G4AntiLambda::AntiLambda(), G4AntiNeutron::AntiNeutron(), G4AntiOmegaMinus::AntiOmegaMinus(), G4AntiProton::AntiProton(), G4AntiSigmaMinus::AntiSigmaMinus(), G4AntiSigmaPlus::AntiSigmaPlus(), G4AntiSigmaZero::AntiSigmaZero(), G4AntiXiMinus::AntiXiMinus(), G4AntiXiZero::AntiXiZero(), G4Deuteron::Deuteron(), G4Gamma::Gamma(), G4He3::He3(), G4KaonMinus::KaonMinus(), G4KaonPlus::KaonPlus(), G4KaonZeroLong::KaonZeroLong(), G4KaonZeroShort::KaonZeroShort(), G4Lambda::Lambda(), G4Neutron::Neutron(), G4OmegaMinus::OmegaMinus(), G4PionMinus::PionMinus(), G4PionPlus::PionPlus(), G4PionZero::PionZero(), G4Proton::Proton(), G4SigmaMinus::SigmaMinus(), G4SigmaPlus::SigmaPlus(), G4SigmaZero::SigmaZero(), G4Triton::Triton(), G4XiMinus::XiMinus(), and G4XiZero::XiZero().

232 // fUpperLimit(100000*GeV),
233  fLowerLimit(10.*MeV),// fLowerLimit(3*GeV),
234  fRadiusConst(1.08*fermi), // 1.1, 1.3 ?
235  fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0), fProductionXsc(0.0),
236  fDiffractionXsc(0.0)
237 // , fHadronNucleonXsc(0.0)
238 {
239  theGamma = G4Gamma::Gamma();
240  theProton = G4Proton::Proton();
241  theNeutron = G4Neutron::Neutron();
242  theAProton = G4AntiProton::AntiProton();
243  theANeutron = G4AntiNeutron::AntiNeutron();
244  thePiPlus = G4PionPlus::PionPlus();
245  thePiMinus = G4PionMinus::PionMinus();
246  thePiZero = G4PionZero::PionZero();
247  theKPlus = G4KaonPlus::KaonPlus();
248  theKMinus = G4KaonMinus::KaonMinus();
250  theK0L = G4KaonZeroLong::KaonZeroLong();
251  theL = G4Lambda::Lambda();
252  theAntiL = G4AntiLambda::AntiLambda();
253  theSPlus = G4SigmaPlus::SigmaPlus();
254  theASPlus = G4AntiSigmaPlus::AntiSigmaPlus();
255  theSMinus = G4SigmaMinus::SigmaMinus();
256  theASMinus = G4AntiSigmaMinus::AntiSigmaMinus();
257  theS0 = G4SigmaZero::SigmaZero();
259  theXiMinus = G4XiMinus::XiMinus();
260  theXi0 = G4XiZero::XiZero();
261  theAXiMinus = G4AntiXiMinus::AntiXiMinus();
262  theAXi0 = G4AntiXiZero::AntiXiZero();
263  theOmega = G4OmegaMinus::OmegaMinus();
264  theAOmega = G4AntiOmegaMinus::AntiOmegaMinus();
265  theD = G4Deuteron::Deuteron();
266  theT = G4Triton::Triton();
267  theA = G4Alpha::Alpha();
268  theHe3 = G4He3::He3();
269 
270  hnXsc = new G4HadronNucleonXsc();
271 }
static G4AntiOmegaMinus * AntiOmegaMinus()
G4VCrossSectionDataSet(const G4String &nam="")
static G4OmegaMinus * OmegaMinus()
static G4KaonZeroLong * KaonZeroLong()
static G4AntiSigmaPlus * AntiSigmaPlus()
static G4SigmaZero * SigmaZero()
Definition: G4SigmaZero.cc:99
static G4KaonMinus * KaonMinus()
Definition: G4KaonMinus.cc:113
static G4AntiSigmaMinus * AntiSigmaMinus()
static G4XiZero * XiZero()
Definition: G4XiZero.cc:106
static G4KaonZeroShort * KaonZeroShort()
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
static G4XiMinus * XiMinus()
Definition: G4XiMinus.cc:106
static G4AntiXiMinus * AntiXiMinus()
static G4Triton * Triton()
Definition: G4Triton.cc:95
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4PionPlus * PionPlus()
Definition: G4PionPlus.cc:98
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
static G4PionZero * PionZero()
Definition: G4PionZero.cc:104
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
static G4SigmaMinus * SigmaMinus()
static G4AntiLambda * AntiLambda()
static G4PionMinus * PionMinus()
Definition: G4PionMinus.cc:98
static G4AntiSigmaZero * AntiSigmaZero()
static G4AntiXiZero * AntiXiZero()
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
static G4SigmaPlus * SigmaPlus()
Definition: G4SigmaPlus.cc:108
static G4Lambda * Lambda()
Definition: G4Lambda.cc:108
static G4KaonPlus * KaonPlus()
Definition: G4KaonPlus.cc:113
static G4He3 * He3()
Definition: G4He3.cc:94
static G4AntiNeutron * AntiNeutron()
G4GlauberGribovCrossSection::~G4GlauberGribovCrossSection ( )
virtual

Definition at line 277 of file G4GlauberGribovCrossSection.cc.

278 {
279  if (hnXsc) delete hnXsc;
280 }

Member Function Documentation

G4double G4GlauberGribovCrossSection::CalcMandelstamS ( const G4double  mp,
const G4double  mt,
const G4double  Plab 
)

Definition at line 1509 of file G4GlauberGribovCrossSection.cc.

Referenced by GetHadronNucleonXsc(), GetHadronNucleonXscNS(), and GetHadronNucleonXscPDG().

1512 {
1513  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1514  G4double sMand = mp*mp + mt*mt + 2*Elab*mt ;
1515 
1516  return sMand;
1517 }
double G4double
Definition: G4Types.hh:76
G4double G4GlauberGribovCrossSection::CalculateEcmValue ( const G4double  mp,
const G4double  mt,
const G4double  Plab 
)

Definition at line 1493 of file G4GlauberGribovCrossSection.cc.

1496 {
1497  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1498  G4double Ecm = std::sqrt ( mp * mp + mt * mt + 2 * Elab * mt );
1499  // G4double Pcm = Plab * mt / Ecm;
1500  // G4double KEcm = std::sqrt ( Pcm * Pcm + mp * mp ) - mp;
1501 
1502  return Ecm ; // KEcm;
1503 }
double G4double
Definition: G4Types.hh:76
void G4GlauberGribovCrossSection::CrossSectionDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 1523 of file G4GlauberGribovCrossSection.cc.

1524 {
1525  outFile << "G4GlauberGribovCrossSection calculates total, inelastic and\n"
1526  << "elastic cross sections for hadron-nucleus cross sections using\n"
1527  << "the Glauber model with Gribov corrections. It is valid for all\n"
1528  << "targets except hydrogen, and for incident p, pbar, n, sigma-,\n"
1529  << "pi+, pi-, K+, K- and gammas with energies above 3 GeV. This is\n"
1530  << "a cross section component which is to be used to build a cross\n"
1531  << "data set.\n";
1532 }
std::ofstream outFile
Definition: GammaRayTel.cc:68
static const char* G4GlauberGribovCrossSection::Default_Name ( )
inlinestatic

Definition at line 62 of file G4GlauberGribovCrossSection.hh.

Referenced by G4BGGNucleonInelasticXS::BuildPhysicsTable(), and G4BGGNucleonElasticXS::BuildPhysicsTable().

62 {return "Glauber-Gribov";}
G4double G4GlauberGribovCrossSection::GetDiffractionGlauberGribovXsc ( )
inline

Definition at line 107 of file G4GlauberGribovCrossSection.hh.

107 { return fDiffractionXsc; };
G4double G4GlauberGribovCrossSection::GetElasticGlauberGribov ( const G4DynamicParticle dp,
G4int  Z,
G4int  A 
)
inline

Definition at line 177 of file G4GlauberGribovCrossSection.hh.

References GetIsoCrossSection().

Referenced by G4BGGPionElasticXS::BuildPhysicsTable(), G4BGGNucleonElasticXS::BuildPhysicsTable(), G4BGGPionElasticXS::GetElementCrossSection(), and G4BGGNucleonElasticXS::GetElementCrossSection().

179 {
180  GetIsoCrossSection(dp, Z, A);
181  return fElasticXsc;
182 }
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
G4double G4GlauberGribovCrossSection::GetElasticGlauberGribovXsc ( )
inline

Definition at line 104 of file G4GlauberGribovCrossSection.hh.

Referenced by G4NeutronElasticXS::GetElementCrossSection().

104 { return fElasticXsc; };
G4double G4GlauberGribovCrossSection::GetHadronNucleonXsc ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 524 of file G4GlauberGribovCrossSection.cc.

References G4lrint(), G4Element::GetN(), and G4Element::GetZ().

526 {
527  G4int At = G4lrint(anElement->GetN()); // number of nucleons
528  G4int Zt = G4lrint(anElement->GetZ()); // number of protons
529 
530  return GetHadronNucleonXsc(aParticle, At, Zt);
531 }
G4double GetN() const
Definition: G4Element.hh:134
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
G4double G4GlauberGribovCrossSection::GetHadronNucleonXsc ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 541 of file G4GlauberGribovCrossSection.cc.

References CalcMandelstamS(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentum(), python.hepunit::GeV, CLHEP::Hep3Vector::mag(), and python.hepunit::millibarn.

543 {
544  G4double xsection;
545 
546  //G4double targ_mass = G4NucleiProperties::GetNuclearMass(At, Zt);
547 
548  G4double targ_mass = 0.939*GeV; // ~mean neutron and proton ???
549 
550  G4double proj_mass = aParticle->GetMass();
551  G4double proj_momentum = aParticle->GetMomentum().mag();
552  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
553 
554  sMand /= GeV*GeV; // in GeV for parametrisation
555  proj_momentum /= GeV;
556 
557  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
558 
559  G4double aa = At;
560 
561  if(theParticle == theGamma)
562  {
563  xsection = aa*(0.0677*std::pow(sMand,0.0808) + 0.129*std::pow(sMand,-0.4525));
564  }
565  else if(theParticle == theNeutron) // as proton ???
566  {
567  xsection = aa*(21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
568  }
569  else if(theParticle == theProton)
570  {
571  xsection = aa*(21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
572  // xsection = At*( 49.51*std::pow(sMand,-0.097) + 0.314*std::log(sMand)*std::log(sMand) );
573  // xsection = At*( 38.4 + 0.85*std::abs(std::pow(log(sMand),1.47)) );
574  }
575  else if(theParticle == theAProton)
576  {
577  xsection = aa*( 21.70*std::pow(sMand,0.0808) + 98.39*std::pow(sMand,-0.4525));
578  }
579  else if(theParticle == thePiPlus)
580  {
581  xsection = aa*(13.63*std::pow(sMand,0.0808) + 27.56*std::pow(sMand,-0.4525));
582  }
583  else if(theParticle == thePiMinus)
584  {
585  // xsection = At*( 55.2*std::pow(sMand,-0.255) + 0.346*std::log(sMand)*std::log(sMand) );
586  xsection = aa*(13.63*std::pow(sMand,0.0808) + 36.02*std::pow(sMand,-0.4525));
587  }
588  else if(theParticle == theKPlus)
589  {
590  xsection = aa*(11.82*std::pow(sMand,0.0808) + 8.15*std::pow(sMand,-0.4525));
591  }
592  else if(theParticle == theKMinus)
593  {
594  xsection = aa*(11.82*std::pow(sMand,0.0808) + 26.36*std::pow(sMand,-0.4525));
595  }
596  else // as proton ???
597  {
598  xsection = aa*(21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
599  }
600  xsection *= millibarn;
601  return xsection;
602 }
G4ParticleDefinition * GetDefinition() const
int millibarn
Definition: hepunit.py:40
G4double GetMass() const
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
double G4double
Definition: G4Types.hh:76
double mag() const
G4ThreeVector GetMomentum() const
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscNS ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 745 of file G4GlauberGribovCrossSection.cc.

References G4lrint(), G4Element::GetN(), and G4Element::GetZ().

Referenced by GetHNinelasticXsc(), GetIsoCrossSection(), GetRatioQE(), and GetRatioSD().

747 {
748  G4int At = G4lrint(anElement->GetN()); // number of nucleons
749  G4int Zt = G4lrint(anElement->GetZ()); // number of protons
750 
751  return GetHadronNucleonXscNS(aParticle, At, Zt);
752 }
G4double GetN() const
Definition: G4Element.hh:134
G4double GetZ() const
Definition: G4Element.hh:131
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4Element *)
int G4int
Definition: G4Types.hh:78
int G4lrint(double ad)
Definition: templates.hh:163
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscNS ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 763 of file G4GlauberGribovCrossSection.cc.

References CalcMandelstamS(), G4DynamicParticle::GetDefinition(), GetHadronNucleonXscPDG(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentum(), G4ParticleTable::GetParticleTable(), G4DynamicParticle::GetTotalEnergy(), python.hepunit::GeV, CLHEP::Hep3Vector::mag(), python.hepunit::millibarn, G4InuclParticleNames::nn, and G4InuclParticleNames::s0.

765 {
766  G4double xsection(0);
767  // G4double Delta; DHW 19 May 2011: variable set but not used
768  G4double A0, B0;
769  G4double hpXscv(0);
770  G4double hnXscv(0);
771 
772  G4int Nt = At-Zt; // number of neutrons
773  if (Nt < 0) Nt = 0;
774 
775  G4double aa = At;
776  G4double zz = Zt;
777  G4double nn = Nt;
778 
780  GetIonTable()->GetIonMass(Zt, At);
781 
782  targ_mass = 0.939*GeV; // ~mean neutron and proton ???
783 
784  G4double proj_mass = aParticle->GetMass();
785  G4double proj_energy = aParticle->GetTotalEnergy();
786  G4double proj_momentum = aParticle->GetMomentum().mag();
787 
788  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
789 
790  sMand /= GeV*GeV; // in GeV for parametrisation
791  proj_momentum /= GeV;
792  proj_energy /= GeV;
793  proj_mass /= GeV;
794 
795  // General PDG fit constants
796 
797  G4double s0 = 5.38*5.38; // in Gev^2
798  G4double eta1 = 0.458;
799  G4double eta2 = 0.458;
800  G4double B = 0.308;
801 
802 
803  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
804 
805 
806  if(theParticle == theNeutron)
807  {
808  if( proj_momentum >= 373.)
809  {
810  return GetHadronNucleonXscPDG(aParticle,At,Zt);
811  }
812  else if( proj_momentum >= 10.)
813  // if( proj_momentum >= 2.)
814  {
815  // Delta = 1.; // DHW 19 May 2011: variable set but not used
816  // if( proj_energy < 40. ) Delta = 0.916+0.0021*proj_energy;
817 
818  if(proj_momentum >= 10.)
819  {
820  B0 = 7.5;
821  A0 = 100. - B0*std::log(3.0e7);
822 
823  xsection = A0 + B0*std::log(proj_energy) - 11
824  + 103*std::pow(2*0.93827*proj_energy + proj_mass*proj_mass+
825  0.93827*0.93827,-0.165); // mb
826  }
827  xsection *= zz + nn;
828  }
829  else
830  {
831  // nn to be pp
832 
833  if( proj_momentum < 0.73 )
834  {
835  hnXscv = 23 + 50*( std::pow( std::log(0.73/proj_momentum), 3.5 ) );
836  }
837  else if( proj_momentum < 1.05 )
838  {
839  hnXscv = 23 + 40*(std::log(proj_momentum/0.73))*
840  (std::log(proj_momentum/0.73));
841  }
842  else // if( proj_momentum < 10. )
843  {
844  hnXscv = 39.0+
845  75*(proj_momentum - 1.2)/(std::pow(proj_momentum,3.0) + 0.15);
846  }
847  // pn to be np
848 
849  if( proj_momentum < 0.8 )
850  {
851  hpXscv = 33+30*std::pow(std::log(proj_momentum/1.3),4.0);
852  }
853  else if( proj_momentum < 1.4 )
854  {
855  hpXscv = 33+30*std::pow(std::log(proj_momentum/0.95),2.0);
856  }
857  else // if( proj_momentum < 10. )
858  {
859  hpXscv = 33.3+
860  20.8*(std::pow(proj_momentum,2.0)-1.35)/
861  (std::pow(proj_momentum,2.50)+0.95);
862  }
863  xsection = hpXscv*zz + hnXscv*nn;
864  }
865  }
866  else if(theParticle == theProton)
867  {
868  if( proj_momentum >= 373.)
869  {
870  return GetHadronNucleonXscPDG(aParticle,At,Zt);
871  }
872  else if( proj_momentum >= 10.)
873  // if( proj_momentum >= 2.)
874  {
875  // Delta = 1.; DHW 19 May 2011: variable set but not used
876  // if( proj_energy < 40. ) Delta = 0.916+0.0021*proj_energy;
877 
878  if(proj_momentum >= 10.)
879  {
880  B0 = 7.5;
881  A0 = 100. - B0*std::log(3.0e7);
882 
883  xsection = A0 + B0*std::log(proj_energy) - 11
884  + 103*std::pow(2*0.93827*proj_energy + proj_mass*proj_mass+
885  0.93827*0.93827,-0.165); // mb
886  }
887  xsection *= zz + nn;
888  }
889  else
890  {
891  // pp
892 
893  if( proj_momentum < 0.73 )
894  {
895  hpXscv = 23 + 50*( std::pow( std::log(0.73/proj_momentum), 3.5 ) );
896  }
897  else if( proj_momentum < 1.05 )
898  {
899  hpXscv = 23 + 40*(std::log(proj_momentum/0.73))*
900  (std::log(proj_momentum/0.73));
901  }
902  else // if( proj_momentum < 10. )
903  {
904  hpXscv = 39.0+
905  75*(proj_momentum - 1.2)/(std::pow(proj_momentum,3.0) + 0.15);
906  }
907  // pn to be np
908 
909  if( proj_momentum < 0.8 )
910  {
911  hnXscv = 33+30*std::pow(std::log(proj_momentum/1.3),4.0);
912  }
913  else if( proj_momentum < 1.4 )
914  {
915  hnXscv = 33+30*std::pow(std::log(proj_momentum/0.95),2.0);
916  }
917  else // if( proj_momentum < 10. )
918  {
919  hnXscv = 33.3+
920  20.8*(std::pow(proj_momentum,2.0)-1.35)/
921  (std::pow(proj_momentum,2.50)+0.95);
922  }
923  xsection = hpXscv*zz + hnXscv*nn;
924  // xsection = hpXscv*(Zt + Nt);
925  // xsection = hnXscv*(Zt + Nt);
926  }
927  // xsection *= 0.95;
928  }
929  else if( theParticle == theAProton )
930  {
931  // xsection = Zt*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
932  // + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2));
933 
934  // xsection += Nt*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
935  // + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2));
936 
937  G4double logP = std::log(proj_momentum);
938 
939  if( proj_momentum <= 1.0 )
940  {
941  xsection = zz*(65.55 + 53.84/(proj_momentum+1.e-6) );
942  }
943  else
944  {
945  xsection = zz*( 41.1 + 77.2*std::pow( proj_momentum, -0.68)
946  + 0.293*logP*logP - 1.82*logP );
947  }
948  if ( nn > 0.)
949  {
950  xsection += nn*( 41.9 + 96.2*std::pow( proj_momentum, -0.99) - 0.154*logP);
951  }
952  else // H
953  {
954  fInelasticXsc = 38.0 + 38.0*std::pow( proj_momentum, -0.96)
955  - 0.169*logP*logP;
956  fInelasticXsc *= millibarn;
957  }
958  }
959  else if( theParticle == thePiPlus )
960  {
961  if(proj_momentum < 0.4)
962  {
963  G4double Ex3 = 180*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.085/0.085);
964  hpXscv = Ex3+20.0;
965  }
966  else if( proj_momentum < 1.15 )
967  {
968  G4double Ex4 = 88*(std::log(proj_momentum/0.75))*(std::log(proj_momentum/0.75));
969  hpXscv = Ex4+14.0;
970  }
971  else if(proj_momentum < 3.5)
972  {
973  G4double Ex1 = 3.2*std::exp(-(proj_momentum-2.55)*(proj_momentum-2.55)/0.55/0.55);
974  G4double Ex2 = 12*std::exp(-(proj_momentum-1.47)*(proj_momentum-1.47)/0.225/0.225);
975  hpXscv = Ex1+Ex2+27.5;
976  }
977  else // if(proj_momentum > 3.5) // mb
978  {
979  hpXscv = 10.6+2.*std::log(proj_energy)+25*std::pow(proj_energy,-0.43);
980  }
981  // pi+n = pi-p??
982 
983  if(proj_momentum < 0.37)
984  {
985  hnXscv = 28.0 + 40*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.07/0.07);
986  }
987  else if(proj_momentum<0.65)
988  {
989  hnXscv = 26+110*(std::log(proj_momentum/0.48))*(std::log(proj_momentum/0.48));
990  }
991  else if(proj_momentum<1.3)
992  {
993  hnXscv = 36.1+
994  10*std::exp(-(proj_momentum-0.72)*(proj_momentum-0.72)/0.06/0.06)+
995  24*std::exp(-(proj_momentum-1.015)*(proj_momentum-1.015)/0.075/0.075);
996  }
997  else if(proj_momentum<3.0)
998  {
999  hnXscv = 36.1+0.079-4.313*std::log(proj_momentum)+
1000  3*std::exp(-(proj_momentum-2.1)*(proj_momentum-2.1)/0.4/0.4)+
1001  1.5*std::exp(-(proj_momentum-1.4)*(proj_momentum-1.4)/0.12/0.12);
1002  }
1003  else // mb
1004  {
1005  hnXscv = 10.6+2*std::log(proj_energy)+30*std::pow(proj_energy,-0.43);
1006  }
1007  xsection = hpXscv*zz + hnXscv*nn;
1008  }
1009  else if(theParticle == thePiMinus)
1010  {
1011  // pi-n = pi+p??
1012 
1013  if(proj_momentum < 0.4)
1014  {
1015  G4double Ex3 = 180*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.085/0.085);
1016  hnXscv = Ex3+20.0;
1017  }
1018  else if(proj_momentum < 1.15)
1019  {
1020  G4double Ex4 = 88*(std::log(proj_momentum/0.75))*(std::log(proj_momentum/0.75));
1021  hnXscv = Ex4+14.0;
1022  }
1023  else if(proj_momentum < 3.5)
1024  {
1025  G4double Ex1 = 3.2*std::exp(-(proj_momentum-2.55)*(proj_momentum-2.55)/0.55/0.55);
1026  G4double Ex2 = 12*std::exp(-(proj_momentum-1.47)*(proj_momentum-1.47)/0.225/0.225);
1027  hnXscv = Ex1+Ex2+27.5;
1028  }
1029  else // if(proj_momentum > 3.5) // mb
1030  {
1031  hnXscv = 10.6+2.*std::log(proj_energy)+25*std::pow(proj_energy,-0.43);
1032  }
1033  // pi-p
1034 
1035  if(proj_momentum < 0.37)
1036  {
1037  hpXscv = 28.0 + 40*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.07/0.07);
1038  }
1039  else if(proj_momentum<0.65)
1040  {
1041  hpXscv = 26+110*(std::log(proj_momentum/0.48))*(std::log(proj_momentum/0.48));
1042  }
1043  else if(proj_momentum<1.3)
1044  {
1045  hpXscv = 36.1+
1046  10*std::exp(-(proj_momentum-0.72)*(proj_momentum-0.72)/0.06/0.06)+
1047  24*std::exp(-(proj_momentum-1.015)*(proj_momentum-1.015)/0.075/0.075);
1048  }
1049  else if(proj_momentum<3.0)
1050  {
1051  hpXscv = 36.1+0.079-4.313*std::log(proj_momentum)+
1052  3*std::exp(-(proj_momentum-2.1)*(proj_momentum-2.1)/0.4/0.4)+
1053  1.5*std::exp(-(proj_momentum-1.4)*(proj_momentum-1.4)/0.12/0.12);
1054  }
1055  else // mb
1056  {
1057  hpXscv = 10.6+2*std::log(proj_energy)+30*std::pow(proj_energy,-0.43);
1058  }
1059  xsection = hpXscv*zz + hnXscv*nn;
1060  }
1061  else if(theParticle == theKPlus)
1062  {
1063  xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.)
1064  + 7.14*std::pow(sMand,-eta1) - 13.45*std::pow(sMand,-eta2));
1065 
1066  xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.)
1067  + 5.17*std::pow(sMand,-eta1) - 7.23*std::pow(sMand,-eta2));
1068  }
1069  else if(theParticle == theKMinus)
1070  {
1071  xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.)
1072  + 7.14*std::pow(sMand,-eta1) + 13.45*std::pow(sMand,-eta2));
1073 
1074  xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.)
1075  + 5.17*std::pow(sMand,-eta1) + 7.23*std::pow(sMand,-eta2));
1076  }
1077  else if(theParticle == theSMinus)
1078  {
1079  xsection = aa*( 35.20 + B*std::pow(std::log(sMand/s0),2.)
1080  - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2));
1081  }
1082  else if(theParticle == theGamma) // modify later on
1083  {
1084  xsection = aa*( 0.0 + B*std::pow(std::log(sMand/s0),2.)
1085  + 0.032*std::pow(sMand,-eta1) - 0.0*std::pow(sMand,-eta2));
1086 
1087  }
1088  else // as proton ???
1089  {
1090  xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
1091  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2));
1092 
1093  xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
1094  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
1095  }
1096  xsection *= millibarn; // parametrised in mb
1097  return xsection;
1098 }
G4double GetTotalEnergy() const
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
int millibarn
Definition: hepunit.py:40
G4double GetMass() const
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4Element *)
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
static G4ParticleTable * GetParticleTable()
double G4double
Definition: G4Types.hh:76
double mag() const
G4ThreeVector GetMomentum() const
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscPDG ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 611 of file G4GlauberGribovCrossSection.cc.

References G4lrint(), G4Element::GetN(), and G4Element::GetZ().

Referenced by GetHadronNucleonXscNS(), and GetKaonNucleonXscVector().

613 {
614  G4int At = G4lrint(anElement->GetN()); // number of nucleons
615  G4int Zt = G4lrint(anElement->GetZ()); // number of protons
616 
617  return GetHadronNucleonXscPDG(aParticle, At, Zt);
618 }
G4double GetN() const
Definition: G4Element.hh:134
G4double GetZ() const
Definition: G4Element.hh:131
int G4int
Definition: G4Types.hh:78
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4Element *)
int G4lrint(double ad)
Definition: templates.hh:163
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscPDG ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 630 of file G4GlauberGribovCrossSection.cc.

References CalcMandelstamS(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentum(), G4ParticleTable::GetParticleTable(), python.hepunit::GeV, CLHEP::Hep3Vector::mag(), python.hepunit::millibarn, G4InuclParticleNames::nn, and G4InuclParticleNames::s0.

632 {
633  G4double xsection;
634 
635  G4int Nt = At-Zt; // number of neutrons
636  if (Nt < 0) Nt = 0;
637 
638  G4double zz = Zt;
639  G4double aa = At;
640  G4double nn = Nt;
641 
643  GetIonTable()->GetIonMass(Zt, At);
644 
645  targ_mass = 0.939*GeV; // ~mean neutron and proton ???
646 
647  G4double proj_mass = aParticle->GetMass();
648  G4double proj_momentum = aParticle->GetMomentum().mag();
649 
650  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
651 
652  sMand /= GeV*GeV; // in GeV for parametrisation
653 
654  // General PDG fit constants
655 
656  G4double s0 = 5.38*5.38; // in Gev^2
657  G4double eta1 = 0.458;
658  G4double eta2 = 0.458;
659  G4double B = 0.308;
660 
661 
662  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
663 
664 
665  if(theParticle == theNeutron) // proton-neutron fit
666  {
667  xsection = zz*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
668  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
669  xsection += nn*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
670  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2)); // pp for nn
671  }
672  else if(theParticle == theProton)
673  {
674 
675  xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
676  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2));
677 
678  xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
679  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
680  }
681  else if(theParticle == theAProton)
682  {
683  xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
684  + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2));
685 
686  xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
687  + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2));
688  }
689  else if(theParticle == thePiPlus)
690  {
691  xsection = aa*( 20.86 + B*std::pow(std::log(sMand/s0),2.)
692  + 19.24*std::pow(sMand,-eta1) - 6.03*std::pow(sMand,-eta2));
693  }
694  else if(theParticle == thePiMinus)
695  {
696  xsection = aa*( 20.86 + B*std::pow(std::log(sMand/s0),2.)
697  + 19.24*std::pow(sMand,-eta1) + 6.03*std::pow(sMand,-eta2));
698  }
699  else if(theParticle == theKPlus || theParticle == theK0L )
700  {
701  xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.)
702  + 7.14*std::pow(sMand,-eta1) - 13.45*std::pow(sMand,-eta2));
703 
704  xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.)
705  + 5.17*std::pow(sMand,-eta1) - 7.23*std::pow(sMand,-eta2));
706  }
707  else if(theParticle == theKMinus || theParticle == theK0S )
708  {
709  xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.)
710  + 7.14*std::pow(sMand,-eta1) + 13.45*std::pow(sMand,-eta2));
711 
712  xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.)
713  + 5.17*std::pow(sMand,-eta1) + 7.23*std::pow(sMand,-eta2));
714  }
715  else if(theParticle == theSMinus)
716  {
717  xsection = aa*( 35.20 + B*std::pow(std::log(sMand/s0),2.)
718  - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2));
719  }
720  else if(theParticle == theGamma) // modify later on
721  {
722  xsection = aa*( 0.0 + B*std::pow(std::log(sMand/s0),2.)
723  + 0.032*std::pow(sMand,-eta1) - 0.0*std::pow(sMand,-eta2));
724 
725  }
726  else // as proton ???
727  {
728  xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
729  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2));
730 
731  xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
732  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
733  }
734  xsection *= millibarn; // parametrised in mb
735  return xsection;
736 }
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
int millibarn
Definition: hepunit.py:40
G4double GetMass() const
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
static G4ParticleTable * GetParticleTable()
double G4double
Definition: G4Types.hh:76
double mag() const
G4ThreeVector GetMomentum() const
G4double G4GlauberGribovCrossSection::GetHNinelasticXsc ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 1139 of file G4GlauberGribovCrossSection.cc.

References G4lrint(), G4Element::GetN(), and G4Element::GetZ().

Referenced by GetIsoCrossSection(), and GetRatioQE().

1141 {
1142  G4int At = G4lrint(anElement->GetN()); // number of nucleons
1143  G4int Zt = G4lrint(anElement->GetZ()); // number of protons
1144 
1145  return GetHNinelasticXsc(aParticle, At, Zt);
1146 }
G4double GetN() const
Definition: G4Element.hh:134
G4double GetZ() const
Definition: G4Element.hh:131
int G4int
Definition: G4Types.hh:78
int G4lrint(double ad)
Definition: templates.hh:163
G4double GetHNinelasticXsc(const G4DynamicParticle *, const G4Element *)
G4double G4GlauberGribovCrossSection::GetHNinelasticXsc ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 1153 of file G4GlauberGribovCrossSection.cc.

References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscNS(), and GetHNinelasticXscVU().

1155 {
1156  G4ParticleDefinition* hadron = aParticle->GetDefinition();
1157  G4double sumInelastic;
1158  G4int Nt = At - Zt;
1159  if(Nt < 0) Nt = 0;
1160 
1161  if( hadron == theKPlus )
1162  {
1163  sumInelastic = GetHNinelasticXscVU(aParticle, At, Zt);
1164  }
1165  else
1166  {
1167  //sumInelastic = Zt*GetHadronNucleonXscMK(aParticle, theProton);
1168  // sumInelastic += Nt*GetHadronNucleonXscMK(aParticle, theNeutron);
1169  sumInelastic = G4double(Zt)*GetHadronNucleonXscNS(aParticle, 1, 1);
1170  sumInelastic += G4double(Nt)*GetHadronNucleonXscNS(aParticle, 1, 0);
1171  }
1172  return sumInelastic;
1173 }
G4double GetHNinelasticXscVU(const G4DynamicParticle *, G4int At, G4int Zt)
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4Element *)
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
double G4double
Definition: G4Types.hh:76
G4double G4GlauberGribovCrossSection::GetHNinelasticXscVU ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 1181 of file G4GlauberGribovCrossSection.cc.

References G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetMomentum(), G4ParticleDefinition::GetPDGEncoding(), G4DynamicParticle::GetTotalEnergy(), python.hepunit::GeV, CLHEP::Hep3Vector::mag(), and python.hepunit::millibarn.

Referenced by GetHNinelasticXsc().

1183 {
1184  G4int PDGcode = aParticle->GetDefinition()->GetPDGEncoding();
1185  G4int absPDGcode = std::abs(PDGcode);
1186 
1187  G4double Elab = aParticle->GetTotalEnergy();
1188  // (s - 2*0.88*GeV*GeV)/(2*0.939*GeV)/GeV;
1189  G4double Plab = aParticle->GetMomentum().mag();
1190  // std::sqrt(Elab * Elab - 0.88);
1191 
1192  Elab /= GeV;
1193  Plab /= GeV;
1194 
1195  G4double LogPlab = std::log( Plab );
1196  G4double sqrLogPlab = LogPlab * LogPlab;
1197 
1198  //G4cout<<"Plab = "<<Plab<<G4endl;
1199 
1200  G4double NumberOfTargetProtons = G4double(Zt);
1201  G4double NumberOfTargetNucleons = G4double(At);
1202  G4double NumberOfTargetNeutrons = NumberOfTargetNucleons - NumberOfTargetProtons;
1203 
1204  if(NumberOfTargetNeutrons < 0.0) NumberOfTargetNeutrons = 0.0;
1205 
1206  G4double Xtotal, Xelastic, Xinelastic;
1207 
1208  if( absPDGcode > 1000 ) //------Projectile is baryon --------
1209  {
1210  G4double XtotPP = 48.0 + 0. *std::pow(Plab, 0. ) +
1211  0.522*sqrLogPlab - 4.51*LogPlab;
1212 
1213  G4double XtotPN = 47.3 + 0. *std::pow(Plab, 0. ) +
1214  0.513*sqrLogPlab - 4.27*LogPlab;
1215 
1216  G4double XelPP = 11.9 + 26.9*std::pow(Plab,-1.21) +
1217  0.169*sqrLogPlab - 1.85*LogPlab;
1218 
1219  G4double XelPN = 11.9 + 26.9*std::pow(Plab,-1.21) +
1220  0.169*sqrLogPlab - 1.85*LogPlab;
1221 
1222  Xtotal = (NumberOfTargetProtons * XtotPP +
1223  NumberOfTargetNeutrons * XtotPN);
1224 
1225  Xelastic = (NumberOfTargetProtons * XelPP +
1226  NumberOfTargetNeutrons * XelPN);
1227  }
1228  else if( PDGcode == 211 ) //------Projectile is PionPlus -------
1229  {
1230  G4double XtotPiP = 16.4 + 19.3 *std::pow(Plab,-0.42) +
1231  0.19 *sqrLogPlab - 0.0 *LogPlab;
1232 
1233  G4double XtotPiN = 33.0 + 14.0 *std::pow(Plab,-1.36) +
1234  0.456*sqrLogPlab - 4.03*LogPlab;
1235 
1236  G4double XelPiP = 0.0 + 11.4*std::pow(Plab,-0.40) +
1237  0.079*sqrLogPlab - 0.0 *LogPlab;
1238 
1239  G4double XelPiN = 1.76 + 11.2*std::pow(Plab,-0.64) +
1240  0.043*sqrLogPlab - 0.0 *LogPlab;
1241 
1242  Xtotal = ( NumberOfTargetProtons * XtotPiP +
1243  NumberOfTargetNeutrons * XtotPiN );
1244 
1245  Xelastic = ( NumberOfTargetProtons * XelPiP +
1246  NumberOfTargetNeutrons * XelPiN );
1247  }
1248  else if( PDGcode == -211 ) //------Projectile is PionMinus -------
1249  {
1250  G4double XtotPiP = 33.0 + 14.0 *std::pow(Plab,-1.36) +
1251  0.456*sqrLogPlab - 4.03*LogPlab;
1252 
1253  G4double XtotPiN = 16.4 + 19.3 *std::pow(Plab,-0.42) +
1254  0.19 *sqrLogPlab - 0.0 *LogPlab;
1255 
1256  G4double XelPiP = 1.76 + 11.2*std::pow(Plab,-0.64) +
1257  0.043*sqrLogPlab - 0.0 *LogPlab;
1258 
1259  G4double XelPiN = 0.0 + 11.4*std::pow(Plab,-0.40) +
1260  0.079*sqrLogPlab - 0.0 *LogPlab;
1261 
1262  Xtotal = ( NumberOfTargetProtons * XtotPiP +
1263  NumberOfTargetNeutrons * XtotPiN );
1264 
1265  Xelastic = ( NumberOfTargetProtons * XelPiP +
1266  NumberOfTargetNeutrons * XelPiN );
1267  }
1268  else if( PDGcode == 111 ) //------Projectile is PionZero -------
1269  {
1270  G4double XtotPiP =(16.4 + 19.3 *std::pow(Plab,-0.42) +
1271  0.19 *sqrLogPlab - 0.0 *LogPlab + //Pi+
1272  33.0 + 14.0 *std::pow(Plab,-1.36) +
1273  0.456*sqrLogPlab - 4.03*LogPlab)/2; //Pi-
1274 
1275  G4double XtotPiN =(33.0 + 14.0 *std::pow(Plab,-1.36) +
1276  0.456*sqrLogPlab - 4.03*LogPlab + //Pi+
1277  16.4 + 19.3 *std::pow(Plab,-0.42) +
1278  0.19 *sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1279 
1280  G4double XelPiP =( 0.0 + 11.4*std::pow(Plab,-0.40) +
1281  0.079*sqrLogPlab - 0.0 *LogPlab + //Pi+
1282  1.76 + 11.2*std::pow(Plab,-0.64) +
1283  0.043*sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1284 
1285  G4double XelPiN =( 1.76 + 11.2*std::pow(Plab,-0.64) +
1286  0.043*sqrLogPlab - 0.0 *LogPlab + //Pi+
1287  0.0 + 11.4*std::pow(Plab,-0.40) +
1288  0.079*sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1289 
1290  Xtotal = ( NumberOfTargetProtons * XtotPiP +
1291  NumberOfTargetNeutrons * XtotPiN );
1292 
1293  Xelastic = ( NumberOfTargetProtons * XelPiP +
1294  NumberOfTargetNeutrons * XelPiN );
1295  }
1296  else if( PDGcode == 321 ) //------Projectile is KaonPlus -------
1297  {
1298  G4double XtotKP = 18.1 + 0. *std::pow(Plab, 0. ) +
1299  0.26 *sqrLogPlab - 1.0 *LogPlab;
1300  G4double XtotKN = 18.7 + 0. *std::pow(Plab, 0. ) +
1301  0.21 *sqrLogPlab - 0.89*LogPlab;
1302 
1303  G4double XelKP = 5.0 + 8.1*std::pow(Plab,-1.8 ) +
1304  0.16 *sqrLogPlab - 1.3 *LogPlab;
1305 
1306  G4double XelKN = 7.3 + 0. *std::pow(Plab,-0. ) +
1307  0.29 *sqrLogPlab - 2.4 *LogPlab;
1308 
1309  Xtotal = ( NumberOfTargetProtons * XtotKP +
1310  NumberOfTargetNeutrons * XtotKN );
1311 
1312  Xelastic = ( NumberOfTargetProtons * XelKP +
1313  NumberOfTargetNeutrons * XelKN );
1314  }
1315  else if( PDGcode ==-321 ) //------Projectile is KaonMinus ------
1316  {
1317  G4double XtotKP = 32.1 + 0. *std::pow(Plab, 0. ) +
1318  0.66 *sqrLogPlab - 5.6 *LogPlab;
1319  G4double XtotKN = 25.2 + 0. *std::pow(Plab, 0. ) +
1320  0.38 *sqrLogPlab - 2.9 *LogPlab;
1321 
1322  G4double XelKP = 7.3 + 0. *std::pow(Plab,-0. ) +
1323  0.29 *sqrLogPlab - 2.4 *LogPlab;
1324 
1325  G4double XelKN = 5.0 + 8.1*std::pow(Plab,-1.8 ) +
1326  0.16 *sqrLogPlab - 1.3 *LogPlab;
1327 
1328  Xtotal = ( NumberOfTargetProtons * XtotKP +
1329  NumberOfTargetNeutrons * XtotKN );
1330 
1331  Xelastic = ( NumberOfTargetProtons * XelKP +
1332  NumberOfTargetNeutrons * XelKN );
1333  }
1334  else if( PDGcode == 311 ) //------Projectile is KaonZero ------
1335  {
1336  G4double XtotKP = ( 18.1 + 0. *std::pow(Plab, 0. ) +
1337  0.26 *sqrLogPlab - 1.0 *LogPlab + //K+
1338  32.1 + 0. *std::pow(Plab, 0. ) +
1339  0.66 *sqrLogPlab - 5.6 *LogPlab)/2; //K-
1340 
1341  G4double XtotKN = ( 18.7 + 0. *std::pow(Plab, 0. ) +
1342  0.21 *sqrLogPlab - 0.89*LogPlab + //K+
1343  25.2 + 0. *std::pow(Plab, 0. ) +
1344  0.38 *sqrLogPlab - 2.9 *LogPlab)/2; //K-
1345 
1346  G4double XelKP = ( 5.0 + 8.1*std::pow(Plab,-1.8 )
1347  + 0.16 *sqrLogPlab - 1.3 *LogPlab + //K+
1348  7.3 + 0. *std::pow(Plab,-0. ) +
1349  0.29 *sqrLogPlab - 2.4 *LogPlab)/2; //K-
1350 
1351  G4double XelKN = ( 7.3 + 0. *std::pow(Plab,-0. ) +
1352  0.29 *sqrLogPlab - 2.4 *LogPlab + //K+
1353  5.0 + 8.1*std::pow(Plab,-1.8 ) +
1354  0.16 *sqrLogPlab - 1.3 *LogPlab)/2; //K-
1355 
1356  Xtotal = ( NumberOfTargetProtons * XtotKP +
1357  NumberOfTargetNeutrons * XtotKN );
1358 
1359  Xelastic = ( NumberOfTargetProtons * XelKP +
1360  NumberOfTargetNeutrons * XelKN );
1361  }
1362  else //------Projectile is undefined, Nucleon assumed
1363  {
1364  G4double XtotPP = 48.0 + 0. *std::pow(Plab, 0. ) +
1365  0.522*sqrLogPlab - 4.51*LogPlab;
1366 
1367  G4double XtotPN = 47.3 + 0. *std::pow(Plab, 0. ) +
1368  0.513*sqrLogPlab - 4.27*LogPlab;
1369 
1370  G4double XelPP = 11.9 + 26.9*std::pow(Plab,-1.21) +
1371  0.169*sqrLogPlab - 1.85*LogPlab;
1372  G4double XelPN = 11.9 + 26.9*std::pow(Plab,-1.21) +
1373  0.169*sqrLogPlab - 1.85*LogPlab;
1374 
1375  Xtotal = ( NumberOfTargetProtons * XtotPP +
1376  NumberOfTargetNeutrons * XtotPN );
1377 
1378  Xelastic = ( NumberOfTargetProtons * XelPP +
1379  NumberOfTargetNeutrons * XelPN );
1380  }
1381  Xinelastic = Xtotal - Xelastic;
1382 
1383  if( Xinelastic < 0.) Xinelastic = 0.;
1384 
1385  return Xinelastic*= millibarn;
1386 }
G4double GetTotalEnergy() const
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
int millibarn
Definition: hepunit.py:40
double G4double
Definition: G4Types.hh:76
double mag() const
G4ThreeVector GetMomentum() const
G4double G4GlauberGribovCrossSection::GetInelasticGlauberGribov ( const G4DynamicParticle dp,
G4int  Z,
G4int  A 
)
inline

Definition at line 188 of file G4GlauberGribovCrossSection.hh.

References GetIsoCrossSection().

Referenced by G4CrossSectionPairGG::BuildPhysicsTable(), G4BGGNucleonInelasticXS::BuildPhysicsTable(), G4BGGPionInelasticXS::BuildPhysicsTable(), G4CrossSectionPairGG::GetElementCrossSection(), G4BGGNucleonInelasticXS::GetElementCrossSection(), and G4BGGPionInelasticXS::GetElementCrossSection().

190 {
191  GetIsoCrossSection(dp, Z, A);
192  return fInelasticXsc;
193 }
virtual G4double GetIsoCrossSection(const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
G4double G4GlauberGribovCrossSection::GetInelasticGlauberGribovXsc ( )
inline

Definition at line 105 of file G4GlauberGribovCrossSection.hh.

Referenced by G4NeutronInelasticXS::GetElementCrossSection().

105 { return fInelasticXsc; };
G4double G4GlauberGribovCrossSection::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 321 of file G4GlauberGribovCrossSection.cc.

References python.hepunit::fermi, G4DynamicParticle::GetDefinition(), G4HadronNucleonXsc::GetHadronNucleonXscNS(), GetHadronNucleonXscNS(), GetHNinelasticXsc(), G4HadronNucleonXsc::GetInelasticHadronNucleonXsc(), GetKaonNucleonXscVector(), GetNucleusRadius(), GetParticleBarCorIn(), GetParticleBarCorTot(), N, and python.hepunit::pi.

Referenced by GetElasticGlauberGribov(), G4NeutronElasticXS::GetElementCrossSection(), G4NeutronInelasticXS::GetElementCrossSection(), and GetInelasticGlauberGribov().

326 {
327  G4double xsection, sigma, cofInelastic, cofTotal, nucleusSquare, ratio;
328  G4double hpInXsc(0.), hnInXsc(0.);
329  G4double R = GetNucleusRadius(A);
330 
331  G4int N = A - Z; // number of neutrons
332  if (N < 0) N = 0;
333 
334  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
335 
336  if( theParticle == theProton ||
337  theParticle == theNeutron ||
338  theParticle == thePiPlus ||
339  theParticle == thePiMinus )
340  {
341  // sigma = GetHadronNucleonXscNS(aParticle, A, Z);
342 
343  sigma = Z*hnXsc->GetHadronNucleonXscNS(aParticle, theProton);
344 
345  hpInXsc = hnXsc->GetInelasticHadronNucleonXsc();
346 
347  sigma += N*hnXsc->GetHadronNucleonXscNS(aParticle, theNeutron);
348 
349  hnInXsc = hnXsc->GetInelasticHadronNucleonXsc();
350 
351  cofInelastic = 2.4;
352  cofTotal = 2.0;
353  }
354  else if( theParticle == theKPlus ||
355  theParticle == theKMinus ||
356  theParticle == theK0S ||
357  theParticle == theK0L )
358  {
359  sigma = GetKaonNucleonXscVector(aParticle, A, Z);
360  cofInelastic = 2.2;
361  cofTotal = 2.0;
362  R = 1.3*fermi;
363  R *= std::pow(G4double(A), 0.3333);
364  }
365  else
366  {
367  sigma = GetHadronNucleonXscNS(aParticle, A, Z);
368  cofInelastic = 2.2;
369  cofTotal = 2.0;
370  }
371  // cofInelastic = 2.0;
372 
373  if( A > 1 )
374  {
375  nucleusSquare = cofTotal*pi*R*R; // basically 2piRR
376  ratio = sigma/nucleusSquare;
377 
378  xsection = nucleusSquare*std::log( 1. + ratio );
379 
380  xsection *= GetParticleBarCorTot(theParticle, Z);
381 
382  fTotalXsc = xsection;
383 
384 
385 
386  fInelasticXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic;
387 
388  fInelasticXsc *= GetParticleBarCorIn(theParticle, Z);
389 
390  fElasticXsc = fTotalXsc - fInelasticXsc;
391 
392  if(fElasticXsc < 0.) fElasticXsc = 0.;
393 
394  G4double difratio = ratio/(1.+ratio);
395 
396  fDiffractionXsc = 0.5*nucleusSquare*( difratio - std::log( 1. + difratio ) );
397 
398 
399  // sigma = GetHNinelasticXsc(aParticle, A, Z);
400 
401  sigma = Z*hpInXsc + N*hnInXsc;
402 
403  ratio = sigma/nucleusSquare;
404 
405  fProductionXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic;
406 
407  if (fElasticXsc < 0.) fElasticXsc = 0.;
408  }
409  else // H
410  {
411  fTotalXsc = sigma;
412  xsection = sigma;
413 
414  if ( theParticle != theAProton )
415  {
416  sigma = GetHNinelasticXsc(aParticle, A, Z);
417  fInelasticXsc = sigma;
418  fElasticXsc = fTotalXsc - fInelasticXsc;
419  }
420  else
421  {
422  fElasticXsc = fTotalXsc - fInelasticXsc;
423  }
424  if (fElasticXsc < 0.) fElasticXsc = 0.;
425 
426  }
427  return xsection;
428 }
G4double GetParticleBarCorTot(const G4ParticleDefinition *theParticle, G4int Z)
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4Element *)
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4ParticleDefinition *)
G4double GetKaonNucleonXscVector(const G4DynamicParticle *, G4int At, G4int Zt)
G4double GetParticleBarCorIn(const G4ParticleDefinition *theParticle, G4int Z)
G4double GetNucleusRadius(const G4DynamicParticle *, const G4Element *)
**D E S C R I P T I O N
double G4double
Definition: G4Types.hh:76
G4double GetHNinelasticXsc(const G4DynamicParticle *, const G4Element *)
G4double GetInelasticHadronNucleonXsc()
G4double G4GlauberGribovCrossSection::GetKaonNucleonXscVector ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 1101 of file G4GlauberGribovCrossSection.cc.

References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscPDG(), G4DynamicParticle::GetKineticEnergy(), G4HadronNucleonXsc::GetKmNeutronTotXscVector(), G4HadronNucleonXsc::GetKmProtonTotXscVector(), G4HadronNucleonXsc::GetKpNeutronTotXscVector(), G4HadronNucleonXsc::GetKpProtonTotXscVector(), and python.hepunit::GeV.

Referenced by GetIsoCrossSection().

1103 {
1104  G4double Tkin, logTkin, xsc, xscP, xscN;
1105  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
1106 
1107  G4int Nt = At-Zt; // number of neutrons
1108  if (Nt < 0) Nt = 0;
1109 
1110  Tkin = aParticle->GetKineticEnergy(); // Tkin in MeV
1111 
1112  if( Tkin > 70*GeV ) return GetHadronNucleonXscPDG(aParticle,At,Zt);
1113 
1114  logTkin = std::log(Tkin); // Tkin in MeV!!!
1115 
1116  if( theParticle == theKPlus )
1117  {
1118  xscP = hnXsc->GetKpProtonTotXscVector(logTkin);
1119  xscN = hnXsc->GetKpNeutronTotXscVector(logTkin);
1120  }
1121  else if( theParticle == theKMinus )
1122  {
1123  xscP = hnXsc->GetKmProtonTotXscVector(logTkin);
1124  xscN = hnXsc->GetKmNeutronTotXscVector(logTkin);
1125  }
1126  else // K-zero as half of K+ and K-
1127  {
1128  xscP = (hnXsc->GetKpProtonTotXscVector(logTkin)+hnXsc->GetKmProtonTotXscVector(logTkin))*0.5;
1129  xscN = (hnXsc->GetKpNeutronTotXscVector(logTkin)+hnXsc->GetKmNeutronTotXscVector(logTkin))*0.5;
1130  }
1131  xsc = xscP*Zt + xscN*Nt;
1132  return xsc;
1133 }
G4double GetKmNeutronTotXscVector(G4double logEnergy)
G4double GetKpProtonTotXscVector(G4double logEnergy)
G4double GetKineticEnergy() const
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4Element *)
G4double GetKmProtonTotXscVector(G4double logEnergy)
double G4double
Definition: G4Types.hh:76
G4double GetKpNeutronTotXscVector(G4double logEnergy)
G4double G4GlauberGribovCrossSection::GetNucleusRadius ( const G4DynamicParticle ,
const G4Element anElement 
)

Definition at line 1393 of file G4GlauberGribovCrossSection.cc.

References G4lrint(), G4Element::GetN(), and G4INCL::Math::oneThird.

Referenced by GetIsoCrossSection(), GetRatioQE(), and GetRatioSD().

1395 {
1396  G4int At = G4lrint(anElement->GetN());
1397  G4double oneThird = 1.0/3.0;
1398  G4double cubicrAt = std::pow(G4double(At), oneThird);
1399 
1400  G4double R; // = fRadiusConst*cubicrAt;
1401  /*
1402  G4double tmp = std::pow( cubicrAt-1., 3.);
1403  tmp += At;
1404  tmp *= 0.5;
1405 
1406  if (At > 20.) // 20.
1407  {
1408  R = fRadiusConst*std::pow (tmp, oneThird);
1409  }
1410  else
1411  {
1412  R = fRadiusConst*cubicrAt;
1413  }
1414  */
1415 
1416  R = fRadiusConst*cubicrAt;
1417 
1418  G4double meanA = 21.;
1419 
1420  G4double tauA1 = 40.;
1421  G4double tauA2 = 10.;
1422  G4double tauA3 = 5.;
1423 
1424  G4double a1 = 0.85;
1425  G4double b1 = 1. - a1;
1426 
1427  G4double b2 = 0.3;
1428  G4double b3 = 4.;
1429 
1430  if (At > 20) // 20.
1431  {
1432  R *= ( a1 + b1*std::exp( -(At - meanA)/tauA1) );
1433  }
1434  else if (At > 3)
1435  {
1436  R *= ( 1.0 + b2*( 1. - std::exp( (At - meanA)/tauA2) ) );
1437  }
1438  else
1439  {
1440  R *= ( 1.0 + b3*( 1. - std::exp( (At - meanA)/tauA3) ) );
1441  }
1442  return R;
1443 
1444 }
G4double GetN() const
Definition: G4Element.hh:134
int G4int
Definition: G4Types.hh:78
int G4lrint(double ad)
Definition: templates.hh:163
double G4double
Definition: G4Types.hh:76
const G4double oneThird
G4double G4GlauberGribovCrossSection::GetNucleusRadius ( G4int  At)

Definition at line 1450 of file G4GlauberGribovCrossSection.cc.

References G4INCL::Math::oneThird.

1451 {
1452  G4double oneThird = 1.0/3.0;
1453  G4double cubicrAt = std::pow(G4double(At), oneThird);
1454 
1455  G4double R; // = fRadiusConst*cubicrAt;
1456 
1457  /*
1458  G4double tmp = std::pow( cubicrAt-1., 3.);
1459  tmp += At;
1460  tmp *= 0.5;
1461 
1462  if (At > 20.)
1463  {
1464  R = fRadiusConst*std::pow (tmp, oneThird);
1465  }
1466  else
1467  {
1468  R = fRadiusConst*cubicrAt;
1469  }
1470  */
1471 
1472  R = fRadiusConst*cubicrAt;
1473 
1474  G4double meanA = 20.;
1475  G4double tauA = 20.;
1476 
1477  if (At > 20) // 20.
1478  {
1479  R *= ( 0.8 + 0.2*std::exp( -(G4double(At) - meanA)/tauA) );
1480  }
1481  else
1482  {
1483  R *= ( 1.0 + 0.1*( 1. - std::exp( (G4double(At) - meanA)/tauA) ) );
1484  }
1485 
1486  return R;
1487 }
double G4double
Definition: G4Types.hh:76
const G4double oneThird
G4double G4GlauberGribovCrossSection::GetParticleBarCorIn ( const G4ParticleDefinition theParticle,
G4int  Z 
)
inline

Definition at line 221 of file G4GlauberGribovCrossSection.hh.

Referenced by GetIsoCrossSection().

223 {
224  if(Z >= 2 && Z <= 92)
225  {
226  if( theParticle == theProton ) return fProtonBarCorrectionIn[Z];
227  else if( theParticle == theNeutron) return fNeutronBarCorrectionIn[Z];
228  else if( theParticle == thePiPlus ) return fPionPlusBarCorrectionIn[Z];
229  else if( theParticle == thePiMinus) return fPionMinusBarCorrectionIn[Z];
230  else return 1.0;
231  }
232  else return 1.0;
233 }
G4double G4GlauberGribovCrossSection::GetParticleBarCorTot ( const G4ParticleDefinition theParticle,
G4int  Z 
)
inline

Definition at line 201 of file G4GlauberGribovCrossSection.hh.

Referenced by GetIsoCrossSection().

203 {
204  if(Z >= 2 && Z <= 92)
205  {
206  if( theParticle == theProton ) return fProtonBarCorrectionTot[Z];
207  else if( theParticle == theNeutron) return fNeutronBarCorrectionTot[Z];
208  else if( theParticle == thePiPlus ) return fPionPlusBarCorrectionTot[Z];
209  else if( theParticle == thePiMinus) return fPionMinusBarCorrectionTot[Z];
210  else return 1.0;
211  }
212  else return 1.0;
213 }
G4double G4GlauberGribovCrossSection::GetProductionGlauberGribovXsc ( )
inline

Definition at line 106 of file G4GlauberGribovCrossSection.hh.

106 { return fProductionXsc; };
G4double G4GlauberGribovCrossSection::GetRadiusConst ( )
inline

Definition at line 108 of file G4GlauberGribovCrossSection.hh.

108 { return fRadiusConst; };
G4double G4GlauberGribovCrossSection::GetRatioQE ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 477 of file G4GlauberGribovCrossSection.cc.

References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscNS(), GetHNinelasticXsc(), GetNucleusRadius(), and python.hepunit::pi.

478 {
479  G4double sigma, cofInelastic, cofTotal, nucleusSquare, ratio;
480  G4double R = GetNucleusRadius(A);
481 
482  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
483 
484  if( theParticle == theProton ||
485  theParticle == theNeutron ||
486  theParticle == thePiPlus ||
487  theParticle == thePiMinus )
488  {
489  sigma = GetHadronNucleonXscNS(aParticle, A, Z);
490  cofInelastic = 2.4;
491  cofTotal = 2.0;
492  }
493  else
494  {
495  sigma = GetHadronNucleonXscNS(aParticle, A, Z);
496  cofInelastic = 2.2;
497  cofTotal = 2.0;
498  }
499  nucleusSquare = cofTotal*pi*R*R; // basically 2piRR
500  ratio = sigma/nucleusSquare;
501 
502  fInelasticXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic;
503 
504  sigma = GetHNinelasticXsc(aParticle, A, Z);
505  ratio = sigma/nucleusSquare;
506 
507  fProductionXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic;
508 
509  if (fInelasticXsc > fProductionXsc) ratio = (fInelasticXsc-fProductionXsc)/fInelasticXsc;
510  else ratio = 0.;
511  if ( ratio < 0. ) ratio = 0.;
512 
513  return ratio;
514 }
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4Element *)
G4ParticleDefinition * GetDefinition() const
G4double GetNucleusRadius(const G4DynamicParticle *, const G4Element *)
double G4double
Definition: G4Types.hh:76
G4double GetHNinelasticXsc(const G4DynamicParticle *, const G4Element *)
G4double G4GlauberGribovCrossSection::GetRatioSD ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 435 of file G4GlauberGribovCrossSection.cc.

References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscNS(), GetNucleusRadius(), and python.hepunit::pi.

436 {
437  G4double sigma, cofInelastic, cofTotal, nucleusSquare, ratio;
438  G4double R = GetNucleusRadius(A);
439 
440  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
441 
442  if( theParticle == theProton ||
443  theParticle == theNeutron ||
444  theParticle == thePiPlus ||
445  theParticle == thePiMinus )
446  {
447  sigma = GetHadronNucleonXscNS(aParticle, A, Z);
448  cofInelastic = 2.4;
449  cofTotal = 2.0;
450  }
451  else
452  {
453  sigma = GetHadronNucleonXscNS(aParticle, A, Z);
454  cofInelastic = 2.2;
455  cofTotal = 2.0;
456  }
457  nucleusSquare = cofTotal*pi*R*R; // basically 2piRR
458  ratio = sigma/nucleusSquare;
459 
460  fInelasticXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic;
461 
462  G4double difratio = ratio/(1.+ratio);
463 
464  fDiffractionXsc = 0.5*nucleusSquare*( difratio - std::log( 1. + difratio ) );
465 
466  if (fInelasticXsc > 0.) ratio = fDiffractionXsc/fInelasticXsc;
467  else ratio = 0.;
468 
469  return ratio;
470 }
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4Element *)
G4ParticleDefinition * GetDefinition() const
G4double GetNucleusRadius(const G4DynamicParticle *, const G4Element *)
double G4double
Definition: G4Types.hh:76
G4double G4GlauberGribovCrossSection::GetTotalGlauberGribovXsc ( )
inline

Definition at line 103 of file G4GlauberGribovCrossSection.hh.

103 { return fTotalXsc; };
G4bool G4GlauberGribovCrossSection::IsIsoApplicable ( const G4DynamicParticle aDP,
G4int  Z,
G4int  A,
const G4Element elm = 0,
const G4Material mat = 0 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 285 of file G4GlauberGribovCrossSection.cc.

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

289 {
290  G4bool applicable = false;
291  // G4int baryonNumber = aDP->GetDefinition()->GetBaryonNumber();
292  G4double kineticEnergy = aDP->GetKineticEnergy();
293 
294  const G4ParticleDefinition* theParticle = aDP->GetDefinition();
295 
296  if ( ( kineticEnergy >= fLowerLimit &&
297  Z > 1 && // >= He
298  ( theParticle == theAProton ||
299  theParticle == theGamma ||
300  theParticle == theKPlus ||
301  theParticle == theKMinus ||
302  theParticle == theK0L ||
303  theParticle == theK0S ||
304  theParticle == theSMinus ||
305  theParticle == theProton ||
306  theParticle == theNeutron ||
307  theParticle == thePiPlus ||
308  theParticle == thePiMinus ) ) ) applicable = true;
309 
310  return applicable;
311 }
G4double GetKineticEnergy() const
G4ParticleDefinition * GetDefinition() const
bool G4bool
Definition: G4Types.hh:79
double G4double
Definition: G4Types.hh:76
void G4GlauberGribovCrossSection::SetEnergyLowerLimit ( G4double  E)
inline

Definition at line 113 of file G4GlauberGribovCrossSection.hh.

113 {fLowerLimit=E;};

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