#ifndef __gxframework_matrix4_h__
	#define __gxframework_matrix4_h__

#include "shared.h"
#include "gxmath.h"

#include "vector4.h"

namespace Gx
{
	class Matrix4;

	static float det( const Gx::Matrix4& m );
	static void inverse( const Gx::Matrix4& m, Gx::Matrix4& inv );

	class Matrix4
	{
	public:
		Matrix4() { }

		Matrix4(const Vector4& ncol0, const Vector4& ncol1, const Vector4& ncol2, const Vector4& ncol3)
			: column0(ncol0), column1(ncol1), column2(ncol2), column3(ncol3) { }

		Matrix4(const Vector3& ncol0, const Vector3& ncol1, const Vector3& ncol2, const Vector3& ncol3)
			: column0(ncol0, 0.0f), column1(ncol1, 0.0f), column2(ncol2, 0.0f), column3(ncol3, 1.0f) { }

		Matrix4(float values[])
			: column0(values[0],values[1],values[2], values[3]), column1(values[4],values[5],values[6], values[7]),
			column2(values[8],values[9],values[10], values[11]), column3(values[12], values[13], values[14], values[15]) { }

		Matrix4(const Vector4& diagonal)
			: column0(diagonal.x,0.0f,0.0f,0.0f), column1(0.0f,diagonal.y,0.0f,0.0f),
			column2(0.0f,0.0f,diagonal.z,0.0f), column3(0.0f,0.0f,0.0f,diagonal.w) { }

		Matrix4(const Matrix4& m)
			: column0(m.column0), column1(m.column1), column2(m.column2), column3(m.column3) { }

		Matrix4& operator=(const Matrix4& m) {
			column0 = m.column0;
			column1 = m.column1;
			column2 = m.column2;
			column3 = m.column3;
			return *this;
		}

		static const Matrix4 Zero;
		static const Matrix4 Identity;

		static Matrix4 createDiagonal(const Vector4& d) {
			return Matrix4(
				Vector4(d.x, 0.0f, 0.0f, 0.0f),
				Vector4(0.0f, d.y, 0.0f, 0.0f),
				Vector4(0.0f, 0.0f, d.z, 0.0f),
				Vector4(0.0f, 0.0f, 0.0f, d.w));
		}

		static Matrix4 createTranslation(float x, float y, float z) {
			Matrix4 translation = Matrix4::Identity;
			translation.column3 = Vector4(x, y, z, 1.0f);
			return translation;
		}

		static Matrix4 createTranslation(const Vector3& v) {
			return createTranslation(v.x, v.y, v.z);
		}

		static Matrix4 createScale(float x, float y, float z) {
			return createDiagonal(Vector4(x, y, z, 1.0f));
		}

		static Matrix4 createScale(float s) {
			return createScale(s, s, s);
			// return createDiagonal(Vector4(1.0f, 1.0f, 1.0f, 1.0f / s));
		}

		static Matrix4 createScale(const Vector3& v) {
			return createScale(v.x, v.y, v.z);
		}

		static Matrix4 createRotationX(float angleRadians) {
			const float c = cos(angleRadians);
			const float s = sin(angleRadians);

			Matrix4 rot = createDiagonal(Vector4(1.0f, c, c, 1.0f));
			rot.column1.z = s; rot.column2.y = -s;

			return rot;
		}

		static Matrix4 createRotationY(float angleRadians) {
			const float c = cos(angleRadians);
			const float s = sin(angleRadians);

			Matrix4 rot = createDiagonal(Vector4(c, 1.0f, c, 1.0f));
			rot.column0.z = -s; rot.column2.x = s;

			return rot;
		}

		static Matrix4 createRotationZ(float angleRadians) {
			const float c = cos(angleRadians);
			const float s = sin(angleRadians);

			Matrix4 rot = createDiagonal(Vector4(c, c, 1.0f, 1.0f));
			rot.column0.y = s; rot.column1.x = -s;

			return rot;
		}

