// rtray.cpp : Defines the entry point for the application.
//



#define RT_DEBUG_MODE

#include "stdafx.h"

#include <windows.h>
#include <stdlib.h>
#include <math.h>
#include <stdio.h>
#include <mmsystem.h>

#include "text.h"

#include "rt.h"
#include "rt_scanline.h"


/*


  - Optimisations
    . Sphere screen space bounding boxes
      . Calculate based on cube
      . Make list of possible spheres at the start of each scanline
    . Work out differential by 4 pixels for a sphere

  # 3D point to screen point function


2004-06-20  First stab at bounding box function for sphere

*/



//---- rt_view_create

RT_VIEW* rt_view_create( int w, int h )
{
  RT_VIEW* view = (RT_VIEW *)malloc( sizeof( RT_VIEW ) );

  //- Initialise view struct

  memset( (void *)view, 0, sizeof( RT_VIEW ) );

  view->w = w;
  view->h = h;
  view->size = w*h;

  view->buf = (int *)malloc( sizeof( int ) * view->size );
  view->obj = (unsigned int *)malloc( sizeof( unsigned int ) * view->size );
  view->sob = (unsigned int *)malloc( sizeof( unsigned int ) * view->size );

  //- Useful precalculated numbers

  view->bx1 = -1.0f*w/h;
  view->bx2 =  1.0f*w/h;
  view->bdx = 2.0f/h;
  //view->dx2 = dx * 2.0f;

  view->by1 =  1.0f;
  view->by2 = -1.0f;
  view->bdy = -2.0f/h;
  //view->dy2 = dy * 2.0f;

  //- Precalculate various values

  view->w2 = w*2;
  view->w3 = w*3;
  view->w4 = w*4;
  view->wm4 = w-4;

  //- Done

  return view;
}



//---- rt_view_destroy

void rt_view_destroy( RT_VIEW *view )
{
  //- Delete screen buffer then view structure

  free( (void *)view->obj );
  free( (void *)view->buf );
  free( (void *)view );
}



//---- SCENE FUNCTIONS ---------------------------------------------------------------------------



//---- rt_scene_create

RT_SCENE* rt_scene_create( )
{
  RT_SCENE* scene = (RT_SCENE *)malloc( sizeof( RT_SCENE ) );

  //- Initialise scene struct

  memset( (void *)scene, 0, sizeof( RT_SCENE ) );

  scene->numSpheres = 1;  //  Sphere 0 is a black background sphere

  //- Done

  return scene;
}



//---- rt_scene_destroy

void rt_scene_destroy( RT_SCENE *scene )
{
  free( (void *)scene );
}



//---- rt_scene_sphere

int rt_scene_sphere( RT_SCENE *scene, float x, float y, float z, float r, int col )
{
  int num = scene->numSpheres;

  //- ERROR: Cannot add sphere to scene - reached maximum (RT_SCENE_MAXSPHERES)

  if( scene->numSpheres >= RT_SCENE_MAXSPHERES ) return false;

  //- Initialise sphere struct

  memset( (void *)(scene->sphere + num), 0, sizeof( RT_SPHERE ) );

  scene->sphere[ num ].ax = x;
  scene->sphere[ num ].ay = y;
  scene->sphere[ num ].az = z;
  scene->sphere[ num ].col = col;
  scene->sphere[ num ].rad = r;

  scene->sphere[ num ].r2 = r*r;

  //- Increase number of spheres in scene

  scene->numSpheres ++;

  //- DONE: Returns an ID for the sphere

  return num;
}



//---- rt_scene_light

int rt_scene_light( RT_SCENE *scene, float x, float y, float z, int col )
{
  //- Initialise light struct

  memset( (void *)(&scene->light), 0, sizeof( RT_LIGHT ) );

  scene->light.ax = scene->light.x = x;
  scene->light.ay = scene->light.y = y;
  scene->light.az = scene->light.z = z;
  scene->light.col = col;

  return 0;
}



//---- rt_scene_setRotation

void rt_scene_setRotation( RT_SCENE *scene, float rx, float ry )
{
  int i;

  //- Per-rotation constants

  float xs = (float)sin(3.14159*rx);
  float xc = (float)cos(3.14159*rx);
  float ys = (float)sin(3.14159*ry);
  float yc = (float)cos(3.14159*ry);

  scene->rx = rx;
  scene->ry = ry;
  scene->xs = xs;
  scene->xc = xc;
  scene->ys = ys;
  scene->yc = yc;

  //- Transform spheres

  RT_SPHERE *sphere = scene->sphere;

  float tx, ty, tz;

  for(i=0; i<scene->numSpheres; i++)
  {
    tx = sphere[i].ax * yc + sphere[i].az * ys;
    tz = sphere[i].az * yc - sphere[i].ax * ys;
    ty = sphere[i].ay;

    sphere[i].x =      tx;
    sphere[i].y =      ty * xc + tz * xs;
    sphere[i].z = 10 + tz * xc - ty * xs;
  }

  //- Work out constants for first ray (assumed to be from origin)

  for(i=0; i<scene->numSpheres; i++)
  {
    sphere[i].c = sphere[i].x*sphere[i].x + sphere[i].y*sphere[i].y + sphere[i].z*sphere[i].z - sphere[i].r2;
  }

}



