My Project
gsl_matrix.hh
Go to the documentation of this file.
1 /* -*- mia-c++ -*-
2  *
3  * This file is part of MIA - a toolbox for medical image analysis
4  * Copyright (c) Leipzig, Madrid 1999-2017 Gert Wollny
5  *
6  * MIA is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 3 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with MIA; if not, see <http://www.gnu.org/licenses/>.
18  *
19  */
20 
21 #ifndef GSLPP_MATRIX_HH
22 #define GSLPP_MATRIX_HH
23 
24 #include <cassert>
25 #include <iterator>
26 #include <gsl/gsl_matrix.h>
27 #include <mia/core/gsl_vector.hh>
28 
29 namespace gsl
30 {
31 
33 {
34  friend class const_matrix_iterator;
35 public:
36  matrix_iterator(gsl_matrix *m, bool begin):
37  m_matrix(m),
38  m_current(begin ? m->data : m->data + m->size1 * m->tda),
39  m_current_column(0),
40  m_row_jump(m->tda - m->size2)
41  {
42  }
43 
45  m_matrix(nullptr),
46  m_current(nullptr),
47  m_current_column(0)
48  {
49  }
50 
51 
52  matrix_iterator(const matrix_iterator& other) = default;
53 
54  double& operator *()
55  {
56  assert(m_current);
57  return *m_current;
58  }
59 
60  matrix_iterator& operator ++()
61  {
62  assert(m_current);
63  ++m_current;
64  ++m_current_column;
65 
66  if (m_current_column == m_matrix->size2) {
67  m_current += m_row_jump;
68  m_current_column = 0;
69  }
70 
71  return *this;
72  }
73 
74 
75  matrix_iterator operator ++(int)
76  {
77  matrix_iterator result(*this);
78  ++(*this);
79  return result;
80  }
81 
82  friend bool operator == (const matrix_iterator& lhs, const matrix_iterator& rhs);
83  friend bool operator != (const matrix_iterator& lhs, const matrix_iterator& rhs);
84 private:
85  gsl_matrix *m_matrix;
86  double *m_current;
87  size_t m_current_column;
88  size_t m_row_jump;
89 };
90 
91 bool EXPORT_GSL operator == (const matrix_iterator& lhs, const matrix_iterator& rhs);
92 bool EXPORT_GSL operator != (const matrix_iterator& lhs, const matrix_iterator& rhs);
93 
94 
96 {
97 public:
98  const_matrix_iterator(const gsl_matrix *m, bool begin):
99  m_matrix(m),
100  m_current(begin ? m->data : m->data + m->size1 * m->tda),
101  m_current_column(0),
102  m_row_jump(m->tda - m->size2)
103  {
104  }
105 
107  m_matrix(other.m_matrix),
108  m_current(other.m_current),
109  m_current_column(other.m_current_column),
110  m_row_jump(other.m_row_jump)
111  {
112  }
113 
115  m_matrix(nullptr),
116  m_current(nullptr),
117  m_current_column(0)
118  {
119  }
120 
121 
122  const_matrix_iterator(const const_matrix_iterator& other) = default;
123 
124  double operator *() const
125  {
126  assert(m_current);
127  return *m_current;
128  }
129 
130  const_matrix_iterator& operator ++()
131  {
132  assert(m_current);
133  ++m_current;
134  ++m_current_column;
135 
136  if (m_current_column == m_matrix->size2) {
137  m_current += m_row_jump;
138  m_current_column = 0;
139  }
140 
141  return *this;
142  }
143 
144  const_matrix_iterator operator ++(int)
145  {
146  const_matrix_iterator result(*this);
147  ++(*this);
148  return result;
149  }
150 
151  friend bool operator == (const const_matrix_iterator& lhs, const const_matrix_iterator& rhs);
152  friend bool operator != (const const_matrix_iterator& lhs, const const_matrix_iterator& rhs);
153 private:
154  const gsl_matrix *m_matrix;
155  const double *m_current;
156  size_t m_current_column;
157  size_t m_row_jump;
158 };
159 
160 bool EXPORT_GSL operator == (const const_matrix_iterator& lhs, const const_matrix_iterator& rhs);
161 bool EXPORT_GSL operator != (const const_matrix_iterator& lhs, const const_matrix_iterator& rhs);
162 
163 
164 
165 
166 
172 {
173 public:
176 
177  Matrix();
178 
185  Matrix(size_t rows, size_t columns, bool clean);
186 
193  Matrix(size_t rows, size_t columns, double init);
194 
202  Matrix(size_t rows, size_t columns, const double *init);
203 
207  Matrix(const Matrix& other);
208 
209 
213  Matrix(Matrix&& other);
214 
215 
220  Matrix(gsl_matrix *m);
221 
226  Matrix(const gsl_matrix *m);
227 
232  Matrix transposed() const;
233 
234 
238  Matrix& operator =(const Matrix& other);
239 
243  Matrix& operator =(Matrix&& other);
244 
245 
252  void reset(size_t rows, size_t columns, bool clean);
253 
260  void reset(size_t rows, size_t columns, double init);
261 
262  ~Matrix();
263 
264 
268  size_t rows()const;
269 
273  size_t cols()const;
274 
281  void set(size_t i, size_t j, double x);
282 
290  double operator ()(size_t i, size_t j) const;
291 
295  operator gsl_matrix *();
296 
300  operator const gsl_matrix *() const;
301 
306  Matrix column_covariance() const;
307 
312  Matrix row_covariance() const;
313 
319  matrix_iterator begin();
320 
326  matrix_iterator end();
327 
333  const_matrix_iterator begin() const;
334 
340  const_matrix_iterator end() const;
341 
348  void set_row(int r, const Vector& row);
349 
355  VectorView get_row(int r);
356 
362  ConstVectorView get_row(int r) const;
363 
370  void set_column(int c, const Vector& col);
371 
377  VectorView get_column(int c);
378 
384  ConstVectorView get_column(int c) const;
385 
393  double dot_row(int r, const Vector& v) const;
394 
403  double dot_column(int c, const Vector& col) const;
404 
409  void print(std::ostream& os) const;
410 
416  Matrix& operator -=(const Matrix& rhs)
417  {
418  gsl_matrix_sub(*this, rhs);
419  return *this;
420  }
421 
427  Matrix& operator +=(const Matrix& rhs)
428  {
429  gsl_matrix_add(*this, rhs);
430  return *this;
431  }
432 
438  Matrix& operator *=(double rhs)
439  {
440  gsl_matrix_scale(*this, rhs);
441  return *this;
442  }
443 
444  bool is_valid() const;
445 
446  bool is_writable() const;
447 
448 private:
449  gsl_matrix *m_matrix;
450  const gsl_matrix *m_const_matrix;
451  bool m_owner;
452 };
453 
454 
455 Matrix EXPORT_GSL operator * (const Matrix& lhs, const Matrix& rhs);
456 Matrix EXPORT_GSL operator + (const Matrix& lhs, const Matrix& rhs);
457 Matrix EXPORT_GSL operator - (const Matrix& lhs, const Matrix& rhs);
458 
459 inline std::ostream& operator << (std::ostream& os, const Matrix& m)
460 {
461  m.print(os);
462  return os;
463 }
464 
470 
473 };
474 
480 
481 
482 } // end namespace
483 
484 namespace std
485 {
486 
487 template <>
488 class iterator_traits< gsl::const_matrix_iterator >
489 {
490 public:
491  typedef size_t difference_type;
492  typedef double value_type;
493  typedef const double *pointer;
494  typedef const double& reference;
495  typedef forward_iterator_tag iterator_category;
496 };
497 
498 template <>
499 class iterator_traits< gsl::matrix_iterator >
500 {
501 public:
502  typedef size_t difference_type;
503  typedef double value_type;
504  typedef double *pointer;
505  typedef double& reference;
506  typedef forward_iterator_tag iterator_category;
507 };
508 
509 }
510 
511 
512 #endif
std::iterator_traits< gsl::const_matrix_iterator >::value_type
double value_type
Definition: gsl_matrix.hh:492
gsl::operator==
bool EXPORT_GSL operator==(const matrix_iterator &lhs, const matrix_iterator &rhs)
gsl::ConstVectorView
Definition: gsl_vector.hh:243
gsl::const_matrix_iterator::const_matrix_iterator
const_matrix_iterator()
Definition: gsl_matrix.hh:114
gsl::Matrix::iterator
matrix_iterator iterator
Definition: gsl_matrix.hh:174
gsl::operator-
vector_iterator operator-(const vector_iterator &it, int dist)
Definition: gsl_iterator.hh:155
gsl::matrix_iterator
Definition: gsl_matrix.hh:32
gsl::operator+
vector_iterator operator+(const vector_iterator &it, int dist)
Definition: gsl_iterator.hh:148
gsl::operator!=
bool EXPORT_GSL operator!=(const matrix_iterator &lhs, const matrix_iterator &rhs)
operator+=
EXPORT_2D C2DFVectorfield & operator+=(C2DFVectorfield &a, const C2DFVectorfield &b)
std::iterator_traits< gsl::matrix_iterator >::iterator_category
forward_iterator_tag iterator_category
Definition: gsl_matrix.hh:506
gsl::Matrix::print
void print(std::ostream &os) const
gsl::VectorView
Definition: gsl_vector.hh:227
gsl::CSymmvEvalEvec::eval
Vector eval
Definition: gsl_matrix.hh:472
gsl::Vector
Definition: gsl_vector.hh:38
std::iterator_traits< gsl::matrix_iterator >::difference_type
size_t difference_type
Definition: gsl_matrix.hh:502
gsl::Matrix
Definition: gsl_matrix.hh:171
gsl::CSymmvEvalEvec
Definition: gsl_matrix.hh:468
std::iterator_traits< gsl::matrix_iterator >::value_type
double value_type
Definition: gsl_matrix.hh:503
gsl::CSymmvEvalEvec::evec
Matrix evec
Definition: gsl_matrix.hh:471
std::iterator_traits< gsl::const_matrix_iterator >::iterator_category
forward_iterator_tag iterator_category
Definition: gsl_matrix.hh:495
std::iterator_traits< gsl::const_matrix_iterator >::difference_type
size_t difference_type
Definition: gsl_matrix.hh:491
std
Definition: gsl_iterator.hh:323
std::iterator_traits< gsl::matrix_iterator >::pointer
double * pointer
Definition: gsl_matrix.hh:504
std::iterator_traits< gsl::matrix_iterator >::reference
double & reference
Definition: gsl_matrix.hh:505
gsl::matrix_iterator::matrix_iterator
matrix_iterator(gsl_matrix *m, bool begin)
Definition: gsl_matrix.hh:36
gsl::matrix_iterator::matrix_iterator
matrix_iterator()
Definition: gsl_matrix.hh:44
gsl::Matrix::const_iterator
const_matrix_iterator const_iterator
Definition: gsl_matrix.hh:175
EXPORT_GSL
#define EXPORT_GSL
Definition: gsl_defines.hh:38
gsl_vector.hh
gsl::const_matrix_iterator::const_matrix_iterator
const_matrix_iterator(const gsl_matrix *m, bool begin)
Definition: gsl_matrix.hh:98
gsl::const_matrix_iterator::const_matrix_iterator
const_matrix_iterator(const matrix_iterator &other)
Definition: gsl_matrix.hh:106
gsl::const_matrix_iterator
Definition: gsl_matrix.hh:95
gsl
Definition: gsl_iterator.hh:28
gsl::operator*
Matrix EXPORT_GSL operator*(const Matrix &lhs, const Matrix &rhs)
gsl::operator<<
std::ostream & operator<<(std::ostream &os, const Matrix &m)
Definition: gsl_matrix.hh:459
std::iterator_traits< gsl::const_matrix_iterator >::pointer
const typedef double * pointer
Definition: gsl_matrix.hh:493
std::iterator_traits< gsl::const_matrix_iterator >::reference
const typedef double & reference
Definition: gsl_matrix.hh:494
gsl::matrix_inv_sqrt
void EXPORT_GSL matrix_inv_sqrt(Matrix &m)