/*  Elias' Computed Pictures  *
 *  By Elias Broms 1995       */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#include "../../tools/vector/linear.h"
#include "../../tools/ppm/ppm.h"
#include "ecp.h"


/*             *
 *    Data     *
 *             */

/* Amount of ambient light for all objects */
scalar globalambient;

int xlights, xspheres, xwalls, xtubes, xpipes, xcones;

light *lights;
sphere *spheres;
wall *walls;
tube *tubes;
pipe *pipes;
cone *cones;

intersection *sphereis, *wallis, *tubeis, *pipeis, *coneis;

/* This is where we will paint */
intcolor **map;


/*                  *
 *    Functions     *
 *                  */

/* Unused function, but it must be defined! (see linear.c) */
void drawline(int px, int py, int qx, int qy)
{
}

/* Scalar square function */
scalar sqr(scalar x)
{
  return (scalar)(x*x);
}

/* Returns the smallest solution to the quadratic ax^2+bx+c=0. If there
   is no solution, 0.0 is returned. */
scalar quadratic(scalar a, scalar b, scalar c)
{
  scalar fourac, sqrb, inv2a, alfa, beta,t;

  fourac=4.0*a*c;
  sqrb=sqr(b);
  if (fourac>sqrb)
    return 0.0;
  inv2a=0.5/a;
  alfa=-inv2a*b;
  beta=inv2a*sqrt(sqrb-fourac);
  t=min(alfa-beta, alfa+beta);
  return t;
}

/* Returns the closest intersection of the line li with the sphere sph. */
intersection sphintersect(line li, sphere sph)
{
  scalar t;
  vector r;
  intersection is;

  is.side=NOINTERSECTION;

  r=vsub(li.p, sph.c);
  t=quadratic(vdotv(li.k), 2.0*dot(li.k, r), vdotv(r)-sqr(sph.r));

  if (t>EPS)
  {
    is.p=vadd(smul(t, li.k), li.p);
    is.side=0;
  }
  return is;
}

/* Returns the closest intersection of the line li with the wall wa. */
intersection wallintersect(line li, wall wa)
{
  vector w1, w2, u, v;
  scalar t;
  intersection is;

  is.side=NOINTERSECTION;

  is.p=planeintersect(li, wa.pl);

  if (is.p.y>EPS)
    is.side=0;

  return is;
}

/* Returns the closest intersection of the line li with the tube tu. */
intersection tubeintersect(line li, tube tu)
{
  vector w1, w2, u, v;
  scalar t, xqTxq;
  intersection is0, is1, is2;
  plane pl1, pl2;
  
  /* Intersection with side 0 */

  is0.side=NOINTERSECTION;
  u=vsub(li.p, tu.li.p);
  w1=vsub(li.k, smul(dot(li.k, tu.li.k)/vdotv(tu.li.k), tu.li.k));
  w2=smul(dot(u, tu.li.k)/vdotv(tu.li.k), tu.li.k);
  v=vsub(u, w2);
  t=quadratic(vdotv(w1), 2.0*dot(w1, v), vdotv(v)-sqr(tu.r));
  if (t>EPS)
    {
      scalar mTxq, mTm;
      
      is0.p=vadd(smul(t, li.k), li.p);
      mTxq=dot(tu.li.k, vsub(is0.p, tu.li.p));
      mTm=vdotv(tu.li.k);
      
      if (mTxq>0 && mTm-mTxq>0)
	is0.side=0;
    }
  
  /* Intersection with side 1 */
  
  is1.side=NOINTERSECTION;
  pl1.n=vneg(tu.li.k);
  pl1.p=tu.li.p;
  t=-dot(pl1.n, vsub(li.p, pl1.p))/dot(pl1.n, li.k);
  if (t>EPS)
  {
    is1.p=vadd(smul(t, li.k), li.p);
    if (vdotv(vsub(is1.p, tu.li.p)) < sqr(tu.r))
      is1.side=1;
  }

  /* Intersection with side 2 */
  
  is2.side=NOINTERSECTION;
  pl2.n=tu.li.k;
  pl2.p=vadd(tu.li.p, tu.li.k);
  t=-dot(pl2.n, vsub(li.p, pl2.p))/dot(pl2.n, li.k);
  if (t>EPS)
  {
    is2.p=vadd(smul(t, li.k), li.p);
    if (vdotv(vsub(vsub(is2.p, tu.li.p), tu.li.k)) < sqr(tu.r))
      is2.side=2;
  }
  
  /* Return the closest of is0, is1 and is2 */

  if (is0.side==NOINTERSECTION)
    if (is1.side==NOINTERSECTION)
      if (is2.side==NOINTERSECTION)
	return is0; /* Could be is1 or is2 too, it doesn't matter */
      else
	return is2;
    else
      if (is2.side==NOINTERSECTION)
	return is1;
      else
	if (vdotv(is1.p) < vdotv(is2.p))
	  return is1;
	else
	  return is2;
  else
    if (is1.side==NOINTERSECTION)
      if (is2.side==NOINTERSECTION)
	return is0;
      else
	if (vdotv(is0.p) < vdotv(is2.p))
	  return is0;
	else
	  return is2;
    else
      if (is2.side==NOINTERSECTION)
	if (vdotv(is0.p) < vdotv(is1.p))
	  return is0;
	else
	  return is1;
      else
	printf("Internal error in tubeintersect()\n");
}

