CN104715487A - 一种基于伪Zernike矩的亚像素边缘检测方法 - Google Patents

Abstract

本发明公开了一种基于伪Zernike矩的亚像素边缘检测方法，其特征在于包括以下步骤：S1：对输入图像进行去除噪处理；S2：将完成去噪处理的图像进行像素级边缘检测：S3：采用伪Zernike矩方法对待处理图像进行亚像素边缘检测：S4：对待处理图像进行边缘位置的误差补偿：S5：获得亚像素边缘检测的修正实际边缘，对待处理图像的所有像素按照S4的方式进行处理、完成图像的亚像素边缘检测。本发明所提的方法对噪声不敏感，提高了亚像素边缘的精度，并且减少了检测边缘所需的计算复杂度。

Description

一种基于伪Zernike矩的亚像素边缘检测方法

技术领域

本发明涉及图像处理技术领域，尤其涉及一种基于伪Zernike矩的亚像素边缘检测方法。

背景技术

随着我国经济的快速发展，工业领域的机器化、智能化的程度将飞快的提高，图像的边缘检测技术的发展空间也是相当可观。基于计算机视觉的边缘检测技术具有精度高、速度快、非接触、自动化程度高等优点。由于图像边缘包含了大量的有用信息，边缘检测结果的精确与否，对后续图像的处理，如物体配准、物体的尺寸测量、物体的检测与识别等，也有着重要的影响，因此准确地提取图像边缘在基于计算机视觉的检测***中占有重要地位。在计算机视觉***中，***的精度与边缘检测的精度成一定性的正比，即如果目标物体的边缘检测的位置精度高，物体的特征信息提取的信息量多，后续的相关的处理的结果及精度将高。提高***精度最直接、有效的方法是提高***的硬件的分辨率，但提高硬件分辨率所需的成本较高，例如，将摄像头256×256的分辨率提高到1024×1024的分辨率，***所需的成本要相差几十倍的价格，用计算机视觉技术提高分辨率，在提高检测精度的同时，可有效降低***成本。

传统的边缘检测算子有：Sobel边缘检测算子，Roberts边缘检测算子，Prewitt边缘检测算子，Log边缘检测算子，Canny边缘检测算子等，这些边缘检测的精度都是像素级，即，定位的精度都是一个像素。随着对检测精度的要求不断提高，像素级的边缘检测已经不能满足实际工业生产的需求，人们提出了亚像素边缘检测技术，例如当检测的精度是0.2像素时，相当于***分辨率提高5倍。

现有技术中专利号为CN101477685公开了一种具有景深零件加工质量的亚像素级图像检测方法。该发明首先对机器视觉***进行分层标定，其次，对原始图像进行插值计算，通过粗精两步法实现零件边缘的精确定位，最后利用建立的图像各层面与零件各层面映射关系计算具有景深零件的形状，关键尺寸参数，通过比较分析获得质量数据。虽然该技术有着较快的检测速度，但是该技术是利用插值的方法来进行亚像素边缘检测的，由于插值技术本身很容易受噪声的影响，那么该技术也很容易受噪声的干扰，这样会造成边缘检测的精度降低，进而影响后续的图像处理性能。

另外2012年，Kaur等在文献“Sub-pixel edge detection using pseudoZernike moment”中提出利用伪Zernike来进行亚像素边缘检测。该技术首先提取获取图像，进行伪Zernike矩运算，其次，获取边缘方向参数分布和边缘方向参数差分值然后，根据预设阈值T来判断是否边缘像素。但是该技术通过伪Zernike矩进行亚像素边缘检测，虽然该技术对噪声不敏感，即克服了噪声的影响，但是该技术使用伪Zernike的计算方法，由于伪Zernike矩计算复杂度较大，那么这样便影响了计算的速度。

发明内容

本发明针对现有的亚像素边缘检测方法精度不高，对噪声敏感，计算复杂度高等问题，提出了一种新的基于伪Zernike矩的亚像素边缘检测方法，具体方案是：

一种基于伪Zernike矩的亚像素边缘检测方法，包括以下步骤：

S1：对输入图像进行去除噪处理；

S2：将完成去噪处理的图像进行像素级边缘检测：以待处理图像中待处理像素为中心，对该像素点的周围四个方向的像素的灰度进行加权运算，从水平和垂直的方向进行边缘检测，按照上述方式对待处理图像的所有像素进行像素级边缘检测；

S3：采用伪Zernike矩方法对待处理图像进行亚像素边缘检测：建立边缘检测模型，对待处理图像的所有像素进行如下处理：计算像素的正交复多项式、利用正交复多形式的结果计算像素相关矩的系数，利用相关矩的系数计算相关矩的大小，利用相关矩的结果来计算出边缘的参数，利用边缘参数计算像素真实的边缘位置；

