dd/d07/pyhistfu_8f_source.html

C*********************************************************************


      SUBROUTINE pyhistfu(KF,X,Q2,XPQ,JBT)


C                       *******JBT specifies beam or target of the particle

C...Gives proton and pi+ parton structure functions according to a few

C...different parametrizations. Note that what is coded is x times the

C...probability distribution, i.e. xq(x,Q2) etc.

      common/hiparnt/hipr1(100),ihpr2(50),hint1(100),ihnt2(50)

      SAVE  /hiparnt/

      common/hijcrdn/yp(3,300),yt(3,300)

      SAVE  /hijcrdn/

C                       ********COMMON BLOCK FROM HIJING

      common/ludat1/mstu(200),paru(200),mstj(200),parj(200)

      SAVE  /ludat1/

      common/ludat2/kchg(500,3),pmas(500,4),parf(2000),vckm(4,4)

      SAVE  /ludat2/

      common/pyhipars/mstp(200),parp(200),msti(200),pari(200)

      SAVE  /pyhipars/

      common/pyhiint1/mint(400),vint(400)

      SAVE  /pyhiint1/

      dimension xpq(-6:6),xq(6),tx(6),tt(6),ts(6),nehlq(8,2),

     &cehlq(6,6,2,8,2),cdo(3,6,5,2),cow(3,5,4,2)


      SAVE xq, tx, tt, ts, nehlq, cehlq, cdo, cow             ! Uzhi


C...The following data lines are coefficients needed in the

C...Eichten, Hinchliffe, Lane, Quigg proton structure function

C...parametrizations, see below.

C...Powers of 1-x in different cases.

      DATA nehlq/3,4,7,5,7,7,7,7,3,4,7,6,7,7,7,7/

C...Expansion coefficients for up valence quark distribution.

      DATA (((cehlq(ix,it,nx,1,1),ix=1,6),it=1,6),nx=1,2)/

     1 7.677e-01,-2.087e-01,-3.303e-01,-2.517e-02,-1.570e-02,-1.000e-04,

     2-5.326e-01,-2.661e-01, 3.201e-01, 1.192e-01, 2.434e-02, 7.620e-03,

     3 2.162e-01, 1.881e-01,-8.375e-02,-6.515e-02,-1.743e-02,-5.040e-03,

     4-9.211e-02,-9.952e-02, 1.373e-02, 2.506e-02, 8.770e-03, 2.550e-03,

     5 3.670e-02, 4.409e-02, 9.600e-04,-7.960e-03,-3.420e-03,-1.050e-03,

     6-1.549e-02,-2.026e-02,-3.060e-03, 2.220e-03, 1.240e-03, 4.100e-04,

     1 2.395e-01, 2.905e-01, 9.778e-02, 2.149e-02, 3.440e-03, 5.000e-04,

     2 1.751e-02,-6.090e-03,-2.687e-02,-1.916e-02,-7.970e-03,-2.750e-03,

     3-5.760e-03,-5.040e-03, 1.080e-03, 2.490e-03, 1.530e-03, 7.500e-04,

     4 1.740e-03, 1.960e-03, 3.000e-04,-3.400e-04,-2.900e-04,-1.800e-04,

     5-5.300e-04,-6.400e-04,-1.700e-04, 4.000e-05, 6.000e-05, 4.000e-05,

     6 1.700e-04, 2.200e-04, 8.000e-05, 1.000e-05,-1.000e-05,-1.000e-05/

      DATA (((cehlq(ix,it,nx,1,2),ix=1,6),it=1,6),nx=1,2)/

     1 7.237e-01,-2.189e-01,-2.995e-01,-1.909e-02,-1.477e-02, 2.500e-04,

     2-5.314e-01,-2.425e-01, 3.283e-01, 1.119e-01, 2.223e-02, 7.070e-03,

     3 2.289e-01, 1.890e-01,-9.859e-02,-6.900e-02,-1.747e-02,-5.080e-03,

     4-1.041e-01,-1.084e-01, 2.108e-02, 2.975e-02, 9.830e-03, 2.830e-03,

     5 4.394e-02, 5.116e-02,-1.410e-03,-1.055e-02,-4.230e-03,-1.270e-03,

     6-1.991e-02,-2.539e-02,-2.780e-03, 3.430e-03, 1.720e-03, 5.500e-04,

     1 2.410e-01, 2.884e-01, 9.369e-02, 1.900e-02, 2.530e-03, 2.400e-04,

     2 1.765e-02,-9.220e-03,-3.037e-02,-2.085e-02,-8.440e-03,-2.810e-03,

     3-6.450e-03,-5.260e-03, 1.720e-03, 3.110e-03, 1.830e-03, 8.700e-04,

     4 2.120e-03, 2.320e-03, 2.600e-04,-4.900e-04,-3.900e-04,-2.300e-04,

     5-6.900e-04,-8.200e-04,-2.000e-04, 7.000e-05, 9.000e-05, 6.000e-05,

     6 2.400e-04, 3.100e-04, 1.100e-04, 0.000e+00,-2.000e-05,-2.000e-05/

