CN104463876B - 一种基于自适应滤波的复杂背景下图像多圆快速检测方法 - Google Patents

Abstract

本发明公开了一种基于自适应滤波的复杂背景下图像多圆快速检测方法,属于计算机视觉在线检测技术领域。应用Sage‑Husa自适应卡尔曼滤波算法对各圆特征的中心坐标及半径进行估计，进而对特征有效区域进行估计；在各圆特征有效区域内，图像局部增强后采用Canny边缘检测算法及最小二乘椭圆拟合法对圆特征的中心坐标及半径进行求解；应用结果有效判据对定位结果是否有效进行判定；若判定结果有效，则将特征提取值作为检测结果输出；若判定结果无效，则将特征参数估计值作为检测结果输出。本发明在复杂背景多圆定位过程中彻底抛弃了非常耗时的Hough变换，输出测量结果前先对定位结果是否有效进行判定，并分别采取不同的对策，在保证精度的同时实现了圆检测速度大幅提升，能满足采用高分辨率图像传感器的工业监测系统在线监测的要求。

Description

一种基于自适应滤波的复杂背景下图像多圆快速检测方法

技术领域

本发明涉及计算机视觉在线检测技术领域,特别是一种基于自适应滤波的复杂背景下多圆快速检测方法。

背景技术

计算机视觉技术具有非接触性、经济性、灵活性和集成性等优点，在工业测试与在线检测领域具有广泛的应用前景。圆检测是许多计算机视觉测量系统的关键步骤之一，近年出现的一些测试系统中，要求计算机视觉算法在高分辨率复杂背景中实现快速精确多圆提取。

传统的圆检测方法包括模版法、重心法、矩估计法、曲线拟合法和Hough变换法等，其中，只有Hough变换法能在复杂背景下实现多圆定位，但实时性差、精度低。随后出现的随机Hough变换、随机圆检测等算法可在一定程度上提高检测速度，但实时性和精度无法同时达到较高的水平，不能满足工业在线检测的要求。

近年，国内外学者对复杂背景下多圆定位进行了深入的研究，提出了许多新的方法，多圆定位的精度和实时性得到了一定程度的提高。现有方法已能够满足采用低分辨率图像传感器的工业监测系统在线监测的要求。但是上述改进仍以速度很慢的Hough变换为基础，随着硬件技术的发展，工业检测图像的分辨率越来越高，这些方法用于高分辨率检测图像时耗时将迅速增加，实时性差仍然是一个突出的问题。

发明内容

本发明所要解决的技术问题是，针对上述现有技术的不足，提供一种基于自适应滤波的复杂背景下图像多圆快速检测方法，对于背景复杂、分辨率高的检测图像快速稳定地实现多圆检测。

为解决上述技术问题，本发明所采用的技术方案是：一种基于自适应滤波的复杂背景下图像多圆快速检测方法，包括以下步骤：

1)利用Sage-Husa自适应卡尔曼滤波方法对各圆特征的中心坐标及半径进行估计，得到特征参数估计值，进而对特征有效区域进行估计；

2)在各圆的特征有效区域内，图像局部增强后采用Canny边缘检测方法及最小二乘椭圆拟合法对各圆特征的中心坐标及半径进行求解，得到初次定位结果；

3)应用优化的特征定位结果有效判据对上述初次定位结果是否有效进行判定；

4)若步骤3)判定结果有效，则将步骤2)的初次定位结果作为检测结果输出；若步骤3)判定结果无效，则将步骤1)的特征参数估计值作为检测结果输出。

所述步骤1)中，特征有效区域为以下列四个点为顶点的正方形：

其中，分别为中心坐标估计值的横坐标和纵坐标；为半径估计值；8≤w≤20。

所述步骤2)中，采用形态学开闭双边滤波方法进行图像局部增强处理。

所述步骤2)中，最小二乘椭圆拟合法对各圆特征的中心坐标求解采用的目标函数为：

