source: HiSusy/trunk/Pythia8/pythia8170/src/SigmaOnia.cc @ 1

Last change on this file since 1 was 1, checked in by zerwas, 11 years ago

first import of structure, PYTHIA8 and DELPHES
File size: 20.1 KB
Line 
1// SigmaOnia.cc is a part of the PYTHIA event generator.
2// Copyright (C) 2012 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) ? settingsPtr->parm("Charmonium:OJpsi3S11")
29    : settingsPtr->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 (rndmPtr->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) ? settingsPtr->parm("Charmonium:Ochic03P01")
94    : settingsPtr->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 (rndmPtr->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) ? settingsPtr->parm("Charmonium:Ochic03P01")
189    : settingsPtr->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) ? settingsPtr->parm("Charmonium:Ochic03P01")
266    : settingsPtr->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) ? settingsPtr->parm("Charmonium:OJpsi3S18")
329      : settingsPtr->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) ? settingsPtr->parm("Charmonium:OJpsi1S08")
335      : settingsPtr->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) ? settingsPtr->parm("Charmonium:OJpsi3P08")
341      : settingsPtr->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 * rndmPtr->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 (rndmPtr->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) ? settingsPtr->parm("Charmonium:OJpsi3S18")
455      : settingsPtr->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) ? settingsPtr->parm("Charmonium:OJpsi1S08")
461      : settingsPtr->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) ? settingsPtr->parm("Charmonium:OJpsi3P08")
467      : settingsPtr->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 * rndmPtr->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) ? settingsPtr->parm("Charmonium:OJpsi3S18")
551      : settingsPtr->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) ? settingsPtr->parm("Charmonium:OJpsi1S08")
557      : settingsPtr->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) ? settingsPtr->parm("Charmonium:OJpsi3P08")
563      : settingsPtr->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 * rndmPtr->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
Note: See TracBrowser for help on using the repository browser.