C++ Interface to Tauola
f3pi_rcht.f
1C JAK = 5
2C
3 COMPLEX FUNCTION f3pi_rcht(IFORM,QQ,SA,SB)
4 IMPLICIT NONE
5 INTEGER IFORM
6 REAL QQ,SA,SB
7C.......................................................................
8C.
9C. F3PI - RchT version of the hadronic curent used in TAUOLA
10C. References: [1] arXiv:0911.4436 (hep-ph) P. Roig et al.
11C. eqs (5)-(9), (32) gives the main part of the current
12C (the part without the sigma contribution)
13C. [2] arXiv:1203.3955,
14C. the manual for the 3pi current without the sigma contribution
15C. [3] http://annapurna.ifj.edu.pl/~wasm/RChL/RChLa1.pdf
16C eq (3) the sigma meson contribution to the 3 pi current
17C. Inputs : QQ,SA,SB - invariant masses**2 [GeV**2]
18C. : IFORM formfactor no.
19C. Outputs : F3PI_RCHT formfactor value
20C.
21C. COMMON : RCHT_3PI content is defined in this routine
22C.
23C. Calls : functions from file ./funct_rpt.f, ./fa1rchl.f
24C. Called : from file ../formf.f
25C. Author : O.S
26C. Created :
27C. Modified : 1. Feb 2011
28C a part with scalars added on the 1st May 2012
29C.......................................................................
30 include '../parameter.inc'
31 include '../funct_declar.inc'
32
33 COMMON / decpar / gfermi,gv,ga,ccabib,scabib,gamel
34 real*4 gfermi,gv,ga,ccabib,scabib,gamel
35 COMMON / taukle / bra1,brk0,brk0b,brks
36 real*4 bra1,brk0,brk0b,brks
37 INTEGER J3PI
38 REAL wid_a1_fit
39*
40 CHARACTER*(*) CRNAME
41 parameter( crname = 'F3PI' )
42*
43 INTEGER IDK
44 DOUBLE PRECISION M1,M2,M3,M1SQ,M2SQ,M3SQ
45 REAL GGMA1
46 REAL S1,S2,S3,FACT_RPT
47*
48 REAL R_RPT,KAP_RPT,FACT_ADD_RPT
49 COMPLEX ALP21_RPT,ALP22_RPT,ALP11_RPT,ALP12_RPT
50 $ ,BETA1_RPT,BETA2_RPT,BETA1_RPT_RHO1,ALP11_RPT_RHO1
51 $ ,ALP21_RPT_RHO1,ALP22_RPT_RHO1
52 COMPLEX FA1RCHL
53 INTEGER IFIRST
54 DATA ifirst/0/
55
56C. GENERAL INITIALIZATION
57C. ======================
58
59 IF (ifirst.EQ.0) THEN
60 ifirst = 1
61 print *,' In F3pi: chpt + 1 resonance + 2 resonance (RchT)'
62 END IF
63
64C******************************************
65C Initilisation of the mass of the particles
66C*****************************************
67 call rchl_parameters(5)
68
69C. We impose isospin symmetry requesting that charged and neutral pion mass
70C. are equal. This may need to be changed MMPI_AV = (2.*MPIC+MPIZ)/3.
71
72 m1 = mmpi_av
73 m2 = mmpi_av
74 m3 = mmpi_av
75
76
77C. normalization factor to compensate on different
78C. convention for normalization constant used in TAUOLA and
79C. TAUOLA documentation on one side and paper [1] on other.
80 fact_add_rpt = 1./fpi_rpt
81
82
83
84C. FUNCTION VALUE, GENERAL CASE
85C. ============================
86
87C Function value set to 0
88 f3pi_rcht = cmplx(0.,0.)
89
90C. First determine whether we are doing pi-2pi0 or 3pi.
91C. we CH3PI GET information (eg from phase space geneator of tauola.f)
92C. whether it is 3 prong J3pi=1 or 1 prong J3pi=2 final state of 3 pion.
93 CALL ch3piget(j3pi)
94
95 IF (j3pi.EQ.2) THEN ! it is pi 2pi0
96 idk = 1
97 r_rpt = -1.
98 kap_rpt = 0.5
99 ELSE IF(j3pi.EQ.1) THEN ! it is 3pi
100 idk = 2
101 r_rpt = 1.
102 kap_rpt = 1.
103 END IF
104
105 m1sq = m1*m1
106 m2sq = m2*m2
107 m3sq = m3*m3
108
109
110C. Calculation, IFORM = 1 or 2.
111C. VECTOR 3 PION FORM FACTORS
112C. ============================
113 IF (iform.EQ.1.OR.iform.EQ.2) THEN
114 s1 = sa ! t variable in [2] (!!! vec1 = v2 !!!)
115 s2 = sb ! s variable in [2]
116 s3 = qq-sa-sb+m1sq+m2sq+m3sq
117
118 IF (s3.LE.0..OR.s2.LE.0.) RETURN
119
120C. FUNCTIONS BETA_RPT, ALP1_RPT are coded in ./funct_rpt.f
121C. they are defined in Eq (6) of [2]
122 beta1_rpt = beta_rpt(qq,s1,s2,m1sq,m2sq,m3sq,mro,mrho1,grho1)
123 alp11_rpt = alp1_rpt(qq,s1,s2,m1sq,m2sq,m3sq,mro,mrho1,grho1)
124
125 f3pi_rcht = - 2.*sqrt(2.)/(3.*fpi_rpt) -
126 $ sqrt(2.)*fv_rpt*gv_rpt/(3.*fpi_rpt**3)*alp11_rpt +
127 $ 4.*fa_rpt*gv_rpt/(3.*fpi_rpt**3)*beta1_rpt*qq
128 $ *fa1rchl(qq)
129
130C. FA1RCHL(QQ) is the a1 propagator, is it coded in ./fa1rchl.f
131
132
133 f3pi_rcht = f3pi_rcht * r_rpt
134
135C.
136C. The contribution from the scalar (sigma) resonance
137C.
138 IF(ff3piscal.EQ.1) THEN
139C. This parametrization does not fit
140 IF(j3pi.EQ.1) THEN
141 f3pi_rcht = f3pi_rcht + (
142 & sqrt(2.)*(r0scal_3pi(qq,s1)+r0scal_3pi(qq,s2)) +
143 & (r2scal_3pi(qq,s1)+r2scal_3pi(qq,s2))
144 & )
145 ELSE IF(j3pi.EQ.2) THEN
146 f3pi_rcht = f3pi_rcht - (
147 & (r0scal_3pi(qq,s1)+r0scal_3pi(qq,s2)) -
148 & sqrt(2.)*(r2scal_3pi(qq,s1)+r2scal_3pi(qq,s2))
149 & )
150 ENDIF
151 ELSE IF(ff3piscal.EQ.2) THEN
152C. A new parametrization 10.02.2013 for the scalar contribution
153C. analytical formulae in [3], eqs (3)-(6),
154C. functions BWsig(Mm,Gg,Qx), FFsig(QQ,Qx) coded in ./funct_rpt.f
155C.
156 IF(j3pi.EQ.1) THEN
157 f3pi_rcht = f3pi_rcht +
158 & sqrt(2.)*fv_rpt*gv_rpt/(3.*fpi_rpt**3)*
159 & (alpsig*bwsig(msig,gsig,s1)*ffsig(qq,s1) +
160 & betasig*bwsig(msig,gsig,s2)*ffsig(qq,s2))
161 & +4.*fa_rpt*gv_rpt/(3.*fpi_rpt**3)*
162 & (gamsig*bwsig(msig,gsig,s1)*ffsig(qq,s1) +
163 & delsig*bwsig(msig,gsig,s2)*ffsig(qq,s2))
164 $ *fa1rchl(qq)
165 ELSE IF(j3pi.EQ.2) THEN
166 f3pi_rcht = f3pi_rcht - (
167 & sqrt(2.)*fv_rpt*gv_rpt/(3.*fpi_rpt**3)*
168 & alpsig*bwsig(msig,gsig,s3)*ffsig(qq,s3)
169 & +4.*fa_rpt*gv_rpt/(3.*fpi_rpt**3)*
170 & gamsig*bwsig(msig,gsig,s3)*ffsig(qq,s3)
171 $ *fa1rchl(qq) )
172 ENDIF
173
174C.
175C. The Coulomb interaction effects for the final pions
176C. only for pi-pi+pi- should be included
177C. Functions fattcoul(m1,m2,s),frepcoul(m1,m2,s) are coded in ./funct_rpt.f
178C.
179 IF (fcoul.EQ.1.AND.(j3pi.EQ.1)) THEN ! it is 3pi
180 f3pi_rcht = f3pi_rcht*sqrt(fattcoul(m2,m3,s1)
181 & *fattcoul(m1,m3,s2)*frepcoul(m1,m2,s3))
182 END IF
183 ENDIF
184
185
186C.
187C. The factor 1/FACT_ADD_RPT = FPI_RPT comes to compensate an additional factor
188C. in the hadronic current F3pi_rcht above compare with eq (6) in [2]
189
190
191 f3pi_rcht = f3pi_rcht/fact_add_rpt
192
193C. Calculation, for IFORM = 3 is not needed
194C. =======================
195C. F3PI_RCHT is set to zero in ../formf.f.
196
197C. Calculation, IFORM = 4.
198C. PSEUDOSCALAR 3 PION FORM FACTOR
199C. =======================
200 ELSEIF (iform.EQ.4) THEN
201 s1 = sa
202 s2 = sb
203 s3 = qq-sa-sb+m1sq+m2sq+m3sq
204
205 IF (s3.LE.0..OR.s2.LE.0.) RETURN
206
207C. Functions ALP21_RPT, ALP22_RPT are Eq(10) in [2]
208C. GRHO_RCHT, GRHO1_RCHT s1 or s2 dependent widths of rho or rho1
209C. coded in ./funct_rpt.f
210 alp21_rpt = 3.*gv_rpt*s1*m1sq*(s3-s2)/
211 & (fv_rpt*qq*(qq-m1sq)*(1.+beta_rho))*
212 & (1./(s1-mro**2+i*mro*grho_rcht(s1,mro))+
213 & beta_rho/(s1-mrho1**2+i*mrho1*grho1_rcht(s1,mrho1,grho1)))
214 alp22_rpt = 3.*gv_rpt*s2*m1sq*(s3-s1)/
215 & (fv_rpt*qq*(qq-m1sq)*(1.+beta_rho))*
216 & (1./(s2-mro**2+i*mro*grho_rcht(s2,mro))+
217 & beta_rho/(s2-mrho1**2+i*mrho1*grho1_rcht(s2,mrho1,grho1)))
218
219C. PSEUDOSCALAR 3 PION FORM FACTOR. Eqs (9) in [2]
220 f3pi_rcht = sqrt(2.)/(3.*fpi_rpt*qq*(qq-m1sq))*
221 & m1sq*(3.*(s3-m3sq)-qq*(1.+2.*kap_rpt*r_rpt))
222
223 f3pi_rcht = f3pi_rcht - sqrt(2.)*fv_rpt*gv_rpt/(3.*fpi_rpt**3)*
224 & ( alp21_rpt + alp22_rpt)
225
226C.
227C. The factor 1/FACT_ADD_RPT = FPI_RPT comes to compensate an additional factor
228C. in the hadronic current F3pi_rcht above compare with eq (9) in [2]
229 f3pi_rcht = r_rpt * f3pi_rcht/fact_add_rpt
230
231C.
232C. The Coulomb interaction effects for the final pions
233C. only for pi-pi+pi- should be included
234C. Functions fattcoul(m1,m2,s),frepcoul(m1,m2,s) are coded in ./funct_rpt.f
235C.
236 IF (fcoul.EQ.1.AND.(j3pi.EQ.1)) THEN ! it is 3pi
237 f3pi_rcht = f3pi_rcht*sqrt(fattcoul(m2,m3,s1)
238 & *fattcoul(m1,m3,s2)*frepcoul(m1,m2,s3))
239 END IF
240
241 END IF
242
243 RETURN
244 END
245