CN105069796A - 基于小波散射网络的sar图像分割方法 - Google Patents

Abstract

本发明公开了一种基于小波散射网络的SAR图像分割方法，解决了传统基于纹理的SAR图像分割低效耗时的技术问题。实现步骤为：对SAR图像进行去噪和归一化的预处理；设定散射变换路径，选定的小波函数和窗函数，产生散射传播算子和散射算子；对SAR图像进行散射变换，得到SAR图像像素点的散射系数；对降维后散射纹理特征进行K-Means聚类，获得初步分割结果；寻找不同类别像素点所处的仿射空间，构成仿射分类器；利用仿射分类器对初步分割结果进行滑窗更正，完成SAR图像分割。本发明无需对SAR图像进行分块操作，且提取了能够降低相同纹理间差异性、增加不同纹理间差异性的散射纹理特征，能够对SAR图像进行精确、快速的分割，用于单幅SAR图像的快速分割。

Description

基于小波散射网络的SAR图像分割方法

技术领域

本发明属于图像处理技术领域，涉及一种SAR图像分割，具体是一种基于小波散射网络的SAR图像分割方法。可应用于单幅SAR图像的分割。

背景技术

合成孔径雷达(SAR)是一种主动式微波对地观测设备，通过发射和接收雷达波而获得地面信息，具有全天时、全天候工作的特点，且雷达波对地面植被具有一定的穿透能力。SAR图像已广泛应用于军事和国民经济领域，包括飞机目标检测，土地监测，导航，自然灾害监测等。对SAR图像的研究主要包括，SAR图像的去噪，融合，分割，目标检测，变化检测等，其中SAR图像分割是一项基础却又重要的信息提取技术。传统的SAR图像分割算法包括基于像素和基于纹理的分割方法，其中基于像素的分割方法常常需要使用滤波效果好的算法对SAR图像进行去噪，同时此方法对于阈值的选取比较敏感，从而增加了分割方法的计算量；而基于纹理的分割算法则将重点放在对纹理特征的提取上。目前较为普遍的纹理特征提取方法是对SAR图像提取小图像块，然后对小图像块提取纹理特征，其常用的纹理特征包括灰度共生矩阵(GLCM)，小波纹理等。由于这些纹理特征提取算法是基于小图像块，其算法时间复杂度非常大。同时，马尔科夫随机场(MRF)也常被用来描述像素点之间的结构关系。但MRF的迭代优化算法也需要很大的计算量，不利于SAR图像的快速分割。

2011年S.Mallat提出了一种散射变换算子，它能够恢复低频信号中丢失的高频成分，从而更好的保持图像在仿射变换和弹性形变时的不变特性，参见Mallat,Stéphane."GroupInvariantScattering."CommunicationsonPure&AppliedMathematics65.10(2012):1331–1398。2013年S.Mallat通过级联小波变换的非线性取模值操作和均值操作，构成了一种小波散射网络，能够在保持平移不变和对弹性形变稳定的同时保留有助于解决分类问题的高频信息，其在纹理图像的分类问题上取得了非常好的效果，参见Bruna,Joan,andS.Mallat."InvariantScatteringConvolutionNetworks."IEEETransPatternAnalMachIntell35.8(2013):1872-1886。吴华娟等人直接将此小波散射网络应用于纹理分割问题中，其做法是对纹理图像进行取块操作，然后对每个小图像块应用小波散射卷积网络提取散射系数，作为此图像块的纹理特征，参见吴华娟,张明新,and郑金龙."基于小波散射卷积网络的纹理分割方法."微电子学与计算机30.5(2013)。与之前基于纹理的SAR图像分割算法一样，吴等人的做法需要大量的计算量，从而降低SAR图像分割算法的效率。

以上各种基于纹理的SAR图像分割算法都是在图像块上进行纹理特征的提取操作，其算法复杂度随着图像大小的增大而快速增加。对于SAR图像中的每一个像素点，上述算法都需要取定一个以此像素点为中心的像素块，然后在此像素块上进行纹理特征提取。对于距离小于像素块边长的两个像素点，以这两个像素点为中心的两个像素块之间具有重叠的像素点，而这些重叠的像素点使得对图像块进行纹理特征提取时产生重复计算。两个像素点之间的距离越小，对应的两个图像块之间的重叠像素点越多，纹理特征提取时重复计算就越多，算法效率就越低。

发明内容

针对上述已有技术存在的缺陷，本发明目的在于提出一种快速提取SAR图像纹理特征并进行分割的技术方案，即基于小波散射网络的SAR图像分割方法。该方法实现了SAR图像的高效分割。

本发明是一种基于小波散射网络的SAR图像分割方法，是基于SAR图像纹理的快速特征提取，包括有如下步骤：

1)选择需要进行分割的原始SAR图像，大小为[m,n]，并对其进行预处理，包括Lee滤波去噪和数据归一化操作。

2)设定散射变换最大尺度J，在对SAR图像进行分割时，每个像素点的特征由一个大小为[2*2^J,2*2^J]的图像块确定；当处理SAR图像边界的像素点时，需要对SAR图像边界进行对称拓展，得到拓展后的SAR图像f，大小为[m+2*2^J,n+2*2^J]。

3)选定散射变换中的小波函数ψ_λ(x)，其中x为需要进行小波变换的图像，λ为参数；根据拓展后图像f的大小和散射变换尺度J，设定散射变换路径p＝{λ₀,λ₁,...,λ_L-1}，其中L为路径总长度，并生成二维散射传播算子U(p)f。

