CN104680181A - 基于似然比特征的sar图像超像素分割方法 - Google Patents

Abstract

本发明公开了一种基于似然比特征的SAR图像超像素分割方法，主要解决现有的SAR图像超像素分割方法分割结果不够准确的问题。其实现步骤为：(1)输入原始SAR图像，初始化聚类中心并对聚类中心标号；(2)对每个像素点聚类并标号；(3)对整幅SAR图像生成备选超像素；(4)去掉备选超像素中像素点个数小于阈值的无效超像素；(5)对整幅SAR图像生成超像素；(6)对每个超像素重新计算聚类中心；(7)重复步骤(2)～(6)，直到迭代次数达到设定值，得到标号图像。本发明能有效保持不同区域的边界，分割结果更准确，可用于SAR目标检测、识别及分类。

Description

基于似然比特征的SAR图像超像素分割方法

技术领域

本发明属于图像处理技术领域，特别涉及一种SAR图像超像素分割方法，可用于SAR目标检测、识别及分类。

背景技术

合成孔径雷达SAR是一种利用微波进行感知的主动传感器，它不受天气、光照等条件的限制，可对感兴趣的目标进行全天候、全天时的侦查，成为目前对地观测和军事侦察的重要手段。因而，SAR图像的处理和识别成为雷达领域研究的热点。

目前，对SAR图像的处理大多以像素为单位，用二维矩阵来表示一幅SAR图像，并未考虑像素之间的空间组织关系，这使得算法处理效率过低。超像素，是指具有相似纹理、亮度等特征的相邻像素构成的图像块。超像素利用像素之间特征的相似程度对像素进行分组，可以很大程度上降低后续SAR图像处理任务的复杂度。

SAR图像超像素分割是指按照一定的相似性准则将SAR图像划分成不同区域，每个区域中的像素具有相似的特征，这些区域被称为超像素。现有的SAR图像超像素分割算法PILS算法是一种基于SAR图像像素强度和位置相似性进行图像分割的方法。在分割SAR图像不同区域时，该算法能够比传统Turbo pixel，Normalized-Cuts，Meanshift，Quick shift图像分割算法更好地保持区域的边缘，对噪声也更加稳健。但是对于复杂场景的SAR图像，该算法分割性能有所下降。

发明内容

本发明的目的在于提出一种基于似然比特征的SAR图像超像素分割方法，以提高对于复杂场景的SAR图像分割质量。

本发明的技术方案是这样实现的：

一.技术思路：

将分割分为三个阶段：初始化聚类中心阶段，划分超像素阶段，迭代划分超像素阶段。其中：

在初始化聚类中心阶段，完成对SAR图像划分、聚类，并对所有聚类中心标号。

在划分超像素阶段，完成对SAR图像的所有像素点聚类和标号，并对SAR图像生成超像素，重新计算每一个超像素聚类中心。

在迭代划分超像素阶段，完成对新的聚类中心执行“划分超像素阶段”的操作，直到达到设定的迭代次数。

二.实现步骤

本发明基于似然比特征的SAR图像超像素分割方法，包括如下步骤：

A.初始化聚类中心步骤：

A1)输入原始SAR图像，将其划分为n个S×S的矩形块T₁,T₂,...T_n，并取矩形块T₁,T₂,...T_n的几何中心c₁,c₂,...c_n作为初始化聚类中心，其中，S为设定的矩形块边长，n＝MN/S²，M,N分别为SAR图像的行数和列数；

A2)将初始化聚类中心c₁,c₂,...c_n分别标号为1,2,...,n；

B.划分超像素步骤：

B1)对原始SAR图像的第i个像素点，将处在以像素i为中心，2S×2S为边长的矩形块内的聚类中心c_i1,c_i2,...c_ik作为像素点i的备选聚类中心；分别计算像素点i与备选聚类中心c_i1,c_i2,...c_ik的差异值D(i,c_i1),D(i,c_i2),....D(i,c_ik)，找到与像素点i差异值最小的备选聚类中心c_ij，将该备选聚类中心c_ij的标号设置为像素点i的标号；

