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)
73 typedef MKL_INT integer;
75 void 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)
88 typedef long int integer;
92 extern "C" integer dsyev_(
char *jobz,
char *uplo, integer *n,
double *a, integer *lda,
93 double *w,
double *WORK, integer *lwork, integer *info);
95 void 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))
112 const unsigned int vpMatrix::m_lapack_min_size_default = 0;
113 unsigned 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)
271 vpMatrix::vpMatrix(
unsigned int nrows,
unsigned int ncols,
const std::initializer_list<double> &list)
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"));
6168 "Cannot compute eigen values on a non square matrix (%dx%d)",
6174 for (
unsigned int i = 0; i <
rowNum; i++) {
6175 for (
unsigned int j = 0; j <
rowNum; j++) {
6177 if (std::fabs(At_A[i][j]) > std::numeric_limits<double>::epsilon()) {
6187 #if defined(VISP_HAVE_LAPACK)
6188 #if defined(VISP_HAVE_GSL)
6190 gsl_vector *eval = gsl_vector_alloc(
rowNum);
6193 gsl_eigen_symmv_workspace *w = gsl_eigen_symmv_alloc(
rowNum);
6196 unsigned int Atda = (
unsigned int)m->tda;
6197 for (
unsigned int i = 0; i <
rowNum; i++) {
6198 unsigned int k = i * Atda;
6199 for (
unsigned int j = 0; j <
colNum; j++)
6200 m->data[k + j] = (*
this)[i][j];
6202 gsl_eigen_symmv(m, eval, evec, w);
6204 gsl_eigen_symmv_sort(eval, evec, GSL_EIGEN_SORT_ABS_ASC);
6206 for (
unsigned int i = 0; i <
rowNum; i++) {
6207 evalue[i] = gsl_vector_get(eval, i);
6209 Atda = (
unsigned int)evec->tda;
6210 for (
unsigned int i = 0; i <
rowNum; i++) {
6211 unsigned int k = i * Atda;
6212 for (
unsigned int j = 0; j <
rowNum; j++) {
6213 evector[i][j] = evec->
data[k + j];
6217 gsl_eigen_symmv_free(w);
6218 gsl_vector_free(eval);
6220 gsl_matrix_free(evec);
6224 const char jobz =
'V';
6225 const char uplo =
'U';
6232 lwork =
static_cast<int>(wkopt);
6241 "Eigen values computation is not implemented. "
6242 "You should install Lapack 3rd party"));
6267 unsigned int nbline =
getRows();
6268 unsigned int nbcol =
getCols();
6273 V.
resize(nbcol, nbcol,
false);
6280 U.
resize(nbcol, nbcol,
true);
6282 U.
resize(nbline, nbcol,
false);
6290 for (
unsigned int i = 0; i < nbcol; i++) {
6291 if (sv[i] > maxsv) {
6296 unsigned int rank = 0;
6297 for (
unsigned int i = 0; i < nbcol; i++) {
6298 if (sv[i] > maxsv * svThreshold) {
6303 kerAt.
resize(nbcol - rank, nbcol);
6304 if (rank != nbcol) {
6305 for (
unsigned int j = 0, k = 0; j < nbcol; j++) {
6307 if ((sv[j] <= maxsv * svThreshold) &&
6308 (std::fabs(V.
getCol(j).
sumSquare()) > std::numeric_limits<double>::epsilon())) {
6309 for (
unsigned int i = 0; i < V.
getRows(); i++) {
6310 kerAt[k][i] = V[i][j];
6338 unsigned int nbrow =
getRows();
6339 unsigned int nbcol =
getCols();
6344 V.
resize(nbcol, nbcol,
false);
6351 U.
resize(nbcol, nbcol,
true);
6353 U.
resize(nbrow, nbcol,
false);
6360 double maxsv = sv[0];
6362 unsigned int rank = 0;
6363 for (
unsigned int i = 0; i < nbcol; i++) {
6364 if (sv[i] > maxsv * svThreshold) {
6369 kerA.
resize(nbcol, nbcol - rank);
6370 if (rank != nbcol) {
6371 for (
unsigned int j = 0, k = 0; j < nbcol; j++) {
6373 if (sv[j] <= maxsv * svThreshold) {
6374 for (
unsigned int i = 0; i < nbcol; i++) {
6375 kerA[i][k] = V[i][j];
6382 return (nbcol - rank);
6403 unsigned int nbrow =
getRows();
6404 unsigned int nbcol =
getCols();
6405 unsigned int dim_ =
static_cast<unsigned int>(dim);
6410 V.
resize(nbcol, nbcol,
false);
6417 U.
resize(nbcol, nbcol,
true);
6419 U.
resize(nbrow, nbcol,
false);
6425 kerA.
resize(nbcol, dim_);
6427 unsigned int rank = nbcol - dim_;
6428 for (
unsigned int k = 0; k < dim_; k++) {
6429 unsigned int j = k + rank;
6430 for (
unsigned int i = 0; i < nbcol; i++) {
6431 kerA[i][k] = V[i][j];
6436 double maxsv = sv[0];
6437 unsigned int rank = 0;
6438 for (
unsigned int i = 0; i < nbcol; i++) {
6439 if (sv[i] > maxsv * 1e-6) {
6443 return (nbcol - rank);
6501 #ifdef VISP_HAVE_GSL
6503 double *b =
new double[size_];
6504 for (
size_t i = 0; i < size_; i++)
6506 gsl_matrix_view m = gsl_matrix_view_array(this->
data,
rowNum, colNum);
6507 gsl_matrix_view em = gsl_matrix_view_array(b,
rowNum,
colNum);
6508 gsl_linalg_exponential_ss(&m.matrix, &em.matrix, 0);
6512 memcpy(expA.
data, b, size_ *
sizeof(
double));
6533 for (
unsigned int i = 0; i <
rowNum; i++) {
6535 for (
unsigned int j = 0; j <
colNum; j++) {
6536 sum += fabs((*
this)[i][j]);
6538 if (
sum > nA || i == 0) {
6547 double sca = 1.0 / pow(2.0, s);
6550 _expE = c * exp + _eye;
6551 _expD = -c * exp + _eye;
6552 for (
int k = 2; k <= q; k++) {
6553 c = c * ((double)(q - k + 1)) / ((
double)(k * (2 * q - k + 1)));
6554 _expcX = exp * _expX;
6557 _expE = _expE + _expcX;
6559 _expD = _expD + _expcX;
6561 _expD = _expD - _expcX;
6565 exp = _expX * _expE;
6566 for (
int k = 1; k <= s; k++) {
6599 for (
unsigned int i = 0; i < col; i++) {
6600 for (
unsigned int j = 0; j < row; j++)
6601 M_comp[i][j] = M[i][j];
6602 for (
unsigned int j = row + 1; j < M.
getRows(); j++)
6603 M_comp[i][j - 1] = M[i][j];
6605 for (
unsigned int i = col + 1; i < M.
getCols(); i++) {
6606 for (
unsigned int j = 0; j < row; j++)
6607 M_comp[i - 1][j] = M[i][j];
6608 for (
unsigned int j = row + 1; j < M.
getRows(); j++)
6609 M_comp[i - 1][j - 1] = M[i][j];
6625 unsigned int nbline =
getRows();
6626 unsigned int nbcol =
getCols();
6631 V.
resize(nbcol, nbcol,
false);
6638 U.
resize(nbcol, nbcol,
true);
6640 U.
resize(nbline, nbcol,
false);
6648 for (
unsigned int i = 0; i < nbcol; i++) {
6649 if (sv[i] > maxsv) {
6655 unsigned int rank = 0;
6656 for (
unsigned int i = 0; i < nbcol; i++) {
6657 if (sv[i] > maxsv * svThreshold) {
6663 double minsv = maxsv;
6664 for (
unsigned int i = 0; i < rank; i++) {
6665 if (sv[i] < minsv) {
6670 if (std::fabs(minsv) > std::numeric_limits<double>::epsilon()) {
6671 return maxsv / minsv;
6674 return std::numeric_limits<double>::infinity();
6692 for (
unsigned int i = 0; i < H.
getCols(); i++) {
6693 HLM[i][i] += alpha * H[i][i];
6707 for (
unsigned int i = 0; i <
dsize; i++) {
6708 double x = *(
data + i);
6725 if(this->
dsize != 0){
6734 unsigned int maxRank = std::min(this->
getCols(), this->
getRows());
6737 unsigned int boundary = std::min(maxRank, w.
size());
6742 for (
unsigned int i = 0; i < boundary; i++) {
6767 for (
unsigned int i = 0; i <
rowNum; i++) {
6769 for (
unsigned int j = 0; j <
colNum; j++) {
6770 x += fabs(*(*(
rowPtrs + i) + j));
6787 double sum_square = 0.0;
6790 for (
unsigned int i = 0; i <
rowNum; i++) {
6791 for (
unsigned int j = 0; j <
colNum; j++) {
6793 sum_square += x * x;
6835 if (mode ==
"valid") {
6842 if (mode ==
"full" || mode ==
"same") {
6843 const unsigned int pad_x =
kernel.getCols()-1;
6844 const unsigned int pad_y =
kernel.getRows()-1;
6846 M_padded.
insert(M, pad_y, pad_x);
6848 if (mode ==
"same") {
6854 }
else if (mode ==
"valid") {
6861 if (mode ==
"same") {
6862 for (
unsigned int i = 0; i < res_same.
getRows(); i++) {
6863 for (
unsigned int j = 0; j < res_same.
getCols(); j++) {
6864 for (
unsigned int k = 0; k <
kernel.getRows(); k++) {
6865 for (
unsigned int l = 0; l <
kernel.getCols(); l++) {
6866 res_same[i][j] += M_padded[i+k][j+l] *
kernel[
kernel.getRows()-k-1][
kernel.getCols()-l-1];
6872 const unsigned int start_i =
kernel.getRows()/2;
6873 const unsigned int start_j =
kernel.getCols()/2;
6874 for (
unsigned int i = 0; i < M.
getRows(); i++) {
6878 for (
unsigned int i = 0; i < res.
getRows(); i++) {
6879 for (
unsigned int j = 0; j < res.
getCols(); j++) {
6880 for (
unsigned int k = 0; k <
kernel.getRows(); k++) {
6881 for (
unsigned int l = 0; l <
kernel.getCols(); l++) {
6882 res[i][j] += M_padded[i+k][j+l] *
kernel[
kernel.getRows()-k-1][
kernel.getCols()-l-1];
6890 #if defined(VISP_BUILD_DEPRECATED_FUNCTIONS)
6941 for (
unsigned int j = 0; j <
getCols(); j++)
6942 c[j] = (*
this)[i - 1][j];
6967 for (
unsigned int i = 0; i <
getRows(); i++)
6968 c[i] = (*
this)[i][j - 1];
6980 for (
unsigned int i = 0; i <
rowNum; i++)
6981 for (
unsigned int j = 0; j <
colNum; j++)
6983 (*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.