CN104361589A - 一种基于尺度间映射的高分辨率遥感影像分割方法 - Google Patents

Abstract

本发明公开一种基于尺度间映射的高分辨率遥感影像分割方法，针对对象级变化检测中涉及的地理对象提取，提出了基于小波变换与改进JSEG算法的高分辨率遥感影像多尺度分割方法；针对传统JSEG算法在高分辨率遥感影像分割中影响分割精度的关键问题，采取了相应的改进策略并取得了良好的效果。引入了小波变换作为多尺度分析工具，放弃了传统JSEG算法过于粗糙的颜色量化，以尽量保留原始影像中的细节信息。进而建立了新颖的尺度间分割结果映射机制，使当前尺度不仅基于上一尺度的分割结果进行图像分割，并且实现了对上一尺度分割结果的修正，从而有效减少了尺度间分割误差的积累。最后，提出了改进的多尺度分割策略及区域合并策略，有效减少了过分割现象及误合并现象。

Description

一种基于尺度间映射的高分辨率遥感影像分割方法

技术领域

本发明涉及一种基于尺度间映射的高分辨率遥感影像分割方法，属于遥感影像分割技术领域。

背景技术

与普通图像相比，遥感影像具有多波段、地物种类众多、纹理特征丰富、多尺度、覆盖范围宽广等特点。首先，遥感影像通常包含了多个波段的数据，使得传统单波段的影像分割方法很难直接应用于多光谱或高光谱遥感影像分割中；另外，遥感影像通常包含了丰富的纹理信息，能够综合反映各种地物复杂的空间结构特征，使得有效的抽取和表述对象的纹理特征更加困难；再者，遥感影像通常覆盖范围广泛的区域，尺寸大，云层遮盖、地物阴影等干扰因素众多，因此需要更加高效的影像分割，有时还需要引入先验知识提高分割精度；最后，遥感影像中的地物通常表现出多尺度特征，即多尺度的影像分割能够提取对象在不同尺度下的空间结构特征。这些原因都造成了成功分割遥感影像的困难，也使遥感影像分割领域始终充满了发展的动力。

近三十年来，学者们对遥感影像分割及其应用展开了广泛而深入的研究。例如，Laprade等人提出了基于分裂及合并的分割算法，其主要思想为利用F检验判断区域的均质程度，并成功应用于针对航空遥感影像的图像分割。Pan等提出将分水岭分割算法与Mean Shift聚类方法相结合，并利用小波变换进行多尺度分割，在卫星遥感影像中取得了良好的效果。Dong等人利用相邻像素间的光谱相似性，提出了基于高斯——马尔科夫模型(GMRF，Gauss-MRF)的SAR影像分割算法，并与Gamma-MRF模型的Radar分割算法进行了比较。需要指出的是，尽管目前遥感影像分割算法很多，但针对高分辨率遥感影像的图像分割研究依然在针对性及系统性方面存在诸多不足。

空间分辨率的提高带来了更加丰富的光谱信息及纹理、形状等空间信息，同时也带来了同种地物可分性增加，而不同种类地物可分性降低等问题，即更加突出的“同谱异物”、“同物异谱”现象，这些都给高分辨率遥感影像分割带来了新的挑战。Deng等提出的JSEG算法综合考虑了影像的光谱信息与纹理特征，对区域彩色纹理特征的同质性具有强大的检测能力。JSEG算法对各种自然场景具有很强的鲁棒性，并已成功应用于遥感影像分割领域。但JSEG算法也存在一些固有的局限，例如结果容易产生过分割；分割结果中不同目标间分界不明显及容易受到阴影等干扰因素的影响等。

针对这些问题，学者已提出了一些改进的JSEG方法：Zheng等利用模糊控制对影像量化过程进行了优化，所有的像素根据隶属度函数被划分为不同的类，更好的保持了颜色信息的空间分布特征，减弱了过分割现象。Chang等提出了利用影像中原始像素的彩色信息代替量化影像的IC-JSEG算法，在自然场景中获得了更好的分割结果。Komati等提出了三种改进版本的JSEG算法：Fractal-JSEG、Fractal-only和Local Fractal Dimension。前两种方法主要提高了分割结果中不同目标间边界区域的检测精度，第三种方法主要针对影像中前景目标与背景区域具有相似光谱信息与纹理特征的应用场合。上述方法均取得了不同程度的改进效果，但算法仅适用于普通图像的分割。在遥感影像分割领域，尽管改进的JSEG算法在中低分辨率遥感影像中已取得了较为理想的效果，但鲜见JSEG算法在高分辨率遥感影像分割中的应用研究。

鉴于JSEG算法固有的局限性及高分辨率遥感影像尤其在城市场景中多变的生态环境和各种人造目标给影像分割造成的困难，直接采用JSEG算法很难取得满意的效果。

发明内容

发明目的：针对现有技术中存在的问题，在影像量化、区域分割以及区域合并三个方面对JSEG算进行了改进，本发明提出了一种基于尺度间映射的高分辨率遥感影像分割算法。

技术方案：一种基于尺度间映射的高分辨率遥感影像分割方法，主要分为三个步骤：量化与小波变换、多尺度分割、区域合并。

量化与小波变换

首先对高分辨率遥感影像进行量化与小波变换。在影像量化过程将采用模糊C均值聚类(FCM，Fuzzy C-Means)方法，FCM的实现过程如下：

设其中，R^P是P维向量空间，N为样本总数。设u_ik为第k个样本属于第i类的隶属度，c为分类的总数，且满足0≤u_ik≤1，0≤k≤N，1≤i≤c，令分类后第i类的聚类中心为v_i，定义目标函数如下：

$J(U,V)=\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\underset{i=1}{\overset{c}{\mathrm{\Σ}}}{\left(u_{\mathrm{ik}}\right)}^m{\left(d_{\mathrm{ik}}\right)}^2---\left(3.2\right)$