4)选取窗函数，结合步骤3)中的二维散射传播算子U(p)f，共同生成二维散射算子S_J(p)f。

5)根据二维散射算子S_J(p)f和路径p，对拓展后的SAR图像f进行路径长度为L的散射变换，级联各层输出系数，得到拓展后的SAR图像f的散射系数图。

6)对SAR图像的散射系数图的上下左右边界进行剪裁，使得散射系数图的大小等于原SAR图像的大小，即[m,n]，于是，SAR图像中每个像素点的散射系数s_ori对应一个长度为N_ori的向量。

7)用离散余弦变换(DCT)对SAR图像的散射系数进行降维，保留前50％的低频信息作为降维后的散射纹理特征s：

7a)对SAR图像中所有像素点的散射系数进行DCT变换，得到了DCT系数；

7b)分别选取DCT系数中前50％的低频系数，作为对应每个像素点散射系数降维后的散射纹理特征，记为s，长度为N；也就是说，用每个像素点降维后的散射纹理特征来表示此像素点。

8)设定SAR图像要被分割成的类别数C，对降维后像素点的散射纹理特征s进行K-Means聚类，得到SAR图像的初步分割结果Y＝{y_i|y_i＝1,2,...C；i＝1,2,...,m*n}，其中y_i对应第i个像素点的类标；至此，以上所有步骤均无需对图像进行取块操作，而实现SAR图像的初步分割；极大地节省了SAR图像分割的运行时间。

9)由初步分割结果Y生成有标签样本集合其中K为有标签样本的总数，表示有标签样本集合T中第k个有类标像素点的散射纹理特征，表示有标签样本集合T中第k个有类标像素点的类标：

9a)根据SAR图像的初步分割结果，从SAR图像中选取所有属于第c类的像素点；

9b)对步骤9a)选出的属于第c类的像素点，求取其平均值作为第c类的类中心o_c；

9c)对所有初步分割结果中类标为c的像素点，计算它们到o_c的距离，选取前20％距离最近的像素点作为属于第c类置信度较高的像素点；

9d)对于C种不同的类别，分别重复步骤9a)～9c)，得到所有对应类别置信度较高的像素点，构成含有C种不同类别的有标签样本集T。

10)在有标签样本集合T上进行训练，寻找不同类别像素点所处的仿射空间A_c，其中c＝1,2,...,C，仿射空间A_c作为第c类像素点所处的子空间，由C个仿射空间组成一个仿射分类器。

11)逐步滑窗对SAR图像初步分割结果Y进行遍历，使用步骤10)中的仿射分类器对初步分割结果Y进行更正，完成SAR图像分割，并输出最终分割结果。

本发明的纹理提取操作是作用于整幅SAR图像的，省去了对SAR图像分块提取纹理特征的步骤，其算法复杂度受图像大小的影响远小于传统算法，能够大大提高SAR图像的分割效率。同时本发明在保持平移不变性和对弹性形变稳定性同时保留了有助于分类的高频信息，使得本发明能够在提高分割效率的同时保证足够的分割精度。

本发明与现有的技术相比具有以下优点：

1、本发明针对传统SAR图像纹理特征提取需进行分块处理的缺点，提出了一种无需进行分块操作的SAR图像纹理提取方法，极大地提高了提取纹理特征的效率；

2、本发明采用一种有效的散射变换算子，其在保持纹理平移不变性和对弹性形变稳定性的同时保留了有助于分类的高频信息，在有效地减小相同纹理间差异性的同时增加不同纹理间的差异性，从而提高SAR图像的分割精度；

3、本发明在SAR图像的初步分割基础上构造了一个简单有效的仿射分类器，对初步分割结果进行更正，避免了人工选取有标签训练样本的过程，同时进一步提高了SAR图像的分割精度。

仿真结果表明，本发明在保证SAR图像分割精度的同时极大地提高了SAR图像分割的效率。

附图说明

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

图2是本发明中待分割的SAR图像；

图3是本发明与传统分割算法的整体分割效果比较图。

具体实施方式

随着合成孔径雷达的快速发展，SAR图像已广泛应用于军事和国民经济领域，对海量SAR图像的研究越来越重要。对SAR图像的研究主要包括，SAR图像的去噪，融合，分割，目标检测，变化检测等，其中SAR图像分割是一项基础却又重要的信息提取技术。现有基于纹理特征的SAR图像分割技术大多是基于图像块进行的，由此造成的重复计算大大降低了SAR图像的分割效率。

针对该技术问题，本发明展开了研究与探讨，提供了一种基于小波散射网络的SAR图像分割方法。

实施例1

本发明的基于小波散射网络的SAR图像分割方法是基于SAR图像纹理的快速特征提取。参见图1，包括有如下步骤：

1)选择需要进行分割的原始SAR图像，大小为[256,256]，参见图2；选定Lee滤波窗口大小5*5对其进行Lee滤波去噪，同时对SAR图像数据除以255，作为数据归一化操作。

2)设定散射变换最大尺度J＝2，在对SAR图像进行分割时，每个像素点的特征由一个大小为[8,8]的图像块确定；当处理SAR图像边界的像素点时，需要对SAR图像边界进行对称拓展，得到拓展后的SAR图像f，大小为[264,264]。

