pythiaeRHIC/pythia/pyoper.f

0001
0002 C*********************************************************************
0003
0004 C...PYOPER
0005 C...Performs operations between histograms.
0006
0007       SUBROUTINE PYOPER(ID1,OPER,ID2,ID3,F1,F2)
0008
0009 C...Double precision declaration.
0010       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
0011       IMPLICIT INTEGER(I-N)
0012 C...Commonblock.
0013       COMMON/PYBINS/IHIST(4),INDX(1000),BIN(20000)
0014       SAVE /PYBINS/
0015 C...Character variable.
0016       CHARACTER OPER*(*)
0017
0018 C...Find initial addresses in memory, and histogram size.
0019       IF(ID1.LE.0.OR.ID1.GT.IHIST(1)) CALL PYERRM(28,
0020      &'(PYFACT:) not allowed histogram number')
0021       IS1=INDX(ID1)
0022       IS2=INDX(MIN(IHIST(1),MAX(1,ID2)))
0023       IS3=INDX(MIN(IHIST(1),MAX(1,ID3)))
0024       NX=NINT(BIN(IS3+1))
0025       IF(OPER.EQ.'M'.AND.ID3.EQ.0) NX=NINT(BIN(IS2+1))
0026
0027 C...Update info on number of histogram entries.
0028       IF(OPER.EQ.'+'.OR.OPER.EQ.'-'.OR.OPER.EQ.'*'.OR.OPER.EQ.'/') THEN
0029         BIN(IS3+5)=BIN(IS1+5)+BIN(IS2+5)
0030       ELSEIF(OPER.EQ.'A'.OR.OPER.EQ.'S'.OR.OPER.EQ.'L') THEN
0031         BIN(IS3+5)=BIN(IS1+5)
0032       ENDIF
0033
0034 C...Operations on pair of histograms: addition, subtraction,
0035 C...multiplication, division.
0036       IF(OPER.EQ.'+') THEN
0037         DO 100 IX=6,8+NX
0038           BIN(IS3+IX)=F1*BIN(IS1+IX)+F2*BIN(IS2+IX)
0039   100   CONTINUE
0040       ELSEIF(OPER.EQ.'-') THEN
0041         DO 110 IX=6,8+NX
0042           BIN(IS3+IX)=F1*BIN(IS1+IX)-F2*BIN(IS2+IX)
0043   110   CONTINUE
0044       ELSEIF(OPER.EQ.'*') THEN
0045         DO 120 IX=6,8+NX
0046           BIN(IS3+IX)=F1*BIN(IS1+IX)*F2*BIN(IS2+IX)
0047   120   CONTINUE
0048       ELSEIF(OPER.EQ.'/') THEN
0049         DO 130 IX=6,8+NX
0050           FA2=F2*BIN(IS2+IX)
0051           IF(ABS(FA2).LE.1D-20) THEN
0052             BIN(IS3+IX)=0D0
0053           ELSE
0054             BIN(IS3+IX)=F1*BIN(IS1+IX)/FA2
0055           ENDIF
0056   130   CONTINUE
0057
0058 C...Operations on single histogram: multiplication+addition,
0059 C...square root+addition, logarithm+addition.
0060       ELSEIF(OPER.EQ.'A') THEN
0061         DO 140 IX=6,8+NX
0062           BIN(IS3+IX)=F1*BIN(IS1+IX)+F2
0063   140   CONTINUE
0064       ELSEIF(OPER.EQ.'S') THEN
0065         DO 150 IX=6,8+NX
0066           BIN(IS3+IX)=F1*SQRT(MAX(0D0,BIN(IS1+IX)))+F2
0067   150   CONTINUE
0068       ELSEIF(OPER.EQ.'L') THEN
0069         ZMIN=1D20
0070         DO 160 IX=9,8+NX
0071           IF(BIN(IS1+IX).LT.ZMIN.AND.BIN(IS1+IX).GT.1D-20)
0072      &    ZMIN=0.8D0*BIN(IS1+IX)
0073   160   CONTINUE
0074         DO 170 IX=6,8+NX
0075           BIN(IS3+IX)=F1*LOG10(MAX(ZMIN,BIN(IS1+IX)))+F2
0076   170   CONTINUE
0077
0078 C...Operation on two or three histograms: average and
0079 C...standard deviation.
0080       ELSEIF(OPER.EQ.'M') THEN
0081         DO 180 IX=6,8+NX
0082           IF(ABS(BIN(IS1+IX)).LE.1D-20) THEN
0083             BIN(IS2+IX)=0D0
0084           ELSE
0085             BIN(IS2+IX)=BIN(IS2+IX)/BIN(IS1+IX)
0086           ENDIF
0087           IF(ID3.NE.0) THEN
0088             IF(ABS(BIN(IS1+IX)).LE.1D-20) THEN
0089               BIN(IS3+IX)=0D0
0090             ELSE
0091               BIN(IS3+IX)=SQRT(MAX(0D0,BIN(IS3+IX)/BIN(IS1+IX)-
0092      &        BIN(IS2+IX)**2))
0093             ENDIF
0094           ENDIF
0095           BIN(IS1+IX)=F1*BIN(IS1+IX)
0096   180   CONTINUE
0097       ENDIF
0098
0099       RETURN
0100       END