东莞市盛裕绒艺玩具有限公司

东莞市盛裕绒艺玩具有限公司

乐动官方客服

15830799767
联系方式
全国服务热线: 15830799767

咨询热线:15858372790
联系人:王辉
地址:泰州市海陵区工业园

亚像素数值极值检测算法总结

来源:乐动官方客服   发布时间:2019-07-02   点击量:263

动机

在计算机视觉领域,经常需要检测极值位置,比如SIFT关键点检测、模板匹配获得最大响应位置、统计直方图峰值位置、边缘检测等等,有时只需要像素精度就可以,有时则需要亚像素精度。本文尝试总结几种常用的一维离散数据极值检测方法,几个算法主要来自论文《A Comparison of Algorithms for Subpixel Peak Detection》,加上自己的理解和推导。

问题定义

给定如下离散值,求其极值位置。可知125为观察极值。

[[60, 80, 100, 120, 125, 105, 70, 55]]

如果这些离散值是从某个分布(f)中等间距采样获得,其真正的极值位置应位于120和125之间。

下面给出形式化的定义:给定一组离散值,令(x)为观测到的极值点位置,其值为(f(x)),其左右相邻位置的值为(f(x-1))(f(x+1)),真正的极值点位置为(x+delta),令(hat{delta})(delta)的估计值。

算法

假设(x)的邻域可通过某个模型进行近似,如高斯近似、抛物线近似,则可以利用(x)的邻域信息根据模型估计出极值。使用的模型不同就有不同的算法,具体如下。

高斯近似

一维高斯函数如下:

[y = y_{max} cdot exp(-frac{(x-mu)^2}{2sigma^2})](y_{max}=frac{1}{sqrt{2sigma}pi})时为标准高斯函数,形如

假设(x)的邻域可用高斯近似,用((x, f(x)))((x-1, f(x-1)))((x+1, f(x+1)))三点对高斯函数进行拟合,获得模型参数(mu)即为峰值位置,(hat{delta}=mu - x)。将三点带入上面的高斯函数两边同时取对数求得:

[hat{delta} = frac{1}{2} frac{ln(f(x-1)) - ln(f(x+1))}{ln(f(x-1)) - 2ln(f(x)) + ln(f(x+1))}]

下面可以看到,高斯近似相当于取对数后的抛物线近似

抛物线近似

使用抛物线近似(x)的局部,可以将((x, f(x)))((x-1, f(x-1)))((x+1, f(x+1)))三点带入(y=a(x-b)^2+c)求参数(b)即为估计的极值位置,也可采用泰勒展开(牛顿法)来求极值。泰勒公式实际上是一种利用高阶导数通过多项式近似函数的方法,下面的图示可直观理解这种近似,图示为通过泰勒公式近似原点附近的正弦曲线:

泰勒近似(x)附近,如只取到二阶则为抛物线近似。假设高阶可导,极值为(f(x+delta)),则根据泰勒公式,

