2 C*********************************************************************
6 C...Generates quantities characterizing the high-pT scattering at the
7 C...parton level according to the matrix elements. Chooses incoming,
8 C...reacting partons, their momentum fractions and one of the possible
10 COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
11 COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4)
12 COMMON/PYSUBS/MSEL,MSUB(200),KFIN(2,-40:40),CKIN(200)
13 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
14 COMMON/PYINT1/MINT(400),VINT(400)
15 COMMON/PYINT2/ISET(200),KFPR(200,2),COEF(200,20),ICOL(40,4,2)
16 COMMON/PYINT3/XSFX(2,-40:40),ISIG(1000,3),SIGH(1000)
17 COMMON/PYINT4/WIDP(21:40,0:40),WIDE(21:40,0:40),WIDS(21:40,3)
18 COMMON/PYINT5/NGEN(0:200,3),XSEC(0:200,3)
19 COMMON/PYINT7/SIGT(0:6,0:6,0:5)
20 COMMON/PYINT9/DXSEC(0:200)
21 DOUBLE PRECISION DXSEC
22 COMMON/PYUPPR/NUP,KUP(20,7),PUP(20,5),NFUP,IFUP(10,2),Q2UP(0:10)
23 SAVE /LUDAT1/,/LUDAT2/
24 SAVE /PYSUBS/,/PYPARS/,/PYINT1/,/PYINT2/,/PYINT3/,/PYINT4/,
25 &/PYINT5/,/PYINT7/,/PYINT9/,/PYUPPR/
26 DIMENSION XPQ(-25:25),PMM(2),PDIF(4),BHAD(4)
28 C...Parameters and data used in elastic/diffractive treatment.
29 DATA EPS/0.0808/, ALP/0.25/, CRES/2./, PMRC/1.062/, SMP/0.880/
30 DATA BHAD/2.3,1.4,1.4,0.23/
32 C...Initial values, specifically for (first) semihard interaction.
39 IF(MSTP(171).EQ.1.AND.MSTP(172).EQ.2) MFAIL=1
45 C...Choice of process type - first event of pileup.
46 IF(MINT(82).EQ.1.AND.(ISUB.LE.90.OR.ISUB.GT.96)) THEN
48 C...For gamma-p or gamma-gamma first pick between alternatives.
49 IF(MINT(121).GT.1) CALL PYSAVE(4,IGA)
52 C...For gamma + gamma with different nature, flip at random.
53 IF(MINT(11).EQ.22.AND.MINT(12).EQ.22.AND.MINT(123).GE.4.AND.
64 IF(MINT(47).EQ.2.OR.MINT(47).EQ.3) MINT(47)=5-MINT(47)
67 C...Pick process type.
70 IF(MSUB(I).NE.1) GOTO 110
73 IF(RSUB.LE.0.) GOTO 120
75 120 IF(ISUB.EQ.95) ISUB=96
76 IF(ISUB.EQ.96) CALL PYMULT(2)
78 C...Choice of inclusive process type - pileup events.
79 ELSEIF(MINT(82).GE.2.AND.ISUB.EQ.0) THEN
82 IF(RSUB.GT.SIGT(0,0,5)) ISUB=94
83 IF(RSUB.GT.SIGT(0,0,5)+SIGT(0,0,4)) ISUB=93
84 IF(RSUB.GT.SIGT(0,0,5)+SIGT(0,0,4)+SIGT(0,0,3)) ISUB=92
85 IF(RSUB.GT.SIGT(0,0,5)+SIGT(0,0,4)+SIGT(0,0,3)+SIGT(0,0,2))
87 IF(ISUB.EQ.96) CALL PYMULT(2)
89 IF(MINT(82).EQ.1) NGEN(0,1)=NGEN(0,1)+1
90 IF(MINT(82).EQ.1) NGEN(ISUB,1)=NGEN(ISUB,1)+1
91 IF(ISUB.EQ.96.AND.LOOP.EQ.1.AND.MINT(82).EQ.1)
92 &NGEN(97,1)=NGEN(97,1)+1
96 C...Find resonances (explicit or implicit in cross-section).
99 IF(ISTSB.EQ.1.OR.ISTSB.EQ.3.OR.ISTSB.EQ.5) THEN
101 ELSEIF(ISUB.EQ.24.OR.ISUB.EQ.25.OR.ISUB.EQ.110.OR.ISUB.EQ.165.OR.
102 &ISUB.EQ.171.OR.ISUB.EQ.176) THEN
104 ELSEIF(ISUB.EQ.23.OR.ISUB.EQ.26.OR.ISUB.EQ.166.OR.ISUB.EQ.172.OR.
107 ELSEIF(ISUB.GE.71.AND.ISUB.LE.77) THEN
109 IF(MSTP(46).EQ.5) THEN
112 PMAS(30,2)=PARP(45)**3/(96.*PARU(1)*PARP(47)**2)
116 IF(CKMX.LE.0.) CKMX=VINT(1)
118 IF(CKIN(1).GT.PMAS(KFR1,1)+20.*PMAS(KFR1,2).OR.
119 & CKMX.LT.PMAS(KFR1,1)-20.*PMAS(KFR1,2)) KFR1=0
122 TAUR1=PMAS(KFR1,1)**2/VINT(2)
123 GAMR1=PMAS(KFR1,1)*PMAS(KFR1,2)/VINT(2)
131 TAUR2=PMAS(KFR2,1)**2/VINT(2)
132 GAMR2=PMAS(KFR2,1)*PMAS(KFR2,2)/VINT(2)
133 IF(CKIN(1).GT.PMAS(KFR2,1)+20.*PMAS(KFR2,2).OR.
134 & CKMX.LT.PMAS(KFR2,1)-20.*PMAS(KFR2,2)) KFR2=0
135 IF(KFR2.NE.0.AND.KFR1.NE.0) THEN
140 ELSEIF(KFR2.NE.0) THEN
151 C...Find product masses and minimum pT of process,
152 C...optionally with broadening according to a truncated Breit-Wigner.
157 IF(MINT(82).GE.2) VINT(71)=0.
159 IF(ISTSB.EQ.2.OR.ISTSB.EQ.4) THEN
162 IF(KFPR(ISUB,I).EQ.0) THEN
163 ELSEIF(MSTP(42).LE.0.OR.PMAS(LUCOMP(KFPR(ISUB,I)),2).LT.
165 VINT(62+I)=PMAS(LUCOMP(KFPR(ISUB,I)),1)**2
171 CALL PYOFSH(4,0,KFPR(ISUB,1),KFPR(ISUB,2),0.,PQM3,PQM4)
172 IF(MINT(51).EQ.1) THEN
173 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
183 IF(MIN(VINT(63),VINT(64)).LT.CKIN(6)**2) MINT(71)=1
184 IF(MINT(71).EQ.1) VINT(71)=MAX(CKIN(3),CKIN(5))
185 ELSEIF(ISTSB.EQ.6) THEN
186 CALL PYOFSH(6,0,KFPR(ISUB,1),KFPR(ISUB,2),0.,PQM3,PQM4)
187 IF(MINT(51).EQ.1) THEN
188 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
199 C...Prepare for additional variable choices in 2 -> 3.
202 IF(KFPR(ISUB,2).GT.0) VINT(201)=PMAS(KFPR(ISUB,2),1)
205 IF(ISUB.EQ.124) VINT(204)=PMAS(24,1)
206 IF(ISUB.EQ.121.OR.ISUB.EQ.122.OR.ISUB.EQ.181.OR.ISUB.EQ.182.OR.
207 & ISUB.EQ.186.OR.ISUB.EQ.187) VINT(204)=VINT(201)
211 C...Select incoming VDM particle (rho/omega/phi/J/psi).
212 IF(ISTSB.NE.0.AND.(MINT(101).GE.2.OR.MINT(102).GE.2).AND.
213 &(MINT(123).EQ.2.OR.MINT(123).EQ.5.OR.MINT(123).EQ.7)) THEN
214 VRN=RLU(0)*SIGT(0,0,5)
215 IF(MINT(101).LE.1) THEN
222 IF(MINT(102).LE.1) THEN
233 VRN=VRN-SIGT(I1,I2,5)
234 IF(VRN.LE.0.) GOTO 160
237 160 IF(MINT(101).GE.2) MINT(103)=KFV1
238 IF(MINT(102).GE.2) MINT(104)=KFV2
242 C...Elastic scattering or single or double diffractive scattering.
244 C...Select incoming particle (rho/omega/phi/J/psi for VDM) and mass.
249 IF(MINT(101).GE.2.OR.MINT(102).GE.2) THEN
251 VRN=RLU(0)*SIGT(0,0,JJ)
252 IF(MINT(101).LE.1) THEN
259 IF(MINT(102).LE.1) THEN
270 VRN=VRN-SIGT(I1,I2,JJ)
271 IF(VRN.LE.0.) GOTO 190
274 190 IF(MINT(101).GE.2) THEN
278 IF(MINT(102).GE.2) THEN
284 C...Side/sides of diffractive system.
287 IF(ISUB.EQ.92.OR.ISUB.EQ.94) MINT(17)=1
288 IF(ISUB.EQ.93.OR.ISUB.EQ.94) MINT(18)=1
290 C...Find masses of particles and minimal masses of diffractive states.
294 IF(MINT(16+JT).EQ.1) PDIF(JT)=PDIF(JT)+PARP(102)
301 SMRES1=(PMM(1)+PMRC)**2
302 SMRES2=(PMM(2)+PMRC)**2
304 C...Find elastic slope and lower limit diffractive slope.
305 IHA=MAX(2,IABS(MINT(103))/110)
307 IHB=MAX(2,IABS(MINT(104))/110)
310 BMN=2.*BHAD(IHA)+2.*BHAD(IHB)+4.*SH**EPS-4.2
311 ELSEIF(ISUB.EQ.92) THEN
312 BMN=MAX(2.,2.*BHAD(IHB))
313 ELSEIF(ISUB.EQ.93) THEN
314 BMN=MAX(2.,2.*BHAD(IHA))
315 ELSEIF(ISUB.EQ.94) THEN
319 C...Determine maximum possible t range and coefficient of generation.
320 SQLA12=(SH-SQM1-SQM2)**2-4.*SQM1*SQM2
321 SQLA34=(SH-SQM3-SQM4)**2-4.*SQM3*SQM4
322 THA=SH-(SQM1+SQM2+SQM3+SQM4)+(SQM1-SQM2)*(SQM3-SQM4)/SH
323 THB=SQRT(MAX(0.,SQLA12))*SQRT(MAX(0.,SQLA34))/SH
324 THC=(SQM3-SQM1)*(SQM4-SQM2)+(SQM1+SQM4-SQM2-SQM3)*
325 & (SQM1*SQM4-SQM2*SQM3)/SH
328 THRND=EXP(MAX(-50.,BMN*(THL-THU)))-1.
330 C...Select diffractive mass/masses according to dm^2/m^2.
332 IF(MINT(16+JT).EQ.0) THEN
336 PMMAX=MAX(VINT(2+JT),VINT(1)-PDIF(3-JT))
337 PDIF(2+JT)=PMMIN*(PMMAX/PMMIN)**RLU(0)
343 C..Additional mass factors, including resonance enhancement.
344 IF(PDIF(3)+PDIF(4).GE.VINT(1)) GOTO 210
346 FSD=(1.-SQM3/SH)*(1.+CRES*SMRES1/(SMRES1+SQM3))
347 IF(FSD.LT.RLU(0)*(1.+CRES)) GOTO 210
348 ELSEIF(ISUB.EQ.93) THEN
349 FSD=(1.-SQM4/SH)*(1.+CRES*SMRES2/(SMRES2+SQM4))
350 IF(FSD.LT.RLU(0)*(1.+CRES)) GOTO 210
351 ELSEIF(ISUB.EQ.94) THEN
352 FDD=(1.-(PDIF(3)+PDIF(4))**2/SH)*(SH*SMP/(SH*SMP+SQM3*SQM4))*
353 & (1.+CRES*SMRES1/(SMRES1+SQM3))*(1.+CRES*SMRES2/(SMRES2+SQM4))
354 IF(FDD.LT.RLU(0)*(1.+CRES)**2) GOTO 210
357 C...Select t according to exp(Bmn*t) and correct to right slope.
358 TH=THU+LOG(1.+THRND*RLU(0))/BMN
361 BADD=2.*ALP*LOG(SH/SQM3)
362 IF(BHAD(IHB).LT.1.) BADD=MAX(0.,BADD+2.*BHAD(IHB)-2.)
363 ELSEIF(ISUB.EQ.93) THEN
364 BADD=2.*ALP*LOG(SH/SQM4)
365 IF(BHAD(IHA).LT.1.) BADD=MAX(0.,BADD+2.*BHAD(IHA)-2.)
366 ELSEIF(ISUB.EQ.94) THEN
367 BADD=2.*ALP*(LOG(EXP(4.)+SH/(ALP*SQM3*SQM4))-4.)
369 IF(EXP(MAX(-50.,BADD*(TH-THU))).LT.RLU(0)) GOTO 210
372 C...Check whether m^2 and t choices are consistent.
373 SQLA34=(SH-SQM3-SQM4)**2-4.*SQM3*SQM4
374 THA=SH-(SQM1+SQM2+SQM3+SQM4)+(SQM1-SQM2)*(SQM3-SQM4)/SH
375 THB=SQRT(MAX(0.,SQLA12))*SQRT(MAX(0.,SQLA34))/SH
376 IF(THB.LE.1E-8) GOTO 210
377 THC=(SQM3-SQM1)*(SQM4-SQM2)+(SQM1+SQM4-SQM2-SQM3)*
378 & (SQM1*SQM4-SQM2*SQM3)/SH
381 IF(TH.LT.THLM.OR.TH.GT.THUM) GOTO 210
383 C...Information to output.
386 VINT(23)=MIN(1.,MAX(-1.,(THA+2.*TH)/THB))
388 VINT(59)=2.*SQRT(MAX(0.,-(THC+THA*TH+TH**2)))/THB
392 C...Note: in the following, by In is meant the integral over the
393 C...quantity multiplying coefficient cn.
394 C...Choose tau according to h1(tau)/tau, where
395 C...h1(tau) = c1 + I1/I2*c2*1/tau + I1/I3*c3*1/(tau+tau_R) +
396 C...I1/I4*c4*tau/((s*tau-m^2)^2+(m*Gamma)^2) +
397 C...I1/I5*c5*1/(tau+tau_R') +
398 C...I1/I6*c6*tau/((s*tau-m'^2)^2+(m'*Gamma')^2) +
399 C...I1/I7*c7*tau/(1.-tau), and
400 C...c1 + c2 + c3 + c4 + c5 + c6 + c7 = 1.
401 ELSEIF(ISTSB.GE.1.AND.ISTSB.LE.6) THEN
403 IF(MINT(51).NE.0) THEN
404 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
413 IF(RTAU.GT.COEF(ISUB,1)) MTAU=2
414 IF(RTAU.GT.COEF(ISUB,1)+COEF(ISUB,2)) MTAU=3
415 IF(RTAU.GT.COEF(ISUB,1)+COEF(ISUB,2)+COEF(ISUB,3)) MTAU=4
416 IF(RTAU.GT.COEF(ISUB,1)+COEF(ISUB,2)+COEF(ISUB,3)+COEF(ISUB,4))
418 IF(RTAU.GT.COEF(ISUB,1)+COEF(ISUB,2)+COEF(ISUB,3)+COEF(ISUB,4)+
419 & COEF(ISUB,5)) MTAU=6
420 IF(RTAU.GT.COEF(ISUB,1)+COEF(ISUB,2)+COEF(ISUB,3)+COEF(ISUB,4)+
421 & COEF(ISUB,5)+COEF(ISUB,6)) MTAU=7
422 CALL PYKMAP(1,MTAU,RLU(0))
424 C...2 -> 3, 4 processes:
425 C...Choose tau' according to h4(tau,tau')/tau', where
426 C...h4(tau,tau') = c1 + I1/I2*c2*(1 - tau/tau')^3/tau' +
427 C...I1/I3*c3*1/(1 - tau'), and c1 + c2 + c3 = 1.
428 IF(ISTSB.GE.3.AND.ISTSB.LE.5) THEN
430 IF(MINT(51).NE.0) THEN
431 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
440 IF(RTAUP.GT.COEF(ISUB,18)) MTAUP=2
441 IF(RTAUP.GT.COEF(ISUB,18)+COEF(ISUB,19)) MTAUP=3
442 CALL PYKMAP(4,MTAUP,RLU(0))
445 C...Choose y* according to h2(y*), where
446 C...h2(y*) = I0/I1*c1*(y*-y*min) + I0/I2*c2*(y*max-y*) +
447 C...I0/I3*c3*1/cosh(y*) + I0/I4*c4*1/(1-exp(y*-y*max)) +
448 C...I0/I5*c5*1/(1-exp(-y*-y*min)), I0 = y*max-y*min,
449 C...and c1 + c2 + c3 + c4 + c5 = 1.
451 IF(MINT(51).NE.0) THEN
452 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
461 IF(RYST.GT.COEF(ISUB,8)) MYST=2
462 IF(RYST.GT.COEF(ISUB,8)+COEF(ISUB,9)) MYST=3
463 IF(RYST.GT.COEF(ISUB,8)+COEF(ISUB,9)+COEF(ISUB,10)) MYST=4
464 IF(RYST.GT.COEF(ISUB,8)+COEF(ISUB,9)+COEF(ISUB,10)+
465 & COEF(ISUB,11)) MYST=5
466 CALL PYKMAP(2,MYST,RLU(0))
468 C...2 -> 2 processes:
469 C...Choose cos(theta-hat) (cth) according to h3(cth), where
470 C...h3(cth) = c0 + I0/I1*c1*1/(A - cth) + I0/I2*c2*1/(A + cth) +
471 C...I0/I3*c3*1/(A - cth)^2 + I0/I4*c4*1/(A + cth)^2,
472 C...A = 1 + 2*(m3*m4/sh)^2 (= 1 for massless products),
473 C...and c0 + c1 + c2 + c3 + c4 = 1.
475 IF(MINT(51).NE.0) THEN
476 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
483 IF(ISTSB.EQ.2.OR.ISTSB.EQ.4.OR.ISTSB.EQ.6) THEN
486 IF(RCTH.GT.COEF(ISUB,13)) MCTH=2
487 IF(RCTH.GT.COEF(ISUB,13)+COEF(ISUB,14)) MCTH=3
488 IF(RCTH.GT.COEF(ISUB,13)+COEF(ISUB,14)+COEF(ISUB,15)) MCTH=4
489 IF(RCTH.GT.COEF(ISUB,13)+COEF(ISUB,14)+COEF(ISUB,15)+
490 & COEF(ISUB,16)) MCTH=5
491 CALL PYKMAP(3,MCTH,RLU(0))
494 C...2 -> 3 : select pT1, phi1, pT2, phi2, y3 for 3 outgoing.
497 IF(MINT(51).NE.0) THEN
498 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
507 C...Low-pT or multiple interactions (first semihard interaction).
508 ELSEIF(ISTSB.EQ.9) THEN
512 C...Generate user-defined process: kinematics plus weight.
513 ELSEIF(ISTSB.EQ.11) THEN
515 CALL PYUPEV(ISUB,SIGS)
519 IF(MINT(82).EQ.1) THEN
520 NGEN(0,1)=NGEN(0,1)-1
521 NGEN(0,2)=NGEN(0,2)-1
522 NGEN(ISUB,1)=NGEN(ISUB,1)-1
524 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
528 C...Construct 'trivial' kinematical variables needed.
531 VINT(41)=2.*PUP(1,4)/VINT(1)
532 VINT(42)=2.*PUP(2,4)/VINT(1)
533 VINT(21)=VINT(41)*VINT(42)
534 VINT(22)=0.5*LOG(VINT(41)/VINT(42))
535 VINT(44)=VINT(21)*VINT(2)
536 VINT(43)=SQRT(MAX(0.,VINT(44)))
538 VINT(55)=SQRT(MAX(0.,VINT(56)))
540 C...Construct other kinematical variables needed (approximately).
543 VINT(45)=-0.5*VINT(44)
544 VINT(46)=-0.5*VINT(44)
554 IF(KUP(IUP,1).EQ.1) VINT(25)=VINT(25)+2.*(PUP(IUP,5)**2+
555 & PUP(IUP,1)**2+PUP(IUP,2)**2)/VINT(1)
556 IF(KUP(IUP,1).EQ.1) VINT(48)=VINT(48)+0.5*(PUP(IUP,1)**2+
559 VINT(47)=SQRT(VINT(48))
561 C...Calculate structure function weights.
562 IF(MINT(47).GE.2) THEN
563 DO 250 I=3-MIN(2,MINT(45)),MIN(2,MINT(46))
564 MINT(105)=MINT(102+I)
565 MINT(109)=MINT(106+I)
566 IF(MSTP(57).LE.1) THEN
567 CALL PYSTFU(MINT(10+I),VINT(40+I),Q2UP(0),XPQ)
569 CALL PYSTFL(MINT(10+I),VINT(40+I),Q2UP(0),XPQ)
578 C...Choose azimuthal angle.
579 VINT(24)=PARU(2)*RLU(0)
581 C...Check against user cuts on kinematics at parton level.
583 IF((ISUB.LE.90.OR.ISUB.GT.100).AND.ISTSB.LE.10) CALL PYKLIM(0)
584 IF(MINT(51).NE.0) THEN
585 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
592 IF(MINT(82).EQ.1.AND.MSTP(141).GE.1.AND.ISTSB.LE.10) THEN
594 IF(MSUB(91)+MSUB(92)+MSUB(93)+MSUB(94)+MSUB(95).EQ.0)
597 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
606 C...Calculate differential cross-section for different subprocesses.
607 IF(ISTSB.LE.10) CALL PYSIGH(NCHN,SIGS)
611 C...Multiply cross-section by user-defined weights.
612 IF(MSTP(173).EQ.1) THEN
615 SIGH(ICHN)=PARP(173)*SIGH(ICHN)
617 SIGLPT=PARP(173)*SIGLPT
623 IF(MINT(82).EQ.1.AND.MSTP(142).GE.1) THEN
624 IF(ISUB.NE.96.AND.MSUB(91)+MSUB(92)+MSUB(93)+MSUB(94)+
625 & MSUB(95).EQ.0) CALL PYEVWT(WTXS)
628 IF(MSTP(142).EQ.1) VINT(100)=1./WTXS
631 C...Calculations for Monte Carlo estimate of all cross-sections.
632 IF(MINT(82).EQ.1.AND.ISUB.LE.90.OR.ISUB.GE.96) THEN
633 IF(MSTP(142).LE.1) THEN
634 XSEC(ISUB,2)=XSEC(ISUB,2)+SIGS
635 DXSEC(ISUB)=DXSEC(ISUB)+SIGS
637 XSEC(ISUB,2)=XSEC(ISUB,2)+SIGSWT
638 DXSEC(ISUB)=DXSEC(ISUB)+SIGSWT
640 ELSEIF(MINT(82).EQ.1) THEN
641 XSEC(ISUB,2)=XSEC(ISUB,2)+SIGS
642 DXSEC(ISUB)=DXSEC(ISUB)+SIGS
644 IF((ISUB.EQ.95.OR.ISUB.EQ.96).AND.LOOP.EQ.1.AND.MINT(82).EQ.1)
646 XSEC(97,2)=XSEC(97,2)+SIGLPT
647 DXSEC(97)=DXSEC(97)+SIGLPT
650 C...Multiple interactions: store results of cross-section calculation.
651 IF(MINT(50).EQ.1.AND.MSTP(82).GE.3) THEN
656 C...Check that weight not negative.
657 VIOL=SIGSWT/XSEC(ISUB,1)
658 IF(ISUB.EQ.96.AND.MSTP(173).EQ.1) VIOL=VIOL/PARP(174)
659 IF(MSTP(123).LE.0) THEN
660 IF(VIOL.LT.-1E-3) THEN
661 WRITE(MSTU(11),5000) VIOL,NGEN(0,3)+1
662 WRITE(MSTU(11),5100) ISUB,VINT(21),VINT(22),VINT(23),VINT(26)
666 IF(VIOL.LT.MIN(-1E-3,VINT(109))) THEN
668 WRITE(MSTU(11),5200) VIOL,NGEN(0,3)+1
669 WRITE(MSTU(11),5100) ISUB,VINT(21),VINT(22),VINT(23),VINT(26)
673 C...Weighting using estimate of maximum of differential cross-section.
675 IF(VIOL.LT.RLU(0)) THEN
676 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
679 ELSEIF(ISUB.NE.95.AND.ISUB.NE.96) THEN
680 IF(VIOL.LT.RLU(0)) THEN
682 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
686 RATND=SIGLPT/XSEC(95,1)
687 IF(LOOP.EQ.1.AND.RATND.LT.RLU(0)) THEN
689 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
693 IF(VIOL.LT.RLU(0)) THEN
694 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
699 C...Check for possible violation of estimated maximum of differential
700 C...cross-section used in weighting.
701 IF(MSTP(123).LE.0) THEN
703 WRITE(MSTU(11),5300) VIOL,NGEN(0,3)+1
704 WRITE(MSTU(11),5100) ISUB,VINT(21),VINT(22),VINT(23),VINT(26)
707 ELSEIF(MSTP(123).EQ.1) THEN
708 IF(VIOL.GT.VINT(108)) THEN
712 WRITE(MSTU(11),5400) VIOL,NGEN(0,3)+1
713 WRITE(MSTU(11),5100) ISUB,VINT(21),VINT(22),VINT(23),
717 ELSEIF(VIOL.GT.VINT(108)) THEN
721 XDIF=XSEC(ISUB,1)*(VIOL-1.)
722 XSEC(ISUB,1)=XSEC(ISUB,1)+XDIF
723 IF(MSUB(ISUB).EQ.1.AND.(ISUB.LE.90.OR.ISUB.GT.96))
724 & XSEC(0,1)=XSEC(0,1)+XDIF
725 WRITE(MSTU(11),5400) VIOL,NGEN(0,3)+1
726 WRITE(MSTU(11),5100) ISUB,VINT(21),VINT(22),VINT(23),VINT(26)
728 WRITE(MSTU(11),5500) ISUB,XSEC(ISUB,1)
729 ELSEIF(ISUB.LE.99) THEN
730 WRITE(MSTU(11),5600) ISUB,XSEC(ISUB,1)
732 WRITE(MSTU(11),5700) ISUB,XSEC(ISUB,1)
738 C...Multiple interactions: choose impact parameter.
740 IF(MINT(50).EQ.1.AND.(ISUB.LE.90.OR.ISUB.GE.96).AND.MSTP(82).GE.3)
743 IF(VINT(150).LT.RLU(0)) THEN
744 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
752 IF(MINT(82).EQ.1) NGEN(0,2)=NGEN(0,2)+1
753 IF(MINT(82).EQ.1.AND.MSUB(95).EQ.1) THEN
754 IF(ISUB.LE.90.OR.ISUB.GE.95) NGEN(95,1)=NGEN(95,1)+1
755 IF(ISUB.LE.90.OR.ISUB.GE.96) NGEN(96,2)=NGEN(96,2)+1
757 IF(ISUB.LE.90.OR.ISUB.GE.96) MINT(31)=MINT(31)+1
759 C...Choose flavour of reacting partons (and subprocess).
760 IF(ISTSB.GE.11) GOTO 280
763 RQQBAR=PARP(87)*(1.-(QT2/(QT2+(PARP(88)*PARP(82))**2))**2)
764 IF(ISUB.NE.95.AND.(ISUB.NE.96.OR.MSTP(82).LE.1.OR.
765 &RLU(0).GT.RQQBAR)) THEN
770 RSIGS=RSIGS-SIGH(ICHN)
771 IF(RSIGS.LE.0.) GOTO 280
774 C...Multiple interactions: choose qq~ preferentially at small pT.
775 ELSEIF(ISUB.EQ.96) THEN
778 CALL PYSPLI(MINT(11),21,KFL1,KFLDUM)
781 CALL PYSPLI(MINT(12),21,KFL2,KFLDUM)
784 IF(KFL1.EQ.KFL2.AND.RLU(0).LT.0.5) MINT(2)=2
786 C...Low-pT: choose string drawing configuration.
792 IF(RSIGS.GT.1.) MINT(2)=2
793 IF(RSIGS.GT.2.) MINT(2)=3
796 C...Reassign QCD process. Partons before initial state radiation.
797 280 IF(MINT(2).GT.10) THEN
799 MINT(2)=MOD(MINT(2),10)
801 IF(MINT(82).EQ.1.AND.MSTP(111).GE.0) NGEN(MINT(1),2)=
812 C...Calculate x value of photon for parton inside photon inside e.
817 IF(JT.EQ.1.AND.MINT(43).LE.2) MSPLI=1
818 IF(JT.EQ.2.AND.MOD(MINT(43),2).EQ.1) MSPLI=1
819 IF(IABS(MINT(14+JT)).LE.8.OR.MINT(14+JT).EQ.21) MSPLI=MSPLI+1
824 MINT(105)=MINT(102+JT)
825 MINT(109)=MINT(106+JT)
826 IF(MSTP(57).LE.1) THEN
827 CALL PYSTFU(22,XHRD,Q2HRD,XPQ)
829 CALL PYSTFL(22,XHRD,Q2HRD,XPQ)
832 IF(MSTP(13).EQ.2) THEN
833 Q2PMS=Q2HRD/PMAS(11,1)**2
834 WTMX=WTMX*LOG(MAX(2.,Q2PMS*(1.-XHRD)/XHRD**2))
837 XG=MIN(0.999999,XHRD/XE)
838 IF(MSTP(57).LE.1) THEN
839 CALL PYSTFU(22,XG,Q2HRD,XPQ)
841 CALL PYSTFL(22,XG,Q2HRD,XPQ)
843 WT=(1.+(1.-XE)**2)*XPQ(KFLH)
844 IF(MSTP(13).EQ.2) WT=WT*LOG(MAX(2.,Q2PMS*(1.-XE)/XE**2))
845 IF(WT.LT.RLU(0)*WTMX) GOTO 290
849 XSFX(JT,KFLS)=XPQ(KFLS)
854 C...Pick scale where photon is resolved.
855 IF(MINT(107).EQ.3) VINT(283)=PARP(15)**2*
856 &(VINT(54)/PARP(15)**2)**RLU(0)
857 IF(MINT(108).EQ.3) VINT(284)=PARP(15)**2*
858 &(VINT(54)/PARP(15)**2)**RLU(0)
859 IF(MINT(121).GT.1) CALL PYSAVE(2,IGA)
861 C...Format statements for differential cross-section maximum violations.
862 5000 FORMAT(1X,'Error: negative cross-section fraction',1P,E11.3,1X,
863 &'in event',1X,I7,'.'/1X,'Execution stopped!')
864 5100 FORMAT(1X,'ISUB = ',I3,'; Point of violation:'/1X,'tau =',1P,
865 &E11.3,', y* =',E11.3,', cthe = ',0P,F11.7,', tau'' =',1P,E11.3)
866 5200 FORMAT(1X,'Warning: negative cross-section fraction',1P,E11.3,1X,
868 5300 FORMAT(1X,'Error: maximum violated by',1P,E11.3,1X,
869 &'in event',1X,I7,'.'/1X,'Execution stopped!')
870 5400 FORMAT(1X,'Warning: maximum violated by',1P,E11.3,1X,
872 5500 FORMAT(1X,'XSEC(',I1,',1) increased to',1P,E11.3)
873 5600 FORMAT(1X,'XSEC(',I2,',1) increased to',1P,E11.3)
874 5700 FORMAT(1X,'XSEC(',I3,',1) increased to',1P,E11.3)