Harris角点检测算法概述

Harris角点检测算法概述

角点特征的提取一般是取某点为中心的窗口计算响应值,然后筛选出响应值显著的点为角点。一般认为,如果一个点是角点,则以角点为中心的窗口向任何方向滑动时,窗口内灰度分布变化均应该较大;如果一个点是边缘点,则滑动窗口时在边缘走向上灰度分布基本不变而垂直边缘的方向上变化剧烈。如图所示,经典的角点提取算子有Moravec算子和Harris算子。Moravec算子分别统计当前点所在窗口向水平、竖直、对角线和反对角线四个方向滑动一个像素时的灰度变化平方和,也就是统计四个方向的灰度梯度平方和,取这四个值中的最小值作为该点的兴趣值:若这一指标大于某阈值,则说明在四个方向上灰度都变化大,该点很可能为角点,经过非极大值抑制后保留真正的角点。这一算法计算简便,但是容易受噪声影响。

Moravec算子计算的四个方向

Harris算子与Moravec算子出发点类似,讨论了角点和非角点的窗口灰度分布特点。而Harris算子考虑了噪声影响,在计算中引入了使用高斯卷积核(下式),且更进一步导出了新的判断测度。

w u , v = e − u 2 + v 2 2 σ 2 w_{u,v}=e^{-\frac{u^2+v^2}{2\sigma^2}} wu,v=e2σ2u2+v2

如下式所示, I ( x , y ) I(x,y) I(x,y)表示图像灰度函数(连续的二元函数),窗口向任意方向滑动 ( d x , d y ) (dx,dy) (dx,dy)时,窗口内的灰度变化情况为加权的差分平方和;窗口内的像素点 ( u , v ) (u,v) (u,v)处作一阶泰勒展开后舍去二次小项,并进一步整理得到一个二次型,其的性质由矩阵 M M M反映。

E d x , d y = ∑ ( u , v ) ∈ N b w u , v [ I d x + u , d y + v − I u , v ] 2 = ∑ ( u , v ) ∈ N b w u , v [ ∂ I ∂ x d x + ∂ I ∂ y d y + o ( d x 2 , d y 2 ) ] 2 = ∑ ( u , v ) ∈ N b ( w u , v ∂ I ∂ x 2 ) ( d x ) 2 + 2 ∑ ( u , v ) ∈ N b w u , v ∂ I ∂ x ∂ I ∂ y d x d y + ∑ ( u , v ) ∈ N b ( w u , v ∂ I ∂ y 2 ) ( d y ) 2 = A d x 2 + B d y 2 + 2 C d x d y = [ d x d y ] M [ d x d y ] \begin{aligned} E_{dx,dy} &= \sum_{(u,v)\in{N_{b}}}{w_{u,v}[I_{dx+u,dy+v}-I_{u,v}]^2}\\ &= \sum_{(u,v)\in{N_{b}}}{w_{u,v}[\frac{\partial I}{\partial x}dx+\frac{\partial I}{\partial y}dy+o(dx^2,dy^2)]^2}\\ &= \sum_{(u,v)\in{N_{b}}}{(w_{u,v}\frac{\partial I}{\partial x}^2)}(dx)^2+2\sum_{(u,v)\in{N_{b}}}{w_{u,v}\frac{\partial I}{\partial x}\frac{\partial I}{\partial y}}dxdy+\sum_{(u,v)\in{N_{b}}}{(w_{u,v}\frac{\partial I}{\partial y}^2)}(dy)^2\\ &= Adx^2+Bdy^2+2Cdxdy\\ &= \begin{bmatrix}dx&dy\end{bmatrix} M \begin{bmatrix}dx\\dy\end{bmatrix} \end{aligned} Edx,dy=(u,v)Nbwu,v[Idx+u,dy+vIu,v]2=(u,v)Nbwu,v[xIdx+yIdy+o(dx2,dy2)]2=(u,v)Nb(wu,vxI2)(dx)2+2(u,v)Nbwu,vxIyIdxdy+(u,v)Nb(wu,vyI2)(dy)2=Adx2+Bdy2+2Cdxdy=[dxdy]M[dxdy]