[f(x+delta) = f(x) + f"(x)delta + frac{1}{2} f""(x)delta^2 + O(delta^3)]

极值处导数为0,这里(x)为常数(delta)为变量,两边同时对(delta)求导,忽略高阶项可得

[f"(x+hat{delta}) = f"(x) + f""(x)hat{delta} = 0]

使用一阶微分和二阶微分近似(f"(x))(f""(x))

[hat{delta} = - frac{f"(x)}{f""(x)} = - frac{(f(x+1)-f(x-1))/2}{(f(x+1)-f(x))-(f(x) - f(x-1))}= frac{1}{2}frac{f(x-1)-f(x+1)}{f(x+1)-2f(x)+ f(x-1)}]与带入抛物线求参数的结果是一致的,加上对数则与高斯近似一致。

质心算法

In physics, the center of mass of a distribution of mass in space is the unique point where the weighted relative position of the distributed mass sums to zero, or the point where if a force is applied it moves in the direction of the force without rotating.——Center of mass wiki

若将(x)(x-1)(x+1)看成质点,将(f(x))(f(x-1))(f(x+1))看成质点的质量,则可以把质心作为极值的估计。根据质点相对质心位置的质量加权和为零,可求得质心位置。令(R)为质心坐标,(m)(r)分别为质点质量和坐标,则(n)个质点的质心满足

[sum_{i=1}^n m_i(r_i - R) = 0]

(M = sum_{i=1}^n m_i),质心坐标为

[R = frac{1}{M} sum_{i=1}^n m_ir_i]

带入得

[x + hat{delta} = frac{(x-1)f(x-1)+xf(x)+(x+1)f(x+1)}{f(x-1)+f(x)+f(x+1)}]

[hat{delta} = frac{f(x+1)-f(x-1)}{f(x-1)+f(x)+f(x+1)}]

以上考虑的是3质点系统的质心,还可考虑5质点、7质点等,甚至考虑所有点。

线性插值

这个模型假设在极值两侧是线性增长和线性下降的,且上升和下降的速度相同,即(y=kx+b),上升侧(k>0),下降侧(k<0),两者绝对值相同,可以利用这个性质求解极值位置。

(f(x+1)>f(x-1))则极值位于((x, x+1))之间,可列等式

[frac{f(x) - f(x-1)}{x-(x-1)} = frac{f(x+delta)-f(x)}{x+delta - x} = frac{f(x+delta)-f(x+1)}{x+1-(x+delta)}]

解得

[hat{delta}=frac{1}{2}frac{f(x+1)-f(x-1)}{f(x)-f(x-1)}]

同理,若(f(x-1)>f(x+1))求得

[hat{delta}=frac{1}{2}frac{f(x+1)-f(x-1)}{f(x)-f(x+1)}]

数值微分滤波

这个方法是利用极值处导数为0的性质,在微分滤波结果上插值得到导数为0的位置,因已知极值点在(x)附近,因此只需在(x)附近做微分和插值即可。插值时取极值点两侧正负值连线的过零点作为极值点的估计,如下图所示

论文Real-time numerical peak detector中定义了4阶和8阶线性滤波器([1, 1, 0, -1, -1])([1,1,1,1,0,-1,-1,-1,-1]),对应的函数形式为

[g_4(x)=f(x-2)+f(x-1)-f(x+1)-f(x+2)]

[g_8(x)=f(x-4)+f(x-3)+f(x-2)+f(x-1) \-f(x+1)-f(x+2)-f(x+3)-f(x+4)]

2阶形式为(g_2(x) = f(x-1) -f(x+1)),这些滤波器的表现与数值微分滤波器相似。

(f(x+1)>f(x-1))时,极值点位于((x, x+1))之间,(g(x)<0)(g(x+1)>0),极值点位置为(g(x))(g(x+1))连线的过零点,通过斜率求得

[hat{delta} = frac{g(x)}{g(x+1)-g(x)}]

(f(x-1)>f(x+1)),则

[hat{delta} = frac{g(x-1)}{g(x-1)-g(x)} - 1]

总结

这些数值极值检测方法均是先获取观测极值(x)及其邻域信息,然后综合邻域信息在各自的模型假设下通过插值估计出极值位置。若能知道数值来自的真实分布,则直接拟合真实分布然后求极值即可,但往往我们并不知道真实的分布是什么,即使知道真实分布,有时为了快速计算,也会采取插值的方式来估计极值,毕竟偏差可接受效果足够好就可以了。应用时,为了抗噪可对数据先平滑然后求极值,具体采用何种方法可在准确和速度间权衡——所用模型与真实分布越相近自然越准确,如果实在不知道怎么选,就实践对比吧(因为我也不知道),毕竟伟大领袖教导过我们——实践是检验真理的唯一标准

参考

A Comparison of Algorithms for Subpixel Peak DetectionReal-time numerical peak detector

个人博客地址:亚像素数值极值检测算法总结

, 1, 0, 9);

相关产品

COPYRIGHTS©2017 乐动官方客服 ALL RIGHTS RESERVED 备案号:263