Harris角點檢測學習

1.角點的定義與性質

      角點是一種局部特徵，具備旋轉不變性和不隨光照條件變化而變化的特色，通常將圖像中曲率足夠高或者曲率變化明顯的點做爲角點。檢測獲得的角點特徵一般用於圖像匹配、目標跟蹤、運動估計等方面。

 

2.Harris角點

1）定性描述

      該算法中，將圖像分爲平坦區域、邊緣、角點三部分。平坦區域中像素灰度在各個方向上變換都很小，邊緣上的像素灰度在某個方向變化很大，可是在另外一些方向變化很小；位於角點上的像素灰度則在各個方向上的變化都比較大。這是經過人眼觀察獲得的直觀感覺。

2）量化的數學表達

      不一樣區域像素的變化與類似性是相反的關係：變化越大，類似性越小，反之亦然。因此咱們可使用類似性度量其變化程度。圖像I(x,y)在中心(x,y)處區域W(x,y)與相對中心(x,y)移動了△x/△y以後的區域相關性能夠經過自相關函數給出

      w(u,v)是加權函數，可選取常數或者高斯加權函數（通常選取後者）

      看到delta標誌就會想到泰勒展開式簡化，簡化流程以下：

 

      因此以前的類似性度量能夠近似表達爲一個二次項函數

 

      其本質是一個橢圓函數，橢圓的主方向就是圖像變化率最大的主要方向。而橢圓的大小又是由M(x,y)的特徵值決定的。

 

      套用這套理論到平坦區域、邊緣、角點上，則表示爲：

（1）邊緣；橢圓狹長，一個特徵值大、另外一個特徵值小

（2）平坦區域；兩個特徵值都很小

（3）角點；兩個特徵值都很大

 

      因此求出矩陣M(x,y)的兩個特徵值能夠用於檢測角點，可是計算量太大。Harris角點檢測中避免直接求取特徵值，而是利用M(x,y)的性質來計算兩個特徵值的相對大小

 

      只有兩個特徵值都比較大的狀況下R取值才比較大

 

3.算法流程

      根據前面的描述，只要計算出圖像中各個像素對應的R，並根據一個設定的閾值二值化處理便可獲得角點。爲了濾除局部區域的非極大值，可使用非極大值抑制法作進一步的處理。算法整理流程以下：

 

4.opencv實現

       我的使用的開發環境是是opencv_2.4.13+vs2012，現將代碼貼出

 1 #include <iostream>
 2 #include <core/core.hpp>
 3 #include <highgui/highgui.hpp>
 4 #include <features2d/features2d.hpp>
 5 #include <imgproc/imgproc.hpp>
 6 
 7 using namespace std;
 8 using namespace cv;
 9 
10 int main(int argc, char* argv[])
11 {
12     /* 1.以bmp格式讀取圖片，注意該圖片此時應該放置於.sln所在目錄 */
13     Mat image = imread("../church01.jpg", 0);
14     if(!image.data)
15         return 0;
16 
17     namedWindow("originalImage");
18     imshow("originalImage", image);
19     
20     /* 2.計算圖像沿x/y方向一階導數 */
21     Mat xKernel = (Mat_<double>(1,3) << -1, 0, 1);
22     Mat yKernel = xKernel.t();
23 
24     Mat Ix, Iy;
25     filter2D(image, Ix, CV_64F, xKernel);
26     filter2D(image, Iy, CV_64F, yKernel);
27     
28     /* 3.計算矩陣M中的各項 */
29     Mat Ix2, Iy2, Ixy;
30     Ix2 = Ix.mul(Ix);
31     Iy2 = Iy.mul(Iy);
32     Ixy = Ix.mul(Iy);
33 
34     /* 注意這裏使用的是一階高斯濾波而非二階，還沒有
35      * 理解其緣由
36      */
37     Mat gaussKernel = getGaussianKernel(7,1);
38     filter2D(Ix2, Ix2, CV_64F, gaussKernel);
39     filter2D(Iy2, Iy2, CV_64F, gaussKernel);
40     filter2D(Ixy, Ixy, CV_64F, gaussKernel);
41 
42     /* 4.根據公式計算R矩陣 */
43     Mat cornerStrength(image.size(), CV_64F);
44     float alpha = 0.1f;
45     int rows, cols;
46     rows = image.rows;
47     cols = image.cols;
48     for(int i = 0; i < rows; i++)
49     {
50         for(int j = 0; j < cols; j++)
51         {
52             double det_m = Ix2.at<double>(i,j) * Iy2.at<double>(i,j) - Ixy.at<double>(i,j) * Ixy.at<double>(i,j);
53             double trace_m = Ix2.at<double>(i,j) + Iy2.at<double>(i,j);
54             cornerStrength.at<double>(i,j) = det_m - alpha*trace_m*trace_m;
55         }
56     }
57 
58     /* 5.閾值化並進行非極大值抑制 */
59     double maxStrength, minStrength;
60     minMaxLoc(cornerStrength, &minStrength, &maxStrength);
61     Mat dilated;
62     Mat localMax;
63     dilate(cornerStrength, dilated, Mat());
64     compare(cornerStrength, dilated, localMax, CMP_EQ);
65 
66     double threshold = 0.01 * maxStrength;
67     Mat cornerMap;
68     cornerMap = cornerStrength > threshold;
69     bitwise_and(cornerMap, localMax, cornerMap);
70     namedWindow("cornerMap");
71     imshow("cornerMap", cornerMap);
72 
73     /* 6.在原圖上繪製角點 */
74     image = imread("../church01.jpg", 0);
75     for(int x = 0;x < cornerMap.cols; x++)
76     {
77         for(int y = 0; y < cornerMap.rows; y++)
78         {
79             if(cornerMap.at<uchar>(y,x))
80             {
81                 circle(image, Point(x,y), 3, Scalar(255), 1);
82 
83             }
84         }
85     }
86 
87     namedWindow("Harris Corner");
88     imshow("Harris Corner", image);
89     waitKey();
90     return 0;
91 }

View Code

      代碼執行效果以下，originalImage是讀入的原圖，cornerMap是計算獲得的角點（白色顯示），Harris Corner是在原圖上繪製角點的效果。

5.參考資料

 [1]《圖像局部不變性特徵與描述》

 [2]http://www.cnblogs.com/polly333/p/5416172.htmlHarris角點算法