refactor(*): move Ctain*.h to folder camera_models

This commit is contained in:
Messier
2019-09-05 10:34:04 +08:00
parent 050e7d6771
commit f8d47e6f3a
6 changed files with 1 additions and 1 deletions

View File

@@ -0,0 +1,39 @@
//
// Created by 顾涵彬 on 2019-08-28.
//
#include <cstdint>
#include "Matrix.h"
#include "SquareMatrix.h"
#include "MatrixSolver.h"
#include "Quaternion.h"
#ifndef MATRIX_CTAIN_H
#define MATRIX_CTAIN_H
namespace Ctain {
typedef SMatrix<double> Matrixd;
typedef Matrix<double> MatrixXd;
typedef Matrix<double> Matrix23d;
typedef SMatrix<double> Matrix3d;
typedef SMatrix<double> Matrix4d;
typedef SMatrix<float> Matrixf;
typedef Matrixf Matrix2f;
typedef Matrixf Matrix3f;
typedef Matrix<float> Vectorf;
typedef Vectorf Vector2f;
typedef Vectorf Vector3f;
typedef Matrix<double> Vectord;
typedef Matrix<double> Vector2d;
typedef Matrix<double> Vector3d;
typedef Matrix<double> MatrixXcd;
typedef Quaternion<double> Quaterniond;
} // end namespace Ctain
#endif // Ctain_CtainBASE_H

View File

