*
* $Id: hwcfor.F,v 1.1.1.1 1996/03/08 17:02:11 mclareni Exp $
*
* $Log: hwcfor.F,v $
* Revision 1.1.1.1  1996/03/08 17:02:11  mclareni
* Herwig58
*
*
*CMZ :          29/08/94  11.51.46  by  Unknown
*-- Author :
CDECK  ID>, HWCFOR.
*CMZ :-        -26/04/91  14.15.56  by  Federico Carminati
*-- Author :    Bryan Webber
C------------------------------------------------------------------------
      SUBROUTINE HWCFOR
C     CONVERTS COLOUR-CONNECTED QUARK-ANTIQUARK PAIRS INTO CLUSTERS
C
C     IF CLUSTER IS BELOW THRESHOLD FOR TWO-HADRON DECAY
C     IT IS SET EQUAL TO THE LOWEST-MASS HADRON OF THAT
C     FLAVOUR AND 4-MOMENTUM IS TRANSFERRED TO ANOTHER
C     CLUSTER
C
C     CHANGED 29/1/92 BY BRW: 1-HADRON DECAYS OF B CLUSTERS
C     ALLOWED UP TO (1+B1LIM)*THRESHOLD (LINEAR DECREASE)
C
C     CHANGED  5/2/92 BY BRW: INTRODUCED CLUSTER SPLITTING
C     PARAMETER CLPOW.  CLUSTERS OF MASS MCL MADE FROM
C     PARTONS ID1 AND ID2 ARE SPLIT IF
C
C     MCL**CLPOW.GT.CLMAX**CLPOW+(RMASS(ID1)+RMASS(ID2))**CLPOW
C
C     THUS PREVIOUS VERSIONS CORRESPONDED TO CLPOW=2.
C
C     CHANGED  5/11/92 BY BRW: DIQUARK-ANTIDIQUARK CLUSTERS BELOW
C     THRESHOLD FOR BARYON-ANTIBARYON DECAY ARE SHIFTED ABOVE
C------------------------------------------------------------------------
#include "herwig58/herwig58.inc"
      LOGICAL BMCLUS,TGCLUS,BTCLUS,SPLIT,HWRLOG
      INTEGER IHEP,IBHEP,IPASS,IBCL,JCL,JHEP,ID1,ID2,ID3,KHEP,
     &        LCL,MCL,L,I
      DOUBLE PRECISION HWUPCM,EMCLP,EM0,EM1,EM2,PC0,PC1,PCL(5)
      IF (IERROR.NE.0) RETURN
C---SPLIT GLUONS
      CALL HWCGSP
C---FIND COLOUR PARTNERS AFTER BARYON NUMBER VIOLATING EVENT
      IF (HVFCEN) CALL HVCBVI
      IF (IERROR.NE.0) RETURN
C---LOOK FOR PARTONS TO CLUSTER
      DO 10 IHEP=1,NHEP
      IF (ISTHEP(IHEP).GE.150.AND.ISTHEP(IHEP).LE.154) GO TO 11
   10 CONTINUE
      IBCL=1
      GO TO 26
   11 IBHEP=IHEP
C---FIRST PASS CHOPS UP BEAM/TARGET CLUSTERS
C   SECOND PASS CHOPS UP MASSIVE PAIRS
C   THIRD PASS CREATES CLUSTERS
      IPASS=0
C---LOOK FOR QUARK OR ANTIDIQUARK
   14 IPASS=IPASS+1
      IF (IPASS.GT.3) GO TO 25
      IBCL=NHEP+1
      JCL=NHEP
      JHEP=IBHEP-1
   15 JHEP=JHEP+1
      IF (JHEP.GT.NHEP) GO TO 14
      IF (ISTHEP(JHEP).LT.150.OR.ISTHEP(JHEP).GT.154) GO TO 15
      ID1=IDHW(JHEP)
      IF (ID1.GT.6.AND.ID1.LT.115) GO TO 15
      IF (ID1.GT.120) GO TO 15