/* Returns the closest intersection of the line li with the pipe pi. */
intersection pipeintersect(line li, pipe pi)
{
  vector l, p, m, q, a, b, u, v;
  scalar r, t;
  intersection is;

  is.side=NOINTERSECTION;

  l=li.k;
  p=li.p;
  m=pi.li.k;
  q=pi.li.p;
  r=pi.r;

  u=vsub(p, q);
  a=vsub(l, smul(dot(l, m)/vdotv(m), m));
  b=smul(dot(u, m)/vdotv(m), m);
  v=vsub(u, b);
  t=quadratic(vdotv(a), 2.0*dot(a, v), vdotv(v)-sqr(r));

  if (t>EPS)
  {
    scalar mTxq, mTm;

    /* Does it intersect between the ends? */
    is.p=vadd(smul(t, l), p);
    mTxq=dot(m, vsub(is.p, q));
    mTm=vdotv(m);

    if (mTxq>0 && mTm-mTxq>0)
      is.side=0;
  }
  return is;
}

/* Returns the closest intersection of the line li with the cone co. */
intersection coneintersect(line li, cone co)
{
  vector l, p, m, q, a, b, u, v;
  scalar r, r2, t, mTm, invmTm, r2dmTm2, lTm, uTm;
  intersection is0, is1;
  plane pl;
  scalar aa, bb, cc, fourac, sqrb, inv2a, alfa, beta, t1, t2;
  vector x1, x2;
  scalar mTx1q, mTx2q;

  l=li.k;
  p=li.p;
  m=co.li.k;
  q=co.li.p;
  r=co.r;

  /* Side 0 intersection */

  is0.side=NOINTERSECTION;
  r2=sqr(r);
  u=vsub(p, q);
  mTm=vdotv(m);
  invmTm=1.0/mTm;
  r2dmTm2=r2*sqr(invmTm);
  lTm=dot(l, m);
  uTm=dot(u, m);
  a=vsub(l, smul(lTm*invmTm, m));
  b=smul(uTm*invmTm, m);
  v=vsub(u, b);
  /* quadratic() can't be used because we need both solutions. */
  aa=vdotv(a) - r2dmTm2*sqr(lTm);
  bb=2.0*(dot(a, v) - r2dmTm2*lTm*uTm);
  cc=vdotv(v) - r2dmTm2*sqr(uTm);
  fourac=4.0*aa*cc;
  sqrb=sqr(bb);
  if (fourac<=sqrb) /* If real solutions */
  {
    inv2a=0.5/aa;
    alfa=-inv2a*bb;
    beta=inv2a*sqrt(sqrb-fourac);
    /* The solutions of the quadratic aa*t^2 + bb*t + cc = 0 */
    t1=alfa-beta;
    t2=alfa+beta;
    if (t1>EPS || t2>EPS)
      {
	/* The intersection points with the quadratic surface */
	x1=vadd(smul(t1, l), p);
	x2=vadd(smul(t2, l), p);
	mTx1q=dot(m, vsub(x1, q));
	mTx2q=dot(m, vsub(x2, q));
	
	/* Get the closest of x1 and x2 */
	is0.side=0;
	if (mTx1q>0 && mTm-mTx1q>0) /* If x1 is on the cone */
	  if (mTx2q>0 && mTm-mTx2q>0) /* If x2 is on the cone */
	    if (vdotv(x1) < vdotv(x2)) /* If x1 is closer than x2 */
	      is0.p=x1;
	    else
	      is0.p=x2;
	  else
	    is0.p=x1;
	else
	  if (mTx2q>0 && mTm-mTx2q>0) /* If x2 is on the cone */
	    is0.p=x2;
	  else
	    is0.side=NOINTERSECTION; /* Neither x1 nor x2 are on the cone */
      }
/*      else
      {
	x1=vadd(smul(t1, l), p);
	mTx1q=dot(m, vsub(x1, q));
	if (mTx1q>0 && mTm-mTx1q>0)
	{
	  is0.p=x1;
	  is0.side=1;
	}
      }
    else
      if (t2>EPS)
      {
	x2=vadd(smul(t2, l), p);
	mTx2q=dot(m, vsub(x2, q));
	if (mTx2q>0 && mTm-mTx2q>0)
	{
	  is0.p=x2;
	  is0.side=2;
	}
      }*/
  }
  
  /* Side 1 intersection */

  is1.side=NOINTERSECTION;
  pl.n=m;
  pl.p=vadd(q, m);
  t=-dot(pl.n, vsub(li.p, pl.p))/dot(pl.n, li.k);
  if (t>EPS)
  {
    is1.p=vadd(smul(t, li.k), li.p);
    if (vdotv(vsub(is1.p, pl.p)) < r2)
      is1.side=1;
  }

  /* Return the closest of is0 and is1 */

  if (is0.side==NOINTERSECTION)
  {
    if (is1.side==NOINTERSECTION)
      return is0; /* We could have is1 here, it doesn't matter. */
    else
      return is1;
  }
  else
  {
    if (is1.side==NOINTERSECTION)
      return is0;
    else
    {         /* Both side 0 and side 1 were intersected */
      if (vdotv(is0.p) < vdotv(is1.p))
	return is0;
      else
	return is1;
    }
  }
}

