Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
pysgtc.f
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file pysgtc.f
1 
2 C*********************************************************************
3 
4 C...PYSGTC
5 C...Subprocess cross sections for Technicolor processes.
6 C...Auxiliary to PYSIGH.
7 
8  SUBROUTINE pysgtc(NCHN,SIGS)
9 
10 C...Double precision and integer declarations
11  IMPLICIT DOUBLE PRECISION(a-h, o-z)
12  IMPLICIT INTEGER(i-n)
13  INTEGER pyk,pychge,pycomp
14 C...Parameter statement to help give large particle numbers.
15  parameter(ksusy1=1000000,ksusy2=2000000,ktechn=3000000,
16  &kexcit=4000000,kdimen=5000000)
17 C...Commonblocks
18  common/pydat1/mstu(200),paru(200),mstj(200),parj(200)
19  common/pydat2/kchg(500,4),pmas(500,4),parf(2000),vckm(4,4)
20  common/pydat3/mdcy(500,3),mdme(8000,2),brat(8000),kfdp(8000,5)
21  common/pypars/mstp(200),parp(200),msti(200),pari(200)
22  common/pyint1/mint(400),vint(400)
23  common/pyint2/iset(500),kfpr(500,2),coef(500,20),icol(40,4,2)
24  common/pyint3/xsfx(2,-40:40),isig(1000,3),sigh(1000)
25  common/pyint4/mwid(500),wids(500,5)
26  common/pytcsm/itcm(0:99),rtcm(0:99)
27  common/pysgcm/isub,isubsv,mmin1,mmax1,mmin2,mmax2,mmina,mmaxa,
28  &kfac(2,-40:40),comfac,fack,faca,sh,th,uh,sh2,th2,uh2,sqm3,sqm4,
29  &shr,sqpth,taup,be34,cth,x(2),sqmz,sqmw,gmmz,gmmw,
30  &aem,as,xw,xw1,xwc,xwv,poll,polr,polll,polrr
31  SAVE /pydat1/,/pydat2/,/pydat3/,/pypars/,/pyint1/,/pyint2/,
32  &/pyint3/,/pyint4/,/pytcsm/,/pysgcm/
33 C...Local arrays and complex variables
34  dimension wdtp(0:400),wdte(0:400,0:5)
35  COMPLEX*16 ssmz,ssmr,ssmo,detd,f2l,f2r,darho,dzrho,daome,dzome
36  COMPLEX*16 daa,dzz,daz,dww,dwrho
37  COMPLEX*16 ztc(6,6),ytc(6,6),dggs,dggt,dggu,dgvs,dgvt,dgvu
38  COMPLEX*16 dqqs,dqqt,dqqu,dqts,dqgs,dtgs
39  COMPLEX*16 dvvs,dvvt,dvvu
40  INTEGER indx(6)
41 
42 C...Combinations of weak mixing angle.
43  tanw=sqrt(xw/xw1)
44  ct2w=(1d0-2d0*xw)/(2d0*xw/tanw)
45 
46 C...Convert almost equivalent technicolor processes into
47 C...a few basic processes, and set distinguishing parameters.
48  IF(isub.GE.361.AND.isub.LE.379) THEN
49  sqtv=rtcm(12)**2
50  sqta=rtcm(13)**2
51  sn2w=2d0*sqrt(paru(102)*(1d0-paru(102)))
52  cs2w=1d0-2d0*paru(102)
53  tanw=sqrt(paru(102)/(1d0-paru(102)))
54  ct2w=cs2w/sn2w
55  csxi=cos(asin(rtcm(3)))
56  csxip=cos(asin(rtcm(4)))
57  qupd=2d0*rtcm(2)-1d0
58  q2ud=rtcm(2)**2+(rtcm(2)-1d0)**2
59 C... rho_tc0 -> W_L W_L
60  IF(isub.EQ.361) THEN
61  kfa=24
62  kfb=24
63  cab2=rtcm(3)**4
64 C... rho_tc0 -> W_L pi_tc-
65  ELSEIF(isub.EQ.362) THEN
66  kfa=24
67  kfb=ktechn+211
68  isub=361
69  cab2=rtcm(3)**2*(1d0-rtcm(3)**2)
70 C... pi_tc pi_tc
71  ELSEIF(isub.EQ.363) THEN
72  kfa=ktechn+211
73  kfb=ktechn+211
74  isub=361
75  cab2=(1d0-rtcm(3)**2)**2
76 C... rho_tc0/omega_tc -> gamma pi_tc
77  ELSEIF(isub.EQ.364) THEN
78  kfa=22
79  kfb=ktechn+111
80  vogp=csxi/rtcm(12)
81 C..........!!!
82  vrgp=vogp*qupd
83  aogp=0d0
84  argp=0d0
85  vagp=2d0*qupd*csxi
86  vzgp=qupd*csxi*(1d0-4d0*paru(102))/sn2w
87 C... gamma pi_tc'
88  ELSEIF(isub.EQ.365) THEN
89  kfa=22
90  kfb=ktechn+221
91  isub=364
92  vrgp=csxip/rtcm(12)
93 C..........!!!!
94  vogp=vrgp*qupd
95  aogp=0d0
96  argp=0d0
97  vagp=2d0*q2ud*csxip
98  vzgp=csxip/sn2w*(1d0-4d0*paru(102)*q2ud)
99 C... Z pi_tc
100  ELSEIF(isub.EQ.366) THEN
101  kfa=23
102  kfb=ktechn+111
103  isub=364
104  vogp=csxi*ct2w/rtcm(12)
105  vrgp=-qupd*csxi*tanw/rtcm(12)
106  aogp=0d0
107  argp=0d0
108  vagp=qupd*csxi*(1d0-4d0*paru(102))/sn2w
109  vzgp=-qupd*csxi*cs2w/(1d0-paru(102))
110 C... Z pi_tc'
111  ELSEIF(isub.EQ.367) THEN
112  kfa=23
113  kfb=ktechn+221
114  isub=364
115  vrgp=csxip*ct2w/rtcm(12)
116  vogp=-qupd*csxip*tanw/rtcm(12)
117  aogp=0d0
118  argp=0d0
119  vagp=csxip*(1d0-4d0*q2ud*paru(102))/sn2w
120  vzgp=2d0*csxip*(cs2w+4d0*q2ud*paru(102)**2)/sn2w**2
121 C... W_T pi_tc
122  ELSEIF(isub.EQ.368) THEN
123  kfa=24
124  kfb=ktechn+211
125  isub=364
126  vogp=csxi/(2d0*sqrt(paru(102)))/rtcm(12)
127  vrgp=0d0
128  aogp=0d0
129 C..........!!!!
130  argp=-csxi/(2d0*sqrt(paru(102)))/rtcm(13)
131  vagp=qupd*csxi/(2d0*sqrt(paru(102)))
132  vzgp=-qupd*csxi/(2d0*sqrt(1d0-paru(102)))
133 C... rho_tc+ -> W_L Z_L
134  ELSEIF(isub.EQ.370) THEN
135  kfa=24
136  kfb=23
137  cab2=rtcm(3)**4
138 C... W_L pi_tc0
139  ELSEIF(isub.EQ.371) THEN
140  kfa=24
141  kfb=ktechn+111
142  isub=370
143  cab2=rtcm(3)**2*(1d0-rtcm(3)**2)
144 C... Z_L pi_tc+
145  ELSEIF(isub.EQ.372) THEN
146  kfa=ktechn+211
147  kfb=23
148  isub=370
149  cab2=rtcm(3)**2*(1d0-rtcm(3)**2)
150 C... pi_tc+ pi_tc0
151  ELSEIF(isub.EQ.373) THEN
152  kfa=ktechn+211
153  kfb=ktechn+111
154  isub=370
155  cab2=(1d0-rtcm(3)**2)**2
156 C... gamma pi_tc+
157  ELSEIF(isub.EQ.374) THEN
158  kfa=ktechn+211
159  kfb=22
160  vrgp=qupd*csxi
161  argp=0d0
162  vwgp=qupd*csxi/(2d0*sqrt(paru(102)))
163 C... Z_T pi_tc+
164  ELSEIF(isub.EQ.375) THEN
165  kfa=ktechn+211
166  kfb=23
167  isub=374
168  vrgp=-qupd*csxi*tanw
169  argp=csxi/(2d0*sqrt(paru(102)*(1d0-paru(102))))
170  vwgp=-qupd*csxi/(2d0*sqrt(1d0-paru(102)))
171 C... W_T pi_tc0
172  ELSEIF(isub.EQ.376) THEN
173  kfa=24
174  kfb=ktechn+111
175  isub=374
176  vrgp=0d0
177  argp=-csxi/(2d0*sqrt(paru(102)))
178  vwgp=0d0
179 C... W_T pi_tc0'
180  ELSEIF(isub.EQ.377) THEN
181  kfa=24
182  kfb=ktechn+221
183  isub=374
184  argp=0d0
185  vrgp=csxip/(2d0*sqrt(paru(102)))
186  vwgp=csxip/(2d0*paru(102))
187  ENDIF
188  ENDIF
189 
190 C...QCD 2 -> 2 processes: corrections from virtual technicolor exchange.
191  IF(isub.GE.381.AND.isub.LE.388) THEN
192  IF(itcm(5).LE.4) THEN
193  sqdqqs=1d0/sh2
194  sqdqqt=1d0/th2
195  sqdqqu=1d0/uh2
196  sqdggs=sqdqqs
197  sqdggt=sqdqqt
198  sqdggu=sqdqqu
199  redggs=1d0/sh
200  redggt=1d0/th
201  redggu=1d0/uh
202  redgtu=1d0/uh/th
203  redgsu=1d0/sh/uh
204  redgst=1d0/sh/th
205  redqst=1d0/sh/th
206  redqtu=1d0/uh/th
207  sqdlgs=0d0
208  sqdlgt=0d0
209  sqdqts=sqdqqs
210  ELSEIF(itcm(5).EQ.5) THEN
211  tant3=rtcm(21)
212  IF(itcm(2).EQ.0) THEN
213  imdl=1
214  ELSE
215  imdl=2
216  ENDIF
217  alprht=2.91d0*(3d0/itcm(1))
218  sin2t=2d0*tant3/(tant3**2+1d0)
219  sint3=tant3/sqrt(tant3**2+1d0)
220  xig=sqrt(pyalps(sh)/alprht)
221  x12=(rtcm(29)*sqrt(1d0-rtcm(29)**2)*cos(rtcm(30))+
222  & rtcm(31)*sqrt(1d0-rtcm(31)**2)*cos(rtcm(32)))/sqrt(2d0)/sin2t
223  x21=(rtcm(29)*sqrt(1d0-rtcm(29)**2)*sin(rtcm(30))+
224  & rtcm(31)*sqrt(1d0-rtcm(31)**2)*sin(rtcm(32)))/sqrt(2d0)/sin2t
225  x11=(.25d0*(rtcm(29)**2+rtcm(31)**2+2d0)-
226  & sint3**2)*2d0/sin2t
227  x22=(.25d0*(2d0-rtcm(29)**2-rtcm(31)**2)-
228  & sint3**2)*2d0/sin2t
229 
230  sm1122=.5d0*(2d0-rtcm(29)**2-rtcm(31)**2)*rtcm(28)**2
231  sm1112=x12*rtcm(28)**2*sin2t
232  sm1121=-x21*rtcm(28)**2*sin2t
233  sm2212=-sm1112
234  sm2221=-sm1121
235  sm1221=-.5d0*((1d0-rtcm(29)**2)*sin(2d0*rtcm(30))+
236  & (1d0-rtcm(31)**2)*sin(2d0*rtcm(32)))*rtcm(28)**2
237 
238 C.........SH LOOP
239  ztc(1,1)=dcmplx(sh,0d0)
240  CALL pywidt(3100021,sh,wdtp,wdte)
241  IF(wdtp(0).GT.rtcm(33)*shr) wdtp(0)=rtcm(33)*shr
242  ztc(2,2)=dcmplx(sh-pmas(pycomp(3100021),1)**2,-shr*wdtp(0))
243  CALL pywidt(3100113,sh,wdtp,wdte)
244  ztc(3,3)=dcmplx(sh-pmas(pycomp(3100113),1)**2,-shr*wdtp(0))
245  CALL pywidt(3400113,sh,wdtp,wdte)
246  ztc(4,4)=dcmplx(sh-pmas(pycomp(3400113),1)**2,-shr*wdtp(0))
247  CALL pywidt(3200113,sh,wdtp,wdte)
248  ztc(5,5)=dcmplx(sh-pmas(pycomp(3200113),1)**2,-shr*wdtp(0))
249  CALL pywidt(3300113,sh,wdtp,wdte)
250  ztc(6,6)=dcmplx(sh-pmas(pycomp(3300113),1)**2,-shr*wdtp(0))
251  ztc(1,2)=(0d0,0d0)
252  ztc(1,3)=dcmplx(sh*xig,0d0)
253  ztc(1,4)=ztc(1,3)
254  ztc(1,5)=ztc(1,2)
255  ztc(1,6)=ztc(1,2)
256  ztc(2,3)=dcmplx(sh*xig*x11,0d0)
257  ztc(2,4)=dcmplx(sh*xig*x22,0d0)
258  ztc(2,5)=dcmplx(sh*xig*x12,0d0)
259  ztc(2,6)=dcmplx(sh*xig*x21,0d0)
260  ztc(3,4)=-sm1122
261  ztc(3,5)=-sm1112
262  ztc(3,6)=-sm1121
263  ztc(4,5)=-sm2212
264  ztc(4,6)=-sm2221
265  ztc(5,6)=-sm1221
266 
267  DO 110 i=1,5
268  DO 100 j=i+1,6
269  ztc(j,i)=ztc(i,j)
270  100 CONTINUE
271  110 CONTINUE
272  CALL pyldcm(ztc,6,6,indx,d)
273  DO 130 i=1,6
274  DO 120 j=1,6
275  ytc(i,j)=(0d0,0d0)
276  IF(i.EQ.j) ytc(i,j)=(1d0,0d0)
277  120 CONTINUE
278  130 CONTINUE
279 
280  DO 140 i=1,6
281  CALL pybksb(ztc,6,6,indx,ytc(1,i))
282  140 CONTINUE
283  dggs=ytc(1,1)
284  dvvs=ytc(2,2)
285  dgvs=ytc(1,2)
286 
287  xig=sqrt(pyalps(-th)/alprht)
288 C.........TH LOOP
289  ztc(1,1)=dcmplx(th)
290  ztc(2,2)=dcmplx(th-pmas(pycomp(3100021),1)**2)
291  ztc(3,3)=dcmplx(th-pmas(pycomp(3100113),1)**2)
292  ztc(4,4)=dcmplx(th-pmas(pycomp(3400113),1)**2)
293  ztc(5,5)=dcmplx(th-pmas(pycomp(3200113),1)**2)
294  ztc(6,6)=dcmplx(th-pmas(pycomp(3300113),1)**2)
295  ztc(1,2)=(0d0,0d0)
296  ztc(1,3)=dcmplx(th*xig,0d0)
297  ztc(1,4)=ztc(1,3)
298  ztc(1,5)=ztc(1,2)
299  ztc(1,6)=ztc(1,2)
300  ztc(2,3)=dcmplx(th*xig*x11,0d0)
301  ztc(2,4)=dcmplx(th*xig*x22,0d0)
302  ztc(2,5)=dcmplx(th*xig*x12,0d0)
303  ztc(2,6)=dcmplx(th*xig*x21,0d0)
304  ztc(3,4)=-sm1122
305  ztc(3,5)=-sm1112
306  ztc(3,6)=-sm1121
307  ztc(4,5)=-sm2212
308  ztc(4,6)=-sm2221
309  ztc(5,6)=-sm1221
310  DO 160 i=1,5
311  DO 150 j=i+1,6
312  ztc(j,i)=ztc(i,j)
313  150 CONTINUE
314  160 CONTINUE
315  CALL pyldcm(ztc,6,6,indx,d)
316  DO 180 i=1,6
317  DO 170 j=1,6
318  ytc(i,j)=(0d0,0d0)
319  IF(i.EQ.j) ytc(i,j)=(1d0,0d0)
320  170 CONTINUE
321  180 CONTINUE
322  DO 190 i=1,6
323  CALL pybksb(ztc,6,6,indx,ytc(1,i))
324  190 CONTINUE
325  dggt=ytc(1,1)
326  dvvt=ytc(2,2)
327  dgvt=ytc(1,2)
328 
329  xig=sqrt(pyalps(-uh)/alprht)
330 C.........UH LOOP
331  ztc(1,1)=dcmplx(uh,0d0)
332  ztc(2,2)=dcmplx(uh-pmas(pycomp(3100021),1)**2)
333  ztc(3,3)=dcmplx(uh-pmas(pycomp(3100113),1)**2)
334  ztc(4,4)=dcmplx(uh-pmas(pycomp(3400113),1)**2)
335  ztc(5,5)=dcmplx(uh-pmas(pycomp(3200113),1)**2)
336  ztc(6,6)=dcmplx(uh-pmas(pycomp(3300113),1)**2)
337  ztc(1,2)=(0d0,0d0)
338  ztc(1,3)=dcmplx(uh*xig,0d0)
339  ztc(1,4)=ztc(1,3)
340  ztc(1,5)=ztc(1,2)
341  ztc(1,6)=ztc(1,2)
342  ztc(2,3)=dcmplx(uh*xig*x11,0d0)
343  ztc(2,4)=dcmplx(uh*xig*x22,0d0)
344  ztc(2,5)=dcmplx(uh*xig*x12,0d0)
345  ztc(2,6)=dcmplx(uh*xig*x21,0d0)
346  ztc(3,4)=-sm1122
347  ztc(3,5)=-sm1112
348  ztc(3,6)=-sm1121
349  ztc(4,5)=-sm2212
350  ztc(4,6)=-sm2221
351  ztc(5,6)=-sm1221
352  DO 210 i=1,5
353  DO 200 j=i+1,6
354  ztc(j,i)=ztc(i,j)
355  200 CONTINUE
356  210 CONTINUE
357  CALL pyldcm(ztc,6,6,indx,d)
358  DO 230 i=1,6
359  DO 220 j=1,6
360  ytc(i,j)=(0d0,0d0)
361  IF(i.EQ.j) ytc(i,j)=(1d0,0d0)
362  220 CONTINUE
363  230 CONTINUE
364  DO 240 i=1,6
365  CALL pybksb(ztc,6,6,indx,ytc(1,i))
366  240 CONTINUE
367  dggu=ytc(1,1)
368  dvvu=ytc(2,2)
369  dgvu=ytc(1,2)
370 
371  IF(imdl.EQ.1) THEN
372  dqqs=dggs+dvvs*dcmplx(tant3**2)-dgvs*dcmplx(2d0*tant3)
373  dqqt=dggt+dvvt*dcmplx(tant3**2)-dgvt*dcmplx(2d0*tant3)
374  dqqu=dggu+dvvu*dcmplx(tant3**2)-dgvu*dcmplx(2d0*tant3)
375  dqts=dggs-dvvs-dgvs*dcmplx(tant3-1d0/tant3)
376  dqgs=dggs-dgvs*dcmplx(tant3)
377  dtgs=dggs+dgvs*dcmplx(1d0/tant3)
378  ELSE
379  dqqs=dggs+dvvs*dcmplx(1d0/tant3**2)+dgvs*dcmplx(2d0/tant3)
380  dqqt=dggt+dvvt*dcmplx(1d0/tant3**2)+dgvt*dcmplx(2d0/tant3)
381  dqqu=dggu+dvvu*dcmplx(1d0/tant3**2)+dgvu*dcmplx(2d0/tant3)
382  dqts=dggs+dvvs*dcmplx(1d0/tant3**2)+dgvs*dcmplx(2d0/tant3)
383  dqgs=dggs+dgvs*dcmplx(1d0/tant3)
384  dtgs=dggs+dgvs*dcmplx(1d0/tant3)
385  ENDIF
386 
387  sqdqts=abs(dqts)**2
388  sqdqqs=abs(dqqs)**2
389  sqdqqt=abs(dqqt)**2
390  sqdqqu=abs(dqqu)**2
391  sqdlgs=abs(dcmplx(sh)*dqgs-dcmplx(1d0))**2
392  redlgs=dble(dqgs)
393  sqdhgs=abs(dcmplx(sh)*dtgs-dcmplx(1d0))**2
394  redhgs=dble(dtgs)
395  sqdlgt=abs(dcmplx(th)*dggt-dcmplx(1d0))**2
396 
397  sqdggs=abs(dggs)**2
398  sqdggt=abs(dggt)**2
399  sqdggu=abs(dggu)**2
400  redggs=dble(dggs)
401  redggt=dble(dggt)
402  redggu=dble(dggu)
403  redgtu=dble(dggu*dconjg(dggt))
404  redgsu=dble(dggu*dconjg(dggs))
405  redgst=dble(dggs*dconjg(dggt))
406  redqst=dble(dqqs*dconjg(dqqt))
407  redqtu=dble(dqqt*dconjg(dqqu))
408  ENDIF
409  ENDIF
410 
411 
412 C...Differential cross section expressions.
413 
414  IF(isub.LE.190) THEN
415  IF(isub.EQ.149) THEN
416 C...g + g -> eta_tc
417  kctc=pycomp(ktechn+331)
418  CALL pywidt(ktechn+331,sh,wdtp,wdte)
419  hs=shr*wdtp(0)
420  facbw=comfac*0.5d0/((sh-pmas(kctc,1)**2)**2+hs**2)
421  IF(abs(shr-pmas(kctc,1)).GT.parp(48)*pmas(kctc,2)) facbw=0d0
422  hp=sh
423  IF(kfac(1,21)*kfac(2,21).EQ.0) goto 250
424  hi=hp*wdtp(3)
425  hf=shr*(wdte(0,1)+wdte(0,2)+wdte(0,4))
426  nchn=nchn+1
427  isig(nchn,1)=21
428  isig(nchn,2)=21
429  isig(nchn,3)=1
430  sigh(nchn)=hi*facbw*hf
431  250 CONTINUE
432 
433  ELSEIF(isub.EQ.165) THEN
434 C...q + qbar -> l+ + l- (including contact term for compositeness)
435  zratr=xwc*sh*(sh-sqmz)/((sh-sqmz)**2+gmmz**2)
436  zrati=xwc*sh*gmmz/((sh-sqmz)**2+gmmz**2)
437  kff=iabs(kfpr(isub,1))
438  ef=kchg(kff,1)/3d0
439  af=sign(1d0,ef+0.1d0)
440  vf=af-4d0*ef*xwv
441  valf=vf+af
442  varf=vf-af
443  fcof=1d0
444  IF(kff.LE.10) fcof=3d0
445  wid2=1d0
446  IF(kff.EQ.6) wid2=wids(6,1)
447  IF(kff.EQ.7.OR.kff.EQ.8) wid2=wids(kff,1)
448  IF(kff.EQ.17.OR.kff.EQ.18) wid2=wids(kff,1)
449  DO 260 i=mmina,mmaxa
450  IF(i.EQ.0.OR.kfac(1,i)*kfac(2,-i).EQ.0) goto 260
451  ei=kchg(iabs(i),1)/3d0
452  ai=sign(1d0,ei+0.1d0)
453  vi=ai-4d0*ei*xwv
454  vali=vi+ai
455  vari=vi-ai
456  fcoi=1d0
457  IF(iabs(i).LE.10) fcoi=faca/3d0
458  IF((itcm(5).EQ.1.AND.iabs(i).LE.2).OR.itcm(5).EQ.2) THEN
459  fgza=(ei*ef+vali*valf*zratr+rtcm(42)*sh/
460  & (aem*rtcm(41)**2))**2+(vali*valf*zrati)**2+
461  & (ei*ef+vari*varf*zratr)**2+(vari*varf*zrati)**2
462  ELSE
463  fgza=(ei*ef+vali*valf*zratr)**2+(vali*valf*zrati)**2+
464  & (ei*ef+vari*varf*zratr)**2+(vari*varf*zrati)**2
465  ENDIF
466  fgzb=(ei*ef+vali*varf*zratr)**2+(vali*varf*zrati)**2+
467  & (ei*ef+vari*valf*zratr)**2+(vari*valf*zrati)**2
468  fgzab=aem**2*(fgza*uh2/sh2+fgzb*th2/sh2)
469  IF((itcm(5).EQ.3.AND.iabs(i).EQ.2).OR.(itcm(5).EQ.4.AND.
470  & mod(iabs(i),2).EQ.0)) fgzab=fgzab+sh2/(2d0*rtcm(41)**4)
471  nchn=nchn+1
472  isig(nchn,1)=i
473  isig(nchn,2)=-i
474  isig(nchn,3)=1
475  sigh(nchn)=comfac*fcoi*fcof*fgzab*wid2
476  260 CONTINUE
477 
478  ELSEIF(isub.EQ.166) THEN
479 C...q + q'bar -> l + nu_l (including contact term for compositeness)
480  wfac=(1d0/4d0)*(aem/xw)**2*uh2/((sh-sqmw)**2+gmmw**2)
481  wcifac=wfac+sh2/(4d0*rtcm(41)**4)
482  kff=iabs(kfpr(isub,1))
483  fcof=1d0
484  IF(kff.LE.10) fcof=3d0
485  DO 280 i=mmin1,mmax1
486  IF(i.EQ.0.OR.kfac(1,i).EQ.0) goto 280
487  ia=iabs(i)
488  DO 270 j=mmin2,mmax2
489  IF(j.EQ.0.OR.kfac(2,j).EQ.0) goto 270
490  ja=iabs(j)
491  IF(i*j.GT.0.OR.mod(ia+ja,2).EQ.0) goto 270
492  IF((ia.LE.10.AND.ja.GT.10).OR.(ia.GT.10.AND.ja.LE.10))
493  & goto 270
494  fcoi=1d0
495  IF(ia.LE.10) fcoi=vckm((ia+1)/2,(ja+1)/2)*faca/3d0
496  wid2=1d0
497  IF((i.GT.0.AND.mod(i,2).EQ.0).OR.(j.GT.0.AND.
498  & mod(j,2).EQ.0)) THEN
499  IF(kff.EQ.5) wid2=wids(6,2)
500  IF(kff.EQ.7) wid2=wids(8,2)*wids(7,3)
501  IF(kff.EQ.17) wid2=wids(18,2)*wids(17,3)
502  ELSE
503  IF(kff.EQ.5) wid2=wids(6,3)
504  IF(kff.EQ.7) wid2=wids(8,3)*wids(7,2)
505  IF(kff.EQ.17) wid2=wids(18,3)*wids(17,2)
506  ENDIF
507  nchn=nchn+1
508  isig(nchn,1)=i
509  isig(nchn,2)=j
510  isig(nchn,3)=1
511  sigh(nchn)=comfac*fcoi*fcof*wfac*wid2
512  IF((itcm(5).EQ.3.AND.ia.LE.2.AND.ja.LE.2).OR.itcm(5).EQ.4)
513  & sigh(nchn)=comfac*fcoi*fcof*wcifac*wid2
514  270 CONTINUE
515  280 CONTINUE
516  ENDIF
517 
518  ELSEIF(isub.LE.200) THEN
519  IF(isub.EQ.191) THEN
520 C...q + qbar -> rho_tc0.
521  kctc=pycomp(ktechn+113)
522  sqmrht=pmas(kctc,1)**2
523  CALL pywidt(ktechn+113,sh,wdtp,wdte)
524  hs=shr*wdtp(0)
525  facbw=12d0*comfac/((sh-sqmrht)**2+hs**2)
526  IF(abs(shr-pmas(kctc,1)).GT.parp(48)*pmas(kctc,2)) facbw=0d0
527  hf=shr*(wdte(0,1)+wdte(0,2)+wdte(0,4))
528  alprht=2.91d0*(3d0/itcm(1))
529  hp=(1d0/6d0)*(aem**2/alprht)*(sqmrht**2/sh)
530  xwrht=(1d0-2d0*xw)/(4d0*xw*(1d0-xw))
531  bwzr=xwrht*sh*(sh-sqmz)/((sh-sqmz)**2+gmmz**2)
532  bwzi=xwrht*sh*gmmz/((sh-sqmz)**2+gmmz**2)
533  DO 290 i=mmina,mmaxa
534  IF(i.EQ.0.OR.kfac(1,i)*kfac(2,-i).EQ.0) goto 290
535  ia=iabs(i)
536  ei=kchg(iabs(i),1)/3d0
537  ai=sign(1d0,ei+0.1d0)
538  vi=ai-4d0*ei*xwv
539  vali=0.5d0*(vi+ai)
540  vari=0.5d0*(vi-ai)
541  hi=hp*((ei+vali*bwzr)**2+(vali*bwzi)**2+
542  & (ei+vari*bwzr)**2+(vari*bwzi)**2)
543  IF(ia.LE.10) hi=hi*faca/3d0
544  nchn=nchn+1
545  isig(nchn,1)=i
546  isig(nchn,2)=-i
547  isig(nchn,3)=1
548  sigh(nchn)=hi*facbw*hf
549  290 CONTINUE
550 
551  ELSEIF(isub.EQ.192) THEN
552 C...q + qbar' -> rho_tc+/-.
553  kctc=pycomp(ktechn+213)
554  sqmrht=pmas(kctc,1)**2
555  CALL pywidt(ktechn+213,sh,wdtp,wdte)
556  hs=shr*wdtp(0)
557  facbw=12d0*comfac/((sh-sqmrht)**2+hs**2)
558  IF(abs(shr-pmas(kctc,1)).GT.parp(48)*pmas(kctc,2)) facbw=0d0
559  alprht=2.91d0*(3d0/itcm(1))
560  hp=(1d0/6d0)*(aem**2/alprht)*(sqmrht**2/sh)*
561  & (0.25d0/xw**2)*sh**2/((sh-sqmw)**2+gmmw**2)
562  DO 310 i=mmin1,mmax1
563  IF(i.EQ.0.OR.kfac(1,i).EQ.0) goto 310
564  ia=iabs(i)
565  DO 300 j=mmin2,mmax2
566  IF(j.EQ.0.OR.kfac(2,j).EQ.0) goto 300
567  ja=iabs(j)
568  IF(i*j.GT.0.OR.mod(ia+ja,2).EQ.0) goto 300
569  IF((ia.LE.10.AND.ja.GT.10).OR.(ia.GT.10.AND.ja.LE.10))
570  & goto 300
571  kchr=(kchg(ia,1)*isign(1,i)+kchg(ja,1)*isign(1,j))/3
572  hf=shr*(wdte(0,1)+wdte(0,(5-kchr)/2)+wdte(0,4))
573  hi=hp
574  IF(ia.LE.10) hi=hi*vckm((ia+1)/2,(ja+1)/2)*faca/3d0
575  nchn=nchn+1
576  isig(nchn,1)=i
577  isig(nchn,2)=j
578  isig(nchn,3)=1
579  sigh(nchn)=hi*facbw*hf
580  300 CONTINUE
581  310 CONTINUE
582 
583  ELSEIF(isub.EQ.193) THEN
584 C...q + qbar -> omega_tc0.
585  kctc=pycomp(ktechn+223)
586  sqmomt=pmas(kctc,1)**2
587  CALL pywidt(ktechn+223,sh,wdtp,wdte)
588  hs=shr*wdtp(0)
589  facbw=12d0*comfac/((sh-sqmomt)**2+hs**2)
590  IF(abs(shr-pmas(kctc,1)).GT.parp(48)*pmas(kctc,2)) facbw=0d0
591  hf=shr*(wdte(0,1)+wdte(0,2)+wdte(0,4))
592  alprht=2.91d0*(3d0/itcm(1))
593  hp=(1d0/6d0)*(aem**2/alprht)*(sqmomt**2/sh)*
594  & (2d0*rtcm(2)-1d0)**2
595  bwzr=(0.5d0/(1d0-xw))*sh*(sh-sqmz)/((sh-sqmz)**2+gmmz**2)
596  bwzi=(0.5d0/(1d0-xw))*sh*gmmz/((sh-sqmz)**2+gmmz**2)
597  DO 320 i=mmina,mmaxa
598  IF(i.EQ.0.OR.kfac(1,i)*kfac(2,-i).EQ.0) goto 320
599  ia=iabs(i)
600  ei=kchg(iabs(i),1)/3d0
601  ai=sign(1d0,ei+0.1d0)
602  vi=ai-4d0*ei*xwv
603  vali=0.5d0*(vi+ai)
604  vari=0.5d0*(vi-ai)
605  hi=hp*((ei-vali*bwzr)**2+(vali*bwzi)**2+
606  & (ei-vari*bwzr)**2+(vari*bwzi)**2)
607  IF(ia.LE.10) hi=hi*faca/3d0
608  nchn=nchn+1
609  isig(nchn,1)=i
610  isig(nchn,2)=-i
611  isig(nchn,3)=1
612  sigh(nchn)=hi*facbw*hf
613  320 CONTINUE
614 
615  ELSEIF(isub.EQ.194) THEN
616 C...f + fbar -> f' + fbar' via s-channel rho_tc and omega_tc.
617  kfa=kfpr(isubsv,1)
618  alprht=2.91d0*(3d0/itcm(1))
619  hp=aem**2*comfac
620  tanw=sqrt(paru(102)/(1d0-paru(102)))
621  ct2w=(1d0-2d0*paru(102))/(2d0*paru(102)/tanw)
622 
623  qupd=2d0*rtcm(2)-1d0
624  far=sqrt(aem/alprht)
625  fao=far*qupd
626  fzr=far*ct2w
627  fzo=-fao*tanw
628  sfar=far**2
629  sfao=fao**2
630  sfzr=fzr**2
631  sfzo=fzo**2
632  CALL pywidt(23,sh,wdtp,wdte)
633  ssmz=dcmplx(1d0-pmas(23,1)**2/sh,wdtp(0)/shr)
634  CALL pywidt(ktechn+113,sh,wdtp,wdte)
635  ssmr=dcmplx(1d0-pmas(pycomp(ktechn+113),1)**2/sh,wdtp(0)/shr)
636  CALL pywidt(ktechn+223,sh,wdtp,wdte)
637  ssmo=dcmplx(1d0-pmas(pycomp(ktechn+223),1)**2/sh,wdtp(0)/shr)
638  detd=(far*fzo-fao*fzr)**2+ssmz*ssmr*ssmo-sfzr*ssmo-
639  $ sfzo*ssmr-sfar*ssmo*ssmz-sfao*ssmr*ssmz
640  daa=(-sfzr*ssmo - sfzo*ssmr + ssmo*ssmr*ssmz)/detd/sh
641  dzz=(-sfar*ssmo - sfao*ssmr + ssmo*ssmr)/detd/sh
642  daz=(far*fzr*ssmo + fao*fzo*ssmr)/detd/sh
643 
644  xwrht=1d0/(4d0*xw*(1d0-xw))
645  kff=iabs(kfpr(isub,1))
646  ef=kchg(kff,1)/3d0
647  af=sign(1d0,ef+0.1d0)
648  vf=af-4d0*ef*xwv
649  valf=0.5d0*(vf+af)
650  varf=0.5d0*(vf-af)
651  fcof=1d0
652  IF(kff.LE.10) fcof=3d0
653 
654  wid2=1d0
655  IF(kff.GE.6.AND.kff.LE.8) wid2=wids(kff,1)
656  IF(kff.EQ.17.OR.kff.EQ.18) wid2=wids(kff,1)
657  dzz=dzz*dcmplx(xwrht,0d0)
658  daz=daz*dcmplx(sqrt(xwrht),0d0)
659 
660  DO 330 i=mmina,mmaxa
661  IF(i.EQ.0.OR.kfac(1,i)*kfac(2,-i).EQ.0) goto 330
662  ei=kchg(iabs(i),1)/3d0
663  ai=sign(1d0,ei+0.1d0)
664  vi=ai-4d0*ei*xwv
665  vali=0.5d0*(vi+ai)
666  vari=0.5d0*(vi-ai)
667  fcoi=fcof
668  IF(iabs(i).LE.10) fcoi=fcoi/3d0
669  difll=abs(ei*ef*daa+vali*valf*dzz+daz*(ei*valf+ef*vali))**2
670  difrr=abs(ei*ef*daa+vari*varf*dzz+daz*(ei*varf+ef*vari))**2
671  diflr=abs(ei*ef*daa+vali*varf*dzz+daz*(ei*varf+ef*vali))**2
672  difrl=abs(ei*ef*daa+vari*valf*dzz+daz*(ei*valf+ef*vari))**2
673  facsig=(difll+difrr)*((uh-sqm4)**2+sh*sqm4)+
674  & (diflr+difrl)*((th-sqm3)**2+sh*sqm3)
675  nchn=nchn+1
676  isig(nchn,1)=i
677  isig(nchn,2)=-i
678  isig(nchn,3)=1
679  sigh(nchn)=hp*fcoi*facsig*wid2
680  330 CONTINUE
681 
682  ELSEIF(isub.EQ.195) THEN
683 C...f + fbar' -> f'' + fbar''' via s-channel rho_tc+
684  kfa=kfpr(isubsv,1)
685  kfb=kfa+1
686  alprht=2.91d0*(3d0/itcm(1))
687  factc=comfac*(aem**2/12d0/xw**2)*(uh-sqm3)*(uh-sqm4)*3d0
688 
689  fwr=sqrt(aem/alprht)/(2d0*sqrt(xw))
690  CALL pywidt(24,sh,wdtp,wdte)
691  ssmz=dcmplx(1d0-pmas(24,1)**2/sh,wdtp(0)/shr)
692  CALL pywidt(ktechn+213,sh,wdtp,wdte)
693  ssmr=dcmplx(1d0-pmas(pycomp(ktechn+213),1)**2/sh,wdtp(0)/shr)
694 
695  fcof=1d0
696  IF(kfa.LE.8) fcof=3d0
697  detd=ssmz*ssmr-dcmplx(fwr**2,0d0)
698  hp=factc*abs(ssmr/detd)**2/sh**2*fcof
699 
700  DO 350 i=mmin1,mmax1
701  IF(i.EQ.0.OR.kfac(1,i).EQ.0) goto 350
702  ia=iabs(i)
703  DO 340 j=mmin2,mmax2
704  IF(j.EQ.0.OR.kfac(2,j).EQ.0) goto 340
705  ja=iabs(j)
706  IF(i*j.GT.0.OR.mod(ia+ja,2).EQ.0) goto 340
707  IF((ia.LE.10.AND.ja.GT.10).OR.(ia.GT.10.AND.ja.LE.10))
708  & goto 340
709  kchr=(kchg(ia,1)*isign(1,i)+kchg(ja,1)*isign(1,j))/3
710  hi=hp
711  IF(ia.LE.10) hi=hi*vckm((ia+1)/2,(ja+1)/2)/3d0
712  nchn=nchn+1
713  isig(nchn,1)=i
714  isig(nchn,2)=j
715  isig(nchn,3)=1
716  sigh(nchn)=hi*wids(kfa,(5-kchr)/2)*wids(kfb,(5+kchr)/2)
717  340 CONTINUE
718  350 CONTINUE
719  ENDIF
720 
721  ELSEIF(isub.LE.380) THEN
722  IF(isub.EQ.361) THEN
723 C...f + fbar -> W_L W_L, W_L pi_tc, pi_tc pi_tc
724  faca=(sh**2*be34**2-(th-uh)**2)
725  alprht=2.91d0*(3d0/itcm(1))
726  hp=(1d0/12d0)*aem**2*cab2*comfac*faca*3d0
727  far=sqrt(aem/alprht)
728  fao=far*qupd
729  fzr=far*ct2w
730  fzo=-fao*tanw
731  sfar=far**2
732  sfao=fao**2
733  sfzr=fzr**2
734  sfzo=fzo**2
735  CALL pywidt(23,sh,wdtp,wdte)
736  ssmz=dcmplx(1d0-pmas(23,1)**2/sh,wdtp(0)/shr)
737  CALL pywidt(ktechn+113,sh,wdtp,wdte)
738  ssmr=dcmplx(1d0-pmas(pycomp(ktechn+113),1)**2/sh,wdtp(0)/shr)
739  CALL pywidt(ktechn+223,sh,wdtp,wdte)
740  ssmo=dcmplx(1d0-pmas(pycomp(ktechn+223),1)**2/sh,wdtp(0)/shr)
741  detd=(far*fzo-fao*fzr)**2+ssmz*ssmr*ssmo-sfzr*ssmo-
742  $ sfzo*ssmr-sfar*ssmo*ssmz-sfao*ssmr*ssmz
743  darho=-(-far*sfzo+fao*fzo*fzr+far*ssmo*ssmz)/detd/sh
744  dzrho=-(-fzr*sfao+fao*fzo*far+fzr*ssmo)/detd/sh
745  daa=-(sfzo*ssmr+sfzr*ssmo-ssmo*ssmr*ssmz)/detd/sh
746  dzz=-(sfao*ssmr+sfar*ssmo-ssmo*ssmr)/detd/sh
747  daz=(far*fzr*ssmo+fao*fzo*ssmr)/detd/sh
748 
749  DO 360 i=mmina,mmaxa
750  IF(i.EQ.0.OR.kfac(1,i)*kfac(2,-i).EQ.0) goto 360
751  ia=iabs(i)
752  ei=kchg(iabs(i),1)/3d0
753  ai=sign(1d0,ei+0.1d0)
754  vi=ai-4d0*ei*xwv
755  vali=0.25d0*(vi+ai)
756  vari=0.25d0*(vi-ai)
757  f2l=ei*(darho/far+daa+ct2w*daz)+
758  $ vali*(ct2w*dzrho/fzr+ct2w*dzz+daz)/sqrt(xw*xw1)
759  f2r=ei*(darho/far+daa+ct2w*daz)+
760  $ vari*(ct2w*dzrho/fzr+ct2w*dzz+daz)/sqrt(xw*xw1)
761  hi=abs(f2l)**2+abs(f2r)**2
762  IF(ia.LE.10) hi=hi/3d0
763  nchn=nchn+1
764  isig(nchn,1)=i
765  isig(nchn,2)=-i
766  isig(nchn,3)=1
767  IF(kfa.EQ.kfb) THEN
768  sigh(nchn)=hi*hp*wids(pycomp(kfa),1)
769  ELSE
770  sigh(nchn)=hi*hp*wids(pycomp(kfa),2)*wids(pycomp(kfb),3)
771  nchn=nchn+1
772  isig(nchn,1)=i
773  isig(nchn,2)=-i
774  isig(nchn,3)=2
775  sigh(nchn)=hi*hp*wids(pycomp(kfa),3)*wids(pycomp(kfb),2)
776  ENDIF
777  360 CONTINUE
778 
779  ELSEIF(isub.EQ.364) THEN
780 C...f + fbar -> gamma pi_tc, gamma pi_tc', Z pi_tc, Z pi_tc',
781 C...W pi_tc
782  vfac=(th**2+uh**2-2d0*sqm3*sqm4)
783  afac=(th**2+uh**2-2d0*sqm3*sqm4+4d0*sh*sqm3)
784  fanom=sqrt(paru(1)*aem)*itcm(1)/paru(2)**2/rtcm(1)
785 
786  alprht=2.91d0*(3d0/itcm(1))
787  hp=(1d0/24d0)*aem**2*comfac*3d0*sh
788  far=sqrt(aem/alprht)
789  fao=far*qupd
790  fzr=far*ct2w
791  fzo=-fao*tanw
792  sfar=far**2
793  sfao=fao**2
794  sfzr=fzr**2
795  sfzo=fzo**2
796  CALL pywidt(23,sh,wdtp,wdte)
797  ssmz=dcmplx(1d0-pmas(23,1)**2/sh,wdtp(0)/shr)
798  CALL pywidt(ktechn+113,sh,wdtp,wdte)
799  ssmr=dcmplx(1d0-pmas(pycomp(ktechn+113),1)**2/sh,wdtp(0)/shr)
800  CALL pywidt(ktechn+223,sh,wdtp,wdte)
801  ssmo=dcmplx(1d0-pmas(pycomp(ktechn+223),1)**2/sh,wdtp(0)/shr)
802  detd=(far*fzo-fao*fzr)**2+ssmz*ssmr*ssmo-sfzr*ssmo-
803  $ sfzo*ssmr-sfar*ssmo*ssmz-sfao*ssmr*ssmz
804  darho=(-far*sfzo+fao*fzo*fzr+far*ssmo*ssmz)/detd/sh
805  dzrho=(-fzr*sfao+fao*fzo*far+fzr*ssmo)/detd/sh
806  daome=(-fao*sfzr+far*fzo*fzr+fao*ssmr*ssmz)/detd/sh
807  dzome=(-fzo*sfar+far*fao*fzr+fzo*ssmr)/detd/sh
808  daa=(sfzo*ssmr+sfzr*ssmo-ssmo*ssmr*ssmz)/detd/sh
809  dzz=(sfao*ssmr+sfar*ssmo-ssmo*ssmr)/detd/sh
810  daz=(far*fzr*ssmo+fao*fzo*ssmr)/detd/sh
811 
812  DO 370 i=mmina,mmaxa
813  IF(i.EQ.0.OR.kfac(1,i)*kfac(2,-i).EQ.0) goto 370
814  ia=iabs(i)
815  ei=kchg(iabs(i),1)/3d0
816  ai=sign(1d0,ei+0.1d0)
817  vi=ai-4d0*ei*xwv
818  vali=0.25d0*(vi+ai)
819  vari=0.25d0*(vi-ai)
820 C...........Add in anomaly contribution
821  f2l=(ei*darho+vali*dzrho/sqrt(xw*xw1))*vrgp
822  f2l=f2l+(ei*daome+vali*dzome/sqrt(xw*xw1))*vogp
823  f2l=f2l+fanom*(vagp*(ei*daa+vali*daz/sqrt(xw*xw1))+
824  $ vzgp*(ei*daz+vali*dzz/sqrt(xw*xw1)))
825  f2r=(ei*darho+vari*dzrho/sqrt(xw*xw1))*vrgp
826  f2r=f2r+(ei*daome+vari*dzome/sqrt(xw*xw1))*vogp
827  f2r=f2r+fanom*(vagp*(ei*daa+vari*daz/sqrt(xw*xw1))+
828  $ vzgp*(ei*daz+vari*dzz/sqrt(xw*xw1)))
829  hi=(abs(f2l)**2+abs(f2r)**2)*vfac
830  f2l=(ei*darho+vali*dzrho/sqrt(xw*xw1))*argp
831  f2l=f2l+(ei*daome+vali*dzome/sqrt(xw*xw1))*aogp
832  f2r=(ei*darho+vari*dzrho/sqrt(xw*xw1))*argp
833  f2r=f2r+(ei*daome+vari*dzome/sqrt(xw*xw1))*aogp
834  hj=(abs(f2l)**2+abs(f2r)**2)*afac
835  hi=hi+hj
836  IF(ia.LE.10) hi=hi/3d0
837  nchn=nchn+1
838  isig(nchn,1)=i
839  isig(nchn,2)=-i
840  isig(nchn,3)=1
841  IF(isubsv.NE.368) THEN
842  sigh(nchn)=hi*hp*wids(pycomp(kfa),2)*wids(pycomp(kfb),2)
843  ELSE
844  sigh(nchn)=hi*hp*wids(pycomp(kfa),2)*wids(pycomp(kfb),3)
845  nchn=nchn+1
846  isig(nchn,1)=i
847  isig(nchn,2)=-i
848  isig(nchn,3)=2
849  sigh(nchn)=hi*hp*wids(pycomp(kfa),3)*wids(pycomp(kfb),2)
850  ENDIF
851  370 CONTINUE
852 
853  ELSEIF(isub.EQ.370) THEN
854 C...f + fbar' -> W_L Z_L, W_L pi_tc, Z_L pi_tc, pi_tc pi_tc
855 
856  faca=(sh**2*be34**2-(th-uh)**2)
857  alprht=2.91d0*(3d0/itcm(1))
858  hp=(1d0/96d0)*aem**2*cab2*comfac*faca*3d0/xw**2
859  fwr=sqrt(aem/alprht)/(2d0*sqrt(xw))
860  CALL pywidt(24,sh,wdtp,wdte)
861  ssmz=dcmplx(1d0-pmas(24,1)**2/sh,wdtp(0)/shr)
862  CALL pywidt(ktechn+213,sh,wdtp,wdte)
863  ssmr=dcmplx(1d0-pmas(pycomp(ktechn+213),1)**2/sh,wdtp(0)/shr)
864  detd=ssmz*ssmr-dcmplx(fwr**2,0d0)
865  dww=ssmr/detd/sh
866  dwrho=-1d0/detd/sh
867  hp=hp*abs(dww+dwrho)**2
868  DO 390 i=mmin1,mmax1
869  IF(i.EQ.0.OR.kfac(1,i).EQ.0) goto 390
870  ia=iabs(i)
871  DO 380 j=mmin2,mmax2
872  IF(j.EQ.0.OR.kfac(2,j).EQ.0) goto 380
873  ja=iabs(j)
874  IF(i*j.GT.0.OR.mod(ia+ja,2).EQ.0) goto 380
875  IF((ia.LE.10.AND.ja.GT.10).OR.(ia.GT.10.AND.ja.LE.10))
876  & goto 380
877  kchr=(kchg(ia,1)*isign(1,i)+kchg(ja,1)*isign(1,j))/3
878  hi=hp
879  IF(ia.LE.10) hi=hi*vckm((ia+1)/2,(ja+1)/2)/3d0
880  nchn=nchn+1
881  isig(nchn,1)=i
882  isig(nchn,2)=j
883  isig(nchn,3)=1
884  sigh(nchn)=hi*wids(pycomp(kfa),(5-kchr)/2)*
885  & wids(pycomp(kfb),2)
886  380 CONTINUE
887  390 CONTINUE
888 
889  ELSEIF(isub.EQ.374) THEN
890 C...f + fbar' -> gamma pi_tc
891  fanom=sqrt(aem)*itcm(1)/2d0/paru(2)/rtcm(1)
892  vfac=(th**2+uh**2-2d0*sqm3*sqm4)
893  afac=(th**2+uh**2-2d0*sqm3*sqm4+4d0*sh*sqm3)/sqta*argp**2
894  alprht=2.91d0*(3d0/itcm(1))
895  hp=(1d0/48d0)*aem**2/xw*comfac*3d0*sh
896  fwr=sqrt(aem/alprht)/(2d0*sqrt(xw))
897  CALL pywidt(24,sh,wdtp,wdte)
898  ssmz=dcmplx(1d0-pmas(24,1)**2/sh,wdtp(0)/shr)
899  CALL pywidt(ktechn+213,sh,wdtp,wdte)
900  ssmr=dcmplx(1d0-pmas(pycomp(ktechn+213),1)**2/sh,wdtp(0)/shr)
901  detd=ssmz*ssmr-dcmplx(fwr**2,0d0)
902  dww=ssmr/detd/sh
903  dwrho=-dcmplx(fwr,0d0)/detd/sh
904  hp=hp*(afac*abs(dwrho)**2+
905  $ vfac*abs(fanom*dww*vwgp+dwrho*vrgp/sqrt(sqtv))**2)
906  DO 410 i=mmin1,mmax1
907  IF(i.EQ.0.OR.kfac(1,i).EQ.0) goto 410
908  ia=iabs(i)
909  DO 400 j=mmin2,mmax2
910  IF(j.EQ.0.OR.kfac(2,j).EQ.0) goto 400
911  ja=iabs(j)
912  IF(i*j.GT.0.OR.mod(ia+ja,2).EQ.0) goto 400
913  IF((ia.LE.10.AND.ja.GT.10).OR.(ia.GT.10.AND.ja.LE.10))
914  & goto 400
915  kchr=(kchg(ia,1)*isign(1,i)+kchg(ja,1)*isign(1,j))/3
916  hi=hp
917  IF(ia.LE.10) hi=hi*vckm((ia+1)/2,(ja+1)/2)/3d0
918  nchn=nchn+1
919  isig(nchn,1)=i
920  isig(nchn,2)=j
921  isig(nchn,3)=1
922  sigh(nchn)=hi*wids(pycomp(kfa),(5-kchr)/2)*
923  & wids(pycomp(kfb),2)
924  400 CONTINUE
925  410 CONTINUE
926  ENDIF
927 
928  ELSEIF(isub.LE.390) THEN
929  IF(isub.EQ.381) THEN
930 C...f + f' -> f + f' (g exchange)
931  facqq1=comfac*as**2*4d0/9d0*(sh2+uh2)*sqdqqt
932  facqqb=comfac*as**2*4d0/9d0*((sh2+uh2)*sqdqqt*faca-
933  & mstp(34)*2d0/3d0*uh2*redqst)
934  facqq2=comfac*as**2*4d0/9d0*(sh2+th2)*sqdqqu
935  facqqi=-comfac*as**2*4d0/9d0*mstp(34)*2d0/3d0*sh2/(th*uh)
936  ratqqi=(facqq1+facqq2+facqqi)/(facqq1+facqq2)
937  IF(itcm(5).GE.1.AND.itcm(5).LE.4) THEN
938 C...Modifications from contact interactions (compositeness)
939  facci1=facqq1+comfac*(sh2/rtcm(41)**4)
940  faccib=facqqb+comfac*(8d0/9d0)*(as*rtcm(42)/rtcm(41)**2)*
941  & (uh2/th+uh2/sh)+comfac*(5d0/3d0)*(uh2/rtcm(41)**4)
942  facci2=facqq2+comfac*(8d0/9d0)*(as*rtcm(42)/rtcm(41)**2)*
943  & (sh2/th+sh2/uh)+comfac*(5d0/3d0)*(sh2/rtcm(41)**4)
944  facci3=facqq1+comfac*(uh2/rtcm(41)**4)
945  ratcii=(facci1+facci2+facqqi)/(facci1+facci2)
946  ELSEIF(itcm(5).EQ.5) THEN
947  facci1=facqq1
948  faccib=facqqb
949  facci2=facqq2
950  facci3=facqq1
951 CSM.......Check this change from
952 CSM RATCII=1D0
953  ratcii=ratqqi
954  ENDIF
955  DO 430 i=mmin1,mmax1
956  ia=iabs(i)
957  IF(i.EQ.0.OR.ia.GT.mstp(58).OR.kfac(1,i).EQ.0) goto 430
958  DO 420 j=mmin2,mmax2
959  ja=iabs(j)
960  IF(j.EQ.0.OR.ja.GT.mstp(58).OR.kfac(2,j).EQ.0) goto 420
961  nchn=nchn+1
962  isig(nchn,1)=i
963  isig(nchn,2)=j
964  isig(nchn,3)=1
965  IF(itcm(5).LE.0.OR.(itcm(5).EQ.1.AND.(ia.GE.3.OR.
966  & ja.GE.3))) THEN
967  sigh(nchn)=facqq1
968  IF(i.EQ.-j) sigh(nchn)=facqqb
969  ELSE
970  sigh(nchn)=facci1
971  IF(i*j.LT.0) sigh(nchn)=facci3
972  IF(i.EQ.-j) sigh(nchn)=faccib
973  ENDIF
974  IF(i.EQ.j) THEN
975  nchn=nchn+1
976  isig(nchn,1)=i
977  isig(nchn,2)=j
978  isig(nchn,3)=2
979  IF(itcm(5).LE.0.OR.(itcm(5).EQ.1.AND.ia.GE.3)) THEN
980  sigh(nchn-1)=0.5d0*facqq1*ratqqi
981  sigh(nchn)=0.5d0*facqq2*ratqqi
982  ELSE
983  sigh(nchn-1)=0.5d0*facci1*ratcii
984  sigh(nchn)=0.5d0*facci2*ratcii
985  ENDIF
986  ENDIF
987  420 CONTINUE
988  430 CONTINUE
989 
990  ELSEIF(isub.EQ.382) THEN
991 C...f + fbar -> f' + fbar' (q + qbar -> q' + qbar' only)
992  CALL pywidt(21,sh,wdtp,wdte)
993  facqqf=comfac*as**2*4d0/9d0*(th2+uh2)
994  facqqb=facqqf*sqdqqs*(wdte(0,1)+wdte(0,2)+wdte(0,4))
995  IF(itcm(5).EQ.1) THEN
996 C...Modifications from contact interactions (compositeness)
997  faccib=facqqb
998  DO 440 i=1,2
999  faccib=faccib+comfac*(uh2/rtcm(41)**4)*(wdte(i,1)+
1000  & wdte(i,2)+wdte(i,4))
1001  440 CONTINUE
1002  ELSEIF(itcm(5).GE.2.AND.itcm(5).LE.4) THEN
1003  faccib=facqqb+comfac*(uh2/rtcm(41)**4)*
1004  & (wdte(0,1)+wdte(0,2)+wdte(0,4))
1005  ELSEIF(itcm(5).EQ.5) THEN
1006  facqqb=facqqf*sqdqqs*(wdte(0,1)+wdte(0,2)+wdte(0,4)-
1007  & wdte(5,1)-wdte(5,2)-wdte(5,4))
1008  faccib=facqqf*sqdqts*(wdte(5,1)+wdte(5,2)+wdte(5,4))
1009  ENDIF
1010  DO 450 i=mmina,mmaxa
1011  IF(i.EQ.0.OR.iabs(i).GT.mstp(58).OR.
1012  & kfac(1,i)*kfac(2,-i).EQ.0) goto 450
1013  nchn=nchn+1
1014  isig(nchn,1)=i
1015  isig(nchn,2)=-i
1016  isig(nchn,3)=1
1017  IF(itcm(5).LE.0.OR.(itcm(5).EQ.1.AND.iabs(i).GE.3)) THEN
1018  sigh(nchn)=facqqb
1019  ELSEIF(itcm(5).EQ.5) THEN
1020  sigh(nchn)=facqqb
1021  nchn=nchn+1
1022  isig(nchn,1)=i
1023  isig(nchn,2)=-i
1024  isig(nchn,3)=2
1025  sigh(nchn)=faccib
1026  ELSE
1027  sigh(nchn)=faccib
1028  ENDIF
1029  450 CONTINUE
1030 
1031  ELSEIF(isub.EQ.383) THEN
1032 C...f + fbar -> g + g (q + qbar -> g + g only)
1033  facgg1=comfac*as**2*32d0/27d0*(uh/th-(2d0+mstp(34)*1d0/4d0)*
1034  & uh2/sh2+9d0/4d0*th*uh/sh2*sqdlgs)
1035  facgg2=comfac*as**2*32d0/27d0*(th/uh-(2d0+mstp(34)*1d0/4d0)*
1036  & th2/sh2+9d0/4d0*th*uh/sh2*sqdlgs)
1037  IF(itcm(5).EQ.5) THEN
1038  facgg3=comfac*as**2*32d0/27d0*(uh/th-(2d0+mstp(34)*1d0/4d0)*
1039  & uh2/sh2+9d0/4d0*th*uh/sh2*sqdhgs)
1040  facgg4=comfac*as**2*32d0/27d0*(th/uh-(2d0+mstp(34)*1d0/4d0)*
1041  & th2/sh2+9d0/4d0*th*uh/sh2*sqdhgs)
1042  ENDIF
1043  DO 460 i=mmina,mmaxa
1044  IF(i.EQ.0.OR.iabs(i).GT.mstp(58).OR.
1045  & kfac(1,i)*kfac(2,-i).EQ.0) goto 460
1046  nchn=nchn+1
1047  isig(nchn,1)=i
1048  isig(nchn,2)=-i
1049  isig(nchn,3)=1
1050  sigh(nchn)=0.5d0*facgg1
1051  IF(itcm(5).EQ.5.AND.iabs(i).EQ.5) sigh(nchn)=0.5d0*facgg3
1052  nchn=nchn+1
1053  isig(nchn,1)=i
1054  isig(nchn,2)=-i
1055  isig(nchn,3)=2
1056  sigh(nchn)=0.5d0*facgg2
1057  IF(itcm(5).EQ.5.AND.iabs(i).EQ.5) sigh(nchn)=0.5d0*facgg4
1058  460 CONTINUE
1059 
1060  ELSEIF(isub.EQ.384) THEN
1061 C...f + g -> f + g (q + g -> q + g only)
1062  facqg1=comfac*as**2*4d0/9d0*((2d0+mstp(34)*1d0/4d0)*uh2/th2-
1063  & uh/sh-9d0/4d0*sh*uh/th2*sqdlgt)*faca
1064  facqg2=comfac*as**2*4d0/9d0*((2d0+mstp(34)*1d0/4d0)*sh2/th2-
1065  & sh/uh-9d0/4d0*sh*uh/th2*sqdlgt)
1066  DO 480 i=mmina,mmaxa
1067  IF(i.EQ.0.OR.iabs(i).GT.10) goto 480
1068  DO 470 isde=1,2
1069  IF(isde.EQ.1.AND.kfac(1,i)*kfac(2,21).EQ.0) goto 470
1070  IF(isde.EQ.2.AND.kfac(1,21)*kfac(2,i).EQ.0) goto 470
1071  nchn=nchn+1
1072  isig(nchn,isde)=i
1073  isig(nchn,3-isde)=21
1074  isig(nchn,3)=1
1075  sigh(nchn)=facqg1
1076  nchn=nchn+1
1077  isig(nchn,isde)=i
1078  isig(nchn,3-isde)=21
1079  isig(nchn,3)=2
1080  sigh(nchn)=facqg2
1081  470 CONTINUE
1082  480 CONTINUE
1083 
1084  ELSEIF(isub.EQ.385) THEN
1085 C...g + g -> f + fbar (g + g -> q + qbar only)
1086  IF(kfac(1,21)*kfac(2,21).EQ.0) goto 500
1087  idc0=mdcy(21,2)-1
1088 C...Begin by d, u, s flavours.
1089  flavwt=0d0
1090  IF(mdme(idc0+1,1).GE.1) flavwt=flavwt+
1091  & sqrt(max(0d0,1d0-4d0*pmas(1,1)**2/sh))
1092  IF(mdme(idc0+2,1).GE.1) flavwt=flavwt+
1093  & sqrt(max(0d0,1d0-4d0*pmas(2,1)**2/sh))
1094  IF(mdme(idc0+3,1).GE.1) flavwt=flavwt+
1095  & sqrt(max(0d0,1d0-4d0*pmas(3,1)**2/sh))
1096  facqq1=comfac*as**2*1d0/6d0*(uh/th-(2d0+mstp(34)*1d0/4d0)*
1097  & uh2/sh2+9d0/4d0*th*uh/sh2*sqdlgs)*flavwt*faca
1098  facqq2=comfac*as**2*1d0/6d0*(th/uh-(2d0+mstp(34)*1d0/4d0)*
1099  & th2/sh2+9d0/4d0*th*uh/sh2*sqdlgs)*flavwt*faca
1100  nchn=nchn+1
1101  isig(nchn,1)=21
1102  isig(nchn,2)=21
1103  isig(nchn,3)=1
1104  sigh(nchn)=facqq1
1105  nchn=nchn+1
1106  isig(nchn,1)=21
1107  isig(nchn,2)=21
1108  isig(nchn,3)=2
1109  sigh(nchn)=facqq2
1110 C...Next c and b flavours: modified that and uhat for fixed
1111 C...cos(theta-hat).
1112  DO 490 ifl=4,5
1113  sqmavg=pmas(ifl,1)**2
1114  IF(mdme(idc0+ifl,1).GE.1.AND.sh.GT.4.04d0*sqmavg) THEN
1115  be34=sqrt(1d0-4d0*sqmavg/sh)
1116  thq=-0.5d0*sh*(1d0-be34*cth)
1117  uhq=-0.5d0*sh*(1d0+be34*cth)
1118  thuhq=thq*uhq-sqmavg*sh
1119  IF(mstp(34).EQ.0) THEN
1120  facqq1=uhq/thq-2d0*uhq**2/sh2+4d0*(sqmavg/sh)*thuhq/thq**2
1121  facqq2=thq/uhq-2d0*thq**2/sh2+4d0*(sqmavg/sh)*thuhq/uhq**2
1122  ELSE
1123  facqq1=uhq/thq-2.25d0*uhq**2/sh2+4.5d0*(sqmavg/sh)*thuhq/
1124  & thq**2+0.5d0*sqmavg*(thq+sqmavg)/thq**2-sqmavg**2/(sh*thq)
1125  facqq2=thq/uhq-2.25d0*thq**2/sh2+4.5d0*(sqmavg/sh)*thuhq/
1126  & uhq**2+0.5d0*sqmavg*(uhq+sqmavg)/uhq**2-sqmavg**2/(sh*uhq)
1127  ENDIF
1128  IF(itcm(5).GE.5) THEN
1129  IF(ifl.EQ.4) THEN
1130  facqq1=facqq1+2.25d0*sqmavg*(thq-uhq)/(sh*thq)*redlgs+
1131  & 2.25d0*thq*uhq/sh2*sqdlgs
1132  facqq2=facqq2+2.25d0*sqmavg*(uhq-thq)/(sh*uhq)*redlgs+
1133  & 2.25d0*thq*uhq/sh2*sqdlgs
1134  ELSE
1135  facqq1=facqq1+2.25d0*sqmavg*(thq-uhq)/(sh*thq)*redhgs+
1136  & 2.25d0*thq*uhq/sh2*sqdhgs
1137  facqq2=facqq2+2.25d0*sqmavg*(uhq-thq)/(sh*uhq)*redhgs+
1138  & 2.25d0*thq*uhq/sh2*sqdhgs
1139  ENDIF
1140  ENDIF
1141  facqq1=comfac*faca*as**2*(1d0/6d0)*facqq1*be34
1142  facqq2=comfac*faca*as**2*(1d0/6d0)*facqq2*be34
1143  nchn=nchn+1
1144  isig(nchn,1)=21
1145  isig(nchn,2)=21
1146  isig(nchn,3)=1+2*(ifl-3)
1147  sigh(nchn)=facqq1
1148  nchn=nchn+1
1149  isig(nchn,1)=21
1150  isig(nchn,2)=21
1151  isig(nchn,3)=2+2*(ifl-3)
1152  sigh(nchn)=facqq2
1153  ENDIF
1154  490 CONTINUE
1155  500 CONTINUE
1156 
1157  ELSEIF(isub.EQ.386) THEN
1158 C...g + g -> g + g
1159  IF(itcm(5).LE.4) THEN
1160  facgg1=comfac*as**2*9d0/4d0*(sh2/th2+2d0*sh/th+3d0+
1161  & 2d0*th/sh+th2/sh2)*faca
1162  facgg2=comfac*as**2*9d0/4d0*(uh2/sh2+2d0*uh/sh+3d0+
1163  & 2d0*sh/uh+sh2/uh2)*faca
1164  facgg3=comfac*as**2*9d0/4d0*(th2/uh2+2d0*th/uh+3d0+
1165  & 2d0*uh/th+uh2/th2)
1166  ELSE
1167  gst= (12d0 + 40d0*th/sh + 56d0*th2/sh2 + 32d0*th**3/sh**3 +
1168  & 16d0*th**4/sh**4 + sqdggs*(4d0*sh2 + 16d0*sh*th + 16d0*th2)+
1169  & 4d0*redgst*(sh + 2d0*th)*
1170  & (2d0*sh**3 - 3d0*sh2*th - 2d0*sh*th2 + 2d0*th**3)/sh2 +
1171  & 2d0*redggs*(2d0*sh - 12d0*th2/sh - 8d0*th**3/sh2) +
1172  & 2d0*redggt*(4d0*sh - 22d0*th - 68d0*th2/sh - 60d0*th**3/sh2-
1173  & 32d0*th**4/sh**3 - 16d0*th**5/sh**4) +
1174  & sqdggt*(16d0*sh2 + 16d0*sh*th + 68d0*th2 + 144d0*th**3/sh +
1175  & 96d0*th**4/sh2 + 32d0*th**5/sh**3 + 16d0*th**6/sh**4))/16d0
1176  gsu= (12d0 + 40d0*uh/sh + 56d0*uh2/sh2 + 32d0*uh**3/sh**3 +
1177  & 16d0*uh**4/sh**4 + sqdggs*(4d0*sh2 + 16d0*sh*uh + 16d0*uh2)+
1178  & 4d0*redgsu*(sh + 2d0*uh)*
1179  & (2d0*sh**3 - 3d0*sh2*uh - 2d0*sh*uh2 + 2d0*uh**3)/sh2 +
1180  & 2d0*redggs*(2d0*sh - 12d0*uh2/sh - 8d0*uh**3/sh2) +
1181  & 2d0*redggu*(4d0*sh - 22d0*uh - 68d0*uh2/sh - 60d0*uh**3/sh2-
1182  & 32d0*uh**4/sh**3 - 16d0*uh**5/sh**4) +
1183  & sqdggu*(16d0*sh2 + 16d0*sh*uh + 68d0*uh2 + 144d0*uh**3/sh +
1184  & 96d0*uh**4/sh2 + 32d0*uh**5/sh**3 + 16d0*uh**6/sh**4))/16d0
1185  gut= (12d0 - 16d0*th*(th - uh)**2*uh/sh**4 +
1186  & 4d0*redggu*(2d0*th**5 - 15d0*th**4*uh - 48d0*th**3*uh2 -
1187  & 58d0*th2*uh**3 - 10d0*th*uh**4 + uh**5)/sh**4 +
1188  & 4d0*redggt*(th**5 - 10d0*th**4*uh - 58d0*th**3*uh2 -
1189  & 48d0*th2*uh**3 - 15d0*th*uh**4 + 2d0*uh**5)/sh**4 +
1190  & 4d0*sqdggu*(4d0*th**6 + 20d0*th**5*uh + 57d0*th**4*uh2 +
1191  & 72d0*th**3*uh**3+ 38d0*th2*uh**4+4d0*th*uh**5 +uh**6)/sh**4+
1192  & 4d0*sqdggt*(4d0*uh**6 + 4d0*th**5*uh + 38d0*th**4*uh2 +
1193  & 72d0*th**3*uh**3 +57d0*th2*uh**4+20d0*th*uh**5+th**6)/sh**4+
1194  & 2d0*redgtu*((th - uh)**2* (th**4 + 20d0*th**3*uh +
1195  & 30d0*th2*uh2 + 20d0*th*uh**3 + uh**4) +
1196  & sh2*(7d0*th**4 + 52d0*th**3*uh + 274d0*th2*uh2 +
1197  & 52d0*th*uh**3 + 7d0*uh**4))/(2d0*sh**4))/16d0
1198  facgg1=comfac*as**2*9d0/4d0*gst*faca
1199  facgg2=comfac*as**2*9d0/4d0*gsu*faca
1200  facgg3=comfac*as**2*9d0/4d0*gut
1201  ENDIF
1202  IF(kfac(1,21)*kfac(2,21).EQ.0) goto 510
1203  nchn=nchn+1
1204  isig(nchn,1)=21
1205  isig(nchn,2)=21
1206  isig(nchn,3)=1
1207  sigh(nchn)=0.5d0*facgg1
1208  nchn=nchn+1
1209  isig(nchn,1)=21
1210  isig(nchn,2)=21
1211  isig(nchn,3)=2
1212  sigh(nchn)=0.5d0*facgg2
1213  nchn=nchn+1
1214  isig(nchn,1)=21
1215  isig(nchn,2)=21
1216  isig(nchn,3)=3
1217  sigh(nchn)=0.5d0*facgg3
1218  510 CONTINUE
1219 
1220  ELSEIF(isub.EQ.387) THEN
1221 C...q + qbar -> Q + Qbar
1222  sqmavg=0.5d0*(sqm3+sqm4)-0.25d0*(sqm3-sqm4)**2/sh
1223  thq=-0.5d0*sh*(1d0-be34*cth)
1224  uhq=-0.5d0*sh*(1d0+be34*cth)
1225  facqqb=comfac*as**2*4d0/9d0*((thq**2+uhq**2)/sh2+
1226  & 2d0*sqmavg/sh)
1227  IF(itcm(5).GE.5) THEN
1228  IF(mint(55).EQ.5.OR.mint(55).EQ.6) THEN
1229  facqqb=facqqb*sh2*sqdqts
1230  ELSE
1231  facqqb=facqqb*sh2*sqdqqs
1232  ENDIF
1233  ENDIF
1234  IF(mstp(35).GE.1) facqqb=facqqb*pyhfth(sh,sqmavg,0d0)
1235  wid2=1d0
1236  IF(mint(55).EQ.6) wid2=wids(6,1)
1237  IF(mint(55).EQ.7.OR.mint(55).EQ.8) wid2=wids(mint(55),1)
1238  facqqb=facqqb*wid2
1239  DO 520 i=mmina,mmaxa
1240  IF(i.EQ.0.OR.iabs(i).GT.mstp(58).OR.
1241  & kfac(1,i)*kfac(2,-i).EQ.0) goto 520
1242  nchn=nchn+1
1243  isig(nchn,1)=i
1244  isig(nchn,2)=-i
1245  isig(nchn,3)=1
1246  sigh(nchn)=facqqb
1247  520 CONTINUE
1248 
1249  ELSEIF(isub.EQ.388) THEN
1250 C...g + g -> Q + Qbar
1251  sqmavg=0.5d0*(sqm3+sqm4)-0.25d0*(sqm3-sqm4)**2/sh
1252  thq=-0.5d0*sh*(1d0-be34*cth)
1253  uhq=-0.5d0*sh*(1d0+be34*cth)
1254  thuhq=thq*uhq-sqmavg*sh
1255  IF(mstp(34).EQ.0) THEN
1256  facqq1=uhq/thq-2d0*uhq**2/sh2+4d0*(sqmavg/sh)*thuhq/thq**2
1257  facqq2=thq/uhq-2d0*thq**2/sh2+4d0*(sqmavg/sh)*thuhq/uhq**2
1258  ELSE
1259  facqq1=uhq/thq-2.25d0*uhq**2/sh2+4.5d0*(sqmavg/sh)*thuhq/
1260  & thq**2+0.5d0*sqmavg*(thq+sqmavg)/thq**2-sqmavg**2/(sh*thq)
1261  facqq2=thq/uhq-2.25d0*thq**2/sh2+4.5d0*(sqmavg/sh)*thuhq/
1262  & uhq**2+0.5d0*sqmavg*(uhq+sqmavg)/uhq**2-sqmavg**2/(sh*uhq)
1263  ENDIF
1264  IF(itcm(5).GE.5) THEN
1265  IF(mint(55).EQ.5.OR.mint(55).EQ.6) THEN
1266  facqq1=facqq1+2.25d0*sqmavg*(thq-uhq)/(sh*thq)*redhgs+
1267  & 2.25d0*thq*uhq/sh2*sqdhgs
1268  facqq2=facqq2+2.25d0*sqmavg*(uhq-thq)/(sh*uhq)*redhgs+
1269  & 2.25d0*thq*uhq/sh2*sqdhgs
1270  ELSE
1271  facqq1=facqq1+2.25d0*sqmavg*(thq-uhq)/(sh*thq)*redlgs+
1272  & 2.25d0*thq*uhq/sh2*sqdlgs
1273  facqq2=facqq2+2.25d0*sqmavg*(uhq-thq)/(sh*uhq)*redlgs+
1274  & 2.25d0*thq*uhq/sh2*sqdlgs
1275  ENDIF
1276  ENDIF
1277  facqq1=comfac*faca*as**2*(1d0/6d0)*facqq1
1278  facqq2=comfac*faca*as**2*(1d0/6d0)*facqq2
1279  IF(mstp(35).GE.1) THEN
1280  fatre=pyhfth(sh,sqmavg,2d0/7d0)
1281  facqq1=facqq1*fatre
1282  facqq2=facqq2*fatre
1283  ENDIF
1284  wid2=1d0
1285  IF(mint(55).EQ.6) wid2=wids(6,1)
1286  IF(mint(55).EQ.7.OR.mint(55).EQ.8) wid2=wids(mint(55),1)
1287  facqq1=facqq1*wid2
1288  facqq2=facqq2*wid2
1289  IF(kfac(1,21)*kfac(2,21).EQ.0) goto 530
1290  nchn=nchn+1
1291  isig(nchn,1)=21
1292  isig(nchn,2)=21
1293  isig(nchn,3)=1
1294  sigh(nchn)=facqq1
1295  nchn=nchn+1
1296  isig(nchn,1)=21
1297  isig(nchn,2)=21
1298  isig(nchn,3)=2
1299  sigh(nchn)=facqq2
1300  530 CONTINUE
1301  ENDIF
1302  ENDIF
1303 
1304 CMRENNA--
1305 
1306  RETURN
1307  END