CN115359221A - 一种基于数字高程模型的内陆湖泊流域提取方法及系统 - Google Patents

一种基于数字高程模型的内陆湖泊流域提取方法及系统 Download PDF

Info

Publication number
CN115359221A
CN115359221A CN202210992072.3A CN202210992072A CN115359221A CN 115359221 A CN115359221 A CN 115359221A CN 202210992072 A CN202210992072 A CN 202210992072A CN 115359221 A CN115359221 A CN 115359221A
Authority
CN
China
Prior art keywords
pixel
lake
flow direction
depression
pixels
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.)
Pending
Application number
CN202210992072.3A
Other languages
English (en)
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.)
Sun Yat Sen University
Original Assignee
Sun Yat Sen 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 Sun Yat Sen University filed Critical Sun Yat Sen University
Priority to CN202210992072.3A priority Critical patent/CN115359221A/zh
Publication of CN115359221A publication Critical patent/CN115359221A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T19/00—Manipulating 3D models or images for computer graphics
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F16/00—Information retrieval; Database structures therefor; File system structures therefor
    • G06F16/20—Information retrieval; Database structures therefor; File system structures therefor of structured data, e.g. relational data
    • G06F16/29—Geographical information databases

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Databases & Information Systems (AREA)
  • Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Remote Sensing (AREA)
  • Data Mining & Analysis (AREA)
  • Computer Graphics (AREA)
  • Computer Hardware Design (AREA)
  • Software Systems (AREA)
  • Image Processing (AREA)

Abstract

本发明涉及地理信息分析技术领域，提出一种基于数字高程模型的内陆湖泊流域提取方法及系统，其中：输入目标地区的数字高程模型和湖泊分布数据；计算目标地区的数字高程模型中各像元的初始栅格流向，并选取无法确定流向的像元组成洼地像元集合；识别所述洼地像元组合中的陆面洼地和湖泊洼地，得到陆面洼地像元集合和湖泊洼地像元集合；对所述数字高程模型中相应的陆面洼地像元进行有选择地填充后计算栅格流向，得到栅格流向矩阵；对所述栅格流向矩阵进行解析并构建流向多叉树；对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线相应的边缘像元，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，提取得到湖泊流域边界。

Description

