StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
WeakShowerMEs.cc
1 // WeakShowerMEs.cc is a part of the PYTHIA event generator.
2 // Copyright (C) 2018 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, 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 // WeakShowerMEs class.
8 
9 #include "Pythia8/WeakShowerMEs.h"
10 
11 namespace Pythia8 {
12 
13 //==========================================================================
14 
15 // The WeakShowerMEs class.
16 // The 2 -> 2 MEs contain the correct kinematics, but for some of
17 // the 2 -> 3 MEs some couplings have been switched off.
18 // This class is used for ME corrections in the weak shower and for merging
19 // whenever weak reclusterings are allowed.
20 
21 //--------------------------------------------------------------------------
22 
23 // Calculate the 2->2 ME qg -> qg, omitting an overall factor of g_s^4 / 9.
24 
25 double WeakShowerMEs::getMEqg2qg(double sH, double tH, double uH) {
26 
27  double sH2 = sH * sH;
28  double tH2 = tH * tH;
29  double uH2 = uH * uH;
30  return (sH2 + uH2) * (9. / tH2 - 4. / (sH * uH));
31 
32 }
33 
34 //--------------------------------------------------------------------------
35 
36 // Calculate the 2->2 ME qq' -> qq', omitting an overall factor of g_s^4 / 9s.
37 
38 double WeakShowerMEs::getMEqq2qq(double sH, double tH, double uH,
39  bool sameID) {
40 
41  double sH2 = sH * sH;
42  double tH2 = tH * tH;
43  double uH2 = uH * uH;
44  if (sameID) return 2. * ( (sH2 + uH2) / tH2 + (sH2 + tH2) / uH2
45  - 2. * sH2 / (3. * tH * uH) );
46  else return 4. * (sH2 + uH2) / tH2;
47 
48 }
49 
50 //--------------------------------------------------------------------------
51 
52 // Calculate the 2->2 ME gg -> gg, omitting an overall factor of g_s^4 / 9.
53 
54 double WeakShowerMEs::getMEgg2gg(double sH, double tH, double uH) {
55 
56  double sH2 = sH * sH;
57  double tH2 = tH * tH;
58  double uH2 = uH * uH;
59  return (81. / 8.) * ( (tH2 + uH2) / sH2 + (sH2 + uH2) / tH2
60  + (sH2 + tH2) / uH2 + 3. );
61 
62 }
63 
64 //--------------------------------------------------------------------------
65 
66 // Calculate the 2->2 ME gg -> qqbar, omitting an overall factor of g_s^4 / 9.
67 
68 double WeakShowerMEs::getMEgg2qqbar(double sH, double tH, double uH) {
69 
70  double sH2 = sH * sH;
71  double tH2 = tH * tH;
72  double uH2 = uH * uH;
73  return (tH2 + uH2) * ( 3. / (2. * tH * uH) - 27. / (8. * sH2) );
74 
75 }
76 
77 //--------------------------------------------------------------------------
78 
79 // Calculate the 2->2 ME qqbar -> gg, omitting an overall factor of g_s^4 / 9.
80 
81 double WeakShowerMEs::getMEqqbar2gg(double sH, double tH, double uH) {
82 
83  double sH2 = sH * sH;
84  double tH2 = tH * tH;
85  double uH2 = uH * uH;
86  return (tH2 + uH2) * ( 16. / (3. * tH * uH) - 12. / sH2 );
87 
88 }
89 
90 //--------------------------------------------------------------------------
91 
92 // Calculate the 2->2 ME qqbar -> qqbar, omitting an overall factor of
93 // g_s^4 / 9.
94 
95 double WeakShowerMEs::getMEqqbar2qqbar(double sH, double tH, double uH,
96  bool sameID) {
97 
98  double sH2 = sH * sH;
99  double tH2 = tH * tH;
100  double uH2 = uH * uH;
101  if (sameID) return 4. * (tH2 + uH2) / sH2 - (8./3.) * uH2 / (sH * tH)
102  + 4. * (sH2 + uH2) / tH2;
103  else return 4. * (tH2 + uH2) / sH2;
104 
105 }
106 
107 //--------------------------------------------------------------------------
108 
109 // Calculate the 2->3 ME ug -> ugZ, omitting an overall factor of
110 // \frac{(g_s^4 * 4\pi \alpha_{em} * (lU^2 + rU^2))}
111 // {9 * \cos^2(\theta_W) \sin^2(\theta_W)}.
112 // p1 = incoming quark, p2 = incoming gluon, p3 = outgoing gluon,
113 // p4 = outgoing Z, p5 = outgoing quark.
114 
115 double WeakShowerMEs::getMEqg2qgZ(Vec4 p1, Vec4 p2, Vec4 p3, Vec4 p4,
116  Vec4 p5) {
117 
118  double p12 = p1*p2;
119  double p13 = p1*p3;
120  double p14 = p1*p4;
121  double p23 = p2*p3;
122  double p24 = p2*p4;
123  double p44 = p4*p4;
124  double p12P = (p1 + p2).m2Calc();
125  double p13M = (p1 - p3).m2Calc();
126  double p14M = (p1 - p4).m2Calc();
127  double p23M = (p2 - p3).m2Calc();
128  double p25M = (p2 - p5).m2Calc();
129  double p35P = (p3 + p5).m2Calc();
130  double p45P = (p4 + p5).m2Calc();
131  double p44S = p44*p44;
132 
133  double dia1 = (-4*p12*(p12* (p44 - 2*p24) -p23 * (p44 - 2*p24) +2*p13*p24))
134  / pow2(p12P * p45P);
135 
136  double dia2 = -(((p12 - p13 - p23) * (p12 * (p44 - 2*p14)
137  - p13 * (p44 - 2*p14) + 2*p14*p23)) / (p12P * p13M * pow2(p45P)));
138 
139  double dia3 = (8*(p44 - 2*p12) * p12 * (p12 - p23 - p24))
140  / (pow2(p12P) * p35P * p45P);
141 
142  double dia4 = -(8*pow3(p12) + 2*pow2(p12) * (p44 - 6*p13 - 6*p14 - 6*p23
143  - 4*p24) + 4*pow2(p13) * (p44 - 2*p14 - p24) - 2*p14*p23
144  * (-p44 + 2*p14 + 2*p23 + 4*p24) + 2*p13 * (-4*pow2(p14)
145  + 2*p14 * (p44 - 3*p23 - 3*p24) + 2*p23 * (p44 - p24) + p44*p24)
146  + p12* (p44S + 4*pow2(p13) + 4*pow2(p14) + 4*pow2(p23) - 4*p14
147  * (p44 - 4*p23 - 3*p24) - 4*p13* (p44 - 4*p14 - 4*p23 - 3*p24)
148  - 2*p44*p24 + 4*p23*p24)) / (4. * p12P * p13M * p25M * p45P);
149 
150  double dia5 = (-9*(2*pow3(p12) + (p13 + p23) * (-(p14* p23)
151  + p13 * (p44 - 2*p14 - p24)) + pow2(p12)* (2*p44 - 4*p13 - 3*p14
152  - 4*p23 - 3*p24) + p12 * (2*pow2(p13) + p23 * (-p44 + 4*p14
153  + 2*p23 + 3*p24) + p13 * (-2*p44 + 5*p14 + 4*p23 + 4*p24))))
154  / (p12P * p23M* pow2(p45P));
155 
156  double dia6 = (-4*(2*pow3(p12) + pow2(p12)* (7*p44 - 4*p13 - 6*p14
157  - 2*p23 - 4*p24) + p12 * (2*pow2(p13) + 4*pow2(p14) - 2*p44*p23
158  - 3*p44*p24 + 2*p23* p24 + 2*pow2(p24) + p14 * (-8*p44 + 6*p23 + 2*p24)
159  + p13* (-7*p44 + 10*p14 + 6*p23 + 6*p24)) + 2*(pow2(p13) * (p44
160  - 2*p14 - p24) + p14* (p14* (p44 - 2*p23) + p23 * (p44 - p24) + p44*p24)
161  + p13 * (-2*pow2(p14) + p23 * (p44 - 2*p24) + p14* (2*p44 - 3*p23
162  - 2*p24) + (p44 - p24)* p24)))) / (p12P * p14M * p35P * p45P);
163 
164  double dia7 = ((p12 - p13 - p14) * (p44 - 2*p23) * (p12 - p13 - p14
165  - p23 - p24)) / (p12P * p14M * p25M * p45P);
166 
167  double dia8 = (9*(2*pow2((p1*p2))* (p44 - 10*p14 - 16*p23 + 4*p24)
168  + p12* (p44S + 20*pow2(p14) - 6*p44*p23 + 16*pow2(p23)
169  + p13* (-18*p44 + 28*p14 + 32*p23 - 8*p24) + 2*p44*p24
170  + 16*p23* p24 - 8*pow2(p24) + p14* (-8*p44 + 52*p23 + 4*p24))
171  + 2*(4*pow2(p13)* (p44 - p14 - 2*p23) - p14*p23* (-3*p44 + 10*p14
172  + 8*p23 + 8*p24) + p13* (-4*pow2(p14) - 8*pow2(p23) + 4*p23* (p44 - p24)
173  + 2*p14* (2*p44 - 8*p23 + p24) + p24* (p44 + 4*p24)))))
174  / (8. * p12P * p23M * p14M * p45P);
175 
176  double dia9 = (-4*p13* (2*pow2(p12) + p12* (p44 - 4*p13 - 2*p14
177  - 4*p23 - 2*p24) + 2*(p13 + p23)* (p13 + p14 + p23 + p24)))
178  / (pow2(p13M) * pow2(p45P));
179 
180  double dia10 = (4*pow3(p12) - 4*pow2(p12) * (3*p13 + 2*p14 + 2*p23 + p24)
181  + p12* (p44S + 8*pow2(p13) + 4*pow2(p14) + 4*pow2(p23)
182  - 4*p14* (p44 - 3*p23 - 3*p24) - 2*p44*p24 + 4*p23*p24
183  + p13* (-6*p44 + 20*p14 + 4*p23 + 12*p24)) + 2*(p14*p23* (p44
184  - 2*p14 - 2*p23 - 4*p24) + 2*pow2(p13)* (p44 - 2*p14 - 2*p24)
185  + p13* (-4*pow2(p14) + 2*p14* (p44 - 2*p23 - 3*p24)
186  + 2*p23* (p44 - p24) + p44*p24))) / (4. * p12P * p13M * p35P* p45P);
187 
188  double dia11 = (4*p13* (p44 + 2*p13)* (p44 + 2*p12 - 2*p14 - 2*p24))
189  / (pow2(p13M)* p25M * p45P);
190 
191  double dia12 = (-9*(-2*pow3(p12) + 2*pow3(p13) - 2*p14* pow2(p23)
192  + pow2(p12)* (p44 + 6*p13 - 2*p14 + 4*p23 + 2*p24) + p13*p23
193  * (-p44 - 2*p14 + 2*p23 + 6*p24) + pow2(p13)* (p44 + 4*p23 + 6*p24)
194  - p12* (6*pow2(p13) + p23* (p44 - 4*p14 + 2*p23 + 2*p24)
195  + p13* (-2*p14 + 8*(p23 + p24))))) / (2.* p13M * p23M * pow2(p45P));
196 
197  double dia13 = -((p12 - p13 - p14)* (p44 - 2*p23)* (p44 - 2*p24))
198  / (2. * p13M * p14M * p35P * p45P);
199 
200  double dia14 = (2*(-8*pow2(p13)* (p44 - p14) - 2*pow2(p12)
201  * (p44 - 2*p14 + 2*p23 - 2*p24) + 2*p14* (-p44S - 2*pow2(p23)
202  + p23 * (p44 + 2*p14 - 2*p24) + 2*p44*p24) + p12* (p44S
203  - 4*pow2(p14) + 4*pow2(p23) + 2*p13* (3*p44 - 6*p14 - 2*p23 - 2*p24)
204  - 4*pow2(p24)) + p13* (-2*p44S + 8*pow2(p14) + 2*(p44 + 2*p23)* p24
205  + 4*pow2(p24) + p14* (-8*p44 + 8*p23 + 4*p24))))
206  / (p13M * p14M * p25M * p45P);
207 
208  double dia15 = (-9*(8*pow3(p12) + 2*pow2(p12)* (p44 - 8*p13 - 2*p14
209  - 4*p23 - 8*p24) + p12* (3*p44S + 8*pow2(p13) - 4*pow2(p14)
210  - 6*p44*p23 - 10*p44*p24 + 24*p23* p24 + 8*pow2(p24) + 4*p14
211  * (p23 + 3*p24) + 2*p13 * (3*p44 - 2*p14 + 8*p23 + 12*p24))
212  + 2*(p14*p23 * (p44 + 2*p14 - 8*p24) + 4*pow2(p13)* (p14
213  + p23 - p24) + p13* (4*pow2(p14) - 4*pow2(p23)
214  + 4*p23* (p44 - 4*p24) - 10*p14* p24 + (3*p44 - 4*p24)* p24))))
215  / (8. * p13M * p23M * p14M * p45P);
216 
217  double dia16 = (4*p12*p23* (p44 + 2*p12 - 2*p14 - 2*p24))
218  / (pow2(p12P) * pow2(p35P));
219 
220  double dia17 = ((p12 - p13 - p14)* (p44 - 2*p24)
221  * (p12 - p13 - p14 - p23 - p24)) / (p12P * p13M * p25M * p35P);
222 
223  double dia18 = (9*(-20*pow3(p12) + 4*pow2(p12)* (4*p44 + p13 + 2*p14
224  + p23) + 2*p14*p23* (p44 + 2*p14 - 4*p23 - 2*p24) + 8*pow2(p13)
225  * (p44 - p14 - p24) + p13* (-8*pow2(p14) + 8*p44*p23
226  + 4*p14* (2*p44 - 2*p23 - 5*p24) + 6*p44*p24 - 4*pow2(p24))
227  + p12* (p44S - 4*pow2(p14) - 14*p44*p23 + 8*pow2(p23)
228  - 8*p14* (p44 - p23 - 2*p24) - 12*p44*p24 + 12*p23* p24
229  + 4*pow2(p24) + p13 * (-22*p44 + 20*p14 - 8*p23 + 24*p24))))
230  / (8. * p12P * p23M * p35P * p45P);
231 
232  double dia19 = (-8*p13* (p44 + 2*p12 - 2*p14 - 2*p24)*(p12 - p14 - p24))
233  / (p12P * p14M * pow2(p35P));
234 
235  double dia20 = (-8*pow3(p12) - p12* (4*pow2(p13) + 4*pow2(p14)
236  - 4*p14* (p44 - 3*p23 - 4*p24) - 4*p13* (p44 - 4*p14 - 3*p23 - 4*p24)
237  + (p44 - 2*p24)* (p44 - 2*p23 - 2*p24)) - 2*pow2(p12)
238  * (p44 - 6*p13 - 6*p14 - 4*p23 - 6*p24) + 2*p14*p23
239  * (-p44 + 2*p14 + 2*p24) + pow2(p13) * (-4*p44 + 8*p14 + 4*p24)
240  + p13* (8*pow2(p14) - 4*p14* (p44 - 3*p23 - 3*p24) - 2*(p44 - 2*p24)
241  * (2*p23 + p24))) / (4. * p12P * p14M * p25M* p35P);
242 
243  double dia21 = (-9*(14*pow3(p12) + pow2(p12)* (7*p44 - 10*p13 - 30*p14
244  - 10*p23 - 22*p24) + 2*p12 * (2*pow2(p13) + 12*pow2(p14) - p44*p23
245  - 3*p44*p24 + 2*p23* p24 + 4*pow2(p24) + p13* (-4*p44 + 10*p14
246  + 6*p23 + 5*p24) + p14* (-5*p44 + 7*p23 + 16*p24)) + 2*(2*pow2(p13)
247  * (p44 - 2*p14 - p24) + p13* (-4*pow2(p14) + 2*p23* (p44 - 2*p24)
248  + 2*p14* (p44 - 3*p23 - 2*p24) + p44*p24) + p14* (-4*pow2(p14)
249  + 2*p14* (p44 - p23 - 3*p24) + p23* (p44 - 2*p24)
250  + 2*(p44 - p24)* p24)))) / (4. * p12P * p23M * p14M* p35P);
251 
252  double dia22 = (-8*p13*p23* (-p12 + p23 + p24))
253  / (pow2(p13M)* pow2(p25M));
254 
255  double dia23 = (-9*(2*pow3(p12) - pow2(p12)* (p44 + 10*p13 + 6*p14)
256  + p12* (p44S + 8*pow2(p13) + 4*pow2(p14) + 3*p44*p23 - 2*pow2(p23)
257  - p44*p24 - 4*p23 * p24 - 2*pow2(p24) + 2*p13* (6*p14 + 7*p23 + 2*p24)
258  + p14* (-2*p44 + 4*p23 + 6*p24)) - 2*(4*pow3(p13) + 2*pow2(p13)
259  * (p44 + 2*p14 + p23 - p24) + p14*p23* (2*p14 - p23 + p24)
260  + p13 * (p44S + 4*pow2(p14) + (-2*p44 + p23)* p24 - pow2(p24)
261  + p14* (-2*p44 + 4*p23 + 2*p24))))) / (4. * p13M * p23M * p25M * p45P);
262 
263  double dia24 = (-2*pow2(p12)* (p44 + p13 + p23 - p24) + p12* (2*pow2(p13)
264  - p13* (p44 - 6*p14) + (p44 + 2*p14 + 2*p23 - 2*p24)* (p23 + p24))
265  + 2*(pow2(p13)* (p44 - 2*p14 - p24) - p14*p23 * (p23 + p24)
266  + p13* (-2*pow2(p14) + p14 * (p44 - p23 - 2*p24) + p24* (p23 + p24))))
267  / (2. * p13M * p14M * p25M* p35P);
268 
269  double dia25 = (8*p12* (p44 + 2*p12 - 2*p23 - 2*p24)*(p12 - p23 - p24))
270  / (p13M * p14M * pow2(p25M));
271 
272  double dia26 = (9*(4*pow3(p12) + 2*pow2(p13)* (p44 - 2*p14 - 3*p24)
273  - 2*pow2(p12)* (p13 + 6*p23 + 2*p24) + p12* (p44S + 6*pow2(p13)
274  - 2*p44*p14 + 4*pow2(p14) - 2*p44*p23 + 8*pow2(p23) + 3*p13
275  * (p44 + 2*p14 - 2*p23 - 2*p24) - 2*p44*p24 + 8*p23*p24)
276  + p14* (-p44S + 2*p14* (p44 - 2*p23 - 2*p24) + 4*p23* p24
277  + 4*pow2(p24)) - p13 * (p44S + 4*pow2(p14) - 4*pow2(p23)
278  + 2*p23 * (p44 - 6*p24) + 2*p44*p24 - 8*pow2(p24)
279  + p14* (-4*p44 + 2*p23 + 8*p24)))) / (4. * p13M * p23M * p14M * p25M);
280 
281  double dia27 = (-9*(3*pow3(p12) - pow3(p13) + p14 * pow2(p23)
282  + pow2(p12)* (p44 - 7*p13 - 3*p14 - 6*p23 - 2*p24) - p13*p23
283  * (2*p44 + p23 - 2*p24) - pow2(p13)* (p14 + 2*p23 - 2*p24)
284  + p12 * (5*pow2(p13) + p13* (p44 + 4*p14 + 8*p23)
285  + p23 * (p44 + 2*p14 + 3*p23 + 2*p24)))) / (pow2(p23M) * pow2(p45P));
286 
287  double dia28 = (9*(2*pow3(p12) - pow2(p12)* (p44 + 4*p13 + 10*p14
288  + 8*p23 - 2*p24) + 2*pow2(p13)* (p44 - 2*p14 - p24) + p14* (p44S
289  + 2*p44*p23 + 2*pow2(p23) - 2*p14* (p44 + 6*p23) - 4*p44*p24)
290  + p12* (-p44S + 2*pow2(p13) + 8*pow2(p14) + p13 * (-5*p44 + 14*p14)
291  - 3*p44*p23 - 2*pow2(p23) + 2*p14* (p44 + 10*p23) + 6*p44*p24
292  + 2*p23* p24 - 4*pow2(p24)) + p13* (-4*pow2(p14) + pow2(p44 - 2*p24)
293  + 2*p23* (p44 + p24) + p14* (-6*p23 + 4*p24))))
294  / (4. * p23M * p14M * p35P * p45P);
295 
296  double dia29 = (-9*(6*pow3(p12) - 2*pow2(p13)* (p44 - 2*p14 + 3*p24)
297  - pow2(p12)* (p44 + 12*p13 + 2*p14 + 8*p23 + 10*p24) + p14* (p44S
298  - 2*p14* (p44 - 6*p23) - 6*p44*p23 - 2*pow2(p23) - 4*p44*p24)
299  + p13* (p44S + 4*pow2(p14) - 8*pow2(p23) - 2*p44* p24 - 4*pow2(p24)
300  - 2*p14* (2*p44 - 5*p23 + 4*p24) - 2*p23* (2*p44 + 5*p24))
301  + p12* (6* pow2(p13) - 4*pow2(p14) + 3*p44*p23 + 2*pow2(p23)
302  + 6*p23* p24 + 4*pow2(p24) + p14 * (6*p44 - 4*p23 + 4*p24)
303  + p13* (7*p44 - 2*p14 + 16*p23 + 16*p24))))
304  / (4. * p23M * p14M * p25M * p45P);
305 
306  double dia30 = (-9*(4*pow3(p12) + 2*pow2(p12)* (p44 - 4*p13 + 2*p14
307  - 6*p24) + p12* (p44S + 4*pow2(p13) - 8*pow2(p14) - 6*p44*p23
308  + 12*pow2(p23) - 6*p44 *p24 + 4*p23* p24 + 8*pow2(p24) + 2*p14
309  * (p44 - 4*p23 + 2*p24) + 8*p13* (p44 - p14 + 2*p23 + 2*p24))
310  - 2*(-2*p14* (p44 + 2*p14 - 3*p23)* p23 + pow2(p13)* (p44 - 2*p14
311  + 2*p24) + p13 * (-2*pow2(p14) + 8*pow2(p23) - 3*p23* (p44 - 2*p24)
312  + p24* (-p44 + 4*p24) + p14* (p44 + 6*p24)))))
313  / (2. * pow2(p23M) * p14M * p45P);
314 
315  double dia31 = (-2*(p13* (p44 - 2*p14) + p14* (p44 + 2*p12 - 2*p14
316  - 2*p23 - 2*p24))* (p44 + 2*p12 - 2*p14 - 2*p24))
317  / (pow2(p14M) * pow2(p35P));
318 
319  double dia32 = -((p44 - 2*p14)* (p12* (p44 - 2*p14) - p13* (p44 - 2*p14)
320  + 2*p14*p23)) / (2. * pow2(p14M) * p25M * p35P);
321 
322  double dia33 = (-9*(4*pow3(p14) + 2*pow2(p12)* (p44 + 2*p14)
323  - 4*pow2(p14)* (p44 + p23 - 3*p24) - 2*p44*p13* p24
324  + p14* (p44S + 2*p44*p23 + (-6*p44 + 4*p13)* p24 + 8*pow2(p24))
325  + p12 * (2*p13* (p44 - 2*p14) + p44*(3*p44 + 2*p23 - 2*p24)
326  - 2*p14* (3*p44 + 2*p23 + 6*p24)))) / (4. * p23M * pow2(p14M) * p35P);
327 
328  double dia34 = (-4*(p12 - p23 - p24)* (p44*p12-2*p14*p24))
329  / (pow2(p14M)* pow2(p25M));
330 
331  double dia35 = (-9*(2*pow2(p12)* (p44 + 2*p14) - 2*p13 * (p44 - 2*p14)
332  * (p44 - 2*p14 + 2*p23 + p24) + 2*p14* (4*pow2(p23) + 8*p23* p24
333  + p24* (-p44 + 2*p14 + 4*p24)) + p12* (2*p13* (p44 - 2*p14)
334  + p44*(p44 - 2*p23 - 2*p24) - 2*p14* (p44 + 6*p23 + 6*p24))))
335  / (4. * p23M * pow2(p14M) * p25M);
336 
337  double dia36 = (-9*(4*pow2(p12)* (p44 + 2*p14) - p13 * (p44 - 2*p14)
338  * (p44 - 2*p14 + 8*p23 + 4*p24) + p14* (-p44S - 4*pow2(p14)
339  + 16*pow2(p23) - 4*p23 * (p44 - 4*p24) - 8*p44*p24 + 16*pow2(p24)
340  + 4*p14* (p44 + 2*p23 + 4*p24)) + p12* (4*p13 * (p44 - 2*p14)
341  - 4*pow2(p14) + p44*(3*p44 + 4*p23 - 4*p24) - 4*p14* (p44 + 6*p23
342  + 6*p24)))) / (4. * pow2(p23M) * pow2(p14M));
343 
344  return dia1 + dia2 + dia3 + dia4 + dia5 + dia6 + dia7 + dia8 + dia9
345  + dia10 + dia11 + dia12 + dia13 + dia14 + dia15 + dia16 + dia17 + dia18
346  + dia19 + dia20 + dia21 + dia22 + dia23 + dia24 + dia25 + dia26 + dia27
347  + dia28 + dia29 + dia30 + dia31 + dia32 + dia33 + dia34 + dia35 + dia36;
348 
349 }
350 
351 //--------------------------------------------------------------------------
352 
353 // Calculate the 2->3 ME ud -> udZ, with the coupling between Z and d
354 // set to zero, and omitting an overall factor of
355 // \frac{(g_s^4 * 4\pi \alpha_{em} * (lU^2 + rU^2))}
356 // {9 * \cos^2(\theta_W) \sin^2(\theta_W)}..
357 // p1 = incoming u, p2 = incoming d, p3 = outgoing Z,
358 // p4 = outgoing d, p5 = outgoing u.
359 
360 double WeakShowerMEs::getMEqq2qqZ(Vec4 p1, Vec4 p2, Vec4 p3, Vec4 p4,
361  Vec4 p5) {
362 
363  double p12 = p1*p2;
364  double p13 = p1*p3;
365  double p14 = p1*p4;
366  double p23 = p2*p3;
367  double p24 = p2*p4;
368  double p33 = p3*p3;
369  double p13M = (p1 - p3).m2Calc();
370  double p24M = (p2 - p4).m2Calc();
371  double p35P = (p3 + p5).m2Calc();
372  double p33S = p33*p33;
373 
374  double dia1 = (-4*(2*pow3(p12) + pow2(p12) * (p33 - 2*p13 - 4*p14
375  - 2*p23 - 4*p24) + p14* (2*p14* p23 - (p33 - 2*p23)* p24)
376  + p12 * (2*pow2(p14) + 2*p24* (p13 + p23 + p24)
377  + p14* (p33 + 2*p13 + 4*p24)))) / pow2(p24M* p35P);
378 
379  double dia2 = (-2*(4*pow3(p12) + 4*pow2(p12)
380  * (p33 - 2*p14 - 3*p23) + p12 * (p33S - 4*pow2(p13) + 4*pow2(p14)
381  - 6*p33*p23 + 8*pow2(p23) + 4*p13 * (p23 - p24) - 4*p33*p24 + 4*p23
382  * p24 + 4*pow2(p24) + 4*p14* (p33 + 4*p23 + 4*p24)) + 2*(-2*pow2(p14)
383  * p23 + p13* (p33 + 2*p13 - 2*p24) * p24 + p14* (-4*pow2(p23) + p23
384  * (p33 - 6*p13 - 6*p24) + 2*(p33 - p13 - 2*p24)* p24))))
385  / (p13M * pow2(p24M) * p35P);
386 
387  double dia3 = (-2*(2*pow2(p12)* (p33 + 2*p13)
388  - 2*p33*p14* (p23 + p24) + 4*pow2(p13)* (2*p23 + p24) + p12
389  * (-4*pow2(p13) + p33*(p33 + 2*p14 - 2*p23) - 4*p13* (p14 + 3*p23
390  + 2*p24)) + p13* (8*pow2(p23) + 2*p24 * (-p33 + 2*p14 + 2*p24) + p23
391  * (-4*p33 + 4*p14 + 8*p24)))) / pow2(p13M * p24M);
392 
393  return dia1 + dia2 + dia3;
394 
395 }
396 
397 //==========================================================================
398 
399 } // end namespace Pythia8
Generated by   doxygen 1.8.5