B2)对原始SAR图像所有像素点进行B1)操作，得到SAR图像对应的标号图像L；

B3)将SAR图像中相同标号且相互连通的像素点，集合为一个备选超像素l，得到备选超像素集合l₁,l₂,...l_p，并对备选超像素集合l₁,l₂,...l_p中的超像素分别标号为1,2,...p；

B4)将备选超像素集合l₁,l₂,...l_p中像素个数小于t_min的超像素称为无效超像素，其余为有效超像素，其中t_min＝S²/h，h为设定的用于控制最小超像素大小的参数，将每个无效超像素的标号设为其邻近某个有效超像素的标号；

B5)将SAR图像中每个像素点的标号设置为其所在的备选超像素的标号，得到新的标号图像L'；

B6)将SAR图像中相同标号且相互连通的像素点集合作为一个超像素l'，得到超像素集合l₁',l₂',...l_j',...l_q'，并将超像素集合l₁',l₂',...l_j',...l_q'中每一个超像素分别标号为1,2,...j,...q；

B7)对每一个超像素l_j'得到其聚类中心c_j的坐标为(mean(x_j),mean(y_j))，将聚类中心c_j的标号设为其所在超像素l_j'的标号，其中，mean(·)表示对向量·求均值，x_j,y_j分别表示超像素l_j'包含的像素的横、纵坐标构成的列向量；

C.迭代划分超像素步骤：

C1)重复步骤B，直到迭代次数达到设定迭代次数It，取值为5。

本发明的特点和进一步改进在于：

在步骤B1)中，计算像素点i与备选聚类中心c_im的差异值D(i,c_im)，按如下步骤进行：

B11)设(x_i,y_i),分别为像素点i与备选聚类中心c_im的坐标，计算像素点i与备选聚类中心c_im的距离d(i,c_im)：

$d(i,c_{\mathrm{im}})=\sqrt{{(x_i-x_{c_{\mathrm{im}}})}^2+{(y_i-y_{c_{\mathrm{im}}})}^2};$

B12)令v₁,v₂分别为以像素点i与备选聚类中心c_im为中心的5×5的矩形块，计算v₁与v₂的强度差异值δ(v₁,v₂)：

当SAR图像强度服从视数为1的gamma分布，即指数分布时，则δ(v₁,v₂)的计算公式为：

$\mathrm{\δ}(v_1,v_2)=2M*\log \frac{\frac{1}{2}\left|\underset{k=1}{\overset{M}{\mathrm{\Σ}}}\right[I_{v_1}\left(k\right)+I_{v_2}\left(k\right)\left]\right|}{\sqrt{\left|\underset{k=1}{\overset{M}{\mathrm{\Σ}}}{I}_{v_1}\left(k\right)\right|\left|\underset{k=1}{\overset{M}{\mathrm{\Σ}}}{I}_{v_2}\left(k\right)\right|}},$

其中，M为v₁矩形块中像素的个数，为矩形块v_i中的像素的强度构成的列向量；

当SAR图像幅度服从lognormal分布时，则δ(v₁,v₂)的计算公式为

$\mathrm{\δ}(v_1,v_2)=-M*\log \left(\sqrt{2\mathrm{\π}}\frac{{\mathrm{\σ}}_1{\mathrm{\σ}}_2}{{{\mathrm{\σ}}_{12}}^2}\right),$

其中，M为v₁矩形块中像素的个数，(μ₁,σ₁ ²)为矩形块v₁像素的幅度分布参数的最大似然估计，(μ₂,σ₂ ²)为矩形块v₂像素的幅度分布参数的最大似然估计，v₁₂表示矩形块v₁,v₂中的所有的像素点，(μ₁₂,σ₁₂ ²)为v₁₂像素的幅度分布参数的最大似然估计；

