Geant4-11
Public Member Functions | Data Fields | Private Member Functions | Private Attributes
G4FTFParameters Class Reference

#include <G4FTFParameters.hh>

Public Member Functions

 G4FTFParameters ()
 
G4double GammaElastic (const G4double impactsquare)
 
G4double GetAvaragePt2ofElasticScattering ()
 
G4double GetAveragePt2 ()
 
G4double GetCofNuclearDestruction ()
 
G4double GetCofNuclearDestructionPr ()
 
G4double GetDeltaProbAtQuarkExchange ()
 
G4double GetDofNuclearDestruction ()
 
G4double GetElasticCrossSection ()
 
G4double GetExcitationEnergyPerWoundedNucleon ()
 
G4double GetInelasticCrossSection ()
 
G4double GetInelasticProbability (const G4double impactsquare)
 
G4double GetMaxNumberOfCollisions ()
 
G4double GetMaxPt2ofNuclearDestruction ()
 
G4double GetProbabilityOfAnnihilation ()
 
G4double GetProbabilityOfElasticScatt ()
 
G4double GetProbabilityOfInteraction (const G4double impactsquare)
 
G4double GetProbLogDistr ()
 
G4double GetProbLogDistrPrD ()
 
G4double GetProbOfInteraction ()
 
G4double GetProbOfSameQuarkExchange ()
 
G4double GetProcProb (const G4int ProcN, const G4double y)
 
G4double GetProjMinDiffMass ()
 
G4double GetProjMinNonDiffMass ()
 
G4double GetPt2Kink ()
 
G4double GetPt2ofNuclearDestruction ()
 
std::vector< G4doubleGetQuarkProbabilitiesAtGluonSplitUp ()
 
G4double GetR2ofNuclearDestruction ()
 
G4double GetSlope ()
 
G4double GetTarMinDiffMass ()
 
G4double GetTarMinNonDiffMass ()
 
G4double GetTotalCrossSection ()
 
void InitForInteraction (const G4ParticleDefinition *, G4int theA, G4int theZ, G4double s)
 
void SetAvaragePt2ofElasticScattering (const G4double aPt2)
 
void SetAveragePt2 (const G4double aValue)
 
void SetCofNuclearDestruction (const G4double aValue)
 
void SetCofNuclearDestructionPr (const G4double aValue)
 
void SetDeltaProbAtQuarkExchange (const G4double aValue)
 
void SetDofNuclearDestruction (const G4double aValue)
 
void SetElastisCrossSection (const G4double Xelastic)
 
void SetExcitationEnergyPerWoundedNucleon (const G4double aValue)
 
void SetGamma0 (const G4double Gamma0)
 
void SethNcmsEnergy (const G4double s)
 
void SetInelasticCrossSection (const G4double Xinelastic)
 
void SetMaxNumberOfCollisions (const G4double aValue, const G4double bValue)
 
void SetMaxPt2ofNuclearDestruction (const G4double aValue)
 
void SetParams (const G4int ProcN, const G4double A1, const G4double B1, const G4double A2, const G4double B2, const G4double A3, const G4double Atop, const G4double Ymin)
 
void SetProbabilityOfAnnihilation (const G4double aValue)
 
void SetProbabilityOfElasticScatt (const G4double aValue)
 
void SetProbabilityOfElasticScatt (const G4double Xtotal, const G4double Xelastic)
 
void SetProbLogDistr (const G4double aValue)
 
void SetProbLogDistrPrD (const G4double aValue)
 
void SetProbOfInteraction (const G4double aValue)
 
void SetProbOfSameQuarkExchange (const G4double aValue)
 
void SetProjMinDiffMass (const G4double aValue)
 
void SetProjMinNonDiffMass (const G4double aValue)
 
void SetPt2Kink (const G4double aValue)
 
void SetPt2ofNuclearDestruction (const G4double aValue)
 
void SetQuarkProbabilitiesAtGluonSplitUp (const G4double Puubar, const G4double Pddbar, const G4double Pssbar)
 
void SetR2ofNuclearDestruction (const G4double aValue)
 
void SetRadiusOfHNinteractions2 (const G4double Radius2)
 
void SetSlope (const G4double Slope)
 
void SetTarMinDiffMass (const G4double aValue)
 
void SetTarMinNonDiffMass (const G4double aValue)
 
void SetTotalCrossSection (const G4double Xtotal)
 
 ~G4FTFParameters ()
 

Data Fields

G4double AvaragePt2ofElasticScattering
 
G4double AveragePt2
 
G4double CofNuclearDestruction
 
G4double CofNuclearDestructionPr
 
G4double DeltaProbAtQuarkExchange
 
G4double DofNuclearDestruction
 
G4double ExcitationEnergyPerWoundedNucleon
 
G4double FTFGamma0
 
G4double FTFhNcmsEnergy
 
G4double FTFSlope
 
G4double FTFXannihilation
 
G4double FTFXelastic
 
G4double FTFXinelastic
 
G4double FTFXtotal
 
G4double MaxNumberOfCollisions
 
G4double MaxPt2ofNuclearDestruction
 
G4double ProbabilityOfAnnihilation
 
G4double ProbabilityOfElasticScatt
 
G4double ProbLogDistr
 
G4double ProbLogDistrPrD
 
G4double ProbOfInelInteraction
 
G4double ProbOfSameQuarkExchange
 
G4double ProcParams [5][7]
 
G4double ProjMinDiffMass
 
G4double ProjMinNonDiffMass
 
G4double Pt2kink
 
G4double Pt2ofNuclearDestruction
 
std::vector< G4doubleQuarkProbabilitiesAtGluonSplitUp
 
G4double R2ofNuclearDestruction
 
G4double RadiusOfHNinteractions2
 
G4double TarMinDiffMass
 
G4double TarMinNonDiffMass
 

Private Member Functions

G4double GetMinMass (const G4ParticleDefinition *aParticle)
 
void Reset ()
 

Private Attributes

G4VComponentCrossSectioncsGGinstance
 
G4FTFParamCollBaryonProj fParCollBaryonProj
 
G4FTFParamCollMesonProj fParCollMesonProj
 
G4FTFParamCollPionProj fParCollPionProj
 
G4LundStringFragmentationStringMass
 

Detailed Description

Definition at line 302 of file G4FTFParameters.hh.

Constructor & Destructor Documentation

◆ G4FTFParameters()

G4FTFParameters::G4FTFParameters ( )

Definition at line 62 of file G4FTFParameters.cc.

63{
64 StringMass = new G4LundStringFragmentation; // for estimation of min. mass of diffr. states
65 Reset();
68 if (!csGGinstance) {
70 }
71
72 // Set parameters of a string kink
73 SetPt2Kink( 0.0*GeV*GeV ); // To switch off kinky strings (bad results obtained with 6.0*GeV*GeV)
74 G4double Puubar( 1.0/3.0 ), Pddbar( 1.0/3.0 ), Pssbar( 1.0/3.0 ); // SU(3) symmetry
75 //G4double Puubar( 0.41 ), Pddbar( 0.41 ), Pssbar( 0.18 ); // Broken SU(3) symmetry
76 SetQuarkProbabilitiesAtGluonSplitUp( Puubar, Pddbar, Pssbar );
77}
static constexpr double GeV
Definition: G4SIunits.hh:203
double G4double
Definition: G4Types.hh:83
G4VComponentCrossSection * GetComponentCrossSection(const G4String &name)
static G4CrossSectionDataSetRegistry * Instance()
void SetQuarkProbabilitiesAtGluonSplitUp(const G4double Puubar, const G4double Pddbar, const G4double Pssbar)
void SetPt2Kink(const G4double aValue)
G4VComponentCrossSection * csGGinstance
G4LundStringFragmentation * StringMass

References csGGinstance, G4CrossSectionDataSetRegistry::GetComponentCrossSection(), GeV, G4CrossSectionDataSetRegistry::Instance(), Reset(), SetPt2Kink(), SetQuarkProbabilitiesAtGluonSplitUp(), and StringMass.

◆ ~G4FTFParameters()

G4FTFParameters::~G4FTFParameters ( )

Definition at line 1525 of file G4FTFParameters.cc.

1525 {
1526 if ( StringMass ) delete StringMass;
1527}

References StringMass.

Member Function Documentation

◆ GammaElastic()

G4double G4FTFParameters::GammaElastic ( const G4double  impactsquare)
inline

Definition at line 485 of file G4FTFParameters.hh.

485 {
486 return ( FTFGamma0 * G4Exp( -FTFSlope * impactsquare ) );
487}
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:179

References FTFGamma0, FTFSlope, and G4Exp().

Referenced by GetInelasticProbability().

◆ GetAvaragePt2ofElasticScattering()

G4double G4FTFParameters::GetAvaragePt2ofElasticScattering ( )
inline

◆ GetAveragePt2()

G4double G4FTFParameters::GetAveragePt2 ( )
inline

Definition at line 718 of file G4FTFParameters.hh.

718 {
719 return AveragePt2;
720}

References AveragePt2.

Referenced by G4DiffractiveExcitation::ExciteParticipants_doNonDiffraction(), and InitForInteraction().

◆ GetCofNuclearDestruction()

G4double G4FTFParameters::GetCofNuclearDestruction ( )
inline

Definition at line 754 of file G4FTFParameters.hh.

754 {
756}
G4double CofNuclearDestruction

References CofNuclearDestruction.

Referenced by InitForInteraction(), and G4FTFModel::ReggeonCascade().

◆ GetCofNuclearDestructionPr()

G4double G4FTFParameters::GetCofNuclearDestructionPr ( )
inline

Definition at line 750 of file G4FTFParameters.hh.

750 {
752}
G4double CofNuclearDestructionPr

References CofNuclearDestructionPr.

Referenced by InitForInteraction(), and G4FTFModel::ReggeonCascade().

◆ GetDeltaProbAtQuarkExchange()

G4double G4FTFParameters::GetDeltaProbAtQuarkExchange ( )
inline

◆ GetDofNuclearDestruction()

G4double G4FTFParameters::GetDofNuclearDestruction ( )
inline

◆ GetElasticCrossSection()

G4double G4FTFParameters::GetElasticCrossSection ( )
inline

Definition at line 654 of file G4FTFParameters.hh.

654 {
655 return FTFXelastic;
656}

References FTFXelastic.

◆ GetExcitationEnergyPerWoundedNucleon()

G4double G4FTFParameters::GetExcitationEnergyPerWoundedNucleon ( )
inline

◆ GetInelasticCrossSection()

G4double G4FTFParameters::GetInelasticCrossSection ( )
inline

Definition at line 658 of file G4FTFParameters.hh.

658 {
659 return FTFXinelastic;
660}

References FTFXinelastic.

◆ GetInelasticProbability()

G4double G4FTFParameters::GetInelasticProbability ( const G4double  impactsquare)
inline

Definition at line 678 of file G4FTFParameters.hh.