C...Expansion coefficients for down valence quark distribution.

      DATA (((cehlq(ix,it,nx,2,1),ix=1,6),it=1,6),nx=1,2)/

     1 3.813e-01,-8.090e-02,-1.634e-01,-2.185e-02,-8.430e-03,-6.200e-04,

     2-2.948e-01,-1.435e-01, 1.665e-01, 6.638e-02, 1.473e-02, 4.080e-03,

     3 1.252e-01, 1.042e-01,-4.722e-02,-3.683e-02,-1.038e-02,-2.860e-03,

     4-5.478e-02,-5.678e-02, 8.900e-03, 1.484e-02, 5.340e-03, 1.520e-03,

     5 2.220e-02, 2.567e-02,-3.000e-05,-4.970e-03,-2.160e-03,-6.500e-04,

     6-9.530e-03,-1.204e-02,-1.510e-03, 1.510e-03, 8.300e-04, 2.700e-04,

     1 1.261e-01, 1.354e-01, 3.958e-02, 8.240e-03, 1.660e-03, 4.500e-04,

     2 3.890e-03,-1.159e-02,-1.625e-02,-9.610e-03,-3.710e-03,-1.260e-03,

     3-1.910e-03,-5.600e-04, 1.590e-03, 1.590e-03, 8.400e-04, 3.900e-04,

     4 6.400e-04, 4.900e-04,-1.500e-04,-2.900e-04,-1.800e-04,-1.000e-04,

     5-2.000e-04,-1.900e-04, 0.000e+00, 6.000e-05, 4.000e-05, 3.000e-05,

     6 7.000e-05, 8.000e-05, 2.000e-05,-1.000e-05,-1.000e-05,-1.000e-05/

      DATA (((cehlq(ix,it,nx,2,2),ix=1,6),it=1,6),nx=1,2)/

     1 3.578e-01,-8.622e-02,-1.480e-01,-1.840e-02,-7.820e-03,-4.500e-04,

     2-2.925e-01,-1.304e-01, 1.696e-01, 6.243e-02, 1.353e-02, 3.750e-03,

     3 1.318e-01, 1.041e-01,-5.486e-02,-3.872e-02,-1.038e-02,-2.850e-03,

     4-6.162e-02,-6.143e-02, 1.303e-02, 1.740e-02, 5.940e-03, 1.670e-03,

     5 2.643e-02, 2.957e-02,-1.490e-03,-6.450e-03,-2.630e-03,-7.700e-04,

     6-1.218e-02,-1.497e-02,-1.260e-03, 2.240e-03, 1.120e-03, 3.500e-04,

     1 1.263e-01, 1.334e-01, 3.732e-02, 7.070e-03, 1.260e-03, 3.400e-04,

     2 3.660e-03,-1.357e-02,-1.795e-02,-1.031e-02,-3.880e-03,-1.280e-03,

     3-2.100e-03,-3.600e-04, 2.050e-03, 1.920e-03, 9.800e-04, 4.400e-04,

     4 7.700e-04, 5.400e-04,-2.400e-04,-3.900e-04,-2.400e-04,-1.300e-04,

     5-2.600e-04,-2.300e-04, 2.000e-05, 9.000e-05, 6.000e-05, 4.000e-05,

     6 9.000e-05, 1.000e-04, 2.000e-05,-2.000e-05,-2.000e-05,-1.000e-05/

C...Expansion coefficients for up and down sea quark distributions.

      DATA (((cehlq(ix,it,nx,3,1),ix=1,6),it=1,6),nx=1,2)/

     1 6.870e-02,-6.861e-02, 2.973e-02,-5.400e-03, 3.780e-03,-9.700e-04,

     2-1.802e-02, 1.400e-04, 6.490e-03,-8.540e-03, 1.220e-03,-1.750e-03,

     3-4.650e-03, 1.480e-03,-5.930e-03, 6.000e-04,-1.030e-03,-8.000e-05,

     4 6.440e-03, 2.570e-03, 2.830e-03, 1.150e-03, 7.100e-04, 3.300e-04,

     5-3.930e-03,-2.540e-03,-1.160e-03,-7.700e-04,-3.600e-04,-1.900e-04,

     6 2.340e-03, 1.930e-03, 5.300e-04, 3.700e-04, 1.600e-04, 9.000e-05,

     1 1.014e+00,-1.106e+00, 3.374e-01,-7.444e-02, 8.850e-03,-8.700e-04,

     2 9.233e-01,-1.285e+00, 4.475e-01,-9.786e-02, 1.419e-02,-1.120e-03,

     3 4.888e-02,-1.271e-01, 8.606e-02,-2.608e-02, 4.780e-03,-6.000e-04,

     4-2.691e-02, 4.887e-02,-1.771e-02, 1.620e-03, 2.500e-04,-6.000e-05,

     5 7.040e-03,-1.113e-02, 1.590e-03, 7.000e-04,-2.000e-04, 0.000e+00,

     6-1.710e-03, 2.290e-03, 3.800e-04,-3.500e-04, 4.000e-05, 1.000e-05/

      DATA (((cehlq(ix,it,nx,3,2),ix=1,6),it=1,6),nx=1,2)/

     1 1.008e-01,-7.100e-02, 1.973e-02,-5.710e-03, 2.930e-03,-9.900e-04,

     2-5.271e-02,-1.823e-02, 1.792e-02,-6.580e-03, 1.750e-03,-1.550e-03,

     3 1.220e-02, 1.763e-02,-8.690e-03,-8.800e-04,-1.160e-03,-2.100e-04,

     4-1.190e-03,-7.180e-03, 2.360e-03, 1.890e-03, 7.700e-04, 4.100e-04,

     5-9.100e-04, 2.040e-03,-3.100e-04,-1.050e-03,-4.000e-04,-2.400e-04,

     6 1.190e-03,-1.700e-04,-2.000e-04, 4.200e-04, 1.700e-04, 1.000e-04,

     1 1.081e+00,-1.189e+00, 3.868e-01,-8.617e-02, 1.115e-02,-1.180e-03,

     2 9.917e-01,-1.396e+00, 4.998e-01,-1.159e-01, 1.674e-02,-1.720e-03,

     3 5.099e-02,-1.338e-01, 9.173e-02,-2.885e-02, 5.890e-03,-6.500e-04,

     4-3.178e-02, 5.703e-02,-2.070e-02, 2.440e-03, 1.100e-04,-9.000e-05,

     5 8.970e-03,-1.392e-02, 2.050e-03, 6.500e-04,-2.300e-04, 2.000e-05,

     6-2.340e-03, 3.010e-03, 5.000e-04,-3.900e-04, 6.000e-05, 1.000e-05/