//---- RAY PRECALC FUNCTIONS -------------------------------------------------//

//  wha?

#define RY_RAY_FIRST_PRECALC( ray )     \
  ray.dx2 = ray.dx*ray.dx;              \
  ray.dy2 = ray.dy*ray.dy;              \
  ray.dz2 = 1.0f;                       \
	ray.len2 = ray.dx2 + ray.dy2 + 1.0f;  \
	ray.leno22 = 0.5f / ray.len2;         \
	ray.len42 = 4.0f * ray.len2;

#define RY_RAY_PRECALC( ray )              \
  ray.dx2 = ray.dx*ray.dx;                 \
  ray.dy2 = ray.dy*ray.dy;                 \
  ray.dz2 = ray.dy*ray.dy;                 \
	ray.len2 = ray.dx2 + ray.dy2 + ray.dz2;  \
	ray.leno22 = 0.5f / ray.len2;            \
	ray.len42 = 4.0f * ray.len2;



inline void rt_ray_first_precalc(RT_RAY *ray)
{
  ray->dx2 = ray->dx*ray->dx;
  ray->dy2 = ray->dy*ray->dy;
  ray->dz2 = 1.0f;
	ray->len2 = ray->dx2 + ray->dy2 + 1.0f;
	ray->leno22 = 0.5f / ray->len2;
	ray->len42 = 4.0f * ray->len2;
}


inline void rt_ray_precalc(RT_RAY *ray)
{
  ray->dx2 = ray->dx*ray->dx;
  ray->dy2 = ray->dy*ray->dy;
  ray->dz2 = ray->dz*ray->dz;
	ray->len2 = ray->dx2 + ray->dy2 + ray->dz2;
	ray->leno22 = 0.5f / ray->len2;
	ray->len42 = 4.0f * ray->len2;
}



//---- SHADOW ROUTINES -------------------------------------------------------//



//---- rt_intersect_shadowRay

//  Returns the index of a sphere intersected by the shadow ray

//  Takes the ray's intersection point as the start of the ray,
//  and the scene's light as the end of the ray.

int rt_intersect_shadowRay( RT_SCENE *scene, RT_RAY *ray )
{
  RT_RAY shadow_ray;

  //  Start point intersection point of ray - ray->ix etc

  shadow_ray.ox = ray->ix;
  shadow_ray.oy = ray->iy;
  shadow_ray.oz = ray->iz;

  //  End point scene->light

  shadow_ray.dx = scene->light.x - ray->ix;
  shadow_ray.dy = scene->light.y - ray->iy;
  shadow_ray.dz = scene->light.z - ray->iz;

  //  Check object normal

  float d = ray->inx * shadow_ray.dx + ray->iny * shadow_ray.dy + ray->inz * shadow_ray.dz;
  if(d < 0) return ray->isphere;

  //  Perform intersection

//  int k = rt_intersect( scene, shadow_ray );

  return 0;
}



//---- GENERAL (NON-FIRST HIT) INTERSECTION ROUTINES ------------------------------------//



//---- rt_intersect_sphere

float rt_intersect_sphere( RT_RAY *ray, RT_SPHERE *sphere )
{
  float x = sphere->x;
  float y = sphere->y;
  float z = sphere->z;
  float r = sphere->rad;

//  float a = ray->dx*ray->dx + ray->dy*ray->dy + ray->dz*ray->dz;
  float b = 2*(ray->dx*(ray->ox-x) + ray->dy*(ray->oy-y) + ray->dz*(ray->oz-z));
  float c = (ray->ox-x)*(ray->ox-x) + (ray->oy-y)*(ray->oy-y) + (ray->oz-z)*(ray->oz-z) - r*r;
//  float d = b*b - 4*a*c;

//  float b = -2.0f*(ray->dx*(ray->ox-sphere->x) + ray->dy*(ray->oy-sphere->y) + ray->dz*(ray->oz-sphere->z));

  float d = b*b - ray->len42*c;

  if( d > 0)
  {
    d = sqrtf(d);

    return (-b - d) * ray->leno22;  //  distance along ray of nearest point on sphere (entry point on ray)
//    return (-b + d) * ray->leno22;  //  distance along ray of FARTHEST point on sphere (exit point of ray)
  }

  return -1;
}



