| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
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 H
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 "G4INCLNNbarToNNbar2piChannel.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 NNbarToNNbar2piChannel::NNbarToNNbar2piChann
50 : particle1(p1), particle2(p2)
51 {}
52
53 NNbarToNNbar2piChannel::~NNbarToNNbar2piChan
54
55 void NNbarToNNbar2piChannel::fillFinalState(
56
57 //brief ppbar
58 // p pbar -> p pbar pi+ pi- (BFMM 167)
59 // p pbar -> p nbar pi- pi0 (same as B
60 // p pbar -> n pbar pi+ pi0 (same as B
61 // p pbar -> n nbar pi+ pi- (BFMM 198)
62 //
63 //brief npbar
64 // n pbar -> p pbar pi- pi0 (BFMM 490)
65 // n pbar -> p nbar pi- pi- (BFMM 492)
66 // n pbar -> n nbar pi- pi0 (same as B
67 // n pbar -> n pbar pi+ pi- (BFMM 494)
68 //
69 //brief nnbar
70 // n nbar -> n nbar pi+ pi- (same as B
71 // n nbar -> p nbar pi- pi0 (same as B
72 // n nbar -> n pbar pi+ pi0 (same as B
73 // n nbar -> p pbar pi+ pi- (same as B
74 //
75 //brief pnbar
76 // p nbar -> p pbar pi+ pi0 (same as B
77 // p nbar -> n pbar pi+ pi+ (same as B
78 // p nbar -> n nbar pi+ pi0 (same as B
79 // p nbar -> p nbar pi+ pi- (same as B
80
81 Particle *nucleon;
82 Particle *antinucleon;
83
84 if(particle1->isNucleon()){
85 nucleon = particle1;
86 antinucleon = particle2;
87 }
88 else{
89 nucleon = particle2;
90 antinucleon = particle1;
91 }
92
93 const G4double plab = 0.001*KinematicsUtil
94 const G4double sqrtS = KinematicsUtils::to
95 const G4double rdm = Random::shoot();
96
97 const std::vector<G4double> BFMM167 =
98 //const G4double Eth_PPbar_PPbar_pip_p
99 const std::vector<G4double> BFMM198 =
100 //const G4double Eth_PPbar_NNbar_pip_p
101 const std::vector<G4double> BFMM490 =
102 //const G4double Eth_PNbar_PPbar_pim_p
103 const std::vector<G4double> BFMM492 =
104 //const G4double Eth_PNbar_NPbar_pim_p
105 const std::vector<G4double> BFMM494 =
106 //const G4double Eth_NPbar_NPbar_pip_p
107
108 // pnbar total is same as for npbar
109 // ppbar total is same as for nnbar
110 const G4double totalppbar = KinematicsUtil
111 const G4double totalpnbar = KinematicsUtil
112 //totalnnbar == totalppbar;
113 //totalpnbar == totalnpbar;
114 ParticleType Pion1;
115 ParticleType Pion2;
116
117 //setting types of new particles
118 if(nucleon->getType()==Proton){
119 if(antinucleon->getType()==antiProton){
120 if(rdm*totalppbar < KinematicsUtils::c
121 Pion1 = PiMinus;
122 Pion2 = PiPlus;
123 if(rdm<0.5){
124 nucleon->setType(Proton);
125 antinucleon->setType(antiProton);
126 }
127 else{
128 nucleon->setType(antiProton);
129 antinucleon->setType(Proton);
130 }
131 }
132 else if(rdm*totalppbar < KinematicsUti
133 Pion1 = PiMinus;
134 Pion2 = PiZero;
135 if(rdm<0.5){
136 nucleon->setType(Proton);
137 antinucleon->setType(antiNeutron);
138 }
139 else{
140 nucleon->setType(antiNeutron);
141 antinucleon->setType(Proton);
142 }
143 }
144 else if(rdm*totalppbar < KinematicsUti
145 Pion1 = PiPlus;
146 Pion2 = PiZero;
147 if(rdm<0.5){
148 nucleon->setType(Neutron);
149 antinucleon->setType(antiProton);
150 }
151 else{
152 nucleon->setType(antiProton);
153 antinucleon->setType(Neutron);
154 }
155 }
156 else{ // n nbar pi+ pi- case case
157 Pion1 = PiMinus;
158 Pion2 = PiPlus;
159 if(rdm<0.5){
160 nucleon->setType(Neutron);
161 antinucleon->setType(antiNeutron);
162 }
163 else{
164 nucleon->setType(antiNeutron);
165 antinucleon->setType(Neutron);
166 }
167 }
168 }
169 else{ //antiNeutron (pnbar case)
170 if(rdm*totalpnbar < KinematicsUtils::c
171 Pion1 = PiZero;
172 Pion2 = PiPlus;
173 if(rdm<0.5){
174 nucleon->setType(Proton);
175 antinucleon->setType(antiProton);
176 }
177 else{
178 nucleon->setType(antiProton);
179 antinucleon->setType(Proton);
180 }
181 }
182 else if(rdm*totalppbar < KinematicsUti
183 Pion1 = PiPlus;
184 Pion2 = PiPlus;
185 if(rdm<0.5){
186 nucleon->setType(Neutron);
187 antinucleon->setType(antiProton);
188 }
189 else{
190 nucleon->setType(antiProton);
191 antinucleon->setType(Neutron);
192 }
193 }
194 else if(rdm*totalppbar < 2*KinematicsU
195 Pion1 = PiZero;
196 Pion2 = PiPlus;
197 if(rdm<0.5){
198 nucleon->setType(Neutron);
199 antinucleon->setType(antiNeutron);
200 }
201 else{
202 nucleon->setType(antiNeutron);
203 antinucleon->setType(Neutron);
204 }
205 }
206 else{ // p nbar pi+ pi- case
207 Pion1 = PiMinus;
208 Pion2 = PiPlus;
209 if(rdm<0.5){
210 nucleon->setType(Proton);
211 antinucleon->setType(antiNeutron);
212 }
213 else{
214 nucleon->setType(antiNeutron);
215 antinucleon->setType(Proton);
216 }
217 }
218 }
219 }
220 else{ // neutron
221 if(antinucleon->getType()==antiProton){
222 if(rdm*totalpnbar < KinematicsUtils::c
223 Pion1 = PiZero;
224 Pion2 = PiMinus;
225 if(rdm<0.5){
226 nucleon->setType(Proton);
227 antinucleon->setType(antiProton);
228 }
229 else{
230 nucleon->setType(antiProton);
231 antinucleon->setType(Proton);
232 }
233 }
234 else if(rdm*totalppbar < KinematicsUti
235 Pion1 = PiMinus;
236 Pion2 = PiMinus;
237 if(rdm<0.5){
238 nucleon->setType(Proton);
239 antinucleon->setType(antiNeutron);
240 }
241 else{
242 nucleon->setType(antiNeutron);
243 antinucleon->setType(Proton);
244 }
245 }
246 else if(rdm*totalppbar < 2*KinematicsU
247 Pion1 = PiZero;
248 Pion2 = PiMinus;
249 if(rdm<0.5){
250 nucleon->setType(Neutron);
251 antinucleon->setType(antiNeutron);
252 }
253 else{
254 nucleon->setType(antiNeutron);
255 antinucleon->setType(Neutron);
256 }
257 }
258 else{ // n pbar pi+ pi- case
259 Pion1 = PiMinus;
260 Pion2 = PiPlus;
261 if(rdm<0.5){
262 nucleon->setType(Neutron);
263 antinucleon->setType(antiProton);
264 }
265 else{
266 nucleon->setType(antiProton);
267 antinucleon->setType(Neutron);
268 }
269 }
270 }
271 else{ //antiNeutron (nnbar case)
272 if(rdm*totalppbar < KinematicsUtils::c
273 Pion1 = PiMinus;
274 Pion2 = PiPlus;
275 if(rdm<0.5){
276 nucleon->setType(Neutron);
277 antinucleon->setType(antiNeutron);
278 }
279 else{
280 nucleon->setType(antiNeutron);
281 antinucleon->setType(Neutron);
282 }
283 }
284 else if(rdm*totalppbar < KinematicsUti
285 Pion1 = PiMinus;
286 Pion2 = PiZero;
287 if(rdm<0.5){
288 nucleon->setType(Proton);
289 antinucleon->setType(antiNeutron);
290 }
291 else{
292 nucleon->setType(antiNeutron);
293 antinucleon->setType(Proton);
294 }
295 }
296 else if(rdm*totalppbar < KinematicsUti
297 Pion1 = PiPlus;
298 Pion2 = PiZero;
299 if(rdm<0.5){
300 nucleon->setType(Neutron);
301 antinucleon->setType(antiProton);
302 }
303 else{
304 nucleon->setType(antiProton);
305 antinucleon->setType(Neutron);
306 }
307 }
308 else{ // p pbar pi+ pi- case
309 Pion1 = PiMinus;
310 Pion2 = PiPlus;
311 if(rdm<0.5){
312 nucleon->setType(Proton);
313 antinucleon->setType(antiProton);
314 }
315 else{
316 nucleon->setType(antiProton);
317 antinucleon->setType(Proton);
318 }
319 }
320 }
321 }
322
323 ParticleList list;
324 list.push_back(nucleon);
325 list.push_back(antinucleon);
326 const ThreeVector &rcol = nucleon->getPosi
327 const ThreeVector zero;
328
329 Particle *pion2 = new Particle(Pion1,zero,
330 Particle *pion1 = new Particle(Pion2,zero,
331 if(rdm < 0.5){
332 pion1->setType(Pion1);
333 pion2->setType(Pion2);
334 }
335
336 list.push_back(pion1);
337 list.push_back(pion2);
338
339 PhaseSpaceGenerator::generate(sqrtS, list)
340
341 fs->addModifiedParticle(nucleon);
342 fs->addModifiedParticle(antinucleon);
343 fs->addCreatedParticle(pion1);
344 fs->addCreatedParticle(pion2);
345
346 }
347 }
348