其中，d_ik为样本k与v_i之间的欧式距离，m为加权指数，取m＝2。进而FCM算法转化为求公式(3.2)的极小值min{J(U,V)}。基于公式(3.2)，利用拉格朗日条件极值法更新隶属度矩阵和聚类中心，公式如下：

更新隶属度矩阵：

$u_{\mathrm{ik}}=\frac{1}{\underset{j=1}{\overset{c}{\mathrm{\Σ}}}\left[{(d_{\mathrm{ik}}/d_{\mathrm{jk}})}^{2/m-1}\right]}---\left(3.3\right)$

更新聚类中心：

$v_i=\frac{\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{u}_{\mathrm{ik}}^m{x}_k}{\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{u}_{\mathrm{ij}}^m}---\left(3.4\right)$

具体聚类过程如下：

Step1：设定分类总数c及停止迭代的误差阈值T_ε，且T_ε＞0。

Step2：取随机数s∈[0,1]，初始化隶属度矩阵，定位为U⁽⁰⁾，且矩阵中元素满足 $\underset{i=1}{\overset{c}{\mathrm{\Σ}}}{u}_{\mathrm{ik}}=1.$

Step3：定义迭代步数t＝1,2,...，并进行迭代计算。

Step4：利用公式(3.3)更新隶属度U^(t+1)。

Step5：利用式(3.4)更新聚类中心

Step6：利用公式(3.2)计算J(D^(t),C^(t))和J(D^(t+1),C^(t+1))，当满足|J(D^(t+1),C^(t+1))-J(D^(t),C^(t))|＜T_ε时停止迭代，将此时的V^(t+1)作为FCM的最优聚类中心，以V^(t+1)V^(t+1)和V^(t+1)作为最优FCM划分和最优聚类中心。反之，t＝t+1，重复Step4。

为了尽可能的保持图像的细节信息，将原始影像压缩为256个灰度级的8bit量化影像，即令公式(3.2)中c＝256。

在此基础上，首先采用小波变换对量化影像进行分解，并对分解后获得的低频影像序列计算多尺度J-image影像，实现由粗到精的多尺度分割。分解次数N可根据图像尺寸及具体应用场合人工设定。Haar小波尺度函数的定义如下：

${\mathrm{\ψ}}_i^j\left(x\right)=\mathrm{\ψ}(2^{j}x-i),i=0,...,(2^j-1)---\left(\text{3.6}\right)$

公式(3.6)中，i为平移参数，控制函数沿x轴平移。j为尺度参数，控制图像的缩小或放大。

多尺度分割

“多尺度”是指采用同一特定尺寸为M×M像素的窗口分别计算原始影像及全部低频图像对应的J-image。定义低频影像序列S_j(j＝1,2...N)，原始影像为S₀。多尺度分割具体步骤如下：

Step1：设定特定窗口尺寸为M×M像素，基于该窗口利用公式(2.4)计算S₀及S_j(j＝1,2...N)对应的J-image。基于J-image的定义可以看出，尽管采用了相同尺寸的窗口，但在空间分辨率更低的低频影像对应的J-image中，某一像素的J-value反映了原始影像中更大范围内的色彩分布。多尺度分割在空间分辨率最低的S_N中开始。

Step2：在S_N对应的J-image中确定种子区域。首先设定阈值a＝a₁，利用公式(3.1)计算阈值T_N，所有J-value小于T_N的像素作为候选种子点，通过4-connectivity方法连接种子点形成候选种子区域。若候选种子区域大于M×M像素，则构成一个种子区域。对剩余的非种子区域像素，设定阈值a＝a₂，且满足a₂<a₁，利用公式(3.1)计算阈值T′_N，采用相同的方法形成新的候选种子区域。对种子区域的二次提取是为了进一步提取对象的细节特征，减小阈值是为了提高新增候选种子区域的可靠性。种子区域确定后，对剩余非种子区域的像素按照J-value由小到大，逐个并入与其相邻的种子区域，直到所有像素计算完毕，从而获得当前尺度下的分割结果。

Step3：将S_N中获得的分割结果根据坐标直接映射到S_N-1对应的J-image中，并且对映射结果在S_N-1中进行修正。由于空间分辨率的提高，S_N中提取的边界映射到S_N-1中对应的是一个块状边界区域K，定义为“疑似边界区域”。为避免尺度间映射造成的误差，首先将K分割得到的区域L作为一个种子区域。对K中像素作为非种子区域像素，按照像素的J-value从小到大逐个并入与其相邻的种子区域，获得最终的映射结果。

Step4：首先对在S_N-1中由映射结果提取的每一个对象判别其在当前尺度下是否需要分割，以保证已经提取的对象不会被进一步过分割。自适应判别准则为：设R_A为待分割的某一对象，利用公式J＝(S_T-S_W)/S_W计算对象R_A对应的J-value为J_A，S_T为窗口中所有像素的总体方差，S_W为窗口中属于同一灰度级像素的方差的和。取a＝a₁，利用公式(3.1)计算阈值T_A。若满足J_A<T_A，则认为区域R_A在当前尺度下内部均质度较高，可视为一个完整目标，无需分割，直接映射到下一尺度。反之，则需要进行分割。对所有需要分割的对象重复Step2的操作，直到当前尺度分割完毕。

Step5：重复Step3到Step4的操作，直至S₀计算完毕。需要指出的是，为了避免欠分割现象，在S₀中对映射获得的所有的对象都进行分割，因此对某一对象R_A不再计算J_A，其他均与Step 4相同。

区域合并

对于分割后产生的碎片区域，首先采用彩色直方图准则进行判别，对符合判别准则的区域采用颜色标准差准则进行二次判别：

首先利用彩色直方图描述量化影像中的各个区域。对任意两个相邻区域R_A和R_B，计算两个直方图间的欧式距离D_H。公式如下：

D_H(A,B)＝||p_A-p_B||   (3.7)

