TAPsMatrix.cpp

Go to the documentation of this file.

/******************************************************************************
TAPsMatrix.cpp

Matrix class is a class for any dimension Matrix.

SUKITTI PUNAK   (01/19/2006)
UPDATE          (01/19/2006)
******************************************************************************/
#include "TAPsMatrix.hpp"
// Using Inclusion Model (i.e. definitions are included in declarations)
//                       (this name.cpp is included in name.hpp)
// Each friend is defined directly inside its declaration.

BEGIN_NAMESPACE_TAPs
//=============================================================================
// Constructors and Destructor
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::Matrix ()
{
    // Identity Matrix
    MakeIdentity();
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::Matrix ( Matrix<T,R,C> const &M )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] = M.e[i];
}
/*
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::Matrix ( Matrix<T,C,R> const &M )
{
    for ( int r = 0; r < R; ++r ) {
        for ( int c = 0; c < C; ++c ) {
            e[r*C + c] = M.e[c*R + r];
        }
    }
}
//*/
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::Matrix ( T val )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] = val;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::Matrix ( T const a[] )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] = a[i];
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::~Matrix ()
{}
//-----------------------------------------------------------------------------
//=============================================================================
// Member Access
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
inline T & Matrix<T,R,C>::operator[] ( int i )
{   return e[i];    }
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
inline T const & Matrix<T,R,C>::operator[] ( int i ) const
{   return e[i];    }
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
inline T & Matrix<T,R,C>::operator() ( int r, int c )
{   return e[r*C + c];  }
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
inline T const & Matrix<T,R,C>::operator() ( int r, int c ) const
{   return e[r*C + c];  }
//-----------------------------------------------------------------------------
//=============================================================================
// Convert to one dimension array
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::operator const T * () const
{   return e;   }
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C>::operator T * ()
{   return e;   }
//-----------------------------------------------------------------------------
//=============================================================================
// Useful Functions
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
void Matrix<T,R,C>::SetAllElements (    T t )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] = t;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
void Matrix<T,R,C>::SetAllElements ( T const a[] )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] = a[i];
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
void Matrix<T,R,C>::MakeIdentity ()
{
    // Identity Matrix
    for ( int r = 0; r < R*C; r+=C )
        for ( int c = 0; c < C; ++c )
            e[r + c] = Math<T>::ZERO;
    if ( R <= C ) {
        for ( int r = 0; r < R; ++r )
            e[r*C + r] = Math<T>::ONE;
    }
    else {
        for ( int r = 0; r < C; ++r )
            e[r*C + r] = Math<T>::ONE;
    }
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
void Matrix<T,R,C>::MakeDiagonal ( T d )
{
    // Diagonal Matrix
    for ( int r = 0; r < R*C; r+=C )
        for ( int c = 0; c < C; ++c )
            e[r + c] = Math<T>::ZERO;
    for ( int r = 0; r < R; ++r )
            e[r*C + r] = d;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
void Matrix<T,R,C>::MakeDiagonal ( T const a[] )
{
    // Diagonal Matrix
    for ( int r = 0; r < R*C; r+=C )
        for ( int c = 0; c < C; ++c )
            e[r + c] = Math<T>::ZERO;
    for ( int r = 0; r < R; ++r )
            e[r*C + r] = a[r];
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
void Matrix<T,R,C>::MakeZero ()
{
    for ( int i = 0; i < R*C; ++i )
        e[i] = Math<T>::ZERO;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
bool Matrix<T,R,C>::IsIdentity () const
{
    int i;
    for ( int r = 0; r < R; ++r ) {
        for ( int c = 0; c < C; ++c ) {
            i = r*C;
            if ( r != c ) {
                if ( e[i + c] != Math<T>::ZERO ) {
                    return false;
                }
            }
            else {
                if ( e[i + c] != Math<T>::ONE ) {
                    return false;
                }
            }
        }
    }
    return true;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
bool Matrix<T,R,C>::IsSymmetric () const
{
    if ( R != C )   return false;   // Not an N-by-N matrix
    for ( int r = 0; r < R; ++r )
        for ( int c = r; c < C; ++c )
            if ( e[r*C + c] != e[c*R + r] )
                return false;
    return true;
}
//=============================================================================
// Matrix Operations
//-----------------------------------------------------------------------------
/*
// Transpose this matrix
template <typename T, int R, int C>
Matrix<T,C,R> & Matrix<T,R,C>::Transposed ()
{
    T temp; int i, j;
    if ( R >= C ) {     // If row >= col
        for ( int r = 0; r < R; ++r ) {
            for ( int c = r; c < C; ++c ) {
                i = r*C + c;
                j = c*R + r;
                e[i] = temp;
                e[i] = e[j];
                e[j] = temp;
            }
        }
    }
    else {              // If col > row
        for ( int c = 0; c < C; ++c ) {
            for ( int r = c; r < R; ++r ) {
                i = r*C + c;
                j = c*R + r;
                e[i] = temp;
                e[i] = e[j];
                e[j] = temp;
            }
        }
    }
    int r = R;
    R = C;
    C = r;
    return *this;
}
//*/
//-----------------------------------------------------------------------------
// Get the transpose matrix of this matrix
template <typename T, int R, int C>
Matrix<T,C,R> Matrix<T,R,C>::GetTranspose() const
{
    Matrix<T,C,R> O;
    int i;
    for ( int r = 0; r < R; ++r ) {
        i = r*C;
        for ( int c = 0; c < C; ++c ) {
            O[c*R + r] = e[i + c];
        }
    }
    return O;
}
/*
//-----------------------------------------------------------------------------
// Inverse the matrix
template <typename T, int R, int C>
inline Matrix<T,R,C> & Matrix<T,R,C>::Inversed ()
{
    *this = (*this).GetInverse();
    return *this;
}
//-----------------------------------------------------------------------------
// Get the matrix inverse
template <typename T, int R, int C>
Matrix<T,R,C> Matrix<T,R,C>::GetInverse () const
{
    T det = GetDeterminant();

    if ( -Math<T>::EPSILON < det && det < Math<T>::EPSILON )
        return Matrix( Math<T>::INFINITY );

    return Matrix<T,R,C>(    ( e[4]*e[8] - e[7]*e[5] ) / det,
                            -( e[1]*e[8] - e[7]*e[2] ) / det,
                             ( e[1]*e[5] - e[4]*e[2] ) / det,
                            -( e[3]*e[8] - e[6]*e[5] ) / det,
                             ( e[0]*e[8] - e[6]*e[2] ) / det,
                            -( e[0]*e[5] - e[3]*e[2] ) / det,
                             ( e[3]*e[7] - e[6]*e[4] ) / det,
                            -( e[0]*e[7] - e[6]*e[1] ) / det,
                             ( e[0]*e[4] - e[3]*e[1] ) / det    );
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
inline T Matrix<T,R,C>::GetDeterminant () const
{
    return      e[0] * ( e[4]*e[8] - e[5]*e[7] ) 
            -   e[3] * ( e[1]*e[8] - e[2]*e[7] )
            +   e[6] * ( e[1]*e[5] - e[2]*e[4] );
}
//*/
//-----------------------------------------------------------------------------
//=============================================================================
// Assignment Overloaded Operator
//--------------------------------------------------------------------------
// Assignment Operator
template <typename T, int R, int C>
inline Matrix<T,R,C> & Matrix<T,R,C>::operator= ( Matrix<T,R,C> const &M )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] = M.e[i];
    return *this;
}
//-----------------------------------------------------------------------------
//=============================================================================
// Unary Overloaded Operators
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> Matrix<T,R,C>::operator- ()
{
    Matrix<T,R,C> O;
    for ( int i = 0; i < R*C; ++i )
        O.e[i] = -e[i];
    return O;
}
//=============================================================================
// Assign Overloaded Operators
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> & Matrix<T,R,C>::operator+= ( Matrix<T,R,C> const &M )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] += M.e[i];
    return *this;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> & Matrix<T,R,C>::operator-= ( Matrix<T,R,C> const &M )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] -= M.e[i];
    return *this;
}
//-----------------------------------------------------------------------------
/*
template <typename T, int R, int C>
Matrix<T,R,C> & Matrix<T,R,C>::operator*= ( Matrix<T,C,R> const &M )
{
    //*this = (*this) * M;
    return *this;
}
//*/
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> & Matrix<T,R,C>::operator*= ( T s )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] *= s;
    return *this;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> & Matrix<T,R,C>::operator/= ( T s )
{
    for ( int i = 0; i < R*C; ++i )
        e[i] /= s;
    return *this;
}
//-----------------------------------------------------------------------------
//=============================================================================
// Binary Overloaded Operators
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> Matrix<T,R,C>::operator+ ( Matrix<T,R,C> const &M ) const
{
    Matrix<T,R,C> O;
    for ( int i = 0; i < R*C; ++i )
        O.e[i] = e[i] + M.e[i];
    return O;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> Matrix<T,R,C>::operator- ( Matrix<T,R,C> const &M ) const
{
    Matrix<T,R,C> O;
    for ( int i = 0; i < R*C; ++i )
        O.e[i] = e[i] - M.e[i];
    return O;
}
/*
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> Matrix<T,R,C>::operator* ( Matrix<T,C,R> const &M ) const
{
    Matrix<T,R,C> O;
    return O;
}
//*/
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> Matrix<T,R,C>::operator* ( T s ) const
{
    Matrix<T,R,C> O;
    for ( int i = 0; i < R*C; ++i )
        O.e[i] = e[i] * s;
    return O;
}
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
Matrix<T,R,C> Matrix<T,R,C>::operator/ ( T s ) const
{
    Matrix<T,R,C> O;
    for ( int i = 0; i < R*C; ++i )
        O.e[i] = e[i] / s;
    return O;
}
/*
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
inline void Matrix<T,R,C>::MultLeft ( Matrix<T,R,C> const &M )
{
    *this = M * (*this);
}
//*/
/*
//-----------------------------------------------------------------------------
template <typename T, int R, int C>
inline void Matrix<T,R,C>::MultRight ( Matrix<T,R,C> const &M )
{
    *this *= M;
}
//*/
//-----------------------------------------------------------------------------
// Matrix * Vector
template <typename T, int R, int C>
Vector<T,R> Matrix<T,R,C>::operator* ( Vector<T,C> const &V ) const
{
    Vector<T,R> O;
    for ( int r = 0, i = 0; r < R; ++r, i+=C ) {
        for ( int c = 0; c < C; ++c ) {
            O[r] += e[i + c] * V[c];
        }
    }
    return O;
}
//-----------------------------------------------------------------------------
//*/
//=============================================================================
END_NAMESPACE_TAPs
//345678901234567890123456789012345678901234567890123456789012345678901234567890
//--+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8