3)选定散射变换中的Morlet小波函数其中为二维母小波函数，x为上述拓展后的SAR图像f，r为旋转参数，分别选取4个旋转方向0°，45°，90°和135°，j为尺度参数，取0～J之间的整数，ψ_λ的频率中心为λ＝2^-jr；根据拓展后SAR图像f的大小和散射变换尺度J，设定散射变换路径其中L＝3为路径总长度；同时由上述Morlet小波函数，产生一组二维滤波器其中 ${\mathrm{\λ}}_0=r,{\mathrm{\λ}}_1=2^{-1}r,{\mathrm{\λ}}_2=2^{-2}r.$

4)取定窗函数其中为二维高斯函数，按照路径p＝{r,2^-1r,2^-2r}和步骤3)中的散射传播算子U(p)f，产生散射变换的框架结构scat，即在每层传播参数之后增加一个窗函数滤波，以保持窗内纹理的空间不变性；对SAR图像进行散射变换的二维散射算子为换句话说，每层散射变换的输出是由小波滤波、取模值和窗函数滤波三个操作产生的。

5)按照步骤4)散射变换的框架结构scat，首先用第一层的小波函数对SAR图像进行滤波，取模值，得到散射变换第一层中的低频成分再用窗函数滤波，得到散射变换的第一层输出用第二层的小波函数对进行滤波，取模值，得到散射变换第二层中的低频成分再用窗函数滤波，得到散射变换的第二层输出重复上述操作，直至完成整个路径p；级联散射变换的各层输出，作为散射变换的散射系数s_ori，即其长度为N_ori＝64。

6)由步骤5)中获得SAR图像的散射变换系数，得到SAR图像的三维纹理特征图，大小为[264,264,64]，其中[264,264]为步骤2)中SAR图像边界拓展后的大小，64为散射系数s_ori的长度；将此三维纹理特征图的边界进行裁剪，上下左右各裁剪2^J＝4行/列的像素点，得到与原SAR图像像素点一一对应的散射系数，即三维纹理特征图的大小被裁剪为[264,264,64]；于是，SAR图像中的每个像素点的散射系数为一个长度为64的向量。

7)对SAR图像中的所有像素点的散射系数s_ori进行余弦离散变换(DCT)，得到变换系数s_DCT；在变换系数s_DCT中选取对应于前50％低频成分的部分，作为像素点降维后的特征散射纹理特征，记为s，长度为N＝32，即用对应的散射纹理特征来表示此像素点。

8)取定SAR图像将要被分割成的类别数C＝3，对步骤7)中获得的降维后的像素点的散射纹理特征s进行K-Means聚类，聚为同一类的像素点被标记为同一种类标；于是得到SAR图像的初步分割结果Y＝{(s_i,y_i)|y_i＝1,2,...C；i＝1,2,...,m*n}，其中y_i对应第i个像素点s_i的类标；至此，以上所有步骤均无需对图像进行取块操作，而实现SAR图像的初步分割；极大地节省了SAR图像分割的运行时间。

9)根据SAR图像的初步分割结果Y，选取SAR图像中类标为的像素点其中c＝1,2,...,C表示类标，表示已选出的有类标像素点的个数；对选出的类标为c的样本集取其所有样本平均值作为此类样本的中心，记为o_c，即其中K_c表示S_c中样本的数目；计算S_c中每个样本到o_c的距离，并按距离的大小升序排列，记为其中对应S_c中第n_k个样本到o_c的距离；选取前20％距离最小的像素点作为属于第c类置信度较高的像素点(即对应的像素点其中符号表示向下取整)，组成大小为k_c的有标签样本集对于C种不同的类别，重复步骤9a)～9c)，得到不同类别中置信度较高的像素点，从而构成有标签样本集 $T=\left\{(s_k^l,y_k^l)\right|y_k^l=1,2,...,C;k=1,2,...,K\},$ 即T＝T₁∪T₂∪...∪T_C。

10)对于步骤9)中得到的有标签集合计算其中样本的协方差矩阵D_c(s)，计算公式如下其中表示T_c中样本的平均值，即其中k_c为集合T_c的大小；对D_c(s)进行特征值分解，得到特征值和对应的特征向量{v₁,v₂,...,v_N}，其中特征值按降序排列；选取前个最大特征值所对应的特征向量{v₁,v₂,...,v_d}，由其线性组合生成d＝6维子空间V_c，即都有下式成立v＝a₁v₁+a₂v₂+...+a_dv_d，其中a₁,a₂,...,a_N为线性组合系数；取有标签样本集T_c中样本的均值E_c(s)，即由线性空间V_c和均值E_c(s)，组合得到第c类样本所处的仿射空间A_c，即A_c＝E_c(s)+V_c，其中都有下式成立其中对于C种不同的类别，重复上述步骤，获得不同类别样本所处仿射空间A_c,c＝1,2,...,C，构成一个仿射分类器。

11)在SAR图像上按步长step＝1滑动取窗，窗口大小为5*5，统计窗内像素点在初步分割结果中的类别，其不同的类标组成集合Y_win；根据Y_win选取对应的A_c，其中c∈Y_win，产生对应于此窗口的仿射分类器，用此仿射分类器对初步分割结果进行更正；滑完整幅图像，完成图像分割过程。

