/*
 * Input description uses these data types:
 *	UINT	An unsigned integer, a la "%u".
 *	FLOAT	A floating-point number, a la "%g".
 *	XYZ	Three whitespace-separated FLOATs: x, y, z.
 *	RGB	Three whitespace-separated FLOATs: r, g, b.
 *	STR	Text string.  This consists of all non-whitespace
 *		  characters or quoted characters up to the next
 *		  unquoted whitespace (", ', and \ quoting can be
 *		  used).
 *	CHAR	Single nonblank character
 *
 * Fixed header:
 *	UINT	Output image X size
 *	UINT	Output image Y size
 *	XYZ	Camera position
 *	XYZ	View direction vector
 *	XYZ	Vertical vector
 *	FLOAT	Number of pixels per 45° of subtended angle
 *	RGB	Colour of the distant background
 *	RGB	Ambient light colour
 *	FLOAT	Brightness multiplier
 *
 * Further sections have a type character, then the rest of the
 *  section.
 *
 * Material definitions:
 *	Type character is 'M'
 *	STR	Material name (cannot be zero-length)
 *	RGB	Self-luminosity colour of the material
 *	RGB	Diffuse reflective colour of the material
 *	FLOAT	Specularity exponent
 *	FLOAT	Specularity multiplier
 *
 * Light source definition:
 *	Type character is 'L'
 *	XYZ	Direction of light source
 *	RGB	Colour of the light
 *
 * Point definition:
 *	Type character is 'P'
 *	UINT	Point number
 *	XYZ	Point location
 *	CHAR	Point type:
 *		's'	Smooth surface internal point, normal follows
 *		'c'	Corner/edge point, use triangle's normal
 *	If type is 's', then
 *	XYZ	Surface normal vector
 *
 * Triangle:
 *	Type character is 'T' or 'D'; 'T' indicates an ordinary surface
 *	  triangle, 'D' a double-sided triangle.  A double-sided
 *	  triangle is equivalent to an ordinary triangle occurring
 *	  twice, once with the vertices in each order.  This is done by
 *	  not doing backface removal on double-sided triangles, and
 *	  negating their normals if necessary so they point towards the
 *	  eye.  Using double-sided triangles with smooth-surface corner
 *	  points may produce surprising results when the triangle is
 *	  seen nearly edge-on.
 *	STR	Material name
 *	UINT	Point number of corner 1
 *	UINT	Point number of corner 2
 *	UINT	Point number of corner 3
 *  Triangles are rendered when they are encountered and do not support
 *  forward references.  The points must be given in clockwise order as
 *  seen from the visible side of the triangle.
 *
 * Output is a PPM file.
 */

#include <math.h>
#include <stdio.h>
#include <ctype.h>
#include <stdarg.h>
#include <stdlib.h>
#include <strings.h>

extern const char *__progname;

typedef struct xyz XYZ;
typedef struct rgb RGB;
typedef struct material MATERIAL;
typedef struct pix PIX;
typedef struct point POINT;
typedef struct light LIGHT;

struct xyz {
  double x;
  double y;
  double z;
  } ;

struct rgb {
  double r;
  double g;
  double b;
  } ;

struct material {
  MATERIAL *link;
  char *name;
  RGB selflum;
  RGB reflect;
  double spec_exp;
  double spec_mul;
  } ;

struct light {
  LIGHT *link;
  XYZ dir;
  RGB colour;
  } ;

struct point {
  XYZ loc;
  XYZ scr;
  unsigned char type;
#define PT_SURFACE 's'
#define PT_CORNER 'c'
  XYZ normal;
  } ;

struct pix {
  unsigned char r;
  unsigned char g;
  unsigned char b;
  double z;
  } ;

static unsigned int oxsize;
static unsigned int oysize;
static XYZ eye;
static XYZ fwd;
static XYZ up;
static XYZ rt;
static double scale;
static RGB bg;
static RGB ambient;
static double bright;
static MATERIAL *materials;
static LIGHT *lights;
static int npoints;
static int apoints;
static POINT **points;
static PIX *screen;
static double xoff;
static double yoff;
static RGB minrgb;
static RGB maxrgb;

