c  program extract
c  mwebring, 6/83
	dimension za(20)
	character*56 ifile,ofile,ifmt
	character*1 charst(80)
	character record*80,id*8,cfmt*16,rfmt*16
	equivalence (record,charst)
	ic=0
	jc=0
	ifmt='(22x,f2.0,f4.2,f4.0,f4.2,t1,a80)'
	cfmt='(80a1)'
	rfmt='(a80)'
	write(*,10)
10	format(' enter filename :'$)
	read(5,37) ifile
37	format(a56)
	write(*,15)
15	format(' enter output filename :'$)
	read(5,37) ofile
20	write(*,30)
30	format(' post=1, bin xyz=2, form xyz=3, std grv=4 :'$)
	read(5,*,err=35) ifunc
	if(ifunc.ge.1 .and. ifunc.le.4) go to 50
35	write(*,40)
40	format(' enter file type :',/,' binary post = 1',/,
     1 ' binary xyz = 2',/,' formatted xyz = 3',/,
     1 ' formatted standard gravity = 4')
	go to 20
50	if(ifunc.le.2) then
	 open(10,file=ifile,status='old',form='unformatted',
     1 share='denywr')
	 nch = lenz(10)
	 nch = nch - 2
	 open(11,file=ofile,status='unknown',form='unformatted')
	else
	 open(10,file=ifile,status='old',form='formatted',
     1 share='denywr')
	 open(11,file=ofile,status='unknown',form='formatted')
	endif
c
	if(ifunc.eq.3) then
	 write(*,70)
70	 format(' enter format :'$)
	 read(5,37) ifmt
	else if(ifunc.eq.4) then
	 do 80 i=1,3
	 read(10,rfmt,end=999) record
80	 write(11,rfmt) record
	endif
c
	if(ifunc.ne.4) then
	write(*,100)
100	format(' enter min x, max x :'$)
	read(5,*) xmin,xmax
	write(*,101)
101	format(' enter min y, max y :'$)
	read(5,*) ymin,ymax
	else
	write(*,102)
102	format(' enter left and right side longitude d,m, d,m :'$)
	read(5,*) a,b,c,d
	xmin = a + sign(b,a) / 60.
	xmax = c + sign(d,c) / 60.
	write(*,103)
103	format(' enter lower and upper latitude d,m, d,m :'$)
	read(5,*) a,b,c,d
	ymin = a + sign(b,a) / 60.
	ymax = c + sign(d,c) / 60.
	endif
	write(*,*)' xmin,xmax =', xmin, xmax
	write(*,*)' ymin,ymax =', ymin, ymax
	if (xmin.ge.xmax) stop ' xmin>xmax'
	if (ymin.ge.ymax) stop ' ymin>ymax'
c
	go to (110,120,130,140) ifunc
110	call io1(id,x,y,za,nch,iend)
	if (iend.eq.1) go to 300
	go to 150
120	read(10,end=300) x,y,z
	go to 150
130	read(10,ifmt,end=300) x,y,z
	go to 150
140	read(10,ifmt,end=300) yd,ym,xd,xm,record
	x = xd + xm*sign(1.666667e-2,xd)
	y = yd + ym*1.666667e-2
c
150	ic=ic+1
151	if (x.gt.xmax .or. x.lt.xmin) go to 260
152	if (y.gt.ymax .or. y.lt.ymin) go to 260
c
200	jc=jc+1
	go to (210,220,230,240) ifunc
210	call io2(id,x,y,za,nch)
	go to 260
220	write(11) x,y,z
	go to 260
230	write(11,ifmt) x,y,z
	go to 260
240	do 245 lenth=80,1,-1
	ich=ichar(charst(lenth))-32
	if(ich) 246,245,246
245	continue
246	write(11,cfmt) (charst(i),i=1,lenth)
260	go to (110,120,130,140) ifunc
c
c
300	write(*,350)ic,jc
350	format(1x,i8,' records input ',i8,' records output ')
	close(10)
	close(11)
	stop
999	close(10)
	stop ' premature end of file'
	end
c*******************************************************************************
	subroutine io1(id,x,y,za,nch,iend)
	dimension za(nch)
	character id*8
	iend = 0
	read(10,end=300) id,x,y,za
	return
300	iend = 1
	return
	end
c*******************************************************************************
	subroutine io2(id,x,y,za,nch)
	dimension za(nch)
	character id*8
	write(11) id,x,y,za
	return
	end
c******************************************************************************
	function lenz(iswt)
c  find the number of channels in a post record
	dimension id(2),z(50)
	do 15 lenz = 50, 1, -1
	  read(iswt,err=10) id,(z(i),i=1,lenz)
	  go to 20
10	  rewind iswt
15	continue
20	rewind iswt
	return
	end