本发明对SAR图像进行纹理特征提取时，省去了对图像的分块操作，从而大大提高了分割效率。同时，本发明采用一种有效的散射变换算子，在保持纹理平移不变性和对弹性形变稳定性的同时保留了有助于分类的高频信息，提高了SAR图像的分割精度。再者，本发明在SAR图像的初步分割基础上构造了一个简单有效的仿射分类器，对初步分割结果进行更正，避免了人工选取有标签训练样本的过程，进一步提高了SAR图像的分割精度。

实施例2

基于小波散射网络的SAR图像分割方法同实施例1，其中步骤3)所述的生成二维散射传播算子U(p)f，具体步骤如下：

3a)选取Morlet二维母小波函数于是散射变换中的小波函数组为ψ_λ(x)＝2^-2jψ(2^-jr^-1x)_，其中r为旋转参数，分别选取4个旋转方向0°，45°，90°和135°，j为尺度参数，取0～J之间的整数，ψ_λ的频率中心为λ＝2^-jr。

3b)根据拓展后SAR图像大小[264,264]和散射变换尺度J＝2设定散射变换路径p＝{λ₀,λ₁,...,λ_L-1}＝{r,2^-1r,2^-2r}，再由步骤3a)中的Morlet小波函数，产生一组二维滤波器其中λ₀＝r,λ₁＝2^-1r,λ₂＝2^-2r。

3c)首先用步骤3b)中的第一个滤波器对图像f进行滤波，得到对应的小波系数然后对其做取模值操作，即作为散射传播到第一层的参数；接着再用下一个滤波器对进行滤波，从而恢复出其中的高频成分，最后再做取模值操作，即作为散射传播到第二层的参数取第三个滤波器重复上述操作得到散射传播到第三层的参数。换句话说，根据二维滤波器组filters构成散射传播算子为 $U\left(p\right)f=\left|\right|\left|\right|f*{\mathrm{\ψ}}_{{\mathrm{\λ}}_0}|*{\mathrm{\ψ}}_{{\mathrm{\λ}}_1}|...|*{\mathrm{\ψ}}_{{\mathrm{\λ}}_l}|,l\≤L-1.$

本实施例中产生的散射传播算子U(p)f，能够将图像f纹理特征的高频成分传播到各层参数之中，增加了不同纹理间的差异性。

实施例3

基于小波散射网络的SAR图像分割方法同实施例1-2，其中步骤4)所述的产生二维散射算子S_J(p)f，具体步骤如下：

4a)取定窗函数其中为二维高斯函数，即

4b)按照路径p＝{r,2^-1r,2^-2r}和步骤3)中的散射传播算子U(p)f，产生散射变换的框架结构scat，即在每层传播参数之后增加一个窗函数滤波，以保持窗内纹理的空间不变性。

4c)对SAR图像进行散射变换的二维散射算子为其中路径长度为L＝3。换句话说，每层散射变换的输出是由小波滤波、取模值和窗函数滤波三个操作产生的。

本实施例中产生的二维散射算子在保持了纹理特征平移不变性和对弹性形变稳定性的同时保留了有助于分类的高频信息，能够在有效地减小相同纹理间差异性的同时增加不同纹理间的差异性，从而提高SAR图像的分割精度。

实施例4

基于小波散射网络的SAR图像分割方法同实施例1-3，其中步骤10)中所述的寻找不同类别样本的仿射空间组成仿射分类器，具体步骤如下：

10a)对于步骤9)中得到的有标签集合计算其中样本的协方差矩阵D_c(s)，计算公式如下：

$D_c\left(s\right)=\frac{\underset{k=1}{\overset{k_c}{\Σ}}(s_k^l-\stackrel{\‾}{s_c})}{k_c-1}$

其中表示T_c中样本的平均值，即其中k_c为集合T_c的大小。

10b)对D_c(s)进行特征值分解，得到特征值和对应的特征向量{v₁,v₂,...,v_N}，其中特征值按降序排列；选取前个最大特征值所对应的特征向量{v₁,v₂,...,v_d}，由其线性组合生成d维子空间V_c，即都有下式成立

v＝a₁v₁+a₂v₂+...+a_dv_d，

其中a₁,a₂,...,a_N为线性组合系数。

10c)取步骤10a)中有标签样本集T_c中样本的均值E_c(s)，即

10d)由线性空间V_c和均值E_c(s)，组合得到第c类样本所处的仿射空间A_c，即

A_c＝E_c(s)+V_c，

其中 $\∀s\∈A_c,$ 都有下式成立 $s=E_c\left(s\right)+s_{V_c},$ 其中 $s_{V_c}\∈V_c.$

10e)对于C种不同的类别，重复步骤10a)～10d)，获得不同类别样本所处仿射空间A_c,c＝1,2,...,C，构成一个仿射分类器。

本实施例在SAR图像的初步分割基础上构造了一个简单有效的仿射分类器，对初步分割结果进行更正，避免了人工选取有标签训练样本的过程，同时进一步提高了SAR图像的分割精度。

实施例5

基于小波散射网络的SAR图像分割方法同实施例1-4，其中步骤10e)所述的仿射分类器，工作方式如下：

给定一个像素点的散射纹理特征s，将其在不同类别的仿射空间A_c,c＝1,2,...,C进行投影，得到计算s和各个投影的差值在所有求使得差值最小的c作为s的类标，记为即以上就是此分类器的工作方式。

实施例6

