CN101329402A - 基于改进Wedgelet的多尺度SAR图像边缘检测方法 - Google Patents

Abstract

本发明公开了一种基于改进多尺度Wedgelet变换的SAR图像边缘检测方法，它涉及图像处理领域。其目的在于解决传统边缘检测方法很难同时兼顾噪声抑制和对边缘完整准确的定位问题。其过程为：1.将SAR图像转化为灰度图像并将其分成多个子图像；2.对每个子图像，自适应地按照自下向上的尺度进行多尺度Wedgelet分解，用四叉树存储各个尺度的分解系数；3.根据分解得到的系数四叉树，对每个子图像进行多尺度Wedgelet逼近，最终得到整幅图像的多尺度Wedgelet逼近；4.根据得到的多尺度Wedgelet逼近图像，通过边缘提取梯度算子检测逼近图像的边缘；5.去除边缘图像中长度小于设定阈值的杂边缘对边缘图像的干扰。本发明具有运行时间快和边缘检测效果好的优点，可用于SAR图像的边缘检测。

Description

基于改进Wedgelet的多尺度SAR图像边缘检测方法

技术领域

本发明属于图像处理技术领域，涉及一种图像边缘检测方法，可用于对合成孔径雷达SAR图像的边缘检测。

背景技术

合成孔径雷达SAR是一种高分辨率成像雷达，具有全天候，多极化，多视角，多俯角数据获取能力，广泛应用于军事侦察及民用遥感探测领域。由于成像雷达发射的是纯相干波，因而SAR图像受到严重的相干斑噪声影响，图像的信噪比较低。边缘信息对于图像理解有着重要的意义，它既是图像分割的有效途径，又是特征提取的有效手段，在图像处理中已逐渐发展成专门的研究领域，也是SAR图像处理中研究的重点。但SAR图像相干斑噪声的存在给边缘检测造成了困难。

目前，常用的传统边缘检测方法主要有：(1)基于微分算子的边缘检测方法如Canny，Sobel算子，这类算子具有算法简单，运行速度快等优点，但是对图像噪声非常敏感。(2)基于假设检验的方法，如Ratio算法，它们能较好的克服噪声的影响，降低检测错误率，但是边缘定位精度较差。(3)基于小波理论的边缘检测方法，利用小波时频分析的优越性，检测图像在多个尺度下的边缘信息，获得多尺度分析下的图像边缘。但是由于SAR图像受强烈的相干斑噪声的影响，传统算法对图像的边缘检测很难同时兼顾噪声抑制和对边缘完整准确定位，对SAR图像的检测效果并不理想。

边缘检测的不确定性表明边缘检测算子的抑噪能力和定位精度是一对矛盾，小尺度算子虽然有利于边缘定位，但对噪声极为敏感；大尺度算子虽然抑噪能力强，但边缘定位精度差，甚至会丢失某些局部细节。多尺度边缘检测为解决此问题的提供了有效途径，例如，在多尺度几何分析工具中，美国学者Donoho等人提出的Wedgelet变换，具有良好的“线”和“面”的特性，但是由于受SAR图像斑点噪声的影响，该变换用于SAR图像时，逼近图像容易产生伪边缘，且运行时间过长，影响检测的效果和效率。

发明的内容

本发明的目的在于：克服现有SAR图像边缘检测中，传统算法很难同时兼顾噪声抑制和对边缘完整准确定位的缺点，而现有的Wedgelet逼近得到的图像具有容易产生伪边缘且运行时间过长的缺点，提出一种基于改进Wedgelet变换的SAR图像边缘检测方法，使检测后得到的边缘图像具有良好的抗噪性和边缘完整连续性的特点。

实现本发明目的的技术方案是：对传统的Wedgelet变换作改进，利用改进后的Wedgelet变换代价函数能够很好地捕捉图像中的“线”和“面”的特性，对输入图像进行多尺度Wedgelet变换，得到输入图像的多尺度Wedgelet逼近图像，再对逼近图像进行边缘检测和后处理，去除细小杂边缘的干扰，得到最终的边缘图像。具体实现过程如下：

(1)输入待检测的SAR图像，将其转化为灰度图像，并将得到的灰度图像分成多个子图像；

(2)用改进的Wedgelet代价函数对每个子图像自适应地按照自下向上的尺度进行多尺度Wedgelet分解，用四叉树保存各个尺度的分解系数，该Wedgelet变换的代价函数为： $\min \{{\left|\right|I\left(S\right)-\hat{I}\left(S\right)\left|\right|}_2^2-\mathrm{\η}*[\left|\widetilde{\mathrm{Ave}}\right\{y\left|P^{\′}\right\}-\widetilde{\mathrm{Ave}}\left\{y\right|P^{\′\′}\left\}\right|*\left({(n\×n)}^2\right)\left]\right\},$