其中，n为Canny边缘检测后的离散边缘点的个数，x_i、y_i为第i个离散边缘点的x和y坐标，约束条件设为A+C＝1，未知参数求解准则为：

圆的中心坐标(x_c,y_c)为：

所述步骤2)中，求解的半径为：

其中，

所述步骤3)中，优化的特征定位结果有效判据为：

[Av(k)]^T[Av(k)]≤γtr{[H(k)p(k|k-1)H^T(k)+R(k)]A²}；

当上式成立时，判定结果有效，当上式不成立时，判定结果无效；

其中，v(k)为k时刻圆特征定位结果后计算的新息向量，Z(k)＝[x(k) y(k) r(k)]^T，x(k)、y(k)、r(k)分别为k时刻圆特征的中心横坐标、纵坐标和半径，Φ(k,k-1)为由k-1时刻到k时刻的圆特征状态转移矩阵，T为k-1时刻到k时刻检测的时间间隔，为k-1时刻圆特征状态参数；p(k|k-1)为一步预测均方差矩阵，p(k|k-1)＝Φ(k,k-1)p(k-1)Φ^T(k,k-1)+Γ(k,k-1)Q(k-1)Γ^T(k,k-1)，p(k-1)为k-1时刻均方差矩阵，Γ(k,k-1)为过程噪声驱动矩阵，Γ(k,k-1)为单位阵；Q(k-1)为k-1时刻过程噪声协方差矩阵；R(k)为观测噪声协方差矩阵，R(k)＝(1-d_k)R(k-1)+d_k{[1-H(k)K(k-1)]v(k)v^T(k)[1-H(k)K(k-1)]^T+H(k)p(k-1)H^T(k)}，d_k＝(1-b)/(1-b^k+1),b为遗忘因子，0.972≤b≤0.993；K(k-1)为k-1时刻滤波增益矩阵，K(k-1)＝p(k-1|k-2)H^T(k-1)[H(k-1)p(k-1|k-2)H^T(k-1)+R(k-1)]^-1；权重矩阵a1、a2取值为1，2≤a3≤5；2≤γ≤10；观测矩阵

与现有技术相比，本发明所具有的有益效果为：本发明将Sage-Husa滤波算法和曲线拟合法结合，在复杂背景多圆定位过程中彻底抛弃了非常耗时的Hough变换，从而实现圆检测速度大幅提升；针对圆检测的特点，引入权重矩阵，给出了优化的特征定位结果有效性判据，使参数估计过程中滤波器能更好的对自身状态做出判断；当判定圆特征提取结果无效时，本发明直接将新的测量值排除在滤波记忆及检测结果以外，而将基于Sage-Husa滤波算法和历史测量数据得到的圆特征状态估计值作为检测结果输出，在工作环境恶劣的工业现场，可有效的抑制测量野值对检测系统的影响，提升系统的鲁棒性；本发明对于分辨率为1280像素x 720像素的高清测量图像，圆特征中心定位rms误差为2.8×10^-6mm,平均耗时小于215.3ms，本发明精度较高，实时性好，鲁棒性强，能满足采用高分辨率图像传感器的工业监测系统在线监测的要求。本发明还可推广到对速度要求高的其他计算机视觉在线检测系统。

附图说明

图1为本发明一实施例方法流程图；图中，1、圆特征区域估计，2、特征精确提取；

图2为本发明一实施例依据某圆特征的中心坐标及半径估计结果，设定特征有效估计区域过程示意图；

图3(a)为本发明一实施例多圆特征第5帧估计结果；图3(b)为本发明一实施例多圆特征第10帧估计结果；图3(c)为本发明一实施例多圆特征第15帧估计结果；图3(d)为本发明一实施例多圆特征第20帧估计结果；图3(e)为本发明一实施例多圆特征第25帧估计结果；图3(f)为本发明一实施例多圆特征第30帧估计结果；图3(g)为本发明一实施例多圆特征第35帧估计结果；图3(h)为本发明一实施例多圆特征第40帧估计结果；

