c   test isometric routines
	program tstiso
	character*1 back
c
	isom=0
	call at_bint('Number of som (no contour) plots',isom,back)
	do ii=1,isom
	  call somplot
	enddo
c
	isoc=0
	call at_bint('Number of soc (contour) plots',isoc,back)
	do ii=1,isoc
	  call socplot
	enddo
c
#ifdef __VMS
	call at_bint('Number of ufile plots',iuf2,back)	
	do ii=1,iuf2		
	  call uplot		
	enddo		
	call at_bint('Number of matrix plots',imat,back)	
	do ii=1,imat		
	  call mplot		
	enddo		
#endif
c
	end
c
	subroutine somplot
	USE SGLIB_MOD
c
	parameter (nxp=10)
	parameter (nyp=10)
	real xdata(nxp)
	real ydata(nyp)
	real zdata(nxp,nyp),zbase(nxp,nyp)
c
	real xmin, xmax 
	real zmin, zmax	!zmin><0
	real zscale
	integer linc
	character*12 xuns	
	character*12 yuns	
	character*12 zuns
	data xmin /10./, xmax /100./
	data zmin /-20.0/, zmax /60.0/	!zmin><0
	data zscale /1./
	data linc /1/
	data xuns	/'X data$'/
	data yuns	/'Y data$'/
	data zuns	/'Z data$'/
	logical zlog, lhide
	character*1 back
	data (xdata(i),i=1,nxp) /10.,20.,30.,40.,50.,60.,70.,80.,90.,100./
	data ((zbase(i,j),i=1,nxp),j=1,nyp) /
     >		60.,40.,20.,00.,-20.,-20.,00.,20.,40.,60.,
     >		60.,40.,20.,00.,-20.,-20.,00.,20.,40.,60.,
     >		62.,42.,22.,02.,-18.,-18.,02.,22.,42.,62.,
     >		64.,44.,24.,04.,-16.,-16.,04.,24.,44.,64.,
     >		66.,46.,26.,06.,-14.,-14.,06.,26.,46.,66.,
     >		66.,46.,26.,06.,-14.,-14.,06.,26.,46.,66.,
     >		64.,44.,24.,04.,-16.,-16.,04.,24.,44.,64.,
     >		62.,42.,22.,02.,-18.,-18.,02.,22.,42.,62.,
     >		60.,40.,20.,00.,-20.,-20.,00.,20.,40.,60.,
     >		60.,40.,20.,00.,-20.,-20.,00.,20.,40.,60./
	yfact=1.
	call at_breal('Y axis scale',yfact,back)
	do jj=1,nxp		!remote exponent on y axis label
	  ydata(jj) = jj *.5 * abs(yfact)
	enddo
	call at_breal('Z scale',zscale,back)
	zmin = 0.
	zmax = 0.
	do jj=1,nyp
	  do kk=1,nxp
	    zdata(kk,jj) = zbase(kk,jj)*zscale
	    zmin = min(zmin,zdata(kk,jj))
	    zmax = max(zmax,zdata(kk,jj))
	  enddo
	enddo
	zlog=.false.
	call at_log ('Z axis as log',zlog)
	if(zlog .and. zmin.lt.0) zmin = 1.
	lhide=.true.
	call at_log('Remove hidden lines',lhide)
	call at_int('Line color',linc)
	theta=45.
	call at_real('Y angle',theta)
	call initt(120)
	call aginit
	if(linc.ne.1) call agccurv(linc)
c
	call agidlx(xmin,xmax)
cx	call agidly(ymin,ymax)
	call agidlz(zmin,zmax)
	call agixnt(0)
	call agiynt(0)
	call agiznt(0)
	if(theta.gt.90) call agislx(430,870)
	call agisly(theta,0)
	call agiskt(.false.,.true.)
	call agihid(lhide)
	if(zlog) call agizty(2)
c
	call agisom(xdata,ydata,zdata,nxp,nxp,nyp,
     >		xuns,yuns,zuns)
	call upause(0)
	call finitt(0,0)
	end
c
	subroutine socplot
c
	USE SGLIB_MOD
c
	parameter (nxp=10)
	parameter (nyp=10)
	real xdata(nxp)
	real ydata(nyp)
	real zdata(nxp,nyp),zbase(nxp,nyp)
	integer m(2,nxp,nyp)	!only for agisoc
	real    c(nyp)
	logical lhide
	real xmin, xmax
	real zmin, zmax
	real zscale
	integer linc
	character*12 xuns	
	character*12 yuns	
	character*12 zuns	

	data xmin /10./, xmax /100./
	data zmin /-60.0/, zmax /20.0/	!zmin><0
	data zscale /1./
	data linc /1/
	data xuns	/'X data$'/
	data yuns	/'Y data$'/
	data zuns	/'Z data$'/
	data (xdata(i),i=1,nxp) /10.,20.,30.,40.,50.,60.,70.,80.,90.,100./
	data ((zbase(i,j),i=1,nxp),j=1,nyp) /
     >		20.,00.,-20.,-40.,-60.,-60.,-40.,-20.,00.,20.,
     >		20.,00.,-20.,-40.,-60.,-60.,-40.,-20.,00.,20.,
     >		22.,02.,-18.,-38.,-58.,-58.,-38.,-18.,02.,22.,
     >		24.,04.,-16.,-36.,-56.,-56.,-36.,-16.,04.,24.,
     >		26.,06.,-14.,-34.,-54.,-54.,-34.,-14.,06.,26.,
     >		26.,06.,-14.,-34.,-54.,-54.,-34.,-14.,06.,26.,
     >		24.,04.,-16.,-36.,-56.,-56.,-36.,-16.,04.,24.,
     >		22.,02.,-18.,-38.,-58.,-58.,-38.,-18.,02.,22.,
     >		20.,00.,-20.,-40.,-60.,-60.,-40.,-20.,00.,20.,
     >		20.,00.,-20.,-40.,-60.,-60.,-40.,-20.,00.,20./
