← Back to branch summary

~maddevelopers/mg5amcnlo/3.0.2-alpha0

« back to all changes in this revision

Viewing changes to HELAS/upsxxx.F

Added Template and HELAS into bzr

Show diffs side-by-side

added added

removed removed

Lines of Context:
HELAS/upsxxx.F
 
1
      subroutine upsxxx(t2,sc,gt, xm,xw,t1)
 
2
 
 
3
c  Subroutines for graviton phase space integration
 
4
c  KEK 2009.11
 
5
c
 
6
      implicit none
 
7
      double complex t1(18), t2(18), sc(3), vertex, tc2(18)
 
8
      double complex gt(2)
 
9
      double precision xm,xw,xmass
 
10
 
 
11
      double complex ft(6,4),ft2(6,4)
 
12
      double precision ps1(4), pt2(4),PT(4),PTT2,PG(0:3)
 
13
      integer i
 
14
      double complex cZero
 
15
      double precision rZero, rTwo,Pi
 
16
      parameter( rZero = 0.0d0, rTwo = 2.0d0 )
 
17
      parameter( cZero = ( 0.0d0, 0.0d0 ) )
 
18
      double precision PADD
 
19
      double precision L_ADD,NADD,MGLOW,MGUP
 
20
      external txxxxx
 
21
 
 
22
      Pi=dacos(-1.0d0)
 
23
      L_ADD=dimag(gt(1))
 
24
      NADD=dreal(gt(1))
 
25
      MGUP=dimag(gt(2))
 
26
      MGLOW=dreal(gt(2))
 
27
 
 
28
 
 
29
      ps1(1) = dreal(sc(2))
 
30
      ps1(2) = dreal(sc(3))
 
31
      ps1(3) = dimag(sc(3))
 
32
      ps1(4) = dimag(sc(2))
 
33
 
 
34
      pt2(1) = -dreal(t2(17))
 
35
      pt2(2) = -dreal(t2(18))
 
36
      pt2(3) = -dimag(t2(18))
 
37
      pt2(4) = -dimag(t2(17))
 
38
 
 
39
      PG(0)=ps1(1)-pt2(1)
 
40
      PG(1)=ps1(2)-pt2(2)
 
41
      PG(2)=ps1(3)-pt2(3)
 
42
      PG(3)=ps1(4)-pt2(4)
 
43
 
 
44
      PTT2=PG(0)**2-PG(1)**2-PG(2)**2-PG(3)**2
 
45
      xmass=dsqrt(PTT2)
 
46
 
 
47
       t1(17) = dcmplx(PG(0),PG(3))
 
48
       t1(18) = dcmplx(PG(1), PG(2))
 
49
 
 
50
      if(xmass.lt.MGLOW.or.xmass.gt.MGUP) then
 
51
      do i=1,16
 
52
      t1(i)=dcmplx(0.0d0,0.0d0)
 
53
      enddo
 
54
      return
 
55
      endif
 
56
 
 
57
      CALL txxxxx(PG,xmass,INT(t2(1)),+1 , t1)
 
58
 
 
59
 
 
60
       if(INT(NADD).eq.2) then
 
61
         PADD=2.0d0*Pi
 
62
        elseif(INT(NADD).eq.3) then
 
63
         PADD=4.0d0*Pi
 
64
        elseif(INT(NADD).eq.4) then
 
65
         PADD=2.0d0*Pi**2
 
66
        elseif(INT(NADD).eq.5) then
 
67
          PADD=8.0d0/3.0d0*Pi**2
 
68
        elseif(INT(NADD).eq.6) then
 
69
           PADD=Pi**3
 
70
        else
 
71
        print *, "OUT CASE"
 
72
        stop
 
73
        endif 
 
74
 
 
75
 
 
76
        do i=1,16
 
77
        t1(i)=-1.0d0*t1(i)*
 
78
     & dsqrt( 
 
79
     & 2.0d0*Pi*8.0d0*Pi  ! to compensate the decay phase factor
 
80
     &* PADD/L_ADD**NADD*xmass**(NADD-1)  ! density factor for d=4 case
 
81
     &/2.0d0/xmass)   ! dm=dm^2/2/m    
 
82
        enddo
 
83
 
 
84
      return
 
85
      end