CN108257142A - Dem中斜坡单元提取方法 - Google Patents
Dem中斜坡单元提取方法 Download PDFInfo
- Publication number
- CN108257142A CN108257142A CN201810027926.8A CN201810027926A CN108257142A CN 108257142 A CN108257142 A CN 108257142A CN 201810027926 A CN201810027926 A CN 201810027926A CN 108257142 A CN108257142 A CN 108257142A
- Authority
- CN
- China
- Prior art keywords
- ridge
- zonule
- region
- mountain valley
- closure
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 84
- 238000000605 extraction Methods 0.000 claims abstract description 39
- 230000000877 morphologic effect Effects 0.000 claims abstract description 26
- 210000004027 cell Anatomy 0.000 claims abstract description 23
- 238000004458 analytical method Methods 0.000 claims abstract description 14
- 230000035772 mutation Effects 0.000 claims abstract description 8
- 238000012876 topography Methods 0.000 claims abstract description 7
- 230000008569 process Effects 0.000 claims description 12
- 239000000284 extract Substances 0.000 claims description 7
- 238000000513 principal component analysis Methods 0.000 claims description 5
- 238000000638 solvent extraction Methods 0.000 claims description 4
- 238000006243 chemical reaction Methods 0.000 claims description 2
- 230000003252 repetitive effect Effects 0.000 claims 1
- 230000007547 defect Effects 0.000 abstract description 4
- 238000012937 correction Methods 0.000 abstract description 2
- 238000010327 methods by industry Methods 0.000 abstract 1
- 230000000052 comparative effect Effects 0.000 description 18
- 238000010586 diagram Methods 0.000 description 10
- 238000012545 processing Methods 0.000 description 9
- 230000015572 biosynthetic process Effects 0.000 description 3
- 210000000988 bone and bone Anatomy 0.000 description 3
- 230000008859 change Effects 0.000 description 3
- 238000013461 design Methods 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 230000007306 turnover Effects 0.000 description 3
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 3
- 238000013459 approach Methods 0.000 description 2
- 238000005192 partition Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 230000011218 segmentation Effects 0.000 description 2
- 239000000725 suspension Substances 0.000 description 2
- 241001269238 Data Species 0.000 description 1
- 240000001439 Opuntia Species 0.000 description 1
- 230000003321 amplification Effects 0.000 description 1
- 230000004888 barrier function Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000007123 defense Effects 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 230000003628 erosive effect Effects 0.000 description 1
- 238000000265 homogenisation Methods 0.000 description 1
- 238000003703 image analysis method Methods 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 238000001556 precipitation Methods 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 230000004044 response Effects 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
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/187—Segmentation; Edge detection involving region growing; involving region merging; involving connected component labelling
-
- 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/181—Segmentation; Edge detection involving edge growing; involving edge linking
-
- 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
-
- 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/64—Analysis of geometric attributes of convexity or concavity
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Geometry (AREA)
- Processing Or Creating Images (AREA)
- Pit Excavations, Shoring, Fill Or Stabilisation Of Slopes (AREA)
Abstract
针对现有技术提取的斜坡单元存在坡度突变缺陷的缺陷，本发明提供一种DEM中斜坡单元提取方法。方法首先确定山脊与山谷区域，其次在山脊山谷区域中沿着地势高低变化突变最集中的点勾勒具有均一性特征的闭合小区域的边界，最后合并具有相同法向量的具有均一性特征的闭合小区域形成斜坡单元。具体，首先利用栅格中心点平均曲率确定栅格单元表面形状，从而确定山脊山谷区域；其次提取山脊山谷形态骨架闭合网络中的高像素部分得到具有均一性特征的闭合小区域边界线；最后对小区域进行平面拟合与合并得到斜坡单元。本发明方法原理可靠，过程科学简便，不需预先完成DEM填洼与后期人工修正。提取的斜坡单元具有均一性特征，能够满足滑坡稳定分析均一性假定。
Description
技术领域
本发明涉及一种图像分析方法，特别是涉及一种从DEM中提取斜坡单元的方法，属于一般的图像数据处理或产生领域、环境地质灾害防治技术领域。
背景技术
预报单元的选取是区域滑坡预报技术的首要前提。恰当的预报单元选定不仅能够提高预报的精度而且有助于指导采取正确的地质灾害应急响应措施。在现有的物理预报、统计预报与成因预报三种预报技术中，预报单元的准确选取对于物理预报模型尤为重要，这是因为物理预报模型需要以预报单元作为计算单元，依据其形态特征选取相匹配的力学分析方法，方能确保模型的物理意义。
在区域滑坡预报的物理预报技术中，预报单元的准确选取存在技术障碍。主要原因在于，形状规则的预报单元是计算机进行地理图像信息分析的理想对象，但对于自然界中斜坡而言，由于坡度、坡高千差万别，因而若机械地以规则边界选取预报单元(通常是栅格单元)则所得预报单元难以反映出实际的斜坡形态，进而导致分析运算对象不明确，预报结果无法准确定位滑坡险情，造成防御工作量增大。更为重要的是，栅格单元无法反映实际的地形特征，与滑坡发生的地质环境边界没有直接关系，导致获得的结果缺乏明确的物理意义。
现有的斜坡单元提取方法基于地表水文过程，斜坡单元被定义为流域内每个子流域的左右两岸。通过正反向DEM等一系列操作，子流域两岸被分割为两个斜坡单元(10)。但是，这种方法提取的斜坡单元识别不出水流方向以外的坡度变化，导致提取的斜坡单元的坡度实际并非均一。例如，D8算法提取的水流方向仅能反映出坡度最陡的特征，而无法反映水流方向以外，尤其是坡脚或台地的坡度变化甚至是骤变。由于区域滑坡评估预报物理模型(无论基于一维无限斜坡模型还是二维极限平衡理论)均建立在条件均一化假定基础上，即假定斜坡单元应当具有均一的坡度和坡向，同时坡体为均质体。因而提取出的斜坡单元的实际非均一性与区域滑坡稳定性分析的均一性根本假定相矛盾，也就动摇了区域滑坡物理预报技术的基础。
发明内容
本发明的目的就是针对现有技术的不足，提供一种能够满足滑坡稳定分析均一性假定的斜坡单元提取方法。
为实现上述目的，本发明的技术方案如下：
一种DEM中斜坡单元提取方法，用于从原始DEM中提取斜坡单元，其特征在于：首先确定山脊与山谷区域，其次在山脊与山谷区域中沿着地势高低变化突变最集中的点勾勒出具有均一性特征的闭合小区域的边界，最后合并具有相同法向量的具有均一性特征的闭合小区域形成斜坡单元。
上述DEM中斜坡单元提取方法，是用于从原始DEM中提取斜坡单元，所获得的斜坡单元的实际地形特征能够满足滑坡稳定分析均一性假定条件。其主要技术原理在于：(1)将满足均一性假定条件的斜坡单元定义为三维空间中的一个连续、均质且闭合的区域。该区域具有均一的坡度与坡向属性，区域内部的地形特征不存在起伏变化、无坡度突变现象，同时区域边界确定为实际地形特征出现微小起伏变化之处，即坡度或坡向发生微小转折之处。这样，在DEM数据中寻找出坡度或坡向发生微小转折之处就能够定义出斜坡单元的边界所在。(2)提取斜坡单元的基础是搜寻出具有均一性特征的闭合小区域。在搜寻得到具有均一性特征的闭合小区域的结果上，对每个小区域包含的空间点集进行平面拟合，以使数目众多的均质小区域合并形成斜坡单元。具体借助法向量指标合并具有相同法向量的相邻小区域，便可最终形成所需的斜坡单元。 (3)山脊区域与山谷区域是坡度、坡向转折或突变的集中区域，故具有均一性特征的闭合小区域的边界线即存在于这两个区域中。利用DEM数据可以确定山脊与山谷区域。在已确定的山脊与山谷区域中，沿着地势高低变化突变最集中的点即可勾勒出具有均一性特征的闭合小区域的边界。
上述DEM中斜坡单元提取方法中，在利用DEM数据确定山脊与山谷区域时，实际地形的微小起伏特征可通过栅格单元的凸凹表征。因此可将地形被离散为栅格中心点集，计算每个栅格点的平均曲率k，利用栅格中心点平均曲率k判断每个栅格单元的凸凹特征。对于每一栅格单元，k＞0的栅格单元表面为凸形(对应山脊区域)、k＜0的栅格单元表面为凹形(对应山谷区域)、k＝0的栅格单元为平面(对应平地区域)。这样就借用栅格中心点平均曲率k的物理意义将实际地形被分为山脊区域和山谷区域两类。每个栅格点的平均曲率k可依式1计算确定：
式中：k为平均曲率，kpl为平面曲率，kpr为剖面曲率。
上述DEM中斜坡单元提取方法中，在山脊与山谷区域中沿着地势高低变化突变最集中的点勾勒出具有均一性特征的闭合小区域的边界时，可以采用形态学方法在已确定的山脊与山谷区域中提取小区域边界。形态学将DEM作为灰度图像处理，实际地形起伏由灰度像素反映，灰度像素值越大，代表所处地势越高。具体可利用形态学骨架线。由于形态学骨架线能够反映区域内部结构与区域轮廓的几何特征，因而构建山脊与山谷区域的形态骨架可以反映实际微小地形的起伏特征。在此基础上，再将山脊与山谷形态骨架连接成闭合网络后，网络中每个区域涵盖的实际地形就具有了均一的坡度与坡向，网络线即可近似等效于具有均一性特征的闭合小区域的边界线。但在骨架线连接成闭合网络后，网络与其内部的栅格点仍然是各自独立的(可以分别理解为GIS中的线矢量图层与点矢量图层)，尚需进一步提取。
在骨架线连接成闭合网络基础上，进一步提取山脊山谷闭合网络中高像素部分，即得到具有均一性特征的闭合小区域的边界线。具体可采用形态学分水岭分割算法提取每个小区域的边界，并对每个小区域赋给唯一值，每个唯一值代表一个坡度坡向均一的闭合区域，得到具有唯一值的整型栅格数据，闭合区域边界线即为具有均一性特征的闭合小区域的边界线。分水岭分割算法一般用于提取二值或灰度图像边界，进行区域图像标注。其基本原理是将图像中像素较高的部分视为高海拔区，像素低的部分视为低海拔区，则图像中每个像素局部极小值影响区域视作一个“集水盆”，“集水盆”边界即为待提取的“分水岭”。在形成的山谷山脊闭合网络二值图像中，像素较高的部分即为待提取的小区域边界。这样，采用分水岭分割后，得到一种具有唯一值的整型栅格数据，每个唯一值代表一个坡度坡向均一的闭合区域。
上述DEM中斜坡单元提取方法中，在合并具有相同法向量的具有均一性特征的闭合小区域形成斜坡单元时，可以采用主成分分析法得到每个小区域相应的拟合平面与法向量，比较相邻的两个小区域的法向量相似度，若满足向量相似度准则即进行合并。合并过程迭代进行。由此小，区域逐渐合并为大区域，最终得到斜坡单元。
与现有技术相比，本发明的有益效果是：本发明方法自原始DEM中提取的斜坡单元具有更为均一的坡度和坡向，克服了现有方法提取的斜坡单元存在的坡度突变缺陷。在方法实现中，只需载入原始DEM数据，即可直接输出斜坡单元栅格文件，不需预先进行DEM填洼与后期人工修正。提取所得的斜坡单元具有均一性特征，能够满足滑坡稳定分析均一性假定。
附图说明
图1是蒋家沟区域地理位置与DEM。
图2是DEM中斜坡单元提取方法流程框图。
图3是斜坡单元划分过程主要步骤处理局部效果图(Ⅰ：DEM，Ⅱ：平均曲率计算，Ⅲ：山脊山谷区域二值图像(黑色为山脊区域，白色为山谷区域)，Ⅳ：山脊山谷线形态学骨架(黑色线条为山脊线骨架，灰色线条为山谷线骨架)，Ⅴ：搜索模式举例，Ⅵ：闭合的山谷山脊骨架(灰色线条为形态骨架网络，黑色线条为小区域合并后的结果))
图4是蒋家沟流域斜坡单元划分结果；图4(a)是对比例一方法提取结果，图4(b)是实施例一方法提取结果。
图5是栅格点簇处理过程示意图(n为法向量；α为坡度；为坡向； kpl为平面曲率，kpr为剖面曲率)；图5(a)是实际地形离散为栅格点集，图5(b)是离散点曲率和法向量计算。
图6是最大圆盘法提取形态骨架示意图；图6(a)示最大圆盘的不同位置，图6(b)示完整的骨架。
图7是蒋家沟流域山脊山谷区域形态学骨架图；图7(a)示山谷区域形态学骨架，图7(b)示山脊区域形态学骨架。
图8所示为蒋家沟流域分水岭算法处理结果(闭合区域数目：7686)。
图9合并小区域方法示意图；图(a)是小区域理想拟合平面，图9(b) 是拟合平面局部放大效果图。
图10是地表水文过程分析方法提取斜坡单元示意图。
图11是实施例一方法与对比例一方法结果对比图；图11(a)是对比例一方法坡度标准差，图11(b)是实施例一方法坡度标准差，图11(c)是对比例一方法坡向标准差，图11(d)是实施例一方法坡向标准差。
具体实施方式
下面结合附图，对本发明的优选实施例作进一步的描述。
实施例一
如图1～图9示，选取中国云南省蒋家沟流域作为研究对象，采用本发明方法提取该流域斜坡单元。
图1是蒋家沟区域地理位置与DEM。蒋家沟位于昆明市东北部，东经103°6′～103°13′，北纬26°13′～26°17′。流域面积47.1km2，长度约12.1km。蒋家沟流域年平均降水900mm，其中雨季(5～10月) 降雨量较大，属于典型的浅层活跃滑坡与泥石流多发地区。
图2是DEM中斜坡单元提取方法流程框图，图3是斜坡单元划分过程主要步骤处理局部效果图(Ⅰ：DEM，Ⅱ：平均曲率计算，Ⅲ：山脊山谷区域二值图像(黑色为山脊区域，白色为山谷区域)，Ⅳ：山脊山谷线形态学骨架(黑色线条为山脊线骨架，灰色线条为山谷线骨架)，Ⅴ：搜索模式举例，Ⅵ：闭合的山谷山脊骨架(灰色线条为形态骨架网络，黑色线条为小区域合并后的结果))。
首先提取斜坡DEM数据。提取过程基于7m精度的蒋家沟流域DEM，栅格数量共计208个。提取结果如图4b所示。图4是蒋家沟流域斜坡单元划分结果，图4(b)是实施例一方法提取结果(单元数：208)。以DEM 数据为基础，依如下步骤实施：
步骤S1、确定DEM栅格单元表面形状
依式1计算每个栅格中心点的平均曲率k。图5是栅格点簇处理过程示意图(n为法向量；α为坡度；为坡向；kpl为平面曲率，kpr为剖面曲率)。图5(a)是实际地形离散为栅格点集，图5(b)是离散点曲率和法向量计算。部分计算结果如表1所示。
表1栅格中心点的平均曲率k(部分数据)
根据k值判断栅格单元表面形状，判断条件是k＞0的栅格单元表面为凸形、k＜0的栅格单元表面为凹形、k＝0的栅格单元为平面。
将曲率结果进行二值化处理。结果如主要步骤处理局部效果图Ⅲ所示，黑色部分为山脊和平地区域，白色部分为山谷区域。
步骤S2、提取具有均一性特征的闭合小区域的边界线
步骤S21、绘制山脊山谷形态骨架闭合网络
分别提取山脊区域、山谷区域形态骨架线，将形态骨架连接起来，构成山脊山谷形态骨架闭合网络。
首先利用形态学腐蚀运算与开运算操作，结合MATLAB提供的图像处理工具箱，分别提取山脊山谷区域形态骨架线，以形态骨架反映山脊山谷区域的几何特征。本实施方式采用最大圆盘法分别提取山脊、山谷区域形态骨架线。图6是最大圆盘法提取形态骨架示意图，图6(a)示最大圆盘的不同位置，图6(b)示完整的骨架。最大圆盘法是用一系列直径可变化的圆盘与图像轮廓相切，如果图像内部找不到更大的圆能够完全包含这些圆，连接这些圆的圆心，即构成了骨架线。主要步骤处理局部效果图Ⅳ所示为处理结果效果图。图7是蒋家沟流域山脊山谷区域形态学骨架图，图7(a)示山谷区域形态学骨架，图7(b)示山脊区域形态学骨架。提取结果显示提取的山脊山谷区域形态骨架可以反映实际微小地形的起伏特征。如此采用一定搜索方式将山脊山谷形态骨架连接成闭合网络后，网络中每个小区域涵盖的实际地形就具有了均一的坡度与坡向。本实施例提取结果图局部如主要步骤处理局部效果图Ⅳ所示。
而后，将形态骨架采用一定方式连接起来，构成山脊山谷形态骨架闭合网络。具体采用迭代搜索的方法构建骨架闭合网络，操作如图3Ⅴ所示，例如，P点山脊骨架线的一个端点，以P点方向向量±45°范围内做圆弧，圆心角Sθ，半径Sr，构成弧形搜索区域。a1、a2、a3为搜索区域中山谷线的中间节点，b1为山谷线的一个端点，寻找a1、a2、a3、b1四点中距离P最近的点。如图3Ⅴ所示，最近点为a1，连接P与a1，即在山谷线和山脊线之间建立了一个连接。如果搜索区域内没有可供连接的点，则扩大搜索范围继续搜索，直到达到设置的搜索范围阈值(对应的中止条件1)。本实施方式中，初始搜索区域是圆心角Sθ ＝90o，搜索半径 Sr为35m(所用DEM栅格分辨率的5倍)。如果初始搜索区域内没有可供连接的点，则以70m为步长(所用DEM栅格分辨率的10倍)增加搜索半径，直到达到设定的搜索半径阈值。本实施方式中，搜索半径阈值设为DEM栅格分辨率的400倍。如此对所有的端点进行循环，即可构建山谷山脊骨架线闭合网络。运算完毕得到山脊、山谷闭合网络二值图像，处理效果如图3Ⅵ所示。
步骤S22、提取具有均一性特征的闭合小区域的边界线
采用形态学分水岭分割算法提取山脊、山谷闭合网络二值图像中高像素部分，得到具有均一性特征的闭合小区域的边界线。
在上一步骤形成的山谷山脊闭合网络二值图像中，采用分水岭算法提取山脊山谷闭合网络中高像素部分，得到每个小区域的边界。再对每个小区域赋给唯一值，分割后得到一种具有唯一值的整型栅格数据，每个唯一值代表一个坡度坡向均一的闭合区域。图8所示为蒋家沟流域分水岭算法处理结果(闭合区域数目：7686)。
步骤S3、形成斜坡单元
采用主成分分析法(PCA)得到每个具有均一性特征的闭合小区域的法向量，比较相邻的两个小区域的法向量相似度，若满足向量相似度准则即进行合并，合并过程迭代进行，最终得到斜坡单元。
首先，采用主成分分析方法提取各个具有均一性特征的闭合小区域所包含的空间点集的拟合平面，以小区域包含的离散点到平面的垂直距离作为拟合平面质量的控制标准。得到每个小区域相应的拟合平面和法向量。然后，通过向量相似度准则，比较相邻小区域的拟合平面的法向量相似度，进行相邻区域的合并。
本实施方式中，将相邻的两个拟合平面的法向量归一化(转换为单位向量)，然后计算其向量相似度。如最后得到的向量相似度在设定阈值范围内(一般接近1)，则说明向量相似程度足以达到进行合并的要求，即可以进行小区域的合并。这一过程迭代进行，使满足要求的小区域逐渐合并为大区域。当合并后的斜坡单元面积达到设置的面积阈值(依绘图所需的斜坡单元精度确定)时，停止合并(对应图2的中止条件2)，最终得到斜坡单元。如图9(a)是小区域理想拟合平面，图9(b)是拟合平面局部放大效果图。
实际操作时，整体依上述提取过程，采用MATLAB编制斜坡单元自动提取程序。程序读入tiff格式的原始DEM栅格数据后，经计算可自动输出tiff格式的斜坡单元栅格数据(图2末尾框)。再将其导入ARCGIS，采用工具箱中的“栅格转面”等工具进行转换，即可得到斜坡单元面矢量数据。
对比例一
如图10所示，选取中国云南省蒋家沟流域作为研究对象，采用基于地表水文过程分析的提取方法提取该流域斜坡单元。
图10是地表水文过程分析方法提取斜坡单元示意图(图中数字1、 2、3代表集水区，LSU代表左侧斜坡单元、RSU代表右侧斜坡单元)。
提取过程基于7m精度的蒋家沟流域DEM，提取的单元数目与实例一接近。图4(a)是对比例一方法提取结果(单元数：203)。
提取方法具体采用Turel和Frost(2011)提供的地表水文过程分析法。 (见参考文件1)
实施例一与对比例一结果对比
图4显示对比例一方法提取出的大部分单元均出现坡度值失真问题，在DEM中部和下部表现的尤为明显。所谓坡度值失真，指一个单元内存在多种明显不同的坡度，说明单元包含的实际地形存在坡度的突变。实施例一方法提取的单元具有相对均一的坡度，很好地克服了坡度值失真问题。
图11给出了实施例一方法与对比例一方法下坡度、坡向标准差的对比，以直方图形式表示。其中图11(a)是对比例一方法坡度标准差，图11(b) 是实施例一方法坡度标准差，图11(c)是对比例一方法坡向标准差，图11(d) 是实施例一方法坡向标准差。根据坡度标准差范围进行了坡度失真程度的分级，如表2所示。斜坡单元的坡度、坡向标准差值越小，说明该单元内部坡度、坡向变化幅度越小，越不易出现坡度、坡向的失真问题。
首先由图11(a)、(b)、(c)、(d)可知，对比例一提取坡度标准差范围0°～18°，坡向标准差范围0°～180°；实施例一方法提取坡度标准差范围0°～13.5°，坡向标准差范围0°～157.6°。说明对比例一方法总体存在着单元内部坡度、坡向波动较大的现象。其次，如表2所示，对比例一方法坡度标准差在(8°,14°)的范围内的单元数所占比例为51.25％，实施例一方法坡度标准差在(8°,14°)的范围内的单元数所占比例为20.19％。对比例一方法远大于实施例一方法相应标准差范围内的结果。对比结果表明，对比例一方法大部分单元(比例57.16％) 出现失真问题，且部分单元(比例5.91％)出现坡度值严重失真问题。而实施例一方法中，大部分单元(接近80％)没有出现坡度值失真问题。
表2坡度失真程度分级
基于上述的提取结果对比分析，可以发现：实施例一方法提取的斜坡单元具有更为均一的坡度和坡向，克服了基于对比例一方法提取的斜坡单元存在的坡度突变缺陷。程序只需载入原始DEM，即可直接输出斜坡单元栅格文件，不需预先进行DEM填洼与后期人工修正。实施例一方法提取的斜坡单元可为区域滑坡分析模型计算单元提供一种更好的选择。
参考文件1：Turel M,Frost JD(2011)Delineation of slope profiles fromdigital elevation models for landslide hazard analysis.American Society ofCivil Engineers(224),829-836.DOI:10.1061/41183(418)87。
Claims (10)
1.DEM中斜坡单元提取方法，用于从原始DEM中提取斜坡单元，其特征在于：首先确定山脊与山谷区域，其次在山脊与山谷区域中沿着地势高低变化突变最集中的点勾勒出具有均一性特征的闭合小区域的边界，最后合并具有相同法向量的具有均一性特征的闭合小区域形成斜坡单元。
2.根据权利要求1所述的方法，其特征在于：在所述确定山脊与山谷区域时，计算每个栅格中心点的平均曲率k，根据k值判断栅格单元表面形状，判断条件是k＞0的栅格单元表面为山脊区域、k＜0的栅格单元表面为山谷区域、k＝0的栅格单元为平地区域。
3.根据权利要求1所述的方法，其特征在于：在所述在山脊与山谷区域中沿着地势高低变化突变最集中的点勾勒出具有均一性特征的闭合小区域的边界时，首先构建山脊与山谷区域的形态骨架，并将山脊与山谷形态骨架连接成闭合网络；再提取山脊山谷闭合网络中高像素部分，即得到具有均一性特征的闭合小区域的边界线。
4.根据权利要求3所述的方法，其特征在于：首先将图像二值化处理，再构建山脊与山谷区域的形态骨架。
5.根据权利要求4所述的方法，其特征在于：在所述提取山脊山谷闭合网络中高像素部分时，采用形态学分水岭分割算法提取每个具有均一性特征的闭合小区域的边界，再对每个小区域赋给唯一值，得到具有唯一值的整型栅格数据，每个唯一值代表一个坡度坡向均一的闭合区域。
6.根据权利要求1所述的方法，其特征在于：在所述合并具有相同法向量的具有均一性特征的闭合小区域形成斜坡单元时，采用主成分分析法得到每个具有均一性特征的闭合小区域的法向量，比较相邻的两个小区域的法向量相似度，若满足向量相似度准则即进行合并得到斜坡单元。
7.根据权利要求1所述的方法，首先提取斜坡DEM数据，其特征在于：然后依如下步骤实施：
步骤S1、确定DEM栅格单元表面形状
计算每个栅格中心点的平均曲率k，根据k值判断栅格单元表面形状，判断条件是k＞0的栅格单元表面为凸形、k＜0的栅格单元表面为凹形、k＝0的栅格单元为平面；
将曲率结果进行二值化处理，黑色部分为山脊和平地区域，白色部分为山谷区域；
步骤S2、提取具有均一性特征的闭合小区域的边界线
步骤S21、绘制山脊山谷形态骨架闭合网络
分别提取山脊区域、山谷区域的形态骨架线，将形态骨架连接起来，构成山脊山谷形态骨架闭合网络；
步骤S22、提取具有均一性特征的闭合小区域的边界线
用形态学分水岭分割算法提取山脊山谷闭合网络二值图像中高像素部分，得到具有均一性特征的闭合小区域的边界线；
步骤S3、形成斜坡单元
采用主成分分析法确定每个具有均一性特征的闭合小区域的法向量，比较相邻的两个小区域的法向量相似度，若满足向量相似度准则即进行合并；合并过程迭代进行，最终得到斜坡单元。
8.根据权利要求7所述的方法，其特征在于：所述步骤11中，依式1计算平均曲率：
式中：k为平均曲率，kpl为平面曲率，kpr为剖面曲率。
9.根据权利要求7所述的方法，其特征在于：所述步骤S21中，采用迭代搜索的方法构建骨架闭合网络，具体是首先确定山脊骨架线的一个端点P，以P点方向向量±45°范围内做圆弧，圆心角Sθ、半径Sr，构成扇形搜索区域；寻找扇形搜索区域中山谷线的中间节点与山谷线端点中距离P点最近的点a，连接P与a在山谷线与山脊线之间建立连接；重复操作直至达到设置的搜索范围阈值；所述搜索范围阈值依绘图所需精度确定。
10.根据权利要求7所述的方法，其特征在于：所述步骤S3中，采用主成分分析方法提取各个小区域包含的空间点集的拟合平面，以小区域包含的离散点到平面的垂直距离作为拟合平面质量的控制标准；再通过向量相似度准则比较相邻小区域的拟合平面的法向量相似度，进行相邻区域的合并形成斜坡单元。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810027926.8A CN108257142B (zh) | 2018-01-11 | 2018-01-11 | Dem中斜坡单元提取方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810027926.8A CN108257142B (zh) | 2018-01-11 | 2018-01-11 | Dem中斜坡单元提取方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN108257142A true CN108257142A (zh) | 2018-07-06 |
| CN108257142B CN108257142B (zh) | 2022-02-18 |
Family
ID=62726443
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810027926.8A Expired - Fee Related CN108257142B (zh) | 2018-01-11 | 2018-01-11 | Dem中斜坡单元提取方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN108257142B (zh) |
Cited By (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110442932A (zh) * | 2019-07-19 | 2019-11-12 | 河海大学 | 一种结合切面曲率和Dinf方法的DEM水流方向模拟方法 |
| CN111191624A (zh) * | 2020-01-04 | 2020-05-22 | 圣点世纪科技股份有限公司 | 一种基于法向量信息的指静脉特征提取算法 |
| CN111739259A (zh) * | 2020-06-15 | 2020-10-02 | 中国科学院、水利部成都山地灾害与环境研究所 | 斜坡单元局部聚集性破坏判定方法、区域滑坡预警方法 |
| CN112560267A (zh) * | 2020-12-17 | 2021-03-26 | 杭州鲁尔物联科技有限公司 | 一种斜坡单元的划分方法、装置、设备和存储介质 |
| CN113012399A (zh) * | 2021-02-25 | 2021-06-22 | 中原工学院 | 一种降雨型滑坡预警方法及系统 |
| CN113554750A (zh) * | 2021-06-09 | 2021-10-26 | 贵州师范学院 | 基于gis的喀斯特槽谷顺逆向坡自动识别方法 |
| CN114329663A (zh) * | 2021-12-27 | 2022-04-12 | 中国自然资源航空物探遥感中心 | 基于历史地质灾害规模的斜坡单元划分方法 |
| CN115457407A (zh) * | 2022-09-29 | 2022-12-09 | 云南大学 | 一种基于YoloV4算法的滑坡及形态识别方法 |
| CN115661004A (zh) * | 2022-12-26 | 2023-01-31 | 武汉天际航信息科技股份有限公司 | 三维地形模型及道路dem的更新方法、设备以及存储介质 |
| CN117828003B (zh) * | 2024-03-05 | 2024-06-04 | 中国铁路设计集团有限公司 | 一种植被覆盖下坡面泥石流灾害判别方法及系统 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103745191A (zh) * | 2013-11-15 | 2014-04-23 | 中国科学院遥感与数字地球研究所 | 一种基于地形分析的黄土地区塬梁峁自动识别方法 |
| CN106981092A (zh) * | 2017-03-28 | 2017-07-25 | 南京师范大学 | 基于Priority‑Flood的内流流域提取方法 |
-
2018
- 2018-01-11 CN CN201810027926.8A patent/CN108257142B/zh not_active Expired - Fee Related
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103745191A (zh) * | 2013-11-15 | 2014-04-23 | 中国科学院遥感与数字地球研究所 | 一种基于地形分析的黄土地区塬梁峁自动识别方法 |
| CN106981092A (zh) * | 2017-03-28 | 2017-07-25 | 南京师范大学 | 基于Priority‑Flood的内流流域提取方法 |
Non-Patent Citations (1)
| Title |
|---|
| 张路路: "基于DEM的古滑坡体提取研究", 《万方数据知识服务平台》 * |
Cited By (15)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110442932B (zh) * | 2019-07-19 | 2022-09-23 | 河海大学 | 一种结合切面曲率和Dinf方法的DEM水流方向模拟方法 |
| CN110442932A (zh) * | 2019-07-19 | 2019-11-12 | 河海大学 | 一种结合切面曲率和Dinf方法的DEM水流方向模拟方法 |
| CN111191624A (zh) * | 2020-01-04 | 2020-05-22 | 圣点世纪科技股份有限公司 | 一种基于法向量信息的指静脉特征提取算法 |
| CN111191624B (zh) * | 2020-01-04 | 2021-02-19 | 圣点世纪科技股份有限公司 | 一种基于法向量信息的指静脉特征提取方法 |
| CN111739259A (zh) * | 2020-06-15 | 2020-10-02 | 中国科学院、水利部成都山地灾害与环境研究所 | 斜坡单元局部聚集性破坏判定方法、区域滑坡预警方法 |
| CN111739259B (zh) * | 2020-06-15 | 2022-04-05 | 中国科学院、水利部成都山地灾害与环境研究所 | 斜坡单元局部聚集性破坏判定方法、区域滑坡预警方法 |
| CN112560267A (zh) * | 2020-12-17 | 2021-03-26 | 杭州鲁尔物联科技有限公司 | 一种斜坡单元的划分方法、装置、设备和存储介质 |
| CN112560267B (zh) * | 2020-12-17 | 2024-04-09 | 杭州鲁尔物联科技有限公司 | 一种斜坡单元的划分方法、装置、设备和存储介质 |
| CN113012399A (zh) * | 2021-02-25 | 2021-06-22 | 中原工学院 | 一种降雨型滑坡预警方法及系统 |
| CN113554750A (zh) * | 2021-06-09 | 2021-10-26 | 贵州师范学院 | 基于gis的喀斯特槽谷顺逆向坡自动识别方法 |
| CN114329663B (zh) * | 2021-12-27 | 2022-08-05 | 中国自然资源航空物探遥感中心 | 基于历史地质灾害规模的斜坡单元划分方法 |
| CN114329663A (zh) * | 2021-12-27 | 2022-04-12 | 中国自然资源航空物探遥感中心 | 基于历史地质灾害规模的斜坡单元划分方法 |
| CN115457407A (zh) * | 2022-09-29 | 2022-12-09 | 云南大学 | 一种基于YoloV4算法的滑坡及形态识别方法 |
| CN115661004A (zh) * | 2022-12-26 | 2023-01-31 | 武汉天际航信息科技股份有限公司 | 三维地形模型及道路dem的更新方法、设备以及存储介质 |
| CN117828003B (zh) * | 2024-03-05 | 2024-06-04 | 中国铁路设计集团有限公司 | 一种植被覆盖下坡面泥石流灾害判别方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN108257142B (zh) | 2022-02-18 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108257142A (zh) | Dem中斜坡单元提取方法 | |
| Rottensteiner et al. | A new method for building extraction in urban areas from high-resolution LIDAR data | |
| CN102013015B (zh) | 一种面向对象的遥感影像海岸线提取方法 | |
| Lohani et al. | Application of airborne scanning laser altimetry to the study of tidal channel geomorphology | |
| CN104008553B (zh) | 融合影像梯度信息和分水岭方法的裂缝检测方法 | |
| CN111696123A (zh) | 超像素分类识别的遥感图像水域分割提取方法 | |
| CN109447160B (zh) | 一种影像和矢量道路交叉点自动匹配的方法 | |
| CN106920235B (zh) | 基于矢量底图匹配的星载光学遥感影像自动校正方法 | |
| Wang et al. | A new slope unit extraction method for regional landslide analysis based on morphological image analysis | |
| CN115641327B (zh) | 一种基于大数据的建筑工程质量监理和预警系统 | |
| CN103106658A (zh) | 一种海岛、礁岸线快速提取方法 | |
| Chen et al. | Analyzing inundation extent in small reservoirs: A combined use of topography, bathymetry and a 3D dam model | |
| CN111046884B (zh) | 一种多特征辅助分水岭算法的斜坡地质灾害提取方法 | |
| CN112509134A (zh) | 一种潮滩数字高程模型构建方法及系统 | |
| Straumann et al. | Delineation of valleys and valley floors | |
| Omidalizarandi et al. | Segmentation and classification of point clouds from dense aerial image matching | |
| CN115937676A (zh) | 一种基于dem的大尺度沙丘地貌精细分类方法 | |
| Mandlburger et al. | Derivation of a countrywide river network based on Airborne Laser Scanning DEMs-results of a pilot study | |
| Moravej et al. | Comparison of automated and manual landform delineation in semi detailed soil survey procedure | |
| Mitchell et al. | Little tools, big job: The periglacial conveyor system in cryoplanated uplands | |
| Šamanović et al. | Influence of pit removal algorithms on surface runoff simulation | |
| CN116681711B (zh) | 一种分区指导下的高分辨率遥感影像多尺度分割方法 | |
| Safari | A Novel Transformation Watershed Image Segmentation Model in Digital Elevation Maps Processing | |
| Ding et al. | An object-based method for Chinese landform types classification | |
| CN112686192B (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 | ||
| CF01 | Termination of patent right due to non-payment of annual fee | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20220218 |