matrix.hpp

一个gps小工具包

#pragma ident "$Id: Matrix.hpp 70 2006-08-01 18:36:21Z ehagen $"



/**
 * @file Matrix.hpp
 * Basic Matrix algorithms
 */
 
#ifndef GPSTK_MATRIX_HPP
#define GPSTK_MATRIX_HPP

//============================================================================
//
//  This file is part of GPSTk, the GPS Toolkit.
//
//  The GPSTk is free software; you can redistribute it and/or modify
//  it under the terms of the GNU Lesser General Public License as published
//  by the Free Software Foundation; either version 2.1 of the License, or
//  any later version.
//
//  The GPSTk is distributed in the hope that it will be useful,
//  but WITHOUT ANY WARRANTY; without even the implied warranty of
//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//  GNU Lesser General Public License for more details.
//
//  You should have received a copy of the GNU Lesser General Public
//  License along with GPSTk; if not, write to the Free Software Foundation,
//  Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
//  
//  Copyright 2004, The University of Texas at Austin
//
//============================================================================

#include "Vector.hpp"
#include "MatrixBase.hpp"

namespace gpstk
{

   /** @defgroup VectorGroup Vector and Matrix mathematics */
   //@{
 
// forward declarations
   template <class T> class MatrixRowSlice;
   template <class T> class ConstMatrixRowSlice;
   template <class T> class MatrixColSlice;
   template <class T> class ConstMatrixColSlice;


/**
 * An implementation of a matrix class using Vector<T> as its internal basis.
 * This class is STL compliant with the iterator proceeding in row major order.
 * Operators +=, -=, *= and /= are implemented in RefMatrixBase.
 * 
 * @sa matvectest.cpp for examples
 */
   template <class T>
   class Matrix : public RefMatrixBase<T, Matrix<T> >
   {
   public:
         /// STL value_type
      typedef typename Vector<T>::value_type value_type;
         /// STL reference type
      typedef typename Vector<T>::reference reference;
         /// STL const reference type
      typedef typename Vector<T>::const_reference const_reference;
         /// STL iterator type
      typedef typename Vector<T>::iterator iterator;
         /// STL const iterator type
      typedef typename Vector<T>::const_iterator const_iterator;  

         /// default constructor
      Matrix();
         /// constructor given an initial size
      Matrix(size_t rows, size_t cols);
         /// constructor for an initial size and value
      Matrix(size_t rows, size_t cols, T initialValue);
         /// copies out the contents of vec to initialize the matrix
      Matrix(size_t rows, size_t cols, const T* vec);
         /// copies out the contents of vec to initialize the matrix
      template <class BaseClass>
      Matrix(size_t rows, size_t cols, const ConstVectorBase<T, BaseClass>& vec)
         : v(vec), r(rows), c(cols), s(rows * cols)
      {}
         /// constructor for a ConstMatrixBase object
      template <class BaseClass>
      Matrix(const ConstMatrixBase<T, BaseClass>& mat) 
            : v(mat.size()), r(mat.rows()), c(mat.cols()), s(mat.size())
         {
            size_t i,j;
            for(i = 0; i < r; i++)
               for(j = 0; j < c; j++)
                  (*this)(i,j) = mat(i, j);
         }

         /// submatrix constructor
      template <class BaseClass>
      Matrix(const ConstMatrixBase<T, BaseClass>& mat, size_t topRow, 
          size_t topCol, size_t numRows, size_t numCols) 
            : v((size_t)0), r(0), c(0), s(0)
         {
               // sanity checks...
            if ( (topCol > mat.cols()) || 
                 (topRow > mat.rows()) ||
                 ((topRow + numRows) > mat.rows()) ||
                 ((topCol + numCols) > mat.cols()) )
            {
               MatrixException e("Invalid dimensions or size for Matrix(MatrixBase)");
               GPSTK_THROW(e);
            }
         
               // seems ok - make the valarray and copy column by column
            r = numRows;
            c = numCols;
            s = r * c;
            v.resize(r * c);
            size_t i, j;
            for(i = 0; i < r; i++)
               for(j = 0; j < c; j++)
                  (*this)(i,j) = mat(topRow + i, topCol + j);
         }

         /// STL begin
      iterator begin() { return v.begin(); }
         /// STL const begin
      const_iterator begin() const { return v.begin(); }
         /// STL end
      iterator end() { return v.end(); }
         /// STL const end
      const_iterator end() const { return v.end(); }
         /// STL front
      value_type front() { return v.front(); }
         /// STL const front
      const_reference front() const { return v.front();}
         /// STL empty
      bool empty() const { return s == 0; }
         /// STL size
      size_t size() const {return s; }
         /// STL max size
      size_t max_size() const { return s; }

