2. 程式人生 > >Harris角點檢測演算法詳解

Harris角點檢測演算法詳解

阿新 • • 發佈：2018-12-18

Harris角點演算法

特徵點檢測廣泛應用到目標匹配、目標跟蹤、三維重建等應用中，在進行目標建模時會對影象進行目標特徵的提取，常用的有顏色、角點、特徵點、輪廓、紋理等特徵。現在開始講解常用的特徵點檢測，其中Harris角點檢測是特徵點檢測的基礎，提出了應用鄰近畫素點灰度差值概念，從而進行判斷是否為角點、邊緣、平滑區域。Harris角點檢測原理是利用移動的視窗在影象中計算灰度變化值，其中關鍵流程包括轉化為灰度影象、計算差分影象、高斯平滑、計算區域性極值、確認角點。

一、基礎知識

影象的變化型別：

在特徵點檢測中經常提出尺度不變、旋轉不變、抗噪聲影響等，這些是判斷特徵點是否穩定的指標。

效能較好的角點：

  1. 檢測出影象中“真實”的角點
  2. 準確的定位效能
  3. 很高的重複檢測率
  4. 噪聲的魯棒性
  5. 較高的計算效率

角點的型別：

基於影象灰度的方法通過計算點的曲率及梯度來檢測角點，避免了第一類方法存在的缺陷，此類方法主要有Moravec運算元、Forstner運算元、Harris運算元、SUSAN運算元等。

二、Harris角點原理

1、演算法思想

角點原理來源於人對角點的感性判斷，即影象在各個方向灰度有明顯變化。演算法的核心是利用區域性視窗在影象上進行移動判斷灰度發生較大的變化，所以此視窗用於計算影象的灰度變化為：[-1,0,1;-1,0,1;-1,0,1][-1,-1,-1;0,0,0;1,1,1]。人各個方向上移動這個特徵的小視窗，如圖3中視窗內區域的灰度發生了較大的變化，那麼就認為在視窗內遇到了角點。如圖1中，視窗內影象的灰度沒有發生變化，那麼視窗內就不存在角點；如果視窗在某一個方向移動時，視窗內影象的灰度發生了較大的變化，而在另一些方向上沒有發生變化，那麼，視窗內的影象可能就是一條直線的線段。

2、數學模型

根據演算法思想，構建數學模型，計算移動視窗的的灰度差值。

為了減小計算量，利用泰勒級數進行簡化公式：

上圖中W函式表示視窗函式，M矩陣為偏導數矩陣。對於矩陣可以進行對稱矩陣的變化，假設利用兩個特徵值進行替代，其幾何含義類似下圖中的表達。在幾何

模型中通過判斷兩個特徵值的大小，來判定畫素的屬性。

M為梯度的協方差矩陣 ，在實際應用中為了能夠應用更好的程式設計，定義了角點響應函式R，通過判定R大小來判斷畫素是否為角點。

R取決於M的特徵值，對於角點|R|很大，平坦的區域|R|很小，邊緣的R為負值。

 

3、演算法流程

1.利用水平，豎直差分運算元對影象的每個畫素進行濾波以求得Ix,Iy,進而求得M中的四個元素的值。

演算法原始碼：

程式碼中如果array為-1,0,1,-1,0,1,-1,0,1}則是求解X方向的，如果為{-1,-1,-1,0,0,0,1,1,1}為Y方向的，則Ix和Iy求解結束

複製程式碼 
//C程式碼
//這裡的size/2是為了不把影象邊界算進去。
//array也就是3*3的視窗函式為：double dx[9]={-1,0,1,-1,0,1,-1,0,1};或者 //定義垂直方向差分運算元並求Iy
    double dy[9]={-1,-1,-1,0,0,0,1,1,1};

