]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | *CMZ : 17/07/98 15.44.33 by Federico Carminati |
2 | *-- Author : | |
3 | C********************************************************************* | |
4 | ||
5 | SUBROUTINE LUSHOW(IP1,IP2,QMAX) | |
6 | ||
7 | C...Purpose: to generate timelike parton showers from given partons. | |
8 | IMPLICIT DOUBLE PRECISION(D) | |
9 | *KEEP,LUJETS. | |
10 | COMMON /LUJETS/ N,K(200000,5),P(200000,5),V(200000,5) | |
11 | SAVE /LUJETS/ | |
12 | *KEEP,LUDAT1. | |
13 | COMMON /LUDAT1/ MSTU(200),PARU(200),MSTJ(200),PARJ(200) | |
14 | SAVE /LUDAT1/ | |
15 | *KEEP,LUDAT2. | |
16 | COMMON /LUDAT2/ KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4) | |
17 | SAVE /LUDAT2/ | |
18 | *KEND. | |
19 | DIMENSION PMTH(5,40),PS(5),PMA(4),PMSD(4),IEP(4),IPA(4), | |
20 | &KFLA(4),KFLD(4),KFL(4),ITRY(4),ISI(4),ISL(4),DP(4),DPT(5,4) | |
21 | ||
22 | C...Initialization of cutoff masses etc. | |
23 | IF(MSTJ(41).LE.0.OR.(MSTJ(41).EQ.1.AND.QMAX.LE.PARJ(82)).OR. | |
24 | &QMAX.LE.MIN(PARJ(82),PARJ(83)).OR.MSTJ(41).GE.3) RETURN | |
25 | PMTH(1,21)=ULMASS(21) | |
26 | PMTH(2,21)=SQRT(PMTH(1,21)**2+0.25*PARJ(82)**2) | |
27 | PMTH(3,21)=2.*PMTH(2,21) | |
28 | PMTH(4,21)=PMTH(3,21) | |
29 | PMTH(5,21)=PMTH(3,21) | |
30 | PMTH(1,22)=ULMASS(22) | |
31 | PMTH(2,22)=SQRT(PMTH(1,22)**2+0.25*PARJ(83)**2) | |
32 | PMTH(3,22)=2.*PMTH(2,22) | |
33 | PMTH(4,22)=PMTH(3,22) | |
34 | PMTH(5,22)=PMTH(3,22) | |
35 | PMQTH1=PARJ(82) | |
36 | IF(MSTJ(41).EQ.2) PMQTH1=MIN(PARJ(82),PARJ(83)) | |
37 | PMQTH2=PMTH(2,21) | |
38 | IF(MSTJ(41).EQ.2) PMQTH2=MIN(PMTH(2,21),PMTH(2,22)) | |
39 | DO 100 IF=1,8 | |
40 | PMTH(1,IF)=ULMASS(IF) | |
41 | PMTH(2,IF)=SQRT(PMTH(1,IF)**2+0.25*PMQTH1**2) | |
42 | PMTH(3,IF)=PMTH(2,IF)+PMQTH2 | |
43 | PMTH(4,IF)=SQRT(PMTH(1,IF)**2+0.25*PARJ(82)**2)+PMTH(2,21) | |
44 | 100 PMTH(5,IF)=SQRT(PMTH(1,IF)**2+0.25*PARJ(83)**2)+PMTH(2,22) | |
45 | PT2MIN=MAX(0.5*PARJ(82),1.1*PARJ(81))**2 | |
46 | ALAMS=PARJ(81)**2 | |
47 | ALFM=LOG(PT2MIN/ALAMS) | |
48 | ||
49 | C...Store positions of shower initiating partons. | |
50 | M3JC=0 | |
51 | IF(IP1.GT.0.AND.IP1.LE.MIN(N,MSTU(4)-MSTU(32)).AND.IP2.EQ.0) THEN | |
52 | NPA=1 | |
53 | IPA(1)=IP1 | |
54 | ELSEIF(MIN(IP1,IP2).GT.0.AND.MAX(IP1,IP2).LE.MIN(N,MSTU(4)- | |
55 | &MSTU(32))) THEN | |
56 | NPA=2 | |
57 | IPA(1)=IP1 | |
58 | IPA(2)=IP2 | |
59 | ELSEIF(IP1.GT.0.AND.IP1.LE.MIN(N,MSTU(4)-MSTU(32)).AND.IP2.LT.0. | |
60 | &AND.IP2.GE.-3) THEN | |
61 | NPA=IABS(IP2) | |
62 | DO 110 I=1,NPA | |
63 | 110 IPA(I)=IP1+I-1 | |
64 | ELSE | |
65 | CALL LUERRM(12, | |
66 | & '(LUSHOW:) failed to reconstruct showering system') | |
67 | IF(MSTU(21).GE.1) RETURN | |
68 | ENDIF | |
69 | ||
70 | C...Check on phase space available for emission. | |
71 | IREJ=0 | |
72 | DO 120 J=1,5 | |
73 | 120 PS(J)=0. | |
74 | PM=0. | |
75 | DO 130 I=1,NPA | |
76 | KFLA(I)=IABS(K(IPA(I),2)) | |
77 | PMA(I)=P(IPA(I),5) | |
78 | IF(KFLA(I).NE.0.AND.(KFLA(I).LE.8.OR.KFLA(I).EQ.21)) | |
79 | &PMA(I)=PMTH(3,KFLA(I)) | |
80 | PM=PM+PMA(I) | |
81 | IF(KFLA(I).EQ.0.OR.(KFLA(I).GT.8.AND.KFLA(I).NE.21).OR. | |
82 | &PMA(I).GT.QMAX) IREJ=IREJ+1 | |
83 | DO 130 J=1,4 | |
84 | 130 PS(J)=PS(J)+P(IPA(I),J) | |
85 | IF(IREJ.EQ.NPA) RETURN | |
86 | PS(5)=SQRT(MAX(0.,PS(4)**2-PS(1)**2-PS(2)**2-PS(3)**2)) | |
87 | IF(NPA.EQ.1) PS(5)=PS(4) | |
88 | IF(PS(5).LE.PM+PMQTH1) RETURN | |
89 | IF(NPA.EQ.2.AND.MSTJ(47).GE.1) THEN | |
90 | IF(KFLA(1).GE.1.AND.KFLA(1).LE.8.AND.KFLA(2).GE.1.AND. | |
91 | & KFLA(2).LE.8) M3JC=1 | |
92 | IF(MSTJ(47).GE.2) M3JC=1 | |
93 | ENDIF | |
94 | ||
95 | C...Define imagined single initiator of shower for parton system. | |
96 | NS=N | |
97 | IF(N.GT.MSTU(4)-MSTU(32)-5) THEN | |
98 | CALL LUERRM(11,'(LUSHOW:) no more memory left in LUJETS') | |
99 | IF(MSTU(21).GE.1) RETURN | |
100 | ENDIF | |
101 | IF(NPA.GE.2) THEN | |
102 | K(N+1,1)=11 | |
103 | K(N+1,2)=21 | |
104 | K(N+1,3)=0 | |
105 | K(N+1,4)=0 | |
106 | K(N+1,5)=0 | |
107 | P(N+1,1)=0. | |
108 | P(N+1,2)=0. | |
109 | P(N+1,3)=0. | |
110 | P(N+1,4)=PS(5) | |
111 | P(N+1,5)=PS(5) | |
112 | V(N+1,5)=PS(5)**2 | |
113 | N=N+1 | |
114 | ENDIF | |
115 | ||
116 | C...Loop over partons that may branch. | |
117 | NEP=NPA | |
118 | IM=NS | |
119 | IF(NPA.EQ.1) IM=NS-1 | |
120 | 140 IM=IM+1 | |
121 | IF(N.GT.NS) THEN | |
122 | IF(IM.GT.N) GOTO 380 | |
123 | KFLM=IABS(K(IM,2)) | |
124 | IF(KFLM.EQ.0.OR.(KFLM.GT.8.AND.KFLM.NE.21)) GOTO 140 | |
125 | IF(P(IM,5).LT.PMTH(2,KFLM)) GOTO 140 | |
126 | IGM=K(IM,3) | |
127 | ELSE | |
128 | IGM=-1 | |
129 | ENDIF | |
130 | IF(N+NEP.GT.MSTU(4)-MSTU(32)-5) THEN | |
131 | CALL LUERRM(11,'(LUSHOW:) no more memory left in LUJETS') | |
132 | IF(MSTU(21).GE.1) RETURN | |
133 | ENDIF | |
134 | ||
135 | C...Position of aunt (sister to branching parton). | |
136 | C...Origin and flavour of daughters. | |
137 | IAU=0 | |
138 | IF(IGM.GT.0) THEN | |
139 | IF(K(IM-1,3).EQ.IGM) IAU=IM-1 | |
140 | IF(N.GE.IM+1.AND.K(IM+1,3).EQ.IGM) IAU=IM+1 | |
141 | ENDIF | |
142 | IF(IGM.GE.0) THEN | |
143 | K(IM,4)=N+1 | |
144 | DO 150 I=1,NEP | |
145 | 150 K(N+I,3)=IM | |
146 | ELSE | |
147 | K(N+1,3)=IPA(1) | |
148 | ENDIF | |
149 | IF(IGM.LE.0) THEN | |
150 | DO 160 I=1,NEP | |
151 | 160 K(N+I,2)=K(IPA(I),2) | |
152 | ELSEIF(KFLM.NE.21) THEN | |
153 | K(N+1,2)=K(IM,2) | |
154 | K(N+2,2)=K(IM,5) | |
155 | ELSEIF(K(IM,5).EQ.21) THEN | |
156 | K(N+1,2)=21 | |
157 | K(N+2,2)=21 | |
158 | ELSE | |
159 | K(N+1,2)=K(IM,5) | |
160 | K(N+2,2)=-K(IM,5) | |
161 | ENDIF | |
162 | ||
163 | C...Reset flags on daughers and tries made. | |
164 | DO 170 IP=1,NEP | |
165 | K(N+IP,1)=3 | |
166 | K(N+IP,4)=0 | |
167 | K(N+IP,5)=0 | |
168 | KFLD(IP)=IABS(K(N+IP,2)) | |
169 | ITRY(IP)=0 | |
170 | ISL(IP)=0 | |
171 | ISI(IP)=0 | |
172 | 170 IF(KFLD(IP).GT.0.AND.(KFLD(IP).LE.8.OR.KFLD(IP).EQ.21)) ISI(IP)=1 | |
173 | ISLM=0 | |
174 | ||
175 | C...Maximum virtuality of daughters. | |
176 | IF(IGM.LE.0) THEN | |
177 | DO 180 I=1,NPA | |
178 | IF(NPA.GE.3) P(N+I,4)=(PS(4)*P(IPA(I),4)-PS(1)*P(IPA(I),1)- | |
179 | & PS(2)*P(IPA(I),2)-PS(3)*P(IPA(I),3))/PS(5) | |
180 | P(N+I,5)=MIN(QMAX,PS(5)) | |
181 | IF(NPA.GE.3) P(N+I,5)=MIN(P(N+I,5),P(N+I,4)) | |
182 | 180 IF(ISI(I).EQ.0) P(N+I,5)=P(IPA(I),5) | |
183 | ELSE | |
184 | IF(MSTJ(43).LE.2) PEM=V(IM,2) | |
185 | IF(MSTJ(43).GE.3) PEM=P(IM,4) | |
186 | P(N+1,5)=MIN(P(IM,5),V(IM,1)*PEM) | |
187 | P(N+2,5)=MIN(P(IM,5),(1.-V(IM,1))*PEM) | |
188 | IF(K(N+2,2).EQ.22) P(N+2,5)=PMTH(1,22) | |
189 | ENDIF | |
190 | DO 190 I=1,NEP | |
191 | PMSD(I)=P(N+I,5) | |
192 | IF(ISI(I).EQ.1) THEN | |
193 | IF(P(N+I,5).LE.PMTH(3,KFLD(I))) P(N+I,5)=PMTH(1,KFLD(I)) | |
194 | ENDIF | |
195 | 190 V(N+I,5)=P(N+I,5)**2 | |
196 | ||
197 | C...Choose one of the daughters for evolution. | |
198 | 200 INUM=0 | |
199 | IF(NEP.EQ.1) INUM=1 | |
200 | DO 210 I=1,NEP | |
201 | 210 IF(INUM.EQ.0.AND.ISL(I).EQ.1) INUM=I | |
202 | DO 220 I=1,NEP | |
203 | IF(INUM.EQ.0.AND.ITRY(I).EQ.0.AND.ISI(I).EQ.1) THEN | |
204 | IF(P(N+I,5).GE.PMTH(2,KFLD(I))) INUM=I | |
205 | ENDIF | |
206 | 220 CONTINUE | |
207 | IF(INUM.EQ.0) THEN | |
208 | RMAX=0. | |
209 | DO 230 I=1,NEP | |
210 | IF(ISI(I).EQ.1.AND.PMSD(I).GE.PMQTH2) THEN | |
211 | RPM=P(N+I,5)/PMSD(I) | |
212 | IF(RPM.GT.RMAX.AND.P(N+I,5).GE.PMTH(2,KFLD(I))) THEN | |
213 | RMAX=RPM | |
214 | INUM=I | |
215 | ENDIF | |
216 | ENDIF | |
217 | 230 CONTINUE | |
218 | ENDIF | |
219 | ||
220 | C...Store information on choice of evolving daughter. | |
221 | INUM=MAX(1,INUM) | |
222 | IEP(1)=N+INUM | |
223 | DO 240 I=2,NEP | |
224 | IEP(I)=IEP(I-1)+1 | |
225 | 240 IF(IEP(I).GT.N+NEP) IEP(I)=N+1 | |
226 | DO 250 I=1,NEP | |
227 | 250 KFL(I)=IABS(K(IEP(I),2)) | |
228 | ITRY(INUM)=ITRY(INUM)+1 | |
229 | IF(ITRY(INUM).GT.200) THEN | |
230 | CALL LUERRM(14,'(LUSHOW:) caught in infinite loop') | |
231 | IF(MSTU(21).GE.1) RETURN | |
232 | ENDIF | |
233 | Z=0.5 | |
234 | IF(KFL(1).EQ.0.OR.(KFL(1).GT.8.AND.KFL(1).NE.21)) GOTO 300 | |
235 | IF(P(IEP(1),5).LT.PMTH(2,KFL(1))) GOTO 300 | |
236 | ||
237 | C...Calculate allowed z range. | |
238 | IF(NEP.EQ.1) THEN | |
239 | PMED=PS(4) | |
240 | ELSEIF(IGM.EQ.0.OR.MSTJ(43).LE.2) THEN | |
241 | PMED=P(IM,5) | |
242 | ELSE | |
243 | IF(INUM.EQ.1) PMED=V(IM,1)*PEM | |
244 | IF(INUM.EQ.2) PMED=(1.-V(IM,1))*PEM | |
245 | ENDIF | |
246 | IF(MOD(MSTJ(43),2).EQ.1) THEN | |
247 | ZC=PMTH(2,21)/PMED | |
248 | ZCE=PMTH(2,22)/PMED | |
249 | ELSE | |
250 | ZC=0.5*(1.-SQRT(MAX(0.,1.-(2.*PMTH(2,21)/PMED)**2))) | |
251 | IF(ZC.LT.1E-4) ZC=(PMTH(2,21)/PMED)**2 | |
252 | ZCE=0.5*(1.-SQRT(MAX(0.,1.-(2.*PMTH(2,22)/PMED)**2))) | |
253 | IF(ZCE.LT.1E-4) ZCE=(PMTH(2,22)/PMED)**2 | |
254 | ENDIF | |
255 | ZC=MIN(ZC,0.491) | |
256 | ZCE=MIN(ZCE,0.491) | |
257 | IF((MSTJ(41).EQ.1.AND.ZC.GT.0.49).OR.(MSTJ(41).EQ.2.AND. | |
258 | &MIN(ZC,ZCE).GT.0.49)) THEN | |
259 | P(IEP(1),5)=PMTH(1,KFL(1)) | |
260 | V(IEP(1),5)=P(IEP(1),5)**2 | |
261 | GOTO 300 | |
262 | ENDIF | |
263 | ||
264 | C...Integral of Altarelli-Parisi z kernel for QCD. | |
265 | IF(MSTJ(49).EQ.0.AND.KFL(1).EQ.21) THEN | |
266 | FBR=6.*LOG((1.-ZC)/ZC)+MSTJ(45)*(0.5-ZC) | |
267 | ELSEIF(MSTJ(49).EQ.0) THEN | |
268 | FBR=(8./3.)*LOG((1.-ZC)/ZC) | |
269 | ||
270 | C...Integral of Altarelli-Parisi z kernel for scalar gluon. | |
271 | ELSEIF(MSTJ(49).EQ.1.AND.KFL(1).EQ.21) THEN | |
272 | FBR=(PARJ(87)+MSTJ(45)*PARJ(88))*(1.-2.*ZC) | |
273 | ELSEIF(MSTJ(49).EQ.1) THEN | |
274 | FBR=(1.-2.*ZC)/3. | |
275 | IF(IGM.EQ.0.AND.M3JC.EQ.1) FBR=4.*FBR | |
276 | ||
277 | C...Integral of Altarelli-Parisi z kernel for Abelian vector gluon. | |
278 | ELSEIF(KFL(1).EQ.21) THEN | |
279 | FBR=6.*MSTJ(45)*(0.5-ZC) | |
280 | ELSE | |
281 | FBR=2.*LOG((1.-ZC)/ZC) | |
282 | ENDIF | |
283 | ||
284 | C...Integral of Altarelli-Parisi kernel for photon emission. | |
285 | IF(MSTJ(41).EQ.2.AND.KFL(1).GE.1.AND.KFL(1).LE.8) | |
286 | &FBRE=(KCHG(KFL(1),1)/3.)**2*2.*LOG((1.-ZCE)/ZCE) | |
287 | ||
288 | C...Inner veto algorithm starts. Find maximum mass for evolution. | |
289 | 260 PMS=V(IEP(1),5) | |
290 | IF(IGM.GE.0) THEN | |
291 | PM2=0. | |
292 | DO 270 I=2,NEP | |
293 | PM=P(IEP(I),5) | |
294 | IF(KFL(I).GT.0.AND.(KFL(I).LE.8.OR.KFL(I).EQ.21)) PM= | |
295 | & PMTH(2,KFL(I)) | |
296 | 270 PM2=PM2+PM | |
297 | PMS=MIN(PMS,(P(IM,5)-PM2)**2) | |
298 | ENDIF | |
299 | ||
300 | C...Select mass for daughter in QCD evolution. | |
301 | B0=27./6. | |
302 | DO 280 IF=4,MSTJ(45) | |
303 | 280 IF(PMS.GT.4.*PMTH(2,IF)**2) B0=(33.-2.*IF)/6. | |
304 | IF(MSTJ(44).LE.0) THEN | |
305 | PMSQCD=PMS*EXP(MAX(-80.,LOG(RLU(0))*PARU(2)/(PARU(111)*FBR))) | |
306 | ELSEIF(MSTJ(44).EQ.1) THEN | |
307 | PMSQCD=4.*ALAMS*(0.25*PMS/ALAMS)**(RLU(0)**(B0/FBR)) | |
308 | ELSE | |
309 | PMSQCD=PMS*RLU(0)**(ALFM*B0/FBR) | |
310 | ENDIF | |
311 | IF(ZC.GT.0.49.OR.PMSQCD.LE.PMTH(4,KFL(1))**2) PMSQCD= | |
312 | &PMTH(2,KFL(1))**2 | |
313 | V(IEP(1),5)=PMSQCD | |
314 | MCE=1 | |
315 | ||
316 | C...Select mass for daughter in QED evolution. | |
317 | IF(MSTJ(41).EQ.2.AND.KFL(1).GE.1.AND.KFL(1).LE.8) THEN | |
318 | PMSQED=PMS*EXP(MAX(-80.,LOG(RLU(0))*PARU(2)/(PARU(101)*FBRE))) | |
319 | IF(ZCE.GT.0.49.OR.PMSQED.LE.PMTH(5,KFL(1))**2) PMSQED= | |
320 | & PMTH(2,KFL(1))**2 | |
321 | IF(PMSQED.GT.PMSQCD) THEN | |
322 | V(IEP(1),5)=PMSQED | |
323 | MCE=2 | |
324 | ENDIF | |
325 | ENDIF | |
326 | ||
327 | C...Check whether daughter mass below cutoff. | |
328 | P(IEP(1),5)=SQRT(V(IEP(1),5)) | |
329 | IF(P(IEP(1),5).LE.PMTH(3,KFL(1))) THEN | |
330 | P(IEP(1),5)=PMTH(1,KFL(1)) | |
331 | V(IEP(1),5)=P(IEP(1),5)**2 | |
332 | GOTO 300 | |
333 | ENDIF | |
334 | ||
335 | C...Select z value of branching: q -> qgamma. | |
336 | IF(MCE.EQ.2) THEN | |
337 | Z=1.-(1.-ZCE)*(ZCE/(1.-ZCE))**RLU(0) | |
338 | IF(1.+Z**2.LT.2.*RLU(0)) GOTO 260 | |
339 | K(IEP(1),5)=22 | |
340 | ||
341 | C...Select z value of branching: q -> qg, g -> gg, g -> qqbar. | |
342 | ELSEIF(MSTJ(49).NE.1.AND.KFL(1).NE.21) THEN | |
343 | Z=1.-(1.-ZC)*(ZC/(1.-ZC))**RLU(0) | |
344 | IF(1.+Z**2.LT.2.*RLU(0)) GOTO 260 | |
345 | K(IEP(1),5)=21 | |
346 | ELSEIF(MSTJ(49).EQ.0.AND.MSTJ(45)*(0.5-ZC).LT.RLU(0)*FBR) THEN | |
347 | Z=(1.-ZC)*(ZC/(1.-ZC))**RLU(0) | |
348 | IF(RLU(0).GT.0.5) Z=1.-Z | |
349 | IF((1.-Z*(1.-Z))**2.LT.RLU(0)) GOTO 260 | |
350 | K(IEP(1),5)=21 | |
351 | ELSEIF(MSTJ(49).NE.1) THEN | |
352 | Z=ZC+(1.-2.*ZC)*RLU(0) | |
353 | IF(Z**2+(1.-Z)**2.LT.RLU(0)) GOTO 260 | |
354 | KFLB=1+INT(MSTJ(45)*RLU(0)) | |
355 | PMQ=4.*PMTH(2,KFLB)**2/V(IEP(1),5) | |
356 | IF(PMQ.GE.1.) GOTO 260 | |
357 | PMQ0=4.*PMTH(2,21)**2/V(IEP(1),5) | |
358 | IF(MOD(MSTJ(43),2).EQ.0.AND.(1.+0.5*PMQ)*SQRT(1.-PMQ).LT. | |
359 | & RLU(0)*(1.+0.5*PMQ0)*SQRT(1.-PMQ0)) GOTO 260 | |
360 | K(IEP(1),5)=KFLB | |
361 | ||
362 | C...Ditto for scalar gluon model. | |
363 | ELSEIF(KFL(1).NE.21) THEN | |
364 | Z=1.-SQRT(ZC**2+RLU(0)*(1.-2.*ZC)) | |
365 | K(IEP(1),5)=21 | |
366 | ELSEIF(RLU(0)*(PARJ(87)+MSTJ(45)*PARJ(88)).LE.PARJ(87)) THEN | |
367 | Z=ZC+(1.-2.*ZC)*RLU(0) | |
368 | K(IEP(1),5)=21 | |
369 | ELSE | |
370 | Z=ZC+(1.-2.*ZC)*RLU(0) | |
371 | KFLB=1+INT(MSTJ(45)*RLU(0)) | |
372 | PMQ=4.*PMTH(2,KFLB)**2/V(IEP(1),5) | |
373 | IF(PMQ.GE.1.) GOTO 260 | |
374 | K(IEP(1),5)=KFLB | |
375 | ENDIF | |
376 | IF(MCE.EQ.1.AND.MSTJ(44).GE.2) THEN | |
377 | IF(Z*(1.-Z)*V(IEP(1),5).LT.PT2MIN) GOTO 260 | |
378 | IF(ALFM/LOG(V(IEP(1),5)*Z*(1.-Z)/ALAMS).LT.RLU(0)) GOTO 260 | |
379 | ENDIF | |
380 | ||
381 | C...Check if z consistent with chosen m. | |
382 | IF(KFL(1).EQ.21) THEN | |
383 | KFLGD1=IABS(K(IEP(1),5)) | |
384 | KFLGD2=KFLGD1 | |
385 | ELSE | |
386 | KFLGD1=KFL(1) | |
387 | KFLGD2=IABS(K(IEP(1),5)) | |
388 | ENDIF | |
389 | IF(NEP.EQ.1) THEN | |
390 | PED=PS(4) | |
391 | ELSEIF(NEP.GE.3) THEN | |
392 | PED=P(IEP(1),4) | |
393 | ELSEIF(IGM.EQ.0.OR.MSTJ(43).LE.2) THEN | |
394 | PED=0.5*(V(IM,5)+V(IEP(1),5)-PM2**2)/P(IM,5) | |
395 | ELSE | |
396 | IF(IEP(1).EQ.N+1) PED=V(IM,1)*PEM | |
397 | IF(IEP(1).EQ.N+2) PED=(1.-V(IM,1))*PEM | |
398 | ENDIF | |
399 | IF(MOD(MSTJ(43),2).EQ.1) THEN | |
400 | PMQTH3=0.5*PARJ(82) | |
401 | IF(KFLGD2.EQ.22) PMQTH3=0.5*PARJ(83) | |
402 | PMQ1=(PMTH(1,KFLGD1)**2+PMQTH3**2)/V(IEP(1),5) | |
403 | PMQ2=(PMTH(1,KFLGD2)**2+PMQTH3**2)/V(IEP(1),5) | |
404 | ZD=SQRT(MAX(0.,(1.-V(IEP(1),5)/PED**2)*((1.-PMQ1-PMQ2)**2- | |
405 | & 4.*PMQ1*PMQ2))) | |
406 | ZH=1.+PMQ1-PMQ2 | |
407 | ELSE | |
408 | ZD=SQRT(MAX(0.,1.-V(IEP(1),5)/PED**2)) | |
409 | ZH=1. | |
410 | ENDIF | |
411 | ZL=0.5*(ZH-ZD) | |
412 | ZU=0.5*(ZH+ZD) | |
413 | IF(Z.LT.ZL.OR.Z.GT.ZU) GOTO 260 | |
414 | IF(KFL(1).EQ.21) V(IEP(1),3)=LOG(ZU*(1.-ZL)/MAX(1E-20,ZL* | |
415 | &(1.-ZU))) | |
416 | IF(KFL(1).NE.21) V(IEP(1),3)=LOG((1.-ZL)/MAX(1E-10,1.-ZU)) | |
417 | ||
418 | C...Three-jet matrix element correction. | |
419 | IF(IGM.EQ.0.AND.M3JC.EQ.1) THEN | |
420 | X1=Z*(1.+V(IEP(1),5)/V(NS+1,5)) | |
421 | X2=1.-V(IEP(1),5)/V(NS+1,5) | |
422 | X3=(1.-X1)+(1.-X2) | |
423 | IF(MCE.EQ.2) THEN | |
424 | KI1=K(IPA(INUM),2) | |
425 | KI2=K(IPA(3-INUM),2) | |
426 | QF1=KCHG(IABS(KI1),1)*ISIGN(1,KI1)/3. | |
427 | QF2=KCHG(IABS(KI2),1)*ISIGN(1,KI2)/3. | |
428 | WSHOW=QF1**2*(1.-X1)/X3*(1.+(X1/(2.-X2))**2)+ | |
429 | & QF2**2*(1.-X2)/X3*(1.+(X2/(2.-X1))**2) | |
430 | WME=(QF1*(1.-X1)/X3-QF2*(1.-X2)/X3)**2*(X1**2+X2**2) | |
431 | ELSEIF(MSTJ(49).NE.1) THEN | |
432 | WSHOW=1.+(1.-X1)/X3*(X1/(2.-X2))**2+ | |
433 | & (1.-X2)/X3*(X2/(2.-X1))**2 | |
434 | WME=X1**2+X2**2 | |
435 | ELSE | |
436 | WSHOW=4.*X3*((1.-X1)/(2.-X2)**2+(1.-X2)/(2.-X1)**2) | |
437 | WME=X3**2 | |
438 | ENDIF | |
439 | IF(WME.LT.RLU(0)*WSHOW) GOTO 260 | |
440 | ||
441 | C...Impose angular ordering by rejection of nonordered emission. | |
442 | ELSEIF(MCE.EQ.1.AND.IGM.GT.0.AND.MSTJ(42).GE.2) THEN | |
443 | MAOM=1 | |
444 | ZM=V(IM,1) | |
445 | IF(IEP(1).EQ.N+2) ZM=1.-V(IM,1) | |
446 | THE2ID=Z*(1.-Z)*(ZM*P(IM,4))**2/V(IEP(1),5) | |
447 | IAOM=IM | |
448 | 290 IF(K(IAOM,5).EQ.22) THEN | |
449 | IAOM=K(IAOM,3) | |
450 | IF(K(IAOM,3).LE.NS) MAOM=0 | |
451 | IF(MAOM.EQ.1) GOTO 290 | |
452 | ENDIF | |
453 | IF(MAOM.EQ.1) THEN | |
454 | THE2IM=V(IAOM,1)*(1.-V(IAOM,1))*P(IAOM,4)**2/V(IAOM,5) | |
455 | IF(THE2ID.LT.THE2IM) GOTO 260 | |
456 | ENDIF | |
457 | ENDIF | |
458 | ||
459 | C...Impose user-defined maximum angle at first branching. | |
460 | IF(MSTJ(48).EQ.1) THEN | |
461 | IF(NEP.EQ.1.AND.IM.EQ.NS) THEN | |
462 | THE2ID=Z*(1.-Z)*PS(4)**2/V(IEP(1),5) | |
463 | IF(THE2ID.LT.1./PARJ(85)**2) GOTO 260 | |
464 | ELSEIF(NEP.EQ.2.AND.IEP(1).EQ.NS+2) THEN | |
465 | THE2ID=Z*(1.-Z)*(0.5*P(IM,4))**2/V(IEP(1),5) | |
466 | IF(THE2ID.LT.1./PARJ(85)**2) GOTO 260 | |
467 | ELSEIF(NEP.EQ.2.AND.IEP(1).EQ.NS+3) THEN | |
468 | THE2ID=Z*(1.-Z)*(0.5*P(IM,4))**2/V(IEP(1),5) | |
469 | IF(THE2ID.LT.1./PARJ(86)**2) GOTO 260 | |
470 | ENDIF | |
471 | ENDIF | |
472 | ||
473 | C...End of inner veto algorithm. Check if only one leg evolved so far. | |
474 | 300 V(IEP(1),1)=Z | |
475 | ISL(1)=0 | |
476 | ISL(2)=0 | |
477 | IF(NEP.EQ.1) GOTO 330 | |
478 | IF(NEP.EQ.2.AND.P(IEP(1),5)+P(IEP(2),5).GE.P(IM,5)) GOTO 200 | |
479 | DO 310 I=1,NEP | |
480 | IF(ITRY(I).EQ.0.AND.KFLD(I).GT.0.AND.(KFLD(I).LE.8.OR.KFLD(I).EQ. | |
481 | &21)) THEN | |
482 | IF(P(N+I,5).GE.PMTH(2,KFLD(I))) GOTO 200 | |
483 | ENDIF | |
484 | 310 CONTINUE | |
485 | ||
486 | C...Check if chosen multiplet m1,m2,z1,z2 is physical. | |
487 | IF(NEP.EQ.3) THEN | |
488 | PA1S=(P(N+1,4)+P(N+1,5))*(P(N+1,4)-P(N+1,5)) | |
489 | PA2S=(P(N+2,4)+P(N+2,5))*(P(N+2,4)-P(N+2,5)) | |
490 | PA3S=(P(N+3,4)+P(N+3,5))*(P(N+3,4)-P(N+3,5)) | |
491 | PTS=0.25*(2.*PA1S*PA2S+2.*PA1S*PA3S+2.*PA2S*PA3S- | |
492 | & PA1S**2-PA2S**2-PA3S**2)/PA1S | |
493 | IF(PTS.LE.0.) GOTO 200 | |
494 | ELSEIF(IGM.EQ.0.OR.MSTJ(43).LE.2.OR.MOD(MSTJ(43),2).EQ.0) THEN | |
495 | DO 320 I1=N+1,N+2 | |
496 | KFLDA=IABS(K(I1,2)) | |
497 | IF(KFLDA.EQ.0.OR.(KFLDA.GT.8.AND.KFLDA.NE.21)) GOTO 320 | |
498 | IF(P(I1,5).LT.PMTH(2,KFLDA)) GOTO 320 | |
499 | IF(KFLDA.EQ.21) THEN | |
500 | KFLGD1=IABS(K(I1,5)) | |
501 | KFLGD2=KFLGD1 | |
502 | ELSE | |
503 | KFLGD1=KFLDA | |
504 | KFLGD2=IABS(K(I1,5)) | |
505 | ENDIF | |
506 | I2=2*N+3-I1 | |
507 | IF(IGM.EQ.0.OR.MSTJ(43).LE.2) THEN | |
508 | PED=0.5*(V(IM,5)+V(I1,5)-V(I2,5))/P(IM,5) | |
509 | ELSE | |
510 | IF(I1.EQ.N+1) ZM=V(IM,1) | |
511 | IF(I1.EQ.N+2) ZM=1.-V(IM,1) | |
512 | PML=SQRT((V(IM,5)-V(N+1,5)-V(N+2,5))**2- | |
513 | & 4.*V(N+1,5)*V(N+2,5)) | |
514 | PED=PEM*(0.5*(V(IM,5)-PML+V(I1,5)-V(I2,5))+PML*ZM)/V(IM,5) | |
515 | ENDIF | |
516 | IF(MOD(MSTJ(43),2).EQ.1) THEN | |
517 | PMQTH3=0.5*PARJ(82) | |
518 | IF(KFLGD2.EQ.22) PMQTH3=0.5*PARJ(83) | |
519 | PMQ1=(PMTH(1,KFLGD1)**2+PMQTH3**2)/V(I1,5) | |
520 | PMQ2=(PMTH(1,KFLGD2)**2+PMQTH3**2)/V(I1,5) | |
521 | ZD=SQRT(MAX(0.,(1.-V(I1,5)/PED**2)*((1.-PMQ1-PMQ2)**2- | |
522 | & 4.*PMQ1*PMQ2))) | |
523 | ZH=1.+PMQ1-PMQ2 | |
524 | ELSE | |
525 | ZD=SQRT(MAX(0.,1.-V(I1,5)/PED**2)) | |
526 | ZH=1. | |
527 | ENDIF | |
528 | ZL=0.5*(ZH-ZD) | |
529 | ZU=0.5*(ZH+ZD) | |
530 | IF(I1.EQ.N+1.AND.(V(I1,1).LT.ZL.OR.V(I1,1).GT.ZU)) ISL(1)=1 | |
531 | IF(I1.EQ.N+2.AND.(V(I1,1).LT.ZL.OR.V(I1,1).GT.ZU)) ISL(2)=1 | |
532 | IF(KFLDA.EQ.21) V(I1,4)=LOG(ZU*(1.-ZL)/MAX(1E-20,ZL*(1.-ZU))) | |
533 | IF(KFLDA.NE.21) V(I1,4)=LOG((1.-ZL)/MAX(1E-10,1.-ZU)) | |
534 | 320 CONTINUE | |
535 | IF(ISL(1).EQ.1.AND.ISL(2).EQ.1.AND.ISLM.NE.0) THEN | |
536 | ISL(3-ISLM)=0 | |
537 | ISLM=3-ISLM | |
538 | ELSEIF(ISL(1).EQ.1.AND.ISL(2).EQ.1) THEN | |
539 | ZDR1=MAX(0.,V(N+1,3)/V(N+1,4)-1.) | |
540 | ZDR2=MAX(0.,V(N+2,3)/V(N+2,4)-1.) | |
541 | IF(ZDR2.GT.RLU(0)*(ZDR1+ZDR2)) ISL(1)=0 | |
542 | IF(ISL(1).EQ.1) ISL(2)=0 | |
543 | IF(ISL(1).EQ.0) ISLM=1 | |
544 | IF(ISL(2).EQ.0) ISLM=2 | |
545 | ENDIF | |
546 | IF(ISL(1).EQ.1.OR.ISL(2).EQ.1) GOTO 200 | |
547 | ENDIF | |
548 | IF(IGM.GT.0.AND.MOD(MSTJ(43),2).EQ.1.AND.(P(N+1,5).GE. | |
549 | &PMTH(2,KFLD(1)).OR.P(N+2,5).GE.PMTH(2,KFLD(2)))) THEN | |
550 | PMQ1=V(N+1,5)/V(IM,5) | |
551 | PMQ2=V(N+2,5)/V(IM,5) | |
552 | ZD=SQRT(MAX(0.,(1.-V(IM,5)/PEM**2)*((1.-PMQ1-PMQ2)**2- | |
553 | & 4.*PMQ1*PMQ2))) | |
554 | ZH=1.+PMQ1-PMQ2 | |
555 | ZL=0.5*(ZH-ZD) | |
556 | ZU=0.5*(ZH+ZD) | |
557 | IF(V(IM,1).LT.ZL.OR.V(IM,1).GT.ZU) GOTO 200 | |
558 | ENDIF | |
559 | ||
560 | C...Accepted branch. Construct four-momentum for initial partons. | |
561 | 330 MAZIP=0 | |
562 | MAZIC=0 | |
563 | IF(NEP.EQ.1) THEN | |
564 | P(N+1,1)=0. | |
565 | P(N+1,2)=0. | |
566 | P(N+1,3)=SQRT(MAX(0.,(P(IPA(1),4)+P(N+1,5))*(P(IPA(1),4)- | |
567 | & P(N+1,5)))) | |
568 | P(N+1,4)=P(IPA(1),4) | |
569 | V(N+1,2)=P(N+1,4) | |
570 | ELSEIF(IGM.EQ.0.AND.NEP.EQ.2) THEN | |
571 | PED1=0.5*(V(IM,5)+V(N+1,5)-V(N+2,5))/P(IM,5) | |
572 | P(N+1,1)=0. | |
573 | P(N+1,2)=0. | |
574 | P(N+1,3)=SQRT(MAX(0.,(PED1+P(N+1,5))*(PED1-P(N+1,5)))) | |
575 | P(N+1,4)=PED1 | |
576 | P(N+2,1)=0. | |
577 | P(N+2,2)=0. | |
578 | P(N+2,3)=-P(N+1,3) | |
579 | P(N+2,4)=P(IM,5)-PED1 | |
580 | V(N+1,2)=P(N+1,4) | |
581 | V(N+2,2)=P(N+2,4) | |
582 | ELSEIF(NEP.EQ.3) THEN | |
583 | P(N+1,1)=0. | |
584 | P(N+1,2)=0. | |
585 | P(N+1,3)=SQRT(MAX(0.,PA1S)) | |
586 | P(N+2,1)=SQRT(PTS) | |
587 | P(N+2,2)=0. | |
588 | P(N+2,3)=0.5*(PA3S-PA2S-PA1S)/P(N+1,3) | |
589 | P(N+3,1)=-P(N+2,1) | |
590 | P(N+3,2)=0. | |
591 | P(N+3,3)=-(P(N+1,3)+P(N+2,3)) | |
592 | V(N+1,2)=P(N+1,4) | |
593 | V(N+2,2)=P(N+2,4) | |
594 | V(N+3,2)=P(N+3,4) | |
595 | ||
596 | C...Construct transverse momentum for ordinary branching in shower. | |
597 | ELSE | |
598 | ZM=V(IM,1) | |
599 | PZM=SQRT(MAX(0.,(PEM+P(IM,5))*(PEM-P(IM,5)))) | |
600 | PMLS=(V(IM,5)-V(N+1,5)-V(N+2,5))**2-4.*V(N+1,5)*V(N+2,5) | |
601 | IF(PZM.LE.0.) THEN | |
602 | PTS=0. | |
603 | ELSEIF(MOD(MSTJ(43),2).EQ.1) THEN | |
604 | PTS=(PEM**2*(ZM*(1.-ZM)*V(IM,5)-(1.-ZM)*V(N+1,5)- | |
605 | & ZM*V(N+2,5))-0.25*PMLS)/PZM**2 | |
606 | ELSE | |
607 | PTS=PMLS*(ZM*(1.-ZM)*PEM**2/V(IM,5)-0.25)/PZM**2 | |
608 | ENDIF | |
609 | PT=SQRT(MAX(0.,PTS)) | |
610 | ||
611 | C...Find coefficient of azimuthal asymmetry due to gluon polarization. | |
612 | HAZIP=0. | |
613 | IF(MSTJ(49).NE.1.AND.MOD(MSTJ(46),2).EQ.1.AND.K(IM,2).EQ.21. | |
614 | & AND.IAU.NE.0) THEN | |
615 | IF(K(IGM,3).NE.0) MAZIP=1 | |
616 | ZAU=V(IGM,1) | |
617 | IF(IAU.EQ.IM+1) ZAU=1.-V(IGM,1) | |
618 | IF(MAZIP.EQ.0) ZAU=0. | |
619 | IF(K(IGM,2).NE.21) THEN | |
620 | HAZIP=2.*ZAU/(1.+ZAU**2) | |
621 | ELSE | |
622 | HAZIP=(ZAU/(1.-ZAU*(1.-ZAU)))**2 | |
623 | ENDIF | |
624 | IF(K(N+1,2).NE.21) THEN | |
625 | HAZIP=HAZIP*(-2.*ZM*(1.-ZM))/(1.-2.*ZM*(1.-ZM)) | |
626 | ELSE | |
627 | HAZIP=HAZIP*(ZM*(1.-ZM)/(1.-ZM*(1.-ZM)))**2 | |
628 | ENDIF | |
629 | ENDIF | |
630 | ||
631 | C...Find coefficient of azimuthal asymmetry due to soft gluon | |
632 | C...interference. | |
633 | HAZIC=0. | |
634 | IF(MSTJ(49).NE.2.AND.MSTJ(46).GE.2.AND.(K(N+1,2).EQ.21.OR. | |
635 | & K(N+2,2).EQ.21).AND.IAU.NE.0) THEN | |
636 | IF(K(IGM,3).NE.0) MAZIC=N+1 | |
637 | IF(K(IGM,3).NE.0.AND.K(N+1,2).NE.21) MAZIC=N+2 | |
638 | IF(K(IGM,3).NE.0.AND.K(N+1,2).EQ.21.AND.K(N+2,2).EQ.21.AND. | |
639 | & ZM.GT.0.5) MAZIC=N+2 | |
640 | IF(K(IAU,2).EQ.22) MAZIC=0 | |
641 | ZS=ZM | |
642 | IF(MAZIC.EQ.N+2) ZS=1.-ZM | |
643 | ZGM=V(IGM,1) | |
644 | IF(IAU.EQ.IM-1) ZGM=1.-V(IGM,1) | |
645 | IF(MAZIC.EQ.0) ZGM=1. | |
646 | HAZIC=(P(IM,5)/P(IGM,5))*SQRT((1.-ZS)*(1.-ZGM)/(ZS*ZGM)) | |
647 | HAZIC=MIN(0.95,HAZIC) | |
648 | ENDIF | |
649 | ENDIF | |
650 | ||
651 | C...Construct kinematics for ordinary branching in shower. | |
652 | 340 IF(NEP.EQ.2.AND.IGM.GT.0) THEN | |
653 | IF(MOD(MSTJ(43),2).EQ.1) THEN | |
654 | P(N+1,4)=PEM*V(IM,1) | |
655 | ELSE | |
656 | P(N+1,4)=PEM*(0.5*(V(IM,5)-SQRT(PMLS)+V(N+1,5)-V(N+2,5))+ | |
657 | & SQRT(PMLS)*ZM)/V(IM,5) | |
658 | ENDIF | |
659 | PHI=PARU(2)*RLU(0) | |
660 | P(N+1,1)=PT*COS(PHI) | |
661 | P(N+1,2)=PT*SIN(PHI) | |
662 | IF(PZM.GT.0.) THEN | |
663 | 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 | ELSE | |
665 | P(N+1,3)=0. | |
666 | ENDIF | |
667 | P(N+2,1)=-P(N+1,1) | |
668 | P(N+2,2)=-P(N+1,2) | |
669 | P(N+2,3)=PZM-P(N+1,3) | |
670 | P(N+2,4)=PEM-P(N+1,4) | |
671 | IF(MSTJ(43).LE.2) THEN | |
672 | V(N+1,2)=(PEM*P(N+1,4)-PZM*P(N+1,3))/P(IM,5) | |
673 | V(N+2,2)=(PEM*P(N+2,4)-PZM*P(N+2,3))/P(IM,5) | |
674 | ENDIF | |
675 | ENDIF | |
676 | ||
677 | C...Rotate and boost daughters. | |
678 | IF(IGM.GT.0) THEN | |
679 | IF(MSTJ(43).LE.2) THEN | |
680 | BEX=P(IGM,1)/P(IGM,4) | |
681 | BEY=P(IGM,2)/P(IGM,4) | |
682 | BEZ=P(IGM,3)/P(IGM,4) | |
683 | GA=P(IGM,4)/P(IGM,5) | |
684 | GABEP=GA*(GA*(BEX*P(IM,1)+BEY*P(IM,2)+BEZ*P(IM,3))/(1.+GA)- | |
685 | & P(IM,4)) | |
686 | ELSE | |
687 | BEX=0. | |
688 | BEY=0. | |
689 | BEZ=0. | |
690 | GA=1. | |
691 | GABEP=0. | |
692 | ENDIF | |
693 | THE=ULANGL(P(IM,3)+GABEP*BEZ,SQRT((P(IM,1)+GABEP*BEX)**2+ | |
694 | & (P(IM,2)+GABEP*BEY)**2)) | |
695 | PHI=ULANGL(P(IM,1)+GABEP*BEX,P(IM,2)+GABEP*BEY) | |
696 | DO 350 I=N+1,N+2 | |
697 | DP(1)=COS(THE)*COS(PHI)*P(I,1)-SIN(PHI)*P(I,2)+ | |
698 | & SIN(THE)*COS(PHI)*P(I,3) | |
699 | DP(2)=COS(THE)*SIN(PHI)*P(I,1)+COS(PHI)*P(I,2)+ | |
700 | & SIN(THE)*SIN(PHI)*P(I,3) | |
701 | DP(3)=-SIN(THE)*P(I,1)+COS(THE)*P(I,3) | |
702 | DP(4)=P(I,4) | |
703 | DBP=BEX*DP(1)+BEY*DP(2)+BEZ*DP(3) | |
704 | DGABP=GA*(GA*DBP/(1D0+GA)+DP(4)) | |
705 | P(I,1)=DP(1)+DGABP*BEX | |
706 | P(I,2)=DP(2)+DGABP*BEY | |
707 | P(I,3)=DP(3)+DGABP*BEZ | |
708 | 350 P(I,4)=GA*(DP(4)+DBP) | |
709 | ENDIF | |
710 | ||
711 | C...Weight with azimuthal distribution, if required. | |
712 | IF(MAZIP.NE.0.OR.MAZIC.NE.0) THEN | |
713 | DO 360 J=1,3 | |
714 | DPT(1,J)=P(IM,J) | |
715 | DPT(2,J)=P(IAU,J) | |
716 | 360 DPT(3,J)=P(N+1,J) | |
717 | DPMA=DPT(1,1)*DPT(2,1)+DPT(1,2)*DPT(2,2)+DPT(1,3)*DPT(2,3) | |
718 | DPMD=DPT(1,1)*DPT(3,1)+DPT(1,2)*DPT(3,2)+DPT(1,3)*DPT(3,3) | |
719 | DPMM=DPT(1,1)**2+DPT(1,2)**2+DPT(1,3)**2 | |
720 | DO 370 J=1,3 | |
721 | DPT(4,J)=DPT(2,J)-DPMA*DPT(1,J)/DPMM | |
722 | 370 DPT(5,J)=DPT(3,J)-DPMD*DPT(1,J)/DPMM | |
723 | DPT(4,4)=SQRT(DPT(4,1)**2+DPT(4,2)**2+DPT(4,3)**2) | |
724 | DPT(5,4)=SQRT(DPT(5,1)**2+DPT(5,2)**2+DPT(5,3)**2) | |
725 | IF(MIN(DPT(4,4),DPT(5,4)).GT.0.1*PARJ(82)) THEN | |
726 | CAD=(DPT(4,1)*DPT(5,1)+DPT(4,2)*DPT(5,2)+ | |
727 | & DPT(4,3)*DPT(5,3))/(DPT(4,4)*DPT(5,4)) | |
728 | IF(MAZIP.NE.0) THEN | |
729 | IF(1.+HAZIP*(2.*CAD**2-1.).LT.RLU(0)*(1.+ABS(HAZIP))) | |
730 | & GOTO 340 | |
731 | ENDIF | |
732 | IF(MAZIC.NE.0) THEN | |
733 | IF(MAZIC.EQ.N+2) CAD=-CAD | |
734 | IF((1.-HAZIC)*(1.-HAZIC*CAD)/(1.+HAZIC**2-2.*HAZIC*CAD). | |
735 | & LT.RLU(0)) GOTO 340 | |
736 | ENDIF | |
737 | ENDIF | |
738 | ENDIF | |
739 | ||
740 | C...Continue loop over partons that may branch, until none left. | |
741 | IF(IGM.GE.0) K(IM,1)=14 | |
742 | N=N+NEP | |
743 | NEP=2 | |
744 | IF(N.GT.MSTU(4)-MSTU(32)-5) THEN | |
745 | CALL LUERRM(11,'(LUSHOW:) no more memory left in LUJETS') | |
746 | IF(MSTU(21).GE.1) N=NS | |
747 | IF(MSTU(21).GE.1) RETURN | |
748 | ENDIF | |
749 | GOTO 140 | |
750 | ||
751 | C...Set information on imagined shower initiator. | |
752 | 380 IF(NPA.GE.2) THEN | |
753 | K(NS+1,1)=11 | |
754 | K(NS+1,2)=94 | |
755 | K(NS+1,3)=IP1 | |
756 | IF(IP2.GT.0.AND.IP2.LT.IP1) K(NS+1,3)=IP2 | |
757 | K(NS+1,4)=NS+2 | |
758 | K(NS+1,5)=NS+1+NPA | |
759 | IIM=1 | |
760 | ELSE | |
761 | IIM=0 | |
762 | ENDIF | |
763 | ||
764 | C...Reconstruct string drawing information. | |
765 | DO 390 I=NS+1+IIM,N | |
766 | IF(K(I,1).LE.10.AND.K(I,2).EQ.22) THEN | |
767 | K(I,1)=1 | |
768 | ELSEIF(K(I,1).LE.10) THEN | |
769 | K(I,4)=MSTU(5)*(K(I,4)/MSTU(5)) | |
770 | K(I,5)=MSTU(5)*(K(I,5)/MSTU(5)) | |
771 | ELSEIF(K(MOD(K(I,4),MSTU(5))+1,2).NE.22) THEN | |
772 | ID1=MOD(K(I,4),MSTU(5)) | |
773 | IF(K(I,2).GE.1.AND.K(I,2).LE.8) ID1=MOD(K(I,4),MSTU(5))+1 | |
774 | ID2=2*MOD(K(I,4),MSTU(5))+1-ID1 | |
775 | K(I,4)=MSTU(5)*(K(I,4)/MSTU(5))+ID1 | |
776 | K(I,5)=MSTU(5)*(K(I,5)/MSTU(5))+ID2 | |
777 | K(ID1,4)=K(ID1,4)+MSTU(5)*I | |
778 | K(ID1,5)=K(ID1,5)+MSTU(5)*ID2 | |
779 | K(ID2,4)=K(ID2,4)+MSTU(5)*ID1 | |
780 | K(ID2,5)=K(ID2,5)+MSTU(5)*I | |
781 | ELSE | |
782 | ID1=MOD(K(I,4),MSTU(5)) | |
783 | ID2=ID1+1 | |
784 | K(I,4)=MSTU(5)*(K(I,4)/MSTU(5))+ID1 | |
785 | K(I,5)=MSTU(5)*(K(I,5)/MSTU(5))+ID1 | |
786 | K(ID1,4)=K(ID1,4)+MSTU(5)*I | |
787 | K(ID1,5)=K(ID1,5)+MSTU(5)*I | |
788 | K(ID2,4)=0 | |
789 | K(ID2,5)=0 | |
790 | ENDIF | |
791 | 390 CONTINUE | |
792 | ||
793 | C...Transformation from CM frame. | |
794 | IF(NPA.GE.2) THEN | |
795 | BEX=PS(1)/PS(4) | |
796 | BEY=PS(2)/PS(4) | |
797 | BEZ=PS(3)/PS(4) | |
798 | GA=PS(4)/PS(5) | |
799 | GABEP=GA*(GA*(BEX*P(IPA(1),1)+BEY*P(IPA(1),2)+BEZ*P(IPA(1),3)) | |
800 | & /(1.+GA)-P(IPA(1),4)) | |
801 | ELSE | |
802 | BEX=0. | |
803 | BEY=0. | |
804 | BEZ=0. | |
805 | GABEP=0. | |
806 | ENDIF | |
807 | THE=ULANGL(P(IPA(1),3)+GABEP*BEZ,SQRT((P(IPA(1),1) | |
808 | &+GABEP*BEX)**2+(P(IPA(1),2)+GABEP*BEY)**2)) | |
809 | PHI=ULANGL(P(IPA(1),1)+GABEP*BEX,P(IPA(1),2)+GABEP*BEY) | |
810 | IF(NPA.EQ.3) THEN | |
811 | CHI=ULANGL(COS(THE)*COS(PHI)*(P(IPA(2),1)+GABEP*BEX)+COS(THE)* | |
812 | & SIN(PHI)*(P(IPA(2),2)+GABEP*BEY)-SIN(THE)*(P(IPA(2),3)+GABEP* | |
813 | & BEZ),-SIN(PHI)*(P(IPA(2),1)+GABEP*BEX)+COS(PHI)*(P(IPA(2),2)+ | |
814 | & GABEP*BEY)) | |
815 | MSTU(33)=1 | |
816 | CALL LUDBRB(NS+1,N,0.,CHI,0D0,0D0,0D0) | |
817 | ENDIF | |
818 | DBEX=DBLE(BEX) | |
819 | DBEY=DBLE(BEY) | |
820 | DBEZ=DBLE(BEZ) | |
821 | MSTU(33)=1 | |
822 | CALL LUDBRB(NS+1,N,THE,PHI,DBEX,DBEY,DBEZ) | |
823 | ||
824 | C...Decay vertex of shower. | |
825 | DO 400 I=NS+1,N | |
826 | DO 400 J=1,5 | |
827 | 400 V(I,J)=V(IP1,J) | |
828 | ||
829 | C...Delete trivial shower, else connect initiators. | |
830 | IF(N.EQ.NS+NPA+IIM) THEN | |
831 | N=NS | |
832 | ELSE | |
833 | DO 410 IP=1,NPA | |
834 | K(IPA(IP),1)=14 | |
835 | K(IPA(IP),4)=K(IPA(IP),4)+NS+IIM+IP | |
836 | K(IPA(IP),5)=K(IPA(IP),5)+NS+IIM+IP | |
837 | K(NS+IIM+IP,3)=IPA(IP) | |
838 | IF(IIM.EQ.1.AND.MSTU(16).NE.2) K(NS+IIM+IP,3)=NS+1 | |
839 | K(NS+IIM+IP,4)=MSTU(5)*IPA(IP)+K(NS+IIM+IP,4) | |
840 | 410 K(NS+IIM+IP,5)=MSTU(5)*IPA(IP)+K(NS+IIM+IP,5) | |
841 | ENDIF | |
842 | ||
843 | RETURN | |
844 | END |