678 {
679 G4double Gamma = GammaElastic( impactsquare );
680 return 2*Gamma - Gamma*Gamma;
681}
G4double GammaElastic(const G4double impactsquare)

References GammaElastic().

◆ GetMaxNumberOfCollisions()

G4double G4FTFParameters::GetMaxNumberOfCollisions ( )
inline

Definition at line 742 of file G4FTFParameters.hh.

742 {
744}
G4double MaxNumberOfCollisions

References MaxNumberOfCollisions.

Referenced by G4FTFModel::ExciteParticipants().

◆ GetMaxPt2ofNuclearDestruction()

G4double G4FTFParameters::GetMaxPt2ofNuclearDestruction ( )
inline

◆ GetMinMass()

G4double G4FTFParameters::GetMinMass ( const G4ParticleDefinition aParticle)
private

Definition at line 824 of file G4FTFParameters.cc.

824 {
825 // The code is used for estimating the minimal string mass produced in diffraction dissociation.
826 // The indices used for minMassQDiQStr must be between 1 and 5, corresponding to the 5 considered
827 // quarks: d, u, s, c and b; enforcing this explicitly avoids compilation errors.
828 G4double EstimatedMass = 0.0;
829 G4int partID = std::abs(aParticle->GetPDGEncoding());
830 G4int Qleft = std::max( partID/100, 1 );
831 G4int Qright = std::max( (partID/ 10)%10, 1 );
832 if ( Qleft < 6 && Qright < 6 ) { // Q-Qbar string
833 EstimatedMass = StringMass->minMassQQbarStr[Qleft-1][Qright-1];
834 } else if ( Qleft < 6 && Qright > 6 ) { // Q - DiQ string
835 G4int q1 = std::max( std::min( Qright/10, 5 ), 1 );
836 G4int q2 = std::max( std::min( Qright%10, 5 ), 1 );
837 EstimatedMass = StringMass->minMassQDiQStr[Qleft-1][q1-1][q2-1];
838 } else if ( Qleft > 6 && Qright < 6 ) { // DiQ - Q string
839 G4int q1 = std::max( std::min( Qleft/10, 5 ), 1 );
840 G4int q2 = std::max( std::min( Qleft%10, 5 ), 1 );
841 EstimatedMass = StringMass->minMassQDiQStr[Qright-1][q1-1][q2-1];
842 }
843 return EstimatedMass;
844}
int G4int
Definition: G4Types.hh:85
T max(const T t1, const T t2)
brief Return the largest of the two arguments
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

References G4ParticleDefinition::GetPDGEncoding(), G4INCL::Math::max(), G4INCL::Math::min(), G4VLongitudinalStringDecay::minMassQDiQStr, G4VLongitudinalStringDecay::minMassQQbarStr, and StringMass.

Referenced by InitForInteraction().

◆ GetProbabilityOfAnnihilation()

G4double G4FTFParameters::GetProbabilityOfAnnihilation ( )
inline

Definition at line 683 of file G4FTFParameters.hh.

683 {
685}
G4double ProbabilityOfAnnihilation

References ProbabilityOfAnnihilation.

Referenced by G4FTFModel::ExciteParticipants().

◆ GetProbabilityOfElasticScatt()

G4double G4FTFParameters::GetProbabilityOfElasticScatt ( )
inline

Definition at line 674 of file G4FTFParameters.hh.

674 {
676}
G4double ProbabilityOfElasticScatt

References ProbabilityOfElasticScatt.

Referenced by G4FTFModel::ExciteParticipants().

◆ GetProbabilityOfInteraction()

G4double G4FTFParameters::GetProbabilityOfInteraction ( const G4double  impactsquare)
inline

Definition at line 666 of file G4FTFParameters.hh.

666 {
667 if ( RadiusOfHNinteractions2 > impactsquare ) {
668 return 1.0;
669 } else {
670 return 0.0;
671 }
672}
G4double RadiusOfHNinteractions2

References RadiusOfHNinteractions2.

Referenced by G4FTFParticipants::GetList().

◆ GetProbLogDistr()

G4double G4FTFParameters::GetProbLogDistr ( )
inline

Definition at line 726 of file G4FTFParameters.hh.

726 {
727 return ProbLogDistr;
728}

References ProbLogDistr.

Referenced by G4DiffractiveExcitation::ExciteParticipants_doNonDiffraction(), and InitForInteraction().

◆ GetProbLogDistrPrD()

G4double G4FTFParameters::GetProbLogDistrPrD ( )
inline

Definition at line 722 of file G4FTFParameters.hh.

722 {
723 return ProbLogDistrPrD;
724}
G4double ProbLogDistrPrD

References ProbLogDistrPrD.

Referenced by G4DiffractiveExcitation::ExciteParticipants_doNonDiffraction(), and InitForInteraction().

◆ GetProbOfInteraction()

G4double G4FTFParameters::GetProbOfInteraction ( )
inline

Definition at line 746 of file G4FTFParameters.hh.

746 {
748}
G4double ProbOfInelInteraction

References ProbOfInelInteraction.

◆ GetProbOfSameQuarkExchange()

G4double G4FTFParameters::GetProbOfSameQuarkExchange ( )
inline

◆ GetProcProb()

G4double G4FTFParameters::GetProcProb ( const G4int  ProcN,
const G4double  y 
)

Definition at line 848 of file G4FTFParameters.cc.

848 {
849 G4double Prob( 0.0 );
850 if ( y < ProcParams[ProcN][6] ) {
851 Prob = ProcParams[ProcN][5];
852 if (Prob < 0.) Prob=0.;
853 return Prob;
854 }
855 Prob = ProcParams[ProcN][0] * G4Exp( -ProcParams[ProcN][1]*y ) +
856 ProcParams[ProcN][2] * G4Exp( -ProcParams[ProcN][3]*y ) +
857 ProcParams[ProcN][4];
858 if (Prob < 0.) Prob=0.;
859 return Prob;
860}
G4double ProcParams[5][7]

References G4Exp(), and ProcParams.

Referenced by G4DiffractiveExcitation::ExciteParticipants().

◆ GetProjMinDiffMass()

G4double G4FTFParameters::GetProjMinDiffMass ( )
inline

◆ GetProjMinNonDiffMass()

G4double G4FTFParameters::GetProjMinNonDiffMass ( )
inline

Definition at line 706 of file G4FTFParameters.hh.

706 {
707 return ProjMinNonDiffMass;
708}
G4double ProjMinNonDiffMass

References ProjMinNonDiffMass.

Referenced by G4DiffractiveExcitation::ExciteParticipants(), and InitForInteraction().

◆ GetPt2Kink()

G4double G4FTFParameters::GetPt2Kink ( )
inline

Definition at line 732 of file G4FTFParameters.hh.

732 {
733 return Pt2kink;
734}

References Pt2kink.

Referenced by G4DiffractiveExcitation::CreateStrings().

◆ GetPt2ofNuclearDestruction()

G4double G4FTFParameters::GetPt2ofNuclearDestruction ( )
inline

◆ GetQuarkProbabilitiesAtGluonSplitUp()

std::vector< G4double > G4FTFParameters::GetQuarkProbabilitiesAtGluonSplitUp ( )
inline

Definition at line 736 of file G4FTFParameters.hh.

736 {
738}
std::vector< G4double > QuarkProbabilitiesAtGluonSplitUp

References QuarkProbabilitiesAtGluonSplitUp.

Referenced by G4DiffractiveExcitation::CreateStrings().

◆ GetR2ofNuclearDestruction()

G4double G4FTFParameters::GetR2ofNuclearDestruction ( )
inline

Definition at line 758 of file G4FTFParameters.hh.

758 {
760}
G4double R2ofNuclearDestruction

References R2ofNuclearDestruction.

Referenced by InitForInteraction(), and G4FTFModel::ReggeonCascade().

◆ GetSlope()

G4double G4FTFParameters::GetSlope ( )
inline

Definition at line 662 of file G4FTFParameters.hh.

662 {
663 return FTFSlope;
664}

References FTFSlope.

Referenced by InitForInteraction().

◆ GetTarMinDiffMass()

G4double G4FTFParameters::GetTarMinDiffMass ( )
inline

◆ GetTarMinNonDiffMass()

G4double G4FTFParameters::GetTarMinNonDiffMass ( )
inline

Definition at line 714 of file G4FTFParameters.hh.

714 {
715 return TarMinNonDiffMass;
716}
G4double TarMinNonDiffMass

References TarMinNonDiffMass.

Referenced by G4DiffractiveExcitation::ExciteParticipants(), and InitForInteraction().

◆ GetTotalCrossSection()

G4double G4FTFParameters::GetTotalCrossSection ( )
inline

Definition at line 650 of file G4FTFParameters.hh.

650 {
651 return FTFXtotal;
652}

References FTFXtotal.

◆ InitForInteraction()

void G4FTFParameters::InitForInteraction ( const G4ParticleDefinition particle,
G4int  theA,
G4int  theZ,
G4double  s 
)

Definition at line 81 of file G4FTFParameters.cc.

