PHOS/shaker/lux3jt.f

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