/**********************************************************

  Line/Object intersection routines.  Routine returns
  TRUE if object is hit, along with the nearest of the
  possibly multiple intersections.  First parameter is
  line, second is object, and third is pointer to parameter
  on line for intersection (filled by routine).  Functions
  are pointed at by entries in ObjData structure.

 **********************************************************/

#include "qrt.h"

/* #define INTERSECTDEBUG */

/**********************************************************

          Line/bounding box intersection test.

   Looks bad, but its pretty fast.  Many times we don't
   even have to go through the whole routine if a miss
   can be detected early.

 **********************************************************/

int LineBbox(line, bbox, t)
  OBJ_PTR line, bbox;
  float *t;
{
  register float tminx, tmaxx, tminy, tmaxy,
                 tminz, tmaxz, tmin, tmax, t1,t2;

  *t=10;

# ifdef ROBUST
    if (line->type!=LINE) Error(INTERNAL_ERROR,301);
    if (bbox->type!=BBOX) Error(INTERNAL_ERROR,302);
# endif

  if (fabs(line->vect1.x) < SMALL) {
    if ((bbox->lower.x < line->loc.x) &&
        (bbox->upper.x > line->loc.x)) {
          tminx = -3e30; tmaxx=3e30;
        } else return(FALSE);

  } else {
    t1 = (bbox->lower.x-line->loc.x)/line->vect1.x;
    t2 = (bbox->upper.x-line->loc.x)/line->vect1.x;
    tminx = MIN(t1,t2);
    tmaxx = MAX(t1,t2);
    if (tmaxx<0) return(FALSE);
  }

  if (fabs(line->vect1.y) < SMALL) {
    if ((bbox->lower.y < line->loc.y) &&
        (bbox->upper.y > line->loc.y)) {
          tminy = -3e30; tmaxy=3e30;
        } else return(FALSE);

  } else {
    t1 = (bbox->lower.y-line->loc.y)/line->vect1.y;
    t2 = (bbox->upper.y-line->loc.y)/line->vect1.y;
    tminy = MIN(t1,t2);
    tmaxy = MAX(t1,t2);
    if (tmaxy<0) return(FALSE);
  }

  if (fabs(line->vect1.z) < SMALL) {
    if ((bbox->lower.z < line->loc.z) &&
        (bbox->upper.z > line->loc.z)) {
          tminz = -3e30; tmaxz=3e30;
        } else return(FALSE);

  } else {
    t1 = (bbox->lower.z-line->loc.z)/line->vect1.z;
    t2 = (bbox->upper.z-line->loc.z)/line->vect1.z;
    tminz = MIN(t1,t2);
    tmaxz = MAX(t1,t2);
    if (tmaxz<0) return(FALSE);
  }

  tmin = MAX(MAX(tminx,tminy),tminz);
  tmax = MIN(MIN(tmaxx,tmaxy),tmaxz);

# ifdef INTERSECTDEBUG
    printf("LINEBBOX: dir = %f %f %f\n",line->vect1.x,
                                        line->vect1.y,
                                        line->vect1.z);

    printf("  tminx=%7.2f, tmaxx=%7.2f, tminy=%7.2f, tmaxy=%7.2f\n",
            tminx, tmaxx, tminy, tmaxy);
    printf("  tminz=%7.2f, tmaxz=%7.2f, tmin=%7.2f,  tmax=%f\n",
            tminz, tmaxz, tmin, tmax);
# endif

  if (tmax<0) return(FALSE);

  if (tmax<tmin) return(FALSE);

# ifdef INTERSECTDEBUG
    printf("  HIT:\n");
# endif

  return(TRUE);
}

/**********************************************************

               Line/ring intersection test

     Similar to, but slower than parallelogram test

 **********************************************************/

int LineRing(line, ring, t)
  OBJ_PTR line, ring;
  float *t;
{
  VECTOR         delta, loc;
  register float dot, rad, pos1, pos2;

# ifdef ROBUST
    if (line->type!=LINE) Error(INTERNAL_ERROR,303);
    if (ring->type!=RING) Error(INTERNAL_ERROR,304);
# endif

  dot  = DotProd((ring->precomp.norm),line->vect1);

  if (fabs(dot)<SMALL) return(FALSE);

  pos1 = ring->precomp.n1;
  pos2 = DotProd((ring->precomp.norm),line->loc);

  *t=(pos1-pos2)/dot;

# ifdef INTERSECTDEBUG
    printf("LINERING: t=%f\n", *t);
# endif

  FindPos(&loc,line,*t);
  VecSubtract(&delta,&loc,&(ring->loc));

  pos1 = DotProd(delta,ring->vect1);
  pos2 = DotProd(delta,ring->vect2);

  rad = sqrt(sqr(pos1)+sqr(pos2));

# ifdef INTERSECTDEBUG
    printf("LINERING  pos1,2 = %f %f\n",pos1,pos2);
    printf("          radius = %f\n",rad);
    printf("          R1,2   = %f %f\n",ring->vect3.x,ring->vect3.y);
# endif

  if (rad<(ring->vect3.x) || rad>(ring->vect3.y))
     return(FALSE);

# ifdef INTERSECTDEBUG
    printf("  HIT:\n");
# endif

  return(TRUE);
}


