mlMatrix5.h

Go to the documentation of this file.

/*************************************************************************************
**
** Copyright 2007, MeVis Medical Solutions AG
**
** The user may use this file in accordance with the license agreement provided with
** the Software or, alternatively, in accordance with the terms contained in a
** written agreement between the user and MeVis Medical Solutions AG.
**
** For further information use the contact form at https://www.mevislab.de/contact
**
**************************************************************************************/
 
#ifndef ML_MATRIX5_H
#define ML_MATRIX5_H
 
 
// Include system independent file and project settings.
#include "mlLinearAlgebraSystem.h"
#include "mlLinearAlgebraDefs.h"
 
 
#include "mlFloatingPointMatrix.h"
 
#include "mlVector5.h"
 
#include "mlMatrix4.h"
 
#include <mlErrorMacros.h>
#include <mlErrorOutput.h>
 
// All declarations of this header will be in the ML_LINEAR_ALGEBRA namespace.
ML_LA_START_NAMESPACE
//--------------------------------------------------------------------
// This is necessary because we do not known whether vector or
// matrix header is included first and we cannot move template
// code into C++ file.
//--------------------------------------------------------------------
template <class DT> class Tvec5;
template <class DT> class Tmat4;
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
class Tmat5 : public FloatingPointMatrix<Tvec5<DT>, 5>
{
public:
 
  typedef DT ComponentType;
 
  // Builds a 5x5 matrix from 25 zero elements.
  Tmat5();
 
  // Builds a matrix that has the argument at the diagonal values, zero otherwise
  Tmat5(const DT diagValue);
 
  // Builds the matrix of the five row vectors row0, ..., row5.
  Tmat5(const Tvec5<DT> &row0, const Tvec5<DT> &row1, const Tvec5<DT> &row2,
        const Tvec5<DT> &row3, const Tvec5<DT> &row4);
 
  // Copy constructor from the Tmat5 mat.
  Tmat5(const Tmat5<DT> &mat);
 
  // Constructor from 25 floats given as array mat, row by row.
  Tmat5(const float mat[25]);
 
  // Constructor from 25 doubles given as array mat, row by row.
  Tmat5(const double mat[25]);
 
  // Initializes all matrix elements explicitly with scalars,
  // filling it row by row.
  Tmat5(const double in00, const double in01, const double in02, const double in03, const double in04,
        const double in10, const double in11, const double in12, const double in13, const double in14,
        const double in20, const double in21, const double in22, const double in23, const double in24,
        const double in30, const double in31, const double in32, const double in33, const double in34,
        const double in40, const double in41, const double in42, const double in43, const double in44);
 
  // Copies the contents from float array m into *this. m must point to an
  // array with at least 25 valid entries, row by row.
  void setValues(const float m[25]);
 
  // Copies the contents of *this into float array m. Note that range
  // and precision of the float values may not be sufficient for
  // higher precision matrix contents. m must point to an array with
  // at least 25 accessible entries, row by row.
  void getValues(float m[25]) const;
 
  // Copies the contents from double array m into *this. m must point to an
  // array with at least 25 valid entries, row by row.
  void setValues(const double m[25]);
 
  // Copies the contents of *this into double array m. m must point to an
  // array with at least 25 accessible entries, row by row.
  void getValues(double m[25]) const;
 
  // Sets the diagonal matrix with scale on diagonal.
  void setScaleMatrix(const DT scale);
 
  // Returns a matrix filled with values val.
  static Tmat5<DT> getMat(const double val);
 
  // Sets all values to val.
  void set(const DT val);
 
  bool operator<(const Tmat5<DT> &) const { return false; }
 
  // Assigns from a Tmat5.
  const Tmat5<DT>& operator=(const Tmat5<DT> &m);
 
  // Adds component wise by a Tmat5.
  const Tmat5<DT>& operator+=(const Tmat5<DT> &m);
 
  // Subtracts component wise by a Tmat5.
  const Tmat5<DT>& operator-=(const Tmat5<DT> &m);
 
  // Multiplies by a scalar constant
  const Tmat5<DT>& operator*=(const DT d);
 
  // Divides by a scalar constant.
  // Division by zero is not handled and must be avoided by caller.
  const Tmat5<DT>& operator/=(const DT d);
 
  //-------------------------------------------------
  // Special functions
  //-------------------------------------------------
  // Determines the (sub)determinant of columns given by col1, col2, col3 and col4.
  DT determinantLower4(const int col1, const int col2, const int col3, const int col4) const;
 
  // Determinant.
  DT det() const;
 
  // Returns the inverse Gauss-Jordan elimination with partial pivoting.
  // If a non-NULL Boolean pointer is passed to isInvertible
  // then true is returned in *isInvertible in the case of a
  // successful inversion or false if the inversion is not possible
  // (function return is the identity then).
  // If a NULL pointer is passed as isInvertible the matrix must
  // be invertible, otherwise errors will occur.
  Tmat5<DT> inverse(bool* isInvertible=nullptr) const;
 
  // Returns the transposed of this matrix.
  Tmat5 transpose() const;
 
  // Returns the identity matrix.
  static Tmat5 getIdentity();
 
  // Applies the function fct to each component.
  const Tmat5<DT> & apply(MLDblFuncPtr fct);
 
}; // end of class *Tmat5*
 
 
 
 
//---------------------------------------------------------------------
// Constructors
//---------------------------------------------------------------------
template <class DT>
inline Tmat5<DT>::Tmat5()
{
  this->v[0] = this->v[1] = this->v[2] = this->v[3] = this->v[4] = Tvec5<DT>(0);
}
 
// Constructs a matrix that has the argument at the diagonal values, zero otherwise
template <class DT>
inline Tmat5<DT>::Tmat5(const DT diagValue)
{
  this->v[0][0] = this->v[1][1] = this->v[2][2] = this->v[3][3] = this->v[4][4] = diagValue;
  this->v[1][0] = this->v[2][0] = this->v[3][0] = this->v[4][0] = this->v[2][1] = this->v[3][1] = this->v[4][1] = this->v[3][2] = this->v[4][2] = this->v[4][3] = 0;
  this->v[0][1] = this->v[0][2] = this->v[0][3] = this->v[0][4] = this->v[1][2] = this->v[1][3] = this->v[1][4] = this->v[2][3] = this->v[2][4] = this->v[3][4] = 0;
}
 
 
template <class DT>
inline Tmat5<DT>::Tmat5(const Tvec5<DT>& row0, const Tvec5<DT>& row1, const Tvec5<DT>& row2,
                        const Tvec5<DT>& row3, const Tvec5<DT>& row4)
{
  this->v[0] = row0;
  this->v[1] = row1;
  this->v[2] = row2;
  this->v[3] = row3;
  this->v[4] = row4;
}
 
template <class DT>
inline Tmat5<DT>::Tmat5(const Tmat5<DT> &mat)
{
  this->v[0] = mat.v[0];
  this->v[1] = mat.v[1];
  this->v[2] = mat.v[2];
  this->v[3] = mat.v[3];
  this->v[4] = mat.v[4];
}
 
template <class DT>
inline Tmat5<DT>::Tmat5(const float mat[25])
{
  setValues(mat);
}
 
template <class DT>
inline Tmat5<DT>::Tmat5(const double mat[25])
{
  setValues(mat);
}
 
template <class DT>
inline Tmat5<DT> Tmat5<DT>::getMat(const double val)
{
  return Tmat5<DT>(val, val, val, val, val,
                   val, val, val, val, val,
                   val, val, val, val, val,
                   val, val, val, val, val,
                   val, val, val, val, val);
}
 
template <class DT>
inline void Tmat5<DT>::set(const DT val)
{
  this->v[0] = this->v[1] = this->v[2] = this->v[3] = this->v[4] = Tvec5<DT>(val);
}
 
template <class DT>
inline const Tmat5<DT> & Tmat5<DT>::operator=(const Tmat5<DT> &m)
{
  if (&m != this){
    this->v[0] = m.v[0];
    this->v[1] = m.v[1];
    this->v[2] = m.v[2];
    this->v[3] = m.v[3];
    this->v[4] = m.v[4];
  }
 
  return *this;
}
 
template <class DT>
inline const Tmat5<DT> & Tmat5<DT>::operator+=(const Tmat5<DT> &m)
{
  this->v[0] += m.v[0];
  this->v[1] += m.v[1];
  this->v[2] += m.v[2];
  this->v[3] += m.v[3];
  this->v[4] += m.v[4];
 
  return *this;
}
 
template <class DT>
inline const Tmat5<DT> & Tmat5<DT>::operator-=(const Tmat5<DT> &m)
{
  this->v[0] -= m.v[0];
  this->v[1] -= m.v[1];
  this->v[2] -= m.v[2];
  this->v[3] -= m.v[3];
  this->v[4] -= m.v[4];
 
  return *this;
}
 
template <class DT>
inline const Tmat5<DT> & Tmat5<DT>::operator*=(const DT d)
{
  this->v[0] *= d;
  this->v[1] *= d;
  this->v[2] *= d;
  this->v[3] *= d;
  this->v[4] *= d;
 
  return *this;
}
 
template <class DT>
inline const Tmat5<DT> & Tmat5<DT>::operator/=(const DT d)
{
  ML_CHECK_FLOAT_THROW(d);
  this->v[0] /= d;
  this->v[1] /= d;
  this->v[2] /= d;
  this->v[3] /= d;
  this->v[4] /= d;
  return *this;
}
 
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
inline Tmat5<DT> operator-(const Tmat5<DT> &a)
{
  return Tmat5<DT>(a) *= static_cast<DT>(-1.0);
}
 
template <class DT>
inline Tmat5<DT> operator+(const Tmat5<DT> &a, const Tmat5<DT> &b)
{
  return Tmat5<DT>(a) += b;
}
 
template <class DT>
inline Tmat5<DT> operator-(const Tmat5<DT> &a, const Tmat5<DT> &b)
{
  return Tmat5<DT>(a) -= b;
}
 
template <class DT>
inline Tmat5<DT> operator*(const Tmat5<DT> &a, const DT d)
{
  return Tmat5<DT>(a) *= d;
}
 
template <class DT>
inline Tmat5<DT> operator*(const DT d, const Tmat5<DT> &a)
{
  return Tmat5<DT>(a) *= d;
}
 
template <class DT>
inline Tmat5<DT> operator/(const Tmat5<DT> &a, const DT d)
{
  Tmat5<DT> divided(a);
  ML_CHECK_FLOAT_THROW(d);
  divided /= d;
  return divided;
}
 
//-------------------------------------------------------------------
//
//      Tmat5 member functions
//
//-------------------------------------------------------------------
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
Tmat5<DT>::Tmat5(const double in00, const double in01, const double in02, const double in03, const double in04,
                 const double in10, const double in11, const double in12, const double in13, const double in14,
                 const double in20, const double in21, const double in22, const double in23, const double in24,
                 const double in30, const double in31, const double in32, const double in33, const double in34,
                 const double in40, const double in41, const double in42, const double in43, const double in44)
{
  this->v[0][0]=static_cast<DT>(in00); this->v[0][1]=static_cast<DT>(in01); this->v[0][2]=static_cast<DT>(in02);  this->v[0][3]=static_cast<DT>(in03); this->v[0][4]=static_cast<DT>(in04);
  this->v[1][0]=static_cast<DT>(in10); this->v[1][1]=static_cast<DT>(in11); this->v[1][2]=static_cast<DT>(in12);  this->v[1][3]=static_cast<DT>(in13); this->v[1][4]=static_cast<DT>(in14);
  this->v[2][0]=static_cast<DT>(in20); this->v[2][1]=static_cast<DT>(in21); this->v[2][2]=static_cast<DT>(in22);  this->v[2][3]=static_cast<DT>(in23); this->v[2][4]=static_cast<DT>(in24);
  this->v[3][0]=static_cast<DT>(in30); this->v[3][1]=static_cast<DT>(in31); this->v[3][2]=static_cast<DT>(in32);  this->v[3][3]=static_cast<DT>(in33); this->v[3][4]=static_cast<DT>(in34);
  this->v[4][0]=static_cast<DT>(in40); this->v[4][1]=static_cast<DT>(in41); this->v[4][2]=static_cast<DT>(in42);  this->v[4][3]=static_cast<DT>(in43); this->v[4][4]=static_cast<DT>(in44);
}
 
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
void Tmat5<DT>::setValues(const float m[25])
{
  this->v[0][0] = static_cast<DT>(m[ 0]);
  this->v[0][1] = static_cast<DT>(m[ 1]);
  this->v[0][2] = static_cast<DT>(m[ 2]);
  this->v[0][3] = static_cast<DT>(m[ 3]);
  this->v[0][4] = static_cast<DT>(m[ 4]);
 
  this->v[1][0] = static_cast<DT>(m[ 5]);
  this->v[1][1] = static_cast<DT>(m[ 6]);
  this->v[1][2] = static_cast<DT>(m[ 7]);
  this->v[1][3] = static_cast<DT>(m[ 8]);
  this->v[1][4] = static_cast<DT>(m[ 9]);
 
  this->v[2][0] = static_cast<DT>(m[10]);
  this->v[2][1] = static_cast<DT>(m[11]);
  this->v[2][2] = static_cast<DT>(m[12]);
  this->v[2][3] = static_cast<DT>(m[13]);
  this->v[2][4] = static_cast<DT>(m[14]);
 
  this->v[3][0] = static_cast<DT>(m[15]);
  this->v[3][1] = static_cast<DT>(m[16]);
  this->v[3][2] = static_cast<DT>(m[17]);
  this->v[3][3] = static_cast<DT>(m[18]);
  this->v[3][4] = static_cast<DT>(m[19]);
 
  this->v[4][0] = static_cast<DT>(m[20]);
  this->v[4][1] = static_cast<DT>(m[21]);
  this->v[4][2] = static_cast<DT>(m[22]);
  this->v[4][3] = static_cast<DT>(m[23]);
  this->v[4][4] = static_cast<DT>(m[24]);
}
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
void Tmat5<DT>::getValues(float m[25]) const
{
  m[ 0] = static_cast<float>(this->v[0][0]);
  m[ 1] = static_cast<float>(this->v[0][1]);
  m[ 2] = static_cast<float>(this->v[0][2]);
  m[ 3] = static_cast<float>(this->v[0][3]);
  m[ 4] = static_cast<float>(this->v[0][4]);
 
  m[ 5] = static_cast<float>(this->v[1][0]);
  m[ 6] = static_cast<float>(this->v[1][1]);
  m[ 7] = static_cast<float>(this->v[1][2]);
  m[ 8] = static_cast<float>(this->v[1][3]);
  m[ 9] = static_cast<float>(this->v[1][4]);
 
  m[10] = static_cast<float>(this->v[2][0]);
  m[11] = static_cast<float>(this->v[2][1]);
  m[12] = static_cast<float>(this->v[2][2]);
  m[13] = static_cast<float>(this->v[2][3]);
  m[14] = static_cast<float>(this->v[2][4]);
 
  m[15] = static_cast<float>(this->v[3][0]);
  m[16] = static_cast<float>(this->v[3][1]);
  m[17] = static_cast<float>(this->v[3][2]);
  m[18] = static_cast<float>(this->v[3][3]);
  m[19] = static_cast<float>(this->v[3][4]);
 
  m[20] = static_cast<float>(this->v[4][0]);
  m[21] = static_cast<float>(this->v[4][1]);
  m[22] = static_cast<float>(this->v[4][2]);
  m[23] = static_cast<float>(this->v[4][3]);
  m[24] = static_cast<float>(this->v[4][4]);
}
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
void Tmat5<DT>::setValues(const double m[25])
{
  this->v[0][0] = static_cast<DT>(m[ 0]);
  this->v[0][1] = static_cast<DT>(m[ 1]);
  this->v[0][2] = static_cast<DT>(m[ 2]);
  this->v[0][3] = static_cast<DT>(m[ 3]);
  this->v[0][4] = static_cast<DT>(m[ 4]);
 
  this->v[1][0] = static_cast<DT>(m[ 5]);
  this->v[1][1] = static_cast<DT>(m[ 6]);
  this->v[1][2] = static_cast<DT>(m[ 7]);
  this->v[1][3] = static_cast<DT>(m[ 8]);
  this->v[1][4] = static_cast<DT>(m[ 9]);
 
  this->v[2][0] = static_cast<DT>(m[10]);
  this->v[2][1] = static_cast<DT>(m[11]);
  this->v[2][2] = static_cast<DT>(m[12]);
  this->v[2][3] = static_cast<DT>(m[13]);
  this->v[2][4] = static_cast<DT>(m[14]);
 
  this->v[3][0] = static_cast<DT>(m[15]);
  this->v[3][1] = static_cast<DT>(m[16]);
  this->v[3][2] = static_cast<DT>(m[17]);
  this->v[3][3] = static_cast<DT>(m[18]);
  this->v[3][4] = static_cast<DT>(m[19]);
 
  this->v[4][0] = static_cast<DT>(m[20]);
  this->v[4][1] = static_cast<DT>(m[21]);
  this->v[4][2] = static_cast<DT>(m[22]);
  this->v[4][3] = static_cast<DT>(m[23]);
  this->v[4][4] = static_cast<DT>(m[24]);
}
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
void Tmat5<DT>::getValues(double m[25]) const
{
  m[ 0] = static_cast<double>(this->v[0][0]);
  m[ 1] = static_cast<double>(this->v[0][1]);
  m[ 2] = static_cast<double>(this->v[0][2]);
  m[ 3] = static_cast<double>(this->v[0][3]);
  m[ 4] = static_cast<double>(this->v[0][4]);
 
  m[ 5] = static_cast<double>(this->v[1][0]);
  m[ 6] = static_cast<double>(this->v[1][1]);
  m[ 7] = static_cast<double>(this->v[1][2]);
  m[ 8] = static_cast<double>(this->v[1][3]);
  m[ 9] = static_cast<double>(this->v[1][4]);
 
  m[10] = static_cast<double>(this->v[2][0]);
  m[11] = static_cast<double>(this->v[2][1]);
  m[12] = static_cast<double>(this->v[2][2]);
  m[13] = static_cast<double>(this->v[2][3]);
  m[14] = static_cast<double>(this->v[2][4]);
 
  m[15] = static_cast<double>(this->v[3][0]);
  m[16] = static_cast<double>(this->v[3][1]);
  m[17] = static_cast<double>(this->v[3][2]);
  m[18] = static_cast<double>(this->v[3][3]);
  m[19] = static_cast<double>(this->v[3][4]);
 
  m[20] = static_cast<double>(this->v[4][0]);
  m[21] = static_cast<double>(this->v[4][1]);
  m[22] = static_cast<double>(this->v[4][2]);
  m[23] = static_cast<double>(this->v[4][3]);
  m[24] = static_cast<double>(this->v[4][4]);
}
 
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
void Tmat5<DT>::setScaleMatrix(const DT scale)
{
  this->v[0][0] = scale;   this->v[0][1] =     0;   this->v[0][2] =     0;  this->v[0][3] =     0;   this->v[0][4] =     0;
  this->v[1][0] =     0;   this->v[1][1] = scale;   this->v[1][2] =     0;  this->v[1][3] =     0;   this->v[1][4] =     0;
  this->v[2][0] =     0;   this->v[2][1] =     0;   this->v[2][2] = scale;  this->v[2][3] =     0;   this->v[2][4] =     0;
  this->v[3][0] =     0;   this->v[3][1] =     0;   this->v[3][2] =     0;  this->v[3][3] = scale;   this->v[3][4] =     0;
  this->v[4][0] =     0;   this->v[4][1] =     0;   this->v[4][2] =     0;  this->v[4][3] =     0;   this->v[4][4] = scale;
}
 
 
 
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
inline DT Tmat5<DT>::determinantLower4(const int col1, const int col2, const int col3, const int col4) const
{
  /* Determinant of 4x4 matrix with given entries */
  return Tmat4<DT>(this->v[1][col1], this->v[1][col2], this->v[1][col3], this->v[1][col4],
                   this->v[2][col1], this->v[2][col2], this->v[2][col3], this->v[2][col4],
                   this->v[3][col1], this->v[3][col2], this->v[3][col3], this->v[3][col4],
                   this->v[4][col1], this->v[4][col2], this->v[4][col3], this->v[4][col4]).det();
}
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
DT Tmat5<DT>::det() const
{
  return (  this->v[0][0] * determinantLower4(1, 2, 3, 4)
          - this->v[0][1] * determinantLower4(0, 2, 3, 4)
          + this->v[0][2] * determinantLower4(0, 1, 3, 4)
          - this->v[0][3] * determinantLower4(0, 1, 2, 4)
          + this->v[0][4] * determinantLower4(0, 1, 2, 3));
}
 
//--------------------------------------------------------------------
//--------------------------------------------------------------------
template <class DT>
Tmat5<DT> Tmat5<DT>::inverse(bool* isInvertible) const
{
  // Epsilon for comparison with 0 in inversion process.
  static const DT Epsilon = static_cast<DT>(10e-13);
 
    // Use helper function from tools to invert the matrix.
  return MLInverseMatHelper(*this,
                            isInvertible,
                            Epsilon,
                            "Tmat5<DT> Tmat5<DT>::inverse(bool* isInvertible) const, matrix not invertable",
                            getIdentity(),
                            5);
}
 
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
inline Tmat5<DT> Tmat5<DT>::transpose() const
{
  return Tmat5<DT>(Tvec5<DT>(this->v[0][0], this->v[1][0], this->v[2][0], this->v[3][0], this->v[4][0]),
                   Tvec5<DT>(this->v[0][1], this->v[1][1], this->v[2][1], this->v[3][1], this->v[4][1]),
                   Tvec5<DT>(this->v[0][2], this->v[1][2], this->v[2][2], this->v[3][2], this->v[4][2]),
                   Tvec5<DT>(this->v[0][3], this->v[1][3], this->v[2][3], this->v[3][3], this->v[4][3]),
                   Tvec5<DT>(this->v[0][4], this->v[1][4], this->v[2][4], this->v[3][4], this->v[4][4]));
}
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
inline Tmat5<DT> Tmat5<DT>::getIdentity()
{
  return Tmat5<DT>(Tvec5<DT>(1,0,0,0,0),
                   Tvec5<DT>(0,1,0,0,0),
                   Tvec5<DT>(0,0,1,0,0),
                   Tvec5<DT>(0,0,0,1,0),
                   Tvec5<DT>(0,0,0,0,1));
}
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
const Tmat5<DT> & Tmat5<DT>::apply(MLDblFuncPtr fct)
{
  this->v[0].apply(fct);
  this->v[1].apply(fct);
  this->v[2].apply(fct);
  this->v[3].apply(fct);
  this->v[4].apply(fct);
  return *this;
}
 
 
#define _ML_MAT5_RC(i, j) a[i][0]*b[0][j] + a[i][1]*b[1][j] + a[i][2]*b[2][j] + \
a[i][3]*b[3][j] + a[i][4]*b[4][j]
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
Tmat5<DT> operator*(const Tmat5<DT> &a, const Tmat5<DT> &b)
{
  return Tmat5<DT>(Tvec5<DT>(_ML_MAT5_RC(0,0), _ML_MAT5_RC(0,1), _ML_MAT5_RC(0,2), _ML_MAT5_RC(0,3), _ML_MAT5_RC(0,4)),
                   Tvec5<DT>(_ML_MAT5_RC(1,0), _ML_MAT5_RC(1,1), _ML_MAT5_RC(1,2), _ML_MAT5_RC(1,3), _ML_MAT5_RC(1,4)),
                   Tvec5<DT>(_ML_MAT5_RC(2,0), _ML_MAT5_RC(2,1), _ML_MAT5_RC(2,2), _ML_MAT5_RC(2,3), _ML_MAT5_RC(2,4)),
                   Tvec5<DT>(_ML_MAT5_RC(3,0), _ML_MAT5_RC(3,1), _ML_MAT5_RC(3,2), _ML_MAT5_RC(3,3), _ML_MAT5_RC(3,4)),
                   Tvec5<DT>(_ML_MAT5_RC(4,0), _ML_MAT5_RC(4,1), _ML_MAT5_RC(4,2), _ML_MAT5_RC(4,3), _ML_MAT5_RC(4,4)));
}
#undef _ML_MAT5_RC
 
 
// Note: Multiplications with Vector5 are implemented with/in class Tvec5.
 
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
inline bool operator==(const Tmat5<DT> &a, const Tmat5<DT> &b)
{
  return (a[0] == b[0]) &&
  (a[1] == b[1]) &&
  (a[2] == b[2]) &&
  (a[3] == b[3]) &&
  (a[4] == b[4]);
}
 
//-------------------------------------------------------------------
//-------------------------------------------------------------------
template <class DT>
inline bool operator!=(const Tmat5<DT> &a, const Tmat5<DT> &b)
{
  return !(a == b);
}
 
 
//-----------------------------------------------------------------------------------
//-----------------------------------------------------------------------------------
 
typedef Tmat5<MLfloat>   Matrix5f;
typedef Tmat5<MLdouble>  Matrix5d;
typedef Tmat5<MLldouble> Matrix5ld;
typedef Tmat5<MLdouble>  Matrix5;
 
 
ML_LA_END_NAMESPACE
 
namespace std
{
  //-------------------------------------------------------------------
  //-------------------------------------------------------------------
  template <class DT>
  inline std::ostream& operator<<(std::ostream& os, const ML_LA_NAMESPACE::Tmat5<DT> & m)
  {
    return os << m[0] << '\n' << m[1] << '\n' << m[2] << '\n' << m[3] << '\n' << m[4];
  }
 
  //-------------------------------------------------------------------
  //-------------------------------------------------------------------
  template <class DT>
  inline std::istream& operator>>(std::istream& is, ML_LA_NAMESPACE::Tmat5<DT> & m)
  {
    ML_LA_NAMESPACE::Tmat5<DT> m_tmp;
 
    is >> m_tmp[0] >> m_tmp[1] >> m_tmp[2] >> m_tmp[3] >> m_tmp[4];
    if (is){ m = m_tmp; }
    return is;
  }
}
 
 
#endif // __mlMatrix5_H