CN102637300A - 改进的Zernike矩边缘检测方法 - Google Patents

Abstract

改进的Zernike矩边缘检测方法，涉及一种图像处理方法。提供一种针对交点附近像素点亚像素定位的改进的Zernike矩边缘检测方法。采用Zernike矩方法对工业CT图像进行亚像素边缘提取，并细化边缘，再采用8链码跟踪方法，获取边缘方向参数φ分布，计算边缘方向参数差分值Δφ，选取阈值T，判断边缘点是否为交点附近的像素点，若不是交点附近边缘点，则进行最小二乘拟合，得到直线方程，再通过直线方程计算交点附近待修正边缘点的亚像素坐标。改进Zernike矩边缘检测方法很好地弥补了Ghosal算法在交点附近像素点定位误差较大的不足，对高精度面积测量和高精度三维测量有重要意义。

Description

改进的Zernike矩边缘检测方法

技术领域

本发明涉及一种图像处理方法，尤其是涉及一种针对交点亚像素定位的改进的Zernike矩边缘检测方法。

背景技术

工业CT(Industry Computerized Tomogragh)广泛应用于航空、航天、汽车等等工业中精密工件内部结构的测量及其缺陷检测领域。工业CT图像尺寸测量方法是一种以工业CT图像为信息载体并从中提取定量数据的测量方法，它以工业CT断层图像为研究对象，计算目标的长度、宽度、高度、壁厚等各种几何参数。与其它测量方法相比，其既能测得实体复杂的内部几何形状(尤其是封闭内腔)的参数，又不破坏实体。边缘提取是工业CT图像测量的基本任务，它对测量结果的精度有很大的影响。传统的边缘检测算子如Sobel，Canny等只能提取像素级边缘不能满足定位的高精度需求。Hueckel首先提出亚像素边缘检测算法，到目前为止，亚像素级边缘检测方法主要分为四大类：插值法、拟合法、基于小波变换的方法、基于矩的方法、基于活动轮廓模型的方法。其中基于矩的方法具有快速、高精度且有一定抗噪能力的优势。

Ghosal和Mehrotal首次提出了利用Zernike正交矩来检测亚像素边缘（参见Ghosal S,Mehrotra R.Orthogonal moment operators for sub-pixel edge detection.PatternRecognition,1993,26(2):295-306），在他们的算法中建立了理想的阶跃灰度模型，通过图像的三个不同阶次Zernike矩计算模型的4个参数，以这4个参数为判断边缘的依据并确定图像中物体边缘的位置。分析发现，Ghosal提出的Zernike矩方法是基于直线边缘的，当检测目标为交点附近像素点时，算法将每个模板内的折线边缘作为直线处理，亚像素检测效果较差。

因此，需要一种改进的Zernike矩边缘检测算法，解决Ghosal算法在利用Zernike矩进行亚像素边缘检测时交点附近的边缘检测误差大的问题。

发明内容

有鉴于此，本发明所要解决的技术问题是提供一种改进的Zernike矩边缘检测方法，对Ghosal算法在利用Zernike矩进行亚像素边缘检测时交点附近的亚像素边缘误差大的不足进行判断和修正。

本发明的目的是这样实现的：

本发明提供的改进的Zernike矩边缘检测方法，包括以下步骤：

1）获取图像，提取图像的边缘；

2）对边缘进行细化处理；

3）获取边缘方向参数φ分布和边缘方向参数差分值Δφ；

4）根据预设阈值T来判断边缘像素点是否为交点附近像素点，如果是，则该像素点为交点附近像素点；

5）如果否，则该像素点为非交点附近像素点，判断为待修正像素点，并利用非交点附近像素点进行最小二乘拟合，得到直线方程P；

6）利用直线方程P来计算待修正像素点的坐标。

进一步，所述步骤2）中边缘点细化是对图像中的所有边缘点做边缘梯度方向上边缘距离参数l的非极小值抑制来细化边缘。

进一步，所述步骤3）中对所有边缘点采用8链码跟踪获取边缘方向参数φ分布。

进一步，所述步骤3）中采用以下公式来计算检测到边缘点的边缘方向参数φ的向前差分值Δφ：

Δφ＝φ(n+1)-φ(n),n＝1,2,3...(NUM-1)，