S4：对待处理图像进行边缘位置的误差补偿：建立误差补偿边缘模型，利用该模型求像素的伪Zernike矩，采用误差补偿边缘求得像素的边缘位置估值、实际的边缘误差和理论边缘误差；

S5：获得亚像素边缘检测的修正实际边缘，对待处理图像的所有像素按照S4的方式进行处理、完成图像的亚像素边缘检测。

S1中对输入图像进行去除噪处理采用如下方式：

S11：以待处理图像中待处理像素为中心，分别计算该像素周围四个窗口所对应的灰度方差；

S12：找出上述灰度方差最小所对应的窗口，并计算其灰度均值；

S13：将计算出的灰度均值代替中心像素的灰度值；对待处理图像的所有像素进行上述操作完成去噪处理。

S3中采用伪Zernike矩方法对图像进行像素边缘检测时采用如下算法：

${\mathrm{PZ}}_{\mathrm{nm}}=\frac{n+1}{\mathrm{\π}}\underset{X}{\mathrm{\Σ}}\underset{Y}{\mathrm{\Σ}}f(x,y)V_{\mathrm{nm}}^{*}(\mathrm{\ρ},\mathrm{\θ})$

其中：(n+1)/π是归一化参数，符号“*”表示复数的共轭计算，θ是边缘与x方向的夹角，ρ是中心到直线的距离，即像素边缘所在的位置，V_nm(ρ,θ)是正交积分内核函数，上述参数利用公式表示为：V_nm(ρ,θ)＝R_nm(ρ)e^imθ，极坐标下，复多项式R_nm(ρ)定义为

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1-s)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}$

其中，0≤|m|≤n，arctan()是反正切函数。

S4中对图像进行边缘位置的误差补偿具体采用如下方式：根据建立的误差补偿边缘模型采用公式(28)、(29)进行像素边缘矩的修正：