/* Computes the intensity in the intersection is of the object with objinfo oi.
   The vector n is the normal in is. If shadow=TRUE then only ambient 
   light is used.
*/
color compute_intensity(objinfo oi, vector n, intersection is, bool shadow)
{
    color c;
    scalar cosalfa,reflterm;
    vector r;

    if (!shadow)
    {
	/* Point light source's intensity times object's intensity */
	c=vassign(lights[0].col.x*oi.col[is.side].x,
		  lights[0].col.y*oi.col[is.side].y,
		  lights[0].col.z*oi.col[is.side].z);
	
	/* Diffuse-reflection */
	c=smul(oi.diff*cosangle(n, lights[0].dir), c);
    }
    else c=zero();
    
    /* Ambient light */
    c=vadd(smul(globalambient*oi.amb, oi.col[is.side]), c);
    
    if (!shadow)
    {
	/* Reflection vector R=2N(N*L)-L */
	r=vsub(smul(2.0*cosangle(n,lights[0].dir)/norm(n),n),
	       smul(1.0/norm(lights[0].dir),lights[0].dir));
	
	/* Angle between reflection of light source and view point */
	cosalfa=-dot(r,is.p)/(norm(r)*norm(is.p));
	if (cosalfa<0.0)cosalfa=0.0;
	
	/* Specular-reflection coefficient */
	reflterm=oi.fspre*exp(oi.refl*log(cosalfa));
	c=vadd(smul(reflterm,lights[0].col),c);
    }
    return c;
}

/* Returns TRUE if the line li intersects with any object other than the 
   object of type objtype with number obj. 
   Used for shadows. */
bool isshadow(line li, int objtype, int obj)
{
  bool shadow=FALSE;
  intersection shadowis;
  int i;

  for (i=0; i<xspheres; i++)
  {
    if (!(objtype==SPHERE && i==obj))
    {
      shadowis=sphintersect(li, spheres[i]);
      if (shadowis.side!=NOINTERSECTION)
	shadow=TRUE;
    }
  }
  for (i=0; i<xtubes; i++)
  {
    if (!(objtype==TUBE && i==obj))
    {
      shadowis=tubeintersect(li, tubes[i]);
      if (shadowis.side!=NOINTERSECTION)
	shadow=TRUE;
    }
  }
  for (i=0; i<xpipes; i++)
  {
    if (!(objtype==PIPE && i==obj))
    {
      shadowis=pipeintersect(li, pipes[i]);
      if (shadowis.side!=NOINTERSECTION)
	shadow=TRUE;
    }
  }
  for (i=0; i<xcones; i++)
  {
    if (!(objtype==CONE && i==obj))
    {
      shadowis=coneintersect(li, cones[i]);
      if (shadowis.side!=NOINTERSECTION)
	shadow=TRUE;
    }
  }
  return shadow;
}

vector spherenormal(sphere sph, intersection is)
{
  return vsub(is.p, sph.c);
}

vector wallnormal(wall wa, intersection is)
{
  return wa.pl.n;
}

vector tubenormal(tube tu, intersection is)
{
  vector w, m;
  scalar wTmdmTm;

  switch (is.side)
  {
  case 0:
    w=vsub(tu.li.p, is.p);
    m=tu.li.k;
    wTmdmTm=dot(w, m)/vdotv(m);
    return vsub(smul(wTmdmTm, m), w);
    break;
  case 1:
    return vneg(tu.li.k);
    break;
  case 2:
    return tu.li.k;
    break;
  }
  fprintf(stderr, "Internal error in tubenormal()\n");
}