         /// The number of rows in the matrix
      inline size_t rows() const { return r; }
         /// The number of columns in the matrix
      inline size_t cols() const { return c; }
         /// A reference slice of a row with a given std::slice
      inline MatrixRowSlice<T> rowRef(size_t rowNum, const std::slice& s);
         /// A reference slice of a row with a starting column (i.e. sub-row)
      inline MatrixRowSlice<T> rowRef(size_t rowNum, size_t colNum = 0);
         /// A const reference slice of a row with a given std::slice
      inline ConstMatrixRowSlice<T> row(size_t rowNum, const std::slice& s) const;
         /// A const reference slice of a row with a starting column (i.e. sub-row)
      inline ConstMatrixRowSlice<T> row(size_t rowNum, size_t colNum = 0) const;

         /// A reference column with a given slice
      inline MatrixColSlice<T> colRef(size_t colNum, const std::slice& s);
         /// A reference column with a starting row number (i.e. sub-column)
      inline MatrixColSlice<T> colRef(size_t colNum, size_t rowNum = 0);
         /// A const reference column with a given slice
      inline ConstMatrixColSlice<T> col(size_t colNum, const std::slice& s) const;
         /// A const reference column with a starting row number (i.e. sub-column)
      inline ConstMatrixColSlice<T> col(size_t colNum, size_t rowNum = 0) const;

         /// Non-const matrix operator(row,col)
      inline T& operator() (size_t rowNum, size_t colNum)
         { return v(rowNum + colNum * r); }
         /// Const matrix operator(row,col)
      inline T operator() (size_t rowNum, size_t colNum) const
         { return v(rowNum + colNum * r); }
         /// operator[] that returns a row slice
      inline MatrixRowSlice<T> operator[] (size_t row)
         { return rowRef(row); }
         /// const operator[] that returns a const row slice
      inline ConstMatrixRowSlice<T> operator[] (size_t rowNum) const 
         { return row(rowNum);}

         /// Resizes the matrix to rows*cols.
         /// @warning YOUR DATA MAY NOT BE RETAINED!!!
      inline Matrix& resize(size_t rows, size_t cols);

      inline Matrix& resize(size_t rows, size_t cols, 
                         const T initialValue);

         /**
          * Assigns this matrix to a T* in column major order.
          * @warning be careful that array is as large as the matrix is!
          */
      inline Matrix& operator=(const T* array)
         { return assignFrom(array); }
         /// Assigns the contents of this matrix to those in array in column
         /// major order.
      inline Matrix& operator=(const std::valarray<T> array)
         { return assignFrom(array); }
         /// Assigns all elements of the matrix to \c t.
      inline Matrix& operator=(const T t)
         { return assignFrom(t); }
         /// Copies the other matrix.
      inline Matrix& operator=(const Matrix& mat)
         { v = mat.v; r = mat.r; c = mat.c; s = mat.s; return *this; }
         /// Copies from any matrix.
      template <class BaseClass>
      inline Matrix& operator=(const ConstMatrixBase<T, BaseClass>& mat)
         { 
            v.resize(mat.size()); 
            r=mat.rows(); 
            c=mat.cols(); 
            s=mat.size();
            return assignFrom(mat);
         }
         /// Copies from any vector.
      template <class BaseClass>
      inline Matrix& operator=(const ConstVectorBase<T, BaseClass>& mat)
         { return assignFrom(mat); }

   private:
         /// the matrix stored in column major order
      Vector<T> v;
      size_t r,  ///< the number of rows
         c,  ///< the number of columns
         s;  ///< the overall size
   };

/**
 * An assignable slice of a matrix.
 */
   template <class T>
   class MatrixSlice : public RefMatrixSliceBase<T, MatrixSlice<T> >
   {
   public:
         /// default constructor
      MatrixSlice() : m(NULL), rSlice(std::slice(0,0,0)), 
         cSlice(std::slice(0,0,0)), s(0)
         {}

         /// Makes a slice of the whole matrix.
      MatrixSlice(Matrix<T>& mat)
            : m(&mat), rSlice(std::slice(0, mat.rows(), 1)),
              cSlice(std::slice(0,mat.cols(), 1)), s(mat.size())
         {
            matSliceCheck(mat.rows(), mat.cols());
         }