$\begin{array}{c}{\mathrm{PZ}}_{20}+4{\mathrm{PZ}}_{10}=\∫{\∫}_{x^2+y^2}{f}^{\′}(x,y)(10x^2+10y^2-5)\mathrm{dxdy}\\ ={\∫}_{-1}^{{\mathrm{\ρ}}_1}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(10x^2+10y^2-5)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_1}^{{\mathrm{\ρ}}_2}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)(10x^2+10y^2-5)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_2}^1{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)(10x^2+10y^2-5)\mathrm{dxdy}\\ =\frac{10\mathrm{\Δk}}{3}[{\mathrm{\ρ}}_1{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]+\frac{10k}{3}{\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}\end{array}---\left(28\right)$

$\begin{array}{c}{\mathrm{PZ}}_{11}=\∫{\∫}_{x^2+y^2}{f}^{\′}(x,y)(x-\mathrm{jy})\mathrm{dxdy}\\ ={\∫}_{-1}^{{\mathrm{\ρ}}_1}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(x-\mathrm{jy})\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_1}^{{\mathrm{\ρ}}_2}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(x-jy)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_2}^1{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)(x-jy)\mathrm{dxdy}\\ =\frac{2\mathrm{\Δk}}{3}[{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]+\frac{2k}{3}{(1-{\mathrm{\ρ}}_2^2)}^{3/2}\end{array}---\left(29\right)$

其中，f′(x,y)是旋转后图像的灰度值，

利用误差补偿边缘采用公式(30)求得像素的边缘位置ρ，采用公式(31)求解像素的实际的边缘和理论边缘误差E：

$\mathrm{\ρ}=\frac{{\mathrm{PZ}}_{20}+4{\mathrm{PZ}}_{10}}{5{\mathrm{PZ}}_{11}}=\frac{(1-\mathrm{\λ}){\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}+\mathrm{\λ}{\mathrm{\ρ}}_1{(1-{\mathrm{\ρ}}_1^2)}^{3/2}}{\mathrm{\λ}{(1-{\mathrm{\ρ}}_1^2)}^{3/2}+(1-\mathrm{\λ}){(1-{\mathrm{\ρ}}_2^2)}^{3/2}}---\left(30\right)$

则像素的实际边缘误差和理论的边缘的误差E为：

$E=\mathrm{\ρ}{-\mathrm{\ρ}}_R=\frac{\mathrm{\λ}(\mathrm{\λ}-1)({\mathrm{\ρ}}_1-{\mathrm{\ρ}}_2)[{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]}{\mathrm{\λ}{(1-{\mathrm{\ρ}}_1^2)}^{3/2}+(1-\mathrm{\λ}){(1-{\mathrm{\ρ}}_2^2)}^{3/2}}---\left(31\right)$

求解修正像素的实际边缘位置采用如下公式

ρ′_R＝ρ-E

其中ρ′_R是修正后像素的实际边缘位置，ρ是实际的边缘误差，E是边缘的误差值。

本发明公开的一种基于伪Zernike矩的亚像素边缘检测方法，首先利用Kuwahara滤波器进行去噪的处理，其次利用Sobel算子进行像素级边缘检测，然后利用伪Zernike矩进行像素边缘检测，最后利用误差补偿对边缘位置进行补偿，本发明所提的方法对噪声不敏感，提高了亚像素边缘的精度，并且减少了检测边缘所需的计算复杂度。

附图说明

为了更清楚地说明本申请实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本申请中记载的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明基于伪Zernike矩的亚像素边缘检测方法的流程图；

图2为本发明中对图像像素进行去噪处理的示意图；

图3(a)为二维阶跃边缘模型中旋转前的阶跃边缘模型示意图；

图3(b)为二维阶跃边缘模型中旋转后的阶跃边缘模型示意图；

图4为本发明中误差补偿边缘模型的示意图；

图5为本发明中不同角度的合成图像的亚像素边缘检测结果示意图；

图6为本发明中实际图像的亚像素边缘检测效果的示意图。

具体实施方式

为使本发明的技术方案和优点更加清楚，下面结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚完整的描述：

本发明首先利用Kuwahara滤波器对图像进行去噪处理，其次利用Sobel算子进行像素级边缘检测，然后利用伪Zernike矩进行像素边缘检测，最后利用误差补偿对图像边缘位置进行补偿，其中利用伪Zernike矩进行亚像素边缘检测的方法和利用边缘位置误差补偿来进行边缘的修正的方法是本发明的发明点。

如图1所示的一种基于伪Zernike矩的亚像素边缘检测方法，包括以下步骤：

S1：对输入图像进行去除噪处理。

输入图像若含有噪声，将影响像素级边缘检测结果，故先应用Kuwahara滤波器对其去噪。Kuwahara滤波器的基本思想是选取如图2所示的四个窗口，对每个窗口分别计算其灰度方差，然后用方差最小的窗口所对应的灰度均值代替中心像素的灰度值，其具体的实现步骤如下：

S11：以待处理图像中待处理像素为中心，分别计算该像素周围四个窗口所对应的灰度方差；

S12：找出上述灰度方差最小所对应的窗口，并计算其灰度均值；

S13：将计算出的灰度均值代替中心像素的灰度值；对待处理图像的所有像素进行上述操作完成去噪处理。

S2：将完成去噪处理的图像进行像素级边缘检测：以待处理图像中待处理像素为中心，对该像素点的周围四个方向的像素的灰度进行加权运算，从水平和垂直的方向进行边缘检测，按照上述方式对待处理图像的所有像素进行像素级边缘检测。

典型的Sobel算子是以处理的像素为中心，然后对像素点的周围4个方向的像素的灰度加权运算，从水平和垂直的方向来检测边缘，具体计算公式如下。

f′_x(x,y)＝f(x-1,y+1)+2f(x,y+1)+f(x+1,y+1)

                                                 (1)

-f(x-1,y-1)-2f(x,y-1)-f(x+1,y-1)

f′_y(x,y)＝f(x-1,y-1)+2f(x-1,y)+f(x-1,y+1)

                                                 (2)

-f(x+1,y-1)-2f(x+1,y)-f(x+1,y+1)

G[f′(x,y)]＝|f′_x(x,y)|+|f′_y(x,y)|    (3)

其中，f_x’(x,y)、f_y’(x,y)分别是x(水平)方向和y(垂直)方向的一阶微分，G[f(x,y)]为Sobel算子的梯度总和，f(x,y)是输入图像的灰度值。

在求出梯度G[f(x,y)]后，可以假给定设一个阈值常数T，进行图像二值化，即当G[f(x,y)]>T时，对应像素点设定为0，否则设定为255或者1，设当调整常数T的大小达到最佳效果，在本发明中T取值为10～30。

S3：采用伪Zernike矩方法对待处理图像进行亚像素边缘检测：建立边缘检测模型，对待处理图像的所有像素进行如下处理：计算像素的正交复多项式、利用正交复多形式的结果计算像素相关矩的系数，利用相关矩的系数计算相关矩的大小，利用相关矩的结果来计算出边缘的参数，利用边缘参数计算像素真实的边缘位置。

自Ghosal等在1993年提出利用Zernike矩进行亚像素边缘检测以来，经过二十多年的研究，学者们已经对其进行了深入的研究，并取得的很好的边缘检测效果。由于伪Zernike矩在具有Zernike矩的旋转不变性的同时，伪Zernike矩比Zernike矩能提供更多的特征向量，并且伪Zernike矩比传统的Zernike矩更对图像噪声不敏感，所以本发明利用伪Zernike矩进行像素边缘检测。

伪Zernike矩原理：在单位圆中，对于数字图像，Zernike矩可以定义为，

${\mathrm{PZ}}_{\mathrm{nm}}=\frac{n+1}{\mathrm{\π}}\underset{X}{\mathrm{\Σ}}\underset{Y}{\mathrm{\Σ}}f(x,y)V_{\mathrm{nm}}^{*}(\mathrm{\ρ},\mathrm{\θ})---\left(5\right)$

其中，(n+1)/π是归一化参数，*表示共轭，θ是边缘与x方向的夹角，ρ是中心到直线的距离，也是边缘所在的位置，V_nm(ρ,θ)是正交积分内核函数，上述参数利用公式表示为，

$\mathrm{\ρ}=\sqrt{x^2+y^2}---\left(6\right)$

$\mathrm{\θ}=\arctan \frac{y}{x}---\left(7\right)$

V_nm(ρ,θ)＝R_nm(ρ)e^imθ    (8)

极坐标下，复多项式R_nm(ρ)定义为，

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1-s)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}---\left(9\right)$