		// view matrix
		static Matrix4 createLookAt(Vector3 cameraPosition, Vector3 cameraTarget, Vector3 cameraUpVector) {
			//Matrix4 proj = Matrix4::Identity;

			Vector3 vectorZ = (cameraPosition - cameraTarget);
			vectorZ.normalize();

			Vector3 vectorX = cameraUpVector.cross(vectorZ);
			vectorX.normalize();

			Vector3 vectorY = vectorZ.cross(vectorX);

			Matrix4 proj = Matrix4(
				vectorX,
				vectorY,
				vectorZ,
				Vector3(
					-vectorX.dot(cameraPosition),
					-vectorY.dot(cameraPosition),
					-vectorZ.dot(cameraPosition)));

			return proj;
		}

		// projection matrix
		static Matrix4 createOrthographic(float xMin, float xMax, float yMin, float yMax, float zMin, float zMax) {
			// Optimize it!
			Matrix4 proj = Matrix4::Identity;

			// TIPS: createDiagonal?
			proj.column0.x = 2.0f / (xMax - xMin);
			proj.column1.y = 2.0f / (yMax - yMin);
			proj.column2.y = -2.0f / (zMax - zMin);

			proj.column3 = Vector4(
				-((xMax + xMin)/(xMax - xMin)),
				-((yMax + yMin)/(yMax - yMin)),
				-((zMax + zMin)/(zMax - zMin)),
				1.0f
				);
		}

		static Matrix4 createPerspective(float fov, float aspect, float zMin, float zMax) {
			Matrix4 proj = Matrix4::Identity;

			const float yMax = zMin * tanf(fov * 0.5f);
			const float yMin = -yMax;
			const float xMin = yMin * aspect;
			const float xMax = -xMin;

			proj.column0.x = (2.0f * zMin) / (xMax - xMin);
			proj.column1.y = (2.0f * zMin) / (yMax - yMin);

			proj.column2 = Vector4(
				(xMax + xMin) / (xMax - xMin),
				(yMax + yMin) / (yMax - yMin),
				-((zMax + zMin) / (zMax - zMin)),
				-1.0f
				);

			proj.column3.z = -((2.0f * (zMax*zMin))/(zMax - zMin));
			proj.column3.w = 0.0f;

			return proj;
		}

		static Matrix4 createOrthographic(float width, float heigth, float near, float far) {
			// Not implemented yet!
			return Matrix4::Identity;
		}

		Matrix4 getTranspose() const {
			return Matrix4(
				Vector4(column0.x, column1.x, column2.x, column3.x),
				Vector4(column0.y, column1.y, column2.y, column3.y),
				Vector4(column0.z, column1.z, column2.z, column3.z),
				Vector4(column0.w, column1.w, column2.w, column3.w));
		}


		Matrix4 getInverse() const {
			Matrix4 inv;
			Gx::inverse(*this, inv);
			return inv;
		}

		float getDeterminant() const {
			return det(*this);
		}


		float getMinor(unsigned int row, unsigned int col) const {

		}

		void transpose() {
			// Rewrite it!
			*this = getTranspose();
		}

		void inverse() {
			// Rewrite it!
			*this = getInverse();
		}

		Matrix4 operator-() const {
			return Matrix4(-column0, -column1, -column2, -column3);
		}

		Matrix4& operator+=(const Matrix4& m) {
			column0 += m.column0;
			column1 += m.column1;
			column2 += m.column2;
			column3 += m.column3;
			return *this;
		}

		Matrix4& operator-=(const Matrix4& m) {
			column0 -= m.column0;
			column1 -= m.column1;
			column2 -= m.column2;
			column3 -= m.column3;
			return *this;
		}

		/*return Matrix4(
			m1.transform(m2.column0),
			m1.transform(m2.column1),
			m1.transform(m2.column2),
			m1.transform(m2.column3));*/

		Matrix4& operator*=(const Matrix4& m) {
			// Not implemented yet!
			return *this;
		}

		Matrix4& operator*=(const float scalar) {
			column0 *= scalar;
			column1 *= scalar;
			column2 *= scalar;
			column3 *= scalar;
			return *this;
		}

