E(u,v)=∑(x,y)∈W [I(x+u,y+v)?I(x,y)]2(w是我们的窗口,[u,v]是我们的shift,也就是移动)
我们想了解微小移动对E到底有何影响
于是我们对I进行一阶泰勒展开
I(x+u,y+v)=I(x,y)+Ixu+Iyv+higherorder terms
≈I(x,y)+Ixu+ Iyv
=I(x,y)+[Ix Iy ][u v]T
再代入下面的公式
E(u,v)=∑[I(x+u,y+v)?I(x,y)]2(x,y)∈W
因此,我们可以写成
其中,M是一个二阶矩矩阵,可以由我们原始图片的差分得到
det = Ix2*Iy2 - Ixy^2 trace = Ix2 + Iy2[还有另一种经验公式,,我自己的实现就是采用了这一种]
1. 用sobel算子分别计算出水平梯度和垂直梯度
-1 0 1
-2 0 2
-1 0 1
-1 -2 -1
0 0 0
1 2 1
2. 计算高斯二阶矩阵,就是一个窗口里所有差分和
3.计算响应函数R
4.设置阀值
5.非极大值抑制,判断这个点的响应函数R是不是周围最大的
强调一下,harris对仿射变换只有部分不变性质。平移和旋转具有covariant,但是scaling不具有,如下图
初识API
|
对每个像素 计算 协方差矩阵 over a neighborhood. 用的是harris响应函数
#include "opencv2/highgui/highgui.hpp" #include "opencv2/imgproc/imgproc.hpp" #include <iostream> #include <stdio.h> #include <stdlib.h> using namespace cv; using namespace std; /// Global variables Mat src, src_gray; int thresh = 200; int max_thresh = 255; char* source_window = "Source image"; char* corners_window = "Corners detected"; /// Function header void cornerHarris_demo( int, void* ); /** @function main */ int main( int argc, char** argv ) { /// Load source image and convert it to gray src = imread( argv[1], 1 ); cvtColor( src, src_gray, CV_BGR2GRAY ); /// Create a window and a trackbar namedWindow( source_window, CV_WINDOW_AUTOSIZE ); createTrackbar( "Threshold: ", source_window, &thresh, max_thresh, cornerHarris_demo ); imshow( source_window, src ); cornerHarris_demo( 0, 0 ); waitKey(0); return(0); } /** @function cornerHarris_demo */ void cornerHarris_demo( int, void* ) { Mat dst, dst_norm, dst_norm_scaled; dst = Mat::zeros( src.size(), CV_32FC1 ); /// Detector parameters int blockSize = 2; int apertureSize = 3; double k = 0.04; /// Detecting corners cornerHarris( src_gray, dst, blockSize, apertureSize, k, BORDER_DEFAULT ); /// Normalizing normalize( dst, dst_norm, 0, 255, NORM_MINMAX, CV_32FC1, Mat() ); convertScaleAbs( dst_norm, dst_norm_scaled ); /// Drawing a circle around corners for( int j = 0; j < dst_norm.rows ; j++ ) { for( int i = 0; i < dst_norm.cols; i++ ) { if( (int) dst_norm.at<float>(j,i) > thresh ) { circle( dst_norm_scaled, Point( i, j ), 5, Scalar(0), 2, 8, 0 ); } } } /// Showing the result namedWindow( corners_window, CV_WINDOW_AUTOSIZE ); imshow( corners_window, dst_norm_scaled ); }
#include "opencv2/highgui/highgui.hpp" #include "opencv2/imgproc/imgproc.hpp" #include <cmath> #include <iostream> using namespace cv; Mat harris(Mat &im, double sigma, int thresh, int radius){ Mat dx,dy,Ix,Iy,Ix2,Iy2,Ixy,cim; Sobel( im, Ix, CV_64F, 1, 0, 3); //算法第一步,计算水平垂直差分 Sobel( im, Iy, CV_64F, 0, 1, 3); int ksize = max(1, (int)(6*sigma)); if(ksize % 2 == 0) ksize++; GaussianBlur(Ix.mul(Ix), Ix2, Size(ksize, ksize), sigma); //算法第二步,计算二阶高斯差分矩阵 GaussianBlur(Iy.mul(Iy), Iy2, Size(ksize, ksize), sigma); GaussianBlur(Ix.mul(Iy), Ixy, Size(ksize, ksize), sigma); //Harris corner measure //cim = (Ix2.*Iy2 - Ixy.^2)./(Ix2 + Iy2); cim = (Ix2.mul(Iy2) - Ixy.mul(Ixy)) / (Ix2+Iy2); //算法第三步,计算响应函数,我使用了另外一种 Mat structedElement(radius, radius, CV_8U, Scalar(1)); Mat mx,norm_cim; normalize( cim, norm_cim, 0, 255, NORM_MINMAX, CV_8U, Mat() ); dilate(norm_cim, mx, structedElement); norm_cim = ( norm_cim == mx) & (norm_cim>thresh); //算法第4第5步融合,非极大值抑制和阀值检测 return norm_cim; } int main( int, char** argv ) { Mat src,gray; src = imread( argv[1] ); cvtColor( src, gray, CV_RGB2GRAY ); Mat corners = harris(gray, 1.5, 30, 2); for( int j = 0; j < corners.rows ; j++ ) { for( int i = 0; i < corners.cols; i++ ) { if( corners.at<unsigned char>(j,i) > 0) { circle( gray, Point( i, j ), 3, Scalar(0), 2, 8, 0 ); } } } namedWindow("result", 1); imshow("result", gray); waitKey(); return 0; }
openCV2马拉松第19圈——Harris角点检测(自己实现),布布扣,bubuko.com
openCV2马拉松第19圈——Harris角点检测(自己实现)
原文:http://blog.csdn.net/abcd1992719g/article/details/26824529