C++InterfacetoTauola
VBF_functions.f
1 C###############################################################################
2 C
3 C Copyright (c) 2010 The ALOHA Development team and Contributors
4 C
5 C This file is a part of the MadGraph 5 project, an application which
6 C automatically generates Feynman diagrams and matrix elements for arbitrary
7 C high-energy processes in the Standard Model and beyond.
8 C
9 C It is subject to the ALOHA license which should accompany this
10 C distribution.
11 C
12 C###############################################################################
13  subroutine ixxxxx(p, fmass, nhel, nsf ,fi)
14 c
15 c This subroutine computes a fermion wavefunction with the flowing-IN
16 c fermion number.
17 c
18 c input:
19 c real p(0:3) : four-momentum of fermion
20 c real fmass : mass of fermion
21 c integer nhel = -1 or 1 : helicity of fermion
22 c integer nsf = -1 or 1 : +1 for particle, -1 for anti-particle
23 c
24 c output:
25 c complex fi(6) : fermion wavefunction |fi>
26 c
27  implicit none
28  double complex fi(6),chi(2)
29  double precision p(0:3),sf(2),sfomeg(2),omega(2),fmass,
30  & pp,pp3,sqp0p3,sqm(0:1)
31  integer nhel,nsf,ip,im,nh
32 
33  double precision rzero, rhalf, rtwo
34  parameter( rzero = 0.0d0, rhalf = 0.5d0, rtwo = 2.0d0 )
35 
36 c#ifdef HELAS_CHECK
37 c double precision p2
38 c double precision epsi
39 c parameter( epsi = 2.0d-5 )
40 c integer stdo
41 c parameter( stdo = 6 )
42 c#endif
43 c
44 c#ifdef HELAS_CHECK
45 c pp = sqrt(p(1)**2+p(2)**2+p(3)**2)
46 c if ( abs(p(0))+pp.eq.rZero ) then
47 c write(stdo,*)
48 c & ' helas-error : p(0:3) in ixxxxx is zero momentum'
49 c endif
50 c if ( p(0).le.rZero ) then
51 c write(stdo,*)
52 c & ' helas-error : p(0:3) in ixxxxx has non-positive energy'
53 c write(stdo,*)
54 c & ' : p(0) = ',p(0)
55 c endif
56 c p2 = (p(0)-pp)*(p(0)+pp)
57 c if ( abs(p2-fmass**2).gt.p(0)**2*epsi ) then
58 c write(stdo,*)
59 c & ' helas-error : p(0:3) in ixxxxx has inappropriate mass'
60 c write(stdo,*)
61 c & ' : p**2 = ',p2,' : fmass**2 = ',fmass**2
62 c endif
63 c if (abs(nhel).ne.1) then
64 c write(stdo,*) ' helas-error : nhel in ixxxxx is not -1,1'
65 c write(stdo,*) ' : nhel = ',nhel
66 c endif
67 c if (abs(nsf).ne.1) then
68 c write(stdo,*) ' helas-error : nsf in ixxxxx is not -1,1'
69 c write(stdo,*) ' : nsf = ',nsf
70 c endif
71 c#endif
72 
73  fi(1) = dcmplx(p(0),p(3))*nsf*(-1)
74  fi(2) = dcmplx(p(1),p(2))*nsf*(-1)
75 
76  nh = nhel*nsf
77 
78  if ( fmass.ne.rzero ) then
79 
80  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
81 
82  if ( pp.eq.rzero ) then
83 
84  sqm(0) = dsqrt(abs(fmass)) ! possibility of negative fermion masses
85  sqm(1) = sign(sqm(0),fmass) ! possibility of negative fermion masses
86  ip = (1+nh)/2
87  im = (1-nh)/2
88 
89  fi(3) = ip * sqm(ip)
90  fi(4) = im*nsf * sqm(ip)
91  fi(5) = ip*nsf * sqm(im)
92  fi(6) = im * sqm(im)
93 
94  else
95 
96  sf(1) = dble(1+nsf+(1-nsf)*nh)*rhalf
97  sf(2) = dble(1+nsf-(1-nsf)*nh)*rhalf
98  omega(1) = dsqrt(p(0)+pp)
99  omega(2) = fmass/omega(1)
100  ip = (3+nh)/2
101  im = (3-nh)/2
102  sfomeg(1) = sf(1)*omega(ip)
103  sfomeg(2) = sf(2)*omega(im)
104  pp3 = max(pp+p(3),rzero)
105  chi(1) = dcmplx( dsqrt(pp3*rhalf/pp) )
106  if ( pp3.eq.rzero ) then
107  chi(2) = dcmplx(-nh )
108  else
109  chi(2) = dcmplx( nh*p(1) , p(2) )/dsqrt(rtwo*pp*pp3)
110  endif
111 
112  fi(3) = sfomeg(1)*chi(im)
113  fi(4) = sfomeg(1)*chi(ip)
114  fi(5) = sfomeg(2)*chi(im)
115  fi(6) = sfomeg(2)*chi(ip)
116 
117  endif
118 
119  else
120 
121  if(p(1).eq.0d0.and.p(2).eq.0d0.and.p(3).lt.0d0) then
122  sqp0p3 = 0d0
123  else
124  sqp0p3 = dsqrt(max(p(0)+p(3),rzero))*nsf
125  end if
126  chi(1) = dcmplx( sqp0p3 )
127  if ( sqp0p3.eq.rzero ) then
128  chi(2) = dcmplx(-nhel )*dsqrt(rtwo*p(0))
129  else
130  chi(2) = dcmplx( nh*p(1), p(2) )/sqp0p3
131  endif
132  if ( nh.eq.1 ) then
133  fi(3) = dcmplx( rzero )
134  fi(4) = dcmplx( rzero )
135  fi(5) = chi(1)
136  fi(6) = chi(2)
137  else
138  fi(3) = chi(2)
139  fi(4) = chi(1)
140  fi(5) = dcmplx( rzero )
141  fi(6) = dcmplx( rzero )
142  endif
143  endif
144 c
145  return
146  end
147 
148 
149  subroutine oxxxxx(p,fmass,nhel,nsf , fo)
150 c
151 c This subroutine computes a fermion wavefunction with the flowing-OUT
152 c fermion number.
153 c
154 c input:
155 c real p(0:3) : four-momentum of fermion
156 c real fmass : mass of fermion
157 c integer nhel = -1 or 1 : helicity of fermion
158 c integer nsf = -1 or 1 : +1 for particle, -1 for anti-particle
159 c
160 c output:
161 c complex fo(6) : fermion wavefunction <fo|
162 c
163  implicit none
164  double complex fo(6),chi(2)
165  double precision p(0:3),sf(2),sfomeg(2),omega(2),fmass,
166  & pp,pp3,sqp0p3,sqm(0:1)
167  integer nhel,nsf,nh,ip,im
168 
169  double precision rzero, rhalf, rtwo
170  parameter( rzero = 0.0d0, rhalf = 0.5d0, rtwo = 2.0d0 )
171 
172 c#ifdef HELAS_CHECK
173 c double precision p2
174 c double precision epsi
175 c parameter( epsi = 2.0d-5 )
176 c integer stdo
177 c parameter( stdo = 6 )
178 c#endif
179 c
180 c#ifdef HELAS_CHECK
181 c pp = sqrt(p(1)**2+p(2)**2+p(3)**2)
182 c if ( abs(p(0))+pp.eq.rZero ) then
183 c write(stdo,*)
184 c & ' helas-error : p(0:3) in oxxxxx is zero momentum'
185 c endif
186 c if ( p(0).le.rZero ) then
187 c write(stdo,*)
188 c & ' helas-error : p(0:3) in oxxxxx has non-positive energy'
189 c write(stdo,*)
190 c & ' : p(0) = ',p(0)
191 c endif
192 c p2 = (p(0)-pp)*(p(0)+pp)
193 c if ( abs(p2-fmass**2).gt.p(0)**2*epsi ) then
194 c write(stdo,*)
195 c & ' helas-error : p(0:3) in oxxxxx has inappropriate mass'
196 c write(stdo,*)
197 c & ' : p**2 = ',p2,' : fmass**2 = ',fmass**2
198 c endif
199 c if ( abs(nhel).ne.1 ) then
200 c write(stdo,*) ' helas-error : nhel in oxxxxx is not -1,1'
201 c write(stdo,*) ' : nhel = ',nhel
202 c endif
203 c if ( abs(nsf).ne.1 ) then
204 c write(stdo,*) ' helas-error : nsf in oxxxxx is not -1,1'
205 c write(stdo,*) ' : nsf = ',nsf
206 c endif
207 c#endif
208 
209  fo(1) = dcmplx(p(0),p(3))*nsf
210  fo(2) = dcmplx(p(1),p(2))*nsf
211 
212  nh = nhel*nsf
213 
214  if ( fmass.ne.rzero ) then
215 
216  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
217 
218  if ( pp.eq.rzero ) then
219 
220  sqm(0) = dsqrt(abs(fmass)) ! possibility of negative fermion masses
221  sqm(1) = sign(sqm(0),fmass) ! possibility of negative fermion masses
222  im = nhel * (1+nh)/2
223  ip = nhel * (-1) * ((1-nh)/2)
224  fo(3) = im * sqm(abs(ip))
225  fo(4) = ip*nsf * sqm(abs(ip))
226  fo(5) = im*nsf * sqm(abs(im))
227  fo(6) = ip * sqm(abs(im))
228  else
229 
230  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
231  sf(1) = dble(1+nsf+(1-nsf)*nh)*rhalf
232  sf(2) = dble(1+nsf-(1-nsf)*nh)*rhalf
233  omega(1) = dsqrt(p(0)+pp)
234  omega(2) = fmass/omega(1)
235  ip = (3+nh)/2
236  im = (3-nh)/2
237  sfomeg(1) = sf(1)*omega(ip)
238  sfomeg(2) = sf(2)*omega(im)
239  pp3 = max(pp+p(3),rzero)
240  chi(1) = dcmplx( dsqrt(pp3*rhalf/pp) )
241  if ( pp3.eq.rzero ) then
242  chi(2) = dcmplx(-nh )
243  else
244  chi(2) = dcmplx( nh*p(1) , -p(2) )/dsqrt(rtwo*pp*pp3)
245  endif
246 
247  fo(3) = sfomeg(2)*chi(im)
248  fo(4) = sfomeg(2)*chi(ip)
249  fo(5) = sfomeg(1)*chi(im)
250  fo(6) = sfomeg(1)*chi(ip)
251 
252  endif
253 
254  else
255 
256  if(p(1).eq.0d0.and.p(2).eq.0d0.and.p(3).lt.0d0) then
257  sqp0p3 = 0d0
258  else
259  sqp0p3 = dsqrt(max(p(0)+p(3),rzero))*nsf
260  end if
261  chi(1) = dcmplx( sqp0p3 )
262  if ( sqp0p3.eq.rzero ) then
263  chi(2) = dcmplx(-nhel )*dsqrt(rtwo*p(0))
264  else
265  chi(2) = dcmplx( nh*p(1), -p(2) )/sqp0p3
266  endif
267  if ( nh.eq.1 ) then
268  fo(3) = chi(1)
269  fo(4) = chi(2)
270  fo(5) = dcmplx( rzero )
271  fo(6) = dcmplx( rzero )
272  else
273  fo(3) = dcmplx( rzero )
274  fo(4) = dcmplx( rzero )
275  fo(5) = chi(2)
276  fo(6) = chi(1)
277  endif
278 
279  endif
280 c
281  return
282  end
283 
284  subroutine pxxxxx(p,tmass,nhel,nst , tc)
285 
286 c CP3 2009.NOV
287 
288 c This subroutine computes a PSEUDOR wavefunction.
289 c
290 c input:
291 c real p(0:3) : four-momentum of tensor boson
292 c real tmass : mass of tensor boson
293 c integer nhel : helicity of tensor boson
294 c = -2,-1,0,1,2 : (0 is forbidden if tmass=0.0)
295 c integer nst = -1 or 1 : +1 for final, -1 for initial
296 c
297 c output:
298 c complex tc(18) : PSEUDOR wavefunction epsilon^mu^nu(t)
299 c
300  implicit none
301  double precision p(0:3), tmass
302  integer nhel, nst
303  double complex tc(18)
304 
305  double complex ft(6,4), ep(4), em(4), e0(4)
306  double precision pt, pt2, pp, pzpt, emp, sqh, sqs
307  integer i, j
308 
309  double precision rzero, rhalf, rone, rtwo
310  parameter( rzero = 0.0d0, rhalf = 0.5d0 )
311  parameter( rone = 1.0d0, rtwo = 2.0d0 )
312 
313 
314  tc(3)=nhel
315  tc(1) = dcmplx(p(0),p(3))*nst
316  tc(2) = dcmplx(p(1),p(2))*nst
317 
318  return
319  end
320 
321  subroutine sxxxxx(p,nss , sc)
322 c
323 c This subroutine computes a complex SCALAR wavefunction.
324 c
325 c input:
326 c real p(0:3) : four-momentum of scalar boson
327 c integer nss = -1 or 1 : +1 for final, -1 for initial
328 c
329 c output:
330 c complex sc(3) : scalar wavefunction s
331 c
332  implicit none
333  double complex sc(3)
334  double precision p(0:3)
335  integer nss
336 
337  double precision rone
338  parameter( rone = 1.0d0 )
339 
340 c#ifdef HELAS_CHECK
341 c double precision p2
342 c double precision epsi
343 c parameter( epsi = 2.0d-5 )
344 c double precision rZero
345 c parameter( rZero = 0.0d0 )
346 c integer stdo
347 c parameter( stdo = 6 )
348 c#endif
349 c
350 c#ifdef HELAS_CHECK
351 c if ( abs(p(0))+abs(p(1))+abs(p(2))+abs(p(3)).eq.rZero ) then
352 c write(stdo,*)
353 c & ' helas-error : p(0:3) in sxxxxx is zero momentum'
354 c endif
355 c if ( p(0).le.rZero ) then
356 c write(stdo,*)
357 c & ' helas-error : p(0:3) in sxxxxx has non-positive energy'
358 c write(stdo,*)
359 c & ' : p(0) = ',p(0)
360 c endif
361 c p2 = p(0)**2-(p(1)**2+p(2)**2+p(3)**2)
362 c if ( p2.lt.-p(0)**2*epsi ) then
363 c write(stdo,*) ' helas-error : p(0:3) in sxxxxx is spacelike'
364 c write(stdo,*) ' : p**2 = ',p2
365 c endif
366 c if ( abs(nss).ne.1 ) then
367 c write(stdo,*) ' helas-error : nss in sxxxxx is not -1,1'
368 c write(stdo,*) ' : nss = ',nss
369 c endif
370 c#endif
371 
372  sc(3) = dcmplx( rone )
373  sc(1) = dcmplx(p(0),p(3))*nss
374  sc(2) = dcmplx(p(1),p(2))*nss
375 c
376  return
377  end
378 
379  subroutine txxxxx(p,tmass,nhel,nst , tc)
380 c
381 c This subroutine computes a TENSOR wavefunction.
382 c
383 c input:
384 c real p(0:3) : four-momentum of tensor boson
385 c real tmass : mass of tensor boson
386 c integer nhel : helicity of tensor boson
387 c = -2,-1,0,1,2 : (0 is forbidden if tmass=0.0)
388 c integer nst = -1 or 1 : +1 for final, -1 for initial
389 c
390 c output:
391 c complex tc(18) : tensor wavefunction epsilon^mu^nu(t)
392 c
393  implicit none
394  double precision p(0:3), tmass
395  integer nhel, nst
396  double complex tc(18)
397 
398  double complex ft(6,4), ep(4), em(4), e0(4)
399  double precision pt, pt2, pp, pzpt, emp, sqh, sqs
400  integer i, j
401 
402  double precision rzero, rhalf, rone, rtwo
403  parameter( rzero = 0.0d0, rhalf = 0.5d0 )
404  parameter( rone = 1.0d0, rtwo = 2.0d0 )
405 
406  integer stdo
407  parameter( stdo = 6 )
408 
409  sqh = sqrt(rhalf)
410  sqs = sqrt(rhalf/3.d0)
411 
412  pt2 = p(1)**2 + p(2)**2
413  pp = min(p(0),sqrt(pt2+p(3)**2))
414  pt = min(pp,sqrt(pt2))
415 
416  ft(5,1) = dcmplx(p(0),p(3))*nst
417  ft(6,1) = dcmplx(p(1),p(2))*nst
418 
419  if ( nhel.ge.0 ) then
420 c construct eps+
421  if ( pp.eq.rzero ) then
422  ep(1) = dcmplx( rzero )
423  ep(2) = dcmplx( -sqh )
424  ep(3) = dcmplx( rzero , nst*sqh )
425  ep(4) = dcmplx( rzero )
426  else
427  ep(1) = dcmplx( rzero )
428  ep(4) = dcmplx( pt/pp*sqh )
429  if ( pt.ne.rzero ) then
430  pzpt = p(3)/(pp*pt)*sqh
431  ep(2) = dcmplx( -p(1)*pzpt , -nst*p(2)/pt*sqh )
432  ep(3) = dcmplx( -p(2)*pzpt , nst*p(1)/pt*sqh )
433  else
434  ep(2) = dcmplx( -sqh )
435  ep(3) = dcmplx( rzero , nst*sign(sqh,p(3)) )
436  endif
437  endif
438  end if
439 
440  if ( nhel.le.0 ) then
441 c construct eps-
442  if ( pp.eq.rzero ) then
443  em(1) = dcmplx( rzero )
444  em(2) = dcmplx( sqh )
445  em(3) = dcmplx( rzero , nst*sqh )
446  em(4) = dcmplx( rzero )
447  else
448  em(1) = dcmplx( rzero )
449  em(4) = dcmplx( -pt/pp*sqh )
450  if ( pt.ne.rzero ) then
451  pzpt = -p(3)/(pp*pt)*sqh
452  em(2) = dcmplx( -p(1)*pzpt , -nst*p(2)/pt*sqh )
453  em(3) = dcmplx( -p(2)*pzpt , nst*p(1)/pt*sqh )
454  else
455  em(2) = dcmplx( sqh )
456  em(3) = dcmplx( rzero , nst*sign(sqh,p(3)) )
457  endif
458  endif
459  end if
460 
461  if ( abs(nhel).le.1 ) then
462 c construct eps0
463  if ( pp.eq.rzero ) then
464  e0(1) = dcmplx( rzero )
465  e0(2) = dcmplx( rzero )
466  e0(3) = dcmplx( rzero )
467  e0(4) = dcmplx( rone )
468  else
469  emp = p(0)/(tmass*pp)
470  e0(1) = dcmplx( pp/tmass )
471  e0(4) = dcmplx( p(3)*emp )
472  if ( pt.ne.rzero ) then
473  e0(2) = dcmplx( p(1)*emp )
474  e0(3) = dcmplx( p(2)*emp )
475  else
476  e0(2) = dcmplx( rzero )
477  e0(3) = dcmplx( rzero )
478  endif
479  end if
480  end if
481 
482  if ( nhel.eq.2 ) then
483  do j = 1,4
484  do i = 1,4
485  ft(i,j) = ep(i)*ep(j)
486  end do
487  end do
488  else if ( nhel.eq.-2 ) then
489  do j = 1,4
490  do i = 1,4
491  ft(i,j) = em(i)*em(j)
492  end do
493  end do
494  else if (tmass.eq.0) then
495  do j = 1,4
496  do i = 1,4
497  ft(i,j) = 0
498  end do
499  end do
500  else if (tmass.ne.0) then
501  if ( nhel.eq.1 ) then
502  do j = 1,4
503  do i = 1,4
504  ft(i,j) = sqh*( ep(i)*e0(j) + e0(i)*ep(j) )
505  end do
506  end do
507  else if ( nhel.eq.0 ) then
508  do j = 1,4
509  do i = 1,4
510  ft(i,j) = sqs*( ep(i)*em(j) + em(i)*ep(j)
511  & + rtwo*e0(i)*e0(j) )
512  end do
513  end do
514  else if ( nhel.eq.-1 ) then
515  do j = 1,4
516  do i = 1,4
517  ft(i,j) = sqh*( em(i)*e0(j) + e0(i)*em(j) )
518  end do
519  end do
520  else
521  write(stdo,*) 'invalid helicity in TXXXXX'
522  stop
523  end if
524  end if
525 
526  tc(3) = ft(1,1)
527  tc(4) = ft(1,2)
528  tc(5) = ft(1,3)
529  tc(6) = ft(1,4)
530  tc(7) = ft(2,1)
531  tc(8) = ft(2,2)
532  tc(9) = ft(2,3)
533  tc(10) = ft(2,4)
534  tc(11) = ft(3,1)
535  tc(12) = ft(3,2)
536  tc(13) = ft(3,3)
537  tc(14) = ft(3,4)
538  tc(15) = ft(4,1)
539  tc(16) = ft(4,2)
540  tc(17) = ft(4,3)
541  tc(18) = ft(4,4)
542  tc(1) = ft(5,1)
543  tc(2) = ft(6,1)
544 
545  return
546  end
547 
548 
549  subroutine vxxxxx(p,vmass,nhel,nsv , vc)
550 c
551 c This subroutine computes a VECTOR wavefunction.
552 c
553 c input:
554 c real p(0:3) : four-momentum of vector boson
555 c real vmass : mass of vector boson
556 c integer nhel = -1, 0, 1: helicity of vector boson
557 c (0 is forbidden if vmass=0.0)
558 c integer nsv = -1 or 1 : +1 for final, -1 for initial
559 c
560 c output:
561 c complex vc(6) : vector wavefunction epsilon^mu(v)
562 c
563  implicit none
564  double complex vc(6)
565  double precision p(0:3),vmass,hel,hel0,pt,pt2,pp,pzpt,emp,sqh
566  integer nhel,nsv,nsvahl
567 
568  double precision rzero, rhalf, rone, rtwo
569  parameter( rzero = 0.0d0, rhalf = 0.5d0 )
570  parameter( rone = 1.0d0, rtwo = 2.0d0 )
571 
572 c#ifdef HELAS_CHECK
573 c double precision p2
574 c double precision epsi
575 c parameter( epsi = 2.0d-5 )
576 c integer stdo
577 c parameter( stdo = 6 )
578 c#endif
579 c
580 c#ifdef HELAS_CHECK
581 c pp = sqrt(p(1)**2+p(2)**2+p(3)**2)
582 c if ( abs(p(0))+pp.eq.rZero ) then
583 c write(stdo,*)
584 c & ' helas-error : p(0:3) in vxxxxx is zero momentum'
585 c endif
586 c if ( p(0).le.rZero ) then
587 c write(stdo,*)
588 c & ' helas-error : p(0:3) in vxxxxx has non-positive energy'
589 c write(stdo,*)
590 c & ' : p(0) = ',p(0)
591 c endif
592 c p2 = (p(0)+pp)*(p(0)-pp)
593 c if ( abs(p2-vmass**2).gt.p(0)**2*2.e-5 ) then
594 c write(stdo,*)
595 c & ' helas-error : p(0:3) in vxxxxx has inappropriate mass'
596 c write(stdo,*)
597 c & ' : p**2 = ',p2,' : vmass**2 = ',vmass**2
598 c endif
599 c if ( vmass.ne.rZero ) then
600 c if ( abs(nhel).gt.1 ) then
601 c write(stdo,*) ' helas-error : nhel in vxxxxx is not -1,0,1'
602 c write(stdo,*) ' : nhel = ',nhel
603 c endif
604 c else
605 c if ( abs(nhel).ne.1 ) then
606 c write(stdo,*) ' helas-error : nhel in vxxxxx is not -1,1'
607 c write(stdo,*) ' : nhel = ',nhel
608 c endif
609 c endif
610 c if ( abs(nsv).ne.1 ) then
611 c write(stdo,*) ' helas-error : nsv in vmxxxx is not -1,1'
612 c write(stdo,*) ' : nsv = ',nsv
613 c endif
614 c#endif
615 
616  sqh = dsqrt(rhalf)
617  hel = dble(nhel)
618  nsvahl = nsv*dabs(hel)
619  pt2 = p(1)**2+p(2)**2
620  pp = min(p(0),dsqrt(pt2+p(3)**2))
621  pt = min(pp,dsqrt(pt2))
622 
623  vc(1) = dcmplx(p(0),p(3))*nsv
624  vc(2) = dcmplx(p(1),p(2))*nsv
625 
626 c#ifdef HELAS_CHECK
627 c nhel=4 option for scalar polarization
628 c if( nhel.eq.4 ) then
629 c if( vmass.eq.rZero ) then
630 c vc(1) = rOne
631 c vc(2) = p(1)/p(0)
632 c vc(3) = p(2)/p(0)
633 c vc(4) = p(3)/p(0)
634 c else
635 c vc(1) = p(0)/vmass
636 c vc(2) = p(1)/vmass
637 c vc(3) = p(2)/vmass
638 c vc(4) = p(3)/vmass
639 c endif
640 c return
641 c endif
642 c#endif
643 
644  if ( vmass.ne.rzero ) then
645 
646  hel0 = rone-dabs(hel)
647 
648  if ( pp.eq.rzero ) then
649 
650  vc(3) = dcmplx( rzero )
651  vc(4) = dcmplx(-hel*sqh )
652  vc(5) = dcmplx( rzero , nsvahl*sqh )
653  vc(6) = dcmplx( hel0 )
654 
655  else
656 
657  emp = p(0)/(vmass*pp)
658  vc(3) = dcmplx( hel0*pp/vmass )
659  vc(6) = dcmplx( hel0*p(3)*emp+hel*pt/pp*sqh )
660  if ( pt.ne.rzero ) then
661  pzpt = p(3)/(pp*pt)*sqh*hel
662  vc(4) = dcmplx( hel0*p(1)*emp-p(1)*pzpt ,
663  & -nsvahl*p(2)/pt*sqh )
664  vc(5) = dcmplx( hel0*p(2)*emp-p(2)*pzpt ,
665  & nsvahl*p(1)/pt*sqh )
666  else
667  vc(4) = dcmplx( -hel*sqh )
668  vc(5) = dcmplx( rzero , nsvahl*sign(sqh,p(3)) )
669  endif
670 
671  endif
672 
673  else
674 
675  pp = p(0)
676  pt = sqrt(p(1)**2+p(2)**2)
677  vc(3) = dcmplx( rzero )
678  vc(6) = dcmplx( hel*pt/pp*sqh )
679  if ( pt.ne.rzero ) then
680  pzpt = p(3)/(pp*pt)*sqh*hel
681  vc(4) = dcmplx( -p(1)*pzpt , -nsv*p(2)/pt*sqh )
682  vc(5) = dcmplx( -p(2)*pzpt , nsv*p(1)/pt*sqh )
683  else
684  vc(4) = dcmplx( -hel*sqh )
685  vc(5) = dcmplx( rzero , nsv*sign(sqh,p(3)) )
686  endif
687 
688  endif
689 c
690  return
691  end
692 
693  subroutine boostx(p,q , pboost)
694 c
695 c This subroutine performs the Lorentz boost of a four-momentum. The
696 c momentum p is assumed to be given in the rest frame of q. pboost is
697 c the momentum p boosted to the frame in which q is given. q must be a
698 c timelike momentum.
699 c
700 c input:
701 c real p(0:3) : four-momentum p in the q rest frame
702 c real q(0:3) : four-momentum q in the boosted frame
703 c
704 c output:
705 c real pboost(0:3) : four-momentum p in the boosted frame
706 c
707  implicit none
708  double precision p(0:3),q(0:3),pboost(0:3),pq,qq,m,lf
709 
710  double precision rzero
711  parameter( rzero = 0.0d0 )
712 
713 c#ifdef HELAS_CHECK
714 c integer stdo
715 c parameter( stdo = 6 )
716 c double precision pp
717 c#endif
718 c
719  qq = q(1)**2+q(2)**2+q(3)**2
720 
721 c#ifdef HELAS_CHECK
722 c if (abs(p(0))+abs(p(1))+abs(p(2))+abs(p(3)).eq.rZero) then
723 c write(stdo,*)
724 c & ' helas-error : p(0:3) in boostx is zero momentum'
725 c endif
726 c if (abs(q(0))+qq.eq.rZero) then
727 c write(stdo,*)
728 c & ' helas-error : q(0:3) in boostx is zero momentum'
729 c endif
730 c if (p(0).le.rZero) then
731 c write(stdo,*)
732 c & ' helas-warn : p(0:3) in boostx has not positive energy'
733 c write(stdo,*)
734 c & ' : p(0) = ',p(0)
735 c endif
736 c if (q(0).le.rZero) then
737 c write(stdo,*)
738 c & ' helas-error : q(0:3) in boostx has not positive energy'
739 c write(stdo,*)
740 c & ' : q(0) = ',q(0)
741 c endif
742 c pp=p(0)**2-p(1)**2-p(2)**2-p(3)**2
743 c if (pp.lt.rZero) then
744 c write(stdo,*)
745 c & ' helas-warn : p(0:3) in boostx is spacelike'
746 c write(stdo,*)
747 c & ' : p**2 = ',pp
748 c endif
749 c if (q(0)**2-qq.le.rZero) then
750 c write(stdo,*)
751 c & ' helas-error : q(0:3) in boostx is not timelike'
752 c write(stdo,*)
753 c & ' : q**2 = ',q(0)**2-qq
754 c endif
755 c if (qq.eq.rZero) then
756 c write(stdo,*)
757 c & ' helas-warn : q(0:3) in boostx has zero spacial components'
758 c endif
759 c#endif
760 
761  if ( qq.ne.rzero ) then
762  pq = p(1)*q(1)+p(2)*q(2)+p(3)*q(3)
763  m = sqrt(max(q(0)**2-qq,1d-99))
764  lf = ((q(0)-m)*pq/qq+p(0))/m
765  pboost(0) = (p(0)*q(0)+pq)/m
766  pboost(1) = p(1)+q(1)*lf
767  pboost(2) = p(2)+q(2)*lf
768  pboost(3) = p(3)+q(3)*lf
769  else
770  pboost(0) = p(0)
771  pboost(1) = p(1)
772  pboost(2) = p(2)
773  pboost(3) = p(3)
774  endif
775 c
776  return
777  end
778 
779  subroutine momntx(energy,mass,costh,phi , p)
780 c
781 c This subroutine sets up a four-momentum from the four inputs.
782 c
783 c input:
784 c real energy : energy
785 c real mass : mass
786 c real costh : cos(theta)
787 c real phi : azimuthal angle
788 c
789 c output:
790 c real p(0:3) : four-momentum
791 c
792  implicit none
793  double precision p(0:3),energy,mass,costh,phi,pp,sinth
794 
795  double precision rzero, rone
796  parameter( rzero = 0.0d0, rone = 1.0d0 )
797 
798 c#ifdef HELAS_CHECK
799 c double precision rPi, rTwo
800 c parameter( rPi = 3.14159265358979323846d0, rTwo = 2.d0 )
801 c integer stdo
802 c parameter( stdo = 6 )
803 c#endif
804 c
805 c#ifdef HELAS_CHECK
806 c if (energy.lt.mass) then
807 c write(stdo,*)
808 c & ' helas-error : energy in momntx is less than mass'
809 c write(stdo,*)
810 c & ' : energy = ',energy,' : mass = ',mass
811 c endif
812 c if (mass.lt.rZero) then
813 c write(stdo,*) ' helas-error : mass in momntx is negative'
814 c write(stdo,*) ' : mass = ',mass
815 c endif
816 c if (abs(costh).gt.rOne) then
817 c write(stdo,*) ' helas-error : costh in momntx is improper'
818 c write(stdo,*) ' : costh = ',costh
819 c endif
820 c if (phi.lt.rZero .or. phi.gt.rTwo*rPi) then
821 c write(stdo,*)
822 c & ' helas-warn : phi in momntx does not lie on 0.0 thru 2.0*pi'
823 c write(stdo,*)
824 c & ' : phi = ',phi
825 c endif
826 c#endif
827 
828  p(0) = energy
829 
830  if ( energy.eq.mass ) then
831 
832  p(1) = rzero
833  p(2) = rzero
834  p(3) = rzero
835 
836  else
837 
838  pp = sqrt((energy-mass)*(energy+mass))
839  sinth = sqrt((rone-costh)*(rone+costh))
840  p(3) = pp*costh
841  if ( phi.eq.rzero ) then
842  p(1) = pp*sinth
843  p(2) = rzero
844  else
845  p(1) = pp*sinth*cos(phi)
846  p(2) = pp*sinth*sin(phi)
847  endif
848 
849  endif
850 c
851  return
852  end
853  subroutine rotxxx(p,q , prot)
854 c
855 c This subroutine performs the spacial rotation of a four-momentum.
856 c the momentum p is assumed to be given in the frame where the spacial
857 c component of q points the positive z-axis. prot is the momentum p
858 c rotated to the frame where q is given.
859 c
860 c input:
861 c real p(0:3) : four-momentum p in q(1)=q(2)=0 frame
862 c real q(0:3) : four-momentum q in the rotated frame
863 c
864 c output:
865 c real prot(0:3) : four-momentum p in the rotated frame
866 c
867  implicit none
868  double precision p(0:3),q(0:3),prot(0:3),qt2,qt,psgn,qq,p1
869 
870  double precision rzero, rone
871  parameter( rzero = 0.0d0, rone = 1.0d0 )
872 
873 c#ifdef HELAS_CHECK
874 c integer stdo
875 c parameter( stdo = 6 )
876 c#endif
877 c
878  prot(0) = p(0)
879 
880  qt2 = q(1)**2 + q(2)**2
881 
882 c#ifdef HELAS_CHECK
883 c if (abs(p(0))+abs(p(1))+abs(p(2))+abs(p(3)).eq.rZero) then
884 c write(stdo,*)
885 c & ' helas-error : p(0:3) in rotxxx is zero momentum'
886 c endif
887 c if (abs(q(0))+abs(q(3))+qt2.eq.rZero) then
888 c write(stdo,*)
889 c & ' helas-error : q(0:3) in rotxxx is zero momentum'
890 c endif
891 c if (qt2+abs(q(3)).eq.rZero) then
892 c write(stdo,*)
893 c & ' helas-warn : q(0:3) in rotxxx has zero spacial momentum'
894 c endif
895 c#endif
896 
897  if ( qt2.eq.rzero ) then
898  if ( q(3).eq.rzero ) then
899  prot(1) = p(1)
900  prot(2) = p(2)
901  prot(3) = p(3)
902  else
903  psgn = dsign(rone,q(3))
904  prot(1) = p(1)*psgn
905  prot(2) = p(2)*psgn
906  prot(3) = p(3)*psgn
907  endif
908  else
909  qq = sqrt(qt2+q(3)**2)
910  qt = sqrt(qt2)
911  p1 = p(1)
912  prot(1) = q(1)*q(3)/qq/qt*p1 -q(2)/qt*p(2) +q(1)/qq*p(3)
913  prot(2) = q(2)*q(3)/qq/qt*p1 +q(1)/qt*p(2) +q(2)/qq*p(3)
914  prot(3) = -qt/qq*p1 +q(3)/qq*p(3)
915  endif
916 c
917  return
918  end
919 
920  subroutine mom2cx(esum,mass1,mass2,costh1,phi1 , p1,p2)
921 c
922 c This subroutine sets up two four-momenta in the two particle rest
923 c frame.
924 c
925 c input:
926 c real esum : energy sum of particle 1 and 2
927 c real mass1 : mass of particle 1
928 c real mass2 : mass of particle 2
929 c real costh1 : cos(theta) of particle 1
930 c real phi1 : azimuthal angle of particle 1
931 c
932 c output:
933 c real p1(0:3) : four-momentum of particle 1
934 c real p2(0:3) : four-momentum of particle 2
935 c
936  implicit none
937  double precision p1(0:3),p2(0:3),
938  & esum,mass1,mass2,costh1,phi1,md2,ed,pp,sinth1
939 
940  double precision rzero, rhalf, rone, rtwo
941  parameter( rzero = 0.0d0, rhalf = 0.5d0 )
942  parameter( rone = 1.0d0, rtwo = 2.0d0 )
943 
944 c#ifdef HELAS_CHECK
945 c double precision rPi
946 c parameter( rPi = 3.14159265358979323846d0 )
947 c integer stdo
948 c parameter( stdo = 6 )
949 c#endif
950 cc
951 c#ifdef HELAS_CHECK
952 c if (esum.lt.mass1+mass2) then
953 c write(stdo,*)
954 c & ' helas-error : esum in mom2cx is less than mass1+mass2'
955 c write(stdo,*)
956 c & ' : esum = ',esum,
957 c & ' : mass1+mass2 = ',mass1,mass2
958 c endif
959 c if (mass1.lt.rZero) then
960 c write(stdo,*) ' helas-error : mass1 in mom2cx is negative'
961 c write(stdo,*) ' : mass1 = ',mass1
962 c endif
963 c if (mass2.lt.rZero) then
964 c write(stdo,*) ' helas-error : mass2 in mom2cx is negative'
965 c write(stdo,*) ' : mass2 = ',mass2
966 c endif
967 c if (abs(costh1).gt.1.) then
968 c write(stdo,*) ' helas-error : costh1 in mom2cx is improper'
969 c write(stdo,*) ' : costh1 = ',costh1
970 c endif
971 c if (phi1.lt.rZero .or. phi1.gt.rTwo*rPi) then
972 c write(stdo,*)
973 c & ' helas-warn : phi1 in mom2cx does not lie on 0.0 thru 2.0*pi'
974 c write(stdo,*)
975 c & ' : phi1 = ',phi1
976 c endif
977 c#endif
978 
979  md2 = (mass1-mass2)*(mass1+mass2)
980  ed = md2/esum
981  if ( mass1*mass2.eq.rzero ) then
982  pp = (esum-abs(ed))*rhalf
983  else
984  pp = sqrt(max((md2/esum)**2-rtwo*(mass1**2+mass2**2)+esum**2,1d-99))*rhalf
985  endif
986  sinth1 = sqrt((rone-costh1)*(rone+costh1))
987 
988  p1(0) = max((esum+ed)*rhalf,rzero)
989  p1(1) = pp*sinth1*cos(phi1)
990  p1(2) = pp*sinth1*sin(phi1)
991  p1(3) = pp*costh1
992 
993  p2(0) = max((esum-ed)*rhalf,rzero)
994  p2(1) = -p1(1)
995  p2(2) = -p1(2)
996  p2(3) = -p1(3)
997 c
998  return
999  end
1000  subroutine irxxxx(p,rmass,nhel,nsr , ri)
1001 c
1002 c This subroutine computes a Rarita-Schwinger wavefunction of spin-3/2
1003 c fermion with the flowing-IN fermion number.
1004 c
1005 c input:
1006 c real p(0:3) : four-momentum of RS fermion
1007 c real rmass : mass of RS fermion
1008 c integer nhel = -3,-1,1,3 : helicity of RS fermion
1009 c (1- and 1 is forbidden if rmass = 0)
1010 c integer nsr = -1 or 1 : +1 for particle, -1 for anti-particle
1011 c
1012 c output:
1013 c complex ri(18) : RS fermion wavefunction |ri>v
1014 c
1015 c- by K.Mawatari - 2008/02/26
1016 c
1017  implicit none
1018  double precision p(0:3),rmass
1019  integer nhel,nsr
1020  double complex ri(18)
1021 
1022  double complex rc(6,4),ep(4),em(4),e0(4),fip(4),fim(4),chi(2)
1023  double precision pp,pt2,pt,pzpt,emp, sf(2),sfomeg(2),omega(2),pp3,
1024  & sqp0p3,sqm
1025  integer i,j,nsv,ip,im,nh
1026 
1027  double precision rzero, rhalf, rone, rtwo, rthree, sqh,sq2,sq3
1028  parameter( rzero = 0.0d0, rhalf = 0.5d0 )
1029  parameter( rone = 1.0d0, rtwo = 2.0d0, rthree = 3.0d0 )
1030 
1031 c#ifdef HELAS_CHECK
1032 c double precision p2
1033 c double precision epsi
1034 c parameter( epsi = 2.0d-5 )
1035 c integer stdo
1036 c parameter( stdo = 6 )
1037 c#endif
1038 c
1039 c#ifdef HELAS_CHECK
1040 c pp = sqrt(p(1)**2+p(2)**2+p(3)**2)
1041 c if ( abs(p(0))+pp.eq.rZero ) then
1042 c write(stdo,*)
1043 c & ' helas-error : p(0:3) in irxxxx is zero momentum'
1044 c endif
1045 c if ( p(0).le.rZero ) then
1046 c write(stdo,*)
1047 c & ' helas-error : p(0:3) in irxxxx has non-positive energy'
1048 c write(stdo,*)
1049 c & ' : p(0) = ',p(0)
1050 c endif
1051 c p2 = (p(0)-pp)*(p(0)+pp)
1052 c if ( abs(p2-rmass**2).gt.p(0)**2*epsi ) then
1053 c write(stdo,*)
1054 c & ' helas-error : p(0:3) in irxxxx has inappropriate mass'
1055 c write(stdo,*)
1056 c & ' : p**2 = ',p2,' : rmass**2 = ',rmass**2
1057 c endif
1058 c if (abs(nhel).gt.3 .or. abs(nhel).eq.2 .or. abs(nhel).eq.0 ) then
1059 c write(stdo,*) ' helas-error : nhel in irxxxx is not -3,-1,1,3'
1060 c write(stdo,*) ' : nhel = ',nhel
1061 c endif
1062 c if (abs(nsr).ne.1) then
1063 c write(stdo,*) ' helas-error : nsr in irxxxx is not -1,1'
1064 c write(stdo,*) ' : nsr = ',nsr
1065 c endif
1066 c#endif
1067 
1068  sqh = sqrt(rhalf)
1069  sq2 = sqrt(rtwo)
1070  sq3 = sqrt(rthree)
1071 
1072  pt2 = p(1)**2 + p(2)**2
1073  pp = min(p(0),sqrt(pt2+p(3)**2))
1074  pt = min(pp,sqrt(pt2))
1075 
1076  rc(5,1) = -1*dcmplx(p(0),p(3))*nsr
1077  rc(6,1) = -1*dcmplx(p(1),p(2))*nsr
1078 
1079  nsv = -nsr ! nsv=+1 for final, -1 for initial
1080 
1081  if ( nhel.ge.1 ) then
1082 c construct eps+
1083  if ( pp.eq.rzero ) then
1084  ep(1) = dcmplx( rzero )
1085  ep(2) = dcmplx( -sqh )
1086  ep(3) = dcmplx( rzero , nsv*sqh )
1087  ep(4) = dcmplx( rzero )
1088  else
1089  ep(1) = dcmplx( rzero )
1090  ep(4) = dcmplx( pt/pp*sqh )
1091  if ( pt.ne.rzero ) then
1092  pzpt = p(3)/(pp*pt)*sqh
1093  ep(2) = dcmplx( -p(1)*pzpt , -nsv*p(2)/pt*sqh )
1094  ep(3) = dcmplx( -p(2)*pzpt , nsv*p(1)/pt*sqh )
1095  else
1096  ep(2) = dcmplx( -sqh )
1097  ep(3) = dcmplx( rzero , nsv*sign(sqh,p(3)) )
1098  endif
1099  endif
1100  end if
1101 
1102  if ( nhel.le.-1 ) then
1103 c construct eps-
1104  if ( pp.eq.rzero ) then
1105  em(1) = dcmplx( rzero )
1106  em(2) = dcmplx( sqh )
1107  em(3) = dcmplx( rzero , nsv*sqh )
1108  em(4) = dcmplx( rzero )
1109  else
1110  em(1) = dcmplx( rzero )
1111  em(4) = dcmplx( -pt/pp*sqh )
1112  if ( pt.ne.rzero ) then
1113  pzpt = -p(3)/(pp*pt)*sqh
1114  em(2) = dcmplx( -p(1)*pzpt , -nsv*p(2)/pt*sqh )
1115  em(3) = dcmplx( -p(2)*pzpt , nsv*p(1)/pt*sqh )
1116  else
1117  em(2) = dcmplx( sqh )
1118  em(3) = dcmplx( rzero , nsv*sign(sqh,p(3)) )
1119  endif
1120  endif
1121  end if
1122 
1123  if ( abs(nhel).le.1 ) then
1124 c construct eps0
1125  if ( pp.eq.rzero ) then
1126  e0(1) = dcmplx( rzero )
1127  e0(2) = dcmplx( rzero )
1128  e0(3) = dcmplx( rzero )
1129  e0(4) = dcmplx( rone )
1130  else
1131  emp = p(0)/(rmass*pp)
1132  e0(1) = dcmplx( pp/rmass )
1133  e0(4) = dcmplx( p(3)*emp )
1134  if ( pt.ne.rzero ) then
1135  e0(2) = dcmplx( p(1)*emp )
1136  e0(3) = dcmplx( p(2)*emp )
1137  else
1138  e0(2) = dcmplx( rzero )
1139  e0(3) = dcmplx( rzero )
1140  endif
1141  end if
1142  end if
1143 
1144  if ( nhel.ge.-1 ) then
1145 c constract spinor+
1146  nh = nsr
1147  if ( rmass.ne.rzero ) then
1148  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
1149  if ( pp.eq.rzero ) then
1150  sqm = dsqrt(rmass)
1151  ip = (1+nh)/2
1152  im = (1-nh)/2
1153  fip(1) = ip * sqm
1154  fip(2) = im*nsr * sqm
1155  fip(3) = ip*nsr * sqm
1156  fip(4) = im * sqm
1157  else
1158  sf(1) = dble(1+nsr+(1-nsr)*nh)*rhalf
1159  sf(2) = dble(1+nsr-(1-nsr)*nh)*rhalf
1160  omega(1) = dsqrt(p(0)+pp)
1161  omega(2) = rmass/omega(1)
1162  ip = (3+nh)/2
1163  im = (3-nh)/2
1164  sfomeg(1) = sf(1)*omega(ip)
1165  sfomeg(2) = sf(2)*omega(im)
1166  pp3 = max(pp+p(3),rzero)
1167  chi(1) = dcmplx( dsqrt(pp3*rhalf/pp) )
1168  if ( pp3.eq.rzero ) then
1169  chi(2) = dcmplx(-nh )
1170  else
1171  chi(2) = dcmplx( nh*p(1) , p(2) )/dsqrt(rtwo*pp*pp3)
1172  endif
1173  fip(1) = sfomeg(1)*chi(im)
1174  fip(2) = sfomeg(1)*chi(ip)
1175  fip(3) = sfomeg(2)*chi(im)
1176  fip(4) = sfomeg(2)*chi(ip)
1177  endif
1178  else
1179  sqp0p3 = dsqrt(max(p(0)+p(3),rzero))*nsr
1180  chi(1) = dcmplx( sqp0p3 )
1181  if ( sqp0p3.eq.rzero ) then
1182  chi(2) = dcmplx(-nhel )*dsqrt(rtwo*p(0))
1183  else
1184  chi(2) = dcmplx( nh*p(1), p(2) )/sqp0p3
1185  endif
1186  if ( nh.eq.1 ) then
1187  fip(1) = dcmplx( rzero )
1188  fip(2) = dcmplx( rzero )
1189  fip(3) = chi(1)
1190  fip(4) = chi(2)
1191  else
1192  fip(1) = chi(2)
1193  fip(2) = chi(1)
1194  fip(3) = dcmplx( rzero )
1195  fip(4) = dcmplx( rzero )
1196  endif
1197  endif
1198  end if
1199 
1200  if ( nhel.le.1 ) then
1201 c constract spinor-
1202  nh = -nsr
1203  if ( rmass.ne.rzero ) then
1204  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
1205  if ( pp.eq.rzero ) then
1206  sqm = dsqrt(rmass)
1207  ip = (1+nh)/2
1208  im = (1-nh)/2
1209  fim(1) = ip * sqm
1210  fim(2) = im*nsr * sqm
1211  fim(3) = ip*nsr * sqm
1212  fim(4) = im * sqm
1213  else
1214  sf(1) = dble(1+nsr+(1-nsr)*nh)*rhalf
1215  sf(2) = dble(1+nsr-(1-nsr)*nh)*rhalf
1216  omega(1) = dsqrt(p(0)+pp)
1217  omega(2) = rmass/omega(1)
1218  ip = (3+nh)/2
1219  im = (3-nh)/2
1220  sfomeg(1) = sf(1)*omega(ip)
1221  sfomeg(2) = sf(2)*omega(im)
1222  pp3 = max(pp+p(3),rzero)
1223  chi(1) = dcmplx( dsqrt(pp3*rhalf/pp) )
1224  if ( pp3.eq.rzero ) then
1225  chi(2) = dcmplx(-nh )
1226  else
1227  chi(2) = dcmplx( nh*p(1) , p(2) )/dsqrt(rtwo*pp*pp3)
1228  endif
1229  fim(1) = sfomeg(1)*chi(im)
1230  fim(2) = sfomeg(1)*chi(ip)
1231  fim(3) = sfomeg(2)*chi(im)
1232  fim(4) = sfomeg(2)*chi(ip)
1233  endif
1234  else
1235  sqp0p3 = dsqrt(max(p(0)+p(3),rzero))*nsr
1236  chi(1) = dcmplx( sqp0p3 )
1237  if ( sqp0p3.eq.rzero ) then
1238  chi(2) = dcmplx(-nhel )*dsqrt(rtwo*p(0))
1239  else
1240  chi(2) = dcmplx( nh*p(1), p(2) )/sqp0p3
1241  endif
1242  if ( nh.eq.1 ) then
1243  fim(1) = dcmplx( rzero )
1244  fim(2) = dcmplx( rzero )
1245  fim(3) = chi(1)
1246  fim(4) = chi(2)
1247  else
1248  fim(1) = chi(2)
1249  fim(2) = chi(1)
1250  fim(3) = dcmplx( rzero )
1251  fim(4) = dcmplx( rzero )
1252  endif
1253  endif
1254  end if
1255 
1256 c spin-3/2 fermion wavefunction
1257  if ( nhel.eq.3 ) then
1258  do j = 1,4
1259  do i = 1,4
1260  rc(i,j) = ep(i)*fip(j)
1261  end do
1262  end do
1263  else if ( nhel.eq.1 ) then
1264  do j = 1,4
1265  do i = 1,4
1266  if ( pt.eq.rzero .and. p(3).ge.0d0 ) then
1267  rc(i,j) = sq2/sq3*e0(i)*fip(j) +rone/sq3*ep(i)*fim(j)
1268  elseif ( pt.eq.rzero .and. p(3).lt.0d0 ) then
1269  rc(i,j) = sq2/sq3*e0(i)*fip(j) -rone/sq3*ep(i)*fim(j)
1270  else
1271  rc(i,j) = sq2/sq3*e0(i)*fip(j)
1272  & +rone/sq3*ep(i)*fim(j) *dcmplx(p(1),nsr*p(2))/pt
1273  endif
1274  end do
1275  end do
1276  else if ( nhel.eq.-1 ) then
1277  do j = 1,4
1278  do i = 1,4
1279  if ( pt.eq.rzero .and.p(3).ge.0d0 ) then
1280  rc(i,j) = rone/sq3*em(i)*fip(j) +sq2/sq3*e0(i)*fim(j)
1281  elseif ( pt.eq.rzero .and.p(3).lt.0d0 ) then
1282  rc(i,j) = rone/sq3*em(i)*fip(j) -sq2/sq3*e0(i)*fim(j)
1283  else
1284  rc(i,j) = rone/sq3*em(i)*fip(j)
1285  & + sq2/sq3*e0(i)*fim(j) *dcmplx(p(1),nsr*p(2))/pt
1286  endif
1287  end do
1288  end do
1289  else
1290  do j = 1,4
1291  do i = 1,4
1292  if ( pt.eq.rzero .and. p(3).ge.0d0 ) then
1293  rc(i,j) = em(i)*fim(j)
1294  elseif ( pt.eq.rzero .and. p(3).lt.0d0 ) then
1295  rc(i,j) = -em(i)*fim(j)
1296  else
1297  rc(i,j) = em(i)*fim(j) *dcmplx(p(1),nsr*p(2))/pt
1298  endif
1299  end do
1300  end do
1301  end if
1302 
1303  ri(3) = rc(1,1)
1304  ri(4) = rc(1,2)
1305  ri(5) = rc(1,3)
1306  ri(6) = rc(1,4)
1307  ri(7) = rc(2,1)
1308  ri(8) = rc(2,2)
1309  ri(9) = rc(2,3)
1310  ri(10) = rc(2,4)
1311  ri(11) = rc(3,1)
1312  ri(12) = rc(3,2)
1313  ri(13) = rc(3,3)
1314  ri(14) = rc(3,4)
1315  ri(15) = rc(4,1)
1316  ri(16) = rc(4,2)
1317  ri(17) = rc(4,3)
1318  ri(18) = rc(4,4)
1319  ri(1) = rc(5,1)
1320  ri(2) = rc(6,1)
1321 
1322  return
1323  end
1324  subroutine orxxxx(p,rmass,nhel,nsr , ro)
1325 c
1326 c This subroutine computes a Rarita-Schwinger wavefunction of spin-3/2
1327 c fermion with the flowing-IN fermion number.
1328 c
1329 c input:
1330 c real p(0:3) : four-momentum of RS fermion
1331 c real rmass : mass of RS fermion
1332 c integer nhel = -3,-1,1,3 : helicity of RS fermion
1333 c (1- and 1 is forbidden if rmass = 0)
1334 c integer nsr = -1 or 1 : +1 for particle, -1 for anti-particle
1335 c
1336 c output:
1337 c complex ro(18) : RS fermion wavefunction |ro>v
1338 c
1339 c- by Y.Takaesu - 2011/01/11
1340 c
1341  implicit none
1342  double precision p(0:3),rmass
1343  integer nhel,nsr
1344  double complex ro(18),fipp(4),fimm(4)
1345 
1346  double complex rc(6,4),ep(4),em(4),e0(4),fop(4),fom(4),chi(2)
1347  double precision pp,pt2,pt,pzpt,emp, sf(2),sfomeg(2),omega(2),pp3,
1348  & sqp0p3,sqm(0:1)
1349  integer i,j,nsv,ip,im,nh
1350 
1351  double precision rzero, rhalf, rone, rtwo, rthree, sqh,sq2,sq3
1352  parameter( rzero = 0.0d0, rhalf = 0.5d0 )
1353  parameter( rone = 1.0d0, rtwo = 2.0d0, rthree = 3.0d0 )
1354 
1355 c#ifdef HELAS_CHECK
1356 c double precision p2
1357 c double precision epsi
1358 c parameter( epsi = 2.0d-5 )
1359 c integer stdo
1360 c parameter( stdo = 6 )
1361 c#endif
1362 c
1363 c#ifdef HELAS_CHECK
1364 c pp = sqrt(p(1)**2+p(2)**2+p(3)**2)
1365 c if ( abs(p(0))+pp.eq.rZero ) then
1366 c write(stdo,*)
1367 c & ' helas-error : p(0:3) in orxxxx is zero momentum'
1368 c endif
1369 c if ( p(0).le.rZero ) then
1370 c write(stdo,*)
1371 c & ' helas-error : p(0:3) in orxxxx has non-positive energy'
1372 c write(stdo,*)
1373 c & ' : p(0) = ',p(0)
1374 c endif
1375 c p2 = (p(0)-pp)*(p(0)+pp)
1376 c if ( abs(p2-rmass**2).gt.p(0)**2*epsi ) then
1377 c write(stdo,*)
1378 c & ' helas-error : p(0:3) in orxxxx has inappropriate mass'
1379 c write(stdo,*)
1380 c & ' : p**2 = ',p2,' : rmass**2 = ',rmass**2
1381 c endif
1382 c if (abs(nhel).gt.3 .or. abs(nhel).eq.2 .or. abs(nhel).eq.0 ) then
1383 c write(stdo,*) ' helas-error : nhel in orxxxx is not -3,-1,1,3'
1384 c write(stdo,*) ' : nhel = ',nhel
1385 c endif
1386 c if (abs(nsr).ne.1) then
1387 c write(stdo,*) ' helas-error : nsr in orxxxx is not -1,1'
1388 c write(stdo,*) ' : nsr = ',nsr
1389 c endif
1390 c#endif
1391 
1392  sqh = sqrt(rhalf)
1393  sq2 = sqrt(rtwo)
1394  sq3 = sqrt(rthree)
1395 
1396  pt2 = p(1)**2 + p(2)**2
1397  pp = min(p(0),sqrt(pt2+p(3)**2))
1398  pt = min(pp,sqrt(pt2))
1399 
1400  rc(5,1) = dcmplx(p(0),p(3))*nsr
1401  rc(6,1) = dcmplx(p(1),p(2))*nsr
1402 
1403  nsv = nsr ! nsv=+1 for final, -1 for initial
1404 
1405  if ( nhel.ge.1 ) then
1406 c construct eps+
1407  if ( pp.eq.rzero ) then
1408  ep(1) = dcmplx( rzero )
1409  ep(2) = dcmplx( -sqh )
1410  ep(3) = dcmplx( rzero , nsv*sqh )
1411  ep(4) = dcmplx( rzero )
1412  else
1413  ep(1) = dcmplx( rzero )
1414  ep(4) = dcmplx( pt/pp*sqh )
1415  if ( pt.ne.rzero ) then
1416  pzpt = p(3)/(pp*pt)*sqh
1417  ep(2) = dcmplx( -p(1)*pzpt , -nsv*p(2)/pt*sqh )
1418  ep(3) = dcmplx( -p(2)*pzpt , nsv*p(1)/pt*sqh )
1419  else
1420  ep(2) = dcmplx( -sqh )
1421  ep(3) = dcmplx( rzero , nsv*sign(sqh,p(3)) )
1422  endif
1423  endif
1424  end if
1425 
1426  if ( nhel.le.-1 ) then
1427 c construct eps-
1428  if ( pp.eq.rzero ) then
1429  em(1) = dcmplx( rzero )
1430  em(2) = dcmplx( sqh )
1431  em(3) = dcmplx( rzero , nsv*sqh )
1432  em(4) = dcmplx( rzero )
1433  else
1434  em(1) = dcmplx( rzero )
1435  em(4) = dcmplx( -pt/pp*sqh )
1436  if ( pt.ne.rzero ) then
1437  pzpt = -p(3)/(pp*pt)*sqh
1438  em(2) = dcmplx( -p(1)*pzpt , -nsv*p(2)/pt*sqh )
1439  em(3) = dcmplx( -p(2)*pzpt , nsv*p(1)/pt*sqh )
1440  else
1441  em(2) = dcmplx( sqh )
1442  em(3) = dcmplx( rzero , nsv*sign(sqh,p(3)) )
1443  endif
1444  endif
1445  end if
1446 
1447  if ( abs(nhel).le.1 ) then
1448 c construct eps0
1449  if ( pp.eq.rzero ) then
1450  e0(1) = dcmplx( rzero )
1451  e0(2) = dcmplx( rzero )
1452  e0(3) = dcmplx( rzero )
1453  e0(4) = dcmplx( rone )
1454  else
1455  emp = p(0)/(rmass*pp)
1456  e0(1) = dcmplx( pp/rmass )
1457  e0(4) = dcmplx( p(3)*emp )
1458  if ( pt.ne.rzero ) then
1459  e0(2) = dcmplx( p(1)*emp )
1460  e0(3) = dcmplx( p(2)*emp )
1461  else
1462  e0(2) = dcmplx( rzero )
1463  e0(3) = dcmplx( rzero )
1464  endif
1465  end if
1466  end if
1467 
1468  if ( nhel.ge.-1 ) then
1469 c constract spinor+
1470  nh = nsr
1471 
1472  if ( rmass.ne.rzero ) then
1473 
1474  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
1475 
1476  if ( pp.eq.rzero ) then
1477 
1478  sqm(0) = dsqrt(abs(rmass)) ! possibility of negative fermion masses
1479  sqm(1) = sign(sqm(0),rmass) ! possibility of negative fermion masses
1480  ip = -((1+nh)/2)
1481  im = (1-nh)/2
1482 
1483  fop(1) = im * sqm(im)
1484  fop(2) = ip*nsr * sqm(im)
1485  fop(3) = im*nsr * sqm(-ip)
1486  fop(4) = ip * sqm(-ip)
1487 
1488  else
1489 
1490  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
1491  sf(1) = dble(1+nsr+(1-nsr)*nh)*rhalf
1492  sf(2) = dble(1+nsr-(1-nsr)*nh)*rhalf
1493  omega(1) = dsqrt(p(0)+pp)
1494  omega(2) = rmass/omega(1)
1495  ip = (3+nh)/2
1496  im = (3-nh)/2
1497  sfomeg(1) = sf(1)*omega(ip)
1498  sfomeg(2) = sf(2)*omega(im)
1499  pp3 = max(pp+p(3),rzero)
1500  chi(1) = dcmplx( dsqrt(pp3*rhalf/pp) )
1501  if ( pp3.eq.rzero ) then
1502  chi(2) = dcmplx(-nh )
1503  else
1504  chi(2) = dcmplx( nh*p(1) , -p(2) )/dsqrt(rtwo*pp*pp3)
1505  endif
1506 
1507  fop(1) = sfomeg(2)*chi(im)
1508  fop(2) = sfomeg(2)*chi(ip)
1509  fop(3) = sfomeg(1)*chi(im)
1510  fop(4) = sfomeg(1)*chi(ip)
1511 
1512  endif
1513 
1514  else
1515 
1516  if(p(1).eq.0d0.and.p(2).eq.0d0.and.p(3).lt.0d0) then
1517  sqp0p3 = 0d0
1518  else
1519  sqp0p3 = dsqrt(max(p(0)+p(3),rzero))*nsr
1520  end if
1521  chi(1) = dcmplx( sqp0p3 )
1522  if ( sqp0p3.eq.rzero ) then
1523  chi(2) = dcmplx(-nhel )*dsqrt(rtwo*p(0))
1524  else
1525  chi(2) = dcmplx( nh*p(1), -p(2) )/sqp0p3
1526  endif
1527  if ( nh.eq.1 ) then
1528  fop(1) = chi(1)
1529  fop(2) = chi(2)
1530  fop(3) = dcmplx( rzero )
1531  fop(4) = dcmplx( rzero )
1532  else
1533  fop(1) = dcmplx( rzero )
1534  fop(2) = dcmplx( rzero )
1535  fop(3) = chi(2)
1536  fop(4) = chi(1)
1537  endif
1538  endif
1539  endif
1540 
1541  if ( nhel.le.1 ) then
1542 c constract spinor+
1543  nh = -nsr
1544 
1545  if ( rmass.ne.rzero ) then
1546 
1547  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
1548 
1549  if ( pp.eq.rzero ) then
1550 
1551  sqm(0) = dsqrt(abs(rmass)) ! possibility of negative fermion masses
1552  sqm(1) = sign(sqm(0),rmass) ! possibility of negative fermion masses
1553  ip = -((1+nh)/2)
1554  im = (1-nh)/2
1555 
1556  fom(1) = im * sqm(im)
1557  fom(2) = ip*nsr * sqm(im)
1558  fom(3) = im*nsr * sqm(-ip)
1559  fom(4) = ip * sqm(-ip)
1560 
1561  else
1562 
1563  pp = min(p(0),dsqrt(p(1)**2+p(2)**2+p(3)**2))
1564  sf(1) = dble(1+nsr+(1-nsr)*nh)*rhalf
1565  sf(2) = dble(1+nsr-(1-nsr)*nh)*rhalf
1566  omega(1) = dsqrt(p(0)+pp)
1567  omega(2) = rmass/omega(1)
1568  ip = (3+nh)/2
1569  im = (3-nh)/2
1570  sfomeg(1) = sf(1)*omega(ip)
1571  sfomeg(2) = sf(2)*omega(im)
1572  pp3 = max(pp+p(3),rzero)
1573  chi(1) = dcmplx( dsqrt(pp3*rhalf/pp) )
1574  if ( pp3.eq.rzero ) then
1575  chi(2) = dcmplx(-nh )
1576  else
1577  chi(2) = dcmplx( nh*p(1) , -p(2) )/dsqrt(rtwo*pp*pp3)
1578  endif
1579 
1580  fom(1) = sfomeg(2)*chi(im)
1581  fom(2) = sfomeg(2)*chi(ip)
1582  fom(3) = sfomeg(1)*chi(im)
1583  fom(4) = sfomeg(1)*chi(ip)
1584 
1585  endif
1586 
1587  else
1588 
1589  if(p(1).eq.0d0.and.p(2).eq.0d0.and.p(3).lt.0d0) then
1590  sqp0p3 = 0d0
1591  else
1592  sqp0p3 = dsqrt(max(p(0)+p(3),rzero))*nsr
1593  end if
1594  chi(1) = dcmplx( sqp0p3 )
1595  if ( sqp0p3.eq.rzero ) then
1596  chi(2) = dcmplx(-nhel )*dsqrt(rtwo*p(0))
1597  else
1598  chi(2) = dcmplx( nh*p(1), -p(2) )/sqp0p3
1599  endif
1600  if ( nh.eq.1 ) then
1601  fom(1) = chi(1)
1602  fom(2) = chi(2)
1603  fom(3) = dcmplx( rzero )
1604  fom(4) = dcmplx( rzero )
1605  else
1606  fom(1) = dcmplx( rzero )
1607  fom(2) = dcmplx( rzero )
1608  fom(3) = chi(2)
1609  fom(4) = chi(1)
1610  endif
1611  endif
1612  endif
1613 
1614 c spin-3/2 fermion wavefunction
1615  if ( nhel.eq.3 ) then
1616  do j = 1,4
1617  do i = 1,4
1618  rc(i,j) = ep(i)*fop(j)
1619  end do
1620  end do
1621  else if ( nhel.eq.1 ) then
1622  do j = 1,4
1623  do i = 1,4
1624  if ( pt.eq.rzero .and. p(3).ge.0d0 ) then
1625  rc(i,j) = sq2/sq3*e0(i)*fop(j)
1626  & +rone/sq3*ep(i)*fom(j)
1627  elseif ( pt.eq.rzero .and. p(3).lt.0d0 ) then
1628  rc(i,j) = sq2/sq3*e0(i)*fop(j)
1629  & -rone/sq3*ep(i)*fom(j)
1630  else
1631  rc(i,j) = sq2/sq3*e0(i)*fop(j)
1632  & +rone/sq3*ep(i)*fom(j)
1633  & *dcmplx(p(1),-nsr*p(2))/pt
1634  endif
1635  end do
1636  end do
1637  else if ( nhel.eq.-1 ) then
1638  do j = 1,4
1639  do i = 1,4
1640  if ( pt.eq.rzero .and.p(3).ge.0d0 ) then
1641  rc(i,j) = rone/sq3*em(i)*fop(j)
1642  & +sq2/sq3*e0(i)*fom(j)
1643  elseif ( pt.eq.rzero .and.p(3).lt.0d0 ) then
1644  rc(i,j) = rone/sq3*em(i)*fop(j)
1645  & -sq2/sq3*e0(i)*fom(j)
1646  else
1647  rc(i,j) = rone/sq3*em(i)*fop(j)
1648  & + sq2/sq3*e0(i)*fom(j)
1649  & *dcmplx(p(1),-nsr*p(2))/pt
1650  endif
1651  end do
1652  end do
1653  else
1654  do j = 1,4
1655  do i = 1,4
1656  if ( pt.eq.rzero .and. p(3).ge.0d0 ) then
1657  rc(i,j) = em(i)*fom(j)
1658  elseif ( pt.eq.rzero .and. p(3).lt.0d0 ) then
1659  rc(i,j) = -em(i)*fom(j)
1660  else
1661  rc(i,j) = em(i)*fom(j)*dcmplx(p(1),-nsr*p(2))/pt
1662  endif
1663  end do
1664  end do
1665  end if
1666 
1667  ro(3) = rc(1,1)
1668  ro(4) = rc(1,2)
1669  ro(5) = rc(1,3)
1670  ro(6) = rc(1,4)
1671  ro(7) = rc(2,1)
1672  ro(8) = rc(2,2)
1673  ro(9) = rc(2,3)
1674  ro(10) = rc(2,4)
1675  ro(11) = rc(3,1)
1676  ro(12) = rc(3,2)
1677  ro(13) = rc(3,3)
1678  ro(14) = rc(3,4)
1679  ro(15) = rc(4,1)
1680  ro(16) = rc(4,2)
1681  ro(17) = rc(4,3)
1682  ro(18) = rc(4,4)
1683  ro(1) = rc(5,1)
1684  ro(2) = rc(6,1)
1685 
1686  return
1687  end
1688 
1689 
1690 
1691 C This File is Automatically generated by ALOHA
1692 C The process calculated in this file is:
1693 C Gamma(3,2,1)
1694 C
1695  SUBROUTINE ffv1_0(F1, F2, V3, COUP,VERTEX)
1696  IMPLICIT NONE
1697  COMPLEX*16 ci
1698  parameter(ci=(0d0,1d0))
1699  COMPLEX*16 tmp15
1700  COMPLEX*16 v3(*)
1701  COMPLEX*16 f1(*)
1702  COMPLEX*16 f2(*)
1703  COMPLEX*16 vertex
1704  COMPLEX*16 coup
1705  tmp15 = (f1(3)*(f2(5)*(v3(3)+v3(6))+f2(6)*(v3(4)+ci*(v3(5))))
1706  $ +(f1(4)*(f2(5)*(v3(4)-ci*(v3(5)))+f2(6)*(v3(3)-v3(6)))
1707  $ +(f1(5)*(f2(3)*(v3(3)-v3(6))-f2(4)*(v3(4)+ci*(v3(5))))
1708  $ +f1(6)*(f2(3)*(+ci*(v3(5))-v3(4))+f2(4)*(v3(3)+v3(6))))))
1709  vertex = coup*(-ci) * tmp15
1710  END
1711 
1712 C This File is Automatically generated by ALOHA
1713 C The process calculated in this file is:
1714 C Gamma(3,2,1)
1715 C
1716  SUBROUTINE ffv1_1(F2, V3, COUP, M1, W1,F1)
1717  IMPLICIT NONE
1718  COMPLEX*16 ci
1719  parameter(ci=(0d0,1d0))
1720  COMPLEX*16 f2(*)
1721  COMPLEX*16 v3(*)
1722  REAL*8 p1(0:3)
1723  REAL*8 m1
1724  REAL*8 w1
1725  COMPLEX*16 f1(6)
1726  COMPLEX*16 denom
1727  COMPLEX*16 coup
1728  f1(1) = +f2(1)+v3(1)
1729  f1(2) = +f2(2)+v3(2)
1730  p1(0) = -dble(f1(1))
1731  p1(1) = -dble(f1(2))
1732  p1(2) = -dimag(f1(2))
1733  p1(3) = -dimag(f1(1))
1734  denom = coup/(p1(0)**2-p1(1)**2-p1(2)**2-p1(3)**2 - m1 * (m1
1735  $ -ci* w1))
1736  f1(3)= denom*ci*(f2(3)*(p1(0)*(v3(6)-v3(3))+(p1(1)*(v3(4)
1737  $ -ci*(v3(5)))+(p1(2)*(v3(5)+ci*(v3(4)))+p1(3)*(v3(6)-v3(3)))))
1738  $ +(f2(4)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*(-1d0)*(v3(3)+v3(6))
1739  $ +(p1(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p1(3)*(v3(4)+ci*(v3(5))))))
1740  $ +m1*(f2(5)*(v3(3)+v3(6))+f2(6)*(v3(4)+ci*(v3(5))))))
1741  f1(4)= denom*(-ci)*(f2(3)*(p1(0)*(+ci*(v3(5))-v3(4))+(p1(1)*(v3(3)
1742  $ -v3(6))+(p1(2)*(-ci*(v3(3))+ci*(v3(6)))+p1(3)*(v3(4)-ci
1743  $ *(v3(5))))))+(f2(4)*(p1(0)*(v3(3)+v3(6))+(p1(1)*(-1d0)*(v3(4)
1744  $ +ci*(v3(5)))+(p1(2)*(+ci*(v3(4))-v3(5))-p1(3)*(v3(3)+v3(6)))))
1745  $ +m1*(f2(5)*(+ci*(v3(5))-v3(4))+f2(6)*(v3(6)-v3(3)))))
1746  f1(5)= denom*(-ci)*(f2(5)*(p1(0)*(v3(3)+v3(6))+(p1(1)*(+ci*(v3(5))
1747  $ -v3(4))+(p1(2)*(-1d0)*(v3(5)+ci*(v3(4)))-p1(3)*(v3(3)+v3(6)))))
1748  $ +(f2(6)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*(v3(6)-v3(3))
1749  $ +(p1(2)*(-ci*(v3(3))+ci*(v3(6)))-p1(3)*(v3(4)+ci*(v3(5))))))
1750  $ +m1*(f2(3)*(v3(6)-v3(3))+f2(4)*(v3(4)+ci*(v3(5))))))
1751  f1(6)= denom*ci*(f2(5)*(p1(0)*(+ci*(v3(5))-v3(4))+(p1(1)*(v3(3)
1752  $ +v3(6))+(p1(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p1(3)*(+ci*(v3(5))
1753  $ -v3(4)))))+(f2(6)*(p1(0)*(v3(6)-v3(3))+(p1(1)*(v3(4)+ci
1754  $ *(v3(5)))+(p1(2)*(v3(5)-ci*(v3(4)))+p1(3)*(v3(6)-v3(3)))))
1755  $ +m1*(f2(3)*(+ci*(v3(5))-v3(4))+f2(4)*(v3(3)+v3(6)))))
1756  END
1757 
1758 
1759 
1760 C This File is Automatically generated by ALOHA
1761 C The process calculated in this file is:
1762 C Gamma(3,2,1)
1763 C
1764  SUBROUTINE ffv1_2(F1, V3, COUP, M2, W2,F2)
1765  IMPLICIT NONE
1766  COMPLEX*16 ci
1767  parameter(ci=(0d0,1d0))
1768  COMPLEX*16 f2(6)
1769  COMPLEX*16 v3(*)
1770  REAL*8 p2(0:3)
1771  REAL*8 w2
1772  COMPLEX*16 f1(*)
1773  REAL*8 m2
1774  COMPLEX*16 denom
1775  COMPLEX*16 coup
1776  f2(1) = +f1(1)+v3(1)
1777  f2(2) = +f1(2)+v3(2)
1778  p2(0) = -dble(f2(1))
1779  p2(1) = -dble(f2(2))
1780  p2(2) = -dimag(f2(2))
1781  p2(3) = -dimag(f2(1))
1782  denom = coup/(p2(0)**2-p2(1)**2-p2(2)**2-p2(3)**2 - m2 * (m2
1783  $ -ci* w2))
1784  f2(3)= denom*ci*(f1(3)*(p2(0)*(v3(3)+v3(6))+(p2(1)*(-1d0)*(v3(4)
1785  $ +ci*(v3(5)))+(p2(2)*(+ci*(v3(4))-v3(5))-p2(3)*(v3(3)+v3(6)))))
1786  $ +(f1(4)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*(v3(6)-v3(3))
1787  $ +(p2(2)*(-ci*(v3(6))+ci*(v3(3)))+p2(3)*(+ci*(v3(5))-v3(4)))))
1788  $ +m2*(f1(5)*(v3(3)-v3(6))+f1(6)*(+ci*(v3(5))-v3(4)))))
1789  f2(4)= denom*(-ci)*(f1(3)*(p2(0)*(-1d0)*(v3(4)+ci*(v3(5)))+(p2(1)
1790  $ *(v3(3)+v3(6))+(p2(2)*(+ci*(v3(3)+v3(6)))-p2(3)*(v3(4)
1791  $ +ci*(v3(5))))))+(f1(4)*(p2(0)*(v3(6)-v3(3))+(p2(1)*(v3(4)
1792  $ -ci*(v3(5)))+(p2(2)*(v3(5)+ci*(v3(4)))+p2(3)*(v3(6)-v3(3)))))
1793  $ +m2*(f1(5)*(v3(4)+ci*(v3(5)))-f1(6)*(v3(3)+v3(6)))))
1794  f2(5)= denom*(-ci)*(f1(5)*(p2(0)*(v3(6)-v3(3))+(p2(1)*(v3(4)
1795  $ +ci*(v3(5)))+(p2(2)*(v3(5)-ci*(v3(4)))+p2(3)*(v3(6)-v3(3)))))
1796  $ +(f1(6)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*(-1d0)*(v3(3)+v3(6))
1797  $ +(p2(2)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)-ci*(v3(5))))))
1798  $ +m2*(f1(3)*(-1d0)*(v3(3)+v3(6))+f1(4)*(+ci*(v3(5))-v3(4)))))
1799  f2(6)= denom*ci*(f1(5)*(p2(0)*(-1d0)*(v3(4)+ci*(v3(5)))+(p2(1)
1800  $ *(v3(3)-v3(6))+(p2(2)*(-ci*(v3(6))+ci*(v3(3)))+p2(3)*(v3(4)
1801  $ +ci*(v3(5))))))+(f1(6)*(p2(0)*(v3(3)+v3(6))+(p2(1)*(+ci*(v3(5))
1802  $ -v3(4))+(p2(2)*(-1d0)*(v3(5)+ci*(v3(4)))-p2(3)*(v3(3)+v3(6)))))
1803  $ +m2*(f1(3)*(v3(4)+ci*(v3(5)))+f1(4)*(v3(3)-v3(6)))))
1804  END
1805 
1806 C This File is Automatically generated by ALOHA
1807 C The process calculated in this file is:
1808 C Gamma(3,2,1)
1809 C
1810  SUBROUTINE ffv1p0_3(F1, F2, COUP, M3, W3,V3)
1811  IMPLICIT NONE
1812  COMPLEX*16 ci
1813  parameter(ci=(0d0,1d0))
1814  COMPLEX*16 f2(*)
1815  COMPLEX*16 v3(6)
1816  REAL*8 w3
1817  REAL*8 p3(0:3)
1818  REAL*8 m3
1819  COMPLEX*16 f1(*)
1820  COMPLEX*16 denom
1821  COMPLEX*16 coup
1822  v3(1) = +f1(1)+f2(1)
1823  v3(2) = +f1(2)+f2(2)
1824  p3(0) = -dble(v3(1))
1825  p3(1) = -dble(v3(2))
1826  p3(2) = -dimag(v3(2))
1827  p3(3) = -dimag(v3(1))
1828  denom = coup/(p3(0)**2-p3(1)**2-p3(2)**2-p3(3)**2 - m3 * (m3
1829  $ -ci* w3))
1830  v3(3)= denom*(-ci)*(f2(5)*f1(3)+f2(6)*f1(4)+f2(3)*f1(5)+f2(4)
1831  $ *f1(6))
1832  v3(4)= denom*(-ci)*(f2(4)*f1(5)+f2(3)*f1(6)-f2(6)*f1(3)-f2(5)
1833  $ *f1(4))
1834  v3(5)= denom*(-ci)*(-ci*(f2(6)*f1(3)+f2(3)*f1(6))+ci*(f2(5)*f1(4)
1835  $ +f2(4)*f1(5)))
1836  v3(6)= denom*(-ci)*(f2(6)*f1(4)+f2(3)*f1(5)-f2(5)*f1(3)-f2(4)
1837  $ *f1(6))
1838  END
1839 
1840 C This File is Automatically generated by ALOHA
1841 C The process calculated in this file is:
1842 C Gamma(3,2,-1)*ProjM(-1,1)
1843 C
1844  SUBROUTINE ffv2_0(F1, F2, V3, COUP,VERTEX)
1845  IMPLICIT NONE
1846  COMPLEX*16 ci
1847  parameter(ci=(0d0,1d0))
1848  COMPLEX*16 v3(*)
1849  COMPLEX*16 f1(*)
1850  COMPLEX*16 f2(*)
1851  COMPLEX*16 vertex
1852  COMPLEX*16 coup
1853  COMPLEX*16 tmp13
1854  tmp13 = (f1(3)*(f2(5)*(v3(3)+v3(6))+f2(6)*(v3(4)+ci*(v3(5))))
1855  $ +f1(4)*(f2(5)*(v3(4)-ci*(v3(5)))+f2(6)*(v3(3)-v3(6))))
1856  vertex = coup*(-ci) * tmp13
1857  END
1858 
1859 
1860 C This File is Automatically generated by ALOHA
1861 C The process calculated in this file is:
1862 C Gamma(3,2,-1)*ProjM(-1,1)
1863 C
1864  SUBROUTINE ffv2_5_0(F1, F2, V3, COUP1, COUP2,VERTEX)
1865  IMPLICIT NONE
1866  COMPLEX*16 ci
1867  parameter(ci=(0d0,1d0))
1868  COMPLEX*16 coup2
1869  COMPLEX*16 v3(*)
1870  COMPLEX*16 f1(*)
1871  COMPLEX*16 coup1
1872  COMPLEX*16 f2(*)
1873  COMPLEX*16 vertex
1874  COMPLEX*16 tmp
1875  CALL ffv2_0(f1,f2,v3,coup1,vertex)
1876  CALL ffv5_0(f1,f2,v3,coup2,tmp)
1877  vertex = vertex + tmp
1878  END
1879 
1880 
1881 C This File is Automatically generated by ALOHA
1882 C The process calculated in this file is:
1883 C Gamma(3,2,-1)*ProjM(-1,1)
1884 C
1885  SUBROUTINE ffv2_3_0(F1, F2, V3, COUP1, COUP2,VERTEX)
1886  IMPLICIT NONE
1887  COMPLEX*16 ci
1888  parameter(ci=(0d0,1d0))
1889  COMPLEX*16 coup2
1890  COMPLEX*16 v3(*)
1891  COMPLEX*16 f1(*)
1892  COMPLEX*16 coup1
1893  COMPLEX*16 f2(*)
1894  COMPLEX*16 vertex
1895  COMPLEX*16 tmp
1896  CALL ffv2_0(f1,f2,v3,coup1,vertex)
1897  CALL ffv3_0(f1,f2,v3,coup2,tmp)
1898  vertex = vertex + tmp
1899  END
1900 
1901 
1902 C This File is Automatically generated by ALOHA
1903 C The process calculated in this file is:
1904 C Gamma(3,2,-1)*ProjM(-1,1)
1905 C
1906  SUBROUTINE ffv2_4_0(F1, F2, V3, COUP1, COUP2,VERTEX)
1907  IMPLICIT NONE
1908  COMPLEX*16 ci
1909  parameter(ci=(0d0,1d0))
1910  COMPLEX*16 coup2
1911  COMPLEX*16 v3(*)
1912  COMPLEX*16 f1(*)
1913  COMPLEX*16 coup1
1914  COMPLEX*16 f2(*)
1915  COMPLEX*16 vertex
1916  COMPLEX*16 tmp
1917  CALL ffv2_0(f1,f2,v3,coup1,vertex)
1918  CALL ffv4_0(f1,f2,v3,coup2,tmp)
1919  vertex = vertex + tmp
1920  END
1921 
1922 C This File is Automatically generated by ALOHA
1923 C The process calculated in this file is:
1924 C Gamma(3,2,-1)*ProjM(-1,1)
1925 C
1926  SUBROUTINE ffv2_1(F2, V3, COUP, M1, W1,F1)
1927  IMPLICIT NONE
1928  COMPLEX*16 ci
1929  parameter(ci=(0d0,1d0))
1930  COMPLEX*16 f2(*)
1931  COMPLEX*16 v3(*)
1932  REAL*8 p1(0:3)
1933  REAL*8 m1
1934  REAL*8 w1
1935  COMPLEX*16 f1(6)
1936  COMPLEX*16 denom
1937  COMPLEX*16 coup
1938  f1(1) = +f2(1)+v3(1)
1939  f1(2) = +f2(2)+v3(2)
1940  p1(0) = -dble(f1(1))
1941  p1(1) = -dble(f1(2))
1942  p1(2) = -dimag(f1(2))
1943  p1(3) = -dimag(f1(1))
1944  denom = coup/(p1(0)**2-p1(1)**2-p1(2)**2-p1(3)**2 - m1 * (m1
1945  $ -ci* w1))
1946  f1(3)= denom*ci * m1*(f2(5)*(v3(3)+v3(6))+f2(6)*(v3(4)+ci
1947  $ *(v3(5))))
1948  f1(4)= denom*(-ci) * m1*(f2(5)*(+ci*(v3(5))-v3(4))+f2(6)*(v3(6)
1949  $ -v3(3)))
1950  f1(5)= denom*(-ci)*(f2(5)*(p1(0)*(v3(3)+v3(6))+(p1(1)*(+ci*(v3(5))
1951  $ -v3(4))+(p1(2)*(-1d0)*(v3(5)+ci*(v3(4)))-p1(3)*(v3(3)+v3(6)))))
1952  $ +f2(6)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*(v3(6)-v3(3))+(p1(2)*(
1953  $ -ci*(v3(3))+ci*(v3(6)))-p1(3)*(v3(4)+ci*(v3(5)))))))
1954  f1(6)= denom*(-ci)*(f2(5)*(p1(0)*(v3(4)-ci*(v3(5)))+(p1(1)*
1955  $ (-1d0)*(v3(3)+v3(6))+(p1(2)*(+ci*(v3(3)+v3(6)))+p1(3)*(v3(4)
1956  $ -ci*(v3(5))))))+f2(6)*(p1(0)*(v3(3)-v3(6))+(p1(1)*(-1d0)*(v3(4)
1957  $ +ci*(v3(5)))+(p1(2)*(+ci*(v3(4))-v3(5))+p1(3)*(v3(3)-v3(6))))))
1958  END
1959 
1960 
1961 C This File is Automatically generated by ALOHA
1962 C The process calculated in this file is:
1963 C Gamma(3,2,-1)*ProjM(-1,1)
1964 C
1965  SUBROUTINE ffv2_3_1(F2, V3, COUP1, COUP2, M1, W1,F1)
1966  IMPLICIT NONE
1967  COMPLEX*16 ci
1968  parameter(ci=(0d0,1d0))
1969  COMPLEX*16 f2(*)
1970  COMPLEX*16 v3(*)
1971  REAL*8 p1(0:3)
1972  REAL*8 m1
1973  REAL*8 w1
1974  COMPLEX*16 f1(6)
1975  COMPLEX*16 coup1
1976  COMPLEX*16 denom
1977  COMPLEX*16 coup2
1978  INTEGER*4 i
1979  COMPLEX*16 ftmp(6)
1980  CALL ffv2_1(f2,v3,coup1,m1,w1,f1)
1981  CALL ffv3_1(f2,v3,coup2,m1,w1,ftmp)
1982  DO i = 3, 6
1983  f1(i) = f1(i) + ftmp(i)
1984  ENDDO
1985  END
1986 
1987 
1988 C This File is Automatically generated by ALOHA
1989 C The process calculated in this file is:
1990 C Gamma(3,2,-1)*ProjM(-1,1)
1991 C
1992  SUBROUTINE ffv2_4_1(F2, V3, COUP1, COUP2, M1, W1,F1)
1993  IMPLICIT NONE
1994  COMPLEX*16 ci
1995  parameter(ci=(0d0,1d0))
1996  COMPLEX*16 f2(*)
1997  COMPLEX*16 v3(*)
1998  REAL*8 p1(0:3)
1999  REAL*8 m1
2000  REAL*8 w1
2001  COMPLEX*16 f1(6)
2002  COMPLEX*16 coup1
2003  COMPLEX*16 denom
2004  COMPLEX*16 coup2
2005  INTEGER*4 i
2006  COMPLEX*16 ftmp(6)
2007  CALL ffv2_1(f2,v3,coup1,m1,w1,f1)
2008  CALL ffv4_1(f2,v3,coup2,m1,w1,ftmp)
2009  DO i = 3, 6
2010  f1(i) = f1(i) + ftmp(i)
2011  ENDDO
2012  END
2013 
2014 
2015 C This File is Automatically generated by ALOHA
2016 C The process calculated in this file is:
2017 C Gamma(3,2,-1)*ProjM(-1,1)
2018 C
2019  SUBROUTINE ffv2_2(F1, V3, COUP, M2, W2,F2)
2020  IMPLICIT NONE
2021  COMPLEX*16 ci
2022  parameter(ci=(0d0,1d0))
2023  COMPLEX*16 f2(6)
2024  COMPLEX*16 v3(*)
2025  REAL*8 p2(0:3)
2026  REAL*8 w2
2027  COMPLEX*16 f1(*)
2028  REAL*8 m2
2029  COMPLEX*16 denom
2030  COMPLEX*16 coup
2031  f2(1) = +f1(1)+v3(1)
2032  f2(2) = +f1(2)+v3(2)
2033  p2(0) = -dble(f2(1))
2034  p2(1) = -dble(f2(2))
2035  p2(2) = -dimag(f2(2))
2036  p2(3) = -dimag(f2(1))
2037  denom = coup/(p2(0)**2-p2(1)**2-p2(2)**2-p2(3)**2 - m2 * (m2
2038  $ -ci* w2))
2039  f2(3)= denom*ci*(f1(3)*(p2(0)*(v3(3)+v3(6))+(p2(1)*(-1d0)*(v3(4)
2040  $ +ci*(v3(5)))+(p2(2)*(+ci*(v3(4))-v3(5))-p2(3)*(v3(3)+v3(6)))))
2041  $ +f1(4)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*(v3(6)-v3(3))+(p2(2)*(
2042  $ -ci*(v3(6))+ci*(v3(3)))+p2(3)*(+ci*(v3(5))-v3(4))))))
2043  f2(4)= denom*ci*(f1(3)*(p2(0)*(v3(4)+ci*(v3(5)))+(p2(1)*
2044  $ (-1d0)*(v3(3)+v3(6))+(p2(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)
2045  $ +ci*(v3(5))))))+f1(4)*(p2(0)*(v3(3)-v3(6))+(p2(1)*(+ci*(v3(5))
2046  $ -v3(4))+(p2(2)*(-1d0)*(v3(5)+ci*(v3(4)))+p2(3)*(v3(3)-v3(6))))))
2047  f2(5)= denom*(-ci) * m2*(f1(3)*(-1d0)*(v3(3)+v3(6))+f1(4)*(
2048  $ +ci*(v3(5))-v3(4)))
2049  f2(6)= denom*ci * m2*(f1(3)*(v3(4)+ci*(v3(5)))+f1(4)*(v3(3)
2050  $ -v3(6)))
2051  END
2052 
2053 
2054 C This File is Automatically generated by ALOHA
2055 C The process calculated in this file is:
2056 C Gamma(3,2,-1)*ProjM(-1,1)
2057 C
2058  SUBROUTINE ffv2_5_2(F1, V3, COUP1, COUP2, M2, W2,F2)
2059  IMPLICIT NONE
2060  COMPLEX*16 ci
2061  parameter(ci=(0d0,1d0))
2062  COMPLEX*16 f2(6)
2063  COMPLEX*16 v3(*)
2064  COMPLEX*16 ftmp(6)
2065  COMPLEX*16 coup2
2066  REAL*8 p2(0:3)
2067  REAL*8 w2
2068  COMPLEX*16 f1(*)
2069  REAL*8 m2
2070  COMPLEX*16 denom
2071  COMPLEX*16 coup1
2072  INTEGER*4 i
2073  CALL ffv2_2(f1,v3,coup1,m2,w2,f2)
2074  CALL ffv5_2(f1,v3,coup2,m2,w2,ftmp)
2075  DO i = 3, 6
2076  f2(i) = f2(i) + ftmp(i)
2077  ENDDO
2078  END
2079 
2080 
2081 C This File is Automatically generated by ALOHA
2082 C The process calculated in this file is:
2083 C Gamma(3,2,-1)*ProjM(-1,1)
2084 C
2085  SUBROUTINE ffv2_4_2(F1, V3, COUP1, COUP2, M2, W2,F2)
2086  IMPLICIT NONE
2087  COMPLEX*16 ci
2088  parameter(ci=(0d0,1d0))
2089  COMPLEX*16 f2(6)
2090  COMPLEX*16 v3(*)
2091  COMPLEX*16 ftmp(6)
2092  COMPLEX*16 coup2
2093  REAL*8 p2(0:3)
2094  REAL*8 w2
2095  COMPLEX*16 f1(*)
2096  REAL*8 m2
2097  COMPLEX*16 denom
2098  COMPLEX*16 coup1
2099  INTEGER*4 i
2100  CALL ffv2_2(f1,v3,coup1,m2,w2,f2)
2101  CALL ffv4_2(f1,v3,coup2,m2,w2,ftmp)
2102  DO i = 3, 6
2103  f2(i) = f2(i) + ftmp(i)
2104  ENDDO
2105  END
2106 
2107 
2108 C This File is Automatically generated by ALOHA
2109 C The process calculated in this file is:
2110 C Gamma(3,2,-1)*ProjM(-1,1)
2111 C
2112  SUBROUTINE ffv2_3_2(F1, V3, COUP1, COUP2, M2, W2,F2)
2113  IMPLICIT NONE
2114  COMPLEX*16 ci
2115  parameter(ci=(0d0,1d0))
2116  COMPLEX*16 f2(6)
2117  COMPLEX*16 v3(*)
2118  COMPLEX*16 ftmp(6)
2119  COMPLEX*16 coup2
2120  REAL*8 p2(0:3)
2121  REAL*8 w2
2122  COMPLEX*16 f1(*)
2123  REAL*8 m2
2124  COMPLEX*16 denom
2125  COMPLEX*16 coup1
2126  INTEGER*4 i
2127  CALL ffv2_2(f1,v3,coup1,m2,w2,f2)
2128  CALL ffv3_2(f1,v3,coup2,m2,w2,ftmp)
2129  DO i = 3, 6
2130  f2(i) = f2(i) + ftmp(i)
2131  ENDDO
2132  END
2133 
2134 
2135 C This File is Automatically generated by ALOHA
2136 C The process calculated in this file is:
2137 C Gamma(3,2,-1)*ProjM(-1,1)
2138 C
2139  SUBROUTINE ffv2_3(F1, F2, COUP, M3, W3,V3)
2140  IMPLICIT NONE
2141  COMPLEX*16 ci
2142  parameter(ci=(0d0,1d0))
2143  COMPLEX*16 denom
2144  COMPLEX*16 v3(6)
2145  REAL*8 w3
2146  COMPLEX*16 tmp0
2147  REAL*8 p3(0:3)
2148  REAL*8 m3
2149  COMPLEX*16 f1(*)
2150  COMPLEX*16 f2(*)
2151  REAL*8 om3
2152  COMPLEX*16 coup
2153  om3 = 0d0
2154  IF (m3.NE.0d0) om3=1d0/m3**2
2155  v3(1) = +f1(1)+f2(1)
2156  v3(2) = +f1(2)+f2(2)
2157  p3(0) = -dble(v3(1))
2158  p3(1) = -dble(v3(2))
2159  p3(2) = -dimag(v3(2))
2160  p3(3) = -dimag(v3(1))
2161  tmp0 = (f1(3)*(f2(5)*(p3(0)+p3(3))+f2(6)*(p3(1)+ci*(p3(2))))
2162  $ +f1(4)*(f2(5)*(p3(1)-ci*(p3(2)))+f2(6)*(p3(0)-p3(3))))
2163  denom = coup/(p3(0)**2-p3(1)**2-p3(2)**2-p3(3)**2 - m3 * (m3
2164  $ -ci* w3))
2165  v3(3)= denom*(-ci)*(f2(5)*f1(3)+f2(6)*f1(4)-p3(0)*om3*tmp0)
2166  v3(4)= denom*(-ci)*(-f2(6)*f1(3)-f2(5)*f1(4)-p3(1)*om3*tmp0)
2167  v3(5)= denom*(-ci)*(-ci*(f2(6)*f1(3))+ci*(f2(5)*f1(4))-p3(2)*om3
2168  $ *tmp0)
2169  v3(6)= denom*(-ci)*(f2(6)*f1(4)-f2(5)*f1(3)-p3(3)*om3*tmp0)
2170  END
2171 
2172 
2173 C This File is Automatically generated by ALOHA
2174 C The process calculated in this file is:
2175 C Gamma(3,2,-1)*ProjM(-1,1)
2176 C
2177  SUBROUTINE ffv2_3_3(F1, F2, COUP1, COUP2, M3, W3,V3)
2178  IMPLICIT NONE
2179  COMPLEX*16 ci
2180  parameter(ci=(0d0,1d0))
2181  COMPLEX*16 denom
2182  COMPLEX*16 v3(6)
2183  REAL*8 w3
2184  REAL*8 p3(0:3)
2185  REAL*8 m3
2186  COMPLEX*16 f1(*)
2187  COMPLEX*16 coup1
2188  COMPLEX*16 f2(*)
2189  COMPLEX*16 coup2
2190  REAL*8 om3
2191  INTEGER*4 i
2192  COMPLEX*16 vtmp(6)
2193  CALL ffv2_3(f1,f2,coup1,m3,w3,v3)
2194  CALL ffv3_3(f1,f2,coup2,m3,w3,vtmp)
2195  DO i = 3, 6
2196  v3(i) = v3(i) + vtmp(i)
2197  ENDDO
2198  END
2199 
2200 
2201 C This File is Automatically generated by ALOHA
2202 C The process calculated in this file is:
2203 C Gamma(3,2,-1)*ProjM(-1,1)
2204 C
2205  SUBROUTINE ffv2_5_3(F1, F2, COUP1, COUP2, M3, W3,V3)
2206  IMPLICIT NONE
2207  COMPLEX*16 ci
2208  parameter(ci=(0d0,1d0))
2209  COMPLEX*16 denom
2210  COMPLEX*16 v3(6)
2211  REAL*8 w3
2212  REAL*8 p3(0:3)
2213  REAL*8 m3
2214  COMPLEX*16 f1(*)
2215  COMPLEX*16 coup1
2216  COMPLEX*16 f2(*)
2217  COMPLEX*16 coup2
2218  REAL*8 om3
2219  INTEGER*4 i
2220  COMPLEX*16 vtmp(6)
2221  CALL ffv2_3(f1,f2,coup1,m3,w3,v3)
2222  CALL ffv5_3(f1,f2,coup2,m3,w3,vtmp)
2223  DO i = 3, 6
2224  v3(i) = v3(i) + vtmp(i)
2225  ENDDO
2226  END
2227 
2228 
2229 C This File is Automatically generated by ALOHA
2230 C The process calculated in this file is:
2231 C Gamma(3,2,-1)*ProjM(-1,1)
2232 C
2233  SUBROUTINE ffv2_4_3(F1, F2, COUP1, COUP2, M3, W3,V3)
2234  IMPLICIT NONE
2235  COMPLEX*16 ci
2236  parameter(ci=(0d0,1d0))
2237  COMPLEX*16 denom
2238  COMPLEX*16 v3(6)
2239  REAL*8 w3
2240  REAL*8 p3(0:3)
2241  REAL*8 m3
2242  COMPLEX*16 f1(*)
2243  COMPLEX*16 coup1
2244  COMPLEX*16 f2(*)
2245  COMPLEX*16 coup2
2246  REAL*8 om3
2247  INTEGER*4 i
2248  COMPLEX*16 vtmp(6)
2249  CALL ffv2_3(f1,f2,coup1,m3,w3,v3)
2250  CALL ffv4_3(f1,f2,coup2,m3,w3,vtmp)
2251  DO i = 3, 6
2252  v3(i) = v3(i) + vtmp(i)
2253  ENDDO
2254  END
2255 
2256 
2257 C This File is Automatically generated by ALOHA
2258 C The process calculated in this file is:
2259 C Gamma(3,2,-1)*ProjM(-1,1) - 2*Gamma(3,2,-1)*ProjP(-1,1)
2260 C
2261  SUBROUTINE ffv3_0(F1, F2, V3, COUP,VERTEX)
2262  IMPLICIT NONE
2263  COMPLEX*16 ci
2264  parameter(ci=(0d0,1d0))
2265  COMPLEX*16 v3(*)
2266  COMPLEX*16 f1(*)
2267  COMPLEX*16 f2(*)
2268  COMPLEX*16 tmp14
2269  COMPLEX*16 vertex
2270  COMPLEX*16 coup
2271  COMPLEX*16 tmp13
2272  tmp14 = (f1(5)*(f2(3)*(v3(3)-v3(6))-f2(4)*(v3(4)+ci*(v3(5))))
2273  $ +f1(6)*(f2(3)*(+ci*(v3(5))-v3(4))+f2(4)*(v3(3)+v3(6))))
2274  tmp13 = (f1(3)*(f2(5)*(v3(3)+v3(6))+f2(6)*(v3(4)+ci*(v3(5))))
2275  $ +f1(4)*(f2(5)*(v3(4)-ci*(v3(5)))+f2(6)*(v3(3)-v3(6))))
2276  vertex = coup*(-ci*(tmp13)+2d0 * ci*(tmp14))
2277  END
2278 
2279 
2280 C This File is Automatically generated by ALOHA
2281 C The process calculated in this file is:
2282 C Gamma(3,2,-1)*ProjM(-1,1) - 2*Gamma(3,2,-1)*ProjP(-1,1)
2283 C
2284  SUBROUTINE ffv3_1(F2, V3, COUP, M1, W1,F1)
2285  IMPLICIT NONE
2286  COMPLEX*16 ci
2287  parameter(ci=(0d0,1d0))
2288  COMPLEX*16 f2(*)
2289  COMPLEX*16 v3(*)
2290  REAL*8 p1(0:3)
2291  REAL*8 m1
2292  REAL*8 w1
2293  COMPLEX*16 f1(6)
2294  COMPLEX*16 denom
2295  COMPLEX*16 coup
2296  f1(1) = +f2(1)+v3(1)
2297  f1(2) = +f2(2)+v3(2)
2298  p1(0) = -dble(f1(1))
2299  p1(1) = -dble(f1(2))
2300  p1(2) = -dimag(f1(2))
2301  p1(3) = -dimag(f1(1))
2302  denom = coup/(p1(0)**2-p1(1)**2-p1(2)**2-p1(3)**2 - m1 * (m1
2303  $ -ci* w1))
2304  f1(3)= denom*(-2d0) * ci*(f2(3)*(p1(0)*(v3(6)-v3(3))+(p1(1)*(v3(4)
2305  $ -ci*(v3(5)))+(p1(2)*(v3(5)+ci*(v3(4)))+p1(3)*(v3(6)-v3(3)))))+(
2306  $ +1d0/2d0*(m1*(+2d0*(f2(5)*(-1d0)/2d0*(v3(3)+v3(6)))-f2(6)*(v3(4)
2307  $ +ci*(v3(5)))))+f2(4)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*
2308  $ (-1d0)*(v3(3)+v3(6))+(p1(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p1(3)*(v3(4)
2309  $ +ci*(v3(5))))))))
2310  f1(4)= denom*(-2d0) * ci*(f2(3)*(p1(0)*(v3(4)-ci*(v3(5)))
2311  $ +(p1(1)*(v3(6)-v3(3))+(p1(2)*(-ci*(v3(6))+ci*(v3(3)))+p1(3)*(
2312  $ +ci*(v3(5))-v3(4)))))+(+1d0/2d0*(m1*(f2(6)*(v3(6)-v3(3))
2313  $ +2d0*(f2(5)*1d0/2d0*(+ci*(v3(5))-v3(4)))))+f2(4)*(p1(0)*
2314  $ (-1d0)*(v3(3)+v3(6))+(p1(1)*(v3(4)+ci*(v3(5)))+(p1(2)*(v3(5)
2315  $ -ci*(v3(4)))+p1(3)*(v3(3)+v3(6)))))))
2316  f1(5)= denom*ci*(f2(5)*(p1(0)*(-1d0)*(v3(3)+v3(6))+(p1(1)*(v3(4)
2317  $ -ci*(v3(5)))+(p1(2)*(v3(5)+ci*(v3(4)))+p1(3)*(v3(3)+v3(6)))))
2318  $ +(f2(6)*(p1(0)*(-1d0)*(v3(4)+ci*(v3(5)))+(p1(1)*(v3(3)-v3(6))
2319  $ +(p1(2)*(-ci*(v3(6))+ci*(v3(3)))+p1(3)*(v3(4)+ci*(v3(5))))))
2320  $ +m1*(f2(3)*2d0*(v3(6)-v3(3))+2d0*(f2(4)*(v3(4)+ci*(v3(5)))))))
2321  f1(6)= denom*(-ci)*(f2(5)*(p1(0)*(v3(4)-ci*(v3(5)))+(p1(1)*
2322  $ (-1d0)*(v3(3)+v3(6))+(p1(2)*(+ci*(v3(3)+v3(6)))+p1(3)*(v3(4)
2323  $ -ci*(v3(5))))))+(f2(6)*(p1(0)*(v3(3)-v3(6))+(p1(1)*(-1d0)*(v3(4)
2324  $ +ci*(v3(5)))+(p1(2)*(+ci*(v3(4))-v3(5))+p1(3)*(v3(3)-v3(6)))))
2325  $ +m1*(f2(3)*2d0*(+ci*(v3(5))-v3(4))+2d0*(f2(4)*(v3(3)+v3(6))))))
2326  END
2327 
2328 
2329 C This File is Automatically generated by ALOHA
2330 C The process calculated in this file is:
2331 C Gamma(3,2,-1)*ProjM(-1,1) - 2*Gamma(3,2,-1)*ProjP(-1,1)
2332 C
2333  SUBROUTINE ffv3_2(F1, V3, COUP, M2, W2,F2)
2334  IMPLICIT NONE
2335  COMPLEX*16 ci
2336  parameter(ci=(0d0,1d0))
2337  COMPLEX*16 f2(6)
2338  COMPLEX*16 v3(*)
2339  REAL*8 p2(0:3)
2340  REAL*8 w2
2341  COMPLEX*16 f1(*)
2342  REAL*8 m2
2343  COMPLEX*16 denom
2344  COMPLEX*16 coup
2345  f2(1) = +f1(1)+v3(1)
2346  f2(2) = +f1(2)+v3(2)
2347  p2(0) = -dble(f2(1))
2348  p2(1) = -dble(f2(2))
2349  p2(2) = -dimag(f2(2))
2350  p2(3) = -dimag(f2(1))
2351  denom = coup/(p2(0)**2-p2(1)**2-p2(2)**2-p2(3)**2 - m2 * (m2
2352  $ -ci* w2))
2353  f2(3)= denom*ci*(f1(3)*(p2(0)*(v3(3)+v3(6))+(p2(1)*(-1d0)*(v3(4)
2354  $ +ci*(v3(5)))+(p2(2)*(+ci*(v3(4))-v3(5))-p2(3)*(v3(3)+v3(6)))))
2355  $ +(f1(4)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*(v3(6)-v3(3))
2356  $ +(p2(2)*(-ci*(v3(6))+ci*(v3(3)))+p2(3)*(+ci*(v3(5))-v3(4)))))
2357  $ +m2*(f1(5)*2d0*(v3(6)-v3(3))+2d0*(f1(6)*(v3(4)-ci*(v3(5)))))))
2358  f2(4)= denom*ci*(f1(3)*(p2(0)*(v3(4)+ci*(v3(5)))+(p2(1)*
2359  $ (-1d0)*(v3(3)+v3(6))+(p2(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)
2360  $ +ci*(v3(5))))))+(f1(4)*(p2(0)*(v3(3)-v3(6))+(p2(1)*(+ci*(v3(5))
2361  $ -v3(4))+(p2(2)*(-1d0)*(v3(5)+ci*(v3(4)))+p2(3)*(v3(3)-v3(6)))))
2362  $ +m2*(f1(5)*2d0*(v3(4)+ci*(v3(5)))-2d0*(f1(6)*(v3(3)+v3(6))))))
2363  f2(5)= denom*2d0 * ci*(f1(5)*(p2(0)*(v3(6)-v3(3))+(p2(1)*(v3(4)
2364  $ +ci*(v3(5)))+(p2(2)*(v3(5)-ci*(v3(4)))+p2(3)*(v3(6)-v3(3)))))+(
2365  $ +1d0/2d0*(m2*(f1(4)*(v3(4)-ci*(v3(5)))+2d0*(f1(3)*1d0/2d0
2366  $ *(v3(3)+v3(6)))))+f1(6)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*
2367  $ (-1d0)*(v3(3)+v3(6))+(p2(2)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)
2368  $ -ci*(v3(5))))))))
2369  f2(6)= denom*2d0 * ci*(f1(5)*(p2(0)*(v3(4)+ci*(v3(5)))+(p2(1)
2370  $ *(v3(6)-v3(3))+(p2(2)*(-ci*(v3(3))+ci*(v3(6)))-p2(3)*(v3(4)
2371  $ +ci*(v3(5))))))+(+1d0/2d0*(m2*(f1(4)*(v3(3)-v3(6))+2d0*(f1(3)
2372  $ *1d0/2d0*(v3(4)+ci*(v3(5))))))+f1(6)*(p2(0)*(-1d0)*(v3(3)+v3(6))
2373  $ +(p2(1)*(v3(4)-ci*(v3(5)))+(p2(2)*(v3(5)+ci*(v3(4)))+p2(3)
2374  $ *(v3(3)+v3(6)))))))
2375  END
2376 
2377 
2378 C This File is Automatically generated by ALOHA
2379 C The process calculated in this file is:
2380 C Gamma(3,2,-1)*ProjM(-1,1) - 2*Gamma(3,2,-1)*ProjP(-1,1)
2381 C
2382  SUBROUTINE ffv3_3(F1, F2, COUP, M3, W3,V3)
2383  IMPLICIT NONE
2384  COMPLEX*16 ci
2385  parameter(ci=(0d0,1d0))
2386  COMPLEX*16 denom
2387  COMPLEX*16 v3(6)
2388  COMPLEX*16 tmp1
2389  REAL*8 w3
2390  COMPLEX*16 tmp0
2391  REAL*8 p3(0:3)
2392  REAL*8 m3
2393  COMPLEX*16 f1(*)
2394  COMPLEX*16 f2(*)
2395  REAL*8 om3
2396  COMPLEX*16 coup
2397  om3 = 0d0
2398  IF (m3.NE.0d0) om3=1d0/m3**2
2399  v3(1) = +f1(1)+f2(1)
2400  v3(2) = +f1(2)+f2(2)
2401  p3(0) = -dble(v3(1))
2402  p3(1) = -dble(v3(2))
2403  p3(2) = -dimag(v3(2))
2404  p3(3) = -dimag(v3(1))
2405  tmp1 = (f1(5)*(f2(3)*(p3(0)-p3(3))-f2(4)*(p3(1)+ci*(p3(2))))
2406  $ +f1(6)*(f2(3)*(+ci*(p3(2))-p3(1))+f2(4)*(p3(0)+p3(3))))
2407  tmp0 = (f1(3)*(f2(5)*(p3(0)+p3(3))+f2(6)*(p3(1)+ci*(p3(2))))
2408  $ +f1(4)*(f2(5)*(p3(1)-ci*(p3(2)))+f2(6)*(p3(0)-p3(3))))
2409  denom = coup/(p3(0)**2-p3(1)**2-p3(2)**2-p3(3)**2 - m3 * (m3
2410  $ -ci* w3))
2411  v3(3)= denom*2d0 * ci*(om3*1d0/2d0 * p3(0)*(tmp0-2d0*(tmp1))
2412  $ +(-1d0/2d0*(f2(5)*f1(3)+f2(6)*f1(4))+f2(3)*f1(5)+f2(4)*f1(6)))
2413  v3(4)= denom*2d0 * ci*(om3*1d0/2d0 * p3(1)*(tmp0-2d0*(tmp1))+(
2414  $ +1d0/2d0*(f2(6)*f1(3)+f2(5)*f1(4))+f2(4)*f1(5)+f2(3)*f1(6)))
2415  v3(5)= denom*(-2d0) * ci*(om3*1d0/2d0 * p3(2)*(+2d0*(tmp1)-tmp0)
2416  $ +(-1d0/2d0 * ci*(f2(6)*f1(3))+1d0/2d0 * ci*(f2(5)*f1(4))
2417  $ -ci*(f2(4)*f1(5))+ci*(f2(3)*f1(6))))
2418  v3(6)= denom*(-2d0) * ci*(om3*1d0/2d0 * p3(3)*(+2d0*(tmp1)-tmp0)
2419  $ +(-1d0/2d0*(f2(5)*f1(3))+1d0/2d0*(f2(6)*f1(4))-f2(3)*f1(5)
2420  $ +f2(4)*f1(6)))
2421  END
2422 
2423 
2424 C This File is Automatically generated by ALOHA
2425 C The process calculated in this file is:
2426 C Gamma(3,2,-1)*ProjM(-1,1) + 2*Gamma(3,2,-1)*ProjP(-1,1)
2427 C
2428  SUBROUTINE ffv4_0(F1, F2, V3, COUP,VERTEX)
2429  IMPLICIT NONE
2430  COMPLEX*16 ci
2431  parameter(ci=(0d0,1d0))
2432  COMPLEX*16 v3(*)
2433  COMPLEX*16 f1(*)
2434  COMPLEX*16 f2(*)
2435  COMPLEX*16 tmp14
2436  COMPLEX*16 vertex
2437  COMPLEX*16 coup
2438  COMPLEX*16 tmp13
2439  tmp14 = (f1(5)*(f2(3)*(v3(3)-v3(6))-f2(4)*(v3(4)+ci*(v3(5))))
2440  $ +f1(6)*(f2(3)*(+ci*(v3(5))-v3(4))+f2(4)*(v3(3)+v3(6))))
2441  tmp13 = (f1(3)*(f2(5)*(v3(3)+v3(6))+f2(6)*(v3(4)+ci*(v3(5))))
2442  $ +f1(4)*(f2(5)*(v3(4)-ci*(v3(5)))+f2(6)*(v3(3)-v3(6))))
2443  vertex = coup*(-1d0)*(+ci*(tmp13)+2d0 * ci*(tmp14))
2444  END
2445 
2446 
2447 C This File is Automatically generated by ALOHA
2448 C The process calculated in this file is:
2449 C Gamma(3,2,-1)*ProjM(-1,1) + 2*Gamma(3,2,-1)*ProjP(-1,1)
2450 C
2451  SUBROUTINE ffv4_1(F2, V3, COUP, M1, W1,F1)
2452  IMPLICIT NONE
2453  COMPLEX*16 ci
2454  parameter(ci=(0d0,1d0))
2455  COMPLEX*16 f2(*)
2456  COMPLEX*16 v3(*)
2457  REAL*8 p1(0:3)
2458  REAL*8 m1
2459  REAL*8 w1
2460  COMPLEX*16 f1(6)
2461  COMPLEX*16 denom
2462  COMPLEX*16 coup
2463  f1(1) = +f2(1)+v3(1)
2464  f1(2) = +f2(2)+v3(2)
2465  p1(0) = -dble(f1(1))
2466  p1(1) = -dble(f1(2))
2467  p1(2) = -dimag(f1(2))
2468  p1(3) = -dimag(f1(1))
2469  denom = coup/(p1(0)**2-p1(1)**2-p1(2)**2-p1(3)**2 - m1 * (m1
2470  $ -ci* w1))
2471  f1(3)= denom*2d0 * ci*(f2(3)*(p1(0)*(v3(6)-v3(3))+(p1(1)*(v3(4)
2472  $ -ci*(v3(5)))+(p1(2)*(v3(5)+ci*(v3(4)))+p1(3)*(v3(6)-v3(3)))))+(
2473  $ +1d0/2d0*(m1*(f2(6)*(v3(4)+ci*(v3(5)))+2d0*(f2(5)*1d0/2d0
2474  $ *(v3(3)+v3(6)))))+f2(4)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*
2475  $ (-1d0)*(v3(3)+v3(6))+(p1(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p1(3)*(v3(4)
2476  $ +ci*(v3(5))))))))
2477  f1(4)= denom*2d0 * ci*(f2(3)*(p1(0)*(v3(4)-ci*(v3(5)))+(p1(1)
2478  $ *(v3(6)-v3(3))+(p1(2)*(-ci*(v3(6))+ci*(v3(3)))+p1(3)*(
2479  $ +ci*(v3(5))-v3(4)))))+(+1d0/2d0*(m1*(f2(6)*(v3(3)-v3(6))
2480  $ +2d0*(f2(5)*1d0/2d0*(v3(4)-ci*(v3(5))))))+f2(4)*(p1(0)*
2481  $ (-1d0)*(v3(3)+v3(6))+(p1(1)*(v3(4)+ci*(v3(5)))+(p1(2)*(v3(5)
2482  $ -ci*(v3(4)))+p1(3)*(v3(3)+v3(6)))))))
2483  f1(5)= denom*(-ci)*(f2(5)*(p1(0)*(v3(3)+v3(6))+(p1(1)*(+ci*(v3(5))
2484  $ -v3(4))+(p1(2)*(-1d0)*(v3(5)+ci*(v3(4)))-p1(3)*(v3(3)+v3(6)))))
2485  $ +(f2(6)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*(v3(6)-v3(3))
2486  $ +(p1(2)*(-ci*(v3(3))+ci*(v3(6)))-p1(3)*(v3(4)+ci*(v3(5))))))
2487  $ +m1*(f2(3)*2d0*(v3(6)-v3(3))+2d0*(f2(4)*(v3(4)+ci*(v3(5)))))))
2488  f1(6)= denom*ci*(f2(5)*(p1(0)*(+ci*(v3(5))-v3(4))+(p1(1)*(v3(3)
2489  $ +v3(6))+(p1(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p1(3)*(+ci*(v3(5))
2490  $ -v3(4)))))+(f2(6)*(p1(0)*(v3(6)-v3(3))+(p1(1)*(v3(4)+ci
2491  $ *(v3(5)))+(p1(2)*(v3(5)-ci*(v3(4)))+p1(3)*(v3(6)-v3(3)))))
2492  $ +m1*(f2(3)*2d0*(+ci*(v3(5))-v3(4))+2d0*(f2(4)*(v3(3)+v3(6))))))
2493  END
2494 
2495 
2496 C This File is Automatically generated by ALOHA
2497 C The process calculated in this file is:
2498 C Gamma(3,2,-1)*ProjM(-1,1) + 2*Gamma(3,2,-1)*ProjP(-1,1)
2499 C
2500  SUBROUTINE ffv4_2(F1, V3, COUP, M2, W2,F2)
2501  IMPLICIT NONE
2502  COMPLEX*16 ci
2503  parameter(ci=(0d0,1d0))
2504  COMPLEX*16 f2(6)
2505  COMPLEX*16 v3(*)
2506  REAL*8 p2(0:3)
2507  REAL*8 w2
2508  COMPLEX*16 f1(*)
2509  REAL*8 m2
2510  COMPLEX*16 denom
2511  COMPLEX*16 coup
2512  f2(1) = +f1(1)+v3(1)
2513  f2(2) = +f1(2)+v3(2)
2514  p2(0) = -dble(f2(1))
2515  p2(1) = -dble(f2(2))
2516  p2(2) = -dimag(f2(2))
2517  p2(3) = -dimag(f2(1))
2518  denom = coup/(p2(0)**2-p2(1)**2-p2(2)**2-p2(3)**2 - m2 * (m2
2519  $ -ci* w2))
2520  f2(3)= denom*ci*(f1(3)*(p2(0)*(v3(3)+v3(6))+(p2(1)*(-1d0)*(v3(4)
2521  $ +ci*(v3(5)))+(p2(2)*(+ci*(v3(4))-v3(5))-p2(3)*(v3(3)+v3(6)))))
2522  $ +(f1(4)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*(v3(6)-v3(3))
2523  $ +(p2(2)*(-ci*(v3(6))+ci*(v3(3)))+p2(3)*(+ci*(v3(5))-v3(4)))))
2524  $ +m2*(f1(5)*2d0*(v3(3)-v3(6))+2d0*(f1(6)*(+ci*(v3(5))-v3(4))))))
2525  f2(4)= denom*ci*(f1(3)*(p2(0)*(v3(4)+ci*(v3(5)))+(p2(1)*
2526  $ (-1d0)*(v3(3)+v3(6))+(p2(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)
2527  $ +ci*(v3(5))))))+(f1(4)*(p2(0)*(v3(3)-v3(6))+(p2(1)*(+ci*(v3(5))
2528  $ -v3(4))+(p2(2)*(-1d0)*(v3(5)+ci*(v3(4)))+p2(3)*(v3(3)-v3(6)))))
2529  $ +m2*(f1(5)*(-2d0)*(v3(4)+ci*(v3(5)))+2d0*(f1(6)*(v3(3)+v3(6))))))
2530  f2(5)= denom*(-2d0) * ci*(f1(5)*(p2(0)*(v3(6)-v3(3))+(p2(1)*(v3(4)
2531  $ +ci*(v3(5)))+(p2(2)*(v3(5)-ci*(v3(4)))+p2(3)*(v3(6)-v3(3)))))+(
2532  $ +1d0/2d0*(m2*(f1(4)*(+ci*(v3(5))-v3(4))+2d0*(f1(3)*(-1d0)/2d0
2533  $ *(v3(3)+v3(6)))))+f1(6)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*
2534  $ (-1d0)*(v3(3)+v3(6))+(p2(2)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)
2535  $ -ci*(v3(5))))))))
2536  f2(6)= denom*(-2d0) * ci*(f1(5)*(p2(0)*(v3(4)+ci*(v3(5)))
2537  $ +(p2(1)*(v3(6)-v3(3))+(p2(2)*(-ci*(v3(3))+ci*(v3(6)))-p2(3)
2538  $ *(v3(4)+ci*(v3(5))))))+(+1d0/2d0*(m2*(f1(4)*(v3(6)-v3(3))
2539  $ +2d0*(f1(3)*(-1d0)/2d0*(v3(4)+ci*(v3(5))))))+f1(6)*(p2(0)*
2540  $ (-1d0)*(v3(3)+v3(6))+(p2(1)*(v3(4)-ci*(v3(5)))+(p2(2)*(v3(5)
2541  $ +ci*(v3(4)))+p2(3)*(v3(3)+v3(6)))))))
2542  END
2543 
2544 
2545 C This File is Automatically generated by ALOHA
2546 C The process calculated in this file is:
2547 C Gamma(3,2,-1)*ProjM(-1,1) + 2*Gamma(3,2,-1)*ProjP(-1,1)
2548 C
2549  SUBROUTINE ffv4_3(F1, F2, COUP, M3, W3,V3)
2550  IMPLICIT NONE
2551  COMPLEX*16 ci
2552  parameter(ci=(0d0,1d0))
2553  COMPLEX*16 denom
2554  COMPLEX*16 v3(6)
2555  COMPLEX*16 tmp1
2556  REAL*8 w3
2557  COMPLEX*16 tmp0
2558  REAL*8 p3(0:3)
2559  REAL*8 m3
2560  COMPLEX*16 f1(*)
2561  COMPLEX*16 f2(*)
2562  REAL*8 om3
2563  COMPLEX*16 coup
2564  om3 = 0d0
2565  IF (m3.NE.0d0) om3=1d0/m3**2
2566  v3(1) = +f1(1)+f2(1)
2567  v3(2) = +f1(2)+f2(2)
2568  p3(0) = -dble(v3(1))
2569  p3(1) = -dble(v3(2))
2570  p3(2) = -dimag(v3(2))
2571  p3(3) = -dimag(v3(1))
2572  tmp1 = (f1(5)*(f2(3)*(p3(0)-p3(3))-f2(4)*(p3(1)+ci*(p3(2))))
2573  $ +f1(6)*(f2(3)*(+ci*(p3(2))-p3(1))+f2(4)*(p3(0)+p3(3))))
2574  tmp0 = (f1(3)*(f2(5)*(p3(0)+p3(3))+f2(6)*(p3(1)+ci*(p3(2))))
2575  $ +f1(4)*(f2(5)*(p3(1)-ci*(p3(2)))+f2(6)*(p3(0)-p3(3))))
2576  denom = coup/(p3(0)**2-p3(1)**2-p3(2)**2-p3(3)**2 - m3 * (m3
2577  $ -ci* w3))
2578  v3(3)= denom*(-2d0) * ci*(om3*(-1d0)/2d0 * p3(0)*(tmp0+2d0*(tmp1))+(
2579  $ +1d0/2d0*(f2(5)*f1(3)+f2(6)*f1(4))+f2(3)*f1(5)+f2(4)*f1(6)))
2580  v3(4)= denom*(-2d0) * ci*(om3*(-1d0)/2d0 * p3(1)*(tmp0+2d0*(tmp1))
2581  $ +(-1d0/2d0*(f2(6)*f1(3)+f2(5)*f1(4))+f2(4)*f1(5)+f2(3)*f1(6)))
2582  v3(5)= denom*2d0 * ci*(om3*1d0/2d0 * p3(2)*(tmp0+2d0*(tmp1))+(
2583  $ +1d0/2d0 * ci*(f2(6)*f1(3))-1d0/2d0 * ci*(f2(5)*f1(4))
2584  $ -ci*(f2(4)*f1(5))+ci*(f2(3)*f1(6))))
2585  v3(6)= denom*2d0 * ci*(om3*1d0/2d0 * p3(3)*(tmp0+2d0*(tmp1))+(
2586  $ +1d0/2d0*(f2(5)*f1(3))-1d0/2d0*(f2(6)*f1(4))-f2(3)*f1(5)
2587  $ +f2(4)*f1(6)))
2588  END
2589 
2590 
2591 C This File is Automatically generated by ALOHA
2592 C The process calculated in this file is:
2593 C Gamma(3,2,-1)*ProjM(-1,1) + 4*Gamma(3,2,-1)*ProjP(-1,1)
2594 C
2595  SUBROUTINE ffv5_0(F1, F2, V3, COUP,VERTEX)
2596  IMPLICIT NONE
2597  COMPLEX*16 ci
2598  parameter(ci=(0d0,1d0))
2599  COMPLEX*16 v3(*)
2600  COMPLEX*16 f1(*)
2601  COMPLEX*16 f2(*)
2602  COMPLEX*16 tmp14
2603  COMPLEX*16 vertex
2604  COMPLEX*16 coup
2605  COMPLEX*16 tmp13
2606  tmp14 = (f1(5)*(f2(3)*(v3(3)-v3(6))-f2(4)*(v3(4)+ci*(v3(5))))
2607  $ +f1(6)*(f2(3)*(+ci*(v3(5))-v3(4))+f2(4)*(v3(3)+v3(6))))
2608  tmp13 = (f1(3)*(f2(5)*(v3(3)+v3(6))+f2(6)*(v3(4)+ci*(v3(5))))
2609  $ +f1(4)*(f2(5)*(v3(4)-ci*(v3(5)))+f2(6)*(v3(3)-v3(6))))
2610  vertex = coup*(-1d0)*(+ci*(tmp13)+4d0 * ci*(tmp14))
2611  END
2612 
2613 C This File is Automatically generated by ALOHA
2614 C The process calculated in this file is:
2615 C Gamma(3,2,-1)*ProjM(-1,1) + 4*Gamma(3,2,-1)*ProjP(-1,1)
2616 C
2617  SUBROUTINE ffv5_2(F1, V3, COUP, M2, W2,F2)
2618  IMPLICIT NONE
2619  COMPLEX*16 ci
2620  parameter(ci=(0d0,1d0))
2621  COMPLEX*16 f2(6)
2622  COMPLEX*16 v3(*)
2623  REAL*8 p2(0:3)
2624  REAL*8 w2
2625  COMPLEX*16 f1(*)
2626  REAL*8 m2
2627  COMPLEX*16 denom
2628  COMPLEX*16 coup
2629  f2(1) = +f1(1)+v3(1)
2630  f2(2) = +f1(2)+v3(2)
2631  p2(0) = -dble(f2(1))
2632  p2(1) = -dble(f2(2))
2633  p2(2) = -dimag(f2(2))
2634  p2(3) = -dimag(f2(1))
2635  denom = coup/(p2(0)**2-p2(1)**2-p2(2)**2-p2(3)**2 - m2 * (m2
2636  $ -ci* w2))
2637  f2(3)= denom*ci*(f1(3)*(p2(0)*(v3(3)+v3(6))+(p2(1)*(-1d0)*(v3(4)
2638  $ +ci*(v3(5)))+(p2(2)*(+ci*(v3(4))-v3(5))-p2(3)*(v3(3)+v3(6)))))
2639  $ +(f1(4)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*(v3(6)-v3(3))
2640  $ +(p2(2)*(-ci*(v3(6))+ci*(v3(3)))+p2(3)*(+ci*(v3(5))-v3(4)))))
2641  $ +m2*(f1(5)*4d0*(v3(3)-v3(6))+4d0*(f1(6)*(+ci*(v3(5))-v3(4))))))
2642  f2(4)= denom*ci*(f1(3)*(p2(0)*(v3(4)+ci*(v3(5)))+(p2(1)*
2643  $ (-1d0)*(v3(3)+v3(6))+(p2(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)
2644  $ +ci*(v3(5))))))+(f1(4)*(p2(0)*(v3(3)-v3(6))+(p2(1)*(+ci*(v3(5))
2645  $ -v3(4))+(p2(2)*(-1d0)*(v3(5)+ci*(v3(4)))+p2(3)*(v3(3)-v3(6)))))
2646  $ +m2*(f1(5)*(-4d0)*(v3(4)+ci*(v3(5)))+4d0*(f1(6)*(v3(3)+v3(6))))))
2647  f2(5)= denom*(-4d0) * ci*(f1(5)*(p2(0)*(v3(6)-v3(3))+(p2(1)*(v3(4)
2648  $ +ci*(v3(5)))+(p2(2)*(v3(5)-ci*(v3(4)))+p2(3)*(v3(6)-v3(3)))))+(
2649  $ +1d0/4d0*(m2*(f1(4)*(+ci*(v3(5))-v3(4))+4d0*(f1(3)*(-1d0)/4d0
2650  $ *(v3(3)+v3(6)))))+f1(6)*(p2(0)*(v3(4)-ci*(v3(5)))+(p2(1)*
2651  $ (-1d0)*(v3(3)+v3(6))+(p2(2)*(+ci*(v3(3)+v3(6)))+p2(3)*(v3(4)
2652  $ -ci*(v3(5))))))))
2653  f2(6)= denom*(-4d0) * ci*(f1(5)*(p2(0)*(v3(4)+ci*(v3(5)))
2654  $ +(p2(1)*(v3(6)-v3(3))+(p2(2)*(-ci*(v3(3))+ci*(v3(6)))-p2(3)
2655  $ *(v3(4)+ci*(v3(5))))))+(+1d0/4d0*(m2*(f1(4)*(v3(6)-v3(3))
2656  $ +4d0*(f1(3)*(-1d0)/4d0*(v3(4)+ci*(v3(5))))))+f1(6)*(p2(0)*
2657  $ (-1d0)*(v3(3)+v3(6))+(p2(1)*(v3(4)-ci*(v3(5)))+(p2(2)*(v3(5)
2658  $ +ci*(v3(4)))+p2(3)*(v3(3)+v3(6)))))))
2659  END
2660 
2661 
2662 C This File is Automatically generated by ALOHA
2663 C The process calculated in this file is:
2664 C Gamma(3,2,-1)*ProjM(-1,1) + 4*Gamma(3,2,-1)*ProjP(-1,1)
2665 C
2666  SUBROUTINE ffv5_3(F1, F2, COUP, M3, W3,V3)
2667  IMPLICIT NONE
2668  COMPLEX*16 ci
2669  parameter(ci=(0d0,1d0))
2670  COMPLEX*16 denom
2671  COMPLEX*16 v3(6)
2672  COMPLEX*16 tmp1
2673  REAL*8 w3
2674  COMPLEX*16 tmp0
2675  REAL*8 p3(0:3)
2676  REAL*8 m3
2677  COMPLEX*16 f1(*)
2678  COMPLEX*16 f2(*)
2679  REAL*8 om3
2680  COMPLEX*16 coup
2681  om3 = 0d0
2682  IF (m3.NE.0d0) om3=1d0/m3**2
2683  v3(1) = +f1(1)+f2(1)
2684  v3(2) = +f1(2)+f2(2)
2685  p3(0) = -dble(v3(1))
2686  p3(1) = -dble(v3(2))
2687  p3(2) = -dimag(v3(2))
2688  p3(3) = -dimag(v3(1))
2689  tmp1 = (f1(5)*(f2(3)*(p3(0)-p3(3))-f2(4)*(p3(1)+ci*(p3(2))))
2690  $ +f1(6)*(f2(3)*(+ci*(p3(2))-p3(1))+f2(4)*(p3(0)+p3(3))))
2691  tmp0 = (f1(3)*(f2(5)*(p3(0)+p3(3))+f2(6)*(p3(1)+ci*(p3(2))))
2692  $ +f1(4)*(f2(5)*(p3(1)-ci*(p3(2)))+f2(6)*(p3(0)-p3(3))))
2693  denom = coup/(p3(0)**2-p3(1)**2-p3(2)**2-p3(3)**2 - m3 * (m3
2694  $ -ci* w3))
2695  v3(3)= denom*(-4d0) * ci*(om3*(-1d0)/4d0 * p3(0)*(tmp0+4d0*(tmp1))+(
2696  $ +1d0/4d0*(f2(5)*f1(3)+f2(6)*f1(4))+f2(3)*f1(5)+f2(4)*f1(6)))
2697  v3(4)= denom*(-4d0) * ci*(om3*(-1d0)/4d0 * p3(1)*(tmp0+4d0*(tmp1))
2698  $ +(-1d0/4d0*(f2(6)*f1(3)+f2(5)*f1(4))+f2(4)*f1(5)+f2(3)*f1(6)))
2699  v3(5)= denom*4d0 * ci*(om3*1d0/4d0 * p3(2)*(tmp0+4d0*(tmp1))+(
2700  $ +1d0/4d0 * ci*(f2(6)*f1(3))-1d0/4d0 * ci*(f2(5)*f1(4))
2701  $ -ci*(f2(4)*f1(5))+ci*(f2(3)*f1(6))))
2702  v3(6)= denom*4d0 * ci*(om3*1d0/4d0 * p3(3)*(tmp0+4d0*(tmp1))+(
2703  $ +1d0/4d0*(f2(5)*f1(3))-1d0/4d0*(f2(6)*f1(4))-f2(3)*f1(5)
2704  $ +f2(4)*f1(6)))
2705  END
2706 
2707 
2708 C This File is Automatically generated by ALOHA
2709 C The process calculated in this file is:
2710 C P(3,1)*Metric(1,2) - P(3,2)*Metric(1,2) - P(2,1)*Metric(1,3) +
2711 C P(2,3)*Metric(1,3) + P(1,2)*Metric(2,3) - P(1,3)*Metric(2,3)
2712 C
2713  SUBROUTINE vvv1_0(V1, V2, V3, COUP,VERTEX)
2714  IMPLICIT NONE
2715  COMPLEX*16 ci
2716  parameter(ci=(0d0,1d0))
2717  COMPLEX*16 v2(*)
2718  COMPLEX*16 tmp2
2719  COMPLEX*16 tmp12
2720  COMPLEX*16 v3(*)
2721  COMPLEX*16 tmp11
2722  REAL*8 p1(0:3)
2723  COMPLEX*16 tmp10
2724  REAL*8 p2(0:3)
2725  COMPLEX*16 tmp7
2726  REAL*8 p3(0:3)
2727  COMPLEX*16 tmp6
2728  COMPLEX*16 tmp5
2729  COMPLEX*16 vertex
2730  COMPLEX*16 coup
2731  COMPLEX*16 tmp9
2732  COMPLEX*16 v1(*)
2733  COMPLEX*16 tmp8
2734  p1(0) = dble(v1(1))
2735  p1(1) = dble(v1(2))
2736  p1(2) = dimag(v1(2))
2737  p1(3) = dimag(v1(1))
2738  p2(0) = dble(v2(1))
2739  p2(1) = dble(v2(2))
2740  p2(2) = dimag(v2(2))
2741  p2(3) = dimag(v2(1))
2742  p3(0) = dble(v3(1))
2743  p3(1) = dble(v3(2))
2744  p3(2) = dimag(v3(2))
2745  p3(3) = dimag(v3(1))
2746  tmp9 = (p3(0)*v2(3)-p3(1)*v2(4)-p3(2)*v2(5)-p3(3)*v2(6))
2747  tmp8 = (v2(3)*p1(0)-v2(4)*p1(1)-v2(5)*p1(2)-v2(6)*p1(3))
2748  tmp2 = (v2(3)*v1(3)-v2(4)*v1(4)-v2(5)*v1(5)-v2(6)*v1(6))
2749  tmp5 = (v3(3)*p1(0)-v3(4)*p1(1)-v3(5)*p1(2)-v3(6)*p1(3))
2750  tmp7 = (v1(3)*v3(3)-v1(4)*v3(4)-v1(5)*v3(5)-v1(6)*v3(6))
2751  tmp6 = (v3(3)*p2(0)-v3(4)*p2(1)-v3(5)*p2(2)-v3(6)*p2(3))
2752  tmp11 = (v1(3)*p2(0)-v1(4)*p2(1)-v1(5)*p2(2)-v1(6)*p2(3))
2753  tmp10 = (v2(3)*v3(3)-v2(4)*v3(4)-v2(5)*v3(5)-v2(6)*v3(6))
2754  tmp12 = (p3(0)*v1(3)-p3(1)*v1(4)-p3(2)*v1(5)-p3(3)*v1(6))
2755  vertex = coup*(tmp10*(-ci*(tmp11)+ci*(tmp12))+(tmp2*(-ci*(tmp5)
2756  $ +ci*(tmp6))+tmp7*(-ci*(tmp9)+ci*(tmp8))))
2757  END
2758 
2759 
2760 C This File is Automatically generated by ALOHA
2761 C The process calculated in this file is:
2762 C ProjM(2,1) + ProjP(2,1)
2763 C
2764  SUBROUTINE ffs4_3(F1, F2, COUP, M3, W3,S3)
2765  IMPLICIT NONE
2766  COMPLEX*16 ci
2767  parameter(ci=(0d0,1d0))
2768  COMPLEX*16 denom
2769  COMPLEX*16 s3(3)
2770  REAL*8 w3
2771  REAL*8 p3(0:3)
2772  REAL*8 m3
2773  COMPLEX*16 f1(*)
2774  COMPLEX*16 f2(*)
2775  COMPLEX*16 tmp4
2776  COMPLEX*16 coup
2777  COMPLEX*16 tmp3
2778  s3(1) = +f1(1)+f2(1)
2779  s3(2) = +f1(2)+f2(2)
2780  p3(0) = -dble(s3(1))
2781  p3(1) = -dble(s3(2))
2782  p3(2) = -dimag(s3(2))
2783  p3(3) = -dimag(s3(1))
2784  tmp4 = (f2(5)*f1(5)+f2(6)*f1(6))
2785  tmp3 = (f2(3)*f1(3)+f2(4)*f1(4))
2786  denom = coup/(p3(0)**2-p3(1)**2-p3(2)**2-p3(3)**2 - m3 * (m3
2787  $ -ci* w3))
2788  s3(3)= denom*(+ci*(tmp3+tmp4))
2789  END
2790 
2791 
2792 C This File is Automatically generated by ALOHA
2793 C The process calculated in this file is:
2794 C Metric(1,2)
2795 C
2796  SUBROUTINE vvs1_0(V1, V2, S3, COUP,VERTEX)
2797  IMPLICIT NONE
2798  COMPLEX*16 ci
2799  parameter(ci=(0d0,1d0))
2800  COMPLEX*16 v2(*)
2801  COMPLEX*16 tmp2
2802  COMPLEX*16 s3(*)
2803  COMPLEX*16 vertex
2804  COMPLEX*16 coup
2805  COMPLEX*16 v1(*)
2806  tmp2 = (v2(3)*v1(3)-v2(4)*v1(4)-v2(5)*v1(5)-v2(6)*v1(6))
2807  vertex = coup*(-ci) * tmp2*s3(3)
2808  END
2809 
2810 C This File is Automatically generated by ALOHA
2811 C The process calculated in this file is:
2812 C Gamma(3,2,-1)*ProjM(-1,1)
2813 C
2814  SUBROUTINE ffv2_5_1(F2, V3, COUP1, COUP2, M1, W1,F1)
2815  IMPLICIT NONE
2816  COMPLEX*16 ci
2817  parameter(ci=(0d0,1d0))
2818  COMPLEX*16 f2(*)
2819  COMPLEX*16 v3(*)
2820  REAL*8 p1(0:3)
2821  REAL*8 m1
2822  REAL*8 w1
2823  COMPLEX*16 f1(6)
2824  COMPLEX*16 coup1
2825  COMPLEX*16 denom
2826  COMPLEX*16 coup2
2827  INTEGER*4 i
2828  COMPLEX*16 ftmp(6)
2829  CALL ffv2_1(f2,v3,coup1,m1,w1,f1)
2830  CALL ffv5_1(f2,v3,coup2,m1,w1,ftmp)
2831  DO i = 3, 6
2832  f1(i) = f1(i) + ftmp(i)
2833  ENDDO
2834  END
2835 
2836 C This File is Automatically generated by ALOHA
2837 C The process calculated in this file is:
2838 C Gamma(3,2,-1)*ProjM(-1,1) + 4*Gamma(3,2,-1)*ProjP(-1,1)
2839 C
2840  SUBROUTINE ffv5_1(F2, V3, COUP, M1, W1,F1)
2841  IMPLICIT NONE
2842  COMPLEX*16 ci
2843  parameter(ci=(0d0,1d0))
2844  COMPLEX*16 f2(*)
2845  COMPLEX*16 v3(*)
2846  REAL*8 p1(0:3)
2847  REAL*8 m1
2848  REAL*8 w1
2849  COMPLEX*16 f1(6)
2850  COMPLEX*16 denom
2851  COMPLEX*16 coup
2852  f1(1) = +f2(1)+v3(1)
2853  f1(2) = +f2(2)+v3(2)
2854  p1(0) = -dble(f1(1))
2855  p1(1) = -dble(f1(2))
2856  p1(2) = -dimag(f1(2))
2857  p1(3) = -dimag(f1(1))
2858  denom = coup/(p1(0)**2-p1(1)**2-p1(2)**2-p1(3)**2 - m1 * (m1
2859  $ -ci* w1))
2860  f1(3)= denom*4d0 * ci*(f2(3)*(p1(0)*(v3(6)-v3(3))+(p1(1)*(v3(4)
2861  $ -ci*(v3(5)))+(p1(2)*(v3(5)+ci*(v3(4)))+p1(3)*(v3(6)-v3(3)))))+(
2862  $ +1d0/4d0*(m1*(f2(6)*(v3(4)+ci*(v3(5)))+4d0*(f2(5)*1d0/4d0
2863  $ *(v3(3)+v3(6)))))+f2(4)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*
2864  $ (-1d0)*(v3(3)+v3(6))+(p1(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p1(3)*(v3(4)
2865  $ +ci*(v3(5))))))))
2866  f1(4)= denom*4d0 * ci*(f2(3)*(p1(0)*(v3(4)-ci*(v3(5)))+(p1(1)
2867  $ *(v3(6)-v3(3))+(p1(2)*(-ci*(v3(6))+ci*(v3(3)))+p1(3)*(
2868  $ +ci*(v3(5))-v3(4)))))+(+1d0/4d0*(m1*(f2(6)*(v3(3)-v3(6))
2869  $ +4d0*(f2(5)*1d0/4d0*(v3(4)-ci*(v3(5))))))+f2(4)*(p1(0)*
2870  $ (-1d0)*(v3(3)+v3(6))+(p1(1)*(v3(4)+ci*(v3(5)))+(p1(2)*(v3(5)
2871  $ -ci*(v3(4)))+p1(3)*(v3(3)+v3(6)))))))
2872  f1(5)= denom*(-ci)*(f2(5)*(p1(0)*(v3(3)+v3(6))+(p1(1)*(+ci*(v3(5))
2873  $ -v3(4))+(p1(2)*(-1d0)*(v3(5)+ci*(v3(4)))-p1(3)*(v3(3)+v3(6)))))
2874  $ +(f2(6)*(p1(0)*(v3(4)+ci*(v3(5)))+(p1(1)*(v3(6)-v3(3))
2875  $ +(p1(2)*(-ci*(v3(3))+ci*(v3(6)))-p1(3)*(v3(4)+ci*(v3(5))))))
2876  $ +m1*(f2(3)*4d0*(v3(6)-v3(3))+4d0*(f2(4)*(v3(4)+ci*(v3(5)))))))
2877  f1(6)= denom*ci*(f2(5)*(p1(0)*(+ci*(v3(5))-v3(4))+(p1(1)*(v3(3)
2878  $ +v3(6))+(p1(2)*(-1d0)*(+ci*(v3(3)+v3(6)))+p1(3)*(+ci*(v3(5))
2879  $ -v3(4)))))+(f2(6)*(p1(0)*(v3(6)-v3(3))+(p1(1)*(v3(4)+ci
2880  $ *(v3(5)))+(p1(2)*(v3(5)-ci*(v3(4)))+p1(3)*(v3(6)-v3(3)))))
2881  $ +m1*(f2(3)*4d0*(+ci*(v3(5))-v3(4))+4d0*(f2(4)*(v3(3)+v3(6))))))
2882  END
2883 
2884 C This File is Automatically generated by ALOHA
2885 C The process calculated in this file is:
2886 C P(3,1)*Metric(1,2) - P(3,2)*Metric(1,2) - P(2,1)*Metric(1,3) +
2887 C P(2,3)*Metric(1,3) + P(1,2)*Metric(2,3) - P(1,3)*Metric(2,3)
2888 C
2889  SUBROUTINE vvv1p0_1(V2, V3, COUP, M1, W1,V1)
2890  IMPLICIT NONE
2891  COMPLEX*16 ci
2892  parameter(ci=(0d0,1d0))
2893  COMPLEX*16 v2(*)
2894  COMPLEX*16 tmp12
2895  COMPLEX*16 v3(*)
2896  REAL*8 p1(0:3)
2897  REAL*8 m1
2898  COMPLEX*16 tmp10
2899  REAL*8 p2(0:3)
2900  REAL*8 p3(0:3)
2901  REAL*8 w1
2902  COMPLEX*16 denom
2903  COMPLEX*16 tmp14
2904  COMPLEX*16 coup
2905  COMPLEX*16 tmp9
2906  COMPLEX*16 v1(6)
2907  COMPLEX*16 tmp8
2908  p2(0) = dble(v2(1))
2909  p2(1) = dble(v2(2))
2910  p2(2) = dimag(v2(2))
2911  p2(3) = dimag(v2(1))
2912  p3(0) = dble(v3(1))
2913  p3(1) = dble(v3(2))
2914  p3(2) = dimag(v3(2))
2915  p3(3) = dimag(v3(1))
2916  v1(1) = +v2(1)+v3(1)
2917  v1(2) = +v2(2)+v3(2)
2918  p1(0) = -dble(v1(1))
2919  p1(1) = -dble(v1(2))
2920  p1(2) = -dimag(v1(2))
2921  p1(3) = -dimag(v1(1))
2922  tmp14 = (v2(3)*v3(3)-v2(4)*v3(4)-v2(5)*v3(5)-v2(6)*v3(6))
2923  tmp12 = (p3(0)*v2(3)-p3(1)*v2(4)-p3(2)*v2(5)-p3(3)*v2(6))
2924  tmp9 = (v3(3)*p2(0)-v3(4)*p2(1)-v3(5)*p2(2)-v3(6)*p2(3))
2925  tmp8 = (v3(3)*p1(0)-v3(4)*p1(1)-v3(5)*p1(2)-v3(6)*p1(3))
2926  tmp10 = (v2(3)*p1(0)-v2(4)*p1(1)-v2(5)*p1(2)-v2(6)*p1(3))
2927  denom = coup/(p1(0)**2-p1(1)**2-p1(2)**2-p1(3)**2 - m1 * (m1
2928  $ -ci* w1))
2929  v1(3)= denom*(tmp14*(-ci*(p2(0))+ci*(p3(0)))+(v2(3)*(-ci*(tmp8)
2930  $ +ci*(tmp9))+v3(3)*(-ci*(tmp12)+ci*(tmp10))))
2931  v1(4)= denom*(tmp14*(-ci*(p2(1))+ci*(p3(1)))+(v2(4)*(-ci*(tmp8)
2932  $ +ci*(tmp9))+v3(4)*(-ci*(tmp12)+ci*(tmp10))))
2933  v1(5)= denom*(tmp14*(-ci*(p2(2))+ci*(p3(2)))+(v2(5)*(-ci*(tmp8)
2934  $ +ci*(tmp9))+v3(5)*(-ci*(tmp12)+ci*(tmp10))))
2935  v1(6)= denom*(tmp14*(-ci*(p2(3))+ci*(p3(3)))+(v2(6)*(-ci*(tmp8)
2936  $ +ci*(tmp9))+v3(6)*(-ci*(tmp12)+ci*(tmp10))))
2937  END
Generated on Tue Mar 12 2024 22:57:51 for C++InterfacetoTauola by   doxygen 1.8.7