@@ -0,0 +1,415 @@
//
// Created by 顾涵彬 on 2019-08-28.
//
#ifndef MATRIX_MATRIX_H
#define MATRIX_MATRIX_H
#include <cstring>
#include <iostream>
#include <cmath>
namespace Ctain {
template<typename _Scalar>
class Matrix {
public:
Matrix(int Rows, int Cols) :
_Rows(Rows), _Cols(Cols), _isSub(0), input_id(0) {
_startRow = 0;
_startCol = 0;
_Rows_raw = Rows;
_Cols_raw = Cols;
data = new _Scalar [_Rows * _Cols];
memset(data, 0, _Rows * _Cols * sizeof(_Scalar));
}
Matrix() :
_Rows(0), _Cols(0), _isSub(0), input_id(0) {
_startRow = 0;
_startCol = 0;
_Rows_raw = 0;
_Cols_raw = 0;
}
Matrix(_Scalar _data[], int Rows, int Cols) :
_Rows(Rows), _Cols(Cols), _isSub(0), input_id(0) {
_startRow = 0;
_startCol = 0;
_Rows_raw = Rows;
_Cols_raw = Cols;
data = new _Scalar [_Rows * _Cols];
memcpy(data, _data, _Rows * _Cols * sizeof(_Scalar));
}
Matrix(_Scalar **_data, int Rows, int Cols) :
_Rows(Rows), _Cols(Cols), _isSub(0), input_id(0) {
_startRow = 0;
_startCol = 0;
_Rows_raw = Rows;
_Cols_raw = Cols;
data = new _Scalar [_Rows * _Cols];
for (int i = 0; i < _Rows; ++i) {
memcpy(data + i * _Cols, *(_data + i), _Cols * sizeof(_Scalar));
}
}
template<typename T>
Matrix<T> cast() {
Matrix<T> res(_Rows, _Cols);
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
res(i,j) = cData(i, j);
}
}
return res;
}
void setIdentity() {
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
if(i == j) {
Data(i,j) = 1;
}else {
Data(i,j) = 0;
}
}
}
}
void setZero() {
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
Data(i, j) = 0;
}
}
}
inline int cols() const { return _Cols; }
inline int rows() const { return _Rows; }
inline int size() const { return cols() * rows(); }
_Scalar * addr() {
return data;
}
friend Matrix<_Scalar> &operator <<(Matrix<_Scalar> &m,_Scalar val) {
m.Data(m.input_id++) = val;
return m;
}
friend std::ostream &operator <<(std::ostream &os,const Matrix<_Scalar> &m) {
// os << std::endl;
for (int i = 0; i < m._Rows; i++) {
for (int j = 0; j < m._Cols - 1; j++) {
std::cout.width(10);
std::cout.setf(std::ios::left);
std::cout.precision(6);
os << m.cData(i,j) ;
}
std::cout.width(9);
std::cout.setf(std::ios::left);
std::cout.precision(8);
os << m.cData(i,m._Cols - 1) << std::endl;
}
return os;
}
friend Matrix<_Scalar> operator *(
double a, const Matrix<_Scalar> &m) {
Matrix<_Scalar> res;
res = m;
for (int i = 0; i < m._Rows; i++) {
for (int j = 0; j < m._Cols; j++) {
res.Data(i,j) *= a;
}
}
return res;
}
friend Matrix<_Scalar> operator *(
const Matrix<_Scalar> &m, double a) {
Matrix<_Scalar> res;
res = m;
for (int i = 0; i < m._Rows; i++) {
for (int j = 0; j < m._Cols; j++) {
res.Data(i,j) *= a;
}
}
return res;
}
friend Matrix<_Scalar> operator -(
const Matrix<_Scalar> &m) {
Matrix<_Scalar> res;
res = m;
for (int i = 0; i < m._Rows; i++) {
for (int j = 0; j < m._Cols; j++) {
res.Data(i,j) *= -1;
}
}
return res;
}
friend Matrix<_Scalar> operator +(
double a, const Matrix<_Scalar> &m) {
Matrix<_Scalar> res;
res = m;
for (int i = 0; i < m._Rows; i++) {
for (int j = 0; j < m._Cols; j++) {
res.Data(i,j) += a;
}
}
return res;
}
void operator =(Matrix<_Scalar> m);
// void operator =(Matrix<_Scalar> &m);
Matrix<_Scalar> operator +(const Matrix<_Scalar> &m) const;
Matrix<_Scalar> operator -(const Matrix<_Scalar> &m) const;
Matrix<_Scalar> operator *(const Matrix<_Scalar> &m) const;
Matrix<_Scalar> operator /(double m) const;
_Scalar &operator()(int i, int j) {
return Data(i,j);
}
_Scalar &operator()(int id) {
return Data(id);
}
_Scalar operator()(int id) const {
return cData(id);
}
Matrix<_Scalar> transpose() const;
Matrix<_Scalar> col(int Col) {
return block(0, Col, _Rows, 1);
}
Matrix<_Scalar> row(int Row) {
return block(Row, 0, 1, _Cols);
}
Matrix<_Scalar> block(int sRow, int sCol, int Rows, int Cols) {
Matrix<_Scalar> sub;
sub = *this;
sub.setSub(sRow, sCol, Rows, Cols, data);
return sub;
}
template<int Rows, int Cols>
Matrix<_Scalar> topLeftCorner() const {
Matrix<_Scalar> sub;
sub = *this;
sub.setSub(0, 0, Rows, Cols, data);
return sub;
}
template<int Rows, int Cols>
Matrix<_Scalar> topRightCorner() const {
Matrix<_Scalar> sub;
sub = *this;
sub.setSub(0, _Cols-Cols, Rows, Cols, data);
return sub;
}
void setSub(int sRow, int sCol, int Rows, int Cols, _Scalar *Data) {
_isSub = true;
_Rows_raw = _Rows;
_Cols_raw = _Cols;
_Rows = Rows;
_Cols = Cols;
_startRow = sRow;
_startCol = sCol;
data = Data;
}
void normalize();
double norm() const;
virtual ~Matrix() {
if(!data)
delete[] data;
}
// template<int _Rows, int _Cols>
// inline Matrix<_Scalar, _Rows, _Cols> block<_Rows, _Cols>(Index i,Index j) {
//
// }
inline _Scalar *Data() {
return data;
}
protected:
_Scalar *data;
int _Rows;
int _Cols;
bool _isSub;
int _startRow;
int _startCol;
int _Rows_raw;
int _Cols_raw;
int input_id ;
inline int id(int i, int j) const {
if(_isSub)
return (i + _startRow) * _Cols_raw + j + _startCol;
else
return i * _Cols + j;
}
inline _Scalar &Data(int i,int j) {
return data[id(i,j)];
}
inline _Scalar &Data(int id) {
int i = id / _Cols;
int j = id % _Cols;
int index;
if(_isSub)
index = (i + _startRow) * _Cols_raw + j + _startCol;
else
index = i * _Cols + j;
return data[index];
}
inline _Scalar cData(int id) const{
int i = id / _Cols;
int j = id % _Cols;
int index;
if(_isSub)
index = (i + _startRow) * _Cols_raw + j + _startCol;
else
index = i * _Cols + j;
return data[index];
}
inline _Scalar cData (int i,int j) const{
return data[id(i,j)];
}
}; //class Matrix end
template<typename _Scalar>
void Matrix<_Scalar>::operator =(Matrix<_Scalar> m) {
if(m._isSub) {
if(_isSub) {
for (int i = 0; i < m._Rows; i++) {
for (int j = 0; j < m._Cols; j++) {
Data(i,j) = m.cData(i, j);
}
}
return;
}
_isSub = true;
_Rows = m._Rows;
_Cols = m._Cols;
_Rows_raw = m._Rows_raw;
_Cols_raw = m._Cols_raw;
_startRow = m._startRow;
_startCol = m._startCol;
data = m.Data();
return;
}
if(!_isSub) {
if(size() != m.size()) {
if(size() > 0) {
delete[] data;
}
_Rows = m._Rows;
_Cols = m._Cols;
data = new _Scalar[_Rows * _Cols];
}else {
_Rows = m._Rows;
_Cols = m._Cols;
}
}
for (int i = 0; i < m._Rows; i++) {
for (int j = 0; j < m._Cols; j++) {
Data(i,j) = m.cData(i, j);
}
}
}
template<typename _Scalar>
Matrix<_Scalar> Matrix<_Scalar>::operator +(const Matrix<_Scalar> &m) const{
Matrix<_Scalar> sum;
sum = *this;
for (int i = 0; i < _Rows * _Cols; i++) {
sum.data[i] += m.data[i];
}
return sum;
}
template<typename _Scalar>
Matrix<_Scalar> Matrix<_Scalar>::operator -(const Matrix<_Scalar> &m) const{
Matrix<_Scalar> sum;
sum = *this;
for (int i = 0; i < _Rows * _Cols; i++) {
sum.data[i] -= m.data[i];
}
return sum;
}
template<typename _Scalar>
Matrix<_Scalar> Matrix<_Scalar>::transpose() const
{
Matrix<_Scalar> res(_Cols, _Rows);
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
res.Data(j, i) = cData(i, j);
}
}
return res;
}
template<typename _Scalar>
Matrix<_Scalar> Matrix<_Scalar>::operator *(const Matrix<_Scalar> &m) const {
if(_Cols != m._Rows) {
//todo:output err
return m;
}
Matrix<_Scalar> res(_Rows, m._Cols);
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < m._Cols; j++) {
_Scalar sum = 0;
for (int k = 0; k < _Cols; k++) {
sum += cData(i, k) * m.cData(k, j);
}
res.Data(i,j) = sum;
}
}
return res;
}
template<typename _Scalar>
Matrix<_Scalar> Matrix<_Scalar>::operator /(double m) const {
Matrix<_Scalar> res;
res = *this;
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
res.Data(i,j) /= m;
}
}
return res;
}
template<typename _Scalar>
void Matrix<_Scalar>::normalize() {
double sum = 0;
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
sum += Matrix::cData(i, j);
}
}
sum = sqrt(sum);
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
Matrix::Data(i, j) /= sum;
}
}
}
template<typename _Scalar>
double Matrix<_Scalar>::norm() const{
double sum = 0;
for (int i = 0; i < _Rows; i++) {
for (int j = 0; j < _Cols; j++) {
sum += Matrix::cData(i, j) * Matrix::cData(i, j);
}
}
sum = sqrt(sum);
return sum;
}
} //namespace Ctain end
#endif //MATRIX_MATRIX_H

