昨日

OpenCV統計應用-Mahalanobis距離

2009-02-03T18:59:00.015+08:00

Mahalanobis距離是一個可以準確找出資料分布上面極端值(Outliers)的統計方法,使用線性迴歸的概念,也就是說他使用的是共變數矩陣以及該資料分布的平均數來找尋極端值的產生,而可以讓一群資料系統具有穩健性(Robust),去除不必要的雜訊訊息,這邊拿前面共變數矩陣的資料為例,並且新增了兩個點座標向量來做Mahalanobis距離的比較

加入兩個座標點找尋極端值

Mahalanobis距離實作
#include <cv.h>
#include <stdio.h>
#include <stdlib.h>

float Coordinates[20]={1.5,2.3,
                                 3.0,1.7,
                                 1.2,2.9,
                                 2.1,2.2,
                                 3.1,3.1,
                                 1.3,2.7,
                                 2.0,1.7,
                                 1.0,2.0,
                                 0.5,0.6,
                                 1.0,0.9};

float Coordinates2[2]={1.3,1.5};
float Coordinates3[2]={200,100};

int main()
{
    CvMat *Vector[1];
    CvMat *Vector1;
    CvMat *CovarMatrix;
    CvMat *InvertCovarMatrix;
    CvMat *AvgVector;
    CvMat *Vector2;
    CvMat *Vector3;

    Vector1=cvCreateMat(10,2,CV_32FC1);
    cvSetData(Vector1,Coordinates1,Vector1->step);
    Vector[0]=Vector1;
    CovarMatrix=cvCreateMat(2,2,CV_32FC1);
    InvertCovarMatrix=cvCreateMat(2,2,CV_32FC1);

    AvgVector=cvCreateMat(1,2,CV_32FC1);
    Vector2=cvCreateMat(1,2,CV_32FC1);
    Vector3=cvCreateMat(1,2,CV_32FC1);
    cvSetData(Vector2,Coordinates2,Vector2->step);
    cvSetData(Vector3,Coordinates3,Vector3->step);

    cvCalcCovarMatrix((const CvArr **)Vector,10,CovarMatrix,AvgVector,CV_COVAR_SCALE+CV_COVAR_NORMAL+CV_COVAR_ROWS);
    cvInvert(CovarMatrix,InvertCovarMatrix,CV_SVD_SYM);

    printf("\nVector2 Mahalanobis Distance\n");
    printf("%f\n",cvMahalanobis(Vector2,AvgVector,InvertCovarMatrix));

    printf("\nVector3 Mahalanobis Distance\n");
    printf("%f\n",cvMahalanobis(Vector3,AvgVector,InvertCovarMatrix));

    system("pause");
}

由上面可以看得出來第一筆座標向量Vector2,他的座標位在(1.3,1.5)的位置,第二筆座標向量Vector2他的座標在(200,100),對於原始資料分布來說,第二筆座標向量它已經遠遠離開了這些座標分佈,因此可以從執行結果看的出來Vector3的Mahalanobis距離遠大於Vector2,而在做Mahalanobis距離的同時,也必須要將共變數矩陣做反矩陣運算,cvMahalanobis()第一個引數為目標要判斷是否是極端值的向量CvMat資料結構,第二個引數為該資料分佈的平均數向量或是該資料分佈的某筆資料的CvMat資料結構,第三個引數為該資料分佈共變數矩陣的反矩陣

而Mahalanobis距離的原理如下

當資料分佈使的共變數矩陣為單位向量I的話,就會退化成歐機理得距離(Euclidean Distance)

而他的計算方式如下

由上面可知,在一個簡單線性迴歸的資料模型下,第二筆資料Vector3為這個線性迴歸的極端值

OpenCV統計應用-PCA主成分分析

2009-01-31T13:25:00.025+08:00

在圖形識別方面,主成分分析(Principal Comonents Analysis,PCA)算是比較快速而且又準確的方式之一,它可以對抗圖形平移旋轉的事件發生,並且藉由主要特徵(主成分)投影過後的資料做資料的比對,在多個特徵資訊裡面,取最主要的K個,做為它的特徵依據,在這邊拿前面共變數矩陣的數據來做沿用,主成分分析使用的方法為計算共變數矩陣,在加上計算共變數矩陣的特徵值及特徵向量,將特徵值以及所對應的特徵向量排序之後,取前面主要K個特徵向量當做主要特徵,而OpenCV也可以對高維度的向量進行主成分分析的計算

數據原始的分佈情況,紅點代表著它的平均數

將座標系位移,以紅點為主要的原點

計算共變數矩陣以及共變數的特徵值以及特徵向量,將特徵向量排序後投影回原始資料的結果的結果,也就是說,對照上面的圖片,EigenVector的作用是找到主軸後,將原本的座標系做旋轉了

再來就是對它做投影,也就是降低維度的動作,將Y軸的數據全部歸零,投影在X軸上

投影完之後,在將它轉回原本的座標系

PCA主成分分析,與線性迴歸有異曲同工之妙,也就是說,這條投影過後的直線,可以稱做迴歸線,當它在做主軸旋轉的時候,所投影的結果為最小均方誤,在將它轉置回來的時候,就形成了一條迴歸直線了

OpenCV的PCA輸入必須要是單通道32位元浮點數格式或是單通道64位元浮點數格式的,參數為CV_32FC1或是CV_64FC1,程式寫法如下

PCA程式實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
#include <stdlib.h>

float Coordinates[20]={1.5,2.3,
                                 3.0,1.7,
                                 1.2,2.9,
                                 2.1,2.2,
                                 3.1,3.1,
                                 1.3,2.7,
                                 2.0,1.7,
                                 1.0,2.0,
                                 0.5,0.6,
                                 1.0,0.9};

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

int main()
{
    CvMat *Vector1;
    CvMat *AvgVector;
    CvMat *EigenValue_Row;
    CvMat *EigenVector;

    Vector1=cvCreateMat(10,2,CV_32FC1);
    cvSetData(Vector1,Coordinates,Vector1->step);
    AvgVector=cvCreateMat(1,2,CV_32FC1);
    EigenValue_Row=cvCreateMat(2,1,CV_32FC1);
    EigenVector=cvCreateMat(2,2,CV_32FC1);

    cvCalcPCA(Vector1,AvgVector,EigenValue_Row,EigenVector,CV_PCA_DATA_AS_ROW);

    printf("Original Data:\n");
    PrintMatrix(Vector1,10,2);

    printf("==========\n");
    PrintMatrix(AvgVector,1,2);

    printf("\nEigne Value:\n");
    PrintMatrix(EigenValue_Row,2,1);

    printf("\nEigne Vector:\n");
    PrintMatrix(EigenVector,2,2);

    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");
    }
}

執行結果:

這部份是把平均數,共變數矩陣,以及特徵值及特徵向量都計算出來了,全部都包在cvCalcPCA()的函式裡面,因此可以不必特地的用cvCalcCovarMatrix()求得共變數矩陣,也不需要再由共變數矩陣套用cvEigenVV()求得它的EigenValue以及EigenVector了,所以說,cvCalcPCA()=cvCalcCovarMatrix()+cvEigenVV(),不僅如此,cvCalcPCA()使用上更是靈活,當向量的維度數目比輸入的資料那的時候(例如Eigenface),它的共變數矩陣就會自動轉成CV_COVAR_SCRAMBLED,而當輸入資料量比向量維度大的時候,它亦會自動轉成CV_COVAR_NORMAL的形態,而OpenCV也提供了計算投影量cvProjectPCA(),以及反向投影的函式cvBackProjectPCA(),cvCalcPCA()的計算結果如下

詳細計算方法可以參考"OpenCV統計應用-共變數矩陣"以及"OpenCV線性代數-cvEigenVV實作"這兩篇,cvCalcPCA()第一個引數為輸入目標要計算的向量,整合在CvMat資料結構裡,第二個引數為空的平均數向量,第三個引數為輸出排序後的EigenValue,以列(Rows)為主的數值,第四個引數為排序後的EigenVector,第五個引數為cvCalcPCA()的參數,它的參數公式如下

#define CV_PCA_DATA_AS_ROW 0
#define CV_PCA_DATA_AS_COL 1
#define CV_PCA_USE_AVG 2

分別代表以列為主,以欄為主的參數設定,以及使用自己定義的平均數,CV_PCA_DATA_AS_ROW與CV_PCA_DATA_AS_COL的參數不可以同時使用,而對於主成分分析的EigenVector對原始資料投影的程式範例如下

PCA特徵向量投影
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
#include <stdlib.h>

float Coordinates[20]={1.5,2.3,
                                 3.0,1.7,
                                 1.2,2.9,
                                 2.1,2.2,
                                 3.1,3.1,
                                 1.3,2.7,
                                 2.0,1.7,
                                 1.0,2.0,
                                 0.5,0.6,
                                 1.0,0.9};

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

int main()
{
    CvMat *Vector1;
    CvMat *AvgVector;
    CvMat *EigenValue_Row;
    CvMat *EigenVector;

    Vector1=cvCreateMat(10,2,CV_32FC1);
    cvSetData(Vector1,Coordinates,Vector1->step);
    AvgVector=cvCreateMat(1,2,CV_32FC1);
    EigenValue_Row=cvCreateMat(2,1,CV_32FC1);
    EigenVector=cvCreateMat(2,2,CV_32FC1);

    cvCalcPCA(Vector1,AvgVector,EigenValue_Row,EigenVector,CV_PCA_DATA_AS_ROW);
    cvProjectPCA(Vector1,AvgVector,EigenVector,Vector1);

    printf("Project Original Data:\n");
    PrintMatrix(Vector1,10,2);

    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");
    }
}

執行結果:

cvProjectPCA(),它的公式定義如下

因此,它所呈現的結果就是將座標旋轉後的特徵向量投影,cvProjectPCA()的地一個引數為輸入原始向量資料,第二個引數為輸入空的(或是以設定好的)平均數向量,第三個引數為輸入以排序好的特徵向量EigenVector,第四個引數為輸出目標投影向量,而投影之外,OpenCV還可以給他做反向投影回原始的資料

PCA反向投影
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
#include <stdlib.h>

float Coordinates[20]={1.5,2.3,
                                 3.0,1.7,
                                 1.2,2.9,
                                 2.1,2.2,
                                 3.1,3.1,
                                 1.3,2.7,
                                 2.0,1.7,
                                 1.0,2.0,
                                 0.5,0.6,
                                 1.0,0.9};

int main()
{
    CvMat *Vector1;
    CvMat *AvgVector;
    IplImage *Image1=cvCreateImage(cvSize(450,450),IPL_DEPTH_8U,3);
    Image1->origin=1;

    Vector1=cvCreateMat(10,2,CV_32FC1);
    cvSetData(Vector1,Coordinates,Vector1->step);

    AvgVector=cvCreateMat(1,2,CV_32FC1);

    CvMat *EigenValue_Row=cvCreateMat(2,1,CV_32FC1);
    CvMat *EigenVector=cvCreateMat(2,2,CV_32FC1);

    cvCalcPCA(Vector1,AvgVector,EigenValue_Row,EigenVector,CV_PCA_DATA_AS_ROW);
    cvProjectPCA(Vector1,AvgVector,EigenVector,Vector1);
    cvBackProjectPCA(Vector1,AvgVector,EigenVector,Vector1);

    printf("Back Project Original Data:\n");
    for(int i=0;i<10;i++)
    {
        printf("%.2f ",cvGetReal2D(Vector1,i,0));
       printf("%.2f ",cvGetReal2D(Vector1,i,1));

        cvCircle(Image1,cvPoint((int)(cvGetReal2D(Vector1,i,0)*100),(int)(cvGetReal2D(Vector1,i,1)*100)),0,CV_RGB(0,0,255),10,CV_AA,0);

        printf("\n");
    }
    cvCircle(Image1,cvPoint((int)((cvGetReal2D(AvgVector,0,0))*100),(int)((cvGetReal2D(AvgVector,0,1))*100)),0,CV_RGB(255,0,0),10,CV_AA,0);

    printf("==========\n");
    printf("%.2f ",cvGetReal2D(AvgVector,0,0));
    printf("%.2f ",cvGetReal2D(AvgVector,0,1));

    cvNamedWindow("Coordinates",1);
    cvShowImage("Coordinates",Image1);
    cvWaitKey(0);
}

執行結果:

而反向投影,只不過是將投影向量還原成原始資料罷了,可以用來做為新進資料反向投影後用來比對的步驟,以下是它的計算公式推導

cvBackProjectPCA()的引數輸入與cvProjectPCA()是一樣的,只不過是裡面的計算公式不同,其實也只是cvProjectPCA()的反運算

cvCalcPCA()
計算多筆高維度資料的主要特徵值及特徵向量,先自動判斷維度大於資料量或維度小於資料量來選擇共變數矩陣的樣式,在由共變數矩陣求得已排序後的特徵值及特徵向量,cvCalcPCA()函式的參數為CV_PCA_DATA_AS_ROW以列為主的資料排列,CV_PCA_DATA_AS_COL以行為主的資料排列,CV_PCA_USE_AVG自行定義平均數數據,CV_PCA_DATA_AS_ROW以及CV_PCA_DATA_AS_COL不可以同時使用,cvCalcPCA()的第一個引數為輸入CvMat維度向量資料結構,第二個引數為輸入空的CvMat平均數向量資料結構(或是自行定義平均數向量),第三個引數為輸入以列為主的特徵值CvMat資料結構,第四個引數為輸入特徵向量CvMat資料結構,第四個引數為輸入cvCalcPCA()函式的參數
cvCalcPCA(輸入目標向量資料CvMat資料結構,輸入或輸出向量平均數CvMat資料結構,輸入以列為主EigenValue的CvMat資料結構,輸入EigenVector的CvMat資料結構,目標參數或代號)

cvProjectPCA()
將計算主成分分析的結果做投影的運算,主要作用是將最具意義的k個EgienVector與位移後的原始資料做矩陣乘法,投影過後的結果將會降低維度,cvProjectPCA()第一個引數為輸入CvMat資料結構原始向量資料數據,第二個引數為輸入CvMat資料結構平均數向量,第三個引數為輸入要降維的特徵向量,第四個引數為輸出CvMat資料結構原始資料投影的結果
cvProjectPCA(輸入CvMat原始向量數據資料結構,輸入CvMat原始向量平均數資料結構,輸入CvMat降低維度的特徵向量資料結構,輸出CvMat投影向量資料結構)

cvBackProjectPCA()
將投影向量轉回原始向量維度座標系,cvProjectPCA()第一個引數為輸入CvMat資料結構投影向量數據,第二個引數為輸入CvMat資料結構平均數向量,第三個引數為輸入特徵向量,第四個引數為輸出CvMat資料結構反向投影的結果
cvProjectPCA(輸入CvMat投影向量資料結構,輸入CvMat原始向量平均數資料結構,輸入CvMat特徵向量資料結構,輸出CvMat反向投影資料結構)

OpenCV統計應用-共變數矩陣

2009-01-23T19:15:00.031+08:00

共變數(Covariance),為兩個隨機變數的離均差除以母體個數,可以判斷兩個事件隨機變數的相依性,而單一變數的共變數,那就是變異數(Variance)了,共變數以及變異數簡單的定義如下

變異數

共變數

而當兩隨機變數為獨立事件的時候

再來提到的是共變數矩陣(Covariance matrix),代表著隨機變數所有共變數的種類,以兩組隨機變數來說所產生的共變數矩陣如下

而三組,四組以上隨機變數的共變數矩陣又是不同情形了,在這邊,隨機變數的個數代表著向量的維度,而他的數對則是同時發生的情形,下面就以座標的簡單例子來示範

這是個維度為2的座標數組,代表的是一個座標平面的點集合,而下面就是給它跑cvCalcCovarMatrix()共變數矩陣的計算方法

共變數矩陣1
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
#include <stdlib.h>

float Coordinates[20]={1.5,2.3,
                                 3.0,1.7,
                                 1.2,2.9,
                                 2.1,2.2,
                                 3.1,3.1,
                                 1.3,2.7,
                                 2.0,1.7,
                                 1.0,2.0,
                                 0.5,0.6,
                                 1.0,0.9};

int main()
{
    CvMat *Vector[10];
    CvMat *CovarMatrix;
    CvMat *AvgMatrix;
    IplImage *Image1=cvCreateImage(cvSize(450,450),IPL_DEPTH_8U,3);
    Image1->origin=1;
    for(int i=0;i<10;i++)
    {
        Vector[i]=cvCreateMat(1,2,CV_32FC1);
        cvSetReal1D(Vector[i],0,Coordinates[i*2]);
        cvSetReal1D(Vector[i],1,Coordinates[i*2+1]);
        cvCircle(Image1,cvPoint((int)(Coordinates[i*2]*100),(int)(Coordinates[i*2+1]*100)),0,CV_RGB(0,0,255),10,CV_AA,0);
    }

    CovarMatrix=cvCreateMat(2,2,CV_32FC1);
    AvgMatrix=cvCreateMat(1,2,CV_32FC1);
    cvCalcCovarMatrix((const CvArr **)Vector,10,CovarMatrix,AvgMatrix,CV_COVAR_SCALE+CV_COVAR_NORMAL);

    for(int i=0;i<2;i++)
    {
        for(int j=0;j<2;j++)
        {
            printf("%f ",cvGetReal2D(CovarMatrix,i,j));
        }
        printf("\n");
    }
    cvNamedWindow("Coordinates",1);
    cvShowImage("Coordinates",Image1);
    cvWaitKey(0);
}

執行結果:

計算方法:

上面的方法是用個別的向量來實作,向量的維度為2,而圖片顯示的結果為它們二維座標的分佈情況,cvCalcCovarMatrix()必須給它空的平均值向量矩陣來計算,而cvCalcCovarMatrix()它的參數被定義如下

#define CV_COVAR_SCRAMBLED 0
#define CV_COVAR_NORMAL 1
#define CV_COVAR_USE_AVG 2
#define CV_COVAR_SCALE 4
#define CV_COVAR_ROWS 8
#define CV_COVAR_COLS 16

由上面可以知道,它是由2的N次方所表達的,所以可以用合成參數的方式來表達cvCalcCovarMatrix()的共變數矩陣函式,也就是說,可以用CV_COVAR_SCRAMBLED+CV_COVAR_ROWS的組合,也可以用CV_COVAR_NORMAL+CV_COVAR_USE_AVG+CV_COVAR_ROWS之類的組合,而它所代表的含意分別是

CV_COVAR_SCRAMBLED
一種共變數矩陣的計算方式,不可以與CV_COVAR_NORMAL合用,表達方式如下

這計算方式在OpenCV說明文件提到為Eigenface的計算方式

CV_COVAR_NORMAL
共變數矩陣計算方式,不可與CV_COVAR_SCRAMBLED合用,表達方式如下

這個則是為一般共變數矩陣的計算方式

CV_COVAR_USE_AVG
不用共變數矩陣內建計算平均數的函式,而是自己給予平均數值,可與任何參數共用

CV_COVAR_SCALE
對共變數矩陣的數據做向量個數總和的純量積,用的是除法計算如下的共變數矩陣

可與任何參數共用,而如果沒這參數則是無除以N的計算

CV_COVAR_ROWS
將共變數矩陣的輸入值用矩陣的方式表達,而不是用向量的方式,矩陣表達方式以列(Rows)為主,並且參數不可與CV_COVAR_COLS合用

CV_COVAR_COLS
將共變數矩陣的輸入值用矩陣的方式表達,而不是用向量的方式,矩陣表達方式以欄(Columns)為主,並且參數不可與CV_COVAR_ROWS共用

