7 IMPLICIT DOUBLE PRECISION(d)
8 common/lujets/
n,
k(9000,5),
p(9000,5),
v(9000,5)
10 common/ludat1/mstu(200),paru(200),mstj(200),parj(200)
12 common/ludat2/kchg(500,3),pmas(500,4),parf(2000),vckm(4,4)
14 dimension pmth(5,40),
ps(5),pma(4),pmsd(4),iep(4),ipa(4),
15 &kfla(4),kfld(4),kfl(4),itry(4),isi(4),isl(4),dp(4),dpt(5,4)
18 IF(mstj(41).LE.0.OR.(mstj(41).EQ.1.AND.qmax.LE.parj(82)).OR.
19 &qmax.LE.min(parj(82),parj(83)).OR.mstj(41).GE.3)
RETURN
21 pmth(2,21)=sqrt(pmth(1,21)**2+0.25*parj(82)**2)
22 pmth(3,21)=2.*pmth(2,21)
26 pmth(2,22)=sqrt(pmth(1,22)**2+0.25*parj(83)**2)
27 pmth(3,22)=2.*pmth(2,22)
31 IF(mstj(41).EQ.2) pmqth1=min(parj(82),parj(83))
33 IF(mstj(41).EQ.2) pmqth2=min(pmth(2,21),pmth(2,22))
36 pmth(2,
if)=sqrt(pmth(1,
if)**2+0.25*pmqth1**2)
37 pmth(3,
if)=pmth(2,
if)+pmqth2
38 pmth(4,
if)=sqrt(pmth(1,
if)**2+0.25*parj(82)**2)+pmth(2,21)
39 100 pmth(5,
if)=sqrt(pmth(1,
if)**2+0.25*parj(83)**2)+pmth(2,22)
40 pt2min=
max(0.5*parj(82),1.1*parj(81))**2
42 alfm=
log(pt2min/alams)
46 IF(ip1.GT.0.AND.ip1.LE.min(
n,mstu(4)-mstu(32)).AND.ip2.EQ.0)
THEN
49 ELSEIF(min(ip1,ip2).GT.0.AND.
max(ip1,ip2).LE.min(
n,mstu(4)-
54 ELSEIF(ip1.GT.0.AND.ip1.LE.min(
n,mstu(4)-mstu(32)).AND.ip2.LT.0.
61 &
'(LUSHOW:) failed to reconstruct showering system')
62 IF(mstu(21).GE.1)
RETURN
71 kfla(
i)=iabs(
k(ipa(
i),2))
73 IF(kfla(
i).NE.0.AND.(kfla(
i).LE.8.OR.kfla(
i).EQ.21))
74 &pma(
i)=pmth(3,kfla(
i))
76 IF(kfla(
i).EQ.0.OR.(kfla(
i).GT.8.AND.kfla(
i).NE.21).OR.
77 &pma(
i).GT.qmax) irej=irej+1
80 IF(irej.EQ.npa)
RETURN
82 IF(npa.EQ.1)
ps(5)=
ps(4)
83 IF(
ps(5).LE.pm+pmqth1)
RETURN
84 IF(npa.EQ.2.AND.mstj(47).GE.1)
THEN
85 IF(kfla(1).GE.1.AND.kfla(1).LE.8.AND.kfla(2).GE.1.AND.
86 & kfla(2).LE.8) m3jc=1
87 IF(mstj(47).GE.2) m3jc=1
92 IF(
n.GT.mstu(4)-mstu(32)-5)
THEN
93 CALL
luerrm(11,
'(LUSHOW:) no more memory left in LUJETS')
94 IF(mstu(21).GE.1)
RETURN
119 IF(kflm.EQ.0.OR.(kflm.GT.8.AND.kflm.NE.21)) goto 140
120 IF(
p(im,5).LT.pmth(2,kflm)) goto 140
125 IF(
n+nep.GT.mstu(4)-mstu(32)-5)
THEN
126 CALL
luerrm(11,
'(LUSHOW:) no more memory left in LUJETS')
127 IF(mstu(21).GE.1)
RETURN
134 IF(
k(im-1,3).EQ.igm) iau=im-1
135 IF(
n.GE.im+1.AND.
k(im+1,3).EQ.igm) iau=im+1
146 160
k(
n+
i,2)=
k(ipa(
i),2)
147 ELSEIF(kflm.NE.21)
THEN
150 ELSEIF(
k(im,5).EQ.21)
THEN
163 kfld(ip)=iabs(
k(
n+ip,2))
167 170
IF(kfld(ip).GT.0.AND.(kfld(ip).LE.8.OR.kfld(ip).EQ.21)) isi(ip)=1
173 IF(npa.GE.3)
p(
n+
i,4)=(
ps(4)*
p(ipa(
i),4)-
ps(1)*
p(ipa(
i),1)-
175 p(
n+
i,5)=min(qmax,
ps(5))
176 IF(npa.GE.3)
p(
n+
i,5)=min(
p(
n+
i,5),
p(
n+
i,4))
177 180
IF(isi(
i).EQ.0)
p(
n+
i,5)=
p(ipa(
i),5)
179 IF(mstj(43).LE.2) pem=
v(im,2)
180 IF(mstj(43).GE.3) pem=
p(im,4)
181 p(
n+1,5)=min(
p(im,5),
v(im,1)*pem)
182 p(
n+2,5)=min(
p(im,5),(1.-
v(im,1))*pem)
183 IF(
k(
n+2,2).EQ.22)
p(
n+2,5)=pmth(1,22)
188 IF(
p(
n+
i,5).LE.pmth(3,kfld(
i)))
p(
n+
i,5)=pmth(1,kfld(
i))
196 210
IF(inum.EQ.0.AND.isl(
i).EQ.1) inum=
i
198 IF(inum.EQ.0.AND.itry(
i).EQ.0.AND.isi(
i).EQ.1)
THEN
199 IF(
p(
n+
i,5).GE.pmth(2,kfld(
i))) inum=
i
205 IF(isi(
i).EQ.1.AND.pmsd(
i).GE.pmqth2)
THEN
207 IF(rpm.GT.rmax.AND.
p(
n+
i,5).GE.pmth(2,kfld(
i)))
THEN
220 240
IF(iep(
i).GT.
n+nep) iep(
i)=
n+1
222 250 kfl(
i)=iabs(
k(iep(
i),2))
223 itry(inum)=itry(inum)+1
224 IF(itry(inum).GT.200)
THEN
225 CALL
luerrm(14,
'(LUSHOW:) caught in infinite loop')
226 IF(mstu(21).GE.1)
RETURN
229 IF(kfl(1).EQ.0.OR.(kfl(1).GT.8.AND.kfl(1).NE.21)) goto 300
230 IF(
p(iep(1),5).LT.pmth(2,kfl(1))) goto 300
235 ELSEIF(igm.EQ.0.OR.mstj(43).LE.2)
THEN
238 IF(inum.EQ.1) pmed=
v(im,1)*pem
239 IF(inum.EQ.2) pmed=(1.-
v(im,1))*pem
241 IF(mod(mstj(43),2).EQ.1)
THEN
245 zc=0.5*(1.-sqrt(
max(0.,1.-(2.*pmth(2,21)/pmed)**2)))
246 IF(zc.LT.1
e-4) zc=(pmth(2,21)/pmed)**2
247 zce=0.5*(1.-sqrt(
max(0.,1.-(2.*pmth(2,22)/pmed)**2)))
248 IF(zce.LT.1
e-4) zce=(pmth(2,22)/pmed)**2
252 IF((mstj(41).EQ.1.AND.zc.GT.0.49).OR.(mstj(41).EQ.2.AND.
253 &min(zc,zce).GT.0.49))
THEN
254 p(iep(1),5)=pmth(1,kfl(1))
255 v(iep(1),5)=
p(iep(1),5)**2
260 IF(mstj(49).EQ.0.AND.kfl(1).EQ.21)
THEN
261 fbr=6.*
log((1.-zc)/zc)+mstj(45)*(0.5-zc)
262 ELSEIF(mstj(49).EQ.0)
THEN
263 fbr=(8./3.)*
log((1.-zc)/zc)
266 ELSEIF(mstj(49).EQ.1.AND.kfl(1).EQ.21)
THEN
267 fbr=(parj(87)+mstj(45)*parj(88))*(1.-2.*zc)
268 ELSEIF(mstj(49).EQ.1)
THEN
270 IF(igm.EQ.0.AND.m3jc.EQ.1) fbr=4.*fbr
273 ELSEIF(kfl(1).EQ.21)
THEN
274 fbr=6.*mstj(45)*(0.5-zc)
276 fbr=2.*
log((1.-zc)/zc)
280 IF(mstj(41).EQ.2.AND.kfl(1).GE.1.AND.kfl(1).LE.8)
281 &fbre=(kchg(kfl(1),1)/3.)**2*2.*
log((1.-zce)/zce)
289 IF(kfl(
i).GT.0.AND.(kfl(
i).LE.8.OR.kfl(
i).EQ.21)) pm=
292 pms=min(pms,(
p(im,5)-pm2)**2)
298 280
IF(pms.GT.4.*pmth(2,
if)**2) b0=(33.-2.*
if)/6.
299 IF(mstj(44).LE.0)
THEN
300 pmsqcd=pms*exp(
max(-100.,
log(
rlu(0))*paru(2)/(paru(111)*fbr)))
301 ELSEIF(mstj(44).EQ.1)
THEN
302 pmsqcd=4.*alams*(0.25*pms/alams)**(
rlu(0)**(b0/fbr))
304 pmsqcd=pms*
rlu(0)**(alfm*b0/fbr)
306 IF(zc.GT.0.49.OR.pmsqcd.LE.pmth(4,kfl(1))**2) pmsqcd=
312 IF(mstj(41).EQ.2.AND.kfl(1).GE.1.AND.kfl(1).LE.8)
THEN
313 pmsqed=pms*exp(
max(-100.,
log(
rlu(0))*paru(2)/(paru(101)*fbre)))
314 IF(zce.GT.0.49.OR.pmsqed.LE.pmth(5,kfl(1))**2) pmsqed=
316 IF(pmsqed.GT.pmsqcd)
THEN
323 p(iep(1),5)=sqrt(
v(iep(1),5))
324 IF(
p(iep(1),5).LE.pmth(3,kfl(1)))
THEN
325 p(iep(1),5)=pmth(1,kfl(1))
326 v(iep(1),5)=
p(iep(1),5)**2
332 z=1.-(1.-zce)*(zce/(1.-zce))**
rlu(0)
333 IF(1.+
z**2.LT.2.*
rlu(0)) goto 260
337 ELSEIF(mstj(49).NE.1.AND.kfl(1).NE.21)
THEN
338 z=1.-(1.-zc)*(zc/(1.-zc))**
rlu(0)
339 IF(1.+
z**2.LT.2.*
rlu(0)) goto 260
341 ELSEIF(mstj(49).EQ.0.AND.mstj(45)*(0.5-zc).LT.
rlu(0)*fbr)
THEN
342 z=(1.-zc)*(zc/(1.-zc))**
rlu(0)
343 IF(
rlu(0).GT.0.5)
z=1.-
z
344 IF((1.-
z*(1.-
z))**2.LT.
rlu(0)) goto 260
346 ELSEIF(mstj(49).NE.1)
THEN
347 z=zc+(1.-2.*zc)*
rlu(0)
348 IF(
z**2+(1.-
z)**2.LT.
rlu(0)) goto 260
349 kflb=1+int(mstj(45)*
rlu(0))
350 pmq=4.*pmth(2,kflb)**2/
v(iep(1),5)
351 IF(pmq.GE.1.) goto 260
352 pmq0=4.*pmth(2,21)**2/
v(iep(1),5)
353 IF(mod(mstj(43),2).EQ.0.AND.(1.+0.5*pmq)*sqrt(1.-pmq).LT.
354 &
rlu(0)*(1.+0.5*pmq0)*sqrt(1.-pmq0)) goto 260
358 ELSEIF(kfl(1).NE.21)
THEN
359 z=1.-sqrt(zc**2+
rlu(0)*(1.-2.*zc))
361 ELSEIF(
rlu(0)*(parj(87)+mstj(45)*parj(88)).LE.parj(87))
THEN
362 z=zc+(1.-2.*zc)*
rlu(0)
365 z=zc+(1.-2.*zc)*
rlu(0)
366 kflb=1+int(mstj(45)*
rlu(0))
367 pmq=4.*pmth(2,kflb)**2/
v(iep(1),5)
368 IF(pmq.GE.1.) goto 260
371 IF(mce.EQ.1.AND.mstj(44).GE.2)
THEN
372 IF(
z*(1.-
z)*
v(iep(1),5).LT.pt2min) goto 260
373 IF(alfm/
log(
v(iep(1),5)*
z*(1.-
z)/alams).LT.
rlu(0)) goto 260
377 IF(kfl(1).EQ.21)
THEN
378 kflgd1=iabs(
k(iep(1),5))
382 kflgd2=iabs(
k(iep(1),5))
386 ELSEIF(nep.GE.3)
THEN
388 ELSEIF(igm.EQ.0.OR.mstj(43).LE.2)
THEN
389 ped=0.5*(
v(im,5)+
v(iep(1),5)-pm2**2)/
p(im,5)
391 IF(iep(1).EQ.
n+1) ped=
v(im,1)*pem
392 IF(iep(1).EQ.
n+2) ped=(1.-
v(im,1))*pem
394 IF(mod(mstj(43),2).EQ.1)
THEN
396 IF(kflgd2.EQ.22) pmqth3=0.5*parj(83)
397 pmq1=(pmth(1,kflgd1)**2+pmqth3**2)/
v(iep(1),5)
398 pmq2=(pmth(1,kflgd2)**2+pmqth3**2)/
v(iep(1),5)
399 zd=sqrt(
max(0.,(1.-
v(iep(1),5)/ped**2)*((1.-pmq1-pmq2)**2-
403 zd=sqrt(
max(0.,1.-
v(iep(1),5)/ped**2))
408 IF(
z.LT.zl.OR.
z.GT.zu) goto 260
409 IF(kfl(1).EQ.21)
v(iep(1),3)=
log(zu*(1.-zl)/
max(1
e-20,zl*
411 IF(kfl(1).NE.21)
v(iep(1),3)=
log((1.-zl)/
max(1
e-10,1.-zu))
414 IF(igm.EQ.0.AND.m3jc.EQ.1)
THEN
415 x1=
z*(1.+
v(iep(1),5)/
v(ns+1,5))
416 x2=1.-
v(iep(1),5)/
v(ns+1,5)
421 qf1=kchg(iabs(ki1),1)*isign(1,ki1)/3.
422 qf2=kchg(iabs(ki2),1)*isign(1,ki2)/3.
423 wshow=qf1**2*(1.-x1)/x3*(1.+(x1/(2.-x2))**2)+
424 & qf2**2*(1.-x2)/x3*(1.+(x2/(2.-x1))**2)
425 wme=(qf1*(1.-x1)/x3-qf2*(1.-x2)/x3)**2*(x1**2+x2**2)
426 ELSEIF(mstj(49).NE.1)
THEN
427 wshow=1.+(1.-x1)/x3*(x1/(2.-x2))**2+
428 & (1.-x2)/x3*(x2/(2.-x1))**2
431 wshow=4.*x3*((1.-x1)/(2.-x2)**2+(1.-x2)/(2.-x1)**2)
434 IF(wme.LT.
rlu(0)*wshow) goto 260
437 ELSEIF(mce.EQ.1.AND.igm.GT.0.AND.mstj(42).GE.2)
THEN
440 IF(iep(1).EQ.
n+2) zm=1.-
v(im,1)
441 the2id=
z*(1.-
z)*(zm*
p(im,4))**2/
v(iep(1),5)
443 290
IF(
k(iaom,5).EQ.22)
THEN
445 IF(
k(iaom,3).LE.ns) maom=0
446 IF(maom.EQ.1) goto 290
449 the2im=
v(iaom,1)*(1.-
v(iaom,1))*
p(iaom,4)**2/
v(iaom,5)
450 IF(the2id.LT.the2im) goto 260
455 IF(mstj(48).EQ.1)
THEN
456 IF(nep.EQ.1.AND.im.EQ.ns)
THEN
457 the2id=
z*(1.-
z)*
ps(4)**2/
v(iep(1),5)
458 IF(the2id.LT.1./parj(85)**2) goto 260
459 ELSEIF(nep.EQ.2.AND.iep(1).EQ.ns+2)
THEN
460 the2id=
z*(1.-
z)*(0.5*
p(im,4))**2/
v(iep(1),5)
461 IF(the2id.LT.1./parj(85)**2) goto 260
462 ELSEIF(nep.EQ.2.AND.iep(1).EQ.ns+3)
THEN
463 the2id=
z*(1.-
z)*(0.5*
p(im,4))**2/
v(iep(1),5)
464 IF(the2id.LT.1./parj(86)**2) goto 260
472 IF(nep.EQ.1) goto 330
473 IF(nep.EQ.2.AND.
p(iep(1),5)+
p(iep(2),5).GE.
p(im,5)) goto 200
475 IF(itry(
i).EQ.0.AND.kfld(
i).GT.0.AND.(kfld(
i).LE.8.OR.kfld(
i).EQ.
477 IF(
p(
n+
i,5).GE.pmth(2,kfld(
i))) goto 200
483 pa1s=(
p(
n+1,4)+
p(
n+1,5))*(
p(
n+1,4)-
p(
n+1,5))
484 pa2s=(
p(
n+2,4)+
p(
n+2,5))*(
p(
n+2,4)-
p(
n+2,5))
485 pa3s=(
p(
n+3,4)+
p(
n+3,5))*(
p(
n+3,4)-
p(
n+3,5))
486 pts=0.25*(2.*pa1s*pa2s+2.*pa1s*pa3s+2.*pa2s*pa3s-
487 & pa1s**2-pa2s**2-pa3s**2)/pa1s
488 IF(pts.LE.0.) goto 200
489 ELSEIF(igm.EQ.0.OR.mstj(43).LE.2.OR.mod(mstj(43),2).EQ.0)
THEN
492 IF(kflda.EQ.0.OR.(kflda.GT.8.AND.kflda.NE.21)) goto 320
493 IF(
p(
i1,5).LT.pmth(2,kflda)) goto 320
502 IF(igm.EQ.0.OR.mstj(43).LE.2)
THEN
503 ped=0.5*(
v(im,5)+
v(
i1,5)-
v(
i2,5))/
p(im,5)
505 IF(
i1.EQ.
n+1) zm=
v(im,1)
506 IF(
i1.EQ.
n+2) zm=1.-
v(im,1)
507 pml=sqrt((
v(im,5)-
v(
n+1,5)-
v(
n+2,5))**2-
508 & 4.*
v(
n+1,5)*
v(
n+2,5))
509 ped=pem*(0.5*(
v(im,5)-pml+
v(
i1,5)-
v(
i2,5))+pml*zm)/
v(im,5)
511 IF(mod(mstj(43),2).EQ.1)
THEN
513 IF(kflgd2.EQ.22) pmqth3=0.5*parj(83)
514 pmq1=(pmth(1,kflgd1)**2+pmqth3**2)/
v(
i1,5)
515 pmq2=(pmth(1,kflgd2)**2+pmqth3**2)/
v(
i1,5)
516 zd=sqrt(
max(0.,(1.-
v(
i1,5)/ped**2)*((1.-pmq1-pmq2)**2-
520 zd=sqrt(
max(0.,1.-
v(
i1,5)/ped**2))
525 IF(
i1.EQ.
n+1.AND.(
v(
i1,1).LT.zl.OR.
v(
i1,1).GT.zu)) isl(1)=1
526 IF(
i1.EQ.
n+2.AND.(
v(
i1,1).LT.zl.OR.
v(
i1,1).GT.zu)) isl(2)=1
527 IF(kflda.EQ.21)
v(
i1,4)=
log(zu*(1.-zl)/
max(1
e-20,zl*(1.-zu)))
528 IF(kflda.NE.21)
v(
i1,4)=
log((1.-zl)/
max(1
e-10,1.-zu))
530 IF(isl(1).EQ.1.AND.isl(2).EQ.1.AND.islm.NE.0)
THEN
533 ELSEIF(isl(1).EQ.1.AND.isl(2).EQ.1)
THEN
534 zdr1=
max(0.,
v(
n+1,3)/
v(
n+1,4)-1.)
535 zdr2=
max(0.,
v(
n+2,3)/
v(
n+2,4)-1.)
536 IF(zdr2.GT.
rlu(0)*(zdr1+zdr2)) isl(1)=0
537 IF(isl(1).EQ.1) isl(2)=0
538 IF(isl(1).EQ.0) islm=1
539 IF(isl(2).EQ.0) islm=2
541 IF(isl(1).EQ.1.OR.isl(2).EQ.1) goto 200
543 IF(igm.GT.0.AND.mod(mstj(43),2).EQ.1.AND.(
p(
n+1,5).GE.
544 &pmth(2,kfld(1)).OR.
p(
n+2,5).GE.pmth(2,kfld(2))))
THEN
545 pmq1=
v(
n+1,5)/
v(im,5)
546 pmq2=
v(
n+2,5)/
v(im,5)
547 zd=sqrt(
max(0.,(1.-
v(im,5)/pem**2)*((1.-pmq1-pmq2)**2-
552 IF(
v(im,1).LT.zl.OR.
v(im,1).GT.zu) goto 200
561 p(
n+1,3)=sqrt(
max(0.,(
p(ipa(1),4)+
p(
n+1,5))*(
p(ipa(1),4)-
565 ELSEIF(igm.EQ.0.AND.nep.EQ.2)
THEN
566 ped1=0.5*(
v(im,5)+
v(
n+1,5)-
v(
n+2,5))/
p(im,5)
569 p(
n+1,3)=sqrt(
max(0.,(ped1+
p(
n+1,5))*(ped1-
p(
n+1,5))))
574 p(
n+2,4)=
p(im,5)-ped1
577 ELSEIF(nep.EQ.3)
THEN
580 p(
n+1,3)=sqrt(
max(0.,pa1s))
583 p(
n+2,3)=0.5*(pa3s-pa2s-pa1s)/
p(
n+1,3)
586 p(
n+3,3)=-(
p(
n+1,3)+
p(
n+2,3))
594 pzm=sqrt(
max(0.,(pem+
p(im,5))*(pem-
p(im,5))))
595 pmls=(
v(im,5)-
v(
n+1,5)-
v(
n+2,5))**2-4.*
v(
n+1,5)*
v(
n+2,5)
598 ELSEIF(mod(mstj(43),2).EQ.1)
THEN
599 pts=(pem**2*(zm*(1.-zm)*
v(im,5)-(1.-zm)*
v(
n+1,5)-
600 & zm*
v(
n+2,5))-0.25*pmls)/pzm**2
602 pts=pmls*(zm*(1.-zm)*pem**2/
v(im,5)-0.25)/pzm**2
608 IF(mstj(49).NE.1.AND.mod(mstj(46),2).EQ.1.AND.
k(im,2).EQ.21.
610 IF(
k(igm,3).NE.0) mazip=1
612 IF(iau.EQ.im+1) zau=1.-
v(igm,1)
613 IF(mazip.EQ.0) zau=0.
614 IF(
k(igm,2).NE.21)
THEN
615 hazip=2.*zau/(1.+zau**2)
617 hazip=(zau/(1.-zau*(1.-zau)))**2
619 IF(
k(
n+1,2).NE.21)
THEN
620 hazip=hazip*(-2.*zm*(1.-zm))/(1.-2.*zm*(1.-zm))
622 hazip=hazip*(zm*(1.-zm)/(1.-zm*(1.-zm)))**2
629 IF(mstj(46).GE.2.AND.(
k(
n+1,2).EQ.21.OR.
k(
n+2,2).EQ.21).
631 IF(
k(igm,3).NE.0) mazic=
n+1
632 IF(
k(igm,3).NE.0.AND.
k(
n+1,2).NE.21) mazic=
n+2
633 IF(
k(igm,3).NE.0.AND.
k(
n+1,2).EQ.21.AND.
k(
n+2,2).EQ.21.AND.
634 & zm.GT.0.5) mazic=
n+2
635 IF(
k(iau,2).EQ.22) mazic=0
637 IF(mazic.EQ.
n+2) zs=1.-zm
639 IF(iau.EQ.im-1) zgm=1.-
v(igm,1)
640 IF(mazic.EQ.0) zgm=1.
641 hazic=(
p(im,5)/
p(igm,5))*sqrt((1.-zs)*(1.-zgm)/(zs*zgm))
642 hazic=min(0.95,hazic)
647 340
IF(nep.EQ.2.AND.igm.GT.0)
THEN
648 IF(mod(mstj(43),2).EQ.1)
THEN
651 p(
n+1,4)=pem*(0.5*(
v(im,5)-sqrt(pmls)+
v(
n+1,5)-
v(
n+2,5))+
652 & sqrt(pmls)*zm)/
v(im,5)
658 p(
n+1,3)=0.5*(
v(
n+2,5)-
v(
n+1,5)-
v(im,5)+2.*pem*
p(
n+1,4))/pzm
664 p(
n+2,3)=pzm-
p(
n+1,3)
665 p(
n+2,4)=pem-
p(
n+1,4)
666 IF(mstj(43).LE.2)
THEN
667 v(
n+1,2)=(pem*
p(
n+1,4)-pzm*
p(
n+1,3))/
p(im,5)
668 v(
n+2,2)=(pem*
p(
n+2,4)-pzm*
p(
n+2,3))/
p(im,5)
674 IF(mstj(43).LE.2)
THEN
675 bex=
p(igm,1)/
p(igm,4)
676 bey=
p(igm,2)/
p(igm,4)
677 bez=
p(igm,3)/
p(igm,4)
679 gabep=ga*(ga*(bex*
p(im,1)+bey*
p(im,2)+bez*
p(im,3))/(1.+ga)-
688 the=
ulangl(
p(im,3)+gabep*bez,sqrt((
p(im,1)+gabep*bex)**2+
689 & (
p(im,2)+gabep*bey)**2))
692 dp(1)=cos(the)*cos(
phi)*
p(
i,1)-sin(
phi)*
p(
i,2)+
693 & sin(the)*cos(
phi)*
p(
i,3)
694 dp(2)=cos(the)*sin(
phi)*
p(
i,1)+cos(
phi)*
p(
i,2)+
695 & sin(the)*sin(
phi)*
p(
i,3)
696 dp(3)=-sin(the)*
p(
i,1)+cos(the)*
p(
i,3)
698 dbp=bex*dp(1)+bey*dp(2)+bez*dp(3)
699 dgabp=ga*(ga*dbp/(1d0+ga)+dp(4))
700 p(
i,1)=dp(1)+dgabp*bex
701 p(
i,2)=dp(2)+dgabp*bey
702 p(
i,3)=dp(3)+dgabp*bez
703 350
p(
i,4)=ga*(dp(4)+dbp)
707 IF(mazip.NE.0.OR.mazic.NE.0)
THEN
711 360 dpt(3,
j)=
p(
n+1,
j)
712 dpma=dpt(1,1)*dpt(2,1)+dpt(1,2)*dpt(2,2)+dpt(1,3)*dpt(2,3)
713 dpmd=dpt(1,1)*dpt(3,1)+dpt(1,2)*dpt(3,2)+dpt(1,3)*dpt(3,3)
714 dpmm=dpt(1,1)**2+dpt(1,2)**2+dpt(1,3)**2
716 dpt(4,
j)=dpt(2,
j)-dpma*dpt(1,
j)/dpmm
717 370 dpt(5,
j)=dpt(3,
j)-dpmd*dpt(1,
j)/dpmm
718 dpt(4,4)=sqrt(dpt(4,1)**2+dpt(4,2)**2+dpt(4,3)**2)
719 dpt(5,4)=sqrt(dpt(5,1)**2+dpt(5,2)**2+dpt(5,3)**2)
720 IF(min(dpt(4,4),dpt(5,4)).GT.0.1*parj(82))
THEN
721 cad=(dpt(4,1)*dpt(5,1)+dpt(4,2)*dpt(5,2)+
722 & dpt(4,3)*dpt(5,3))/(dpt(4,4)*dpt(5,4))
724 IF(1.+hazip*(2.*cad**2-1.).LT.
rlu(0)*(1.+abs(hazip)))
728 IF(mazic.EQ.
n+2) cad=-cad
729 IF((1.-hazic)*(1.-hazic*cad)/(1.+hazic**2-2.*hazic*cad).
730 & lt.
rlu(0)) goto 340
736 IF(igm.GE.0)
k(im,1)=14
739 IF(
n.GT.mstu(4)-mstu(32)-5)
THEN
740 CALL
luerrm(11,
'(LUSHOW:) no more memory left in LUJETS')
741 IF(mstu(21).GE.1)
n=ns
742 IF(mstu(21).GE.1)
RETURN
747 380
IF(npa.GE.2)
THEN
751 IF(ip2.GT.0.AND.ip2.LT.ip1)
k(ns+1,3)=ip2
761 IF(
k(
i,1).LE.10.AND.
k(
i,2).EQ.22)
THEN
763 ELSEIF(
k(
i,1).LE.10)
THEN
764 k(
i,4)=mstu(5)*(
k(
i,4)/mstu(5))
765 k(
i,5)=mstu(5)*(
k(
i,5)/mstu(5))
766 ELSEIF(
k(mod(
k(
i,4),mstu(5))+1,2).NE.22)
THEN
767 id1=mod(
k(
i,4),mstu(5))
768 IF(
k(
i,2).GE.1.AND.
k(
i,2).LE.8) id1=mod(
k(
i,4),mstu(5))+1
769 id2=2*mod(
k(
i,4),mstu(5))+1-id1
770 k(
i,4)=mstu(5)*(
k(
i,4)/mstu(5))+id1
771 k(
i,5)=mstu(5)*(
k(
i,5)/mstu(5))+id2
772 k(id1,4)=
k(id1,4)+mstu(5)*
i
773 k(id1,5)=
k(id1,5)+mstu(5)*id2
774 k(id2,4)=
k(id2,4)+mstu(5)*id1
775 k(id2,5)=
k(id2,5)+mstu(5)*
i
777 id1=mod(
k(
i,4),mstu(5))
779 k(
i,4)=mstu(5)*(
k(
i,4)/mstu(5))+id1
780 k(
i,5)=mstu(5)*(
k(
i,5)/mstu(5))+id1
781 k(id1,4)=
k(id1,4)+mstu(5)*
i
782 k(id1,5)=
k(id1,5)+mstu(5)*
i
794 gabep=ga*(ga*(bex*
p(ipa(1),1)+bey*
p(ipa(1),2)+bez*
p(ipa(1),3))
795 & /(1.+ga)-
p(ipa(1),4))
802 the=
ulangl(
p(ipa(1),3)+gabep*bez,sqrt((
p(ipa(1),1)
803 &+gabep*bex)**2+(
p(ipa(1),2)+gabep*bey)**2))
804 phi=
ulangl(
p(ipa(1),1)+gabep*bex,
p(ipa(1),2)+gabep*bey)
806 chi=
ulangl(cos(the)*cos(
phi)*(
p(ipa(2),1)+gabep*bex)+cos(the)*
807 & sin(
phi)*(
p(ipa(2),2)+gabep*bey)-sin(the)*(
p(ipa(2),3)+gabep*
808 & bez),-sin(
phi)*(
p(ipa(2),1)+gabep*bex)+cos(
phi)*(
p(ipa(2),2)+
810 CALL ludbrb(ns+1,
n,0.,chi,0d0,0d0,0d0)
815 CALL ludbrb(ns+1,
n,the,
phi,dbex,dbey,dbez)
823 IF(
n.EQ.ns+npa+iim)
THEN
828 k(ipa(ip),4)=
k(ipa(ip),4)+ns+iim+ip
829 k(ipa(ip),5)=
k(ipa(ip),5)+ns+iim+ip
830 k(ns+iim+ip,3)=ipa(ip)
831 IF(iim.EQ.1.AND.mstu(16).NE.2)
k(ns+iim+ip,3)=ns+1
832 k(ns+iim+ip,4)=mstu(5)*ipa(ip)+
k(ns+iim+ip,4)
833 410
k(ns+iim+ip,5)=mstu(5)*ipa(ip)+
k(ns+iim+ip,5)