// @(#)root/smatrix:$Name: $:$Id: SMatrix.h,v 1.10 2005/12/13 18:28:09 moneta Exp $
// Authors: T. Glebe, L. Moneta 2005
#ifndef ROOT_Math_SMatrix
#define ROOT_Math_SMatrix
// ********************************************************************
//
// source:
//
// type: source code
//
// created: 20. Mar 2001
//
// author: Thorsten Glebe
// HERA-B Collaboration
// Max-Planck-Institut fuer Kernphysik
// Saupfercheckweg 1
// 69117 Heidelberg
// Germany
// E-mail: T.Glebe@mpi-hd.mpg.de
//
// Description: A fixed size two dimensional Matrix class
//
// changes:
// 20 Mar 2001 (TG) creation
// 21 Mar 2001 (TG) added operators +=, -=, *=, /=
// 26 Mar 2001 (TG) place_in_row(), place_in_col() added
// 02 Apr 2001 (TG) non-const Array() added
// 03 Apr 2001 (TG) invert() added
// 07 Apr 2001 (TG) CTOR from SVertex (dyadic product) added
// 09 Apr 2001 (TG) CTOR from array added
// 11 Apr 2001 (TG) rows(), cols(), size() replaced by rows, cols, size
// 25 Mai 2001 (TG) row(), col() added
// 04 Sep 2001 (TG) moved inlined functions to .icc file
// 11 Jan 2002 (TG) added operator==(), operator!=()
// 14 Jan 2002 (TG) added more operator==(), operator!=(), operator>(), operator<()
//
// ********************************************************************
// for platform specific configurations
#include "Math/MConfig.h"
#include <iosfwd>
// expression engine
#include "Math/Expression.h"
//doxygen tag
/**
@defgroup SMatrix Matrix and Vector classes
*/
namespace ROOT {
namespace Math {
template <class T, unsigned int D> class SVector;
/**
SMatrix: a generic fixed size n x m Matrix class.
The class is template on the scalar type and on the matrix sizes:
D1 = number of rows and D2 = number of columns.
@ingroup SMatrix
@memo SMatrix
@author T. Glebe
*/
//==============================================================================
// SMatrix: column-wise storage
//==============================================================================
template <class T, unsigned int D1, unsigned int D2 = D1>
class SMatrix {
public:
/** @name --- Typedefs --- */
///
typedef T value_type;
/** STL iterator interface. */
typedef T* iterator;
/** STL const_iterator interface. */
typedef const T* const_iterator;
/** @name --- Constructors --- */
/**
Default constructor:
*/
SMatrix();
///
SMatrix(const SMatrix<T,D1,D2>& rhs);
///
template <class A>
SMatrix(const Expr<A,T,D1,D2>& rhs);
// new constructs using STL iterator interface
/**
* Constructor with STL iterator interface. The data will be copied into the matrix
* Size of the matrix must match size of the iterators
*/
template<class InputIterator>
SMatrix(InputIterator begin, InputIterator end);
/**
* Constructor with STL iterator interface. The data will be copied into the matrix
* In this case the value passed size must be equal to the matrix size (D1*D2)
*/
template<class InputIterator>
SMatrix(InputIterator begin, unsigned int size);
// skip this methods (they are too ambigous)
#ifdef OLD_IMPL
/// 2nd arg: set only diagonal?
SMatrix(const T& rhs, bool diagonal=false);
/// constructor via dyadic product
SMatrix(const SVector<T,D1>& rhs);
/// constructor via dyadic product
template <class A>
SMatrix(const Expr<A,T,D1>& rhs);
/** constructor via array, triag=true: array contains only upper/lower
triangular part of a symmetric matrix, len: length of array */
template <class T1>
//SMatrix(const T1* a, bool triang=false, unsigned int len=D1*D2);
// to avoid clash with TRootIOconst
SMatrix(const T1* a, bool triang, unsigned int len=D1*D2);
/// assign from a scalar value
SMatrix<T,D1,D2>& operator=(const T& rhs);
#endif
/**
construct a symmetric matrix from a SVector containing the upper(lower)
part of a triangular matrix
*/
#ifndef UNSUPPORTED_TEMPLATE_EXPRESSION
SMatrix(const SVector<T, D1*(D2+1)/2> & v, bool lower = false );
#else
template<unsigned int N>
SMatrix(const SVector<T,N> & v, bool lower = false );
#endif
///
template <class A>
SMatrix<T,D1,D2>& operator=(const Expr<A,T,D1,D2>& rhs);
#ifdef OLD_IMPL
/// return no. of matrix rows
static const unsigned int kRows = D1;
/// return no. of matrix columns
static const unsigned int kCols = D2;
/// return no of elements: rows*columns
static const unsigned int kSize = D1*D2;
#else
enum {
/// return no. of matrix rows
kRows = D1,
/// return no. of matrix columns
kCols = D2,
/// return no of elements: rows*columns
kSize = D1*D2
};
#endif
/** @name --- Access functions --- */
/// access the parse tree
T apply(unsigned int i) const;
/// return read-only pointer to internal array
const T* Array() const;
/// return pointer to internal array
T* Array();
// STL interface
/** STL iterator interface. */
iterator begin();
/** STL iterator interface. */
iterator end();
/** STL const_iterator interface. */
const_iterator begin() const;
/** STL const_iterator interface. */
const_iterator end() const;
/** @name --- Operators --- */
/// element wise comparison
bool operator==(const T& rhs) const;
/// element wise comparison
bool operator!=(const T& rhs) const;
/// element wise comparison
bool operator==(const SMatrix<T,D1,D2>& rhs) const;
/// element wise comparison
bool operator!=(const SMatrix<T,D1,D2>& rhs) const;
/// element wise comparison
template <class A>
bool operator==(const Expr<A,T,D1,D2>& rhs) const;
/// element wise comparison
template <class A>
bool operator!=(const Expr<A,T,D1,D2>& rhs) const;
/// element wise comparison
bool operator>(const T& rhs) const;
/// element wise comparison
bool operator<(const T& rhs) const;
/// element wise comparison
bool operator>(const SMatrix<T,D1,D2>& rhs) const;
/// element wise comparison
bool operator<(const SMatrix<T,D1,D2>& rhs) const;
/// element wise comparison
template <class A>
bool operator>(const Expr<A,T,D1,D2>& rhs) const;
/// element wise comparison
template <class A>
bool operator<(const Expr<A,T,D1,D2>& rhs) const;
/// read-only access
const T& operator()(unsigned int i, unsigned int j) const;
/// read/write access
T& operator()(unsigned int i, unsigned int j);
///
SMatrix<T,D1,D2>& operator+=(const SMatrix<T,D1,D2>& rhs);
SMatrix<T,D1,D2>& operator-=(const SMatrix<T,D1,D2>& rhs);
#ifdef OLD_IMPL
// this operations are not well defines
// in th eold impl they were implemented not as matrix - matrix multiplication, but as
// m(i,j)*m(i,j) multiplication
SMatrix<T,D1,D2>& operator*=(const SMatrix<T,D1,D2>& rhs);
SMatrix<T,D1,D2>& operator/=(const SMatrix<T,D1,D2>& rhs);
#endif
#ifndef __CINT__
///
template <class A>
SMatrix<T,D1,D2>& operator+=(const Expr<A,T,D1,D2>& rhs);
///
///
template <class A>
SMatrix<T,D1,D2>& operator-=(const Expr<A,T,D1,D2>& rhs);
///
///
#ifdef OLD_IMPL
template <class A>
SMatrix<T,D1,D2>& operator*=(const Expr<A,T,D1,D2>& rhs);
///
///
template <class A>
SMatrix<T,D1,D2>& operator/=(const Expr<A,T,D1,D2>& rhs);
#endif
#endif
/** @name --- Expert functions --- */
/**
invert symmetric, pos. def. Matrix via Dsinv.
This method change the current matrix
*/
bool Sinvert();
/**
invert symmetric, pos. def. Matrix via Dsinv.
This method returns a new matrix. In case the inversion fails
the current matrix is returned
*/
SMatrix<T,D1,D2> Sinverse() const;
/** determinant of symmetrc, pos. def. Matrix via Dsfact. \textbf{Note:} this
will destroy the contents of the Matrix!
*/
bool Sdet(T& det);
/** determinant of symmetrc, pos. def. Matrix via Dsfact. \textbf{Note:}
this method will preserve the contents of the Matrix!
*/
bool Sdet2(T& det) const;
/**
invert square Matrix via Dinv.
This method change the current matrix
*/
bool Invert();
/**
invert square Matrix via Dinv.
This method returns a new matrix. In case the inversion fails
the current matrix is returned
*/
SMatrix<T,D1,D2> Inverse() const;
/**
determinant of square Matrix via Dfact. \textbf{Note:} this will destroy
the contents of the Matrix!
*/
bool Det(T& det);
/**
determinant of square Matrix via Dfact. \textbf{Note:} this will preserve
the content of the Matrix!
*/
bool Det2(T& det) const;
/// place a vector in a Matrix row
template <unsigned int D>
SMatrix<T,D1,D2>& Place_in_row(const SVector<T,D>& rhs,
unsigned int row,
unsigned int col);
/// place a vector expression in a Matrix row
template <class A, unsigned int D>
SMatrix<T,D1,D2>& Place_in_row(const Expr<A,T,D>& rhs,
unsigned int row,
unsigned int col);
/// place a vector in a Matrix column
template <unsigned int D>
SMatrix<T,D1,D2>& Place_in_col(const SVector<T,D>& rhs,
unsigned int row,
unsigned int col);
/// place a vector expression in a Matrix column
template <class A, unsigned int D>
SMatrix<T,D1,D2>& Place_in_col(const Expr<A,T,D>& rhs,
unsigned int row,
unsigned int col);
/// place a matrix in this matrix
template <unsigned int D3, unsigned int D4>
SMatrix<T,D1,D2>& Place_at(const SMatrix<T,D3,D4>& rhs,
unsigned int row,
unsigned int col);
/// place a matrix expression in this matrix
template <class A, unsigned int D3, unsigned int D4>
SMatrix<T,D1,D2>& Place_at(const Expr<A,T,D3,D4>& rhs,
unsigned int row,
unsigned int col);
/**
return a full Matrix row as a vector (copy the content in a new vector)
*/
SVector<T,D2> Row(unsigned int therow) const;
/**
return a full Matrix column as a vector (copy the content in a new vector)
*/
SVector<T,D1> Col(unsigned int thecol) const;
/**
return a slice of therow as a vector starting at the colum value col0 until col0+N.
Condition col0+N <= D2
*/
template <unsigned int N>
SVector<T,N> SubRow(unsigned int therow, unsigned int col0 = 0 ) const;
/**
return a slice of the column as a vector starting at the row value row0 until row0+Dsub.
Condition row0+N <= D1
*/
template <unsigned int N>
SVector<T,N> SubCol(unsigned int thecol, unsigned int row0 = 0) const;
/**
return a submatrix with the upper left corner at the values (row0, col0) and with sizes N1, N2
Condition row0+N1 <= D1 && col0+N2 <=D2
*/
template <unsigned int N1, unsigned int N2 >
SMatrix<T,N1,N2> Sub(unsigned int row0, unsigned int col0) const;
/**
return diagonal elements of a matrix as a Vector.
It works only for squared matrices D1 == D2, otherwise it will produce a compile error
*/
SVector<T,D1> Diagonal() const;
/**
return the upper Triangular block of the matrices (including the diagonal) as
a vector of sizes N = D1 * (D1 + 1)/2.
It works only for square matrices with D1==D2, otherwise it will produce a compile error
*/
#ifndef UNSUPPORTED_TEMPLATE_EXPRESSION
SVector<T, D1 * (D2 +1)/2> UpperBlock() const;
#else
template<unsigned int N>
SVector<T,N> UpperBlock() const;
#endif
/**
return the lower Triangular block of the matrices (including the diagonal) as
a vector of sizes N = D1 * (D1 + 1)/2.
It works only for square matrices with D1==D2, otherwise it will produce a compile error
*/
#ifndef UNSUPPORTED_TEMPLATE_EXPRESSION
SVector<T, D1 * (D2 +1)/2> LowerBlock() const;
#else
template<unsigned int N>
SVector<T,N> LowerBlock() const;
#endif
// submatrices
/// Print: used by operator<<()
std::ostream& Print(std::ostream& os) const;
private:
T fArray[D1*D2];
}; // end of class SMatrix
//==============================================================================
// operator<<
//==============================================================================
template <class T, unsigned int D1, unsigned int D2>
inline std::ostream& operator<<(std::ostream& os, const ROOT::Math::SMatrix<T,D1,D2>& rhs) {
return rhs.Print(os);
}
} // namespace Math
} // namespace ROOT
#ifndef __CINT__
#include "Math/SMatrix.icc"
// include Matrix-Vector multiplication
#include "Math/MatrixFunctions.h"
#endif //__CINT__
#endif /* ROOT_Math_SMatrix */
ROOT page - Class index - Class Hierarchy - Top of the page
This page has been automatically generated. If you have any comments or suggestions about the page layout send a mail to ROOT support, or contact the developers with any questions or problems regarding ROOT.