diff --git a/source/processes/hadronic/cross_sections/include/G4ComponentAntiNuclNuclearXS.hh b/source/processes/hadronic/cross_sections/include/G4ComponentAntiNuclNuclearXS.hh index 596add3dcd11ff9c124c8b4e11926ad121da8a90..4686628295c995afbe892b8648f902cb92ef382c 100644 --- a/source/processes/hadronic/cross_sections/include/G4ComponentAntiNuclNuclearXS.hh +++ b/source/processes/hadronic/cross_sections/include/G4ComponentAntiNuclNuclearXS.hh @@ -52,84 +52,52 @@ #include "G4VComponentCrossSection.hh" -class G4ParticleDefinition; - -class G4ComponentAntiNuclNuclearXS : public G4VComponentCrossSection -{ -public: - - G4ComponentAntiNuclNuclearXS (); - virtual ~G4ComponentAntiNuclNuclearXS (); - -virtual - G4double GetTotalIsotopeCrossSection(const G4ParticleDefinition* aParticle, - G4double kinEnergy, - G4int Z, G4int A); - -virtual - G4double GetTotalElementCrossSection(const G4ParticleDefinition* aParticle, - G4double kinEnergy, - G4int Z, G4double A); - -virtual - G4double GetInelasticIsotopeCrossSection(const G4ParticleDefinition* aParticle, - G4double kinEnergy, - G4int Z, G4int A); - -virtual - G4double GetInelasticElementCrossSection(const G4ParticleDefinition* aParticle, - G4double kinEnergy, - G4int Z, G4double A); - -virtual - G4double GetElasticElementCrossSection(const G4ParticleDefinition* aParticle, - G4double kinEnergy, - G4int Z, G4double A); - -virtual - G4double GetElasticIsotopeCrossSection(const G4ParticleDefinition* aParticle, - G4double kinEnergy, - G4int Z, G4int A); - -virtual - void BuildPhysicsTable(const G4ParticleDefinition&) - {} -virtual - void DumpPhysicsTable(const G4ParticleDefinition&) - {} - -virtual void CrossSectionDescription(std::ostream&) const; - -// Method for calculation of Anti-Hadron Nucleon Total Cross-section -G4double GetAntiHadronNucleonTotCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy); - - -// Method for calculation of Anti-Hadron Nucleon Elastic Cross-section -G4double GetAntiHadronNucleonElCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy); - -// Method for calculation of scaling factor for cross-section -G4double GetScalingFactorCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy); - -private: +class G4ParticleDefinition; -// const G4double fUpperLimit; -// const G4double fLowerLimit; - G4double fRadiusEff; // Effective Radius for AntiNucleus - G4double fRadiusNN2; // Sqr of radius of NN collision - G4double fTotalXsc, fElasticXsc, fInelasticXsc; - G4double fAntiHadronNucleonTotXsc, fAntiHadronNucleonElXsc; - G4double Elab, S, SqrtS ; - G4double Mn, b0, b2, SqrtS0, S0, R0; //parameters for AntiHadron-Nucleon Xsc - G4double fScalingFactorXsc; // scaling factor for cross-sections - - G4ParticleDefinition* theAProton; - G4ParticleDefinition* theANeutron; - G4ParticleDefinition* theADeuteron; - G4ParticleDefinition* theATriton; - G4ParticleDefinition* theAAlpha; - G4ParticleDefinition* theAHe3; +class G4ComponentAntiNuclNuclearXS : public G4VComponentCrossSection { + + public: + G4ComponentAntiNuclNuclearXS (); + virtual ~G4ComponentAntiNuclNuclearXS (); + virtual G4double GetTotalIsotopeCrossSection(const G4ParticleDefinition* aParticle, + G4double kinEnergy, G4int Z, G4int A); + virtual G4double GetTotalElementCrossSection(const G4ParticleDefinition* aParticle, + G4double kinEnergy, G4int Z, G4double A); + virtual G4double GetInelasticIsotopeCrossSection(const G4ParticleDefinition* aParticle, + G4double kinEnergy, G4int Z, G4int A); + virtual G4double GetInelasticElementCrossSection(const G4ParticleDefinition* aParticle, + G4double kinEnergy, G4int Z, G4double A); + virtual G4double GetElasticElementCrossSection(const G4ParticleDefinition* aParticle, + G4double kinEnergy, G4int Z, G4double A); + virtual G4double GetElasticIsotopeCrossSection(const G4ParticleDefinition* aParticle, + G4double kinEnergy, G4int Z, G4int A); + virtual void BuildPhysicsTable(const G4ParticleDefinition&) {} + virtual void DumpPhysicsTable(const G4ParticleDefinition&) {} + virtual void CrossSectionDescription(std::ostream&) const; + // Method for calculation of Anti-Hadron Nucleon Total Cross-section + G4double GetAntiHadronNucleonTotCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy); + // Method for calculation of Anti-Hadron Nucleon Elastic Cross-section + G4double GetAntiHadronNucleonElCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy); + // Method for calculation of scaling factor for cross-section + G4double GetScalingFactorCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy); + + private: + G4double fRadiusEff; // Effective Radius for AntiNucleus + G4double fRadiusNN2; // Sqr of radius of NN collision + G4double fTotalXsc, fElasticXsc, fInelasticXsc; + G4double fAntiHadronNucleonTotXsc, fAntiHadronNucleonElXsc; + G4double Elab, S, SqrtS ; + G4double Mn, b0, b2, SqrtS0, S0, R0; // Parameters for AntiHadron-Nucleon Xsc + G4double fScalingFactorXsc; // scaling factor for cross-sections + + G4ParticleDefinition* theAProton; + G4ParticleDefinition* theANeutron; + G4ParticleDefinition* theADeuteron; + G4ParticleDefinition* theATriton; + G4ParticleDefinition* theAAlpha; + G4ParticleDefinition* theAHe3; }; diff --git a/source/processes/hadronic/cross_sections/src/G4ComponentAntiNuclNuclearXS.cc b/source/processes/hadronic/cross_sections/src/G4ComponentAntiNuclNuclearXS.cc index be774aa8d7efa87f392fd71d1181f0c5219112ac..f448dd10a4f881a941fc81f795dab55676716302 100644 --- a/source/processes/hadronic/cross_sections/src/G4ComponentAntiNuclNuclearXS.cc +++ b/source/processes/hadronic/cross_sections/src/G4ComponentAntiNuclNuclearXS.cc @@ -45,266 +45,228 @@ #include <dlfcn.h> -/////////////////////////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////////////////////// G4ComponentAntiNuclNuclearXS::G4ComponentAntiNuclNuclearXS() : G4VComponentCrossSection("AntiAGlauber"), -// fUpperLimit(10000*GeV), fLowerLimit(10*MeV), fRadiusEff(0.0), fRadiusNN2(0.0), fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0), fAntiHadronNucleonTotXsc(0.0), fAntiHadronNucleonElXsc(0.0), Elab(0.0), S(0.0), SqrtS(0), fScalingFactorXsc(1.0) { - theAProton = G4AntiProton::AntiProton(); - theANeutron = G4AntiNeutron::AntiNeutron(); - theADeuteron = G4AntiDeuteron::AntiDeuteron(); - theATriton = G4AntiTriton::AntiTriton(); - theAAlpha = G4AntiAlpha::AntiAlpha(); - theAHe3 = G4AntiHe3::AntiHe3(); - - Mn = 0.93827231; // GeV - b0 = 11.92; // GeV^(-2) - b2 = 0.3036; // GeV^(-2) - SqrtS0 = 20.74; // GeV - S0 = 33.0625; // GeV^2 - R0 = 1.0; // default value (V.Ivanchenko) + theAProton = G4AntiProton::AntiProton(); + theANeutron = G4AntiNeutron::AntiNeutron(); + theADeuteron = G4AntiDeuteron::AntiDeuteron(); + theATriton = G4AntiTriton::AntiTriton(); + theAAlpha = G4AntiAlpha::AntiAlpha(); + theAHe3 = G4AntiHe3::AntiHe3(); + Mn = 0.93827231; // GeV + b0 = 11.92; // GeV^(-2) + b2 = 0.3036; // GeV^(-2) + SqrtS0 = 20.74; // GeV + S0 = 33.0625; // GeV^2 + R0 = 1.0; // default value (V.Ivanchenko) } -/////////////////////////////////////////////////////////////////////////////////////// -// -// + +///////////////////////////////////////////////////////////////////////////// G4ComponentAntiNuclNuclearXS::~G4ComponentAntiNuclNuclearXS() { } -//////////////////////////////////////////////////////////////////////////////// - - -//////////////////////////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////////////////////// // // Calculation of total CrossSection of Anti-Nucleus - Nucleus - G4double G4ComponentAntiNuclNuclearXS::GetTotalElementCrossSection (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A) { - G4double xsection, sigmaTotal, sigmaElastic; - - const G4ParticleDefinition* theParticle = aParticle; - - sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy); - sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy); - -// calculation of squared radius of NN-collision - fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi) ; //fm^2 + G4double xsection, sigmaTotal, sigmaElastic; + const G4ParticleDefinition* theParticle = aParticle; + sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy); + sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy); -// calculation of effective nuclear radius for Pbar and Nbar interactions (can be changed) + // calculation of squared radius of NN-collision + fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi) ; //fm^2 + // calculation of effective nuclear radius for Pbar and Nbar interactions (can be changed) //A.R. 29-Jan-2013 : use antiprotons/antineutrons as the default case, // to be used for instance, as first approximation // without validation, for anti-hyperons. - //if ( (theParticle == theAProton) || (theParticle == theANeutron) ) { - if(A==1) - { fTotalXsc = sigmaTotal * millibarn; - return fTotalXsc; } - - fRadiusEff = 1.34*G4Pow::GetInstance()->powA(A,0.23)+1.35/G4Pow::GetInstance()->powA(A,1./3.); //fm - - if( (Z==1) && (A==2) ) fRadiusEff = 3.800; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 3.300; - if( (Z==2) && (A==3) ) fRadiusEff = 3.300; - if( (Z==2) && (A==4) ) fRadiusEff = 2.376; - //} - -//calculation of effective nuclear radius for AntiDeuteron interaction (can be changed) - if (theParticle == theADeuteron) - { fRadiusEff = 1.46 * G4Pow::GetInstance()->powA(A,0.21) + 1.45 / G4Pow::GetInstance()->powA(A,1./3.); - - if( (Z==1) && (A==2) ) fRadiusEff = 3.238; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 3.144; - if( (Z==2) && (A==3) ) fRadiusEff = 3.144; - if( (Z==2) && (A==4) ) fRadiusEff = 2.544; - } -// calculation of effective nuclear radius for AntiHe3 interaction (can be changed) - - if( (theParticle ==theAHe3) || (theParticle ==theATriton) ) - { fRadiusEff = 1.40* G4Pow::GetInstance()->powA(A,0.21)+1.63/G4Pow::GetInstance()->powA(A,1./3.); - - if( (Z==1) && (A==2) ) fRadiusEff = 3.144; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 3.075; - if( (Z==2) && (A==3) ) fRadiusEff = 3.075; - if( (Z==2) && (A==4) ) fRadiusEff = 2.589; + if (A==1) { + fTotalXsc = sigmaTotal * millibarn; + return fTotalXsc; + } + fRadiusEff = 1.34*G4Pow::GetInstance()->powA(A,0.23)+1.35/G4Pow::GetInstance()->powA(A,1./3.); //fm + if ( (Z==1) && (A==2) ) fRadiusEff = 3.800; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 3.300; + if ( (Z==2) && (A==3) ) fRadiusEff = 3.300; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.376; + + // calculation of effective nuclear radius for AntiDeuteron interaction (can be changed) + if (theParticle == theADeuteron) { + fRadiusEff = 1.46 * G4Pow::GetInstance()->powA(A,0.21) + 1.45 / G4Pow::GetInstance()->powA(A,1./3.); + if ( (Z==1) && (A==2) ) fRadiusEff = 3.238; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 3.144; + if ( (Z==2) && (A==3) ) fRadiusEff = 3.144; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.544; } -//calculation of effective nuclear radius for AntiAlpha interaction (can be changed) + // calculation of effective nuclear radius for AntiHe3 interaction (can be changed) + if ( (theParticle ==theAHe3) || (theParticle ==theATriton) ) { + fRadiusEff = 1.40* G4Pow::GetInstance()->powA(A,0.21)+1.63/G4Pow::GetInstance()->powA(A,1./3.); + if ( (Z==1) && (A==2) ) fRadiusEff = 3.144; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 3.075; + if ( (Z==2) && (A==3) ) fRadiusEff = 3.075; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.589; + } - if (theParticle == theAAlpha) - { - fRadiusEff = 1.35* G4Pow::GetInstance()->powA(A,0.21)+1.1/G4Pow::GetInstance()->powA(A,1./3.); - - if( (Z==1) && (A==2) ) fRadiusEff = 2.544; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 2.589; - if( (Z==2) && (A==3) ) fRadiusEff = 2.589; - if( (Z==2) && (A==4) ) fRadiusEff = 2.241; - - } + // calculation of effective nuclear radius for AntiAlpha interaction (can be changed) + if (theParticle == theAAlpha) { + fRadiusEff = 1.35* G4Pow::GetInstance()->powA(A,0.21)+1.1/G4Pow::GetInstance()->powA(A,1./3.); + if ( (Z==1) && (A==2) ) fRadiusEff = 2.544; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 2.589; + if ( (Z==2) && (A==3) ) fRadiusEff = 2.589; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.241; + } - G4double R2 = fRadiusEff*fRadiusEff; - G4double REf2 = R2+fRadiusNN2; - G4double ApAt = std::abs(theParticle->GetBaryonNumber()) * A; + G4double R2 = fRadiusEff*fRadiusEff; + G4double REf2 = R2+fRadiusNN2; + G4double ApAt = std::abs(theParticle->GetBaryonNumber()) * A; - xsection = 2*pi*REf2*10.*G4Log(1+(ApAt*sigmaTotal/(2*pi*REf2*10.))); //mb - xsection =xsection *millibarn; - fTotalXsc = xsection; + xsection = 2*pi*REf2*10.*G4Log(1+(ApAt*sigmaTotal/(2*pi*REf2*10.))); //mb + xsection = xsection *millibarn; + fTotalXsc = xsection; - return fTotalXsc; + return fTotalXsc; } -//////////////////////////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////////////////////// // // Calculation of total CrossSection of Anti-Nucleus - Nucleus -////////////////////////////////////////////////////////////////////////////// + G4double G4ComponentAntiNuclNuclearXS::GetTotalIsotopeCrossSection (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A ) -{ return GetTotalElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); } +{ + return GetTotalElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); +} + -//////////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////////////////////// // Calculation of inelastic CrossSection of Anti-Nucleus - Nucleus -//////////////////////////////////////////////////////////////// G4double G4ComponentAntiNuclNuclearXS::GetInelasticElementCrossSection (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A) { G4double inelxsection, sigmaTotal, sigmaElastic; - const G4ParticleDefinition* theParticle = aParticle; - - sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy); - sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy); + sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy); + sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy); -// calculation of sqr of radius NN-collision - fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi); // fm^2 - - -// calculation of effective nuclear radius for Pbar and Nbar interaction (can be changed) + // calculation of sqr of radius NN-collision + fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi); // fm^2 + // calculation of effective nuclear radius for Pbar and Nbar interaction (can be changed) //A.R. 29-Jan-2013 : use antiprotons/antineutrons as the default case, // to be used for instance, as first approximation // without validation, for anti-hyperons. - //if ( (theParticle == theAProton) || (theParticle == theANeutron) ) { - if (A==1) - { fInelasticXsc = (sigmaTotal - sigmaElastic)*millibarn; - return fInelasticXsc; - } - - - fRadiusEff = 1.31*G4Pow::GetInstance()->powA(A, 0.22)+0.9/G4Pow::GetInstance()->powA(A, 1./3.); //fm - - if( (Z==1) && (A==2) ) fRadiusEff = 3.582; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 3.105; - if( (Z==2) && (A==3) ) fRadiusEff = 3.105; - if( (Z==2) && (A==4) ) fRadiusEff = 2.209; - //} - -//calculation of effective nuclear radius for AntiDeuteron interaction (can be changed) - - if (theParticle ==theADeuteron) -{ - fRadiusEff = 1.38*G4Pow::GetInstance()->powA(A, 0.21)+1.55/G4Pow::GetInstance()->powA(A, 1./3.); - - if( (Z==1) && (A==2) ) fRadiusEff = 3.169; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 3.066; - if( (Z==2) && (A==3) ) fRadiusEff = 3.066; - if( (Z==2) && (A==4) ) fRadiusEff = 2.498; - } + if (A==1) { + fInelasticXsc = (sigmaTotal - sigmaElastic) * millibarn; + return fInelasticXsc; + } + fRadiusEff = 1.31*G4Pow::GetInstance()->powA(A, 0.22)+0.9/G4Pow::GetInstance()->powA(A, 1./3.); //fm + if ( (Z==1) && (A==2) ) fRadiusEff = 3.582; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 3.105; + if ( (Z==2) && (A==3) ) fRadiusEff = 3.105; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.209; + + // calculation of effective nuclear radius for AntiDeuteron interaction (can be changed) + if (theParticle ==theADeuteron) { + fRadiusEff = 1.38*G4Pow::GetInstance()->powA(A, 0.21)+1.55/G4Pow::GetInstance()->powA(A, 1./3.); + if ( (Z==1) && (A==2) ) fRadiusEff = 3.169; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 3.066; + if ( (Z==2) && (A==3) ) fRadiusEff = 3.066; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.498; + } -//calculation of effective nuclear radius for AntiHe3 interaction (can be changed) + // calculation of effective nuclear radius for AntiHe3 interaction (can be changed) + if ( (theParticle ==theAHe3) || (theParticle ==theATriton) ) { + fRadiusEff = 1.34 * G4Pow::GetInstance()->powA(A, 0.21)+1.51/G4Pow::GetInstance()->powA(A, 1./3.); + if ( (Z==1) && (A==2) ) fRadiusEff = 3.066; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 2.973; + if ( (Z==2) && (A==3) ) fRadiusEff = 2.973; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.508; + } - if( (theParticle ==theAHe3) || (theParticle ==theATriton) ) - { - fRadiusEff = 1.34 * G4Pow::GetInstance()->powA(A, 0.21)+1.51/G4Pow::GetInstance()->powA(A, 1./3.); - - if( (Z==1) && (A==2) ) fRadiusEff = 3.066; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 2.973; - if( (Z==2) && (A==3) ) fRadiusEff = 2.973; - if( (Z==2) && (A==4) ) fRadiusEff = 2.508; - - } + // calculation of effective nuclear radius for AntiAlpha interaction (can be changed) + if (theParticle == theAAlpha) { + fRadiusEff = 1.3*G4Pow::GetInstance()->powA(A, 0.21)+1.05/G4Pow::GetInstance()->powA(A, 1./3.); + if ( (Z==1) && (A==2) ) fRadiusEff = 2.498; //fm + if ( (Z==1) && (A==3) ) fRadiusEff = 2.508; + if ( (Z==2) && (A==3) ) fRadiusEff = 2.508; + if ( (Z==2) && (A==4) ) fRadiusEff = 2.158; + } -//calculation of effective nuclear radius for AntiAlpha interaction (can be changed) + G4double R2 = fRadiusEff*fRadiusEff; + G4double REf2 = R2+fRadiusNN2; + G4double ApAt = std::abs(theParticle->GetBaryonNumber()) * A; - if (theParticle == theAAlpha) - { - fRadiusEff = 1.3*G4Pow::GetInstance()->powA(A, 0.21)+1.05/G4Pow::GetInstance()->powA(A, 1./3.); - - if( (Z==1) && (A==2) ) fRadiusEff = 2.498; //fm - if( (Z==1) && (A==3) ) fRadiusEff = 2.508; - if( (Z==2) && (A==3) ) fRadiusEff = 2.508; - if( (Z==2) && (A==4) ) fRadiusEff = 2.158; - } - G4double R2 = fRadiusEff*fRadiusEff; - G4double REf2 = R2+fRadiusNN2; - G4double ApAt= std::abs(theParticle->GetBaryonNumber()) * A; - - inelxsection = pi*REf2 *10* G4Log(1+(ApAt*sigmaTotal/(pi*REf2*10.))); //mb - inelxsection = inelxsection * millibarn; - fInelasticXsc = inelxsection; - - // scale inelastic cross-section with a momentum-dependent factor - // if (theParticle == theADeuteron || theParticle == theAProton) { - fScalingFactorXsc = GetScalingFactorCrSc(aParticle,kinEnergy); - fInelasticXsc = inelxsection*fScalingFactorXsc; - // } + inelxsection = pi*REf2 *10* G4Log(1+(ApAt*sigmaTotal/(pi*REf2*10.))); //mb + inelxsection = inelxsection * millibarn; + fInelasticXsc = inelxsection; - return fInelasticXsc; + // scale inelastic cross-section with a momentum-dependent factor + fScalingFactorXsc = GetScalingFactorCrSc(aParticle,kinEnergy); + fInelasticXsc = inelxsection*fScalingFactorXsc; + + return fInelasticXsc; } -/////////////////////////////////////////////////////////////////////////////// + +///////////////////////////////////////////////////////////////////////////// // // Calculates Inelastic Anti-nucleus-Nucleus cross-section -// + G4double G4ComponentAntiNuclNuclearXS::GetInelasticIsotopeCrossSection (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A) -{return GetInelasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); } - +{ + return GetInelasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); +} -/////////////////////////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////////////////////// // // Calculates elastic Anti-nucleus-Nucleus cross-section as Total - Inelastic -// + G4double G4ComponentAntiNuclNuclearXS::GetElasticElementCrossSection (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A) { - fElasticXsc = GetTotalElementCrossSection(aParticle, kinEnergy, Z, A)- - GetInelasticElementCrossSection(aParticle, kinEnergy, Z, A); - - // keep elastic cross-section unchanged - // if (aParticle == theADeuteron || aParticle == theAProton) { - fScalingFactorXsc = GetScalingFactorCrSc(aParticle,kinEnergy); - - fElasticXsc = GetTotalElementCrossSection(aParticle, kinEnergy, Z, A)- - GetInelasticElementCrossSection(aParticle, kinEnergy, Z, A)/fScalingFactorXsc; - // } + // keep elastic cross-section unchanged + fScalingFactorXsc = GetScalingFactorCrSc(aParticle,kinEnergy); - if (fElasticXsc < 0.) fElasticXsc = 0.; + fElasticXsc = GetTotalElementCrossSection(aParticle, kinEnergy, Z, A)- + GetInelasticElementCrossSection(aParticle, kinEnergy, Z, A)/fScalingFactorXsc; - return fElasticXsc; + if (fElasticXsc < 0.) fElasticXsc = 0.; + return fElasticXsc; } + -/////////////////////////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////////////////////// // // Calculates elastic Anti-nucleus-Nucleus cross-section -// + G4double G4ComponentAntiNuclNuclearXS::GetElasticIsotopeCrossSection (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A) -{ return GetElasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); } +{ + return GetElasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); +} -/////////////////////////////////////////////////////////////////////////////////// +///////////////////////////////////////////////////////////////////////////// // Calculation of Antihadron - hadron Total Cross-section G4double G4ComponentAntiNuclNuclearXS::GetAntiHadronNucleonTotCrSc @@ -317,97 +279,61 @@ G4double G4ComponentAntiNuclNuclearXS::GetAntiHadronNucleonTotCrSc momentum=std::sqrt(Energy*Energy-Pmass*Pmass)/std::abs(theParticle->GetBaryonNumber()); G4double Plab = momentum / GeV; - G4double B, SigAss; - G4double C, d1, d2, d3 ; - - Elab = std::sqrt(Mn*Mn + Plab*Plab); // GeV - S = 2.*Mn*Mn + 2. *Mn*Elab; // GeV^2 - SqrtS = std::sqrt(S); // GeV - - B = b0+b2*G4Log(SqrtS/SqrtS0)*G4Log(SqrtS/SqrtS0); //GeV^(-2) - SigAss = 36.04 +0.304*G4Log(S/S0)*G4Log(S/S0); //mb - R0 = std::sqrt(0.40874044*SigAss - B); //GeV^(-2) - - C = 13.55; - d1 = -4.47; - d2 = 12.38; - d3 = -12.43; - xsection = SigAss*(1 + 1./(std::sqrt(S-4.*Mn*Mn)) / (G4Pow::GetInstance()->powA(R0, 3.)) - *C* (1+d1/SqrtS+d2/(G4Pow::GetInstance()->powA(SqrtS,2.))+d3/(G4Pow::GetInstance()->powA(SqrtS,3.)) )); - -// xsection *= millibarn; - + G4double B, SigAss; + G4double C, d1, d2, d3; + Elab = std::sqrt(Mn*Mn + Plab*Plab); // GeV + S = 2.*Mn*Mn + 2. *Mn*Elab; // GeV^2 + SqrtS = std::sqrt(S); // GeV + B = b0+b2*G4Log(SqrtS/SqrtS0)*G4Log(SqrtS/SqrtS0); //GeV^(-2) + SigAss = 36.04 +0.304*G4Log(S/S0)*G4Log(S/S0); //mb + R0 = std::sqrt(0.40874044*SigAss - B); //GeV^(-2) + C = 13.55; + d1 = -4.47; + d2 = 12.38; + d3 = -12.43; + + xsection = SigAss * ( 1 + 1./(std::sqrt(S-4.*Mn*Mn)) / (G4Pow::GetInstance()->powA(R0, 3.)) + * C * ( 1 + d1/SqrtS + d2/(G4Pow::GetInstance()->powA(SqrtS,2.)) + + d3/(G4Pow::GetInstance()->powA(SqrtS,3.)) ) ); + + //xsection *= millibarn; fAntiHadronNucleonTotXsc = xsection; + return fAntiHadronNucleonTotXsc; } -// -// ///////////////////////////////////////////////////////////////////////////////// +// ////////////////////////////////////////////////////////////////////////// // Calculation of Antihadron - hadron Elastic Cross-section G4double G4ComponentAntiNuclNuclearXS :: GetAntiHadronNucleonElCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy) { - G4double xsection; - - G4double SigAss; - G4double C, d1, d2, d3 ; - - GetAntiHadronNucleonTotCrSc(aParticle,kinEnergy); - - SigAss = 4.5 + 0.101*G4Log(S/S0)*G4Log(S/S0); //mb - - C = 59.27; - d1 = -6.95; - d2 = 23.54; - d3 = -25.34; - - xsection = SigAss* (1 + 1. / (std::sqrt(S-4.*Mn*Mn)) / (G4Pow::GetInstance()->powA(R0, 3.)) - *C* ( 1+d1/SqrtS+d2/(G4Pow::GetInstance()->powA(SqrtS,2.))+d3/(G4Pow::GetInstance()->powA(SqrtS,3.)) )); - -// xsection *= millibarn; - + G4double xsection; + G4double SigAss; + G4double C, d1, d2, d3; + GetAntiHadronNucleonTotCrSc(aParticle,kinEnergy); + SigAss = 4.5 + 0.101*G4Log(S/S0)*G4Log(S/S0); //mb + C = 59.27; + d1 = -6.95; + d2 = 23.54; + d3 = -25.34; + + xsection = SigAss * ( 1 + 1. / (std::sqrt(S-4.*Mn*Mn)) / (G4Pow::GetInstance()->powA(R0, 3.)) + * C * ( 1 + d1/SqrtS + d2/(G4Pow::GetInstance()->powA(SqrtS,2.)) + + d3/(G4Pow::GetInstance()->powA(SqrtS,3.)) ) ); + + //xsection *= millibarn; fAntiHadronNucleonElXsc = xsection; + return fAntiHadronNucleonElXsc; } + +//////////////////////////////////////////////////////////////////////////// // Calculation of scaling factor for Cross-section G4double G4ComponentAntiNuclNuclearXS::GetScalingFactorCrSc(const G4ParticleDefinition* aParticle, G4double kinEnergy){ -// const G4ParticleDefinition* theParticle = aParticle; -// G4double mass = theParticle->GetPDGMass(); -// G4double p = std::sqrt(kinEnergy*kinEnergy + 2.0*kinEnergy*mass); // p in MeV -// p = p/1000.; // make p in GeV -// -// static void* lib = dlopen("./CustomScalingFunction.so", RTLD_NOW); -// static void* fun_ad = dlsym(lib, "CustomScalingFunctionAntiDeuteron"); -// static void* fun_ap = dlsym(lib, "CustomScalingFunctionAntiProton"); -// -// if (!lib){ -// return 1.0; -// } -// -// if ((theParticle == theADeuteron && !fun_ad) || -// (theParticle == theAProton && !fun_ap)){ -// return 1.0; -// } -// -// static double (*scaling_ap)(double) = (double (*)(double)) fun_ap; -// static double (*scaling_ad)(double) = (double (*)(double)) fun_ad; -// -// if (theParticle == theADeuteron){ -// return scaling_ad(p); -// } -// -// if (theParticle == theAProton){ -// return scaling_ap(p); -// } -// -// return 1.0; - - // generic implementation from Max - const G4ParticleDefinition* theParticle = aParticle; G4double mass = theParticle->GetPDGMass(); G4double p = std::sqrt(kinEnergy*kinEnergy + 2.0*kinEnergy*mass); // p in MeV @@ -429,13 +355,15 @@ G4double G4ComponentAntiNuclNuclearXS::GetScalingFactorCrSc(const G4ParticleDefi } +///////////////////////////////////////////////////////////////////////////// + void G4ComponentAntiNuclNuclearXS::CrossSectionDescription(std::ostream& outFile) const { outFile << "The G4ComponentAntiNuclNuclearXS calculates total,\n" << "inelastic, elastic cross sections of anti-nucleons and light \n" - << "anti-nucleus interactions with nuclei using Glauber's approach.\n" + << "anti-nucleus interactions with nuclei using Glauber's approach.\n" << "It uses parametrizations of antiproton-proton total and elastic \n" - << "cross sections and Wood-Saxon distribution of nuclear density.\n" - << "The lower limit is 10 MeV, the upper limit is 10 TeV. \n" + << "cross sections and Wood-Saxon distribution of nuclear density.\n" << "See details in Phys.Lett. B705 (2011) 235. \n"; } +