#include <CGAL/Homogeneous.h>
#include <CGAL/leda_real.h>
#include <CGAL/Point_3.h>
#include <CGAL/predicates_on_points_3.h>
#include <iostream>

using std::cout;
using std::cin;

typedef CGAL::Homogeneous<leda_real> Rep_class;
typedef CGAL::Point_3<Rep_class> Point_3D;


void error(char *s)
{
	cerr << "error: " << s << endl;
	exit(1);
}

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

// Definition of the class Point

class Point {
	Point_3D p;		// The point in 3D

	public:
		void set();
		void set(leda_real,leda_real,leda_real);
		void set(Point_3D);
		Point_3D get_point();
		void print();
};

void Point::set()
{
	leda_real x,y,z;
	
	cout << "Enter coordinates for the point" << endl;
	cout << ">>" << std::flush;
	cin >> x >> y >> z;

	Point_3D q(x,y,z);
	p=q;
}

void Point::set(leda_real x,leda_real y, leda_real z)
{
	Point_3D q(x,y,z);
	p=q;
}

void Point::set(Point_3D q)
{
	p=q;
}

Point_3D Point::get_point()
{
	return p;
}


void Point::print()
{
	cout << p << endl;
}


// End of definition for class Point 


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


// Definition for class Line

class Line {
	Point p,q;
	leda_real plucker[6];
	void calculate_plucker_coeffs();

	public:
		void set();
		void set(Point, Point);
		Point point(int);
		leda_real& operator[](int);
		void print_coords();
		void print_plucker();
};	

void Line::set()
{
	cout << " Enter two points for the Line " << endl;
	p.set();
	q.set();
	cout << endl;
	calculate_plucker_coeffs();
}

void Line::set(Point x, Point y)
{
	p=x;
	q=y;
	calculate_plucker_coeffs();
	
}

void Line::calculate_plucker_coeffs()
{
	Point_3D p1,q1;
	p1 = p.get_point();
	q1 = q.get_point();

	plucker[0]= p1.hx()*q1.hy()-q1.hx()*p1.hy();
	plucker[1]= p1.hx()*q1.hz()-q1.hx()*p1.hz();
	plucker[2]= p1.hx()-q1.hx();
	plucker[3]= p1.hy()*q1.hz()-q1.hy()*p1.hz();
	plucker[4]= p1.hz()-q1.hz();
	plucker[5]= q1.hy()-p1.hy();
}

Point Line::point(int i)
{
	if(i<0 || i>1) error(" Invalid parameter for Line ");	
	if(i==0) return p;
	else  return q;
}

leda_real& Line::operator[](int i)
{
	if(i<0 || i>5) error(" Plucker index out of range"); 
	return plucker[i];
}

void Line::print_coords()
{
	p.print();
	q.print();
	cout << endl;	
}

void Line::print_plucker()
{
	int i;
	for(i=0;i<=5;i++) cout << plucker[i] << "  ";
	cout << endl;
}

// End of definition for class Line


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

// Definition of class triangle

class Triangle {
	Point a,b,c;
	Line l1,l2,l3;

	public:
		void set();
		void set(Point,Point,Point);
		Point vertice(int);
		Line edge(int);
		Point vertice(int,int);
		void print();
};

void Triangle::set()
{
	cout << " Enter three points for the triangle " << endl;
	a.set();
	b.set();
	c.set();
	cout << endl << endl;	
		
	l1.set(a,b);
	l2.set(b,c);
	l3.set(c,a);
}

void Triangle::set(Point x, Point y, Point z)
{
	a=x;b=y;c=z;
	l1.set(a,b);
	l2.set(b,c);
	l3.set(c,a);
}
	
Point Triangle::vertice(int i)
{
	if(i<0 || i>2) error(" Parameter out of range for triangle ");
	if(i==0) return a;
	else if(i==1) return b;
	else return c;
}

Line Triangle::edge(int i)
{
	if(i<0 || i>2) error(" Parameter out of range for triangle ");
	if(i==0) return l1;
	else if(i==1) return l2;
	else return l3;
}

Point Triangle::vertice(int m, int n) //Returns the vertex between edge m and n
{
	if(m==0 && n==1) return b;
	else if(m==0 && n==2) return a;
	else return c;
}

void Triangle::print()
{
	a.print();
	b.print();
	c.print();
	cout << endl;
}

