JHUGen MELA  JHUGen v7.5.6, MELA v2.4.2
Matrix element calculations as used in JHUGen.
Functions/Subroutines
modhiggs Module Reference

Functions/Subroutines

subroutine, public evalamp_gg_h_vv (p, MY_IDUP, res)
 
subroutine, public calchelamp (ordering, VVMode, MY_IDUP, p, i1, i2, i3, i4, A)
 
subroutine gghzzampl (VVMode, p, sp, res)
 
subroutine, public evalamp_h_vv (p, MY_IDUP, res)
 
subroutine, public calchelamp2 (ordering, VVMode, MY_IDUP, p, i3, i4, A)
 
subroutine hzzampl (VVMode, p, sp, res)
 
subroutine, public evalamp_h_ff (pin, mass_F, res)
 
subroutine, public evalamp_h_tt_decay (pin, mass_F, ga_F, res)
 
subroutine getdecay_couplings_spinors_props (VVMode, idordered, pordered, h3
 
subroutine getdecay_vvmode_ordering (MY_IDUP, VVMode, ordering, VVMode_swap
 

Function/Subroutine Documentation

◆ calchelamp()

subroutine, public modhiggs::calchelamp ( integer, dimension(1:4)  ordering,
integer  VVMode,
integer, dimension(6: real(dp) :: p(1:4,1:6) complex(dp) :: propg real(dp) :: pin(4,4) complex(dp) :: a(1:1), sp(4,4) l1=ordering(1) l2=ordering(2) l3=ordering(3) l4=ordering(4) !print *,"l1=",l1 !print *,"l2=",l2 !print *,"l3=",l3 !print *,"l4=",l4 !print *,"id(6)=",convertlhe(my_idup(l1)) !print *,"id(7)=",convertlhe(my_idup(l2)) !print *,"id(8)=",convertlhe(my_idup(l3)) !print *,"id(9)=",convertlhe(my_idup(l4)) !print *,"p(6)=",(p(:,l1)) !print *,"p(7)=",(p(:,l2)) !print *,"p(8)=",(p(:,l3)) !print *,"p(9)=",(p(:,l4)) propg = one/dcmplx(2d0*scr(p(:,1),p(:,2)) - m_reso**2,m_reso*ga_re pin(1,:) = p(:,1) pin(2,:) = p(:,2) sp(1,:) = pol_mless2(dcmplx(p(:,1)),-3+2*i1,'in')  MY_IDUP,
  p,
integer  i1,
integer  i2,
integer  i3,
integer  i4,
  A 
)

Definition at line 168 of file mod_Higgs.F90.

168  implicit none
169  integer :: ordering(1:4),VVMode,i1,i2,i3,i4,l1,l2,l3,l4,MY_IDUP(6:9)
170  real(dp) :: p(1:4,1:6)
171  complex(dp) :: propG
172  real(dp) :: pin(4,4)
173  complex(dp) :: A(1:1), sp(4,4)
174 
175  l1=ordering(1)
176  l2=ordering(2)
177  l3=ordering(3)
178  l4=ordering(4)
179 
180  !print *,"l1=",l1
181  !print *,"l2=",l2

◆ calchelamp2()

subroutine, public modhiggs::calchelamp2 ( integer, dimension(1:4)  ordering,
integer  VVMode,
integer, dimension(6:9)  MY_IDUP,
real(dp), dimension(1:4,1:6)  p,
integer  i3,
integer  i4,
complex(dp), dimension(1:1)  A 
)

Definition at line 482 of file mod_Higgs.F90.

482 
483 
484 ! MADGRAPH CHECK
485 ! call coupsm(0)
486 ! if( (MY_IDUP(6).eq.MY_IDUP(8)) .and. (MY_IDUP(7).eq.MY_IDUP(9)) ) then
487 ! call SH_EMEPEMEP((/-P(1:4,1)-P(1:4,2),P(1:4,3),P(1:4,4),P(1:4,5),P(1:4,6)/)/GeV,res2)
488 ! else
489 ! call SH_EMEPEMEP_NOINT((/-P(1:4,1)-P(1:4,2),P(1:4,3),P(1:4,4),P(1:4,5),P(1:4,6)/)/GeV,res2)
490 ! endif
491 ! res2=res2* cdabs( (0d0,1d0)/dcmplx(2d0*scr(p(:,1),p(:,2))-M_Reso**2,M_Reso*Ga_Reso) * dconjg((0d0,1d0)/dcmplx(2d0*scr(p(:,1),p(:,2))-M_Reso**2,M_Reso*Ga_Reso)) )
492 ! res2=res2*GeV**4
493 ! pause
494 ! res=res2; RETURN
495 
496 ! res= res*prefactor
497 ! print *, "checker 1",res
498 ! print *, "checker 2",res2
499 ! print *, "checker 1/2",res/res2
500 ! pause
501 
502 

◆ evalamp_gg_h_vv()

subroutine, public modhiggs::evalamp_gg_h_vv ( real(dp), dimension(4,6), intent(in)  p,
integer, dimension(6:9), intent(in)  MY_IDUP,
real(dp), intent(out)  res 
)

Definition at line 21 of file mod_Higgs.F90.

21  implicit none
22  real(dp), intent(out) :: res
23  real(dp), intent(in) :: p(4,6)
24  integer, intent(in) :: MY_IDUP(6:9)
25  complex(dp) :: A_VV(1:18), A0_VV(1:2)
26  integer :: i1,i2,i3,i4,VVMode,VVmode_swap
27  real(dp) :: prefactor!,res2
28  real(dp) :: intcolfac
29  integer :: ordering(1:4),ordering_swap(1:4)
30  logical :: doInterference
31 
32  if(isaquark(my_idup(6)) .and. isaquark(my_idup(8))) then
33  intcolfac=1.0_dp/3.0_dp
34  else
35  intcolfac=1.0_dp
36  endif
37 
38  call getdecay_vvmode_ordering(my_idup(6:9),vvmode,ordering,vvmode_swap,ordering_swap)
39 
40 ! Global normalization
41  if( vvmode.eq.zzmode ) then! Z decay
42  prefactor = 8d0*overallcouplvffsq**2
43  elseif( vvmode.eq.wwmode ) then ! W decay
44  prefactor = 8d0*overallcouplvffsq**2
45  elseif( vvmode.eq.zgmode ) then ! Z+photon "decay"
46  prefactor = 8d0*overallcouplvffsq ! Only single powers
47  elseif( vvmode.eq.ggmode ) then ! photon "decay"
48  prefactor = 8d0
49  else
50  prefactor = 0d0
51  endif
52  prefactor = prefactor * (alphas/(3.0_dp*pi*vev))**2
53 
54 
55 ! ! MADGRAPH CHECK
56 ! res=0d0
57 ! if (MY_IDUP(6).ne.MY_IDUP(8) ) return
58 ! if (MY_IDUP(7).ne.MY_IDUP(9) ) return
59 ! if (MY_IDUP(6).eq.MY_IDUP(8) ) return
60 ! if (MY_IDUP(7).eq.MY_IDUP(9) ) return
61 ! print *, "MY COUPL",al1*dsqrt(overallCouplVffsq),ar1*dsqrt(overallCouplVffsq)
62 ! ar1=-0.158480099490745d0! this is MadGraphs gzl(R) , the MG couplings differ by a global minus sign, relative differences are because of different input parameters
63 ! al1=+0.210150647402957d0! this is MadGraphs gzl(L)
64 ! al2=al1
65 ! ar2=ar1
66 ! print *, "MG COUPL",al1,ar1
67 ! pause
68 
69  res = zero
70  a_vv(:) = 0d0
71  dointerference = includeinterference .and. ( &
72  ((vvmode.eq.zzmode) .and. (vvmode_swap.eq.zzmode)) &
73  )
74  if ( includevprime .and. .not.(vvmode.eq.zzmode .or. vvmode.eq.zgmode .or. vvmode.eq.wwmode) ) then
75  call error("Contact terms only for WW, ZZ or Zg!")
76  endif
77  do i1=1,2; do i2=1,2; do i3=1,2; do i4=1,2! sum over helicities
78  call calchelamp(ordering,vvmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(1))
79  if( vvmode.eq.zzmode ) then
80  if( includegammastar ) then
81  call calchelamp(ordering,zgsmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(3))
82  call calchelamp(ordering,gszmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(5))
83  call calchelamp(ordering,gsgsmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(7))
84  endif
85  if( includevprime ) then
86  call calchelamp(ordering,zzpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(9))
87  call calchelamp(ordering,zpzmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(11))
88  call calchelamp(ordering,zpzpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(13))
89  endif
90  if( includegammastar .and. includevprime ) then
91  call calchelamp(ordering,gszpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(15))
92  call calchelamp(ordering,zpgsmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(17))
93  endif
94  elseif( vvmode.eq.zgmode ) then
95  if(includegammastar) then
96  call calchelamp(ordering,gsgmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(3))
97  endif
98  if( includevprime ) then
99  call calchelamp(ordering,zpgmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(5))
100  endif
101  elseif( vvmode.eq.wwmode .and. includevprime ) then
102  call calchelamp(ordering,wwpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(9))
103  call calchelamp(ordering,wpwmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(11))
104  call calchelamp(ordering,wpwpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(13))
105  endif
106 
107  if( dointerference ) then
108  call calchelamp(ordering_swap,vvmode_swap,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(2))
109  if( includegammastar ) then
110  call calchelamp(ordering_swap,zgsmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(4))
111  call calchelamp(ordering_swap,gszmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(6))
112  call calchelamp(ordering_swap,gsgsmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(8))
113  endif
114  if( includevprime ) then
115  call calchelamp(ordering_swap,zzpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(10))
116  call calchelamp(ordering_swap,zpzmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(12))
117  call calchelamp(ordering_swap,zpzpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(14))
118  endif
119  if( includegammastar .and. includevprime ) then
120  call calchelamp(ordering_swap,gszpmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(16))
121  call calchelamp(ordering_swap,zpgsmode,my_idup,p(1:4,1:6),i1,i2,i3,i4,a_vv(18))
122  endif
123  endif
124 
125  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
126  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
127  res = res + dreal(a0_vv(1)*dconjg(a0_vv(1)))
128  res = res + dreal(a0_vv(2)*dconjg(a0_vv(2)))
129  if( dointerference .and. (i3.eq.i4) ) then! interfere the 3456 with 5436 pieces
130  res = res - 2d0*intcolfac*dreal( a0_vv(1)*dconjg(a0_vv(2)) ) ! minus from Fermi statistics
131  endif
132  enddo; enddo; enddo; enddo
133 
134 
135 ! MADGRAPH CHECK
136 ! call coupsm(0)
137 ! if( (MY_IDUP(6).eq.MY_IDUP(8)) .and. (MY_IDUP(7).eq.MY_IDUP(9)) ) then
138 ! call SH_EMEPEMEP((/-P(1:4,1)-P(1:4,2),P(1:4,3),P(1:4,4),P(1:4,5),P(1:4,6)/)/GeV,res2)
139 ! else
140 ! call SH_EMEPEMEP_NOINT((/-P(1:4,1)-P(1:4,2),P(1:4,3),P(1:4,4),P(1:4,5),P(1:4,6)/)/GeV,res2)
141 ! call SH_TAMTAPTAMTAP_NOINT((/-P(1:4,1)-P(1:4,2),P(1:4,3),P(1:4,4),P(1:4,5),P(1:4,6)/)/GeV,res2)
142 ! endif
143 ! res2=res2* cdabs( (0d0,1d0)/dcmplx(2d0*scr(p(:,1),p(:,2))-M_Reso**2,M_Reso*Ga_Reso) * dconjg((0d0,1d0)/dcmplx(2d0*scr(p(:,1),p(:,2))-M_Reso**2,M_Reso*Ga_Reso)) )
144 ! res2=res2*GeV**4
145 ! pause
146 ! res=res2; RETURN
147 
148 ! res= res*prefactor
149 ! print *, "checker 1",res
150 ! print *, "checker 2",res2
151 ! print *, "checker 1/2",res/res2
152 ! pause
153 
154 
155  res = res*prefactor
156  if( (vvmode.eq.zzmode) .and. dointerference ) res = res * symmfac
157 
158  !print *,"VVmode:",VVmode
159  !print *,"VVmode_swap:",VVmode_swap
160  !print *,"ids:",MY_IDUP
161  !print *,"res:",res
162  !pause
163 
164  RETURN

◆ evalamp_h_ff()

subroutine, public modhiggs::evalamp_h_ff ( real(dp), dimension(1:4,1:2), intent(in)  pin,
real(dp), intent(in)  mass_F,
real(dp), intent(out)  res 
)

Definition at line 641 of file mod_Higgs.F90.

641  yyy3 = ahz3
642  elseif( (vvmode.eq.zgmode) .OR. (vvmode.eq.gszmode) .OR. (vvmode.eq.zgsmode) ) then! decay (Z-photon) OR (gamma*-Z) OR (Z-gamma*)
643  yyy1 = ahz1
644  yyy2 = -2*ahz1*q_q/(q_q-q3_q3)
645  yyy3 = ahz3
646  else
647  yyy1=czero
648  yyy2=czero
649  yyy3=czero
650  print *,"VVMode",vvmode,"not implemented in generate_as"
651  endif
652  endif

◆ evalamp_h_tt_decay()

subroutine, public modhiggs::evalamp_h_tt_decay ( real(dp), dimension(1:4,1:6), intent(in)  pin,
real(dp), intent(in)  mass_F,
real(dp), intent(in)  ga_F,
real(dp), intent(out)  res 
)

Definition at line 664 of file mod_Higgs.F90.

664 ! or H --> tbar(p1) + top(p2)
665 ! with tau/top controlled by value of mass_F
666 ! Since this is an isotropic scalar decay, just provide a matrix element SQUARED
667 ! R.Rontsch July 2015
668  SUBROUTINE evalamp_h_ff(pin,mass_F,res)
669  implicit none
670  real(dp), intent(out) :: res
671  complex(dp) :: amp2
672  real(dp), intent(in) :: pin(1:4,1:2),mass_F
673  real(dp) :: s12
674 
675  s12=2d0*(pin(1,1)*pin(1,2)-pin(2,1)*pin(2,2)-pin(3,1)*pin(3,2)-pin(4,1)*pin(4,2)) + 2d0*mass_f**2
676  amp2 = 2d0*s12*(kappa_tilde**2 + kappa**2) - 8d0*mass_f**2*kappa**2
677  amp2 = amp2*mass_f**2/vev**2
678  res = cdabs(amp2)
679 
680  RETURN
681  END SUBROUTINE
682 
683 
684 ! Decay amplitude H --> tau^-(-->l^-(p1)+nubar(p2)+nutau(p3)) + tau^+(-->nu(p4)+l^+(p5)+nutaubar(p6))
685 ! or H --> tbar^-(-->l^-(p1)+nubar(p2)+bbar(p3))+top(-->nu(p4)+l^+(p5)+b(p6))
686 ! Note: 1. Value of m allows to change between H->tau^+tau^- and H->ttbar
687 ! 2. Full kinematics of decaying particles -- no NWA.
688 ! 3. Final state b quark is massless.
689 ! 4. At present, no width in the tau/top propagator
690 ! R. Rontsch July 2015
691  SUBROUTINE evalamp_h_tt_decay(pin,mass_F,ga_F,res)
692  implicit none
693  real(dp), intent(out) :: res
694  real(dp), intent(in) :: pin(1:4,1:6),mass_F,ga_F

◆ evalamp_h_vv()

subroutine, public modhiggs::evalamp_h_vv ( real(dp), dimension(4,6), intent(in)  p,
integer, dimension(6:9), intent(in)  MY_IDUP,
real(dp), intent(out)  res 
)

Definition at line 341 of file mod_Higgs.F90.

341  yyy2 = -2*ahz1 !ahz2 ! gauge invariance fixes ahz2 in this case
342  yyy3 = ahz3
343  elseif( (vvmode.eq.zgmode) .OR. (vvmode.eq.gszmode) .OR. (vvmode.eq.zgsmode) ) then! decay (Z-photon) OR (gamma*-Z) OR (Z-gamma*)
344  yyy1 = ahz1
345  yyy2 = -2*ahz1*q_q/(q_q-q3_q3)
346  yyy3 = ahz3
347  else
348  yyy1=czero
349  yyy2=czero
350  yyy3=czero
351  print *,"VVMode",vvmode,"not implemented in generate_as"
352  endif
353  endif
354 
355  res= e1_e2*e3_e4*q_q**2*yyy1*xxx1 &
356  + e1_e2*e3_q4*e4_q3*q_q*yyy2*xxx1 &
357  + et1(e1,e2,q1,q2)*e3_e4*q_q*yyy1*xxx3 &
358  + et1(e1,e2,q1,q2)*e3_q4*e4_q3*yyy2*xxx3 &
359  + et1(e1,e2,q1,q2)*et1(e3,e4,q3,q4)*yyy3*xxx3 &
360  + et1(e3,e4,q3,q4)*e1_e2*q_q*yyy3*xxx1
361 
362  END SUBROUTINE gghzzampl
363 
364 
365 !----- a subroutinefor H -> (Z+gamma*)(Z+gamma*)/WW/gammagamma
366 !----- all outgoing convention and the following momentum assignment
367 !----- 0 -> Higgs(p1) + e-(p3) + e+(p4) +mu-(p5) +mu+(p6)
368  subroutine evalamp_h_vv(p,MY_IDUP,res)
369  implicit none
370  real(dp), intent(out) :: res
371  real(dp), intent(in) :: p(4,6)
372  integer, intent(in) :: MY_IDUP(6:9)
373  complex(dp) :: A_VV(1:18),A0_VV(1:2)
374  integer :: i3,i4,VVMode,VVmode_swap
375  real(dp) :: prefactor
376  real(dp) :: intcolfac
377  integer :: ordering(1:4),ordering_swap(1:4)
378  logical :: doInterference
379 
380  if(isaquark(my_idup(6)) .and. isaquark(my_idup(8))) then
381  intcolfac=1.0_dp/3.0_dp
382  else
383  intcolfac=1.0_dp
384  endif
385 
386  call getdecay_vvmode_ordering(my_idup(6:9),vvmode,ordering,vvmode_swap,ordering_swap)
387 
388 ! Global normalization
389  if( vvmode.eq.zzmode ) then! Z decay
390  prefactor = overallcouplvffsq**2
391  elseif( vvmode.eq.wwmode ) then ! W decay
392  prefactor = overallcouplvffsq**2
393  elseif( vvmode.eq.zgmode ) then ! Z+photon "decay"
394  prefactor = overallcouplvffsq ! Only single powers
395  elseif( vvmode.eq.ggmode ) then ! photon "decay"
396  prefactor = 1d0
397  else
398  prefactor = 0d0
399  endif
400 
401 ! ! MADGRAPH CHECK
402 ! res=0d0
403 ! if (MY_IDUP(6).ne.MY_IDUP(8) ) return
404 ! if (MY_IDUP(7).ne.MY_IDUP(9) ) return
405 ! if (MY_IDUP(6).eq.MY_IDUP(8) ) return
406 ! if (MY_IDUP(7).eq.MY_IDUP(9) ) return
407 ! print *, "MY COUPL",al1*dsqrt(overallCouplVffsq),ar1*dsqrt(overallCouplVffsq)
408 ! ar1=-0.158480099490745d0! this is MadGraphs gzl(R) , the MG couplings differ by a global minus sign, relative differences are because of different input parameters
409 ! al1=+0.210150647402957d0! this is MadGraphs gzl(L)
410 ! al2=al1
411 ! ar2=ar1
412 ! print *, "MG COUPL",al1,ar1
413 ! pause
414 
415 
416  res = zero
417  a_vv(:) = 0d0
418  dointerference = includeinterference .and. ( &
419  ((vvmode.eq.zzmode) .and. (vvmode_swap.eq.zzmode)) &
420  )
421  if ( includevprime .and. .not.(vvmode.eq.zzmode .or. vvmode.eq.zgmode .or. vvmode.eq.wwmode) ) then
422  call error("Contact terms only for WW, ZZ or Zg!")
423  endif
424  do i3=1,2; do i4=1,2! sum over helicities
425  call calchelamp2(ordering,vvmode,my_idup,p(1:4,1:6),i3,i4,a_vv(1))
426  if( vvmode.eq.zzmode ) then
427  if( includegammastar ) then
428  call calchelamp2(ordering,zgsmode,my_idup,p(1:4,1:6),i3,i4,a_vv(3))
429  call calchelamp2(ordering,gszmode,my_idup,p(1:4,1:6),i3,i4,a_vv(5))
430  call calchelamp2(ordering,gsgsmode,my_idup,p(1:4,1:6),i3,i4,a_vv(7))
431  endif
432  if( includevprime ) then
433  call calchelamp2(ordering,zzpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(9))
434  call calchelamp2(ordering,zpzmode,my_idup,p(1:4,1:6),i3,i4,a_vv(11))
435  call calchelamp2(ordering,zpzpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(13))
436  endif
437  if( includegammastar .and. includevprime ) then
438  call calchelamp2(ordering,gszpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(15))
439  call calchelamp2(ordering,zpgsmode,my_idup,p(1:4,1:6),i3,i4,a_vv(17))
440  endif
441  elseif( vvmode.eq.zgmode ) then
442  if(includegammastar) then
443  call calchelamp2(ordering,gsgmode,my_idup,p(1:4,1:6),i3,i4,a_vv(3))
444  endif
445  if( includevprime ) then
446  call calchelamp2(ordering,zpgmode,my_idup,p(1:4,1:6),i3,i4,a_vv(5))
447  endif
448  elseif( vvmode.eq.wwmode .and. includevprime ) then
449  call calchelamp2(ordering,wwpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(9))
450  call calchelamp2(ordering,wpwmode,my_idup,p(1:4,1:6),i3,i4,a_vv(11))
451  call calchelamp2(ordering,wpwpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(13))
452  endif
453 
454  if( dointerference ) then
455  call calchelamp2(ordering_swap,vvmode,my_idup,p(1:4,1:6),i3,i4,a_vv(2))
456  if( includegammastar ) then
457  call calchelamp2(ordering_swap,zgsmode,my_idup,p(1:4,1:6),i3,i4,a_vv(4))
458  call calchelamp2(ordering_swap,gszmode,my_idup,p(1:4,1:6),i3,i4,a_vv(6))
459  call calchelamp2(ordering_swap,gsgsmode,my_idup,p(1:4,1:6),i3,i4,a_vv(8))
460  endif
461  if( includevprime ) then
462  call calchelamp2(ordering_swap,zzpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(10))
463  call calchelamp2(ordering_swap,zpzmode,my_idup,p(1:4,1:6),i3,i4,a_vv(12))
464  call calchelamp2(ordering_swap,zpzpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(14))
465  endif
466  if( includegammastar .and. includevprime ) then
467  call calchelamp2(ordering_swap,gszpmode,my_idup,p(1:4,1:6),i3,i4,a_vv(16))
468  call calchelamp2(ordering_swap,zpgsmode,my_idup,p(1:4,1:6),i3,i4,a_vv(18))
469  endif
470  endif
471 
472  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
473  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
474  res = res + dreal(a0_vv(1)*dconjg(a0_vv(1)))
475  res = res + dreal(a0_vv(2)*dconjg(a0_vv(2)))
476  if( dointerference .and. (i3.eq.i4) ) then! interfere the 3456 with 5436 pieces
477  res = res - 2d0*intcolfac*dreal( a0_vv(1)*dconjg(a0_vv(2)) ) ! minus from Fermi statistics
478  endif

◆ getdecay_couplings_spinors_props()

subroutine modhiggs::getdecay_couplings_spinors_props ( integer, intent(in)  VVMode,
integer, dimension(6:9), intent(in)  idordered,
real(dp), dimension(1:4,6:9), intent(in)  pordered,
integer, intent(in)  h3 
)
private

Definition at line 701 of file mod_Higgs.F90.

701  call convert_to_mcfm(pin(1:4,j),p(j,1:4))
702  enddo
703 ! call spinoru(6,p,za,zb,s)
704  call spinoru(p,za,zb,s)
705 
706  s12=s(1,2)
707  s45=s(4,5)
708  s123=s(1,2)+s(1,3)+s(2,3)
709  s456=s(4,5)+s(4,6)+s(5,6)
710 
711  kl=-mass_f/vev*( kappa -(0d0,1d0)*kappa_tilde )
712  kr=-mass_f/vev*( kappa +(0d0,1d0)*kappa_tilde )
713 
714  amp = + kr * ( za(1,3)*za(1,4)*zb(1,2)*zb(5,6)- za(1,3)*za(3,4)*zb(2,3)*zb(5,6)) &
715  + kl * (- za(1,3)*za(4,5)*zb(2,5)*zb(5,6)- za(1,3)*za(4,6)*zb(2,6)*zb(5,6))
716 
717 ! overall factors and propagators
718  amp=amp/(s123-mass_f**2+ci*mass_f*ga_f)/(s456-mass_f**2+ci*mass_f*ga_f)/(s12-m_w**2+ci*m_w*ga_w)/(s45-m_w**2+ci*m_w*ga_w)
719  amp=amp*16d0*ci*mass_f*gwsq**2
720  res = cdabs(amp)**2
721 
722  RETURN
723  END SUBROUTINE
724 
725 
726 
727 
728 subroutine getdecay_couplings_spinors_props(VVMode,idordered,pordered,h3,h4, sp,pV)
729  implicit none
730  integer, intent(in) :: VVMode,idordered(6:9),h3,h4
731  real(dp), intent(in) :: pordered(1:4,6:9)
732  complex(dp), intent(out) :: sp(3:4,1:4)
733  real(dp), intent(out) :: pV(3:4,1:4)
734  real(dp) :: s
735  complex(dp) :: propV(1:2), aL1,aR1,aL2,aR2
736 
737  ! h3/h4 helicities: -1 == left, 1 == right
738  if( vvmode.eq.zzmode ) then
739  ! ZZ DECAYS
740  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
741  al1=al_lep * dsqrt(scale_alpha_z_ll)
742  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
743  elseif( abs(idordered(6)).eq.abs(tam_) ) then
744  al1=al_lep * dsqrt(scale_alpha_z_tt)
745  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
746  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
747  al1=al_neu * dsqrt(scale_alpha_z_nn)
748  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
749  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
750  al1=al_qup * dsqrt(scale_alpha_z_uu)
751  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
752  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
753  al1=al_qdn * dsqrt(scale_alpha_z_dd)
754  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
755  else
756  al1=0d0
757  ar1=0d0
758  endif
759  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
760  al2=al_lep * dsqrt(scale_alpha_z_ll)
761  ar2=ar_lep * dsqrt(scale_alpha_z_ll)
762  elseif( abs(idordered(8)).eq.abs(tam_) ) then
763  al2=al_lep * dsqrt(scale_alpha_z_tt)
764  ar2=ar_lep * dsqrt(scale_alpha_z_tt)
765  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
766  al2=al_neu * dsqrt(scale_alpha_z_nn)
767  ar2=ar_neu * dsqrt(scale_alpha_z_nn)
768  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
769  al2=al_qup * dsqrt(scale_alpha_z_uu)
770  ar2=ar_qup * dsqrt(scale_alpha_z_uu)
771  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
772  al2=al_qdn * dsqrt(scale_alpha_z_dd)
773  ar2=ar_qdn * dsqrt(scale_alpha_z_dd)
774  else
775  al2=0d0
776  ar2=0d0
777  endif
778  pv(3,:) = pordered(:,6)+pordered(:,7)
779  pv(4,:) = pordered(:,8)+pordered(:,9)
780  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
781  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
782  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
783  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
784  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
785  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
786  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
787  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
788 
789 
790 
791  elseif( vvmode.eq.wwmode ) then
792  ! WW DECAYS
793  if( isaquark(idordered(6)) ) then
794  al1 = bl * dsqrt(scale_alpha_w_ud)
795  ar1 = br * dsqrt(scale_alpha_w_ud)! = 0
796  elseif( &
797  (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. &
798  (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_)) &
799  ) then
800  al1 = bl * dsqrt(scale_alpha_w_ln)
801  ar1 = br * dsqrt(scale_alpha_w_ln)! = 0
802  elseif( &
803  (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_)) &
804  ) then
805  al1 = bl * dsqrt(scale_alpha_w_tn)
806  ar1 = br * dsqrt(scale_alpha_w_tn)! = 0
807  else
808  al1=0d0
809  ar1=0d0
810  endif
811  if( isaquark(idordered(8)) ) then
812  al2 = bl * dsqrt(scale_alpha_w_ud)
813  ar2 = br * dsqrt(scale_alpha_w_ud)! = 0
814  elseif( &
815  (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. &
816  (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_)) &
817  ) then
818  al2 = bl * dsqrt(scale_alpha_w_ln)
819  ar2 = br * dsqrt(scale_alpha_w_ln)! = 0
820  elseif( &
821  (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_)) &
822  ) then
823  al2 = bl * dsqrt(scale_alpha_w_tn)
824  ar2 = br * dsqrt(scale_alpha_w_tn)! = 0
825  else
826  al2=0d0
827  ar2=0d0
828  endif
829  pv(3,:) = pordered(:,6)+pordered(:,7)
830  pv(4,:) = pordered(:,8)+pordered(:,9)
831  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
832  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
833  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
834  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
835  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
836  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
837  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
838  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
839 
840 
841  elseif( vvmode.eq.zgmode ) then
842  ! Zgamma DECAYS
843  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
844  al1=al_lep * dsqrt(scale_alpha_z_ll)
845  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
846  elseif( abs(idordered(6)).eq.abs(tam_) ) then
847  al1=al_lep * dsqrt(scale_alpha_z_tt)
848  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
849  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
850  al1=al_neu * dsqrt(scale_alpha_z_nn)
851  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
852  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
853  al1=al_qup * dsqrt(scale_alpha_z_uu)
854  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
855  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
856  al1=al_qdn * dsqrt(scale_alpha_z_dd)
857  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
858  else
859  al1=0d0
860  ar1=0d0
861  endif
862  al2=1d0
863  ar2=1d0
864  pv(3,:) = pordered(:,6)+pordered(:,7)
865  pv(4,:) = pordered(:,8)
866  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
867  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
868  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
869  !sp(4,1:4)=pV(4,1:4); print *, "this checks gauge invariance"
870  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
871  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
872  propv(2)=1d0
873 
874 
875  elseif( vvmode.eq.ggmode ) then
876  ! gamma gamma DECAYS
877  al1=1d0
878  ar1=1d0
879  al2=1d0
880  ar2=1d0
881  pv(3,:) = pordered(:,6)
882  pv(4,:) = pordered(:,8)
883  sp(3,:) = pol_mless2(dcmplx(pordered(:,6)),h3,'out') ! photon
884  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
885  !sp(3,1:4)=pV(3,1:4); print *, "this checks gauge invariance"
886  !sp(4,1:4)=pV(4,1:4)
887  propv(1)=1d0
888  propv(2)=1d0
889 
890 
891  elseif( vvmode.eq.gsgmode ) then
892  ! gamma* gamma DECAYS
893  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
894  al1=cl_lep * dsqrt(scale_alpha_z_ll)
895  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
896  elseif( abs(idordered(6)).eq.abs(tam_) ) then
897  al1=cl_lep * dsqrt(scale_alpha_z_tt)
898  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
899  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
900  al1=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
901  ar1=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
902  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
903  al1=cl_qup * dsqrt(scale_alpha_z_uu)
904  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
905  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
906  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
907  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
908  else
909  al1=0d0
910  ar1=0d0
911  endif
912  al2=1d0
913  ar2=1d0
914  pv(3,:) = pordered(:,6)+pordered(:,7)
915  pv(4,:) = pordered(:,8)
916  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
917  sp(3,:) = -sp(3,:) ! photon propagator
918  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
919  !sp(4,1:4)=pV(4,1:4); print *, "this checks gauge invariance"
920  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
921  propv(1) = 1d0
922  propv(2) = 1d0
923  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
924 
925 
926  elseif( vvmode.eq.gszmode ) then
927  ! gamma* Z DECAYS
928  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
929  al1=cl_lep * dsqrt(scale_alpha_z_ll)
930  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
931  elseif( abs(idordered(6)).eq.abs(tam_) ) then
932  al1=cl_lep * dsqrt(scale_alpha_z_tt)
933  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
934  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
935  al1=cl_neu * dsqrt(scale_alpha_z_nn)
936  ar1=cr_neu * dsqrt(scale_alpha_z_nn)
937  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
938  al1=cl_qup * dsqrt(scale_alpha_z_uu)
939  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
940  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
941  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
942  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
943  else
944  al1=0d0
945  ar1=0d0
946  endif
947  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
948  al2=al_lep * dsqrt(scale_alpha_z_ll)
949  ar2=ar_lep * dsqrt(scale_alpha_z_ll)
950  elseif( abs(idordered(8)).eq.abs(tam_) ) then
951  al2=al_lep * dsqrt(scale_alpha_z_tt)
952  ar2=ar_lep * dsqrt(scale_alpha_z_tt)
953  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
954  al2=al_neu * dsqrt(scale_alpha_z_nn)
955  ar2=ar_neu * dsqrt(scale_alpha_z_nn)
956  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
957  al2=al_qup * dsqrt(scale_alpha_z_uu)
958  ar2=ar_qup * dsqrt(scale_alpha_z_uu)
959  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
960  al2=al_qdn * dsqrt(scale_alpha_z_dd)
961  ar2=ar_qdn * dsqrt(scale_alpha_z_dd)
962  else
963  al2=0d0
964  ar2=0d0
965  endif
966  pv(3,:) = pordered(:,6)+pordered(:,7)
967  pv(4,:) = pordered(:,8)+pordered(:,9)
968  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
969  sp(3,:) = -sp(3,:)
970  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
971  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
972  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
973  propv(1) = 1d0! = s/dcmplx(s)
974  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
975  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
976  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
977 
978 
979  elseif( vvmode.eq.zgsmode ) then
980  ! Z gamma* DECAYS
981  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
982  al1=al_lep * dsqrt(scale_alpha_z_ll)
983  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
984  elseif( abs(idordered(6)).eq.abs(tam_) ) then
985  al1=al_lep * dsqrt(scale_alpha_z_tt)
986  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
987  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
988  al1=al_neu * dsqrt(scale_alpha_z_nn)
989  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
990  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
991  al1=al_qup * dsqrt(scale_alpha_z_uu)
992  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
993  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
994  al1=al_qdn * dsqrt(scale_alpha_z_dd)
995  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
996  else
997  al1=0d0
998  ar1=0d0
999  endif
1000  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
1001  al2=cl_lep * dsqrt(scale_alpha_z_ll)
1002  ar2=cr_lep * dsqrt(scale_alpha_z_ll)
1003  elseif( abs(idordered(8)).eq.abs(tam_) ) then
1004  al2=cl_lep * dsqrt(scale_alpha_z_tt)
1005  ar2=cr_lep * dsqrt(scale_alpha_z_tt)
1006  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
1007  al2=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
1008  ar2=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
1009  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
1010  al2=cl_qup * dsqrt(scale_alpha_z_uu)
1011  ar2=cr_qup * dsqrt(scale_alpha_z_uu)
1012  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
1013  al2=cl_qdn * dsqrt(scale_alpha_z_dd)
1014  ar2=cr_qdn * dsqrt(scale_alpha_z_dd)
1015  else
1016  al2=0d0
1017  ar2=0d0
1018  endif
1019  pv(3,:) = pordered(:,6)+pordered(:,7)
1020  pv(4,:) = pordered(:,8)+pordered(:,9)
1021  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1022  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1023  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1024  sp(4,:) = -sp(4,:)
1025  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1026  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
1027  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1028  propv(2) = 1d0 ! = s/dcmplx(s)
1029  if( s.lt.mphotoncutoff**2 ) propv(2)=czero
1030 
1031 
1032  elseif( vvmode.eq.gsgsmode ) then
1033  ! gamma* gamma* DECAYS
1034  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1035  al1=cl_lep * dsqrt(scale_alpha_z_ll)
1036  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
1037  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1038  al1=cl_lep * dsqrt(scale_alpha_z_tt)
1039  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
1040  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1041  al1=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
1042  ar1=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
1043  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1044  al1=cl_qup * dsqrt(scale_alpha_z_uu)
1045  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
1046  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1047  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
1048  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
1049  else
1050  al1=0d0
1051  ar1=0d0
1052  endif
1053  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
1054  al2=cl_lep * dsqrt(scale_alpha_z_ll)
1055  ar2=cr_lep * dsqrt(scale_alpha_z_ll)
1056  elseif( abs(idordered(8)).eq.abs(tam_) ) then
1057  al2=cl_lep * dsqrt(scale_alpha_z_tt)
1058  ar2=cr_lep * dsqrt(scale_alpha_z_tt)
1059  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
1060  al2=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
1061  ar2=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
1062  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
1063  al2=cl_qup * dsqrt(scale_alpha_z_uu)
1064  ar2=cr_qup * dsqrt(scale_alpha_z_uu)
1065  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
1066  al2=cl_qdn * dsqrt(scale_alpha_z_dd)
1067  ar2=cr_qdn * dsqrt(scale_alpha_z_dd)
1068  else
1069  al2=0d0
1070  ar2=0d0
1071  endif
1072  pv(3,:) = pordered(:,6)+pordered(:,7)
1073  pv(4,:) = pordered(:,8)+pordered(:,9)
1074  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1075  sp(3,:) = -sp(3,:)
1076  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1077  sp(4,:) = -sp(4,:)
1078  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1079  propv(1) = 1d0 ! = s/dcmplx(s)
1080  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
1081  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1082  propv(2) = 1d0 ! = s/dcmplx(s)
1083  if( s.lt.mphotoncutoff**2 ) propv(2)=czero
1084 
1085 
1086  elseif( vvmode.eq.zpzmode ) then
1087  ! Z'Z DECAYS
1088  al1 = vpffcoupling(idordered(6),-1,.false.)
1089  ar1 = vpffcoupling(idordered(6),+1,.false.)
1090  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
1091  al2=al_lep * dsqrt(scale_alpha_z_ll)
1092  ar2=ar_lep * dsqrt(scale_alpha_z_ll)
1093  elseif( abs(idordered(8)).eq.abs(tam_) ) then
1094  al2=al_lep * dsqrt(scale_alpha_z_tt)
1095  ar2=ar_lep * dsqrt(scale_alpha_z_tt)
1096  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
1097  al2=al_neu * dsqrt(scale_alpha_z_nn)
1098  ar2=ar_neu * dsqrt(scale_alpha_z_nn)
1099  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
1100  al2=al_qup * dsqrt(scale_alpha_z_uu)
1101  ar2=ar_qup * dsqrt(scale_alpha_z_uu)
1102  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
1103  al2=al_qdn * dsqrt(scale_alpha_z_dd)
1104  ar2=ar_qdn * dsqrt(scale_alpha_z_dd)
1105  else
1106  al2=0d0
1107  ar2=0d0
1108  endif
1109 
1110  pv(3,:) = pordered(:,6)+pordered(:,7)
1111  pv(4,:) = pordered(:,8)+pordered(:,9)
1112  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1113  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1114  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1115  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1116  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1117  if( m_vprime .gt. 0d0 ) then
1118  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1119  elseif( m_vprime .eq. 0d0 ) then
1120  propv(1) = 1d0
1121  else
1122  propv(1) = s/m_z**2
1123  endif
1124  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1125  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
1126 
1127 
1128  elseif( vvmode.eq.zzpmode ) then
1129  ! ZZ' DECAYS
1130  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1131  al1=al_lep * dsqrt(scale_alpha_z_ll)
1132  ar1=ar_lep * dsqrt(scale_alpha_z_ll)
1133  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1134  al1=al_lep * dsqrt(scale_alpha_z_tt)
1135  ar1=ar_lep * dsqrt(scale_alpha_z_tt)
1136  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1137  al1=al_neu * dsqrt(scale_alpha_z_nn)
1138  ar1=ar_neu * dsqrt(scale_alpha_z_nn)
1139  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1140  al1=al_qup * dsqrt(scale_alpha_z_uu)
1141  ar1=ar_qup * dsqrt(scale_alpha_z_uu)
1142  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1143  al1=al_qdn * dsqrt(scale_alpha_z_dd)
1144  ar1=ar_qdn * dsqrt(scale_alpha_z_dd)
1145  else
1146  al1=0d0
1147  ar1=0d0
1148  endif
1149  al2 = vpffcoupling(idordered(8),-1,.false.)
1150  ar2 = vpffcoupling(idordered(8),+1,.false.)
1151 
1152  pv(3,:) = pordered(:,6)+pordered(:,7)
1153  pv(4,:) = pordered(:,8)+pordered(:,9)
1154  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1155  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1156  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1157  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1158  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1159  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
1160  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1161  if( m_vprime .gt. 0d0 ) then
1162  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1163  elseif( m_vprime .eq. 0d0 ) then
1164  propv(2) = 1d0
1165  else
1166  propv(2) = s/m_z**2
1167  endif
1168 
1169 
1170  elseif( vvmode.eq.zpzpmode ) then
1171  ! Z'Z' DECAYS
1172  al1 = vpffcoupling(idordered(6),-1,.false.)
1173  ar1 = vpffcoupling(idordered(6),+1,.false.)
1174  al2 = vpffcoupling(idordered(8),-1,.false.)
1175  ar2 = vpffcoupling(idordered(8),+1,.false.)
1176 
1177  pv(3,:) = pordered(:,6)+pordered(:,7)
1178  pv(4,:) = pordered(:,8)+pordered(:,9)
1179  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1180  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1181  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1182  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1183  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1184  if( m_vprime .gt. 0d0 ) then
1185  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1186  elseif( m_vprime .eq. 0d0 ) then
1187  propv(1) = 1d0
1188  else
1189  propv(1) = s/m_z**2
1190  endif
1191  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1192  if( m_vprime .gt. 0d0 ) then
1193  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1194  elseif( m_vprime .eq. 0d0 ) then
1195  propv(2) = 1d0
1196  else
1197  propv(2) = s/m_z**2
1198  endif
1199 
1200 
1201  elseif( vvmode.eq.zpgsmode ) then
1202  ! Z' gamma* DECAYS
1203  al1 = vpffcoupling(idordered(6),-1,.false.)
1204  ar1 = vpffcoupling(idordered(6),+1,.false.)
1205  if( abs(idordered(8)).eq.abs(elm_) .or. abs(idordered(8)).eq.abs(mum_) ) then
1206  al2=cl_lep * dsqrt(scale_alpha_z_ll)
1207  ar2=cr_lep * dsqrt(scale_alpha_z_ll)
1208  elseif( abs(idordered(8)).eq.abs(tam_) ) then
1209  al2=cl_lep * dsqrt(scale_alpha_z_tt)
1210  ar2=cr_lep * dsqrt(scale_alpha_z_tt)
1211  elseif( abs(idordered(8)).eq.abs(nue_) .or. abs(idordered(8)).eq.abs(num_) .or. abs(idordered(8)).eq.abs(nut_) ) then
1212  al2=cl_neu * dsqrt(scale_alpha_z_nn)! = 0
1213  ar2=cr_neu * dsqrt(scale_alpha_z_nn)! = 0
1214  elseif( abs(idordered(8)).eq.abs(up_) .or. abs(idordered(8)).eq.abs(chm_) ) then
1215  al2=cl_qup * dsqrt(scale_alpha_z_uu)
1216  ar2=cr_qup * dsqrt(scale_alpha_z_uu)
1217  elseif( abs(idordered(8)).eq.abs(dn_) .or. abs(idordered(8)).eq.abs(str_) .or. abs(idordered(8)).eq.abs(bot_) ) then
1218  al2=cl_qdn * dsqrt(scale_alpha_z_dd)
1219  ar2=cr_qdn * dsqrt(scale_alpha_z_dd)
1220  else
1221  al2=0d0
1222  ar2=0d0
1223  endif
1224 
1225  pv(3,:) = pordered(:,6)+pordered(:,7)
1226  pv(4,:) = pordered(:,8)+pordered(:,9)
1227  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1228  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1229  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1230  sp(4,:) = -sp(4,:)
1231  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1232  if( m_vprime .gt. 0d0 ) then
1233  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1234  elseif( m_vprime .eq. 0d0 ) then
1235  propv(1) = 1d0
1236  else
1237  propv(1) = s/m_z**2
1238  endif
1239  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1240  propv(2) = 1d0 ! = s/dcmplx(s)
1241  if( s.lt.mphotoncutoff**2 ) propv(2)=czero
1242 
1243  elseif( vvmode.eq.gszpmode ) then
1244  ! gamma* Z' DECAYS
1245  if( abs(idordered(6)).eq.abs(elm_) .or. abs(idordered(6)).eq.abs(mum_) ) then
1246  al1=cl_lep * dsqrt(scale_alpha_z_ll)
1247  ar1=cr_lep * dsqrt(scale_alpha_z_ll)
1248  elseif( abs(idordered(6)).eq.abs(tam_) ) then
1249  al1=cl_lep * dsqrt(scale_alpha_z_tt)
1250  ar1=cr_lep * dsqrt(scale_alpha_z_tt)
1251  elseif( abs(idordered(6)).eq.abs(nue_) .or. abs(idordered(6)).eq.abs(num_) .or. abs(idordered(6)).eq.abs(nut_) ) then
1252  al1=cl_neu * dsqrt(scale_alpha_z_nn)
1253  ar1=cr_neu * dsqrt(scale_alpha_z_nn)
1254  elseif( abs(idordered(6)).eq.abs(up_) .or. abs(idordered(6)).eq.abs(chm_) ) then
1255  al1=cl_qup * dsqrt(scale_alpha_z_uu)
1256  ar1=cr_qup * dsqrt(scale_alpha_z_uu)
1257  elseif( abs(idordered(6)).eq.abs(dn_) .or. abs(idordered(6)).eq.abs(str_) .or. abs(idordered(6)).eq.abs(bot_) ) then
1258  al1=cl_qdn * dsqrt(scale_alpha_z_dd)
1259  ar1=cr_qdn * dsqrt(scale_alpha_z_dd)
1260  else
1261  al1=0d0
1262  ar1=0d0
1263  endif
1264  al2 = vpffcoupling(idordered(8),-1,.false.)
1265  ar2 = vpffcoupling(idordered(8),+1,.false.)
1266 
1267  pv(3,:) = pordered(:,6)+pordered(:,7)
1268  pv(4,:) = pordered(:,8)+pordered(:,9)
1269  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1270  sp(3,:) = -sp(3,:)
1271  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1272  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1273  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1274  propv(1) = 1d0! = s/dcmplx(s)
1275  if( s.lt.mphotoncutoff**2 ) propv(1)=czero
1276  s = scr(pordered(:,8)+pordered(:,9),pordered(:,8)+pordered(:,9))
1277  if( m_vprime .gt. 0d0 ) then
1278  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1279  elseif( m_vprime .eq. 0d0 ) then
1280  propv(2) = 1d0
1281  else
1282  propv(2) = s/m_z**2
1283  endif
1284 
1285 
1286  elseif( vvmode.eq.zpgmode ) then
1287  ! Z' gamma DECAYS
1288  al1 = vpffcoupling(idordered(6),-1,.false.)
1289  ar1 = vpffcoupling(idordered(6),+1,.false.)
1290  al2=1d0
1291  ar2=1d0
1292  pv(3,:) = pordered(:,6)+pordered(:,7)
1293  pv(4,:) = pordered(:,8)
1294  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1295  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1296  sp(4,:) = pol_mless2(dcmplx(pordered(:,8)),h4,'out') ! photon
1297  s = scr(pordered(:,6)+pordered(:,7),pordered(:,6)+pordered(:,7))
1298  if( m_vprime .gt. 0d0 ) then
1299  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1300  elseif( m_vprime .eq. 0d0 ) then
1301  propv(1) = 1d0
1302  else
1303  propv(1) = s/m_z**2
1304  endif
1305  propv(2)=1d0
1306 
1307 
1308  elseif( vvmode.eq.wpwmode ) then
1309  ! W'W DECAYS
1310  al1 = vpffcoupling(idordered(6),-1,.true.)
1311  ar1 = vpffcoupling(idordered(6),+1,.true.)
1312  if( isaquark(idordered(8)) ) then
1313  al2 = bl * dsqrt(scale_alpha_w_ud)
1314  ar2 = br * dsqrt(scale_alpha_w_ud)! = 0
1315  elseif( &
1316  (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. &
1317  (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_)) &
1318  ) then
1319  al2 = bl * dsqrt(scale_alpha_w_ln)
1320  ar2 = br * dsqrt(scale_alpha_w_ln)! = 0
1321  elseif( &
1322  (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_)) &
1323  ) then
1324  al2 = bl * dsqrt(scale_alpha_w_tn)
1325  ar2 = br * dsqrt(scale_alpha_w_tn)! = 0
1326  else
1327  al2=0d0
1328  ar2=0d0
1329  endif
1330  pv(3,:) = pordered(:,6)+pordered(:,7)
1331  pv(4,:) = pordered(:,8)+pordered(:,9)
1332  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1333  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1334  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1335  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1336  s = scr(pv(3,:),pv(3,:))
1337  if( m_vprime .gt. 0d0 ) then
1338  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1339  elseif( m_vprime .eq. 0d0 ) then
1340  propv(1) = 1d0
1341  else
1342  propv(1) = s/m_w**2
1343  endif
1344  s = scr(pv(4,:),pv(4,:))
1345  propv(2) = s/dcmplx(s - m_v**2,m_v*ga_v)
1346 
1347 
1348  elseif( vvmode.eq.wwpmode ) then
1349  ! WW' DECAYS
1350  if( isaquark(idordered(6)) ) then
1351  al1 = bl * dsqrt(scale_alpha_w_ud)
1352  ar1 = br * dsqrt(scale_alpha_w_ud)! = 0
1353  elseif( &
1354  (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. &
1355  (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_)) &
1356  ) then
1357  al1 = bl * dsqrt(scale_alpha_w_ln)
1358  ar1 = br * dsqrt(scale_alpha_w_ln)! = 0
1359  elseif( &
1360  (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_)) &
1361  ) then
1362  al1 = bl * dsqrt(scale_alpha_w_tn)
1363  ar1 = br * dsqrt(scale_alpha_w_tn)! = 0
1364  else
1365  al1=0d0
1366  ar1=0d0
1367  endif
1368  al2 = vpffcoupling(idordered(8),-1,.true.)
1369  ar2 = vpffcoupling(idordered(8),+1,.true.)
1370 
1371  pv(3,:) = pordered(:,6)+pordered(:,7)
1372  pv(4,:) = pordered(:,8)+pordered(:,9)
1373  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1374  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1375  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1376  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1377  s = scr(pv(3,:),pv(3,:))
1378  propv(1) = s/dcmplx(s - m_v**2,m_v*ga_v)
1379  s = scr(pv(4,:),pv(4,:))
1380  if( m_vprime .gt. 0d0 ) then
1381  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1382  elseif( m_vprime .eq. 0d0 ) then
1383  propv(2) = 1d0
1384  else
1385  propv(2) = s/m_w**2
1386  endif
1387 
1388 
1389  elseif( vvmode.eq.wpwpmode ) then
1390  ! W'W' DECAYS
1391  al1 = vpffcoupling(idordered(6),-1,.true.)
1392  ar1 = vpffcoupling(idordered(6),+1,.true.)
1393  al2 = vpffcoupling(idordered(8),-1,.true.)
1394  ar2 = vpffcoupling(idordered(8),+1,.true.)
1395 
1396  pv(3,:) = pordered(:,6)+pordered(:,7)
1397  pv(4,:) = pordered(:,8)+pordered(:,9)
1398  sp(3,:) = pol_dk2mom(dcmplx(pordered(:,6)),dcmplx(pordered(:,7)),h3) ! ubar(l1), v(l2)
1399  sp(3,:) = -sp(3,:) + pv(3,:)*( sc(sp(3,:),dcmplx(pv(3,:))) )/scr(pv(3,:),pv(3,:))! full propagator numerator
1400  sp(4,:) = pol_dk2mom(dcmplx(pordered(:,8)),dcmplx(pordered(:,9)),h4) ! ubar(l3), v(l4)
1401  sp(4,:) = -sp(4,:) + pv(4,:)*( sc(sp(4,:),dcmplx(pv(4,:))) )/scr(pv(4,:),pv(4,:))! full propagator numerator
1402  s = scr(pv(3,:),pv(3,:))
1403  if( m_vprime .gt. 0d0 ) then
1404  propv(1) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1405  elseif( m_vprime .eq. 0d0 ) then
1406  propv(1) = 1d0
1407  else
1408  propv(1) = s/m_w**2
1409  endif
1410  s = scr(pv(4,:),pv(4,:))
1411  if( m_vprime .gt. 0d0 ) then
1412  propv(2) = s/dcmplx(s - m_vprime**2,m_vprime*ga_vprime)
1413  elseif( m_vprime .eq. 0d0 ) then
1414  propv(2) = 1d0
1415  else
1416  propv(2) = s/m_w**2

