pythiaeRHIC/pythia/pyplot.f

0001
0002 C*********************************************************************
0003
0004 C...PYPLOT
0005 C...Prints a histogram (but does not reset it).
0006
0007       SUBROUTINE PYPLOT(ID)
0008
0009 C...Double precision declaration.
0010       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
0011       IMPLICIT INTEGER(I-N)
0012 C...Commonblocks.
0013       COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0014       COMMON/PYBINS/IHIST(4),INDX(1000),BIN(20000)
0015       SAVE /PYDAT1/,/PYBINS/
0016 C...Local arrays and character variables.
0017       DIMENSION IDATI(6), IROW(100), IFRA(100), DYAC(10)
0018       CHARACTER TITLE*60, OUT*100, CHA(0:11)*1
0019
0020 C...Steps in histogram scale. Character sequence.
0021       DATA DYAC/.04,.05,.06,.08,.10,.12,.15,.20,.25,.30/
0022       DATA CHA/'0','1','2','3','4','5','6','7','8','9','X','-'/
0023
0024 C...Find initial address in memory; skip if empty histogram.
0025       IF(ID.LE.0.OR.ID.GT.IHIST(1)) RETURN
0026       IS=INDX(ID)
0027       IF(IS.EQ.0) RETURN
0028       IF(NINT(BIN(IS+5)).LE.0) THEN
0029         WRITE(MSTU(11),5000) ID
0030         RETURN
0031       ENDIF
0032
0033 C...Number of histogram lines and x bins.
0034       LIN=IHIST(3)-18
0035       NX=NINT(BIN(IS+1))
0036
0037 C...Extract title by conversion from double precision via integer.
0038       DO 100 IT=1,20
0039         IEQ=NINT(BIN(IS+8+NX+IT))
0040         TITLE(3*IT-2:3*IT)=CHAR(IEQ/256**2)//CHAR(MOD(IEQ,256**2)/256)
0041      &  //CHAR(MOD(IEQ,256))
0042   100 CONTINUE
0043
0044 C...Find time; print title.
0045       CALL PYTIME(IDATI)
0046       IF(IDATI(1).GT.0) THEN
0047         WRITE(MSTU(11),5100) ID, TITLE, (IDATI(J),J=1,5)
0048       ELSE
0049         WRITE(MSTU(11),5200) ID, TITLE
0050       ENDIF
0051
0052 C...Find minimum and maximum bin content.
0053       YMIN=BIN(IS+9)
0054       YMAX=BIN(IS+9)
0055       DO 110 IX=IS+10,IS+8+NX
0056         IF(BIN(IX).LT.YMIN) YMIN=BIN(IX)
0057         IF(BIN(IX).GT.YMAX) YMAX=BIN(IX)
0058   110 CONTINUE
0059
0060 C...Determine scale and step size for y axis.
0061       IF(YMAX-YMIN.GT.LIN*DYAC(1)*1D-9) THEN
0062         IF(YMIN.GT.0D0.AND.YMIN.LT.0.1D0*YMAX) YMIN=0D0
0063         IF(YMAX.LT.0D0.AND.YMAX.GT.0.1D0*YMIN) YMAX=0D0
0064         IPOT=INT(LOG10(YMAX-YMIN)+10D0)-10
0065         IF(YMAX-YMIN.LT.LIN*DYAC(1)*10D0**IPOT) IPOT=IPOT-1
0066         IF(YMAX-YMIN.GT.LIN*DYAC(10)*10D0**IPOT) IPOT=IPOT+1
0067         DELY=DYAC(1)
0068         DO 120 IDEL=1,9
0069           IF(YMAX-YMIN.GE.LIN*DYAC(IDEL)*10D0**IPOT) DELY=DYAC(IDEL+1)
0070   120   CONTINUE
0071         DY=DELY*10D0**IPOT
0072
0073 C...Convert bin contents to integer form; fractional fill in top row.
0074         DO 130 IX=1,NX
0075           CTA=ABS(BIN(IS+8+IX))/DY
0076           IROW(IX)=SIGN(CTA+0.95D0,BIN(IS+8+IX))
0077           IFRA(IX)=10D0*(CTA+1.05D0-DBLE(INT(CTA+0.95D0)))
0078   130   CONTINUE
0079         IRMI=SIGN(ABS(YMIN)/DY+0.95D0,YMIN)
0080         IRMA=SIGN(ABS(YMAX)/DY+0.95D0,YMAX)
0081
0082 C...Print histogram row by row.
0083         DO 150 IR=IRMA,IRMI,-1
0084           IF(IR.EQ.0) GOTO 150
0085           OUT=' '
0086           DO 140 IX=1,NX
0087             IF(IR.EQ.IROW(IX)) OUT(IX:IX)=CHA(IFRA(IX))
0088             IF(IR*(IROW(IX)-IR).GT.0) OUT(IX:IX)=CHA(10)
0089   140     CONTINUE
0090           WRITE(MSTU(11),5300) IR*DELY, IPOT, OUT
0091   150   CONTINUE
0092
0093 C...Print sign and value of bin contents.
0094         IPOT=INT(LOG10(MAX(YMAX,-YMIN))+10.0001D0)-10
0095         OUT=' '
0096         DO 160 IX=1,NX
0097           IF(BIN(IS+8+IX).LT.-10D0**(IPOT-4)) OUT(IX:IX)=CHA(11)
0098           IROW(IX)=NINT(10D0**(3-IPOT)*ABS(BIN(IS+8+IX)))
0099   160   CONTINUE
0100         WRITE(MSTU(11),5400) OUT
0101         DO 180 IR=4,1,-1
0102           DO 170 IX=1,NX
0103             OUT(IX:IX)=CHA(MOD(IROW(IX),10**IR)/10**(IR-1))
0104   170     CONTINUE
0105           WRITE(MSTU(11),5500) IPOT+IR-4, OUT
0106   180   CONTINUE
0107
0108 C...Print sign and value of lower bin edge.
0109         IPOT=INT(LOG10(MAX(-BIN(IS+2),BIN(IS+3)-BIN(IS+4)))+
0110      &  10.0001D0)-10
0111         OUT=' '
0112         DO 190 IX=1,NX
0113           IF(BIN(IS+2)+(IX-1)*BIN(IS+4).LT.-10D0**(IPOT-3))
0114      &    OUT(IX:IX)=CHA(11)
0115           IROW(IX)=NINT(10D0**(2-IPOT)*ABS(BIN(IS+2)+(IX-1)*BIN(IS+4)))
0116   190   CONTINUE
0117         WRITE(MSTU(11),5600) OUT
0118         DO 210 IR=3,1,-1
0119           DO 200 IX=1,NX
0120             OUT(IX:IX)=CHA(MOD(IROW(IX),10**IR)/10**(IR-1))
0121   200     CONTINUE
0122           WRITE(MSTU(11),5500) IPOT+IR-3, OUT
0123   210   CONTINUE
0124       ENDIF
0125
0126 C...Calculate and print statistics.
0127       CSUM=0D0
0128       CXSUM=0D0
0129       CXXSUM=0D0
0130       DO 220 IX=1,NX
0131         CTA=ABS(BIN(IS+8+IX))
0132         X=BIN(IS+2)+(IX-0.5D0)*BIN(IS+4)
0133         CSUM=CSUM+CTA
0134         CXSUM=CXSUM+CTA*X
0135         CXXSUM=CXXSUM+CTA*X**2
0136   220 CONTINUE
0137       XMEAN=CXSUM/MAX(CSUM,1D-20)
0138       XRMS=SQRT(MAX(0D0,CXXSUM/MAX(CSUM,1D-20)-XMEAN**2))
0139       WRITE(MSTU(11),5700) NINT(BIN(IS+5)),XMEAN,BIN(IS+6),
0140      &BIN(IS+2),BIN(IS+7),XRMS,BIN(IS+8),BIN(IS+3)
0141
0142 C...Formats for output.
0143  5000 FORMAT(/5X,'Histogram no',I5,' : no entries')
0144  5100 FORMAT('1'/5X,'Histogram no',I5,6X,A60,5X,I4,'-',I2,'-',I2,1X,
0145      &I2,':',I2/)
0146  5200 FORMAT('1'/5X,'Histogram no',I5,6X,A60/)
0147  5300 FORMAT(2X,F7.2,'*10**',I2,3X,A100)
0148  5400 FORMAT(/8X,'Contents',3X,A100)
0149  5500 FORMAT(9X,'*10**',I2,3X,A100)
0150  5600 FORMAT(/8X,'Low edge',3X,A100)
0151  5700 FORMAT(/5X,'Entries  =',I12,1P,6X,'Mean =',D12.4,6X,'Underflow ='
0152      &,D12.4,6X,'Low edge  =',D12.4/5X,'All chan =',D12.4,6X,
0153      &'Rms  =',D12.4,6X,'Overflow  =',D12.4,6X,'High edge =',D12.4)
0154
0155       RETURN
0156       END