图4为本发明一实施例中根据本发明优化的特征定位结果有效判据计算得到的γ与标准判据计算得到的γ值对比曲线；

图5为本发明实施例与现有技术耗时对比曲线。

具体实施方式

如图1所示，本发明一实施例方法步骤为：

1)圆特征区域估计：应用Sage-Husa自适应卡尔曼滤波算法对各圆特征的中心坐标及半径进行估计，进而对特征有效区域进行估计；

2)特征精确提取：在各圆特征的有效区域内，图像局部增强后采用Canny边缘检测算法及最小二乘椭圆拟合法对圆特征的中心坐标及半径进行求解；

3)提取结果有效性判定：在得到提取结果后，应用一种优化的特征定位结果有效判据对初次定位结果是否有效进行判定；

4)输出检测结果：若步骤3判定结果有效，则将步骤2的特征提取值作为检测结果输出；若步骤3判定结果无效，则将步骤1的特征参数估计值作为检测结果输出。

所述步骤1)中，圆特征参数估计是基于Sage-Husa自适应卡尔曼滤波对圆特征的参数进行估计实现的。由卡尔曼原理可知，线性离散系统状态方程为：

X(k)＝Φ(k,k-1)X(k-1)+Γ(k,k-1)W(k-1)

其中，X(k)为系统的状态向量，Φ(k,k-1)为由k-1时刻到k时刻的系统状态转移矩阵，Γ(k,k-1)过程噪声驱动矩阵，W(k)为过程噪声向量。

选取中心坐标(x、y)、中心坐标对时间的一阶导数和半径r为状态变量，则k时刻滤波器的状态向量X(k)为：

依据卡尔曼滤波原理，以平面内的匀加速直线运动来近似特征的运动规律，滤波器k-1时刻到k时刻的状态转移矩阵Φ(k,k-1)为：

其中，T为两次检测的时间间隔。圆特征提取检测中，k时刻过程噪声与k-1刻过程噪声无确定的转移关系，过程噪声驱动矩阵Γ(k,k-1)为单位阵。

圆特征提取检测中，直接观测量为中心坐标(x、y)和半径r，则k时刻滤波器的观测向量Z(k)为：

Z(k)＝[x(k) y(k) r(k)]^T

滤波器的观测矩阵H(k)为：

将圆特征运动系统视为线性高斯系统，过程噪声W(k)和观测噪声V(k)可处理成具有时变均值和时变协方差阵的不相关高斯随机过程,即：

E(W(k))＝q(k)

E{[W(k)-q(k)][W(j)-q(j)]^T}＝Q(k)δ_kj

E(V(k))＝r(k)

E{[V(k)-r(k)][V(j)-r(j)]^T}＝R(k)δ_kj

其中，δ_kj为克罗内克函数，假设过程噪声和观测噪声均是无偏的,即：q(k)＝r(k)＝0。过程和观测协方差矩阵Q(k)和R(k)阵需要依据离线实验误差方差的历史统计数据进行初始化。若历史数据估计圆中心坐标(x、y)、中心坐标对时间的一阶导数和半径r的系统模型误差方差分别为D(x)、D(y)、D(r)，则Q(k)初始化为：

若历史数据估计圆中心坐标(x、y)和半径r的观测误差方差分别为D₁(x)、D₁(y)、D₁(r)，则R(k)初始化为：

依据Sage-Husa自适应滤波算法和历史检测数据对圆的中心坐标和半径进行迭代估计。迭代估计过程为：

由圆特征k-1时刻的状态预测圆特征k时刻状态：

得到k时刻圆特征定位结果后计算新息向量v(k)：

更新滤波增益矩阵K(k)：

K(k)＝p(k|k-1)H^T(k)[H(k)p(k|k-1)H^T(k)+R(k)]^-1

更新一步预测均方差矩阵p(k|k-1)：

