]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PYTHIA8/pythia8130/src/SigmaOnia.cxx
using option '-treename HLTesdTree' for EsdCollector, adding default parameter for...
[u/mrichter/AliRoot.git] / PYTHIA8 / pythia8130 / src / SigmaOnia.cxx
CommitLineData
5ad4eb21 1// SigmaOnia.cc is a part of the PYTHIA event generator.
2// Copyright (C) 2008 Torbjorn Sjostrand.
3// PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
4// Please respect the MCnet Guidelines, see GUIDELINES for details.
5
6// Function definitions (not found in the header) for the
7// charmonia/bottomonia simulation classes.
8
9#include "SigmaOnia.h"
10
11namespace Pythia8 {
12
13//**************************************************************************
14
15// Sigma2gg2QQbar3S11g class.
16// Cross section g g -> QQbar[3S1(1)] g (Q = c or b).
17
18//*********
19
20// Initialize process.
21
22void Sigma2gg2QQbar3S11g::initProc() {
23
24 // Produced state. Process name. Onium matrix element.
25 idHad = (idNew == 4) ? 443 : 553;
26 nameSave = (idNew == 4) ? "g g -> ccbar[3S1(1)] g"
27 : "g g -> bbbar[3S1(1)] g";
28 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi3S11")
29 : Settings::parm("Bottomonium:OUpsilon3S11");
30
31}
32
33//*********
34
35// Evaluate d(sigmaHat)/d(tHat); no explicit flavour dependence.
36
37void Sigma2gg2QQbar3S11g::sigmaKin() {
38
39 // Calculate kinematics dependence.
40 double stH = sH + tH;
41 double tuH = tH + uH;
42 double usH = uH + sH;
43 double sig = (10. * M_PI / 81.) * m3 * ( pow2(sH * tuH)
44 + pow2(tH * usH) + pow2(uH * stH) ) / pow2( stH * tuH * usH );
45
46 // Answer.
47 sigma = (M_PI/sH2) * pow3(alpS) * oniumME * sig;
48
49}
50
51//*********
52
53// Select identity, colour and anticolour.
54
55void Sigma2gg2QQbar3S11g::setIdColAcol() {
56
57 // Flavours are trivial.
58 setId( id1, id2, idHad, 21);
59
60 // Two orientations of colour flow..
61 setColAcol( 1, 2, 2, 3, 0, 0, 1, 3);
62 if (Rndm::flat() > 0.5) swapColAcol();
63
64}
65
66//**************************************************************************
67
68// Sigma2gg2QQbar3PJ1g class.
69// Cross section g g -> QQbar[3PJ(1)] g (Q = c or b, J = 0, 1 or 2).
70
71//*********
72
73// Initialize process.
74
75void Sigma2gg2QQbar3PJ1g::initProc() {
76
77 // Produced state. Process name. Onium matrix element.
78 idHad = 0;
79 nameSave = "illegal process";
80 if (jSave == 0) {
81 idHad = (idNew == 4) ? 10441 : 10551;
82 nameSave = (idNew == 4) ? "g g -> ccbar[3P0(1)] g"
83 : "g g -> bbbar[3P0(1)] g";
84 } else if (jSave == 1) {
85 idHad = (idNew == 4) ? 20443 : 20553;
86 nameSave = (idNew == 4) ? "g g -> ccbar[3P1(1)] g"
87 : "g g -> bbbar[3P1(1)] g";
88 } else if (jSave == 2) {
89 idHad = (idNew == 4) ? 445 : 555;
90 nameSave = (idNew == 4) ? "g g -> ccbar[3P2(1)] g"
91 : "g g -> bbbar[3P2(1)] g";
92 }
93 oniumME = (idNew == 4) ? Settings::parm("Charmonium:Ochic03P01")
94 : Settings::parm("Bottomonium:Ochib03P01");
95
96}
97
98//*********
99
100// Evaluate d(sigmaHat)/d(tHat); no explicit flavour dependence.
101
102void Sigma2gg2QQbar3PJ1g::sigmaKin() {
103
104 // Useful derived kinematics quantities.
105 double pRat = (sH * uH + uH * tH + tH * sH)/ sH2;
106 double qRat = tH * uH / sH2;
107 double rRat = s3 / sH;
108 double pRat2 = pRat * pRat;
109 double pRat3 = pRat2 * pRat;
110 double pRat4 = pRat3 * pRat;
111 double qRat2 = qRat * qRat;
112 double qRat3 = qRat2 * qRat;
113 double qRat4 = qRat3 * qRat;
114 double rRat2 = rRat * rRat;
115 double rRat3 = rRat2 * rRat;
116 double rRat4 = rRat3 * rRat;
117
118 // Calculate kinematics dependence.
119 double sig = 0.;
120 if (jSave == 0) {
121 sig = (8. * M_PI / (9. * m3 * sH))
122 * ( 9. * rRat2 * pRat4 * (rRat4 - 2. * rRat2 * pRat + pRat2)
123 - 6. * rRat * pRat3 * qRat * (2. * rRat4 - 5. * rRat2 * pRat
124 + pRat2) - pRat2 * qRat2 * (rRat4 + 2. * rRat2 * pRat - pRat2)
125 + 2. * rRat * pRat * qRat3 * (rRat2 - pRat) + 6. * rRat2 * qRat4)
126 / (qRat * pow4(qRat - rRat * pRat));
127 } else if (jSave == 1) {
128 sig = (8. * M_PI / (3.* m3 * sH)) * pRat2
129 * (rRat * pRat2 * (rRat2 - 4. * pRat)
130 + 2. * qRat * (-rRat4 + 5. * rRat2 * pRat + pRat2)
131 - 15. * rRat * qRat2) / pow4(qRat - rRat * pRat);
132 } else if (jSave == 2) {
133 sig = (8. * M_PI / (9. * m3 * sH))
134 * (12. * rRat2 * pRat4 * (rRat4 - 2. * rRat2 * pRat + pRat2)
135 - 3. * rRat * pRat3 * qRat * (8. * rRat4 - rRat2 * pRat + 4. * pRat2)
136 + 2. * pRat2 * qRat2 * (-7. * rRat4 + 43. * rRat2 * pRat + pRat2)
137 + rRat * pRat * qRat3 * (16. * rRat2 - 61. * pRat)
138 + 12. * rRat2 * qRat4) / (qRat * pow4(qRat-rRat * pRat));
139 }
140
141 // Answer.
142 sigma = (M_PI/sH2) * pow3(alpS) * oniumME * sig;
143
144}
145
146//*********
147
148// Select identity, colour and anticolour.
149
150void Sigma2gg2QQbar3PJ1g::setIdColAcol() {
151
152 // Flavours are trivial.
153 setId( id1, id2, idHad, 21);
154
155 // Two orientations of colour flow.
156 setColAcol( 1, 2, 2, 3, 0, 0, 1, 3);
157 if (Rndm::flat() > 0.5) swapColAcol();
158
159}
160
161//**************************************************************************
162
163// Sigma2qg2QQbar3PJ1q class.
164// Cross section q g -> QQbar[3PJ(1)] q (Q = c or b, J = 0, 1 or 2).
165
166//*********
167
168// Initialize process.
169
170void Sigma2qg2QQbar3PJ1q::initProc() {
171
172 // Produced state. Process name. Onium matrix element.
173 idHad = 0;
174 nameSave = "illegal process";
175 if (jSave == 0) {
176 idHad = (idNew == 4) ? 10441 : 10551;
177 nameSave = (idNew == 4) ? "q g -> ccbar[3P0(1)] q"
178 : "q g -> bbbar[3P0(1)] q";
179 } else if (jSave == 1) {
180 idHad = (idNew == 4) ? 20443 : 20553;
181 nameSave = (idNew == 4) ? "q g -> ccbar[3P1(1)] q"
182 : "q g -> bbbar[3P1(1)] q";
183 } else if (jSave == 2) {
184 idHad = (idNew == 4) ? 445 : 555;
185 nameSave = (idNew == 4) ? "q g -> ccbar[3P2(1)] q"
186 : "q g -> bbbar[3P2(1)] q";
187 }
188 oniumME = (idNew == 4) ? Settings::parm("Charmonium:Ochic03P01")
189 : Settings::parm("Bottomonium:Ochib03P01");
190
191}
192
193//*********
194
195// Evaluate d(sigmaHat)/d(tHat); no explicit flavour dependence.
196
197void Sigma2qg2QQbar3PJ1q::sigmaKin() {
198
199 // Calculate kinematics dependence.
200 double usH = uH + sH;
201 double sig = 0.;
202 if (jSave == 0) {
203 sig = - (16. * M_PI / 81.) * pow2(tH - 3. * s3) * (sH2 + uH2)
204 / (m3 * tH * pow4(usH));
205 } else if (jSave == 1) {
206 sig = - (32. * M_PI / 27.) * (4. * s3 * sH * uH + tH * (sH2 + uH2))
207 / (m3 * pow4(usH));
208 } else if (jSave == 2) {
209 sig = - (32. *M_PI / 81.) * ( (6. * s3*s3 + tH2) * pow2(usH)
210 - 2. * sH * uH * (tH2 + 6. * s3 * usH)) / (m3 * tH * pow4(usH));
211 }
212
213 // Answer.
214 sigma = (M_PI/sH2) * pow3(alpS) * oniumME * sig;
215
216}
217
218//*********
219
220// Select identity, colour and anticolour.
221
222void Sigma2qg2QQbar3PJ1q::setIdColAcol() {
223
224 // Flavours are trivial.
225 int idq = (id2 == 21) ? id1 : id2;
226 setId( id1, id2, idHad, idq);
227
228 // tH defined between q_in and q_out: must swap tHat <-> uHat if q g in.
229 swapTU = (id2 == 21);
230
231 // Colour flow topologies. Swap when antiquarks.
232 if (id2 == 21) setColAcol( 1, 0, 2, 1, 0, 0, 2, 0);
233 else setColAcol( 2, 1, 1, 0, 0, 0, 2, 0);
234 if (idq < 0) swapColAcol();
235
236}
237
238//**************************************************************************
239
240// Sigma2qqbar2QQbar3PJ1g class.
241// Cross section q qbar -> QQbar[3PJ(1)] g (Q = c or b, J = 0, 1 or 2).
242
243//*********
244
245// Initialize process.
246
247void Sigma2qqbar2QQbar3PJ1g::initProc() {
248
249 // Produced state. Process name. Onium matrix element.
250 idHad = 0;
251 nameSave = "illegal process";
252 if (jSave == 0) {
253 idHad = (idNew == 4) ? 10441 : 10551;
254 nameSave = (idNew == 4) ? "q qbar -> ccbar[3P0(1)] g"
255 : "q qbar -> bbbar[3P0(1)] g";
256 } else if (jSave == 1) {
257 idHad = (idNew == 4) ? 20443 : 20553;
258 nameSave = (idNew == 4) ? "q qbar -> ccbar[3P1(1)] g"
259 : "q qbar -> bbbar[3P1(1)] g";
260 } else if (jSave == 2) {
261 idHad = (idNew == 4) ? 445 : 555;
262 nameSave = (idNew == 4) ? "q qbar -> ccbar[3P2(1)] g"
263 : "q qbar -> bbbar[3P2(1)] g";
264 }
265 oniumME = (idNew == 4) ? Settings::parm("Charmonium:Ochic03P01")
266 : Settings::parm("Bottomonium:Ochib03P01");
267
268}
269
270//*********
271
272// Evaluate d(sigmaHat)/d(tHat); no explicit flavour dependence.
273
274void Sigma2qqbar2QQbar3PJ1g::sigmaKin() {
275
276 // Calculate kinematics dependence.
277 double tuH = tH + uH;
278 double sig = 0.;
279 if (jSave == 0) {
280 sig =(128. * M_PI / 243.) * pow2(sH - 3. * s3) * (tH2 + uH2)
281 / (m3 * sH * pow4(tuH));
282 } else if (jSave == 1) {
283 sig = (256. * M_PI / 81.) * (4. * s3 * tH * uH + sH * (tH2 + uH2))
284 / (m3 * pow4(tuH));
285 } else if (jSave == 2) {
286 sig = (256. * M_PI / 243.) * ( (6. * s3*s3 + sH2) * pow2(tuH)
287 - 2. * tH * uH * (sH2 + 6. * s3 * tuH) )/ (m3 * sH * pow4(tuH));
288 }
289
290 // Answer.
291 sigma = (M_PI/sH2) * pow3(alpS) * oniumME * sig;
292
293}
294
295//*********
296
297// Select identity, colour and anticolour.
298
299void Sigma2qqbar2QQbar3PJ1g::setIdColAcol() {
300
301 // Flavours are trivial.
302 setId( id1, id2, idHad, 21);
303
304 // Colour flow topologies. Swap when antiquarks.
305 setColAcol( 1, 0, 0, 2, 0, 0, 1, 2);
306 if (id1 < 0) swapColAcol();
307
308}
309
310//**************************************************************************
311
312// Sigma2gg2QQbarX8g class.
313// Cross section g g -> QQbar[X(8)] g (Q = c or b, X = 3S1, 1S0 or 3PJ).
314
315//*********
316
317// Initialize process.
318
319void Sigma2gg2QQbarX8g::initProc() {
320
321 // Produced state. Process name. Onium matrix element.
322 idHad = 0;
323 nameSave = "illegal process";
324 if (stateSave == 0) {
325 idHad = (idNew == 4) ? 9900443 : 9900553;
326 nameSave = (idNew == 4) ? "g g -> ccbar[3S1(8)] g"
327 : "g g -> bbbar[3S1(8)] g";
328 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi3S18")
329 : Settings::parm("Bottomonium:OUpsilon3S18");
330 } else if (stateSave == 1) {
331 idHad = (idNew == 4) ? 9900441 : 9900551;
332 nameSave = (idNew == 4) ? "g g -> ccbar[1S0(8)] g"
333 : "g g -> bbbar[1S0(8)] g";
334 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi1S08")
335 : Settings::parm("Bottomonium:OUpsilon1S08");
336 } else if (stateSave == 2) {
337 idHad = (idNew == 4) ? 9910441 : 9910551;
338 nameSave = (idNew == 4) ? "g g -> ccbar[3PJ(8)] g"
339 : "g g -> bbbar[3PJ(8)] g";
340 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi3P08")
341 : Settings::parm("Bottomonium:OUpsilon3P08");
342 }
343
344}
345
346//*********
347
348// Evaluate d(sigmaHat)/d(tHat); no explicit flavour dependence.
349
350void Sigma2gg2QQbarX8g::sigmaKin() {
351
352 // Calculate kinematics dependence.
353 double stH = sH + tH;
354 double tuH = tH + uH;
355 double usH = uH + sH;
356 double sig = 0.;
357 if (stateSave == 0) {
358 sig = (M_PI / 72.) * m3 * ( 27. * (pow2(stH) + pow2(tuH)
359 + pow2(usH)) / (s3*s3) - 16. ) * ( pow2(sH * tuH)
360 + pow2(tH * usH) + pow2(uH * stH) ) / pow2( stH * tuH * usH );
361 } else if (stateSave == 1) {
362 sig = (5. * M_PI / 16.) * m3 * ( pow2(uH / (tuH * usH))
363 + pow2(sH / (stH * usH)) + pow2(tH / (stH * tuH)) ) * ( 12.
364 + (pow4(stH) + pow4(tuH) + pow4(usH)) / (s3 * sH * tH * uH) );
365 } else if (stateSave == 2) {
366 double sH3 = sH2 * sH;
367 double sH4 = sH3 * sH;
368 double sH5 = sH4 * sH;
369 double sH6 = sH5 * sH;
370 double sH7 = sH6 * sH;
371 double sH8 = sH7 * sH;
372 double tH3 = tH2 * tH;
373 double tH4 = tH3 * tH;
374 double tH5 = tH4 * tH;
375 double tH6 = tH5 * tH;
376 double tH7 = tH6 * tH;
377 double tH8 = tH7 * tH;
378 double ssttH = sH * sH + sH * tH + tH * tH;
379 sig = 5. * M_PI * (3. * sH * tH * stH * pow4(ssttH)
380 - s3 * pow2(ssttH) * (7. * sH6 + 36. * sH5 * tH + 45. * sH4 * tH2
381 + 28. * sH3 * tH3 + 45. * sH2 * tH4 + 36. * sH * tH5 + 7. * tH6)
382 + pow2(s3) * stH * (35. *sH8 + 169. * sH7 * tH + 299. * sH6 * tH2
383 + 401. * sH5 * tH3 + 418. * sH4 * tH4 + 401. * sH3 * tH5
384 + 299. * sH2 * tH6 + 169. * sH * tH7 + 35. * tH8)
385 - pow3(s3) * (84. *sH8+432. *sH7*tH+905. *sH6*tH2
386 + 1287. * sH5 * tH3 + 1436. * sH4 * tH4 +1287. * sH3 * tH5
387 + 905. * sH2 * tH6 + 432. * sH * tH7 + 84. * tH8)
388 + pow4(s3) * stH * (126. * sH6 + 451. * sH5 * tH +677. * sH4 * tH2
389 + 836. * sH3 * tH3 + 677. * sH2 * tH4 + 451. * sH * tH5
390 + 126. * tH6)
391 - pow5(s3) * 3. * (42. * sH6 + 171. * sH5 * tH + 304. * sH4 * tH2
392 + 362. * sH3 * tH3 + 304. * sH2 * tH4 + 171. * sH * tH5
393 + 42. * tH6)
394 + pow3(s3 * s3) * 2. * stH * (42. * sH4 + 106. * sH3 * tH
395 + 119. * sH2 * tH2 + 106. * sH * tH3 + 42. * tH4)
396 - pow4(s3) * pow3(s3) * (35. * sH4 + 99. * sH3 * tH
397 + 120. * sH2 * tH2 + 99. * sH * tH3 + 35. * tH4)
398 + pow4(s3 * s3) * 7. * stH * ssttH)
399 / (sH * tH * uH * s3 * m3 * pow3(stH * tuH * usH));
400 }
401
402 // Answer.
403 sigma = (M_PI/sH2) * pow3(alpS) * oniumME * sig;
404
405}
406
407//*********
408
409// Select identity, colour and anticolour.
410
411void Sigma2gg2QQbarX8g::setIdColAcol() {
412
413 // Flavours are trivial.
414 setId( id1, id2, idHad, 21);
415
416 // Split total contribution into different colour flows just like in
417 // g g -> g g (with kinematics recalculated for massless partons).
418 double sHr = - (tH + uH);
419 double sH2r = sHr * sHr;
420 double sigTS = tH2/sH2r + 2.*tH/sHr + 3. + 2.*sHr/tH + sH2r/tH2;
421 double sigUS = uH2/sH2r + 2.*uH/sHr + 3. + 2.*sHr/uH + sH2r/uH2;
422 double sigTU = tH2/uH2 + 2.*tH/uH + 3. + 2.*uH/tH + uH2/tH2;
423 double sigSum = sigTS + sigUS + sigTU;
424
425 // Three colour flow topologies, each with two orientations.
426 double sigRand = sigSum * Rndm::flat();
427 if (sigRand < sigTS) setColAcol( 1, 2, 2, 3, 1, 4, 4, 3);
428 else if (sigRand < sigTS + sigUS)
429 setColAcol( 1, 2, 3, 1, 3, 4, 4, 2);
430 else setColAcol( 1, 2, 3, 4, 1, 4, 3, 2);
431 if (Rndm::flat() > 0.5) swapColAcol();
432
433
434}
435
436//**************************************************************************
437
438// Sigma2qg2QQbarX8q class.
439// Cross section q g -> QQbar[X(8)] q (Q = c or b, X = 3S1, 1S0 or 3PJ).
440
441//*********
442
443// Initialize process.
444
445void Sigma2qg2QQbarX8q::initProc() {
446
447 // Produced state. Process name. Onium matrix element.
448 idHad = 0;
449 nameSave = "illegal process";
450 if (stateSave == 0) {
451 idHad = (idNew == 4) ? 9900443 : 9900553;
452 nameSave = (idNew == 4) ? "q g -> ccbar[3S1(8)] q"
453 : "q g -> bbbar[3S1(8)] q";
454 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi3S18")
455 : Settings::parm("Bottomonium:OUpsilon3S18");
456 } else if (stateSave == 1) {
457 idHad = (idNew == 4) ? 9900441 : 9900551;
458 nameSave = (idNew == 4) ? "q g -> ccbar[1S0(8)] q"
459 : "q g -> bbbar[1S0(8)] q";
460 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi1S08")
461 : Settings::parm("Bottomonium:OUpsilon1S08");
462 } else if (stateSave == 2) {
463 idHad = (idNew == 4) ? 9910441 : 9910551;
464 nameSave = (idNew == 4) ? "q g -> ccbar[3PJ(8)] q"
465 : "q g -> bbbar[3PJ(8)] q";
466 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi3P08")
467 : Settings::parm("Bottomonium:OUpsilon3P08");
468 }
469
470}
471
472//*********
473
474// Evaluate d(sigmaHat)/d(tHat); no explicit flavour dependence.
475
476void Sigma2qg2QQbarX8q::sigmaKin() {
477
478 // Calculate kinematics dependence.
479 double stH = sH + tH;
480 double tuH = tH + uH;
481 double usH = uH + sH;
482 double stH2 = stH * stH;
483 double tuH2 = tuH * tuH;
484 double usH2 = usH * usH;
485 double sig = 0.;
486 if (stateSave == 0) {
487 sig = - (M_PI / 27.)* (4. * (sH2 + uH2) - sH * uH) * (stH2 +tuH2)
488 / (s3 * m3 * sH * uH * usH2);
489 } else if (stateSave == 1) {
490 sig = - (5. * M_PI / 18.) * (sH2 + uH2) / (m3 * tH * usH2);
491 } else if (stateSave == 2) {
492 sig = - (10. * M_PI / 9.) * ( (7. * usH + 8. * tH) * (sH2 + uH2)
493 + 4. * tH * (2. * pow2(s3) - stH2 - tuH2) )
494 / (s3 * m3 * tH * usH2 * usH);
495 }
496
497 // Answer.
498 sigma = (M_PI/sH2) * pow3(alpS) * oniumME * sig;
499
500}
501
502//*********
503
504// Select identity, colour and anticolour.
505
506void Sigma2qg2QQbarX8q::setIdColAcol() {
507
508 // Flavours are trivial.
509 int idq = (id2 == 21) ? id1 : id2;
510 setId( id1, id2, idHad, idq);
511
512 // tH defined between q_in and q_out: must swap tHat <-> uHat if q g in.
513 swapTU = (id2 == 21);
514
515 // Split total contribution into different colour flows just like in
516 // q g -> q g (with kinematics recalculated for massless partons).
517 double sHr = - (tH + uH);
518 double sH2r = sHr * sHr;
519 double sigTS = uH2/tH2 - (4./9.) * uH/sHr;
520 double sigTU = sH2r/tH2 - (4./9.) * sHr/uH;
521 double sigSum = sigTS + sigTU;
522
523 // Two colour flow topologies. Swap if first is gluon, or when antiquark.
524 double sigRand = sigSum * Rndm::flat();
525 if (sigRand < sigTS) setColAcol( 1, 0, 2, 1, 2, 3, 3, 0);
526 else setColAcol( 1, 0, 2, 3, 1, 3, 2, 0);
527 if (id1 == 21) swapCol12();
528 if (idq < 0) swapColAcol();
529
530}
531
532//**************************************************************************
533
534// Sigma2qqbar2QQbarX8g class.
535// Cross section q qbar -> QQbar[X(8)] g (Q = c or b, X = 3S1, 1S0 or 3PJ).
536
537//*********
538
539// Initialize process.
540
541void Sigma2qqbar2QQbarX8g::initProc() {
542
543 // Produced state. Process name. Onium matrix element.
544 idHad = 0;
545 nameSave = "illegal process";
546 if (stateSave == 0) {
547 idHad = (idNew == 4) ? 9900443 : 9900553;
548 nameSave = (idNew == 4) ? "q qbar -> ccbar[3S1(8)] g"
549 : "q qbar -> bbbar[3S1(8)] g";
550 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi3S18")
551 : Settings::parm("Bottomonium:OUpsilon3S18");
552 } else if (stateSave == 1) {
553 idHad = (idNew == 4) ? 9900441 : 9900551;
554 nameSave = (idNew == 4) ? "q qbar -> ccbar[1S0(8)] g"
555 : "q qbar -> bbbar[1S0(8)] g";
556 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi1S08")
557 : Settings::parm("Bottomonium:OUpsilon1S08");
558 } else if (stateSave == 2) {
559 idHad = (idNew == 4) ? 9910441 : 9910551;
560 nameSave = (idNew == 4) ? "q qbar -> ccbar[3PJ(8)] g"
561 : "q qbar -> bbbar[3PJ(8)] g";
562 oniumME = (idNew == 4) ? Settings::parm("Charmonium:OJpsi3P08")
563 : Settings::parm("Bottomonium:OUpsilon3P08");
564 }
565
566}
567
568//*********
569
570// Evaluate d(sigmaHat)/d(tHat); no explicit flavour dependence.
571
572void Sigma2qqbar2QQbarX8g::sigmaKin() {
573
574 // Calculate kinematics dependence.
575 double stH = sH + tH;
576 double tuH = tH + uH;
577 double usH = uH + sH;
578 double stH2 = stH * stH;
579 double tuH2 = tuH * tuH;
580 double usH2 = usH * usH;
581 double sig = 0.;
582 if (stateSave == 0) {
583 sig = (8. * M_PI / 81.) * (4. * (tH2 + uH2) - tH * uH)
584 * (stH2 + usH2) / (s3 * m3 * tH * uH * tuH2);
585 } else if (stateSave == 1) {
586 sig = (20. * M_PI / 27.) * (tH2 + uH2) / (m3 * sH * tuH2);
587 } else if (stateSave == 2) {
588 sig = (80. * M_PI / 27.) * ( (7. * tuH + 8. * sH) * (tH2 + uH2)
589 + 4. * sH * (2. * pow2(s3) - stH2 -usH2) )
590 / (s3 * m3 * sH * tuH2 * tuH);
591 }
592
593 // Answer.
594 sigma = (M_PI/sH2) * pow3(alpS) * oniumME * sig;
595
596}
597
598//*********
599
600// Select identity, colour and anticolour.
601
602void Sigma2qqbar2QQbarX8g::setIdColAcol() {
603
604 // Flavours are trivial.
605 setId( id1, id2, idHad, 21);
606
607 // Split total contribution into different colour flows just like in
608 // q qbar -> g g (with kinematics recalculated for massless partons).
609 double sHr = - (tH + uH);
610 double sH2r = sHr * sHr;
611 double sigTS = (4. / 9.) * uH / tH - uH2 / sH2r;
612 double sigUS = (4. / 9.) * tH / uH - tH2 / sH2r;
613 double sigSum = sigTS + sigUS;
614
615 // Two colour flow topologies. Swap if first is antiquark.
616 double sigRand = sigSum * Rndm::flat();
617 if (sigRand < sigTS) setColAcol( 1, 0, 0, 2, 1, 3, 3, 2);
618 else setColAcol( 1, 0, 0, 2, 3, 2, 1, 3);
619 if (id1 < 0) swapColAcol();
620
621}
622
623//**************************************************************************
624
625} // end namespace Pythia8
626