◆ getdecay_vvmode_ordering()

subroutine modhiggs::getdecay_vvmode_ordering ( integer, dimension(6:9), intent(in)  MY_IDUP,
integer, intent(out)  VVMode,
integer, dimension(1:4), intent(out)  ordering,
integer, intent(out)  VVMode_swap 
)
private

Definition at line 1420 of file mod_Higgs.F90.

1420  else
1421  call error("Unsupported decay modes")
1422  endif
1423 
1424  !print *,"sp(3)=",sp(3,:)
1425  !print *,"propV(1)=",propV(1)
1426  !print *,"aL1,aR1=",aL1,aR1
1427  !print *,"sp(4)=",sp(4,:)
1428  !print *,"propV(2)=",propV(2)
1429  !print *,"aL2,aR2=",aL2,aR2
1430 
1431  sp(3,:) = sp(3,:)*propv(1)
1432  sp(4,:) = sp(4,:)*propv(2)
1433  if (h3.eq.-1) then
1434  sp(3,:) = al1 * sp(3,:)
1435  elseif(h3.eq.1) then
1436  sp(3,:) = ar1 * sp(3,:)
1437  endif
1438  if (h4.eq.-1) then
1439  sp(4,:) = al2 * sp(4,:)
1440  elseif(h4.eq.1) then
1441  sp(4,:) = ar2 * sp(4,:)
1442  endif
1443 
1444  return
1445 end subroutine
1446 
1447 subroutine getdecay_vvmode_ordering(MY_IDUP, VVMode,ordering,VVMode_swap,ordering_swap)
1448  implicit none
1449  integer, intent(in) :: MY_IDUP(6:9)
1450  integer, intent(out) :: VVMode,ordering(1:4),VVMode_swap,ordering_swap(1:4)
1451  integer :: idV(1:2),idV_swap(1:2)
1452 
1453  ordering=(/3,4,5,6/)
1454  idv(1)=coupledvertex(my_idup(6:7),-1)
1455  idv(2)=coupledvertex(my_idup(8:9),-1)
1456  if(my_idup(6).eq.pho_ .or. my_idup(7).eq.pho_) idv(1)=pho_
1457  if(my_idup(8).eq.pho_ .or. my_idup(9).eq.pho_) idv(2)=pho_
1458  if(convertlhe(my_idup(6)).lt.0 .or. my_idup(6).eq.not_a_particle_) then
1459  call swap(ordering(1),ordering(2))
1460  endif
1461  if(convertlhe(my_idup(8)).lt.0 .or. my_idup(8).eq.not_a_particle_) then
1462  call swap(ordering(3),ordering(4))
1463  endif
1464  if( &
1465  (idv(1).eq.wm_ .and. idv(2).eq.wp_) .or. &
1466  (idv(2).eq.z0_ .and. idv(1).eq.pho_) &
1467  ) then
1468  call swap(ordering(1),ordering(3))
1469  call swap(ordering(2),ordering(4))
1470  call swap(idv(1),idv(2))
1471  endif
1472  ordering_swap(:)=ordering(:)
1473  call swap(ordering_swap(1),ordering_swap(3))
1474 

