// Implementation file for the Matrix class
// Created by:		Michael Edwards
// Date Started:	6/07/04
// Date Finished:
// Compiler Used:	Microsoft Visual ++, Version 6.0 Enterprise Edition

#include "Matrix.h"

Matrix::Matrix(): column(0), data(0), row(0)		{}

Matrix::Matrix(unsigned userRow, unsigned userColumn): column(userColumn), row(userRow)
// Creates a matrix of row and column with all elements = 0
// PRECONDITION:  none.
// POSTCONDITION: a matrix is created with size (row, column) 
{
	data = new float*[row];
	for (unsigned index = 0; index < row; ++index)
	{
		data[index] = new float[column];
		memset(data[index], 0, sizeof(float) * column);
	}
}

Matrix::Matrix(unsigned userRow, unsigned userColumn, float** userData): column(userColumn), row(userRow)
// Converts a float matrix to a matrix object
// PRECONDITION:  for all (row, column), userData[row][column] is valid
// POSTCONDITION: a matrix is created with size (row, column) and filled
//				  with userData
{
	data = CreateMatrix(userRow, userColumn);
	for (unsigned r = 0; r < userRow; ++r)
		memcpy(data[r], &userData[r][0], sizeof(float) * column);
}

Matrix::Matrix(unsigned userRow, unsigned userColumn, float* userData): column(userColumn), row(userRow)
// Fills a Matrix with an array of data, userData
// PRECONDITION:  length of userData >= row * column
// POSTCONDITION: a matrix is created with size (row, column) and filled
//				  with userData
{
	data = CreateMatrix(userRow, userColumn);
	unsigned index = 0;
	for (unsigned r = 0; r < userRow; ++r, index += column)
		memcpy(data[r], &userData[index], sizeof(float) * column);
}

Matrix::Matrix(const Matrix& rhs): row(rhs.row), column(rhs.column)
// Copy constructor
{
	data = CreateMatrix(row, column);
	for (unsigned r = 0; r < row; ++r)
		memcpy(data[r], &rhs.data[r][0], sizeof(float) * column);
}

Matrix::~Matrix()
// Default destructor
{
	DestroyMatrix(row, data);
}

Matrix Matrix::operator + (const Matrix &rhs) const
// returns the matrix addition of this and rhs
// PRECONDITION:  matrices must be of same dimension.
// POSTCONDITION: returns this + rhs
{
	assert( (row == rhs.row) && (column == rhs.column) );
	Matrix temp(row, column);
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			temp.data[r][c] = data[r][c] + rhs.data[r][c];
	return temp;
}

Matrix Matrix::operator - (const Matrix &rhs) const
// returns the matrix substraction of this and rhs
// PRECONDITION:  matrices must be of same dimension.
// POSTCONDITION: returns this - rhs
{
	assert( (row == rhs.row) && (column == rhs.column) );
	Matrix temp(row, column);
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			temp.data[r][c] = data[r][c] - rhs.data[r][c];
	return temp;
}

Matrix Matrix::operator * (const Matrix &rhs) const
// returns the cross product of this and rhs
// NOTE: this throws an exception if this->column != rhs.row
// PRECONDITION:  this->column == rhs.row
// POSTCONDITION: returns the cross product of this and rhs
{
	assert(column == rhs.row);
	Matrix temp(row, column);
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			for (unsigned i = 0; i < column; ++i)
				temp.data[r][c] += data[r][i] * rhs.data[i][c];
	return temp;
}

Matrix Matrix::operator * (const float scalar) const
// does the scalar multiplication of this and a scalar, scalar
// PRECONDITION:  none.
// POSTCONDITION: returns the scalar multiplication of this and rhs
{
	Matrix temp(row, column);
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			temp.data[r][c] = data[r][c] * scalar;
	return temp;
}

Matrix operator * (const float scalar, const Matrix& rhs)
{
	Matrix temp(rhs.row, rhs.column);
	for (unsigned r = 0; r < rhs.row; ++r)
		for (unsigned c = 0; c < rhs.column; ++c)
			temp.data[r][c] = rhs.data[r][c] * scalar;
	return temp;
}

Matrix Matrix::operator / (const float scalar) const
// does the scalar division of this and a scalar, scalar
// PRECONDITION:  scalar != 0.
// POSTCONDITION: returns the scalar division of this and rhs
{
	assert(scalar != 0);
	Matrix temp(row, column);
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			temp.data[r][c] = data[r][c] / scalar;
	return temp;
}

Matrix& Matrix::operator = (const Matrix &rhs)
// assignment operator for the matrix class
// PRECONDITION:  none.
// POSTCONDITION: this = rhs
{
	if ((row != rhs.row) || (column != rhs.column))
	{
		DestroyMatrix(row, data);
		data = CreateMatrix(rhs.row, rhs.column);
		row = rhs.row;
		column = rhs.column;
	}

	for (unsigned r = 0; r < row; ++r)
		memcpy(data[r], rhs.data[r], sizeof(float) * column);
			
	return *this;
}