p(k|k-1)＝Φ(k,k-1)p(k-1)Φ^T(k,k-1)+Γ(k,k-1)Q(k-1)Γ^T(k,k-1)

更新估计均方差矩阵p(k)：

p(k)＝[1-K(k)H(k)]p(k|k-1)[1-K(k)H(k)]^T+K(k)R(k-1)K^T(k)

修正k时刻圆特征状态参数：

更新观测噪声协方差矩阵R(k)：

R(k)＝(1-d_k)R(k-1)+d_k{[1-H(k)K(k-1)]v(k)v^T(k)[1-H(k)K(k-1)]^T+H(k)p(k-1)H^T(k)}

其中，d_k＝(1-b)/(1-b^k+1),b为遗忘因子，0.972≤b≤0.993。

迭代估计开始前，需获取前两个位置圆特征的状态参数对滤波器的状态向量初始化，该过状态参数可通过Hough圆变换获取。

若由两个位置圆特征状态参数构建的观测向量分别为Z'(1)＝{x'(1)y'(1)r'(1)}和Z'(2)＝{x'(2)y'(2)r'(2)}，则滤波器状态向量初始化为：

滤波器观测向量初始化为：

Z(0)＝Z'(2)

所述步骤1)中，若依据k-1时刻的观测对k时刻某一特征的基本参数估计结果为：特征中心坐标半径为确保圆特征全部在有效区域内，设置余量参数w，w取值范围为8≤w≤20，则特征估计区域为以下四个点为顶点的正方形：

所述步骤2)中，作为优选方案，图像局部增强的方法为形态学开闭双边滤波。

所述步骤2)中，最小二乘椭圆拟合法对圆特征的中心坐标及半径进行求解的过程为：应用椭圆的标准方程对有效区域内边缘检测后的离散点进行最小二乘处理，可得目标函数：

其中，n为边缘检测后的离散边缘点的个数，x_i、y_i为第i个离散点的x和y坐标。为了避免零解,并将解的任何整数倍都视为对同一椭圆的表述,对参数做一些限制,约束条件设为：

A+C＝1

各系数按照目标函数取最小值准则确定，由极值原理可知,欲使函数f值为最小,必有：

由上式可得到一个线性方程组，结合约束条件可求解得到方程各系数的值，椭圆中心(x_c,y_c)坐标为：

椭圆的长短半轴a、b分别为：

输出的圆半径为得到的椭圆长短半轴的平均值：

所述步骤3)中，以卡尔曼滤波发散判据为基础，为更好的符合实际情况，对状态向量中不同元素赋予不同权重，定义权重矩阵A：

其中，a1、a2、a3分别为圆特征中心x坐标、y坐标及半径在发散判定过程中所占的权重。a1、a2取为1，a3取值范围限定为2≤a3≤5，所述优化的特征定位结果有效判据为：

[Av(k)]^T[Av(k)]≤γtr{[H(k)p(k|k-1)H^T(k)+R(k)]A²}

其中，γ为储备系数，γ取值范围限定为2≤γ≤10。

本发明实施例检测图像为分辨率1280像素x 720像素的高清图像，如图1所示，检测过程如下：

由于实施例中检测图像包含4个圆特征，先建立4个滤波器对图像中4个圆特征分别进行参数估计，过程噪声和观测噪声的协方差矩阵依据先验统计数据分别初始化为：

本实施例中，为遗忘因子b取值为0.99。依据Hough圆变换获取本发明方法实施以前两个位置圆特征的状态参数，对滤波器的状态向量初始化。依据Sage-Husa自适应滤波算法和历史检测数据对圆的中心坐标和半径进行迭代估计。

得到圆的中心坐标和半径估计值后，本发明通过确定圆特征区域的4个顶点设定圆特征估计区域，过程参见附图2，本实施例中，余量常数w取值为7.5像素，实施例执行过程中特征区域的预测结果参见附图3，可见，本发明方法对所有位置的特征区域均做出了有效预测，某些时刻，滤波器对特征参数的预测出现了明显的偏差(如第20和第25帧)，导致特征没有出现在预测区域的中心，但由于余量常数w取值恰当，特征并没有出现在预测区域之外，且由于Sage-Husa滤波方法能够对观测噪声进行在线自整定，后续时刻特征又逐渐回到了预测区域的中心位置。