C...Expansion coefficients for gluon distribution.

      DATA (((cehlq(ix,it,nx,4,1),ix=1,6),it=1,6),nx=1,2)/

     1 9.482e-01,-9.578e-01, 1.009e-01,-1.051e-01, 3.456e-02,-3.054e-02,

     2-9.627e-01, 5.379e-01, 3.368e-01,-9.525e-02, 1.488e-02,-2.051e-02,

     3 4.300e-01,-8.306e-02,-3.372e-01, 4.902e-02,-9.160e-03, 1.041e-02,

     4-1.925e-01,-1.790e-02, 2.183e-01, 7.490e-03, 4.140e-03,-1.860e-03,

     5 8.183e-02, 1.926e-02,-1.072e-01,-1.944e-02,-2.770e-03,-5.200e-04,

     6-3.884e-02,-1.234e-02, 5.410e-02, 1.879e-02, 3.350e-03, 1.040e-03,

     1 2.948e+01,-3.902e+01, 1.464e+01,-3.335e+00, 5.054e-01,-5.915e-02,

     2 2.559e+01,-3.955e+01, 1.661e+01,-4.299e+00, 6.904e-01,-8.243e-02,

     3-1.663e+00, 1.176e+00, 1.118e+00,-7.099e-01, 1.948e-01,-2.404e-02,

     4-2.168e-01, 8.170e-01,-7.169e-01, 1.851e-01,-1.924e-02,-3.250e-03,

     5 2.088e-01,-4.355e-01, 2.239e-01,-2.446e-02,-3.620e-03, 1.910e-03,

     6-9.097e-02, 1.601e-01,-5.681e-02,-2.500e-03, 2.580e-03,-4.700e-04/

      DATA (((cehlq(ix,it,nx,4,2),ix=1,6),it=1,6),nx=1,2)/

     1 2.367e+00, 4.453e-01, 3.660e-01, 9.467e-02, 1.341e-01, 1.661e-02,

     2-3.170e+00,-1.795e+00, 3.313e-02,-2.874e-01,-9.827e-02,-7.119e-02,

     3 1.823e+00, 1.457e+00,-2.465e-01, 3.739e-02, 6.090e-03, 1.814e-02,

     4-1.033e+00,-9.827e-01, 2.136e-01, 1.169e-01, 5.001e-02, 1.684e-02,

     5 5.133e-01, 5.259e-01,-1.173e-01,-1.139e-01,-4.988e-02,-2.021e-02,

     6-2.881e-01,-3.145e-01, 5.667e-02, 9.161e-02, 4.568e-02, 1.951e-02,

     1 3.036e+01,-4.062e+01, 1.578e+01,-3.699e+00, 6.020e-01,-7.031e-02,

     2 2.700e+01,-4.167e+01, 1.770e+01,-4.804e+00, 7.862e-01,-1.060e-01,

     3-1.909e+00, 1.357e+00, 1.127e+00,-7.181e-01, 2.232e-01,-2.481e-02,

     4-2.488e-01, 9.781e-01,-8.127e-01, 2.094e-01,-2.997e-02,-4.710e-03,

     5 2.506e-01,-5.427e-01, 2.672e-01,-3.103e-02,-1.800e-03, 2.870e-03,

     6-1.128e-01, 2.087e-01,-6.972e-02,-2.480e-03, 2.630e-03,-8.400e-04/

C...Expansion coefficients for strange sea quark distribution.

      DATA (((cehlq(ix,it,nx,5,1),ix=1,6),it=1,6),nx=1,2)/

     1 4.968e-02,-4.173e-02, 2.102e-02,-3.270e-03, 3.240e-03,-6.700e-04,

     2-6.150e-03,-1.294e-02, 6.740e-03,-6.890e-03, 9.000e-04,-1.510e-03,

     3-8.580e-03, 5.050e-03,-4.900e-03,-1.600e-04,-9.400e-04,-1.500e-04,

     4 7.840e-03, 1.510e-03, 2.220e-03, 1.400e-03, 7.000e-04, 3.500e-04,

     5-4.410e-03,-2.220e-03,-8.900e-04,-8.500e-04,-3.600e-04,-2.000e-04,

     6 2.520e-03, 1.840e-03, 4.100e-04, 3.900e-04, 1.600e-04, 9.000e-05,

     1 9.235e-01,-1.085e+00, 3.464e-01,-7.210e-02, 9.140e-03,-9.100e-04,

     2 9.315e-01,-1.274e+00, 4.512e-01,-9.775e-02, 1.380e-02,-1.310e-03,

     3 4.739e-02,-1.296e-01, 8.482e-02,-2.642e-02, 4.760e-03,-5.700e-04,

     4-2.653e-02, 4.953e-02,-1.735e-02, 1.750e-03, 2.800e-04,-6.000e-05,

     5 6.940e-03,-1.132e-02, 1.480e-03, 6.500e-04,-2.100e-04, 0.000e+00,

     6-1.680e-03, 2.340e-03, 4.200e-04,-3.400e-04, 5.000e-05, 1.000e-05/

      DATA (((cehlq(ix,it,nx,5,2),ix=1,6),it=1,6),nx=1,2)/

     1 6.478e-02,-4.537e-02, 1.643e-02,-3.490e-03, 2.710e-03,-6.700e-04,

     2-2.223e-02,-2.126e-02, 1.247e-02,-6.290e-03, 1.120e-03,-1.440e-03,

     3-1.340e-03, 1.362e-02,-6.130e-03,-7.900e-04,-9.000e-04,-2.000e-04,

     4 5.080e-03,-3.610e-03, 1.700e-03, 1.830e-03, 6.800e-04, 4.000e-04,

     5-3.580e-03, 6.000e-05,-2.600e-04,-1.050e-03,-3.800e-04,-2.300e-04,

     6 2.420e-03, 9.300e-04,-1.000e-04, 4.500e-04, 1.700e-04, 1.100e-04,

     1 9.868e-01,-1.171e+00, 3.940e-01,-8.459e-02, 1.124e-02,-1.250e-03,

     2 1.001e+00,-1.383e+00, 5.044e-01,-1.152e-01, 1.658e-02,-1.830e-03,

     3 4.928e-02,-1.368e-01, 9.021e-02,-2.935e-02, 5.800e-03,-6.600e-04,

     4-3.133e-02, 5.785e-02,-2.023e-02, 2.630e-03, 1.600e-04,-8.000e-05,

     5 8.840e-03,-1.416e-02, 1.900e-03, 5.800e-04,-2.500e-04, 1.000e-05,

     6-2.300e-03, 3.080e-03, 5.500e-04,-3.700e-04, 7.000e-05, 1.000e-05/

