using option '-treename HLTesdTree' for EsdCollector, adding default parameter for...
[u/mrichter/AliRoot.git] / PYTHIA8 / pythia8130 / src / SigmaExtraDim.cxx
CommitLineData
5ad4eb21 1// SigmaExtraDim.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// extra-dimensional simulation classes.
8
9#include "SigmaExtraDim.h"
10
11namespace Pythia8 {
12
13//**************************************************************************
14
15// Sigma1gg2GravitonStar class.
16// Cross section for g g -> G* (excited graviton state).
17
18//*********
19
20// Initialize process.
21
22void Sigma1gg2GravitonStar::initProc() {
23
24 // Store G* mass and width for propagator.
25 idGstar = 5000039;
26 mRes = ParticleDataTable::m0(idGstar);
27 GammaRes = ParticleDataTable::mWidth(idGstar);
28 m2Res = mRes*mRes;
29 GamMRat = GammaRes / mRes;
30
31 // Overall coupling strength kappa * m_G*.
32 kappaMG = Settings::parm("ExtraDimensionsG*:kappaMG");
33
34 // Set pointer to particle properties and decay table.
35 gStarPtr = ParticleDataTable::particleDataPtr(idGstar);
36
37}
38
39//*********
40
41// Evaluate sigmaHat(sHat), part independent of incoming flavour.
42
43void Sigma1gg2GravitonStar::sigmaKin() {
44
45 // Incoming width for gluons.
46 double widthIn = pow2(kappaMG) * mH / (160. * M_PI);
47
48 // Set up Breit-Wigner. Width out only includes open channels.
49 double sigBW = 5. * M_PI/ ( pow2(sH - m2Res) + pow2(sH * GamMRat) );
50 double widthOut = gStarPtr->resWidthOpen(idGstar, mH);
51
52 // Modify cross section in wings of peak. Done.
53 sigma = widthIn * sigBW * widthOut * pow2(sH / m2Res);
54
55}
56
57//*********
58
59// Select identity, colour and anticolour.
60
61void Sigma1gg2GravitonStar::setIdColAcol() {
62
63 // Flavours trivial.
64 setId( 21, 21, idGstar);
65
66 // Colour flow topology.
67 setColAcol( 1, 2, 2, 1, 0, 0);
68
69}
70
71//*********
72
73// Evaluate weight for G* decay angle.
74
75double Sigma1gg2GravitonStar::weightDecay( Event& process, int iResBeg,
76 int iResEnd) {
77
78 // Identity of mother of decaying reseonance(s).
79 int idMother = process[process[iResBeg].mother1()].idAbs();
80
81 // For top decay hand over to standard routine.
82 if (idMother == 6)
83 return weightTopDecay( process, iResBeg, iResEnd);
84
85 // G* should sit in entry 5.
86 if (iResBeg != 5 || iResEnd != 5) return 1.;
87
88 // Phase space factors. Reconstruct decay angle.
89 double mr1 = pow2(process[6].m()) / sH;
90 double mr2 = pow2(process[7].m()) / sH;
91 double betaf = sqrtpos( pow2(1. - mr1 - mr2) - 4. * mr1 * mr2);
92 double cosThe = (process[3].p() - process[4].p())
93 * (process[7].p() - process[6].p()) / (sH * betaf);
94
95 // Default is isotropic decay.
96 double wt = 1.;
97
98 // Angular weight for g + g -> G* -> f + fbar.
99 if (process[6].idAbs() < 19) wt = 1. - pow4(cosThe);
100
101 // Angular weight for g + g -> G* -> g + g or gamma + gamma.
102 else if (process[6].id() == 21 || process[6].id() == 22)
103 wt = (1. + 6. * pow2(cosThe) + pow4(cosThe)) / 8.;
104
105 // Done.
106 return wt;
107
108}
109
110//**************************************************************************
111
112// Sigma1ffbar2GravitonStar class.
113// Cross section for f fbar -> G* (excited graviton state).
114
115//*********
116
117// Initialize process.
118
119void Sigma1ffbar2GravitonStar::initProc() {
120
121 // Store G* mass and width for propagator.
122 idGstar = 5000039;
123 mRes = ParticleDataTable::m0(idGstar);
124 GammaRes = ParticleDataTable::mWidth(idGstar);
125 m2Res = mRes*mRes;
126 GamMRat = GammaRes / mRes;
127
128 // Overall coupling strength kappa * m_G*.
129 kappaMG = Settings::parm("ExtraDimensionsG*:kappaMG");
130
131 // Set pointer to particle properties and decay table.
132 gStarPtr = ParticleDataTable::particleDataPtr(idGstar);
133
134}
135
136//*********
137
138// Evaluate sigmaHat(sHat), part independent of incoming flavour.
139
140void Sigma1ffbar2GravitonStar::sigmaKin() {
141
142 // Incoming width for fermions, disregarding colour factor.
143 double widthIn = pow2(kappaMG) * mH / (80. * M_PI);
144
145 // Set up Breit-Wigner. Width out only includes open channels.
146 double sigBW = 5. * M_PI/ ( pow2(sH - m2Res) + pow2(sH * GamMRat) );
147 double widthOut = gStarPtr->resWidthOpen(idGstar, mH);
148
149 // Modify cross section in wings of peak. Done.
150 sigma0 = widthIn * sigBW * widthOut * pow2(sH / m2Res);
151
152}
153
154//*********
155
156// Select identity, colour and anticolour.
157
158void Sigma1ffbar2GravitonStar::setIdColAcol() {
159
160 // Flavours trivial.
161 setId( id1, id2, idGstar);
162
163 // Colour flow topologies. Swap when antiquarks.
164 if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0);
165 else setColAcol( 0, 0, 0, 0, 0, 0);
166 if (id1 < 0) swapColAcol();
167
168}
169
170//*********
171
172// Evaluate weight for G* decay angle.
173
174double Sigma1ffbar2GravitonStar::weightDecay( Event& process, int iResBeg,
175 int iResEnd) {
176
177 // Identity of mother of decaying reseonance(s).
178 int idMother = process[process[iResBeg].mother1()].idAbs();
179
180 // For top decay hand over to standard routine.
181 if (idMother == 6)
182 return weightTopDecay( process, iResBeg, iResEnd);
183
184 // G* should sit in entry 5.
185 if (iResBeg != 5 || iResEnd != 5) return 1.;
186
187 // Phase space factors. Reconstruct decay angle.
188 double mr1 = pow2(process[6].m()) / sH;
189 double mr2 = pow2(process[7].m()) / sH;
190 double betaf = sqrtpos( pow2(1. - mr1 - mr2) - 4. * mr1 * mr2);
191 double cosThe = (process[3].p() - process[4].p())
192 * (process[7].p() - process[6].p()) / (sH * betaf);
193
194 // Default is isotropic decay.
195 double wt = 1.;
196
197 // Angular weight for f + fbar -> G* -> f + fbar.
198 if (process[6].idAbs() < 19)
199 wt = (1. - 3. * pow2(cosThe) + 4. * pow4(cosThe)) / 2.;
200
201 // Angular weight for f + fbar -> G* -> g + g or gamma + gamma.
202 else if (process[6].id() == 21 || process[6].id() == 22)
203 wt = 1. - pow4(cosThe);
204
205 // Done.
206 return wt;
207
208}
209
210//**************************************************************************
211
212// Sigma2gg2GravitonStarg class.
213// Cross section for g g -> G* g (excited graviton state).
214
215//*********
216
217// Initialize process.
218
219void Sigma2gg2GravitonStarg::initProc() {
220
221 // Store G* mass and width for propagator.
222 idGstar = 5000039;
223 mRes = ParticleDataTable::m0(idGstar);
224 GammaRes = ParticleDataTable::mWidth(idGstar);
225 m2Res = mRes*mRes;
226 GamMRat = GammaRes / mRes;
227
228 // Overall coupling strength kappa * m_G*.
229 kappaMG = Settings::parm("ExtraDimensionsG*:kappaMG");
230
231 // Secondary open width fraction.
232 openFrac = ParticleDataTable::resOpenFrac(idGstar);
233
234}
235
236//*********
237
238// Evaluate sigmaHat(sHat), part independent of incoming flavour.
239
240void Sigma2gg2GravitonStarg::sigmaKin() {
241
242 // Evaluate cross section. Secondary width for G*.
243 sigma = (3. * pow2(kappaMG) * alpS) / (32. * sH * s3)
244 * ( pow2(tH2 + tH * uH + uH2) / (sH2 * tH * uH)
245 + 2. * (tH2 / uH + uH2 / tH) / sH + 3. * (tH / uH + uH / tH)
246 + 2. * (sH / uH + sH/tH) + sH2 / (tH * uH) );
247 sigma *= openFrac;
248
249}
250
251//*********
252
253// Select identity, colour and anticolour.
254
255void Sigma2gg2GravitonStarg::setIdColAcol() {
256
257 // Flavours trivial.
258 setId( 21, 21, idGstar, 21);
259
260 // Colour flow topologies: random choice between two mirrors.
261 if (Rndm::flat() < 0.5) setColAcol( 1, 2, 2, 3, 0, 0, 1, 3);
262 else setColAcol( 1, 2, 3, 1, 0, 0, 3, 2);
263
264}
265
266//*********
267
268// Evaluate weight for decay angles: currently G* assumed isotropic.
269
270double Sigma2gg2GravitonStarg::weightDecay( Event& process, int iResBeg,
271 int iResEnd) {
272
273 // Identity of mother of decaying reseonance(s).
274 int idMother = process[process[iResBeg].mother1()].idAbs();
275
276 // For top decay hand over to standard routine.
277 if (idMother == 6)
278 return weightTopDecay( process, iResBeg, iResEnd);
279
280 // No equations for G* decay so assume isotropic.
281 return 1.;
282
283}
284
285//**************************************************************************
286
287// Sigma2qg2GravitonStarq class.
288// Cross section for q g -> G* q (excited graviton state).
289
290//*********
291
292// Initialize process.
293
294void Sigma2qg2GravitonStarq::initProc() {
295
296 // Store G* mass and width for propagator.
297 idGstar = 5000039;
298 mRes = ParticleDataTable::m0(idGstar);
299 GammaRes = ParticleDataTable::mWidth(idGstar);
300 m2Res = mRes*mRes;
301 GamMRat = GammaRes / mRes;
302
303 // Overall coupling strength kappa * m_G*.
304 kappaMG = Settings::parm("ExtraDimensionsG*:kappaMG");
305
306 // Secondary open width fraction.
307 openFrac = ParticleDataTable::resOpenFrac(idGstar);
308
309}
310
311//*********
312
313// Evaluate sigmaHat(sHat), part independent of incoming flavour.
314
315void Sigma2qg2GravitonStarq::sigmaKin() {
316
317 // Evaluate cross section. Secondary width for G*.
318 sigma = -(pow2(kappaMG) * alpS) / (192. * sH * s3)
319 * ( 4. * (sH2 + uH2) / (tH * sH) + 9. * (sH + uH) / sH + sH / uH
320 + uH2 / sH2 + 3. * tH * (4. + sH / uH + uH / sH) / sH
321 + 4. * tH2 * (1. / uH + 1. / sH) / sH + 2. * tH2 * tH / (uH * sH2) );
322 sigma *= openFrac;
323
324}
325
326//*********
327
328// Select identity, colour and anticolour.
329
330void Sigma2qg2GravitonStarq::setIdColAcol() {
331
332 // Flavour set up for q g -> H q.
333 int idq = (id2 == 21) ? id1 : id2;
334 setId( id1, id2, idGstar, idq);
335
336 // tH defined between f and f': must swap tHat <-> uHat if q g in.
337 swapTU = (id2 == 21);
338
339 // Colour flow topologies. Swap when antiquarks.
340 if (id2 == 21) setColAcol( 1, 0, 2, 1, 0, 0, 2, 0);
341 else setColAcol( 2, 1, 1, 0, 0, 0, 2, 0);
342 if (idq < 0) swapColAcol();
343
344}
345
346//*********
347
348// Evaluate weight for decay angles: currently G* assumed isotropic.
349
350double Sigma2qg2GravitonStarq::weightDecay( Event& process, int iResBeg,
351 int iResEnd) {
352
353 // Identity of mother of decaying reseonance(s).
354 int idMother = process[process[iResBeg].mother1()].idAbs();
355
356 // For top decay hand over to standard routine.
357 if (idMother == 6)
358 return weightTopDecay( process, iResBeg, iResEnd);
359
360 // No equations for G* decay so assume isotropic.
361 return 1.;
362
363}
364
365//**************************************************************************
366
367// Sigma2qqbar2GravitonStarg class.
368// Cross section for q qbar -> G* g (excited graviton state).
369
370//*********
371
372// Initialize process.
373
374void Sigma2qqbar2GravitonStarg::initProc() {
375
376 // Store G* mass and width for propagator.
377 idGstar = 5000039;
378 mRes = ParticleDataTable::m0(idGstar);
379 GammaRes = ParticleDataTable::mWidth(idGstar);
380 m2Res = mRes*mRes;
381 GamMRat = GammaRes / mRes;
382
383 // Overall coupling strength kappa * m_G*.
384 kappaMG = Settings::parm("ExtraDimensionsG*:kappaMG");
385
386 // Secondary open width fraction.
387 openFrac = ParticleDataTable::resOpenFrac(idGstar);
388
389}
390
391//*********
392
393// Evaluate sigmaHat(sHat), part independent of incoming flavour.
394
395void Sigma2qqbar2GravitonStarg::sigmaKin() {
396
397 // Evaluate cross section. Secondary width for G*.
398 sigma = (pow2(kappaMG) * alpS) / (72. * sH * s3)
399 * ( 4. * (tH2 + uH2) / sH2 + 9. * (tH + uH) / sH
400 + (tH2 / uH + uH2 / tH) / sH + 3. * (4. + tH / uH + uH/ tH)
401 + 4. * (sH / uH + sH / tH) + 2. * sH2 / (tH * uH) );
402 sigma *= openFrac;
403
404}
405
406//*********
407
408// Select identity, colour and anticolour.
409
410void Sigma2qqbar2GravitonStarg::setIdColAcol() {
411
412 // Flavours trivial.
413 setId( id1, id2, idGstar, 21);
414
415 // Colour flow topologies. Swap when antiquarks.
416 setColAcol( 1, 0, 0, 2, 0, 0, 1, 2);
417 if (id1 < 0) swapColAcol();
418
419}
420
421//*********
422
423// Evaluate weight for decay angles: currently G* assumed isotropic.
424
425double Sigma2qqbar2GravitonStarg::weightDecay( Event& process, int iResBeg,
426 int iResEnd) {
427
428 // Identity of mother of decaying reseonance(s).
429 int idMother = process[process[iResBeg].mother1()].idAbs();
430
431 // For top decay hand over to standard routine.
432 if (idMother == 6)
433 return weightTopDecay( process, iResBeg, iResEnd);
434
435 // No equations for G* decay so assume isotropic.
436 return 1.;
437
438}
439
440
441//**************************************************************************
442
443} // end namespace Pythia8