B13)根据B11)与B12)的计算结果，则像素点i与备选聚类中心c_im的差异值D(i,c_im)计算公式为

D(i,c_im)＝δ(v₁,v₂)+λ*d(i,c_im)

其中λ为设定的用于调节距离d(i,c_im)在总差异值D(i,c_im)中所占比重的参数。

本发明与现有技术相比，具有如下优点：

1.本发明在SAR图像强度服从gamma分布时，以似然比特征衡量像素点与聚类中心的差异，使本发明采用gamma分布特征与现有的PILS算法相比在分割超像素时能更好的保持图像的边界，以及同一个超像素中像素点的同质性，分割效果优于现有技术；

2.本发明将似然比特征扩展至lognormal分布数据，在SAR图像强度服从lognormal分布时，采用本发明采用lognormal分布特征分割比采用现有gamma形式分割超像素时能更好的保持图像的边界。

附图说明

图1为本发明的实现流程图；

图2为本发明和现有技术使用的仿真图像；

图3为用本发明采用gamma分布特征和现有PILS技术采用像素相似性特征，分别对强度服从gamma分布的图像进行仿真分割结果图；

图4为用本发明采用gamma分布特征和现有PILS技术采用像素相似性特征，分别对强度服从gamma分布的图像进行仿真分割的性能图；

图5为用本发明采用lognormal分布特征和gamma分布特征，分别对幅度服从lognormal分布的图像进行仿真分割结果图；

图6为用本发明采用lognormal分布特征和gamma分布特征，分别对幅度服从lognormal分布的图像进行仿真分割性能图；

图7为用本发明采用gamma分布特征和现有PILS技术采用像素相似性特征，分别对真实SAR图像进行分割结果图。

具体实施方式

下面结合附图对本发明的实施步骤和效果作进一步说明。

参照图1，本发明基于似然比特征的SAR图像超像素分割方法，其实现步骤如下：

步骤1，初始化聚类中心：

1a)输入原始SAR图像，将其划分为n个S×S的矩形块T₁,T₂,...T_n，并取矩形块T₁,T₂,...T_n的几何中心c₁,c₂,...c_n作为初始化聚类中心，其中，S为设定的矩形块边长，n＝MN/S²，M,N分别为SAR图像的行数和列数。

1b)将初始化聚类中心c₁,c₂,...c_n分别标号为1,2,...,n。

步骤2，划分超像素：

2a)对原始SAR图像的第i个像素点，将处在以像素i为中心，2S×2S为边长的矩形块内的聚类中心c_i1,c_i2,...c_ik作为像素点i的备选聚类中心；

2b)分别计算像素点i与备选聚类中心c_i1,c_i2,...c_im,...c_ik的差异值D(i,c_i1),D(i,c_i2),....D(i,c_im),...D(i,c_ik)，其中m为1,2,...k：

2b1)设(x_i,y_i),分别为像素点i与备选聚类中心c_im的坐标，计算像素点i与备选聚类中心c_im的距离d(i,c_im)：

$d(i,c_{\mathrm{im}})=\sqrt{{(x_i-x_{c_{\mathrm{im}}})}^2+{(y_i-y_{c_{\mathrm{im}}})}^2};$

2b2)令v₁,v₂分别为以像素点i与备选聚类中心c_im为中心的5×5的矩形块，v₁₂表示矩形块v₁,v₂中的所有的像素点，计算v₁与v₂的强度差异值δ(v₁,v₂)：

$\mathrm{\δ}(v_1,v_2)=-\log \left(\frac{{\sup }_{{\mathrm{\θ}}_{12}}\left[p\left(v_1\right|{\mathrm{\θ}}_{12})p\left(v_2\right|{\mathrm{\θ}}_{12})\right]}{{\sup }_{{\mathrm{\θ}}_1}\left[p\left(v_1\right|{\mathrm{\θ}}_1)\right]{\sup }_{{\mathrm{\θ}}_2}\left[p\left(v_2\right|{\mathrm{\θ}}_2)\right]}\right),$

