c  SHODVALS
c  This program roughly shows nodata areas of a grid on the terminal.
c  Designed to give a quick check of merging or extracting operations
c  when working with irregular grids.
c  Tien Grauch
c
      character infile*50,istar*1,idot*1,ofile*50
      character*1 prnt(80)
      dimension id(14),pgm(2),f(5000),index(132),in(3,5000),iplace(3)
c if increase no. of cols from 5000, need to check 1st dimensiion
c of array in and array iplace because they are the alog10 of ncol
      data dval/1.0e30/
      istar='*'
      idot='.'
      print 820
820   format(' This program roughly shows no-data areas of a grid on',
     1 ' the terminal or printer'/' Output is scaled to page size ',
     2 'although you can specify vertical scale'/' Note that listing '
     3 'starts with row 1, which is the southernmost row,'/' so the '
     4 'grid appears inverted.')
11    print 822
822   format(/' Enter input file: ',$)
      read 804, infile
804   format(a50)
      open(10,file=infile,status='old',form='unformatted',mode='read',
     &err=299)
       read(10) id,pgm,ncol,nrow,nz,x0,dx,y0,dy
      if(nz.eq.1.and.ncol.le.5000) go to 1
      print*, 'nz must=1 and ncol must <=5000'
      go to 300
1	print*,'Scale vertical to how many lines?'
	print 814, nrow
814	format(3x,'(Enter 0 for default, ',i10,' for all rows): ',$)
	read*,jrow
	if(jrow.lt.0) go to 1
	write(6,810)
810	format(' Output file for the printer (132-char lines)',
     1 ' (optional):')
	read 804, ofile
	iunit=6
	imax=72
	jmax=30
	if(ofile.eq.' ') go to 5
	iunit=12
	imax=126
	jmax=48
	open(12,file=ofile,form='formatted',status='new')
c
5      	if(jrow.ne.0) jmax=jrow
c
	write(iunit,803) infile
803   format(1x,a50/'  . = data points,  * = dvals',/,5x,'Grids come out
     1 upside-down (but not backwards); that is,',/,'    south comes
     2 out on top and east comes out on the left.')
      write(iunit,802)
802   format(1x)
	write(iunit,808)
808	format(' Actual'/' Row No.')
      scalex=1.0
      scaley=1.0
	if(nrow.gt.jmax) then
		scaley=float(jmax)/float(nrow)
		jend=jmax
	else
		jend=nrow
	endif
      if(ncol.gt.imax) then
		scalex=imax/float(ncol)
		iend=imax
	else
	      iend=ncol
	endif
	jndex2=0
      index2=0
	incount=0
      do 40 i=1,ncol
      ind=ifix(scalex*float(i))
	if(ind.lt.1) go to 40
      if(ind.eq.index2.or.ind.gt.iend) go to 40
	incount=incount+1
	index(incount)=i
      index2=ind
40     continue
      do 50 j=1,nrow
      read(10,end=50) dum,(f(i),i=1,ncol)
      jndex=ifix(scaley*float(j))
      if(jndex.gt.jend.or.jndex.eq.jndex2) go to 50
      jndex2=jndex
       do 35 k=1,imax
35    prnt(k)=idot
      do 45 i=1,incount
	ii=index(i)
      if(f(ii).ge.dval) prnt(i)=istar
45     continue
	if(iunit.eq.6) then
           write(iunit,801) j,(prnt(iend+1-i),i=1,iend)
801        format(1x,i5,1x,72a1)
	else
           write(iunit,811) j,(prnt(iend+1-i),i=1,iend)
811        format(1x,i5,1x,126a1)
	endif
50    continue
c figure out how to label columns, first break index array into tens pieces
	iorder= ifix(alog10(float(ncol)))
	do 65 i=1,incount
	call places(index(i),iorder,iplace)
	do 60 k=1,iorder+1
	in(k,i)=iplace(k)
60 	continue
65	continue
	do 75 k=1,iorder+1
	write(iunit,830) (in(k,incount+1-i),i=1,incount)
830	format(7x,132i1)
75	continue
	write(iunit,831)
831	format(15x,'Actual column numbers written vertically')
       go to 300
299   print*,'File incompatible or not found.  Try again'
      go to 11
300   close(10)
	if(iunit.eq.12) close(12)
       stop
      end
c *******************
c subroutine places breaks an integer down into its units places
c
c  integer= the integer to break down
c  iorder= the maximum order that the integer could be
c  iplace= array of dimension iorder that has the breakdown,
c           where the first integer in the array is the ones place
c
	subroutine places(integer,iorder,iplace)
	dimension iplace(iorder+1)
	tenth=10.0**(-iorder)
	iplace(1)= ifix(integer*tenth)
	if(iorder.eq.0) return
	ten=10.0**iorder
	number=integer
	do 10 k=2,iorder+1
	tenth=tenth*10.0
	number= number - iplace(k-1)* ifix(ten)
	iplace(k)= ifix( number * tenth)
	ten=ten*0.1
10	continue
	return
	end