CvMat *mbys(CvMat *mat,int xwidth,int ywidth,double *array,int size1,int size2)//size 
{
    int i,j;
    int i1,j1;
    int px,py;
    int m;
    CvMat *mat1;
    mat1=cvCloneMat(mat);
//為了去除邊界，從框體一半開始遍歷
    for(i=size1/2;i<ywidth-size1/2;i++)
        for(j=size2/2;j<xwidth-size2/2;j++)                         
        {
            m=0;
            for(i1=0;i1<size1;i1++)
                for(j1=0;j1<size2;j1++)
                {
                    px=i-size1/2+i1;
                    py=j-size2/2+j1;
                    //CV_MAT_ELEM訪問矩陣元素
                  //每個元素和窗體函式遍歷相加
                    m+=CV_MAT_ELEM(*mat,double,px,py)*array[i1*size1+j1];            
                }
               //賦值
                CV_MAT_ELEM(*mat1,double,i,j)=m;
        }
        return mat1;
}
複製程式碼

求解IX2相對比較簡單，畫素相乘即可。下面為基於opencv的C++版本，很是簡單

複製程式碼 
//構建模板
Mat xKernel = (Mat_<double>(1,3) << -1, 0, 1);
Mat yKernel = xKernel.t();
//計算IX和IY
Mat Ix,Iy;
//可參考filter2d的定義
filter2D(gray, Ix, CV_64F, xKernel);
filter2D(gray, Iy, CV_64F, yKernel);
//計算其平方
Mat Ix2,Iy2,Ixy;
Ix2 = Ix.mul(Ix);
Iy2 = Iy.mul(Iy);
Ixy = Ix.mul(Iy);
複製程式碼

2.對M的四個元素進行高斯平滑濾波，為的是消除一些不必要的孤立點和凸起，得到新的矩陣M。

複製程式碼 
//本例中使用5×5的高斯模板
    //計算模板引數
    //int gausswidth=5;
    //double sigma=0.8;
    double *g=new double[gausswidth*gausswidth];
    for(i=0;i<gausswidth;i++)//定義模板
        for(j=0;j<gausswidth;j++)
            g[i*gausswidth+j]=exp(-((i-int(gausswidth/2))*(i-int(gausswidth/2))+(j-int(gausswidth/2))*(j-int(gausswidth/2)))/(2*sigma));

    //歸一化：使模板引數之和為1（其實此步可以省略）
    double total=0;
    for(i=0;i<gausswidth*gausswidth;i++)
        total+=g[i];
    for(i=0;i<gausswidth;i++)
        for(j=0;j<gausswidth;j++)
            g[i*gausswidth+j]/=total;

    //進行高斯平滑
    mat_Ix2=mbys(mat_Ix2,cxDIB,cyDIB,g,gausswidth,gausswidth);
    mat_Iy2=mbys(mat_Iy2,cxDIB,cyDIB,g,gausswidth,gausswidth);
    mat_Ixy=mbys(mat_Ixy,cxDIB,cyDIB,g,gausswidth,gausswidth);
複製程式碼

利用opencv介面則是很簡單

複製程式碼 
//構建高斯函式
Mat gaussKernel = getGaussianKernel(7, 1);
//卷積計算
filter2D(Ix2, Ix2, CV_64F, gaussKernel);
filter2D(Iy2, Iy2, CV_64F, gaussKernel);
filter2D(Ixy, Ixy, CV_64F, gaussKernel);
複製程式碼

 

3.接下來利用M計算對應每個畫素的角點響應函式R，即：

也可以使用改進的R：
R=[Ix^2*Iy^2-(Ix*Iy)^2]/(Ix^2+Iy^2);裡面沒有隨意給定的引數k，取值應當比第一個令人滿意。

複製程式碼 
//計算cim：即cornerness of image，我們把它稱做‘角點量’
    CvMat *mat_cim;
    mat_cim=mbcim(mat_Ix2,mat_Iy2,mat_Ixy,cxDIB,cyDIB);