C...Expansion coefficients for charm sea quark distribution.

      DATA (((cehlq(ix,it,nx,6,1),ix=1,6),it=1,6),nx=1,2)/

     1 9.270e-03,-1.817e-02, 9.590e-03,-6.390e-03, 1.690e-03,-1.540e-03,

     2 5.710e-03,-1.188e-02, 6.090e-03,-4.650e-03, 1.240e-03,-1.310e-03,

     3-3.960e-03, 7.100e-03,-3.590e-03, 1.840e-03,-3.900e-04, 3.400e-04,

     4 1.120e-03,-1.960e-03, 1.120e-03,-4.800e-04, 1.000e-04,-4.000e-05,

     5 4.000e-05,-3.000e-05,-1.800e-04, 9.000e-05,-5.000e-05,-2.000e-05,

     6-4.200e-04, 7.300e-04,-1.600e-04, 5.000e-05, 5.000e-05, 5.000e-05,

     1 8.098e-01,-1.042e+00, 3.398e-01,-6.824e-02, 8.760e-03,-9.000e-04,

     2 8.961e-01,-1.217e+00, 4.339e-01,-9.287e-02, 1.304e-02,-1.290e-03,

     3 3.058e-02,-1.040e-01, 7.604e-02,-2.415e-02, 4.600e-03,-5.000e-04,

     4-2.451e-02, 4.432e-02,-1.651e-02, 1.430e-03, 1.200e-04,-1.000e-04,

     5 1.122e-02,-1.457e-02, 2.680e-03, 5.800e-04,-1.200e-04, 3.000e-05,

     6-7.730e-03, 7.330e-03,-7.600e-04,-2.400e-04, 1.000e-05, 0.000e+00/

      DATA (((cehlq(ix,it,nx,6,2),ix=1,6),it=1,6),nx=1,2)/

     1 9.980e-03,-1.945e-02, 1.055e-02,-6.870e-03, 1.860e-03,-1.560e-03,

     2 5.700e-03,-1.203e-02, 6.250e-03,-4.860e-03, 1.310e-03,-1.370e-03,

     3-4.490e-03, 7.990e-03,-4.170e-03, 2.050e-03,-4.400e-04, 3.300e-04,

     4 1.470e-03,-2.480e-03, 1.460e-03,-5.700e-04, 1.200e-04,-1.000e-05,

     5-9.000e-05, 1.500e-04,-3.200e-04, 1.200e-04,-6.000e-05,-4.000e-05,

     6-4.200e-04, 7.600e-04,-1.400e-04, 4.000e-05, 7.000e-05, 5.000e-05,

     1 8.698e-01,-1.131e+00, 3.836e-01,-8.111e-02, 1.048e-02,-1.300e-03,

     2 9.626e-01,-1.321e+00, 4.854e-01,-1.091e-01, 1.583e-02,-1.700e-03,

     3 3.057e-02,-1.088e-01, 8.022e-02,-2.676e-02, 5.590e-03,-5.600e-04,

     4-2.845e-02, 5.164e-02,-1.918e-02, 2.210e-03,-4.000e-05,-1.500e-04,

     5 1.311e-02,-1.751e-02, 3.310e-03, 5.100e-04,-1.200e-04, 5.000e-05,

     6-8.590e-03, 8.380e-03,-9.200e-04,-2.600e-04, 1.000e-05,-1.000e-05/

C...Expansion coefficients for bottom sea quark distribution.

      DATA (((cehlq(ix,it,nx,7,1),ix=1,6),it=1,6),nx=1,2)/

     1 9.010e-03,-1.401e-02, 7.150e-03,-4.130e-03, 1.260e-03,-1.040e-03,

     2 6.280e-03,-9.320e-03, 4.780e-03,-2.890e-03, 9.100e-04,-8.200e-04,

     3-2.930e-03, 4.090e-03,-1.890e-03, 7.600e-04,-2.300e-04, 1.400e-04,

     4 3.900e-04,-1.200e-03, 4.400e-04,-2.500e-04, 2.000e-05,-2.000e-05,

     5 2.600e-04, 1.400e-04,-8.000e-05, 1.000e-04, 1.000e-05, 1.000e-05,

     6-2.600e-04, 3.200e-04, 1.000e-05,-1.000e-05, 1.000e-05,-1.000e-05,

     1 8.029e-01,-1.075e+00, 3.792e-01,-7.843e-02, 1.007e-02,-1.090e-03,

     2 7.903e-01,-1.099e+00, 4.153e-01,-9.301e-02, 1.317e-02,-1.410e-03,

     3-1.704e-02,-1.130e-02, 2.882e-02,-1.341e-02, 3.040e-03,-3.600e-04,

     4-7.200e-04, 7.230e-03,-5.160e-03, 1.080e-03,-5.000e-05,-4.000e-05,

     5 3.050e-03,-4.610e-03, 1.660e-03,-1.300e-04,-1.000e-05, 1.000e-05,

     6-4.360e-03, 5.230e-03,-1.610e-03, 2.000e-04,-2.000e-05, 0.000e+00/

      DATA (((cehlq(ix,it,nx,7,2),ix=1,6),it=1,6),nx=1,2)/

     1 8.980e-03,-1.459e-02, 7.510e-03,-4.410e-03, 1.310e-03,-1.070e-03,

     2 5.970e-03,-9.440e-03, 4.800e-03,-3.020e-03, 9.100e-04,-8.500e-04,

     3-3.050e-03, 4.440e-03,-2.100e-03, 8.500e-04,-2.400e-04, 1.400e-04,

     4 5.300e-04,-1.300e-03, 5.600e-04,-2.700e-04, 3.000e-05,-2.000e-05,

     5 2.000e-04, 1.400e-04,-1.100e-04, 1.000e-04, 0.000e+00, 0.000e+00,

     6-2.600e-04, 3.200e-04, 0.000e+00,-3.000e-05, 1.000e-05,-1.000e-05,

     1 8.672e-01,-1.174e+00, 4.265e-01,-9.252e-02, 1.244e-02,-1.460e-03,

     2 8.500e-01,-1.194e+00, 4.630e-01,-1.083e-01, 1.614e-02,-1.830e-03,

     3-2.241e-02,-5.630e-03, 2.815e-02,-1.425e-02, 3.520e-03,-4.300e-04,

     4-7.300e-04, 8.030e-03,-5.780e-03, 1.380e-03,-1.300e-04,-4.000e-05,

     5 3.460e-03,-5.380e-03, 1.960e-03,-2.100e-04, 1.000e-05, 1.000e-05,

     6-4.850e-03, 5.950e-03,-1.890e-03, 2.600e-04,-3.000e-05, 0.000e+00/

