~maddevelopers/mg5amcnlo/3.0.2-alpha0

« back to all changes in this revision

Viewing changes to HELAS/txxxxx.f

Added Template and HELAS into bzr

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
      subroutine txxxxx(p,tmass,nhel,nst , tc)
 
2
c
 
3
c This subroutine computes a TENSOR wavefunction.
 
4
c
 
5
c input:
 
6
c       real    p(0:3)         : four-momentum of tensor boson
 
7
c       real    tmass          : mass          of tensor boson
 
8
c       integer nhel           : helicity      of tensor boson
 
9
c                = -2,-1,0,1,2 : (0 is forbidden if tmass=0.0)
 
10
c       integer nst  = -1 or 1 : +1 for final, -1 for initial
 
11
c
 
12
c output:
 
13
c       complex tc(18)         : tensor wavefunction    epsilon^mu^nu(t)
 
14
c     
 
15
      implicit none
 
16
      double precision p(0:3), tmass
 
17
      integer nhel, nst
 
18
      double complex tc(18)
 
19
 
 
20
      double complex ft(6,4), ep(4), em(4), e0(4)
 
21
      double precision pt, pt2, pp, pzpt, emp, sqh, sqs
 
22
      integer i, j
 
23
 
 
24
      double precision rZero, rHalf, rOne, rTwo
 
25
      parameter( rZero = 0.0d0, rHalf = 0.5d0 )
 
26
      parameter( rOne = 1.0d0, rTwo = 2.0d0 )
 
27
 
 
28
      integer stdo
 
29
      parameter( stdo = 6 )
 
30
 
 
31
      sqh = sqrt(rHalf)
 
32
      sqs = sqrt(rHalf/3.d0)
 
33
 
 
34
      pt2 = p(1)**2 + p(2)**2
 
35
      pp = min(p(0),sqrt(pt2+p(3)**2))
 
36
      pt = min(pp,sqrt(pt2))
 
37
 
 
38
      ft(5,1) = dcmplx(p(0),p(3))*nst
 
39
      ft(6,1) = dcmplx(p(1),p(2))*nst
 
40
 
 
41
      if ( nhel.ge.0 ) then 
 
42
c construct eps+
 
43
         if ( pp.eq.rZero ) then
 
44
            ep(1) = dcmplx( rZero )
 
45
            ep(2) = dcmplx( -sqh )
 
46
            ep(3) = dcmplx( rZero , nst*sqh )
 
47
            ep(4) = dcmplx( rZero )
 
48
         else
 
49
            ep(1) = dcmplx( rZero )
 
50
            ep(4) = dcmplx( pt/pp*sqh )
 
51
            if ( pt.ne.rZero ) then
 
52
               pzpt = p(3)/(pp*pt)*sqh
 
53
               ep(2) = dcmplx( -p(1)*pzpt , -nst*p(2)/pt*sqh )
 
54
               ep(3) = dcmplx( -p(2)*pzpt ,  nst*p(1)/pt*sqh )
 
55
            else
 
56
               ep(2) = dcmplx( -sqh )
 
57
               ep(3) = dcmplx( rZero , nst*sign(sqh,p(3)) )
 
58
            endif
 
59
         endif
 
60
      end if
 
61
 
 
62
      if ( nhel.le.0 ) then 
 
63
c construct eps-
 
64
         if ( pp.eq.rZero ) then
 
65
            em(1) = dcmplx( rZero )
 
66
            em(2) = dcmplx( sqh )
 
67
            em(3) = dcmplx( rZero , nst*sqh )
 
68
            em(4) = dcmplx( rZero )
 
69
         else
 
70
            em(1) = dcmplx( rZero )
 
71
            em(4) = dcmplx( -pt/pp*sqh )
 
72
            if ( pt.ne.rZero ) then
 
73
               pzpt = -p(3)/(pp*pt)*sqh
 
74
               em(2) = dcmplx( -p(1)*pzpt , -nst*p(2)/pt*sqh )
 