C---FIND CONNECTED ANTIQUARK OR DIQUARK
      KHEP=JMOHEP(2,JHEP)
      BMCLUS=ISTHEP(JHEP).EQ.153 .OR.ISTHEP(KHEP).EQ.153
      TGCLUS=ISTHEP(JHEP).EQ.154 .OR.ISTHEP(KHEP).EQ.154
      IF (IPASS.EQ.1) THEN
        BTCLUS=(BMCLUS.OR.TGCLUS).AND.
     &       ISTHEP(JHEP).NE.151.AND.ISTHEP(KHEP).NE.151
C---CHOP UP BEAM, TARGET CLUSTERS
        IF (BTCLUS) THEN
          CALL HWVSUM(4,PHEP(1,JHEP),PHEP(1,KHEP),PCL)
          CALL HWUMAS(PCL)
          CALL HWCCUT(JHEP,KHEP,PCL,BTCLUS,SPLIT)
        ENDIF
      ELSE
C---Q-QBAR PAIR FOUND. CHECK MASS/MAKE CLUSTER
        ID2=IDHW(KHEP)
        CALL HWVSUM(4,PHEP(1,JHEP),PHEP(1,KHEP),PCL)
        CALL HWUMAS(PCL)
        IF (IPASS.EQ.2) THEN
          EMCLP=SIGN(ABS(PCL(5))**CLPOW,PCL(5))
          IF (EMCLP.GT.CLMAX**CLPOW+
     &       (RMASS(ID1)+RMASS(ID2))**CLPOW) THEN
            CALL HWCCUT(JHEP,KHEP,PCL,.FALSE.,SPLIT)
            IF (SPLIT) THEN
              JHEP=JHEP-1
              GO TO 15
            ENDIF
          ENDIF
        ELSE
          JCL=JCL+1
          IDHW(JCL)=19
          IDHEP(JCL)=91
          CALL HWVEQU(5,PCL,PHEP(1,JCL))
          IF (BMCLUS) THEN
            ISTHEP(JCL)=164
          ELSE IF (TGCLUS) THEN
            ISTHEP(JCL)=165
          ELSE
            ISTHEP(JCL)=163
          ENDIF
          JMOHEP(1,JCL)=JHEP
          JMOHEP(2,JCL)=KHEP
          JDAHEP(1,JCL)=0
          JDAHEP(2,JCL)=0
          JDAHEP(1,JHEP)=JCL
          JDAHEP(1,KHEP)=JCL
          ISTHEP(JHEP)=ISTHEP(JHEP)+8
          ISTHEP(KHEP)=ISTHEP(KHEP)+8
        ENDIF
      ENDIF
      GO TO 15
   25 NHEP=JCL
C---FIX UP MOMENTA FOR SINGLE-HADRON CLUSTERS
   26 JCL=IBCL-1
   30 JCL=JCL+1
      IF (JCL.GT.NHEP) GO TO 45
      IF (ISTHEP(JCL).LT.163.OR.ISTHEP(JCL).GT.165) GO TO 30
C---DON'T HADRONIZE BEAM/TARGET CLUSTERS
      IF (ISTHEP(JCL).NE.163.AND.GENSOF) GO TO 30
C---SET UP FLAVOURS FOR CLUSTER DECAY
      CALL HWCFLA(IDHW(JMOHEP(1,JCL)),IDHW(JMOHEP(2,JCL)),ID1,ID3)
      EM0=PHEP(5,JCL)
      IF ((B1LIM.EQ.ZERO).OR.(ID1.NE.11.AND.ID3.NE.11))THEN
        IF (EM0.GT.RMIN(ID1,2)+RMIN(2,ID3)) GO TO 30
      ELSE
C---SPECIAL FOR B CLUSTERS: ALLOW 1-HADRON DECAY ABOVE THRESHOLD
        IF (HWRLOG((EM0/(RMIN(ID1,2)+RMIN(2,ID3))
     &     -ONE)/B1LIM)) GO TO 30
      ENDIF
      EM1=RMIN(ID1,ID3)
C---DIQUARK-ANTIDIQUARK CLUSTERS WILL RETURN EM1=0
C   SHIFT THEM ABOVE THRESHOLD FOR BARYON-ANTIBARYON
      IF (EM1.EQ.ZERO) EM1=RMIN(ID1,2)+RMIN(2,ID3)+1D-2
      IF (ABS(EM0-EM1).LT.1D-3) GO TO 30
