#pragma once
#include "../../stdafx.h"
#include "glLibMathMatrix.h"

namespace glLibMathVec {
	void glLibInternal_math_init(void);
	void glLibInternal_math_deinit(void);


	template <typename type_float> type_float length(type_float vec[], int vec_length) {
		type_float squared_distance = 0.0;
		for (int i=0;i<vec_length;++i) { squared_distance += vec[i]*vec[i]; }
		return pow(squared_distance,(type_float)0.5);
	}


	template <typename type_float> type_float dot(type_float vec1[], type_float vec2[], int vec_length) {
		type_float result = 0.0;
		for (int i=0;i<vec_length;++i) {
			result += vec1[i]*vec2[i];
		}
		return result;
	}
	template <typename type_float> void cross(type_float vec1[3], type_float vec2[3], type_float result[3]) {
		result[0] = vec1[1]*vec2[2]-vec1[2]*vec2[1];
		result[1] = vec1[2]*vec2[0]-vec1[0]*vec2[2];
		result[2] = vec1[0]*vec2[1]-vec1[1]*vec2[0];
	}


	template <typename type_float> type_float min(type_float vec[], int vec_length) {
		type_float minimum = vec[0];
		int index = 1;
		while (index<vec_length) {
			minimum = std::min(minimum,vec[index++]);
		}
		return minimum;
	}
	template <typename type_float> type_float max(type_float vec[], int vec_length) {
		type_float maximum = vec[0];
		int index = 1;
		while (index<vec_length) {
			maximum = std::max(maximum,vec[index++]);
		}
		return maximum;
	}

	
	void rotate_rad(float  point[3], float  result[3], float  vec[3],         float  theta);
	void rotate_rad(double point[3], double result[3], double vec[3],         double theta);
	void rotate_deg(float  point[3], float  result[3], float  vec[3],         float  theta);
	void rotate_deg(double point[3], double result[3], double vec[3],         double theta);
	void rotate(float  point[3], float  result[3], glLibMathMatrix< float,char>* rot_matrix);
	void rotate(double point[3], double result[3], glLibMathMatrix<double,char>* rot_matrix);


	template <typename type_float> bool normalize(type_float vec[], int vec_length) {
		type_float l = length(vec,vec_length);
		if (l!=0.0) {
			for (int i=0;i<vec_length;++i) {
				vec[i] = vec[i]/l;
			}
			return true;
		}
		for (int i=0;i<vec_length;++i) { vec[i]=0.0; }
		return false;
	}
	template <typename type_float> bool normalize(type_float vec[], type_float result[], int vec_length) {
		type_float l = length(vec,vec_length);
		if (l!=0.0) {
			for (int i=0;i<vec_length;++i) {
				result[i] = vec[i]/l;
			}
			return true;
		}
		for (int i=0;i<vec_length;++i) { result[i]=0.0; }
		return false;
	}


	template <typename type_num> void scale(type_num sc, type_num vec[], int vec_length) {
		for (int i=0;i<vec_length;++i) { vec[i]=sc*vec[i]; }
	}
	template <typename type_num> void scale(type_num sc, type_num vec[], type_num result[], int vec_length) {
		for (int i=0;i<vec_length;++i) { result[i]=sc*vec[i]; }
	}
	template <typename type_num> void negate(type_num vec[], int vec_length) {
		for (int i=0;i<vec_length;++i) { vec[i]=-vec[i]; }
	}
	template <typename type_num> void negate(type_num vec[], type_num result[], int vec_length) {
		for (int i=0;i<vec_length;++i) { result[i]=-vec[i]; }
	}
	template <typename type_num> void add(type_num vec1[], type_num vec2[], type_num result[], int vec_length) {
		for (int i=0;i<vec_length;++i) { result[i]=vec1[i]+vec2[i]; }
	}
	template <typename type_num> void subtract(type_num vec1[], type_num vec2[], type_num result[], int vec_length) {
		for (int i=0;i<vec_length;++i) { result[i]=vec1[i]-vec2[i]; }
	}
	template <typename type_num1, typename type_num2, typename type_num3> void multiply(type_num1 vec1[], type_num2 vec2[], type_num3 result[], int vec_length) {
		for (int i=0;i<vec_length;++i) { result[i]=vec1[i]*vec2[i]; }
	}
	template <typename type_num1, typename type_num2, typename type_num3> void divide(type_num1 vec1[], type_num2 vec2[], type_num3 result[], int vec_length) {
		for (int i=0;i<vec_length;++i) { result[i]=vec1[i]/vec2[i]; }
	}


	template <typename type_num> void zero(type_num vec[], int vec_length) {
		for (int i=0;i<vec_length;++i) { vec[i]=0; }
	}


	template <typename type_num1, typename type_num2> void copy(type_num1 src[], type_num2 dst[], int vec_length) {
		for (int i=0;i<vec_length;++i) { dst[i] = src[i]; }
	}


	template <typename type_num1, typename type_num2> bool equals(type_num1 vec1[], type_num2 vec2[], int vec_length) {
		for (int i=0;i<vec_length;++i) { if (vec1[i]!=vec2[i]) return false; }
		return true;
	}


	template <typename type_float> void print(type_float vec[], int vec_length, bool add_newline=true) {
		printf("[");
		for (int i=0;i<vec_length-1;++i) {
			printf("%f, ",vec[i]);
		}
		printf("%f",vec[vec_length-1]);
		(add_newline) ? printf("]\n") : printf("]");
	}


	//TODO reflect
	template <typename type_float> void refract(type_float* incident, type_float* normal, type_float* result, type_float eta, int vec_length) {
		//printf("incident = <%f,%f,%f>\n",incident[0],incident[1],incident[2]);
		//printf("normal = <%f,%f,%f>\n",normal[0],normal[1],normal[2]);
		//printf("eta = %f\n",eta);
	
		type_float n_dot_i = dot(normal,incident,vec_length);
		//printf("n * i = %f\n",n_dot_i);
		type_float k = 1.0 - eta * eta * (1.0 - n_dot_i * n_dot_i);
		//printf("k = %f\n",k);
		if (k<0.0) return;
		type_float* eta_i = new type_float[vec_length];
		sc_vec(eta,incident,eta_i,vec_length);
	
		//printf("eta_i = <%f,%f,%f>\n",eta_i[0],eta_i[1],eta_i[2]);
		sc_vec( eta*n_dot_i+sqrt(k), normal, result, vec_length);
		//printf("result = <%f,%f,%f>\n",result[0],result[1],result[2]);
		vec_subt(eta_i,result,result,vec_length);
		//printf("result = <%f,%f,%f>\n\n",result[0],result[1],result[2]);
		delete eta_i;
	}
};