// End of definition for class triangle

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


// Computes the side operator of two lines 

leda_real side_product(Line l1, Line l2)
{
	leda_real result = 0;
	
	result = l1[0]*l2[4]+l1[1]*l2[5]+l1[2]*l2[3]+l1[3]*l2[2]+l1[4]*l2[0] \
			+ l1[5]*l2[1];
	return result;
}

// Computes the 2D case when the triangle and the line lie in the same plane 

void compute_2D_intersection(Triangle T, Line L)
{
	Point p[4];	// Four trial points
	Point x;	// The selected point : This point does not lie in the same plane
	Triangle T1;	// Trial triangle constructed to test coplanarity
	int edge[] = {0,0,0},
	    vertex[] ={0,0,0}; 	// Markers for intersection

	int i,intersection=0;
	leda_real s1,s2,s3;
	
	p[0].set(0,0,0);
	p[1].set(0,1,0);
	p[2].set(0,0,1);
	p[3].set(1,0,0);
	
	for(i=0;i<4;i++)
	{
		T1.set(T.vertice(0),T.vertice(1),p[i]);
		s1 = side_product(L,T1.edge(0));
		s2 = side_product(L,T1.edge(1));
		s3 = side_product(L,T1.edge(2));
	
		if(!(s1==0 && s2==0 && s3==0)) break;
	}
	
	x = p[i]; // This plane is not in plane containing the Triangle T  and the Line L
	
	for(i=0;i<3;i++)
	{
		T1.set(T.edge(i).point(0),T.edge(i).point(1),x);
		s2 = side_product(L,T1.edge(1));
		s3 = side_product(L,T1.edge(2));
	
		if(s2*s3<0) ;
		else if(s2*s3>0) 
			 { edge[i]++;
			   intersection++; }
		else if(s2==0 && s3==0) 
			 {  edge[i] +=2;
			    intersection += 2;
			    break; }
		else if(s2==0)
			{  vertex[(i+1) -3*((i+1)/3)]++;
			   intersection++; }
		else if(s3==0)
			{  vertex[i]++;
			   intersection++; }
	} 

	// Printing the result 
	
	if(intersection==0) cout << " But the Line and Triangle do not intersect \n";
	else {
		for(i=0;i<3;i++) {
			if (edge[i]==1) { cout << " The Line intersects the edge \n ";
					  T.edge(i).print_coords(); }
			else if(edge[i]==2) { cout << "The Line touches the edge \n ";
					      T.edge(i).print_coords(); break; }
		}
		
		for(i=0;i<3 && intersection <=3;i++) {
			if(vertex[i]!=0) { cout << " The Line passes through the vertex ";
					   T.vertice(i).print(); }
		}
	}
}


// Computes the intersection between a triangle and line

void compute_intersection( Triangle T, Line L)
{
	leda_real s1,s2,s3;
	
	s1 = side_product(L,T.edge(0));
	s2 = side_product(L,T.edge(1));
	s3 = side_product(L,T.edge(2));

	if(s1==0 && s2==0 && s3==0)  
		 {  cout << " Line and Triangle are coplanar \n\n"; 
	            compute_2D_intersection(T,L); }
	else if(( s1>0 && s2>0 && s3>0)||(s1<0 && s2<0 && s3<0)) 
		cout << " Line passes through the triangle \n";
	else if(s1==0 && s2*s3>0)
		 { cout << " Line passes through the edge  \n" ;
		 T.edge(0).print_coords(); }
	else if(s2==0 && s1*s3>0)
		 { cout << " Line passes through the edge \n" ;
		 T.edge(1).print_coords(); }
	else if(s3==0 && s1*s2>0)
		 { cout << " Line passes through the edge \n" ;
		 T.edge(2).print_coords(); }
	else if(s1==0 && s2==0)
		 { cout << " Line passes through the vertex " ;
		 T.vertice(1).print(); }
	else if(s1==0 && s3==0)
		 { cout << " Line passes through the vertex " ; 
		 T.vertice(0).print(); }
	else if(s2==0 && s3==0)
		 { cout << " Line passes through the vertex " ;
		 T.vertice(2).print(); }
	else cout << " The Line and triangle do not intersect "<< endl;

}


int main()
{
	Triangle T;
	Line L;

	cout << endl;	
	T.set();
	L.set();

	compute_intersection(T,L);	

	return 0;
}