CN107657618A - 基于遥感影像和地形数据的区域尺度侵蚀沟自动提取方法 - Google Patents

Abstract

本发明公开一种基于遥感影像和地形数据的区域尺度侵蚀沟自动提取方法，属于区域土壤侵蚀领域。针对现有侵蚀沟提取方法无法应用于较大范围的技术瓶颈，本发明提出了一种融合地形骨架信息和影像特征的算法思想，在传统基于影像特征进行侵蚀沟提取的基础上，通过构建与侵蚀沟分布相适配的地形骨架信息，对提取结果进行修正，有效提高了区域尺度侵蚀沟提取精度和效率。该方法可采用中分辨率的地形数据和免费获取的Google Earth影像，扩展了方法的适用范围，为区域尺度土壤侵蚀调查与评价、水土保持决策等提供了重要方法支撑。

Description

基于遥感影像和地形数据的区域尺度侵蚀沟自动提取方法

技术领域

本发明属于区域土壤侵蚀领域，特别涉及一种基于遥感影像和地形数据的区域尺度侵蚀沟自动提取方法。

背景技术

沟蚀是一种地表径流冲刷并破坏土壤及其母质，形成切入地表以下沟壑的土壤侵蚀形式(唐克丽,2004)。相比较面蚀和细沟侵蚀，沟蚀对农业生产、生态环境的危害要更大(Valentin et al.,2005)，沟蚀也被认为是土地退化的重要形式(UNEP，1994)。

侵蚀沟是沟蚀作用所形成的沟渠，其规模介于细沟和河道之间(Posen et al.,2003)。侵蚀沟一方面被沟蚀作用所塑造，另一方面也是侵蚀产沙的重要场所，因此对侵蚀沟形态特征的研究可为沟蚀评价、沟蚀监测以及沟蚀过程与机理研究提供重要的基础数据(郑粉莉等，2016)。侵蚀沟提取可理解为侵蚀沟制图，是沟蚀研究的重要内容，主要包括提取侵蚀沟的空间分布及其二维地形参数(沟长、沟面积、沟道密度等)。根据研究尺度，侵蚀沟的提取大致可分为小区尺度、坡面尺度、小流域尺度和区域尺度四种。小区尺度和坡面尺度的侵蚀沟提取主要和侵蚀学研究相结合，探究侵蚀机理与构建侵蚀模型，这两类尺度主要采用了人工测量方法提取侵蚀沟信息，包含了径流小区法(肖培青等.,2008)、测针法(张新和,2007)、摄影测量法(Marzolff and Poesen,2009)、GPS及三维激光扫描技术等(Perroy et al.,2010)。小流域尺度侵蚀沟提取方法则更多基于数字地形分析方法(Evansand Lindsay.,2010；Castillo et al.,2014)和遥感影像分析方法(Shruthi et al.,2011；d’Oleire-Oltmanns et al.,2014)，其数据源主要包括了高精度卫星影像数据、三维激光扫描数据和无人机摄影测量数据。

区域尺度侵蚀沟提取主要是指在较大空间范围上进行侵蚀沟形态信息的获取，主要服务于区域土壤侵蚀评价以及区域水土保持治理。现有研究中主要方法包含了以下两种：(1)目视解译方法，即通过人工判读，基于影像特征，手动勾绘出侵蚀沟的范围。该方法对解译人员具有一定专业度要求，同时效率较低，但是其精度可以得到较好的保证(闫业超等，2006；Zhang et al.,2015)。(2)基于影像数据的自动提取方法，其中面向像元的分析方法较为普遍,该方法以单个栅格为分析单元，通过设定规则逐栅格判断以实现整个研究区的侵蚀沟制图(Knightet al.,2007)。随着高空间分辨率影像的普及，综合光谱信息、几何信息和结构信息的面向对象分析方法得到发展，面向对象方法通过将具有较高同质性的栅格组合成对象，以对象为分析单元进行特征提取分析，其主要优点在于充分利用了几何、结构和光谱信息，效率高，分类效果优。在侵蚀沟提取研究中，部分学者也开始使用面向对象方法(Shruthiet al.,2014)。

总体而言，面向大区域范围的侵蚀沟提取方法研究相对较少，主要有三点原因：第一，覆盖大区域的高质量数据获取成本较高，这也限制学者在大区域范围进行侵蚀沟提取方法的研究；第二，研究区域的扩展导致研究数据量的激增，对计算与存储资源提出了较高要求；第三，随着研究区的扩大，侵蚀沟的区域差异性特征会凸显，因此，如何设计侵蚀沟提取规则并确定其适用域是一大难点。当前，随着数据获取能力、计算机分析能力的不断增强，区域尺度的土壤侵蚀研究日益增加，相关的研究成果也与全球变化、生态系统紧密联系，在此背景下，对区域尺度侵蚀沟提取方法的研究，可为区域土壤侵蚀及相关研究提供重要的方法和技术支撑。

