spdf/sasano.F

   1 #include "pdf/pilot.h"
   2       SUBROUTINE SASANO(KF,X,Q2,P2,ALAM,XPGA,VXPGA)
   3 C...Purpose: to evaluate the parton distributions of the anomalous
   4 C...photon, inhomogeneously evolved from a scale P2 (where it vanishes)
   5 C...to Q2.
   6 C...KF=0 gives the sum over (up to) 5 flavours,
   7 C...KF<0 limits to flavours up to abs(KF),
   8 C...KF>0 is for flavour KF only.
   9 C...ALAM is the 4-flavour Lambda, which is automatically converted
  10 C...to 3- and 5-flavour equivalents as needed.
  11       DIMENSION XPGA(-6:6), VXPGA(-6:6), ALAMSQ(3:5)
  12       DATA PMC/1.3/, PMB/4.6/, AEM/0.007297/, AEM2PI/0.0011614/
  13
  14 C...Reset output.
  15       DO 100 KFL=-6,6
  16       XPGA(KFL)=0.
  17       VXPGA(KFL)=0.
  18   100 CONTINUE
  19       IF(Q2.LE.P2) RETURN
  20       KFA=IABS(KF)
  21
  22 C...Calculate Lambda; protect against unphysical Q2 and P2 input.
  23       ALAMSQ(3)=(ALAM*(PMC/ALAM)**(2./27.))**2
  24       ALAMSQ(4)=ALAM**2
  25       ALAMSQ(5)=(ALAM*(ALAM/PMB)**(2./23.))**2
  26       P2EFF=MAX(P2,1.2*ALAMSQ(3))
  27       IF(KF.EQ.4) P2EFF=MAX(P2EFF,PMC**2)
  28       IF(KF.EQ.5) P2EFF=MAX(P2EFF,PMB**2)
  29       Q2EFF=MAX(Q2,P2EFF)
  30       XL=-LOG(X)
  31
  32 C...Find number of flavours at lower and upper scale.
  33       NFP=4
  34       IF(P2EFF.LT.PMC**2) NFP=3
  35       IF(P2EFF.GT.PMB**2) NFP=5
  36       NFQ=4
  37       IF(Q2EFF.LT.PMC**2) NFQ=3
  38       IF(Q2EFF.GT.PMB**2) NFQ=5
  39
  40 C...Define range of flavour loop.
  41       IF(KF.EQ.0) THEN
  42         KFLMN=1
  43         KFLMX=5
  44       ELSEIF(KF.LT.0) THEN
  45         KFLMN=1
  46         KFLMX=KFA
  47       ELSE
  48         KFLMN=KFA
  49         KFLMX=KFA
  50       ENDIF
  51
  52 C...Loop over flavours the photon can branch into.
  53       DO 110 KFL=KFLMN,KFLMX
  54
  55 C...Light flavours: calculate t range and (approximate) s range.
  56       IF(KFL.LE.3.AND.(KFL.EQ.1.OR.KFL.EQ.KF)) THEN
  57         TDIFF=LOG(Q2EFF/P2EFF)
  58         S=(6./(33.-2.*NFQ))*LOG(LOG(Q2EFF/ALAMSQ(NFQ))/
  59      &  LOG(P2EFF/ALAMSQ(NFQ)))
  60         IF(NFQ.GT.NFP) THEN
  61           Q2DIV=PMB**2
  62           IF(NFQ.EQ.4) Q2DIV=PMC**2
  63           SNFQ=(6./(33.-2.*NFQ))*LOG(LOG(Q2DIV/ALAMSQ(NFQ))/
  64      &    LOG(P2EFF/ALAMSQ(NFQ)))
  65           SNFP=(6./(33.-2.*(NFQ-1)))*LOG(LOG(Q2DIV/ALAMSQ(NFQ-1))/
  66      &    LOG(P2EFF/ALAMSQ(NFQ-1)))
  67           S=S+(LOG(Q2DIV/P2EFF)/LOG(Q2EFF/P2EFF))*(SNFP-SNFQ)
  68         ENDIF
  69         IF(NFQ.EQ.5.AND.NFP.EQ.3) THEN
  70           Q2DIV=PMC**2
  71           SNF4=(6./(33.-2.*4))*LOG(LOG(Q2DIV/ALAMSQ(4))/
  72      &    LOG(P2EFF/ALAMSQ(4)))
  73           SNF3=(6./(33.-2.*3))*LOG(LOG(Q2DIV/ALAMSQ(3))/
  74      &    LOG(P2EFF/ALAMSQ(3)))
  75           S=S+(LOG(Q2DIV/P2EFF)/LOG(Q2EFF/P2EFF))*(SNF3-SNF4)
  76         ENDIF
  77
  78 C...u and s quark do not need a separate treatment when d has been done.
  79       ELSEIF(KFL.EQ.2.OR.KFL.EQ.3) THEN
  80
  81 C...Charm: as above, but only include range above c threshold.
  82       ELSEIF(KFL.EQ.4) THEN
  83         IF(Q2.LE.PMC**2) GOTO 110
  84         P2EFF=MAX(P2EFF,PMC**2)
  85         Q2EFF=MAX(Q2EFF,P2EFF)
  86         TDIFF=LOG(Q2EFF/P2EFF)
  87         S=(6./(33.-2.*NFQ))*LOG(LOG(Q2EFF/ALAMSQ(NFQ))/
  88      &  LOG(P2EFF/ALAMSQ(NFQ)))
  89         IF(NFQ.EQ.5.AND.NFP.EQ.4) THEN
  90           Q2DIV=PMB**2
  91           SNFQ=(6./(33.-2.*NFQ))*LOG(LOG(Q2DIV/ALAMSQ(NFQ))/
  92      &    LOG(P2EFF/ALAMSQ(NFQ)))
  93           SNFP=(6./(33.-2.*(NFQ-1)))*LOG(LOG(Q2DIV/ALAMSQ(NFQ-1))/
  94      &    LOG(P2EFF/ALAMSQ(NFQ-1)))
  95           S=S+(LOG(Q2DIV/P2EFF)/LOG(Q2EFF/P2EFF))*(SNFP-SNFQ)
  96         ENDIF
  97
  98 C...Bottom: as above, but only include range above b threshold.
  99       ELSEIF(KFL.EQ.5) THEN
 100         IF(Q2.LE.PMB**2) GOTO 110
 101         P2EFF=MAX(P2EFF,PMB**2)
 102         Q2EFF=MAX(Q2,P2EFF)
 103         TDIFF=LOG(Q2EFF/P2EFF)
 104         S=(6./(33.-2.*NFQ))*LOG(LOG(Q2EFF/ALAMSQ(NFQ))/
 105      &  LOG(P2EFF/ALAMSQ(NFQ)))
 106       ENDIF
 107
 108 C...Evaluate flavour-dependent prefactor (charge^2 etc.).
 109       CHSQ=1./9.
 110       IF(KFL.EQ.2.OR.KFL.EQ.4) CHSQ=4./9.
 111       FAC=AEM2PI*2.*CHSQ*TDIFF
 112
 113 C...Evaluate parton distributions (normalized to unit momentum sum).
 114       IF(KFL.EQ.1.OR.KFL.EQ.4.OR.KFL.EQ.5.OR.KFL.EQ.KF) THEN
 115         XVAL= ((1.5+2.49*S+26.9*S**2)/(1.+32.3*S**2)*X**2 +
 116      &  (1.5-0.49*S+7.83*S**2)/(1.+7.68*S**2)*(1.-X)**2 +
 117      &  1.5*S/(1.-3.2*S+7.*S**2)*X*(1.-X)) *
 118      &  X**(1./(1.+0.58*S)) * (1.-X**2)**(2.5*S/(1.+10.*S))
 119         XGLU= 2.*S/(1.+4.*S+7.*S**2) *
 120      &  X**(-1.67*S/(1.+2.*S)) * (1.-X**2)**(1.2*S) *
 121      &  ((4.*X**2+7.*X+4.)*(1.-X)/3. - 2.*X*(1.+X)*XL)
 122         XSEA= 0.333*S**2/(1.+4.90*S+4.69*S**2+21.4*S**3) *
 123      &  X**(-1.18*S/(1.+1.22*S)) * (1.-X)**(1.2*S) *
 124      &  ((8.-73.*X+62.*X**2)*(1.-X)/9. + (3.-8.*X**2/3.)*X*XL +
 125      &  (2.*X-1.)*X*XL**2)
 126
 127 C...Threshold factors for c and b sea.
 128         SLL=LOG(LOG(Q2EFF/ALAM**2)/LOG(P2EFF/ALAM**2))
 129         XCHM=0.
 130         IF(Q2.GT.PMC**2.AND.Q2.GT.1.001*P2EFF) THEN
 131           SCH=MAX(0.,LOG(LOG(PMC**2/ALAM**2)/LOG(P2EFF/ALAM**2)))
 132           XCHM=XSEA*(1.-(SCH/SLL)**3)
 133         ENDIF
 134         XBOT=0.
 135         IF(Q2.GT.PMB**2.AND.Q2.GT.1.001*P2EFF) THEN
 136           SBT=MAX(0.,LOG(LOG(PMB**2/ALAM**2)/LOG(P2EFF/ALAM**2)))
 137           XBOT=XSEA*(1.-(SBT/SLL)**3)
 138         ENDIF
 139       ENDIF
 140
 141 C...Add contribution of each valence flavour.
 142       XPGA(0)=XPGA(0)+FAC*XGLU
 143       XPGA(1)=XPGA(1)+FAC*XSEA
 144       XPGA(2)=XPGA(2)+FAC*XSEA
 145       XPGA(3)=XPGA(3)+FAC*XSEA
 146       XPGA(4)=XPGA(4)+FAC*XCHM
 147       XPGA(5)=XPGA(5)+FAC*XBOT
 148       XPGA(KFL)=XPGA(KFL)+FAC*XVAL
 149       VXPGA(KFL)=VXPGA(KFL)+FAC*XVAL
 150   110 CONTINUE
 151       DO 120 KFL=1,5
 152       XPGA(-KFL)=XPGA(KFL)
 153       VXPGA(-KFL)=VXPGA(KFL)
 154   120 CONTINUE
 155
 156       RETURN
 157       END