subroutine wq(val,eval)

      implicit none

      integer i,n
      parameter(n=37)
      double precision q(n), w(n)
      double precision val,eval,slop,y_int

      data q/ 0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0,
     &     1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0,
     &     2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3.0,
     &     3.1, 3.2, 3.3, 3.4, 3.5, 3.6/             ! (meV)
      data w/  0.000000, 0.235735, 0.484567, 0.739947, 0.969134,
     &     1.198321, 1.407863, 1.538827, 1.656695, 1.722177, 1.774562,
     &     1.807303, 1.800755, 1.768014, 1.695984, 1.597761, 1.466797,
     &     1.342381, 1.211417, 1.139387, 1.172128, 1.309640, 1.525730,
     &     1.774562, 2.029942, 2.154358, 2.226388, 2.265677, 2.291870,
     &     2.318063, 2.337707, 2.357352, 2.376997, 2.396641, 2.403189,
     &     2.409738, 2.416286/ ! no unit (normalization factor)

      do i=1,n-1
         if(val.lt.q(1).or.val.gt.q(n)) then
            write(*,*) i, val,q(1),q(n),
     +           ' The value is out of range ... wq code'
            if(val.lt.q(1)) eval=w(1)
            if(val.gt.q(n)) eval = w(n)
         elseif(val.eq.q(i)) then
            eval = w(i)
         elseif(val.eq.q(i+1)) then
            eval = w(i+1)
         elseif(val.gt.q(i).and.val.lt.q(i+1)) then
            slop = (w(i+1)-w(i))/(q(i+1)-q(i))
            y_int = w(i) - slop*q(i)
            eval = slop*val+y_int
         endif
      enddo

      return
      end