Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
pyglui.f
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file pyglui.f
1 
2 C*********************************************************************
3 
4 C...PYGLUI
5 C...Calculates gluino decay modes.
6 
7  SUBROUTINE pyglui(KFIN,XLAM,IDLAM,IKNT)
8 
9 C...Double precision and integer declarations.
10  IMPLICIT DOUBLE PRECISION(a-h, o-z)
11  IMPLICIT INTEGER(i-n)
12  INTEGER pyk,pychge,pycomp
13 C...Parameter statement to help give large particle numbers.
14  parameter(ksusy1=1000000,ksusy2=2000000,ktechn=3000000,
15  &kexcit=4000000,kdimen=5000000)
16 C...Commonblocks.
17  common/pydat1/mstu(200),paru(200),mstj(200),parj(200)
18  common/pydat2/kchg(500,4),pmas(500,4),parf(2000),vckm(4,4)
19  common/pymssm/imss(0:99),rmss(0:99)
20  common/pyssmt/zmix(4,4),umix(2,2),vmix(2,2),smz(4),smw(2),
21  &sfmix(16,4),zmixi(4,4),umixi(2,2),vmixi(2,2)
22 CC &SFMIX(16,4),
23 C COMMON/PYINTS/XXM(20)
24  COMPLEX*16 cxc
25  common/pyintc/xxc(10),cxc(8)
26  SAVE /pydat1/,/pydat2/,/pymssm/,/pyssmt/,/pyintc/
27 
28 C...Local variables
29  COMPLEX*16 zmixc(4,4),vmixc(2,2),umixc(2,2),olpp,orpp,glij,grij
30  DOUBLE PRECISION xmi,xmj,xmf,axmj,axmi
31  DOUBLE PRECISION xmi2,xmi3,xma2,xmb2,xmfp
32  DOUBLE PRECISION pylamf,xl
33  DOUBLE PRECISION tanw,xw,aem,c1,as,s12max,s12min
34  DOUBLE PRECISION ca,cb,al,ar,bl,br
35  DOUBLE PRECISION xlam(0:400)
36  INTEGER idlam(400,3)
37  INTEGER lknt,ix,ilr,i,iknt,ifl
38  DOUBLE PRECISION sr2
39  DOUBLE PRECISION gam
40  DOUBLE PRECISION pyalem,pi,pyalps,ei,t3i
41  EXTERNAL pygaus,pyxxz6
42  DOUBLE PRECISION pygaus,pyxxz6
43  DOUBLE PRECISION prec
44  INTEGER kfnchi(4),kfcchi(2)
45  DATA pi/3.141592654d0/
46  DATA sr2/1.4142136d0/
47  DATA prec/1d-2/
48  DATA kfnchi/1000022,1000023,1000025,1000035/
49  DATA kfcchi/1000024,1000037/
50 
51 C...COUNT THE NUMBER OF DECAY MODES
52  lknt=0
53  IF(kfin.NE.ksusy1+21) RETURN
54  kcin=pycomp(kfin)
55 
56  xw=paru(102)
57  tanw = sqrt(xw/(1d0-xw))
58 
59  xmi=pmas(kcin,1)
60  axmi=abs(xmi)
61  xmi2=xmi**2
62  aem=pyalem(xmi2)
63  as =pyalps(xmi2)
64  c1=aem/xw
65  xmi3=axmi**3
66 
67  xmi=sign(xmi,rmss(3))
68 
69 C...2-BODY DECAYS OF GLUINO -> GRAVITINO GLUON
70 
71  IF(imss(11).EQ.1) THEN
72  xmp=rmss(29)
73  idg=39+ksusy1
74  xmgr=pmas(pycomp(idg),1)
75  xfac=(xmi2/(xmp*xmgr))**2*axmi/48d0/pi
76  IF(axmi.GT.xmgr) THEN
77  lknt=lknt+1
78  idlam(lknt,1)=idg
79  idlam(lknt,2)=21
80  idlam(lknt,3)=0
81  xlam(lknt)=xfac
82  ENDIF
83  ENDIF
84 
85 C...2-BODY DECAYS OF GLUINO -> QUARK SQUARK
86 
87  DO 110 ifl=1,6
88  DO 100 ilr=1,2
89  xmj=pmas(pycomp(ilr*ksusy1+ifl),1)
90  axmj=abs(xmj)
91  xmf=pmas(ifl,1)
92  IF(axmi.GE.axmj+xmf) THEN
93 C...Minus sign difference from gluino-quark-squark feynman rules
94  al=sfmix(ifl,1)
95  bl=-sfmix(ifl,3)
96  ar=sfmix(ifl,2)
97  br=-sfmix(ifl,4)
98 C...F1 -> F CHI
99  IF(ilr.EQ.1) THEN
100  ca=al
101  cb=bl
102 C...F2 -> F CHI
103  ELSE
104  ca=ar
105  cb=br
106  ENDIF
107  lknt=lknt+1
108  xma2=xmj**2
109  xmb2=xmf**2
110  xl=pylamf(xmi2,xma2,xmb2)
111  xlam(lknt)=4d0/8d0*as/4d0/xmi3*sqrt(xl)*((xmi2+xmb2-xma2)*
112  & (ca**2+cb**2)-4d0*ca*cb*xmi*xmf)
113  idlam(lknt,1)=ilr*ksusy1+ifl
114  idlam(lknt,2)=-ifl
115  idlam(lknt,3)=0
116  lknt=lknt+1
117  xlam(lknt)=xlam(lknt-1)
118  idlam(lknt,1)=-idlam(lknt-1,1)
119  idlam(lknt,2)=-idlam(lknt-1,2)
120  idlam(lknt,3)=0
121  ENDIF
122  100 CONTINUE
123  110 CONTINUE
124 
125 C...3-BODY DECAYS TO GAUGINO FERMION-FERMION
126 C...GLUINO -> NI Q QBAR
127  DO 170 ix=1,4
128  xmj=smz(ix)
129  axmj=abs(xmj)
130  IF(axmi.GE.axmj) THEN
131  DO 120 i=1,4
132  zmixc(ix,i)=dcmplx(zmix(ix,i),zmixi(ix,i))
133  120 CONTINUE
134  olpp=dcmplx(cos(rmss(32)),sin(rmss(32)))/sr2
135  orpp=dconjg(olpp)
136  xxc(1)=0d0
137  xxc(2)=xmj
138  xxc(3)=0d0
139  xxc(4)=xmi
140  ia=1
141  xxc(5)=pmas(pycomp(ksusy1+ia),1)
142  xxc(6)=pmas(pycomp(ksusy2+ia),1)
143  xxc(7)=xxc(5)
144  xxc(8)=xxc(6)
145  xxc(9)=1d6
146  xxc(10)=0d0
147  ei=kchg(ia,1)/3d0
148  t3i=sign(1d0,ei+1d-6)/2d0
149  glij=(t3i*zmixc(ix,2)-tanw*(t3i-ei)*zmixc(ix,1))*olpp
150  grij=zmixc(ix,1)*(ei*tanw)*orpp
151  cxc(1)=0d0
152  cxc(2)=-glij
153  cxc(3)=0d0
154  cxc(4)=dconjg(glij)
155  cxc(5)=0d0
156  cxc(6)=grij
157  cxc(7)=0d0
158  cxc(8)=-dconjg(grij)
159  s12min=0d0
160  s12max=(axmi-axmj)**2
161  IF( xxc(5).LT.axmi .OR. xxc(6).LT.axmi ) goto 130
162  IF(axmi.GE.axmj+2d0*pmas(1,1)) THEN
163  lknt=lknt+1
164  xlam(lknt)=c1*as/xmi3/(16d0*pi)*
165  & pygaus(pyxxz6,s12min,s12max,1d-2)
166  idlam(lknt,1)=kfnchi(ix)
167  idlam(lknt,2)=1
168  idlam(lknt,3)=-1
169  ENDIF
170  IF(axmi.GE.axmj+2d0*pmas(3,1)) THEN
171  lknt=lknt+1
172  xlam(lknt)=xlam(lknt-1)
173  idlam(lknt,1)=kfnchi(ix)
174  idlam(lknt,2)=3
175  idlam(lknt,3)=-3
176  ENDIF
177  130 CONTINUE
178  IF(axmi.GE.axmj+2d0*pmas(5,1)) THEN
179  pmold=pmas(pycomp(ksusy1+5),1)
180  IF(axmi.GT.pmas(pycomp(ksusy2+5),1)+pmas(5,1)) THEN
181  goto 140
182  ELSEIF(axmi.GT.pmas(pycomp(ksusy1+5),1)+pmas(5,1)) THEN
183  pmas(pycomp(ksusy1+5),1)=100d0*xmi
184  ENDIF
185  CALL pytbbn(ix,100,-1d0/3d0,xmi,gam)
186  lknt=lknt+1
187  xlam(lknt)=gam
188  idlam(lknt,1)=kfnchi(ix)
189  idlam(lknt,2)=5
190  idlam(lknt,3)=-5
191  pmas(pycomp(ksusy1+5),1)=pmold
192  ENDIF
193 C...U-TYPE QUARKS
194  140 CONTINUE
195  ia=2
196  xxc(5)=pmas(pycomp(ksusy1+ia),1)
197  xxc(6)=pmas(pycomp(ksusy2+ia),1)
198 C IF( XXC(5).LT.AXMI .OR. XXC(6).LT.AXMI ) GOTO 290
199  xxc(7)=xxc(5)
200  xxc(8)=xxc(6)
201  ei=kchg(ia,1)/3d0
202  t3i=sign(1d0,ei+1d-6)/2d0
203  glij=(t3i*zmixc(ix,2)-tanw*(t3i-ei)*zmixc(ix,1))*olpp
204  grij=zmixc(ix,1)*(ei*tanw)*orpp
205  cxc(2)=-glij
206  cxc(4)=dconjg(glij)
207  cxc(6)=grij
208  cxc(8)=-dconjg(grij)
209  IF( xxc(5).LT.axmi .OR. xxc(6).LT.axmi ) goto 150
210  IF(axmi.GE.axmj+2d0*pmas(2,1)) THEN
211  lknt=lknt+1
212  xlam(lknt)=c1*as/xmi3/(16d0*pi)*
213  & pygaus(pyxxz6,s12min,s12max,1d-2)
214  idlam(lknt,1)=kfnchi(ix)
215  idlam(lknt,2)=2
216  idlam(lknt,3)=-2
217  ENDIF
218  IF(axmi.GE.axmj+2d0*pmas(4,1)) THEN
219  lknt=lknt+1
220  xlam(lknt)=xlam(lknt-1)
221  idlam(lknt,1)=kfnchi(ix)
222  idlam(lknt,2)=4
223  idlam(lknt,3)=-4
224  ENDIF
225  150 CONTINUE
226 C...INCLUDE THE DECAY GLUINO -> NJ + T + T~
227 C...IF THE DECAY GLUINO -> ST + T CANNOT OCCUR
228  xmf=pmas(6,1)
229  IF(axmi.GE.axmj+2d0*xmf) THEN
230  pmold=pmas(pycomp(ksusy1+6),1)
231  IF(axmi.GT.pmas(pycomp(ksusy2+6),1)+xmf) THEN
232  goto 160
233  ELSEIF(axmi.GT.pmas(pycomp(ksusy1+6),1)+xmf) THEN
234  pmas(pycomp(ksusy1+6),1)=100d0*xmi
235  ENDIF
236  CALL pytbbn(ix,100,2d0/3d0,xmi,gam)
237  lknt=lknt+1
238  xlam(lknt)=gam
239  idlam(lknt,1)=kfnchi(ix)
240  idlam(lknt,2)=6
241  idlam(lknt,3)=-6
242  pmas(pycomp(ksusy1+6),1)=pmold
243  ENDIF
244  160 CONTINUE
245  ENDIF
246  170 CONTINUE
247 
248 C...GLUINO -> CI Q QBAR'
249  DO 210 ix=1,2
250  xmj=smw(ix)
251  axmj=abs(xmj)
252  IF(axmi.GE.axmj) THEN
253  DO 180 i=1,2
254  vmixc(ix,i)=dcmplx(vmix(ix,i),vmixi(ix,i))
255  umixc(ix,i)=dcmplx(umix(ix,i),umixi(ix,i))
256  180 CONTINUE
257  s12min=0d0
258  s12max=(axmi-axmj)**2
259  xxc(1)=0d0
260  xxc(2)=xmj
261  xxc(3)=0d0
262  xxc(4)=xmi
263  xxc(5)=pmas(pycomp(ksusy1+1),1)
264  xxc(6)=pmas(pycomp(ksusy1+2),1)
265  xxc(9)=1d6
266  xxc(10)=0d0
267  olpp=dcmplx(cos(rmss(32)),sin(rmss(32)))
268  orpp=dconjg(olpp)
269  cxc(1)=dcmplx(0d0,0d0)
270  cxc(3)=dcmplx(0d0,0d0)
271  cxc(5)=dcmplx(0d0,0d0)
272  cxc(7)=dcmplx(0d0,0d0)
273  cxc(2)=umixc(ix,1)*olpp/sr2
274  cxc(4)=-dconjg(vmixc(ix,1))*orpp/sr2
275  cxc(6)=dcmplx(0d0,0d0)
276  cxc(8)=dcmplx(0d0,0d0)
277  IF(xxc(5).LT.axmi) THEN
278  xxc(5)=1d6
279  ELSEIF(xxc(6).LT.axmi) THEN
280  xxc(6)=1d6
281  ENDIF
282  xxc(7)=xxc(6)
283  xxc(8)=xxc(5)
284  IF( xxc(5).LT.axmi .OR. xxc(6).LT.axmi ) goto 190
285  IF(axmi.GE.axmj+pmas(1,1)+pmas(2,1)) THEN
286  lknt=lknt+1
287  xlam(lknt)=0.5d0*c1*as/xmi3/(16d0*pi)*
288  & pygaus(pyxxz6,s12min,s12max,prec)
289  idlam(lknt,1)=kfcchi(ix)
290  idlam(lknt,2)=1
291  idlam(lknt,3)=-2
292  lknt=lknt+1
293  xlam(lknt)=xlam(lknt-1)
294  idlam(lknt,1)=-idlam(lknt-1,1)
295  idlam(lknt,2)=-idlam(lknt-1,2)
296  idlam(lknt,3)=-idlam(lknt-1,3)
297  ENDIF
298  IF(axmi.GE.axmj+pmas(3,1)+pmas(4,1)) THEN
299  lknt=lknt+1
300  xlam(lknt)=xlam(lknt-1)
301  idlam(lknt,1)=kfcchi(ix)
302  idlam(lknt,2)=3
303  idlam(lknt,3)=-4
304  lknt=lknt+1
305  xlam(lknt)=xlam(lknt-1)
306  idlam(lknt,1)=-idlam(lknt-1,1)
307  idlam(lknt,2)=-idlam(lknt-1,2)
308  idlam(lknt,3)=-idlam(lknt-1,3)
309  ENDIF
310  190 CONTINUE
311 
312  xmf=pmas(6,1)
313  xmfp=pmas(5,1)
314  IF(axmi.GE.axmj+xmf+xmfp) THEN
315  IF(xmi.GT.min(pmas(pycomp(ksusy1+5),1)+xmfp,
316  $ pmas(pycomp(ksusy2+6),1)+xmf)) goto 200
317  pmolt2=pmas(pycomp(ksusy2+6),1)
318  pmolb2=pmas(pycomp(ksusy2+5),1)
319  pmolt1=pmas(pycomp(ksusy1+6),1)
320  pmolb1=pmas(pycomp(ksusy1+5),1)
321  IF(xmi.GT.pmolt2+xmf) pmas(pycomp(ksusy2+6),1)=100d0*axmi
322  IF(xmi.GT.pmolt1+xmf) pmas(pycomp(ksusy1+6),1)=100d0*axmi
323  IF(xmi.GT.pmolb2+xmfp) pmas(pycomp(ksusy2+5),1)=100d0*axmi
324  IF(xmi.GT.pmolb1+xmfp) pmas(pycomp(ksusy1+5),1)=100d0*axmi
325  CALL pytbbc(ix,100,xmi,gam)
326  lknt=lknt+1
327  xlam(lknt)=gam
328  idlam(lknt,1)=kfcchi(ix)
329  idlam(lknt,2)=5
330  idlam(lknt,3)=-6
331  lknt=lknt+1
332  xlam(lknt)=xlam(lknt-1)
333  idlam(lknt,1)=-idlam(lknt-1,1)
334  idlam(lknt,2)=-idlam(lknt-1,2)
335  idlam(lknt,3)=-idlam(lknt-1,3)
336  pmas(pycomp(ksusy2+6),1)=pmolt2
337  pmas(pycomp(ksusy2+5),1)=pmolb2
338  pmas(pycomp(ksusy1+6),1)=pmolt1
339  pmas(pycomp(ksusy1+5),1)=pmolb1
340  ENDIF
341  200 CONTINUE
342  ENDIF
343  210 CONTINUE
344 
345 C...R-parity violating (3-body) decays.
346  CALL pyrvgl(kfin,xlam,idlam,lknt)
347 
348  iknt=lknt
349  xlam(0)=0d0
350  DO 220 i=1,iknt
351  IF(xlam(i).LT.0d0) xlam(i)=0d0
352  xlam(0)=xlam(0)+xlam(i)
353  220 CONTINUE
354  IF(xlam(0).EQ.0d0) xlam(0)=1d-6
355 
356  RETURN
357  END