#include <G4GlauberGribovCrossSection.hh>
Inheritance diagram for G4GlauberGribovCrossSection:
Definition at line 55 of file G4GlauberGribovCrossSection.hh.
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().
00231 : G4VCrossSectionDataSet(Default_Name()), 00232 fUpperLimit(100000*GeV), fLowerLimit(10.*MeV),// fLowerLimit(3*GeV), 00233 fRadiusConst(1.08*fermi), // 1.1, 1.3 ? 00234 fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0), fProductionXsc(0.0), 00235 fDiffractionXsc(0.0), fHadronNucleonXsc(0.0) 00236 { 00237 theGamma = G4Gamma::Gamma(); 00238 theProton = G4Proton::Proton(); 00239 theNeutron = G4Neutron::Neutron(); 00240 theAProton = G4AntiProton::AntiProton(); 00241 theANeutron = G4AntiNeutron::AntiNeutron(); 00242 thePiPlus = G4PionPlus::PionPlus(); 00243 thePiMinus = G4PionMinus::PionMinus(); 00244 thePiZero = G4PionZero::PionZero(); 00245 theKPlus = G4KaonPlus::KaonPlus(); 00246 theKMinus = G4KaonMinus::KaonMinus(); 00247 theK0S = G4KaonZeroShort::KaonZeroShort(); 00248 theK0L = G4KaonZeroLong::KaonZeroLong(); 00249 theL = G4Lambda::Lambda(); 00250 theAntiL = G4AntiLambda::AntiLambda(); 00251 theSPlus = G4SigmaPlus::SigmaPlus(); 00252 theASPlus = G4AntiSigmaPlus::AntiSigmaPlus(); 00253 theSMinus = G4SigmaMinus::SigmaMinus(); 00254 theASMinus = G4AntiSigmaMinus::AntiSigmaMinus(); 00255 theS0 = G4SigmaZero::SigmaZero(); 00256 theAS0 = G4AntiSigmaZero::AntiSigmaZero(); 00257 theXiMinus = G4XiMinus::XiMinus(); 00258 theXi0 = G4XiZero::XiZero(); 00259 theAXiMinus = G4AntiXiMinus::AntiXiMinus(); 00260 theAXi0 = G4AntiXiZero::AntiXiZero(); 00261 theOmega = G4OmegaMinus::OmegaMinus(); 00262 theAOmega = G4AntiOmegaMinus::AntiOmegaMinus(); 00263 theD = G4Deuteron::Deuteron(); 00264 theT = G4Triton::Triton(); 00265 theA = G4Alpha::Alpha(); 00266 theHe3 = G4He3::He3(); 00267 00268 hnXsc = new G4HadronNucleonXsc(); 00269 }
G4GlauberGribovCrossSection::~G4GlauberGribovCrossSection | ( | ) | [virtual] |
G4double G4GlauberGribovCrossSection::CalcMandelstamS | ( | const | G4double, | |
const | G4double, | |||
const | G4double | |||
) |
Definition at line 1507 of file G4GlauberGribovCrossSection.cc.
Referenced by GetHadronNucleonXsc(), GetHadronNucleonXscNS(), and GetHadronNucleonXscPDG().
01510 { 01511 G4double Elab = std::sqrt ( mp * mp + Plab * Plab ); 01512 G4double sMand = mp*mp + mt*mt + 2*Elab*mt ; 01513 01514 return sMand; 01515 }
G4double G4GlauberGribovCrossSection::CalculateEcmValue | ( | const | G4double, | |
const | G4double, | |||
const | G4double | |||
) |
Definition at line 1491 of file G4GlauberGribovCrossSection.cc.
01494 { 01495 G4double Elab = std::sqrt ( mp * mp + Plab * Plab ); 01496 G4double Ecm = std::sqrt ( mp * mp + mt * mt + 2 * Elab * mt ); 01497 // G4double Pcm = Plab * mt / Ecm; 01498 // G4double KEcm = std::sqrt ( Pcm * Pcm + mp * mp ) - mp; 01499 01500 return Ecm ; // KEcm; 01501 }
void G4GlauberGribovCrossSection::CrossSectionDescription | ( | std::ostream & | ) | const [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 1521 of file G4GlauberGribovCrossSection.cc.
01522 { 01523 outFile << "G4GlauberGribovCrossSection calculates total, inelastic and\n" 01524 << "elastic cross sections for hadron-nucleus cross sections using\n" 01525 << "the Glauber model with Gribov corrections. It is valid for all\n" 01526 << "targets except hydrogen, and for incident p, pbar, n, sigma-,\n" 01527 << "pi+, pi-, K+, K- and gammas with energies above 3 GeV. This is\n" 01528 << "a cross section component which is to be used to build a cross\n" 01529 << "data set.\n"; 01530 }
static const char* G4GlauberGribovCrossSection::Default_Name | ( | ) | [inline, static] |
Definition at line 62 of file G4GlauberGribovCrossSection.hh.
Referenced by G4BGGNucleonInelasticXS::BuildPhysicsTable(), and G4BGGNucleonElasticXS::BuildPhysicsTable().
G4double G4GlauberGribovCrossSection::GetDiffractionGlauberGribovXsc | ( | ) | [inline] |
G4double G4GlauberGribovCrossSection::GetElasticGlauberGribov | ( | const G4DynamicParticle * | , | |
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().
00179 { 00180 GetIsoCrossSection(dp, Z, A); 00181 return fElasticXsc; 00182 }
G4double G4GlauberGribovCrossSection::GetElasticGlauberGribovXsc | ( | ) | [inline] |
Definition at line 104 of file G4GlauberGribovCrossSection.hh.
Referenced by G4NeutronElasticXS::GetElementCrossSection().
G4double G4GlauberGribovCrossSection::GetHadronNucleonXsc | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 539 of file G4GlauberGribovCrossSection.cc.
References CalcMandelstamS(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetMass(), and G4DynamicParticle::GetMomentum().
00541 { 00542 G4double xsection; 00543 00544 //G4double targ_mass = G4NucleiProperties::GetNuclearMass(At, Zt); 00545 00546 G4double targ_mass = 0.939*GeV; // ~mean neutron and proton ??? 00547 00548 G4double proj_mass = aParticle->GetMass(); 00549 G4double proj_momentum = aParticle->GetMomentum().mag(); 00550 G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum ); 00551 00552 sMand /= GeV*GeV; // in GeV for parametrisation 00553 proj_momentum /= GeV; 00554 00555 const G4ParticleDefinition* theParticle = aParticle->GetDefinition(); 00556 00557 G4double aa = At; 00558 00559 if(theParticle == theGamma) 00560 { 00561 xsection = aa*(0.0677*std::pow(sMand,0.0808) + 0.129*std::pow(sMand,-0.4525)); 00562 } 00563 else if(theParticle == theNeutron) // as proton ??? 00564 { 00565 xsection = aa*(21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525)); 00566 } 00567 else if(theParticle == theProton) 00568 { 00569 xsection = aa*(21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525)); 00570 // xsection = At*( 49.51*std::pow(sMand,-0.097) + 0.314*std::log(sMand)*std::log(sMand) ); 00571 // xsection = At*( 38.4 + 0.85*std::abs(std::pow(log(sMand),1.47)) ); 00572 } 00573 else if(theParticle == theAProton) 00574 { 00575 xsection = aa*( 21.70*std::pow(sMand,0.0808) + 98.39*std::pow(sMand,-0.4525)); 00576 } 00577 else if(theParticle == thePiPlus) 00578 { 00579 xsection = aa*(13.63*std::pow(sMand,0.0808) + 27.56*std::pow(sMand,-0.4525)); 00580 } 00581 else if(theParticle == thePiMinus) 00582 { 00583 // xsection = At*( 55.2*std::pow(sMand,-0.255) + 0.346*std::log(sMand)*std::log(sMand) ); 00584 xsection = aa*(13.63*std::pow(sMand,0.0808) + 36.02*std::pow(sMand,-0.4525)); 00585 } 00586 else if(theParticle == theKPlus) 00587 { 00588 xsection = aa*(11.82*std::pow(sMand,0.0808) + 8.15*std::pow(sMand,-0.4525)); 00589 } 00590 else if(theParticle == theKMinus) 00591 { 00592 xsection = aa*(11.82*std::pow(sMand,0.0808) + 26.36*std::pow(sMand,-0.4525)); 00593 } 00594 else // as proton ??? 00595 { 00596 xsection = aa*(21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525)); 00597 } 00598 xsection *= millibarn; 00599 return xsection; 00600 }
G4double G4GlauberGribovCrossSection::GetHadronNucleonXsc | ( | const G4DynamicParticle * | , | |
const G4Element * | ||||
) |
Definition at line 522 of file G4GlauberGribovCrossSection.cc.
References G4lrint(), G4Element::GetN(), and G4Element::GetZ().
00524 { 00525 G4int At = G4lrint(anElement->GetN()); // number of nucleons 00526 G4int Zt = G4lrint(anElement->GetZ()); // number of protons 00527 00528 return GetHadronNucleonXsc(aParticle, At, Zt); 00529 }
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscNS | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 761 of file G4GlauberGribovCrossSection.cc.
References CalcMandelstamS(), G4DynamicParticle::GetDefinition(), GetHadronNucleonXscPDG(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentum(), G4ParticleTable::GetParticleTable(), G4DynamicParticle::GetTotalEnergy(), G4InuclParticleNames::nn, and G4InuclParticleNames::s0.
00763 { 00764 G4double xsection(0); 00765 // G4double Delta; DHW 19 May 2011: variable set but not used 00766 G4double A0, B0; 00767 G4double hpXscv(0); 00768 G4double hnXscv(0); 00769 00770 G4int Nt = At-Zt; // number of neutrons 00771 if (Nt < 0) Nt = 0; 00772 00773 G4double aa = At; 00774 G4double zz = Zt; 00775 G4double nn = Nt; 00776 00777 G4double targ_mass = G4ParticleTable::GetParticleTable()-> 00778 GetIonTable()->GetIonMass(Zt, At); 00779 00780 targ_mass = 0.939*GeV; // ~mean neutron and proton ??? 00781 00782 G4double proj_mass = aParticle->GetMass(); 00783 G4double proj_energy = aParticle->GetTotalEnergy(); 00784 G4double proj_momentum = aParticle->GetMomentum().mag(); 00785 00786 G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum ); 00787 00788 sMand /= GeV*GeV; // in GeV for parametrisation 00789 proj_momentum /= GeV; 00790 proj_energy /= GeV; 00791 proj_mass /= GeV; 00792 00793 // General PDG fit constants 00794 00795 G4double s0 = 5.38*5.38; // in Gev^2 00796 G4double eta1 = 0.458; 00797 G4double eta2 = 0.458; 00798 G4double B = 0.308; 00799 00800 00801 const G4ParticleDefinition* theParticle = aParticle->GetDefinition(); 00802 00803 00804 if(theParticle == theNeutron) 00805 { 00806 if( proj_momentum >= 373.) 00807 { 00808 return GetHadronNucleonXscPDG(aParticle,At,Zt); 00809 } 00810 else if( proj_momentum >= 10.) 00811 // if( proj_momentum >= 2.) 00812 { 00813 // Delta = 1.; // DHW 19 May 2011: variable set but not used 00814 // if( proj_energy < 40. ) Delta = 0.916+0.0021*proj_energy; 00815 00816 if(proj_momentum >= 10.) 00817 { 00818 B0 = 7.5; 00819 A0 = 100. - B0*std::log(3.0e7); 00820 00821 xsection = A0 + B0*std::log(proj_energy) - 11 00822 + 103*std::pow(2*0.93827*proj_energy + proj_mass*proj_mass+ 00823 0.93827*0.93827,-0.165); // mb 00824 } 00825 xsection *= zz + nn; 00826 } 00827 else 00828 { 00829 // nn to be pp 00830 00831 if( proj_momentum < 0.73 ) 00832 { 00833 hnXscv = 23 + 50*( std::pow( std::log(0.73/proj_momentum), 3.5 ) ); 00834 } 00835 else if( proj_momentum < 1.05 ) 00836 { 00837 hnXscv = 23 + 40*(std::log(proj_momentum/0.73))* 00838 (std::log(proj_momentum/0.73)); 00839 } 00840 else // if( proj_momentum < 10. ) 00841 { 00842 hnXscv = 39.0+ 00843 75*(proj_momentum - 1.2)/(std::pow(proj_momentum,3.0) + 0.15); 00844 } 00845 // pn to be np 00846 00847 if( proj_momentum < 0.8 ) 00848 { 00849 hpXscv = 33+30*std::pow(std::log(proj_momentum/1.3),4.0); 00850 } 00851 else if( proj_momentum < 1.4 ) 00852 { 00853 hpXscv = 33+30*std::pow(std::log(proj_momentum/0.95),2.0); 00854 } 00855 else // if( proj_momentum < 10. ) 00856 { 00857 hpXscv = 33.3+ 00858 20.8*(std::pow(proj_momentum,2.0)-1.35)/ 00859 (std::pow(proj_momentum,2.50)+0.95); 00860 } 00861 xsection = hpXscv*zz + hnXscv*nn; 00862 } 00863 } 00864 else if(theParticle == theProton) 00865 { 00866 if( proj_momentum >= 373.) 00867 { 00868 return GetHadronNucleonXscPDG(aParticle,At,Zt); 00869 } 00870 else if( proj_momentum >= 10.) 00871 // if( proj_momentum >= 2.) 00872 { 00873 // Delta = 1.; DHW 19 May 2011: variable set but not used 00874 // if( proj_energy < 40. ) Delta = 0.916+0.0021*proj_energy; 00875 00876 if(proj_momentum >= 10.) 00877 { 00878 B0 = 7.5; 00879 A0 = 100. - B0*std::log(3.0e7); 00880 00881 xsection = A0 + B0*std::log(proj_energy) - 11 00882 + 103*std::pow(2*0.93827*proj_energy + proj_mass*proj_mass+ 00883 0.93827*0.93827,-0.165); // mb 00884 } 00885 xsection *= zz + nn; 00886 } 00887 else 00888 { 00889 // pp 00890 00891 if( proj_momentum < 0.73 ) 00892 { 00893 hpXscv = 23 + 50*( std::pow( std::log(0.73/proj_momentum), 3.5 ) ); 00894 } 00895 else if( proj_momentum < 1.05 ) 00896 { 00897 hpXscv = 23 + 40*(std::log(proj_momentum/0.73))* 00898 (std::log(proj_momentum/0.73)); 00899 } 00900 else // if( proj_momentum < 10. ) 00901 { 00902 hpXscv = 39.0+ 00903 75*(proj_momentum - 1.2)/(std::pow(proj_momentum,3.0) + 0.15); 00904 } 00905 // pn to be np 00906 00907 if( proj_momentum < 0.8 ) 00908 { 00909 hnXscv = 33+30*std::pow(std::log(proj_momentum/1.3),4.0); 00910 } 00911 else if( proj_momentum < 1.4 ) 00912 { 00913 hnXscv = 33+30*std::pow(std::log(proj_momentum/0.95),2.0); 00914 } 00915 else // if( proj_momentum < 10. ) 00916 { 00917 hnXscv = 33.3+ 00918 20.8*(std::pow(proj_momentum,2.0)-1.35)/ 00919 (std::pow(proj_momentum,2.50)+0.95); 00920 } 00921 xsection = hpXscv*zz + hnXscv*nn; 00922 // xsection = hpXscv*(Zt + Nt); 00923 // xsection = hnXscv*(Zt + Nt); 00924 } 00925 // xsection *= 0.95; 00926 } 00927 else if( theParticle == theAProton ) 00928 { 00929 // xsection = Zt*( 35.45 + B*std::pow(std::log(sMand/s0),2.) 00930 // + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2)); 00931 00932 // xsection += Nt*( 35.80 + B*std::pow(std::log(sMand/s0),2.) 00933 // + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2)); 00934 00935 G4double logP = std::log(proj_momentum); 00936 00937 if( proj_momentum <= 1.0 ) 00938 { 00939 xsection = zz*(65.55 + 53.84/(proj_momentum+1.e-6) ); 00940 } 00941 else 00942 { 00943 xsection = zz*( 41.1 + 77.2*std::pow( proj_momentum, -0.68) 00944 + 0.293*logP*logP - 1.82*logP ); 00945 } 00946 if ( nn > 0.) 00947 { 00948 xsection += nn*( 41.9 + 96.2*std::pow( proj_momentum, -0.99) - 0.154*logP); 00949 } 00950 else // H 00951 { 00952 fInelasticXsc = 38.0 + 38.0*std::pow( proj_momentum, -0.96) 00953 - 0.169*logP*logP; 00954 fInelasticXsc *= millibarn; 00955 } 00956 } 00957 else if( theParticle == thePiPlus ) 00958 { 00959 if(proj_momentum < 0.4) 00960 { 00961 G4double Ex3 = 180*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.085/0.085); 00962 hpXscv = Ex3+20.0; 00963 } 00964 else if( proj_momentum < 1.15 ) 00965 { 00966 G4double Ex4 = 88*(std::log(proj_momentum/0.75))*(std::log(proj_momentum/0.75)); 00967 hpXscv = Ex4+14.0; 00968 } 00969 else if(proj_momentum < 3.5) 00970 { 00971 G4double Ex1 = 3.2*std::exp(-(proj_momentum-2.55)*(proj_momentum-2.55)/0.55/0.55); 00972 G4double Ex2 = 12*std::exp(-(proj_momentum-1.47)*(proj_momentum-1.47)/0.225/0.225); 00973 hpXscv = Ex1+Ex2+27.5; 00974 } 00975 else // if(proj_momentum > 3.5) // mb 00976 { 00977 hpXscv = 10.6+2.*std::log(proj_energy)+25*std::pow(proj_energy,-0.43); 00978 } 00979 // pi+n = pi-p?? 00980 00981 if(proj_momentum < 0.37) 00982 { 00983 hnXscv = 28.0 + 40*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.07/0.07); 00984 } 00985 else if(proj_momentum<0.65) 00986 { 00987 hnXscv = 26+110*(std::log(proj_momentum/0.48))*(std::log(proj_momentum/0.48)); 00988 } 00989 else if(proj_momentum<1.3) 00990 { 00991 hnXscv = 36.1+ 00992 10*std::exp(-(proj_momentum-0.72)*(proj_momentum-0.72)/0.06/0.06)+ 00993 24*std::exp(-(proj_momentum-1.015)*(proj_momentum-1.015)/0.075/0.075); 00994 } 00995 else if(proj_momentum<3.0) 00996 { 00997 hnXscv = 36.1+0.079-4.313*std::log(proj_momentum)+ 00998 3*std::exp(-(proj_momentum-2.1)*(proj_momentum-2.1)/0.4/0.4)+ 00999 1.5*std::exp(-(proj_momentum-1.4)*(proj_momentum-1.4)/0.12/0.12); 01000 } 01001 else // mb 01002 { 01003 hnXscv = 10.6+2*std::log(proj_energy)+30*std::pow(proj_energy,-0.43); 01004 } 01005 xsection = hpXscv*zz + hnXscv*nn; 01006 } 01007 else if(theParticle == thePiMinus) 01008 { 01009 // pi-n = pi+p?? 01010 01011 if(proj_momentum < 0.4) 01012 { 01013 G4double Ex3 = 180*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.085/0.085); 01014 hnXscv = Ex3+20.0; 01015 } 01016 else if(proj_momentum < 1.15) 01017 { 01018 G4double Ex4 = 88*(std::log(proj_momentum/0.75))*(std::log(proj_momentum/0.75)); 01019 hnXscv = Ex4+14.0; 01020 } 01021 else if(proj_momentum < 3.5) 01022 { 01023 G4double Ex1 = 3.2*std::exp(-(proj_momentum-2.55)*(proj_momentum-2.55)/0.55/0.55); 01024 G4double Ex2 = 12*std::exp(-(proj_momentum-1.47)*(proj_momentum-1.47)/0.225/0.225); 01025 hnXscv = Ex1+Ex2+27.5; 01026 } 01027 else // if(proj_momentum > 3.5) // mb 01028 { 01029 hnXscv = 10.6+2.*std::log(proj_energy)+25*std::pow(proj_energy,-0.43); 01030 } 01031 // pi-p 01032 01033 if(proj_momentum < 0.37) 01034 { 01035 hpXscv = 28.0 + 40*std::exp(-(proj_momentum-0.29)*(proj_momentum-0.29)/0.07/0.07); 01036 } 01037 else if(proj_momentum<0.65) 01038 { 01039 hpXscv = 26+110*(std::log(proj_momentum/0.48))*(std::log(proj_momentum/0.48)); 01040 } 01041 else if(proj_momentum<1.3) 01042 { 01043 hpXscv = 36.1+ 01044 10*std::exp(-(proj_momentum-0.72)*(proj_momentum-0.72)/0.06/0.06)+ 01045 24*std::exp(-(proj_momentum-1.015)*(proj_momentum-1.015)/0.075/0.075); 01046 } 01047 else if(proj_momentum<3.0) 01048 { 01049 hpXscv = 36.1+0.079-4.313*std::log(proj_momentum)+ 01050 3*std::exp(-(proj_momentum-2.1)*(proj_momentum-2.1)/0.4/0.4)+ 01051 1.5*std::exp(-(proj_momentum-1.4)*(proj_momentum-1.4)/0.12/0.12); 01052 } 01053 else // mb 01054 { 01055 hpXscv = 10.6+2*std::log(proj_energy)+30*std::pow(proj_energy,-0.43); 01056 } 01057 xsection = hpXscv*zz + hnXscv*nn; 01058 } 01059 else if(theParticle == theKPlus) 01060 { 01061 xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.) 01062 + 7.14*std::pow(sMand,-eta1) - 13.45*std::pow(sMand,-eta2)); 01063 01064 xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.) 01065 + 5.17*std::pow(sMand,-eta1) - 7.23*std::pow(sMand,-eta2)); 01066 } 01067 else if(theParticle == theKMinus) 01068 { 01069 xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.) 01070 + 7.14*std::pow(sMand,-eta1) + 13.45*std::pow(sMand,-eta2)); 01071 01072 xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.) 01073 + 5.17*std::pow(sMand,-eta1) + 7.23*std::pow(sMand,-eta2)); 01074 } 01075 else if(theParticle == theSMinus) 01076 { 01077 xsection = aa*( 35.20 + B*std::pow(std::log(sMand/s0),2.) 01078 - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2)); 01079 } 01080 else if(theParticle == theGamma) // modify later on 01081 { 01082 xsection = aa*( 0.0 + B*std::pow(std::log(sMand/s0),2.) 01083 + 0.032*std::pow(sMand,-eta1) - 0.0*std::pow(sMand,-eta2)); 01084 01085 } 01086 else // as proton ??? 01087 { 01088 xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.) 01089 + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2)); 01090 01091 xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.) 01092 + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2)); 01093 } 01094 xsection *= millibarn; // parametrised in mb 01095 return xsection; 01096 }
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscNS | ( | const G4DynamicParticle * | , | |
const G4Element * | ||||
) |
Definition at line 743 of file G4GlauberGribovCrossSection.cc.
References G4lrint(), G4Element::GetN(), and G4Element::GetZ().
Referenced by GetHNinelasticXsc(), GetIsoCrossSection(), GetRatioQE(), and GetRatioSD().
00745 { 00746 G4int At = G4lrint(anElement->GetN()); // number of nucleons 00747 G4int Zt = G4lrint(anElement->GetZ()); // number of protons 00748 00749 return GetHadronNucleonXscNS(aParticle, At, Zt); 00750 }
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscPDG | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 628 of file G4GlauberGribovCrossSection.cc.
References CalcMandelstamS(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentum(), G4ParticleTable::GetParticleTable(), G4InuclParticleNames::nn, and G4InuclParticleNames::s0.
00630 { 00631 G4double xsection; 00632 00633 G4int Nt = At-Zt; // number of neutrons 00634 if (Nt < 0) Nt = 0; 00635 00636 G4double zz = Zt; 00637 G4double aa = At; 00638 G4double nn = Nt; 00639 00640 G4double targ_mass = G4ParticleTable::GetParticleTable()-> 00641 GetIonTable()->GetIonMass(Zt, At); 00642 00643 targ_mass = 0.939*GeV; // ~mean neutron and proton ??? 00644 00645 G4double proj_mass = aParticle->GetMass(); 00646 G4double proj_momentum = aParticle->GetMomentum().mag(); 00647 00648 G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum ); 00649 00650 sMand /= GeV*GeV; // in GeV for parametrisation 00651 00652 // General PDG fit constants 00653 00654 G4double s0 = 5.38*5.38; // in Gev^2 00655 G4double eta1 = 0.458; 00656 G4double eta2 = 0.458; 00657 G4double B = 0.308; 00658 00659 00660 const G4ParticleDefinition* theParticle = aParticle->GetDefinition(); 00661 00662 00663 if(theParticle == theNeutron) // proton-neutron fit 00664 { 00665 xsection = zz*( 35.80 + B*std::pow(std::log(sMand/s0),2.) 00666 + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2)); 00667 xsection += nn*( 35.45 + B*std::pow(std::log(sMand/s0),2.) 00668 + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2)); // pp for nn 00669 } 00670 else if(theParticle == theProton) 00671 { 00672 00673 xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.) 00674 + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2)); 00675 00676 xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.) 00677 + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2)); 00678 } 00679 else if(theParticle == theAProton) 00680 { 00681 xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.) 00682 + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2)); 00683 00684 xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.) 00685 + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2)); 00686 } 00687 else if(theParticle == thePiPlus) 00688 { 00689 xsection = aa*( 20.86 + B*std::pow(std::log(sMand/s0),2.) 00690 + 19.24*std::pow(sMand,-eta1) - 6.03*std::pow(sMand,-eta2)); 00691 } 00692 else if(theParticle == thePiMinus) 00693 { 00694 xsection = aa*( 20.86 + B*std::pow(std::log(sMand/s0),2.) 00695 + 19.24*std::pow(sMand,-eta1) + 6.03*std::pow(sMand,-eta2)); 00696 } 00697 else if(theParticle == theKPlus || theParticle == theK0L ) 00698 { 00699 xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.) 00700 + 7.14*std::pow(sMand,-eta1) - 13.45*std::pow(sMand,-eta2)); 00701 00702 xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.) 00703 + 5.17*std::pow(sMand,-eta1) - 7.23*std::pow(sMand,-eta2)); 00704 } 00705 else if(theParticle == theKMinus || theParticle == theK0S ) 00706 { 00707 xsection = zz*( 17.91 + B*std::pow(std::log(sMand/s0),2.) 00708 + 7.14*std::pow(sMand,-eta1) + 13.45*std::pow(sMand,-eta2)); 00709 00710 xsection += nn*( 17.87 + B*std::pow(std::log(sMand/s0),2.) 00711 + 5.17*std::pow(sMand,-eta1) + 7.23*std::pow(sMand,-eta2)); 00712 } 00713 else if(theParticle == theSMinus) 00714 { 00715 xsection = aa*( 35.20 + B*std::pow(std::log(sMand/s0),2.) 00716 - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2)); 00717 } 00718 else if(theParticle == theGamma) // modify later on 00719 { 00720 xsection = aa*( 0.0 + B*std::pow(std::log(sMand/s0),2.) 00721 + 0.032*std::pow(sMand,-eta1) - 0.0*std::pow(sMand,-eta2)); 00722 00723 } 00724 else // as proton ??? 00725 { 00726 xsection = zz*( 35.45 + B*std::pow(std::log(sMand/s0),2.) 00727 + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2)); 00728 00729 xsection += nn*( 35.80 + B*std::pow(std::log(sMand/s0),2.) 00730 + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2)); 00731 } 00732 xsection *= millibarn; // parametrised in mb 00733 return xsection; 00734 }
G4double G4GlauberGribovCrossSection::GetHadronNucleonXscPDG | ( | const G4DynamicParticle * | , | |
const G4Element * | ||||
) |
Definition at line 609 of file G4GlauberGribovCrossSection.cc.
References G4lrint(), G4Element::GetN(), and G4Element::GetZ().
Referenced by GetHadronNucleonXscNS(), and GetKaonNucleonXscVector().
00611 { 00612 G4int At = G4lrint(anElement->GetN()); // number of nucleons 00613 G4int Zt = G4lrint(anElement->GetZ()); // number of protons 00614 00615 return GetHadronNucleonXscPDG(aParticle, At, Zt); 00616 }
G4double G4GlauberGribovCrossSection::GetHNinelasticXsc | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 1151 of file G4GlauberGribovCrossSection.cc.
References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscNS(), and GetHNinelasticXscVU().
01153 { 01154 G4ParticleDefinition* hadron = aParticle->GetDefinition(); 01155 G4double sumInelastic; 01156 G4int Nt = At - Zt; 01157 if(Nt < 0) Nt = 0; 01158 01159 if( hadron == theKPlus ) 01160 { 01161 sumInelastic = GetHNinelasticXscVU(aParticle, At, Zt); 01162 } 01163 else 01164 { 01165 //sumInelastic = Zt*GetHadronNucleonXscMK(aParticle, theProton); 01166 // sumInelastic += Nt*GetHadronNucleonXscMK(aParticle, theNeutron); 01167 sumInelastic = G4double(Zt)*GetHadronNucleonXscNS(aParticle, 1, 1); 01168 sumInelastic += G4double(Nt)*GetHadronNucleonXscNS(aParticle, 1, 0); 01169 } 01170 return sumInelastic; 01171 }
G4double G4GlauberGribovCrossSection::GetHNinelasticXsc | ( | const G4DynamicParticle * | , | |
const G4Element * | ||||
) |
Definition at line 1137 of file G4GlauberGribovCrossSection.cc.
References G4lrint(), G4Element::GetN(), and G4Element::GetZ().
Referenced by GetIsoCrossSection(), and GetRatioQE().
01139 { 01140 G4int At = G4lrint(anElement->GetN()); // number of nucleons 01141 G4int Zt = G4lrint(anElement->GetZ()); // number of protons 01142 01143 return GetHNinelasticXsc(aParticle, At, Zt); 01144 }
G4double G4GlauberGribovCrossSection::GetHNinelasticXscVU | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 1179 of file G4GlauberGribovCrossSection.cc.
References G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetMomentum(), G4ParticleDefinition::GetPDGEncoding(), and G4DynamicParticle::GetTotalEnergy().
Referenced by GetHNinelasticXsc().
01181 { 01182 G4int PDGcode = aParticle->GetDefinition()->GetPDGEncoding(); 01183 G4int absPDGcode = std::abs(PDGcode); 01184 01185 G4double Elab = aParticle->GetTotalEnergy(); 01186 // (s - 2*0.88*GeV*GeV)/(2*0.939*GeV)/GeV; 01187 G4double Plab = aParticle->GetMomentum().mag(); 01188 // std::sqrt(Elab * Elab - 0.88); 01189 01190 Elab /= GeV; 01191 Plab /= GeV; 01192 01193 G4double LogPlab = std::log( Plab ); 01194 G4double sqrLogPlab = LogPlab * LogPlab; 01195 01196 //G4cout<<"Plab = "<<Plab<<G4endl; 01197 01198 G4double NumberOfTargetProtons = G4double(Zt); 01199 G4double NumberOfTargetNucleons = G4double(At); 01200 G4double NumberOfTargetNeutrons = NumberOfTargetNucleons - NumberOfTargetProtons; 01201 01202 if(NumberOfTargetNeutrons < 0.0) NumberOfTargetNeutrons = 0.0; 01203 01204 G4double Xtotal, Xelastic, Xinelastic; 01205 01206 if( absPDGcode > 1000 ) //------Projectile is baryon -------- 01207 { 01208 G4double XtotPP = 48.0 + 0. *std::pow(Plab, 0. ) + 01209 0.522*sqrLogPlab - 4.51*LogPlab; 01210 01211 G4double XtotPN = 47.3 + 0. *std::pow(Plab, 0. ) + 01212 0.513*sqrLogPlab - 4.27*LogPlab; 01213 01214 G4double XelPP = 11.9 + 26.9*std::pow(Plab,-1.21) + 01215 0.169*sqrLogPlab - 1.85*LogPlab; 01216 01217 G4double XelPN = 11.9 + 26.9*std::pow(Plab,-1.21) + 01218 0.169*sqrLogPlab - 1.85*LogPlab; 01219 01220 Xtotal = (NumberOfTargetProtons * XtotPP + 01221 NumberOfTargetNeutrons * XtotPN); 01222 01223 Xelastic = (NumberOfTargetProtons * XelPP + 01224 NumberOfTargetNeutrons * XelPN); 01225 } 01226 else if( PDGcode == 211 ) //------Projectile is PionPlus ------- 01227 { 01228 G4double XtotPiP = 16.4 + 19.3 *std::pow(Plab,-0.42) + 01229 0.19 *sqrLogPlab - 0.0 *LogPlab; 01230 01231 G4double XtotPiN = 33.0 + 14.0 *std::pow(Plab,-1.36) + 01232 0.456*sqrLogPlab - 4.03*LogPlab; 01233 01234 G4double XelPiP = 0.0 + 11.4*std::pow(Plab,-0.40) + 01235 0.079*sqrLogPlab - 0.0 *LogPlab; 01236 01237 G4double XelPiN = 1.76 + 11.2*std::pow(Plab,-0.64) + 01238 0.043*sqrLogPlab - 0.0 *LogPlab; 01239 01240 Xtotal = ( NumberOfTargetProtons * XtotPiP + 01241 NumberOfTargetNeutrons * XtotPiN ); 01242 01243 Xelastic = ( NumberOfTargetProtons * XelPiP + 01244 NumberOfTargetNeutrons * XelPiN ); 01245 } 01246 else if( PDGcode == -211 ) //------Projectile is PionMinus ------- 01247 { 01248 G4double XtotPiP = 33.0 + 14.0 *std::pow(Plab,-1.36) + 01249 0.456*sqrLogPlab - 4.03*LogPlab; 01250 01251 G4double XtotPiN = 16.4 + 19.3 *std::pow(Plab,-0.42) + 01252 0.19 *sqrLogPlab - 0.0 *LogPlab; 01253 01254 G4double XelPiP = 1.76 + 11.2*std::pow(Plab,-0.64) + 01255 0.043*sqrLogPlab - 0.0 *LogPlab; 01256 01257 G4double XelPiN = 0.0 + 11.4*std::pow(Plab,-0.40) + 01258 0.079*sqrLogPlab - 0.0 *LogPlab; 01259 01260 Xtotal = ( NumberOfTargetProtons * XtotPiP + 01261 NumberOfTargetNeutrons * XtotPiN ); 01262 01263 Xelastic = ( NumberOfTargetProtons * XelPiP + 01264 NumberOfTargetNeutrons * XelPiN ); 01265 } 01266 else if( PDGcode == 111 ) //------Projectile is PionZero ------- 01267 { 01268 G4double XtotPiP =(16.4 + 19.3 *std::pow(Plab,-0.42) + 01269 0.19 *sqrLogPlab - 0.0 *LogPlab + //Pi+ 01270 33.0 + 14.0 *std::pow(Plab,-1.36) + 01271 0.456*sqrLogPlab - 4.03*LogPlab)/2; //Pi- 01272 01273 G4double XtotPiN =(33.0 + 14.0 *std::pow(Plab,-1.36) + 01274 0.456*sqrLogPlab - 4.03*LogPlab + //Pi+ 01275 16.4 + 19.3 *std::pow(Plab,-0.42) + 01276 0.19 *sqrLogPlab - 0.0 *LogPlab)/2; //Pi- 01277 01278 G4double XelPiP =( 0.0 + 11.4*std::pow(Plab,-0.40) + 01279 0.079*sqrLogPlab - 0.0 *LogPlab + //Pi+ 01280 1.76 + 11.2*std::pow(Plab,-0.64) + 01281 0.043*sqrLogPlab - 0.0 *LogPlab)/2; //Pi- 01282 01283 G4double XelPiN =( 1.76 + 11.2*std::pow(Plab,-0.64) + 01284 0.043*sqrLogPlab - 0.0 *LogPlab + //Pi+ 01285 0.0 + 11.4*std::pow(Plab,-0.40) + 01286 0.079*sqrLogPlab - 0.0 *LogPlab)/2; //Pi- 01287 01288 Xtotal = ( NumberOfTargetProtons * XtotPiP + 01289 NumberOfTargetNeutrons * XtotPiN ); 01290 01291 Xelastic = ( NumberOfTargetProtons * XelPiP + 01292 NumberOfTargetNeutrons * XelPiN ); 01293 } 01294 else if( PDGcode == 321 ) //------Projectile is KaonPlus ------- 01295 { 01296 G4double XtotKP = 18.1 + 0. *std::pow(Plab, 0. ) + 01297 0.26 *sqrLogPlab - 1.0 *LogPlab; 01298 G4double XtotKN = 18.7 + 0. *std::pow(Plab, 0. ) + 01299 0.21 *sqrLogPlab - 0.89*LogPlab; 01300 01301 G4double XelKP = 5.0 + 8.1*std::pow(Plab,-1.8 ) + 01302 0.16 *sqrLogPlab - 1.3 *LogPlab; 01303 01304 G4double XelKN = 7.3 + 0. *std::pow(Plab,-0. ) + 01305 0.29 *sqrLogPlab - 2.4 *LogPlab; 01306 01307 Xtotal = ( NumberOfTargetProtons * XtotKP + 01308 NumberOfTargetNeutrons * XtotKN ); 01309 01310 Xelastic = ( NumberOfTargetProtons * XelKP + 01311 NumberOfTargetNeutrons * XelKN ); 01312 } 01313 else if( PDGcode ==-321 ) //------Projectile is KaonMinus ------ 01314 { 01315 G4double XtotKP = 32.1 + 0. *std::pow(Plab, 0. ) + 01316 0.66 *sqrLogPlab - 5.6 *LogPlab; 01317 G4double XtotKN = 25.2 + 0. *std::pow(Plab, 0. ) + 01318 0.38 *sqrLogPlab - 2.9 *LogPlab; 01319 01320 G4double XelKP = 7.3 + 0. *std::pow(Plab,-0. ) + 01321 0.29 *sqrLogPlab - 2.4 *LogPlab; 01322 01323 G4double XelKN = 5.0 + 8.1*std::pow(Plab,-1.8 ) + 01324 0.16 *sqrLogPlab - 1.3 *LogPlab; 01325 01326 Xtotal = ( NumberOfTargetProtons * XtotKP + 01327 NumberOfTargetNeutrons * XtotKN ); 01328 01329 Xelastic = ( NumberOfTargetProtons * XelKP + 01330 NumberOfTargetNeutrons * XelKN ); 01331 } 01332 else if( PDGcode == 311 ) //------Projectile is KaonZero ------ 01333 { 01334 G4double XtotKP = ( 18.1 + 0. *std::pow(Plab, 0. ) + 01335 0.26 *sqrLogPlab - 1.0 *LogPlab + //K+ 01336 32.1 + 0. *std::pow(Plab, 0. ) + 01337 0.66 *sqrLogPlab - 5.6 *LogPlab)/2; //K- 01338 01339 G4double XtotKN = ( 18.7 + 0. *std::pow(Plab, 0. ) + 01340 0.21 *sqrLogPlab - 0.89*LogPlab + //K+ 01341 25.2 + 0. *std::pow(Plab, 0. ) + 01342 0.38 *sqrLogPlab - 2.9 *LogPlab)/2; //K- 01343 01344 G4double XelKP = ( 5.0 + 8.1*std::pow(Plab,-1.8 ) 01345 + 0.16 *sqrLogPlab - 1.3 *LogPlab + //K+ 01346 7.3 + 0. *std::pow(Plab,-0. ) + 01347 0.29 *sqrLogPlab - 2.4 *LogPlab)/2; //K- 01348 01349 G4double XelKN = ( 7.3 + 0. *std::pow(Plab,-0. ) + 01350 0.29 *sqrLogPlab - 2.4 *LogPlab + //K+ 01351 5.0 + 8.1*std::pow(Plab,-1.8 ) + 01352 0.16 *sqrLogPlab - 1.3 *LogPlab)/2; //K- 01353 01354 Xtotal = ( NumberOfTargetProtons * XtotKP + 01355 NumberOfTargetNeutrons * XtotKN ); 01356 01357 Xelastic = ( NumberOfTargetProtons * XelKP + 01358 NumberOfTargetNeutrons * XelKN ); 01359 } 01360 else //------Projectile is undefined, Nucleon assumed 01361 { 01362 G4double XtotPP = 48.0 + 0. *std::pow(Plab, 0. ) + 01363 0.522*sqrLogPlab - 4.51*LogPlab; 01364 01365 G4double XtotPN = 47.3 + 0. *std::pow(Plab, 0. ) + 01366 0.513*sqrLogPlab - 4.27*LogPlab; 01367 01368 G4double XelPP = 11.9 + 26.9*std::pow(Plab,-1.21) + 01369 0.169*sqrLogPlab - 1.85*LogPlab; 01370 G4double XelPN = 11.9 + 26.9*std::pow(Plab,-1.21) + 01371 0.169*sqrLogPlab - 1.85*LogPlab; 01372 01373 Xtotal = ( NumberOfTargetProtons * XtotPP + 01374 NumberOfTargetNeutrons * XtotPN ); 01375 01376 Xelastic = ( NumberOfTargetProtons * XelPP + 01377 NumberOfTargetNeutrons * XelPN ); 01378 } 01379 Xinelastic = Xtotal - Xelastic; 01380 01381 if( Xinelastic < 0.) Xinelastic = 0.; 01382 01383 return Xinelastic*= millibarn; 01384 }
G4double G4GlauberGribovCrossSection::GetInelasticGlauberGribov | ( | const G4DynamicParticle * | , | |
G4int | Z, | |||
G4int | A | |||
) | [inline] |
Definition at line 188 of file G4GlauberGribovCrossSection.hh.
References GetIsoCrossSection().
Referenced by G4CrossSectionPairGG::BuildPhysicsTable(), G4BGGPionInelasticXS::BuildPhysicsTable(), G4BGGNucleonInelasticXS::BuildPhysicsTable(), G4CrossSectionPairGG::GetElementCrossSection(), G4BGGPionInelasticXS::GetElementCrossSection(), and G4BGGNucleonInelasticXS::GetElementCrossSection().
00190 { 00191 GetIsoCrossSection(dp, Z, A); 00192 return fInelasticXsc; 00193 }
G4double G4GlauberGribovCrossSection::GetInelasticGlauberGribovXsc | ( | ) | [inline] |
Definition at line 105 of file G4GlauberGribovCrossSection.hh.
Referenced by G4NeutronInelasticXS::GetElementCrossSection().
G4double G4GlauberGribovCrossSection::GetIsoCrossSection | ( | const G4DynamicParticle * | , | |
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.
References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscNS(), G4HadronNucleonXsc::GetHadronNucleonXscNS(), GetHNinelasticXsc(), G4HadronNucleonXsc::GetInelasticHadronNucleonXsc(), GetKaonNucleonXscVector(), GetNucleusRadius(), GetParticleBarCorIn(), GetParticleBarCorTot(), and G4INCL::Math::pi.
Referenced by GetElasticGlauberGribov(), G4NeutronInelasticXS::GetElementCrossSection(), G4NeutronElasticXS::GetElementCrossSection(), and GetInelasticGlauberGribov().
00324 { 00325 G4double xsection, sigma, cofInelastic, cofTotal, nucleusSquare, ratio; 00326 G4double hpInXsc(0.), hnInXsc(0.); 00327 G4double R = GetNucleusRadius(A); 00328 00329 G4int N = A - Z; // number of neutrons 00330 if (N < 0) N = 0; 00331 00332 const G4ParticleDefinition* theParticle = aParticle->GetDefinition(); 00333 00334 if( theParticle == theProton || 00335 theParticle == theNeutron || 00336 theParticle == thePiPlus || 00337 theParticle == thePiMinus ) 00338 { 00339 // sigma = GetHadronNucleonXscNS(aParticle, A, Z); 00340 00341 sigma = Z*hnXsc->GetHadronNucleonXscNS(aParticle, theProton); 00342 00343 hpInXsc = hnXsc->GetInelasticHadronNucleonXsc(); 00344 00345 sigma += N*hnXsc->GetHadronNucleonXscNS(aParticle, theNeutron); 00346 00347 hnInXsc = hnXsc->GetInelasticHadronNucleonXsc(); 00348 00349 cofInelastic = 2.4; 00350 cofTotal = 2.0; 00351 } 00352 else if( theParticle == theKPlus || 00353 theParticle == theKMinus || 00354 theParticle == theK0S || 00355 theParticle == theK0L ) 00356 { 00357 sigma = GetKaonNucleonXscVector(aParticle, A, Z); 00358 cofInelastic = 2.2; 00359 cofTotal = 2.0; 00360 R = 1.3*fermi; 00361 R *= std::pow(G4double(A), 0.3333); 00362 } 00363 else 00364 { 00365 sigma = GetHadronNucleonXscNS(aParticle, A, Z); 00366 cofInelastic = 2.2; 00367 cofTotal = 2.0; 00368 } 00369 // cofInelastic = 2.0; 00370 00371 if( A > 1 ) 00372 { 00373 nucleusSquare = cofTotal*pi*R*R; // basically 2piRR 00374 ratio = sigma/nucleusSquare; 00375 00376 xsection = nucleusSquare*std::log( 1. + ratio ); 00377 00378 xsection *= GetParticleBarCorTot(theParticle, Z); 00379 00380 fTotalXsc = xsection; 00381 00382 00383 00384 fInelasticXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic; 00385 00386 fInelasticXsc *= GetParticleBarCorIn(theParticle, Z); 00387 00388 fElasticXsc = fTotalXsc - fInelasticXsc; 00389 00390 if(fElasticXsc < 0.) fElasticXsc = 0.; 00391 00392 G4double difratio = ratio/(1.+ratio); 00393 00394 fDiffractionXsc = 0.5*nucleusSquare*( difratio - std::log( 1. + difratio ) ); 00395 00396 00397 // sigma = GetHNinelasticXsc(aParticle, A, Z); 00398 00399 sigma = Z*hpInXsc + N*hnInXsc; 00400 00401 ratio = sigma/nucleusSquare; 00402 00403 fProductionXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic; 00404 00405 if (fElasticXsc < 0.) fElasticXsc = 0.; 00406 } 00407 else // H 00408 { 00409 fTotalXsc = sigma; 00410 xsection = sigma; 00411 00412 if ( theParticle != theAProton ) 00413 { 00414 sigma = GetHNinelasticXsc(aParticle, A, Z); 00415 fInelasticXsc = sigma; 00416 fElasticXsc = fTotalXsc - fInelasticXsc; 00417 } 00418 else 00419 { 00420 fElasticXsc = fTotalXsc - fInelasticXsc; 00421 } 00422 if (fElasticXsc < 0.) fElasticXsc = 0.; 00423 00424 } 00425 return xsection; 00426 }
G4double G4GlauberGribovCrossSection::GetKaonNucleonXscVector | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 1099 of file G4GlauberGribovCrossSection.cc.
References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscPDG(), G4DynamicParticle::GetKineticEnergy(), G4HadronNucleonXsc::GetKmNeutronTotXscVector(), G4HadronNucleonXsc::GetKmProtonTotXscVector(), G4HadronNucleonXsc::GetKpNeutronTotXscVector(), and G4HadronNucleonXsc::GetKpProtonTotXscVector().
Referenced by GetIsoCrossSection().
01101 { 01102 G4double Tkin, logTkin, xsc, xscP, xscN; 01103 const G4ParticleDefinition* theParticle = aParticle->GetDefinition(); 01104 01105 G4int Nt = At-Zt; // number of neutrons 01106 if (Nt < 0) Nt = 0; 01107 01108 Tkin = aParticle->GetKineticEnergy(); // Tkin in MeV 01109 01110 if( Tkin > 70*GeV ) return GetHadronNucleonXscPDG(aParticle,At,Zt); 01111 01112 logTkin = std::log(Tkin); // Tkin in MeV!!! 01113 01114 if( theParticle == theKPlus ) 01115 { 01116 xscP = hnXsc->GetKpProtonTotXscVector(logTkin); 01117 xscN = hnXsc->GetKpNeutronTotXscVector(logTkin); 01118 } 01119 else if( theParticle == theKMinus ) 01120 { 01121 xscP = hnXsc->GetKmProtonTotXscVector(logTkin); 01122 xscN = hnXsc->GetKmNeutronTotXscVector(logTkin); 01123 } 01124 else // K-zero as half of K+ and K- 01125 { 01126 xscP = (hnXsc->GetKpProtonTotXscVector(logTkin)+hnXsc->GetKmProtonTotXscVector(logTkin))*0.5; 01127 xscN = (hnXsc->GetKpNeutronTotXscVector(logTkin)+hnXsc->GetKmNeutronTotXscVector(logTkin))*0.5; 01128 } 01129 xsc = xscP*Zt + xscN*Nt; 01130 return xsc; 01131 }
Definition at line 1448 of file G4GlauberGribovCrossSection.cc.
References G4INCL::Math::oneThird.
01449 { 01450 G4double oneThird = 1.0/3.0; 01451 G4double cubicrAt = std::pow(G4double(At), oneThird); 01452 01453 G4double R; // = fRadiusConst*cubicrAt; 01454 01455 /* 01456 G4double tmp = std::pow( cubicrAt-1., 3.); 01457 tmp += At; 01458 tmp *= 0.5; 01459 01460 if (At > 20.) 01461 { 01462 R = fRadiusConst*std::pow (tmp, oneThird); 01463 } 01464 else 01465 { 01466 R = fRadiusConst*cubicrAt; 01467 } 01468 */ 01469 01470 R = fRadiusConst*cubicrAt; 01471 01472 G4double meanA = 20.; 01473 G4double tauA = 20.; 01474 01475 if (At > 20) // 20. 01476 { 01477 R *= ( 0.8 + 0.2*std::exp( -(G4double(At) - meanA)/tauA) ); 01478 } 01479 else 01480 { 01481 R *= ( 1.0 + 0.1*( 1. - std::exp( (G4double(At) - meanA)/tauA) ) ); 01482 } 01483 01484 return R; 01485 }
G4double G4GlauberGribovCrossSection::GetNucleusRadius | ( | const G4DynamicParticle * | , | |
const G4Element * | ||||
) |
Definition at line 1391 of file G4GlauberGribovCrossSection.cc.
References G4lrint(), G4Element::GetN(), and G4INCL::Math::oneThird.
Referenced by GetIsoCrossSection(), GetRatioQE(), and GetRatioSD().
01393 { 01394 G4int At = G4lrint(anElement->GetN()); 01395 G4double oneThird = 1.0/3.0; 01396 G4double cubicrAt = std::pow(G4double(At), oneThird); 01397 01398 G4double R; // = fRadiusConst*cubicrAt; 01399 /* 01400 G4double tmp = std::pow( cubicrAt-1., 3.); 01401 tmp += At; 01402 tmp *= 0.5; 01403 01404 if (At > 20.) // 20. 01405 { 01406 R = fRadiusConst*std::pow (tmp, oneThird); 01407 } 01408 else 01409 { 01410 R = fRadiusConst*cubicrAt; 01411 } 01412 */ 01413 01414 R = fRadiusConst*cubicrAt; 01415 01416 G4double meanA = 21.; 01417 01418 G4double tauA1 = 40.; 01419 G4double tauA2 = 10.; 01420 G4double tauA3 = 5.; 01421 01422 G4double a1 = 0.85; 01423 G4double b1 = 1. - a1; 01424 01425 G4double b2 = 0.3; 01426 G4double b3 = 4.; 01427 01428 if (At > 20) // 20. 01429 { 01430 R *= ( a1 + b1*std::exp( -(At - meanA)/tauA1) ); 01431 } 01432 else if (At > 3) 01433 { 01434 R *= ( 1.0 + b2*( 1. - std::exp( (At - meanA)/tauA2) ) ); 01435 } 01436 else 01437 { 01438 R *= ( 1.0 + b3*( 1. - std::exp( (At - meanA)/tauA3) ) ); 01439 } 01440 return R; 01441 01442 }
G4double G4GlauberGribovCrossSection::GetParticleBarCorIn | ( | const G4ParticleDefinition * | theParticle, | |
G4int | Z | |||
) | [inline] |
Definition at line 221 of file G4GlauberGribovCrossSection.hh.
Referenced by GetIsoCrossSection().
00223 { 00224 if(Z >= 2 && Z <= 92) 00225 { 00226 if( theParticle == theProton ) return fProtonBarCorrectionIn[Z]; 00227 else if( theParticle == theNeutron) return fNeutronBarCorrectionIn[Z]; 00228 else if( theParticle == thePiPlus ) return fPionPlusBarCorrectionIn[Z]; 00229 else if( theParticle == thePiMinus) return fPionMinusBarCorrectionIn[Z]; 00230 else return 1.0; 00231 } 00232 else return 1.0; 00233 }
G4double G4GlauberGribovCrossSection::GetParticleBarCorTot | ( | const G4ParticleDefinition * | theParticle, | |
G4int | Z | |||
) | [inline] |
Definition at line 201 of file G4GlauberGribovCrossSection.hh.
Referenced by GetIsoCrossSection().
00203 { 00204 if(Z >= 2 && Z <= 92) 00205 { 00206 if( theParticle == theProton ) return fProtonBarCorrectionTot[Z]; 00207 else if( theParticle == theNeutron) return fNeutronBarCorrectionTot[Z]; 00208 else if( theParticle == thePiPlus ) return fPionPlusBarCorrectionTot[Z]; 00209 else if( theParticle == thePiMinus) return fPionMinusBarCorrectionTot[Z]; 00210 else return 1.0; 00211 } 00212 else return 1.0; 00213 }
G4double G4GlauberGribovCrossSection::GetProductionGlauberGribovXsc | ( | ) | [inline] |
G4double G4GlauberGribovCrossSection::GetRadiusConst | ( | ) | [inline] |
G4double G4GlauberGribovCrossSection::GetRatioQE | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 475 of file G4GlauberGribovCrossSection.cc.
References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscNS(), GetHNinelasticXsc(), GetNucleusRadius(), and G4INCL::Math::pi.
00476 { 00477 G4double sigma, cofInelastic, cofTotal, nucleusSquare, ratio; 00478 G4double R = GetNucleusRadius(A); 00479 00480 const G4ParticleDefinition* theParticle = aParticle->GetDefinition(); 00481 00482 if( theParticle == theProton || 00483 theParticle == theNeutron || 00484 theParticle == thePiPlus || 00485 theParticle == thePiMinus ) 00486 { 00487 sigma = GetHadronNucleonXscNS(aParticle, A, Z); 00488 cofInelastic = 2.4; 00489 cofTotal = 2.0; 00490 } 00491 else 00492 { 00493 sigma = GetHadronNucleonXscNS(aParticle, A, Z); 00494 cofInelastic = 2.2; 00495 cofTotal = 2.0; 00496 } 00497 nucleusSquare = cofTotal*pi*R*R; // basically 2piRR 00498 ratio = sigma/nucleusSquare; 00499 00500 fInelasticXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic; 00501 00502 sigma = GetHNinelasticXsc(aParticle, A, Z); 00503 ratio = sigma/nucleusSquare; 00504 00505 fProductionXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic; 00506 00507 if (fInelasticXsc > fProductionXsc) ratio = (fInelasticXsc-fProductionXsc)/fInelasticXsc; 00508 else ratio = 0.; 00509 if ( ratio < 0. ) ratio = 0.; 00510 00511 return ratio; 00512 }
G4double G4GlauberGribovCrossSection::GetRatioSD | ( | const G4DynamicParticle * | , | |
G4int | At, | |||
G4int | Zt | |||
) |
Definition at line 433 of file G4GlauberGribovCrossSection.cc.
References G4DynamicParticle::GetDefinition(), GetHadronNucleonXscNS(), GetNucleusRadius(), and G4INCL::Math::pi.
00434 { 00435 G4double sigma, cofInelastic, cofTotal, nucleusSquare, ratio; 00436 G4double R = GetNucleusRadius(A); 00437 00438 const G4ParticleDefinition* theParticle = aParticle->GetDefinition(); 00439 00440 if( theParticle == theProton || 00441 theParticle == theNeutron || 00442 theParticle == thePiPlus || 00443 theParticle == thePiMinus ) 00444 { 00445 sigma = GetHadronNucleonXscNS(aParticle, A, Z); 00446 cofInelastic = 2.4; 00447 cofTotal = 2.0; 00448 } 00449 else 00450 { 00451 sigma = GetHadronNucleonXscNS(aParticle, A, Z); 00452 cofInelastic = 2.2; 00453 cofTotal = 2.0; 00454 } 00455 nucleusSquare = cofTotal*pi*R*R; // basically 2piRR 00456 ratio = sigma/nucleusSquare; 00457 00458 fInelasticXsc = nucleusSquare*std::log( 1. + cofInelastic*ratio )/cofInelastic; 00459 00460 G4double difratio = ratio/(1.+ratio); 00461 00462 fDiffractionXsc = 0.5*nucleusSquare*( difratio - std::log( 1. + difratio ) ); 00463 00464 if (fInelasticXsc > 0.) ratio = fDiffractionXsc/fInelasticXsc; 00465 else ratio = 0.; 00466 00467 return ratio; 00468 }
G4double G4GlauberGribovCrossSection::GetTotalGlauberGribovXsc | ( | ) | [inline] |
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.
References G4DynamicParticle::GetDefinition(), and G4DynamicParticle::GetKineticEnergy().
00287 { 00288 G4bool applicable = false; 00289 // G4int baryonNumber = aDP->GetDefinition()->GetBaryonNumber(); 00290 G4double kineticEnergy = aDP->GetKineticEnergy(); 00291 00292 const G4ParticleDefinition* theParticle = aDP->GetDefinition(); 00293 00294 if ( ( kineticEnergy >= fLowerLimit && 00295 Z > 1 && // >= He 00296 ( theParticle == theAProton || 00297 theParticle == theGamma || 00298 theParticle == theKPlus || 00299 theParticle == theKMinus || 00300 theParticle == theK0L || 00301 theParticle == theK0S || 00302 theParticle == theSMinus || 00303 theParticle == theProton || 00304 theParticle == theNeutron || 00305 theParticle == thePiPlus || 00306 theParticle == thePiMinus ) ) ) applicable = true; 00307 00308 return applicable; 00309 }
void G4GlauberGribovCrossSection::SetEnergyLowerLimit | ( | G4double | E | ) | [inline] |