/*----------------------------------------------------------------------+
| Spring NASA R&D Ray Tracing Project					|
| 3dda.c								|
|    Identify all cubes (seads) through which a ray passes.		|
+----------------------------------------------------------------------*/

#include <math.h>
#include "rtd.h"
#include "extern.h"
#define SIGN(x) (x > 0) ? 1 : ((x < 0) ? -1 : 0)
  
int box_x,box_y,box_z;			/* current seads coordinates 	*/
int sign_x,sign_y,sign_z;		/* seads step directions	*/
float ctrl_x, ctrl_y, ctrl_z;		/* control terms		*/
float slope_x,slope_y,slope_z;		/* unit step in driving axis	*/
int index;
int axis[3];

struct obj_link *three_dda(r, x, y, z)
/*----------------------------------------------------------------------+
| Sets up the control terms and step data which will be used by		|
| three_next(). 							|
+----------------------------------------------------------------------*/
struct ray *r;
int *x, *y, *z;

{
float org_x,org_y,org_z;		/* origin of the ray		*/
float seads_x,seads_y,seads_z;		/* seads that contains the org	*/
float abs_x,abs_y,abs_z;		/* slope terms			*/
struct obj_link *temp, *three_next();	

org_x  = r->org.x;			/* get the origin of the ray	*/
org_y  = r->org.y;
org_z  = r->org.z;
					/* The next section(s) of code  */
					/* clip the origin to start	*/
					/* within the seads space, if	*/
					/* possible.			*/

if (r->dir.x != 0.0) {			/* clip to x direction		*/
  if (org_x < x_world_start) {
    org_y = org_y + (r->dir.y / r->dir.x) * (x_world_start - org_x);
    org_z = org_z + (r->dir.z / r->dir.x) * (x_world_start - org_x);
    org_x = x_world_start;
    }
  else if (org_x > x_world_end) {
    org_y = org_y + (r->dir.y / r->dir.x) * (x_world_end - org_x);
    org_z = org_z + (r->dir.z / r->dir.x) * (x_world_end - org_x);
    org_x = x_world_start;
    }
  }
if (r->dir.y != 0.0) {			/* clip to y direction		*/
  if (org_y < y_world_start) {
    org_x = org_x + (r->dir.x / r->dir.y) * (y_world_start - org_y);
    org_z = org_z + (r->dir.z / r->dir.y) * (y_world_start - org_y);
    org_y = y_world_start;
    }
  else if (org_y > y_world_end) {
    org_x = org_x + (r->dir.x / r->dir.y) * (y_world_end - org_y);
    org_z = org_z + (r->dir.z / r->dir.y) * (y_world_end - org_y);
    org_y = y_world_end;
    }
  }
if (r->dir.z != 0.0) {			/* clip to z direction		*/
  if (org_z < z_world_start) {
    org_x = org_x + (r->dir.x / r->dir.z) * (z_world_start - org_z);
    org_y = org_y + (r->dir.y / r->dir.z) * (z_world_start - org_z);
    org_z = z_world_start;
    }
  else if (org_z > z_world_end) {
    org_x = org_x + (r->dir.x / r->dir.z) * (z_world_start - org_z);
    org_y = org_y + (r->dir.y / r->dir.z) * (z_world_start - org_z);
    org_z = z_world_end;
    }
  }

seads_x = (org_x - x_world_start) / num_pts_per_sead;	/* Identify starting */
seads_y = (org_y - y_world_start) / num_pts_per_sead;	/* seads box	     */
seads_z = (org_z - z_world_start) / num_pts_per_sead;
box_x  = seads_x;
box_y  = seads_y;
box_z  = seads_z;
sign_x = SIGN(r->dir.x); abs_x = fabs(r->dir.x);	/* Determine step    */
sign_y = SIGN(r->dir.y); abs_y = fabs(r->dir.y);	/* factor for each   */
sign_z = SIGN(r->dir.z); abs_z = fabs(r->dir.z);	/* direction	     */

slope_x = slope_y = slope_z = 0.0;
ctrl_x = ctrl_y = ctrl_z = 0.0;
					/* If the ray travels most in the   */
					/* x direction, then let the x axis */
					/* be the driving axis and compute  */
					/* the y and z control terms.	    */
if (abs_x >= abs_y && abs_x >= abs_z) {		
  if (sign_y != 0) {
    slope_y = abs_x / abs_y;
    ctrl_y = seads_x + (box_y - seads_y) * r->dir.x / r->dir.y - box_x;
    ctrl_y = (sign_x > 0) ? -ctrl_y : ctrl_y - 1;
    if (sign_y > 0)
      ctrl_y= ctrl_y - slope_y;
    }
  if (sign_z != 0) {  
    slope_z= abs_x / abs_z;
    ctrl_z = seads_x + (box_z - seads_z) * r->dir.x / r->dir.z - box_x;
    ctrl_z = (sign_x > 0) ? -ctrl_z : ctrl_z - 1;
    if (sign_z > 0)
      ctrl_z= ctrl_z - slope_z;
    }
  }
else if (abs_y >= abs_x && abs_y >= abs_z) {	/* same as above, except */
  if (sign_x != 0) {				/* the y axis is the	 */
    slope_x= abs_y / abs_x;			/* driving axis...	 */
    ctrl_x = seads_y + (box_x - seads_x) * r->dir.y / r->dir.x - box_y;
    ctrl_x = (sign_y > 0) ? -ctrl_x : ctrl_x - 1;
    if (sign_x > 0)
      ctrl_x = ctrl_x - slope_x;
    }
  if (sign_z != 0) {
    slope_z= abs_y / abs_z;
    ctrl_z = seads_y + (box_z - seads_z) * r->dir.y / r->dir.z - box_y;
    ctrl_z = (sign_y > 0) ? -ctrl_z : ctrl_z - 1;
    if (sign_z > 0)
      ctrl_z = ctrl_z - slope_z;
    }
  }
else {						/* same as above with	 */
  if (sign_x != 0) {				/* z as the driving axis */
    slope_x= abs_z / abs_x;
    ctrl_x = seads_z + (box_x - seads_x) * r->dir.z / r->dir.x - box_z;
    ctrl_x = (sign_z > 0) ? -ctrl_x : ctrl_x - 1;
    if (sign_x > 0)
      ctrl_x = ctrl_x - slope_x;
    }
  if (sign_y != 0) {
    slope_y= abs_z / abs_y;
    ctrl_y = seads_z + (box_y - seads_y) * r->dir.z / r->dir.y - box_z;
    ctrl_y = (sign_z > 0) ? -ctrl_y : ctrl_y - 1;
    if (sign_y > 0)
      ctrl_y = ctrl_y - slope_y;
    }
  }
index = 0;					/* clear the index counter */
if (box_x < 0 || box_y < 0 || box_z < 0 ||
    box_x > NUM_SEADS || box_y > NUM_SEADS || box_z > NUM_SEADS)
  return(NULL);					/* the ray is outside seads */
else 
  if (seads[box_x][box_y][box_z] != NULL)
    {
    *x = box_x;
    *y = box_y;
    *z = box_z;
    return(seads[box_x][box_y][box_z]);		/* the origin sead contains */
    }
						/* an object		    */
  else 
  {
  temp = three_next(x, y, z);		/* find the next sead which */
  return(temp);					/* contains an object and   */
						/* return it.		    */
  }
}
  