其中，Δφ表示边缘方向参数向前差分值，φ(n+1)表示后一个边缘点的边缘方向参数，φ(n)表示当前边缘点的边缘方向参数，NUM为检测到的像素点总数。

进一步，所述步骤4）中判断像素点是否为交点附近像素点包括以下步骤：

如果边缘方向参数差分值满足关系式|Δφ|≥T时，则边缘像素点为交点附近待修正像素点，如果边缘方向参数差分值满足关系式|Δφ|＜T，则边缘像素点为非交点附近像素点。

进一步，所述步骤5）中最小二乘拟合直线方程根据以下公式进行：

$\min \underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}^2(X_i,P)$

其中，F(X,P)＝Ax+By+C＝0表示拟合的直线方程，P表示直线参数P＝[A,B,C]，X＝[x,y,1]^T表示坐标向量，对于直线上任意边缘点

则F(X,P)表示该点到F(X,P)＝0的代数距离，根据目标函数求解直线参数P。

进一步，所述步骤6）中交点附近的像素点的精确坐标是通过直线参数P来确定的。

本发明的优点在于：本发明与传统的Zernike矩边缘检测方法相比，针对交点附近像素点亚像素定位误差较大进行改进，采用Zernike矩方法对工业CT图像进行亚像素边缘提取，并细化边缘，然后采用8链码跟踪方法获取边缘方向参数φ分布，并计算边缘方向参数差分值Δφ，再选取阈值T，判断边缘点是否为交点附近的像素点，若为非交点附近的边缘点，进行最小二乘拟合从而得到直线方程，再通过直线方程计算交点附近待修正边缘点的亚像素坐标。改进Zernike矩边缘检测方法很好地弥补了Ghosal算法在交点附近像素点定位误差较大的不足，对高精度面积测量和高精度三维测量有重要意义。

本发明的其它优点、目标和特征在某种程度上将在随后的说明书中进行阐述，并且在某种程度上，基于对下文的考察研究对本领域技术人员而言将是显而易见的，或者可以从本发明的实践中得到教导。本发明的目标和其它优点可以通过下面的说明书，权利要求书，以及附图中所特别指出的结构来实现和获得。

附图说明

为了使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明作进一步的详细描述，其中：

图1为本发明实施例提供的改进的Zernike矩边缘检测方法流程图；

图2为本发明实施例的矩形数字模拟二值图；

图3为本发明实施例的边缘方向参数φ分布图；

图4为本发明实施例的边缘方向参数差分值Δφ分布图；

图5为本发明实施例的亚像素边缘定位结果图。

具体实施方式

以下将结合附图，对本发明的优选实施例进行详细的描述；应当理解，优选实施例仅为了说明本发明，而不是为了限制本发明的保护范围。

图1为本发明实施例提供的改进的Zernike矩边缘检测方法流程图，

实施例1

本发明提供的改进的Zernike矩边缘检测方法，包括以下步骤：

1）获取图像，提取图像的边缘；本实施例中采用Zernike矩方法对图像进行边缘提取，本实施例中的图像为工业CT图像。

2）对边缘进行细化处理；本实施例的细化边缘是通过对图像中的每一个边缘点在梯度最大方向上做边缘距离参数l的非极小值抑制来进行的。

3）获取边缘方向参数φ分布和边缘方向参数差分值Δφ；本实施例采取8链码跟踪获取边缘方向参数φ分布；本实施例采用以下公式来计算检测到边缘点的边缘方向参数φ的向前差分值Δφ：

Δφ＝φ(n+1)-φ(n),n＝1,2,3...(NUM-1)，

其中，Δφ表示边缘方向参数向前差分值，φ(n+1)表示后一个边缘点的边缘方向参数，φ(n)表示当前边缘点的边缘方向参数，NUM为检测到的像素点总数。

4）根据预设阈值T来判断边缘像素点是否为交点附近像素点；判断像素点是否为交点附近像素点如下步骤：

如果边缘方向参数差分值满足关系式|Δφ|≥T时，则边缘像素点为交点像素点，如果边缘方向参数差分值满足关系式|Δφ|＜T，则边缘像素点为待修正像素点，如果是，则该像素点为交点附近像素点；

5）如果否，则该像素点为非交点附近像素点，判断为待修正像素点，并利用非交点附近像素点进行最小二乘拟合，得到直线方程P；最小二乘拟合直线方程根据以下公式进行：

