HIJING/hijing1_36/hijsft.F

   1 * $Id$
   2 C
   3 C
   4 C*******************************************************************
   5 C                                                                  *
   6 C               Subroutine HIJSFT                                  *
   7 C                                                                  *
   8 C  Scatter two excited strings, JP from proj and JT from target    *
   9 C*******************************************************************
  10         SUBROUTINE HIJSFT(JP,JT,JOUT,IERROR)
  11 #define BLANKET_SAVE
  12 #include "hijcrdn.inc"
  13 #include "hiparnt.inc"
  14 #include "hijdat.inc"
  15 #include "hijjet1.inc"
  16 #include "hijjet2.inc"
  17 #include "histrng.inc"
  18 #include "dpmcom1.inc"
  19 #include "dpmcom2.inc"
  20         SAVE
  21 C*******************************************************************
  22 C       JOUT-> the number
  23 C       of hard scatterings preceding this soft collision.
  24 C       IHNT2(13)-> 1=
  25 C       double diffrac 2=single diffrac, 3=non-single diffrac.
  26 C*******************************************************************
  27         IERROR=0
  28         JJP=JP
  29         JJT=JT
  30         NDPM=0
  31         IOPMAIN=0
  32         IF(JP.GT.IHNT2(1) .OR. JT.GT.IHNT2(3)) RETURN
  33
  34         EPP=PP(JP,4)+PP(JP,3)
  35         EPM=PP(JP,4)-PP(JP,3)
  36         ETP=PT(JT,4)+PT(JT,3)
  37         ETM=PT(JT,4)-PT(JT,3)
  38
  39         WP=EPP+ETP
  40         WM=EPM+ETM
  41         SW=WP*WM
  42 C               ********total W+,W- and center-of-mass energy
  43
  44         IF(WP.LT.0.0 .OR. WM.LT.0.0) GO TO 1000
  45
  46         IF(JOUT.EQ.0) THEN
  47                 IF(EPP.LT.0.0) GO TO 1000
  48                 IF(EPM.LT.0.0) GO TO 1000
  49                 IF(ETP.LT.0.0) GO TO 1000
  50                 IF(ETM.LT.0.0) GO TO 1000
  51                 IF(EPP/(EPM+0.01).LE.ETP/(ETM+0.01)) RETURN
  52         ENDIF
  53 C               ********For strings which does not follow a jet-prod,
  54 C                       scatter only if Ycm(JP)>Ycm(JT). When jets
  55 C                       are produced just before this collision
  56 C                       this requirement has already be enforced
  57 C                       (see SUBROUTINE HIJHRD)
  58         IHNT2(11)=JP
  59         IHNT2(12)=JT
  60 C
  61 C
  62 C
  63         MISS=0
  64         PKC1=0.0
  65         PKC2=0.0
  66         PKC11=0.0
  67         PKC12=0.0
  68         PKC21=0.0
  69         PKC22=0.0
  70         DPKC11=0.0
  71         DPKC12=0.0
  72         DPKC21=0.0
  73         DPKC22=0.0
  74         IF(NFP(JP,10).EQ.1.OR.NFT(JT,10).EQ.1) THEN
  75            IF(NFP(JP,10).EQ.1) THEN
  76               PHI1=ULANGL_HIJING(PP(JP,10),PP(JP,11))
  77               PPJET=SQRT(PP(JP,10)**2+PP(JP,11)**2)
  78               PKC1=PPJET
  79               PKC11=PP(JP,10)
  80               PKC12=PP(JP,11)
  81            ENDIF
  82            IF(NFT(JT,10).EQ.1) THEN
  83               PHI2=ULANGL_HIJING(PT(JT,10),PT(JT,11))
  84               PTJET=SQRT(PT(JT,10)**2+PT(JT,11)**2)
  85               PKC2=PTJET
  86               PKC21=PT(JT,10)
  87               PKC22=PT(JT,11)
  88            ENDIF
  89 C
  90 C Quenching requested ?
  91 C
  92            IF(IHPR2(4).GT.0.AND.IHNT2(1).GT.1.AND.IHNT2(3).GT.1) THEN
  93               IF(NFP(JP,10).EQ.0) THEN
  94                  PHI=-PHI2
  95               ELSE IF(NFT(JT,10).EQ.0) THEN
  96                  PHI=PHI1
  97               ELSE
  98                  IF (PHI1.LT.0. .AND. PHI2.GT.0) THEN ! PH AM
  99                     PHI=(PHI1+PHI2-HIPR1(40))/2.0     ! PH AM
 100                  ELSE                                 ! PH AM
 101                     PHI=(PHI1+PHI2+HIPR1(40))/2.0     ! PH AM
 102                  ENDIF                                ! PH AM
 103               ENDIF
 104               BX=HINT1(19)*COS(HINT1(20))
 105               BY=HINT1(19)*SIN(HINT1(20))
 106               XP0=YP(1,JP)
 107               YP0=YP(2,JP)
 108               XT0=YT(1,JT)+BX
 109               YT0=YT(2,JT)+BY
 110               R1=MAX(1.2*IHNT2(1)**0.3333333,
 111      &               SQRT(XP0**2+YP0**2))
 112               R2=MAX(1.2*IHNT2(3)**0.3333333,
 113      &               SQRT((XT0-BX)**2+(YT0-BY)**2))
 114               IF(ABS(COS(PHI)).LT.1.0E-5) THEN
 115 C AM             DD1=R1
 116 C AM             DD2=R1
 117                  DD1 = ABS(SQRT(R1**2-(XT0-BX)**2) - YT0)
 118                  DD2 = ABS(SQRT(R1**2-(XT0-BX)**2) + YT0)
 119 C PH             DD3=ABS(BY+SQRT(R2**2-(XP0-BX)**2)-YP0)
 120 C PH             DD4=ABS(BY-SQRT(R2**2-(XP0-BX)**2)-YP0)
 121                  DD3=ABS(BY+SQRT(R2**2-(XT0-BX)**2)-YP0)
 122                  DD4=ABS(BY-SQRT(R2**2-(XT0-BX)**2)-YP0)
 123                  GO TO 5
 124               ENDIF
 125               BB=2.0*SIN(PHI)*(COS(PHI)*YP0-SIN(PHI)*XP0)
 126               CC=(YP0**2-R1**2)*COS(PHI)**2+XP0*SIN(PHI)*(
 127      &                          XP0*SIN(PHI)-2.0*YP0*COS(PHI))
 128               DD=BB**2-4.0*CC
 129               IF(DD.LT.0.0) GO TO 10
 130               XX1=(-BB+SQRT(DD))/2.0
 131               XX2=(-BB-SQRT(DD))/2.0
 132               DD1=ABS((XX1-XP0)/COS(PHI))
 133               DD2=ABS((XX2-XP0)/COS(PHI))
 134 C
 135               BB=2.0*SIN(PHI)*(COS(PHI)*(YT0-BY)-SIN(PHI)*XT0)-2.0*BX
 136               CC=(BX**2+(YT0-BY)**2-R2**2)*COS(PHI)**2+XT0*SIN(PHI)
 137      &           *(XT0*SIN(PHI)-2.0*COS(PHI)*(YT0-BY))
 138      &           -2.0*BX*SIN(PHI)*(COS(PHI)*(YT0-BY)-SIN(PHI)*XT0)
 139               DD=BB**2-4.0*CC
 140               IF(DD.LT.0.0) GO TO 10
 141               XX1=(-BB+SQRT(DD))/2.0
 142               XX2=(-BB-SQRT(DD))/2.0
 143               DD3=ABS((XX1-XT0)/COS(PHI))
 144               DD4=ABS((XX2-XT0)/COS(PHI))
 145 C
 146  5            DD1=MIN(DD1,DD3)
 147               DD2=MIN(DD2,DD4)
 148               IF(DD1.LT.HIPR1(13)) DD1=0.0
 149               IF(DD2.LT.HIPR1(13)) DD2=0.0
 150               IF(NFP(JP,10).EQ.1.AND.PPJET.GT.HIPR1(11)) THEN
 151                  IF (IHPR2(50) .EQ. 1) THEN
 152                     DEDX0 = HIPR1(14)*LOG10(PPJET)/LOG10(5.)
 153                  ELSE
 154                     DEDX0 = HIPR1(14)
 155                  ENDIF
 156                  DP1=DD1*DEDX0/2.0
 157                  DP1=MIN(DP1,PPJET-HIPR1(11))
 158                  PKC1=PPJET-DP1
 159                  DPX1=COS(PHI1)*DP1
 160                  DPY1=SIN(PHI1)*DP1
 161                  PKC11=PP(JP,10)-DPX1
 162                  PKC12=PP(JP,11)-DPY1
 163                  IF(DP1.GT.0.0) THEN
 164                     CTHEP=PP(JP,12)/SQRT(PP(JP,12)**2+PPJET**2)
 165                     DPZ1=DP1*CTHEP/SQRT(1.0-CTHEP**2)
 166                     DPE1=SQRT(DPX1**2+DPY1**2+DPZ1**2)
 167                     EPPPRM=PP(JP,4)+PP(JP,3)-DPE1-DPZ1
 168                     EPMPRM=PP(JP,4)-PP(JP,3)-DPE1+DPZ1
 169                     IF(EPPPRM.LE.0.0.OR.EPMPRM.LE.0.0) GO TO 15
 170                     EPP=EPPPRM
 171                     EPM=EPMPRM
 172                     PP(JP,10)=PKC11
 173                     PP(JP,11)=PKC12
 174                     NPJ(JP)=NPJ(JP)+1
 175                     KFPJ(JP,NPJ(JP))=21
 176                     PJPX(JP,NPJ(JP))=DPX1
 177                     PJPY(JP,NPJ(JP))=DPY1
 178                     PJPZ(JP,NPJ(JP))=DPZ1
 179                     PJPE(JP,NPJ(JP))=DPE1
 180                     PJPM(JP,NPJ(JP))=0.0
 181                     PP(JP,3)=PP(JP,3)-DPZ1
 182                     PP(JP,4)=PP(JP,4)-DPE1
 183                  ENDIF
 184               ENDIF
 185  15           IF(NFT(JT,10).EQ.1.AND.PTJET.GT.HIPR1(11)) THEN
 186                  IF (IHPR2(50) .EQ. 1) THEN
 187                     DEDX0 = HIPR1(14)*LOG10(PTJET)/LOG10(5.)
 188                  ELSE
 189                     DEDX0 = HIPR1(14)
 190                  ENDIF
 191                  DP2=DD2*DEDX0/2.0
 192                  DP2=MIN(DP2,PTJET-HIPR1(11))
 193                  PKC2=PTJET-DP2
 194                  DPX2=COS(PHI2)*DP2
 195                  DPY2=SIN(PHI2)*DP2
 196                  PKC21=PT(JT,10)-DPX2
 197                  PKC22=PT(JT,11)-DPY2
 198                  IF(DP2.GT.0.0) THEN
 199                     CTHET=PT(JT,12)/SQRT(PT(JT,12)**2+PTJET**2)
 200                     DPZ2=DP2*CTHET/SQRT(1.0-CTHET**2)
 201                     DPE2=SQRT(DPX2**2+DPY2**2+DPZ2**2)
 202                     ETPPRM=PT(JT,4)+PT(JT,3)-DPE2-DPZ2
 203                     ETMPRM=PT(JT,4)-PT(JT,3)-DPE2+DPZ2
 204                     IF(ETPPRM.LE.0.0.OR.ETMPRM.LE.0.0) GO TO 16
 205                     ETP=ETPPRM
 206                     ETM=ETMPRM
 207                     PT(JT,10)=PKC21
 208                     PT(JT,11)=PKC22
 209                     NTJ(JT)=NTJ(JT)+1
 210                     KFTJ(JT,NTJ(JT))=21
 211                     PJTX(JT,NTJ(JT))=DPX2
 212                     PJTY(JT,NTJ(JT))=DPY2
 213                     PJTZ(JT,NTJ(JT))=DPZ2
 214                     PJTE(JT,NTJ(JT))=DPE2
 215                     PJTM(JT,NTJ(JT))=0.0
 216                     PT(JT,3)=PT(JT,3)-DPZ2
 217                     PT(JT,4)=PT(JT,4)-DPE2
 218                  ENDIF
 219               ENDIF
 220  16           CONTINUE
 221               PKC11=0.
 222               PKC12=0.
 223               PKC21=0.
 224               PKC22=0.
 225               DPKC11=0.
 226               DPKC12=0.
 227               DPKC21=0.
 228               DPKC22=0.
 229
 230               WP=EPP+ETP
 231               WM=EPM+ETM
 232               SW=WP*WM
 233            ENDIF
 234         ENDIF
 235
 236         IF(NFP(JP,10).EQ.1) THEN
 237            PPJET=SQRT(PP(JP,10)**2+PP(JP,11)**2)
 238            PKC1=PPJET
 239         ENDIF
 240         IF(NFT(JT,10).EQ.1) THEN
 241            PTJET=SQRT(PT(JT,10)**2+PT(JT,11)**2)
 242            PKC2=PTJET
 243         ENDIF
 244 C               ********If jet is quenched the pt from valence quark
 245 C                       hard scattering has to reduced by d*kapa
 246 C
 247 C
 248 10      PTP02=PP(JP,1)**2+PP(JP,2)**2
 249         PTT02=PT(JT,1)**2+PT(JT,2)**2
 250 C
 251         AMQ=MAX(PP(JP,14)+PP(JP,15),PT(JT,14)+PT(JT,15))
 252         AMX=HIPR1(1)+AMQ
 253 C               ********consider mass cut-off for strings which
 254 C                       must also include quark's mass
 255         AMP0=AMX
 256         DPM0=AMX
 257         NFDP=0
 258         IF(NFP(JP,5).LE.2.AND.NFP(JP,3).NE.0) THEN
 259                 AMP0=ULMASS_HIJING(NFP(JP,3))
 260                 NFDP=NFP(JP,3)+2*NFP(JP,3)/ABS(NFP(JP,3))
 261                 DPM0=ULMASS_HIJING(NFDP)
 262                 IF(DPM0.LE.0.0) THEN
 263                         NFDP=NFDP-2*NFDP/ABS(NFDP)
 264                         DPM0=ULMASS_HIJING(NFDP)
 265                 ENDIF
 266         ENDIF
 267         AMT0=AMX
 268         DTM0=AMX
 269         NFDT=0
 270         IF(NFT(JT,5).LE.2.AND.NFT(JT,3).NE.0) THEN
 271                 AMT0=ULMASS_HIJING(NFT(JT,3))
 272                 NFDT=NFT(JT,3)+2*NFT(JT,3)/ABS(NFT(JT,3))
 273                 DTM0=ULMASS_HIJING(NFDT)
 274                 IF(DTM0.LE.0.0) THEN
 275                         NFDT=NFDT-2*NFDT/ABS(NFDT)
 276                         DTM0=ULMASS_HIJING(NFDT)
 277                 ENDIF
 278         ENDIF
 279 C
 280         AMPN=SQRT(AMP0**2+PTP02)
 281         AMTN=SQRT(AMT0**2+PTT02)
 282         SNN=(AMPN+AMTN)**2+0.001
 283 C
 284         IF(SW.LT.SNN+0.001) GO TO 4000
 285 C               ********Scatter only if SW>SNN
 286 C*****give some PT kick to the two exited strings******************
 287 20      SWPTN=4.0*(MAX(AMP0,AMT0)**2+MAX(PTP02,PTT02))
 288         SWPTD=4.0*(MAX(DPM0,DTM0)**2+MAX(PTP02,PTT02))
 289         SWPTX=4.0*(AMX**2+MAX(PTP02,PTT02))
 290         IF(SW.LE.SWPTN) THEN
 291                 PKCMX=0.0
 292         ELSE IF(SW.GT.SWPTN .AND. SW.LE.SWPTD
 293      &          .AND.NPJ(JP).EQ.0.AND.NTJ(JT).EQ.0) THEN
 294            PKCMX=SQRT(SW/4.0-MAX(AMP0,AMT0)**2)
 295      &           -SQRT(MAX(PTP02,PTT02))
 296         ELSE IF(SW.GT.SWPTD .AND. SW.LE.SWPTX
 297      &          .AND.NPJ(JP).EQ.0.AND.NTJ(JT).EQ.0) THEN
 298            PKCMX=SQRT(SW/4.0-MAX(DPM0,DTM0)**2)
 299      &           -SQRT(MAX(PTP02,PTT02))
 300         ELSE IF(SW.GT.SWPTX) THEN
 301            PKCMX=SQRT(SW/4.0-AMX**2)-SQRT(MAX(PTP02,PTT02))
 302         ENDIF
 303 C               ********maximun PT kick
 304 C*********************************************************
 305         IF(NFP(JP,10).EQ.1.OR.NFT(JT,10).EQ.1) THEN
 306                 IF(PKC1.GT.PKCMX) THEN
 307                         PKC1=PKCMX
 308                         PKC11=PKC1*COS(PHI1)
 309                         PKC12=PKC1*SIN(PHI1)
 310                         DPKC11=-(PP(JP,10)-PKC11)/2.0
 311                         DPKC12=-(PP(JP,11)-PKC12)/2.0
 312                 ENDIF
 313                 IF(PKC2.GT.PKCMX) THEN
 314                         PKC2=PKCMX
 315                         PKC21=PKC2*COS(PHI2)
 316                         PKC22=PKC2*SIN(PHI2)
 317                         DPKC21=-(PT(JT,10)-PKC21)/2.0
 318                         DPKC22=-(PT(JT,11)-PKC22)/2.0
 319                 ENDIF
 320                 DPKC1=DPKC11+DPKC21
 321                 DPKC2=DPKC12+DPKC22
 322                 NFP(JP,10)=-NFP(JP,10)
 323                 NFT(JT,10)=-NFT(JT,10)
 324                 GO TO 40
 325         ENDIF
 326 C               ********If the valence quarks had a hard-collision
 327 C                       the pt kick is the pt from hard-collision.
 328         I_SNG=0
 329         IF(IHPR2(13).NE.0 .AND. RLU_HIJING(0).LE.HIDAT(4)) I_SNG=1
 330         IF((NFP(JP,5).EQ.3 .OR.NFT(JT,5).EQ.3).OR.
 331      &          (NPJ(JP).NE.0.OR.NFP(JP,10).NE.0).OR.
 332      &          (NTJ(JT).NE.0.OR.NFT(JT,10).NE.0)) I_SNG=0
 333 C
 334 C               ********decite whether to have single-diffractive
 335         IF(IHPR2(5).EQ.0) THEN
 336                 PKC=HIPR1(2)*SQRT(-ALOG(1.0-RLU_HIJING(0)
 337      &                  *(1.0-EXP(-PKCMX**2/HIPR1(2)**2))))
 338                 GO TO 30
 339         ENDIF
 340         PKC=HIRND2(3,0.0,PKCMX**2)
 341         PKC=SQRT(PKC)
 342         IF(PKC.GT.HIPR1(20))
 343      &     PKC=HIPR1(2)*SQRT(-ALOG(EXP(-HIPR1(20)**2/HIPR1(2)**2)
 344      &         -RLU_HIJING(0)*(EXP(-HIPR1(20)**2/HIPR1(2)**2)-
 345      &         EXP(-PKCMX**2/HIPR1(2)**2))))
 346 C
 347         IF(I_SNG.EQ.1) PKC=0.65*SQRT(
 348      &          -ALOG(1.0-RLU_HIJING(0)*(1.0-EXP(-PKCMX**2/0.65**2))))
 349 C                       ********select PT kick
 350 30      PHI0=2.0*HIPR1(40)*RLU_HIJING(0)
 351         PKC11=PKC*SIN(PHI0)
 352         PKC12=PKC*COS(PHI0)
 353         PKC21=-PKC11
 354         PKC22=-PKC12
 355         DPKC1=0.0
 356         DPKC2=0.0
 357 40      PP11=PP(JP,1)+PKC11-DPKC1
 358         PP12=PP(JP,2)+PKC12-DPKC2
 359         PT11=PT(JT,1)+PKC21-DPKC1
 360         PT12=PT(JT,2)+PKC22-DPKC2
 361         PTP2=PP11**2+PP12**2
 362         PTT2=PT11**2+PT12**2
 363 C
 364         AMPN=SQRT(AMP0**2+PTP2)
 365         AMTN=SQRT(AMT0**2+PTT2)
 366         SNN=(AMPN+AMTN)**2+0.001
 367 C***************************************
 368         WP=EPP+ETP
 369         WM=EPM+ETM
 370         SW=WP*WM
 371 C****************************************
 372         IF(SW.LT.SNN) THEN
 373            MISS=MISS+1
 374            IF(MISS.LE.100) then
 375               PKC=0.0
 376               GO TO 30
 377            ENDIF
 378            IF(IHPR2(10).NE.0)
 379      &       WRITE(6,*) 'Error occured in Pt kick section of HIJSFT'
 380            GO TO 4000
 381         ENDIF
 382 C******************************************************************
 383         AMPD=SQRT(DPM0**2+PTP2)
 384         AMTD=SQRT(DTM0**2+PTT2)
 385
 386         AMPX=SQRT(AMX**2+PTP2)
 387         AMTX=SQRT(AMX**2+PTT2)
 388
 389         DPN=AMPN**2/SW
 390         DTN=AMTN**2/SW
 391         DPD=AMPD**2/SW
 392         DTD=AMTD**2/SW
 393         DPX=AMPX**2/SW
 394         DTX=AMTX**2/SW
 395 C
 396         SPNTD=(AMPN+AMTD)**2
 397         SPNTX=(AMPN+AMTX)**2
 398 C                       ********CM energy if proj=N,targ=N*
 399         SPDTN=(AMPD+AMTN)**2
 400         SPXTN=(AMPX+AMTN)**2
 401 C                       ********CM energy if proj=N*,targ=N
 402         SPDTX=(AMPD+AMTX)**2
 403         SPXTD=(AMPX+AMTD)**2
 404         SDD=(AMPD+AMTD)**2
 405         SXX=(AMPX+AMTX)**2
 406
 407 C
 408 C
 409 C               ********CM energy if proj=delta, targ=delta
 410 C****************There are many different cases**********
 411 c       IF(IHPR2(15).EQ.1) GO TO 500
 412 C
 413 C               ********to have DPM type soft interactions
 414 C
 415  45     CONTINUE
 416         IF(SW.GT.SXX+0.001) THEN
 417            IF(I_SNG.EQ.0) THEN
 418               D1=DPX
 419               D2=DTX
 420               NFP3=0
 421               NFT3=0
 422               GO TO 400
 423            ELSE
 424 c**** 5/30/1998 this is identical to the above statement. Added to
 425 c**** avoid questional branching to block.
 426               IF((NFP(JP,5).EQ.3 .AND.NFT(JT,5).EQ.3).OR.
 427      &           (NPJ(JP).NE.0.OR.NFP(JP,10).NE.0).OR.
 428      &           (NTJ(JT).NE.0.OR.NFT(JT,10).NE.0)) THEN
 429                  D1=DPX
 430                  D2=DTX
 431                  NFP3=0
 432                  NFT3=0
 433                  GO TO 400
 434               ENDIF
 435 C               ********do not allow excited strings to have
 436 C                       single-diffr
 437               IF(RLU_HIJING(0).GT.0.5.OR.(NFT(JT,5).GT.2.OR.
 438      &                NTJ(JT).NE.0.OR.NFT(JT,10).NE.0)) THEN
 439                  D1=DPN
 440                  D2=DTX
 441                  NFP3=NFP(JP,3)
 442                  NFT3=0
 443                  GO TO 220
 444               ELSE
 445                  D1=DPX
 446                  D2=DTN
 447                  NFP3=0
 448                  NFT3=NFT(JT,3)
 449                  GO TO 240
 450               ENDIF
 451 C               ********have single diffractive collision
 452            ENDIF
 453         ELSE IF(SW.GT.MAX(SPDTX,SPXTD)+0.001 .AND.
 454      &                          SW.LE.SXX+0.001) THEN
 455            IF(((NPJ(JP).EQ.0.AND.NTJ(JT).EQ.0.AND.
 456      &         RLU_HIJING(0).GT.0.5).OR.(NPJ(JP).EQ.0
 457      &         .AND.NTJ(JT).NE.0)).AND.NFP(JP,5).LE.2) THEN
 458               D1=DPD
 459               D2=DTX
 460               NFP3=NFDP
 461               NFT3=0
 462               GO TO 220
 463            ELSE IF(NTJ(JT).EQ.0.AND.NFT(JT,5).LE.2) THEN
 464               D1=DPX
 465               D2=DTD
 466               NFP3=0
 467               NFT3=NFDT
 468               GO TO 240
 469            ENDIF
 470            GO TO 4000
 471         ELSE IF(SW.GT.MIN(SPDTX,SPXTD)+0.001.AND.
 472      &                  SW.LE.MAX(SPDTX,SPXTD)+0.001) THEN
 473            IF(SPDTX.LE.SPXTD.AND.NPJ(JP).EQ.0
 474      &                       .AND.NFP(JP,5).LE.2) THEN
 475               D1=DPD
 476               D2=DTX
 477               NFP3=NFDP
 478               NFT3=0
 479               GO TO 220
 480            ELSE IF(SPDTX.GT.SPXTD.AND.NTJ(JT).EQ.0
 481      &                       .AND.NFT(JT,5).LE.2) THEN
 482               D1=DPX
 483               D2=DTD
 484               NFP3=0
 485               NFT3=NFDT
 486               GO TO 240
 487            ENDIF
 488 c*** 5/30/1998 added to avoid questional branching to another block
 489 c*** this is identical to the statement following the next ELSE IF
 490            IF(((NPJ(JP).EQ.0.AND.NTJ(JT).EQ.0
 491      &       .AND.RLU_HIJING(0).GT.0.5).OR.(NPJ(JP).EQ.0
 492      &        .AND.NTJ(JT).NE.0)).AND.NFP(JP,5).LE.2) THEN
 493               D1=DPN
 494               D2=DTX
 495               NFP3=NFP(JP,3)
 496               NFT3=0
 497               GO TO 220
 498            ELSE IF(NTJ(JT).EQ.0.AND.NFT(JT,5).LE.2) THEN
 499               D1=DPX
 500               D2=DTN
 501               NFP3=0
 502               NFT3=NFT(JT,3)
 503               GO TO 240
 504            ENDIF
 505            GO TO 4000
 506         ELSE IF(SW.GT.MAX(SPNTX,SPXTN)+0.001 .AND.
 507      &                  SW.LE.MIN(SPDTX,SPXTD)+0.001) THEN
 508            IF(((NPJ(JP).EQ.0.AND.NTJ(JT).EQ.0
 509      &       .AND.RLU_HIJING(0).GT.0.5).OR.(NPJ(JP).EQ.0
 510      &        .AND.NTJ(JT).NE.0)).AND.NFP(JP,5).LE.2) THEN
 511               D1=DPN
 512               D2=DTX
 513               NFP3=NFP(JP,3)
 514               NFT3=0
 515               GO TO 220
 516            ELSE IF(NTJ(JT).EQ.0.AND.NFT(JT,5).LE.2) THEN
 517               D1=DPX
 518               D2=DTN
 519               NFP3=0
 520               NFT3=NFT(JT,3)
 521               GO TO 240
 522            ENDIF
 523            GO TO 4000
 524         ELSE IF(SW.GT.MIN(SPNTX,SPXTN)+0.001 .AND.
 525      &                  SW.LE.MAX(SPNTX,SPXTN)+0.001) THEN
 526            IF(SPNTX.LE.SPXTN.AND.NPJ(JP).EQ.0
 527      &                           .AND.NFP(JP,5).LE.2) THEN
 528               D1=DPN
 529               D2=DTX
 530               NFP3=NFP(JP,3)
 531               NFT3=0
 532               GO TO 220
 533            ELSEIF(SPNTX.GT.SPXTN.AND.NTJ(JT).EQ.0
 534      &                           .AND.NFT(JT,5).LE.2) THEN
 535               D1=DPX
 536               D2=DTN
 537               NFP3=0
 538               NFT3=NFT(JT,3)
 539               GO TO 240
 540            ENDIF
 541            GO TO 4000
 542         ELSE IF(SW.LE.MIN(SPNTX,SPXTN)+0.001 .AND.
 543      &                  (NPJ(JP).NE.0 .OR.NTJ(JT).NE.0)) THEN
 544            GO TO 4000
 545         ELSE IF(SW.LE.MIN(SPNTX,SPXTN)+0.001 .AND.
 546      &          NFP(JP,5).GT.2.AND.NFT(JT,5).GT.2) THEN
 547            GO TO 4000
 548         ELSE IF(SW.GT.SDD+0.001.AND.SW.LE.
 549      &                     MIN(SPNTX,SPXTN)+0.001) THEN
 550            D1=DPD
 551            D2=DTD
 552            NFP3=NFDP
 553            NFT3=NFDT
 554            GO TO 100
 555         ELSE IF(SW.GT.MAX(SPNTD,SPDTN)+0.001
 556      &                      .AND. SW.LE.SDD+0.001) THEN
 557            IF(RLU_HIJING(0).GT.0.5) THEN
 558               D1=DPD
 559               D2=DTN
 560               NFP3=NFDP
 561               NFT3=NFT(JT,3)
 562               GO TO 100
 563            ELSE
 564               D1=DPN
 565               D2=DTD
 566               NFP3=NFP(JP,3)
 567               NFT3=NFDT
 568               GO TO 100
 569            ENDIF
 570         ELSE IF(SW.GT.MIN(SPNTD,SPDTN)+0.001
 571      &          .AND. SW.LE.MAX(SPNTD,SPDTN)+0.001) THEN
 572            IF(SPNTD.GT.SPDTN) THEN
 573               D1=DPD
 574               D2=DTN
 575               NFP3=NFDP
 576               NFT3=NFT(JT,3)
 577               GO TO 100
 578            ELSE
 579               D1=DPN
 580               D2=DTD
 581               NFP3=NFP(JP,3)
 582               NFT3=NFDT
 583               GO TO 100
 584            ENDIF
 585         ELSE IF(SW.LE.MIN(SPNTD,SPDTN)+0.001) THEN
 586            D1=DPN
 587            D2=DTN
 588            NFP3=NFP(JP,3)
 589            NFT3=NFT(JT,3)
 590            GO TO 100
 591         ENDIF
 592         WRITE(6,*) ' Error in HIJSFT: There is no path to here'
 593         RETURN
 594 C
 595 C***************  elastic scattering ***************
 596 C       this is like elastic, both proj and targ mass
 597 C       must be fixed
 598 C***************************************************
 599 100     NFP5=MAX(2,NFP(JP,5))
 600         NFT5=MAX(2,NFT(JT,5))
 601         BB1=1.0+D1-D2
 602         BB2=1.0+D2-D1
 603         IF(BB1**2.LT.4.0*D1 .OR. BB2**2.LT.4.0*D2) THEN
 604                 MISS=MISS+1
 605                 IF(MISS.GT.100.OR.PKC.EQ.0.0) GO TO 3000
 606                 PKC=PKC*0.5
 607                 GO TO 30
 608         ENDIF
 609         IF(RLU_HIJING(0).LT.0.5) THEN
 610                 X1=(BB1-SQRT(BB1**2-4.0*D1))/2.0
 611                 X2=(BB2-SQRT(BB2**2-4.0*D2))/2.0
 612         ELSE
 613                 X1=(BB1+SQRT(BB1**2-4.0*D1))/2.0
 614                 X2=(BB2+SQRT(BB2**2-4.0*D2))/2.0
 615         ENDIF
 616         IHNT2(13)=2
 617         GO TO 600
 618 C
 619 C********** Single diffractive ***********************
 620 C either proj or targ's mass is fixed
 621 C*****************************************************
 622 220     NFP5=MAX(2,NFP(JP,5))
 623         NFT5=3
 624         IF(NFP3.EQ.0) NFP5=3
 625         BB2=1.0+D2-D1
 626         IF(BB2**2.LT.4.0*D2) THEN
 627                 MISS=MISS+1
 628                 IF(MISS.GT.100.OR.PKC.EQ.0.0) GO TO 3000
 629                 PKC=PKC*0.5
 630                 GO TO 30
 631         ENDIF
 632         XMIN=(BB2-SQRT(BB2**2-4.0*D2))/2.0
 633         XMAX=(BB2+SQRT(BB2**2-4.0*D2))/2.0
 634         MISS4=0
 635 222     X2=HIRND2(6,XMIN,XMAX)
 636         X1=D1/(1.0-X2)
 637         IF(X2*(1.0-X1).LT.(D2+1.E-4/SW)) THEN
 638                 MISS4=MISS4+1
 639                 IF(MISS4.LE.1000) GO TO 222
 640                 GO TO 5000
 641         ENDIF
 642         IHNT2(13)=2
 643         GO TO 600
 644 C                       ********Fix proj mass*********
 645 240     NFP5=3
 646         NFT5=MAX(2,NFT(JT,5))
 647         IF(NFT3.EQ.0) NFT5=3
 648         BB1=1.0+D1-D2
 649         IF(BB1**2.LT.4.0*D1) THEN
 650                 MISS=MISS+1
 651                 IF(MISS.GT.100.OR.PKC.EQ.0.0) GO TO 3000
 652                 PKC=PKC*0.5
 653                 GO TO 30
 654         ENDIF
 655         XMIN=(BB1-SQRT(BB1**2-4.0*D1))/2.0
 656         XMAX=(BB1+SQRT(BB1**2-4.0*D1))/2.0
 657         MISS4=0
 658 242     X1=HIRND2(6,XMIN,XMAX)
 659         X2=D2/(1.0-X1)
 660         IF(X1*(1.0-X2).LT.(D1+1.E-4/SW)) THEN
 661                 MISS4=MISS4+1
 662                 IF(MISS4.LE.1000) GO TO 242
 663                 GO TO 5000
 664         ENDIF
 665         IHNT2(13)=2
 666         GO TO 600
 667 C                       ********Fix targ mass*********
 668 C
 669 C*************non-single diffractive**********************
 670 C       both proj and targ may not be fixed in mass
 671 C*********************************************************
 672 C
 673 400     NFP5=3
 674         NFT5=3
 675         BB1=1.0+D1-D2
 676         BB2=1.0+D2-D1
 677         IF(BB1**2.LT.4.0*D1 .OR. BB2**2.LT.4.0*D2) THEN
 678                 MISS=MISS+1
 679                 IF(MISS.GT.100.OR.PKC.EQ.0.0) GO TO 3000
 680                 PKC=PKC*0.5
 681                 GO TO 30
 682         ENDIF
 683         XMIN1=(BB1-SQRT(BB1**2-4.0*D1))/2.0
 684         XMAX1=(BB1+SQRT(BB1**2-4.0*D1))/2.0
 685         XMIN2=(BB2-SQRT(BB2**2-4.0*D2))/2.0
 686         XMAX2=(BB2+SQRT(BB2**2-4.0*D2))/2.0
 687         MISS4=0
 688 410     X1=HIRND2(4,XMIN1,XMAX1)
 689         X2=HIRND2(4,XMIN2,XMAX2)
 690         IF(NFP(JP,5).EQ.3.OR.NFT(JT,5).EQ.3) THEN
 691                 X1=HIRND2(6,XMIN1,XMAX1)
 692                 X2=HIRND2(6,XMIN2,XMAX2)
 693         ENDIF
 694 C                       ********
 695         IF(ABS(NFP(JP,1)*NFP(JP,2)).GT.1000000.OR.
 696      &                  ABS(NFP(JP,1)*NFP(JP,2)).LT.100) THEN
 697                 X1=HIRND2(5,XMIN1,XMAX1)
 698         ENDIF
 699         IF(ABS(NFT(JT,1)*NFT(JT,2)).GT.1000000.OR.
 700      &                  ABS(NFT(JT,1)*NFT(JT,2)).LT.100) THEN
 701                 X2=HIRND2(5,XMIN2,XMAX2)
 702         ENDIF
 703 c       IF(IOPMAIN.EQ.3) X1=HIRND2(6,XMIN1,XMAX1)
 704 c       IF(IOPMAIN.EQ.2) X2=HIRND2(6,XMIN2,XMAX2)
 705 C       ********For q-qbar or (qq)-(qq)bar system use symetric
 706 C               distribution, for q-(qq) or qbar-(qq)bar use
 707 C               unsymetrical distribution
 708 C
 709         IF(ABS(NFP(JP,1)*NFP(JP,2)).GT.1000000) X1=1.0-X1
 710         XXP=X1*(1.0-X2)
 711         XXT=X2*(1.0-X1)
 712         IF(XXP.LT.(D1+1.E-4/SW) .OR. XXT.LT.(D2+1.E-4/SW)) THEN
 713                 MISS4=MISS4+1
 714                 IF(MISS4.LE.1000) GO TO 410
 715                 GO TO 5000
 716         ENDIF
 717         IHNT2(13)=3
 718 C***************************************************
 719 C***************************************************
 720 600     CONTINUE
 721         IF(X1*(1.0-X2).LT.(AMPN**2-1.E-4)/SW.OR.
 722      &                  X2*(1.0-X1).LT.(AMTN**2-1.E-4)/SW) THEN
 723                 MISS=MISS+1
 724                 IF(MISS.GT.100.OR.PKC.EQ.0.0) GO TO 2000
 725                 PKC=0.0
 726                 GO TO 30
 727         ENDIF
 728 C
 729         EPP=(1.0-X2)*WP
 730         EPM=X1*WM
 731         ETP=X2*WP
 732         ETM=(1.0-X1)*WM
 733         PP(JP,3)=(EPP-EPM)/2.0
 734         PP(JP,4)=(EPP+EPM)/2.0
 735         IF(EPP*EPM-PTP2.LT.0.0) GO TO 6000
 736         PP(JP,5)=SQRT(EPP*EPM-PTP2)
 737         NFP(JP,3)=NFP3
 738         NFP(JP,5)=NFP5
 739
 740         PT(JT,3)=(ETP-ETM)/2.0
 741         PT(JT,4)=(ETP+ETM)/2.0
 742         IF(ETP*ETM-PTT2.LT.0.0) GO TO 6000
 743         PT(JT,5)=SQRT(ETP*ETM-PTT2)
 744         NFT(JT,3)=NFT3
 745         NFT(JT,5)=NFT5
 746 C*****recoil PT from hard-inter is shared by two end-partons
 747 C       so that pt=p1+p2
 748         PP(JP,1)=PP11-PKC11
 749         PP(JP,2)=PP12-PKC12
 750
 751         KICKDIP=1
 752         KICKDIT=1
 753         IF(ABS(NFP(JP,1)*NFP(JP,2)).GT.1000000.OR.
 754      &                  ABS(NFP(JP,1)*NFP(JP,2)).LT.100) THEN
 755                 KICKDIP=0
 756         ENDIF
 757         IF(ABS(NFT(JT,1)*NFT(JT,2)).GT.1000000.OR.
 758      &                  ABS(NFT(JT,1)*NFT(JT,2)).LT.100) THEN
 759                 KICKDIT=0
 760         ENDIF
 761         IF((KICKDIP.EQ.0.AND.RLU_HIJING(0).LT.0.5)
 762      &     .OR.(KICKDIP.NE.0.AND.RLU_HIJING(0)
 763      &     .LT.0.5/(1.0+(PKC11**2+PKC12**2)/HIPR1(22)**2))) THEN
 764            PP(JP,6)=(PP(JP,1)-PP(JP,6)-PP(JP,8)-DPKC1)/2.0+PP(JP,6)
 765            PP(JP,7)=(PP(JP,2)-PP(JP,7)-PP(JP,9)-DPKC2)/2.0+PP(JP,7)
 766            PP(JP,8)=(PP(JP,1)-PP(JP,6)-PP(JP,8)-DPKC1)/2.0
 767      &              +PP(JP,8)+PKC11
 768            PP(JP,9)=(PP(JP,2)-PP(JP,7)-PP(JP,9)-DPKC2)/2.0
 769      &              +PP(JP,9)+PKC12
 770         ELSE
 771            PP(JP,8)=(PP(JP,1)-PP(JP,6)-PP(JP,8)-DPKC1)/2.0+PP(JP,8)
 772            PP(JP,9)=(PP(JP,2)-PP(JP,7)-PP(JP,9)-DPKC2)/2.0+PP(JP,9)
 773            PP(JP,6)=(PP(JP,1)-PP(JP,6)-PP(JP,8)-DPKC1)/2.0
 774      &              +PP(JP,6)+PKC11
 775            PP(JP,7)=(PP(JP,2)-PP(JP,7)-PP(JP,9)-DPKC2)/2.0
 776      &              +PP(JP,7)+PKC12
 777         ENDIF
 778         PP(JP,1)=PP(JP,6)+PP(JP,8)
 779         PP(JP,2)=PP(JP,7)+PP(JP,9)
 780 C                               ********pt kick for proj
 781         PT(JT,1)=PT11-PKC21
 782         PT(JT,2)=PT12-PKC22
 783         IF((KICKDIT.EQ.0.AND.RLU_HIJING(0).LT.0.5)
 784      &     .OR.(KICKDIT.NE.0.AND.RLU_HIJING(0)
 785      &     .LT.0.5/(1.0+(PKC21**2+PKC22**2)/HIPR1(22)**2))) THEN
 786            PT(JT,6)=(PT(JT,1)-PT(JT,6)-PT(JT,8)-DPKC1)/2.0+PT(JT,6)
 787            PT(JT,7)=(PT(JT,2)-PT(JT,7)-PT(JT,9)-DPKC2)/2.0+PT(JT,7)
 788            PT(JT,8)=(PT(JT,1)-PT(JT,6)-PT(JT,8)-DPKC1)/2.0
 789      &              +PT(JT,8)+PKC21
 790            PT(JT,9)=(PT(JT,2)-PT(JT,7)-PT(JT,9)-DPKC2)/2.0
 791      &              +PT(JT,9)+PKC22
 792         ELSE
 793            PT(JT,8)=(PT(JT,1)-PT(JT,6)-PT(JT,8)-DPKC1)/2.0+PT(JT,8)
 794            PT(JT,9)=(PT(JT,2)-PT(JT,7)-PT(JT,9)-DPKC2)/2.0+PT(JT,9)
 795            PT(JT,6)=(PT(JT,1)-PT(JT,6)-PT(JT,8)-DPKC1)/2.0
 796      &              +PT(JT,6)+PKC21
 797            PT(JT,7)=(PT(JT,2)-PT(JT,7)-PT(JT,9)-DPKC2)/2.0
 798      &              +PT(JT,7)+PKC22
 799         ENDIF
 800         PT(JT,1)=PT(JT,6)+PT(JT,8)
 801         PT(JT,2)=PT(JT,7)+PT(JT,9)
 802 C                       ********pt kick for targ
 803
 804         IF(NPJ(JP).NE.0) NFP(JP,5)=3
 805         IF(NTJ(JT).NE.0) NFT(JT,5)=3
 806 C                       ********jets must be connected to string
 807         IF(EPP/(EPM+0.0001).LT.ETP/(ETM+0.0001).AND.
 808      &                  ABS(NFP(JP,1)*NFP(JP,2)).LT.1000000)THEN
 809                 DO 620 JSB=1,15
 810                 PSB=PP(JP,JSB)
 811                 PP(JP,JSB)=PT(JT,JSB)
 812                 PT(JT,JSB)=PSB
 813                 NSB=NFP(JP,JSB)
 814                 NFP(JP,JSB)=NFT(JT,JSB)
 815                 NFT(JT,JSB)=NSB
 816 620             CONTINUE
 817 C               ********when Ycm(JP)<Ycm(JT) after the collision
 818 C                       exchange the positions of the two
 819         ENDIF
 820 C
 821         RETURN
 822 C**************************************************
 823 C**************************************************
 824 1000    IERROR=1
 825         IF(IHPR2(10).EQ.0) RETURN
 826         WRITE(6,*) '    Fatal HIJSFT start error,abandon this event'
 827         WRITE(6,*) '    PROJ E+,E-,W+',EPP,EPM,WP
 828         WRITE(6,*) '    TARG E+,E-,W-',ETP,ETM,WM
 829         WRITE(6,*) '    W+*W-, (APN+ATN)^2',SW,SNN
 830         RETURN
 831 2000    IERROR=0
 832         IF(IHPR2(10).EQ.0) RETURN
 833         WRITE(6,*) '    (2)energy partition fail,'
 834         WRITE(6,*) '    HIJSFT not performed, but continue'
 835         WRITE(6,*) '    MP1,MPN',X1*(1.0-X2)*SW,AMPN**2
 836         WRITE(6,*) '    MT2,MTN',X2*(1.0-X1)*SW,AMTN**2
 837         RETURN
 838 3000    IERROR=0
 839         IF(IHPR2(10).EQ.0) RETURN
 840         WRITE(6,*) '    (3)something is wrong with the pt kick, '
 841         WRITE(6,*) '    HIJSFT not performed, but continue'
 842         WRITE(6,*) '    D1=',D1,' D2=',D2,' SW=',SW
 843         WRITE(6,*) '    HISTORY NFP5=',NFP(JP,5),' NFT5=',NFT(JT,5)
 844         WRITE(6,*) '    THIS COLLISON NFP5=',NFP5, ' NFT5=',NFT5
 845         WRITE(6,*) '    # OF JET IN PROJ',NPJ(JP),' IN TARG',NTJ(JT)
 846         RETURN
 847 4000    IERROR=0
 848         IF(IHPR2(10).EQ.0) RETURN
 849         WRITE(6,*) '    (4)unable to choose process, but not harmful'
 850         WRITE(6,*) '    HIJSFT not performed, but continue'
 851         WRITE(6,*) '    PTP=',SQRT(PTP2),' PTT=',SQRT(PTT2),' SW=',SW
 852         WRITE(6,*) '    AMCUT=',AMX,' JP=',JP,' JT=',JT
 853         WRITE(6,*) '    HISTORY NFP5=',NFP(JP,5),' NFT5=',NFT(JT,5)
 854         RETURN
 855 5000    IERROR=0
 856         IF(IHPR2(10).EQ.0) RETURN
 857         WRITE(6,*) '    energy partition failed(5),for limited try'
 858         WRITE(6,*) '    HIJSFT not performed, but continue'
 859         WRITE(6,*) '    NFP5=',NFP5,' NFT5=',NFT5
 860         WRITE(6,*) '    D1',D1,' X1(1-X2)',X1*(1.0-X2)
 861         WRITE(6,*) '    D2',D2,' X2(1-X1)',X2*(1.0-X1)
 862         RETURN
 863 6000    PKC=0.0
 864         MISS=MISS+1
 865         IF(MISS.LT.100) GO TO 30
 866         IERROR=1
 867         IF(IHPR2(10).EQ.0) RETURN
 868         WRITE(6,*) ' ERROR OCCURED, HIJSFT NOT PERFORMED'
 869         WRITE(6,*) ' Abort this event'
 870         WRITE(6,*) 'MTP,PTP2',EPP*EPM,PTP2,'  MTT,PTT2',ETP*ETM,PTT2
 871         RETURN
 872         END