struct obj_link *three_next(x, y, z)
int *x, *y, *z;
/*----------------------------------------------------------------------+
| Step through the seads along the ray initialized by three_dda() 	|
| above. When we come across a sead with some objects in it, return	|
| a pointer to a list of those objects. Otherwise, keep stepping until  |
| we leave the mesh and return NULL.					|
+----------------------------------------------------------------------*/
{
int temp;

while (1) {
  if (index == 0) {
    axis[1]= (ctrl_x < ctrl_y) ? ((ctrl_z < ctrl_x) ? 2 : 0)
			       : ((ctrl_z < ctrl_y) ? 2 : 1);
    axis[2]= (ctrl_x > ctrl_y) ? ((ctrl_z > ctrl_x) ? 2 : 0)
			       : ((ctrl_z > ctrl_y) ? 2 : 1);
    axis[0]= 3 - axis[2] - axis[1];
    }
  temp = axis[index];
  index= (index + 1) % 3;
  switch (temp) {
    case 0: if (++ctrl_x >= 0.0 && sign_x != 0) {	/* step along x axis */
    	      ctrl_x= ctrl_x - slope_x;
    	      box_x = box_x + sign_x;
              *x = box_x;
              *y = box_y;
              *z = box_z;
    	      if (box_x < 0 || box_x > NUM_SEADS)
    	        return(NULL);
    	      if (seads[box_x][box_y][box_z] != NULL)
    	        return(seads[box_x][box_y][box_z]);
    	      }
    	    break;
    case 1: if (++ctrl_y >= 0.0 && sign_y != 0) {	/* step along y axis */
              ctrl_y= ctrl_y - slope_y;
              box_y = box_y + sign_y;
              *x = box_x;
              *y = box_y;
              *z = box_z;
              if (box_y < 0 || box_y > NUM_SEADS)
                return(NULL);
              if (seads[box_x][box_y][box_z] != NULL)
                return(seads[box_x][box_y][box_z]);
              }
            break;
    case 2: if (++ctrl_z >= 0.0 && sign_z != 0) {	/* step along z axis */
              ctrl_z= ctrl_z - slope_z;
              box_z = box_z + sign_z;
              *x = box_x;
              *y = box_y;
              *z = box_z;
              if (box_z < 0 || box_z > NUM_SEADS)
                return(NULL);
              if (seads[box_x][box_y][box_z] != NULL)
                return(seads[box_x][box_y][box_z]);
              }
            break;
    }
  }
}