[u/mrichter/AliRoot.git] / ISAJET / code / setdky.F

#include "isajet/pilot.h"
      SUBROUTINE SETDKY(LPRINT)
C
C          Read in decay table from tape ITDKY and set up /DKYTAB/.
C          Then append forced decay modes and set LOOK to negative
C          number pointing to LOOK2, which points to table.
C          Forced decays for antiparticles are stored in conjugated 
C          form so that DECAY can always conjugate them.
C
C          Logical flag LPRINT controls printing of table.
C
C          Ver 7.41: Check version of decay table. Also read matrix 
C          element flags and save in MELEM:
C          MELEM=0: Phase space
C          MELEM=1: Dalitz decay
C          MELEM=2: omega/phi decay
C          MELEM=3: V-A
C          MELEM=4: V-A plus W propagator (for top)
C          MELEM=5: tau -> ell nu nu
C          MELEM=6: tau -> nu pi/K
C          MELEM=7: tau -> nu rho/a1
C
#if defined(CERNLIB_IMPNONE)
      IMPLICIT NONE
#endif
#include "isajet/itapes.inc"
#include "isajet/force.inc"
#include "isajet/dkytab.inc"
#include "isajet/nodcay.inc"
#include "isajet/ssmode.inc"
#include "isajet/sstype.inc"
#include "isajet/xmssm.inc"
#include "isajet/keys.inc"
C
      INTEGER IMODE(6),LOOP,IOLD,I,IRES,ITYPE,K,J,IPOINT
      INTEGER IFL1,IFL2,IFL3,JSPIN,INDEX,ID1,IDANTI,KTYPE,IRES2
      REAL    BR
      CHARACTER*8 LABEL,LMODE(6),LRES
      CHARACTER*8 IBLANK,LREAD(10),IQUIT
      LOGICAL LPRINT
      INTEGER NOUT,NTHAD
      PARAMETER (NOUT=33)
      PARAMETER (NTHAD=12)
      INTEGER IDOUT(NOUT),ITHAD(NTHAD),IDUMMY(5),MEOUT
      REAL SUMBR,SUMBR2,SUMGAM
      CHARACTER*40,V,VOLD,VISAJE
C
      DATA IDOUT/
     $IDTP,ISGL,ISUPL,ISDNL,ISSTL,ISCHL,ISBT1,ISTP1,ISUPR,ISDNR,
     $ISSTR,ISCHR,ISBT2,ISTP2,ISEL,ISMUL,ISTAU1,ISNEL,ISNML,ISNTL,
     $ISER,ISMUR,ISTAU2,ISZ1,ISZ2,ISZ3,ISZ4,ISW1,ISW2,
     $ISHL,ISHH,ISHA,ISHC/
      DATA IQUIT/'////'/,IBLANK/' '/
      DATA ITHAD/-160,-260,-360,
     $  1160,1260,2260,2160,1360,2360,3160,3260,3360/
C
C          Print header for table.
C
      IF(LPRINT) WRITE(ITLIS,10)
10    FORMAT('1',30('*')/' *',28X,'*'/
     1' *',5X,'ISAJET DECAY TABLE',5X,'*'/
     2' *',28X,'*'/' ',30('*')//
     33X,'PART',16X,'DECAY MODE',16X,'CUM BR',10X,'IDENT',18X,
     4'DECAY IDENT'/)
C
C          Initialize. LOOP is the decay mode counter.
C
      LOOP=0
      IOLD=0
      DO 100 I=1,MXLOOK
        LOOK(I)=0
100   CONTINUE
      DO 110 I=1,MXFORC
        LOOK2(1,I)=0
        LOOK2(2,I)=0
110   CONTINUE
C
C          Read in table, checking for valid version.
C
      IF(NODCAY.OR.ITDKY.EQ.0) RETURN
      REWIND ITDKY
C
      VOLD=VISAJE()
      READ(ITDKY,*) V
      IF(V.NE.VOLD) THEN
        WRITE(ITLIS,2000) V,VOLD