一种基于数字高程模型的内陆湖泊流域提取方法及系统
技术领域
本发明涉及地理信息分析技术领域，更具体地，涉及一种基于数字高程模型的内陆湖泊流域提取方法及系统。
背景技术
内陆湖泊是指内陆洼地中相对静止、有一定面积，不与海洋发生直接联系的水体。部分内陆湖泊没有河流出流，部分内陆湖泊的河流出流至下一个湖泊或消失于内陆。研究内陆湖泊的流域边界对了解内陆盆区的水资源分布特征、湖泊对气候变化与农业生产的影响有重要意义。但目前应用广泛的基于数字高程模型的水文信息提取程序适用于与海洋有水力联系的外流流域，而对内陆湖泊这种内流的水体则不适用，需要特殊的地形处理步骤。在实际工作中，内陆湖泊的流域边界通常依靠人为地描绘山脊线确定，当湖泊数量过多时容易出错。
发明内容
本发明为克服上述现有技术所述的缺少对内陆湖泊水文信息提取的方法，在内陆湖泊的流域边界确定时出错率高的缺陷，提供一种基于数字高程模型的内陆湖泊流域提取方法及系统。
为解决上述技术问题，本发明的技术方案如下：
一种基于数字高程模型的内陆湖泊流域提取方法，包括以下步骤：
S1、输入目标地区的数字高程模型和湖泊分布数据；
S2、计算目标地区的数字高程模型中各像元的初始栅格流向，并根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合；
S3、根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地，得到陆面洼地像元集合和湖泊洼地像元集合；
S4、对所述数字高程模型中相应的陆面洼地像元进行有选择地填充，然后计算栅格流向，得到栅格流向矩阵；
S5、对所述栅格流向矩阵进行解析并构建流向多叉树；
S6、对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线相应的边缘像元，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，提取得到湖泊流域边界。
本技术方案面向内陆盆区的地形处理与湖泊流域提取，根据内陆盆区地形特点对洼地像元进行解析和填充处理，再对栅格流向进行解析并构建流向多叉树，将树状数据结构遍历算法应用于湖泊流域提取，实现内陆湖泊等内流水体的流域边界提取。
进一步地，本发明还提出一种基于数字高程模型的内陆湖泊流域提取系统，应用上述技术方案提出的内陆湖泊流域提取方法。其中包括数据采集模块，洼地识别模块，洼地填充模块，流向多叉树构建模块，和湖泊流域边界提取模块。
其中，数据采集模块用于获取目标地区的数字高程模型和湖泊分布数据；
洼地识别模块用于通过计算目标地区的数字高程模型中各像元的初始栅格流向，并根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合；以及用于根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地，输出陆面洼地像元集合和湖泊洼地像元集合；
洼地填充模块用于对所述数字高程模型中相应的陆面洼地像元进行有选择地填充，然后计算栅格流向，输出栅格流向矩阵；
流向多叉树构建模块，用于对所述栅格流向矩阵进行解析并构建流向多叉树；
湖泊流域边界提取模块用于对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线相应的边缘像元，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，输出提取得到的湖泊流域边界数据。
进一步地，本发明还提出一种基于数字高程模型的内陆湖泊流域提取系统，包括存储器和处理器，所述存储器存储有计算机程序，其中所述处理器执行所述计算机程序时实现本发明提出的内陆湖泊流域提取方法的步骤。
与现有技术相比，本发明技术方案的有益效果是：本发明通过结合内陆湖泊等内流水体的地形特征，对目标地区的数字高程模型和湖泊分布数据进行处理和解析，根据各像元的栅格流向，对数字高程模型中像元是否为内陆湖泊进行判断，并对洼地像元进行有选择地填充，有效反映内陆湖泊流域真实的水文信息，进一步用于构建流向多叉树，将树状数据结构遍历算法应用于湖泊流域提取，实现高效率且高准确率的内陆湖泊等内流水体的流域边界提取。
附图说明
图1为本发明实施例1的基于数字高程模型的内陆湖泊流域提取方法的流程图。
图2为本发明实施例2的羌塘高原区数字高程模型与湖泊分布特征图。
图3为本发明实施例2的羌塘高原区的洼地分布示意图。
图4为本发明实施例2的羌塘高原区的湖泊洼地的分布示意图。
图5为本发明实施例2的陆面洼地填充原理示意图。
图6为本发明实施例2的羌塘高原区的湖泊流域边界示意图。
图7为本发明实施例3的内陆湖泊流域提取系统的架构图。
具体实施方式
附图仅用于示例性说明，不能理解为对本专利的限制；
对于本领域技术人员来说，附图中某些公知结构及其说明可能省略是可以理解的。
下面结合附图和实施例对本发明的技术方案做进一步的说明。
实施例1
本实施例为解决流域管理中内陆盆区地形复杂、洼地遍布，用常规水文处理程序无法提取湖泊流域边界的问题，提出一种基于数字高程模型的内陆湖泊流域提取方法。
请参阅图1，本实施例提出的内陆湖泊流域提取方法的步骤如下：
S1、输入目标地区的数字高程模型和湖泊分布数据。
S2、计算目标地区的数字高程模型中各像元的初始栅格流向，并根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合。
S3、根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地，得到陆面洼地像元集合和湖泊洼地像元集合。
S4、对所述数字高程模型中相应的陆面洼地像元进行有选择地填充，然后计算栅格流向，得到栅格流向矩阵。
S5、对所述栅格流向矩阵进行解析并构建流向多叉树。
S6、对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线相应的边缘像元，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，提取得到湖泊流域边界。
本实施例中，首先计算初始栅格流向，并根据栅格流向从中提取无法确定流向的洼地像元集合，然后将湖泊分布数据转化为掩膜数组，识别分辨洼地集合中的陆面洼地与湖泊洼地，便于后续的洼地处理步骤。对数字高程模型中属于湖泊的洼地进行识别后，可有选择地对陆面洼地进行填充，对湖泊洼地予以保留，然后计算栅格流向数据，构建流向多叉树，用于后续的湖泊流域提取步骤。最后对湖泊数据集进行迭代遍历，提取湖泊矢量的外围边界线，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，提取得到湖泊流域边界。
在一可选实施例中，S2步骤中计算目标地区的数字高程模型中各像元的初始栅格流向的步骤包括：
S2.1、计算目标地区的数字高程模型中各像元与周围8个像元的高程差值。
S2.2、将各像元与相邻像元的高程差值除以像元之间的距离，计算得到像元间的斜率。
S2.3、根据各像元间的斜率，设定沿斜率最大的方向作为初始栅格流向。
其中，本实施例基于水沿着斜率最大的方向进行流动的假设得到栅格流向。
进一步地，根据各像元的初始栅格流向，对于数字高程模型中的洼地像元，其高程低于周围像元，8个方向的高程差值均为负值，水无法从该像元流往其他像元，其流向待定。
由此，根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合的步骤包括：
遍历所述数字高程模型中各个像元并进行判断：
若当前像元的高程值低于周围8个像元，或当前像元与周围8个像元的高程差值均为负数，则判断当前像元为无法确定流向的像元，将当前像元加入洼地像元集合中，并对当前像元填充预设的洼地像元填充值。
在一具体实施过程中，使用正整数1、2、3、4、5、6、7、8分别代表东、东北、北、西北、西、西南、南、东南8个流动方向，对于洼地像元，其流动方向待定，使用填充值-32768表示。
其中，该栅格流向数据用于辅助洼地识别，不参与后续的内陆湖泊流域提取。
在一可选实施例中，S3步骤中根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地的步骤包括：
S3.1、根据所述洼地像元集合提取相应的像元编号，组成洼地像元编号集合sinktotal；其表达式如下：
sinktotal＝{pixeli,j|flowdiri,j＝fillvalue}
pixeli,j＝i*rows+j
Figure BDA0003802769240000051
式中，pixeli,j表示第i行第j列的像元编号，flowdiri,j表示第i行第j列像元的流向，fillvalue为预设的洼地像元填充值；rows为数字高程模型的行数，cols为数字高程模型的列数。
本步骤中，将流向为填充值的像元作为洼地像元，得到洼地像元编号集合。
S3.2、根据所述湖泊分布数据，将第k个湖泊从矢量转化为与原数字高程模型形状相同、分辨率一致的栅格，其中被湖泊覆盖的像元为1，不被湖泊覆盖的像元为0，得到第k个湖泊的掩膜矩阵masklake,k。
S3.3、从所述洼地像元编号集合sinktotal中提取被湖泊所覆盖的洼地，遍历所有湖泊分布数据后，得到湖泊洼地像元集合和陆面洼地像元集合；其表达式如下：
Figure BDA0003802769240000052
式中，sinklake,k表示第k个湖泊所覆盖的洼地像元集合，其中，sinklake,k取0值表示该像元属于陆面洼地像元，sinklake,k表示取1值表示该像元属于湖泊洼地像元；sinkt表示洼地像元编号集合sinktotal中的第t个像元；
最后根据sinklake,k的取值划分得到陆面洼地像元集合和湖泊洼地像元集合。
本实施例中，通过提取洼地像元的编号集合，便于后续逐个洼地进行判断其属于湖泊洼地还是陆面洼地并进行分别处理。
此外，考虑到在洼地像元集合中，一部分洼地是由于数字高程模型摄影成像、数据滤波处理等过程产生的，并不代表真实地形；一部分洼地反映真实地形，需要在数字高程模型处理程序中予以保留，本实施例的S3.2步骤和S3.3步骤用于识别陆面洼地和湖泊洼地，利用将湖泊分布数据转化的湖泊掩膜矩阵覆盖在像元上，从洼地像元编号集合中提取该湖泊所覆盖的洼地像元编号，所覆盖的洼地像元编号即为湖泊像元编号。
而内陆湖泊的水位随季节发生变化，部分洼地在丰水年被湖泊覆盖，在平水年或枯水年将暴露在湖泊外成为陆面洼地。因此需要辨识湖泊洼地中的最低洼地，其余则视作陆面洼地，这样的处理不影响湖泊流域的提取。
由此，进一步地，本实施例的S3步骤中还包括以下步骤：
识别第k个湖泊所覆盖的洼地像元集合后，通过数字高程模型查询各像元的高程，并保留高程值最小的像元；其表达式如下：
Figure BDA0003802769240000061
sink1,sink2,...,sinkn∈sinklake,k
式中，sink'表示第k个湖泊所覆盖的最低洼地，
Figure BDA0003802769240000062
表示sinkt所对应的高程。
对所有湖泊重复上述操作后，得到所有湖泊最低洼地集合sinklowest，再通过差值运算得到陆面洼地集合sinkland＝sinktotal–sinklowest。
进一步地，在一可选实施例中，S4步骤的具体步骤包括：
S4.1、将所述湖泊洼地像元集合转化为与原数字高程模型形状相同、分辨率一致的掩膜数组masksink；掩膜数组masksink中，将属于湖泊洼地的像元标记为1，将不属于湖泊洼地的像元标记为0。
S4.2、根据所述掩膜数组masksink，对标记为0的洼地进行填平去除。
S4.3、对经过填平后的掩膜数组masksink进行流向计算，对数字高程模型中的得到栅格流向矩阵fdir。
本步骤通过对数字高程模型中属于湖泊的洼地进行识别后，有选择地对陆面洼地进行填充，对湖泊洼地予以保留，然后计算栅格流向数据，用于后续的湖泊流域提取步骤。
其中，首先将湖泊洼地像元编号的集合转化为与原数字高程模型形状相同、分辨率一致的掩膜数组masksink，该掩膜数组标记了需要保留的洼地，不被该数组标记的洼地将被填平去除。
在一具体实施过程中，采用TauDEM程序中的pitremove函数，将掩膜数组masksink传入depmask参数，对不被掩膜数组masksink标记的洼地进行填平去除。
TauDEM程序在进行流向计算时，默认将所有洼地进行填平处理，在消除伪洼地的影响时，也会将真实的洼地去掉，不能反映内陆湖泊流域真实的水文信息。而depmask参数可标记数字高程模型中的真实洼地，被标记的这部分洼地不会被填平，被标记洼地对应的流域像元流向仍可以指向该洼地。
由此，在一具体实施过程中，通过调用TauDEM程序中的flowdir流向计算程序，计算得到栅格流向矩阵fdir。
进一步地，在一可选实施例中，S5步骤中的具体步骤包括：
S5.1、对所述栅格流向矩阵fdir进行解析，构建上游像元与下游像元的对应关系；其表达式如下：
dir(pixeli,j)＝pixeli,j+offset(fdir(i,j))
pixeli,j＝i*rows+j
Figure BDA0003802769240000071
式中，pixeli,j表示第i行第j列的像元编号，rows为数字高程模型的行数，cols为数字高程模型的列数；fdir(i,j)表示第i行第j列的像元的栅格流向，offset(·)表示上游像元与下游像元的编号偏移量，dir(·)表示上游像元与下游像元的对应关系。
S5.2、对所述栅格流向矩阵fdir进行解析，构建下游像元与上游像元的关系，保存为字典；其表达式如下：
dirinverse(pixeli,j)＝{pixelm,n|dir(pixelm,n)＝pixeli,j}
Figure BDA0003802769240000072
Figure BDA0003802769240000073
式中，dirinverse(·)表示下游像元与上游像元的对应关系，pixelm,n表示满足下游像元是pixeli,j这一条件的上游像元。
S5.3、根据上游像元与下游像元的对应关系，以及下游像元与上游像元的关系，构建得到流向多叉树。
本实施例中，对于8种不同流向，上游像元与下游像元的编号偏移关系的表达式如下：
Figure BDA0003802769240000074
式中，offset(·)表示上游像元与下游像元的编号偏移量，y表示上游像元的栅格流向，其中1、2、3、4、5、6、7、8分别代表东、东北、北、西北、西、西南、南、东南8个流动方向。
一个上游像元仅对应一个下游像元，但一个下游像元与多个上游像元对应。本实施例中的流域湖泊提取的过程是根据栅格流向从下游向上游迭代搜索的过程需要对上游像元与下游像元的对应关系dir进行解析，构建下游像元与上游像元的关系。最后以字典的形式保存流向多叉树，用于后续湖泊流域提取。
进一步地，在一可选实施例中，S6步骤中的具体步骤包括：
S6.1、对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线，并将所述外围边界线由矢量数据转为栅格数据以提取湖泊对应的边缘像元。
S6.2、以湖泊对应的边缘像元作为起始像元，利用所述流向多叉树进行子树层遍历；其中，查询起始像元的上游像元后，再次将当前上游像元作为起始像元，往复迭代提取，得到湖泊流域边界。
在一具体实施过程中，对湖泊数据集进行迭代遍历，调用python第三方库shapely中的exterior方法提取湖泊矢量的外围边界线，利用第三方库rasterio中的geometry_mask方法将边界线由矢量数据转为栅格数据以提取湖泊对应的边缘像元，以湖泊边缘像元为起始像元，利用流向多叉树进行子树层序遍历，查询起始像元的上游像元后，再次将上游像元作为起始像元，往复迭代提取湖泊流域边界。
本实施例中，通过结合内陆湖泊等内流水体的地形特征，对目标地区的数字高程模型和湖泊分布数据进行处理和解析，根据各像元的栅格流向，对数字高程模型中像元是否为内陆湖泊进行判断，并对洼地像元进行有选择地填充，有效反映内陆湖泊流域真实的水文信息，进一步用于构建流向多叉树，将树状数据结构遍历算法应用于湖泊流域提取，实现高效率且高准确率的内陆湖泊等内流水体的流域边界提取。
实施例2
请参阅图2～图6，本实施例应用实施例1提出的内陆湖泊流域提取方法对羌塘高原区的内陆湖泊流域进行提取。
步骤1：本实施例将羌塘高原区湖泊流域提取所需要的数字高程模型、湖泊位置数据储存在相同的文件夹中，然后分别定义变量dem、lake、储存两个原始数据的文件路径，通过pygeoc中的d8flowdir函数计算初始流向。
步骤2：使用numpy中的where函数识别羌塘高原区的洼地，得到如图2所示的羌塘高原区数字高程模型与湖泊分布特征图。羌塘高原区的洼地分布情况如图3所示，其中，灰度颜色的深浅代表洼地的深度。可见在洼地当中一部分洼地被湖泊所覆盖，但仍有大部分洼地属于陆面洼地。
步骤3：对羌塘高原区的洼地进行区分，识别当中的湖泊洼地与陆面洼地。利用rasterio中的geometry_mask函数将湖泊矢量位置数据转化为栅格数据，与步骤2中识别得到的洼地进行叠加分析，对于各个湖泊所覆盖的洼地，仅保留当中高程最低的洼地作为湖泊洼地，其他洼地视作陆面洼地，得到如图4所示的湖泊洼地的分布图，其中点状图案代表湖泊所覆盖的最低洼地。
步骤4：利用TauDEM程序中的Pitremove功能对陆面洼地进行填充，湖泊洼地予以保留。该方法的示意图如图5所示，倒三角形代表洼地的最低点，箭头代表水流方向。图5中的左图(a)为处理前的洼地示意图，图5中的右图(b)为处理后的洼地示意图，其中陆面洼地被填充至与洼地边缘等高，湖泊洼地则予以保留。
步骤5：利用TauDEM程序中的flowdir功能对处理后的羌塘高原区数字高程模型进行栅格流向计算，然后根据不同流向的上游像元与下游像元编号的偏移关系解析流向数据，用python中的dict数据类型保存上游像元与下游像元之间的映射关系。
步骤6：对步骤5中的上游像元与下游像元之间的映射关系进行解析，利用python中的tuple数据类型保存各下游像元的上游像元，建立下游像元与上游像元的映射关系，构建流向多叉树。
步骤7：利用python中的for循环对湖泊数据集进行遍历，对于任意一个湖泊，利用shapely中的exterior方法提取湖泊边界，然后将湖泊边界转化为栅格像元，利用流向多叉树对这些边缘像元进行子树遍历，提取该湖泊的流域边界，得到shapely中的Polygon对象。使用geopandas中的to_file方法保存各湖泊流域。最终得到的羌塘高原区的湖泊流域边界示意图如图6所示。
实施例3
本实施例提出一种基于数字高程模型的内陆湖泊流域提取系统，应用实施例1提出的内陆湖泊流域提取方法。请参阅图7，为本实施例的内陆湖泊流域提取系统的架构图。
本实施例提出的内陆湖泊流域提取系统中，包括数据采集模块，洼地识别模块，洼地填充模块，流向多叉树构建模块，和湖泊流域边界提取模块。
其中：
数据采集模块，用于获取目标地区的数字高程模型和湖泊分布数据。
洼地识别模块，用于通过计算目标地区的数字高程模型中各像元的初始栅格流向，并根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合；以及用于根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地，输出陆面洼地像元集合和湖泊洼地像元集合。
洼地填充模块，用于对所述数字高程模型中相应的陆面洼地像元进行有选择地填充，然后计算栅格流向，输出栅格流向矩阵。
流向多叉树构建模块，用于对所述栅格流向矩阵进行解析并构建流向多叉树。
湖泊流域边界提取模块，用于对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线相应的边缘像元，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，输出提取得到的湖泊流域边界数据。
在一具体实施过程中，所述数据采集模块将获取的目标地区的数字高程模型和湖泊分布数据传输至洼地识别模块中进行洼地识别处理。
洼地识别模块根据输入的数字高程模型数据，计算目标地区的数字高程模型中各像元与周围8个像元的高程差值，将各像元与相邻像元的高程差值除以像元之间的距离，计算得到像元间的斜率，并根据各像元间的斜率，设定沿斜率最大的方向作为初始栅格流向。然后根据所述初始栅格流向，遍历所述数字高程模型中各个像元并进行判断：若当前像元的高程值低于周围8个像元，或当前像元与周围8个像元的高程差值均为负数，则判断当前像元为无法确定流向的像元，将当前像元加入洼地像元集合中，并对当前像元填充预设的洼地像元填充值。
进一步地，洼地识别模块根据所述洼地像元集合提取相应的像元编号，组成洼地像元编号集合，并根据所述湖泊分布数据，将第k个湖泊从矢量转化为与原数字高程模型形状相同、分辨率一致的栅格，其中被湖泊覆盖的像元为1，不被湖泊覆盖的像元为0，得到第k个湖泊的掩膜矩阵masklake,k；最后从所述洼地像元编号集合sinktotal中提取被湖泊所覆盖的洼地，遍历所有湖泊分布数据后，得到湖泊洼地像元集合和陆面洼地像元集合。洼地识别模块将湖泊洼地像元集合和陆面洼地像元集合发送至洼地填充模块中。
洼地填充模块将所述湖泊洼地像元集合转化为与原数字高程模型形状相同、分辨率一致的掩膜数组masksink；掩膜数组masksink中，将属于湖泊洼地的像元标记为1，将不属于湖泊洼地的像元标记为0，然后根据所述掩膜数组masksink，对标记为0的洼地进行填平去除，并对经过填平后的掩膜数组masksink进行流向计算，对数字高程模型中的得到栅格流向矩阵fdir。洼地填充模块将其计算生成的栅格流向矩阵fdir发送至流向多叉树构建模块。
流向多叉树构建模块对输入的栅格流向矩阵fdir进行解析，构建上游像元与下游像元的对应关系，以及对所述栅格流向矩阵fdir进行解析，构建下游像元与上游像元的关系，最后根据上游像元与下游像元的对应关系，以及下游像元与上游像元的关系，构建得到流向多叉树。流向多叉树构建模块将其生成的流向多叉树发送至湖泊流域边界提取模块中。
湖泊流域边界提取模块对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线，并将所述外围边界线由矢量数据转为栅格数据以提取湖泊对应的边缘像元，然后以湖泊对应的边缘像元作为起始像元，利用所述流向多叉树进行子树层遍历；其中，查询起始像元的上游像元后，再次将当前上游像元作为起始像元，往复迭代提取，得到湖泊流域边界并输出。
进一步地，本实施例还提出一种基于数字高程模型的内陆湖泊流域提取系统，包括存储器和处理器，所述存储器存储有计算机程序，所述处理器执行所述计算机程序时实现实施例1提出的内陆湖泊流域提取方法的步骤。
具体地，通过定义类class以及定义函数def()的方式，将S1～S6步骤中的洼地地形处理与湖泊流域提取过程封装成多个类及函数，将代码保存为.py格式文件。在需要调用时，使用import函数导入类及函数后，即可使用。
附图中描述位置关系的用语仅用于示例性说明，不能理解为对本专利的限制；
显然，本发明的上述实施例仅仅是为清楚地说明本发明所作的举例，而并非是对本发明的实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或变动。这里无需也无法对所有的实施方式予以穷举。凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明权利要求的保护范围之内。