static unsigned int getuint(const char *, ...)
	__attribute__((__format__(__printf__,1,2)));
static unsigned int getuint(const char *fmt, ...)
{
 unsigned int v;
 va_list ap;

 if (scanf("%u",&v) != 1)
  { fprintf(stderr,"%s: can't read ",__progname);
    va_start(ap,fmt);
    vfprintf(stderr,fmt,ap);
    va_end(ap);
    fprintf(stderr,"\n");
    exit(1);
  }
 return(v);
}

static XYZ getxyz(const char *, ...)
	__attribute__((__format__(__printf__,1,2)));
static XYZ getxyz(const char *fmt, ...)
{
 XYZ v;
 va_list ap;

 if (scanf("%lg%lg%lg",&v.x,&v.y,&v.z) != 3)
  { fprintf(stderr,"%s: can't read ",__progname);
    va_start(ap,fmt);
    vfprintf(stderr,fmt,ap);
    va_end(ap);
    fprintf(stderr,"\n");
    exit(1);
  }
 return(v);
}

static double getfloat(const char *, ...)
	__attribute__((__format__(__printf__,1,2)));
static double getfloat(const char *fmt, ...)
{
 double v;
 va_list ap;

 if (scanf("%lg",&v) != 1)
  { fprintf(stderr,"%s: can't read ",__progname);
    va_start(ap,fmt);
    vfprintf(stderr,fmt,ap);
    va_end(ap);
    fprintf(stderr,"\n");
    exit(1);
  }
 return(v);
}

static RGB getrgb(const char *, ...)
	__attribute__((__format__(__printf__,1,2)));
static RGB getrgb(const char *fmt, ...)
{
 double r;
 double g;
 double b;
 va_list ap;

 if (scanf("%lg%lg%lg",&r,&g,&b) != 3)
  { fprintf(stderr,"%s: can't read ",__progname);
    va_start(ap,fmt);
    vfprintf(stderr,fmt,ap);
    va_end(ap);
    fprintf(stderr,"\n");
    exit(1);
  }
 if ((r < 0) || (r > 1) || (g < 0) || (g > 1) || (b < 0) || (b > 1))
  { fprintf(stderr,"%s: out-of-range value in ",__progname);
    va_start(ap,fmt);
    vfprintf(stderr,fmt,ap);
    va_end(ap);
    fprintf(stderr,"\n");
    exit(1);
  }
 return((RGB){.r=r,.g=g,.b=b});
}

static char *getstr(const char *, ...)
	__attribute__((__format__(__printf__,1,2)));
static char *getstr(const char *fmt, ...)
{
 char *s;
 int l;
 int a;
 int q;
 int bq;
 int c;
 int any;
 va_list ap;

 static void savechar(int ch)
  { if (l >= a) s = realloc(s,a=l+8);
    s[l++] = ch;
  }

 static void fail(const char *how)
  { fprintf(stderr,"%s: %s ",__progname,how);
    vfprintf(stderr,fmt,ap);
    fprintf(stderr,"\n");
    exit(1);
  }

 va_start(ap,fmt);
 s = 0;
 l = 0;
 a = 0;
 q = 0;
 bq = 0;
 any = 0;
 while (1)
  { c = getchar();
    switch (c)
     { case '\0':
	  fail("NUL in");
	  break;
       case EOF:
	  fail("unexpected EOF in");
	  break;
     }
    if (!q && !bq && isspace(c))
     { if (any) break;
       continue;
     }
    any = 1;
    if (bq)
     { savechar(c);
       bq = 0;
     }
    else if (c == '\\')
     { bq = 1;
     }
    else if (c == q)
     { q = 0;
     }
    else if ((c == '"') || (c == '\''))
     { q = c;
     }
    else
     { savechar(c);
     }
  }
 savechar(0);
 va_end(ap);
 return(s);
}

