]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ISAJET/code/zjj.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Using AliITSgeomTGeo
[u/mrichter/AliRoot.git] / ISAJET / code / zjj.F
CommitLineData
0795afa3 1#include "isajet/pilot.h"
2 SUBROUTINE ZJJ
3C-----------------------------------------------------------------------
4C
5C Use MadGraph/Helas to generate Z + 2 jets after setup by
6C ZJJ0 using cross section routines from MadGraph:
7C ZJJ1: q1 q1b -> Z q2 q2b, q1 != q2
8C ZJJ2: g g -> Z q2 q2b
9C ZJJ3: q1 q1b -> Z g g
10C ZJJ4: q1 q1b -> Z q1 q1b
11C ZJJ5: q1 q2 -> Z q1 q2
12C ZJJ6: q1 q1 -> Z q1 q1
13C ZJJ7: g q -> Z g q
14C
15C Note: The Z is always jet1, but the other two jets are
16C symmetrized so a symmetry factor of 1/2 is needed for every
17C subprocess. This is included by MadGraph for identical
18C particles!
19C
20C-----------------------------------------------------------------------
21#if defined(CERNLIB_IMPNONE)
22 IMPLICIT NONE
23#endif
24#include "isajet/const.inc"
25#include "isajet/q1q2.inc"
26#include "isajet/itapes.inc"
27#include "isajet/jetlim.inc"
28#include "isajet/jetpar.inc"
29#include "isajet/jetset.inc"
30#include "isajet/partcl.inc"
31#include "isajet/pinits.inc"
32#include "isajet/pjets.inc"
33#include "isajet/primar.inc"
34#include "isajet/sstype.inc"
35#include "isajet/totals.inc"
36#include "isajet/mgkin.inc"
37#include "isajet/mgcoms.inc"
38#include "isajet/mgsigs.inc"
39C
40 INTEGER IMAD(6)
41 REAL*8 P1(0:3),P2(0:3),P3(0:3),P4(0:3),P5(0:3)
42 EQUIVALENCE (P1(0),PJETS8(0,1))
43 EQUIVALENCE (P2(0),PJETS8(0,2))
44 EQUIVALENCE (P3(0),PJETS8(0,3))
45 EQUIVALENCE (P4(0),PJETS8(0,4))
46 EQUIVALENCE (P5(0),PJETS8(0,5))
47 REAL QFCN,XX,QQ,RND,RANF,SIG,FJAC,STRUC,ALQCD
48 REAL*8 SZJJ1,WT8,TERM,SUM,SZJJ2,SZJJ3,SZJJ4,SIG8,SIGI8
49 REAL*8 SZJJ5,SZJJ6,SZJJ7
50 INTEGER IQ,IH,ISIG8,IFL1,IFL2,IM1,IM2,IQ1,IQ2,NTRY,I,II,K,IWT8
51C
52C Map Jettype/2 to MadGraph
53 DATA IMAD/3,4,8,7,12,11/
54C
55C Parton distributions
56 QFCN(XX,IQ,IH)=STRUC(XX,QQ,IQ,IDIN(IH))/XX
57C
58C Begin
59C
60 NTRY=0
61 NJSET=0
62 NPTCL=0
63C
64C Select process
65C
66 RND=RANF()
67 ISIG8=0
68 SIG=0
69 DO 10 I=1,NSIG8
70 SIG=SIG+WTSUM8(I)/NWT8(I)
7110 CONTINUE
72 SUM=0
73 DO 20 I=1,NSIG8
74 II=ISORT8(NSIG8+1-I)
75 SUM=SUM+WTSUM8(II)/NWT8(II)
76 IF(SUM.GE.RND*SIG) THEN
77 ISIG8=II
78 GO TO 100
79 ENDIF
8020 CONTINUE
81 WRITE(ITLIS,*) 'ERROR IN ZJJ: NO MODE FOUND'
82 STOP99
83C
84100 CONTINUE
85 SIG8=0
86 FJAC=UNITS/SCM
87 NTRY=NTRY+1
88 IF(NTRY.GT.NTRIES) THEN
89 WRITE(ITLIS,*) 'ERROR IN ZJJ: NTRY = ',NTRY
90 WRITE(ITLIS,*) 'PROCESS WAS ',(IDENT8(K,ISIG8),K=1,5)
91 SIGI8=WTSUM8(ISIG8)/NWT8(ISIG8)
92 WRITE(ITLIS,*) 'PROCESS SIGMA/MAX = ',SIGI8,WTMAX8(ISIG8)
93 WRITE(ITLIS,*) 'CHECK YOUR LIMITS OR INCREASE NTRIES'
94 STOP99
95 ENDIF
96C
97C Cases 1,4: q1 q1b -> z q2 q2b
98C
99 IF(IFUNC8(ISIG8).EQ.1.OR.IFUNC8(ISIG8).EQ.4) THEN
100 AMJET8(3)=ZMASS
101 IFL1=IABS(IDENT8(1,ISIG8))
102 IM1=IMAD(IFL1)
103 IQ1=2*IFL1
104 IQ2=IQ1+1
105 AMJET8(1)=FMASS(IM1)
106 AMJET8(2)=FMASS(IM1)
107 IFL2=IABS(IDENT8(4,ISIG8))
108 IM2=IMAD(IFL2)
109 AMJET8(4)=FMASS(IM2)
110 AMJET8(5)=FMASS(IM2)
111 DO 210 I=1,NTRIES
112 IWT8=I
113 CALL MULJET(WT8)
114 IF(WT8.GT.0) GO TO 220
115210 CONTINUE
116 WRITE(ITLIS,*) 'ERROR IN ZJJ: NO PHASE SPACE POINT IN ',
117 $ NTRIES,' TRIES'
118 STOP99
119220 NWTTOT=NWTTOT+IWT8-1
120 NWT8(ISIG8)=NWT8(ISIG8)+IWT8-1
121 X1=(P1(0)+P1(3))/ECM
122 X2=(P2(0)-P2(3))/ECM
123 QQ=P3(1)**2+P3(2)**2+P4(1)**2+P4(2)**2+P5(1)**2+
124 $ P5(2)**2+AMJET8(3)**2+AMJET8(4)**2+AMJET8(5)**2
125C
126C Subcases
127C
128 IF(IDENT8(1,ISIG8).GT.0.AND.IDENT8(4,ISIG8).GT.0) THEN
129 IF(IFUNC8(ISIG8).EQ.1) THEN
130 TERM=SZJJ1(P1,P2,P3,P4,P5,IM1,IM2)
131 ELSE
132 TERM=SZJJ4(P1,P2,P3,P4,P5,IM1)
133 ENDIF
134 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
135 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
136 ELSEIF(IDENT8(1,ISIG8).GT.0.AND.IDENT8(4,ISIG8).LT.0) THEN
137 IF(IFUNC8(ISIG8).EQ.1) THEN
138 TERM=SZJJ1(P1,P2,P3,P5,P4,IM1,IM2)
139 ELSE
140 TERM=SZJJ4(P1,P2,P3,P5,P4,IM1)
141 ENDIF
142 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
143 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
144 ELSEIF(IDENT8(1,ISIG8).LT.0.AND.IDENT8(4,ISIG8).GT.0) THEN
145 IF(IFUNC8(ISIG8).EQ.1) THEN
146 TERM=SZJJ1(P1,P2,P3,P4,P5,IM1,IM2)
147 ELSE
148 TERM=SZJJ4(P1,P2,P3,P4,P5,IM1)
149 ENDIF
150 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
151 TERM=TERM*WT8*FJAC*QFCN(X1,IQ2,1)*QFCN(X2,IQ1,2)
152 ELSEIF(IDENT8(1,ISIG8).LT.0.AND.IDENT8(4,ISIG8).LT.0) THEN
153 IF(IFUNC8(ISIG8).EQ.1) THEN
154 TERM=SZJJ1(P1,P2,P3,P5,P4,IM1,IM2)
155 ELSE
156 TERM=SZJJ4(P1,P2,P3,P5,P4,IM1)
157 ENDIF
158 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
159 TERM=TERM*WT8*FJAC*QFCN(X1,IQ2,1)*QFCN(X2,IQ1,2)
160 ELSE
161 WRITE(ITLIS,*) 'ERROR IN ZJJ...INVALID FLAVOR FOR ZJJ1'
162 STOP99
163 ENDIF
164 SIG8=0.5*TERM
165 GO TO 900
166 ENDIF
167C
168C Case 2: g g -> z q2 q2b
169C
170 IF(IFUNC8(ISIG8).EQ.2) THEN
171 AMJET8(3)=ZMASS
172 IFL1=IABS(IDENT8(1,ISIG8))
173 AMJET8(1)=0
174 AMJET8(2)=0
175 IFL2=IABS(IDENT8(4,ISIG8))
176 IM2=IMAD(IFL2)
177 AMJET8(4)=FMASS(IM2)
178 AMJET8(5)=FMASS(IM2)
179 DO 310 I=1,NTRIES
180 IWT8=I
181 CALL MULJET(WT8)
182 IF(WT8.GT.0) GO TO 320
183310 CONTINUE
184 WRITE(ITLIS,*) 'ERROR IN ZJJ: NO PHASE SPACE POINT IN ',
185 $ NTRIES,' TRIES'
186 STOP99
187320 NWTTOT=NWTTOT+IWT8-1
188 NWT8(ISIG8)=NWT8(ISIG8)+IWT8-1
189 X1=(P1(0)+P1(3))/ECM
190 X2=(P2(0)-P2(3))/ECM
191 QQ=P3(1)**2+P3(2)**2+P4(1)**2+P4(2)**2+P5(1)**2+
192 $ P5(2)**2+AMJET8(3)**2+AMJET8(4)**2+AMJET8(5)**2
193C
194C Subcases
195C
196 IF(IDENT8(4,ISIG8).GT.0) THEN
197 TERM=SZJJ2(P1,P2,P3,P4,P5,IM2)
198 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
199 TERM=TERM*WT8*FJAC*QFCN(X1,1,1)*QFCN(X2,1,2)
200 ELSEIF(IDENT8(4,ISIG8).LT.0) THEN
201 TERM=SZJJ2(P1,P2,P3,P5,P4,IM2)
202 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
203 TERM=TERM*WT8*FJAC*QFCN(X1,1,1)*QFCN(X2,1,2)
204 ELSE
205 WRITE(ITLIS,*) 'ERROR IN ZJJ...INVALID FLAVOR FOR ZJJ2'
206 STOP99
207 ENDIF
208 SIG8=0.5*TERM
209 GO TO 900
210 ENDIF
211C
212C Case 3: q1 q1b -> z g g
213C
214 IF(IFUNC8(ISIG8).EQ.3) THEN
215 AMJET8(3)=ZMASS
216 IFL1=IABS(IDENT8(1,ISIG8))
217 IQ1=2*IFL1
218 IQ2=IQ1+1
219 IM1=IMAD(IFL1)
220 AMJET8(1)=FMASS(IM1)
221 AMJET8(2)=FMASS(IM1)
222 IFL2=9
223 AMJET8(4)=0
224 AMJET8(5)=0
225 DO 410 I=1,NTRIES
226 IWT8=I
227 CALL MULJET(WT8)
228 IF(WT8.GT.0) GO TO 420
229410 CONTINUE
230 WRITE(ITLIS,*) 'ERROR IN ZJJ: NO PHASE SPACE POINT IN ',
231 $ NTRIES,' TRIES'
232 STOP99
233420 NWTTOT=NWTTOT+IWT8-1
234 NWT8(ISIG8)=NWT8(ISIG8)+IWT8-1
235 X1=(P1(0)+P1(3))/ECM
236 X2=(P2(0)-P2(3))/ECM
237 QQ=P3(1)**2+P3(2)**2+P4(1)**2+P4(2)**2+P5(1)**2+
238 $ P5(2)**2+AMJET8(3)**2+AMJET8(4)**2+AMJET8(5)**2
239C
240C Subcases
241C
242 IF(IDENT8(1,ISIG8).GT.0) THEN
243 TERM=SZJJ3(P1,P2,P3,P4,P5,IM1)
244 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
245 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
246 ELSEIF(IDENT8(1,ISIG8).LT.0) THEN
247 TERM=SZJJ3(P2,P1,P3,P4,P5,IM1)
248 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
249 TERM=TERM*WT8*FJAC*QFCN(X1,IQ2,1)*QFCN(X2,IQ1,2)
250 ELSE
251 WRITE(ITLIS,*) 'ERROR IN ZJJ...INVALID FLAVOR FOR ZJJ3'
252 STOP99
253 ENDIF
254 SIG8=TERM
255 GO TO 900
256 ENDIF
257C
258C Cases 5,6: q1 q2 -> z q1 q2
259C
260 IF(IFUNC8(ISIG8).EQ.5.OR.IFUNC8(ISIG8).EQ.6) THEN
261 IFL1=IABS(IDENT8(1,ISIG8))
262 IM1=IMAD(IFL1)
263 IFL2=IABS(IDENT8(2,ISIG8))
264 IM2=IMAD(IFL2)
265 IQ1=2*IFL1
266 IQ2=2*IFL2
267 IF(IDENT8(1,ISIG8).LT.0) IQ1=IQ1+1
268 IF(IDENT8(2,ISIG8).LT.0) IQ2=IQ2+1
269 AMJET8(1)=FMASS(IM1)
270 AMJET8(2)=FMASS(IM2)
271 AMJET8(3)=ZMASS
272 AMJET8(4)=FMASS(IM1)
273 AMJET8(5)=FMASS(IM2)
274 DO 510 I=1,NTRIES
275 IWT8=I
276 CALL MULJET(WT8)
277 IF(WT8.GT.0) GO TO 520
278510 CONTINUE
279 WRITE(ITLIS,*) 'ERROR IN ZJJ: NO PHASE SPACE POINT IN ',
280 $ NTRIES,' TRIES'
281 STOP99
282520 NWTTOT=NWTTOT+IWT8-1
283 NWT8(ISIG8)=NWT8(ISIG8)+IWT8-1
284 X1=(P1(0)+P1(3))/ECM
285 X2=(P2(0)-P2(3))/ECM
286 QQ=P3(1)**2+P3(2)**2+P4(1)**2+P4(2)**2+P5(1)**2+
287 $ P5(2)**2+AMJET8(3)**2+AMJET8(4)**2+AMJET8(5)**2
288C
289C Subcases
290C
291 IF(IDENT8(1,ISIG8).EQ.IDENT8(4,ISIG8)) THEN
292 IF(IFUNC8(ISIG8).EQ.5) THEN
293 TERM=SZJJ5(P1,P2,P3,P4,P5,IM1,IM2)
294 ELSE
295 TERM=SZJJ6(P1,P2,P3,P4,P5,IM1)
296 ENDIF
297 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
298 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
299 ELSEIF(IDENT8(1,ISIG8).EQ.IDENT8(5,ISIG8)) THEN
300 TERM=SZJJ5(P1,P2,P3,P5,P4,IM1,IM2)
301 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
302 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
303 ELSE
304 WRITE(ITLIS,*) 'ERROR IN ZJJ...INVALID FLAVOR FOR ZJJ1'
305 STOP99
306 ENDIF
307 SIG8=TERM
308 IF(IFL1.NE.IFL2) SIG8=0.5*SIG8
309 GO TO 900
310 ENDIF
311C
312C Case 7: g q -> z g q
313C
314 IF(IFUNC8(ISIG8).EQ.7) THEN
315 IF(IDENT8(1,ISIG8).EQ.9) THEN
316 IFL1=IABS(IDENT8(2,ISIG8))
317 IM1=IMAD(IFL1)
318 AMJET8(1)=0
319 AMJET8(2)=FMASS(IM1)
320 IQ1=1
321 IQ2=2*IFL1
322 IF(IDENT8(2,ISIG8).LT.0) IQ2=IQ2+1
323 ELSE
324 IFL1=IABS(IDENT8(1,ISIG8))
325 IM1=IMAD(IFL1)
326 AMJET8(1)=FMASS(IM1)
327 AMJET8(2)=0
328 IQ2=1
329 IQ1=2*IFL1
330 IF(IDENT8(1,ISIG8).LT.0) IQ1=IQ1+1
331 ENDIF
332 AMJET8(3)=ZMASS
333 IF(IDENT8(4,ISIG8).EQ.9) THEN
334 AMJET8(4)=0
335 AMJET8(5)=FMASS(IM1)
336 ELSE
337 AMJET8(4)=FMASS(IM1)
338 AMJET8(5)=0
339 ENDIF
340 DO 610 I=1,NTRIES
341 IWT8=I
342 CALL MULJET(WT8)
343 IF(WT8.GT.0) GO TO 620
344610 CONTINUE
345 WRITE(ITLIS,*) 'ERROR IN ZJJ: NO PHASE SPACE POINT IN ',
346 $ NTRIES,' TRIES'
347 STOP99
348620 NWTTOT=NWTTOT+IWT8-1
349 NWT8(ISIG8)=NWT8(ISIG8)+IWT8-1
350 X1=(P1(0)+P1(3))/ECM
351 X2=(P2(0)-P2(3))/ECM
352 QQ=P3(1)**2+P3(2)**2+P4(1)**2+P4(2)**2+P5(1)**2+
353 $ P5(2)**2+AMJET8(3)**2+AMJET8(4)**2+AMJET8(5)**2
354C
355C Subcases
356C
357 IF(IDENT8(1,ISIG8).EQ.9.AND.IDENT8(4,ISIG8).EQ.9) THEN
358 TERM=SZJJ7(P1,P2,P3,P4,P5,IM1)
359 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
360 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
361 ELSEIF(IDENT8(2,ISIG8).EQ.9.AND.IDENT8(4,ISIG8).EQ.9) THEN
362 TERM=SZJJ7(P2,P1,P3,P4,P5,IM1)
363 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
364 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
365 ELSEIF(IDENT8(1,ISIG8).EQ.9.AND.IDENT8(5,ISIG8).EQ.9) THEN
366 TERM=SZJJ7(P1,P2,P3,P5,P4,IM1)
367 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
368 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
369 ELSEIF(IDENT8(2,ISIG8).EQ.9.AND.IDENT8(5,ISIG8).EQ.9) THEN
370 TERM=SZJJ7(P2,P1,P3,P5,P4,IM1)
371 TERM=TERM*(4*PI*ALQCD(REAL(QQ)))**2
372 TERM=TERM*WT8*FJAC*QFCN(X1,IQ1,1)*QFCN(X2,IQ2,2)
373 ELSE
374 WRITE(ITLIS,*) 'ERROR IN ZJJ...INVALID FLAVOR FOR ZJJ1'
375 STOP99
376 ENDIF
377 SIG8=0.5*TERM
378 GO TO 900
379 ENDIF
380C
381C Increment totals and test
382C
383900 WTTOT8=WTTOT8+SIG8
384 NWTTOT=NWTTOT+1
385 WTSUM8(ISIG8)=WTSUM8(ISIG8)+SIG8
386 WTMAX8(ISIG8)=MAX(WTMAX8(ISIG8),SIG8)
387 NWT8(ISIG8)=NWT8(ISIG8)+1
388 IF(SIG8.LT.RANF()*WTMAX8(ISIG8)) GO TO 100
389C
390C Good event
391C
392 DO 910 I=1,3
393 DO 911 K=1,3
394 PJETS(K,I)=PJETS8(K,I+2)
395911 CONTINUE
396 PJETS(4,I)=PJETS8(0,I+2)
397 PJETS(5,I)=AMJET8(I+2)
398 IDJETS(I)=IDENT8(I+2,ISIG8)
399910 CONTINUE
400 DO 920 I=1,2
401 DO 921 K=1,3
402 PINITS(K,I)=PJETS8(K,I)
403921 CONTINUE
404 PINITS(4,I)=PJETS8(0,I)
405 PINITS(5,I)=AMJET8(I)
406 IDINIT(I)=IDENT8(I,ISIG8)
407920 CONTINUE
408C
409 QSQ=QQ
410 SHAT=(P1(0)+P2(0))**2-(P1(3)+P2(3))**2
411 PBEAM(1)=(1.-X1)*HALFE
412 PBEAM(2)=(1.-X2)*HALFE
413C
414C Set /TOTALS/
415C
416 NKINPT=NWTTOT
417 SUMWT=WTTOT8
418C
419 RETURN
420 END