其中 A = ∑ ( u , v ) ∈ N b ( w u , v ∂ I ∂ x 2 ) A=\sum_{(u,v)\in{N_{b}}}{(w_{u,v}\frac{\partial I}{\partial x}^2)} A=(u,v)Nb(wu,vxI2) B = ∑ ( u , v ) ∈ N b ( w u , v ∂ I ∂ y 2 ) B=\sum_{(u,v)\in{N_{b}}}{(w_{u,v}\frac{\partial I}{\partial y}^2)} B=(u,v)Nb(wu,vyI2) C = ∑ ( u , v ) ∈ N b w u , v ∂ I ∂ x ∂ I ∂ y C=\sum_{(u,v)\in{N_{b}}}{w_{u,v}\frac{\partial I}{\partial x}\frac{\partial I}{\partial y}} C=(u,v)Nbwu,vxIyI M = [ A C C B ] M=\begin{bmatrix}A&C\\C&B\end{bmatrix} M=[ACCB]

M M M的具体计算即下式所示:

M = ∑ ( u , v ) ∈ N b w u , v [ I x 2 I x I y I x I y I y 2 ] M=\sum_{(u,v)\in{N_{b}}}w_{u,v} \begin{bmatrix} I_{x}^{2} & I_{x}I_{y}\\ I_{x}I_{y} & I_{y}^{2} \end{bmatrix} M=(u,v)Nbwu,v[Ix2IxIyIxIyIy2]

M M M可以特征分解为以 M = A T P A M=A^{T}PA M=ATPA的形式, A A A是正交矩阵; P P P是对角阵,主对角线元素为特征值,所以上面式子可以分解为这一形式:

E x , y ≈ [ d x ′ d y ′ ] [ λ 1 0 0 λ 2 ] [ d x ′ d y ′ ] E_{x,y}\approx \begin{bmatrix}dx'&dy'\end{bmatrix} \begin{bmatrix} \lambda_{1} & 0\\ 0 & \lambda_{2} \end{bmatrix} \begin{bmatrix}dx'\\ dy'\end{bmatrix} Ex,y[dxdy][λ100λ2][dxdy]

[ d x ′ d y ′ ] = A [ d x d y ] \begin{bmatrix}dx'\\ dy'\end{bmatrix}=A\begin{bmatrix}dx\\ dy\end{bmatrix} [dxdy]=A[dxdy]
是正交变换后的窗口平移滑动矢量。如果两个特征值 λ 1 , λ 2 \lambda_1,\lambda_2 λ1,λ2都比较大,则说明窗口向任何方向滑动都导致较大的灰度变化,当前点很可能为角点;如果只有其中一个比较大,则说明可能是边缘点。因此,将两个特征值的大小情况作为角点判据是合理的。
实际上借助式\eqref{eq:Rindex}构造了指标 R = d e t ( M ) − k ⋅ [ t r ( M ) ] 2 R=det(M)-k·[tr(M)]^{2} R=det(M)k[tr(M)]2避免进行特征分解的不必要计算。 R R R称为响应值函数,它综合反映了两个特征值的情况;R为正值时,值越大,越可能是角点。R的取值与两特征值的情况见图~\ref{fig:Rindex}

d e t ( M ) = λ 1 ∗ λ 2 t r ( M ) = λ 1 + λ 2 \begin{aligned} det(M)&=\lambda_{1}*\lambda_{2}\\ tr(M)&=\lambda_{1}+\lambda_{2} \end{aligned} det(M)tr(M)=λ1λ2=λ1+λ2
R函数等值线和两特征值的关系示意图

算法的具体步骤如下:

  • 计算窗口内每个像素的梯度;
  • 按照\eqref{eq:Mexpansion}计算梯度自相关矩阵的带权和矩阵 M M M
  • 从矩阵 M M M计算性能指标R,如果大于给定阈值则保留;
  • 全部像素处理后进行非极大值抑制。
  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值