p_A和p_B分别是区域R_A、R_B的彩色直方图向量。设定的阈值T_H，当D_H≤T_H时，进行下一步判别。否则，不进行合并。定义了一种颜色标准差准则如下：定义D_color：

$D_{\mathrm{color}}=\underset{b=1}{\overset{q}{\mathrm{\&Sigma;}}}\left[\frac{{\mathrm{\&sigma;}}_{b-\mathrm{merge}}}{\frac{n_D}{n_D+n_E}{\mathrm{\&sigma;}}_{\mathrm{Db}}+\frac{n_E}{n_D+n_E}{\mathrm{\&sigma;}}_{\mathrm{Eb}}}\right]---\left(\text{3.8}\right)$

R_D、R_E为符合彩色直方图准则的两个相邻待合并区域，q为遥感影像的波段总数，b为波段号，n_D和n_E分别为区域R_D、R_E所包含的像素数，σ_Db、σ_Eb、σ_b-merge分别为区域R_D、R_E以及合并后区域的标准差。单一波段下的标准差定义如下：

$\mathrm{\&sigma;}=\underset{i=0}{\overset{255}{\mathrm{\&Sigma;}}}{[i-\mathrm{\&mu;}]}^2{p}_i---\left(3.9\right)$

其中，μ为区域像素均值，p_i为像素值为i的像素所占的比重。通过公式(3.8)可以看出，与彩色直方图准则相比，D_color进一步计算了原始影像中所有波段的颜色信息。因此，不同对象即使具有相似的彩色直方图分布，但其颜色标准差距离会很大，从而有效避免误合并现象。设定区间T_C，若D_color∈T_C，则合并区域R_D、R_E。区域合并采用区域邻接图(Region Adjacency Graphics，RAG)方法。合并策略如下：

Step1：从图像的分割结果中生成Region Adjacency Graphics。

Step2：对于每个对象R_A，从RAG中确定所有R_B，若R_A与R_B邻接。

Step3：对于所有R_A与R_B，分别利用上文的区域合并准则进行判别。若满足条件，则进行合并。

Step4：更新区域邻接图，重复Step1至Step3，直到没有合并的区域为止，从而得到最终分割结果。

本发明采用上述技术方案，具有以下有益效果：针对对象级变化检测中涉及的地理对象提取，提出了基于小波变换与改进JSEG算法的高分辨率遥感影像多尺度分割方法。针对传统JSEG算法在高分辨率遥感影像分割中影响分割精度的关键问题，采取了相应的改进策略并取得了良好的效果。引入了小波变换作为多尺度分析工具，放弃了传统JSEG算法过于粗糙的颜色量化，以尽量保留原始影像中的细节信息。进而建立了新颖的尺度间分割结果映射机制，使当前尺度不仅基于上一尺度的分割结果进行图像分割，并且实现了对上一尺度分割结果的修正，从而有效减少了尺度间分割误差的积累。最后，提出了改进的多尺度分割策略及区域合并策略，有效减少了过分割现象及误合并现象。通过对多光谱及全色高分辨率遥感影像的分割实验可以得出以下结论：在分割结果中的区域总数相当时，本发明方法较JSEG算法及eCognition提取目标的边缘更加准确，保持目标轮廓更加完整，同时能有效识别场景中的地物遮挡、阴影区域，显著提高了分割精度，同时采用本发明方法较eCognition更有利于提高对象级变化检测算法的透明性与鲁棒性。

附图说明

图1为SPOT 5融合影像；

图2为12灰度级量化影像；

图3为本发明实施例的方法流程图；

图4为直接映射结果；

图5为最终映射结果；

图6为eCognition分割流程；

图7为本发明方法区域合并前；

图8为本发明方法区域合并后；

图9为JSEG分割结果；

图10为eCognition分割结果；

图11为IKONOS 1m全色影像；

图12为本发明方法分割结果；

图13为JSEG分割结果；

图14为eCognition分割结果。

具体实施方式

下面结合具体实施例，进一步阐明本发明，应理解这些实施例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。

本实施例首先介绍了JSEG算法分割的基本流程，进而详细分析了JSEG算法在高分辨率遥感影像分割中影响分割精度的关键问题，并针对这些问题阐述了本发明方法的改进策略与具体实现。实验采用高分辨率IKONOS全色影像以及SPOT 5多光谱影像，并分别对JSEG算法、国际知名商业软件eCognition以及本发明方法的分割结果进行了比较和分析。

JSEG算法及存在的关键问题

JSEG算法实现

传统JSEG算法将影像分割划分为两个步骤：颜色量化与空间分割。在颜色量化中，首先将图像转换到LUV颜色空间，利用对等组滤波(PGF，Peer GroupFiltering)对影像进行平滑去噪，最后采用经典的硬聚类算法(HCM，HardC-Means)方法获得量化影像。

在空间分割中，首先计算多尺度J-image影像序列，J-value体现了以像素为z中心的特定尺寸窗口中灰度分布的均质程度，即J-value值越大，窗口中地物的类间差异更大，而类内差异更小。在某一尺度J-image中，首先根据阈值T_J确定种子区域，分割出影像中对象的中心区域，T_J的定义如下：

T_J＝μ_J+aσ_J   (3.1)

其中，μ_J和σ_J分别为像素z所在窗口内J-value的均值和标准差。a为预先设定的阈值，通常取a∈[-0.4,0.4]。当前尺度J-image中，小于阈值T_J的点采用4-connectivity方法获得联通区域，即从种子点开始，在上下左右四个方向上移动，扩大区域面积。若该区域尺寸区域大于当前尺度下种子区域的最小尺寸时，即可作为一个种子区域，该最小尺寸由人工预先设定。