其中，表示v_i的无噪声图像参数为θ_i时，v_i中所有元素强度的联合概率密度，其中θ_i通过最大似然估计求得，v_i分别取v₁,v₂,v₁₂；

2b3)推导当SAR图像分别服从gamma分布和lognormal分布时，v₁与v₂的强度差异值δ(v₁,v₂)表达式：

当SAR图像强度服从视数为1的gamma分布，即指数分布时，图像中的像素强度x的概率分布为：

$f\left(x\right)=\frac{1}{\mathrm{\&beta;}}\exp (-\frac{x}{\mathrm{\&beta;}}),x\&GreaterEqual;0,---<1>$

其中，β为指数分布参数，β的最大似然估计为β＝Ε(x)；

将<1>式代入2b2)计算出的v₁与v₂的强度差异值δ(v₁,v₂)

$\begin{array}{c}\mathrm{\&delta;}(v_1,v_2)=-\log \frac{\underset{k=1}{\overset{2M}{\mathrm{\&Pi;}}}{f}_{12}\left[I_{v_{12}}\left(k\right)\right]}{\underset{k=1}{\overset{M}{\mathrm{\&Pi;}}}{f}_1\left[I_{v_1}\left(k\right)\right]*\underset{k=1}{\overset{m}{\mathrm{\&Pi;}}}{f}_2\left[I_{v_2}\left(k\right)\right]}\\ =-\log \frac{\underset{k=1}{\overset{2M}{\mathrm{\&Pi;}}}\frac{1}{E\left(I_{v_{12}}\right)}\exp (-\frac{I_{v_{12}}\left(k\right)}{E\left(I_{v_{12}}\right)})}{\underset{k=1}{\overset{M}{\mathrm{\&Pi;}}}\frac{1}{E\left(I_{v_1}\right)}\exp (-\frac{I_{v_1}\left(k\right)}{E\left(I_{v_1}\right)})*\underset{k=1}{\overset{M}{\mathrm{\&Pi;}}}\frac{1}{E\left(I_{v_2}\right)}\exp (-\frac{I_{v_2}\left(k\right)}{E\left(I_{v_2}\right)})}\\ =2M*\log \frac{\left|E\left(I_{v_{12}}\right)\right|}{\sqrt{\left|E\left(I_{v_1}\right)\right|\left|E\left(I_{v_2}\right)\right|}}\end{array}$

由于 $E\left(I_{v_1}\right)=\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{I}_{v_1}\left(k\right),E\left(I_{v_2}\right)=\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{I}_{v_2}\left(k\right),E\left(I_{v_{12}}\right)=\frac{1}{2M}\underset{k=1}{\overset{2M}{\mathrm{\&Sigma;}}}{I}_{v_{12}}\left(k\right),$ 于是v₁与v₂的强度差异值δ(v₁,v₂)可以写为

$\mathrm{\&delta;}(v_1,v_2)=2M*\log \frac{\frac{1}{2}\left|\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}\right[I_{v_1}\left(k\right)+I_{v_2}\left(k\right)\left]\right|}{\sqrt{\left|\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{I}_{v_1}\left(k\right)\right|\left|\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{I}_{v_2}\left(k\right)\right|}},$

其中，M为v₁矩形块中像素的个数，为矩形块v₁中的像素强度构成的列向量，为矩形块v₂中的像素强度构成的列向量，v₁₂表示矩形块v₁,v₂中的所有的像素点，为v₁₂中的像素强度构成的列向量；

当SAR图像幅度服从lognormal分布时，即对图像中的像素幅度x，其概率分布为

