LinearAlgebra.H File Reference

---

Detailed Description

Definition in file LinearAlgebra.H.

#include "Image/LinearAlgebraFlags.H"
#include "Util/Promotions.H"

Include dependency graph for LinearAlgebra.H:

This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Functions
void	svd (const Image< double > &M, Image< double > &U, Image< double > &S, Image< double > &V, SvdAlgo algo, const SvdFlag flags=0)
	Compute the singular decomposition of M: M = U * S * V^T.
void	svdf (const Image< float > &M, Image< float > &U, Image< float > &S, Image< float > &V, SvdAlgo algo, const SvdFlag flags=0)
	Single-precision counterpart of svd().
Image< double >	svdPseudoInv (const Image< double > &x, const SvdAlgo algo, int *const rank=0, const double tolfactor=1.0e-8)
	Compute the pseudo-inverse using the method described in svd().
Image< float >	svdPseudoInvf (const Image< float > &x, const SvdAlgo algo, int *const rank=0, const float tolfactor=1.0e-8f)
	Single-precision counterpart of svdPseudoInv().
template<class T>
Image< typename promote_trait< T, float >::TP >	naiveUnstablePseudoInv (const Image< T > &M)
	NOTE that this procedure is called "naive, unstable" for a reason!

---

Function Documentation

template<class T> Image<typename promote_trait<T, float>::TP> naiveUnstablePseudoInv (  const Image< T > &  M  )

NOTE that this procedure is called "naive, unstable" for a reason! This procedure computes the pseudo-inverse using the direct formula pinv(M) = (M^T * M)^-1 * M^T. The formula is analytically correct, but in practice it can be numerically unstable. Indeed it might work just fine for your matrix, but to be safe you'll need to check the identity pinv(M)*M ~= I -- for example, do RMSerr(matrixMult(M,pinv(M)), eye<float>(M.getHeight())) and make sure that it's small (say, less than 0.01 or so). If you do get a large RMS error, then you can consider using svdPseudoInv() instead, which offers more numerical stability. Definition at line 217 of file LinearAlgebra.C. References GVX_TRACE, matrixInv(), matrixMult(), and transpose(). Referenced by Image_xx_svd_gsl_xx_1(), Image_xx_svd_lapack_xx_1(), and svdIdent().

void svd (  const Image< double > &  M, Image< double > &  U, Image< double > &  S, Image< double > &  V, SvdAlgo  algo, const SvdFlag  flags = 0 )

Compute the singular decomposition of M: M = U * S * V^T. The U, S, V matrices have special properties: U is column-orthogonal, such that (U^T * U) = I S is square, diagonal, so that S is easily invertible V is square and orthogonal, so (V^T * V) = (V * V^T) = I Note that these properties mean that the svd can be used to compute the pseudo-inverse of a general matrix M: M = U*S*V^T U^T*M = S*V^T (since U is column-orthogonal) S^-1*U^T*M = V^T (since S is diagonal and invertible) V*S^-1*U^T*M = I (since V is square and orthogonal) by definition, pinv(M)*M = I, therefore V*S^-1*U^T = pinv(M) (Note ^T is used to designate the matrix transpose, and ^-1 is used to designate the matrix inverse). You can use the routine svdPseudoInv to compute the pseudo-inverse using just this method.

void svdf (  const Image< float > &  M, Image< float > &  U, Image< float > &  S, Image< float > &  V, SvdAlgo  algo, const SvdFlag  flags = 0 )

Single-precision counterpart of svd().

Image<double> svdPseudoInv (  const Image< double > &  x, const SvdAlgo  algo, int *const   rank = 0, const double  tolfactor = 1.0e-8 )

Compute the pseudo-inverse using the method described in svd(). Optionally throw out small values from the diagonal S matrix -- this helps numerical stability. The number of eigenvalues retained is controlled by tolfactor. If tolfactor is less than 1.0, then it serves as a multiplier of the largest eigenvalue, below which any smaller eigenvalues are thrown out. If tolfactor is greater than 1.0, then the integer portion of tolfactor specifies the number of eigenvalues to be retained. If rank is a non-null pointer, then the rank of the diagonal matrix (after setting small values to zero) will be returned through that pointer.

Image<float> svdPseudoInvf (  const Image< float > &  x, const SvdAlgo  algo, int *const   rank = 0, const float  tolfactor = 1.0e-8f )

Single-precision counterpart of svdPseudoInv().

---