PDF/spdf/sfwhi6.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.4  2000/09/18 10:02:36  hristov
   6 * Makefile added to PDF8
   7 *
   8 * Revision 1.1.1.2  1996/10/30 08:30:46  cernlib
   9 * Version 7.04
  10 *
  11 * Revision 1.1.1.1  1996/04/12 15:29:46  plothow
  12 * Version 7.01
  13 *
  14 *
  15 #include "pdf/pilot.h"
  16 c-------------------------------------------------------
  17       subroutine SFWHI6(ZX,ZQ,ZUV,ZDV,ZUB,ZDB,ZSB,ZCB,ZGL)
  18 c-------------------------------------------------------
  19 c     WHIT6 parton distribution in the photon
  20 c
  21 c     INPUT:  integer ic  : if ic=0 then qc=0
  22 c                           else qc is calculated
  23 c             DOUBLE PRECISION  Q2  : energy scale Q^2 (GeV^2)
  24 c             DOUBLE PRECISION  x   : energy fraction
  25 c
  26 c     OUTPUT: DOUBLE PRECISION  qu  : up-quark dist.
  27 c             DOUBLE PRECISION  qd  : down- or strange-quark dist.
  28 c             DOUBLE PRECISION  qc  : charm-quark dist.
  29 c             DOUBLE PRECISION  g   : gluon dist.
  30 c-------------------------------------------------------
  31 c     Modified by M.Tanaka on July 22, 1994.
  32 c     The bug pointed out by M.Drees is fixed.
  33 c-------------------------------------------------------
  34 c     Modified by I.Watanabe on July 22, 1994.
  35 c-------------------------------------------------------
  36       implicit none
  37       external WHIT6G
  38 #if defined(CERNLIB_DOUBLE)
  39       DOUBLE PRECISION
  40 #endif
  41 #if defined(CERNLIB_SINGLE)
  42       REAL
  43 #endif
  44      +       ZX,ZQ,ZUV,ZDV,ZUB,ZDB,ZSB,ZCB,ZGL
  45 c arg
  46       integer ic
  47       DOUBLE PRECISION Q2,x
  48       DOUBLE PRECISION qu,qd,qc,g
  49 c const
  50       DOUBLE PRECISION q42it,q52it,lam42,lam52
  51       DOUBLE PRECISION alinv,mc,PI
  52 c local
  53       DOUBLE PRECISION qv,qsea,cv,cs,dcv,dcs
  54       DOUBLE PRECISION A0val,A1val,A2val,Bval,Cval,
  55      $       A0sea,B0sea,BB0sea,C0sea
  56       DOUBLE PRECISION A0dcv,A1dcv,A2dcv,A3dcv,Bdcv,Cdcv
  57       DOUBLE PRECISION Adcs, B0dcs, B1dcs, Cdcs
  58       DOUBLE PRECISION x1,x2,mc2q2
  59       DOUBLE PRECISION s,s2,s3,s4,prsccf,alstpi
  60       DOUBLE PRECISION WHIT6G
  61 c parameters
  62       parameter(lam42=0.16d0, lam52=0.091411319d0)
  63       parameter(Q42IT=4.0d0, Q52IT=100.0d0)
  64       parameter(alinv=137.036d0, mc=1.5d0)
  65       parameter(pi=3.14159265358979323846d0)
  66       common /scale/ s,s2,s3,s4,prsccf
  67 c
  68 c begin
  69       x=ZX
  70       Q2=ZQ*ZQ
  71       ic=1
  72 c
  73       x1=1.0d0-x
  74       x2=x**2
  75       mc2q2=mc**2/Q2
  76 c
  77       if(Q2.lt.100.0d0) then
  78 c  under 100 GeV^2
  79 c
  80 c  set  scale s
  81          if(Q2.lt.4.0d0) then
  82 cccc  for under 4GeV^2 prescription
  83             s=  0.0d0
  84             prsccf =  log(Q2/LAM42)/ log(Q42IT/LAM42)
  85             alstpi = 6.0d0/25.0d0/ log(Q42IT/LAM42)
  86          else
  87             s=   log(  log(Q2/LAM42)/ log(Q42IT/LAM42))
  88             prsccf = 1.0d0
  89             alstpi = 6.0d0/25.0d0/ log(Q2/LAM42)
  90          endif
  91             s2=s**2
  92             s3=s2*s
  93             s4=s2**2
  94 c
  95 cccccc   WHIT6 quark (U100)
  96 c
  97       A0val= 2.540000d+00+s*( 2.000000d+00)+s2*( 7.180000d-01)
  98       A1val= 6.230000d-02+s*(-7.010000d+00)+s2*( 1.251000d-01)
  99       A2val=-1.642000d-01+s*(-4.360000d-01)+s2*( 1.048000d+01)
 100      $           +s3*(-5.200000d+00)
 101       Bval = 6.990000d-01+s*(-2.796000d-02)+s2*(-3.650000d-03)
 102       Cval = 4.420000d-01+s*(-1.255000d+00)+s2*( 1.941000d+00)
 103      $           +s3*(-9.950000d-01)
 104       A0sea= 3.180000d+00+s*( 8.690000d+00)+s2*(-2.287000d+01)
 105      $           +s3*( 1.896000d+01)+s4*(-5.140000d+00)
 106       B0sea=-1.003000d-01+s*( 1.603000d-01)+s2*(-1.037000d+00)
 107      $           +s3*( 9.440000d-01)+s4*(-2.915000d-01)
 108       BB0sea=5.690000d+00+s*( 1.867000d+01)+s2*(-4.670000d+01)
 109      $           +s3*( 5.050000d+01)+s4*(-1.835000d+01)
 110       C0sea= 2.149000d+01+s*(-5.650000d+01)+s2*( 1.293000d+02)
 111      $           +s3*(-1.459000d+02)+s4*( 5.750000d+01)
 112 c
 113          qv  = prsccf/alinv/x*
 114      $         (A0val+A1val*x+A2val*x2) * x**Bval * x1**Cval
 115          qsea= prsccf/alinv/x*
 116      $         A0sea * x**(B0sea+BB0sea*x) * x1**C0sea
 117 c
 118          qu  = qv/3.0d0  + qsea/6.0d0
 119          qu  = qu*x
 120          ZUV=qu
 121          ZUB=qu
 122          qd  = qv/12.0d0 + qsea/6.0d0
 123          qd  = qd*x
 124          ZDV=qd
 125          ZDB=qd
 126          ZSB=qd
 127 c
 128          if((ic.ne.0) .and. (x*(1.0d0+4.0d0*mc2q2).lt.1.0d0)) then
 129             call WHIT6Q(x,mc2q2,cv,cs)
 130             qc = cv/alinv/2.0d0/PI + cs*alstpi
 131             qc  = qc*x
 132             ZCB=qc
 133          else
 134             qc = 0.0d0
 135             ZCB=qc
 136          endif
 137 c
 138          g   = WHIT6G(x,Q2)
 139          g   = g*x
 140          ZGL=g
 141 c
 142       else
 143 c over 100 GeV^2
 144 c
 145 c set scale s
 146          s=   log(  log(Q2/LAM52)/ log(Q52IT/LAM52))
 147          prsccf = 1.0d0
 148          alstpi = 6.0d0/23.0d0/ log(Q2/LAM52)
 149             s2=s**2
 150             s3=s2*s
 151             s4=s2**2
 152 c
 153 cccccc   WHIT6 quark (O100)
 154 c
 155       A0val= 4.270000d+00+s*( 3.096000d+00)+s2*( 1.621000d+00)
 156       A1val=-4.740000d+00+s*(-6.900000d+00)+s2*(-2.439000d+00)
 157       A2val= 2.837000d+00+s*( 6.460000d+00)+s2*( 4.100000d+00)
 158       Bval = 6.780000d-01+s*(-3.940000d-02)+s2*( 1.758000d-02)
 159       Cval = 1.728000d-01+s*(-2.493000d-02)+s2*( 1.451000d-01)
 160       A0sea= 3.340000d+00+s*(-5.610000d+00)+s2*( 5.000000d+01)
 161      $           +s3*(-2.207000d+02)+s4*( 3.028000d+02)
 162       B0sea=-2.402000d-01+s*(-4.090000d-01)+s2*( 2.263000d+00)
 163      $           +s3*(-1.050000d+01)+s4*( 1.487000d+01)
 164       BB0sea=8.790000d+00+s*(-8.860000d+00)+s2*( 1.640000d+02)
 165      $           +s3*(-7.120000d+02)+s4*( 9.730000d+02)
 166       C0sea= 9.160000d+00+s*( 9.290000d+00)+s2*(-2.784000d+02)
 167      $           +s3*( 1.175000d+03)+s4*(-1.592000d+03)
 168 c
 169          qv  = 1.0d0/alinv/x*
 170      $         (A0val+A1val*x+A2val*x2) * x**Bval * x1**Cval
 171          qsea= 1.0d0/alinv/x*
 172      $         A0sea * x**(B0sea+BB0sea*x) * x1**C0sea
 173 c
 174          qu  = qv/3.0d0  + qsea/6.0d0
 175          qu  = qu*x
 176          ZUV=qu
 177          ZUB=qu
 178          qd  = qv/12.0d0 + qsea/6.0d0
 179          qd  = qd*x
 180          ZDV=qd
 181          ZDB=qd
 182          ZSB=qd
 183          g   = WHIT6G(x,Q2)
 184          g   = g*x
 185          ZGL=g
 186 c
 187          if((ic.ne.0) .and. (x*(1.0d0+4.0d0*mc2q2).lt.1.0d0)) then
 188       A0dcv=              s*( 1.219000d-01)+s2*( 6.200000d+00)
 189      $           +s3*(-2.504000d+01)+s4*( 3.098000d+01)
 190       A1dcv=              s*( 1.913000d+00)+s2*(-7.690000d+01)
 191      $           +s3*( 3.180000d+02)+s4*(-3.920000d+02)
 192       A2dcv=              s*(-7.160000d+00)+s2*( 2.503000d+02)
 193      $           +s3*(-1.062000d+03)+s4*( 1.308000d+03)
 194       A3dcv=              s*( 3.190000d+00)+s2*(-2.301000d+02)
 195      $           +s3*( 1.012000d+03)+s4*(-1.250000d+03)
 196       Bdcv = 4.990000d-01+s*( 3.470000d+00)+s2*(-1.526000d+01)
 197      $           +s3*( 1.967000d+01)
 198       Cdcv = 3.290000d-01+s*( 8.240000d+00)+s2*(-3.800000d+01)
 199      $           +s3*( 4.630000d+01)
 200       Adcs =              s*(-4.990000d-02)+s2*( 1.026000d-01)
 201      $           +s3*(-7.870000d-02)
 202       B0dcs=-3.610000d-01+s*(-5.760000d-01)+s2*( 2.257000d-01)
 203       B1dcs= 7.680000d+00+s*(-8.830000d+00)+s2*( 3.880000d+00)
 204       Cdcs = 2.548000d+00+s*( 6.910000d-01)+s2*(-8.700000d+00)
 205      $           +s3*( 1.065000d+01)
 206 c
 207          dcv = 1.0d0/alinv/x*
 208      $         (A0dcv+x*A1dcv+x2*A2dcv+x2*x*A3dcv) * x**Bdcv * x1**Cdcv
 209          dcs = 1.0d0/alinv/x*
 210      $         Adcs * x**(B0dcs+B1dcs*x) * x1**Cdcs
 211 c
 212            call WHIT6Q(x,mc*mc/Q2,cv,cs)
 213            qc = cv/alinv/2.0d0/PI + cs*alstpi + dcs + dcv
 214            qc  = qc*x
 215            ZCB=qc
 216          else
 217            qc = 0.0d0
 218            ZCB=qc
 219          endif
 220       endif
 221 c
 222       return
 223       end