pyspli.f

Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file pyspli.f

 
C*********************************************************************
 
C...PYSPLI
C...Splits a hadron remnant into two (partons or hadron + parton)
C...in case it is more complicated than just a quark or a diquark.
 
      SUBROUTINE pyspli(KF,KFLIN,KFLCH,KFLSP)
 
C...Double precision and integer declarations.
      IMPLICIT DOUBLE PRECISION(a-h, o-z)
      IMPLICIT INTEGER(i-n)
      INTEGER pyk,pychge,pycomp
C...Commonblocks. PYDAT1 temporary
      common/pypars/mstp(200),parp(200),msti(200),pari(200)
      common/pyint1/mint(400),vint(400)
      common/pydat1/mstu(200),paru(200),mstj(200),parj(200)
      SAVE /pypars/,/pyint1/,/pydat1/
C...Local array.
      dimension kfl(3)
 
C...Preliminaries. Parton composition.
      kfa=iabs(kf)
      kfs=isign(1,kf)
      kfl(1)=mod(kfa/1000,10)
      kfl(2)=mod(kfa/100,10)
      kfl(3)=mod(kfa/10,10)
      IF(kfa.EQ.22.AND.mint(109).EQ.2) THEN
        kfl(2)=int(1.5d0+pyr(0))
        IF(mint(105).EQ.333) kfl(2)=3
        IF(mint(105).EQ.443) kfl(2)=4
        kfl(3)=kfl(2)
      ELSEIF((kfa.EQ.111.OR.kfa.EQ.113).AND.pyr(0).GT.0.5d0) THEN
        kfl(2)=2
        kfl(3)=2
      ELSEIF(kfa.EQ.223.AND.pyr(0).GT.0.5d0) THEN
        kfl(2)=1
        kfl(3)=1
      ELSEIF((kfa.EQ.130.OR.kfa.EQ.310).AND.pyr(0).GT.0.5d0) THEN
        kfl(2)=mod(kfa/10,10)
        kfl(3)=mod(kfa/100,10)
      ENDIF
      IF(kflin.NE.21.AND.kflin.NE.22.AND.kflin.NE.23) THEN
        kflr=kflin*kfs
      ELSE
        kflr=kflin
      ENDIF
      kflch=0
 
C...Subdivide lepton.
      IF(kfa.GE.11.AND.kfa.LE.18) THEN
        IF(kflr.EQ.kfa) THEN
          kflsp=kfs*22
        ELSEIF(kflr.EQ.22) THEN
          kflsp=kfa
        ELSEIF(kflr.EQ.-24.AND.mod(kfa,2).EQ.1) THEN
          kflsp=kfa+1
        ELSEIF(kflr.EQ.24.AND.mod(kfa,2).EQ.0) THEN
          kflsp=kfa-1
        ELSEIF(kflr.EQ.21) THEN
          kflsp=kfa
          kflch=kfs*21
        ELSE
          kflsp=kfa
          kflch=-kflr
        ENDIF
 
C...Subdivide photon.
      ELSEIF(kfa.EQ.22.AND.mint(109).NE.2) THEN
        IF(kflr.NE.21) THEN
          kflsp=-kflr
        ELSE
          ragr=0.75d0*pyr(0)
          kflsp=1
          IF(ragr.GT.0.125d0) kflsp=2
          IF(ragr.GT.0.625d0) kflsp=3
          IF(pyr(0).GT.0.5d0) kflsp=-kflsp
          kflch=-kflsp
        ENDIF
 
C...Subdivide Reggeon or Pomeron.
      ELSEIF(kfa.EQ.110.OR.kfa.EQ.990) THEN
        IF(kflin.EQ.21) THEN
          kflsp=kfs*21
        ELSE
          kflsp=-kflin
        ENDIF
 