c
	do jj=1,nxp
	  ydata(jj) = jj *.5
	enddo
	do jj=1,nyp
	  do kk=1,nxp
	    zdata(kk,jj) = zbase(kk,jj)*zscale
	  enddo
	enddo
	lhide=.true.
	call at_log('Remove hidden lines',lhide)
	nc=-1
	call at_int('Contours (-1 for uniform, 0 for none)',nc)
	call at_int('Line color',linc)
	theta=45.
	call at_real('Y angle',theta)
	call initt(120)
	call aginit
	if(linc.ne.1) call agccurv(linc)
c
	call agidlx(xmin,xmax)
cx	call agidly(ymin,ymax)
	call agidlz(zmin,zmax)
	call agixnt(0)
	call agiynt(0)
	call agiznt(0)
	call agiskt(.false.,.true.)
	if(theta.gt.90) call agislx(430,870)
	call agisly(theta,0)
	call agihid(lhide)
	call agiszc(nc,c)
c
	call agisoc(xdata,ydata,zdata,nxp,nxp,nyp,
     >		xuns,yuns,zuns,m)
	call upause(0)
	call finitt(0,0)
	end
c
#ifdef __VMS
	subroutine uplot 
	parameter (mxd=2050)	!2050 
	parameter (myd=640) 
	real xdata(mxd) 
	real ydata(myd) 
	real zdata(mxd,myd) 
c 
	logical lhide, zlog, xlog, ylog 
	character fspec*48,fpref*24,fextn*4,fdisk*24,farea*24,fshot*5 
	character xuns*30,yuns*30,zuns*30 
	character tdev*4,tdate*8,suns*30 
	data fspec /'res$:[tftr.fm]a44008.fma'/ 
	!res$:[tftr.ym]t49120.ecm 
c 
	call at_char('File spec',fspec) 
	call ufspec(fspec,fpref,fextn,fdisk,farea,fshot,jshot) 
c 
	call ufsetr(8,fpref,fextn,fdisk,farea) 
	call ufoprd(8,jshot,ier) 
	call uf2drd(8,tdev,tdate,zdata,mxd,xdata,mxd,ydata,myd,iproc, 
     >	  nx,ny,0,zuns,xuns,yuns,1,nsc,scval,suns,ier) 
	call ufclos(8) 
	if(ier.ne.0) then 
	  type *,'cant open file',ier 
	  return 
	endif 
	call at_msgi('  x points',nx) 
	call at_msgi('  y points',ny) 
	kxst=max(nx/20,1) 
	call at_int('X step',kxst) 
	kyst=max(ny/20,1) 
	call at_int('Y step',kyst) 
	call at_real('Y window, low',ylo) 
	call at_real('Y window, high',yhi) 
	call at_real('X window, low',xlo) 
	call at_real('X window, high',xhi) 
	zlog=.false. 
	call at_log ('Z axis as log',zlog) 
	xlog=.false. 
	call at_log ('X axis as log',xlog) 
	ylog=.false. 
	call at_log ('Y ayis as log',ylog) 
	theta=45. 
	call at_real('Y angle',theta) 
	lhide=.true. 
	call at_log('Remove hidden lines',lhide) 
	iside=1 
	call at_int('0 for top&bottom, 1 for top only',iside) 
	neat=1 
	call at_int('0 to extend to data, 1 to minor ticks',neat) 
c 
	call initt(120) 
	call aginit 
c 
	call agixnt(neat) 
	call agiynt(neat) 
	call agiznt(neat) 
	call agistx(kxst) 
	call agisty(kyst) 
	if(ylo.ne.yhi) call agidly(ylo,yhi) 
	if(xlo.ne.xhi) call agidlx(xlo,xhi) 
	call agisid(iside) 
	call agiskt(.false.,.true.) 
	if(zlog) call agizty(2) 
	if(xlog) call agixty(2) 
	if(ylog) call agiyty(2) 
	call agihid(lhide) 
	call agislx(110,510)	!match RPLOT/UGRAF2 
	call agislz(090,350) 
	if(theta.eq.45.) then 
	  call agisly(theta,550) 
	elseif(theta.gt.90.) then 
	  call agislx(430,870) 
	  call agisly(theta,0) 
	else 
	  call agisly(theta,0) 
	endif 