//---- rt_intersect

//  Intersects the ray with the scene and returns the index of the sphere it hits

int rt_intersect( RT_SCENE *scene, RT_RAY *ray )
{
  return 0;
}



//---- FIRST-HIT INTERSECTION ROUTINES -------------------------------------------------//



//---- rt_intersect_first_sphere

//  Given a sphere and a ray returns the distance along the ray that the sphere intersects
//  If the ray does not intersect the sphere, returns < 0

inline float rt_intersect_first_sphere( RT_RAY *ray, RT_SPHERE *sphere )
{
  float b = -2.0f*(ray->dx*sphere->x + ray->dy*sphere->y + sphere->z);

  float d = b*b - ray->len42*sphere->c;

  if( d > 0 )
  {
    d = (float)sqrtf(d);

    return (-b - d) * ray->leno22;
  }

  return -1.0f;
}



//---- rt_intersect_first

int rt_intersect_first( RT_SCENE *scene, RT_RAY *ray )
{
  //- Find which sphere this ray intersects

  int i, foundi = 0;
  float dd, d = RT_FAR_PLANE;
/*
  // old, simple way of doing it -- check every sphere!

  for(i=1; i<scene->numSpheres; i++)
  {
    dd = rt_intersect_first_sphere( ray, scene->sphere+i );
    if(dd > RT_NEAR_PLANE && dd < d) { d = dd; foundi = i; }
  }
*/

  //  New shiny way of doing it - just check spheres whose bounding box we are within
  //  They are stored in a linked list in scene->sphere

  float current_z = 99999;

  for( i = scene->sphere[0].in_next; i; i = scene->sphere[i].in_next )
  {
    //  v0013: check near bounding plane of sphere before slower test??

    //  If the sphere's nearest Z is further away than the current Z it's not worth doing slower test

    //  0.77 fps -> 0.85 fps after adding this line :)  - 10000 spheres at 320x240

    if( scene->sphere[i].z - scene->sphere[i].rad > current_z ) continue;

    //  Get intersection point of sphere

    dd = rt_intersect_first_sphere( ray, scene->sphere+i );

    //  If it's nearer than our best so far then this is our new nearest hit

    if( dd > RT_NEAR_PLANE && dd < d ) { d = dd; foundi = i; current_z = ray->dz * d; }
  }

  //  For intersected sphere, store intersection info in ray:
  //    Sphere index, intersection co-ords, surface normal @ intersection

  float a = 1.0f/scene->sphere[foundi].rad;

  ray->ix = /*ray->ox +*/ ray->dx * d;
  ray->inx = a*(ray->ix - scene->sphere[foundi].x );

  ray->iy = /*ray->oy +*/ ray->dy * d;
  ray->iny = a*(ray->iy - scene->sphere[foundi].y );

  ray->iz = /*ray->oz +*/ ray->dz * d;
  ray->inz = a*(ray->iz - scene->sphere[foundi].z );

  ray->id = d;
  ray->isphere = foundi;

  return foundi;
}



//---- SCREEN SPACE & 3D FUNCTIONS -------------------------------------------//



//---- rt_rotatePoint

//  Transforms the point by the scene's object roataion and z offset

void rt_rotatePoint( RT_SCENE *scene, float x, float y, float z, float *ax, float *ay, float *az )
{
  float tx, ty, tz;

  tx = x * scene->yc + z * scene->ys;
  tz = z * scene->yc - x * scene->ys;
  ty = y;

  *ax =      tx;
  *ay =      ty * scene->xc + tz * scene->xs;
  *az = 10 + tz * scene->xc - ty * scene->xs;
}



//---- rt_screenPoint

//  Calculates the screen co-ords of the given 3d point

//  Screen point returned is -1..1

void rt_viewPoint( RT_SCENE *scene, RT_VIEW *view, float x, float y, float z, float *vx, float *vy )
{
  //- Transform (x,y,z) by the scene's transformation

  rt_rotatePoint( scene, x, y, z, &x, &y, &z );

  //- Screen point ( ranges - X: -w/h to w/h, Y: -1 to 1 ) 

  float oz = 1.0f / z;

  x *= oz;
  y *= -oz;

  //  Calculate view point by scaling by view buffer size etc

  *vx = 0.5f * view->w + 0.5f * x * view->h;
  *vy = 0.5f * (y+1) * view->h;
}