其中，P′，P″分别表示图像子块S被划分成的两个楔形区域，

图像块S的大小为n×n，

分别表示区域P′，P″上像素的均值，

η为 $\left[\right|\widetilde{\mathrm{Av}}e\left\{y\right|P^{\′}\}-\widetilde{\mathrm{Av}}e\{y\left|P^{\′\′}\right\}|*\left({(n\×n)}^2\right)]$ 的权值；

(3)根据分解得到的Wedgelet系数，进行多尺度Wedgelet重构，得到每个子图像的多尺度Wedgelet逼近，对所有子图像都逼近后，得到输入图像的多尺度Wedgelet逼近；

(4)对得到的多尺度Wedgelet逼近图像，通过边缘提取梯度算子检测该逼近图像的边缘；

(5)搜索检测得到的边缘图像，存储每一条边缘的长度，去除边缘长度小于设定阈值的细小边缘，得到最终的边缘检测结果。

本发明由于采用多尺度Wedgelet逼近，利用Wedgelet能够良好捕捉图像中“线”和“面”的特性，因而抗噪性好和边缘检测结果完整连续；同时由于本发明对现有多尺度Wedgelet逼近方法进行了改进，并且对图像进行分块处理，进一步提高了其抗噪性和运行速度，降低了伪边缘的检测率，并且使其具有可应用性的特点；此外由于本发明采用了后处理操作，去除了杂边缘的干扰，得到了更好的边缘检测效果。

附图说明

图1是本发明的流程示意图；

图2是本发明采用的Wedgelet字典结构图；

图3是现有Wedgelet基函数示意图；

图4是本发明对一幅简单图像的Wedgelet逼近示意图；

图5是本发明采用的多尺度Wedgelet逼近图像边缘检测梯度算子结构图；

图6是本发明对一幅SAR图像Wedgelet逼近的结果图；

图7是本发明用于河流的SAR图像与现有方法的边缘检测结果比较图；

图8是本发明用于机场SAR图像与现有方法的边缘检测结果比较图；

图9是本发明用于机场跑道SAR图像与现有方法的边缘检测结果比较图；

图10是本发明用于港口SAR图像与现有方法的边缘检测结果比较图；

图11是本发明边缘检测结果与SAR图像原图叠加后的结果图。

具体实施方式

参照图1，本发明的具体实现过程如下：

步骤1，输入待检测的SAR图像，将其转化为灰度图像，并将得到的灰度图像分成多个子图像。

每个子图像大小为：M×M，假设输入图像大小为N×N，令R＝N/M，则得到R×R个子块图像。

设S(k₁，k₂，j)＝{(x₁，x₂)：[k₁/2^j，k₁+1/2^j]×[k₂/2^j，k₂+1/2^j]}(其中j≥0为一个整数，0≤k₁，k₂≤2^j)为二进正方形中的点的集合，任意两个不在同一条边上的顶点之间的连线就构成了一条edgelet线，每一条edgelet线和所在的二进正方形的边构成的楔形区域上的特征函数，称为非退化的Wedgelet，而一个没有被划分的正方形上的特征函数，则称为退化的Wedgelet。

步骤2，用改进的Wedgelet代价函数对每个子图像自适应地按照自下向上的尺度进行多尺度Wedgelet分解，并用四叉树存储各个尺度的分解系数。具体实施步骤如下：

1.初始化深度J＝log₂(M)的四叉数的最底层，其对应尺度为j＝J；

2.令j＝j-1，计算j尺度每个节点的代价函数，分以下三种情况：

(1)对于退化的叶节点，将其所对应的图像子块看作一个均匀的图像块，代价函数为：

$\cos \mathrm{tE}={\left|\right|I\left(S\right)-\widetilde{\mathrm{Ave}}\left\{S\right\}\left|\right|}^2+\mathrm{\λ}\#\left|P\right|---\left(1\right)$

式中，

表示图像块S上所有像素点的均值，P表示图像的ED-RDP，即逼近图像的Wedgelet的集合，#|P|表示集合P中的Wedgelet数目，λ为复杂度惩罚参数；

(2)对于非退化的叶节点，在其对应的图像块内，按图2所示Wedgelet字典结构图，根据本发明的最优Wedgelet函数搜索出最优的一条edgelet线，该函数为：