$\min \underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}^2(X_i,P),$

其中，F(X,P)＝Ax+By+C＝0表示拟合的直线方程，P表示直线参数P＝[A,B,C]，X＝[x,y,1]^T表示坐标向量，对于直线上任意边缘点

则F(X,P)表示该点到F(X,P)＝0的代数距离，根据目标函数求解直线参数P。

6）利用直线方程P来计算待修正像素点的坐标。

实施例2

本发明提供的实施例2采用如下方式进行改进的Zernike矩边缘检测方法。

步骤一、Zernike矩方法提取边缘轮廓，Zernike矩亚像素边缘检测的具体实现步骤如下：

(1)计算模板Re[M₁₁]，Im[M₁₁]和M₂₀；

其中M_pq表示p阶q次Zernike矩模板，M₁₁表示一阶一次Zernike矩模板，M₂₀表示二阶零次Zernike矩模板，Re[M₁₁]表示M₁₁模板实部，Im[M₁₁]表示M₁₁模板虚部。

(2)利用模板Re[M₁₁]，Im[M₁₁]分别与图像卷积运算得到对应的Zernike矩实部Re[A₁₁]，和虚部Im[A₁₁]，M₂₀与图像卷积运算得到二阶零次Zernike矩A₂₀。p阶q次Zernike矩模板M_pq与图像卷积运算后得到p阶q次Zernike矩A_pq。

(3)根据公式

$\mathrm{\φ}={\tan }^{-1}\left(\frac{\mathrm{Im}\left[A_{11}\right]}{\mathrm{Re}\left[A_{11}\right]}\right),$

$l=\frac{A_{20}}{A_{11}^{\′}},$

$k=\frac{{3A}_{11}^{\′}}{2{(1-l^2)}^{3/2}},$

$h=\frac{A_{00}-\frac{\mathrm{k\π}}{2}+k{\sin }^{-1}\left(l\right)+\mathrm{kl}\sqrt{(1-l^2)}}{\mathrm{\π}},$

计算边缘模型参数φ，l，k，h；

其中，φ表示边缘中心与对应像素中心的连线和水平方向的夹角，l表示边缘中心到像素中心距离，k表示边缘相对于背景的灰度阶跃值，h表示背景灰度值，A′_pq为A_pq旋转φ角后的Zernike矩。

(4)选取阈值δ和τ，对每个像素点进行边缘点判定，若该点的参数满足条件l≤δ，k≥τ，则认为此点边缘点，利用式子

$\left[\begin{array}{c}{x}_s\\ y_s\end{array}\right]=\left[\begin{array}{c}x\\ y\end{array}\right]+\frac{\mathrm{Nl}}{2}\left[\begin{array}{c}\cos \left(\mathrm{\φ}\right)\\ \sin \left(\mathrm{\φ}\right)\end{array}\right],$

计算亚像素边缘点坐标，N为模版的大小，(x_s,y_s)表示亚像素位置，(x,y)表示中心像素位置。

步骤二、边缘细化

由于l≤δ,k≥τ的判定具有全局性，会产生提取边缘较粗的现象，所以进行了边缘细化处理，其步骤如下：

(1)Zernike矩算法亚像素边缘检测中保留每个点的边缘方向φ和边缘距离l；

(2)对图像中的所有边缘点做边缘距离l的非极小值抑制。所谓的l的非极小值抑制，即将边缘点的边缘方向φ分为四个范围0~π/4，π/4~π/2，π/2~3π/4，π/2~π，再比较每个边缘点的边缘距离l₀在边缘方向上相邻的两个像素点的边缘距离值l₁,l₂，判断l₀是否比l₁,l₂小，若不是，则认为此点不是边缘点。

图2为本发明实施例的矩形数字模拟二值图；图3为本发明实施例的边缘方向参数φ分布图；图4为本发明实施例的边缘方向参数差分值Δφ分布图；图3中横坐标为边缘像素点序号，纵坐标为像素点的边缘方向参数φ；图4中横坐标为边缘像素点序号，纵坐标为像素点的边缘方向参数差分值Δφ。

步骤三、判断边缘点是否为交点附近边缘点，其步骤如下：

(1)边缘跟踪，获取边缘方向参数φ分布

对检测结果中亚像素边缘进行8链码跟踪，获取边缘参数φ分布。

