#ifndef MATRIX_H
#define MATRIX_H


#include <cmath>

#include "vector.h"


class Matrix {
public:
	union {
		Real m[3][3];
		struct {
			Real m00, m01, m02;
			Real m10, m11, m12;
			Real m20, m21, m22;
		};
	};

	Matrix(
		Real m00, Real m01, Real m02,
		Real m10, Real m11, Real m12,
		Real m20, Real m21, Real m22
	):
		m00(m00), m01(m01), m02(m02),
		m10(m10), m11(m11), m12(m12),
		m20(m20), m21(m21), m22(m22)
	{ }

	Matrix(): Matrix(
		1, 0, 0,
		0, 1, 0,
		0, 0, 1 ) { }

	Matrix(
		const Vector &axis_x,
		const Vector &axis_y,
		const Vector &offset = Vector()
	): Matrix(
		axis_x.x, axis_x.y, 0,
		axis_y.x, axis_y.y, 0,
		offset.x, offset.y, 1 ) { }

	const Vector& axis(int index) const { return *(const Vector*)(m[index]); }
	const Vector& axis_x() const { return axis(0); }
	const Vector& axis_y() const { return axis(1); }
	const Vector& offset() const { return axis(2); }

	Vector& axis(int index) { return *(Vector*)(m[index]); }
	Vector& axis_x() { return axis(0); }
	Vector& axis_y() { return axis(1); }
	Vector& offset() { return axis(2); }

	Matrix& set_identity()
		{ return *this = Matrix(); }
	bool is_identity() const
		{ return *this == Matrix(); }

	Matrix& set_scale(const Real &x, const Real &y) {
		return *this = Matrix( x, 0, 0,
							   0, y, 0,
							   0, 0, 1 );
	}
	Matrix& set_scale(const Real &s)
		{ return set_scale(s, s); }
	Matrix& set_scale(const Vector &s)
		{ return set_scale(s.x, s.y); }

	Matrix& set_rotate(const Real &angle) {
		Real s = sin(angle), c = cos(angle);
		return *this = Matrix( c, s, 0,
							  -s, c, 0,
							   0, 0, 1 );
	}

	Matrix& set_translate(const Real &x, const Real &y) {
		return *this = Matrix( 1, 0, 0,
							   0, 1, 0,
							   x, y, 1 );
	}
	Matrix& set_translate(const Vector &t)
		{ return set_translate(t.x, t.y); }

	void transform(
		Real &out_x, Real &out_y, Real &out_z,
		const Real &x, const Real &y, const Real &z ) const
	{
		out_x = x*m00 + y*m10 + z*m20;
		out_y = x*m01 + y*m11 + z*m21;
		out_z = x*m02 + y*m12 + z*m22;
	}
	void transform(Real &out_x, Real &out_y, const Real &x, const Real &y, bool translate = true) const
		{ Real z; transform(out_x, out_y, z, x, y, translate ? 1.0 : 0.0); }
	Vector transform(const Vector &v, bool translate = true)const
		{ Vector vv; transform(vv.x, vv.y, v.x, v.y, translate); return vv; }

	bool operator==(const Matrix &other) const {
		return equal(m00, other.m00) && equal(m01, other.m01) && equal(m02, other.m02)
			&& equal(m10, other.m10) && equal(m11, other.m11) && equal(m12, other.m12)
			&& equal(m20, other.m20) && equal(m21, other.m21) && equal(m22, other.m22);
	}
	bool operator!=(const Matrix &other) const
		{ return !(*this == other); }

	Matrix& operator*=(const Matrix &other)
		{ return *this = *this * other; }

	Matrix operator*(const Matrix &other) const {
		Matrix m;
		transform(m.m00, m.m01, m.m02, other.m00, other.m01, other.m02);
		transform(m.m10, m.m11, m.m12, other.m10, other.m11, other.m12);
		transform(m.m20, m.m21, m.m22, other.m20, other.m21, other.m22);
		return m;
	}

	Real det() const {
		return m00*(m11*m22 - m12*m21)
			 - m01*(m10*m22 - m12*m20)
			 + m02*(m10*m21 - m11*m20);
	}
	
	bool is_invertible() const
		{ return !equal(det(), 0); }
	
	Matrix get_inverted() const;

	Matrix& invert()
		{ return *this = get_inverted(); }
};

#endif