static MATERIAL *material_find(const char *n)
{
 MATERIAL *m;
 MATERIAL **mp;

 if (! materials) return(0);
 if (! strcmp(materials->name,n)) return(materials);
 mp = &materials->link;
 while ((m = *mp))
  { if (!strcmp(m->name,n))
     { *mp = m->link;
       m->link = materials;
       materials = m;
       return(m);
     }
    mp = &m->link;
  }
 return(0);
}

static double dot(XYZ a, XYZ b)
{
 return((a.x*b.x)+(a.y*b.y)+(a.z*b.z));
}

static XYZ smul3(double f, XYZ v)
{
 return((XYZ){.x=v.x*f,.y=v.y*f,.z=v.z*f});
}

static double norm(XYZ v)
{
 return(sqrt((v.x*v.x)+(v.y*v.y)+(v.z*v.z)));
}

static XYZ unit(XYZ v)
{
 return(smul3(1/norm(v),v));
}

static XYZ cross(XYZ a, XYZ b)
{
 return((XYZ){ .x = (a.y*b.z) - (a.z*b.y),
	       .y = (a.z*b.x) - (a.x*b.z),
	       .z = (a.x*b.y) - (a.y*b.x) });
}

static XYZ sub3(XYZ a, XYZ b)
{
 return((XYZ){.x=a.x-b.x,.y=a.y-b.y,.z=a.z-b.z});
}

static XYZ interp3(XYZ a, double alpha, XYZ b)
{
 return((XYZ){ .x = (alpha*b.x) + ((1-alpha)*a.x),
	       .y = (alpha*b.y) + ((1-alpha)*a.y),
	       .z = (alpha*b.z) + ((1-alpha)*a.z) });
}

static POINT interpp(POINT a, double alpha, POINT b)
{
 return((POINT){ .loc = interp3(a.loc,alpha,b.loc),
		 .scr = interp3(a.scr,alpha,b.scr),
		 .type = PT_SURFACE,
		 .normal = unit(interp3(a.normal,alpha,b.normal)) });
}

static double interp1(double a, double alpha, double b)
{
 return((alpha*b)+((1-alpha)*a));
}

static XYZ screen_project(XYZ p)
{
 double d;
 double u;
 double r;

 p = sub3(p,eye);
 d = dot(p,fwd);
 u = dot(p,up);
 r = dot(p,rt);
 return((XYZ){.x=r*scale/d,.y=u*scale/d,.z=d});
}

static RGB compute_colour(MATERIAL *m, int ds, XYZ loc, XYZ normal)
{
 LIGHT *l;
 RGB v;
 double rp;
 double dc;
 XYZ e;

 static void addv(double r, double g, double b)
  { v.r += r;
    v.g += g;
    v.b += b;
  }

 v = m->selflum;
 addv(ambient.r*m->reflect.r,ambient.g*m->reflect.g,ambient.b*m->reflect.b);
 e = unit(sub3(eye,loc));
 if (ds && (dot(normal,e) < 0)) normal = smul3(-1,normal);
 for (l=lights;l;l=l->link)
  { double t1;
    XYZ t2;
    dc = dot(normal,l->dir);
    if (dc < 0) continue;
    t2 = smul3(2*dc,normal);
    t2 = sub3(t2,l->dir);
    t1 = dot(t2,e);
    rp = (t1 < 0) ? 0 : (m->spec_mul > 0) ? m->spec_mul*pow(t1,m->spec_exp) : 0;
    addv( (rp+(m->reflect.r*dc)) * l->colour.r,
	  (rp+(m->reflect.g*dc)) * l->colour.g,
	  (rp+(m->reflect.b*dc)) * l->colour.b );
  }
 return((RGB){ .r = v.r*bright, .g = v.g*bright, .b = v.b*bright });
}

