#! /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 shell archive."
# Contents:  Makefile SCCS make_rotation.c quadtofile.c torusknot.c
#   tuber.c tubeseg.c vec4.h xform.c
# Wrapped by gunn@bolyai on Tue Jul 17 11:26:36 1990
PATH=/bin:/usr/bin:/usr/ucb ; export PATH
if test -f 'Makefile' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'Makefile'\"
else
echo shar: Extracting \"'Makefile'\" \(766 characters\)
sed "s/^X//" >'Makefile' <<'END_OF_FILE'
XHACKS = -g
XPICIO = /usr/pixar/host/lib/libpicio.a
XINCLUDES = -I/u/gunn/include -I/usr/pixar/include -I/usr/pixar/chreyes/include -I/usr/pixar/chreyes/host/src/lib/libchreyes
X#INCLUDES = -I/u/gunn/include -I/usr/pixar/include -I/usr/pixar/chreyes/include 
XCFLAGS = ${HACKS} ${INCLUDES} -yVcosd
X
Xtuber: tuber.o quadtofile.o make_rotation.o 
X	cc ${CFLAGS} tuber.o quadtofile.o make_rotation.o ${PICIO} /usr/pixar/host/lib/libaa.a /usr/pixar/host/lib/libaa.i.a -lm -o tuber
X	mv tuber ..
X
Xtubeseg: tubeseg.o quadtofile.o make_rotation.o 
X	cc tubeseg.o quadtofile.o make_rotation.o ${PICIO} /usr/pixar/host/lib/libaa.a -lm ${CFLAGS} -o tubeseg
X
Xtorusknot: torusknot.o
X	cc torusknot.o -lm -o torusknot
X	mv torusknot ..
X
Xxform: xform.o
X	cc xform.o -lm -o xform
X	mv xform ..
END_OF_FILE
if test 766 -ne `wc -c <'Makefile'`; then
    echo shar: \"'Makefile'\" unpacked with wrong size!
fi
# end of 'Makefile'
fi
if test ! -d 'SCCS' ; then
    echo shar: Creating directory \"'SCCS'\"
    mkdir 'SCCS'
fi
if test -f 'make_rotation.c' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'make_rotation.c'\"
else
echo shar: Extracting \"'make_rotation.c'\" \(4280 characters\)
sed "s/^X//" >'make_rotation.c' <<'END_OF_FILE'
X#include <math.h>
X#include "vec4.h"
X
X	static matrix4 last_mat = {1,0,0,0,0,1,0,0,0,0,1,0,0,0,01};;
X
Xmake_tform(p1, p2, p3, m)
Xpoint4 p1, p2, p3;
Xmatrix4 m;
X/* 	Generate a Euclidean isometry which moves (0,0,0) -> p1.
X   Furthermore, 
X	the vector (1,0,0) will go a vector tx based at p1 and ending at p2.
X	the vector (0,1,0) will go a vector ty based at p1, perpendicular to tx,
X   and lying in the plane of p1-p2-p3, in the "direction" of p3.	
X	the vector (0,0,1) will go a vector tz based at p1, perpendicular to tx and ty.
X	The rotational part is computed by make_rotation. See comments there */
X{
X    register int i;
X    point4 nv1, nv2;
X
X    for (i=0; i<3; ++i)
X	{
X	nv1[i] = p2[i] - p1[i];
X	nv2[i] = p3[i] - p2[i];
X	}
X    nv1[3] = 1.0; nv2[3] = 1.0;
X    if (VDOT3(nv1,nv1) == 0.0)
X	{
X	return(-1);
X	}
X    if (VDOT3(nv2,nv2) == 0.0)
X	{
X	if (nv1[0] != 0.0 || nv1[1] != 0.0)
X	    {
X	    nv2[0] = nv1[1];
X 	    nv2[1] = -nv1[0];
X	    nv2[2] = 0;
X	    }
X	else
X	    {
X	    nv2[1] = nv1[2];
X	    nv2[2] = -nv1[1];
X	    nv2[0] = 0;
X	    }	
X	}
X	
X    normalize3(nv1);
X    normalize3(nv2);
X    if ( fabs(VDOT3(nv1,nv2) > .99999))
X	{
X	if (nv1[0] != 0.0 || nv1[1] != 0.0)
X	    {
X	    nv2[0] = nv1[1];
X 	    nv2[1] = -nv1[0];
X	    nv2[2] = 0;
X	    }
X	else
X	    {
X	    nv2[1] = nv1[2];
X	    nv2[2] = -nv1[1];
X	    nv2[0] = 0;
X	    }	
X	}
X    if ( !make_rotation(nv1, nv2, m))	/* if collinear, use most recent rotation */
X	COPYMAT4(last_mat, m)
X    else 
X	COPYMAT4(m, last_mat)
X    /* compute translation to take p1->(0,0,0) */
X    for (i=0; i<3; ++i)
X	{
X	m[i][3] = p1[i];
X	m[3][i] = 0.0;
X	}
X    m[3][3] = 1.0;
X    }
X
Xnormalize3(v)
Xpoint4 v;
X{
X    register int i; 
X    double tsum;
X    for (i=0, tsum=0; i<3; ++i)
X	tsum += v[i]*v[i];
X    tsum = sqrt(tsum);
X    if (tsum == 0.0)	tsum = 1.0;
X    for (i=0; i<3; ++i)
X	v[i] = v[i]/tsum;
X}
X
X/* make rotation that takes (1,0,0) to v1; and (0,1,0) to the orthog. proj of v2 onto
X   the perpendicular subspace of v1.  Return (0) if v1 and v2 are collinear. 
X*/
Xmake_rotation(v1, v2, m)
Xpoint4 v1, v2;
Xmatrix4 m;
X{
X    int i;
X    double a;
X    point4 t1, t2;
X
X    a = VDOT3(v1, v2);
X    if (a > .9999 || a < -.9999)  return(0); 
X    VSCALE3(v1, a, t1)
X    VSUB3(v2, t1, t2) 		/* now t2 is orthogonal to v1 */
X    normalize3(t2);
X    m[0][0] = v1[0];	m[1][0] = v1[1]; 	m[2][0] = v1[2];
X    m[0][1] = t2[0];	m[1][1] = t2[1]; 	m[2][1] = t2[2];
X    /* last row is cross product of the first two rows */
X    m[0][2] = v1[1]*t2[2] - v1[2]*t2[1];
X    m[1][2] = v1[2]*t2[0] - v1[0]*t2[2];
X    m[2][2] = v1[0]*t2[1] - v1[1]*t2[0];
X
X    for (i=0; i<3; ++i)
X	m[i][3] = m[3][i] = 0.0;
X    m[3][3] = 1.0;
X    return(1);
X/* 
X    double x,y,z,tt,tmp;
X    x = v[0];
X    y = v[1];
X    z = v[2];
X    tt = y*y + z*z;
X    if (tt == 0.0)
X	{
X	makeid4(m);
X	return;
X	}
X    m[0][0] = x;
X    m[1][0] = y;
X    m[2][0] = z;
X    m[0][1] = -y;
X    m[0][2] = -z;
X    m[1][1] = (z*z + y*y*x)/tt;
X    m[2][2] = (y*y + z*z*x)/tt;
X    m[2][1] = tmp = (x*y*z - y*z)/tt;
X    m[1][2] = -tmp;
X*/
X
X}
X
Xmakeid4(m)
Xmatrix4 m;
X{
X    register int i, j;
X    for (i=0; i<4; ++i)
X        for (j=0; j<4; ++j)	m[i][j] = (i == j);
X}
X
Xcopymat4(m1, m2)
Xmatrix4 m1, m2;
X{
X    register int i, j;
X    for (i=0; i<4; ++i)
X        for (j=0; j<4; ++j)	m2[i][j] = m1[i][j];
X}
X
Xmatmatmul4(m1,m2,m3)
Xmatrix4 m1, m2, m3;
X{
X    matrix4 mtmp;
X    double tmp;
X    register int i,j,k;
X
X    for (i=0; i<4; ++i)
X	for (j=0;j<4;++j)
X	    {
X	    for (tmp = 0, k=0; k<4; ++k)
X		tmp += m1[i][k]*m2[k][j];
X	    mtmp[i][j] = tmp;
X	    }	
X    for (i=0; i<4; ++i)
X	for (j=0;j<4;++j)	m3[i][j] = mtmp[i][j];
X}
X
Xaxis_rotation(axis, angle, mat)
Xdouble angle;
Xmatrix4 mat;
X{
X    register int i,j, k;
X    double c,s;
X
X    switch (axis)	{
X	case 0:
X	    i = 0;
X	    j = 1;
X	    k = 2;
X	    break;	
X	case 1:
X	    i = 1;
X	    j = 0;
X	    k = 2;
X	    break;	
X	case 2:
X	    i = 2;
X	    j = 0;
X	    k = 1;
X	    break;	
X	}
X
X    c = cos(angle);
X    s = sin(angle);
X
X    makeid4(mat);
X    mat[j][j] = mat[k][k] = c;
X    mat[k][j] = -s;
X    mat[j][k] = s;
X}
Xmatvecmul4(m,v1,v2)
Xmatrix4 m;
Xpoint4 v1,v2;
X{
X    register int i,j,k;
X    double tmp[4], *tv1, *tv2;
X 
X    tv1 = (double *) v1; 
X    tv2 = (double *) v2;
X    for (i=0; i<4; ++i)
X	for (tmp[i]=0, j = 0; j<4; ++j)
X	    tmp[i] += m[i][j] * tv1[j];
X    for (i=0; i<4; ++i)
X	tv2[i] = tmp[i];
X}
X
END_OF_FILE
if test 4280 -ne `wc -c <'make_rotation.c'`; then
    echo shar: \"'make_rotation.c'\" unpacked with wrong size!