其中，0≤|m|≤n。

基于伪Zernike矩的图像像素边缘检测方法具体步骤是：

为了进行边缘检测，需建立理想的边缘检测模型，下式和图3进行详细的介绍，

$f(x,y,\mathrm{\ρ},\mathrm{\θ})=\left\{\begin{array}{cc}h,& x\cos \mathrm{\θ}+y\sin \mathrm{\θ}\≤\mathrm{\ρ}\\ h+k,& x\cos \mathrm{\θ}+y\sin \mathrm{\θ}>\mathrm{\ρ}\end{array}\right.---\left(10\right)$

其中，h和h+k分别是直线左右两侧的灰度值，θ是边缘与x方向的夹角，ρ是中心到直线的距离，也是边缘所在的位置。

如图3(a)旋转前的阶跃边缘模型，如图3(b)旋转后的阶跃边缘模型

如图3所示为了计算边缘的参数，PZ₀₀、PZ₁₀、PZ₁₁、PZ₂₀四个矩需要计算，本发明仅利用了PZ₁₀，PZ₁₁，PZ₂₀三个矩进行计算。

根据图3的阶跃边缘模型，本发明的详细步骤如下所示

(1)计算正交复多项式，其计算过程详见下式：

$R_{10}=\underset{s=0}{\overset{1}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1-s)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}=-2+3\mathrm{\ρ}---\left(11\right)$

$R_{11}=\underset{s=0}{\overset{0}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1-s)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}=\mathrm{\ρ}---\left(12\right)$

$R_{20}=\underset{s=0}{\overset{2}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1-s)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}=3-12\mathrm{\ρ}+10{\mathrm{\ρ}}^2---\left(13\right)$

为了方便计算，将R₂₀和R₁₀进行合并运算，其计算结果是：

R₂₀+4R₁₀＝-5+10ρ²    (14)

(2)计算PZ₁₁，PZ₂₀+4PZ₁₀的系数其计算过程详见下式，

${\mathrm{CPZ}}_{\mathrm{nm}}=\frac{n+1}{\mathrm{\π}}{\∫\∫}_{x^2+y^2\≤1}{V}_{\mathrm{nm}}^{*}(\mathrm{\ρ},\mathrm{\θ})\mathrm{dxdy}---\left(15\right)$

${\mathrm{CPZ}}_{11R}={\∫\∫}_{x^2+y^2}\left(x\right)\mathrm{dxdy}---\left(16\right)$

${\mathrm{CPZ}}_{11I}={\∫\∫}_{x^2+y^2}\left(y\right)\mathrm{dxdy}---\left(17\right)$

${\mathrm{CPZ}}_{20}+4{\mathrm{PZ}}_{10}={\∫\∫}_{x^2+y^2}(10x^2+10y^2-5)\mathrm{dxdy}---\left(18\right)$

本发明用5×5的模板经过计算得到的PZ₁₁，PZ₂₀+4PZ₁₀三个矩的系数分别为如下表所示：

表1PZ₁₁矩的实系数CPZ_11R

0.0452	0.0226	0.0	-0.0226	-0.0452
					0.0453	0.0226	0.0	-0.0226	-0.0453
0.0452	0.0226	0.0	-0.0226	-0.0452
					0.0453	0.0226	0.0	-0.0226	-0.0453
0.0452	0.0226	0.0	-0.0226	-0.0452

表2 PZ₁₁矩的虚系数CPZ_11I

0.0452	0.0453	0.0452	0.0453	0.0452
					0.0226	0.0226	0.0226	0.0226	0.0226
0.0	0.0	0.0	0.0	0.0
					-0.0226	-0.0226	-0.0226	-0.0226	-0.0226
