]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PYTHIA/jetset/lutabu.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merging the VirtualMC branch to the main development branch (HEAD)
[u/mrichter/AliRoot.git] / PYTHIA / jetset / lutabu.F
CommitLineData
fe4da5cc 1
2C*********************************************************************
3
4 SUBROUTINE LUTABU(MTABU)
5
6C...Purpose: to evaluate various properties of an event, with
7C...statistics accumulated during the course of the run and
8C...printed at the end.
9 COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5)
10 COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
11 COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4)
12 COMMON/LUDAT3/MDCY(500,3),MDME(2000,2),BRAT(2000),KFDP(2000,5)
13 SAVE /LUJETS/,/LUDAT1/,/LUDAT2/,/LUDAT3/
14 DIMENSION KFIS(100,2),NPIS(100,0:10),KFFS(400),NPFS(400,4),
15 &FEVFM(10,4),FM1FM(3,10,4),FM2FM(3,10,4),FMOMA(4),FMOMS(4),
16 &FEVEE(50),FE1EC(50),FE2EC(50),FE1EA(25),FE2EA(25),
17 &KFDM(8),KFDC(200,0:8),NPDC(200)
18 SAVE NEVIS,NKFIS,KFIS,NPIS,NEVFS,NPRFS,NFIFS,NCHFS,NKFFS,
19 &KFFS,NPFS,NEVFM,NMUFM,FM1FM,FM2FM,NEVEE,FE1EC,FE2EC,FE1EA,
20 &FE2EA,NEVDC,NKFDC,NREDC,KFDC,NPDC
21 CHARACTER CHAU*16,CHIS(2)*12,CHDC(8)*12
22 DATA NEVIS/0/,NKFIS/0/,NEVFS/0/,NPRFS/0/,NFIFS/0/,NCHFS/0/,
23 &NKFFS/0/,NEVFM/0/,NMUFM/0/,FM1FM/120*0./,FM2FM/120*0./,
24 &NEVEE/0/,FE1EC/50*0./,FE2EC/50*0./,FE1EA/25*0./,FE2EA/25*0./,
25 &NEVDC/0/,NKFDC/0/,NREDC/0/
26
27C...Reset statistics on initial parton state.
28 IF(MTABU.EQ.10) THEN
29 NEVIS=0
30 NKFIS=0
31
32C...Identify and order flavour content of initial state.
33 ELSEIF(MTABU.EQ.11) THEN
34 NEVIS=NEVIS+1
35 KFM1=2*IABS(MSTU(161))
36 IF(MSTU(161).GT.0) KFM1=KFM1-1
37 KFM2=2*IABS(MSTU(162))
38 IF(MSTU(162).GT.0) KFM2=KFM2-1
39 KFMN=MIN(KFM1,KFM2)
40 KFMX=MAX(KFM1,KFM2)
41 DO 100 I=1,NKFIS
42 IF(KFMN.EQ.KFIS(I,1).AND.KFMX.EQ.KFIS(I,2)) THEN
43 IKFIS=-I
44 GOTO 110
45 ELSEIF(KFMN.LT.KFIS(I,1).OR.(KFMN.EQ.KFIS(I,1).AND.
46 & KFMX.LT.KFIS(I,2))) THEN
47 IKFIS=I
48 GOTO 110
49 ENDIF
50 100 CONTINUE
51 IKFIS=NKFIS+1
52 110 IF(IKFIS.LT.0) THEN
53 IKFIS=-IKFIS
54 ELSE
55 IF(NKFIS.GE.100) RETURN
56 DO 130 I=NKFIS,IKFIS,-1
57 KFIS(I+1,1)=KFIS(I,1)
58 KFIS(I+1,2)=KFIS(I,2)
59 DO 120 J=0,10
60 NPIS(I+1,J)=NPIS(I,J)
61 120 CONTINUE
62 130 CONTINUE
63 NKFIS=NKFIS+1
64 KFIS(IKFIS,1)=KFMN
65 KFIS(IKFIS,2)=KFMX
66 DO 140 J=0,10
67 NPIS(IKFIS,J)=0
68 140 CONTINUE
69 ENDIF
70 NPIS(IKFIS,0)=NPIS(IKFIS,0)+1
71
72C...Count number of partons in initial state.
73 NP=0
74 DO 160 I=1,N
75 IF(K(I,1).LE.0.OR.K(I,1).GT.12) THEN
76 ELSEIF(IABS(K(I,2)).GT.80.AND.IABS(K(I,2)).LE.100) THEN
77 ELSEIF(IABS(K(I,2)).GT.100.AND.MOD(IABS(K(I,2))/10,10).NE.0)
78 & THEN
79 ELSE
80 IM=I
81 150 IM=K(IM,3)
82 IF(IM.LE.0.OR.IM.GT.N) THEN
83 NP=NP+1
84 ELSEIF(K(IM,1).LE.0.OR.K(IM,1).GT.20) THEN
85 NP=NP+1
86 ELSEIF(IABS(K(IM,2)).GT.80.AND.IABS(K(IM,2)).LE.100) THEN
87 ELSEIF(IABS(K(IM,2)).GT.100.AND.MOD(IABS(K(IM,2))/10,10).NE.0)
88 & THEN
89 ELSE
90 GOTO 150
91 ENDIF
92 ENDIF
93 160 CONTINUE
94 NPCO=MAX(NP,1)
95 IF(NP.GE.6) NPCO=6
96 IF(NP.GE.8) NPCO=7
97 IF(NP.GE.11) NPCO=8
98 IF(NP.GE.16) NPCO=9
99 IF(NP.GE.26) NPCO=10
100 NPIS(IKFIS,NPCO)=NPIS(IKFIS,NPCO)+1
101 MSTU(62)=NP
102
103C...Write statistics on initial parton state.
104 ELSEIF(MTABU.EQ.12) THEN
105 FAC=1./MAX(1,NEVIS)
106 WRITE(MSTU(11),5000) NEVIS
107 DO 170 I=1,NKFIS
108 KFMN=KFIS(I,1)
109 IF(KFMN.EQ.0) KFMN=KFIS(I,2)
110 KFM1=(KFMN+1)/2
111 IF(2*KFM1.EQ.KFMN) KFM1=-KFM1
112 CALL LUNAME(KFM1,CHAU)
113 CHIS(1)=CHAU(1:12)
114 IF(CHAU(13:13).NE.' ') CHIS(1)(12:12)='?'
115 KFMX=KFIS(I,2)
116 IF(KFIS(I,1).EQ.0) KFMX=0
117 KFM2=(KFMX+1)/2
118 IF(2*KFM2.EQ.KFMX) KFM2=-KFM2
119 CALL LUNAME(KFM2,CHAU)
120 CHIS(2)=CHAU(1:12)
121 IF(CHAU(13:13).NE.' ') CHIS(2)(12:12)='?'
122 WRITE(MSTU(11),5100) CHIS(1),CHIS(2),FAC*NPIS(I,0),
123 & (NPIS(I,J)/FLOAT(NPIS(I,0)),J=1,10)
124 170 CONTINUE
125
126C...Copy statistics on initial parton state into /LUJETS/.
127 ELSEIF(MTABU.EQ.13) THEN
128 FAC=1./MAX(1,NEVIS)
129 DO 190 I=1,NKFIS
130 KFMN=KFIS(I,1)
131 IF(KFMN.EQ.0) KFMN=KFIS(I,2)
132 KFM1=(KFMN+1)/2
133 IF(2*KFM1.EQ.KFMN) KFM1=-KFM1
134 KFMX=KFIS(I,2)
135 IF(KFIS(I,1).EQ.0) KFMX=0
136 KFM2=(KFMX+1)/2
137 IF(2*KFM2.EQ.KFMX) KFM2=-KFM2
138 K(I,1)=32
139 K(I,2)=99
140 K(I,3)=KFM1
141 K(I,4)=KFM2
142 K(I,5)=NPIS(I,0)
143 DO 180 J=1,5
144 P(I,J)=FAC*NPIS(I,J)
145 V(I,J)=FAC*NPIS(I,J+5)
146 180 CONTINUE
147 190 CONTINUE
148 N=NKFIS
149 DO 200 J=1,5
150 K(N+1,J)=0
151 P(N+1,J)=0.
152 V(N+1,J)=0.
153 200 CONTINUE
154 K(N+1,1)=32
155 K(N+1,2)=99
156 K(N+1,5)=NEVIS
157 MSTU(3)=1
158
159C...Reset statistics on number of particles/partons.
160 ELSEIF(MTABU.EQ.20) THEN
161 NEVFS=0
162 NPRFS=0
163 NFIFS=0
164 NCHFS=0
165 NKFFS=0
166
167C...Identify whether particle/parton is primary or not.
168 ELSEIF(MTABU.EQ.21) THEN
169 NEVFS=NEVFS+1
170 MSTU(62)=0
171 DO 260 I=1,N
172 IF(K(I,1).LE.0.OR.K(I,1).GT.20.OR.K(I,1).EQ.13) GOTO 260
173 MSTU(62)=MSTU(62)+1
174 KC=LUCOMP(K(I,2))
175 MPRI=0
176 IF(K(I,3).LE.0.OR.K(I,3).GT.N) THEN
177 MPRI=1
178 ELSEIF(K(K(I,3),1).LE.0.OR.K(K(I,3),1).GT.20) THEN
179 MPRI=1
180 ELSEIF(K(K(I,3),2).GE.91.AND.K(K(I,3),2).LE.93) THEN
181 MPRI=1
182 ELSEIF(KC.EQ.0) THEN
183 ELSEIF(K(K(I,3),1).EQ.13) THEN
184 IM=K(K(I,3),3)
185 IF(IM.LE.0.OR.IM.GT.N) THEN
186 MPRI=1
187 ELSEIF(K(IM,1).LE.0.OR.K(IM,1).GT.20) THEN
188 MPRI=1
189 ENDIF
190 ELSEIF(KCHG(KC,2).EQ.0) THEN
191 KCM=LUCOMP(K(K(I,3),2))
192 IF(KCM.NE.0) THEN
193 IF(KCHG(KCM,2).NE.0) MPRI=1
194 ENDIF
195 ENDIF
196 IF(KC.NE.0.AND.MPRI.EQ.1) THEN
197 IF(KCHG(KC,2).EQ.0) NPRFS=NPRFS+1
198 ENDIF
199 IF(K(I,1).LE.10) THEN
200 NFIFS=NFIFS+1
201 IF(LUCHGE(K(I,2)).NE.0) NCHFS=NCHFS+1
202 ENDIF
203
204C...Fill statistics on number of particles/partons in event.
205 KFA=IABS(K(I,2))
206 KFS=3-ISIGN(1,K(I,2))-MPRI
207 DO 210 IP=1,NKFFS
208 IF(KFA.EQ.KFFS(IP)) THEN
209 IKFFS=-IP
210 GOTO 220
211 ELSEIF(KFA.LT.KFFS(IP)) THEN
212 IKFFS=IP
213 GOTO 220
214 ENDIF
215 210 CONTINUE
216 IKFFS=NKFFS+1
217 220 IF(IKFFS.LT.0) THEN
218 IKFFS=-IKFFS
219 ELSE
220 IF(NKFFS.GE.400) RETURN
221 DO 240 IP=NKFFS,IKFFS,-1
222 KFFS(IP+1)=KFFS(IP)
223 DO 230 J=1,4
224 NPFS(IP+1,J)=NPFS(IP,J)
225 230 CONTINUE
226 240 CONTINUE
227 NKFFS=NKFFS+1
228 KFFS(IKFFS)=KFA
229 DO 250 J=1,4
230 NPFS(IKFFS,J)=0
231 250 CONTINUE
232 ENDIF
233 NPFS(IKFFS,KFS)=NPFS(IKFFS,KFS)+1
234 260 CONTINUE
235
236C...Write statistics on particle/parton composition of events.
237 ELSEIF(MTABU.EQ.22) THEN
238 FAC=1./MAX(1,NEVFS)
239 WRITE(MSTU(11),5200) NEVFS,FAC*NPRFS,FAC*NFIFS,FAC*NCHFS
240 DO 270 I=1,NKFFS
241 CALL LUNAME(KFFS(I),CHAU)
242 KC=LUCOMP(KFFS(I))
243 MDCYF=0
244 IF(KC.NE.0) MDCYF=MDCY(KC,1)
245 WRITE(MSTU(11),5300) KFFS(I),CHAU,MDCYF,(FAC*NPFS(I,J),J=1,4),
246 & FAC*(NPFS(I,1)+NPFS(I,2)+NPFS(I,3)+NPFS(I,4))
247 270 CONTINUE
248
249C...Copy particle/parton composition information into /LUJETS/.
250 ELSEIF(MTABU.EQ.23) THEN
251 FAC=1./MAX(1,NEVFS)
252 DO 290 I=1,NKFFS
253 K(I,1)=32
254 K(I,2)=99
255 K(I,3)=KFFS(I)
256 K(I,4)=0
257 K(I,5)=NPFS(I,1)+NPFS(I,2)+NPFS(I,3)+NPFS(I,4)
258 DO 280 J=1,4
259 P(I,J)=FAC*NPFS(I,J)
260 V(I,J)=0.
261 280 CONTINUE
262 P(I,5)=FAC*K(I,5)
263 V(I,5)=0.
264 290 CONTINUE
265 N=NKFFS
266 DO 300 J=1,5
267 K(N+1,J)=0
268 P(N+1,J)=0.
269 V(N+1,J)=0.
270 300 CONTINUE
271 K(N+1,1)=32
272 K(N+1,2)=99
273 K(N+1,5)=NEVFS
274 P(N+1,1)=FAC*NPRFS
275 P(N+1,2)=FAC*NFIFS
276 P(N+1,3)=FAC*NCHFS
277 MSTU(3)=1
278
279C...Reset factorial moments statistics.
280 ELSEIF(MTABU.EQ.30) THEN
281 NEVFM=0
282 NMUFM=0
283 DO 330 IM=1,3
284 DO 320 IB=1,10
285 DO 310 IP=1,4
286 FM1FM(IM,IB,IP)=0.
287 FM2FM(IM,IB,IP)=0.
288 310 CONTINUE
289 320 CONTINUE
290 330 CONTINUE
291
292C...Find particles to include, with (pion,pseudo)rapidity and azimuth.
293 ELSEIF(MTABU.EQ.31) THEN
294 NEVFM=NEVFM+1
295 NLOW=N+MSTU(3)
296 NUPP=NLOW
297 DO 410 I=1,N
298 IF(K(I,1).LE.0.OR.K(I,1).GT.10) GOTO 410
299 IF(MSTU(41).GE.2) THEN
300 KC=LUCOMP(K(I,2))
301 IF(KC.EQ.0.OR.KC.EQ.12.OR.KC.EQ.14.OR.KC.EQ.16.OR.
302 & KC.EQ.18) GOTO 410
303 IF(MSTU(41).GE.3.AND.KCHG(KC,2).EQ.0.AND.LUCHGE(K(I,2)).EQ.0)
304 & GOTO 410
305 ENDIF
306 PMR=0.
307 IF(MSTU(42).EQ.1.AND.K(I,2).NE.22) PMR=ULMASS(211)
308 IF(MSTU(42).GE.2) PMR=P(I,5)
309 PR=MAX(1E-20,PMR**2+P(I,1)**2+P(I,2)**2)
310 YETA=SIGN(LOG(MIN((SQRT(PR+P(I,3)**2)+ABS(P(I,3)))/SQRT(PR),
311 & 1E20)),P(I,3))
312 IF(ABS(YETA).GT.PARU(57)) GOTO 410
313 PHI=ULANGL(P(I,1),P(I,2))
314 IYETA=512.*(YETA+PARU(57))/(2.*PARU(57))
315 IYETA=MAX(0,MIN(511,IYETA))
316 IPHI=512.*(PHI+PARU(1))/PARU(2)
317 IPHI=MAX(0,MIN(511,IPHI))
318 IYEP=0
319 DO 340 IB=0,9
320 IYEP=IYEP+4**IB*(2*MOD(IYETA/2**IB,2)+MOD(IPHI/2**IB,2))
321 340 CONTINUE
322
323C...Order particles in (pseudo)rapidity and/or azimuth.
324 IF(NUPP.GT.MSTU(4)-5-MSTU(32)) THEN
325 CALL LUERRM(11,'(LUTABU:) no more memory left in LUJETS')
326 RETURN
327 ENDIF
328 NUPP=NUPP+1
329 IF(NUPP.EQ.NLOW+1) THEN
330 K(NUPP,1)=IYETA
331 K(NUPP,2)=IPHI
332 K(NUPP,3)=IYEP
333 ELSE
334 DO 350 I1=NUPP-1,NLOW+1,-1
335 IF(IYETA.GE.K(I1,1)) GOTO 360
336 K(I1+1,1)=K(I1,1)
337 350 CONTINUE
338 360 K(I1+1,1)=IYETA
339 DO 370 I1=NUPP-1,NLOW+1,-1
340 IF(IPHI.GE.K(I1,2)) GOTO 380
341 K(I1+1,2)=K(I1,2)
342 370 CONTINUE
343 380 K(I1+1,2)=IPHI
344 DO 390 I1=NUPP-1,NLOW+1,-1
345 IF(IYEP.GE.K(I1,3)) GOTO 400
346 K(I1+1,3)=K(I1,3)
347 390 CONTINUE
348 400 K(I1+1,3)=IYEP
349 ENDIF
350 410 CONTINUE
351 K(NUPP+1,1)=2**10
352 K(NUPP+1,2)=2**10
353 K(NUPP+1,3)=4**10
354
355C...Calculate sum of factorial moments in event.
356 DO 480 IM=1,3
357 DO 430 IB=1,10
358 DO 420 IP=1,4
359 FEVFM(IB,IP)=0.
360 420 CONTINUE
361 430 CONTINUE
362 DO 450 IB=1,10
363 IF(IM.LE.2) IBIN=2**(10-IB)
364 IF(IM.EQ.3) IBIN=4**(10-IB)
365 IAGR=K(NLOW+1,IM)/IBIN
366 NAGR=1
367 DO 440 I=NLOW+2,NUPP+1
368 ICUT=K(I,IM)/IBIN
369 IF(ICUT.EQ.IAGR) THEN
370 NAGR=NAGR+1
371 ELSE
372 IF(NAGR.EQ.1) THEN
373 ELSEIF(NAGR.EQ.2) THEN
374 FEVFM(IB,1)=FEVFM(IB,1)+2.
375 ELSEIF(NAGR.EQ.3) THEN
376 FEVFM(IB,1)=FEVFM(IB,1)+6.
377 FEVFM(IB,2)=FEVFM(IB,2)+6.
378 ELSEIF(NAGR.EQ.4) THEN
379 FEVFM(IB,1)=FEVFM(IB,1)+12.
380 FEVFM(IB,2)=FEVFM(IB,2)+24.
381 FEVFM(IB,3)=FEVFM(IB,3)+24.
382 ELSE
383 FEVFM(IB,1)=FEVFM(IB,1)+NAGR*(NAGR-1.)
384 FEVFM(IB,2)=FEVFM(IB,2)+NAGR*(NAGR-1.)*(NAGR-2.)
385 FEVFM(IB,3)=FEVFM(IB,3)+NAGR*(NAGR-1.)*(NAGR-2.)*(NAGR-3.)
386 FEVFM(IB,4)=FEVFM(IB,4)+NAGR*(NAGR-1.)*(NAGR-2.)*(NAGR-3.)*
387 & (NAGR-4.)
388 ENDIF
389 IAGR=ICUT
390 NAGR=1
391 ENDIF
392 440 CONTINUE
393 450 CONTINUE
394
395C...Add results to total statistics.
396 DO 470 IB=10,1,-1
397 DO 460 IP=1,4
398 IF(FEVFM(1,IP).LT.0.5) THEN
399 FEVFM(IB,IP)=0.
400 ELSEIF(IM.LE.2) THEN
401 FEVFM(IB,IP)=2.**((IB-1)*IP)*FEVFM(IB,IP)/FEVFM(1,IP)
402 ELSE
403 FEVFM(IB,IP)=4.**((IB-1)*IP)*FEVFM(IB,IP)/FEVFM(1,IP)
404 ENDIF
405 FM1FM(IM,IB,IP)=FM1FM(IM,IB,IP)+FEVFM(IB,IP)
406 FM2FM(IM,IB,IP)=FM2FM(IM,IB,IP)+FEVFM(IB,IP)**2
407 460 CONTINUE
408 470 CONTINUE
409 480 CONTINUE
410 NMUFM=NMUFM+(NUPP-NLOW)
411 MSTU(62)=NUPP-NLOW
412
413C...Write accumulated statistics on factorial moments.
414 ELSEIF(MTABU.EQ.32) THEN
415 FAC=1./MAX(1,NEVFM)
416 IF(MSTU(42).LE.0) WRITE(MSTU(11),5400) NEVFM,'eta'
417 IF(MSTU(42).EQ.1) WRITE(MSTU(11),5400) NEVFM,'ypi'
418 IF(MSTU(42).GE.2) WRITE(MSTU(11),5400) NEVFM,'y '
419 DO 510 IM=1,3
420 WRITE(MSTU(11),5500)
421 DO 500 IB=1,10
422 BYETA=2.*PARU(57)
423 IF(IM.NE.2) BYETA=BYETA/2**(IB-1)
424 BPHI=PARU(2)
425 IF(IM.NE.1) BPHI=BPHI/2**(IB-1)
426 IF(IM.LE.2) BNAVE=FAC*NMUFM/FLOAT(2**(IB-1))
427 IF(IM.EQ.3) BNAVE=FAC*NMUFM/FLOAT(4**(IB-1))
428 DO 490 IP=1,4
429 FMOMA(IP)=FAC*FM1FM(IM,IB,IP)
430 FMOMS(IP)=SQRT(MAX(0.,FAC*(FAC*FM2FM(IM,IB,IP)-FMOMA(IP)**2)))
431 490 CONTINUE
432 WRITE(MSTU(11),5600) BYETA,BPHI,BNAVE,(FMOMA(IP),FMOMS(IP),
433 & IP=1,4)
434 500 CONTINUE
435 510 CONTINUE
436
437C...Copy statistics on factorial moments into /LUJETS/.
438 ELSEIF(MTABU.EQ.33) THEN
439 FAC=1./MAX(1,NEVFM)
440 DO 540 IM=1,3
441 DO 530 IB=1,10
442 I=10*(IM-1)+IB
443 K(I,1)=32
444 K(I,2)=99
445 K(I,3)=1
446 IF(IM.NE.2) K(I,3)=2**(IB-1)
447 K(I,4)=1
448 IF(IM.NE.1) K(I,4)=2**(IB-1)
449 K(I,5)=0
450 P(I,1)=2.*PARU(57)/K(I,3)
451 V(I,1)=PARU(2)/K(I,4)
452 DO 520 IP=1,4
453 P(I,IP+1)=FAC*FM1FM(IM,IB,IP)
454 V(I,IP+1)=SQRT(MAX(0.,FAC*(FAC*FM2FM(IM,IB,IP)-P(I,IP+1)**2)))
455 520 CONTINUE
456 530 CONTINUE
457 540 CONTINUE
458 N=30
459 DO 550 J=1,5
460 K(N+1,J)=0
461 P(N+1,J)=0.
462 V(N+1,J)=0.
463 550 CONTINUE
464 K(N+1,1)=32
465 K(N+1,2)=99
466 K(N+1,5)=NEVFM
467 MSTU(3)=1
468
469C...Reset statistics on Energy-Energy Correlation.
470 ELSEIF(MTABU.EQ.40) THEN
471 NEVEE=0
472 DO 560 J=1,25
473 FE1EC(J)=0.
474 FE2EC(J)=0.
475 FE1EC(51-J)=0.
476 FE2EC(51-J)=0.
477 FE1EA(J)=0.
478 FE2EA(J)=0.
479 560 CONTINUE
480
481C...Find particles to include, with proper assumed mass.
482 ELSEIF(MTABU.EQ.41) THEN
483 NEVEE=NEVEE+1
484 NLOW=N+MSTU(3)
485 NUPP=NLOW
486 ECM=0.
487 DO 570 I=1,N
488 IF(K(I,1).LE.0.OR.K(I,1).GT.10) GOTO 570
489 IF(MSTU(41).GE.2) THEN
490 KC=LUCOMP(K(I,2))
491 IF(KC.EQ.0.OR.KC.EQ.12.OR.KC.EQ.14.OR.KC.EQ.16.OR.
492 & KC.EQ.18) GOTO 570
493 IF(MSTU(41).GE.3.AND.KCHG(KC,2).EQ.0.AND.LUCHGE(K(I,2)).EQ.0)
494 & GOTO 570
495 ENDIF
496 PMR=0.
497 IF(MSTU(42).EQ.1.AND.K(I,2).NE.22) PMR=ULMASS(211)
498 IF(MSTU(42).GE.2) PMR=P(I,5)
499 IF(NUPP.GT.MSTU(4)-5-MSTU(32)) THEN
500 CALL LUERRM(11,'(LUTABU:) no more memory left in LUJETS')
501 RETURN
502 ENDIF
503 NUPP=NUPP+1
504 P(NUPP,1)=P(I,1)
505 P(NUPP,2)=P(I,2)
506 P(NUPP,3)=P(I,3)
507 P(NUPP,4)=SQRT(PMR**2+P(I,1)**2+P(I,2)**2+P(I,3)**2)
508 P(NUPP,5)=MAX(1E-10,SQRT(P(I,1)**2+P(I,2)**2+P(I,3)**2))
509 ECM=ECM+P(NUPP,4)
510 570 CONTINUE
511 IF(NUPP.EQ.NLOW) RETURN
512
513C...Analyze Energy-Energy Correlation in event.
514 FAC=(2./ECM**2)*50./PARU(1)
515 DO 580 J=1,50
516 FEVEE(J)=0.
517 580 CONTINUE
518 DO 600 I1=NLOW+2,NUPP
519 DO 590 I2=NLOW+1,I1-1
520 CTHE=(P(I1,1)*P(I2,1)+P(I1,2)*P(I2,2)+P(I1,3)*P(I2,3))/
521 & (P(I1,5)*P(I2,5))
522 THE=ACOS(MAX(-1.,MIN(1.,CTHE)))
523 ITHE=MAX(1,MIN(50,1+INT(50.*THE/PARU(1))))
524 FEVEE(ITHE)=FEVEE(ITHE)+FAC*P(I1,4)*P(I2,4)
525 590 CONTINUE
526 600 CONTINUE
527 DO 610 J=1,25
528 FE1EC(J)=FE1EC(J)+FEVEE(J)
529 FE2EC(J)=FE2EC(J)+FEVEE(J)**2
530 FE1EC(51-J)=FE1EC(51-J)+FEVEE(51-J)
531 FE2EC(51-J)=FE2EC(51-J)+FEVEE(51-J)**2
532 FE1EA(J)=FE1EA(J)+(FEVEE(51-J)-FEVEE(J))
533 FE2EA(J)=FE2EA(J)+(FEVEE(51-J)-FEVEE(J))**2
534 610 CONTINUE
535 MSTU(62)=NUPP-NLOW
536
537C...Write statistics on Energy-Energy Correlation.
538 ELSEIF(MTABU.EQ.42) THEN
539 FAC=1./MAX(1,NEVEE)
540 WRITE(MSTU(11),5700) NEVEE
541 DO 620 J=1,25
542 FEEC1=FAC*FE1EC(J)
543 FEES1=SQRT(MAX(0.,FAC*(FAC*FE2EC(J)-FEEC1**2)))
544 FEEC2=FAC*FE1EC(51-J)
545 FEES2=SQRT(MAX(0.,FAC*(FAC*FE2EC(51-J)-FEEC2**2)))
546 FEECA=FAC*FE1EA(J)
547 FEESA=SQRT(MAX(0.,FAC*(FAC*FE2EA(J)-FEECA**2)))
548 WRITE(MSTU(11),5800) 3.6*(J-1),3.6*J,FEEC1,FEES1,FEEC2,FEES2,
549 & FEECA,FEESA
550 620 CONTINUE
551
552C...Copy statistics on Energy-Energy Correlation into /LUJETS/.
553 ELSEIF(MTABU.EQ.43) THEN
554 FAC=1./MAX(1,NEVEE)
555 DO 630 I=1,25
556 K(I,1)=32
557 K(I,2)=99
558 K(I,3)=0
559 K(I,4)=0
560 K(I,5)=0
561 P(I,1)=FAC*FE1EC(I)
562 V(I,1)=SQRT(MAX(0.,FAC*(FAC*FE2EC(I)-P(I,1)**2)))
563 P(I,2)=FAC*FE1EC(51-I)
564 V(I,2)=SQRT(MAX(0.,FAC*(FAC*FE2EC(51-I)-P(I,2)**2)))
565 P(I,3)=FAC*FE1EA(I)
566 V(I,3)=SQRT(MAX(0.,FAC*(FAC*FE2EA(I)-P(I,3)**2)))
567 P(I,4)=PARU(1)*(I-1)/50.
568 P(I,5)=PARU(1)*I/50.
569 V(I,4)=3.6*(I-1)
570 V(I,5)=3.6*I
571 630 CONTINUE
572 N=25
573 DO 640 J=1,5
574 K(N+1,J)=0
575 P(N+1,J)=0.
576 V(N+1,J)=0.
577 640 CONTINUE
578 K(N+1,1)=32
579 K(N+1,2)=99
580 K(N+1,5)=NEVEE
581 MSTU(3)=1
582
583C...Reset statistics on decay channels.
584 ELSEIF(MTABU.EQ.50) THEN
585 NEVDC=0
586 NKFDC=0
587 NREDC=0
588
589C...Identify and order flavour content of final state.
590 ELSEIF(MTABU.EQ.51) THEN
591 NEVDC=NEVDC+1
592 NDS=0
593 DO 670 I=1,N
594 IF(K(I,1).LE.0.OR.K(I,1).GE.6) GOTO 670
595 NDS=NDS+1
596 IF(NDS.GT.8) THEN
597 NREDC=NREDC+1
598 RETURN
599 ENDIF
600 KFM=2*IABS(K(I,2))
601 IF(K(I,2).LT.0) KFM=KFM-1
602 DO 650 IDS=NDS-1,1,-1
603 IIN=IDS+1
604 IF(KFM.LT.KFDM(IDS)) GOTO 660
605 KFDM(IDS+1)=KFDM(IDS)
606 650 CONTINUE
607 IIN=1
608 660 KFDM(IIN)=KFM
609 670 CONTINUE
610
611C...Find whether old or new final state.
612 DO 690 IDC=1,NKFDC
613 IF(NDS.LT.KFDC(IDC,0)) THEN
614 IKFDC=IDC
615 GOTO 700
616 ELSEIF(NDS.EQ.KFDC(IDC,0)) THEN
617 DO 680 I=1,NDS
618 IF(KFDM(I).LT.KFDC(IDC,I)) THEN
619 IKFDC=IDC
620 GOTO 700
621 ELSEIF(KFDM(I).GT.KFDC(IDC,I)) THEN
622 GOTO 690
623 ENDIF
624 680 CONTINUE
625 IKFDC=-IDC
626 GOTO 700
627 ENDIF
628 690 CONTINUE
629 IKFDC=NKFDC+1
630 700 IF(IKFDC.LT.0) THEN
631 IKFDC=-IKFDC
632 ELSEIF(NKFDC.GE.200) THEN
633 NREDC=NREDC+1
634 RETURN
635 ELSE
636 DO 720 IDC=NKFDC,IKFDC,-1
637 NPDC(IDC+1)=NPDC(IDC)
638 DO 710 I=0,8
639 KFDC(IDC+1,I)=KFDC(IDC,I)
640 710 CONTINUE
641 720 CONTINUE
642 NKFDC=NKFDC+1
643 KFDC(IKFDC,0)=NDS
644 DO 730 I=1,NDS
645 KFDC(IKFDC,I)=KFDM(I)
646 730 CONTINUE
647 NPDC(IKFDC)=0
648 ENDIF
649 NPDC(IKFDC)=NPDC(IKFDC)+1
650
651C...Write statistics on decay channels.
652 ELSEIF(MTABU.EQ.52) THEN
653 FAC=1./MAX(1,NEVDC)
654 WRITE(MSTU(11),5900) NEVDC
655 DO 750 IDC=1,NKFDC
656 DO 740 I=1,KFDC(IDC,0)
657 KFM=KFDC(IDC,I)
658 KF=(KFM+1)/2
659 IF(2*KF.NE.KFM) KF=-KF
660 CALL LUNAME(KF,CHAU)
661 CHDC(I)=CHAU(1:12)
662 IF(CHAU(13:13).NE.' ') CHDC(I)(12:12)='?'
663 740 CONTINUE
664 WRITE(MSTU(11),6000) FAC*NPDC(IDC),(CHDC(I),I=1,KFDC(IDC,0))
665 750 CONTINUE
666 IF(NREDC.NE.0) WRITE(MSTU(11),6100) FAC*NREDC
667
668C...Copy statistics on decay channels into /LUJETS/.
669 ELSEIF(MTABU.EQ.53) THEN
670 FAC=1./MAX(1,NEVDC)
671 DO 780 IDC=1,NKFDC
672 K(IDC,1)=32
673 K(IDC,2)=99
674 K(IDC,3)=0
675 K(IDC,4)=0
676 K(IDC,5)=KFDC(IDC,0)
677 DO 760 J=1,5
678 P(IDC,J)=0.
679 V(IDC,J)=0.
680 760 CONTINUE
681 DO 770 I=1,KFDC(IDC,0)
682 KFM=KFDC(IDC,I)
683 KF=(KFM+1)/2
684 IF(2*KF.NE.KFM) KF=-KF
685 IF(I.LE.5) P(IDC,I)=KF
686 IF(I.GE.6) V(IDC,I-5)=KF
687 770 CONTINUE
688 V(IDC,5)=FAC*NPDC(IDC)
689 780 CONTINUE
690 N=NKFDC
691 DO 790 J=1,5
692 K(N+1,J)=0
693 P(N+1,J)=0.
694 V(N+1,J)=0.
695 790 CONTINUE
696 K(N+1,1)=32
697 K(N+1,2)=99
698 K(N+1,5)=NEVDC
699 V(N+1,5)=FAC*NREDC
700 MSTU(3)=1
701 ENDIF
702
703C...Format statements for output on unit MSTU(11) (default 6).
704 5000 FORMAT(///20X,'Event statistics - initial state'/
705 &20X,'based on an analysis of ',I6,' events'//
706 &3X,'Main flavours after',8X,'Fraction',4X,'Subfractions ',
707 &'according to fragmenting system multiplicity'/
708 &4X,'hard interaction',24X,'1',7X,'2',7X,'3',7X,'4',7X,'5',
709 &6X,'6-7',5X,'8-10',3X,'11-15',3X,'16-25',4X,'>25'/)
710 5100 FORMAT(3X,A12,1X,A12,F10.5,1X,10F8.4)
711 5200 FORMAT(///20X,'Event statistics - final state'/
712 &20X,'based on an analysis of ',I7,' events'//
713 &5X,'Mean primary multiplicity =',F10.4/
714 &5X,'Mean final multiplicity =',F10.4/
715 &5X,'Mean charged multiplicity =',F10.4//
716 &5X,'Number of particles produced per event (directly and via ',
717 &'decays/branchings)'/
718 &5X,'KF Particle/jet MDCY',10X,'Particles',13X,'Antiparticles',
719 &8X,'Total'/35X,'prim seco prim seco'/)
720 5300 FORMAT(1X,I6,4X,A16,I2,5(1X,F11.6))
721 5400 FORMAT(///20X,'Factorial moments analysis of multiplicity'/
722 &20X,'based on an analysis of ',I6,' events'//
723 &3X,'delta-',A3,' delta-phi <n>/bin',10X,'<F2>',18X,'<F3>',
724 &18X,'<F4>',18X,'<F5>'/35X,4(' value error '))
725 5500 FORMAT(10X)
726 5600 FORMAT(2X,2F10.4,F12.4,4(F12.4,F10.4))
727 5700 FORMAT(///20X,'Energy-Energy Correlation and Asymmetry'/
728 &20X,'based on an analysis of ',I6,' events'//
729 &2X,'theta range',8X,'EEC(theta)',8X,'EEC(180-theta)',7X,
730 &'EECA(theta)'/2X,'in degrees ',3(' value error')/)
731 5800 FORMAT(2X,F4.1,' - ',F4.1,3(F11.4,F9.4))
732 5900 FORMAT(///20X,'Decay channel analysis - final state'/
733 &20X,'based on an analysis of ',I6,' events'//
734 &2X,'Probability',10X,'Complete final state'/)
735 6000 FORMAT(2X,F9.5,5X,8(A12,1X))
736 6100 FORMAT(2X,F9.5,5X,'into other channels (more than 8 particles ',
737 &'or table overflow)')
738
739 RETURN
740 END