hijing/src/luzdis.f

0001
0002 C*********************************************************************
0003
0004       SUBROUTINE LUZDIS(KFL1,KFL2,PR,Z)
0005
0006 C...Purpose: to generate the longitudinal splitting variable z.
0007       COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0008       SAVE /LUDAT1/
0009
0010 C...Check if heavy flavour fragmentation.
0011       KFLA=IABS(KFL1)
0012       KFLB=IABS(KFL2)
0013       KFLH=KFLA
0014       IF(KFLA.GE.10) KFLH=MOD(KFLA/1000,10)
0015
0016 C...Lund symmetric scaling function: determine parameters of shape.
0017       IF(MSTJ(11).EQ.1.OR.(MSTJ(11).EQ.3.AND.KFLH.LE.3)) THEN
0018         FA=PARJ(41)
0019         IF(MSTJ(91).EQ.1) FA=PARJ(43)
0020         IF(KFLB.GE.10) FA=FA+PARJ(45)
0021         FB=PARJ(42)*PR
0022         IF(MSTJ(91).EQ.1) FB=PARJ(44)*PR
0023         FC=1.
0024         IF(KFLA.GE.10) FC=FC-PARJ(45)
0025         IF(KFLB.GE.10) FC=FC+PARJ(45)
0026         MC=1
0027         IF(ABS(FC-1.).GT.0.01) MC=2
0028
0029 C...Determine position of maximum. Special cases for a = 0 or a = c.
0030         IF(FA.LT.0.02) THEN
0031           MA=1
0032           ZMAX=1.
0033           IF(FC.GT.FB) ZMAX=FB/FC
0034         ELSEIF(ABS(FC-FA).LT.0.01) THEN
0035           MA=2
0036           ZMAX=FB/(FB+FC)
0037         ELSE
0038           MA=3
0039           ZMAX=0.5*(FB+FC-SQRT((FB-FC)**2+4.*FA*FB))/(FC-FA)
0040           IF(ZMAX.GT.0.99.AND.FB.GT.100.) ZMAX=1.-FA/FB
0041         ENDIF
0042
0043 C...Subdivide z range if distribution very peaked near endpoint.
0044         MMAX=2
0045         IF(ZMAX.LT.0.1) THEN
0046           MMAX=1
0047           ZDIV=2.75*ZMAX
0048           IF(MC.EQ.1) THEN
0049             FINT=1.-LOG(ZDIV)
0050           ELSE
0051             ZDIVC=ZDIV**(1.-FC)
0052             FINT=1.+(1.-1./ZDIVC)/(FC-1.)
0053           ENDIF
0054         ELSEIF(ZMAX.GT.0.85.AND.FB.GT.1.) THEN
0055           MMAX=3
0056           FSCB=SQRT(4.+(FC/FB)**2)
0057           ZDIV=FSCB-1./ZMAX-(FC/FB)*LOG(ZMAX*0.5*(FSCB+FC/FB))
0058           IF(MA.GE.2) ZDIV=ZDIV+(FA/FB)*LOG(1.-ZMAX)
0059           ZDIV=MIN(ZMAX,MAX(0.,ZDIV))
0060           FINT=1.+FB*(1.-ZDIV)
0061         ENDIF
0062
0063 C...Choice of z, preweighted for peaks at low or high z.
0064   100   Z=RLU(0)
0065         FPRE=1.
0066         IF(MMAX.EQ.1) THEN
0067           IF(FINT*RLU(0).LE.1.) THEN
0068             Z=ZDIV*Z
0069           ELSEIF(MC.EQ.1) THEN
0070             Z=ZDIV**Z
0071             FPRE=ZDIV/Z
0072           ELSE
0073             Z=1./(ZDIVC+Z*(1.-ZDIVC))**(1./(1.-FC))
0074             FPRE=(ZDIV/Z)**FC
0075           ENDIF
0076         ELSEIF(MMAX.EQ.3) THEN
0077           IF(FINT*RLU(0).LE.1.) THEN
0078             Z=ZDIV+LOG(Z)/FB
0079             FPRE=EXP(FB*(Z-ZDIV))
0080           ELSE
0081             Z=ZDIV+Z*(1.-ZDIV)
0082           ENDIF
0083         ENDIF
0084
0085 C...Weighting according to correct formula.
0086         IF(Z.LE.FB/(50.+FB).OR.Z.GE.1.) GOTO 100
0087         FVAL=(ZMAX/Z)**FC*EXP(FB*(1./ZMAX-1./Z))
0088         IF(MA.GE.2) FVAL=((1.-Z)/(1.-ZMAX))**FA*FVAL
0089         IF(FVAL.LT.RLU(0)*FPRE) GOTO 100
0090
0091 C...Generate z according to Field-Feynman, SLAC, (1-z)**c OR z**c.
0092       ELSE
0093         FC=PARJ(50+MAX(1,KFLH))
0094         IF(MSTJ(91).EQ.1) FC=PARJ(59)
0095   110   Z=RLU(0)
0096         IF(FC.GE.0..AND.FC.LE.1.) THEN
0097           IF(FC.GT.RLU(0)) Z=1.-Z**(1./3.)
0098         ELSEIF(FC.GT.-1.) THEN
0099           IF(-4.*FC*Z*(1.-Z)**2.LT.RLU(0)*((1.-Z)**2-FC*Z)**2) GOTO 110
0100         ELSE
0101           IF(FC.GT.0.) Z=1.-Z**(1./FC)
0102           IF(FC.LT.0.) Z=Z**(-1./FC)
0103         ENDIF
0104       ENDIF
0105
0106       RETURN
0107       END