-0.0452	-0.0453	-0.0452	-0.0453	-0.0452

表3 PZ₂₀+4PZ₁₀矩的系数CPZ₂₀+4CPZ₁₀

0.1227	-0.0692	-0.1333	-0.0692	0.1227
					-0.0692	-0.2613	-0.3255	-0.2613	-0.0692
-0.1333	-0.3255	-0.3893	-0.3255	-0.1573
					-0.0692	-0.2613	-0.3255	-0.2613	-0.0692
0.1227	-0.0692	-0.1573	-0.0692	0.1227

(3)然后利用以下的公式，分别求PZ₁₁，PZ₂₀+4PZ₁₀，

${\mathrm{PZ}}_{11R}=\underset{i=-2}{\overset{2}{\mathrm{\Σ}}}\underset{j=-2}{\overset{2}{\mathrm{\Σ}}}f(i+m,j+n){\mathrm{CPZ}}_{11R}---\left(19\right)$

${\mathrm{PZ}}_{11I}=\underset{i=-2}{\overset{2}{\mathrm{\Σ}}}\underset{j=-2}{\overset{2}{\mathrm{\Σ}}}f(i+m,j+n){\mathrm{CPZ}}_{11I}---\left(20\right)$

${\mathrm{PZ}}_{20}+4{\mathrm{PZ}}_{10}=\underset{i=-2}{\overset{2}{\mathrm{\Σ}}}\underset{j=-2}{\overset{2}{\mathrm{\Σ}}}f(i+m,j+n)({\mathrm{CPZ}}_{20}+4{\mathrm{CPZ}}_{10})---\left(21\right)$

其中，f(n,m)是像素边缘检测的位置的灰度值。

(4)边缘参数的求解，其计算过程详细下式，

$\mathrm{\θ}=\arctan \frac{\mathrm{IM}\left[{\mathrm{PZ}}_{11}\right]}{\mathrm{RE}\left[{\mathrm{PZ}}_{11}\right]}---\left(22\right)$

$\mathrm{\ρ}=\frac{{\mathrm{PZ}}_{20}+4{\mathrm{PZ}}_{10}}{5{\mathrm{PZ}}_{11}}---\left(23\right)$

$k=\frac{{3\mathrm{PZ}}_{11}}{{(1-{\mathrm{\ρ}}^2)}^{3/2}}---\left(24\right)$

(5)利用式(22～24)就能将边缘的参数(ρ,θ)求出来，那么真实的边缘的计算过程为,

$\left[\begin{array}{c}{x}_{\mathrm{sub}}\\ y_{\mathrm{sub}}\end{array}\right]=\left[\begin{array}{c}x\\ y\end{array}\right]+\frac{N}{2}\mathrm{\ρ}\begin{array}{}\end{array}\left[\begin{array}{c}\cos \mathrm{\θ}\\ \sin \mathrm{\θ}\end{array}\right]---\left(25\right)$

其中，x，y是像素级边缘检测出来的位置，N＝25。

S4：对待处理图像进行边缘位置的误差补偿：建立误差补偿边缘模型，利用该模型求像素的伪Zernike矩，采用误差补偿边缘求得像素的边缘位置、实际的边缘误差和理论边缘误差。

现有技术中的边缘模型是基于阶跃边缘模型的，而实际的边缘都不会是阶跃模型，这是因为实际的图像经过基于CCD的摄像头的凸镜都会有一定的模糊，退一步是即使经过凸镜以后不会模糊，但是数据经过数据A/D采集变换，量化等数字化操作，边缘也不可能是阶跃模型，所以经过现有技术和上文所介绍的伪Zernike矩进行亚像素边缘检测得到数据还是不准确，那么在此，为了改善上面所介绍的问题，本发明将使用一种边缘误差技术来减少因为模型误差所造成的影响。

为了说明边缘的误差补偿原理，本发明建立一个误差补偿边缘模型，具有模型如图4所示：图4中的h是背景灰度值，h+Δk是过度带的灰度值，h+k是目标对象的灰度值，ρ₁和ρ₂是单位圆的中心到两个边缘的距离。由图4的中各个参数的关系，可推出，

$\frac{\mathrm{\Δk}}{k}=\frac{{\mathrm{\ρ}}_2-{\mathrm{\ρ}}_R}{{\mathrm{\ρ}}_2-{\mathrm{\ρ}}_1}---\left(26\right)$

其中，ρ_R是实际的边缘。

如果令Δk/k＝λ，那么由上式可得

ρ_R＝ρ₂-λ(ρ₂-ρ₁)    (27)

利用上面的补偿边缘模型来求伪Zernike矩，那么修正后的各个矩详见下式：