		// Element access, mathematical way!
		float operator()(unsigned int row, unsigned int col) const {
			return (*this)[col][row];
		}

		// Element access, mathematical way!
		float& operator()(unsigned int row, unsigned int col) {
			return (*this)[col][row];
		}

		// Transform vector by matrix, equal to v' = M*v
		Vector4 transform(const Vector4& v) const {
			return column0 * v.x + column1 * v.y + column2 * v.z + column3 * v.w;
		}

		// Transform vector by matrix, equal to v' = M*v
		Vector3 transform(const Vector3& v) const {
			return transform(Vector4(v, 1.0f)).getXyz();
		}

		// Rotate vector by matrix, equal to v' = M*v
		Vector4 rotate(const Vector4& v) const {
			return column0 * v.x + column1 * v.y + column2 * v.z; //+ column3 * 0;
		}

		// Rotate vector by matrix, equal to v' = M*v
		Vector3 rotate(const Vector3& v) const {
			return rotate(Vector4(v, 1.0f)).getXyz();
		}

		Vector3 getBasis(int num) const {
			ASSERT(num >= 0 && num <= 3);
			return (&column0)[num].getXyz();
		}

		Vector3 getPosition() const {
			return column3.getXyz();
		}

		void setPosition(const Vector3& position) {
			column3.x = position.x;
			column3.y = position.y;
			column3.z = position.z;
		}

		const float* begin() const {
			return &column0.x;
		}

		Vector4& operator[](unsigned int col) {
			return (&column0)[col];
		}

		const Vector4& operator[](unsigned int col) const {
			return (&column0)[col];
		}

		Vector4 column0, column1, column2, column3;
	};

	__declspec(selectany) const Matrix4 Matrix4::Zero = Matrix4(Vector4(0.0f), Vector4(0.0f), Vector4(0.0f), Vector4(0.0f));

	__declspec(selectany) const Matrix4 Matrix4::Identity = Matrix4(
		Vector4(1.0f, 0.0f, 0.0f, 0.0f),
		Vector4(0.0f, 1.0f, 0.0f, 0.0f),
		Vector4(0.0f, 0.0f, 1.0f, 0.0f),
		Vector4(0.0f, 0.0f, 0.0f, 1.0f));

	static Matrix4 operator+(const Matrix4& m1, const Matrix4& m2)
	{
		return Matrix4(
			m1.column0 + m2.column0,
			m1.column1 + m2.column1,
			m1.column2 + m2.column2,
			m1.column3 + m2.column3);
	}

	static Matrix4 operator-(const Matrix4& m1, const Matrix4& m2)
	{
		return Matrix4(
			m1.column0 - m2.column0,
			m1.column1 - m2.column1,
			m1.column2 - m2.column2,
			m1.column3 - m2.column3);
	}

	static Matrix4 operator*(const Matrix4& m, const float scalar)
	{
		return Matrix4(m.column0 * scalar, m.column1 * scalar, m.column2 * scalar, m.column3 * scalar);
	}

	static Matrix4 operator*(const Matrix4& m1, const Matrix4& m2)
	{
		// Rows from this <dot> columns from other
		// column0 = transform(other.column0) etc
		return Matrix4(
			m1.transform(m2.column0),
			m1.transform(m2.column1),
			m1.transform(m2.column2),
			m1.transform(m2.column3));
	}