         /// Makes a partial slice of a matrix.
      MatrixSlice(Matrix<T>& mat, const std::slice& rowSlice,
               const std::slice& colSlice)
            : m(&mat), rSlice(rowSlice), cSlice(colSlice),
              s(rSlice.size() * cSlice.size())
         {
            matSliceCheck(mat.rows(), mat.cols());
         }

         /// Submatrix slice.
      MatrixSlice(Matrix<T>& mat, size_t topRow, size_t topCol, 
               size_t numRows, size_t numCols)
            : m(&mat), rSlice(std::slice(topRow, numRows, 1)),
              cSlice(std::slice(topCol, numCols, 1)),
              s(rSlice.size() * cSlice.size())
         {
            matSliceCheck(mat.rows(), mat.cols());
         }
      
         /// Copies from x to (*this).
      template <class V>
      MatrixSlice& operator=(const ConstMatrixBase<T, V>& x)
         { return assignFrom(x); }

         /// Copies from x to (*this).
      template <class V>
      MatrixSlice& operator=(const ConstVectorBase<T, V>& x)
         { return assignFrom(x); }

         /// Copies from x to (*this).
      MatrixSlice& operator=(const std::valarray<T>& x)
         { return assignFrom(x); }
         /// Copies from x to (*this).
      MatrixSlice& operator=(const T x)
         { return assignFrom(x); }
         /// Copies from x to (*this).
      MatrixSlice& operator=(const T* x)
         { return assignFrom(x); }

         /// returns the size of this slice
      size_t size() const { return s; }
         /// returns the number of columns in the slice
      size_t cols() const { return colSize(); }
         /// returns the number of rows in the slice
      size_t rows() const { return rowSize(); }
         /// returns the (i,j) element of the slice.
      T& operator() (size_t i, size_t j)
         { return (*m)(i * rowStride() + rowStart(), 
                       j * colStride() + colStart()); }
         /// returns the (i,j) element of the slice, const version.
      T operator() (size_t i, size_t j) const
         { return (*m)(i * rowStride() + rowStart(), 
                       j * colStride() + colStart()); }


         /// returns the number of rows in this slice
      size_t rowSize() const { return rSlice.size(); }
         /// returns the starting row in the base matrix of this slice
      size_t rowStart() const{ return rSlice.start(); }
         /// returns the number of elements between the i'th and i+1'th row
      size_t rowStride() const { return rSlice.stride(); }
         /// returns the number of columns in this slice
      size_t colSize() const { return cSlice.size(); }
         /// returns the starting row in the base matrix of this slice
      size_t colStart() const { return cSlice.start(); }
         /// returns the number of elements between the i'th and i+1'th row
      size_t colStride() const { return cSlice.stride(); }

   private:
         /// The matrix this slice refers to
      Matrix<T>* m;
      std::slice rSlice, ///< a row slice
         cSlice; ///< a column slice
      size_t s; ///< the overall size
   };

/**
 * An unmodifiable matrix slice.
 */
   template <class T>
   class ConstMatrixSlice : public ConstMatrixSliceBase<T, ConstMatrixSlice<T> >
   {
   public:
         /// default constructor
      ConstMatrixSlice(void) : m(NULL), rSlice(std::slice(0,0,0)), 
         cSlice(std::slice(0,0,0)), s(0)
         {}

         /// makes a const slice of the whole matrix
      ConstMatrixSlice(const Matrix<T>& mat)
            : m(&mat), rSlice(std::slice(0, mat.rows(), 1)),
              cSlice(std::slice(0,mat.cols(), 1)), s(mat.size())
         {
            matSliceCheck(mat.rows(), mat.cols());
         }

         /// makes a slice given std::slices for rows and columns
      ConstMatrixSlice(const Matrix<T>& mat, const std::slice& rowSlice,
               const std::slice& colSlice)
            : m(&mat), rSlice(rowSlice), cSlice(colSlice),
              s(rSlice.size() * cSlice.size())
         {
            matSliceCheck(mat.rows(), mat.cols());
         }

         /// submatrix slice
      ConstMatrixSlice(const Matrix<T>& mat, size_t topRow, size_t topCol, 
               size_t numRows, size_t numCols)
            : m(&mat), rSlice(std::slice(topRow, numRows, 1)),
              cSlice(std::slice(topCol, numCols, 1)),
              s(rSlice.size() * cSlice.size())
         {
            matSliceCheck(mat.rows(), mat.cols());
         }