]>
Commit | Line | Data |
---|---|---|
9419eeef | 1 | // SigmaCompositeness.cc is a part of the PYTHIA event generator. |
2 | // Copyright (C) 2010 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 | // compositeness simulation classes. | |
8 | ||
9 | #include "SigmaCompositeness.h" | |
10 | ||
11 | namespace Pythia8 { | |
12 | ||
13 | //========================================================================== | |
14 | ||
15 | // Sigma1qg2qStar class. | |
16 | // Cross section for q g -> q^* (excited quark state). | |
17 | // Note: for simplicity decay is assumed isotropic. | |
18 | ||
19 | //-------------------------------------------------------------------------- | |
20 | ||
21 | // Initialize process. | |
22 | ||
23 | void Sigma1qg2qStar::initProc() { | |
24 | ||
25 | // Set up process properties from the chosen quark flavour. | |
26 | idRes = 4000000 + idq; | |
27 | codeSave = 4000 + idq; | |
28 | if (idq == 1) nameSave = "d g -> d^*"; | |
29 | else if (idq == 2) nameSave = "u g -> u^*"; | |
30 | else if (idq == 3) nameSave = "s g -> s^*"; | |
31 | else if (idq == 4) nameSave = "c g -> c^*"; | |
32 | else nameSave = "b g -> b^*"; | |
33 | ||
34 | // Store q* mass and width for propagator. | |
35 | mRes = particleDataPtr->m0(idRes); | |
36 | GammaRes = particleDataPtr->mWidth(idRes); | |
37 | m2Res = mRes*mRes; | |
38 | GamMRat = GammaRes / mRes; | |
39 | ||
40 | // Locally stored properties and couplings. | |
41 | Lambda = settingsPtr->parm("ExcitedFermion:Lambda"); | |
42 | coupFcol = settingsPtr->parm("ExcitedFermion:coupFcol"); | |
43 | ||
44 | // Set pointer to particle properties and decay table. | |
45 | qStarPtr = particleDataPtr->particleDataEntryPtr(idRes); | |
46 | ||
47 | } | |
48 | ||
49 | //-------------------------------------------------------------------------- | |
50 | ||
51 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
52 | ||
53 | void Sigma1qg2qStar::sigmaKin() { | |
54 | ||
55 | // Incoming width for correct quark. | |
56 | widthIn = pow3(mH) * alpS * pow2(coupFcol) / (3. * pow2(Lambda)); | |
57 | ||
58 | // Set up Breit-Wigner. | |
59 | sigBW = M_PI/ ( pow2(sH - m2Res) + pow2(sH * GamMRat) ); | |
60 | ||
61 | } | |
62 | ||
63 | //-------------------------------------------------------------------------- | |
64 | ||
65 | // Evaluate sigmaHat(sHat) for specific incoming flavours. | |
66 | ||
67 | double Sigma1qg2qStar::sigmaHat() { | |
68 | ||
69 | // Identify whether correct incoming flavours. | |
70 | int idqNow = (id2 == 21) ? id1 : id2; | |
71 | if (abs(idqNow) != idq) return 0.; | |
72 | ||
73 | // Outgoing width and total sigma. Done. | |
74 | return widthIn * sigBW * qStarPtr->resWidthOpen(idqNow, mH); | |
75 | ||
76 | } | |
77 | ||
78 | //-------------------------------------------------------------------------- | |
79 | ||
80 | // Select identity, colour and anticolour. | |
81 | ||
82 | void Sigma1qg2qStar::setIdColAcol() { | |
83 | ||
84 | // Flavours. | |
85 | int idqNow = (id2 == 21) ? id1 : id2; | |
86 | int idqStar = (idqNow > 0) ? idRes : -idRes; | |
87 | setId( id1, id2, idqStar); | |
88 | ||
89 | // Colour flow topology. | |
90 | if (id1 == idqNow) setColAcol( 1, 0, 2, 1, 2, 0); | |
91 | else setColAcol( 2, 1, 1, 0, 2, 0); | |
92 | if (idqNow < 0) swapColAcol(); | |
93 | ||
94 | } | |
95 | ||
96 | //-------------------------------------------------------------------------- | |
97 | ||
98 | // Evaluate weight for q* decay angle. | |
99 | ||
100 | double Sigma1qg2qStar::weightDecay( Event& process, int iResBeg, | |
101 | int iResEnd) { | |
102 | ||
103 | // q* should sit in entry 5. Sequential Z/W decay assumed isotropic. | |
104 | if (iResBeg != 5 || iResEnd != 5) return 1.; | |
105 | ||
106 | // Sign of asymmetry. | |
107 | int sideIn = (process[3].idAbs() < 20) ? 1 : 2; | |
108 | int sideOut = (process[6].idAbs() < 20) ? 1 : 2; | |
109 | double eps = (sideIn == sideOut) ? 1. : -1.; | |
110 | ||
111 | // Phase space factors. | |
112 | double mr1 = pow2(process[6].m()) / sH; | |
113 | double mr2 = pow2(process[7].m()) / sH; | |
114 | double betaf = sqrtpos( pow2(1. - mr1 - mr2) - 4. * mr1 * mr2); | |
115 | ||
116 | // Reconstruct decay angle. Default isotropic decay. | |
117 | double cosThe = (process[3].p() - process[4].p()) | |
118 | * (process[7].p() - process[6].p()) / (sH * betaf); | |
119 | double wt = 1.; | |
120 | double wtMax = 1.; | |
121 | ||
122 | // Decay q* -> q (g/gamma) or q (Z^0/W^+-). | |
123 | int idBoson = (sideOut == 1) ? process[7].idAbs() : process[6].idAbs(); | |
124 | if (idBoson == 21 || idBoson == 22) { | |
125 | wt = 1. + eps * cosThe; | |
126 | wtMax = 2.; | |
127 | } else if (idBoson == 23 || idBoson == 24) { | |
128 | double mrB = (sideOut == 1) ? mr2 : mr1; | |
129 | double ratB = (1. - 0.5 * mrB) / (1 + 0.5 * mrB); | |
130 | wt = 1. + eps * cosThe * ratB; | |
131 | wtMax = 1. + ratB; | |
132 | } | |
133 | ||
134 | // Done. | |
135 | return (wt / wtMax); | |
136 | ||
137 | } | |
138 | ||
139 | //========================================================================== | |
140 | ||
141 | // Sigma1lgm2lStar class. | |
142 | // Cross section for l gamma -> l^* (excited lepton state). | |
143 | // Note: for simplicity decay is assumed isotropic. | |
144 | ||
145 | //-------------------------------------------------------------------------- | |
146 | ||
147 | // Initialize process. | |
148 | ||
149 | void Sigma1lgm2lStar::initProc() { | |
150 | ||
151 | // Set up process properties from the chosen lepton flavour. | |
152 | idRes = 4000000 + idl; | |
153 | codeSave = 4000 + idl; | |
154 | if (idl == 11) nameSave = "e gamma -> e^*"; | |
155 | else if (idl == 13) nameSave = "mu gamma -> mu^*"; | |
156 | else nameSave = "tau gamma -> tau^*"; | |
157 | ||
158 | // Store l* mass and width for propagator. | |
159 | mRes = particleDataPtr->m0(idRes); | |
160 | GammaRes = particleDataPtr->mWidth(idRes); | |
161 | m2Res = mRes*mRes; | |
162 | GamMRat = GammaRes / mRes; | |
163 | ||
164 | // Locally stored properties and couplings. | |
165 | Lambda = settingsPtr->parm("ExcitedFermion:Lambda"); | |
166 | double coupF = settingsPtr->parm("ExcitedFermion:coupF"); | |
167 | double coupFp = settingsPtr->parm("ExcitedFermion:coupFprime"); | |
168 | coupChg = -0.5 * coupF - 0.5 * coupFp; | |
169 | ||
170 | // Set pointer to particle properties and decay table. | |
171 | qStarPtr = particleDataPtr->particleDataEntryPtr(idRes); | |
172 | ||
173 | } | |
174 | ||
175 | //-------------------------------------------------------------------------- | |
176 | ||
177 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
178 | ||
179 | void Sigma1lgm2lStar::sigmaKin() { | |
180 | ||
181 | // Incoming width for correct lepton. | |
182 | widthIn = pow3(mH) * alpEM * pow2(coupChg) / pow2(Lambda); | |
183 | ||
184 | // Set up Breit-Wigner. | |
185 | sigBW = M_PI/ ( pow2(sH - m2Res) + pow2(sH * GamMRat) ); | |
186 | ||
187 | } | |
188 | ||
189 | //-------------------------------------------------------------------------- | |
190 | ||
191 | // Evaluate sigmaHat(sHat) for specific incoming flavours. | |
192 | ||
193 | double Sigma1lgm2lStar::sigmaHat() { | |
194 | ||
195 | // Identify whether correct incoming flavours. | |
196 | int idlNow = (id2 == 22) ? id1 : id2; | |
197 | if (abs(idlNow) != idl) return 0.; | |
198 | ||
199 | // Outgoing width and total sigma. Done. | |
200 | return widthIn * sigBW * qStarPtr->resWidthOpen(idlNow, mH); | |
201 | ||
202 | } | |
203 | ||
204 | //-------------------------------------------------------------------------- | |
205 | ||
206 | // Select identity, colour and anticolour. | |
207 | ||
208 | void Sigma1lgm2lStar::setIdColAcol() { | |
209 | ||
210 | // Flavours. | |
211 | int idlNow = (id2 == 22) ? id1 : id2; | |
212 | int idlStar = (idlNow > 0) ? idRes : -idRes; | |
213 | setId( id1, id2, idlStar); | |
214 | ||
215 | // No colour flow. | |
216 | setColAcol( 0, 0, 0, 0, 0, 0); | |
217 | ||
218 | } | |
219 | ||
220 | //-------------------------------------------------------------------------- | |
221 | ||
222 | // Evaluate weight for l* decay angle. | |
223 | ||
224 | double Sigma1lgm2lStar::weightDecay( Event& process, int iResBeg, | |
225 | int iResEnd) { | |
226 | ||
227 | // l* should sit in entry 5. Sequential Z/W decay assumed isotropic. | |
228 | if (iResBeg != 5 || iResEnd != 5) return 1.; | |
229 | ||
230 | // Sign of asymmetry. | |
231 | int sideIn = (process[3].idAbs() < 20) ? 1 : 2; | |
232 | int sideOut = (process[6].idAbs() < 20) ? 1 : 2; | |
233 | double eps = (sideIn == sideOut) ? 1. : -1.; | |
234 | ||
235 | // Phase space factors. | |
236 | double mr1 = pow2(process[6].m()) / sH; | |
237 | double mr2 = pow2(process[7].m()) / sH; | |
238 | double betaf = sqrtpos( pow2(1. - mr1 - mr2) - 4. * mr1 * mr2); | |
239 | ||
240 | // Reconstruct decay angle. Default isotropic decay. | |
241 | double cosThe = (process[3].p() - process[4].p()) | |
242 | * (process[7].p() - process[6].p()) / (sH * betaf); | |
243 | double wt = 1.; | |
244 | double wtMax = 1.; | |
245 | ||
246 | // Decay l* -> l gamma or l (Z^0/W^+-). | |
247 | int idBoson = (sideOut == 1) ? process[7].idAbs() : process[6].idAbs(); | |
248 | if (idBoson == 22) { | |
249 | wt = 1. + eps * cosThe; | |
250 | wtMax = 2.; | |
251 | } else if (idBoson == 23 || idBoson == 24) { | |
252 | double mrB = (sideOut == 1) ? mr2 : mr1; | |
253 | double ratB = (1. - 0.5 * mrB) / (1 + 0.5 * mrB); | |
254 | wt = 1. + eps * cosThe * ratB; | |
255 | wtMax = 1. + ratB; | |
256 | } | |
257 | ||
258 | // Done. | |
259 | return (wt / wtMax); | |
260 | ||
261 | } | |
262 | ||
263 | //========================================================================== | |
264 | ||
265 | // Sigma2qq2qStarq class. | |
266 | // Cross section for q q' -> q^* q' (excited quark state). | |
267 | // Note: for simplicity decay is assumed isotropic. | |
268 | ||
269 | //-------------------------------------------------------------------------- | |
270 | ||
271 | // Initialize process. | |
272 | ||
273 | void Sigma2qq2qStarq::initProc() { | |
274 | ||
275 | // Set up process properties from the chosen quark flavour. | |
276 | idRes = 4000000 + idq; | |
277 | codeSave = 4020 + idq; | |
278 | if (idq == 1) nameSave = "q q -> d^* q"; | |
279 | else if (idq == 2) nameSave = "q q -> u^* q"; | |
280 | else if (idq == 3) nameSave = "q q -> s^* q"; | |
281 | else if (idq == 4) nameSave = "q q -> c^* q"; | |
282 | else nameSave = "q q -> b^* q"; | |
283 | ||
284 | // Locally stored properties and couplings. | |
285 | Lambda = settingsPtr->parm("ExcitedFermion:Lambda"); | |
286 | preFac = M_PI / pow4(Lambda); | |
287 | ||
288 | // Secondary open width fractions. | |
289 | openFracPos = particleDataPtr->resOpenFrac( idRes); | |
290 | openFracNeg = particleDataPtr->resOpenFrac(-idRes); | |
291 | ||
292 | } | |
293 | ||
294 | //-------------------------------------------------------------------------- | |
295 | ||
296 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
297 | ||
298 | void Sigma2qq2qStarq::sigmaKin() { | |
299 | ||
300 | // Two possible expressions, for like or unlike sign. | |
301 | sigmaA = preFac * (1. - s3 / sH); | |
302 | sigmaB = preFac * (-uH) * (sH + tH) / sH2; | |
303 | ||
304 | } | |
305 | ||
306 | //-------------------------------------------------------------------------- | |
307 | ||
308 | // Evaluate sigmaHat(sHat) for specific incoming flavours. | |
309 | ||
310 | double Sigma2qq2qStarq::sigmaHat() { | |
311 | ||
312 | // Identify different allowed incoming flavour combinations. | |
313 | int id1Abs = abs(id1); | |
314 | int id2Abs = abs(id2); | |
315 | double open1 = (id1 > 0) ? openFracPos : openFracNeg; | |
316 | double open2 = (id2 > 0) ? openFracPos : openFracNeg; | |
317 | double sigma = 0.; | |
318 | if (id1 * id2 > 0) { | |
319 | if (id1Abs == idq) sigma += (4./3.) * sigmaA * open1; | |
320 | if (id2Abs == idq) sigma += (4./3.) * sigmaA * open2; | |
321 | } else if (id1Abs == idq && id2 == -id1) | |
322 | sigma = (8./3.) * sigmaB * (open1 + open2); | |
323 | else if (id2 == -id1) sigma = sigmaB * (open1 + open2); | |
324 | else if (id1Abs == idq) sigma = sigmaB * open1; | |
325 | else if (id2Abs == idq) sigma = sigmaB * open2; | |
326 | ||
327 | // Done. | |
328 | return sigma; | |
329 | ||
330 | } | |
331 | ||
332 | //-------------------------------------------------------------------------- | |
333 | ||
334 | // Select identity, colour and anticolour. | |
335 | ||
336 | void Sigma2qq2qStarq::setIdColAcol() { | |
337 | ||
338 | // Flavours: either side may have been excited. | |
339 | double open1 = 0.; | |
340 | double open2 = 0.; | |
341 | if (abs(id1) == idq) open1 = (id1 > 0) ? openFracPos : openFracNeg; | |
342 | if (abs(id2) == idq) open2 = (id2 > 0) ? openFracPos : openFracNeg; | |
343 | if (open1 == 0. && open2 == 0.) { | |
344 | open1 = (id1 > 0) ? openFracPos : openFracNeg; | |
345 | open2 = (id2 > 0) ? openFracPos : openFracNeg; | |
346 | } | |
347 | bool excite1 = (open1 > 0.); | |
348 | if (open1 > 0. && open2 > 0.) excite1 | |
349 | = (rndmPtr->flat() * (open1 + open2) < open1); | |
350 | ||
351 | // Always excited quark in slot 3 so colour flow flipped or not. | |
352 | if (excite1) { | |
353 | id3 = (id1 > 0) ? idq : -idq; | |
354 | id4 = id2; | |
355 | if (id1 * id2 > 0) setColAcol( 1, 0, 2, 0, 1, 0, 2, 0); | |
356 | else setColAcol( 1, 0, 0, 2, 1, 0, 0, 2); | |
357 | if (id1 < 0) swapColAcol(); | |
358 | } else { | |
359 | id3 = (id2 > 0) ? idq : -idq; | |
360 | id4 = id1; | |
361 | swapTU = true; | |
362 | if (id1 * id2 > 0) setColAcol( 1, 0, 2, 0, 2, 0, 1, 0); | |
363 | else setColAcol( 1, 0, 0, 2, 0, 2, 1, 0); | |
364 | if (id1 < 0) swapColAcol(); | |
365 | } | |
366 | setId( id1, id2, id3, id4); | |
367 | ||
368 | } | |
369 | ||
370 | //========================================================================== | |
371 | ||
372 | // Sigma2qqbar2lStarlbar class. | |
373 | // Cross section for q qbar -> l^* lbar (excited lepton state). | |
374 | // Note: for simplicity decay is assumed isotropic. | |
375 | ||
376 | //-------------------------------------------------------------------------- | |
377 | ||
378 | // Initialize process. | |
379 | ||
380 | void Sigma2qqbar2lStarlbar::initProc() { | |
381 | ||
382 | // Set up process properties from the chosen lepton flavour. | |
383 | idRes = 4000000 + idl; | |
384 | codeSave = 4020 + idl; | |
385 | if (idl == 11) nameSave = "q qbar -> e^*+- e^-+"; | |
386 | else if (idl == 12) nameSave = "q qbar -> nu_e^* nu_ebar"; | |
387 | else if (idl == 13) nameSave = "q qbar -> mu^*+- mu^-+"; | |
388 | else if (idl == 14) nameSave = "q qbar -> nu_mu^* nu_mubar"; | |
389 | else if (idl == 15) nameSave = "q qbar -> tau^*+- tau^-+"; | |
390 | else nameSave = "q qbar -> nu_tau^* nu_taubar"; | |
391 | ||
392 | // Secondary open width fractions. | |
393 | openFracPos = particleDataPtr->resOpenFrac( idRes); | |
394 | openFracNeg = particleDataPtr->resOpenFrac(-idRes); | |
395 | ||
396 | // Locally stored properties and couplings. | |
397 | Lambda = settingsPtr->parm("ExcitedFermion:Lambda"); | |
398 | preFac = (M_PI / pow4(Lambda)) * (openFracPos + openFracNeg) / 3.; | |
399 | ||
400 | } | |
401 | ||
402 | //-------------------------------------------------------------------------- | |
403 | ||
404 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
405 | ||
406 | void Sigma2qqbar2lStarlbar::sigmaKin() { | |
407 | ||
408 | // Only one possible expressions | |
409 | sigma = preFac * (-uH) * (sH + tH) / sH2; | |
410 | ||
411 | } | |
412 | ||
413 | //-------------------------------------------------------------------------- | |
414 | ||
415 | // Select identity, colour and anticolour. | |
416 | ||
417 | void Sigma2qqbar2lStarlbar::setIdColAcol() { | |
418 | ||
419 | // Flavours: either lepton or antilepton may be excited. | |
420 | if (rndmPtr->flat() * (openFracPos + openFracNeg) < openFracPos) { | |
421 | setId( id1, id2, idRes, -idl); | |
422 | if (id1 < 0) swapTU = true; | |
423 | } else { | |
424 | setId( id1, id2, -idRes, idl); | |
425 | if (id1 > 0) swapTU = true; | |
426 | } | |
427 | ||
428 | // Colour flow trivial. | |
429 | if (id1 > 0) setColAcol( 1, 0, 0, 1, 0, 0, 0, 0); | |
430 | else setColAcol( 0, 1, 1, 0, 0, 0, 0, 0); | |
431 | ||
432 | } | |
433 | ||
434 | //========================================================================== | |
435 | ||
436 | // Sigma2QCqq2qq class. | |
437 | // Cross section for q q -> q q (quark contact interactions). | |
438 | ||
439 | //-------------------------------------------------------------------------- | |
440 | ||
441 | // Initialize process. | |
442 | ||
443 | void Sigma2QCqq2qq::initProc() { | |
444 | ||
445 | m_Lambda2 = settingsPtr->parm("ContactInteractions:Lambda"); | |
446 | m_etaLL = settingsPtr->mode("ContactInteractions:etaLL");; | |
447 | m_etaRR = settingsPtr->mode("ContactInteractions:etaRR");; | |
448 | m_etaLR = settingsPtr->mode("ContactInteractions:etaLR");; | |
449 | m_Lambda2 *= m_Lambda2; | |
450 | ||
451 | } | |
452 | ||
453 | //-------------------------------------------------------------------------- | |
454 | ||
455 | // Evaluate d(sigmaHat)/d(tHat), part independent of incoming flavour. | |
456 | ||
457 | void Sigma2QCqq2qq::sigmaKin() { | |
458 | ||
459 | // Calculate kinematics dependence for different terms. | |
460 | sigT = (4./9.) * (sH2 + uH2) / tH2; | |
461 | sigU = (4./9.) * (sH2 + tH2) / uH2; | |
462 | sigTU = - (8./27.) * sH2 / (tH * uH); | |
463 | sigST = - (8./27.) * uH2 / (sH * tH); | |
464 | ||
465 | sigQCSTU = sH2 * (1 / tH + 1 / uH); | |
466 | sigQCUTS = uH2 * (1 / tH + 1 / sH); | |
467 | ||
468 | } | |
469 | ||
470 | //-------------------------------------------------------------------------- | |
471 | ||
472 | // Evaluate d(sigmaHat)/d(tHat), including incoming flavour dependence. | |
473 | ||
474 | double Sigma2QCqq2qq::sigmaHat() { | |
475 | ||
476 | // Terms from QC contact interactions. | |
477 | double sigQCLL = 0; | |
478 | double sigQCRR = 0; | |
479 | double sigQCLR = 0; | |
480 | ||
481 | // Combine cross section terms; factor 1/2 when identical quarks. | |
482 | // q q -> q q | |
483 | if (id2 == id1) { | |
484 | ||
485 | sigSum = 0.5 * (sigT + sigU + sigTU); // SM terms. | |
486 | ||
487 | // Contact terms. | |
488 | sigQCLL = (8./9.) * alpS * (m_etaLL/m_Lambda2) * sigQCSTU | |
489 | + (8./3.) * pow2(m_etaLL/m_Lambda2) * sH2; | |
490 | sigQCRR = (8./9.) * alpS * (m_etaRR/m_Lambda2) * sigQCSTU | |
491 | + (8./3.) * pow2(m_etaRR/m_Lambda2) * sH2; | |
492 | sigQCLR = 2. * (uH2 + tH2) * pow2(m_etaLR/m_Lambda2); | |
493 | ||
494 | sigQCLL /= 2; | |
495 | sigQCRR /= 2; | |
496 | sigQCLR /= 2; | |
497 | ||
498 | // q qbar -> q qbar, without pure s-channel term. | |
499 | } else if (id2 == -id1) { | |
500 | ||
501 | sigSum = sigT + sigST; // SM terms. | |
502 | ||
503 | // Contact terms, minus the terms included in qqbar2qqbar. | |
504 | sigQCLL = (8./9.) * alpS * (m_etaLL/m_Lambda2) * sigQCUTS | |
505 | + (5./3.) * pow2(m_etaLL/m_Lambda2) * uH2; | |
506 | sigQCRR = (8./9.) * alpS * (m_etaRR/m_Lambda2) * sigQCUTS | |
507 | + (5./3.) * pow2(m_etaRR/m_Lambda2) * uH2; | |
508 | sigQCLR = 2. * sH2 * pow2(m_etaLR/m_Lambda2); | |
509 | ||
510 | // q q' -> q q' or q qbar' -> q qbar' | |
511 | } else { | |
512 | ||
513 | sigSum = sigT; // SM terms. | |
514 | ||
515 | // Contact terms. | |
516 | if (id1 * id2 > 0) { | |
517 | sigQCLL = pow2(m_etaLL/m_Lambda2) * sH2; | |
518 | sigQCRR = pow2(m_etaRR/m_Lambda2) * sH2; | |
519 | sigQCLR = 2 * pow2(m_etaLR/m_Lambda2) * uH2; | |
520 | } else { | |
521 | sigQCLL = pow2(m_etaLL/m_Lambda2) * uH2; | |
522 | sigQCRR = pow2(m_etaRR/m_Lambda2) * uH2; | |
523 | sigQCLR = 2 * pow2(m_etaLR/m_Lambda2) * sH2; | |
524 | } | |
525 | } | |
526 | ||
527 | // Answer. | |
528 | double sigma = (M_PI/sH2) * ( pow2(alpS) * sigSum | |
529 | + sigQCLL + sigQCRR + sigQCLR ); | |
530 | return sigma; | |
531 | ||
532 | } | |
533 | ||
534 | //-------------------------------------------------------------------------- | |
535 | ||
536 | // Select identity, colour and anticolour. | |
537 | ||
538 | void Sigma2QCqq2qq::setIdColAcol() { | |
539 | ||
540 | // Outgoing = incoming flavours. | |
541 | setId( id1, id2, id1, id2); | |
542 | ||
543 | // Colour flow topologies. Swap when antiquarks. | |
544 | if (id1 * id2 > 0) setColAcol( 1, 0, 2, 0, 2, 0, 1, 0); | |
545 | else setColAcol( 1, 0, 0, 1, 2, 0, 0, 2); | |
546 | if (id2 == id1 && (sigT + sigU) * rndmPtr->flat() > sigT) | |
547 | setColAcol( 1, 0, 2, 0, 1, 0, 2, 0); | |
548 | if (id1 < 0) swapColAcol(); | |
549 | ||
550 | } | |
551 | ||
552 | //========================================================================== | |
553 | ||
554 | // Sigma2QCqqbar2qqbar class. | |
555 | // Cross section for q qbar -> q' qbar' (quark contact interactions). | |
556 | ||
557 | //-------------------------------------------------------------------------- | |
558 | ||
559 | // Initialize process. | |
560 | ||
561 | void Sigma2QCqqbar2qqbar::initProc() { | |
562 | ||
563 | m_nQuarkNew = settingsPtr->mode("ContactInteractions:nQuarkNew"); | |
564 | m_Lambda2 = settingsPtr->parm("ContactInteractions:Lambda"); | |
565 | m_etaLL = settingsPtr->mode("ContactInteractions:etaLL"); | |
566 | m_etaRR = settingsPtr->mode("ContactInteractions:etaRR"); | |
567 | m_etaLR = settingsPtr->mode("ContactInteractions:etaLR"); | |
568 | m_Lambda2 *= m_Lambda2; | |
569 | ||
570 | } | |
571 | ||
572 | //-------------------------------------------------------------------------- | |
573 | ||
574 | // Evaluate d(sigmaHat)/d(tHat) - no incoming flavour dependence. | |
575 | ||
576 | void Sigma2QCqqbar2qqbar::sigmaKin() { | |
577 | ||
578 | // Pick new flavour. | |
579 | idNew = 1 + int( m_nQuarkNew * rndmPtr->flat() ); | |
580 | mNew = particleDataPtr->m0(idNew); | |
581 | m2New = mNew*mNew; | |
582 | ||
583 | // Calculate kinematics dependence. | |
584 | double sigQC = 0.; | |
585 | sigS = 0.; | |
586 | if (sH > 4. * m2New) { | |
587 | sigS = (4./9.) * (tH2 + uH2) / sH2; | |
588 | sigQC = pow2(m_etaLL/m_Lambda2) * uH2 | |
589 | + pow2(m_etaRR/m_Lambda2) * uH2 | |
590 | + 2 * pow2(m_etaLR/m_Lambda2) * tH2; | |
591 | } | |
592 | ||
593 | // Answer is proportional to number of outgoing flavours. | |
594 | sigma = (M_PI / sH2) * m_nQuarkNew * ( pow2(alpS) * sigS + sigQC); | |
595 | ||
596 | } | |
597 | ||
598 | //-------------------------------------------------------------------------- | |
599 | ||
600 | // Select identity, colour and anticolour. | |
601 | ||
602 | void Sigma2QCqqbar2qqbar::setIdColAcol() { | |
603 | ||
604 | // Set outgoing flavours ones. | |
605 | id3 = (id1 > 0) ? idNew : -idNew; | |
606 | setId( id1, id2, id3, -id3); | |
607 | ||
608 | // Colour flow topologies. Swap when antiquarks. | |
609 | setColAcol( 1, 0, 0, 2, 1, 0, 0, 2); | |
610 | if (id1 < 0) swapColAcol(); | |
611 | ||
612 | } | |
613 | ||
614 | //========================================================================== | |
615 | ||
616 | } // end namespace Pythia8 |