共變數矩陣有幾個規則,也就是,輸入一定要是方陣,平均數的長度要是向量的維度,而平均數的長度也一定要是向量的大小,然後一定要是用單通道CV_32FC1或是CV_64FC1做為輸入,cvCalcCovarMatrix()函式,第一個引數則必須要用(const CvArr **)強制型別轉換,第二個引數為輸入向量的數目,第三個引數為空的或是非空平均數向量,第四個引數為cvCalcCovarMatrix()這函式要輸入的參數,而下面,則是使用矩陣方式表達共變數矩陣的範例

共變數矩陣2
#include <cv.h>
#include <stdio.h>
#include <stdlib.h>

float Coordinates[20]={1.5,2.3,
                                 3.0,1.7,
                                 1.2,2.9,
                                 2.1,2.2,
                                 3.1,3.1,
                                 1.3,2.7,
                                 2.0,1.7,
                                 1.0,2.0,
                                 0.5,0.6,
                                 1.0,0.9};

int main()
{
    CvMat *Vector[1];
    CvMat *Vector1;
    CvMat *CovarMatrix;
    CvMat *avg;

    Vector1=cvCreateMat(10,2,CV_32FC1);
    cvSetData(Vector1,Coordinates,Vector1->step);
    Vector[0]=Vector1;
    CovarMatrix=cvCreateMat(2,2,CV_32FC1);
    avg=cvCreateMat(1,2,CV_32FC1);

    cvCalcCovarMatrix((const CvArr **)Vector,10,CovarMatrix,avg,CV_COVAR_SCALE+CV_COVAR_NORMAL+CV_COVAR_ROWS);

    for(int i=0;i<2;i++)
    {
        for(int j=0;j<2;j++)
        {
            printf("%f ",cvGetReal2D(CovarMatrix,i,j));
        }
        printf("\n");
    }
    system("pause");
}

執行結果:

上面的程式碼,如果想用矩陣來表達向量的方式計算共變數矩陣,就必須要將第一個引數設為二維陣列當做輸入,輸入的方式就如上面程式的寫法,而其他地方則是沒什麼差異.

共變數矩陣在OpenCV內,不但可以做到計算主成分分析(Principal Cmponents Analysis,PCA),以及另一個Mahalanobis距離的計算

cvCalcCovarMatrix()
計算共變數矩陣,輸入可為多個IplImage或CvMat資料結構,可輸入高維度資料的向量,有多種輸入方式,在資料輸入方面,可以用單一CvMat資料結構,以列(Rows)為主或以行(Columns)為主,使用CV_COVAR_ROWS以及CV_COVAR_COLS的參數輸入,而要做多個IplImage或CvMat資料結構輸入則不需提供CV_COVAR_ROWS,CV_COVAR_COLS的參數,在計算方面,又分為以維度為主的共變數矩陣參數輸入CV_COVAR_NORMAL以及以個數為主的共變數矩陣CV_COVAR_SCRAMBLED,而他也可一自行定義平均值CV_COVAR_USE_AVG,以及除以是否除以共變數矩陣的共變數個數CV_COVAR_SCALE,而這些參數可以自行組合運算,第一個引數為輸入目標向量,第二個引數為輸出目標共變數矩陣,第三個引數為輸入或輸出平均數向量,第四個引數為cvCalcCovarMatrix()共變數計算函式的參數輸入
cvCalcCovarMatrix(輸入多個IplImage或CvMat資料結構,輸出目標共變數矩陣,輸入/出平均數向量,共變數矩陣參數或代號)

OpenCV統計應用-直方圖比較

2009-01-02T21:45:00.023+08:00

cvCompareHist(),是比較兩個統計直方圖的分布,總共有四個方法,被定義如下：

#define CV_COMP_CORREL 0
#define CV_COMP_CHISQR 1
#define CV_COMP_INTERSECT 2
#define CV_COMP_BHATTACHARYYA 3

而這些方法分別為相關係數,卡方,交集法以及在做常態分布比對的Bhattacharyya距離,這些方法都是用來做統計直方圖的相似度比較的方法,而且,都是根據統計學的概念,這邊就簡單的拿來用灰階統計直方圖來比較,而這部份的比較方式,是由圖形的色彩結構來著手,下面就簡單的用三種情況來分析它們距離比較的方式

直方圖比較實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
#include <stdlib.h>

int HistogramBins = 256;
float HistogramRange1[2]={0,255};
float *HistogramRange[1]={&HistogramRange1[0]};

int main()
{
    IplImage *Image1=cvLoadImage("RiverBank.jpg",0);
    IplImage *Image2=cvLoadImage("DarkClouds.jpg",0);

    CvHistogram *Histogram1=cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    CvHistogram *Histogram2=cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);

    cvCalcHist(&Image1,Histogram1);
    cvCalcHist(&Image2,Histogram2);

    cvNormalizeHist(Histogram1,1);
    cvNormalizeHist(Histogram2,1);

    printf("CV_COMP_CORREL : %.4f\n",cvCompareHist(Histogram1,Histogram2,CV_COMP_CORREL));
    printf("CV_COMP_CHISQR : %.4f\n",cvCompareHist(Histogram1,Histogram2,CV_COMP_CHISQR));
    printf("CV_COMP_INTERSECT : %.4f\n",cvCompareHist(Histogram1,Histogram2,CV_COMP_INTERSECT));
    printf("CV_COMP_BHATTACHARYYA : %.4f\n",cvCompareHist(Histogram1,Histogram2,CV_COMP_BHATTACHARYYA));

    cvNamedWindow("Image1",1);
    cvNamedWindow("Image2",1);
    cvShowImage("Image1",Image1);
    cvShowImage("Image2",Image2);
    cvWaitKey(0);
}

原始圖片:

執行結果:

(1)RiverBank.jpg & DarkClouds.jpg

Output:
CV_COMP_CORREL : -0.1407
CV_COMP_CHISQR : 0.6690
CV_COMP_INTERSECT : 0.4757
CV_COMP_BHATTACHARYYA : 0.4490

(2)RiverBank.jpg & RiverBank.jpg

Output:
CV_COMP_CORREL : 1
CV_COMP_CHISQR : 0
CV_COMP_INTERSECT : 1
CV_COMP_BHATTACHARYYA : 0

(3)Black.jpg & White.jpg

Output:
CV_COMP_CORREL : 1
CV_COMP_CHISQR : 1
CV_COMP_INTERSECT : 0
CV_COMP_BHATTACHARYYA : 1

這邊的直方圖比較,則是將它們的統計直方圖用cvNormalizeHist()正規化成1,在由正規化的統計分佈來做直方圖的比較,從上面的輸出結果可以推測出,卡方法以及Bhattacharyya是數值越小圖形越相似,而相關係數則是看圖形的分佈程度,因此第三個Black.jpg&White.jpg所顯示相關係數的結果才會是1,這邊的圖形比較用的是統計學的方法,而一般的比較方式還有歐幾里德距離的方式,在這個cvCompareHist()的函式,也可以實作出多通道的多維度直方圖比較,而且支援CV_HIST_ARRAY及CV_HIST_SPARSE這兩種直方圖資料結構的格式,而這些比對方法的相關公式如下

相關係數法

而它的公式推導如下

再來是卡方的方式

這邊的卡方法跟一般的卡方檢定不太一樣,下面是一般適合度檢定(Goodness of fit test)的公式

o為觀察者次數,e為期望值次數

交集的方式就比較簡單了,兩個直方圖取最小的做累加

再來就是常態分配比對的Bhattacharyya距離

OpenCV統計應用-直方圖反向投影

2008-12-29T12:51:00.015+08:00

影像處理的統計直方圖,可以知道一張圖片在該色彩空間的數據分布狀況,而這邊,就要介紹到直方圖反向投影的函式,直方圖反向投影,也就是將數據分布的狀況依照Look-up table的方式對應回去,實際上,這個函式是跟前面介紹到的cvLUT()是一樣的,只不過,差別是差異在cvLUT()的第三個引數改變成CvHistogram資料結構的輸入,直方圖反向投影,cvCalcBackProject()的第一個引數是輸入單通道IplImage資料結構,第二個引數是輸出單通道IplImage反向投影圖形資料結構,第三個引數是選定要被反向投影的CvHistogram直方圖資料結構,而cvCalcBackProject()把前面提到的Look-up table的計算方式包在cvCalcBackProject()函式的底層,因此,它可以整合CvHistogram這個資料結構做更多的應用,下面這個就是修改前面的範例"OpenCV統計應用-CvHistogram直方圖資料結構",來做直方圖反向投影的程式

灰階直方圖反向投影
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int HistogramBins = 50;
int HistogramBinWidth;
float HistogramRange1[2]={0,255};
float *HistogramRange[1]={&HistogramRange1[0]};

CvPoint Point1;
CvPoint Point2;

int main()
{
    IplImage *GrayImage1;
    IplImage *Image1;
    IplImage *Image2;
    IplImage *BackProjectImage;
    CvHistogram *Histogram1;
    IplImage *HistogramImage1;

    Image1=cvLoadImage("Riverbank.jpg",1);
    Image2=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,3);
    GrayImage1=cvLoadImage("Riverbank.jpg",0);
    BackProjectImage=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);
    Histogram1 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    HistogramImage1 = cvCreateImage(cvSize(256,300),8,3);

    cvSetZero(HistogramImage1);
    HistogramImage1->origin=1;
    HistogramBinWidth=256/HistogramBins;

    cvCalcHist(&GrayImage1,Histogram1);
    cvNormalizeHist(Histogram1,5000);

    cvThreshHist(Histogram1,50);
    cvCalcBackProject(&GrayImage1,BackProjectImage,Histogram1);
    cvCopy(Image1,Image2,BackProjectImage);

    for(int i=0;i<HistogramBins;i++)
    {
        printf("%f\n",cvQueryHistValue_1D(Histogram1,i));
        Point1=cvPoint(i*HistogramBinWidth,0);
        Point2=cvPoint((i+1)*HistogramBinWidth,(int)        cvQueryHistValue_1D(Histogram1,i));

        cvRectangle(HistogramImage1,Point1,Point2,CV_RGB(127,127,127));
    }

    cvNamedWindow("Histogram1",1);
    cvNamedWindow("Riverbank",1);
    cvNamedWindow("Back Project RiverBank",1);
    cvShowImage("Riverbank",Image1);
    cvShowImage("Back Project RiverBank",Image2);
    cvShowImage("Histogram1",HistogramImage1);
    cvWaitKey(0);
}

執行結果:

這邊就是拿前面灰階去除較小直方圖區塊的程式碼做修改,然後將前面去除最小區塊的部份,對應到彩色的圖片去了,顯示的結果會是,只要是影像裡面灰階值分佈數量比較少的數據,全部都被對應成黑色的像素,也就是全部都變成0,而cvCalcBackProject(),就跟cvLUT()一樣直接拿直方圖的數據去對應,假設一個直方圖從頭開始的數據為254,129,80,70....那麼只要是圖片內像素值為1的數據就會對應到254,1的圖片裡面像素值的數據就會直接變成254,像素值為2的就會對應到129,2的數據就會直接變成129,以此類推,因此,用cvThreshHist()去除小於50的直方圖區段,讓小於50的全部歸0,再來,對他做一個圖片的反向投影,所對應出來的結果雖然不是0或255,可是它卻可以直接拿來當做是遮罩,也就是說,直接拿來給cvCopy()來做對應,因此,反向投影的結果就出來啦,在遮罩的部份就要參考"資料結構操作與運算-圖形的Mask遮罩實作",而這段程式碼則是用到"OpenCV統計應用-CvHistogram資料結構操作"裡面的部份.

在OpenCV Documentation的部份有提到cvCalcBackProject()可以對HSV色彩空間的Hue做反向投影,到底是怎麼實作出來呢？在OpenCV的Sample Code裡面有一個camshift.c的程式,就是用到這個反向投影的函式,而它的反向投影的簡單範例,原理就如下所示

HSV色彩空間反向投影
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

IplImage *Image1,*Image2;
IplImage *HSVImage;
IplImage *HueImage;
IplImage *BackProjectHueImage,*BackProjectImage;
CvHistogram *Histogram1;
IplImage *HistogramImage1;
CvPoint Point1,Point2;

int HueValue=0;

int HistogramBins = 180;
int HistogramBinWidth;
float HistogramRange1[2]={0,180};
float *HistogramRange[1]={&HistogramRange1[0]};
void onTrackbar(int position);

int main()
{

    Image1 =cvLoadImage("Riverbank.jpg",1);
    HSVImage = cvCreateImage( cvGetSize(Image1),8,3);
    HueImage = cvCreateImage( cvGetSize(Image1),8,1);

    BackProjectHueImage = cvCreateImage( cvGetSize(Image1),8,1);
    BackProjectImage = cvCreateImage( cvGetSize(Image1),8,3);

    Histogram1 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    HistogramImage1 = cvCreateImage(cvSize(180,300),8,3);
    HistogramImage1->origin=1;

    cvCvtColor( Image1, HSVImage, CV_BGR2HSV );
    cvSplit(HSVImage,HueImage,0,0,0);
    cvCalcHist( &HueImage, Histogram1);
    cvNormalizeHist(Histogram1,5000);
    cvZero( HistogramImage1 );
    cvNot(HistogramImage1,HistogramImage1);
    HistogramBinWidth = HistogramImage1->width/HistogramBins;
    for(int i=0;i<HistogramBins;i++)
    {

        Point1=cvPoint(i,0);
        Point2=cvPoint(i,(int)cvQueryHistValue_1D(Histogram1,i));
        printf("%d\n",(int)cvQueryHistValue_1D(Histogram1,i));
        cvLine(HistogramImage1,Point1,Point2,CV_RGB(127,127,127));
    }
    cvNamedWindow("Riverbank",1 );
    cvNamedWindow("Hue Histogram",1);
    cvCreateTrackbar("Hue Thresh","Riverbank",&HueValue,250,onTrackbar);
    cvShowImage("Riverbank",Image1);
    cvShowImage("Hue Histogram",HistogramImage1);
    cvWaitKey(0);

}

void onTrackbar(int position)
{
    IplImage *Image2=cvCreateImage( cvGetSize(Image1),8,3);
    CvHistogram *Histogram2= cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    cvCopyHist(Histogram1,&Histogram2);

    cvThreshHist(Histogram2,position);
    cvCalcBackProject(&HueImage, BackProjectHueImage, Histogram2);
    cvCopy(Image1,Image2,BackProjectHueImage);

    cvZero( HistogramImage1 );
    cvNot(HistogramImage1,HistogramImage1);
    HistogramBinWidth = HistogramImage1->width/HistogramBins;
    for(int i=0;i<HistogramBins;i++)
    {

        Point1=cvPoint(i,0);
        Point2=cvPoint(i,(int)cvQueryHistValue_1D(Histogram2,i));
        printf("%d\n",(int)cvQueryHistValue_1D(Histogram2,i));
        cvLine(HistogramImage1,Point1,Point2,CV_RGB(127,127,127));
    }
    cvShowImage("Hue Histogram",HistogramImage1);
    cvShowImage("Riverbank",Image2);
}

執行結果:

在OpenCV裡面,HSV色彩空間,色調(Hue)值的範圍在0~180,飽和度(Saturation)的範圍在0~255,亮度(Value)的範圍在0~255,而這邊就只取色調(Hue)值在做反向投影,開啟了一個Track bar的功能,並且利用cvCvtColor()將BGR的色彩空間轉換成HSV,並且用cvSplit()通道分割取Hue通道的圖片,計算Hue值的直方圖,在onTrackbar()的部份,則是用Trackbar來調整去除cvThreshHist()的最小區塊的臨界值,去除之後在反向投影到原始的Hue的圖片,在由反向投影的結果當做遮罩,直接跟彩色圖片做對應.

而cvCalcBackProject()不單單只有這樣的功能,它可以對多維度空間的色彩直方圖做對應,cvCalcBackProject()提供了一個當CvHistogram資料結構維度為3的時候的一個反向投影,下面的這個例子就以HSV的色彩空間為例,建構一個三維的CvHistogram資料結構

HSV三維直方圖反向投影
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

IplImage *Image1,*Image2;
IplImage *HSV;
IplImage *HueImage,*SaturationImage,*ValueImage;
IplImage *ImageArray[3];
IplImage *BackProjectImage;
CvHistogram *Histogram1;
IplImage *HistogramImage1;
CvPoint Point1,Point2;

int HueValue=0;

int HistogramBins[3] ={180,256,256};
int HistogramBinWidth;
float HistogramRange1[6]={0,180,0,255,0,255};
float *HistogramRange[3]={&HistogramRange1[0],&HistogramRange1[2],&HistogramRange1[4]};
void onTrackbar(int position);

int main()
{

    Image1 =cvLoadImage("Riverbank.jpg",1);
    HSV = cvCreateImage( cvGetSize(Image1),8,3 );
    HueImage = cvCreateImage( cvGetSize(Image1),8,1 );
    SaturationImage = cvCreateImage( cvGetSize(Image1),8,1 );
    ValueImage = cvCreateImage( cvGetSize(Image1),8,1);
    ImageArray[0]=HueImage;
    ImageArray[1]=SaturationImage;
    ImageArray[2]=ValueImage;

    BackProjectImage = cvCreateImage( cvGetSize(Image1),8,3 );

    Histogram1 = cvCreateHist(3,HistogramBins,CV_HIST_SPARSE,HistogramRange);

    cvCvtColor( Image1, HSV, CV_BGR2HSV );
    cvSplit(HSV,HueImage,SaturationImage,ValueImage,0);

    cvCalcHist( ImageArray, Histogram1);

    cvNamedWindow("Riverbank",1 );
    cvCreateTrackbar("Hue Thresh","Riverbank",&HueValue,200,onTrackbar);
    cvShowImage("Riverbank",Image1);
    cvWaitKey(0);

}

void onTrackbar(int position)
{
    CvHistogram *Histogram2= cvCreateHist(3,HistogramBins,CV_HIST_SPARSE,HistogramRange);
    IplImage *Image2=cvCreateImage( cvGetSize(Image1),8,3 );
    IplImage *BackProjectImage = cvCreateImage( cvGetSize(Image1),8,1 );

    cvCopyHist(Histogram1,&Histogram2);

    cvThreshHist(Histogram2,position);
    cvCalcBackProject( ImageArray, BackProjectImage, Histogram2);
    cvCopy(Image1,Image2,BackProjectImage);

    cvShowImage("Riverbank",Image2);
}

執行結果:

在三維空間的作法上面,就要參考到前面"OpenCV統計應用-CvHistogram直方圖資料結構"關於三維空間製作的部份,除了用cvCvtColor()將色彩空間轉換,用cvSplit()將通道做分割,還要做個圖形陣列(ImageArray)來讓cvCalcHist()這個函式做運算,計算出來的結果為一個CvHistogram的三維空間稀疏矩陣直方圖,而在onTrackbar()的部份,cvCalcBackProject()直方圖反向投影亦是同樣要用ImageArray做輸入,而輸出則是一個單通道的圖形,在稀疏矩陣裡面,由於維度為三維,所以他所形成的統計直方圖數值都是極小,所以門檻值只要一點點就快要全部都分佈了,而這個三維空間的反向投影可以如此建構,是基於Look-up table的功能來實現,只不過他的缺點是,每一個維度的Look-up table範圍是0~255,因此如果是像Hue值的範圍0~180,它的數值就會被模糊化,也就是數據會被些許位移,但是仍不會影響它出來結果的精確度