2000    FORMAT(//
     $  '  ***WARNING: DECAY TABLE DOES NOT MATCH ISAJET VERSION'/
     $  '  ***DECAY VERSION  : ',A40/
     $  '  ***PROGRAM VERSION: ',A40)
      ENDIF
C
200   LOOP=LOOP+1
      IF(LOOP.GT.MXDKY) GO TO 9999
220   DO 210 I=1,5
        IMODE(I)=0
        LMODE(I)=IBLANK
210   CONTINUE
      READ(ITDKY,*) IRES,ITYPE,BR,IMODE
C
      IF(IRES.NE.0) THEN
        IF(NOPI0.AND.IRES.EQ.110) GO TO 220
        IF(NOETA.AND.IRES.EQ.220) GO TO 220
        IF(IRES.NE.IOLD) THEN
          CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
          LOOK(INDEX)=LOOP
        ENDIF
        IOLD=IRES
        CBR(LOOP)=BR
        MELEM(LOOP)=ITYPE
        DO 240 I=1,5
          MODE(I,LOOP)=IMODE(I)
          IF(IMODE(I).NE.0) LMODE(I)=LABEL(IMODE(I))
240     CONTINUE
        LRES=LABEL(IRES)
        IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
     1  BR,IRES,(IMODE(K),K=1,5)
20      FORMAT(3X,A5,4X,5(A5,2X),F8.5,10X,I5,4X,5(I7,2X))
        GO TO 200
      ENDIF
C
C          Add HIGGS FOR WHIGGS
C
      IF(KEYS(10).AND..NOT.GOMSSM) THEN
          SUMGAM=0
          SUMBR=0
          DO 244 J=1,NSSMOD
            IF(ISSMOD(J).EQ.81.AND.GSSMOD(J).GT.0) THEN
              SUMGAM=SUMGAM+GSSMOD(J)
            ENDIF
244       CONTINUE
          DO 245 J=1,NSSMOD
            IF(ISSMOD(J).EQ.81.AND.GSSMOD(J).GT.0) THEN
              BSSMOD(J)=GSSMOD(J)/SUMGAM
            ENDIF
245       CONTINUE
          DO 246 J=1,NSSMOD
            IF(ISSMOD(J).EQ.81.AND.BSSMOD(J).GT.0) THEN
              SUMBR=SUMBR+BSSMOD(J)
            ENDIF
246       CONTINUE
C          If modes exist, add them
          IF(SUMBR.LE.0) GO TO 249
          IRES=81
          LRES=LABEL(IRES)
          CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
          LOOK(INDEX)=LOOP+1
          SUMBR2=0
          DO 247 J=1,NSSMOD
            IF(ISSMOD(J).EQ.81.AND.BSSMOD(J).GT.0) THEN
              LOOP=LOOP+1
              SUMBR2=SUMBR2+BSSMOD(J)
              BR=SUMBR2/SUMBR
              CBR(LOOP)=BR
              MELEM(LOOP)=0
              DO 248 K=1,5
                MODE(K,LOOP)=JSSMOD(K,J)
                LMODE(K)=LABEL(MODE(K,LOOP))
248           CONTINUE
              IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
     $        BR,IRES,(MODE(K,LOOP),K=1,5)
            ENDIF
247       CONTINUE
249     CONTINUE
      END IF
C
C          Add MSSM decay modes if applicable, OR H_SM FOR WHIGGS
C
      IF(GOMSSM) THEN
        DO 250 I=1,NOUT
C          Check for modes
          SUMBR=0
          DO 251 J=1,NSSMOD
            IF(ISSMOD(J).EQ.IDOUT(I).AND.BSSMOD(J).GT.0) THEN
              SUMBR=SUMBR+BSSMOD(J)
            ENDIF