Matrix& Matrix::operator -= (const Matrix &rhs)
// returns the matrix substraction of this and rhs
// NOTE: this and rhs must be the same size
// PRECONDITION:  matrices must be of same dimension.
// POSTCONDITION: returns this - rhs
{
	assert( (row == rhs.row) && (column == rhs.column) );
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			data[r][c] -= rhs.data[r][c];
	return *this;
}

Matrix& Matrix::operator += (const Matrix &rhs)
// returns the matrix addition of this and rhs
// NOTE: this and rhs must be the same size
// PRECONDITION:  matrices must be of same dimension.
// POSTCONDITION: returns this + rhs
{
	assert( (row == rhs.row) && (column == rhs.column) );
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			data[r][c] += rhs.data[r][c];
	return *this;
}

Matrix& Matrix::operator *= (const float scalar)
// does the scalar multiplication of this and a scalar, scalar
// PRECONDITION:  none.
// POSTCONDITION: returns the scalar multiplication of this and rhs
{
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			data[r][c] *= scalar;
	return *this;
}

Matrix& Matrix::operator /= (const float scalar)
// does the scalar division of this and a scalar, scalar
// PRECONDITION:  scalar != 0.
// POSTCONDITION: returns the scalar division of this and rhs
{
	assert(scalar != 0);
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			data[r][c] /= scalar;
	return *this;
}

float** Matrix::CreateMatrix(unsigned row, unsigned column) const
// creates a matrix that is row by column
// PRECONDITION:  none.
// POSTCONDITION: the matrix is returned
{
	float** temp;
	temp = new float*[row];
	for (unsigned index = 0; index < row; ++index)
		temp[index] = new float[column];
	return temp;
}

void Matrix::DestroyMatrix(unsigned row, float** myMatrix)
// destroys a matrix pointed to by myMatrix
// PRECONDITION:
// POSTCONDITION:
{
	for (unsigned r = 0; r < row; ++r)
		delete [] myMatrix[r];
	delete [] myMatrix;
}
/*
float Matrix::Determinant() const
// Finds the determinant of the matrix
// PRECONDITION:  row == column
// POSTCONDITION:
{
	assert(row == column);

	return 0;
}
*/
Matrix Matrix::GenerateIdentity(unsigned userRow, unsigned userColumn) const
// creates an identity matrix of size [userRow, userColumn]
// PRECONDITION:
// POSTCONDITION:
{
	Matrix temp(userRow, userColumn);
	for (unsigned i = 0; i < ((userRow < userColumn) ? userRow: userColumn); ++i)
		temp.data[i][i] = 1;
	return temp;
}

Matrix Matrix::Inverse() const
// returns an inverted form of the matrix
// PRECONDITION:  the current matrix is square
// POSTCONDITION: returns the inverted form of the current matrix
{
	assert(row == column);
	Matrix temp(*this);
	Matrix inverse = GenerateIdentity(3, 3);
	
	for (unsigned c = 0; c < column; ++c)
	{
		Matrix lhs = GenerateIdentity(row, column);
		if (temp.data[c][c] != 0)
			lhs.data[c][c] = 1 / temp.data[c][c];
		else 
			return Matrix(row, column);
		temp	= lhs * temp;
		inverse = lhs * inverse;

		for (unsigned r = 0; r < row; ++r)
			if (r != c)
			{
				Matrix lhs = GenerateIdentity(row, column);
				lhs.data[r][c] = -1 * temp.data[r][c];
				temp	= lhs * temp;
				inverse = lhs * inverse;
			}
	}
	return inverse;
}

Matrix Matrix::Reduce() const
// returns the Matrix in Row-Echelon Form 
// PRECONDITION:  none.
// POSTCONDITION: returns the reduced matrix of this
{
	Matrix temp(*this);

	for (unsigned c = 0; (c < column) && (c < row); ++c)
	{
		Matrix lhs = GenerateIdentity(row, column);
		if (temp.data[c][c] != 0)
			lhs.data[c][c] = 1 / temp.data[c][c];
		else 
			return temp;
		temp = lhs * temp;

		for (unsigned r = 0; r < row; ++r)
			if (r != c)
			{
				Matrix lhs = GenerateIdentity(row, column);
				lhs.data[r][c] = -1 * temp.data[r][c];
				temp = lhs * temp;
			}
	}
	return temp;
}

Matrix Matrix::Transpose() const
// returns the transposed matrix of this
// PRECONDITION:  none.
// POSTCONDITION: returns the transposition of the current matrix
{
	Matrix temp(column, row);
	for (unsigned r = 0; r < row; ++r)
		for (unsigned c = 0; c < column; ++c)
			temp.data[c][r] = data[r][c];
	return temp;
}