//---- BOUNDING BOX CODE -----------------------------------------------------//

/*
  Given a sphere's position and radius, we work out the bounding rect onscreen
  Pretty clever actually, took 4 pages of scribbles and help from Derive.

  Assumes view point is at (0,0,0)
  Assumes screen plane is (u,v,1) - u and v are our screen_x and screen_y positions
  
  Returns false if sphere offscreen

  float *bounds is a pointer to an array with 4 elements ( x1, y1, x2, y2 )
*/

bool rt_sphere_boundingRect( float x, float y, float z, float r, float *bounds )
{
  //  Main constants, based on sphere info

  float const_d2 = x*x + y*y + z*z;
  float const_d  = sqrtf( const_d2 );  //  Distance from view point to

  float const_l  = const_d2 - r*r;

  if( const_l <= 0 ) return false;  //  Sphere offscreen

  float const_c  = sqrtf( const_l ) / const_d;
  float const_c2 = const_c * const_c;

  float k = const_c2 * const_d2;

  //  The elipse equation is au^2 + bu + cuv + dv + ev^2 = 0

  float a = x*x - k;
  float b = 2*x*z;
  float c = 2*x*y;
  float d = 2*y*z;
  float e = y*y - k;
  float f = z*z - k;

  //  The elipse equation can be expressed as a quadratic of u or v:

  //    [a] u^2 + [b+cv] u + [ev^2+dv+f] = 0
  //    [e] v^2 + [d+cu] v + [au^2+bu+f] = 0

  //  Now, from school days, a quadratic like this:

  //      A x^2 +      B x +           C = 0

  //  Has one solution when B^2 - 4AC = 0

  float alpha, beta, gamma, g;

  //  A magic constant which appears in derive's solutions for both u and v

  float brian = ( a*( d*d - 4*e*f ) + b*(b*e - c*d) + c*c*f);

  //  So we can find the v vaules at which the two u values converge (screen min and max v )
  //  by throwing Solve(   (b + c·v)^2 - 4·a·(d·v + e·v^2 + f) = 0  ,   v ) into Derive

  //    v = ( alpha +/- 2 sqrt( gamma ) ) / beta

  alpha = 2*a*d - b*c;
  beta  = c*c - 4*a*e;
  gamma = a * brian;

  if( gamma <= 0 ) return false;  //  No roots - means sphere is offscreen

  g = 2 * sqrtf( gamma );

  float v1 = ( alpha - g ) / beta;
  float v2 = ( alpha + g ) / beta;

  //  Now find the u vaules at which the two v values converge (screen min and max u )
  //  by throwing Solve(   (d + c·u)^2 - 4·e·(b·u + a·u^2 + f) = 0  ,   u ) into Derive

  //    u = ( alpha +/- 2 sqrt( gamma ) ) / beta

  alpha = 2*e*b - d*c;
//  beta  = c*c - 4*a*e;  //  Unchanged from our u calculation
  gamma = e * brian;

  if( gamma <= 0 ) return false;  //  No roots - means sphere is offscreen

  g = 2 * sqrtf( gamma );

  float u1 = ( alpha - g ) / beta;
  float u2 = ( alpha + g ) / beta;

  //  Done!!  Now store in user-supplied array

  bounds[0] = u1<u2 ? u1 : u2;
  bounds[1] = v1<v2 ? v2 : v1;
  bounds[2] = u1<u2 ? u2 : u1;
  bounds[3] = v1<v2 ? v1 : v2;

  return true;
}



//---- BASIC PIXEL ROUTINES --------------------------------------------------//



#define RT_COL_MIDDLE( c1, c2 )    ( ((c1>>1)&0x7F7F7F7F) + ((c2>>1)&7F7F7F7F) )




int rt_col_blend( int c1, int c2, float a )
{
  int r1 = (c1>>16)&255;
  int g1 = (c1>> 8)&255;
  int b1 = (c1    )&255;

  int r2 = (c2>>16)&255;
  int g2 = (c2>> 8)&255;
  int b2 = (c2    )&255;

  r1 += (r2-r1)*a;
  g1 += (g2-g1)*a;
  b1 += (b2-b1)*a;

  if( r1 < 0 ) r1 = 0;
  if( g1 < 0 ) g1 = 0;
  if( b1 < 0 ) b1 = 0;
  if( r1 > 255 ) r1 = 255;
  if( g1 > 255 ) g1 = 255;
  if( b1 > 255 ) b1 = 255;

  return (r1<<16) + (g1<<8) + b1;
}



