]>
Commit | Line | Data |
---|---|---|
1 | #include "isajet/pilot.h" | |
2 | SUBROUTINE SIGDY | |
3 | C | |
4 | C Compute the Drell-Yan and Drell-Yan plus jet cross sections | |
5 | C d(sigma)/d(qmw**2)d(qtw**2)d(yw)d(yj) | |
6 | C | |
7 | C SIGMA = cross section summed over quark types allowed by | |
8 | C JETTYPE3 and WTYPE cards. | |
9 | C SIGS(I) = partial cross section for I1 + I2 --> I3 + I4. | |
10 | C INOUT(I) = IOPAK**3*I4 + IOPAK**2*I3 + IOPAK*I2 + I1 | |
11 | C using JETTYPE code. | |
12 | C | |
13 | C QT cutoff for W+JET taken from Parisi and Petronzio, | |
14 | C Nucl Phys B154, 427 | |
15 | C qk + gl --> qk + w suppressed at low QTW by extra factor | |
16 | C of qtw**2/(qtw**2+qt2cut(qmw)) | |
17 | C | |
18 | C Ver 7.17: include top mass for gb --> Wt and gt --> Zt | |
19 | C with no extra qt suppression factor. Note we do NOT include | |
20 | C gt --> Wb; while this process makes sense for qt >> m_t, | |
21 | C it has a pole in the physical region at low qt from the | |
22 | C on-shell decay t --> Wb. We let Q**2 --> Q**2 + m_t**2 | |
23 | C in the scale for the parton distributions. | |
24 | C | |
25 | C Ver 7.32: Rewrite AJLWT for gb --> Wt, etc., in terms of | |
26 | C scaled variables, and restore SWT**5 later to avoid | |
27 | C floating errors on VMS. | |
28 | C | |
29 | C Ver 7.41: Recalculate COUT for each mass(!). | |
30 | C | |
31 | #if defined(CERNLIB_IMPNONE) | |
32 | IMPLICIT NONE | |
33 | #endif | |
34 | #include "isajet/itapes.inc" | |
35 | #include "isajet/qcdpar.inc" | |
36 | #include "isajet/jetpar.inc" | |
37 | #include "isajet/primar.inc" | |
38 | #include "isajet/q1q2.inc" | |
39 | #include "isajet/jetsig.inc" | |
40 | #include "isajet/qsave.inc" | |
41 | #include "isajet/wcon.inc" | |
42 | #include "isajet/const.inc" | |
43 | #include "isajet/nodcay.inc" | |
44 | C | |
45 | REAL X(2) | |
46 | REAL Z,S,T,U,QMW2,QZW,EHAT,Q2SAVE,YHAT,EY,P3Z,P1,P2,AMASS,ANEFF, | |
47 | $SIG0,DENOM,QT2CUT,SIGT,SIGU,FAC,PROP,FACTOR,SIG,AMT,AMT2,SWT, | |
48 | $P1WT,P2WT,X1WT,X2WT,TWT,UWT,Q2,QFCN,STRUC,XX,ACOSH,ATANH,P2M,P1M | |
49 | REAL AMI2,AMF2,EFWT | |
50 | REAL AJLWT,AJLZT1,AJLZT2,A2,A2B2,QQ,TM2 | |
51 | INTEGER I,IQ,IH,IQ1,IFL,IQ2,IW | |
52 | INTEGER NZERO(4) | |
53 | REAL AMFAC(13) | |
54 | INTEGER NUTYP(25) | |
55 | INTEGER IFL1,IFL2 | |
56 | REAL TERM | |
57 | EQUIVALENCE (S,SHAT),(T,THAT),(U,UHAT),(X(1),X1) | |
58 | C | |
59 | DATA NZERO/11,9,9,11/ | |
60 | DATA AMFAC/11*0.,2*1./ | |
61 | DATA NUTYP/13*0,1,1,0,0,1,1,0,0,1,1,0,0/ | |
62 | C | |
63 | C Functions | |
64 | ACOSH(Z)=ALOG(Z+SQRT(Z**2-1.)) | |
65 | ATANH(Z)=.5*ALOG((1.+Z)/(1.-Z)) | |
66 | PROP(I)=1./((QMW2-WMASS(I)**2)**2+(WMASS(I)*WGAM(I))**2) | |
67 | C Qt cutoff function | |
68 | QT2CUT(QMW)=CUTOFF*QMW**CUTPOW | |
69 | C Parton distributions | |
70 | QFCN(XX,IQ,IH)=STRUC(XX,QSQ+AMT2,IQ,IDIN(IH))/XX | |
71 | C Integrated matrix elements JLint from FORM | |
72 | AJLWT(S,T,QQ,TM2)= | |
73 | $ - 32*QQ**3*S*T + 32*QQ**3*S*TM2 + 32*QQ**2*S**2*T | |
74 | $ + 32*QQ**2*S*T**2 - 16*QQ**2*S*T*TM2 - 16*QQ**2*S*TM2**2 | |
75 | $ - 16*QQ*S**3*T + 16*QQ*S**3*TM2 - 16*QQ*S**2*T*TM2 | |
76 | $ - 16*QQ*S*T**3 + 32*QQ*S*T**2*TM2 - 16*QQ*S*T*TM2**2 | |
77 | $ - 8*S**3*T*TM2 + 8*S**3*TM2**2 - 16*S**2*T**2*TM2 | |
78 | $ + 16*S**2*T*TM2**2 - 16*S**2*TM2**3 - 8*S*T**3*TM2 | |
79 | $ + 8*S*T**2*TM2**2 - 8*S*T*TM2**3 + 8*S*TM2**4 | |
80 | C | |
81 | AJLZT1(S,T,QQ,TM2)= | |
82 | $ + A2 * ( - 96*QQ**2*S*T*TM2 + 96*QQ**2*S*TM2**2 | |
83 | $ + 96*QQ**2*T*TM2**2 - 96*QQ**2*TM2**3 + 96*QQ*S**2*T*TM2 | |
84 | $ + 96*QQ*S*T**2*TM2 - 192*QQ*S*T*TM2**2 - 96*QQ*S*TM2**3 | |
85 | $ - 96*QQ*T*TM2**3 + 192*QQ*TM2**4 + 16*S**3*T*TM2 | |
86 | $ - 16*S**3*TM2**2 + 32*S**2*T**2*TM2 - 112*S**2*T*TM2**2 | |
87 | $ + 80*S**2*TM2**3 + 16*S*T**3*TM2 - 112*S*T**2*TM2**2 | |
88 | $ + 224*S*T*TM2**3 - 128*S*TM2**4 - 16*T**3*TM2**2 | |
89 | $ + 80*T**2*TM2**3 - 128*T*TM2**4 + 64*TM2**5 ) | |
90 | AJLZT2(S,T,QQ,TM2)= | |
91 | $ + A2B2 * ( - 16*QQ**3*S*T + 16*QQ**3*S*TM2 + 16*QQ**3*T*TM2 | |
92 | $ - 16*QQ**3*TM2**2 + 16*QQ**2*S**2*T + 16*QQ**2*S*T**2 | |
93 | $ + 32*QQ**2*S*T*TM2 - 80*QQ**2*S*TM2**2 - 80*QQ**2*T*TM2**2 | |
94 | $ + 96*QQ**2*TM2**3 - 8*QQ*S**3*T + 8*QQ*S**3*TM2 | |
95 | $ - 40*QQ*S**2*T*TM2 - 24*QQ*S**2*TM2**2 - 8*QQ*S*T**3 | |
96 | $ - 40*QQ*S*T**2*TM2 + 80*QQ*S*T*TM2**2 + 96*QQ*S*TM2**3 | |
97 | $ + 8*QQ*T**3*TM2 - 24*QQ*T**2*TM2**2 + 96*QQ*T*TM2**3 | |
98 | $ - 144*QQ*TM2**4 - 16*S**3*T*TM2 + 16*S**3*TM2**2 | |
99 | $ - 32*S**2*T**2*TM2 + 112*S**2*T*TM2**2 - 80*S**2*TM2**3 | |
100 | $ - 16*S*T**3*TM2 + 112*S*T**2*TM2**2 - 224*S*T*TM2**3 | |
101 | $ + 128*S*TM2**4 + 16*T**3*TM2**2 - 80*T**2*TM2**3 | |
102 | $ + 128*T*TM2**4 - 64*TM2**5 ) | |
103 | C | |
104 | C Kinematics | |
105 | C | |
106 | QMW2=QMW**2 | |
107 | QTMW=SQRT(QMW2+QTW**2) | |
108 | Q0W=QTMW*COSH(YW) | |
109 | QZW=QTMW*SINH(YW) | |
110 | QW=SQRT(QZW**2+QTW**2) | |
111 | C Protect against errors | |
112 | IF(QW.NE.0.) THEN | |
113 | CTHW=QZW/QW | |
114 | STHW=QTW/QW | |
115 | IF(ABS(CTHW).LT.1.) THEN | |
116 | THW=ACOS(CTHW) | |
117 | ELSE | |
118 | CTHW=0. | |
119 | STHW=1. | |
120 | THW=.5*PI | |
121 | ENDIF | |
122 | ELSE | |
123 | CTHW=0. | |
124 | STHW=1. | |
125 | THW=.5*PI | |
126 | ENDIF | |
127 | C | |
128 | IF(STDDY) THEN | |
129 | C Kinematics for standard Drell-Yan | |
130 | EHAT=QMW | |
131 | SHAT=QMW**2 | |
132 | QSQ=SHAT | |
133 | Q2SAVE=QSQ | |
134 | YHAT=YW | |
135 | EY=EXP(YHAT) | |
136 | X1=EHAT/ECM*EY | |
137 | X2=EHAT/(ECM*EY) | |
138 | ELSE | |
139 | C Kinematics for Drell-Yan plus jet | |
140 | P3Z=P(3)*CTH(3) | |
141 | SHAT=QMW2+2.*Q0W*P(3)-2.*QZW*P3Z+2.*PT(3)**2 | |
142 | P1=.5*(P(3)+P3Z+Q0W+QZW) | |
143 | P2=.5*(P(3)-P3Z+Q0W-QZW) | |
144 | X1=P1/HALFE | |
145 | X2=P2/HALFE | |
146 | THAT=-2.*P1*(P(3)-P3Z) | |
147 | UHAT=-2.*P2*(P(3)+P3Z) | |
148 | QSQ=QTW**2 | |
149 | QSQ=AMAX1(QSQ,4.) | |
150 | ANEFF=4.+QSQ/(QSQ+AMASS(5)**2)+QSQ/(QSQ+AMASS(6)**2) | |
151 | ALFQSQ=12.*PI/((33.-2.*ANEFF)*ALOG(QSQ/ALAM2)) | |
152 | Q2SAVE=QSQ | |
153 | QSQ=SHAT | |
154 | ENDIF | |
155 | C | |
156 | C Initialize | |
157 | C | |
158 | SIGMA=0. | |
159 | NSIGS=0 | |
160 | DO 100 I=1,MXSIGS | |
161 | SIGS(I)=0. | |
162 | 100 CONTINUE | |
163 | IF(X1.GE.1..OR.X2.GE.1.) RETURN | |
164 | C | |
165 | C Compute structure functions | |
166 | C | |
167 | DO 110 IH=1,2 | |
168 | DO 120 IQ=1,11 | |
169 | QSAVE(IQ,IH)=STRUC(X(IH),QSQ,IQ,IDIN(IH))/X(IH) | |
170 | 120 CONTINUE | |
171 | QSAVE(12,IH)=0 | |
172 | QSAVE(13,IH)=0 | |
173 | 110 CONTINUE | |
174 | QSQ=Q2SAVE | |
175 | C | |
176 | C Recompute COUT for this mass | |
177 | C | |
178 | DO 130 IW=1,4 | |
179 | COUT(IW)=0. | |
180 | IF(.NOT.GODY(IW)) GO TO 130 | |
181 | DO 140 IQ1=2,25 | |
182 | IQ2=MATCH(IQ1,IW) | |
183 | IF(IQ2.EQ.0) GO TO 140 | |
184 | IF(.NOT.(GOQ(IQ1,1).AND.GOQ(IQ2,2))) GO TO 140 | |
185 | IF(NUTYP(IQ1)*NUTYP(IQ2).EQ.1.AND.NONUNU) GO TO 140 | |
186 | IFL1=IQ1/2 | |
187 | IFL2=IQ2/2 | |
188 | IF(AMASS(IFL1)+AMASS(IFL2).GE.QMW) GO TO 140 | |
189 | TERM=.5*(AQ(IFL1,IW)**2+BQ(IFL1,IW)**2) | |
190 | IF(IQ1.LE.13) TERM=3.*TERM | |
191 | COUT(IW)=COUT(IW)+TERM | |
192 | 140 CONTINUE | |
193 | 130 CONTINUE | |
194 | C | |
195 | IF(STDDY) GO TO 400 | |
196 | C | |
197 | C Compute cross section for types allowed by WTYPE and | |
198 | C JETTYPE cards. | |
199 | C | |
200 | C qk + gl --> qk + W | |
201 | C | |
202 | SIG0=ALFA**2*ALFQSQ*QMW2/(9.*SCM*S)*UNITS | |
203 | DENOM=S**2*EXP(.5*ALOG(QTW**4+QT2CUT(QMW)**2)) | |
204 | SIGT=SIG0*(S**2+U**2+2.*T*QMW2)*(-T)/DENOM | |
205 | SIGU=SIG0*(S**2+T**2+2.*U*QMW2)*(-U)/DENOM | |
206 | DO 200 IW=1,4 | |
207 | IF(.NOT.GODY(IW)) GO TO 200 | |
208 | FAC=COUT(IW)*PROP(IW) | |
209 | DO 210 IQ=2,NZERO(IW) | |
210 | IF(.NOT.GOQ(IQ,3)) GO TO 210 | |
211 | IQ1=MATCH(IQ,4) | |
212 | IQ1=MATCH(IQ1,IW) | |
213 | IF(IQ1.EQ.0) GO TO 210 | |
214 | IFL=IQ/2 | |
215 | FACTOR=FAC*(AQ(IFL,IW)**2+BQ(IFL,IW)**2) | |
216 | $ *QTW**2/(QTW**2+QT2CUT(QMW)) | |
217 | SIG=FACTOR*SIGT*QSAVE(IQ1,1)*QSAVE(1,2) | |
218 | CALL SIGFIL(SIG,IQ1,1,IW,IQ) | |
219 | SIG=FACTOR*SIGU*QSAVE(IQ1,2)*QSAVE(1,1) | |
220 | CALL SIGFIL(SIG,1,IQ1,IW,IQ) | |
221 | 210 CONTINUE | |
222 | 200 CONTINUE | |
223 | C | |
224 | C bt,tp + gl -> bt,tp + W,Z | |
225 | C | |
226 | AMT=AMASS(6) | |
227 | AMT2=AMT**2 | |
228 | Q2=QMW2 | |
229 | DO 220 IW=2,4 | |
230 | IF(.NOT.GODY(IW)) GO TO 220 | |
231 | DO 230 IQ=NZERO(IW)+1,13 | |
232 | IF(.NOT.GOQ(IQ,3)) GO TO 230 | |
233 | IQ1=MATCH(IQ,4) | |
234 | IQ1=MATCH(IQ1,IW) | |
235 | IF(IQ1.EQ.0) GO TO 230 | |
236 | IF(IQ1.GE.12.AND.IW.NE.4) GO TO 230 | |
237 | C Assign zero or top masses for initial/final quarks | |
238 | AMF2=AMT2*AMFAC(IQ) | |
239 | AMI2=AMT2*AMFAC(IQ1) | |
240 | EFWT=SQRT(P(3)**2+AMF2) | |
241 | SWT=QMW2+AMF2+2.*Q0W*EFWT-2.*QZW*P3Z+2.*PT(3)**2 | |
242 | C | |
243 | C qk + gl initial state | |
244 | C Do kinematics using p(small) = 0 for gluon | |
245 | C | |
246 | P1WT=EFWT+P3Z+Q0W+QZW | |
247 | P1M=AMI2/P1WT | |
248 | P2WT=EFWT-P3Z+Q0W-QZW-P1M | |
249 | X1WT=.5*P1WT/HALFE | |
250 | X2WT=.5*P2WT/HALFE | |
251 | TWT=-P1WT*(EFWT-P3Z)-P1M*(EFWT+P3Z)+AMI2+AMF2 | |
252 | UWT=-P2WT*(EFWT+P3Z)+AMF2 | |
253 | IF(X1WT.LT.0.OR.X1WT.GT.1.OR.X2WT.LT.0.OR.X2WT.GT.1) | |
254 | $ GO TO 240 | |
255 | C Cross sections | |
256 | IF(IW.EQ.2.OR.IW.EQ.3) THEN | |
257 | SIG0=ALFA**2*ALFQSQ/(144*SCM*SWT)*UNITS | |
258 | SIG0=SIG0*(AQ(5,IW)**2+BQ(5,IW)**2)*COUT(IW)*PROP(IW) | |
259 | SIGU=SIG0*AJLWT(SWT/SWT,UWT/SWT,Q2/SWT,AMT2/SWT)*SWT* | |
260 | $ (SWT/(SWT-AMI2))**2*(SWT/(UWT-AMF2))**2 | |
261 | SIG=SIGU*QFCN(X1WT,IQ1,1)*QFCN(X2WT,1,2) | |
262 | CALL SIGFIL(SIG,IQ1,1,IW,IQ) | |
263 | ELSEIF(IW.EQ.4) THEN | |
264 | SIG0=ALFA**2*ALFQSQ/(144*SCM*SWT)*UNITS | |
265 | SIG0=SIG0*COUT(IW)*PROP(IW) | |
266 | A2=AQ(6,IW)**2 | |
267 | A2B2=AQ(6,IW)**2+BQ(6,IW)**2 | |
268 | SIGU=SIG0*(AJLZT1(SWT/SWT,UWT/SWT,Q2/SWT,AMT2/SWT)+ | |
269 | $ AJLZT2(SWT/SWT,UWT/SWT,Q2/SWT,AMT2/SWT))*SWT* | |
270 | $ (SWT/(SWT-AMI2))**2*(SWT/(UWT-AMF2))**2 | |
271 | SIG=SIGU*QFCN(X1WT,IQ1,1)*QFCN(X2WT,1,2) | |
272 | CALL SIGFIL(SIG,IQ1,1,IW,IQ) | |
273 | ENDIF | |
274 | C | |
275 | C gl + qk initial state | |
276 | C Do kinematics using p(small) = 0 for gluon | |
277 | C | |
278 | 240 P2WT=EFWT-P3Z+Q0W-QZW | |
279 | P2M=AMI2/P2WT | |
280 | P1WT=EFWT+P3Z+Q0W+QZW-P2M | |
281 | X1WT=.5*P1WT/HALFE | |
282 | X2WT=.5*P2WT/HALFE | |
283 | TWT=-P1WT*(EFWT-P3Z)+AMF2 | |
284 | UWT=-P2WT*(EFWT+P3Z)-P2M*(EFWT-P3Z)+AMI2+AMF2 | |
285 | IF(X1WT.LT.0.OR.X1WT.GT.1.OR.X2WT.LT.0.OR.X2WT.GT.1) | |
286 | $ GO TO 230 | |
287 | C Cross sections | |
288 | IF(IW.EQ.2.OR.IW.EQ.3) THEN | |
289 | SIG0=ALFA**2*ALFQSQ/(144*SCM*SWT)*UNITS | |
290 | SIG0=SIG0*(AQ(5,IW)**2+BQ(5,IW)**2)*COUT(IW)*PROP(IW) | |
291 | SIGT=SIG0*AJLWT(SWT/SWT,TWT/SWT,Q2/SWT,AMT2/SWT)*SWT* | |
292 | $ (SWT/(SWT-AMI2))**2*(SWT/(TWT-AMF2)**2) | |
293 | SIG=SIGT*QFCN(X1WT,1,1)*QFCN(X2WT,IQ1,2) | |
294 | CALL SIGFIL(SIG,1,IQ1,IW,IQ) | |
295 | ELSEIF(IW.EQ.4) THEN | |
296 | SIG0=ALFA**2*ALFQSQ/(144*SCM*SWT)*UNITS | |
297 | SIG0=SIG0*COUT(IW)*PROP(IW) | |
298 | A2=AQ(6,IW)**2 | |
299 | A2B2=AQ(6,IW)**2+BQ(6,IW)**2 | |
300 | SIGU=SIG0*(AJLZT1(SWT/SWT,TWT/SWT,Q2/SWT,AMT2/SWT)+ | |
301 | $ AJLZT2(SWT/SWT,TWT/SWT,Q2/SWT,AMT2/SWT))*SWT* | |
302 | $ (SWT/(SWT-AMI2))**2*(SWT/(TWT-AMF2))**2 | |
303 | SIG=SIGU*QFCN(X1WT,1,1)*QFCN(X2WT,IQ1,2) | |
304 | CALL SIGFIL(SIG,1,IQ1,IW,IQ) | |
305 | ENDIF | |
306 | 230 CONTINUE | |
307 | 220 CONTINUE | |
308 | C | |
309 | C qk + qb --> gl + W | |
310 | C | |
311 | IF(.NOT.GOQ(1,3)) RETURN | |
312 | SIG0=8.*ALFA**2*ALFQSQ*QMW2/(27.*SCM*S)*UNITS | |
313 | DENOM=S*EXP(.5*ALOG(QTW**4+QT2CUT(QMW)**2)) | |
314 | SIG0=SIG0*(T**2+U**2+2.*S*QMW2)/DENOM | |
315 | DO 300 IW=1,4 | |
316 | IF(.NOT.GODY(IW)) GO TO 300 | |
317 | FAC=COUT(IW)*PROP(IW) | |
318 | DO 310 IQ1=2,11 | |
319 | IQ2=MATCH(IQ1,IW) | |
320 | IF(IQ2.EQ.0) GO TO 310 | |
321 | IFL=IQ1/2 | |
322 | SIG=FAC*SIG0*(AQ(IFL,IW)**2+BQ(IFL,IW)**2) | |
323 | $ *QSAVE(IQ1,1)*QSAVE(IQ2,2) | |
324 | CALL SIGFIL(SIG,IQ1,IQ2,IW,1) | |
325 | 310 CONTINUE | |
326 | 300 CONTINUE | |
327 | RETURN | |
328 | C | |
329 | C Standard Drell-Yan for QT=0. | |
330 | C | |
331 | 400 CONTINUE | |
332 | SIG0=4.*PI*ALFA**2*QMW2/(9.*SCM)*UNITS | |
333 | DO 410 IW=1,4 | |
334 | IF(.NOT.GODY(IW)) GO TO 410 | |
335 | FAC=COUT(IW)*PROP(IW) | |
336 | DO 420 IQ1=2,13 | |
337 | IQ2=MATCH(IQ1,IW) | |
338 | IF(IQ2.EQ.0) GO TO 420 | |
339 | IFL=IQ1/2 | |
340 | SIG=FAC*SIG0*(AQ(IFL,IW)**2+BQ(IFL,IW)**2) | |
341 | $ *QSAVE(IQ1,1)*QSAVE(IQ2,2) | |
342 | CALL SIGFIL(SIG,IQ1,IQ2,IW,0) | |
343 | 420 CONTINUE | |
344 | 410 CONTINUE | |
345 | C | |
346 | RETURN | |
347 | END |