//用來求得響應度
CvMat *mbcim(CvMat *mat1,CvMat *mat2,CvMat *mat3,int xwidth,int ywidth)
{
    int i,j;
    CvMat *mat;
    mat=cvCloneMat(mat1);
    for(i = 0; i <ywidth; i++)
    {
        for(j = 0; j < xwidth; j++)
        {
            //注意：要在分母中加入一個極小量以防止除數為零溢位
            CV_MAT_ELEM(*mat,double,i,j)=(CV_MAT_ELEM(*mat1,double,i,j)*CV_MAT_ELEM(*mat2,double,i,j)-
                CV_MAT_ELEM(*mat3,double,i,j)*CV_MAT_ELEM(*mat3,double,i,j))/
                (CV_MAT_ELEM(*mat1,double,i,j)+CV_MAT_ELEM(*mat2,double,i,j)+0.000001);
        }
    }
    return mat;
}
複製程式碼 複製程式碼 
//opencv程式碼
Mat cornerStrength(gray.size(), gray.type());
    for (int i = 0; i < gray.rows; i++)
    {
        for (int j = 0; j < gray.cols; j++)
        {
            double det_m = Ix2.at<double>(i,j) * Iy2.at<double>(i,j) - Ixy.at<double>(i,j) * Ixy.at<double>(i,j);
            double trace_m = Ix2.at<double>(i,j) + Iy2.at<double>(i,j);
            cornerStrength.at<double>(i,j) = det_m - alpha * trace_m *trace_m;
        }
    }
複製程式碼

4、區域性極大值抑制，同時選取其極大值

複製程式碼 
//--------------------------------------------------------------------------
    //                 第四步：進行區域性非極大值抑制
    //--------------------------------------------------------------------------
    CvMat *mat_locmax;
    //const int size=7;
    mat_locmax=mblocmax(mat_cim,cxDIB,cyDIB,size);
//     cout<<CV_MAT_ELEM(*mat_locmax,double,cyDIB-1,cxDIB-1)<<endl;
//用來求得區域性極大值
CvMat *mblocmax(CvMat *mat1,int xwidth,int ywidth,int size)
{
    int i,j;
    double max=-1000;
    int i1,j1;
    int px,py;
    CvMat *mat;
    mat=cvCloneMat(mat1);
    for(i=size/2;i<ywidth-size/2;i++)
        for(j=size/2;j<xwidth-size/2;j++)
        {
            max=-10000;
            for(i1=0;i1<size;i1++)
                for(j1=0;j1<size;j1++)
                {
                    px=i-size/2+i1;
                    py=j-size/2+j1;
                    if(CV_MAT_ELEM(*mat1,double,px,py)>max)
                        max=CV_MAT_ELEM(*mat1,double,px,py);

                }
                if(max>0)
                    CV_MAT_ELEM(*mat,double,i,j)=max;
                else
                    CV_MAT_ELEM(*mat,double,i,j)=0;
        }
        return mat;
}
複製程式碼

在opencv中應用maxminloc函式

複製程式碼 
// threshold
double maxStrength;
minMaxLoc(cornerStrength, NULL, &maxStrength, NULL, NULL);
Mat dilated;
Mat localMax;
//預設3*3核膨脹，膨脹之後，除了區域性最大值點和原來相同，其它非區域性最大值點被  
 //3*3鄰域內的最大值點取代
dilate(cornerStrength, dilated, Mat());
//與原圖相比，只剩下和原圖值相同的點，這些點都是區域性最大值點，儲存到localMax 
compare(cornerStrength, dilated, localMax, CMP_EQ);
複製程式碼

5.在矩陣R中，同時滿足R(i,j)大於一定閾值threshold和R(i,j)是某領域內的區域性極大值，則被認為是角點。

複製程式碼 
cv::Mat cornerMap;  
  // 根據角點響應最大值計算閾值  
  threshold= qualityLevel*maxStrength;  
  cv::threshold(cornerStrength,cornerTh,  
   threshold,255,cv::THRESH_BINARY);  
  // 轉為8-bit圖 
 cornerTh.convertTo(cornerMap,CV_8U);// 和區域性最大值圖與，剩下角點區域性最大值圖，即：完成非最大值抑制  
  cv::bitwise_and(cornerMap,localMax,cornerMap);
imgDst = cornerMap.clone();
複製程式碼 
for( int y = 0; y < cornerMap.rows; y++ ) 
{  
        const uchar* rowPtr = cornerMap.ptr<uchar>(y);  
            for( int x = 0; x < cornerMap.cols; x++ ) 
    {  
               // 非零點就是角點  
              if (rowPtr[x]) 
    {  
                        points.push_back(cv::Point(x,y));  
                 }  
                 }  
 }