fi
# end of 'make_rotation.c'
fi
if test -f 'quadtofile.c' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'quadtofile.c'\"
else
echo shar: Extracting \"'quadtofile.c'\" \(5243 characters\)
sed "s/^X//" >'quadtofile.c' <<'END_OF_FILE'
X#include <gobject.h>
X#include <picio.h>
X#include <pixeldef.h>
X#include <conversion.h>
X
Xstatic struct { int type; char *pre, *lc_pre; } quadtypes[] = {
X	{ QUAD,	 "",	"" },
X	{ NQUAD, "N",	"n" },
X	{ CQUAD, "C",	"c" },
X	{ CNQUAD,"CN",	"cn"},
X	{ 0,  "NOT-A-", "not-a-"}
X};
Xint
XQuadtoFile(type,numquad,v, outfile)
Xint	type;
Xint	numquad;
Xfloat	*v;		/* float v[numquad][4][3,6,7 or 10]; */
Xchar 	*outfile;
X{
X	register int		i, off;
X	register RGBAPixelType	*b;
X	RGBAPixelType		buf[1024];
X	int			rows;
X	PFILE			*pf;
X	int 			debug=1;
X
X
X	rows = (numquad/ ((type & GPM_NORMALS) ? 32 : 64)) + 1;
X	if (debug)	printf("quads:	%d	rows:	%d\n",numquad,rows);
X	PicSetPsize(1024, rows);
X	PicSetTsize(1024,rows);
X	PicSetPformat(PF_RGBA);
X	PicSetPstorage(PS_12DUMP);
X	PicSetPmatting(PM_NONE);
X	PicSetLabel("Created by QuadtoFile\n");
X
X	pf = PicCreat(outfile,0755);
X	picPreEncodeScanline(pf);
X
X	buf[0].Red	=numquad;
X	buf[0].Green	=(numquad << 8) + (numquad >> 8);
X	buf[0].Blue 	= rows;
X	buf[0].Alpha	=(rows << 8) + (rows >> 8);
X	
X	b = buf + ((type & GPM_NORMALS) ? 4*(4+4) : 4*4); /* Offset one quad */
X	off = (type & GPM_NORMALS) ? 3 : 0;
X
X	for (i = 0; i<numquad*4; ++i) {
X		FLTOPIX16(*v, b); ++v; ++b;
X		FLTOPIX16(*v, b); ++v; ++b;
X		FLTOPIX16(*v, b); ++v; ++b;
X
X		if (type & GPM_COLORS) {
X			FLTOPIXPK11(v[off+0],v[off+1],v[off+2],v[off+3], b);
X		} else {
X			FLTOPIXINT((float)0.0, b);
X		}
X		++b;
X
X		if (type & GPM_NORMALS) {
X			FLTOPIX14(*v, b); ++v; ++b;
X			FLTOPIX14(*v, b); ++v; ++b;
X			FLTOPIX14(*v, b); ++v; ++b;
X			FLTOPIXINT((float)0.0, b); ++b;
X		}
X		if (type & GPM_COLORS) v += 4;
X	
X		if (b >= buf+1024) {
X			picEncodeScanline(pf, buf);
X			bzero(buf, sizeof(buf));
X			b = buf;
X		}
X	}
X
X	/* send remaining data */
X	if (b > buf) {
X		picEncodeScanline(pf, buf);
X		bzero(buf, sizeof(buf));
X	}
X
X	picPostEncodeScanline(pf);
X	PicClose(pf);
X
X}
X
Xint
XBbptoFile(type,numpatch,v, outfile)
Xint	type;
Xint	numpatch;
Xfloat	*v;		/* float v[numpatch]([16][3] {+ [8]}); */
Xchar  	*outfile;
X{
X	register int		i,j,k;
X	RGBAPixelType		*b, buf[1024];
X	int			rows;
X	PFILE			*pf;
X	int 			debug=1;
X
X	rows = numpatch/ 15 + 1;
X	if (debug)	printf("patches:  %d	rows:	%d\n",numpatch,rows);
X	PicSetPsize(1024, rows);
X	PicSetTsize(1024,rows);
X	PicSetPformat(PF_RGBA);
X	PicSetPstorage(PS_12DUMP);
X	PicSetPmatting(PM_NONE);
X	PicSetLabel("Created by BbptoFile\n");
X
X	pf = PicCreat(outfile,0755);
X	picPreEncodeScanline(pf);
X
X	bzero(buf, sizeof(buf));
X	buf[0].Red	= numpatch;
X	buf[0].Green	=(numpatch << 8) + (numpatch >> 8);
X	buf[0].Blue 	= rows;
X	buf[0].Alpha	=(rows << 8) + (rows >> 8);
X
X	b = buf + 68;		/* skip over one patch space */
X	for (i=1; i<=numpatch; ++i)	/* Skip one space for count */
X		{
X  		convert_patch (type, &v, &b);
X			
X		if ((i % 15) == 14) {	/* 15 patches done */
X			picEncodeScanline(pf, buf);
X			bzero(buf, sizeof(buf));
X			b = buf;
X			}
X		}
X
X	if ((i % 15) != 0) /* flush the buffer */
X		{
X		picEncodeScanline(pf, buf);
X		}
X
X	picPostEncodeScanline(pf);
X	PicClose(pf);
X}
X
X/*  48 (56) floats at location "patch" converted to pixar fb format,
X    stored in buffer, and patch id number stored in its place too */
Xstatic
Xconvert_patch (type, patch, buf)
X    int type;
X    float **patch;
X    RGBAPixelType **buf;
X    {
X    register int j,k;
X
X        for (j=0;j<16;++j) {
X            for (k=0;k<3;++k) {
X                FLTOPIX16(**patch, *buf);
X                ++*patch; ++*buf;
X            }
X	    FLTOPIX16(1.0, *buf);
X            ++*buf;
X        }
X        if (type & GPM_ST) {
X	    for (k = 0; k < 4; ++k) {
X	        PACK_ST_COORDS((*patch)[0], (*patch)[1], *buf);
X		*patch += 2; ++*buf;
X	    }
X        } else {
X	    PACK_ST_COORDS(0.0, 0.0, *buf); ++*buf;
X	    PACK_ST_COORDS(1.0, 0.0, *buf); ++*buf;
X	    PACK_ST_COORDS(0.0, 1.0, *buf); ++*buf;
X	    PACK_ST_COORDS(1.0, 1.0, *buf); ++*buf;
X        }
X    }
X
X#ifdef FOODACIOUS
X	if (pr_ascii) {
X	    register int i, j;
X	    char fn[80];
X	    FILE *f;
X
X	    for(i=0; quadtypes[i].type!=0; ++i)
X		if (type==quadtypes[i].type) break;
X	    sprintf(fn, "%s.%s.%squad", pr_aprefix, name, quadtypes[i].lc_pre);
X	    f = fopen(fn, "w");
X
X	    fprintf(f, "%sQUAD\n", quadtypes[i].pre);
X	    for(j=0; j<numquad*4; ++j) {
X		fprintf(f, "%g ", *v++);
X		fprintf(f, "%g ", *v++);
X		fprintf(f, "%g ", *v++);
X		if (type&GPM_NORMALS) {
X		    fprintf(f, "%g ", *v++);
X		    fprintf(f, "%g ", *v++);
X		    fprintf(f, "%g ", *v++);
X		}
X		if (type&GPM_COLORS) {
X		    fprintf(f, "%g ", *v++);
X		    fprintf(f, "%g ", *v++);
X		    fprintf(f, "%g ", *v++);
X		    fprintf(f, "%g ", *v++);
X		}
X	    }
X
X	    fclose(f);
X	    fprintf(pr_ascii, "gprim %s %squad %s nosave\n",
X			name, quadtypes[i].lc_pre, fn);
X	    return 1;
X	}
X
X	if (pr_ascii) {
X		int n;
X		char fn[80];
X		FILE *f;
X
X		for(n=0; bbptypes[n].type!=0; ++n)
X			if (type == bbptypes[n].type) break;
X		sprintf(fn, "%s.%s.%sbbp", pr_aprefix,name,bbptypes[n].lc_pre);
X		f = fopen(fn, "w");
X
X		fprintf(f, "%sBBP\n", bbptypes[n].pre);
X		for (i = 1; i <= numpatch; ++i) {
X			for(j=0; j<16; ++j) {
X				for(k=0; k<3; ++k)
X					fprintf(f, "%g ", *v++);
X				fprintf(f, "\n");
X			}
X			if (type & GPM_ST)
X				for(k=0; k<8; ++k)
X					fprintf(f, "%g ", *v++);
X			fprintf(f, "\n");
X		}
X		fclose(f);
X
X		fprintf(pr_ascii, "gprim %s %sbbp %s nosave\n",
X				  name, bbptypes[n].lc_pre, fn);
X		return 1;
X	}
X
X#endif
END_OF_FILE
if test 5243 -ne `wc -c <'quadtofile.c'`; then
    echo shar: \"'quadtofile.c'\" unpacked with wrong size!