在确定种子区域后，首先计算非种子区域像素的J-value均值，小于均值的像素构成增长区域。若该增长区域仅与一个种子区域相邻，则纳入该种子区域。对于剩下的非种子区域像素，利用更小尺寸的窗口计算J-image。重复上一步，生成新的增长区域直到最小尺度。对剩余的像素根据J-value从小到大逐个纳入相邻的种子区域。最后，为解决过分割问题，根据区域颜色相似性信息根据欧式距离最小进行区域合并，获得最终的分割结果。

JSEG在高分辨率遥感影像分割中的关键问题

通过分析JSEG的基本原理及相关文献，总结JSEG算法在高分辨率遥感影像分割中影响分割精度的关键问题主要包括：

问题1：在颜色量化阶段，JSEG算法量化后的影像通常只包含十几个灰度级^[97-101]。这样的影像量化尽管能够有效减少计算量，但同时会造成高分辨率遥感影像中大量细节信息的损失，降低了类间可分性，从而难以准确定位对象的边缘。以2004年中国上海地区空间分辨为2.5m的全色—多光谱融合SPOT 5影像为例，如图1所示。若直接采用JSEG中的方法量化影像结果如图2所示，量化影像中只包含12个灰度级数。

比较图1、图2可以看出，量化影像中道路、建筑以及其他人造目标轮廓信息损失严重，植被等地物与周边环境具有相似的光谱特征，此时，单纯为了减少运算量而采用的过于粗糙的量化已变得没有意义。

问题2：在空间分割阶段，基于多尺度J-image影像的分割通过改变特定窗口尺寸能够提取不同范围内的颜色分布特征，但本质上只利用了原始影像在单一空间分辨率下的空间信息。根据专家的知识，单一空间分辨率下获得的结果不一定可靠，且更容易受到噪声的干扰。

问题3：JSEG算法的区域增长准则中，某一尺度下的种子区域一旦确定，就不会再生成新的种子区域，非种子区域的像素在区域合并中只会作为增长区域最终纳入原先的种子区域。这样的区域增长准则很容易造成分割结果中图像细节的缺失和欠分割现象。

问题4：较大尺寸的窗口可能包含了多个对象，获得的J-image也更为粗糙，因而大尺度上分割结果很难准确的定位目标的边缘。如果窗口过大，甚至造成分割的结果失去意义。JSEG算法中大尺度上获得的粗糙分割结果在最终区域合并前得到保留，下一个更小尺度的分割直接在上一尺度粗分割基础上进行。因而，大尺度上获得的粗糙边界无法在下一尺度得到修正，这样的误差积累将导致最终分割结果难以准确定位对象的边缘。

问题5：JSEG算法在分割过程中，当前尺度的分割建立在上一尺度分割结果的基础之上。因此无论在上一尺度是否已经完整分割出了某一对象的完整轮廓，JSEG算法都会用一个更小的尺度对这一对象再次进行判别。但在高分辨率遥感影像中，对象内部的纹理特征更加复杂，因此JSEG算法的分割策略很容易造成过分割现象。即使JSEG算法的合并策略能够恢复部分对象的整体轮廓，但仍然造成了不必要的运算。

问题6：JSEG算法采用R.O.Duda等提出的方法进行区域合并，以基于区域彩色直方图的欧式距离作为合并准则。尽管彩色直方图能够较好的表示区域的光谱信息与纹理特征，但在高分辨率遥感影像中，地物种类众多，影像中很多不同种类的地物可能具有相似的直方图分布，原有合并策略很容易产生错误的合并结果。

所提出方法实现流程

针对上述所提出的关键问题，基于尺度间映射的高分辨率遥感影像分割方法，主要分为三个步骤：量化与小波变换、多尺度分割、区域合并。算法流程如图3所示。

量化与小波变换

首先对高分辨率遥感影像进行量化与小波变换。针对上述“问题1”，在影像量化过程将采用模糊C均值聚类(FCM，Fuzzy C-Means)方法代替HCM方法。

FCM算法是最初由Bezkek提出，是HCM方法的自然推广，它采用模糊划分，通过优化目标函数得到每个样本点对所有类中心的隶属度，从而决定样本点的类属。FCM的实现过程如下：

设其中，R^P是P维向量空间，N为样本总数。设u_ik为第k个样本属于第i类的隶属度，c为分类的总数，且满足0≤u_ik≤1，0≤k≤N，1≤i≤c，令分类后第i类的聚类中心为v_i，定义目标函数如下：

$J(U,V)=\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}\underset{i=1}{\overset{c}{\mathrm{\&Sigma;}}}{\left(u_{\mathrm{ik}}\right)}^m{\left(d_{\mathrm{ik}}\right)}^2---\left(3.2\right)$

其中，d_ik为样本k与v_i之间的欧式距离，m为加权指数，取m＝2。进而FCM算法转化为求公式(3.2)的极小值min{J(U,V)}。基于公式(3.2)，利用拉格朗日条件极值法更新隶属度矩阵和聚类中心，公式如下：

更新隶属度矩阵：

$u_{\mathrm{ik}}=\frac{1}{\underset{j=1}{\overset{c}{\mathrm{\&Sigma;}}}\left[{(d_{\mathrm{ik}}/d_{\mathrm{jk}})}^{2/m-1}\right]}---\left(3.3\right)$

更新聚类中心：

$v_i=\frac{\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{u}_{\mathrm{ik}}^m{x}_k}{\underset{j=1}{\overset{N}{\mathrm{\&Sigma;}}}{u}_{\mathrm{ij}}^m}---\left(3.4\right)$

具体聚类过程如下：

Step1：设定分类总数c及停止迭代的误差阈值T_ε，且T_ε＞0。

Step2：取随机数s∈[0,1]，初始化隶属度矩阵，定位为U⁽⁰⁾，且矩阵中元素满足 $\underset{i=1}{\overset{c}{\mathrm{\&Sigma;}}}{u}_{\mathrm{ik}}=1.$

Step3：定义迭代步数t＝1,2,...，并进行迭代计算。

Step4：利用公式(3.3)更新隶属度U^(t+1)。

