c     mediplug.for - dos version
c     J.Phillips, USGS, Denver, Colorado
c     changes from vms version are indicated by cv
c
c     fills holes using the median of surrounding values
c
      dimension a(1000), b(1000, 3), c(1000), id(14), pgm(2)
      character ifile*50, ifil1*50, ifil2*50, ifilt*50
      data ifil1/'mediplug.tmp'/
      data ddval / 1.0e+36 /
      ndold=9999
    1 write(*, 200) 
  200 format(14h enter ifile: ,$)
      read(5, 100) ifile
  100 format(a50)
cv    open(10,file=ifile,status='old',form='unformatted')
      open(10, file=ifile, status='old', form='unformatted', err=1) 
      write(*, 201) 
  201 format(14h enter ofile: ,$)
      read(5, 100) ifile
      ifil2=ifile
cv    open(11,file=ifile,status='new',form='unformatted')
      open(11, file=ifil1, status='unknown', form='unformatted')
9     read(10) id, pgm, nc, nr, nz, xo, dx, yo, dy
      write(11) id, pgm, nc, nr, nz, xo, dx, yo, dy
      read(10) dlt, (b(i,3),i = 1, nc)
      read(10) dlt, (b(i,2),i = 1, nc)
      nd = 0
      if(ndold.eq.9999) ndold=nc*nr
      do 20 i=1,nc
      if(b(i,3).gt.ddval) b(i,3)=b(i,2)
      if(b(i,3).gt.ddval) nd = nd + 1
   20 continue
      write(11) dlt, (b(i,3),i = 1, nc)
      do 14 ii = 3, nr
        read(10) dlt, (b(i,1),i = 1, nc)
        c(1) = b(1,2)
        if(c(1).gt.ddval) c(1)=b(2,2)
        c(nc) = b(nc,2)
        if(c(nc).gt.ddval) c(nc)=b(nc-1,2)
        do 12 i = 2, nc - 1
          c(i) = b(i,2)
          if (c(i) .gt. ddval) then
            j = 0
            do 11 k = 1, 3
              do 10 l = i - 1, i + 1
                if (b(l,k) .lt. ddval) then
                  j = j + 1
                  a(j) = b(l,k)
                end if
10            continue
11          continue
            if (j .gt. 1) call median(j, a)
            if (j .ge. 1) c(i) = a(1)
          end if
12      continue
        do 22 i=1,nc
22      if(c(i).gt.ddval) nd = nd + 1
        write(11) dlt, (c(i),i = 1, nc)
        do 13 i = 1, nc
          b(i,3) = b(i,2)
c          b(i,3) = c(i)
          b(i,2) = b(i,1)
13      continue
14    continue
      do 21 i=1,nc
      if(b(i,1).gt.ddval) b(i,1)=b(i,3)
      if(b(i,1).gt.ddval) nd = nd + 1
   21 continue
      write(11) dlt, (b(i,1),i = 1, nc)
      write(*, 101) nd,ndold
  101 format(1h+,i10,13h dvals remain)
      if(nd.lt.ndold) then
        ndold=nd
        close(10)
        close(11)
        ifilt=ifil1
        ifil1=ifil2
        ifil2=ifilt
        open(10, file=ifil2, form='unformatted', status='old')
        open(11, file=ifil1, form='unformatted', status='unknown')
        go to 9
      endif
      close(10)
      if(ifil1.eq.ifile) go to 7
      rewind(11)
      open(10,file=ifile,form='unformatted',status='unknown')
      read(11) id, pgm, nc, nr, nz, xo, dx, yo, dy
      write(10) id, pgm, nc, nr, nz, xo, dx, yo, dy
      do 8 j=1,nr
      read(11) dlt, (b(i,3),i = 1, nc)
      write(10) dlt, (b(i,3),i = 1, nc)
    8 continue
    7 stop
      end
c
      subroutine median(n, a)
      dimension a(n)
   10 if (n .lt. 2) return
      if(n.eq.2) then
       a(1)=(a(1)+a(2))/2.
       return
      endif
      ax = a(1)
      an = a(1)
      iax = 1
      ian = 1
      do 11 i = 2, n
        if (a(i) .gt. ax) then
          ax = a(i)
          iax = i
        end if
        if (a(i) .lt. an) then
          an = a(i)
          ian = i
        end if
11    continue
      if (an .eq. ax) then
        a(1) = an
        return 
      end if
      j = 0
      do 20 i = 1, n
      temp = a(i)
      if ((i .eq. iax) .or. (i .eq. ian)) goto 20
      j = j + 1
      a(j) = temp
   20 continue
      n = j
      goto 10
      end
