FairRoot/PandaRoot
matrix.h
Go to the documentation of this file.
1 /*
2 Copyright 2011. All rights reserved.
3 Institute of Measurement and Control Systems
4 Karlsruhe Institute of Technology, Germany
5 
6 Authors: Andreas Geiger
7 
8 matrix is free software; you can redistribute it and/or modify it under the
9 terms of the GNU General Public License as published by the Free Software
10 Foundation; either version 2 of the License, or any later version.
11 
12 matrix is distributed in the hope that it will be useful, but WITHOUT ANY
13 WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
14 PARTICULAR PURPOSE. See the GNU General Public License for more details.
15 
16 You should have received a copy of the GNU General Public License along with
17 matrix; if not, write to the Free Software Foundation, Inc., 51 Franklin
18 Street, Fifth Floor, Boston, MA 02110-1301, USA
19  */
20 
21 #ifndef MATRIX_H
22 #define MATRIX_H
23 
24 #include <stdio.h>
25 #include <string.h>
26 #include <stdlib.h>
27 #include <iostream>
28 #include <vector>
29 
30 #ifndef _MSC_VER
31 //#include <stdint.h>
32 typedef int int32_t;
33 #else
34 typedef __int8 int8_t;
35 typedef __int16 int16_t;
36 typedef __int32 int32_t;
37 typedef __int64 int64_t;
38 typedef unsigned __int8 uint8_t;
39 typedef unsigned __int16 uint16_t;
40 typedef unsigned __int32 uint32_t;
41 typedef unsigned __int64 uint64_t;
42 #endif
43 
44 //#define endll endl << endl // double end line definition
45 
46 //typedef long double FLOAT; //quad precision
47 typedef double FLOAT; // double precision
48 //typedef float FLOAT; // single precision
49 
50 class Matrix {
51 
52 public:
53 
54  // constructor / deconstructor
55  Matrix (); // init empty 0x0 matrix
56  Matrix (const int32_t m,const int32_t n); // init empty mxn matrix
57  Matrix (const int32_t m,const int32_t n,const FLOAT* val_); // init mxn matrix with values from array 'val'
58  Matrix (const Matrix &M); // creates deepcopy of M
59  ~Matrix ();
60 
61  // assignment operator, copies contents of M
62  Matrix& operator= (const Matrix &M);
63 
64  // copies submatrix of M into array 'val', default values copy whole row/column/matrix
65  void getData(FLOAT* val_,int32_t i1=0,int32_t j1=0,int32_t i2=-1,int32_t j2=-1);
66 
67  // set or get submatrices of current matrix
68  Matrix getMat(int32_t i1,int32_t j1,int32_t i2=-1,int32_t j2=-1);
69  void setMat(const Matrix &M,const int32_t i,const int32_t j);
70 
71  //set mxn matrix with values from array 'val'
72  void setVal(const int32_t m,const int32_t n,const FLOAT* val_);
73 
74  // set sub-matrix to scalar (default 0), -1 as end replaces whole row/column/matrix
75  void setVal(FLOAT s,int32_t i1=0,int32_t j1=0,int32_t i2=-1,int32_t j2=-1);
76 
77  // set (part of) diagonal to scalar, -1 as end replaces whole diagonal
78  void setDiag(FLOAT s,int32_t i1=0,int32_t i2=-1);
79 
80  // clear matrix
81  void zero();
82 
83  // extract columns with given index
84  Matrix extractCols (std::vector<int> idx);
85 
86  // create identity matrix
87  static Matrix eye (const int32_t m);
88  void eye ();
89 
90  // create matrix with ones
91  static Matrix ones(const int32_t m,const int32_t n);
92 
93  // create diagonal matrix with nx1 or 1xn matrix M as elements
94  static Matrix diag(const Matrix &M);
95 
96  // returns the m-by-n matrix whose elements are taken column-wise from M
97  static Matrix reshape(const Matrix &M,int32_t m,int32_t n);
98 
99  // create 3x3 rotation matrices (convention: http://en.wikipedia.org/wiki/Rotation_matrix)
100  static Matrix rotMatX(const FLOAT &angle);
101  static Matrix rotMatY(const FLOAT &angle);
102  static Matrix rotMatZ(const FLOAT &angle);
103 
104  //homogenize 3x3 matrix or 3 vector to 4x4 matrix or 4 vector
105  static Matrix homogenize(const Matrix &M);
106  //de-homogenize a 4 vector
107  static Matrix dehomogenize(const Matrix &M);
108  //homogenize 3 vector to 4x4 translation matrix
109  static Matrix homogenizeTranslation(const Matrix &M);
110  //homogenize rotation and translation simultaniously
111  static Matrix homogenizeRotTrans(const Matrix &R, const Matrix &t);
112 
113  // simple arithmetic operations
114  Matrix operator+ (const Matrix &M); // add matrix
115  Matrix operator- (const Matrix &M); // subtract matrix
116  Matrix operator* (const Matrix &M); // multiply with matrix
117  Matrix operator* (const FLOAT &s); // multiply with scalar
118  Matrix operator/ (const Matrix &M); // divide elementwise by matrix (or vector)
119  Matrix operator/ (const FLOAT &s); // divide by scalar
120  Matrix operator- (); // negative matrix
121  Matrix operator~ (); // transpose
122  FLOAT l2norm (); // euclidean norm (vectors) / frobenius norm (matrices)
123  FLOAT mean (); // mean of all elements in matrix
124 
125  // complex arithmetic operations
126  static Matrix cross (const Matrix &a, const Matrix &b); // cross product of two vectors
127  static Matrix inv (const Matrix &M); // invert matrix M
128  bool inv (); // invert this matrix
129  FLOAT det (); // returns determinant of matrix
130  bool solve (const Matrix &M,FLOAT eps=1e-20); // solve linear system M*x=B, replaces *this and M
131  bool lu(int32_t *idx, FLOAT &d, FLOAT eps=1e-20); // replace *this by lower upper decomposition
132  void svd(Matrix &U,Matrix &W,Matrix &V); // singular value decomposition *this = U*diag(W)*V^T
133 
134  // print matrix to stream
135  friend std::ostream& operator<< (std::ostream& out,const Matrix& M);
136 
137  // direct data access
140 
141 private:
142 
143  void allocateMemory (const int32_t m_,const int32_t n_);
144  void releaseMemory ();
145  inline FLOAT pythag(FLOAT a,FLOAT b);
146 
147 };
148 
149 #endif // MATRIX_H
bool solve(const Matrix &M, FLOAT eps=1e-20)
void setMat(const Matrix &M, const int32_t i, const int32_t j)
void zero()
Matrix operator-()
static Matrix ones(const int32_t m, const int32_t n)
bool lu(int32_t *idx, FLOAT &d, FLOAT eps=1e-20)
TObjArray * d
static Matrix diag(const Matrix &M)
Int_t i
Definition: run_full.C:25
TTree * b
FLOAT det()
Matrix operator~()
TLorentzVector s
Definition: Pnd2DStar.C:50
static Matrix homogenizeTranslation(const Matrix &M)
static Matrix homogenize(const Matrix &M)
Matrix operator*(const Matrix &M)
double FLOAT
Definition: matrix.h:47
Matrix operator/(const Matrix &M)
void allocateMemory(const int32_t m_, const int32_t n_)
Matrix & operator=(const Matrix &M)
static Matrix reshape(const Matrix &M, int32_t m, int32_t n)
void setVal(const int32_t m, const int32_t n, const FLOAT *val_)
int int32_t
Definition: matrix.h:32
int idx[MAX]
Definition: autocutx.C:38
void setDiag(FLOAT s, int32_t i1=0, int32_t i2=-1)
Int_t a
Definition: anaLmdDigi.C:126
static Matrix cross(const Matrix &a, const Matrix &b)
#define W
Definition: createSTT.C:76
FLOAT pythag(FLOAT a, FLOAT b)
FLOAT mean()
double eps(TVector3 v1, TVector3 v2)
static Matrix rotMatX(const FLOAT &angle)
int32_t m
Definition: matrix.h:139
TFile * out
Definition: reco_muo.C:20
Definition: matrix.h:50
bool inv()
int32_t n
Definition: matrix.h:139
static Matrix homogenizeRotTrans(const Matrix &R, const Matrix &t)
void releaseMemory()
void getData(FLOAT *val_, int32_t i1=0, int32_t j1=0, int32_t i2=-1, int32_t j2=-1)
friend std::ostream & operator<<(std::ostream &out, const Matrix &M)
TTree * t
Definition: bump_analys.C:13
Matrix operator+(const Matrix &M)
static Matrix rotMatZ(const FLOAT &angle)
static Matrix dehomogenize(const Matrix &M)
Double_t angle
void svd(Matrix &U, Matrix &W, Matrix &V)
FLOAT l2norm()
Matrix extractCols(std::vector< int > idx)
Double_t R
Definition: checkhelixhit.C:61
void eye()
Matrix getMat(int32_t i1, int32_t j1, int32_t i2=-1, int32_t j2=-1)
FLOAT ** val
Definition: matrix.h:138
static Matrix rotMatY(const FLOAT &angle)