C...Expansion coefficients for top sea quark distribution.

      DATA (((cehlq(ix,it,nx,8,1),ix=1,6),it=1,6),nx=1,2)/

     1 4.410e-03,-7.480e-03, 3.770e-03,-2.580e-03, 7.300e-04,-7.100e-04,

     2 3.840e-03,-6.050e-03, 3.030e-03,-2.030e-03, 5.800e-04,-5.900e-04,

     3-8.800e-04, 1.660e-03,-7.500e-04, 4.700e-04,-1.000e-04, 1.000e-04,

     4-8.000e-05,-1.500e-04, 1.200e-04,-9.000e-05, 3.000e-05, 0.000e+00,

     5 1.300e-04,-2.200e-04,-2.000e-05,-2.000e-05,-2.000e-05,-2.000e-05,

     6-7.000e-05, 1.900e-04,-4.000e-05, 2.000e-05, 0.000e+00, 0.000e+00,

     1 6.623e-01,-9.248e-01, 3.519e-01,-7.930e-02, 1.110e-02,-1.180e-03,

     2 6.380e-01,-9.062e-01, 3.582e-01,-8.479e-02, 1.265e-02,-1.390e-03,

     3-2.581e-02, 2.125e-02, 4.190e-03,-4.980e-03, 1.490e-03,-2.100e-04,

     4 7.100e-04, 5.300e-04,-1.270e-03, 3.900e-04,-5.000e-05,-1.000e-05,

     5 3.850e-03,-5.060e-03, 1.860e-03,-3.500e-04, 4.000e-05, 0.000e+00,

     6-3.530e-03, 4.460e-03,-1.500e-03, 2.700e-04,-3.000e-05, 0.000e+00/

      DATA (((cehlq(ix,it,nx,8,2),ix=1,6),it=1,6),nx=1,2)/

     1 4.260e-03,-7.530e-03, 3.830e-03,-2.680e-03, 7.600e-04,-7.300e-04,

     2 3.640e-03,-6.050e-03, 3.030e-03,-2.090e-03, 5.900e-04,-6.000e-04,

     3-9.200e-04, 1.710e-03,-8.200e-04, 5.000e-04,-1.200e-04, 1.000e-04,

     4-5.000e-05,-1.600e-04, 1.300e-04,-9.000e-05, 3.000e-05, 0.000e+00,

     5 1.300e-04,-2.100e-04,-1.000e-05,-2.000e-05,-2.000e-05,-1.000e-05,

     6-8.000e-05, 1.800e-04,-5.000e-05, 2.000e-05, 0.000e+00, 0.000e+00,

     1 7.146e-01,-1.007e+00, 3.932e-01,-9.246e-02, 1.366e-02,-1.540e-03,

     2 6.856e-01,-9.828e-01, 3.977e-01,-9.795e-02, 1.540e-02,-1.790e-03,

     3-3.053e-02, 2.758e-02, 2.150e-03,-4.880e-03, 1.640e-03,-2.500e-04,

     4 9.200e-04, 4.200e-04,-1.340e-03, 4.600e-04,-8.000e-05,-1.000e-05,

     5 4.230e-03,-5.660e-03, 2.140e-03,-4.300e-04, 6.000e-05, 0.000e+00,

     6-3.890e-03, 5.000e-03,-1.740e-03, 3.300e-04,-4.000e-05, 0.000e+00/


C...The following data lines are coefficients needed in the

C...Duke, Owens proton structure function parametrizations, see below.

C...Expansion coefficients for (up+down) valence quark distribution.

      DATA ((cdo(ip,is,1,1),is=1,6),ip=1,3)/

     1 4.190e-01, 3.460e+00, 4.400e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2 4.000e-03, 7.240e-01,-4.860e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     3-7.000e-03,-6.600e-02, 1.330e+00, 0.000e+00, 0.000e+00, 0.000e+00/

      DATA ((cdo(ip,is,1,2),is=1,6),ip=1,3)/

     1 3.740e-01, 3.330e+00, 6.030e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2 1.400e-02, 7.530e-01,-6.220e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     3 0.000e+00,-7.600e-02, 1.560e+00, 0.000e+00, 0.000e+00, 0.000e+00/

C...Expansion coefficients for down valence quark distribution.

      DATA ((cdo(ip,is,2,1),is=1,6),ip=1,3)/

     1 7.630e-01, 4.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2-2.370e-01, 6.270e-01,-4.210e-01, 0.000e+00, 0.000e+00, 0.000e+00,

     3 2.600e-02,-1.900e-02, 3.300e-02, 0.000e+00, 0.000e+00, 0.000e+00/

      DATA ((cdo(ip,is,2,2),is=1,6),ip=1,3)/

     1 7.610e-01, 3.830e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2-2.320e-01, 6.270e-01,-4.180e-01, 0.000e+00, 0.000e+00, 0.000e+00,

     3 2.300e-02,-1.900e-02, 3.600e-02, 0.000e+00, 0.000e+00, 0.000e+00/

C...Expansion coefficients for (up+down+strange) sea quark distribution.

      DATA ((cdo(ip,is,3,1),is=1,6),ip=1,3)/

     1 1.265e+00, 0.000e+00, 8.050e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2-1.132e+00,-3.720e-01, 1.590e+00, 6.310e+00,-1.050e+01, 1.470e+01,

     3 2.930e-01,-2.900e-02,-1.530e-01,-2.730e-01,-3.170e+00, 9.800e+00/

      DATA ((cdo(ip,is,3,2),is=1,6),ip=1,3)/

     1 1.670e+00, 0.000e+00, 9.150e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2-1.920e+00,-2.730e-01, 5.300e-01, 1.570e+01,-1.010e+02, 2.230e+02,

     3 5.820e-01,-1.640e-01,-7.630e-01,-2.830e+00, 4.470e+01,-1.170e+02/