int rt_col_add( int c1, int c2, float a )
{
  int r1 = (c1>>16)&255;
  int g1 = (c1>> 8)&255;
  int b1 = (c1    )&255;

  int r2 = (c2>>16)&255;
  int g2 = (c2>> 8)&255;
  int b2 = (c2    )&255;

  r1 += r2*a;
  g1 += g2*a;
  b1 += b2*a;

  if( r1 < 0 ) r1 = 0;
  if( g1 < 0 ) g1 = 0;
  if( b1 < 0 ) b1 = 0;
  if( r1 > 255 ) r1 = 255;
  if( g1 > 255 ) g1 = 255;
  if( b1 > 255 ) b1 = 255;

  return (r1<<16) + (g1<<8) + b1;
}





//---- RENDERING ROUTINES ----------------------------------------------------//



//---- rt_scene_pixel_sphere

//  Returns the index of the sphere at this point
//  Inefficient for repeated calling

int rt_scene_pixel_sphere( RT_VIEW *view, RT_SCENE *scene, int x, int y )
{
  //- Set up ray

  RT_RAY ray;

  ray.ox = 0;
  ray.oy = 0;
  ray.oz = 0;

  ray.dx = view->bx1 + x*view->bdx;
  ray.dy = view->by1 + y*view->bdy;
  ray.dz = 1;

  //- Precalculate some ray values

  rt_ray_first_precalc( &ray );

  //- Find out what this ray intersects

  return rt_intersect_first( scene, &ray );
}



//---- rt_scene_pixel

//  Returns the colour of the pixel

int rt_scene_pixel( RT_VIEW *view, RT_SCENE *scene, int x, int y )
{
  //- Set up ray

  RT_RAY ray;

  ray.ox = 0;
  ray.oy = 0;
  ray.oz = 0;

  ray.dx = view->bx1 + x*view->bdx;
  ray.dy = view->by1 + y*view->bdy;
  ray.dz = 1;

  //- Precalculate some ray values

  rt_ray_first_precalc( &ray );

  //- Find out what this ray intersects

  int obj = rt_intersect_first( scene, &ray );
  int col = scene->sphere[ obj ].col;
  int sob = rt_intersect_shadowRay( scene, &ray );

  //  v0012: work out lighting

  if( !sob )
  {
    //  unit vector towards light

    double lx = scene->light.x - ray.ix;
    double ly = scene->light.y - ray.iy;
    double lz = scene->light.z - ray.iz;

    double d = lx*lx + ly*ly + lz*lz;
    if( d > 0 )
    {
      d = 1.0 / sqrtf( d );
    }
    else
    {
      d = 0.0;
    }

    lx *= d;
    ly *= d;
    lz *= d;

    //  dot product with unit normal vector

    float light = lx*ray.inx + ly*ray.iny + lz*ray.inz;

    light *= light*light;
    light *= light*light;

    col = rt_col_add( col, 0xFFFFFFFF, light );
  }

  //  darken colour if there's a shadow on this pixel

//  if( sob ) col = (col>>1) & 0x7F7F7F7F;

  return col;
}



//---- rt_render_area

//  Just a basic (slow) render

void rt_render_area( RT_SCENE *scene, int x1, int y1, int x2, int y2, RT_VIEW *view )
{
  int i, j;

  if( x1 < 0 ) x1 = 0;
  if( x1 >= view->w ) return;
  if( x2 > view->w ) x2 = view->w;
  if( x2 <= 0 ) return;

  if( y1 < 0 ) y1 = 0;
  if( y1 >= view->h ) return;
  if( y2 > view->h ) y2 = view->h;
  if( y2 <= 0 ) return;

  for( j=y1; j<y2; j++ )
    for( i=x1; i<x2; i++ )
      view->buf[ i + j * view->w ] = rt_scene_pixel( view, scene, i, j );
}



//---- ACCELERATED RENDERING FUNCTIONS ---------------------------------------//

/*

  Renders pixels in top corner of 4x4 squares

*/


//---- rt_render_sq_interp

void rt_render_sq_interp( RT_VIEW *view, int offset )
{
  int *buf = view->buf;
  int w4 = view->w4;
  int w = view->w;

  int c00 = buf[ offset ];
  int c01 = buf[ offset + 4 ];
  int c10 = buf[ offset + view->w4 ];
  int c11 = buf[ offset + view->w4 + 4];

  int i, j;

  for( j=0; j<4; j++ )
  {
    for( i=0; i<4; i++ )
    {
      int c0 = rt_col_blend( c00, c01, 0.25*i );
      int c1 = rt_col_blend( c10, c11, 0.25*i );

      int col = rt_col_blend( c0, c1, 0.25*j );

      buf[ offset + i + j*w ] = col;
    }
  }

  /*
  int w2 = view->w2;

  buf[ offset + 2 ] = RT_COL_MIDDLE( c00, c01 );
  buf[ offset + w2 ] = RT_COL_MIDDLE( c00, c10 );
  buf[ offset + w2 + 4 ] = RT_COL_MIDDLE( c01, c11 );
  buf[ offset + w2 + 2 ] = RT_COL_MIDDLE( c00, c11 );
*/
  /*


    *.o.+
    .....
    o.o..
    .....
    +...+

  */

}