Step5：利用式(3.4)更新聚类中心

Step6：利用公式(3.2)计算J(D^(t),C^(t))和J(D^(t+1),C^(t+1))，当满足|J(D^(t+1),C^(t+1))-J(D^(t),C^(t))|＜T_ε时停止迭代，将此时的V^(t+1)作为FCM的最优聚类中心，以V^(t+1)V^(t+1)和V^(t+1)作为最优FCM划分和最优聚类中心。反之，t＝t+1，重复Step4。

为了尽可能的保持图像的细节信息，将原始影像压缩为256个灰度级的8bit量化影像，即令公式(3.2)中c＝256。

在此基础上，针对“问题2”首先采用小波变换对量化影像进行分解，并对分解后获得的低频影像序列计算多尺度J-image影像，实现由粗到精的多尺度分割。Alfréd Haar等人提出的Haar小波是目前唯一既具有对称性又是有限支撑的正交小波。同时Haar小波仅取+1和－1，计算简单。因此采用Haar小波，分解次数N可根据图像尺寸及具体应用场合人工设定。Haar小波尺度函数的定义如下：

${\mathrm{\&psi;}}_i^j\left(x\right)=\mathrm{\&psi;}(2^{j}x-i),i=0,...,(2^j-1)---\left(\text{3.6}\right)$

公式(3.6)中，i为平移参数，控制函数沿x轴平移。j为尺度参数，控制图像的缩小或放大。

多尺度分割

在多尺度分割阶段，鉴于小波分解后的各尺度低频影像都包含了下一尺度低频及高频影像的全部信息，因此，影像分割只在各尺度低频图像序列及原始影像中进行。在JSEG算法中，“多尺度”是利用不同特定尺寸窗口对原始影像计算获得的J-image序列。而本发明方法中的“多尺度”是指采用同一特定尺寸为M×M像素的窗口分别计算原始影像及全部低频图像对应的J-image。定义低频影像序列S_j(j＝1,2...N)，原始影像为S₀。多尺度分割具体步骤如下：

Step1：设定特定窗口尺寸为M×M像素，基于该窗口利用公式(2.4)计算S₀及S_j(j＝1,2...N)对应的J-image。基于J-image的定义可以看出，尽管采用了相同尺寸的窗口，但在空间分辨率更低的低频影像对应的J-image中，某一像素的J-value反映了原始影像中更大范围内的色彩分布。多尺度分割在空间分辨率最低的S_N中开始。

Step2：在S_N对应的J-image中确定种子区域。首先设定阈值a＝a₁，利用公式(3.1)计算阈值T_N，所有J-value小于T_N的像素作为候选种子点，通过4-connectivity方法连接种子点形成候选种子区域。若候选种子区域大于M×M像素，则构成一个种子区域。针对“问题3”，对剩余的非种子区域像素，设定阈值a＝a₂，且满足a₂<a₁，利用公式(3.1)计算阈值T′_N，采用相同的方法形成新的候选种子区域。对种子区域的二次提取是为了进一步提取对象的细节特征，减小阈值是为了提高新增候选种子区域的可靠性。种子区域确定后，对剩余非种子区域的像素按照J-value由小到大，逐个并入与其相邻的种子区域，直到所有像素计算完毕，从而获得当前尺度下的分割结果。

Step3：针对“问题4”，将S_N中获得的分割结果根据坐标直接映射到S_N-1对应的J-image中，并且对映射结果在S_N-1中进行修正。由于空间分辨率的提高，S_N中提取的边界映射到S_N-1中对应的是一个块状边界区域K，定义为“疑似边界区域”，如图4所示。为避免尺度间映射造成的误差，首先将K分割得到的区域L作为一个种子区域。对K中像素作为非种子区域像素，按照像素的J-value从小到大逐个并入与其相邻的种子区域，获得最终的映射结果，如图5所示。其中，图4及图5为实验部分第一组实验过程中产生的中间图像。

Step4：针对“问题5”，首先对在S_N-1中由映射结果提取的每一个对象判别其在当前尺度下是否需要分割，以保证已经提取的对象不会被进一步过分割。自适应判别准则为：设R_A为待分割的某一对象，如图5所示，利用公式J＝(S_T-S_W)/S_W计算对象R_A对应的J-value为J_A，S_T为窗口中所有像素的总体方差，S_W为窗口中属于同一灰度级像素的方差的和。取a＝a₁，利用公式(3.1)计算阈值T_A。若满足J_A<T_A，则认为区域R_A在当前尺度下内部均质度较高，可视为一个完整目标，无需分割，直接映射到下一尺度。反之，则需要进行分割。对所有需要分割的对象重复Step2的操作，直到当前尺度分割完毕。

Step5：重复Step3到Step4的操作，直至S₀计算完毕。需要指出的是，为了避免欠分割现象，在S₀中对映射获得的所有的对象都进行分割，因此对某一对象R_A不再计算J_A，其他均与Step 4相同。

区域合并

针对“问题6”，对于分割后产生的碎片区域，首先采用彩色直方图准则进行判别，对符合判别准则的区域采用颜色标准差准则进行二次判别：

首先利用彩色直方图描述量化影像中的各个区域。对任意两个相邻区域R_A和R_B，计算两个直方图间的欧式距离D_H。公式如下：

D_H(A,B)＝||p_A-p_B||   (3.7)

p_A和p_B分别是区域R_A、R_B的彩色直方图向量。设定的阈值T_H，当D_H≤T_H时，进行下一步判别。否则，不进行合并。定义了一种颜色标准差准则如下：定义D_color：

$D_{\mathrm{color}}=\underset{b=1}{\overset{q}{\mathrm{\&Sigma;}}}\left[\frac{{\mathrm{\&sigma;}}_{b-\mathrm{merge}}}{\frac{n_D}{n_D+n_E}{\mathrm{\&sigma;}}_{\mathrm{Db}}+\frac{n_E}{n_D+n_E}{\mathrm{\&sigma;}}_{\mathrm{Eb}}}\right]---\left(\text{3.8}\right)$