C...Expansion coefficients for charm sea quark distribution.

      DATA ((cdo(ip,is,4,1),is=1,6),ip=1,3)/

     1 0.000e+00,-3.600e-02, 6.350e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2 1.350e-01,-2.220e-01, 3.260e+00,-3.030e+00, 1.740e+01,-1.790e+01,

     3-7.500e-02,-5.800e-02,-9.090e-01, 1.500e+00,-1.130e+01, 1.560e+01/

       DATA ((cdo(ip,is,4,2),is=1,6),ip=1,3)/

     1 0.000e+00,-1.200e-01, 3.510e+00, 0.000e+00, 0.000e+00, 0.000e+00,

     2 6.700e-02,-2.330e-01, 3.660e+00,-4.740e-01, 9.500e+00,-1.660e+01,

     3-3.100e-02,-2.300e-02,-4.530e-01, 3.580e-01,-5.430e+00, 1.550e+01/

C...Expansion coefficients for gluon distribution.

      DATA ((cdo(ip,is,5,1),is=1,6),ip=1,3)/

     1 1.560e+00, 0.000e+00, 6.000e+00, 9.000e+00, 0.000e+00, 0.000e+00,

     2-1.710e+00,-9.490e-01, 1.440e+00,-7.190e+00,-1.650e+01, 1.530e+01,

     3 6.380e-01, 3.250e-01,-1.050e+00, 2.550e-01, 1.090e+01,-1.010e+01/

      DATA ((cdo(ip,is,5,2),is=1,6),ip=1,3)/

     1 8.790e-01, 0.000e+00, 4.000e+00, 9.000e+00, 0.000e+00, 0.000e+00,

     2-9.710e-01,-1.160e+00, 1.230e+00,-5.640e+00,-7.540e+00,-5.960e-01,

     3 4.340e-01, 4.760e-01,-2.540e-01,-8.170e-01, 5.500e+00, 1.260e-01/


C...The following data lines are coefficients needed in the

C...Owens pion structure function parametrizations, see below.

C...Expansion coefficients for up and down valence quark distributions.

      DATA ((cow(ip,is,1,1),is=1,5),ip=1,3)/

     1  4.0000e-01,  7.0000e-01,  0.0000e+00,  0.0000e+00,  0.0000e+00,

     2 -6.2120e-02,  6.4780e-01,  0.0000e+00,  0.0000e+00,  0.0000e+00,

     3 -7.1090e-03,  1.3350e-02,  0.0000e+00,  0.0000e+00,  0.0000e+00/

      DATA ((cow(ip,is,1,2),is=1,5),ip=1,3)/

     1  4.0000e-01,  6.2800e-01,  0.0000e+00,  0.0000e+00,  0.0000e+00,

     2 -5.9090e-02,  6.4360e-01,  0.0000e+00,  0.0000e+00,  0.0000e+00,

     3 -6.5240e-03,  1.4510e-02,  0.0000e+00,  0.0000e+00,  0.0000e+00/

C...Expansion coefficients for gluon distribution.

      DATA ((cow(ip,is,2,1),is=1,5),ip=1,3)/

     1  8.8800e-01,  0.0000e+00,  3.1100e+00,  6.0000e+00,  0.0000e+00,

     2 -1.8020e+00, -1.5760e+00, -1.3170e-01,  2.8010e+00, -1.7280e+01,

     3  1.8120e+00,  1.2000e+00,  5.0680e-01, -1.2160e+01,  2.0490e+01/

      DATA ((cow(ip,is,2,2),is=1,5),ip=1,3)/

     1  7.9400e-01,  0.0000e+00,  2.8900e+00,  6.0000e+00,  0.0000e+00,

     2 -9.1440e-01, -1.2370e+00,  5.9660e-01, -3.6710e+00, -8.1910e+00,

     3  5.9660e-01,  6.5820e-01, -2.5500e-01, -2.3040e+00,  7.7580e+00/

C...Expansion coefficients for (up+down+strange) quark sea distribution.

      DATA ((cow(ip,is,3,1),is=1,5),ip=1,3)/

     1  9.0000e-01,  0.0000e+00,  5.0000e+00,  0.0000e+00,  0.0000e+00,

     2 -2.4280e-01, -2.1200e-01,  8.6730e-01,  1.2660e+00,  2.3820e+00,

     3  1.3860e-01,  3.6710e-03,  4.7470e-02, -2.2150e+00,  3.4820e-01/

      DATA ((cow(ip,is,3,2),is=1,5),ip=1,3)/

     1  9.0000e-01,  0.0000e+00,  5.0000e+00,  0.0000e+00,  0.0000e+00,

     2 -1.4170e-01, -1.6970e-01, -2.4740e+00, -2.5340e+00,  5.6210e-01,

     3 -1.7400e-01, -9.6230e-02,  1.5750e+00,  1.3780e+00, -2.7010e-01/

C...Expansion coefficients for charm quark sea distribution.

      DATA ((cow(ip,is,4,1),is=1,5),ip=1,3)/

     1  0.0000e+00, -2.2120e-02,  2.8940e+00,  0.0000e+00,  0.0000e+00,

     2  7.9280e-02, -3.7850e-01,  9.4330e+00,  5.2480e+00,  8.3880e+00,

     3 -6.1340e-02, -1.0880e-01, -1.0852e+01, -7.1870e+00, -1.1610e+01/

      DATA ((cow(ip,is,4,2),is=1,5),ip=1,3)/

     1  0.0000e+00, -8.8200e-02,  1.9240e+00,  0.0000e+00,  0.0000e+00,

     2  6.2290e-02, -2.8920e-01,  2.4240e-01, -4.4630e+00, -8.3670e-01,

     3 -4.0990e-02, -1.0820e-01,  2.0360e+00,  5.2090e+00, -4.8400e-02/


C...Euler's beta function, requires ordinary Gamma function

      eulbet(x,y)=pyhigamm(x)*pyhigamm(y)/pyhigamm(x+y)