83{
84 Reset();
85
86 G4int ProjectilePDGcode = particle->GetPDGEncoding();
87 G4int ProjectileabsPDGcode = std::abs( ProjectilePDGcode );
88 G4double ProjectileMass = particle->GetPDGMass();
89 G4double ProjectileMass2 = ProjectileMass * ProjectileMass;
90
91 G4int ProjectileBaryonNumber( 0 ), AbsProjectileBaryonNumber( 0 ), AbsProjectileCharge( 0 );
92 G4bool ProjectileIsNucleus = false;
93
94 if ( std::abs( particle->GetBaryonNumber() ) > 1 ) { // The projectile is a nucleus
95 ProjectileIsNucleus = true;
96 ProjectileBaryonNumber = particle->GetBaryonNumber();
97 AbsProjectileBaryonNumber = std::abs( ProjectileBaryonNumber );
98 AbsProjectileCharge = std::abs( G4int( particle->GetPDGCharge() ) );
99 if ( ProjectileBaryonNumber > 1 ) {
100 ProjectilePDGcode = 2212; ProjectileabsPDGcode = 2212; // Proton
101 } else {
102 ProjectilePDGcode = -2212; ProjectileabsPDGcode = 2212; // Anti-Proton
103 }
104 ProjectileMass = G4Proton::Proton()->GetPDGMass();
105 ProjectileMass2 = sqr( ProjectileMass );
106 }
107
108 G4double TargetMass = G4Proton::Proton()->GetPDGMass();
109 G4double TargetMass2 = TargetMass * TargetMass;
110
111 G4double Plab = PlabPerParticle;
112 G4double Elab = std::sqrt( Plab*Plab + ProjectileMass2 );
113 G4double KineticEnergy = Elab - ProjectileMass;
114
115 G4double S = ProjectileMass2 + TargetMass2 + 2.0*TargetMass*Elab;
116
117 #ifdef debugFTFparams
118 G4cout << "--------- FTF Parameters --------------" << G4endl << "Proj Plab "
119 << ProjectilePDGcode << " " << Plab << G4endl << "Mass KinE " << ProjectileMass
120 << " " << KineticEnergy << G4endl << " A Z " << theA << " " << theZ << G4endl;
121 #endif
122
123 G4double Ylab, Xtotal( 0.0 ), Xelastic( 0.0 ), Xannihilation( 0.0 );
124 G4int NumberOfTargetNucleons;
125
126 Ylab = 0.5 * G4Log( (Elab + Plab)/(Elab - Plab) );
127
128 G4double ECMSsqr = S/GeV/GeV;
129 G4double SqrtS = std::sqrt( S )/GeV;
130
131 #ifdef debugFTFparams
132 G4cout << "Sqrt(s) " << SqrtS << G4endl;
133 #endif
134
135 TargetMass /= GeV; TargetMass2 /= (GeV*GeV);
136 ProjectileMass /= GeV; ProjectileMass2 /= (GeV*GeV);
137
138 Plab /= GeV;
139 G4double Xftf = 0.0;
140
141 G4int NumberOfTargetProtons = theZ;
142 G4int NumberOfTargetNeutrons = theA - theZ;
143 NumberOfTargetNucleons = NumberOfTargetProtons + NumberOfTargetNeutrons;
144
145 // ---------- hadron projectile ----------------
146 if ( AbsProjectileBaryonNumber <= 1 ) { // Projectile is hadron or baryon
147
148 // Interaction on P
149 G4double xTtP = csGGinstance->GetTotalIsotopeCrossSection( particle, KineticEnergy, 1, 1);
150 G4double xElP = csGGinstance->GetElasticIsotopeCrossSection(particle, KineticEnergy, 1, 1);
151
152 // Interaction on N
153 G4double xTtN = csGGinstance->GetTotalIsotopeCrossSection( particle, KineticEnergy, 0, 1);
154 G4double xElN = csGGinstance->GetElasticIsotopeCrossSection(particle, KineticEnergy, 0, 1);
155
156 // Average properties of h+N interactions
157 Xtotal = ( NumberOfTargetProtons * xTtP + NumberOfTargetNeutrons * xTtN ) / NumberOfTargetNucleons;
158 Xelastic = ( NumberOfTargetProtons * xElP + NumberOfTargetNeutrons * xElN ) / NumberOfTargetNucleons;
159 Xannihilation = 0.0;
160
161 Xtotal /= millibarn;
162 Xelastic /= millibarn;
163
164 #ifdef debugFTFparams
165 G4cout<<"Estimated cross sections (total and elastic) of h+N interactions "<<Xtotal<<" "<<Xelastic<<" (mb)"<<G4endl;
166 #endif
167 }
168
169 // ---------- nucleus projectile ----------------
170 if ( ProjectileIsNucleus && ProjectileBaryonNumber > 1 ) {
171
172 #ifdef debugFTFparams
173 G4cout<<"Projectile is a nucleus: A and Z - "<<ProjectileBaryonNumber<<" "<<ProjectileCharge<<G4endl;
174 #endif
175
177 // Interaction on P
178 G4double XtotPP = csGGinstance->GetTotalIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
179 G4double XelPP = csGGinstance->GetElasticIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
180
182 // Interaction on N
183 G4double XtotPN = csGGinstance->GetTotalIsotopeCrossSection(Neutron, KineticEnergy, 0, 1);
184 G4double XelPN = csGGinstance->GetElasticIsotopeCrossSection(Neutron, KineticEnergy, 0, 1);
185
186 #ifdef debugFTFparams
187 G4cout << "XsPP (total and elastic) " << XtotPP/millibarn << " " << XelPP/millibarn <<" (mb)"<< G4endl
188 << "XsPN (total and elastic) " << XtotPN/millibarn << " " << XelPN/millibarn <<" (mb)"<< G4endl;
189 #endif
190
191 Xtotal = (
192 AbsProjectileCharge * NumberOfTargetProtons * XtotPP +
193 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
194 NumberOfTargetNeutrons * XtotPP
195 +
196 ( AbsProjectileCharge * NumberOfTargetNeutrons +
197 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
198 NumberOfTargetProtons ) * XtotPN
199 ) / ( AbsProjectileBaryonNumber * NumberOfTargetNucleons );
200 Xelastic= (
201 AbsProjectileCharge * NumberOfTargetProtons * XelPP +
202 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
203 NumberOfTargetNeutrons * XelPP
204 +
205 ( AbsProjectileCharge * NumberOfTargetNeutrons +
206 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
207 NumberOfTargetProtons ) * XelPN
208 ) / ( AbsProjectileBaryonNumber * NumberOfTargetNucleons );
209
210 Xannihilation = 0.0;
211 Xtotal /= millibarn;
212 Xelastic /= millibarn;
213 }
214
215 // ---------- The projectile is anti-baryon or anti-nucleus ----------------
216 // anti Sigma^0_c anti Delta^-
217 if ( ProjectilePDGcode >= -4112 && ProjectilePDGcode <= -1114 ) {
218 // Only non-strange and strange baryons are considered
219
220 #ifdef debugFTFparams
221 G4cout<<"Projectile is a anti-baryon or anti-nucleus - "<<ProjectileBaryonNumber<<" "<<ProjectileCharge<<G4endl;
222 G4cout<<"(Only non-strange and strange baryons are considered)"<<G4endl;
223 #endif
224
225 G4double X_a( 0.0 ), X_b( 0.0 ), X_c( 0.0 ), X_d( 0.0 );
226 G4double MesonProdThreshold = ProjectileMass + TargetMass +
227 ( 2.0 * 0.14 + 0.016 ); // 2 Mpi + DeltaE;
228
229 if ( PlabPerParticle < 40.0*MeV ) { // Low energy limits. Projectile at rest.
230 Xtotal = 1512.9; // mb
231 Xelastic = 473.2; // mb
232 X_a = 625.1; // mb
233 X_b = 9.780; // mb
234 X_c = 49.989; // mb
235 X_d = 6.614; // mb
236 } else { // Total and elastic cross section of PbarP interactions a'la Arkhipov
237 G4double LogS = G4Log( ECMSsqr / 33.0625 );
238 G4double Xasmpt = 36.04 + 0.304*LogS*LogS; // mb
239 LogS = G4Log( SqrtS / 20.74 );
240 G4double Basmpt = 11.92 + 0.3036*LogS*LogS; // GeV^(-2)
241 G4double R0 = std::sqrt( 0.40874044*Xasmpt - Basmpt ); // GeV^(-1)
242
243 G4double FlowF = SqrtS / std::sqrt( ECMSsqr*ECMSsqr + ProjectileMass2*ProjectileMass2 +
244 TargetMass2*TargetMass2 - 2.0*ECMSsqr*ProjectileMass2
245 - 2.0*ECMSsqr*TargetMass2
246 - 2.0*ProjectileMass2*TargetMass2 );
247
248 Xtotal = Xasmpt * ( 1.0 + 13.55*FlowF/R0/R0/R0*
249 (1.0 - 4.47/SqrtS + 12.38/ECMSsqr - 12.43/SqrtS/ECMSsqr) ); // mb
250
251 Xasmpt = 4.4 + 0.101*LogS*LogS; // mb
252 Xelastic = Xasmpt * ( 1.0 + 59.27*FlowF/R0/R0/R0*
253 (1.0 - 6.95/SqrtS + 23.54/ECMSsqr - 25.34/SqrtS/ECMSsqr ) ); // mb
254
255 X_a = 25.0*FlowF; // mb, 3-shirts diagram
256
257 if ( SqrtS < MesonProdThreshold ) {
258 X_b = 3.13 + 140.0*G4Pow::GetInstance()->powA( MesonProdThreshold - SqrtS, 2.5 ); // mb anti-quark-quark annihilation
259 Xelastic -= 3.0*X_b; // Xel-X(PbarP->NNbar)
260 } else {
261 X_b = 6.8/SqrtS; // mb anti-quark-quark annihilation
262 Xelastic -= 3.0*X_b; // Xel-X(PbarP->NNbar)
263 }
264
265 X_c = 2.0*FlowF*sqr( ProjectileMass + TargetMass )/ECMSsqr; // mb rearrangement
266
267 X_d = 23.3/ECMSsqr; // mb anti-quark-quark string creation
268 }
269
270 G4double Xann_on_P( 0.0), Xann_on_N( 0.0 );
271
272 if ( ProjectilePDGcode == -2212 ) { // Pbar+P/N
273 Xann_on_P = X_a + X_b*5.0 + X_c*5.0 + X_d*6.0;
274 Xann_on_N = X_a + X_b*4.0 + X_c*4.0 + X_d*4.0;
275 } else if ( ProjectilePDGcode == -2112 ) { // NeutrBar+P/N
276 Xann_on_P = X_a + X_b*4.0 + X_c*4.0 + X_d*4.0;
277 Xann_on_N = X_a + X_b*5.0 + X_c*5.0 + X_d*6.0;
278 } else if ( ProjectilePDGcode == -3122 ) { // LambdaBar+P/N
279 Xann_on_P = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
280 Xann_on_N = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
281 } else if ( ProjectilePDGcode == -3112 ) { // Sigma-Bar+P/N
282 Xann_on_P = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
283 Xann_on_N = X_a + X_b*4.0 + X_c*4.0 + X_d*2.0;
284 } else if ( ProjectilePDGcode == -3212 ) { // Sigma0Bar+P/N
285 Xann_on_P = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
286 Xann_on_N = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
287 } else if ( ProjectilePDGcode == -3222 ) { // Sigma+Bar+P/N
288 Xann_on_P = X_a + X_b*4.0 + X_c*4.0 + X_d*2.0;
289 Xann_on_N = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
290 } else if ( ProjectilePDGcode == -3312 ) { // Xi-Bar+P/N
291 Xann_on_P = X_a + X_b*1.0 + X_c*1.0 + X_d*0.0;
292 Xann_on_N = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
293 } else if ( ProjectilePDGcode == -3322 ) { // Xi0Bar+P/N
294 Xann_on_P = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
295 Xann_on_N = X_a + X_b*1.0 + X_c*1.0 + X_d*0.0;
296 } else if ( ProjectilePDGcode == -3334 ) { // Omega-Bar+P/N
297 Xann_on_P = X_a + X_b*0.0 + X_c*0.0 + X_d*0.0;
298 Xann_on_N = X_a + X_b*0.0 + X_c*0.0 + X_d*0.0;
299 } else {
300 G4cout << "Unknown anti-baryon for FTF annihilation" << G4endl;
301 }
302
303 //G4cout << "Sum " << Xann_on_P << G4endl;
304
305 if ( ! ProjectileIsNucleus ) { // Projectile is anti-baryon
306 Xannihilation = ( NumberOfTargetProtons * Xann_on_P + NumberOfTargetNeutrons * Xann_on_N )
307 / NumberOfTargetNucleons;
308 } else { // Projectile is a nucleus
309 Xannihilation = (
310 ( AbsProjectileCharge * NumberOfTargetProtons +
311 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
312 NumberOfTargetNeutrons ) * Xann_on_P
313 +
314 ( AbsProjectileCharge * NumberOfTargetNeutrons +
315 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
316 NumberOfTargetProtons ) * Xann_on_N
317 ) / ( AbsProjectileBaryonNumber * NumberOfTargetNucleons );
318 }
319
320 //G4double Xftf = 0.0;
321 MesonProdThreshold = ProjectileMass + TargetMass + (0.14 + 0.08); // Mpi + DeltaE
322 if ( SqrtS > MesonProdThreshold ) {
323 Xftf = 36.0 * ( 1.0 - MesonProdThreshold/SqrtS );
324 }
325
326 Xtotal = Xelastic + Xannihilation + Xftf;
327
328 #ifdef debugFTFparams
329 G4cout << "Plab Xtotal, Xelastic Xinel Xftf " << Plab << " " << Xtotal << " " << Xelastic
330 << " " << Xtotal - Xelastic << " " << Xtotal - Xelastic - Xannihilation << " (mb)"<< G4endl
331 << "Plab Xelastic/Xtotal, Xann/Xin " << Plab << " " << Xelastic/Xtotal << " "
332 << Xannihilation/(Xtotal - Xelastic) << G4endl;
333 #endif
334
335 }
336
337 if ( Xtotal == 0.0 ) { // Projectile is undefined, Nucleon assumed
338
340 // Interaction on P
341 G4double XtotPP = csGGinstance->GetTotalIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
342 G4double XelPP = csGGinstance->GetElasticIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
343
344 // Interaction on N
345 G4double XtotPN = csGGinstance->GetTotalIsotopeCrossSection(Proton, KineticEnergy, 0, 1);
346 G4double XelPN = csGGinstance->GetElasticIsotopeCrossSection(Proton, KineticEnergy, 0, 1);
347
348 Xtotal = ( NumberOfTargetProtons * XtotPP + NumberOfTargetNeutrons * XtotPN )
349 / NumberOfTargetNucleons;
350 Xelastic = ( NumberOfTargetProtons * XelPP + NumberOfTargetNeutrons * XelPN )
351 / NumberOfTargetNucleons;
352 Xannihilation = 0.0;
353 Xtotal /= millibarn;
354 Xelastic /= millibarn;
355 };
356
357 // Geometrical parameters
358 SetTotalCrossSection( Xtotal );
359 SetElastisCrossSection( Xelastic );
360 SetInelasticCrossSection( Xtotal - Xelastic );
361
362 // Interactions with elastic and inelastic collisions
363 SetProbabilityOfElasticScatt( Xtotal, Xelastic );
364
365 SetRadiusOfHNinteractions2( Xtotal/pi/10.0 );
366
367 if ( ( Xtotal - Xelastic ) == 0.0 ) {
369 } else {
370 SetProbabilityOfAnnihilation( Xannihilation / (Xtotal - Xelastic) );
371 }
372
373 if(Xelastic > 0.0) {
374 SetSlope( Xtotal*Xtotal/16.0/pi/Xelastic/0.3894 );// Slope parameter of elastic scattering
375 // (GeV/c)^(-2))
376 // Parameters of elastic scattering
377 // Gaussian parametrization of elastic scattering amplitude assumed
378 SetAvaragePt2ofElasticScattering( 1.0/( Xtotal*Xtotal/16.0/pi/Xelastic/0.3894 )*GeV*GeV );
379 } else {
380 SetSlope(1.0);
382 }
383 SetGamma0( GetSlope()*Xtotal/10.0/2.0/pi );
384
385 G4double Xinel = Xtotal - Xelastic;
386
387 #ifdef debugFTFparams
388 G4cout<< "Slope of hN elastic scattering" << GetSlope() << G4endl;
389 G4cout << "AvaragePt2ofElasticScattering " << GetAvaragePt2ofElasticScattering() << G4endl;
390 G4cout<<"Parameters of excitation for projectile "<<ProjectilePDGcode<< G4endl;
391 #endif
392
393 if ( (ProjectilePDGcode == 2212) || (ProjectilePDGcode == 2112) ) { // Projectile is proton or neutron
394
395 // A process probability is parameterized as Prob = A_1*exp(-A_2*y) + A_3*exp(-A_4*y) + A_top
396 // y is a rapidity of a partcle in the target nucleus. Ymin is a minimal rapidity below it X=0
397
398 // Proc# A1 B1 A2 B2 A3 Atop Ymin
399 /* original hadr-string-diff-V10-03-07 (similar to 10.3.x)
400 SetParams( 0, 13.71, 1.75, -214.5, 4.25, 0.0, 0.5 , 1.1 ); // Qexchange without Exc.
401 SetParams( 1, 25.0, 1.0, -50.34, 1.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc.
402 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Projectile diffraction
403 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Target diffraction
404 SetParams( 4, 1.0, 0.0 , -2.01 , 0.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc. Additional multiplier
405 */
406
407 // Proc#
411 fParCollBaryonProj.GetProc0Ymin() ); // Qexchange without Exc.
415 fParCollBaryonProj.GetProc1Ymin() ); // Qexchange with Exc.
416 if ( Xinel > 0.0 ) {
417 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93 ); // Projectile diffraction
418 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93 ); // Target diffraction
419
423 fParCollBaryonProj.GetProc4Ymin() ); // Qexchange with Exc. Additional multiplier
424 } else { // if Xinel=0., zero everything out (obviously)
425 SetParams( 2, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
426 SetParams( 3, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
427 SetParams( 4, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
428 }
429
430 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
431 //
432 // It is not decided what to do with diffraction dissociation in Had-Nucl and Nucl-Nucl interactions
433 // For the moment both ProjDiffDisso & TgtDiffDisso for A > 10 are set to false,
434 // so both projectile and target diffraction are turned OFF
435 //
437 SetParams( 2, 0.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
439 SetParams( 3, 0.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
440 }
441
443
444 if ( NumberOfTargetNucleons > 26 ) {
446 } else {
448 }
449
452
455
459
460 } else if ( ProjectilePDGcode == -2212 || ProjectilePDGcode == -2112 ) { // Projectile is anti_proton or anti_neutron
461
462 // Proc# A1 B1 A2 B2 A3 Atop Ymin
463 SetParams( 0, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 1000.0 ); // Qexchange without Exc.
464 SetParams( 1, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 1000.0 ); // Qexchange with Exc.
465 if ( Xinel > 0.) {
466 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93 ); // Projectile diffraction
467 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93 ); // Target diffraction
468 SetParams( 4, 1.0, 0.0 , 0.0, 0.0 , 0.0, 0.0 , 0.93 ); // Qexchange with Exc. Additional multiply
469 } else {
470 SetParams( 2, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0 );
471 SetParams( 3, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0 );
472 SetParams( 4, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0 );
473 }
474
475 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
476 //
477 // It is not decided what to do with diffraction dissociation in Had-Nucl and Nucl-Nucl interactions
478 // For the moment both ProjDiffDisso & TgtDiffDisso are set to false,
479 // so both projectile and target diffraction are turned OFF
480 //
482 SetParams( 2, 0.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
484 SetParams( 3, 0.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
485 }
486
489 SetProjMinDiffMass( ProjectileMass + 0.22 ); // GeV
490 SetProjMinNonDiffMass( ProjectileMass + 0.22 ); // GeV
491 SetTarMinDiffMass( TargetMass + 0.22 ); // GeV
492 SetTarMinNonDiffMass( TargetMass + 0.22 ); // GeV
493 SetAveragePt2( 0.3 ); // GeV^2
494 SetProbLogDistrPrD( 0.55 );
495 SetProbLogDistr( 0.55 );
496
497 } else if ( ProjectileabsPDGcode == 211 || ProjectilePDGcode == 111 ) { // Projectile is Pion
498
499 // Proc# A1 B1 A2 B2 A3 Atop Ymin
500 /* --> original code
501 SetParams( 0, 150.0, 1.8 , -247.3, 2.3, 0., 1. , 2.3 );
502 SetParams( 1, 5.77, 0.6 , -5.77, 0.8, 0., 0. , 0.0 );
503 SetParams( 2, 2.27, 0.5 , -98052.0, 4.0, 0., 0. , 3.0 );
504 SetParams( 3, 7.0, 0.9, -85.28, 1.9, 0.08, 0. , 2.2 );
505 SetParams( 4, 1.0, 0.0 , -11.02, 1.0, 0.0, 0. , 2.4 ); // Qexchange with Exc. Additional multiply
506 */
507 // Proc#
511 fParCollPionProj.GetProc0Ymin() ); // Qexchange without Exc.
515 fParCollPionProj.GetProc1Ymin() ); // Qexchange with Exc.
519 fParCollPionProj.GetProc2Ymin() ); // Projectile diffraction
523 fParCollPionProj.GetProc3Ymin() ); // Target diffraction
527 fParCollPionProj.GetProc4Ymin() ); // Qexchange with Exc. Additional multiply
528
529 // NOTE: how can it be |ProjectileBaryonNumber| > 10 if projectile is a pion ???
530 //
531 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
533 SetParams( 2, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
535 SetParams( 3, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
536 }
537
538 /* original code -->
539 SetDeltaProbAtQuarkExchange( 0.56 );
540 SetProjMinDiffMass( 1.0 ); // GeV
541 SetProjMinNonDiffMass( 1.0 ); // GeV
542 SetTarMinDiffMass( 1.16 ); // GeV
543 SetTarMinNonDiffMass( 1.16 ); // GeV
544 SetAveragePt2( 0.3 ); // GeV^2
545 SetProbLogDistrPrD( 0.55 );
546 SetProbLogDistr( 0.55 );
547 */
548
549 // JVY update, Aug.8, 2018 --> Feb.14, 2019
550 //
559
560 // ---> end update
561
562 } else if ( ProjectileabsPDGcode == 321 || ProjectileabsPDGcode == 311 ||
563 ProjectilePDGcode == 130 || ProjectilePDGcode == 310 ) { // Projectile is Kaon
564
565 // Proc# A1 B1 A2 B2 A3 Atop Ymin
566 SetParams( 0, 60.0 , 2.5 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Qexchange without Exc.
567 SetParams( 1, 6.0 , 1.0 , -24.33 , 2.0 , 0.0, 0.0 , 1.40 ); // Qexchange with Exc.
568 SetParams( 2, 2.76, 1.2 , -22.5 , 2.7 ,0.04, 0.0 , 1.40 ); // Projectile diffraction
569 SetParams( 3, 1.09, 0.5 , -8.88 , 2. ,0.05, 0.0 , 1.40 ); // Target diffraction
570 SetParams( 4, 1.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 0.93 ); // Qexchange with Exc. Additional multiply
571 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
572 SetParams( 2, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
573 SetParams( 3, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
574 }
575
577 SetProjMinDiffMass( 0.7 ); // GeV
578 SetProjMinNonDiffMass( 0.7 ); // GeV
579 SetTarMinDiffMass( 1.16 ); // GeV
580 SetTarMinNonDiffMass( 1.16 ); // GeV
581 SetAveragePt2( 0.3 ); // GeV^2
582 SetProbLogDistrPrD( 0.55 );
583 SetProbLogDistr( 0.55 );
584
585 } else { // Projectile is not p, n, Pi0, Pi+, Pi-, K+, K-, K0 or their anti-particles
586
587 if ( ProjectileabsPDGcode > 1000 ) { // The projectile is a baryon as P or N
588 // Proc# A1 B1 A2 B2 A3 Atop Ymin
589 SetParams( 0, 13.71, 1.75, -30.69, 3.0 , 0.0, 1.0 , 0.93 ); // Qexchange without Exc.
590 SetParams( 1, 25.0, 1.0, -50.34, 1.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc.
591 if ( Xinel > 0.) {
592 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Projectile diffraction
593 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Target diffraction
594 SetParams( 4, 1.0, 0.0 , -2.01 , 0.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc. Additional multiply
595 } else {
596 SetParams( 2, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
597 SetParams( 3, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
598 SetParams( 4, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
599 }
600
601 } else { // The projectile is a meson as K+-0
602 // Proc# A1 B1 A2 B2 A3 Atop Ymin
603 SetParams( 0, 60.0 , 2.5 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Qexchange without Exc.
604 SetParams( 1, 6.0 , 1.0 , -24.33 , 2.0 , 0.0, 0.0 , 1.40 ); // Qexchange with Exc.
605 SetParams( 2, 2.76, 1.2 , -22.5 , 2.7 ,0.04, 0.0 , 1.40 ); // Projectile diffraction
606 SetParams( 3, 1.09, 0.5 , -8.88 , 2. ,0.05, 0.0 , 1.40 ); // Target diffraction
607 SetParams( 4, 1.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 0.93 ); // Qexchange with Exc. Additional multiply
608 }
609
610 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
611 SetParams( 2, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
612 SetParams( 3, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
613 }
614
617
618 SetProjMinDiffMass( GetMinMass(particle)/GeV );
620
624
625 SetAveragePt2( 0.3 ); // GeV^2
626 SetProbLogDistrPrD( 0.55 );
627 SetProbLogDistr( 0.55 );
628
629 }
630
631 #ifdef debugFTFparams
632 G4cout<<"DeltaProbAtQuarkExchange "<< GetDeltaProbAtQuarkExchange() << G4endl;
633 G4cout<<"ProbOfSameQuarkExchange "<< GetProbOfSameQuarkExchange() << G4endl;
634 G4cout<<"ProjMinDiffMass "<< GetProjMinDiffMass()/GeV <<" GeV"<< G4endl;
635 G4cout<<"ProjMinNonDiffMass "<< GetProjMinNonDiffMass() <<" GeV"<< G4endl;
636 G4cout<<"TarMinDiffMass "<< GetTarMinDiffMass() <<" GeV"<< G4endl;
637 G4cout<<"TarMinNonDiffMass "<< GetTarMinNonDiffMass() <<" GeV"<< G4endl;
638 G4cout<<"AveragePt2 "<< GetAveragePt2() <<" GeV^2"<< G4endl;
639 G4cout<<"ProbLogDistrPrD "<< GetProbLogDistrPrD() << G4endl;
640 G4cout<<"ProbLogDistrTrD "<< GetProbLogDistr() << G4endl;
641 #endif
642
643 // Set parameters of nuclear destruction
644
645 if ( ProjectileabsPDGcode < 1000 ) { // Meson projectile
646
647 SetMaxNumberOfCollisions( Plab, 2.0 ); // 3.0 )
648 //
649 // target destruction
650 //
651 /* original code --->
652 SetCofNuclearDestruction( 0.00481*G4double(NumberOfTargetNucleons)*
653 G4Exp( 4.0*(Ylab - 2.1) )/( 1.0 + G4Exp( 4.0*(Ylab - 2.1) ) ) );
654
655 SetR2ofNuclearDestruction( 1.5*fermi*fermi );
656 SetDofNuclearDestruction( 0.3 );
657 SetPt2ofNuclearDestruction( ( 0.035 + 0.04*G4Exp( 4.0*(Ylab - 2.5) )/
658 ( 1.0 + G4Exp( 4.0*(Ylab - 2.5) ) ) )*GeV*GeV );
659 SetMaxPt2ofNuclearDestruction( 1.0*GeV*GeV );
660 SetExcitationEnergyPerWoundedNucleon( 40.0*MeV );
661 */
663 //
664 // NOTE (JVY): Set this switch to false/true on line 138
665 //
667 {
668 coeff *= G4double(NumberOfTargetNucleons);
669 }
672 coeff *= exfactor;
673 coeff /= ( 1.+ exfactor );
674
676
682 coeff *= exfactor;
683 coeff /= ( 1. + exfactor );
685
688
689 } else if ( ProjectilePDGcode == -2212 || ProjectilePDGcode == -2112 ) { // for anti-baryon projectile
690
691 SetMaxNumberOfCollisions( Plab, 2.0 );
692
693 SetCofNuclearDestruction( 0.00481*G4double(NumberOfTargetNucleons)*
694 G4Exp( 4.0*(Ylab - 2.1) )/( 1.0 + G4Exp( 4.0*(Ylab - 2.1) ) ) );
697 SetPt2ofNuclearDestruction( ( 0.035 + 0.04*G4Exp( 4.0*(Ylab - 2.5) )/
698 ( 1.0 + G4Exp( 4.0*(Ylab - 2.5) ) ) )*GeV*GeV );
701 if ( Plab < 2.0 ) { // 2 GeV/c
702 // For slow anti-baryon we have to garanty putting on mass-shell
704 SetR2ofNuclearDestruction( 1.5*fermi*fermi ); // this is equivalent to setting a few line above
705 // is it even necessary ?
709 }
710
711 } else { // Projectile baryon assumed
712
713 // NOTE (JVY) FIXME !!! Will decide later how/if to make this one configurable...
714 //
715 SetMaxNumberOfCollisions( Plab, 2.0 );
716
717 // projectile destruction - does NOT really matter for particle projectile, only for a nucleus projectile
718 //
719 double coeff = 0.;
721 //
722 // NOTE (JVY): Set this switch to false/true on line 136
723 //
725 {
726 coeff *= G4double(AbsProjectileBaryonNumber);
727 }
729 coeff *= exfactor;
730 coeff /= ( 1.+ exfactor );
732
733 // target desctruction
734 //
736 //
737 // NOTE (JVY): Set this switch to false/true on line 138
738 //
740 {
741 coeff *= G4double(NumberOfTargetNucleons);
742 }
744 coeff *= exfactor;
745 coeff /= ( 1.+ exfactor );
747
750
753 coeff *= exfactor;
754 coeff /= ( 1. + exfactor );
756
759
760 }
761
762 #ifdef debugFTFparams
763 G4cout<<"CofNuclearDestructionPr "<< GetCofNuclearDestructionPr() << G4endl;
764 G4cout<<"CofNuclearDestructionTr "<< GetCofNuclearDestruction() << G4endl;
765 G4cout<<"R2ofNuclearDestruction "<< GetR2ofNuclearDestruction()/fermi/fermi <<" fermi^2"<< G4endl;
766 G4cout<<"DofNuclearDestruction "<< GetDofNuclearDestruction() << G4endl;
767 G4cout<<"Pt2ofNuclearDestruction "<< GetPt2ofNuclearDestruction()/GeV/GeV <<" GeV^2"<< G4endl;
768 G4cout<<"ExcitationEnergyPerWoundedNucleon "<< GetExcitationEnergyPerWoundedNucleon() <<" MeV"<< G4endl;
769 #endif
770
771 //SetCofNuclearDestruction( 0.47*G4Exp( 2.0*(Ylab - 2.5) )/( 1.0 + G4Exp( 2.0*(Ylab - 2.5) ) ) );
772 //SetPt2ofNuclearDestruction( ( 0.035 + 0.1*G4Exp( 4.0*(Ylab - 3.0) )/( 1.0 + G4Exp( 4.0*(Ylab - 3.0) ) ) )*GeV*GeV );
773
774 //SetMagQuarkExchange( 120.0 ); // 210.0 PipP
775 //SetSlopeQuarkExchange( 2.0 );
776 //SetDeltaProbAtQuarkExchange( 0.6 );
777 //SetProjMinDiffMass( 0.7 ); // GeV 1.1
778 //SetProjMinNonDiffMass( 0.7 ); // GeV
779 //SetProbabilityOfProjDiff( 0.0); // 0.85*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.5 ) ); // 40/32 X-dif/X-inel
780 //SetTarMinDiffMass( 1.1 ); // GeV
781 //SetTarMinNonDiffMass( 1.1 ); // GeV
782 //SetProbabilityOfTarDiff( 0.0 ); // 0.85*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.5 ) ); // 40/32 X-dif/X-inel
783
784 //SetAveragePt2( 0.0 ); // GeV^2 0.3
785 //------------------------------------
786 //SetProbabilityOfElasticScatt( 1.0, 1.0); //(Xtotal, Xelastic);
787 //SetProbabilityOfProjDiff( 1.0*0.62*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.51 ) ); // 0->1
788 //SetProbabilityOfTarDiff( 4.0*0.62*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.51 ) ); // 2->4
789 //SetAveragePt2( 0.3 ); // (0.15)
790 //SetAvaragePt2ofElasticScattering( 0.0 );
791
792 //SetMaxNumberOfCollisions( Plab, 6.0 ); //(4.0*(Plab + 0.01), Plab); // 6.0 );
793 //SetAveragePt2( 0.15 );
794 //SetCofNuclearDestruction(-1.);//( 0.75 ); // (0.25)
795 //SetExcitationEnergyPerWoundedNucleon(0.);//( 30.0*MeV ); // (75.0*MeV)
796 //SetDofNuclearDestruction(0.);//( 0.2 ); //0.4 // 0.3 0.5
797
798 /*
799 SetAveragePt2(0.3);
800 SetCofNuclearDestructionPr(0.);
801 SetCofNuclearDestruction(0.); //( 0.5 ); (0.25)
802 SetExcitationEnergyPerWoundedNucleon(0.); // 30.0*MeV; (75.0*MeV)
803 SetDofNuclearDestruction(0.); // 0.2; 0.4; 0.3; 0.5
804 SetPt2ofNuclearDestruction(0.); //(2.*0.075*GeV*GeV); ( 0.3*GeV*GeV ); (0.168*GeV*GeV)
805 */
806
807 //SetExcitationEnergyPerWoundedNucleon(0.001);
808 //SetPt2Kink( 0.0*GeV*GeV );
809
810 //SetRadiusOfHNinteractions2( Xtotal/pi/10.0 /2.);
811 //SetRadiusOfHNinteractions2( (Xtotal - Xelastic)/pi/10.0 );
812 //SetProbabilityOfElasticScatt( 1.0, 0.0);
813 /*
814 G4cout << "Pt2 " << GetAveragePt2()<<" "<<GetAveragePt2()/GeV/GeV<<G4endl;
815 G4cout << "Cnd " << GetCofNuclearDestruction() << G4endl;
816 G4cout << "Dnd " << GetDofNuclearDestruction() << G4endl;
817 G4cout << "Pt2 " << GetPt2ofNuclearDestruction()/GeV/GeV << G4endl;
818 */
819
820}
G4double S(G4double temp)
G4double G4Log(G4double x)
Definition: G4Log.hh:226
static constexpr double millibarn
Definition: G4SIunits.hh:86
static constexpr double fermi
Definition: G4SIunits.hh:83
static constexpr double MeV
Definition: G4SIunits.hh:200
static constexpr double pi
Definition: G4SIunits.hh:55
bool G4bool
Definition: G4Types.hh:86
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
double GetProc3B2() const
double GetProc4A3() const
double GetNuclearTgtDestructP3() const
double GetProc3B1() const
bool IsNuclearProjDestructP1_NBRNDEP() const
double GetProc0A2() const
double GetExciEnergyPerWoundedNucleon() const
double GetProc1Atop() const
double GetProc2A3() const
double GetNuclearProjDestructP2() const
double GetProbLogDistrPrD() const
double GetPt2NuclearDestructP4() const
double GetPt2NuclearDestructP3() const
double GetNuclearProjDestructP3() const
double GetProc4A2() const
double GetPt2NuclearDestructP1() const
double GetProc2Ymin() const
double GetProjMinNonDiffMass() const
double GetProc3Ymin() const
double GetProc2A2() const
double GetProc4B1() const
double GetProc1A1() const
double GetProc0B2() const
double GetProbLogDistr() const
double GetProc1B2() const
double GetProc3Atop() const
double GetProc1A3() const
double GetProc1B1() const
double GetProc1Ymin() const
double GetProc2B1() const
double GetNuclearTgtDestructP2() const
double GetPt2NuclearDestructP2() const
double GetR2ofNuclearDestruct() const
double GetDeltaProbAtQuarkExchange() const
double GetDofNuclearDestruct() const
double GetProc0A3() const
bool IsProjDiffDissociation() const
double GetProc0Atop() const
double GetProc2B2() const
double GetProc2A1() const
double GetProc2Atop() const
double GetProc1A2() const
double GetProc4Atop() const
double GetProc3A1() const
double GetAveragePt2() const
double GetProc4Ymin() const
double GetProc0B1() const
double GetTgtMinNonDiffMass() const
bool IsTgtDiffDissociation() const
double GetProbOfSameQuarkExchange() const
double GetProc4B2() const
double GetProc0Ymin() const
double GetProc3A2() const
double GetProc3A3() const
bool IsNuclearTgtDestructP1_ADEP() const
double GetMaxPt2ofNuclearDestruct() const
double GetProc4A1() const
double GetTgtMinDiffMass() const
double GetNuclearTgtDestructP1() const
double GetProjMinDiffMass() const
double GetProc0A1() const
double GetNuclearProjDestructP1() const
void SetTarMinNonDiffMass(const G4double aValue)
void SetProjMinDiffMass(const G4double aValue)
G4double GetProbLogDistrPrD()
G4double GetCofNuclearDestructionPr()
void SetTotalCrossSection(const G4double Xtotal)
void SetExcitationEnergyPerWoundedNucleon(const G4double aValue)
G4double GetAveragePt2()
void SetProbabilityOfElasticScatt(const G4double Xtotal, const G4double Xelastic)
void SetElastisCrossSection(const G4double Xelastic)
G4double GetPt2ofNuclearDestruction()
G4double GetProbLogDistr()
G4FTFParamCollBaryonProj fParCollBaryonProj
void SetMaxPt2ofNuclearDestruction(const G4double aValue)
G4double GetProjMinNonDiffMass()
G4double GetAvaragePt2ofElasticScattering()
void SetSlope(const G4double Slope)
G4double GetTarMinDiffMass()
void SetDeltaProbAtQuarkExchange(const G4double aValue)
G4double GetTarMinNonDiffMass()
G4double GetMinMass(const G4ParticleDefinition *aParticle)
void SetAvaragePt2ofElasticScattering(const G4double aPt2)
void SetTarMinDiffMass(const G4double aValue)
void SetProjMinNonDiffMass(const G4double aValue)
void SetCofNuclearDestructionPr(const G4double aValue)
void SetR2ofNuclearDestruction(const G4double aValue)
void SetProbLogDistrPrD(const G4double aValue)
void SetGamma0(const G4double Gamma0)
void SetAveragePt2(const G4double aValue)
G4double GetDeltaProbAtQuarkExchange()
G4double GetExcitationEnergyPerWoundedNucleon()
void SetProbabilityOfAnnihilation(const G4double aValue)
void SetDofNuclearDestruction(const G4double aValue)
G4double GetProjMinDiffMass()
G4double GetProbOfSameQuarkExchange()
G4double GetDofNuclearDestruction()
void SetRadiusOfHNinteractions2(const G4double Radius2)
void SetParams(const G4int ProcN, const G4double A1, const G4double B1, const G4double A2, const G4double B2, const G4double A3, const G4double Atop, const G4double Ymin)
void SetPt2ofNuclearDestruction(const G4double aValue)
G4double GetR2ofNuclearDestruction()
G4double GetCofNuclearDestruction()
void SetProbLogDistr(const G4double aValue)
void SetCofNuclearDestruction(const G4double aValue)
void SetMaxNumberOfCollisions(const G4double aValue, const G4double bValue)
void SetProbOfSameQuarkExchange(const G4double aValue)
G4FTFParamCollPionProj fParCollPionProj
G4FTFParamCollMesonProj fParCollMesonProj
void SetInelasticCrossSection(const G4double Xinelastic)
static G4Neutron * Neutron()
Definition: G4Neutron.cc:103
G4double GetPDGCharge() const
static G4Pow * GetInstance()
Definition: G4Pow.cc:41
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:230
static G4Proton * Proton()
Definition: G4Proton.cc:92
virtual G4double GetTotalIsotopeCrossSection(const G4ParticleDefinition *, G4double kinEnergy, G4int, G4int)=0
virtual G4double GetElasticIsotopeCrossSection(const G4ParticleDefinition *, G4double kinEnergy, G4int, G4int)=0
static constexpr double GeV
T sqr(const T &x)
Definition: templates.hh:128

References csGGinstance, fermi, fParCollBaryonProj, fParCollMesonProj, fParCollPionProj, G4cout, G4endl, G4Exp(), G4Log(), GetAvaragePt2ofElasticScattering(), GetAveragePt2(), G4FTFParamCollection::GetAveragePt2(), G4ParticleDefinition::GetBaryonNumber(), GetCofNuclearDestruction(), GetCofNuclearDestructionPr(), GetDeltaProbAtQuarkExchange(), G4FTFParamCollection::GetDeltaProbAtQuarkExchange(), G4FTFParamCollection::GetDofNuclearDestruct(), GetDofNuclearDestruction(), G4VComponentCrossSection::GetElasticIsotopeCrossSection(), G4FTFParamCollection::GetExciEnergyPerWoundedNucleon(), GetExcitationEnergyPerWoundedNucleon(), G4Pow::GetInstance(), G4FTFParamCollection::GetMaxPt2ofNuclearDestruct(), GetMinMass(), G4FTFParamCollection::GetNuclearProjDestructP1(), G4FTFParamCollection::GetNuclearProjDestructP2(), G4FTFParamCollection::GetNuclearProjDestructP3(), G4FTFParamCollection::GetNuclearTgtDestructP1(), G4FTFParamCollection::GetNuclearTgtDestructP2(), G4FTFParamCollection::GetNuclearTgtDestructP3(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGEncoding(), G4ParticleDefinition::GetPDGMass(), GetProbLogDistr(), G4FTFParamCollection::GetProbLogDistr(), GetProbLogDistrPrD(), G4FTFParamCollection::GetProbLogDistrPrD(), GetProbOfSameQuarkExchange(), G4FTFParamCollection::GetProbOfSameQuarkExchange(), G4FTFParamCollection::GetProc0A1(), G4FTFParamCollection::GetProc0A2(), G4FTFParamCollection::GetProc0A3(), G4FTFParamCollection::GetProc0Atop(), G4FTFParamCollection::GetProc0B1(), G4FTFParamCollection::GetProc0B2(), G4FTFParamCollection::GetProc0Ymin(), G4FTFParamCollection::GetProc1A1(), G4FTFParamCollection::GetProc1A2(), G4FTFParamCollection::GetProc1A3(), G4FTFParamCollection::GetProc1Atop(), G4FTFParamCollection::GetProc1B1(), G4FTFParamCollection::GetProc1B2(), G4FTFParamCollection::GetProc1Ymin(), G4FTFParamCollection::GetProc2A1(), G4FTFParamCollection::GetProc2A2(), G4FTFParamCollection::GetProc2A3(), G4FTFParamCollection::GetProc2Atop(), G4FTFParamCollection::GetProc2B1(), G4FTFParamCollection::GetProc2B2(), G4FTFParamCollection::GetProc2Ymin(), G4FTFParamCollection::GetProc3A1(), G4FTFParamCollection::GetProc3A2(), G4FTFParamCollection::GetProc3A3(), G4FTFParamCollection::GetProc3Atop(), G4FTFParamCollection::GetProc3B1(), G4FTFParamCollection::GetProc3B2(), G4FTFParamCollection::GetProc3Ymin(), G4FTFParamCollection::GetProc4A1(), G4FTFParamCollection::GetProc4A2(), G4FTFParamCollection::GetProc4A3(), G4FTFParamCollection::GetProc4Atop(), G4FTFParamCollection::GetProc4B1(), G4FTFParamCollection::GetProc4B2(), G4FTFParamCollection::GetProc4Ymin(), GetProjMinDiffMass(), G4FTFParamCollection::GetProjMinDiffMass(), GetProjMinNonDiffMass(), G4FTFParamCollection::GetProjMinNonDiffMass(), G4FTFParamCollection::GetPt2NuclearDestructP1(), G4FTFParamCollection::GetPt2NuclearDestructP2(), G4FTFParamCollection::GetPt2NuclearDestructP3(), G4FTFParamCollection::GetPt2NuclearDestructP4(), GetPt2ofNuclearDestruction(), G4FTFParamCollection::GetR2ofNuclearDestruct(), GetR2ofNuclearDestruction(), GetSlope(), GetTarMinDiffMass(), GetTarMinNonDiffMass(), G4FTFParamCollection::GetTgtMinDiffMass(), G4FTFParamCollection::GetTgtMinNonDiffMass(), G4VComponentCrossSection::GetTotalIsotopeCrossSection(), CLHEP::GeV, GeV, G4FTFParamCollection::IsNuclearProjDestructP1_NBRNDEP(), G4FTFParamCollection::IsNuclearTgtDestructP1_ADEP(), G4FTFParamCollection::IsProjDiffDissociation(), G4FTFParamCollection::IsTgtDiffDissociation(), MeV, millibarn, G4Neutron::Neutron(), Neutron, pi, G4Pow::powA(), G4Proton::Proton(), Proton, Reset(), S(), SetAvaragePt2ofElasticScattering(), SetAveragePt2(), SetCofNuclearDestruction(), SetCofNuclearDestructionPr(), SetDeltaProbAtQuarkExchange(), SetDofNuclearDestruction(), SetElastisCrossSection(), SetExcitationEnergyPerWoundedNucleon(), SetGamma0(), SetInelasticCrossSection(), SetMaxNumberOfCollisions(), SetMaxPt2ofNuclearDestruction(), SetParams(), SetProbabilityOfAnnihilation(), SetProbabilityOfElasticScatt(), SetProbLogDistr(), SetProbLogDistrPrD(), SetProbOfSameQuarkExchange(), SetProjMinDiffMass(), SetProjMinNonDiffMass(), SetPt2ofNuclearDestruction(), SetR2ofNuclearDestruction(), SetRadiusOfHNinteractions2(), SetSlope(), SetTarMinDiffMass(), SetTarMinNonDiffMass(), SetTotalCrossSection(), and sqr().

Referenced by G4FTFModel::Init().

◆ Reset()

void G4FTFParameters::Reset ( )
private

Definition at line 1531 of file G4FTFParameters.cc.

1532{
1533 FTFhNcmsEnergy = 0.0;
1534 FTFXtotal = 0.0;
1535 FTFXelastic = 0.0;
1536 FTFXinelastic = 0.0;
1537 FTFXannihilation = 0.0;
1541 FTFSlope = 0.0;
1543 FTFGamma0 = 0.0;
1546 ProjMinDiffMass = 0.0;
1547 ProjMinNonDiffMass = 0.0;
1548 ProbLogDistrPrD = 0.0;
1549 TarMinDiffMass = 0.0;
1550 TarMinNonDiffMass = 0.0;
1551 AveragePt2 = 0.0;
1552 ProbLogDistr = 0.0;
1553 Pt2kink = 0.0;
1554 MaxNumberOfCollisions = 0.0;
1555 ProbOfInelInteraction = 0.0;
1560 DofNuclearDestruction = 0.0;
1563
1564 for ( G4int i = 0; i < 4; i++ ) {
1565 for ( G4int j = 0; j < 7; j++ ) {
1566 ProcParams[i][j] = 0.0;
1567 }
1568 }
1569
1570 return;
1571}
G4double FTFXannihilation

References AvaragePt2ofElasticScattering, AveragePt2, CofNuclearDestruction, CofNuclearDestructionPr, DeltaProbAtQuarkExchange, DofNuclearDestruction, ExcitationEnergyPerWoundedNucleon, FTFGamma0, FTFhNcmsEnergy, FTFSlope, FTFXannihilation, FTFXelastic, FTFXinelastic, FTFXtotal, MaxNumberOfCollisions, MaxPt2ofNuclearDestruction, ProbabilityOfAnnihilation, ProbabilityOfElasticScatt, ProbLogDistr, ProbLogDistrPrD, ProbOfInelInteraction, ProbOfSameQuarkExchange, ProcParams, ProjMinDiffMass, ProjMinNonDiffMass, Pt2kink, Pt2ofNuclearDestruction, R2ofNuclearDestruction, RadiusOfHNinteractions2, TarMinDiffMass, and TarMinNonDiffMass.

Referenced by G4FTFParameters(), and InitForInteraction().

◆ SetAvaragePt2ofElasticScattering()

void G4FTFParameters::SetAvaragePt2ofElasticScattering ( const G4double  aPt2)
inline

Definition at line 537 of file G4FTFParameters.hh.

537 {
539}

References AvaragePt2ofElasticScattering.

Referenced by InitForInteraction().

◆ SetAveragePt2()

void G4FTFParameters::SetAveragePt2 ( const G4double  aValue)
inline

Definition at line 577 of file G4FTFParameters.hh.

577 {
579}

References AveragePt2, and CLHEP::GeV.

Referenced by InitForInteraction().

◆ SetCofNuclearDestruction()

void G4FTFParameters::SetCofNuclearDestruction ( const G4double  aValue)
inline

Definition at line 625 of file G4FTFParameters.hh.

625 {
626 CofNuclearDestruction = aValue;
627}

References CofNuclearDestruction.

Referenced by InitForInteraction().

◆ SetCofNuclearDestructionPr()

void G4FTFParameters::SetCofNuclearDestructionPr ( const G4double  aValue)
inline

Definition at line 621 of file G4FTFParameters.hh.

621 {
623}

References CofNuclearDestructionPr.

Referenced by InitForInteraction().

◆ SetDeltaProbAtQuarkExchange()

void G4FTFParameters::SetDeltaProbAtQuarkExchange ( const G4double  aValue)
inline

Definition at line 553 of file G4FTFParameters.hh.

553 {
555}

References DeltaProbAtQuarkExchange.

Referenced by InitForInteraction().

◆ SetDofNuclearDestruction()

void G4FTFParameters::SetDofNuclearDestruction ( const G4double  aValue)
inline

Definition at line 637 of file G4FTFParameters.hh.

637 {
638 DofNuclearDestruction = aValue;
639}

References DofNuclearDestruction.

Referenced by InitForInteraction().

◆ SetElastisCrossSection()

void G4FTFParameters::SetElastisCrossSection ( const G4double  Xelastic)
inline

Definition at line 499 of file G4FTFParameters.hh.

499 {
500 FTFXelastic = Xelastic;
501}

References FTFXelastic.

Referenced by InitForInteraction().

◆ SetExcitationEnergyPerWoundedNucleon()

void G4FTFParameters::SetExcitationEnergyPerWoundedNucleon ( const G4double  aValue)
inline

Definition at line 633 of file G4FTFParameters.hh.

633 {
635}

References ExcitationEnergyPerWoundedNucleon.

Referenced by InitForInteraction().

◆ SetGamma0()

void G4FTFParameters::SetGamma0 ( const G4double  Gamma0)
inline

Definition at line 532 of file G4FTFParameters.hh.

532 {
533 FTFGamma0 = Gamma0;
534}

References FTFGamma0.

Referenced by InitForInteraction().

◆ SethNcmsEnergy()

void G4FTFParameters::SethNcmsEnergy ( const G4double  s)
inline

Definition at line 489 of file G4FTFParameters.hh.

489 {
491}

References FTFhNcmsEnergy, and S().

◆ SetInelasticCrossSection()

void G4FTFParameters::SetInelasticCrossSection ( const G4double  Xinelastic)
inline

Definition at line 503 of file G4FTFParameters.hh.

503 {
504 FTFXinelastic = Xinelastic;
505}

References FTFXinelastic.

Referenced by InitForInteraction().

◆ SetMaxNumberOfCollisions()

void G4FTFParameters::SetMaxNumberOfCollisions ( const G4double  aValue,
const G4double  bValue 
)
inline

Definition at line 604 of file G4FTFParameters.hh.

605 {
606 if ( Plab > Pbound ) {
607 MaxNumberOfCollisions = Plab/Pbound;
608 SetProbOfInteraction( -1.0 );
609 } else {
610 //MaxNumberOfCollisions = -1.0;
611 //SetProbOfInteraction( G4Exp( 0.25*(Plab-Pbound) ) );
613 SetProbOfInteraction( -1.0 );
614 }
615}
void SetProbOfInteraction(const G4double aValue)

References MaxNumberOfCollisions, and SetProbOfInteraction().

Referenced by InitForInteraction().

◆ SetMaxPt2ofNuclearDestruction()

void G4FTFParameters::SetMaxPt2ofNuclearDestruction ( const G4double  aValue)
inline

Definition at line 645 of file G4FTFParameters.hh.

645 {
647}

References MaxPt2ofNuclearDestruction.

Referenced by InitForInteraction().

◆ SetParams()

void G4FTFParameters::SetParams ( const G4int  ProcN,
const G4double  A1,
const G4double  B1,
const G4double  A2,
const G4double  B2,
const G4double  A3,
const G4double  Atop,
const G4double  Ymin 
)
inline

Definition at line 543 of file G4FTFParameters.hh.

546 {
547 ProcParams[ProcN][0] = A1; ProcParams[ProcN][1] = B1;
548 ProcParams[ProcN][2] = A2; ProcParams[ProcN][3] = B2;
549 ProcParams[ProcN][4] = A3;
550 ProcParams[ProcN][5] = Atop; ProcParams[ProcN][6] = Ymin;
551}

References ProcParams.

Referenced by InitForInteraction().

◆ SetProbabilityOfAnnihilation()

void G4FTFParameters::SetProbabilityOfAnnihilation ( const G4double  aValue)
inline

Definition at line 520 of file G4FTFParameters.hh.

520 {
522}

References ProbabilityOfAnnihilation.

Referenced by InitForInteraction().

◆ SetProbabilityOfElasticScatt() [1/2]

void G4FTFParameters::SetProbabilityOfElasticScatt ( const G4double  aValue)
inline

Definition at line 516 of file G4FTFParameters.hh.

516 {
518}

References ProbabilityOfElasticScatt.

◆ SetProbabilityOfElasticScatt() [2/2]

void G4FTFParameters::SetProbabilityOfElasticScatt ( const G4double  Xtotal,
const G4double  Xelastic 
)
inline

Definition at line 507 of file G4FTFParameters.hh.

508 {
509 if ( Xtotal == 0.0 ) {
511 } else {
512 ProbabilityOfElasticScatt = Xelastic / Xtotal;
513 }
514}

References ProbabilityOfElasticScatt.

Referenced by G4FTFModel::Init(), and InitForInteraction().

◆ SetProbLogDistr()

void G4FTFParameters::SetProbLogDistr ( const G4double  aValue)
inline

Definition at line 585 of file G4FTFParameters.hh.

585 {
586 ProbLogDistr = aValue;
587}

References ProbLogDistr.

Referenced by InitForInteraction().

◆ SetProbLogDistrPrD()

void G4FTFParameters::SetProbLogDistrPrD ( const G4double  aValue)
inline

Definition at line 581 of file G4FTFParameters.hh.

581 {
582 ProbLogDistrPrD = aValue;
583}

References ProbLogDistrPrD.

Referenced by InitForInteraction().

◆ SetProbOfInteraction()

void G4FTFParameters::SetProbOfInteraction ( const G4double  aValue)
inline

Definition at line 617 of file G4FTFParameters.hh.

617 {
618 ProbOfInelInteraction = aValue;
619}

References ProbOfInelInteraction.

Referenced by SetMaxNumberOfCollisions().

◆ SetProbOfSameQuarkExchange()

void G4FTFParameters::SetProbOfSameQuarkExchange ( const G4double  aValue)
inline

Definition at line 557 of file G4FTFParameters.hh.

557 {
559}

References ProbOfSameQuarkExchange.

Referenced by InitForInteraction().

◆ SetProjMinDiffMass()

void G4FTFParameters::SetProjMinDiffMass ( const G4double  aValue)
inline

◆ SetProjMinNonDiffMass()

void G4FTFParameters::SetProjMinNonDiffMass ( const G4double  aValue)
inline

◆ SetPt2Kink()

void G4FTFParameters::SetPt2Kink ( const G4double  aValue)
inline

Definition at line 591 of file G4FTFParameters.hh.

591 {
592 Pt2kink = aValue;
593}

References Pt2kink.

Referenced by G4FTFParameters().

◆ SetPt2ofNuclearDestruction()

void G4FTFParameters::SetPt2ofNuclearDestruction ( const G4double  aValue)
inline

Definition at line 641 of file G4FTFParameters.hh.

641 {
643}

References Pt2ofNuclearDestruction.

Referenced by InitForInteraction().

◆ SetQuarkProbabilitiesAtGluonSplitUp()

void G4FTFParameters::SetQuarkProbabilitiesAtGluonSplitUp ( const G4double  Puubar,
const G4double  Pddbar,
const G4double  Pssbar 
)
inline

Definition at line 595 of file G4FTFParameters.hh.

597 {
598 QuarkProbabilitiesAtGluonSplitUp.push_back( Puubar );
599 QuarkProbabilitiesAtGluonSplitUp.push_back( Puubar + Pddbar );
600 QuarkProbabilitiesAtGluonSplitUp.push_back( Puubar + Pddbar + Pssbar );
601}

References QuarkProbabilitiesAtGluonSplitUp.

Referenced by G4FTFParameters().

◆ SetR2ofNuclearDestruction()

void G4FTFParameters::SetR2ofNuclearDestruction ( const G4double  aValue)
inline

Definition at line 629 of file G4FTFParameters.hh.

629 {
630 R2ofNuclearDestruction = aValue;
631}

References R2ofNuclearDestruction.

Referenced by InitForInteraction().

◆ SetRadiusOfHNinteractions2()

void G4FTFParameters::SetRadiusOfHNinteractions2 ( const G4double  Radius2)
inline

Definition at line 524 of file G4FTFParameters.hh.

524 {
525 RadiusOfHNinteractions2 = Radius2;
526}

References RadiusOfHNinteractions2.

Referenced by InitForInteraction().

◆ SetSlope()

void G4FTFParameters::SetSlope ( const G4double  Slope)
inline

Definition at line 528 of file G4FTFParameters.hh.

528 {
529 FTFSlope = 12.84 / Slope; // Slope is in GeV^-2, FTFSlope in fm^-2
530}

References FTFSlope.

Referenced by InitForInteraction().

◆ SetTarMinDiffMass()

void G4FTFParameters::SetTarMinDiffMass ( const G4double  aValue)
inline

◆ SetTarMinNonDiffMass()

void G4FTFParameters::SetTarMinNonDiffMass ( const G4double  aValue)
inline

◆ SetTotalCrossSection()

void G4FTFParameters::SetTotalCrossSection ( const G4double  Xtotal)
inline

Definition at line 495 of file G4FTFParameters.hh.

495 {
496 FTFXtotal = Xtotal;
497}

References FTFXtotal.

Referenced by InitForInteraction().

Field Documentation

◆ AvaragePt2ofElasticScattering

G4double G4FTFParameters::AvaragePt2ofElasticScattering

◆ AveragePt2

G4double G4FTFParameters::AveragePt2

Definition at line 445 of file G4FTFParameters.hh.

Referenced by GetAveragePt2(), Reset(), and SetAveragePt2().

◆ CofNuclearDestruction

G4double G4FTFParameters::CofNuclearDestruction

Definition at line 457 of file G4FTFParameters.hh.

Referenced by GetCofNuclearDestruction(), Reset(), and SetCofNuclearDestruction().

◆ CofNuclearDestructionPr

G4double G4FTFParameters::CofNuclearDestructionPr

◆ csGGinstance

G4VComponentCrossSection* G4FTFParameters::csGGinstance
private

Definition at line 481 of file G4FTFParameters.hh.

Referenced by G4FTFParameters(), and InitForInteraction().

◆ DeltaProbAtQuarkExchange

G4double G4FTFParameters::DeltaProbAtQuarkExchange

◆ DofNuclearDestruction

G4double G4FTFParameters::DofNuclearDestruction

Definition at line 462 of file G4FTFParameters.hh.

Referenced by GetDofNuclearDestruction(), Reset(), and SetDofNuclearDestruction().

◆ ExcitationEnergyPerWoundedNucleon

G4double G4FTFParameters::ExcitationEnergyPerWoundedNucleon

◆ fParCollBaryonProj

G4FTFParamCollBaryonProj G4FTFParameters::fParCollBaryonProj
private

Definition at line 474 of file G4FTFParameters.hh.

Referenced by InitForInteraction().

◆ fParCollMesonProj

G4FTFParamCollMesonProj G4FTFParameters::fParCollMesonProj
private

Definition at line 477 of file G4FTFParameters.hh.

Referenced by InitForInteraction().

◆ fParCollPionProj

G4FTFParamCollPionProj G4FTFParameters::fParCollPionProj
private

Definition at line 478 of file G4FTFParameters.hh.

Referenced by InitForInteraction().

◆ FTFGamma0

G4double G4FTFParameters::FTFGamma0

Definition at line 431 of file G4FTFParameters.hh.

Referenced by GammaElastic(), Reset(), and SetGamma0().

◆ FTFhNcmsEnergy

G4double G4FTFParameters::FTFhNcmsEnergy

Definition at line 419 of file G4FTFParameters.hh.

Referenced by Reset(), and SethNcmsEnergy().

◆ FTFSlope

G4double G4FTFParameters::FTFSlope

Definition at line 429 of file G4FTFParameters.hh.

Referenced by GammaElastic(), GetSlope(), Reset(), and SetSlope().

◆ FTFXannihilation

G4double G4FTFParameters::FTFXannihilation

Definition at line 425 of file G4FTFParameters.hh.

Referenced by Reset().

◆ FTFXelastic

G4double G4FTFParameters::FTFXelastic

Definition at line 423 of file G4FTFParameters.hh.

Referenced by GetElasticCrossSection(), Reset(), and SetElastisCrossSection().

◆ FTFXinelastic

G4double G4FTFParameters::FTFXinelastic

Definition at line 424 of file G4FTFParameters.hh.

Referenced by GetInelasticCrossSection(), Reset(), and SetInelasticCrossSection().

◆ FTFXtotal

G4double G4FTFParameters::FTFXtotal

Definition at line 422 of file G4FTFParameters.hh.

Referenced by GetTotalCrossSection(), Reset(), and SetTotalCrossSection().

◆ MaxNumberOfCollisions

G4double G4FTFParameters::MaxNumberOfCollisions

Definition at line 453 of file G4FTFParameters.hh.

Referenced by GetMaxNumberOfCollisions(), Reset(), and SetMaxNumberOfCollisions().

◆ MaxPt2ofNuclearDestruction

G4double G4FTFParameters::MaxPt2ofNuclearDestruction

◆ ProbabilityOfAnnihilation

G4double G4FTFParameters::ProbabilityOfAnnihilation

◆ ProbabilityOfElasticScatt

G4double G4FTFParameters::ProbabilityOfElasticScatt

◆ ProbLogDistr

G4double G4FTFParameters::ProbLogDistr

Definition at line 446 of file G4FTFParameters.hh.

Referenced by GetProbLogDistr(), Reset(), and SetProbLogDistr().

◆ ProbLogDistrPrD

G4double G4FTFParameters::ProbLogDistrPrD

Definition at line 441 of file G4FTFParameters.hh.

Referenced by GetProbLogDistrPrD(), Reset(), and SetProbLogDistrPrD().

◆ ProbOfInelInteraction

G4double G4FTFParameters::ProbOfInelInteraction

Definition at line 454 of file G4FTFParameters.hh.

Referenced by GetProbOfInteraction(), Reset(), and SetProbOfInteraction().

◆ ProbOfSameQuarkExchange

G4double G4FTFParameters::ProbOfSameQuarkExchange

◆ ProcParams

G4double G4FTFParameters::ProcParams[5][7]

Definition at line 434 of file G4FTFParameters.hh.

Referenced by GetProcProb(), Reset(), and SetParams().

◆ ProjMinDiffMass

G4double G4FTFParameters::ProjMinDiffMass

Definition at line 439 of file G4FTFParameters.hh.

Referenced by GetProjMinDiffMass(), Reset(), and SetProjMinDiffMass().

◆ ProjMinNonDiffMass

G4double G4FTFParameters::ProjMinNonDiffMass

Definition at line 440 of file G4FTFParameters.hh.

Referenced by GetProjMinNonDiffMass(), Reset(), and SetProjMinNonDiffMass().

◆ Pt2kink

G4double G4FTFParameters::Pt2kink

Definition at line 449 of file G4FTFParameters.hh.

Referenced by GetPt2Kink(), Reset(), and SetPt2Kink().

◆ Pt2ofNuclearDestruction

G4double G4FTFParameters::Pt2ofNuclearDestruction

◆ QuarkProbabilitiesAtGluonSplitUp

std::vector< G4double > G4FTFParameters::QuarkProbabilitiesAtGluonSplitUp

◆ R2ofNuclearDestruction

G4double G4FTFParameters::R2ofNuclearDestruction

Definition at line 458 of file G4FTFParameters.hh.

Referenced by GetR2ofNuclearDestruction(), Reset(), and SetR2ofNuclearDestruction().

◆ RadiusOfHNinteractions2

G4double G4FTFParameters::RadiusOfHNinteractions2

◆ StringMass

G4LundStringFragmentation* G4FTFParameters::StringMass
private

Definition at line 467 of file G4FTFParameters.hh.

Referenced by G4FTFParameters(), GetMinMass(), and ~G4FTFParameters().

◆ TarMinDiffMass

G4double G4FTFParameters::TarMinDiffMass

Definition at line 442 of file G4FTFParameters.hh.

Referenced by GetTarMinDiffMass(), Reset(), and SetTarMinDiffMass().

◆ TarMinNonDiffMass

G4double G4FTFParameters::TarMinNonDiffMass

Definition at line 443 of file G4FTFParameters.hh.

Referenced by GetTarMinNonDiffMass(), Reset(), and SetTarMinNonDiffMass().


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