vector pipenormal(pipe pi, intersection is)
{
  vector w, m;

  w=vsub(pi.li.p, is.p);
  m=pi.li.k;

  return vsub(smul(dot(w, m)/vdotv(m), m), w);
  fprintf(stderr, "Internal error in pipenormal()\n");
}

vector conenormal(cone co, intersection is)
{
  vector m, w;

  m=co.li.k;

  switch(is.side)
  {
  case 0:
    w=vsub(is.p, co.li.p);
    return vsub(w, smul(vdotv(w)/dot(w, m), m));
    break;
  case 1:
    return m;
    break;
  }
  fprintf(stderr, "Internal error in conenormal()\n");
}

/* Returns the color of a certain pixel. The vector v is supposed to be 
a point on the projection plane */
color pixel(vector v, bool shadowsoff)
{
  intersection objis[OBJECTTYPES];
  line projli, shadowli;
  vector normal;
  int i, csphere=0, cwall=0, tube=0, pipe=0, cone=0, objtype=0, obj;
  bool hitsphere=FALSE, hittube=FALSE, hitpipe=FALSE, hitcone=FALSE;
  bool shadow=FALSE;
  objinfo oi;

  /* The projection line through the origin */
  projli.k=vneg(v);
  projli.p=zero();

  /* Compute the line's closest intersection point with the spheres */
  for (i=0; i<xspheres; i++)
  {
    sphereis[i]=sphintersect(projli, spheres[i]);
    if (sphereis[i].side != NOINTERSECTION)
      csphere=i;
  }
  for (i=0; i<xspheres; i++)
    if (sphereis[i].side != NOINTERSECTION)
    {
      if (vdotv(sphereis[i].p) < vdotv(sphereis[csphere].p))
	csphere=i;
      hitsphere=TRUE;
    }

  /* Compute the line's closest intersection point with the tubes */
  for (i=0; i<xtubes; i++)
  {
    tubeis[i]=tubeintersect(projli, tubes[i]);
    if (tubeis[i].side != NOINTERSECTION)
      tube=i;
  }
  for (i=0; i<xtubes; i++)
    if (tubeis[i].side != NOINTERSECTION)
    {
      if (vdotv(tubeis[i].p) < vdotv(tubeis[tube].p))
	tube=i;
      hittube=TRUE;
    }

  /* Compute the line's closest intersection point with the pipes */
  for (i=0; i<xpipes; i++)
  {
    pipeis[i]=pipeintersect(projli, pipes[i]);
    if (pipeis[i].side != NOINTERSECTION)
      pipe=i;
  }
  for (i=0; i<xpipes; i++)
    if (pipeis[i].side != NOINTERSECTION)
    {
      if (vdotv(pipeis[i].p) < vdotv(pipeis[pipe].p))
	pipe=i;
      hitpipe=TRUE;
    }

  /* Compute the line's closest intersection point with the cones */
  for (i=0; i<xcones; i++)
  {
    coneis[i]=coneintersect(projli, cones[i]);
    if (coneis[i].side != NOINTERSECTION)
      cone=i;
  }
  for (i=0; i<xcones; i++)
    if (coneis[i].side != NOINTERSECTION)
    {
      if (vdotv(coneis[i].p) < vdotv(coneis[cone].p))
	cone=i;
      hitcone=TRUE;
    }

  /* Compute the lines' closest intersection point with the walls */
  for (i=0; i<xwalls; i++)
  {
      wallis[i]=wallintersect(projli, walls[i]);
      if (wallis[i].side != NOINTERSECTION)
	  cwall=i;
  }
  for (i=0; i<xwalls; i++)
    if (wallis[i].side != NOINTERSECTION)
    {
      if (vdotv(wallis[i].p)<vdotv(wallis[cwall].p))
	cwall=i;
    }  


  /* Now we have one (at most) intersection point for every object type.
     Which of them is closest to the viewer? */

  if (xspheres)
    objis[SPHERE]=sphereis[csphere];
  if (xwalls)
    objis[WALL]=wallis[cwall];
  if (xtubes)
    objis[TUBE]=tubeis[tube];
  if (xpipes)
    objis[PIPE]=pipeis[pipe];
  if (xcones)
    objis[CONE]=coneis[cone];

  if (xspheres && objis[SPHERE].side != NOINTERSECTION)
    objtype=SPHERE;
  else
  if (xwalls && objis[WALL].side !=NOINTERSECTION)
    objtype=WALL;
  else
  if (xtubes && objis[TUBE].side !=NOINTERSECTION)
    objtype=TUBE;
  else
  if (xpipes && objis[PIPE].side !=NOINTERSECTION)
    objtype=PIPE;
  else
  if (xcones && objis[CONE].side !=NOINTERSECTION)
    objtype=CONE;

  if (xspheres && objis[SPHERE].side != NOINTERSECTION)
    if (vdotv(objis[SPHERE].p) < vdotv(objis[objtype].p))
      objtype=SPHERE;
  if (xwalls && objis[WALL].side != NOINTERSECTION)
    if (vdotv(objis[WALL].p) < vdotv(objis[objtype].p))
      objtype=WALL;
  if (xtubes && objis[TUBE].side != NOINTERSECTION)
    if (vdotv(objis[TUBE].p) < vdotv(objis[objtype].p))
      objtype=TUBE;
  if (xpipes && objis[PIPE].side != NOINTERSECTION)
    if (vdotv(objis[PIPE].p) < vdotv(objis[objtype].p))
      objtype=PIPE;
  if (xcones && objis[CONE].side != NOINTERSECTION)
    if (vdotv(objis[CONE].p) < vdotv(objis[objtype].p))
      objtype=CONE;

  /* objtype now contains the type of the object whose closest object
     is closest among *all* objects, i. e. objis[objtype] is the
     intersection point that will be drawn in this pixel. */
  

  switch(objtype)
  {

  case SPHERE:
    oi=spheres[csphere].oi;
    obj=csphere;
    normal=spherenormal(spheres[csphere], objis[SPHERE]);
    break;

  case WALL:
    oi=walls[cwall].oi;
    obj=cwall;
    normal=wallnormal(walls[cwall], objis[WALL]);
    break;

  case TUBE:
    oi=tubes[tube].oi;
    obj=tube;
    normal=tubenormal(tubes[tube], objis[TUBE]);
    break;

  case PIPE:
    oi=pipes[pipe].oi;
    obj=pipe;
    normal=pipenormal(pipes[pipe], objis[PIPE]);
    break;

  case CONE:
    oi=cones[cone].oi;
    obj=cone;
    normal=conenormal(cones[cone], objis[CONE]);
    break;

  default:
    fprintf(stderr, "Internal error in pixel()\n");
    break;
  }

  /* Get the line from the intersection point to the light source */
  shadowli.p=objis[objtype].p;
  shadowli.k=lights[0].dir;

  if (shadowsoff)
    return compute_intensity(oi, normal, objis[objtype], FALSE);
  else
    return compute_intensity(oi, normal, objis[objtype], 
			     isshadow(shadowli, objtype, obj));
}