C...Reset structure functions, check x and hadron flavour.

      alam=0.

      DO 100 kfl=-6,6

  100 xpq(kfl)=0.

      IF(x.LT.0..OR.x.GT.1.) THEN

        WRITE(mstu(11),1000) x

        RETURN

      ENDIF

      kfa=iabs(kf)

      IF(kfa.NE.211.AND.kfa.NE.2212.AND.kfa.NE.2112) THEN

        WRITE(mstu(11),1100) kf

        RETURN

      ENDIF


C...Call user-supplied structure function. Select proton/neutron/pion.

      IF(mstp(51).EQ.0.OR.mstp(52).GE.2) THEN

        kfe=kfa

        IF(kfa.EQ.2112) kfe=2212

        CALL pyhistfe(kfe,x,q2,xpq)

        goto 230

      ENDIF

      IF(kfa.EQ.211) goto 200


      IF(mstp(51).EQ.1.OR.mstp(51).EQ.2) THEN

C...Proton structure functions from Eichten, Hinchliffe, Lane, Quigg.

C...Allowed variable range: 5 GeV2 < Q2 < 1E8 GeV2; 1E-4 < x < 1


C...Determine set, Lamdba and x and t expansion variables.

        nset=mstp(51)

        IF(nset.EQ.1) alam=0.2

        IF(nset.EQ.2) alam=0.29

        tmin=log(5./alam**2)

        tmax=log(1e8/alam**2)

        IF(mstp(52).EQ.0) THEN

          t=tmin

        ELSE

          t=log(q2/alam**2)

        ENDIF

        vt=max(-1.,min(1.,(2.*t-tmax-tmin)/(tmax-tmin)))

        nx=1

        IF(x.LE.0.1) nx=2

        IF(nx.EQ.1) vx=(2.*x-1.1)/0.9

        IF(nx.EQ.2) vx=max(-1.,(2.*log(x)+11.51293)/6.90776)

        cxs=1.

        IF(x.LT.1e-4.AND.abs(parp(51)-1.).GT.0.01) cxs=

     &  (1e-4/x)**(parp(51)-1.)


C...Chebyshev polynomials for x and t expansion.

        tx(1)=1.

        tx(2)=vx

        tx(3)=2.*vx**2-1.

        tx(4)=4.*vx**3-3.*vx

        tx(5)=8.*vx**4-8.*vx**2+1.

        tx(6)=16.*vx**5-20.*vx**3+5.*vx

        tt(1)=1.

        tt(2)=vt

        tt(3)=2.*vt**2-1.

        tt(4)=4.*vt**3-3.*vt

        tt(5)=8.*vt**4-8.*vt**2+1.

        tt(6)=16.*vt**5-20.*vt**3+5.*vt


C...Calculate structure functions.

        DO 120 kfl=1,6

        xqsum=0.

        DO 110 it=1,6

        DO 110 ix=1,6

  110   xqsum=xqsum+cehlq(ix,it,nx,kfl,nset)*tx(ix)*tt(it)

  120   xq(kfl)=xqsum*(1.-x)**nehlq(kfl,nset)*cxs


C...Put into output array.

        xpq(0)=xq(4)

        xpq(1)=xq(2)+xq(3)

        xpq(2)=xq(1)+xq(3)

        xpq(3)=xq(5)

        xpq(4)=xq(6)

        xpq(-1)=xq(3)

        xpq(-2)=xq(3)

        xpq(-3)=xq(5)

        xpq(-4)=xq(6)


C...Special expansion for bottom (threshold effects).

        IF(mstp(54).GE.5) THEN

          IF(nset.EQ.1) tmin=8.1905

          IF(nset.EQ.2) tmin=7.4474

          IF(t.LE.tmin) goto 140

          vt=max(-1.,min(1.,(2.*t-tmax-tmin)/(tmax-tmin)))

          tt(1)=1.

          tt(2)=vt

          tt(3)=2.*vt**2-1.

          tt(4)=4.*vt**3-3.*vt

          tt(5)=8.*vt**4-8.*vt**2+1.

          tt(6)=16.*vt**5-20.*vt**3+5.*vt

          xqsum=0.

          DO 130 it=1,6

          DO 130 ix=1,6

  130     xqsum=xqsum+cehlq(ix,it,nx,7,nset)*tx(ix)*tt(it)

          xpq(5)=xqsum*(1.-x)**nehlq(7,nset)

          xpq(-5)=xpq(5)

  140     CONTINUE

        ENDIF


C...Special expansion for top (threshold effects).

        IF(mstp(54).GE.6) THEN

          IF(nset.EQ.1) tmin=11.5528

          IF(nset.EQ.2) tmin=10.8097

          tmin=tmin+2.*log(pmas(6,1)/30.)

          tmax=tmax+2.*log(pmas(6,1)/30.)

          IF(t.LE.tmin) goto 160

          vt=max(-1.,min(1.,(2.*t-tmax-tmin)/(tmax-tmin)))

          tt(1)=1.

          tt(2)=vt

          tt(3)=2.*vt**2-1.

          tt(4)=4.*vt**3-3.*vt

          tt(5)=8.*vt**4-8.*vt**2+1.

          tt(6)=16.*vt**5-20.*vt**3+5.*vt

          xqsum=0.

          DO 150 it=1,6

          DO 150 ix=1,6

  150     xqsum=xqsum+cehlq(ix,it,nx,8,nset)*tx(ix)*tt(it)

          xpq(6)=xqsum*(1.-x)**nehlq(8,nset)

          xpq(-6)=xpq(6)

  160     CONTINUE

        ENDIF


      ELSEIF(mstp(51).EQ.3.OR.mstp(51).EQ.4) THEN

C...Proton structure functions from Duke, Owens.

C...Allowed variable range: 4 GeV2 < Q2 < approx 1E6 GeV2.


C...Determine set, Lambda and s expansion parameter.

        nset=mstp(51)-2

        IF(nset.EQ.1) alam=0.2

        IF(nset.EQ.2) alam=0.4

        IF(mstp(52).LE.0) THEN

          sd=0.

        ELSE

          sd=log(log(max(q2,4.)/alam**2)/log(4./alam**2))

        ENDIF