$\begin{array}{c}{\mathrm{PZ}}_{20}+4{\mathrm{PZ}}_{10}=\∫{\∫}_{x^2+y^2}{f}^{\′}(x,y)(10x^2+10y^2-5)\mathrm{dxdy}\\ ={\∫}_{-1}^{{\mathrm{\ρ}}_1}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(10x^2+10y^2-5)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_1}^{{\mathrm{\ρ}}_2}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)(10x^2+10y^2-5)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_2}^1{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)(10x^2+10y^2-5)\mathrm{dxdy}\\ =\frac{10\mathrm{\Δk}}{3}[{\mathrm{\ρ}}_1{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]+\frac{10k}{3}{\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}\end{array}---\left(28\right)$

$\begin{array}{c}{\mathrm{PZ}}_{11}=\∫{\∫}_{x^2+y^2}{f}^{\′}(x,y)(x-\mathrm{jy})\mathrm{dxdy}\\ ={\∫}_{-1}^{{\mathrm{\ρ}}_1}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(x-\mathrm{jy})\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_1}^{{\mathrm{\ρ}}_2}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(x-jy)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_2}^1{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)(x-jy)\mathrm{dxdy}\\ =\frac{2\mathrm{\Δk}}{3}[{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]+\frac{2k}{3}{(1-{\mathrm{\ρ}}_2^2)}^{3/2}\end{array}---\left(29\right)$

其中，f′(x,y)是旋转后图像的灰度值。

利用误差补偿边缘求得边缘位置ρ为下式所示，

$\mathrm{\ρ}=\frac{{\mathrm{PZ}}_{20}+4{\mathrm{PZ}}_{10}}{5{\mathrm{PZ}}_{11}}=\frac{(1-\mathrm{\λ}){\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}+\mathrm{\λ}{\mathrm{\ρ}}_1{(1-{\mathrm{\ρ}}_1^2)}^{3/2}}{\mathrm{\λ}{(1-{\mathrm{\ρ}}_1^2)}^{3/2}+(1-\mathrm{\λ}){(1-{\mathrm{\ρ}}_2^2)}^{3/2}}---\left(30\right)$

实际的边缘误差和理论的边缘的误差E为

$E=\mathrm{\ρ}{-\mathrm{\ρ}}_R=\frac{\mathrm{\λ}(\mathrm{\λ}-1)({\mathrm{\ρ}}_1-{\mathrm{\ρ}}_2)[{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]}{\mathrm{\λ}{(1-{\mathrm{\ρ}}_1^2)}^{3/2}+(1-\mathrm{\λ}){(1-{\mathrm{\ρ}}_2^2)}^{3/2}}---\left(31\right)$

由于本发明选取的模板大小是5×5，故λ可取0.5，ρ₁，ρ₂分别取-0.2,0.2。基于误差补偿的亚像素边缘检测的修正实际边缘为

ρ′_R＝ρ-E    (32)

S5：获得亚像素边缘检测的修正实际边缘，对待处理图像的所有像素按照S4的方式进行处理完成图像的亚像素边缘检测。

实施例：

本发明的有益效果在于：为了验证本发明的有效性，进行了计算机仿真实验。在实验中，实验参数为CPU IntelR CoreTM i3 2.4GHz,2G内存，显卡是ATIMobility Radeon HD 5470，***为Window7家庭版，软件编程环境为Matlab2010b。本发明实验的图像是利用人工合成的图像和实际图像，对于人工合成的图片的大小是400像素×400像素，而对于实际的图片选择是512像素×512像素。

本发明参考两个文献，即文献[2]、文献[5]，文献[6]，将发明方法与文献[2]、文献[5]，文献[6]进行了仿真对比实验，具体的仿真结果如图5及表5、表6和表7所示。

表5不同方法的边缘位置的误差

表6不同方法的边缘方向的误差

表7不同方法的边缘检测时间

从表5可以看出，在不同角度直线的实验情况下，文献[6]的位置检测的精度比文献[5]的检测精度高，而文献[2]比文献[6]的检测精度高，而本文对文献[2]进行了改善，其检测的精度又有了进一步的提高；从表6可以看出，对不同角度直线的检测，不同的方法有着不同的结果，文献[6]的方向检测的精度比文献[5]和文献[2]的检测精度高，这是因为文献[6]利用了迭代的方法，这在一定程度上改善了方向检测精度，而文献[2]比文献[5]的检测精度高，而本文对文献[2]进行了改善，其检测的精度又有了进一步的提高，并且其检测精度最高。从表7可以看出本文所提的方法对文献[2]计算复杂度进行很大的改善，与文献[5]相当。

如图5所示：图5不同角度的合成图像的亚像素边缘检测结果：(a)～(f)没有加白噪声的边缘检测结果；(a)角度为零的检测结果；(b)角度为零的检测结果的局部放大图；(c)角度为-45度的检测结果；(d)角度为-45度的检测结果的局部放大图；(e)角度为45度的检测结果；(f)角度为45度的检测结果的局部放大图；(g)～(l)加入噪声的边缘检测结果；(g)角度为零的检测结果；(h)角度为零的检测结果的局部放大图；(i)角度为-45度的检测结果；(j)角度为-45度的检测结果的局部放大图；(k)角度为45度的检测结果；(l)角度为45度的检测结果的局部放大图。

图6实际图像的亚像素边缘检测效果：(a)～(c)具有丰富边缘的图像检测结果；(a)b图的局部放大图；(b)边缘的检测结果；(c)b图的局部放大图；(d)～(g)简单圆环工件的检测图像检测结果；(d)e图的局部放大图；(e)圆环的检测结果；(g)e图的局部放大图。

从图5可以看出，本文基于一维灰度矩的方法对合成的直线的灰度图像有着很好的检测精度，并且有很好的抗噪能力。而从图6可以看出，本文的方法在提取实际的物体边缘也有着较好的效果，图6(c)显示对于有着较多边缘的玩具虎，其提取的边缘精度较高，图6(g)表明对于简单的环形工件，同样也有着很好的检测结果。

本发明的去噪技术是利用Kuwahara滤波器进行平滑去噪，其他可以代替的平滑滤波可以为，均值滤波器，中值滤波器，高斯滤波器，方向平滑滤波器等。

本发明的像素级边缘检测采用的是Sobel边缘检测，目前已经很多像素级边缘检测技术，比如Robert边缘检测，Laplace边缘检测，Log边缘检测，Canny边缘等，可以用这些方法代替本发明的像素级边缘检测方法。

本发明中伪Zernike矩采用模板的尺寸是5×5的大小，而其他尺寸的大小比如3×3，7×7，9×9等都可以代替本发明中模板的尺寸。

参考文献(如专利/论文/标准)

[1]曹建福，史博，汪霖，张家良.具有景深零件加工质量的亚像素级图像检测方法.中国，公开号：CN101477685B.2011.06.01.

[2]Kaur A,Singh C.Sub-pixel edge detection using pseudo Zernike moment[J].Int.J.Signal Process.Image Process.Pattern Recognit.2011,4(2):107-118.

[3]Jain A K.Fundamentals of digital image processing[M].New York:Prentice-Hall,1989.

[4]Sobel I.Camera models and machine perception[D].Stanford:StanfordUniversity,1970.

[5]Sun Q,Hou Y,Tan Q,et al.A robust edge detection method with sub-pixelaccuracy[J].Optik International Journal for Light and Electron Optics,2014,125(14):3449-3453.

[6]Chen P,Chen F,Han Y,et al.Sub-pixel dimensional measurement with Logisticedge model[J].Optik-International Journal for Light and Electron Optics,2014,125(9):2076-2080.

以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，根据本发明的技术方案及其发明构思加以等同替换或改变，都应涵盖在本发明的保护范围之内。

Claims (5)

1.一种基于伪Zernike矩的亚像素边缘检测方法，其特征在于包括以下步骤：

S1：对输入图像进行去除噪处理；

S2：将完成去噪处理的图像进行像素级边缘检测：以待处理图像中待处理像素为中心，对该像素点的周围四个方向的像素的灰度进行加权运算，从水平和垂直的方向进行边缘检测，按照上述方式对待处理图像的所有像素进行像素级边缘检测；

S3：采用伪Zernike矩方法对待处理图像进行亚像素边缘检测：建立边缘检测模型，对待处理图像的所有像素进行如下处理：计算像素的正交复多项式、利用正交复多形式的结果计算像素相关矩的系数，利用相关矩的系数计算相关矩的大小，利用相关矩的结果来计算出边缘的参数，利用边缘参数计算像素真实的边缘位置；

S4：对待处理图像进行边缘位置的误差补偿：建立误差补偿边缘模型，利用该模型求像素的伪Zernike矩，采用误差补偿边缘求得像素的边缘位置估值、实际的边缘误差和理论边缘误差；

S5：获得亚像素边缘检测的修正实际边缘，对待处理图像的所有像素按照S4的方式进行处理、完成图像的亚像素边缘检测。

2.根据权利要求1所述的一种基于伪Zernike矩的亚像素边缘检测方法，其特征还在于：S1中对输入图像进行去除噪处理采用如下方式：

S11：以待处理图像中待处理像素为中心，分别计算该像素周围四个窗口所对应的灰度方差；

S12：找出上述灰度方差最小所对应的窗口，并计算其灰度均值；

S13：将计算出的灰度均值代替中心像素的灰度值；对待处理图像的所有像素进行上述操作完成去噪处理。

3.根据权利要求1所述的一种基于伪Zernike矩的亚像素边缘检测方法，其特征还在于：S3中采用伪Zernike矩方法对图像进行像素边缘检测时采用如下算法：

$P{Z}_{\mathrm{nm}}=\frac{n+1}{\mathrm{\π}}\underset{X}{\mathrm{\Σ}}\underset{Y}{\mathrm{\Σ}}f(x,y)V_{\mathrm{nm}}^{*}(\mathrm{\ρ},\mathrm{\θ})$

其中：(n+1)/π是归一化参数，符号“*”表示复数的共轭计算，θ是边缘与x方向的夹角，ρ是中心到直线的距离，即像素边缘所在的位置，V_nm(ρ,θ)是正交积分内核函数，上述参数利用公式表示为：V_nm(ρ,θ)＝R_nm(ρ)e^imθ，极坐标下，复多项式R_nm(ρ)定义为

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1-s)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}$

