#include "geometry_global.h"
#include "vertex.h"
#include <CGAL/bounding_box.h>
#include <CGAL/linear_least_squares_fitting_3.h>

#include <cmath>



namespace Global {

	//////////////////////////////////////////////////////////////////////////
	//-----------------------   global variables

	float	   light_position[4] = { 1.0f, 1.0f, 1.0f, 1.0f };
	float      light_specular = 1.0f;
	float      light_diffuse = 1.0f; 

	// colors for selected and highlighted objects
	Colorf	   color_selected		    = Colorf(1.0f, 0.0f, 0.0f, 1.0f);   
	Colorf	   color_highlighted	    = Colorf(0.0f, 1.0f, 1.0f, 1.0f);    

	// drawing parameters for OpenGL
	//Colorf	   wire_frame_color 	    = Colorf(0.8f, 0.35f, 0.0f, 1.0f); 
	Colorf	   wire_frame_color 	    = Colorf(0.0f, 0.0f, 1.0f, 1.0f); 
	float	   facet_transparency		= 0.7f;
	Colorf	   facet_color				= Colorf(1.0f, 1.0f, 0.0f, facet_transparency); 

	float      line_width			    = 1.0f;
	float      highlighted_line_width	= line_width + 2.0f;

	float      point_size			    = 2.0f;
	float      highlighted_point_size	= 4.0f;

	//////////////////////////////////////////////////////////////////////////

	Bbox2f bbox_of(const std::vector<Point2f>& points) {
		// axis-aligned bounding box of 2D points

		return CGAL::bounding_box(points.begin(), points.end()).bbox(); 
	}

	Bbox3f bbox_of(const std::vector<Point3f>& points) {
		// axis-aligned bounding box of 3D points

		return CGAL::bounding_box(points.begin(), points.end()).bbox(); 
	}

	Bbox3f bbox_of(const std::vector<Vertex*> pts) {
		// axis-aligned bounding box of 3D points		

		unsigned int size = pts.size();
		assert(size > 0);		
		unsigned int mini = 100000;
		unsigned int partial = 1;

		if (size > mini)
			partial = size / 100;	

		std::vector<Point3f> points;	
		for (unsigned int i=0; i<pts.size(); i+=partial) { // I don't need the accurate one
			points.push_back(pts[i]->point()) ;
		}		

		return Global::bbox_of(points);
	}

	float normalize(Vector3f& v) {
		float len = std::sqrt(v.squared_length());
		if (len != 0.0f)
			v = v / len;
		return len;
	} 

	float length_of(const Vector3f& v) {
		return std::sqrt(v.squared_length());
	}

	float normalize(Vector2f& v) {
		float len = std::sqrt(v.squared_length());
		if (len != 0.0f)
			v = v / len;
		return len;
	} 

	float length_of(const Vector2f& v) {
		return std::sqrt(v.squared_length());
	}


	float fit_plane(const std::vector<Point3f>& vertices, Plane3f& plane) {
		unsigned int num = vertices.size();

		if (num < 3) 
			return 0;

		float q = CGAL::linear_least_squares_fitting_3(
			vertices.begin(), 
			vertices.end(), 
			plane,
			CGAL::Dimension_tag<0>()
			);

		return q; 
	}	


	float fit_plane(const std::vector<Vertex*>& pts, Plane3f& plane) {
		unsigned int num = pts.size();

		if (num < 3) 
			return 0;

		std::vector<Point3f>  vertices(num);
		for (unsigned int i=0; i < num; ++i) {
			const Point3f& p = pts[i]->point();
			vertices[i] = Point3f(p.x(), p.y(), p.z());
		}

		return fit_plane(vertices, plane); 
	}	


	Vector3f normal_of(const Plane3f& plane) {
		Vector3f n = plane.orthogonal_vector();
		normalize(n);
		return n;
	} 


	float radian_to_degree(float rad) {
		return (rad * 180.0f / M_PI);
	}

	float degree_to_radian(float deg) {
		return (M_PI * deg / 180.0f); 
	}

	Point3f rotate_in_radian(const Point3f& p, const Point3f& axis_pos, const Vector3f& axis_dir, float angle) {
		//if (Numeric::is_nan(angle))	{
		//	std::cout << "angle is not a number" << std::endl;
		//}

		Line3f axis_line(axis_pos, axis_dir);
		Point3f proj = axis_line.projection(p);

		// compute the translation
		AffTransformation3f mat1(CGAL::Translation(), CGAL::ORIGIN - proj);

		// compute the rotation
		Vector3f u = axis_dir;
		Global::normalize(u);

		float cos_angle = std::cos(angle);
		float sin_angle = std::sin(angle);

		float rot[3][3];
		rot[0][0] = u.x() * u.x() * (1.0f - cos_angle) + cos_angle;
		rot[0][1] = u.x() * u.y() * (1.0f - cos_angle) - u.z() * sin_angle;
		rot[0][2] = u.x() * u.z() * (1.0f - cos_angle) + u.y() * sin_angle;
		rot[1][0] = u.y() * u.x() * (1.0f - cos_angle) + u.z() * sin_angle;
		rot[1][1] = u.y() * u.y() * (1.0f - cos_angle) + cos_angle;
		rot[1][2] = u.y() * u.z() * (1.0f - cos_angle) - u.x() * sin_angle;
		rot[2][0] = u.z() * u.x() * (1.0f - cos_angle) - u.y() * sin_angle;
		rot[2][1] = u.z() * u.y() * (1.0f - cos_angle) + u.x() * sin_angle;
		rot[2][2] = u.z() * u.z() * (1.0f - cos_angle) + cos_angle;

		AffTransformation3f mat2( 
			rot[0][0], rot[0][1], rot[0][2],  
			rot[1][0], rot[1][1], rot[1][2], 
			rot[2][0], rot[2][1], rot[2][2] );

		AffTransformation3f mat3( CGAL::Translation(), proj - CGAL::ORIGIN);

		AffTransformation3f mat = mat3 * mat2 * mat1;

		Point3f new_p = mat.transform(p);
		//if (Numeric::is_nan(p.x()) || Numeric::is_nan(p.y()) || Numeric::is_nan(p.z()))	{
		//	std::cout << "not a number. How can this happen? "
		//		<< "point: " << p << std::endl
		//		<< "rot pos: " << axis_pos << std::endl
		//		<< "rot dir: " << axis_dir << std::endl;
		//}
		return new_p;
	}

	Point3f rotate_in_degree(const Point3f& p, const Point3f& axis_pos, const Vector3f& axis_dir, float angle) {
		float rad = degree_to_radian(angle);  
		return rotate_in_radian(p, axis_pos, axis_dir, rad);
	}


	bool is_point_in_polygon(const Point2f& P, const std::vector<Point2f>& polygon)
	{
		bool bIsInside = false;

		int N = (int)polygon.size();

		int i = 0;
		int j = N-1;

		for (; i<N; j=i, ++i)
		{
			const Point2f& U0 = polygon[i];
			const Point2f& U1 = polygon[j];  // current edge

			if (((U0.y() <= P.y()) && (P.y() < U1.y())) ||  // U1 is above the ray, U0 is on or below the ray
				((U1.y() <= P.y()) && (P.y() < U0.y())))    // U0 is above the ray, U1 is on or below the ray
			{
				// find x-intersection of current edge with the ray. 
				// Only consider edge crossings on the ray to the right of P.
				float x = U0.x() + (P.y() - U0.y()) * (U1.x() - U0.x()) / (U1.y() - U0.y());
				if (x > P.x())
					bIsInside = !bIsInside;
			}
		}

		return bIsInside;
	}

}