| 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 "G4INCLNKbToNKbChannel.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 NKbToNKbChannel::NKbToNKbChannel(Particle *p
50 : particle1(p1), particle2(p2)
51 {}
52
53 NKbToNKbChannel::~NKbToNKbChannel(){}
54
55 void NKbToNKbChannel::fillFinalState(FinalSt
56
57 Particle *nucleon;
58 Particle *kaon;
59
60 if(particle1->isNucleon()){
61 nucleon = particle1;
62 kaon = particle2;
63 }
64 else{
65 nucleon = particle2;
66 kaon = particle1;
67 }
68
69 // assert((ParticleTable::getIsospin(nucleon->
70
71 ThreeVector mom_kaon = KaonMomentum(kaon,n
72
73 if(kaon->getType() == KZeroBar){
74 nucleon->setType(Proton);
75 kaon->setType(KMinus);
76 }
77 else{
78 nucleon->setType(Neutron);
79 kaon->setType(KZeroBar);
80 }
81
82 G4double norm = KinematicsUtils::momentumI
83
84 kaon->setMomentum(mom_kaon*norm);
85 nucleon->setMomentum(-mom_kaon*norm);
86
87 kaon->adjustEnergyFromMomentum();
88 nucleon->adjustEnergyFromMomentum();
89
90
91 fs->addModifiedParticle(nucleon);
92 fs->addModifiedParticle(kaon);
93
94 }
95
96 ThreeVector NKbToNKbChannel::KaonMomentum(Pa
97
98 const G4double pLab = KinematicsUtils::mom
99
100 if(pLab < 235.) return Random::normVector(
101
102 G4double cos_theta = 1.;
103 G4double sin_theta = 0.;
104 const G4double cos_phi = std::cos(Random::
105 const G4double sin_phi = std::sqrt(1-cos_p
106
107 const G4double x = kaon->getMomentum().get
108 const G4double y = kaon->getMomentum().get
109 const G4double z = kaon->getMomentum().get
110
111 const G4double r = std::sqrt(x*x+y*y+z*z);
112 const G4double rho = std::sqrt(x*x+y*y);
113
114 if(pLab >= 1355.){
115 const G4double b = 12. * pLab/2375.; //
116 cos_theta = std::log(Random::shoot()*(st
117 sin_theta = std::sqrt(1-cos_theta*cos_th
118
119 }
120 else{
121 const G4double Legendre_coef[225][9] = {
122 {235,-0.30755,0,-0.04859,-0.01348,-9e-
123 {240,-0.3041,0,-0.03549,-0.01096,-8e-0
124 {245,-0.30451,0,-0.02232,-0.00843,-8e-
125 {250,-0.31028,0,-0.00899,-0.00585,-7e-
126 {255,-0.32259,0,0.00462,-0.00319,-6e-0
127 {260,-0.34185,0,0.01859,-0.00041,-5e-0
128 {265,-0.36746,0,0.03305,0.00251,-4e-05
129 {270,-0.39766,0,0.04809,0.00563,-3e-05
130 {275,-0.42975,0,0.0638,0.00898,-2e-05,
131 {280,-0.46069,1e-05,0.08023,0.0126,0,1
132 {285,-0.48736,5e-05,0.0973,0.0165,1e-0
133 {290,-0.50729,2e-04,0.11487,0.02069,4e
134 {295,-0.51855,0.00063,0.13282,0.02519,
135 {300,-0.52024,0.00134,0.15096,0.02998,
136 {305,-0.51207,0.00124,0.16892,0.03503,
137 {310,-0.49468,-0.002,0.18618,0.04031,0
138 {315,-0.46892,-0.00933,0.20212,0.04568
139 {320,-0.43609,-0.01586,0.21604,0.05101
140 {325,-0.3973,-0.01303,0.22726,0.05617,
141 {330,-0.35379,0.00167,0.23511,0.06107,
142 {335,-0.3067,0.0219,0.23911,0.06576,0.
143 {340,-0.25739,0.04678,0.23878,0.07033,
144 {345,-0.20755,0.09094,0.23367,0.07488,
145 {350,-0.15908,0.17482,0.22334,0.0795,0
146 {355,-0.11446,0.30682,0.20753,0.08421,
147 {360,-0.07625,0.48075,0.1861,0.08906,0
148 {365,-0.04747,0.68345,0.1593,0.09421,0
149 {370,-0.03064,0.90388,0.1275,0.09983,0
150 {375,-0.02812,1.13397,0.09104,0.10608,
151 {380,-0.04087,1.3636,0.05039,0.11311,0
152 {385,-0.06879,1.57657,0.00644,0.12093,
153 {390,-0.10972,1.75346,-0.03978,0.1295,
154 {395,-0.16022,1.87981,-0.08707,0.13859
155 {400,-0.21553,1.9509,-0.1342,0.14795,0
156 {405,-0.27043,1.97287,-0.17995,0.15732
157 {410,-0.32009,1.95458,-0.22323,0.16643
158 {415,-0.36095,1.89945,-0.26354,0.1749,
159 {420,-0.3912,1.80489,-0.30052,0.18232,
160 {425,-0.4106,1.68228,-0.3338,0.18828,-
161 {430,-0.42021,1.56358,-0.36315,0.19246
162 {435,-0.42128,1.47378,-0.38888,0.19478
163 {440,-0.415,1.42688,-0.41138,0.19524,-
164 {445,-0.40201,1.39998,-0.43103,0.19376
165 {450,-0.38247,1.36801,-0.44826,0.19036
166 {455,-0.35589,1.33185,-0.46368,0.18542
167 {460,-0.3216,1.29215,-0.47787,0.17941,
168 {465,-0.27917,1.2499,-0.49102,0.1727,-
169 {470,-0.22973,1.20597,-0.50321,0.16567
170 {475,-0.17713,1.16167,-0.5145,0.15865,
171 {480,-0.12841,1.11813,-0.52493,0.15196
172 {485,-0.09149,1.07666,-0.53436,0.14575
173 {490,-0.07028,1.03839,-0.5426,0.14013,
174 {495,-0.06296,1.00433,-0.54947,0.13518
175 {500,-0.06366,0.97503,-0.55477,0.13101
176 {505,-0.06664,0.95085,-0.55833,0.12773
177 {510,-0.06771,0.93169,-0.55996,0.12544
178 {515,-0.06476,0.91724,-0.55947,0.12423
179 {520,-0.05705,0.90697,-0.55703,0.12437
180 {525,-0.04488,0.90021,-0.55318,0.12594
181 {530,-0.02927,0.8963,-0.54843,0.12893,
182 {535,-0.01183,0.8945,-0.5433,0.13334,-
183 {540,0.00596,0.89424,-0.53831,0.13906,
184 {545,0.02306,0.89492,-0.5339,0.1458,-0
185 {550,0.039,0.89617,-0.53048,0.15329,-0
186 {555,0.05391,0.89758,-0.52847,0.16125,
187 {560,0.06806,0.89897,-0.52824,0.16947,
188 {565,0.08166,0.90012,-0.52974,0.17794,
189 {570,0.09478,0.90105,-0.53283,0.18651,
190 {575,0.10746,0.9017,-0.53736,0.19494,-
191 {580,0.11985,0.90216,-0.54318,0.20297,
192 {585,0.13218,0.90248,-0.55002,0.21039,
193 {590,0.14481,0.90282,-0.55752,0.21698,
194 {595,0.15769,0.90328,-0.56531,0.22255,
195 {600,0.17046,0.90397,-0.57302,0.2269,-
196 {605,0.18205,0.90496,-0.58026,0.23009,
197 {610,0.19103,0.90635,-0.58673,0.23229,
198 {615,0.19589,0.90814,-0.59209,0.2337,-
199 {620,0.19545,0.91035,-0.59603,0.2345,0
200 {625,0.18873,0.91297,-0.59864,0.23485,
201 {630,0.17514,0.91598,-0.60046,0.2349,0
202 {635,0.15461,0.91935,-0.60188,0.23486,
203 {640,0.12791,0.92305,-0.60283,0.23506,
204 {645,0.09694,0.92707,-0.60322,0.23576,
205 {650,0.06435,0.9314,-0.60304,0.23717,0
206 {655,0.03278,0.93603,-0.60224,0.23951,
207 {660,0.00403,0.94096,-0.60081,0.24298,
208 {665,-0.02138,0.9462,-0.59882,0.24774,
209 {670,-0.04399,0.95174,-0.59638,0.25384
210 {675,-0.06489,0.95757,-0.59359,0.26128
211 {680,-0.08497,0.96367,-0.59053,0.27011
212 {685,-0.10471,0.97,-0.58712,0.2802,0.0
213 {690,-0.12371,0.97649,-0.58315,0.29136
214 {695,-0.14097,0.98305,-0.57828,0.3031,
215 {700,-0.15469,0.98957,-0.57209,0.31472
216 {705,-0.16271,0.99591,-0.56443,0.3256,
217 {710,-0.16255,1.00194,-0.55574,0.33535
218 {715,-0.15199,1.00752,-0.54657,0.34363
219 {720,-0.12951,1.01251,-0.53743,0.3499,
220 {725,-0.09529,1.0168,-0.52879,0.35376,
221 {730,-0.05159,1.0203,-0.5212,0.35562,0
222 {735,-0.00245,1.02293,-0.51519,0.35595
223 {740,0.04772,1.02468,-0.51129,0.35522,
224 {745,0.09584,1.02553,-0.51005,0.35387,
225 {750,0.14086,1.02552,-0.51202,0.35229,
226 {755,0.18338,1.02466,-0.51778,0.35084,
227 {760,0.22444,1.023,-0.52787,0.34965,0.
228 {765,0.2645,1.0206,-0.54268,0.34874,-0
229 {770,0.30276,1.01748,-0.56211,0.34808,
230 {775,0.33729,1.01369,-0.58582,0.34762,
231 {780,0.36532,1.00926,-0.61343,0.34741,
232 {785,0.38445,1.00421,-0.64453,0.34758,
233 {790,0.39323,0.99857,-0.67869,0.34828,
234 {795,0.39129,0.99243,-0.71549,0.34967,
235 {800,0.37889,0.98587,-0.75456,0.35269,
236 {805,0.35624,0.97903,-0.79561,0.35848,
237 {810,0.32292,0.97211,-0.83834,0.36717,
238 {815,0.27791,0.96533,-0.88246,0.37885,
239 {820,0.22029,0.95894,-0.92768,0.39359,
240 {825,0.15048,0.9532,-0.9737,0.41145,-0
241 {830,0.07108,0.94826,-1.02023,0.43252,
242 {835,-0.0136,0.94426,-1.06697,0.45687,
243 {840,-0.09855,0.94118,-1.11361,0.48455
244 {845,-0.17951,0.9389,-1.15976,0.51525,
245 {850,-0.25366,0.93727,-1.20484,0.54811
246 {855,-0.31965,0.93606,-1.24827,0.58227
247 {860,-0.37707,0.93506,-1.28948,0.61685
248 {865,-0.42612,0.93404,-1.3279,0.651,-0
249 {870,-0.46711,0.93281,-1.36298,0.68418
250 {875,-0.50042,0.93129,-1.39424,0.7161,
251 {880,-0.52649,0.92936,-1.42118,0.74648
252 {885,-0.54572,0.92699,-1.44339,0.77492
253 {890,-0.55824,0.92417,-1.46054,0.801,-
254 {895,-0.56399,0.92098,-1.47259,0.82492
255 {900,-0.56271,0.91754,-1.47975,0.84731
256 {905,-0.55466,0.91396,-1.48218,0.86857
257 {910,-0.54086,0.91044,-1.48007,0.88907
258 {915,-0.52332,0.90717,-1.47358,0.90921
259 {920,-0.5045,0.90438,-1.46288,0.92935,
260 {925,-0.48654,0.9023,-1.4482,0.94975,-
261 {930,-0.47053,0.90114,-1.42984,0.97058
262 {935,-0.45659,0.90112,-1.40806,0.99201
263 {940,-0.44423,0.90246,-1.38314,1.01421
264 {945,-0.43287,0.90525,-1.35536,1.0371,
265 {950,-0.42193,0.90965,-1.32498,1.06041
266 {955,-0.41093,0.91577,-1.29226,1.08389
267 {960,-0.39924,0.92369,-1.25746,1.10725
268 {965,-0.38604,0.93346,-1.22095,1.13008
269 {970,-0.37027,0.94514,-1.18303,1.15197
270 {975,-0.35077,0.95876,-1.14385,1.17271
271 {980,-0.32717,0.97435,-1.10358,1.19206
272 {985,-0.3004,0.99188,-1.0624,1.20979,-
273 {990,-0.27312,1.01135,-1.02046,1.2257,
274 {995,-0.24834,1.03272,-0.97795,1.23954
275 {1000,-0.22769,1.0559,-0.93501,1.2511,
276 {1005,-0.21032,1.0808,-0.89184,1.26017
277 {1010,-0.19369,1.10728,-0.84866,1.2666
278 {1015,-0.17436,1.13513,-0.80576,1.2706
279 {1020,-0.1494,1.16411,-0.7634,1.27208,
280 {1025,-0.11749,1.19391,-0.72184,1.2709
281 {1030,-0.081,1.22417,-0.68137,1.26732,
282 {1035,-0.04533,1.25445,-0.64225,1.2611
283 {1040,-0.01646,1.28424,-0.60476,1.2524
284 {1045,0.00361,1.31304,-0.56915,1.24121
285 {1050,0.01552,1.34024,-0.53567,1.22751
286 {1055,0.02103,1.36529,-0.50435,1.21159
287 {1060,0.02215,1.38771,-0.47519,1.19376
288 {1065,0.02047,1.40712,-0.4482,1.1743,-
289 {1070,0.01735,1.4233,-0.42337,1.15353,
290 {1075,0.0137,1.43585,-0.40071,1.13173,
291 {1080,0.00983,1.44489,-0.38022,1.10922
292 {1085,0.00516,1.45066,-0.3619,1.08628,
293 {1090,-0.00091,1.45292,-0.34574,1.0632
294 {1095,-0.00921,1.4522,-0.33171,1.04013
295 {1100,-0.01996,1.44884,-0.31977,1.0172
296 {1105,-0.03239,1.44326,-0.30987,0.9945
297 {1110,-0.0447,1.43582,-0.30197,0.97238
298 {1115,-0.05508,1.42673,-0.29604,0.9507
299 {1120,-0.06264,1.41667,-0.29202,0.9297
300 {1125,-0.06741,1.40585,-0.28988,0.9096
301 {1130,-0.06996,1.39476,-0.28958,0.8905
302 {1135,-0.07105,1.38376,-0.29107,0.8724
303 {1140,-0.0711,1.3732,-0.29431,0.85567,
304 {1145,-0.07085,1.36335,-0.29927,0.8402
305 {1150,-0.07116,1.35452,-0.30589,0.8263
306 {1155,-0.07304,1.34687,-0.31414,0.8140
307 {1160,-0.07771,1.34061,-0.32397,0.8035
308 {1165,-0.08664,1.3359,-0.33535,0.79499
309 {1170,-0.10173,1.33284,-0.34823,0.7884
310 {1175,-0.12509,1.33145,-0.36258,0.7835
311 {1180,-0.15682,1.33177,-0.37838,0.7801
312 {1185,-0.19464,1.33374,-0.39559,0.7778
313 {1190,-0.23352,1.33727,-0.41419,0.7763
314 {1195,-0.26876,1.34221,-0.43415,0.7752
315 {1200,-0.299,1.34841,-0.45544,0.77437,
316 {1205,-0.32631,1.35562,-0.47804,0.7733
317 {1210,-0.35505,1.36358,-0.50188,0.7718
318 {1215,-0.38961,1.372,-0.52686,0.7699,-
319 {1220,-0.43176,1.38055,-0.55282,0.7674
320 {1225,-0.47936,1.38889,-0.57962,0.7644
321 {1230,-0.52635,1.39666,-0.60711,0.7609
322 {1235,-0.56649,1.40351,-0.63516,0.7568
323 {1240,-0.59709,1.40914,-0.66361,0.7522
324 {1245,-0.61879,1.41323,-0.69232,0.7469
325 {1250,-0.63452,1.41555,-0.72116,0.7411
326 {1255,-0.6472,1.4159,-0.74996,0.73469,
327 {1260,-0.65841,1.4142,-0.77859,0.7276,
328 {1265,-0.66819,1.41043,-0.80691,0.7198
329 {1270,-0.67506,1.40465,-0.83476,0.7114
330 {1275,-0.67695,1.39703,-0.86202,0.7023
331 {1280,-0.67215,1.38781,-0.88852,0.6925
332 {1285,-0.65883,1.37728,-0.91412,0.6819
333 {1290,-0.63639,1.36579,-0.93873,0.6707
334 {1295,-0.60463,1.35371,-0.96237,0.6589
335 {1300,-0.56527,1.34141,-0.98514,0.6466
336 {1305,-0.52229,1.32923,-1.0071,0.63407
337 {1310,-0.4806,1.31751,-1.02834,0.62132
338 {1315,-0.44384,1.3065,-1.04893,0.60855
339 {1320,-0.41298,1.29641,-1.06896,0.5959
340 {1325,-0.38708,1.28739,-1.0885,0.58352
341 {1330,-0.36429,1.27952,-1.10764,0.5715
342 {1335,-0.34347,1.27284,-1.12645,0.5601
343 {1340,-0.32464,1.26732,-1.14501,0.5493
344 {1345,-0.30854,1.2629,-1.1634,0.5395,-
345 {1350,-0.29584,1.25949,-1.18169,0.5305
346 {1355,-0.28655,1.25697,-1.19998,0.5228
347
348 const G4int coef_ener = G4int((pLab-Lege
349 const G4double sup_ener = pLab/5. - coef
350
351 // assert(pLab >= Legendre_coef[coef_ener][0]
352
353 // Legendre coefficient normalized
354 const G4double A0 = 1.;
355 const G4double A1 = (1-sup_ener)*Legendr
356 const G4double A2 = (1-sup_ener)*Legendr
357 const G4double A3 = (1-sup_ener)*Legendr
358 const G4double A4 = (1-sup_ener)*Legendr
359 const G4double A5 = (1-sup_ener)*Legendr
360 const G4double A6 = (1-sup_ener)*Legendr
361 const G4double A7 = (1-sup_ener)*Legendr
362 const G4double A8 = (1-sup_ener)*Legendr
363
364 // Theoritical max if all Ai > 0 (often
365 const G4double A = std::fabs(A0) + std::
366
367 G4bool success = false;
368 G4int maxloop = 0;
369
370 while(!success && maxloop < 1000){
371
372 cos_theta = Random::shoot()*2-1.; // n
373
374 // Legendre Polynomial
375 G4double P0 = A0;
376 G4double P1 = A1*cos_theta;
377 G4double P2 = A2/2.*(3*std::pow(cos_th
378 G4double P3 = A3/2.*(5*std::pow(cos_th
379 G4double P4 = A4/8.*(35*std::pow(cos_t
380 G4double P5 = A5/8.*(63*std::pow(cos_t
381 G4double P6 = A6/16.*(231*std::pow(cos
382 G4double P7 = A7/16.*(429*std::pow(cos
383 G4double P8 = A8/128.*(6435*std::pow(c
384
385 G4double P = (P0 + P1 + P2 + P3 + P4 +
386
387 if(Random::shoot()*A < P) success = tr
388 maxloop +=1 ;
389 if(maxloop==1000) cos_theta = std::log
390 }
391 sin_theta = std::sqrt(1-cos_theta*cos_th
392 }
393
394 if(rho == 0) return ThreeVector(sin_theta*
395 // Rotation in the direction of the incide
396 const G4double px = x/r*cos_theta - y/rho*
397 const G4double py = y/r*cos_theta + x/rho*
398 const G4double pz = z/r*cos_theta - rho/r*
399
400 return ThreeVector(px,py,pz);
401 }
402 }
403