複製程式碼 複製程式碼

 

 

 

 

 

 

 

 

 

 

 

三、演算法原始碼

1、C版本

裡面函式最基礎的程式碼解釋和註釋

  C原始碼

 

2、基於opencv原始碼

  opencv+C++推薦

 

3、opencv中Harris介面

OpenCV的Hairrs角點檢測的函式為cornerHairrs()，但是它的輸出是一幅浮點值影象，浮點值越高，表明越可能是特徵角點，我們需要對影象進行閾值化。

C++: void cornerHarris(InputArray src, OutputArray dst, int blockSize, int apertureSize, double k, int borderType = BORDER_DEFAULT);

 

  • src – 輸入的單通道8-bit或浮點影象。
  • dst – 儲存著Harris角點響應的影象矩陣，大小與輸入影象大小相同，是一個浮點型矩陣。
  • blockSize – 鄰域大小。
  • apertureSize – 擴充套件的微分運算元大。
  • k – 響應公式中的，引數$\alpha$。
  • boderType – 邊界處理的型別。
複製程式碼 
int main()
{
  Mat image = imread("../buliding.png");
  Mat gray;
  cvtColor(image, gray, CV_BGR2GRAY);

  Mat cornerStrength;
  cornerHarris(gray, cornerStrength, 3, 3, 0.01);
  threshold(cornerStrength, cornerStrength, 0.001, 255, THRESH_BINARY);
  return 0;
}
複製程式碼

 

 

 

 

 

 

 

從上面上間一幅影象我們可以看到，有很多角點都是粘連在一起的，我們下面通過加入非極大值抑制來進一步去除一些粘在一起的角點。

非極大值抑制原理是，在一個視窗內，如果有多個角點則用值最大的那個角點，其他的角點都刪除，視窗大小這裡我們用3*3，程式中通過影象的膨脹運算來達到檢測極大值的目的，因為預設引數的膨脹運算就是用視窗內的最大值替代當前的灰度值。

  View Code

 

四、Harris角點性質

 

1、閾值決定檢測點數量

增大$\alpha$的值，將減小角點響應值$R$，降低角點檢測的靈性，減少被檢測角點的數量；減小$\alpha$值，將增大角點響應值$R$，增加角點檢測的靈敏性，增加被檢測角點的數量
 

2、Harris角點檢測運算元對亮度和對比度的變化不敏感

這是因為在進行Harris角點檢測時，使用了微分運算元對影象進行微分運算，而微分運算對影象密度的拉昇或收縮和對亮度的擡高或下降不敏感。換言之，對亮度和對比度的仿射變換並不改變Harris響應的極值點出現的位置，但是，由於閾值的選擇，可能會影響角點檢測的數量。

 

3. Harris角點檢測運算元具有旋轉不變性

Harris角點檢測運算元使用的是角點附近的區域灰度二階矩矩陣。而二階矩矩陣可以表示成一個橢圓，橢圓的長短軸正是二階矩矩陣特徵值平方根的倒數。當特徵橢圓轉動時，特徵值並不發生變化，所以判斷角點響應值也不發生變化，由此說明Harris角點檢測運算元具有旋轉不變性。

4. Harris角點檢測運算元不具有尺度不變性

如下圖所示，當右圖被縮小時，在檢測視窗尺寸不變的前提下，在視窗內所包含影象的內容是完全不同的。左側的影象可能被檢測為邊緣或曲線，而右側的影象則可能被檢測為一個角點。

 

五、函式備註

1、compare

功能：兩個陣列之間或者一個數組和一個常數之間的比較

結構：

void compare(InputArray src1, InputArray src2, OutputArray dst, int cmpop)


src1 ：第一個陣列或者標量，如果是陣列，必須是單通道陣列。
src2 ：第二個陣列或者標量，如果是陣列，必須是單通道陣列。
dst ：輸出陣列，和輸入陣列有同樣的size和type=CV_8UC1
cmpop ：

標誌指明瞭元素之間的對比關係
CMP_EQ src1 相等 src2.

