C++ Interface to Tauola
frho_pi_belle.f
1  Complex*16 FUNCTION fpibel(W,flpar)
2 c****************************************************************************
3 c pion form factor from belle
4 c
5 c rho(770)+rho(1450)+rho(1700) model is used as in:
6 c m.fujikawa et al., phys.rev.d 78 (2008) 072006
7 c -------------------
8 c * flpar: = 0 - (all free) fit result for par(1...11)
9 c = 1 - (par(1)=f_pi(0)=1-fixed) fit result for par(2...11)
10 c
11 c * parameter : par(1)=overall norm. factor(=f_pi(0))
12 c : par(2)=rho(770) mass
13 c : par(3)=rho(770) width
14 c : par(4)=rho(1450) mass
15 c : par(5)=rho(1450) width
16 c : par(6)=rho(1450) admixture abs. value(|beta|)
17 c : par(7)=rho(1450) admixture phase(phi=arg(beta))
18 c : par(8)=rho(1700) mass
19 c : par(9)=rho(1700) width
20 c : par(10)=rho(1700) admixture abs. value(|gamma|)
21 c : par(11)=rho(1700) admixture phase(phi3=arg(gamma))
22 c
23 c * x:=s=(m_pipi0)^2=w^2
24 c
25 c * notice: the following code was extracted(and checked) directly
26 c from the belle fit function with minor cosmetic changes
27 c(this is internal comment of belle)
28 c
29 c * notice: the following code was extracted(and checked) directly
30 c from tauola_with_belle_fpi.f sent by d. epifanov and
31 c it is an identical copy of
32 c Complex*16 FUNCTION fpibel(W,flpar)
33 c used in tauola_with_belle_fpi.f
34 c
35 c * called: curr_pipi0 if (ff2pirho = 2 or ff2pirho =3)
36 c in ../rchl-currents/value_parameter.f
37 c ff2pirho = 2 is for flpar = 0
38 c ff2pirho = 3 is for flpar = 1
39 c
40 c****************************************************************************
41  implicit none
42  integer flpar
43  real W
44  real*8 x,par(11),pi,pimas,pi0mas,taumas
45  real*8 beta,berho1,berho2,berho3
46  real*8 ps,prho1,prho2,prho3
47  real*8 gamma1,gamma2,gamma3
48  real*8 hs,h1,h2,h3,dhds1,dhds2,dhds3
49  real*8 d1,d2,d3,fs1,fs2,fs3
50  real*8 a1,b1,c1,bw_re1,bw_im1
51  real*8 a2,b2,c2,bw_re2,bw_im2
52  real*8 a3,b3,c3,bw_re3,bw_im3
53  real*8 sinphi,cosphi,sinphi3,cosphi3
54  complex*16 mbw1,mbw2,mbw3,mff,mephi,mephi3
55  parameter(pi=3.141592653589)
56 
57  x=dble(w**2)
58 c------------- parameters --------------------
59 c particle mass(unit: gev)
60 
61  pimas = 0.1395702
62  pi0mas = 0.13498
63  taumas = 1.77699
64 
65 c---------------------------------------------
66 c opt. par. from table vii of prd78(2008) 072006
67 
68  if(flpar.eq.0) then
69  par(1) = 1.02
70  par(2) = 0.7749
71  par(3) = 0.1486
72  par(4) = 1.428
73  par(5) = 0.413
74  par(6) = 0.13
75  par(7) = 197
76  par(8) = 1.694
77  par(9) = 0.135
78  par(10)= 0.028
79  par(11)= -3
80  else
81  par(1) = 1.0
82  par(2) = 0.7746
83  par(3) = 0.1481
84  par(4) = 1.446
85  par(5) = 0.434
86  par(6) = 0.15
87  par(7) = 202
88  par(8) = 1.728
89  par(9) = 0.164
90  par(10)= 0.037
91  par(11)= 24
92  endif
93 
94 c========= beta(s, rho(770), rho(1450), rho(1700)) ==================
95 
96  beta =sqrt( (1.0 - (pimas - pi0mas)**2/x)
97  + *(1.0 - (pimas + pi0mas)**2/x) )
98 
99 c ----- rho(770) -----
100  berho1=sqrt( ( 1.0 - (pimas - pi0mas)**2/(par(2)*par(2)) )
101  + *( 1.0 - (pimas + pi0mas)**2/(par(2)*par(2)) ) )
102 
103 c ----- rho(1450) -----
104  berho2=sqrt( ( 1.0 - (pimas - pi0mas)**2/(par(4)*par(4)) )
105  + *( 1.0 - (pimas + pi0mas)**2/(par(4)*par(4)) ) )
106 
107 c ----- rho(1700) -----
108  berho3=sqrt( ( 1.0 - (pimas - pi0mas)**2/(par(8)*par(8)) )
109  + *( 1.0 - (pimas + pi0mas)**2/(par(8)*par(8)) ) )
110 
111 c========= momentum(s, rho(770), rho(1450), rho(1700)) ==============
112 c
113 c ps : momentum of pi at s
114 c prho1 : momentum of pi at s=rho(770)**2
115 c prho2 : momentum of pi at s=rho(1450)**2
116 c prho3 : momentum of pi at s=rho(1700)**2
117 
118  ps = 0.5 * sqrt(x) * beta
119 c ----- rho(770) -----
120  prho1 = 0.5* par(2) * berho1
121 c ----- rho(1450) -----
122  prho2 = 0.5* par(4) * berho2
123 c ----- rho(1700) -----
124  prho3 = 0.5* par(8) * berho3
125 
126 c========= width(rho(770), rho(1450), rho(1700)) ====================
127 
128 c ----- rho(770) -----
129  gamma1 = par(3) * ( par(2)*par(2)/x ) * (ps/prho1)**3
130 c ----- rho(1450) -----
131  gamma2 = par(5) * ( par(4)*par(4)/x ) * (ps/prho2)**3
132 c ----- rho(1700) -----
133  gamma3 = par(9) * ( par(8)*par(8)/x ) * (ps/prho3)**3
134 
135 c============== h(rho(770), rho(1450), rho(1700)) ====================
136 c
137 c hs : h(s)...s(s=x)
138 c h1 : h(s) for s=rho(770)**2
139 c h2 : h(s) for s=rho(1450)**2
140 c h3 : h(s) for s=rho(1700)**2
141 
142  hs = (2.0/pi)*(ps/sqrt(x))*
143  + log( (sqrt(x) + 2.0*ps)/(2.0*pimas) )
144 c ----- rho(770) -----
145  h1 = (2.0/pi)*(prho1/par(2))*
146  + log( (par(2) + 2.0*prho1)/(2.0*pimas) )
147 c ----- rho(1450) -----
148  h2 = (2.0/pi)*(prho2/par(4))*
149  + log( (par(4) + 2.0*prho2)/(2.0*pimas) )
150 c ----- rho(1700) -----
151  h3 = (2.0/pi)*(prho3/par(8))*
152  + log( (par(8) + 2.0*prho3)/(2.0*pimas) )
153 
154 c============ dhds(rho(770), rho(1450), rho(1700)) ===================
155 c
156 c dhds = dh/ds = h(m_rho)[1/(8*p(m_rho)^2) - 1/(2m_rho^2)
157 c + 1/(2pi*m_rho^2)
158 c dhds1 : dh/ds for m_rho=rho(770)
159 c dhds2 : dh/ds for m_rho=rho(1450)
160 c dhds3 : dh/ds for m_rho=rho(1700)
161 
162 c ----- rho(770) -----
163  dhds1 = h1*( 1.0/(8.0*prho1**2) - 1.0/(2.0*par(2)**2) )
164  + + 1.0/(2.0*pi*par(2)**2)
165 c ----- rho(1450) -----
166  dhds2 = h2*( 1.0/(8.0*prho2**2) - 1.0/(2.0*par(4)**2) )
167  + + 1.0/(2.0*pi*par(4)**2)
168 c ----- rho(1700) -----
169  dhds3 = h3*( 1.0/(8.0*prho3**2) - 1.0/(2.0*par(8)**2) )
170  + + 1.0/(2.0*pi*par(8)**2)
171 
172 c============ d(rho(770), rho(1450), rho(1700)) ======================
173 c
174 c d = 3/pi * (m_pi^2/p(m_rho)^2) *
175 c ln {(m_rho+2*p(m_rho))/2m_pi} + m_rho/(2pi*p(m_rho))
176 c - (m_pi^2 * m_rho)/(pi*p(m_rho)^3)
177 c d1 : d for m_rho=rho(770)
178 c d2 : d for m_rho=rho(1450)
179 c d2 : d for m_rho=rho(1700)
180 
181 c ----- rho(770) -----
182  d1 = (3.0/pi)*(pimas**2/prho1**2)
183  + * log((par(2) + 2.0*prho1)/(2.0*pimas))
184  + + par(2)/(2.0*pi*prho1)
185  + - ((pimas**2)*par(2))/(pi*prho1**3)
186 c ----- rho(1450) -----
187  d2 = (3.0/pi)*(pimas**2/prho2**2)
188  + * log((par(4) + 2.0*prho2)/(2.0*pimas))
189  + + par(4)/(2.0*pi*prho2)
190  + - ((pimas**2)*par(4))/(pi*prho2**3)
191 c ----- rho(1700) -----
192  d3 = (3.0/pi)*(pimas**2/prho3**2)
193  + * log((par(8) + 2.0*prho3)/(2.0*pimas))
194  + + par(8)/(2.0*pi*prho3)
195  + - ((pimas**2)*par(8))/(pi*prho3**3)
196 
197 c================ f(s) (rho(770), rho(1450), rho(1700) ) ==========
198 c
199 c f(s) = gamma-rho * m_rhp**2/p(m_rho)[ p(s)^2(h(s)-h(m_rho)
200 c + (m_rho^2 -s)p(m_rho)^2 * dh/ds]
201 c fs1 : f(s) for m_rho=rho(770)
202 c fs2 : f(s) for m_rho=rho(1450)
203 c fs3 : f(s) for m_rho=rho(1700)
204 c
205 c ----- rho(770) -----
206  fs1 = par(3) * (par(2)**2/prho1**3) *
207  + ( ps**2 *(hs - h1) + (par(2)**2 -x)*prho1**2 * dhds1 )
208 c ----- rho(1450) -----
209  fs2 = par(5) * (par(4)**2/prho2**3) *
210  + ( ps**2 *(hs - h2) + (par(4)**2 -x)*prho2**2 * dhds2 )
211 c ----- rho(1700) -----
212  fs3 = par(9) * (par(8)**2/prho3**3) *
213  + ( ps**2 *(hs - h3) + (par(8)**2 -x)*prho3**2 * dhds3 )
214 
215 c====== bw form g&s model(rho(770) + rho(1450) + rho(1700)) ======
216 c bw_re - real part ; bw_im - imaginary part
217 
218 c ----- rho(770) -----
219  a1=par(2)*par(2) - x + fs1
220  b1=sqrt(x)*gamma1
221  c1=par(2)*par(2) + d1*par(2)*par(3)
222 
223  bw_re1= (a1*c1)/(a1**2 + b1**2)
224  bw_im1= (b1*c1)/(a1**2 + b1**2)
225  mbw1=dcmplx(bw_re1,bw_im1)
226 
227 c ----- rho(1450) -----
228  a2=par(4)*par(4) - x + fs2
229  b2=sqrt(x)*gamma2
230  c2=par(4)*par(4)+ d2*par(4)*par(5)
231 
232  bw_re2= (a2*c2)/(a2**2 + b2**2)
233  bw_im2= (b2*c2)/(a2**2 + b2**2)
234  mbw2=dcmplx(bw_re2,bw_im2)
235 
236 c ----- rho(1700) -----
237  a3=par(8)*par(8) - x + fs3
238  b3=sqrt(x)*gamma3
239  c3=par(8)*par(8)+ d3*par(8)*par(9)
240 
241  bw_re3= (a3*c3)/(a3**2 + b3**2)
242  bw_im3= (b3*c3)/(a3**2 + b3**2)
243  mbw3=dcmplx(bw_re3,bw_im3)
244 
245 c====== admixtures of rho(1450) and rho(1700) ======
246 
247 c --------- exp(i*arg(beta)) ----------
248  sinphi=sin(par(7)*(pi/180.0)) ! use degree
249  cosphi=cos(par(7)*(pi/180.0)) ! use degree
250  mephi=dcmplx(cosphi,sinphi)
251 
252 c --------- exp(i*arg(gamma)) ---------
253  sinphi3=sin(par(11)*(pi/180.0)) ! use degree
254  cosphi3=cos(par(11)*(pi/180.0)) ! use degree
255  mephi3=dcmplx(cosphi3,sinphi3)
256 
257 c================= form factor =====================
258 
259  mff = 1.0/(1.0+par(6)*mephi+par(10)*mephi3)
260  fpibel=par(1)*mff*(mbw1+par(6)*mephi*mbw2+par(10)*mephi3*mbw3)
261  RETURN
262  END
Generated by doxygen 1.9.1