iRoCS Toolbox  1.1.0
ATBLinAlg.hh
Go to the documentation of this file.
1 /**************************************************************************
2  *
3  * Copyright (C) 2015 Thorsten Falk
4  *
5  * Image Analysis Lab, University of Freiburg, Germany
6  *
7  * This program is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation; either version 3 of the License, or
10  * (at your option) any later version.
11  *
12  * This program is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program; if not, write to the Free Software Foundation,
19  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  *
21  **************************************************************************/
22 
23 #ifndef ATBLINALG_HH
24 #define ATBLINALG_HH
25 
26 #if defined HAVE_CONFIG_H
27 #include <config.hh>
28 #endif
29 
30 #include <gsl/gsl_linalg.h>
31 #include <gsl/gsl_eigen.h>
32 #include <gsl/gsl_permute_vector.h>
33 #include <gsl/gsl_sort_vector.h>
34 #include <gsl/gsl_sort.h>
35 #include <gsl/gsl_permutation.h>
36 
37 #include "TypeTraits.hh"
38 #include "ATBGSLWrapper.hh"
39 #include "Quaternion.hh"
40 
41 #ifdef HAVE_BLITZ_V9
42 #include <blitz/tinyvec-et.h>
43 #endif
44 
45 namespace atb
46 {
47 
49 
50 /*======================================================================*/
63 /*======================================================================*/
64  double factorialFraction(int n, int k);
65 
66 /*======================================================================*/
81 /*======================================================================*/
82  ptrdiff_t modulo(ptrdiff_t a, ptrdiff_t b);
83 
84 /*======================================================================*/
100 /*======================================================================*/
101  template<int Dim>
102  blitz::TinyVector<ptrdiff_t,Dim> modulo(
103  blitz::TinyVector<ptrdiff_t,Dim> const &a,
104  blitz::TinyVector<ptrdiff_t,Dim> const &b);
105 
106 /*======================================================================*/
114 /*======================================================================*/
115  template<typename DataT,int Dim>
116  blitz::TinyVector<DataT,Dim+1>
117  euclideanToHomogeneous(blitz::TinyVector<DataT,Dim> const &pos);
118 
119 /*======================================================================*/
127 /*======================================================================*/
128  template<typename DataT,int Dim>
129  blitz::TinyVector<DataT,Dim-1>
130  homogeneousToEuclidean(blitz::TinyVector<DataT,Dim> const &pos);
131 
132 /*======================================================================*/
148 /*======================================================================*/
149  blitz::TinyVector<double,3>
150  euclideanToSpherical(blitz::TinyVector<double,3> const &eucl,
151  blitz::Array<double,2> const &trafo);
152 
153 /*======================================================================*/
164 /*======================================================================*/
165  blitz::TinyVector<double,3>
166  euclideanToSpherical(blitz::TinyVector<double,3> const &eucl);
167 
168 /*======================================================================*/
185 /*======================================================================*/
186  blitz::TinyVector<double,3>
187  sphericalToEuclidean(blitz::TinyVector<double,3> const &spherical,
188  blitz::Array<double,2> const &trafo);
189 
190 /*======================================================================*/
201 /*======================================================================*/
202  blitz::TinyVector<double,3>
203  sphericalToEuclidean(blitz::TinyVector<double,3> const &spherical);
204 
205 /*======================================================================*/
216 /*======================================================================*/
217  template<typename DataT>
218  blitz::TinyVector<DataT,3>
219  rotate(blitz::TinyVector<DataT,3> const &vec,
220  blitz::TinyVector<DataT,3> const &axis,
221  DataT angle);
222 
223 /*======================================================================*/
233 /*======================================================================*/
234  template<typename DataT>
235  double dot(blitz::Array<DataT,1> const &vec1,
236  blitz::Array<DataT,1> const &vec2);
237 
238 /*======================================================================*/
247 /*======================================================================*/
248  template<typename DataT, int Dim1, int Dim2>
249  blitz::TinyMatrix<DataT,Dim1,Dim2>
250  outerProduct(blitz::TinyVector<DataT,Dim1> const &vec1,
251  blitz::TinyVector<DataT,Dim2> const &vec2);
252 
253 /*======================================================================*/
261 /*======================================================================*/
262  template<typename DataT, int Dim1, int Dim2>
263  blitz::TinyMatrix<DataT,Dim2,Dim1>
264  transpose(blitz::TinyMatrix<DataT,Dim1,Dim2> const &in);
265 
266 /*======================================================================*/
276 /*======================================================================*/
277  template<typename DataT>
278  blitz::Array<DataT,2> transpose(blitz::Array<DataT,2> const &in);
279 
280 /*======================================================================*/
291 /*======================================================================*/
292  template<typename MatrixT, typename BaseT, int Dim>
293  blitz::TinyVector<BaseT,Dim>
294  mvMult(MatrixT const &m, blitz::TinyVector<BaseT,Dim> const &v);
295 
296 /*======================================================================*/
305 /*======================================================================*/
306  template<typename BaseT, int Dim>
307  blitz::TinyVector<BaseT,Dim>
308  mvMult(blitz::Array<BaseT,2> const &m, blitz::TinyVector<BaseT,Dim> const &v);
309 
310 /*======================================================================*/
321 /*======================================================================*/
322  template<typename BaseT>
323  blitz::Array<BaseT,1>
324  mvMult(blitz::Array<BaseT,2> const &m, blitz::Array<BaseT,1> const &v);
325 
326 /*======================================================================*/
337 /*======================================================================*/
338  template<typename DataT>
339  blitz::Array<DataT,2>
340  mmMult(blitz::Array<DataT,2> const &A, blitz::Array<DataT,2> const &B);
341 
342 /*======================================================================*/
351 /*======================================================================*/
352  template<typename DataT, int Dim1, int Dim2, int Dim3>
353  blitz::TinyMatrix<DataT,Dim1,Dim3>
354  mmMult(blitz::TinyMatrix<DataT,Dim1,Dim2> const &A,
355  blitz::TinyMatrix<DataT,Dim2,Dim3> const &B);
356 
357 /*======================================================================*/
364 /*======================================================================*/
365  void invert(blitz::Array<double,2> const &A, blitz::Array<double,2> &Ainv);
366 
367 /*======================================================================*/
376 /*======================================================================*/
377  blitz::Array<double,2>
378  invert(blitz::Array<double,2> const &A);
379 
380 /*======================================================================*/
387 /*======================================================================*/
388  template<int Dim>
389  blitz::TinyMatrix<double,Dim,Dim>
390  invert(blitz::TinyMatrix<double,Dim,Dim> const &A);
391 
392 /*======================================================================*/
400 /*======================================================================*/
401  template<int Dim>
402  double determinant(blitz::TinyMatrix<double,Dim,Dim> const &A);
403 
404 /*======================================================================*/
414 /*======================================================================*/
415  template<int Dim>
417  blitz::TinyMatrix<double,Dim,Dim> const &A,
418  blitz::TinyVector<double,Dim> &lambda,
419  SortingType sort = Descending);
420 
421 /*======================================================================*/
431 /*======================================================================*/
432  template<int Dim>
433  void computeEigenvectors(
434  blitz::TinyMatrix<double,Dim,Dim> const &A,
435  blitz::TinyMatrix<double,Dim,Dim> &U,
436  blitz::TinyVector<double,Dim> const &lambda);
437 
438 /*======================================================================*/
450 /*======================================================================*/
451  template<int Dim>
453  blitz::TinyMatrix<double,Dim,Dim> const &A,
454  blitz::TinyMatrix<double,Dim,Dim> &U,
455  blitz::TinyVector<double,Dim> &lambda,
456  SortingType sort = Descending);
457 
458 /*======================================================================*/
466 /*======================================================================*/
467  template<typename DataT>
468  void sortElements(blitz::TinyVector<DataT,3> &v, SortingType sort);
469 
470 /*======================================================================*/
482 /*======================================================================*/
484  blitz::Array<double,2> const &A,
485  blitz::Array<double,2> &U, blitz::Array<double,1> &lambda,
486  SortingType sort = Descending);
487 
488 }
489 
490 #include "ATBLinAlg.icc"
491 
492 #endif
ptrdiff_t modulo(ptrdiff_t a, ptrdiff_t b)
Strictly positive modulo, to compute valid indices for Array access modulo Array size.
blitz::TinyVector< DataT, 3 > rotate(blitz::TinyVector< DataT, 3 > const &vec, blitz::TinyVector< DataT, 3 > const &axis, DataT angle)
Rotate the given 3-D vector by a specified angle angle around the given axis.
void eigenvalueDecompositionRealSymmetric(blitz::TinyMatrix< double, Dim, Dim > const &A, blitz::TinyVector< double, Dim > &lambda, SortingType sort=Descending)
Computation of the eigen value decomposition of real symmetric matrices.
blitz::TinyVector< double, 3 > euclideanToSpherical(blitz::TinyVector< double, 3 > const &eucl, blitz::Array< double, 2 > const &trafo)
Get the spherical (theta, phi, r) coordinates of the given euclidean point.
double determinant(blitz::TinyMatrix< double, Dim, Dim > const &A)
Compute the determinant of matrix A.
Query specific information about different data types.
ATBGSLWrapper.hh provides wrappers and setup routines to use GSL functions with ArrayToolbox data str...
SortingType
Definition: ATBLinAlg.hh:48
blitz::TinyVector< double, 3 > sphericalToEuclidean(blitz::TinyVector< double, 3 > const &spherical, blitz::Array< double, 2 > const &trafo)
Get the euclidean (z, y, x) coordinates of the given point on the 2-sphere.
void invert(blitz::Array< double, 2 > const &A, blitz::Array< double, 2 > &Ainv)
Compute the inverse of matrix A.
blitz::TinyMatrix< DataT, Dim1, Dim2 > outerProduct(blitz::TinyVector< DataT, Dim1 > const &vec1, blitz::TinyVector< DataT, Dim2 > const &vec2)
Calculation of the outer product of two vectors.
blitz::TinyVector< BaseT, Dim > mvMult(MatrixT const &m, blitz::TinyVector< BaseT, Dim > const &v)
Quadratic Matrix * Vector product (generic).
blitz::TinyMatrix< DataT, Dim2, Dim1 > transpose(blitz::TinyMatrix< DataT, Dim1, Dim2 > const &in)
Compute the transpose of the given Matrix.
double factorialFraction(int n, int k)
Compute the value of expression .
blitz::Array< DataT, 2 > mmMult(blitz::Array< DataT, 2 > const &A, blitz::Array< DataT, 2 > const &B)
Naive matrix multiplication A * B = C.
blitz::TinyVector< DataT, Dim+1 > euclideanToHomogeneous(blitz::TinyVector< DataT, Dim > const &pos)
Convert the given vector to homogeneous coordinates.
blitz::TinyVector< DataT, Dim-1 > homogeneousToEuclidean(blitz::TinyVector< DataT, Dim > const &pos)
Convert the given homogeneous vector to euclidean coordinates.
double dot(blitz::Array< DataT, 1 > const &vec1, blitz::Array< DataT, 1 > const &vec2)
Calculation of the inner product (dot product) of the given vectors.
void computeEigenvectors(blitz::TinyMatrix< double, Dim, Dim > const &A, blitz::TinyMatrix< double, Dim, Dim > &U, blitz::TinyVector< double, Dim > const &lambda)
Computates the eigenvectors of the real symmetric matrix A given the eigenvalues. ...
void sortElements(blitz::TinyVector< DataT, 3 > &v, SortingType sort)
Sort the values of the given TinyVector according to the given SortingType.