CMP_GT src1 大於 src2.
CMP_GE src1 大於或等於 src2.
CMP_LT src1 小於 src2.
CMP_LE src1 小於或等於 src2.
CMP_NE src1 不等於 src2.

如果對比結果為true，那麼輸出陣列對應元素的值為255，否則為0

複製程式碼 
//獲取角點圖
    Mat getCornerMap(double qualityLevel) {
      Mat cornerMap;
      // 根據角點響應最大值計算閾值
      thresholdvalue= qualityLevel*maxStrength;
      threshold(cornerStrength,cornerTh,
      thresholdvalue,255,cv::THRESH_BINARY);
      // 轉為8-bit圖
      cornerTh.convertTo(cornerMap,CV_8U);
      // 和區域性最大值圖與，剩下角點區域性最大值圖，即：完成非最大值抑制
      bitwise_and(cornerMap,localMax,cornerMap);
      return cornerMap;
    }
複製程式碼

2、bitwise_and（位運算函式）

功能：計算兩個陣列或陣列和常量之間與的關係

結構：

void bitwise_and(InputArray src1, InputArray src2, OutputArray dst, InputArray mask=noArray())


src1 ：第一個輸入的陣列或常量
src2 ：第二個輸入的陣列或常量
dst ：輸出陣列，和輸入陣列有同樣的size和type
mask ：可選擇的操作掩碼，為8位單通道陣列，指定了輸出陣列哪些元素可以被改變，哪些不可以

操作過程為：

如果為多通道陣列，每個通道單獨處理

 

 

 

3、filter2D

OpenCV中的內部函式filter2D可以直接用來做影象卷積操作

void filter2D(InputArray src, OutputArray dst, int ddepth, InputArray kernel, Point anchor=Point(-1,-1), double delta=0, int borderType=BORDER_DEFAULT )

六、參考文章 

 

1、Opencv學習筆記（五）Harris角點檢測

2、尺度空間理論

3、影象特徵檢測：Harris

4、影象特徵ppt經典版

5、影象特徵ppt經典版

相關推薦

Harris角點檢測演算法詳解

Harris角點演算法 特徵點檢測廣泛應用到目標匹配、目標跟蹤、三維重建等應用中，在進行目標建模時會對影象進行目標特徵的提取，常用的有顏色、角點、特徵點、輪廓、紋理等特徵。現在開始講解常用的特徵點檢測，其中Harris角點檢測是特徵點檢測的基礎，提出了應用鄰近畫素

Harris角點檢測原理詳解外加原始碼分析

 網上有很多Harris角點檢測的資料，但是或多或少有一些缺陷，有的講清楚了原理，卻沒有給出程式碼的實現，有的給出了程式碼的實現，卻沒給出原理性的解釋。為此，我找出了幾篇部落格，將他們綜合整理在一起，有一套理論到程式碼的完整實現過程，這才是學一門演算法或運算元的正確方法。

Harris角點檢測實現詳解

本篇文章是對Harris角點檢測基於OpenCV 3.2.0的C++實現總結與記錄，程式程式碼力圖實現最基本的演算法步驟而無關優化，目的在於理解演算法，關於Harris角點的基本原理，本篇不予贅述，網上搜索即可，下面只列出出自Computer Vision:Al

Harris角點檢測演算法

Harris演算法是使用比較多的角點演算法之一，opencv中實現的角點檢測演算法就是Harris和它的改進。 Harris角點檢測的原理： Harris 角點最初由Harris在A combined corner and edge detector一文中提出,該文實

OpenCV計算機視覺學習（13）——影象特徵點檢測（Harris角點檢測，sift演算法）

如果需要處理的原圖及程式碼，請移步小編的GitHub地址 　　傳送門：請點選我 　　如果點選有誤：https://github.com/LeBron-Jian/ComputerVisionPractice 前言 　　特徵點檢測廣泛應用到目標匹配，目標跟蹤，三維重建等應用中，在進行目標建模時會對影象進行目標特徵

第十一節、Harris角點檢測原理

