标签：orm [] continue maps strong tran tail get 消元
最近在搞图像处理，碰到了透视变换的问题。
同事给我一些代码，里边有误差，挺严重，让我帮他想想哪里出错了。捣鼓了很久，我猜测肯定是透视变换矩阵求错了，然后我的透视变换之旅就开始了。
后来问题解决了，发现他的矩阵和我求得矩阵一摸一样，他的代码并没有错误，是他采用的图片在做广角变换的时候有误差，导致程序结果误差。
首先感谢xiaowei_cqu，是她的一篇博客（http://blog.csdn.net/xiaowei_cqu/article/details/26471527）教会了我变换原理。
后来去用opencv来验证，发现函数调来调去太麻烦了，代码量也不小。
可是我只想要变换矩阵，自己写个不就好了！然后根据opencv源码里的原理，配上高斯消元自己写了个qwq。
talk is cheap,show me the code!!!
//Calculates coefficients of perspective transformation
/* Calculates coefficients of perspective transformation
* which maps (xi,yi) to (ui,vi), (i=1,2,3,4):
*
* c00*xi + c01*yi + c02
* ui = ---------------------
* c20*xi + c21*yi + c22
*
* c10*xi + c11*yi + c12
* vi = ---------------------
* c20*xi + c21*yi + c22
*
* Coefficients are calculated by solving linear system:
* / x0 y0 1 0 0 0 -x0*u0 -y0*u0 \ /c00\ /u0* | x1 y1 1 0 0 0 -x1*u1 -y1*u1 | |c01| |u1|
* | x2 y2 1 0 0 0 -x2*u2 -y2*u2 | |c02| |u2|
* | x3 y3 1 0 0 0 -x3*u3 -y3*u3 |.|c10|=|u3|,
* | 0 0 0 x0 y0 1 -x0*v0 -y0*v0 | |c11| |v0|
* | 0 0 0 x1 y1 1 -x1*v1 -y1*v1 | |c12| |v1|
* | 0 0 0 x2 y2 1 -x2*v2 -y2*v2 | |c20| |v2|
* \ 0 0 0 x3 y3 1 -x3*v3 -y3*v3 / \c21/ \v3/
*
* where:
* cij - matrix coefficients, c22 = 1
*/
void Gauss(double A[][9], int equ, int var, double* ans) { //epu:A‘s row var:A‘s col-1
int row, col;
for (row = 0, col = 0; col<var&&row<equ; col++, row++) {
int max_r = row;
for (int i = row + 1; i<equ; i++) {
if ((1e-12)<fabs(A[i][col]) - fabs(A[max_r][col])) {
max_r = i;
}
}
if (max_r != row)
for (int j = 0; j<var + 1; j++)
swap(A[row][j], A[max_r][j]);
for (int i = row + 1; i<equ; i++) {
if (fabs(A[i][col])<(1e-12))
continue;
double tmp = -A[i][col] / A[row][col];
for (int j = col; j<var + 1; j++) {
A[i][j] += tmp*A[row][j];
}
}
}
for (int i = var - 1; i >= 0; i--) { //计算唯一解。
double tmp = 0;
for (int j = i + 1; j<var; j++) {
tmp += A[i][j] * (*(ans + j));
}
ans[i] = (A[i][var] - tmp) / A[i][i];
}
}
void byx_getPerspectiveTransform(PointD * src, PointD * dst, double* ret){
double x0 = src[0].x, x1 = src[1].x, x2 = src[3].x, x3 = src[2].x;
double y0 = src[0].y, y1 = src[1].y, y2 = src[3].y, y3 = src[2].y;
double u0 = dst[0].x, u1 = dst[1].x, u2 = dst[3].x, u3 = dst[2].x;
double v0 = dst[0].y, v1 = dst[1].y, v2 = dst[3].y, v3 = dst[2].y;
double A[8][9] = {
{ x0, y0, 1, 0, 0, 0, -x0*u0, -y0*u0, u0 },
{ x1, y1, 1, 0, 0, 0, -x1*u1, -y1*u1, u1 },
{ x2, y2, 1, 0, 0, 0, -x2*u2, -y2*u2, u2 },
{ x3, y3, 1, 0, 0, 0, -x3*u3, -y3*u3, u3 },
{ 0, 0, 0, x0, y0, 1, -x0*v0, -y0*v0, v0 },
{ 0, 0, 0, x1, y1, 1, -x1*v1, -y1*v1, v1 },
{ 0, 0, 0, x2, y2, 1, -x2*v2, -y2*v2, v2 },
{ 0, 0, 0, x3, y3, 1, -x3*v3, -y3*v3, v3 },
};
double ans[8] = { 0 };
Gauss(A, 8, 8, ans);
*(ret + 0) = ans[0]; *(ret + 1) = ans[1]; *(ret + 2) = 0; *(ret + 3) = ans[2];
*(ret + 4) = ans[3]; *(ret + 5) = ans[4]; *(ret + 6) = 0; *(ret + 7) = ans[5];
*(ret + 8) = 0; *(ret + 9) = 0; *(ret + 10) = 1; *(ret + 11) = 0;
*(ret + 12) = ans[6]; *(ret + 13) = ans[7]; *(ret + 14) = 0; *(ret + 15) = 1;
}
PointD byx_Transform(PointD p,double * mat){
PointD ret;
double D = p.x*mat[12]+p.y*mat[13]+mat[15];
ret.x=(p.x*mat[0]+p.y*mat[1]+mat[3])/D;
ret.y=(p.x*mat[4]+p.y*mat[5]+mat[7])/D;
return ret;
}
ps:变换矩阵是3*3的，我同事需要4*4的（第三行和第三列自行忽略）
标签：orm [] continue maps strong tran tail get 消元
原文地址：http://www.cnblogs.com/baoyinhang/p/7202200.html
