]>
Commit | Line | Data |
---|---|---|
0795afa3 | 1 | #include "isajet/pilot.h" |
2 | SUBROUTINE WHIGGS | |
3 | C | |
4 | C Finish generation of whiggs events started bY TWOJET. | |
5 | C Select W decay modes as allowed by WMODE1, WMODE2. | |
6 | C Generate W decay angles and put vectors in PPAIR. | |
7 | C | |
8 | C | |
9 | #if defined(CERNLIB_IMPNONE) | |
10 | IMPLICIT NONE | |
11 | #endif | |
12 | #include "isajet/itapes.inc" | |
13 | #include "isajet/qcdpar.inc" | |
14 | #include "isajet/jetpar.inc" | |
15 | #include "isajet/primar.inc" | |
16 | #include "isajet/q1q2.inc" | |
17 | #include "isajet/jetsig.inc" | |
18 | #include "isajet/const.inc" | |
19 | #include "isajet/qsave.inc" | |
20 | #include "isajet/wcon.inc" | |
21 | #include "isajet/pjets.inc" | |
22 | #include "isajet/pinits.inc" | |
23 | #include "isajet/keys.inc" | |
24 | #include "isajet/hcon.inc" | |
25 | #include "isajet/wwpar.inc" | |
26 | #include "isajet/xmssm.inc" | |
27 | C | |
28 | DIMENSION X(2),LIST(25),P1WCM(4),P2WCM(4),P1LAB(4),P2LAB(4) | |
29 | 1,PBOOST(4) | |
30 | EQUIVALENCE (X(1),X1) | |
31 | DIMENSION JWWTYP(2) | |
32 | REAL GVQ(2),GAQ(2),GVL(2),GAL(2) | |
33 | REAL X,RND,RANF,CBRWW,AMASS,AM0,AM1,AM2, | |
34 | $E1CM,E2CM,P12CM,CTHCM,STHCM,PHICM,CPHICM,SPHICM,P1WCM,P2WCM, | |
35 | $PBOOST,P1LAB,P2LAB,ZHSIG,ZHMAX | |
36 | INTEGER JWWTYP,JET,JWT,JQ,IQ1,IQ2,LIST,NREJ,NJ0,K | |
37 | REAL BRANCH(12),SUMBR,BETAWH,GAMWH,PZWHCM,CTHD,WHSIG | |
38 | INTEGER IDABS,IDABS1,IDABSJ,IDIABS | |
39 | C | |
40 | DATA LIST/9,1,-1,2,-2,3,-3,4,-4,5,-5,6,-6, | |
41 | $11,-11,12,-12,13,-13,14,-14,15,-15,16,-16/ | |
42 | C | |
43 | GVQ(1)=.25-2*SIN2W/3. | |
44 | GVQ(2)=-.25+SIN2W/3. | |
45 | GAQ(1)=-.25 | |
46 | GAQ(2)=.25 | |
47 | GVL(1)=.25 | |
48 | GVL(2)=-.25+SIN2W | |
49 | GAL(1)=-.25 | |
50 | GAL(2)=.25 | |
51 | NPAIR=0 | |
52 | IF(KEYS(10).AND..NOT.GOMSSM) THEN | |
53 | JWWTYP(1)=JETTYP(1)-25 | |
54 | JWWTYP(2)=JETTYP(2)-25 | |
55 | ELSEIF(KEYS(10).AND.GOMSSM) THEN | |
56 | JWWTYP(1)=JETTYP(1)-76 | |
57 | JWWTYP(2)=JETTYP(2)-76 | |
58 | ENDIF | |
59 | C | |
60 | C Select W decay modes and put in /JETSET/. First particle | |
61 | C is always the fermion. | |
62 | ||
63 | DO 200 JET=1,2 | |
64 | IDABS=IABS(IDJETS(JET)) | |
65 | IF(IDABS.NE.80.AND.IDABS.NE.90) GO TO 200 | |
66 | RND=RANF() | |
67 | JWT=JWWTYP(JET) | |
68 | C Must only consider allowed decays for this mass | |
69 | SUMBR=0. | |
70 | DO 201 JQ=1,12 | |
71 | IQ1=2*JQ | |
72 | IQ2=MATCH(IQ1,JWT) | |
73 | IF(IQ2.EQ.0) THEN | |
74 | BRANCH(JQ)=0. | |
75 | GO TO 201 | |
76 | ENDIF | |
77 | AM1=AMASS(LIST(IQ1)) | |
78 | AM2=AMASS(LIST(IQ2)) | |
79 | IF(AM1+AM2.LT.PJETS(5,JET)) THEN | |
80 | BRANCH(JQ)=RBRWW(JQ,JWT,JET) | |
81 | SUMBR=SUMBR+BRANCH(JQ) | |
82 | ELSE | |
83 | BRANCH(JQ)=0. | |
84 | ENDIF | |
85 | 201 CONTINUE | |
86 | IF(SUMBR.LE.0.) GO TO 998 | |
87 | DO 202 JQ=1,12 | |
88 | 202 BRANCH(JQ)=BRANCH(JQ)/SUMBR | |
89 | C | |
90 | CBRWW=0. | |
91 | DO 210 JQ=1,12 | |
92 | CBRWW=CBRWW+BRANCH(JQ) | |
93 | IF(RND.GT.CBRWW) GO TO 210 | |
94 | IQ1=2*JQ | |
95 | IQ2=MATCH(IQ1,JWT) | |
96 | IDPAIR(NPAIR+1)=LIST(IQ1) | |
97 | IDPAIR(NPAIR+2)=LIST(IQ2) | |
98 | PPAIR(5,NPAIR+1)=AMASS(LIST(IQ1)) | |
99 | PPAIR(5,NPAIR+2)=AMASS(LIST(IQ2)) | |
100 | JPAIR(NPAIR+1)=JET | |
101 | JPAIR(NPAIR+2)=JET | |
102 | NPAIR=NPAIR+2 | |
103 | JQWW(JET)=JQ | |
104 | GO TO 200 | |
105 | 210 CONTINUE | |
106 | 200 CONTINUE | |
107 | C | |
108 | C Generate decay uniformly in angle and put in PPAIR. | |
109 | C Will check cross section later. | |
110 | C | |
111 | NREJ=0 | |
112 | 300 NJ0=2 | |
113 | DO 310 JET=1,2 | |
114 | IDABS=IABS(IDJETS(JET)) | |
115 | IF(IDABS.NE.80.AND.IDABS.NE.90) GO TO 310 | |
116 | C Construct W com momenta. | |
117 | IDABSJ=IDABS | |
118 | AM0=PJETS(5,JET) | |
119 | AM1=PPAIR(5,NJ0-1) | |
120 | AM2=PPAIR(5,NJ0) | |
121 | E1CM=(AM0**2+AM1**2-AM2**2)/(2.*AM0) | |
122 | E2CM=(AM0**2+AM2**2-AM1**2)/(2.*AM0) | |
123 | P12CM=(AM0**2-AM1**2-AM2**2)**2-4.*(AM1*AM2)**2 | |
124 | P12CM=SQRT(P12CM)/(2.*AM0) | |
125 | CTHCM=2.*RANF()-1. | |
126 | STHCM=SQRT(1.-CTHCM**2) | |
127 | PHICM=2.*PI*RANF() | |
128 | CPHICM=COS(PHICM) | |
129 | SPHICM=SIN(PHICM) | |
130 | P1WCM(1)=P12CM*STHCM*CPHICM | |
131 | P2WCM(1)=-P1WCM(1) | |
132 | P1WCM(2)=P12CM*STHCM*SPHICM | |
133 | P2WCM(2)=-P1WCM(2) | |
134 | P1WCM(3)=P12CM*CTHCM | |
135 | P2WCM(3)=-P1WCM(3) | |
136 | P1WCM(4)=E1CM | |
137 | P2WCM(4)=E2CM | |
138 | C Boost to lab frame. | |
139 | DO 320 K=1,3 | |
140 | 320 PBOOST(K)=-PJETS(K,JET) | |
141 | PBOOST(4)=PJETS(4,JET) | |
142 | CALL LBOOST(PBOOST,1,P1WCM,P1LAB) | |
143 | CALL LBOOST(PBOOST,1,P2WCM,P2LAB) | |
144 | DO 330 K=1,4 | |
145 | PPAIR(K,NJ0-1)=P1LAB(K) | |
146 | PPAIR(K,NJ0)=P2LAB(K) | |
147 | 330 CONTINUE | |
148 | NJ0=NJ0+2 | |
149 | 310 CONTINUE | |
150 | C | |
151 | C Impose simple (1+-cos(theta))**2 decay distribution for WH | |
152 | C Must use P1 in WH CoM frame | |
153 | IF (IDABSJ.NE.80.AND.IDABSJ.NE.90) GO TO 400 | |
154 | BETAWH=(PJETS(3,1)+PJETS(3,2))/(PJETS(4,1)+PJETS(4,2)) | |
155 | GAMWH=1./SQRT(1.-BETAWH**2) | |
156 | PZWHCM=GAMWH*(P1LAB(3)-BETAWH*P1LAB(4)) | |
157 | CTHD=PZWHCM/SQRT(P1LAB(1)**2+P1LAB(2)**2+PZWHCM**2) | |
158 | IF (IDINIT(1).LT.0) CTHD=-CTHD | |
159 | IDIABS=IABS(IDINIT(1)) | |
160 | IDABS1=IABS(IDPAIR(1)) | |
161 | IF (IDABSJ.EQ.80) THEN | |
162 | WHSIG=(1.+CTHD)**2 | |
163 | IF(WHSIG.GT.4*RANF()) GO TO 400 | |
164 | END IF | |
165 | IF (IDABSJ.EQ.90) THEN | |
166 | IF (IDIABS.EQ.1.OR.IDIABS.EQ.4) THEN | |
167 | IF (IDABS1.EQ.1.OR.IDABS1.EQ.4) THEN | |
168 | ZHSIG=(GVQ(1)**2+GAQ(1)**2)**2*(1.+CTHD**2) | |
169 | $ +8*GVQ(1)*GAQ(1)*GVQ(1)*GAQ(1)*CTHD | |
170 | ZHMAX=2*(GVQ(1)**2+GAQ(1)**2)**2 | |
171 | $ +8*GVQ(1)*GAQ(1)*GVQ(1)*GAQ(1) | |
172 | ELSEIF (IDABS1.EQ.2.OR.IDABS1.EQ.3.OR.IDABS1.EQ.5) THEN | |
173 | ZHSIG=(GVQ(1)**2+GAQ(1))*(GVQ(2)**2+GAQ(2)**2)*(1.+CTHD**2) | |
174 | $ +8*GVQ(1)*GAQ(1)*GVQ(2)*GAQ(2)*CTHD | |
175 | ZHMAX=(GVQ(1)**2+GAQ(1))*(GVQ(2)**2+GAQ(2)**2)*2 | |
176 | $ +8*GVQ(1)*GAQ(1)*GVQ(2)*GAQ(2) | |
177 | ELSEIF (IDABS1.EQ.11.OR.IDABS1.EQ.13.OR.IDABS1.EQ.15) THEN | |
178 | ZHSIG=(GVQ(1)**2+GAQ(1))*(GVL(1)**2+GAL(1)**2)*(1.+CTHD**2) | |
179 | $ +8*GVQ(1)*GAQ(1)*GVL(1)*GAL(1)*CTHD | |
180 | ZHMAX=(GVQ(1)**2+GAQ(1))*(GVL(1)**2+GAL(1)**2)*2 | |
181 | $ +8*GVQ(1)*GAQ(1)*GVL(1)*GAL(1) | |
182 | ELSEIF (IDABS1.EQ.12.OR.IDABS1.EQ.14.OR.IDABS1.EQ.16) THEN | |
183 | ZHSIG=(GVQ(1)**2+GAQ(1))*(GVL(2)**2+GAL(2)**2)*(1.+CTHD**2) | |
184 | $ +8*GVQ(1)*GAQ(1)*GVL(2)*GAL(2)*CTHD | |
185 | ZHMAX=(GVQ(1)**2+GAQ(1))*(GVL(2)**2+GAL(2)**2)*2 | |
186 | $ +8*GVQ(1)*GAQ(1)*GVL(2)*GAL(2) | |
187 | END IF | |
188 | ELSE IF (IDIABS.EQ.2.OR.IDIABS.EQ.3.OR.IDIABS.EQ.5) THEN | |
189 | IF (IDABS1.EQ.1.OR.IDABS1.EQ.4) THEN | |
190 | ZHSIG=(GVQ(2)**2+GAQ(2)**2)**2*(1.+CTHD**2) | |
191 | $ +8*GVQ(2)*GAQ(2)*GVQ(1)*GAQ(1)*CTHD | |
192 | ZHMAX=(GVQ(2)**2+GAQ(2)**2)**2*2 | |
193 | $ +8*GVQ(2)*GAQ(2)*GVQ(1)*GAQ(1) | |
194 | ELSEIF (IDABS1.EQ.2.OR.IDABS1.EQ.3.OR.IDABS1.EQ.5) THEN | |
195 | ZHSIG=(GVQ(2)**2+GAQ(2))*(GVQ(2)**2+GAQ(2)**2)*(1.+CTHD**2) | |
196 | $ +8*GVQ(2)*GAQ(2)*GVQ(2)*GAQ(2)*CTHD | |
197 | ZHMAX=(GVQ(2)**2+GAQ(2))*(GVQ(2)**2+GAQ(2)**2)*2 | |
198 | $ +8*GVQ(2)*GAQ(2)*GVQ(2)*GAQ(2) | |
199 | ELSEIF (IDABS1.EQ.11.OR.IDABS1.EQ.13.OR.IDABS1.EQ.15) THEN | |
200 | ZHSIG=(GVQ(2)**2+GAQ(2))*(GVL(1)**2+GAL(1)**2)*(1.+CTHD**2) | |
201 | $ +8*GVQ(2)*GAQ(2)*GVL(1)*GAL(1)*CTHD | |
202 | ZHMAX=(GVQ(2)**2+GAQ(2))*(GVL(1)**2+GAL(1)**2)*2 | |
203 | $ +8*GVQ(2)*GAQ(2)*GVL(1)*GAL(1) | |
204 | ELSEIF (IDABS1.EQ.12.OR.IDABS1.EQ.14.OR.IDABS1.EQ.16) THEN | |
205 | ZHSIG=(GVQ(2)**2+GAQ(2))*(GVL(2)**2+GAL(2)**2)*(1.+CTHD**2) | |
206 | $ +8*GVQ(2)*GAQ(2)*GVL(2)*GAL(2)*CTHD | |
207 | ZHMAX=(GVQ(2)**2+GAQ(2))*(GVL(2)**2+GAL(2)**2)*2 | |
208 | $ +8*GVQ(2)*GAQ(2)*GVL(2)*GAL(2) | |
209 | END IF | |
210 | END IF | |
211 | IF(ZHSIG.GT.RANF()*ZHMAX) GO TO 400 | |
212 | END IF | |
213 | NREJ=NREJ+1 | |
214 | IF(NREJ.LT.NTRIES) GO TO 300 | |
215 | GO TO 999 | |
216 | C | |
217 | C Good event | |
218 | C | |
219 | 400 CONTINUE | |
220 | RETURN | |
221 | C | |
222 | 999 CALL PRTEVT(0) | |
223 | WRITE(ITLIS,9991) NREJ | |
224 | 9991 FORMAT(//' IT IS TAKING MORE THAN',I5,' TRIES TO GENERATE ', | |
225 | 1'A GOOD WHIGGS EVENT.'/' CHECK LIMITS OR INCREASE NTRIES.') | |
226 | STOP 99 | |
227 | 998 CALL PRTEVT(0) | |
228 | WRITE(ITLIS,9981) JET | |
229 | 9981 FORMAT(//' ERROR IN WHIGGS ... NO DECAY POSSIBLE FOR JET',I3) | |
230 | STOP 99 | |
231 | END |