Line data Source code
1 : // SigmaLeftRightSym.cc is a part of the PYTHIA event generator.
2 : // Copyright (C) 2015 Torbjorn Sjostrand.
3 : // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
4 : // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 :
6 : // Function definitions (not found in the header) for the
7 : // left-right-symmetry simulation classes.
8 :
9 : #include "Pythia8/SigmaLeftRightSym.h"
10 :
11 : namespace Pythia8 {
12 :
13 : //==========================================================================
14 :
15 : // Sigma1ffbar2ZRight class.
16 : // Cross section for f fbar -> Z_R^0 (righthanded gauge boson).
17 :
18 : //--------------------------------------------------------------------------
19 :
20 : // Initialize process.
21 :
22 : void Sigma1ffbar2ZRight::initProc() {
23 :
24 : // Store Z_R mass and width for propagator.
25 0 : idZR = 9900023;
26 0 : mRes = particleDataPtr->m0(idZR);
27 0 : GammaRes = particleDataPtr->mWidth(idZR);
28 0 : m2Res = mRes*mRes;
29 0 : GamMRat = GammaRes / mRes;
30 0 : sin2tW = couplingsPtr->sin2thetaW();
31 :
32 : // Set pointer to particle properties and decay table.
33 0 : ZRPtr = particleDataPtr->particleDataEntryPtr(idZR);
34 :
35 0 : }
36 :
37 : //--------------------------------------------------------------------------
38 :
39 : // Evaluate sigmaHat(sHat), part independent of incoming flavour.
40 :
41 : void Sigma1ffbar2ZRight::sigmaKin() {
42 :
43 : // Set up Breit-Wigner. Width out only includes open channels.
44 0 : double sigBW = 12. * M_PI/ ( pow2(sH - m2Res) + pow2(sH * GamMRat) );
45 0 : double widthOut = ZRPtr->resWidthOpen(idZR, mH);
46 :
47 : // Prefactor for incoming widths. Combine. Done.
48 0 : double preFac = alpEM * mH / ( 48. * sin2tW * (1. - sin2tW)
49 0 : * (1. - 2. * sin2tW) );
50 0 : sigma0 = preFac * sigBW * widthOut;
51 :
52 0 : }
53 :
54 : //--------------------------------------------------------------------------
55 :
56 : // Evaluate sigmaHat(sHat), including incoming flavour dependence.
57 :
58 : double Sigma1ffbar2ZRight::sigmaHat() {
59 :
60 : // Vector and axial couplings of incoming fermion pair.
61 0 : int idAbs = abs(id1);
62 : double af = 0.;
63 : double vf = 0.;
64 0 : if (idAbs < 9 && idAbs%2 == 1) {
65 0 : af = -1. + 2. * sin2tW;
66 0 : vf = -1. + 4. * sin2tW / 3.;
67 0 : } else if (idAbs < 9) {
68 0 : af = 1. - 2. * sin2tW;
69 0 : vf = 1. - 8. * sin2tW / 3.;
70 0 : } else if (idAbs < 19 && idAbs%2 == 1) {
71 0 : af = -1. + 2. * sin2tW;
72 0 : vf = -1. + 4. * sin2tW;
73 0 : }
74 :
75 : // Colour factor. Answer.
76 0 : double sigma = (vf*vf + af*af) * sigma0;
77 0 : if (idAbs < 9) sigma /= 3.;
78 0 : return sigma;
79 :
80 : }
81 :
82 : //--------------------------------------------------------------------------
83 :
84 : // Select identity, colour and anticolour.
85 :
86 : void Sigma1ffbar2ZRight::setIdColAcol() {
87 :
88 : // Flavours trivial.
89 0 : setId( id1, id2, idZR);
90 :
91 : // Colour flow topologies. Swap when antiquarks.
92 0 : if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0);
93 0 : else setColAcol( 0, 0, 0, 0, 0, 0);
94 0 : if (id1 < 0) swapColAcol();
95 :
96 0 : }
97 :
98 : //--------------------------------------------------------------------------
99 :
100 : // Evaluate weight for Z_R decay angle.
101 :
102 : double Sigma1ffbar2ZRight::weightDecay( Event& process, int iResBeg,
103 : int iResEnd) {
104 :
105 : // Identity of mother of decaying reseonance(s).
106 0 : int idMother = process[process[iResBeg].mother1()].idAbs();
107 :
108 : // For top decay hand over to standard routine.
109 0 : if (idMother == 6)
110 0 : return weightTopDecay( process, iResBeg, iResEnd);
111 :
112 : // Z_R should sit in entry 5.
113 0 : if (iResBeg != 5 || iResEnd != 5) return 1.;
114 :
115 : // Couplings for in- and out-flavours.
116 : double ai, vi, af, vf;
117 0 : int idInAbs = process[3].idAbs();
118 0 : if (idInAbs < 9 && idInAbs%2 == 1) {
119 0 : ai = -1. + 2. * sin2tW;
120 0 : vi = -1. + 4. * sin2tW / 3.;
121 0 : } else if (idInAbs < 9) {
122 0 : ai = 1. - 2. * sin2tW;
123 0 : vi = 1. - 8. * sin2tW / 3.;
124 0 : } else {
125 0 : ai = -1. + 2. * sin2tW;
126 0 : vi = -1. + 4. * sin2tW;
127 : }
128 0 : int idOutAbs = process[6].idAbs();
129 0 : if (idOutAbs < 9 && idOutAbs%2 == 1) {
130 0 : af = -1. + 2. * sin2tW;
131 0 : vf = -1. + 4. * sin2tW / 3.;
132 0 : } else if (idOutAbs < 9) {
133 0 : af = 1. - 2. * sin2tW;
134 0 : vf = 1. - 8. * sin2tW / 3.;
135 0 : } else {
136 0 : af = -1. + 2. * sin2tW;
137 0 : vf = -1. + 4. * sin2tW;
138 : }
139 :
140 : // Phase space factors. Reconstruct decay angle.
141 0 : double mr1 = pow2(process[6].m()) / sH;
142 0 : double mr2 = pow2(process[7].m()) / sH;
143 0 : double betaf = sqrtpos( pow2(1. - mr1 - mr2) - 4. * mr1 * mr2);
144 0 : double cosThe = (process[3].p() - process[4].p())
145 0 : * (process[7].p() - process[6].p()) / (sH * betaf);
146 :
147 : // Angular weight and its maximum.
148 0 : double wt1 = (vi*vi + ai*ai) * (vf*vf + af*af * betaf*betaf);
149 0 : double wt2 = (1. - betaf*betaf) * (vi*vi + ai*ai) * vf*vf;
150 0 : double wt3 = betaf * 4. * vi * ai * vf * af;
151 0 : if (process[3].id() * process[6].id() < 0) wt3 = -wt3;
152 0 : double wt = wt1 * (1. + cosThe*cosThe) + wt2 * (1. - cosThe*cosThe)
153 0 : + 2. * wt3 * cosThe;
154 0 : double wtMax = 2. * (wt1 + abs(wt3));
155 :
156 : // Done.
157 0 : return wt / wtMax;
158 :
159 0 : }
160 :
161 : //==========================================================================
162 :
163 : // Sigma1ffbar2WRight class.
164 : // Cross section for f fbar' -> W_R^+- (righthanded gauge boson).
165 :
166 : //--------------------------------------------------------------------------
167 :
168 : // Initialize process.
169 :
170 : void Sigma1ffbar2WRight::initProc() {
171 :
172 : // Store W_R^+- mass and width for propagator.
173 0 : idWR = 9900024;
174 0 : mRes = particleDataPtr->m0(idWR);
175 0 : GammaRes = particleDataPtr->mWidth(idWR);
176 0 : m2Res = mRes*mRes;
177 0 : GamMRat = GammaRes / mRes;
178 0 : thetaWRat = 1. / (12. * couplingsPtr->sin2thetaW());
179 :
180 : // Set pointer to particle properties and decay table.
181 0 : particlePtr = particleDataPtr->particleDataEntryPtr(idWR);
182 :
183 0 : }
184 :
185 : //--------------------------------------------------------------------------
186 :
187 : // Evaluate sigmaHat(sHat), part independent of incoming flavour.
188 :
189 : void Sigma1ffbar2WRight::sigmaKin() {
190 :
191 : // Common coupling factors.
192 0 : double colQ = 3. * (1. + alpS / M_PI);
193 :
194 : // Reset quantities to sum. Declare variables inside loop.
195 : double widOutPos = 0.;
196 : double widOutNeg = 0.;
197 : int id1Now, id2Now, id1Abs, id2Abs, id1Neg, id2Neg, onMode;
198 : double widNow, widSecPos, widSecNeg, mf1, mf2, mr1, mr2, kinFac;
199 :
200 : // Loop over all W_R^+- decay channels.
201 0 : for (int i = 0; i < particlePtr->sizeChannels(); ++i) {
202 0 : id1Now = particlePtr->channel(i).product(0);
203 0 : id2Now = particlePtr->channel(i).product(1);
204 0 : id1Abs = abs(id1Now);
205 0 : id2Abs = abs(id2Now);
206 :
207 : // Check that above threshold. Phase space.
208 0 : mf1 = particleDataPtr->m0(id1Abs);
209 0 : mf2 = particleDataPtr->m0(id2Abs);
210 0 : if (mH > mf1 + mf2 + MASSMARGIN) {
211 0 : mr1 = pow2(mf1 / mH);
212 0 : mr2 = pow2(mf2 / mH);
213 0 : kinFac = (1. - 0.5 * (mr1 + mr2) - 0.5 * pow2(mr1 - mr2))
214 0 : * sqrtpos( pow2(1. - mr1 - mr2) - 4. * mr1 * mr2 );
215 :
216 : // Combine kinematics with colour factor and CKM couplings.
217 : widNow = kinFac;
218 0 : if (id1Abs < 9) widNow *= colQ * couplingsPtr->V2CKMid(id1Abs, id2Abs);
219 :
220 : // Secondary width from top and righthanded neutrino decay.
221 0 : id1Neg = (id1Abs < 19) ? -id1Now : id1Abs;
222 0 : id2Neg = (id2Abs < 19) ? -id2Now : id2Abs;
223 0 : widSecPos = particleDataPtr->resOpenFrac(id1Now, id2Now);
224 0 : widSecNeg = particleDataPtr->resOpenFrac(id1Neg, id2Neg);
225 :
226 : // Add weight for channels on for all, W_R^+ and W_R^-, respectively.
227 0 : onMode = particlePtr->channel(i).onMode();
228 0 : if (onMode == 1 || onMode == 2) widOutPos += widNow * widSecPos;
229 0 : if (onMode == 1 || onMode == 3) widOutNeg += widNow * widSecNeg;
230 :
231 : // End loop over fermions.
232 : }
233 : }
234 :
235 : // Set up Breit-Wigner. Cross section for W_R^+ and W_R^- separately.
236 0 : double sigBW = 12. * M_PI * pow2(alpEM * thetaWRat) * sH
237 0 : / ( pow2(sH - m2Res) + pow2(sH * GamMRat) );
238 0 : sigma0Pos = sigBW * widOutPos;
239 0 : sigma0Neg = sigBW * widOutNeg;
240 :
241 0 : }
242 :
243 : //--------------------------------------------------------------------------
244 :
245 : // Evaluate sigmaHat(sHat), including incoming flavour dependence.
246 :
247 : double Sigma1ffbar2WRight::sigmaHat() {
248 :
249 : // Secondary width for W_R^+ or W_R^-. CKM and colour factors.
250 0 : int idUp = (abs(id1)%2 == 0) ? id1 : id2;
251 0 : double sigma = (idUp > 0) ? sigma0Pos : sigma0Neg;
252 0 : if (abs(id1) < 9) sigma *= couplingsPtr->V2CKMid(abs(id1), abs(id2)) / 3.;
253 :
254 : // Answer.
255 0 : return sigma;
256 :
257 : }
258 :
259 : //--------------------------------------------------------------------------
260 :
261 : // Select identity, colour and anticolour.
262 :
263 : void Sigma1ffbar2WRight::setIdColAcol() {
264 :
265 : // Sign of outgoing W_R.
266 0 : int sign = (abs(id1)%2 == 0) ? 1 : -1;
267 0 : if (id1 < 0) sign = -sign;
268 0 : setId( id1, id2, idWR * sign);
269 :
270 : // Colour flow topologies. Swap when antiquarks.
271 0 : if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0);
272 0 : else setColAcol( 0, 0, 0, 0, 0, 0);
273 0 : if (id1 < 0) swapColAcol();
274 :
275 0 : }
276 :
277 : //--------------------------------------------------------------------------
278 :
279 : // Evaluate weight for W_R decay angle.
280 :
281 : double Sigma1ffbar2WRight::weightDecay( Event& process, int iResBeg,
282 : int iResEnd) {
283 :
284 : // Identity of mother of decaying reseonance(s).
285 0 : int idMother = process[process[iResBeg].mother1()].idAbs();
286 :
287 : // For top decay hand over to standard routine.
288 0 : if (idMother == 6)
289 0 : return weightTopDecay( process, iResBeg, iResEnd);
290 :
291 : // W_R should sit in entry 5.
292 0 : if (iResBeg != 5 || iResEnd != 5) return 1.;
293 :
294 : // Phase space factors.
295 0 : double mr1 = pow2(process[6].m()) / sH;
296 0 : double mr2 = pow2(process[7].m()) / sH;
297 0 : double betaf = sqrtpos( pow2(1. - mr1 - mr2) - 4. * mr1 * mr2);
298 :
299 : // Sign of asymmetry.
300 0 : double eps = (process[3].id() * process[6].id() > 0) ? 1. : -1.;
301 :
302 : // Reconstruct decay angle and weight for it.
303 0 : double cosThe = (process[3].p() - process[4].p())
304 0 : * (process[7].p() - process[6].p()) / (sH * betaf);
305 : double wtMax = 4.;
306 0 : double wt = pow2(1. + betaf * eps * cosThe) - pow2(mr1 - mr2);
307 :
308 : // Done.
309 0 : return (wt / wtMax);
310 :
311 0 : }
312 :
313 : //==========================================================================
314 :
315 : // Sigma1ll2Hchgchg class.
316 : // Cross section for l l -> H_L^++-- or H_R^++-- (doubly charged Higgs).
317 :
318 : //--------------------------------------------------------------------------
319 :
320 : // Initialize process.
321 :
322 : void Sigma1ll2Hchgchg::initProc() {
323 :
324 : // Set process properties: H_L^++-- or H_R^++--.
325 0 : if (leftRight == 1) {
326 0 : idHLR = 9900041;
327 0 : codeSave = 3121;
328 0 : nameSave = "l l -> H_L^++--";
329 0 : } else {
330 0 : idHLR = 9900042;
331 0 : codeSave = 3141;
332 0 : nameSave = "l l -> H_R^++--";
333 : }
334 :
335 : // Read in Yukawa matrix for couplings to a lepton pair.
336 0 : yukawa[1][1] = settingsPtr->parm("LeftRightSymmmetry:coupHee");
337 0 : yukawa[2][1] = settingsPtr->parm("LeftRightSymmmetry:coupHmue");
338 0 : yukawa[2][2] = settingsPtr->parm("LeftRightSymmmetry:coupHmumu");
339 0 : yukawa[3][1] = settingsPtr->parm("LeftRightSymmmetry:coupHtaue");
340 0 : yukawa[3][2] = settingsPtr->parm("LeftRightSymmmetry:coupHtaumu");
341 0 : yukawa[3][3] = settingsPtr->parm("LeftRightSymmmetry:coupHtautau");
342 :
343 : // Store H_L/R mass and width for propagator.
344 0 : mRes = particleDataPtr->m0(idHLR);
345 0 : GammaRes = particleDataPtr->mWidth(idHLR);
346 0 : m2Res = mRes*mRes;
347 0 : GamMRat = GammaRes / mRes;
348 :
349 : // Set pointer to particle properties and decay table.
350 0 : particlePtr = particleDataPtr->particleDataEntryPtr(idHLR);
351 :
352 0 : }
353 :
354 : //--------------------------------------------------------------------------
355 :
356 : // Evaluate sigmaHat(sHat), including incoming flavour dependence.
357 :
358 : double Sigma1ll2Hchgchg::sigmaHat() {
359 :
360 : // Initial state must consist of two identical-sign leptons.
361 0 : if (id1 * id2 < 0) return 0.;
362 0 : int id1Abs = abs(id1);
363 0 : int id2Abs = abs(id2);
364 0 : if (id1Abs != 11 && id1Abs != 13 && id1Abs != 15) return 0.;
365 0 : if (id2Abs != 11 && id2Abs != 13 && id2Abs != 15) return 0.;
366 :
367 : // Set up Breit-Wigner, inwidth and outwidth.
368 0 : double sigBW = 8. * M_PI / ( pow2(sH - m2Res) + pow2(sH * GamMRat) );
369 0 : double widIn = pow2(yukawa[(id1Abs-9)/2][(id2Abs-9)/2])
370 0 : * mH / (8. * M_PI);
371 0 : int idSgn = (id1 < 0) ? idHLR : -idHLR;
372 0 : double widOut = particlePtr->resWidthOpen( idSgn, mH);
373 :
374 : // Answer.
375 0 : return widIn * sigBW * widOut;
376 :
377 0 : }
378 :
379 : //--------------------------------------------------------------------------
380 :
381 : // Select identity, colour and anticolour.
382 :
383 : void Sigma1ll2Hchgchg::setIdColAcol() {
384 :
385 : // Sign of outgoing H_L/R.
386 0 : int idSgn = (id1 < 0) ? idHLR : -idHLR;
387 0 : setId( id1, id2, idSgn);
388 :
389 : // No colours whatsoever.
390 0 : setColAcol( 0, 0, 0, 0, 0, 0);
391 :
392 0 : }
393 :
394 : //--------------------------------------------------------------------------
395 :
396 : // Evaluate weight for H_L/R sequential decay angles.
397 :
398 : double Sigma1ll2Hchgchg::weightDecay( Event& process, int iResBeg,
399 : int iResEnd) {
400 :
401 : // Identity of mother of decaying reseonance(s).
402 0 : int idMother = process[process[iResBeg].mother1()].idAbs();
403 :
404 : // For top decay hand over to standard routine.
405 0 : if (idMother == 6)
406 0 : return weightTopDecay( process, iResBeg, iResEnd);
407 :
408 : // Else isotropic decay.
409 0 : return 1.;
410 :
411 0 : }
412 :
413 : //==========================================================================
414 :
415 : // Sigma2lgm2Hchgchgl class.
416 : // Cross section for l gamma -> H_L^++-- l or H_R^++-- l
417 : // (doubly charged Higgs).
418 :
419 : //--------------------------------------------------------------------------
420 :
421 : // Initialize process.
422 :
423 : void Sigma2lgm2Hchgchgl::initProc() {
424 :
425 : // Set process properties: H_L^++-- or H_R^++-- and e/mu/tau.
426 0 : idHLR = (leftRight == 1) ? 9900041 : 9900042;
427 0 : codeSave = (leftRight == 1) ? 3122 : 3142;
428 0 : if (idLep == 13) codeSave += 2;
429 0 : if (idLep == 15) codeSave += 4;
430 0 : if (codeSave == 3122) nameSave = "l^+- gamma -> H_L^++-- e^-+";
431 0 : else if (codeSave == 3123) nameSave = "l^+- gamma -> H_L^++-- mu^-+";
432 0 : else if (codeSave == 3124) nameSave = "l^+- gamma -> H_L^++-- tau^-+";
433 0 : else if (codeSave == 3142) nameSave = "l^+- gamma -> H_R^++-- e^-+";
434 0 : else if (codeSave == 3143) nameSave = "l^+- gamma -> H_R^++-- mu^-+";
435 0 : else nameSave = "l^+- gamma -> H_R^++-- tau^-+";
436 :
437 : // Read in relevantYukawa matrix for couplings to a lepton pair.
438 0 : if (idLep == 11) {
439 0 : yukawa[1] = settingsPtr->parm("LeftRightSymmmetry:coupHee");
440 0 : yukawa[2] = settingsPtr->parm("LeftRightSymmmetry:coupHmue");
441 0 : yukawa[3] = settingsPtr->parm("LeftRightSymmmetry:coupHtaue");
442 0 : } else if (idLep == 13) {
443 0 : yukawa[1] = settingsPtr->parm("LeftRightSymmmetry:coupHmue");
444 0 : yukawa[2] = settingsPtr->parm("LeftRightSymmmetry:coupHmumu");
445 0 : yukawa[3] = settingsPtr->parm("LeftRightSymmmetry:coupHtaumu");
446 0 : } else {
447 0 : yukawa[1] = settingsPtr->parm("LeftRightSymmmetry:coupHtaue");
448 0 : yukawa[2] = settingsPtr->parm("LeftRightSymmmetry:coupHtaumu");
449 0 : yukawa[3] = settingsPtr->parm("LeftRightSymmmetry:coupHtautau");
450 : }
451 :
452 : // Secondary open width fractions.
453 0 : openFracPos = particleDataPtr->resOpenFrac( idHLR);
454 0 : openFracNeg = particleDataPtr->resOpenFrac(-idHLR);
455 :
456 0 : }
457 :
458 : //--------------------------------------------------------------------------
459 :
460 : // Evaluate sigmaHat(sHat), including incoming flavour dependence.
461 :
462 : double Sigma2lgm2Hchgchgl::sigmaHat() {
463 :
464 : // Initial state must consist of a lepton and a photon.
465 0 : int idIn = (id2 == 22) ? id1 : id2;
466 0 : int idInAbs = abs(idIn);
467 0 : if (idInAbs != 11 && idInAbs != 13 && idInAbs != 15) return 0.;
468 :
469 : // Incoming squared lepton mass.
470 0 : double s1 = pow2( particleDataPtr->m0(idInAbs) );
471 :
472 : // Kinematical expressions.
473 0 : double smm1 = 8. * (sH + tH - s3) * (sH + tH - 2. * s3 - s1 - s4)
474 0 : / pow2(uH - s3);
475 0 : double smm2 = 2. * ( (2. * s3 - 3. * s1) * s4 + (s1 - 2. * s4) * tH
476 0 : - (tH - s4) * sH ) / pow2(tH - s4);
477 0 : double smm3 = 2. * ( (2. * s3 - 3. * s4 + tH) * s1
478 0 : - (2. * s1 - s4 + tH) * sH ) / pow2(sH - s1);
479 0 : double smm12 = 4. * ( (2. * s1 - s4 - 2. * s3 + tH) * sH
480 0 : + (tH - 3. * s3 - 3. * s4) * tH + (2. * s3 - 2. * s1
481 0 : + 3. * s4) * s3 ) / ( (uH - s3) * (tH - s4) );
482 0 : double smm13 = -4. * ( (tH + s1 - 2. * s4) * tH - (s3 + 3. * s1 - 2. * s4)
483 0 : * s3 + (s3 + 3. * s1 + tH) * sH - pow2(tH - s3 + sH) )
484 0 : / ( (uH - s3) * (sH - s1) );
485 0 : double smm23 = -4. * ( (s1 - s4 + s3) * tH - s3*s3 + s3 * (s1 + s4)
486 0 : - 3. * s1 * s4 - (s1 - s4 - s3 + tH) * sH)
487 0 : / ( (sH - s1) * (tH - s4) );
488 0 : double sigma = alpEM * pow2(sH / (sH - s1) ) * (smm1 + smm2 + smm3
489 0 : + smm12 + smm13 + smm23) / (4. * sH2);
490 :
491 : // Lepton Yukawa and secondary widths.
492 0 : sigma *= pow2(yukawa[(idInAbs-9)/2]);
493 0 : sigma *= (idIn < 0) ? openFracPos : openFracNeg;
494 :
495 : // Answer.
496 : return sigma;
497 :
498 0 : }
499 :
500 : //--------------------------------------------------------------------------
501 :
502 : // Select identity, colour and anticolour.
503 :
504 : void Sigma2lgm2Hchgchgl::setIdColAcol() {
505 :
506 : // Sign of outgoing H_L/R.
507 0 : int idIn = (id2 == 22) ? id1 : id2;
508 0 : int idSgn = (idIn < 0) ? idHLR : -idHLR;
509 0 : int idOut = (idIn < 0) ? idLep : -idLep;
510 0 : setId( id1, id2, idSgn, idOut);
511 :
512 : // tHat is defined between incoming lepton and outgoing Higgs.
513 0 : if (id1 == 22) swapTU = true;
514 :
515 : // No colours whatsoever.
516 0 : setColAcol( 0, 0, 0, 0, 0, 0, 0, 0);
517 :
518 0 : }
519 :
520 : //--------------------------------------------------------------------------
521 :
522 : // Evaluate weight for H_L/R sequential decay angles.
523 :
524 : double Sigma2lgm2Hchgchgl::weightDecay( Event& process, int iResBeg,
525 : int iResEnd) {
526 :
527 : // Identity of mother of decaying reseonance(s).
528 0 : int idMother = process[process[iResBeg].mother1()].idAbs();
529 :
530 : // For top decay hand over to standard routine.
531 0 : if (idMother == 6)
532 0 : return weightTopDecay( process, iResBeg, iResEnd);
533 :
534 : // Else isotropic decay.
535 0 : return 1.;
536 :
537 0 : }
538 :
539 : //==========================================================================
540 :
541 : // Sigma3ff2HchgchgfftWW class.
542 : // Cross section for f_1 f_2 -> H_(L/R)^++-- f_3 f_4 (W+- W+- fusion).
543 :
544 : //--------------------------------------------------------------------------
545 :
546 : // Initialize process.
547 :
548 : void Sigma3ff2HchgchgfftWW::initProc() {
549 :
550 : // Set process properties: H_L^++-- or H_R^++--.
551 0 : if (leftRight == 1) {
552 0 : idHLR = 9900041;
553 0 : codeSave = 3125;
554 0 : nameSave = "f_1 f_2 -> H_L^++-- f_3 f_4 (W+- W+- fusion)";
555 0 : } else {
556 0 : idHLR = 9900042;
557 0 : codeSave = 3145;
558 0 : nameSave = "f_1 f_2 -> H_R^++-- f_3 f_4 (W+- W+- fusion)";
559 : }
560 :
561 : // Common fixed mass and coupling factor.
562 0 : double mW = particleDataPtr->m0(24);
563 0 : double mWR = particleDataPtr->m0(9900024);
564 0 : mWS = (leftRight == 1) ? pow2(mW) : pow2(mWR);
565 0 : double gL = settingsPtr->parm("LeftRightSymmmetry:gL");
566 0 : double gR = settingsPtr->parm("LeftRightSymmmetry:gR");
567 0 : double vL = settingsPtr->parm("LeftRightSymmmetry:vL");
568 0 : prefac = (leftRight == 1) ? pow2(pow4(gL) * vL)
569 0 : : 2. * pow2(pow3(gR) * mWR);
570 : // Secondary open width fractions.
571 0 : openFracPos = particleDataPtr->resOpenFrac( idHLR);
572 0 : openFracNeg = particleDataPtr->resOpenFrac(-idHLR);
573 :
574 0 : }
575 :
576 : //--------------------------------------------------------------------------
577 :
578 : // Evaluate sigmaHat(sHat), part independent of incoming flavour.
579 :
580 : void Sigma3ff2HchgchgfftWW::sigmaKin() {
581 :
582 : // Required four-vector products.
583 0 : double pp12 = 0.5 * sH;
584 0 : double pp14 = 0.5 * mH * p4cm.pNeg();
585 0 : double pp15 = 0.5 * mH * p5cm.pNeg();
586 0 : double pp24 = 0.5 * mH * p4cm.pPos();
587 0 : double pp25 = 0.5 * mH * p5cm.pPos();
588 0 : double pp45 = p4cm * p5cm;
589 :
590 : // Cross section: kinematics part. Combine with couplings.
591 0 : double propT = 1. / ( (2. * pp14 + mWS) * (2. * pp25 + mWS) );
592 0 : double propU = 1. / ( (2. * pp24 + mWS) * (2. * pp15 + mWS) );
593 0 : sigma0TU = prefac * pp12 * pp45 * pow2(propT + propU);
594 0 : sigma0T = prefac * pp12 * pp45 * 2. * pow2(propT);
595 :
596 0 : }
597 :
598 : //--------------------------------------------------------------------------
599 :
600 : // Evaluate sigmaHat(sHat), including incoming flavour dependence.
601 :
602 : double Sigma3ff2HchgchgfftWW::sigmaHat() {
603 :
604 : // Do not allow creation of righthanded neutrinos for H_R.
605 0 : int id1Abs = abs(id1);
606 0 : int id2Abs = abs(id2);
607 0 : if ( leftRight == 2 && (id1Abs > 10 || id2Abs > 10) ) return 0.;
608 :
609 : // Many flavour combinations not possible because of charge.
610 0 : int chg1 = (( id1Abs%2 == 0 && id1 > 0)
611 0 : || (id1Abs%2 == 1 && id1 < 0) ) ? 1 : -1;
612 0 : int chg2 = (( id2Abs%2 == 0 && id2 > 0)
613 0 : || (id2Abs%2 == 1 && id2 < 0) ) ? 1 : -1;
614 0 : if (abs(chg1 + chg2) != 2) return 0.;
615 :
616 : // Basic cross section. CKM factors for final states.
617 0 : double sigma = (id2 == id1 && id1Abs > 10) ? sigma0TU : sigma0T;
618 0 : sigma *= couplingsPtr->V2CKMsum(id1Abs)
619 0 : * couplingsPtr->V2CKMsum(id2Abs);
620 :
621 : // Secondary width for H0.
622 0 : sigma *= (chg1 + chg2 == 2) ? openFracPos : openFracNeg;
623 :
624 : // Spin-state extra factor 2 per incoming neutrino.
625 0 : if (id1Abs == 12 || id1Abs == 14 || id1Abs == 16) sigma *= 2.;
626 0 : if (id2Abs == 12 || id2Abs == 14 || id2Abs == 16) sigma *= 2.;
627 :
628 : // Answer.
629 : return sigma;
630 :
631 0 : }
632 :
633 : //--------------------------------------------------------------------------
634 :
635 : // Select identity, colour and anticolour.
636 :
637 : void Sigma3ff2HchgchgfftWW::setIdColAcol() {
638 :
639 : // Pick out-flavours by relative CKM weights.
640 0 : int id1Abs = abs(id1);
641 0 : int id2Abs = abs(id2);
642 0 : id4 = couplingsPtr->V2CKMpick(id1);
643 0 : id5 = couplingsPtr->V2CKMpick(id2);
644 :
645 : // Find charge of Higgs .
646 0 : id3 = (( id1Abs%2 == 0 && id1 > 0) || (id1Abs%2 == 1 && id1 < 0) )
647 0 : ? idHLR : -idHLR;
648 0 : setId( id1, id2, id3, id4, id5);
649 :
650 : // Colour flow topologies. Swap when antiquarks.
651 0 : if (id1Abs < 9 && id2Abs < 9 && id1*id2 > 0)
652 0 : setColAcol( 1, 0, 2, 0, 0, 0, 1, 0, 2, 0);
653 0 : else if (id1Abs < 9 && id2Abs < 9)
654 0 : setColAcol( 1, 0, 0, 2, 0, 0, 1, 0, 0, 2);
655 0 : else if (id1Abs < 9) setColAcol( 1, 0, 0, 0, 0, 0, 1, 0, 0, 0);
656 0 : else if (id2Abs < 9) setColAcol( 0, 0, 1, 0, 0, 0, 0, 0, 1, 0);
657 0 : else setColAcol( 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
658 0 : if ( (id1Abs < 9 && id1 < 0) || (id1Abs > 10 && id2 < 0) )
659 0 : swapColAcol();
660 :
661 0 : }
662 :
663 : //--------------------------------------------------------------------------
664 :
665 : // Evaluate weight for decay angles.
666 :
667 : double Sigma3ff2HchgchgfftWW::weightDecay( Event& process, int iResBeg,
668 : int iResEnd) {
669 :
670 : // Identity of mother of decaying reseonance(s).
671 0 : int idMother = process[process[iResBeg].mother1()].idAbs();
672 :
673 : // For top decay hand over to standard routine.
674 0 : if (idMother == 6)
675 0 : return weightTopDecay( process, iResBeg, iResEnd);
676 :
677 : // Else done.
678 0 : return 1.;
679 :
680 0 : }
681 :
682 : //==========================================================================
683 :
684 : // Sigma2ffbar2HchgchgHchgchg class.
685 : // Cross section for f fbar -> H_(L/R)^++ H_(L/R)^-- (doubly charged Higgs).
686 :
687 : //--------------------------------------------------------------------------
688 :
689 : // Initialize process.
690 :
691 : void Sigma2ffbar2HchgchgHchgchg::initProc() {
692 :
693 : // Set process properties: H_L^++ H_L^-- or H_R^++ H_R^--.
694 0 : if (leftRight == 1) {
695 0 : idHLR = 9900041;
696 0 : codeSave = 3126;
697 0 : nameSave = "f fbar -> H_L^++ H_L^--";
698 0 : } else {
699 0 : idHLR = 9900042;
700 0 : codeSave = 3146;
701 0 : nameSave = "f fbar -> H_R^++ H_R^--";
702 : }
703 :
704 : // Read in Yukawa matrix for couplings to a lepton pair.
705 0 : yukawa[1][1] = settingsPtr->parm("LeftRightSymmmetry:coupHee");
706 0 : yukawa[2][1] = settingsPtr->parm("LeftRightSymmmetry:coupHmue");
707 0 : yukawa[2][2] = settingsPtr->parm("LeftRightSymmmetry:coupHmumu");
708 0 : yukawa[3][1] = settingsPtr->parm("LeftRightSymmmetry:coupHtaue");
709 0 : yukawa[3][2] = settingsPtr->parm("LeftRightSymmmetry:coupHtaumu");
710 0 : yukawa[3][3] = settingsPtr->parm("LeftRightSymmmetry:coupHtautau");
711 :
712 : // Electroweak parameters.
713 0 : mRes = particleDataPtr->m0(23);
714 0 : GammaRes = particleDataPtr->mWidth(23);
715 0 : m2Res = mRes*mRes;
716 0 : GamMRat = GammaRes / mRes;
717 0 : sin2tW = couplingsPtr->sin2thetaW();
718 0 : preFac = (1. - 2. * sin2tW) / ( 8. * sin2tW * (1. - sin2tW) );
719 :
720 : // Open fraction from secondary widths.
721 0 : openFrac = particleDataPtr->resOpenFrac( idHLR, -idHLR);
722 :
723 0 : }
724 :
725 : //--------------------------------------------------------------------------
726 :
727 : // Evaluate sigmaHat(sHat), including incoming flavour dependence.
728 :
729 : double Sigma2ffbar2HchgchgHchgchg::sigmaHat() {
730 :
731 : // Electroweak couplings to gamma^*/Z^0.
732 0 : int idAbs = abs(id1);
733 0 : double ei = couplingsPtr->ef(idAbs);
734 0 : double vi = couplingsPtr->vf(idAbs);
735 0 : double ai = couplingsPtr->af(idAbs);
736 :
737 : // Part via gamma^*/Z^0 propagator. No Z^0 coupling to H_R.
738 0 : double resProp = 1. / ( pow2(sH - m2Res) + pow2(sH * GamMRat) );
739 0 : double sigma = 8. * pow2(alpEM) * ei*ei / sH2;
740 0 : if (leftRight == 1) sigma += 8. * pow2(alpEM)
741 0 : * (2. * ei * vi * preFac * (sH - m2Res) * resProp / sH
742 0 : + (vi * vi + ai * ai) * pow2(preFac) * resProp);
743 :
744 : // Part via t-channel lepton + interference; sum over possibilities.
745 0 : if (idAbs == 11 || idAbs == 13 || idAbs == 15) {
746 0 : double yuk2Sum;
747 0 : if (idAbs == 11) yuk2Sum
748 0 : = pow2(yukawa[1][1]) + pow2(yukawa[2][1]) + pow2(yukawa[3][1]);
749 0 : else if (idAbs == 13) yuk2Sum
750 0 : = pow2(yukawa[2][1]) + pow2(yukawa[2][2]) + pow2(yukawa[3][2]);
751 : else yuk2Sum
752 0 : = pow2(yukawa[3][1]) + pow2(yukawa[3][2]) + pow2(yukawa[3][3]);
753 0 : yuk2Sum /= 4. * M_PI;
754 0 : sigma += 8. * alpEM * ei * yuk2Sum / (sH * tH)
755 0 : + 4. * pow2(yuk2Sum) / tH2;
756 0 : if (leftRight == 1) sigma += 8. * alpEM * (vi + ai) * yuk2Sum
757 0 : * preFac * (sH - m2Res) * resProp / tH;
758 0 : }
759 :
760 : // Common kinematical factor. Colour factor.
761 0 : sigma *= M_PI * (tH * uH - s3 * s4) / sH2;
762 0 : if (idAbs < 9) sigma /= 3.;
763 :
764 : // Answer.
765 0 : return sigma;
766 :
767 : }
768 :
769 : //--------------------------------------------------------------------------
770 :
771 : // Select identity, colour and anticolour.
772 :
773 : void Sigma2ffbar2HchgchgHchgchg::setIdColAcol() {
774 :
775 : // Outgoing flavours trivial.
776 0 : setId( id1, id2, idHLR, -idHLR);
777 :
778 : // tHat is defined between incoming fermion and outgoing H--.
779 0 : if (id1 > 0) swapTU = true;
780 :
781 : // No colours at all or one flow topology. Swap if first is antiquark.
782 0 : if (abs(id1) < 9) setColAcol( 1, 0, 0, 1, 0, 0, 0, 0);
783 0 : else setColAcol( 0, 0, 0, 0, 0, 0, 0, 0);
784 0 : if (id1 < 0) swapColAcol();
785 :
786 0 : }
787 :
788 : //--------------------------------------------------------------------------
789 :
790 : // Evaluate weight for H_L/R sequential decay angles.
791 :
792 : double Sigma2ffbar2HchgchgHchgchg::weightDecay( Event& process,
793 : int iResBeg, int iResEnd) {
794 :
795 : // Identity of mother of decaying reseonance(s).
796 0 : int idMother = process[process[iResBeg].mother1()].idAbs();
797 :
798 : // For top decay hand over to standard routine.
799 0 : if (idMother == 6)
800 0 : return weightTopDecay( process, iResBeg, iResEnd);
801 :
802 : // Else isotropic decay.
803 0 : return 1.;
804 :
805 0 : }
806 :
807 : //==========================================================================
808 :
809 : } // end namespace Pythia8
|