精华内容
下载资源
问答
  • harris角点检测优化算法
    千次阅读
    2015-06-28 16:56:07

    Harris角点检测算法优化

    一、综述

    用 Harris 算法进行检测,有三点不足:

    (1 )该算法不具有尺度不变性;

    (2 )该算法提取的角点是像素级的;

    (3 )该算法检测时间不是很令人满意。

    基于以上认识,我们主要针对第(3)点对Harris 角点检测算法提出了改进。


    二、改进 Harris 算法原理

          在介绍方法之前,我们先提出如下概念:图像区域像素的相似度

    我们知道, Harris角点检测是基于图像像素灰度值变化梯度的, 灰度值图像的角点附近,是其像素灰度值变化非常大的区域,其梯度也非常大。

    换句话说,在非角点位置邻域里,各点的像素值变化不大,甚至几乎相等,其梯度相对也比较小。从这个角度着眼,于是提出了图像区域像素的相似度的概念,它是指检测窗口中心点灰度值与其周围n 邻域内其他像素点灰度值的相似程度,这种相似程度是用其灰度值之差来描述的。如果邻域内点的灰度值与中心点Image (i,j) 的灰度值之差 的绝对值 在一个阈值t 范围 内,那就认为这个点与中心点是相似的。

    与此同时,属于该Image (i,j) 点的相似点计数器nlike(i,j) 也随之加一。在Image (i,j) 点的n 邻域 全部被遍历一边之后,就能得到在这个邻域范围内 与中心点相似的 点个数的统计值nlike(i,j) 。根据nlike(i,j)的大小,就可以判断这个中心点是否可能为角点。

    由于我们选择3*3 的检测窗口,所以,对于中心像素点 , 在下面的讨论中只考虑其8 邻域内像素点的相似度。计算该范围的像素点与中心像素点的灰度值之差的绝对值 ( 记为 Δ ) , 如果该值小于等于设定的阈值( 记为 t), 则认为该像素点与目标像素点相似 。


    nlike( i , j ) = sum( R ( i+ x , j + y ) )   ,( - 1 ≤ x ≤ 1 , -  1 ≤ y≤ 1 , 且 x ≠ 0 , y ≠ 0),

    其中 : 1 , Δ( i + x, j + y)  ≤ t  ,R(i+x, j+y)= 0 ,  Δ( i + x , j + y)  >  t


    从定义中可以看出 : 0 ≤nlike ( i , j ) ≤ 8 。 


    现在讨论 nlike( i , j) 值的含义 :

    (1)nlike (i , j) = 8 , 表示当前中心像素点的 8 邻域范围内都是与之相似的像素点 , 所以该像素点邻域范围内的梯度不会很大 , 因此角点检测时 , 应该排除此类像素点,不将其作为角点的候选点 。


    (2)nlike (i , j) = 0 , 表示当前中心像素点的 8 邻域范围内没有与之相似的像素点 , 所以该像素点为孤立像素点或者是噪声点 , 因此角点检测时 , 也应该排除此类像素点 。


    (3)nlike (i , j) = 7 , 可以归结为以下的两者情况 , 其他情形都可以通过旋转来得到( 图中黑色区域仅表示与中心像素相似 , 而两个黑色区域像素可能是相似的, 也可能不相似 ) 。 对于图 1 (a) 中 , 可能的角点应该是中心像素点的正上方的那个像素点 , 1(b) 图中可能的角点应该是中心像素点右上方的那个像素点 , 故这种情况下, 中心像素点不应该作为角点的候选点 。

     


    (4)nlike (i , j) = 1 ,可以归结为图2中的两种情况 ( 图中白色区域仅表示与中心像素不相似 , 而两个白色区域像素可能是相似的 , 也可能不相似 ) , 在这两种情况下 , 中心像素点也不可能为角点 。

      


    (5)  2 ≤ nlike ( i , j) ≤ 6 , 情况比较复杂 , 无法确认像素点准确的性质。我们采取的方法是先将其列入候选角点之列,对其进行计算CRF 等后续操作。 









    更多相关内容
  • 针对这一情况,文章在传统Harris算法基础上提出一种新的检测方法,采用多阈值的圆形非极大值抑制法提取角点, 以此降低算法检测时间并增强图像旋转不变性,再借鉴SUSAN思想消去大部分伪角点。通过实验对比,该算法...
  • 摘要:针对Harris角点检测算法中提取出较多的伪角点和计算量大的问题,提出了一种基于Harris角点检测的改进算法. 为抑制 Harris 角点检测中的伪角点
  • 基于matlab的改进的harris角点检测算法,与harris算法可以对比。
  • 在研究Harris角点检测算法时,发现检测出的角点常常会受到噪声的影响。通过优化Har-ris角点响应函数,避免了角点响应函数中k值的影响,提高了每个目标像素的响应值精度。把该方法运用到人脸特征的检测中,实验结果表明...
  • harris角点检测原理

    千次阅读 2021-07-23 22:23:36
    基于图像灰度的方法通过计算的曲率及梯度来检测角点。 2、数学知识 1)泰勒展开 泰勒展开公式是一种统一的形式,非常完美。 一维泰勒展开公式: 二维泰勒展开公式: 2)矩阵的特征值和特征向量 ...

    目录

    1、角点概述

    2、数学知识

    3、Harris角点检测基本原理

    4、优化改进


    1、角点概述

            如果一个点在任意方向的一个微小变动都会引起灰度很大的变化,那么我们就把它称之为角点,也就是一阶导数(即灰度图的梯度)中的局部最大所对应的像素点就是角点。在现实世界中,角点对应于物体的拐角,道路的十字路口、丁字路口等。基于图像灰度的方法通过计算点的曲率及梯度来检测角点。

    2、数学知识

    1)泰勒展开

            泰勒展开公式是一种统一的形式,非常完美。

            一维泰勒展开公式:

            二维泰勒展开公式:

    2)矩阵的特征值和特征向量

     https://zhuanlan.zhihu.com/p/42490675

    3、Harris角点检测基本原理

            人眼对角点的识别通常是在一个局部的小区域或小窗口完成的。如果在各个方向上移动这个特征的小窗口,窗口内区域的灰度发生了较大的变化,那么就认为在窗口内遇到了角点。如果这个特定的窗口在图像各个方向上移动时,窗口内图像的灰度没有发生变化,那么窗口内就不存在角点;如果窗口在某一个方向移动时,窗口内图像的灰度发生了较大的变化,而在另一些方向上没有发生变化,那么,窗口内的图像可能就是一条直线的线段。

            Harris 检测器具有旋转不变性,但不具有尺度不变性,也就是说尺度变化可能会导致角点变为边缘,如下图所示:

            想要尺度不变特性的话,可以关注SIFT特征。

    Harris 角点检测算法分为以下三步:

            1、当窗口同时向 x 和 y 两个方向移动时,计算窗口内部的像素值变化量E(u,v);

            2、对于每个窗口,都计算其对应的一个角点响应函数R;

            3、然后对该函数进行阈值处理,如果R > threshold,表示该窗口对应一个角点特征。

    E(u,v)推导过程:

            首先,将图像窗口平移[u,v]产生灰度变化的自相关函数如下:

            其中窗口函数(权重矩阵)可以是平坦的,也可以是高斯的,是一个二维的滤波器。对于一个角点来说, E(u,v)会非常大。因此,我们可以最大化上面这个函数来得到图像中的角点。用上面的函数计算会非常慢。因此,我们使用泰勒展开式(只有一阶)来得到这个公式的近似形式。

            将平移后的式子进行泰勒展开如下:

             其中Ix和Iy是I的偏微分,在图像中就是在x和y方向的梯度图(可以通过cv2.Sobel()来得到):

            接下来继续推导:

             把u和v拿出来,得到最终的形式:

            其中矩阵M为:

             最后是把实对称矩阵对角化处理后的结果,可以把R看成旋转因子,其不影响两个正交方向的变化分量。经对角化处理后,将两个正交方向的变化分量提取出来,就是 λ1 和 λ2(特征值)。

            对于图像的每一个像素点(x,y),对应一个以该像素为中心的窗口w(x,y),然后该像素平移(u,v)得到新的像素点(x+u,y+v),而E(u,v)就是窗口中所有像素的加权和乘以不同位置像素的灰度差值。

            矩阵M又称为Harris矩阵。w(x,y)的宽度决定了在像素x 周围的感兴趣区域。

    计算响应函数R:

            得到E(u,v)的最终形式,我们的目的是要找到会引起较大的灰度值变化的那些窗口。灰度值变化的大小则取决于矩阵M,那么如何找到这些窗口,我们可以使用矩阵的特征值来实现。

            忽略余项之后的表达式为一个二项式函数,然而二项式函数的本质上就是一个椭圆函数,椭圆的扁率和尺寸是由M(x,y)的特征值λ1、λ2决定的,椭圆的方向是由M(x,y)的特征矢量决定的,如下图所示,椭圆方程为:

             椭圆函数特征值与图像中的角点、直线(边缘)和平面之间的关系如下图所示。共可分为三种情况:

            a)图像中的直线。一个特征值大,另一个特征值小,λ1>λ2或λ2>λ1。自相关函数值在某一方向上大,在其他方向上小。

            b)图像中的平面。两个特征值都小,且近似相等;自相关函数数值在各个方向上都小。

            c)图像中的角点。两个特征值都大,且近似相等,自相关函数在所有方向都增大。

            通过M的两个特征值λ1和λ2的大小对图像点进行分类:

             如果λ1和λ2都很小,图像窗口在所有方向上移动都无明显灰度变化。由于我们是通过M的两个特征值的大小对图像进行分类,所以,定义角点相应函数R:

              其中k为经验常数,一般取k=0.04~0.06。为了去除加权常数κ,我们通常使用商数detM/(traceM)2作为指示器。所以,上图可以转化为:

             因为特征值λ1和λ2决定了R的值,R 只与M的特征值有关,所以我们可以用特征值来决定一个窗口是平面、边缘还是角点。

            平面:该窗口在平坦区域上滑动,窗口内的灰度值基本不会发生变化,所以|R|值非常小,在水平和竖直方向的变化量均较小,即Ix和Iy都较小,那么λ1和λ2都较小;

            边缘:R值为负数,仅在水平或竖直方向有较大的变化量,即Ix和Iy只有一个较大,也就是λ1>>λ2或λ2>>λ1;

            角点:R值很大,在水平、竖直两个方向上变化均较大的点,即Ix和Iy都较大,也就是λ1和λ2都很大。

    最优角点判别:

            根据R的值,将这个窗口所在的区域划分为平面、边缘或角点。为了得到最优的角点,我们还可以使用非极大值抑制。        

            Harris角点检测的结果是带有这些分数R的灰度图像,设定一个阈值,R > threshold,分数大于这个阈值的像素就对应角点。

    4、优化改进

    1)由于Harris角点检测算法的稳定性和k值有关,而k是个经验值,不好设定最佳值。

            Shi-Tomasi发现,角点的稳定性其实和矩阵M的较小特征值有关,于是直接用较小的那个特征值作为分数。这样就不用调整k值了。

            所以Shi-Tomasi将分数公式改为如下形式:

            和Harris一样,如果该分数大于设定的阈值,我们就认为它是一个角点。

    2)Harris和Shi-Tomasi都是基于梯度计算的角点检测方法,Shi-Tomasi的效果要好一些。基于梯度的检测方法有一些缺点: 计算复杂度高,图像中的噪声可以阻碍梯度计算。

            想要提高检测速度的话,可以考虑基于模板的方法:FAST角点检测算法。该算法原理比较简单,但实时性很强。

    3)Harris 检测器具有旋转不变性,但不具有尺度不变性,也就是说尺度变化可能会导致角点变为边缘,如下图所示:

            想要尺度不变特性的话,可以关注SIFT特征。

    相关链接:

    1、harris角点检测算法实现

    2、SHI-TOMASI角点检测

    展开全文
  • Harris 角点检测是基于图像像素灰度值变化梯度的, 灰度值图像的角点附近,是其像素灰度值变化非常大的区域,其梯度也非常大。换句话说,在非角点位置邻域里,各的像素值变化不大,甚至几乎相等,其梯度相对也比较...
  • 2、Harris角点检测算法

    千次阅读 2018-07-26 16:23:39
    (1)Harris算子用高斯函数代替Moravec算子的二值窗口函数,如图1所示,窗口函数应对离中心越的像素赋予越大的权重,以减少噪声影响; 图1 对于图像I(x,y),当在(x,y)处平移(u,v)后的自相似性,可以通过自相关...

    1、基本原理

    Harris算子是对Moravec算子的改进,包括:

    (1)Harris算子用高斯函数代替Moravec算子的二值窗口函数,如图1所示,窗口函数应对离中心点越的像素赋予越大的权重,以减少噪声影响;

    图1

    对于图像I(x,y),当在点(x,y)处平移(u,v)后的自相似性,可以通过自相关函数给出:

    E(u,v)=\sum_{x,y}\omega (x,y)[I(x+u,y+v)-I(x,y)]^2   (1.1)

    其中,w(x,y)是以点(x,y)为中心的窗口加权函数,它既可是常数(图1左),也可以是高斯加权函数(图1右)。高斯函数的表达式如下:

    \omega (x,y)=\frac{1}{2\pi \sigma ^2}e^{-(x^2+y^2)/2\sigma ^2}    (1.2)

    (2)Moravec算子只考虑每隔45度方向,Harris用Taylor展开去近似任意方向;

    利用泰勒级数将I(x+u,y+v)展开得:

    I(x+u,y+v)=I(x,y)+I_xu+I_yv+o(u^2,v^2)     (1.3)

    则E(u,v)表达式可以更新为:

    E(u,v)\cong \sum_{x,y}\omega (x,y)[I_xu+I_yv]^2       (1.4)

    其中,Ix为x方向的差分(一阶微分近似),Iy为y方向的差分,w(x,y)为高斯函数。通过推导,E(u,v)以矩阵的形式表示为:

    E(u,v)\cong [u,v]M[u,v]^T   (1.5)

    矩阵M为实对称矩阵,且

    M=\omega (x,y)\bigotimes \begin{bmatrix} I_x^2 & I_xI_y\\ I_xI_y& I_y^2 \end{bmatrix}=\begin{bmatrix} A & C\\ C& B \end{bmatrix}      (1.6)

    则E(u,v)化简为二次项函数:

    E(u,v)=Au^2+Bv^2+2Cuv          (1.7)

    二次项函数本质上就是一个椭圆函数。椭圆的扁率和尺寸是由M的特征值\lambda _1,\lambda _2决定的,椭圆的方向是由M的特征向量决定的,如图2所示。

    图2

    椭圆函数特征值与图像中的角点、直线(边缘)和平面之间的关系如图3所示。共可分为三种情况:

    1)图像上的直线。一个特征值比较大,另外一个特征值比较小;自相关函数在某一个方向上大,在另一个方向上小。

    2)图像中的平滑区域。两个特征值都小,且近似相等;自相关函数值子各个方向上都小。

    3)图像中的角点。两个特征值都比较大,且近似相等,自相关函数在所有方向上都增大。

    图3

    (3)角点响应函数

    根据二次项函数特征值的计算公式,我们可以求M矩阵的特征值。但是Harris给出的角点差别方法并不需要计算具体的特征值,而是计算一个角点响应值R来判断角点。R的计算公式为:

    R=detM-k(traceM)^2       (1.8)

    其中,detM为矩阵M的行列式,traceM为矩阵M的迹,k为常数,取值范围为0.04~0.06。事实上,特征是隐含在detM和traceM中,因为

    detM=\lambda _1\lambda _2=AB-C^2        (1.9)

    traceM=\lambda _1+\lambda _2=A+B        (1.10)

    其实,角点量R的计算方式可以自由发挥,只要能反应角点的特征即可。例如,Nobel于1988年提出利用如下公式计算角点的响应值,无需设定参数k的值:

    cim=\frac{I_x^2I_y^2-(I_xI_y)^2}{I_x^2+I_y^2}             (1.11)

    采用上述公式计算角点的CRF值,从而避免的参数k对角点选取的影响,在实际应用中,通常选用这个改进的Harris角点检测算法进行检测:当cim值大于预定的阈值,则该点为角点候选点,通过非极大值抑制挑选出最终的角点。

    此外,或如式(1.12)计算R值,也无需考虑参数k。

    R=\frac{detM}{traceM}           (1.12)

    (4)最后设定R的阈值进行角点判断,以及角点的极大值抑制等。

    二、算法步骤

    输入:源单通道图像,参数k

    输出:角点检测图

    具体步骤:

    (1)计算图像I(x,y)在X和Y方向的梯度Ix和Iy;

    I_x=\frac{\partial I}{\partial x}=I\bigotimes (-1,0,1)

    I_y=\frac{\partial I}{\partial x}=I\bigotimes (-1,0,1)^T

    (2)计算梯度方向的乘积;

    I_x^2=I_x\cdot I_x

    I_y^2=I_y\cdot I_y

    I_x_y=I_x\cdot I_y

    (3)使用高斯核进行加权,计算矩阵M的元素A,B,C;

    A=g(I_x^2)=I_x^2\bigotimes \omega

    B=g(I_y^2)=I_y^2\bigotimes \omega

    C=g(I_x_y)=I_x_y\bigotimes \omega

    (4)计算角点响应函数R,并设定阈值,当R小于阈值时,不是候选角点;

    R=(AB-C^2)-k(A+B)^2

    (5)进行局部极大值抑制。

    算法结束。

    三、优缺点

    Harris角点检测算法有诸多优点,但也有不完善的地方。

    (1)Harris角点检测算子具有旋转不变性

    Harris角点检测算子使用的是角点附近的区域灰度二阶矩矩阵。而二阶矩矩阵可以表示成一个椭圆,椭圆的长短轴正是二阶矩矩阵特征值平方根的倒数。当特征椭圆转动时,特征值并不发生变化,所以判断角点响应值R也不发生变化,由此说明Harris角点检测算子具有旋转不变性。

    (2)Harris角点检测算子对灰度平移和灰度尺度变化不敏感

    这是因为在进行Harris角点检测时,使用了微分算子对图像进行微分运算,而微分运算对图像亮度的抬高或下降(I=I+a)、密度的拉升或收缩(I=bI)不敏感。换言之,对亮度和对比度的仿射变换并不改变Harris响应的极值点出现的位置,但是,由于阈值的选择,可能会影响角点检测的数量。

    (3)Harris角点检测算子不具有尺度不变性

    如图4所示,当图像被缩小时,在检测窗口尺寸不变的前提下,窗口内所包含图像的内容可能是完全不同的。左侧的图像可能被检测为边缘或曲线,而右侧的图像则可能被检测为一个角点。或者说如果图像尺度发生变化,原来是角点的点在新的尺度可能就不是角点了。

    图 4

    注:尺度不变性问题可通过图像金字塔解决,例如,在运算的开始先将图像转化到尺度空间表示,即将原图像进行尺度变换,而尺度变换的方式就是源图像与尺度核函数做卷积运算:

    L(x,y,\sigma )=g(X,Y,\sigma )\bigotimes I(x,y)

    其中,sigma表示尺度。然后使用L代替原图像去进行运算,尺度为运算的参数。关于图像金字塔的具体内容回继续学习!

    Harris角点本身就不受光照,旋转的影响,现在又使其满足尺度不变性,至此,Harris角点可以成为一个优秀的特征了。

    四、代码分析

    /**********************************************************************************
    *函数 Mat detectHarrisCorners(const Mat& imgSrc, double alpha)
    *输入:
          *imgSrc : 源单通道图像
          *alpha  : Harris响应函数参数
    *输出
         *imgDst : 提取到角点的图像
    ***************************************************************************************/
    
    #include <iostream> 
    #include <vector>
    #include "opencv2/opencv.hpp"
    #include "opencv2/imgproc/imgproc.hpp"  
    #include "opencv2/highgui/highgui.hpp"   
    using namespace cv;
    using namespace std;
    
    Mat detectHarrisCorners(const Mat& imgSrc, double alpha)
    {
    	Mat gray;
    	gray = imgSrc.clone();
        gray.convertTo(gray, CV_64F);
    
    	Mat xKernel = (Mat_<double>(1, 3) << -1, 0, 1);//水平方向模板计算Ix 
    	Mat yKernel = xKernel.t();
    	Mat Ix, Iy;
    	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);
    
    	Mat gaussKernel = getGaussianKernel(5, 1);//获得高斯核size=5,sigma=1
    	filter2D(Ix2, Ix2, CV_64F, gaussKernel);
    	filter2D(Iy2, Iy2, CV_64F, gaussKernel);
    	filter2D(Ixy, Ixy, CV_64F, gaussKernel);
    
    	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;//响应函数值R
    		}
    	}
    	
    	double maxStrength;
    	minMaxLoc(cornerStrength, NULL, &maxStrength, NULL, NULL);
    	double qualityLevel = 0.1;
    	double thresh = qualityLevel * maxStrength;// 设置threshold
     
    	Mat dilated, localMax;
    	//默认3 * 3核膨胀,膨胀之后,除了局部最大值点和原来相同,其它非局部最大值点被3*3邻域内的最大值点取代
    	dilate(cornerStrength, dilated, Mat()); 
    	//与原图相比,只剩下和原图值相同的点,这些点都是局部最大值点,保存到localMax 
    	compare(cornerStrength, dilated, localMax, CMP_EQ);
    	//和局部最大值图与,剩下角点局部最大值图,即:完成非最大值抑制
    	Mat cornerMap;
    	cornerMap = cornerStrength > thresh;
    	bitwise_and(cornerMap, localMax, cornerMap);
    
    
    	vector<Point> points;
    	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(Point(x,y));
    		}
    	}
    	//画角点
    	Mat gray1 = imgSrc.clone();
    	vector<Point>::const_iterator it = points.begin();
    	while (it != points.end())
    	{
    		circle(gray1, *it, 3, Scalar(255, 255, 255), 1);
    		++it;
    	}
    
    	return gray1;
    }
    
    int main()
    {
    	Mat SrcImage = imread("lena.png");
    	if (!SrcImage.data)
    		return -1;
        Mat gray;
    	if (SrcImage.channels() == 3)
    	{
    		cvtColor(SrcImage, gray, CV_BGR2GRAY);  //RGB转换成灰度图像
    	}
    	else
    	{
    		gray = SrcImage.clone();
    	}
    
    	double time0 = static_cast<double>(getTickCount());
    	Mat HarrisImage = detectHarrisCorners(gray, 0.05);
    
    	time0 = ((double)getTickCount() - time0) / getTickFrequency();
    	cout << "runtime :" << time0 << "s" << endl;
    	imshow("HarrisImage", HarrisImage);
    	waitKey(0);
    	return 0;
    }
    

    实验结果图如下,代码是拼凑的并未优化,运行时间还是比较长的。阈值设置的比较大,检测到的角点数量也比较少。

    图5

     

    展开全文
  • 针对传统Otsu算法只用于单阈值分割的不足,将Otsu算法推广到多阈值彩色图像分割中,提出先在众多极大值中寻找有意义峰值,根据峰值将直方图划分成多个待分割区间,再在每个区间进行阈值选取的方法;并且综合运用了...
  • harris角点检测的简要总结

    千次阅读 2019-04-13 23:32:30
    简要介绍了harris角点检测的原理与具体的实现过程。

    1. 概述相关

    harris角点检测是一种特征提取的方法,而特征提取正是计算机视觉的一种重要手段。尽管它看起来很复杂,其实也是基于数学原理和简单的图像处理来实现的。
    本文之前可以参看笔者写的几篇图像处理的文章,将会有助于更深入了解harris角点检测的实现。

    1. 图像的卷积(滤波)运算(一)——图像梯度
    2. 图像的卷积(滤波)运算(二)——高斯滤波
    3. 图像的膨胀与腐蚀——OpenCV与C++的具体实现

    2. 原理详解

    1) 算法思想

    为了判断图像的角点,可以利用卷积窗口滑动的思想,让以该点为中心的窗口在附近滑动。如下图是所有描述角点文章的初始图例,它表征的正是这一特性:当滑动窗口在所有方向移动时,窗口内的像素灰度出现了较大的变化,就可能是角点。
    在这里插入图片描述

    2) 数学模型

    根据上述的算法思想,可以构建数学模型,图像窗口平移[u,v]产生灰度变化E(u,v)为:
    在这里插入图片描述
    其中w(x,y)是一种加权函数,几乎所有的应用都把它设为高斯函数。由上述公式,进行推导如下:
    在这里插入图片描述
    最后得到的公式(6),在几何意义上表征的是一个椭圆。椭圆的长短轴分别沿着矩阵M的两个特征向量的方向,而两个与之对应的特征值分别是半长轴和半短轴的长度的平方的倒数。
    在这里插入图片描述
    那么根据矩阵M的两个特征值λ1和λ2,可以将图像上的像素点分类成直线、平面与角点:当λ1和λ2 都比较大,且近似相等时,可以认为是角点。如下图所示:
    在这里插入图片描述

    3) 优化推导

    而上述表达不太方便使用,又定义了一个角点响应函数R,通过R的大小来判断像素是否为角点:
    在这里插入图片描述
    式中,detM为矩阵M的行列式,traceM为矩阵M的直迹。α为经常常数,取值范围为0.04~0.06。对于R公式,有推导如下:
    在这里插入图片描述
    可以知道,角点响应值R仍然表征了矩阵M两个特征值λ1和λ2,同样可以进行上述分类:当R为大数值正数的时候,表示为角点。如下图所示:
    在这里插入图片描述

    3. 具体实现

    在OpenCV中,已经提供了Harris角点检测函数cornerHarris()。为了更好地理解Harris角点提取的原理,这里参考了网上代码,自己实现了其算法,不过也调用了OpenCV中一些基本函数。
    根据上述原理,Harris图像角点检测算法的关键是计算M矩阵,M矩阵是图像I(x,y)的偏导数矩阵,也就是要先求出图像的梯度。

    1) 详细步骤

    1.计算图像I(x,y)在X,Y方向的梯度。在这里是通过卷积函数filter2D实现的,具体原理可以看(1)中提到的相关文章。

    Mat gray;
    imgSrc.convertTo(gray, CV_64F);
    
    Mat xKernel = (Mat_<double>(1, 3) << -1, 0, 1);
    Mat yKernel = xKernel.t();
    
    Mat Ix, Iy;
    filter2D(gray, Ix, CV_64F, xKernel);
    filter2D(gray, Iy, CV_64F, yKernel);
    

    2.计算图像两个方向梯度的乘积。

    Mat Ix2, Iy2, Ixy;
    Ix2 = Ix.mul(Ix);
    Iy2 = Iy.mul(Iy);
    Ixy = Ix.mul(Iy);
    

    3.对Ix2、Iy2和Ixy进行高斯滤波,生成矩阵M的元素A、B和C。

    Mat gaussKernel = getGaussianKernel(7, 1);
    filter2D(Ix2, Ix2, CV_64F, gaussKernel);
    filter2D(Iy2, Iy2, CV_64F, gaussKernel);
    filter2D(Ixy, Ixy, CV_64F, gaussKernel);
    

    4.根据公式计算每个像素的Harris响应值R,得到图像对应的响应值矩阵。

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

    5.在3×3的邻域内进行非最大值抑制,找到局部最大值点,即为图像中的角点。在这里非最大值抑制是通过图像膨胀的实现的。比较膨胀前后的响应值矩阵,得到局部最大值。

    //在3×3的邻域内进行非最大值抑制,找到局部最大值点,即为图像中的角点
    double maxStrength;
    minMaxLoc(cornerStrength, NULL, &maxStrength, NULL, NULL);
    Mat dilated;
    Mat localMax;
    dilate(cornerStrength, dilated, Mat());				//膨胀
    compare(cornerStrength, dilated, localMax, CMP_EQ);      //比较保留最大值的点
    
    //得到角点的位置
    Mat cornerMap;
    double qualityLevel = 0.01;
    double thresh = qualityLevel * maxStrength;   
    cornerMap = cornerStrength > thresh;				//小于阈值t的R置为零。
    bitwise_and(cornerMap, localMax, cornerMap);			//位与运算,有0则为0, 全为1则为1
    
    imgDst = cornerMap.clone();
    

    2) 最终实现

    合并以上步骤,传入参数,最终的实现代码:

    #include <iostream>
    #include <algorithm>
    #include <opencv2\opencv.hpp>
    
    using namespace cv;
    using namespace std;
    
    void detectHarrisCorners(const Mat& imgSrc, Mat& imgDst, double alpha)
    {
    	//
    	Mat gray;
    	imgSrc.convertTo(gray, CV_64F);
    
    	//计算图像I(x,y)在X,Y方向的梯度
    	Mat xKernel = (Mat_<double>(1, 3) << -1, 0, 1);
    	Mat yKernel = xKernel.t();
    
    	Mat Ix, Iy;
    	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);
    
    	//对Ix2、Iy2和Ixy进行高斯滤波,生成矩阵M的元素A、B和C。
    	Mat gaussKernel = getGaussianKernel(7, 1);
    	filter2D(Ix2, Ix2, CV_64F, gaussKernel);
    	filter2D(Iy2, Iy2, CV_64F, gaussKernel);
    	filter2D(Ixy, Ixy, CV_64F, gaussKernel);
    
    	//根据公式计算每个像素的Harris响应值R,得到图像对应的响应值矩阵。
    	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;
    		}
    	}
    
    	//在3×3的邻域内进行非最大值抑制,找到局部最大值点,即为图像中的角点
    	double maxStrength;
    	minMaxLoc(cornerStrength, NULL, &maxStrength, NULL, NULL);
    	Mat dilated;
    	Mat localMax;
    	dilate(cornerStrength, dilated, Mat());				//膨胀
    	compare(cornerStrength, dilated, localMax, CMP_EQ);      //比较保留最大值的点
    	
    	//得到角点的位置
    	Mat cornerMap;
    	double qualityLevel = 0.01;
    	double thresh = qualityLevel * maxStrength;   
    	cornerMap = cornerStrength > thresh;				//小于阈值t的R置为零。
    	bitwise_and(cornerMap, localMax, cornerMap);			//位与运算,有0则为0, 全为1则为1
    
    	imgDst = cornerMap.clone();
    }
    
    //在角点位置绘制标记
    void drawCornerOnImage(Mat& image, const Mat&binary)
    {
    	Mat_<uchar>::const_iterator it = binary.begin<uchar>();
    	Mat_<uchar>::const_iterator itd = binary.end<uchar>();
    	for (int i = 0; it != itd; it++, i++)
    	{
    		if (*it)
    			circle(image, Point(i%image.cols, i / image.cols), 3, Scalar(0, 255, 0), 1);
    	}
    }
    
    int main()
    {
    	//从文件中读取成灰度图像
    	const char* imagename = "D:\\Data\\imgDemo\\whdx.jpg";
    	Mat img = imread(imagename, IMREAD_GRAYSCALE);
    	if (img.empty())
    	{
    		fprintf(stderr, "Can not load image %s\n", imagename);
    		return -1;
    	}
    
    	//
    	Mat imgDst;
    	double alpha = 0.05;	
    	detectHarrisCorners(img, imgDst, alpha);
    	
    	//在角点位置绘制标记
    	drawCornerOnImage(img, imgDst);
    
    	//
    	imshow("Harris角点检测", img);
    	waitKey();
    	
    	return 0;
    }
    

    其运行结果为:
    在这里插入图片描述

    4. 参考文献

    1. Harris角点
    2. Harris角点算法
    3. 矩阵特征值和椭圆长短轴的关系?-Eathen的回答
    展开全文
  • Harris 算子是一种计算简单、应用广泛的角点检测算子,只使用灰度的一阶差分和滤波,可以定量地提取特征角点并且提取的角点特征均匀。由于计算过程用到图像的一阶导数,故对存在灰度变化、图像旋转、视点变换和噪声...
  • 角点检测算法的原始思想: 我们可以从角点具有的特征出发: 即选取一个局部窗口,将这个窗口沿着各个方向移动,计算移动前后窗口内像素的差异的多少进而判断窗口对应的区域是否是角点。 有了基本的对角点的描述思想...
  • 文章目录一、Harris角点检测1.Harris角点检测算法2.何为角点3.对比Moravec角点检测算子二、不同场景下的实验分析1.纹理平坦场景2.1.1原图2.1.2光线变暗2.1.3侧面2.1.4变近2.1.5旋转2.1.6实验结果分析2.垂直或水平...
  • Harris角点检测

    2019-03-08 22:57:03
    在图像上进行Harris角点检测,对于每个像素 ( x , y ) (x,y) ( x , y ) ,计算2x2的梯度的协方差矩阵 M M M over blockSize x blockSize neighborhood, 然后计算: d s t ( x , y ) = d e t ( M ) − k ( T r a c e ...
  • Harris角点检测(python-opencv) 1. 图像特征 图像的特征我认为可以分为基于一幅图像的平面,边缘和角点。 2.Harris角点检测 图像是N*M维度的矩阵,选择一个合适的窗口滑过图像的每个区域,然后计算所有差值的总和...
  • 针对传统Harris角点检测算法和目前一些改进算法应用在图像拼接时,仍然可能存在只可在单一尺度上检测角点位置不准确、伪检和对噪声敏感致使检测率不高等缺点,提出一种基于AP聚类角点提取优化的双边滤波(BF)角点...
  • Harris 角点检测算法是由 ​C.Harris 和 J.Stephens 在 1988 年提出,此算法的原理是建立在 ​Moravec 算法的基础上进行大幅改进的,相比之前的算法更具明显的优势。​Moravec如上文中的介绍说明,相比 Harris 检测...
  • Harris角点检测python实现及基于opencv实现,(角点检测、非极大值抑制)
  • 基于曲率尺度空间的角点检测图技术,优化设计图像匹配算法,基于曲率尺度空间的角点检测算法进行图像特征的提取,归一化处理特征,有助于提高图像匹配精度。利用该算法最终实现图像匹配需求,验证了算法的有效性...
  • 目录 @[TOC] # 一、Harris角点检测算法 ## 1.1角点是什么
  • Harris角点检测原理分析

    千次阅读 2019-04-01 16:22:14
    转自:Harris角点检测原理分析 Harris角点检测算子是于1988年由CHris Harris & Mike Stephens提出来的。在具体展开之前,不得不提一下Moravec早在1981就提出来的Moravec角点检测算子。 1.Moravec角点检测算子 ...
  • 局部图像描述子——图像特征匹配特征 特征 在这里插入代码片
  • 一、Harris角点检测算法 1.1 角点是什么 角点具有的特征: <1>轮廓之间的交点; <2>局部窗口沿任意方向移动,均产生明显变化的; <3>图像局部曲率突变的; <4>对于同一场景,即使...
  • 在摄像机自标定过程中,可根据Harris检测算法提取对角点。该算法简单有效,非常稳定。在图像旋转、灰度、噪声影响和视点变换的条件下,与其他算子相比,是最稳定的一种特征提取算子。为了获得亚像素级的角点坐标...
  • 1. Harris角点检测的原理 2. Harris角点检测的详细源码分析

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 2,006
精华内容 802
关键字:

harris角点检测优化算法