Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4INCLCrossSectionsAntiparticles.cc
Go to the documentation of this file.
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26// INCL++ intra-nuclear cascade model
27// Alain Boudard, CEA-Saclay, France
28// Joseph Cugnon, University of Liege, Belgium
29// Jean-Christophe David, CEA-Saclay, France
30// Pekka Kaitaniemi, CEA-Saclay, France, and Helsinki Institute of Physics, Finland
31// Sylvie Leray, CEA-Saclay, France
32// Davide Mancusi, CEA-Saclay, France
33//
34#define INCLXX_IN_GEANT4_MODE 1
35
36#include "globals.hh"
37
38/** \file G4INCLCrossSectionsAntiparticles.cc
39 * \brief Multipion, mesonic Resonances, strange cross sections and antinucleon as projectile
40 *
41 * \date 31st March 2023
42 * \author Demid Zharenov
43 */
44
45#include "G4INCLCrossSectionsAntiparticles.hh"
46#include "G4INCLKinematicsUtils.hh"
47#include "G4INCLParticleTable.hh"
48// #include <cassert>
49
50namespace G4INCL {
51
52 template<G4int N>
53 struct BystrickyEvaluator {
54 static G4double eval(const G4double pLab, const G4double oneOverThreshold, HornerCoefficients<N> const &coeffs) {
55 const G4double pMeV = pLab*1E3;
56 const G4double ekin=std::sqrt(ParticleTable::effectiveNucleonMass2+pMeV*pMeV)-ParticleTable::effectiveNucleonMass;
57 const G4double xrat=ekin*oneOverThreshold;
58 const G4double x=std::log(xrat);
59 return HornerEvaluator<N>::eval(x, coeffs) * x * std::exp(-0.5*x);
60 }
61 };
62
63 const G4int CrossSectionsAntiparticles::nMaxPiNN = 4;
64 const G4int CrossSectionsAntiparticles::nMaxPiPiN = 4;
65
66 CrossSectionsAntiparticles::CrossSectionsAntiparticles() :
67 s11pzHC(-2.228000000000294018,8.7560000000005723725,-0.61000000000023239325,-5.4139999999999780324,3.3338333333333348023,-0.75835000000000022049,0.060623611111111114688),
68 s01ppHC(2.0570000000126518344,-6.029000000012135826,36.768500000002462784,-45.275666666666553533,25.112666666666611953,-7.2174166666666639187,1.0478875000000000275,-0.060804365079365080846),
69 s01pzHC(0.18030000000000441851,7.8700999999999953598,-4.0548999999999990425,0.555199999999999959),
70 s11pmHC(0.20590000000000031866,3.3450999999999993936,-1.4401999999999997825,0.17076666666666664973),
71 s12pmHC(-0.77235999999999901328,4.2626599999999991117,-1.9008899999999997323,0.30192266666666663379,-0.012270833333333331986),
72 s12ppHC(-0.75724999999999975664,2.0934399999999998565,-0.3803099999999999814),
73 s12zzHC(-0.89599999999996965072,7.882999999999978632,-7.1049999999999961928,1.884333333333333089),
74 s02pzHC(-1.0579999999999967036,11.113999999999994089,-8.5259999999999990196,2.0051666666666666525),
75 s02pmHC(2.4009000000012553286,-7.7680000000013376183,20.619000000000433505,-16.429666666666723928,5.2525708333333363472,-0.58969166666666670206),
76 s12mzHC(-0.21858699999999976269,1.9148999999999999722,-0.31727500000000001065,-0.027695000000000000486)
77 {
78 }
79
80 /// \brief redefining previous cross sections
81
82 G4double CrossSectionsAntiparticles::total(Particle const * const p1, Particle const * const p2) {
83 G4double inelastic;
84 if ((p1->isNucleon() && p2->isAntiNucleon()) || (p1->isAntiNucleon() && p2->isNucleon()))
85 inelastic = NNbarCEX(p1, p2) + NNbarToNNbarpi(p1, p2) + NNbarToNNbar2pi(p1, p2) + NNbarToNNbar3pi(p1, p2) + NNbarToAnnihilation(p1, p2) + NNbarToLLbar(p1, p2);
86 else if(p1->isNucleon() && p2->isNucleon()) {
87 return CrossSectionsMultiPions::NNTot(p1, p2);
88 } else if((p1->isNucleon() && p2->isDelta()) ||
89 (p1->isDelta() && p2->isNucleon())) {
90 inelastic = CrossSectionsMultiPions::NDeltaToNN(p1, p2) + NDeltaToNLK(p1, p2) + NDeltaToNSK(p1, p2) + NDeltaToDeltaLK(p1, p2) + NDeltaToDeltaSK(p1, p2) + NDeltaToNNKKb(p1, p2);
91 } else if((p1->isNucleon() && p2->isPion()) ||
92 (p1->isPion() && p2->isNucleon())) {
93 return CrossSectionsMultiPions::piNTot(p1,p2);
94 } else if((p1->isNucleon() && p2->isEta()) ||
95 (p1->isEta() && p2->isNucleon())) {
96 inelastic = CrossSectionsMultiPionsAndResonances::etaNToPiN(p1,p2) + CrossSectionsMultiPionsAndResonances::etaNToPiPiN(p1,p2);
97 } else if((p1->isNucleon() && p2->isOmega()) ||
98 (p1->isOmega() && p2->isNucleon())) {
99 inelastic = CrossSectionsMultiPionsAndResonances::omegaNInelastic(p1,p2);
100 } else if((p1->isNucleon() && p2->isEtaPrime()) ||
101 (p1->isEtaPrime() && p2->isNucleon())) {
102 inelastic = CrossSectionsMultiPionsAndResonances::etaPrimeNToPiN(p1,p2);
103 } else if((p1->isNucleon() && p2->isLambda()) ||
104 (p1->isLambda() && p2->isNucleon())) {
105 inelastic = CrossSectionsStrangeness::NLToNS(p1,p2);
106 } else if((p1->isNucleon() && p2->isSigma()) ||
107 (p1->isSigma() && p2->isNucleon())) {
108 inelastic = CrossSectionsStrangeness::NSToNL(p1,p2) + CrossSectionsStrangeness::NSToNS(p1,p2);
109 } else if((p1->isNucleon() && p2->isKaon()) ||
110 (p1->isKaon() && p2->isNucleon())) {
111 inelastic = CrossSectionsStrangeness::NKToNK(p1,p2) + CrossSectionsStrangeness::NKToNKpi(p1,p2) + CrossSectionsStrangeness::NKToNK2pi(p1,p2);
112 } else if((p1->isNucleon() && p2->isAntiKaon()) ||
113 (p1->isAntiKaon() && p2->isNucleon())) {
114 inelastic = CrossSectionsStrangeness::NKbToLpi(p1,p2)
115 + CrossSectionsStrangeness::NKbToSpi(p1,p2) + CrossSectionsStrangeness::NKbToL2pi(p1,p2)
116 + CrossSectionsStrangeness::NKbToS2pi(p1,p2) + CrossSectionsStrangeness::NKbToNKb(p1,p2)
117 + CrossSectionsStrangeness::NKbToNKbpi(p1,p2) + CrossSectionsStrangeness::NKbToNKb2pi(p1,p2);
118 } else {
119 inelastic = 0.;
120 }
121 return inelastic + elastic(p1, p2);
122 }
123
124 // without NNbar!
125 G4double CrossSectionsAntiparticles::elastic(Particle const * const p1, Particle const * const p2) {
126 if ((p1->isNucleon() && p2->isAntiNucleon()) || (p1->isAntiNucleon() && p2->isNucleon()))
127 return NNbarElastic(p1, p2);
128 if((p1->isNucleon()||p1->isDelta()) && (p2->isNucleon()||p2->isDelta())){ // N-N, N-Delta, Delta-Delta
129 return CrossSectionsMultiPions::elastic(p1, p2);
130 }
131 else if ((p1->isNucleon() && p2->isPion()) || (p2->isNucleon() && p1->isPion())){
132 return CrossSectionsMultiPions::elastic(p1, p2);
133 }
134 else if ((p1->isNucleon() && p2->isEta()) || (p2->isNucleon() && p1->isEta())){
135 return CrossSectionsMultiPionsAndResonances::etaNElastic(p1, p2);
136 }
137 else if ((p1->isNucleon() && p2->isHyperon()) || (p2->isNucleon() && p1->isHyperon())){
138 return CrossSectionsStrangeness::NYelastic(p1, p2);
139 }
140 else if ((p1->isNucleon() && p2->isKaon()) || (p2->isNucleon() && p1->isKaon())){
141 return CrossSectionsStrangeness::NKelastic(p1, p2);
142 }
143 else if ((p1->isNucleon() && p2->isAntiKaon()) || (p2->isNucleon() && p1->isAntiKaon())){
144 return CrossSectionsStrangeness::NKbelastic(p1, p2);
145 }
146 else {
147 return 0.0;
148 }
149 }
150
151 G4double CrossSectionsAntiparticles::NNbarCEX(Particle const * const p1, Particle const * const p2) {
152 //brief ppbar
153 // p pbar -> n nbar (BFMM 204)
154 //
155 //brief nnbar
156 // n nbar -> p pbar (same as BFMM 204, but no threshold)
157 //
158
159// assert((p1->isAntiNucleon() && p2->isNucleon()) || (p1->isNucleon() && p2->isAntiNucleon()));
160
161 G4double sigma=0.;
162 const G4int iso=ParticleTable::getIsospin(p1->getType()) + ParticleTable::getIsospin(p2->getType());
163 // iso == 2 || iso == -2 (n pbar or p nbar)
164
165 const std::vector<G4double> BFMM204 = {7.549, -0.041, -2.959, -6.835, 1.629, 0.114};
166 //{6.875, 0.590, -0.003, -6.629, 1.532, 0.114}
167 //const G4double Eth_PPbar_NNbar = 0.114;
168 const std::vector<G4double> BFMM204nn = {7.549, -0.041, -2.959, -6.835, 1.629};
169 //const G4double Eth_NNbar_PPbar = 0.0;
170
171 const Particle *antinucleon;
172 const Particle *nucleon;
173
174 if (p1->isAntiNucleon()) {
175 antinucleon = p1;
176 nucleon = p2;
177 }
178 else {
179 antinucleon = p2;
180 nucleon = p1;
181 }
182
183 const G4double pLab = 0.001*KinematicsUtils::momentumInLab(antinucleon, nucleon); // GeV
184
185 if(iso == 2 || iso == -2){ //npbar or pnbar
186 sigma = 0.0;
187 return sigma;
188 }
189 else{ // ppbar or nnbar
190 if(p1->getType()==antiProton || p1->getType()==Proton)
191 sigma = KinematicsUtils::compute_xs(BFMM204, pLab); // ppbar case
192 else
193 sigma = KinematicsUtils::compute_xs(BFMM204nn, pLab); // nnbar case
194 return sigma;
195 }
196 }
197
198 G4double CrossSectionsAntiparticles::NNbarElastic(Particle const * const p1, Particle const * const p2) {
199 //brief ppbar
200 // p pbar -> p pbar (BFMM 2)
201 //
202 //brief npbar
203 // n pbar -> n pbar (BFMM 472)
204 //
205 //brief nnbar
206 // n nbar -> n nbar (same as BFMM 2)
207 //
208 //brief pnbar
209 // p nbar -> p nbar (same as BFMM 472)
210 //
211
212// assert((p1->isAntiNucleon() && p2->isNucleon()) || (p1->isNucleon() && p2->isAntiNucleon()));
213
214 G4double sigma=0.;
215 const G4int iso=ParticleTable::getIsospin(p1->getType()) + ParticleTable::getIsospin(p2->getType());
216 // iso == 2 || iso == -2 (n pbar or p nbar)
217
218 const std::vector<G4double> BFMM2 = {110.496, -65.605, -0.198, -34.813, 4.317};
219 //elastic ppbar;
220 const std::vector<G4double> BFMM472 = {14.625, 23.413, -0.288, -9.002, 1.084};
221 //elastic pnbar;
222
223 const Particle *antinucleon;
224 const Particle *nucleon;
225
226 if (p1->isAntiNucleon()) {
227 antinucleon = p1;
228 nucleon = p2;
229 }
230 else {
231 antinucleon = p2;
232 nucleon = p1;
233 }
234
235 const G4double pLab = 0.001*KinematicsUtils::momentumInLab(antinucleon, nucleon); // GeV
236
237 if(iso == 2 || iso == -2){ //npbar or pnbar
238 sigma = KinematicsUtils::compute_xs(BFMM472, pLab);
239 return sigma;
240 }
241 else{ // ppbar or nnbar
242 if(p1->getType()==antiProton || p1->getType()==Proton)
243 sigma = KinematicsUtils::compute_xs(BFMM2, pLab); // ppbar case
244 else
245 sigma = KinematicsUtils::compute_xs(BFMM2, pLab); // nnbar case
246 return sigma;
247 }
248 }
249
250 G4double CrossSectionsAntiparticles::NNbarToLLbar(Particle const * const p1, Particle const * const p2) {
251 // this channel includes all states with lambdas, sigmas and xis and their antiparticles
252
253 //brief ppbar
254 // p pbar -> l lbar (BFMM 121)
255 // ppbar -> l lbar pi0 (BFMM 113)
256 // ppbar -> splus pim lbar || sminusbar pim l (BFMM 136)
257 // ppbar -> sminus pip lbar || splusbar l pip (BFMM 146)
258 // ppbar -> sp spbar (BFMM 139)
259 // ppbar -> sm smbar (BFMM 149)
260 // ppbar -> szero szerobar (BFMM 144)
261 // ppbar -> ximinus ximinusbar (BFMM 101)
262 // ppbar -> szero lbar || szerobar l (BFMM 143)
263 //
264 //
265 //brief npbar
266 // n pbar -> l lbar pi- (BFMM 487)
267 // n pbar -> l sbarplus || lbar sminus (BFMM 488)
268 //
269 //
270 //brief nnbar
271 // all same as for ppbar
272 //
273 //
274 //brief pnbar
275 // p nbar -> l lbar pi+ (same as BFMM 487)
276 // p nbar -> l sbarminus || lbar splus (same as BFMM 488)
277 //
278
279 const std::vector<G4double> BFMM121 = {2.379, -2.738, -1.260, -1.915, 0.430, 1.437};
280 //const G4double Eth_PPbar_LLbar = 1.437;
281 const std::vector<G4double> BFMM113 = {-0.105, 0.000, -5.099, 0.188, -0.050, 1.820};
282 //const G4double Eth_PPbar_LLbar_pi0 = 1.820;
283 const std::vector<G4double> BFMM139 = {0.142, -0.291, -1.702, -0.058, 0.001, 1.851};
284 //const G4double Eth_PPbar_SpSpbar = 1.851;
285 const std::vector<G4double> BFMM149 = {1.855, -2.238, -1.002, -1.279, 0.252, 1.896};
286 //const G4double Eth_PPbar_SmSmbar = 1.896;
287 const std::vector<G4double> BFMM136 = {1.749, -2.506, -1.222, -1.262, 0.274, 2.042};
288 //const G4double Eth_PPbar_SpLbar_pim = 2.042;
289 const std::vector<G4double> BFMM146 = {1.037, -1.437, -1.155, -0.709, 0.138, 2.065};
290 //const G4double Eth_PPbar_SmLbar_pip = 2.065;
291 const std::vector<G4double> BFMM143 = {0.652, -1.006, -1.805, -0.537, 0.121, 1.653};
292 //const G4double Eth_PPbar_Szero_Lbar = 1.653;
293
294// assert((p1->isAntiNucleon() && p2->isNucleon()) || (p1->isNucleon() && p2->isAntiNucleon()));
295
296 G4double sigma=0.;
297 const G4int iso=ParticleTable::getIsospin(p1->getType()) + ParticleTable::getIsospin(p2->getType());
298 // iso == 2 || iso == -2 (n pbar or p nbar)
299
300 const Particle *antinucleon;
301 const Particle *nucleon;
302
303 if (p1->isAntiNucleon()) {
304 antinucleon = p1;
305 nucleon = p2;
306 }
307 else {
308 antinucleon = p2;
309 nucleon = p1;
310 }
311
312 const G4double pLab = 0.001*KinematicsUtils::momentumInLab(antinucleon, nucleon); // GeV
313
314 //fixed due to limited data
315 G4double BFMM144;
316 if(pLab > 1.868) BFMM144 = 0.008; //sigmazero sigmazerobar
317 else BFMM144 = 0.0;
318 G4double BFMM101;
319 if(pLab > 1.868) BFMM101 = 0.002; //xizero xizerobar
320 else BFMM101 = 0.0;
321
322 // npbar cross sections (fixed due to limited data)
323 G4double BFMM487;
324 if(pLab > 2.1) BFMM487 = 0.048; //llbar piminus
325 else BFMM487 = 0.0;
326 G4double BFMM488;
327 if(pLab > 2.0) BFMM488 = 0.139; //lsigmaminus +cc
328 else BFMM488 = 0.0;
329
330 if(iso == 2 || iso == -2){ //npbar or pnbar
331 sigma = BFMM487 + BFMM488;
332 return sigma;
333 }
334 else{ // ppbar or nnbar
335 sigma = KinematicsUtils::compute_xs(BFMM113, pLab)
336 +KinematicsUtils::compute_xs(BFMM139, pLab) +KinematicsUtils::compute_xs(BFMM136, pLab)
337 +KinematicsUtils::compute_xs(BFMM146, pLab)+KinematicsUtils::compute_xs(BFMM143, pLab)
338 +KinematicsUtils::compute_xs(BFMM121, pLab)+KinematicsUtils::compute_xs(BFMM149, pLab)
339 +BFMM144 +BFMM101; // nnbar case totally same as ppbar
340 return sigma;
341 }
342 }
343
344 G4double CrossSectionsAntiparticles::NNbarToNNbarpi(Particle const * const p1, Particle const * const p2) {
345 //brief ppbar
346 // p pbar -> p pbar pi0 (BFMM 185)
347 // p pbar -> p nbar pi- (BFMM 188)
348 // p pbar -> n pbar pi+ (BFMM 199)
349 // p pbar -> n nbar pi0 (no data)
350 //
351 //brief npbar
352 // n pbar -> p pbar pi- (BFMM 491)
353 // n pbar -> p nbar pion (impossible)
354 // n pbar -> n pbar pi0 (BFMM 495)
355 // n pbar -> n nbar pi- (same as BFMM 188)
356 //
357 //brief nnbar
358 // n nbar -> n nbar pi0 (same as BFMM 185)
359 // n nbar -> p nbar pi- (same as BFMM 188)
360 // n nbar -> n pbar pi+ (same as BFMM 199)
361 // n nbar -> p pbar pi0 (no data)
362 //
363 //brief pnbar
364 // p nbar -> p pbar pi+ (same as BFMM 491)
365 // p nbar -> n pbar pion (impossible)
366 // p nbar -> p nbar pi0 (BFMM 495)
367 // p nbar -> n nbar pi- (same as BFMM 188)
368 //
369 //
370 // BFMM 188,199 are very close in value, 491 is larger
371
372// assert((p1->isAntiNucleon() && p2->isNucleon()) || (p1->isNucleon() && p2->isAntiNucleon()));
373
374 G4double sigma=0.;
375 const G4int iso=ParticleTable::getIsospin(p1->getType()) + ParticleTable::getIsospin(p2->getType());
376 // iso == 2 || iso == -2 (n pbar or p nbar)
377
378 const std::vector<G4double> BFMM185 = {-0.734, 0.841, 0.905, 3.415, -2.316, 0.775};
379 //{22.781, -22.602, -0.752, -11.036, 1.548, 0.775}
380 //const G4double Eth_PPbar_PPbar_pi0 = 0.775;
381 const std::vector<G4double> BFMM188 = { -0.442, 0.501, 0.002, 3.434, -1.201, 0.798};
382 //const G4double Eth_PPbar_PNbar_pim = 0.798;
383 const std::vector<G4double> BFMM199 = {-2.025, 2.055, -2.355, 6.064, -2.004, 0.798};
384 //const G4double Eth_PPbar_NPbar_pip = 0.798;
385 const std::vector<G4double> BFMM491 = {24.125, -20.669, -1.534, -19.573, 4.493, 0.787};
386 //const G4double Eth_NPbar_PPbar_pim = 0.787;
387 const std::vector<G4double> BFMM495 = {-0.650, -0.140, -0.058, 5.166, -1.705, 0.777};
388 //const G4double Eth_NPbar_NPbar_pi0 = 0.777;
389
390 const Particle *antinucleon;
391 const Particle *nucleon;
392
393 if (p1->isAntiNucleon()) {
394 antinucleon = p1;
395 nucleon = p2;
396 }
397 else {
398 antinucleon = p2;
399 nucleon = p1;
400 }
401
402 const G4double pLab = 0.001*KinematicsUtils::momentumInLab(antinucleon, nucleon); // GeV
403
404 if(iso == 2 || iso == -2){ //npbar or pnbar
405 sigma = KinematicsUtils::compute_xs(BFMM491, pLab) + KinematicsUtils::compute_xs(BFMM185, pLab) + KinematicsUtils::compute_xs(BFMM188, pLab);
406 return sigma;
407 }
408 else{ // ppbar or nnbar
409 sigma = KinematicsUtils::compute_xs(BFMM199, pLab) + KinematicsUtils::compute_xs(BFMM185, pLab) + KinematicsUtils::compute_xs(BFMM188, pLab);
410 return sigma;
411 }
412 }
413
414 G4double CrossSectionsAntiparticles::NNbarToNNbar2pi(Particle const * const p1, Particle const * const p2) {
415 //brief ppbar
416 // p pbar -> p pbar pi+ pi- (BFMM 167)
417 // p pbar -> p nbar pi- pi0 (same as BFMM 490)
418 // p pbar -> n pbar pi+ pi0 (same as BFMM 490)
419 // p pbar -> n nbar pi+ pi- (BFMM 198)
420 //
421 //brief npbar
422 // n pbar -> p pbar pi- pi0 (BFMM 490)
423 // n pbar -> p nbar pi- pi- (BFMM 492)
424 // n pbar -> n pbar pi+ pi- (BFMM 494)
425 // n pbar -> n nbar pi- pi0 (same as BFMM 490)
426 //
427 //brief nnbar
428 // n nbar -> n nbar pi+ pi- (same as BFMM 167)
429 // n nbar -> p nbar pi- pi0 (same as BFMM 490)
430 // n nbar -> n pbar pi+ pi0 (same as BFMM 490)
431 // n nbar -> p pbar pi+ pi- (same as BFMM 198)
432 //
433 //brief pnbar
434 // p nbar -> p pbar pi+ pi0 (same as BFMM 490)
435 // p nbar -> n pbar pi+ pi+ (same as BFMM 492)
436 // p nbar -> p nbar pi+ pi- (same as BFMM 494)
437 // p nbar -> n nbar pi+ pi0 (same as BFMM 490)
438 //
439 //
440 // BFMM 188,199 are very close in value, 491 is larger
441
442// assert((p1->isAntiNucleon() && p2->isNucleon()) || (p1->isNucleon() && p2->isAntiNucleon()));
443
444 G4double sigma=0.;
445 const G4int iso=ParticleTable::getIsospin(p1->getType()) + ParticleTable::getIsospin(p2->getType());
446 // iso == 2 || iso == -2 (n pbar or p nbar)
447
448 const std::vector<G4double> BFMM167 = {-6.885, 0.476, 1.206, 13.857, -5.728, 1.220};
449 //const G4double Eth_PPbar_PPbar_pip_pim = 1.220;
450 const std::vector<G4double> BFMM198 = {1.857, -21.213, -3.448, 0.827, -0.390, 1.231};
451 //const G4double Eth_PPbar_NNbar_pip_pim = 1.231;
452 const std::vector<G4double> BFMM490 = {-3.594, 0.811, 0.306, 5.108, -1.625, 1.201};
453 //const G4double Eth_PNbar_PPbar_pim_pi0 = 1.201;
454 const std::vector<G4double> BFMM492 = {-5.443, 7.254, -2.936, 8.441, -2.588, 1.221};
455 //const G4double Eth_PNbar_NPbar_pim_pim = 1.221;
456 const std::vector<G4double> BFMM494 = {21.688, -38.709, -2.062, -17.783, 3.895, 1.221};
457 //const G4double Eth_NPbar_NPbar_pip_pim = 1.221;
458
459 const Particle *antinucleon;
460 const Particle *nucleon;
461
462 if (p1->isAntiNucleon()) {
463 antinucleon = p1;
464 nucleon = p2;
465 }
466 else {
467 antinucleon = p2;
468 nucleon = p1;
469 }
470
471 const G4double pLab = 0.001*KinematicsUtils::momentumInLab(antinucleon, nucleon); // GeV
472
473 if(iso == 2 || iso == -2){ // pnbar or npbar
474 sigma = KinematicsUtils::compute_xs(BFMM490, pLab) + KinematicsUtils::compute_xs(BFMM490, pLab) + KinematicsUtils::compute_xs(BFMM167, pLab) + KinematicsUtils::compute_xs(BFMM198, pLab);
475 return sigma;
476 }
477 else{ // ppbar or nnbar
478 sigma = KinematicsUtils::compute_xs(BFMM490, pLab) + KinematicsUtils::compute_xs(BFMM490, pLab) + KinematicsUtils::compute_xs(BFMM492, pLab) + KinematicsUtils::compute_xs(BFMM494, pLab);
479 return sigma;
480 }
481 }
482
483 G4double CrossSectionsAntiparticles::NNbarToNNbar3pi(Particle const * const p1, Particle const * const p2) {
484 //brief ppbar
485 // p pbar -> p pbar pi+ pi- pi0 (BFMM 161)
486 // p pbar -> p nbar 2pi- pi+ (BFMM 169)
487 // p pbar -> n pbar 2pi+ pi- (BFMM 201)
488 // p pbar -> n nbar pi+ pi- pi0 (BFMM 197)
489 //
490 //brief npbar
491 // n pbar -> p pbar 2pi- pi+ (same as BFMM 169)
492 // n pbar -> p nbar 2pi- pi0 (same as BFMM 197)
493 // n pbar -> n pbar pi+ pi- pi0 (same as BFMM 161)
494 // n pbar -> n nbar 2pi- pi+ (same as BFMM 169)
495 //
496 //brief nnbar
497 // n nbar -> n nbar pi+ pi- pi0 (same as BFMM 161)
498 // n nbar -> p nbar 2pi- pi+ (same as BFMM 169)
499 // n nbar -> n pbar 2pi+ pi- (same as BFMM 201)
500 // n nbar -> p pbar pi+ pi- pi0 (same as BFMM 197)
501 //
502 //brief pnbar
503 // p nbar -> p pbar 2pi+ pi- (same as BFMM 169)
504 // p nbar -> n pbar 2pi+ pi0 (same as BFMM 197)
505 // p nbar -> p nbar pi+ pi- pi0 (same as BFMM 161)
506 // p nbar -> n nbar 2pi+ pi- (same as BFMM 169)
507 //
508
509// assert((p1->isAntiNucleon() && p2->isNucleon()) || (p1->isNucleon() && p2->isAntiNucleon()));
510
511 G4double sigma=0.;
512 const G4int iso=ParticleTable::getIsospin(p1->getType()) + ParticleTable::getIsospin(p2->getType());
513 // iso == 2 || iso == -2 (n pbar or p nbar)
514
515 const std::vector<G4double> BFMM161 = {-6.434, 1.351, -5.185, 7.754, -1.692, 1.604};
516 //const G4double Eth_PPbar_PPbar_pip_pim_pi0 = 1.604;
517 const std::vector<G4double> BFMM169 = {3.696, -5.356, -0.053, 1.941, -0.432, 1.624};
518 //const G4double Eth_PPbar_PNbar_2pim_pip = 1.624;
519 const std::vector<G4double> BFMM201 = {-1.070, -0.636, -0.009, 2.335, -0.499, 1.624};
520 //const G4double Eth_PPbar_NPbar_2pip_pim = 1.624;
521 const std::vector<G4double> BFMM197 = {1.857, -21.213, -3.448, 0.827, -0.390, 1.616};
522 //const G4double Eth_PPbar_NNbar_pip_pim_pi0 = 1.616;
523
524 const Particle *antinucleon;
525 const Particle *nucleon;
526
527 if (p1->isAntiNucleon()) {
528 antinucleon = p1;
529 nucleon = p2;
530 }
531 else {
532 antinucleon = p2;
533 nucleon = p1;
534 }
535
536 const G4double pLab = 0.001*KinematicsUtils::momentumInLab(antinucleon, nucleon); // GeV
537
538 if(iso == 2 || iso == -2){ // pnbar or npbar
539 sigma = KinematicsUtils::compute_xs(BFMM169, pLab) + KinematicsUtils::compute_xs(BFMM169, pLab) + KinematicsUtils::compute_xs(BFMM197, pLab) + KinematicsUtils::compute_xs(BFMM161, pLab);
540 return sigma;
541 }
542 else{ // ppbar or nnbar
543 sigma = KinematicsUtils::compute_xs(BFMM161, pLab) + KinematicsUtils::compute_xs(BFMM169, pLab) + KinematicsUtils::compute_xs(BFMM197, pLab) + KinematicsUtils::compute_xs(BFMM201, pLab);
544 return sigma;
545 }
546 }
547
548 G4double CrossSectionsAntiparticles::NNbarToAnnihilation(Particle const * const p1, Particle const * const p2) {
549 //brief ppbar
550 /*
551 This part only contains total annihilation xs, the choice of a particular final state
552 will be done in the channel file.
553 As long as we only have good data for ppbar, we assume that for npbar, pnbar and nnbar the xs
554 will be the same, but in order to compensate for the Coulombic effect the ppbar annihilation xs
555 is multiplied by the pnbar total xs and divided by the ppbar total xs.
556 */
557
558// assert((p1->isAntiNucleon() && p2->isNucleon()) || (p1->isNucleon() && p2->isAntiNucleon()));
559
560 G4double sigma=0.;
561 const G4int iso=ParticleTable::getIsospin(p1->getType()) + ParticleTable::getIsospin(p2->getType());
562 // iso == 2 || iso == -2 (n pbar or p nbar)
563
564 const std::vector<G4double> BFMM6 = {66.098, 0.153, -4.576, -38.319, 6.625}; //ppbar annihilation xs
565 const std::vector<G4double> BFMM1 = {119.066, 6.251, -0.006, -60.046, 11.958}; //ppbar total xs
566 const std::vector<G4double> BFMM471 = {108.104, 15.708, 0.832, -54.632, -6.958}; //npbar total xs
567
568 const Particle *antinucleon;
569 const Particle *nucleon;
570
571 if (p1->isAntiNucleon()) {
572 antinucleon = p1;
573 nucleon = p2;
574 }
575 else {
576 antinucleon = p2;
577 nucleon = p1;
578 }
579
580 const G4double pLab = 0.001*KinematicsUtils::momentumInLab(antinucleon, nucleon); // GeV
581
582 if(iso == 2 || iso == -2){ // pnbar or npbar
583 sigma = KinematicsUtils::compute_xs(BFMM6, pLab)*KinematicsUtils::compute_xs(BFMM471, pLab)/KinematicsUtils::compute_xs(BFMM1, pLab);
584 return sigma;
585 }
586 else if(p1->getType()==antiProton || p2->getType()==Proton){ // ppbar case
587 sigma = KinematicsUtils::compute_xs(BFMM6, pLab);
588 return sigma;
589 }
590 else{ // nnbar case
591 sigma = KinematicsUtils::compute_xs(BFMM6, pLab)*KinematicsUtils::compute_xs(BFMM471, pLab)/KinematicsUtils::compute_xs(BFMM1, pLab);
592 return sigma;
593 }
594 }
595
596} // namespace G4INCL
597
G4INCLCrossSectionsAntiparticles.hh
Multipion, mesonic Resonances, strange cross sections and antinucleon as projectile.
G4double
double G4double
Definition G4Types.hh:83
G4int
int G4int
Definition G4Types.hh:85
G4INCL::CrossSectionsAntiparticles::NNbarElastic
virtual G4double NNbarElastic(Particle const *const p1, Particle const *const p2)
Nucleon-AntiNucleon to Nucleon-AntiNucleon cross sections.
Definition G4INCLCrossSectionsAntiparticles.cc:198
G4INCL::CrossSectionsAntiparticles::NNbarCEX
virtual G4double NNbarCEX(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsAntiparticles.cc:151
G4INCL::CrossSectionsAntiparticles::NNbarToAnnihilation
virtual G4double NNbarToAnnihilation(Particle const *const p1, Particle const *const p2)
Nucleon-AntiNucleon total annihilation cross sections.
Definition G4INCLCrossSectionsAntiparticles.cc:548
G4INCL::CrossSectionsAntiparticles::elastic
virtual G4double elastic(Particle const *const p1, Particle const *const p2)
old elastic particle-particle cross section
Definition G4INCLCrossSectionsAntiparticles.cc:125
G4INCL::CrossSectionsAntiparticles::NNbarToLLbar
virtual G4double NNbarToLLbar(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsAntiparticles.cc:250
G4INCL::CrossSectionsAntiparticles::NNbarToNNbarpi
virtual G4double NNbarToNNbarpi(Particle const *const p1, Particle const *const p2)
Nucleon-AntiNucleon to Nucleon-AntiNucleon + pions cross sections.
Definition G4INCLCrossSectionsAntiparticles.cc:344
G4INCL::CrossSectionsAntiparticles::nMaxPiNN
static const G4int nMaxPiNN
Maximum number of outgoing pions in NN collisions.
Definition G4INCLCrossSectionsAntiparticles.hh:82
G4INCL::CrossSectionsAntiparticles::NNbarToNNbar3pi
virtual G4double NNbarToNNbar3pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsAntiparticles.cc:483
G4INCL::CrossSectionsAntiparticles::NNbarToNNbar2pi
virtual G4double NNbarToNNbar2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsAntiparticles.cc:414
G4INCL::CrossSectionsAntiparticles::nMaxPiPiN
static const G4int nMaxPiPiN
Maximum number of outgoing pions in piN collisions.
Definition G4INCLCrossSectionsAntiparticles.hh:85
G4INCL::CrossSectionsAntiparticles::total
virtual G4double total(Particle const *const p1, Particle const *const p2)
second new total particle-particle cross section
Definition G4INCLCrossSectionsAntiparticles.cc:82
G4INCL::CrossSectionsMultiPionsAndResonances::etaNElastic
virtual G4double etaNElastic(Particle const *const p1, Particle const *const p2)
Cross sections for mesonic resonance absorption on nucleon - elastic Channel.
Definition G4INCLCrossSectionsMultiPionsAndResonances.cc:340
G4INCL::CrossSectionsMultiPionsAndResonances::etaNToPiPiN
virtual G4double etaNToPiPiN(Particle const *const p1, Particle const *const p2)
Cross sections for mesonic resonance absorption on nucleon - pipiN Channel.
Definition G4INCLCrossSectionsMultiPionsAndResonances.cc:303
G4INCL::CrossSectionsMultiPionsAndResonances::omegaNInelastic
virtual G4double omegaNInelastic(Particle const *const p1, Particle const *const p2)
Cross sections for mesonic resonance absorption on nucleon - inelastic Channel.
Definition G4INCLCrossSectionsMultiPionsAndResonances.cc:377
G4INCL::CrossSectionsMultiPionsAndResonances::etaPrimeNToPiN
virtual G4double etaPrimeNToPiN(Particle const *const p1, Particle const *const p2)
Cross section for EtaPrimeN->PiN.
Definition G4INCLCrossSectionsMultiPionsAndResonances.cc:493
G4INCL::CrossSectionsMultiPionsAndResonances::etaNToPiN
virtual G4double etaNToPiN(Particle const *const p1, Particle const *const p2)
Cross sections for mesonic resonance absorption on nucleon - piN Channel.
Definition G4INCLCrossSectionsMultiPionsAndResonances.cc:256
G4INCL::CrossSectionsMultiPions::piNTot
G4double piNTot(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsMultiPions.cc:712
G4INCL::CrossSectionsMultiPions::elastic
virtual G4double elastic(Particle const *const p1, Particle const *const p2)
Elastic particle-particle cross section.
Definition G4INCLCrossSectionsMultiPions.cc:807
G4INCL::CrossSectionsMultiPions::NNTot
G4double NNTot(Particle const *const part1, Particle const *const part2)
Internal implementation of the NN total cross section.
Definition G4INCLCrossSectionsMultiPions.cc:164
G4INCL::CrossSectionsMultiPions::NDeltaToNN
virtual G4double NDeltaToNN(Particle const *const p1, Particle const *const p2)
Cross section for NDelta->NN.
Definition G4INCLCrossSectionsMultiPions.cc:749
G4INCL::CrossSectionsStrangeness::NKbToNKb2pi
virtual G4double NKbToNKb2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1957
G4INCL::CrossSectionsStrangeness::NKbToNKbpi
virtual G4double NKbToNKbpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1911
G4INCL::CrossSectionsStrangeness::NKbToL2pi
virtual G4double NKbToL2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1865
G4INCL::CrossSectionsStrangeness::NKbelastic
virtual G4double NKbelastic(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:372
G4INCL::CrossSectionsStrangeness::NLToNS
virtual G4double NLToNS(Particle const *const p1, Particle const *const p2)
Nucleon-Hyperon quasi-elastic cross sections.
Definition G4INCLCrossSectionsStrangeness.cc:1423
G4INCL::CrossSectionsStrangeness::NKToNK2pi
virtual G4double NKToNK2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1607
G4INCL::CrossSectionsStrangeness::NDeltaToDeltaSK
virtual G4double NDeltaToDeltaSK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:875
G4INCL::CrossSectionsStrangeness::NKToNKpi
virtual G4double NKToNKpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1563
G4INCL::CrossSectionsStrangeness::NDeltaToNSK
virtual G4double NDeltaToNSK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:778
G4INCL::CrossSectionsStrangeness::NKelastic
virtual G4double NKelastic(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:335
G4INCL::CrossSectionsStrangeness::NKbToLpi
virtual G4double NKbToLpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1764
G4INCL::CrossSectionsStrangeness::NSToNS
virtual G4double NSToNS(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1497
G4INCL::CrossSectionsStrangeness::NSToNL
virtual G4double NSToNL(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1458
G4INCL::CrossSectionsStrangeness::NKbToS2pi
virtual G4double NKbToS2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1813
G4INCL::CrossSectionsStrangeness::NKbToSpi
virtual G4double NKbToSpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:1715
G4INCL::CrossSectionsStrangeness::NKToNK
virtual G4double NKToNK(Particle const *const p1, Particle const *const p2)
Nucleon-Kaon cross sections.
Definition G4INCLCrossSectionsStrangeness.cc:1526
G4INCL::CrossSectionsStrangeness::NDeltaToNNKKb
virtual G4double NDeltaToNNKKb(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:937
G4INCL::CrossSectionsStrangeness::NDeltaToDeltaLK
virtual G4double NDeltaToDeltaLK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSectionsStrangeness.cc:831
G4INCL::CrossSectionsStrangeness::NYelastic
virtual G4double NYelastic(Particle const *const p1, Particle const *const p2)
elastic scattering for Nucleon-Strange Particles cross sections
Definition G4INCLCrossSectionsStrangeness.cc:300
G4INCL::CrossSectionsStrangeness::NDeltaToNLK
virtual G4double NDeltaToNLK(Particle const *const p1, Particle const *const p2)
Nucleon-Delta to Stange particles cross sections.
Definition G4INCLCrossSectionsStrangeness.cc:732
G4INCL::CrossSectionsStrangeness::NKbToNKb
virtual G4double NKbToNKb(Particle const *const p1, Particle const *const p2)
Nucleon-antiKaon cross sections.
Definition G4INCLCrossSectionsStrangeness.cc:1663
G4INCL::Particle::isLambda
G4bool isLambda() const
Is this a Lambda?
Definition G4INCLParticle.hh:434
G4INCL::Particle::isEtaPrime
G4bool isEtaPrime() const
Is this an etaprime?
Definition G4INCLParticle.hh:411
G4INCL::Particle::isOmega
G4bool isOmega() const
Is this an omega?
Definition G4INCLParticle.hh:408
G4INCL::Particle::isSigma
G4bool isSigma() const
Is this a Sigma?
Definition G4INCLParticle.hh:425
G4INCL::Particle::isAntiNucleon
G4bool isAntiNucleon() const
Is this an antinucleon?
Definition G4INCLParticle.hh:455
G4INCL::Particle::isEta
G4bool isEta() const
Is this an eta?
Definition G4INCLParticle.hh:405
G4INCL::Particle::isPion
G4bool isPion() const
Is this a pion?
Definition G4INCLParticle.hh:402
G4INCL::Particle::isAntiKaon
G4bool isAntiKaon() const
Is this an antiKaon?
Definition G4INCLParticle.hh:431
G4INCL::Particle::getType
G4INCL::ParticleType getType() const
Definition G4INCLParticle.hh:187
G4INCL::Particle::isKaon
G4bool isKaon() const
Is this a Kaon?
Definition G4INCLParticle.hh:428
G4INCL::Particle::isDelta
G4bool isDelta() const
Is it a Delta?
Definition G4INCLParticle.hh:420
G4INCL::Particle::isHyperon
G4bool isHyperon() const
Is this an Hyperon?
Definition G4INCLParticle.hh:440
G4INCL::Particle::isNucleon
G4bool isNucleon() const
Definition G4INCLParticle.hh:362
G4INCL::KinematicsUtils::compute_xs
G4double compute_xs(const std::vector< G4double > coefficients, const G4double pLab)
Definition G4INCLKinematicsUtils.cc:49
G4INCL::KinematicsUtils::momentumInLab
G4double momentumInLab(Particle const *const p1, Particle const *const p2)
gives the momentum in the lab frame of two particles.
Definition G4INCLKinematicsUtils.cc:162
G4INCL::ParticleTable::effectiveNucleonMass2
const G4double effectiveNucleonMass2
Definition G4INCLParticleTable.hh:70
G4INCL::ParticleTable::getIsospin
G4int getIsospin(const ParticleType t)
Get the isospin of a particle.
Definition G4INCLParticleTable.cc:531
G4INCL::ParticleTable::effectiveNucleonMass
const G4double effectiveNucleonMass
Definition G4INCLParticleTable.hh:69
G4INCL::BystrickyEvaluator::eval
static G4double eval(const G4double pLab, const G4double oneOverThreshold, HornerCoefficients< N > const &coeffs)
Definition G4INCLCrossSectionsAntiparticles.cc:54
G4INCL::HornerEvaluator::eval
static G4double eval(const G4double x, HornerCoefficients< N > const &coeffs)
Definition G4INCLHornerFormEvaluator.hh:182