str 物體 per 權重 模式 windows www http 特定 OpenCV可以檢測圖像的主要特征，然後提取這些特征、使其成為圖像描述符，這類似於人的眼睛和大腦。這些圖像特征可作為圖像搜索的數據庫。此外，人們可以利用這些關鍵點將圖像拼接起來，組成一個更大的圖像，比

opencv-角點檢測之Harris角點檢測

trunc get data 圖像旋轉 ror 協方差矩陣 -a double 特定 轉自：https://blog.csdn.net/poem_qianmo/article/details/29356187 先看看程序運行截圖： 一、引言：關於興趣點（i

Opencv Harris角點檢測

demo one pla fault end ges for using opencv #include <iostream>#include <opencv2/opencv.hpp> using namespace std;using namesp

Moravec角點檢測演算法

一、Moravec角點檢測演算法原理     Moravec角點檢測演算法的思想是： 在影象中設計一個區域性檢測視窗，當該視窗沿各個方向作微小移動時 ， 考查視窗的平均能量變化，當該能量變化值超過設定的閾值時，就將視窗的中心畫素點提取為角點。    

OpenCv-C++-Harris角點檢測

首先，感謝賈志剛老師的課程教學。 理論部分： 我現在將程式碼放入： #include<opencv2/opencv.hpp> #include<iostream> #include<math.h> using namespace cv; using

OpenCV-Python系列·第二十五集：Harris角點檢測

Tip:Harris角點檢測 # -*- coding: utf-8 -*- """ Created on Sat Sep 22 21:18:41 2018 @author: Administrat

Harris角點檢測原理

1. 首先，我們不禁要問什麼是harris角點？        對於角點，到目前為止還沒有明確的數學定義。但是你可以認為角點就是極值點，即在某方面屬性特別突出的點。一般的角點檢測都是對有具體定義的、或者是能夠具體檢測出來的興趣點的檢測。這意味著興趣點可以是角點，是在某些屬性

OpenCV入門：Harris角點檢測

Harris角點檢測 #include <opencv2/opencv.hpp> #include <opencv2/imgproc/imgproc.hpp> using namespace cv; /* Harris角點檢測:void cornerH

Harris角點提取演算法及實現

角點：最直觀的印象就是在水平、豎直兩個方向上變化均較大的點，即Ix、Iy都較大 邊緣：僅在水平、或者僅在豎直方向有較大的變化量，即Ix和Iy只有其一較大 平坦地區：在水平、豎直方向的變化量均較小，即Ix、Iy都較小 2 strong eigenv

MATLAB——Harris角點檢測

為完成對棋盤標板的角點檢測，特此記錄。 角點：最直接的定義為在水平、垂直兩個方向上變化均較大的點 最常用的角點檢測演算法為Harris角點檢測方法，原理如下： 計算影象I(x,y)在X和Y兩個方向的梯度Ix、Iy 計算影象兩個方向梯度的乘積 使用高斯函式對Ix^2、

Python 實現Harris 角點檢測

具體原理見Feature extraction & image processing for computer vision (3rd edition) 4.4.1.4小節 # Harris Corner Detector # Code - 4.18

Harris角點檢測+sift描述匹配

    最近用到：Harris檢測影象中邊緣處的特徵點，並且需要兩張圖特徵點匹配起來。Harris只是一個角點檢測的演算法，最終只得到特徵點座標，想要匹配需要描述，而Harris是單尺度的，自己寫描述函式又有些麻煩。找到MATLAB和opencv都有整合的函式：

【OpenCV3經典程式設計100例】（24）2D特徵：Harris角點檢測cornerHarris（）

一、2d特徵相關知識Learn about how to use the feature points detectors, descriptors and matching inside OpenCV

Harris角點檢測原理分析

http://blog.csdn.net/newthinker_wei/article/details/45603583 http://www.360doc.com/content/15/1212/23/20007814_519967668.shtml 本文將該文拷

harris角點檢測計算詳細步驟

說明 某人說，演算法看不懂，公式不明白，舉個栗子咯，能看到具體數值的那種。 好吧，舉個最簡單的角點檢測演算法吧。以下程式碼為matlab。 步驟 一、生成影象 構造一個棋盤影