/**********************************************************

        Line/Parallelogram intersection test
        Returns Parameter T for intersection

 **********************************************************/

int LineParallelogram(line, para, t)
  OBJ_PTR line, para;
  float *t;
{
  VECTOR         delta, loc;
  register float dot, in1, in2;

# ifdef ROBUST
    if (line->type!=LINE)          Error(INTERNAL_ERROR,305);
    if (para->type!=PARALLELOGRAM) Error(INTERNAL_ERROR,306);
# endif

  dot = DotProd((para->precomp.norm),line->vect1);

  if (fabs(dot)<SMALL) return(FALSE);

  in1 = para->precomp.n1;
  in2 = DotProd((para->precomp.norm),line->loc);

  *t=(in1-in2)/dot;

# ifdef INTERSECTDEBUG
    printf("LINEPARALL: t=%f\n", *t);
# endif

  FindPos(&loc,line,*t);
  VecSubtract(&delta,&loc,&(para->loc));

  in1  = DotProd(delta,para->vect1)/(para->precomp.len1);
  in2  = DotProd(delta,para->vect2)/(para->precomp.len2);

# ifdef INTERSECTDEBUG
    printf("LINEPARALL: in1,2 = %f %f\n",in1,in2);
# endif

  if (!((in1>=0) && (in2>=0) && (in1<=1) && (in2<=1)))
    return(FALSE);

# ifdef INTERSECTDEBUG
    printf("  HIT:\n");
# endif

  return(TRUE);
}


/**********************************************************

        Line/Triangle intersection test
        Returns Parameter T for intersection

 **********************************************************/

int LineTriangle(line, obj, t)
  OBJ_PTR line, obj;
  float *t;
{
  VECTOR         delta, loc;
  register float dot, in1, in2;

# ifdef ROBUST
    if (line->type != LINE)      Error(INTERNAL_ERROR,307);
    if (obj->type  != TRIANGLE)  Error(INTERNAL_ERROR,308);
# endif

  dot = DotProd((obj->precomp.norm),line->vect1);

  if (fabs(dot)<SMALL) return(FALSE);

  in1 = obj->precomp.n1;
  in2 = DotProd((obj->precomp.norm),line->loc);

  *t=(in1-in2)/dot;

# ifdef INTERSECTDEBUG
    printf("LINEPARALL: t=%f\n", *t);
# endif

  FindPos(&loc,line,*t);
  VecSubtract(&delta,&loc,&(obj->loc));

  in1  = DotProd(delta,obj->vect1)/(obj->precomp.len1);
  in2  = DotProd(delta,obj->vect2)/(obj->precomp.len2);

# ifdef INTERSECTDEBUG
    printf("LINETRIANGLE: in1,2 = %f %f\n",in1,in2);
# endif

  if (!((in1>=0) && (in2>=0) && (in1+in2<=1)))
    return(FALSE);

# ifdef INTERSECTDEBUG
    printf("  HIT:\n");
# endif

  return(TRUE);
}


/**********************************************************

           Line/sphere intersection test
        Returns parameter T for intersection

 **********************************************************/

int LineSphere(line, sph, t)
  OBJ_PTR line, sph;
  float *t;
{
  register float a,b,c,d,t1, tmpx,tmpy,tmpz;

# ifdef ROBUST
    if (line->type!=LINE) Error(INTERNAL_ERROR,309);
    if (!(sph->type==SPHERE || sph->type==LAMP))
      Error(INTERNAL_ERROR,310);
# endif

  tmpx = sph->loc.x-line->loc.x;
  tmpy = sph->loc.y-line->loc.y;
  tmpz = sph->loc.z-line->loc.z;

  c = sqr(tmpx)+ sqr(tmpy)+ sqr(tmpz) - (sph->precomp.n1);

  b = -2*(line->vect1.x*tmpx+                      /* find b */
          line->vect1.y*tmpy+
          line->vect1.z*tmpz);

  a = sqr(line->vect1.x)+                          /* find a */
      sqr(line->vect1.y)+
      sqr(line->vect1.z);

  d = sqr(b)-4.0*a*c;

# ifdef INTERSECTDEBUG
    printf("LINESPHERE: a=%f, b=%f, c=%f, d=%f\n",a,b,c,d);
# endif

  if (d<=0) return(FALSE);                          /* does sphere hit? */

  d=sqrt(d); *t=(-b+d)/(a+a);
             t1=(-b-d)/(a+a);

  if (t1<*t && t1>SMALL) *t=t1;                     /* find 1st collision */

  if (*t > SMALL) {

#   ifdef INTERSECTDEBUG
      printf("LINESPHERE: collision @ t=%f\n",*t);
#   endif

    return(TRUE);
  }

  return(FALSE);
}

/**********************************************************

        Line/quadratic intersection test
        Returns parameter T for intersection

  newline is the input line translated and rotated so that
  the quadratic is @ 0,0,0 and pointed up.

 **********************************************************/

