C++InterfacetoTauola
tauola-F/curr_cleo.F
1 *AJW CLEO version of CURR from KORALB.
2  SUBROUTINE curr_cleo(MNUM,PIM1,PIM2,PIM3,PIM4,HADCUR)
3 C ==================================================================
4 C AJW, 11/97 - based on original CURR from TAUOLA:
5 C hadronic current for 4 pi final state
6 C R. Fisher, J. Wess and F. Wagner Z. Phys C3 (1980) 313
7 C R. Decker Z. Phys C36 (1987) 487.
8 C M. Gell-Mann, D. Sharp, W. Wagner Phys. Rev. Lett 8 (1962) 261.
9 C BUT, rewritten to be more general and less "theoretical",
10 C using parameters tuned by Vasia and DSC.
11 C ==================================================================
12 
13  COMMON / parmas / amtau,amnuta,amel,amnue,ammu,amnumu
14  * ,ampiz,ampi,amro,gamro,ama1,gama1
15  * ,amk,amkz,amkst,gamkst
16 C
17  REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
18  * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
19  * ,AMK,AMKZ,AMKST,GAMKST
20 C
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
38 C
39  bwign(a,xm,xg)=1.0/cmplx(a-xm**2,xm*xg)
40 C*******************************************************************************
41 C
42 C --- 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 
53 C 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
58 C
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)
65 C
66 C 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)))
72 C Amplitudes for (pi0pi0pi0pi-) -> PS, rho-.
73  ampl(6) = cmplx(pkorb(3,36)*coef1)
74  ampl(7) = cmplx(pkorb(3,37)*coef1)
75 C
76 C 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)
81 C 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)
86 C 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)
90 C
91  END IF
92 C**************************************************
93 C
94 C --- 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)
104 C
105  IF (mnum.EQ.1) THEN
106 C ===================================================================
107 C pi- pi- p0 pi+ case ====
108 C ===================================================================
109  qq=paa(4)**2-paa(3)**2-paa(2)**2-paa(1)**2
110 
111 C 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 
116 C --- loop over five contributions of the rho-pi-pi
117  DO 201 k1=1,3
118  DO 201 k2=3,4
119 C
120  IF (k2.EQ.k1) THEN
121  GOTO 201
122  ELSEIF (k2.EQ.3) THEN
123 C rho-
124  ampr = ampl(3)
125  ampa = ampiz
126  ELSEIF (k1.EQ.3) THEN
127 C rho+
128  ampr = ampl(4)
129  ampa = ampiz
130  ELSE
131 C rho0
132  ampr = ampl(2)
133  ampa = ampi
134  END IF
135 C
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 
139 C -- definition of AA matrix
140 C -- 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
145 C ... 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
159 C
160 C --- lets add something to HADCURR
161 C FORM1= FPIKM(SQRT(SK),AMPI,AMPI) *FPIKM(SQRT(QQ),AMPI,AMPI)
162 C 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
168 C
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
173 C --- end of the rho-pi-pi current (5 possibilities)
174  201 CONTINUE
175 C
176 C ===================================================================
177 C Now modify the coefficient for the omega-pi current: =
178 C ===================================================================
179  IF (ampl(5).EQ.cmplx(0.,0.)) GOTO 311
180 
181 C Overall rho+rhoprime for the 4pi system:
182 C FORM2=AMPL(5)*(BWIGN(QQ,AMRO,GAMRO)+ELPHA*BWIGN(QQ,AMROP,GAMROP))
183 C 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))
187 C
188 C --- there are two possibilities for omega current
189 C --- 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
195 C --- 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
220 C
221 C omega -> rho pi for the 3pi system:
222 C FORM3=BWIGN(QQA,AMOM,GAMOM)*(BWIGN(SS23,AMRO,GAMRO)+
223 C $ BWIGN(SS24,AMRO,GAMRO)+BWIGN(SS34,AMRO,GAMRO))
224 C No omega -> rho pi; just straight omega:
225  form3=bwign(qqa,amom,gamom)
226 C
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
235 C
236  ELSE
237 C ===================================================================
238 C pi0 pi0 p0 pi- case ====
239 C ===================================================================
240  qq=paa(4)**2-paa(3)**2-paa(2)**2-paa(1)**2
241 
242 C --- 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 
247 C -- definition of AA matrix
248 C -- 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
253 C
254 C ... 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 
269 C --- lets add something to HADCURR
270 C FORM1= FPIKM(SQRT(SK),AMPI,AMPI) *FPIKMD(SQRT(QQ),AMPI,AMPI)
271 CCCCCCCCCCCCC FORM1=WIGFOR(SK,AMRO,GAMRO) (tests)
272 C 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
279 C --- end of the non omega current (3 possibilities)
280  101 CONTINUE
281 
282  ENDIF
283  END
284 
285 
286