基于小波散射网络的SAR图像分割方法同实施例1-5，其中步骤11)所述的对初步SAR图像的初步分割结果Y进行更正，具体步骤如下：

11a)取窗口大小为5*5，窗口滑动的步长为step＝1。

11b)在SAR图像上取大小为5*5的窗口，统计窗内像素点在初步分割结果中的类别，其不同的类标组成集合Y_win，组成仿射分类器对窗口内部的像素点进行更正。

11c)按照步骤11a)中设定的步长，逐步在SAR图像上滑动窗口，重复步骤11a)和11b)，直至窗口滑完整幅图像，完成图像分割过程。

本实施例利用了一种简单有效的仿射分类器，在结合了像素点空间信息的基础上对步骤8)中获得的初步分割结果进行了更正，进一步提高分割精度。

实施例7

基于小波散射网络的SAR图像分割方法同实施例1-6，其中步骤11b)所述的组成仿射分类器对窗口内部的像素点进行更正，工作方式如下：

11b1)若窗内所有像素点的类别相同，即Y_win中只有一种类标，则不做任何操作；

11b2)若窗内位于四周的像素点拥有相同的类别c，则将窗内所有像素点的类别数更正为c；

11b3)若窗内像素点有两种不同的类别，即Y_win中有两种不同的类标，则认为此窗口对应图像分割的边缘，不做任何操作；

11b4)若窗内有三种或三种以上不同类别的像素点，即Y_win中有三种或三种以上不同的类标，则用步骤10)中的仿射分类器对窗内像素点进行重新分类，更新此窗口内像素点的类标；设窗口内有设窗口内有M种不同的类别，记为Y_win＝{c₁,c₂,...,c_M}，其中M≥3，则此时的仿射分类器仅由与此M个对应的仿射空间构成，即给定像素点s，其分类结果由下式得到：

$\hat{c}=\underset{c\∈Y_{win}}{\arg min}\left|\right|s-P_{A_c}\left(s\right)\left|\right|.$

下面给出一个详尽的实施流程过程：

实施例8

基于小波散射网络的SAR图像分割方法同实施例1-7，参照图1，本发明的具体实现步骤如下：

步骤1，选择待分割的原始SAR图像，并对其进行包括Lee滤波去噪和数据归一化在内的预处理。

该步骤的具体实现如下：

1a)选取待分割的SAR图像，大小为[256,256]；

1b)选择Lee滤波窗口大小为5*5，根据此窗口大小对SAR图像的边界进行对称拓展，并对其进行Lee滤波，产生跟原图大小相同的去噪后的SAR图像；

1c)通过对SAR图像数据除以255的操作，将在图像数据归一化到0～1之间。

步骤2，对预处理之后的SAR图像边界进行拓展，为后续的散射变换做准备。

该步骤的具体实现如下：

确定散射变换的最大尺度J＝2，则SAR图像中每个像素点的特征由其周围大小为[8,8]的图像块决定；据此对SAR图像的边界进行对称拓展，拓展后的图像记为f，大小为[264,264]。

步骤3，设定散射变换的小波函数和变换路径，生成散射变换所需的二维散射传播算子U(p)f。

该步骤的具体实现如下：

3a)选取Morlet二维母小波函数于是散射变换中的小波函数组为ψ_λ(x)＝2^-2jψ(2^-jr^-1x)，其中r为旋转参数，分别选取4个旋转方向0°，45°，90°和135°，j为尺度参数，取0～J之间的整数，ψ_λ的频率中心为λ＝2^-jr；

3b)根据图像大小和散射变换尺度J设定散射变换路径p＝{λ₀,λ₁,...,λ_L-1}＝{r,2^-1r,2^-2r}，其中L＝3，再由步骤3a)中的Morlet小波函数，产生一组二维滤波器其中λ₀＝r,λ₁＝2^-1r,λ₂＝2^-2r；

3c)首先用步骤3b)中的第一个滤波器对图像f进行滤波，得到对应的小波系数然后对其做取模值操作，即作为散射传播到第一层的参数；接着再用下一个滤波器对进行滤波，从而恢复出其中的高频成分，最后再做取模值操作，即作为散射传播到第二层的参数取第三个滤波器重复上述操作得到散射传播到第三层的参数。换句话说，根据二维滤波器组filters构成散射传播算子为其中L＝3，将图像f纹理特征的高频成分传播到各层参数之中。

步骤4，产生二维散射算子。

该步骤的具体实现如下：

4a)取定窗函数其中为二维高斯函数，即

4b)按照路径p＝{r,2^-1r,2^-2r}和步骤3)中的散射传播算子U(p)f，产生散射变换的框架结构scat，即在每层传播参数之后增加一个窗函数滤波，以保持窗内纹理的空间不变性；

4c)对SAR图像进行散射变换的二维散射算子为换句话说，每层散射变换的输出是由小波滤波、取模值和窗函数滤波三个操作产生的。

步骤5，计算SAR图像的散射系数s_ori。

该步骤的具体实现如下：

按照步骤4)散射变换的框架结构scat，首先用第一层的小波函数对SAR图像进行滤波，取模值，得到散射变换第一层中的低频成分再用窗函数滤波，得到散射变换的第一层输出用第二层的小波函数对进行滤波，取模值，得到散射变换第二层中的低频成分再用窗函数滤波，得到散射变换的第二层输出重复上述操作，直至完成整个路径p；级联散射变换的各层输出，作为散射变换的散射系数s_ori，即其长度为N_ori＝64。

