v15i024: Graphics Gems, Part 2/5
Craig Kolb
craig at weedeater.math.yale.edu
Mon Oct 15 11:12:16 AEST 1990
Posting-number: Volume 15, Issue 24
Submitted-by: Craig Kolb <craig at weedeater.math.yale.edu>
Archive-name: ggems/part02
#! /bin/sh
# This is a shell archive. Remove anything before this line, then unpack
# it by saving it into a file and typing "sh file". To overwrite existing
# files, type "sh file -c". You can also feed this as standard input via
# unshar, or by typing "sh <file", e.g.. If this archive is complete, you
# will see the following message at the end:
# "End of archive 2 (of 5)."
# Contents: 2DClip/bio.c 2DClip/clip.c BoxSphere.c Label.c Makefile
# MatrixPost.c Median.c OrderDither.c PntOnLine.c
# PolyScan/fancytest.c PolyScan/poly.c PolyScan/poly.h
# PolyScan/scantest.c Quaternions.c SeedFill.c SquareRoot.c
# Sturm/main.c TriPoints.c
# Wrapped by craig at weedeater on Fri Oct 12 15:53:12 1990
PATH=/bin:/usr/bin:/usr/ucb ; export PATH
if test -f '2DClip/bio.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'2DClip/bio.c'\"
else
echo shar: Extracting \"'2DClip/bio.c'\" \(3146 characters\)
sed "s/^X//" >'2DClip/bio.c' <<'END_OF_FILE'
X
X/*
X * file bio.c
X * contains the basic operations
X *
X */
X#include <stdio.h>
X#include "GraphicsGems.h"
X#include "line.h"
X
X/*
X * def_contour
X *
X * Purpose:
X * add a contour to the list
X * NOTE: coordinates must already be converted into longs!
X *
X * x x coordinates of the end points of segments
X * y y coordinates of the end points of segments
X * n number of coordinate pairs
X * no contour number (id), does no have to be unique!
X * type type of clip operation desired CLIP_NORMAL means
X * clip everything inside the contour
X */
Xdef_contour(x, y, n, no, type)
Xlong x[], y[];
Xint n, no, type;
X{
X short i;
X long dx1, dx2, dy1, dy2;
X long minx, miny, maxx, maxy;
X CONTOUR *cp;
X SEGMENT *sp, *spp;
X
X if((cp = CL) == (CONTOUR *)NULL) {
X cp = CL = NEWTYPE(CONTOUR);
X }
X else {
X while(cp->_next != (CONTOUR *)NULL)
X cp = cp->_next;
X i = cp->_no;
X cp = cp->_next = NEWTYPE(CONTOUR);
X }
X
X cp->_next = (CONTOUR *)NULL;
X cp->_no = no;
X SET_ON(cp);
X if(type == CLIP_NORMAL)
X SET_INVERSE(cp);
X else
X SET_NORMAL(cp);
X minx = miny = 1000000;
X maxx = maxy = -1;
X for(i=0; i<n; i++) {
X if(i == 0) {
X cp->_s = sp = NEWTYPE(SEGMENT);
X sp->_from._x = x[0];
X sp->_from._y = y[0];
X sp->_to._x = x[1];
X sp->_to._y = y[1];
X }
X else {
X /*
X * if necessary stretch the contour
X * and skip the point
X */
X dx1 = sp->_to._x - sp->_from._x;
X dy1 = sp->_to._y - sp->_from._y;
X dx2 = x[(i == n-1) ? 0 : i+1] - sp->_to._x;
X dy2 = y[(i == n-1) ? 0 : i+1] - sp->_to._y;
X if(dy2*dx1 == dy1*dx2) {
X sp->_to._x = x[(i == n-1) ? 0 : i+1];
X sp->_to._y = y[(i == n-1) ? 0 : i+1];
X }
X else {
X spp = sp;
X sp = sp->_next = NEWTYPE(SEGMENT);
X sp->_prev = spp;
X sp->_from._x = x[i];
X sp->_from._y = y[i];
X sp->_to._x = x[(i == n-1) ? 0 : i+1];
X sp->_to._y = y[(i == n-1) ? 0 : i+1];
X }
X }
X
X/*
X * calculate the enclosing box
X */
X if(x[i] < minx)
X minx = x[i];
X if(x[i] > maxx)
X maxx = x[i];
X if(y[i] < miny)
X miny = y[i];
X if(y[i] > maxy)
X maxy = y[i];
X
X }
X cp->_minx = minx;
X cp->_maxx = maxx;
X cp->_miny = miny;
X cp->_maxy = maxy;
X sp->_next = cp->_s;
X cp->_s->_prev = sp;
X cp->_next = (CONTOUR *)NULL;
X}
X
X/*
X * get_contour_ptr
X *
X * PURPOSE
X * get the pointer to a contour given its id
X * with multiple id's first fit algorithm is
X * used. Returns NULL in case of error.
X *
X * no id of contour
X */
XCONTOUR *get_contour_ptr(no)
Xint no;
X{
X CONTOUR *cp;
X
X if((cp = CL) == (CONTOUR *)NULL)
X return((CONTOUR *)NULL);
X else {
X while(cp != (CONTOUR *)NULL) {
X if(cp->_no == no)
X return(cp);
X else
X cp = cp->_next;
X }
X return((CONTOUR *)NULL);
X }
X}
X
X
X/*
X * del_contour
X *
X * PURPOSE
X * delete contour(s) from the list with id
X * no
X */
Xdel_contour(no)
Xint no;
X{
X CONTOUR *cp, *cpp;
X CONTOUR *qp = (CONTOUR *)NULL;
X SEGMENT *sp, *spp;
X
X if((cp = CL) == (CONTOUR *)NULL)
X return;
X while(cp != (CONTOUR *)NULL) {
X if(cp->_no == no) {
X sp = cp->_s;
X do {
X spp = sp->_next;
X free(sp);
X sp = spp;
X } while(sp != cp->_s);
X cpp = cp->_next;
X free(cp);
X if(qp)
X qp->_next = cpp;
X else
X CL = cpp;
X cp = cpp;
X }
X else {
X qp = cp;
X cp = cp->_next;
X }
X }
X}
X
X
END_OF_FILE
if test 3146 -ne `wc -c <'2DClip/bio.c'`; then
echo shar: \"'2DClip/bio.c'\" unpacked with wrong size!
fi
# end of '2DClip/bio.c'
fi
if test -f '2DClip/clip.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'2DClip/clip.c'\"
else
echo shar: Extracting \"'2DClip/clip.c'\" \(2904 characters\)
sed "s/^X//" >'2DClip/clip.c' <<'END_OF_FILE'
X/*
X * file clip.c
X * contains the actual clipping routines
X */
X#include <stdio.h>
X#include "GraphicsGems.h"
X#include "line.h"
X
X/*
X * vis_vector
X *
X * PURPOSE
X * actual user interface. Draws a clipped line
X * NOTE: coordinates are given in converted LONGS!
X *
X * xf, yf from coordinates of vector to be drawn
X * xt, yt to coordinates of vector to be drawn
X */
Xvis_vector(xf, yf, xt, yt)
Xlong xf, yf, xt, yt;
X{
X SEGMENT l;
X
X if(xf == xt && yf == yt)
X return;
X l._from._x = xf;
X l._from._y = yf;
X l._to._x = xt;
X l._to._y = yt;
X/*
X * start at top of list
X */
X clip(CL, &l);
X}
X
X/*
X * clip
X *
X * PURPOSE
X *
X * p pointer to polygon
X * l pointer to line segment
X */
Xclip(p, l)
XCONTOUR *p;
XSEGMENT *l;
X{
X SEGMENT *sp, ss;
X CLIST *sol;
X POINT pt;
X boolean up, delay, inside, p_inside(), disjunct();
X int i;
X short nsol, nsmax = 2;
X
X
X/*
X * list exhausted do what you like
X * we want to plot
X */
X if(p == (CONTOUR *)NULL) {
X move((l->_from._x), (l->_from._y));
X cont((l->_to._x), (l->_to._y));
X return;
X }
X/*
X * polygon is switched off
X * take next one
X */
X if(!IS_ON(p)) {
X clip(p->_next, l);
X return;
X }
X/*
X * comparison on basis of the
X * enclosing rectangle
X */
X if(disjunct(p, l)) {
X if(!IS_NORMAL(p)) {
X clip(p->_next, l);
X }
X return;
X }
X/*
X * calculate possible intersections
X */
X sol = (CLIST *) calloc(2, sizeof(CLIST));
X sol[0]._p._x = l->_from._x;
X sol[0]._p._y = l->_from._y;
X sol[0]._type = STD;
X sol[1]._p._x = l->_to._x;
X sol[1]._p._y = l->_to._y;
X sol[1]._type = STD;
X nsol = 2;
X cross_calc(p, l, &sol, &nsol, nsmax);
X pt._x = sol[0]._p._x;
X pt._y = sol[0]._p._y;
X
X/*
X * determine status of first point
X */
X inside = p_inside(p, &pt);
X if((!inside && IS_NORMAL(p)) || (inside && !IS_NORMAL(p)))
X up = TRUE;
X else
X up = FALSE;
X delay = FALSE;
X/*
X * process list of intersections
X */
X for(i=1; i<nsol; i++) {
X if(!up) {
X ss._from._x = sol[i-1]._p._x;
X ss._from._y = sol[i-1]._p._y;
X ss._to._x = sol[i]._p._x;
X ss._to._y = sol[i]._p._y;
X clip(p->_next, &ss);
X }
X if(!delay) {
X if(sol[i]._type != DELAY)
X up = (up) ? FALSE : TRUE;
X else
X delay = TRUE;
X }
X else {
X up = (up) ? FALSE : TRUE;
X delay = FALSE;
X }
X }
X free(sol);
X}
X
X/*
X * disjunct
X *
X * PURPOSE
X * determine if the box enclosing the polygon
X * stored in p and the box enclosing the line
X * segment stored in l are disjunct.
X * Return TRUE if disjunct else FALSE
X *
X * p points to the polygon structure
X * l points to the linesegment structure
X *
X */
Xboolean disjunct(p, l)
XCONTOUR *p;
XSEGMENT *l;
X
X{
X if((max(l->_from._x, l->_to._x) < p->_minx) ||
X (min(l->_from._x, l->_to._x) > p->_maxx) ||
X (max(l->_from._y, l->_to._y) < p->_miny) ||
X (min(l->_from._y, l->_to._y) > p->_maxy) )
X return(TRUE);
X else
X return(FALSE);
X}
X
X#define DEBUG
X#ifdef DEBUG
Xmove(x, y)
Xlong x, y;
X{
X printf("(%d,%d) ->", x, y);
X}
X
Xcont(x, y)
Xlong x, y;
X{
X printf("(%d,%d)\n", x, y);
X}
X
X#endif
X
X
X
X
END_OF_FILE
if test 2904 -ne `wc -c <'2DClip/clip.c'`; then
echo shar: \"'2DClip/clip.c'\" unpacked with wrong size!
fi
# end of '2DClip/clip.c'
fi
if test -f 'BoxSphere.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'BoxSphere.c'\"
else
echo shar: Extracting \"'BoxSphere.c'\" \(3760 characters\)
sed "s/^X//" >'BoxSphere.c' <<'END_OF_FILE'
X/*
XA Simple Method for Box-Sphere Intersection Testing
Xby Jim Arvo
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
X#include "GraphicsGems.h"
X
X/*
X * This routine tests for intersection between an n-dimensional
X * axis-aligned box and an n-dimensional sphere. The mode argument
X * indicates whether the objects are to be regarded as surfaces or
X * solids. The values are:
X *
X * mode
X *
X * 0 Hollow Box, Hollow Sphere
X * 1 Hollow Box, Solid Sphere
X * 2 Solid Box, Hollow Sphere
X * 3 Solid Box, Solid Sphere
X *
X*/
Xint Box_Sphere_Intersect( n, Bmin, Bmax, C, r, mode )
Xint n; /* The dimension of the space. */
Xfloat Bmin[]; /* The minium of the box for each axis. */
Xfloat Bmax[]; /* The maximum of the box for each axis. */
Xfloat C[]; /* The sphere center in n-space. */
Xfloat r; /* The radius of the sphere. */
Xint mode; /* Selects hollow or solid. */
X {
X float a, b;
X float dmin, dmax;
X float r2 = SQR( r );
X int i, face;
X
X
X switch( mode )
X {
X case 0: /* Hollow Box - Hollow Sphere */
X dmin = 0;
X dmax = 0;
X face = FALSE;
X for( i = 0; i < n; i++ ) {
X a = SQR( C[i] - Bmin[i] );
X b = SQR( C[i] - Bmax[i] );
X dmax += MAX( a, b );
X if( C[i] < Bmin[i] ) {
X face = TRUE;
X dmin += a;
X }
X else if( C[i] > Bmax[i] ) {
X face = TRUE;
X dmin += b;
X }
X else if( MIN( a, b ) <= r2 ) face = TRUE;
X }
X if(face && ( dmin <= r2 ) && ( r2 <= dmax)) return(TRUE);
X break;
X
X case 1: /* Hollow Box - Solid Sphere */
X dmin = 0;
X face = FALSE;
X for( i = 0; i < n; i++ ) {
X if( C[i] < Bmin[i] ) {
X face = TRUE;
X dmin += SQR( C[i] - Bmin[i] );
X }
X else if( C[i] > Bmax[i] ) {
X face = TRUE;
X dmin += SQR( C[i] - Bmax[i] );
X }
X else if( C[i] - Bmin[i] <= r ) face = TRUE;
X else if( Bmax[i] - C[i] <= r ) face = TRUE;
X }
X if( face && ( dmin <= r2 ) ) return( TRUE );
X break;
X
X
X case 2: /* Solid Box - Hollow Sphere */
X dmax = 0;
X dmin = 0;
X for( i = 0; i < n; i++ ) {
X a = SQR( C[i] - Bmin[i] );
X b = SQR( C[i] - Bmax[i] );
X dmax += MAX( a, b );
X if( C[i] < Bmin[i] ) dmin += a; else
X if( C[i] > Bmax[i] ) dmin += b;
X }
X if( dmin <= r2 && r2 <= dmax ) return( TRUE );
X break;
X
X case 3: /* Solid Box - Solid Sphere */
X dmin = 0;
X for( i = 0; i < n; i++ ) {
X if( C[i] < Bmin[i] ) dmin += SQR(C[i] - Bmin[i] ); else
X if( C[i] > Bmax[i] ) dmin += SQR( C[i] - Bmax[i] );
X }
X if( dmin <= r2 ) return( TRUE );
X break;
X
X } /* end switch */
X
X return( FALSE );
X }
END_OF_FILE
if test 3760 -ne `wc -c <'BoxSphere.c'`; then
echo shar: \"'BoxSphere.c'\" unpacked with wrong size!
fi
# end of 'BoxSphere.c'
fi
if test -f 'Label.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'Label.c'\"
else
echo shar: Extracting \"'Label.c'\" \(2291 characters\)
sed "s/^X//" >'Label.c' <<'END_OF_FILE'
X/*
X * Nice Numbers for Graph Labels
X * by Paul Heckbert
X * from "Graphics Gems", Academic Press, 1990
X */
X
X/*
X * label.c: demonstrate nice graph labeling
X *
X * Paul Heckbert 2 Dec 88
X */
X
X#include <stdio.h>
X#include <math.h>
X#include "GraphicsGems.h"
X
Xdouble nicenum();
X
X/* expt(a,n)=a^n for integer n */
X
X#ifdef POW_NOT_TRUSTWORTHY
X/* if roundoff errors in pow cause problems, use this: */
X
Xdouble expt(a, n)
Xdouble a;
Xregister int n;
X{
X double x;
X
X x = 1.;
X if (n>0) for (; n>0; n--) x *= a;
X else for (; n<0; n++) x /= a;
X return x;
X}
X
X#else
X# define expt(a, n) pow(a, (double)(n))
X#endif
X
X#define NTICK 5 /* desired number of tick marks */
X
Xmain(ac, av)
Xint ac;
Xchar **av;
X{
X double min, max;
X
X if (ac!=3) {
X fprintf(stderr, "Usage: label <min> <max>\n");
X exit(1);
X }
X min = atof(av[1]);
X max = atof(av[2]);
X
X loose_label(min, max);
X}
X
X/*
X * loose_label: demonstrate loose labeling of data range from min to max.
X * (tight method is similar)
X */
X
Xloose_label(min, max)
Xdouble min, max;
X{
X char str[6], temp[20];
X int nfrac;
X double d; /* tick mark spacing */
X double graphmin, graphmax; /* graph range min and max */
X double range, x;
X
X /* we expect min!=max */
X range = nicenum(max-min, 0);
X d = nicenum(range/(NTICK-1), 1);
X graphmin = floor(min/d)*d;
X graphmax = ceil(max/d)*d;
X nfrac = MAX(-floor(log10(d)), 0); /* # of fractional digits to show */
X sprintf(str, "%%.%df", nfrac); /* simplest axis labels */
X
X printf("graphmin=%g graphmax=%g increment=%g\n", graphmin, graphmax, d);
X for (x=graphmin; x<graphmax+.5*d; x+=d) {
X sprintf(temp, str, x);
X printf("(%s)\n", temp);
X }
X}
X
X/*
X * nicenum: find a "nice" number approximately equal to x.
X * Round the number if round=1, take ceiling if round=0
X */
X
Xstatic double nicenum(x, round)
Xdouble x;
Xint round;
X{
X int exp; /* exponent of x */
X double f; /* fractional part of x */
X double nf; /* nice, rounded fraction */
X
X exp = floor(log10(x));
X f = x/expt(10., exp); /* between 1 and 10 */
X if (round)
X if (f<1.5) nf = 1.;
X else if (f<3.) nf = 2.;
X else if (f<7.) nf = 5.;
X else nf = 10.;
X else
X if (f<=1.) nf = 1.;
X else if (f<=2.) nf = 2.;
X else if (f<=5.) nf = 5.;
X else nf = 10.;
X return nf*expt(10., exp);
X}
END_OF_FILE
if test 2291 -ne `wc -c <'Label.c'`; then
echo shar: \"'Label.c'\" unpacked with wrong size!
fi
# end of 'Label.c'
fi
if test -f 'Makefile' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'Makefile'\"
else
echo shar: Extracting \"'Makefile'\" \(2894 characters\)
sed "s/^X//" >'Makefile' <<'END_OF_FILE'
X#
X# Makefile for Graphics Gems source
X#
X# Craig Kolb, 8/90
X#
X# This make file will build "gemslib.a" and a number of executables.
X# Gemslib is built solely for debugging purposes -- it is not intended
X# to be used as a library.
X#
X# Note that some of the gems need additional macros, functions, tables
X# driving routines, etc. before they will compile or run properly.
X# These include:
X#
X# AALines
X# Needs variables set in the Makefile.
X# Dissolve
X# Needs a table filled.
X# FitCurves
X# Needs a functioning DrawBezierCurve().
X# MatrixOrtho
X# Needs a number of matrix routines.
X# PolyScan
X# Needs driving routines and other additional code.
X# RayPolygon
X# Needs surrounding function structure, declaration of
X# variable types, etc.
X
X#
X# C compiler flags
X#
XCFLAGS = -g
X#
X# Location of Graphics Gems library
X#
XLIBFILE = gemslib.a
X
X#
X# Graphics Gems Vector Library
X#
XVECLIB = GGVecLib.o
X
XMFLAGS = "LIBFILE = ../$(LIBFILE)" "GENCFLAGS = $(CFLAGS)"
X
XFTPDIR = /usr/spool/uucppublic/ftp/pub/GraphicsGems/src
XZOOFILE = Gems.zoo
X
XSHELL = /bin/sh
X
XOFILES = PntOnLine.o ViewTrans.o AAPolyScan.o Albers.o \
X Interleave.o BoundSphere.o BoxSphere.o CircleRect.o \
X ConcaveScan.o Roots3And4.o Dissolve.o DigitalLine.o \
X FastJitter.o FixedTrig.o HSLtoRGB.o HypotApprox.o LineEdge.o \
X MatrixInvert.o MatrixPost.o Median.o PixelInteger.o \
X TriPoints.o Quaternions.o RGBTo4Bits.o RayBox.o \
X SeedFill.o SquareRoot.o DoubleLine.o TransBox.o
X
XDIRS = 2DClip PolyScan Sturm
X
XALL = Hash3D FitCurves Forms NearestPoint Label \
X OrderDither BinRec $(LIBFILE)
X
Xall: $(ALL)
X @for d in $(DIRS) ; do \
X (cd $$d ; $(MAKE) $(MFLAGS)) ;\
X done
X
X$(LIBFILE): $(OFILES) $(VECLIB)
X ar rcs $(LIBFILE) $(OFILES) $(VECLIB)
X
XHash3D: Hash3D.o
X $(CC) $(CFLAGS) -o $@ Hash3D.o
X
XFitCurves: FitCurves.o $(VECLIB)
X $(CC) $(CFLAGS) -o $@ FitCurves.o $(VECLIB) -lm
X
XForms: Forms.o
X $(CC) $(CFLAGS) -o $@ Forms.o
X
XNearestPoint: NearestPoint.o $(VECLIB)
X $(CC) $(CFLAGS) -o $@ NearestPoint.o $(VECLIB) -lm
X
XLabel: Label.o
X $(CC) $(CFLAGS) -o $@ Label.o -lm
X
XOrderDither: OrderDither.o
X $(CC) $(CFLAGS) -o $@ OrderDither.o
X
XBinRec: BinRec.o
X $(CC) $(CFLAGS) -o $@ BinRec.o
X
Xftpreadme:
X /bin/rm -f $(FTPDIR)/README
X sed -e "s/__DATE__/`date`/g" < README.ftp > $(FTPDIR)/README
X
Xreadme:
X /bin/rm -rf README
X sed -e "s/__DATE__/`date`/g" < README.dist > README
X
X
Xclean:
X @for d in $(DIRS) ; do \
X (cd $$d ; $(MAKE) $(MFLAGS) clean) ;\
X done
X /bin/rm -f $(OFILES) $(VECLIB)
X /bin/rm -f Hash3D.o FitCurves.o \
X Forms.o NearestPoint.o Label.o OrderDither.o BinRec.o \
X Hash3D FitCurves Forms NearestPoint Label OrderDither BinRec \
X bugs a.out core Part??.Z Part?? $(ZOOFILE)
X
Xkit: readme
X /bin/rm -f Part??.Z Part??
X (makekit -iPACKING_LIST -oMANIFEST)
X
Xzoo: readme
X /bin/rm -f Gems.zoo
X cut -d" " -f2 PACKING_LIST | zoo aI $(ZOOFILE)
X
Xftp: kit zoo ftpreadme
X compress Part??
X /bin/cp Part??.Z $(ZOOFILE) $(FTPDIR)
X
X$(ALL): GraphicsGems.h
END_OF_FILE
if test 2894 -ne `wc -c <'Makefile'`; then
echo shar: \"'Makefile'\" unpacked with wrong size!
fi
# end of 'Makefile'
fi
if test -f 'MatrixPost.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'MatrixPost.c'\"
else
echo shar: Extracting \"'MatrixPost.c'\" \(3390 characters\)
sed "s/^X//" >'MatrixPost.c' <<'END_OF_FILE'
X/*
XEfficient Post-Concatenation of Transformation Matrics
Xby Joseph M. Cychosz
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
X#include <math.h>
X#include "GraphicsGems.h"
X
X
X/* M4xform.c - Basic 4x4 transformation package.
X *
X * Description:
X * M4xform.c contains a collection of routines used to perform
X * direct post-concatenated transformation operations.
X *
X *
X * Contents:
X * M4RotateX Post-concatenate a x-axis rotation.
X * M4RotateY Post-concatenate a y-axis rotation.
X * M4RotateZ Post-concatenate a z-axis rotation.
X * M4Scale Post-concatenate a scaling.
X * M4Translate Post-concatenate a translation.
X * M4ZPerspective Post-concatenate a z-axis perspective
X * transformation.
X *
X * Externals:
X * cos, sin.
X */
X
X
X
X/* M4RotateX - Post-concatenate a x-axis rotation matrix. */
X
XMatrix4 *M4RotateX (m,a)
X Matrix4 *m; /* Current transformation matrix*/
X double a; /* Rotation angle */
X{
X double c, s;
X double t;
X int i;
X
X
X c = cos (a); s = sin (a);
X
X for (i = 0 ; i < 4 ; i++) {
X t = m->element[i][1];
X m->element[i][1] = t*c - m->element[i][2]*s;
X m->element[i][2] = t*s + m->element[i][2]*c;
X }
X return (m);
X}
X
X
X/* M4RotateY - Post-concatenate a y-axis rotation matrix. */
X
XMatrix4 *M4RotateY (m,a)
X Matrix4 *m; /* Current transformation matrix*/
X double a; /* Rotation angle */
X{
X double c, s;
X double t;
X int i;
X
X c = cos (a); s = sin (a);
X
X for (i = 0 ; i < 4 ; i++) {
X t = m->element[i][0];
X m->element[i][0] = t*c + m->element[i][2]*s;
X m->element[i][2] = m->element[i][2]*c - t*s;
X }
X return (m);
X}
X
X
X/* M4RotateZ - Post-concatenate a z-axis rotation matrix. */
X
XMatrix4 *M4RotateZ (m,a)
X Matrix4 *m; /* Current transformation matrix*/
X double a; /* Rotation angle */
X{
X double c, s;
X double t;
X int i;
X
X c = cos (a); s = sin (a);
X
X for (i = 0 ; i < 4 ; i++) {
X t = m->element[i][0];
X m->element[i][0] = t*c - m->element[i][1]*s;
X m->element[i][1] = t*s + m->element[i][1]*c;
X }
X return (m);
X}
X
X
X
X/* M4Scale - Post-concatenate a scaling. */
X
XMatrix4 *M4Scale (m,sx,sy,sz)
X Matrix4 *m; /* Current transformation matrix */
X double sx, sy, sz; /* Scale factors about x,y,z */
X{
X int i;
X
X for (i = 0 ; i < 4 ; i++) {
X m->element[i][0] *= sx;
X m->element[i][1] *= sy;
X m->element[i][2] *= sz;
X }
X return (m);
X}
X
X
X/* M4Translate - Post-concatenate a translation. */
X
XMatrix4 *M4Translate (m,tx,ty,tz)
X Matrix4 *m; /* Current transformation matrix */
X double tx, ty, tz; /* Translation distance */
X{
X int i;
X
X for (i = 0 ; i < 4 ; i++) {
X m->element[i][0] += m->element[i][3]*tx;
X m->element[i][1] += m->element[i][3]*ty;
X m->element[i][2] += m->element[i][3]*tz;
X }
X return (m);
X}
X
X
X /* M4ZPerspective Post-concatenate a perspective */ /*transformation. */
X /* */
X /* Perspective is along the z-axis with the eye at +z. */
X
X
XMatrix4 *M4ZPerspective (m,d)
X Matrix4 *m; /* Current transformation matrix */
X double d; /* Perspective distance */
X{
X int i;
X double f = 1. / d;
X
X for (i = 0 ; i < 4 ; i++) {
X m->element[i][3] += m->element[i][2]*f;
X m->element[i][2] = 0.;
X }
X return (m);
X}
X
X
X
X
END_OF_FILE
if test 3390 -ne `wc -c <'MatrixPost.c'`; then
echo shar: \"'MatrixPost.c'\" unpacked with wrong size!
fi
# end of 'MatrixPost.c'
fi
if test -f 'Median.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'Median.c'\"
else
echo shar: Extracting \"'Median.c'\" \(2941 characters\)
sed "s/^X//" >'Median.c' <<'END_OF_FILE'
X/*
XMedian Finding on a 3-by-3 Grid
Xby Alan Paeth
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
X#define s2(a,b) {register int t; if ((t=b-a)<0) {a+=t; b-=t;}}
X#define mn3(a,b,c) s2(a,b); s2(a,c);
X#define mx3(a,b,c) s2(b,c); s2(a,c);
X#define mnmx3(a,b,c) mx3(a,b,c); s2(a,b);
X#define mnmx4(a,b,c,d) s2(a,b); s2(c,d); s2(a,c); s2(b,d);
X#define mnmx5(a,b,c,d,e) s2(a,b); s2(c,d); mn3(a,c,e); mx3(b,d,e);
X#define mnmx6(a,b,c,d,e,f) s2(a,d); s2(b,e); s2(c,f);\
X mn3(a,b,c); mx3(d,e,f);
Xmed3x3(b1, b2, b3)
X int *b1, *b2, *b3;
X/*
X * Find median on a 3x3 input box of integers.
X * b1, b2, b3 are pointers to the left-hand edge of three
X * parallel scan-lines to form a 3x3 spatial median.
X * Rewriting b2 and b3 as b1 yields code which forms median
X * on input presented as a linear array of nine elements.
X */
X {
X register int r1, r2, r3, r4, r5, r6;
X r1 = *b1++; r2 = *b1++; r3 = *b1++;
X r4 = *b2++; r5 = *b2++; r6 = *b2++;
X mnmx6(r1, r2, r3, r4, r5, r6);
X r1 = *b3++;
X mnmx5(r1, r2, r3, r4, r5);
X r1 = *b3++;
X mnmx4(r1, r2, r3, r4);
X r1 = *b3++;
X mnmx3(r1, r2, r3);
X return(r2);
X }
X
X
X/* t2(i,j) transposes elements in A[] such that A[i] <= A[j] */
X
X#define t2(i, j) s2(A[i-1], A[j-1])
X
X
Xint median25(A)
X int A[25];
X {
X/*
X * median25(A) partitions the array A[0..24] such that element
X * A[12] is the median and subarrays A[0..11] and A[13..24]
X * are partitions containing elements of smaller and larger
X * value (rank), respectively.
X *
X * The exchange table lists element indices on the range 1..25,
X * this accounts for the "-1" offsets in the macro t2 and in
X * the final return value used to adjust subscripts to C-code
X * conventions (array indices begin at zero).
X */
X t2( 1, 2); t2( 4, 5); t2( 3, 5); t2( 3, 4); t2( 7, 8);
X t2( 6, 8); t2( 6, 7); t2(10,11); t2( 9,11); t2( 9,10);
X t2(13,14); t2(12,14); t2(12,13); t2(16,17); t2(15,17);
X t2(15,16); t2(19,20); t2(18,20); t2(18,19); t2(22,23);
X t2(21,23); t2(21,22); t2(24,25); t2( 3, 6); t2( 4, 7);
X t2( 1, 7); t2( 1, 4); t2( 5, 8); t2( 2, 8); t2( 2, 5);
X t2(12,15); t2( 9,15); t2( 9,12); t2(13,16); t2(10,16);
X t2(10,13); t2(14,17); t2(11,17); t2(11,14); t2(21,24);
X t2(18,24); t2(18,21); t2(22,25); t2(19,25); t2(19,22);
X t2(20,23); t2( 9,18); t2(10,19); t2( 1,19); t2( 1,10);
X t2(11,20); t2( 2,20); t2( 2,11); t2(12,21); t2( 3,21);
X t2( 3,12); t2(13,22); t2( 4,22); t2( 4,13); t2(14,23);
X t2( 5,23); t2( 5,14); t2(15,24); t2( 6,24); t2( 6,15);
X t2(16,25); t2( 7,25); t2( 7,16); t2( 8,17); t2( 8,20);
X t2(14,22); t2(16,24); t2( 8,14); t2( 8,16); t2( 2,10);
X t2( 4,12); t2( 6,18); t2(12,18); t2(10,18); t2( 5,11);
X t2( 7,13); t2( 8,15); t2( 5, 7); t2( 5, 8); t2(13,15);
X t2(11,15); t2( 7, 8); t2(11,13); t2( 7,11); t2( 7,18);
X t2(13,18); t2( 8,18); t2( 8,11); t2(13,19); t2( 8,13);
X t2(11,19); t2(13,21); t2(11,21); t2(11,13);
X return(A[13-1]);
X }
END_OF_FILE
if test 2941 -ne `wc -c <'Median.c'`; then
echo shar: \"'Median.c'\" unpacked with wrong size!
fi
# end of 'Median.c'
fi
if test -f 'OrderDither.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'OrderDither.c'\"
else
echo shar: Extracting \"'OrderDither.c'\" \(2483 characters\)
sed "s/^X//" >'OrderDither.c' <<'END_OF_FILE'
X/*
XOrdered Dithering
Xby Stephen Hawley
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
X/* Program to generate dithering matrices.
X * written by Jim Blandy, Oberlin College, jimb at occs.oberlin.edu
X * Gifted to, documented and revised by Stephen Hawley,
X * sdh at flash.bellcore.com
X *
X * Generates a dithering matrix from the command line arguments.
X * The first argument, size, determines the dimensions of the
X * matrix: 2^size by 2^size
X * The optional range argument is the range of values to be
X * dithered over. By default, it is (2^size)^2, or simply the
X * total number of elements in the matrix.
X * The final output is suitable for inclusion in a C program.
X * A typical dithering function is something like this:
X * extern int dm[], size;
X *
X * int
X * dither(x,y, level)
X * register int x,y, level;
X * {
X * return(level > dm[(x % size) + size * (y % size)];
X * }
X */
X
Xmain(argc, argv)
Xint argc;
Xchar **argv;
X{
X register int size, range;
X
X if (argc >= 2) size = atoi(argv[1]);
X else size = 2;
X
X if (argc == 3) range = atoi(argv[2]);
X else range = (1 << size) * (1 << size);
X
X printdither (size, range);
X}
X
X
Xprintdither (size, range)
Xregister int size, range;
X{
X register int l = (1 << size), i;
X /*
X * print a dithering matrix.
X * l is the length on a side.
X */
X range = range / (l * l);
X puts("int dm[] = {");
X for (i=0; i < l*l; i++) {
X if (i % l == 0) /* tab in 4 spaces per row */
X printf(" ");
X /* print the dither value for this location
X * scaled to the given range
X */
X printf("%4d", range * dithervalue(i / l, i % l, size));
X
X /* commas after all but the last */
X if (i + 1 < l * l)
X putchar(',');
X /* newline at the end of the row */
X if ((i + 1) % l == 0)
X putchar('\n');
X }
X puts("\n}; ");
X}
Xdithervalue(x, y, size)
Xregister int x, y, size;
X{
X register int d;
X /*
X * calculate the dither value at a particular
X * (x, y) over the size of the matrix.
X */
X d=0;
X while (size-->0) {
X /* Think of d as the density. At every iteration,
X * d is shifted left one and a new bit is put in the
X * low bit based on x and y. If x is odd and y is even,
X * or x is even and y is odd, a bit is put in. This
X * generates the checkerboard seen in dithering.
X * This quantity is shifted left again and the low bit of
X * y is added in.
X * This whole thing interleaves a checkerboard bit pattern
X * and y's bits, which is the value you want.
X */
X d = (d <<1 | (x&1 ^ y&1))<<1 | y&1;
X x >>= 1;
X y >>= 1;
X }
X return(d);
X}
X
END_OF_FILE
if test 2483 -ne `wc -c <'OrderDither.c'`; then
echo shar: \"'OrderDither.c'\" unpacked with wrong size!
fi
# end of 'OrderDither.c'
fi
if test -f 'PntOnLine.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'PntOnLine.c'\"
else
echo shar: Extracting \"'PntOnLine.c'\" \(2024 characters\)
sed "s/^X//" >'PntOnLine.c' <<'END_OF_FILE'
X/*
XA Fast 2D Point-On-Line Test
Xby Alan Paeth
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
X#include "GraphicsGems.h"
X
Xint PntOnLine(px,py,qx,qy,tx,ty)
X int px, py, qx, qy, tx, ty;
X {
X/*
X * given a line through P:(px,py) Q:(qx,qy) and T:(tx,ty)
X * return 0 if T is not on the line through <--P--Q-->
X * 1 if T is on the open ray ending at P: <--P
X * 2 if T is on the closed interior along: P--Q
X * 3 if T is on the open ray beginning at Q: Q-->
X *
X * Example: consider the line P = (3,2), Q = (17,7). A plot
X * of the test points T(x,y) (with 0 mapped onto '.') yields:
X *
X * 8| . . . . . . . . . . . . . . . . . 3 3
X * Y 7| . . . . . . . . . . . . . . 2 2 Q 3 3 Q = 2
X * 6| . . . . . . . . . . . 2 2 2 2 2 . . .
X * a 5| . . . . . . . . 2 2 2 2 2 2 . . . . .
X * x 4| . . . . . 2 2 2 2 2 2 . . . . . . . .
X * i 3| . . . 2 2 2 2 2 . . . . . . . . . . .
X * s 2| 1 1 P 2 2 . . . . . . . . . . . . . . P = 2
X * 1| 1 1 . . . . . . . . . . . . . . . . .
X * +--------------------------------------
X * 1 2 3 4 5 X-axis 10 15 19
X *
X * Point-Line distance is normalized with the Infinity Norm
X * avoiding square-root code and tightening the test vs the
X * Manhattan Norm. All math is done on the field of integers.
X * The latter replaces the initial ">= MAX(...)" test with
X * "> (ABS(qx-px) + ABS(qy-py))" loosening both inequality
X * and norm, yielding a broader target line for selection.
X * The tightest test is employed here for best discrimination
X * in merging collinear (to grid coordinates) vertex chains
X * into a larger, spanning vectors within the Lemming editor.
X */
X
X
X if ( ABS((qy-py)*(tx-px)-(ty-py)*(qx-px)) >=
X (MAX(ABS(qx-px), ABS(qy-py)))) return(0);
X if (((qx<px)&&(px<tx)) || ((qy<py)&&(py<ty))) return(1);
X if (((tx<px)&&(px<qx)) || ((ty<py)&&(py<qy))) return(1);
X if (((px<qx)&&(qx<tx)) || ((py<qy)&&(qy<ty))) return(3);
X if (((tx<qx)&&(qx<px)) || ((ty<qy)&&(qy<py))) return(3);
X return(2);
X }
END_OF_FILE
if test 2024 -ne `wc -c <'PntOnLine.c'`; then
echo shar: \"'PntOnLine.c'\" unpacked with wrong size!
fi
# end of 'PntOnLine.c'
fi
if test -f 'PolyScan/fancytest.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'PolyScan/fancytest.c'\"
else
echo shar: Extracting \"'PolyScan/fancytest.c'\" \(3298 characters\)
sed "s/^X//" >'PolyScan/fancytest.c' <<'END_OF_FILE'
X/*
X * fancytest.c: subroutine illustrating the use of poly_clip and poly_scan
X * for Phong-shading and texture mapping.
X *
X * Note: lines enclosed in angle brackets '<', '>' should be replaced
X * with the code described.
X * Makes calls to hypothetical packages "shade", "image", "texture", "zbuffer".
X *
X * Paul Heckbert Dec 1989
X */
X
X#include <stdio.h>
X#include <math.h>
X#include "poly.h"
X
X#define XMAX 1280 /* hypothetical image width */
X#define YMAX 1024 /* hypothetical image height */
X#define LIGHT_INTEN 255 /* light source intensity */
X
Xvoid pixelproc();
X
Xfancytest()
X{
X int i;
X double WS[4][4]; /* world space to screen space transform */
X Poly p; /* a polygon */
X Poly_vert *v;
X
X static Poly_box box = {0, XMAX, 0, YMAX, -32678, 32767.99};
X /* 3-D screen clipping box, continuous coordinates */
X
X static Window win = {0, 0, XMAX-1, YMAX-1};
X /* 2-D screen clipping window, discrete coordinates */
X
X <initialize world space position (x,y,z), normal (nx,ny,nz), and texture
X position (u,v) at each vertex of p; set p.n>
X <set WS to world-to-screen transform>
X
X /* transform vertices from world space to homogeneous screen space */
X for (i=0; i<p.n; i++) {
X v = &p.vert[i];
X mx4_transform(v->x, v->y, v->z, 1., WS, &v->sx, &v->sy, &v->sz, &v->sw);
X }
X
X /* interpolate sx, sy, sz, sw, nx, ny, nz, u, v in poly_clip */
X p.mask = POLY_MASK(sx) | POLY_MASK(sy) | POLY_MASK(sz) | POLY_MASK(sw) |
X POLY_MASK(nx) | POLY_MASK(ny) | POLY_MASK(nz) |
X POLY_MASK(u) | POLY_MASK(v);
X
X poly_print("before clipping", &p);
X if (poly_clip_to_box(&p, &box) == POLY_CLIP_OUT) /* clip polygon */
X return; /* quit if off-screen */
X
X /* do homogeneous division of screen position and texture position */
X for (i=0; i<p.n; i++) {
X v = &p.vert[i];
X v->sx /= v->sw;
X v->sy /= v->sw;
X v->sz /= v->sw;
X v->u /= v->sw;
X v->v /= v->sw;
X v->q = 1./v->sw;
X }
X /*
X * previously we ignored q (assumed q=1), henceforth ignore sw (assume sw=1)
X * Interpolate sx, sy, sz, nx, ny, nz, u, v, q in poly_scan
X */
X p.mask &= ~POLY_MASK(sw);
X p.mask |= POLY_MASK(q);
X
X poly_print("scan converting the polygon", &p);
X poly_scan(&p, &win, pixelproc); /* scan convert! */
X}
X
Xstatic void pixelproc(x, y, pt) /* called at each pixel by poly_scan */
Xint x, y;
XPoly_vert *pt;
X{
X int sz, u, v, inten;
X double len, nor[3], diff, spec;
X
X sz = pt->sz;
X if (sz < zbuffer_read(x, y)) {
X len = sqrt(pt->nx*pt->nx + pt->ny*pt->ny + pt->nz*pt->nz);
X nor[0] = pt->nx/len; /* unitize the normal vector */
X nor[1] = pt->ny/len;
X nor[2] = pt->nz/len;
X shade(nor, &diff, &spec); /* compute specular and diffuse coeffs*/
X u = pt->u/pt->q; /* do homogeneous div. of texture pos */
X v = pt->v/pt->q;
X inten = texture_read(u, v)*diff + LIGHT_INTEN*spec;
X image_write(x, y, inten);
X zbuffer_write(x, y, sz);
X }
X}
X
X/* mx4_transform: transform 4-vector p by matrix m yielding q: q = p*m */
X
Xmx4_transform(px, py, pz, pw, m, qxp, qyp, qzp, qwp)
Xdouble px, py, pz, pw, m[4][4], *qxp, *qyp, *qzp, *qwp;
X{
X *qxp = px*m[0][0] + py*m[1][0] + pz*m[2][0] + pw*m[3][0];
X *qyp = px*m[0][1] + py*m[1][1] + pz*m[2][1] + pw*m[3][1];
X *qzp = px*m[0][2] + py*m[1][2] + pz*m[2][2] + pw*m[3][2];
X *qwp = px*m[0][3] + py*m[1][3] + pz*m[2][3] + pw*m[3][3];
X}
END_OF_FILE
if test 3298 -ne `wc -c <'PolyScan/fancytest.c'`; then
echo shar: \"'PolyScan/fancytest.c'\" unpacked with wrong size!
fi
# end of 'PolyScan/fancytest.c'
fi
if test -f 'PolyScan/poly.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'PolyScan/poly.c'\"
else
echo shar: Extracting \"'PolyScan/poly.c'\" \(2310 characters\)
sed "s/^X//" >'PolyScan/poly.c' <<'END_OF_FILE'
X/*
X * poly.c: simple utilities for polygon data structures
X */
X
X#include "poly.h"
X
XPoly_vert *poly_dummy; /* used superficially by POLY_MASK macro */
X
X/*
X * poly_print: print Poly p to stdout, prefixed by the label str
X */
X
Xvoid poly_print(str, p)
Xchar *str;
XPoly *p;
X{
X int i;
X
X printf("%s: %d sides\n", str, p->n);
X poly_vert_label(" ", p->mask);
X for (i=0; i<p->n; i++) {
X printf(" v[%d] ", i);
X poly_vert_print("", &p->vert[i], p->mask);
X }
X}
X
Xvoid poly_vert_label(str, mask)
Xchar *str;
Xint mask;
X{
X printf("%s", str);
X if (mask&POLY_MASK(sx)) printf(" sx ");
X if (mask&POLY_MASK(sy)) printf(" sy ");
X if (mask&POLY_MASK(sz)) printf(" sz ");
X if (mask&POLY_MASK(sw)) printf(" sw ");
X if (mask&POLY_MASK(x)) printf(" x ");
X if (mask&POLY_MASK(y)) printf(" y ");
X if (mask&POLY_MASK(z)) printf(" z ");
X if (mask&POLY_MASK(u)) printf(" u ");
X if (mask&POLY_MASK(v)) printf(" v ");
X if (mask&POLY_MASK(q)) printf(" q ");
X if (mask&POLY_MASK(r)) printf(" r ");
X if (mask&POLY_MASK(g)) printf(" g ");
X if (mask&POLY_MASK(b)) printf(" b ");
X if (mask&POLY_MASK(nx)) printf(" nx ");
X if (mask&POLY_MASK(ny)) printf(" ny ");
X if (mask&POLY_MASK(nz)) printf(" nz ");
X printf("\n");
X}
X
Xvoid poly_vert_print(str, v, mask)
Xchar *str;
XPoly_vert *v;
Xint mask;
X{
X printf("%s", str);
X if (mask&POLY_MASK(sx)) printf(" %6.1f", v->sx);
X if (mask&POLY_MASK(sy)) printf(" %6.1f", v->sy);
X if (mask&POLY_MASK(sz)) printf(" %6.2f", v->sz);
X if (mask&POLY_MASK(sw)) printf(" %6.2f", v->sw);
X if (mask&POLY_MASK(x)) printf(" %6.2f", v->x);
X if (mask&POLY_MASK(y)) printf(" %6.2f", v->y);
X if (mask&POLY_MASK(z)) printf(" %6.2f", v->z);
X if (mask&POLY_MASK(u)) printf(" %6.2f", v->u);
X if (mask&POLY_MASK(v)) printf(" %6.2f", v->v);
X if (mask&POLY_MASK(q)) printf(" %6.2f", v->q);
X if (mask&POLY_MASK(r)) printf(" %6.4f", v->r);
X if (mask&POLY_MASK(g)) printf(" %6.4f", v->g);
X if (mask&POLY_MASK(b)) printf(" %6.4f", v->b);
X if (mask&POLY_MASK(nx)) printf(" %6.3f", v->nx);
X if (mask&POLY_MASK(ny)) printf(" %6.3f", v->ny);
X if (mask&POLY_MASK(nz)) printf(" %6.3f", v->nz);
X printf("\n");
X}
END_OF_FILE
if test 2310 -ne `wc -c <'PolyScan/poly.c'`; then
echo shar: \"'PolyScan/poly.c'\" unpacked with wrong size!
fi
# end of 'PolyScan/poly.c'
fi
if test -f 'PolyScan/poly.h' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'PolyScan/poly.h'\"
else
echo shar: Extracting \"'PolyScan/poly.h'\" \(1914 characters\)
sed "s/^X//" >'PolyScan/poly.h' <<'END_OF_FILE'
X/* poly.h: definitions for polygon package */
X
X#ifndef POLY_HDR
X#define POLY_HDR
X
X#define POLY_NMAX 8 /* max #sides to a polygon; change if needed */
X
Xtypedef struct { /* A POLYGON VERTEX */
X double sx, sy, sz, sw; /* screen space position (sometimes homo.) */
X double x, y, z; /* world space position */
X double u, v, q; /* texture position (sometimes homogeneous) */
X double r, g, b; /* (red,green,blue) color */
X double nx, ny, nz; /* world space normal vector */
X} Poly_vert;
X/* update poly.c if you change this structure */
X
Xtypedef struct { /* A POLYGON */
X int n; /* number of sides */
X int mask; /* interpolation mask for vertex elems */
X Poly_vert vert[POLY_NMAX]; /* vertices */
X} Poly;
X/*
X * mask is an interpolation mask whose kth bit indicates whether the kth
X * double in a Poly_vert is relevant.
X * For example, if the valid attributes are sx, sy, and sz, then set
X * mask = POLY_MASK(sx) | POLY_MASK(sy) | POLY_MASK(sz);
X */
X
Xtypedef struct { /* A BOX (TYPICALLY IN SCREEN SPACE) */
X double x0, x1; /* left and right */
X double y0, y1; /* top and bottom */
X double z0, z1; /* near and far */
X} Poly_box;
X
Xtypedef struct { /* WINDOW: A DISCRETE 2-D RECTANGLE */
X int x0, y0; /* xmin and ymin */
X int x1, y1; /* xmax and ymax (inclusive) */
X} Window;
X
X#define POLY_MASK(elem) (1 << (&poly_dummy->elem - (double *)poly_dummy))
X
X#define POLY_CLIP_OUT 0 /* polygon entirely outside box */
X#define POLY_CLIP_PARTIAL 1 /* polygon partially inside */
X#define POLY_CLIP_IN 2 /* polygon entirely inside box */
X
Xextern Poly_vert *poly_dummy; /* used superficially by POLY_MASK macro */
X
Xvoid poly_print(/* str, p */);
Xvoid poly_vert_label(/* str, mask */);
Xvoid poly_vert_print(/* str, v, mask */);
Xint poly_clip_to_box(/* p1, box */);
Xvoid poly_clip_to_halfspace(/* p, q, index, sign, k, name */);
Xvoid poly_scan(/* p, win, pixelproc */);
X
X#endif
END_OF_FILE
if test 1914 -ne `wc -c <'PolyScan/poly.h'`; then
echo shar: \"'PolyScan/poly.h'\" unpacked with wrong size!
fi
# end of 'PolyScan/poly.h'
fi
if test -f 'PolyScan/scantest.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'PolyScan/scantest.c'\"
else
echo shar: Extracting \"'PolyScan/scantest.c'\" \(1854 characters\)
sed "s/^X//" >'PolyScan/scantest.c' <<'END_OF_FILE'
X/*
X * scantest.c: use poly_scan() for Gouraud shading and z-buffer demo.
X * Given the screen space X, Y, and Z of N-gon on command line,
X * print out all pixels during scan conversion.
X * This code could easily be modified to actually read and write pixels.
X *
X * Paul Heckbert Dec 1989
X */
X
X#include <stdio.h>
X#include <math.h>
X#include "poly.h"
X
X#define XMAX 1280 /* hypothetical image width */
X#define YMAX 1024 /* hypothetical image height */
X
X#define FRANDOM() ((rand()&32767)/32767.) /* random number between 0 and 1 */
X
Xstatic void pixelproc();
X
Xmain(ac, av)
Xint ac;
Xchar **av;
X{
X int i;
X Poly p;
X static Window win = {0, 0, XMAX-1, YMAX-1}; /* screen clipping window */
X
X if (ac<2 || ac != 2+3*(p.n = atoi(av[1]))) {
X fprintf(stderr, "Usage: scantest N X1 Y1 Z1 X2 Y2 Z2 ... XN YN ZN\n");
X exit(1);
X }
X for (i=0; i<p.n; i++) {
X p.vert[i].sx = atof(av[2+3*i]); /* set screen space x,y,z */
X p.vert[i].sy = atof(av[3+3*i]);
X p.vert[i].sz = atof(av[4+3*i]);
X p.vert[i].r = FRANDOM(); /* random vertex colors, for kicks */
X p.vert[i].g = FRANDOM();
X p.vert[i].b = FRANDOM();
X }
X /* interpolate sx, sy, sz, r, g, and b in poly_scan */
X p.mask = POLY_MASK(sx) | POLY_MASK(sy) | POLY_MASK(sz) |
X POLY_MASK(r) | POLY_MASK(g) | POLY_MASK(b);
X
X poly_print("scan converting the polygon", &p);
X
X poly_scan(&p, &win, pixelproc); /* scan convert! */
X}
X
Xstatic void pixelproc(x, y, point) /* called at each pixel by poly_scan */
Xint x, y;
XPoly_vert *point;
X{
X printf("pixel (%d,%d) screenz=%g rgb=(%g,%g,%g)\n",
X x, y, point->sz, point->r, point->g, point->b);
X
X /*
X * in real graphics program you could read and write pixels, e.g.:
X *
X * if (point->sz < zbuffer_read(x, y)) {
X * image_write_rgb(x, y, point->r, point->g, point->b);
X * zbuffer_write(x, y, point->sz);
X * }
X */
X}
END_OF_FILE
if test 1854 -ne `wc -c <'PolyScan/scantest.c'`; then
echo shar: \"'PolyScan/scantest.c'\" unpacked with wrong size!
fi
# end of 'PolyScan/scantest.c'
fi
if test -f 'Quaternions.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'Quaternions.c'\"
else
echo shar: Extracting \"'Quaternions.c'\" \(2216 characters\)
sed "s/^X//" >'Quaternions.c' <<'END_OF_FILE'
X/*
XUsing Quaternions for Coding 3D Transformations
XPatrick-Gilles Maillot
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
Xextern double P[3], Q[4], M[4][4];
X
Xset_obs_position(x,y,z)
Xfloat x, y, z;
X{
Xint i;
X
X/*
X * Set the values of the eye's position.
X * The position here represents the position of the orthonormal base
X * in respect to the observer.
X */
X P[0] = -x;
X P[1] = -y;
X P[2] = -z;
X/*
X * Set the visualization to be in the decreasing x axis
X */
X Q[0] = 1.;
X for (i = 1; i < 4; i++) Q[i] = 0.;
X}
X
Xtranslate_quaternion(x,i,w)
Xfloat x;
Xint i, w;
X{
Xint j, k;
Xfloat A, B, D, E, F;
X
X if (w < 0) {
X/*
X * The observer moves in respect to the scene.
X */
X P[i - 1] -= x;
X } else {
X
X/*
X * The scene moves in respect to the observer.
X * Compute the successor axis of i [1,2,3];
X * and then the successor axis of j [1,2,3];
X */
X if ((j = i + 1) > 3) j = 1;
X if ((k = j + 1) > 3) k = 1;
X A = Q[j]; B = Q[k]; F = Q[0]; E = Q[i];
X P[i - 1] += x * (E * E + F * F - A * A - B * B);
X D = x + x;
X P[j - 1] += D * (E * A + F * B);
X P[k - 1] += D * (E * B + F * A);
X }
X}
X
Xrotate_quaternion(x,y,i,w)
Xfloat x, y;
Xint i, w;
X{
Xint j, k;
Xfloat E, F, R1;
X/*
X * Compute the successor axis of i [1,2,3] and j [1,2,3];
X */
X if ((j = i + 1) > 3) j = 1;
X if ((k = j + 1) > 3) k = 1;
X E = Q[i];
X Q[i] = E * x + w * y * Q[0];
X Q[0] = Q[0] * x - w * y * E;
X E = Q[j];
X Q[j] = E * x + y * Q[k];
X Q[k] = Q[k] * x - y * E;
X if (w < 0) {
X/* Compute a new position if the observer moves in respect to the scene. */
X j -= 1; k -= 1;
X R1 = x * x - y * y;
X F = 2. * x * y;
X E = P[j];
X P[j] = E * R1 + F * P[k];
X P[k] = P[k] * R1 - F * E;
X }
X}
X
X
XEvaluate_matrix()
X{
Xfloat e, f, r[4];
Xint i, j, k;
X/*
X * We will need some square values!
X */
X for (i = 0; i < 4; i++) r[i] = Q[i] * Q[i];
X/*
X * Compute each element of the matrix.
X * j is the successor of i (in 1,2,3), while k is the successor of j.
X */
X for (i = 1; i < 4; i++) {
X if ((j = i + 1) > 3) j = 1;
X if ((k = j + 1) > 3) k = 1;
X e = 2. * Q[i] * Q[j];
X f = 2. * Q[k] * Q[0];
X M[j][i] = e - f;
X M[i][j] = e + f;
X M[i][i] = r[i] + r[0] - r[j] - r[k];
X M[0][i] = P[i - 1];
X M[i][0] = 0.;
X }
X M[0][0] = 1.;
X}
X
X
X
X
END_OF_FILE
if test 2216 -ne `wc -c <'Quaternions.c'`; then
echo shar: \"'Quaternions.c'\" unpacked with wrong size!
fi
# end of 'Quaternions.c'
fi
if test -f 'SeedFill.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'SeedFill.c'\"
else
echo shar: Extracting \"'SeedFill.c'\" \(2371 characters\)
sed "s/^X//" >'SeedFill.c' <<'END_OF_FILE'
X/*
X * A Seed Fill Algorithm
X * by Paul Heckbert
X * from "Graphics Gems", Academic Press, 1990
X */
X
X/*
X * fill.c : simple seed fill program
X * Calls pixelread() to read pixels, pixelwrite() to write pixels.
X *
X * Paul Heckbert 13 Sept 1982, 28 Jan 1987
X */
X
Xtypedef struct { /* window: a discrete 2-D rectangle */
X int x0, y0; /* xmin and ymin */
X int x1, y1; /* xmax and ymax (inclusive) */
X} Window;
X
Xtypedef int Pixel; /* 1-channel frame buffer assumed */
X
XPixel pixelread();
X
Xtypedef struct {short y, xl, xr, dy;} Segment;
X/*
X * Filled horizontal segment of scanline y for xl<=x<=xr.
X * Parent segment was on line y-dy. dy=1 or -1
X */
X
X#define MAX 10000 /* max depth of stack */
X
X#define PUSH(Y, XL, XR, DY) /* push new segment on stack */ \
X if (sp<stack+MAX && Y+(DY)>=win->y0 && Y+(DY)<=win->y1) \
X {sp->y = Y; sp->xl = XL; sp->xr = XR; sp->dy = DY; sp++;}
X
X#define POP(Y, XL, XR, DY) /* pop segment off stack */ \
X {sp--; Y = sp->y+(DY = sp->dy); XL = sp->xl; XR = sp->xr;}
X
X/*
X * fill: set the pixel at (x,y) and all of its 4-connected neighbors
X * with the same pixel value to the new pixel value nv.
X * A 4-connected neighbor is a pixel above, below, left, or right of a pixel.
X */
X
Xfill(x, y, win, nv)
Xint x, y; /* seed point */
XWindow *win; /* screen window */
XPixel nv; /* new pixel value */
X{
X int l, x1, x2, dy;
X Pixel ov; /* old pixel value */
X Segment stack[MAX], *sp = stack; /* stack of filled segments */
X
X ov = pixelread(x, y); /* read pv at seed point */
X if (ov==nv || x<win->x0 || x>win->x1 || y<win->y0 || y>win->y1) return;
X PUSH(y, x, x, 1); /* needed in some cases */
X PUSH(y+1, x, x, -1); /* seed segment (popped 1st) */
X
X while (sp>stack) {
X /* pop segment off stack and fill a neighboring scan line */
X POP(y, x1, x2, dy);
X /*
X * segment of scan line y-dy for x1<=x<=x2 was previously filled,
X * now explore adjacent pixels in scan line y
X */
X for (x=x1; x>=win->x0 && pixelread(x, y)==ov; x--)
X pixelwrite(x, y, nv);
X if (x>=x1) goto skip;
X l = x+1;
X if (l<x1) PUSH(y, l, x1-1, -dy); /* leak on left? */
X x = x1+1;
X do {
X for (; x<=win->x1 && pixelread(x, y)==ov; x++)
X pixelwrite(x, y, nv);
X PUSH(y, l, x-1, dy);
X if (x>x2+1) PUSH(y, x2+1, x-1, -dy); /* leak on right? */
Xskip: for (x++; x<=x2 && pixelread(x, y)!=ov; x++);
X l = x;
X } while (x<=x2);
X }
X}
END_OF_FILE
if test 2371 -ne `wc -c <'SeedFill.c'`; then
echo shar: \"'SeedFill.c'\" unpacked with wrong size!
fi
# end of 'SeedFill.c'
fi
if test -f 'SquareRoot.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'SquareRoot.c'\"
else
echo shar: Extracting \"'SquareRoot.c'\" \(2048 characters\)
sed "s/^X//" >'SquareRoot.c' <<'END_OF_FILE'
X/*
X * A High Speed, Low Precision Square Root
X * by Paul Lalonde and Robert Dawson
X * from "Graphics Gems", Academic Press, 1990
X */
X
X
X/*
X * SPARC implementation of a fast square root by table
X * lookup.
X * SPARC floating point format is as follows:
X *
X * BIT 31 30 23 22 0
X * sign exponent mantissa
X */
X
X#include <math.h>
X
Xstatic short sqrttab[0x100]; /* declare table of square roots */
X
Xvoid
Xbuild_table()
X{
X unsigned short i;
X float f;
X unsigned int *fi = &f; /* to access the bits of a float in */
X /* C quickly we must misuse pointers */
X
X for(i=0; i<= 0x7f; i++) {
X *fi = 0;
X
X /*
X * Build a float with the bit pattern i as mantissa
X * and an exponent of 0, stored as 127
X */
X
X *fi = (i << 16) | (127 << 23);
X f = sqrt(f);
X
X /*
X * Take the square root then strip the first 7 bits of
X * the mantissa into the table
X */
X
X sqrttab[i] = (*fi & 0x7fffff) >> 16;
X
X /*
X * Repeat the process, this time with an exponent of
X * 1, stored as 128
X */
X
X *fi = 0;
X *fi = (i << 16) | (128 << 23);
X f = sqrt(f);
X sqrttab[i+0x80] = (*fi & 0x7fffff) >> 16;
X }
X}
X
X/*
X * fsqrt - fast square root by table lookup
X */
Xfloat
Xfsqrt(float n)
X{
X unsigned int *num = &n; /* to access the bits of a float in C
X * we must misuse pointers */
X
X short e; /* the exponent */
X if (n == 0) return (0); /* check for square root of 0 */
X e = (*num >> 23) - 127; /* get the exponent - on a SPARC the */
X /* exponent is stored with 127 added */
X *num & = 0x7fffff; /* leave only the mantissa */
X if (e & 0x01) *num | = 0x800000;
X /* the exponent is odd so we have to */
X /* look it up in the second half of */
X /* the lookup table, so we set the */
X /* high bit */
X e >>= 1; /* divide the exponent by two */
X /* note that in C the shift */
X /* operators are sign preserving */
X /* for signed operands */
X /* Do the table lookup, based on the quaternary mantissa,
X * then reconstruct the result back into a float
X */
X *num = ((sqrttab[*num >> 16]) << 16) | ((e + 127) << 23);
X return(n);
X}
END_OF_FILE
if test 2048 -ne `wc -c <'SquareRoot.c'`; then
echo shar: \"'SquareRoot.c'\" unpacked with wrong size!
fi
# end of 'SquareRoot.c'
fi
if test -f 'Sturm/main.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'Sturm/main.c'\"
else
echo shar: Extracting \"'Sturm/main.c'\" \(2173 characters\)
sed "s/^X//" >'Sturm/main.c' <<'END_OF_FILE'
X/*
XUsing Sturm Sequences to Bracket Real Roots of Polynomial Equations
Xby D.G. Hook and P.R. McAree
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
X/*
X * main.c
X *
X * a sample driver program.
X */
X#include <stdio.h>
X#include <math.h>
X#include "solve.h"
X
X/*
X * a driver program for a root solver.
X */
Xmain()
X{
X poly sseq[MAX_ORDER];
X double min, max, roots[MAX_ORDER];
X int i, j, order, nroots, nchanges, np, atmin, atmax;
X
X /*
X * get the details...
X */
X
X printf("Please enter order of polynomial: ");
X scanf("%d", &order);
X
X printf("\n");
X
X for (i = order; i >= 0; i--) {
X printf("Please enter coefficient number %d: ", i);
X scanf("%lf", &sseq[0].coef[i]);
X }
X
X printf("\n");
X
X /*
X * build the Sturm sequence
X */
X np = buildsturm(order, sseq);
X
X printf("Sturm sequence for:\n");
X
X for (i = order; i >= 0; i--)
X printf("%lf ", sseq[0].coef[i]);
X
X printf("\n\n");
X
X for (i = 0; i <= np; i++) {
X for (j = sseq[i].ord; j >= 0; j--)
X printf("%lf ", sseq[i].coef[j]);
X printf("\n");
X }
X
X printf("\n");
X
X
X /*
X * get the number of real roots
X */
X nroots = numroots(np, sseq, &atmin, &atmax);
X
X if (nroots == 0) {
X printf("solve: no real roots\n");
X exit(0);
X }
X
X /*
X * calculate the bracket that the roots live in
X */
X min = -1.0;
X nchanges = numchanges(np, sseq, min);
X for (i = 0; nchanges != atmin && i != MAXPOW; i++) {
X min *= 10.0;
X nchanges = numchanges(np, sseq, min);
X }
X
X if (nchanges != atmin) {
X printf("solve: unable to bracket all negetive roots\n");
X atmin = nchanges;
X }
X
X max = 1.0;
X nchanges = numchanges(np, sseq, max);
X for (i = 0; nchanges != atmax && i != MAXPOW; i++) {
X max *= 10.0;
X nchanges = numchanges(np, sseq, max);
X }
X
X if (nchanges != atmax) {
X printf("solve: unable to bracket all positive roots\n");
X atmax = nchanges;
X }
X
X nroots = atmin - atmax;
X
X /*
X * perform the bisection.
X */
X sbisect(np, sseq, min, max, atmin, atmax, roots);
X
X
X /*
X * write out the roots...
X */
X if (nroots == 1) {
X printf("\n1 distinct real root at x = %f\n", roots[0]);
X } else {
X printf("\n%d distinct real roots for x: ", nroots);
X
X for (i = 0; i != nroots; i++)
X printf("%f ", roots[i]);
X printf("\n");
X }
X}
END_OF_FILE
if test 2173 -ne `wc -c <'Sturm/main.c'`; then
echo shar: \"'Sturm/main.c'\" unpacked with wrong size!
fi
# end of 'Sturm/main.c'
fi
if test -f 'TriPoints.c' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'TriPoints.c'\"
else
echo shar: Extracting \"'TriPoints.c'\" \(2746 characters\)
sed "s/^X//" >'TriPoints.c' <<'END_OF_FILE'
X/*
XGenerating Random Points in Triangles
Xby Greg Turk
Xfrom "Graphics Gems", Academic Press, 1990
X*/
X
X#include <math.h>
X#include "GraphicsGems.h"
X
X/*****************************************************************
XCompute relative areas of sub-triangles that form a convex polygon.
XThere are vcount-2 sub-triangles, each defined by the first point
Xin the polygon and two other adjacent points.
X
XThis procedure should be called once before using
Xsquare_to_polygon().
X
XEntry:
X vertices - list of the vertices of a convex polygon
X vcount - number of vertices of polygon
XExit:
X areas - relative areas of sub-triangles of polygon
X*******************************************************************/
X
Xtriangle_areas (vertices, vcount, areas)
X Point3 vertices[];
X int vcount;
X float areas[];
X{
X int i;
X float area_sum = 0;
X Vector3 v1,v2,v3;
X
X /* compute relative areas of the sub-triangles of polygon */
X
X for (i = 0; i < vcount - 2; i++) {
X V3Sub(&vertices[i+1], &vertices[0], &v1);
X V3Sub(&vertices[i+2], &vertices[0], &v2);
X V3Cross(&v1, &v2, &v3);
X areas[i] = LengthVector3(&v3);
X area_sum += areas[i];
X }
X
X /* normalize areas so that the sum of all sub-triangles is one */
X
X for (i = 0; i < vcount - 2; i++)
X areas[i] /= area_sum;
X}
X
X
X/******************************************************************
XMap a point from the square [0,1] x [0,1] into a convex polygon.
XUniform random points in the square will generate uniform random
Xpoints in the polygon.
X
XThe procedure triangle_areas() must be called once to compute
X'areas', and then this procedure can be called repeatedly.
X
XEntry:
X vertices - list of the vertices of a convex polygon
X vcount - number of vertices of polygon
X areas - relative areas of sub-triangles of polygon
X s,t - position in the square [0,1] x [0,1]
XExit:
X p - position in polygon
X*********************************************************************/
X
Xsquare_to_polygon (vertices, vcount, areas, s, t, p)
X Point3 vertices[];
X int vcount;
X float areas[];
X float s,t;
X Point3 *p;
X{
X int i;
X float area_sum = 0;
X float a,b,c;
X
X /* use 's' to pick one sub-triangle, weighted by relative */
X /* area of triangles */
X
X for (i = 0; i < vcount - 2; i++) {
X area_sum += areas[i];
X if (area_sum >= s)
X break;
X }
X
X /* map 's' into the interval [0,1] */
X
X s = (s - area_sum + areas[i]) / areas[i];
X
X /* map (s,t) to a point in that sub-triangle */
X
X t = sqrt(t);
X
X a = 1 - t;
X b = (1 - s) * t;
X c = s * t;
X
X p->x = a * vertices[0].x + b * vertices[i+1].x + c * vertices[i+2].x;
X p->y = a * vertices[0].y + b * vertices[i+1].y + c * vertices[i+2].y;
X p->z = a * vertices[0].z + b * vertices[i+1].z + c * vertices[i+2].z;
X}
END_OF_FILE
if test 2746 -ne `wc -c <'TriPoints.c'`; then
echo shar: \"'TriPoints.c'\" unpacked with wrong size!
fi
# end of 'TriPoints.c'
fi
echo shar: End of archive 2 \(of 5\).
cp /dev/null ark2isdone
MISSING=""
for I in 1 2 3 4 5 ; do
if test ! -f ark${I}isdone ; then
MISSING="${MISSING} ${I}"
fi
done
if test "${MISSING}" = "" ; then
echo You have unpacked all 5 archives.
rm -f ark[1-9]isdone
else
echo You still need to unpack the following archives:
echo " " ${MISSING}
fi
## End of shell archive.
exit 0
More information about the Comp.sources.misc
mailing list