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