View File

@@ -0,0 +1,250 @@
//
// Created by 顾涵彬 on 2019-08-30.
//
#ifndef MATRIX_MATRIXSOLVER_H
#define MATRIX_MATRIXSOLVER_H
#include <cmath>
#include <complex>
static bool Matrix_EigenValue(double *K1,int n,int LoopNumber,double Error1,double *Ret);
static void Matrix_Hessenberg(double *A1,int n,double *ret);
namespace Ctain {
class EigenSolver {
public:
EigenSolver(SMatrix<double> s) {
double *A = new double[s.rows()*2];
double *B = new double[s.size()];
for (int i = 0; i < s.size(); i++)
B[i] = s(i);
memset(A, 0, sizeof(s.rows()*2));
Matrix_EigenValue(B, s.rows(),1000,1e-10,A);
Matrix<double> tt(A, s.rows(), 2);
t=tt;
std::cout<<"s:"<<s;
SMatrix<double> s2(A, s.rows());
std::cout<<"tt:"<<tt;
std::cout<<"s2:"<<s2;
delete []A;
delete []B;
}
Matrix<double> eigenvalues() {
return t;
}
private:
Matrix<double> t;
};
} //namespace Ctain end
static void Matrix_Hessenberg(double *A1,int n,double *ret)
{
int MaxNumber;
double temp,*A;
A=new double[n*n];
memset(A, 0, sizeof(double)*n*n);
for (int i=0;i<n;i++) {
int k=i*n;
for (int j=0;j<n;j++)
{
A[k+j]=A1[k+j];
}
}
for (int k=1;k<n-1;k++) {
int i=k-1;
MaxNumber=k;
temp=fabs(A[k*n+i]);
for (int j=k+1;j<n;j++) {
if (fabs(A[j*n+i])>temp) {
temp=fabs(A[j*n+i]);
MaxNumber=j;
}
}
ret[0]=A[MaxNumber*n+i];
if (ret[0]!=0) {
if (MaxNumber!=k) {
for (int j=k-1;j<n;j++) {
temp=A[i*n+j];
A[i*n+j]=A[k*n+j];
A[k*n+j]=temp;
}
for (int j=0;j<n;j++) {
temp=A[j*n+i];
A[j*n+i]=A[j*n+k];
A[j*n+k]=temp;
}
}
for (int i=k+1;i<n;i++) {
temp=A[i*n+k-1]/ret[0];
A[i*n+k-1]=0;
for (int j=k;j<n;j++) {
A[i*n+j]-=temp*A[k*n+j];
}
for (int j=0;j<n;j++) {
A[j*n+k]+=temp*A[j*n+i];
}
}
}
}
for (int i=0;i<n;i++) {
int k=i*n;
for (int j=0;j<n;j++) {
ret[k+j]=A[k+j];
}
}
delete []A;
}
static bool Matrix_EigenValue(double *K1,int n,int LoopNumber,double Error1,double *Ret)
{
int i,j,k,t,m,Loop1;
double b,c,d,g,xy,p,q,r,x,s,e,f,z,y,temp,*A;
A=new double[n*n];
memset(A, 0, sizeof(double)*n*n);
Matrix_Hessenberg(K1,n,A);
m=n;
Loop1=LoopNumber;
while(m!=0) {
t=m-1;
while(t>0) {
temp=fabs(A[(t-1)*n+t-1]);
temp+=fabs(A[t*n+t]);
temp=temp*Error1;
if (fabs(A[t*n+t-1])>temp) {
t--;
}
else {
break;
}
}
if (t==m-1) {
Ret[(m-1)*2]=A[(m-1)*n+m-1];
Ret[(m-1)*2+1]=0;
m-=1;
Loop1=LoopNumber;
}
else if(t==m-2) {
b=-A[(m-1)*n+m-1]-A[(m-2)*n+m-2];
c=A[(m-1)*n+m-1]*A[(m-2)*n+m-2]-A[(m-1)*n+m-2]*A[(m-2)*n+m-1];
d=b*b-4*c;
y=sqrt(fabs(d));
if (d>0) {
xy=1;
if (b<0) {
xy=-1;
}
Ret[(m-1)*2]=-(b+xy*y)/2;
Ret[(m-1)*2+1]=0;
Ret[(m-2)*2]=c/Ret[(m-1)*2];
Ret[(m-2)*2+1]=0;
}
else {
Ret[(m-1)*2]=-b/2;
Ret[(m-2)*2]=-b/2;
Ret[(m-1)*2+1]=y/2;
Ret[(m-2)*2+1]=-y/2;
}
m-=2;
Loop1=LoopNumber;
}
else {
if (Loop1<1) {
return false;
}
Loop1--;
j=t+2;
while (j<m) {
A[j*n+j-2]=0;
j++;
}
j=t+3;
while (j<m) {
A[j*n+j-3]=0;
j++;
}
k=t;
while (k<m-1) {
if (k!=t) {
p=A[k*n+k-1];
q=A[(k+1)*n+k-1];
if (k!=m-2) {
r=A[(k+2)*n+k-1];
}
else {
r=0;
}
}
else {
b=A[(m-1)*n+m-1];
c=A[(m-2)*n+m-2];
x=b+c;
y=b*c-A[(m-2)*n+m-1]*A[(m-1)*n+m-2];
p=A[t*n+t]*(A[t*n+t]-x)+A[t*n+t+1]*A[(t+1)*n+t]+y;
q=A[(t+1)*n+t]*(A[t*n+t]+A[(t+1)*n+t+1]-x);
r=A[(t+1)*n+t]*A[(t+2)*n+t+1];
}
if (p!=0 || q!=0 || r!=0) {
if (p<0) {
xy=-1;
}
else {
xy=1;
}
s=xy*sqrt(p*p+q*q+r*r);
if (k!=t) {
A[k*n+k-1]=-s;
}
e=-q/s;
f=-r/s;
x=-p/s;
y=-x-f*r/(p+s);
g=e*r/(p+s);
z=-x-e*q/(p+s);
for (j=k;j<m;j++) {
b=A[k*n+j];
c=A[(k+1)*n+j];
p=x*b+e*c;
q=e*b+y*c;
r=f*b+g*c;
if (k!=m-2) {
b=A[(k+2)*n+j];
p+=f*b;
q+=g*b;
r+=z*b;
A[(k+2)*n+j]=r;
}
A[(k+1)*n+j]=q;
A[k*n+j]=p;
}
j=k+3;
if (j>m-2) {
j=m-1;
}
for (i=t;i<j+1;i++) {
b=A[i*n+k];
c=A[i*n+k+1];
p=x*b+e*c;
q=e*b+y*c;
r=f*b+g*c;
if (k!=m-2) {
b=A[i*n+k+2];
p+=f*b;
q+=g*b;
r+=z*b;
A[i*n+k+2]=r;
}
A[i*n+k+1]=q;
A[i*n+k]=p;
}
}
k++;
}
}
}
delete []A;
return true;
}
#endif //MATRIX_MATRIXSOLVER_H

