CN110457771B - 一种基于高程偏差传递的dem水流方向计算方法 - Google Patents

一种基于高程偏差传递的dem水流方向计算方法 Download PDF

Info

Publication number
CN110457771B
CN110457771B CN201910655040.2A CN201910655040A CN110457771B CN 110457771 B CN110457771 B CN 110457771B CN 201910655040 A CN201910655040 A CN 201910655040A CN 110457771 B CN110457771 B CN 110457771B
Authority
CN
China
Prior art keywords
unit
elevation
dem
water flow
array
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.)
Active
Application number
CN201910655040.2A
Other languages
English (en)
Other versions
CN110457771A (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.)
Hohai University HHU
Original Assignee
Hohai University HHU
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 Hohai University HHU filed Critical Hohai University HHU
Priority to CN201910655040.2A priority Critical patent/CN110457771B/zh
Publication of CN110457771A publication Critical patent/CN110457771A/zh
Application granted granted Critical
Publication of CN110457771B publication Critical patent/CN110457771B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • 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

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Software Systems (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Computational Mathematics (AREA)
  • Geometry (AREA)
  • Operations Research (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Computer Graphics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明提供了一种基于高程偏差传递的DEM水流方向计算方法，属于数字地形分析技术领域。其技术方案为：一种基于高程偏差传递的DEM水流方向计算方法，采用水流方向组成的水流路径接近水流在地表的真实移动路径；对数值模拟水和污染物输移过程的水文物理模型使用了DEM获取更高精度水流路径的方法。本发明的有益效果为：本发明通过使用水流到达下游单元的真实高程和下游单元中心高程的高度偏差作为传递参数校正水流方向的方法；相较于D8方法，本发明提供了水流方向组成的水流路径更接近水流在地表的真实移动路径。

Description

一种基于高程偏差传递的DEM水流方向计算方法
技术领域
本发明涉及数字地形分析技术领域，具体涉及一种通过使用水流到达下游单元的真实高程和下游单元中心高程的高度偏差作为传递参数校正水流方向的方法。
背景技术
水流方向是一个重要的水文地貌参数，广泛应用于许多对水流和污染物进行输移过程模拟的数值模型。目前的水流方向计算主要基于数字高程模型(Digital ElevationModel,DEM)进行。DEM是对地形的一种概化描述，将地形简化为大量正方向栅格单元，每个单元赋予一个高程值视为该单元平均高度。在目前的水流方向计算的研究中，主要将每个单元的高程作为该单元中心点的高度。
目前最经典的水流方向计算方法是O’Callaghan和Mark(1984)提出的D8方法，该方法将DEM栅格单元的水流方向指向相邻单元中与中心单元坡度最大者，但是由于这种方法只考虑了局部坡度，而且水流的真实方向可能并非8个允许方向之一，因此从整体上看D8方法提供的水流方向连接成的完整水流路径存在较大偏差。尽管Quinn等(1991)提出了允许一个单元具备多个水流方向的方法，但该方法计算更为复杂且不适用于流域划分、流动距离计算等领域。目前广泛接受的是允许水流方向从0°到360°之间任意方向的无穷流向方法(Tarboton，1997)，但是由于每个单元的流向都起始于单元中心，如果水流方向不是8方向之一，当水流沿着该方向流动后就没有了进一步的移动方向，因此该方法在实现上仍是一种多流向算法，同样存在前面所说的问题。
因此，为了更好地服务于流域划分、流动距离计算，需要有一种改进的八方向方法。考虑到水流实际应该是沿0°到360°之间任意方向流动，其实际到达位置的高度也并非DEM栅格单元中心的高程，因此使用这部分高度偏差对下游的水流方向进行修正是一个很好的想法。
发明内容
本发明要解决的技术问题在于，解决现有的应用于DEM地形的单流向方法中只考虑局部坡度，忽略的地形整体变化，水流路径偏差较大的问题，而提供一种基于高程偏差传递的DEM水流方向计算方法。
本发明是通过如下措施实现的：一种基于高程偏差传递的DEM水流方向计算方法，主要包括以下步骤：
步骤S1：加载DEM数据，DEM为二维数组格式，构建3个与DEM数据尺寸相同的数组A、B、C，数组A用于保存最低相邻单元方向、数组B用于保存第二可能流动方向、数组C用于保存最终确定的流向，并且3个数组中所有位置的初始值都设为0，再构建一个能保存栅格单元三维坐标信息，并根据栅格高程从低到高依次排列的优先队列Q，以及另一个从高到低排列栅格单元三维坐标信息的优先队列T，两个优先队列中同一高程的不同单元按插入先后顺序先后排列，新插入的单元排在同一高程单元的最后；
步骤S2：扫描DEM，将DEM中的有效地形单元加入优先队列T排序，同时将位于有效地形边缘的单元加入优先队列Q进行排序，并同时在数组A、B、C中赋予相同的流向，该流向指向DEM有效地形外侧；
步骤S3：不断取出优先队列Q头部的单元，检索该头部单元的8个相邻单元中的未在数组A中赋予流向的有效单元，为其在数组A中赋予指向此头部单元的方向并加入队列Q，若该未赋流向值的单元高程值不高于P，还需为其在数组B与C中赋予指向此头部单元的方向；如果该单元在数组B中对应方向值为初始值0，则比较该单元邻近其最低相邻单元的另外两个相邻单元的高程，将其中较低者所在方向作为第二可能流向并存入数组B，重复此步骤至队列Q为空队列；
步骤S4：读取优先队列T的头部单元K，若单元K在数组C中对应的流向为初始值0，则从单元K开始在最低下坡单元和第二可能方向间使用传递的高程偏差判断从K开始到某个已经具有最终流向的单元的流动路径上各单元的水流方向，然后将K从队列T中移除，重新读取队列T的头部单元以重复此步骤，直至队列T为空队列。
作为本发明提供的一种基于高程偏差传递的DEM水流方向计算方法进一步优化方案，所述步骤S4中，从单元K开始在最低下坡单元方向和第二可能方向间使用传递的高程偏差判断从K开始到某个已经具有最终流向的单元的流动路径上各单元的水流方向，其具体包括如下步骤：
步骤S4-1：初始化高程偏差传递的参数，以单元K为当前计算单元，以K的高程为h0，初始高程偏差d为0，K的最低下坡单元方向和第二可能方向指向的两个单元相邻，二者一个在K的正方向，一个在K的斜方向，设其中正方向单元高程为h1，斜方向单元高程为h2；
步骤S4-2：构建h0、h1、h2这3个点组成的三角形平面，确定如果水流从点h0出发到达平面底部边界时目标点的高度h3；如果三角形平面的坡度方向超出三角面边缘，水流沿与坡度方向最邻近的三角面边缘方向，否则就沿着坡度方向直线移动，移动线路跟h1、h2连接线的交点所在高程即为h3，具体计算方程如下：
Figure GDA0003753481900000031
步骤S4-3：由于单流向方法是从上游单元中心指向下游单元中心，上游单元水流到达下游单元时到达点的高度并不一定是下游单元中心点的高程值(即该单元对应的DEM栅格值)，水流从自上游到达本单元时的真实高度出发，下降高度与自本单元中心流至下游单元时下降高度h3相同，如此到达的高度作为最终下降高度，比较最低相邻单元与第二可能流动方向对应的单元高程，选取其中高程与最终下降高度最接近的单元作为下游单元，由于在地形走向发生旋转时的水流方向过于复杂，只在上下游单元的最低相邻单元方向和第二可能方向都相同时如此考虑，如果不同则认为水流到达最低相邻单元，到达高程为最低相邻单元中心点高程；
步骤S4-4：若最终流向指向的单元L在数组C中对应的流向不是初始值0，则步骤结束，否则以该被指向的单元L为当前计算单元，以L在数组A、B中对应方向指向的两个单元中，位于正向的单元高程作为新的h1值，斜向的单元高程作为新的h2值，以下文所述的d*的值作为d的值，返回步骤4-2重新进行。
作为本发明提供的一种基于高程偏差传递的DEM水流方向计算方法进一步优化方案，所述步骤S2具体内容为：从左到右、从上到下逐个检索DEM栅格单元，如果某一单元的高程不是无效值，则将其加入优先队列T，并根据两个原则判断该有效值单元是否为边缘单元；
步骤S2-1：若该单元位于整个矩形DEM数组的第一行、最后一行、第一列、最后一列，则将该单元加入优先队列Q，并将其在数组A、B和C中对应的单元的流向指向DEM外；
步骤S2-2：若该单元位于DEM范围的上下左右边缘志之外，就按中从P1到P8的顺序检索当前单元P0周围8个相邻单元，若存在某个相邻单元高程为无效地形值，则将该单元加入优先队列Q，并将该单元的流向指向无效值相邻单元，同时记入数组A、B和C对应位置，如果存在多个相邻单元为无效值，则最终流向指向第一个检索到的单元。
作为本发明提供的一种基于高程偏差传递的DEM水流方向计算方法进一步优化方案，所述步骤S3具体内容为：读取优先队列Q头部的栅格单元P，检索其周边的8个相邻栅格，比较与其在数组A中的方向相邻的两个方向对应的相邻单元的高程，将其中高程较低者所在方向作为第二可能水流方向存入P在数组B中对应的单元，并将所有高程值不是无效值且在数组A中对应的流向是初始值0的相邻单元加入队列Q，同时将他们在数组A中对应的流向指向P，随后将P从优先队列Q中移除，之后重新读取Q头部的单元并重复此步骤直至队列Q中不再存有单元，此时DEM中的有效地形单元都在数组A和B中被赋予了指向最低相邻单元的流向和第二可能流向；在DEM中存在部分洼地或者平原，这类地形指的是部分DEM栅格的8个相邻单元的高程都不低于中心单元，赋予这类地形向外流的水流路径，对这类地形不做基于偏差传递的处理，直接以此处的方向作为最终方向，即当P的邻近单元对应的流向是初始值0且高程值不高于P时，同时将他们在数组A、B、C中对应的流向都指向P。
作为本发明提供的一种基于高程偏差传递的DEM水流方向计算方法进一步优化方案，所述步骤S4-3中，最终下降高度与下游单元中心高程间的差值d进行数学表示，比较h3加上传递高程偏差d的综合高度与h1和h2高度中点的比值，若
Figure GDA0003753481900000041
则最终流向指向h1并计入数组C，新的传递偏差为
Figure GDA0003753481900000042
若
Figure GDA0003753481900000043
则最终流向指向h2并计入数组C，新的传递偏差为
Figure GDA0003753481900000044
所述DEM流向即为数字高程模型水流方向。
本发明的有益效果为：本发明通过使用水流到达下游单元的真实高程和下游单元中心高程的高度偏差作为传递参数校正水流方向的方法；相较于D8方法，本发明提供了水流方向组成的水流路径更接近水流在地表的真实移动路径；对数值模拟水和污染物输移过程的水文物理模型使用了DEM获取更高精度水流路径的方法，此方法同样以8个方向作为输出，具备简单的数据处理结构；具体是使用水流到达下游单元的真实高程与下游单元中心高程间的高度偏差作为传递参数校正水流方向，本方法提供的水流方向连接成的完整水流路径更为合理，本发明为水文地貌模型，尤其是只允许单个DEM栅格单元具备唯一流向的模型提供了获取更高精度水流方向的方法。
附图说明
图1为本发明的整体流程图。
图2为本发明的中心单元及其周围8个单元位置示意图。
图3为本发明的8个方向数字指代示意图，其中1代表正北，其他方向顺时针旋转。
图4为本发明从中心单元中点出发的三角形平面上的水流示意图。
图5为本发明实施例1使用的DEM示意图。
图6为本发明实施例1处理步骤S3的示意图A、B、C数组结果示意图。
图7为本发明实施例1处理步骤S4的第一条流动路径最终流向和处理完步骤S4的最终结果示意图。
图8(a)、图8(b)为本发明的倾斜平面和锥形面三维图。
图9(a)、图9(b)、图9(c)和图9(d)为可视化后的等高线和D8方法及本发明得到的最终水流路径示意图。
具体实施方式
为能清楚说明本方案的技术特点，下面通过具体实施方式，对本方案进行阐述。
参见图1至图8，本发明是：一种基于高程偏差传递的DEM水流方向计算方法，主要包括以下步骤：
步骤S1：加载DEM数据，DEM表现为二维数组格式，构建3个与DEM数据尺寸相同的数组A、B、C，数组A用于保存最低相邻单元方向、数组B用于保存第二可能流动方向、数组C用于保存最终确定的流向，并且3个数组中所有位置的初始值都设为0，再构建一个能保存栅格单元(径度X，纬度Y，高程Z)三维坐标信息，并根据栅格高程从低到高依次排列的优先队列Q，以及另一个从高到低排列栅格单元三维坐标信息的优先队列T，两个优先队列中同一高程的不同单元都按插入先后顺序先后排列，新插入的单元排在同一高程单元的最后；
步骤S2：扫描DEM，将DEM中的有效地形单元加入优先队列T排序，同时将位于有效地形边缘的单元加入优先队列Q进行排序，并同时在数组A、B、C中赋予相同的流向，该流向指向DEM有效地形外侧，具体为：
从左到右、从上到下逐个检索DEM栅格单元，如果某一单元的高程不是无效值，则将其加入优先队列T，并根据两个原则判断该有效值单元是否为边缘单元：
(1)若该单元位于整个矩形DEM范围的上下左右边缘(即DEM数组的第一行、最后一行、第一列、最后一列)，则将该单元加入优先队列Q，并将其在数组A、B和C中对应的单元的流向指向DEM外；
(2)若该单元不位于DEM范围的上下左右边缘，就按图2中从P1到P8的顺序检索当前单元P0周围8个相邻单元，若存在某个相邻单元高程为无效地形值(一般DEM中设为-9999)，则将该单元加入优先队列Q，并将该单元的流向指向无效值相邻单元，同时记入数组A、B和C对应位置，如果存在多个相邻单元为无效值，则最终流向指向第一个检索到的单元；
上述过程中所指的流向指8个相邻单元所在的东、南、西、北4个正方向和东南、西南、西北、东北4个斜方向，这里以图3所示的数字1-8指代8个方向记录入数组；
步骤S3、不断取出优先队列Q头部的单元，检索该头部单元的8个相邻单元中的未在数组A中赋予流向的有效单元，为其在数组A中赋予指向此头部单元的方向并加入队列Q，若该未赋流向值的单元高程值不高于P，还需为其在数组B与C中赋予指向此头部单元的方向；此外，如果该单元在数组B中对应方向值为初始值0，则比较该单元邻近其最低相邻单元的另外两个相邻单元的高程，将其中较低者所在方向作为第二可能流向并存入数组B，重复此步骤至队列Q为空队列，具体为：
读取优先队列Q头部的栅格单元P，检索其周边的8个相邻栅格，比较与其在数组A中的方向相邻的两个方向对应的相邻单元的高程，将其中高程较低者所在方向作为第二可能水流方向存入P在数组B中对应的单元，并将所有高程值不是无效值且在数组A中对应的流向是初始值0的相邻单元加入队列Q，同时将他们在数组A中对应的流向指向P，随后将P从优先队列Q中移除，之后重新读取Q头部的单元并重复此步骤直至队列Q中不再存有单元，此时DEM中的有效地形单元都在数组A和B中被赋予了指向最低相邻单元的流向和第二可能流向；
值得注意的是，在DEM中存在部分洼地或者平原的现象，这类地形指的是部分DEM栅格的8个相邻单元的高程都不低于中心单元，这类地形无法简单得到向外流的水流方向，因此需要强行赋予其向外流的水流路径，因此对这类地形不做基于偏差传递的处理，直接以此处的方向作为最终方向，即当P的邻近单元对应的流向是初始值0且高程值不高于P时，同时将他们在数组A、B、C中对应的流向都指向P；
步骤S4，读取优先队列T的头部单元K，若单元K在数组C中对应的流向为初始值0，则从单元K开始在最低下坡单元和第二可能方向间使用传递的高程偏差判断从K开始到某个已经具有最终流向的单元的流动路径上各单元的水流方向，然后将K从队列T中移除，重新读取队列T的头部单元以重复此步骤，直至队列T为空队列。
进一步地，所述从单元K开始在最低下坡单元方向和第二可能方向间使用传递的高程偏差判断从K开始到某个已经具有最终流向的单元的流动路径上各单元的水流方向的步骤包括：
步骤S4-1：初始化高程偏差传递的参数，以单元K为当前计算单元，以K的高程为h0，初始高程偏差d为0，易知K的最低下坡单元方向和第二可能方向指向的两个单元相邻，也就是二者一个在K的正方向，一个在K的斜方向，设其中正方向单元高程为h1，斜方向单元高程为h2；
步骤S4-2：如图4，构建h0、h1、h2这3个点组成的三角形平面，确定如果水流从点h0出发到达平面底部边界时目标点的高度h3。这里如果三角形平面的坡度方向超出三角面边缘则认为水流沿与坡度方向最邻近的三角面边缘方向，否则就沿着坡度方向直线移动，移动线路跟h1、h2连接线的交点所在高程即为h3，具体计算方程如下：
Figure GDA0003753481900000071
步骤S4-3：由于单流向方法只能从上游单元中心指向下游单元中心，但上游单元水流到达下游单元时到达点的高度并不一定是下游单元中心点的高程值也就是该单元对应的DEM栅格值，这里认为水流从自上游到达本单元时的高度出发，下降高度与自本单元中心流至下游单元时下降高度相同，如此到达的高度作为最终下降高度，比较最低相邻单元与第二可能流动方向对应的单元高程，选取其中高程与最终下降高度最接近的单元作为下游单元，由于在地形走向发生旋转时的水流方向过于复杂，这里只在上下游单元的最低相邻单元方向和第二可能方向都相同时如此考虑，如果不同则认为水流到达最低相邻单元，到达高程为最低相邻单元中心点高程，这里使用最终下降高度与下游单元中心高程间的差值d进行数学表示，具体为：
比较h3加上传递高程偏差d的综合高度与h1和h2高度中点的比值，若
Figure GDA0003753481900000072
则最终流向指向h1并计入数组C，新的传递偏差为
Figure GDA0003753481900000073
若
Figure GDA0003753481900000074
则最终流向指向h2并计入数组C，新的传递偏差为
Figure GDA0003753481900000075
步骤S4-4：若最终流向指向的单元L在数组C中对应的流向不是初始值0，则步骤结束，否则以该被指向的单元L为当前计算单元，以L在数组A、B中对应方向指向的两个单元中，位于正向的单元高程作为新的h1值，斜向的单元高程作为新的h2值，以d*的值作为d的值，返回步骤S4-2重新进行。
为了更好地实现本发明的目的，发明进行水流方向计算的过程以及发明的优越性，下面以实例进行说明。
使用图5所示的DEM，其中A5为无效值单元，这里用null表示无效值，根据步骤S1创建数组A、B、C和优先队列Q和T。根据步骤S2，将位于DEM边缘和有效地形边缘的非无效值单元加入队列Q按高程排序，同时赋予指向无效值单元或者区域外的初始流向值，并同时存入数组A、B、C，如图6中(a)所示，图中在此步骤中已处理单元被以箭头覆盖，箭头方向是数组A的数值决定的流向。
位于队列Q最前端的是单元D1，检索D1的相邻单元，其中C2和D2在数组A中对应值为0，于是为其赋予指向D1的流向值存入数组A。由于D1在数组B中对应值不为0，故不再判断第二方向。将单元D1从队列Q中移除，并将C2和D2插入队列Q中(图6中(b))。此时位于队列Q最前端的单元为D5，于是开始检索D5的相邻单元，重复步骤3，直至队列Q里不再包含任何单元。此后的数组A、B、C如图6中(d)、6中(e)、6中(f)所示，图6中(c)为数组A的方向表示。
队列T是根据栅格高程进行从高到低的排序，这里未做图示，总之，步骤S4从最高的单元C3开始，数组A中的最低邻近单元方向为6，B中的第二方向为5，即最低邻近单元方向指向D2，第二方向指向D3，故h0为7，h1为4，h2为2，且d＝0；可知
Figure GDA0003753481900000081
所以最终流向指向h1即D3方向，在数组C的C3处记录该方向值5，根据
Figure GDA0003753481900000082
传递偏差值为1.333，由于下游单元D3的最低邻近单元方向和第二方向值与C3相同，故D3从C3继承传递偏差，D3的d值依然为1.333。随后计算D3，同样对应A中最低邻近单元方向为6，B中的第二方向为5，即最低邻近单元方向指向E2，第二方向指向E3，故h0为4，h1为3，h2为2，且d＝1.333，于是
Figure GDA0003753481900000083
故D3的最终方向指向h2所在E2单元，在数组C的C3处记录该方向值6。由于E2在数组C中已有最终流向，故从C3单元开始的独立水流路径已结束，此时的最终流向数组C见图7中(a)。从队列T中移除C3，再次检索队列T头部单元，即DEM中第二高的单元B2，从步骤S4-2开始重复循环，如此往复直至T中不再包含任何单元。此时的数组C(图7中(b))即为使用本发明得到的最终水流方向成果。
以图8(a)中的倾斜平面和图8(b)中的锥形面和为例，倾斜平面的理想流动路径应该是垂直于等高线的直线，锥形面的理想流动路径应该是发源于圆锥顶点的直线并垂直于等高线，使用1m分辨率的DEM，图9(a)是使用传统D8方法得到的流向结果可视化出的倾斜平面流动路径，图9(b)是使用本发明得到的流向结果可视化出的倾斜平面流动路径，图9(c)是使用传统D8方法得到的流向结果可视化出的锥形面流动路径，图9(d)是使用本发明得到的流向结果可视化出的锥形面流动路径。可见本发明得到的流动路径在两种地形DEM上大致与等高线垂直，较传统D8方法更为准确。
本发明未经描述的技术特征可以通过或采用现有技术实现，在此不再赘述，当然，上述说明并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的普通技术人员在本发明的实质范围内所做出的变化、改型、添加或替换，也应属于本发明的保护范围。

Claims (5)

1.一种基于高程偏差传递的DEM水流方向计算方法，其特征在于，包括以下步骤：
步骤S1：加载DEM数据，DEM为二维数组格式，构建3个与DEM数据尺寸相同的数组A、B、C，数组A用于保存最低相邻单元方向、数组B用于保存第二可能流动方向、数组C用于保存最终确定的流向，并且3个数组中所有位置的初始值都设为0，再构建一个能保存栅格单元三维坐标信息，并根据栅格高程从低到高依次排列的优先队列Q，以及另一个从高到低排列栅格单元三维坐标信息的优先队列T，两个优先队列中同一高程的不同单元按插入先后顺序先后排列，新插入的单元排在同一高程单元的最后；
步骤S2：扫描DEM，将DEM中的有效地形单元加入优先队列T排序，同时将位于有效地形边缘的单元加入优先队列Q进行排序，并同时在数组A、B、C中赋予相同的流向，该流向指向DEM有效地形外侧；
步骤S3：不断取出优先队列Q头部的单元，检索该头部单元的8个相邻单元中的未在数组A中赋予流向的有效单元，为其在数组A中赋予指向此头部单元的方向并加入队列Q，若未赋流向值的单元高程值不高于P，还需为其在数组B与C中赋予指向此头部单元的方向；如果该单元在数组B中对应方向值为初始值0，则比较该单元邻近其最低相邻单元的另外两个相邻单元的高程，将其中较低者所在方向作为第二可能流向并存入数组B，重复此步骤至队列Q为空队列；
步骤S4：读取优先队列T的头部单元K，若单元K在数组C中对应的流向为初始值0，则从单元K开始在最低下坡单元和第二可能方向间使用传递的高程偏差判断从K开始到某个已经具有最终流向的单元的流动路径上各单元的水流方向，然后将K从队列T中移除，重新读取队列T的头部单元以重复此步骤，直至队列T为空队列。
2.根据权利要求1所述的基于高程偏差传递的DEM水流方向计算方法，其特征在于，所述步骤S4中，从单元K开始在最低下坡单元方向和第二可能方向间使用传递的高程偏差判断从K开始到某个已经具有最终流向的单元的流动路径上各单元的水流方向，其具体包括如下步骤：
步骤S4-1：初始化高程偏差传递的参数，以单元K为当前计算单元，以K的高程为h0，初始高程偏差d为0，K的最低下坡单元方向和第二可能方向指向的两个单元相邻，二者一个在K的正方向，一个在K的斜方向，设其中正方向单元高程为h1，斜方向单元高程为h2；
步骤S4-2：构建h0、h1、h2这3个点组成的三角形平面，确定如果水流从点h0出发到达平面底部边界时目标点的高度h3；如果三角形平面的坡度方向超出三角面边缘，水流沿与坡度方向最邻近的三角面边缘方向，否则就沿着坡度方向直线移动，移动线路跟h1、h2连接线的交点所在高程即为h3，具体计算方程如下：
Figure FDA0003734601290000021
步骤S4-3：由于单流向方法是从上游单元中心指向下游单元中心，上游单元水流到达下游单元时到达点的高度并不一定是下游单元中心点的高程值该单元对应的DEM栅格值，水流从自上游到达本单元时的真实高度出发，下降高度与自本单元中心流至下游单元时下降高度h3相同，如此到达的高度作为最终下降高度，比较最低相邻单元与第二可能流动方向对应的单元高程，选取其中高程与最终下降高度最接近的单元作为下游单元，由于在地形走向发生旋转时的水流方向过于复杂，只在上下游单元的最低相邻单元方向和第二可能流动方向都相同时，如果不同则认为水流到达单元为最低相邻单元，到达高程为最低相邻单元的中心点高程；
步骤S4-4：若步骤S4-3中确定的最终流向指向的单元L在数组C中对应的流向不是初始值0，则步骤结束，否则以该被指向的单元L为当前计算单元，以L在数组A、B中对应方向指向的两个单元中，位于正向的单元高程作为新的h1值，斜向的单元高程作为新的h2值，以d*的值作为d的值，返回步骤S4-2重新进行。
3.根据权利要求1所述的基于高程偏差传递的DEM水流方向计算方法，其特征在于，所述步骤S2具体内容为：从左到右、从上到下逐个检索DEM栅格单元，如果某一单元的高程不是无效值，则将其加入优先队列T，并根据两个原则判断有效值单元是否为边缘单元；
步骤S2-1：若该单元位于整个矩形DEM数组的第一行、最后一行、第一列、最后一列，则将该单元加入优先队列Q，并将其在数组A、B和C中对应的单元的流向指向DEM外；
步骤S2-2：若该单元位于DEM范围的上下左右边缘志之外，就按中从P1到P8的顺序检索当前单元P0周围8个相邻单元，若存在某个相邻单元高程为无效地形值，则将该单元加入优先队列Q，并将该单元的流向指向无效值相邻单元，同时记入数组A、B和C对应位置，如果存在多个相邻单元为无效值，则最终流向指向第一个检索到的单元。
4.根据权利要求1所述的基于高程偏差传递的DEM水流方向计算方法，其特征在于，所述步骤S3具体内容为：读取优先队列Q头部的栅格单元P，检索其周边的8个相邻栅格，比较与其在数组A中的方向相邻的两个方向对应的相邻单元的高程，将其中高程较低者所在方向作为第二可能水流方向存入P在数组B中对应的单元，并将所有高程值不是无效值且在数组A中对应的流向是初始值0的相邻单元加入队列Q，同时将在数组A中对应的流向指向P，随后将P从优先队列Q中移除，之后重新读取Q头部的单元并重复此步骤直至队列Q中不再存有单元，此时DEM中的有效地形单元都在数组A和B中被赋予了指向最低相邻单元的流向和第二可能流向；在DEM中存在部分洼地或者平原，这类地形指的是部分DEM栅格的8个相邻单元的高程都不低于中心单元，赋予这类地形向外流的水流路径，对这类地形不做基于偏差传递的处理，直接以此处的方向作为最终方向，即当P的邻近单元对应的流向是初始值0且高程值不高于P时，同时将在数组A、B、C中对应的流向都指向P。
5.根据权利要求2所述的基于高程偏差传递的DEM水流方向计算方法，其特征在于，所述步骤S4-3中，最终下降高度与下游单元中心高程间的差值d进行数学表示，比较h3加上传递高程偏差d的综合高度与h1和h2高度中点的比值，若
Figure FDA0003734601290000031
则最终流向指向h1并计入数组C，新的传递偏差为
Figure FDA0003734601290000032
若
Figure FDA0003734601290000033
则最终流向指向h2并计入数组C，新的传递偏差为
Figure FDA0003734601290000034
CN201910655040.2A 2019-07-19 2019-07-19 一种基于高程偏差传递的dem水流方向计算方法 Active CN110457771B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910655040.2A CN110457771B (zh) 2019-07-19 2019-07-19 一种基于高程偏差传递的dem水流方向计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910655040.2A CN110457771B (zh) 2019-07-19 2019-07-19 一种基于高程偏差传递的dem水流方向计算方法

Publications (2)

Publication Number Publication Date
CN110457771A CN110457771A (zh) 2019-11-15
CN110457771B true CN110457771B (zh) 2022-09-23

Family

ID=68481569

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910655040.2A Active CN110457771B (zh) 2019-07-19 2019-07-19 一种基于高程偏差传递的dem水流方向计算方法

Country Status (1)

Country Link
CN (1) CN110457771B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112989639B (zh) * 2019-12-12 2022-09-16 河海大学 一种基于平均化处理的dem栅格局部排水方向确定方法
CN111125893B (zh) * 2019-12-12 2022-11-04 河海大学 一种基于dem和流量汇集的非弥散水流路径模拟方法
CN114547531B (zh) * 2022-02-25 2022-09-06 中国水利水电科学研究院 一种城市不透水面有效性量化方法
CN117172142B (zh) * 2023-06-27 2024-05-17 长江水利委员会水文局 考虑地形影响的水文模型水流沿程再分配方法
CN117629147B (zh) * 2024-01-25 2024-03-26 北京易控智驾科技有限公司 一种障碍物的检测方法、云控平台及无人车

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105279317A (zh) * 2015-09-30 2016-01-27 西北农林科技大学 一种基于dem的平地河网水流方向估算方法
CN105303612A (zh) * 2014-12-03 2016-02-03 河南理工大学 一种基于不规则三角网模型的数字河网提取方法
CN106981092A (zh) * 2017-03-28 2017-07-25 南京师范大学 基于Priority‑Flood的内流流域提取方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105303612A (zh) * 2014-12-03 2016-02-03 河南理工大学 一种基于不规则三角网模型的数字河网提取方法
CN105279317A (zh) * 2015-09-30 2016-01-27 西北农林科技大学 一种基于dem的平地河网水流方向估算方法
CN106981092A (zh) * 2017-03-28 2017-07-25 南京师范大学 基于Priority‑Flood的内流流域提取方法

Also Published As

Publication number Publication date
CN110457771A (zh) 2019-11-15

Similar Documents

Publication Publication Date Title
CN110457771B (zh) 一种基于高程偏差传递的dem水流方向计算方法
CN111504325B (zh) 一种基于扩大搜索邻域的加权a*算法的全局路径规划方法
CN110487279A (zh) 一种基于改进a*算法的路径规划方法
CN105825510B (zh) 一种兴趣点与道路网的自动配准方法
CN110515094B (zh) 基于改进rrt*的机器人点云地图路径规划方法及系统
CN106525047A (zh) 一种基于floyd算法的无人机路径规划方法
CN107180450A (zh) 一种基于dem的河谷横断面形态的算法
CN106981092B (zh) 基于Priority-Flood的内流流域提取方法
CN102890703B (zh) 一种网络异质多维标度方法
CN101251592A (zh) 一种无线传感器网络的节点定位方法
CN102496187A (zh) 一种基于三角形网格的追踪等值线至边界及断层的方法
CN111858810A (zh) 一种面向道路dem构建的建模高程点筛选方法
CN112965485A (zh) 一种基于二次区域划分的机器人全覆盖路径规划方法
CN101630366A (zh) 大量分块地形数据的动态淹没区提取方法、装置及系统
Ben-Dor et al. On constructing radiation hybrid maps
CN109459052A (zh) 一种扫地机全覆盖路径规划方法
CN105787998B (zh) 一种离散元仿真中多球颗粒的两层网格搜索接触检测方法
CN111177917B (zh) 一种基于srtm的坡长提取方法
CN108875936A (zh) 求解三维空间内任意两个多面体间的最近距离的方法
CN105279317A (zh) 一种基于dem的平地河网水流方向估算方法
CN110555189B (zh) 一种基于反向计算思维的空间插值方法
CN116415318B (zh) 一种基于数学形态学的内流区湖泊水文连通性建模方法
CN102663761A (zh) 用于影像地图的线状矢量与遥感影像自动配准方法
CN115469673B (zh) 一种基于空地信息协同的无人车路径规划方法
CN111858809A (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