R_D、R_E为符合彩色直方图准则的两个相邻待合并区域，q为遥感影像的波段总数，b为波段号，n_D和n_E分别为区域R_D、R_E所包含的像素数，σ_Db、σ_Eb、σ_b-merge分别为区域R_D、R_E以及合并后区域的标准差。单一波段下的标准差定义如下：

$\mathrm{\&sigma;}=\underset{i=0}{\overset{255}{\mathrm{\&Sigma;}}}{[i-\mathrm{\&mu;}]}^2{p}_i---\left(3.9\right)$

其中，μ为区域像素均值，p_i为像素值为i的像素所占的比重。通过公式(3.8)可以看出，与彩色直方图准则相比，D_color进一步计算了原始影像中所有波段的颜色信息。因此，不同对象即使具有相似的彩色直方图分布，但其颜色标准差距离会很大，从而有效避免误合并现象。设定区间T_C，若D_color∈T_C，则合并区域R_D、R_E。区域合并采用区域邻接图(Region Adjacency Graphics，RAG)方法。合并策略如下：

Step1：从图像的分割结果中生成Region Adjacency Graphics。

Step2：对于每个对象R_A，从RAG中确定所有R_B，若R_A与R_B邻接。

Step3：对于所有R_A与R_B，分别利用上文的区域合并准则进行判别。若满足条件，则进行合并。

Step4：更新区域邻接图，重复Step1至Step3，直到没有合并的区域为止，从而得到最终分割结果。

实验结果与分析

实验采用两组不同传感器类型、不同分辨率的高分辨率遥感影像，分别检验本发明方法对多光谱及全色影像的分割效果。实验结果除与传统JSEG算法比较外，还将与德国Definiens Imaging公司开发的商业软件eCognition进行比较。下文首先阐述了eCognition分割基本原理及关键参数，进而对两组实验结果进行了分析和比较。

eCognition影像分割

eCognition是目前公认最好的面向高分辨率遥感影像分割及分类的商业软件，其采用基于分形网络演化方法(Fractal Net Evolution Approach，FNEA)的多分辨率分割，通过模拟人类感知系统，同时考虑了光谱特征与分割单元的形状特征。实现流程如图6所示。

如图6所示，eCognition用多尺度方法对影像进行多次分割，构建网络层次结构：当前尺度分割是基于上一个尺度中获得的对象实现的，并对分割结果中满足合并条件的相邻区域进行合并；另一方面，当前尺度的影像分割也受到下一个尺度分割结果所提取对象边界的约束，构成拓扑关系的网络结构。分割过程主要包括三个参数需要人工设定：

(1)尺度参数，决定了影像分割能够提取的最小对象的尺寸大小，是影像分割结果的最重要参数。参数值越大，对象的总数越少，而对象内部的平均面积越大。尺度参数直接决定了对象提取的精度，合理的尺度参数应能够使分割后的多边形准确定位对象的边缘，同时不能破坏对象完整性，避免造成大量的碎片区域。

(2)形状参数，用于保持分割结果中对象轮廓的完整性。

(3)紧致度参数，用于提高具有相似光谱特征的对象的可分性。

在实验中，通过改变以上三个参数对实验结果进行了比较，从而确定分割效果最佳的参数组合。

数据集1实验结果与目视分析

实验中，第一组数据采用2003年上海地区SPOT 5全色—多光谱融合影像，全色影像空间分辨率为2.5m，多光谱波段包括红色、绿色及近红外波段，空间分辨率为10m，影像尺寸为420×420像素。实验中，本发明方法采用两次Haar小波分解，组成包括原始量化影像及两幅低频影像的图像序列，设定特定尺寸窗口大小为M×M＝5×5像素，其他参设定为a₁＝-0.2，a₂＝-0.4，T_H＝0.18，T_C＝[2.5,3.5]。JSEG算法采用的特定窗口尺寸为20×20像素、10×10像素、5×5像素，阈值a＝0.2。在eCongition分割中，设定尺度参数为80，形状参数为0.5，紧致度参数为0.5。实验结果如图7～图10所示。为便于分析本发明方法所提出的区域合并策略对分割结果的影响，图7，图8分别为区域合并前后的分割结果。图9、图10分别为JSEG算法及eCognition算法分割结果。

由于实验影像空间分辨率为2.5m，可首先根据目视检验分割结果的准确性，为便于分析，在图8、图9、图10中选择部分典型地物进行了标号。对比三种算法的分割结果可以看出，JSEG算法的分割精度明显低于本发明方法和eCognition，例如图中的区域具有规则的形状，JSEG没有完整的提取对象轮廓，而本发明方法和eCognition则效果较好。因此实验分析着重比较本发明方法和eCognition分割结果：两种算法都较好的分割出图像中的湖泊区域水体区域本发明方法明显对湖泊的边缘的定位更加准确，如而eCognition对水体边界的分割比本发明方法更为准确，如对于典型的人造目标如 本发明方法对位于和的目标的提取明显好于eCognition。对于对象内部色彩同质性较高的区域如eCognition存在过分割及误分割现象，而本发明方法对区域的整体性保持更好。而对于内部光谱特征存在较大差异的斑块如本发明方法提取的轮廓信息也较eCognition更加准确。最后，通过对比图7和图8，可以看出本发明方法提出的合并规则能够明显改善过分割现象如同时基本没有产生误合并现象。

数据集2实验结果与目视分析

