Commit 21c0c84a authored by Manuel Günther's avatar Manuel Günther

Removed the last occurrences of ..._ and ..._nocheck functions

parent fe06cddc
Pipeline #10434 passed with stages
in 7 minutes and 19 seconds
...@@ -57,16 +57,6 @@ void bob::math::eig(const blitz::Array<double,2>& A, ...@@ -57,16 +57,6 @@ void bob::math::eig(const blitz::Array<double,2>& A,
bob::core::array::assertSameShape(V,shape2); bob::core::array::assertSameShape(V,shape2);
bob::core::array::assertSameShape(D,shape1); bob::core::array::assertSameShape(D,shape1);
bob::math::eig_(A, V, D);
}
void bob::math::eig_(const blitz::Array<double,2>& A,
blitz::Array<std::complex<double>,2>& V,
blitz::Array<std::complex<double>,1>& D)
{
// Size variable
const int N = A.extent(0);
// Prepares to call LAPACK function // Prepares to call LAPACK function
// Initialises LAPACK variables // Initialises LAPACK variables
const char jobvl = 'N'; // Do NOT compute left eigen-vectors const char jobvl = 'N'; // Do NOT compute left eigen-vectors
...@@ -147,15 +137,6 @@ void bob::math::eigSym(const blitz::Array<double,2>& A, ...@@ -147,15 +137,6 @@ void bob::math::eigSym(const blitz::Array<double,2>& A,
bob::core::array::assertSameShape(V,shape2); bob::core::array::assertSameShape(V,shape2);
bob::core::array::assertSameShape(D,shape1); bob::core::array::assertSameShape(D,shape1);
bob::math::eigSym_(A, V, D);
}
void bob::math::eigSym_(const blitz::Array<double,2>& A,
blitz::Array<double,2>& V, blitz::Array<double,1>& D)
{
// Size variable
const int N = A.extent(0);
// Prepares to call LAPACK function // Prepares to call LAPACK function
// Initialises LAPACK variables // Initialises LAPACK variables
const char jobz = 'V'; // Get both the eigenvalues and the eigenvectors const char jobz = 'V'; // Get both the eigenvalues and the eigenvectors
...@@ -234,15 +215,6 @@ void bob::math::eigSym(const blitz::Array<double,2>& A, const blitz::Array<doubl ...@@ -234,15 +215,6 @@ void bob::math::eigSym(const blitz::Array<double,2>& A, const blitz::Array<doubl
bob::core::array::assertSameShape(V,shape2); bob::core::array::assertSameShape(V,shape2);
bob::core::array::assertSameShape(D,shape1); bob::core::array::assertSameShape(D,shape1);
bob::math::eigSym_(A, B, V, D);
}
void bob::math::eigSym_(const blitz::Array<double,2>& A, const blitz::Array<double,2>& B,
blitz::Array<double,2>& V, blitz::Array<double,1>& D)
{
// Size variable
const int N = A.extent(0);
// Prepares to call LAPACK function // Prepares to call LAPACK function
// Initialises LAPACK variables // Initialises LAPACK variables
const int itype = 1; const int itype = 1;
......
...@@ -18,6 +18,7 @@ void bob::math::pinv(const blitz::Array<double,2>& A, ...@@ -18,6 +18,7 @@ void bob::math::pinv(const blitz::Array<double,2>& A,
// Size variables // Size variables
const int M = A.extent(0); const int M = A.extent(0);
const int N = A.extent(1); const int N = A.extent(1);
const int nb_singular = std::min(M,N);
// Checks zero base // Checks zero base
bob::core::array::assertZeroBase(A); bob::core::array::assertZeroBase(A);
...@@ -26,18 +27,6 @@ void bob::math::pinv(const blitz::Array<double,2>& A, ...@@ -26,18 +27,6 @@ void bob::math::pinv(const blitz::Array<double,2>& A,
bob::core::array::assertSameDimensionLength(B.extent(0), N); bob::core::array::assertSameDimensionLength(B.extent(0), N);
bob::core::array::assertSameDimensionLength(B.extent(1), M); bob::core::array::assertSameDimensionLength(B.extent(1), M);
// Calls pinv_()
bob::math::pinv_(A, B, rcond);
}
void bob::math::pinv_(const blitz::Array<double,2>& A,
blitz::Array<double,2>& B, const double rcond)
{
// Size variables
const int M = A.extent(0);
const int N = A.extent(1);
const int nb_singular = std::min(M,N);
// Allocates arrays for the SVD // Allocates arrays for the SVD
blitz::Array<double,2> U(M,M); blitz::Array<double,2> U(M,M);
blitz::Array<double,1> sigma(nb_singular); blitz::Array<double,1> sigma(nb_singular);
...@@ -46,7 +35,7 @@ void bob::math::pinv_(const blitz::Array<double,2>& A, ...@@ -46,7 +35,7 @@ void bob::math::pinv_(const blitz::Array<double,2>& A,
blitz::Array<double,2> sigmai_Ut(N, M); blitz::Array<double,2> sigmai_Ut(N, M);
blitz::Array<double,1> Ut_sum(M); blitz::Array<double,1> Ut_sum(M);
// Computes the SVD // Computes the SVD
bob::math::svd_(A, U, sigma, Vt); bob::math::svd(A, U, sigma, Vt);
// Cuts off small sigmas for the inverse similar to Numpy: // Cuts off small sigmas for the inverse similar to Numpy:
// cf. https://github.com/numpy/numpy/blob/maintenance/1.7.x/numpy/linalg/linalg.py#L1577 // cf. https://github.com/numpy/numpy/blob/maintenance/1.7.x/numpy/linalg/linalg.py#L1577
...@@ -68,4 +57,3 @@ void bob::math::pinv_(const blitz::Array<double,2>& A, ...@@ -68,4 +57,3 @@ void bob::math::pinv_(const blitz::Array<double,2>& A,
blitz::Array<double,2> V = Vt.transpose(1,0); blitz::Array<double,2> V = Vt.transpose(1,0);
bob::math::prod(V, sigmai_Ut, B); bob::math::prod(V, sigmai_Ut, B);
} }
...@@ -34,7 +34,7 @@ void bob::math::sqrtSymReal(const blitz::Array<double,2>& A, ...@@ -34,7 +34,7 @@ void bob::math::sqrtSymReal(const blitz::Array<double,2>& A,
blitz::Array<double,2> Vt = V.transpose(1,0); blitz::Array<double,2> Vt = V.transpose(1,0);
blitz::Array<double,1> D(N); blitz::Array<double,1> D(N);
blitz::Array<double,2> tmp(N,N); // Cache for multiplication blitz::Array<double,2> tmp(N,N); // Cache for multiplication
bob::math::eigSym_(A,V,D); bob::math::eigSym(A,V,D);
// 2/ Updates the diagonal matrix D, such that D=sqrt(|D|) // 2/ Updates the diagonal matrix D, such that D=sqrt(|D|)
// |.| is used to deal with values close to zero (-epsilon) // |.| is used to deal with values close to zero (-epsilon)
......
...@@ -79,10 +79,10 @@ static void svd_lapack( const char jobz, const int M, const int N, ...@@ -79,10 +79,10 @@ static void svd_lapack( const char jobz, const int M, const int N,
if (info != 0) if (info != 0)
throw std::runtime_error("The LAPACK dgesdd function returned a non-zero value. You may consider using LAPACK dgsevd instead (see #171) by enabling the 'safe' option."); throw std::runtime_error("The LAPACK dgesdd function returned a non-zero value. You may consider using LAPACK dgsevd instead (see #171) by enabling the 'safe' option.");
} }
// Defining the sign of the eigenvectors // Defining the sign of the eigenvectors
// Approch extracted from page 8 - http://prod.sandia.gov/techlib/access-control.cgi/2007/076422.pdf // Approch extracted from page 8 - http://prod.sandia.gov/techlib/access-control.cgi/2007/076422.pdf
if(U[0] < 0){ if(U[0] < 0){
int ucol=0; ucol= (jobz=='A')? M : std::min(M,N); int ucol=0; ucol= (jobz=='A')? M : std::min(M,N);
for (int i=0; i<ldu*ucol; i++){ for (int i=0; i<ldu*ucol; i++){
U[i] = -1*U[i]; U[i] = -1*U[i];
...@@ -114,17 +114,6 @@ void bob::math::svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U, ...@@ -114,17 +114,6 @@ void bob::math::svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
bob::core::array::assertSameDimensionLength(Vt.extent(0), N); bob::core::array::assertSameDimensionLength(Vt.extent(0), N);
bob::core::array::assertSameDimensionLength(Vt.extent(1), N); bob::core::array::assertSameDimensionLength(Vt.extent(1), N);
bob::math::svd_(A, U, sigma, Vt, safe);
}
void bob::math::svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
blitz::Array<double,1>& sigma, blitz::Array<double,2>& Vt, bool safe)
{
// Size variables
const int M = A.extent(0);
const int N = A.extent(1);
const int nb_singular = std::min(M,N);
// Prepares to call LAPACK function: // Prepares to call LAPACK function:
// We will decompose A^T rather than A to reduce the required number of copy // We will decompose A^T rather than A to reduce the required number of copy
// We recall that FORTRAN/LAPACK is column-major order whereas blitz arrays // We recall that FORTRAN/LAPACK is column-major order whereas blitz arrays
...@@ -161,7 +150,7 @@ void bob::math::svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U, ...@@ -161,7 +150,7 @@ void bob::math::svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
// Call the LAPACK function // Call the LAPACK function
svd_lapack(jobz, N, M, A_lapack, lda, S_lapack, U_lapack, ldu, svd_lapack(jobz, N, M, A_lapack, lda, S_lapack, U_lapack, ldu,
VT_lapack, ldvt, safe); VT_lapack, ldvt, safe);
// Copy singular vectors back to U, V and sigma if required // Copy singular vectors back to U, V and sigma if required
...@@ -188,17 +177,6 @@ void bob::math::svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U, ...@@ -188,17 +177,6 @@ void bob::math::svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
bob::core::array::assertSameDimensionLength(U.extent(1), nb_singular); bob::core::array::assertSameDimensionLength(U.extent(1), nb_singular);
bob::core::array::assertSameDimensionLength(sigma.extent(0), nb_singular); bob::core::array::assertSameDimensionLength(sigma.extent(0), nb_singular);
bob::math::svd_(A, U, sigma, safe);
}
void bob::math::svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
blitz::Array<double,1>& sigma, bool safe)
{
// Size variables
const int M = A.extent(0);
const int N = A.extent(1);
const int nb_singular = std::min(M,N);
// Prepares to call LAPACK function // Prepares to call LAPACK function
// Initialises LAPACK variables // Initialises LAPACK variables
...@@ -249,16 +227,6 @@ void bob::math::svd(const blitz::Array<double,2>& A, blitz::Array<double,1>& sig ...@@ -249,16 +227,6 @@ void bob::math::svd(const blitz::Array<double,2>& A, blitz::Array<double,1>& sig
// Checks and resizes if required // Checks and resizes if required
bob::core::array::assertSameDimensionLength(sigma.extent(0), nb_singular); bob::core::array::assertSameDimensionLength(sigma.extent(0), nb_singular);
bob::math::svd_(A, sigma, safe);
}
void bob::math::svd_(const blitz::Array<double,2>& A, blitz::Array<double,1>& sigma, bool safe)
{
// Size variables
const int M = A.extent(0);
const int N = A.extent(1);
const int nb_singular = std::min(M,N);
// Prepares to call LAPACK function // Prepares to call LAPACK function
// Initialises LAPACK variables // Initialises LAPACK variables
......
...@@ -32,23 +32,6 @@ namespace bob { namespace math { ...@@ -32,23 +32,6 @@ namespace bob { namespace math {
blitz::Array<std::complex<double>,2>& V, blitz::Array<std::complex<double>,2>& V,
blitz::Array<std::complex<double>,1>& D); blitz::Array<std::complex<double>,1>& D);
/**
* @brief Function which performs an eigenvalue decomposition of a real
* matrix, using the LAPACK function <code>dgeev</code>. This version does
* <b>NOT</b> perform any checks on the input data and should be, therefore,
* faster.
*
* @param A The A matrix to decompose (size NxN)
*
* @param V The V matrix of eigenvectors (size NxN) stored in columns
*
* @param D The vector of eigenvalues (size N) (in ascending order). Note
* that this is a 1D array rather than a 2D diagonal matrix!
*/
void eig_(const blitz::Array<double,2>& A,
blitz::Array<std::complex<double>,2>& V,
blitz::Array<std::complex<double>,1>& D);
/** /**
* @brief Function which performs an eigenvalue decomposition of a real * @brief Function which performs an eigenvalue decomposition of a real
* symmetric matrix, using the LAPACK function <code>dsyevd</code>. * symmetric matrix, using the LAPACK function <code>dsyevd</code>.
...@@ -65,24 +48,6 @@ namespace bob { namespace math { ...@@ -65,24 +48,6 @@ namespace bob { namespace math {
void eigSym(const blitz::Array<double,2>& A, blitz::Array<double,2>& V, void eigSym(const blitz::Array<double,2>& A, blitz::Array<double,2>& V,
blitz::Array<double,1>& D); blitz::Array<double,1>& D);
/**
* @brief Function which performs an eigenvalue decomposition of a real
* symmetric matrix, using the LAPACK function <code>dsyevd</code>. This
* version does <b>NOT</b> perform any checks on the input data and should
* be, therefore, faster.
*
* @warning The input matrix should be symmetric.
*
* @param A The A matrix to decompose (size NxN)
*
* @param V The V matrix of eigenvectors (size NxN) stored in columns
*
* @param D The vector of eigenvalues (size N) (in ascending order). Note
* that this is a 1D array rather than a 2D diagonal matrix!
*/
void eigSym_(const blitz::Array<double,2>& A, blitz::Array<double,2>& V,
blitz::Array<double,1>& D);
/** /**
* @brief Computes all the eigenvalues and the eigenvectors of a real * @brief Computes all the eigenvalues and the eigenvectors of a real
* generalized symmetric-definite eigenproblem, of the form: * generalized symmetric-definite eigenproblem, of the form:
...@@ -103,28 +68,6 @@ namespace bob { namespace math { ...@@ -103,28 +68,6 @@ namespace bob { namespace math {
void eigSym(const blitz::Array<double,2>& A, const blitz::Array<double,2>& B, void eigSym(const blitz::Array<double,2>& A, const blitz::Array<double,2>& B,
blitz::Array<double,2>& V, blitz::Array<double,1>& D); blitz::Array<double,2>& V, blitz::Array<double,1>& D);
/**
* @brief Computes all the eigenvalues and the eigenvectors of a real
* generalized symmetric-definite eigenproblem, of the form:
* A*x=(lambda)*B*x, using the LAPACK function <code>dsygvd</code>. This
* version does <b>NOT</b> perform any checks on the input data and should
* be, therefore, faster.
*
* @warning The input matrices A and B are assumed to be symmetric and B is
* also positive definite.
*
* @param A The A input matrix (size NxN) of the problem
*
* @param B The B input matrix (size NxN) of the problem
*
* @param V The V matrix of eigenvectors (size NxN) stored in columns
*
* @param D The vector of eigenvalues (size N) (in ascending order). Note
* that this is a 1D array rather than a 2D diagonal matrix!
*/
void eigSym_(const blitz::Array<double,2>& A, const blitz::Array<double,2>& B,
blitz::Array<double,2>& V, blitz::Array<double,1>& D);
}} }}
#endif /* BOB_MATH_EIG_H */ #endif /* BOB_MATH_EIG_H */
...@@ -27,20 +27,7 @@ namespace bob { namespace math { ...@@ -27,20 +27,7 @@ namespace bob { namespace math {
*/ */
void pinv(const blitz::Array<double,2>& A, blitz::Array<double,2>& B, void pinv(const blitz::Array<double,2>& A, blitz::Array<double,2>& B,
const double rcond=1e-15); const double rcond=1e-15);
/**
* @brief Function which computes the pseudo-inverse using the SVD method.
* @warning The output blitz::array B should have the correct
* size, with zero base index. Checks are NOT performed.
* @param A The A matrix to decompose (size MxN)
* @param B The pseudo-inverse of the matrix A (size NxM)
* @param rcond Cutoff for small singular values. Singular values smaller
* (in modulus) than rcond * largest_singular_value (again, in modulus)
* are set to zero.
*/
void pinv_(const blitz::Array<double,2>& A, blitz::Array<double,2>& B,
const double rcond=1e-15);
}} }}
#endif /* BOB_MATH_PINV_H */ #endif /* BOB_MATH_PINV_H */
...@@ -23,15 +23,20 @@ namespace bob { namespace math { ...@@ -23,15 +23,20 @@ namespace bob { namespace math {
* organized row-wise (each sample is a row, each feature is a column). * organized row-wise (each sample is a row, each feature is a column).
* Outputs the sample mean M and the scatter matrix S. * Outputs the sample mean M and the scatter matrix S.
* *
* @warning No checks are performed on the array sizes and is recommended * The input and output data have their sizes checked and this method will
* only in scenarios where you have previously checked conformity and is * raise an appropriate exception if that is not cased.
* focused only on speed.
* *
* This version of the method also returns the sample mean of the array. * This version of the method also returns the sample mean of the array.
*/ */
template<typename T> template<typename T>
void scatter_(const blitz::Array<T,2>& A, blitz::Array<T,2>& S, void scatter(const blitz::Array<T,2>& A, blitz::Array<T,2>& S,
blitz::Array<T,1>& M) { blitz::Array<T,1>& M) {
// Check output
bob::core::array::assertSameDimensionLength(A.extent(1), M.extent(0));
bob::core::array::assertSameDimensionLength(A.extent(1), S.extent(0));
bob::core::array::assertSameDimensionLength(A.extent(1), S.extent(1));
blitz::firstIndex i; blitz::firstIndex i;
blitz::secondIndex j; blitz::secondIndex j;
blitz::Range a = blitz::Range::all(); blitz::Range a = blitz::Range::all();
...@@ -46,52 +51,13 @@ namespace bob { namespace math { ...@@ -46,52 +51,13 @@ namespace bob { namespace math {
} }
} }
/**
* @brief Computes the scatter matrix of a 2D array considering data is
* organized row-wise (each sample is a row, each feature is a column).
* Outputs the sample mean M and the scatter matrix S.
*
* The input and output data have their sizes checked and this method will
* raise an appropriate exception if that is not cased. If you know that
* the input and output matrices conform, use the scatter_() variant.
*
* This version of the method also returns the sample mean of the array.
*/
template<typename T>
void scatter(const blitz::Array<T,2>& A, blitz::Array<T,2>& S,
blitz::Array<T,1>& M) {
// Check output
bob::core::array::assertSameDimensionLength(A.extent(1), M.extent(0));
bob::core::array::assertSameDimensionLength(A.extent(1), S.extent(0));
bob::core::array::assertSameDimensionLength(A.extent(1), S.extent(1));
scatter_<T>(A, S, M);
}
/**
* @brief Computes the scatter matrix of a 2D array considering data is
* organized row-wise (each sample is a row, each feature is a column).
* Outputs the sample scatter matrix S.
*
* @warning No checks are performed on the array sizes and is recommended
* only in scenarios where you have previously checked conformity and is
* focused only on speed.
*/
template<typename T>
void scatter_(const blitz::Array<T,2>& A, blitz::Array<T,2>& S) {
blitz::Array<T,1> M;
scatter_<T>(A, S, M);
}
/** /**
* @brief Computes the scatter matrix of a 2D array considering data is * @brief Computes the scatter matrix of a 2D array considering data is
* organized row-wise (each sample is a row, each feature is a column). * organized row-wise (each sample is a row, each feature is a column).
* Outputs the sample scatter matrix S. * Outputs the sample scatter matrix S.
* *
* The input and output data have their sizes checked and this method will * The input and output data have their sizes checked and this method will
* raise an appropriate exception if that is not cased. If you know that * raise an appropriate exception if that is not cased.
* the input and output matrices conform, use the scatter_() variant.
*/ */
template<typename T> template<typename T>
void scatter(const blitz::Array<T,2>& A, blitz::Array<T,2>& S) { void scatter(const blitz::Array<T,2>& A, blitz::Array<T,2>& S) {
...@@ -149,16 +115,20 @@ namespace bob { namespace math { ...@@ -149,16 +115,20 @@ namespace bob { namespace math {
* *
* This method was designed based on the previous design at * This method was designed based on the previous design at
* torch3Vision 2.1, by SM. * torch3Vision 2.1, by SM.
*
* @warning No checks are performed on the array sizes and is recommended
* only in scenarios where you have previously checked conformity and is
* focused only on speed.
*/ */
template <typename T> template <typename T>
void scatters_(const std::vector<blitz::Array<T,2> >& data, void scatters(const std::vector<blitz::Array<T,2> >& data,
blitz::Array<T,2>& Sw, blitz::Array<T,2>& Sb, blitz::Array<T,2>& Sw, blitz::Array<T,2>& Sb,
blitz::Array<T,1>& m) blitz::Array<T,1>& m)
{ {
// Check output
for (size_t i=0; i<data.size(); ++i)
bob::core::array::assertSameDimensionLength(data[i].extent(1), m.extent(0));
bob::core::array::assertSameDimensionLength(m.extent(0), Sw.extent(0));
bob::core::array::assertSameDimensionLength(m.extent(0), Sw.extent(1));
bob::core::array::assertSameDimensionLength(m.extent(0), Sb.extent(0));
bob::core::array::assertSameDimensionLength(m.extent(0), Sb.extent(1));
// checks for data shape should have been done before... // checks for data shape should have been done before...
const int n_features = data[0].extent(1); const int n_features = data[0].extent(1);
...@@ -193,74 +163,6 @@ namespace bob { namespace math { ...@@ -193,74 +163,6 @@ namespace bob { namespace math {
} }
} }
/**
* @brief Calculates the within and between class scatter matrices Sw and
* Sb. Returns those matrices and the overall means vector (m).
*
* Strategy implemented:
* 1. Evaluate the overall mean (m), class means (m_k) and the total class
* counts (N).
* 2. Evaluate Sw and Sb using normal loops.
*
* Note that Sw and Sb, in this implementation, will be normalized by N-1
* (number of samples) and K (number of classes). This procedure makes
* the eigen values scaled by (N-1)/K, effectively increasing their
* values. The main motivation for this normalization are numerical
* precision concerns with the increasing number of samples causing a
* rather large Sw matrix. A normalization strategy mitigates this
* problem. The eigen vectors will see no effect on this normalization as
* they are normalized in the euclidean sense (||a|| = 1) so that does
* not change those.
*
* This method was designed based on the previous design at
* torch3Vision 2.1, by SM.
*/
template <typename T>
void scatters(const std::vector<blitz::Array<T,2> >& data,
blitz::Array<T,2>& Sw, blitz::Array<T,2>& Sb,
blitz::Array<T,1>& m)
{
// Check output
for (size_t i=0; i<data.size(); ++i)
bob::core::array::assertSameDimensionLength(data[i].extent(1), m.extent(0));
bob::core::array::assertSameDimensionLength(m.extent(0), Sw.extent(0));
bob::core::array::assertSameDimensionLength(m.extent(0), Sw.extent(1));
bob::core::array::assertSameDimensionLength(m.extent(0), Sb.extent(0));
bob::core::array::assertSameDimensionLength(m.extent(0), Sb.extent(1));
scatters_<T>(data, Sw, Sb, m);
}
/**
* @brief Calculates the within and between class scatter matrices Sw and
* Sb. Returns those matrices.
*
* Strategy implemented:
* 1. Evaluate the overall mean (m), class means (m_k) and the total class
* counts (N).
* 2. Evaluate Sw and Sb using normal loops.
*
* Note that Sw and Sb, in this implementation, will be normalized by N-1
* (number of samples) and K (number of classes). This procedure makes
* the eigen values scaled by (N-1)/K, effectively increasing their
* values. The main motivation for this normalization are numerical
* precision concerns with the increasing number of samples causing a
* rather large Sw matrix. A normalization strategy mitigates this
* problem. The eigen vectors will see no effect on this normalization as
* they are normalized in the euclidean sense (||a|| = 1) so that does
* not change those.
*
* This method was designed based on the previous design at
* torch3Vision 2.1, by SM.
*/
template<typename T>
void scatters_(const std::vector<blitz::Array<T,2> >& data,
blitz::Array<T,2>& Sw, blitz::Array<T,2>& Sb)
{
blitz::Array<T,1> M(data[0].extent(1));
scatters_<T>(data, Sw, Sb, M);
}
/** /**
* @brief Calculates the within and between class scatter matrices Sw and * @brief Calculates the within and between class scatter matrices Sw and
* Sb. Returns those matrices. * Sb. Returns those matrices.
......
...@@ -22,7 +22,7 @@ namespace bob { namespace math { ...@@ -22,7 +22,7 @@ namespace bob { namespace math {
/** /**
* @brief Function which performs a 'full' Singular Value Decomposition * @brief Function which performs a 'full' Singular Value Decomposition
* using the divide and conquer routine dgesdd of LAPACK. * using the divide and conquer routine dgesdd of LAPACK.
* @warning The output blitz::array U, sigma and Vt should have the correct * @warning The output blitz::array U, sigma and Vt should have the correct
* size, with zero base index. Checks are performed. * size, with zero base index. Checks are performed.
* @param A The A matrix to decompose (size MxN) * @param A The A matrix to decompose (size MxN)
* @param U The U matrix of left singular vectors (size MxM) * @param U The U matrix of left singular vectors (size MxM)
...@@ -33,28 +33,13 @@ namespace bob { namespace math { ...@@ -33,28 +33,13 @@ namespace bob { namespace math {
*/ */
void svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U, void svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
blitz::Array<double,1>& sigma, blitz::Array<double,2>& Vt, bool safe=false); blitz::Array<double,1>& sigma, blitz::Array<double,2>& Vt, bool safe=false);
/**
* @brief Function which performs a 'full' Singular Value Decomposition
* using the divide and conquer routine dgesdd of LAPACK.
* @warning The output blitz::array U, sigma and Vt should have the correct
* size, with zero base index. Checks are NOT performed.
* @param A The A matrix to decompose (size MxN)
* @param U The U matrix of left singular vectors (size MxM)
* @param sigma The vector of singular values (size min(M,N))
* Please note that this is a 1D array rather than a 2D diagonal matrix!
* @param Vt The V^T matrix of right singular vectors (size NxN)
* @param safe If enabled, use LAPACK dgesvd instead of dgesdd
*/
void svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
blitz::Array<double,1>& sigma, blitz::Array<double,2>& Vt, bool safe=false);
/** /**
* @brief Function which performs a 'partial' Singular Value Decomposition * @brief Function which performs a 'partial' Singular Value Decomposition
* using the 'simple' driver routine dgesvd of LAPACK. It only returns * using the 'simple' driver routine dgesvd of LAPACK. It only returns
* the first min(M,N) columns of U, and is somehow similar to the * the first min(M,N) columns of U, and is somehow similar to the
* 'economical' SVD variant of matlab (except that it does not return V). * 'economical' SVD variant of matlab (except that it does not return V).
* @warning The output blitz::array U and sigma should have the correct * @warning The output blitz::array U and sigma should have the correct
* size, with zero base index. Checks are performed. * size, with zero base index. Checks are performed.
* @param A The A matrix to decompose (size MxN) * @param A The A matrix to decompose (size MxN)
* @param U The U matrix of left singular vectors (size M x min(M,N)) * @param U The U matrix of left singular vectors (size M x min(M,N))
...@@ -64,28 +49,12 @@ void svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U, ...@@ -64,28 +49,12 @@ void svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
*/ */
void svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U, void svd(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
blitz::Array<double,1>& sigma, bool safe=false); blitz::Array<double,1>& sigma, bool safe=false);
/**
* @brief Function which performs a 'partial' Singular Value Decomposition
* using the 'simple' driver routine dgesvd of LAPACK. It only returns
* the first min(M,N) columns of U, and is somehow similar to the
* 'economical' SVD variant of matlab (except that it does not return V).
* @warning The output blitz::array U and sigma should have the correct
* size, with zero base index. Checks are NOT performed.
* @param A The A matrix to decompose (size MxN)
* @param U The U matrix of left singular vectors (size M x min(M,N))
* @param sigma The vector of singular values (size min(M,N))
* Please note that this is a 1D array rather than a 2D diagonal matrix!
* @param safe If enabled, use LAPACK dgesvd instead of dgesdd
*/
void svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
blitz::Array<double,1>& sigma, bool safe=false);
/** /**
* @brief Function which performs a 'partial' Singular Value Decomposition * @brief Function which performs a 'partial' Singular Value Decomposition
* using the 'simple' driver routine dgesvd of LAPACK. It only returns * using the 'simple' driver routine dgesvd of LAPACK. It only returns
* the singular values. * the singular values.
* @warning The output blitz::array sigma should have the correct * @warning The output blitz::array sigma should have the correct
* size, with zero base index. Checks are performed. * size, with zero base index. Checks are performed.
* @param A The A matrix to decompose (size MxN) * @param A The A matrix to decompose (size MxN)
* @param sigma The vector of singular values (size min(M,N)) * @param sigma The vector of singular values (size min(M,N))
...@@ -94,19 +63,6 @@ void svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U, ...@@ -94,19 +63,6 @@ void svd_(const blitz::Array<double,2>& A, blitz::Array<double,2>& U,
*/ */
void svd(const blitz::Array<double,2>& A, blitz::Array<double,1>& sigma, void svd(const blitz::Array<double,2>& A, blitz::Array<double,1>& sigma,
bool safe=false); bool safe=false);
/**
* @brief Function which performs a 'partial' Singular Value Decomposition
* using the 'simple' driver routine dgesvd of LAPACK. It only returns
* the singular values.
* @warning The output blitz::array sigma should have the correct
* size, with zero base index. Checks are NOT performed.
* @param A The A matrix to decompose (size MxN)
* @param sigma The vector of singular values (size min(M,N))
* Please note that this is a 1D array rather than a 2D diagonal matrix!
* @param safe If enabled, use LAPACK dgesvd instead of dgesdd
*/
void svd_(const blitz::Array<double,2>& A, blitz::Array<double,1>& sigma,
bool safe=false);
/** /**
* @} * @}
......
...@@ -169,23 +169,6 @@ static bob::extension::FunctionDoc s_scatter = bob::extension::FunctionDoc( ...@@ -169,23 +169,6 @@ static bob::extension::FunctionDoc s_scatter = bob::extension::FunctionDoc(
.add_return("m", "array_like (float, 1D)", "The mean matrix, with with the row means of ``a``") .add_return("m", "array_like (float, 1D)", "The mean matrix, with with the row means of ``a``")
; ;
static bob::extension::FunctionDoc s_scatter_nocheck = bob::extension::FunctionDoc(
"scatter_",
"Computes scatter matrix of a 2D array.",
".. warning:: This variant does not perform any checks on the input matrices and is faster then :py:func:`scatter`."
"Use it when you are sure your input matrices sizes match.\n\n"
"Computes the scatter matrix of a 2D array *considering data is organized row-wise* (each sample is a row, each feature is a column). "
"The resulting array ``s`` is squared with extents equal to the number of columns in ``a``. "
"The resulting array ``m`` is a 1D array with the row means of ``a``. "
"This function supports many calling modes, but you should provide, at least, the input data matrix ``a``. "
"All non-provided arguments will be allocated internally and returned."
)
.add_prototype("a, s, m")
.add_parameter("a", "array_like (float, 2D)", "The sample matrix, *considering data is organized row-wise* (each sample is a row, each feature is a column)")
.add_parameter("s", "array_like (float, 2D)", "The scatter matrix, squared with extents equal to the number of columns in ``a``")
.add_parameter("m", "array_like (float,1D)", "The mean matrix, with with the row means of ``a``")
;
static bob::extension::FunctionDoc s_scatters = bob::extension::FunctionDoc(