fi
# end of 'quadtofile.c'
fi
if test -f 'torusknot.c' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'torusknot.c'\"
else
echo shar: Extracting \"'torusknot.c'\" \(833 characters\)
sed "s/^X//" >'torusknot.c' <<'END_OF_FILE'
X/*  torusknot.c:	generate torus knots */
X#include <math.h>
X#define	TWOPI	6.2831837
X	int bignum, littlenum, numsegs;
X	double phi, theta, bigrad, littlerad;
X
Xmain(argc, argv)	char **argv;
X{
X
X	register int i,k;
X	double ct,st,cp,sp, x,y,z;
X
X	if (argc != 6)	
X	    {
X	    printf("Usage: bignum, littlenum, numsegs, bigrad, littlerad\n");
X	    exit(1);
X	    }
X	bignum = atoi(*++argv);
X	littlenum = atoi(*++argv);
X	numsegs = atoi(*++argv);
X	bigrad = atof(*++argv);
X	littlerad = atof(*++argv);
X
X	for (i=0; i<=(numsegs+1); ++i)	
X	    {
X	    theta = (bignum*i*TWOPI)/numsegs;
X	    phi = (littlenum*i*TWOPI)/numsegs;
X	    cp = cos(phi);
X	    sp = sin(phi);
X	    ct = cos(theta);
X	    st = sin(theta);
X	    x = ct*(bigrad + littlerad*cp); 
X	    y = st*(bigrad + littlerad*cp); 
X	    z = littlerad*sp; 
X	    printf("%f\t%f\t%f\n",x,y,z);
X	    }
X}
END_OF_FILE
if test 833 -ne `wc -c <'torusknot.c'`; then
    echo shar: \"'torusknot.c'\" unpacked with wrong size!