◆ gghzzampl()

subroutine modhiggs::gghzzampl ( integer, intent(in)  VVMode,
real(dp), dimension(4,4), intent(in)  p,
complex(dp), dimension(4,4), intent(in)  sp,
complex(dp), intent(out)  res 
)
private

Definition at line 185 of file mod_Higgs.F90.

185  !print *,"id(7)=",convertLHE(MY_IDUP(l2))
186  !print *,"id(8)=",convertLHE(MY_IDUP(l3))
187  !print *,"id(9)=",convertLHE(MY_IDUP(l4))
188  !print *,"p(6)=",(p(:,l1))
189  !print *,"p(7)=",(p(:,l2))
190  !print *,"p(8)=",(p(:,l3))
191  !print *,"p(9)=",(p(:,l4))
192 
193  propg = one/dcmplx(2d0*scr(p(:,1),p(:,2)) - m_reso**2,m_reso*ga_reso)
194 
195  pin(1,:) = p(:,1)
196  pin(2,:) = p(:,2)
197  sp(1,:) = pol_mless2(dcmplx(p(:,1)),-3+2*i1,'in') ! gluon
198  sp(2,:) = pol_mless2(dcmplx(p(:,2)),-3+2*i2,'in') ! gluon
199  !sp(1,1:4)=pin(1,1:4);print *, "this checks IS gauge invariance"
200  !sp(2,1:4)=pin(2,1:4);print *, "this checks IS gauge invariance"
201  call getdecay_couplings_spinors_props( &
202  vvmode, &
203  (/my_idup(l1+3),my_idup(l2+3),my_idup(l3+3),my_idup(l4+3)/), &
204  (/p(:,l1),p(:,l2),p(:,l3),p(:,l4)/), &
205  -3+2*i3,-3+2*i4, &
206  sp(3:4,:),pin(3:4,:) &
207  )
208  call gghzzampl(vvmode,pin,sp,a(1))
209  a(1) = a(1) * propg
210  end subroutine
211 
212  SUBROUTINE gghzzampl(VVMode,p,sp,res)
213  implicit none
214  integer, intent(in) :: VVMode
215  real(dp), intent(in) :: p(4,4)
216  complex(dp), intent(in) :: sp(4,4)
217  complex(dp), intent(out) :: res
218  complex(dp) :: e1_e2, e1_e3, e1_e4
219  complex(dp) :: e2_e3, e2_e4
220  complex(dp) :: e3_e4
221  complex(dp) :: q1_e3,q1_e4,q2_e3,q2_e4
222  complex(dp) :: e1_q3,e1_q4,e2_q3,e2_q4
223  complex(dp) :: e3_q4,e4_q3
224  complex(dp) :: q1(4),q2(4),q3(4),q4(4),q(4)
225  complex(dp) :: e1(4),e2(4),e3(4),e4(4)
226  complex(dp) :: xxx1,xxx2,xxx3,yyy1,yyy2,yyy3,yyy4
227  complex(dp) :: ghz1_dyn,ghz2_dyn,ghz3_dyn,ghz4_dyn
228  real(dp) :: q34
229  real(dp) :: q_q, q3_q3, q4_q4
230 
231 
232  res = 0d0
233  q1 = dcmplx(p(1,:),0d0)
234  q2 = dcmplx(p(2,:),0d0)
235  q3 = dcmplx(p(3,:),0d0)
236  q4 = dcmplx(p(4,:),0d0)
237 
238  e1 = sp(1,:)
239  e2 = sp(2,:)
240  e3 = sp(3,:)
241  e4 = sp(4,:)
242 
243  q = -q1-q2
244 
245  q_q =sc(q,q)
246  q3_q3 = sc(q3,q3)
247  q4_q4 = sc(q4,q4)
248 
249  e1_e2 = sc(e1,e2)
250  e1_e3 = sc(e1,e3)
251  e1_e4 = sc(e1,e4)
252  e2_e3 = sc(e2,e3)
253  e2_e4 = sc(e2,e4)
254  e3_e4 = sc(e3,e4)
255 
256  q1_e3 = sc(q1,e3)
257  q1_e4 = sc(q1,e4)
258  q2_e3 = sc(q2,e3)
259  q2_e4 = sc(q2,e4)
260  e1_q3 = sc(e1,q3)
261  e1_q4 = sc(e1,q4)
262  e2_q3 = sc(e2,q3)
263  e2_q4 = sc(e2,q4)
264  e3_q4 = sc(e3,q4)
265  e4_q3 = sc(e4,q3)
266 
267 
268  if (q_q.lt.-0.1d0 .or. q3_q3.lt.-0.1d0 .or. q4_q4.lt.-0.1d0) return ! if negative invariant masses return zero
269 
270 !---- data that defines couplings
271  if( (vvmode.eq.zzmode) .or. (vvmode.eq.wwmode) ) then! decay ZZ's or WW's
272  ghz1_dyn = hvvspinzerodynamiccoupling(1,q3_q3,q4_q4,q_q)
273  ghz2_dyn = hvvspinzerodynamiccoupling(2,q3_q3,q4_q4,q_q)
274  ghz3_dyn = hvvspinzerodynamiccoupling(3,q3_q3,q4_q4,q_q)
275  ghz4_dyn = hvvspinzerodynamiccoupling(4,q3_q3,q4_q4,q_q)
276  elseif( (vvmode.eq.gszmode) ) then
277  ghz1_dyn = hvvspinzerodynamiccoupling(5,0d0,q3_q3,q_q)
278  ghz2_dyn = hvvspinzerodynamiccoupling(6,0d0,q3_q3,q_q)
279  ghz3_dyn = hvvspinzerodynamiccoupling(7,0d0,q3_q3,q_q)
280  ghz4_dyn = hvvspinzerodynamiccoupling(8,0d0,q3_q3,q_q)
281  elseif( (vvmode.eq.zgmode) .OR. (vvmode.eq.zgsmode) ) then
282  ghz1_dyn = hvvspinzerodynamiccoupling(5,0d0,q4_q4,q_q)
283  ghz2_dyn = hvvspinzerodynamiccoupling(6,0d0,q4_q4,q_q)
284  ghz3_dyn = hvvspinzerodynamiccoupling(7,0d0,q4_q4,q_q)
285  ghz4_dyn = hvvspinzerodynamiccoupling(8,0d0,q4_q4,q_q)
286  elseif( (vvmode.eq.ggmode) .or. (vvmode.eq.gsgsmode) .or. (vvmode.eq.gsgmode) ) then
287  ghz1_dyn = czero
288  ghz2_dyn = hvvspinzerodynamiccoupling(9,q3_q3,q4_q4,q_q)
289  ghz3_dyn = hvvspinzerodynamiccoupling(10,q3_q3,q4_q4,q_q)
290  ghz4_dyn = hvvspinzerodynamiccoupling(11,q3_q3,q4_q4,q_q)
291  elseif( (vvmode.eq.zzpmode) .or. (vvmode.eq.wwpmode) ) then
292  ghz1_dyn = hvvspinzerodynamiccoupling(12,q3_q3,q4_q4,q_q)
293  ghz2_dyn = hvvspinzerodynamiccoupling(13,q3_q3,q4_q4,q_q)
294  ghz3_dyn = hvvspinzerodynamiccoupling(14,q3_q3,q4_q4,q_q)
295  ghz4_dyn = hvvspinzerodynamiccoupling(15,q3_q3,q4_q4,q_q)
296  elseif( (vvmode.eq.zpzmode) .or. (vvmode.eq.wpwmode) ) then
297  ghz1_dyn = hvvspinzerodynamiccoupling(12,q4_q4,q3_q3,q_q)
298  ghz2_dyn = hvvspinzerodynamiccoupling(13,q4_q4,q3_q3,q_q)
299  ghz3_dyn = hvvspinzerodynamiccoupling(14,q4_q4,q3_q3,q_q)
300  ghz4_dyn = hvvspinzerodynamiccoupling(15,q4_q4,q3_q3,q_q)
301  elseif( (vvmode.eq.zpzpmode) .or. (vvmode.eq.wpwpmode) ) then
302  ghz1_dyn = hvvspinzerodynamiccoupling(16,q3_q3,q4_q4,q_q)
303  ghz2_dyn = hvvspinzerodynamiccoupling(17,q3_q3,q4_q4,q_q)
304  ghz3_dyn = hvvspinzerodynamiccoupling(18,q3_q3,q4_q4,q_q)
305  ghz4_dyn = hvvspinzerodynamiccoupling(19,q3_q3,q4_q4,q_q)
306  elseif( (vvmode.eq.gszpmode) ) then
307  ghz1_dyn = hvvspinzerodynamiccoupling(20,0d0,q3_q3,q_q)
308  ghz2_dyn = hvvspinzerodynamiccoupling(21,0d0,q3_q3,q_q)
309  ghz3_dyn = hvvspinzerodynamiccoupling(22,0d0,q3_q3,q_q)
310  ghz4_dyn = hvvspinzerodynamiccoupling(23,0d0,q3_q3,q_q)
311  elseif( (vvmode.eq.zpgmode) .OR. (vvmode.eq.zpgsmode) ) then
312  ghz1_dyn = hvvspinzerodynamiccoupling(20,0d0,q4_q4,q_q)
313  ghz2_dyn = hvvspinzerodynamiccoupling(21,0d0,q4_q4,q_q)
314  ghz3_dyn = hvvspinzerodynamiccoupling(22,0d0,q4_q4,q_q)
315  ghz4_dyn = hvvspinzerodynamiccoupling(23,0d0,q4_q4,q_q)
316  else
317  ghz1_dyn = czero
318  ghz2_dyn = czero
319  ghz3_dyn = czero
320  ghz4_dyn = czero
321  print *,"VVMode",vvmode,"not implemented"
322  endif
323 
324  if( .not. generate_as ) then
325  xxx1 = ghg2+ghg3/4d0/lambda**2*q_q
326  xxx3 = -2d0*ghg4
327  yyy1 = ghz1_dyn*m_v**2/q_q & ! in this line M_V is indeed correct, not a misprint
328  + ghz2_dyn*(q_q-q3_q3-q4_q4)/q_q &
329  + ghz3_dyn/lambda**2*(q_q-q3_q3-q4_q4)*(q_q-q4_q4-q3_q3)/4d0/q_q
330  yyy2 = -2d0*ghz2_dyn-ghz3_dyn/2d0/lambda**2*(q_q-q3_q3-q4_q4)
331  yyy3 = -2d0*ghz4_dyn
332  else
333  xxx1 = ahg1