/* Calls pixel() for all pixels in the picture.
   If antialiasing is used, then pixel() is called 16 times per pixel and
   the pixel gets the mean value of the 16 colors. */
void render_the_whole_picture(bool antialias, bool preview, bool shadowsoff)
{
  scalar sx, sy, px, py, d=1.0;
  color c;
  
  if (preview)
    d=3.0;

  for (sy=0.0; sy<(scalar)HEIGHT-EPS; sy+=d)
  {
    if (((int)sy)%10 == 0)
      printf("Line %d\n",(int)sy);
    for (sx=0.0; sx<(scalar)WIDTH-EPS; sx+=d)
    {
      if (antialias)
      {
	c=zero();
	if (preview)
	{
	  for (px=.75; px<3.0-EPS; px+=1.5)
	    for (py=.75; py<3.0-EPS; py+=1.5)
	      c=vadd(pixel(vassign((scalar)CENTRX-(sx+px),
				   -ZOOM,
				   sy+py-(scalar)CENTRY), shadowsoff),c);
	  c=smul(.25, c);
	}
	else
	{
	  for (px=.125; px<1.0-EPS; px+=.25)
	    for (py=.125; py<1.0-EPS; py+=.25)
	      c=vadd(pixel(vassign((scalar)CENTRX-(sx+px),
				   -ZOOM,
				   sy+py-(scalar)CENTRY), shadowsoff),c);
	  c=smul(.0625, c);
	}
      }
      else
      {
	px=py=0.5;
	c=pixel(vassign((scalar)CENTRX-(sx+px),
			-ZOOM,
			sy+py-(scalar)CENTRY), shadowsoff);
      }
      if (preview)
      {
	long dx, dy, x, y;

	for (dx=-1; dx<=1; dx++)
	  for (dy=-1; dy<=1; dy++)
	  {
	    x=dx+(long)sx;
	    y=dy+(long)sy;
	    if (x>=0 && x<WIDTH && y>=0 && y<HEIGHT)
	      colorplot(c, x, y, map);
	  }
      }
      else
	colorplot(c, (long)sx, (long)sy, map);
    }
  }
}

