C++InterfacetoTauola
gfact.f
1  DOUBLE PRECISION FUNCTION gfact(QQ)
2  IMPLICIT NONE
3 C factor G to be used as inteligent retabulation as in paper
4 C Kuhn Santamaria
5  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
6  * ,ampiz,ampi,amro,gamro,ama1,gama1
7  * ,amk,amkz,amkst,gamkst
8 C
9  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
10  * ,ampiz,ampi,amro,gamro,ama1,gama1
11  * ,amk,amkz,amkst,gamkst
12  DOUBLE PRECISION dgamqq,qq,thr,thu,th0,lam
13  INTEGER ii,init,inum
14  DOUBLE PRECISION af4,af5,aa1,aa2,aa4
15  DOUBLE PRECISION f1,f2,f3,f4,f5
16  DOUBLE PRECISION aa,bb,cc,dd,a,b,c,d,x
17  SAVE a,b,c,d,af4,af5,aa,bb,cc,dd,thr,thu,th0
18 
19  CALL ifgfact(1,ii,init)
20  IF(init.EQ.0) THEN
21  thr=(ampi+amro)**2
22  th0=9*ampi**2
23  thu=thr
24  CALL getff3piscal(inum) ! we switch off part of the contr to 3 pi mode
25  CALL setff3piscal(0)
26 C below THR , we calculate at distance 1/4 1/2 and 1 between
27 C minimum and maximum for this range
28  lam=(thr-th0)/4
29  aa1=dgamqq(th0+ lam)
30  aa2=dgamqq(th0+2*lam)
31  aa4=dgamqq(th0+4*lam)
32 C above THR we calculate at THR times 1 2 3 and 4 for higher range
33  f1=dgamqq(thu)
34  f2=dgamqq(1.5*thu)
35  f3=dgamqq(2.*thu)
36  f4=dgamqq(3*thu)
37  f5=dgamqq(3.5*thu)
38  CALL setff3piscal(inum) ! we switch back part of the contr to 3 pi mode
39 
40 C we calculate coefs for expansion ( polynomial of order 2 once
41 C X**2 factorized out, X=QQ-TH0 )
42 
43  aa1 = aa1/lam**2 ! RESCALING due to factorized X**2
44  aa2 = aa2/(2.*lam)**2 !
45  aa4 = aa4/(4.*lam)**2 !
46 
47  bb=1./8.*(10.*aa2-aa4-16.*aa1)
48  aa=(8.*aa1-aa2-4.*bb)/6.
49  cc=aa1-aa-bb
50  aa=aa/lam
51  bb=bb/lam**2
52  cc=cc/lam**3
53 C calulate coefs, assuming it is polynomial order 3 note negative pwrs
54  d=-9*(-4*f3-f1+4*f2+f4)
55  c=3*(f4+f1-2*f3-11./18.*d)
56  a=f3-f1+0.5*c+0.75*d
57  b=f1-a-c-d
58  a=a/thu
59  b=b
60  c=c*thu
61  d=d*thu**2
62  af4=f4
63  af5=f5
64 c write(*,*) "A=",AA,"B=",-BB,"C=",CC,"D=",A,"E=",-B,"F=",C
65 c write(*,*) "G=",-D,"H=",AF4,"P=",AF5-AF4
66  ENDIF
67  IF (qq.GT.3*thu) THEN
68  gfact=af4+(af5-af4)*(qq-3*thu)*2/thu
69  ELSEIF (qq.GT.thr) THEN
70  gfact=a*qq+b+c/qq+d/qq**2
71  ELSEIF(qq.LE.th0) THEN
72  gfact=0.0
73  ELSE
74  x=qq-th0
75  gfact=aa*x+bb*x**2+cc*x**3
76  gfact=x**2*gfact
77  ENDIF
78  END
79 
80 
81  DOUBLE PRECISION FUNCTION dgamqq(XQQB)
82 C **************************************************************
83 C calculates \tau^- -> pi^- pi^- pi^+ nu width as function of QQ (XQQB)
84 C formulas (19) of ref [2a] integration over S1
85 C limit of integration (21) of ref [2a] see also [4]
86 C called from main function
87 C **************************************************************
88  IMPLICIT NONE
89  common/precint/ epssq,abs1
90  DOUBLE PRECISION epssq,abs1
91  COMMON /EXTERNAL/ xqqa
92  DOUBLE PRECISION xqqa
93  EXTERNAL gaus,dgamqqs1
94  DOUBLE PRECISION gaus,dgamqqs1
95  DOUBLE PRECISION xqqb,eps,ups1,downs1
96  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
97  * ,ampiz,ampi,amro,gamro,ama1,gama1
98  * ,amk,amkz,amkst,gamkst
99 C
100  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
101  * ,ampiz,ampi,amro,gamro,ama1,gama1
102  * ,amk,amkz,amkst,gamkst
103 
104  xqqa = xqqb
105  eps = epssq/3.d0
106  ups1=(dsqrt(xqqb) - ampi)**2-abs1 ! limits on S1
107  downs1=4.d0*ampi**2+abs1
108 
109  dgamqq = gaus(dgamqqs1,downs1,ups1,eps)
110 
111  RETURN
112  END
113 
114 
115 
116  DOUBLE PRECISION FUNCTION dgamqqs1(S1)
117  IMPLICIT NONE
118  DOUBLE PRECISION s1
119 C **************************************************************
120 C calculates \tau^- -> pi^- pi^- pi^+ nu width
121 C as function of QQ,S1
122 C GAUS integrant in DGAMQQ (XQQA- hidden argument)
123 C calculates tau^- -> pi^- pi^- pi^+ nu spectrum as function of S1
124 C formulas (19) of ref [2a] see also [4]
125 C limit of integration (21) of ref [2a]
126 C **************************************************************
127  EXTERNAL gaus2,ffwid3pi,dgamqqs1s3
128  DOUBLE PRECISION gaus2,ffwid3pi,dgamqqs1s3
129  COMMON /internal/ s1a
130  DOUBLE PRECISION s1a
131  COMMON /EXTERNAL/ xqqa
132  DOUBLE PRECISION xqqa
133  common/precint/ epssq,abs1
134  DOUBLE PRECISION epssq,abs1
135  DOUBLE PRECISION eps,ups3,downs3
136  DOUBLE PRECISION xlam,x,y,z
137  DOUBLE PRECISION xampi2
138  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
139  * ,ampiz,ampi,amro,gamro,ama1,gama1
140  * ,amk,amkz,amkst,gamkst
141 C
142  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
143  * ,ampiz,ampi,amro,gamro,ama1,gama1
144  * ,amk,amkz,amkst,gamkst
145 
146  xlam(x,y,z) = sqrt(abs((x-y-z)**2 - 4.*y*z))
147  s1a = s1
148  eps = epssq/9.d0
149  xampi2=ampi**2
150 
151  ups3 = (xqqa - ampi**2)**2 - ! limits on S3
152  & ( xlam(xqqa,s1,xampi2)
153  & - xlam(s1,xampi2,xampi2) )**2
154  downs3 = (xqqa - ampi**2)**2 -
155  & (xlam(xqqa,s1,xampi2)
156  & + xlam(s1,xampi2,xampi2) )**2
157 
158  ups3 = ups3/4./s1
159  downs3 = downs3/4./s1
160 
161  dgamqqs1 = gaus2(dgamqqs1s3,downs3,ups3,eps)
162 
163  RETURN
164  END
165 
166 
167  DOUBLE PRECISION FUNCTION dgamqqs1s3(XS3)
168  IMPLICIT NONE
169 C **************************************************************
170 C calculates \tau^- -> pi^- pi^- pi^+ nu width
171 C as function of QQ,S1,S3
172 C **************************************************************
173 C EXTERNAL GAUS2,FFWID3PI
174  DOUBLE PRECISION ffwid3pi,xs3
175  COMMON /internal/ xs1a
176  DOUBLE PRECISION xs1a
177  COMMON /EXTERNAL/ xqqa
178  DOUBLE PRECISION xqqa
179  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
180  * ,ampiz,ampi,amro,gamro,ama1,gama1
181  * ,amk,amkz,amkst,gamkst
182 C
183  REAL*4 amtau,amnuta,amel,amnue,ammu,amnumu
184  * ,ampiz,ampi,amro,gamro,ama1,gama1
185  * ,amk,amkz,amkst,gamkst
186 
187 
188  dgamqqs1s3 = ffwid3pi(xqqa,xs1a,xs3)
189 
190  RETURN
191 
192  END