55#include <visp3/core/vpConfig.h>
56#include <visp3/core/vpCPUFeatures.h>
57#include <visp3/core/vpColVector.h>
58#include <visp3/core/vpDebug.h>
59#include <visp3/core/vpException.h>
60#include <visp3/core/vpMath.h>
61#include <visp3/core/vpMatrix.h>
62#include <visp3/core/vpTranslationVector.h>
64#include <Simd/SimdLib.hpp>
66#ifdef VISP_HAVE_LAPACK
68# include <gsl/gsl_eigen.h>
69# include <gsl/gsl_math.h>
70# include <gsl/gsl_linalg.h>
71# elif defined(VISP_HAVE_MKL)
73typedef MKL_INT integer;
75void vpMatrix::blas_dsyev(
char jobz,
char uplo,
unsigned int n_,
double *a_data,
unsigned int lda_,
76 double *w_data,
double *work_data,
int lwork_,
int &info_)
78 MKL_INT n =
static_cast<MKL_INT
>(n_);
79 MKL_INT lda =
static_cast<MKL_INT
>(lda_);
80 MKL_INT lwork =
static_cast<MKL_INT
>(lwork_);
81 MKL_INT info =
static_cast<MKL_INT
>(info_);
83 dsyev(&jobz, &uplo, &n, a_data, &lda, w_data, work_data, &lwork, &info);
87# if defined(VISP_HAVE_LAPACK_BUILT_IN)
88typedef long int integer;
92extern "C" integer dsyev_(
char *jobz,
char *uplo, integer *n,
double *a, integer *lda,
93 double *w,
double *WORK, integer *lwork, integer *info);
95void vpMatrix::blas_dsyev(
char jobz,
char uplo,
unsigned int n_,
double *a_data,
unsigned int lda_,
96 double *w_data,
double *work_data,
int lwork_,
int &info_)
98 integer n =
static_cast<integer
>(n_);
99 integer lda =
static_cast<integer
>(lda_);
100 integer lwork =
static_cast<integer
>(lwork_);
101 integer info =
static_cast<integer
>(info_);
103 dsyev_(&jobz, &uplo, &n, a_data, &lda, w_data, work_data, &lwork, &info);
105 lwork_ =
static_cast<int>(lwork);
106 info_ =
static_cast<int>(info);
111#if !defined(VISP_USE_MSVC) || (defined(VISP_USE_MSVC) && !defined(VISP_BUILD_SHARED_LIBS))
112const unsigned int vpMatrix::m_lapack_min_size_default = 0;
113unsigned int vpMatrix::m_lapack_min_size = vpMatrix::m_lapack_min_size_default;
120 unsigned int ncols,
double svThreshold,
vpMatrix &Ap,
int &rank_out,
123 Ap.
resize(ncols, nrows,
true,
false);
126 double maxsv = sv[0];
130 for (
unsigned int i = 0; i < ncols; i++) {
131 if (sv[i] > maxsv * svThreshold) {
136 unsigned int rank =
static_cast<unsigned int>(rank_out);
138 rank =
static_cast<unsigned int>(*rank_in);
141 for (
unsigned int i = 0; i < ncols; i++) {
142 for (
unsigned int j = 0; j < nrows; j++) {
143 for (
unsigned int k = 0; k < rank; k++) {
144 Ap[i][j] += V[i][k] * U[j][k] / sv[k];
153 for (
unsigned int i = 0; i < nrows; i++) {
154 for (
unsigned int j = 0; j < rank; j++) {
155 (*imA)[i][j] = U[i][j];
162 imAt->
resize(ncols, rank);
163 for (
unsigned int i = 0; i < ncols; i++) {
164 for (
unsigned int j = 0; j < rank; j++) {
165 (*imAt)[i][j] = V[i][j];
172 kerAt->
resize(ncols - rank, ncols);
174 for (
unsigned int k = 0; k < (ncols - rank); k++) {
175 unsigned j = k + rank;
176 for (
unsigned int i = 0; i < V.
getRows(); i++) {
177 (*kerAt)[k][i] = V[i][j];
192 if (((r + nrows) > M.
rowNum) || ((c + ncols) > M.
colNum)) {
194 "Cannot construct a sub matrix (%dx%d) starting at "
195 "position (%d,%d) that is not contained in the "
196 "original matrix (%dx%d)",
200 init(M, r, c, nrows, ncols);
203#if (VISP_CXX_STANDARD >= VISP_CXX_STANDARD_11)
272 :
vpArray2D<double>(nrows, ncols, list) {}
348 unsigned int rnrows = r + nrows;
349 unsigned int cncols = c + ncols;
357 resize(nrows, ncols,
false,
false);
361 for (
unsigned int i = 0; i < nrows; i++) {
362 memcpy((*
this)[i], &M[i + r][c], ncols *
sizeof(
double));
409 unsigned int rnrows = r + nrows;
410 unsigned int cncols = c + ncols;
420 M.
resize(nrows, ncols,
false,
false);
421 for (
unsigned int i = 0; i < nrows; i++) {
422 memcpy(M[i], &(*
this)[i + r][c], ncols *
sizeof(
double));
451 for (
unsigned int i = 0; i <
rowNum; i++) {
452 for (
unsigned int j = 0; j <
colNum; j++) {
491 for (
unsigned int i = 0; i <
rowNum; i++) {
492 for (
unsigned int j = 0; j <
colNum; j++) {
493 At[j][i] = (*this)[i][j];
533 bool useLapack = (
rowNum > vpMatrix::m_lapack_min_size ||
colNum > vpMatrix::m_lapack_min_size);
534#if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
539#if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
540 const double alpha = 1.0;
541 const double beta = 0.0;
542 const char transa =
't';
543 const char transb =
'n';
545 vpMatrix::blas_dgemm(transa, transb,
rowNum,
rowNum,
colNum, alpha,
data,
colNum,
data,
colNum, beta, B.
data,
rowNum);
550 for (
unsigned int i = 0; i <
rowNum; i++) {
551 for (
unsigned int j = i; j <
rowNum; j++) {
557 for (
unsigned int k = 0; k <
colNum; k++)
558 ssum += *(pi++) * *(pj++);
585 bool useLapack = (
rowNum > vpMatrix::m_lapack_min_size ||
colNum > vpMatrix::m_lapack_min_size);
586#if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
591#if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
592 const double alpha = 1.0;
593 const double beta = 0.0;
594 const char transa =
'n';
595 const char transb =
't';
597 vpMatrix::blas_dgemm(transa, transb,
colNum,
colNum,
rowNum, alpha,
data,
colNum,
data,
colNum, beta, B.
data,
colNum);
601 for (
unsigned int i = 0; i <
colNum; i++) {
603 for (
unsigned int j = 0; j < i; j++) {
606 for (
unsigned int k = 0; k <
rowNum; k++) {
607 s += (*(ptr + i)) * (*(ptr + j));
615 for (
unsigned int k = 0; k <
rowNum; k++) {
616 s += (*(ptr + i)) * (*(ptr + i));
665#if (VISP_CXX_STANDARD >= VISP_CXX_STANDARD_11)
679 if (
this != &other) {
691 other.rowPtrs = NULL;
725 if (
dsize !=
static_cast<unsigned int>(list.size())) {
726 resize(1,
static_cast<unsigned int>(list.size()),
false,
false);
729 std::copy(list.begin(), list.end(),
data);
759 unsigned int nrows =
static_cast<unsigned int>(lists.size()), ncols = 0;
760 for (
auto& l : lists) {
761 if (
static_cast<unsigned int>(l.size()) > ncols) {
762 ncols =
static_cast<unsigned int>(l.size());
766 resize(nrows, ncols,
false,
false);
767 auto it = lists.begin();
768 for (
unsigned int i = 0; i <
rowNum; i++, ++it) {
769 std::copy(it->begin(), it->end(),
rowPtrs[i]);
790 for (
unsigned int i = 0; i <
rowNum; i++) {
791 for (
unsigned int j = 0; j <
colNum; j++) {
800 resize(1, 1,
false,
false);
849 unsigned int rows = A.
getRows();
853 for (
unsigned int i = 0; i < rows; i++)
854 (*
this)[i][i] = A[i];
891 for (
unsigned int i = 0; i < min_; i++)
908 unsigned int rows = A.
getRows();
909 DA.
resize(rows, rows,
true);
911 for (
unsigned int i = 0; i < rows; i++)
928 for (
unsigned int j = 0; j < 3; j++)
931 for (
unsigned int j = 0; j < 3; j++) {
933 for (
unsigned int i = 0; i < 3; i++) {
934 t_out[i] +=
rowPtrs[i][j] * tj;
970 bool useLapack = (A.
rowNum > vpMatrix::m_lapack_min_size || A.
colNum > vpMatrix::m_lapack_min_size);
971#if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
976#if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
982 vpMatrix::blas_dgemv(trans, A.
colNum, A.
rowNum, alpha, A.
data, A.
colNum, v.
data, incr, beta, w.
data, incr);
987 for (
unsigned int j = 0; j < A.
colNum; j++) {
989 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1020 bool useLapack = (A.
rowNum > vpMatrix::m_lapack_min_size || A.
colNum > vpMatrix::m_lapack_min_size || B.
colNum > vpMatrix::m_lapack_min_size);
1021#if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1026#if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1027 const double alpha = 1.0;
1028 const double beta = 0.0;
1029 const char trans =
'n';
1030 vpMatrix::blas_dgemm(trans, trans, B.
colNum, A.
rowNum, A.
colNum, alpha, B.
data, B.
colNum, A.
data, A.
colNum, beta,
1036 const unsigned int BcolNum = B.
colNum;
1037 const unsigned int BrowNum = B.
rowNum;
1038 double **BrowPtrs = B.
rowPtrs;
1039 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1040 const double *rowptri = A.
rowPtrs[i];
1042 for (
unsigned int j = 0; j < BcolNum; j++) {
1044 for (
unsigned int k = 0; k < BrowNum; k++)
1045 s += rowptri[k] * BrowPtrs[k][j];
1069 "Cannot multiply (%dx%d) matrix by (%dx%d) matrix as a "
1074 const unsigned int BcolNum = B.
colNum;
1075 const unsigned int BrowNum = B.
rowNum;
1076 double **BrowPtrs = B.
rowPtrs;
1077 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1078 const double *rowptri = A.
rowPtrs[i];
1080 for (
unsigned int j = 0; j < BcolNum; j++) {
1082 for (
unsigned int k = 0; k < BrowNum; k++)
1083 s += rowptri[k] * BrowPtrs[k][j];
1106 "Cannot multiply (%dx%d) matrix by (%dx%d) matrix as a "
1115 bool useLapack = (A.
rowNum > vpMatrix::m_lapack_min_size || A.
colNum > vpMatrix::m_lapack_min_size || B.
colNum > vpMatrix::m_lapack_min_size);
1116#if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1121#if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1122 const double alpha = 1.0;
1123 const double beta = 0.0;
1124 const char trans =
'n';
1125 vpMatrix::blas_dgemm(trans, trans, B.
colNum, A.
rowNum, A.
colNum, alpha, B.
data, B.
colNum, A.
data, A.
colNum, beta,
1131 const unsigned int BcolNum = B.
colNum;
1132 const unsigned int BrowNum = B.
rowNum;
1133 double **BrowPtrs = B.
rowPtrs;
1134 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1135 const double *rowptri = A.
rowPtrs[i];
1137 for (
unsigned int j = 0; j < BcolNum; j++) {
1139 for (
unsigned int k = 0; k < BrowNum; k++)
1140 s += rowptri[k] * BrowPtrs[k][j];
1190 unsigned int RcolNum = R.
getCols();
1191 unsigned int RrowNum = R.
getRows();
1192 for (
unsigned int i = 0; i <
rowNum; i++) {
1195 for (
unsigned int j = 0; j < RcolNum; j++) {
1197 for (
unsigned int k = 0; k < RrowNum; k++)
1198 s += rowptri[k] * R[k][j];
1219 const unsigned int McolNum = M.
getCols();
1220 const unsigned int MrowNum = M.
getRows();
1221 for (
unsigned int i = 0; i <
rowNum; i++) {
1222 const double *rowptri =
rowPtrs[i];
1224 for (
unsigned int j = 0; j < McolNum; j++) {
1226 for (
unsigned int k = 0; k < MrowNum; k++)
1227 s += rowptri[k] * M[k][j];
1253 bool useLapack = (
rowNum > vpMatrix::m_lapack_min_size ||
colNum > vpMatrix::m_lapack_min_size || V.
colNum > vpMatrix::m_lapack_min_size);
1254#if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1259#if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1260 const double alpha = 1.0;
1261 const double beta = 0.0;
1262 const char trans =
'n';
1263 vpMatrix::blas_dgemm(trans, trans, V.
colNum,
rowNum,
colNum, alpha, V.
data, V.
colNum,
data,
colNum, beta,
1292 bool useLapack = (
rowNum > vpMatrix::m_lapack_min_size ||
colNum > vpMatrix::m_lapack_min_size || V.
getCols() > vpMatrix::m_lapack_min_size);
1293#if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1298#if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1299 const double alpha = 1.0;
1300 const double beta = 0.0;
1301 const char trans =
'n';
1302 vpMatrix::blas_dgemm(trans, trans, V.
getCols(),
rowNum,
colNum, alpha, V.
data, V.
getCols(),
data,
colNum, beta,
1334 double **ArowPtrs = A.
rowPtrs;
1335 double **BrowPtrs = B.
rowPtrs;
1336 double **CrowPtrs = C.
rowPtrs;
1338 for (
unsigned int i = 0; i < A.
rowNum; i++)
1339 for (
unsigned int j = 0; j < A.
colNum; j++)
1340 CrowPtrs[i][j] = wB * BrowPtrs[i][j] + wA * ArowPtrs[i][j];
1362 double **ArowPtrs = A.
rowPtrs;
1363 double **BrowPtrs = B.
rowPtrs;
1364 double **CrowPtrs = C.
rowPtrs;
1366 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1367 for (
unsigned int j = 0; j < A.
colNum; j++) {
1368 CrowPtrs[i][j] = BrowPtrs[i][j] + ArowPtrs[i][j];
1395 double **ArowPtrs = A.
rowPtrs;
1396 double **BrowPtrs = B.
rowPtrs;
1397 double **CrowPtrs = C.
rowPtrs;
1399 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1400 for (
unsigned int j = 0; j < A.
colNum; j++) {
1401 CrowPtrs[i][j] = BrowPtrs[i][j] + ArowPtrs[i][j];
1442 double **ArowPtrs = A.
rowPtrs;
1443 double **BrowPtrs = B.
rowPtrs;
1444 double **CrowPtrs = C.
rowPtrs;
1446 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1447 for (
unsigned int j = 0; j < A.
colNum; j++) {
1448 CrowPtrs[i][j] = ArowPtrs[i][j] - BrowPtrs[i][j];
1475 double **ArowPtrs = A.
rowPtrs;
1476 double **BrowPtrs = B.
rowPtrs;
1477 double **CrowPtrs = C.
rowPtrs;
1479 for (
unsigned int i = 0; i < A.
rowNum; i++) {
1480 for (
unsigned int j = 0; j < A.
colNum; j++) {
1481 CrowPtrs[i][j] = ArowPtrs[i][j] - BrowPtrs[i][j];
1506 double **BrowPtrs = B.
rowPtrs;
1508 for (
unsigned int i = 0; i <
rowNum; i++)
1509 for (
unsigned int j = 0; j <
colNum; j++)
1510 rowPtrs[i][j] += BrowPtrs[i][j];
1523 double **BrowPtrs = B.
rowPtrs;
1524 for (
unsigned int i = 0; i <
rowNum; i++)
1525 for (
unsigned int j = 0; j <
colNum; j++)
1526 rowPtrs[i][j] -= BrowPtrs[i][j];
1545 double **ArowPtrs = A.
rowPtrs;
1546 double **CrowPtrs = C.
rowPtrs;
1549 for (
unsigned int i = 0; i < A.
rowNum; i++)
1550 for (
unsigned int j = 0; j < A.
colNum; j++)
1551 CrowPtrs[i][j] = -ArowPtrs[i][j];
1568 for (
unsigned int i = 0; i <
rowNum; i++) {
1569 for (
unsigned int j = 0; j <
colNum; j++) {
1591 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1595 unsigned int Brow = B.
getRows();
1596 unsigned int Bcol = B.
getCols();
1599 C.
resize(Brow, Bcol,
false,
false);
1601 for (
unsigned int i = 0; i < Brow; i++)
1602 for (
unsigned int j = 0; j < Bcol; j++)
1603 C[i][j] = B[i][j] * x;
1614 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1621 for (
unsigned int i = 0; i <
rowNum; i++)
1622 for (
unsigned int j = 0; j <
colNum; j++)
1631 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1635 if (std::fabs(x) < std::numeric_limits<double>::epsilon()) {
1642 double xinv = 1 / x;
1644 for (
unsigned int i = 0; i <
rowNum; i++)
1645 for (
unsigned int j = 0; j <
colNum; j++)
1646 C[i][j] =
rowPtrs[i][j] * xinv;
1654 for (
unsigned int i = 0; i <
rowNum; i++)
1655 for (
unsigned int j = 0; j <
colNum; j++)
1664 for (
unsigned int i = 0; i <
rowNum; i++)
1665 for (
unsigned int j = 0; j <
colNum; j++)
1677 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1681 for (
unsigned int i = 0; i <
rowNum; i++)
1682 for (
unsigned int j = 0; j <
colNum; j++)
1691 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1695 if (std::fabs(x) < std::numeric_limits<double>::epsilon())
1698 double xinv = 1 / x;
1700 for (
unsigned int i = 0; i <
rowNum; i++)
1701 for (
unsigned int j = 0; j <
colNum; j++)
1720 double *optr = out.
data;
1721 for (
unsigned int j = 0; j <
colNum; j++) {
1722 for (
unsigned int i = 0; i <
rowNum; i++) {
1789 unsigned int r1 = m1.
getRows();
1790 unsigned int c1 = m1.
getCols();
1791 unsigned int r2 = m2.
getRows();
1792 unsigned int c2 = m2.
getCols();
1794 out.
resize(r1*r2, c1*c2,
false,
false);
1796 for (
unsigned int r = 0; r < r1; r++) {
1797 for (
unsigned int c = 0; c < c1; c++) {
1798 double alpha = m1[r][c];
1799 double *m2ptr = m2[0];
1800 unsigned int roffset = r * r2;
1801 unsigned int coffset = c * c2;
1802 for (
unsigned int rr = 0; rr < r2; rr++) {
1803 for (
unsigned int cc = 0; cc < c2; cc++) {
1804 out[roffset + rr][coffset + cc] = alpha * *(m2ptr++);
1827 unsigned int r1 = m1.
getRows();
1828 unsigned int c1 = m1.
getCols();
1829 unsigned int r2 = m2.
getRows();
1830 unsigned int c2 = m2.
getCols();
1833 out.
resize(r1 * r2, c1 * c2,
false,
false);
1835 for (
unsigned int r = 0; r < r1; r++) {
1836 for (
unsigned int c = 0; c < c1; c++) {
1837 double alpha = m1[r][c];
1838 double *m2ptr = m2[0];
1839 unsigned int roffset = r * r2;
1840 unsigned int coffset = c * c2;
1841 for (
unsigned int rr = 0; rr < r2; rr++) {
1842 for (
unsigned int cc = 0; cc < c2; cc++) {
1843 out[roffset + rr][coffset + cc] = alpha * *(m2ptr++);
2032#if defined(VISP_HAVE_LAPACK)
2034#elif defined(VISP_HAVE_EIGEN3)
2036#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
2101#if defined(VISP_HAVE_LAPACK)
2103#elif defined(VISP_HAVE_EIGEN3)
2105#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
2111 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2177#if defined(VISP_HAVE_LAPACK)
2179#elif defined(VISP_HAVE_EIGEN3)
2181#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
2187 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2243#if defined(VISP_HAVE_LAPACK)
2245#elif defined(VISP_HAVE_EIGEN3)
2247#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
2252 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2308#if defined(VISP_HAVE_LAPACK)
2310#elif defined(VISP_HAVE_EIGEN3)
2312#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
2317 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2321#ifndef DOXYGEN_SHOULD_SKIP_THIS
2322#if defined(VISP_HAVE_LAPACK)
2361 unsigned int nrows =
getRows();
2362 unsigned int ncols =
getCols();
2368 Ap.
resize(ncols, nrows,
false,
false);
2370 if (nrows < ncols) {
2371 U.
resize(ncols, ncols,
true);
2374 U.
resize(nrows, ncols,
false);
2381 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
2428 unsigned int nrows =
getRows();
2429 unsigned int ncols =
getCols();
2435 Ap.
resize(ncols, nrows,
false,
false);
2437 if (nrows < ncols) {
2438 U.
resize(ncols, ncols,
true);
2441 U.
resize(nrows, ncols,
false);
2448 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
2450 return static_cast<unsigned int>(rank_out);
2501 unsigned int nrows =
getRows();
2502 unsigned int ncols =
getCols();
2508 Ap.
resize(ncols, nrows,
false,
false);
2510 if (nrows < ncols) {
2511 U.
resize(ncols, ncols,
true);
2514 U.
resize(nrows, ncols,
false);
2521 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
2527 return static_cast<unsigned int>(rank_out);
2639 unsigned int nrows =
getRows();
2640 unsigned int ncols =
getCols();
2645 if (nrows < ncols) {
2646 U.
resize(ncols, ncols,
true);
2649 U.
resize(nrows, ncols,
false);
2656 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, &imA, &imAt, &kerAt);
2662 return static_cast<unsigned int>(rank_out);
2703 unsigned int nrows =
getRows();
2704 unsigned int ncols =
getCols();
2706 double svThreshold = 1e-26;
2711 Ap.
resize(ncols, nrows,
false,
false);
2713 if (nrows < ncols) {
2714 U.
resize(ncols, ncols,
true);
2717 U.
resize(nrows, ncols,
false);
2724 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
2777 unsigned int nrows =
getRows();
2778 unsigned int ncols =
getCols();
2780 double svThreshold = 1e-26;
2785 Ap.
resize(ncols, nrows,
false,
false);
2787 if (nrows < ncols) {
2788 U.
resize(ncols, ncols,
true);
2791 U.
resize(nrows, ncols,
false);
2798 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
2858 unsigned int nrows =
getRows();
2859 unsigned int ncols =
getCols();
2861 double svThreshold = 1e-26;
2866 Ap.
resize(ncols, nrows,
false,
false);
2868 if (nrows < ncols) {
2869 U.
resize(ncols, ncols,
true);
2872 U.
resize(nrows, ncols,
false);
2879 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
3002 unsigned int nrows =
getRows();
3003 unsigned int ncols =
getCols();
3005 double svThreshold = 1e-26;
3009 if (nrows < ncols) {
3010 U.
resize(ncols, ncols,
true);
3013 U.
resize(nrows, ncols,
false);
3020 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, &imA, &imAt, &kerAt);
3030#if defined(VISP_HAVE_EIGEN3)
3069 unsigned int nrows =
getRows();
3070 unsigned int ncols =
getCols();
3076 Ap.
resize(ncols, nrows,
false,
false);
3078 if (nrows < ncols) {
3079 U.
resize(ncols, ncols,
true);
3082 U.
resize(nrows, ncols,
false);
3089 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
3136 unsigned int nrows =
getRows();
3137 unsigned int ncols =
getCols();
3143 Ap.
resize(ncols, nrows,
false,
false);
3145 if (nrows < ncols) {
3146 U.
resize(ncols, ncols,
true);
3149 U.
resize(nrows, ncols,
false);
3156 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
3158 return static_cast<unsigned int>(rank_out);
3209 unsigned int nrows =
getRows();
3210 unsigned int ncols =
getCols();
3216 Ap.
resize(ncols, nrows,
false,
false);
3218 if (nrows < ncols) {
3219 U.
resize(ncols, ncols,
true);
3222 U.
resize(nrows, ncols,
false);
3229 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
3235 return static_cast<unsigned int>(rank_out);
3347 unsigned int nrows =
getRows();
3348 unsigned int ncols =
getCols();
3353 if (nrows < ncols) {
3354 U.
resize(ncols, ncols,
true);
3357 U.
resize(nrows, ncols,
false);
3364 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, &imA, &imAt, &kerAt);
3370 return static_cast<unsigned int>(rank_out);
3411 unsigned int nrows =
getRows();
3412 unsigned int ncols =
getCols();
3414 double svThreshold = 1e-26;
3419 Ap.
resize(ncols, nrows,
false,
false);
3421 if (nrows < ncols) {
3422 U.
resize(ncols, ncols,
true);
3425 U.
resize(nrows, ncols,
false);
3432 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
3485 unsigned int nrows =
getRows();
3486 unsigned int ncols =
getCols();
3488 double svThreshold = 1e-26;
3493 Ap.
resize(ncols, nrows,
false,
false);
3495 if (nrows < ncols) {
3496 U.
resize(ncols, ncols,
true);
3499 U.
resize(nrows, ncols,
false);
3506 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
3566 unsigned int nrows =
getRows();
3567 unsigned int ncols =
getCols();
3569 double svThreshold = 1e-26;
3574 Ap.
resize(ncols, nrows,
false,
false);
3576 if (nrows < ncols) {
3577 U.
resize(ncols, ncols,
true);
3580 U.
resize(nrows, ncols,
false);
3587 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
3710 unsigned int nrows =
getRows();
3711 unsigned int ncols =
getCols();
3713 double svThreshold = 1e-26;
3717 if (nrows < ncols) {
3718 U.
resize(ncols, ncols,
true);
3721 U.
resize(nrows, ncols,
false);
3728 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, &imA, &imAt, &kerAt);
3738#if (VISP_HAVE_OPENCV_VERSION >= 0x020101)
3777 unsigned int nrows =
getRows();
3778 unsigned int ncols =
getCols();
3784 Ap.
resize(ncols, nrows,
false,
false);
3786 if (nrows < ncols) {
3787 U.
resize(ncols, ncols,
true);
3790 U.
resize(nrows, ncols,
false);
3797 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
3844 unsigned int nrows =
getRows();
3845 unsigned int ncols =
getCols();
3851 Ap.
resize(ncols, nrows,
false,
false);
3853 if (nrows < ncols) {
3854 U.
resize(ncols, ncols,
true);
3857 U.
resize(nrows, ncols,
false);
3864 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
3866 return static_cast<unsigned int>(rank_out);
3917 unsigned int nrows =
getRows();
3918 unsigned int ncols =
getCols();
3924 Ap.
resize(ncols, nrows,
false,
false);
3926 if (nrows < ncols) {
3927 U.
resize(ncols, ncols,
true);
3930 U.
resize(nrows, ncols,
false);
3937 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, NULL, NULL, NULL);
3943 return static_cast<unsigned int>(rank_out);
4055 unsigned int nrows =
getRows();
4056 unsigned int ncols =
getCols();
4061 if (nrows < ncols) {
4062 U.
resize(ncols, ncols,
true);
4065 U.
resize(nrows, ncols,
false);
4072 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, NULL, &imA, &imAt, &kerAt);
4078 return static_cast<unsigned int>(rank_out);
4119 unsigned int nrows =
getRows();
4120 unsigned int ncols =
getCols();
4122 double svThreshold = 1e-26;
4127 Ap.
resize(ncols, nrows,
false,
false);
4129 if (nrows < ncols) {
4130 U.
resize(ncols, ncols,
true);
4133 U.
resize(nrows, ncols,
false);
4140 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
4193 unsigned int nrows =
getRows();
4194 unsigned int ncols =
getCols();
4196 double svThreshold = 1e-26;
4201 Ap.
resize(ncols, nrows,
false,
false);
4203 if (nrows < ncols) {
4204 U.
resize(ncols, ncols,
true);
4207 U.
resize(nrows, ncols,
false);
4214 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
4274 unsigned int nrows =
getRows();
4275 unsigned int ncols =
getCols();
4277 double svThreshold = 1e-26;
4282 Ap.
resize(ncols, nrows,
false,
false);
4284 if (nrows < ncols) {
4285 U.
resize(ncols, ncols,
true);
4288 U.
resize(nrows, ncols,
false);
4295 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, NULL, NULL, NULL);
4418 unsigned int nrows =
getRows();
4419 unsigned int ncols =
getCols();
4421 double svThreshold = 1e-26;
4425 if (nrows < ncols) {
4426 U.
resize(ncols, ncols,
true);
4429 U.
resize(nrows, ncols,
false);
4436 compute_pseudo_inverse(U, sv_, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, &imA, &imAt, &kerAt);
4511#if defined(VISP_HAVE_LAPACK)
4513#elif defined(VISP_HAVE_EIGEN3)
4515#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
4522 "Install Lapack, Eigen3 or OpenCV 3rd party"));
4594#if defined(VISP_HAVE_LAPACK)
4596#elif defined(VISP_HAVE_EIGEN3)
4598#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
4605 "Install Lapack, Eigen3 or OpenCV 3rd party"));
4685 return pseudoInverse(Ap, sv, svThreshold, imA, imAt, kerAt);
4908#if defined(VISP_HAVE_LAPACK)
4910#elif defined(VISP_HAVE_EIGEN3)
4912#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
4922 "Install Lapack, Eigen3 or OpenCV 3rd party"));
5067#if defined(VISP_HAVE_LAPACK)
5069#elif defined(VISP_HAVE_EIGEN3)
5071#elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
5081 "Install Lapack, Eigen3 or OpenCV 3rd party"));
5131 for (
unsigned int i = 0; i < column_size; i++)
5132 c[i] = (*
this)[i_begin + i][j];
5265 if (r.
data != NULL &&
data != NULL) {
5266 memcpy(r.
data, (*
this)[i] + j_begin, row_size*
sizeof(
double));
5313 diag.resize(min_size,
false);
5315 for (
unsigned int i = 0; i < min_size; i++) {
5316 diag[i] = (*this)[i][i];
5355 unsigned int nra = A.
getRows();
5356 unsigned int nrb = B.
getRows();
5366 std::cerr <<
"A and C must be two different objects!" << std::endl;
5371 std::cerr <<
"B and C must be two different objects!" << std::endl;
5377 if (C.
data != NULL && A.
data != NULL && A.
size() > 0) {
5382 if (C.
data != NULL && B.
data != NULL && B.
size() > 0) {
5419 std::cerr <<
"A and C must be two different objects!" << std::endl;
5458 std::cerr <<
"A and C must be two different objects!" << std::endl;
5505 for (
unsigned int i = 0; i < A.
getRows(); i++) {
5506 for (
unsigned int j = 0; j < A.
getCols(); j++) {
5507 if (i >= r && i < (r + B.
getRows()) && j >= c && j < (c + B.
getCols())) {
5508 C[i][j] = B[i - r][j - c];
5554 unsigned int nca = A.
getCols();
5555 unsigned int ncb = B.
getCols();
5564 if (B.
getRows() == 0 || nca + ncb == 0) {
5565 std::cerr <<
"B.getRows() == 0 || nca+ncb == 0" << std::endl;
5600 typedef std::string::size_type size_type;
5605 std::vector<std::string> values(m * n);
5606 std::ostringstream oss;
5607 std::ostringstream ossFixed;
5608 std::ios_base::fmtflags original_flags = oss.flags();
5611 ossFixed.setf(std::ios::fixed, std::ios::floatfield);
5613 size_type maxBefore = 0;
5614 size_type maxAfter = 0;
5616 for (
unsigned int i = 0; i < m; ++i) {
5617 for (
unsigned int j = 0; j < n; ++j) {
5619 oss << (*this)[i][j];
5620 if (oss.str().find(
"e") != std::string::npos) {
5622 ossFixed << (*this)[i][j];
5623 oss.str(ossFixed.str());
5626 values[i * n + j] = oss.str();
5627 size_type thislen = values[i * n + j].size();
5628 size_type p = values[i * n + j].find(
'.');
5630 if (p == std::string::npos) {
5640 size_type totalLength = length;
5644 maxAfter = (std::min)(maxAfter, totalLength - maxBefore);
5651 if (! intro.empty())
5653 s <<
"[" << m <<
"," << n <<
"]=\n";
5655 for (
unsigned int i = 0; i < m; i++) {
5657 for (
unsigned int j = 0; j < n; j++) {
5658 size_type p = values[i * n + j].find(
'.');
5659 s.setf(std::ios::right, std::ios::adjustfield);
5660 s.width((std::streamsize)maxBefore);
5661 s << values[i * n + j].substr(0, p).c_str();
5664 s.setf(std::ios::left, std::ios::adjustfield);
5665 if (p != std::string::npos) {
5666 s.width((std::streamsize)maxAfter);
5667 s << values[i * n + j].substr(p, maxAfter).c_str();
5669 assert(maxAfter > 1);
5670 s.width((std::streamsize)maxAfter);
5680 s.flags(original_flags);
5682 return (
int)(maxBefore + maxAfter);
5724 for (
unsigned int i = 0; i < this->
getRows(); ++i) {
5725 for (
unsigned int j = 0; j < this->
getCols(); ++j) {
5726 os << (*this)[i][j] <<
", ";
5728 if (this->
getRows() != i + 1) {
5729 os <<
";" << std::endl;
5731 os <<
"]" << std::endl;
5767 os <<
"([ " << std::endl;
5768 for (
unsigned int i = 0; i < this->
getRows(); ++i) {
5770 for (
unsigned int j = 0; j < this->
getCols(); ++j) {
5771 os << (*this)[i][j] <<
", ";
5773 os <<
"]," << std::endl;
5775 os <<
"])" << std::endl;
5808 for (
unsigned int i = 0; i < this->
getRows(); ++i) {
5809 for (
unsigned int j = 0; j < this->
getCols(); ++j) {
5810 os << (*this)[i][j];
5811 if (!(j == (this->
getCols() - 1)))
5857 os <<
"vpMatrix " << matrixName <<
" (" << this->
getRows() <<
", " << this->
getCols() <<
"); " << std::endl;
5859 for (
unsigned int i = 0; i < this->
getRows(); ++i) {
5860 for (
unsigned int j = 0; j < this->
getCols(); ++j) {
5862 os << matrixName <<
"[" << i <<
"][" << j <<
"] = " << (*this)[i][j] <<
"; " << std::endl;
5864 for (
unsigned int k = 0; k <
sizeof(double); ++k) {
5865 os <<
"((unsigned char*)&(" << matrixName <<
"[" << i <<
"][" << j <<
"]) )[" << k <<
"] = 0x" << std::hex
5866 << (
unsigned int)((
unsigned char *)&((*this)[i][j]))[k] <<
"; " << std::endl;
5889 unsigned int rowNumOld =
rowNum;
5920 if (r.
size() == 0) {
5924 unsigned int oldSize =
size();
5929 memcpy(
data + oldSize, r.
data,
sizeof(
double) * r.
size());
5960 if (c.
size() == 0) {
5965 unsigned int oldColNum =
colNum;
5970 for (
unsigned int i = 0; i <
rowNum; i++) {
5971 memcpy(
data + i*
colNum, tmp.
data + i*oldColNum,
sizeof(
double) * oldColNum);
5994 for (
unsigned int i = r; i < (r + A.
getRows()); i++) {
6046 "Cannot compute eigen values on a non square matrix (%dx%d)",
6052 for (
unsigned int i = 0; i <
rowNum; i++) {
6053 for (
unsigned int j = 0; j <
rowNum; j++) {
6055 if (std::fabs(At_A[i][j]) > std::numeric_limits<double>::epsilon()) {
6061#if defined(VISP_HAVE_LAPACK)
6062#if defined(VISP_HAVE_GSL)
6064 gsl_vector *eval = gsl_vector_alloc(
rowNum);
6067 gsl_eigen_symmv_workspace *w = gsl_eigen_symmv_alloc(
rowNum);
6070 unsigned int Atda = (
unsigned int)m->tda;
6071 for (
unsigned int i = 0; i <
rowNum; i++) {
6072 unsigned int k = i * Atda;
6073 for (
unsigned int j = 0; j <
colNum; j++)
6074 m->data[k + j] = (*
this)[i][j];
6076 gsl_eigen_symmv(m, eval, evec, w);
6078 gsl_eigen_symmv_sort(eval, evec, GSL_EIGEN_SORT_ABS_ASC);
6080 for (
unsigned int i = 0; i <
rowNum; i++) {
6081 evalue[i] = gsl_vector_get(eval, i);
6084 gsl_eigen_symmv_free(w);
6085 gsl_vector_free(eval);
6087 gsl_matrix_free(evec);
6091 const char jobz =
'N';
6092 const char uplo =
'U';
6099 lwork =
static_cast<int>(wkopt);
6107 "Eigen values computation is not implemented. "
6108 "You should install Lapack 3rd party"));
6167 "Cannot compute eigen values on a non square matrix (%dx%d)",
6173 for (
unsigned int i = 0; i <
rowNum; i++) {
6174 for (
unsigned int j = 0; j <
rowNum; j++) {
6176 if (std::fabs(At_A[i][j]) > std::numeric_limits<double>::epsilon()) {
6186#if defined(VISP_HAVE_LAPACK)
6187#if defined(VISP_HAVE_GSL)
6189 gsl_vector *eval = gsl_vector_alloc(
rowNum);
6192 gsl_eigen_symmv_workspace *w = gsl_eigen_symmv_alloc(
rowNum);
6195 unsigned int Atda = (
unsigned int)m->tda;
6196 for (
unsigned int i = 0; i <
rowNum; i++) {
6197 unsigned int k = i * Atda;
6198 for (
unsigned int j = 0; j <
colNum; j++)
6199 m->data[k + j] = (*
this)[i][j];
6201 gsl_eigen_symmv(m, eval, evec, w);
6203 gsl_eigen_symmv_sort(eval, evec, GSL_EIGEN_SORT_ABS_ASC);
6205 for (
unsigned int i = 0; i <
rowNum; i++) {
6206 evalue[i] = gsl_vector_get(eval, i);
6208 Atda = (
unsigned int)evec->tda;
6209 for (
unsigned int i = 0; i <
rowNum; i++) {
6210 unsigned int k = i * Atda;
6211 for (
unsigned int j = 0; j <
rowNum; j++) {
6212 evector[i][j] = evec->
data[k + j];
6216 gsl_eigen_symmv_free(w);
6217 gsl_vector_free(eval);
6219 gsl_matrix_free(evec);
6223 const char jobz =
'V';
6224 const char uplo =
'U';
6231 lwork =
static_cast<int>(wkopt);
6240 "Eigen values computation is not implemented. "
6241 "You should install Lapack 3rd party"));
6266 unsigned int nbline =
getRows();
6267 unsigned int nbcol =
getCols();
6272 V.
resize(nbcol, nbcol,
false);
6279 U.
resize(nbcol, nbcol,
true);
6281 U.
resize(nbline, nbcol,
false);
6289 for (
unsigned int i = 0; i < nbcol; i++) {
6290 if (sv[i] > maxsv) {
6295 unsigned int rank = 0;
6296 for (
unsigned int i = 0; i < nbcol; i++) {
6297 if (sv[i] > maxsv * svThreshold) {
6302 kerAt.
resize(nbcol - rank, nbcol);
6303 if (rank != nbcol) {
6304 for (
unsigned int j = 0, k = 0; j < nbcol; j++) {
6306 if ((sv[j] <= maxsv * svThreshold) &&
6307 (std::fabs(V.
getCol(j).
sumSquare()) > std::numeric_limits<double>::epsilon())) {
6308 for (
unsigned int i = 0; i < V.
getRows(); i++) {
6309 kerAt[k][i] = V[i][j];
6337 unsigned int nbrow =
getRows();
6338 unsigned int nbcol =
getCols();
6343 V.
resize(nbcol, nbcol,
false);
6350 U.
resize(nbcol, nbcol,
true);
6352 U.
resize(nbrow, nbcol,
false);
6359 double maxsv = sv[0];
6361 unsigned int rank = 0;
6362 for (
unsigned int i = 0; i < nbcol; i++) {
6363 if (sv[i] > maxsv * svThreshold) {
6368 kerA.
resize(nbcol, nbcol - rank);
6369 if (rank != nbcol) {
6370 for (
unsigned int j = 0, k = 0; j < nbcol; j++) {
6372 if (sv[j] <= maxsv * svThreshold) {
6373 for (
unsigned int i = 0; i < nbcol; i++) {
6374 kerA[i][k] = V[i][j];
6381 return (nbcol - rank);
6402 unsigned int nbrow =
getRows();
6403 unsigned int nbcol =
getCols();
6404 unsigned int dim_ =
static_cast<unsigned int>(dim);
6409 V.
resize(nbcol, nbcol,
false);
6416 U.
resize(nbcol, nbcol,
true);
6418 U.
resize(nbrow, nbcol,
false);
6424 kerA.
resize(nbcol, dim_);
6426 unsigned int rank = nbcol - dim_;
6427 for (
unsigned int k = 0; k < dim_; k++) {
6428 unsigned int j = k + rank;
6429 for (
unsigned int i = 0; i < nbcol; i++) {
6430 kerA[i][k] = V[i][j];
6435 double maxsv = sv[0];
6436 unsigned int rank = 0;
6437 for (
unsigned int i = 0; i < nbcol; i++) {
6438 if (sv[i] > maxsv * 1e-6) {
6442 return (nbcol - rank);
6502 double *b =
new double[size_];
6503 for (
size_t i = 0; i < size_; i++)
6505 gsl_matrix_view m = gsl_matrix_view_array(this->
data,
rowNum, colNum);
6506 gsl_matrix_view em = gsl_matrix_view_array(b,
rowNum,
colNum);
6507 gsl_linalg_exponential_ss(&m.matrix, &em.matrix, 0);
6511 memcpy(expA.
data, b, size_ *
sizeof(
double));
6532 for (
unsigned int i = 0; i <
rowNum; i++) {
6534 for (
unsigned int j = 0; j <
colNum; j++) {
6535 sum += fabs((*
this)[i][j]);
6537 if (
sum > nA || i == 0) {
6546 double sca = 1.0 / pow(2.0, s);
6549 _expE = c * exp + _eye;
6550 _expD = -c * exp + _eye;
6551 for (
int k = 2; k <= q; k++) {
6552 c = c * ((double)(q - k + 1)) / ((
double)(k * (2 * q - k + 1)));
6553 _expcX = exp * _expX;
6556 _expE = _expE + _expcX;
6558 _expD = _expD + _expcX;
6560 _expD = _expD - _expcX;
6564 exp = _expX * _expE;
6565 for (
int k = 1; k <= s; k++) {
6598 for (
unsigned int i = 0; i < col; i++) {
6599 for (
unsigned int j = 0; j < row; j++)
6600 M_comp[i][j] = M[i][j];
6601 for (
unsigned int j = row + 1; j < M.
getRows(); j++)
6602 M_comp[i][j - 1] = M[i][j];
6604 for (
unsigned int i = col + 1; i < M.
getCols(); i++) {
6605 for (
unsigned int j = 0; j < row; j++)
6606 M_comp[i - 1][j] = M[i][j];
6607 for (
unsigned int j = row + 1; j < M.
getRows(); j++)
6608 M_comp[i - 1][j - 1] = M[i][j];
6624 unsigned int nbline =
getRows();
6625 unsigned int nbcol =
getCols();
6630 V.
resize(nbcol, nbcol,
false);
6637 U.
resize(nbcol, nbcol,
true);
6639 U.
resize(nbline, nbcol,
false);
6647 for (
unsigned int i = 0; i < nbcol; i++) {
6648 if (sv[i] > maxsv) {
6654 unsigned int rank = 0;
6655 for (
unsigned int i = 0; i < nbcol; i++) {
6656 if (sv[i] > maxsv * svThreshold) {
6662 double minsv = maxsv;
6663 for (
unsigned int i = 0; i < rank; i++) {
6664 if (sv[i] < minsv) {
6669 if (std::fabs(minsv) > std::numeric_limits<double>::epsilon()) {
6670 return maxsv / minsv;
6673 return std::numeric_limits<double>::infinity();
6691 for (
unsigned int i = 0; i < H.
getCols(); i++) {
6692 HLM[i][i] += alpha * H[i][i];
6706 for (
unsigned int i = 0; i <
dsize; i++) {
6707 double x = *(
data + i);
6724 if(this->
dsize != 0){
6733 unsigned int maxRank = std::min(this->
getCols(), this->
getRows());
6736 unsigned int boundary = std::min(maxRank, w.
size());
6741 for (
unsigned int i = 0; i < boundary; i++) {
6766 for (
unsigned int i = 0; i <
rowNum; i++) {
6768 for (
unsigned int j = 0; j <
colNum; j++) {
6769 x += fabs(*(*(
rowPtrs + i) + j));
6786 double sum_square = 0.0;
6789 for (
unsigned int i = 0; i <
rowNum; i++) {
6790 for (
unsigned int j = 0; j <
colNum; j++) {
6792 sum_square += x * x;
6834 if (mode ==
"valid") {
6841 if (mode ==
"full" || mode ==
"same") {
6842 const unsigned int pad_x =
kernel.getCols()-1;
6843 const unsigned int pad_y =
kernel.getRows()-1;
6845 M_padded.
insert(M, pad_y, pad_x);
6847 if (mode ==
"same") {
6853 }
else if (mode ==
"valid") {
6860 if (mode ==
"same") {
6861 for (
unsigned int i = 0; i < res_same.
getRows(); i++) {
6862 for (
unsigned int j = 0; j < res_same.
getCols(); j++) {
6863 for (
unsigned int k = 0; k <
kernel.getRows(); k++) {
6864 for (
unsigned int l = 0; l <
kernel.getCols(); l++) {
6865 res_same[i][j] += M_padded[i+k][j+l] *
kernel[
kernel.getRows()-k-1][
kernel.getCols()-l-1];
6871 const unsigned int start_i =
kernel.getRows()/2;
6872 const unsigned int start_j =
kernel.getCols()/2;
6873 for (
unsigned int i = 0; i < M.
getRows(); i++) {
6877 for (
unsigned int i = 0; i < res.
getRows(); i++) {
6878 for (
unsigned int j = 0; j < res.
getCols(); j++) {
6879 for (
unsigned int k = 0; k <
kernel.getRows(); k++) {
6880 for (
unsigned int l = 0; l <
kernel.getCols(); l++) {
6881 res[i][j] += M_padded[i+k][j+l] *
kernel[
kernel.getRows()-k-1][
kernel.getCols()-l-1];
6889#if defined(VISP_BUILD_DEPRECATED_FUNCTIONS)
6940 for (
unsigned int j = 0; j <
getCols(); j++)
6941 c[j] = (*
this)[i - 1][j];
6966 for (
unsigned int i = 0; i <
getRows(); i++)
6967 c[i] = (*
this)[i][j - 1];
6979 for (
unsigned int i = 0; i <
rowNum; i++)
6980 for (
unsigned int j = 0; j <
colNum; j++)
6982 (*this)[i][j] = val;
Implementation of a generic 2D array used as base class for matrices and vectors.
unsigned int getCols() const
double * data
Address of the first element of the data array.
double ** rowPtrs
Address of the first element of each rows.
void resize(unsigned int nrows, unsigned int ncols, bool flagNullify=true, bool recopy_=true)
unsigned int rowNum
Number of rows in the array.
unsigned int dsize
Current array size (rowNum * colNum)
unsigned int size() const
Return the number of elements of the 2D array.
unsigned int getRows() const
unsigned int colNum
Number of columns in the array.
Implementation of column vector and the associated operations.
void resize(unsigned int i, bool flagNullify=true)
error that can be emited by ViSP classes.
@ functionNotImplementedError
Function not implemented.
@ dimensionError
Bad dimension.
@ divideByZeroError
Division by zero.
Implementation of an homogeneous matrix and operations on such kind of matrices.
static Type maximum(const Type &a, const Type &b)
Implementation of a matrix and operations on matrices.
void svdLapack(vpColVector &w, vpMatrix &V)
static vpMatrix conv2(const vpMatrix &M, const vpMatrix &kernel, const std::string &mode="full")
vpColVector eigenValues() const
vpMatrix pseudoInverseOpenCV(double svThreshold=1e-6) const
std::ostream & cppPrint(std::ostream &os, const std::string &matrixName="A", bool octet=false) const
vpMatrix & operator<<(double *)
double cond(double svThreshold=1e-6) const
vpMatrix pseudoInverseEigen3(double svThreshold=1e-6) const
vpMatrix hadamard(const vpMatrix &m) const
vpMatrix operator-() const
vp_deprecated void setIdentity(const double &val=1.0)
vpMatrix & operator/=(double x)
Divide all the element of the matrix by x : Aij = Aij / x.
int print(std::ostream &s, unsigned int length, const std::string &intro="") const
std::ostream & maplePrint(std::ostream &os) const
unsigned int kernel(vpMatrix &kerAt, double svThreshold=1e-6) const
vpMatrix inverseByLU() const
vpMatrix operator*(const vpMatrix &B) const
vpMatrix operator/(double x) const
Cij = Aij / x (A is unchanged)
void stack(const vpMatrix &A)
vp_deprecated void stackMatrices(const vpMatrix &A)
vpMatrix & operator-=(const vpMatrix &B)
Operation A = A - B.
vpMatrix & operator*=(double x)
Multiply all the element of the matrix by x : Aij = Aij * x.
vpMatrix operator*(const double &x, const vpMatrix &B)
static vpMatrix juxtaposeMatrices(const vpMatrix &A, const vpMatrix &B)
vpColVector stackColumns()
vp_deprecated vpColVector column(unsigned int j)
static void sub2Matrices(const vpMatrix &A, const vpMatrix &B, vpMatrix &C)
void svd(vpColVector &w, vpMatrix &V)
vpMatrix operator+(const vpMatrix &B) const
void diag(const double &val=1.0)
vpRowVector getRow(unsigned int i) const
void svdOpenCV(vpColVector &w, vpMatrix &V)
vpMatrix & operator+=(const vpMatrix &B)
Operation A = A + B.
static void add2WeightedMatrices(const vpMatrix &A, const double &wA, const vpMatrix &B, const double &wB, vpMatrix &C)
double inducedL2Norm() const
double infinityNorm() const
vpMatrix & operator=(const vpArray2D< double > &A)
double frobeniusNorm() const
vpColVector getDiag() const
static void createDiagonalMatrix(const vpColVector &A, vpMatrix &DA)
void solveBySVD(const vpColVector &B, vpColVector &x) const
static void computeHLM(const vpMatrix &H, const double &alpha, vpMatrix &HLM)
vpMatrix pseudoInverseLapack(double svThreshold=1e-6) const
vpMatrix pseudoInverse(double svThreshold=1e-6) const
vpMatrix & operator,(double val)
static void multMatrixVector(const vpMatrix &A, const vpColVector &v, vpColVector &w)
vp_deprecated vpRowVector row(unsigned int i)
vpColVector getCol(unsigned int j) const
double det(vpDetMethod method=LU_DECOMPOSITION) const
vp_deprecated void init()
static void add2Matrices(const vpMatrix &A, const vpMatrix &B, vpMatrix &C)
void insert(const vpMatrix &A, unsigned int r, unsigned int c)
void svdEigen3(vpColVector &w, vpMatrix &V)
void kron(const vpMatrix &m1, vpMatrix &out) const
vpMatrix transpose() const
static void mult2Matrices(const vpMatrix &A, const vpMatrix &B, vpMatrix &C)
static void negateMatrix(const vpMatrix &A, vpMatrix &C)
std::ostream & csvPrint(std::ostream &os) const
unsigned int nullSpace(vpMatrix &kerA, double svThreshold=1e-6) const
std::ostream & matlabPrint(std::ostream &os) const
double euclideanNorm() const
vpMatrix extract(unsigned int r, unsigned int c, unsigned int nrows, unsigned int ncols) const
Implementation of a rotation matrix and operations on such kind of matrices.
Implementation of row vector and the associated operations.
void resize(unsigned int i, bool flagNullify=true)
Class that consider the case of a translation vector.