参考文献：

[1]Castillo C,Taguas E V,Zarco‐Tejada P,et al.The normalizedtopographic method:an automated procedure for gully mapping using GIS[J].Earth Surface Processes and Landforms,2014,39(15):2002-2015.

[2]d’Oleire-Oltmanns S,Marzolff I,Tiede D,et al.Detection of Gully-Affected Areas by Applying Object-Based Image Analysis(OBIA)in the Region ofTaroudannt,Morocco[J].Remote Sensing,2014,6(9):8287-8309.

[3]Evans M,Lindsay J.High resolution quantification of gully erosionin upland peatlands at the landscape scale[J].Earth Surface Processes andLandforms,2010,35(8):876-886.

[4]Knight J,Spencer J,Brooks A,et al.Large-area,high-resolutionremote sensing based mapping of alluvial gully erosion in Australia’stropical rivers[C]//Proceedings of the 5th Australian Stream ManagementConference.Charles Sturt University,2007:199-204.

[5]Marzolff I,Poesen J.The potential of 3D gully monitoring with GISusing high-resolution aerial photography and a digital photogrammetry system[J].Geomorphology,2009,111(1):48-60.

[6]Perroy R L,Bookhagen B,Asner G P,et al.Comparison of gully erosionestimates using airborne and ground-based LiDAR on Santa Cruz Island,California[J].Geomorphology,2010,118(3):288-300.

[7]Poesen J,Nachtergaele J,Verstraeten G,et al.Gully erosion andenvironmental change:importance and research needs[J].Catena,2003,50(2):91-133.

[8]UNEP,1994.United Nations Convention to CombatDesertification.United Nations Environmental Programme,Geneva.

[9]Valentin C,Poesen J,Li Y.Gully erosion:impacts,factors and control[J].Catena,2005,63(2):132-153.

[10]Zhang S W,Li F,Li T Q,Yang J C,Bu K,Chang L P,Wang W J,Yan YC.Remote sensing monitoring of gullies on a regional scale:a case study ofKebai region in Heilongjiang Province,China[J].Chinese Geographical Science,2015,25,602-611.

[11]肖培青,郑粉莉,汪晓勇,等.黄土坡面侵蚀方式演变与侵蚀产沙过程试验研究[J].水土保持学报,2008,22(1):24-27.

[12]唐克丽.中国水土保持.科学出版社,北京,2004.

[13]闫业超,张树文,李晓燕,等.黑龙江克拜黑土区50多年来侵蚀沟时空变化[J].地理学报,2006,60(6):1015-1020.

[14]闫业超,张树文,岳书平.基于Corona和Spot影像的近40年黑土典型区侵蚀沟动态变化[J].资源科学,2006,28(6):154-160.

[15]郑粉莉,徐锡蒙,覃超.沟蚀过程研究进展[J].农业机械学报,2016,47(8):48-59.

[16]张新和.黄土坡面片蚀-细沟侵蚀-切沟侵蚀演变与侵蚀产沙过程研究[J].陕西杨凌:中国科学院水土保持研究所,2007.

发明内容

为解决现有区域尺度侵蚀沟提取方法效率低、精度差等问题，本发明提出一种融合地形骨架信息和影像特征的侵蚀沟提取方法，基于可免费获取的Google Earth影像数据和Aster地形数据，实现对区域尺度侵蚀沟的自动提取。

针对上述发明目的，本发明所采用的技术方案如下：

一种基于遥感影像和地形数据的区域尺度侵蚀沟自动提取方法，包括以下步骤：

步骤1，基于多层次流域单元的数据划分：

获取研究区地形数据，基于侵蚀沟形态的空间异质性确定流域划分阈值，将研究区划分为若干个样区单元，对每个样区单元确定一个小流域，用来生成训练数据；流域是天然的地理边界，具有严格的地理含义，同时其本身也具有多层次的特点，因此本发明中将不同层次的流域作为单元划分的主要依据，将将研究区划分为若干个样区单元，每一个样区单元内部采用相同的侵蚀沟提取规则；

步骤2，根据划分的样区单元进行影像数据的下载，并对影像数据进行预处理；

步骤3，基于预处理影像，对步骤1中每个样区单元确定的小流域的侵蚀沟进行目视解译，获取训练数据；