bool read_scene(char *filename)
{
  FILE *fh;
  int i, n=0;
  
  if (fh=fopen(filename, "r"))
    { 
      n+=fscanf(fh, "globalambient %lf", &globalambient);
      n+=fscanf(fh, "\nlights %d", &xlights);
      if (xlights)
	if (lights=malloc(sizeof(light)*xlights))
	  {
	    for (i=0; i<xlights; i++)
	      {
		n+=fscanf(fh, "\tdirection %lf %lf %lf",
			  &lights[i].dir.x, &lights[i].dir.y, &lights[i].dir.z);
		n+=fscanf(fh, "\tcolor %lf %lf %lf",
			  &lights[i].col.x, &lights[i].col.y, &lights[i].col.z);
	      }
	    printf("\nNumber of lights:  %d\n", xlights);
	  }
	else
	  {
	    fprintf(stderr, "Out of memory!\n");
	    fclose(fh);
	    return FALSE;
	  }
      if (n != 2 + xlights*6)
	{
	  fprintf(stderr, "Error in lightsource data!\n");
	  fclose(fh);
	  return FALSE;
	}
      
      n+=fscanf(fh, "\nspheres %d", &xspheres);
      if (feof(fh))
	n++;
      if (xspheres)
	if ((spheres=malloc(sizeof(sphere)*xspheres)) &&
	    (sphereis=malloc(sizeof(intersection)*xspheres)))
	  {
	    for (i=0; i<xspheres; i++)
	      {
		n+=fscanf(fh, "\tpos %lf %lf %lf",
			  &spheres[i].c.x, &spheres[i].c.y, &spheres[i].c.z);
		n+=fscanf(fh, "\tradius %lf", &spheres[i].r);
		n+=fscanf(fh, "\tcolor0 %lf %lf %lf",
			  &spheres[i].oi.col[0].x, 
			  &spheres[i].oi.col[0].y, 
			  &spheres[i].oi.col[0].z);
		n+=fscanf(fh, "\tambient %lf", &spheres[i].oi.amb);
		n+=fscanf(fh, "\tdiffusion %lf", &spheres[i].oi.diff);
		n+=fscanf(fh, "\tfracspref %lf", &spheres[i].oi.fspre);
		n+=fscanf(fh, "\tsprefexp %lf", &spheres[i].oi.refl);
	      }
	    n--;
	  }
	else
	  {
	    fprintf(stderr, "Out of memory!\n");
	    fclose(fh);
	    return FALSE;
	  }
      if (n != 2 + xlights*6 + xspheres*11)
	{
	  fprintf(stderr, "Error in sphere data!\n");
	  fclose(fh);
	  return FALSE;
	}
      printf("Number of spheres: %d\n", xspheres);

      n+=fscanf(fh, "\nwalls %d", &xwalls);
      if (feof(fh))
	n++;
      if (xwalls)
	if ((walls=malloc(sizeof(wall)*xwalls)) &&
	    (wallis=malloc(sizeof(intersection)*xwalls)))
	  {
	    for (i=0; i<xwalls; i++)
	      {
		n+=fscanf(fh, "\tplane %lf %lf %lf, %lf %lf %lf",
			  &walls[i].pl.n.x, &walls[i].pl.n.y, &walls[i].pl.n.z,
			  &walls[i].pl.p.x, &walls[i].pl.p.y, &walls[i].pl.p.z);
		n+=fscanf(fh, "\tcolor0 %lf %lf %lf",
			  &walls[i].oi.col[0].x, 
			  &walls[i].oi.col[0].y, 
			  &walls[i].oi.col[0].z);
		n+=fscanf(fh, "\tambient %lf", &walls[i].oi.amb);
		n+=fscanf(fh, "\tdiffusion %lf", &walls[i].oi.diff);
		n+=fscanf(fh, "\tfracspref %lf", &walls[i].oi.fspre);
		n+=fscanf(fh, "\tsprefexp %lf", &walls[i].oi.refl);
	      }
	    n--;
	  }
	else
	  {
	    fprintf(stderr, "Out of memory!\n");
	    fclose(fh);
	    return FALSE;
	  }
      if (n != 2 + xlights*6 + xspheres*11 + xwalls*13)
	{
	  fprintf(stderr, "Error in wall data!\n");
	  fclose(fh);
	  return FALSE;
	}
      printf("Number of walls:   %d\n", xwalls);
      
      n+=fscanf(fh, "\ntubes %d", &xtubes);
      if (feof(fh))
	n++;
      if (xtubes)
	if ((tubes=malloc(sizeof(tube)*xtubes)) &&
	    (tubeis=malloc(sizeof(intersection)*xtubes)))
	  {
	    for (i=0; i<xtubes; i++)
	      {
	      n+=fscanf(fh, "\tline %lf %lf %lf, %lf %lf %lf",
			&tubes[i].li.k.x, &tubes[i].li.k.y, &tubes[i].li.k.z,
			&tubes[i].li.p.x, &tubes[i].li.p.y, &tubes[i].li.p.z);
	      n+=fscanf(fh, "\tradius %lf", &tubes[i].r);
	      n+=fscanf(fh, "\tcolor0 %lf %lf %lf",
			&tubes[i].oi.col[0].x, 
			&tubes[i].oi.col[0].y, 
			&tubes[i].oi.col[0].z);
	      n+=fscanf(fh, "\tcolor1 %lf %lf %lf",
			&tubes[i].oi.col[1].x, 
			&tubes[i].oi.col[1].y, 
			&tubes[i].oi.col[1].z);
	      n+=fscanf(fh, "\tcolor2 %lf %lf %lf",
			&tubes[i].oi.col[2].x, 
			&tubes[i].oi.col[2].y, 
			&tubes[i].oi.col[2].z);
	      n+=fscanf(fh, "\tambient %lf",&tubes[i].oi.amb);
	      n+=fscanf(fh, "\tdiffusion %lf",&tubes[i].oi.diff);
	      n+=fscanf(fh, "\tfracspref %lf",&tubes[i].oi.fspre);
	      n+=fscanf(fh, "\tsprefexp %lf",&tubes[i].oi.refl);
	      }
	    n--;
	  }
	else
	  {
	    fprintf(stderr, "Out of memory!\n");
	    fclose(fh);
	    return FALSE;
	  }
      if (n != 2 + xlights*6 + xspheres*11 + xwalls*13 + xtubes*20)
	{
	  fprintf(stderr, "Error in tube data!\n");
	  fclose(fh);
	  return FALSE;
	}
      printf("Number of tubes:   %d\n", xtubes);
    
      n+=fscanf(fh, "\npipes %d", &xpipes);
      if (feof(fh))
	n++;
      if (xpipes)
	if ((pipes=malloc(sizeof(pipe)*xpipes)) &&
	    (pipeis=malloc(sizeof(intersection)*xpipes)))
	  {
	    for (i=0; i<xpipes; i++)
	      {
		n+=fscanf(fh, "\tline %lf %lf %lf, %lf %lf %lf",
			  &pipes[i].li.k.x, &pipes[i].li.k.y, &pipes[i].li.k.z,
			  &pipes[i].li.p.x, &pipes[i].li.p.y, &pipes[i].li.p.z);
		n+=fscanf(fh, "\tradius %lf", &pipes[i].r);
		n+=fscanf(fh, "\tcolor0 %lf %lf %lf", 
			  &pipes[i].oi.col[0].x, 
			  &pipes[i].oi.col[0].y, 
			  &pipes[i].oi.col[0].z);
		n+=fscanf(fh, "\tambient %lf", &pipes[i].oi.amb);
		n+=fscanf(fh, "\tdiffusion %lf", &pipes[i].oi.diff);
		n+=fscanf(fh, "\tfracspref %lf", &pipes[i].oi.fspre);
		n+=fscanf(fh, "\tsprefexp %lf", &pipes[i].oi.refl);
	      }
	    n--;
	  }
	else
	  {
	    fprintf(stderr, "Out of memory!\n");
	    fclose(fh);
	    return FALSE;
	  }
      if (n != 2 + xlights*6 + xspheres*11 + xwalls*13 + xtubes*20 + 
	  xpipes*14)
	{
	  fprintf(stderr, "Error in pipe data!\n");
	  fclose(fh);
	  return FALSE;
	}
      printf("Number of pipes:   %d\n", xpipes);
      
      n+=fscanf(fh, "\ncones %d", &xcones);
      if (feof(fh))
	n++;
      if (xcones)
	if ((cones=malloc(sizeof(cone)*xcones)) &&
	    (coneis=malloc(sizeof(intersection)*xcones)))
	  {
	    for (i=0; i<xcones; i++)
	      {
		n+=fscanf(fh, "\tline %lf %lf %lf, %lf %lf %lf",
			  &cones[i].li.k.x, &cones[i].li.k.y, &cones[i].li.k.z,
			  &cones[i].li.p.x, &cones[i].li.p.y, &cones[i].li.p.z);
		n+=fscanf(fh, "\tradius %lf",&cones[i].r);
		n+=fscanf(fh, "\tcolor0 %lf %lf %lf", 
			  &cones[i].oi.col[0].x, 
			  &cones[i].oi.col[0].y, 
			  &cones[i].oi.col[0].z);
		n+=fscanf(fh, "\tcolor1 %lf %lf %lf", 
			  &cones[i].oi.col[1].x, 
			  &cones[i].oi.col[1].y, 
			  &cones[i].oi.col[1].z);
		n+=fscanf(fh, "\tambient %lf", &cones[i].oi.amb);
		n+=fscanf(fh, "\tdiffusion %lf", &cones[i].oi.diff);
		n+=fscanf(fh, "\tfracspref %lf", &cones[i].oi.fspre);
		n+=fscanf(fh, "\tsprefexp %lf", &cones[i].oi.refl);
	      }
	    n--;
	  }
	else
	  {
	    fprintf(stderr, "Out of memory!\n");
	    fclose(fh);
	    return FALSE;
	  }
      if (n != 2 + xlights*6 + xspheres*11 + xwalls*13 + xtubes*20 + 
	  xpipes*14 + xcones*17)
	{
	  fprintf(stderr, "Error in cone data!\n");
	  fclose(fh);
	  return FALSE;
	}
      printf("Number of cones:   %d\n", xcones);
      
      fclose(fh);
      return TRUE;
    }
  else
    {
      fprintf(stderr, "Could not open file \"%s\"!\n", filename);
      return FALSE;
    }
}

