LHAPDF/lhapdf5.3.1/wraph1.f

   1       subroutine H1evolve(xin,qin,pdf)
   2       implicit real*8 (a-h,o-z)
   3 *******************************************************
   4 c  done on  13/07/04 at  10.04.39
   5 c evolution has been made starting at q2_input =      4.000
   6 c  available for :       1.500 <= q2 <= 1000000.000
   7 c            and :  0.000057 <= x <=  0.906052
   8 c*
   9 c   for x outside limits, the closest limit
  10 c        is assumed : f2(x>xmax,q2)=f2(xmax,q2)
  11 c                     f2(x<xmin,q2)=f2(xmin,q2)
  12 c  for q2 outside limits, the closest limit
  13 c        is assumed : f2(x,q2>q2max)=f2(x,q2max)
  14 c                     f2(x,q2<q2min)=f2(x,q2min)
  15 c*
  16 c  comments, etc... to C. Pascaud or F. Zomer
  17 ********************************************************
  18       include 'parmsetup.inc'
  19       PARAMETER(n_bin_q2=86)
  20       PARAMETER(n_bin_x=100)
  21       REAL*4 xl_bin(n_bin_x),q2l_bin(n_bin_q2)
  22       PARAMETER(ngrid=20)
  23       REAL*4 f(0:ngrid,8,n_bin_x,n_bin_q2),val(8)
  24       real*8 pdf(-6:6)
  25       real*4 q2in,x,y
  26       character*16 name(nmxset)
  27       integer nmem(nmxset),ndef(nmxset),mmem
  28       common/NAME/name,nmem,ndef,mmem
  29       integer nset,iset
  30       save
  31 c
  32 c      enddo
  33       call getnset(iset)
  34       call getnmem(iset,imem)
  35 c
  36       q2in = qin*qin
  37       x=log(xin)
  38       y=log(q2in)
  39       DO i=2,n_bin_x
  40         IF(x.LT.xl_bin(i))  goto 1
  41         IF(xl_bin(i).ge.0.)  goto 1
  42       ENDDO
  43       i=n_bin_x
  44     1 i=i-1
  45       DO j=2,n_bin_q2
  46         IF(y.LT.q2l_bin(j))  GOTO 2
  47       ENDDO
  48       j=n_bin_q2
  49     2 j=j-1
  50       dx=xl_bin(i+1)-xl_bin(i)
  51       xd=(x-xl_bin(i))/dx
  52       dy=q2l_bin(j+1)-q2l_bin(j)
  53       yd=(y-q2l_bin(j))/dy
  54 c
  55       do k=1,8
  56       val(k)=f(imem,k,i,j)+xd*(f(imem,k,i+1,j)-f(imem,k,i,j))
  57      &+yd*(f(imem,k,i,j+1)-f(imem,k,i,j))
  58      &+xd*yd*(f(imem,k,i+1,j+1)+f(imem,k,i,j)
  59      &-f(imem,k,i+1,j)-f(imem,k,i,j+1))
  60       enddo
  61       pdf(-6) = 0.0d0
  62        pdf(6) = 0.0d0
  63       pdf(-5) = val(7)
  64        pdf(5) = val(7)
  65       pdf(-4) = val(6)
  66        pdf(4) = val(6)
  67       pdf(-3) = val(5)
  68        pdf(3) = val(5)
  69       pdf(-2) = val(4)
  70        pdf(2) = val(3)+val(4)
  71       pdf(-1) = val(2)
  72        pdf(1) = val(1)+val(2)
  73        pdf(0) = val(8)
  74       return
  75 c
  76       entry H1read(nset)
  77 c
  78       read(1,*)nmem(nset),ndef(nset)
  79       read(1,1000)xl_bin
  80       read(1,1000)q2l_bin
  81         do i=1,n_bin_q2
  82           q2l_bin(i)=log(q2l_bin(i))
  83         enddo
  84       do nm = 0,nmem(nset)
  85         do jval = 1,8
  86       read(1,1000)((f(nm,jval,nx,nq2),nx=1,n_bin_x),nq2=1,n_bin_q2)
  87         enddo
  88       enddo
  89       return
  90 c
  91       entry H1alfa(alfas,qalfa)
  92       call alphah1(alfas,Qalfa)
  93       return
  94 c
  95       entry H1init(Eorder,Q2fit)
  96       return
  97 c
  98       entry H1pdf(mem)
  99       call getnset(iset)
 100       call setnmem(iset,mem)
 101 c      imem = mem
 102       return
 103 c
 104  1000 format(5e13.5)
 105       END
 106 *
 107       subroutine alphah1(alpha,Qin)
 108       implicit real*8 (a-h,o-z)
 109         call getnset(nset)
 110         call GetOrderAsM(nset,iord)
 111         if(iord.eq.1) then
 112           call alphah1nlo(alpha,Qin)
 113         elseif(iord.eq.0) then
 114           call alphah1lo(alpha,Qin)
 115         else
 116           print *,'iord = ',iord
 117           stop
 118         endif
 119       return
 120       end
 121 *
 122       subroutine alphah1nlo(alpha,Qin)
 123       implicit real*8 (a-h,o-z)
 124 ****************************************************
 125 c  done on  13/07/04 at  09.10.39
 126 c evolution has been made starting at q2_input =      4.000
 127 c  available for :       1.500 <= q2 <= 1000000.000
 128 c  for q2 outside limits, the closest limit
 129 c        is assumed : f2(x,q2>q2max)=f2(x,q2max)
 130 c                     f2(x,q2<q2min)=f2(x,q2min)
 131 c*
 132 c  comments, etc... to C. Pascaud or F. Zomer
 133 ***************************************************
 134       PARAMETER(n_bin_q2=102)
 135       dimension q2l_bin(n_bin_q2)
 136       dimension f(n_bin_q2)
 137       data q2l_bin/
 138      +1.500000E+00,1.600000E+00,1.700000E+00,1.800000E+00,1.900000E+00,
 139      +1.959902E+00,1.960000E+00,1.960098E+00,2.000000E+00,2.100000E+00,
 140      +2.200000E+00,2.300000E+00,2.400000E+00,2.500000E+00,3.000000E+00,
 141      +3.500000E+00,4.000000E+00,4.500000E+00,5.000000E+00,6.000000E+00,
 142      +7.000000E+00,8.000000E+00,9.000000E+00,1.000000E+01,1.500000E+01,
 143      +2.000000E+01,2.024899E+01,2.025000E+01,2.025101E+01,2.500000E+01,
 144      +3.000000E+01,3.500000E+01,4.000000E+01,4.500000E+01,5.000000E+01,
 145      +5.500000E+01,6.000000E+01,6.500000E+01,7.000000E+01,7.500000E+01,
 146      +8.000000E+01,8.500000E+01,9.000000E+01,9.500000E+01,1.000000E+02,
 147      +1.500000E+02,2.000000E+02,2.500000E+02,3.000000E+02,3.500000E+02,
 148      +4.000000E+02,4.500000E+02,5.000000E+02,5.500000E+02,6.000000E+02,
 149      +6.500000E+02,7.000000E+02,7.500000E+02,8.000000E+02,8.500000E+02,
 150      +9.000000E+02,9.500000E+02,1.000000E+03,1.500000E+03,2.000000E+03,
 151      +2.500000E+03,3.000000E+03,3.500000E+03,4.000000E+03,4.500000E+03,
 152      +5.000000E+03,5.500000E+03,6.000000E+03,6.500000E+03,7.000000E+03,
 153      +7.500000E+03,8.000000E+03,8.500000E+03,9.000000E+03,9.500000E+03,
 154      +1.000000E+04,1.500000E+04,2.000000E+04,2.500000E+04,3.000000E+04,
 155      +3.500000E+04,4.000000E+04,4.500000E+04,5.000000E+04,5.500000E+04,
 156      +6.000000E+04,6.500000E+04,7.000000E+04,7.500000E+04,8.000000E+04,
 157      +8.500000E+04,9.000000E+04,9.500000E+04,1.000000E+05,1.500000E+05,
 158      +2.000000E+05,1.000000E+06/
 159       data f/
 160      +3.935326E-01,3.849873E-01,3.773198E-01,3.703884E-01,3.640814E-01,
 161      +3.605647E-01,3.605591E-01,3.605540E-01,3.585220E-01,3.537029E-01,
 162      +3.492356E-01,3.450786E-01,3.411967E-01,3.375601E-01,3.222805E-01,
 163      +3.104663E-01,3.009529E-01,2.930611E-01,2.863649E-01,2.755124E-01,
 164      +2.669931E-01,2.600508E-01,2.542359E-01,2.492618E-01,2.318890E-01,
 165      +2.210261E-01,2.205828E-01,2.205810E-01,2.205794E-01,2.139838E-01,
 166      +2.085977E-01,2.042587E-01,2.006485E-01,1.975722E-01,1.949020E-01,
 167      +1.925502E-01,1.904538E-01,1.885667E-01,1.868536E-01,1.852875E-01,
 168      +1.838468E-01,1.825145E-01,1.812765E-01,1.801214E-01,1.790394E-01,
 169      +1.709373E-01,1.656324E-01,1.617456E-01,1.587067E-01,1.562275E-01,
 170      +1.541435E-01,1.523521E-01,1.507856E-01,1.493968E-01,1.481517E-01,
 171      +1.470250E-01,1.459974E-01,1.450540E-01,1.441827E-01,1.433740E-01,
 172      +1.426200E-01,1.419142E-01,1.412513E-01,1.362263E-01,1.328778E-01,
 173      +1.303943E-01,1.284347E-01,1.268244E-01,1.254625E-01,1.242858E-01,
 174      +1.232522E-01,1.223323E-01,1.215046E-01,1.207533E-01,1.200661E-01,
 175      +1.194337E-01,1.188480E-01,1.183032E-01,1.177942E-01,1.173169E-01,
 176      +1.168677E-01,1.134367E-01,1.111245E-01,1.093963E-01,1.080244E-01,
 177      +1.068916E-01,1.059298E-01,1.050959E-01,1.043613E-01,1.037057E-01,
 178      +1.031144E-01,1.025766E-01,1.020837E-01,1.016292E-01,1.012077E-01,
 179      +1.008150E-01,1.004476E-01,1.001026E-01,9.977747E-02,9.728142E-02,
 180      +9.558620E-02,8.711093E-02/
 181       data init/0/
 182       if(init.eq.0) then
 183         do i=1,n_bin_q2
 184          q2l_bin(i) = log(q2l_bin(i))
 185         enddo
 186         init=1
 187       endif
 188 *
 189       q2in = qin*qin
 190 *
 191       y = log(q2in)
 192       do j = 2, n_bin_q2
 193         if (y.lt.q2l_bin(j))  goto 2
 194       enddo
 195       j=n_bin_q2
 196  2    j=j-1
 197 *
 198       dy = q2l_bin(j+1) - q2l_bin(j)
 199       yd = (y - q2l_bin(j)) / dy
 200       alpha = f(j) + yd*(f(j+1)-f(j))
 201 *
 202       return
 203       END
 204
 205 *
 206       subroutine alphah1lo(alpha,Qin)
 207       implicit real*8 (a-h,o-z)
 208 c  done on  13/07/04 at  11.32.26
 209 c evolution has been made starting at q2_input =      4.000
 210 c  available for :       1.500 <= q2 <= 1000000.000
 211 c  for q2 outside limits, the closest limit
 212 c        is assumed : f2(x,q2>q2max)=f2(x,q2max)
 213 c                     f2(x,q2<q2min)=f2(x,q2min)
 214 c*
 215 c  comments, etc... to C. Pascaud or F. Zomer
 216 ***************************************************
 217       PARAMETER(n_bin_q2=102)
 218       dimension q2l_bin(n_bin_q2)
 219       dimension f(n_bin_q2)
 220       data q2l_bin/
 221      +1.500000E+00,1.600000E+00,1.700000E+00,1.800000E+00,1.900000E+00,
 222      +1.959902E+00,1.960000E+00,1.960098E+00,2.000000E+00,2.100000E+00,
 223      +2.200000E+00,2.300000E+00,2.400000E+00,2.500000E+00,3.000000E+00,
 224      +3.500000E+00,4.000000E+00,4.500000E+00,5.000000E+00,6.000000E+00,
 225      +7.000000E+00,8.000000E+00,9.000000E+00,1.000000E+01,1.500000E+01,
 226      +2.000000E+01,2.024899E+01,2.025000E+01,2.025101E+01,2.500000E+01,
 227      +3.000000E+01,3.500000E+01,4.000000E+01,4.500000E+01,5.000000E+01,
 228      +5.500000E+01,6.000000E+01,6.500000E+01,7.000000E+01,7.500000E+01,
 229      +8.000000E+01,8.500000E+01,9.000000E+01,9.500000E+01,1.000000E+02,
 230      +1.500000E+02,2.000000E+02,2.500000E+02,3.000000E+02,3.500000E+02,
 231      +4.000000E+02,4.500000E+02,5.000000E+02,5.500000E+02,6.000000E+02,
 232      +6.500000E+02,7.000000E+02,7.500000E+02,8.000000E+02,8.500000E+02,
 233      +9.000000E+02,9.500000E+02,1.000000E+03,1.500000E+03,2.000000E+03,
 234      +2.500000E+03,3.000000E+03,3.500000E+03,4.000000E+03,4.500000E+03,
 235      +5.000000E+03,5.500000E+03,6.000000E+03,6.500000E+03,7.000000E+03,
 236      +7.500000E+03,8.000000E+03,8.500000E+03,9.000000E+03,9.500000E+03,
 237      +1.000000E+04,1.500000E+04,2.000000E+04,2.500000E+04,3.000000E+04,
 238      +3.500000E+04,4.000000E+04,4.500000E+04,5.000000E+04,5.500000E+04,
 239      +6.000000E+04,6.500000E+04,7.000000E+04,7.500000E+04,8.000000E+04,
 240      +8.500000E+04,9.000000E+04,9.500000E+04,1.000000E+05,1.500000E+05,
 241      +2.000000E+05,1.000000E+06/
 242       data f/
 243      +4.395646E-01,4.308115E-01,4.229009E-01,4.157042E-01,4.091185E-01,
 244      +4.054310E-01,4.054251E-01,4.054197E-01,4.032349E-01,3.980418E-01,
 245      +3.932134E-01,3.887078E-01,3.844897E-01,3.805290E-01,3.637916E-01,
 246      +3.507479E-01,3.401822E-01,3.313772E-01,3.238784E-01,3.116737E-01,
 247      +3.020502E-01,2.941818E-01,2.875740E-01,2.819097E-01,2.620464E-01,
 248      +2.495700E-01,2.490599E-01,2.490579E-01,2.490560E-01,2.413308E-01,
 249      +2.350218E-01,2.299395E-01,2.257114E-01,2.221089E-01,2.189825E-01,
 250      +2.162291E-01,2.137753E-01,2.115667E-01,2.095621E-01,2.077297E-01,
 251      +2.060444E-01,2.044861E-01,2.030382E-01,2.016875E-01,2.004225E-01,
 252      +1.909551E-01,1.847628E-01,1.802294E-01,1.766873E-01,1.737993E-01,
 253      +1.713728E-01,1.692881E-01,1.674658E-01,1.658508E-01,1.644033E-01,
 254      +1.630939E-01,1.619001E-01,1.608043E-01,1.597925E-01,1.588537E-01,
 255      +1.579785E-01,1.571596E-01,1.563904E-01,1.505657E-01,1.466892E-01,
 256      +1.438172E-01,1.415527E-01,1.396930E-01,1.381211E-01,1.367637E-01,
 257      +1.355719E-01,1.345115E-01,1.335578E-01,1.326924E-01,1.319011E-01,
 258      +1.311728E-01,1.304988E-01,1.298719E-01,1.292864E-01,1.287374E-01,
 259      +1.282208E-01,1.242789E-01,1.216259E-01,1.196449E-01,1.180735E-01,
 260      +1.167767E-01,1.156763E-01,1.147226E-01,1.138828E-01,1.131336E-01,
 261      +1.124583E-01,1.118440E-01,1.112813E-01,1.107625E-01,1.102815E-01,
 262      +1.098335E-01,1.094145E-01,1.090210E-01,1.086503E-01,1.058065E-01,
 263      +1.038775E-01,9.426288E-02/
 264       data init/0/
 265       if(init.eq.0) then
 266         do i=1,n_bin_q2
 267          q2l_bin(i) = log(q2l_bin(i))
 268         enddo
 269         init=1
 270       endif
 271
 272       q2in = qin*qin
 273 *
 274       y = log(q2in)
 275       do j = 2, n_bin_q2
 276         if (y.lt.q2l_bin(j))  goto 2
 277       enddo
 278       j=n_bin_q2
 279  2    j=j-1
 280 *
 281       dy = q2l_bin(j+1) - q2l_bin(j)
 282       yd = (y - q2l_bin(j)) / dy
 283       alpha = f(j) + yd*(f(j+1)-f(j))
 284 *
 285       return
 286       END
 287 *
 288