步骤6，对SAR图像散射系数图的边界进行裁剪。

该步骤的具体实现如下：

由步骤5)中获得SAR图像的散射变换系数，得到SAR图像的三维纹理特征图，大小为[264,264,64]，其中[264,264]为步骤2)中SAR图像边界拓展后的大小，64为散射系数s_ori的长度；将此三维纹理特征图的边界进行裁剪，上下左右各裁剪2^J＝4行/列的像素点，得到与原SAR图像像素点一一对应的散射系数，即三维纹理特征图的大小被裁剪为[264,264,64]；于是，SAR图像中的每个像素点的散射系数为一个长度为64的向量。

步骤7，对SAR图像中像素点的散射系数进行数据将维。

该步骤的具体实现如下：

7a)对SAR图像中的所有像素点的散射系数s_ori进行余弦离散变换(DCT)，得到变换系数s_DCT；

7b)在变换系数s_DCT中选取对应于前50％低频成分的部分，作为像素点散射系数降维后的散射纹理特征，记为s，长度为N＝32，用对应的散射纹理特征来表示此像素点。

步骤8，对SAR图像进行初步分割。

该步骤的具体实现如下：

取定SAR图像将要被分割成的类别数C＝3，对步骤7)中获得的降维后的像素点散射系数s进行K-Means聚类，聚为同一类的像素点被标记为同一种类标；于是得到SAR图像的初步分割结果Y＝{(s_i,y_i)|y_i＝1,2,...C；i＝1,2,...,m*n}，其中y_i对应第i个像素点s_i的类标。

步骤9，产生有标签训练样本集T。

该步骤的具体实现如下：

9a)根据SAR图像的初步分割结果Y，选取SAR图像中类标为的像素点其中c＝1,2,...,C表示类标，表示已选出的有类标像素点的个数；

9b)对选出的类标为c的样本集取其所有样本平均值作为此类样本的中心，记为o_c，即

$o_c=\frac{1}{K_c}\underset{(s_k^l,y_k^l)\∈S_c}{\Σ}{s}_k^l$

其中K_c表示S_c中样本的数目；

9c)计算S_c中每个样本到o_c的距离，并按距离的大小升序排列，记为其中对应S_c中第n_k个样本到o_c的距离；选取前20％距离最小的像素点作为属于第c类置信度较高的像素点(即对应的像素点其中符号表示向下取整)，组成大小为k_c的有标签样本集 $T_c=\{(s_{n_1},y_{n_1}),(s_{n_2},y_{n_2}),...,(s_{n_{k_c}},y_{n_{k_c}})\};$

9d)对于C种不同的类别，重复步骤9a)～9c)，得到不同类别中置信度较高的像素点，从而构成有标签样本集即T＝T₁∪T₂∪...∪T_C；

步骤10，根据由标签样本集T，生成一个仿射分类器。

该步骤的具体实现如下：

10a)对于步骤9)中得到的有标签集合计算其中样本的协方差矩阵D_c(s)，计算公式如下

$D_c\left(s\right)=\frac{\underset{k=1}{\overset{k_c}{\Σ}}(s_k^l-\stackrel{\‾}{s_c})}{k_c-1}$

其中表示T_c中样本的平均值，即其中k_c为集合T_c的大小；

10b)对D_c(s)进行特征值分解，得到特征值和对应的特征向量{v₁,v₂,...,v_N}，其中特征值按降序排列；选取前个最大特征值所对应的特征向量{v₁,v₂,...,v_d}，由其线性组合生成d维子空间V_c，即都有下式成立

v＝a₁v₁+a₂v₂+...+a_dv_d，

其中a₁,a₂,...,a_N为线性组合系数；

10c)取步骤10a)中有标签样本集T_c中样本的均值E_c(s)，即

10d)由线性空间V_c和均值E_c(s)，组合得到第c类样本所处的仿射空间A_c，即

A_c＝E_c(s)+V_c，

其中 $\∀s\∈A_c,$ 都有下式成立 $s=E_c\left(s\right)+s_{V_c},$ 其中 $s_{V_c}\∈V_c;$

10e)对于C种不同的类别，重复步骤10a)～10d)，获得不同类别样本所处仿射空间A_c,c＝1,2,...,C，构成一个仿射分类器，此分类器的工作方式如下：

给定一个像素点的散射系数s，将其在不同类别的仿射空间A_c,c＝1,2,...,C进行投影，得到计算s和各个投影的差值在所有 $\left\{\right||s-P_{A_c}\left(s\right)||,c=\mathrm{1,2},...,C\}$ 求使得差值最小的c作为s的类标，记为即

步骤11，通过滑窗的方式对SAR图像的分割结果进行更正。

该步骤的具体实现如下：

11a)取窗口大小为5*5，窗口滑动的步长为step＝1；

11b)在SAR图像上取大小为5*5的窗口，统计窗内像素点在初步分割结果中的类别，其不同的类标组成集合Y_win，于是根据如下规则进行更正：

11b1)若窗内所有像素点的类别相同，即Y_win中只有一种类标，则不做任何操作；

11b2)若窗内位于四周的像素点拥有相同的类别c，则将窗内所有像素点的类别数更正为c；

11b3)若窗内像素点有两种不同的类别，即Y_win中有两种不同的类标，则认为此窗口对应图像分割的边缘，不做任何操作；