为进一步验证本发明方法算法对单波段影像的分割性能，选择南京地区2008年采集的IKONOS全色影像进行实验，影像空间分辨率为1m，尺寸为256×256像素，如图11所示。实验中，所提出本发明方法算法的参数设置为：Haar小波分解层数为2，局部窗口尺寸设置为M×M＝5×5像素，a₁＝0.2，a₂＝-0.2，T_H＝0.18，T_C＝[2.5,3.5]。JSEG算法采用的特定窗口尺寸为20×20像素、10×10像素、5×5像素，阈值a＝0.2。采用eCognition分割时尺度参数设置为30，形状参数为0.8，紧致度参数为0.5。其他参数设定与数据集1实验相同。由于是全色影像，在本发明方法中不需对图像进行量化，而直接进行影像分割。三种算法的实验结果如图12～图14所示：

通过目视分析可以看出，场景中主要包含了桥梁，水体，人造目标以及其他纹理特征丰富的地物信息。在桥梁轮廓不明显的区域只有本发明方法准确的分割出了正确的结果。而eCognition在两个较大分割单元之间存在一些狭长的虚假单元如三种算法虽然都分割出了的流区域，但本发明方法提取的边界信息最准确，同时较为完整的提取了水体与堤岸之间的边坡区域如JSEG则没有区分开水体与边坡区域，eCognition则出现了明显的误分割。本发明方法分割结果中的桥梁以及水体区域的轮廓非常完整，有效避免了对光谱同质性较高区域的过分割现象。JSEG尽管进行了区域合并，依然存在严重的过分割现象，eCognition也存在同样的问题。图像中最为典型的人造目标是位于 的建筑，三种算法虽然都分割出了楼房建筑的轮廓特征，但效果有所不同。以位于位置的建筑为例，楼房与楼房的阴影光谱特征差异明显，三种算法都能够区分出楼房区域与阴影区域，而本发明方法与eCognition进一步区分出了建筑的楼顶区域。而位于位置的建筑，建筑的阴影区域较小，且楼顶与阴影以及相邻地面区域具有相似的光谱特征。只有本发明方法准确分割出了楼顶、阴影与地面区域，而eCognition对建筑与相邻地面区域的分割不够准确，JSEG算法则没有分割出阴影与楼顶区域。

实验结果定量分析

以上检验均通过目视分析获得，为进一步验证算法的效果，采用Deng等提出的分割算法精度评价方法(文献Deng Y,Manjunath B S.Unsupervisedsegmentation of color-texture regions in images and video[J].Pattern Analysis andMachine Intelligence,IEEE Transactions on,2001,23(8):800-810.)对两组实验进行定量分析：

定义精度评价指标如下所示：

$\stackrel{\&OverBar;}{J}=\frac{1}{R}{\mathrm{\&Sigma;}}_{t=1}^T{S}_t{J}_t---\left(3.10\right)$

其中，R为图像中像素总数，T为分割结果中区域的总数，S_t为第t个区域像素总数，J_t为第t个区域的J-value，可通过公式J＝(S_T-S_W)/S_W计算得到。当分割结果中区域的个数一定时，越小，则分割结果中各区域的平均内部同质性越高，分割效果越好。

表1 分割结果与精度评价

如表1所示，尽管三种算法中需要设置的参数不尽相同，但在对同一幅图像的分割结果中区域总数相当，因此这样的实验结果是具有可比性的，参数设置的差异对比较各算法分割效果的影响不大。另外，两组实验中本发明方法算法在区域合并后指标均明显降低，证明了所提出区域合并准则的有效性。

Claims (1)

1.一种基于尺度间映射的高分辨率遥感影像分割方法，其特征在于，主要分为三个步骤：量化与小波变换、多尺度分割、区域合并；

量化与小波变换

首先对高分辨率遥感影像进行量化与小波变换；在影像量化过程将采用模糊C均值聚类(FCM，Fuzzy C-Means)方法，FCM的实现过程如下：

设其中，R^P是P维向量空间，N为样本总数；设u_ik为第k个样本属于第i类的隶属度，c为分类的总数，且满足0≤u_ik≤1， 0≤k≤N，1≤i≤c，令分类后第i类的聚类中心为v_i，定义目标函数如下：

$J(U,V)=\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}\underset{i=1}{\overset{c}{\mathrm{\&Sigma;}}}{\left(u_{\mathrm{ik}}\right)}^m{\left(d_{\mathrm{ik}}\right)}^2---\left(3.2\right)$

其中，d_ik为样本k与v_i之间的欧式距离，m为加权指数，取m＝2；进而FCM算法转化为求公式(3.2)的极小值min{J(U,V)}；基于公式(3.2)，利用拉格朗日条件极值法更新隶属度矩阵和聚类中心，公式如下：

更新隶属度矩阵：

$u_{\mathrm{ik}}=\frac{1}{\underset{j=1}{\overset{c}{\mathrm{\&Sigma;}}}{[(d_{\mathrm{ik}}/d_{\mathrm{jk}})}^{2/m-1}]}---\left(3.3\right)$

更新聚类中心：

$v_i=\frac{\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{u}_{\mathrm{ik}}^m{x}_k}{\underset{j=1}{\overset{N}{\mathrm{\&Sigma;}}}{u}_{\mathrm{ij}}^m}---\left(3.4\right)$

具体聚类过程如下：

Step1：设定分类总数c及停止迭代的误差阈值T_ε，且T_ε＞0；

Step2：取随机数s∈[0,1]，初始化隶属度矩阵，定位为U⁽⁰⁾，且矩阵中元素满足 $\underset{i=1}{\overset{c}{\mathrm{\&Sigma;}}}{u}_{\mathrm{ik}}=1;$

Step3：定义迭代步数t＝1,2,...，并进行迭代计算；

Step4：利用公式(3.3)更新隶属度U^(t+1)；

Step5：利用式(3.4)更新聚类中心

Step6：利用公式(3.2)计算J(D^(t),C^(t))和J(D^(t+1),C^(t+1))，当满足|J(D^(t+1),C^(t+1))-J(D^(t),C^(t))|＜T_ε时停止迭代，将此时的V^(t+1)作为FCM的最优聚类中心，以V^(t+1)V^(t+1)和V^(t+1)作为最优FCM划分和最优聚类中心；反之，t＝t+1，重复Step4；