//---- rt_render_sq

//  Renders a 4x4 square in the given colour at the offset in view->buf

//  NOTE: Do this inline for maximum efficiency

void rt_render_sq( RT_VIEW *view, int offset, int col )
{
  int wm4 = view->w - 4;

  offset++;//view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  offset += wm4;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  offset += wm4;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  offset += wm4;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
  view->buf[ offset++ ] = col;
}



//---- rt_render_h4

void rt_render_h4( RT_SCENE *scene, RT_VIEW *view, int offset, int i, int j )
{
  unsigned int *obj = view->obj;
  int *buf = view->buf;
/*
  int obj10 = obj[ offset + 1 ] = rt_scene_pixel_sphere( view, scene, i+1, j );
  int obj20 = obj[ offset + 2 ] = rt_scene_pixel_sphere( view, scene, i+2, j );
  int obj30 = obj[ offset + 3 ] = rt_scene_pixel_sphere( view, scene, i+3, j );

  RT_SPHERE *sphere = scene->sphere;

  int col10 = buf[ offset + 1 ] = sphere[ obj10 ].col;
  int col20 = buf[ offset + 2 ] = sphere[ obj20 ].col;
  int col30 = buf[ offset + 3 ] = sphere[ obj30 ].col;
*/

#if RT_SUPERPIXEL_SIZE == 4

  buf[ offset + 1 ] = rt_scene_pixel( view, scene, i+1, j );
  buf[ offset + 2 ] = rt_scene_pixel( view, scene, i+2, j );
  buf[ offset + 3 ] = rt_scene_pixel( view, scene, i+3, j );

#endif
}



//---- rt_render_v4

void rt_render_v4( RT_SCENE *scene, RT_VIEW *view, int offset, int i, int j )
{
  unsigned int *obj = view->obj;
  int *buf = view->buf;
/*
  int obj10 = obj[ offset + view->w ] = rt_scene_pixel_sphere( view, scene, i, j+1 );
  int obj20 = obj[ offset + view->w2 ] = rt_scene_pixel_sphere( view, scene, i, j+2 );
  int obj30 = obj[ offset + view->w3 ] = rt_scene_pixel_sphere( view, scene, i, j+3 );

  RT_SPHERE *sphere = scene->sphere;

  int col10 = buf[ offset + view->w ] = sphere[ obj10 ].col;
  int col20 = buf[ offset + view->w2 ] = sphere[ obj20 ].col;
  int col30 = buf[ offset + view->w3 ] = sphere[ obj30 ].col;
*/

#if RT_SUPERPIXEL_SIZE == 4
  buf[ offset + view->w  ] = rt_scene_pixel( view, scene, i, j+1 );
  buf[ offset + view->w2 ] = rt_scene_pixel( view, scene, i, j+2 );
  buf[ offset + view->w3 ] = rt_scene_pixel( view, scene, i, j+3 );
#endif

}



//---- rt_render

//  Renders full scene using optimisations for realtime performance