251       CONTINUE
C          If modes exist, add them
          IF(SUMBR.LE.0) GO TO 250
          IRES=IDOUT(I)
          LRES=LABEL(IRES)
          CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
          LOOK(INDEX)=LOOP+1
          SUMBR2=0
          DO 252 J=1,NSSMOD
            IF(ISSMOD(J).EQ.IDOUT(I).AND.BSSMOD(J).GT.0) THEN
              LOOP=LOOP+1
              SUMBR2=SUMBR2+BSSMOD(J)
              BR=SUMBR2/SUMBR
              CBR(LOOP)=BR
              MELEM(LOOP)=MSSMOD(J)
              DO 253 K=1,5
                MODE(K,LOOP)=JSSMOD(K,J)
                LMODE(K)=LABEL(MODE(K,LOOP))
253           CONTINUE
              IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
     $        BR,IRES,(MODE(K,LOOP),K=1,5)
            ENDIF
252       CONTINUE
250     CONTINUE
C
C          Top hadron decays
C
        DO 260 I=1,NTHAD
C          Check for modes
          SUMBR=0
          DO 261 J=1,NSSMOD
            IF(ISSMOD(J).EQ.6.AND.BSSMOD(J).GT.0) THEN
              SUMBR=SUMBR+BSSMOD(J)
            ENDIF
261       CONTINUE
C          If modes exist, add them -- conjugate for antimesons
          IF(SUMBR.LE.0) GO TO 260
          IRES=IABS(ITHAD(I))
          LRES=LABEL(IRES)
          CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
          LOOK(INDEX)=LOOP+1
          SUMBR2=0
          DO 262 J=1,NSSMOD
            IF(ISSMOD(J).EQ.6.AND.BSSMOD(J).GT.0) THEN
              LOOP=LOOP+1
              SUMBR2=SUMBR2+BSSMOD(J)
              BR=SUMBR2/SUMBR
              CBR(LOOP)=BR
              IF(IABS(JSSMOD(1,J)).LT.20.AND.IABS(JSSMOD(2,J)).LT.20
     $        .AND.IABS(JSSMOD(3,J)).LT.20.AND.IABS(JSSMOD(4,J)).LT.20
     $        .AND.IABS(JSSMOD(5,J)).LT.20) THEN
                MELEM(LOOP)=4
              ELSE
                MELEM(LOOP)=0
              ENDIF
              DO 263 K=1,5
                IF(ITHAD(I).GT.0) THEN
                  MODE(K,LOOP)=JSSMOD(K,J)
                ELSE
                  MODE(K,LOOP)=IDANTI(JSSMOD(K,J))
                ENDIF
                LMODE(K)=LABEL(MODE(K,LOOP))
263           CONTINUE
              IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
     $        BR,IRES,(MODE(K,LOOP),K=1,5)
            ENDIF
262       CONTINUE
260     CONTINUE
      ENDIF
C
C          Set forced decay modes.
C          LOOK(INDEX) = -IRES, where LOOK2(K,IRES) points to entries in
C          decay table for IDENT>0 and IDENT<0.
C          LOOKST(IRES) = standard LOOK value.
C
      IF(NFORCE.EQ.0) GO TO 400
C          Append each forced decay to table
      IRES=0
      DO 310 I=1,NFORCE
        IF(IFORCE(I).EQ.0) GO TO 310
        LOOP=LOOP+1
        IF(LOOP.GT.MXDKY) GO TO 9999
        CALL FLAVOR(IFORCE(I),IFL1,IFL2,IFL3,JSPIN,INDEX)
        IF(IFORCE(I).GT.0) THEN
          KTYPE=1
        ELSE
          KTYPE=2
        ENDIF
C
        IF(LOOK(INDEX).GE.0) THEN
          IRES=IRES+1
          IF(IRES.GT.MXFORC) GO TO 9998
          LOOKST(IRES)=LOOK(INDEX)
          LOOK2(KTYPE,IRES)=LOOP
          LOOK2(3-KTYPE,IRES)=LOOKST(IRES)
          LOOK(INDEX)=-IRES
        ELSE
          IRES2=-LOOK(INDEX)
          IF(IRES2.GT.MXFORC) GO TO 9998
          LOOK2(KTYPE,IRES2)=LOOP
        ENDIF