步骤4，基于影像数据，采用面向对象的方法，经过对象分割、分割对象特征计算和侵蚀沟提取模型构建，获取侵蚀沟的初始提取结果；

步骤5，基于地形数据，生成与侵蚀沟分布相适配的汇水网络，并结合河道数据构建地形骨架，基于所构建的地形骨架信息对侵蚀沟的初始提取结果进行修正；

步骤6，对修正后的提取结果进行精度分析并将结果合并输出。

本发明的方法，所述步骤1还包括，当单个样区单元的数据量超过计算机处理能力时，采用更小的流域划分阈值，将单个样区单元划分为若干个处理单元；当选取的小流域面积小于处理单元时，则将其隶属的处理单元选中，采用适合的流域划分阈值，提取出采样单元。

所述步骤2中，影像数据选用春秋季节影像。春秋季节的影像受植被及天气影响较小，可减少天气因素对处理结果精度的影响；进一步的，当涉及到多幅影像时，对所述多幅影像进行影像融合和匀色处理。多幅影像的情况下，对影像进行影像融合和匀色处理可减小单个样区单元内部影像数据之间的差异。

所述步骤4中，对象分割采用多尺度分割算法，结合训练样区的样本确定每个样区的分割参数，每一个样区单元采用一组分割参数；对于进一步划分了处理单元的样区单元，对样区单元所包含的处理单元进行批处理；分割对象特征计算的计算参数包含光谱特征、纹理特征和形状特征；侵蚀沟提取模型采用随机森林算法，基于训练数据，构建预测模型，并将其运用到全部区域。

所述步骤5中，基于地形数据构建地形骨架的方法如下：

根据地形数据，确定与侵蚀沟的相适配的汇水网络起点和终点，其起点阈值的设定时保证大多数的汇水网络的源点处在沟头区域的分割对象内部，即相交但不越界；其终点阈值的设定根据研究区河道和侵蚀沟的界限确定。地形骨架信息包含汇水网络和河道两部分，针对侵蚀沟提取结果的修正，合理的汇水网络需要在整体上不超过侵蚀沟的沟沿线，同时尽可能多的表达出侵蚀沟的支沟发育情况。因此，对于不同样区而言，汇水网络的起点阈值是需要根据侵蚀沟的发育特征来确定。汇水网络的终点即河道的起点，其阈值确定需要考虑当地河道发育的实际情况，在方法简化情况下，可以采用50平方公里的经验值。

基于地形骨架信息的结果对侵蚀沟初始提取结果的修正方法如下：

将表征地形骨架的线对象和分割后的面对象进行空间关联分析，实现对提取结果的修正，其修正规则如下：

(a)预测为非侵蚀沟的对象，与任意一二级汇水网络相交则标记为侵蚀沟区域；

(b)预测为非侵蚀沟的对象，与任意三级及以上汇水网络中点相交则标记为侵蚀沟区域；

(c)预测为侵蚀沟的对象，与任意汇水网络均不相交则标记为非侵蚀沟区域；

(d)与河网相交的任何对象均标记为非侵蚀沟区域。

基于地形骨架信息的结果修正主要解决基于影像数据提取结果中错分和漏分问题。

流水侵蚀是侵蚀沟发育主要动力，侵蚀沟的分布必然受到地表汇流的影响，同时，侵蚀沟的发育客观上也形成或改变了汇流结构。因此本发明中将基于DEM提取的汇水网络和河道作为辅助侵蚀沟提取的地形骨架信息。地形骨架信息一方面可以框定侵蚀沟的主体结构特征，另一方面，其对数据分辨率的包容度较大。可以在保证精度的同时降低对原始数据精度的要求。

本发明的方法，地形数据和影像数据精度为中分辨率；地形数据来源为AsterGDEM数据，影像数据来源为Google Earth影像数据。Aster GDEM数据和Google Earth影像数据可免费获取，来源方便。

本发明具有以下两点优势：

(1)本发明提出了一种可用于区域尺度的侵蚀沟自动提取方法，为大范围的侵蚀沟信息普查、区域土壤侵蚀定量评价等提供了技术支撑。

(2)本发明基于融合地形骨架和影像特征的思路，方法对数据的依赖程度较低，基于可获取性强的中分辨率地形数据和Google Earth影像数据即可得到较高的提取精度。

附图说明

图1本发明实施例提供的样区图；

图2本发明算法流程图；

图3本发明采用的多层次流域划分策略；

图4本发明实施例提供的不同汇流阈值情况下的汇水网络；