◆ hzzampl()

subroutine modhiggs::hzzampl ( integer, intent(in)  VVMode,
real(dp), dimension(4,4), intent(in)  p,
complex(dp), dimension(3:4,4), intent(in)  sp,
complex(dp), intent(out)  res 
)
private

Definition at line 506 of file mod_Higgs.F90.

506  RETURN
507  END SUBROUTINE
508 
509  subroutine calchelamp2(ordering,VVMode,MY_IDUP,p,i3,i4,A)
510  implicit none
511  integer :: ordering(1:4),VVMode,i3,i4,l1,l2,l3,l4,MY_IDUP(6:9)
512  real(dp) :: p(1:4,1:6)
513  real(dp) :: pin(4,4)
514  complex(dp) :: A(1:1), sp(3:4,4)
515 
516  l1=ordering(1)
517  l2=ordering(2)
518  l3=ordering(3)
519  l4=ordering(4)
520 
521  pin(1,:) = p(:,1)
522 
523  call getdecay_couplings_spinors_props( &
524  vvmode, &
525  (/my_idup(l1+3),my_idup(l2+3),my_idup(l3+3),my_idup(l4+3)/), &
526  (/p(:,l1),p(:,l2),p(:,l3),p(:,l4)/), &
527  -3+2*i3,-3+2*i4, &
528  sp(3:4,:),pin(3:4,:) &
529  )
530  call hzzampl(vvmode,pin,sp,a(1))
531  end subroutine
532 
533  subroutine hzzampl(VVMode,p,sp,res)
534  implicit none
535  integer, intent(in) :: VVMode
536  real(dp), intent(in) :: p(4,4)
537  complex(dp), intent(in) :: sp(3:4,4)
538  complex(dp), intent(out) :: res
539  complex(dp) :: e3_e4
540  complex(dp) :: e3_q4,e4_q3
541  complex(dp) :: q1(4),q3(4),q4(4),q(4)
542  complex(dp) :: e3(4),e4(4)
543  complex(dp) :: yyy1,yyy2,yyy3,yyy4
544  real(dp) :: q34
545  real(dp) :: q_q, q3_q3, q4_q4
546  complex(dp) :: ghz1_dyn,ghz2_dyn,ghz3_dyn,ghz4_dyn
547 
548 
549 
550  res = 0d0
551  q1 = dcmplx(p(1,:),0d0)
552  q3 = dcmplx(p(3,:),0d0)
553  q4 = dcmplx(p(4,:),0d0)
554 
555 
556  e3 = sp(3,:)
557  e4 = sp(4,:)
558 
559  q = -q1
560 
561  q_q =sc(q,q)
562  q3_q3 = sc(q3,q3)
563  q4_q4 = sc(q4,q4)
564 
565  e3_e4 = sc(e3,e4)
566  e3_q4 = sc(e3,q4)
567  e4_q3 = sc(e4,q3)
568 
569 
570  if ((q_q).lt.-0.1d0 .or. (q3_q3).lt.-0.1d0 .or. (q4_q4).lt.-0.1d0) return ! if negative invariant masses return zero
571 
572 
573 !---- data that defines couplings
574  if( (vvmode.eq.zzmode) .or. (vvmode.eq.wwmode) ) then! decay ZZ's or WW's
575  ghz1_dyn = hvvspinzerodynamiccoupling(1,q3_q3,q4_q4,q_q)
576  ghz2_dyn = hvvspinzerodynamiccoupling(2,q3_q3,q4_q4,q_q)
577  ghz3_dyn = hvvspinzerodynamiccoupling(3,q3_q3,q4_q4,q_q)
578  ghz4_dyn = hvvspinzerodynamiccoupling(4,q3_q3,q4_q4,q_q)
579  elseif( (vvmode.eq.gszmode) ) then
580  ghz1_dyn = hvvspinzerodynamiccoupling(5,0d0,q3_q3,q_q)
581  ghz2_dyn = hvvspinzerodynamiccoupling(6,0d0,q3_q3,q_q)
582  ghz3_dyn = hvvspinzerodynamiccoupling(7,0d0,q3_q3,q_q)
583  ghz4_dyn = hvvspinzerodynamiccoupling(8,0d0,q3_q3,q_q)
584  elseif( (vvmode.eq.zgmode) .OR. (vvmode.eq.zgsmode) ) then
585  ghz1_dyn = hvvspinzerodynamiccoupling(5,0d0,q4_q4,q_q)
586  ghz2_dyn = hvvspinzerodynamiccoupling(6,0d0,q4_q4,q_q)
587  ghz3_dyn = hvvspinzerodynamiccoupling(7,0d0,q4_q4,q_q)
588  ghz4_dyn = hvvspinzerodynamiccoupling(8,0d0,q4_q4,q_q)
589  elseif( (vvmode.eq.ggmode) .or. (vvmode.eq.gsgsmode) .or. (vvmode.eq.gsgmode) ) then
590  ghz1_dyn = czero
591  ghz2_dyn = hvvspinzerodynamiccoupling(9,q3_q3,q4_q4,q_q)
592  ghz3_dyn = hvvspinzerodynamiccoupling(10,q3_q3,q4_q4,q_q)
593  ghz4_dyn = hvvspinzerodynamiccoupling(11,q3_q3,q4_q4,q_q)
594  elseif( (vvmode.eq.zzpmode) .or. (vvmode.eq.wwpmode) ) then
595  ghz1_dyn = hvvspinzerodynamiccoupling(12,q3_q3,q4_q4,q_q)
596  ghz2_dyn = hvvspinzerodynamiccoupling(13,q3_q3,q4_q4,q_q)
597  ghz3_dyn = hvvspinzerodynamiccoupling(14,q3_q3,q4_q4,q_q)
598  ghz4_dyn = hvvspinzerodynamiccoupling(15,q3_q3,q4_q4,q_q)
599  elseif( (vvmode.eq.zpzmode) .or. (vvmode.eq.wpwmode) ) then
600  ghz1_dyn = hvvspinzerodynamiccoupling(12,q4_q4,q3_q3,q_q)
601  ghz2_dyn = hvvspinzerodynamiccoupling(13,q4_q4,q3_q3,q_q)
602  ghz3_dyn = hvvspinzerodynamiccoupling(14,q4_q4,q3_q3,q_q)
603  ghz4_dyn = hvvspinzerodynamiccoupling(15,q4_q4,q3_q3,q_q)
604  elseif( (vvmode.eq.zpzpmode) .or. (vvmode.eq.wpwpmode) ) then
605  ghz1_dyn = hvvspinzerodynamiccoupling(16,q3_q3,q4_q4,q_q)
606  ghz2_dyn = hvvspinzerodynamiccoupling(17,q3_q3,q4_q4,q_q)
607  ghz3_dyn = hvvspinzerodynamiccoupling(18,q3_q3,q4_q4,q_q)
608  ghz4_dyn = hvvspinzerodynamiccoupling(19,q3_q3,q4_q4,q_q)
609  elseif( (vvmode.eq.gszpmode) ) then
610  ghz1_dyn = hvvspinzerodynamiccoupling(20,0d0,q3_q3,q_q)
611  ghz2_dyn = hvvspinzerodynamiccoupling(21,0d0,q3_q3,q_q)
612  ghz3_dyn = hvvspinzerodynamiccoupling(22,0d0,q3_q3,q_q)
613  ghz4_dyn = hvvspinzerodynamiccoupling(23,0d0,q3_q3,q_q)
614  elseif( (vvmode.eq.zpgmode) .OR. (vvmode.eq.zpgsmode) ) then
615  ghz1_dyn = hvvspinzerodynamiccoupling(20,0d0,q4_q4,q_q)
616  ghz2_dyn = hvvspinzerodynamiccoupling(21,0d0,q4_q4,q_q)
617  ghz3_dyn = hvvspinzerodynamiccoupling(22,0d0,q4_q4,q_q)
618  ghz4_dyn = hvvspinzerodynamiccoupling(23,0d0,q4_q4,q_q)
619  else
620  ghz1_dyn = czero
621  ghz2_dyn = czero
622  ghz3_dyn = czero
623  ghz4_dyn = czero
624  print *,"VVMode",vvmode,"not implemented"
625  endif
626 
627  if( .not. generate_as ) then
628  yyy1 = ghz1_dyn*m_v**2/q_q & ! in this line M_V is indeed correct, not a misprint
629  + ghz2_dyn*(q_q-q3_q3-q4_q4)/q_q &
630  + ghz3_dyn/lambda**2*(q_q-q3_q3-q4_q4)*(q_q-q4_q4-q3_q3)/4d0/q_q
hto_units::one
real *8, parameter one
Definition: CALLING_cpHTO.f:184
hto_units::ci
real *8, dimension(1:2) ci
Definition: CALLING_cpHTO.f:189
hto_units::zero
real *8, parameter zero
Definition: CALLING_cpHTO.f:185