subroutine intrp(wltab,absop,wlgrid,abgrd,nfr,nfgrid)
c     =====================================================
c
c     a more efficient interpolation routine - using bisection
c
      INCLUDE 'PARAMS.FOR'
      dimension wltab(1),absop(1),wlgrid(1),abgrd(1)
      dimension yint(mfgrid),jint(mfgrid)
c
c      set up interpolation coefficients for an interpolation
c      by bisection
c
         fr1=wltab(1)
         fr2=wltab(nfr)
         do ij=1,nfgrid
            xint=wlgrid(ij)
            jl=0
            ju=nfr+1
   10       continue
            if(ju-jl.gt.1) then
               jm=(ju+jl)/2
               if((fr2.gt.fr1).eqv.(xint.gt.wltab(jm))) then
                  jl=jm
                else
                  ju=jm
               end if
               go to 10
            end if
            j=jl
            if(j.eq.nfr) j=j-1
            if(j.eq.0) j=j+1
            jint(ij)=j
c           yint(ij)=un/log10(wltab(j+1)/wltab(j))
            yint(ij)=1./(wltab(j+1)-wltab(j))
         end do
c
         do ij=1,nfgrid
            j=jint(ij)
            rc=(absop(j+1)-absop(j))*yint(ij)
c           abgrd(ij)=rc*log10(wlgrid(ij)/wltab(j))+absop(j)
            abgrd(ij)=rc*(wlgrid(ij)-wltab(j))+absop(j)
         end do
      return
      end