$\min \{{\left|\right|I\left(S\right)-\hat{I}\left(S\right)\left|\right|}_2^2-\mathrm{\η}*[\left|\widetilde{\mathrm{Ave}}\right\{y\left|P^{\′}\right\}-\widetilde{\mathrm{Ave}}\left\{y\right|P^{\′\′}\left\}\right|*\left({(n\×n)}^2\right)\left]\right\}---\left(2\right)$

式中，P′，P″分别为图像子块S被一条edgelet线划分得到的两个楔形区域，图像块S的大小为n×n，(n×n)²用来调节均值差在不同尺度的图像块中所占的比例的要求，η为均值差部分 $\left[\right|\widetilde{\mathrm{Ave}}\left\{y\right|P^{\′}\}-\widetilde{\mathrm{Ave}}\{y\left|P^{\′\′}\right\}|*\left({(n\×n)}^2\right)]$ 的权值；通过调节参数η，可以改变均值差部分的信息对代价函数的影响，η越大，表示均值差部分的信息对代价函数影响越大，本发明取η＝0.3。

将搜索到的最优的edgelet线对应的尺度、位置、方向和特征函数信息进行存储，并计算该edgelet线对应的最优Wedgelet表示的代价函数：

$\cos \mathrm{tP}={\left|\right|I\left(S\right)-\widetilde{\mathrm{Ave}}\left\{y\right|P\left\}\right||}^2+\mathrm{\λ}\#\left|P\right|---\left(3\right)$

参照图3，在j＝2时，将在单位正方形每条边上间隔δ＝2^-j设一个顶点，该单位正方形的四条边上共有4×2^j＝16个顶点，在二进正方形中的edgelet线左侧区域P′，构成了值为常数c_a的Wedgelet基函数，而在二进正方形中的edgelet线右侧区域P″，构成了值为常数c_b的Wedgelet基函数；c_a，c_b的值可由下式得到：

c_a＝Ave(I(S_j，k)|R_a)

                        (4)

c_b＝Ave(I(S_j，k)|R_b)

式中，I为该二进正方形[0，1]²，并且 $S_{j,,k}\⋐{\left[\mathrm{0,1}\right]}^2$ (j，k分别表示图像块S尺度和位置)，c_a，c_b即为相应Wedgelet区域内各像素点的均值。

(3)对于中间节点，由其4个子节点在下一层更细的尺度上表示，其代价函数为四个子节点的代价函数之和，即：

costD＝costS₁+costS₂+costS₃+costS₄    (5)

其中，S₁，S₂，S₃，S₄分别表示该节点对应的四个子节点。

3.比较每个节点三种代价函数值，取最小值作为当前节点代价函数值：

cost_optS＝min{costE，costP，costD}    (6)

4.如果cost_optS＝costE或costP，标记当前节点为叶子节点，剪枝掉以当前节点为根节点的所有子节点，并在分解四叉树中对应位置存储当前节点的最优Wedgelet表示信息，否则记录当前节点为中间节点；

5.循环执行(2)～(4)，满足停机条件j＝0，即对所有尺度执行完后停止，得到每个子图像完整的分解系数四叉树。

步骤3，根据分解得到的Wedgelet系数，进行多尺度Wedgelet重构，得到每个子图像的多尺度Wedgelet逼近，对所有子图像都逼近后，得到输入图像的多尺度Wedgelet逼近，具体实施步骤如下：

1.将分解系数四叉树，自顶向下遍历搜索整个四叉树；

2.如果搜索到叶节点，根据存储在该节点中的其对应二进正方形的Wedgelet表示信息，用Wedgelet逼近该节点所对应的图像子块；

3.遍历搜索完整个四叉树后，将所有叶节点对应的图像子块用Wedgelet逼近，最终得到整幅图像的Wedgelet逼近，如图4所示。

步骤4，根据得到的多尺度Wedgelet逼近图像，通过边缘检测梯度提取算子检测逼近图像的边缘，具体实施步骤如下：

1.对得到的多尺度Wedgelet逼近图像，本发明利用缘检测梯度公式Y(i，j)＝max{Sx，Sy}对Wedgelet逼近图像进行边缘检测，得到为当前像素点的梯度值，其中，Sx和Sy分别对应水平和垂直方向的梯度值。该公式的检测模板如图5所示，图5中，左上角的像素点为当前像素点，

2.根据设定的边缘阈值Th，对当前像素点的梯度值进行判断，如果Y(i，j)＞Th，则该像素点的值为1，即该像素点为边缘点，否则为0，即该像素点为非边缘点。

步骤5，根据检测到的边缘图像，去除长度小于某一阈值的细小边缘杂点，得到最终的边缘检测结果，具体实施步骤如下：

1.对步骤4所得到的边缘图像，按照8连通方式扫描边缘图像中的连通曲线，分别标记每条连通曲线；