static PIX makepix(RGB c, double z)
{
 PIX p;
 int i;

 if (c.r < minrgb.r) minrgb.r = c.r;
 if (c.g < minrgb.g) minrgb.g = c.g;
 if (c.b < minrgb.b) minrgb.b = c.b;
 if (c.r > maxrgb.r) maxrgb.r = c.r;
 if (c.g > maxrgb.g) maxrgb.g = c.g;
 if (c.b > maxrgb.b) maxrgb.b = c.b;
 i = c.r * 256;
 p.r = (i<0) ? 0 : (i>255) ? 255 : i;
 i = c.g * 256;
 p.g = (i<0) ? 0 : (i>255) ? 255 : i;
 i = c.b * 256;
 p.b = (i<0) ? 0 : (i>255) ? 255 : i;
 p.z = z;
 return(p);
}

static void render_rows(MATERIAL *m, int ds, POINT t1, POINT t2, POINT b1, POINT b2)
{
 int y1;
 int y2;
 int y;
 POINT e1;
 POINT e2;
 int x1;
 int x2;
 int x;
 double xi;
 double z;
 PIX *pp;
 RGB c;

 y1 = ceil(t1.scr.y+yoff);
 y2 = floor(b1.scr.y+yoff);
 if (y1 < 0) y1 = 0;
 if (y2 >= (int)oysize) y2 = oysize - 1;
 if (y1 > y2) return;
 for (y=y1;y<=y2;y++)
  { e1 = interpp(t1,(y-yoff-t1.scr.y)/(b1.scr.y-t1.scr.y),b1);
    e2 = interpp(t2,(y-yoff-t1.scr.y)/(b1.scr.y-t1.scr.y),b2);
    if (e1.scr.x == e2.scr.x) continue;
    x1 = ceil(e1.scr.x+xoff);
    x2 = floor(e2.scr.x+xoff);
    if (x1 < 0) x1 = 0;
    if (x2 >= (int)oxsize) x2 = oxsize - 1;
    if (x1 > x2) continue;
    pp = &screen[((oysize-1-y)*oxsize)+x1];
    for (x=x1;x<=x2;x++,pp++)
     { xi = (x - xoff - e1.scr.x) / (e2.scr.x - e1.scr.x);
       z = interp1(e1.scr.z,xi,e2.scr.z);
       if ((z < 0) || (z > pp->z)) continue;
       c = compute_colour(m,ds,interp3(e1.loc,xi,e2.loc),interp3(e1.normal,xi,e2.normal));
       *pp = makepix(c,z);
     }
  }
}

static void render_t_up(MATERIAL *m, int ds, POINT top, POINT b1, POINT b2)
{
 if (b2.scr.x < b1.scr.x)
  { render_rows(m,ds,top,top,b2,b1);
  }
 else
  { render_rows(m,ds,top,top,b1,b2);
  }
}

static void render_t_dn(MATERIAL *m, int ds, POINT t1, POINT t2, POINT bot)
{
 if (t2.scr.x < t1.scr.x)
  { render_rows(m,ds,t2,t1,bot,bot);
  }
 else
  { render_rows(m,ds,t1,t2,bot,bot);
  }
}

