C********************************************************************* C...PYTABU C...Evaluates various properties of an event, with statistics C...accumulated during the course of the run and C...printed at the end. SUBROUTINE PYTABU(MTABU) C...Double precision and integer declarations. IMPLICIT DOUBLE PRECISION(A-H, O-Z) IMPLICIT INTEGER(I-N) INTEGER PYK,PYCHGE,PYCOMP C...Commonblocks. COMMON/PYJETS/N,NPAD,K(4000,5),P(4000,5),V(4000,5) COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200) COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4) COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5) SAVE /PYJETS/,/PYDAT1/,/PYDAT2/,/PYDAT3/ C...Local arrays, character variables, saved variables and data. DIMENSION KFIS(100,2),NPIS(100,0:10),KFFS(400),NPFS(400,4), &FEVFM(10,4),FM1FM(3,10,4),FM2FM(3,10,4),FMOMA(4),FMOMS(4), &FEVEE(50),FE1EC(50),FE2EC(50),FE1EA(25),FE2EA(25), &KFDM(8),KFDC(200,0:8),NPDC(200) SAVE NEVIS,NKFIS,KFIS,NPIS,NEVFS,NPRFS,NFIFS,NCHFS,NKFFS, &KFFS,NPFS,NEVFM,NMUFM,FM1FM,FM2FM,NEVEE,FE1EC,FE2EC,FE1EA, &FE2EA,NEVDC,NKFDC,NREDC,KFDC,NPDC CHARACTER CHAU*16,CHIS(2)*12,CHDC(8)*12 DATA NEVIS/0/,NKFIS/0/,NEVFS/0/,NPRFS/0/,NFIFS/0/,NCHFS/0/, &NKFFS/0/,NEVFM/0/,NMUFM/0/,FM1FM/120*0D0/,FM2FM/120*0D0/, &NEVEE/0/,FE1EC/50*0D0/,FE2EC/50*0D0/,FE1EA/25*0D0/,FE2EA/25*0D0/, &NEVDC/0/,NKFDC/0/,NREDC/0/ C...Reset statistics on initial parton state. IF(MTABU.EQ.10) THEN NEVIS=0 NKFIS=0 C...Identify and order flavour content of initial state. ELSEIF(MTABU.EQ.11) THEN NEVIS=NEVIS+1 KFM1=2*IABS(MSTU(161)) IF(MSTU(161).GT.0) KFM1=KFM1-1 KFM2=2*IABS(MSTU(162)) IF(MSTU(162).GT.0) KFM2=KFM2-1 KFMN=MIN(KFM1,KFM2) KFMX=MAX(KFM1,KFM2) DO 100 I=1,NKFIS IF(KFMN.EQ.KFIS(I,1).AND.KFMX.EQ.KFIS(I,2)) THEN IKFIS=-I GOTO 110 ELSEIF(KFMN.LT.KFIS(I,1).OR.(KFMN.EQ.KFIS(I,1).AND. & KFMX.LT.KFIS(I,2))) THEN IKFIS=I GOTO 110 ENDIF 100 CONTINUE IKFIS=NKFIS+1 110 IF(IKFIS.LT.0) THEN IKFIS=-IKFIS ELSE IF(NKFIS.GE.100) RETURN DO 130 I=NKFIS,IKFIS,-1 KFIS(I+1,1)=KFIS(I,1) KFIS(I+1,2)=KFIS(I,2) DO 120 J=0,10 NPIS(I+1,J)=NPIS(I,J) 120 CONTINUE 130 CONTINUE NKFIS=NKFIS+1 KFIS(IKFIS,1)=KFMN KFIS(IKFIS,2)=KFMX DO 140 J=0,10 NPIS(IKFIS,J)=0 140 CONTINUE ENDIF NPIS(IKFIS,0)=NPIS(IKFIS,0)+1 C...Count number of partons in initial state. NP=0 DO 160 I=1,N IF(K(I,1).LE.0.OR.K(I,1).GT.12) THEN ELSEIF(IABS(K(I,2)).GT.80.AND.IABS(K(I,2)).LE.100) THEN ELSEIF(IABS(K(I,2)).GT.100.AND.MOD(IABS(K(I,2))/10,10).NE.0) & THEN ELSE IM=I 150 IM=K(IM,3) IF(IM.LE.0.OR.IM.GT.N) THEN NP=NP+1 ELSEIF(K(IM,1).LE.0.OR.K(IM,1).GT.20) THEN NP=NP+1 ELSEIF(IABS(K(IM,2)).GT.80.AND.IABS(K(IM,2)).LE.100) THEN ELSEIF(IABS(K(IM,2)).GT.100.AND.MOD(IABS(K(IM,2))/10,10) & .NE.0) THEN ELSE GOTO 150 ENDIF ENDIF 160 CONTINUE NPCO=MAX(NP,1) IF(NP.GE.6) NPCO=6 IF(NP.GE.8) NPCO=7 IF(NP.GE.11) NPCO=8 IF(NP.GE.16) NPCO=9 IF(NP.GE.26) NPCO=10 NPIS(IKFIS,NPCO)=NPIS(IKFIS,NPCO)+1 MSTU(62)=NP C...Write statistics on initial parton state. ELSEIF(MTABU.EQ.12) THEN FAC=1D0/MAX(1,NEVIS) WRITE(MSTU(11),5000) NEVIS DO 170 I=1,NKFIS KFMN=KFIS(I,1) IF(KFMN.EQ.0) KFMN=KFIS(I,2) KFM1=(KFMN+1)/2 IF(2*KFM1.EQ.KFMN) KFM1=-KFM1 CALL PYNAME(KFM1,CHAU) CHIS(1)=CHAU(1:12) IF(CHAU(13:13).NE.' ') CHIS(1)(12:12)='?' KFMX=KFIS(I,2) IF(KFIS(I,1).EQ.0) KFMX=0 KFM2=(KFMX+1)/2 IF(2*KFM2.EQ.KFMX) KFM2=-KFM2 CALL PYNAME(KFM2,CHAU) CHIS(2)=CHAU(1:12) IF(CHAU(13:13).NE.' ') CHIS(2)(12:12)='?' WRITE(MSTU(11),5100) CHIS(1),CHIS(2),FAC*NPIS(I,0), & (NPIS(I,J)/DBLE(NPIS(I,0)),J=1,10) 170 CONTINUE C...Copy statistics on initial parton state into /PYJETS/. ELSEIF(MTABU.EQ.13) THEN FAC=1D0/MAX(1,NEVIS) DO 190 I=1,NKFIS KFMN=KFIS(I,1) IF(KFMN.EQ.0) KFMN=KFIS(I,2) KFM1=(KFMN+1)/2 IF(2*KFM1.EQ.KFMN) KFM1=-KFM1 KFMX=KFIS(I,2) IF(KFIS(I,1).EQ.0) KFMX=0 KFM2=(KFMX+1)/2 IF(2*KFM2.EQ.KFMX) KFM2=-KFM2 K(I,1)=32 K(I,2)=99 K(I,3)=KFM1 K(I,4)=KFM2 K(I,5)=NPIS(I,0) DO 180 J=1,5 P(I,J)=FAC*NPIS(I,J) V(I,J)=FAC*NPIS(I,J+5) 180 CONTINUE 190 CONTINUE N=NKFIS DO 200 J=1,5 K(N+1,J)=0 P(N+1,J)=0D0 V(N+1,J)=0D0 200 CONTINUE K(N+1,1)=32 K(N+1,2)=99 K(N+1,5)=NEVIS MSTU(3)=1 C...Reset statistics on number of particles/partons. ELSEIF(MTABU.EQ.20) THEN NEVFS=0 NPRFS=0 NFIFS=0 NCHFS=0 NKFFS=0 C...Identify whether particle/parton is primary or not. ELSEIF(MTABU.EQ.21) THEN NEVFS=NEVFS+1 MSTU(62)=0 DO 260 I=1,N IF(K(I,1).LE.0.OR.K(I,1).GT.20.OR.K(I,1).EQ.13) GOTO 260 MSTU(62)=MSTU(62)+1 KC=PYCOMP(K(I,2)) MPRI=0 IF(K(I,3).LE.0.OR.K(I,3).GT.N) THEN MPRI=1 ELSEIF(K(K(I,3),1).LE.0.OR.K(K(I,3),1).GT.20) THEN MPRI=1 ELSEIF(K(K(I,3),2).GE.91.AND.K(K(I,3),2).LE.93) THEN MPRI=1 ELSEIF(KC.EQ.0) THEN ELSEIF(K(K(I,3),1).EQ.13) THEN IM=K(K(I,3),3) IF(IM.LE.0.OR.IM.GT.N) THEN MPRI=1 ELSEIF(K(IM,1).LE.0.OR.K(IM,1).GT.20) THEN MPRI=1 ENDIF ELSEIF(KCHG(KC,2).EQ.0) THEN KCM=PYCOMP(K(K(I,3),2)) IF(KCM.NE.0) THEN IF(KCHG(KCM,2).NE.0) MPRI=1 ENDIF ENDIF IF(KC.NE.0.AND.MPRI.EQ.1) THEN IF(KCHG(KC,2).EQ.0) NPRFS=NPRFS+1 ENDIF IF(K(I,1).LE.10) THEN NFIFS=NFIFS+1 IF(PYCHGE(K(I,2)).NE.0) NCHFS=NCHFS+1 ENDIF C...Fill statistics on number of particles/partons in event. KFA=IABS(K(I,2)) KFS=3-ISIGN(1,K(I,2))-MPRI DO 210 IP=1,NKFFS IF(KFA.EQ.KFFS(IP)) THEN IKFFS=-IP GOTO 220 ELSEIF(KFA.LT.KFFS(IP)) THEN IKFFS=IP GOTO 220 ENDIF 210 CONTINUE IKFFS=NKFFS+1 220 IF(IKFFS.LT.0) THEN IKFFS=-IKFFS ELSE IF(NKFFS.GE.400) RETURN DO 240 IP=NKFFS,IKFFS,-1 KFFS(IP+1)=KFFS(IP) DO 230 J=1,4 NPFS(IP+1,J)=NPFS(IP,J) 230 CONTINUE 240 CONTINUE NKFFS=NKFFS+1 KFFS(IKFFS)=KFA DO 250 J=1,4 NPFS(IKFFS,J)=0 250 CONTINUE ENDIF NPFS(IKFFS,KFS)=NPFS(IKFFS,KFS)+1 260 CONTINUE C...Write statistics on particle/parton composition of events. ELSEIF(MTABU.EQ.22) THEN FAC=1D0/MAX(1,NEVFS) WRITE(MSTU(11),5200) NEVFS,FAC*NPRFS,FAC*NFIFS,FAC*NCHFS DO 270 I=1,NKFFS CALL PYNAME(KFFS(I),CHAU) KC=PYCOMP(KFFS(I)) MDCYF=0 IF(KC.NE.0) MDCYF=MDCY(KC,1) WRITE(MSTU(11),5300) KFFS(I),CHAU,MDCYF,(FAC*NPFS(I,J),J=1,4), & FAC*(NPFS(I,1)+NPFS(I,2)+NPFS(I,3)+NPFS(I,4)) 270 CONTINUE C...Copy particle/parton composition information into /PYJETS/. ELSEIF(MTABU.EQ.23) THEN FAC=1D0/MAX(1,NEVFS) DO 290 I=1,NKFFS K(I,1)=32 K(I,2)=99 K(I,3)=KFFS(I) K(I,4)=0 K(I,5)=NPFS(I,1)+NPFS(I,2)+NPFS(I,3)+NPFS(I,4) DO 280 J=1,4 P(I,J)=FAC*NPFS(I,J) V(I,J)=0D0 280 CONTINUE P(I,5)=FAC*K(I,5) V(I,5)=0D0 290 CONTINUE N=NKFFS DO 300 J=1,5 K(N+1,J)=0 P(N+1,J)=0D0 V(N+1,J)=0D0 300 CONTINUE K(N+1,1)=32 K(N+1,2)=99 K(N+1,5)=NEVFS P(N+1,1)=FAC*NPRFS P(N+1,2)=FAC*NFIFS P(N+1,3)=FAC*NCHFS MSTU(3)=1 C...Reset factorial moments statistics. ELSEIF(MTABU.EQ.30) THEN NEVFM=0 NMUFM=0 DO 330 IM=1,3 DO 320 IB=1,10 DO 310 IP=1,4 FM1FM(IM,IB,IP)=0D0 FM2FM(IM,IB,IP)=0D0 310 CONTINUE 320 CONTINUE 330 CONTINUE C...Find particles to include, with (pion,pseudo)rapidity and azimuth. ELSEIF(MTABU.EQ.31) THEN NEVFM=NEVFM+1 NLOW=N+MSTU(3) NUPP=NLOW DO 410 I=1,N IF(K(I,1).LE.0.OR.K(I,1).GT.10) GOTO 410 IF(MSTU(41).GE.2) THEN KC=PYCOMP(K(I,2)) IF(KC.EQ.0.OR.KC.EQ.12.OR.KC.EQ.14.OR.KC.EQ.16.OR. & KC.EQ.18) GOTO 410 IF(MSTU(41).GE.3.AND.KCHG(KC,2).EQ.0.AND. & PYCHGE(K(I,2)).EQ.0) GOTO 410 ENDIF PMR=0D0 IF(MSTU(42).EQ.1.AND.K(I,2).NE.22) PMR=PYMASS(211) IF(MSTU(42).GE.2) PMR=P(I,5) PR=MAX(1D-20,PMR**2+P(I,1)**2+P(I,2)**2) YETA=SIGN(LOG(MIN((SQRT(PR+P(I,3)**2)+ABS(P(I,3)))/SQRT(PR), & 1D20)),P(I,3)) IF(ABS(YETA).GT.PARU(57)) GOTO 410 PHI=PYANGL(P(I,1),P(I,2)) IYETA=512D0*(YETA+PARU(57))/(2D0*PARU(57)) IYETA=MAX(0,MIN(511,IYETA)) IPHI=512D0*(PHI+PARU(1))/PARU(2) IPHI=MAX(0,MIN(511,IPHI)) IYEP=0 DO 340 IB=0,9 IYEP=IYEP+4**IB*(2*MOD(IYETA/2**IB,2)+MOD(IPHI/2**IB,2)) 340 CONTINUE C...Order particles in (pseudo)rapidity and/or azimuth. IF(NUPP.GT.MSTU(4)-5-MSTU(32)) THEN CALL PYERRM(11,'(PYTABU:) no more memory left in PYJETS') RETURN ENDIF NUPP=NUPP+1 IF(NUPP.EQ.NLOW+1) THEN K(NUPP,1)=IYETA K(NUPP,2)=IPHI K(NUPP,3)=IYEP ELSE DO 350 I1=NUPP-1,NLOW+1,-1 IF(IYETA.GE.K(I1,1)) GOTO 360 K(I1+1,1)=K(I1,1) 350 CONTINUE 360 K(I1+1,1)=IYETA DO 370 I1=NUPP-1,NLOW+1,-1 IF(IPHI.GE.K(I1,2)) GOTO 380 K(I1+1,2)=K(I1,2) 370 CONTINUE 380 K(I1+1,2)=IPHI DO 390 I1=NUPP-1,NLOW+1,-1 IF(IYEP.GE.K(I1,3)) GOTO 400 K(I1+1,3)=K(I1,3) 390 CONTINUE 400 K(I1+1,3)=IYEP ENDIF 410 CONTINUE K(NUPP+1,1)=2**10 K(NUPP+1,2)=2**10 K(NUPP+1,3)=4**10 C...Calculate sum of factorial moments in event. DO 480 IM=1,3 DO 430 IB=1,10 DO 420 IP=1,4 FEVFM(IB,IP)=0D0 420 CONTINUE 430 CONTINUE DO 450 IB=1,10 IF(IM.LE.2) IBIN=2**(10-IB) IF(IM.EQ.3) IBIN=4**(10-IB) IAGR=K(NLOW+1,IM)/IBIN NAGR=1 DO 440 I=NLOW+2,NUPP+1 ICUT=K(I,IM)/IBIN IF(ICUT.EQ.IAGR) THEN NAGR=NAGR+1 ELSE IF(NAGR.EQ.1) THEN ELSEIF(NAGR.EQ.2) THEN FEVFM(IB,1)=FEVFM(IB,1)+2D0 ELSEIF(NAGR.EQ.3) THEN FEVFM(IB,1)=FEVFM(IB,1)+6D0 FEVFM(IB,2)=FEVFM(IB,2)+6D0 ELSEIF(NAGR.EQ.4) THEN FEVFM(IB,1)=FEVFM(IB,1)+12D0 FEVFM(IB,2)=FEVFM(IB,2)+24D0 FEVFM(IB,3)=FEVFM(IB,3)+24D0 ELSE FEVFM(IB,1)=FEVFM(IB,1)+NAGR*(NAGR-1D0) FEVFM(IB,2)=FEVFM(IB,2)+NAGR*(NAGR-1D0)*(NAGR-2D0) FEVFM(IB,3)=FEVFM(IB,3)+NAGR*(NAGR-1D0)*(NAGR-2D0)* & (NAGR-3D0) FEVFM(IB,4)=FEVFM(IB,4)+NAGR*(NAGR-1D0)*(NAGR-2D0)* & (NAGR-3D0)*(NAGR-4D0) ENDIF IAGR=ICUT NAGR=1 ENDIF 440 CONTINUE 450 CONTINUE C...Add results to total statistics. DO 470 IB=10,1,-1 DO 460 IP=1,4 IF(FEVFM(1,IP).LT.0.5D0) THEN FEVFM(IB,IP)=0D0 ELSEIF(IM.LE.2) THEN FEVFM(IB,IP)=2D0**((IB-1)*IP)*FEVFM(IB,IP)/FEVFM(1,IP) ELSE FEVFM(IB,IP)=4D0**((IB-1)*IP)*FEVFM(IB,IP)/FEVFM(1,IP) ENDIF FM1FM(IM,IB,IP)=FM1FM(IM,IB,IP)+FEVFM(IB,IP) FM2FM(IM,IB,IP)=FM2FM(IM,IB,IP)+FEVFM(IB,IP)**2 460 CONTINUE 470 CONTINUE 480 CONTINUE NMUFM=NMUFM+(NUPP-NLOW) MSTU(62)=NUPP-NLOW C...Write accumulated statistics on factorial moments. ELSEIF(MTABU.EQ.32) THEN FAC=1D0/MAX(1,NEVFM) IF(MSTU(42).LE.0) WRITE(MSTU(11),5400) NEVFM,'eta' IF(MSTU(42).EQ.1) WRITE(MSTU(11),5400) NEVFM,'ypi' IF(MSTU(42).GE.2) WRITE(MSTU(11),5400) NEVFM,'y ' DO 510 IM=1,3 WRITE(MSTU(11),5500) DO 500 IB=1,10 BYETA=2D0*PARU(57) IF(IM.NE.2) BYETA=BYETA/2**(IB-1) BPHI=PARU(2) IF(IM.NE.1) BPHI=BPHI/2**(IB-1) IF(IM.LE.2) BNAVE=FAC*NMUFM/DBLE(2**(IB-1)) IF(IM.EQ.3) BNAVE=FAC*NMUFM/DBLE(4**(IB-1)) DO 490 IP=1,4 FMOMA(IP)=FAC*FM1FM(IM,IB,IP) FMOMS(IP)=SQRT(MAX(0D0,FAC*(FAC*FM2FM(IM,IB,IP)- & FMOMA(IP)**2))) 490 CONTINUE WRITE(MSTU(11),5600) BYETA,BPHI,BNAVE,(FMOMA(IP),FMOMS(IP), & IP=1,4) 500 CONTINUE 510 CONTINUE C...Copy statistics on factorial moments into /PYJETS/. ELSEIF(MTABU.EQ.33) THEN FAC=1D0/MAX(1,NEVFM) DO 540 IM=1,3 DO 530 IB=1,10 I=10*(IM-1)+IB K(I,1)=32 K(I,2)=99 K(I,3)=1 IF(IM.NE.2) K(I,3)=2**(IB-1) K(I,4)=1 IF(IM.NE.1) K(I,4)=2**(IB-1) K(I,5)=0 P(I,1)=2D0*PARU(57)/K(I,3) V(I,1)=PARU(2)/K(I,4) DO 520 IP=1,4 P(I,IP+1)=FAC*FM1FM(IM,IB,IP) V(I,IP+1)=SQRT(MAX(0D0,FAC*(FAC*FM2FM(IM,IB,IP)- & P(I,IP+1)**2))) 520 CONTINUE 530 CONTINUE 540 CONTINUE N=30 DO 550 J=1,5 K(N+1,J)=0 P(N+1,J)=0D0 V(N+1,J)=0D0 550 CONTINUE K(N+1,1)=32 K(N+1,2)=99 K(N+1,5)=NEVFM MSTU(3)=1 C...Reset statistics on Energy-Energy Correlation. ELSEIF(MTABU.EQ.40) THEN NEVEE=0 DO 560 J=1,25 FE1EC(J)=0D0 FE2EC(J)=0D0 FE1EC(51-J)=0D0 FE2EC(51-J)=0D0 FE1EA(J)=0D0 FE2EA(J)=0D0 560 CONTINUE C...Find particles to include, with proper assumed mass. ELSEIF(MTABU.EQ.41) THEN NEVEE=NEVEE+1 NLOW=N+MSTU(3) NUPP=NLOW ECM=0D0 DO 570 I=1,N IF(K(I,1).LE.0.OR.K(I,1).GT.10) GOTO 570 IF(MSTU(41).GE.2) THEN KC=PYCOMP(K(I,2)) IF(KC.EQ.0.OR.KC.EQ.12.OR.KC.EQ.14.OR.KC.EQ.16.OR. & KC.EQ.18) GOTO 570 IF(MSTU(41).GE.3.AND.KCHG(KC,2).EQ.0.AND. & PYCHGE(K(I,2)).EQ.0) GOTO 570 ENDIF PMR=0D0 IF(MSTU(42).EQ.1.AND.K(I,2).NE.22) PMR=PYMASS(211) IF(MSTU(42).GE.2) PMR=P(I,5) IF(NUPP.GT.MSTU(4)-5-MSTU(32)) THEN CALL PYERRM(11,'(PYTABU:) no more memory left in PYJETS') RETURN ENDIF NUPP=NUPP+1 P(NUPP,1)=P(I,1) P(NUPP,2)=P(I,2) P(NUPP,3)=P(I,3) P(NUPP,4)=SQRT(PMR**2+P(I,1)**2+P(I,2)**2+P(I,3)**2) P(NUPP,5)=MAX(1D-10,SQRT(P(I,1)**2+P(I,2)**2+P(I,3)**2)) ECM=ECM+P(NUPP,4) 570 CONTINUE IF(NUPP.EQ.NLOW) RETURN C...Analyze Energy-Energy Correlation in event. FAC=(2D0/ECM**2)*50D0/PARU(1) DO 580 J=1,50 FEVEE(J)=0D0 580 CONTINUE DO 600 I1=NLOW+2,NUPP DO 590 I2=NLOW+1,I1-1 CTHE=(P(I1,1)*P(I2,1)+P(I1,2)*P(I2,2)+P(I1,3)*P(I2,3))/ & (P(I1,5)*P(I2,5)) THE=ACOS(MAX(-1D0,MIN(1D0,CTHE))) ITHE=MAX(1,MIN(50,1+INT(50D0*THE/PARU(1)))) FEVEE(ITHE)=FEVEE(ITHE)+FAC*P(I1,4)*P(I2,4) 590 CONTINUE 600 CONTINUE DO 610 J=1,25 FE1EC(J)=FE1EC(J)+FEVEE(J) FE2EC(J)=FE2EC(J)+FEVEE(J)**2 FE1EC(51-J)=FE1EC(51-J)+FEVEE(51-J) FE2EC(51-J)=FE2EC(51-J)+FEVEE(51-J)**2 FE1EA(J)=FE1EA(J)+(FEVEE(51-J)-FEVEE(J)) FE2EA(J)=FE2EA(J)+(FEVEE(51-J)-FEVEE(J))**2 610 CONTINUE MSTU(62)=NUPP-NLOW C...Write statistics on Energy-Energy Correlation. ELSEIF(MTABU.EQ.42) THEN FAC=1D0/MAX(1,NEVEE) WRITE(MSTU(11),5700) NEVEE DO 620 J=1,25 FEEC1=FAC*FE1EC(J) FEES1=SQRT(MAX(0D0,FAC*(FAC*FE2EC(J)-FEEC1**2))) FEEC2=FAC*FE1EC(51-J) FEES2=SQRT(MAX(0D0,FAC*(FAC*FE2EC(51-J)-FEEC2**2))) FEECA=FAC*FE1EA(J) FEESA=SQRT(MAX(0D0,FAC*(FAC*FE2EA(J)-FEECA**2))) WRITE(MSTU(11),5800) 3.6D0*(J-1),3.6D0*J,FEEC1,FEES1, & FEEC2,FEES2,FEECA,FEESA 620 CONTINUE C...Copy statistics on Energy-Energy Correlation into /PYJETS/. ELSEIF(MTABU.EQ.43) THEN FAC=1D0/MAX(1,NEVEE) DO 630 I=1,25 K(I,1)=32 K(I,2)=99 K(I,3)=0 K(I,4)=0 K(I,5)=0 P(I,1)=FAC*FE1EC(I) V(I,1)=SQRT(MAX(0D0,FAC*(FAC*FE2EC(I)-P(I,1)**2))) P(I,2)=FAC*FE1EC(51-I) V(I,2)=SQRT(MAX(0D0,FAC*(FAC*FE2EC(51-I)-P(I,2)**2))) P(I,3)=FAC*FE1EA(I) V(I,3)=SQRT(MAX(0D0,FAC*(FAC*FE2EA(I)-P(I,3)**2))) P(I,4)=PARU(1)*(I-1)/50D0 P(I,5)=PARU(1)*I/50D0 V(I,4)=3.6D0*(I-1) V(I,5)=3.6D0*I 630 CONTINUE N=25 DO 640 J=1,5 K(N+1,J)=0 P(N+1,J)=0D0 V(N+1,J)=0D0 640 CONTINUE K(N+1,1)=32 K(N+1,2)=99 K(N+1,5)=NEVEE MSTU(3)=1 C...Reset statistics on decay channels. ELSEIF(MTABU.EQ.50) THEN NEVDC=0 NKFDC=0 NREDC=0 C...Identify and order flavour content of final state. ELSEIF(MTABU.EQ.51) THEN NEVDC=NEVDC+1 NDS=0 DO 670 I=1,N IF(K(I,1).LE.0.OR.K(I,1).GE.6) GOTO 670 NDS=NDS+1 IF(NDS.GT.8) THEN NREDC=NREDC+1 RETURN ENDIF KFM=2*IABS(K(I,2)) IF(K(I,2).LT.0) KFM=KFM-1 DO 650 IDS=NDS-1,1,-1 IIN=IDS+1 IF(KFM.LT.KFDM(IDS)) GOTO 660 KFDM(IDS+1)=KFDM(IDS) 650 CONTINUE IIN=1 660 KFDM(IIN)=KFM 670 CONTINUE C...Find whether old or new final state. DO 690 IDC=1,NKFDC IF(NDS.LT.KFDC(IDC,0)) THEN IKFDC=IDC GOTO 700 ELSEIF(NDS.EQ.KFDC(IDC,0)) THEN DO 680 I=1,NDS IF(KFDM(I).LT.KFDC(IDC,I)) THEN IKFDC=IDC GOTO 700 ELSEIF(KFDM(I).GT.KFDC(IDC,I)) THEN GOTO 690 ENDIF 680 CONTINUE IKFDC=-IDC GOTO 700 ENDIF 690 CONTINUE IKFDC=NKFDC+1 700 IF(IKFDC.LT.0) THEN IKFDC=-IKFDC ELSEIF(NKFDC.GE.200) THEN NREDC=NREDC+1 RETURN ELSE DO 720 IDC=NKFDC,IKFDC,-1 NPDC(IDC+1)=NPDC(IDC) DO 710 I=0,8 KFDC(IDC+1,I)=KFDC(IDC,I) 710 CONTINUE 720 CONTINUE NKFDC=NKFDC+1 KFDC(IKFDC,0)=NDS DO 730 I=1,NDS KFDC(IKFDC,I)=KFDM(I) 730 CONTINUE NPDC(IKFDC)=0 ENDIF NPDC(IKFDC)=NPDC(IKFDC)+1 C...Write statistics on decay channels. ELSEIF(MTABU.EQ.52) THEN FAC=1D0/MAX(1,NEVDC) WRITE(MSTU(11),5900) NEVDC DO 750 IDC=1,NKFDC DO 740 I=1,KFDC(IDC,0) KFM=KFDC(IDC,I) KF=(KFM+1)/2 IF(2*KF.NE.KFM) KF=-KF CALL PYNAME(KF,CHAU) CHDC(I)=CHAU(1:12) IF(CHAU(13:13).NE.' ') CHDC(I)(12:12)='?' 740 CONTINUE WRITE(MSTU(11),6000) FAC*NPDC(IDC),(CHDC(I),I=1,KFDC(IDC,0)) 750 CONTINUE IF(NREDC.NE.0) WRITE(MSTU(11),6100) FAC*NREDC C...Copy statistics on decay channels into /PYJETS/. ELSEIF(MTABU.EQ.53) THEN FAC=1D0/MAX(1,NEVDC) DO 780 IDC=1,NKFDC K(IDC,1)=32 K(IDC,2)=99 K(IDC,3)=0 K(IDC,4)=0 K(IDC,5)=KFDC(IDC,0) DO 760 J=1,5 P(IDC,J)=0D0 V(IDC,J)=0D0 760 CONTINUE DO 770 I=1,KFDC(IDC,0) KFM=KFDC(IDC,I) KF=(KFM+1)/2 IF(2*KF.NE.KFM) KF=-KF IF(I.LE.5) P(IDC,I)=KF IF(I.GE.6) V(IDC,I-5)=KF 770 CONTINUE V(IDC,5)=FAC*NPDC(IDC) 780 CONTINUE N=NKFDC DO 790 J=1,5 K(N+1,J)=0 P(N+1,J)=0D0 V(N+1,J)=0D0 790 CONTINUE K(N+1,1)=32 K(N+1,2)=99 K(N+1,5)=NEVDC V(N+1,5)=FAC*NREDC MSTU(3)=1 ENDIF C...Format statements for output on unit MSTU(11) (default 6). 5000 FORMAT(///20X,'Event statistics - initial state'/ &20X,'based on an analysis of ',I6,' events'// &3X,'Main flavours after',8X,'Fraction',4X,'Subfractions ', &'according to fragmenting system multiplicity'/ &4X,'hard interaction',24X,'1',7X,'2',7X,'3',7X,'4',7X,'5', &6X,'6-7',5X,'8-10',3X,'11-15',3X,'16-25',4X,'>25'/) 5100 FORMAT(3X,A12,1X,A12,F10.5,1X,10F8.4) 5200 FORMAT(///20X,'Event statistics - final state'/ &20X,'based on an analysis of ',I7,' events'// &5X,'Mean primary multiplicity =',F10.4/ &5X,'Mean final multiplicity =',F10.4/ &5X,'Mean charged multiplicity =',F10.4// &5X,'Number of particles produced per event (directly and via ', &'decays/branchings)'/ &8X,'KF Particle/jet MDCY',10X,'Particles',13X,'Antiparticles', &8X,'Total'/35X,'prim seco prim seco'/) 5300 FORMAT(1X,I9,4X,A16,I2,5(1X,F11.6)) 5400 FORMAT(///20X,'Factorial moments analysis of multiplicity'/ &20X,'based on an analysis of ',I6,' events'// &3X,'delta-',A3,' delta-phi /bin',10X,'',18X,'', &18X,'',18X,''/35X,4(' value error ')) 5500 FORMAT(10X) 5600 FORMAT(2X,2F10.4,F12.4,4(F12.4,F10.4)) 5700 FORMAT(///20X,'Energy-Energy Correlation and Asymmetry'/ &20X,'based on an analysis of ',I6,' events'// &2X,'theta range',8X,'EEC(theta)',8X,'EEC(180-theta)',7X, &'EECA(theta)'/2X,'in degrees ',3(' value error')/) 5800 FORMAT(2X,F4.1,' - ',F4.1,3(F11.4,F9.4)) 5900 FORMAT(///20X,'Decay channel analysis - final state'/ &20X,'based on an analysis of ',I6,' events'// &2X,'Probability',10X,'Complete final state'/) 6000 FORMAT(2X,F9.5,5X,8(A12,1X)) 6100 FORMAT(2X,F9.5,5X,'into other channels (more than 8 particles ', &'or table overflow)') RETURN END