8 IMPLICIT DOUBLE PRECISION(d)
9 common/lujets/
n,
k(9000,5),
p(9000,5),
v(9000,5)
11 common/ludat1/mstu(200),paru(200),mstj(200),parj(200)
13 common/ludat2/kchg(500,3),pmas(500,4),parf(2000),vckm(4,4)
15 common/ludat3/mdcy(500,3),mdme(2000,2),brat(2000),kfdp(2000,5)
17 dimension dps(5),dpc(5),ue(3)
23 IF(
k(
i,1).NE.3) goto 120
27 IF(kq.EQ.0.OR.(mqgst.EQ.1.AND.kq.EQ.2)) goto 120
31 IF(kq*isign(1,
k(
i,2)).LT.0) kcs=5
36 CALL
luerrm(14,
'(LUPREP:) caught in infinite loop')
42 IF(
i1.GE.mstu(4)-mstu(32)-5)
THEN
43 CALL
luerrm(11,
'(LUPREP:) no more memory left in LUJETS')
48 IF(nstp.GE.2.AND.iabs(
k(ia,2)).NE.21)
k(
i1,1)=1
57 IF(
k(
i1,1).EQ.1) goto 120
62 IF(mod(
k(ib,kcs)/mstu(5)**2,2).EQ.0.AND.mod(
k(ib,kcs),mstu(5)).
64 ia=mod(
k(ib,kcs),mstu(5))
65 k(ib,kcs)=
k(ib,kcs)+mstu(5)**2
68 IF(
k(ib,kcs).GE.2*mstu(5)**2.OR.mod(
k(ib,kcs)/mstu(5),mstu(5)).
70 ia=mod(
k(ib,kcs)/mstu(5),mstu(5))
71 k(ib,kcs)=
k(ib,kcs)+2*mstu(5)**2
74 IF(ia.LE.0.OR.ia.GT.
n)
THEN
75 CALL
luerrm(12,
'(LUPREP:) colour rearrangement failed')
78 IF(mod(
k(ia,4)/mstu(5),mstu(5)).EQ.ib.OR.mod(
k(ia,5)/mstu(5),
80 IF(mrev.EQ.1) kcs=9-kcs
81 IF(mod(
k(ia,kcs)/mstu(5),mstu(5)).NE.ib) kcs=9-kcs
82 k(ia,kcs)=
k(ia,kcs)+2*mstu(5)**2
84 IF(mrev.EQ.0) kcs=9-kcs
85 IF(mod(
k(ia,kcs),mstu(5)).NE.ib) kcs=9-kcs
86 k(ia,kcs)=
k(ia,kcs)+mstu(5)**2
95 IF(mstj(14).LE.0) goto 320
100 DO 190
i=
max(1,ip),ns
101 IF(
k(
i,1).NE.1.AND.
k(
i,1).NE.2)
THEN
102 ELSEIF(
k(
i,1).EQ.2.AND.ic.EQ.0)
THEN
109 ELSEIF(
k(
i,1).EQ.2)
THEN
111 160 dps(
j)=dps(
j)+
p(
i,
j)
112 ELSEIF(ic.NE.0.AND.kchg(
lucomp(
k(
i,2)),2).NE.0)
THEN
114 170 dps(
j)=dps(
j)+
p(
i,
j)
117 pd=sqrt(
max(0d0,dps(4)**2-dps(1)**2-dps(2)**2-dps(3)**2))-dps(5)
130 IF(pdm.GE.parj(32)) goto 320
134 pecm=sqrt(
max(0d0,dpc(4)**2-dpc(1)**2-dpc(2)**2-dpc(3)**2))
149 IF(mstu(16).NE.2)
THEN
160 IF(iabs(
k(ic1,2)).NE.21)
THEN
163 IF(kc1.EQ.0.OR.kc2.EQ.0) goto 320
164 kq1=kchg(kc1,2)*isign(1,
k(ic1,2))
165 kq2=kchg(kc2,2)*isign(1,
k(ic2,2))
166 IF(kq1+kq2.NE.0) goto 320
167 200 CALL
lukfdi(
k(ic1,2),0,kfln,
k(
n+2,2))
168 CALL
lukfdi(
k(ic2,2),-kfln,kfldmp,
k(
n+3,2))
169 IF(
k(
n+2,2).EQ.0.OR.
k(
n+3,2).EQ.0) goto 200
171 IF(iabs(
k(ic2,2)).NE.21) goto 320
172 210 CALL
lukfdi(1+int((2.+parj(2))*
rlu(0)),0,kfln,kfdmp)
174 CALL
lukfdi(-kfln,-kflm,kfldmp,
k(
n+3,2))
175 IF(
k(
n+2,2).EQ.0.OR.
k(
n+3,2).EQ.0) goto 210
179 IF(
p(
n+2,5)+
p(
n+3,5)+parj(64).GE.pecm.AND.nsin.EQ.1) goto 320
180 IF(
p(
n+2,5)+
p(
n+3,5)+parj(64).GE.pecm) goto 260
183 IF(pecm.GE.0.02*dpc(4))
THEN
184 pa=sqrt((pecm**2-(
p(
n+2,5)+
p(
n+3,5))**2)*(pecm**2-
185 & (
p(
n+2,5)-
p(
n+3,5))**2))/(2.*pecm)
188 ue(1)=sqrt(1.-ue(3)**2)*cos(
phi)
189 ue(2)=sqrt(1.-ue(3)**2)*sin(
phi)
192 220
p(
n+3,
j)=-pa*ue(
j)
193 p(
n+2,4)=sqrt(pa**2+
p(
n+2,5)**2)
194 p(
n+3,4)=sqrt(pa**2+
p(
n+3,5)**2)
195 CALL ludbrb(
n+2,
n+3,0.,0.,dpc(1)/dpc(4),dpc(2)/dpc(4),
200 230
IF(
k(
i,1).EQ.1.OR.
k(
i,1).EQ.2)
np=
np+1
201 ha=
p(ic1,4)*
p(ic2,4)-
p(ic1,1)*
p(ic2,1)-
p(ic1,2)*
p(ic2,2)-
203 IF(
np.GE.3.OR.ha.LE.1.25*
p(ic1,5)*
p(ic2,5)) goto 260
204 hd1=0.5*(
p(
n+2,5)**2-
p(ic1,5)**2)
205 hd2=0.5*(
p(
n+3,5)**2-
p(ic2,5)**2)
206 hr=sqrt(
max(0.,((ha-hd1-hd2)**2-(
p(
n+2,5)*
p(
n+3,5))**2)/
207 & (ha**2-(
p(ic1,5)*
p(ic2,5))**2)))-1.
208 hc=
p(ic1,5)**2+2.*ha+
p(ic2,5)**2
209 hk1=((
p(ic2,5)**2+ha)*hr+hd1-hd2)/hc
210 hk2=((
p(ic1,5)**2+ha)*hr+hd2-hd1)/hc
212 p(
n+2,
j)=(1.+hk1)*
p(ic1,
j)-hk2*
p(ic2,
j)
213 240
p(
n+3,
j)=(1.+hk2)*
p(ic2,
j)-hk1*
p(ic1,
j)
218 250
v(
n+3,
j)=
v(ic2,
j)
227 IF(iabs(
k(ic1,2)).GT.100.AND.iabs(
k(ic2,2)).GT.100)
THEN
229 ELSEIF(iabs(
k(ic1,2)).NE.21)
THEN
230 CALL
lukfdi(
k(ic1,2),
k(ic2,2),kfldmp,
k(
n+2,2))
232 kfln=1+int((2.+parj(2))*
rlu(0))
233 CALL
lukfdi(kfln,-kfln,kfldmp,
k(
n+2,2))
235 IF(
k(
n+2,2).EQ.0) goto 260
244 IF(
k(
i,1).LE.0.OR.
k(
i,1).GT.10.OR.(
i.GE.ic1.AND.
i.LE.ic2.
245 &and.
k(
i,1).GE.1.AND.
k(
i,1).LE.2)) goto 270
247 IF(mcomb.EQ.1.AND.kci.EQ.0) goto 270
248 IF(mcomb.EQ.1.AND.kchg(kci,2).EQ.0.AND.
i.LE.ns) goto 270
249 IF(mcomb.EQ.2.AND.iabs(
k(
i,2)).GT.10.AND.iabs(
k(
i,2)).LE.100)
251 hcr=dpc(4)*
p(
i,4)-dpc(1)*
p(
i,1)-dpc(2)*
p(
i,2)-dpc(3)*
p(
i,3)
265 IF(ha**2-(pecm*
p(ir,5))**2.EQ.0.0.OR.hb+hd.EQ.0.0)
go to 285
267 hk2=0.5*(hb*sqrt(((hb+hc)**2-4.*(hb+hd)*
p(
n+2,5)**2)/
268 &(ha**2-(pecm*
p(ir,5))**2))-(hb+hc))/(hb+hd)
269 285 hk1=(0.5*(
p(
n+2,5)**2-pecm**2)+hd*hk2)/hb
271 p(
n+2,
j)=(1.+hk1)*dpc(
j)-hk2*
p(ir,
j)
272 p(ir,
j)=(1.+hk2)*
p(ir,
j)-hk1*dpc(
j)
274 290
v(
n+2,
j)=
v(ic1,
j)
281 IF((
k(
i,1).EQ.1.OR.
k(
i,1).EQ.2).AND.kchg(
lucomp(
k(
i,2)),2).NE.0)
284 IF(mstu(16).NE.2)
THEN
293 IF(
n.LT.mstu(4)-mstu(32)-5) goto 140
302 IF(
k(
i,1).LE.0.OR.
k(
i,1).GT.10) goto 360
305 kq=kchg(kc,2)*isign(1,
k(
i,2))
315 340 dps(
j)=dps(
j)+
p(
i,
j)
317 IF(
np.NE.1.AND.(kfn.EQ.1.OR.kfn.GE.3.OR.kqs.NE.0)) CALL
318 &
luerrm(2,
'(LUPREP:) unphysical flavour combination')
319 IF(
np.NE.1.AND.dps(4)**2-dps(1)**2-dps(2)**2-dps(3)**2.LT.
320 & (0.9*parj(32)+dps(5))**2) CALL
luerrm(3,
321 &
'(LUPREP:) too small mass in jet system')