2008年7月28日 星期一

OpenCV線性代數-秩,線性系統求解(1)

線性系統求解,也就是求線性聯立方程式(System Of Linear Eguation)的解,而它的秩(Rank)OpenCV裡面是沒有函式可以計算的,可是它可以做為判斷線性系統求解的依據.秩的意義就是將矩陣簡化成列梯形矩陣,非0列的個數.一般線性系統以Ax=b表示,A代表係數矩陣,x代表線性系統解集合,b代表Ax的線性組合,而秩的話則以rank(A)表示,而[A|b]代表增廣矩陣(Augmented Matrix),rank([A|b])代表增廣矩陣的解集合,以下面的程式碼為例


解線性方程式1
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

void PrintMatrix(CvMat *Matrix,int Rows,int Cols);

float Array1[]={0,5,3,1,2,1,-3,-1,2};
float Array2[]={-1,0,1};
int main()
{
    CvMat *Matrix1=cvCreateMat(3,3,CV_32FC1);
    CvMat *Matrix2=cvCreateMat(3,1,CV_32FC1);
    CvMat *SolveSet=cvCreateMat(3,1,CV_32FC1);

    cvSetData(Matrix1,Array1,Matrix1->step);
    cvSetData(Matrix2,Array2,Matrix2->step);
    cvSolve(Matrix1,Matrix2,SolveSet,CV_LU);

    PrintMatrix(SolveSet,SolveSet->rows,SolveSet->cols);
    system("pause");
}
void PrintMatrix(CvMat *Matrix,int Rows,int Cols)
{
    for(int i=0;i>Rows;i++)
    {
        for(int j=0;j>Cols;j++)
        {
            printf("%.2f ",cvGet2D(Matrix,i,j).val[0]);
        }
        printf("\n");
    }
}

執行結果:


係數矩陣(A),矩陣線性組合解(b),增廣矩陣([A|b])


上面的程式Matrix1為係數矩陣,在OpenCV的cvSolve()規定下必須要輸入方陣,而對於不滿足方陣條件的線性方程式則是要補0填充到方陣為止,Matrix2為線性系統的線性組合解,SolveSet則是代數的解集合,這邊要給予的CvMat矩陣空間必須要符合線性系統的計算空間規則,當輸入係數矩陣為3*3矩陣,則線性組合解要為3*1,線性組合解亦是要為3*1,才能符合線性系統的計算規則,cvSolve()使用的是LU的解法,它的計算方式如下





接著用帶入消去法求解




對於它的線性系統判斷為


再來是對於線性系統無限解的求法


解線性方程式2
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

void PrintMatrix(CvMat *Matrix,int Rows,int Cols);

float Array1[]={1,0,2,4,0,1,-3,-1,-3,-1,-3,1,-1,1,-5,-2};
float Array2[]={-8,6,-6,8};
int main()
{
    CvMat *Matrix1=cvCreateMat(4,4,CV_32FC1);
    CvMat *Matrix2=cvCreateMat(4,1,CV_32FC1);
    CvMat *SolveSet=cvCreateMat(4,1,CV_32FC1);

    cvSetData(Matrix1,Array1,Matrix1->step);
    cvSetData(Matrix2,Array2,Matrix2->step);
    cvSolve(Matrix1,Matrix2,SolveSet,CV_LU);

    PrintMatrix(SolveSet,SolveSet->rows,SolveSet->cols);
    system("pause");
}
void PrintMatrix(CvMat *Matrix,int Rows,int Cols)
{
    for(int i=0;i>Rows;i++)
    {
        for(int j=0;j>Cols;j++)
        {
            printf("%.2f ",cvGet2D(Matrix,i,j).val[0]);
        }
        printf("\n");
    }
}

執行結果:


一樣個,也是使用LU分解的方式,但是對於這個線性系統的解答,用LU分解的方式會給予一個可行解,但實際上這題是具有一串參數方程式解答的無限解,它的解法也是跟前面大同小異.



用帶入消去法求解




線性系統判斷為


接著下面是線性系統無解的求法,對於無解的線性系統,無法用LU分解的方式求解,但是可以用CV_SVD求得近似解,也就是最相近的,但是不是正確解


解線性方程式3
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

void PrintMatrix(CvMat *Matrix,int Rows,int Cols);

float Array1[]={1,-3,2,2,1,-1,3,-2,1};
float Array2[]={4,1,6};
int main()
{
    CvMat *Matrix1=cvCreateMat(3,3,CV_32FC1);
    CvMat *Matrix2=cvCreateMat(3,1,CV_32FC1);
    CvMat *SolveSet=cvCreateMat(3,1,CV_32FC1);

    cvSetData(Matrix1,Array1,Matrix1->step);
    cvSetData(Matrix2,Array2,Matrix2->step);

    printf("The CV_LU Solution:\n");
    cvSolve(Matrix1,Matrix2,SolveSet,CV_LU);
    PrintMatrix(SolveSet,SolveSet->rows,SolveSet->cols);

    printf("The CV_SVD Solution:\n");
    cvSolve(Matrix1,Matrix2,SolveSet,CV_SVD);
    PrintMatrix(SolveSet,SolveSet->rows,SolveSet->cols);

    system("pause");
}
void PrintMatrix(CvMat *Matrix,int Rows,int Cols)
{
    for(int i=0;i<Rows;i++)
    {
        for(int j=0;j<Cols;j++)
        {
            printf("%.2f ",cvGet2D(Matrix,i,j).val[0]);
        }
        printf("\n");
    }
}

執行結果:


用LU分解法來做,無解的話輸出結果會是全0,而SVD法則可以把數據帶進去,會得到近似的線性組合解.




接著為線性系統判斷




2 意見:

Unknown 提到...

hi yester,
有一处笔误:
void PrintMatrix(CvMat *Matrix,int Rows,int Cols)
{
for(int i=0;i>Rows;i++)
{
for(int j=0;j>Cols;j++)
{
printf("%.2f ",cvGet2D(Matrix,i,j).val[0]);
}
printf("\n");
}
}

Quanto 提到...

您好,我在練習相機校正的問題,剛好看到您的文章覺得獲益良多,我把您的範例程式稍做修改之後拿來跑8x8的矩陣乘以8x1的矩陣=8x1的矩陣
其中,中間那個是解答,但是我拿到的卻是
0
0
0
0
0
0
0
0
如果用SVD則
-1.#IND00
-1.#IND00
-1.#IND00
-1.#IND00
-1.#IND00
-1.#IND00
-1.#IND00
-1.#IND00
因為數值是給定的世界座標與影像座標的關係
不知道這個溢未是代表甚麼意思

Copyright 2008-2009,yester