void free_scene(void)
{
	if (lights) free(lights);
	if (spheres) free(spheres);
	if (walls) free(walls);
	if (tubes) free(tubes);
	if (pipes) free(pipes);
	if (cones) free(cones);
	if (sphereis) free(sphereis);
	if (wallis) free(wallis);
	if (tubeis) free(tubeis);
	if (pipeis) free(pipeis);
	if (coneis) free(coneis);
}

/* Tells the user about how to use this program and then exits */

bool parse_commandline(int argc, char *argv[], 
		       char **ecpfile, char **ppmfile,
		       bool *antialias, bool *preview, bool *shadowsoff)
{
  int i;

  *antialias=FALSE;
  *preview=FALSE;
  *shadowsoff=FALSE;

  for (i=1; argc>1 && i<argc; i++)
    if (argv[i][0]=='-')
      switch(argv[i][1])
	{
	case 'a':
	  *antialias=TRUE;
	  break;
	case 'p':
	  *preview=TRUE;
	  break;
	case 's':
	  *shadowsoff=TRUE;
	  break;
	case 'h':
	case '?':
	  printf("Usage: %s [-a | -p | -s | -h] ecpfile ppmfile\n", argv[0]);
	  printf("Flags:\n");
	  printf("-a      switches antialiasing on\n");
	  printf("-p      preview mode; half resolution\n");
	  printf("-s      turns shadows off\n");
	  printf("-h, -?  shows this message\n");
	  return FALSE;
	  break;
	default:
	  fprintf(stderr, "Unknown flag -%c.\n", argv[i][1]);
	  return FALSE;
	  break;
	}
  
  if (argc==1)
    {
      fprintf(stderr, "Missing arguments.\n");
      return FALSE;
    }

  for (i=1; i<argc && argv[i][0]=='-' && 
       (argv[i][1]=='a' || argv[i][1]=='p' || argv[i][1]=='s'); i++);
  if (i >= argc)
    {
      fprintf(stderr, "Missing filenames.\n");
      return FALSE;
    }

  *ecpfile=argv[i++];

  for (; i<argc && argv[i][0]=='-' && 
       (argv[i][1]=='a' || argv[i][1]=='p' || argv[i][1]=='s'); i++);
  if (i >= argc)
    {
      fprintf(stderr, "Missing output file.\n");
      return FALSE;
    }

  *ppmfile=argv[i];
  
  return TRUE;
}

