* * $Id: hwdhig.F,v 1.1.1.1 1996/03/08 17:02:11 mclareni Exp $ * * $Log: hwdhig.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>, HWDHIG. *CMZ :- -24/04/92 14.23.44 by Mike Seymour *-- Author : Mike Seymour C----------------------------------------------------------------------- SUBROUTINE HWDHIG(GAMINP) C HIGGS DECAY ROUTINE C A) FOR GAMinp=0 FIND AND DECAY HIGGS C B) FOR GAMinp>0 CALCULATE TOTAL HIGGS WIDTH C FOR EMH=GAMINP. STORE RESULT IN GAMINP. C----------------------------------------------------------------------- #include "herwig58/herwig58.inc" DOUBLE PRECISION GAMINP,EMH,EMF,COLFAC,ENF,K1,K0,BET0,BET1,GAM0, & GAM1,SCLOG,CFAC,XF,EM,GAMLIM,GAM,XW,EMW,XZ,EMZ,YW,YZ,EMI, & TAUT,TAUW,WIDHIG,VECDEC,EMB,GAMB,TMIN,TMAX1,EM1,TMAX2,EM2,X1,X2, & PROB,PCM,SUMR,SUMI,TAUTR,TAUTI,TAUWR,TAUWI,GFACTR, & HWDHGF,HWRGEN,HWRUNI,HWUSQR,HWUPCM INTEGER IHIG,I,IFERM,NLOOK,I1,I2,IPART,IMODE,IDEC,MMAX,HWRINT LOGICAL HWRLOG PARAMETER (NLOOK=100) DIMENSION VECDEC(2,0:NLOOK) SAVE GAM,EM,VECDEC EQUIVALENCE (EMW,RMASS(198)) EQUIVALENCE (EMZ,RMASS(200)) DATA GAMLIM,GAM,EM/10D0,2*0D0/ C---IF DECAY, FIND HIGGS (HWDHAD WILL HAVE GIVEN IT STATUS=1) IF (GAMINP.EQ.0) THEN IHIG=0 DO 10 I=1,NHEP 10 IF (IHIG.EQ.0.AND.IDHW(I).EQ.201.AND.ISTHEP(I).EQ.1) IHIG=I IF (IHIG.EQ.0) CALL HWWARN('HWDHIG',101,*999) EMH=PHEP(5,IHIG) IF (EMH.LE.0) CALL HWWARN('HWDHIG',102,*999) EMSCA=EMH ELSE EMH=GAMINP IF (EMH.LE.0) THEN GAMINP=0 RETURN ENDIF ENDIF C---CALCULATE BRANCHING FRACTIONS C---FERMIONS C---NLL CORRECTION TO QUARK DECAY RATE (HHG eq 2.6-9) ENF=0 DO 1 I=1,6 1 IF (2*RMASS(I).LT.EMH) ENF=ENF+1 K1=5/PIFAC**2 K0=3/(4*PIFAC**2) BET0=(11*CAFAC-2*ENF)/3 BET1=(34*CAFAC**2-(10*CAFAC+6*CFFAC)*ENF)/3 GAM0=-8 GAM1=-404./3+40*ENF/9 SCLOG=LOG(EMH**2/QCDLAM**2) CFAC=1 + ( K1/K0 - 2*GAM0 + GAM0*BET1/BET0**2*LOG(SCLOG) & + (GAM0*BET1-GAM1*BET0)/BET0**2) / (BET0*SCLOG) DO 100 IFERM=1,9 IF (IFERM.LE.6) THEN EMF=RMASS(IFERM) XF=(EMF/EMH)**2 COLFAC=FLOAT(NCOLO) IF (EMF.GT.QCDLAM) & EMF=EMF*(LOG(EMH/QCDLAM)/LOG(EMF/QCDLAM))**(GAM0/(2*BET0)) ELSE EMF=RMASS(107+IFERM*2) XF=(EMF/EMH)**2 COLFAC=1 CFAC=1 ENDIF IF (XF.LT.0.25) THEN GFACTR=ALPHEM/(8.*SWEIN*EMW**2) BRHIG(IFERM)=COLFAC*GFACTR*EMH*EMF**2 * (1-4*XF)**1.5 * CFAC ELSE BRHIG(IFERM)=0 ENDIF 100 CONTINUE C---W*W*/Z*Z* IF (ABS(EM-EMH).GE.GAMLIM*GAM) THEN C---OFF EDGE OF LOOK-UP TABLE XW=(EMW/EMH)**2 XZ=(EMZ/EMH)**2 YW=EMW*GAMW/EMH**2 YZ=EMZ*GAMZ/EMH**2 BRHIG(10)=.50*GFACTR * EMH**3 * HWDHGF(XW,YW) BRHIG(11)=.25*GFACTR * EMH**3 * HWDHGF(XZ,YZ) ELSE C---LOOK IT UP EMI=((EMH-EM)/(GAM*GAMLIM)+1)*NLOOK/2.0 I1=INT(EMI) I2=INT(EMI+1) BRHIG(10)=.50*GFACTR * EMH**3 * ( VECDEC(1,I1)*(I2-EMI) + & VECDEC(1,I2)*(EMI-I1) ) BRHIG(11)=.25*GFACTR * EMH**3 * ( VECDEC(2,I1)*(I2-EMI) + & VECDEC(2,I2)*(EMI-I1) ) ENDIF C---GAMMAGAMMA TAUT=(2*RMASS(6)/EMH)**2 TAUW=(2*EMW/EMH)**2 CALL HWDHGC(TAUT,TAUTR,TAUTI) CALL HWDHGC(TAUW,TAUWR,TAUWI) SUMR=4./3*( - 2*TAUT*( 1 + (1-TAUT)*TAUTR ) ) * ENHANC(6) & +(2 + 3*TAUW*( 1 + (2-TAUW)*TAUWR ) ) * ENHANC(10) SUMI=4./3*( - 2*TAUT*( (1-TAUT)*TAUTI ) ) * ENHANC(6) & +( 3*TAUW*( (2-TAUW)*TAUWI ) ) * ENHANC(10) BRHIG(12)=GFACTR*.03125*(ALPHEM/PIFAC)**2 & *EMH**3 * (SUMR**2 + SUMI**2) WIDHIG=0 DO 200 IPART=1, 12 IF (IPART.LT.12) BRHIG(IPART)=BRHIG(IPART)*ENHANC(IPART)**2 200 WIDHIG=WIDHIG+BRHIG(IPART) IF (WIDHIG.EQ.0) CALL HWWARN('HWDHIG',103,*999) DO 300 IPART=1, 12 300 BRHIG(IPART)=BRHIG(IPART)/WIDHIG IF (EM.NE.RMASS(201)) THEN C---SET UP W*W*/Z*Z* LOOKUP TABLES EM=EMH GAM=WIDHIG GAMLIM=MAX(GAMLIM,GAMMAX) DO 400 I=0,NLOOK EMH=(I*2.0/NLOOK-1)*GAM*GAMLIM+EM XW=(EMW/EMH)**2 XZ=(EMZ/EMH)**2 YW=EMW*GAMW/EMH**2 YZ=EMZ*GAMZ/EMH**2 VECDEC(1,I)=HWDHGF(XW,YW) VECDEC(2,I)=HWDHGF(XZ,YZ) 400 CONTINUE EMH=EM ENDIF IF (GAMINP.GT.0) THEN GAMINP=WIDHIG RETURN ENDIF C---SEE IF USER SPECIFIED A DECAY MODE IMODE=MOD(IPROC,100) C---IF NOT, CHOOSE ONE IF (IMODE.LT.1.OR.IMODE.GT.12) THEN MMAX=12 IF (IMODE.LT.1) MMAX=6 500 IMODE=HWRINT(1,MMAX) IF (BRHIG(IMODE).LT.HWRGEN(0)) GOTO 500 ENDIF C---SEE IF SPECIFIED DECAY IS POSSIBLE IF (BRHIG(IMODE).EQ.0) CALL HWWARN('HWDHIG',104,*999) IF (IMODE.LE.6) THEN IDEC=IMODE ELSEIF (IMODE.LE.9) THEN IDEC=107+IMODE*2 ELSEIF (IMODE.EQ.10) THEN IDEC=198 ELSEIF (IMODE.EQ.11) THEN IDEC=200 ELSEIF (IMODE.EQ.12) THEN IDEC=59 ENDIF C---STATUS, IDs AND POINTERS ISTHEP(IHIG)=195 DO 600 I=1,2 ISTHEP(NHEP+I)=193 IDHW(NHEP+I)=IDEC IDHEP(NHEP+I)=IDPDG(IDEC) JDAHEP(I,IHIG)=NHEP+I JMOHEP(1,NHEP+I)=IHIG JMOHEP(2,NHEP+I)=NHEP+(3-I) JDAHEP(2,NHEP+I)=NHEP+(3-I) PHEP(5,NHEP+I)=RMASS(IDEC) IDEC=IDEC+6 IF (IDEC.EQ.204) IDEC=199 IF (IDEC.EQ.206) IDEC=200 IF (IDEC.EQ. 65) IDEC= 59 600 CONTINUE C---ALLOW W/Z TO BE OFF-SHELL IF (IMODE.EQ.10.OR.IMODE.EQ.11) THEN IF (IMODE.EQ.10) THEN EMB=EMW GAMB=GAMW ELSE EMB=EMZ GAMB=GAMZ ENDIF C---STANDARD MASS DISTRIBUTION 700 TMIN=ATAN(-EMB/GAMB) TMAX1=ATAN((EMH**2/EMB-EMB)/GAMB) EM1=HWUSQR(EMB*(GAMB*TAN(HWRUNI(0,TMIN,TMAX1))+EMB)) TMAX2=ATAN(((EMH-EM1)**2/EMB-EMB)/GAMB) EM2=HWUSQR(EMB*(GAMB*TAN(HWRUNI(0,TMIN,TMAX2))+EMB)) X1=(EM1/EMH)**2 X2=(EM2/EMH)**2 C---CORRECT MASS DISTRIBUTION PROB=HWUSQR(1+X1**2+X2**2-2*X1-2*X2-2*X1*X2) & * ((X1+X2-1)**2 + 8*X1*X2) IF (.NOT.HWRLOG(PROB)) GOTO 700 C---CALCULATE SPIN DENSITY MATRIX RHOHEP(1,NHEP+1)=4*X1*X2 / (8*X1*X2 + (X1+X2-1)**2) RHOHEP(2,NHEP+1)=(X1+X2-1)**2 / (8*X1*X2 + (X1+X2-1)**2) RHOHEP(3,NHEP+1)=RHOHEP(1,NHEP+1) C---SYMMETRIZE DISTRIBUTIONS IN PARTICLES 1,2 IF (HWRLOG(HALF)) THEN PHEP(5,NHEP+1)=EM1 PHEP(5,NHEP+2)=EM2 ELSE PHEP(5,NHEP+1)=EM2 PHEP(5,NHEP+2)=EM1 ENDIF ENDIF C---DO DECAY PCM=HWUPCM(EMH,PHEP(5,NHEP+1),PHEP(5,NHEP+2)) IF (PCM.LT.0) CALL HWWARN('HWDHIG',105,*999) CALL HWDTWO(PHEP(1,IHIG),PHEP(1,NHEP+1),PHEP(1,NHEP+2), & PCM,TWO,.TRUE.) NHEP=NHEP+2 C---IF QUARK DECAY, HADRONIZE IF (IMODE.LE.6) THEN ISTHEP(NHEP-1)=113 ISTHEP(NHEP)=114 CALL HWBGEN CALL HWDHQK CALL HWCFOR CALL HWCDEC ENDIF 999 END