static void render_triangle(MATERIAL *m, POINT *pa, POINT *pb, POINT *pc, int doublesided)
{
 POINT *t;
 XYZ n;
 int have_n;

 static XYZ normal(void)
  { if (! have_n)
     { n = unit(cross(sub3(pa->loc,pb->loc),sub3(pc->loc,pb->loc)));
       have_n = 1;
     }
    return(n);
  }

printf("# render_triangle: (%g,%g,%g) - (%g,%g,%g) - (%g,%g,%g)\n",
pa->loc.x,pa->loc.y,pa->loc.z,
pb->loc.x,pb->loc.y,pb->loc.z,
pc->loc.x,pc->loc.y,pc->loc.z);
 if (!doublesided && (dot(fwd,cross(sub3(pc->loc,pb->loc),sub3(pa->loc,pb->loc))) < 0)) return;
 have_n = 0;
 if (pa->type == PT_CORNER) pa->normal = normal();
 if (pb->type == PT_CORNER) pb->normal = normal();
 if (pc->type == PT_CORNER) pc->normal = normal();
 if ((pb->scr.y <= pa->scr.y) && (pb->scr.y <= pc->scr.y))
  { t = pa; pa = pb; pb = t; }
 else if ((pc->scr.y <= pa->scr.y) && (pc->scr.y <= pb->scr.y))
  { t = pa; pa = pc; pc = t; }
 if ((pb->scr.y >= pa->scr.y) && (pb->scr.y >= pc->scr.y))
  { t = pc; pc = pb; pb = t; }
 else if ((pa->scr.y >= pc->scr.y) && (pa->scr.y >= pb->scr.y))
  { t = pc; pc = pa; pa = t; }
 if (pa->scr.y < pb->scr.y)
  { render_t_up(m,doublesided,
		*pa,
		interpp( *pa,
			 (pb->scr.y-pa->scr.y)/(pc->scr.y-pa->scr.y),
			 *pc ),
		*pb);
  }
 if (pc->scr.y > pb->scr.y)
  { render_t_dn(m,doublesided,
		interpp( *pa,
			 (pb->scr.y-pa->scr.y)/(pc->scr.y-pa->scr.y),
			 *pc ),
		*pb,
		*pc);
  }
}

static void init(void)
{
 materials = 0;
 points = 0;
 apoints = 0;
 npoints = 0;
}

static void read_header(void)
{
 oxsize = getuint("output X size");
 oysize = getuint("output Y size");
 eye = getxyz("camera position");
 fwd = getxyz("view direction");
 up = getxyz("vertical vector");
 scale = getfloat("pixel size");
 bg = getrgb("background colour");
 ambient = getrgb("ambient light");
 bright = getfloat("brightness compensation");
}

static void init_render(void)
{
 int n;
 int i;
 int j;

 printf("P3\n%u %u\n255\n",oxsize,oysize);
 printf("# eye (%g,%g,%g)\n",eye.x,eye.y,eye.z);
 printf("# fwd (%g,%g,%g)\n",fwd.x,fwd.y,fwd.z);
 printf("# up (%g,%g,%g)\n",up.x,up.y,up.z);
 printf("# scale %g\n",scale);
 n = oxsize * oysize;
 screen = malloc(n*sizeof(PIX));
 screen[0] = makepix(bg,1e30);
 j = 1;
 for (i=2;i<=n;j=i,i<<=1) bcopy(screen,screen+j,j*sizeof(PIX));
 if (j < n) bcopy(screen,screen+j,(n-j)*sizeof(PIX));
 fwd = unit(fwd);
 up = unit(sub3(up,smul3(dot(up,fwd),fwd)));
 rt = cross(fwd,up);
 printf("# ----\n");
 printf("# fwd (%g,%g,%g)\n",fwd.x,fwd.y,fwd.z);
 printf("# up (%g,%g,%g)\n",up.x,up.y,up.z);
 printf("# rt (%g,%g,%g)\n",rt.x,rt.y,rt.z);
 xoff = oxsize / 2.0;
 yoff = oysize / 2.0;
 minrgb = (RGB) { .r = 1, .g = 1, .b = 1 };
 maxrgb = (RGB) { .r = 0, .g = 0, .b = 0 };
}

static void read_comment(void)
{
 int c;

 do
  { c = getchar();
  } while ((c != '\n') && (c != EOF));
}

static void read_material(void)
{
 char *s;
 RGB lum;
 RGB ref;
 double s_e;
 double s_m;
 MATERIAL *m;

 s = getstr("material name");
 lum = getrgb("material `%s' self-luminant colour",s);
 ref = getrgb("material `%s' reflective colour",s);
 s_e = getfloat("material `%s' specularity exponent",s);
 s_m = getfloat("material `%s' specularity multiplier",s);
 m = malloc(sizeof(MATERIAL));
 m->name = s;
 m->selflum = lum;
 m->reflect = ref;
 m->spec_exp = s_e;
 m->spec_mul = s_m;
 m->link = materials;
 materials = m;
}