(2)计算边缘方向参数差分值Δφ

计算检测到边缘点的边缘方向参数φ的向前差分值Δφ：

Δφ＝φ(n+1)-φ(n),n＝1,2,3...(NUM-1)，

其中，Δφ表示边缘方向参数向前差分值，φ(n+1)表示后一个边缘点的边缘方向参数，φ(n)表示当前边缘点的边缘方向参数，NUM为检测到的像素点总数。

(3)选取阈值T，判断像素点是否为交点附近的像素点

选取阈值T，如果当差分值|Δφ|≥T时，则边缘像素点为交点附近像素点，如果|Δφ|＜T，则边缘像素点为非交点附近像素点。

步骤四、直线拟合

直线的一般方程为：

F(X,P)＝Ax+By+C＝0，

P＝[A,B,C]表示直线参数，X＝[x,y,1]^T表示坐标向量。对于直线上任意边缘点

则F(X,P)表示该点到F(X,P)＝0的代数距离，根据最小二乘原理，目标函数为：

$\min \underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}^2(X_i,P),$

求解出直线参数P。

步骤五、求交点附近像素点的亚像素坐标，

利用直线参数P，求取交点附近亚像素坐标。

图5为实施例中模拟二值图片的亚像素边缘定位结果，其中横坐标为列数，纵坐标为行数。

最后应说明的是，以上实施例仅用以说明本发明技术方案而非限制，尽管参照较佳实施例对本发明进行了说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者同等替换，而不脱离本发明方案的精神和范围。本发明就每个步骤所涉及的技术来说，都是可以用本技术领域的常识做到的，但这些步骤组合起来形成一种基于改进Zernike矩边缘技术的工业CT图像边缘定位方法。

Claims (7)

1.改进的Zernike矩边缘检测方法，其特点在于，包括以下步骤：

1）获取图像并提取图像的边缘；

2）对边缘进行细化处理；

3）获取边缘方向参数φ分布和边缘方向参数差分值Δφ；

4）根据预设阈值T来判断边缘像素点是否为交点附近像素点，如果是，则该像素点为交点附近像素点；

5）如果否，则该像素点为非交点附近像素点，判断为待修正像素点，并利用非交点附近像素点进行最小二乘拟合，得到直线方程P；

6）利用直线方程P来计算待修正像素点的坐标。

2.根据权利要求1所述改进的Zernike矩边缘检测方法，其特点在于，所述步骤2）中边缘点细化是对图像中的所有边缘点做边缘方向上边缘距离参数l的非极小值抑制来细化边缘。

3.根据权利要求1所述改进的Zernike矩边缘检测方法，其特点在于，所述步骤3）中对所有边缘点采用8链码跟踪获取边缘方向参数φ分布。

4.根据权利要求1所述改进的Zernike矩边缘检测方法，其特点在于，所述步骤3）中采用以下公式来计算检测到边缘点的边缘方向参数φ的向前差分值Δφ：

Δφ＝φ(n+1)-φ(n),n＝1,2,3...(NUM-1)，

其中，Δφ表示边缘方向参数向前差分值，φ(n+1)表示后一个边缘点的边缘方向参数，φ(n)表示当前边缘点的边缘方向参数，NUM为检测到的像素点总数。

5.根据权利要求1所述改进的Zernike矩边缘检测方法，其特点在于，所述步骤4）中判断像素点是否为交点附近像素点包括以下步骤：

如果边缘方向参数差分值满足关系式|Δφ|≥T时，则边缘像素点为交点附近待修正像素点，如果边缘方向参数差分值满足关系式|Δφ|＜T，则边缘像素点为非交点附近像素点。

6.根据权利要求1所述改进的Zernike矩边缘检测方法，其特点在于，所述步骤5）中最小二乘拟合直线方程根据以下公式进行：

$\min \underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}^2(X_i,P),$

其中，F(X,P)＝Ax+By+C＝0表示拟合的直线方程，P表示直线参数P＝[A,B,C]，X＝[x,y,1]^T表示坐标向量，对于直线上任意边缘点则F(X,P)表示该点到F(X,P)＝0的代数距离，根据目标函数求解直线参数P。

7.根据权利要求1所述改进的Zernike矩边缘检测方法，其特点在于，所述步骤6）中交点附近的像素点的精确坐标是通过直线参数P来确定的。

Images