C++ Interface to Tauola
tauola-factory/Standart_Tauola/curr_cleo.F
1*AJW CLEO version of CURR from KORALB.
2 SUBROUTINE curr_cleo(MNUM,PIM1,PIM2,PIM3,PIM4,HADCUR)
3C ==================================================================
4C AJW, 11/97 - based on original CURR from TAUOLA:
5C hadronic current for 4 pi final state
6C R. Fisher, J. Wess and F. Wagner Z. Phys C3 (1980) 313
7C R. Decker Z. Phys C36 (1987) 487.
8C M. Gell-Mann, D. Sharp, W. Wagner Phys. Rev. Lett 8 (1962) 261.
9C BUT, rewritten to be more general and less "theoretical",
10C using parameters tuned by Vasia and DSC.
11C ==================================================================
12
13 COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
14 * ,ampiz,ampi,amro,gamro,ama1,gama1
15 * ,amk,amkz,amkst,gamkst
16C
17 real*4 amtau,amnuta,amel,amnue,ammu,amnumu
18 * ,ampiz,ampi,amro,gamro,ama1,gama1
19 * ,amk,amkz,amkst,gamkst
20C
21 REAL PIM1(4),PIM2(4),PIM3(4),PIM4(4)
22 COMPLEX HADCUR(4)
23
24 INTEGER K,L,MNUM,K1,K2,IRO,I,J,KK
25 REAL PA(4),PB(4),PAA(4)
26 REAL AA(4,4),PP(4,4)
27 REAL A,XM,XG,G1,G2,G,AMRO2,GAMRO2,AMRO3,GAMRO3,AMOM,GAMOM
28 REAL FRO,COEF1,FPI,COEF2,QQ,SK,DENOM,SIG,QQA,SS23,SS24,SS34,QP1P2
29 REAL QP1P3,QP1P4,P1P2,P1P3,P1P4,SIGN
30 REAL PKORB,AMPA
31 COMPLEX ALF0,ALF1,ALF2,ALF3
32 COMPLEX LAM0,LAM1,LAM2,LAM3
33 COMPLEX BET1,BET2,BET3
34 COMPLEX FORM1,FORM2,FORM3,FORM4,FORM2PI
35 COMPLEX BWIGM,WIGFOR,FPIKM,FPIKMD
36 COMPLEX AMPL(7),AMPR
37 COMPLEX BWIGN
38C
39 bwign(a,xm,xg)=1.0/cmplx(a-xm**2,xm*xg)
40C*******************************************************************************
41C
42C --- masses and constants
43 IF (g1.NE.12.924) THEN
44 g1=12.924
45 g2=1475.98
46 fpi=93.3e-3
47 g =g1*g2
48 fro=0.266*amro**2
49 coef1=2.0*sqrt(3.0)/fpi**2
50 coef2=fro*g ! overall constant for the omega current
51 coef2= coef2*0.56 ! factor 0.56 reduces contribution of omega from 68% to 40 %
52
53C masses and widths for for rho-prim and rho-bis:
54 amro2 = 1.465
55 gamro2= 0.310
56 amro3=1.700
57 gamro3=0.235
58C
59 amom = pkorb(1,14)
60 gamom = pkorb(2,14)
61 amro2 = pkorb(1,21)
62 gamro2= pkorb(2,21)
63 amro3 = pkorb(1,22)
64 gamro3= pkorb(2,22)
65C
66C Amplitudes for (pi-pi-pi0pi+) -> PS, rho0, rho-, rho+, omega.
67 ampl(1) = cmplx(pkorb(3,31)*coef1,0.)
68 ampl(2) = cmplx(pkorb(3,32)*coef1,0.)*cexp(cmplx(0.,pkorb(3,42)))
69 ampl(3) = cmplx(pkorb(3,33)*coef1,0.)*cexp(cmplx(0.,pkorb(3,43)))
70 ampl(4) = cmplx(pkorb(3,34)*coef1,0.)*cexp(cmplx(0.,pkorb(3,44)))
71 ampl(5) = cmplx(pkorb(3,35)*coef2,0.)*cexp(cmplx(0.,pkorb(3,45)))
72C Amplitudes for (pi0pi0pi0pi-) -> PS, rho-.
73 ampl(6) = cmplx(pkorb(3,36)*coef1)
74 ampl(7) = cmplx(pkorb(3,37)*coef1)
75C
76C rho' contributions to rho' -> pi-omega:
77 alf0 = cmplx(pkorb(3,51),0.0)
78 alf1 = cmplx(pkorb(3,52)*amro**2,0.0)
79 alf2 = cmplx(pkorb(3,53)*amro2**2,0.0)
80 alf3 = cmplx(pkorb(3,54)*amro3**2,0.0)
81C rho' contribtions to rho' -> rhopipi:
82 lam0 = cmplx(pkorb(3,55),0.0)
83 lam1 = cmplx(pkorb(3,56)*amro**2,0.0)
84 lam2 = cmplx(pkorb(3,57)*amro2**2,0.0)
85 lam3 = cmplx(pkorb(3,58)*amro3**2,0.0)
86C rho contributions to rhopipi, rho -> 2pi:
87 bet1 = cmplx(pkorb(3,59)*amro**2,0.0)
88 bet2 = cmplx(pkorb(3,60)*amro2**2,0.0)
89 bet3 = cmplx(pkorb(3,61)*amro3**2,0.0)
90C
91 END IF
92C**************************************************
93C
94C --- initialization of four vectors
95 DO 7 k=1,4
96 DO 8 l=1,4
97 8 aa(k,l)=0.0
98 hadcur(k)=cmplx(0.0)
99 paa(k)=pim1(k)+pim2(k)+pim3(k)+pim4(k)
100 pp(1,k)=pim1(k)
101 pp(2,k)=pim2(k)
102 pp(3,k)=pim3(k)
103 7 pp(4,k)=pim4(k)
104C
105 IF (mnum.EQ.1) THEN
106C ===================================================================
107C pi- pi- p0 pi+ case ====
108C ===================================================================
109 qq=paa(4)**2-paa(3)**2-paa(2)**2-paa(1)**2
110
111C Add M(4pi)-dependence to rhopipi channels:
112 form4= lam0+lam1*bwign(qq,amro,gamro)
113 * +lam2*bwign(qq,amro2,gamro2)
114 * +lam3*bwign(qq,amro3,gamro3)
115
116C --- loop over five contributions of the rho-pi-pi
117 DO 201 k1=1,3
118 DO 201 k2=3,4
119C
120 IF (k2.EQ.k1) THEN
121 GOTO 201
122 ELSEIF (k2.EQ.3) THEN
123C rho-
124 ampr = ampl(3)
125 ampa = ampiz
126 ELSEIF (k1.EQ.3) THEN
127C rho+
128 ampr = ampl(4)
129 ampa = ampiz
130 ELSE
131C rho0
132 ampr = ampl(2)
133 ampa = ampi
134 END IF
135C
136 sk=(pp(k1,4)+pp(k2,4))**2-(pp(k1,3)+pp(k2,3))**2
137 $ -(pp(k1,2)+pp(k2,2))**2-(pp(k1,1)+pp(k2,1))**2
138
139C -- definition of AA matrix
140C -- cronecker delta
141 DO 202 i=1,4
142 DO 203 j=1,4
143 203 aa(i,j)=0.0
144 202 aa(i,i)=1.0
145C ... and the rest ...
146 DO 204 l=1,4
147 IF (l.NE.k1.AND.l.NE.k2) THEN
148 denom=(paa(4)-pp(l,4))**2-(paa(3)-pp(l,3))**2
149 $ -(paa(2)-pp(l,2))**2-(paa(1)-pp(l,1))**2
150 DO 205 i=1,4
151 DO 205 j=1,4
152 sig= 1.0
153 IF(j.NE.4) sig=-sig
154 aa(i,j)=aa(i,j)
155 $ -sig*(paa(i)-2.0*pp(l,i))*(paa(j)-pp(l,j))/denom
156 205 CONTINUE
157 ENDIF
158 204 CONTINUE
159C
160C --- lets add something to HADCURR
161C FORM1= FPIKM(SQRT(SK),AMPI,AMPI) *FPIKM(SQRT(QQ),AMPI,AMPI)
162C FORM1= AMPL(1)+AMPR*FPIKM(SQRT(SK),AMPI,AMPI)
163
164 form2pi= bet1*bwigm(sk,amro,gamro,ampa,ampi)
165 1 +bet2*bwigm(sk,amro2,gamro2,ampa,ampi)
166 2 +bet3*bwigm(sk,amro3,gamro3,ampa,ampi)
167 form1= ampl(1)+ampr*form2pi
168C
169 DO 206 i=1,4
170 DO 206 j=1,4
171 hadcur(i)=hadcur(i)+form1*form4*aa(i,j)*(pp(k1,j)-pp(k2,j))
172 206 CONTINUE
173C --- end of the rho-pi-pi current (5 possibilities)
174 201 CONTINUE
175C
176C ===================================================================
177C Now modify the coefficient for the omega-pi current: =
178C ===================================================================
179 IF (ampl(5).EQ.cmplx(0.,0.)) GOTO 311
180
181C Overall rho+rhoprime for the 4pi system:
182C FORM2=AMPL(5)*(BWIGN(QQ,AMRO,GAMRO)+ELPHA*BWIGN(QQ,AMROP,GAMROP))
183C Modified M(4pi)-dependence:
184 form2=ampl(5)*(alf0+alf1*bwign(qq,amro,gamro)
185 * +alf2*bwign(qq,amro2,gamro2)
186 * +alf3*bwign(qq,amro3,gamro3))
187C
188C --- there are two possibilities for omega current
189C --- PA PB are corresponding first and second pi-s
190 DO 301 kk=1,2
191 DO 302 i=1,4
192 pa(i)=pp(kk,i)
193 pb(i)=pp(3-kk,i)
194 302 CONTINUE
195C --- lorentz invariants
196 qqa=0.0
197 ss23=0.0
198 ss24=0.0
199 ss34=0.0
200 qp1p2=0.0
201 qp1p3=0.0
202 qp1p4=0.0
203 p1p2 =0.0
204 p1p3 =0.0
205 p1p4 =0.0
206 DO 303 k=1,4
207 sign=-1.0
208 IF (k.EQ.4) sign= 1.0
209 qqa=qqa+sign*(paa(k)-pa(k))**2
210 ss23=ss23+sign*(pb(k) +pim3(k))**2
211 ss24=ss24+sign*(pb(k) +pim4(k))**2
212 ss34=ss34+sign*(pim3(k)+pim4(k))**2
213 qp1p2=qp1p2+sign*(paa(k)-pa(k))*pb(k)
214 qp1p3=qp1p3+sign*(paa(k)-pa(k))*pim3(k)
215 qp1p4=qp1p4+sign*(paa(k)-pa(k))*pim4(k)
216 p1p2=p1p2+sign*pa(k)*pb(k)
217 p1p3=p1p3+sign*pa(k)*pim3(k)
218 p1p4=p1p4+sign*pa(k)*pim4(k)
219 303 CONTINUE
220C
221C omega -> rho pi for the 3pi system:
222C FORM3=BWIGN(QQA,AMOM,GAMOM)*(BWIGN(SS23,AMRO,GAMRO)+
223C $ BWIGN(SS24,AMRO,GAMRO)+BWIGN(SS34,AMRO,GAMRO))
224C No omega -> rho pi; just straight omega:
225 form3=bwign(qqa,amom,gamom)
226C
227 DO 304 k=1,4
228 hadcur(k)=hadcur(k)+form2*form3*(
229 $ pb(k)*(qp1p3*p1p4-qp1p4*p1p3)
230 $ +pim3(k)*(qp1p4*p1p2-qp1p2*p1p4)
231 $ +pim4(k)*(qp1p2*p1p3-qp1p3*p1p2) )
232 304 CONTINUE
233 301 CONTINUE
234 311 CONTINUE
235C
236 ELSE
237C ===================================================================
238C pi0 pi0 p0 pi- case ====
239C ===================================================================
240 qq=paa(4)**2-paa(3)**2-paa(2)**2-paa(1)**2
241
242C --- loop over three contribution of the non-omega current
243 DO 101 k=1,3
244 sk=(pp(k,4)+pim4(4))**2-(pp(k,3)+pim4(3))**2
245 $ -(pp(k,2)+pim4(2))**2-(pp(k,1)+pim4(1))**2
246
247C -- definition of AA matrix
248C -- cronecker delta
249 DO 102 i=1,4
250 DO 103 j=1,4
251 103 aa(i,j)=0.0
252 102 aa(i,i)=1.0
253C
254C ... and the rest ...
255 DO 104 l=1,3
256 IF (l.NE.k) THEN
257 denom=(paa(4)-pp(l,4))**2-(paa(3)-pp(l,3))**2
258 $ -(paa(2)-pp(l,2))**2-(paa(1)-pp(l,1))**2
259 DO 105 i=1,4
260 DO 105 j=1,4
261 sig=1.0
262 IF(j.NE.4) sig=-sig
263 aa(i,j)=aa(i,j)
264 $ -sig*(paa(i)-2.0*pp(l,i))*(paa(j)-pp(l,j))/denom
265 105 CONTINUE
266 ENDIF
267 104 CONTINUE
268
269C --- lets add something to HADCURR
270C FORM1= FPIKM(SQRT(SK),AMPI,AMPI) *FPIKMD(SQRT(QQ),AMPI,AMPI)
271CCCCCCCCCCCCC FORM1=WIGFOR(SK,AMRO,GAMRO) (tests)
272C FORM1= FPIKM(SQRT(SK),AMPI,AMPI) *FPIKM(SQRT(QQ),AMPI,AMPI)
273 form1 = ampl(6)+ampl(7)*fpikm(sqrt(sk),ampi,ampi)
274
275 DO 106 i=1,4
276 DO 106 j=1,4
277 hadcur(i)=hadcur(i)+form1*aa(i,j)*(pp(k,j)-pp(4,j))
278 106 CONTINUE
279C --- end of the non omega current (3 possibilities)
280 101 CONTINUE
281
282 ENDIF
283 END
284
285
286