View File

@@ -0,0 +1,57 @@
//
// Created by 顾涵彬 on 2019-08-30.
//
#ifndef MATRIX_QUATERNION_H
#define MATRIX_QUATERNION_H
#include "SquareMatrix.h"
#include <cmath>
namespace Ctain {
// using SMatrix<double>;
template<typename T>
class Quaternion {
public:
Quaternion(){}
Quaternion(SMatrix<double> m){
//double f = 1e-10;
_w = sqrt(m(0)+m(4)+m(8)+1)/2;
_x = (m(5)-m(7))/(4*_w);
_y = (m(6)-m(2))/(4*_w);
_z = (m(1)-m(3))/(4*_w);
}
Quaternion(T X,T Y,T Z,T W) : _x(X),_y(Y),_z(Z),_w(W){}
void normalize() {
double len;
len = sqrt(_x*_x+_y*_y+_z*_z+_w*_w);
_x = _x / len;
_y = _y / len;
_z = _z / len;
_w = _w / len;
}
T x(){return _x;}
T y(){return _y;}
T z(){return _z;}
T w(){return _w;}
SMatrix<double> toRotationMatrix() const {
SMatrix<double> r(3);
double q0=_w,q1=_x,q2=_y,q3=_z;
r(0) = 1 - 2*q2*q2-2*q3*q3;
r(1) = 2*q1*q2+2*q0*q3;
r(2) = 2*q1*q3-2*q0*q2;
r(3) = 2*q1*q2 - 2*q0*q3;
r(4) = 1-2*q1*q1 - 2*q3*q3;
r(5) = 2*q2*q3 +2*q0*q1;
r(6) = 2*q1*q3+2*q0*q2;
r(7) = 2*q2*q3 - 2*q0*q1;
r(8) = 1-2*q1*q1 -2*q2*q2;
return r;
}
private:
T _x;
T _y;
T _z;
T _w;
};
}
#endif //MATRIX_QUATERNION_H