图5本发明实施例提供的陕西延安样区侵蚀沟提取结果；

图6本发明实施例提供的山西河曲样区侵蚀沟提取结果；

图7本发明实施例提供的甘肃华池样区侵蚀沟提取结果；

图8本发明实施例提供的宁夏泾源样区侵蚀沟提取结果。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。以下实施案例用于说明本发明，但不用来限制本发明的范围。

本发明的实施例以黄土高原水土重点流失区为例，对本发明的方法作进一步阐述。

如图1所示，黄土高原水土重点流失区主要包含了晋西—陕北—陇东一线，总面积约为15.72万平方公里。影像数据采用Google Earth影像数据，地形数据则采用Aster GDEM数据。

如图2所示，是本发明的流程图，本实施例包括以下步骤：

步骤1、获取研究区地形数据，基于侵蚀沟形态的空间异质性确定流域划分阈值，将研究区划分为若干个样区单元；如图3所示，结合野外调查和相关专家知识，基于2000平方公里的流域面积阈值，将研究区划分为50个样区单元。基于第一个层次数据划分后，所得到的流域单元虽确保了地貌、侵蚀沟形态以及数据特征的一致性，但是对于算法分析而言，整个单元统一处理对计算资源要求过高。特别是部分流域影像数据量超过10GB，在进行多尺度分割等计算复杂度较高的分析时，计算时间过长。因此，进一步采用100平方公里流域面积阈值，将研究区划分为952个处理单元。本发明中提取规则的确定采用了监督分类的思路，因此对于每一个样区单元需要确定一个训练区域。训练样区的确定一方面要保证样本的代表性，使得基于该区域构建的侵蚀沟提取模型能够在整个样区单元内推广。同时，由于训练数据的生成需要人工解译出侵蚀沟的范围，如采用第二层次的流域单元，则人工解译的工作量过大。因此本研究中，进行了第三层次的流域划分，采用10平方公里的流域面积阈值进一步划分出采样单元。

步骤2、根据划分的样区单元进行影像数据的下载，并对影像数据进行预处理；对于单个样区单元，当涉及到多幅影像时，进行影像融合和匀色处理。

步骤3、基于预处理影像，对每一个样区单元，确定其对应的采样单元，通过野外考察和室内目视解译，得到训练数据。

步骤4、基于Google Earth影像，采用面向对象的方法，经过对象分割、分割对象特征计算和侵蚀沟提取模型构建，获取侵蚀沟的初始提取结果；本实施例中具体处理方法如下：

(1)采用多尺度分割算法，确定分割参数，对每一个样区单元内部的处理单元进行批处理，实现对象分割；

(2)对分割后的对象计算其特征值；

分割对象特征计算的特征列表如表1所示；

表1分割对象特征计算列表

(3)采用随机森林算法，实现侵蚀沟的初步提取。

步骤5、基于DEM数据，针对每个样区单元侵蚀沟发育特征，确定汇流累积阈值，生成与分割对象相适配的汇水网络。汇水网络与侵蚀沟的适配程度决定了地形骨架信息的价值，也决定了最终的修正效果。当汇流累计值阈值较小时，汇水网络可能会出现超过侵蚀沟沟头区域的现象，此时，虽然对于沟头部分漏分区域得到修正，但同时会造成新的错分区域，这种情况可认为是一种过修正方案。而如果将汇流累计值调整到相对较大时，会出现侵蚀沟沟头区域表达不足，甚至部分小的支沟没法表达的现象，此时沟头区域虽然不会出现错分现象，但会造成大量的漏分问题，这种情况可认为是一种欠修正方案。合理的情况是采用一种较为折中的方案，整体上保证大多数的汇水网络源点处在沟头区域的分割对象内部，即相交但不越界。如图4所示，基于不同汇流累积阈值生成的汇水网络在整体上均表达了较好的侵蚀沟的结构特征。当阈值为100时，放大区域主沟道能够被汇水网络有效表达，然而其左侧支沟区域并无对应的汇水网络；而当阈值为20时，放大图区域主沟道与支沟均能通过汇水网络予以表达，但同时出现了汇水网络超过沟沿线现象；相比较之下阈值为50时，结果汇水网络较为理想，即能表达出主沟道和主沟，同时也没有出现越界，基于该汇水网络，通过与分割对象进行空间连接可实现对侵蚀沟提取结果的修正。