11b4)若窗内有三种或三种以上不同类别的像素点，即Y_win中有三种或三种以上不同的类标，则用步骤10)中的仿射分类器对窗内像素点进行重新分类，更新此窗口内像素点的类标；设窗口内有设窗口内有M种不同的类别，记为Y_win＝{c₁,c₂,...,c_M}，其中M≥3，则此时的仿射分类器仅由与此M个对应的仿射空间构成，即给定像素点s，其分类结果由下式得到：

$\hat{c}=\underset{c\∈Y_{win}}{\arg min}\left|\right|s-P_{A_c}\left(s\right)\left|\right|;$

11c)按照步骤11a)中设定的步长，逐步在SAR图像上滑动窗口，重复步骤11a)和11b)，直至窗口滑完整幅图像，完成图像分割过程。

本发明的效果可通过仿真进一步说明：

实施例9

基于小波散射网络的SAR图像分割方法同实施例1-8，其仿真内容如下。

1.仿真条件

本实例仿真条件如下：COREi53.2GHzPCWindows7操作系统，Matlab2010运行平台。

2.仿真内容及结果

选择一副大小为256*256的SAR图像，参见图2，用本发明和传统基于纹理的方法分别对此SAR图像进行分割仿真；所用的传统纹理特征提取方法有：灰度共生矩阵(GLCM)和小波纹理。

所选SAR图像如图2所示，仿真结果如图3所示，其中图3(a)表示基于灰度共生矩阵纹理的分割结果，图3(b)表示基于小波纹理的分割结果，图3(c)表示本发明的分割结果。在各个分割结果中，灰色表示道路，灰色表示农田，白色表示城区。从图3(a)中可以看出，基于GLCM纹理的SAR图像分割方法不能有效的区分道路和农田，其分割结果非常杂乱，各个类别之间的区分度很差；由图3(b)可知，基于小波纹理的SAR图像分割方法虽然能够较好的区分道路和农田，但仍未能保持道路和农田内部的区域一致性；参见图3(c)，本发明的分割结果不仅能够较好的区分不同的区域(道路、农田和城区)，同时能够保持各类别内部的区域一致性。宏观来说，本发明分割结果各类别的区分度和区域一致性均优于传统基于纹理的分割结果。

从分割效率角度来说，本发明也优于传统基于纹理的分割结果，具体参见表1。

表1.各个分割算法的运行时间

方法	运行时间(s)
		灰度共生矩阵	205.42
小波纹理	190.15
		本发明	5.76

表1列出了各个方法的运行时间。表中数据说明，本发明在保证分割精度的同时，节省了大量的运行时间。

综上，针对SAR图像分割问题，本发明不论是分割精度还是分割效率，均优于传统基于纹理的分割方法。尤其是本发明的分割效率，在本仿真中大小为256*256的SAR图像分割中，提高了40倍左右。随着合成孔径雷达的快速发展，需要对海量的SAR图像进行高效的分割，而本发明分割方法的高效高精度性使之成为可能。

Claims (7)

1.一种基于小波散射网络的SAR图像分割方法，其特征在于，是基于SAR图像纹理的快速特征提取，包括有如下步骤：

1)选择需要进行分割的原始SAR图像，大小为[m,n]，并对其进行预处理，包括Lee滤波去噪和数据归一化操作；

2)设定散射变换最大尺度J，在对SAR图像进行分割时，每个像素点的特征由一个大小为[2*2^J,2*2^J]的图像块确定；当处理SAR图像边界的像素点时，需要对SAR图像边界进行对称拓展，得到拓展后的SAR图像f，大小为[m+2*2^J,n+2*2^J]；

3)选定散射变换中的小波函数ψ_λ(x)，其中x为需要进行小波变换的图像，λ为参数；根据拓展后图像f的大小和散射变换尺度J，设定散射变换路径p＝{λ₀,λ₁,...,λ_L-1}，其中L为路径总长度，并生成二维散射传播算子U(p)f；

4)选取窗函数结合二维散射传播算子U(p)f，共同生成二维散射算子S_J(p)f；

5)根据二维散射算子S_J(p)f和路径p，对拓展后的SAR图像f进行路径长度为L的散射变换，级联各层输出系数，得到拓展后的SAR图像f的散射系数图；

6)对SAR图像的散射系数图的上下左右边界进行剪裁，使得散射系数图的大小等于原SAR图像的大小，即[m,n]，于是，SAR图像中每个像素点的散射系数s_ori对应一个长度为N_ori的向量；

7)用离散余弦变换(DCT)对SAR图像的散射系数进行降维，保留前50％的低频信息作为降维后的散射纹理特征s：

7a)对SAR图像中所有像素点的散射系数进行DCT变换，得到了DCT系数；

7b)分别选取DCT系数中前50％的低频系数，作为对应每个像素点散射系数降维后的散射纹理特征，记为s，长度为N；

8)设定SAR图像要被分割成的类别数C，对降维后像素点的散射纹理特征s进行K-Means聚类，得到SAR图像的初步分割结果Y＝{y_i|y_i＝1,2,...C；i＝1,2,...,m*n}，其中y_i对应第i个像素点的类标；至此，以上所有步骤均无需对图像进行取块操作，而实现SAR图像的初步分割；

9)根据初步分割结果Y，生成有标签样本集合其中K为有标签样本的总数，表示有标签样本集合T中第k个有类标像素点的散射纹理特征，表示有标签样本集合T中第k个有类标像素点的类标：