View File

@@ -0,0 +1,91 @@
//
// Created by 顾涵彬 on 2019-08-29.
//
#ifndef MATRIX_SQUAREMATRIX_H
#define MATRIX_SQUAREMATRIX_H
#include "Matrix.h"
namespace Ctain {
#define Matrix Matrix<_Scalar>
template<typename _Scalar>
class SMatrix: public Matrix{
public:
SMatrix(int D) : Matrix(D, D) {}
SMatrix() : Matrix(0, 0) {}
SMatrix(_Scalar _data[], int D) :
Matrix(_data, D, D) {}
SMatrix(_Scalar **_data, int D) :
Matrix(_data, D, D) {}
SMatrix(Matrix m) :
Matrix(m) {}
// void operator =(const Matrix &m){
// }
_Scalar determinant();
_Scalar M(int m, int n);
SMatrix<_Scalar> inverse() {
SMatrix<_Scalar> res(Matrix::_Rows);
_Scalar d = determinant();
for (int i = 0; i < Matrix::_Rows; i++) {
for (int j = 0; j < Matrix::_Cols; j++) {
res.Data(j, i) = 1.0*M(i, j)/d;
}
}
return res;
}
};//class Matrix end
template<typename _Scalar>
_Scalar SMatrix<_Scalar>::determinant() {
int r, c, m;
int lop = 0;
int n = Matrix::_Rows;
_Scalar result = 0;
_Scalar mid = 1;
if (n != 1) {
lop = (n == 2) ? 1 : n;
for (m = 0; m < lop; m++) {
mid = 1;
for (r = 0, c = m; r < n; r++, c++) {
mid = mid * (*(Matrix::data+r*n+c%n));
}
result += mid;
}
for (m = 0; m < lop; m++) {
mid = 1;
for (r = 0, c = n-1-m+n; r < n; r++, c--) {
mid = mid * (*(Matrix::data+r*n+c%n));
}
result -= mid;
}
}
else
result = Matrix::data[0];
return result;
}
template<typename _Scalar>
_Scalar SMatrix<_Scalar>::M(int m, int n) {
float mid_result = 0;
int sign = 1;
int k = Matrix::_Rows;
SMatrix mid(k-1);
int c = 0;
for (int i = 0; i < k; i++) {
for (int j = 0; j < k; j++) {
if (i != m && j != n)
{
mid.Data(c++) = Matrix::cData(i,j);
}
}
}
sign = (m+n)%2 == 0 ? 1 : -1;
mid_result = (float)sign*mid.determinant();
return mid_result;
}
#undef Matrix
}//namespace Ctain end
#endif //MATRIX_SQUAREMATRIX_H

View File

@@ -16,7 +16,7 @@
#include <vector>
#include <memory>
#include "Ctain/CtainBase.h"
#include "CtainBase.h"
#include <opencv2/core/core.hpp>
#include "mynteye/mynteye.h"