Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4INCLNNToMultiPionsChannel.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#include "G4INCLNNToMultiPionsChannel.hh"
39#include "G4INCLKinematicsUtils.hh"
40#include "G4INCLBinaryCollisionAvatar.hh"
41#include "G4INCLRandom.hh"
42#include "G4INCLGlobals.hh"
43#include "G4INCLLogger.hh"
44#include <algorithm>
45#include "G4INCLPhaseSpaceGenerator.hh"
46
47namespace G4INCL {
48
49 const G4double NNToMultiPionsChannel::angularSlope = 6.;
50
51 NNToMultiPionsChannel::NNToMultiPionsChannel(const G4int npi, Particle *p1, Particle *p2)
52 : npion(npi),
53 iso1(0),
54 iso2(0),
55 particle1(p1),
56 particle2(p2)
57 {
58 std::fill(isosp, isosp+4, 0);
59 }
60
61 NNToMultiPionsChannel::~NNToMultiPionsChannel(){
62
63 }
64
65 void NNToMultiPionsChannel::fillFinalState(FinalState *fs) {
66// assert(npion > 0 && npion < 5);
67
68 iso1=ParticleTable::getIsospin(particle1->getType());
69 iso2=ParticleTable::getIsospin(particle2->getType());
70
71 ParticleList list;
72 list.push_back(particle1);
73 list.push_back(particle2);
74 fs->addModifiedParticle(particle1);
75 fs->addModifiedParticle(particle2);
76
77 isospinRepartition();
78
79 const ParticleType tn1=ParticleTable::getNucleonType(iso1);
80 particle1->setType(tn1);
81 const ParticleType tn2=ParticleTable::getNucleonType(iso2);
82 particle2->setType(tn2);
83 const ThreeVector &rcolnucleon1 = particle1->getPosition();
84 const ThreeVector &rcolnucleon2 = particle2->getPosition();
85 const ThreeVector rcol = (rcolnucleon1+rcolnucleon2)*0.5;
86 const ThreeVector zero;
87 for(G4int i=0; i<npion; ++i) {
88 const ParticleType pionType=ParticleTable::getPionType(isosp[i]);
89 Particle *pion = new Particle(pionType,zero,rcol);
90 list.push_back(pion);
91 fs->addCreatedParticle(pion);
92 }
93
94 const G4double sqrtS = KinematicsUtils::totalEnergyInCM(particle1, particle2);
95 G4int biasIndex = ((Random::shoot()<0.5) ? 0 : 1);
96 PhaseSpaceGenerator::generateBiased(sqrtS, list, biasIndex, angularSlope);
97
98 }
99
100 void NNToMultiPionsChannel::isospinRepartition() {
101 const G4double rjcd=Random::shoot();
102 G4double p;
103 const G4int itot=iso1+iso2;
104
105 if (npion == 1) {
106 p=3.*rjcd;
107 if (p < 1.) pn_ppPim();
108 else if (p < 2.) pn_pnPi0();
109 else pn_nnPip();
110 }
111 else if (npion == 2) {
112 if (itot == 2) {
113 p=20.*rjcd;
114 if (p >= 14.) pp_nnPipPip();
115 else if (p >= 11.) pp_pnPipPi0();
116 else if (p >= 7.) pp_ppPi0Pi0();
117 else pp_ppPipPim();
118 }
119 else if (itot == -2) {
120 p=20.*rjcd;
121 if (p >= 14.) nn_ppPimPim();
122 else if (p >= 11.) nn_pnPimPi0();
123 else if (p >= 7.) nn_nnPi0Pi0();
124 else nn_nnPipPim();
125 }
126 else {
127 G4double pp=Random::shoot();
128 if (pp > 0.5) {
129 p=3.*rjcd;
130 if (p < 2.) {
131 pn_pnPipPim();
132 }
133 else {
134 pn_pnPi0Pi0();
135 }
136 }
137 else {
138 p=60.*rjcd;
139 if (p >= 51.) pn_nnPipPi0();
140 else if (p >= 33.) pn_pnPi0Pi0();
141 else if (p >= 9.) pn_pnPipPim();
142 else pn_ppPimPi0();
143 }
144 }
145 }
146 else if (npion == 3) {
147 p=60.*rjcd;
148 if (itot == 2) {
149 if (p >= 42.) pp_nnPipPipPi0();
150 else if (p >= 39.) pp_pnPipPi0Pi0();
151 else if (p >= 33.) pp_pnPipPipPim();
152 else if (p >= 22.) pp_ppPi0Pi0Pi0();
153 else pp_ppPipPimPi0();
154 }
155 else if (itot == -2) {
156 if (p >= 42.) nn_ppPimPimPi0();
157 else if (p >= 39.) nn_pnPimPi0Pi0();
158 else if (p >= 33.) nn_pnPipPimPim();
159 else if (p >= 22.) nn_nnPi0Pi0Pi0();
160 else nn_nnPipPimPi0();
161 }
162 else {
163 if (p >= 57.) pn_nnPipPi0Pi0();
164 else if (p >= 51.) pn_nnPipPipPim();
165 else if (p >= 37.) pn_pnPi0Pi0Pi0();
166 else if (p >= 9.) pn_pnPi0PipPim();
167 else if (p >= 6.) pn_ppPimPi0Pi0();
168 else pn_ppPimPimPip();
169
170 }
171 }
172 else if (npion == 4) {
173 p=60.*rjcd;
174 if (itot == 2) {
175 if (p >= 48.) pp_nnPipPipPipPim();
176 else if (p >= 42.) pp_nnPipPipPi0Pi0();
177 else if (p >= 36.) pp_pnPipPipPi0Pim();
178 else if (p >= 33.) pp_pnPipPi0Pi0Pi0();
179 else if (p >= 19.) pp_ppPipPipPimPim();
180 else if (p >= 4.) pp_ppPipPi0Pi0Pim();
181 else pp_ppPi0Pi0Pi0Pi0();
182 }
183 else if (itot == -2) {
184 if (p >= 48.) nn_ppPipPimPimPim();
185 else if (p >= 42.) nn_ppPi0Pi0PimPim();
186 else if (p >= 36.) nn_pnPipPi0PimPim();
187 else if (p >= 33.) nn_pnPi0Pi0Pi0Pim();
188 else if (p >= 19.) nn_nnPipPipPimPim();
189 else if (p >= 4.) nn_nnPipPi0Pi0Pim();
190 else nn_nnPi0Pi0Pi0Pi0();
191 }
192 else {
193 G4double pp=Random::shoot();
194 if (pp > 0.5) {
195 p=9.*rjcd;
196 if (p < 1.) pn_pnPi0Pi0Pi0Pi0();
197 else if (p < 5.) pn_pnPipPi0Pi0Pim();
198 else pn_pnPipPipPimPim();
199 }
200 else {
201 if (p < 3.) pn_ppPi0Pi0Pi0Pim();
202 else if (p < 9.) pn_ppPipPi0PimPim();
203 else if (p < 15.) pn_pnPi0Pi0Pi0Pi0();
204 else if (p < 35.) pn_pnPipPi0Pi0Pim();
205 else if (p < 51.) pn_pnPipPipPimPim();
206 else if (p < 54.) pn_nnPipPi0Pi0Pi0();
207 else pn_nnPipPipPi0Pim();
208 }
209 }
210 }
211
212 std::shuffle(isosp,isosp+npion,Random::getAdapter());
213 inter2Part(0.5);
214 }
215
216
217 void NNToMultiPionsChannel::pn_ppPim() {
218 isosp[0]=-2;
219 iso1=1;
220 iso2=1;
221 }
222 void NNToMultiPionsChannel::pn_pnPi0() {
223 isosp[0]=0;
224 }
225 void NNToMultiPionsChannel::pn_nnPip() {
226 isosp[0]=2;
227 iso1=-1;
228 iso2=-1;
229 }
230 void NNToMultiPionsChannel::pp_nnPipPip() {
231 isosp[0]=2;
232 isosp[1]=2;
233 iso1=-1;
234 iso2=-1;
235 }
236 void NNToMultiPionsChannel::nn_ppPimPim() {
237 isosp[0]=-2;
238 isosp[1]=-2;
239 iso1=1;
240 iso2=1;
241 }
242 void NNToMultiPionsChannel::pn_pnPipPim() {
243 isosp[0]=2;
244 isosp[1]=-2;
245 }
246 void NNToMultiPionsChannel::pn_pnPi0Pi0() {
247 isosp[0]=0;
248 isosp[1]=0;
249 }
250 void NNToMultiPionsChannel::pp_ppPipPim() {
251 isosp[0]=2;
252 isosp[1]=-2;
253 }
254 void NNToMultiPionsChannel::nn_nnPipPim() {
255 isosp[0]=2;
256 isosp[1]=-2;
257 }
258 void NNToMultiPionsChannel::pp_ppPi0Pi0() {
259 isosp[0]=0;
260 isosp[1]=0;
261 }
262 void NNToMultiPionsChannel::nn_nnPi0Pi0() {
263 isosp[0]=0;
264 isosp[1]=0;
265 }
266 void NNToMultiPionsChannel::pp_pnPipPi0() {
267 isosp[0]=2;
268 isosp[1]=0;
269 iso1=1;
270 iso2=-1;
271 }
272 void NNToMultiPionsChannel::pn_ppPimPi0() {
273 isosp[0]=-2;
274 isosp[1]=0;
275 iso1=1;
276 iso2=1;
277 }
278 void NNToMultiPionsChannel::pn_nnPipPi0() {
279 isosp[0]=2;
280 isosp[1]=0;
281 iso1=-1;
282 iso2=-1;
283 }
284 void NNToMultiPionsChannel::nn_pnPimPi0() {
285 isosp[0]=-2;
286 isosp[1]=0;
287 iso1=1;
288 iso2=-1;
289 }
290 void NNToMultiPionsChannel::pp_pnPipPi0Pi0() {
291 isosp[0]=2;
292 isosp[1]=0;
293 isosp[2]=0;
294 iso1=1;
295 iso2=-1;
296 }
297 void NNToMultiPionsChannel::nn_pnPimPi0Pi0() {
298 isosp[0]=-2;
299 isosp[1]=0;
300 isosp[2]=0;
301 iso1=1;
302 iso2=-1;
303 }
304 void NNToMultiPionsChannel::pn_nnPipPi0Pi0() {
305 isosp[0]=2;
306 isosp[1]=0;
307 isosp[2]=0;
308 iso1=-1;
309 iso2=-1;
310 }
311 void NNToMultiPionsChannel::pp_ppPipPimPi0() {
312 isosp[0]=2;
313 isosp[1]=-2;
314 isosp[2]=0;
315 }
316 void NNToMultiPionsChannel::nn_nnPipPimPi0() {
317 isosp[0]=2;
318 isosp[1]=-2;
319 isosp[2]=0;
320 }
321 void NNToMultiPionsChannel::pp_ppPi0Pi0Pi0() {
322 isosp[0]=0;
323 isosp[1]=0;
324 isosp[2]=0;
325 }
326 void NNToMultiPionsChannel::nn_nnPi0Pi0Pi0() {
327 isosp[0]=0;
328 isosp[1]=0;
329 isosp[2]=0;
330 }
331 void NNToMultiPionsChannel::pp_pnPipPipPim() {
332 isosp[0]=2;
333 isosp[1]=2;
334 isosp[2]=-2;
335 iso1=1;
336 iso2=-1;
337 }
338 void NNToMultiPionsChannel::pp_nnPipPipPi0() {
339 isosp[0]=2;
340 isosp[1]=2;
341 isosp[2]=0;
342 iso1=-1;
343 iso2=-1;
344 }
345 void NNToMultiPionsChannel::pn_ppPimPi0Pi0() {
346 isosp[0]=-2;
347 isosp[1]=0;
348 isosp[2]=0;
349 iso1=1;
350 iso2=1;
351 }
352 void NNToMultiPionsChannel::pn_ppPimPimPip() {
353 isosp[0]=-2;
354 isosp[1]=-2;
355 isosp[2]=2;
356 iso1=1;
357 iso2=1;
358 }
359 void NNToMultiPionsChannel::pn_pnPi0PipPim() {
360 isosp[0]=0;
361 isosp[1]=2;
362 isosp[2]=-2;
363 }
364 void NNToMultiPionsChannel::pn_pnPi0Pi0Pi0() {
365 isosp[0]=0;
366 isosp[1]=0;
367 isosp[2]=0;
368 }
369 void NNToMultiPionsChannel::pn_nnPipPipPim() {
370 isosp[0]=2;
371 isosp[1]=2;
372 isosp[2]=-2;
373 iso1=-1;
374 iso2=-1;
375 }
376 void NNToMultiPionsChannel::nn_pnPipPimPim() {
377 isosp[0]=2;
378 isosp[1]=-2;
379 isosp[2]=-2;
380 iso1=1;
381 iso2=-1;
382 }
383 void NNToMultiPionsChannel::nn_ppPimPimPi0() {
384 isosp[0]=-2;
385 isosp[1]=-2;
386 isosp[2]=0;
387 iso1=1;
388 iso2=1;
389 }
390 void NNToMultiPionsChannel::pp_nnPipPipPi0Pi0() {
391 isosp[0]=2;
392 isosp[1]=2;
393 isosp[2]=0;
394 isosp[3]=0;
395 iso1=-1;
396 iso2=-1;
397 }
398 void NNToMultiPionsChannel::pp_nnPipPipPipPim() {
399 isosp[0]=2;
400 isosp[1]=2;
401 isosp[2]=2;
402 isosp[3]=-2;
403 iso1=-1;
404 iso2=-1;
405 }
406 void NNToMultiPionsChannel::nn_ppPi0Pi0PimPim() {
407 isosp[0]=0;
408 isosp[1]=0;
409 isosp[2]=-2;
410 isosp[3]=-2;
411 iso1=1;
412 iso2=1;
413 }
414 void NNToMultiPionsChannel::nn_ppPipPimPimPim() {
415 isosp[0]=2;
416 isosp[1]=-2;
417 isosp[2]=-2;
418 isosp[3]=-2;
419 iso1=1;
420 iso2=1;
421 }
422 void NNToMultiPionsChannel::pp_ppPi0Pi0Pi0Pi0() {
423 isosp[0]=0;
424 isosp[1]=0;
425 isosp[2]=0;
426 isosp[3]=0;
427 }
428 void NNToMultiPionsChannel::nn_nnPi0Pi0Pi0Pi0() {
429 isosp[0]=0;
430 isosp[1]=0;
431 isosp[2]=0;
432 isosp[3]=0;
433 }
434 void NNToMultiPionsChannel::pn_pnPi0Pi0Pi0Pi0() {
435 isosp[0]=0;
436 isosp[1]=0;
437 isosp[2]=0;
438 isosp[3]=0;
439 }
440 void NNToMultiPionsChannel::pp_ppPipPi0Pi0Pim() {
441 isosp[0]=2;
442 isosp[1]=0;
443 isosp[2]=0;
444 isosp[3]=-2;
445 }
446 void NNToMultiPionsChannel::nn_nnPipPi0Pi0Pim() {
447 isosp[0]=2;
448 isosp[1]=0;
449 isosp[2]=0;
450 isosp[3]=-2;
451 }
452 void NNToMultiPionsChannel::pn_pnPipPi0Pi0Pim() {
453 isosp[0]=2;
454 isosp[1]=0;
455 isosp[2]=0;
456 isosp[3]=-2;
457 }
458 void NNToMultiPionsChannel::pp_ppPipPipPimPim() {
459 isosp[0]=2;
460 isosp[1]=2;
461 isosp[2]=-2;
462 isosp[3]=-2;
463 }
464 void NNToMultiPionsChannel::nn_nnPipPipPimPim() {
465 isosp[0]=2;
466 isosp[1]=2;
467 isosp[2]=-2;
468 isosp[3]=-2;
469 }
470 void NNToMultiPionsChannel::pn_pnPipPipPimPim() {
471 isosp[0]=2;
472 isosp[1]=2;
473 isosp[2]=-2;
474 isosp[3]=-2;
475 }
476 void NNToMultiPionsChannel::pp_pnPipPi0Pi0Pi0() {
477 isosp[0]=2;
478 isosp[1]=0;
479 isosp[2]=0;
480 isosp[3]=0;
481 iso1=1;
482 iso2=-1;
483 }
484 void NNToMultiPionsChannel::pn_nnPipPi0Pi0Pi0() {
485 isosp[0]=2;
486 isosp[1]=0;
487 isosp[2]=0;
488 isosp[3]=0;
489 iso1=-1;
490 iso2=-1;
491 }
492 void NNToMultiPionsChannel::pp_nnPipPi0Pi0Pi0() {
493 isosp[0]=2;
494 isosp[1]=0;
495 isosp[2]=0;
496 isosp[3]=0;
497 iso1=-1;
498 iso2=-1;
499 }
500 void NNToMultiPionsChannel::pp_pnPipPipPi0Pim() {
501 isosp[0]=2;
502 isosp[1]=2;
503 isosp[2]=0;
504 isosp[3]=-2;
505 iso1=1;
506 iso2=-1;
507 }
508 void NNToMultiPionsChannel::pn_nnPipPipPi0Pim() {
509 isosp[0]=2;
510 isosp[1]=2;
511 isosp[2]=0;
512 isosp[3]=-2;
513 iso1=-1;
514 iso2=-1;
515 }
516 void NNToMultiPionsChannel::pp_nnPipPipPi0Pim() {
517 isosp[0]=2;
518 isosp[1]=2;
519 isosp[2]=0;
520 isosp[3]=-2;
521 iso1=-1;
522 iso2=-1;
523 }
524 void NNToMultiPionsChannel::nn_pnPi0Pi0Pi0Pim() {
525 isosp[0]=0;
526 isosp[1]=0;
527 isosp[2]=0;
528 isosp[3]=-2;
529 iso1=1;
530 iso2=-1;
531 }
532 void NNToMultiPionsChannel::pn_ppPi0Pi0Pi0Pim() {
533 isosp[0]=0;
534 isosp[1]=0;
535 isosp[2]=0;
536 isosp[3]=-2;
537 iso1=1;
538 iso2=1;
539 }
540 void NNToMultiPionsChannel::nn_pnPipPi0PimPim() {
541 isosp[0]=2;
542 isosp[1]=0;
543 isosp[2]=-2;
544 isosp[3]=-2;
545 iso1=1;
546 iso2=-1;
547 }
548 void NNToMultiPionsChannel::pn_ppPipPi0PimPim() {
549 isosp[0]=2;
550 isosp[1]=0;
551 isosp[2]=-2;
552 isosp[3]=-2;
553 iso1=1;
554 iso2=1;
555 }
556
557 void NNToMultiPionsChannel::inter2Part(const G4double p) {
558
559 if (Random::shoot() < p) std::swap(iso1,iso2);
560
561 }
562
563
564}
G4double
double G4double
Definition: G4Types.hh:83
G4int
int G4int
Definition: G4Types.hh:85
G4INCL::FinalState::addModifiedParticle
void addModifiedParticle(Particle *p)
Definition: G4INCLFinalState.cc:60
G4INCL::FinalState::addCreatedParticle
void addCreatedParticle(Particle *p)
Definition: G4INCLFinalState.cc:75
G4INCL::NNToMultiPionsChannel::NNToMultiPionsChannel
NNToMultiPionsChannel(const G4int, Particle *, Particle *)
Definition: G4INCLNNToMultiPionsChannel.cc:51
G4INCL::Particle::getPosition
const G4INCL::ThreeVector & getPosition() const
Definition: G4INCLParticle.hh:820
G4INCL::Particle::getType
G4INCL::ParticleType getType() const
Definition: G4INCLParticle.hh:182
G4INCL::Particle::setType
void setType(ParticleType t)
Definition: G4INCLParticle.hh:191
G4INCL::KinematicsUtils::totalEnergyInCM
G4double totalEnergyInCM(Particle const *const p1, Particle const *const p2)
Definition: G4INCLKinematicsUtils.cc:94
G4INCL::ParticleTable::getIsospin
G4int getIsospin(const ParticleType t)
Get the isospin of a particle.
Definition: G4INCLParticleTable.cc:478
G4INCL::ParticleTable::getNucleonType
ParticleType getNucleonType(const G4int isosp)
Get the type of nucleon.
Definition: G4INCLParticleTable.cc:1253
G4INCL::ParticleTable::getPionType
ParticleType getPionType(const G4int isosp)
Get the type of pion.
Definition: G4INCLParticleTable.cc:1240
G4INCL::PhaseSpaceGenerator::generateBiased
void generateBiased(const G4double sqrtS, ParticleList &particles, const size_t index, const G4double slope)
Generate a biased event in the CM system.
Definition: G4INCLPhaseSpaceGenerator.cc:98
G4INCL::Random::getAdapter
Adapter const & getAdapter()
Definition: G4INCLRandom.cc:257
G4INCL::Random::shoot
G4double shoot()
Definition: G4INCLRandom.cc:93