在各圆特征的有效估计区域内采用形态学开闭双边滤波方法对检测图像进行局部增强，并通过Canny算法实现局部图像的边缘二值化，最后，最小二乘椭圆拟合法对圆特征的中心坐标及半径进行求解。

得到新的圆特征的中心坐标及半径计算值后，依据本发明优化的特征定位结果有效判据对新的计算值是否有效进行判定，判定过程中，权重矩阵中参数a1取值1、a2取值1、a3取值3，贮备系数γ取值为3.5。实施例中，位于图像左上方的圆特征依据标准滤波收敛判据和本发明优化的有效判据计算得到的贮备系数γ值对比参见附图4，采用本发明的有效判据,在定位检测结果有效时γ值的波动范围更小，而当出现无效检测结果时(第17帧检测图像)，本发明判据计算得到的γ值增加更显著，本发明判据对于圆检测具有更好的适用性。

最后，输出检测结果，本发明的圆检测精度与发明《噪声环境下激光束中心高效精确检测方法》中公开的方法精度完全一致，耗时对比曲线参见附图5，对于高清测量图像，本发明方法定位耗时在110ms上下波动，而现有技术特征定位平均耗时为632ms。本发明用Sage-Husa自适应滤波算法替代耗时的Hough变换实现对特征区域的有效定位，特征提取速度有显著的提升。

本发明将Sage-Husa滤波算法和曲线拟合法结合，在复杂背景多圆定位过程中彻底抛弃了非常耗时的Hough变换，从而实现圆检测速度大幅提升；针对圆检测的特点，引入权重矩阵，给出了优化的特征定位结果有效性判据，使参数估计过程中滤波器能更好的对自身状态做出判断；当判定圆特征提取结果无效时，本发明直接将新的测量值排除在滤波记忆及检测结果以外，而将基于Sage-Husa滤波算法和历史测量数据得到的圆特征状态估计值作为检测结果输出，在工作环境恶劣的工业现场，可有效的抑制测量野值对检测系统的影响，提升系统的鲁棒性；本发明对于分辨率为1280像素x 720像素的高清测量图像，圆特征中心定位rms误差为2.8×10^-6mm,平均耗时小于215.3ms，本发明精度较高，实时性好，鲁棒性强，能满足采用高分辨率图像传感器的工业监测系统在线监测的要求。本发明还可推广到对速度要求高的其他计算机视觉在线检测系统。

Claims (5)

1.一种基于自适应滤波的复杂背景下图像多圆快速检测方法，其特征在于，包括以下步骤：

1)利用Sage-Husa自适应卡尔曼滤波方法对各圆特征的中心坐标及半径进行估计，得到特征参数估计值，进而对特征有效区域进行估计；

2)在各圆的特征有效区域内，图像局部增强后采用Canny边缘检测方法及最小二乘椭圆拟合法对各圆特征的中心坐标及半径进行求解，得到初次定位结果；

3)应用优化的特征定位结果有效判据对上述初次定位结果是否有效进行判定；优化的特征定位结果有效判据为：[Av(k)]^T[Av(k)]≤γtr{[H(k)p(k|k-1)H^T(k)+R(k)]A²}；

