33 #define INCLXX_IN_GEANT4_MODE 1
110 if(plab > 2000.)
goto sel13;
115 sel1:
if (i == 0)
goto sel2;
116 sel3:
if (plab < 800.)
goto sel4;
117 if (plab > 2000.)
goto sel13;
118 sel=(1250./(50.+p1)-4.*std::pow(p1-1.3, 2))*scale;
121 sel4:
if (plab < 440.) {
122 sel=34.*std::pow(p1/0.4, (-2.104))*scale;
124 sel=(23.5+1000.*std::pow(p1-0.7, 4))*scale;
128 sel13: sel=77./(p1+1.5)*scale;
131 sel2:
if (plab < 800.)
goto sel11;
132 if (plab > 2000.)
goto sel13;
133 sel=31./std::sqrt(p1)*scale;
136 sel11:
if (plab < 450.) {
138 sel=6.3555*std::exp(-3.2481*alp-0.377*alp*alp)*scale;
140 sel=(33.+196.*std::sqrt(std::pow(std::abs(p1-0.95),5)))*scale;
150 const G4double momentumGeV = 0.001 * pLab;
155 if(isospin==2 || isospin==-2) {
157 xs = (41.0 + (60.0*momentumGeV - 54.0)*std::exp(-1.2*momentumGeV) - 77.0/(momentumGeV + 1.5));
158 }
else if(pLab >= 1500.0 && pLab < 2000.0) {
159 xs = (41.0 + 60.0*(momentumGeV - 0.9)*std::exp(-1.2*momentumGeV) - 1250.0/(momentumGeV+50.0)+ 4.0*std::pow(momentumGeV - 1.3, 2));
160 }
else if(pLab < 1500.0) {
161 xs = (23.5 + 24.6/(1.0 + std::exp(-10.0*momentumGeV + 12.0))
162 -1250.0/(momentumGeV +50.0)+4.0*std::pow(momentumGeV - 1.3,2));
164 }
else if(isospin==0) {
166 xs = (42.0 - 77.0/(momentumGeV + 1.5));
167 }
else if(pLab >= 1000.0 && pLab < 2000.0) {
168 xs = (24.2 + 8.9*momentumGeV - 31.1/std::sqrt(momentumGeV));
169 }
else if(pLab < 1000.0) {
170 xs = (33.0 + 196.0*std::sqrt(std::pow(std::abs(momentumGeV - 0.95),5))
171 -31.1/std::sqrt(momentumGeV));
175 if(xs < 0.0)
return 0.0;
182 return -2.33730e-06*std::pow(x, 3)+1.13819e-02*std::pow(x,2)
183 -1.83993e+01*x+9893.4;
184 }
else if(x > 1750.0 && x <= 2175.0) {
185 return 1.13531e-06*std::pow(x, 3)-6.91694e-03*std::pow(x, 2)
186 +1.39907e+01*x-9360.76;
188 return -3.18087*std::log(x)+52.9784;
195 return 0.00120683*(x-1372.52)*(x-1372.52)+26.2058;
196 }
else if(x > 1475.0 && x <= 1565.0) {
197 return 1.15873e-05*x*x+49965.6/((x-1519.59)*(x-1519.59)+2372.55);
198 }
else if(x > 1565.0 && x <= 2400.0) {
199 return 34.0248+43262.2/((x-1681.65)*(x-1681.65)+1689.35);
200 }
else if(x > 2400.0 && x <= 7500.0) {
201 return 3.3e-7*(x-7500.0)*(x-7500.0)+24.5;
221 return inelastic +
elastic(p1, p2);
234 if(x>10000.)
return 0.0;
242 }
else if(particle2->
isPion()) {
253 G4double q2=(y-1076.0*1076.0)*(y-800.0*800.0)/y/4.0;
257 G4double q3 = std::pow(std::sqrt(q2),3);
259 G4double spnResult = 326.5/(std::pow((x-1215.0-ramass)*2.0/(110.0-ramass), 2)+1.0);
260 spnResult = spnResult*(1.0-5.0*ramass/1215.0);
262 spnResult = spnResult*f3*cg/6.0;
264 if(x < 1200.0 && spnResult < 5.0) {
270 if((ind2t3 == 1 && ipit3 == 2) || (ind2t3 == -1 && ipit3 == -2))
272 else if((ind2t3 == 1 && ipit3 == -2) || (ind2t3 == -1 && ipit3 == 2))
276 INCL_ERROR(
"Unknown configuration!" << std::endl);
285 if(isospin==4 || isospin==-4)
return 0.0;
299 if(Ecm <= 938.3 + deltaMass) {
303 if(Ecm < 938.3 + deltaMass + 2.0) {
304 Ecm = 938.3 + deltaMass + 2.0;
316 result *= 3.*(32.0 + isospin * isospin * (deltaIsospin * deltaIsospin - 5))/64.0;
317 result /= 1.0 + 0.25 * isospin * isospin;
346 return 5.5e-6 * x/(7.7 + x);
348 return (5.34 + 0.67*(x - 2.0)) * 1.0e-6;
352 G4double b = (7.16 - 1.63*x) * 1.0e-6;
353 return b/(1.0 + std::exp(-(x - 0.45)/0.05));
354 }
else if(pl < 1100.0) {
355 return (9.87 - 4.88 * x) * 1.0e-6;
357 return (3.68 + 0.76*x) * 1.0e-6;
G4double getMass() const
Get the cached particle mass.
virtual G4double deltaProduction(Particle const *const p1, Particle const *const p2)
Cross section for NN->NDelta.
G4double elasticNNLegacy(Particle const *const part1, Particle const *const part2)
Internal implementation of the elastic cross section.
G4double squareTotalEnergyInCM(Particle const *const p1, Particle const *const p2)
virtual G4double calculateNNAngularSlope(G4double energyCM, G4int iso)
Calculate the slope of the NN DDXS.
virtual G4double recombination(Particle const *const p1, Particle const *const p2)
Cross section for NDelta->NN.
G4bool isDelta() const
Is it a Delta?
G4double spnPiMinusPHE(const G4double x)
G4double spnPiPlusPHE(const G4double x)
const G4double effectivePionMass
virtual G4double total(Particle const *const p1, Particle const *const p2)
Total (elastic+inelastic) particle-particle cross section.
virtual G4double pionNucleon(Particle const *const p1, Particle const *const p2)
Total (elastic+inelastic) pion-nucleon cross section.
G4int getIsospin(const ParticleType t)
Get the isospin of a particle.
const G4double effectiveNucleonMass2
G4INCL::ParticleType getType() const
Get the particle type.
G4bool isNucleon() const
Is this a nucleon?
virtual G4double elastic(Particle const *const p1, Particle const *const p2)
Elastic particle-particle cross section.
G4double momentumInLab(Particle const *const p1, Particle const *const p2)
gives the momentum in the lab frame of two particles.
G4double totalEnergyInCM(Particle const *const p1, Particle const *const p2)
G4bool isPion() const
Is this a pion?
Cross sections used in INCL4.6.
const G4double effectiveNucleonMass