$f\left(x\right)=\frac{1}{\sqrt{2\mathrm{\&pi;}}\mathrm{\&sigma;}}\exp (-\frac{{(\ln x-\mathrm{\&mu;})}^2}{2{\mathrm{\&sigma;}}^2}),x>0,---<2>$

其中，(μ,σ²)为lognormal分布参数，(μ,σ²)最大似然估计为μ＝E(lnx)， $\mathrm{\&sigma;}=\sqrt{D\left(\ln x\right)}.$

将<2>式代入2b2)计算出的v₁与v₂的强度差异值δ(v₁,v₂)：

$\begin{array}{c}\mathrm{\&delta;}(v_1,v_2)=-\log \frac{\underset{k=1}{\overset{2M}{\mathrm{\&Pi;}}}{f}_{12}\left[I_{v_{12}}\left(k\right)\right]}{\underset{k=1}{\overset{M}{\mathrm{\&Pi;}}}{f}_1\left[I_{v_1}\left(k\right)\right]*\underset{k=1}{\overset{M}{\mathrm{\&Pi;}}}{f}_2\left[I_{v_2}\left(k\right)\right]}\\ =-M*\log \left(\sqrt{2\mathrm{\&pi;}}\frac{{\mathrm{\&sigma;}}_1{\mathrm{\&sigma;}}_2}{{{\mathrm{\&sigma;}}_{12}}^2}\right)+\underset{k=1}{\overset{2M}{\mathrm{\&Sigma;}}}\frac{{\left[\ln \right[I_{v_{12}}\left(k\right)]-{\mathrm{\&mu;}}_{12}]}^2}{2{{\mathrm{\&sigma;}}_{12}}^2}-\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\left[\ln \right[I_{v_1}\left(k\right)]-{\mathrm{\&mu;}}_1]}^2}{2{{\mathrm{\&sigma;}}_1}^2}-\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\left[\ln \right[I_{v_2}\left(k\right)]-{\mathrm{\&mu;}}_2]}^2}{2{{\mathrm{\&sigma;}}_2}^2}\end{array},---<3>$

将 ${{\mathrm{\&sigma;}}_i}^2=\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{\left(\ln \right[I_{v_i}\left(k\right)]-{\mathrm{\&mu;}}_i)}^2$ 和 $\begin{array}{c}\underset{k=1}{\overset{2M}{\mathrm{\&Sigma;}}}\frac{{\left[\ln \right[I_{v_{12}}\left(k\right)]-{\mathrm{\&mu;}}_{12}]}^2}{2{{\mathrm{\&sigma;}}_{12}}^2}-\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\left[\ln \right[I_{v_1}\left(k\right)]-{\mathrm{\&mu;}}_1]}^2}{2{{\mathrm{\&sigma;}}_1}^2}-\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\left[\ln \right[I_{v_2}\left(k\right)]-{\mathrm{\&mu;}}_2]}^2}{2{{\mathrm{\&sigma;}}_2}^2}\\ =\frac{2M*{{\mathrm{\&sigma;}}_{12}}^2}{2{{\mathrm{\&sigma;}}_{12}}^2}-\frac{M*{{\mathrm{\&sigma;}}_1}^2}{2{{\mathrm{\&sigma;}}_1}^2}-\frac{M*{{\mathrm{\&sigma;}}_2}^2}{2{{\mathrm{\&sigma;}}_2}^2}=0\end{array}$

代入<3>式，于是v₁与v₂的强度差异值δ(v₁,v₂)可以写为：

$\mathrm{\&delta;}(v_1,v_2)=-M*\log \left(\sqrt{2\mathrm{\&pi;}}\frac{{\mathrm{\&sigma;}}_1{\mathrm{\&sigma;}}_2}{{{\mathrm{\&sigma;}}_{12}}^2}\right),$