2.分别统计每条连通曲线的长度，当其长度小于设定阈值Th_post时，在边缘图像中去除该曲线。

本发明的效果可以通过以下仿真结果进一步说明。

仿真条件：

本发明选用了4幅公共测试SAR图像进行仿真实验，计算机硬件配置为Pentium(R)4，主频3.00GHz。现有Wedgelet逼近技术Beamlab200的软件平台为Matlab7.0，本发明的软件平台为Matlab7.0和VC2005，选取的测试图像为Cameraman图像，大小256×256，运行时间取3次运行的平均结果。

仿真结果：

1.本发明的逼近时间与现有技术的逼近时间如表1所示：

表1Beamlab200工具箱和本发明方法运行时间对比

	现有技术Beamlab200方法(运行平台Matlab7.0)	本发明方法(运行平台Matlab7.0)	本发明方法(运行平台VC2005)
				测试图像	Cameraman图像	Cameraman图像	Cameraman图像
Wedgelet逼近时间	1806秒	178秒	15秒

从表1可见，在相同软件平台下，本发明的逼近时间仅为现有技术运行时间的1/10，而在VC2005的软件平台下，运行时间低于现有技术的1/100，可以看出，与现有技术相比，本发明更具有可应用性。

2.本发明的逼近效果如图6所示。其中，图6(a)是SAR图像原图，图6(b)是现有Donoho方法的逼近结果，图6(c)本发明方法的逼近结果。从逼近结果可以看出，本发明方法的逼近结果产生的伪边缘较少，抗噪性较好。

3.本发明的检测结果如图7、图8、图9、图10所示，其中：

图7(a)为河流的SAR图像原图，图7(b)和7(b1)分别为Canny算子后处理前后边缘检测结果，图7(c)和7(c1)分别为Sobel算子后处理前后边缘检测结果，图7(d)和7(d1)分别为Ratio算法后处理前后边缘检测结果，图7(e)和7(e1)分别为小波模极大值方法后处理前后边缘检测结果，图7(f)和7(f1)分别为本发明方法后处理前后边缘检测结果。

图8(a)为机场的SAR图像原图，图8(b)和8(b1)分别为Canny算子后处理前后边缘检测结果，图8(c)和8(c1)分别为Sobel算子后处理前后边缘检测结果，图8(d)和8(d1)分别为Ratio算法后处理前后边缘检测结果，图8(e)和8(e1)分别为小波模极大值方法后处理前后边缘检测结果，图8(f)和8(f1)分别为本发明方法后处理前后边缘检测结果。

图9(a)为机场跑道SAR图像原图，图9(b)和9(b1)分别为Canny算子后处理前后边缘检测结果，图9(c)和9(c1)分别为Sobel算子后处理前后边缘检测结果，图9(d)和9(d1)分别为Ratio算法后处理前后边缘检测结果，图9(e)和9(e1)分别为小波模极大值方法后处理前后边缘检测结果，图9(f)和9(f1)分别为本发明方法后处理前后边缘检测结果。

图10(a)为港口SAR图像原图，图10(b)和10(b1)分别为Canny算子后处理前后边缘检测结果，图10(c)和10(c1)分别为Sobel算子后处理前后边缘检测结果，图10(d)和10(d1)分别为Ratio算法后处理前后边缘检测结果，图10(e)和10(e1)分别为小波模极大值方法后处理前后边缘检测结果，图10(f)和10(f1)分别为本发明方法后处理前后边缘检测结果。

从图7、图8、图9、图10检测结果可以看出，传统的边缘检测方法对SAR图像的相干斑噪声抑制能力不够强且边缘完整连续性不够好，主要表现如下：

Canny方法的检测到的边缘比较精细，对边缘的提取比较完整，但抑噪能力较差，体现在检测到的错误边缘较多。

Sobel方法虽然检测到的伪边缘较Canny算子少，但边缘检出率不及Canny算子，且抑噪能力不理想。这主要是因为经典的边缘检测算子本质上都为高通滤波器，而SAR图像的边缘和噪声都为高频信息，所以该类算子并不能很好地解决SAR图像的边缘检测问题。

Ratio方法对SAR图像有较高的检出率，但检测到的边缘线较粗，难以判别真实边缘的位置。

小波模极大值边缘检测方法相对经典梯度算子方法和Ratio方法，其抑噪能力和检测效果有所改善，但从对其后处理的结果可以看出，其边缘存在许多间断现象。

本发明利用Wedgelet能够较好地捕捉图像中的“线”和“面”的特性，进一步对Wedgelet逼近过程，提高了其抑噪能力和运行速度，并在抑噪能力、边缘检出率和主轮廓完整性方面都较传统方法有明显改善。