static void read_light(void)
{
 XYZ dir;
 RGB col;
 LIGHT *l;

 dir = getxyz("light direction");
 col = getrgb("light colour");
 l = malloc(sizeof(LIGHT));
 l->dir = unit(dir);
 l->colour = col;
 l->link = lights;
 lights = l;
}

static void read_point(void)
{
 int n;
 XYZ at;
 POINT *p;
 char type;
 XYZ normal;

 n = getuint("point number");
 if (n < 0)
  { fprintf(stderr,"%s: negative point number %d\n",__progname,n);
    exit(1);
  }
 at = getxyz("point location");
 if (scanf(" %c",&type) != 1)
  { fprintf(stderr,"%s: can't read point type\n",__progname);
    exit(1);
  }
 switch (type)
  { case 's':
       type = PT_SURFACE;
       normal = unit(getxyz("point surface normal"));
       break;
    case 'c':
       type = PT_CORNER;
       break;
    default:
       fprintf(stderr,"%s: bad point type `%c'\n",__progname,type);
       exit(1);
       break;
  }
 p = malloc(sizeof(POINT));
 p->loc = at;
 p->scr = screen_project(at);
 p->type = type;
 switch (type)
  { case PT_SURFACE:
       p->normal = normal;
       break;
    case PT_CORNER:
       break;
  }
 if (n >= apoints) points = realloc(points,(apoints=n+16)*sizeof(POINT *));
 for (;npoints<=n;npoints++) points[npoints] = 0;
 free(points[n]);
 points[n] = p;
}

static void read_triangle(int doublesided)
{
 int i;
 char *mn;
 MATERIAL *m;
 unsigned int v[3];

 mn = getstr("triangle material name");
 m = material_find(mn);
 if (m == 0)
  { fprintf(stderr,"%s: undefined material %s\n",__progname,mn);
    exit(1);
  }
 free(mn);
 for (i=0;i<3;i++)
  { v[i] = getuint("triangle vertex %d",i+1);
    if ((v[i] >= npoints) || !points[v[i]])
     { fprintf(stderr,"%s: undefined triangle point %d\n",__progname,v[i]);
       exit(1);
     }
  }
 printf("# triangle %s %d %d %d\n",m->name,v[0],v[1],v[2]);
 render_triangle(m,points[v[0]],points[v[1]],points[v[2]],doublesided);
}

static void read_body(void)
{
 int type;

 while (1)
  { scanf(" ");
    type = getchar();
    if (type == EOF) break;
    switch (type)
     { case '#': read_comment();   break;
       case 'M': read_material();  break;
       case 'L': read_light();     break;
       case 'P': read_point();     break;
       case 'T': read_triangle(0); break;
       case 'D': read_triangle(1); break;
       default:
	  fprintf(stderr,"%s: bad body type character `%c'\n",__progname,type);
	  exit(1);
     }
  }
}

static void dump_output(void)
{
 int i;
 PIX *p;

 p = screen;
 for (i=oxsize*oysize;i>0;i--,p++) printf("%d %d %d\n",p->r,p->g,p->b);
 if ( (minrgb.r < 0) || (minrgb.g < 0) || (minrgb.b < 0) ||
      (maxrgb.r > 1) || (maxrgb.g > 1) || (maxrgb.b > 1) )
  { fprintf(stderr,"RGB ranges [%g..%g], [%g..%g], [%g..%g]\n",minrgb.r,maxrgb.r,minrgb.g,maxrgb.g,minrgb.b,maxrgb.b);
  }
}

int main(void);
int main(void)
{
 init();
 read_header();
 init_render();
 read_body();
 dump_output();
 return(0);
}