其中，M为v₁矩形块中像素的个数，(μ₁,σ₁ ²)为矩形块v₁像素的幅度分布参数的最大似然估计，(μ₂,σ₂ ²)为矩形块v₂像素的幅度分布参数的最大似然估计，v₁₂表示矩形块v₁,v₂中的所有的像素点，(μ₁₂,σ₁₂ ²)为v₁₂像素幅度分布参数的最大似然估计；

2b4)根据上面计算出的像素点i与备选聚类中心c_im的距离d(i,c_im)及v₁与v₂的强度差异值δ(v₁,v₂)，计算像素点i与备选聚类中心c_im的差异值D(i,c_im)

D(i,c_im)＝δ(v₁,v₂)+λ*d(i,c_im)，

其中λ为设定的用于调节距离d(i,c_im)在总差异值D(i,c_im)中所占比重的参数；

2c)找到与像素点i差异值最小的备选聚类中心c_ij，将该备选聚类中心c_ij的标号设置为像素点i的标号；

2d)对每个像素点完成标号后，得到SAR图像对应的标号图像L；

2e)将SAR图像中相同标号且相互连通的像素点集合作为一个备选超像素l，得到备选超像素l₁,l₂,...l_p，并对每一个备选超像素l₁,l₂,...l_p分别标号为1,2,...p；

2f)将备选超像素l₁,l₂,...l_p中像素个数小于t_min的超像素称为无效超像素，其余为有效超像素，其中t_min＝S²/h，h为设定的用于控制最小有效超像素大小的参数，对每个无效超像素，将该无效超像素的标号设为其临近某个有效超像素的标号；

2g)将SAR图像中每个像素点的标号设置为其所在的备选超像素的标号，得到新的标号图像L'；

2h)将SAR图像中相同标号且相互连通的像素点集合作为一个超像素l'，得到超像素集合l₁',l₂',...l_j',...l_q'，并将超像素集合l₁',l₂',...l_j',...l_q'中每一个超像素分别标号为1,2,...j,...q；

2i)对每一个超像素l_j得到其聚类中心c_j的坐标为(mean(x_j),mean(y_j))，将聚类中心c_j的标号设为其所在超像素l_j的标号，其中，mean(·)表示对向量·求均值，x_j,y_j分别表示超像素l_j包含的像素的横、纵坐标拉成的列向量。

步骤3，重复步骤2，直到迭代次数达到设定迭代次数It＝5为止，得到标号图像L'，即为本发明的超像素分割结果。

本发明的效果通过以下仿真实验说明：

1.仿真数据：

为了验证本发明对区域分割的性能，按照图2(a)设置一个含有多种类型区域的图像，其中不同强度表示不同类型的区域。图2(b)为对图2(a)的每一个同类型区域用同参数gamma分布数据的填充图；图2(c)为对图2(a)的每一个同类型区域用同参数lognormal分布数据的填充图。

实验所用的实测数据为美国Sandia实验室公开的miniSAR数据集，其分辨率为0.1m×0.1m，图像大小为2510×1638。

2.仿真实验内容：

仿真1，用本发明采用的gamma分布特征和现有PILS技术采用的像素相似性特征对图2(b)进行分割，分割结果为图3。其中，图3(a)为用现有PILS技术采用像素相似性特征对图2(b)数据的分割结果，图3(b)为用本发明采用gamma分布特征对图2(b)数据的分割结果。

从图3(a)和图3(b)可以看出，对于服从gamma分布的数据，本发明分割结果比现有PILS分割结果更符合真实区域，且错误的区域边界较少。

仿真2，用本发明采用的gamma分布特征和现有PILS技术采用的像素相似性特征在不同超像素大小下对图2(b)进行分割，分割性能曲线为图4。其中，图4(a)为本发明分割性能指标under segmentation error曲线与现有算法PILS分割性能指标undersegmentation error曲线对比，图4(b)为本发明分割性能指标boundary recall曲线与现有算法PILS分割性能指标boundary recall曲线对比。

