PYTHIA/jetset/lux3jt.F

   1
   2 C*********************************************************************
   3
   4       SUBROUTINE LUX3JT(NJET,CUT,KFL,ECM,X1,X2)
   5
   6 C...Purpose: to select the kinematical variables of three-jet events.
   7       COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
   8       SAVE /LUDAT1/
   9       DIMENSION ZHUP(5,12)
  10
  11 C...Coefficients of Zhu second order parametrization.
  12       DATA ((ZHUP(IC1,IC2),IC2=1,12),IC1=1,5)/
  13      &    18.29,    89.56,    4.541,   -52.09,   -109.8,    24.90,
  14      &    11.63,    3.683,    17.50, 0.002440,   -1.362,  -0.3537,
  15      &    11.42,    6.299,   -22.55,   -8.915,    59.25,   -5.855,
  16      &   -32.85,   -1.054,   -16.90, 0.006489,  -0.8156,  0.01095,
  17      &    7.847,   -3.964,   -35.83,    1.178,    29.39,   0.2806,
  18      &    47.82,   -12.36,   -56.72,  0.04054,  -0.4365,   0.6062,
  19      &    5.441,   -56.89,   -50.27,    15.13,    114.3,   -18.19,
  20      &    97.05,   -1.890,   -139.9,  0.08153,  -0.4984,   0.9439,
  21      &   -17.65,    51.44,   -58.32,    70.95,   -255.7,   -78.99,
  22      &    476.9,    29.65,   -239.3,   0.4745,   -1.174,    6.081/
  23
  24 C...Dilogarithm of x for x<0.5 (x>0.5 obtained by analytic trick).
  25       DILOG(X)=X+X**2/4.+X**3/9.+X**4/16.+X**5/25.+X**6/36.+X**7/49.
  26
  27 C...Event type. Mass effect factors and other common constants.
  28       MSTJ(120)=2
  29       MSTJ(121)=0
  30       PMQ=ULMASS(KFL)
  31       QME=(2.*PMQ/ECM)**2
  32       IF(MSTJ(109).NE.1) THEN
  33         CUTL=LOG(CUT)
  34         CUTD=LOG(1./CUT-2.)
  35         IF(MSTJ(109).EQ.0) THEN
  36           CF=4./3.
  37           CN=3.
  38           TR=2.
  39           WTMX=MIN(20.,37.-6.*CUTD)
  40           IF(MSTJ(110).EQ.2) WTMX=2.*(7.5+80.*CUT)
  41         ELSE
  42           CF=1.
  43           CN=0.
  44           TR=12.
  45           WTMX=0.
  46         ENDIF
  47
  48 C...Alpha_strong and effects of optimized Q^2 scale. Maximum weight.
  49         ALS2PI=PARU(118)/PARU(2)
  50         WTOPT=0.
  51         IF(MSTJ(111).EQ.1) WTOPT=(33.-2.*MSTU(112))/6.*LOG(PARJ(169))*
  52      &  ALS2PI
  53         WTMAX=MAX(0.,1.+WTOPT+ALS2PI*WTMX)
  54
  55 C...Choose three-jet events in allowed region.
  56   100   NJET=3
  57   110   Y13L=CUTL+CUTD*RLU(0)
  58         Y23L=CUTL+CUTD*RLU(0)
  59         Y13=EXP(Y13L)
  60         Y23=EXP(Y23L)
  61         Y12=1.-Y13-Y23
  62         IF(Y12.LE.CUT) GOTO 110
  63         IF(Y13**2+Y23**2+2.*Y12.LE.2.*RLU(0)) GOTO 110
  64
  65 C...Second order corrections.
  66         IF(MSTJ(101).EQ.2.AND.MSTJ(110).LE.1) THEN
  67           Y12L=LOG(Y12)
  68           Y13M=LOG(1.-Y13)
  69           Y23M=LOG(1.-Y23)
  70           Y12M=LOG(1.-Y12)
  71           IF(Y13.LE.0.5) Y13I=DILOG(Y13)
  72           IF(Y13.GE.0.5) Y13I=1.644934-Y13L*Y13M-DILOG(1.-Y13)
  73           IF(Y23.LE.0.5) Y23I=DILOG(Y23)
  74           IF(Y23.GE.0.5) Y23I=1.644934-Y23L*Y23M-DILOG(1.-Y23)
  75           IF(Y12.LE.0.5) Y12I=DILOG(Y12)
  76           IF(Y12.GE.0.5) Y12I=1.644934-Y12L*Y12M-DILOG(1.-Y12)
  77           WT1=(Y13**2+Y23**2+2.*Y12)/(Y13*Y23)
  78           WT2=CF*(-2.*(CUTL-Y12L)**2-3.*CUTL-1.+3.289868+
  79      &    2.*(2.*CUTL-Y12L)*CUT/Y12)+
  80      &    CN*((CUTL-Y12L)**2-(CUTL-Y13L)**2-(CUTL-Y23L)**2-11.*CUTL/6.+
  81      &    67./18.+1.644934-(2.*CUTL-Y12L)*CUT/Y12+(2.*CUTL-Y13L)*
  82      &    CUT/Y13+(2.*CUTL-Y23L)*CUT/Y23)+
  83      &    TR*(2.*CUTL/3.-10./9.)+
  84      &    CF*(Y12/(Y12+Y13)+Y12/(Y12+Y23)+(Y12+Y23)/Y13+(Y12+Y13)/Y23+
  85      &    Y13L*(4.*Y12**2+2.*Y12*Y13+4.*Y12*Y23+Y13*Y23)/(Y12+Y23)**2+
  86      &    Y23L*(4.*Y12**2+2.*Y12*Y23+4.*Y12*Y13+Y13*Y23)/(Y12+Y13)**2)/
  87      &    WT1+
  88      &    CN*(Y13L*Y13/(Y12+Y23)+Y23L*Y23/(Y12+Y13))/WT1+
  89      &    (CN-2.*CF)*((Y12**2+(Y12+Y13)**2)*(Y12L*Y23L-Y12L*Y12M-Y23L*
  90      &    Y23M+1.644934-Y12I-Y23I)/(Y13*Y23)+(Y12**2+(Y12+Y23)**2)*
  91      &    (Y12L*Y13L-Y12L*Y12M-Y13L*Y13M+1.644934-Y12I-Y13I)/
  92      &    (Y13*Y23)+(Y13**2+Y23**2)/(Y13*Y23*(Y13+Y23))-
  93      &    2.*Y12L*Y12**2/(Y13+Y23)**2-4.*Y12L*Y12/(Y13+Y23))/WT1-
  94      &    CN*(Y13L*Y23L-Y13L*Y13M-Y23L*Y23M+1.644934-Y13I-Y23I)
  95           IF(1.+WTOPT+ALS2PI*WT2.LE.0.) MSTJ(121)=1
  96           IF(1.+WTOPT+ALS2PI*WT2.LE.WTMAX*RLU(0)) GOTO 110
  97           PARJ(156)=(WTOPT+ALS2PI*WT2)/(1.+WTOPT+ALS2PI*WT2)
  98
  99         ELSEIF(MSTJ(101).EQ.2.AND.MSTJ(110).EQ.2) THEN
 100 C...Second order corrections; Zhu parametrization of ERT.
 101           ZX=(Y23-Y13)**2
 102           ZY=1.-Y12
 103           IZA=0
 104           DO 120 IY=1,5
 105           IF(ABS(CUT-0.01*IY).LT.0.0001) IZA=IY
 106   120     CONTINUE
 107           IF(IZA.NE.0) THEN
 108             IZ=IZA
 109             WT2=ZHUP(IZ,1)+ZHUP(IZ,2)*ZX+ZHUP(IZ,3)*ZX**2+(ZHUP(IZ,4)+
 110      &      ZHUP(IZ,5)*ZX)*ZY+(ZHUP(IZ,6)+ZHUP(IZ,7)*ZX)*ZY**2+
 111      &      (ZHUP(IZ,8)+ZHUP(IZ,9)*ZX)*ZY**3+ZHUP(IZ,10)/(ZX-ZY**2)+
 112      &      ZHUP(IZ,11)/(1.-ZY)+ZHUP(IZ,12)/ZY
 113           ELSE
 114             IZ=100.*CUT
 115             WTL=ZHUP(IZ,1)+ZHUP(IZ,2)*ZX+ZHUP(IZ,3)*ZX**2+(ZHUP(IZ,4)+
 116      &      ZHUP(IZ,5)*ZX)*ZY+(ZHUP(IZ,6)+ZHUP(IZ,7)*ZX)*ZY**2+
 117      &      (ZHUP(IZ,8)+ZHUP(IZ,9)*ZX)*ZY**3+ZHUP(IZ,10)/(ZX-ZY**2)+
 118      &      ZHUP(IZ,11)/(1.-ZY)+ZHUP(IZ,12)/ZY
 119             IZ=IZ+1
 120             WTU=ZHUP(IZ,1)+ZHUP(IZ,2)*ZX+ZHUP(IZ,3)*ZX**2+(ZHUP(IZ,4)+
 121      &      ZHUP(IZ,5)*ZX)*ZY+(ZHUP(IZ,6)+ZHUP(IZ,7)*ZX)*ZY**2+
 122      &      (ZHUP(IZ,8)+ZHUP(IZ,9)*ZX)*ZY**3+ZHUP(IZ,10)/(ZX-ZY**2)+
 123      &      ZHUP(IZ,11)/(1.-ZY)+ZHUP(IZ,12)/ZY
 124             WT2=WTL+(WTU-WTL)*(100.*CUT+1.-IZ)
 125           ENDIF
 126           IF(1.+WTOPT+2.*ALS2PI*WT2.LE.0.) MSTJ(121)=1
 127           IF(1.+WTOPT+2.*ALS2PI*WT2.LE.WTMAX*RLU(0)) GOTO 110
 128           PARJ(156)=(WTOPT+2.*ALS2PI*WT2)/(1.+WTOPT+2.*ALS2PI*WT2)
 129         ENDIF
 130
 131 C...Impose mass cuts (gives two jets). For fixed jet number new try.
 132         X1=1.-Y23
 133         X2=1.-Y13
 134         X3=1.-Y12
 135         IF(4.*Y23*Y13*Y12/X3**2.LE.QME) NJET=2
 136         IF(MOD(MSTJ(103),4).GE.2.AND.IABS(MSTJ(101)).LE.1.AND.QME*X3+
 137      &  0.5*QME**2+(0.5*QME+0.25*QME**2)*((1.-X2)/(1.-X1)+
 138      &  (1.-X1)/(1.-X2)).GT.(X1**2+X2**2)*RLU(0)) NJET=2
 139         IF(MSTJ(101).EQ.-1.AND.NJET.EQ.2) GOTO 100
 140
 141 C...Scalar gluon model (first order only, no mass effects).
 142       ELSE
 143   130   NJET=3
 144   140   X3=SQRT(4.*CUT**2+RLU(0)*((1.-CUT)**2-4.*CUT**2))
 145         IF(LOG((X3-CUT)/CUT).LE.RLU(0)*LOG((1.-2.*CUT)/CUT)) GOTO 140
 146         YD=SIGN(2.*CUT*((X3-CUT)/CUT)**RLU(0)-X3,RLU(0)-0.5)
 147         X1=1.-0.5*(X3+YD)
 148         X2=1.-0.5*(X3-YD)
 149         IF(4.*(1.-X1)*(1.-X2)*(1.-X3)/X3**2.LE.QME) NJET=2
 150         IF(MSTJ(102).GE.2) THEN
 151           IF(X3**2-2.*(1.+X3)*(1.-X1)*(1.-X2)*PARJ(171).LT.
 152      &    X3**2*RLU(0)) NJET=2
 153         ENDIF
 154         IF(MSTJ(101).EQ.-1.AND.NJET.EQ.2) GOTO 130
 155       ENDIF
 156
 157       RETURN
 158       END