C          Set forced decay mode - conjugate if necessary
        IF(KTYPE.EQ.1) THEN
          DO 320 K=1,5
320       MODE(K,LOOP)=MFORCE(K,I)
        ELSE
          DO 330 K=1,5
330       MODE(K,LOOP)=IDANTI(MFORCE(K,I))
        ENDIF
        CBR(LOOP)=1.
C          Set matrix element flag
        CALL ORDER(IFORCE(I),MFORCE(1,I),IDUMMY,MEOUT)
        MELEM(LOOP)=MEOUT
        MEFORC(I)=MEOUT
310   CONTINUE
C
400   RETURN
C
C          Errors
C
9999  WRITE(ITLIS,3001) LOOP
3001  FORMAT(//' ***** ERROR IN SETDKY ... DECAY COUNTER LOOP = ',
     $I6,' *****')
      STOP 99
9998  WRITE(ITLIS,3002) IRES
3002  FORMAT(//' ***** ERROR IN SETDKY ... FORCE COUNTER IRES = ',
     $I6,' *****')
      STOP 99
      END
Commit	Line	Data
0795afa3	1	#include "isajet/pilot.h"
	2	SUBROUTINE SETDKY(LPRINT)
	3	C
	4	C Read in decay table from tape ITDKY and set up /DKYTAB/.
	5	C Then append forced decay modes and set LOOK to negative
	6	C number pointing to LOOK2, which points to table.
	7	C Forced decays for antiparticles are stored in conjugated
	8	C form so that DECAY can always conjugate them.
	9	C
	10	C Logical flag LPRINT controls printing of table.
	11	C
	12	C Ver 7.41: Check version of decay table. Also read matrix
	13	C element flags and save in MELEM:
	14	C MELEM=0: Phase space
	15	C MELEM=1: Dalitz decay
	16	C MELEM=2: omega/phi decay
	17	C MELEM=3: V-A
	18	C MELEM=4: V-A plus W propagator (for top)
	19	C MELEM=5: tau -> ell nu nu
	20	C MELEM=6: tau -> nu pi/K
	21	C MELEM=7: tau -> nu rho/a1
	22	C
	23	#if defined(CERNLIB_IMPNONE)
	24	IMPLICIT NONE
	25	#endif
	26	#include "isajet/itapes.inc"
	27	#include "isajet/force.inc"
	28	#include "isajet/dkytab.inc"
	29	#include "isajet/nodcay.inc"
	30	#include "isajet/ssmode.inc"
	31	#include "isajet/sstype.inc"
	32	#include "isajet/xmssm.inc"
	33	#include "isajet/keys.inc"
	34	C
	35	INTEGER IMODE(6),LOOP,IOLD,I,IRES,ITYPE,K,J,IPOINT
	36	INTEGER IFL1,IFL2,IFL3,JSPIN,INDEX,ID1,IDANTI,KTYPE,IRES2
	37	REAL BR
	38	CHARACTER*8 LABEL,LMODE(6),LRES
	39	CHARACTER*8 IBLANK,LREAD(10),IQUIT
	40	LOGICAL LPRINT
	41	INTEGER NOUT,NTHAD
	42	PARAMETER (NOUT=33)
	43	PARAMETER (NTHAD=12)
	44	INTEGER IDOUT(NOUT),ITHAD(NTHAD),IDUMMY(5),MEOUT
	45	REAL SUMBR,SUMBR2,SUMGAM
	46	CHARACTER*40,V,VOLD,VISAJE
	47	C
	48	DATA IDOUT/
	49	$IDTP,ISGL,ISUPL,ISDNL,ISSTL,ISCHL,ISBT1,ISTP1,ISUPR,ISDNR,
	50	$ISSTR,ISCHR,ISBT2,ISTP2,ISEL,ISMUL,ISTAU1,ISNEL,ISNML,ISNTL,
	51	$ISER,ISMUR,ISTAU2,ISZ1,ISZ2,ISZ3,ISZ4,ISW1,ISW2,
	52	$ISHL,ISHH,ISHA,ISHC/
	53	DATA IQUIT/'////'/,IBLANK/' '/
	54	DATA ITHAD/-160,-260,-360,
	55	$ 1160,1260,2260,2160,1360,2360,3160,3260,3360/
	56	C
	57	C Print header for table.
	58	C
	59	IF(LPRINT) WRITE(ITLIS,10)
	60	10 FORMAT('1',30('')/' ',28X,'*'/
	61	1' ',5X,'ISAJET DECAY TABLE',5X,''/
	62	2' ',28X,''/' ',30('*')//
	63	33X,'PART',16X,'DECAY MODE',16X,'CUM BR',10X,'IDENT',18X,
	64	4'DECAY IDENT'/)