关于修正规则，并非对汇水网络和分割对象之间直接判定是否相交即可。在实验中我们发现，对于沟头区域而言，一般为单一侵蚀沟个体，反应在汇水网络上即保持同一方向的直线，曲折现象较少。而随着汇流现象的发生，汇水网络会发生不断曲折的现象，此时很容易出现汇水网络与侵蚀沟两侧坡面上的对象相交，造成侵蚀沟范围被扩大化。基于这种情况，本发明中对判定规则进行了一定优化，对于第一级和第二级汇水网络，采用较为宽松的Intersect算子，即接触就判定为相交；而对于汇水网络的其他部分，采用较为严格的Have center in算子，即线段的中心点落入对象才判定为相交。

本实施例中采用Python语言实现了基于地形骨架的侵蚀沟提取结果修正，核心代码如表2所示。

表2地形骨架信息修正核心代码

步骤6、对每个样区样区单元完成侵蚀沟提取后，进行精度评价，并将结果合并输出，完成对整体研究区侵蚀沟的自动提取。

图5、图6、图7、图8分别是采用本发明的方法在陕西延川、山西河曲、甘肃华池和宁夏泾源样区的侵蚀沟提取结果。

Claims (10)

1.一种基于遥感影像和地形数据的区域尺度侵蚀沟自动提取方法，其特征在于，包括以下过程：

步骤1，基于多层次流域单元的数据划分：

获取研究区地形数据，基于侵蚀沟形态的空间异质性确定流域划分阈值，将研究区划分为若干个样区单元，对每个样区单元确定一个小流域，用来生成训练数据；

步骤2，根据划分的样区单元进行影像数据的下载，并对影像数据进行预处理；

步骤3，基于预处理影像，对步骤1中每个样区单元确定的小流域的侵蚀沟进行目视解译，获取训练数据；

步骤4，基于影像数据，采用面向对象的方法，经过对象分割、分割对象特征计算和侵蚀沟提取模型构建，获取侵蚀沟的初始提取结果；

步骤5，基于地形数据，生成与侵蚀沟分布相适配的汇水网络，并结合河道数据构建地形骨架，基于所构建的地形骨架信息对侵蚀沟的初始提取结果进行修正；

步骤6，对修正后的提取结果进行精度分析并将结果合并输出。

2.根据权利要求1所述的方法，其特征在于，所述步骤1还包括，当单个样区单元的数据量超过计算机处理能力时，采用更小的流域划分阈值，将单个样区单元划分为若干个处理单元；当选取的小流域面积小于处理单元时，则将其隶属的处理单元选中，采用适合的流域划分阈值，提取出采样单元。

3.根据权利要求1所述的方法，其特征在于，所述步骤2中，影像数据选用春秋季节影像。

4.根据权利要求1所述的方法，其特征在于，所述步骤2中，当涉及到多幅影像时，对所述多幅影像进行影像融合和匀色处理。

5.根据权利要求1或2所述的方法，其特征在于，所述步骤4中，对象分割采用多尺度分割算法，每一个样区单元采用一组分割参数；对于进一步划分了处理单元的样区单元，对样区单元所包含的处理单元进行批处理。

6.根据权利要求1所述的方法，其特征在于，所述步骤4中，分割对象特征计算的计算参数包含光谱特征、纹理特征和形状特征。

7.根据权利要求1所述的方法，其特征在于，所述步骤4中，侵蚀沟提取模型采用随机森林算法，基于训练数据，构建预测模型，并将其运用到全部区域。

8.根据权利要求1所述的方法，其特征在于，所述步骤5中，基于地形数据构建地形骨架的方法如下：

根据地形数据，确定与侵蚀沟的相适配的汇水网络起点和终点，其起点阈值的设定时保证大多数的汇水网络的源点处在沟头区域的分割对象内部，即相交但不越界；其终点阈值的设定根据研究区河道和侵蚀沟的界限确定。

9.根据权利要求1所述的方法，其特征在于，所述步骤5中，基于地形骨架信息的结果对侵蚀沟初始提取结果的修正方法如下：

将表征地形骨架的线对象和分割后的面对象进行空间关联分析，实现对提取结果的修正，其修正规则如下：

（a）预测为非侵蚀沟的对象，与任意一二级汇水网络相交则标记为侵蚀沟区域；

（b）预测为非侵蚀沟的对象，与任意三级及以上汇水网络中点相交则标记为侵蚀沟区域；

（c）预测为侵蚀沟的对象，与任意汇水网络均不相交则标记为非侵蚀沟区域；

（d）与河网相交的任何对象均标记为非侵蚀沟区域。

10. 根据权利要求1所述的方法，其特征在于，地形数据来源为Aster GDEM数据，影像数据来源为Google Earth影像数据；地形数据和影像数据精度为中分辨率。