其中，Under segmentation error是一种衡量分割结果与真实区域差异性的数值，当已知真实区域g₁,g₂,...g_M和超像素分割结果s₁,s₂,...s_L，其计算方式为

$U=\frac{1}{N}[\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}\left(\underset{\left[s_j\right|s_j\&cap;g_i>B]}{\mathrm{\&Sigma;}}\right|s_j\left|\right)-N]$

其中，|s_j|表示超像素s_j中像素的个数，N是图像中像素的总数，s_j∩g_i表示超像素s_j中泄漏到真实区域g_i的像素点个数，B设定为|s_j|的5％；

Boundary recall是一种衡量分割结果的边界与真实边界相似性的数值，即在真实边界中位于距离分割结果边界两个像素点以内的部分，该部分像素点个数占真实边界总像素个数的比例即Boundary recall。

从图4(a)和图4(b)可以看出，对于服从gamma分布的数据，在不同的超像素大小下，本发明分割结果的under segmentation error都低于现有算法PILS，本发明分割结果的boundary recall都高于现有PILS算法，且比现有PILS算法性能更加稳定。

仿真3，用本发明采用lognormal分布特征和采用gamma分布特征对图2(c)进行分割，分割结果为图5。其中，图5(a)为用本发明采用gamma分布特征对图2(c)数据的分割结果，5(b)为用本发明采用lognormal分布特征对图2(c)数据的分割结果。

从图5(a)、图5(b)可以看出，对于服从lognormal分布的数据，采用lognormal分布特征分割结果比用gamma分布特征分割结果更符合真实区域，且错误的区域边缘较少。

仿真4，用本发明采用lognormal分布特征和采用gamma分布特征在不同超像素大小下对图2(c)进行分割，分割性能曲线为图6。其中，图6(a)为用本发明技术分别采用gamma分布特征的分割性能指标under segmentation error曲线与采用lognormal分布特征的分割性能指标under segmentation error曲线对比。图6(b)为用本发明技术分别采用gamma分布特征的分割性能指标boundary recall曲线与采用lognormal分布特征的分割性能指标boundary recall曲线对比。

从图6(a)和图6(b)可以看出，对于服从lognormal分布的数据，在不同的超像素大小下，采用lognormal分布特征的分割结果的under segmentation error更低，且分割结果的boundary recall更高。

仿真5，用本发明采用gamma分布特征和现有PILS技术采用像素相似性特征对真实SAR图像进行分割，分割结果为图7。其中，图7(a)为用现有PILS技术采用像素相似性特征对真实SAR图像的分割结果。图7(b)为用本发明采用gamma分布特征对真实SAR图像的分割结果。

从图7(a)和图7(b)可以看出，对于真实SAR图像，本发明分割结果比现有PILS技术分割结果更符合真实区域，且错误的区域边缘大大减少。

Claims (2)

1.一种基于似然比特征的SAR图像超像素分割方法，包括：

A.初始化聚类中心步骤：

A1)输入原始SAR图像，将其划分为n个S×S的矩形块T₁,T₂,...T_n，并取矩形块T₁,T₂,...T_n的几何中心c₁,c₂,...c_n作为初始化聚类中心，其中，S为设定的矩形块边长，n＝MN/S²，M,N分别为SAR图像的行数和列数；

A2)将初始化聚类中心c₁,c₂,...c_n分别标号为1,2,...,n；

B.划分超像素步骤：

B1)对原始SAR图像的第i个像素点，将处在以像素i为中心，2S×2S为边长的矩形块内的聚类中心c_i1,c_i2,...c_ik作为像素点i的备选聚类中心；分别计算像素点i与备选聚类中心c_i1,c_i2,...c_ik的差异值D(i,c_i1),D(i,c_i2),....D(i,c_ik)，找到与像素点i差异值最小的备选聚类中心c_ij，将该备选聚类中心c_ij的标号设置为像素点i的标号；

B2)对原始SAR图像所有像素点进行B1)操作，得到SAR图像对应的标号图像L；

