	subroutine expand_iy0(a,u,nx,nxmax,ny,xfill,ixml,yinc,y)
c
c	puts a small non-equally spaced (in y) array into a bigger one 
c	and inverts in y
c
	dimension a(1),u(nxmax,ny),y(ny)
	nxnew=(nx-1)*ixml+1
	nynew=(y(ny)-y(1))/yinc+1.1
	print*,nx,ny
	print*,y(1),y(ny),yinc
	print*,nxnew,nynew
	call init_array(a,nxnew,nynew,1,xfill)
	do 100 iy=1,ny
	do 100 ix=1,nx
	ly=(y(iy)-y(1))/yinc+1.1
	ky=nynew-ly
	k=(ky)*nxnew+1+(ix-1)*ixml
100	a(k)=u(ix,iy)
	nx=nxnew
	ny=nynew
	return
	end
	subroutine expand_iy(a,u,nx,nxmax,ny,xfill,ixml,iyml)
c
c	puts a small equally spaced array into a bigger one and inverts in y
c
	dimension a(1),u(nxmax,ny)
	nxnew=(nx-1)*ixml+1
	nynew=(ny-1)*iyml+1
	call init_array(a,nxnew,nynew,1,xfill)
	do 100 iy=1,ny
	do 100 ix=1,nx
	ky=ny-iy
	k=(ky)*nxnew*iyml+1+(ix-1)*ixml
100	a(k)=u(ix,iy)
	nx=nxnew
	ny=nynew
	return
	end
	subroutine expand(a,u,nx,nxmax,ny,xfill,ixml,iyml)
c
c	puts a small array into a big one
c
	dimension a(1),u(nxmax,ny)
	nxnew=(nx-1)*ixml+1
	nynew=(ny-1)*iyml+1
	call init_array(a,nxnew,nynew,1,xfill)
	do 100 iy=1,ny
	do 100 ix=1,nx
	k=(iy-1)*nxnew*iyml+1+(ix-1)*ixml
100	a(k)=u(ix,iy)
	nx=nxnew
	ny=nynew
	return
	end
      SUBROUTINE INTY(A,NX,NY,XFILL)
      DIMENSION A(NX,NY)
      DO 100 IX=1,NX
	NGOOD=0
	do 5 iy=1,ny
5	if(a(ix,iy).ne.xfill) ngood=ngood+1
	IF(NGOOD.LT.2.or.NGOOD.GE.NY) GO TO 100
C	FIND FIRST GOOD POINT
	L1=0
10	L1=L1+1
	IF(L1.GT.NY) GO TO 100
	IF(A(IX,L1).EQ.XFILL) GO TO 10
C	FIND NEXT GOOD POINT
	L2=L1
	DO 50 N=2,NGOOD
11	L2=L2+1
	IF(L2.GT.NY) GO TO 100
	IF(A(IX,L2).EQ.XFILL) GO TO 11
	A1=A(IX,L1)
	A2=A(IX,L2)
      DY=L2-L1
      DADY=(A2-A1)/DY 
C
C	EXTRAPOLATE LOWER DEPTHS
C
	IF(N.NE.2.OR.L1.EQ.1) GO TO 21
	LTIME=L1-1
c
c	do not extrapolate more than 25% of values
c
	if(ltime.gt. 0.25*ny) go to 21
	DO 20 L=1,LTIME
	DYL=L-L1
20      A(IX,L)=A1+DADY*DYL
C
C	INTERPOLATE
C
21	LTIME=L2-L1-1
	IF(LTIME.LE.0) GO TO 31
	L11=L1+1
	L21=L2-1
	DO 30 L=L11,L21
	DYL=L-L1
30      A(IX,L)=A1+DADY*DYL
C
C	EXTRAPOLATE UPPER DEPTHS
C
31	IF(N.NE.NGOOD.OR.L2.EQ.NY) GO TO 50
	LTIME=NY-L2
c
c	do not extrapolate more than 25% of values
c
	if(ltime.gt. 0.25*ny) go to 50
	L21=L2+1
	DO 40 L=L21,NY
	DYL=L-L1
40      A(IX,L)=A1+DADY*DYL
50	L1=L2
100   CONTINUE
      RETURN
	END
      SUBROUTINE INTX(A,NX,NY,XFILL)
      DIMENSION A(NX,NY)
C
	DO 100 IY=1,NY
C	FIND FIRST GOOD COLUMN
	IX=0
10	IX=IX+1
	IF(IX.GT.NX) GO TO 100
	IF(A(IX,IY).EQ.XFILL) GO TO 10
	IX2=IX
20	IX2=IX2+1
	IF(IX2.GT.NX) GO TO 100
	IF(A(IX2,IY).EQ.XFILL) GO TO 20
	A1=A(IX,IY)
	A2=A(IX2,IY)
	DA=A2-A1
	DX=IX2-IX
	DADX=DA/DX
	DO 30 N=IX+1,IX2-1
30	A(N,IY)=A1+DADX*(N-IX)	
	IX=IX2
	GO TO 20
100	CONTINUE
	RETURN
	END
        subroutine dmpmat(a,nx,nxm,ny)
        dimension a(nxm,ny)
      NC1=20
C
C     WRITE OUT DATA
C
      DO 20 K=1,NX,NC1
      L1 = K
      L2 = L1 + NC1 -1
      IF (L2.GT.NX) L2 = NX
      WRITE (16,4) (N,N=L1,L2)
4     FORMAT (/11X,20i6)
      DO 20 I=1,NY
        WRITE (16, 6 ) I,(A(L,I),L=L1,L2)
    6 FORMAT (6X,i5,20F6.0)
20      CONTINUE
        return
        end
	subroutine init_array0(ii,nx,ny,nv,ifill)
	integer ii(nx,ny,nv)
	do 10 ix=1,nx
	do 10 iy=1,ny
	do 10 iv=1,nv
10	ii(ix,iy,iv)=ifill
	return
	end