在這個直方圖反向投影的部份,也可以結合連通成分來做去除某一門檻值的連通分量

cvCalcBackProject()
將統計直方圖的分布數據根據Look-up table對應回去,也就是說,當今天CvHistogram資料結構內的數據分佈是243,110,0,60...則使用cvCalcBackProject()函式單通道的圖片像素值會是,當遇到像素值為1的時候變243,像素值為2的時候變110,依此類推,cvCalcBackProject()直方圖反向投影可以根據多維度設計,而cvCalcBackProject()第一個引數為輸入單通道IplImage或CvMat資料結構,第二個引數為輸入單通道反向投影IplImage或CvMat資料結構,第三個引數為輸入CvHistogram資料結構
cvCalcBackProject(輸入單通道IplImage或CvMat資料結構,輸入單通道反向投影IplImage或CvMat資料結構,輸入CvHistogram資料結構)

OpenCV統計應用-影像增強,亮度/對比實作

2008-11-12T13:08:00.020+08:00

在一般顯示螢幕以及圖形處理的應用軟體上,都會有一個亮度/對比的色彩(Brightness/Contrast)調整,它是屬於影像增強的部份,在OpenCV裡面的Sample Code裡面就有這樣的灰階程式的實作,在這邊就修改了OpenCV的Sample Code,來做色彩增強的亮度/對比的程式,而在一般的亮度/對比來講亮度(Brightness)的範圍為0~200而對比(Contrast)亦是0~200,它們由一條線性函數的公式所定義,對比所代表的是斜率,亮度則是偏移量,這條線性公式代表的是Look-up table的對應,它的數學式定義如下

原始的亮度對比數值範圍為-100~100之間,C代表對比,B代表亮度

對於對比率(Contrast ratio)來講,delta範圍應該落在0~255,這邊將對比率的公式做重新的調整

對比率代表著斜率的α值,而亮度則是決定線性公式位移的情況,也就是β值,而Y=αX+β這個線性公式它所表達的情況如下

α值的範圍落在0~255之間,而它的情況如下

再來下面是用虛擬碼的方式表達亮度/對比的演算法

下面就是亮度/對比的程式了

亮度/對比實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int BrightnessPosition = 100;
int ContrastPosition = 100;

int HistogramBins = 64;
int HistogramBinWidth;
float HistogramRange1[2]={0,256};
float *HistogramRange[1]={&HistogramRange1[0]};

IplImage *Image1,*Image2;
CvHistogram *Histogram1;
IplImage *HistogramImage;

uchar LookupTableData[256];
CvMat *LookupTableMatrix;
IplImage *LookupTableImage;
CvPoint Point1,Point2;

void OnTrackbar(int Position)
{
    int Brightness=BrightnessPosition-100;
    int Contrast=ContrastPosition -100;
    double Delta;
    double a,b;
    int y;

    //Brightness/Contrast Formula
    if(Contrast>0)
    {
        Delta=127*Contrast/100;
        a=255/(255-Delta*2);
        b=a*(Brightness-Delta);

        for(int x=0;x<256;x++)
        {
            y=(int)(a*x+b);

            if(y<0)
                y=0;
            if(y>255)
                y=255;

            LookupTableData[x]=(uchar)y;
        }
    }
    else
    {
        Delta=-128*Contrast/100;
        a=(256-Delta*2)/255;
        b=a*Brightness+Delta;

        for(int x=0;x<256;x++)
        {
            y=(int)(a*x+b);

            if(y<0)
                y=0;
            if(y>255)
                y=255;

            LookupTableData[x]=(uchar)y;
        }
    }
    //End

    //Look up table sketch
    cvSetZero(LookupTableImage);
    cvNot(LookupTableImage,LookupTableImage);
    Point2=cvPoint(0,LookupTableData[0]);
    for(int i=0;i<256;i++)
    {
        Point1=cvPoint(i,LookupTableData[i]);
        cvLine(LookupTableImage,Point1,Point2,CV_RGB(0,0,0),3);
        Point2=Point1;
    }
    cvLUT(Image1,Image2,LookupTableMatrix);
    //End

    //Gray Level Histogram
    cvCalcHist(&Image2,Histogram1);
    cvNormalizeHist(Histogram1,3000);

    cvSetZero(HistogramImage);
    cvNot(HistogramImage,HistogramImage);
    HistogramBinWidth = HistogramImage->width/HistogramBins;
    for(int i=0;i<HistogramBins;i++)
    {
        Point1=cvPoint(i*HistogramBinWidth,0);
        Point2=cvPoint((i+1)*HistogramBinWidth,(int)cvQueryHistValue_1D(Histogram1,i));
        cvRectangle(HistogramImage,Point1,Point2,CV_RGB(0,0,0),CV_FILLED);
    }
    //End

    cvShowImage("Gray Level Histogram",HistogramImage);
    cvShowImage("Brightness/Contrast",Image2);
    cvShowImage("Image Enhance",LookupTableImage);
    cvZero(Image2);
}

int main()
{
    Image1=cvLoadImage("DarkClouds.jpg",0);
    Image2=cvCloneImage(Image1);

    Histogram1=cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    HistogramImage = cvCreateImage(cvSize(320,200),8,1);

    LookupTableImage=cvCreateImage(cvSize(256,256),8,3);
    LookupTableMatrix=cvCreateMatHeader(1,256,CV_8UC1);
    cvSetData(LookupTableMatrix,LookupTableData,0);

    LookupTableImage->origin=1;
    HistogramImage->origin=1;

    cvNamedWindow("Brightness/Contrast",1);
    cvNamedWindow("Gray Level Histogram",1);
    cvNamedWindow("Image Enhance",1);

    cvCreateTrackbar("brightness","Brightness/Contrast",&BrightnessPosition,200,OnTrackbar);
    cvCreateTrackbar("contrast","Brightness/Contrast",&ContrastPosition,200,OnTrackbar);

    OnTrackbar(0);

    cvWaitKey(0);
}

執行結果:

這隻程式同樣也是用到CvHistogram資料結構,使用到兩個拉軸(Trackbar),以及Look-up table的應用,在//Brightness/Contrast Formula的註解內所包的就是亮度/對比演算法虛擬碼的實作,再來就是把它的線性系統化出來,也就是Y=αX+β的函數方程式,這個方程式,當然同等於Look up table,而之後,在把他們灰階直方圖的分布畫出來,在main()裡面,當然是先讀取目標圖片轉成灰階,初始化繪製直方圖與線性系統圖片的空間,創立三個視窗介面,設立兩個拉軸,並且將拉軸的事件函式設定成同一個的副程式的名稱.而對於影像增強(Image Enhance)這個視窗介面,它所代表的含意如下

X軸代表為是原始灰階的輸入值,而Y軸代表的是灰階值所對應的結果,而X軸跟Y軸的範圍都是0~255,而這條直線公式也會受到斜率(α)以及平移(β)的結果改面灰階值輸入以及輸出的對應,它是將一張原始灰階圖片的每一個像素值做線性函式的對應,使得每個灰階值對應出來的結果產生了變化,由下面可以知道它(LUT)對應的關係

(a)亮度條為0因此小於100的灰階值都為0而灰階值方圖也像左偏移

(b)亮度條為100,因此大於156的灰階值都為255,而灰階值方圖也都向右偏移

(c)對比為0,這個時候斜率α為0,因此輸入的0~255的灰階值輸出都固定為128,因此整張圖片都是灰階值128的影像,而灰階直方圖則是所有數據都集中在128

(d)對比為100,這個時候斜率為255,而這樣的圖片又可以叫做二值化圖片,因為輸出結果非黑即白,而移動亮度則是在平移二值化的門檻值,由灰階值方圖可以得知,所有數據都被分開到0跟255兩邊

上面所表達的其實就是一種Look-up table的表達方式,藉由一個輸入灰階值的矩陣,對應岀另一個不同的灰階值數據,因此改變了原始灰階值的數據,而整張圖片也因此產生了變化

OpenCV統計應用-直方圖等化

2008-11-07T17:44:00.029+08:00

直方圖等化(Histogram Equalization)為一種使用統計方法的影像處理程式設計,它的功能為將統計直方圖的色彩分布平均的打散在直方圖裡,也就是說,讓一張圖的直方圖分布均勻化,同樣的也是使用到LUT(Look-up Table)的方法

而在設計直方圖等化不可或缺的就是需要先知道統計學的機率密度函數(Probability Density Function,PDF)以及累積分配函數(Cumlative Distribution Function,CDF)相關的基本概念,在這邊會簡單的說明機率密度函數及累積分配函數的概念還有直方圖等化的程式設計實作,而OpenCV也提供了直接使用直方圖等化的函式,cvEqualizeHist(),下面就簡單的製作直方圖等化的程式設計

RGB結構直方圖等化
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int main()
{
    IplImage *Image1;
    IplImage *Image2;
    IplImage *RedImage;
    IplImage *GreenImage;
    IplImage *BlueImage;

    Image1=cvLoadImage("DarkClouds.jpg",1);
    Image2=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,3);
    RedImage=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);
    GreenImage=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);
    BlueImage=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);

    cvSplit(Image1,BlueImage,GreenImage,RedImage,0);

    cvEqualizeHist(BlueImage,BlueImage);
    cvEqualizeHist(GreenImage,GreenImage);
    cvEqualizeHist(RedImage,RedImage);

    cvMerge(BlueImage,GreenImage,RedImage,0,Image2);

    cvNamedWindow("DarkClouds",1);
    cvShowImage("DarkClouds",Image1);
    cvNamedWindow("Equalize DarkClouds",1);
    cvShowImage("Equalize DarkClouds",Image2);

    cvWaitKey(0);

}

原始圖片:

執行結果:

當然,使用了函式之後就完全不用知道直方圖等化在跑什麼啦,cvEqualizeHist()只支援單通道uchar型別的圖片資料結構,因此,就需要用cvSplit()分割,而這裡所做的是將RGB三原色直接做直方圖等化,這個直方圖等化的結果是彩色的,將個別的RGB直分割後再用cvMerge()合併,通道分割與合併的部份就要參考"資料結構操作與運算-通道的分割,合併與混合"這個地方了.

再來就是比較複雜的直方圖等化演算法,要知道直方圖等化如何計算就要知道以下的步驟

1.資料結構初始化
2.機率密度函數(PDF)的計算
3.累積分配函數(CDF)的計算
4.LUT對應

在機率密度函數的部分,在統計直方圖來講,它的X軸範圍會落在0~255,而它的Y軸範圍代表著機率,由0~255對應的數據代表著它發生的機率,也就是說,這代表著0到255數值出現的機率分佈,0~255代表著所有的機率的範圍,那從0加到255的總和會為1!它所代表的公式如下

PDF機率分佈

而累積分配函數就更簡單了,它只是將機率密度函數做累積的計算,也就是說,當數值為255的時候它的機率會為1,而且曲線會是由下往上的成長,表達的方式就是將機率密度函數做累加,它的公式如下

CDF機率分佈

再來是LUT的部份了,它只是開一個0~255的陣列做對應,將累積分配函數的數值乘以255,將它存放在0~255的陣列裡,這樣會出現許多重複的部份,因此,將Look-up table對應原來的圖形的灰階值,就可以將原本的統計直方圖分布打散掉了,而且直方圖等化的結果也會發現很多地方數值是空心的,也是因為如此,直方圖的分布就會向外推擠了

累積分配直方圖乘255存在LUT陣列的情況

再來,就是直方圖等化的程式了

灰階直方圖等化
#include <stdio.h>
#include <cv.h>
#include <highgui.h>

int HistogramBins = 256;
float HistogramRange1[2]={0,255};
float *HistogramRange[1]={&HistogramRange1[0]};

float CumulativeNumber;
float CDFArray[256];
uchar LookupTableData[256];

void DrawHistogramImage(CvHistogram *Histogram,IplImage *HistogramImage,int HistogramBins);

int main()
{
    IplImage *Image1;
    IplImage *Image2;
    CvHistogram *Histogram1;
    CvHistogram *Histogram2;
    CvMat *LookupTableMatrix;
    IplImage *GrayHistogramImage;
    IplImage *EqualizeHistogramImage;

    Image1=cvLoadImage("DarkClouds.jpg",0);
    Image2=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);
    Histogram1 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    Histogram2 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);

    LookupTableMatrix = cvCreateMat(1,256,CV_8UC1);
    GrayHistogramImage=cvCreateImage(cvSize(256,250),IPL_DEPTH_8U,3);
    EqualizeHistogramImage=cvCreateImage(cvSize(256,250),IPL_DEPTH_8U,3);

    GrayHistogramImage->origin=1;
    EqualizeHistogramImage->origin=1;

    cvCalcHist(&Image1,Histogram1);
    DrawHistogramImage(Histogram1,GrayHistogramImage,HistogramBins);

    //Probability Density Function (PDF)
    cvNormalizeHist(Histogram1,1);
    //End

    //Cumulative Distribution Function (CDF)
    CumulativeNumber=0;
    for(int i=0;i<HistogramBins;i++)
    {
        CumulativeNumber=CumulativeNumber+cvQueryHistValue_1D(Histogram1,i);
        CDFArray[i]=CumulativeNumber;
    }
    //End

    //Make Look-up Table
    printf("Look-up Table Number:\n");
    for(int i=0;i<HistogramBins;i++)
    {
        LookupTableData[i]=(uchar)(255*CDFArray[i]);
        printf("%f\n",(255*CDFArray[i]));
    }
    //End

    cvSetData(LookupTableMatrix,LookupTableData,CV_AUTOSTEP);
    cvLUT(Image1,Image2,LookupTableMatrix);

    cvCalcHist(&Image2,Histogram2);
    DrawHistogramImage(Histogram2,EqualizeHistogramImage,HistogramBins);

    cvNamedWindow("DarkClouds",1);
    cvNamedWindow("Equalize DarkClouds",1);
    cvNamedWindow("Gray Histogram",1);
    cvNamedWindow("Equalize Histogram",1);
    cvShowImage("DarkClouds",Image1);
    cvShowImage("Equalize DarkClouds",Image2);
    cvShowImage("Gray Histogram",GrayHistogramImage);
    cvShowImage("Equalize Histogram",EqualizeHistogramImage);

    cvWaitKey(0);

}

void DrawHistogramImage(CvHistogram *Histogram,IplImage *HistogramImage,int HistogramBins)
{
    CvPoint Point1;
    CvPoint Point2;

    for(int i=0;i<HistogramBins;i++)
    {
        Point1=cvPoint(i,(int)(cvQueryHistValue_1D(Histogram,i)/20));
        Point2=cvPoint(i,0);
        cvLine(HistogramImage,Point1,Point2,CV_RGB(127,127,127));

    }
}

執行結果:

前面資料結構初始化的部份,分別製造了兩個IplImage圖形資料結構,兩個CvHistogram直方圖資料結構,一個Look-up table對應矩陣結構,兩個直方圖圖形輸出的資料結構,並且將它做基本的設定,在用cvCalcHist()放入圖形的資料,繪製出灰階直方圖的圖形,再來,利用cvNormalizeHist()計算機率密度函數(PDF),機率密度函數有一條規則,所有機率的總和為1,因此,利用cvNormalizeHist()很快的就可以把機率密度函數算出,而在做累積分配直方圖,則是用for迴圈慢慢的累加,製造Look-up Table的時候,則是讓它乘以255並且給它用uchar型別轉換,將轉換後的結果存到CvMatrix資料結構裡,在用cvLUT()函式做對應,因此均化的灰階圖片就這樣被製造出來啦.

cvEqualizeHist()
將單通道8bits uchar型別的圖片做直方圖等化的演算法,輸入為單通道uchar型別的IplImage資料結構,輸出為直方圖等化後單通道uchar型別IplImage資料結構
cvEqualizeHist(輸入單通道uchar型別IplImage資料結構,輸出單通道uchar型別資料結構)

OpenCV統計應用-CvHistogram資料結構操作

2008-11-01T16:28:00.020+08:00

CvHistogram資料結構有提供一些基本的功能,有許多可以不必去直接處理CvHistogram資料結構內部資料的函數,接著就修改前面灰階直方圖的部分做簡單操作的實作

不同範圍灰階直方圖
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int HistogramBins = 50;
int HistogramBinWidth;
float HistogramRange1[2]={0,255};
float *HistogramRange[1]={&HistogramRange1[0]};
float MinValue,MaxValue;
int MinIndex,MaxIndex;
int main()
{
    IplImage *Image1;
    CvHistogram *Histogram1;
    IplImage *HistogramImage1;
    IplImage *HistogramImage2;

    CvPoint Point1;
    CvPoint Point2;

    Image1=cvLoadImage("Riverbank.jpg",0);

    Histogram1 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    HistogramImage1 = cvCreateImage(cvSize(256,300),8,3);
    HistogramImage2 = cvCreateImage(cvSize(256,300),8,3);

    cvSetZero(HistogramImage1);
    HistogramImage1->origin=1;
    HistogramImage2->origin=1;
    HistogramBinWidth=256/HistogramBins;

    cvCalcHist(&Image1,Histogram1);
    cvNormalizeHist(Histogram1,5000);

    cvGetMinMaxHistValue(Histogram1,&MinValue,&MaxValue,&MinIndex,&MaxIndex);
    printf("Histogram1\n");
    printf("The Minimum value : %f\tIndex :%d\n",MinValue,MinIndex);
    printf("The Maximum value : %f\tIndex :%d\n",MaxValue,MaxIndex);

    for(int i=0;i<HistogramBins;i++)
    {

        Point1=cvPoint(i*HistogramBinWidth,0);
        Point2=cvPoint((i+1)*HistogramBinWidth,(int)cvQueryHistValue_1D(Histogram1,i));

        cvRectangle(HistogramImage1,Point1,Point2,CV_RGB(127,127,127));
    }

    cvClearHist(Histogram1);

    HistogramRange1[0]=30;
    HistogramRange1[1]=255;
    HistogramRange[0]=&HistogramRange1[0];

    cvSetHistBinRanges(Histogram1,HistogramRange);
    cvCalcHist(&Image1,Histogram1);

    cvNormalizeHist(Histogram1,5000);

    cvGetMinMaxHistValue(Histogram1,&MinValue,&MaxValue,&MinIndex,&MaxIndex);
    printf("\nHistogram2\n");
    printf("The Minimum value : %f\tIndex :%d\n",MinValue,MinIndex);
    printf("The Maximum value : %f\tIndex :%d\n",MaxValue,MaxIndex);

    for(int i=0;i<HistogramBins;i++)
    {
        Point1=cvPoint(i*HistogramBinWidth,0);
        Point2=cvPoint((i+1)*HistogramBinWidth,(int)cvQueryHistValue_1D(Histogram1,i));

        cvRectangle(HistogramImage2,Point1,Point2,CV_RGB(127,127,127));
    }

    cvReleaseHist(&Histogram1);
    cvNamedWindow("Histogram1",1);
    cvNamedWindow("Histogram2",1);
    cvNamedWindow("Riverbank",1);
    cvShowImage("Riverbank",Image1);
    cvShowImage("Histogram1",HistogramImage1);
    cvShowImage("Histogram2",HistogramImage2);
    cvWaitKey(0);
}

執行結果:

