基于尺度优选的遥感影像分割方法与流程
本发明属于遥感影像处理领域，具体涉及一种基于对象复杂度的遥感影像分割方法，该方法可以自适应地选取最优尺度进行对象分割。
背景技术：
高(空间)分辨率遥感影像记录了丰富的精细的地物空间结构、纹理特征和场景格局，使得面向对象的遥感分析(object-orientedanalysis,ooa)成为当前高分辨率遥感应用的主要技术手段。相对于传统的基于像元(per-pixel)的遥感分析技术，面向对象的影像分析方法在利用影像光谱特征的基础上，能够从影像分析层面出发、以更类似人工解译的方式，充分挖掘高分辨率影像丰富的几何特征、纹理特征、空间格局特征等，实现更高精度和更高效率的影像分析；并且该方法能够与gis空间分析相结合，进一步融入社会经济、空间模型等高层次信息，为快速高效、高精度的遥感影像分析提供新思路。相应的面向对象遥感分析的参考文献包括，周亚男,骆剑承,程熙等.多特征融入的自适应遥感影像多尺度分割[j].武汉大学学报：信息科学版,2013,38(1):19-22、周成虎,骆剑承等.高分辨率卫星遥感影像地学计算[m].北京:科学出版社,2009、blaschke,t.,2010.objectbasedimageanalysisforremotesensing.isprsjournalofphotogrammetryandremotesensing,65(1):2-16、myint,s.w.,gober,p.,brazel,a.,grossman-clarke,s.,weng,q.,2011.per-pixelvs.object-basedclassificationofurbanlandcoverextractionusinghighspatialresolutionimagery.remotesens.environ.115(5):1145–1161等。
遥感影像分割是实现从像元到对象转变的首要步骤，是面向对象遥感分析的基础，并直接影响到后续的地物分类和目标识别的精度和效率。与一般图像相比，遥感影像具有多尺度、多目标、宽覆盖、地物类型多样等特点：(1)一般图像往往是具有“目标-背景”模式的单目标场景，目标的表达尺度被固定、目标/背景相对明确，而遥感影像中几乎不存在背景地物，每一个像元/对象都是需要分类识别的目标；(2)一般图像的覆盖范围较小、目标类型单一，而遥感影像幅宽较大、场景中地物类型多样，不同类型不同尺度的地物、同类型不同尺度的地物、甚至同一地物的多尺度形态都在充斥于整景影像。因此常规的图像分割方法难以处理遥感影像的多目标、多尺度问题；不适宜的尺度参数所导致的过分割和欠分割错误也将由此传播到后续的对象分类、目标识别和变化监测等过程中。另一方面，地理现象的发生演变都有其固有的尺度特征。在地理学研究中，需要自适应地以地理事物本身内在的时间或空间(多)尺度去认识和分析它，而不是将人为规定的时空尺度框架强加于地理事物。相应的多尺度分割参考文献包括，baatzm,
在地理学中，尺度是描述地理事物或现象的范围和大小的抽象概念；在遥感影像分割模型或方法中并没有具体的一个参数直接与地物尺度相对应，而是利用待分割影像的空间分辨率、成图比例尺、最小分割对象的区域大小、对象合并的迭代次数、对象合并的同质性或异质性阈值等参数来间接表达；分割对象之间的合并与分裂操作也恰当地模拟了小尺度与大尺度之间的转换过程。例如在自底向上的均值漂移(meanshift)多尺度合并方法中，小于最小区域阈值的分割区域将被合并到相邻的较大分割对象中，较小分割对象之间的合并也就生成了较大的分割对象，即大尺度分割结果。鉴于影像分割尺度的抽象特点以及尺度参数对影像分割与面向对象分析的重要影响，遥感学者提出了一些尺度选择方法；按照所选最优尺度的数目，可以分为单尺度优选和多尺度优选两种方法，如表1所示。
表1.遥感影像分割的尺度优选方法
单尺度优选方法能够为每一待分割影像选择一个最优尺度；但在遥感影像中一方面单一尺度难以有效处理地物场景中的多类型、多目标、多尺度问题；另一方面有时又需要从多个尺度来综合表达地物的全面特性；例如从小尺度来看，单体建筑物或道路往往与其周边环境关系紧密，而从更大尺度来看还与其位于城市或农村、所在的气候区位等格局信息相关。多尺度优选方法能够建立待分割影像的尺度序列(包含从小到大的多个分割尺度)，并试图更准确地描述影像场景中多种类型、多个目标、多个尺度的特性；相对于单尺度优选方法，多尺度是更为理想的影像分割尺度优选方法，也更符合实际应用的需求。然而无论是单尺度优选还是多尺度优选方法，仍存在几个方面的问题：(1)大多数方法是基于试错(trialanderror)或者后评估(post-evaluation)的方法，不但计算复杂，并且严重依靠经验判断、难以实现自动化与流程化；(2)一些方法强烈依赖先期的土地利用/覆盖分类的专题知识，难以运用在先验知识缺失的影像区域中；(3)方法均是为整个或者局部的影像场景选择最优尺度，并没有确认哪个分割对象是最优尺度对象；(4)影像分割生成的多尺度对象分布在多个尺度图层中，未能有效地融合到一个图层中，难以用于后续的面向对象分析。相应的分割尺度优选参考文献包括，mingd,cit,caih,etal.semivariogram-basedspatialbandwidthselectionforremotesensingimagesegmentationwithmean-shiftalgorithm[j].geoscienceandremotesensingletters,ieee,2012,9(5):813-817、yangj,lip,hey.amulti-bandapproachtounsupervisedscaleparameterselectionformulti-scaleimagesegmentation[j].isprsjournalofphotogrammetryandremotesensing,2014,94:13-24、dragut,l.,tiede,d.,levick,s.r.,2010.esp:atooltoestimatescaleparameterformultiresolutionimagesegmentationofremotelysenseddata.int.j.geograph.inform.sci.24(6):859–871、mingd,lij,wangj,etal.scaleparameterselectionbyspatialstatisticsforgeobia:usingmean-shiftbasedmulti-scalesegmentationasanexample[j].isprsjournalofphotogrammetryandremotesensing,2015,106:28-41、李秦,高锡章,张涛,等.最优分割尺度下的多层次遥感地物分类实验分析[j].地球信息科学学报,2011,13(3):409-417.等。
技术实现要素：
发明目的：针对现有技术中存在的问题，本发明公开了一种基于尺度优选的遥感影像分割方法，该方法可以自适应地选取最优尺度进行对象分割。
技术方案：本发明采用如下技术方案：
一种基于尺度优选的遥感影像分割方法，包括如下步骤：
步骤1、选定一种基于区域生长与合并的遥感影像的分割算法，并设置除分割尺度外的其他分割参数；
步骤2、依照原始待分割遥感影像的大小，预先构建以q为比率的分割尺度序列；并从大到小排序；
步骤3、选择值最大的分割尺度作为当前分割尺度；
步骤4、用当前分割尺度对遥感影像中的某个对象进行分割，生成多个较小尺度的分割子对象，计算每个分割子对象的复杂度fc；
步骤5、根据复杂度判断每个分割子对象是否需要继续分割，如果需要继续分割，选取分割尺度序列中当前分割尺度的下一个值作为当前分割尺度，跳转到步骤(4)继续分割；否则，将所述分割子对象标定为分割最优对象；
步骤6、对遥感影像中所有对象按步骤(3)-(5)进行分割，直到所有分割子对象均被标定为分割最优对象，或者当前分割尺度是分割尺度序列中的最小值时，完成分割；
步骤7、将分割出的所有对象融合到一个分割对象图层中。
优选地，步骤5中判断分割子对象是否需要继续分割的判断条件为：将步骤(4)得到的每个分割子对象复杂度fc与预设的对象复杂度阈值tfc比较，如果fc>tfc，需要继续分割，否则，不需要继续分割。
具体地，步骤4中分割子对象的复杂度fc为光谱复杂度、纹理复杂度、形状复杂度或角点复杂度中的一种；
所述光谱复杂度fcspc的计算公式为：fcspc＝d+u+h；
其中
所述纹理复杂度fct的计算公式为：fct＝j+g+s+q+cv；
其中纹理能量j：
所述形状复杂度fcs的计算公式为：fcs＝gs+gc；
其中对象光滑度gs为：gs＝po/pc；对象紧凑度gc为：
所述角点复杂度fcc的计算公式为：fcc＝re＝pe/n；
其中re为对象区域内的边缘点比率，pe为对象区域内边缘点的个数，n为对象内像元的个数。
作为一种优选，步骤4中分割子对象的复杂度fc为光谱复杂度、纹理复杂度、形状复杂度和角点复杂度中影响因素的加权和，计算公式为：
fc＝wd·d+wu·u+wh·h+
wj·j+wg·g+ws·s+wq·q+wcv·cv+
wgs·gs+wgc·gc+
wre·re
其中d为光谱标准差；u为对象光谱一致度；h为对象信息熵；j为纹理能量；g为纹理反差；s为纹理熵；q为纹理逆差；cv为纹理相关度；gs为对象光滑度；gc为对象紧凑度；re为对象区域内的边缘点比率；wd、wu、wh、wj、wg、ws、wq、wcv、wgs、wgc和wre均为加权系数。
作为另一种优选，步骤5中判断分割子对象是否需要继续分割的判断条件为：cseg＝(fc>tfc)∨(np>tnp∧rp<trp)
或：cseg＝(fc>tfc)∨(nt>tnt∧rt<trt)；
或：cseg＝(fc>tfc)∨(np>tnp∧rp<trp)∨(nt>tnt∧rt<trt)；
即按上述三式之一计算cseg的值，当cseg为真时，需要继续分割，否则，不需要继续分割；
其中tfc为预设的对象复杂度阈值；np为土地专题图上的地块个数；tnp为地块个数阈值；rp为分割对象范围内最大的地块面积与分割对象的面积之比；trp为最大地块面积占比阈值；nt为类比个数；tnt为类比个数阈值；rt为分割对象范围内最多类别的面积与分割对象的面积之比，trt为最大类别面积占比阈值；∨为或运算；∧为与运算。
作为另一种优选，分割尺度序列中元素比率q取值为2。
作为另一种优选，分割尺度序列中元素最小值为待分割影像像元数目的四百分之一。
作为另一种优选，分割算法为分水岭分割算法、均值漂移分割算法、多分辨率分割算法中的一种。
有益效果：与现有技术相比，本发明具有以下优点：(1)在自顶向下的迭代分割流程中，能够有效结合分割对象的属性特征和前期的土地专题图来分析分割对象的复杂度，并据此自适应地选取最优的分割尺度；(2)本发明公开的分割方法不仅适用于均值漂移算法，也适合于其他区域生长与合并算法，如分水岭分割算法、多分辨率分割算法等；(3)本发明公开的遥感影像分割方法将分割最优对象融合到一个分割对象图层中，便于后续的分析。
附图说明
图1是遥感影像分割尺度的优选方法示意图；
图2是遥感影像分割中尺度自适应优选的实现流程图；
图3是实施例中遥感影像分割与分割对象特征提取的流程图；
图4是分割对象是否继续下一尺度分割的判断示意图；
图5是最优尺度分割对象集合的层次组织示意图；
图6是某地zy-3融合影像的最优尺度分割效果示意图；
图7是某地zy-3融合影像的最优分割对象的尺度分布示意图。
具体实施方式
下面结合附图和具体实施方式，进一步阐明本发明。
如图1所示，是遥感影像分割尺度的优选方法示意图，其中在影像分割与特征提取步骤中，由预设分割尺度序列中从大到小地选取一个分割尺度，并利用选定的分割算法实施影像/区域分割与分割对象的特征提取；在分割对象复杂度的分析步骤中，根据所提取的对象特征来分析分割对象的复杂度，并可以通过结合前期的土地专题图判断该分割对象是否需要在下一尺度继续分割；如果需要下一尺度的分割，则由分割序列中选取下一个较小的尺度继续影像分割；以此实现影像区域或分割对象的最优尺度选择。图2是遥感影像分割尺度自适应优选的具体实现流程图，包含有影像分割与特征提取、分割尺度序列构建、分割对象复杂度分析、最优分割尺度判定、最优分割对象的层次组织等5个实现单元。
本发明公开了一种基于尺度优选的遥感影像分割方法，包括到如下步骤：
步骤1、选定一种基于区域生长与合并的遥感影像的分割算法，并设置除分割尺度外的其他分割参数；
本发明公开的遥感影像分割方法适用于多种基于区域生长与合并的遥感影像的分割算法，包括分水岭分割算法、均值漂移分割算法、多分辨率分割算法等。本实施例以均值漂移分割算法为例来说明。
如图3所示，基于均值漂移的遥感影像分割主要包含影像滤波和区域合并两个步骤；其中均值漂移算法运用于滤波过程，目的在于寻找影像中的局部极值点并生成均值的影像区域；合并过程则是寻找均值块的连通区域并形成最终的分割对象。在滤波过程中，联合考虑影像的空间位置和光谱特征，组成一个d＝p+2维的联合向量x＝(xs,xr)，其中xs表示影像中网格/像素的位置坐标，xr表示影像中网格/像素的p维向量特征。均值漂移算法依据核函数的概率密度估计，在上述向量空间中寻找局部极值点。其中多维核函数定义为：
其中k(x)为空间和光谱域中的核函数，hs,hr分别为空间域和光谱域的核带宽，c为标准化常数。则均值漂移迭代函数为：
其中xt为经过t次迭代后模点(mode)的位置，wi为x领域中像素点i的作用权重。迭代中x所经过的位置，即序列{x,m(x),m(m(x)),...}为x的轨迹；均值漂移总是指向具有最大局部密度的地方，在密度函数极大值处漂移量趋于零，则迭代结束。均值滤波通过计算特征空间中局部(邻域内)样本点的加权和来实现漂移；一方面对影像中地物内部的像元进行了平滑处理，另一方面又保留了地物的边界特征。
在区域合并中首先利用传递闭包(transitiveclosure)算法迭代地将在空间上相邻的两个(多个)具有相似模式(mode)和较小边界强度(edgestrength)的均质区域合并，以生成更大的分割区域。然后将较小的分割区域，如像素数目少于sc，合并到相邻的分割区域中，生成最终的分割对象。
在基于均值漂移的影像分割算法中有三个重要的参数(hs,hr,sc)控制着分割对象的尺寸(尺度)。在影像滤波中，空间带宽hs指示滤波移动窗口的大小，光谱带宽hr表示滤波所能允许的光谱差异；在区域合并中像素数目阈值sc决定着分割对象的最小尺寸。而在实际应用中，通常将空间带宽和光谱带宽参数固定，并通过调整sc的值来生成不同尺度的分割对象。本实施例中取hs＝7.0,hr＝6.5。
步骤2、依照原始待分割遥感影像的大小，预先构建以q为比率的分割尺度序列；并从大到小排序；
设尺度序列s＝(s0,s1,…,sn)。本发明参照影像金字塔模型中相邻层之间的抽取比率，选用比率2来生成影像分割的尺度序列，即si＝2×si+1。s0的值为经验值，需依具体的影像场景调整，本实施例中取值为待分割影像像元数目的四百分之一，然后设置s1为s0的一半、s2为s1的一半，直至sn小于20。例如一景1000×1000的遥感影像对应的分割尺度序列为(2560,1280,640,320,160,80,40,20)。
步骤3、选择值最大的分割尺度作为当前分割尺度；即初始分割尺度为s0；
步骤4、用当前分割尺度对遥感影像中的某个对象进行分割，生成多个较小尺度的分割子对象，计算每个分割子对象的复杂度fc；
采用均值漂移分割对遥感影像中的某个对象进行分割后，提取分割子对象的多种属性特征。具体包括各波段的光谱均质、光谱标准差、光谱一致度、光谱信息熵，基于灰度共生矩阵的各波段的纹理能量、纹理反差、纹理熵、纹理逆差、纹理相关度，形状光滑度、形状紧致度等多维特征。
影像分割区域即分割对象具有颜色、形状、纹理、几何、边缘、空间结构等多种特征，相应的对象复杂度也表现在光谱复杂度、纹理复杂度、角点复杂度、形状复杂度等多个方面。
(1)光谱复杂度；光谱复杂度反映对象区域内光谱值的变化程度，可以用对象区域内像元的光谱标准差d、对象光谱一致度u和对象信息熵h等进行表征：
h＝∑nj/n·log(nj/n)
其中gi为像元i的光谱值，
(2)纹理复杂度；灰度共生矩阵是一种有效描述影像纹理的方法，表达了影像灰度空间的分布形态和影像整体复杂度。基于灰度共生矩阵，提取纹理的能量j、反差g、熵s、逆差q、相关度cv等描述因子：
其中p(i,j)为影像灰度共生矩阵(i,j)处元素，μx,μy,σx,σy分别为px和py的均值和标准差，px和py分别是共生矩阵每列与每行元素和；多波段遥感影像的纹理复杂度可以用来表示所有波段复杂度的总和。
(3)形状复杂度；形状复杂度是对分割对象的空间分布和几何形状信息的描述，可以选用对象光滑度gs和紧凑度gc之和来定量表达：
gs＝po/pc
其中ao、po分别表示对象的面积和周长，pc为与对象面积相同圆的周长。
(4)角点复杂度；影像边缘点的疏密反映了影像的复杂度，利用对象区域内的边缘点比率re表达对象的角点复杂度：
re＝pe/n
其中pe为对象区域内边缘点的个数，n为对象内像元的个数。本实施例利用sobel算子提取待分割影像的角点强度图，并在分割对象区域内统计角点个数。
步骤5、根据复杂度判断每个分割子对象是否需要继续分割，如果需要继续分割，选取分割尺度序列中当前分割尺度的下一个值作为当前分割尺度，跳转到步骤(4)继续分割；否则，将所述分割子对象标定为分割最优对象；
分割对象的复杂度是影像分割中尺度内在的本质的反映。在自顶向下的影像多尺度分割流程中，如果一个分割对象具有较高的复杂度，那么该对象区域需要继续下一(更小)尺度的分割。例如在影像的水库或者湖泊区域，经过四次迭代分割生成的分割对象就足够简单，实现最优的尺度分割；而在影像的城市区域，可能需要十次迭代分割才能实现最优的尺度分割。
由此判断是否继续分割的条件可以是将步骤4得到的每个分割子对象复杂度fc与预设的对象复杂度阈值tfc比较，如果fc>tfc，需要继续分割，否则，不需要继续分割。
前期的土地专题图保存着地物的土地利用/覆盖属性和几何边界信息，能够为遥感影像的分割提供先验的指导知识。例如在城市区域往往预示着较小的地物尺度和分割对象，而在森林、草原等区域往往对应较大的地物尺度和分割对象。本发明将进一步分析分割对象局部区域内的土地专题属性和几何边界，并优化分割对象的复杂度。
(1)地块个数；在土地专题图上，对应到分割对象空间范围内所含有的土地利用单元(地块)的个数np。地块越多，则该区域的尺度越小，需要下一尺度的继续分割；反之，不需要进一步的分割。
(2)类别个数；在土地专题图上，对应到分割对象范围内所含有的土地利用单元(地块)的类别数nt。类别越多，则该区域的尺度越小，需要下一尺度的继续分割；反之，不需要进一步的分割。
(3)最大地块面积占比；在土地专题图上，对应到分割对象范围内最大的地块面积与分割对象的面积之比rp。占比越大，则该对象的复杂度越低，下一尺度分割的可能就越小；反之，继续分割的可能性越大。
(4)最大类别面积占比；在土地专题图上，对应到分割对象范围内最多类别的面积与分割对象的面积之比rt。占比越大，则该对象的复杂度越低，下一尺度分割的可能就越小；反之，继续分割的可能性越大。
如图4所示，本实施例从分割对象特征和前期土地专题图两个层面来分析分割对象的复杂度；由于这两种复杂度在信息来源、提取方法、归一化方法、单位量纲等的差异，本实施例设计了一个多层停止判定准则来选取最优的分割尺度。
首先计算分割对象的光谱复杂度、纹理复杂度、形状复杂度和角点复杂度加权和，用此加权和来表征分割对象的特征复杂度fc，计算公式如下：
fc＝wd·d+wu·u+wh·h+
wj·j+wg·g+ws·s+wq·q+wcv·cv+
wgs·gs+wgc·gc+
wre·re
其中wd、wu、wh、wj、wg、ws、wq、wcv、wgs、wgc和wre均为加权系数。
然后选择多个阈值来判断分割对象是否需要继续下一尺度的分割。当分割对象的特征复杂度fc大于tfc，或者地块个数np大于tnp并且最大地块面积占比rp小于trp，或者类比个数nt大于tnt并且最大类别面积占比rt小于trt时，该对象需要继续下一尺度分割；反之，不需要下一尺度的分割。即判断是否进行下一尺度分割的条件为：
cseg＝(fc>tfc)∨(np>tnp∧rp<trp)∨(nt>tnt∧rt<trt)
当cseg为真时，需要继续分割，否则，不需要继续分割。
即依据分割对象的复杂度，实现了最优分割尺度的选择。每个阈值可以通过多次实验中评估得到的经验值，本实施例中设置tfc为180、tnp为4、trp为0.5、tnt为3、trt为0.5。
步骤6、对遥感影像中所有对象按步骤(3)-(5)进行分割，直到所有分割子对象均被标定为分割最优对象，或者当前分割尺度是分割尺度序列中的最小值时，完成分割；
步骤7、将分割出的所有对象融合到一个分割对象图层中。
如图5所示，自顶向下的分割过程中，每次迭代分割都将生成多个分割对象，其中有的分割对象较为复杂，即对象复杂度大于阈值，被称之为分割中间对象，并在下一次迭代中继续分割；有的分割对象足够简单，被称之为分割最优对象，并停止向下的迭代分割。当整个自顶向下的分割完成时，就形成了以整个影像为根、以分割中间对象为中间节点、以分割最优对象为叶子的分割对象树结构。在分割对象的组织中，将全部分割最优对象向下映射到一个图层，并记录每个对象的迭代分割次数与最优分割尺度，生成最终的多尺度分割结果。
图6示意了利用本发明中的方法实施的某地zy-3融合影像的最优尺度分割效果，图7示意了该zy-3影像中最优分割尺度的分布。从整体分割效果来看，本发明方法在影像中的水体、林地等区域的分割对象的尺度较大，而在影像中的耕地区域的分割对象的尺度较小，而在影像中的建设用地、水浇地、道路等区域的分割对象的尺度最小。这恰好同我们对影像空间中地物最优尺度的判断和认识相一致。
本发明的实例是在pc平台上实现的。经实验证明，本发明能够找到影像局部区域最优的分割尺度，并生成更优的多尺度分割结果；并且在后续的面向对象分类中，较常规方法有着较大程度的精度改进，本实施例中影像的土地利用分类精度提高约7.4％。本发明所提及的方法可广泛应用于高分辨率遥感影像的面向对象的影像分析、分类和识别等过程中，如第三次全国农业普查、国土资源调查等大型应用。
- 还没有人留言评论。精彩留言会获得点赞！