pythiaeRHIC/pythia/pypdfl.f

0001
0002 C*********************************************************************
0003
0004 C...PYPDFL
0005 C...Gives proton parton distribution at small x and/or Q^2 according to
0006 C...correct limiting behaviour.
0007
0008       SUBROUTINE PYPDFL(KF,X,Q2,XPQ)
0009
0010 C...Double precision and integer declarations.
0011       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
0012       IMPLICIT INTEGER(I-N)
0013       INTEGER PYK,PYCHGE,PYCOMP
0014 C...Commonblocks.
0015       COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0016       COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
0017       COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
0018       COMMON/PYINT1/MINT(400),VINT(400)
0019       SAVE /PYDAT1/,/PYDAT2/,/PYPARS/,/PYINT1/
0020 C...Local arrays.
0021       DIMENSION XPQ(-25:25),XPA(-25:25),XPB(-25:25),WTSB(-3:3)
0022       DATA RMR/0.92D0/,RMP/0.38D0/,WTSB/0.5D0,1D0,1D0,5D0,1D0,1D0,0.5D0/
0023
0024 C...Send everything but protons/neutrons/VMD pions directly to PYPDFU.
0025       MINT(92)=0
0026       KFA=IABS(KF)
0027       IACC=0
0028       IF((KFA.EQ.2212.OR.KFA.EQ.2112).AND.MSTP(57).GE.2) IACC=1
0029       IF(KFA.EQ.211.AND.MSTP(57).GE.3) IACC=1
0030       IF(KFA.EQ.22.AND.MINT(109).EQ.2.AND.MSTP(57).GE.3) IACC=1
0031       IF(IACC.EQ.0) THEN
0032         CALL PYPDFU(KF,X,Q2,XPQ)
0033         RETURN
0034       ENDIF
0035
0036 C...Reset. Check x.
0037       DO 100 KFL=-25,25
0038         XPQ(KFL)=0D0
0039   100 CONTINUE
0040       IF(X.LE.0D0.OR.X.GE.1D0) THEN
0041         WRITE(MSTU(11),5000) X
0042         RETURN
0043       ENDIF
0044
0045 C...Define valence content.
0046       KFC=KF
0047       NV1=2
0048       NV2=1
0049       IF(KF.EQ.2212) THEN
0050         KFV1=2
0051         KFV2=1
0052       ELSEIF(KF.EQ.-2212) THEN
0053         KFV1=-2
0054         KFV2=-1
0055       ELSEIF(KF.EQ.2112) THEN
0056         KFV1=1
0057         KFV2=2
0058       ELSEIF(KF.EQ.-2112) THEN
0059         KFV1=-1
0060         KFV2=-2
0061       ELSEIF(KF.EQ.211) THEN
0062         NV1=1
0063         KFV1=2
0064         KFV2=-1
0065       ELSEIF(KF.EQ.-211) THEN
0066         NV1=1
0067         KFV1=-2
0068         KFV2=1
0069       ELSEIF(MINT(105).LE.223) THEN
0070         KFV1=1
0071         WTV1=0.2D0
0072         KFV2=2
0073         WTV2=0.8D0
0074       ELSEIF(MINT(105).EQ.333) THEN
0075         KFV1=3
0076         WTV1=1.0D0
0077         KFV2=1
0078         WTV2=0.0D0
0079       ELSEIF(MINT(105).EQ.443) THEN
0080         KFV1=4
0081         WTV1=1.0D0
0082         KFV2=1
0083         WTV2=0.0D0
0084       ENDIF
0085
0086 C...Do naive evaluation and find min Q^2, boundary Q^2 and x_0.
0087       MINT30=MINT(30)
0088       CALL PYPDFU(KFC,X,Q2,XPA)
0089       Q2MN=MAX(3D0,VINT(231))
0090       Q2B=2D0+0.052D0**2*EXP(3.56D0*SQRT(MAX(0D0,-LOG(3D0*X))))
0091       XMN=EXP(-(LOG((Q2MN-2D0)/0.052D0**2)/3.56D0)**2)/3D0
0092
0093 C...Large Q2 and large x: naive call is enough.
0094       IF(Q2.GT.Q2MN.AND.Q2.GT.Q2B) THEN
0095         DO 110 KFL=-25,25
0096           XPQ(KFL)=XPA(KFL)
0097   110   CONTINUE
0098         MINT(92)=1
0099
0100 C...Small Q2 and large x: dampen boundary value.
0101       ELSEIF(X.GT.XMN) THEN
0102
0103 C...Evaluate at boundary and define dampening factors.
0104         MINT(30)=MINT30
0105         CALL PYPDFU(KFC,X,Q2MN,XPA)
0106         FV=(Q2*(Q2MN+RMR)/(Q2MN*(Q2+RMR)))**(0.55D0*(1D0-X)/(1D0-XMN))
0107         FS=(Q2*(Q2MN+RMP)/(Q2MN*(Q2+RMP)))**1.08D0
0108
0109 C...Separate valence and sea parts of parton distribution.
0110         IF(KFA.NE.22) THEN
0111           XFV1=XPA(KFV1)-XPA(-KFV1)
0112           XPA(KFV1)=XPA(-KFV1)
0113           XFV2=XPA(KFV2)-XPA(-KFV2)
0114           XPA(KFV2)=XPA(-KFV2)
0115         ELSE
0116           XPA(KFV1)=XPA(KFV1)-WTV1*VINT(232)
0117           XPA(-KFV1)=XPA(-KFV1)-WTV1*VINT(232)
0118           XPA(KFV2)=XPA(KFV2)-WTV2*VINT(232)
0119           XPA(-KFV2)=XPA(-KFV2)-WTV2*VINT(232)
0120         ENDIF
0121
0122 C...Dampen valence and sea separately. Put back together.
0123         DO 120 KFL=-25,25
0124           XPQ(KFL)=FS*XPA(KFL)
0125   120   CONTINUE
0126         IF(KFA.NE.22) THEN
0127           XPQ(KFV1)=XPQ(KFV1)+FV*XFV1
0128           XPQ(KFV2)=XPQ(KFV2)+FV*XFV2
0129         ELSE
0130           XPQ(KFV1)=XPQ(KFV1)+FV*WTV1*VINT(232)
0131           XPQ(-KFV1)=XPQ(-KFV1)+FV*WTV1*VINT(232)
0132           XPQ(KFV2)=XPQ(KFV2)+FV*WTV2*VINT(232)
0133           XPQ(-KFV2)=XPQ(-KFV2)+FV*WTV2*VINT(232)
0134         ENDIF
0135         MINT(92)=2
0136
0137 C...Large Q2 and small x: interpolate behaviour.
0138       ELSEIF(Q2.GT.Q2MN) THEN
0139
0140 C...Evaluate at extremes and define coefficients for interpolation.
0141         MINT(30)=MINT30
0142         CALL PYPDFU(KFC,XMN,Q2MN,XPA)
0143         VI232A=VINT(232)
0144         MINT(30)=MINT30
0145         CALL PYPDFU(KFC,X,Q2B,XPB)
0146         VI232B=VINT(232)
0147         FLA=LOG(Q2B/Q2)/LOG(Q2B/Q2MN)
0148         FVA=(X/XMN)**0.45D0*FLA
0149         FSA=(X/XMN)**(-0.08D0)*FLA
0150         FB=1D0-FLA
0151
0152 C...Separate valence and sea parts of parton distribution.
0153         IF(KFA.NE.22) THEN
0154           XFVA1=XPA(KFV1)-XPA(-KFV1)
0155           XPA(KFV1)=XPA(-KFV1)
0156           XFVA2=XPA(KFV2)-XPA(-KFV2)
0157           XPA(KFV2)=XPA(-KFV2)
0158           XFVB1=XPB(KFV1)-XPB(-KFV1)
0159           XPB(KFV1)=XPB(-KFV1)
0160           XFVB2=XPB(KFV2)-XPB(-KFV2)
0161           XPB(KFV2)=XPB(-KFV2)
0162         ELSE
0163           XPA(KFV1)=XPA(KFV1)-WTV1*VI232A
0164           XPA(-KFV1)=XPA(-KFV1)-WTV1*VI232A
0165           XPA(KFV2)=XPA(KFV2)-WTV2*VI232A
0166           XPA(-KFV2)=XPA(-KFV2)-WTV2*VI232A
0167           XPB(KFV1)=XPB(KFV1)-WTV1*VI232B
0168           XPB(-KFV1)=XPB(-KFV1)-WTV1*VI232B
0169           XPB(KFV2)=XPB(KFV2)-WTV2*VI232B
0170           XPB(-KFV2)=XPB(-KFV2)-WTV2*VI232B
0171         ENDIF
0172
0173 C...Interpolate for valence and sea. Put back together.
0174         DO 130 KFL=-25,25
0175           XPQ(KFL)=FSA*XPA(KFL)+FB*XPB(KFL)
0176   130   CONTINUE
0177         IF(KFA.NE.22) THEN
0178           XPQ(KFV1)=XPQ(KFV1)+(FVA*XFVA1+FB*XFVB1)
0179           XPQ(KFV2)=XPQ(KFV2)+(FVA*XFVA2+FB*XFVB2)
0180         ELSE
0181           XPQ(KFV1)=XPQ(KFV1)+WTV1*(FVA*VI232A+FB*VI232B)
0182           XPQ(-KFV1)=XPQ(-KFV1)+WTV1*(FVA*VI232A+FB*VI232B)
0183           XPQ(KFV2)=XPQ(KFV2)+WTV2*(FVA*VI232A+FB*VI232B)
0184           XPQ(-KFV2)=XPQ(-KFV2)+WTV2*(FVA*VI232A+FB*VI232B)
0185         ENDIF
0186         MINT(92)=3
0187
0188 C...Small Q2 and small x: dampen boundary value and add term.
0189       ELSE
0190
0191 C...Evaluate at boundary and define dampening factors.
0192         MINT(30)=MINT30
0193         CALL PYPDFU(KFC,XMN,Q2MN,XPA)
0194         FB=(XMN-X)*(Q2MN-Q2)/(XMN*Q2MN)
0195         FA=1D0-FB
0196         FVC=(X/XMN)**0.45D0*(Q2/(Q2+RMR))**0.55D0
0197         FVA=FVC*FA*((Q2MN+RMR)/Q2MN)**0.55D0
0198         FVB=FVC*FB*1.10D0*XMN**0.45D0*0.11D0
0199         FSC=(X/XMN)**(-0.08D0)*(Q2/(Q2+RMP))**1.08D0
0200         FSA=FSC*FA*((Q2MN+RMP)/Q2MN)**1.08D0
0201         FSB=FSC*FB*0.21D0*XMN**(-0.08D0)*0.21D0
0202
0203 C...Separate valence and sea parts of parton distribution.
0204         IF(KFA.NE.22) THEN
0205           XFV1=XPA(KFV1)-XPA(-KFV1)
0206           XPA(KFV1)=XPA(-KFV1)
0207           XFV2=XPA(KFV2)-XPA(-KFV2)
0208           XPA(KFV2)=XPA(-KFV2)
0209         ELSE
0210           XPA(KFV1)=XPA(KFV1)-WTV1*VINT(232)
0211           XPA(-KFV1)=XPA(-KFV1)-WTV1*VINT(232)
0212           XPA(KFV2)=XPA(KFV2)-WTV2*VINT(232)
0213           XPA(-KFV2)=XPA(-KFV2)-WTV2*VINT(232)
0214         ENDIF
0215
0216 C...Dampen valence and sea separately. Add constant terms.
0217 C...Put back together.
0218         DO 140 KFL=-25,25
0219           XPQ(KFL)=FSA*XPA(KFL)
0220   140   CONTINUE
0221         IF(KFA.NE.22) THEN
0222           DO 150 KFL=-3,3
0223             XPQ(KFL)=XPQ(KFL)+FSB*WTSB(KFL)
0224   150     CONTINUE
0225           XPQ(KFV1)=XPQ(KFV1)+(FVA*XFV1+FVB*NV1)
0226           XPQ(KFV2)=XPQ(KFV2)+(FVA*XFV2+FVB*NV2)
0227         ELSE
0228           DO 160 KFL=-3,3
0229             XPQ(KFL)=XPQ(KFL)+VINT(281)*FSB*WTSB(KFL)
0230   160     CONTINUE
0231           XPQ(KFV1)=XPQ(KFV1)+WTV1*(FVA*VINT(232)+FVB*VINT(281))
0232           XPQ(-KFV1)=XPQ(-KFV1)+WTV1*(FVA*VINT(232)+FVB*VINT(281))
0233           XPQ(KFV2)=XPQ(KFV2)+WTV2*(FVA*VINT(232)+FVB*VINT(281))
0234           XPQ(-KFV2)=XPQ(-KFV2)+WTV2*(FVA*VINT(232)+FVB*VINT(281))
0235         ENDIF
0236         XPQ(21)=XPQ(0)
0237         MINT(92)=4
0238       ENDIF
0239
0240 C...Format for error printout.
0241  5000 FORMAT(' Error: x value outside physical range; x =',1P,D12.3)
0242
0243       RETURN
0244       END