這個執行出來得結果可以看得出來,直方圖的範圍改變之後,長方形的數據有向右平移的樣子,而在之前的前面的程式碼cvQueryHistValue_1D()所直接輸出的數據實在是過大了,在繪製直方圖的時候一定會超出座標軸的範圍,所以就在前面的程式碼以全部除以50來代替,這邊,改用另外一個解決的方法,cvNormalizeHist(),這個函數是將直方圖所有區塊加起來會等於5000,也就是直方圖數據的總和會等於5000,因此,用這種量化的方式可以避免直方圖因為數據過大而圖形繪不出來的問題,近而縮小直方圖所有的數據,當然,這樣的結果直方圖數據會變成相對的,而cvNormalizeHist()內部則是使用到了cvScale()這個函式來計算,cvGetMinMaxHistValue(),則是找出直方圖的最大值跟最小值的方法,跟前面的最大最小值不太一樣的地方,這邊找出了直方圖的最大值將會代表的是眾數,在統計上是很具意義的東西,cvGetMinMaxHistValue()同樣內部使用的是cvMinMaxLoc()的函式.再來的部分就是用cvClearHist()清除直方圖資料,從新設定直方圖的範圍為30~255,使用cvSetHistBinRanges()的函數,最後再用cvReleaseHist()釋放CvHistogram的記憶體

接下來的這個的結果也是有趣,它可以將直方圖較小的區塊去除掉,這個函式在可以用來製作去除最小面積的圖形區塊,當然在這前面就必須要建立一些相關的演算法來建構直方圖,這邊就簡單介紹cvThreshHist()的操作

去除較小直方圖區塊
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int HistogramBins = 50;
int HistogramBinWidth;
float HistogramRange1[2]={0,255};
float *HistogramRange[1]={&HistogramRange1[0]};
float MinValue,MaxValue;
int MinIndex,MaxIndex;
int main()
{
    IplImage *Image1;
    CvHistogram *Histogram1;
    CvHistogram *Histogram2;
    IplImage *HistogramImage1;
    IplImage *HistogramImage2;

    CvPoint Point1;
    CvPoint Point2;

    Image1=cvLoadImage("Riverbank.jpg",0);

    Histogram1 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    HistogramImage1 = cvCreateImage(cvSize(256,300),8,3);
    HistogramImage2 = cvCreateImage(cvSize(256,300),8,3);

    cvSetZero(HistogramImage1);
    HistogramImage1->origin=1;
    HistogramImage2->origin=1;
    HistogramBinWidth=256/HistogramBins;

    cvCalcHist(&Image1,Histogram1);
    cvNormalizeHist(Histogram1,5000);

    for(int i=0;i<HistogramBins;i++)
    {

        Point1=cvPoint(i*HistogramBinWidth,0);
        Point2=cvPoint((i+1)*HistogramBinWidth,(int)cvQueryHistValue_1D(Histogram1,i));

        cvRectangle(HistogramImage1,Point1,Point2,CV_RGB(127,127,127));
    }

    Histogram2 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    cvCopyHist(Histogram1,&Histogram2);
    cvThreshHist(Histogram2,50);

    for(int i=0;i<HistogramBins;i++)
    {
        Point1=cvPoint(i*HistogramBinWidth,0);
        Point2=cvPoint((i+1)*HistogramBinWidth,(int)cvQueryHistValue_1D(Histogram2,i));

        cvRectangle(HistogramImage2,Point1,Point2,CV_RGB(127,127,127));
    }

    cvReleaseHist(&Histogram1);
    cvReleaseHist(&Histogram2);
    cvNamedWindow("Histogram1",1);
    cvNamedWindow("Histogram2",1);
    cvNamedWindow("Riverbank",1);
    cvShowImage("Riverbank",Image1);
    cvShowImage("Histogram1",HistogramImage1);
    cvShowImage("Histogram2",HistogramImage2);
    cvWaitKey(0);
}

執行結果:

在這邊同樣建立了一個50個區塊,範圍為0~255的直方圖,這裡開啟了兩個CvHistogram的資料結構,並且初始化用cvCopyHist()複製直方圖裡面的資料,而第二個直方圖資料結構Histogram2則是給定一個區塊門檻值,也就是當小於50的區塊數據就會被剃掉,因此圖形顯示的結果,小於50數據的直方圖都被刪除了.

cvNormalizeHist()
將所有直方圖的數據標準化,也就是將所有數據總和為指定的一個數,第一個引數為輸入CvHistogram資料結構,第二個引數為輸入double型別指定總和數據,cvNormalizeHist()內部使用cvScale()函式
cvNormalizeHist(輸入CvHistogram資料結構,輸入double型別總和數據)

cvGetMinMaxHistValue()
取得直方圖的最大最小值的數據,以及它得位置(Index),對於灰階統計直方圖,可以知道是哪個色彩數據擁有最大值,第一個引數為輸入CvHistogram資料結構,第二個引數為輸出float型別最小值數據,第三個引數為輸出float型別最大值數據,第四個引數為輸出int型別最小值的陣列索引(Index),第五個引數為輸出int型別最大值的陣列索引(Index)
cvGetMinMaxHistValue(輸入CvHistogram資料結構,輸出float最小值數據,輸出float最大值數據,輸出int最小值索引,輸出int最大值索引)

cvClearHist()
清除目標CvHistogram直方圖資料結構內的資料.
cvClearHist(目標CvHistogram資料結構)

cvSetHistBinRanges()
從新設立直方圖的數據範圍,也就是要挑選的直方圖最小值下限以及最大值上限,第一個引數為輸入CvHistogram資料結構,第二個引數為輸入float型別的範圍二維陣列
cvSetHistBinRanges(輸入CvHistogram資料結構,輸入float型別範圍二維陣列)

cvReleaseHist()
釋放CvHistogram直方圖資料結構記憶體位址
cvReleaseHist(目標CvHistogram直方圖資料結構)

cvCopyHist()
複製目標CvHistogram資料結構,包括設定直以及直方圖數據資料,第一個引數為輸入目標CvHistogram資料結構直方圖,第二個引數為輸出CvHistogram資料結構直方圖
cvCopyHist(輸入CvHistogram資料結構,輸出要被複製的CvHistogram資料結構)

cvThreshHist()
去除小於目標數據的直方圖區塊,第一個引數為輸入CvHistogram直方圖資料結構,第二個引數為輸入double型別要被刪除的直方圖區塊數據大小
cvThreshHist(輸入CvHistogram資料結構,輸入小於門檻值double型別的數據)

OpenCV 1.1

2008-10-23T12:28:00.017+08:00

新版的OpenCV 1.1 ,在2008年10月15號發佈,基本上只是優化前一版的內容,刪除了cvcam函式庫,並且提供了一本歐萊里的工具書
http://sourceforge.net/project/showfiles.php?group_id=22870

Notes:

This is OpenCV 1.1pre1, the first intermediate update after 1.0.It features better stereo support, SURF features, better performance of some functions, DirectShow support in highgui etc. It also includes many bug-fixes.

Changes:

>>> New functionality/features: <<<

    - General:

        * Octave bindings have been added. See interfaces/swig/octave
         (for now, Linux only)

        * [Windows] OpenCV is now built with VS2005 with SSE2 and OpenMP support included
         (if you want to rebuild OpenCV using Express or Standard Edition of VS, use
         _make\opencv.vs2005.no_openmp.sln).

        * [Windows] Python bindings have been updated to use Python 2.6

        * [Windows] cvcam has been removed (as videoInput is now supported by highgui)

    - CXCORE, CV, CVAUX:

        * Speeded-up Robust Features (SURF), contributed by Liu Liu.
         see samples/c/find_obj.cpp and the documentation opencvref_cv.htm

        * Many improvements in camera calibration:
            - Added stereo camera calibration: cvStereoCalibrate, cvStereoRectify etc.
            - Single camera calibration now uses Levenberg-Marquardt method and supports extra flags to
             switch on/off optimization of individual camera parameters
            - The optional 3rd radial distortion parameter (k3*r^6) is now supported in every calibration-related
             function

        * 2 stereo correspondence algorithms:
            - very fast block matching method by Kurt Konolige
             (processes the Tsukuba stereo pair in <10ms on Core2Duo laptop)
            - slow but more accurate graph-cut based algorithm by Kolmogorov and Zabin

        * Better homography estimation algorithms (RANSAC and LMEDs)

        * new C++ template image classes contributed by Daniel Filip (Google inc.). see
            opencv/cxcore/include/cvwimage.h

        * Fast approximate nearest neighbor search (by Xavier Delacour)

        * Codebook method for background/foreground segmentation (by Gary Bradski)

        * Sort function (contributed by Shiqi Yu)

        * [OpenCV+IPP] Face Detection (cvHaarDetectObjects) now runs much faster (up to 2x faster) when
            using IPP 5.3 or higher.

        * Much faster (~4x faster) fixed-point variant of cvRemap has been added

    - MLL:

        * Python bindings for MLL have been added. There are no samples yet.

    - HighGUI:

        * [Windows, 32bit] Added support for videoInput library. Hence, cvcam is [almost] not needed anymore

        * [Windows, 32bit] FFMPEG can now be used for video decoding/encoding via ffopencv*.dll

        * [Linux] Added unicap support

        * Improved internal video capturing and video encoding APIs

    - Documentation:

        * OpenCV book has been published (sold separately :) see docs/index.htm)

        - New samples (opencv/samples):

        * Many Octave samples

        * find_obj.cpp (SURF), bgfg_codebook.cpp (Codebook BG/FG segmentation), stereo_calib.cpp
            (Stereo calibration and stereo correspondence)