65	C
66	C Initialize. LOOP is the decay mode counter.
67	C
68	LOOP=0
69	IOLD=0
70	DO 100 I=1,MXLOOK
71	LOOK(I)=0
72	100 CONTINUE
73	DO 110 I=1,MXFORC
74	LOOK2(1,I)=0
75	LOOK2(2,I)=0
76	110 CONTINUE
77	C
78	C Read in table, checking for valid version.
79	C
80	IF(NODCAY.OR.ITDKY.EQ.0) RETURN
81	REWIND ITDKY
82	C
83	VOLD=VISAJE()
84	READ(ITDKY,*) V
85	IF(V.NE.VOLD) THEN
86	WRITE(ITLIS,2000) V,VOLD
87	2000 FORMAT(//
88	$ ' ***WARNING: DECAY TABLE DOES NOT MATCH ISAJET VERSION'/
89	$ ' ***DECAY VERSION : ',A40/
90	$ ' ***PROGRAM VERSION: ',A40)
91	ENDIF
92	C
93	200 LOOP=LOOP+1
94	IF(LOOP.GT.MXDKY) GO TO 9999
95	220 DO 210 I=1,5
96	IMODE(I)=0
97	LMODE(I)=IBLANK
98	210 CONTINUE
99	READ(ITDKY,*) IRES,ITYPE,BR,IMODE
100	C
101	IF(IRES.NE.0) THEN
102	IF(NOPI0.AND.IRES.EQ.110) GO TO 220
103	IF(NOETA.AND.IRES.EQ.220) GO TO 220
104	IF(IRES.NE.IOLD) THEN
105	CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
106	LOOK(INDEX)=LOOP
107	ENDIF
108	IOLD=IRES
109	CBR(LOOP)=BR
110	MELEM(LOOP)=ITYPE
111	DO 240 I=1,5
112	MODE(I,LOOP)=IMODE(I)
113	IF(IMODE(I).NE.0) LMODE(I)=LABEL(IMODE(I))
114	240 CONTINUE
115	LRES=LABEL(IRES)
116	IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
117	1 BR,IRES,(IMODE(K),K=1,5)
118	20 FORMAT(3X,A5,4X,5(A5,2X),F8.5,10X,I5,4X,5(I7,2X))
119	GO TO 200
120	ENDIF
121	C
122	C Add HIGGS FOR WHIGGS
123	C
124	IF(KEYS(10).AND..NOT.GOMSSM) THEN
125	SUMGAM=0
126	SUMBR=0
127	DO 244 J=1,NSSMOD
128	IF(ISSMOD(J).EQ.81.AND.GSSMOD(J).GT.0) THEN
129	SUMGAM=SUMGAM+GSSMOD(J)
130	ENDIF
131	244 CONTINUE
132	DO 245 J=1,NSSMOD
133	IF(ISSMOD(J).EQ.81.AND.GSSMOD(J).GT.0) THEN
134	BSSMOD(J)=GSSMOD(J)/SUMGAM
135	ENDIF
136	245 CONTINUE
137	DO 246 J=1,NSSMOD
138	IF(ISSMOD(J).EQ.81.AND.BSSMOD(J).GT.0) THEN
139	SUMBR=SUMBR+BSSMOD(J)
140	ENDIF
141	246 CONTINUE
142	C If modes exist, add them
143	IF(SUMBR.LE.0) GO TO 249
144	IRES=81
145	LRES=LABEL(IRES)
146	CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
147	LOOK(INDEX)=LOOP+1
148	SUMBR2=0
149	DO 247 J=1,NSSMOD
150	IF(ISSMOD(J).EQ.81.AND.BSSMOD(J).GT.0) THEN
151	LOOP=LOOP+1
152	SUMBR2=SUMBR2+BSSMOD(J)
153	BR=SUMBR2/SUMBR
154	CBR(LOOP)=BR
155	MELEM(LOOP)=0
156	DO 248 K=1,5
157	MODE(K,LOOP)=JSSMOD(K,J)
158	LMODE(K)=LABEL(MODE(K,LOOP))
159	248 CONTINUE
160	IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
161	$ BR,IRES,(MODE(K,LOOP),K=1,5)
162	ENDIF
163	247 CONTINUE
164	249 CONTINUE
165	END IF
166	C
167	C Add MSSM decay modes if applicable, OR H_SM FOR WHIGGS
168	C
169	IF(GOMSSM) THEN
170	DO 250 I=1,NOUT
171	C Check for modes
172	SUMBR=0
173	DO 251 J=1,NSSMOD
174	IF(ISSMOD(J).EQ.IDOUT(I).AND.BSSMOD(J).GT.0) THEN
175	SUMBR=SUMBR+BSSMOD(J)
176	ENDIF
177	251 CONTINUE
178	C If modes exist, add them
179	IF(SUMBR.LE.0) GO TO 250
180	IRES=IDOUT(I)
181	LRES=LABEL(IRES)
182	CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
183	LOOK(INDEX)=LOOP+1
184	SUMBR2=0
185	DO 252 J=1,NSSMOD
186	IF(ISSMOD(J).EQ.IDOUT(I).AND.BSSMOD(J).GT.0) THEN
187	LOOP=LOOP+1
188	SUMBR2=SUMBR2+BSSMOD(J)
189	BR=SUMBR2/SUMBR
190	CBR(LOOP)=BR
191	MELEM(LOOP)=MSSMOD(J)
192	DO 253 K=1,5
193	MODE(K,LOOP)=JSSMOD(K,J)
194	LMODE(K)=LABEL(MODE(K,LOOP))
195	253 CONTINUE
196	IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
197	$ BR,IRES,(MODE(K,LOOP),K=1,5)
198	ENDIF
199	252 CONTINUE
200	250 CONTINUE
201	C
202	C Top hadron decays
203	C
204	DO 260 I=1,NTHAD
205	C Check for modes
206	SUMBR=0
207	DO 261 J=1,NSSMOD
208	IF(ISSMOD(J).EQ.6.AND.BSSMOD(J).GT.0) THEN
209	SUMBR=SUMBR+BSSMOD(J)
210	ENDIF
211	261 CONTINUE
212	C If modes exist, add them -- conjugate for antimesons
213	IF(SUMBR.LE.0) GO TO 260
214	IRES=IABS(ITHAD(I))
215	LRES=LABEL(IRES)
216	CALL FLAVOR(IRES,IFL1,IFL2,IFL3,JSPIN,INDEX)
217	LOOK(INDEX)=LOOP+1
218	SUMBR2=0
219	DO 262 J=1,NSSMOD
220	IF(ISSMOD(J).EQ.6.AND.BSSMOD(J).GT.0) THEN
221	LOOP=LOOP+1
222	SUMBR2=SUMBR2+BSSMOD(J)
223	BR=SUMBR2/SUMBR
224	CBR(LOOP)=BR
225	IF(IABS(JSSMOD(1,J)).LT.20.AND.IABS(JSSMOD(2,J)).LT.20
226	$ .AND.IABS(JSSMOD(3,J)).LT.20.AND.IABS(JSSMOD(4,J)).LT.20
227	$ .AND.IABS(JSSMOD(5,J)).LT.20) THEN
228	MELEM(LOOP)=4
229	ELSE
230	MELEM(LOOP)=0
231	ENDIF
232	DO 263 K=1,5
233	IF(ITHAD(I).GT.0) THEN
234	MODE(K,LOOP)=JSSMOD(K,J)
235	ELSE
236	MODE(K,LOOP)=IDANTI(JSSMOD(K,J))
237	ENDIF
238	LMODE(K)=LABEL(MODE(K,LOOP))
239	263 CONTINUE
240	IF(LPRINT) WRITE(ITLIS,20) LRES,(LMODE(K),K=1,5),
241	$ BR,IRES,(MODE(K,LOOP),K=1,5)
242	ENDIF
243	262 CONTINUE
244	260 CONTINUE
245	ENDIF
246	C
247	C Set forced decay modes.
248	C LOOK(INDEX) = -IRES, where LOOK2(K,IRES) points to entries in
249	C decay table for IDENT>0 and IDENT<0.
250	C LOOKST(IRES) = standard LOOK value.
251	C
252	IF(NFORCE.EQ.0) GO TO 400
253	C Append each forced decay to table
254	IRES=0
255	DO 310 I=1,NFORCE
256	IF(IFORCE(I).EQ.0) GO TO 310
257	LOOP=LOOP+1
258	IF(LOOP.GT.MXDKY) GO TO 9999
259	CALL FLAVOR(IFORCE(I),IFL1,IFL2,IFL3,JSPIN,INDEX)
260	IF(IFORCE(I).GT.0) THEN
261	KTYPE=1
262	ELSE
263	KTYPE=2
264	ENDIF
265	C
266	IF(LOOK(INDEX).GE.0) THEN
267	IRES=IRES+1
268	IF(IRES.GT.MXFORC) GO TO 9998
269	LOOKST(IRES)=LOOK(INDEX)
270	LOOK2(KTYPE,IRES)=LOOP
271	LOOK2(3-KTYPE,IRES)=LOOKST(IRES)
272	LOOK(INDEX)=-IRES
273	ELSE
274	IRES2=-LOOK(INDEX)
275	IF(IRES2.GT.MXFORC) GO TO 9998
276	LOOK2(KTYPE,IRES2)=LOOP
277	ENDIF
278	C Set forced decay mode - conjugate if necessary
279	IF(KTYPE.EQ.1) THEN
280	DO 320 K=1,5
281	320 MODE(K,LOOP)=MFORCE(K,I)
282	ELSE
283	DO 330 K=1,5
284	330 MODE(K,LOOP)=IDANTI(MFORCE(K,I))
285	ENDIF
286	CBR(LOOP)=1.
287	C Set matrix element flag
288	CALL ORDER(IFORCE(I),MFORCE(1,I),IDUMMY,MEOUT)
289	MELEM(LOOP)=MEOUT
290	MEFORC(I)=MEOUT
291	310 CONTINUE
292	C
293	400 RETURN
294	C
295	C Errors
296	C
297	9999 WRITE(ITLIS,3001) LOOP
298	3001 FORMAT(//' ***** ERROR IN SETDKY ... DECAY COUNTER LOOP = ',
299	$I6,' *****')
300	STOP 99
301	9998 WRITE(ITLIS,3002) IRES
302	3002 FORMAT(//' ***** ERROR IN SETDKY ... FORCE COUNTER IRES = ',
303	$I6,' *****')
304	STOP 99
305	END