]>
Commit | Line | Data |
---|---|---|
9419eeef | 1 | // SigmaHiggs.cc is a part of the PYTHIA event generator. |
2 | // Copyright (C) 2010 Torbjorn Sjostrand. | |
3 | // Part of code written by Marc Montull, CERN summer student 2007. | |
4 | // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details. | |
5 | // Please respect the MCnet Guidelines, see GUIDELINES for details. | |
6 | // Function definitions (not found in the header) for the | |
7 | // Higgs simulation classes. | |
8 | ||
9 | #include "SigmaHiggs.h" | |
10 | ||
11 | namespace Pythia8 { | |
12 | ||
13 | //========================================================================== | |
14 | ||
15 | // Sigma1ffbar2H class. | |
16 | // Cross section for f fbar -> H0 , H1, H2 or A3. | |
17 | // (f is quark or lepton, H0 SM Higgs and H1, H2, A3 BSM Higgses ). | |
18 | ||
19 | //-------------------------------------------------------------------------- | |
20 | ||
21 | // Initialize process. | |
22 | ||
23 | void Sigma1ffbar2H::initProc() { | |
24 | ||
25 | // Properties specific to Higgs state. | |
26 | if (higgsType == 0) { | |
27 | nameSave = "f fbar -> H (SM)"; | |
28 | codeSave = 901; | |
29 | idRes = 25; | |
30 | } | |
31 | else if (higgsType == 1) { | |
32 | nameSave = "f fbar -> h0(H1)"; | |
33 | codeSave = 1001; | |
34 | idRes = 25; | |
35 | } | |
36 | else if (higgsType == 2) { | |
37 | nameSave = "f fbar -> H0(H2)"; | |
38 | codeSave = 1021; | |
39 | idRes = 35; | |
40 | } | |
41 | else if (higgsType == 3) { | |
42 | nameSave = "f fbar -> A0(A3)"; | |
43 | codeSave = 1041; | |
44 | idRes = 36; | |
45 | } | |
46 | ||
47 | // Find pointer to H0, H1, H2 or A3 depending on the value of idRes. | |
48 | HResPtr = particleDataPtr->particleDataEntryPtr(idRes); | |
49 | ||
50 | // Store H0, H1, H2 or A3 mass and width for propagator. | |
51 | mRes = HResPtr->m0(); | |
52 | GammaRes = HResPtr->mWidth(); | |
53 | m2Res = mRes*mRes; | |
54 | GamMRat = GammaRes / mRes; | |
55 | ||
56 | } | |
57 | ||
58 | ||
59 | //-------------------------------------------------------------------------- | |
60 | ||
61 | // Evaluate d(sigmaHat)/d(tHat), part independent of incoming flavour. | |
62 | ||
63 | void Sigma1ffbar2H::sigmaKin() { | |
64 | ||
65 | // Set up Breit-Wigner. | |
66 | double width = HResPtr->resWidth(idRes, mH); | |
67 | sigBW = 4. * M_PI/ ( pow2(sH - m2Res) + pow2(mH * width) ); | |
68 | ||
69 | // Width out only includes open channels. | |
70 | widthOut = width * HResPtr->resOpenFrac(idRes); | |
71 | ||
72 | } | |
73 | ||
74 | //-------------------------------------------------------------------------- | |
75 | ||
76 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
77 | ||
78 | double Sigma1ffbar2H::sigmaHat() { | |
79 | ||
80 | // Calculate mass-dependent incoming width, including colour factor. | |
81 | int idAbs = abs(id1); | |
82 | double widthIn = HResPtr->resWidthChan( mH, idAbs, -idAbs); | |
83 | if (idAbs < 9) widthIn /= 9.; | |
84 | ||
85 | // Done. | |
86 | return widthIn * sigBW * widthOut; | |
87 | ||
88 | } | |
89 | ||
90 | //-------------------------------------------------------------------------- | |
91 | ||
92 | // Select identity, colour and anticolour. | |
93 | ||
94 | void Sigma1ffbar2H::setIdColAcol() { | |
95 | ||
96 | // Flavours trivial. | |
97 | setId( id1, id2, idRes); | |
98 | ||
99 | // Colour flow topologies. Swap when antiquarks. | |
100 | if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0); | |
101 | else setColAcol( 0, 0, 0, 0, 0, 0); | |
102 | if (id1 < 0) swapColAcol(); | |
103 | ||
104 | } | |
105 | ||
106 | //-------------------------------------------------------------------------- | |
107 | ||
108 | // Evaluate weight for decay angles. | |
109 | ||
110 | double Sigma1ffbar2H::weightDecay( Event& process, int iResBeg, | |
111 | int iResEnd) { | |
112 | ||
113 | // Identity of mother of decaying reseonance(s). | |
114 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
115 | ||
116 | // For Higgs decay hand over to standard routine. | |
117 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
118 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
119 | ||
120 | // For top decay hand over to standard routine. | |
121 | if (idMother == 6) | |
122 | return weightTopDecay( process, iResBeg, iResEnd); | |
123 | ||
124 | // Else done. | |
125 | return 1.; | |
126 | ||
127 | } | |
128 | ||
129 | //========================================================================== | |
130 | ||
131 | // Sigma1gg2H class. | |
132 | // Cross section for g g -> H0, H1, H2 or A3 (H0 SM Higgs, H1, H2, A3 BSM). | |
133 | ||
134 | //-------------------------------------------------------------------------- | |
135 | ||
136 | // Initialize process. | |
137 | ||
138 | void Sigma1gg2H::initProc() { | |
139 | ||
140 | // Properties specific to Higgs state. | |
141 | if (higgsType == 0) { | |
142 | nameSave = "g g -> H (SM)"; | |
143 | codeSave = 902; | |
144 | idRes = 25; | |
145 | } | |
146 | else if (higgsType == 1) { | |
147 | nameSave = "g g -> h0(H1)"; | |
148 | codeSave = 1002; | |
149 | idRes = 25; | |
150 | } | |
151 | else if (higgsType == 2) { | |
152 | nameSave = "g g -> H0(H2)"; | |
153 | codeSave = 1022; | |
154 | idRes = 35; | |
155 | } | |
156 | else if (higgsType == 3) { | |
157 | nameSave = "g g -> A0(A3)"; | |
158 | codeSave = 1042; | |
159 | idRes = 36; | |
160 | } | |
161 | ||
162 | // Find pointer to H0, H1, H2 or A3 depending on idRes. | |
163 | HResPtr = particleDataPtr->particleDataEntryPtr(idRes); | |
164 | ||
165 | // Store H0, H1, H2 or A3 mass and width for propagator. | |
166 | mRes = HResPtr->m0(); | |
167 | GammaRes = HResPtr->mWidth(); | |
168 | m2Res = mRes*mRes; | |
169 | GamMRat = GammaRes / mRes; | |
170 | ||
171 | } | |
172 | ||
173 | //-------------------------------------------------------------------------- | |
174 | ||
175 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
176 | ||
177 | void Sigma1gg2H::sigmaKin() { | |
178 | ||
179 | // Incoming width for gluons, gives colour factor of 1/8 * 1/8. | |
180 | double widthIn = HResPtr->resWidthChan( mH, 21, 21) / 64.; | |
181 | ||
182 | // Set up Breit-Wigner. | |
183 | double width = HResPtr->resWidth(idRes, mH); | |
184 | double sigBW = 8. * M_PI/ ( pow2(sH - m2Res) + pow2(mH * width) ); | |
185 | ||
186 | // Width out only includes open channels. | |
187 | double widthOut = width * HResPtr->resOpenFrac(idRes); | |
188 | ||
189 | // Done. | |
190 | sigma = widthIn * sigBW * widthOut; | |
191 | ||
192 | } | |
193 | ||
194 | //-------------------------------------------------------------------------- | |
195 | ||
196 | // Select identity, colour and anticolour. | |
197 | ||
198 | void Sigma1gg2H::setIdColAcol() { | |
199 | ||
200 | // Flavours trivial. | |
201 | setId( 21, 21, idRes); | |
202 | ||
203 | // Colour flow topology. | |
204 | setColAcol( 1, 2, 2, 1, 0, 0); | |
205 | ||
206 | } | |
207 | ||
208 | //-------------------------------------------------------------------------- | |
209 | ||
210 | // Evaluate weight for decay angles. | |
211 | ||
212 | double Sigma1gg2H::weightDecay( Event& process, int iResBeg, | |
213 | int iResEnd) { | |
214 | ||
215 | // Identity of mother of decaying reseonance(s). | |
216 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
217 | ||
218 | // For Higgs decay hand over to standard routine. | |
219 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
220 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
221 | ||
222 | // For top decay hand over to standard routine. | |
223 | if (idMother == 6) | |
224 | return weightTopDecay( process, iResBeg, iResEnd); | |
225 | ||
226 | // Else done. | |
227 | return 1.; | |
228 | ||
229 | } | |
230 | ||
231 | //========================================================================== | |
232 | ||
233 | // Sigma1gmgm2H class. | |
234 | // Cross section for gamma gamma -> H0, H1, H2 or H3. | |
235 | // (H0 SM Higgs, H1, H2 and A3 BSM Higgses). | |
236 | ||
237 | //-------------------------------------------------------------------------- | |
238 | ||
239 | // Initialize process. | |
240 | ||
241 | void Sigma1gmgm2H::initProc() { | |
242 | ||
243 | // Properties specific to Higgs state. | |
244 | if (higgsType == 0) { | |
245 | nameSave = "gamma gamma -> H (SM)"; | |
246 | codeSave = 903; | |
247 | idRes = 25; | |
248 | } | |
249 | else if (higgsType == 1) { | |
250 | nameSave = "gamma gamma -> h0(H1)"; | |
251 | codeSave = 1003; | |
252 | idRes = 25; | |
253 | } | |
254 | else if (higgsType == 2) { | |
255 | nameSave = "gamma gamma -> H0(H2)"; | |
256 | codeSave = 1023; | |
257 | idRes = 35; | |
258 | } | |
259 | else if (higgsType == 3) { | |
260 | nameSave = "gamma gamma -> A0(A3)"; | |
261 | codeSave = 1043; | |
262 | idRes = 36; | |
263 | } | |
264 | ||
265 | // Find pointer to H0, H1, H2 or A3. | |
266 | HResPtr = particleDataPtr->particleDataEntryPtr(idRes); | |
267 | ||
268 | // Store H0, H1, H2 or A3 mass and width for propagator. | |
269 | mRes = HResPtr->m0(); | |
270 | GammaRes = HResPtr->mWidth(); | |
271 | m2Res = mRes*mRes; | |
272 | GamMRat = GammaRes / mRes; | |
273 | ||
274 | } | |
275 | ||
276 | //-------------------------------------------------------------------------- | |
277 | ||
278 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
279 | ||
280 | void Sigma1gmgm2H::sigmaKin() { | |
281 | ||
282 | // Incoming width for photons. | |
283 | double widthIn = HResPtr->resWidthChan( mH, 22, 22); | |
284 | ||
285 | // Set up Breit-Wigner. | |
286 | double width = HResPtr->resWidth(idRes, mH); | |
287 | double sigBW = 8. * M_PI/ ( pow2(sH - m2Res) + pow2(mH * width) ); | |
288 | ||
289 | // Width out only includes open channels. | |
290 | double widthOut = width * HResPtr->resOpenFrac(idRes); | |
291 | ||
292 | // Done. | |
293 | sigma = widthIn * sigBW * widthOut; | |
294 | ||
295 | } | |
296 | ||
297 | //-------------------------------------------------------------------------- | |
298 | ||
299 | // Select identity, colour and anticolour. | |
300 | ||
301 | void Sigma1gmgm2H::setIdColAcol() { | |
302 | ||
303 | // Flavours trivial. | |
304 | setId( 22, 22, idRes); | |
305 | ||
306 | // Colour flow trivial. | |
307 | setColAcol( 0, 0, 0, 0, 0, 0); | |
308 | ||
309 | } | |
310 | ||
311 | //-------------------------------------------------------------------------- | |
312 | ||
313 | // Evaluate weight for decay angles. | |
314 | ||
315 | double Sigma1gmgm2H::weightDecay( Event& process, int iResBeg, | |
316 | int iResEnd) { | |
317 | ||
318 | // Identity of mother of decaying reseonance(s). | |
319 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
320 | ||
321 | // For Higgs decay hand over to standard routine. | |
322 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
323 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
324 | ||
325 | // For top decay hand over to standard routine. | |
326 | if (idMother == 6) | |
327 | return weightTopDecay( process, iResBeg, iResEnd); | |
328 | ||
329 | // Else done. | |
330 | return 1.; | |
331 | ||
332 | } | |
333 | ||
334 | //========================================================================== | |
335 | ||
336 | // Sigma2ffbar2HZ class. | |
337 | // Cross section for f fbar -> H0 Z0, H1 Z0, H2 Z0 or A3 Z0. | |
338 | // (H0 SM Higgs, H1, H2 and A3 BSM Higgses). | |
339 | ||
340 | //-------------------------------------------------------------------------- | |
341 | ||
342 | // Initialize process. | |
343 | ||
344 | void Sigma2ffbar2HZ::initProc() { | |
345 | ||
346 | // Properties specific to Higgs state. | |
347 | if (higgsType == 0) { | |
348 | nameSave = "f fbar -> H0 Z0 (SM)"; | |
349 | codeSave = 904; | |
350 | idRes = 25; | |
351 | coup2Z = 1.; | |
352 | } | |
353 | else if (higgsType == 1) { | |
354 | nameSave = "f fbar -> h0(H1) Z0"; | |
355 | codeSave = 1004; | |
356 | idRes = 25; | |
357 | coup2Z = settingsPtr->parm("HiggsH1:coup2Z"); | |
358 | } | |
359 | else if (higgsType == 2) { | |
360 | nameSave = "f fbar -> H0(H2) Z0"; | |
361 | codeSave = 1024; | |
362 | idRes = 35; | |
363 | coup2Z = settingsPtr->parm("HiggsH2:coup2Z"); | |
364 | } | |
365 | else if (higgsType == 3) { | |
366 | nameSave = "f fbar -> A0(A3) ZO"; | |
367 | codeSave = 1044; | |
368 | idRes = 36; | |
369 | coup2Z = settingsPtr->parm("HiggsA3:coup2Z"); | |
370 | } | |
371 | ||
372 | // Store Z0 mass and width for propagator. Common coupling factor. | |
373 | mZ = particleDataPtr->m0(23); | |
374 | widZ = particleDataPtr->mWidth(23); | |
375 | mZS = mZ*mZ; | |
376 | mwZS = pow2(mZ * widZ); | |
377 | thetaWRat = 1. / (16. * couplingsPtr->sin2thetaW() | |
378 | * couplingsPtr->cos2thetaW()); | |
379 | ||
380 | // Secondary open width fraction. | |
381 | openFracPair = particleDataPtr->resOpenFrac(idRes, 23); | |
382 | ||
383 | } | |
384 | ||
385 | //-------------------------------------------------------------------------- | |
386 | ||
387 | // Evaluate d(sigmaHat)/d(tHat), part independent of incoming flavour. | |
388 | ||
389 | void Sigma2ffbar2HZ::sigmaKin() { | |
390 | ||
391 | // Evaluate differential cross section. | |
392 | sigma0 = (M_PI / sH2) * 8. * pow2(alpEM * thetaWRat * coup2Z) | |
393 | * (tH * uH - s3 * s4 + 2. * sH * s4) / (pow2(sH - mZS) + mwZS); | |
394 | ||
395 | } | |
396 | ||
397 | //-------------------------------------------------------------------------- | |
398 | ||
399 | // Evaluate d(sigmaHat)/d(tHat), including incoming flavour dependence. | |
400 | ||
401 | double Sigma2ffbar2HZ::sigmaHat() { | |
402 | ||
403 | // Coupling a_f^2 + v_f^2 to s-channel Z0 and colour factor. | |
404 | int idAbs = abs(id1); | |
405 | double sigma = sigma0 * couplingsPtr->vf2af2(idAbs); | |
406 | if (idAbs < 9) sigma /= 3.; | |
407 | ||
408 | // Secondary width for H0 and Z0 or H1 and Z0 or H2 and Z0 or A3 and Z0. | |
409 | sigma *= openFracPair; | |
410 | ||
411 | // Answer. | |
412 | return sigma; | |
413 | ||
414 | } | |
415 | ||
416 | //-------------------------------------------------------------------------- | |
417 | ||
418 | // Select identity, colour and anticolour. | |
419 | ||
420 | void Sigma2ffbar2HZ::setIdColAcol() { | |
421 | ||
422 | // Flavours trivial. | |
423 | setId( id1, id2, idRes, 23); | |
424 | ||
425 | // Colour flow topologies. Swap when antiquarks. | |
426 | if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0); | |
427 | else setColAcol( 0, 0, 0, 0, 0, 0); | |
428 | if (id1 < 0) swapColAcol(); | |
429 | ||
430 | } | |
431 | ||
432 | //-------------------------------------------------------------------------- | |
433 | ||
434 | // Evaluate weight for decay angles. | |
435 | ||
436 | double Sigma2ffbar2HZ::weightDecay( Event& process, int iResBeg, | |
437 | int iResEnd) { | |
438 | ||
439 | // Identity of mother of decaying reseonance(s). | |
440 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
441 | ||
442 | // For Higgs decay hand over to standard routine. | |
443 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
444 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
445 | ||
446 | // For top decay hand over to standard routine. | |
447 | if (idMother == 6) | |
448 | return weightTopDecay( process, iResBeg, iResEnd); | |
449 | ||
450 | // If not decay of Z0 created along with Higgs then done. | |
451 | if (iResBeg != 5 || iResEnd != 6) return 1.; | |
452 | ||
453 | // Order so that fbar(1) f(2) -> H() f'(3) fbar'(4). | |
454 | int i1 = (process[3].id() < 0) ? 3 : 4; | |
455 | int i2 = 7 - i1; | |
456 | int i3 = process[6].daughter1(); | |
457 | int i4 = process[6].daughter2(); | |
458 | if (process[i3].id() < 0) swap( i3, i4); | |
459 | ||
460 | // Find left- and righthanded couplings of fermion pairs. | |
461 | int idAbs = process[i1].idAbs(); | |
462 | double liS = pow2( couplingsPtr->lf(idAbs) ); | |
463 | double riS = pow2( couplingsPtr->rf(idAbs) ); | |
464 | idAbs = process[i3].idAbs(); | |
465 | double lfS = pow2( couplingsPtr->lf(idAbs) ); | |
466 | double rfS = pow2( couplingsPtr->rf(idAbs) ); | |
467 | ||
468 | // Evaluate relevant four-products. | |
469 | double pp13 = process[i1].p() * process[i3].p(); | |
470 | double pp14 = process[i1].p() * process[i4].p(); | |
471 | double pp23 = process[i2].p() * process[i3].p(); | |
472 | double pp24 = process[i2].p() * process[i4].p(); | |
473 | ||
474 | // Weight and maximum. | |
475 | double wt = (liS * lfS + riS * rfS) * pp13 * pp24 | |
476 | + (liS * rfS + riS * lfS) * pp14 * pp23; | |
477 | double wtMax = (liS + riS) * (lfS + rfS) * (pp13 + pp14) * (pp23 + pp24); | |
478 | ||
479 | // Done. | |
480 | return wt / wtMax; | |
481 | ||
482 | } | |
483 | ||
484 | //========================================================================== | |
485 | ||
486 | // Sigma2ffbar2HW class. | |
487 | // Cross section for f fbar -> H0 W+-, H1 W+-, H2 W+- or A3 W+-. | |
488 | // (H0 SM Higgs, H1, H2 and A3 BSM Higgses). | |
489 | ||
490 | //-------------------------------------------------------------------------- | |
491 | ||
492 | // Initialize process. | |
493 | ||
494 | void Sigma2ffbar2HW::initProc() { | |
495 | ||
496 | // Properties specific to Higgs state. | |
497 | if (higgsType == 0) { | |
498 | nameSave = "f fbar -> H0 W+- (SM)"; | |
499 | codeSave = 905; | |
500 | idRes = 25; | |
501 | coup2W = 1.; | |
502 | } | |
503 | else if (higgsType == 1) { | |
504 | nameSave = "f fbar -> h0(H1) W+-"; | |
505 | codeSave = 1005; | |
506 | idRes = 25; | |
507 | coup2W = settingsPtr->parm("HiggsH1:coup2W"); | |
508 | } | |
509 | else if (higgsType == 2) { | |
510 | nameSave = "f fbar -> H0(H2) W+-"; | |
511 | codeSave = 1025; | |
512 | idRes = 35; | |
513 | coup2W = settingsPtr->parm("HiggsH2:coup2W"); | |
514 | } | |
515 | else if (higgsType == 3) { | |
516 | nameSave = "f fbar -> A0(A3) W+-"; | |
517 | codeSave = 1045; | |
518 | idRes = 36; | |
519 | coup2W = settingsPtr->parm("HiggsA3:coup2W"); | |
520 | } | |
521 | ||
522 | // Store W+- mass and width for propagator. Common coupling factor. | |
523 | mW = particleDataPtr->m0(24); | |
524 | widW = particleDataPtr->mWidth(24); | |
525 | mWS = mW*mW; | |
526 | mwWS = pow2(mW * widW); | |
527 | thetaWRat = 1. / (4. * couplingsPtr->sin2thetaW()); | |
528 | ||
529 | // Secondary open width fractions. | |
530 | openFracPairPos = particleDataPtr->resOpenFrac(idRes, 24); | |
531 | openFracPairNeg = particleDataPtr->resOpenFrac(idRes, -24); | |
532 | ||
533 | } | |
534 | ||
535 | //-------------------------------------------------------------------------- | |
536 | ||
537 | // Evaluate d(sigmaHat)/d(tHat), part independent of incoming flavour. | |
538 | ||
539 | void Sigma2ffbar2HW::sigmaKin() { | |
540 | ||
541 | // Evaluate differential cross section. | |
542 | sigma0 = (M_PI / sH2) * 2. * pow2(alpEM * thetaWRat * coup2W) | |
543 | * (tH * uH - s3 * s4 + 2. * sH * s4) / (pow2(sH - mWS) + mwWS); | |
544 | ||
545 | } | |
546 | ||
547 | //-------------------------------------------------------------------------- | |
548 | ||
549 | // Evaluate d(sigmaHat)/d(tHat), including incoming flavour dependence. | |
550 | ||
551 | double Sigma2ffbar2HW::sigmaHat() { | |
552 | ||
553 | // CKM and colour factors. | |
554 | double sigma = sigma0; | |
555 | if (abs(id1) < 9) sigma *= couplingsPtr->V2CKMid(abs(id1), abs(id2)) / 3.; | |
556 | ||
557 | // Secondary width for H0 and W+-. | |
558 | int idUp = (abs(id1)%2 == 0) ? id1 : id2; | |
559 | sigma *= (idUp > 0) ? openFracPairPos : openFracPairNeg; | |
560 | ||
561 | // Answer. | |
562 | return sigma; | |
563 | ||
564 | } | |
565 | ||
566 | //-------------------------------------------------------------------------- | |
567 | ||
568 | // Select identity, colour and anticolour. | |
569 | ||
570 | void Sigma2ffbar2HW::setIdColAcol() { | |
571 | ||
572 | // Sign of outgoing W. | |
573 | int sign = 1 - 2 * (abs(id1)%2); | |
574 | if (id1 < 0) sign = -sign; | |
575 | setId( id1, id2, idRes, 24 * sign); | |
576 | ||
577 | // Colour flow topologies. Swap when antiquarks. | |
578 | if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0); | |
579 | else setColAcol( 0, 0, 0, 0, 0, 0); | |
580 | if (id1 < 0) swapColAcol(); | |
581 | ||
582 | } | |
583 | ||
584 | //-------------------------------------------------------------------------- | |
585 | ||
586 | // Evaluate weight for decay angles. | |
587 | ||
588 | double Sigma2ffbar2HW::weightDecay( Event& process, int iResBeg, | |
589 | int iResEnd) { | |
590 | ||
591 | // Identity of mother of decaying reseonance(s). | |
592 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
593 | ||
594 | // For Higgs decay hand over to standard routine. | |
595 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
596 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
597 | ||
598 | // For top decay hand over to standard routine. | |
599 | if (idMother == 6) | |
600 | return weightTopDecay( process, iResBeg, iResEnd); | |
601 | ||
602 | // If not decay of W+- created along with Higgs then done. | |
603 | if (iResBeg != 5 || iResEnd != 6) return 1.; | |
604 | ||
605 | // Order so that fbar(1) f(2) -> H() f'(3) fbar'(4). | |
606 | int i1 = (process[3].id() < 0) ? 3 : 4; | |
607 | int i2 = 7 - i1; | |
608 | int i3 = process[6].daughter1(); | |
609 | int i4 = process[6].daughter2(); | |
610 | if (process[i3].id() < 0) swap( i3, i4); | |
611 | ||
612 | // Evaluate relevant four-products. | |
613 | double pp13 = process[i1].p() * process[i3].p(); | |
614 | double pp14 = process[i1].p() * process[i4].p(); | |
615 | double pp23 = process[i2].p() * process[i3].p(); | |
616 | double pp24 = process[i2].p() * process[i4].p(); | |
617 | ||
618 | // Weight and maximum. | |
619 | double wt = pp13 * pp24; | |
620 | double wtMax = (pp13 + pp14) * (pp23 + pp24); | |
621 | ||
622 | // Done. | |
623 | return wt / wtMax; | |
624 | ||
625 | } | |
626 | ||
627 | //========================================================================== | |
628 | ||
629 | // Sigma3ff2HfftZZ class. | |
630 | // Cross section for f f' -> H f f' (Z0 Z0 fusion of SM or BSM Higgs). | |
631 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
632 | ||
633 | //-------------------------------------------------------------------------- | |
634 | ||
635 | // Initialize process. | |
636 | ||
637 | void Sigma3ff2HfftZZ::initProc() { | |
638 | ||
639 | // Properties specific to Higgs state. | |
640 | if (higgsType == 0) { | |
641 | nameSave = "f f' -> H0 f f'(Z0 Z0 fusion) (SM)"; | |
642 | codeSave = 906; | |
643 | idRes = 25; | |
644 | coup2Z = 1.; | |
645 | } | |
646 | else if (higgsType == 1) { | |
647 | nameSave = "f f' -> h0(H1) f f' (Z0 Z0 fusion)"; | |
648 | codeSave = 1006; | |
649 | idRes = 25; | |
650 | coup2Z = settingsPtr->parm("HiggsH1:coup2Z"); | |
651 | } | |
652 | else if (higgsType == 2) { | |
653 | nameSave = "f f' -> H0(H2) f f' (Z0 Z0 fusion)"; | |
654 | codeSave = 1026; | |
655 | idRes = 35; | |
656 | coup2Z = settingsPtr->parm("HiggsH2:coup2Z"); | |
657 | } | |
658 | else if (higgsType == 3) { | |
659 | nameSave = "f f' -> A0(A3) f f' (Z0 Z0 fusion)"; | |
660 | codeSave = 1046; | |
661 | idRes = 36; | |
662 | coup2Z = settingsPtr->parm("HiggsA3:coup2Z"); | |
663 | } | |
664 | ||
665 | // Common fixed mass and coupling factor. | |
666 | mZS = pow2( particleDataPtr->m0(23) ); | |
667 | prefac = 0.25 * mZS * pow3( 4. * M_PI / (couplingsPtr->sin2thetaW() | |
668 | * couplingsPtr->cos2thetaW()) ); | |
669 | ||
670 | // Secondary open width fraction. | |
671 | openFrac = particleDataPtr->resOpenFrac(idRes); | |
672 | ||
673 | } | |
674 | ||
675 | //-------------------------------------------------------------------------- | |
676 | ||
677 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
678 | ||
679 | void Sigma3ff2HfftZZ::sigmaKin() { | |
680 | ||
681 | // Required four-vector products. | |
682 | double pp12 = 0.5 * sH; | |
683 | double pp14 = 0.5 * mH * p4cm.pNeg(); | |
684 | double pp15 = 0.5 * mH * p5cm.pNeg(); | |
685 | double pp24 = 0.5 * mH * p4cm.pPos(); | |
686 | double pp25 = 0.5 * mH * p5cm.pPos(); | |
687 | double pp45 = p4cm * p5cm; | |
688 | ||
689 | // Propagator factors and two possible numerators. | |
690 | double prop = pow2( (2. * pp14 + mZS) * (2. * pp25 + mZS) ); | |
691 | sigma1 = prefac * pp12 * pp45 / prop; | |
692 | sigma2 = prefac * pp15 * pp24 / prop; | |
693 | ||
694 | } | |
695 | ||
696 | //-------------------------------------------------------------------------- | |
697 | ||
698 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
699 | ||
700 | double Sigma3ff2HfftZZ::sigmaHat() { | |
701 | ||
702 | // Flavour-dependent coupling factors for two incoming flavours. | |
703 | int id1Abs = abs(id1); | |
704 | int id2Abs = abs(id2); | |
705 | double lf1S = pow2( couplingsPtr->lf(id1Abs) ); | |
706 | double rf1S = pow2( couplingsPtr->rf(id1Abs) ); | |
707 | double lf2S = pow2( couplingsPtr->lf(id2Abs) ); | |
708 | double rf2S = pow2( couplingsPtr->rf(id2Abs) ); | |
709 | double c1 = lf1S * lf2S + rf1S * rf2S; | |
710 | double c2 = lf1S * rf2S + rf1S * lf2S; | |
711 | ||
712 | // Combine couplings and kinematics factors. | |
713 | double sigma = pow3(alpEM) * (c1 * sigma1 + c2 * sigma2) * pow2(coup2Z); | |
714 | ||
715 | // Secondary width for H0, H1, H2 or A3. | |
716 | sigma *= openFrac; | |
717 | ||
718 | // Answer. | |
719 | return sigma; | |
720 | ||
721 | } | |
722 | ||
723 | //-------------------------------------------------------------------------- | |
724 | ||
725 | // Select identity, colour and anticolour. | |
726 | ||
727 | void Sigma3ff2HfftZZ::setIdColAcol() { | |
728 | ||
729 | // Trivial flavours: out = in. | |
730 | setId( id1, id2, idRes, id1, id2); | |
731 | ||
732 | // Colour flow topologies. Swap when antiquarks. | |
733 | if (abs(id1) < 9 && abs(id2) < 9 && id1*id2 > 0) | |
734 | setColAcol( 1, 0, 2, 0, 0, 0, 1, 0, 2, 0); | |
735 | else if (abs(id1) < 9 && abs(id2) < 9) | |
736 | setColAcol( 1, 0, 0, 2, 0, 0, 1, 0, 0, 2); | |
737 | else if (abs(id1) < 9) setColAcol( 1, 0, 0, 0, 0, 0, 1, 0, 0, 0); | |
738 | else if (abs(id2) < 9) setColAcol( 0, 0, 1, 0, 0, 0, 0, 0, 1, 0); | |
739 | else setColAcol( 0, 0, 0, 0, 0, 0, 0, 0, 0, 0); | |
740 | if ( (abs(id1) < 9 && id1 < 0) || (abs(id1) > 10 && id2 < 0) ) | |
741 | swapColAcol(); | |
742 | ||
743 | } | |
744 | ||
745 | //-------------------------------------------------------------------------- | |
746 | ||
747 | // Evaluate weight for decay angles. | |
748 | ||
749 | double Sigma3ff2HfftZZ::weightDecay( Event& process, int iResBeg, | |
750 | int iResEnd) { | |
751 | ||
752 | // Identity of mother of decaying reseonance(s). | |
753 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
754 | ||
755 | // For Higgs decay hand over to standard routine. | |
756 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
757 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
758 | ||
759 | // For top decay hand over to standard routine. | |
760 | if (idMother == 6) | |
761 | return weightTopDecay( process, iResBeg, iResEnd); | |
762 | ||
763 | // Else done. | |
764 | return 1.; | |
765 | ||
766 | } | |
767 | ||
768 | //========================================================================== | |
769 | ||
770 | // Sigma3ff2HfftWW class. | |
771 | // Cross section for f_1 f_2 -> H0 f_3 f_4 (W+ W- fusion of SM or BSM Higgs). | |
772 | ||
773 | //-------------------------------------------------------------------------- | |
774 | ||
775 | // Initialize process. | |
776 | ||
777 | void Sigma3ff2HfftWW::initProc() { | |
778 | ||
779 | // Properties specific to Higgs state. | |
780 | if (higgsType == 0) { | |
781 | nameSave = "f_1 f_2 -> H0 f_3 f_4 (W+ W- fusion) (SM)"; | |
782 | codeSave = 907; | |
783 | idRes = 25; | |
784 | coup2W = 1.; | |
785 | } | |
786 | else if (higgsType == 1) { | |
787 | nameSave = "f_1 f_2 -> h0(H1) f_3 f_4 (W+ W- fusion)"; | |
788 | codeSave = 1007; | |
789 | idRes = 25; | |
790 | coup2W = settingsPtr->parm("HiggsH1:coup2W"); | |
791 | } | |
792 | else if (higgsType == 2) { | |
793 | nameSave = "f_1 f_2 -> H0(H2) f_3 f_4 (W+ W- fusion)"; | |
794 | codeSave = 1027; | |
795 | idRes = 35; | |
796 | coup2W = settingsPtr->parm("HiggsH2:coup2W"); | |
797 | } | |
798 | else if (higgsType == 3) { | |
799 | nameSave = "f_1 f_2 -> A0(A3) f_3 f_4 (W+ W- fusion)"; | |
800 | codeSave = 1047; | |
801 | idRes = 36; | |
802 | coup2W = settingsPtr->parm("HiggsA3:coup2W"); | |
803 | } | |
804 | ||
805 | // Common fixed mass and coupling factor. | |
806 | mWS = pow2( particleDataPtr->m0(24) ); | |
807 | prefac = mWS * pow3( 4. * M_PI / couplingsPtr->sin2thetaW() ); | |
808 | ||
809 | // Secondary open width fraction. | |
810 | openFrac = particleDataPtr->resOpenFrac(idRes); | |
811 | ||
812 | } | |
813 | ||
814 | //-------------------------------------------------------------------------- | |
815 | ||
816 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
817 | ||
818 | void Sigma3ff2HfftWW::sigmaKin() { | |
819 | ||
820 | // Required four-vector products. | |
821 | double pp12 = 0.5 * sH; | |
822 | double pp14 = 0.5 * mH * p4cm.pNeg(); | |
823 | double pp25 = 0.5 * mH * p5cm.pPos(); | |
824 | double pp45 = p4cm * p5cm; | |
825 | ||
826 | // Cross section: kinematics part. Combine with couplings. | |
827 | double prop = pow2( (2. * pp14 + mWS) * (2. * pp25 + mWS) ); | |
828 | sigma0 = prefac * pp12 * pp45 * pow2(coup2W) / prop; | |
829 | ||
830 | } | |
831 | ||
832 | //-------------------------------------------------------------------------- | |
833 | ||
834 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
835 | ||
836 | double Sigma3ff2HfftWW::sigmaHat() { | |
837 | ||
838 | // Some flavour combinations not possible. | |
839 | int id1Abs = abs(id1); | |
840 | int id2Abs = abs(id2); | |
841 | if ( (id1Abs%2 == id2Abs%2 && id1 * id2 > 0) | |
842 | || (id1Abs%2 != id2Abs%2 && id1 * id2 < 0) ) return 0.; | |
843 | ||
844 | // Basic cross section. CKM factors for final states. | |
845 | double sigma = sigma0 * pow3(alpEM) * couplingsPtr->V2CKMsum(id1Abs) | |
846 | * couplingsPtr->V2CKMsum(id2Abs); | |
847 | ||
848 | // Secondary width for H0, H1, H2 or A3. | |
849 | sigma *= openFrac; | |
850 | ||
851 | // Spin-state extra factor 2 per incoming neutrino. | |
852 | if (id1Abs == 12 || id1Abs == 14 || id1Abs == 16) sigma *= 2.; | |
853 | if (id2Abs == 12 || id2Abs == 14 || id2Abs == 16) sigma *= 2.; | |
854 | ||
855 | // Answer. | |
856 | return sigma; | |
857 | ||
858 | } | |
859 | ||
860 | //-------------------------------------------------------------------------- | |
861 | ||
862 | // Select identity, colour and anticolour. | |
863 | ||
864 | void Sigma3ff2HfftWW::setIdColAcol() { | |
865 | ||
866 | // Pick out-flavours by relative CKM weights. | |
867 | id4 = couplingsPtr->V2CKMpick(id1); | |
868 | id5 = couplingsPtr->V2CKMpick(id2); | |
869 | setId( id1, id2, idRes, id4, id5); | |
870 | ||
871 | // Colour flow topologies. Swap when antiquarks. | |
872 | if (abs(id1) < 9 && abs(id2) < 9 && id1*id2 > 0) | |
873 | setColAcol( 1, 0, 2, 0, 0, 0, 1, 0, 2, 0); | |
874 | else if (abs(id1) < 9 && abs(id2) < 9) | |
875 | setColAcol( 1, 0, 0, 2, 0, 0, 1, 0, 0, 2); | |
876 | else if (abs(id1) < 9) setColAcol( 1, 0, 0, 0, 0, 0, 1, 0, 0, 0); | |
877 | else if (abs(id2) < 9) setColAcol( 0, 0, 1, 0, 0, 0, 0, 0, 1, 0); | |
878 | else setColAcol( 0, 0, 0, 0, 0, 0, 0, 0, 0, 0); | |
879 | if ( (abs(id1) < 9 && id1 < 0) || (abs(id1) > 10 && id2 < 0) ) | |
880 | swapColAcol(); | |
881 | ||
882 | } | |
883 | ||
884 | //-------------------------------------------------------------------------- | |
885 | ||
886 | // Evaluate weight for decay angles. | |
887 | ||
888 | double Sigma3ff2HfftWW::weightDecay( Event& process, int iResBeg, | |
889 | int iResEnd) { | |
890 | ||
891 | // Identity of mother of decaying reseonance(s). | |
892 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
893 | ||
894 | // For Higgs decay hand over to standard routine. | |
895 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
896 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
897 | ||
898 | // For top decay hand over to standard routine. | |
899 | if (idMother == 6) | |
900 | return weightTopDecay( process, iResBeg, iResEnd); | |
901 | ||
902 | // Else done. | |
903 | return 1.; | |
904 | ||
905 | } | |
906 | ||
907 | //========================================================================== | |
908 | ||
909 | // Sigma3gg2HQQbar class. | |
910 | // Cross section for g g -> H0 Q Qbar (Q Qbar fusion of SM or BSM Higgs). | |
911 | ||
912 | //-------------------------------------------------------------------------- | |
913 | ||
914 | // Initialize process. | |
915 | ||
916 | void Sigma3gg2HQQbar::initProc() { | |
917 | ||
918 | // Properties specific to Higgs state for the "g g -> H ttbar" process. | |
919 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
920 | if (higgsType == 0 && idNew == 6) { | |
921 | nameSave = "g g -> H t tbar (SM)"; | |
922 | codeSave = 908; | |
923 | idRes = 25; | |
924 | coup2Q = 1.; | |
925 | } | |
926 | else if (higgsType == 1 && idNew == 6) { | |
927 | nameSave = "g g -> h0(H1) t tbar"; | |
928 | codeSave = 1008; | |
929 | idRes = 25; | |
930 | coup2Q = settingsPtr->parm("HiggsH1:coup2u"); | |
931 | } | |
932 | else if (higgsType == 2 && idNew == 6) { | |
933 | nameSave = "g g -> H0(H2) t tbar"; | |
934 | codeSave = 1028; | |
935 | idRes = 35; | |
936 | coup2Q = settingsPtr->parm("HiggsH2:coup2u"); | |
937 | } | |
938 | else if (higgsType == 3 && idNew == 6) { | |
939 | nameSave = "g g -> A0(A3) t tbar"; | |
940 | codeSave = 1048; | |
941 | idRes = 36; | |
942 | coup2Q = settingsPtr->parm("HiggsA3:coup2u"); | |
943 | } | |
944 | ||
945 | // Properties specific to Higgs state for the "g g -> H b bbar" process. | |
946 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
947 | if (higgsType == 0 && idNew == 5) { | |
948 | nameSave = "g g -> H b bbar (SM)"; | |
949 | codeSave = 912; | |
950 | idRes = 25; | |
951 | coup2Q = 1.; | |
952 | } | |
953 | else if (higgsType == 1 && idNew == 5) { | |
954 | nameSave = "g g -> h0(H1) b bbar"; | |
955 | codeSave = 1012; | |
956 | idRes = 25; | |
957 | coup2Q = settingsPtr->parm("HiggsH1:coup2d"); | |
958 | } | |
959 | else if (higgsType == 2 && idNew == 5) { | |
960 | nameSave = "g g -> H0(H2) b bbar"; | |
961 | codeSave = 1032; | |
962 | idRes = 35; | |
963 | coup2Q = settingsPtr->parm("HiggsH2:coup2d"); | |
964 | } | |
965 | else if (higgsType == 3 && idNew == 5) { | |
966 | nameSave = "g g -> A0(A3) b bbar"; | |
967 | codeSave = 1052; | |
968 | idRes = 36; | |
969 | coup2Q = settingsPtr->parm("HiggsA3:coup2d"); | |
970 | } | |
971 | ||
972 | // Common mass and coupling factors. | |
973 | double mWS = pow2(particleDataPtr->m0(24)); | |
974 | prefac = (4. * M_PI / couplingsPtr->sin2thetaW()) * pow2(4. * M_PI) | |
975 | * 0.25 / mWS; | |
976 | ||
977 | // Secondary open width fraction. | |
978 | openFracTriplet = particleDataPtr->resOpenFrac(idRes, idNew, -idNew); | |
979 | ||
980 | } | |
981 | ||
982 | //-------------------------------------------------------------------------- | |
983 | ||
984 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
985 | ||
986 | void Sigma3gg2HQQbar::sigmaKin() { | |
987 | ||
988 | // Running mass of heavy quark. | |
989 | double mQ2run = pow2( particleDataPtr->mRun(idNew, mH) ); | |
990 | ||
991 | // Linear combination of p_Q and p_Qbar to ensure common mass. | |
992 | double mQ2 = m4 * m5; | |
993 | double epsi = 0.; | |
994 | if (m4 != m5) { | |
995 | double s45 = (p4cm + p5cm).m2Calc(); | |
996 | mQ2 = 0.5 * (s4 + s5) - 0.25 * pow2(s4 - s5) / s45; | |
997 | epsi = 0.5 * (s5 - s4) / s45; | |
998 | } | |
999 | ||
1000 | // Set up kinematics: g(4) g(5) -> H(3) Q(1) Qbar(2) in outgoing sense. | |
1001 | Vec4 pTemp[6]; | |
1002 | pTemp[4] = Vec4( 0., 0., -0.5* mH, -0.5* mH); | |
1003 | pTemp[5] = Vec4( 0., 0., 0.5* mH, -0.5* mH); | |
1004 | pTemp[1] = p4cm + epsi * (p4cm + p5cm); | |
1005 | pTemp[2] = p5cm - epsi * (p4cm + p5cm); | |
1006 | pTemp[3] = p3cm; | |
1007 | ||
1008 | // Four-product combinations. | |
1009 | double z1 = pTemp[1] * pTemp[2]; | |
1010 | double z2 = pTemp[1] * pTemp[3]; | |
1011 | double z3 = pTemp[1] * pTemp[4]; | |
1012 | double z4 = pTemp[1] * pTemp[5]; | |
1013 | double z5 = pTemp[2] * pTemp[3]; | |
1014 | double z6 = pTemp[2] * pTemp[4]; | |
1015 | double z7 = pTemp[2] * pTemp[5]; | |
1016 | double z8 = pTemp[3] * pTemp[4]; | |
1017 | double z9 = pTemp[3] * pTemp[5]; | |
1018 | double z10 = pTemp[4] * pTemp[5]; | |
1019 | ||
1020 | // Powers required as shorthand in matriz elements. | |
1021 | double mQ4 = mQ2 * mQ2; | |
1022 | double mQ6 = mQ2 * mQ4; | |
1023 | double z1S = z1 * z1; | |
1024 | double z2S = z2 * z2; | |
1025 | double z3S = z3 * z3; | |
1026 | double z4S = z4 * z4; | |
1027 | double z5S = z5 * z5; | |
1028 | double z6S = z6 * z6; | |
1029 | double z7S = z7 * z7; | |
1030 | double z8S = z8 * z8; | |
1031 | double z9S = z9 * z9; | |
1032 | double z10S = z10 * z10; | |
1033 | ||
1034 | // Evaluate matriz elements for g + g -> Q + Qbar + H. | |
1035 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1036 | double fm[9][9]; | |
1037 | fm[1][1] = 64*mQ6+16*mQ4*s3+32*mQ4*(z1+2*z2+z4+z9+2* | |
1038 | z7+z5)+8*mQ2*s3*(-z1-z4+2*z7)+16*mQ2*(z2*z9+4*z2* | |
1039 | z7+z2*z5-2*z4*z7-2*z9*z7)+8*s3*z4*z7-16*z2*z9*z7; | |
1040 | fm[1][2] = 16*mQ6+8*mQ4*(-2*z1+z2-2*z3-2*z4-4*z10+z9-z8+2 | |
1041 | *z7-4*z6+z5)+8*mQ2*(-2*z1*z2-2*z2*z4-2*z2*z10+z2*z7-2* | |
1042 | z2*z6-2*z3*z7+2*z4*z7+4*z10*z7-z9*z7-z8*z7)+16*z2*z7*(z4+ | |
1043 | z10); | |
1044 | fm[1][3] = 16*mQ6-4*mQ4*s3+8*mQ4*(-2*z1+2*z2-2*z3-4* | |
1045 | z4-8*z10+z9+z8-2*z7-4*z6+2*z5)-(4*mQ2*s3)*(z1+z4+z10 | |
1046 | +z6)+8*mQ2*(-2*z1*z2-2*z1*z10+z1*z9+z1*z8-2*z1*z5+z2S | |
1047 | -4*z2*z4-5*z2*z10+z2*z8-z2*z7-3*z2*z6+z2*z5+z3*z9+2*z3*z7 | |
1048 | -z3*z5+z4*z8+2*z4*z6-3*z4*z5-5*z10*z5+z9*z8+z9*z6+z9*z5+ | |
1049 | z8*z7-4*z6*z5+z5S)-(16*z2*z5)*(z1+z4+z10+z6); | |
1050 | fm[1][4] = 16*mQ6+4*mQ4*s3+16*mQ4*(-z1+z2-z3-z4+z10- | |
1051 | z9-z8+2*z7+2*z6-z5)+4*mQ2*s3*(z1+z3+z4+z10+2*z7+2*z6 | |
1052 | )+8*mQ2*(4*z1*z10+4*z1*z7+4*z1*z6+2*z2*z10-z2*z9-z2*z8+ | |
1053 | 4*z2*z7+4*z2*z6-z2*z5+4*z10*z5+4*z7*z5+4*z6*z5)-(8*s3* | |
1054 | z1)*(z10+z7+z6)+16*z2*z5*(z10+z7+z6); | |
1055 | fm[1][5] = 8*mQ4*(-2*z1-2*z4+z10-z9)+4*mQ2*(4*z1S-2*z1* | |
1056 | z2+8*z1*z3+6*z1*z10-2*z1*z9+4*z1*z8+4*z1*z7+4*z1*z6+2*z1* | |
1057 | z5+z2*z10+4*z3*z4-z3*z9+2*z3*z7+3*z4*z8-2*z4*z6+2*z4*z5-4 | |
1058 | *z10*z7+3*z10*z5-3*z9*z6+3*z8*z7-4*z7S+4*z7*z5)+8*(z1S | |
1059 | *z9-z1S*z8-z1*z2*z7+z1*z2*z6+z1*z3*z9+z1*z3*z5-z1*z4* | |
1060 | z8-z1*z4*z5+z1*z10*z9+z1*z9*z7+z1*z9*z6-z1*z8*z7-z2*z3*z7 | |
1061 | +z2*z4*z6-z2*z10*z7-z2*z7S+z3*z7*z5-z4*z10*z5-z4*z7*z5- | |
1062 | z4*z6*z5); | |
1063 | fm[1][6] = 16*mQ4*(-4*z1-z4+z9-z7)+4*mQ2*s3*(-2*z1-z4- | |
1064 | z7)+16*mQ2*(-2*z1S-3*z1*z2-2*z1*z4-3*z1*z9-2*z1*z7-3* | |
1065 | z1*z5-2*z2*z4-2*z7*z5)-8*s3*z4*z7+8*(-z1*z2*z9-2*z1*z2 | |
1066 | *z5-z1*z9S-z1*z9*z5+z2S*z7-z2*z4*z5+z2*z9*z7-z2*z7*z5 | |
1067 | +z4*z9*z5+z4*z5S); | |
1068 | fm[1][7] = 8*mQ4*(2*z3+z4+3*z10+z9+2*z8+3*z7+6*z6)+2*mQ2* | |
1069 | s3*(-2*z3-z4+3*z10+3*z7+6*z6)+4*mQ2*(4*z1*z10+4*z1* | |
1070 | z7+8*z1*z6+6*z2*z10+z2*z9+2*z2*z8+6*z2*z7+12*z2*z6-8*z3* | |
1071 | z7+4*z4*z7+4*z4*z6+4*z10*z5+4*z9*z7+4*z9*z6-8*z8*z7+4*z7* | |
1072 | z5+8*z6*z5)+4*s3*(-z1*z10-z1*z7-2*z1*z6+2*z3*z7-z4*z7- | |
1073 | z4*z6)+8*z2*(z10*z5+z9*z7+z9*z6-2*z8*z7+z7*z5+2*z6*z5); | |
1074 | fm[1][8] = 8*mQ4*(2*z3+z4+3*z10+2*z9+z8+3*z7+6*z6)+2*mQ2* | |
1075 | s3*(-2*z3-z4+2*z10+z7+2*z6)+4*mQ2*(4*z1*z10-2*z1*z9+ | |
1076 | 2*z1*z8+4*z1*z7+8*z1*z6+5*z2*z10+2*z2*z9+z2*z8+4*z2*z7+8* | |
1077 | z2*z6-z3*z9-8*z3*z7+2*z3*z5+2*z4*z9-z4*z8+4*z4*z7+4*z4*z6 | |
1078 | +4*z4*z5+5*z10*z5+z9S-z9*z8+2*z9*z7+5*z9*z6+z9*z5-7*z8* | |
1079 | z7+2*z8*z5+2*z7*z5+10*z6*z5)+2*s3*(-z1*z10+z3*z7-2*z4* | |
1080 | z7+z4*z6)+4*(-z1*z9S+z1*z9*z8-2*z1*z9*z5-z1*z8*z5+2*z2* | |
1081 | z10*z5+z2*z9*z7+z2*z9*z6-2*z2*z8*z7+3*z2*z6*z5+z3*z9*z5+ | |
1082 | z3*z5S+z4*z9*z5-2*z4*z8*z5+2*z4*z5S); | |
1083 | fm[2][2] = 16*mQ6+16*mQ4*(-z1+z3-z4-z10+z7-z6)+16*mQ2*( | |
1084 | z3*z10+z3*z7+z3*z6+z4*z7+z10*z7)-16*z3*z10*z7; | |
1085 | fm[2][3] = 16*mQ6+8*mQ4*(-2*z1+z2+2*z3-4*z4-4*z10-z9+z8-2 | |
1086 | *z7-2*z6+z5)+8*mQ2*(-2*z1*z5+4*z3*z10-z3*z9-z3*z8-2*z3* | |
1087 | z7+2*z3*z6+z3*z5-2*z4*z5-2*z10*z5-2*z6*z5)+16*z3*z5*(z10+ | |
1088 | z6); | |
1089 | fm[2][4] = 8*mQ4*(-2*z1-2*z3+z10-z8)+4*mQ2*(4*z1S-2*z1* | |
1090 | z2+8*z1*z4+6*z1*z10+4*z1*z9-2*z1*z8+4*z1*z7+4*z1*z6+2*z1* | |
1091 | z5+z2*z10+4*z3*z4+3*z3*z9-2*z3*z7+2*z3*z5-z4*z8+2*z4*z6-4 | |
1092 | *z10*z6+3*z10*z5+3*z9*z6-3*z8*z7-4*z6S+4*z6*z5)+8*(-z1S | |
1093 | *z9+z1S*z8+z1*z2*z7-z1*z2*z6-z1*z3*z9-z1*z3*z5+z1*z4 | |
1094 | *z8+z1*z4*z5+z1*z10*z8-z1*z9*z6+z1*z8*z7+z1*z8*z6+z2*z3* | |
1095 | z7-z2*z4*z6-z2*z10*z6-z2*z6S-z3*z10*z5-z3*z7*z5-z3*z6* | |
1096 | z5+z4*z6*z5); | |
1097 | fm[2][5] = 16*mQ4*z10+8*mQ2*(2*z1S+2*z1*z3+2*z1*z4+2*z1 | |
1098 | *z10+2*z1*z7+2*z1*z6+z3*z7+z4*z6)+8*(-2*pow3(z1)-2*z1S*z3- | |
1099 | 2*z1S*z4-2*z1S*z10-2*z1S*z7-2*z1S*z6-2*z1*z3*z4- | |
1100 | z1*z3*z10-2*z1*z3*z6-z1*z4*z10-2*z1*z4*z7-z1*z10S-z1* | |
1101 | z10*z7-z1*z10*z6-2*z1*z7*z6+z3S*z7-z3*z4*z7-z3*z4*z6+z3 | |
1102 | *z10*z7+z3*z7S-z3*z7*z6+z4S*z6+z4*z10*z6-z4*z7*z6+z4* | |
1103 | z6S); | |
1104 | fm[2][6] = 8*mQ4*(-2*z1+z10-z9-2*z7)+4*mQ2*(4*z1S+2*z1* | |
1105 | z2+4*z1*z3+4*z1*z4+6*z1*z10-2*z1*z9+4*z1*z8+8*z1*z6-2*z1* | |
1106 | z5+4*z2*z4+3*z2*z10+2*z2*z7-3*z3*z9-2*z3*z7-4*z4S-4*z4* | |
1107 | z10+3*z4*z8+2*z4*z6+z10*z5-z9*z6+3*z8*z7+4*z7*z6)+8*(z1S | |
1108 | *z9-z1S*z8-z1*z2*z7+z1*z2*z6+z1*z3*z9+z1*z3*z5+z1*z4* | |
1109 | z9-z1*z4*z8-z1*z4*z5+z1*z10*z9+z1*z9*z6-z1*z8*z7-z2*z3*z7 | |
1110 | -z2*z4*z7+z2*z4*z6-z2*z10*z7+z3*z7*z5-z4S*z5-z4*z10*z5- | |
1111 | z4*z6*z5); | |
1112 | fm[2][7] = 8*mQ4*(z3+2*z4+3*z10+z7+2*z6)+4*mQ2*(-4*z1*z3- | |
1113 | 2*z1*z4-2*z1*z10+z1*z9-z1*z8-4*z1*z7-2*z1*z6+z2*z3+2*z2* | |
1114 | z4+3*z2*z10+z2*z7+2*z2*z6-6*z3*z4-6*z3*z10-2*z3*z9-2*z3* | |
1115 | z7-4*z3*z6-z3*z5-6*z4S-6*z4*z10-3*z4*z9-z4*z8-4*z4*z7-2 | |
1116 | *z4*z6-2*z4*z5-3*z10*z9-3*z10*z8-6*z10*z7-6*z10*z6+z10*z5 | |
1117 | +z9*z7-2*z8*z7-2*z8*z6-6*z7*z6+z7*z5-6*z6S+2*z6*z5)+4*( | |
1118 | -z1S*z9+z1S*z8-2*z1*z2*z10-3*z1*z2*z7-3*z1*z2*z6+z1* | |
1119 | z3*z9-z1*z3*z5+z1*z4*z9+z1*z4*z8+z1*z4*z5+z1*z10*z9+z1* | |
1120 | z10*z8-z1*z9*z6+z1*z8*z6+z2*z3*z7-3*z2*z4*z7-z2*z4*z6-3* | |
1121 | z2*z10*z7-3*z2*z10*z6-3*z2*z7*z6-3*z2*z6S-2*z3*z4*z5-z3 | |
1122 | *z10*z5-z3*z6*z5-z4S*z5-z4*z10*z5+z4*z6*z5); | |
1123 | fm[2][8] = 8*mQ4*(z3+2*z4+3*z10+z7+2*z6)+4*mQ2*(-4*z1*z3- | |
1124 | 2*z1*z4-2*z1*z10-z1*z9+z1*z8-4*z1*z7-2*z1*z6+z2*z3+2*z2* | |
1125 | z4+z2*z10-z2*z7-2*z2*z6-6*z3*z4-6*z3*z10-2*z3*z9+z3*z8-2* | |
1126 | z3*z7-4*z3*z6+z3*z5-6*z4S-6*z4*z10-2*z4*z9-4*z4*z7-2*z4 | |
1127 | *z6+2*z4*z5-3*z10*z9-3*z10*z8-6*z10*z7-6*z10*z6+3*z10*z5- | |
1128 | z9*z6-2*z8*z7-3*z8*z6-6*z7*z6+z7*z5-6*z6S+2*z6*z5)+4*( | |
1129 | z1S*z9-z1S*z8-z1*z2*z7+z1*z2*z6-3*z1*z3*z5+z1*z4*z9- | |
1130 | z1*z4*z8-3*z1*z4*z5+z1*z10*z9+z1*z10*z8-2*z1*z10*z5+z1*z9 | |
1131 | *z6+z1*z8*z7+z1*z8*z6-z2*z4*z7+z2*z4*z6-z2*z10*z7-z2*z10* | |
1132 | z6-2*z2*z7*z6-z2*z6S-3*z3*z4*z5-3*z3*z10*z5+z3*z7*z5-3* | |
1133 | z3*z6*z5-3*z4S*z5-3*z4*z10*z5-z4*z6*z5); | |
1134 | fm[3][3] = 64*mQ6+16*mQ4*s3+32*mQ4*(z1+z2+2*z3+z8+z6 | |
1135 | +2*z5)+8*mQ2*s3*(-z1+2*z3-z6)+16*mQ2*(z2*z5-2*z3* | |
1136 | z8-2*z3*z6+4*z3*z5+z8*z5)+8*s3*z3*z6-16*z3*z8*z5; | |
1137 | fm[3][4] = 16*mQ4*(-4*z1-z3+z8-z6)+4*mQ2*s3*(-2*z1-z3- | |
1138 | z6)+16*mQ2*(-2*z1S-3*z1*z2-2*z1*z3-3*z1*z8-2*z1*z6-3* | |
1139 | z1*z5-2*z2*z3-2*z6*z5)-8*s3*z3*z6+8*(-z1*z2*z8-2*z1*z2 | |
1140 | *z5-z1*z8S-z1*z8*z5+z2S*z6-z2*z3*z5+z2*z8*z6-z2*z6*z5 | |
1141 | +z3*z8*z5+z3*z5S); | |
1142 | fm[3][5] = 8*mQ4*(-2*z1+z10-z8-2*z6)+4*mQ2*(4*z1S+2*z1* | |
1143 | z2+4*z1*z3+4*z1*z4+6*z1*z10+4*z1*z9-2*z1*z8+8*z1*z7-2*z1* | |
1144 | z5+4*z2*z3+3*z2*z10+2*z2*z6-4*z3S-4*z3*z10+3*z3*z9+2*z3 | |
1145 | *z7-3*z4*z8-2*z4*z6+z10*z5+3*z9*z6-z8*z7+4*z7*z6)+8*(-z1S | |
1146 | *z9+z1S*z8+z1*z2*z7-z1*z2*z6-z1*z3*z9+z1*z3*z8-z1*z3 | |
1147 | *z5+z1*z4*z8+z1*z4*z5+z1*z10*z8-z1*z9*z6+z1*z8*z7+z2*z3* | |
1148 | z7-z2*z3*z6-z2*z4*z6-z2*z10*z6-z3S*z5-z3*z10*z5-z3*z7* | |
1149 | z5+z4*z6*z5); | |
1150 | fm[3][6] = 16*mQ6+4*mQ4*s3+16*mQ4*(-z1-z2+2*z3+2*z4+ | |
1151 | z10-z9-z8-z7-z6+z5)+4*mQ2*s3*(z1+2*z3+2*z4+z10+z7+z6 | |
1152 | )+8*mQ2*(4*z1*z3+4*z1*z4+4*z1*z10+4*z2*z3+4*z2*z4+4*z2* | |
1153 | z10-z2*z5+4*z3*z5+4*z4*z5+2*z10*z5-z9*z5-z8*z5)-(8*s3* | |
1154 | z1)*(z3+z4+z10)+16*z2*z5*(z3+z4+z10); | |
1155 | fm[3][7] = 8*mQ4*(3*z3+6*z4+3*z10+z9+2*z8+2*z7+z6)+2*mQ2* | |
1156 | s3*(z3+2*z4+2*z10-2*z7-z6)+4*mQ2*(4*z1*z3+8*z1*z4+4* | |
1157 | z1*z10+2*z1*z9-2*z1*z8+2*z2*z3+10*z2*z4+5*z2*z10+2*z2*z9+ | |
1158 | z2*z8+2*z2*z7+4*z2*z6-7*z3*z9+2*z3*z8-8*z3*z7+4*z3*z6+4* | |
1159 | z3*z5+5*z4*z8+4*z4*z6+8*z4*z5+5*z10*z5-z9*z8-z9*z6+z9*z5+ | |
1160 | z8S-z8*z7+2*z8*z6+2*z8*z5)+2*s3*(-z1*z10+z3*z7-2*z3* | |
1161 | z6+z4*z6)+4*(-z1*z2*z9-2*z1*z2*z8+z1*z9*z8-z1*z8S+z2S | |
1162 | *z7+2*z2S*z6+3*z2*z4*z5+2*z2*z10*z5-2*z2*z9*z6+z2*z8*z7 | |
1163 | +z2*z8*z6-2*z3*z9*z5+z3*z8*z5+z4*z8*z5); | |
1164 | fm[3][8] = 8*mQ4*(3*z3+6*z4+3*z10+2*z9+z8+2*z7+z6)+2*mQ2* | |
1165 | s3*(3*z3+6*z4+3*z10-2*z7-z6)+4*mQ2*(4*z1*z3+8*z1*z4+ | |
1166 | 4*z1*z10+4*z2*z3+8*z2*z4+4*z2*z10-8*z3*z9+4*z3*z8-8*z3*z7 | |
1167 | +4*z3*z6+6*z3*z5+4*z4*z8+4*z4*z6+12*z4*z5+6*z10*z5+2*z9* | |
1168 | z5+z8*z5)+4*s3*(-z1*z3-2*z1*z4-z1*z10+2*z3*z7-z3*z6-z4 | |
1169 | *z6)+8*z5*(z2*z3+2*z2*z4+z2*z10-2*z3*z9+z3*z8+z4*z8); | |
1170 | fm[4][4] = 64*mQ6+16*mQ4*s3+32*mQ4*(z1+2*z2+z3+z8+2* | |
1171 | z6+z5)+8*mQ2*s3*(-z1-z3+2*z6)+16*mQ2*(z2*z8+4*z2* | |
1172 | z6+z2*z5-2*z3*z6-2*z8*z6)+8*s3*z3*z6-16*z2*z8*z6; | |
1173 | fm[4][5] = 16*mQ6+8*mQ4*(-2*z1+z2-2*z3-2*z4-4*z10-z9+z8-4 | |
1174 | *z7+2*z6+z5)+8*mQ2*(-2*z1*z2-2*z2*z3-2*z2*z10-2*z2*z7+ | |
1175 | z2*z6+2*z3*z6-2*z4*z6+4*z10*z6-z9*z6-z8*z6)+16*z2*z6*(z3+ | |
1176 | z10); | |
1177 | fm[4][6] = 16*mQ6-4*mQ4*s3+8*mQ4*(-2*z1+2*z2-4*z3-2* | |
1178 | z4-8*z10+z9+z8-4*z7-2*z6+2*z5)-(4*mQ2*s3)*(z1+z3+z10 | |
1179 | +z7)+8*mQ2*(-2*z1*z2-2*z1*z10+z1*z9+z1*z8-2*z1*z5+z2S | |
1180 | -4*z2*z3-5*z2*z10+z2*z9-3*z2*z7-z2*z6+z2*z5+z3*z9+2*z3*z7 | |
1181 | -3*z3*z5+z4*z8+2*z4*z6-z4*z5-5*z10*z5+z9*z8+z9*z6+z8*z7+ | |
1182 | z8*z5-4*z7*z5+z5S)-(16*z2*z5)*(z1+z3+z10+z7); | |
1183 | fm[4][7] = 8*mQ4*(-z3-2*z4-3*z10-2*z9-z8-6*z7-3*z6)+2*mQ2 | |
1184 | *s3*(z3+2*z4-3*z10-6*z7-3*z6)+4*mQ2*(-4*z1*z10-8*z1* | |
1185 | z7-4*z1*z6-6*z2*z10-2*z2*z9-z2*z8-12*z2*z7-6*z2*z6-4*z3* | |
1186 | z7-4*z3*z6+8*z4*z6-4*z10*z5+8*z9*z6-4*z8*z7-4*z8*z6-8*z7* | |
1187 | z5-4*z6*z5)+4*s3*(z1*z10+2*z1*z7+z1*z6+z3*z7+z3*z6-2* | |
1188 | z4*z6)+8*z2*(-z10*z5+2*z9*z6-z8*z7-z8*z6-2*z7*z5-z6*z5); | |
1189 | fm[4][8] = 8*mQ4*(-z3-2*z4-3*z10-z9-2*z8-6*z7-3*z6)+2*mQ2 | |
1190 | *s3*(z3+2*z4-2*z10-2*z7-z6)+4*mQ2*(-4*z1*z10-2*z1*z9 | |
1191 | +2*z1*z8-8*z1*z7-4*z1*z6-5*z2*z10-z2*z9-2*z2*z8-8*z2*z7-4 | |
1192 | *z2*z6+z3*z9-2*z3*z8-4*z3*z7-4*z3*z6-4*z3*z5+z4*z8+8*z4* | |
1193 | z6-2*z4*z5-5*z10*z5+z9*z8+7*z9*z6-2*z9*z5-z8S-5*z8*z7-2 | |
1194 | *z8*z6-z8*z5-10*z7*z5-2*z6*z5)+2*s3*(z1*z10-z3*z7+2*z3 | |
1195 | *z6-z4*z6)+4*(-z1*z9*z8+z1*z9*z5+z1*z8S+2*z1*z8*z5-2*z2 | |
1196 | *z10*z5+2*z2*z9*z6-z2*z8*z7-z2*z8*z6-3*z2*z7*z5+2*z3*z9* | |
1197 | z5-z3*z8*z5-2*z3*z5S-z4*z8*z5-z4*z5S); | |
1198 | fm[5][5] = 16*mQ6+16*mQ4*(-z1-z3+z4-z10-z7+z6)+16*mQ2*( | |
1199 | z3*z6+z4*z10+z4*z7+z4*z6+z10*z6)-16*z4*z10*z6; | |
1200 | fm[5][6] = 16*mQ6+8*mQ4*(-2*z1+z2-4*z3+2*z4-4*z10+z9-z8-2 | |
1201 | *z7-2*z6+z5)+8*mQ2*(-2*z1*z5-2*z3*z5+4*z4*z10-z4*z9-z4* | |
1202 | z8+2*z4*z7-2*z4*z6+z4*z5-2*z10*z5-2*z7*z5)+16*z4*z5*(z10+ | |
1203 | z7); | |
1204 | fm[5][7] = 8*mQ4*(-2*z3-z4-3*z10-2*z7-z6)+4*mQ2*(2*z1*z3+ | |
1205 | 4*z1*z4+2*z1*z10+z1*z9-z1*z8+2*z1*z7+4*z1*z6-2*z2*z3-z2* | |
1206 | z4-3*z2*z10-2*z2*z7-z2*z6+6*z3S+6*z3*z4+6*z3*z10+z3*z9+ | |
1207 | 3*z3*z8+2*z3*z7+4*z3*z6+2*z3*z5+6*z4*z10+2*z4*z8+4*z4*z7+ | |
1208 | 2*z4*z6+z4*z5+3*z10*z9+3*z10*z8+6*z10*z7+6*z10*z6-z10*z5+ | |
1209 | 2*z9*z7+2*z9*z6-z8*z6+6*z7S+6*z7*z6-2*z7*z5-z6*z5)+4*(- | |
1210 | z1S*z9+z1S*z8+2*z1*z2*z10+3*z1*z2*z7+3*z1*z2*z6-z1*z3 | |
1211 | *z9-z1*z3*z8-z1*z3*z5-z1*z4*z8+z1*z4*z5-z1*z10*z9-z1*z10* | |
1212 | z8-z1*z9*z7+z1*z8*z7+z2*z3*z7+3*z2*z3*z6-z2*z4*z6+3*z2* | |
1213 | z10*z7+3*z2*z10*z6+3*z2*z7S+3*z2*z7*z6+z3S*z5+2*z3*z4 | |
1214 | *z5+z3*z10*z5-z3*z7*z5+z4*z10*z5+z4*z7*z5); | |
1215 | fm[5][8] = 8*mQ4*(-2*z3-z4-3*z10-2*z7-z6)+4*mQ2*(2*z1*z3+ | |
1216 | 4*z1*z4+2*z1*z10-z1*z9+z1*z8+2*z1*z7+4*z1*z6-2*z2*z3-z2* | |
1217 | z4-z2*z10+2*z2*z7+z2*z6+6*z3S+6*z3*z4+6*z3*z10+2*z3*z8+ | |
1218 | 2*z3*z7+4*z3*z6-2*z3*z5+6*z4*z10-z4*z9+2*z4*z8+4*z4*z7+2* | |
1219 | z4*z6-z4*z5+3*z10*z9+3*z10*z8+6*z10*z7+6*z10*z6-3*z10*z5+ | |
1220 | 3*z9*z7+2*z9*z6+z8*z7+6*z7S+6*z7*z6-2*z7*z5-z6*z5)+4*( | |
1221 | z1S*z9-z1S*z8-z1*z2*z7+z1*z2*z6+z1*z3*z9-z1*z3*z8+3* | |
1222 | z1*z3*z5+3*z1*z4*z5-z1*z10*z9-z1*z10*z8+2*z1*z10*z5-z1*z9 | |
1223 | *z7-z1*z9*z6-z1*z8*z7-z2*z3*z7+z2*z3*z6+z2*z10*z7+z2*z10* | |
1224 | z6+z2*z7S+2*z2*z7*z6+3*z3S*z5+3*z3*z4*z5+3*z3*z10*z5+ | |
1225 | z3*z7*z5+3*z4*z10*z5+3*z4*z7*z5-z4*z6*z5); | |
1226 | fm[6][6] = 64*mQ6+16*mQ4*s3+32*mQ4*(z1+z2+2*z4+z9+z7 | |
1227 | +2*z5)+8*mQ2*s3*(-z1+2*z4-z7)+16*mQ2*(z2*z5-2*z4* | |
1228 | z9-2*z4*z7+4*z4*z5+z9*z5)+8*s3*z4*z7-16*z4*z9*z5; | |
1229 | fm[6][7] = 8*mQ4*(-6*z3-3*z4-3*z10-2*z9-z8-z7-2*z6)+2*mQ2 | |
1230 | *s3*(-2*z3-z4-2*z10+z7+2*z6)+4*mQ2*(-8*z1*z3-4*z1*z4 | |
1231 | -4*z1*z10+2*z1*z9-2*z1*z8-10*z2*z3-2*z2*z4-5*z2*z10-z2*z9 | |
1232 | -2*z2*z8-4*z2*z7-2*z2*z6-5*z3*z9-4*z3*z7-8*z3*z5-2*z4*z9+ | |
1233 | 7*z4*z8-4*z4*z7+8*z4*z6-4*z4*z5-5*z10*z5-z9S+z9*z8-2*z9 | |
1234 | *z7+z9*z6-2*z9*z5+z8*z7-z8*z5)+2*s3*(z1*z10-z3*z7+2*z4 | |
1235 | *z7-z4*z6)+4*(2*z1*z2*z9+z1*z2*z8+z1*z9S-z1*z9*z8-2* | |
1236 | z2S*z7-z2S*z6-3*z2*z3*z5-2*z2*z10*z5-z2*z9*z7-z2*z9*z6+ | |
1237 | 2*z2*z8*z7-z3*z9*z5-z4*z9*z5+2*z4*z8*z5); | |
1238 | fm[6][8] = 8*mQ4*(-6*z3-3*z4-3*z10-z9-2*z8-z7-2*z6)+2*mQ2 | |
1239 | *s3*(-6*z3-3*z4-3*z10+z7+2*z6)+4*mQ2*(-8*z1*z3-4*z1* | |
1240 | z4-4*z1*z10-8*z2*z3-4*z2*z4-4*z2*z10-4*z3*z9-4*z3*z7-12* | |
1241 | z3*z5-4*z4*z9+8*z4*z8-4*z4*z7+8*z4*z6-6*z4*z5-6*z10*z5-z9 | |
1242 | *z5-2*z8*z5)+4*s3*(2*z1*z3+z1*z4+z1*z10+z3*z7+z4*z7-2* | |
1243 | z4*z6)+8*z5*(-2*z2*z3-z2*z4-z2*z10-z3*z9-z4*z9+2*z4*z8); | |
1244 | fm[7][7] = 72*mQ4*z10+18*mQ2*s3*z10+8*mQ2*(z1*z10+9* | |
1245 | z2*z10+7*z3*z7+2*z3*z6+2*z4*z7+7*z4*z6+z10*z5+2*z9*z7+7* | |
1246 | z9*z6+7*z8*z7+2*z8*z6)+2*s3*(-z1*z10-7*z3*z7-2*z3*z6-2 | |
1247 | *z4*z7-7*z4*z6)+4*z2*(z10*z5+2*z9*z7+7*z9*z6+7*z8*z7+2*z8 | |
1248 | *z6); | |
1249 | fm[7][8] = 72*mQ4*z10+2*mQ2*s3*z10+4*mQ2*(2*z1*z10+ | |
1250 | 10*z2*z10+7*z3*z9+2*z3*z8+14*z3*z7+4*z3*z6+2*z4*z9+7*z4* | |
1251 | z8+4*z4*z7+14*z4*z6+10*z10*z5+z9S+7*z9*z8+2*z9*z7+7*z9* | |
1252 | z6+z8S+7*z8*z7+2*z8*z6)+2*s3*(7*z1*z10-7*z3*z7-2*z3* | |
1253 | z6-2*z4*z7-7*z4*z6)+2*(-2*z1*z9S-14*z1*z9*z8-2*z1*z8S | |
1254 | +2*z2*z10*z5+2*z2*z9*z7+7*z2*z9*z6+7*z2*z8*z7+2*z2*z8*z6+ | |
1255 | 7*z3*z9*z5+2*z3*z8*z5+2*z4*z9*z5+7*z4*z8*z5); | |
1256 | fm[8][8] = 72*mQ4*z10+18*mQ2*s3*z10+8*mQ2*(z1*z10+z2 | |
1257 | *z10+7*z3*z9+2*z3*z8+7*z3*z7+2*z3*z6+2*z4*z9+7*z4*z8+2*z4 | |
1258 | *z7+7*z4*z6+9*z10*z5)+2*s3*(-z1*z10-7*z3*z7-2*z3*z6-2* | |
1259 | z4*z7-7*z4*z6)+4*z5*(z2*z10+7*z3*z9+2*z3*z8+2*z4*z9+7*z4* | |
1260 | z8); | |
1261 | double fm99 = -4*mQ4*z10-mQ2*s3*z10+4*mQ2*(-z1*z10-z2*z10+ | |
1262 | z3*z7+z4*z6-z10*z5+z9*z6+z8*z7)+s3*(z1*z10-z3*z7-z4*z6 | |
1263 | )+2*z2*(-z10*z5+z9*z6+z8*z7); | |
1264 | double fm910 = -4*mQ4*z10-mQ2*s3*z10+2*mQ2*(-2*z1*z10-2*z2* | |
1265 | z10+2*z3*z9+2*z3*z7+2*z4*z6-2*z10*z5+z9*z8+2*z8*z7)+s3 | |
1266 | *(z1*z10-z3*z7-z4*z6)+2*(-z1*z9*z8-z2*z10*z5+z2*z8*z7+z3* | |
1267 | z9*z5); | |
1268 | double fmxx = -4*mQ4*z10-mQ2*s3*z10+2*mQ2*(-2*z1*z10-2*z2* | |
1269 | z10+2*z4*z8+2*z4*z6+2*z3*z7-2*z10*z5+z9*z8+2*z9*z6)+s3 | |
1270 | *(z1*z10-z3*z7-z4*z6)+2*(-z1*z9*z8-z2*z10*z5+z2*z9*z6+z4* | |
1271 | z8*z5); | |
1272 | fm910 = 0.5*(fmxx+fm910); | |
1273 | double fm1010 = -4*mQ4*z10-mQ2*s3*z10+4*mQ2*(-z1*z10-z2*z10+ | |
1274 | z3*z7+z4*z6-z10*z5+z9*z3+z8*z4)+s3*(z1*z10-z3*z7-z4*z6 | |
1275 | )+2*z5*(-z10*z2+z9*z3+z8*z4); | |
1276 | fm[7][7] -= 2. * fm99; | |
1277 | fm[7][8] -= 2. * fm910; | |
1278 | fm[8][8] -= 2. * fm1010; | |
1279 | ||
1280 | // Propagators. | |
1281 | double ss1 = (pTemp[1] + pTemp[3]).m2Calc() - mQ2; | |
1282 | double ss2 = (pTemp[1] + pTemp[4]).m2Calc() - mQ2; | |
1283 | double ss3 = (pTemp[1] + pTemp[5]).m2Calc() - mQ2; | |
1284 | double ss4 = (pTemp[2] + pTemp[3]).m2Calc() - mQ2; | |
1285 | double ss5 = (pTemp[2] + pTemp[4]).m2Calc() - mQ2; | |
1286 | double ss6 = (pTemp[2] + pTemp[5]).m2Calc() - mQ2; | |
1287 | double ss7 = sH; | |
1288 | ||
1289 | // Propagator combinations. | |
1290 | double dz[9]; | |
1291 | dz[1] = ss1 * ss6; | |
1292 | dz[2] = ss2 * ss6; | |
1293 | dz[3] = ss2 * ss4; | |
1294 | dz[4] = ss1 * ss5; | |
1295 | dz[5] = ss3 * ss5; | |
1296 | dz[6] = ss3 * ss4; | |
1297 | dz[7] = ss7 * ss1; | |
1298 | dz[8] = ss7 * ss4; | |
1299 | ||
1300 | // Colour factors. | |
1301 | double clr[9][9]; | |
1302 | for (int i = 1; i < 4; ++i) | |
1303 | for (int j = 1; j < 4; ++j) { | |
1304 | clr[i][j] = 16. / 3.; | |
1305 | clr[i][j+3] = -2. / 3.; | |
1306 | clr[i+3][j] = -2. / 3.; | |
1307 | clr[i+3][j+3] = 16. / 3.; | |
1308 | } | |
1309 | for (int i = 1; i < 4; ++i) | |
1310 | for (int j = 1; j < 3; ++j) { | |
1311 | clr[i][j+6] = -6.; | |
1312 | clr[i+3][j+6] = 6.; | |
1313 | clr[j+6][i] = -6.; | |
1314 | clr[j+6][i+3] = 6.; | |
1315 | } | |
1316 | for (int i = 1; i < 3; ++i) | |
1317 | for (int j = 1; j < 3; ++j) | |
1318 | clr[i+6][j+6] = 12.; | |
1319 | ||
1320 | // Produce final result: matrix elements * colours * propagators. | |
1321 | double wtSum = 0.; | |
1322 | for (int i = 1; i < 9; ++i) | |
1323 | for (int j = i; j < 9; ++j) { | |
1324 | double fac = (j == i) ? 4. : 8.; | |
1325 | wtSum += fm[i][j] * fac * clr[i][j] / (dz[i] * dz[j]); | |
1326 | } | |
1327 | wtSum *= -1./256.; | |
1328 | ||
1329 | // Combine factors. | |
1330 | sigma = prefac * alpEM * pow2(alpS) * mQ2run * wtSum *pow2(coup2Q); | |
1331 | ||
1332 | // Secondary width for H, Q and Qbar (latter for top only). | |
1333 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1334 | sigma *= openFracTriplet; | |
1335 | ||
1336 | } | |
1337 | ||
1338 | //-------------------------------------------------------------------------- | |
1339 | ||
1340 | // Select identity, colour and anticolour. | |
1341 | ||
1342 | void Sigma3gg2HQQbar::setIdColAcol() { | |
1343 | ||
1344 | // Pick out-flavours by relative CKM weights. | |
1345 | setId( id1, id2, idRes, idNew, -idNew); | |
1346 | ||
1347 | // Colour flow topologies. | |
1348 | if (rndmPtr->flat() < 0.5) setColAcol( 1, 2, 2, 3, 0, 0, 1, 0, 0, 3); | |
1349 | else setColAcol( 1, 2, 3, 1, 0, 0, 3, 0, 0, 2); | |
1350 | ||
1351 | } | |
1352 | ||
1353 | //-------------------------------------------------------------------------- | |
1354 | ||
1355 | // Evaluate weight for decay angles. | |
1356 | ||
1357 | double Sigma3gg2HQQbar::weightDecay( Event& process, int iResBeg, | |
1358 | int iResEnd) { | |
1359 | ||
1360 | // Identity of mother of decaying reseonance(s). | |
1361 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
1362 | ||
1363 | // For Higgs decay hand over to standard routine. | |
1364 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
1365 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
1366 | ||
1367 | // For top decay hand over to standard routine. | |
1368 | if (idMother == 6) | |
1369 | return weightTopDecay( process, iResBeg, iResEnd); | |
1370 | ||
1371 | // Else done. | |
1372 | return 1.; | |
1373 | ||
1374 | } | |
1375 | ||
1376 | //========================================================================== | |
1377 | ||
1378 | // Sigma3qqbar2HQQbar class. | |
1379 | // Cross section for q qbar -> H0 Q Qbar (Q Qbar fusion of SM Higgs). | |
1380 | // REDUCE output and part of the rest courtesy Z. Kunszt, | |
1381 | // see Z. Kunszt, Nucl. Phys. B247 (1984) 339. | |
1382 | ||
1383 | //-------------------------------------------------------------------------- | |
1384 | ||
1385 | // Initialize process. | |
1386 | ||
1387 | void Sigma3qqbar2HQQbar::initProc() { | |
1388 | ||
1389 | // Properties specific to Higgs state for the "q qbar -> H ttbar" process. | |
1390 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1391 | ||
1392 | if (higgsType == 0 && idNew == 6) { | |
1393 | nameSave = "q qbar -> H t tbar (SM)"; | |
1394 | codeSave = 909; | |
1395 | idRes = 25; | |
1396 | coup2Q = 1.; | |
1397 | } | |
1398 | else if (higgsType == 1 && idNew == 6) { | |
1399 | nameSave = "q qbar -> h0(H1) t tbar"; | |
1400 | codeSave = 1009; | |
1401 | idRes = 25; | |
1402 | coup2Q = settingsPtr->parm("HiggsH1:coup2u"); | |
1403 | } | |
1404 | else if (higgsType == 2 && idNew == 6) { | |
1405 | nameSave = "q qbar -> H0(H2) t tbar"; | |
1406 | codeSave = 1029; | |
1407 | idRes = 35; | |
1408 | coup2Q = settingsPtr->parm("HiggsH2:coup2u"); | |
1409 | } | |
1410 | else if (higgsType == 3 && idNew == 6) { | |
1411 | nameSave = "q qbar -> A0(A3) t tbar"; | |
1412 | codeSave = 1049; | |
1413 | idRes = 36; | |
1414 | coup2Q = settingsPtr->parm("HiggsA3:coup2u"); | |
1415 | } | |
1416 | ||
1417 | // Properties specific to Higgs state for the "q qbar -> H b bbar" process. | |
1418 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1419 | if (higgsType == 0 && idNew == 5) { | |
1420 | nameSave = "q qbar -> H b bbar (SM)"; | |
1421 | codeSave = 913; | |
1422 | idRes = 25; | |
1423 | coup2Q = 1.; | |
1424 | } | |
1425 | else if (higgsType == 1 && idNew == 5) { | |
1426 | nameSave = "q qbar -> h0(H1) b bbar"; | |
1427 | codeSave = 1013; | |
1428 | idRes = 25; | |
1429 | coup2Q = settingsPtr->parm("HiggsH1:coup2d"); | |
1430 | } | |
1431 | else if (higgsType == 2 && idNew == 5) { | |
1432 | nameSave = "q qbar -> H0(H2) b bbar"; | |
1433 | codeSave = 1033; | |
1434 | idRes = 35; | |
1435 | coup2Q = settingsPtr->parm("HiggsH2:coup2d"); | |
1436 | } | |
1437 | else if (higgsType == 3 && idNew == 5) { | |
1438 | nameSave = "q qbar -> A0(A3) b bbar"; | |
1439 | codeSave = 1053; | |
1440 | idRes = 36; | |
1441 | coup2Q = settingsPtr->parm("HiggsA3:coup2d"); | |
1442 | } | |
1443 | ||
1444 | // Common mass and coupling factors. | |
1445 | double mWS = pow2(particleDataPtr->m0(24)); | |
1446 | prefac = (4. * M_PI / couplingsPtr->sin2thetaW()) * pow2(4. * M_PI) | |
1447 | * 0.25 / mWS; | |
1448 | ||
1449 | // Secondary open width fraction. | |
1450 | openFracTriplet = particleDataPtr->resOpenFrac(idRes, idNew, -idNew); | |
1451 | ||
1452 | } | |
1453 | ||
1454 | //-------------------------------------------------------------------------- | |
1455 | ||
1456 | // Evaluate sigma(sHat), part independent of incoming flavour. | |
1457 | ||
1458 | void Sigma3qqbar2HQQbar::sigmaKin() { | |
1459 | ||
1460 | // Running mass of heavy quark. | |
1461 | double mQ2run = pow2( particleDataPtr->mRun(idNew, mH) ); | |
1462 | ||
1463 | // Linear combination of p_Q and p_Qbar to ensure common mass. | |
1464 | double mQ2 = m4 * m5; | |
1465 | double epsi = 0.; | |
1466 | if (m4 != m5) { | |
1467 | double s45 = (p4cm + p5cm).m2Calc(); | |
1468 | mQ2 = 0.5 * (s4 + s5) - 0.25 * pow2(s4 - s5) / s45; | |
1469 | epsi = 0.5 * (s5 - s4) / s45; | |
1470 | } | |
1471 | ||
1472 | // Set up kinematics: q(4) qbar(5) -> H(3) Q(1) Qbar(2) in outgoing sense. | |
1473 | Vec4 pTemp[6]; | |
1474 | pTemp[4] = Vec4( 0., 0., -0.5* mH, -0.5* mH); | |
1475 | pTemp[5] = Vec4( 0., 0., 0.5* mH, -0.5* mH); | |
1476 | pTemp[1] = p4cm + epsi * (p4cm + p5cm); | |
1477 | pTemp[2] = p5cm - epsi * (p4cm + p5cm); | |
1478 | pTemp[3] = p3cm; | |
1479 | ||
1480 | // Four-product combinations. | |
1481 | double z1 = pTemp[1] * pTemp[2]; | |
1482 | double z2 = pTemp[1] * pTemp[3]; | |
1483 | double z3 = pTemp[1] * pTemp[4]; | |
1484 | double z4 = pTemp[1] * pTemp[5]; | |
1485 | double z5 = pTemp[2] * pTemp[3]; | |
1486 | double z6 = pTemp[2] * pTemp[4]; | |
1487 | double z7 = pTemp[2] * pTemp[5]; | |
1488 | double z8 = pTemp[3] * pTemp[4]; | |
1489 | double z9 = pTemp[3] * pTemp[5]; | |
1490 | double z10 = pTemp[4] * pTemp[5]; | |
1491 | ||
1492 | // Powers required as shorthand in matriz elements. | |
1493 | double mQ4 = mQ2 * mQ2; | |
1494 | ||
1495 | // Evaluate matrix elements for q + qbar -> Q + Qbar + H. | |
1496 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1497 | double a11 = -8.*mQ4*z10-2.*mQ2*s3*z10-(8.*mQ2)*(z2*z10+z3 | |
1498 | *z7+z4*z6+z9*z6+z8*z7)+2.*s3*(z3*z7+z4*z6)-(4.*z2)*(z9 | |
1499 | *z6+z8*z7); | |
1500 | double a12 = -8.*mQ4*z10+4.*mQ2*(-z2*z10-z3*z9-2.*z3*z7-z4*z8- | |
1501 | 2.*z4*z6-z10*z5-z9*z8-z9*z6-z8*z7)+2.*s3*(-z1*z10+z3*z7 | |
1502 | +z4*z6)+2.*(2.*z1*z9*z8-z2*z9*z6-z2*z8*z7-z3*z9*z5-z4*z8* | |
1503 | z5); | |
1504 | double a22 = -8.*mQ4*z10-2.*mQ2*s3*z10-(8.*mQ2)*(z3*z9+z3* | |
1505 | z7+z4*z8+z4*z6+z10*z5)+2.*s3*(z3*z7+z4*z6)-(4.*z5)*(z3 | |
1506 | *z9+z4*z8); | |
1507 | ||
1508 | // Propagators and propagator combinations. | |
1509 | double ss1 = (pTemp[1] + pTemp[3]).m2Calc() - mQ2; | |
1510 | double ss4 = (pTemp[2] + pTemp[3]).m2Calc() - mQ2; | |
1511 | double ss7 = sH; | |
1512 | double dz7 = ss7 * ss1; | |
1513 | double dz8 = ss7 * ss4; | |
1514 | ||
1515 | // Produce final result: matrix elements * propagators. | |
1516 | a11 /= (dz7 * dz7); | |
1517 | a12 /= (dz7 * dz8); | |
1518 | a22 /= (dz8 * dz8); | |
1519 | double wtSum = -(a11 + a22 + 2.*a12) * (8./9.); | |
1520 | ||
1521 | // Combine factors. | |
1522 | sigma = prefac * alpEM * pow2(alpS) * mQ2run * wtSum * pow2(coup2Q); | |
1523 | ||
1524 | // Secondary width for H, Q and Qbar (latter for top only). | |
1525 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1526 | sigma *= openFracTriplet; | |
1527 | ||
1528 | } | |
1529 | ||
1530 | //-------------------------------------------------------------------------- | |
1531 | ||
1532 | // Select identity, colour and anticolour. | |
1533 | ||
1534 | void Sigma3qqbar2HQQbar::setIdColAcol() { | |
1535 | ||
1536 | // Pick out-flavours by relative CKM weights. | |
1537 | setId( id1, id2, idRes, idNew, -idNew); | |
1538 | ||
1539 | // Colour flow topologies. | |
1540 | if (id1 > 0) setColAcol( 1, 0, 0, 2, 0, 0, 1, 0, 0, 2); | |
1541 | else setColAcol( 0, 1, 2, 0, 0, 0, 2, 0, 0, 1); | |
1542 | ||
1543 | } | |
1544 | ||
1545 | //-------------------------------------------------------------------------- | |
1546 | ||
1547 | // Evaluate weight for decay angles. | |
1548 | ||
1549 | double Sigma3qqbar2HQQbar::weightDecay( Event& process, int iResBeg, | |
1550 | int iResEnd) { | |
1551 | ||
1552 | // Identity of mother of decaying reseonance(s). | |
1553 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
1554 | ||
1555 | // For Higgs decay hand over to standard routine. | |
1556 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
1557 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
1558 | ||
1559 | // For top decay hand over to standard routine. | |
1560 | if (idMother == 6) | |
1561 | return weightTopDecay( process, iResBeg, iResEnd); | |
1562 | ||
1563 | // Else done. | |
1564 | return 1.; | |
1565 | ||
1566 | } | |
1567 | ||
1568 | //========================================================================== | |
1569 | ||
1570 | // Sigma2qg2Hq class. | |
1571 | // Cross section for q g -> H q. | |
1572 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1573 | ||
1574 | //-------------------------------------------------------------------------- | |
1575 | ||
1576 | // Initialize process. | |
1577 | ||
1578 | void Sigma2qg2Hq::initProc() { | |
1579 | ||
1580 | // Properties specific to Higgs state for the "c g -> H c" process. | |
1581 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1582 | if (higgsType == 0 && idNew == 4) { | |
1583 | nameSave = "c g -> H c (SM)"; | |
1584 | codeSave = 911; | |
1585 | idRes = 25; | |
1586 | } | |
1587 | else if (higgsType == 1 && idNew == 4) { | |
1588 | nameSave = "c g -> h0(H1) c"; | |
1589 | codeSave = 1011; | |
1590 | idRes = 25; | |
1591 | } | |
1592 | else if (higgsType == 2 && idNew == 4) { | |
1593 | nameSave = "c g -> H0(H2) c"; | |
1594 | codeSave = 1031; | |
1595 | idRes = 35; | |
1596 | } | |
1597 | else if (higgsType == 3 && idNew == 4) { | |
1598 | nameSave = "c g -> A0(A3) c"; | |
1599 | codeSave = 1051; | |
1600 | idRes = 36; | |
1601 | } | |
1602 | ||
1603 | // Properties specific to Higgs state for the "b g -> H b" process. | |
1604 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1605 | if (higgsType == 0 && idNew == 5) { | |
1606 | nameSave = "b g -> H b (SM)"; | |
1607 | codeSave = 911; | |
1608 | idRes = 25; | |
1609 | } | |
1610 | else if (higgsType == 1 && idNew == 5) { | |
1611 | nameSave = "b g -> h0(H1) b"; | |
1612 | codeSave = 1011; | |
1613 | idRes = 25; | |
1614 | } | |
1615 | else if (higgsType == 2 && idNew == 5) { | |
1616 | nameSave = "b g -> H0(H2) b"; | |
1617 | codeSave = 1031; | |
1618 | idRes = 35; | |
1619 | } | |
1620 | else if (higgsType == 3 && idNew == 5) { | |
1621 | nameSave = "b g -> A0(A3) b"; | |
1622 | codeSave = 1051; | |
1623 | idRes = 36; | |
1624 | } | |
1625 | ||
1626 | // Standard parameters. | |
1627 | m2W = pow2( particleDataPtr->m0(24) ); | |
1628 | thetaWRat = 1. / (24. * couplingsPtr->sin2thetaW()); | |
1629 | ||
1630 | // Secondary open width fraction. | |
1631 | openFrac = particleDataPtr->resOpenFrac(idRes); | |
1632 | ||
1633 | ||
1634 | } | |
1635 | ||
1636 | //-------------------------------------------------------------------------- | |
1637 | ||
1638 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
1639 | ||
1640 | void Sigma2qg2Hq::sigmaKin() { | |
1641 | ||
1642 | // Running mass provides coupling. | |
1643 | double m2Run = pow2( particleDataPtr->mRun(idNew, mH) ); | |
1644 | ||
1645 | // Cross section, including couplings and kinematics. | |
1646 | sigma = (M_PI / sH2) * alpS * alpEM * thetaWRat * (m2Run/m2W) | |
1647 | * (sH / (s4 - uH) + 2. * s4 * (s3 - uH) / pow2(s4 - uH) | |
1648 | + (s4 - uH) / sH - 2. * s4 / (s4 - uH) | |
1649 | + 2. * (s3 - uH) * (s3 - s4 - sH) / ((s4 - uH) * sH) ); | |
1650 | ||
1651 | // Include secondary width for H0, H1, H2 or A3. Done. | |
1652 | sigma *= openFrac; | |
1653 | ||
1654 | } | |
1655 | ||
1656 | //-------------------------------------------------------------------------- | |
1657 | ||
1658 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
1659 | ||
1660 | double Sigma2qg2Hq::sigmaHat() { | |
1661 | ||
1662 | // Check that specified flavour present. | |
1663 | if (abs(id1) != idNew && abs(id2) != idNew) return 0.; | |
1664 | ||
1665 | // Answer. | |
1666 | return sigma; | |
1667 | ||
1668 | } | |
1669 | ||
1670 | ||
1671 | //-------------------------------------------------------------------------- | |
1672 | ||
1673 | // Select identity, colour and anticolour. | |
1674 | ||
1675 | void Sigma2qg2Hq::setIdColAcol() { | |
1676 | ||
1677 | // Flavour set up for q g -> H0 q. | |
1678 | int idq = (id2 == 21) ? id1 : id2; | |
1679 | setId( id1, id2, idRes, idq); | |
1680 | ||
1681 | // tH defined between f and f': must swap tHat <-> uHat if q g in. | |
1682 | swapTU = (id2 == 21); | |
1683 | ||
1684 | // Colour flow topologies. Swap when antiquarks. | |
1685 | if (id2 == 21) setColAcol( 1, 0, 2, 1, 0, 0, 2, 0); | |
1686 | else setColAcol( 2, 1, 1, 0, 0, 0, 2, 0); | |
1687 | if (idq < 0) swapColAcol(); | |
1688 | ||
1689 | } | |
1690 | ||
1691 | //-------------------------------------------------------------------------- | |
1692 | ||
1693 | // Evaluate weight for decay angles. | |
1694 | ||
1695 | double Sigma2qg2Hq::weightDecay( Event& process, int iResBeg, | |
1696 | int iResEnd) { | |
1697 | ||
1698 | // Identity of mother of decaying reseonance(s). | |
1699 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
1700 | ||
1701 | // For Higgs decay hand over to standard routine. | |
1702 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
1703 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
1704 | ||
1705 | // For top decay hand over to standard routine. | |
1706 | if (idMother == 6) | |
1707 | return weightTopDecay( process, iResBeg, iResEnd); | |
1708 | ||
1709 | // Else done. | |
1710 | return 1.; | |
1711 | ||
1712 | } | |
1713 | ||
1714 | //========================================================================== | |
1715 | ||
1716 | // Sigma2gg2Hglt class. | |
1717 | // Cross section for g g -> H g (H SM Higgs or BSM Higgs) via top loop. | |
1718 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1719 | ||
1720 | //-------------------------------------------------------------------------- | |
1721 | ||
1722 | // Initialize process. | |
1723 | ||
1724 | void Sigma2gg2Hglt::initProc() { | |
1725 | ||
1726 | // Properties specific to Higgs state. | |
1727 | if (higgsType == 0) { | |
1728 | nameSave = "g g -> H g (SM; top loop)"; | |
1729 | codeSave = 914; | |
1730 | idRes = 25; | |
1731 | } | |
1732 | else if (higgsType == 1) { | |
1733 | nameSave = "g g -> h0(H1) g (BSM; top loop)"; | |
1734 | codeSave = 1014; | |
1735 | idRes = 25; | |
1736 | } | |
1737 | else if (higgsType == 2) { | |
1738 | nameSave = "g g -> H0(H2) g (BSM; top loop)"; | |
1739 | codeSave = 1034; | |
1740 | idRes = 35; | |
1741 | } | |
1742 | else if (higgsType == 3) { | |
1743 | nameSave = "g g -> A0(A3) g (BSM; top loop)"; | |
1744 | codeSave = 1054; | |
1745 | idRes = 36; | |
1746 | } | |
1747 | ||
1748 | // Normalization factor by g g -> H partial width. | |
1749 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1750 | double mHiggs = particleDataPtr->m0(idRes); | |
1751 | widHgg = particleDataPtr->resWidthChan(idRes, mHiggs, 21, 21); | |
1752 | ||
1753 | // Secondary open width fraction. | |
1754 | openFrac = particleDataPtr->resOpenFrac(idRes); | |
1755 | ||
1756 | } | |
1757 | ||
1758 | //-------------------------------------------------------------------------- | |
1759 | ||
1760 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
1761 | ||
1762 | void Sigma2gg2Hglt::sigmaKin() { | |
1763 | ||
1764 | // Evaluate cross section. Secondary width for H0, H1, H2 or A3. | |
1765 | sigma = (M_PI / sH2) * (3. / 16.) * alpS * (widHgg / m3) | |
1766 | * (sH2 * sH2 + tH2 * tH2 + uH2 * uH2 + pow4(s3)) | |
1767 | / (sH * tH * uH * s3); | |
1768 | sigma *= openFrac; | |
1769 | ||
1770 | } | |
1771 | ||
1772 | //-------------------------------------------------------------------------- | |
1773 | ||
1774 | // Select identity, colour and anticolour. | |
1775 | ||
1776 | void Sigma2gg2Hglt::setIdColAcol() { | |
1777 | ||
1778 | // Flavour set up for g g -> H g trivial. | |
1779 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1780 | setId( 21, 21, idRes, 21); | |
1781 | ||
1782 | // Colour flow topologies: random choice between two mirrors. | |
1783 | if (rndmPtr->flat() < 0.5) setColAcol( 1, 2, 2, 3, 0, 0, 1, 3); | |
1784 | else setColAcol( 1, 2, 3, 1, 0, 0, 3, 2); | |
1785 | ||
1786 | } | |
1787 | ||
1788 | //-------------------------------------------------------------------------- | |
1789 | ||
1790 | // Evaluate weight for decay angles. | |
1791 | ||
1792 | double Sigma2gg2Hglt::weightDecay( Event& process, int iResBeg, | |
1793 | int iResEnd) { | |
1794 | ||
1795 | // Identity of mother of decaying reseonance(s). | |
1796 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
1797 | ||
1798 | // For Higgs decay hand over to standard routine. | |
1799 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
1800 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
1801 | ||
1802 | // For top decay hand over to standard routine. | |
1803 | if (idMother == 6) | |
1804 | return weightTopDecay( process, iResBeg, iResEnd); | |
1805 | ||
1806 | // Else done. | |
1807 | return 1.; | |
1808 | ||
1809 | } | |
1810 | ||
1811 | //========================================================================== | |
1812 | ||
1813 | // Sigma2qg2Hqlt class. | |
1814 | // Cross section for q g -> H q (H SM or BSM Higgs) via top loop. | |
1815 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1816 | ||
1817 | //-------------------------------------------------------------------------- | |
1818 | ||
1819 | // Initialize process. | |
1820 | ||
1821 | void Sigma2qg2Hqlt::initProc() { | |
1822 | ||
1823 | // Properties specific to Higgs state. | |
1824 | if (higgsType == 0) { | |
1825 | nameSave = "q g -> H q (SM; top loop)"; | |
1826 | codeSave = 915; | |
1827 | idRes = 25; | |
1828 | } | |
1829 | else if (higgsType == 1) { | |
1830 | nameSave = "q g -> h0(H1) q (BSM; top loop)"; | |
1831 | codeSave = 1015; | |
1832 | idRes = 25; | |
1833 | } | |
1834 | else if (higgsType == 2) { | |
1835 | nameSave = "q g -> H0(H2) q (BSM; top loop)"; | |
1836 | codeSave = 1035; | |
1837 | idRes = 35; | |
1838 | } | |
1839 | else if (higgsType == 3) { | |
1840 | nameSave = "q g -> A0(A3) q (BSM; top loop)"; | |
1841 | codeSave = 1055; | |
1842 | idRes = 36; | |
1843 | } | |
1844 | ||
1845 | // Normalization factor by g g -> H partial width. | |
1846 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1847 | double mHiggs = particleDataPtr->m0(idRes); | |
1848 | widHgg = particleDataPtr->resWidthChan(idRes, mHiggs, 21, 21); | |
1849 | ||
1850 | // Secondary open width fraction. | |
1851 | openFrac = particleDataPtr->resOpenFrac(idRes); | |
1852 | ||
1853 | } | |
1854 | ||
1855 | //-------------------------------------------------------------------------- | |
1856 | ||
1857 | // Evaluate sigmaHat(sHat, part independent of incoming flavour). | |
1858 | ||
1859 | void Sigma2qg2Hqlt::sigmaKin() { | |
1860 | ||
1861 | // Evaluate cross section. Secondary width for H0, H1, H2 or A3. | |
1862 | sigma = (M_PI / sH2) * (1. / 12.) * alpS * (widHgg / m3) | |
1863 | * (sH2 + uH2) / (-tH * s3); | |
1864 | sigma *= openFrac; | |
1865 | ||
1866 | } | |
1867 | ||
1868 | //-------------------------------------------------------------------------- | |
1869 | ||
1870 | // Select identity, colour and anticolour. | |
1871 | ||
1872 | void Sigma2qg2Hqlt::setIdColAcol() { | |
1873 | ||
1874 | // Flavour set up for q g -> H q. | |
1875 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1876 | int idq = (id2 == 21) ? id1 : id2; | |
1877 | setId( id1, id2, idRes, idq); | |
1878 | ||
1879 | // tH defined between f and f': must swap tHat <-> uHat if q g in. | |
1880 | swapTU = (id2 == 21); | |
1881 | ||
1882 | // Colour flow topologies. Swap when antiquarks. | |
1883 | if (id2 == 21) setColAcol( 1, 0, 2, 1, 0, 0, 2, 0); | |
1884 | else setColAcol( 2, 1, 1, 0, 0, 0, 2, 0); | |
1885 | if (idq < 0) swapColAcol(); | |
1886 | ||
1887 | } | |
1888 | ||
1889 | //-------------------------------------------------------------------------- | |
1890 | ||
1891 | // Evaluate weight for decay angles. | |
1892 | ||
1893 | double Sigma2qg2Hqlt::weightDecay( Event& process, int iResBeg, | |
1894 | int iResEnd) { | |
1895 | ||
1896 | // Identity of mother of decaying reseonance(s). | |
1897 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
1898 | ||
1899 | // For Higgs decay hand over to standard routine. | |
1900 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
1901 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
1902 | ||
1903 | // For top decay hand over to standard routine. | |
1904 | if (idMother == 6) | |
1905 | return weightTopDecay( process, iResBeg, iResEnd); | |
1906 | ||
1907 | // Else done. | |
1908 | return 1.; | |
1909 | ||
1910 | } | |
1911 | ||
1912 | //========================================================================== | |
1913 | ||
1914 | // Sigma2qqbar2Hglt class. | |
1915 | // Cross section for q qbar -> H g (H SM or BSM Higgs) via top loop. | |
1916 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1917 | ||
1918 | //-------------------------------------------------------------------------- | |
1919 | ||
1920 | // Initialize process. | |
1921 | ||
1922 | void Sigma2qqbar2Hglt::initProc() { | |
1923 | ||
1924 | // Properties specific to Higgs state. | |
1925 | if (higgsType == 0) { | |
1926 | nameSave = "q qbar -> H g (SM; top loop)"; | |
1927 | codeSave = 916; | |
1928 | idRes = 25; | |
1929 | } | |
1930 | else if (higgsType == 1) { | |
1931 | nameSave = "q qbar -> h0(H1) g (BSM; top loop)"; | |
1932 | codeSave = 1016; | |
1933 | idRes = 25; | |
1934 | } | |
1935 | else if (higgsType == 2) { | |
1936 | nameSave = "q qbar -> H0(H2) g (BSM; top loop)"; | |
1937 | codeSave = 1036; | |
1938 | idRes = 35; | |
1939 | } | |
1940 | else if (higgsType == 3) { | |
1941 | nameSave = "q qbar -> A0(A3) g (BSM; top loop)"; | |
1942 | codeSave = 1056; | |
1943 | idRes = 36; | |
1944 | } | |
1945 | ||
1946 | // Normalization factor by g g -> H partial width. | |
1947 | // (H can be H0 SM or H1, H2, A3 from BSM). | |
1948 | double mHiggs = particleDataPtr->m0(idRes); | |
1949 | widHgg = particleDataPtr->resWidthChan(idRes, mHiggs, 21, 21); | |
1950 | ||
1951 | // Secondary open width fraction. | |
1952 | openFrac = particleDataPtr->resOpenFrac(idRes); | |
1953 | ||
1954 | ||
1955 | } | |
1956 | ||
1957 | //-------------------------------------------------------------------------- | |
1958 | ||
1959 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
1960 | ||
1961 | void Sigma2qqbar2Hglt::sigmaKin() { | |
1962 | ||
1963 | // Evaluate cross section. Secondary width for H0, H1, H2 or A3. | |
1964 | sigma = (M_PI / sH2) * (2. / 9.) * alpS * (widHgg / m3) | |
1965 | * (tH2 + uH2) / (sH * s3); | |
1966 | sigma *= openFrac; | |
1967 | ||
1968 | } | |
1969 | ||
1970 | //-------------------------------------------------------------------------- | |
1971 | ||
1972 | // Select identity, colour and anticolour. | |
1973 | ||
1974 | void Sigma2qqbar2Hglt::setIdColAcol() { | |
1975 | ||
1976 | // Flavours trivial. | |
1977 | setId( id1, id2, idRes, 21); | |
1978 | ||
1979 | // Colour flow topologies. Swap when antiquarks. | |
1980 | setColAcol( 1, 0, 0, 2, 0, 0, 1, 2); | |
1981 | if (id1 < 0) swapColAcol(); | |
1982 | ||
1983 | } | |
1984 | ||
1985 | //-------------------------------------------------------------------------- | |
1986 | ||
1987 | // Evaluate weight for decay angles. | |
1988 | ||
1989 | double Sigma2qqbar2Hglt::weightDecay( Event& process, int iResBeg, | |
1990 | int iResEnd) { | |
1991 | ||
1992 | // Identity of mother of decaying reseonance(s). | |
1993 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
1994 | ||
1995 | // For Higgs decay hand over to standard routine. | |
1996 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
1997 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
1998 | ||
1999 | // For top decay hand over to standard routine. | |
2000 | if (idMother == 6) | |
2001 | return weightTopDecay( process, iResBeg, iResEnd); | |
2002 | ||
2003 | // Else done. | |
2004 | return 1.; | |
2005 | ||
2006 | } | |
2007 | ||
2008 | ||
2009 | //========================================================================== | |
2010 | ||
2011 | // Sigma1ffbar2Hchg class. | |
2012 | // Cross section for f fbar -> H+- (f is quark or lepton). | |
2013 | ||
2014 | //-------------------------------------------------------------------------- | |
2015 | ||
2016 | // Initialize process. | |
2017 | ||
2018 | void Sigma1ffbar2Hchg::initProc() { | |
2019 | ||
2020 | // Find pointer to H+-. | |
2021 | HResPtr = particleDataPtr->particleDataEntryPtr(37); | |
2022 | ||
2023 | // Store H+- mass and width for propagator. | |
2024 | mRes = HResPtr->m0(); | |
2025 | GammaRes = HResPtr->mWidth(); | |
2026 | m2Res = mRes*mRes; | |
2027 | GamMRat = GammaRes / mRes; | |
2028 | ||
2029 | // Couplings. | |
2030 | m2W = pow2(particleDataPtr->m0(24)); | |
2031 | thetaWRat = 1. / (8. * couplingsPtr->sin2thetaW()); | |
2032 | tan2Beta = pow2(settingsPtr->parm("HiggsHchg:tanBeta")); | |
2033 | ||
2034 | } | |
2035 | ||
2036 | //-------------------------------------------------------------------------- | |
2037 | ||
2038 | // Evaluate d(sigmaHat)/d(tHat), part independent of incoming flavour. | |
2039 | ||
2040 | void Sigma1ffbar2Hchg::sigmaKin() { | |
2041 | ||
2042 | // Set up Breit-Wigner. Width out only includes open channels. | |
2043 | sigBW = 4. * M_PI / ( pow2(sH - m2Res) + pow2(sH * GamMRat) ); | |
2044 | widthOutPos = HResPtr->resWidthOpen( 37, mH); | |
2045 | widthOutNeg = HResPtr->resWidthOpen(-37, mH); | |
2046 | ||
2047 | } | |
2048 | ||
2049 | //-------------------------------------------------------------------------- | |
2050 | ||
2051 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
2052 | ||
2053 | double Sigma1ffbar2Hchg::sigmaHat() { | |
2054 | ||
2055 | // Only allow generation-diagonal states. | |
2056 | int id1Abs = abs(id1); | |
2057 | int id2Abs = abs(id2); | |
2058 | int idUp = max(id1Abs, id2Abs); | |
2059 | int idDn = min(id1Abs, id2Abs); | |
2060 | if (idUp%2 != 0 || idUp - idDn != 1) return 0.; | |
2061 | ||
2062 | // Calculate mass-dependent incoming width. Total cross section. | |
2063 | double m2RunUp = pow2(particleDataPtr->mRun(idUp, mH)); | |
2064 | double m2RunDn = pow2(particleDataPtr->mRun(idDn, mH)); | |
2065 | double widthIn = alpEM * thetaWRat * (mH/m2W) | |
2066 | * (m2RunDn * tan2Beta + m2RunUp / tan2Beta); | |
2067 | int idUpChg = (id1Abs%2 == 0) ? id1 : id2; | |
2068 | double sigma = (idUpChg > 0) ? widthIn * sigBW * widthOutPos | |
2069 | : widthIn * sigBW * widthOutNeg; | |
2070 | ||
2071 | // Colour factor. Answer. | |
2072 | if (idUp < 9) sigma /= 3.; | |
2073 | return sigma; | |
2074 | ||
2075 | } | |
2076 | ||
2077 | //-------------------------------------------------------------------------- | |
2078 | ||
2079 | // Select identity, colour and anticolour. | |
2080 | ||
2081 | void Sigma1ffbar2Hchg::setIdColAcol() { | |
2082 | ||
2083 | // Charge of Higgs. Fill flavours. | |
2084 | int idUpChg = (abs(id1)%2 == 0) ? id1 : id2; | |
2085 | int idHchg = (idUpChg > 0) ? 37 : -37; | |
2086 | setId( id1, id2, idHchg); | |
2087 | ||
2088 | // Colour flow topologies. Swap when antiquarks. | |
2089 | if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0); | |
2090 | else setColAcol( 0, 0, 0, 0, 0, 0); | |
2091 | if (id1 < 0) swapColAcol(); | |
2092 | ||
2093 | } | |
2094 | ||
2095 | //-------------------------------------------------------------------------- | |
2096 | ||
2097 | // Evaluate weight for decay angles. | |
2098 | ||
2099 | double Sigma1ffbar2Hchg::weightDecay( Event& process, int iResBeg, | |
2100 | int iResEnd) { | |
2101 | ||
2102 | // Identity of mother of decaying reseonance(s). | |
2103 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
2104 | ||
2105 | // For Higgs decay hand over to standard routine. | |
2106 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
2107 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
2108 | ||
2109 | // For top decay hand over to standard routine. | |
2110 | if (idMother == 6) | |
2111 | return weightTopDecay( process, iResBeg, iResEnd); | |
2112 | ||
2113 | // Else done. | |
2114 | return 1.; | |
2115 | ||
2116 | } | |
2117 | ||
2118 | //========================================================================== | |
2119 | ||
2120 | // Sigma2qg2Hq class. | |
2121 | // Cross section for q g -> H+- q'. | |
2122 | ||
2123 | //-------------------------------------------------------------------------- | |
2124 | ||
2125 | // Initialize process. | |
2126 | ||
2127 | void Sigma2qg2Hchgq::initProc() { | |
2128 | ||
2129 | // Standard parameters. | |
2130 | m2W = pow2( particleDataPtr->m0(24) ); | |
2131 | thetaWRat = 1. / (24. * couplingsPtr->sin2thetaW()); | |
2132 | tan2Beta = pow2(settingsPtr->parm("HiggsHchg:tanBeta")); | |
2133 | ||
2134 | // Incoming flavour within same doublet. Uptype and downtype flavours. | |
2135 | idOld = (idNew%2 == 0) ? idNew - 1 : idNew + 1; | |
2136 | idUp = max(idOld, idNew); | |
2137 | idDn = min(idOld, idNew); | |
2138 | ||
2139 | // Secondary open width fraction. | |
2140 | openFracPos = (idOld%2 == 0) ? particleDataPtr->resOpenFrac( 37, idNew) | |
2141 | : particleDataPtr->resOpenFrac(-37, idNew); | |
2142 | openFracNeg = (idOld%2 == 0) ? particleDataPtr->resOpenFrac(-37, -idNew) | |
2143 | : particleDataPtr->resOpenFrac( 37, -idNew); | |
2144 | ||
2145 | } | |
2146 | ||
2147 | //-------------------------------------------------------------------------- | |
2148 | ||
2149 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
2150 | ||
2151 | void Sigma2qg2Hchgq::sigmaKin() { | |
2152 | ||
2153 | // Running masses provides coupling. | |
2154 | double m2RunUp = pow2(particleDataPtr->mRun(idUp, mH)); | |
2155 | double m2RunDn = pow2(particleDataPtr->mRun(idDn, mH)); | |
2156 | ||
2157 | // Cross section, including couplings and kinematics. | |
2158 | sigma = (M_PI / sH2) * alpS * alpEM * thetaWRat | |
2159 | * (m2RunDn * tan2Beta + m2RunUp / tan2Beta) / m2W | |
2160 | * (sH / (s4 - uH) + 2. * s4 * (s3 - uH) / pow2(s4 - uH) | |
2161 | + (s4 - uH) / sH - 2. * s4 / (s4 - uH) | |
2162 | + 2. * (s3 - uH) * (s3 - s4 - sH) / ((s4 - uH) * sH) ); | |
2163 | ||
2164 | } | |
2165 | ||
2166 | //-------------------------------------------------------------------------- | |
2167 | ||
2168 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
2169 | ||
2170 | double Sigma2qg2Hchgq::sigmaHat() { | |
2171 | ||
2172 | // Check that specified flavour present. | |
2173 | if (abs(id1) != idOld && abs(id2) != idOld) return 0.; | |
2174 | ||
2175 | // Answer. | |
2176 | return (id1 == idOld || id2 == idOld) ? sigma * openFracPos | |
2177 | : sigma * openFracNeg; | |
2178 | ||
2179 | } | |
2180 | ||
2181 | //-------------------------------------------------------------------------- | |
2182 | ||
2183 | // Select identity, colour and anticolour. | |
2184 | ||
2185 | void Sigma2qg2Hchgq::setIdColAcol() { | |
2186 | ||
2187 | // Flavour set up for q g -> H+- q'. | |
2188 | int idq = (id2 == 21) ? id1 : id2; | |
2189 | id3 = ( (idq > 0 && idOld%2 == 0) || (idq < 0 && idOld%2 != 0) ) | |
2190 | ? 37 : -37; | |
2191 | id4 = (idq > 0) ? idNew : -idNew; | |
2192 | setId( id1, id2, id3, id4); | |
2193 | ||
2194 | // tH defined between f and f': must swap tHat <-> uHat if q g in. | |
2195 | swapTU = (id2 == 21); | |
2196 | ||
2197 | // Colour flow topologies. Swap when antiquarks. | |
2198 | if (id2 == 21) setColAcol( 1, 0, 2, 1, 0, 0, 2, 0); | |
2199 | else setColAcol( 2, 1, 1, 0, 0, 0, 2, 0); | |
2200 | if (idq < 0) swapColAcol(); | |
2201 | ||
2202 | } | |
2203 | ||
2204 | //-------------------------------------------------------------------------- | |
2205 | ||
2206 | // Evaluate weight for decay angles. | |
2207 | ||
2208 | double Sigma2qg2Hchgq::weightDecay( Event& process, int iResBeg, | |
2209 | int iResEnd) { | |
2210 | ||
2211 | // Identity of mother of decaying reseonance(s). | |
2212 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
2213 | ||
2214 | // For Higgs decay hand over to standard routine. | |
2215 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
2216 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
2217 | ||
2218 | // For top decay hand over to standard routine. | |
2219 | if (idMother == 6) | |
2220 | return weightTopDecay( process, iResBeg, iResEnd); | |
2221 | ||
2222 | // Else done. | |
2223 | return 1.; | |
2224 | ||
2225 | } | |
2226 | ||
2227 | //========================================================================== | |
2228 | ||
2229 | // Sigma2ffbar2A3H12 class. | |
2230 | // Cross section for f fbar -> A0(H_3) h0(H_1) or A0(H_3) H0(H_2). | |
2231 | ||
2232 | //-------------------------------------------------------------------------- | |
2233 | ||
2234 | // Initialize process. | |
2235 | ||
2236 | void Sigma2ffbar2A3H12::initProc() { | |
2237 | ||
2238 | // Set up whether h0(H_1) or H0(H_2). | |
2239 | higgs12 = (higgsType == 1) ? 25 : 35; | |
2240 | codeSave = (higgsType == 1) ? 1081 : 1082; | |
2241 | nameSave = (higgsType == 1) ? "f fbar -> A0(H3) h0(H1)" | |
2242 | : "f fbar -> A0(H3) H0(H2)"; | |
2243 | coupZA3H12 = (higgsType == 1) ? settingsPtr->parm("HiggsA3:coup2H1Z") | |
2244 | : settingsPtr->parm("HiggsA3:coup2H2Z"); | |
2245 | ||
2246 | // Standard parameters. | |
2247 | double mZ = particleDataPtr->m0(23); | |
2248 | double GammaZ = particleDataPtr->mWidth(23); | |
2249 | m2Z = mZ * mZ; | |
2250 | mGammaZ = mZ * GammaZ; | |
2251 | thetaWRat = 1. / (4. * couplingsPtr->sin2thetaW() * couplingsPtr->cos2thetaW()); | |
2252 | ||
2253 | // Secondary open width fraction. | |
2254 | openFrac = particleDataPtr->resOpenFrac(36, higgs12); | |
2255 | ||
2256 | } | |
2257 | ||
2258 | //-------------------------------------------------------------------------- | |
2259 | ||
2260 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
2261 | ||
2262 | void Sigma2ffbar2A3H12::sigmaKin() { | |
2263 | ||
2264 | // Common kinematics factora. | |
2265 | sigma0 = (M_PI / sH2) * pow2(alpEM * thetaWRat * coupZA3H12) | |
2266 | * (uH * tH - s3 * s4) / ( pow2(sH - m2Z) + pow2(mGammaZ) ); | |
2267 | ||
2268 | } | |
2269 | ||
2270 | //-------------------------------------------------------------------------- | |
2271 | ||
2272 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
2273 | ||
2274 | double Sigma2ffbar2A3H12::sigmaHat() { | |
2275 | ||
2276 | // Couplings for incoming flavour. | |
2277 | int idAbs = abs(id1); | |
2278 | double lIn = couplingsPtr->lf(idAbs); | |
2279 | double rIn = couplingsPtr->rf(idAbs); | |
2280 | ||
2281 | // Combine to total cross section. Colour factor. | |
2282 | double sigma = (pow2(lIn) + pow2(rIn)) * sigma0 * openFrac; | |
2283 | if (idAbs < 9) sigma /= 3.; | |
2284 | return sigma; | |
2285 | ||
2286 | } | |
2287 | ||
2288 | //-------------------------------------------------------------------------- | |
2289 | ||
2290 | // Select identity, colour and anticolour. | |
2291 | ||
2292 | void Sigma2ffbar2A3H12::setIdColAcol() { | |
2293 | ||
2294 | // Flavours trivial | |
2295 | setId( id1, id2, 36, higgs12); | |
2296 | ||
2297 | // Colour flow topologies. Swap when antiquarks. | |
2298 | if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0); | |
2299 | else setColAcol( 0, 0, 0, 0, 0, 0); | |
2300 | if (id1 < 0) swapColAcol(); | |
2301 | ||
2302 | } | |
2303 | ||
2304 | //-------------------------------------------------------------------------- | |
2305 | ||
2306 | // Evaluate weight for decay angles. | |
2307 | ||
2308 | double Sigma2ffbar2A3H12::weightDecay( Event& process, int iResBeg, | |
2309 | int iResEnd) { | |
2310 | ||
2311 | // Identity of mother of decaying reseonance(s). | |
2312 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
2313 | ||
2314 | // For Higgs decay hand over to standard routine. | |
2315 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
2316 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
2317 | ||
2318 | // For top decay hand over to standard routine. | |
2319 | if (idMother == 6) | |
2320 | return weightTopDecay( process, iResBeg, iResEnd); | |
2321 | ||
2322 | // Else done. | |
2323 | return 1.; | |
2324 | ||
2325 | } | |
2326 | ||
2327 | //========================================================================== | |
2328 | ||
2329 | // Sigma2ffbar2HchgH12 class. | |
2330 | // Cross section for f fbar -> H+- h0(H_1) or H+- H0(H_2). | |
2331 | ||
2332 | //-------------------------------------------------------------------------- | |
2333 | ||
2334 | // Initialize process. | |
2335 | ||
2336 | void Sigma2ffbar2HchgH12::initProc() { | |
2337 | ||
2338 | // Set up whether h0(H_1) or H0(H_2). | |
2339 | higgs12 = (higgsType == 1) ? 25 : 35; | |
2340 | codeSave = (higgsType == 1) ? 1083 : 1084; | |
2341 | nameSave = (higgsType == 1) ? "f fbar' -> H+- h0(H1)" | |
2342 | : "f fbar' -> H+- H0(H2)"; | |
2343 | coupWHchgH12 = (higgsType == 1) ? settingsPtr->parm("HiggsHchg:coup2H1W") | |
2344 | : settingsPtr->parm("HiggsHchg:coup2H2W"); | |
2345 | ||
2346 | // Standard parameters. | |
2347 | double mW = particleDataPtr->m0(24); | |
2348 | double GammaW = particleDataPtr->mWidth(24); | |
2349 | m2W = mW * mW; | |
2350 | mGammaW = mW * GammaW; | |
2351 | thetaWRat = 1. / (2. * couplingsPtr->sin2thetaW()); | |
2352 | ||
2353 | // Secondary open width fraction. | |
2354 | openFracPos = particleDataPtr->resOpenFrac( 37, higgs12); | |
2355 | openFracNeg = particleDataPtr->resOpenFrac(-37, higgs12); | |
2356 | ||
2357 | } | |
2358 | ||
2359 | //-------------------------------------------------------------------------- | |
2360 | ||
2361 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
2362 | ||
2363 | void Sigma2ffbar2HchgH12::sigmaKin() { | |
2364 | ||
2365 | // Common kinematics factora. | |
2366 | sigma0 = 0.5 * (M_PI / sH2) * pow2(alpEM * thetaWRat * coupWHchgH12) | |
2367 | * (uH * tH - s3 * s4) / ( pow2(sH - m2W) + pow2(mGammaW) ); | |
2368 | ||
2369 | } | |
2370 | ||
2371 | //-------------------------------------------------------------------------- | |
2372 | ||
2373 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
2374 | ||
2375 | double Sigma2ffbar2HchgH12::sigmaHat() { | |
2376 | ||
2377 | // Combine to total cross section. CKM and colour factor. | |
2378 | int idUp = (abs(id1)%2 == 0) ? id1 : id2; | |
2379 | double sigma = (idUp > 0) ? sigma0 * openFracPos : sigma0 * openFracNeg; | |
2380 | if (abs(id1) < 9) sigma *= couplingsPtr->V2CKMid(abs(id1), abs(id2)) / 3.; | |
2381 | return sigma; | |
2382 | ||
2383 | } | |
2384 | ||
2385 | //-------------------------------------------------------------------------- | |
2386 | ||
2387 | // Select identity, colour and anticolour. | |
2388 | ||
2389 | void Sigma2ffbar2HchgH12::setIdColAcol() { | |
2390 | ||
2391 | // Charge of Higgs. Fill flavours. | |
2392 | int idUpChg = (abs(id1)%2 == 0) ? id1 : id2; | |
2393 | int idHchg = (idUpChg > 0) ? 37 : -37; | |
2394 | setId( id1, id2, idHchg, higgs12); | |
2395 | ||
2396 | // Colour flow topologies. Swap when antiquarks. | |
2397 | if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0); | |
2398 | else setColAcol( 0, 0, 0, 0, 0, 0); | |
2399 | if (id1 < 0) swapColAcol(); | |
2400 | ||
2401 | } | |
2402 | ||
2403 | //-------------------------------------------------------------------------- | |
2404 | ||
2405 | // Evaluate weight for decay angles. | |
2406 | ||
2407 | double Sigma2ffbar2HchgH12::weightDecay( Event& process, int iResBeg, | |
2408 | int iResEnd) { | |
2409 | ||
2410 | // Identity of mother of decaying reseonance(s). | |
2411 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
2412 | ||
2413 | // For Higgs decay hand over to standard routine. | |
2414 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
2415 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
2416 | ||
2417 | // For top decay hand over to standard routine. | |
2418 | if (idMother == 6) | |
2419 | return weightTopDecay( process, iResBeg, iResEnd); | |
2420 | ||
2421 | // Else done. | |
2422 | return 1.; | |
2423 | ||
2424 | } | |
2425 | ||
2426 | //========================================================================== | |
2427 | ||
2428 | // Sigma2ffbar2HposHneg class. | |
2429 | // Cross section for q g -> H+- q'. | |
2430 | ||
2431 | //-------------------------------------------------------------------------- | |
2432 | ||
2433 | // Initialize process. | |
2434 | ||
2435 | void Sigma2ffbar2HposHneg::initProc() { | |
2436 | ||
2437 | // Standard parameters. | |
2438 | double mZ = particleDataPtr->m0(23); | |
2439 | double GammaZ = particleDataPtr->mWidth(23); | |
2440 | m2Z = mZ * mZ; | |
2441 | mGammaZ = mZ * GammaZ; | |
2442 | thetaWRat = 1. / (4. * couplingsPtr->sin2thetaW() * couplingsPtr->cos2thetaW()); | |
2443 | ||
2444 | // Charged Higgs coupling to gamma and Z0. | |
2445 | eH = -1.; | |
2446 | lH = -1. + 2. * couplingsPtr->sin2thetaW(); | |
2447 | ||
2448 | // Secondary open width fraction. | |
2449 | openFrac = particleDataPtr->resOpenFrac(37, -37); | |
2450 | ||
2451 | } | |
2452 | ||
2453 | //-------------------------------------------------------------------------- | |
2454 | ||
2455 | // Evaluate sigmaHat(sHat), part independent of incoming flavour. | |
2456 | ||
2457 | void Sigma2ffbar2HposHneg::sigmaKin() { | |
2458 | ||
2459 | // Common kinematics factora. | |
2460 | double preFac = M_PI * pow2(alpEM) * ((uH * tH - s3 * s4) / sH2); | |
2461 | double propZ = 1. / ( pow2(sH - m2Z) + pow2(mGammaZ) ); | |
2462 | ||
2463 | // Separate parts for gamma*, interference and Z0. | |
2464 | gamSig = preFac * 2. * pow2(eH) / sH2; | |
2465 | intSig = preFac * 2. * eH * lH * thetaWRat * propZ * (sH - m2Z) / sH; | |
2466 | resSig = preFac * pow2(lH * thetaWRat) * propZ; | |
2467 | ||
2468 | } | |
2469 | ||
2470 | //-------------------------------------------------------------------------- | |
2471 | ||
2472 | // Evaluate sigmaHat(sHat), including incoming flavour dependence. | |
2473 | ||
2474 | double Sigma2ffbar2HposHneg::sigmaHat() { | |
2475 | ||
2476 | // Couplings for incoming flavour. | |
2477 | int idAbs = int(id1); | |
2478 | double eIn = couplingsPtr->ef(idAbs); | |
2479 | double lIn = couplingsPtr->lf(idAbs); | |
2480 | double rIn = couplingsPtr->rf(idAbs); | |
2481 | ||
2482 | // Combine to total cross section. Colour factor. | |
2483 | double sigma = (pow2(eIn) * gamSig + eIn * (lIn + rIn) * intSig | |
2484 | + (pow2(lIn) + pow2(rIn)) * resSig) * openFrac; | |
2485 | if (idAbs < 9) sigma /= 3.; | |
2486 | return sigma; | |
2487 | ||
2488 | } | |
2489 | ||
2490 | //-------------------------------------------------------------------------- | |
2491 | ||
2492 | // Select identity, colour and anticolour. | |
2493 | ||
2494 | void Sigma2ffbar2HposHneg::setIdColAcol() { | |
2495 | ||
2496 | // Flavours trivial | |
2497 | setId( id1, id2, 37, -37); | |
2498 | ||
2499 | // Colour flow topologies. Swap when antiquarks. | |
2500 | if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0); | |
2501 | else setColAcol( 0, 0, 0, 0, 0, 0); | |
2502 | if (id1 < 0) swapColAcol(); | |
2503 | ||
2504 | } | |
2505 | ||
2506 | //-------------------------------------------------------------------------- | |
2507 | ||
2508 | // Evaluate weight for decay angles. | |
2509 | ||
2510 | double Sigma2ffbar2HposHneg::weightDecay( Event& process, int iResBeg, | |
2511 | int iResEnd) { | |
2512 | ||
2513 | // Identity of mother of decaying reseonance(s). | |
2514 | int idMother = process[process[iResBeg].mother1()].idAbs(); | |
2515 | ||
2516 | // For Higgs decay hand over to standard routine. | |
2517 | if (idMother == 25 || idMother == 35 || idMother == 36) | |
2518 | return weightHiggsDecay( process, iResBeg, iResEnd); | |
2519 | ||
2520 | // For top decay hand over to standard routine. | |
2521 | if (idMother == 6) | |
2522 | return weightTopDecay( process, iResBeg, iResEnd); | |
2523 | ||
2524 | // Else done. | |
2525 | return 1.; | |
2526 | ||
2527 | } | |
2528 | ||
2529 | //========================================================================== | |
2530 | ||
2531 | } // end namespace Pythia8 |