C...Subdivide meson.
      ELSEIF(kfl(1).EQ.0) THEN
        kfl(2)=kfl(2)*(-1)**kfl(2)
        kfl(3)=-kfl(3)*(-1)**iabs(kfl(2))
        IF(kflr.EQ.kfl(2)) THEN
          kflsp=kfl(3)
        ELSEIF(kflr.EQ.kfl(3)) THEN
          kflsp=kfl(2)
        ELSEIF(kflr.EQ.21.AND.pyr(0).GT.0.5d0) THEN
          kflsp=kfl(2)
          kflch=kfl(3)
        ELSEIF(kflr.EQ.21) THEN
          kflsp=kfl(3)
          kflch=kfl(2)
        ELSEIF(kflr*kfl(2).GT.0) THEN
          ntry=0
  100     ntry=ntry+1
          CALL pykfdi(-kflr,kfl(2),kfdump,kflch)
          IF(kflch.EQ.0.AND.ntry.LT.100) THEN
            goto 100
          ELSEIF(kflch.EQ.0) THEN
            CALL pyerrm(14,'(PYSPLI:) caught in infinite loop')
            mint(51)=1
            RETURN
          ENDIF
          kflsp=kfl(3)
        ELSE
          ntry=0
  110     ntry=ntry+1
          CALL pykfdi(-kflr,kfl(3),kfdump,kflch)
          IF(kflch.EQ.0.AND.ntry.LT.100) THEN
            goto 110
          ELSEIF(kflch.EQ.0) THEN
            CALL pyerrm(14,'(PYSPLI:) caught in infinite loop')
            mint(51)=1
            RETURN
          ENDIF
          kflsp=kfl(2)
        ENDIF

C...Special case for extracting photon from baryon without splitting
C...the latter. (Currently only used by external programs.)
      ELSEIF(kflin.EQ.22.AND.mstp(98).EQ.1) then
        kflsp=kfa
        kflch=0
 
C...Subdivide baryon.
      ELSE
        nagr=0
        DO 120 j=1,3
          IF(kflr.EQ.kfl(j)) nagr=nagr+1
  120   CONTINUE
        IF(nagr.GE.1) THEN
          ragr=0.00001d0+(nagr-0.00002d0)*pyr(0)
          iagr=0
          DO 130 j=1,3
            IF(kflr.EQ.kfl(j)) ragr=ragr-1d0
            IF(iagr.EQ.0.AND.ragr.LE.0d0) iagr=j
  130     CONTINUE
        ELSE
          iagr=1.00001d0+2.99998d0*pyr(0)
        ENDIF
        id1=1
        IF(iagr.EQ.1) id1=2
        IF(iagr.EQ.1.AND.kfl(3).GT.kfl(2)) id1=3
        id2=6-iagr-id1
        ksp=3
        IF(mod(kfa,10).EQ.2.AND.kfl(1).EQ.kfl(2)) THEN
          IF(iagr.NE.3.AND.pyr(0).GT.0.25d0) ksp=1
        ELSEIF(mod(kfa,10).EQ.2.AND.kfl(2).GE.kfl(3)) THEN
          IF(iagr.NE.1.AND.pyr(0).GT.0.25d0) ksp=1
        ELSEIF(mod(kfa,10).EQ.2) THEN
          IF(iagr.EQ.1) ksp=1
          IF(iagr.NE.1.AND.pyr(0).GT.0.75d0) ksp=1
        ENDIF
        kflsp=1000*kfl(id1)+100*kfl(id2)+ksp
        IF(kflr.EQ.21) THEN
          kflch=kfl(iagr)
        ELSEIF(nagr.EQ.0.AND.kflr.GT.0) THEN
          ntry=0
  140     ntry=ntry+1
          CALL pykfdi(-kflr,kfl(iagr),kfdump,kflch)
          IF(kflch.EQ.0.AND.ntry.LT.100) THEN
            goto 140
          ELSEIF(kflch.EQ.0) THEN
            CALL pyerrm(14,'(PYSPLI:) caught in infinite loop')
            mint(51)=1
            RETURN
          ENDIF
        ELSEIF(nagr.EQ.0) THEN
          ntry=0
  150     ntry=ntry+1
          CALL pykfdi(10000*kfl(id1)+kflsp,-kflr,kfdump,kflch)
          IF(kflch.EQ.0.AND.ntry.LT.100) THEN
            goto 150
          ELSEIF(kflch.EQ.0) THEN
            CALL pyerrm(14,'(PYSPLI:) caught in infinite loop')
            mint(51)=1
            RETURN
          ENDIF
          kflsp=kfl(iagr)
        ENDIF
      ENDIF
 
C...Add on correct sign for result.
      kflch=kflch*kfs
      kflsp=kflsp*kfs
 
      RETURN
      END