fi
# end of 'torusknot.c'
fi
if test -f 'tuber.c' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'tuber.c'\"
else
echo shar: Extracting \"'tuber.c'\" \(18815 characters\)
sed "s/^X//" >'tuber.c' <<'END_OF_FILE'
X#include "vec4.h"
X#include "gobject.h"
X#include <stdio.h>
X#include <string.h>
X#include <aarg.h>
X#include <math.h>
X#define LSIZE 8192
X#define MAXPOINTS 64
X#define	MAXLINKS	32
X/*
X#define DEBUG
X*/
X
X	float *v;		/* pointer to list of quadrilateral data */
X	int vptr, BlockSize, array_size;
X	int debug = 1;
X
X	static char closed[] = "uclosed";
X	static char open[] = "uopen";
X 	point4 profile[] = { 
X		0.0, 0.5,     0.0, 1.0,      
X		0.0, 0.5,     0.27614 ,  1.0,
X		0.0, 0.27614,      0.5, 1.0, 
X		0.0, 0.0,       0.5,  1.0, 
X		0.0, 0.0,       0.5,  1.0, 
X		0.0, -0.27614,      0.5, 1.0, 
X		0.0, -0.5,     0.27614 ,  1.0,     
X		0.0, -0.5,     0.0, 1.0,      
X		0.0, -0.5,     0.0, 1.0,      
X		0.0, -0.5,     -0.27614 ,  1.0,     
X		0.0, -0.27614 ,     -0.5, 1.0, 
X		0.0, 0.0,       -0.5,  1.0, 
X		0.0, 0.0,       -0.5,  1.0, 
X		0.0, 0.27614 ,     -0.5, 1.0, 
X		0.0, 0.5,     -0.27614 ,  1.0,     
X		0.0, 0.5,     0.0, 1.0 
X
X	} ;
X
X 	point4 *vlist[MAXLINKS], *profile_init, *profilem,
X  		*profile0, *profile1, *profile2, *profile3;
X	double totaltwist = 0.0;
X 	matrix4 lastmat;
X	point4 y_vec = {0,1,0,1};
X	point4 z_vec = {0,0,1,1};
X        char outfile[128], profile_file[64] = {""};
X	char softname[128] = "link";
X	int vsize, profile_size;
X	int make_binary_pixar = 0;
X	int firsttime, wrap = 0, split= 0;
X	int numlinks, linksize[MAXLINKS];
X	FILE *outf, *inf;
X 	int total, totalon = 0, vcnt;
X	int twist, twiston = 0;
X	int softimage = 0;
X	double radius, length[MAXLINKS];
X
Xdouble
Xgetlength(vlist, nn)
Xpoint4 *vlist;
X{
X	int i;
X	double ll;
X	point4 diff0, diff1, diff2;
X	for (ll = 0.0, i=0; i<nn; ++i)
X		{
X		VADD3(vlist[(i+2)%nn],vlist[(i+1)%nn],diff0);	
X		VADD3(vlist[(i+1)%nn],vlist[i],diff1);	
X		VSUB3(diff1, diff0, diff2);
X		VSCALE3(diff2, .5, diff2);
X		ll +=  MAGNITUDE3(diff2);
X		}
X	return(ll);
X}
X
Xremain(argc, argv)	char **argv;
X{
X    
X    int numseg, i, linkno, nn;
X    double atof(), angle;
X    double ol, ll;
X    point4 pbuffer[MAXPOINTS], *profileptr, diff0, diff1, diff2;
X    FILE *fd, *fopen();
X    
X    if (argc == 1) 
X	exit(printf("Need a file name\n"));
X
X    inf = fopen(argv[1],"r");
X    if (inf == NULL)	exit(fprintf(stderr,"Unable to open infile %s\n",argv[1]));
X
X
X    /* process profile curve */
X    if (strlen(profile_file) == 0)	/* use standard circular arc */
X	{
X        profile_size = sizeof(profile)/sizeof(point4);
X	profileptr = profile;
X	}
X    else			/* read profile curve from file */
X	{
X	if ( (fd = fopen(profile_file,"r")) == NULL)
X	    exit(fprintf(stderr,"Unable to open profile curve file\n"));
X	for (i = 0, profileptr = pbuffer; ; ++i )
X	    {
X	    if (fscanf(fd,"%F%F%F\n",
X	        &(profileptr[i][0]), &(profileptr[i][1]),&(profileptr[i][2])) != 3)	
X		break;
X	    profileptr[i][3] = 1.0;
X	    }
X
X	fclose(fd);
X	profileptr = pbuffer;
X	profile_size = i;
X	}
X	
X    /* Allocate space for the various profile curves which we'll need */
X    profile_init = (point4 *) malloc (sizeof(point4) * profile_size);
X    profile0 = (point4 *) malloc (sizeof(point4) * profile_size);
X    profile3 = (point4 *) malloc (sizeof(point4) * profile_size);
X    profilem = (point4 *) malloc (sizeof(point4) * profile_size);
X    profile1 = (point4 *) malloc (sizeof(point4) * profile_size);
X    profile2 = (point4 *) malloc (sizeof(point4) * profile_size);
X
X    for (i=0; i<profile_size; ++i)	
X	{
X	profile_init[i][0] = profileptr[i][0]*radius;
X	profile_init[i][1] = profileptr[i][1]*radius;
X	profile_init[i][2] = profileptr[i][2]*radius;
X	profile_init[i][3] = profileptr[i][3];
X	}
X
X    /* read the input link file */
X    {
X    char type[MAXLINKS];
X    fscanf(inf,"%s",type);
X    if (type[0] <= 'z'  && type[0] >= 'A')		/* new version: type identifier  at top */
X	{
X	if (type[0] != 'L')	
X	    exit(fprintf(stderr,"tuber: no such type: %s\n",type));
X
X	fscanf(inf,"%d",&numlinks);
X	for (i=0; i<numlinks; ++i)
X	    fscanf(inf,"%d", &(linksize[i]) ); 
X	for (linkno=0; linkno<numlinks; ++linkno)
X	    {
X	    if ( (vlist[linkno] = (point4 *) malloc(linksize[linkno] * sizeof(point4))) == NULL)
X		exit(fprintf(stderr,"Tuber: Out of space getting linkno %d\n",linkno));
X	    for (numseg=0; numseg<linksize[linkno]; ++numseg)
X	        {
X	        if ( (fscanf(inf,"%F%F%F",vlist[linkno][numseg], vlist[linkno][numseg]+1, vlist[linkno][numseg]+2) < 3))
X		    exit(fprintf(stderr,"tuber: premature EOF reading input\n"));
X#ifdef DEBUG
X	fprintf(stderr,"%f %f %f\n",vlist[linkno][numseg][0],vlist[linkno][numseg][1],vlist[linkno][numseg][2]);
X#endif
X	        vlist[linkno][numseg][3] = 1.0;
X	        }
X	    length[linkno] = getlength(vlist[linkno], linksize[linkno]);
X	    }
X	}
X    else
X      {
X      vsize = 1;
X      vlist[0] = (point4 *) malloc(LSIZE * sizeof(point4));
X      /* get list of points */
X      for (numseg = 0 ; ; ++numseg)
X	{
X	if (numseg+2 > LSIZE)
X	    {
X	    vsize *= 2;
X	    vlist[0] = (point4 *) realloc(LSIZE * vsize * sizeof(point4));
X	    printf("Growing vertex list: size is %d\n",vsize*LSIZE);
X	    }
X	if (numseg == 0)
X	    {
X	    vlist[0][numseg][0] = atof(type);	/* convert the string we've read */
X	    scanf("%F%F",vlist[0][numseg]+1, vlist[0][numseg]+2); 
X	    }
X	else if ( (scanf("%F%F%F",vlist[0][numseg], vlist[0][numseg]+1, vlist[0][numseg]+2) < 3)) break;
X	vlist[0][numseg][3] = 1.0;
X	}
X      linksize[0] = numseg;
X      numlinks = 1;
X      }
X    }
X
X    /* for single file output, do prep work outside loop */
X    if (!split)
X	{
X	if (softimage)
X		{
X		sprintf(outfile,"%s.def",softname);
X	    	outf = fopen(outfile,"w");
X	    	if (outf == NULL)	exit(fprintf(stderr,"Unable to open outfile %s\n",outfile));
X		fprintf(outf,"SOFTIMAGE 4D Creative Environment   v 1.55  \"ASCII\"\n\n");
X		}
X
X	else	
X	    {
X	    if (!make_binary_pixar)
X	    	{
X		sprintf(outfile,"%s.bbp",argv[1]);
X	    	outf = fopen(outfile,"w");
X	    	if (outf == NULL)	
X		    exit(fprintf(stderr,"Unable to open outfile %s\n",outfile));
X            	fprintf(outf,"STBBP\n"); /* ascii header for Pixar batch files */
X		}
X    	    else
X	    	{
X		sprintf(outfile,"%s.bbp.fbi",argv[1]);
X	    	outf = fopen(outfile,"w");
X	    	if (outf == NULL)	
X		    exit(fprintf(stderr,"Unable to open outfile %s\n",outfile));
X	    	vptr = 0;
X	    	array_size = 1;
X	    	BlockSize = 8192;
X	    	v = (float *) malloc (sizeof(float) * BlockSize );
X	    	}
X	    }
X	}
X	
X    /* loop through each component of link */
X    for (linkno = 0; linkno < numlinks ; ++ linkno)
X	{
X	if (split)	/* compute separate file name for this component */
X	    {
X	    if (softimage)
X		{
X		sprintf(outfile,"%s.%d.def",softname,linkno);
X	    	outf = fopen(outfile,"w");
X	    	if (outf == NULL)	exit(fprintf(stderr,"Unable to open outfile %s\n",outfile));
X		fprintf(outf,"SOFTIMAGE 4D Creative Environment   v 1.55  \"ASCII\"\n\n");
X		}
X	    else {
X		if (!make_binary_pixar)
X		    {
X		    sprintf(outfile,"%s.%d.bbp",argv[1],linkno); 
X	    	    outf = fopen(outfile,"w");
X	    	    if (outf == NULL)	exit(fprintf(stderr,"Unable to open outfile %s\n",outfile));
X    	    	    fprintf(outf,"STBBP\n");	
X		    }
X	    	else
X		    {
X	 	    sprintf(outfile,"%s.%d.bbp.fbi",argv[1],linkno); 
X	    	    outf = fopen(outfile,"w");
X	    	    if (outf == NULL)	exit(fprintf(stderr,"Unable to open outfile %s\n",outfile));
X	    	    vptr = 0;
X	    	    array_size = 1;
X	    	    BlockSize = 8192;
X	    	    v = (float *) malloc (sizeof(float) * BlockSize );
X	    	    }
X		}
X	    }
X
X 	if (softimage)	{
X	    fprintf(outf,"PTCH \"%s.%d\"\n{\n",softname,linkno);
X	    fprintf(outf,"   utype\tBEZIER\n   upoint	%d\n   %s\n   ustep\t10\n   ucurve\t1\n\n", wrap ? 3*(linksize[linkno]) : 3*(linksize[linkno]-2), wrap ? closed : open);
X	    fprintf(outf,"   vtype\tBEZIER\n   vpoint	%d\n   vclosed\n   vstep\t10\n   vcurve\t1\n\n",3*(profile_size/4));
X	    fprintf(outf,"   point\n");
X	    }
X      for (firsttime = 1, angle = 0.0, i=0, ol = 0.0, ll = 0.0; i < (wrap ? linksize[linkno] : linksize[linkno] - 2); ++i)
X	{
X	nn = linksize[linkno];
X	VADD3(vlist[linkno][(i+2)%nn],vlist[linkno][(i+1)%nn],diff0);	
X	VADD3(vlist[linkno][(i+1)%nn],vlist[linkno][i],diff1);	
X	VSUB3(diff1, diff0, diff2);
X	VSCALE3(diff2, .5, diff2);
X	ll += MAGNITUDE3(diff2);
X	/* meat of the matter here: compute one segment of tube */
X   	make_cylinder(vlist[linkno][i], vlist[linkno][(i+1)%linksize[linkno]], 
X		vlist[linkno][(i+2)%linksize[linkno]], 
X		ol/length[linkno], ll/length[linkno]); 
X	ol = ll;
X  	}	
X	if (softimage) fprintf(outf,"}\n");
X 	else if (split && make_binary_pixar) BbptoFile(STBBP, vptr/56, v, outfile);
X	else fclose(outf);
X	}
X
X    if (!split && make_binary_pixar)
X	BbptoFile(STBBP, vptr/56, v, outfile);
X    fclose(outf);
X
X}	
X
X/* make_cylinder creates a bicubic patch tube which goes from the midpoint
X   of the segment v0-v1 and ends at the midpoint of segment v1-v2.
X   It is also "tangent" to these segments at the beginning and end.
X   	The central computation is make_tform, which computes the Euclidean
X   isometry which sends a standard frame at the origin to a frame based
X   at a given point, and the y-z plane to the plane through this
X   point perpendicular to given segment.  
X	(Simultaneously another rotation is prefixed which maintains
X   the continuity of texture parameters from one segment to the next.)
X	The composite transformation is applied
X   to the profile curve to yield beginning and ending control points for the tube:
X   call these curves p0 and p3.
X	A third transformation is computed which yields a profile curve
X   around the point v1, call it pm.  Intuitively, this is the intersection
X   of the two linear cylinders around the lines v0-v1 and v1-v2.
X	Intermediate control point curves p1 and p2 are computed by
X   taking a linear interpolation between p0 and pm and pm and p3 with certain
X   nice properties.  To be exact, 
X	p1 = p0 + f(alpha) * (pm - p0)
X	p2 = p3 + f(alpha) * (pm - p3)
X   where f(alpha) is a quadratic interpolation function of the angle alpha
X   between v0-v1 and v1-v2.  This 'f' has certain nice properties:
X	To be exact, the magic numbers for interpolating stiffness of bends:
X   chosen so that (assuming midpoints of two adjacent segments are
X   normalized to be a unit distance apart ):
X	two collinear segments yield cylinder
X   	two at right angles yield a circular arc of 90 degrees 
X	two parallel yield close to a circular arc of 180 degrees
X
X   	A quadratic polynomial f(x) = a*x*x + b*x + c is interpolated 
X   through these constraints. 
X*/
Xstatic double a = .437902833, b = -.1045695, c = .3333333 ;
X 
Xmake_cylinder(v0,v1,v2, tt0, tt1)
Xpoint4 v0, v1, v2;
Xdouble tt0, tt1;
X{
X    register int i,j,k;
X    point4 m1, m2, m3, t1, t2, p1, p2;
X    point4 ty_vec0, tz_vec0, ty_vec1;
X    plane4 pl;
X    matrix4 mat1, mat2, mat3, mat1n, mat2n, mat3n, tmat;
X    double scale, angle, fraction, sinx, cosx, d;
X
X    /* the cylinder will start and end at the midpoints of these segments */
X    MIDPOINT4(v0,v1,m1)
X    MIDPOINT4(v1,v2,m2)
X
X    /* compute tform which takes canonical frame to one at m1 */
X    make_tform(m1,v1,v2,mat1);
X
X    /* compute amount of twist needed to line up seam with previous segment */
X    if (firsttime == 1)
X	{
X	firsttime = 0;
X	totaltwist = 0.0;
X	}
X    else
X	{
X	matvecmul4(lastmat, y_vec, ty_vec0);
X	matvecmul4(lastmat, z_vec, tz_vec0);
X	VSUB3(ty_vec0, m1, ty_vec0)
X	VSUB3(tz_vec0, m1, tz_vec0)
X	matvecmul4(mat1, y_vec, ty_vec1);
X	VSUB3(ty_vec1, m1, ty_vec1)
X	cosx = VDOT3(ty_vec1, ty_vec0);
X	sinx = VDOT3(ty_vec1, tz_vec0);
X	if (sinx > 1.0) sinx = 1.0; 
X	if (cosx > 1.0) cosx = 1.0; 
X	if (sinx < -1.0) sinx = -1.0; 
X	if (cosx < -1.0) cosx = -1.0; 
X	angle = atan2(sinx, cosx);
X    /*printf("angle is %f\n",angle); */
X	totaltwist += angle;
X	}
X		
X    /* compute rotation in y-z plane to match profile curve at seams */
X    axis_rotation(0, totaltwist, tmat);
X
X    matmatmul4(mat1, tmat, mat1n);
X
X    /* now the hard part: compute the tform which sends the profile curve
X 	onto the plane through v1 which is perpendicular to the plane
X 	of v0-v1-v2 and which bisects the angle v0-v1-v2.  This is the
X	curve which is the intersection of the two linear cylinders
X 	around the two segments v0-v1 and v1-v2 */ 
X    VSUB3(v0,v1,p1)		
X    VSUB3(v2,v1,p2)
X    d = MAGNITUDE3(p2);
X    normalize3(p1);
X    normalize3(p2);
X   
X    VSUB3(p2,p1,t2)
X    VSCALE3(t2, d, t2)
X    VADD3(v1,t2,t2)
X    make_tform(v1,t2,v2,mat3);	
X
X    /* compute the corresponding tform for the ending profile curve */
X    make_tform(m2,v2,v0,mat2);
X    matmatmul4(mat2, tmat, mat2n);
X    copymat4(mat2, lastmat);
X
X    /* correct the scaling in the y-direction for middle curve to reflect stretching */
X    angle = VDOT3(p1,p2);
X    if (angle > 1.0) angle = 1.0;
X    if (angle < 1.0) angle = -1.0;
X    angle = acos(angle);
X    /* compute scaling factor for profile curve in direction of stretching */
X    scale = sin(angle/2.0);
X    if (scale != 0.0) 
X	for (i=0; i<3; ++i)	tmat[1][i] /= scale; 
X
X    matmatmul4(mat3, tmat, mat3n);	/* follow with the twist rotation */
X
X    for (i=0; i<profile_size; ++i)	
X	matvecmul4(mat1n, profile_init[i], profile0[i]);
X    for (i=0; i<profile_size; ++i)	
X	matvecmul4(mat2n, profile_init[i], profile3[i]);
X    for (i=0; i<profile_size; ++i)	
X	matvecmul4(mat3n, profile_init[i], profilem[i]);
X   
X    /* now compute the intermediate profile curves: f(alpha) from the end profiles to
X       the middle one. 
X       Compute how far out the segment to go to place intermediate control points */
X    angle = angle/3.1415926;
X    angle = 1.0-angle;				
X    fraction = a*angle*angle + b * angle + c;	/* our quadratic interpolant */
X    fraction = (2*fraction*scale); 
X    /* Do the linear interpolation between profile curves */
X    for (i=0; i<profile_size; ++i)	
X	{
X	VSUB3(profilem[i], profile0[i], t1)
X 	VSCALE3(t1, fraction, t1)
X	VADD3(t1, profile0[i], profile1[i])
X	}
X    for (i=0; i<profile_size; ++i)	
X	{
X	VSUB3(profilem[i], profile3[i], t1)
X 	VSCALE3(t1, fraction, t1)
X	VADD3(t1, profile3[i], profile2[i])
X	}
X      
X    /* now output the patches */
X    if (softimage)
X	{
X        for (k = 0, i=0; k<profile_size; k++)
X	  {
X          if (k%4 != 3) 
X	    fprintf(outf,"     %d\t%f\t%f\t%f\t1.0\n",
X		i++,profile0[k][0],profile0[k][1], profile0[k][2]);
X	  }
X	  fprintf(outf,"\n");
X        for (k = 0, i=0; k<profile_size; k++)
X	  {
X          if (k%4 != 3) 
X	    fprintf(outf,"     %d\t%f\t%f\t%f\t1.0\n",
X		i++,profile1[k][0],profile1[k][1], profile1[k][2]);
X	  }
X	  fprintf(outf,"\n");
X        for (k = 0, i=0; k<profile_size; k++)
X	  {
X          if (k%4 != 3) 
X	    fprintf(outf,"     %d\t%f\t%f\t%f\t1.0\n",
X		i++,profile2[k][0],profile2[k][1], profile2[k][2]);
X	  }
X	  fprintf(outf,"\n");
X	}
X    else if (!make_binary_pixar)
X    {
X    for (i=0; i<profile_size; i += 4)
X	{
X	for (j=0; j<4; ++j) fprintf(outf,"%f\t%f\t%f\n",
X		profile0[i+j][0],profile0[i+j][1], profile0[i+j][2]);
X	fprintf(outf,"\n");
X	for (j=0; j<4; ++j) fprintf(outf,"%f\t%f\t%f\n",
X		profile1[i+j][0],profile1[i+j][1], profile1[i+j][2]);
X	fprintf(outf,"\n");
X	for (j=0; j<4; ++j) fprintf(outf,"%f\t%f\t%f\n",
X		profile2[i+j][0],profile2[i+j][1], profile2[i+j][2]);
X	fprintf(outf,"\n");
X	for (j=0; j<4; ++j) fprintf(outf,"%f\t%f\t%f\n",
X		profile3[i+j][0],profile3[i+j][1], profile3[i+j][2]);
X	fprintf(outf,"\n");
X	{
X	double s0,s1,dt0,dt1;	
X	s0 = .25+(i)/((double)2*profile_size) ;
X	s1 = .25+(i+4)/((double)2*profile_size) ;
X	if (totalon)
X	    {
X	    dt0 = total * tt0;
X	    dt1 = total * tt1;
X	    dt0 = dt0 - ((int) (dt0));	/* fraction from 0 -> 2 */
X	    dt1 = dt1 - ((int) (dt1));
X	    if (dt1 < dt0) dt1 += 1.0;
X	    dt0 = .01 + dt0/2.0;
X	    dt1 = .01 + dt1/2.0;
X	    }
X	else
X	    {
X	    dt0 = 0.25;
X	    dt1 = .75;
X	    }
X	fprintf(outf,"%f %f %f %f %f %f %f %f \n",
X		s0,dt0,s1,dt0,s0,dt1,s1,dt1);
X	}
X	fprintf(outf,"\n");
X	}
X    }
X    else
X    {
X	if (vptr + (profile_size/4)*56*sizeof(float) > BlockSize*array_size)
X	    {
X	    if (debug)
X		{
X		printf("allocating again: size is now %d\n",array_size*BlockSize);
X	 	}
X	    array_size = array_size*2;
X            v = (float *) realloc(v,array_size*BlockSize*sizeof(float));
X	    }
X	
X    for (i=0; i<profile_size; i += 4)
X	{
X	for (j=0; j<4; ++j) 
X	    for (k=0;k<3;++k)	*(v+vptr++) = profile0[i+j][k];
X	for (j=0; j<4; ++j) 
X	    for (k=0;k<3;++k)	*(v+vptr++) = profile1[i+j][k];
X	for (j=0; j<4; ++j) 
X	    for (k=0;k<3;++k)	*(v+vptr++) = profile2[i+j][k];
X	for (j=0; j<4; ++j) 
X	    for (k=0;k<3;++k)	*(v+vptr++) = profile3[i+j][k];
X	{
X	double s0,s1,dt0,dt1;
X	s0 = .25+(i)/((double)2*profile_size) ;
X	s1 = .25+(i+4)/((double)2*profile_size) ;
X	if (totalon)
X	    {
X	    dt0 = total * tt0;
X	    dt1 = total * tt1;
X	    dt0 = dt0 - ((int) (dt0));	/* fraction from 0 -> 2 */
X	    dt1 = dt1 - ((int) (dt1));
X	    if (dt1 < dt0) dt1 += 1.0;
X	    dt0 = .01 + dt0/2.0;
X	    dt1 = .01 + dt1/2.0;
X	    }
X	else
X	    {
X	    dt0 = 0.25;
X	    dt1 = .75;
X	    }
X	*(v+vptr++) = s0;
X	*(v+vptr++) = dt0;
X	*(v+vptr++) = s1;
X	*(v+vptr++) = dt0;
X	*(v+vptr++) = s0;
X	*(v+vptr++) = dt1;
X	*(v+vptr++) = s1;
X	*(v+vptr++) = dt1;
X	}
X 	}	
X    }
X}
X
Xstatic aaSetSplit () {
X	if (aargv[1][1] == 'n')		/* turn on */
X		split = 1;
X	else split = 0;
X}
X
X
Xstatic aaSetWrap () {
X	if (aargv[1][1] == 'n')		/* turn on */
X		wrap = 1;
X	else wrap = 0;
X}
X
Xstatic aaSetBinary () {
X	if (aargv[1][1] == 'n')		/* turn on */
X		make_binary_pixar = 1;
X	else make_binary_pixar = 0;
X}
X
Xstatic aaSetProfileFile() {
X	strcpy(profile_file, aargv[1]);
X}
X
Xstatic aaSetRadius()
X{
X    radius = atof(aargv[1]);
X}
X
Xstatic aaSetTotal()
X{
X    total = atoi(aargv[1]);
X    totalon = 1;
X}
X
Xstatic aaSetTwist()
X{
X    twist = atoi(aargv[1]);
X    twiston = 1;
X}
X
Xstatic aaSetSoftimage()
X{
X    strcpy(softname,aargv[1]);
X    softimage = 1;
X}
X
Xstruct aarg_cmd cmd[] = {
X    AA_BEGIN
X	remain, "tuber %s",
X/*
X	"Read named inputfile  and put tubes around space curve(s) for pixar rendering.\nOutput file will be created as "inputfile.bbp[.fbi]"\nFirst line is 'LINK'\nNext line is # of components\nEach component consists of # of points, then {xyz} list of points\n"
X*/
X	"Read inputfile and...\n"
X    AA_NEXT
X	aaSetSoftimage, "[-soft %s]",
X	"Set softimage output, use given name to create object"
X    AA_NEXT
X	aaSetSplit, "[-s[plit] %s]",
X	"Set split-link-into-components flag (on/off)  [off]"
X    AA_NEXT
X	aaSetWrap, "[-w[rap] %s]",
X	"Set wraparound flag (on/off)  [off]"
X    AA_NEXT
X	aaSetBinary, "[-b[inary] %s]",
X	"Set binary Pixar output  (on/off) [off]"
X    AA_NEXT
X	aaSetRadius, "[-r[adius] %f]",
X	"Set radius of tubes (in scaled coordinates)[1.0]"
X    AA_NEXT
X	aaSetTotal, "[-t[otal] %d ]",
X	"Set total # of repeats of texture space (parametrized by arc length) "
X    AA_NEXT
X	aaSetTwist, "[-tw[ist] %d ]",
X	"Set total # of twists of profile curve (parametrized by arc length) "
X    AA_NEXT
X	aaSetProfileFile, "[-p[rofile] %s]",
X	"Give name of file to read profile curve (default is circular)"
X    AA_END,
X    AA_STDCMD,
X    AA_LASTCMD
X};
X
Xmain(argc,argv)
X	int argc;
X	char *argv[];
X{
X	aa_flag(cmd, "[-?", AA_UNSEEN);
X	aa_flag(cmd, "[-help", AA_UNSEEN);
X	aa_main(cmd,argc,argv);
X}
X
X
END_OF_FILE
if test 18815 -ne `wc -c <'tuber.c'`; then
    echo shar: \"'tuber.c'\" unpacked with wrong size!