void main(int argc, char *argv[])
{
  char *ecpfile, *ppmfile;
  bool antialias, preview, shadowsoff;
  
  if (parse_commandline(argc, argv, &ecpfile, &ppmfile,
			&antialias, &preview, &shadowsoff))
    {
      /* Allocate memory for the drawing area */
      if (allocpicture(&map, WIDTH, HEIGHT))
	{
	  printf("Opening file \"%s\".\n", ecpfile);
	  if (read_scene(ecpfile))
	    {
	      printf("\n");
	      if (antialias) printf("Using antialiasing.\n");
	      if (preview) printf("Preview mode.\n");
	      if (shadowsoff) printf("Shadows are turned off.\n");

	      /* Clear the drawing area */
	      zeromap(makecol(0,0,0), map);
	      
	      printf("\nComputing:\n");
	      /* Render the whole picture :) */
	      render_the_whole_picture(antialias, preview, shadowsoff);
	      
	      /* Write it to a file */
	      printf("\nSaving \"%s\"...", ppmfile);
	      writeppm(ppmfile, map);
	      printf("done.\n");
	      
	      /* Free memory */
	      free_scene();
	    }
	  /* Free drawing area */
	  freepicture(map, WIDTH);
	}
      else fprintf(stderr, "Out of memory!\n");
    }
}