Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4GlauberGribovCrossSection Class Reference

#include <G4GlauberGribovCrossSection.hh>

+ Inheritance diagram for G4GlauberGribovCrossSection:

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)
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=0, const G4Material *mat=0)
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4double GetIsoCrossSection (const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
 
virtual G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
virtual void CrossSectionDescription (std::ostream &) const
 
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::G4GlauberGribovCrossSection ( )

Definition at line 230 of file G4GlauberGribovCrossSection.cc.

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

◆ ~G4GlauberGribovCrossSection()

G4GlauberGribovCrossSection::~G4GlauberGribovCrossSection ( )
virtual

Definition at line 275 of file G4GlauberGribovCrossSection.cc.

276{
277 if (hnXsc) delete hnXsc;
278}

Member Function Documentation

◆ CalcMandelstamS()

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

Definition at line 1507 of file G4GlauberGribovCrossSection.cc.

1510{
1511 G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1512 G4double sMand = mp*mp + mt*mt + 2*Elab*mt ;
1513
1514 return sMand;
1515}
double G4double
Definition: G4Types.hh:64

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

◆ CalculateEcmValue()

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

Definition at line 1491 of file G4GlauberGribovCrossSection.cc.

1494{
1495 G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1496 G4double Ecm = std::sqrt ( mp * mp + mt * mt + 2 * Elab * mt );
1497 // G4double Pcm = Plab * mt / Ecm;
1498 // G4double KEcm = std::sqrt ( Pcm * Pcm + mp * mp ) - mp;
1499
1500 return Ecm ; // KEcm;
1501}

◆ CrossSectionDescription()

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

Reimplemented from G4VCrossSectionDataSet.

Definition at line 1521 of file G4GlauberGribovCrossSection.cc.

1522{
1523 outFile << "G4GlauberGribovCrossSection calculates total, inelastic and\n"
1524 << "elastic cross sections for hadron-nucleus cross sections using\n"
1525 << "the Glauber model with Gribov corrections. It is valid for all\n"
1526 << "targets except hydrogen, and for incident p, pbar, n, sigma-,\n"
1527 << "pi+, pi-, K+, K- and gammas with energies above 3 GeV. This is\n"
1528 << "a cross section component which is to be used to build a cross\n"
1529 << "data set.\n";
1530}

◆ Default_Name()

static const char * G4GlauberGribovCrossSection::Default_Name ( )
inlinestatic

Definition at line 62 of file G4GlauberGribovCrossSection.hh.

62{return "Glauber-Gribov";}

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

◆ GetDiffractionGlauberGribovXsc()

G4double G4GlauberGribovCrossSection::GetDiffractionGlauberGribovXsc ( )
inline

Definition at line 107 of file G4GlauberGribovCrossSection.hh.

107{ return fDiffractionXsc; };

◆ GetElasticGlauberGribov()

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

Definition at line 177 of file G4GlauberGribovCrossSection.hh.

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)

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

◆ GetElasticGlauberGribovXsc()

G4double G4GlauberGribovCrossSection::GetElasticGlauberGribovXsc ( )
inline

Definition at line 104 of file G4GlauberGribovCrossSection.hh.

104{ return fElasticXsc; };

Referenced by G4NeutronElasticXS::GetElementCrossSection().

◆ GetHadronNucleonXsc() [1/2]

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

Definition at line 522 of file G4GlauberGribovCrossSection.cc.

524{
525 G4int At = G4lrint(anElement->GetN()); // number of nucleons
526 G4int Zt = G4lrint(anElement->GetZ()); // number of protons
527
528 return GetHadronNucleonXsc(aParticle, At, Zt);
529}
int G4int
Definition: G4Types.hh:66
G4double GetZ() const
Definition: G4Element.hh:131
G4double GetN() const
Definition: G4Element.hh:134
G4double GetHadronNucleonXsc(const G4DynamicParticle *, const G4Element *)
int G4lrint(double ad)
Definition: templates.hh:163

Referenced by GetHadronNucleonXsc().

◆ GetHadronNucleonXsc() [2/2]

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

Definition at line 539 of file G4GlauberGribovCrossSection.cc.

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

◆ GetHadronNucleonXscNS() [1/2]

G4double G4GlauberGribovCrossSection::GetHadronNucleonXscNS ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 743 of file G4GlauberGribovCrossSection.cc.

745{
746 G4int At = G4lrint(anElement->GetN()); // number of nucleons
747 G4int Zt = G4lrint(anElement->GetZ()); // number of protons
748
749 return GetHadronNucleonXscNS(aParticle, At, Zt);
750}
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4Element *)

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

◆ GetHadronNucleonXscNS() [2/2]

G4double G4GlauberGribovCrossSection::GetHadronNucleonXscNS ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 761 of file G4GlauberGribovCrossSection.cc.

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

◆ GetHadronNucleonXscPDG() [1/2]

G4double G4GlauberGribovCrossSection::GetHadronNucleonXscPDG ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 609 of file G4GlauberGribovCrossSection.cc.

611{
612 G4int At = G4lrint(anElement->GetN()); // number of nucleons
613 G4int Zt = G4lrint(anElement->GetZ()); // number of protons
614
615 return GetHadronNucleonXscPDG(aParticle, At, Zt);
616}

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

