c
c  return magnetic field Bz in kG for z and x in cm
c
      subroutine bend( z, x, hz)
      implicit none
      integer i, j, iz, ix, init
      real z, x, hz, z0, x0, bz0, bx0
      real dz, dx, delz, delx
      real bz(26, 301)
      real bzdz,bzdx
c
c
      data init/0/, dz/0.5/, dx/1./

      data bz/7826*0.0/
c
      if(init.eq.0) then
         init = 1

         open(12,FILE='bend.dat',form='FORMATTED')

         do i=1, 301
            read(12, *) z0, (bz(j,i),j=1,26)
         end do
      end if

c------------------------------  check z and r
      if(abs(z).ge.150.0.or.abs(x).ge.25.0) then
         hz = 0.0
         return
      end if

         iz = (abs(z)+0.001)/dz
         ix = (abs(x)+0.001)/dx

         iz = iz + 1
         ix = ix + 1

c         write(6,777) z, x, iz,ix
c 777     format(2x,'z=',2f8.2,2i4)
c

      bz0 = bz(ix,iz)

c      write(6,750) bz0
c 750  format(2x,'b0=',f10.2)

      delx = 0.
      if(abs(x).gt.0.0) delx = abs(x) - float(ix-1) * dx

      bzdz = (bz(ix,iz+1)-bz0)/dz
      bzdx = (bz(ix+1,iz)-bz0)/dx

c         write(6,752) bz(ix,iz+1),bz(ix+1,iz)
c 752     format(2x,'br=',2f10.2)

      delz = abs(z) - float(iz-1) * dz

      hz = bz0+bzdz*delz+bzdx*delx

      hz = hz * 0.001

c      write(6,800) delz,delx
c 800  format(2x,'delz,delx=',2f10.2)
c      write(6,880) bzdz,bzdx
c 880  format(2x,'bzdz,bzdx ',2f10.2)
c      write(6,888) hz
c 888  format(2x,'hz=',f10.2)
     
      return
      end