CN113850822A - 一种基于汇流分割的斜坡单元自动划分方法 - Google Patents

一种基于汇流分割的斜坡单元自动划分方法 Download PDF

Info

Publication number
CN113850822A
CN113850822A CN202111100800.7A CN202111100800A CN113850822A CN 113850822 A CN113850822 A CN 113850822A CN 202111100800 A CN202111100800 A CN 202111100800A CN 113850822 A CN113850822 A CN 113850822A
Authority
CN
China
Prior art keywords
slope
unit
basin
slope unit
threshold
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202111100800.7A
Other languages
English (en)
Other versions
CN113850822B (zh
Inventor
魏进兵
杨仲康
邓建辉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Sichuan University
Original Assignee
Sichuan University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Sichuan University filed Critical Sichuan University
Priority to CN202111100800.7A priority Critical patent/CN113850822B/zh
Publication of CN113850822A publication Critical patent/CN113850822A/zh
Application granted granted Critical
Publication of CN113850822B publication Critical patent/CN113850822B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/11—Region-based segmentation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/05—Geographic models
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/60—Analysis of geometric attributes
    • G06T7/62—Analysis of geometric attributes of area, perimeter, diameter or volume
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30181—Earth observation
    • G06T2207/30184—Infrastructure

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Geometry (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Software Systems (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Remote Sensing (AREA)
  • Computer Graphics (AREA)
  • Image Analysis (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种基于汇流分割的斜坡单元自动划分方法，通过给定数字高程模型和预定义的斜坡单元最小面积阈值与圆方差阈值，实现斜坡单元的自动划分；同时遵循数字地貌最大内部均一性和外部异质性原则，基于内部同质性指数和外部异质性指数建立斜坡单元划分质量的客观定量评价方法，实现斜坡单元的质量评定和优选。本发明能够为山区滑坡灾害敏感性分析和水文分析提供高质量的评价单元。

Description

一种基于汇流分割的斜坡单元自动划分方法
技术领域
本发明属于地貌景观数字化研究领域，具体涉及一种基于汇流分割的斜坡单元自动划分方法。
技术背景
斜坡单元是一种重要的地貌单元类型，定义为山谷线和山脊线所围成的区域。通过水文汇流分析，提取山脊线、山谷线、台地和谷底边界等地貌分界线，理论上能够实现地貌景观的细分。将广域而复杂的地理空间区域，甚至一个完整的自然流域自动划分为可重复、地貌学上具有一致性的斜坡单元，能够为地貌学和水文建模、基于遥感影像的山区滑坡灾害风险评价建模提供可靠的评价单元。斜坡单元需要保证其内部地貌特征具有较好的同质性，与外部单元间具有较大异质性。目前，斜坡单元的最佳细分仍然停留在理论探索阶段，在实际操作中存在巨大的挑战。已有斜坡单元划分方法往往将多个斜坡单元甚至是一个流域作为一个斜坡单元，斜坡单元划分结果不够准确，不能满足实际应用的需要。
发明内容
本发明的目的在于针对现有技术的不足，提出一种基于汇流分割的斜坡单元自动划分方法，通过给定数字高程模型(DEM)和预定义的参数，实现斜坡单元的最佳细分，为山区滑坡灾害敏感性分析和水文分析提供高质量的评价单元。
本发明提供的基于汇流分割的斜坡单元自动划分方法，包括以下步骤：
(1)DEM数据和参数准备
准备覆盖待划分的目标区域的数字高程模型(DEM)，通过GIS软件表面分析工具获取三维视角的山体阴影图，依据人工视觉选取典型区域，手动绘制部分斜坡单元作为斜坡单元样本集，统计斜坡单元面积的分布区间，获取斜坡单元最小面积阈值区间和最大限定面积。
(2)生成最佳初始子流域盆地
将DEM数据中每个栅格单元标记为流量单元(此处流量并不是水文意义上的河流流量，实际表示具有高程属性的栅格单元自高向低处汇聚，累加的栅格单元数量)，按照DEM表示的区域地形特征和由高到低的汇流方向，累加每个流量单元周围8个方向内产生的汇聚栅格单元数量，并逐步累加到下一个流量单元，获取区域汇流累计流量。标记为高流量的栅格单元表示集中汇水区域，一般为河道。定义累计流量阈值FT，当栅格单元中累计流量大于累计流量阈值FT时，该栅格单元被标记为汇水点，将所有标记的汇水点连接起来，形成汇水线；流量为0的栅格单元为局部地形高点，标记为排水点，所有标记的排水点连接起来，形成排水线；排水线和汇水线相交形成的闭合环形为流域盆地，根据在汇水线两侧的位置，汇水盆地划分为左右子流域盆地。
将执行斜坡单元划分时，第一次进行分割的子流域盆地称为初始子流域盆地。当累计流量阈值FT取值较低时，产生密集、小面积的子流域盆地，斜坡单元划分形态细碎；当累计流量阈值FT取值较高时，产生数量较少的大规模子流域盆地，斜坡单元划分结果较粗糙。因此，需要根据试算得出的最佳初始子流域盆地累计流量阈值FTOPT，对目标区域进行划分，生成目标区域最佳初始子流域盆地。
(3)基于最佳初始子流域盆地的斜坡单元细分
对步骤(2)确定的最佳初始子流域盆地，在预先定义的斜坡单元最小面积阈值s和圆方差阈值D组合下执行迭代细分，划分为两个以上的下一级子流域盆地；当子流域盆地面积小于斜坡单元最小面积阈值s或子流域盆地圆方差小于圆方差阈值D时，则判断该子流域盆地符合预定义的斜坡单元标准，将其作为候选斜坡单元，不再执行细分；否则该子流域盆地继续作为初始子流域盆地，返回迭代细分环节，继续执行细分；当所有子流域盆地符合预定义的斜坡单元标准时，迭代中止，斜坡单元细分完成。每种斜坡单元最小面积阈值s和圆方差阈值D组合下划分得到的斜坡单元为一个斜坡单元组。
斜坡单元最小面积阈值s的上限、下限从手绘斜坡单元样本集中得到客观参考；斜坡单元圆方差D为地貌特征中每个栅格单元坡向的法向量在平均值上的变化，介于0～1之间，计算公式为：
Figure BDA0003270637020000021
式中：D为圆方差，NV为每个流域盆地栅格单元的数量，R为栅格单元坡向单位向量绝对值的总和。通过坡向均一性标准控制斜坡单元内部地貌一致性特征。坡向的圆方差越小，说明该斜坡单元内部所有栅格单元的坡向越接近，即为具有均一朝向的斜坡面。
(4)斜坡单元清理
该步骤旨在识别和清理不切实际或者难以接受的斜坡单元，并进行合并。例如，在多个排水线汇聚的平坦地形河谷出口可能会产生不切实际的细碎多边形，在地形陡变的区域或由于DEM数据本身的噪声可能会产生狭长畸形斜坡单元。斜坡单元清理方法为：设置斜坡单元清理面积阈值(Cleansize)参数，对候选斜坡单元细分中产生的细碎多边形进行删除，合并填充到相邻较大的候选斜坡单元中；同时依据经验和实际勘测，将狭长畸形斜坡单元定义为最窄处小于100m的斜坡单元，对符合这些特征的斜坡单元进行识别和合并。
(5)斜坡单元质量评价与优选
步骤(3)和(4)已经实现了斜坡单元的精细化划分和清理，但仍可能存在较多的细碎斜坡单元，导致后期计算的冗余。最佳的斜坡单元应该保证斜坡单元内部地貌特征具有最大一致性，而不同斜坡单元之间显示最大异质性。在同质性合并和异质性分割时，同一种分类共享某些相同的像素特征，而不同分类间像素特征显著不同。将DEM栅格单元视为要分割的对象，地貌坡向等同于像素值，则可以建立斜坡单元分类质量的定量评价指标，对斜坡单元进行质量评价和优选。
基于以上假设，定义自相关指数S为内部同质性指数，内部同质性指数S数值越低，表明斜坡单元内部同质性越好。定义自相关指数L为外部异质性指数，外部异质性指数L数值越高，表明斜坡单元间具有更好的异质性：
Figure BDA0003270637020000031
Figure BDA0003270637020000032
式中：S为内部同质性指数；在每个斜坡单元组中将每个斜坡单元标记为i，第n个斜坡单元的面积标记为Sn；第n个斜坡单元的圆方差标记为dn；L是外部异质性指数；αn为第n个斜坡单元的平均坡向；
Figure BDA0003270637020000033
为研究区栅格单元的平均坡向；wnl为空间临近度指标，当编号为n和编号为l的斜坡单元相邻时，该值为1，否则该值为0。
同时，需要特别强调，地貌特征坡向由角度表示，并不能直接用于数值计算，需要进行三角函数转换，转换格式如下：
Figure BDA0003270637020000041
Figure BDA0003270637020000042
斜坡单元细分时，斜坡单元最小面积阈值s和圆方差阈值D的最佳组合是较低的内部同质性指数和较高的外部异质性指数的结合，综合考虑两种细分指标，通过以下函数进行量化：
Figure BDA0003270637020000043
式中：Smax、Smin分别是不同的斜坡单元最小面积阈值s和圆方差阈值D组合得到的所有斜坡单元组内部同质性指数的最大值和最小值，Lmax、Lmin分别是所有斜坡单元组外部异质性指数的最大值和最小值。
当F(S，L)取值最高时的斜坡单元圆方差阈值D和斜坡单元最小面积阈值s参数组合为最佳参数组合，所对应的斜坡单元组为最佳斜坡单元划分结果。
在本发明的上述技术方案中，进一步地，步骤(2)中最佳初始子流域盆地累计流量阈值FTOPT的确定方法：预先定义多个累计流量阈值FT和斜坡单元最小面积阈值s，将不同累计流量阈值FT和斜坡单元最小面积阈值s进行交叉组合；依据不同参数组合对步骤(1)选取的典型区域执行划分，搜索和计算每个子流域盆地的面积和累计流量，当子流域盆地面积小于最小面积阈值s或者累计流量小于累计流量阈值FT时，则停止斜坡单元划分迭代，作为候选斜坡单元，由此获取多个斜坡单元组；将获取的斜坡单元与步骤(1)手动绘制的斜坡单元进行比较，获取的斜坡单元与手动绘制的斜坡单元最为接近时，对应的累计流量阈值确定为最佳初始子流域盆地累计流域阈值FTOPT。
进一步地，获取多组斜坡单元组，以步骤(1)手动绘制的斜坡单元作为参照标准，计算候选斜坡单元划分的总体精度，计算公式如下：
Figure BDA0003270637020000044
式中，OA为总体精度，TP表示正确分类的斜坡单元面积，TN表示正确分类的非斜坡单元面积，FP表示错误分类的斜坡单元面积，FN表示错误分类的非斜坡单元面积。当OA取值最高时，获取的斜坡单元与步骤(1)手动绘制的斜坡单元最为接近，对应的累计流量阈值为最佳初始子流域盆地累计流域阈值FTOPT。
在本发明的上述技术方案中，进一步地，步骤(2)中，在迭代细分环节，将累计流量阈值进行折减，在第i次迭代中，折减的累计流量阈值为：
Figure BDA0003270637020000051
式中：FTi为第i次迭代时的累计流量阈值，FTi-1为上一步迭代时的累计流量阈值，r为折减系数，建议取值10以上。虽然折减系数取值较大时，会导致FTi的缓慢降低和迭代次数的增加，但能够精细化控制子流域盆地的细分规模，在视觉上产生更好的结果。
在本发明的上述技术方案中,进一步地，步骤(4)中斜坡单元清理面积阈值取值建议为1～6万m2；狭长畸形斜坡单元的定义为最窄处小于100m的斜坡单元。
与现有技术相比，本发明具有以下有益效果：
1.本发明弥补了斜坡单元自动精细化划分的技术空白，基于数字高程模型和少量预定义的参数，在广泛复杂地理区域建立斜坡单元自动划分的通用框架，推动了地貌数字化和计算机科学深度融合，能够为山区滑坡灾害敏感性评价、分布式水文建模提供可靠的分析单元。
2.本发明在执行斜坡单元细分的同时，遵循数字地貌最大内部均一性和外部异质性原则，基于内部同质性指数和外部异质性指数，建立斜坡单元划分质量的客观定量评价方法，通过优选实现斜坡单元的最佳细分。
附图说明
图1为本发明所述斜坡单元自动划分方法流程示意图；
图2为子流域盆地划分示意图；
图3为最佳初始子流域盆地确定流程；
图4为基于最佳初始子流域盆地的斜坡单元细分及优选流程；
图5为基于坡向圆方差控制的斜坡单元划分示意图(图5a表示坡向散度为37°，即圆方差为0.2时获取的斜坡单元。图5b表示坡向散度为64°，即圆方差为0.55时获取的斜坡单元)；
图6为实施例的空间位置示意图；
图7为典型区域滑坡和手绘斜坡单元面积统计图；
图8为不同初始子流域盆地下斜坡单元划分结果；
图9为典型区域不同累计流量阈值下斜坡单元划分精度；
图10为不同参数组合下斜坡单元细分结果；
图11为斜坡单元划分质量评价图；
图12为最佳斜坡单元划分结果。
具体实施方式
下面通过实施例对本发明做进一步说明。有必要指出，以下实施例只用于对本发明作进一步说明，不能理解为对本发明保护范围的限制，所属领域技术人员根据上述发明内容，对本发明做出一些非本质的改进和调整进行具体实施，仍属于发明保护的范围。
实施例1
本实施划分的目标区域为青藏高原玉曲河中下游流域。玉曲河是怒江中游左岸一级支流，年径流量约36亿m3，水力资源理论蕴藏量2510MW，在玉曲河干流中下游金达-河口河段规划8级水电站，同时中下游流域内还分布有318国道、左贡县城及部分乡镇。受青藏高原隆升的影响，河流切割强烈，发育典型的高山峡谷地貌，地表破碎裸露，滑坡等地质灾害频发，给区域经济发展和基础设施建设造成了巨大的威胁。在高山峡谷地区，实施斜坡单元细分工作，能够为区域滑坡灾害敏感性分析及风险评价提供高质量的评价单元。
本实施例是在GNU/Linux系统下GRASS GIS软件平台，通过自主编写程序实现，具体包括以下内容：
(1)DEM数据和参数准备
获取覆盖示例实施区精度为12.5m的ALOS DEM数据，总面积5306km2。利用GoogleEarth影像，结合现场调查，对玉曲流域滑坡灾害进行了遥感解译，滑坡灾害数量1313个，滑坡累计面积76.75km2。同时，利用GoogleEarth三维立体影像，依据人工视觉选取典型区域，手动绘制斜坡单元作为斜坡单元样本集，累计斜坡单元面积979.41km2。面积统计如图7所示，滑坡面积区间在1～12万m2区间，主要集中在4～10万m2。手绘斜坡单元面积分布在14～350万m2，主要集中在100～200万m2区间。综合考虑两类数据集，斜坡单元最小面积阈值为15～50万m2，按照5万m2为间隔，最小面积阈值分别设置为15万、20万、25万、30万、35万、40万、45万、50万m2，最大限定面积为350万m2。
(2)生成最佳初始子流域盆地
将DEM数据中每个栅格单元标记为流量单元，按照DEM表示的区域地形特征，和由高到低的汇流方向，累加每个流量单元周围8个方向内产生的汇聚栅格单元数量，并逐步累加到下一个流量单元，获取区域汇流累计流量。示例实施区累计流量分布区间为0.5～3000万m2。采取试算的策略确定最佳初始子流域盆地累计流量阈值，设置累计流量阈值依次为50万、75万、100万、125万、150万和170万(单位m2)，斜坡单元最小面积阈值为15万、20万、25万、30万、35万、40万、45万、50万m2。将累计流量阈值和斜坡单元最小面积阈值进行交叉组合，共记60组不同参数。
采用上述参数组合对典型区域进行斜坡单元划分，生成初始子流域盆地及斜坡单元。搜索和计算每个子流域盆地的面积和累计流量，当子流域盆地面积小于最小面积阈值s或者累计流量小于累计流量阈值FT时，则停止斜坡单元划分迭代，作为候选斜坡单元，由此获取多组斜坡单元组。图8为斜坡单元最小面积阈值30万m2时，不同累计流量阈值FT所对应的初始子流域盆地的斜坡单元划分结果。随着累计流量阈值FT的逐渐增加，所生成的子流域盆地规模越大，盆地数量越来越小。在相同面积条件下，FT＝50万m2时，子流域盆地个数为134，FT＝100万m2时，子流域盆地个数为103。而当FT＝150万m2时，子流域盆地个数仅为63。
不同初始子流域盆地对斜坡单元的划分也会产生深刻的影响。手动绘制的斜坡单元如图8(a)所示，提取不同累计流量阈值对应位置处的斜坡单元划分结果：当累计流量阈值FT较小时，产生了密集的初始子流域盆地，斜坡单元自动划分结果与手动结果相比，表现为细小破碎的多边形，如图8(b)、(c)。随着FT阈值的逐渐增大，自动划分的斜坡单元规模逐渐增大，如图8(d)、(e)，直至超过手动划分的结果，如图8(f)。
为了定量研究不同初始子流域盆地对斜坡单元规模划分的影响，将自动划分的斜坡单元与手动划分的斜坡单元进行比较，采用公式(1)计算自动划分斜坡单元的总体精度OA(即正确划分的斜坡和非斜坡单元的面积与参与划分区域面积的比值)。
Figure BDA0003270637020000081
式中，OA为总体精度，TP表示正确分类的斜坡单元面积，TN表示正确分类的非斜坡单元面积，FP表示错误分类的斜坡单元面积，FN表示错误分类的非斜坡单元面积。当OA取值最高时，获取的斜坡单元规模与实际情形最为接近，对应的累计流量阈值为最佳初始子流域盆地累计流域阈值。图9为典型区域不同累计流量阈值FT下斜坡单元划分精度，可见累计流量阈值FT对斜坡单元划分结果的准确性有主导性影响，当累计流量阈值FT为100万m2时，斜坡单元划分精度最高。因此，确定最佳初始子流域盆地累计流量阈值FTOPT为100万m2，生成整个示例实施区的最佳初始子流域盆地。
(3)基于最佳初始子流域盆地的斜坡单元细分
定义斜坡单元圆方差D为地貌特征中每个栅格单元坡向的法向量在平均值上的变化，介于0～1之间，计算公式为：
Figure BDA0003270637020000082
式中：D为圆方差，NV为每个流域盆地栅格单元的数量，R为栅格单元坡向单位向量绝对值的总和。通过坡向均一性标准控制斜坡单元内部地貌一致性特征。坡向的圆方差越小，说明该斜坡单元内部所有栅格单元的坡向越接近，即为具有均一朝向的斜坡面。
对步骤(2)确定的最佳初始子流域盆地，在预先定义的斜坡单元最小面积阈值s和圆方差阈值D组合下执行迭代细分，划分为两个以上的下一级子流域盆地；当子流域盆地面积小于斜坡单元最小面积阈值s，或子流域盆地圆方差小于圆方差阈值D时，则判断该子流域盆地符合预定义的斜坡单元标准，作为候选斜坡单元，不再执行细分；否则该子流域盆地继续作为初始子流域盆地，返回迭代细分环节，继续执行细分；当所有子流域盆地符合预定义的斜坡单元标准时，迭代中止，斜坡单元细分完成，生成的斜坡单位为一个斜坡单元组。其中圆方差阈值D为0.1、0.15、0.20、0.25、0.30、0.35、0.40、0.45、0.50、0.55、0.60，斜坡单元最小面积阈值s分别设置为15万、20万、25万、30万、35万、40万、45万、50万m2，最大限定面积为350万m2，上述参数进行交叉组合，得到88组参数组合，相应地生成88个斜坡单元组。在迭代细分环节，将累计流量阈值进行折减，在第i次迭代中，折减的累计流量阈值为：
Figure BDA0003270637020000091
式中：FTi为第i次迭代时的累计流量阈值，FTi-1为上一步迭代时的累计流量阈值，r为折减系数，取值为10。
(4)斜坡单元清理
设置斜坡单元清理面积阈值为6万m2，对每组参数组合细分的斜坡单元组进行清理：对候选斜坡单元细分中产生的细碎多边形进行删除，合并填充到相邻较大的候选斜坡单元中；依据经验和实际勘测，将狭长畸形斜坡单元定义为最窄处小于100m的斜坡单元，对符合这些特征的斜坡单元进行识别和合并。
图10展示了使用不同斜坡单元圆方差阈值D和斜坡单元最小面积阈值s组合中的9例细分结果，从左上角到右下角显示了斜坡单元由最粗到最细的斜坡单元划分。显而易见，斜坡单元的面积和规模随着圆方差阈值的增大显著增大，斜坡单元最小面积阈值的影响则没有那么敏感。从视觉显示效果来看，当圆方差在0.35～0.60，斜坡单元最小面积阈值在35万～50万m2时，能够取得较好的斜坡单元细分结果。最佳的斜坡单元细分结果则需要在下一环节进行定量评价和选择。
(5)斜坡单元质量评价与优选
最佳的斜坡单元应该保证斜坡单元内部地貌特征具有最大一致性，而不同斜坡单元之间显示最大异质性。在同质性合并和异质性分割时，同一种分类共享某些相同的像素特征，而不同分类间像素特征显著不同，将DEM栅格单元视为要分割的对象，地貌坡向等同于像素值，则可以建立斜坡单元分类质量的定量评价指标。基于以上假设，定义自相关指数S为内部同质性指数，定义L为外部异质性指数：
Figure BDA0003270637020000092
Figure BDA0003270637020000093
式中：S为内部同质性指数；在每个斜坡单元组中将每个斜坡单元标记为i，第n个斜坡单元的面积标记为Sn；第n个斜坡单元的圆方差标记为dn；L是外部异质性指数；αn为第n个斜坡单元的平均坡向；
Figure BDA0003270637020000104
为研究区栅格单元的平均坡向；wnl为空间临近度指标，当编号为n和编号为l的斜坡单元相邻时，该值为1，否则该值为0。
同时，地貌特征坡向由角度表示，并不能直接用于数值计算，需要进行三角函数转换，转换格式如下：
Figure BDA0003270637020000101
Figure BDA0003270637020000102
输入参数的最佳组合是较低的内部同质性指数和较高的外部异质性指数的结合，综合考虑两种细分指标，通过以下函数进行量化：
Figure BDA0003270637020000103
式中：Smax、Smin分别是所有斜坡单元组同质性指数的最大值和最小值，Lmax、Lmin分别是所有斜坡单元组异质性指数的最大值和最小值。
当F(S，L)取值最高时的斜坡单元圆方差阈值D和斜坡单元最小面积阈值s参数组合为最佳参数组合，所对应的斜坡单元组为最佳斜坡单元划分结果。
按照公式(4)～(8)对88个斜坡单元组的划分质量进行定量评价。评价结果如图11所示，当划分参数较小时(圆方差阈值0.1～0.20，最小面积阈值15～25万m2)，划分的斜坡单元表现出较好的内部均一性，但是不满足相邻单元之间所要求的异质性。当划分参数过大，获取的斜坡单元规模较大，则不能满足所要求的内部一致性。对内部均一性和外部异质性进行综合量化评定，如图11(c)。在本发明示例实施区，最佳划分参数为圆方差阈值0.35、斜坡单元最小面积阈值为30万m2，最佳斜坡单元划分结果如图12所示。
结论
将广域地貌景观自动划分为均一的斜坡单元仍然是一个巨大的挑战。斜坡单元是排水线和汇水线界定的一种重要地形单元类型，但由于缺乏可靠的自动划分方法，在水文和地貌学研究中的使用受到了极大的限制。本发明提出一种创新性框架和方法，给定一个数字高程模型和预定义参数，实现斜坡单元的自动划分。遵循斜坡单元内部最大均一性和不同斜坡间最大异质性原则，进一步提出了斜坡单元质量评价和最佳划分参数的确定方法，通过优选得到最佳斜坡单元。在青藏高原东南缘典型高山峡谷区进行了实例测试，获取了高质量斜坡单元划分结果，为该区域滑坡灾害敏感性分析及风险评价提供了可靠的评价单元。

Claims (6)

1.一种基于汇流分割的斜坡单元自动划分方法，其特征在于，包括以下步骤：
(1)DEM数据和参数准备
准备覆盖待划分的目标区域的数字高程模型DEM，通过GIS软件获取三维山体阴影图，依据人工视觉选取典型区域，手动绘制部分斜坡单元作为斜坡单元样本集，统计斜坡单元面积的分布区间，获取斜坡单元最小面积阈值s区间和最大限定面积；
(2)生成最佳初始子流域盆地
将DEM数据中每个栅格单元标记为流量单元，按照DEM表示的区域地形特征和由高到低的汇流方向，累加每个流量单元周围8个方向内产生的汇聚栅格单元数量，并逐步累加到下一个流量单元，获取区域汇流累计流量；定义累计流量阈值FT，当栅格单元中累计流量大于FT时，该栅格单元被标记为汇水点，将所有标记的汇水点连接起来，形成汇水线；流量为0的栅格单元为局部地形高点，标记为排水点，将所有标记的排水点连接起来，形成排水线；排水线和汇水线相交形成的闭合环形为流域盆地，根据在汇水线两侧的位置，将汇水盆地划分为左右子流域盆地。依据试算得出的最佳初始子流域盆地累计流域阈值FTOPT对目标区域进行划分，生成目标区域最佳初始子流域盆地；
(3)基于最佳初始子流域盆地的斜坡单元细分
定义斜坡单元圆方差D为地貌特征中每个栅格单元坡向的法向量在平均值上的变化，介于0～1之间，计算公式为：
Figure FDA0003270637010000011
式中，D为圆方差，NV为每个流域盆地栅格单元的数量，R为栅格单元坡向单位向量绝对值的总和，通过坡向均一性标准控制斜坡单元内部地貌一致性特征；
在预先定义的斜坡单元最小面积阈值s和圆方差阈值D组合下,对步骤(2)确定的最佳初始子流域盆地执行迭代细分，划分为两个以上的下一级子流域盆地；当子流域盆地面积小于斜坡单元最小面积阈值s，或子流域盆地圆方差小于圆方差阈值D时，则判断该子流域盆地符合预定义的斜坡单元标准，作为候选斜坡单元，不再执行细分；否则该子流域盆地继续作为初始子流域盆地，返回迭代细分环节，继续执行细分；当所有子流域盆地符合预定义的斜坡单元标准时，迭代中止，斜坡单元细分完成；每种斜坡单元最小面积阈值s和圆方差阈值D组合下划分得到的斜坡单元为一个斜坡单元组；
(4)斜坡单元清理
设置斜坡单元清理面积阈值参数，对候选斜坡单元细分中产生的细碎多边形进行删除，合并填充到相邻较大的候选斜坡单元中；依据经验和实际勘测，定义和狭长畸形斜坡单元，对这些斜坡单元进行识别和合并；
(5)斜坡单元质量评价与优选
最佳的斜坡单元应该保证斜坡单元内部地貌特征具有最大一致性，而不同斜坡单元之间显示最大异质性；定义自相关指数S为内部同质性指数，定义自相关指数L为外部异质性指数：
Figure FDA0003270637010000021
Figure FDA0003270637010000022
式中：S为内部同质性指数；在每个斜坡单元组中将每个斜坡单元标记为i，第n个斜坡单元的面积标记为Sn；第n个斜坡单元的圆方差标记为dn；L是外部异质性指数；αn为第n个斜坡单元的平均坡向；
Figure FDA0003270637010000023
为研究区栅格单元的平均坡向；wnl为空间临近度指标，当编号为n和编号为l的斜坡单元相邻时，该值为1，否则该值为0；
对
Figure FDA0003270637010000024
和αn进行三角函数转换：
Figure FDA0003270637010000025
Figure FDA0003270637010000026
通过以下函数进行斜坡单元质量评价与优选：
Figure FDA0003270637010000027
式中：Smax、Smin分别是所有斜坡单元组内部同质性指数的最大值和最小值，Lmax、Lmin分别是所有斜坡单元组外部异质性指数的最大值和最小值；当F(S，L)取值最高时的斜坡单元圆方差阈值D和斜坡单元最小面积阈值s参数组合为最佳参数，所获取斜坡单元划分结果为最佳斜坡单元划分结果。
2.根据权利要求1所述方法，其特征在于，步骤(2)中最佳初始子流域盆地累计流量阈值FTOPT的确定方法如下：
预先定义多个累计流量阈值FT和斜坡单元最小面积阈值s，将不同累计流量阈值FT和斜坡单元最小面积阈值s进行交叉组合；依据不同参数组合对步骤(1)选取的典型区域执行划分，搜索和计算每个子流域盆地的面积和累计流量，当子流域盆地面积小于最小面积阈值s或者累计流量小于累计流量阈值FT时，则停止斜坡单元划分迭代，作为候选斜坡单元，由此获取多个斜坡单元组；将获取的斜坡单元规模与(1)手动绘制的斜坡单元进行对比，获取的斜坡单元与手动绘制的斜坡单元最为接近时，对应的累计流量阈值确定为最佳初始子流域盆地累计流域阈值FTOPT。
3.根据权利要求2所述方法，其特征在于，获取多个斜坡单元组后，以步骤(1)手动绘制的斜坡单元作为参照标准，计算候选斜坡单元划分的总体精度，计算公式如下：
Figure FDA0003270637010000031
式中，OA为总体精度，TP表示正确分类的斜坡单元面积，TN表示正确分类的非斜坡单元面积，FP表示错误分类的斜坡单元面积，FN表示错误分类的非斜坡单元面积；当OA取值最高时，获取的斜坡单元规模与实际情形最为接近，对应的累计流量阈值确定为最佳初始子流域盆地累计流域阈值FTOPT。
4.根据权利要求1所述方法，其特征在于，步骤(3)中，在迭代细分环节，将累计流量阈值进行折减，在第i次迭代中，折减的累计流量阈值为：
Figure FDA0003270637010000032
式中：FTi为第i次迭代时的累计流量阈值，FTi-1为上一步迭代时的累计流量阈值，r为折减系数，取值为10以上。
5.根据权利要求1所述方法，其特征在于，步骤(4)中斜坡单元清理面积阈值取值为1～6万m2。
6.根据权利要求1所述方法，其特征在于，步骤(4)中，狭长畸形斜坡单元的定义为最窄处小于100m的斜坡单元。
CN202111100800.7A 2021-09-18 2021-09-18 一种基于汇流分割的斜坡单元自动划分方法 Active CN113850822B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202111100800.7A CN113850822B (zh) 2021-09-18 2021-09-18 一种基于汇流分割的斜坡单元自动划分方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202111100800.7A CN113850822B (zh) 2021-09-18 2021-09-18 一种基于汇流分割的斜坡单元自动划分方法

Publications (2)

Publication Number Publication Date
CN113850822A true CN113850822A (zh) 2021-12-28
CN113850822B CN113850822B (zh) 2023-04-25

Family

ID=78974645

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111100800.7A Active CN113850822B (zh) 2021-09-18 2021-09-18 一种基于汇流分割的斜坡单元自动划分方法

Country Status (1)

Country Link
CN (1) CN113850822B (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114662950A (zh) * 2022-03-31 2022-06-24 中南大学 一种适用于区域地震滑坡危险性评估的斜坡单元细化方法
CN117828003A (zh) * 2024-03-05 2024-04-05 中国铁路设计集团有限公司 一种植被覆盖下坡面泥石流灾害判别方法及系统
CN118364320A (zh) * 2024-06-14 2024-07-19 中南大学 一种斜坡单元聚类划分方法、装置、终端设备及介质
CN118397378A (zh) * 2024-06-27 2024-07-26 湖南科技大学 一种坡体单元分类与提取方法、系统、设备及存储介质

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103514599A (zh) * 2013-08-30 2014-01-15 中国公路工程咨询集团有限公司 一种基于邻域总变分的影像最佳分割尺度选择方法
CN108776851A (zh) * 2018-06-11 2018-11-09 广州地理研究所 一种暴雨诱发的浅层滑坡灾害预警阈值确定方法
CN110532969A (zh) * 2019-09-02 2019-12-03 中南大学 基于多尺度图像分割的斜坡单元划分方法
CN110927360A (zh) * 2019-11-29 2020-03-27 成都理工大学 基于降雨渗流路径和含水率变化的斜坡稳定动态评价方法
JP2020060078A (ja) * 2018-10-12 2020-04-16 日本信号株式会社 土砂崩落の発生時刻を推定するシステム

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103514599A (zh) * 2013-08-30 2014-01-15 中国公路工程咨询集团有限公司 一种基于邻域总变分的影像最佳分割尺度选择方法
CN108776851A (zh) * 2018-06-11 2018-11-09 广州地理研究所 一种暴雨诱发的浅层滑坡灾害预警阈值确定方法
JP2020060078A (ja) * 2018-10-12 2020-04-16 日本信号株式会社 土砂崩落の発生時刻を推定するシステム
CN110532969A (zh) * 2019-09-02 2019-12-03 中南大学 基于多尺度图像分割的斜坡单元划分方法
CN110927360A (zh) * 2019-11-29 2020-03-27 成都理工大学 基于降雨渗流路径和含水率变化的斜坡稳定动态评价方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
HAO FANG 等: "Semi-automatic terrain slope unit division method based on human–computer interaction" *
李晓靖 等: "基于最优尺度和规则的高分辨率影像分类研究" *
王凯 等: "斜坡单元提取方法研究进展和展望" *
颜阁 等: "基于GIS的斜坡单元划分方法改进与实现" *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114662950A (zh) * 2022-03-31 2022-06-24 中南大学 一种适用于区域地震滑坡危险性评估的斜坡单元细化方法
CN117828003A (zh) * 2024-03-05 2024-04-05 中国铁路设计集团有限公司 一种植被覆盖下坡面泥石流灾害判别方法及系统
CN117828003B (zh) * 2024-03-05 2024-06-04 中国铁路设计集团有限公司 一种植被覆盖下坡面泥石流灾害判别方法及系统
CN118364320A (zh) * 2024-06-14 2024-07-19 中南大学 一种斜坡单元聚类划分方法、装置、终端设备及介质
CN118397378A (zh) * 2024-06-27 2024-07-26 湖南科技大学 一种坡体单元分类与提取方法、系统、设备及存储介质

Also Published As

Publication number Publication date
CN113850822B (zh) 2023-04-25

Similar Documents

Publication Publication Date Title
CN113850822A (zh) 一种基于汇流分割的斜坡单元自动划分方法
Evans Geomorphometry and landform mapping: What is a landform?
Alexander et al. Integrating building footprints and LiDAR elevation data to classify roof structures and visualise buildings
CN100533486C (zh) 一种数字城市全自动生成的方法
CN113283802A (zh) 一种复杂艰险山区滑坡危险性评估方法
CN104574512A (zh) 一种顾及地形语义信息的多尺度dem构建方法
CN109086649A (zh) 卫星遥感图像水体识别方法
Martinello et al. Landform classification: a high-performing mapping unit partitioning tool for landslide susceptibility assessment—a test in the Imera River basin (northern Sicily, Italy)
CN107944089B (zh) 一种基于现状视廊的地块限高分析系统及其分析方法
Arastounia Automatic classification of lidar point clouds in a railway environment
Zhang Spatiotemporal features of the three-dimensional architectural landscape in Qingdao, China
Rau et al. Semi-automatic shallow landslide detection by the integration of airborne imagery and laser scanning data
Firmansyah et al. Building typology classification using convolutional neural networks utilizing multiple ground-level image process for city-scale rapid seismic vulnerability assessment
Shannon et al. LiDAR-based sinkhole detection and mapping in Knox County, Tennessee
CN105279305B (zh) 一种地面三维激光扫描技术测站选取方法
Xiao Detecting changes in trees using multi-temporal airborne LIDAR point clouds
Doumit Digital Terrain Analysis of Lebanon: A Study of Geomorphometry
Silva et al. Hydro-geomorphological conditions for the classification of terrain susceptibility to shallow translational landslides: a geo-hydro ecological approach
Ródenas-López et al. Three-dimensional urban models in complex rural environments. Proposal for automation in the historical centre of Cehegín.
Jiang et al. Study on quantitative terrain analysis based on DEM
US20230298207A1 (en) Information processing apparatus, information processing system, information processing method, and non-transitory computer-executable medium
Rahimi et al. Investigation and forecasting of land use changes in Shahrekord city Using LCM and GIS
Ingensand et al. Automated reconstruction of 3D buildings in historic city centers from LIDAR data and 2D building footprints
Zhong et al. Landform classification based on hillslope units from DEMs
CN117272473A (zh) 基于空间三维量化系统的山地聚落保护更新设计方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant