#ifndef TRIANGLE_HXX
#define TRIANGLE_HXX

#define POINTA 0
#define POINTB 1
#define POINTC 2

#include "Primitive.hxx"

class Triangle : public Primitive
{
	Vec3f point[3], a, b, c, p, origin, normal,
	      b_minus_a, c_minus_b, a_minus_c,
	      bisec[3],
	      ray_delta;
	float t, cos_angle2[3],
	      numerator, denominator;
	int   first, second;
public:
	// This is a optimized version of the triangle intersection
	// which is no longer basing on barycentric coordinates.
	// Instead, we check whether the ray hitpoint is in the fractum
	// spanned by the first angle. We proceed only if this is
	// the case and then quickly check whether we are in front of
	// or behind the limiting edge

	Triangle(Vec3f a, Vec3f b, Vec3f c, Shader *shader)
		: a(a),b(b),c(c)
	{
		this->shader = shader;

		// We may always operate with cosines of angles instead of
		// angles themselves, because of the monotonicity of the
		// cosine in [0,pi] and we only regard comparisons rather
		// real values.

		// Introduce array notation for points, for easy access
		point[POINTA] = a;
		point[POINTB] = b;
		point[POINTC] = c;

		// Normalize the edges of the triangle
		b_minus_a = b - a;
		c_minus_b = c - b;
		a_minus_c = a - c;

		Normalize(b_minus_a);
		Normalize(c_minus_b);
		Normalize(a_minus_c);
    
		// Get the bisec vectors
		bisec[POINTA] = (b_minus_a - a_minus_c) / 2;
		bisec[POINTB] = (c_minus_b - b_minus_a) / 2;
		bisec[POINTC] = (a_minus_c - c_minus_b) / 2;

		Normalize(bisec[POINTA]);
		Normalize(bisec[POINTB]);
		Normalize(bisec[POINTC]);

		// Get the cosines of halved angles
		cos_angle2[POINTA] = Dot(b_minus_a, bisec[POINTA]);
		cos_angle2[POINTB] = Dot(c_minus_b, bisec[POINTB]);
		cos_angle2[POINTC] = Dot(a_minus_c, bisec[POINTC]);

		// Case distinction. The angles should be traversed smallest to greatest,
		// because this ensures that more mismatches are excluded very early.
		for (int i = 0; i < 3; ++i)
		{
			if (cos_angle2[i] > cos_angle2[(i + 1) % 3] && cos_angle2[i] > cos_angle2[(i + 2) % 3])
			{
				first = i;

				if (cos_angle2[(i + 1) % 3] > cos_angle2[(i + 2) % 3])
					second = (i + 1) % 3;
				else
					second = (i + 2) % 3;
			}
		}

		origin      = a;
		normal      = Cross(b - a, c - a);

		Normalize(normal);

		// Precomputation bounds
		box.min.SetMin(a);
		box.min.SetMin(b);
		box.min.SetMin(c);

		box.max.SetMax(a);
		box.max.SetMax(b);
		box.max.SetMax(c);

		box.min -= Vec3f(Epsilon);
		box.max += Vec3f(Epsilon);
	};

	virtual bool Intersect(Ray &ray)
	{
		numerator   = Dot(origin - ray.org, normal);
		denominator = Dot(         ray.dir, normal);

		if (denominator != 0.0)
		{
			t = numerator / denominator;

			if (t < ray.t && t > Epsilon)
			{
				p = ray.org + t * ray.dir;

				ray_delta = p - point[first];
				Normalize(ray_delta);

				if (Dot(bisec[first], ray_delta) < cos_angle2[first])
					return false;

				ray_delta = p - point[second];
				Normalize(ray_delta);

				if (Dot(bisec[second], ray_delta) < cos_angle2[second])
					return false;

				ray.t = t;
				ray.hit = this;
				return true;
			}
		}

		return false;
	};

	virtual Vec3f GetNormal(Ray &ray)
	{
		Vec3f normal = this->normal;
		return normal;
	};
};

#endif
