Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:15:44

0001     
0002 C*********************************************************************  
0003     
0004       SUBROUTINE LUX3JT(NJET,CUT,KFL,ECM,X1,X2) 
0005     
0006 C...Purpose: to select the kinematical variables of three-jet events.   
0007       COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200) 
0008       SAVE /LUDAT1/ 
0009       DIMENSION ZHUP(5,12)  
0010     
0011 C...Coefficients of Zhu second order parametrization.   
0012       DATA ((ZHUP(IC1,IC2),IC2=1,12),IC1=1,5)/  
0013      &    18.29,    89.56,    4.541,   -52.09,   -109.8,    24.90,  
0014      &    11.63,    3.683,    17.50, 0.002440,   -1.362,  -0.3537,  
0015      &    11.42,    6.299,   -22.55,   -8.915,    59.25,   -5.855,  
0016      &   -32.85,   -1.054,   -16.90, 0.006489,  -0.8156,  0.01095,  
0017      &    7.847,   -3.964,   -35.83,    1.178,    29.39,   0.2806,  
0018      &    47.82,   -12.36,   -56.72,  0.04054,  -0.4365,   0.6062,  
0019      &    5.441,   -56.89,   -50.27,    15.13,    114.3,   -18.19,  
0020      &    97.05,   -1.890,   -139.9,  0.08153,  -0.4984,   0.9439,  
0021      &   -17.65,    51.44,   -58.32,    70.95,   -255.7,   -78.99,  
0022      &    476.9,    29.65,   -239.3,   0.4745,   -1.174,    6.081/  
0023     
0024 C...Dilogarithm of x for x<0.5 (x>0.5 obtained by analytic trick).  
0025       DILOG(X)=X+X**2/4.+X**3/9.+X**4/16.+X**5/25.+X**6/36.+X**7/49.    
0026     
0027 C...Event type. Mass effect factors and other common constants. 
0028       MSTJ(120)=2   
0029       MSTJ(121)=0   
0030       PMQ=ULMASS(KFL)   
0031       QME=(2.*PMQ/ECM)**2   
0032       IF(MSTJ(109).NE.1) THEN   
0033         CUTL=LOG(CUT)   
0034         CUTD=LOG(1./CUT-2.) 
0035         IF(MSTJ(109).EQ.0) THEN 
0036           CF=4./3.  
0037           CN=3. 
0038           TR=2. 
0039           WTMX=MIN(20.,37.-6.*CUTD) 
0040           IF(MSTJ(110).EQ.2) WTMX=2.*(7.5+80.*CUT)  
0041         ELSE    
0042           CF=1. 
0043           CN=0. 
0044           TR=12.    
0045           WTMX=0.   
0046         ENDIF   
0047     
0048 C...Alpha_strong and effects of optimized Q^2 scale. Maximum weight.    
0049         ALS2PI=PARU(118)/PARU(2)    
0050         WTOPT=0.    
0051         IF(MSTJ(111).EQ.1) WTOPT=(33.-2.*MSTU(112))/6.*LOG(PARJ(169))*  
0052      &  ALS2PI  
0053         WTMAX=MAX(0.,1.+WTOPT+ALS2PI*WTMX)  
0054     
0055 C...Choose three-jet events in allowed region.  
0056   100   NJET=3  
0057   110   Y13L=CUTL+CUTD*RLU(0)   
0058         Y23L=CUTL+CUTD*RLU(0)   
0059         Y13=EXP(Y13L)   
0060         Y23=EXP(Y23L)   
0061         Y12=1.-Y13-Y23  
0062         IF(Y12.LE.CUT) GOTO 110 
0063         IF(Y13**2+Y23**2+2.*Y12.LE.2.*RLU(0)) GOTO 110  
0064     
0065 C...Second order corrections.   
0066         IF(MSTJ(101).EQ.2.AND.MSTJ(110).LE.1) THEN  
0067           Y12L=LOG(Y12) 
0068           Y13M=LOG(1.-Y13)  
0069           Y23M=LOG(1.-Y23)  
0070           Y12M=LOG(1.-Y12)  
0071           IF(Y13.LE.0.5) Y13I=DILOG(Y13)    
0072           IF(Y13.GE.0.5) Y13I=1.644934-Y13L*Y13M-DILOG(1.-Y13)  
0073           IF(Y23.LE.0.5) Y23I=DILOG(Y23)    
0074           IF(Y23.GE.0.5) Y23I=1.644934-Y23L*Y23M-DILOG(1.-Y23)  
0075           IF(Y12.LE.0.5) Y12I=DILOG(Y12)    
0076           IF(Y12.GE.0.5) Y12I=1.644934-Y12L*Y12M-DILOG(1.-Y12)  
0077           WT1=(Y13**2+Y23**2+2.*Y12)/(Y13*Y23)  
0078           WT2=CF*(-2.*(CUTL-Y12L)**2-3.*CUTL-1.+3.289868+   
0079      &    2.*(2.*CUTL-Y12L)*CUT/Y12)+   
0080      &    CN*((CUTL-Y12L)**2-(CUTL-Y13L)**2-(CUTL-Y23L)**2-11.*CUTL/6.+ 
0081      &    67./18.+1.644934-(2.*CUTL-Y12L)*CUT/Y12+(2.*CUTL-Y13L)*   
0082      &    CUT/Y13+(2.*CUTL-Y23L)*CUT/Y23)+  
0083      &    TR*(2.*CUTL/3.-10./9.)+   
0084      &    CF*(Y12/(Y12+Y13)+Y12/(Y12+Y23)+(Y12+Y23)/Y13+(Y12+Y13)/Y23+  
0085      &    Y13L*(4.*Y12**2+2.*Y12*Y13+4.*Y12*Y23+Y13*Y23)/(Y12+Y23)**2+  
0086      &    Y23L*(4.*Y12**2+2.*Y12*Y23+4.*Y12*Y13+Y13*Y23)/(Y12+Y13)**2)/ 
0087      &    WT1+  
0088      &    CN*(Y13L*Y13/(Y12+Y23)+Y23L*Y23/(Y12+Y13))/WT1+   
0089      &    (CN-2.*CF)*((Y12**2+(Y12+Y13)**2)*(Y12L*Y23L-Y12L*Y12M-Y23L*  
0090      &    Y23M+1.644934-Y12I-Y23I)/(Y13*Y23)+(Y12**2+(Y12+Y23)**2)* 
0091      &    (Y12L*Y13L-Y12L*Y12M-Y13L*Y13M+1.644934-Y12I-Y13I)/   
0092      &    (Y13*Y23)+(Y13**2+Y23**2)/(Y13*Y23*(Y13+Y23))-    
0093      &    2.*Y12L*Y12**2/(Y13+Y23)**2-4.*Y12L*Y12/(Y13+Y23))/WT1-   
0094      &    CN*(Y13L*Y23L-Y13L*Y13M-Y23L*Y23M+1.644934-Y13I-Y23I) 
0095           IF(1.+WTOPT+ALS2PI*WT2.LE.0.) MSTJ(121)=1 
0096           IF(1.+WTOPT+ALS2PI*WT2.LE.WTMAX*RLU(0)) GOTO 110  
0097           PARJ(156)=(WTOPT+ALS2PI*WT2)/(1.+WTOPT+ALS2PI*WT2)    
0098     
0099         ELSEIF(MSTJ(101).EQ.2.AND.MSTJ(110).EQ.2) THEN  
0100 C...Second order corrections; Zhu parametrization of ERT.   
0101           ZX=(Y23-Y13)**2   
0102           ZY=1.-Y12 
0103           IZA=0 
0104           DO 120 IY=1,5 
0105   120     IF(ABS(CUT-0.01*IY).LT.0.0001) IZA=IY 
0106           IF(IZA.NE.0) THEN 
0107             IZ=IZA  
0108             WT2=ZHUP(IZ,1)+ZHUP(IZ,2)*ZX+ZHUP(IZ,3)*ZX**2+(ZHUP(IZ,4)+  
0109      &      ZHUP(IZ,5)*ZX)*ZY+(ZHUP(IZ,6)+ZHUP(IZ,7)*ZX)*ZY**2+ 
0110      &      (ZHUP(IZ,8)+ZHUP(IZ,9)*ZX)*ZY**3+ZHUP(IZ,10)/(ZX-ZY**2)+    
0111      &      ZHUP(IZ,11)/(1.-ZY)+ZHUP(IZ,12)/ZY  
0112           ELSE  
0113             IZ=100.*CUT 
0114             WTL=ZHUP(IZ,1)+ZHUP(IZ,2)*ZX+ZHUP(IZ,3)*ZX**2+(ZHUP(IZ,4)+  
0115      &      ZHUP(IZ,5)*ZX)*ZY+(ZHUP(IZ,6)+ZHUP(IZ,7)*ZX)*ZY**2+ 
0116      &      (ZHUP(IZ,8)+ZHUP(IZ,9)*ZX)*ZY**3+ZHUP(IZ,10)/(ZX-ZY**2)+    
0117      &      ZHUP(IZ,11)/(1.-ZY)+ZHUP(IZ,12)/ZY  
0118             IZ=IZ+1 
0119             WTU=ZHUP(IZ,1)+ZHUP(IZ,2)*ZX+ZHUP(IZ,3)*ZX**2+(ZHUP(IZ,4)+  
0120      &      ZHUP(IZ,5)*ZX)*ZY+(ZHUP(IZ,6)+ZHUP(IZ,7)*ZX)*ZY**2+ 
0121      &      (ZHUP(IZ,8)+ZHUP(IZ,9)*ZX)*ZY**3+ZHUP(IZ,10)/(ZX-ZY**2)+    
0122      &      ZHUP(IZ,11)/(1.-ZY)+ZHUP(IZ,12)/ZY  
0123             WT2=WTL+(WTU-WTL)*(100.*CUT+1.-IZ)  
0124           ENDIF 
0125           IF(1.+WTOPT+2.*ALS2PI*WT2.LE.0.) MSTJ(121)=1  
0126           IF(1.+WTOPT+2.*ALS2PI*WT2.LE.WTMAX*RLU(0)) GOTO 110   
0127           PARJ(156)=(WTOPT+2.*ALS2PI*WT2)/(1.+WTOPT+2.*ALS2PI*WT2)  
0128         ENDIF   
0129     
0130 C...Impose mass cuts (gives two jets). For fixed jet number new try.    
0131         X1=1.-Y23   
0132         X2=1.-Y13   
0133         X3=1.-Y12   
0134         IF(4.*Y23*Y13*Y12/X3**2.LE.QME) NJET=2  
0135         IF(MOD(MSTJ(103),4).GE.2.AND.IABS(MSTJ(101)).LE.1.AND.QME*X3+   
0136      &  0.5*QME**2+(0.5*QME+0.25*QME**2)*((1.-X2)/(1.-X1)+  
0137      &  (1.-X1)/(1.-X2)).GT.(X1**2+X2**2)*RLU(0)) NJET=2    
0138         IF(MSTJ(101).EQ.-1.AND.NJET.EQ.2) GOTO 100  
0139     
0140 C...Scalar gluon model (first order only, no mass effects). 
0141       ELSE  
0142   130   NJET=3  
0143   140   Y12=SQRT(4.*CUT**2+RLU(0)*((1.-CUT)**2-4.*CUT**2))  
0144         IF(LOG((Y12-CUT)/CUT).LE.RLU(0)*LOG((1.-2.*CUT)/CUT)) GOTO 140  
0145         YD=SIGN(2.*CUT*((Y12-CUT)/CUT)**RLU(0)-Y12,RLU(0)-0.5)  
0146         X1=1.-0.5*(Y12+YD)  
0147         X2=1.-0.5*(Y12-YD)  
0148         IF(4.*(1.-X1)*(1.-X2)*Y12/(1.-Y12)**2.LE.QME) NJET=2    
0149         IF(MSTJ(101).EQ.-1.AND.NJET.EQ.2) GOTO 130  
0150       ENDIF 
0151     
0152       RETURN    
0153       END