	// www.euclideanspace.com/maths/algebra/matrix/functions/determinant/fourD/index.htm
	static float det( const Gx::Matrix4& m )
	{
		return m.column0.w * m.column1.z * m.column2.y * m.column3.x - m.column0.z * m.column1.w * m.column2.y * m.column3.x
			- m.column0.w * m.column1.y * m.column2.z * m.column3.x + m.column0.y * m.column1.w * m.column2.z * m.column3.x
			+ m.column0.z * m.column1.y * m.column2.w * m.column3.x - m.column0.y * m.column1.z * m.column2.w * m.column3.x
			- m.column0.w * m.column1.z * m.column2.x * m.column3.y + m.column0.z * m.column1.w * m.column2.x * m.column3.y
			+ m.column0.w * m.column1.x * m.column2.z * m.column3.y - m.column0.x * m.column1.w * m.column2.z * m.column3.y
			- m.column0.z * m.column1.x * m.column2.w * m.column3.y + m.column0.x * m.column1.z * m.column2.w * m.column3.y
			+ m.column0.w * m.column1.y * m.column2.x * m.column3.z - m.column0.y * m.column1.w * m.column2.x * m.column3.z
			- m.column0.w * m.column1.x * m.column2.y * m.column3.z + m.column0.x * m.column1.w * m.column2.y * m.column3.z
			+ m.column0.y * m.column1.x * m.column2.w * m.column3.z - m.column0.x * m.column1.y * m.column2.w * m.column3.z
			- m.column0.z * m.column1.y * m.column2.x * m.column3.w + m.column0.y * m.column1.z * m.column2.x * m.column3.w
			+ m.column0.z * m.column1.x * m.column2.y * m.column3.w - m.column0.x * m.column1.z * m.column2.y * m.column3.w
			- m.column0.y * m.column1.x * m.column2.z * m.column3.w + m.column0.x * m.column1.y * m.column2.z * m.column3.w;
	}

