C++InterfacetoTauola
ffwid3pi.f
1  DOUBLE PRECISION FUNCTION ffwid3pi(QQ,S1,S3)
2  IMPLICIT NONE
3  DOUBLE PRECISION qq,s1,s3
4 C **************************************************************
5 C Input: QQ S1 S3 ! mpi-pi-pi+**2 mpi-pi+**2 mpi-pi-**2
6 C Calls: functions FORM1, FORM2, FORM4
7 C Uses constants: tau mass, pi mass, normalization constant.
8 C Load: Initialized tauola library,
9 C Output: d\Gamma(tau --> 3pi nu)/(dQQ dS1 dS3)
10 C Remark: If QQ S1 S3 are outside of the phase space
11 C function FFWID3PI returns zero.
12 C **************************************************************
13  COMPLEX f1,f2,f4, form1,form2, form4
14  DOUBLE PRECISION v11,v12,v22,ggf2,vud2,abs1,qqmin,
15  & qqmax,s3max,s3min,s1min,s1max
16  REAL xqq,xs1,xs3, xs2,rqq
17  COMMON / decpar / gfermi,gv,ga,ccabib,scabib,gamel
18  REAL*4 gfermi,gv,ga,ccabib,scabib,gamel
19  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
20  * ,ampiz,ampi,amro,gamro,ama1,gama1
21  * ,amk,amkz,amkst,gamkst
22  DOUBLE PRECISION xlam,x,y,z
23  DOUBLE PRECISION xampi2
24  DOUBLE PRECISION getfpirpt
25  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
26  * ,ampiz,ampi,amro,gamro,ama1,gama1
27  * ,amk,amkz,amkst,gamkst
28  DOUBLE PRECISION pi
29  DATA pi /3.141592653589793238462643d0/
30  INTEGER imode,idum,ifrchl
31  REAL rrq,rchlwida1pi
32 
33 
34  xlam(x,y,z)= sqrt(abs((x-y-z)**2 - 4.*y*z))
35 
36  abs1 = 1.d-5
37 
38 
39  ggf2 = gfermi**2
40  vud2 = ccabib**2
41 
42 
43 C TO CHANGE THE VARIABLES TO SINGLE PRECISION
44 C INPUT FOR FORM1,FORM2,FORM4 IS SINGLE PRECISION
45 
46  xqq = qq
47  xs1 = s1
48  xs2 = qq -s1-s3 + 3.*ampi**2
49  xs3 = s3
50  xampi2 = ampi**2
51 C Limits for PHASE SPACE
52 C Limits for QQ
53  qqmin = 9.d0*ampi**2
54  qqmax = (amtau-amnuta)**2
55 
56 C Limits for S1
57  s1max=(dsqrt(qq) - ampi)**2 -abs1
58  s1min=4.d0*ampi**2 +abs1
59 
60 C LIMIT FOR XS3
61  s3max = (qq - ampi**2)**2 -
62  & ( xlam(qq,s1,xampi2)
63  & - xlam(s1,xampi2,xampi2) )**2
64  s3min = (qq - ampi**2)**2 -
65  & (xlam(qq,s1,xampi2)
66  & + xlam(s1,xampi2,xampi2) )**2
67 
68  s3max = s3max/4./s1
69  s3min = s3min/4./s1
70 
71 C Check on PHASE SPACE
72 C
73  IF((xs2.LE.0.) .OR.(s3max.LE.s3min)
74  & .OR.(xs1.LE.s1min).OR.(xs1.GE.s1max)
75  & .OR.(xs3.LE.s3min).OR.(xs3.GE.s3max)
76  & .OR.(qq.LE.qqmin).OR.(qq.GE.qqmax)
77  & ) THEN
78  ffwid3pi = 0.d0
79  RETURN
80  ENDIF
81 
82 C
83  v11 = -xs1+4.d0*ampi**2 -(xs2-xs3)**2/(4.d0*xqq)
84  v22 = -xs2+4.d0*ampi**2 - (xs3-xs1)**2/(4.d0*xqq)
85  v12 = 0.5d0*(xs3-xs1-xs2+4.d0*ampi**2)-0.25d0*(xs3-xs2)*(xs3-xs1)/xqq
86 
87 
88  f1 = form1(0,xqq,xs1,xs2)
89  f2 = form2(0,xqq,xs2,xs1)
90  f4 = form4(0,xqq,xs2,xs1,xs3)
91 
92 C formula 3.46 of [3]
93  ffwid3pi = abs(f1*conjg(f1))*v11+abs(f2*conjg(f2))*v22+
94  $ 2.d0*REAL(f1*conjg(f2))*v12
95 
96  CALL ifgfact(2,imode,idum)
97  CALL inirchlget(ifrchl)
98  IF (imode.EQ.0) THEN
99 C VERSION A: The 3 pion contribution to the a1 width
100 C factor of a1 phase space and zeroing a1 propagator etc.is
101 C done in RCHLWIDA1PI
102  rqq=qq
103  IF (ifrchl.EQ.1) THEN
104 C CASE OF RCHL
105  ffwid3pi = rchlwida1pi(rqq,ffwid3pi)
106  ELSE
107 C NOT READY YET
108  WRITE(*,*) 'FFWID3PI is not ready for non rchl currents'
109  stop
110  ENDIF
111 C to get the total a1 width the contribution from (KKPI)- and K-K0pi0
112 C channels have to be added
113  ELSE
114 C VERSION B: calculation of 3 pion spectra in tau to 3pi nu channel.
115 
116 C factor for phase space of tau to XQQ nu decay and contribution from F4
117 C (formula 3.21 of [3])
118 
119  ffwid3pi = (- ffwid3pi/3.d0*(1.d0+2.d0*xqq/(amtau**2))+
120  $ xqq*abs(f4*conjg(f4)))*(amtau**2/xqq-1.d0)**2
121 
122 C Flux factor and normalization const.
123  ffwid3pi =
124  & ggf2*vud2/(128.d0*(2.d0*pi)**5*amtau)/2.d0
125  & *ffwid3pi
126  IF (ifrchl.EQ.1) THEN
127 C CASE OF RCHL
128 C RChL normalization constant
129  ffwid3pi =ffwid3pi/getfpirpt(1)**2
130  ELSE
131 C NOT READY YET
132  WRITE(*,*) 'FFWID3PI is not ready for non rchl currents'
133  stop
134  ENDIF
135 
136  ENDIF
137 
138  RETURN
139  END
140 
141  REAL FUNCTION rchlwida1pi(RQQ,FFWID3PI)
142 C The 3 pion contribution to the a1 width
143 C (in [3] simple pretabulation is used through formula 3.48)
144 C for calculation of g(QQ) of 3.45 3.46 of [3] in RChL style
145 C a1 propagator has to be taken with the zero width.
146 
147  IMPLICIT NONE
148  COMPLEX fa1rchl
149  REAL rqq
150  DOUBLE PRECISION ffwid3pi
151  COMMON / decpar / gfermi,gv,ga,ccabib,scabib,gamel
152  REAL*4 gfermi,gv,ga,ccabib,scabib,gamel
153  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
154  * ,ampiz,ampi,amro,gamro,ama1,gama1
155  * ,amk,amkz,amkst,gamkst
156  DOUBLE PRECISION xlam,x,y,z
157  DOUBLE PRECISION xampi2
158  DOUBLE PRECISION getfpirpt
159  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
160  * ,ampiz,ampi,amro,gamro,ama1,gama1
161  * ,amk,amkz,amkst,gamkst
162 
163  common/rcht_3pi/ fpi_rpt,fv_rpt,gv_rpt,fa_rpt,beta_rho,fk_rpt
164  & ,fv1_rpt,gv1_rpt
165  DOUBLE PRECISION fpi_rpt,fv_rpt,gv_rpt,fa_rpt,beta_rho,fk_rpt
166  & ,fv1_rpt,gv1_rpt
167  DOUBLE PRECISION pi
168  DATA pi /3.141592653589793238462643d0/
169 
170 C AMA1 should be replaced by variable from the rchl namespace.
171  rchlwida1pi=- 1.0/REAL(fa1rchl(rqq)*conjg(fa1rchl(rqq)))/rqq**2
172  $ /(96.d0*8.d0*pi**3*ama1)/(fa_rpt**2*fpi_rpt**2)
173  $ *ffwid3pi/2.d0
174  END
175 
176  DOUBLE PRECISION FUNCTION getfpirpt(I)
177  IMPLICIT NONE
178  common/rcht_3pi/ fpi_rpt,fv_rpt,gv_rpt,fa_rpt,beta_rho,fk_rpt
179  & ,fv1_rpt,gv1_rpt
180  DOUBLE PRECISION fpi_rpt,fv_rpt,gv_rpt,fa_rpt,beta_rho,fk_rpt
181  & ,fv1_rpt,gv1_rpt
182  INTEGER i
183  getfpirpt=fpi_rpt
184  END