Claims (10)

1.一种基于数字高程模型的内陆湖泊流域提取方法，其特征在于，包括以下步骤：
S1、输入目标地区的数字高程模型和湖泊分布数据；
S2、计算目标地区的数字高程模型中各像元的初始栅格流向，并根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合；
S3、根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地，得到陆面洼地像元集合和湖泊洼地像元集合；
S4、对所述数字高程模型中相应的陆面洼地像元进行有选择地填充，然后计算栅格流向，得到栅格流向矩阵；
S5、对所述栅格流向矩阵进行解析并构建流向多叉树；
S6、对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线相应的边缘像元，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，提取得到湖泊流域边界。
2.根据权利要求1所述的内陆湖泊流域提取方法，其特征在于，所述S2步骤中，计算目标地区的数字高程模型中各像元的初始栅格流向的步骤包括：
S2.1、计算目标地区的数字高程模型中各像元与周围8个像元的高程差值；
S2.2、将各像元与相邻像元的高程差值除以像元之间的距离，计算得到像元间的斜率；
S2.3、根据各像元间的斜率，设定沿斜率最大的方向作为初始栅格流向。
3.根据权利要求2所述的内陆湖泊流域提取方法，其特征在于，所述S2步骤中，根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合的步骤包括：
遍历所述数字高程模型中各个像元并进行判断：若当前像元的高程值低于周围8个像元，或当前像元与周围8个像元的高程差值均为负数，则判断当前像元为无法确定流向的像元，将当前像元加入洼地像元集合中，并对当前像元填充预设的洼地像元填充值。
4.根据权利要求3所述的内陆湖泊流域提取方法，其特征在于，所述S3步骤中，根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地的步骤包括：
S3.1、根据所述洼地像元集合提取相应的像元编号，组成洼地像元编号集合sinktotal；其表达式如下：
sinktotal＝{pixeli,j|flowdiri,j＝fillvalue}
Figure FDA0003802769230000021
式中，pixeli,j表示第i行第j列的像元编号，flowdiri,j表示第i行第j列像元的流向，fillvalue为预设的洼地像元填充值；rows为数字高程模型的行数，cols为数字高程模型的列数；
S3.2、根据所述湖泊分布数据，将第k个湖泊从矢量转化为与原数字高程模型形状相同、分辨率一致的栅格，其中被湖泊覆盖的像元为1，不被湖泊覆盖的像元为0，得到第k个湖泊的掩膜矩阵masklake,k；
S3.3、从所述洼地像元编号集合sinktotal中提取被湖泊所覆盖的洼地，遍历所有湖泊分布数据后，得到湖泊洼地像元集合和陆面洼地像元集合；其表达式如下：
Figure FDA0003802769230000022
式中，sinklake,k表示第k个湖泊所覆盖的洼地像元集合，其中，sinklake,k取0值表示该像元属于陆面洼地像元，sinklake,k表示取1值表示该像元属于湖泊洼地像元；sinkt表示洼地像元编号集合sinktotal中的第t个像元；
根据sinklake,k的取值划分得到陆面洼地像元集合和湖泊洼地像元集合。
5.根据权利要求4所述的内陆湖泊流域提取方法，其特征在于，所述S3步骤中，还包括以下步骤：
识别第k个湖泊所覆盖的洼地像元集合后，通过数字高程模型查询各像元的高程，并保留高程值最小的像元；其表达式如下：
Figure FDA0003802769230000023
sink1,sink2,...,sinkn∈sinklake,k
式中，sink'表示第k个湖泊所覆盖的最低洼地，
Figure FDA0003802769230000024
表示sinkt所对应的高程；
对所有湖泊重复上述操作后，得到所有湖泊最低洼地集合sinklowest，再通过差值运算得到陆面洼地集合sinkland＝sinktotal–sinklowest。
6.根据权利要求1～5任一项所述的内陆湖泊流域提取方法，其特征在于，所述S4步骤中，其具体步骤包括：
S4.1、将所述湖泊洼地像元集合转化为与原数字高程模型形状相同、分辨率一致的掩膜数组masksink；掩膜数组masksink中，将属于湖泊洼地的像元标记为1，将不属于湖泊洼地的像元标记为0；
S4.2、根据所述掩膜数组masksink，对标记为0的洼地进行填平去除；
S4.3、对经过填平后的掩膜数组masksink进行流向计算，对数字高程模型中的得到栅格流向矩阵fdir。
7.根据权利要求6所述的内陆湖泊流域提取方法，其特征在于，所述S5步骤中，其具体步骤包括：
S5.1、对所述栅格流向矩阵fdir进行解析，构建上游像元与下游像元的对应关系；其表达式如下：
dir(pixeli,j)＝pixeli,j+offset(fdir(i,j))
Figure FDA0003802769230000031
式中，pixeli,j表示第i行第j列的像元编号，rows为数字高程模型的行数，cols为数字高程模型的列数；fdir(i,j)表示第i行第j列的像元的栅格流向，offset(·)表示上游像元与下游像元的编号偏移量，dir(·)表示上游像元与下游像元的对应关系；
S5.2、对所述栅格流向矩阵fdir进行解析，构建下游像元与上游像元的关系；其表达式如下：
Figure FDA0003802769230000032
式中，dirinverse(·)表示下游像元与上游像元的对应关系，pixelm,n表示满足下游像元是pixeli,j这一条件的上游像元；
S5.3、根据上游像元与下游像元的对应关系，以及下游像元与上游像元的关系，构建得到流向多叉树。
8.根据权利要求7所述的内陆湖泊流域提取方法，其特征在于，所述S6步骤中，其具体步骤包括：
S6.1、对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线，并将所述外围边界线由矢量数据转为栅格数据以提取湖泊对应的边缘像元；
S6.2、以湖泊对应的边缘像元作为起始像元，利用所述流向多叉树进行子树层遍历；其中，查询起始像元的上游像元后，再次将当前上游像元作为起始像元，往复迭代提取，得到湖泊流域边界。
9.一种基于数字高程模型的内陆湖泊流域提取系统，应用于权利要求1～8任一项所述的内陆湖泊流域提取方法，其特征在于，包括：
数据采集模块，用于获取目标地区的数字高程模型和湖泊分布数据；
洼地识别模块，用于通过计算目标地区的数字高程模型中各像元的初始栅格流向，并根据所述初始栅格流向选取无法确定流向的像元组成洼地像元集合；以及用于根据所述湖泊分布数据识别所述洼地像元组合中的陆面洼地和湖泊洼地，输出陆面洼地像元集合和湖泊洼地像元集合；
洼地填充模块，用于对所述数字高程模型中相应的陆面洼地像元进行有选择地填充，然后计算栅格流向，输出栅格流向矩阵；
流向多叉树构建模块，用于对所述栅格流向矩阵进行解析并构建流向多叉树；
湖泊流域边界提取模块，用于对所述湖泊分布数据进行迭代遍历，提取湖泊矢量的外围边界线相应的边缘像元，然后以湖泊边缘像元为起始像元，利用所述流向多叉树进行子树层遍历，输出提取得到的湖泊流域边界数据。
10.一种基于数字高程模型的内陆湖泊流域提取系统，包括存储器和处理器，所述存储器存储有计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至8中任一项所述的内陆湖泊流域提取方法的步骤。
CN202210992072.3A 2022-08-17 2022-08-17 一种基于数字高程模型的内陆湖泊流域提取方法及系统 Pending CN115359221A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210992072.3A CN115359221A (zh) 2022-08-17 2022-08-17 一种基于数字高程模型的内陆湖泊流域提取方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210992072.3A CN115359221A (zh) 2022-08-17 2022-08-17 一种基于数字高程模型的内陆湖泊流域提取方法及系统