	// www.euclideanspace.com/maths/algebra/matrix/functions/inverse/fourD/index.htm
	static void inverse( const Gx::Matrix4& m, Gx::Matrix4& inv )
	{
		const float det = m.getDeterminant();
		inv = Gx::Matrix4::Identity;

		if ( det != 0 )
		{
			const float invDet = 1.0f / det;

			inv.column0.x = invDet * ( m.column1.z * m.column2.w * m.column3.y - m.column1.w * m.column2.z * m.column3.y + m.column1.w * m.column2.y * m.column3.z - m.column1.y * m.column2.w * m.column3.z - m.column1.z * m.column2.y * m.column3.w + m.column1.y * m.column2.z * m.column3.w );
			inv.column0.y = invDet * ( m.column0.w * m.column2.z * m.column3.y - m.column0.z * m.column2.w * m.column3.y - m.column0.w * m.column2.y * m.column3.z + m.column0.y * m.column2.w * m.column3.z + m.column0.z * m.column2.y * m.column3.w - m.column0.y * m.column2.z * m.column3.w );
			inv.column0.z = invDet * ( m.column0.z * m.column1.w * m.column3.y - m.column0.w * m.column1.z * m.column3.y + m.column0.w * m.column1.y * m.column3.z - m.column0.y * m.column1.w * m.column3.z - m.column0.z * m.column1.y * m.column3.w + m.column0.y * m.column1.z * m.column3.w );
			inv.column0.w = invDet * ( m.column0.w * m.column1.z * m.column2.y - m.column0.z * m.column1.w * m.column2.y - m.column0.w * m.column1.y * m.column2.z + m.column0.y * m.column1.w * m.column2.z + m.column0.z * m.column1.y * m.column2.w - m.column0.y * m.column1.z * m.column2.w );
			inv.column1.x = invDet * ( m.column1.w * m.column2.z * m.column3.x - m.column1.z * m.column2.w * m.column3.x - m.column1.w * m.column2.x * m.column3.z + m.column1.x * m.column2.w * m.column3.z + m.column1.z * m.column2.x * m.column3.w - m.column1.x * m.column2.z * m.column3.w );
			inv.column1.y = invDet * ( m.column0.z * m.column2.w * m.column3.x - m.column0.w * m.column2.z * m.column3.x + m.column0.w * m.column2.x * m.column3.z - m.column0.x * m.column2.w * m.column3.z - m.column0.z * m.column2.x * m.column3.w + m.column0.x * m.column2.z * m.column3.w );
			inv.column1.z = invDet * ( m.column0.w * m.column1.z * m.column3.x - m.column0.z * m.column1.w * m.column3.x - m.column0.w * m.column1.x * m.column3.z + m.column0.x * m.column1.w * m.column3.z + m.column0.z * m.column1.x * m.column3.w - m.column0.x * m.column1.z * m.column3.w );
			inv.column1.w = invDet * ( m.column0.z * m.column1.w * m.column2.x - m.column0.w * m.column1.z * m.column2.x + m.column0.w * m.column1.x * m.column2.z - m.column0.x * m.column1.w * m.column2.z - m.column0.z * m.column1.x * m.column2.w + m.column0.x * m.column1.z * m.column2.w );
			inv.column2.x = invDet * ( m.column1.y * m.column2.w * m.column3.x - m.column1.w * m.column2.y * m.column3.x + m.column1.w * m.column2.x * m.column3.y - m.column1.x * m.column2.w * m.column3.y - m.column1.y * m.column2.x * m.column3.w + m.column1.x * m.column2.y * m.column3.w );
			inv.column2.y = invDet * ( m.column0.w * m.column2.y * m.column3.x - m.column0.y * m.column2.w * m.column3.x - m.column0.w * m.column2.x * m.column3.y + m.column0.x * m.column2.w * m.column3.y + m.column0.y * m.column2.x * m.column3.w - m.column0.x * m.column2.y * m.column3.w );
			inv.column2.z = invDet * ( m.column0.y * m.column1.w * m.column3.x - m.column0.w * m.column1.y * m.column3.x + m.column0.w * m.column1.x * m.column3.y - m.column0.x * m.column1.w * m.column3.y - m.column0.y * m.column1.x * m.column3.w + m.column0.x * m.column1.y * m.column3.w );
			inv.column2.w = invDet * ( m.column0.w * m.column1.y * m.column2.x - m.column0.y * m.column1.w * m.column2.x - m.column0.w * m.column1.x * m.column2.y + m.column0.x * m.column1.w * m.column2.y + m.column0.y * m.column1.x * m.column2.w - m.column0.x * m.column1.y * m.column2.w );
			inv.column3.x = invDet * ( m.column1.z * m.column2.y * m.column3.x - m.column1.y * m.column2.z * m.column3.x - m.column1.z * m.column2.x * m.column3.y + m.column1.x * m.column2.z * m.column3.y + m.column1.y * m.column2.x * m.column3.z - m.column1.x * m.column2.y * m.column3.z );
			inv.column3.y = invDet * ( m.column0.y * m.column2.z * m.column3.x - m.column0.z * m.column2.y * m.column3.x + m.column0.z * m.column2.x * m.column3.y - m.column0.x * m.column2.z * m.column3.y - m.column0.y * m.column2.x * m.column3.z + m.column0.x * m.column2.y * m.column3.z );
			inv.column3.z = invDet * ( m.column0.z * m.column1.y * m.column3.x - m.column0.y * m.column1.z * m.column3.x - m.column0.z * m.column1.x * m.column3.y + m.column0.x * m.column1.z * m.column3.y + m.column0.y * m.column1.x * m.column3.z - m.column0.x * m.column1.y * m.column3.z );
			inv.column3.w = invDet * ( m.column0.y * m.column1.z * m.column2.x - m.column0.z * m.column1.y * m.column2.x + m.column0.z * m.column1.x * m.column2.y - m.column0.x * m.column1.z * m.column2.y - m.column0.y * m.column1.x * m.column2.z + m.column0.x * m.column1.y * m.column2.z );
		}
	}

	static std::ostream& operator<<(std::ostream& stream, const Matrix4& m)
	{
		// TIPS: można transponować i skorzystać z << dla wektorów

		stream << std::endl <<
			"{" << m.column0.x << " " << m.column1.x << " " << m.column2.x << " " << m.column3.x <<"}" << std::endl <<
			"{" << m.column0.y << " " << m.column1.y << " " << m.column2.y << " " << m.column3.y <<"}" << std::endl <<
			"{" << m.column0.z << " " << m.column1.z << " " << m.column2.z << " " << m.column3.z <<"}" << std::endl <<
			"{" << m.column0.w << " " << m.column1.w << " " << m.column2.w << " " << m.column3.w <<"}";

		return stream;
	}
}

#endif /* ~__gxframework_matrix4_h__ */