fi
# end of 'tuber.c'
fi
if test -f 'tubeseg.c' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'tubeseg.c'\"
else
echo shar: Extracting \"'tubeseg.c'\" \(4237 characters\)
sed "s/^X//" >'tubeseg.c' <<'END_OF_FILE'
X#include "vec4.h"
X#include "gobject.h"
X#include <stdio.h>
X#include <aarg.h>
X#include <math.h>
X#define LSIZE 8192
X#define MAXPOINTS 64
X#define CYLFILE	"/n/poincare/usr1/gunn/data/cyl.bbp"
X#define SPHFILE	"/n/poincare/usr1/gunn/data/sphere.bbp"
X
X	float *v;		/* pointer to list of quadrilateral data */
X	int vptr, BlockSize, array_size;
X	int debug = 1;
X
X
X	point4 sphere[16*8], cylinder[16*4], tsph[16*8], tcyl[16*4];;
X
X 	point4 *vlist;
X	double radius = 1.0, scale = 1.0;
X 	matrix4 mat0, mat1;
X        char outfile[64] = {""}, profile_file[64] = {""};
X	int vsize, profile_size;
X	int make_binary_pixar = 0;
X	int firsttime;
X
Xremain(argc, argv)	char **argv;
X{
X    
X    int numseg, j, i, numvert, i0, i1, k;
X    double atof(), angle, d;
X    point4 pbuffer[MAXPOINTS], *profileptr, tv, uv;
X    FILE *fd=NULL, *fopen(), *fp;
X    
X    if ( argc < 2)
X	    exit(fprintf(stderr,"Please specify input file\n"));
X
X    fd = fopen(argv[1],"r");
X    if (fd == NULL)	
X            exit(fprintf(stderr, "Unable to open file, %s\n", argv[1]));
X
X    /* read in sphere, cylinder */
X    if ( (fp=fopen(CYLFILE, "r")) == NULL) exit(fprintf(stderr,"Unable to open cylinder file\n"));
X    for (i=0; i<16*4; ++i)
X	{
X	fscanf(fp,"%F%F%F",cylinder[i]+2,cylinder[i]+1,cylinder[i]);
X	for (j=1; j<3;++j)	cylinder[i][j] *= 2;
X	cylinder[i][3] = 1.0;
X	}
X    fclose(fp);
X	
X    if ( (fp=fopen(SPHFILE, "r")) == NULL) exit(fprintf(stderr,"Unable to open sphere file\n"));
X    for (i=0; i<16*8; ++i)
X	{
X	fscanf(fp,"%F%F%F",sphere[i]+2,sphere[i]+1,sphere[i]);
X	for (j=0; j<3;++j)	sphere[i][j] *= 2;
X	sphere[i][3] = 1.0;
X	}
X    fclose(fp);
X	
X    vsize = 1;
X    vlist = (point4 *) malloc(LSIZE * sizeof(point4));
X
X    /* get list of vertices */
X    fscanf(fd,"%d",&numvert);
X    for (i=0; i<numvert ;++i) 	
X	{
X	fscanf(fd,"%F%F%F",vlist[i]+2, vlist[i]+1, vlist[i]);
X	for (j=0;j<3;++j)	vlist[i][j] *= scale;
X	vlist[i][3] = 1.0;
X	}
X
X    if (!make_binary_pixar) printf("BBP\n");	/* ascii header for Pixar batch files */
X    else
X	    {
X	    vptr = 0;
X	    array_size = 1;
X	    BlockSize = 8192;
X	    v = (float *) malloc (sizeof(float) * BlockSize );
X	    }
X    /* output the spheres */
X    for (k=0; k<numvert ;++k) 	
X	{
X	makeid4(mat0);
X	mat0[0][0] = mat0[1][1] = mat0[2][2] = radius;
X	mat0[0][3] = vlist[k][0];
X	mat0[1][3] = vlist[k][1];
X	mat0[2][3] = vlist[k][2];
X	for (i=0; i<8; ++i)
X	  {
X	  for (j=0; j<16; ++j)
X	    {
X	    matvecmul4(mat0, sphere[16*i+j], tsph[16*i+j]);
X	    if (!make_binary_pixar)
X	   	printf("%f %f %f \n",tsph[16*i+j][2], tsph[16*i+j][1], tsph[16*i+j][0]);
X	    else printf("Not implemented\n");
X	    }
X	  if (!make_binary_pixar)	printf("\n");
X	  }
X	if (!make_binary_pixar)	printf("\n");
X	}
X
X    /* process the edges */
X    while (fscanf(fd,"%d%d",&i0, &i1) != EOF)
X	{
X	VSUB3(vlist[i1], vlist[i0], tv)	
X	d = MAGNITUDE3	(tv);
X	/* first, stretch standard cyl */
X	makeid4(mat0);
X	mat0[2][2] = mat0[1][1] = radius;
X	mat0[0][0] = d;
X	/* then, rotate and translate the result into position */
X	/* this next line is a kludge */
X	VADD3(vlist[i0], vlist[i1], uv)
X	make_tform(vlist[i0], vlist[i1], uv, mat1);
X	matmatmul4(mat1, mat0, mat1);
X	for (i=0; i<4; ++i)
X	  {
X	  for (j=0; j< 16; ++j)
X	    {
X	    matvecmul4(mat1, cylinder[16*i+j], tcyl[16*i+j]);
X	    if (!make_binary_pixar)
X	   	printf("%f %f %f \n",tcyl[16*i+j][2], tcyl[16*i+j][1], tcyl[16*i+j][0]);
X	    else printf("Not implemented\n");
X	    }
X	  if (!make_binary_pixar)	printf("\n");
X	  }
X	if (!make_binary_pixar)	printf("\n");
X	}
X}	
Xstatic aaSetScale()
X{
X    scale = atof(aargv[1]);
X}
X
Xstatic aaSetRadius()
X{
X    radius = atof(aargv[1]);
X}
X
Xstruct aarg_cmd cmd[] = {
X    AA_BEGIN
X	remain, "tubeseg %s",
X	"Put tubes around a list of line segments for pixar rendering.\n\tFirst line of file contains n, # of vertices.\n\tThen follow the list of vertices:\n\tFinally, pairs of integers specifying the edges, indices into vertex list"
X    AA_NEXT
X	aaSetScale, "[-s[cale] %f]",
X	"Set scaling factor for vertex coordinates [1.0]"
X    AA_NEXT
X	aaSetRadius, "[-r[adius] %f]",
X	"Set radius of tubes (in scaled coordinates)[1.0]"
X    AA_END,
X    AA_STDCMD,
X    AA_LASTCMD
X};
X
Xmain(argc,argv)
X	int argc;
X	char *argv[];
X{
X	aa_flag(cmd, "[-?", AA_UNSEEN);
X	aa_flag(cmd, "[-help", AA_UNSEEN);
X	aa_main(cmd,argc,argv);
X}
END_OF_FILE
if test 4237 -ne `wc -c <'tubeseg.c'`; then
    echo shar: \"'tubeseg.c'\" unpacked with wrong size!
