HIJING/hipyset1_35/pymult_hijing.F

   1 * $Id$
   2
   3 C*********************************************************************
   4
   5       SUBROUTINE PYMULT_HIJING(MMUL)
   6
   7 C...Initializes treatment of multiple interactions, selects kinematics
   8 C...of hardest interaction if low-pT physics included in run, and
   9 C...generates all non-hardest interactions.
  10 #include "lujets_hijing.inc"
  11 #include "ludat1_hijing.inc"
  12 #include "ludat2_hijing.inc"
  13 #include "pysubs_hijing.inc"
  14 #include "pypars_hijing.inc"
  15 #include "pyint1_hijing.inc"
  16 #include "pyint2_hijing.inc"
  17 #include "pyint3_hijing.inc"
  18 #include "pyint5_hijing.inc"
  19       DIMENSION NMUL(20),SIGM(20),KSTR(500,2)
  20       SAVE XT2,XT2FAC,XC2,XTS,IRBIN,RBIN,NMUL,SIGM
  21
  22 C...Initialization of multiple interaction treatment.
  23       IF(MMUL.EQ.1) THEN
  24         IF(MSTP(122).GE.1) WRITE(MSTU(11),1000) MSTP(82)
  25         ISUB=96
  26         MINT(1)=96
  27         VINT(63)=0.
  28         VINT(64)=0.
  29         VINT(143)=1.
  30         VINT(144)=1.
  31
  32 C...Loop over phase space points: xT2 choice in 20 bins.
  33   100   SIGSUM=0.
  34         DO 120 IXT2=1,20
  35         NMUL(IXT2)=MSTP(83)
  36         SIGM(IXT2)=0.
  37         DO 110 ITRY=1,MSTP(83)
  38         RSCA=0.05*((21-IXT2)-RLU_HIJING(0))
  39         XT2=VINT(149)*(1.+VINT(149))/(VINT(149)+RSCA)-VINT(149)
  40         XT2=MAX(0.01*VINT(149),XT2)
  41         VINT(25)=XT2
  42
  43 C...Choose tau and y*. Calculate cos(theta-hat).
  44         IF(RLU_HIJING(0).LE.COEF(ISUB,1)) THEN
  45           TAUP=(2.*(1.+SQRT(1.-XT2))/XT2-1.)**RLU_HIJING(0)
  46           TAU=XT2*(1.+TAUP)**2/(4.*TAUP)
  47         ELSE
  48           TAU=XT2*(1.+TAN(RLU_HIJING(0)*ATAN(SQRT(1./XT2-1.)))**2)
  49         ENDIF
  50         VINT(21)=TAU
  51         CALL PYKLIM_HIJING(2)
  52         RYST=RLU_HIJING(0)
  53         MYST=1
  54         IF(RYST.GT.COEF(ISUB,7)) MYST=2
  55         IF(RYST.GT.COEF(ISUB,7)+COEF(ISUB,8)) MYST=3
  56         CALL PYKMAP_HIJING(2,MYST,RLU_HIJING(0))
  57         VINT(23)=SQRT(MAX(0.,1.-XT2/TAU))*(-1)**INT(1.5+RLU_HIJING(0))
  58
  59 C...Calculate differential cross-section.
  60         VINT(71)=0.5*VINT(1)*SQRT(XT2)
  61         CALL PYSIGH_HIJING(NCHN,SIGS)
  62   110   SIGM(IXT2)=SIGM(IXT2)+SIGS
  63   120   SIGSUM=SIGSUM+SIGM(IXT2)
  64         SIGSUM=SIGSUM/(20.*MSTP(83))
  65
  66 C...Reject result if sigma(parton-parton) is smaller than hadronic one.
  67         IF(SIGSUM.LT.1.1*VINT(106)) THEN
  68           IF(MSTP(122).GE.1) WRITE(MSTU(11),1100) PARP(82),SIGSUM
  69           PARP(82)=0.9*PARP(82)
  70           VINT(149)=4.*PARP(82)**2/VINT(2)
  71           GOTO 100
  72         ENDIF
  73         IF(MSTP(122).GE.1) WRITE(MSTU(11),1200) PARP(82), SIGSUM
  74
  75 C...Start iteration to find k factor.
  76         YKE=SIGSUM/VINT(106)
  77         SO=0.5
  78         XI=0.
  79         YI=0.
  80         XK=0.5
  81         IIT=0
  82   130   IF(IIT.EQ.0) THEN
  83           XK=2.*XK
  84         ELSEIF(IIT.EQ.1) THEN
  85           XK=0.5*XK
  86         ELSE
  87           XK=XI+(YKE-YI)*(XF-XI)/(YF-YI)
  88         ENDIF
  89
  90 C...Evaluate overlap integrals.
  91         IF(MSTP(82).EQ.2) THEN
  92           SP=0.5*PARU(1)*(1.-EXP(-XK))
  93           SOP=SP/PARU(1)
  94         ELSE
  95           IF(MSTP(82).EQ.3) DELTAB=0.02
  96           IF(MSTP(82).EQ.4) DELTAB=MIN(0.01,0.05*PARP(84))
  97           SP=0.
  98           SOP=0.
  99           B=-0.5*DELTAB
 100   140     B=B+DELTAB
 101           IF(MSTP(82).EQ.3) THEN
 102             OV=EXP(-B**2)/PARU(2)
 103           ELSE
 104             CQ2=PARP(84)**2
 105             OV=((1.-PARP(83))**2*EXP(-MIN(100.,B**2))+2.*PARP(83)*
 106      &      (1.-PARP(83))*2./(1.+CQ2)*EXP(-MIN(100.,B**2*2./(1.+CQ2)))+
 107      &      PARP(83)**2/CQ2*EXP(-MIN(100.,B**2/CQ2)))/PARU(2)
 108           ENDIF
 109           PACC=1.-EXP(-MIN(100.,PARU(1)*XK*OV))
 110           SP=SP+PARU(2)*B*DELTAB*PACC
 111           SOP=SOP+PARU(2)*B*DELTAB*OV*PACC
 112           IF(B.LT.1..OR.B*PACC.GT.1E-6) GOTO 140
 113         ENDIF
 114         YK=PARU(1)*XK*SO/SP
 115
 116 C...Continue iteration until convergence.
 117         IF(YK.LT.YKE) THEN
 118           XI=XK
 119           YI=YK
 120           IF(IIT.EQ.1) IIT=2
 121         ELSE
 122           XF=XK
 123           YF=YK
 124           IF(IIT.EQ.0) IIT=1
 125         ENDIF
 126         IF(ABS(YK-YKE).GE.1E-5*YKE) GOTO 130
 127
 128 C...Store some results for subsequent use.
 129         VINT(145)=SIGSUM
 130         VINT(146)=SOP/SO
 131         VINT(147)=SOP/SP
 132
 133 C...Initialize iteration in xT2 for hardest interaction.
 134       ELSEIF(MMUL.EQ.2) THEN
 135         IF(MSTP(82).LE.0) THEN
 136         ELSEIF(MSTP(82).EQ.1) THEN
 137           XT2=1.
 138           XT2FAC=XSEC(96,1)/VINT(106)*VINT(149)/(1.-VINT(149))
 139         ELSEIF(MSTP(82).EQ.2) THEN
 140           XT2=1.
 141           XT2FAC=VINT(146)*XSEC(96,1)/VINT(106)*VINT(149)*(1.+VINT(149))
 142         ELSE
 143           XC2=4.*CKIN(3)**2/VINT(2)
 144           IF(CKIN(3).LE.CKIN(5).OR.MINT(82).GE.2) XC2=0.
 145         ENDIF
 146
 147       ELSEIF(MMUL.EQ.3) THEN
 148 C...Low-pT or multiple interactions (first semihard interaction):
 149 C...choose xT2 according to dpT2/pT2**2*exp(-(sigma above pT2)/norm)
 150 C...or (MSTP(82)>=2) dpT2/(pT2+pT0**2)**2*exp(-....).
 151         ISUB=MINT(1)
 152         IF(MSTP(82).LE.0) THEN
 153           XT2=0.
 154         ELSEIF(MSTP(82).EQ.1) THEN
 155           XT2=XT2FAC*XT2/(XT2FAC-XT2*LOG(RLU_HIJING(0)))
 156         ELSEIF(MSTP(82).EQ.2) THEN
 157           IF(XT2.LT.1..AND.EXP(-XT2FAC*XT2/(VINT(149)*(XT2+
 158      &    VINT(149)))).GT.RLU_HIJING(0)) XT2=1.
 159           IF(XT2.GE.1.) THEN
 160             XT2=(1.+VINT(149))*XT2FAC/(XT2FAC-(1.+VINT(149))*LOG(1.-
 161      &            RLU_HIJING(0)*(1.-EXP(-XT2FAC/(VINT(149)*(1.+VINT(149)
 162      $            ))))))-VINT(149)
 163           ELSE
 164             XT2=-XT2FAC/LOG(EXP(-XT2FAC/(XT2+VINT(149)))+RLU_HIJING(0)*
 165      &      (EXP(-XT2FAC/VINT(149))-EXP(-XT2FAC/(XT2+VINT(149)))))-
 166      &      VINT(149)
 167           ENDIF
 168           XT2=MAX(0.01*VINT(149),XT2)
 169         ELSE
 170           XT2=(XC2+VINT(149))*(1.+VINT(149))/(1.+VINT(149)-
 171      &    RLU_HIJING(0)*(1.-XC2))-VINT(149)
 172           XT2=MAX(0.01*VINT(149),XT2)
 173         ENDIF
 174         VINT(25)=XT2
 175
 176 C...Low-pT: choose xT2, tau, y* and cos(theta-hat) fixed.
 177         IF(MSTP(82).LE.1.AND.XT2.LT.VINT(149)) THEN
 178           IF(MINT(82).EQ.1) NGEN(0,1)=NGEN(0,1)-1
 179           IF(MINT(82).EQ.1) NGEN(ISUB,1)=NGEN(ISUB,1)-1
 180           ISUB=95
 181           MINT(1)=ISUB
 182           VINT(21)=0.01*VINT(149)
 183           VINT(22)=0.
 184           VINT(23)=0.
 185           VINT(25)=0.01*VINT(149)
 186
 187         ELSE
 188 C...Multiple interactions (first semihard interaction).
 189 C...Choose tau and y*. Calculate cos(theta-hat).
 190           IF(RLU_HIJING(0).LE.COEF(ISUB,1)) THEN
 191             TAUP=(2.*(1.+SQRT(1.-XT2))/XT2-1.)**RLU_HIJING(0)
 192             TAU=XT2*(1.+TAUP)**2/(4.*TAUP)
 193           ELSE
 194             TAU=XT2*(1.+TAN(RLU_HIJING(0)*ATAN(SQRT(1./XT2-1.)))**2)
 195           ENDIF
 196           VINT(21)=TAU
 197           CALL PYKLIM_HIJING(2)
 198           RYST=RLU_HIJING(0)
 199           MYST=1
 200           IF(RYST.GT.COEF(ISUB,7)) MYST=2
 201           IF(RYST.GT.COEF(ISUB,7)+COEF(ISUB,8)) MYST=3
 202           CALL PYKMAP_HIJING(2,MYST,RLU_HIJING(0))
 203           VINT(23)=SQRT(MAX(0.,1.-XT2/TAU))*(-1)**INT(1.5+RLU_HIJING(0))
 204         ENDIF
 205         VINT(71)=0.5*VINT(1)*SQRT(VINT(25))
 206
 207 C...Store results of cross-section calculation.
 208       ELSEIF(MMUL.EQ.4) THEN
 209         ISUB=MINT(1)
 210         XTS=VINT(25)
 211         IF(ISET(ISUB).EQ.1) XTS=VINT(21)
 212         IF(ISET(ISUB).EQ.2) XTS=(4.*VINT(48)+2.*VINT(63)+2.*VINT(64))/
 213      &  VINT(2)
 214         IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) XTS=VINT(26)
 215         RBIN=MAX(0.000001,MIN(0.999999,XTS*(1.+VINT(149))/
 216      &  (XTS+VINT(149))))
 217         IRBIN=INT(1.+20.*RBIN)
 218         IF(ISUB.EQ.96) NMUL(IRBIN)=NMUL(IRBIN)+1
 219         IF(ISUB.EQ.96) SIGM(IRBIN)=SIGM(IRBIN)+VINT(153)
 220
 221 C...Choose impact parameter.
 222       ELSEIF(MMUL.EQ.5) THEN
 223         IF(MSTP(82).EQ.3) THEN
 224           VINT(148)=RLU_HIJING(0)/(PARU(2)*VINT(147))
 225         ELSE
 226           RTYPE=RLU_HIJING(0)
 227           CQ2=PARP(84)**2
 228           IF(RTYPE.LT.(1.-PARP(83))**2) THEN
 229             B2=-LOG(RLU_HIJING(0))
 230           ELSEIF(RTYPE.LT.1.-PARP(83)**2) THEN
 231             B2=-0.5*(1.+CQ2)*LOG(RLU_HIJING(0))
 232           ELSE
 233             B2=-CQ2*LOG(RLU_HIJING(0))
 234           ENDIF
 235           VINT(148)=((1.-PARP(83))**2*EXP(-MIN(100.,B2))+2.*PARP(83)*
 236      &    (1.-PARP(83))*2./(1.+CQ2)*EXP(-MIN(100.,B2*2./(1.+CQ2)))+
 237      &    PARP(83)**2/CQ2*EXP(-MIN(100.,B2/CQ2)))/(PARU(2)*VINT(147))
 238         ENDIF
 239
 240 C...Multiple interactions (variable impact parameter) : reject with
 241 C...probability exp(-overlap*cross-section above pT/normalization).
 242         RNCOR=(IRBIN-20.*RBIN)*NMUL(IRBIN)
 243         SIGCOR=(IRBIN-20.*RBIN)*SIGM(IRBIN)
 244         DO 150 IBIN=IRBIN+1,20
 245         RNCOR=RNCOR+NMUL(IBIN)
 246   150   SIGCOR=SIGCOR+SIGM(IBIN)
 247         SIGABV=(SIGCOR/RNCOR)*VINT(149)*(1.-XTS)/(XTS+VINT(149))
 248         VINT(150)=EXP(-MIN(100.,VINT(146)*VINT(148)*SIGABV/VINT(106)))
 249
 250 C...Generate additional multiple semihard interactions.
 251       ELSEIF(MMUL.EQ.6) THEN
 252
 253 C...Reconstruct strings in hard scattering.
 254         ISUB=MINT(1)
 255         NMAX=MINT(84)+4
 256         IF(ISET(ISUB).EQ.1) NMAX=MINT(84)+2
 257         NSTR=0
 258         DO 170 I=MINT(84)+1,NMAX
 259         KCS=KCHG(LUCOMP_HIJING(K(I,2)),2)*ISIGN(1,K(I,2))
 260         IF(KCS.EQ.0) GOTO 170
 261         DO 160 J=1,4
 262         IF(KCS.EQ.1.AND.(J.EQ.2.OR.J.EQ.4)) GOTO 160
 263         IF(KCS.EQ.-1.AND.(J.EQ.1.OR.J.EQ.3)) GOTO 160
 264         IF(J.LE.2) THEN
 265           IST=MOD(K(I,J+3)/MSTU(5),MSTU(5))
 266         ELSE
 267           IST=MOD(K(I,J+1),MSTU(5))
 268         ENDIF
 269         IF(IST.LT.MINT(84).OR.IST.GT.I) GOTO 160
 270         IF(KCHG(LUCOMP_HIJING(K(IST,2)),2).EQ.0) GOTO 160
 271         NSTR=NSTR+1
 272         IF(J.EQ.1.OR.J.EQ.4) THEN
 273           KSTR(NSTR,1)=I
 274           KSTR(NSTR,2)=IST
 275         ELSE
 276           KSTR(NSTR,1)=IST
 277           KSTR(NSTR,2)=I
 278         ENDIF
 279   160   CONTINUE
 280   170   CONTINUE
 281
 282 C...Set up starting values for iteration in xT2.
 283         XT2=VINT(25)
 284         IF(ISET(ISUB).EQ.1) XT2=VINT(21)
 285         IF(ISET(ISUB).EQ.2) XT2=(4.*VINT(48)+2.*VINT(63)+2.*VINT(64))/
 286      &  VINT(2)
 287         IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) XT2=VINT(26)
 288         ISUB=96
 289         MINT(1)=96
 290         IF(MSTP(82).LE.1) THEN
 291           XT2FAC=XSEC(ISUB,1)*VINT(149)/((1.-VINT(149))*VINT(106))
 292         ELSE
 293           XT2FAC=VINT(146)*VINT(148)*XSEC(ISUB,1)/VINT(106)*
 294      &    VINT(149)*(1.+VINT(149))
 295         ENDIF
 296         VINT(63)=0.
 297         VINT(64)=0.
 298         VINT(151)=0.
 299         VINT(152)=0.
 300         VINT(143)=1.-VINT(141)
 301         VINT(144)=1.-VINT(142)
 302
 303 C...Iterate downwards in xT2.
 304   180   IF(MSTP(82).LE.1) THEN
 305           XT2=XT2FAC*XT2/(XT2FAC-XT2*LOG(RLU_HIJING(0)))
 306           IF(XT2.LT.VINT(149)) GOTO 220
 307         ELSE
 308           IF(XT2.LE.0.01*VINT(149)) GOTO 220
 309           XT2=XT2FAC*(XT2+VINT(149))/(XT2FAC-(XT2+VINT(149))*
 310      &    LOG(RLU_HIJING(0)))-VINT(149)
 311           IF(XT2.LE.0.) GOTO 220
 312           XT2=MAX(0.01*VINT(149),XT2)
 313         ENDIF
 314         VINT(25)=XT2
 315
 316 C...Choose tau and y*. Calculate cos(theta-hat).
 317         IF(RLU_HIJING(0).LE.COEF(ISUB,1)) THEN
 318           TAUP=(2.*(1.+SQRT(1.-XT2))/XT2-1.)**RLU_HIJING(0)
 319           TAU=XT2*(1.+TAUP)**2/(4.*TAUP)
 320         ELSE
 321           TAU=XT2*(1.+TAN(RLU_HIJING(0)*ATAN(SQRT(1./XT2-1.)))**2)
 322         ENDIF
 323         VINT(21)=TAU
 324         CALL PYKLIM_HIJING(2)
 325         RYST=RLU_HIJING(0)
 326         MYST=1
 327         IF(RYST.GT.COEF(ISUB,7)) MYST=2
 328         IF(RYST.GT.COEF(ISUB,7)+COEF(ISUB,8)) MYST=3
 329         CALL PYKMAP_HIJING(2,MYST,RLU_HIJING(0))
 330         VINT(23)=SQRT(MAX(0.,1.-XT2/TAU))*(-1)**INT(1.5+RLU_HIJING(0))
 331
 332 C...Check that x not used up. Accept or reject kinematical variables.
 333         X1M=SQRT(TAU)*EXP(VINT(22))
 334         X2M=SQRT(TAU)*EXP(-VINT(22))
 335         IF(VINT(143)-X1M.LT.0.01.OR.VINT(144)-X2M.LT.0.01) GOTO 180
 336         VINT(71)=0.5*VINT(1)*SQRT(XT2)
 337         CALL PYSIGH_HIJING(NCHN,SIGS)
 338         IF(SIGS.LT.XSEC(ISUB,1)*RLU_HIJING(0)) GOTO 180
 339
 340 C...Reset K, P and V vectors. Select some variables.
 341         DO 190 I=N+1,N+2
 342         DO 190 J=1,5
 343         K(I,J)=0
 344         P(I,J)=0.
 345   190   V(I,J)=0.
 346         RFLAV=RLU_HIJING(0)
 347         PT=0.5*VINT(1)*SQRT(XT2)
 348         PHI=PARU(2)*RLU_HIJING(0)
 349         CTH=VINT(23)
 350
 351 C...Add first parton to event record.
 352         K(N+1,1)=3
 353         K(N+1,2)=21
 354         IF(RFLAV.GE.MAX(PARP(85),PARP(86))) K(N+1,2)=
 355      &  1+INT((2.+PARJ(2))*RLU_HIJING(0))
 356         P(N+1,1)=PT*COS(PHI)
 357         P(N+1,2)=PT*SIN(PHI)
 358         P(N+1,3)=0.25*VINT(1)*(VINT(41)*(1.+CTH)-VINT(42)*(1.-CTH))
 359         P(N+1,4)=0.25*VINT(1)*(VINT(41)*(1.+CTH)+VINT(42)*(1.-CTH))
 360         P(N+1,5)=0.
 361
 362 C...Add second parton to event record.
 363         K(N+2,1)=3
 364         K(N+2,2)=21
 365         IF(K(N+1,2).NE.21) K(N+2,2)=-K(N+1,2)
 366         P(N+2,1)=-P(N+1,1)
 367         P(N+2,2)=-P(N+1,2)
 368         P(N+2,3)=0.25*VINT(1)*(VINT(41)*(1.-CTH)-VINT(42)*(1.+CTH))
 369         P(N+2,4)=0.25*VINT(1)*(VINT(41)*(1.-CTH)+VINT(42)*(1.+CTH))
 370         P(N+2,5)=0.
 371
 372         IF(RFLAV.LT.PARP(85).AND.NSTR.GE.1) THEN
 373 C....Choose relevant string pieces to place gluons on.
 374           DO 210 I=N+1,N+2
 375           DMIN=1E8
 376           DO 200 ISTR=1,NSTR
 377           I1=KSTR(ISTR,1)
 378           I2=KSTR(ISTR,2)
 379           DIST=(P(I,4)*P(I1,4)-P(I,1)*P(I1,1)-P(I,2)*P(I1,2)-
 380      &    P(I,3)*P(I1,3))*(P(I,4)*P(I2,4)-P(I,1)*P(I2,1)-
 381      &    P(I,2)*P(I2,2)-P(I,3)*P(I2,3))/MAX(1.,P(I1,4)*P(I2,4)-
 382      &    P(I1,1)*P(I2,1)-P(I1,2)*P(I2,2)-P(I1,3)*P(I2,3))
 383           IF(ISTR.EQ.1.OR.DIST.LT.DMIN) THEN
 384             DMIN=DIST
 385             IST1=I1
 386             IST2=I2
 387             ISTM=ISTR
 388           ENDIF
 389   200     CONTINUE
 390
 391 C....Colour flow adjustments, new string pieces.
 392           IF(K(IST1,4)/MSTU(5).EQ.IST2) K(IST1,4)=MSTU(5)*I+
 393      &    MOD(K(IST1,4),MSTU(5))
 394           IF(MOD(K(IST1,5),MSTU(5)).EQ.IST2) K(IST1,5)=
 395      &    MSTU(5)*(K(IST1,5)/MSTU(5))+I
 396           K(I,5)=MSTU(5)*IST1
 397           K(I,4)=MSTU(5)*IST2
 398           IF(K(IST2,5)/MSTU(5).EQ.IST1) K(IST2,5)=MSTU(5)*I+
 399      &    MOD(K(IST2,5),MSTU(5))
 400           IF(MOD(K(IST2,4),MSTU(5)).EQ.IST1) K(IST2,4)=
 401      &    MSTU(5)*(K(IST2,4)/MSTU(5))+I
 402           KSTR(ISTM,2)=I
 403           KSTR(NSTR+1,1)=I
 404           KSTR(NSTR+1,2)=IST2
 405   210     NSTR=NSTR+1
 406
 407 C...String drawing and colour flow for gluon loop.
 408         ELSEIF(K(N+1,2).EQ.21) THEN
 409           K(N+1,4)=MSTU(5)*(N+2)
 410           K(N+1,5)=MSTU(5)*(N+2)
 411           K(N+2,4)=MSTU(5)*(N+1)
 412           K(N+2,5)=MSTU(5)*(N+1)
 413           KSTR(NSTR+1,1)=N+1
 414           KSTR(NSTR+1,2)=N+2
 415           KSTR(NSTR+2,1)=N+2
 416           KSTR(NSTR+2,2)=N+1
 417           NSTR=NSTR+2
 418
 419 C...String drawing and colour flow for q-qbar pair.
 420         ELSE
 421           K(N+1,4)=MSTU(5)*(N+2)
 422           K(N+2,5)=MSTU(5)*(N+1)
 423           KSTR(NSTR+1,1)=N+1
 424           KSTR(NSTR+1,2)=N+2
 425           NSTR=NSTR+1
 426         ENDIF
 427
 428 C...Update remaining energy; iterate.
 429         N=N+2
 430         IF(N.GT.MSTU(4)-MSTU(32)-10) THEN
 431            CALL LUERRM_HIJING(11
 432      $          ,'(PYMULT_HIJING:) no more memory left in LUJETS_HIJING'
 433      $          )
 434           IF(MSTU(21).GE.1) RETURN
 435         ENDIF
 436         MINT(31)=MINT(31)+1
 437         VINT(151)=VINT(151)+VINT(41)
 438         VINT(152)=VINT(152)+VINT(42)
 439         VINT(143)=VINT(143)-VINT(41)
 440         VINT(144)=VINT(144)-VINT(42)
 441         IF(MINT(31).LT.240) GOTO 180
 442   220   CONTINUE
 443       ENDIF
 444
 445 C...Format statements for printout.
 446  1000 FORMAT(/1X
 447      $     ,'****** PYMULT_HIJING: initialization of multiple inter'
 448      $     ,'actions for MSTP(82) =',I2,' ******')
 449  1100 FORMAT(8X,'pT0 =',F5.2,' GeV gives sigma(parton-parton) =',1P,
 450      &E9.2,' mb: rejected')
 451  1200 FORMAT(8X,'pT0 =',F5.2,' GeV gives sigma(parton-parton) =',1P,
 452      &E9.2,' mb: accepted')
 453
 454       RETURN
 455       END