      SUBROUTINE READ_R8SEG(n,xdum,fid,irec)
C     ================================================
C     This subroutine reads a real*8 1d record of length n
C     ================================================
      implicit none
      integer n,fid
      real*8 xdum(n)
      integer irec

      read(fid,rec=irec) xdum

      END

C================================================
C This subroutine reads a real*4 1d record of length n
C================================================

      SUBROUTINE READ_R4SEG(n,xdum,fid,irec)

      implicit none
      integer n,fid
      real*4 xdum(n)
      integer irec

      read(fid,rec=irec) xdum

      END
      
C================================================
C This subroutine writes a real*8 1d record of length n
C================================================

      SUBROUTINE WRITE_R8SEG(n,xdum,fid,irec)

      implicit none
      integer n,fid
      real*8 xdum(n)
      integer irec

      write(fid,rec=irec) xdum

      END

C================================================
C This subroutine writes a real*4 1d record of length n
C================================================

      SUBROUTINE WRITE_R4SEG(n,xdum,fid,irec)

      implicit none
      integer n,fid
      real*4 xdum(n)
      integer irec

      write(fid,rec=irec) xdum

      END

C==================================================
C This subroutine reads a real*4 array from file
C==================================================

      SUBROUTINE READ_R4_FIELD(nx,ny,iSlice,fld,filename)

      implicit none
      integer nx,ny,iSlice
      real fld(nx,ny)
      character*(*) filename

      integer i,j,irec,ioUnit

      call findunit(ioUnit)

      open(unit=ioUnit,file=filename,status='old',
     &     access='direct',recl=nx*4)
     
      do j=1,ny
        irec=j+(iSlice-1)*ny
        call read_r4seg(nx,fld(1,j),ioUnit,irec)
      enddo

      close(ioUnit)
      
      END

C==================================================
C This subroutine reads a real*8 array from file
C==================================================

      SUBROUTINE READ_R8_FIELD(nx,ny,iSlice,fld,filename)

      implicit none
      integer nx,ny,iSlice
      real fld(nx,ny)
      character*(*) filename

      integer i,j,irec,ioUnit

      call findunit(ioUnit)

      open(unit=ioUnit,file=filename,status='old',
     &     access='direct',recl=nx*8)
     
      do j=1,ny
        irec=j+(iSlice-1)*ny
        call read_r8seg(nx,fld(1,j),ioUnit,irec)
      enddo

      close(ioUnit)
      
      END
            
C==================================================
C This subroutine writes out a real*4 array to file
C==================================================

      SUBROUTINE WRITE_R4_FIELD(nx,ny,iSlice,fld,filename)

      implicit none
      integer nx,ny,iSlice
      real fld(nx,ny)
      character*(*) filename

      integer i,j,irec,ioUnit

      call findunit(ioUnit)

      if (iSlice.eq.1) then  !  open file here for first time
        open(unit=ioUnit,file=filename,status='unknown',
     &       access='direct',recl=nx*4)
      else  !  open existing file for appending
        open(unit=ioUnit,file=filename,status='old',
     &       access='direct',recl=nx*4)
      endif

      do j=1,ny
        irec=j+(iSlice-1)*ny
        call write_r4seg(nx,fld(1,j),ioUnit,irec)
      enddo

      close(ioUnit)
      
      END

C==================================================
C This subroutine writes out a real*8 array to file
C==================================================

      SUBROUTINE WRITE_R8_FIELD(nx,ny,iSlice,fld,filename)

      implicit none
      integer nx,ny,iSlice
      real*8 fld(nx,ny)
      character*(*) filename

      integer i,j,irec,ioUnit

      call findunit(ioUnit)

      if (iSlice.eq.1) then  !  open file here for first time
        open(unit=ioUnit,file=filename,status='unknown',
     &       access='direct',recl=nx*8)
      else  !  open existing file for appending
        open(unit=ioUnit,file=filename,status='old',
     &       access='direct',recl=nx*8)
      endif

      do j=1,ny
        irec=j+(iSlice-1)*ny
        call write_r8seg(nx,fld(1,j),ioUnit,irec)
      enddo

      close(ioUnit)
      
      END
      
C=============================================================
C=============================================================
      SUBROUTINE FINDUNIT( iounit )
C OUT:
C     iounit   integer - unit number
C
C MDSFINDUNIT returns a valid, unused unit number for f77 I/O
C The routine stops the program is an error occurs in the process
C of searching the I/O channels.
C
C Created: 03/20/99 adcroft@mit.edu

      implicit none

C Arguments
      integer iounit
C Local
      integer ii
      logical op
      integer ios
C     ------------------------------------------------------------------

C Sweep through a valid range of unit numbers
      iounit=-1
      do ii=9,99
        if (iounit.eq.-1) then
          inquire(unit=ii,iostat=ios,opened=op)
          if (ios.ne.0) then
            stop 'ABNORMAL END: S/R MDSFINDUNIT'
          endif
          if (.NOT. op) then
            iounit=ii
          endif
        endif
      enddo

C Was there an available unit number
      if (iounit.eq.-1) then
        stop ' MDSFINDUNIT: could not find an available unit number!'
      endif

C     ------------------------------------------------------------------
      return
      END