C---DECIDE TO GO BACK OR FORWARD TO TRANSFER 4-MOMENTUM
      L=1
      IF (HWRLOG(HALF)) L=-1
      MCL=NHEP-IBCL+1
      LCL=JCL
      DO 40 I=1,MCL
      LCL=LCL+L
      IF (LCL.LT.IBCL) LCL=LCL+MCL
      IF (LCL.GT.NHEP) LCL=LCL-MCL
      IF (LCL.EQ.JCL) THEN
        IF (EM0.GE.EM1+RMIN(1,1)) GO TO 30
        CALL HWWARN('HWCFOR',101,*999)
      ENDIF
      IF (ISTHEP(LCL).LT.163.OR.ISTHEP(LCL).GT.165) GO TO 40
C---RESCALE MOMENTA IN 2-CLUSTER C.M.
      CALL HWVSUM(4,PHEP(1,JCL),PHEP(1,LCL),PCL)
      CALL HWUMAS(PCL)
      EM2=PHEP(5,LCL)
      PC0=HWUPCM(PCL(5),EM0,EM2)
      PC1=HWUPCM(PCL(5),EM1,EM2)
      IF (PC1.LT.0.) THEN
C---NEED TO RESCALE OTHER MASS AS WELL
        CALL HWCFLA(IDHW(JMOHEP(1,LCL)),IDHW(JMOHEP(2,LCL)),ID1,ID3)
        EM2=RMIN(ID1,ID3)
C---DIQUARK-ANTIDIQUARK CLUSTERS WILL RETURN EM2=0
C   SHIFT THEM ABOVE THRESHOLD FOR BARYON-ANTIBARYON
        IF (EM2.EQ.ZERO) EM2=RMIN(ID1,2)+RMIN(2,ID3)+1D-2
        PC1=HWUPCM(PCL(5),EM1,EM2)
        IF (PC1.LT.0.) GO TO 40
        PHEP(5,LCL)=EM2
      ENDIF
      IF (PC0.GT.0.) THEN
        PC0=PC1/PC0
        CALL HWULOF(PCL,PHEP(1,JCL),PHEP(1,JCL))
        CALL HWVSCA(3,PC0,PHEP(1,JCL),PHEP(1,JCL))
        PHEP(4,JCL)=SQRT(PC1**2+EM1**2)
        PHEP(5,JCL)=EM1
        CALL HWULOB(PCL,PHEP(1,JCL),PHEP(1,JCL))
        CALL HWVDIF(4,PCL,PHEP(1,JCL),PHEP(1,LCL))
        GO TO 30
      ELSEIF (PC0.EQ.0.) THEN
        PHEP(5,JCL)=EM1
        CALL HWDTWO(PCL,PHEP(1,JCL),PHEP(1,LCL),PC1,TWO,.TRUE.)
        GO TO 30
      ELSE
        CALL HWWARN('HWCFOR',103,*999)
      ENDIF
   40 CONTINUE
      CALL HWWARN('HWCFOR',102,*999)
   45 ISTAT=60
C---NON-PARTONS LABELLED AS PARTONS (IE PHOTONS) SHOULD GET COPIED
      DO 50 IHEP=1,NHEP
        IF (ISTHEP(IHEP).EQ.150) THEN
          NHEP=NHEP+1
          JDAHEP(1,IHEP)=NHEP
          ISTHEP(IHEP)=157
          ISTHEP(NHEP)=190
          IDHW(NHEP)=IDHW(IHEP)
          IDHEP(NHEP)=IDPDG(IDHW(IHEP))
          CALL HWVEQU(5,PHEP(1,IHEP),PHEP(1,NHEP))
          JMOHEP(1,NHEP)=IHEP
          JMOHEP(2,NHEP)=JMOHEP(1,IHEP)
          JDAHEP(1,NHEP)=0
          JDAHEP(2,NHEP)=0
        ENDIF
 50   CONTINUE
  999 END
