#include <math.h>
#include "rtd.h"
#include "extern.h"
#include "macros.h"

int d_bug = FALSE;

double  findo (m, s) /* finds where a ray inside the ball exits. */
struct mat *m;
struct sphere_data  *s;
{
                           /* foops id the rotated position vector. */
    struct vector   foops;
    double  t;
    MTV (foops, (*m), s -> cent);

                          /* see if it hits the ball (it better)*/
    t = s -> rad * s -> rad - foops.y * foops.y - foops.z * foops.z;
    if (t > 0)
	t = foops.x + sqrt (t);
    else
	t = 0;
                          /* return how far along the ray you 
                             were when you hit */
    return (t);
}

double  find (m, s)/* finds whether a ray hits a ball*/
struct mat *m;
struct sphere_data  *s;
{
    struct vector   foops;
    double  t;
    MTV (foops, (*m), s -> cent);
    t = s -> rad * s -> rad - foops.y * foops.y - foops.z * foops.z;
    if (t > 0)
	t = foops.x - sqrt (t);
    else
	t = 0;
    return (t);
}

double  finds (m, s)/* finds if a ball is between a point and a 
			lightsource. Returns how obscuring the ball is */
struct mat *m;
struct sphere_data  *s;
{
    struct vector   foops;
    double  t;
    MTV (foops, (*m), s -> cent);
    t = s -> rad - sqrt (foops.y * foops.y + foops.z * foops.z);
    if (t > 0)
	t = t / foops.x;
    else
	t = 0;
    return (t);
}




double  shadow (p,objid)/* finds if a point is in a shadow, or if it is on edge */
struct vector  *p;
int objid;
{
    struct mat  trans;
    struct sphere_data   ss;
    struct vector   d, temp, new;
    struct ray r;
    struct obj_link *ptr, *three_dda(), *three_next();
    struct object_data *c;
    struct vector intersection_pts;
    int     i;
    double  l,
            k,
            x,
            y,
            z,
            find_facet(),
            finds ();
    int sead_x, sead_y, sead_z;
    double bound_x, bound_y, bound_z;

    l = 0.0;
    c = NULL;
    SV (d, light_source.cent, (*p));
    d.l = LEN (d);
    d.xzl = XZL (d);
    mt (&(d), &trans);

    MV(p->x, p->y, p-> z, r.org);
    MV(d.x, d.y, d.z, r.dir);
    ptr = three_dda(&r, &sead_x, &sead_y, &sead_z);
    r.src = objid;
    while (ptr != NULL) {
      while (ptr != NULL) 
        {
        if (ptr->the_object->type == SPHERE_TYPE)
          { 
	  ss.rad = ptr->the_object->object.sphere.rad;
	  SV (ss.cent, ptr->the_object->object.sphere.cent, (*p));
	  if ((k = finds (&trans, &ss)) > l) 
            {
	    c = ptr->the_object;
            l = k;
	    }  /* k > l */
          }   /* end of if SPHERE_TYPE */
          else{
            if((k=find_facet(ptr->the_object,&r,&intersection_pts))>0.0)
               {
               if (0.9>l) {l = 0.9; c = ptr->the_object;}
               /* facets are very occlusive in this code, snicker snicker */
               }
          }
        ptr = ptr->next_obj;
        }     /* end of inner while */
      ptr = three_next(&sead_x, &sead_y, &sead_z);
    }         /* end of outer while */
    if (c == NULL)
	k = 200.0;
    else {
	k = 1.0 - l / ((light_source.rad) / (d.l));
	if (k < 0.0)
	    k = 0.0;
	k *= 200.0;
    }
    return (k);
}

/* find if a vector intersects a triangle  in three space */

double find_facet(ptr, r, intersect_pt)
struct object_data *ptr;
struct ray *r;
struct vector *intersect_pt;

{
  
double x, y, z;
double t, d, x5, y5, z5, divisor;
struct vector v1, v2, vec1, vec2, vec3, norm, new;

  if (ptr->id == r->src) return(0.0);
  /* if the ray began at the current object then we should have ret'd above */
  MV(ptr->object.facet.point[0].x, ptr->object.facet.point[0].y, 
     ptr->object.facet.point[0].z, vec1); 
  MV(ptr->object.facet.point[1].x, ptr->object.facet.point[1].y, 
     ptr->object.facet.point[1].z, vec2); 
  MV(ptr->object.facet.point[2].x, ptr->object.facet.point[2].y, 
     ptr->object.facet.point[2].z, vec3); 


  MV(ptr->object.facet.normal.x, ptr->object.facet.normal.y,
     ptr->object.facet.normal.z, norm);
  norm.xzl = ptr->object.facet.normal.xzl;

  norm.l = ptr->object.facet.normal.l; 
  MV(r->dir.x,r->dir.y,r->dir.z,new);
  new.l = LEN(r->dir);
 if (acos(DOT(norm,new)/(norm.l*new.l)) <= (double)(3.14159265358/2.0)) return(0.0);
 /*if angle twixt ray and normal is greater than pi/2 we should have returned */

  SCMLT((1.0/norm.l), norm);

  divisor = DOT(norm, r->dir);
  if (divisor != 0)
    {
    t = (DOT(norm, vec1) - DOT(norm, r->org)) / divisor;
    MV(r->dir.x, r->dir.y, r->dir.z, new);
    new.l = LEN(new);
    SCMLT(t, new);
    t = new.l;
    AV((*intersect_pt), r->org, new);

    if (point_inclusion(intersect_pt, &vec1, 
         &vec2, &vec3, ptr->object.facet.normal))
           /* we have an intersection */
      {
	if (d_bug)
	   {
    	    printf("FACET INTERSECTION t= %5.2f\n", t);
      	    printf("   i_pt x = %5.2f  y = %5.2f  z = %5.2f\n",
            intersect_pt->x, intersect_pt->y, intersect_pt->z);
	   }
      return(t);
      }
    else
      return(0.0);
    }
  else
    return(0.0);
} /* end of find_facet */


int point_inclusion(intersect_pt, vec1, vec2, vec3, normal)
struct vector *intersect_pt, *vec1, *vec2, *vec3, normal;
{
  struct vector s1, check_vec, new_norm;
  float dir;

  SV(s1, (*vec1), (*vec2))
  CROSS_PRODUCT(new_norm, s1, normal);
  SV(check_vec, (*intersect_pt), (*vec1))
  dir = DOT(check_vec, new_norm);
  if (dir < 0)
    return(FALSE);

  SV(s1, (*vec2), (*vec3))
  CROSS_PRODUCT(new_norm, s1, normal);
  SV(check_vec, (*intersect_pt), (*vec2))
  dir = DOT(check_vec, new_norm);
  if (dir < 0)
    return(FALSE);

  SV(s1, (*vec3), (*vec1))
  CROSS_PRODUCT(new_norm, s1, normal);
  SV(check_vec, (*intersect_pt), (*vec3))
  dir = DOT(check_vec, new_norm);
  if (dir < 0)
    return(FALSE);

  return(TRUE);
}