int LineQuadratic(line, quad, t)
  OBJ_PTR line, quad;
  float *t;
{
  register   float a,b,c,d, t1;
  VECTOR     loc, loc1, tempdir;
  OBJ_STRUCT newline;

# ifdef ROBUST
    if (line->type!=LINE)      Error(INTERNAL_ERROR,311);
    if (quad->type!=QUADRATIC) Error(INTERNAL_ERROR,312);
# endif

  newline.type=LINE;

  /*** translation transformation to newpos for line ***/

  VecSubtract(&(newline.loc),&(line->loc),&(quad->loc));

  if ((quad->vect1.x == 0) &&                  /* no rotation  */
      (quad->vect1.y == 1) &&                  /* if aligned   */
      (quad->vect1.z == 0))   {

      VectEQ(&(newline.vect1),&(line->vect1));

  } else {                                     /* here we must rot */

#   ifdef INTERSECTDEBUG

      printf("LINEQUADRATIC 0\n");
      printf("  inline.loc  = %f %f %f\n",    line->loc.x,
                                              line->loc.y,
                                              line->loc.z);

      printf("  inline.dir  = %f %f %f\n",    line->vect1.x,
                                              line->vect1.y,
                                              line->vect1.z);

      printf("  newline.loc = %f %f %f\n",    newline.loc.x,
                                              newline.loc.y,
                                              newline.loc.z);
#   endif

    Rot12( &(line->vect1),&(newline.vect1),    /* rot view direction */
           quad->precomp.cos1,
           quad->precomp.sin1,
           quad->precomp.cos2,
           quad->precomp.sin2 );

    Rot12( &(newline.loc),&(newline.loc),      /* rotate view location */
           quad->precomp.cos1,
           quad->precomp.sin1,
           quad->precomp.cos2,
           quad->precomp.sin2 );

#   ifdef INTERSECTDEBUG

      printf("LINEQUADRATIC 1\n");
      printf("  inline.loc  = %f %f %f\n",    line->loc.x,
                                              line->loc.y,
                                              line->loc.z);

      printf("  inline.dir  = %f %f %f\n",    line->vect1.x,
                                              line->vect1.y,
                                              line->vect1.z);

      printf("  newline.loc = %f %f %f\n",    newline.loc.x,
                                              newline.loc.y,
                                              newline.loc.z);

      printf("  newline.dir = %f %f %f\n",    newline.vect1.x,
                                              newline.vect1.y,
                                              newline.vect1.z);
#   endif

  }

  c = -(quad->cterm) + quad->vect2.x * sqr(newline.loc.x) +
                       quad->vect2.y * sqr(newline.loc.y) +
                       quad->vect2.z * sqr(newline.loc.z);

  b = 2*( quad->vect2.x * newline.loc.x * newline.vect1.x +
          quad->vect2.y * newline.loc.y * newline.vect1.y +
          quad->vect2.z * newline.loc.z * newline.vect1.z);

  a = quad->vect2.x * sqr(newline.vect1.x) +
      quad->vect2.y * sqr(newline.vect1.y) +
      quad->vect2.z * sqr(newline.vect1.z);

  d = sqr(b)-4.0*a*c;

# ifdef INTERSECTDEBUG
    printf("LINEQUADRATIC 2: a=%f, b=%f, c=%f, d=%f\n",a,b,c,d);

    printf("               newpos = %f %f %f\n",newline.loc.x,
                                                newline.loc.y,
                                                newline.loc.z);

    printf("               newdir = %f %f %f\n",newline.vect1.x,
                                                newline.vect1.y,
                                                newline.vect1.z);
# endif

  if (d<0) return(FALSE);                           /* we missed it */


  d=sqrt(d); *t=(-b+d)/(a+a);
             t1=(-b-d)/(a+a);

  FindPos(&loc,&newline,*t);                        /* find locations */
  FindPos(&loc1,&newline,t1);


  if ((loc.x < quad->lower.x) ||                    /* 1st in range ? */
      (loc.x > quad->upper.x) ||
      (loc.y < quad->lower.y) ||
      (loc.y > quad->upper.y) ||
      (loc.z < quad->lower.z) ||
      (loc.z > quad->upper.z))    *t = -1;

  if ((loc1.x < quad->lower.x) ||                   /* 2nd in range ? */
      (loc1.x > quad->upper.x) ||
      (loc1.y < quad->lower.y) ||
      (loc1.y > quad->upper.y) ||
      (loc1.z < quad->lower.z) ||
      (loc1.z > quad->upper.z))   t1 = -1;

# ifdef INTERSECTDEBUG
    printf("    t,t1 = %f %f\n",*t,t1);
# endif

  if (*t<=SMALL && t1 <=SMALL) return(FALSE);

  if (*t<=SMALL && t1>SMALL) *t=t1;
  if (t1<*t && t1>SMALL)     *t=t1;   /* find 1st collision */

  if (*t>SMALL) {

#   ifdef INTERSECTDEBUG
      printf("LINEQUADRATIC: collision @ t=%f\n",*t);
#   endif

    return(TRUE);
  }

  return(FALSE);
}