void rt_render( RT_SCENE *scene, RT_VIEW *view )
{
  int i, j, k, l;

  //- Precalc

  int w2 = view->w*2;
  int w3 = view->w*3;
  int w4 = view->w*4;
  int wm4 = view->w - 4;
  int hm4 = view->h - 4;

  int w = view->w;
  int h = view->h;

  rt_doScanlinePrecalc( view, scene );

  //  Debug: clear everything

  memset( (void *)view->buf, 0, view->w*view->h*sizeof(int) );
  memset( (void *)view->obj, 0, view->w*view->h*sizeof(int) );

  //  Set up ray

  RT_RAY ray;

  ray.ox = 0;
  ray.oy = 0;
  ray.oz = 0;
  ray.dz = 1;

  ray.dy = view->by1;

  //  Set up vars

  int offset;

  int *buf = view->buf;
  unsigned int *obj = view->obj;
  unsigned int *sob = view->sob;

  unsigned int col;

  //  Do 4x4 subsample first

  offset = 0;

  rt_scanline_firstRow( scene );

  for( j=0; j<h; j+=4 )
  {
    ray.dx = view->bx1;

    //  Do horizontal line as subsample
    //  This will also create linked lists attached to each superpixel of spheres
    //    which have been intersected

    for( i=0; i<w; i+=4 )
    {
      RY_RAY_FIRST_PRECALC( ray )

      view->obj[ offset ] = k = rt_intersect_first( scene, &ray );
      register int sob = view->sob[ offset ] = rt_intersect_shadowRay( scene, &ray );

      //- plot pixel

      col = scene->sphere[ k ].col;

      //  v0012: work out lighting

      //  unit vector towards light

      if( ! sob )
      {
        double lx = scene->light.x - ray.ix;
        double ly = scene->light.y - ray.iy;
        double lz = scene->light.z - ray.iz;

        double d = lx*lx + ly*ly + lz*lz;
        if( d > 0 )
        {
          d = 1.0 / sqrtf( d );
        }
        else
        {
          d = 0.0;
        }

        lx *= d;
        ly *= d;
        lz *= d;

        //  dot product with unit normal vector

        float light = lx*ray.inx + ly*ray.iny + lz*ray.inz;

        light *= light*light;
        light *= light*light;

        col = rt_col_add( col, 0xFFFFFFFF, light );
      }

      //  darken colour if there's a shadow on this pixel

//      if( sob ) col = (col>>1) & 0x7F7F7F7F;



      view->buf[ offset ] = col;




      ray.dx += view->bdx*4.0f;
      offset += 4;

      if( i < w-4 ) rt_scanline_nextColumn( scene, i+4 );
    }

    offset -= w;

    //  Move onto next 4x scanline

    ray.dy += view->bdy*4.0f;
    offset += w4;

    //  Move sphere list on to next row

    if( j < h-4 ) rt_scanline_nextRow( scene, j+4 );
  }

  //  Fill in any which are solid

  RT_SPHERE *sphere = scene->sphere;

  offset = 0;

  rt_scanline_firstRow( scene );

  for( j=0; j<hm4; j+=4 )
  {
    for( i=0; i<wm4; i+=4 )
    {
      //  Do 4x4 square if necessary

      unsigned int obj00 = obj[ offset ];
      unsigned int col00 = buf[ offset ];
      unsigned int sob00 = sob[ offset ];

      //- Render object's colour
      
      if( obj00 == obj[ offset+w4 ] && obj00 == obj[ offset+4 ] && obj00 == obj[ offset+w4+4 ] )
      {
        //  Same object at each corner - render rectangle

        if( sob00 == sob[ offset+w4 ] && sob00 == sob[ offset+4 ] && sob00 == sob[ offset+w4+4 ] )
        {
          //  Same shadow object at all four corners
//          rt_render_sq( view, offset, (col00 >> 1) & 0x7F7F7F7F );//  DEBUG: dark version
          
          rt_render_sq_interp( view, offset );

//          rt_render_sq( view, offset, col00 );
        }
        else
        {
          //  Different shadow object
          rt_render_v4( scene, view, offset, i, j );

          rt_render_h4( scene, view, offset   , i, j   );
          rt_render_h4( scene, view, offset+w , i, j+1 );
          rt_render_h4( scene, view, offset+w2, i, j+2 );
          rt_render_h4( scene, view, offset+w3, i, j+3 );
        }
      }
      else
      {
        rt_render_v4( scene, view, offset, i, j );

        rt_render_h4( scene, view, offset   , i, j   );
        rt_render_h4( scene, view, offset+w , i, j+1 );
        rt_render_h4( scene, view, offset+w2, i, j+2 );
        rt_render_h4( scene, view, offset+w3, i, j+3 );
      }

      offset += 4;

      rt_scanline_nextColumn( scene, i+4 );
    }

    offset += w3+4;

    rt_scanline_nextRow( scene, j+4 );
  }


#ifdef RT_DEBUG_RENDER_BAR_GRAPHS

  //  Render side bar graph

  int count=0;

  for( j=0; j<h>>2; j++ )
  {
    int old_count = count;

    count=0;

    //  Increase count for every sphere we come to
    for( i = y_spheres_start[j].header.first; i != -1; i = y_spheres_start[i].node.next )
    {
      count++;
      rt_setPixelAbs( view, w-50-count, j, 0xFFCC99 );
    }

    count=0;

    //  Decrease count for every sphere we are leaving
    for( i = y_spheres_end  [j].header.first; i != -1; i = y_spheres_end  [i].node.next )
    {
      count--;
      rt_setPixelAbs( view, w-50-count, j, 0x99CCFF );
    }

    //  Draw bar on left
//    for( i=0; i<count; i++ ) rt_setPixelAbs( view, i, j, 0xFFCC00 );

    //  show difference on right
  //  for( i=0; i<count-old_count; i++ ) rt_setPixelAbs( view, view->w-16-i, j, 0xFFFFFF );  // count greater than old_count
  //  for( i=0; i<old_count-count; i++ ) rt_setPixelAbs( view, view->w-16-i, j, 0xFF0000 );  // count less than old_count
  }

  //  Render top bar graph

  count=0;

  for( j=0; j<w>>2; j++ )
  {
    int old_count = count;

    count = 0;
    //  Increase count for every sphere we are leaving
    for( i = x_spheres_start[j].header.first; i != -1; i = x_spheres_start[i].node.next )
    {
      count++;
      rt_setPixelAbs( view, j, h-50-count, 0xFFCC99 );
    }

    count = 0;

    //  Decrease count for every sphere we are leaving
    for( i = x_spheres_end  [j].header.first; i != -1; i = x_spheres_end  [i].node.next )
    {
      count--;
      rt_setPixelAbs( view, j, h-50-count, 0x99CCFF );
    }

    //  Draw totals barchart up top
//    for( i=0; i<count; i++ ) rt_setPixelAbs( view, j, i, 0x6699CC );

    //  show difference on right
//    for( i=0; i<count-old_count; i++ ) rt_setPixelAbs( view, j, view->h-16-i, 0xFFFFFF );  // count greater than old_count
//    for( i=0; i<old_count-count; i++ ) rt_setPixelAbs( view, j, view->h-16-i, 0xFF0000 );  // count less than old_count
  }
// 

#endif

}