当上式成立时，判定结果有效，当上式不成立时，判定结果无效；其中，v(k)为k时刻圆特征定位结果后计算的新息向量，Z(k)＝[x(k) y(k) r(k)]^T，x(k)、y(k)、r(k)分别为k时刻圆特征的中心横坐标、纵坐标和半径，Φ(k,k-1)为由k-1时刻到k时刻的圆特征状态转移矩阵，T为k-1时刻到k时刻检测的时间间隔，为k-1时刻圆特征状态参数；p(k|k-1)为一步预测均方差矩阵，p(k|k-1)＝Φ(k,k-1)p(k-1)Φ^T(k,k-1)+Γ(k,k-1)Q(k-1)Γ^T(k,k-1)，p(k-1)为k-1时刻均方差矩阵，Γ(k,k-1)为过程噪声驱动矩阵，Γ(k,k-1)为单位阵；Q(k-1)为k-1时刻过程噪声协方差矩阵；R(k)为观测噪声协方差矩阵，R(k)＝(1-d_k)R(k-1)+d_k{[1-H(k)K(k-1)]v(k)v^T(k)[1-H(k)K(k-1)]^T+H(k)p(k-1)H^T(k)}，d_k＝(1-b)/(1-b^k+1),b为遗忘因子，0.972≤b≤0.993；K(k-1)为k-1时刻滤波增益矩阵，K(k-1)＝p(k-1|k-2)H^T(k-1)[H(k-1)p(k-1|k-2)H^T(k-1)+R(k-1)]^-1；权重矩阵a1、a2取值为1，2≤a3≤5；2≤γ≤10；观测矩阵

4)若步骤3)判定结果有效，则将步骤2)的初次定位结果作为检测结果输出；若步骤3)判定结果无效，则将步骤1)的特征参数估计值作为检测结果输出。

2.根据权利要求1所述的基于自适应滤波的复杂背景下图像多圆快速检测方法，其特征在于，所述步骤1)中，特征有效区域为以下列四个点为顶点的正方形：

$S\left(\hat{x}\right(k|k-1)-\hat{r}(k|k-1)-w,\hat{y}(k|k-1)+\hat{r}(k|k-1)+w)$

$L\left(\hat{x}\right(k|k-1)+\hat{r}(k|k-1)+w,\hat{y}(k|k-1)+\hat{r}(k|k-1)+w)$

$M\left(\hat{x}\right(k|k-1)+\hat{r}(k|k-1)+w,\hat{y}(k|k-1)-\hat{r}(k|k-1)-w)$

$N\left(\hat{x}\right(k|k-1)-\hat{r}(k|k-1)-w,\hat{y}(k|k-1)-\hat{r}(k|k-1)-w)$

其中，分别为中心坐标估计值的横坐标和纵坐标；为半径估计值；8≤w≤20。

3.根据权利要求1所述的基于自适应滤波的复杂背景下图像多圆快速检测方法，其特征在于，所述步骤2)中，采用形态学开闭双边滤波方法进行图像局部增强处理。

4.根据权利要求2或3所述的基于自适应滤波的复杂背景下图像多圆快速检测方法，其特征在于，所述步骤2)中，最小二乘椭圆拟合法对各圆特征的中心坐标求解采用的目标函数为：

$f(A,B,C,D,E,F)=\underset{i=1}{\overset{n}{\Σ}}{({\mathrm{Ax}}_i^2+{\mathrm{Bx}}_i{y}_i+{{\mathrm{Cy}}_i}^2+{\mathrm{Dx}}_i+{\mathrm{Ey}}_i+F)}^2$

其中，n为Canny边缘检测后的离散边缘点的个数，x_i、y_i为第i个离散边缘点的x和y坐标，约束条件设为A+C＝1，未知参数求解准则为：

$\frac{\∂f}{\∂B}=\frac{\∂f}{\∂C}=\frac{\∂f}{\∂D}=\frac{\∂f}{\∂E}=\frac{\∂f}{\∂F}=0;$

圆的中心坐标(x_c,y_c)为：

$x_c=\frac{BE-2CD}{4AC-B^2},y_c=\frac{BD-2AE}{4AC-B^2}.$

5.根据权利要求4所述的基于自适应滤波的复杂背景下图像多圆快速检测方法，其特征在于，所述步骤2)中，求解的半径为：

$r=\frac{a+b}{2};$

其中，