◆ GetHadronNucleonXscPDG() [2/2]

G4double G4GlauberGribovCrossSection::GetHadronNucleonXscPDG ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 628 of file G4GlauberGribovCrossSection.cc.

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

◆ GetHNinelasticXsc() [1/2]

G4double G4GlauberGribovCrossSection::GetHNinelasticXsc ( const G4DynamicParticle aParticle,
const G4Element anElement 
)

Definition at line 1137 of file G4GlauberGribovCrossSection.cc.

1139{
1140 G4int At = G4lrint(anElement->GetN()); // number of nucleons
1141 G4int Zt = G4lrint(anElement->GetZ()); // number of protons
1142
1143 return GetHNinelasticXsc(aParticle, At, Zt);
1144}
G4double GetHNinelasticXsc(const G4DynamicParticle *, const G4Element *)

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

◆ GetHNinelasticXsc() [2/2]

G4double G4GlauberGribovCrossSection::GetHNinelasticXsc ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 1151 of file G4GlauberGribovCrossSection.cc.

1153{
1154 G4ParticleDefinition* hadron = aParticle->GetDefinition();
1155 G4double sumInelastic;
1156 G4int Nt = At - Zt;
1157 if(Nt < 0) Nt = 0;
1158
1159 if( hadron == theKPlus )
1160 {
1161 sumInelastic = GetHNinelasticXscVU(aParticle, At, Zt);
1162 }
1163 else
1164 {
1165 //sumInelastic = Zt*GetHadronNucleonXscMK(aParticle, theProton);
1166 // sumInelastic += Nt*GetHadronNucleonXscMK(aParticle, theNeutron);
1167 sumInelastic = G4double(Zt)*GetHadronNucleonXscNS(aParticle, 1, 1);
1168 sumInelastic += G4double(Nt)*GetHadronNucleonXscNS(aParticle, 1, 0);
1169 }
1170 return sumInelastic;
1171}
G4double GetHNinelasticXscVU(const G4DynamicParticle *, G4int At, G4int Zt)

◆ GetHNinelasticXscVU()

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

Definition at line 1179 of file G4GlauberGribovCrossSection.cc.

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

Referenced by GetHNinelasticXsc().

◆ GetInelasticGlauberGribov()

◆ GetInelasticGlauberGribovXsc()

G4double G4GlauberGribovCrossSection::GetInelasticGlauberGribovXsc ( )
inline

Definition at line 105 of file G4GlauberGribovCrossSection.hh.

105{ return fInelasticXsc; };

Referenced by G4NeutronInelasticXS::GetElementCrossSection().

◆ GetIsoCrossSection()

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 319 of file G4GlauberGribovCrossSection.cc.

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

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

◆ GetKaonNucleonXscVector()

G4double G4GlauberGribovCrossSection::GetKaonNucleonXscVector ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 1099 of file G4GlauberGribovCrossSection.cc.

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

Referenced by GetIsoCrossSection().

◆ GetNucleusRadius() [1/2]

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

Definition at line 1391 of file G4GlauberGribovCrossSection.cc.

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

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

◆ GetNucleusRadius() [2/2]

G4double G4GlauberGribovCrossSection::GetNucleusRadius ( G4int  At)

Definition at line 1448 of file G4GlauberGribovCrossSection.cc.

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

◆ GetParticleBarCorIn()

G4double G4GlauberGribovCrossSection::GetParticleBarCorIn ( const G4ParticleDefinition theParticle,
G4int  Z 
)
inline

Definition at line 221 of file G4GlauberGribovCrossSection.hh.

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}

Referenced by GetIsoCrossSection().

◆ GetParticleBarCorTot()

G4double G4GlauberGribovCrossSection::GetParticleBarCorTot ( const G4ParticleDefinition theParticle,
G4int  Z 
)
inline

Definition at line 201 of file G4GlauberGribovCrossSection.hh.

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}

Referenced by GetIsoCrossSection().

◆ GetProductionGlauberGribovXsc()

G4double G4GlauberGribovCrossSection::GetProductionGlauberGribovXsc ( )
inline

Definition at line 106 of file G4GlauberGribovCrossSection.hh.

106{ return fProductionXsc; };

◆ GetRadiusConst()

G4double G4GlauberGribovCrossSection::GetRadiusConst ( )
inline

Definition at line 108 of file G4GlauberGribovCrossSection.hh.

108{ return fRadiusConst; };

◆ GetRatioQE()

G4double G4GlauberGribovCrossSection::GetRatioQE ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 474 of file G4GlauberGribovCrossSection.cc.

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

◆ GetRatioSD()

G4double G4GlauberGribovCrossSection::GetRatioSD ( const G4DynamicParticle aParticle,
G4int  At,
G4int  Zt 
)

Definition at line 432 of file G4GlauberGribovCrossSection.cc.

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

◆ GetTotalGlauberGribovXsc()

G4double G4GlauberGribovCrossSection::GetTotalGlauberGribovXsc ( )
inline

Definition at line 103 of file G4GlauberGribovCrossSection.hh.

103{ return fTotalXsc; };

◆ IsIsoApplicable()

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 283 of file G4GlauberGribovCrossSection.cc.

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

◆ SetEnergyLowerLimit()

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: