mitkMatrix.h

Go to the documentation of this file.

/*============================================================================

The Medical Imaging Interaction Toolkit (MITK)

Copyright (c) German Cancer Research Center (DKFZ)
All rights reserved.

Use of this source code is governed by a 3-clause BSD license that can be
found in the LICENSE file.

============================================================================*/

#ifndef MITKMATRIX_H_
#define MITKMATRIX_H_

#include <itkMatrix.h>

#include "mitkArray.h"
#include "mitkEqual.h"
#include "mitkNumericConstants.h"

namespace mitk
{
  template <class T, unsigned int NRows = 3, unsigned int NColumns = 3>
  class Matrix : public itk::Matrix<T, NRows, NColumns>
  {
  public:
    typedef Matrix Self;

    typedef typename itk::Matrix<T, NRows, NColumns>::InternalMatrixType InternalMatrixType;

    explicit Matrix<T, NRows, NColumns>() : itk::Matrix<T, NRows, NColumns>() {}
    explicit Matrix<T, NRows, NColumns>(const Matrix &matrix) : itk::Matrix<T, NRows, NColumns>(matrix) {}
    Matrix<T, NRows, NColumns>(const itk::Matrix<T, NRows, NColumns> &matrix) : itk::Matrix<T, NRows, NColumns>(matrix)
    {
    }
    inline Matrix<T, NRows, NColumns>(const vnl_matrix<T> &matrix) : itk::Matrix<T, NRows, NColumns>(matrix) {}
    inline explicit Matrix<T, NRows, NColumns>(InternalMatrixType &matrix) : itk::Matrix<T, NRows, NColumns>(matrix) {}
    using itk::Matrix<T, NRows, NColumns>::operator=;

    template <typename ArrayType>
    void FillMatrix(const ArrayType &array)
    {
      for (unsigned i = 0; i < NRows; i++)
      {
        for (unsigned j = 0; j < NColumns; j++)
        {
          (*this)[i][j] = array[i][j];
        }
      }
    };

    template <typename MatrixType>
    void ToArray(MatrixType matrix) const
    {
      for (unsigned i = 0; i < NRows; i++)
      {
        for (unsigned j = 0; j < NColumns; j++)
        {
          matrix[i][j] = (*this)[i][j];
        }
      }
    }
  };

  typedef Matrix<ScalarType, 2, 2> Matrix2D;
  typedef Matrix<ScalarType, 3, 3> Matrix3D;
  typedef Matrix<ScalarType, 4, 4> Matrix4D;

  template <typename TCoordRep, unsigned int NRows, unsigned int NCols>
  inline bool MatrixEqualRMS(const vnl_matrix_fixed<TCoordRep, NRows, NCols> &matrix1,
                             const vnl_matrix_fixed<TCoordRep, NRows, NCols> &matrix2,
                             mitk::ScalarType epsilon = mitk::eps)
  {
    if ((matrix1.rows() == matrix2.rows()) && (matrix1.cols() == matrix2.cols()))
    {
      vnl_matrix_fixed<TCoordRep, NRows, NCols> differenceMatrix = matrix1 - matrix2;
      if (differenceMatrix.rms() < epsilon)
      {
        return true;
      }
      else
      {
        return false;
      }
    }
    else
    {
      return false;
    }
  }

  template <typename TCoordRep, unsigned int NRows, unsigned int NCols>
  inline bool MatrixEqualRMS(const itk::Matrix<TCoordRep, NRows, NCols> &matrix1,
                             const itk::Matrix<TCoordRep, NRows, NCols> &matrix2,
                             mitk::ScalarType epsilon = mitk::eps)
  {
    return mitk::MatrixEqualRMS(matrix1.GetVnlMatrix(), matrix2.GetVnlMatrix(), epsilon);
  }

  template <typename TCoordRep, unsigned int NRows, unsigned int NCols>
  inline bool MatrixEqualElementWise(const vnl_matrix_fixed<TCoordRep, NRows, NCols> &matrix1,
                                     const vnl_matrix_fixed<TCoordRep, NRows, NCols> &matrix2,
                                     mitk::ScalarType epsilon = mitk::eps)
  {
    if ((matrix1.rows() == matrix2.rows()) && (matrix1.cols() == matrix2.cols()))
    {
      for (unsigned int r = 0; r < NRows; r++)
      {
        for (unsigned int c = 0; c < NCols; c++)
        {
          TCoordRep difference = matrix1(r, c) - matrix2(r, c);
          if (DifferenceBiggerOrEqualEps(difference, epsilon))
          {
            return false;
          }
        }
      }
      return true;
    }
    else
    {
      return false;
    }
  }

  template <typename TCoordRep, unsigned int NRows, unsigned int NCols>
  inline bool MatrixEqualElementWise(const itk::Matrix<TCoordRep, NRows, NCols> &matrix1,
                                     const itk::Matrix<TCoordRep, NRows, NCols> &matrix2,
                                     mitk::ScalarType epsilon = mitk::eps)
  {
    return mitk::MatrixEqualElementWise(matrix1.GetVnlMatrix(), matrix2.GetVnlMatrix(), epsilon);
  }
}

#endif /* MITKMATRIX_H_ */