//---- DEBUG FUNCTIONS -------------------------------------------------------//


//---- rt_setPixelAbs

//  x and y are [-1..1]

void rt_setPixelAbs( RT_VIEW *view, int x, int y, int col )
{
  if( x < 0 || x >= view->w || y < 0 || y >= view->h ) return;

  view->buf[ x + y * view->w ] = col;
}



//---- rt_setPixel

//  x and y are [-1..1]

void rt_setPixel( RT_VIEW *view, float x, float y, int col )
{
  int i = (x-view->bx1)/view->bdx;
  int j = (y-view->by1)/view->bdy;

  if( i < 0 || i >= view->w || j < 0 || j >= view->h ) return;

  view->buf[ i + j * view->w ] = col;
}



//---- rt_drawFrameRate

void rt_drawFrameRate( RT_VIEW *view )
{
  int i, j;

  //  Frame rate indicator

  static int frame = 0;
  frame ++;

  for(j=0; j<8; j++) for(i=0; i<8; i++)
  {
    view->buf[j*view->w+i] = ((-(((i+j+frame)>>2)&1)) ^ 0xFFFF00 );
  }

  //  FPS based on current frame

  static int lastTime = 0;
  int time = timeGetTime();
  float fps = time > lastTime ? (1000.0f / (time-lastTime)) : 0.0f;
  lastTime = time;

  static int startTime = 0;
  if( startTime == 0 ) startTime = time;

  //- Frame rate based on last 10 frames

  static int times[10] ={0,0,0,0,0,0,0,0,0,0};
  static int timeptr = 0;

  int time10 = time-times[timeptr];

  times[timeptr] = time;
  timeptr = (timeptr+1)%10;

  //- Frame rate based on last second

  static float fps1sec = 0;
  static int fps1secLast = 0;
  static int fps1secFrame = 0;

  if( time > fps1secLast+1000 )
  {
    fps1sec = ( frame-fps1secFrame ) / (0.001f * (time-fps1secLast) );

    fps1secLast  = time;
    fps1secFrame = frame;
  }

  //- Display string

  char s[64] = "hello";

  sprintf( s, "%d %9.03f %9.03f %9.03f %9.03f",
              frame, (time-startTime)*0.001f, fps, fps1sec, 10000.0f / time10 );

  text_draw( view->buf, view->w, view->h ,9, 1, s, 0xFFFFFF );
}



//---- rt_drawHUD

void rt_drawHUD( RT_SCENE *scene, RT_VIEW *view )
{
  int i;

  float sx, sy;
  int x, y;

  //- Sphere centres

  RT_SPHERE *sphere = scene->sphere;

  for( i=0; i<scene->numSpheres; i++ )
  {
    rt_viewPoint( scene, view, sphere->ax, sphere->ay, sphere->az, &sx, &sy );
    x = (int)sx;
    y = (int)sy;

    if( x >= 0 && x < view->w && y >= 0 && y < view->h )
    {
      view->buf[ x + y * view->w ] = 0xFF00FF;
    }

    sphere ++;
  }
}



//---- END OF FILE -----------------------------------------------------------//