c 
	call agisom(xdata,ydata,zdata,mxd,nx,ny, 
     >	  xuns(21:),yuns(21:),zuns(21:)) 
	ll=jc_trim(fspec) 
cx	call agtmhs(400,760,fspec(1:ll+1)//tdate) 
	call agtcgt(760,fspec(1:ll+1)//tdate) 
	call agtmhs(20,740,zuns(1:20)) 
	call upause(0) 
	call finitt(0,0) 
	end 
c 
	subroutine mplot 
	parameter (maxraw=400) 
	real xraw(maxraw),yraw(maxraw),zraw(maxraw) 
	parameter (mxd=100) 
	parameter (myd=100) 
	real xdata(mxd),ydata(myd),zdata(mxd,myd) 
	character*1 back 
	character fname*48 /'kowens$:[toci]bay_avg.dat'/ 
	call at_bchar('File name',fname,back) 
	if(back.ne.' ') return 
	open(8,file=fname, 
     >	  access='sequential', 
     >	  form='formatted', 
     >	  status='old', 
     >	  readonly,shared, 
     >	  iostat=ierr) 
	if(ierr.ne.0) then 
	  call at_msgi('** cant open '//fname,ierr) 
	  return 
	endif 
	ierr=0 
	nraw=0 
	do while (ierr.eq.0) 
	  nraw=nraw+1 
	  read(8,*,iostat=ierr) xraw(nraw),yraw(nraw),zraw(nraw) 
	  if(ierr.ne.0) nraw=nraw-1 
	enddo 
	close(8) 
	nx=mxd 
	call at_int('Number of X points',nx) 
	ny=myd 
	call at_int('Number of Y points',ny) 
	ncell=2 
	call at_int('Furthest cell to use',ncell) 
	wgt=2. 
	call at_real('Weighting factor',wgt) 
	call makmat(xraw,yraw,zraw,nraw,xdata,nx,ydata,ny,zdata,mxd 
     >	  ,ncell,wgt) 
c 
	neat=1 
	lhide=.true. 
	theta=45. 
	call initt(120) 
	call aginit 
c 
	call agixnt(neat) 
	call agiynt(neat) 
	call agiznt(neat) 
cx	call agistx(kxst) 
cx	call agisty(kyst) 
	call agisid(1)	!top only 
	call agiskt(.false.,.true.) 
	call agihid(lhide) 
	call agislx(110,510)	!match RPLOT/UGRAF2 
	call agislz(090,350) 
	if(theta.eq.45.) then 
	  call agisly(theta,550) 
	elseif(theta.gt.90.) then 
	  call agislx(430,870) 
	  call agisly(theta,0) 
	else 
	  call agisly(theta,0) 
	endif 
c 
	call agisom(xdata,ydata,zdata,mxd,nx,ny,'X','Y','Z') 
	call agtmhs(400,760,fname) 
	call upause(0) 
	call finitt(0,0) 
	end 
c 
	subroutine makmat(x,y,z,nxyz,xmat,nx,ymat,ny,zmat,nzd,ncell,wgt) 
	real x(*),y(*),z(*) 
	real xmat(nx),ymat(ny),zmat(nzd,ny) 
	nc=ncell 
	if(nc.eq.0) nc=2 
	w=wgt 
	if(w.le.0.) w=2. 
	call jc_amnmx(x,nxyz,xmin,xmax) 
	call jc_amnmx(y,nxyz,ymin,ymax) 
	xval=xmin 
	xdel=(xmax-xmin)/(nx-1) 
	do i=1,nx 
	  xmat(i)=xval 
	  xval=xval+xdel 
	enddo 
	yval=ymin 
	ydel=(ymax-ymin)/(ny-1) 
	do i=1,ny 
	  ymat(i)=yval 
	  yval=yval+ydel 
	enddo 
	xfact=1./xdel 
	yfact=1./ydel 
	do i=1,nx 
	  do j=1,ny 
	    ztop=0. 
	    zbot=0. 
	    do k=1,nxyz 
	      xd=(x(k)-xmat(i))*xfact 
	      ixd=abs(xd)	!convert to number of cells 
	      yd=(y(k)-ymat(j))*yfact 
	      iyd=abs(yd) 
	      if(ixd.le.nc .and. iyd.le.nc) then 
		dw=sqrt((xd/(nx-1))**2+(yd/(ny-1))**2)**w 
		ztop=ztop+z(k)*dw 
		zbot=zbot+dw 
	      endif 
	    enddo 
	    if(ztop.ne.0) ztop=ztop/zbot 
	    zmat(i,j)=ztop 
	  enddo 
	enddo 
	end 
       subroutine jc_amnmx(a,n,a1,a2)  !find minimum and maximum  !01/15/85
        real a(*)
        t1=a(1)         !faster than using a1/a2 inside loop
        t2=t1
        do i=2,n
          t1=min(t1,a(i))
          t2=max(t2,a(i))
        enddo
        a1=t1
        a2=t2
        end
#endif







