/*
 * Header Begins--Do not remove################################################
 *
 *	C-code to exercise the plotting package by animating a teapot.
 *
 *	(c) Eugene Fiume, University of Toronto, 1995.
 *
 *	You are permitted to use this code for noncommercial, educational
 *	purposes only.  This header must not be removed from this file.
 *
 * Header Ends--Do not remove##################################################
 */

#include <stdio.h>
#include "Plot.h"
#define PI	3.141592654
#define PI4	(PI/4.0)
#define PI3	(PI/3.0)
#define TWOPI	(2.0*PI)
#define FRAMES	30
#define copybits(s,d,n)	memcpy( (d), (s), (n))

typedef double matrix[4][4];
typedef double vector[4];

matrix xmatrix, ymatrix, zmatrix;
matrix *xmatrices, *ymatrices, *zmatrices;


void GetMatrix(mat)
	matrix mat;
{
	int i,j;

	for (i=0; i<4; i++) {
	    for (j=0; j<4; j++) {
		scanf("%lf", &(mat[i][j]));
	    }
	}

} /* GetMatrix */

void PutMatrix(mat)
	matrix mat;
{
	int i,j;

	for (i=0; i<4; i++)
		printf("%f %f %f %f\n",mat[i][0],mat[i][1],mat[i][2],mat[i][3]);
	putchar('\n');

} /* PutMatrix */

/*
 *	Compute [s^3  s^2  s  1]M[t^3 t^2 t 1]+, where M is
 *	BPB+ for basis B and observations P.
 *	s,t are specific values in [0,1].
 */
double poly(s, t, mat)
	double s, t;
	matrix mat;
{
	vector T, S, sums;
	double sum;
	int i,j;

	T[3] = 1.0;
	T[2] = t;
	T[1] = t*t;
	T[0] = T[1]*t;

	S[3] = 1.0;
	S[2] = s;
	S[1] = s*s;
	S[0] = S[1]*s;

	for (i=0; i<4; i++) {			/* first evaluate MT+ */
		sums[i] = 0.0;
	    	for (j=0; j<4; j++)
			sums[i] += T[j]*mat[j][i];
	}

	sum = 0.0;
	for (i=0; i<4; i++)			/* now compute S*sums */
		sum += S[i]*sums[i];

	return(sum);
} /* poly */


double polyx(s,t)
	double s,t;
{
	return(poly(s,t,xmatrix));
} /* polyx */

double polyy(s,t)
	double s,t;
{
	return(poly(s,t,ymatrix));
} /* polyy */

double polyz(s,t)
	double s,t;
{
	return(poly(s,t,zmatrix));
} /* polyz */



main(argc, argv)
	int argc;
	char *argv[];
{
	int i,j, n;
	double angle, dangle;

	printf("Animating 3D Parametric Polynomials\n");
	XOpen(600,600, argc, argv);

	printf("Hit a mouse button to start\n");
	GetEvent();

	scanf("%d", &n);			/* get number of patches */

	xmatrices = (matrix *) malloc(n*sizeof(matrix));	/* allocate */
	ymatrices = (matrix *) malloc(n*sizeof(matrix));	/* patch */
	zmatrices = (matrix *) malloc(n*sizeof(matrix));	/* data  */

	printf("Plotting %d spline segments for each of x,y,z\n", n);

/*
 *	Now get patch data.
 */
	for (i=0; i<n; i++) {
		GetMatrix(xmatrix);
		GetMatrix(ymatrix);
		GetMatrix(zmatrix);
		copybits(xmatrix,xmatrices[i],sizeof(xmatrix));
		copybits(ymatrix,ymatrices[i],sizeof(ymatrix));
		copybits(zmatrix,zmatrices[i],sizeof(zmatrix));
	}

/*
 *	Try to "animate" the data by successively rotating it.
 */
	dangle = TWOPI/FRAMES;
	angle  = 0.0;
	SetScale(75.0,-75.0);
	SetOffset(300,300);
	SetColour(5);

	for (angle=0,j=0; j<=FRAMES; j++,angle += dangle) {
	    ClearScreen();
	    SetRotation(0.0, 0.0, angle);
	    for (i=0; i<n; i++) {
		copybits(xmatrices[i],xmatrix,sizeof(xmatrix));
		copybits(ymatrices[i],ymatrix,sizeof(ymatrix));
		copybits(zmatrices[i],zmatrix,sizeof(zmatrix));
		Plot3D(polyx,polyy,polyz, 0.0,1.0, 8, 0.0,1.0, 8);
	    }
	}

	printf("Hit a mouse button to end\n");
	GetEvent();
} /* main */