9a)根据SAR图像的初步分割结果，从SAR图像中选取所有属于第c类的像素点；

9b)对步骤9a)选出的属于第c类的像素点，求取其平均值作为第c类的类中心o_c；

9c)对所有初步分割结果中类标为c的像素点，计算它们到o_c的距离，选取前20％距离最近的像素点作为属于第c类置信度较高的像素点；

9d)对于C种不同的类别，分别重复步骤9a)～9c)，得到所有对应类别置信度较高的像素点，构成含有C种不同类别的有标签样本集T；

10)在有标签样本集合T上进行训练，寻找不同类别像素点所处的仿射空间A_c，其中c＝1,2,...,C，仿射空间A_c作为第c类像素点所处的子空间，由C个仿射空间组成一个仿射分类器；

11)逐步滑窗对SAR图像初步分割结果Y进行遍历，使用仿射分类器对初步分割结果Y进行更正，完成SAR图像分割，并输出最终分割结果。

2.根据权利要求1所述的基于小波散射网络的SAR图像分割方法，其特征在于，其中步骤3)所述的生成二维散射传播算子U(p)f，按如下步骤进行：

3a)选定Morlet小波ψ_λ(x)＝2^-2jψ(2^-jr^-1x),λ＝2^-jr作为散射变换中的小波函数，其中r为旋转参数，j为尺度参数；

3b)根据拓展后SAR图像f的大小和散射变换尺度J设定散射变换路径p＝{λ₀,λ₁,...,λ_L-1}，其路径总长度为L，同时根据由步骤3a)中选定的Morlet小波，产生二维滤波器组 $filters=\{{\mathrm{\ψ}}_{{\mathrm{\λ}}_0},{\mathrm{\ψ}}_{{\mathrm{\λ}}_1},...,{\mathrm{\ψ}}_{{\mathrm{\λ}}_{L-1}}\};$

3c)由二维滤波器组filters构成散射传播算子

3.根据权利要求1所述的基于小波散射网络的SAR图像分割方法，其中步骤4)中所述的生成二维散射算子S_J(p)f，具体步骤如下：

4a)取定窗函数其中为二维高斯函数；

4b)根据拓展后SAR图像f的大小和路径p，产生散射变换的框架结构scat；

4c)生成二维散射算子

4.根据权利要求1所述的基于小波散射网络的SAR图像分割方法，其中步骤10)中所述的寻找不同类别样本的仿射空间组成仿射分类器，具体步骤如下：

10a)选择有标签训练样本集T中类标为c的像素点，求取此类像素点散射纹理特征s的协方差矩阵D_c(s)；

10b)对D_c(s)进行特征值分解，选取前d个最大特征值所对应的特征向量，生成d维线性子空间V_c，其中d＜N；

10c)求取步骤10b)中所选像素点散射纹理特征s的平均值，记为E_c(s)；

10d)由线性子空间V_c和平均值E_c(s)，得到第c类像素点所处的仿射空间A_c，即A_c＝E_c(s)+V_c；

10e)对于C种不同的类别，重复步骤10a)～10d)，获得不同类别像素点所处仿射空间A_c，其中c＝1,2,...,C，构成一个仿射分类器。

5.根据权利要求4所述的基于小波散射网络的SAR图像分割方法，其中步骤10e)中所述的由仿射空间A_c,c＝1,2,...,C构成一个仿射分类器，通过如下方式对像素点进行分类：

给定一个像素点，用其散射纹理特征s来表示该像素点，此像素点的类标通过如下公式确定：

$\hat{c}=\underset{1\≤c\≤C}{\arg min}\left|\right|s-P_{A_c}\left(s\right)\left|\right|,$

其中，表示像素点的散射系数s在仿射空间A_c上的投影。

6.根据权利要求1所述的基于小波散射网络的SAR图像分割方法，其中步骤11)中所述的对初步SAR图像的初步分割结果Y进行更正，具体步骤如下：

11a)设定窗口大小win*win和滑窗步长step；

11b)在SAR图像上取大小为win*win的窗口，统计窗内像素点在初步分割结果中的类别；若窗内所有像素点的类别相同，则不做任何操作；若窗内位于四周的像素点拥有相同的类别c，则将窗内所有像素点的类别数更正为c；若窗内像素点有两种不同的类别，则认为此窗口对应图像分割的边缘，不做任何操作；若窗内有三种或三种以上不同类别的像素点，则用步骤10)中的仿射分类器对窗内像素点进行重新分类，更新此窗口内像素点的类标；

11c)按照步长step，将窗口在SAR图像上滑动，重复步骤11a)和11b)，直至窗口滑完整幅图像。

7.根据权利要求6所述的基于小波散射网络的SAR图像分割方法，其中步骤11b)中所述的对窗口内像素点进行重新分类所用的仿射分类器，只由部分仿射空间构成；具体来说，设窗口内有M种不同的类别，记为Y_win＝{c₁,c₂,...,c_M}，其中M≥3，则此时的仿射分类器仅由跟此M种类别对应的仿射空间构成，即给定窗内某像素点，用其散射纹理特征s来表示该像素点，该像素点的具体分类结果通过下式得到：

$\hat{c}=\underset{c\∈Y_{win}}{\arg min}\left|\right|s-P_{A_c}\left(s\right)\left|\right|.$