fi
# end of 'tubeseg.c'
fi
if test -f 'vec4.h' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'vec4.h'\"
else
echo shar: Extracting \"'vec4.h'\" \(2852 characters\)
sed "s/^X//" >'vec4.h' <<'END_OF_FILE'
X#include <math.h>
Xtypedef	double  point4[4];
Xtypedef	double  plane4[4];
Xtypedef	double  matrix4[4][4];
X
X#define	ABS(aaa)	(( (aaa) > 0) ? (aaa) : -(aaa))
X#define	VDOT3(v0,v1)	(v0[0]*v1[0] + v0[1]*v1[1] + v0[2]*v1[2])
X#define	INPRO31(v0,v1)	(v0[0]*v1[0] + v0[1]*v1[1] + v0[2]*v1[2] - v0[3]*v1[3])
X#define	INPRO4(v0,v1)	(v0[0]*v1[0] + v0[1]*v1[1] + v0[2]*v1[2] + v0[3]*v1[3])
X
X#define MAGNITUDE3(v0)	(sqrt (VDOT3(v0, v0)))
X#define MAGNITUDE31(v0)	(sqrt (ABS (INPRO31(v0, v0))))
X#define MAGNITUDE4(v0)	(sqrt (INPRO4(v0, v0)))
X
X#define PROJECT4(v0)	\
X	{float d = v0[3];	\
X	if (d)	{	\
X	d = 1.0/d;	\
X	v0[0] = v0[0] * d;	\
X	v0[1] = v0[1] * d;	\
X	v0[2] = v0[2] * d;	\
X	v0[3] = 1.0; }}	
X
X#define NORMALIZE3(v0)	\
X	{double xxx;		\
X	if (xxx = MAGNITUDE3(v0))	\
X	{v0[0] /= xxx;  v0[1] /= xxx;  v0[2] /= xxx; }  }
X
X#define 	NORMALIZE31(v)	\
X 	{double tmp = INPRO31(v,v);		\
X	if (tmp) { tmp = 1.0/sqrt(ABS(tmp)); v[0] *= tmp; v[1] *= tmp; v[2] *= tmp; v[3] *= tmp; }}
X
X#define VPRINT3(name,v)		printf("name\t%f %f %f\n",v[0],v[1],v[2]);
X#define VPRINT4(name,v)		printf("name\t%f %f %f %f %f\n",v[0],v[1],v[2],v[3]);
X#define VADD3(v0,v1,v2)	\
X	{v2[0] = v0[0] + v1[0]; \
X	v2[1] = v0[1] + v1[1]; \
X	v2[2] = v0[2] + v1[2]; \
X	v2[3] = 1.0; }
X	
X#define VADD4(v0,v1,v2)	\
X	{v2[0] = v0[0] + v1[0]; \
X	v2[1] = v0[1] + v1[1]; \
X	v2[2] = v0[2] + v1[2]; \
X	v2[3] = v0[3] + v1[3]; }
X
X#define VSUB3(v0,v1,v2)	\
X	{v2[0] = v0[0] - v1[0]; \
X	v2[1] = v0[1] - v1[1]; \
X	v2[2] = v0[2] - v1[2]; \
X	v2[3] = 1.0;}
X	
X#define VSUB4(v0,v1,v2)	\
X	{v2[0] = v0[0] - v1[0]; \
X	v2[1] = v0[1] - v1[1]; \
X	v2[2] = v0[2] - v1[2]; \
X	v2[3] = v0[3] - v1[3];}
X
X#define MIDPOINT4(v0,v1,v2)	\
X	{v2[0] = (v0[0] + v1[0])/2.0; \
X	v2[1] = (v0[1] + v1[1])/2.0; \
X	v2[2] = (v0[2] + v1[2])/2.0; \
X	v2[3] = (v0[3] + v1[3])/2.0; }
X
X#define VISZERO4(v)	(v[0] == 0 && v[1] == 0 && v[2] == 0 && v[3] == 0)
X
X#define MIDPOINT31(v0,v1,v2)	\
X	{float m00, m11, a,v3[4];	\
X	VSUB4(v0, v1, v3);	\
X	if (VISZERO4(v3))	\
X		VCOPY4(v0,v2)	\
X	else	{		\
X	m00 =  INPRO31(v0,v3);	\
X	m11 = INPRO31(v1,v3);	\
X	v2[0] = m11*v0[0] - m00*v1[0]; \
X	v2[1] = m11*v0[1] - m00*v1[1]; \
X	v2[2] = m11*v0[2] - m00*v1[2]; \
X	v2[3] = m11*v0[3] - m00*v1[3]; }}
X
X#define XPRO3(v0,v1,v2)	\
X	{v2[0] = v0[1]*v1[2] - v0[2]*v1[1];	\
X	v2[1] = v0[2]*v1[0] - v0[0]*v1[2];	\
X	v2[2] = v0[0]*v1[1] - v0[1]*v1[0];  }
X
X#define VSCALE3(v0, s, v1)	\
X	{v1[0] = s * v0[0];	\
X	v1[1] = s * v0[1];	\
X	v1[2] = s * v0[2];	}
X
X#define VSCALE4(v0, s, v1)	\
X	{v1[0] = s * v0[0];	\
X	v1[1] = s * v0[1];	\
X	v1[2] = s * v0[2];	\
X	v1[3] = s * v0[3];	}
X
X#define VCOPY3(v0, v1)	\
X	{v1[0] = v0[0];	\
X	v1[1] = v0[1];	\
X	v1[2] = v0[2];	}
X
X#define VCOPY4(v0, v1)	\
X	{v1[0] = v0[0];	\
X	v1[1] = v0[1];	\
X	v1[2] = v0[2];	\
X	v1[3] = v0[3];}
X
X#define VZERO3(v0)	\
X	{v0[0] = 0;	\
X	v0[1] = 0;	\
X	v0[2] = 0;	\
X	v0[3] = 1.0; }
X
X#define COPYMAT4(m1, m2)	\
X	bcopy(m1, m2, (sizeof(double) * 16));
X
X
X
X
X
END_OF_FILE
if test 2852 -ne `wc -c <'vec4.h'`; then
    echo shar: \"'vec4.h'\" unpacked with wrong size!