其中，0≤|m|≤n，arctan()是反正切函数。

4.根据权利要求1所述的一种基于伪Zernike矩的亚像素边缘检测方法，其特征还在于：S4中对图像进行边缘位置的误差补偿具体采用如下方式：根据建立的误差补偿边缘模型采用公式(28)、(29)进行像素边缘矩的修正：

$\begin{array}{c}{\mathrm{PZ}}_{20}+4P{Z}_{10}=\∫{\∫}_{x^2+y^2}{f}^{\′}(x,y)({10x}^2+{10y}^2-5)\mathrm{dxdy}\\ ={\∫}_{-1}^{{\mathrm{\ρ}}_1}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h({10x}^2+{10y}^2-5)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_1}^{{\mathrm{\ρ}}_2}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)({10x}^2+{10y}^2-5)\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_2}^1{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)({10x}^2+{10y}^2-5)\mathrm{dxdy}\\ =\frac{10\mathrm{\Δk}}{3}[{\mathrm{\ρ}}_1{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]+\frac{10k}{3}{\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}\end{array}---\left(28\right)$

$\begin{array}{c}{\mathrm{PZ}}_{11}={\∫\∫}_{x^2+y^2}{f}^{\′}(x,y)(x-\mathrm{jy})\mathrm{dxdy}\\ ={\∫}_{-1}^{{\mathrm{\ρ}}_1}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(x-\mathrm{jy})\mathrm{dxdy}\\ +{\∫}_{p_1}^{{\mathrm{\ρ}}_2}{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}h(x-\mathrm{jy})\mathrm{dxdy}\\ +{\∫}_{{\mathrm{\ρ}}_2}^1{\∫}_{-\sqrt{1-x^2}}^{\sqrt{1-x^2}}(h+k)(x-\mathrm{jy})\mathrm{dxdy}\\ =\frac{2\mathrm{\Δk}}{3}[{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]+\frac{2k}{3}{(1-{\mathrm{\ρ}}_2^2)}^{3/2}\end{array}---\left(29\right)$

其中，f′(x,y)是旋转后图像的灰度值，

利用误差补偿边缘采用公式(30)求得像素的边缘位置ρ，采用公式(31)求解像素的实际的边缘和理论边缘误差E：

$\mathrm{\ρ}=\frac{{\mathrm{PZ}}_{20}+4P{Z}_{10}}{5{\mathrm{PZ}}_{11}}=\frac{(1-\mathrm{\λ}){\mathrm{\ρ}}_2{(1-{\mathrm{\ρ}}_2^2)}^{3/2}+{\mathrm{\λ\ρ}}_1{(1-{\mathrm{\ρ}}_1^2)}^{3/2}}{\mathrm{\λ}{(1-{\mathrm{\ρ}}_1^2)}^{3/2}+(1-\mathrm{\λ}){(1-{\mathrm{\ρ}}_2^2)}^{3/2}}---\left(30\right)$

则像素的实际边缘误差和理论的边缘的误差E为：

$E=\mathrm{\ρ}-{\mathrm{\ρ}}_R=\frac{\mathrm{\λ}(\mathrm{\λ}-1)({\mathrm{\ρ}}_1-{\mathrm{\ρ}}_2)[{(1-{\mathrm{\ρ}}_1^2)}^{3/2}-{(1-{\mathrm{\ρ}}_2^2)}^{3/2}]}{{\mathrm{\λ}(1-{\mathrm{\ρ}}_1^2)}^{3/2}+(1-\mathrm{\λ}){(1-{\mathrm{\ρ}}_2^2)}^{3/2}}.---\left(31\right)$

5.根据权利要求1所述的一种基于伪Zernike矩的亚像素边缘检测方法，其特征还在于：求解修正像素的实际边缘位置采用如下公式

ρ′_R＝ρ-E

其中ρ′_R是修正后像素的实际边缘位置，ρ是实际的边缘误差，E是边缘的误差值。