>>> Bug fixes: <<<

    Many thanks to everybody who submitted bug reports and/or provided the patches!

    * added dma_unlisten to icvCloseCAM_DC1394 (thanks to Victor Benso)

    * LMEDs algorithm for cvFindFundamentalMat estimation has been fixed

    * Broken debug build of highgui in VS2005/2008 (SF #2056185, SF #1709435)

    * cvFindChessboardCorners memory leak and incorrect size checks
    (SF #1972896, SF #1910323)

    * disabling GTK causes v4l runtime error (SF #2088102)

    * cvsetmousecallback bug (SF #2053529)

    * libhighgui needed deprecated "img_convert" replacement (SF #2011753)

    * Segfault in createsamples caused by uninitialized variable (SF #1977989)

    * Data Alignment Issue in bgfg_gaussmix (SF #1961755)

    * libpng need to be updated (SF #1952793)

    * cvCreateVideoWriter_Win32 - identifier not found (SF #1944254)

    * Bug in cvmorph.cpp (SF #1908844)

    * dilate (cvDilate) works bogus with default kernel (SF #1887130)

    * CvEM non-default constructor declared but not defined (SF #1830346)

    * cvFloodFill (in ver 1.0) Hangs (SF #1790565)

    * double delete in CvImage (SF #1733437)

    * cvFilter2D anchor default value is not working properly (SF #1713779)

    * cvCvtColor - Bug? in converting HSV2RGB (SF #1676344)

    * Invalid selection of the MKL-dll version in cvUseOptimized() (SF #1633017)

>>> Known issues:

    * Borland compiler is not supported (but might work) in this update.

    * 64-bit Windows is not supported (but might work) in this update.

    * SF bug tracker still contains over 50 records of open bugs. Many of them will be addressed in the next
     update.

OpenCV統計應用-CvHistogram直方圖資料結構

2008-10-19T10:06:00.020+08:00

CvHistogram可以提供直方圖的計算,並且可以支援多個維度的直方圖設計,但是在繪製直方圖圖形的時候就需要自己用繪圖函式來繪製,CvHistogram算是一個比較複雜的資料結構,由於它在維度小於二的可以用密集的資料結構,CvMatND來設計,可以當它的直方圖維度大於二的時候,就必須要用到稀疏矩陣,CvSparseMat的方式,下面就是用CvHistogram結構計算一維直方圖的例子

CvHistogram一維設計
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int HistogramBins = 256;
float HistogramRange1[2]={0,255};
float *HistogramRange[1]={&HistogramRange1[0]};

int main()
{
    IplImage *Image1;
    CvHistogram *Histogram1;
    IplImage *HistogramImage1;

    Image1=cvLoadImage("Riverbank.jpg",0);

    Histogram1 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    HistogramImage1 = cvCreateImage(cvSize(256,300),8,3);
    HistogramImage1->origin=1;

    cvCalcHist(&Image1,Histogram1);

    printf("type is : %d\n",Histogram1->type);
    printf("Low Bound is : %.f\n",Histogram1->thresh[0][0]);
    printf("Up Bound is : %.f\n",Histogram1->thresh[0][1]);
    printf("The Bins are : %d\n",((CvMatND *) Histogram1->bins)->dim[0].size);

    printf("\nGray Level Values:\n");
    for(int i=0;i<HistogramBins;i++)
    {
        printf("%.f \n",((CvMatND *) Histogram1->bins)->data.fl[i]);
        cvLine(HistogramImage1,cvPoint(i,0),cvPoint(i,(int)(cvQueryHistValue_1D(Histogram1,i)/10)),CV_RGB(127,127,127));
    }

    cvNamedWindow("Histogram",1);
    cvNamedWindow("Riverbank",1);
    cvShowImage("Riverbank",Image1);
    cvShowImage("Histogram",HistogramImage1);
    cvWaitKey(0);
}

原始圖片:

執行結果:

跟前面直接用矩陣累加設計出的統計直方圖一樣,到最後還是要自己用繪圖函式自己畫,CvHistogram這個資料結構的特色,它可以設定要用幾個區塊來呈現,這邊設立為256個區塊,而實際上,它可以用比較模糊的方式設定直方圖的區域,也就是可以小於256個的區塊數,是利用切割成n等分的方式,然後它可以自行設定上界(Up Bound)以及下界(Low Bound),在灰階圖裡面數據的範圍為0~255,而CvHistogram資料結構可以設立為上界30,下界200這樣的方式縮小範圍,而0~29,201~266這範圍內的數據將不會被計算,下面是CvHistogram資料結構的內容,分別為CV_HIST_ARRAY及CV_HIST_SPARSE兩種

1.

2.

CvHistogram的type固定都為CV_HIST_MAGIC_VAL這的參數,與上面程式cvCreateHist()所設定的CV_HIST_ARRAY無關,這也許是OpenCV裡面對於CvHistogram這個結構設計不良的地方,而對於CvHistogram資料結構的參數,定義如下

#define CV_HIST_ARRAY 0
#define CV_HIST_SPARSE 1
#define CV_HIST_TREE CV_HIST_SPARSE

因此,它的參數只有CV_HIST_ARRAY的多維矩陣CvMatND以及CV_HIST_SPARSE的稀疏矩陣CvSparseMat資料結構,而它的thresh是放上下界的資料,thresh2是放動態的維度上下界資料,而一般的直方圖數據資料都是放在bins裡面.而CvHistogram內的CvMatND結構則是在做快速初始化,bins以及mat都是使用同一個記憶體空間,而將直方圖資料的提取就要用到cvQueryHistValue_1D()這個函式了.

接著是將直方圖的空間分區塊的方式實作,將它分為50塊,並且上界為30,下界為200

CvHistogram結構區塊與上下界
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int HistogramBins = 50;
int HistogramBinWidth;
float HistogramRange1[2]={30,200};
float *HistogramRange[1]={&HistogramRange1[0]};

int main()
{
    IplImage *Image1;
    CvHistogram *Histogram1;
    IplImage *HistogramImage1;
    CvPoint Point1;
    CvPoint Point2;

    Image1=cvLoadImage("Riverbank.jpg",0);

    Histogram1 = cvCreateHist(1,&HistogramBins,CV_HIST_ARRAY,HistogramRange);
    HistogramImage1 = cvCreateImage(cvSize(256,300),8,3);
    cvSetZero(HistogramImage1);
    HistogramImage1->origin=1;
    HistogramBinWidth=256/HistogramBins;
    printf("The Bolck Width is : %d\n",HistogramBinWidth);
    cvCalcHist(&Image1,Histogram1);

    printf("Gray Level Values:\n");
    for(int i=0;i<HistogramBins;i++)
    {
        printf("%.f \n",((CvMatND *) Histogram1->bins)->data.fl[i]);
        Point1=cvPoint(i*HistogramBinWidth,0);
        Point2=cvPoint((i+1)*HistogramBinWidth,(int)cvQueryHistValue_1D(Histogram1,i)/50);

        cvRectangle(HistogramImage1,Point1,Point2,CV_RGB(127,127,127));
    }

    cvNamedWindow("Histogram",1);
    cvNamedWindow("Riverbank",1);
    cvShowImage("Riverbank",Image1);
    cvShowImage("Histogram",HistogramImage1);
    cvWaitKey(0);
}

執行結果:

因此,這邊就用cvRectangle()的方式來表達,對於CvHistogram這個資料結構,為什麼可以任意定義上下界以及區塊的個數呢？因為CvHistogram使用的是Look-up table(LUT)的方式,也就是查表法,開一個256大小空間的陣列,利用比例的縮放,縮放的數據存入Look-up table裡面,在利用索引的方式對應,而使用Look-up table最大的缺點為,無法取得很精確的數據,由於它是被比例縮放過的,因此除了範圍為256的可以求的完整的數值,LUT的用法會在後面使用到.

對於其他維度的直方圖,它的做法如下

三維直方圖擷取
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int Histogram3DBins[3] = {256,256,256};
float HistogramRange1[6]={0,255,0,255,0,255};
float *HistogramRange[3]={&HistogramRange1[0],&HistogramRange1[2],&HistogramRange1[4]};

void Print3DHistogram(CvHistogram *Histogram,int BinSize);

int main()
{
    CvHistogram *Histogram1;
    IplImage *Image1=cvLoadImage("Riverbank.jpg",1);
    IplImage *RedImage=cvCreateImage(cvGetSize(Image1),8,1);
    IplImage *GreenImage=cvCreateImage(cvGetSize(Image1),8,1);
    IplImage *BlueImage=cvCreateImage(cvGetSize(Image1),8,1);
    IplImage *ImageArray[3]={RedImage,GreenImage,BlueImage};

    cvSplit(Image1,BlueImage,GreenImage,RedImage,0);
    Histogram1 = cvCreateHist(3,Histogram3DBins,CV_HIST_SPARSE,HistogramRange);

    cvCalcHist(ImageArray,Histogram1);

    printf("3D Historgram Data\n");
    Print3DHistogram(Histogram1,256);

    cvNamedWindow("Riverbank",1);
    cvShowImage("Riverbank",Image1);

    cvWaitKey(0);
}

void Print3DHistogram(CvHistogram *Histogram,int BinSize)
{
    for(int i=0;i<BinSize;i++)
    {
        for(int j=0;j<BinSize;j++)
        {
            for(int k=0;k<BinSize;k++)
            {
                if(cvQueryHistValue_3D(Histogram,i,j,k)>10)
                {
                    printf("%.f\n",cvQueryHistValue_3D(Histogram,i,j,k));
                }
            }
        }
    }
}

執行結果:

上面是三個維度的統計直方圖,分別為R維度,G維度,B維度,並且分別被量化成256等分的區塊,同樣的,也是使用LUT的方法,由於這個統計直方圖非常的大,它所佔的區域為256*256*256的大小,而實際上累積的分佈只有一點點,所以說,對於高維度的統計直方圖就要用到稀疏矩陣,要不然會浪費大量的記憶體空間,但是對於多維度的圖表呈現,OpenCV最多也只能支援到一維的方式,二維的方式可以用OpenGL(glut)呈現,當然高的維度本來就很難用視覺化的方式呈現.

cvCreateHist()
初始化CvHistogram資料結構,可以選擇密集矩陣(CvMatND)CV_HIST_ARRAY以及稀疏矩陣(CvSparseMat)CV_HIST_SPARSE,第一個引數為維度的選擇,第二個引數為要將直方圖切割為多少區塊,第三個引數為選擇cvCreateHist()函數的參數或代號,第四個引數為每一個維度的上下界
cvCreateHist(輸入int型別直方圖維度,輸入int型別直方圖區塊數,輸入參數或代號,輸入每個維度的上下界數據)

cvCalcHist()
計算直方圖累積數據,第一個引數為輸入目標IplImage圖形陣列資料結構,第二個引數輸出為CvHistogram資料結構
cvCalcHist(輸入目標IplImage圖形陣列資料結構,輸出為CvHistogram資料結構)

cvQueryHistValue_1D()
讀取CvHistogram資料結構一維空間的圖形直方圖資料,它在OpenCV的"cvcompat.h"函式庫被定義為

#define cvQueryHistValue_1D( hist, idx0 ) ((float)cvGetReal1D( (hist)->bins, (idx0)))

因此回傳的是float型別的資料,輸入CvHistogram資料結構,以及輸入一維資料結構的Index索引數據
cvQueryHistValue_1D(輸入CvHistogram資料結構,輸入一維int型別數據索引)

cvQueryHistValue_2D()
讀取CvHistogram資料結構二維空間的圖形直方圖資料,它在OpenCV的"cvcompat.h"函式庫被定義為

#define cvQueryHistValue_2D( hist, idx0, idx1 ) ((float)cvGetReal2D( (hist)->bins, (idx0), (idx1)))

回傳float型別的資料,第一個引數為輸入CvHistogram資料結構,第二個引數為輸入二維資料結構的x軸索引數據,第三個引數為y軸索引數據
cvQueryHistValue_2D(輸入CvHistogram資料結構,輸入二維int型別x軸數據索引,輸入二維int型別y軸數據索引)

cvQueryHistValue_3D()
讀取CvHistogram資料結構三維空間的圖形直方圖資料,它在OpenCV的"cvcompat.h"函式庫被定義為

#define cvQueryHistValue_3D( hist, idx0, idx1, idx2 ) ((float)cvGetReal3D( (hist)->bins, (idx0), (idx1), (idx2)))

回傳float型別的資料,第一個引數為輸入CvHistogram資料結構,第二個引數為輸入三維資料結構的x軸索引數據,第三個引數為輸入y軸索引數據,第四個引數為輸入z軸索引數據
cvQueryHistValue_3D(輸入CvHistogram資料結構,輸入三維int型別x軸數據索引,輸入三維int型別y軸數據索引,輸入三維int型別z軸數據索引)

OpenCV統計應用-極端值,cvReduce

2008-10-11T15:45:00.014+08:00

在統計模型裡面,總是會有一些特別高或特別低的數值,這個數值已經脫離了統計分配的假設,不符合統計的模型,通常會有這種極大或是極小的極端值(outlier)有可能是因為量測誤差,或是在統計模型中那0.0001的機率發生,因此,在統計學裡,挑出極端值來做最精簡的量測也是很重要的,就好比一群低收入戶的人裡面出現了一位身價上億的人而拉高了一群人的平均所得,使得低收入戶的族群變成中高收入的族群,這邊要抓取那個最大最小值就用到了cvMinMaxLoc()的函式啦.在圖片裡面,也許很難有如此極端的事件產生,cvMinMaxLoc()可以處理一般的數據,亦可以處理圖片找出它最大最小值以及它的位置.

cvMinMaxLoc()的實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);

    double MinValue;
    double MaxValue;

    CvPoint MinLocation;
    CvPoint MaxLocation;

    cvSetImageCOI(Image1,1);
    cvMinMaxLoc(Image1,&MinValue,&MaxValue,&MinLocation,&MaxLocation);

    printf("The Min number is : %.f\n",MinValue);
    printf("The position is : ( %d , %d )\n",MinLocation.x,MinLocation.y);
    printf("The Max number is : %.f\n",MaxValue);
    printf("The position is : (%d , %d )\n",MaxLocation.x,MaxLocation.y);

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvWaitKey(0);
}

執行結果:

上面的結果是抓出這張圖片最大最小值的數據,選擇綠色這個通道,雖然這已經不算是極端值的意義了,不過它仍然是可以對一般圖形做處理,cvMinMaxLoc()可以同時找出最大最小值,也可以指出最大值的位置跟最小值的位置,而cvMinMaxLoc()必須對單通道做處理因此必須要用,cvSetImageCOI選定顏色,也可以支援ROI,甚至,cvMinMaxLoc()可以用遮罩的方式實作,使用的方法如下

#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);
    IplImage *MaskImage1=cvLoadImage("grotto_Threshold.bmp",0);
    double MinValue;
    double MaxValue;

    CvPoint MinLocation;
    CvPoint MaxLocation;

    cvSetImageCOI(Image1,1);
    cvMinMaxLoc(Image1,&MinValue,&MaxValue,&MinLocation,&MaxLocation,MaskImage1);

    printf("The Min number is : %.f\n",MinValue);
    printf("The position is : ( %d , %d )\n",MinLocation.x,MinLocation.y);
    printf("The Max number is : %.f\n",MaxValue);
    printf("The position is : (%d , %d )\n",MaxLocation.x,MaxLocation.y);

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvWaitKey(0);
}

執行結果:

上面的程式也只對grotto_Threshold.bmp白色的部份做處理,從白色區域找出它的最大最小值以及它的位置,跟前面的程式差不多.

cvReduce()為將矩陣維度降低為向量的維度,也就勢將該行或該列做加總,平均,找出最大,最小值,cvReduce()可使用的參數或代號如下

#define CV_REDUCE_SUM 0
#define CV_REDUCE_AVG 1
#define CV_REDUCE_MAX 2
#define CV_REDUCE_MIN 3

而它的數據要從圖形縮成行或縮成列,則輸入的向量就必須要跟它的行或列相等長度,而且是一維的形式,下面這個是求出每一列的RGB值的平均數,並且用圖表的方式

cvReduce()平均數向量
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);
    IplImage *AVGImage=cvCreateImage(cvSize(Image1->width,255),IPL_DEPTH_8U,3);
    CvMat *AVGVector=cvCreateMat(1,Image1->width,CV_32FC3);

    cvReduce(Image1,AVGVector,0,CV_REDUCE_AVG);

    for(int i=0;i<Image1->width;i++)
    {
        printf("%f\t%f\t%f\t\n",cvGet1D(AVGVector,i).val[0],cvGet1D(SumVector,i).val[1],cvGet1D(SumVector,i).val[2]);
        if(i!=Image1->width-1)
        {
            cvLine(AVGImage,cvPoint(i,(int)cvGet1D(AVGVector,i).val[0]),cvPoint(i+1,(int)cvGet1D(AVGVector,i+1).val[0]),CV_RGB(0,0,255));
            cvLine(AVGImage,cvPoint(i,(int)cvGet1D(AVGVector,i).val[1]),cvPoint(i+1,(int)cvGet1D(AVGVector,i+1).val[1]),CV_RGB(0,255,0));
            cvLine(AVGImage,cvPoint(i,(int)cvGet1D(AVGVector,i).val[2]),cvPoint(i+1,(int)cvGet1D(AVGVector,i+1).val[2]),CV_RGB(255,0,0));
        }

    }

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvNamedWindow("AVG Image",1);
    cvShowImage("AVG Image",AVGImage);

    cvWaitKey(0);
}

執行結果:

cvReduce()可以支援多通道,而AVGVector則是通道三的向量,然後圖表的呈現以cvLine()實作,cvReduce()第一個引數為目標要被統計的圖片,第二個引數為要被計算出來的向量,而第二個引數有一些特殊的規定,它的使用條件如下

1.通道數必須與輸入圖片的數目相同
2.所需要計算的列或欄的長度必須要與輸入圖片
3.一定要用浮點數型別來輸入,如CV_32FC3,CV64FC3,CV_32FC1等,由於uchar型別在做像是CV_REDUCE_SUM時數據一定會超過,因此cvReduce()會做檢查驗證

而且而CV_REDUCE_SUM則是在二值化的圖片的計算比較具有意義,它可以求出二值化圖片的投影向量,並且以圖表的方式呈現.下面就是以列投影以及欄投影的方式來做圖表的呈現.

二值化圖形投影計算
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int main()
{
    IplImage *Image1=cvLoadImage("grotto_Threshold.bmp",0);
    IplImage *RowProjectionImage=cvCreateImage(cvSize(Image1->width,300),IPL_DEPTH_8U,3);
    IplImage *ColumnProjectionImage=cvCreateImage(cvSize(300,Image1->height),IPL_DEPTH_8U,3);
    CvMat *RowSumVector=cvCreateMat(1,Image1->width,CV_32FC1);
    CvMat *ColumnSumVector=cvCreateMat(Image1->height,1,CV_32FC1);

    RowProjectionImage->origin=1;
    cvReduce(Image1,RowSumVector,0,CV_REDUCE_SUM);
    cvReduce(Image1,ColumnSumVector,1,CV_REDUCE_SUM);

    for(int i=0;i<Image1->width;i++)
    {
        cvLine(RowProjectionImage,cvPoint(i,(int)cvGet1D(RowSumVector,i).val[0]/(2*255)),cvPoint(i,0),CV_RGB(255,255,255));
    }
    for(int i=0;i<Image1->height;i++)
    {
        cvLine(ColumnProjectionImage,cvPoint((int)cvGet1D(ColumnSumVector,i).val[0]/(2*255),i),cvPoint(0,i),CV_RGB(255,255,255));
    }

    cvNamedWindow("grotto Threshold",1);
    cvShowImage("grotto Threshold",Image1);

    cvNamedWindow("Row Projection Image",1);
    cvShowImage("Row Projection Image",RowProjectionImage);

    cvNamedWindow("Column Projection Image",1);
    cvShowImage("Column Projection Image",ColumnProjectionImage);
    cvWaitKey(0);
}

執行結果:

上面的程式,由於是二值化的圖片,只有0跟255,所以必須要將累加的數值除以255才能知道它所累積的個數,在將該數目除以二平移,使他在圖表呈現上不會太大,對於投影的實作,可以藉由選取較高的區間,提出二值化裡面有意義的資料區塊,在人臉辨識可以藉由這種方法提取特定臉部器官的特徵.

cvMinMaxLoc()
找出圖片或一組數據中最大值及最小值的數據,以及最大值及最小值的位置,第一個引數為輸入IplImage資料結構或CvMat資料結構,第二個引數為輸出最小值double型別數據,第三個引數為輸出最大值double型別數據,第四個引數為輸出最小值位置CvPoint資料結構,第五個引數為輸出最大值位置CvPoint資料結構.
cvMinMaxLoc(輸入IplImage或CvMat資料結構,輸出double型別最小值數據,輸出double型別最大值數據,輸出最小值CvPoint資料結構,輸出最大值CvPoint資料結構)

cvReduce()
將二維圖形縮減成列向量或欄向量,cvReduce()具有四個參數,分別為CV_REDUCE_SUM總和,CV_REDUCE_AVG平均,CV_REDUCE_MAX最大值,CV_REDUCE_MIN最小值,可以由多通道的輸入,但輸出的行向量或欄向量一定要符合圖形的長度.第一個引數為輸入目標IplImage資料結構或CvMat資料結構,第二個引數為輸出目標CvMat結構向量,第三個引數為縮減為列向量或欄向量,0為欄向量,1為列向量,第四個引數為cvReduce()的參數或代號.
cvReduce(輸入IplImage或CvMat資料結構,輸出CvMat向量資料結構,輸入縮減欄向量或列向量代號,目標參數或代號)

OpenCV統計應用-平均數與標準差

2008-10-10T19:37:00.013+08:00

平均數跟標準差是統計學裡最基本的東西,cvAvg()為平均數的計算,也支援多通道的使用,而它使用的方法如下

簡單cvAvg()實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);
    CvScalar Scalar1;
    Scalar1 = cvAvg(Image1);

    printf("Blue Channel Avg is : %.f\n",Scalar1.val[0]);
    printf("Green Channel Avg is : %.f\n",Scalar1.val[1]);
    printf("Red Channel Avg is : %.f\n",Scalar1.val[2]);

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvWaitKey(0);
}

執行結果:

與cvSum()相同,而它與cvSum()唯一不同的地方為它支援Mask遮罩的操作,cvAvg()將黑白圖遮罩白色的部份做平均數的統計,使用方式如下

簡單cvAvg()遮罩實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);
    IplImage *MaskImage1=cvLoadImage("grotto_Threshold.jpg",0);
    CvScalar Scalar1;
    Scalar1 = cvAvg(Image1,MaskImage1);

    printf("Blue Channel Avg is : %.f\n",Scalar1.val[0]);
    printf("Green Channel Avg is : %.f\n",Scalar1.val[1]);
    printf("Red Channel Avg is : %.f\n",Scalar1.val[2]);

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvWaitKey(0);
}

執行結果:

將之前洞穴的二值化黑白圖片當遮罩對應出白色部分做平均數的計算,而OpenCV也有另外一個算平均數的函式cvAvgSdv(),它可以同時算出平均數跟標準差,也支援Mask遮罩的實作,使用的方式如下

簡單cvAvgSdv()實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);
    CvScalar MeanScalar;
    CvScalar StandardDeviationScalar;
    cvAvgSdv(Image1,&MeanScalar,&StandardDeviationScalar);

    printf("Blue Channel Avg is : %.f\n",MeanScalar.val[0]);
    printf("Blue Channel Standard Deviation is : %.f\n",StandardDeviationScalar.val[0]);
    printf("Green Channel Avg is : %.f\n",MeanScalar.val[1]);
    printf("Green Channel Standard Deviation is : %.f\n",StandardDeviationScalar.val[1]);
    printf("Red Channel Avg is : %.f\n",MeanScalar.val[2]);
    printf("Red Channel Standard Deviation is : %.f\n",StandardDeviationScalar.val[2]);

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvWaitKey(0);
}

執行結果:

上面的是用cvAvgSdv()計算平均數跟標準差,也同樣是用CvScalar資料結構來做數值的接收,再來,就是Mask的部份

簡單cvAvgSdv()遮罩實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);
    IplImage *MaskImage1=cvLoadImage("grotto_Threshold.jpg",0);
    CvScalar MeanScalar;
    CvScalar StandardDeviationScalar;
    cvAvgSdv(Image1,&MeanScalar,&StandardDeviationScalar,MaskImage1);

    printf("Blue Channel Avg is : %.f\n",MeanScalar.val[0]);
    printf("Blue Channel Standard Deviation is : %.f\n",StandardDeviationScalar.val[0]);
    printf("Green Channel Avg is : %.f\n",MeanScalar.val[1]);
    printf("Green Channel Standard Deviation is : %.f\n",StandardDeviationScalar.val[1]);
    printf("Red Channel Avg is : %.f\n",MeanScalar.val[2]);
    printf("Red Channel Standard Deviation is : %.f\n",StandardDeviationScalar.val[2]);

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvNamedWindow("grotto_Threshold",1);
    cvShowImage("grotto_Threshold",MaskImage1);

    cvWaitKey(0);
}

執行結果:

也是跟cvAvg()一樣用單通道二值化影像做遮罩,計算平均數與標準差,也不限於圖形的輸入,實際上任何的資料數據都可以做輸入,下面就是任意一筆資料來做平均數跟標準差的計算

cvMat資料結構計算平均數標準差
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

uchar Array[10]={1,2,3,4,5,6,7,8,9};
int main()
{
    CvMat *Matrix1=cvCreateMat(1,10,CV_8UC1);
    cvSetData(Matrix1,Array,Matrix1->step);

    CvScalar MeanScalar;
    CvScalar StandardDeviationScalar;
    cvAvgSdv(Matrix1,&MeanScalar,&StandardDeviationScalar);

    printf("Matrix1 Avg is : %.f\n",MeanScalar.val[0]);
    printf("Matrix1 Standard Deviation is : %.f\n",StandardDeviationScalar.val[0]);

    system("pause");
}

執行結果:

實際上,很多的函式IplImage資料結構跟CvMat資料結構都是可以共用的,包括前面的cvCountNonZero(),cvSum(),cvAvg()及cvAvgSdv().

cvAvg()
計算IplImage資料結構或CvMat資料結構的平均數,支援Mask遮罩的運算,第一個引數輸入為IplImage資料結構或CvMat資料結構,第二個引數為IplImage資料結構單通道uchar型別的遮罩對應,輸出為CvScalar資料結構
cvAvg(輸入IplImage資料結構或CvMat資料結構,輸入單通道uchar型別IplImage資料結構)

cvAvgSdv()
計算IplImage資料結構或CvMat資料結構的平均數,標準差,支援Mask通道的計算,第一個引數為輸入IplImage資料結構或CvMat資料結構,第二個引數為輸出CvScalar平均數資料結構,第三個引數為輸出CvScalar標準差資料結構,第四個引數為IplImage資料結構單通道uchar型別遮罩的輸入
cvAvgSdv(輸入IplImage資料結構或CvMat資料結構,輸出平均數CvScalar資料結構,輸出標準差CvScalar資料結構,輸入單通道uchar型別IplImage資料結構)

OpenCV統計應用-總和,計算非零

2008-10-07T21:02:00.008+08:00

再來就是一些計數的函式啦,下面簡單的介紹一張圖形,所有數值總和以及圖形上所有非0的計算,在一些圖形數值的量測方面,稍微具有一些實質的意義.

簡單cvSum()使用
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);
    CvScalar Scalar1;
    Scalar1 = cvSum(Image1);

    printf("Blue Channel Sum is : %.f\n",Scalar1.val[0]);
    printf("Green Channel Sum is : %.f\n",Scalar1.val[1]);
    printf("Red Channel Sum is : %.f\n",Scalar1.val[2]);

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvWaitKey(0);
}

原始圖片:

執行結果:

cvSum()是計算所有圖片每個通道數據的總和,因此在這個通道為3的RGB圖,所有的總和都被計算在CvScalar資料結構上面了,cvSum()輸入為IplImage或CvMat資料結構,輸出為CvScalar的資料結構,裡面存放著RGB的總和數值,接著就是cvContNonZero()的使用.

簡單cvCountNonZero()使用
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
int main()
{
    IplImage *Image1=cvLoadImage("grotto_Threshold.jpg",0);

    printf("Thresholding Non Zero is : %d\n",cvCountNonZero(Image1));

    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);
    cvWaitKey(0);
}

原始圖片：

執行結果:

cvCountNonZero()為計算非0的個數,只能使用在單通道的圖形,而且對於這個函式在黑白圖的使用上會比較具有意義,由於黑白圖非0即1因此計算出非0的個數同時也算出了1的個數,而cvCountNonZero(),以及cvSum()也可以支援COI（Color Of Interesting）以及ROI（Region Of Interesting）的使用,使用方式如下

cvCountNonZero()及cvSum()的ROI,COI
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int main()
{
    IplImage *Image1=cvLoadImage("grotto.jpg",1);

    cvSetImageROI(Image1,cvRect(100,100,100,100));

    cvSetImageCOI(Image1,1);
    printf("Blue Channel Non Zero is : %d\n",cvCountNonZero(Image1));
    printf("Blue Channel Sum is : %.f\n",cvSum(Image1).val[0]);

    cvSetImageCOI(Image1,2);
    printf("\nGreen Channel Non Zero is : %d\n",cvCountNonZero(Image1));
    printf("Green Channel Sum is : %.f\n",cvSum(Image1).val[0]);

    cvSetImageCOI(Image1,3);
    printf("\nRed Channel Non Zero is : %d\n",cvCountNonZero(Image1));
    printf("Red Channel Sum is : %.f\n",cvSum(Image1).val[0]);

    cvNamedWindow("grotto(ROI)",1);
    cvShowImage("grotto(ROI)",Image1);
    cvResetImageROI(Image1);
    cvNamedWindow("grotto",1);
    cvShowImage("grotto",Image1);

    cvWaitKey(0);
}

執行結果：

對於一個被選定COI的他會對於特定顏色做統計,cvCountNonZero()以及cvSum()亦是如此,而被選定ROI他會對於特定的區域做計算,因此上面的結果會是被選定的區域,加上被選定的顏色,他所計算出來的統計結果,而COI選定的顏色後,因為是單通道的,所以cvSum()計算出來的結果會被放在CvScalar資料結構val[0]的位置,這邊用cvSum(Image1).val[0]表示直接存取CvScalar資料結構val[0]的數據.

cvSum()
計算所有通道數據的總和,計算出來的結果會回傳放在CvScalar資料結構上,如果是使用單通道或是COI的話會放在CvScalar資料結構val[0]上面.輸入的引數為IplImage資料結構或CvMat資料結構.
cvSum(輸入IplImage或CvMat資料結構)

cvCountNonZero()
計算所有非0數據的總和,只支援單通道或是COI,輸入為IplImage資料結構或CvMat資料結構,輸出為int變數型別.
cvCountNonZero(輸入IplImage或CvMat資料結構)

OpenCV統計應用-基本統計直方圖

2008-10-07T11:11:00.008+08:00

在統計學裡面,最基本影像處理的應用,就是統計直方圖了,它可以將一張圖片的色彩結構統計起來,對於一張圖片的R值,G值,B值或灰階做統計計算,或者是,轉換到另一個色彩空間,將圖形做量化及統計的計算,在OpenCV裡面,也有不少統計直方圖的計算,並且可以將它圖表化,呈現出統計的概況.

基本統計直方圖計算
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int BlueHistogramData[256]={0};
int GreenHistogramData[256]={0};
int RedHistogramData[256]={0};
int GrayLevelHistogramData[256]={0};

int main()
{
    IplImage *Image1=cvLoadImage("waterfall.jpg",1);
    CvScalar Scalar1;
    int BlueColor,GreenColor,RedColor,GrayLevelColor;

    for(int i=0;i<Image1->height;i++)
    {
        for(int j=0;j<Image1->width;j++)
        {

            Scalar1=cvGet2D(Image1,i,j);
            BlueColor=(int)Scalar1.val[0];
            GreenColor=(int)Scalar1.val[1];
            RedColor=(int)Scalar1.val[2];
            GrayLevelColor=(int) (0.299*RedColor+0.587*GreenColor+0.114*BlueColor);

            BlueHistogramData[BlueColor]++;
            GreenHistogramData[GreenColor]++;
            RedHistogramData[RedColor]++;
            GrayLevelHistogramData[GrayLevelColor]++;

        }
    }
    printf("Red\tGreen\tBlue\tGrayLevel\n");
    for(int i=0;i<256;i++)
    {

        printf("%d\t%d\t%d\t%d\t\n",RedHistogramData[i],GreenHistogramData[i],BlueHistogramData[i],GrayLevelHistogramData[i]);
    }

    cvNamedWindow("waterfall",1);
    cvShowImage("waterfall",Image1);
    cvWaitKey(0);
}

原始圖片:

執行結果:

上面的程式碼,就是將RGB值做統計,將R值,G值,B值累加在0~255的範圍裡,因此就可以知道這張圖的色彩分布的結構,灰階圖亦然,而OpenCV亦可以做到圖表的呈現的方式,利用的是繪圖函式來實作,用cvLine()或是cvRectangle()就可以做到統計值方圖的繪製,OpenCV在統計直方圖製作部分有一個CvHistogram的資料結構可以實作,而下面這個,則是直接將上面的程式碼修改繪製上去的統計值方圖製作.

統計直方圖計算與繪製
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

int BlueHistogramData[256]={0};
int GreenHistogramData[256]={0};
int RedHistogramData[256]={0};
int GrayLevelHistogramData[256]={0};

int main()
{
    IplImage *Image1=cvLoadImage("waterfall.jpg",1);
    IplImage *BlueHistogramImage,*GreenHistogramImage,*RedHistogramImage;
    IplImage *GrayLevelHistogramImage;

    BlueHistogramImage=cvCreateImage(cvSize(256,250),IPL_DEPTH_8U,3);
    GreenHistogramImage=cvCreateImage(cvSize(256,250),IPL_DEPTH_8U,3);
    RedHistogramImage=cvCreateImage(cvSize(256,250),IPL_DEPTH_8U,3);
    GrayLevelHistogramImage=cvCreateImage(cvSize(256,250),IPL_DEPTH_8U,3);

    BlueHistogramImage->origin=1;
    GreenHistogramImage->origin=1;
    RedHistogramImage->origin=1;
    GrayLevelHistogramImage->origin=1;

    CvScalar Scalar1;
    int BlueColor,GreenColor,RedColor,GrayLevelColor;

    for(int i=0;i<Image1->height;i++)
    {
        for(int j=0;j<Image1->width;j++)
        {

            Scalar1=cvGet2D(Image1,i,j);
            BlueColor=(int)Scalar1.val[0];
            GreenColor=(int)Scalar1.val[1];
            RedColor=(int)Scalar1.val[2];
            GrayLevelColor=(int) (0.299*RedColor+0.587*GreenColor+0.114*BlueColor);

            BlueHistogramData[BlueColor]++;
            GreenHistogramData[GreenColor]++;
            RedHistogramData[RedColor]++;
            GrayLevelHistogramData[GrayLevelColor]++;

        }
    }
    for(int i=0;i<255;i++)
    {
        cvLine(BlueHistogramImage,cvPoint(i,BlueHistogramData[i]/10),cvPoint(i,0),CV_RGB(0,0,255));
        cvLine(GreenHistogramImage,cvPoint(i,GreenHistogramData[i]/10),cvPoint(i,0),CV_RGB(0,255,0));
        cvLine(RedHistogramImage,cvPoint(i,RedHistogramData[i]/10),cvPoint(i,0),CV_RGB(255,0,0));
        cvLine(GrayLevelHistogramImage,cvPoint(i,GrayLevelHistogramData[i]/10),cvPoint(i,0),CV_RGB(128,128,128));
    }

    cvNamedWindow("waterfall",1);
    cvShowImage("waterfall",Image1);
    cvNamedWindow("Blue Histogram",1);
    cvShowImage("Blue Histogram",BlueHistogramImage);
    cvNamedWindow("Green Histogram",1);
    cvShowImage("Green Histogram",GreenHistogramImage);
    cvNamedWindow("Red Histogram",1);
    cvShowImage("Red Histogram",RedHistogramImage);
    cvNamedWindow("GrayLevel Histogram",1);
    cvShowImage("GrayLevel Histogram",GrayLevelHistogramImage);

    cvWaitKey(0);
}

執行結果:

在繪製圖形的時候,亦要開啟圖表的空間,IplImage資料結構,在將IplImage資料結構的origin設為1,為Window XP預設圖形起始點,這部份在"資料結構操作與運算-IplImage資料結構"有詳細的介紹,再來就是讀取RGB值做統計,以及使用cvLine()做繪製,cvLine()的繪製方法則是在"OpenCV繪圖的實作-cvLine,cvRectangle"有介紹到,而CvHistogram資料結構也是使用類似上面的原理,在OpenCV裡面,統計直方圖的部分還是需要自己動手畫,沒有很方便的函式可以直接實作出視覺化圖形的呈現.

OpenCV線性代數-cvGEMM通用矩陣乘法

2008-09-14T10:49:00.025+08:00

OpenCV提供了一個通用矩陣乘法的函式,cvGEMM(),代表的是GEneralized Matrix Multiplication,cvGEMM()可以處理線性代數方面許多的乘法運算,cvmMul()矩陣乘法這個函式就是從這裡來的,cvmMul()在"cvcompat.h"及"cxcore.h"這兩個函式庫被定義為

#define cvmMul(src1,src2,dst) cvMatMulAdd(src1,src2,0,dst)
#define cvMatMul(src1,src2,dst) cvMatMulAdd((src1),(src2),NULL,(dst))
#define cvMatMulAdd(src1,src2,src3,dst) cvGEMM((src1),(src2),1.,(src3),1.,(dst),0)

先簡單的介紹cvmMul()跟cvMatMul()的用法

cvmMul()及cvMatMul()的實作
#include <cv.h>
#include <stdio.h>

float array1[]={3,1,2,0,1,5};
float array2[]={4,3,1,1,6,0};

void PrintMatrix(CvMat *Matrix,int Rows,int Cols,int Channels);
int main()
{
    CvMat *A=cvCreateMat(2,3,CV_32FC1);
    CvMat *B=cvCreateMat(3,2,CV_32FC1);
    CvMat *ResultMatrix=cvCreateMat(2,2,CV_32FC1);

    cvSetData(A,array1,A->step);
    cvSetData(B,array2,B->step);

    printf("cvmMul():\n");
    cvmMul(A,B,ResultMatrix);
    PrintMatrix(ResultMatrix,2,2,1);

    printf("\ncvMatMul():\n");
    cvMatMul(A,B,ResultMatrix);
    PrintMatrix(ResultMatrix,2,2,1);

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

執行結果:

所以說cvmMul()跟cvMatMul()是相同的,而他們都源自於cvMatMulAdd()這個函式,cvMatMulAdd()這個函式的前身而是cvGEMM()定義而來的,因此cvGEMM()包含了許多乘法運算的應用,下面這個為cvmMul()跟cvMatMul()的計算方式

因此,cvmMul()跟cvMatMul()是簡單的乘法運算,而cvMatMulAdd()則是如何呢？下面就是cvMatMulAdd()的使用方法

cvMatMulAdd()實作
#include <cv.h>
#include <stdio.h>

float array1[]={3,1,2,0,1,5};
float array2[]={4,3,1,1,6,0};
float array3[]={-10,-2,-20,10};

void PrintMatrix(CvMat *Matrix,int Rows,int Cols,int Channels);
int main()
{
    CvMat *A=cvCreateMat(2,3,CV_32FC1);
    CvMat *B=cvCreateMat(3,2,CV_32FC1);
    CvMat *C=cvCreateMat(2,2,CV_32FC1);
    CvMat *ResultMatrix=cvCreateMat(2,2,CV_32FC1);

    cvSetData(A,array1,A->step);
    cvSetData(B,array2,B->step);
    cvSetData(C,array3,C->step);

    printf("cvMatMulAdd():\n");
    cvMatMulAdd(A,B,C,ResultMatrix);
    PrintMatrix(ResultMatrix,2,2,1);

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

執行結果:

cvMatMulAdd()是簡單的相乘在相加的函式,只要符合矩陣大小的規格,它就可以對三個矩陣做相乘在相加的動作,對於上面的程式碼,它的計算方式如下

這是將矩陣相乘後的結果在加上另外一個矩陣,這跟cvGEMM()什麼關係呢？cvGEMM()函式裡面包含了很多種計算方式,而cvmMul(),cvMatMul(),cvMatMulAdd()則是用#define來擷取它的算法,而cvGEMM()的使用方式如下

簡單cvGEMM()函式實作
#include <cv.h>
#include <stdio.h>

float array1[]={3,1,2,0,1,5};
float array2[]={4,3,1,1,6,0};
float array3[]={-10,-2,-20,10};

void PrintMatrix(CvMat *Matrix,int Rows,int Cols,int Channels);
int main()
{
    CvMat *A=cvCreateMat(2,3,CV_32FC1);
    CvMat *B=cvCreateMat(3,2,CV_32FC1);
    CvMat *C=cvCreateMat(2,2,CV_32FC1);
    CvMat *ResultMatrix=cvCreateMat(2,2,CV_32FC1);
    double alpha =0.5;
    double beta =0.1;

    cvSetData(A,array1,A->step);
    cvSetData(B,array2,B->step);
    cvSetData(C,array3,C->step);

    printf("cvGEMM():\n");
    cvGEMM(A,B,alpha,C,beta,ResultMatrix);
    PrintMatrix(ResultMatrix,2,2,1);

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

執行結果:

因此,cvGEMM()利用#define的方式,將不必要的引數歸零或設為空或是將純量積設為1,所以才有cvmMul()及cvMulAdd()的變化,cvGEMM()的程式碼,它的數學運算方式如下所示

cvGEMM()也可以做矩陣轉置的運算,它將矩陣轉置的定義如下

#define CV_GEMM_A_T 1
#define CV_GEMM_B_T 2
#define CV_GEMM_C_T 4

如果要轉置A矩陣的話就用CV_GEMM_A_T,如果要轉置A跟B矩陣的話就用CV_GEMM_A_T+CV_GEMM_B_T(1+2=3),如果要轉置A,B,C矩陣的話就要用CV_GEMM_A_T+CV_GEMM_B_T+CV_GEMM_C_T(1+2+4=7),它可以用參數來表示,也可以用參數代號表示,只要將參數代號的數值累加即可,參數代號用1,2,4這種方式是用類似linux權限控制規劃的方法,在數學上的意義來講就是2ⁿ的加法可以表達出所有數字的組合,n=0,1,2...,在離散數學裡面可以找到這種計算規則,以下就是這個轉置運算的程式碼.

cvGEMM()轉置實作
#include <cv.h>
#include <stdio.h>

float array1[]={3,1,2,0,1,5};
float array2[]={4,3,1,1,6,0};
float array3[]={-70,-4,15,5,15,5,35,-7,20};

void PrintMatrix(CvMat *Matrix,int Rows,int Cols,int Channels);
int main()
{
    CvMat *A=cvCreateMat(2,3,CV_32FC1);
    CvMat *B=cvCreateMat(3,2,CV_32FC1);
    CvMat *C=cvCreateMat(3,3,CV_32FC1);
    CvMat *ResultMatrix=cvCreateMat(3,3,CV_32FC1);
    double alpha =0.5;
    double beta =0.1;

    cvSetData(A,array1,A->step);
    cvSetData(B,array2,B->step);
    cvSetData(C,array3,C->step);

    printf("cvMatMulAdd():\n");
    cvGEMM(A,B,alpha,C,beta,ResultMatrix,CV_GEMM_A_T+CV_GEMM_B_T+CV_GEMM_C_T);
    PrintMatrix(ResultMatrix,3,3,1);

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

執行結果:

這邊CV_GEMM_A_T+CV_GEMM_B_T+CV_GEMM_C_T的參數可以直接用7來代替,上面的程式碼將A,B,C三個矩陣都轉置了,而他的計算方式如下

cvMatMul()
與cvmMul()相同,做簡單的矩陣乘法運算,輸入與輸出的矩陣大小要符合矩陣運算的規格,第一個引數為輸入CvMat資料結構,第二個引數為輸入CvMat資料結構,第三個引數為輸出CvMat資料結構
cvMatMul(輸入CvMat資料結構,輸入CvMat資料結構,輸出CvMat資料結構)

cvMulAdd()
先相乘後相加的矩陣運算,輸入相乘與相加的矩陣及輸出的矩陣必須要符合矩陣運算規格,第一個引數為輸入CvMat相乘矩陣資料結構,第二個引數為輸入CvMat相乘矩陣資料結構,第三的引數為輸入CvMat相加矩陣資料結構,第四個引數為CvMat輸出資料結構
cvMulAdd(輸入CvMat資料結構,輸入CvMat資料結構,輸入CvMat資料結構,輸出CvMat資料結構)

cvGEMM()
矩陣乘法通用運算,具有純量積,矩陣相乘以及矩陣相加,矩陣轉置的功能,包含所有的矩陣運算,並且有一些函式是由cvGEMM()縮減而來如cvmMul(),cvMatMul(),cvMulAdd()等,可以藉由將一些數值設為空值縮減運算的式子,並且可以指定目標矩陣做轉置運算,轉置的規則則由CV_GEMM_A_T(1),CV_GEMM_B_T(2),CV_GEMM_C_T(4)的參數而來,如果要轉置B跟C矩陣則要用CV_GEMM_B_T+CV_GEMM_C_T參數代號則為6(2+4),其他轉置規則依此類推.第一個引數為輸入相乘CvMat資料結構,第二個引數為輸入相乘CvMat資料結構,第三個引數為輸入double型別相乘結構的純量積,第四個引數為輸入相加CvMat資料結構純量積,第五個引數為輸入double型別相加結構的純量積,第六個引數為輸出CvMat資料結構,第七個引數為輸入cvGEMM()的轉置矩陣參數
cvGEMM(輸入CvMat相乘資料結構,輸入CvMat相乘資料結構,輸入double型別相乘純量積,輸入CvMat相加資料結構,輸入double型別相加純量積,輸出CvMat結果資料結構,目標參數或代號)

OpenCV線性代數-正交投影向量計算

2008-09-11T15:06:00.017+08:00

向量空間,為一個可被多個向量基底線性組合的向量集合,也就是有許多的向量基底它可以用任意數據的線性組合產生另外一個向量,存在的那個向量就被稱為是向量空間的一部分,在這邊,要用一個向量空間來計算另外一個向量的正交投影,正交投影表示的是存在向量空間的一個向量它是另外一個向量的投影向量,而且該向量減去它的投影向量會跟向量空間正交(垂直),存在於向量空間的那個向量就稱為正交投影向量(orthogonal projection vector)

由這張圖可以看出,W為一個向量空間,v為一個向量,而proj_v為v的正交投影向量,v-proj_v則會與W向量空間正交,正交投影向量的計算方式如下

A為W空間向量基底的集合,為一個行向量(column vector)的矩陣,v則是要計算正交投影的目標向量,而它的程式計算,則是下面的方式

正交投影向量計算
#include <cv.h>
#include <stdio.h>

float array1[]={1,1,2,-1,-2,0,1,1};
float array2[]={1,2,3,4};

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

int main()
{
    CvMat *A=cvCreateMat(4,2,CV_32FC1);
    CvMat *v=cvCreateMat(4,1,CV_32FC1);
    CvMat *Transpose_A=cvCreateMat(2,4,CV_32FC1);
    CvMat *MulTransposed_A=cvCreateMat(2,2,CV_32FC1);
    CvMat *tempMatrix1=cvCreateMat(4,2,CV_32FC1);
    CvMat *tempMatrix2=cvCreateMat(2,1,CV_32FC1);
    CvMat *proj_v=cvCreateMat(4,1,CV_32FC1);
    cvSetData(A,array1,A->step);
    cvSetData(v,array2,v->step);

    cvTranspose(A,Transpose_A);
    cvMulTransposed(A,MulTransposed_A,1);
    cvInvert(MulTransposed_A,MulTransposed_A);
    cvmMul(A,MulTransposed_A,tempMatrix1);
    cvmMul(Transpose_A,v,tempMatrix2);
    cvmMul(tempMatrix1,tempMatrix2,proj_v);

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

執行結果:

這裡用到了前面許多運算的函數,cvTranspose(),cvInvert(),cvmMul()以及cvMulTransposed(),cvMulTransposed()轉置相乘的乘積第一個引數為輸入CvMat矩陣向量,第二個引數為輸出CvMat矩陣向量,第三個引數為計算方式,分別為原矩陣乘以原矩陣轉置或是原矩陣轉置在乘以原矩陣的選擇,0為原矩陣與原矩陣轉置相乘,1為原矩陣轉置乘以原矩陣.正交投影的計算方法如下.

cvMulTransposed()
計算原矩陣與原矩陣轉置相乘,或是原矩陣轉置與原矩陣的乘積,因為原矩陣轉置相乘為對稱矩陣,因此可以用一些特殊算法來加快矩陣相乘的乘積,cvMulTransposed()第一個引數為輸入CvMat原矩陣資料結構,第二個引數為輸出CvMat資料結構,第三個引數為cvMulTransposed()這函式的參數,0為原矩陣與原矩陣轉置相乘,1為原矩陣轉置與原矩陣的乘積,第四個引數為輸入CvMat資料結構,將原矩陣與第四個引數的矩陣做相減,再做原矩陣轉置相乘原矩陣或是原矩陣相乘原矩陣轉置的運算
cvMulTransposed(輸入CvMat原矩陣資料結構,輸出CvMat資料結構,0或1的參數,輸入CvMat矩陣相減資料結構)

OpenCV線性代數-向量的表達,內積與外積運算

2008-09-06T12:58:00.027+08:00

接下來就簡單的介紹向量,以及向量的計算了,在OpenCV內,向量的表達有很多種方式,向量是以CvMat資料結構為基準,實際上,矩陣本身就可以看成是一群的向量,它可以切割成行向量,列向量,在OpenCV內,甚至可以以通道為單位做為一個向量的表達,在一般的圖形裡,<R,G,B>格式的三通道本身就可以當做是一種向量的形式.

在下面就簡單的做內積跟外積的運算,並且用不同的方式來表達向量的運作.下面這個是用列向量的形式來做表達

列向量內積與外積的計算
#include <cv.h>
#include <stdio.h>

float Array1[]={1,2,3};
float Array2[]={4,5,6};
void PrintMatrix(CvMat *Matrix,int Rows,int Cols,int Channels);

int main()
{
    CvMat *Matrix1=cvCreateMat(1,3,CV_32FC1);
    CvMat *Matrix2=cvCreateMat(1,3,CV_32FC1);
    CvMat *ResultMatrix=cvCreateMat(1,3,CV_32FC1);

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

    printf("The DotProduct is : %.f\n\n",cvDotProduct(Matrix1,Matrix2));

    printf("The CrossProduct is :\n");
    cvCrossProduct(Matrix1,Matrix2,ResultMatrix);
    PrintMatrix(ResultMatrix,1,3,1);

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

執行結果:

再來就是以行向量的方式來表達,並做內積與外積的運算.

行向量內積與外積的計算
#include <cv.h>
#include <stdio.h>

float Array1[]={1,2,3};
float Array2[]={4,5,6};
void PrintMatrix(CvMat *Matrix,int Rows,int Cols,int Channels);

int main()
{
    CvMat *Matrix1=cvCreateMat(3,1,CV_32FC1);
    CvMat *Matrix2=cvCreateMat(3,1,CV_32FC1);
    CvMat *ResultMatrix=cvCreateMat(3,1,CV_32FC1);

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

    printf("The DotProduct is : %.f\n\n",cvDotProduct(Matrix1,Matrix2));

    printf("The CrossProduct is :\n");
    cvCrossProduct(Matrix1,Matrix2,ResultMatrix);
    PrintMatrix(ResultMatrix,3,1,1);

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

執行結果:

接著是通道的方法,並且將通道的向量來做內外積的運算.

以通道為單位的計算方式
#include <cv.h>
#include <stdio.h>

float Array1[]={1,2,3};
float Array2[]={4,5,6};
void PrintMatrix(CvMat *Matrix,int Rows,int Cols,int Channels);

int main()
{
    CvMat *Matrix1=cvCreateMat(1,1,CV_32FC3);
    CvMat *Matrix2=cvCreateMat(1,1,CV_32FC3);
    CvMat *ResultMatrix=cvCreateMat(1,1,CV_32FC3);

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

    printf("The DotProduct is : %.f\n\n",cvDotProduct(Matrix1,Matrix2));

    printf("The CrossProduct is :\n");
    cvCrossProduct(Matrix1,Matrix2,ResultMatrix);
    PrintMatrix(ResultMatrix,1,1,3);

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

執行結果:

內積與外積的計算方式

從上面可以看的出來,向量的計算方式只有一些些微的數據有變化,在cvCreateMat()上要改變格式,在PrintMatrix()上要改變列印方式,其他的地方都是一樣,而輸出的結果也影響不大,在OpenCV上,向量的實作是實現在矩陣上,而它使用的方式也是很自由的,用這三種方式表達向量都是可行的.

cvDotProduct()
向量內積的計算,第一個引數為輸入CvMat資料結構的向量,第二個引數為輸入CvMat資料結構的向量,輸入的向量每個維度的個數必須要相等,且必須具備為向量輸入的格式,輸出為一個double型別的數據為內積計算的結果.
double cvDotProduct(輸入CvMat資料結構向量,輸入CvMat資料結構向量)

cvCrossProduct()
向量外積的計算,必須要是三維向量,且必須要為向量具備的輸入格式,第一個引數為輸入CvMat三維向量資料結構,第二個引數為輸入CvMat三維向量資料結構第三個引數為輸出CvMat三維向量外積計算結果.
cvCrossProduct(輸入CvMat三維向量資料結構,輸入CvMat三維向量資料結構,輸出CvMat三維向量資料結構)

OpenCV線性代數-cvNorm向量長度,距離計算

2008-09-03T21:21:00.015+08:00

cvNorm()這個函數,可以用來做為特徵值比對的公式,也就是在資訊檢索這個領域內,將兩張圖形的特徵值擷取出來,利用紋理值擷取,統計直方圖,小波轉換,或是輪廓資訊等方法,在用裡面cvNorm()內建的CV_L1,CV_L2,CV_C參數來做圖形特徵值比對,cvNorm()的參數定義如下

#define CV_C 1
#define CV_L1 2
#define CV_L2 4
#define CV_DIFF_C (CV_DIFF CV_C)
#define CV_DIFF_L1 (CV_DIFF CV_L1)
#define CV_DIFF_L2 (CV_DIFF CV_L2)
#define CV_RELATIVE_C (CV_RELATIVE | CV_C)
#define CV_RELATIVE_L1 (CV_RELATIVE | CV_L1)
#define CV_RELATIVE_L2 (CV_RELATIVE | CV_L2)

CV_L1是一般的長度計算,CV_L2則是論文常用的歐幾里德距離,CV_C則是拿絕對值後的最大值,在圖形檢索比對這方面cvNorm()還算是不錯用的函式,cvNorm()還提供了遮罩的功能,可以將不必要的特徵資訊做掩蓋的功能,遮罩的部份就要參考前面"圖形的Mask遮罩實作"這地方,與前面圖形Mask使用原理相同.

在線性代數的部分Norm,為向量長度的意思,最常見的長度計算方式為1-Norm將所有向量數據取絕對值的總和,2-Norm為歐機里德距離,將所有數據平方在開根號,∞-Norm為將所有向量數據絕對值後取對大的那個,當然,這分別代表著CV_L1,CV_L2,CV_C的使用,在下面則用簡單的方式做cvNorm()的長度計算.

cvNorm()長度的計算
#include <cv.h>
#include <stdio.h>

float Array2[]={6,-1,2};

int main()
{
    CvMat *Matrix2=cvCreateMat(1,3,CV_32FC1);
    cvSetData(Matrix2,Array2,Matrix2->step);

    printf("1-norm is : %.1f\n",cvNorm(Matrix2,NULL,CV_L1));
    printf("2-norm is : %.1f\n",cvNorm(Matrix2,NULL,CV_L2));
    printf("Infininte norm is : %.1f\n",cvNorm(Matrix2,NULL,CV_C));

    system("pause");
}

執行結果:

cvNorm()的輸出為一個double型別的數據,上面的程式可以簡單的看出來cvNorm()對長度的計算,要計算Matrix2的長度,則是要把第二個引數設置為空,而第二個引數的輸入則是代表著距離,在下面的程式碼範例會介紹到,對於長度的計算公式則是如下所示.

再來是長度計算Mask使用的範例,利用一個遮罩將不想被計算的數據給取消掉,0代表不想使用的數據,1代表要被計算的數據.

cvNorm()遮罩與長度的計算
#include <cv.h>
#include <stdio.h>

float Array2[]={6,-1,2};
uchar MaskArray[]={1,1,0};

int main()
{
    CvMat *Matrix2=cvCreateMat(1,3,CV_32FC1);
    CvMat *MaskMatrix=cvCreateMat(1,3,CV_8UC1);

    cvSetData(Matrix2,Array2,Matrix2->step);
    cvSetData(MaskMatrix,MaskArray,MaskMatrix->step);

    printf("1-norm is : %.1f\n",cvNorm(Matrix2,NULL,CV_L1,MaskMatrix));
    printf("2-norm is : %.1f\n",cvNorm(Matrix2,NULL,CV_L2,MaskMatrix));
    printf("Infininte norm is : %.1f\n",cvNorm(Matrix2,NULL,CV_C,MaskMatrix));

    system("pause");
}

執行結果:

因此,實際上做長度計算的只有6跟-1這兩個數據,在這邊,遮罩的使用CvMat資料結構必須要用單通道非浮點數的數據型別,如CV_32FC1及CV_64FC1都是不可用來作遮罩矩陣的參數宣告.而矩陣也不可以用float及double的型別.遮罩長度的計算方式如下

再來,就是距離的計算了,它可以比較兩個向量的距離大小,也可以用許多向量的距離來做分數(Score)的參考目標,也就是說,在做圖形比對的時候,特徵值距離越小的兩張圖片某些特殊的特徵就越接近(看是用什麼方法擷取特徵值),因此,在一個龐大的圖形特徵資料庫中,將目標圖形的特徵值對每一筆資料做比對,可以藉由取得最小距離來找出最相似的圖形,距離的簡單例子就如下

cvNorm()距離的計算
#include <cv.h>
#include <stdio.h>

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

    printf("1-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_DIFF_L1));
    printf("2-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_DIFF_L2));
    printf("Infininte norm distance is : %.1f\n\n",cvNorm(Matrix1,Matrix2,CV_DIFF_C));

    printf("Relative 1-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_RELATIVE_L1));
    printf("Relative 2-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_RELATIVE_L2));
    printf("Relative infininte norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_RELATIVE_C));

    system("pause");
}

執行結果:

同樣的,它跟長度的計算有點像,比較不同的地方是將兩個向量個別值做相減在對它做1-Norm,2-Norm或是Infininte-Norm的處理,而2-Norm則是論文常用的歐幾里德距離,CV_RELATIVE_L1,CV_RELATIVE_L2,CV_RELATIVE_C則是取它們相對差異的距離.以下則是它們的計算方式.

再來就是遮罩的部份了,實作方式也跟長度的做法差不多.

cvNorm()遮罩與距離的計算
#include <cv.h>
#include <stdio.h>

float Array1[]={5,3,1};
float Array2[]={6,-1,2};
uchar MaskArray[]={1,1,0};

int main()
{
    CvMat *Matrix1=cvCreateMat(1,3,CV_32FC1);
    CvMat *Matrix2=cvCreateMat(1,3,CV_32FC1);
    CvMat *MaskMatrix=cvCreateMat(1,3,CV_8UC1);

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

    printf("1-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_DIFF_L1,MaskMatrix));
    printf("2-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_DIFF_L2,MaskMatrix));
    printf("Infininte norm distance is : %.1f\n\n",cvNorm(Matrix1,Matrix2,CV_DIFF_C,MaskMatrix));

    printf("Relative 1-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_RELATIVE_L1,MaskMatrix));
    printf("Relative 2-norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_RELATIVE_L2,MaskMatrix));
    printf("Relative infininte norm distance is : %.1f\n",cvNorm(Matrix1,Matrix2,CV_RELATIVE_C,MaskMatrix));

    system("pause");
}

執行結果:

由上面看,有用到的Matrix1向量為<5,3>,Matrix2的向量為<6,-1>,因此Mask為0的向量數據都被忽略了,而Mask為1的被保留,計算的方式如下所示.

cvNorm()
計算向量的長度或距離,分為1-Norm(CV_L1),2-Norm(CV_L2),∞-Norm(CV_C)等三種方法,也可以計算相對距離,分別為CV_RELATIVE_C,CV_RELATIVE_L1,CV_RELATIVE_L2,在做資訊檢索比對上有很大的用處.第一個跟第二個引數為輸入向量的CvMat資料結構,如果要計算長度的話第二個引數要為空,第三個引數為cvNorm()的參數,第四個引數為單通道CvMat資料結構的CV_8UC1參數(非浮點數型別)的遮罩矩陣,cvNorm()的輸出為一個double型別的數據,代表著向量的長度或距離.
cvNorm(輸入CvMat資料結構,輸入CvMat資料結構或為空,參數或代號,輸入CvMat遮罩)

OpenCV線性代數-線性映射,矩陣變換

2008-08-19T12:00:00.021+08:00

線性映射(Linear transformation),是由向量為基礎的映射函式,在OpenCV內標準的向量表達方式為通道的向量,也就是今天為一組<r,g,b>三通道所組成的圖形,則圖形內的<r,g,b>個別值則代表一個向量,cvTransform()則是為此而做的轉換,下面則是簡單的一個通道向量的運算

簡單線性映射實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

float Array1[]={4,1,0,2,1,1,3,2,1,1,4,1.5};
float Array2[]={2,1,2};
int main()
{
    CvMat *A=cvCreateMat(4,3,CV_32FC1);
    CvMat *x=cvCreateMat(1,1,CV_32FC3);
    CvMat *T_x=cvCreateMat(1,1,CV_32FC4);

    cvSetData(A,Array1,A->step);
    cvSetData(x,Array2,x->step);
    cvTransform(x,T_x,A);

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

執行結果:

A矩陣則是叫轉移矩陣(Transition Matrix),可以專門在做基底座標變換的運用,這邊要特別注意的是,cvTransform()它拿通道來代替向量,而不是拿矩陣來代替向量,事實上,矩陣也可以作為向量的表達,但不適用於cvTransform().這邊cvTransform()輸入,x矩陣為通道三,而A則一定要輸入通道一的矩陣,輸出的T_x矩陣為通道四,而對於一個向量而言,此種運算方式叫做線性映射,如果對於一個矩陣圖形的話,則叫做矩陣變換(Matrix Transformation).

矩陣變換的部分,就拿RGB轉YCbCr的公式來實作,跟前面的線性映射不同的是,它直接拿通道三的矩陣來做色彩空間的轉換.

RGB轉YCbCr實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

float Array1[]={0.2989,0.5867,0.1140,-0.1687,-0.3312,0.5,0.5,-0.4183,-0.0816};
float Array2[]={0,128,128};

int main()
{
    CvMat *TransitionMatrix=cvCreateMat(3,3,CV_32FC1);
    CvMat *ShiftVector=cvCreateMat(3,1,CV_32FC1);

    IplImage *Image1=cvLoadImage("MineCar.jpg",1);
    IplImage *GrayImage1=cvLoadImage("MineCar.jpg",0);
    IplImage *YCbCrImage1=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,3);
    IplImage *YImage=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);
    IplImage *CbImage=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);
    IplImage *CrImage=cvCreateImage(cvGetSize(Image1),IPL_DEPTH_8U,1);

    cvSetData(TransitionMatrix,Array1,TransitionMatrix->step);
    cvSetData(ShiftVector,Array2,ShiftVector->step);

    cvTransform(Image1,YCbCrImage1,TransitionMatrix,ShiftVector);

    cvSplit(YCbCrImage1,YImage,CbImage,CrImage,0);

    cvNamedWindow("YImage",1);
    cvShowImage("YImage",YImage);
    cvNamedWindow("GrayImage",1);
    cvShowImage("GrayImage",GrayImage1);

    cvWaitKey(0);
}

原始圖片:

執行結果:

由上面線性映射的基礎來對一個通道三個圖形做矩陣變換,這邊將RGB的圖片藉由公式轉成YCbCr的格式,這個部分就叫做色彩空間的轉換,Y代表的是亮度,也就是灰階值,Cb則是藍色的色差,Cr則為紅色的色差,因此,這邊直接拿cvLoadImage()讀取灰階圖形來跟Y值圖形來做比較.中間有用到通道分割的函式cvSplit()則是要參考前面資料結構"通道的分割,合併與混合"的部份,cvTransform()的第四個引數則是對轉移矩陣計算後的數值結果做平移的運算.

cvTransform()
將圖形或矩陣做矩陣變換,矩陣變換的向量維度是以通道為基礎,也就是轉移矩陣會對通道向量做相乘而完成色彩空間的轉換,第四個引數則是對計算結果的向量做位移的調整,而轉移矩陣的m*n大小可以決定輸出通道的維度,這邊m代表輸入通道數,n代表輸出通道數.
cvTransform(輸入CvMat或IplImage資料結構,輸出CvMat或IplImage資料結構,CvMat轉移矩陣資料結構,CvMat平移向量資料結構)

OpenCV線性代數-cvSVBkSb奇異值倒回替代法

2008-08-16T13:20:00.017+08:00

cvSVBkSb()為解線性系統解專用的,實際上它跟前面的cvSolve()的線性系統求解的奇異值法是一樣的東西,cvSolve()內部的CV_SVD參數它的計算方式實際上就是呼叫cvSVBkSb()這個函式,所以,用起來的方式差不多,在這邊簡單的做cvSVBkSb()的使用方法.

線性系統求解-奇異值倒回替代法
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

float Array1[]={-1,2,-1,3,1,-4,5,-5,-2,6,-5,7,-1,-4,11,-2};
float Array2[]={6,-7,7,7};
int main()
{
    CvMat *A=cvCreateMat(4,4,CV_32FC1);
    CvMat *b=cvCreateMat(4,1,CV_32FC1);
    CvMat *W=cvCreateMat(4,4,CV_32FC1);
    CvMat *U=cvCreateMat(4,4,CV_32FC1);
    CvMat *V=cvCreateMat(4,4,CV_32FC1);
    CvMat *x=cvCreateMat(4,1,CV_32FC1);

    cvSetData(A,Array1,A->step);
    cvSetData(b,Array2,b->step);
    cvSVD(A,W,U,V);
    cvSVBkSb(W,U,V,b,x,0);

    PrintMatrix(x,x->rows,x->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");
    }
}

在這邊,要用cvSVBkSb()之前就要先用到前面的cvSVD()了,在將cvSVD()奇異值分解的數值餵進cvSVBkSb()裡面解線性系統解,這邊奇異值倒回替代法用的是虛反矩陣的觀念,虛反矩陣則是在前面cvInvert()裡面有提到,而它的第四個引數為cvSVBkSb()的參數及代號,與cvSVD()相同,可輸入的參數為

#define CV_SVD_U_T 2
#define CV_SVD_V_T 4

分別代表U矩陣跟V矩陣輸入是否為轉置,如果不需要任何轉置輸入則將參數代號預設為0,而此題的線性系統,驗證及計算的方式如下

cvSVBkSb()
為線性系統求解的奇異值替代倒回法,對於一個無解的線性系統,使用奇異值替代倒回法可以求得近似解,但是使用cvSVBkSb()這個函式時必須要先做cvSVD()奇異值分解的計算,在將分解後的矩陣丟給cvSVBkSb()求線性系統的解,使用到的是虛反矩陣的概念,第一個引數為CvMat結構的W矩陣,第二個引數為CvMat結構的U矩陣,第三個引數為CvMat結構的V矩陣,前面的矩陣結構都必須經由cvSVD()奇異值分解計算而得,第四個引數為Ax的線性組合數據,第五個引數為x的解集合,第六個引數則是cvSVBkSb()的參數輸入,只有CV_SVD_U_T及CV_SVD_V_T或者是0(無參數)這三種.
cvSVBkSb(CvMat資料結構W矩陣輸入,CvMat資料結構U矩陣輸入,CvMat資料結構V矩陣輸入,CvMat資料結構b矩陣輸入,CvMat資料結構線性系統解集合x矩陣,目標參數或代號)

OpenCV線性代數-cvSVD奇異值分解(2)

2008-08-14T07:28:00.010+08:00

這邊就要簡單介紹cvSVD()的所有參數及用法,cvSVD()在"cxcore.h"裡面的參數分別被定義為

#define CV_SVD_MODIFY_A 1
#define CV_SVD_U_T 2
#define CV_SVD_V_T 4

第一個CV_SVD_MODIFY_A是將A矩陣的空間做填充,用空間換取時間的方式做加速的運算,而後面的CV_SVD_U_T,及CV_SVD_V_T則是單純的將U跟V做轉置.接著下面的程式拿非奇異矩陣做奇異值分解並且加上cvSVD()的參數.

非奇異矩陣奇異值分解
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

double Array1[]={1,3,3,-3,-5,-3,3,3,1};
int main()
{
    CvMat *Matrix1=cvCreateMat(3,3,CV_64FC1);
    CvMat *W=cvCreateMat(3,3,CV_64FC1);
    CvMat *V=cvCreateMat(3,3,CV_64FC1);
    CvMat *U=cvCreateMat(3,3,CV_64FC1);
    CvMat *V_T=cvCreateMat(3,3,CV_64FC1);
    CvMat *ResultMatrix=cvCreateMat(3,3,CV_64FC1);

    cvSetData(Matrix1,Array1,Matrix1->step);
    cvSVD(Matrix1,W,U,V,CV_SVD_MODIFY_A);

    printf("\nW\n");
    PrintMatrix(W,W->rows,W->cols);

    printf("\nU\n");
    PrintMatrix(U,U->rows,U->cols);
    printf("\nV\n");
    PrintMatrix(V,V->rows,V->cols);

    printf("\nValid\n");
    cvmMul(U,W,ResultMatrix);
    cvTranspose(V,V_T);
    cvmMul(ResultMatrix,V_T,ResultMatrix);
    PrintMatrix(ResultMatrix,ResultMatrix->rows,ResultMatrix->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");
    }
}

執行結果:

上面使用了CV_SVD_MODIFY_A,由上面的非奇異矩陣得知,無法求得他的Eign Value及Eigen Vector,但是可以求得它的奇異值以及它的分解矩陣,並且藉由矩陣相乘的方式計算回原矩陣.計算及驗證的方式如下

而使用了CV_SVD_V_T這個參數則可以省略掉對V做cvTranspose()轉置的步驟.

cvSVD()使用CV_SVD_V_T參數
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

double Array1[]={1,3,3,-3,-5,-3,3,3,1};
int main()
{
    CvMat *Matrix1=cvCreateMat(3,3,CV_64FC1);
    CvMat *W=cvCreateMat(3,3,CV_64FC1);
    CvMat *U=cvCreateMat(3,3,CV_64FC1);
    CvMat *V_T=cvCreateMat(3,3,CV_64FC1);
    CvMat *ResultMatrix=cvCreateMat(3,3,CV_64FC1);

    cvSetData(Matrix1,Array1,Matrix1->step);
    cvSVD(Matrix1,W,U,V_T,CV_SVD_V_T);

    printf("\nW\n");
    PrintMatrix(W,W->rows,W->cols);

    printf("\nU\n");
    PrintMatrix(U,U->rows,U->cols);
    printf("\nV_T\n");
    PrintMatrix(V_T,V_T->rows,V_T->cols);

    printf("\nValid\n");
    cvmMul(U,W,ResultMatrix);
    cvmMul(ResultMatrix,V_T,ResultMatrix);
    PrintMatrix(ResultMatrix,ResultMatrix->rows,ResultMatrix->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");
    }
}

執行結果:

上面的程式碼跟前面的差別,就是直接在cvSVD()內部將V矩陣轉置,可以直接將奇異值分解做驗證,而CV_SVD_U_T亦是直接將U矩陣在內部轉置.

cvSVD()
將任何的矩陣做奇異值分解,假如要分解的矩陣是m*n,則W矩陣CvMat資料結構輸入格式必須是m*n,U矩陣CvMat資料結構輸入格式必須是m*m,V矩陣CvMat資料結構輸入格式必須是n*n,cvSVD()第一個引數為輸入目標CvMat矩陣資料結構,第二個引數為輸出W矩陣CvMat資料結構,第三個引數為U矩陣CvMat資料結構,第四個引數為V矩陣CvMat資料結構,第五個引數為cvSVD()的參數,分別為CV_SVD_MODIFY_A,CV_SVD_U_T,CV_SVD_V_T.
cvSVD(輸入CvMat資料結構,輸出W矩陣CvMat結構,輸出U矩陣CvMat結構,輸出V矩陣CvMat結構,目標參數或代號)

OpenCV線性代數-cvSVD奇異值分解(1)

2008-08-11T16:08:00.041+08:00

奇異值分解(Singular Value Decomposition,SVD),可以說是對角化應用的特例,它具有與Eigen Value,Eigen Vector相同的特性,分解後向量的矩陣乘積可以還原為原矩陣,而他可以對奇異矩陣(Singular Matrix)做分解,亦可以對非奇異矩陣做分解.奇異矩陣的定義就是不能計算為反矩陣的矩陣,包括方陣,以及非方陣的長方形矩陣.而奇異值分解的計算它的數學定義為

它將奇異矩陣轉置後相乘,藉由這個方法計算出對稱矩陣,在對對稱矩陣做SVD分解,由上面可知,原矩陣大小為m*n,而W(Sigma)矩陣的舉陣大小必須為m*n,U的矩陣大小必須為m*m,而V的矩陣大小為n*n,它可以藉由矩陣矩陣相乘的方式算回原矩陣.

下面這個是拿2*3的奇異矩陣做奇異值分解

奇異矩陣SVD分解
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

double Array1[]={3,2,2,2,3,-2};
int main()
{
    CvMat *Matrix1=cvCreateMat(2,3,CV_64FC1);
    CvMat *W=cvCreateMat(2,3,CV_64FC1);
    CvMat *V=cvCreateMat(3,3,CV_64FC1);
    CvMat *U=cvCreateMat(2,2,CV_64FC1);
    CvMat *V_T=cvCreateMat(3,3,CV_64FC1);
    CvMat *ResultMatrix=cvCreateMat(2,3,CV_64FC1);

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

    cvSVD(Matrix1,W,U,V);

    printf("\nW\n");
    PrintMatrix(W,W->rows,W->cols);
    printf("\nU\n");
    PrintMatrix(U,U->rows,U->cols);
    printf("\nV\n");
    PrintMatrix(V,V->rows,V->cols);

    printf("\nValid\n");
    cvmMul(U,W,ResultMatrix);
    cvTranspose(V,V_T);
    cvmMul(ResultMatrix,V_T,ResultMatrix);
    PrintMatrix(ResultMatrix,ResultMatrix->rows,ResultMatrix->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");
    }
}

執行結果:

上面的程式,利用了前面奇異值分解公式說明的方式,cvSVD()的內部運算,先將矩陣做轉置在相乘,求得對稱矩陣,在將該矩陣做奇異值分解,而它的W,U,V矩陣的大小皆不相同,W將會是原矩陣大小2*3,V則是2*2,U則會是3*3,因此利用公式倒回運算,可以計算出原矩陣,計算以及驗證的公式如下.

而對稱矩陣的奇異值分解,則是不需要再做轉置後相乘,直接做奇異值分解.

對稱矩陣奇異值分解
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

double Array1[]={2,-2,4,-2,2,-4,4,-4,8};
int main()
{
    CvMat *Matrix1=cvCreateMat(3,3,CV_64FC1);
    CvMat *W=cvCreateMat(3,3,CV_64FC1);
    CvMat *V=cvCreateMat(3,3,CV_64FC1);
    CvMat *U=cvCreateMat(3,3,CV_64FC1);
    CvMat *V_T=cvCreateMat(3,3,CV_64FC1);
    CvMat *ResultMatrix=cvCreateMat(3,3,CV_64FC1);

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

    cvSVD(Matrix1,W,U,V);

    printf("\nW\n");
    PrintMatrix(W,W->rows,W->cols);
    printf("\nU\n");
    PrintMatrix(U,U->rows,U->cols);
    printf("\nV\n");
    PrintMatrix(V,V->rows,V->cols);

    printf("\nValid\n");
    cvmMul(U,W,ResultMatrix);
    cvTranspose(V,V_T);
    cvmMul(ResultMatrix,V_T,ResultMatrix);
    PrintMatrix(ResultMatrix,ResultMatrix->rows,ResultMatrix->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");
    }
}

執行結果:

這是一個直接用對稱矩陣做奇異值分解的程式,奇異值分解很特別的地方,它對對稱矩陣做奇異值分解求得的結果為特徵值(Eigen Value)及特徵向量(Eigen Vector),但是相較於cvEigenVV()這兩個使用的是不同的演算法,而一般的矩陣則無法用cvSVD()求特徵值及特徵向量,上面的程式用cvSVD()分別求出它的W,U,V,在代回原公式做驗證,驗證的方法如下

驗證對稱矩陣奇異值分解

OpenCV線性代數-cvEigenVV實作

2008-08-03T05:32:00.026+08:00

cvEigenVV()為計算CvMat方陣的特徵值(Eigen Value)跟特徵向量(Eigen Vector)的函式,但是,OpenCV的Eigen Value跟Eigen Vector計算並不是一般性的用法,無法處理正常的特徵值跟特徵向量的計算,cvEigenVV()是使用到Jacobi Eigenvalue Algorithm的方法,輸入必須要對稱矩陣,輸入的對稱矩陣做Jacobi Transformation的轉換,這部份的資料就要參考數值分析(Numerical Recipes)的相關書籍了,而Jacobi Eigenvalue Algorithm有將Eigen Value重新Sort過,因此特徵值的排序會是由大到小排列,而不適用於很多矩陣對角化的解題應用.

cvEigenVV()實作
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

double Array1[]={2,3,0,3,5,-1,0,-1,2};

int main()
{
    CvMat *Matrix1=cvCreateMat(3,3,CV_64FC1);
    CvMat *EigenValue_Row=cvCreateMat(3,1,CV_64FC1);
    CvMat *EigenValue=cvCreateMat(3,3,CV_64FC1);
    CvMat *EigenVector=cvCreateMat(3,3,CV_64FC1);
    CvMat *EigenVector_Invert=cvCreateMat(3,3,CV_64FC1);
    CvMat *ResultMatrix=cvCreateMat(3,3,CV_64FC1);
    cvSetData(Matrix1,Array1,Matrix1->step);
    cvSetZero(EigenValue);

    cvEigenVV(Matrix1,EigenVector,EigenValue_Row,DBL_EPSILON);

    printf("\nThe EigenValue_Row is:\n");
    PrintMatrix(EigenValue_Row,EigenValue_Row->rows,EigenValue_Row->cols);

    printf("\nThe EigenVector is:\n");
    PrintMatrix(EigenVector,EigenVector->rows,EigenVector->cols);

    printf("\nThe EigenValue is:\n");
    cvSet2D(EigenValue,0,0,cvGet2D(EigenValue_Row,0,0));
    cvSet2D(EigenValue,1,1,cvGet2D(EigenValue_Row,1,0));
    cvSet2D(EigenValue,2,2,cvGet2D(EigenValue_Row,2,0));
    PrintMatrix(EigenValue,EigenValue->rows,EigenValue->cols);

    cvTranspose(EigenVector,EigenVector);
    cvmMul(EigenVector,EigenValue,ResultMatrix);
    cvInvert(EigenVector,EigenVector_Invert,CV_LU);
    cvmMul(ResultMatrix,EigenVector_Invert,ResultMatrix);

    printf("\nTo validate Matrix\n");
    PrintMatrix(ResultMatrix,ResultMatrix->rows,ResultMatrix->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");
    }
}
執行結果:

由上面可以知道Eigen Value的排列方式為

而Eigen Vector的排列方式為

因此,要驗證是否可行,利用

的方法將它還原為對稱矩陣,則必須要將輸出的列矩陣補0讓它成為對角矩陣,並且將Eigen Vector矩陣做轉置,帶入驗證公式,輸出結果為原對稱矩陣,因此可以證明此函式及輸入值無誤

cvEigenVV()第一的引數為特徵向量(Eigen Vector)第二個引數為特徵值(Eigen Value)第三個引數為精確度,DBL_EPSILON為最小double型別精確度,而一般精確度的使用定義為

#define FLT_EPSILON 1.19209290E-07F
#define DBL_EPSILON 2.2204460492503131E-16

這可以來做浮點數的精確度比較運算使用.
而他的使用方式就如同下所述

bool IsEqual(double x,double y)
{
    if(x-y<DBL_EPSILON)
    {
        return true;
    }
    else
    {
        return false;
    }
}

double型別則對照DBL_EPSILON參數,float型別則對照FLT_EPSILON,而cvEigenVV()則可以自行定義它的精確度比較運算的大小.

cvEigenVV()
利用Jacobi Eigenvalue Algorithm法計算Eigen Value,Eigen Vector,輸入只能為對稱矩陣,輸出則是Eigen Value的列舉陣及Eigen Vector的列舉陣,第一個引數為要計算的CvMat資料結構對稱矩陣,第二個引數為CvMat資料結構的Eigen Value列舉陣,第三個引數是CvMat資料結構的Eigen Vector列舉陣,第四個引數為精確度比較運算的大小.
cvEigenVV(輸入CvMat資料結構對稱矩陣,輸出CvMat資料結構的EigenValue列舉陣,輸出CvMat資料結構的EigenVector列舉陣,浮點型別精準度誤差比較數值)

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

2008-08-01T12:15:00.009+08:00

在cvSolve()裡面,輸入的係數矩陣規定要方陣,那假如所求的線性系統不是方陣怎麼辦呢？那當然就是要用補0了方式填充成方陣,使用的方法請看以下程式碼.

解線性方程式4
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
void PrintMatrix(CvMat *Matrix,int Rows,int Cols);

float Array1[]={1,1,0,2,1,0,3,2,0};
float Array2[]={1,3,4};
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_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");
    }
}

執行結果:

線性系統判斷

上面的程式是尾端補0,因為這個線性系統只存在兩個未知數,三個方程式,要給它弄成方陣,則只好假設有第三個未知數存在,但係數為0,而只要是需要補0的線性系統運算,全部都需要用CV_SVD來解,LU是無法解的.再來,下面是整欄補零的範例.

解線性方程式5
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
void PrintMatrix(CvMat *Matrix,int Rows,int Cols);

float Array1[]={1,1,1,1,2,0,0,0,0};
float Array2[]={4,6,0};
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_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");
    }
}

執行結果:

線性系統判斷

這個線性系統是四個未知數,三個方程式,很顯然的是無限解,用SVD的方法的話就會任意給一個可行解,而至於cvSolve()的另外一個參數CV_SVD_SYM,則是線性系統對對稱係數矩陣的加速算法,也是跟前面相關的程式碼同樣的用法,這邊不再累述.

cvSolve()
為解線性方程式的函式,可以對唯一解,無解,無限解,而無限解的問題會給予一個答案,無解的問題則是用CV_SVD參數則可得到近似解,而無限解的問題則用cvSolve()參數可以得到可行解.第一個引數為CvMat係數矩陣,第二個引數為線性方程式的線性組合解,第三個引數為代數的解集合,第四個為解線性方程的參數輸入,分別為CV_LU,CV_SVD,CV_SVD_SYM,與求反矩陣的函式相同參數.
cvSolve(CvMat係數矩陣結構,CvMat結構線性組合解,CvMat結構代數解集合,輸入計算參數或代號)

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

2008-07-28T15:50:00.032+08:00

線性系統求解,也就是求線性聯立方程式(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法則可以把數據帶進去,會得到近似的線性組合解.

接著為線性系統判斷

OpenCV線性代數-行列式,Cramer's rule

2008-07-28T06:56:00.009+08:00

再來就是簡單的介紹行列式的函式cvDet(),行列式的輸入CvMat資料結構必須要為方陣,對矩陣的內部做一些特殊運算,可以用來求反矩陣,判斷矩陣線性解,特徵值及特徵向量等應用,其他行列式的性質則不再這裡詳述.

簡單行列式運算
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

float Array[]={1,3,1,3,4,2,2,1,4,2,2,3,3,1,4,1};
int main()
{
    double Determinant;
    CvMat *Matrix1=cvCreateMat(4,4,CV_32FC1);
    cvSetData(Matrix1,Array,Matrix1->step);

    printf("\nThe Matrix is:\n");
    PrintMatrix(Matrix1,Matrix1->rows,Matrix1->cols);

    printf("\nThe Determinant is:\n");
    Determinant=cvDet(Matrix1);
    printf("%.f\n",Determinant);

    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");
    }
}

執行結果:

這是一個4*4的矩陣,行列式的輸出為一個數據,當行列視為0則為不可逆矩陣(不具反矩陣),行列式也可以解線性系統的解,而用行列式解線性系統的方法叫做Cramer's rule,以下則由此方程式求x,y,z的解

行列式Cramer's rule
#include <cv.h>
#include <highgui.h>
#include <stdio.h>

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

float Array[]={2,3,-1,3,5,2,1,-2,-3};
float Array1[]={1,3,-1,8,5,2,-1,-2,-3};
float Array2[]={2,1,-1,3,8,2,1,-1,-3};
float Array3[]={2,3,1,3,5,8,1,-2,-1};

int main()
{
    double x,y,z;
    double Determinant;
    double Determinant1;
    double Determinant2;
    double Determinant3;

    CvMat *Matrix=cvCreateMat(3,3,CV_32FC1);
    CvMat *Matrix1=cvCreateMat(3,3,CV_32FC1);
    CvMat *Matrix2=cvCreateMat(3,3,CV_32FC1);
    CvMat *Matrix3=cvCreateMat(3,3,CV_32FC1);

    cvSetData(Matrix,Array,Matrix->step);
    cvSetData(Matrix1,Array1,Matrix1->step);
    cvSetData(Matrix2,Array2,Matrix2->step);
    cvSetData(Matrix3,Array3,Matrix3->step);

    Determinant=cvDet(Matrix);
    Determinant1=cvDet(Matrix1);
    Determinant2=cvDet(Matrix2);
    Determinant3=cvDet(Matrix3);

    x=Determinant1/Determinant;
    y=Determinant2/Determinant;
    z=Determinant3/Determinant;

    printf("The equation solution is :\n");
    printf("x = %.f\n",x);
    printf("y = %.f\n",y);
    printf("z = %.f\n",z);

    system("pause");
}

執行結果:

在這裡Matrix為解集合的線性組合,而Matrix1為將線性組合的第一列與解集合做交換,Matrix2為第二列與解集合做交換,Matrix3為第三列與解集合做交換,在將Matrix,Matrix1,Matrix2,Matrix3解行列式,而x,y,z的解就分別是x=det(Matrix1)/det(Matrix),y=det(Matrix2)/det(Matrix),z=det(Matrix3)/det(Matrix).

Cramer's rule在使用上執行效率並不會很好,一般都是以高斯消去法來解聯立的線性系統.

cvDet()
計算CvMat資料結構的行列式值,回傳結果為double型別資料型態,但是在線性代數的定義上行列式的值為純量(Scale)
double cvDet(CvMat資料結構)