fi
# end of 'vec4.h'
fi
if test -f 'xform.c' -a "${1}" != "-c" ; then 
  echo shar: Will not clobber existing file \"'xform.c'\"
else
echo shar: Extracting \"'xform.c'\" \(628 characters\)
sed "s/^X//" >'xform.c' <<'END_OF_FILE'
X#include <math.h>
X
Xmain(argc, argv)	char **argv;
X{ 
X	int numframes, strlength;
X	double scale, c, s, angle;
X	char *name;
X 	register int i;
X
X	if (argc != 4)  exit(printf("Usage: %s numframes scale name\n",argv[0]));
X	numframes = atoi(*++argv);
X	scale= atof(*++argv);
X	name = *++argv;
X	strlength = strlen(name);
X
X	for (i=0;i<numframes;++i)
X		{
X		sprintf(name+strlength,".%d.pic\0",i);	
X		angle = (i*6.28315)/numframes;
X		c = cos(angle);  s = sin(angle);
X		printf("modelx sph %f 0 %f 0 0 %f 0 0 %f 0 %f 0 0 0 0 1\n",
X			c*scale, s*scale, scale, -s*scale, c*scale);
X		printf("doframe -sync\n");
X		printf("save %s\n",name);
X		} 
X		
X}
END_OF_FILE
if test 628 -ne `wc -c <'xform.c'`; then
    echo shar: \"'xform.c'\" unpacked with wrong size!
fi
# end of 'xform.c'
fi
echo shar: End of shell archive.
exit 0
