hijing/src/pyhiklim.f

0001
0002 C***********************************************************************
0003
0004       SUBROUTINE PYHIKLIM(ILIM)
0005
0006 C...Checks generated variables against pre-set kinematical limits;
0007 C...also calculates limits on variables used in generation.
0008       COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0009       SAVE /LUDAT1/
0010       COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4)
0011       SAVE /LUDAT2/
0012       COMMON/LUDAT3/MDCY(500,3),MDME(2000,2),BRAT(2000),KFDP(2000,5)
0013       SAVE /LUDAT3/
0014       COMMON/PYHIPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
0015       SAVE /PYHIPARS/
0016       COMMON/PYHISUBS/MSEL,MSUB(200),KFIN(2,-40:40),CKIN(200)
0017       SAVE /PYHISUBS/
0018       COMMON/PYHIINT1/MINT(400),VINT(400)
0019       SAVE /PYHIINT1/
0020       COMMON/PYHIINT2/ISET(200),KFPR(200,2),COEF(200,20),ICOL(40,4,2)
0021       SAVE /PYHIINT2/
0022
0023 C...Common kinematical expressions.
0024       ISUB=MINT(1)
0025       IF(ISUB.EQ.96) GOTO 110
0026       SQM3=VINT(63)
0027       SQM4=VINT(64)
0028       IF(ILIM.NE.1) THEN
0029         TAU=VINT(21)
0030         RM3=SQM3/(TAU*VINT(2))
0031         RM4=SQM4/(TAU*VINT(2))
0032         BE34=SQRT((1.-RM3-RM4)**2-4.*RM3*RM4)
0033       ENDIF
0034       PTHMIN=CKIN(3)
0035       IF(MIN(SQM3,SQM4).LT.CKIN(6)**2) PTHMIN=MAX(CKIN(3),CKIN(5))
0036
0037       IF(ILIM.EQ.0) THEN
0038 C...Check generated values of tau, y*, cos(theta-hat), and tau' against
0039 C...pre-set kinematical limits.
0040         YST=VINT(22)
0041         CTH=VINT(23)
0042         TAUP=VINT(26)
0043         IF(ISET(ISUB).LE.2) THEN
0044           X1=SQRT(TAU)*EXP(YST)
0045           X2=SQRT(TAU)*EXP(-YST)
0046         ELSE
0047           X1=SQRT(TAUP)*EXP(YST)
0048           X2=SQRT(TAUP)*EXP(-YST)
0049         ENDIF
0050         XF=X1-X2
0051         IF(TAU*VINT(2).LT.CKIN(1)**2) MINT(51)=1
0052         IF(CKIN(2).GE.0..AND.TAU*VINT(2).GT.CKIN(2)**2) MINT(51)=1
0053         IF(X1.LT.CKIN(21).OR.X1.GT.CKIN(22)) MINT(51)=1
0054         IF(X2.LT.CKIN(23).OR.X2.GT.CKIN(24)) MINT(51)=1
0055         IF(XF.LT.CKIN(25).OR.XF.GT.CKIN(26)) MINT(51)=1
0056         IF(YST.LT.CKIN(7).OR.YST.GT.CKIN(8)) MINT(51)=1
0057         IF(ISET(ISUB).EQ.2.OR.ISET(ISUB).EQ.4) THEN
0058           PTH=0.5*BE34*SQRT(TAU*VINT(2)*(1.-CTH**2))
0059           Y3=YST+0.5*LOG((1.+RM3-RM4+BE34*CTH)/(1.+RM3-RM4-BE34*CTH))
0060           Y4=YST+0.5*LOG((1.+RM4-RM3-BE34*CTH)/(1.+RM4-RM3+BE34*CTH))
0061           YLARGE=MAX(Y3,Y4)
0062           YSMALL=MIN(Y3,Y4)
0063           ETALAR=10.
0064           ETASMA=-10.
0065           STH=SQRT(1.-CTH**2)
0066           IF(STH.LT.1.E-6) GOTO 100
0067           EXPET3=((1.+RM3-RM4)*SINH(YST)+BE34*COSH(YST)*CTH+
0068      &    SQRT(((1.+RM3-RM4)*COSH(YST)+BE34*SINH(YST)*CTH)**2-4.*RM3))/
0069      &    (BE34*STH)
0070           EXPET4=((1.-RM3+RM4)*SINH(YST)-BE34*COSH(YST)*CTH+
0071      &    SQRT(((1.-RM3+RM4)*COSH(YST)-BE34*SINH(YST)*CTH)**2-4.*RM4))/
0072      &    (BE34*STH)
0073           ETA3=LOG(MIN(1.E10,MAX(1.E-10,EXPET3)))
0074           ETA4=LOG(MIN(1.E10,MAX(1.E-10,EXPET4)))
0075           ETALAR=MAX(ETA3,ETA4)
0076           ETASMA=MIN(ETA3,ETA4)
0077   100     CTS3=((1.+RM3-RM4)*SINH(YST)+BE34*COSH(YST)*CTH)/
0078      &    SQRT(((1.+RM3-RM4)*COSH(YST)+BE34*SINH(YST)*CTH)**2-4.*RM3)
0079           CTS4=((1.-RM3+RM4)*SINH(YST)-BE34*COSH(YST)*CTH)/
0080      &    SQRT(((1.-RM3+RM4)*COSH(YST)-BE34*SINH(YST)*CTH)**2-4.*RM4)
0081           CTSLAR=MAX(CTS3,CTS4)
0082           CTSSMA=MIN(CTS3,CTS4)
0083           IF(PTH.LT.PTHMIN) MINT(51)=1
0084           IF(CKIN(4).GE.0..AND.PTH.GT.CKIN(4)) MINT(51)=1
0085           IF(YLARGE.LT.CKIN(9).OR.YLARGE.GT.CKIN(10)) MINT(51)=1
0086           IF(YSMALL.LT.CKIN(11).OR.YSMALL.GT.CKIN(12)) MINT(51)=1
0087           IF(ETALAR.LT.CKIN(13).OR.ETALAR.GT.CKIN(14)) MINT(51)=1
0088           IF(ETASMA.LT.CKIN(15).OR.ETASMA.GT.CKIN(16)) MINT(51)=1
0089           IF(CTSLAR.LT.CKIN(17).OR.CTSLAR.GT.CKIN(18)) MINT(51)=1
0090           IF(CTSSMA.LT.CKIN(19).OR.CTSSMA.GT.CKIN(20)) MINT(51)=1
0091           IF(CTH.LT.CKIN(27).OR.CTH.GT.CKIN(28)) MINT(51)=1
0092         ENDIF
0093         IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) THEN
0094           IF(TAUP*VINT(2).LT.CKIN(31)**2) MINT(51)=1
0095           IF(CKIN(32).GE.0..AND.TAUP*VINT(2).GT.CKIN(32)**2) MINT(51)=1
0096         ENDIF
0097
0098       ELSEIF(ILIM.EQ.1) THEN
0099 C...Calculate limits on tau
0100 C...0) due to definition
0101         TAUMN0=0.
0102         TAUMX0=1.
0103 C...1) due to limits on subsystem mass
0104         TAUMN1=CKIN(1)**2/VINT(2)
0105         TAUMX1=1.
0106         IF(CKIN(2).GE.0.) TAUMX1=CKIN(2)**2/VINT(2)
0107 C...2) due to limits on pT-hat (and non-overlapping rapidity intervals)
0108         TM3=SQRT(SQM3+PTHMIN**2)
0109         TM4=SQRT(SQM4+PTHMIN**2)
0110         YDCOSH=1.
0111         IF(CKIN(9).GT.CKIN(12)) YDCOSH=COSH(CKIN(9)-CKIN(12))
0112         TAUMN2=(TM3**2+2.*TM3*TM4*YDCOSH+TM4**2)/VINT(2)
0113         TAUMX2=1.
0114 C...3) due to limits on pT-hat and cos(theta-hat)
0115         CTH2MN=MIN(CKIN(27)**2,CKIN(28)**2)
0116         CTH2MX=MAX(CKIN(27)**2,CKIN(28)**2)
0117         TAUMN3=0.
0118         IF(CKIN(27)*CKIN(28).GT.0.) TAUMN3=
0119      &  (SQRT(SQM3+PTHMIN**2/(1.-CTH2MN))+
0120      &  SQRT(SQM4+PTHMIN**2/(1.-CTH2MN)))**2/VINT(2)
0121         TAUMX3=1.
0122         IF(CKIN(4).GE.0..AND.CTH2MX.LT.1.) TAUMX3=
0123      &  (SQRT(SQM3+CKIN(4)**2/(1.-CTH2MX))+
0124      &  SQRT(SQM4+CKIN(4)**2/(1.-CTH2MX)))**2/VINT(2)
0125 C...4) due to limits on x1 and x2
0126         TAUMN4=CKIN(21)*CKIN(23)
0127         TAUMX4=CKIN(22)*CKIN(24)
0128 C...5) due to limits on xF
0129         TAUMN5=0.
0130         TAUMX5=MAX(1.-CKIN(25),1.+CKIN(26))
0131         VINT(11)=MAX(TAUMN0,TAUMN1,TAUMN2,TAUMN3,TAUMN4,TAUMN5)
0132         VINT(31)=MIN(TAUMX0,TAUMX1,TAUMX2,TAUMX3,TAUMX4,TAUMX5)
0133         IF(MINT(43).EQ.1.AND.(ISET(ISUB).EQ.1.OR.ISET(ISUB).EQ.2)) THEN
0134           VINT(11)=0.99999
0135           VINT(31)=1.00001
0136         ENDIF
0137         IF(VINT(31).LE.VINT(11)) MINT(51)=1
0138
0139       ELSEIF(ILIM.EQ.2) THEN
0140 C...Calculate limits on y*
0141         IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) TAU=VINT(26)
0142         TAURT=SQRT(TAU)
0143 C...0) due to kinematics
0144         YSTMN0=LOG(TAURT)
0145         YSTMX0=-YSTMN0
0146 C...1) due to explicit limits
0147         YSTMN1=CKIN(7)
0148         YSTMX1=CKIN(8)
0149 C...2) due to limits on x1
0150         YSTMN2=LOG(MAX(TAU,CKIN(21))/TAURT)
0151         YSTMX2=LOG(MAX(TAU,CKIN(22))/TAURT)
0152 C...3) due to limits on x2
0153         YSTMN3=-LOG(MAX(TAU,CKIN(24))/TAURT)
0154         YSTMX3=-LOG(MAX(TAU,CKIN(23))/TAURT)
0155 C...4) due to limits on xF
0156         YEPMN4=0.5*ABS(CKIN(25))/TAURT
0157         YSTMN4=SIGN(LOG(SQRT(1.+YEPMN4**2)+YEPMN4),CKIN(25))
0158         YEPMX4=0.5*ABS(CKIN(26))/TAURT
0159         YSTMX4=SIGN(LOG(SQRT(1.+YEPMX4**2)+YEPMX4),CKIN(26))
0160 C...5) due to simultaneous limits on y-large and y-small
0161         YEPSMN=(RM3-RM4)*SINH(CKIN(9)-CKIN(11))
0162         YEPSMX=(RM3-RM4)*SINH(CKIN(10)-CKIN(12))
0163         YDIFMN=ABS(LOG(SQRT(1.+YEPSMN**2)-YEPSMN))
0164         YDIFMX=ABS(LOG(SQRT(1.+YEPSMX**2)-YEPSMX))
0165         YSTMN5=0.5*(CKIN(9)+CKIN(11)-YDIFMN)
0166         YSTMX5=0.5*(CKIN(10)+CKIN(12)+YDIFMX)
0167 C...6) due to simultaneous limits on cos(theta-hat) and y-large or
0168 C...   y-small
0169         CTHLIM=SQRT(1.-4.*PTHMIN**2/(BE34*TAU*VINT(2)))
0170         RZMN=BE34*MAX(CKIN(27),-CTHLIM)
0171         RZMX=BE34*MIN(CKIN(28),CTHLIM)
0172         YEX3MX=(1.+RM3-RM4+RZMX)/MAX(1E-10,1.+RM3-RM4-RZMX)
0173         YEX4MX=(1.+RM4-RM3-RZMN)/MAX(1E-10,1.+RM4-RM3+RZMN)
0174         YEX3MN=MAX(1E-10,1.+RM3-RM4+RZMN)/(1.+RM3-RM4-RZMN)
0175         YEX4MN=MAX(1E-10,1.+RM4-RM3-RZMX)/(1.+RM4-RM3+RZMX)
0176         YSTMN6=CKIN(9)-0.5*LOG(MAX(YEX3MX,YEX4MX))
0177         YSTMX6=CKIN(12)-0.5*LOG(MIN(YEX3MN,YEX4MN))
0178         VINT(12)=MAX(YSTMN0,YSTMN1,YSTMN2,YSTMN3,YSTMN4,YSTMN5,YSTMN6)
0179         VINT(32)=MIN(YSTMX0,YSTMX1,YSTMX2,YSTMX3,YSTMX4,YSTMX5,YSTMX6)
0180         IF(MINT(43).EQ.1) THEN
0181           VINT(12)=-0.00001
0182           VINT(32)=0.00001
0183         ELSEIF(MINT(43).EQ.2) THEN
0184           VINT(12)=0.99999*YSTMX0
0185           VINT(32)=1.00001*YSTMX0
0186         ELSEIF(MINT(43).EQ.3) THEN
0187           VINT(12)=-1.00001*YSTMX0
0188           VINT(32)=-0.99999*YSTMX0
0189         ENDIF
0190         IF(VINT(32).LE.VINT(12)) MINT(51)=1
0191
0192       ELSEIF(ILIM.EQ.3) THEN
0193 C...Calculate limits on cos(theta-hat)
0194         YST=VINT(22)
0195 C...0) due to definition
0196         CTNMN0=-1.
0197         CTNMX0=0.
0198         CTPMN0=0.
0199         CTPMX0=1.
0200 C...1) due to explicit limits
0201         CTNMN1=MIN(0.,CKIN(27))
0202         CTNMX1=MIN(0.,CKIN(28))
0203         CTPMN1=MAX(0.,CKIN(27))
0204         CTPMX1=MAX(0.,CKIN(28))
0205 C...2) due to limits on pT-hat
0206         CTNMN2=-SQRT(1.-4.*PTHMIN**2/(BE34**2*TAU*VINT(2)))
0207         CTPMX2=-CTNMN2
0208         CTNMX2=0.
0209         CTPMN2=0.
0210         IF(CKIN(4).GE.0.) THEN
0211           CTNMX2=-SQRT(MAX(0.,1.-4.*CKIN(4)**2/(BE34**2*TAU*VINT(2))))
0212           CTPMN2=-CTNMX2
0213         ENDIF
0214 C...3) due to limits on y-large and y-small
0215         CTNMN3=MIN(0.,MAX((1.+RM3-RM4)/BE34*TANH(CKIN(11)-YST),
0216      &  -(1.-RM3+RM4)/BE34*TANH(CKIN(10)-YST)))
0217         CTNMX3=MIN(0.,(1.+RM3-RM4)/BE34*TANH(CKIN(12)-YST),
0218      &  -(1.-RM3+RM4)/BE34*TANH(CKIN(9)-YST))
0219         CTPMN3=MAX(0.,(1.+RM3-RM4)/BE34*TANH(CKIN(9)-YST),
0220      &  -(1.-RM3+RM4)/BE34*TANH(CKIN(12)-YST))
0221         CTPMX3=MAX(0.,MIN((1.+RM3-RM4)/BE34*TANH(CKIN(10)-YST),
0222      &  -(1.-RM3+RM4)/BE34*TANH(CKIN(11)-YST)))
0223         VINT(13)=MAX(CTNMN0,CTNMN1,CTNMN2,CTNMN3)
0224         VINT(33)=MIN(CTNMX0,CTNMX1,CTNMX2,CTNMX3)
0225         VINT(14)=MAX(CTPMN0,CTPMN1,CTPMN2,CTPMN3)
0226         VINT(34)=MIN(CTPMX0,CTPMX1,CTPMX2,CTPMX3)
0227         IF(VINT(33).LE.VINT(13).AND.VINT(34).LE.VINT(14)) MINT(51)=1
0228
0229       ELSEIF(ILIM.EQ.4) THEN
0230 C...Calculate limits on tau'
0231 C...0) due to kinematics
0232         TAPMN0=TAU
0233         TAPMX0=1.
0234 C...1) due to explicit limits
0235         TAPMN1=CKIN(31)**2/VINT(2)
0236         TAPMX1=1.
0237         IF(CKIN(32).GE.0.) TAPMX1=CKIN(32)**2/VINT(2)
0238         VINT(16)=MAX(TAPMN0,TAPMN1)
0239         VINT(36)=MIN(TAPMX0,TAPMX1)
0240         IF(MINT(43).EQ.1) THEN
0241           VINT(16)=0.99999
0242           VINT(36)=1.00001
0243         ENDIF
0244         IF(VINT(36).LE.VINT(16)) MINT(51)=1
0245
0246       ENDIF
0247       RETURN
0248
0249 C...Special case for low-pT and multiple interactions:
0250 C...effective kinematical limits for tau, y*, cos(theta-hat).
0251   110 IF(ILIM.EQ.0) THEN
0252       ELSEIF(ILIM.EQ.1) THEN
0253         IF(MSTP(82).LE.1) VINT(11)=4.*PARP(81)**2/VINT(2)
0254         IF(MSTP(82).GE.2) VINT(11)=PARP(82)**2/VINT(2)
0255         VINT(31)=1.
0256       ELSEIF(ILIM.EQ.2) THEN
0257         VINT(12)=0.5*LOG(VINT(21))
0258         VINT(32)=-VINT(12)
0259       ELSEIF(ILIM.EQ.3) THEN
0260         IF(MSTP(82).LE.1) ST2EFF=4.*PARP(81)**2/(VINT(21)*VINT(2))
0261         IF(MSTP(82).GE.2) ST2EFF=0.01*PARP(82)**2/(VINT(21)*VINT(2))
0262         VINT(13)=-SQRT(MAX(0.,1.-ST2EFF))
0263         VINT(33)=0.
0264         VINT(14)=0.
0265         VINT(34)=-VINT(13)
0266       ENDIF
0267
0268       RETURN
0269       END