从每种方法后处理结果可以看出，传统算法都存在主轮廓边缘不连续现象，而本发明方法经后处理操作后，由于其良好的主轮廓边缘连续性，去除了细小杂边缘的干扰，得到了较好的边缘检测结果。图11(a)、图11(b)、图11(c)、图11(d)分别给出了本发明方法对图7(a)、图8(a)、图9(a)、图10(a)边缘检测结果与原图叠加的结果，从结果可以看出，本发明方法对SAR图像的边缘检测是有效的。

Claims (5)

1.一种基于改进Wedgelet的多尺度SAR图像边缘检测方法，包括如下过程：

(1)输入待检测的SAR图像，将其转化为灰度图像，并将得到的灰度图像分成多个子图像；

(2)用改进的Wedgelet代价函数对每个子图像自适应地按照自下向上的尺度进行多尺度Wedgelet分解，用四叉树保存各个尺度的分解系数，该Wedgelet变换的代价函数为：

其中，P′，P″分别表示图像子块S被划分成的两个楔形区域，

图像块S的大小为n×n，

和

分别表示区域P′，P″上像素的均值，

η为

的权值；

(3)根据分解得到的Wedgelet系数，进行多尺度Wedgelet重构，得到每个子图像的多尺度Wedgelet逼近，对所有子图像都逼近后，得到输入图像的多尺度Wedgelet逼近；

(4)对得到的多尺度Wedgelet逼近图像，通过边缘提取梯度算子检测该逼近图像的边缘；

(5)搜索检测得到的边缘图像，存储每一条边缘的长度，去除边缘长度小于设定阈值的细小边缘，得到最终的边缘检测结果。

2.根据权利要求1所述的多尺度SAR图像边缘检测方法，其中步骤(2)按如下过程进行：

(2a)根据代价函数计算公式建立分解系数四叉树的最底层，得到最底层系数，

其中，P表示图像的ED-RDP，即逼近图像的Wedgelet的集合，

λ为复杂度惩罚参数；

(2b)按自下向上的尺度分别搜索各尺度上每个二进正方形的最优Wedgelet表示，该最优表示分为非退化Wedgelet表示、退化Wedgelet表示和该二进正方形为中间节点三种情况，对于非退化的Wedgelet表示，搜索其最优Wedgelet基的代价函数为：

对于退化的Wedgelet表示，其代价函数为：

对于该二进正方形为中间节点的情况，代价函数为其对应四个子节点的代价函数之和；

(2c)取所述的三个代价函数值中的最小者作为当前二进正方形对应节点的代价函数，存储于分解系数四叉树对应的节点中；

(2d)如果最小代价函数对应退化或非退化的Wedgelet表示，标记当前节点为叶节点，剪枝掉以当前节点为根节点的所有子节点，并在分解系数四叉树中的对应位置存储当前节点的最优Wedgelet表示信息，否则记录当前节点为中间节点；

(2e)循环执行(2b)～(2d)，对所有尺度执行完后停止，得到每个子图像完整的分解系数四叉树。

3.根据权利要求1所述的多尺度SAR图像边缘检测方法，其中步骤(3)按如下过程进行：

(3a)将分解系数四叉树，自顶向下遍历搜索整个四叉树；

(3b)对于搜索到每个叶节点，根据存储在该节点中的其对应二进正方形的Wedgelet表示信息，用Wedgelet逼近该节点所对应的图像子块；

(3c)遍历搜索完整个四叉树后，则所有叶节点对应的图像子块用Wedgelet逼近，最终得到整幅图像的Wedgelet逼近。

4.根据权利要求1所述的多尺度SAR图像边缘检测方法，其中步骤(4)按如下过程进行：

(4a)根据多尺度Wedgelet逼近图像的特点，利用边缘检测梯度公式Y(i，j)＝max{Sx，Sy}对Wedgelet逼近图像进行边缘检测，得到当前像素点的梯度值，其中，Sx和Sy分别对应水平和垂直方向的梯度值；

(4b)根据设定的边缘阈值Th，对当前像素点的梯度值进行判断，如果Y(i，j)＞Th，则该像素点的值为1，即该像素点为边缘点，否则为0，即该像素点为非边缘点。

5.根据权利要求1所述的多尺度SAR图像边缘检测方法，其中步骤(5)按如下过程进行：

(1).对步骤4所得到的边缘图像，按照8连通方式扫描边缘图像中的连通曲线，分别标记每条连通曲线；

(2).分别统计每条连通曲线的长度，当其长度小于设定阈值Th_post时，在边缘图像中去除该曲线。

Images