75
               em(3) = dcmplx( -p(2)*pzpt ,  nst*p(1)/pt*sqh )
 
76
            else
 
77
               em(2) = dcmplx( sqh )
 
78
               em(3) = dcmplx( rZero , nst*sign(sqh,p(3)) )
 
79
            endif
 
80
         endif
 
81
      end if
 
82
 
 
83
      if ( abs(nhel).le.1 ) then  
 
84
c construct eps0
 
85
         if ( pp.eq.rZero ) then
 
86
            e0(1) = dcmplx( rZero )
 
87
            e0(2) = dcmplx( rZero )
 
88
            e0(3) = dcmplx( rZero )
 
89
            e0(4) = dcmplx( rOne )
 
90
         else
 
91
            emp = p(0)/(tmass*pp)
 
92
            e0(1) = dcmplx( pp/tmass )
 
93
            e0(4) = dcmplx( p(3)*emp )
 
94
            if ( pt.ne.rZero ) then
 
95
               e0(2) = dcmplx( p(1)*emp )
 
96
               e0(3) = dcmplx( p(2)*emp )
 
97
            else
 
98
               e0(2) = dcmplx( rZero )
 
99
               e0(3) = dcmplx( rZero )
 
100
            endif
 
101
         end if
 
102
      end if
 
103
 
 
104
      if ( nhel.eq.2 ) then
 
105
         do j = 1,4
 
106
            do i = 1,4
 
107
               ft(i,j) = ep(i)*ep(j)
 
108
            end do
 
109
         end do
 
110
      else if ( nhel.eq.-2 ) then
 
111
         do j = 1,4
 
112
            do i = 1,4
 
113
               ft(i,j) = em(i)*em(j)
 
114
            end do
 
115
         end do
 
116
      else if (tmass.eq.0) then
 
117
         do j = 1,4
 
118
            do i = 1,4
 
119
               ft(i,j) = 0
 
120
            end do
 
121
         end do
 
122
      else if (tmass.ne.0) then
 
123
        if  ( nhel.eq.1 ) then
 
124
           do j = 1,4
 
125
              do i = 1,4
 
126
                 ft(i,j) = sqh*( ep(i)*e0(j) + e0(i)*ep(j) )
 
127
              end do
 
128
           end do
 
129
        else if ( nhel.eq.0 ) then
 
130
           do j = 1,4
 
131
              do i = 1,4
 
132
                 ft(i,j) = sqs*( ep(i)*em(j) + em(i)*ep(j)
 
133
     &                                + rTwo*e0(i)*e0(j) )
 
134
              end do
 
135
           end do
 
136
        else if ( nhel.eq.-1 ) then
 
137
           do j = 1,4
 
138
              do i = 1,4
 
139
                 ft(i,j) = sqh*( em(i)*e0(j) + e0(i)*em(j) )
 
140
              end do
 
141
           end do
 
142
        else
 
143
           write(stdo,*) 'invalid helicity in TXXXXX'
 
144
           stop
 
145
        end if
 
146
      end if
 
147
 
 
148
      tc(1) = ft(1,1)
 
149
      tc(2) = ft(1,2)
 
150
      tc(3) = ft(1,3)
 
151
      tc(4) = ft(1,4)
 
152
      tc(5) = ft(2,1)
 
153
      tc(6) = ft(2,2)
 
154
      tc(7) = ft(2,3)
 
155
      tc(8) = ft(2,4)
 
156
      tc(9) = ft(3,1)
 
157
      tc(10) = ft(3,2)
 
158
      tc(11) = ft(3,3)
 
159
      tc(12) = ft(3,4)
 
160
      tc(13) = ft(4,1)
 
161
      tc(14) = ft(4,2)
 
162
      tc(15) = ft(4,3)
 
163
      tc(16) = ft(4,4)
 
164
      tc(17) = ft(5,1)
 
165
      tc(18) = ft(6,1)
 
166
 
 
167
      return
 
168
      end