More volume overlaps corrected
[u/mrichter/AliRoot.git] / ISAJET / code / sigh.F
1 #include "isajet/pilot.h"
2       SUBROUTINE SIGH
3 C
4 C          COMPUTE THE INTEGRATED WEINBERG-SALAM HIGGS CROSS SECTION
5 C          D(SIGMA)/D(QMW**2)D(YW)
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          VER. 7.14: CHECK INITIAL QUARK MASS IS ALLOWED
14 C
15 #include "isajet/itapes.inc"
16 #include "isajet/qcdpar.inc"
17 #include "isajet/jetpar.inc"
18 #include "isajet/primar.inc"
19 #include "isajet/q1q2.inc"
20 #include "isajet/jetsig.inc"
21 #include "isajet/qsave.inc"
22 #include "isajet/wcon.inc"
23 #include "isajet/const.inc"
24 #include "isajet/jetlim.inc"
25 #include "isajet/hcon.inc"
26 C
27       DIMENSION AMQCUR(6),LISTW(4),WTHELI(4),FINT(9)
28       DIMENSION X(2)
29       EQUIVALENCE (S,SHAT),(T,THAT),(U,UHAT),(X(1),X1)
30 #if defined(CERNLIB_DOUBLE)
31       DOUBLE PRECISION C,TERM,SUM,FINT,ZLIM
32 #endif
33       DATA AMQCUR/.005,.009,.175,1.25,4.50,30./
34       DATA LISTW/10,80,-80,90/
35 C          WTHELI ARE WEIGHTS OF HELICITY AMPLITUDES IN SIGMA.
36       DATA WTHELI/1.,2.,2.,4./
37 C
38 C          FUNCTIONS
39       ACOSH(Z)=ALOG(Z+SQRT(Z**2-1.))
40       ATANH(Z)=.5*ALOG((1.+Z)/(1.-Z))
41 C
42 C          KINEMATICS (IDENTICAL TO DRELL-YAN)
43 C
44       AMQCUR(6)=AMASS(6)
45       QMW2=QMW**2
46       QTMW=SQRT(QMW2+QTW**2)
47       Q0W=QTMW*COSH(YW)
48       QZW=QTMW*SINH(YW)
49       QW=SQRT(QZW**2+QTW**2)
50       IF(QW.NE.0.) THEN
51         CTHW=QZW/QW
52         STHW=QTW/QW
53         IF(ABS(CTHW).LT.1.) THEN
54           THW=ACOS(CTHW)
55         ELSE
56           CTHW=0.
57           STHW=1.
58           THW=.5*PI
59         ENDIF
60       ELSE
61         CTHW=0.
62         STHW=1.
63         THW=.5*PI
64       ENDIF
65       EHAT=QMW
66       SHAT=QMW**2
67       QSQ=SHAT
68       ANEFF=4.+QSQ/(QSQ+AMASS(5)**2)+QSQ/(QSQ+AMASS(6)**2)
69       ALFQSQ=12.*PI/((33.-ANEFF)*ALOG(QSQ/ALAM2))
70       Q2SAVE=QSQ
71       YHAT=YW
72       EY=EXP(YHAT)
73       X1=EHAT/ECM*EY
74       X2=EHAT/(ECM*EY)
75 C
76 C          INITIALIZE
77 C
78       SIGMA=0.
79       NSIGS=0
80       DO 100 I=1,MXSIGS
81 100   SIGS(I)=0
82 C
83       IF(X1.GE.1..OR.X2.GE.1.) RETURN
84 C
85 C          COMPUTE STRUCTURE FUNCTIONS
86       DO 110 IH=1,2
87       DO 120 IQ=1,13
88 120   QSAVE(IQ,IH)=STRUC(X(IH),QSQ,IQ,IDIN(IH))/X(IH)
89       DO 130 IQ=14,26
90 130   QSAVE(IQ,IH)=0.
91       DO 140 IW=2,4
92       AMW=AMASS(LISTW(IW))
93       IF(QMW.GT.2.*AMW) THEN
94         QSAVE(25+IW,IH)=STRUCW(X(IH),IW,IDIN(IH))/X(IH)
95       ELSE
96         QSAVE(25+IW,IH)=0.
97       ENDIF
98 140   CONTINUE
99 110   CONTINUE
100 C
101 C          CALCULATE HIGGS-GLUON-GLUON COUPLING FOR GIVEN Q**2
102       ETAR=0.
103       ETAI=0.
104       DO 150 IQ=1,8
105       AMQ=AMASS(IQ)
106       IF(AMQ.LE.0.) GO TO 150
107       RQ=(2.*AMQ/HMASS)**2
108       IF(RQ.GE.1.) THEN
109         ETAR=ETAR+.5*RQ*(1.+(1.-RQ)*ASIN(1./SQRT(RQ))**2)
110       ELSE
111         RQLOG=ALOG((1.+SQRT(1.-RQ))/(1.-SQRT(1.-RQ)))
112         PHIR=.25*(RQLOG**2-PI**2)
113         ETAR=ETAR+.5*RQ*(1.+(RQ-1.)*PHIR)
114         PHII=.5*PI*RQLOG
115         ETAI=ETAI+.5*RQ*(RQ-1.)*PHII
116       ENDIF
117 150   CONTINUE
118       ETAHGG=ETAR**2+ETAI**2
119 C
120 C          GL + GL --> HIGGS
121 C
122       SIG0=GF*ALFQSQ**2/(32.*PI*SQRT2)*ETAHGG*X1*X2*UNITS
123       SIG0=SIG0*S/(PI*HMASS*((S-HMASS**2)**2+(HMASS*HGAM)**2))
124       SIG0=SIG0*QSAVE(1,1)*QSAVE(1,2)
125       DO 160 IQ1=2,29
126       IQ2=MATCHH(IQ1)
127       IF(GOQ(IQ1,1).AND.GOQ(IQ2,2)) THEN
128         SIG=SIG0*HGAMS(IQ1)
129         IF(IQ1.GT.25) SIG=SIG*TBRWW(IQ1-25,1)*TBRWW(IQ2-25,2)
130         CALL SIGFIL(SIG,1,1,IQ1,IQ2)
131       ENDIF
132 160   CONTINUE
133 C
134 C          QK + QB --> HIGGS
135 C
136       SIG0=PI*GF/(3.*SQRT2*HMASS**2)*X1*X2*UNITS
137       SIG0=SIG0*S/(PI*HMASS*((S-HMASS**2)**2+(HMASS*HGAM)**2))
138       DO 210 IQ1=2,13
139       IQ2=MATCHH(IQ1)
140       AMQ=AMQCUR(IQ1/2)
141       IF(QMW.LE.2*AMQ) GO TO 210
142       SIG1=SIG0*AMQ**2*QSAVE(IQ1,1)*QSAVE(IQ2,2)
143       DO 220 IQ3=2,29
144       IQ4=MATCHH(IQ3)
145       IF(GOQ(IQ3,1).AND.GOQ(IQ4,2)) THEN
146         SIG=SIG1*HGAMS(IQ3)
147         IF(IQ3.GT.25) SIG=SIG*TBRWW(IQ3-25,1)*TBRWW(IQ4-25,2)
148         CALL SIGFIL(SIG,IQ1,IQ2,IQ3,IQ4)
149       ENDIF
150 220   CONTINUE
151 210   CONTINUE
152 C
153 C          W+W FUSION AND W+W->W+W IN EFFECTIVE W APPROXIMATION WITH
154 C          ANGULAR DISTRIBUTION CUT OFF BY PTMIN.
155 C          Z0 Z0 FINAL STATE HAS SYMMETRY FACTOR OF .5
156 C
157       IF(QMW.LE.2.*AMASS(80)) GO TO 500
158 C
159 C          W+ W- --> W+ W-
160 C
161       IF(.NOT.((GOQ(27,1).AND.GOQ(28,2)).OR.(GOQ(28,1).AND.GOQ(27,2))))
162      $GO TO 400
163       WM=AMASS(80)
164       PWWCM=.5*SQRT(QMW**2-4.*WM**2)
165       STHLIM=PTMIN(1)/PWWCM
166       IF(STHLIM.LE.1) THEN
167         ZLIM=SQRT(1.-STHLIM**2)
168       ELSE
169         GO TO 400
170       ENDIF
171 C          SET UP AMPLITUDES
172       CALL XWWWW
173 C          SUM CROSS SECTION TERMS. I,J RUN OVER AMPLITUDE TERMS.
174 C          L RUNS OVER HELICITY STATES. N RUNS OVER POWERS.
175 C          REMEMBER THAT L=4 IS MISSING SIN(THETA)/SQRT(2)
176       SUM=0.
177       DO 311 I=1,4
178       DO 311 J=I,4
179       CALL SIGINT(FINT,ZLIM,ADWWWW(1,I),ADWWWW(2,I),ADWWWW(1,J),
180      $ADWWWW(2,J))
181         DO 312 L=1,4
182         TERM=0.
183           DO 313 N=0,6
184           C=0.
185           N1=MAX(N-3,0)
186           N2=MIN(3,N)
187           DO 314 K=N1,N2
188 314       C=C+ANWWWW(K+1,I,L)*ANWWWW(N-K+1,J,L)
189           C=C*WTHELI(L)
190           IF(J.NE.I) C=2.*C
191           IF(L.EQ.4) THEN
192             TERM=TERM+.5*C*FINT(N+1)-.5*C*FINT(N+3)
193           ELSE
194             TERM=TERM+C*FINT(N+1)
195           ENDIF
196 313       CONTINUE
197         SUM=SUM+TERM
198 312     CONTINUE
199 311   CONTINUE
200 C          ADD INTEGRAL OF IMAGINARY PART SQUARED.
201       SUM=SUM+2.*ZLIM*(WTHELI(1)*AIWWWW(1)**2+WTHELI(2)*AIWWWW(2)**2
202      $+WTHELI(3)*AIWWWW(3)**2+WTHELI(4)*AIWWWW(4)**2)
203 C          CROSS SECTION
204       SIG0=SUM/(32.*PI*S*SCM)*UNITS
205       SIG1=.5*SIG0*QSAVE(27,1)*QSAVE(28,2)
206       IF(GOQ(27,1).AND.GOQ(28,2)) THEN
207         SIG=SIG1*TBRWW(2,1)*TBRWW(3,2)
208         CALL SIGFIL(SIG,27,28,27,28)
209       ENDIF
210       IF(GOQ(28,1).AND.GOQ(27,2)) THEN
211         SIG=SIG1*TBRWW(3,1)*TBRWW(2,2)
212         CALL SIGFIL(SIG,27,28,28,27)
213       ENDIF
214       SIG1=.5*SIG0*QSAVE(28,1)*QSAVE(27,2)
215       IF(GOQ(27,1).AND.GOQ(28,2)) THEN
216         SIG=SIG1*TBRWW(2,1)*TBRWW(3,2)
217         CALL SIGFIL(SIG,28,27,27,28)
218       ENDIF
219       IF(GOQ(28,1).AND.GOQ(27,2)) THEN
220         SIG=SIG1*TBRWW(3,1)*TBRWW(2,2)
221         CALL SIGFIL(SIG,28,27,28,27)
222       ENDIF
223 C
224 C          Z0 Z0 --> W+ W-
225 C
226 C          SET UP AMPLITUDES
227       IF(QMW.LE.2.*AMASS(90)) GO TO 500
228       CALL XZZWW
229 C          SUM CROSS SECTION TERMS. I,J RUN OVER AMPLITUDE TERMS.
230 C          L RUNS OVER HELICITY STATES. N RUNS OVER POWERS.
231 C          REMEMBER THAT L=4 IS MISSING SIN(THETA)/SQRT(2)
232       SUM=0.
233       DO 321 I=1,4
234       DO 321 J=I,4
235       CALL SIGINT(FINT,ZLIM,ADWWWW(1,I),ADWWWW(2,I),ADWWWW(1,J),
236      $ADWWWW(2,J))
237         DO 322 L=1,4
238         TERM=0.
239           DO 323 N=0,6
240           C=0.
241           N1=MAX(N-3,0)
242           N2=MIN(3,N)
243           DO 324 K=N1,N2
244 324       C=C+ANWWWW(K+1,I,L)*ANWWWW(N-K+1,J,L)
245           C=C*WTHELI(L)
246           IF(J.NE.I) C=2.*C
247           IF(L.EQ.4) THEN
248             TERM=TERM+.5*C*FINT(N+1)-.5*C*FINT(N+3)
249           ELSE
250             TERM=TERM+C*FINT(N+1)
251           ENDIF
252 323       CONTINUE
253         SUM=SUM+TERM
254 322     CONTINUE
255 321   CONTINUE
256 C          ADD INTEGRAL OF IMAGINARY PART SQUARED.
257       SUM=SUM+2.*ZLIM*(WTHELI(1)*AIWWWW(1)**2+WTHELI(2)*AIWWWW(2)**2
258      $+WTHELI(3)*AIWWWW(3)**2+WTHELI(4)*AIWWWW(4)**2)
259 C          CROSS SECTION
260       SIG0=SUM/(32.*PI*S*SCM)*UNITS
261       SIG1=.5*SIG0*QSAVE(29,1)*QSAVE(29,2)
262       IF(GOQ(27,1).AND.GOQ(28,2)) THEN
263         SIG=SIG1*TBRWW(2,1)*TBRWW(3,2)
264         CALL SIGFIL(SIG,29,29,27,28)
265       ENDIF
266       IF(GOQ(28,1).AND.GOQ(27,2)) THEN
267         SIG=SIG1*TBRWW(3,1)*TBRWW(2,2)
268         CALL SIGFIL(SIG,29,29,28,27)
269       ENDIF
270 C
271 C          W+ W- --> Z0 Z0
272 C
273 400   IF(QMW.LE.2.*AMASS(90)) GO TO 500
274       IF(.NOT.(GOQ(29,1).AND.GOQ(29,2))) GO TO 500
275       WM=AMASS(90)
276       PWWCM=.5*SQRT(QMW**2-4.*WM**2)
277       STHLIM=PTMIN(1)/PWWCM
278       IF(STHLIM.LE.1) THEN
279         ZLIM=SQRT(1.-STHLIM**2)
280       ELSE
281         GO TO 500
282       ENDIF
283 C          SET UP AMPLITUDES
284       CALL XWWZZ
285 C          SUM CROSS SECTION TERMS. I,J RUN OVER AMPLITUDE TERMS.
286 C          L RUNS OVER HELICITY STATES. N RUNS OVER POWERS.
287 C          REMEMBER THAT L=4 IS MISSING SIN(THETA)/SQRT(2)
288       SUM=0.
289       DO 411 I=1,4
290       DO 411 J=I,4
291       CALL SIGINT(FINT,ZLIM,ADWWWW(1,I),ADWWWW(2,I),ADWWWW(1,J),
292      $ADWWWW(2,J))
293         DO 412 L=1,4
294         TERM=0.
295           DO 413 N=0,6
296           C=0.
297           N1=MAX(N-3,0)
298           N2=MIN(3,N)
299           DO 414 K=N1,N2
300 414       C=C+ANWWWW(K+1,I,L)*ANWWWW(N-K+1,J,L)
301           C=C*WTHELI(L)
302           IF(J.NE.I) C=2.*C
303           IF(L.EQ.4) THEN
304             TERM=TERM+.5*C*FINT(N+1)-.5*C*FINT(N+3)
305           ELSE
306             TERM=TERM+C*FINT(N+1)
307           ENDIF
308 413       CONTINUE
309         SUM=SUM+TERM
310 412     CONTINUE
311 411   CONTINUE
312 C          ADD INTEGRAL OF IMAGINARY PART SQUARED.
313       SUM=SUM+2.*ZLIM*(WTHELI(1)*AIWWWW(1)**2+WTHELI(2)*AIWWWW(2)**2
314      $+WTHELI(3)*AIWWWW(3)**2+WTHELI(4)*AIWWWW(4)**2)
315 C          CROSS SECTION
316       SIG0=SUM/(32.*PI*S*SCM)*UNITS
317       SIG0=.5*SIG0
318       SIG0=SIG0*TBRWW(4,1)*TBRWW(4,2)
319       SIG=SIG0*QSAVE(27,1)*QSAVE(28,2)
320       CALL SIGFIL(SIG,27,28,29,29)
321       SIG=SIG0*QSAVE(28,1)*QSAVE(27,2)
322       CALL SIGFIL(SIG,28,27,29,29)
323 C
324 C          Z0 Z0 --> Z0 Z0
325 C
326 C          SET UP AMPLITUDES
327       CALL XZZZZ
328 C          SUM CROSS SECTION TERMS. I,J RUN OVER AMPLITUDE TERMS.
329 C          L RUNS OVER HELICITY STATES. N RUNS OVER POWERS.
330 C          REMEMBER THAT L=4 IS MISSING SIN(THETA)/SQRT(2)
331       SUM=0.
332       DO 421 I=1,4
333       DO 421 J=I,4
334       CALL SIGINT(FINT,ZLIM,ADWWWW(1,I),ADWWWW(2,I),ADWWWW(1,J),
335      $ADWWWW(2,J))
336         DO 422 L=1,4
337         TERM=0.
338           DO 423 N=0,6
339           C=0.
340           N1=MAX(N-3,0)
341           N2=MIN(3,N)
342           DO 424 K=N1,N2
343 424       C=C+ANWWWW(K+1,I,L)*ANWWWW(N-K+1,J,L)
344           C=C*WTHELI(L)
345           IF(J.NE.I) C=2.*C
346           IF(L.EQ.4) THEN
347             TERM=TERM+.5*C*FINT(N+1)-.5*C*FINT(N+3)
348           ELSE
349             TERM=TERM+C*FINT(N+1)
350           ENDIF
351 423       CONTINUE
352         SUM=SUM+TERM
353 422     CONTINUE
354 421   CONTINUE
355 C          ADD INTEGRAL OF IMAGINARY PART SQUARED.
356       SUM=SUM+2.*ZLIM*(WTHELI(1)*AIWWWW(1)**2+WTHELI(2)*AIWWWW(2)**2
357      $+WTHELI(3)*AIWWWW(3)**2+WTHELI(4)*AIWWWW(4)**2)
358 C          CROSS SECTION
359       SIG0=SUM/(32.*PI*S*SCM)*UNITS
360       SIG0=.5*SIG0
361       SIG0=SIG0*TBRWW(4,1)*TBRWW(4,2)
362       SIG=SIG0*QSAVE(29,1)*QSAVE(29,2)
363       CALL SIGFIL(SIG,29,29,29,29)
364 C
365 500   RETURN
366       END