#include <cmath>
#include <float.h>
#include "polymesh.h"

static bool
triIntersect( const ray &r, const vec3 &a, const vec3 &b, const vec3 &c,
              vec3 &bary, float &t, vec3 &normal );

Polymesh::~Polymesh()
{
    for( Materials::iterator i = materials.begin(); i != materials.end(); ++i )
    {
        delete *i;
    }
}

                                // must add vertices, normals, and materials IN ORDER
void
Polymesh::addVertex( const vec3 &v )
{
    vertices.push_back( v );
}

void
Polymesh::addMaterial( Material *m )
{
    materials.push_back( m );
}

void
Polymesh::addNormal( const vec3 &n )
{
    normals.push_back( n );
}

bool
Polymesh::addFace( int a, int b, int c )
// RETURNS false if the vertices a,b,c don't all exist    
{
    int vcnt = vertices.size();

    if( a >= vcnt || b >= vcnt || c >= vcnt )
        return false;

    faces.push_back( Triangle( a,b,c ) );
    return true;
}

char *
Polymesh::doubleCheck()
// Check to make sure that if we have per-vertex materials or normals
// they are the right number.
{
    if( materials.size() && materials.size() != vertices.size() )
        return "Bad Polymesh: Wrong number of materials.";
    if( normals.size() && normals.size() != vertices.size() )
        return "Bad Polymesh: Wrong number of normals.";

    return 0;
}
        
bool
Polymesh::intersect( const ray& r, isect& i ) const
{
    i.setT( DBL_MAX );
    isect ii;
    ii.obj = this;

    bool intersection = false;
    for( Faces::const_iterator fi = faces.begin(); fi != faces.end(); ++fi )
    {
        if( intersect( r, *fi, ii ) && ii.t < i.t )
        {
            i = ii;
            intersection = true;
        }
    }
    
    return intersection;
}

bool
Polymesh::intersect( const ray& r ) const
{
    for( Faces::const_iterator fi = faces.begin(); fi != faces.end(); ++fi )
        if( intersect( r, *fi ) )
            return true;
    return false;
}

bool
Polymesh::intersect( const ray& r, const Triangle &tri, isect& i ) const
{
    vec3 bary;
    float t;
    vec3 normal;
    
    if( !triIntersect( r, vertices[tri[0]], vertices[tri[1]], vertices[tri[2]],
                       bary, t, normal ) )
        return false;

    i.setT( t );
    if( normals.size() )
    {
                                // linearly interpolate vertex normals
        i.setN( (bary[0]*normals[tri[0]]+bary[1]*normals[tri[1]]+bary[2]*normals[tri[2]])
                .normalize() );
    }
    else
    {
        i.setN( normal );       // otherwise use face normal
    }
    if( materials.size() )
    {
        Material *m = new Material( Material::zero );
        for( int jj = 0; jj < 3; ++jj )
            (*m) += bary[jj] * (*materials[ tri[jj] ]);
        i.setMaterial( m );
    }
    
    return true;
}

bool
Polymesh::intersect( const ray& r, const Triangle &tri ) const
{
    vec3 bary;
    float t;
    vec3 normal;
    
    return triIntersect( r, vertices[tri[0]], vertices[tri[1]], vertices[tri[2]],
                         bary, t, normal );
}

static bool
triIntersect( const ray &r, const vec3 &a, const vec3 &b, const vec3 &c,
              vec3 &bary, float &t, vec3 &n )
// Intersect ray r with the triangle abc.  If it hits returns true,
// and put the parameter in t and the barycentric coordinates of the
// intersection in bary.
// Uses the algorithm and notation from _Graphic Gems 5_, p. 232.
//
// Calculates and returns the normal of the triangle too.
{
    vec3 p = r.getPosition();
    vec3 v = r.getDirection();
    
    vec3 ab = b - a;
    vec3 ac = c - a;
    vec3 ap = p - a;
    
    n = (ab^ac).normalize();

    double vdotn = v*n;
    // Wyvern removed this exitpoint.
    // In raycasting, we can safely cull the backfaces.  But in raytracing, we need
    // the backface intersections as well in order to render internal surfaces and
    // reflections.
    // if( -vdotn < RayGlobals::normal_epsilon )
    //     return false;
    
    t = - (ap*n)/vdotn;
    
    if( t < RayGlobals::ray_epsilon )
        return false;

                                // find k where k is the index of the component
                                // of normal vector with greatest absolute value
    float greatestMag = FLT_MIN;
    int k = -1;
    for( int i = 0; i < 3; ++i )
    {
        float val = n[i];
        if( val < 0 )
            val *= -1;
        if( val > greatestMag )
        {
            k = i;
            greatestMag = val;
        }
    }

    vec3 am = ap + t * v;

    bary[1] = (am^ac)[k]/(ab^ac)[k];
    bary[2] = (ab^am)[k]/(ab^ac)[k];
    bary[0] = 1-bary[1]-bary[2];
    if( bary[0] < 0 || bary[1] < 0 || bary[1] > 1 || bary[2] < 0 || bary[2] > 1 )
        return false;
    else
        return true;
}

void
Polymesh::generateNormals()
// Once you've loaded all the verts and faces, we can generate per
// vertex normals by averaging the normals of the neighboring faces.
{
    int cnt = vertices.size();
    normals.resize( cnt );
    int *numFaces = new int[ cnt ]; // the number of faces assoc. with each vertex
    memset( numFaces, 0, sizeof(int)*cnt );
    
    for( Faces::iterator fi = faces.begin(); fi != faces.end(); ++fi )
    {
        vec3 a = vertices[(*fi)[0]];
        vec3 b = vertices[(*fi)[1]];
        vec3 c = vertices[(*fi)[2]];
        
        vec3 faceNormal = ((b-a)^(c-a)).normalize();
        
        for( int i = 0; i < 3; ++i )
        {
            normals[(*fi)[i]] += faceNormal;
            ++numFaces[(*fi)[i]];
        }
    }

    for( int i = 0; i < cnt; ++i )
    {
        if( numFaces[i] )
            normals[i]  /= numFaces[i];
    }

    delete [] numFaces;
}