JHUGen MELA  v2.4.1
Matrix element calculations as used in JHUGen. MELA is an important tool that was used for the Higgs boson discovery and for precise measurements of its structure and interactions. Please see the website https://spin.pha.jhu.edu/ and papers cited there for more details, and kindly cite those papers when using this code.
mod_Graviton.F90
Go to the documentation of this file.
1  module modgraviton
2  use modparameters
3  use modmisc
4  implicit none
5  private
6 
7 
8 !----- notation for subroutines
9  public :: evalamp_gg_g_vv,evalamp_qqb_g_vv,evalamp_g_vv
10 
11  contains
12 
13 
14 !----- a subroutinefor gg -> G -> ZZ/WW
15 !----- all outgoing convention and the following momentum assignment
16 !----- 0 -> g(p1) + g(p2) + e-(p3) + e+(p4) +mu-(p5) +mu+(p6)
17  subroutine evalamp_gg_g_vv(p,MY_IDUP,res)
18  implicit none
19  real(dp), intent(out) :: res
20  real(dp), intent(in) :: p(4,6)
21  integer, intent(in) :: my_idup(6:9)
22  complex(dp) :: a_vv(1:18), a0_vv(1:2)
23  integer :: i1,i2,i3,i4,vvmode,vvmode_swap
24  real(dp) :: prefactor
25  real(dp) :: intcolfac
26  integer :: ordering(1:4),ordering_swap(1:4)
27  logical :: dointerference
28 
29  if(isaquark(my_idup(6)) .and. isaquark(my_idup(8))) then
30  intcolfac=1.0_dp/3.0_dp
31  else
32  intcolfac=1.0_dp
33  endif
34 
35  call getdecay_vvmode_ordering(my_idup(6:9),vvmode,ordering,vvmode_swap,ordering_swap)
36 
37  if( vvmode.eq.zzmode ) then! Z decay
38  prefactor = 8d0*overallcouplvffsq**2
39  elseif( vvmode.eq.wwmode ) then ! W decay
40  prefactor = 8d0*overallcouplvffsq**2
41  elseif( vvmode.eq.zgmode ) then ! Z+photon "decay"
42  prefactor = 8d0*overallcouplvffsq ! Only single powers
43  elseif( vvmode.eq.ggmode ) then ! photon "decay"
44  prefactor = 8d0
45  else
46  prefactor = 0d0
47  endif
48 
49 
50  res = zero
51  a_vv(:) = 0d0
52  dointerference = includeinterference .and. ( &
53  ((vvmode.eq.zzmode) .and. (vvmode_swap.eq.zzmode)) &
54  )
55  if ( includevprime .and. .not.(vvmode.eq.zzmode .or. vvmode.eq.zgmode .or. vvmode.eq.wwmode) ) then
56  call error("Contact terms only for WW, ZZ or Zg!")
57  endif
58  do i1=1,2; do i2=1,2; do i3=1,2; do i4=1,2! sum over helicities
59  call calchelamp_gg(ordering,vvmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(1))
60  if( vvmode.eq.zzmode ) then
61  if( includegammastar ) then
62  call calchelamp_gg(ordering,zgsmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(3))
63  call calchelamp_gg(ordering,gszmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(5))
64  call calchelamp_gg(ordering,gsgsmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(7))
65  endif
66  if( includevprime ) then
67  call calchelamp_gg(ordering,zzpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(9))
68  call calchelamp_gg(ordering,zpzmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(11))
69  call calchelamp_gg(ordering,zpzpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(13))
70  endif
71  if( includegammastar .and. includevprime ) then
72  call calchelamp_gg(ordering,gszpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(15))
73  call calchelamp_gg(ordering,zpgsmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(17))
74  endif
75  elseif( vvmode.eq.zgmode ) then
76  if(includegammastar) then
77  call calchelamp_gg(ordering,gsgmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(3))
78  endif
79  if( includevprime ) then
80  call calchelamp_gg(ordering,zpgmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(5))
81  endif
82  elseif( vvmode.eq.wwmode .and. includevprime ) then
83  call calchelamp_gg(ordering,wwpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(9))
84  call calchelamp_gg(ordering,wpwmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(11))
85  call calchelamp_gg(ordering,wpwpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(13))
86  endif
87 
88  if( dointerference ) then
89  call calchelamp_gg(ordering_swap,vvmode_swap,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(2))
90  if( includegammastar ) then
91  call calchelamp_gg(ordering_swap,zgsmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(4))
92  call calchelamp_gg(ordering_swap,gszmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(6))
93  call calchelamp_gg(ordering_swap,gsgsmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(8))
94  endif
95  if( includevprime ) then
96  call calchelamp_gg(ordering_swap,zzpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(10))
97  call calchelamp_gg(ordering_swap,zpzmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(12))
98  call calchelamp_gg(ordering_swap,zpzpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(14))
99  endif
100  if( includegammastar .and. includevprime ) then
101  call calchelamp_gg(ordering_swap,gszpmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(16))
102  call calchelamp_gg(ordering_swap,zpgsmode,p(1:4,1:6),my_idup,i1,i2,i3,i4,a_vv(18))
103  endif
104  endif
105 
106  a0_vv(1) = a_vv(1)+a_vv(3)+a_vv(5)+a_vv(7)+a_vv(9)+a_vv(11)+a_vv(13)+a_vv(15)+a_vv(17) ! 3456 pieces
107  a0_vv(2) = a_vv(2)+a_vv(4)+a_vv(6)+a_vv(8)+a_vv(10)+a_vv(12)+a_vv(14)+a_vv(16)+a_vv(18) ! 5436 pieces
108  res = res + dreal(a0_vv(1)*dconjg(a0_vv(1)))
109  res = res + dreal(a0_vv(2)*dconjg(a0_vv(2)))
110  if( dointerference .and. (i3.eq.i4) ) then! interfere the 3456 with 5436 pieces
111  res = res - 2d0*intcolfac*dreal( a0_vv(1)*dconjg(a0_vv(2)) ) ! minus from Fermi statistics
112  endif
113  enddo; enddo; enddo; enddo
114 
115  res = res*prefactor
116  if( (vvmode.eq.zzmode) .and. dointerference ) res = res * symmfac
117 
118  end subroutine
119 
120  subroutine calchelamp_gg(ordering,VVMode,p,MY_IDUP,i1,i2,i3,i4,A)
121  implicit none
122  integer :: ordering(1:4),i1,i2,i3,i4,l1,l2,l3,l4,MY_IDUP(6:9),VVMode
123  real(dp) :: p(1:4,1:6)
124  complex(dp) :: propG
125  real(dp) :: pin(4,4)
126  complex(dp) :: A(1:1), sp(4,4)
127 
128 
129  l1=ordering(1)
130  l2=ordering(2)
131  l3=ordering(3)
132  l4=ordering(4)
133 
134  propg = one/dcmplx(2d0*scr(p(:,1),p(:,2)) - m_reso**2,m_reso*ga_reso)
135 
136  pin(1,:) = p(:,1)
137  pin(2,:) = p(:,2)
138  sp(1,:) = pol_mless2(dcmplx(p(:,1)),-3+2*i1,'in') ! gluon
139  sp(2,:) = pol_mless2(dcmplx(p(:,2)),-3+2*i2,'in') ! gluon
140  !sp(1,1:4)=pin(1,1:4);print *, "this checks IS gauge invariance"
141  !sp(2,1:4)=pin(2,1:4);print *, "this checks IS gauge invariance"
142  !careful: for gluon gauge invariance check the terms ~c3,c4 are needed because e1.q2 is not zero for e1-->q1
143  call getdecay_couplings_spinors_props( &
144  vvmode, &
145  (/my_idup(l1+3),my_idup(l2+3),my_idup(l3+3),my_idup(l4+3)/), &
146  (/p(:,l1),p(:,l2),p(:,l3),p(:,l4)/), &
147  -3+2*i3,-3+2*i4, &
148  sp(3:4,:),pin(3:4,:) &
149  )
150  call gggzzampl(vvmode,pin,sp,a(1))
151  a(1) = a(1)*propg
152 
153  end subroutine
154 
155  subroutine gggzzampl(VVMode,p,sp,res)
156  implicit none
157  integer, intent(in) :: VVMode
158  real(dp), intent(in) :: p(4,4)
159  complex(dp), intent(in) :: sp(4,4)
160  complex(dp), intent(out) :: res
161  complex(dp) :: e1_e2, e1_e3, e1_e4
162  complex(dp) :: e2_e3, e2_e4
163  complex(dp) :: e3_e4
164  complex(dp) :: q_q
165  complex(dp) :: q1_q2,q1_q3,q1_q4,q3_q3,q4_q4
166  complex(dp) :: q2_q3,q2_q4
167  complex(dp) :: q3_q4
168  complex(dp) :: q1_e3,q1_e4,q2_e3,q2_e4
169  complex(dp) :: e1_q3,e1_q4,e2_q3,e2_q4, q1_e1,q1_e2,q2_e1,q2_e2,q4_e1,q3_e1
170  complex(dp) :: e3_q4,e4_q3
171  complex(dp) :: q1(4),q2(4),q3(4),q4(4),q(4)
172  complex(dp) :: e1(4),e2(4),e3(4),e4(4)
173  complex(dp) :: xxx1,xxx2,xxx3,xxx4,yyy1,yyy2,yyy3,yyy4,yyy41,yyy42,yyy5,yyy6
174  complex(dp) :: yyy7,abr1
175  complex(dp) :: b_dyn(1:10)
176  real(dp) :: q34,MZ3,MZ4,MG
177  logical :: new
178  real(dp) :: rr_gam, rr
179 
180 
181  new = .true.
182 
183  b_dyn(:)=czero
184 
185  q1 = dcmplx(p(1,:),0d0)
186  q2 = dcmplx(p(2,:),0d0)
187  q3 = dcmplx(p(3,:),0d0)
188  q4 = dcmplx(p(4,:),0d0)
189 
190 
191  e1 = sp(1,:)
192  e2 = sp(2,:)
193  e3 = sp(3,:)
194  e4 = sp(4,:)
195 
196  q = -q1-q2
197  q_q = sc(q,q)
198  q3_q3 = sc(q3,q3)
199  q4_q4 = sc(q4,q4)
200 
201 
202  q1_q2 = sc(q1,q2)
203  q1_q3 = sc(q1,q3)
204  q1_q4 = sc(q1,q4)
205  q2_q3 = sc(q2,q3)
206  q2_q4 = sc(q2,q4)
207  q3_q4 = sc(q3,q4)
208 
209  e1_e2 = sc(e1,e2)
210  e1_e3 = sc(e1,e3)
211  e1_e4 = sc(e1,e4)
212  e2_e3 = sc(e2,e3)
213  e2_e4 = sc(e2,e4)
214  e3_e4 = sc(e3,e4)
215 
216 
217  q1_e3 = sc(q1,e3)
218  q1_e4 = sc(q1,e4)
219  q2_e3 = sc(q2,e3)
220  q2_e4 = sc(q2,e4)
221  e1_q3 = sc(e1,q3)
222  e1_q4 = sc(e1,q4)
223  e2_q3 = sc(e2,q3)
224  e2_q4 = sc(e2,q4)
225  e3_q4 = sc(e3,q4)
226  e4_q3 = sc(e4,q3)
227 
228  mz3=dsqrt(cdabs(q3_q3))
229  mz4=dsqrt(cdabs(q4_q4))
230  if( use_dynamic_mg ) then
231  mg = dsqrt(cdabs(q_q))
232  else
233  mg = m_reso
234  endif
235 
236 !---- define couplings
237  q34 = (mg**2-mz3**2-mz4**2)/2d0! = s = pV1.pV2 = (q_q-MZ3^2-MZ4^2)/2
238 
239  rr_gam = q_q/two/lambda**2! kappa for IS
240  xxx1 = (a1 + a2*rr_gam) ! those a's correspond to g's in eq.(7) for the IS
241  xxx2 = -a1/two + (a3+two*a4)*rr_gam
242  xxx3 = 4d0*a5*rr_gam * mg**2/q_q
243  ! for gluon gauge invariance check the terms ~c3,c4 are needed because e1.q2 is not zero for e1-->q1
244 
245  if (generate_bis) then
246  rr = q34/lambda**2! kappa for FS
247 
248  if( (vvmode.eq.zzmode) .or. (vvmode.eq.wwmode) ) then! decay ZZ's or WW's
249  b_dyn(1)=b1
250  b_dyn(2)=b2
251  b_dyn(3)=b3
252  b_dyn(4)=b4
253  b_dyn(5)=b5
254  b_dyn(6)=b6
255  b_dyn(7)=b7
256  b_dyn(8)=b8
257  b_dyn(9)=b9
258  b_dyn(10)=b10
259  elseif( (vvmode.eq.zgmode) .OR. (vvmode.eq.gszmode) .OR. (vvmode.eq.zgsmode) ) then
260  b_dyn(1)=bzgs1
261  b_dyn(2)=bzgs2
262  b_dyn(3)=bzgs3
263  b_dyn(4)=bzgs4
264  b_dyn(8)=bzgs8
265  elseif( (vvmode.eq.ggmode) .or. (vvmode.eq.gsgsmode) .or. (vvmode.eq.gsgmode) ) then
266  b_dyn(1)=bgsgs1
267  b_dyn(2)=bgsgs2
268  b_dyn(3)=bgsgs3
269  b_dyn(4)=bgsgs4
270  b_dyn(8)=bgsgs8
271  elseif( (vvmode.eq.zzpmode) .or. (vvmode.eq.wwpmode) .or. (vvmode.eq.zpzmode) .or. (vvmode.eq.wpwmode) ) then
272  b_dyn(1)=bzzp1
273  b_dyn(2)=bzzp2
274  b_dyn(3)=bzzp3
275  b_dyn(4)=bzzp4
276  b_dyn(5)=bzzp5
277  b_dyn(6)=bzzp6
278  b_dyn(7)=bzzp7
279  b_dyn(8)=bzzp8
280  b_dyn(9)=bzzp9
281  b_dyn(10)=bzzp10
282  elseif( (vvmode.eq.zpzpmode) .or. (vvmode.eq.wpwpmode) ) then
283  b_dyn(1)=bzpzp1
284  b_dyn(2)=bzpzp2
285  b_dyn(3)=bzpzp3
286  b_dyn(4)=bzpzp4
287  b_dyn(5)=bzpzp5
288  b_dyn(6)=bzpzp6
289  b_dyn(7)=bzpzp7
290  b_dyn(8)=bzpzp8
291  b_dyn(9)=bzpzp9
292  b_dyn(10)=bzpzp10
293  elseif( (vvmode.eq.zpgmode) .OR. (vvmode.eq.gszpmode) .OR. (vvmode.eq.zpgsmode) ) then
294  b_dyn(1)=bzpgs1
295  b_dyn(2)=bzpgs2
296  b_dyn(3)=bzpgs3
297  b_dyn(4)=bzpgs4
298  b_dyn(8)=bzpgs8
299  else
300  print *,"VVMode",vvmode,"not implemented"
301  endif
302 
303  yyy1 = q34*( b_dyn(1) + b_dyn(2)*rr*(one+mz3**2/q34)*(one+mz4**2/q34) ) + b_dyn(5)*m_v**2
304  yyy2 = -b_dyn(1)/two + b_dyn(3)*rr*(1d0-(mz3**2+mz4**2)/(2d0*q34)) + two*b_dyn(4)*rr + b_dyn(7)*rr*m_v**2/q34
305  yyy3 = (-b_dyn(2)/two - b_dyn(3)- two*b_dyn(4))*rr/q34
306  yyy41 = -b_dyn(1) - b_dyn(2)*(q34+mz3**2)/lambda**2 - b_dyn(3)*mz4**2/lambda**2 - 2d0*b_dyn(6)*m_v**2/lambda**2
307  yyy42 = -b_dyn(1) - b_dyn(2)*(q34+mz4**2)/lambda**2 - b_dyn(3)*mz3**2/lambda**2 - 2d0*b_dyn(6)*m_v**2/lambda**2
308  yyy5 = two*b_dyn(8)*rr*mg**2/q34
309  yyy6 = b_dyn(9) * m_v**2/lambda**2
310  yyy7 = b_dyn(10) * mg**2 * m_v**2/lambda**4
311 
312  else
313  yyy1 = q34*c1/2d0
314  yyy2 = c2
315  yyy3 = c3/mg**2
316  yyy41 = c41
317  yyy42 = c42
318  yyy5 = c5
319  yyy6 = czero
320  yyy7 = czero
321  if(vvmode.eq.zzmode .or. vvmode.eq.wwmode) then
322  yyy6 = c6
323  yyy7 = c7
324  endif
325  endif
326 
327  res = czero
328 
329 
330  if (new) then
331 
332 
333 ! res = &
334 ! + 8.D0*q1_e3*q1_e4*e1_e2*yyy1*xxx2 - 8.D0*q1_e3*q1_q3*e1_e2* &
335 ! e4_q3*yyy4*xxx2 + 4.D0*q1_e3*e1_e2*e4_q3*yyy1*xxx2 + 2.D0* &
336 ! q1_e3*e1_e2*e4_q3*MZ4**2*yyy4*xxx2 - 2.D0*q1_e3*e1_e2*e4_q3* &
337 ! MZ3**2*yyy4*xxx2 - 2.D0*q1_e3*e1_e2*e4_q3*MG**2*yyy4*xxx2 + 8.d0 &
338 ! *q1_e4*q1_q3*e1_e2*e3_q4*yyy4*xxx2 + 4.D0*q1_e4*e1_e2*e3_q4 &
339 ! *yyy1*xxx2 - 2.D0*q1_e4*e1_e2*e3_q4*MZ4**2*yyy4*xxx2 + 2.D0* &
340 ! q1_e4*e1_e2*e3_q4*MZ3**2*yyy4*xxx2 + 2.D0*q1_e4*e1_e2*e3_q4* &
341 ! MG**2*yyy4*xxx2 - 8.D0*q1_q3*e1_e2*e3_e4*MZ4**2*yyy2*xxx2 + 8.d0 &
342 ! *q1_q3*e1_e2*e3_e4*MZ3**2*yyy2*xxx2 + 8.D0*q1_q3*e1_e2* &
343 ! e3_e4*MG**2*yyy2*xxx2 - 8.D0*q1_q3*e1_e2*e3_q4*e4_q3*MZ4**2* &
344 ! yyy3*xxx2 + 8.D0*q1_q3*e1_e2*e3_q4*e4_q3*MZ3**2*yyy3*xxx2 + 8.d0 &
345 ! *q1_q3*e1_e2*e3_q4*e4_q3*MG**2*yyy3*xxx2 + 16.D0*q1_q3**2* &
346 ! e1_e2*e3_e4*yyy2*xxx2 + 16.D0*q1_q3**2*e1_e2*e3_q4*e4_q3*yyy3 &
347 ! *xxx2 + 2.D0/3.D0*e1_e2*e3_e4*MZ4**4*yyy2*xxx2 + 1.D0/3.D0* &
348 ! e1_e2*e3_e4*MZ4**4*yyy2*xxx1
349 ! res = res - 4.D0/3.D0*e1_e2*e3_e4*MZ3**2*MZ4**2*yyy2*xxx2 - 2.D0/ &
350 ! 3.D0*e1_e2*e3_e4*MZ3**2*MZ4**2*yyy2*xxx1 + 2.D0/3.D0*e1_e2* &
351 ! e3_e4*MZ3**4*yyy2*xxx2 + 1.D0/3.D0*e1_e2*e3_e4*MZ3**4*yyy2* &
352 ! xxx1 + 2.D0/3.D0*e1_e2*e3_e4*MG**2*yyy1*xxx2 - 2.D0/3.D0* &
353 ! e1_e2*e3_e4*MG**2*yyy1*xxx1 - 4.D0/3.D0*e1_e2*e3_e4*MG**2* &
354 ! MZ4**2*yyy2*xxx2 - 2.D0/3.D0*e1_e2*e3_e4*MG**2*MZ4**2*yyy2* &
355 ! xxx1 + 8.D0/3.D0*e1_e2*e3_e4*MG**2*MZ3**2*yyy2*xxx2 - 2.D0/3.D0 &
356 ! *e1_e2*e3_e4*MG**2*MZ3**2*yyy2*xxx1 + 2.D0/3.D0*e1_e2*e3_e4* &
357 ! MG**4*yyy2*xxx2 + 1.D0/3.D0*e1_e2*e3_e4*MG**4*yyy2*xxx1 + 4.D0 &
358 ! /3.D0*e1_e2*e3_q4*e4_q3*yyy1*xxx2 + 2.D0/3.D0*e1_e2*e3_q4* &
359 ! e4_q3*yyy1*xxx1 + 2.D0/3.D0*e1_e2*e3_q4*e4_q3*MZ4**4*yyy3* &
360 ! xxx2 + 1.D0/3.D0*e1_e2*e3_q4*e4_q3*MZ4**4*yyy3*xxx1 - 4.D0/3.D0 &
361 ! *e1_e2*e3_q4*e4_q3*MZ3**2*MZ4**2*yyy3*xxx2 - 2.D0/3.D0*e1_e2 &
362 ! *e3_q4*e4_q3*MZ3**2*MZ4**2*yyy3*xxx1 + 2.D0/3.D0*e1_e2*e3_q4* &
363 ! e4_q3*MZ3**4*yyy3*xxx2 + 1.D0/3.D0*e1_e2*e3_q4*e4_q3*MZ3**4* &
364 ! yyy3*xxx1
365 ! res = res + 2.D0/3.D0*e1_e2*e3_q4*e4_q3*MG**2*yyy4*xxx2 - 2.D0/3.D0 &
366 ! *e1_e2*e3_q4*e4_q3*MG**2*yyy4*xxx1 - 4.D0/3.D0*e1_e2*e3_q4* &
367 ! e4_q3*MG**2*MZ4**2*yyy3*xxx2 - 2.D0/3.D0*e1_e2*e3_q4*e4_q3* &
368 ! MG**2*MZ4**2*yyy3*xxx1 + 8.D0/3.D0*e1_e2*e3_q4*e4_q3*MG**2* &
369 ! MZ3**2*yyy3*xxx2 - 2.D0/3.D0*e1_e2*e3_q4*e4_q3*MG**2*MZ3**2* &
370 ! yyy3*xxx1 + 2.D0/3.D0*e1_e2*e3_q4*e4_q3*MG**4*yyy3*xxx2 + 1.D0 &
371 ! /3.D0*e1_e2*e3_q4*e4_q3*MG**4*yyy3*xxx1 + e1_e3*e2_e4*MG**2* &
372 ! yyy1*xxx1 - e1_e3*e2_q3*e4_q3*MG**2*yyy4*xxx1 + e1_e4*e2_e3* &
373 ! MG**2*yyy1*xxx1 + e1_e4*e2_q3*e3_q4*MG**2*yyy4*xxx1 - e1_q3* &
374 ! e2_e3*e4_q3*MG**2*yyy4*xxx1 + e1_q3*e2_e4*e3_q4*MG**2*yyy4* &
375 ! xxx1 + 4.D0*e1_q3*e2_q3*e3_e4*MG**2*yyy2*xxx1 + 4.D0*e1_q3* &
376 ! e2_q3*e3_q4*e4_q3*MG**2*yyy3*xxx1 + 8.D0*et1(q1,q2,e1,e2)* &
377 ! q1_e3*q1_e4*MG**(-2)*yyy1*xxx3 - 8.D0*et1(q1,q2,e1,e2)*q1_e3* &
378 ! q1_q3*e4_q3*MG**(-2)*yyy4*xxx3 + 4.D0*et1(q1,q2,e1,e2)*q1_e3* &
379 ! e4_q3*MG**(-2)*yyy1*xxx3 + 2.D0*et1(q1,q2,e1,e2)*q1_e3*e4_q3* &
380 ! MG**(-2)*MZ4**2*yyy4*xxx3
381 ! res = res - 2.D0*et1(q1,q2,e1,e2)*q1_e3*e4_q3*MG**(-2)*MZ3**2* &
382 ! yyy4*xxx3 - 2.D0*et1(q1,q2,e1,e2)*q1_e3*e4_q3*yyy4*xxx3 + 8.D0* &
383 ! et1(q1,q2,e1,e2)*q1_e4*q1_q3*e3_q4*MG**(-2)*yyy4*xxx3 + 4.D0* &
384 ! et1(q1,q2,e1,e2)*q1_e4*e3_q4*MG**(-2)*yyy1*xxx3 - 2.D0*et1(q1 &
385 ! ,q2,e1,e2)*q1_e4*e3_q4*MG**(-2)*MZ4**2*yyy4*xxx3 + 2.D0*et1( &
386 ! q1,q2,e1,e2)*q1_e4*e3_q4*MG**(-2)*MZ3**2*yyy4*xxx3 + 2.D0* &
387 ! et1(q1,q2,e1,e2)*q1_e4*e3_q4*yyy4*xxx3 - 8.D0*et1(q1,q2,e1,e2 &
388 ! )*q1_q3*e3_e4*MG**(-2)*MZ4**2*yyy2*xxx3 + 8.D0*et1(q1,q2,e1, &
389 ! e2)*q1_q3*e3_e4*MG**(-2)*MZ3**2*yyy2*xxx3 + 8.D0*et1(q1,q2,e1 &
390 ! ,e2)*q1_q3*e3_e4*yyy2*xxx3 - 8.D0*et1(q1,q2,e1,e2)*q1_q3* &
391 ! e3_q4*e4_q3*MG**(-2)*MZ4**2*yyy3*xxx3 + 8.D0*et1(q1,q2,e1,e2) &
392 ! *q1_q3*e3_q4*e4_q3*MG**(-2)*MZ3**2*yyy3*xxx3 + 8.D0*et1(q1,q2 &
393 ! ,e1,e2)*q1_q3*e3_q4*e4_q3*yyy3*xxx3 + 16.D0*et1(q1,q2,e1,e2)* &
394 ! q1_q3**2*e3_e4*MG**(-2)*yyy2*xxx3 + 16.D0*et1(q1,q2,e1,e2)* &
395 ! q1_q3**2*e3_q4*e4_q3*MG**(-2)*yyy3*xxx3 + 2.D0/3.D0*et1(q1,q2 &
396 ! ,e1,e2)*e3_e4*MG**(-2)*MZ4**4*yyy2*xxx3
397 ! res = res - 4.D0/3.D0*et1(q1,q2,e1,e2)*e3_e4*MG**(-2)*MZ3**2* &
398 ! MZ4**2*yyy2*xxx3 + 2.D0/3.D0*et1(q1,q2,e1,e2)*e3_e4*MG**(-2)* &
399 ! MZ3**4*yyy2*xxx3 + 2.D0/3.D0*et1(q1,q2,e1,e2)*e3_e4*yyy1*xxx3 &
400 ! - 4.D0/3.D0*et1(q1,q2,e1,e2)*e3_e4*MZ4**2*yyy2*xxx3 + 8.D0/3.D0 &
401 ! *et1(q1,q2,e1,e2)*e3_e4*MZ3**2*yyy2*xxx3 + 2.D0/3.D0*et1(q1 &
402 ! ,q2,e1,e2)*e3_e4*MG**2*yyy2*xxx3 + 4.D0/3.D0*et1(q1,q2,e1,e2) &
403 ! *e3_q4*e4_q3*MG**(-2)*yyy1*xxx3 + 2.D0/3.D0*et1(q1,q2,e1,e2)* &
404 ! e3_q4*e4_q3*MG**(-2)*MZ4**4*yyy3*xxx3 - 4.D0/3.D0*et1(q1,q2, &
405 ! e1,e2)*e3_q4*e4_q3*MG**(-2)*MZ3**2*MZ4**2*yyy3*xxx3 + 2.D0/3.D0 &
406 ! *et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**(-2)*MZ3**4*yyy3*xxx3 + 2.D0 &
407 ! /3.D0*et1(q1,q2,e1,e2)*e3_q4*e4_q3*yyy4*xxx3 - 4.D0/3.D0* &
408 ! et1(q1,q2,e1,e2)*e3_q4*e4_q3*MZ4**2*yyy3*xxx3 + 8.D0/3.D0* &
409 ! et1(q1,q2,e1,e2)*e3_q4*e4_q3*MZ3**2*yyy3*xxx3 + 2.D0/3.D0* &
410 ! et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**2*yyy3*xxx3 - et1(q1,q2,e1, &
411 ! e2)*et1(q1,q,e3,e4)*MG**(-2)*MZ4**2*yyy6*xxx3 + et1(q1,q2,e1, &
412 ! e2)*et1(q1,q,e3,e4)*MG**(-2)*MZ3**2*yyy6*xxx3
413 ! res = res + et1(q1,q2,e1,e2)*et1(q1,q,e3,e4)*yyy6*xxx3 + 4.D0* &
414 ! et1(q1,q2,e1,e2)*et1(q1,q,e3,e4)*q1_q3*MG**(-2)*yyy6*xxx3 - 4.0d0 &
415 ! *et1(q1,q2,e1,e2)*et1(q1,q,e3,q3)*q1_q3*e4_q3*MG**(-4)*yyy7 &
416 ! *xxx3 + et1(q1,q2,e1,e2)*et1(q1,q,e3,q3)*e4_q3*MG**(-4)* &
417 ! MZ4**2*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e3,q3)*e4_q3* &
418 ! MG**(-4)*MZ3**2*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e3,q3)* &
419 ! e4_q3*MG**(-2)*yyy7*xxx3 + 4.D0*et1(q1,q2,e1,e2)*et1(q1,q,e3, &
420 ! q4)*q1_q3*e4_q3*MG**(-4)*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1, &
421 ! q,e3,q4)*e4_q3*MG**(-4)*MZ4**2*yyy7*xxx3 + et1(q1,q2,e1,e2)* &
422 ! et1(q1,q,e3,q4)*e4_q3*MG**(-4)*MZ3**2*yyy7*xxx3 + et1(q1,q2, &
423 ! e1,e2)*et1(q1,q,e3,q4)*e4_q3*MG**(-2)*yyy7*xxx3 - 4.D0*et1(q1 &
424 ! ,q2,e1,e2)*et1(q1,q,e4,q3)*q1_q3*e3_q4*MG**(-4)*yyy7*xxx3 + &
425 ! et1(q1,q2,e1,e2)*et1(q1,q,e4,q3)*e3_q4*MG**(-4)*MZ4**2*yyy7* &
426 ! xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e4,q3)*e3_q4*MG**(-4)*MZ3**2 &
427 ! *yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e4,q3)*e3_q4*MG**(-2)* &
428 ! yyy7*xxx3
429 ! res = res + 4.D0*et1(q1,q2,e1,e2)*et1(q1,q,e4,q4)*q1_q3*e3_q4* &
430 ! MG**(-4)*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e4,q4)*e3_q4* &
431 ! MG**(-4)*MZ4**2*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q1,q,e4,q4)* &
432 ! e3_q4*MG**(-4)*MZ3**2*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q1,q, &
433 ! e4,q4)*e3_q4*MG**(-2)*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q2,q, &
434 ! e3,e4)*MG**(-2)*MZ4**2*yyy6*xxx3 - et1(q1,q2,e1,e2)*et1(q2,q, &
435 ! e3,e4)*MG**(-2)*MZ3**2*yyy6*xxx3 - et1(q1,q2,e1,e2)*et1(q2,q, &
436 ! e3,e4)*yyy6*xxx3 - 4.D0*et1(q1,q2,e1,e2)*et1(q2,q,e3,e4)* &
437 ! q1_q3*MG**(-2)*yyy6*xxx3 + 4.D0*et1(q1,q2,e1,e2)*et1(q2,q,e3, &
438 ! q3)*q1_q3*e4_q3*MG**(-4)*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q2, &
439 ! q,e3,q3)*e4_q3*MG**(-4)*MZ4**2*yyy7*xxx3 + et1(q1,q2,e1,e2)* &
440 ! et1(q2,q,e3,q3)*e4_q3*MG**(-4)*MZ3**2*yyy7*xxx3 + et1(q1,q2, &
441 ! e1,e2)*et1(q2,q,e3,q3)*e4_q3*MG**(-2)*yyy7*xxx3 - 4.D0*et1(q1 &
442 ! ,q2,e1,e2)*et1(q2,q,e3,q4)*q1_q3*e4_q3*MG**(-4)*yyy7*xxx3 + &
443 ! et1(q1,q2,e1,e2)*et1(q2,q,e3,q4)*e4_q3*MG**(-4)*MZ4**2*yyy7* &
444 ! xxx3
445 ! res = res - et1(q1,q2,e1,e2)*et1(q2,q,e3,q4)*e4_q3*MG**(-4)* &
446 ! MZ3**2*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q2,q,e3,q4)*e4_q3* &
447 ! MG**(-2)*yyy7*xxx3 + 4.D0*et1(q1,q2,e1,e2)*et1(q2,q,e4,q3)* &
448 ! q1_q3*e3_q4*MG**(-4)*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q2,q,e4 &
449 ! ,q3)*e3_q4*MG**(-4)*MZ4**2*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1( &
450 ! q2,q,e4,q3)*e3_q4*MG**(-4)*MZ3**2*yyy7*xxx3 + et1(q1,q2,e1,e2 &
451 ! )*et1(q2,q,e4,q3)*e3_q4*MG**(-2)*yyy7*xxx3 - 4.D0*et1(q1,q2, &
452 ! e1,e2)*et1(q2,q,e4,q4)*q1_q3*e3_q4*MG**(-4)*yyy7*xxx3 + et1( &
453 ! q1,q2,e1,e2)*et1(q2,q,e4,q4)*e3_q4*MG**(-4)*MZ4**2*yyy7*xxx3 &
454 ! - et1(q1,q2,e1,e2)*et1(q2,q,e4,q4)*e3_q4*MG**(-4)*MZ3**2* &
455 ! yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q2,q,e4,q4)*e3_q4*MG**(-2)* &
456 ! yyy7*xxx3 - 1.D0/3.D0*et1(q1,q2,e1,e2)*et1(q,e3,e4,q3)*yyy6* &
457 ! xxx3 + 1.D0/3.D0*et1(q1,q2,e1,e2)*et1(q,e3,e4,q4)*yyy6*xxx3 &
458 ! + 2.D0/3.D0*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*MG**(-4)* &
459 ! MZ4**4*yyy5*xxx3 - 4.D0/3.D0*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4 &
460 ! )*MG**(-4)*MZ3**2*MZ4**2*yyy5*xxx3
461 ! res = res + 2.D0/3.D0*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*MG**(-4)* &
462 ! MZ3**4*yyy5*xxx3 - 4.D0/3.D0*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)* &
463 ! MG**(-2)*MZ4**2*yyy5*xxx3 + 8.D0/3.D0*et1(q1,q2,e1,e2)*et1(e3 &
464 ! ,e4,q3,q4)*MG**(-2)*MZ3**2*yyy5*xxx3 + 2.D0/3.D0*et1(q1,q2,e1 &
465 ! ,e2)*et1(e3,e4,q3,q4)*yyy5*xxx3 - 8.D0*et1(q1,q2,e1,e2)*et1( &
466 ! e3,e4,q3,q4)*q1_q3*MG**(-4)*MZ4**2*yyy5*xxx3 + 8.D0*et1(q1,q2 &
467 ! ,e1,e2)*et1(e3,e4,q3,q4)*q1_q3*MG**(-4)*MZ3**2*yyy5*xxx3 + 8.D0 &
468 ! *et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*q1_q3*MG**(-2)*yyy5*xxx3 &
469 ! + 16.D0*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*q1_q3**2*MG**(-4)* &
470 ! yyy5*xxx3 + 4.D0*et1(q1,q,e3,e4)*q1_q3*e1_e2*yyy6*xxx2 - et1( &
471 ! q1,q,e3,e4)*e1_e2*MZ4**2*yyy6*xxx2 + et1(q1,q,e3,e4)*e1_e2* &
472 ! MZ3**2*yyy6*xxx2 + et1(q1,q,e3,e4)*e1_e2*MG**2*yyy6*xxx2 - 4.D0 &
473 ! *et1(q1,q,e3,q3)*q1_q3*e1_e2*e4_q3*MG**(-2)*yyy7*xxx2 + et1( &
474 ! q1,q,e3,q3)*e1_e2*e4_q3*MG**(-2)*MZ4**2*yyy7*xxx2 - et1(q1,q, &
475 ! e3,q3)*e1_e2*e4_q3*MG**(-2)*MZ3**2*yyy7*xxx2 - et1(q1,q,e3,q3 &
476 ! )*e1_e2*e4_q3*yyy7*xxx2
477 ! res = res + 4.D0*et1(q1,q,e3,q4)*q1_q3*e1_e2*e4_q3*MG**(-2)*yyy7* &
478 ! xxx2 - et1(q1,q,e3,q4)*e1_e2*e4_q3*MG**(-2)*MZ4**2*yyy7*xxx2 + &
479 ! et1(q1,q,e3,q4)*e1_e2*e4_q3*MG**(-2)*MZ3**2*yyy7*xxx2 + et1( &
480 ! q1,q,e3,q4)*e1_e2*e4_q3*yyy7*xxx2 - 4.D0*et1(q1,q,e4,q3)* &
481 ! q1_q3*e1_e2*e3_q4*MG**(-2)*yyy7*xxx2 + et1(q1,q,e4,q3)*e1_e2* &
482 ! e3_q4*MG**(-2)*MZ4**2*yyy7*xxx2 - et1(q1,q,e4,q3)*e1_e2*e3_q4 &
483 ! *MG**(-2)*MZ3**2*yyy7*xxx2 - et1(q1,q,e4,q3)*e1_e2*e3_q4*yyy7 &
484 ! *xxx2 + 4.D0*et1(q1,q,e4,q4)*q1_q3*e1_e2*e3_q4*MG**(-2)*yyy7* &
485 ! xxx2 - et1(q1,q,e4,q4)*e1_e2*e3_q4*MG**(-2)*MZ4**2*yyy7*xxx2 &
486 ! + et1(q1,q,e4,q4)*e1_e2*e3_q4*MG**(-2)*MZ3**2*yyy7*xxx2 + &
487 ! et1(q1,q,e4,q4)*e1_e2*e3_q4*yyy7*xxx2 - 4.D0*et1(q2,q,e3,e4)* &
488 ! q1_q3*e1_e2*yyy6*xxx2 + et1(q2,q,e3,e4)*e1_e2*MZ4**2*yyy6* &
489 ! xxx2 - et1(q2,q,e3,e4)*e1_e2*MZ3**2*yyy6*xxx2 - et1(q2,q,e3, &
490 ! e4)*e1_e2*MG**2*yyy6*xxx2 + 4.D0*et1(q2,q,e3,q3)*q1_q3*e1_e2* &
491 ! e4_q3*MG**(-2)*yyy7*xxx2 - et1(q2,q,e3,q3)*e1_e2*e4_q3* &
492 ! MG**(-2)*MZ4**2*yyy7*xxx2
493 ! res = res + et1(q2,q,e3,q3)*e1_e2*e4_q3*MG**(-2)*MZ3**2*yyy7*xxx2 &
494 ! + et1(q2,q,e3,q3)*e1_e2*e4_q3*yyy7*xxx2 - 4.D0*et1(q2,q,e3, &
495 ! q4)*q1_q3*e1_e2*e4_q3*MG**(-2)*yyy7*xxx2 + et1(q2,q,e3,q4)* &
496 ! e1_e2*e4_q3*MG**(-2)*MZ4**2*yyy7*xxx2 - et1(q2,q,e3,q4)*e1_e2 &
497 ! *e4_q3*MG**(-2)*MZ3**2*yyy7*xxx2 - et1(q2,q,e3,q4)*e1_e2* &
498 ! e4_q3*yyy7*xxx2 + 4.D0*et1(q2,q,e4,q3)*q1_q3*e1_e2*e3_q4* &
499 ! MG**(-2)*yyy7*xxx2 - et1(q2,q,e4,q3)*e1_e2*e3_q4*MG**(-2)* &
500 ! MZ4**2*yyy7*xxx2 + et1(q2,q,e4,q3)*e1_e2*e3_q4*MG**(-2)* &
501 ! MZ3**2*yyy7*xxx2 + et1(q2,q,e4,q3)*e1_e2*e3_q4*yyy7*xxx2 - 4.D0 &
502 ! *et1(q2,q,e4,q4)*q1_q3*e1_e2*e3_q4*MG**(-2)*yyy7*xxx2 + et1( &
503 ! q2,q,e4,q4)*e1_e2*e3_q4*MG**(-2)*MZ4**2*yyy7*xxx2 - et1(q2,q, &
504 ! e4,q4)*e1_e2*e3_q4*MG**(-2)*MZ3**2*yyy7*xxx2 - et1(q2,q,e4,q4 &
505 ! )*e1_e2*e3_q4*yyy7*xxx2 - et1(q,e1,e3,e4)*e2_q3*MG**2*yyy6* &
506 ! xxx1 + et1(q,e1,e3,q3)*e2_q3*e4_q3*yyy7*xxx1 - et1(q,e1,e3,q4 &
507 ! )*e2_q3*e4_q3*yyy7*xxx1 + et1(q,e1,e4,q3)*e2_q3*e3_q4*yyy7* &
508 ! xxx1
509 ! res = res - et1(q,e1,e4,q4)*e2_q3*e3_q4*yyy7*xxx1 - et1(q,e2,e3, &
510 ! e4)*e1_q3*MG**2*yyy6*xxx1 + et1(q,e2,e3,q3)*e1_q3*e4_q3*yyy7* &
511 ! xxx1 - et1(q,e2,e3,q4)*e1_q3*e4_q3*yyy7*xxx1 + et1(q,e2,e4,q3 &
512 ! )*e1_q3*e3_q4*yyy7*xxx1 - et1(q,e2,e4,q4)*e1_q3*e3_q4*yyy7* &
513 ! xxx1 - 1.D0/3.D0*et1(q,e3,e4,q3)*e1_e2*MG**2*yyy6*xxx2 + 1.D0/ &
514 ! 3.D0*et1(q,e3,e4,q3)*e1_e2*MG**2*yyy6*xxx1 + 1.D0/3.D0*et1(q, &
515 ! e3,e4,q4)*e1_e2*MG**2*yyy6*xxx2 - 1.D0/3.D0*et1(q,e3,e4,q4)* &
516 ! e1_e2*MG**2*yyy6*xxx1 - 8.D0*et1(e3,e4,q3,q4)*q1_q3*e1_e2* &
517 ! MG**(-2)*MZ4**2*yyy5*xxx2 + 8.D0*et1(e3,e4,q3,q4)*q1_q3*e1_e2 &
518 ! *MG**(-2)*MZ3**2*yyy5*xxx2 + 8.D0*et1(e3,e4,q3,q4)*q1_q3* &
519 ! e1_e2*yyy5*xxx2 + 16.D0*et1(e3,e4,q3,q4)*q1_q3**2*e1_e2* &
520 ! MG**(-2)*yyy5*xxx2 + 2.D0/3.D0*et1(e3,e4,q3,q4)*e1_e2* &
521 ! MG**(-2)*MZ4**4*yyy5*xxx2 + 1.D0/3.D0*et1(e3,e4,q3,q4)*e1_e2* &
522 ! MG**(-2)*MZ4**4*yyy5*xxx1 - 4.D0/3.D0*et1(e3,e4,q3,q4)*e1_e2* &
523 ! MG**(-2)*MZ3**2*MZ4**2*yyy5*xxx2 - 2.D0/3.D0*et1(e3,e4,q3,q4) &
524 ! *e1_e2*MG**(-2)*MZ3**2*MZ4**2*yyy5*xxx1
525 ! res = res + 2.D0/3.D0*et1(e3,e4,q3,q4)*e1_e2*MG**(-2)*MZ3**4*yyy5 &
526 ! *xxx2 + 1.D0/3.D0*et1(e3,e4,q3,q4)*e1_e2*MG**(-2)*MZ3**4*yyy5* &
527 ! xxx1 - 4.D0/3.D0*et1(e3,e4,q3,q4)*e1_e2*MZ4**2*yyy5*xxx2 - 2.D0 &
528 ! /3.D0*et1(e3,e4,q3,q4)*e1_e2*MZ4**2*yyy5*xxx1 + 8.D0/3.D0* &
529 ! et1(e3,e4,q3,q4)*e1_e2*MZ3**2*yyy5*xxx2 - 2.D0/3.D0*et1(e3,e4 &
530 ! ,q3,q4)*e1_e2*MZ3**2*yyy5*xxx1 + 2.D0/3.D0*et1(e3,e4,q3,q4)* &
531 ! e1_e2*MG**2*yyy5*xxx2 + 1.D0/3.D0*et1(e3,e4,q3,q4)*e1_e2* &
532 ! MG**2*yyy5*xxx1 + 4.D0*et1(e3,e4,q3,q4)*e1_q3*e2_q3*yyy5*xxx1
533 !
534 ! print *, "old res GG",res
535 
536 
537 ! this is the new code that includes couplings yyy41 and yyy42 instead of yyy4
538 ! res =&
539 ! & + 8.*q1_e3*q1_e4*e1_e2*yyy1*xxx2 - 8.*q1_e3*q1_q3*e1_e2*e4_q3*&
540 ! & yyy41*xxx2 + 4.*q1_e3*e1_e2*e4_q3*yyy1*xxx2 + 2.*q1_e3*e1_e2*&
541 ! & e4_q3*MZ4**2*yyy41*xxx2 - 2.*q1_e3*e1_e2*e4_q3*MZ3**2*yyy41*&
542 ! & xxx2 - 2.*q1_e3*e1_e2*e4_q3*MG**2*yyy41*xxx2 + 8.*q1_e4*q1_q3&
543 ! & *e1_e2*e3_q4*yyy42*xxx2 + 4.*q1_e4*e1_e2*e3_q4*yyy1*xxx2 - 2.&
544 ! & *q1_e4*e1_e2*e3_q4*MZ4**2*yyy42*xxx2 + 2.*q1_e4*e1_e2*e3_q4*&
545 ! & MZ3**2*yyy42*xxx2 + 2.*q1_e4*e1_e2*e3_q4*MG**2*yyy42*xxx2 - 8.&
546 ! & *q1_q3*e1_e2*e3_e4*MZ4**2*yyy2*xxx2 + 8.*q1_q3*e1_e2*e3_e4*&
547 ! & MZ3**2*yyy2*xxx2 + 8.*q1_q3*e1_e2*e3_e4*MG**2*yyy2*xxx2 + 4.*&
548 ! & q1_q3*e1_e2*e3_q4*e4_q3*yyy42*xxx2 - 4.*q1_q3*e1_e2*e3_q4*&
549 ! & e4_q3*yyy41*xxx2 - 8.*q1_q3*e1_e2*e3_q4*e4_q3*MZ4**2*yyy3*&
550 ! & xxx2 + 8.*q1_q3*e1_e2*e3_q4*e4_q3*MZ3**2*yyy3*xxx2 + 8.*q1_q3&
551 ! & *e1_e2*e3_q4*e4_q3*MG**2*yyy3*xxx2 + 16.*q1_q3**2*e1_e2*e3_e4&
552 ! & *yyy2*xxx2 + 16.*q1_q3**2*e1_e2*e3_q4*e4_q3*yyy3*xxx2 + 2./3.&
553 ! & *e1_e2*e3_e4*MZ4**4*yyy2*xxx2
554 ! res = res + 1./3.*e1_e2*e3_e4*MZ4**4*yyy2*xxx1 - 4./3.*e1_e2*&
555 ! & e3_e4*MZ3**2*MZ4**2*yyy2*xxx2 - 2./3.*e1_e2*e3_e4*MZ3**2*&
556 ! & MZ4**2*yyy2*xxx1 + 2./3.*e1_e2*e3_e4*MZ3**4*yyy2*xxx2 + 1./3.&
557 ! & *e1_e2*e3_e4*MZ3**4*yyy2*xxx1 + 2./3.*e1_e2*e3_e4*MG**2*yyy1*&
558 ! & xxx2 - 2./3.*e1_e2*e3_e4*MG**2*yyy1*xxx1 - 4./3.*e1_e2*e3_e4*&
559 ! & MG**2*MZ4**2*yyy2*xxx2 - 2./3.*e1_e2*e3_e4*MG**2*MZ4**2*yyy2*&
560 ! & xxx1 + 8./3.*e1_e2*e3_e4*MG**2*MZ3**2*yyy2*xxx2 - 2./3.*e1_e2&
561 ! & *e3_e4*MG**2*MZ3**2*yyy2*xxx1 + 2./3.*e1_e2*e3_e4*MG**4*yyy2*&
562 ! & xxx2 + 1./3.*e1_e2*e3_e4*MG**4*yyy2*xxx1 + 4./3.*e1_e2*e3_q4*&
563 ! & e4_q3*yyy1*xxx2 + 2./3.*e1_e2*e3_q4*e4_q3*yyy1*xxx1 - 2./3.*&
564 ! & e1_e2*e3_q4*e4_q3*MZ4**2*yyy42*xxx2 - 1./3.*e1_e2*e3_q4*e4_q3&
565 ! & *MZ4**2*yyy42*xxx1 + 2./3.*e1_e2*e3_q4*e4_q3*MZ4**2*yyy41*&
566 ! & xxx2 + 1./3.*e1_e2*e3_q4*e4_q3*MZ4**2*yyy41*xxx1 + 2./3.*&
567 ! & e1_e2*e3_q4*e4_q3*MZ4**4*yyy3*xxx2 + 1./3.*e1_e2*e3_q4*e4_q3*&
568 ! & MZ4**4*yyy3*xxx1 + 2./3.*e1_e2*e3_q4*e4_q3*MZ3**2*yyy42*xxx2&
569 ! & + 1./3.*e1_e2*e3_q4*e4_q3*MZ3**2*yyy42*xxx1
570 ! res = res - 2./3.*e1_e2*e3_q4*e4_q3*MZ3**2*yyy41*xxx2 - 1./3.*&
571 ! & e1_e2*e3_q4*e4_q3*MZ3**2*yyy41*xxx1 - 4./3.*e1_e2*e3_q4*e4_q3&
572 ! & *MZ3**2*MZ4**2*yyy3*xxx2 - 2./3.*e1_e2*e3_q4*e4_q3*MZ3**2*&
573 ! & MZ4**2*yyy3*xxx1 + 2./3.*e1_e2*e3_q4*e4_q3*MZ3**4*yyy3*xxx2&
574 ! & + 1./3.*e1_e2*e3_q4*e4_q3*MZ3**4*yyy3*xxx1 + 4./3.*e1_e2*&
575 ! & e3_q4*e4_q3*MG**2*yyy42*xxx2 - 1./3.*e1_e2*e3_q4*e4_q3*MG**2*&
576 ! & yyy42*xxx1 - 2./3.*e1_e2*e3_q4*e4_q3*MG**2*yyy41*xxx2 - 1./3.&
577 ! & *e1_e2*e3_q4*e4_q3*MG**2*yyy41*xxx1 - 4./3.*e1_e2*e3_q4*e4_q3&
578 ! & *MG**2*MZ4**2*yyy3*xxx2 - 2./3.*e1_e2*e3_q4*e4_q3*MG**2*&
579 ! & MZ4**2*yyy3*xxx1 + 8./3.*e1_e2*e3_q4*e4_q3*MG**2*MZ3**2*yyy3*&
580 ! & xxx2 - 2./3.*e1_e2*e3_q4*e4_q3*MG**2*MZ3**2*yyy3*xxx1 + 2./3.&
581 ! & *e1_e2*e3_q4*e4_q3*MG**4*yyy3*xxx2 + 1./3.*e1_e2*e3_q4*e4_q3*&
582 ! & MG**4*yyy3*xxx1 + e1_e3*e2_e4*MG**2*yyy1*xxx1 - e1_e3*e2_q3*&
583 ! & e4_q3*MG**2*yyy41*xxx1 + e1_e4*e2_e3*MG**2*yyy1*xxx1 + e1_e4*&
584 ! & e2_q3*e3_q4*MG**2*yyy42*xxx1 - e1_q3*e2_e3*e4_q3*MG**2*yyy41*&
585 ! & xxx1
586 ! res = res + e1_q3*e2_e4*e3_q4*MG**2*yyy42*xxx1 + 4.*e1_q3*e2_q3*&
587 ! & e3_e4*MG**2*yyy2*xxx1 + 4.*e1_q3*e2_q3*e3_q4*e4_q3*MG**2*yyy3&
588 ! & *xxx1 + 8.*et1(q1,q2,e1,e2)*q1_e3*q1_e4*MG**(-2)*yyy1*xxx3 - &
589 ! & 8.*et1(q1,q2,e1,e2)*q1_e3*q1_q3*e4_q3*MG**(-2)*yyy41*xxx3 + 4.&
590 ! & *et1(q1,q2,e1,e2)*q1_e3*e4_q3*MG**(-2)*yyy1*xxx3 + 2.*et1(q1,&
591 ! & q2,e1,e2)*q1_e3*e4_q3*MG**(-2)*MZ4**2*yyy41*xxx3 - 2.*et1(q1,&
592 ! & q2,e1,e2)*q1_e3*e4_q3*MG**(-2)*MZ3**2*yyy41*xxx3 - 2.*et1(q1,&
593 ! & q2,e1,e2)*q1_e3*e4_q3*yyy41*xxx3 + 8.*et1(q1,q2,e1,e2)*q1_e4*&
594 ! & q1_q3*e3_q4*MG**(-2)*yyy42*xxx3 + 4.*et1(q1,q2,e1,e2)*q1_e4*&
595 ! & e3_q4*MG**(-2)*yyy1*xxx3 - 2.*et1(q1,q2,e1,e2)*q1_e4*e3_q4*&
596 ! & MG**(-2)*MZ4**2*yyy42*xxx3 + 2.*et1(q1,q2,e1,e2)*q1_e4*e3_q4*&
597 ! & MG**(-2)*MZ3**2*yyy42*xxx3 + 2.*et1(q1,q2,e1,e2)*q1_e4*e3_q4*&
598 ! & yyy42*xxx3 - 8.*et1(q1,q2,e1,e2)*q1_q3*e3_e4*MG**(-2)*MZ4**2*&
599 ! & yyy2*xxx3 + 8.*et1(q1,q2,e1,e2)*q1_q3*e3_e4*MG**(-2)*MZ3**2*&
600 ! & yyy2*xxx3 + 8.*et1(q1,q2,e1,e2)*q1_q3*e3_e4*yyy2*xxx3 + 4.*&
601 ! & et1(q1,q2,e1,e2)*q1_q3*e3_q4*e4_q3*MG**(-2)*yyy42*xxx3
602 ! res = res - 4.*et1(q1,q2,e1,e2)*q1_q3*e3_q4*e4_q3*MG**(-2)*yyy41*&
603 ! & xxx3 - 8.*et1(q1,q2,e1,e2)*q1_q3*e3_q4*e4_q3*MG**(-2)*MZ4**2*&
604 ! & yyy3*xxx3 + 8.*et1(q1,q2,e1,e2)*q1_q3*e3_q4*e4_q3*MG**(-2)*&
605 ! & MZ3**2*yyy3*xxx3 + 8.*et1(q1,q2,e1,e2)*q1_q3*e3_q4*e4_q3*yyy3&
606 ! & *xxx3 + 16.*et1(q1,q2,e1,e2)*q1_q3**2*e3_e4*MG**(-2)*yyy2*&
607 ! & xxx3 + 16.*et1(q1,q2,e1,e2)*q1_q3**2*e3_q4*e4_q3*MG**(-2)*&
608 ! & yyy3*xxx3 + 2./3.*et1(q1,q2,e1,e2)*e3_e4*MG**(-2)*MZ4**4*yyy2&
609 ! & *xxx3 - 4./3.*et1(q1,q2,e1,e2)*e3_e4*MG**(-2)*MZ3**2*MZ4**2*&
610 ! & yyy2*xxx3 + 2./3.*et1(q1,q2,e1,e2)*e3_e4*MG**(-2)*MZ3**4*yyy2&
611 ! & *xxx3 + 2./3.*et1(q1,q2,e1,e2)*e3_e4*yyy1*xxx3 - 4./3.*et1(q1&
612 ! & ,q2,e1,e2)*e3_e4*MZ4**2*yyy2*xxx3 + 8./3.*et1(q1,q2,e1,e2)*&
613 ! & e3_e4*MZ3**2*yyy2*xxx3 + 2./3.*et1(q1,q2,e1,e2)*e3_e4*MG**2*&
614 ! & yyy2*xxx3 + 4./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**(-2)*yyy1*&
615 ! & xxx3 - 2./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**(-2)*MZ4**2*&
616 ! & yyy42*xxx3 + 2./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**(-2)*&
617 ! & MZ4**2*yyy41*xxx3
618 ! res = res + 2./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**(-2)*MZ4**4*&
619 ! & yyy3*xxx3 + 2./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**(-2)*MZ3**2*&
620 ! & yyy42*xxx3 - 2./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**(-2)*&
621 ! & MZ3**2*yyy41*xxx3 - 4./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*&
622 ! & MG**(-2)*MZ3**2*MZ4**2*yyy3*xxx3 + 2./3.*et1(q1,q2,e1,e2)*&
623 ! & e3_q4*e4_q3*MG**(-2)*MZ3**4*yyy3*xxx3 + 4./3.*et1(q1,q2,e1,e2&
624 ! & )*e3_q4*e4_q3*yyy42*xxx3 - 2./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3&
625 ! & *yyy41*xxx3 - 4./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MZ4**2*yyy3*&
626 ! & xxx3 + 8./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MZ3**2*yyy3*xxx3 + &
627 ! & 2./3.*et1(q1,q2,e1,e2)*e3_q4*e4_q3*MG**2*yyy3*xxx3 - et1(q1,&
628 ! & q2,e1,e2)*et1(q1,q,e3,e4)*MG**(-2)*MZ4**2*yyy6*xxx3 + et1(q1,&
629 ! & q2,e1,e2)*et1(q1,q,e3,e4)*MG**(-2)*MZ3**2*yyy6*xxx3 + et1(q1,&
630 ! & q2,e1,e2)*et1(q1,q,e3,e4)*yyy6*xxx3 + 4.*et1(q1,q2,e1,e2)*&
631 ! & et1(q1,q,e3,e4)*q1_q3*MG**(-2)*yyy6*xxx3 - 4.*et1(q1,q2,e1,e2&
632 ! & )*et1(q1,q,e3,q3)*q1_q3*e4_q3*MG**(-4)*yyy7*xxx3 + et1(q1,q2,&
633 ! & e1,e2)*et1(q1,q,e3,q3)*e4_q3*MG**(-4)*MZ4**2*yyy7*xxx3
634 ! res = res - et1(q1,q2,e1,e2)*et1(q1,q,e3,q3)*e4_q3*MG**(-4)*&
635 ! & MZ3**2*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e3,q3)*e4_q3*&
636 ! & MG**(-2)*yyy7*xxx3 + 4.*et1(q1,q2,e1,e2)*et1(q1,q,e3,q4)*&
637 ! & q1_q3*e4_q3*MG**(-4)*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e3&
638 ! & ,q4)*e4_q3*MG**(-4)*MZ4**2*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(&
639 ! & q1,q,e3,q4)*e4_q3*MG**(-4)*MZ3**2*yyy7*xxx3 + et1(q1,q2,e1,e2&
640 ! & )*et1(q1,q,e3,q4)*e4_q3*MG**(-2)*yyy7*xxx3 - 4.*et1(q1,q2,e1,&
641 ! & e2)*et1(q1,q,e4,q3)*q1_q3*e3_q4*MG**(-4)*yyy7*xxx3 + et1(q1,&
642 ! & q2,e1,e2)*et1(q1,q,e4,q3)*e3_q4*MG**(-4)*MZ4**2*yyy7*xxx3 - &
643 ! & et1(q1,q2,e1,e2)*et1(q1,q,e4,q3)*e3_q4*MG**(-4)*MZ3**2*yyy7*&
644 ! & xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e4,q3)*e3_q4*MG**(-2)*yyy7*&
645 ! & xxx3 + 4.*et1(q1,q2,e1,e2)*et1(q1,q,e4,q4)*q1_q3*e3_q4*&
646 ! & MG**(-4)*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q1,q,e4,q4)*e3_q4*&
647 ! & MG**(-4)*MZ4**2*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q1,q,e4,q4)*&
648 ! & e3_q4*MG**(-4)*MZ3**2*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q1,q,&
649 ! & e4,q4)*e3_q4*MG**(-2)*yyy7*xxx3
650 ! res = res + et1(q1,q2,e1,e2)*et1(q2,q,e3,e4)*MG**(-2)*MZ4**2*yyy6&
651 ! & *xxx3 - et1(q1,q2,e1,e2)*et1(q2,q,e3,e4)*MG**(-2)*MZ3**2*yyy6*&
652 ! & xxx3 - et1(q1,q2,e1,e2)*et1(q2,q,e3,e4)*yyy6*xxx3 - 4.*et1(q1&
653 ! & ,q2,e1,e2)*et1(q2,q,e3,e4)*q1_q3*MG**(-2)*yyy6*xxx3 + 4.*et1(&
654 ! & q1,q2,e1,e2)*et1(q2,q,e3,q3)*q1_q3*e4_q3*MG**(-4)*yyy7*xxx3&
655 ! & - et1(q1,q2,e1,e2)*et1(q2,q,e3,q3)*e4_q3*MG**(-4)*MZ4**2*&
656 ! & yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q2,q,e3,q3)*e4_q3*MG**(-4)*&
657 ! & MZ3**2*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q2,q,e3,q3)*e4_q3*&
658 ! & MG**(-2)*yyy7*xxx3 - 4.*et1(q1,q2,e1,e2)*et1(q2,q,e3,q4)*&
659 ! & q1_q3*e4_q3*MG**(-4)*yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q2,q,e3&
660 ! & ,q4)*e4_q3*MG**(-4)*MZ4**2*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(&
661 ! & q2,q,e3,q4)*e4_q3*MG**(-4)*MZ3**2*yyy7*xxx3 - et1(q1,q2,e1,e2&
662 ! & )*et1(q2,q,e3,q4)*e4_q3*MG**(-2)*yyy7*xxx3 + 4.*et1(q1,q2,e1,&
663 ! & e2)*et1(q2,q,e4,q3)*q1_q3*e3_q4*MG**(-4)*yyy7*xxx3 - et1(q1,&
664 ! & q2,e1,e2)*et1(q2,q,e4,q3)*e3_q4*MG**(-4)*MZ4**2*yyy7*xxx3 + &
665 ! & et1(q1,q2,e1,e2)*et1(q2,q,e4,q3)*e3_q4*MG**(-4)*MZ3**2*yyy7*&
666 ! & xxx3
667 ! res = res + et1(q1,q2,e1,e2)*et1(q2,q,e4,q3)*e3_q4*MG**(-2)*yyy7*&
668 ! & xxx3 - 4.*et1(q1,q2,e1,e2)*et1(q2,q,e4,q4)*q1_q3*e3_q4*MG**(-4)*&
669 ! & yyy7*xxx3 + et1(q1,q2,e1,e2)*et1(q2,q,e4,q4)*e3_q4*MG**(-4)*&
670 ! & MZ4**2*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q2,q,e4,q4)*e3_q4*&
671 ! & MG**(-4)*MZ3**2*yyy7*xxx3 - et1(q1,q2,e1,e2)*et1(q2,q,e4,q4)*&
672 ! & e3_q4*MG**(-2)*yyy7*xxx3 - 1./3.*et1(q1,q2,e1,e2)*et1(q,e3,e4&
673 ! & ,q3)*yyy6*xxx3 + 1./3.*et1(q1,q2,e1,e2)*et1(q,e3,e4,q4)*yyy6*&
674 ! & xxx3 + 2./3.*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*MG**(-4)*&
675 ! & MZ4**4*yyy5*xxx3 - 4./3.*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*&
676 ! & MG**(-4)*MZ3**2*MZ4**2*yyy5*xxx3 + 2./3.*et1(q1,q2,e1,e2)*&
677 ! & et1(e3,e4,q3,q4)*MG**(-4)*MZ3**4*yyy5*xxx3 - 4./3.*et1(q1,q2,&
678 ! & e1,e2)*et1(e3,e4,q3,q4)*MG**(-2)*MZ4**2*yyy5*xxx3 + 8./3.*&
679 ! & et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*MG**(-2)*MZ3**2*yyy5*xxx3&
680 ! & + 2./3.*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*yyy5*xxx3 - 8.*&
681 ! & et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*q1_q3*MG**(-4)*MZ4**2*yyy5*&
682 ! & xxx3
683 ! res = res + 8.*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*q1_q3*MG**(-4)*&
684 ! & MZ3**2*yyy5*xxx3 + 8.*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*q1_q3*&
685 ! & MG**(-2)*yyy5*xxx3 + 16.*et1(q1,q2,e1,e2)*et1(e3,e4,q3,q4)*&
686 ! & q1_q3**2*MG**(-4)*yyy5*xxx3 + 4.*et1(q1,q,e3,e4)*q1_q3*e1_e2*&
687 ! & yyy6*xxx2 - et1(q1,q,e3,e4)*e1_e2*MZ4**2*yyy6*xxx2 + et1(q1,q&
688 ! & ,e3,e4)*e1_e2*MZ3**2*yyy6*xxx2 + et1(q1,q,e3,e4)*e1_e2*MG**2*&
689 ! & yyy6*xxx2 - 4.*et1(q1,q,e3,q3)*q1_q3*e1_e2*e4_q3*MG**(-2)*&
690 ! & yyy7*xxx2 + et1(q1,q,e3,q3)*e1_e2*e4_q3*MG**(-2)*MZ4**2*yyy7*&
691 ! & xxx2 - et1(q1,q,e3,q3)*e1_e2*e4_q3*MG**(-2)*MZ3**2*yyy7*xxx2&
692 ! & - et1(q1,q,e3,q3)*e1_e2*e4_q3*yyy7*xxx2 + 4.*et1(q1,q,e3,q4)&
693 ! & *q1_q3*e1_e2*e4_q3*MG**(-2)*yyy7*xxx2 - et1(q1,q,e3,q4)*e1_e2&
694 ! & *e4_q3*MG**(-2)*MZ4**2*yyy7*xxx2 + et1(q1,q,e3,q4)*e1_e2*&
695 ! & e4_q3*MG**(-2)*MZ3**2*yyy7*xxx2 + et1(q1,q,e3,q4)*e1_e2*e4_q3&
696 ! & *yyy7*xxx2 - 4.*et1(q1,q,e4,q3)*q1_q3*e1_e2*e3_q4*MG**(-2)*&
697 ! & yyy7*xxx2 + et1(q1,q,e4,q3)*e1_e2*e3_q4*MG**(-2)*MZ4**2*yyy7*&
698 ! & xxx2
699 ! res = res - et1(q1,q,e4,q3)*e1_e2*e3_q4*MG**(-2)*MZ3**2*yyy7*xxx2&
700 ! & - et1(q1,q,e4,q3)*e1_e2*e3_q4*yyy7*xxx2 + 4.*et1(q1,q,e4,q4)&
701 ! & *q1_q3*e1_e2*e3_q4*MG**(-2)*yyy7*xxx2 - et1(q1,q,e4,q4)*e1_e2&
702 ! & *e3_q4*MG**(-2)*MZ4**2*yyy7*xxx2 + et1(q1,q,e4,q4)*e1_e2*&
703 ! & e3_q4*MG**(-2)*MZ3**2*yyy7*xxx2 + et1(q1,q,e4,q4)*e1_e2*e3_q4&
704 ! & *yyy7*xxx2 - 4.*et1(q2,q,e3,e4)*q1_q3*e1_e2*yyy6*xxx2 + et1(&
705 ! & q2,q,e3,e4)*e1_e2*MZ4**2*yyy6*xxx2 - et1(q2,q,e3,e4)*e1_e2*&
706 ! & MZ3**2*yyy6*xxx2 - et1(q2,q,e3,e4)*e1_e2*MG**2*yyy6*xxx2 + 4.&
707 ! & *et1(q2,q,e3,q3)*q1_q3*e1_e2*e4_q3*MG**(-2)*yyy7*xxx2 - et1(&
708 ! & q2,q,e3,q3)*e1_e2*e4_q3*MG**(-2)*MZ4**2*yyy7*xxx2 + et1(q2,q,&
709 ! & e3,q3)*e1_e2*e4_q3*MG**(-2)*MZ3**2*yyy7*xxx2 + et1(q2,q,e3,q3&
710 ! & )*e1_e2*e4_q3*yyy7*xxx2 - 4.*et1(q2,q,e3,q4)*q1_q3*e1_e2*&
711 ! & e4_q3*MG**(-2)*yyy7*xxx2 + et1(q2,q,e3,q4)*e1_e2*e4_q3*&
712 ! & MG**(-2)*MZ4**2*yyy7*xxx2 - et1(q2,q,e3,q4)*e1_e2*e4_q3*&
713 ! & MG**(-2)*MZ3**2*yyy7*xxx2 - et1(q2,q,e3,q4)*e1_e2*e4_q3*yyy7*&
714 ! & xxx2
715 ! res = res + 4.*et1(q2,q,e4,q3)*q1_q3*e1_e2*e3_q4*MG**(-2)*yyy7*&
716 ! & xxx2 - et1(q2,q,e4,q3)*e1_e2*e3_q4*MG**(-2)*MZ4**2*yyy7*xxx2 + &
717 ! & et1(q2,q,e4,q3)*e1_e2*e3_q4*MG**(-2)*MZ3**2*yyy7*xxx2 + et1(&
718 ! & q2,q,e4,q3)*e1_e2*e3_q4*yyy7*xxx2 - 4.*et1(q2,q,e4,q4)*q1_q3*&
719 ! & e1_e2*e3_q4*MG**(-2)*yyy7*xxx2 + et1(q2,q,e4,q4)*e1_e2*e3_q4*&
720 ! & MG**(-2)*MZ4**2*yyy7*xxx2 - et1(q2,q,e4,q4)*e1_e2*e3_q4*&
721 ! & MG**(-2)*MZ3**2*yyy7*xxx2 - et1(q2,q,e4,q4)*e1_e2*e3_q4*yyy7*&
722 ! & xxx2 - et1(q,e1,e3,e4)*e2_q3*MG**2*yyy6*xxx1 + et1(q,e1,e3,q3&
723 ! & )*e2_q3*e4_q3*yyy7*xxx1 - et1(q,e1,e3,q4)*e2_q3*e4_q3*yyy7*&
724 ! & xxx1 + et1(q,e1,e4,q3)*e2_q3*e3_q4*yyy7*xxx1 - et1(q,e1,e4,q4&
725 ! & )*e2_q3*e3_q4*yyy7*xxx1 - et1(q,e2,e3,e4)*e1_q3*MG**2*yyy6*&
726 ! & xxx1 + et1(q,e2,e3,q3)*e1_q3*e4_q3*yyy7*xxx1 - et1(q,e2,e3,q4&
727 ! & )*e1_q3*e4_q3*yyy7*xxx1 + et1(q,e2,e4,q3)*e1_q3*e3_q4*yyy7*&
728 ! & xxx1 - et1(q,e2,e4,q4)*e1_q3*e3_q4*yyy7*xxx1 - 1./3.*et1(q,e3&
729 ! & ,e4,q3)*e1_e2*MG**2*yyy6*xxx2 + 1./3.*et1(q,e3,e4,q3)*e1_e2*&
730 ! & MG**2*yyy6*xxx1
731 ! res = res + 1./3.*et1(q,e3,e4,q4)*e1_e2*MG**2*yyy6*xxx2 - 1./3.*&
732 ! & et1(q,e3,e4,q4)*e1_e2*MG**2*yyy6*xxx1 - 8.*et1(e3,e4,q3,q4)*&
733 ! & q1_q3*e1_e2*MG**(-2)*MZ4**2*yyy5*xxx2 + 8.*et1(e3,e4,q3,q4)*&
734 ! & q1_q3*e1_e2*MG**(-2)*MZ3**2*yyy5*xxx2 + 8.*et1(e3,e4,q3,q4)*&
735 ! & q1_q3*e1_e2*yyy5*xxx2 + 16.*et1(e3,e4,q3,q4)*q1_q3**2*e1_e2*&
736 ! & MG**(-2)*yyy5*xxx2 + 2./3.*et1(e3,e4,q3,q4)*e1_e2*MG**(-2)*&
737 ! & MZ4**4*yyy5*xxx2 + 1./3.*et1(e3,e4,q3,q4)*e1_e2*MG**(-2)*&
738 ! & MZ4**4*yyy5*xxx1 - 4./3.*et1(e3,e4,q3,q4)*e1_e2*MG**(-2)*&
739 ! & MZ3**2*MZ4**2*yyy5*xxx2 - 2./3.*et1(e3,e4,q3,q4)*e1_e2*&
740 ! & MG**(-2)*MZ3**2*MZ4**2*yyy5*xxx1 + 2./3.*et1(e3,e4,q3,q4)*&
741 ! & e1_e2*MG**(-2)*MZ3**4*yyy5*xxx2 + 1./3.*et1(e3,e4,q3,q4)*&
742 ! & e1_e2*MG**(-2)*MZ3**4*yyy5*xxx1 - 4./3.*et1(e3,e4,q3,q4)*&
743 ! & e1_e2*MZ4**2*yyy5*xxx2 - 2./3.*et1(e3,e4,q3,q4)*e1_e2*MZ4**2*&
744 ! & yyy5*xxx1 + 8./3.*et1(e3,e4,q3,q4)*e1_e2*MZ3**2*yyy5*xxx2 - 2.&
745 ! & /3.*et1(e3,e4,q3,q4)*e1_e2*MZ3**2*yyy5*xxx1 + 2./3.*et1(e3,e4&
746 ! & ,q3,q4)*e1_e2*MG**2*yyy5*xxx2
747 ! res = res + 1./3.*et1(e3,e4,q3,q4)*e1_e2*MG**2*yyy5*xxx1 + 4.*&
748 ! & et1(e3,e4,q3,q4)*e1_q3*e2_q3*yyy5*xxx1
749 
750 
751 ! print *, "new ",res
752 
753 
754 ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! SAME AS ABOVE BUT SHORTER ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! !
755 
756  abr1 = mg**2 + mz3**2 - mz4**2 + 4*q1_q3
757 
758  res = mg**2*xxx1*(e1_e4*(yyy1*e2_e3 + yyy42*e2_q3*e3_q4) + e1_e3*(yyy1*e2_e4 - yyy41*e2_q3*e4_q3) + &
759  e1_q3*(yyy42*e2_e4*e3_q4 - yyy41*e2_e3*e4_q3 + 4*e2_q3*(yyy2*e3_e4 + yyy3*e3_q4*e4_q3))) + &
760  e1_e2*((xxx1*((mg**4*yyy2 + (mz3**2 - mz4**2)**2*yyy2 - 2*mg**2*(yyy1 + (mz3**2 + mz4**2)*yyy2))*e3_e4 + &
761  (2*yyy1 + mg**4*yyy3 + (mz3 - mz4)*(mz3 + mz4)*((mz3 - mz4)*(mz3 + mz4)*yyy3 - yyy41 + yyy42) - mg**2*(2*(mz3**2 + mz4**2)*yyy3 + yyy41 + yyy42))* &
762  e3_q4*e4_q3))/3d0 + (2*xxx2*(3*q1_e3*((2*yyy1 - abr1*yyy41)*e4_q3 + 4*yyy1*q1_e4) + &
763  e3_e4*(abr1**2*yyy2 + abr1*(yyy1 + 2*mz3**2*yyy2) - (mz3 - mz4)*(mz3 + mz4)*(yyy1 + 2*mz3**2*yyy2) - &
764  4*q1_q3*(yyy1 - abr1*yyy2 + 2*mz3**2*yyy2 + 8*yyy2*q1_q3) + 24*yyy2*q1_q3**2) + &
765  e3_q4*(3*(2*yyy1 + abr1*yyy42)*q1_e4 + e4_q3* &
766  (2*yyy1 + abr1**2*yyy3 - (mz3 - mz4)*(mz3 + mz4)*(2*mz3**2*yyy3 + yyy42) + abr1*(2*mz3**2*yyy3 - yyy41 + 2*yyy42 ) - &
767  2*q1_q3*(-2*abr1*yyy3 + 4*mz3**2*yyy3 + yyy41 + yyy42 + 16*yyy3*q1_q3) + 24*yyy3*q1_q3**2))))/3d0) - &
768  mg**2*xxx1*yyy6*e2_q3*et1(q,e1,e3,e4) + xxx1*yyy7*e2_q3*e4_q3*et1(q,e1,e3,q3) - xxx1*yyy7*e2_q3*e4_q3*et1(q,e1,e3,q4) + &
769  xxx1*yyy7*e2_q3*e3_q4*et1(q,e1,e4,q3) - xxx1*yyy7*e2_q3*e3_q4*et1(q,e1,e4,q4) - mg**2*xxx1*yyy6*e1_q3*et1(q,e2,e3,e4) + &
770  xxx1*yyy7*e1_q3*e4_q3*et1(q,e2,e3,q3) - xxx1*yyy7*e1_q3*e4_q3*et1(q,e2,e3,q4) + xxx1*yyy7*e1_q3*e3_q4*et1(q,e2,e4,q3) - &
771  xxx1*yyy7*e1_q3*e3_q4*et1(q,e2,e4,q4) + et1(q,e3,e4,q3)* &
772  (((mg**2*xxx1*yyy6)/3d0 - (mg**2*xxx2*yyy6)/3d0)*e1_e2 - (xxx3*yyy6*et1(q1,q2,e1,e2))/3d0) + &
773  et1(q,e3,e4,q4)*((-(mg**2*xxx1*yyy6)/3d0 + (mg**2*xxx2*yyy6)/3d0)*e1_e2 + (xxx3*yyy6*et1(q1,q2,e1,e2))/3d0) + &
774  et1(q1,q,e3,e4)*(abr1*xxx2*yyy6*e1_e2 + (abr1*xxx3*yyy6*et1(q1,q2,e1,e2))/mg**2) + &
775  et1(q1,q,e4,q3)*(-((abr1*xxx2*yyy7*e1_e2*e3_q4)/mg**2) - (abr1*xxx3*yyy7*e3_q4*et1(q1,q2,e1,e2))/mg**4) + &
776  et1(q1,q,e4,q4)*((abr1*xxx2*yyy7*e1_e2*e3_q4)/mg**2 + (abr1*xxx3*yyy7*e3_q4*et1(q1,q2,e1,e2))/mg**4) + &
777  et1(q1,q,e3,q3)*(-((abr1*xxx2*yyy7*e1_e2*e4_q3)/mg**2) - (abr1*xxx3*yyy7*e4_q3*et1(q1,q2,e1,e2))/mg**4) + &
778  et1(q1,q,e3,q4)*((abr1*xxx2*yyy7*e1_e2*e4_q3)/mg**2 + (abr1*xxx3*yyy7*e4_q3*et1(q1,q2,e1,e2))/mg**4) + &
779  et1(e3,e4,q3,q4)*(4*xxx1*yyy5*e1_q3*e2_q3 + e1_e2* &
780  (((mg - mz3 - mz4)*(mg + mz3 - mz4)*(mg - mz3 + mz4)*(mg + mz3 + mz4)*xxx1*yyy5)/(3d0*mg**2) + &
781  (2*xxx2*yyy5*(abr1**2 + 2*abr1*mz3**2 - 2*mz3**4 + 2*mz3**2*mz4**2 + 4*(abr1 - 2*mz3**2 - 8*q1_q3)*q1_q3 + 24*q1_q3**2))/(3d0*mg**2)) + &
782  (2*xxx3*yyy5*(abr1**2 + 2*abr1*mz3**2 - 2*mz3**4 + 2*mz3**2*mz4**2 + 4*(abr1 - 2*mz3**2 - 8*q1_q3)*q1_q3 + 24*q1_q3**2)*et1(q1,q2,e1,e2))/ &
783  (3d0*mg**4)) - abr1*xxx2*yyy6*e1_e2*et1(q2,q,e3,e4) + (abr1*xxx2*yyy7*e1_e2*e4_q3*et1(q2,q,e3,q3))/mg**2 - &
784  (abr1*xxx2*yyy7*e1_e2*e4_q3*et1(q2,q,e3,q4))/mg**2 + (abr1*xxx2*yyy7*e1_e2*e3_q4*et1(q2,q,e4,q3))/mg**2 - &
785  (abr1*xxx2*yyy7*e1_e2*e3_q4*et1(q2,q,e4,q4))/mg**2 + &
786  et1(q1,q2,e1,e2)*((2*xxx3*(3*q1_e3*((2*yyy1 - abr1*yyy41)*e4_q3 + 4*yyy1*q1_e4) + &
787  e3_e4*(abr1**2*yyy2 + abr1*(yyy1 + 2*mz3**2*yyy2) - (mz3 - mz4)*(mz3 + mz4)*(yyy1 + 2*mz3**2*yyy2) - &
788  4*q1_q3*(yyy1 - abr1*yyy2 + 2*mz3**2*yyy2 + 8*yyy2*q1_q3) + 24*yyy2*q1_q3**2) + &
789  e3_q4*(3*(2*yyy1 + abr1*yyy42)*q1_e4 + e4_q3* &
790  (2*yyy1 + abr1**2*yyy3 - (mz3 - mz4)*(mz3 + mz4)*(2*mz3**2*yyy3 + yyy42) + abr1*(2*mz3**2*yyy3 - yyy41 + 2*yyy42) - &
791  2*q1_q3*(-2*abr1*yyy3 + 4*mz3**2*yyy3 + yyy41 + yyy42 + 16*yyy3*q1_q3) + 24*yyy3*q1_q3**2))))/(3d0*mg**2) - &
792  (abr1*xxx3*yyy6*et1(q2,q,e3,e4))/mg**2 + (abr1*xxx3*yyy7*e4_q3*et1(q2,q,e3,q3))/mg**4 - (abr1*xxx3*yyy7*e4_q3*et1(q2,q,e3,q4))/mg**4 + &
793  (abr1*xxx3*yyy7*e3_q4*et1(q2,q,e4,q3))/mg**4 - (abr1*xxx3*yyy7*e3_q4*et1(q2,q,e4,q4))/mg**4)
794 
795 
796 ! print *, "newer",res
797 ! pause
798 
799 ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! !
800 
801  else
802 
803 
804 
805  res = q1_e3*q2_e4*e1_e2 * ( m_reso**2 )
806  res = res + q1_e3*q2_q4*e1_e2*e4_q3 * ( - 2.0 )
807  res = res + q1_e4*q2_e3*e1_e2 * ( m_reso**2 )
808  res = res + q1_e4*q2_q3*e1_e2*e3_q4 * ( - 2.0 )
809  res = res + q1_q3*q2_e4*e1_e2*e3_q4 * ( - 2.0 )
810  res = res + q1_q3*q2_q4*e1_e2*e3_e4 * ( 2.0 )
811  res = res + q1_q3*e1_e2*e3_e4 * ( 1./2.*m_reso**2 )
812  res = res + q1_q3*e1_e2*e3_q4*e4_q3 * ( - 1. )
813  res = res + q1_q4*q2_e3*e1_e2*e4_q3 * ( - 2. )
814  res = res + q1_q4*q2_q3*e1_e2*e3_e4 * ( 2 )
815  res = res + q1_q4*e1_e2*e3_e4 * ( 1./2.*m_reso**2 )
816  res = res + q1_q4*e1_e2*e3_q4*e4_q3 * ( - 1. )
817  res = res + q2_q3*e1_e2*e3_e4 * ( 1./2.*m_reso**2 )
818  res = res + q2_q3*e1_e2*e3_q4*e4_q3 * ( - 1 )
819  res = res + q2_q4*e1_e2*e3_e4 * ( 1./2.*m_reso**2 )
820  res = res + q2_q4*e1_e2*e3_q4*e4_q3 * ( - 1. )
821  res = res + e1_e2*e3_e4 * ( m_reso**2*m_v**2 - 1./2.*m_reso**4 )
822  res = res + e1_e2*e3_q4*e4_q3 * ( m_reso**2 )
823  res = res + e1_e3*e2_e4 * ( 1./2.*m_reso**4 )
824  res = res + e1_e3*e2_q4*e4_q3 * ( - m_reso**2 )
825  res = res + e1_e4*e2_e3 * ( 1./2.*m_reso**4 )
826  res = res + e1_e4*e2_q3*e3_q4 * ( - m_reso**2 )
827  res = res + e1_q3*e2_e4*e3_q4 * ( - m_reso**2 )
828  res = res + e1_q3*e2_q4*e3_e4 * ( m_reso**2 )
829  res = res + e1_q4*e2_e3*e4_q3 * ( - m_reso**2 )
830  res = res + e1_q4*e2_q3*e3_e4 * ( m_reso**2 )
831 
832  print *,"this code should no longer be used"; stop 1
833 
834 !print *, "res GG old",res
835 !pause
836 
837 
838  endif
839 
840  end subroutine gggzzampl
841 
842 
843 !----- a subroutine for q qbar -> G -> Z -> lept + Z --> 2 lepts
844 !----- all outgoing convention and the following momentum assignment
845 !----- 0 -> bq(p1) + q(p2) + e-(p3) + e+(p4) +mu-(p5) +mu+(p6)
846  subroutine evalamp_qqb_g_vv(p,MY_IDUP,res)
847  implicit none
848  real(dp), intent(out) :: res
849  real(dp), intent(in) :: p(4,6)
850  integer, intent(in) :: MY_IDUP(6:9)
851  real(dp) :: pin(4,4)
852  complex(dp) :: A_VV(1:18), A0_VV(1:2)
853  integer :: i1,i3,i4,VVMode,VVmode_swap
854  real(dp) :: prefactor
855  real(dp) :: intcolfac
856  integer :: ordering(1:4),ordering_swap(1:4)
857  logical :: doInterference
858 
859  if(isaquark(my_idup(6)) .and. isaquark(my_idup(8))) then
860  intcolfac=1.0_dp/3.0_dp
861  else
862  intcolfac=1.0_dp
863  endif
864 
865  call getdecay_vvmode_ordering(my_idup(6:9),vvmode,ordering,vvmode_swap,ordering_swap)
866 
867 !---- full prefactor; 3 is the color factor
868  if( vvmode.eq.zzmode ) then! Z decay
869  prefactor = 3d0*overallcouplvffsq**2
870  elseif( vvmode.eq.wwmode ) then ! W decay
871  prefactor = 3d0*overallcouplvffsq**2
872  elseif( vvmode.eq.zgmode ) then ! Z+photon "decay"
873  prefactor = 3d0*overallcouplvffsq ! Only single powers
874  elseif( vvmode.eq.ggmode ) then ! photon "decay"
875  prefactor = 3d0
876  else
877  prefactor = 0d0
878  endif
879 
880 
881  res = zero
882  a_vv(:) = 0d0
883  dointerference = includeinterference .and. ( &
884  ((vvmode.eq.zzmode) .and. (vvmode_swap.eq.zzmode)) &
885  )
886  if ( includevprime .and. .not.(vvmode.eq.zzmode .or. vvmode.eq.zgmode .or. vvmode.eq.wwmode) ) then
887  call error("Contact terms only for WW, ZZ or Zg!")
888  endif
889  do i1=1,2; do i3=1,2; do i4=1,2! sum over helicities
890  call calchelamp_qq(ordering,vvmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(1))
891  if( vvmode.eq.zzmode ) then
892  if( includegammastar ) then
893  call calchelamp_qq(ordering,zgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(3))
894  call calchelamp_qq(ordering,gszmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(5))
895  call calchelamp_qq(ordering,gsgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(7))
896  endif
897  if( includevprime ) then
898  call calchelamp_qq(ordering,zzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(9))
899  call calchelamp_qq(ordering,zpzmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(11))
900  call calchelamp_qq(ordering,zpzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(13))
901  endif
902  if( includegammastar .and. includevprime ) then
903  call calchelamp_qq(ordering,gszpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(15))
904  call calchelamp_qq(ordering,zpgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(17))
905  endif
906  elseif( vvmode.eq.zgmode ) then
907  if(includegammastar) then
908  call calchelamp_qq(ordering,gsgmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(3))
909  endif
910  if( includevprime ) then
911  call calchelamp_qq(ordering,zpgmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(5))
912  endif
913  elseif( vvmode.eq.wwmode .and. includevprime ) then
914  call calchelamp_qq(ordering,wwpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(9))
915  call calchelamp_qq(ordering,wpwmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(11))
916  call calchelamp_qq(ordering,wpwpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(13))
917  endif
918 
919  if( dointerference ) then
920  call calchelamp_qq(ordering_swap,vvmode_swap,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(2))
921  if( includegammastar ) then
922  call calchelamp_qq(ordering_swap,zgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(4))
923  call calchelamp_qq(ordering_swap,gszmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(6))
924  call calchelamp_qq(ordering_swap,gsgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(8))
925  endif
926  if( includevprime ) then
927  call calchelamp_qq(ordering_swap,zzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(10))
928  call calchelamp_qq(ordering_swap,zpzmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(12))
929  call calchelamp_qq(ordering_swap,zpzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(14))
930  endif
931  if( includegammastar .and. includevprime ) then
932  call calchelamp_qq(ordering_swap,gszpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(16))
933  call calchelamp_qq(ordering_swap,zpgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(18))
934  endif
935  endif
936 
937  a0_vv(1) = a_vv(1)+a_vv(3)+a_vv(5)+a_vv(7)+a_vv(9)+a_vv(11)+a_vv(13)+a_vv(15)+a_vv(17) ! 3456 pieces
938  a0_vv(2) = a_vv(2)+a_vv(4)+a_vv(6)+a_vv(8)+a_vv(10)+a_vv(12)+a_vv(14)+a_vv(16)+a_vv(18) ! 5436 pieces
939  res = res + dreal(a0_vv(1)*dconjg(a0_vv(1)))
940  res = res + dreal(a0_vv(2)*dconjg(a0_vv(2)))
941  if( dointerference .and. (i3.eq.i4) ) then! interfere the 3456 with 5436 pieces
942  res = res - 2d0*intcolfac*dreal( a0_vv(1)*dconjg(a0_vv(2)) ) ! minus from Fermi statistics
943  endif
944  enddo; enddo; enddo
945 
946  res = res*prefactor
947  if( (vvmode.eq.zzmode) .and. dointerference ) res = res * symmfac
948 
949  end subroutine
950 
951  subroutine calchelamp_qq(ordering,VVMode,p,MY_IDUP,i1,i3,i4,A)
952  implicit none
953  integer :: ordering(1:4),i1,i3,i4,l1,l2,l3,l4,MY_IDUP(6:9),VVMode
954  real(dp) :: p(1:4,1:6)
955  complex(dp) :: propG, propZ1, propZ2
956  real(dp) :: s, pin(4,4)
957  complex(dp) :: A(1:1), sp(4,4)
958 
959  l1=ordering(1)
960  l2=ordering(2)
961  l3=ordering(3)
962  l4=ordering(4)
963 
964  !print *,"p(6)=",(p(:,l1))
965  !print *,"p(7)=",(p(:,l2))
966  !print *,"p(8)=",(p(:,l3))
967  !print *,"p(9)=",(p(:,l4))
968  !pause
969 
970  s = 2d0*scr(p(:,1),p(:,2))
971  propg = s/dcmplx(s - m_reso**2,m_reso*ga_reso)
972  pin(1,:) = p(:,1)
973  pin(2,:) = p(:,2)
974  sp(1,:) = pol_dk2mom(dcmplx(p(:,2)),dcmplx(p(:,1)),-3+2*i1) !qbq
975  sp(2,:) = sp(1,:) !-- the same, isn't really needed but for uniform bookeeping
976  call getdecay_couplings_spinors_props( &
977  vvmode, &
978  (/my_idup(l1+3),my_idup(l2+3),my_idup(l3+3),my_idup(l4+3)/), &
979  (/p(:,l1),p(:,l2),p(:,l3),p(:,l4)/), &
980  -3+2*i3,-3+2*i4, &
981  sp(3:4,:),pin(3:4,:) &
982  )
983  call qqgzzampl(vvmode,pin,sp,a(1))
984 !---- chiral couplings of quarks to gravitons
985  if (i1.eq.1) then
986  a(1) = graviton_qq_left*a(1)
987  elseif(i1.eq.2) then
988  a(1) = graviton_qq_right*a(1)
989  endif
990  a(1) = a(1)*propg
991  end subroutine
992 
993  subroutine qqgzzampl(VVMode,p,sp,res)
994  implicit none
995  integer, intent(in) :: VVMode
996  real(dp), intent(in) :: p(4,4)
997  complex(dp), intent(in) :: sp(4,4)
998  complex(dp), intent(out) :: res
999  complex(dp) :: e1_e2, e1_e3, e1_e4
1000  complex(dp) :: e2_e3, e2_e4
1001  complex(dp) :: e3_e4
1002  complex(dp) :: q_q,q3_q3,q4_q4
1003  complex(dp) :: q1_q2,q1_q3,q1_q4
1004  complex(dp) :: q2_q3,q2_q4
1005  complex(dp) :: q3_q4
1006  complex(dp) :: q1_e3,q1_e4,q2_e3,q2_e4
1007  complex(dp) :: e1_q3,e1_q4,e2_q3,e2_q4
1008  complex(dp) :: e3_q4,e4_q3
1009  complex(dp) :: q1(4),q2(4),q3(4),q4(4),q(4)
1010  complex(dp) :: e1(4),e2(4),e3(4),e4(4),abr1
1011  complex(dp) :: yyy1,yyy2,yyy3,yyy41,yyy42,yyy5,yyy6,yyy7,yyy4
1012  complex(dp) :: b_dyn(1:10)
1013  real(dp) :: q34,MG,MZ3,MZ4
1014  real(dp) :: rr
1015 
1016  b_dyn(:)=czero
1017 
1018  q1 = dcmplx(p(1,:),0d0)
1019  q2 = dcmplx(p(2,:),0d0)
1020  q3 = dcmplx(p(3,:),0d0)
1021  q4 = dcmplx(p(4,:),0d0)
1022 
1023 
1024  e1 = sp(1,:)
1025  e2 = sp(2,:)
1026  e3 = sp(3,:)
1027  e4 = sp(4,:)
1028 
1029  q = -q1-q2
1030 
1031  q_q =sc(q,q)
1032  q3_q3 = sc(q3,q3)
1033  q4_q4 = sc(q4,q4)
1034 
1035 
1036  q1_q2 = sc(q1,q2)
1037  q1_q3 = sc(q1,q3)
1038  q1_q4 = sc(q1,q4)
1039  q2_q3 = sc(q2,q3)
1040  q2_q4 = sc(q2,q4)
1041  q3_q4 = sc(q3,q4)
1042 
1043  e1_e2 = sc(e1,e2)
1044  e1_e3 = sc(e1,e3)
1045  e1_e4 = sc(e1,e4)
1046 
1047  e2_e3 = sc(e2,e3)
1048  e2_e4 = sc(e2,e4)
1049 
1050  e3_e4 = sc(e3,e4)
1051 
1052 
1053  q1_e3 = sc(q1,e3)
1054  q1_e4 = sc(q1,e4)
1055  q2_e3 = sc(q2,e3)
1056  q2_e4 = sc(q2,e4)
1057  e1_q3 = sc(e1,q3)
1058  e1_q4 = sc(e1,q4)
1059  e2_q3 = sc(e2,q3)
1060  e2_q4 = sc(e2,q4)
1061  e3_q4 = sc(e3,q4)
1062  e4_q3 = sc(e4,q3)
1063 
1064  mz3=dsqrt(cdabs(q3_q3))
1065  mz4=dsqrt(cdabs(q4_q4))
1066  if( use_dynamic_mg ) then
1067  mg = dsqrt(cdabs(q_q))
1068  else
1069  mg = m_reso
1070  endif
1071 
1072  q34 = (mg**2-mz3**2-mz4**2)/2d0
1073 
1074  if (generate_bis) then
1075  rr = q34/lambda**2! kappa for FS
1076 
1077  if( (vvmode.eq.zzmode) .or. (vvmode.eq.wwmode) ) then! decay ZZ's or WW's
1078  b_dyn(1)=b1
1079  b_dyn(2)=b2
1080  b_dyn(3)=b3
1081  b_dyn(4)=b4
1082  b_dyn(5)=b5
1083  b_dyn(6)=b6
1084  b_dyn(7)=b7
1085  b_dyn(8)=b8
1086  b_dyn(9)=b9
1087  b_dyn(10)=b10
1088  elseif( (vvmode.eq.zgmode) .OR. (vvmode.eq.gszmode) .OR. (vvmode.eq.zgsmode) ) then
1089  b_dyn(1)=bzgs1
1090  b_dyn(2)=bzgs2
1091  b_dyn(3)=bzgs3
1092  b_dyn(4)=bzgs4
1093  b_dyn(8)=bzgs8
1094  elseif( (vvmode.eq.ggmode) .or. (vvmode.eq.gsgsmode) .or. (vvmode.eq.gsgmode) ) then
1095  b_dyn(1)=bgsgs1
1096  b_dyn(2)=bgsgs2
1097  b_dyn(3)=bgsgs3
1098  b_dyn(4)=bgsgs4
1099  b_dyn(8)=bgsgs8
1100  elseif( (vvmode.eq.zzpmode) .or. (vvmode.eq.wwpmode) .or. (vvmode.eq.zpzmode) .or. (vvmode.eq.wpwmode) ) then
1101  b_dyn(1)=bzzp1
1102  b_dyn(2)=bzzp2
1103  b_dyn(3)=bzzp3
1104  b_dyn(4)=bzzp4
1105  b_dyn(5)=bzzp5
1106  b_dyn(6)=bzzp6
1107  b_dyn(7)=bzzp7
1108  b_dyn(8)=bzzp8
1109  b_dyn(9)=bzzp9
1110  b_dyn(10)=bzzp10
1111  elseif( (vvmode.eq.zpzpmode) .or. (vvmode.eq.wpwpmode) ) then
1112  b_dyn(1)=bzpzp1
1113  b_dyn(2)=bzpzp2
1114  b_dyn(3)=bzpzp3
1115  b_dyn(4)=bzpzp4
1116  b_dyn(5)=bzpzp5
1117  b_dyn(6)=bzpzp6
1118  b_dyn(7)=bzpzp7
1119  b_dyn(8)=bzpzp8
1120  b_dyn(9)=bzpzp9
1121  b_dyn(10)=bzpzp10
1122  elseif( (vvmode.eq.zpgmode) .OR. (vvmode.eq.gszpmode) .OR. (vvmode.eq.zpgsmode) ) then
1123  b_dyn(1)=bzpgs1
1124  b_dyn(2)=bzpgs2
1125  b_dyn(3)=bzpgs3
1126  b_dyn(4)=bzpgs4
1127  b_dyn(8)=bzpgs8
1128  else
1129  print *,"VVMode",vvmode,"not implemented"
1130  endif
1131 
1132  yyy1 = q34*( b_dyn(1) + b_dyn(2)*rr*(one+mz3**2/q34)*(one+mz4**2/q34) ) + b_dyn(5)*m_v**2
1133  yyy2 = -b_dyn(1)/two + b_dyn(3)*rr*(1d0-(mz3**2+mz4**2)/(2d0*q34)) + two*b_dyn(4)*rr + b_dyn(7)*rr*m_v**2/q34
1134  yyy3 = (-b_dyn(2)/two - b_dyn(3)- two*b_dyn(4))*rr/q34
1135  yyy41 = -b_dyn(1) - b_dyn(2)*(q34+mz3**2)/lambda**2 - b_dyn(3)*mz4**2/lambda**2 - 2d0*b_dyn(6)*m_v**2/lambda**2
1136  yyy42 = -b_dyn(1) - b_dyn(2)*(q34+mz4**2)/lambda**2 - b_dyn(3)*mz3**2/lambda**2 - 2d0*b_dyn(6)*m_v**2/lambda**2
1137  yyy5 = two*b_dyn(8)*rr*mg**2/q34
1138  yyy6 = b_dyn(9) * m_v**2/lambda**2
1139  yyy7 = b_dyn(10) * mg**2 * m_v**2/lambda**4
1140  else
1141  yyy1 = q34*c1/2d0
1142  yyy2 = c2
1143  yyy3 = c3/mg**2
1144  yyy41 = c41
1145  yyy42 = c42
1146  yyy5 = c5
1147  yyy6 = czero
1148  yyy7 = czero
1149  if(vvmode.eq.zzmode .or. vvmode.eq.wwmode) then
1150  yyy6 = c6
1151  yyy7 = c7
1152  endif
1153  endif
1154 
1155  res = czero
1156 
1157 ! old code without c41 c42 coupligs
1158 ! res = &
1159 ! q1_e3*e1_e4*yyy1 - q1_e3*e1_q3*e4_q3*yyy4 + q1_e4*e1_e3*yyy1 + &
1160 ! q1_e4*e1_q3*e3_q4*yyy4 - q1_q3*e1_e3*e4_q3*yyy4 + q1_q3*e1_e4* &
1161 ! e3_q4*yyy4 + 4*q1_q3*e1_q3*e3_e4*yyy2 + 4*q1_q3*e1_q3*e3_q4* &
1162 ! e4_q3*yyy3 + 1./2.*e1_e3*e4_q3*yyy1 - 1./4.*e1_e3*e4_q3*MG**2* &
1163 ! yyy4 + 1./2.*e1_e4*e3_q4*yyy1 + 1./4.*e1_e4*e3_q4*MG**2*yyy4 + &
1164 ! e1_q3*e3_e4*MG**2*yyy2 + e1_q3*e3_q4*e4_q3*MG**2*yyy3
1165 !
1166 ! res = res + &
1167 ! 1./2.*et1(e3,e4,q1,q)*e1_q3*yyy6 - 1./2.*et1(e3,e4,q2,q)*e1_q3* &
1168 ! yyy6 + 1./4.*et1(e3,e4,e1,q)*MG**2*yyy6 + et1(e3,e4,e1,q)*q1_q3* &
1169 ! yyy6 + 4*et1(e3,e4,q3,q4)*q1_q3*e1_q3*MG**(-2)*yyy5 + et1(e3,e4, &
1170 ! q3,q4)*e1_q3*yyy5
1171 !
1172 ! res = res + &
1173 ! 1./2.*et1(q1,e3,q,q3)*e1_q3*e4_q3*MG**(-2)*yyy7 - 1./2.*et1(q1, &
1174 ! e3,q,q4)*e1_q3*e4_q3*MG**(-2)*yyy7 + 1./2.*et1(q1,e4,q,q3)*e1_q3 &
1175 ! *e3_q4*MG**(-2)*yyy7 - 1./2.*et1(q1,e4,q,q4)*e1_q3*e3_q4* &
1176 ! MG**(-2)*yyy7 - 1./2.*et1(q2,e3,q,q3)*e1_q3*e4_q3*MG**(-2)*yyy7 &
1177 ! + 1./2.*et1(q2,e3,q,q4)*e1_q3*e4_q3*MG**(-2)*yyy7 - 1./2.*et1( &
1178 ! q2,e4,q,q3)*e1_q3*e3_q4*MG**(-2)*yyy7 + 1./2.*et1(q2,e4,q,q4)* &
1179 ! e1_q3*e3_q4*MG**(-2)*yyy7 + et1(e1,e3,q,q3)*q1_q3*e4_q3*MG**(-2) &
1180 ! *yyy7 + 1./4.*et1(e1,e3,q,q3)*e4_q3*yyy7 - et1(e1,e3,q,q4)*q1_q3 &
1181 ! *e4_q3*MG**(-2)*yyy7 - 1./4.*et1(e1,e3,q,q4)*e4_q3*yyy7 + et1(e1 &
1182 ! ,e4,q,q3)*q1_q3*e3_q4*MG**(-2)*yyy7 + 1./4.*et1(e1,e4,q,q3)* &
1183 ! e3_q4*yyy7 - et1(e1,e4,q,q4)*q1_q3*e3_q4*MG**(-2)*yyy7 - 1./4.* &
1184 ! et1(e1,e4,q,q4)*e3_q4*yyy7
1185 ! print *, "res old QQB",res
1186 
1187 
1188 
1189 ! ! new code with c41 c42 coupligs
1190 ! res =&
1191 ! & q1_e3*e1_e4*yyy1 - q1_e3*e1_q3*e4_q3*yyy41 + q1_e4*e1_e3*yyy1 + &
1192 ! & q1_e4*e1_q3*e3_q4*yyy42 - q1_q3*e1_e3*e4_q3*yyy41 + q1_q3*e1_e4*&
1193 ! & e3_q4*yyy42 + 4.*q1_q3*e1_q3*e3_e4*yyy2 + 4.*q1_q3*e1_q3*e3_q4*&
1194 ! & e4_q3*yyy3 + 1./2.*e1_e3*e4_q3*yyy1 + 1./4.*e1_e3*e4_q3*MZ4**2*&
1195 ! & yyy41 - 1./4.*e1_e3*e4_q3*MZ3**2*yyy41 - 1./4.*e1_e3*e4_q3*MG**2&
1196 ! & *yyy41 + 1./2.*e1_e4*e3_q4*yyy1 - 1./4.*e1_e4*e3_q4*MZ4**2*yyy42&
1197 ! & + 1./4.*e1_e4*e3_q4*MZ3**2*yyy42 + 1./4.*e1_e4*e3_q4*MG**2*&
1198 ! & yyy42 - e1_q3*e3_e4*MZ4**2*yyy2 + e1_q3*e3_e4*MZ3**2*yyy2 + &
1199 ! & e1_q3*e3_e4*MG**2*yyy2 + 1./2.*e1_q3*e3_q4*e4_q3*yyy42 - 1./2.*&
1200 ! & e1_q3*e3_q4*e4_q3*yyy41 - e1_q3*e3_q4*e4_q3*MZ4**2*yyy3 + e1_q3*&
1201 ! & e3_q4*e4_q3*MZ3**2*yyy3 + e1_q3*e3_q4*e4_q3*MG**2*yyy3 + 1./2.*&
1202 ! & et1(q1,e3,q,q3)*e1_q3*e4_q3*MG**(-2)*yyy7 - 1./2.*et1(q1,e3,q,q4&
1203 ! & )*e1_q3*e4_q3*MG**(-2)*yyy7 + 1./2.*et1(q1,e4,q,q3)*e1_q3*e3_q4*&
1204 ! & MG**(-2)*yyy7 - 1./2.*et1(q1,e4,q,q4)*e1_q3*e3_q4*MG**(-2)*yyy7&
1205 ! & - 1./2.*et1(q2,e3,q,q3)*e1_q3*e4_q3*MG**(-2)*yyy7
1206 ! res = res + 1./2.*et1(q2,e3,q,q4)*e1_q3*e4_q3*MG**(-2)*yyy7 - 1./&
1207 ! & 2.*et1(q2,e4,q,q3)*e1_q3*e3_q4*MG**(-2)*yyy7 + 1./2.*et1(q2,e4,q&
1208 ! & ,q4)*e1_q3*e3_q4*MG**(-2)*yyy7 + et1(e1,e3,q,q3)*q1_q3*e4_q3*&
1209 ! & MG**(-2)*yyy7 - 1./4.*et1(e1,e3,q,q3)*e4_q3*MG**(-2)*MZ4**2*yyy7&
1210 ! & + 1./4.*et1(e1,e3,q,q3)*e4_q3*MG**(-2)*MZ3**2*yyy7 + 1./4.*et1(&
1211 ! & e1,e3,q,q3)*e4_q3*yyy7 - et1(e1,e3,q,q4)*q1_q3*e4_q3*MG**(-2)*&
1212 ! & yyy7 + 1./4.*et1(e1,e3,q,q4)*e4_q3*MG**(-2)*MZ4**2*yyy7 - 1./4.*&
1213 ! & et1(e1,e3,q,q4)*e4_q3*MG**(-2)*MZ3**2*yyy7 - 1./4.*et1(e1,e3,q,&
1214 ! & q4)*e4_q3*yyy7 + et1(e1,e4,q,q3)*q1_q3*e3_q4*MG**(-2)*yyy7 - 1./&
1215 ! & 4.*et1(e1,e4,q,q3)*e3_q4*MG**(-2)*MZ4**2*yyy7 + 1./4.*et1(e1,e4,&
1216 ! & q,q3)*e3_q4*MG**(-2)*MZ3**2*yyy7 + 1./4.*et1(e1,e4,q,q3)*e3_q4*&
1217 ! & yyy7 - et1(e1,e4,q,q4)*q1_q3*e3_q4*MG**(-2)*yyy7 + 1./4.*et1(e1,&
1218 ! & e4,q,q4)*e3_q4*MG**(-2)*MZ4**2*yyy7 - 1./4.*et1(e1,e4,q,q4)*&
1219 ! & e3_q4*MG**(-2)*MZ3**2*yyy7 - 1./4.*et1(e1,e4,q,q4)*e3_q4*yyy7 + &
1220 ! & 1./2.*et1(e3,e4,q1,q)*e1_q3*yyy6 - 1./2.*et1(e3,e4,q2,q)*e1_q3*&
1221 ! & yyy6
1222 ! res = res - 1./4.*et1(e3,e4,e1,q)*MZ4**2*yyy6 + 1./4.*et1(e3,e4,&
1223 ! & e1,q)*MZ3**2*yyy6 + 1./4.*et1(e3,e4,e1,q)*MG**2*yyy6 + et1(e3,e4&
1224 ! & ,e1,q)*q1_q3*yyy6 + 4.*et1(e3,e4,q3,q4)*q1_q3*e1_q3*MG**(-2)*&
1225 ! & yyy5 - et1(e3,e4,q3,q4)*e1_q3*MG**(-2)*MZ4**2*yyy5 + et1(e3,e4,&
1226 ! & q3,q4)*e1_q3*MG**(-2)*MZ3**2*yyy5 + et1(e3,e4,q3,q4)*e1_q3*yyy5
1227 
1228 ! print *, "res new QQB",res
1229 
1230 
1231 
1232 ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! SAME AS ABOVE BUT SHORTER ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! ! !
1233 
1234  abr1 = (mg**2 + mz3**2 - mz4**2 + 4*q1_q3)
1235 
1236  res = (e1_e3*(4*yyy1*q1_e4 + e4_q3*(2*yyy1 - (mg**2 + mz3**2 - mz4**2)*yyy41 - 4*yyy41*q1_q3)) + &
1237  e1_e4*(4*yyy1*q1_e3 + e3_q4*(2*yyy1 + (mg**2 + mz3**2 - mz4**2)*yyy42 + 4*yyy42*q1_q3)) + &
1238  2*e1_q3*(-2*yyy41*e4_q3*q1_e3 + 2*yyy2*e3_e4*abr1 + &
1239  e3_q4*(2*yyy42*q1_e4 + e4_q3*(2*(mg**2 + mz3**2 - mz4**2)*yyy3 - yyy41 + yyy42 + 8*yyy3*q1_q3))))/4d0 + &
1240  (yyy6*abr1*et1(e3,e4,e1,q))/4d0 + (yyy6*e1_q3*et1(e3,e4,q1,q))/2d0 &
1241  - (yyy6*e1_q3*et1(e3,e4,q2,q))/2d0 + &
1242  (yyy5*e1_q3*abr1*et1(e3,e4,q3,q4))/mg**2 + &
1243  yyy7*((e4_q3*abr1*et1(e1,e3,q,q3))/(4d0*mg**2) - &
1244  (e4_q3*abr1*et1(e1,e3,q,q4))/(4d0*mg**2) + &
1245  (e3_q4*abr1*et1(e1,e4,q,q3))/(4d0*mg**2) - &
1246  (e3_q4*abr1*et1(e1,e4,q,q4))/(4d0*mg**2) &
1247  + (e1_q3*e4_q3*et1(q1,e3,q,q3))/(2d0*mg**2) - &
1248  (e1_q3*e4_q3*et1(q1,e3,q,q4))/(2d0*mg**2) + (e1_q3*e3_q4*et1(q1,e4,q,q3))/(2d0*mg**2) - &
1249  (e1_q3*e3_q4*et1(q1,e4,q,q4))/(2d0*mg**2) - (e1_q3*e4_q3*et1(q2,e3,q,q3))/(2d0*mg**2) + &
1250  (e1_q3*e4_q3*et1(q2,e3,q,q4))/(2d0*mg**2) - (e1_q3*e3_q4*et1(q2,e4,q,q3))/(2d0*mg**2) + &
1251  (e1_q3*e3_q4*et1(q2,e4,q,q4))/(2d0*mg**2))
1252 
1253 ! print *, "res newer QQB ",res
1254 ! pause
1255 
1256  end subroutine qqgzzampl
1257 
1258 
1259 !----- a subroutine for G -> ZZ/WW/AA
1260 !----- all outgoing convention and the following momentum assignment
1261 !----- 0 -> G(p1) + e-(p3) + e+(p4) +mu-(p5) +mu+(p6)
1262  subroutine evalamp_g_vv(p,MY_IDUP,res)
1263  implicit none
1264  real(dp), intent(out) :: res
1265  real(dp), intent(in) :: p(4,6)
1266  integer, intent(in) :: MY_IDUP(6:9)
1267  complex(dp) :: A_VV(1:18),A0_VV(1:2)
1268  integer :: i1,i3,i4,VVMode,VVmode_swap
1269  real(dp) :: prefactor
1270  real(dp) :: intcolfac
1271  integer :: ordering(1:4),ordering_swap(1:4)
1272  logical :: doInterference
1273 
1274  if(isaquark(my_idup(6)) .and. isaquark(my_idup(8))) then
1275  intcolfac=1.0_dp/3.0_dp
1276  else
1277  intcolfac=1.0_dp
1278  endif
1279 
1280  call getdecay_vvmode_ordering(my_idup(6:9),vvmode,ordering,vvmode_swap,ordering_swap)
1281 
1282 !---- full prefactor
1283  if( vvmode.eq.zzmode ) then! Z decay
1284  prefactor = overallcouplvffsq**2
1285  elseif( vvmode.eq.wwmode ) then ! W decay
1286  prefactor = overallcouplvffsq**2
1287  elseif( vvmode.eq.zgmode ) then ! Z+photon "decay"
1288  prefactor = overallcouplvffsq ! Only single powers
1289  elseif( vvmode.eq.ggmode ) then ! photon "decay"
1290  prefactor = 1d0
1291  else
1292  prefactor = 0d0
1293  endif
1294 
1295  res = zero
1296  a_vv(:) = 0d0
1297  dointerference = includeinterference .and. ( &
1298  ((vvmode.eq.zzmode) .and. (vvmode_swap.eq.zzmode)) &
1299  )
1300  if ( includevprime .and. .not.(vvmode.eq.zzmode .or. vvmode.eq.zgmode .or. vvmode.eq.wwmode) ) then
1301  call error("Contact terms only for WW, ZZ or Zg!")
1302  endif
1303  do i1 =-2,2; do i3=1,2; do i4=1,2! sum over helicities
1304  call calchelamp2(ordering,vvmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(1))
1305  if( vvmode.eq.zzmode ) then
1306  if( includegammastar ) then
1307  call calchelamp2(ordering,zgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(3))
1308  call calchelamp2(ordering,gszmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(5))
1309  call calchelamp2(ordering,gsgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(7))
1310  endif
1311  if( includevprime ) then
1312  call calchelamp2(ordering,zzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(9))
1313  call calchelamp2(ordering,zpzmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(11))
1314  call calchelamp2(ordering,zpzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(13))
1315  endif
1316  if( includegammastar .and. includevprime ) then
1317  call calchelamp2(ordering,gszpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(15))
1318  call calchelamp2(ordering,zpgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(17))
1319  endif
1320  elseif( vvmode.eq.zgmode ) then
1321  if(includegammastar) then
1322  call calchelamp2(ordering,gsgmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(3))
1323  endif
1324  if( includevprime ) then
1325  call calchelamp2(ordering,zpgmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(5))
1326  endif
1327  elseif( vvmode.eq.wwmode .and. includevprime ) then
1328  call calchelamp2(ordering,wwpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(9))
1329  call calchelamp2(ordering,wpwmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(11))
1330  call calchelamp2(ordering,wpwpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(13))
1331  endif
1332 
1333  if( dointerference ) then
1334  call calchelamp2(ordering_swap,vvmode_swap,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(2))
1335  if( includegammastar ) then
1336  call calchelamp2(ordering_swap,zgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(4))
1337  call calchelamp2(ordering_swap,gszmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(6))
1338  call calchelamp2(ordering_swap,gsgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(8))
1339  endif
1340  if( includevprime ) then
1341  call calchelamp2(ordering_swap,zzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(10))
1342  call calchelamp2(ordering_swap,zpzmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(12))
1343  call calchelamp2(ordering_swap,zpzpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(14))
1344  endif
1345  if( includegammastar .and. includevprime ) then
1346  call calchelamp2(ordering_swap,gszpmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(16))
1347  call calchelamp2(ordering_swap,zpgsmode,p(1:4,1:6),my_idup,i1,i3,i4,a_vv(18))
1348  endif
1349  endif
1350 
1351  a0_vv(1) = a_vv(1)+a_vv(3)+a_vv(5)+a_vv(7)+a_vv(9)+a_vv(11)+a_vv(13)+a_vv(15)+a_vv(17) ! 3456 pieces
1352  a0_vv(2) = a_vv(2)+a_vv(4)+a_vv(6)+a_vv(8)+a_vv(10)+a_vv(12)+a_vv(14)+a_vv(16)+a_vv(18) ! 5436 pieces
1353  res = res + dreal(a0_vv(1)*dconjg(a0_vv(1)))
1354  res = res + dreal(a0_vv(2)*dconjg(a0_vv(2)))
1355  if( dointerference .and. (i3.eq.i4) ) then! interfere the 3456 with 5436 pieces
1356  res = res - 2d0*intcolfac*dreal( a0_vv(1)*dconjg(a0_vv(2)) ) ! minus from Fermi statistics
1357  endif
1358  enddo; enddo; enddo
1359 
1360  res = res*prefactor
1361  if( (vvmode.eq.zzmode) .and. dointerference ) res = res * symmfac
1362 
1363  end subroutine
1364 
1365  subroutine calchelamp2(ordering,VVMode,p,MY_IDUP,i1,i3,i4,A)
1366  implicit none
1367  integer :: ordering(1:4),i1,i3,i4,l1,l2,l3,l4,mu,nu,MY_IDUP(6:9),VVMode
1368  real(dp) :: p(1:4,1:6)
1369  real(dp) :: pin(4,4)
1370  complex(dp) :: A(1:1), sp(0:4,1:4)
1371 
1372  l1=ordering(1)
1373  l2=ordering(2)
1374  l3=ordering(3)
1375  l4=ordering(4)
1376 
1377  pin(1,:) = p(:,1)
1378  pin(2,:) = 0d0 ! dummy
1379 
1380  sp(0,1:4) = pol_mass2(dcmplx(p(1:4,1)), 0,'in')
1381  sp(1,1:4) = pol_mass2(dcmplx(p(1:4,1)),-1,'in')
1382  sp(2,1:4) = pol_mass2(dcmplx(p(1:4,1)),+1,'in')
1383  call getdecay_couplings_spinors_props( &
1384  vvmode, &
1385  (/my_idup(l1+3),my_idup(l2+3),my_idup(l3+3),my_idup(l4+3)/), &
1386  (/p(:,l1),p(:,l2),p(:,l3),p(:,l4)/), &
1387  -3+2*i3,-3+2*i4, &
1388  sp(3:4,:),pin(3:4,:) &
1389  )
1390  call gzzampl(vvmode,pin,sp,i1,a(1))
1391  end subroutine
1392 
1393  subroutine gzzampl(VVMode,p,sp,i1,res)
1394  implicit none
1395  integer, intent(in) :: VVMode
1396  real(dp), intent(in) :: p(4,4)
1397  complex(dp), intent(in) :: sp(0:4,4)
1398  integer,intent(in) :: i1
1399  complex(dp), intent(out) :: res
1400  complex(dp) :: e1_e2, e1_e3, e1_e4
1401  complex(dp) :: e2_e3, e2_e4
1402  complex(dp) :: e3_e4
1403  complex(dp) :: q_q,q3_q3,q4_q4
1404  complex(dp) :: q1_q2,q1_q3,q1_q4
1405  complex(dp) :: q2_q3,q2_q4
1406  complex(dp) :: q3_q4
1407  complex(dp) :: q1_e3,q1_e4,q2_e3,q2_e4,e0_e3,e0_e4
1408  complex(dp) :: e1_q3,e1_q4,e2_q3,e2_q4,e0_q3,e0_q4,q_e3,q_e4
1409  complex(dp) :: e3_q4,e4_q3
1410  complex(dp) :: q1(4),q2(4),q3(4),q4(4),q(4)
1411  complex(dp) :: e1(4),e2(4),e3(4),e4(4),e0(4)
1412  complex(dp) :: yyy1,yyy2,yyy3,yyy41,yyy42,yyy5,yyy6,yyy7
1413  complex(dp) :: b_dyn(1:10)
1414  real(dp) :: q34,MG,MZ3,MZ4
1415  real(dp) :: rr
1416  real(dp),parameter :: sqrt6=dsqrt(6d0)
1417 
1418  b_dyn(:)=czero
1419 
1420  q1 = dcmplx(p(1,:),0d0)
1421  q2 = dcmplx(p(2,:),0d0)
1422  q3 = dcmplx(p(3,:),0d0)
1423  q4 = dcmplx(p(4,:),0d0)
1424 
1425 ! requirement: sp(0,:)= 0 pol
1426 ! requirement: sp(1,:)= - pol
1427 ! requirement: sp(2,:)= + pol
1428  if( i1.eq.+2 ) then
1429  e0 = 1d9! dummy
1430  e1 = sp(2,:)! +
1431  e2 = sp(2,:)! +
1432  elseif( i1.eq.+1 ) then
1433  e0 = 1d9! dummy
1434  e1 = sp(2,:)! +
1435  e2 = sp(0,:)! 0
1436  elseif( i1.eq.0 ) then
1437  e0 = sp(0,:)! 0
1438  e1 = sp(1,:)! -
1439  e2 = sp(2,:)! +
1440  elseif( i1.eq.-1 ) then
1441  e0 = 1d9! dummy
1442  e1 = sp(1,:)! -
1443  e2 = sp(0,:)! 0
1444  elseif( i1.eq.-2 ) then
1445  e0 = 1d9! dummy
1446  e1 = sp(1,:)! -
1447  e2 = sp(1,:)! -
1448  endif
1449 
1450  e3 = sp(3,:)
1451  e4 = sp(4,:)
1452 
1453  q = -q1-q2
1454 
1455  q_q =sc(q,q)
1456  q3_q3 = sc(q3,q3)
1457  q4_q4 = sc(q4,q4)
1458 
1459 
1460  q1_q2 = sc(q1,q2)
1461  q1_q3 = sc(q1,q3)
1462  q1_q4 = sc(q1,q4)
1463  q2_q3 = sc(q2,q3)
1464  q2_q4 = sc(q2,q4)
1465  q3_q4 = sc(q3,q4)
1466 
1467  e1_e2 = sc(e1,e2)
1468  e1_e3 = sc(e1,e3)
1469  e1_e4 = sc(e1,e4)
1470 
1471  e2_e3 = sc(e2,e3)
1472  e2_e4 = sc(e2,e4)
1473 
1474  e3_e4 = sc(e3,e4)
1475 
1476 ! new
1477  e0_e3 = sc(e0,e3)
1478  e0_e4 = sc(e0,e4)
1479  e0_q3 = sc(q3,e0)
1480  e0_q4 = sc(q4,e0)
1481  q_e3 = sc(q4,e3)
1482  q_e4 = sc(q3,e4)
1483 
1484  q1_e3 = sc(q1,e3)
1485  q1_e4 = sc(q1,e4)
1486  q2_e3 = sc(q2,e3)
1487  q2_e4 = sc(q2,e4)
1488  e1_q3 = sc(e1,q3)
1489  e1_q4 = sc(e1,q4)
1490  e2_q3 = sc(e2,q3)
1491  e2_q4 = sc(e2,q4)
1492  e3_q4 = sc(e3,q4)
1493  e4_q3 = sc(e4,q3)
1494 
1495  mz3=dsqrt(cdabs(q3_q3))
1496  mz4=dsqrt(cdabs(q4_q4))
1497  if( use_dynamic_mg ) then
1498  mg = dsqrt(cdabs(q_q))
1499  else
1500  mg = m_reso
1501  endif
1502 
1503  q34 = (mg**2-mz3**2-mz4**2)/2d0
1504 
1505  if (generate_bis) then
1506  rr = q34/lambda**2! kappa for FS
1507 
1508  if( (vvmode.eq.zzmode) .or. (vvmode.eq.wwmode) ) then! decay ZZ's or WW's
1509  b_dyn(1)=b1
1510  b_dyn(2)=b2
1511  b_dyn(3)=b3
1512  b_dyn(4)=b4
1513  b_dyn(5)=b5
1514  b_dyn(6)=b6
1515  b_dyn(7)=b7
1516  b_dyn(8)=b8
1517  b_dyn(9)=b9
1518  b_dyn(10)=b10
1519  elseif( (vvmode.eq.zgmode) .OR. (vvmode.eq.gszmode) .OR. (vvmode.eq.zgsmode) ) then
1520  b_dyn(1)=bzgs1
1521  b_dyn(2)=bzgs2
1522  b_dyn(3)=bzgs3
1523  b_dyn(4)=bzgs4
1524  b_dyn(8)=bzgs8
1525  elseif( (vvmode.eq.ggmode) .or. (vvmode.eq.gsgsmode) .or. (vvmode.eq.gsgmode) ) then
1526  b_dyn(1)=bgsgs1
1527  b_dyn(2)=bgsgs2
1528  b_dyn(3)=bgsgs3
1529  b_dyn(4)=bgsgs4
1530  b_dyn(8)=bgsgs8
1531  elseif( (vvmode.eq.zzpmode) .or. (vvmode.eq.wwpmode) .or. (vvmode.eq.zpzmode) .or. (vvmode.eq.wpwmode) ) then
1532  b_dyn(1)=bzzp1
1533  b_dyn(2)=bzzp2
1534  b_dyn(3)=bzzp3
1535  b_dyn(4)=bzzp4
1536  b_dyn(5)=bzzp5
1537  b_dyn(6)=bzzp6
1538  b_dyn(7)=bzzp7
1539  b_dyn(8)=bzzp8
1540  b_dyn(9)=bzzp9
1541  b_dyn(10)=bzzp10
1542  elseif( (vvmode.eq.zpzpmode) .or. (vvmode.eq.wpwpmode) ) then
1543  b_dyn(1)=bzpzp1
1544  b_dyn(2)=bzpzp2
1545  b_dyn(3)=bzpzp3
1546  b_dyn(4)=bzpzp4
1547  b_dyn(5)=bzpzp5
1548  b_dyn(6)=bzpzp6
1549  b_dyn(7)=bzpzp7
1550  b_dyn(8)=bzpzp8
1551  b_dyn(9)=bzpzp9
1552  b_dyn(10)=bzpzp10
1553  elseif( (vvmode.eq.zpgmode) .OR. (vvmode.eq.gszpmode) .OR. (vvmode.eq.zpgsmode) ) then
1554  b_dyn(1)=bzpgs1
1555  b_dyn(2)=bzpgs2
1556  b_dyn(3)=bzpgs3
1557  b_dyn(4)=bzpgs4
1558  b_dyn(8)=bzpgs8
1559  else
1560  print *,"VVMode",vvmode,"not implemented"
1561  endif
1562 
1563  yyy1 = q34*( b_dyn(1) + b_dyn(2)*rr*(one+mz3**2/q34)*(one+mz4**2/q34) ) + b_dyn(5)*m_v**2
1564  yyy2 = -b_dyn(1)/two + b_dyn(3)*rr*(1d0-(mz3**2+mz4**2)/(2d0*q34)) + two*b_dyn(4)*rr + b_dyn(7)*rr*m_v**2/q34
1565  yyy3 = (-b_dyn(2)/two - b_dyn(3)- two*b_dyn(4))*rr/q34
1566  yyy41 = -b_dyn(1) - b_dyn(2)*(q34+mz3**2)/lambda**2 - b_dyn(3)*mz4**2/lambda**2 - 2d0*b_dyn(6)*m_v**2/lambda**2
1567  yyy42 = -b_dyn(1) - b_dyn(2)*(q34+mz4**2)/lambda**2 - b_dyn(3)*mz3**2/lambda**2 - 2d0*b_dyn(6)*m_v**2/lambda**2
1568  yyy5 = two*b_dyn(8)*rr*mg**2/q34
1569  yyy6 = b_dyn(9) * m_v**2/lambda**2
1570  yyy7 = b_dyn(10) * mg**2 * m_v**2/lambda**4
1571 
1572  else
1573  yyy1 = q34*c1/2d0
1574  yyy2 = c2
1575  yyy3 = c3/mg**2
1576  yyy41 = c41
1577  yyy42 = c42
1578  yyy5 = c5
1579  yyy6 = czero
1580  yyy7 = czero
1581  if(vvmode.eq.zzmode .or. vvmode.eq.wwmode) then
1582  yyy6 = c6
1583  yyy7 = c7
1584  endif
1585  endif
1586 
1587 
1588  res = czero
1589 
1590 
1591  if( abs(i1).eq.2 ) then
1592 
1593  res =&
1594  & + yyy7 * ( et1(e1,e3,q,q3)*q_e4*e2_q3*mg**(-2) - et1(e1,e3,q,q3&
1595  & )*q_e4*e2_q4*mg**(-2) - et1(e1,e3,q,q4)*q_e4*e2_q3*mg**(-2)&
1596  & + et1(e1,e3,q,q4)*q_e4*e2_q4*mg**(-2) + et1(e1,e4,q,q3)*q_e3&
1597  & *e2_q3*mg**(-2) - et1(e1,e4,q,q3)*q_e3*e2_q4*mg**(-2) - et1(&
1598  & e1,e4,q,q4)*q_e3*e2_q3*mg**(-2) + et1(e1,e4,q,q4)*q_e3*e2_q4*&
1599  & mg**(-2) )
1600  res = res + yyy6 * ( et1(e3,e4,e1,q)*e2_q3 - et1(e3,e4,e1,q)*&
1601  & e2_q4 )
1602  res = res + yyy5 * ( et1(e3,e4,q3,q4)*e1_q3*e2_q3*mg**(-2) - et1(&
1603  & e3,e4,q3,q4)*e1_q3*e2_q4*mg**(-2) - et1(e3,e4,q3,q4)*e1_q4*&
1604  & e2_q3*mg**(-2) + et1(e3,e4,q3,q4)*e1_q4*e2_q4*mg**(-2) )
1605  res = res + yyy42 * ( e1_e4*e2_q3*e3_q4 + e1_q3*e2_e4*e3_q4 )
1606  res = res + yyy41 * ( e1_e3*e2_q4*e4_q3 + e1_q4*e2_e3*e4_q3 )
1607  res = res + yyy3 * ( e1_q3*e2_q3*e3_q4*e4_q3 - e1_q3*e2_q4*e3_q4*&
1608  & e4_q3 - e1_q4*e2_q3*e3_q4*e4_q3 + e1_q4*e2_q4*e3_q4*e4_q3 )
1609  res = res + yyy2 * ( e1_q3*e2_q3*e3_e4 - e1_q3*e2_q4*e3_e4 - &
1610  & e1_q4*e2_q3*e3_e4 + e1_q4*e2_q4*e3_e4 )
1611  res = res + yyy1 * ( e1_e3*e2_e4 + e1_e4*e2_e3 )
1612 
1613 
1614  elseif( abs(i1).eq.1 ) then
1615 
1616  res =&
1617  & + yyy7*sqrt2**(-1) * ( et1(e1,e3,q,q3)*q_e4*e2_q3*mg**(-2) - &
1618  & et1(e1,e3,q,q3)*q_e4*e2_q4*mg**(-2) - et1(e1,e3,q,q4)*q_e4*&
1619  & e2_q3*mg**(-2) + et1(e1,e3,q,q4)*q_e4*e2_q4*mg**(-2) + et1(e1&
1620  & ,e4,q,q3)*q_e3*e2_q3*mg**(-2) - et1(e1,e4,q,q3)*q_e3*e2_q4*&
1621  & mg**(-2) - et1(e1,e4,q,q4)*q_e3*e2_q3*mg**(-2) + et1(e1,e4,q,&
1622  & q4)*q_e3*e2_q4*mg**(-2) + et1(e2,e3,q,q3)*q_e4*e1_q3*mg**(-2)&
1623  & - et1(e2,e3,q,q3)*q_e4*e1_q4*mg**(-2) - et1(e2,e3,q,q4)*q_e4&
1624  & *e1_q3*mg**(-2) + et1(e2,e3,q,q4)*q_e4*e1_q4*mg**(-2) + et1(&
1625  & e2,e4,q,q3)*q_e3*e1_q3*mg**(-2) - et1(e2,e4,q,q3)*q_e3*e1_q4*&
1626  & mg**(-2) - et1(e2,e4,q,q4)*q_e3*e1_q3*mg**(-2) + et1(e2,e4,q,&
1627  & q4)*q_e3*e1_q4*mg**(-2) )
1628  res = res + yyy6*sqrt2**(-1) * ( et1(e3,e4,e1,q)*e2_q3 - et1(e3,&
1629  & e4,e1,q)*e2_q4 + et1(e3,e4,e2,q)*e1_q3 - et1(e3,e4,e2,q)*&
1630  & e1_q4 )
1631  res = res + yyy5*sqrt2**(-1) * ( 2.0d0*et1(e3,e4,q3,q4)*e1_q3*e2_q3*&
1632  & mg**(-2) - 2.0d0*et1(e3,e4,q3,q4)*e1_q3*e2_q4*mg**(-2) - 2.0d0*et1(&
1633  & e3,e4,q3,q4)*e1_q4*e2_q3*mg**(-2) + 2.0d0*et1(e3,e4,q3,q4)*e1_q4&
1634  & *e2_q4*mg**(-2) )
1635  res = res + yyy42*sqrt2**(-1) * ( 2.0d0*e1_e4*e2_q3*e3_q4 + 2.0d0*e1_q3&
1636  & *e2_e4*e3_q4 )
1637  res = res + yyy41*sqrt2**(-1) * ( 2.0d0*e1_e3*e2_q4*e4_q3 + 2.0d0*e1_q4&
1638  & *e2_e3*e4_q3 )
1639  res = res + yyy3*sqrt2**(-1) * ( 2.0d0*e1_q3*e2_q3*e3_q4*e4_q3 - 2.0d0*&
1640  & e1_q3*e2_q4*e3_q4*e4_q3 - 2.0d0*e1_q4*e2_q3*e3_q4*e4_q3 + 2.0d0*&
1641  & e1_q4*e2_q4*e3_q4*e4_q3 )
1642  res = res + yyy2*sqrt2**(-1) * ( 2.0d0*e1_q3*e2_q3*e3_e4 - 2.0d0*e1_q3*&
1643  & e2_q4*e3_e4 - 2.0d0*e1_q4*e2_q3*e3_e4 + 2.0d0*e1_q4*e2_q4*e3_e4 )
1644  res = res + yyy1*sqrt2**(-1) * ( 2.0d0*e1_e3*e2_e4 + 2.0d0*e1_e4*e2_e3&
1645  & )
1646 
1647 
1648  elseif( abs(i1).eq.0 ) then
1649  res =&
1650  & + yyy7*sqrt6**(-1) * ( - 2.0d0*et1(e0,e3,q,q3)*q_e4*e0_q3*&
1651  & mg**(-2) + 2.0d0*et1(e0,e3,q,q3)*q_e4*e0_q4*mg**(-2) + 2.0d0*et1(e0&
1652  & ,e3,q,q4)*q_e4*e0_q3*mg**(-2) - 2.0d0*et1(e0,e3,q,q4)*q_e4*e0_q4&
1653  & *mg**(-2) - 2.0d0*et1(e0,e4,q,q3)*q_e3*e0_q3*mg**(-2) + 2.0d0*et1(&
1654  & e0,e4,q,q3)*q_e3*e0_q4*mg**(-2) + 2.0d0*et1(e0,e4,q,q4)*q_e3*&
1655  & e0_q3*mg**(-2) - 2.0d0*et1(e0,e4,q,q4)*q_e3*e0_q4*mg**(-2) + &
1656  & et1(e1,e3,q,q3)*q_e4*e2_q3*mg**(-2) - et1(e1,e3,q,q3)*q_e4*&
1657  & e2_q4*mg**(-2) - et1(e1,e3,q,q4)*q_e4*e2_q3*mg**(-2) + et1(e1&
1658  & ,e3,q,q4)*q_e4*e2_q4*mg**(-2) + et1(e1,e4,q,q3)*q_e3*e2_q3*&
1659  & mg**(-2) - et1(e1,e4,q,q3)*q_e3*e2_q4*mg**(-2) - et1(e1,e4,q,&
1660  & q4)*q_e3*e2_q3*mg**(-2) + et1(e1,e4,q,q4)*q_e3*e2_q4*mg**(-2)&
1661  & + et1(e2,e3,q,q3)*q_e4*e1_q3*mg**(-2) - et1(e2,e3,q,q3)*q_e4&
1662  & *e1_q4*mg**(-2) - et1(e2,e3,q,q4)*q_e4*e1_q3*mg**(-2) + et1(&
1663  & e2,e3,q,q4)*q_e4*e1_q4*mg**(-2) + et1(e2,e4,q,q3)*q_e3*e1_q3*&
1664  & mg**(-2) - et1(e2,e4,q,q3)*q_e3*e1_q4*mg**(-2) - et1(e2,e4,q,&
1665  & q4)*q_e3*e1_q3*mg**(-2) )
1666  res = res + yyy7*sqrt6**(-1) * ( et1(e2,e4,q,q4)*q_e3*e1_q4*&
1667  & mg**(-2) )
1668  res = res + yyy6*sqrt6**(-1) * ( - 2.0d0*et1(e3,e4,e0,q)*e0_q3 + 2.0d0&
1669  & *et1(e3,e4,e0,q)*e0_q4 + et1(e3,e4,e1,q)*e2_q3 - et1(e3,e4,e1&
1670  & ,q)*e2_q4 + et1(e3,e4,e2,q)*e1_q3 - et1(e3,e4,e2,q)*e1_q4 )
1671  res = res + yyy5*sqrt6**(-1) * ( 4.0d0*et1(e3,e4,q3,q4)*e0_q3*e0_q4*&
1672  & mg**(-2) - 2.0d0*et1(e3,e4,q3,q4)*e0_q3**2*mg**(-2) - 2.0d0*et1(e3,&
1673  & e4,q3,q4)*e0_q4**2*mg**(-2) + 2.0d0*et1(e3,e4,q3,q4)*e1_q3*e2_q3&
1674  & *mg**(-2) - 2.0d0*et1(e3,e4,q3,q4)*e1_q3*e2_q4*mg**(-2) - 2.0d0*&
1675  & et1(e3,e4,q3,q4)*e1_q4*e2_q3*mg**(-2) + 2.0d0*et1(e3,e4,q3,q4)*&
1676  & e1_q4*e2_q4*mg**(-2) )
1677  res = res + yyy42*sqrt6**(-1) * ( - 4.0d0*e0_e4*e0_q3*e3_q4 + 2.0d0*&
1678  & e1_e4*e2_q3*e3_q4 + 2.0d0*e1_q3*e2_e4*e3_q4 )
1679  res = res + yyy41*sqrt6**(-1) * ( - 4.0d0*e0_e3*e0_q4*e4_q3 + 2.0d0*&
1680  & e1_e3*e2_q4*e4_q3 + 2.0d0*e1_q4*e2_e3*e4_q3 )
1681  res = res + yyy3*sqrt6**(-1) * ( 4.0d0*e0_q3*e0_q4*e3_q4*e4_q3 - 2.0d0*&
1682  & e0_q3**2*e3_q4*e4_q3 - 2.0d0*e0_q4**2*e3_q4*e4_q3 + 2.0d0*e1_q3*&
1683  & e2_q3*e3_q4*e4_q3 - 2.0d0*e1_q3*e2_q4*e3_q4*e4_q3 - 2.0d0*e1_q4*&
1684  & e2_q3*e3_q4*e4_q3 + 2.0d0*e1_q4*e2_q4*e3_q4*e4_q3 )
1685  res = res + yyy2*sqrt6**(-1) * ( 4.0d0*e0_q3*e0_q4*e3_e4 - 2.0d0*&
1686  & e0_q3**2*e3_e4 - 2.0d0*e0_q4**2*e3_e4 + 2.0d0*e1_q3*e2_q3*e3_e4 - 2.0d0&
1687  & *e1_q3*e2_q4*e3_e4 - 2.0d0*e1_q4*e2_q3*e3_e4 + 2.0d0*e1_q4*e2_q4*&
1688  & e3_e4 )
1689  res = res + yyy1*sqrt6**(-1) * ( - 4.0d0*e0_e3*e0_e4 + 2.0d0*e1_e3*&
1690  & e2_e4 + 2.0d0*e1_e4*e2_e3 )
1691 
1692 
1693  endif
1694 
1695  end subroutine gzzampl
1696 
1697 
1698 subroutine getdecay_couplings_spinors_props(VVMode,idordered,pordered,h3,h4, sp,pV)
1699  implicit none
1700  integer, intent(in) :: VVMode,idordered(6:9),h3,h4
1701  real(dp), intent(in) :: pordered(1:4,6:9)
1702  complex(dp), intent(out) :: sp(3:4,1:4)
1703  real(dp), intent(out) :: pV(3:4,1:4)
1704  real(dp) :: s, aL1,aR1,aL2,aR2
1705  complex(dp) :: propV(1:2)
1706 
1707  ! h3/h4 helicities: -1 == left, 1 == right
1708  if( vvmode.eq.zzmode ) then
1709  ! ZZ DECAYS
1710  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1711  al1=al_lep * dsqrt(scale_alpha_z_ll)
1712  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
1713  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1714  al1=al_lep * dsqrt(scale_alpha_z_tt)
1715  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
1716  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1717  al1=al_neu * dsqrt(scale_alpha_z_nn)
1718  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
1719  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1720  al1=al_qup * dsqrt(scale_alpha_z_uu)
1721  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
1722  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1723  al1=al_qdn * dsqrt(scale_alpha_z_dd)
1724  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
1725  else
1726  al1=0d0
1727  ar1=0d0
1728  endif
1729  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
1730  al2=al_lep * dsqrt(scale_alpha_z_ll)
1731  ar2=ar_lep * dsqrt(scale_alpha_z_ll)
1732  elseif( abs(idordered(8)).eq.abs(tam_) ) then
1733  al2=al_lep * dsqrt(scale_alpha_z_tt)
1734  ar2=ar_lep * dsqrt(scale_alpha_z_tt)
1735  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
1736  al2=al_neu * dsqrt(scale_alpha_z_nn)
1737  ar2=ar_neu * dsqrt(scale_alpha_z_nn)
1738  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
1739  al2=al_qup * dsqrt(scale_alpha_z_uu)
1740  ar2=ar_qup * dsqrt(scale_alpha_z_uu)
1741  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
1742  al2=al_qdn * dsqrt(scale_alpha_z_dd)
1743  ar2=ar_qdn * dsqrt(scale_alpha_z_dd)
1744  else
1745  al2=0d0
1746  ar2=0d0
1747  endif
1748  pv(3,:) = pordered(:,6)+pordered(:,7)
1749  pv(4,:) = pordered(:,8)+pordered(:,9)
1750  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1751  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1752  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1753  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1754  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1755  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
1756  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1757  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
1758 
1759 
1760 
1761  elseif( vvmode.eq.wwmode ) then
1762  ! WW DECAYS
1763  if( isaquark(idordered(6)) ) then
1764  al1 = bl * dsqrt(scale_alpha_w_ud)
1765  ar1 = br * dsqrt(scale_alpha_w_ud)! = 0
1766  elseif( &
1767  (abs(idordered(6)).eq.abs(elp_) .and. abs(idordered(7)).eq.abs(nue_)) .or. (abs(idordered(7)).eq.abs(elp_) .and. abs(idordered(6)).eq.abs(nue_)) .or. &
1768  (abs(idordered(6)).eq.abs(mup_) .and. abs(idordered(7)).eq.abs(num_)) .or. (abs(idordered(7)).eq.abs(mup_) .and. abs(idordered(6)).eq.abs(num_)) &
1769  ) then
1770  al1 = bl * dsqrt(scale_alpha_w_ln)
1771  ar1 = br * dsqrt(scale_alpha_w_ln)! = 0
1772  elseif( &
1773  (abs(idordered(6)).eq.abs(tap_) .and. abs(idordered(7)).eq.abs(nut_)) .or. (abs(idordered(7)).eq.abs(tap_) .and. abs(idordered(6)).eq.abs(nut_)) &
1774  ) then
1775  al1 = bl * dsqrt(scale_alpha_w_tn)
1776  ar1 = br * dsqrt(scale_alpha_w_tn)! = 0
1777  else
1778  al1=0d0
1779  ar1=0d0
1780  endif
1781  if( isaquark(idordered(8)) ) then
1782  al2 = bl * dsqrt(scale_alpha_w_ud)
1783  ar2 = br * dsqrt(scale_alpha_w_ud)! = 0
1784  elseif( &
1785  (abs(idordered(8)).eq.abs(elm_) .and. abs(idordered(9)).eq.abs(anue_)) .or. (abs(idordered(9)).eq.abs(elm_) .and. abs(idordered(8)).eq.abs(anue_)) .or. &
1786  (abs(idordered(8)).eq.abs(mum_) .and. abs(idordered(9)).eq.abs(anum_)) .or. (abs(idordered(9)).eq.abs(mum_) .and. abs(idordered(8)).eq.abs(anum_)) &
1787  ) then
1788  al2 = bl * dsqrt(scale_alpha_w_ln)
1789  ar2 = br * dsqrt(scale_alpha_w_ln)! = 0
1790  elseif( &
1791  (abs(idordered(8)).eq.abs(tam_) .and. abs(idordered(9)).eq.abs(anut_)) .or. (abs(idordered(9)).eq.abs(tam_) .and. abs(idordered(8)).eq.abs(anut_)) &
1792  ) then
1793  al2 = bl * dsqrt(scale_alpha_w_tn)
1794  ar2 = br * dsqrt(scale_alpha_w_tn)! = 0
1795  else
1796  al2=0d0
1797  ar2=0d0
1798  endif
1799  pv(3,:) = pordered(:,6)+pordered(:,7)
1800  pv(4,:) = pordered(:,8)+pordered(:,9)
1801  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1802  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1803  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1804  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1805  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1806  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
1807  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1808  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
1809 
1810 
1811  elseif( vvmode.eq.zgmode ) then
1812  ! Zgamma DECAYS
1813  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1814  al1=al_lep * dsqrt(scale_alpha_z_ll)
1815  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
1816  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1817  al1=al_lep * dsqrt(scale_alpha_z_tt)
1818  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
1819  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1820  al1=al_neu * dsqrt(scale_alpha_z_nn)
1821  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
1822  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1823  al1=al_qup * dsqrt(scale_alpha_z_uu)
1824  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
1825  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1826  al1=al_qdn * dsqrt(scale_alpha_z_dd)
1827  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
1828  else
1829  al1=0d0
1830  ar1=0d0
1831  endif
1832  al2=1d0
1833  ar2=1d0
1834  pv(3,:) = pordered(:,6)+pordered(:,7)
1835  pv(4,:) = pordered(:,8)
1836  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1837  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1838  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
1839  !sp(4,1:4)=pV(4,1:4); print *, "this checks gauge invariance"
1840  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1841  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
1842  propv(2)=1d0
1843 
1844 
1845  elseif( vvmode.eq.ggmode ) then
1846  ! gamma gamma DECAYS
1847  al1=1d0
1848  ar1=1d0
1849  al2=1d0
1850  ar2=1d0
1851  pv(3,:) = pordered(:,6)
1852  pv(4,:) = pordered(:,8)
1853  sp(3,:) = pol_mless2(dcmplx(pordered(:,6)),h3,'out') ! photon
1854  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
1855  !sp(3,1:4)=pV(3,1:4); print *, "this checks gauge invariance"
1856  !sp(4,1:4)=pV(4,1:4)
1857  propv(1)=1d0
1858  propv(2)=1d0
1859 
1860 
1861  elseif( vvmode.eq.gsgmode ) then
1862  ! gamma* gamma DECAYS
1863  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1864  al1=cl_lep * dsqrt(scale_alpha_z_ll)
1865  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
1866  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1867  al1=cl_lep * dsqrt(scale_alpha_z_tt)
1868  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
1869  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1870  al1=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
1871  ar1=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
1872  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1873  al1=cl_qup * dsqrt(scale_alpha_z_uu)
1874  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
1875  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1876  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
1877  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
1878  else
1879  al1=0d0
1880  ar1=0d0
1881  endif
1882  al2=1d0
1883  ar2=1d0
1884  pv(3,:) = pordered(:,6)+pordered(:,7)
1885  pv(4,:) = pordered(:,8)
1886  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1887  sp(3,:) = -sp(3,:) ! photon propagator
1888  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
1889  !sp(4,1:4)=pV(4,1:4); print *, "this checks gauge invariance"
1890  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1891  propv(1) = 1d0
1892  propv(2) = 1d0
1893  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
1894 
1895 
1896  elseif( vvmode.eq.gszmode ) then
1897  ! gamma* Z DECAYS
1898  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1899  al1=cl_lep * dsqrt(scale_alpha_z_ll)
1900  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
1901  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1902  al1=cl_lep * dsqrt(scale_alpha_z_tt)
1903  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
1904  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1905  al1=cl_neu * dsqrt(scale_alpha_z_nn)
1906  ar1=cr_neu * dsqrt(scale_alpha_z_nn)
1907  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1908  al1=cl_qup * dsqrt(scale_alpha_z_uu)
1909  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
1910  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1911  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
1912  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
1913  else
1914  al1=0d0
1915  ar1=0d0
1916  endif
1917  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
1918  al2=al_lep * dsqrt(scale_alpha_z_ll)
1919  ar2=ar_lep * dsqrt(scale_alpha_z_ll)
1920  elseif( abs(idordered(8)).eq.abs(tam_) ) then
1921  al2=al_lep * dsqrt(scale_alpha_z_tt)
1922  ar2=ar_lep * dsqrt(scale_alpha_z_tt)
1923  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
1924  al2=al_neu * dsqrt(scale_alpha_z_nn)
1925  ar2=ar_neu * dsqrt(scale_alpha_z_nn)
1926  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
1927  al2=al_qup * dsqrt(scale_alpha_z_uu)
1928  ar2=ar_qup * dsqrt(scale_alpha_z_uu)
1929  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
1930  al2=al_qdn * dsqrt(scale_alpha_z_dd)
1931  ar2=ar_qdn * dsqrt(scale_alpha_z_dd)
1932  else
1933  al2=0d0
1934  ar2=0d0
1935  endif
1936  pv(3,:) = pordered(:,6)+pordered(:,7)
1937  pv(4,:) = pordered(:,8)+pordered(:,9)
1938  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1939  sp(3,:) = -sp(3,:)
1940  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1941  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1942  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1943  propv(1) = 1d0! = s/dcmplx(s)
1944  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
1945  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1946  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
1947 
1948 
1949  elseif( vvmode.eq.zgsmode ) then
1950  ! Z gamma* DECAYS
1951  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1952  al1=al_lep * dsqrt(scale_alpha_z_ll)
1953  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
1954  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1955  al1=al_lep * dsqrt(scale_alpha_z_tt)
1956  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
1957  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1958  al1=al_neu * dsqrt(scale_alpha_z_nn)
1959  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
1960  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1961  al1=al_qup * dsqrt(scale_alpha_z_uu)
1962  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
1963  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1964  al1=al_qdn * dsqrt(scale_alpha_z_dd)
1965  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
1966  else
1967  al1=0d0
1968  ar1=0d0
1969  endif
1970  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
1971  al2=cl_lep * dsqrt(scale_alpha_z_ll)
1972  ar2=cr_lep * dsqrt(scale_alpha_z_ll)
1973  elseif( abs(idordered(8)).eq.abs(tam_) ) then
1974  al2=cl_lep * dsqrt(scale_alpha_z_tt)
1975  ar2=cr_lep * dsqrt(scale_alpha_z_tt)
1976  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
1977  al2=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
1978  ar2=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
1979  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
1980  al2=cl_qup * dsqrt(scale_alpha_z_uu)
1981  ar2=cr_qup * dsqrt(scale_alpha_z_uu)
1982  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
1983  al2=cl_qdn * dsqrt(scale_alpha_z_dd)
1984  ar2=cr_qdn * dsqrt(scale_alpha_z_dd)
1985  else
1986  al2=0d0
1987  ar2=0d0
1988  endif
1989  pv(3,:) = pordered(:,6)+pordered(:,7)
1990  pv(4,:) = pordered(:,8)+pordered(:,9)
1991  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1992  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1993  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1994  sp(4,:) = -sp(4,:)
1995  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1996  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
1997  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1998  propv(2) = 1d0 ! = s/dcmplx(s)
1999  if( s.lt.mphotoncutoff**2 ) propv(2)=czero
2000 
2001 
2002  elseif( vvmode.eq.gsgsmode ) then
2003  ! gamma* gamma* DECAYS
2004  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
2005  al1=cl_lep * dsqrt(scale_alpha_z_ll)
2006  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
2007  elseif( abs(idordered(6)).eq.abs(tam_) ) then
2008  al1=cl_lep * dsqrt(scale_alpha_z_tt)
2009  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
2010  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
2011  al1=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
2012  ar1=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
2013  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
2014  al1=cl_qup * dsqrt(scale_alpha_z_uu)
2015  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
2016  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
2017  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
2018  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
2019  else
2020  al1=0d0
2021  ar1=0d0
2022  endif
2023  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
2024  al2=cl_lep * dsqrt(scale_alpha_z_ll)
2025  ar2=cr_lep * dsqrt(scale_alpha_z_ll)
2026  elseif( abs(idordered(8)).eq.abs(tam_) ) then
2027  al2=cl_lep * dsqrt(scale_alpha_z_tt)
2028  ar2=cr_lep * dsqrt(scale_alpha_z_tt)
2029  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
2030  al2=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
2031  ar2=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
2032  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
2033  al2=cl_qup * dsqrt(scale_alpha_z_uu)
2034  ar2=cr_qup * dsqrt(scale_alpha_z_uu)
2035  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
2036  al2=cl_qdn * dsqrt(scale_alpha_z_dd)
2037  ar2=cr_qdn * dsqrt(scale_alpha_z_dd)
2038  else
2039  al2=0d0
2040  ar2=0d0
2041  endif
2042  pv(3,:) = pordered(:,6)+pordered(:,7)
2043  pv(4,:) = pordered(:,8)+pordered(:,9)
2044  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2045  sp(3,:) = -sp(3,:)
2046  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2047  sp(4,:) = -sp(4,:)
2048  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
2049  propv(1) = 1d0 ! = s/dcmplx(s)
2050  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
2051  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
2052  propv(2) = 1d0 ! = s/dcmplx(s)
2053  if( s.lt.mphotoncutoff**2 ) propv(2)=czero
2054 
2055 
2056  elseif( vvmode.eq.zpzmode ) then
2057  ! Z'Z DECAYS
2058  al1 = vpffcoupling(idordered(6),-1,.false.)
2059  ar1 = vpffcoupling(idordered(6),+1,.false.)
2060  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
2061  al2=al_lep * dsqrt(scale_alpha_z_ll)
2062  ar2=ar_lep * dsqrt(scale_alpha_z_ll)
2063  elseif( abs(idordered(8)).eq.abs(tam_) ) then
2064  al2=al_lep * dsqrt(scale_alpha_z_tt)
2065  ar2=ar_lep * dsqrt(scale_alpha_z_tt)
2066  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
2067  al2=al_neu * dsqrt(scale_alpha_z_nn)
2068  ar2=ar_neu * dsqrt(scale_alpha_z_nn)
2069  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
2070  al2=al_qup * dsqrt(scale_alpha_z_uu)
2071  ar2=ar_qup * dsqrt(scale_alpha_z_uu)
2072  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
2073  al2=al_qdn * dsqrt(scale_alpha_z_dd)
2074  ar2=ar_qdn * dsqrt(scale_alpha_z_dd)
2075  else
2076  al2=0d0
2077  ar2=0d0
2078  endif
2079 
2080  pv(3,:) = pordered(:,6)+pordered(:,7)
2081  pv(4,:) = pordered(:,8)+pordered(:,9)
2082  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2083  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2084  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2085  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
2086  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
2087  if( m_vprime .gt. 0d0 ) then
2088  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2089  elseif( m_vprime .eq. 0d0 ) then
2090  propv(1) = 1d0
2091  else
2092  propv(1) = s/m_z**2
2093  endif
2094  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
2095  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
2096 
2097 
2098  elseif( vvmode.eq.zzpmode ) then
2099  ! ZZ' DECAYS
2100  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
2101  al1=al_lep * dsqrt(scale_alpha_z_ll)
2102  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
2103  elseif( abs(idordered(6)).eq.abs(tam_) ) then
2104  al1=al_lep * dsqrt(scale_alpha_z_tt)
2105  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
2106  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
2107  al1=al_neu * dsqrt(scale_alpha_z_nn)
2108  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
2109  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
2110  al1=al_qup * dsqrt(scale_alpha_z_uu)
2111  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
2112  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
2113  al1=al_qdn * dsqrt(scale_alpha_z_dd)
2114  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
2115  else
2116  al1=0d0
2117  ar1=0d0
2118  endif
2119  al2 = vpffcoupling(idordered(8),-1,.false.)
2120  ar2 = vpffcoupling(idordered(8),+1,.false.)
2121 
2122  pv(3,:) = pordered(:,6)+pordered(:,7)
2123  pv(4,:) = pordered(:,8)+pordered(:,9)
2124  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2125  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2126  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2127  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
2128  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
2129  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
2130  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
2131  if( m_vprime .gt. 0d0 ) then
2132  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2133  elseif( m_vprime .eq. 0d0 ) then
2134  propv(2) = 1d0
2135  else
2136  propv(2) = s/m_z**2
2137  endif
2138 
2139 
2140  elseif( vvmode.eq.zpzpmode ) then
2141  ! Z'Z' DECAYS
2142  al1 = vpffcoupling(idordered(6),-1,.false.)
2143  ar1 = vpffcoupling(idordered(6),+1,.false.)
2144  al2 = vpffcoupling(idordered(8),-1,.false.)
2145  ar2 = vpffcoupling(idordered(8),+1,.false.)
2146 
2147  pv(3,:) = pordered(:,6)+pordered(:,7)
2148  pv(4,:) = pordered(:,8)+pordered(:,9)
2149  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2150  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2151  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2152  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
2153  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
2154  if( m_vprime .gt. 0d0 ) then
2155  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2156  elseif( m_vprime .eq. 0d0 ) then
2157  propv(1) = 1d0
2158  else
2159  propv(1) = s/m_z**2
2160  endif
2161  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
2162  if( m_vprime .gt. 0d0 ) then
2163  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2164  elseif( m_vprime .eq. 0d0 ) then
2165  propv(2) = 1d0
2166  else
2167  propv(2) = s/m_z**2
2168  endif
2169 
2170 
2171  elseif( vvmode.eq.zpgsmode ) then
2172  ! Z' gamma* DECAYS
2173  al1 = vpffcoupling(idordered(6),-1,.false.)
2174  ar1 = vpffcoupling(idordered(6),+1,.false.)
2175  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
2176  al2=cl_lep * dsqrt(scale_alpha_z_ll)
2177  ar2=cr_lep * dsqrt(scale_alpha_z_ll)
2178  elseif( abs(idordered(8)).eq.abs(tam_) ) then
2179  al2=cl_lep * dsqrt(scale_alpha_z_tt)
2180  ar2=cr_lep * dsqrt(scale_alpha_z_tt)
2181  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
2182  al2=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
2183  ar2=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
2184  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
2185  al2=cl_qup * dsqrt(scale_alpha_z_uu)
2186  ar2=cr_qup * dsqrt(scale_alpha_z_uu)
2187  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
2188  al2=cl_qdn * dsqrt(scale_alpha_z_dd)
2189  ar2=cr_qdn * dsqrt(scale_alpha_z_dd)
2190  else
2191  al2=0d0
2192  ar2=0d0
2193  endif
2194 
2195  pv(3,:) = pordered(:,6)+pordered(:,7)
2196  pv(4,:) = pordered(:,8)+pordered(:,9)
2197  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2198  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2199  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2200  sp(4,:) = -sp(4,:)
2201  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
2202  if( m_vprime .gt. 0d0 ) then
2203  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2204  elseif( m_vprime .eq. 0d0 ) then
2205  propv(1) = 1d0
2206  else
2207  propv(1) = s/m_z**2
2208  endif
2209  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
2210  propv(2) = 1d0 ! = s/dcmplx(s)
2211  if( s.lt.mphotoncutoff**2 ) propv(2)=czero
2212 
2213  elseif( vvmode.eq.gszpmode ) then
2214  ! gamma* Z' DECAYS
2215  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
2216  al1=cl_lep * dsqrt(scale_alpha_z_ll)
2217  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
2218  elseif( abs(idordered(6)).eq.abs(tam_) ) then
2219  al1=cl_lep * dsqrt(scale_alpha_z_tt)
2220  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
2221  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
2222  al1=cl_neu * dsqrt(scale_alpha_z_nn)
2223  ar1=cr_neu * dsqrt(scale_alpha_z_nn)
2224  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
2225  al1=cl_qup * dsqrt(scale_alpha_z_uu)
2226  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
2227  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
2228  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
2229  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
2230  else
2231  al1=0d0
2232  ar1=0d0
2233  endif
2234  al2 = vpffcoupling(idordered(8),-1,.false.)
2235  ar2 = vpffcoupling(idordered(8),+1,.false.)
2236 
2237  pv(3,:) = pordered(:,6)+pordered(:,7)
2238  pv(4,:) = pordered(:,8)+pordered(:,9)
2239  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2240  sp(3,:) = -sp(3,:)
2241  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2242  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
2243  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
2244  propv(1) = 1d0! = s/dcmplx(s)
2245  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
2246  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
2247  if( m_vprime .gt. 0d0 ) then
2248  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2249  elseif( m_vprime .eq. 0d0 ) then
2250  propv(2) = 1d0
2251  else
2252  propv(2) = s/m_z**2
2253  endif
2254 
2255 
2256  elseif( vvmode.eq.zpgmode ) then
2257  ! Z' gamma DECAYS
2258  al1 = vpffcoupling(idordered(6),-1,.false.)
2259  ar1 = vpffcoupling(idordered(6),+1,.false.)
2260  al2=1d0
2261  ar2=1d0
2262  pv(3,:) = pordered(:,6)+pordered(:,7)
2263  pv(4,:) = pordered(:,8)
2264  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2265  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2266  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
2267  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
2268  if( m_vprime .gt. 0d0 ) then
2269  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2270  elseif( m_vprime .eq. 0d0 ) then
2271  propv(1) = 1d0
2272  else
2273  propv(1) = s/m_z**2
2274  endif
2275  propv(2)=1d0
2276 
2277 
2278  elseif( vvmode.eq.wpwmode ) then
2279  ! W'W DECAYS
2280  al1 = vpffcoupling(idordered(6),-1,.true.)
2281  ar1 = vpffcoupling(idordered(6),+1,.true.)
2282  if( isaquark(idordered(8)) ) then
2283  al2 = bl * dsqrt(scale_alpha_w_ud)
2284  ar2 = br * dsqrt(scale_alpha_w_ud)! = 0
2285  elseif( &
2286  (abs(idordered(8)).eq.abs(elm_) .and. abs(idordered(9)).eq.abs(anue_)) .or. (abs(idordered(9)).eq.abs(elm_) .and. abs(idordered(8)).eq.abs(anue_)) .or. &
2287  (abs(idordered(8)).eq.abs(mum_) .and. abs(idordered(9)).eq.abs(anum_)) .or. (abs(idordered(9)).eq.abs(mum_) .and. abs(idordered(8)).eq.abs(anum_)) &
2288  ) then
2289  al2 = bl * dsqrt(scale_alpha_w_ln)
2290  ar2 = br * dsqrt(scale_alpha_w_ln)! = 0
2291  elseif( &
2292  (abs(idordered(8)).eq.abs(tam_) .and. abs(idordered(9)).eq.abs(anut_)) .or. (abs(idordered(9)).eq.abs(tam_) .and. abs(idordered(8)).eq.abs(anut_)) &
2293  ) then
2294  al2 = bl * dsqrt(scale_alpha_w_tn)
2295  ar2 = br * dsqrt(scale_alpha_w_tn)! = 0
2296  else
2297  al2=0d0
2298  ar2=0d0
2299  endif
2300  pv(3,:) = pordered(:,6)+pordered(:,7)
2301  pv(4,:) = pordered(:,8)+pordered(:,9)
2302  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2303  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2304  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2305  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
2306  s = scr(pv(3,:),pv(3,:))
2307  if( m_vprime .gt. 0d0 ) then
2308  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2309  elseif( m_vprime .eq. 0d0 ) then
2310  propv(1) = 1d0
2311  else
2312  propv(1) = s/m_w**2
2313  endif
2314  s = scr(pv(4,:),pv(4,:))
2315  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
2316 
2317 
2318  elseif( vvmode.eq.wwpmode ) then
2319  ! WW' DECAYS
2320  if( isaquark(idordered(6)) ) then
2321  al1 = bl * dsqrt(scale_alpha_w_ud)
2322  ar1 = br * dsqrt(scale_alpha_w_ud)! = 0
2323  elseif( &
2324  (abs(idordered(6)).eq.abs(elp_) .and. abs(idordered(7)).eq.abs(nue_)) .or. (abs(idordered(7)).eq.abs(elp_) .and. abs(idordered(6)).eq.abs(nue_)) .or. &
2325  (abs(idordered(6)).eq.abs(mup_) .and. abs(idordered(7)).eq.abs(num_)) .or. (abs(idordered(7)).eq.abs(mup_) .and. abs(idordered(6)).eq.abs(num_)) &
2326  ) then
2327  al1 = bl * dsqrt(scale_alpha_w_ln)
2328  ar1 = br * dsqrt(scale_alpha_w_ln)! = 0
2329  elseif( &
2330  (abs(idordered(6)).eq.abs(tap_) .and. abs(idordered(7)).eq.abs(nut_)) .or. (abs(idordered(7)).eq.abs(tap_) .and. abs(idordered(6)).eq.abs(nut_)) &
2331  ) then
2332  al1 = bl * dsqrt(scale_alpha_w_tn)
2333  ar1 = br * dsqrt(scale_alpha_w_tn)! = 0
2334  else
2335  al1=0d0
2336  ar1=0d0
2337  endif
2338  al2 = vpffcoupling(idordered(8),-1,.true.)
2339  ar2 = vpffcoupling(idordered(8),+1,.true.)
2340 
2341  pv(3,:) = pordered(:,6)+pordered(:,7)
2342  pv(4,:) = pordered(:,8)+pordered(:,9)
2343  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2344  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2345  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2346  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
2347  s = scr(pv(3,:),pv(3,:))
2348  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
2349  s = scr(pv(4,:),pv(4,:))
2350  if( m_vprime .gt. 0d0 ) then
2351  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2352  elseif( m_vprime .eq. 0d0 ) then
2353  propv(2) = 1d0
2354  else
2355  propv(2) = s/m_w**2
2356  endif
2357 
2358 
2359  elseif( vvmode.eq.wpwpmode ) then
2360  ! W'W' DECAYS
2361  al1 = vpffcoupling(idordered(6),-1,.true.)
2362  ar1 = vpffcoupling(idordered(6),+1,.true.)
2363  al2 = vpffcoupling(idordered(8),-1,.true.)
2364  ar2 = vpffcoupling(idordered(8),+1,.true.)
2365 
2366  pv(3,:) = pordered(:,6)+pordered(:,7)
2367  pv(4,:) = pordered(:,8)+pordered(:,9)
2368  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
2369  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
2370  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
2371  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
2372  s = scr(pv(3,:),pv(3,:))
2373  if( m_vprime .gt. 0d0 ) then
2374  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2375  elseif( m_vprime .eq. 0d0 ) then
2376  propv(1) = 1d0
2377  else
2378  propv(1) = s/m_w**2
2379  endif
2380  s = scr(pv(4,:),pv(4,:))
2381  if( m_vprime .gt. 0d0 ) then
2382  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
2383  elseif( m_vprime .eq. 0d0 ) then
2384  propv(2) = 1d0
2385  else
2386  propv(2) = s/m_w**2
2387  endif
2388 
2389 
2390  else
2391  call error("Unsupported decay modes")
2392  endif
2393 
2394  sp(3,:) = sp(3,:)*propv(1)
2395  sp(4,:) = sp(4,:)*propv(2)
2396  if (h3.eq.-1) then
2397  sp(3,:) = al1 * sp(3,:)
2398  elseif(h3.eq.1) then
2399  sp(3,:) = ar1 * sp(3,:)
2400  endif
2401  if (h4.eq.-1) then
2402  sp(4,:) = al2 * sp(4,:)
2403  elseif(h4.eq.1) then
2404  sp(4,:) = ar2 * sp(4,:)
2405  endif
2406 
2407  return
2408 end subroutine
2409 
2410 subroutine getdecay_vvmode_ordering(MY_IDUP, VVMode,ordering,VVMode_swap,ordering_swap)
2411  implicit none
2412  integer, intent(in) :: MY_IDUP(6:9)
2413  integer, intent(out) :: VVMode,ordering(1:4),VVMode_swap,ordering_swap(1:4)
2414  integer :: idV(1:2),idV_swap(1:2)
2415 
2416  ordering=(/3,4,5,6/)
2417  idv(1)=coupledvertex(my_idup(6:7),-1)
2418  idv(2)=coupledvertex(my_idup(8:9),-1)
2419  if(my_idup(6).eq.pho_ .or. my_idup(7).eq.pho_) idv(1)=pho_
2420  if(my_idup(8).eq.pho_ .or. my_idup(9).eq.pho_) idv(2)=pho_
2421  if(convertlhe(my_idup(6)).lt.0 .or. my_idup(6).eq.not_a_particle_) then
2422  call swap(ordering(1),ordering(2))
2423  endif
2424  if(convertlhe(my_idup(8)).lt.0 .or. my_idup(8).eq.not_a_particle_) then
2425  call swap(ordering(3),ordering(4))
2426  endif
2427  if( &
2428  (idv(1).eq.wm_ .and. idv(2).eq.wp_) .or. &
2429  (idv(2).eq.z0_ .and. idv(1).eq.pho_) &
2430  ) then
2431  call swap(ordering(1),ordering(3))
2432  call swap(ordering(2),ordering(4))
2433  call swap(idv(1),idv(2))
2434  endif
2435  ordering_swap(:)=ordering(:)
2436  call swap(ordering_swap(1),ordering_swap(3))
2437 
2438  idv_swap(1) = coupledvertex( (/ my_idup(3+ordering_swap(1)), my_idup(3+ordering_swap(2)) /), -1)
2439  idv_swap(2) = coupledvertex( (/ my_idup(3+ordering_swap(3)), my_idup(3+ordering_swap(4)) /), -1)
2440  if ( (idv_swap(1).eq.wm_) .and. (idv_swap(2).eq.wp_) ) then
2441  call swap(ordering_swap(1),ordering_swap(3))
2442  call swap(ordering_swap(2),ordering_swap(4))
2443  call swap(idv_swap(1),idv_swap(2))
2444  endif
2445 
2446  if(idv(1).eq.z0_ .and. idv(2).eq.z0_) then
2447  vvmode=zzmode
2448  elseif(idv(1).eq.z0_ .and. idv(2).eq.pho_) then
2449  vvmode=zgmode
2450  elseif(idv(1).eq.pho_ .and. idv(2).eq.pho_) then
2451  vvmode=ggmode
2452  elseif(idv(1).eq.wp_ .and. idv(2).eq.wm_) then
2453  vvmode=wwmode
2454  else
2455  print *,"idV=",idv
2456  call error("Unsupported decay mode")
2457  endif
2458 
2459  vvmode_swap=invalidmode
2460  if(idv_swap(1).eq.z0_ .and. idv_swap(2).eq.z0_) then
2461  vvmode_swap=zzmode
2462  elseif(idv_swap(1).eq.wp_ .and. idv_swap(2).eq.wm_) then
2463  vvmode_swap=wwmode
2464  endif
2465  return
2466 end subroutine
2467 
2468 
2469  end module
2470 
2471 
modparameters::b5
complex(8), public b5
Definition: mod_Parameters.F90:940
modparameters::b1
complex(8), public b1
Definition: mod_Parameters.F90:936
modparameters::bzpgs8
complex(8), public bzpgs8
Definition: mod_Parameters.F90:985
modparameters::m_vprime
real(8), public m_vprime
Definition: mod_Parameters.F90:78
modparameters::includegammastar
logical, public includegammastar
Definition: mod_Parameters.F90:213
modparameters::bzpgs3
complex(8), public bzpgs3
Definition: mod_Parameters.F90:983
modparameters::a5
complex(8), public a5
Definition: mod_Parameters.F90:927
modgraviton::getdecay_couplings_spinors_props
subroutine getdecay_couplings_spinors_props(VVMode, idordered, pordered, h3
Definition: mod_Graviton.F90:1699
modparameters::elm_
integer, target, public elm_
Definition: mod_Parameters.F90:1112
modparameters::zpgsmode
integer, parameter, public zpgsmode
Definition: mod_Parameters.F90:16
modparameters::mphotoncutoff
real(8), public mphotoncutoff
Definition: mod_Parameters.F90:215
modparameters::bzpzp1
complex(8), public bzpzp1
Definition: mod_Parameters.F90:958
modparameters::coupledvertex
integer function coupledvertex(id, hel, useAHcoupl)
Definition: mod_Parameters.F90:2467
modparameters::c7
complex(8), parameter, public c7
Definition: mod_Parameters.F90:995
modparameters::bzpzp6
complex(8), public bzpzp6
Definition: mod_Parameters.F90:963
modparameters::zgmode
integer, parameter, public zgmode
Definition: mod_Parameters.F90:13
modparameters::bzgs1
complex(8), public bzgs1
Definition: mod_Parameters.F90:969
modparameters::cl_neu
real(8), public cl_neu
Definition: mod_Parameters.F90:1071
modparameters::b10
complex(8), public b10
Definition: mod_Parameters.F90:945
modparameters::b8
complex(8), public b8
Definition: mod_Parameters.F90:943
modparameters::bzpgs2
complex(8), public bzpgs2
Definition: mod_Parameters.F90:982
modmisc::et1
double complex function et1(e1, e2, e3, e4)
Definition: mod_Misc.F90:109
modparameters::a2
complex(8), public a2
Definition: mod_Parameters.F90:924
modparameters::bzpzp8
complex(8), public bzpzp8
Definition: mod_Parameters.F90:965
modparameters::dn_
integer, target, public dn_
Definition: mod_Parameters.F90:1085
modparameters::bzgs4
complex(8), public bzgs4
Definition: mod_Parameters.F90:972
modparameters::bzzp2
complex(8), public bzzp2
Definition: mod_Parameters.F90:948
modparameters::pol_mless2
complex(dp) function, dimension(4) pol_mless2(p, i, out)
Definition: mod_Parameters.F90:3053
modgraviton::evalamp_gg_g_vv
subroutine, public evalamp_gg_g_vv(p, MY_IDUP, res)
Definition: mod_Graviton.F90:18
modparameters::scale_alpha_w_tn
real(8), public scale_alpha_w_tn
Definition: mod_Parameters.F90:335
modparameters::zpzpmode
integer, parameter, public zpzpmode
Definition: mod_Parameters.F90:15
modgraviton::calchelamp_gg
subroutine calchelamp_gg(ordering, VVMode, p, MY_IDUP, i1, i2, i3, i4, A)
Definition: mod_Graviton.F90:121
modparameters::one
real(8), parameter, public one
Definition: mod_Parameters.F90:83
modparameters::zgsmode
integer, parameter, public zgsmode
Definition: mod_Parameters.F90:13
modparameters::sqrt2
real(dp), parameter, public sqrt2
Definition: mod_Parameters.F90:80
modparameters::scale_alpha_w_ln
real(8), public scale_alpha_w_ln
Definition: mod_Parameters.F90:334
modparameters::bzpzp3
complex(8), public bzpzp3
Definition: mod_Parameters.F90:960
modparameters::two
real(8), parameter, public two
Definition: mod_Parameters.F90:84
modmisc::error
subroutine error(Message, ErrNum)
Definition: mod_Misc.F90:366
modparameters::includevprime
logical, public includevprime
Definition: mod_Parameters.F90:214
modparameters::c6
complex(8), parameter, public c6
Definition: mod_Parameters.F90:994
modparameters::lambda
real(8), parameter, public lambda
Definition: mod_Parameters.F90:354
modparameters::bgsgs1
complex(8), public bgsgs1
Definition: mod_Parameters.F90:975
modparameters::bgsgs8
complex(8), public bgsgs8
Definition: mod_Parameters.F90:979
modparameters::graviton_qq_left
complex(8), public graviton_qq_left
Definition: mod_Parameters.F90:928
modparameters::cr_qdn
real(8), public cr_qdn
Definition: mod_Parameters.F90:1074
modparameters::bzzp10
complex(8), public bzzp10
Definition: mod_Parameters.F90:956
modparameters::bgsgs4
complex(8), public bgsgs4
Definition: mod_Parameters.F90:978
modparameters::c42
complex(8), parameter, public c42
Definition: mod_Parameters.F90:992
modparameters::pol_mass2
complex(dp) function, dimension(4) pol_mass2(p, i, out)
Definition: mod_Parameters.F90:3176
modmisc::sc
double complex function sc(q1, q2)
Definition: mod_Misc.F90:127
modparameters::cr_neu
real(8), public cr_neu
Definition: mod_Parameters.F90:1070
modparameters::wwpmode
integer, parameter, public wwpmode
Definition: mod_Parameters.F90:14
modparameters::wm_
integer, target, public wm_
Definition: mod_Parameters.F90:1115
modparameters::zero
real(8), parameter, public zero
Definition: mod_Parameters.F90:85
modparameters::chm_
integer, target, public chm_
Definition: mod_Parameters.F90:1086
modparameters::up_
integer, target, public up_
Definition: mod_Parameters.F90:1084
modparameters::bot_
integer, target, public bot_
Definition: mod_Parameters.F90:1089
modparameters::cl_lep
real(8), public cl_lep
Definition: mod_Parameters.F90:1069
modparameters::cr_lep
real(8), public cr_lep
Definition: mod_Parameters.F90:1068
modparameters::a3
complex(8), public a3
Definition: mod_Parameters.F90:925
modparameters::ga_reso
real(8), public ga_reso
Definition: mod_Parameters.F90:231
modparameters::bzzp1
complex(8), public bzzp1
Definition: mod_Parameters.F90:947
modparameters::zzmode
integer, parameter, public zzmode
Definition: mod_Parameters.F90:13
modgraviton
Definition: mod_Graviton.F90:1
modparameters::b9
complex(8), public b9
Definition: mod_Parameters.F90:944
modparameters::ar_qdn
real(8), public ar_qdn
Definition: mod_Parameters.F90:1064
modparameters::bzpzp9
complex(8), public bzpzp9
Definition: mod_Parameters.F90:966
modparameters::c41
complex(8), parameter, public c41
Definition: mod_Parameters.F90:991
modparameters::bzgs2
complex(8), public bzgs2
Definition: mod_Parameters.F90:970
modparameters::bzgs8
complex(8), public bzgs8
Definition: mod_Parameters.F90:973
modparameters::gsgmode
integer, parameter, public gsgmode
Definition: mod_Parameters.F90:13
modparameters::scale_alpha_z_ll
real(8), public scale_alpha_z_ll
Definition: mod_Parameters.F90:329
modparameters::bzpzp5
complex(8), public bzpzp5
Definition: mod_Parameters.F90:962
modparameters::bl
real(8), public bl
Definition: mod_Parameters.F90:1066
modparameters::includeinterference
logical, public includeinterference
Definition: mod_Parameters.F90:77
modparameters::generate_bis
logical, parameter, public generate_bis
Definition: mod_Parameters.F90:932
modparameters::ar_neu
real(8), public ar_neu
Definition: mod_Parameters.F90:1060
modparameters::tam_
integer, target, public tam_
Definition: mod_Parameters.F90:1114
modparameters::scale_alpha_w_ud
real(8), public scale_alpha_w_ud
Definition: mod_Parameters.F90:332
modparameters::c5
complex(8), parameter, public c5
Definition: mod_Parameters.F90:993
modparameters::bzzp5
complex(8), public bzzp5
Definition: mod_Parameters.F90:951
modparameters::bzpgs4
complex(8), public bzpgs4
Definition: mod_Parameters.F90:984
modparameters::al_qdn
real(8), public al_qdn
Definition: mod_Parameters.F90:1065
modgraviton::calchelamp_qq
subroutine calchelamp_qq(ordering, VVMode, p, MY_IDUP, i1, i3, i4, A)
Definition: mod_Graviton.F90:952
modparameters::bzzp3
complex(8), public bzzp3
Definition: mod_Parameters.F90:949
modparameters::cl_qdn
real(8), public cl_qdn
Definition: mod_Parameters.F90:1075
modparameters::a4
complex(8), public a4
Definition: mod_Parameters.F90:926
modparameters::nue_
integer, target, public nue_
Definition: mod_Parameters.F90:1097
modparameters::bzgs3
complex(8), public bzgs3
Definition: mod_Parameters.F90:971
modparameters::ar_lep
real(8), public ar_lep
Definition: mod_Parameters.F90:1058
modparameters::b7
complex(8), public b7
Definition: mod_Parameters.F90:942
modgraviton::gzzampl
subroutine gzzampl(VVMode, p, sp, i1, res)
Definition: mod_Graviton.F90:1394
modparameters::bzzp8
complex(8), public bzzp8
Definition: mod_Parameters.F90:954
modparameters::zzpmode
integer, parameter, public zzpmode
Definition: mod_Parameters.F90:15
modparameters::czero
complex(8), parameter, public czero
Definition: mod_Parameters.F90:86
modparameters::gszmode
integer, parameter, public gszmode
Definition: mod_Parameters.F90:13
modparameters::bzpzp4
complex(8), public bzpzp4
Definition: mod_Parameters.F90:961
modparameters::c3
complex(8), parameter, public c3
Definition: mod_Parameters.F90:990
modparameters::elp_
integer, target, public elp_
Definition: mod_Parameters.F90:1090
modparameters::zpzmode
integer, parameter, public zpzmode
Definition: mod_Parameters.F90:15
modparameters::bzzp7
complex(8), public bzzp7
Definition: mod_Parameters.F90:953
modparameters::m_w
real(8), public m_w
Definition: mod_Parameters.F90:228
modparameters::bzpzp7
complex(8), public bzpzp7
Definition: mod_Parameters.F90:964
modparameters::scale_alpha_z_tt
real(8), public scale_alpha_z_tt
Definition: mod_Parameters.F90:330
modparameters::m_v
real(8), public m_v
Definition: mod_Parameters.F90:78
modparameters::anum_
integer, target, public anum_
Definition: mod_Parameters.F90:1117
modparameters::bzzp4
complex(8), public bzzp4
Definition: mod_Parameters.F90:950
modparameters::al_qup
real(8), public al_qup
Definition: mod_Parameters.F90:1063
modparameters::bzzp6
complex(8), public bzzp6
Definition: mod_Parameters.F90:952
modparameters::ga_v
real(8), public ga_v
Definition: mod_Parameters.F90:78
modparameters::anut_
integer, target, public anut_
Definition: mod_Parameters.F90:1118
modparameters::pol_dk2mom
complex(dp) function, dimension(4) pol_dk2mom(plepton, antilepton, i, outgoing)
Definition: mod_Parameters.F90:3068
modparameters::anue_
integer, target, public anue_
Definition: mod_Parameters.F90:1116
modgraviton::calchelamp2
subroutine calchelamp2(ordering, VVMode, p, MY_IDUP, i1, i3, i4, A)
Definition: mod_Graviton.F90:1366
modparameters
Definition: mod_Parameters.F90:1
modparameters::scale_alpha_z_nn
real(8), public scale_alpha_z_nn
Definition: mod_Parameters.F90:331
modparameters::wpwmode
integer, parameter, public wpwmode
Definition: mod_Parameters.F90:14
modmisc
Definition: mod_Misc.F90:1
modgraviton::qqgzzampl
subroutine qqgzzampl(VVMode, p, sp, res)
Definition: mod_Graviton.F90:994
modparameters::bzpzp10
complex(8), public bzpzp10
Definition: mod_Parameters.F90:967
modparameters::convertlhe
integer function convertlhe(Part)
Definition: mod_Parameters.F90:1798
modparameters::graviton_qq_right
complex(8), public graviton_qq_right
Definition: mod_Parameters.F90:929
modparameters::cr_qup
real(8), public cr_qup
Definition: mod_Parameters.F90:1072
modparameters::bzpzp2
complex(8), public bzpzp2
Definition: mod_Parameters.F90:959
modparameters::a1
complex(8), public a1
Definition: mod_Parameters.F90:923
modparameters::scale_alpha_z_dd
real(8), public scale_alpha_z_dd
Definition: mod_Parameters.F90:328
modparameters::symmfac
real(8), parameter, public symmfac
Definition: mod_Parameters.F90:91
modparameters::b2
complex(8), public b2
Definition: mod_Parameters.F90:937
modparameters::vpffcoupling
complex(8) function vpffcoupling(jhuid, hel, useWp)
Definition: mod_Parameters.F90:1370
modparameters::br
real(8), public br
Definition: mod_Parameters.F90:1067
modparameters::m_z
real(8), public m_z
Definition: mod_Parameters.F90:226
modparameters::nut_
integer, target, public nut_
Definition: mod_Parameters.F90:1099
modparameters::wwmode
integer, parameter, public wwmode
Definition: mod_Parameters.F90:13
modparameters::dp
integer, parameter, public dp
Definition: mod_Parameters.F90:11
modparameters::b6
complex(8), public b6
Definition: mod_Parameters.F90:941
modparameters::mum_
integer, target, public mum_
Definition: mod_Parameters.F90:1113
modparameters::b3
complex(8), public b3
Definition: mod_Parameters.F90:938
modparameters::c1
complex(8), parameter, public c1
Definition: mod_Parameters.F90:988
modparameters::wp_
integer, target, public wp_
Definition: mod_Parameters.F90:1096
modmisc::scr
double precision function scr(p1, p2)
Definition: mod_Misc.F90:135
modgraviton::gggzzampl
subroutine gggzzampl(VVMode, p, sp, res)
Definition: mod_Graviton.F90:156
modparameters::m_reso
real(8), public m_reso
Definition: mod_Parameters.F90:230
modparameters::str_
integer, target, public str_
Definition: mod_Parameters.F90:1087
modparameters::isaquark
logical function isaquark(PartType)
Definition: mod_Parameters.F90:2369
modparameters::bzpgs1
complex(8), public bzpgs1
Definition: mod_Parameters.F90:981
modparameters::scale_alpha_z_uu
real(8), public scale_alpha_z_uu
Definition: mod_Parameters.F90:327
modparameters::zpgmode
integer, parameter, public zpgmode
Definition: mod_Parameters.F90:16
modparameters::not_a_particle_
integer, parameter, public not_a_particle_
Definition: mod_Parameters.F90:1121
modparameters::bgsgs2
complex(8), public bgsgs2
Definition: mod_Parameters.F90:976
modparameters::b4
complex(8), public b4
Definition: mod_Parameters.F90:939
modparameters::ggmode
integer, parameter, public ggmode
Definition: mod_Parameters.F90:13
modparameters::al_lep
real(8), public al_lep
Definition: mod_Parameters.F90:1059
modparameters::num_
integer, target, public num_
Definition: mod_Parameters.F90:1098
modparameters::wpwpmode
integer, parameter, public wpwpmode
Definition: mod_Parameters.F90:14
modparameters::use_dynamic_mg
logical, parameter, public use_dynamic_mg
Definition: mod_Parameters.F90:933
modparameters::cl_qup
real(8), public cl_qup
Definition: mod_Parameters.F90:1073
modparameters::gszpmode
integer, parameter, public gszpmode
Definition: mod_Parameters.F90:16
modparameters::tap_
integer, target, public tap_
Definition: mod_Parameters.F90:1092
modparameters::bzzp9
complex(8), public bzzp9
Definition: mod_Parameters.F90:955
modparameters::invalidmode
integer, parameter, public invalidmode
Definition: mod_Parameters.F90:13
modparameters::gsgsmode
integer, parameter, public gsgsmode
Definition: mod_Parameters.F90:13
modparameters::overallcouplvffsq
real(8), public overallcouplvffsq
Definition: mod_Parameters.F90:1057
modparameters::ga_vprime
real(8), public ga_vprime
Definition: mod_Parameters.F90:78
modparameters::pho_
integer, target, public pho_
Definition: mod_Parameters.F90:1094
modgraviton::getdecay_vvmode_ordering
subroutine getdecay_vvmode_ordering(MY_IDUP, VVMode, ordering, VVMode_swap
Definition: mod_Graviton.F90:2411
modparameters::z0_
integer, target, public z0_
Definition: mod_Parameters.F90:1095
modparameters::bgsgs3
complex(8), public bgsgs3
Definition: mod_Parameters.F90:977
modparameters::c2
complex(8), parameter, public c2
Definition: mod_Parameters.F90:989
modmisc::swap
Definition: mod_Misc.F90:5
modparameters::al_neu
real(8), public al_neu
Definition: mod_Parameters.F90:1061
modparameters::mup_
integer, target, public mup_
Definition: mod_Parameters.F90:1091
modparameters::ar_qup
real(8), public ar_qup
Definition: mod_Parameters.F90:1062
Generated by   doxygen 1.8.17