Publications (1)

Publication Number Publication Date
CN115359221A true CN115359221A (zh) 2022-11-18

Family

ID=84003521

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210992072.3A Pending CN115359221A (zh) 2022-08-17 2022-08-17 一种基于数字高程模型的内陆湖泊流域提取方法及系统

Country Status (1)

Country Link
CN (1) CN115359221A (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116415318A (zh) * 2023-03-20 2023-07-11 长江水利委员会长江科学院 一种基于数学形态学的内流区湖泊水文连通性建模方法

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116415318A (zh) * 2023-03-20 2023-07-11 长江水利委员会长江科学院 一种基于数学形态学的内流区湖泊水文连通性建模方法
CN116415318B (zh) * 2023-03-20 2024-03-08 长江水利委员会长江科学院 一种基于数学形态学的内流区湖泊水文连通性建模方法

Similar Documents

Publication Publication Date Title
CN105551028B (zh) 一种基于遥感影像的地理空间数据动态更新的方法及系统
Fairfield et al. Drainage networks from grid digital elevation models
Mackay et al. Extraction and representation of nested catchment areas from digital elevation models in lake‐dominated topography
CN111028255A (zh) 基于先验信息与深度学习的农田区域预筛选方法及装置
CN110415265A (zh) 基于无人机高精度dem数据的梯田自动提取方法
CN114821355B (zh) 海岸线自动化识别方法和装置
CN104715474B (zh) 基于标记分水岭算法的高分辨率合成孔径雷达图像线性建筑物检测方法
CN111881816A (zh) 一种长时序的河湖围埂养殖区域监测方法
CN109614513A (zh) 一种面向区域覆盖的遥感影像数据选取方法和系统
McCormack et al. Feature-based derivation of drainage networks
CN115359221A (zh) 一种基于数字高程模型的内陆湖泊流域提取方法及系统
CN112037244A (zh) 联合指数和轮廓指示器SLIC的Landsat-8图像养殖池提取方法
CN104182976B (zh) 一种野外运动目标精细化提取方法
Lyu et al. A vector-based method for drainage network analysis based on LiDAR data
CN116415318B (zh) 一种基于数学形态学的内流区湖泊水文连通性建模方法
CN112861669B (zh) 地表坡向约束的高分辨率dem地形特征增强提取方法
CN110569733A (zh) 基于遥感大数据平台的湖泊长时序连续水域变化重建方法
CN112966657A (zh) 一种大尺度水体覆盖的遥感自动分类方法
CN109064430B (zh) 一种针对航拍区域含云图的除云方法和系统
CN115223000B (zh) 面向耕地资源监测的深度学习样本制作方法
WO2024020744A1 (zh) 基于高分辨率影像的生态图斑提取方法
Ning et al. Precipitation scale effect of the TRMM satellite in Tianshan, China
CN114494530A (zh) 一种证件照抠图方法及系统
CN111400932A (zh) 水环境仿真中网格的生成方法、装置、存储介质及终端
Sowmya et al. Watershed delineation using SRTM DEM for Bangalore South Region

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