B3)将SAR图像中相同标号且相互连通的像素点，集合为一个备选超像素l，得到备选超像素集合l₁,l₂,...l_p，并对备选超像素集合l₁,l₂,...l_p中的超像素分别标号为1,2,...p；

B4)将备选超像素集合l₁,l₂,...l_p中像素个数小于t_min的超像素称为无效超像素，其余为有效超像素，其中t_min＝S²/h，h为设定的用于控制最小超像素大小的参数，将每个无效超像素的标号设为其邻近某个有效超像素的标号；

B5)将SAR图像中每个像素点的标号设置为其所在的备选超像素的标号，得到新的标号图像L'；

B6)将SAR图像中相同标号且相互连通的像素点集合作为一个超像素l'，得到超像素集合l₁',l₂',...l_j',...l_q'，并将超像素集合l₁',l₂',...l_j',...l_q'中每一个超像素分别标号为1,2,...j,...q；

B7)对每一个超像素l_j'得到其聚类中心c_j的坐标为(mean(x_j),mean(y_j))，将聚类中心c_j的标号设为其所在超像素l_j'的标号，其中，mean(·)表示对向量·求均值，x_j,y_j分别表示超像素l_j'包含的像素的横、纵坐标构成的列向量；

C.迭代划分超像素步骤：

C1)重复步骤B，直到迭代次数达到设定迭代次数It，It取值为5。

2.根据权利要求1所述方法，其中步骤B1所述的计算像素点i与备选聚类中心c_im的差异值D(i,c_im)，按如下步骤进行：

B11)设(x_i,y_i),分别为像素点i与备选聚类中心c_im的坐标，计算像素点i与备选聚类中心c_im的距离d(i,c_im)：

$d(i,c_{\mathrm{im}})=\sqrt{{(x_i-x_{c_{\mathrm{im}}})}^2+{(y_i-y_{c_{\mathrm{im}}})}^2};$

B12)令v₁,v₂分别为以像素点i与备选聚类中心c_im为中心的5×5的矩形块，计算v₁与v₂的强度差异值δ(v₁,v₂)：

当SAR图像强度服从视数为1的gamma分布，即指数分布时，则δ(v₁,v₂)的计算公式为：

$\mathrm{\&delta;}(v_1,v_2)=2M*\log \frac{\frac{1}{2}\left|\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}\right[I_{v_1}\left(k\right)+I_{v_2}\left(k\right)\left]\right|}{\sqrt{\left|\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{I}_{v_1}\left(k\right)\right|\left|\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{I}_{v_2}\left(k\right)\right|}},$

其中，M为v₁矩形块中像素的个数，为矩形块v_i中的像素的强度构成的列向量；

当SAR图像幅度服从lognormal分布时，则δ(v₁,v₂)的计算公式为

$\mathrm{\&delta;}(v_1,v_2)=-M*\log \left(\sqrt{2\mathrm{\&pi;}}\frac{{\mathrm{\&sigma;}}_1{\mathrm{\&sigma;}}_2}{{{\mathrm{\&sigma;}}_{12}}^2}\right),$

其中，M为v₁矩形块中像素的个数，(μ₁,σ₁ ²)为矩形块v₁像素的幅度分布参数的最大似然估计，(μ₂,σ₂ ²)为矩形块v₂像素的幅度分布参数的最大似然估计，v₁₂表示矩形块v₁,v₂中的所有的像素点，(μ₁₂,σ₁₂ ²)为v₁₂像素的幅度分布参数的最大似然估计；

B13)根据B11)与B12)的计算结果，则像素点i与备选聚类中心c_im的差异值D(i,c_im)计算公式为

D(i,c_im)＝δ(v₁,v₂)+λ*d(i,c_im)

其中λ为设定的用于调节距离d(i,c_im)在总差异值D(i,c_im)中所占比重的参数。