/* minimal ray tracer, hybrid version - 888 tokens
 * Paul Heckbert, ucbvax!pixar!ph, 13 Jun 87
 * Using tricks from Darwyn Peachey and Joe Cychosz. */
/* cc minray.c -lc_s -s -lgl_s -lm -o minray */

#define TOL 1e-7
#define UAMBIENT vec U, black, amb
#define SPHERE struct sphere {vec cen, color; double rad, kd, ks, kt, kl, ir;} \
    *s, *best, sph[]
typedef struct { double x,y,z; } vec;
#include "ray.h"
#include <stdio.h>
#include <gl/gl.h>
#include <gl/device.h> 

int yx;
double u, b, tmin, sqrt(), tan();

double vdot(A, B)
vec A, B;
{
    return A.x*B.x + A.y*B.y + A.z*B.z;
}

vec vcomb(a, A, B)	/* aA+B */
double a;
vec A, B;
{
    B.x += a*A.x;
    B.y += a*A.y;
    B.z += a*A.z;
    return B;
}

vec vunit(A)
vec A;
{
    return vcomb(1./sqrt(vdot(A, A)), A, black);
}

struct sphere *intersect(P, D)
vec P, D;
{
    best = 0;
    tmin = 1e30;
    s = sph+NSPHERE;
    while (s-->sph)
	b = vdot(D, U = vcomb(-1., P, s->cen)),
	u = b*b-vdot(U, U)+s->rad*s->rad,
	u = u>0 ? sqrt(u) : 1e31,
	u = b-u>TOL ? b-u : b+u,
	tmin = u>=TOL && u<tmin ?
	    best = s, u : tmin;
    return best;
}

vec trace(level, P, D)
vec P, D;
{
    double d, eta, e;
    vec N, color;
    struct sphere *s, *l;

    if (!level--) return black;
    if ((s = intersect(P, D))!=0);
    else return amb;

    color = amb;
    eta = s->ir;
    d = -vdot(D, N = vunit(vcomb(-1., P = vcomb(tmin, D, P), s->cen)));
    if (d<0)
	N = vcomb(-1., N, black),
	eta = 1/eta,
	d = -d;
    l = sph+NSPHERE;
    while (l-->sph)
	if ((e = l->kl*vdot(N, U = vunit(vcomb(-1., P, l->cen)))) > 0 &&
	    intersect(P, U)==l)
		color = vcomb(e, l->color, color);
    U = s->color;
    color.x *= U.x;
    color.y *= U.y;
    color.z *= U.z;
    e = 1-eta*eta*(1-d*d);
    /* the following is non-portable: we assume right to left arg evaluation.
     * (use U before call to trace, which modifies U) */
    return vcomb(s->kt,
	e>0 ? trace(level, P, vcomb(eta, D, vcomb(eta*d-sqrt(e), N, black)))
	    : black,
	vcomb(s->ks, trace(level, P, vcomb(2*d, N, D)),
	    vcomb(s->kd, color, vcomb(s->kl, U, black))));
}

void uspset(int x,int y,int r,int g,int b)
{
   short linecolvec[3];
   short xy1[2],xy2[2];

   linecolvec[0]=r;
   linecolvec[1]=g;
   linecolvec[2]=b;
   c3s(linecolvec);
   xy1[0]=x;xy1[1]=y;xy2[0]=x;xy2[1]=y;
   bgnline();
      v2s(xy1);v2s(xy2);
   endline();
}

void gend()
{
    short dat;
    while ( 1 ) {
       if ( qread( &dat ) == WINQUIT ) exit(0);
    }
}

void setscreen(int nnx,int mmy)
{
   prefsize(nnx,mmy);
   winopen("draw bezier");
   color(BLACK);
   clear();
   color(WHITE);
}

void main()
{
   int x,y;

/*   printf("%d %d \n", SIZE, SIZE); */
   setscreen(SIZE,SIZE);
   RGBmode();
   gconfig();
   clear(); 

   yx=0;
   for(y=0;y<=SIZE;y++){
     for(x=0;x<SIZE;x++){
	U.x = yx%SIZE-SIZE/2,
	U.z = SIZE/2-yx++/SIZE,
	U.y = SIZE/2/tan(AOV/114.5915590261),	/* 360/PI~=114 */
	U = vcomb(255., trace(DEPTH, black, vunit(U)), black),
        /*printf("%.0f %.0f %.0f \n", U); */ 	/* yowsa! non-portable! */
        /* printf("%f %f %f \n", U); */
        uspset(x,SIZE-y,(int)U.x,(int)U.y,(int)U.z);
      }
   } 
   gend();
   gexit();
   exit(0); 
}