home *** CD-ROM | disk | FTP | other *** search
- Path: uunet!husc6!think!ames!sri-spam!rutgers!clyde!cbosgd!mandrill!hal!ncoast!allbery
- From: edward@pur-ee.UUCP (Edward L Haletky)
- Newsgroups: comp.sources.misc
- Subject: 2D contouring routines
- Message-ID: <6331@ncoast.UUCP>
- Date: 4 Dec 87 01:41:45 GMT
- Sender: allbery@ncoast.UUCP
- Organization: Purdue University Engineering Computer Network
- Lines: 437
- Approved: allbery@ncoast.UUCP
- X-Archive: comp.sources.misc/8712/2
-
-
-
-
- NetLanders;
-
- The contour routine is not big. It is actually quite small.
- The key to using the routines is as follows:
-
- The routines expect the data to be upon some plane of coordinates x
- and y or whatever is being used. The random points are located within this
- 2d mesh of x and y points. The routine requires that 3 2d arrays should be
- present the first array is the location of the x point in mesh coordinates,
- and so on. You can read about it in the Man page supplied. This could be
- used to do random point contouring by masking the unused coordinates within
- the arrays with some value that you know the random values will not touch.
- Either way above the max or way below the minimun value of the function to
- be plotted. The source and man page follow. It uses unix type plot calls
- from a package developed here at Purdue (CRC Graphics). I am sure your
- system has a similiar package....Ed
- ---------------8<------------------------8<--------------------8<-------------
- # This is a shell archive. Save this into a file, edit it
- # and delete all lines above this comment. Then give this
- # file to sh by executing the command "sh file". The files
- # will be extracted into the current directory owned by
- # you with default permissions.
- #
- # The files contained herein are:
- # UContour.f UContour.3G
- #
- echo 'x - UContour.f'
- sed 's/^X//' <<'________This_Is_The_END________' >>UContour.f
- Xc
- Xc UContour was originally written by Dr. W.J. Usab for the
- Xc PEPL graphics package. The routines have been modified to
- Xc use the CRC graphics package by Edward L. Haletky and Tom
- Xc Jentink.
- Xc
- Xc The following modifications were made: 11/13/87
- Xc Use of output files and not direct printing (gf)
- Xc gf options: filename - plot to file filename
- Xc scr - plot to standard output
- Xc lpr - plot to Printronix
- Xc
- Xc added ml option - specifies length of column in 2D matrix
- Xc Call to plots,fname were added to routines
- Xc Added Label option
- Xc Call to plot(0,0,999) added
- Xc Number call changed to reflect C syntax on formats
- Xc Transformed to lower case
- Xc
- Xc To Compile:
- Xc f77 -c UContour.f -lG
- Xc
- Xc
- Xc subroutine contour(x,y,f,n,m,cmin,cmax,inc)
- X subroutine ucontour(x,y,f,ml,in,jmax,cmin,cmax,inc,gf,
- X +lbl,bd)
- X real xmin,xmax,ymin,ymax,gscale
- X common /climit/ xmin,xmax,ymin,ymax,gscale
- X character *8 option
- X character *20 lbl
- X character *3 dev1,dev2,gf
- X integer dev,blank
- X real x(ml,*),y(ml,*),f(ml,*),bd(4),jmax(ml)
- Xc define constants
- X option = ' '
- X dev1 = "scr"
- X dev2 = "lpr"
- X dev = 0
- X if (gf.eq.dev1) dev = 6
- X if (gf.eq.dev2) dev = 16
- X blank = 0
- X xfac = 10.0
- X yfac = 10.0
- X
- Xc define min and max ranges for x,y,f
- X
- X call plots(dev,blank,option)
- X if ((dev.ne.16).and.(dev.ne.6)) call fname(gf)
- X call symbol(2,9.5,0.075,lbl,0.0)
- Xc
- Xc define min and max ranges for x,y,f
- Xc
- X xmin = x(1,1)
- X xmax = x(1,1)
- X ymin = y(1,1)
- X ymax = y(1,1)
- X fmin = f(1,1)
- X fmax = f(1,1)
- Xc
- Xc do 10 j = 1,m
- Xc do 10 i = 1,n
- X do 10 i=1,in
- X m=jmax(i)
- X do 10 j=1,m
- X xmin = min(x(i,j),xmin)
- X xmax = max(x(i,j),xmax)
- X ymin = min(y(i,j),ymin)
- X ymax = max(y(i,j),ymax)
- X fmin = min(f(i,j),fmin)
- X fmax = max(f(i,j),fmax)
- X 10 continue
- Xc
- X xscale = 1.0*xfac/(xmax-xmin)
- X yscale = 1.0*yfac/(ymax-ymin)
- X gscale = min(xscale,yscale)
- Xc
- Xc do 20 j = 1,m-1
- Xc do 20 i = 1,n-1
- X do 20 i=1,in-1
- X m=jmax(i)
- X do 20 j=1,m-1
- X x1 = gscale*(x(i,j)-xmin)
- X y1 = gscale*(y(i,j)-ymin)
- X f1 = f(i,j)
- X x2 = gscale*(x(i,j+1)-xmin)
- X y2 = gscale*(y(i,j+1)-ymin)
- X f2 = f(i,j+1)
- Xc******************************
- X if (j.eq.(m-1)) then
- X x3=gscale*x(i+1,j)
- X y3=gscale*y(i+1,j)
- X f3=f(i+1,j)
- X else
- X x3 = gscale*(x(i+1,j+1)-xmin)
- X y3 = gscale*(y(i+1,j+1)-ymin)
- X f3 = f(i+1,j+1)
- X endif
- Xc******************************
- X x4 = gscale*(x(i+1,j)-xmin)
- X y4 = gscale*(y(i+1,j)-ymin)
- X f4 = f(i+1,j)
- X diag1 = (x3-x1)**2+(y3-y1)**2
- X diag2 = (x4-x2)**2+(y4-y2)**2
- X if(diag1.le.diag2) then
- X call tricon(x1,y1,f1,x2,y2,f2,x3,y3,f3,
- X 1 cmin,cmax,inc)
- X call tricon(x1,y1,f1,x3,y3,f3,x4,y4,f4,
- X 1 cmin,cmax,inc)
- X else
- X call tricon(x1,y1,f1,x2,y2,f2,x4,y4,f4,
- X 1 cmin,cmax,inc)
- X call tricon(x2,y2,f2,x3,y3,f3,x4,y4,f4,
- X 1 cmin,cmax,inc)
- X endif
- Xc
- X 20 continue
- Xc
- Xc plot contour region boundary (clockwise from 0,0)
- Xc
- X x1 = gscale*(x(1,1)-xmin)
- X y1 = gscale*(y(1,1)-ymin)
- X call plot(x1,y1,3)
- Xc
- X m=jmax(1)
- X x1=gscale*x(1,m)
- X y1=gscale*y(1,m)
- X call plot(x1,y1,2)
- Xc
- X
- X x1=gscale*bd(1)
- X y1=gscale*bd(2)
- X call plot(x1,y1,2)
- Xc
- X x1=gscale*bd(3)
- X y1=gscale*bd(4)
- X call plot(x1,y1,2)
- Xc
- X if (bd(3).lt.x(in,1)) then
- X x1=gscale*x(in,1)
- X y1=gscale*bd(4)
- X call plot(x1,y1,2)
- Xc
- X x1 = gscale*x(in,1)
- X y1 = gscale*y(in,1)
- X call plot(x1,y1,2)
- Xc
- X x1 = gscale*(x(1,1)-xmin)
- X y1 = gscale*(y(1,1)-ymin)
- X call plot(x1,y1,2)
- X endif
- X call plot(0,0,999)
- Xc
- Xc
- X return
- X end
- Xc
- X subroutine tricon(xa,ya,fa,xb,yb,fb,xc,yc,fc,cmin,cmax,inc)
- X real x(3),y(3),f(3)
- Xc
- X x(1) = xa
- X x(2) = xb
- X x(3) = xc
- X y(1) = ya
- X y(2) = yb
- X y(3) = yc
- X f(1) = fa
- X f(2) = fb
- X f(3) = fc
- Xc
- Xc sort points in ascending order
- Xc
- X do 10 n=1,2
- X do 10 i=1,3-n
- X if(f(i).gt.f(i+1))then
- X xt = x(i)
- X yt = y(i)
- X ft = f(i)
- X x(i) = x(i+1)
- X y(i) = y(i+1)
- X f(i) = f(i+1)
- X x(i+1) = xt
- X y(i+1) = yt
- X f(i+1) = ft
- X endif
- X 10 continue
- Xc
- Xc locate contour range within triangle
- Xc
- X dc = (cmax-cmin)/float(inc)
- X icmin = int((f(1)-cmin)/dc)
- X icmin = max(icmin,0)
- X icmax = int((f(3)-cmin)/dc)
- X icmax = min(icmax,inc)
- X if(icmax.lt.icmin) return
- Xc
- Xc plot contours within range
- Xc
- X do 20 i = icmin,icmax
- Xc do 20 i=0,inc
- X clev = cmin+dc*float(i)
- Xc if(clev.lt.f(1)) go to 20
- X if (clev.le.f(1)) goto 20
- Xc if(clev.gt.f(3)) go to 20
- X if (clev.ge.f(3)) goto 20
- X if(f(2).gt.clev) then
- X fac = (clev-f(1))/(f(2)-f(1))
- X x1 = x(1)+fac*(x(2)-x(1))
- X y1 = y(1)+fac*(y(2)-y(1))
- X else
- X fac = (clev-f(2))/(f(3)-f(2))
- X x1 = x(2)+fac*(x(3)-x(2))
- X y1 = y(2)+fac*(y(3)-y(2))
- X endif
- X fac = (clev-f(1))/(f(3)-f(1))
- X x2 = x(1)+fac*(x(3)-x(1))
- X y2 = y(1)+fac*(y(3)-y(1))
- Xc
- X call plot(x1,y1,3)
- X call plot(x2,y2,2)
- Xc call labinc(i,clev,x1,y1,x2,y2)
- X 20 continue
- Xc
- X return
- X end
- Xc
- X subroutine labinc(i,clev,x1,y1,x2,y2)
- Xc common /clabinc/ icount
- Xc data icount/0/
- X tdist = abs(x2-x1)+abs(y2-y1)
- X if(tdist.lt.1e-20) then
- X write(6,*)' ** error in labinc: tdist=',tdist
- X return
- X endif
- X interv = 21
- X icount = icount+1
- X if(icount.lt.interv) return
- X icount = 0
- X pi = acos(-1.0)
- X hpi = 0.5*pi
- X theta = atan2((y2-y1),(x2-x1))
- X if(theta.ge.hpi) then
- X theta = theta-pi
- X else if(theta.lt.-hpi) then
- X theta = theta+pi
- X endif
- X thetan = theta+hpi
- X x = 0.5*(x1+x2)+0.05*cos(thetan)-0.225*cos(theta)
- X y = 0.5*(y1+y2)+0.05*sin(thetan)-0.225*sin(theta)
- X thetad = 180.0*theta/pi
- Xcc call number(x,y,.1,i,0.0,'i2')
- Xc call number(x,y,.075,clev,thetad,'f6.3')
- Xc call number(x,y,.055,thetad,"%6.3f",clev)
- X return
- X end
- ________This_Is_The_END________
- echo 'x - UContour.3G'
- sed 's/^X//' <<'________This_Is_The_END________' >>UContour.3G
- X.\" troff -man
- X.EN
- X.TH UContour 3G 11/13/87
- X.SH NAME
- XUContour - EGraphics Contour Routine.
- X.SH AUTHORS
- XDr. W.J.Usab Jr.
- X.br
- XEdward L. Haletky
- X.br
- XTom Jentink
- X.SH SYNOPSIS
- Xcall ucontour(x,y,f,ml,in,jmax,cmin,cmax,inc,gf,lbl,bd)
- X.SH DESCRIPTION
- X.fi
- XUContour will plot a family of curves of number
- X.I inc
- Xwithin the range
- X.I [cmin:cmax].
- XAll the curves will be of a constant value of the function,
- X.I f.
- X.sp
- X.RS
- X.fi
- X.B x
- X- two dimensional array for x-coordinate at grid point (n,m) at which
- Xthe function,
- X.I f,
- Xwas computed.
- X.sp 1
- X.B y
- X- two dimensional array for y-coordinate at grid point (n,m) at which
- Xthe function,
- X.I f,
- Xwas computed.
- X.sp 1
- X.B f
- X- two dimensional array for the function at grid point (n,m).
- X.sp 1
- X.B ml
- X- maximun length of the arrays
- X.I x,
- X.I y,
- X.I f.
- XThe arrays are declared in the following manor:
- X.br
- X.ce 1
- Xreal x(ml,*),y(ml,*),f(ml,*)
- X.br
- X.sp 1
- X.B in
- X- maximun number of grid points in n-direction.
- X.sp 1
- X.B jmax
- X- an array that gives maximun number of grid points in the m-direction,
- Xfor each,
- X.I n
- Xstation.
- X.sp 1
- X.B cmin
- X- minimun value of the constant-value contour of the function,
- X.I f.
- X.sp 1
- X.B cmax
- X- maximun value of the constant-value contour of the function,
- X.I f.
- X.sp 1
- X.B inc
- X- integer number of contours between the range
- X.I [cmin:cmax].
- X.sp 1
- X.B gf
- X- device which to plot
- X.RE
- X.in +10
- X.B scr
- X- plot to stdout.
- X.br
- X.B lpr
- X- plot to Printronix Printer.
- X.br
- X.B gf
- X- plot to file. File can be printed with the
- X.I gplp(1G)
- Xcommand. The
- X.I filename
- Xcan not be more than 3 characters long.
- X.in -10
- X.RS
- X.sp 1
- X.B lbl
- X- label the contour plot.
- X.B lbl
- Xcan be a maximun of 20 characters long.
- X.br
- X.sp 1
- X.B bd
- X- an array that holds the bounding area values.
- X.RE
- X.in 10
- X.B bd(1)
- X- xb1
- X.br
- X.B bd(2)
- X- yb1
- X.br
- X.B bd(3)
- X- xb2
- X.br
- X.B bd(4)
- X- yb2
- X.br
- X.in -10
- X.RS
- X.sp 1
- X.B "To Compile"
- Xf77 prog.f -L/x/edward/lib -lE -lG
- X.RE
- X.SH FILES
- X.RS
- X/x/edward/lib/libE.a
- X.br
- X/x/edward/src/UContour.f
- X.RE
- X.SH DIAGNOSTICS
- X.LP
- XUContour will blank the page before plotting. No CRC graphics options are used.
- XThe Border of the Channel will be printed, it is a crude representation. The
- Xlabeling of the lines will be implemented at a later date. There are a few bugs
- Xstill.
- X.SH "SEE ALSO"
- X.B "CRC Graphics Manaul ,"
- X.B gplp(1G)
- X.SH BUGS
- X.RS
- XPlease Report all bugs to:
- X.br
- Xedward@ga.ecn.purdue.edu
- X.RE
- ________This_Is_The_END________
- exit
- ---------------8<------------------------8<--------------------8<-------------
-
-
- ===================================== ^ =================================
- Edward L. Haletky |E |o| U| "To race against the sky..."
- Usenet: ~!pur-ee!edward |L ^/| |\^ S| "The world of the
- Arpa: edward@ga.ecn.purdue.edu |H / | | \ A| immagination is boundless."
- =====================================/__-|-__\=================================
-