C...Calculate structure functions.

        DO 180 kfl=1,5

        DO 170 is=1,6

  170   ts(is)=cdo(1,is,kfl,nset)+cdo(2,is,kfl,nset)*sd+

     &  cdo(3,is,kfl,nset)*sd**2

        IF(kfl.LE.2) THEN

          xq(kfl)=x**ts(1)*(1.-x)**ts(2)*(1.+ts(3)*x)/(eulbet(ts(1),

     &    ts(2)+1.)*(1.+ts(3)*ts(1)/(ts(1)+ts(2)+1.)))

        ELSE

          xq(kfl)=ts(1)*x**ts(2)*(1.-x)**ts(3)*(1.+ts(4)*x+ts(5)*x**2+

     &    ts(6)*x**3)

        ENDIF

  180   CONTINUE


C...Put into output arrays.

        xpq(0)=xq(5)

        xpq(1)=xq(2)+xq(3)/6.

        xpq(2)=3.*xq(1)-xq(2)+xq(3)/6.

        xpq(3)=xq(3)/6.

        xpq(4)=xq(4)

        xpq(-1)=xq(3)/6.

        xpq(-2)=xq(3)/6.

        xpq(-3)=xq(3)/6.

        xpq(-4)=xq(4)


C...Proton structure functions from Diemoz, Ferroni, Longo, Martinelli.

C...These are accessed via PYHISTFE since the files needed may not always

C...available.

      ELSEIF(mstp(51).GE.11.AND.mstp(51).LE.13) THEN

        CALL pyhistfe(2212,x,q2,xpq)


C...Unknown proton parametrization.

      ELSE

        WRITE(mstu(11),1200) mstp(51)

      ENDIF

      goto 230


  200 IF((mstp(51).GE.1.AND.mstp(51).LE.4).OR.

     &(mstp(51).GE.11.AND.mstp(51).LE.13)) THEN

C...Pion structure functions from Owens.

C...Allowed variable range: 4 GeV2 < Q2 < approx 2000 GeV2.


C...Determine set, Lambda and s expansion variable.

        nset=1

        IF(mstp(51).EQ.2.OR.mstp(51).EQ.4.OR.mstp(51).EQ.13) nset=2

        IF(nset.EQ.1) alam=0.2

        IF(nset.EQ.2) alam=0.4

        IF(mstp(52).LE.0) THEN

          sd=0.

        ELSE

          sd=log(log(max(q2,4.)/alam**2)/log(4./alam**2))

        ENDIF


C...Calculate structure functions.

        DO 220 kfl=1,4

        DO 210 is=1,5

  210   ts(is)=cow(1,is,kfl,nset)+cow(2,is,kfl,nset)*sd+

     &  cow(3,is,kfl,nset)*sd**2

        IF(kfl.EQ.1) THEN

          xq(kfl)=x**ts(1)*(1.-x)**ts(2)/eulbet(ts(1),ts(2)+1.)

        ELSE

          xq(kfl)=ts(1)*x**ts(2)*(1.-x)**ts(3)*(1.+ts(4)*x+ts(5)*x**2)

        ENDIF

  220   CONTINUE


C...Put into output arrays.

        xpq(0)=xq(2)

        xpq(1)=xq(3)/6.

        xpq(2)=xq(1)+xq(3)/6.

        xpq(3)=xq(3)/6.

        xpq(4)=xq(4)

        xpq(-1)=xq(1)+xq(3)/6.

        xpq(-2)=xq(3)/6.

        xpq(-3)=xq(3)/6.

        xpq(-4)=xq(4)


C...Unknown pion parametrization.

      ELSE

        WRITE(mstu(11),1200) mstp(51)

      ENDIF


C...Isospin conjugation for neutron, charge conjugation for antipart.

  230 IF(kfa.EQ.2112) THEN

        xps=xpq(1)

        xpq(1)=xpq(2)

        xpq(2)=xps

        xps=xpq(-1)

        xpq(-1)=xpq(-2)

        xpq(-2)=xps

      ENDIF

      IF(kf.LT.0) THEN

        DO 240 kfl=1,4

        xps=xpq(kfl)

        xpq(kfl)=xpq(-kfl)

  240   xpq(-kfl)=xps

      ENDIF


C...Check positivity and reset above maximum allowed flavour.

      DO 250 kfl=-6,6

      xpq(kfl)=max(0.,xpq(kfl))

  250 IF(iabs(kfl).GT.mstp(54)) xpq(kfl)=0.


C...consider nuclear effect on the structure function

        IF((jbt.NE.1.AND.jbt.NE.2).OR.ihpr2(6).EQ.0

     &                  .OR.ihnt2(16).EQ.1) go to 400

        atnm=ihnt2(2*jbt-1)

        IF(atnm.LE.1.0) go to 400

        IF(jbt.EQ.1) THEN

          bbr2=(yp(1,ihnt2(11))**2+yp(2,ihnt2(11))**2)/1.44/atnm**0.66666

        ELSEIF(jbt.EQ.2) THEN

          bbr2=(yt(1,ihnt2(12))**2+yt(2,ihnt2(12))**2)/1.44/atnm**0.66666

        ENDIF

        bbr2=min(1.0,bbr2)

        abx=(atnm**0.33333333-1.0)

        apx=hipr1(6)*4.0/3.0*abx*sqrt(1.0-bbr2)

        aax=1.192*alog(atnm)**0.1666666

        rrx=aax*(x**3-1.2*x**2+0.21*x)+1.0

     &          -(apx-1.079*abx*sqrt(x)/alog(atnm+1.0))*exp(-x**2.0/0.01)

        DO 300 kfl=-6,6

                xpq(kfl)=xpq(kfl)*rrx

300     CONTINUE

C                       ********consider the nuclear effect on the structure

C                               fucntion which also depends on the impact

C                               parameter of the nuclear reaction


400     CONTINUE

C...Formats for error printouts.

 1000 FORMAT(' Error: x value outside physical range, x =',1p,e12.3)

 1100 FORMAT(' Error: illegal particle code for structure function,',

     &' KF =',i5)

 1200 FORMAT(' Error: bad value of parameter MSTP(51) in PYHISTFU,',

     &' MSTP(51) =',i5)


      RETURN

      END