将原始影像压缩为256个灰度级的8bit量化影像，即令公式(3.2)中c＝256；

在此基础上，首先采用小波变换对量化影像进行分解，并对分解后获得的低频影像序列计算多尺度J-image影像，实现由粗到精的多尺度分割；分解次数N可根据图像尺寸及具体应用场合设定；Haar小波尺度函数的定义如下：

${\mathrm{\&psi;}}_i^j\left(x\right)=\mathrm{\&psi;}(2^{j}x-i),i=0,...,(2^j-1)---\left(3.6\right)$

公式(3.6)中，i为平移参数，控制函数沿x轴平移；j为尺度参数，控制图像的缩小或放大；

多尺度分割

“多尺度”是指采用同一特定尺寸为M×M像素的窗口分别计算原始影像及全部低频图像对应的J-image；定义低频影像序列S_j(j＝1,2...N)，原始影像为S₀；多尺度分割具体步骤如下：

Step1：设定特定窗口尺寸为M×M像素，基于该窗口利用公式J＝(S_T-S_W)/S_W计算S₀及S_j(j＝1,2...N)对应的J-image；多尺度分割在空间分辨率最低的S_N中开始；

Step2：在S_N对应的J-image中确定种子区域；首先设定阈值a＝a₁，利用公式(3.1)计算阈值T_N，所有J-value小于T_N的像素作为候选种子点，通过4-connectivity方法连接种子点形成候选种子区域；若候选种子区域大于M×M像素，则构成一个种子区域；对剩余的非种子区域像素，设定阈值a＝a₂，且满足a₂<a₁，利用公式(3.1)计算阈值T′_N，采用相同的方法形成新的候选种子区域；种子区域确定后，对剩余非种子区域的像素按照J-value由小到大，逐个并入与其相邻的种子区域，直到所有像素计算完毕，从而获得当前尺度下的分割结果；

Step3：将S_N中获得的分割结果根据坐标直接映射到S_N-1对应的J-image中，并且对映射结果在S_N-1中进行修正；由于空间分辨率的提高，S_N中提取的边界映射到S_N-1中对应的是一个块状边界区域K，定义为“疑似边界区域”；首先将K分割得到的区域L作为一个种子区域；对K中像素作为非种子区域像素，按照像素的J-value从小到大逐个并入与其相邻的种子区域，获得最终的映射结果；

Step4：首先对在S_N-1中由映射结果提取的每一个对象判别其在当前尺度下是否需要分割，以保证已经提取的对象不会被进一步过分割；自适应判别准则为：设R_A为待分割的某一对象，利用公式J＝(S_T-S_W)/S_W计算对象R_A对应的J-value为J_A，S_T为窗口中所有像素的总体方差，S_W为窗口中属于同一灰度级像素的方差的和；取a＝a₁，利用公式(3.1)计算阈值T_A；若满足J_A<T_A，则认为区域R_A在当前尺度下内部均质度较高，可视为一个完整目标，无需分割，直接映射到下一尺度；反之，则需要进行分割；对所有需要分割的对象重复Step2的操作，直到当前尺度分割完毕；

Step5：重复Step3到Step4的操作，直至S₀计算完毕；为了避免欠分割现象，在S₀中对映射获得的所有的对象都进行分割，因此对某一对象R_A不再计算J_A，其他均与Step 4相同。

区域合并

对于分割后产生的碎片区域，首先采用彩色直方图准则进行判别，对符合判别准则的区域采用颜色标准差准则进行二次判别：

首先利用彩色直方图描述量化影像中的各个区域；对任意两个相邻区域R_A和R_B，计算两个直方图间的欧式距离D_H；公式如下：

D_H(A,B)＝||p_A-p_B||    (3.7)

p_A和p_B分别是区域R_A、R_B的彩色直方图向量；设定的阈值T_H，当D_H≤T_H时，进行下一步判别；否则，不进行合并；定义了一种颜色标准差准则如下：定义D_color：

$D_{\mathrm{color}}=\underset{b=1}{\overset{q}{\mathrm{\&Sigma;}}}\left[\frac{{\mathrm{\&sigma;}}_{b-\mathrm{merge}}}{\frac{n_D}{n_D+n_E}{\mathrm{\&sigma;}}_{\mathrm{Db}}+\frac{n_E}{n_D+n_E}{\mathrm{\&sigma;}}_{\mathrm{Eb}}}\right]---\left(3.8\right)$

R_D、R_E为符合彩色直方图准则的两个相邻待合并区域，q为遥感影像的波段总数，b为波段号，n_D和n_E分别为区域R_D、R_E所包含的像素数，σ_Db、σ_Eb、σ_b-merge分别为区域R_D、R_E以及合并后区域的标准差；单一波段下的标准差定义如下：

$\mathrm{\&sigma;}=\underset{i=0}{\overset{255}{\mathrm{\&Sigma;}}}{[i-\mathrm{\&mu;}]}^2{p}_i---\left(3.9\right)$

其中，μ为区域像素均值，p_i为像素值为i的像素所占的比重。通过公式(3.8)可以看出，与彩色直方图准则相比，D_color进一步计算了原始影像中所有波段的颜色信息；因此，不同对象即使具有相似的彩色直方图分布，但其颜色标准差距离会很大，从而有效避免误合并现象。设定区间T_C，若D_color∈T_C，则合并区域R_D、R_E；区域合并采用区域邻接图(Region Adjacency Graphics，RAG)方法；合并策略如下：

Step1：从图像的分割结果中生成Region Adjacency Graphics；

Step2：对于每个对象R_A，从RAG中确定所有R_B，若R_A与R_B邻接；

Step3：对于所有R_A与R_B，分别利用上文的区域合并准则进行判别；若满足条件，则进行合并；

Step4：更新区域邻接图，重复Step1至Step3，直到没有合并的区域为止，从而得到最终分割结果。