CN111243091B - 基于分布式系统的海量dem金字塔切片并行构建方法 - Google Patents

基于分布式系统的海量dem金字塔切片并行构建方法 Download PDF

Info

Publication number
CN111243091B
CN111243091B CN202010271605.XA CN202010271605A CN111243091B CN 111243091 B CN111243091 B CN 111243091B CN 202010271605 A CN202010271605 A CN 202010271605A CN 111243091 B CN111243091 B CN 111243091B
Authority
CN
China
Prior art keywords
tile
slice
key
dem
value
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
CN202010271605.XA
Other languages
English (en)
Other versions
CN111243091A (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.)
Speed China Technology Co Ltd
Original Assignee
Speed Space Time Information Technology Co Ltd
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 Speed Space Time Information Technology Co Ltd filed Critical Speed Space Time Information Technology Co Ltd
Priority to CN202010271605.XA priority Critical patent/CN111243091B/zh
Publication of CN111243091A publication Critical patent/CN111243091A/zh
Application granted granted Critical
Publication of CN111243091B publication Critical patent/CN111243091B/zh
Priority to PCT/CN2021/084081 priority patent/WO2021204026A1/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
    • 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
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00—Geometric image transformations in the plane of the image
    • G06T3/40—Scaling of whole images or parts thereof, e.g. expanding or contracting
    • G06T3/4007—Scaling of whole images or parts thereof, e.g. expanding or contracting based on interpolation, e.g. bilinear interpolation

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Geometry (AREA)
  • Software Systems (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Graphics (AREA)
  • Remote Sensing (AREA)
  • Image Processing (AREA)
  • Image Generation (AREA)

Abstract

本发明公开了一种基于分布式系统的海量DEM金字塔切片并行构建方法，具体步骤：获取待切分DEM的元数据，建立DEM局部坐标系，对DEM进行切分形成若干个可并行计算的切片Tile，计算的Maxzoom，设定的输出切片的宽高、坐标系和采样算法，每个计算节点对所分配的输入切片进行切分计算；在每个计算节点上合并相同Key的cell值；合并不同计算节点上相同Key的cell输出多个Key/Value对；利用8邻域算法修正切片边缘的高程值；再以上一级的输出List<Key，Tile>作为输入，循环计算下一级的zoom、Tile的宽高数据，在各个计算节点上切分、采样、合并输出下一级的切片，直到循环终止。

Description

基于分布式系统的海量DEM金字塔切片并行构建方法
技术领域
本发明涉及地理空间数据并行处理技术领域，具体涉及一种基于分布式系统的海量DEM金字塔切片并行构建方法。
背景技术
目前，随着DEM(数字高程模型Digital Elevation Model)的需求日益增大，高精度、高分辨率DEM数据量呈几何级的增长，DEM金字塔切片的理论与方法也趋于成熟，但如何对海量的DEM数据高效切片处理，是DEM广泛应用亟待解决的关键问题。随着并行计算、分布式计算平台的快速发展和趋于成熟，为海量DEM金字塔切片的并行计算提供高效的算力。结合并行计算和成熟的DEM采样算法、基于四叉树的分割算法实现DEM金字塔切片的并行构建方法，可提高DEM数据计算的效率。
近年来，国内外学者开始重点研究海量空间数据的分布式计算和分析，在构建DEM金字塔方面，一般还是集中在线程并行策略，主要是结合高性能计算机下利用多线程来构建海量DEM金字塔切片，这种方式能够充分发挥机器的性能，但是如果DEM数据量较大，一般的计算机内存不能够承受，达不到计算的效果。关于DEM分析计算方面也有一些专利。中国专利文献089098709949.5公开了一种面向并行数字地形分析海量DEM部署与调度方法，解决了海量DEM的地形分析的效率，此专利未涉及DEM金字塔切片的内容。中国专利文献CN201710107308.X公开了一种基于DEM的内流湖集水区划定方法，其包括获取设定区域内高程数字模型的DEM数据，并对其进行切割，得到基础DEM数据，根据基础DEM数据确立设定区域内的山脊线，并对空间单元进行分割计算空间栅格数据。本专利没有涉及并行和DEM金字塔的算法。中国专利文献201410687961.4公布了一种适于规则网格DEM的结构化多尺度重构方法，同样没有涉及金字塔建立和并行方法。中国专利文献201410788404.1公开了一种面向关系数据库扩展的自适应影像金字塔切片方法，提出了自适应影像切片模型，根据自适应影像切片模型自适应的选择影像切片的大小。影像切片涉及影像元的空值，不涉及DEM切片的拼接方面的方法。在软件方面，主要是商业GIS软件具有DEM金字塔切片的功能，这些方法的主要采用不同的采样算法来切分DEM，未实现分布式的DEM金字塔切片等内容。
因此，有必要开发一种基于分布式系统的海量DEM金字塔切片并行构建方法，能够对大规模或海量DEM的金字塔切片进行并行化切分，能够在高性能计算机集群环境下进行分布式计算，减少金字塔构建的时间，提高计算效率，同时可以结合临近平均值方法计算DEM切片边缘的值，以解决相邻切片间的平滑过渡。
发明内容
本发明要解决的技术问题是提供一种基于分布式系统的海量DEM金字塔切片并行构建方法，能够对大规模或海量DEM的金字塔切片进行并行化切分，能够在高性能计算机集群环境下进行分布式计算，减少金字塔构建的时间，提高计算效率，同时可以结合临近平均值方法计算DEM切片边缘的值，以解决相邻切片间的平滑过渡。
为了解决上述技术问题，本发明采用的技术方案是：该基于分布式系统的海量DEM金字塔切片并行构建方法，具体包括以下步骤：
S1获取源DEM的元数据Metadata；
S2建立所述源DEM的局部坐标系，对所述源DEM进行切分形成若干个可并行计算的切片Tile(切片，就是将一幅DEM切成很多大小一致的小块)；并根据所述源DEM的参数计算每一个所述切片Tile的索引值记为Key，以<Key，Tile>键值对为计算单元组成一个保存源DEM的元数据的链表List<Key，Tile>，将<Key，Tile>按照索引值Key分布到不同的计算节点进行计算；
S3根据是否定义了最大层次MaxZoom，进行设置或计算最大层次MaxZoom，再设置需要的输出切片Tile的参数；再根据设定的坐标系和切分规则计算输出切片Tile的参数，所述输出切片Tile的参数包括切片的长、宽、分辨率等；
S4计算节点并行构建所分配的键值对<Key，Tile>，根据索引值Key合并同一计算节点上的高程值Value；分别根据索引值Key按照所述输出切片Tile的长、宽、分辨率参数计算所对应的切片Tile的外包矩形Extend，计算输入切片Tile中每个网格单元cell的坐标，再根据设定的输出坐标系，计算每个所述网格单元cell在输出坐标系的坐标；根据每个所述网格单元cell的坐标计算相应的所述网格单元cell所在的输出切片Tile的索引值Key，并将其加入一个映射Map中，组成一个映射Map<Key，List<cell>>，其中Key表示切片的索引值，List表示链表，cell表示网格单元；
S5将计算节点上相同索引值Key的链表List<cell>合并；根据索引值Key合并所有计算节点上的链表List<cell>，形成一个键值对<Key，List<cell>>；
S6根据输出切片Tile的参数和坐标系，计算输出切片Tile的网格单元cell的坐标，接着根据采样方式从合并的键值对<Key，List<cell>>中计算网格单元cell的高程值Value；
S7利用8邻域算法，修正输出切片Tile边缘的高程值；根据生成的所述链表List<Key，Tile>，循环获取每个切片Tile的四边网格单元cell和与之相邻的8个切片Tile，取相邻切片Tile的网格单元高程值的平均值作为相应网格单元cell的高程值Value，达到切片Tile拼接的平滑；
S8判断当前层次Zoom是否满足终止条件，若满足终止条件，则结束；若不满足终止条件则根据上一级层次计算下一级层次和下一级层次的切片Tile的外包矩形Extend的宽度和高度，计算方法：NZ＝CZ–1，W＝EW/2(NZ+1)，H＝EH/2(NZ+1)，其中NZ表示下一级层次，CZ表示当前层次，EW、EH表示原始DEM的宽度和高度，W、H表示切片Tile的宽度和高度；
S9根据上一级层次切片Tile和采样算法，并行计算生成下一级层次切片Tile和索引值Key；根据索引值Key、下一级的层次和切片Tile的参数计算当前层次切片Tile所对应的下一级层次的索引值Key，利用所述采样算法从当前层次切片Tile采样下一级层次切片的高程值，构成下一级层次切片Tile的子切片和索引值Key，在每个节点上合并相同的索引值Key，构成键值对<Key，List<Tile>>；接着根据索引值Key，合并所有节点上相同的索引值Key的切片Tile，组成下一级层次的切片Tile，再执行步骤S7；
S10循环执行步骤S8～S9，直到判断层次Zoom满足终止条件，所述终止条件为：若设定了最小层次MinZoom，则判断当前层次Zoom是否达到设定最小层次MinZoom，若达到，则终止；若未设定最小层次MinZoom，若当前层次zoom小于0，则终止。
其中Tile表示切片，Key表示切片的索引值，Value表示高程值，List表示链表，cell表示网格单元，Zoom表示层次，Map表示映射，Extend表示源DEM的外包矩形；采用上述技术方案，获取待切分源DEM的元数据Metadata，作为数据切分和投影转换的基础；建立源DEM左上角为原点，向下向右为正方向的坐标系，根据默认的切片Tile的宽度和高度对源DEM图幅进行切分，根据初步切分的Tile的行列号组成索引值Key，构建一个以索引值Key为键，切片Tile为值的键值对；将键值对分布在不同的计算节点进行以下的计算；根据设定的MaxZoom或计算的MaxZoom，以及设定的输出切片Tile的宽高、坐标系和采样算法，每个计算节点对所分配的输入切片Tile进行切分计算，输出索引值Key和一个包含网格单元cell的链表；接着在每个计算节点上合并相同索引值的网格单元；最后合并不同计算节点上相同的索引值Key的cell，输出多个Key/Value对；根据索引值计算相邻的8个索引值Key，获取相邻切片Tile与当前切片Tile的边缘高程值，将相应的边缘值的平均值作为当前切片Tile的边缘值，达到切片Tile拼接时的平滑；接着循环计算下一级的层次Zoom、切片Tile的宽高等数据，以上一级的输出链表List<Key，Tile>作为输入，计算下一级的链表List<Key，Tile>；循环终止条件是当前层次Zoom小于设定的MinZoom或0；其中采样方式为利用采样算法进行采样，采样算法可以是最邻近像元法，双线性内插法，三次立方卷积法等，均为现有技术。
作为本发明的优选技术方案，所述步骤S1中所述元数据Metadata包括原始分辨率、坐标系，以及在源坐标系下源DEM的外包矩形Extent数据。获取源DEM的Metadata的目的是为了后续计算的需要，这些计算包括将源数据进行初步的切分，以便于后续计算的并行化。
作为本发明的优选技术方案，所述步骤S2中建立源DEM的局部坐标系具体为：建立以源DEM左上角为原点，向下向右为正方向的坐标系；根据设定切片Tile的宽和高，将源DEM数据切分为多个切片Tile，为了将源DEM进行初步的切分，以构成链表List<Key，Tile>，其中索引值Key的计算方法：
Key＝C+"_"+R,C＝(x-XMin)/TW,R＝(y-YMin)/TH
其中，C、R表示切片Tile的行列号，x、y表示切片Tile任意一点的切片坐标，XMin表示源DEM的外包矩形Extent的X最小值，YMin表示源DEM的外包矩形Extent的Y最小值，TH是初始切片Tile的高，TW是初始切片Tile的宽。
作为本发明的优选技术方案，所述步骤S3若未设定层次MaxZoom，则需要计算层次MaxZoom，计算过程是根据源DEM的外包矩形Extend的左下点坐标和右上点的坐标，并将左下点坐标和右上点的坐标进行投影转换，转换为经纬度坐标，其中左下点坐标是(SLon，SLat)，右上点的坐标是(ELon，ELat)；接着根据公式计算层次MaxZoom，计算公式如下：
dLat＝((ELat-SLat)/180.0)·π
dLon＝((ELon-SLon)/180.0)·π
lat1＝(SLat/180.0)·π
lat2＝(ELat/180.0)·π
a＝(sin(dLat/2)·sin(dLat/2))+(sin(dLon/2)·sin(dLon/2)·coslat1·coslat2)
Figure GDA0002539814390000051
MaxZoom＝(log(2πR/dist-tilesize)/log(2))
其中，MaxZoom为最大层次，dLat表示纬度差，dLon表示经度差，lat1表示SLat的弧度，lat2表示ELat的弧度，a表示计算的中间值，dist表示两点的距离，tilesize表示Tile设置的最大边长，R表示WGS84椭球的长半轴，R＝637817。计算输出切片Tile的Extend，需要适当的扩大，目的是为了避免输出切片Tile拼接不会出现裂缝等问题。
作为本发明的优选技术方案，所述步骤S4中网格单元cell＝<x,y,Value>，其中x，y是在切片Tile的坐标，Value是高程值。
作为本发明的优选技术方案，所述步骤S4并行计算输入切片Tile，根据输出切片Tile的参数将输入切片Tile的DEM值划分到不同的输出切片Tile中；其中首先根据输入切片Tile的索引值Key计算Tile的外包矩形Extend，计算公式：
MinX＝ExMinX+col·W
MinY＝ExMinY+row·H
MaxX＝ExMinX+(col+1)·W
MaxY＝ExMinY+(row+1)·H
其中，MinX、MinY、MaxX、MaxY分别表示外包矩形Extend的最小最大XY坐标值，col、row是输入切片Tile的行列号，ExMinX、ExMinY表示源DEM的外包矩形Extend的最小的X、Y坐标值，H表示切片Tile的高，W表示切片Tile的宽。
作为本发明的优选技术方案，所述步骤S7利用8邻域算法，修正切片Tile边缘的高程值Value时，对输出的链表List<Key，Tile>并行计算相邻的切片Tile，取四周除交点外的相邻网格单元cell值，并对上下或左右的cell取值，将取到的cell值取平均值作为指定网格单元cell的高程值Value；四个角点的高程值的计算是取切片Tile的相邻的另外的三个切片Tile与相应角点的值，并取平均值作为相应网格单元cell的高程值Value。
与现有技术相比，本发明具有的有益效果为：该基于分布式系统的海量DEM金字塔切片并行构建方法具有以下优点：利用高性能计算集群，实现了对大规模DEM数据构建DEM金字塔切片的并行化，可减少DEM金字塔构建的时间，提高计算效率；对切片的边缘高程值进行8邻域平滑处理，以避免DEM可视化和分析计算的平滑。
附图说明
图1为本发明的基于分布式系统的海量DEM金字塔切片并行构建方法的工作流程图；
图2为本发明的基于分布式系统的海量DEM金字塔切片并行构建方法的局部坐标系切分DEM示意图；
图3为本发明的基于分布式系统的海量DEM金字塔切片并行构建方法的单个计算节点计算示意图；
图4为本发明的基于分布式系统的海量DEM金字塔切片并行构建方法的多个计算节点合并示意图；
图5为本发明的基于分布式系统的海量DEM金字塔切片并行构建方法边缘平滑示意图；
图6为本发明的基于分布式系统的海量DEM金字塔切片并行构建方法可视化图；
其中：图(a)为17级图；图(b)为17、16级镶嵌显示图；图(c)为17级的Wireframe图；图(d)为17、16级镶嵌Wireframe图。
具体实施方式
下面将结合本发明的实施例图中的附图，对本发明实施例中的技术方案进行清楚、完整的描述。
实施例：如图1～6所示，该基于分布式系统的海量DEM金字塔切片并行构建方法，具体包括以下步骤：
其中Tile表示切片，Key表示切片的索引值，Value表示高程值，List表示链表，cell表示网格单元，zoom表示层次，Map表示映射，Extend表示源DEM的外包矩形，<Key，Tile>表示键值对；
S1使用分辨率30米的源DEM数据和一部分高精度源DEM数据，数据格式是tif，设定的坐标系为WGS84；分别研发基于Flink分布式平台的金字塔切片程序；
获取全国30米源DEM的元数据Metadata，包括原始的分辨率、坐标系WGS84以及Extend；
S2对所述输入源DEM数据建立局部坐标系，设定初始划分切片的高和宽均为126(网格)，利用滑动窗口方法对源数据进行切分，接着根据Tile在局部坐标系的行列号计算Key，构建一系列的<Key，Tile>键值对，按照Flink的算法将键值对数据分布到不同的计算节点上进行计算；建立以源DEM左上角为原点，向下向右为正方向的坐标系；其中索引值Key的计算方法：
Key＝C+"_"+R,C＝(x-XMin)/TW,R＝(y-YMin)/TH
其中，C、R表示切片Tile的行列号，x、y表示切片Tile任意一点的切片坐标，XMin表示源DEM的外包矩形Extent的X最小值，YMin表示源DEM的外包矩形Extent的Y最小值，TH是初始切片Tile的高，TW是初始切片Tile的宽；
S3根据是否定义了最大层次MaxZoom，进行设置或计算最大层次MaxZoom，再根据设定的坐标系和切分规则计算输出切片Tile的参数；所述输出切片Tile的参数包括切片的长、宽、分辨率等；若未设定最大层次MaxZoom，则需要计算最大层次MaxZoom，计算过程是根据源DEM的外包矩形Extend的左下点坐标和右上点的坐标，并将左下点坐标和右上点的坐标进行投影转换，转换为经纬度坐标，其中左下点坐标是(SLon，SLat)，右上点的坐标是(ELon，ELat)；接着根据公式计算最大层次MaxZoom，计算公式如下：
dLat＝((ELat-SLat)/180.0)·π
dLon＝((ELon-SLon)/180.0)·π
lat1＝(SLat/180.0)·π
lat2＝(ELat/180.0)·π
a＝(sin(dLat/2)·sin(dLat/2))+(sin(dLon/2)·sin(dLon/2)·coslat1·coslat2)
Figure GDA0002539814390000081
MaxZoom＝(log(2πR/dist-tilesize)/log(2))
其中，MaxZoom为最大层次，dLat表示纬度差，dLon表示经度差，lat1表示SLat的弧度，lat2表示ELat的弧度，a表示计算的中间值，dist表示两点的距离，tilesize表示Tile设置的最大边长，R表示WGS84椭球的长半轴，R＝637817；
计算输出切片Tile的Extend，需要适当的扩大，目的是为了避免输出切片Tile拼接不会出现裂缝等问题；所述步骤S3中设置输出切片Tile的参数为：宽高均为65(网格)、坐标系为WGS84、索引值Key的命名遵循TMS标准；
S4计算节点并行构建所分配的键值对<Key，Tile>，根据索引值Key合并同一计算节点上的高程值Value；分别根据索引值Key按照所述输出切片Tile的长、宽、分辨率参数计算所对应的切片Tile的外包矩形Extend，计算输入切片Tile中每个网格单元cell的坐标，对初步切分的Tile分别计算；首先根据索引值Key计算Tile的Extend；计算输入切片Tile中每个网格单元cell的坐标，其中cell＝<x,y,Value>，其中x，y是在Tile的坐标，Value是高程值的坐标，根据设定的输出坐标系WGS84，计算cell在输出坐标系的坐标，根据计算结果计算此cell所在的输出切片Tile的索引值Key，并将其加入Map中，组成一个Map<Key，List<cell>>；将计算节点上相同索引值Key的List<cell>合并，原理见图3所示；所述步骤S4并行计算输入切片Tile，根据输出切片Tile的参数将输入切片Tile的DEM值划分到不同的输出切片Tile中；其中首先根据输入切片Tile的索引值Key计算Tile的外包矩形Extend，计算公式：
MinX＝ExMinX+col·W
MinY＝ExMinY+row·H
MaxX＝ExMinX+(col+1)·W
MaxY＝ExMinY+(row+1)·H
其中，MinX、MinY、MaxX、MaxY分别表示外包矩形Extend的最小最大XY坐标值，col、row是输入切片Tile的行列号，ExMinX、ExMinY表示源DEM的外包矩形Extend的最小的X、Y坐标值，H表示切片Tile的高，W表示切片Tile的宽；
S5根据索引值Key合并所有计算节点上的List<cell>；由于输入切片Tile分布在不同的计算节点上，而输出切片Tile可能包含不同的输入切片Tile，也就是说，输出切片Tile的数据可能在不同的计算节点，需要将这些数据拉到同一个节点上进行合并，也就是Flink的Reduce操作，原理见图4所示。
S6根据输出切片Tile的参数和坐标系，计算输出切片Tile的网格单元cell的坐标，使用设置的采样方法，从合并的<Key，List<cell>>中计算cell的高程值，本实施例中使用的是最邻近采样算法；
S7利用8邻域算法，修正切片Tile边缘的高程值Value；根据生成的所述链表List<Key，Tile>，循环获取每个切片Tile的四边网格单元cell和与之相邻的8个切片Tile，取相邻切片Tile的网格单元高程值的平均值作为相应网格单元cell的高程值Value，达到切片Tile拼接的平滑；根据输出的键值对<Key，Tile>，根据索引值Key的TMS命名规则计算每个切片Tile相邻的8个切片Tile，取相邻Tile的cell值，并作平均值作为相应cell的值，达到Tile拼接的平滑；原理示意见图5所示；
S8判断当前层次Zoom是否满足终止条件，若满足终止条件，则结束；若不满足终止条件则根据上一级层次计算下一级层次和下一级层次的切片Tile的外包矩形Extend的宽度和高度，计算方法：NZ＝CZ–1，W＝EW/2(NZ+1)，H＝EH/2(NZ+1)，其中NZ表示下一级层次，CZ表示当前层次，EW、EH表示原始DEM的宽度和高度，W、H表示切片Tile的宽度和高度；
S9根据上一级层次切片Tile和采样算法，并行计算生成下一级层次切片Tile和索引值Key；根据索引值Key、下一级的层次和切片Tile的参数计算当前层次切片Tile所对应的下一级层次的索引值Key，利用所述采样算法从当前层次切片Tile采样下一级层次切片的高程值，构成下一级层次切片Tile的子切片和索引值Key，在每个节点上合并相同的索引值Key，构成键值对<Key，List<Tile>>；接着根据索引值Key，合并所有节点上相同的索引值Key的切片Tile，组成下一级层次的切片Tile，再执行步骤S7；本实施例中使用的是最邻近采样算法；
S10循环执行步骤S8～S9，计算下一级的List<Key，Tile>，直到判断当前层次Zoom满足终止条件，所述终止条件为：若设定了最小层次MinZoom，则判断当前层次Zoom是否达到设定最小层次MinZoom，若达到，则终止；若未设定最小层次MinZoom，则当前层次Zoom小于0，则终止。
以上所述仅为本发明的较佳实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种基于分布式系统的海量DEM金字塔切片并行构建方法，其特征在于，具体包括以下步骤：
S1获取源DEM的元数据Metadata；
S2建立所述源DEM的局部坐标系，对所述源DEM进行切分形成若干个可并行计算的切片Tile；并根据所述源DEM的参数计算每一个所述切片Tile的索引值记为Key，以<Key，Tile>键值对为计算单元组成一个保存源DEM的元数据的链表List<Key，Tile>，将<Key，Tile>按照索引值Key分布到不同的计算节点进行计算；
S3根据是否定义了最大层次MaxZoom，进行设置或计算最大层次MaxZoom，设置需要的输出切片Tile的参数；根据设定的坐标系和切分规则计算输出切片Tile的参数；
S4计算节点并行构建所分配的键值对<Key，Tile>，根据索引值Key合并同一计算节点上的高程值Value；分别根据索引值Key按照所述输出切片Tile的长、宽、分辨率参数计算所对应的切片Tile的外包矩形Extend，计算输入切片Tile中每个网格单元cell的坐标，再根据设定的输出坐标系，计算每个所述网格单元cell在输出坐标系的坐标；根据每个所述网格单元cell的坐标计算相应的所述网格单元cell所在的输出切片Tile的索引值Key，并将其加入一个映射Map中，组成一个映射Map<Key，List<cell>>；
S5将计算节点上相同索引值Key的链表List<cell>合并；根据索引值Key合并所有计算节点上的链表List<cell>，形成一个键值对<Key，List<cell>>；
S6根据输出切片Tile的参数和坐标系，计算输出切片Tile的网格单元cell的坐标，接着根据采样方式从合并的键值对<Key，List<cell>>中计算网格单元cell的高程值Value；
S7利用8邻域算法，修正输出切片Tile边缘的高程值；根据生成的所述链表List<Key，Tile>，循环获取每个切片Tile的四边网格单元cell和与之相邻的8个切片Tile，取相邻切片Tile的网格单元高程值的平均值作为相应网格单元cell的高程值Value，达到切片Tile拼接的平滑；
S8判断当前层次Zoom是否满足终止条件，若满足终止条件，则结束；若不满足终止条件则根据上一级层次计算下一级层次和下一级层次的切片Tile的外包矩形Extend的宽度和高度，计算方法：NZ＝CZ–1，W＝EW/2(NZ+1)，H＝EH/2(NZ+1)，其中NZ表示下一级层次，CZ表示当前层次，EW、EH表示原始DEM的宽度和高度，W、H表示切片Tile的宽度和高度；
S9根据上一级层次切片Tile和采样算法，并行计算生成下一级层次切片Tile和索引值Key；根据索引值Key、下一级的层次和切片Tile的参数计算当前层次切片Tile所对应的下一级层次的索引值Key，利用所述采样算法从当前层次切片Tile采样下一级层次切片的高程值，构成下一级层次切片Tile的子切片和索引值Key，在每个节点上合并相同的索引值Key，构成键值对<Key，List<Tile>>；接着根据索引值Key，合并所有节点上相同的索引值Key的切片Tile，组成下一级层次的切片Tile，再执行步骤S7；
S10循环执行步骤S8～S9，直到判断层次Zoom满足终止条件，所述终止条件为：若设定了层次MinZoom，则判断当前层次Zoom是否达到设定层次MinZoom，若达到，则终止；若未设定层次MinZoom，若当前层次zoom小于0，则终止。
2.根据权利要求1所述的基于分布式系统的海量DEM金字塔切片并行构建方法，其特征在于，所述步骤S1中所述元数据Metadata包括原始分辨率、坐标系，以及在原坐标系下DEM的外包矩形Extent数据。
3.根据权利要求1所述的基于分布式系统的海量DEM金字塔切片并行构建方法，其特征在于，所述步骤S2中建立源DEM的局部坐标系具体为：建立以源DEM左上角为原点，向下向右为正方向的坐标系；根据设定切片Tile的宽和高，将源DEM数据切分为多个切片Tile，为了将源DEM进行初步的切分，以构成链表List<Key，Tile>，其中索引值Key的计算方法：
Key＝C+"_"+R,C＝(x-XMin)/TW,R＝(y-YMin)/TH
其中，C、R表示切片Tile的行列号，x、y表示切片Tile任意一点的切片坐标，XMin表示源DEM的外包矩形Extent的X最小值，YMin表示源DEM的外包矩形Extent的Y最小值，TH是初始切片Tile的高，TW是初始切片Tile的宽。
4.根据权利要求3所述的基于分布式系统的海量DEM金字塔切片并行构建方法，其特征在于，所述步骤S3若未设定层次MaxZoom，则需要计算层次MaxZoom，计算过程是根据源DEM的外包矩形Extend的左下点坐标和右上点的坐标，并将左下点坐标和右上点的坐标进行投影转换，转换为经纬度坐标，其中左下点坐标是(SLon，SLat)，右上点的坐标是(ELon，ELat)；接着根据公式计算层次MaxZoom，计算公式如下：
dLat＝((ELat-SLat)/180.0)·π
dLon＝((ELon-SLon)/180.0)·π
lat1＝(SLat/180.0)·π
lat2＝(ELat/180.0)·π
a＝(sin(dLat/2)·sin(dLat/2))+(sin(dLon/2)·sin(dLon/2)·coslat1·coslat2)
Figure FDA0002539814380000031
MaxZoom＝(log(2πR/dist-tilesize)/log(2))
其中，MaxZoom为最大层次，dLat表示纬度差，dLon表示经度差，lat1表示SLat的弧度，lat2表示ELat的弧度，a表示计算的中间值，dist表示两点的距离，tilesize表示Tile设置的最大边长，R表示WGS84椭球的长半轴，R＝637817。
5.根据权利要求3所述的基于分布式系统的海量DEM金字塔切片并行构建方法，其特征在于，所述步骤S4中网格单元cell＝<x,y,Value>，其中x，y是在切片Tile的坐标，Value是高程值。
6.根据权利要求5所述的基于分布式系统的海量DEM金字塔切片并行构建方法，其特征在于，所述步骤S4并行计算输入切片Tile，根据输出切片Tile的参数将输入切片Tile的DEM值划分到不同的输出切片Tile中；其中首先根据输入切片Tile的索引值Key计算Tile的外包矩形Extend，计算公式：
MinX＝ExMinX+col·W
MinY＝ExMinY+row·H
MaxX＝ExMinX+(col+1)·W
MaxY＝ExMinY+(row+1)·H
其中，MinX、MinY、MaxX、MaxY分别表示外包矩形Extend的最小最大XY坐标值，col、row是输入切片Tile的行列号，ExMinX、ExMinY表示源DEM的外包矩形Extend的最小的X、Y坐标值，H表示切片Tile的高，W表示切片Tile的宽。
7.根据权利要求5所述的基于分布式系统的海量DEM金字塔切片并行构建方法，其特征在于，所述步骤S7利用8邻域算法，修正切片Tile边缘的高程值Value时，对输出的链表List<Key，Tile>并行计算相邻的切片Tile，取四周除交点外的相邻网格单元cell值，并对上下或左右的cell取值，将取到的cell值取平均值作为指定网格单元cell的高程值Value；四个角点的高程值的计算是取切片Tile的相邻的另外的三个切片Tile与相应角点的值，并取平均值作为相应网格单元cell的高程值Value。
CN202010271605.XA 2020-04-09 2020-04-09 基于分布式系统的海量dem金字塔切片并行构建方法 Active CN111243091B (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN202010271605.XA CN111243091B (zh) 2020-04-09 2020-04-09 基于分布式系统的海量dem金字塔切片并行构建方法
PCT/CN2021/084081 WO2021204026A1 (zh) 2020-04-09 2021-03-30 基于分布式系统的海量dem金字塔切片并行构建方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010271605.XA CN111243091B (zh) 2020-04-09 2020-04-09 基于分布式系统的海量dem金字塔切片并行构建方法

Publications (2)

Publication Number Publication Date
CN111243091A CN111243091A (zh) 2020-06-05
CN111243091B true CN111243091B (zh) 2020-07-24

Family

ID=70871339

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010271605.XA Active CN111243091B (zh) 2020-04-09 2020-04-09 基于分布式系统的海量dem金字塔切片并行构建方法

Country Status (2)

Country Link
CN (1) CN111243091B (zh)
WO (1) WO2021204026A1 (zh)

Families Citing this family (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111243091B (zh) * 2020-04-09 2020-07-24 速度时空信息科技股份有限公司 基于分布式系统的海量dem金字塔切片并行构建方法
CN112381715B (zh) * 2020-11-16 2024-04-09 航天科工(北京)空间信息应用股份有限公司 海量遥感影像并行生成地图瓦片的方法和装置
CN112883130A (zh) * 2021-01-05 2021-06-01 中国人民解放军国防科技大学 基于空间数据库输出地图高程数据的方法及设备、介质
CN113724229B (zh) * 2021-08-31 2023-07-28 北京英视睿达科技股份有限公司 高程差的确定方法、装置及电子设备
CN113989454B (zh) * 2021-12-30 2022-03-18 山东省地质测绘院 适用于地质数据和地理信息数据的融合方法、装置及系统
CN114357104A (zh) * 2022-02-25 2022-04-15 北京嘉诚瑞杰信息技术有限公司 一种无索引地图切片聚合方法
CN116302488B (zh) * 2023-01-17 2024-01-23 重庆市地理信息和遥感应用中心(重庆市测绘产品质量检验测试中心) 一种地形图数据坐标多进程自动识别转换方法
CN117724850B (zh) * 2024-02-06 2024-04-19 天津先进技术研究院 野外预通行路径可行性评估方法、系统、设备及介质
CN117874301B (zh) * 2024-03-11 2024-05-31 浙江省气象台 基于网格数据金字塔切片的处理、存储、调用的实现方法

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP1819044B1 (en) * 2001-06-26 2009-04-22 Nokia Corporation A multi-level quantizer with current mode DEM switch matrices and separate DEM decision logic for a multibit sigma delta modulator
CN103034999A (zh) * 2012-12-05 2013-04-10 上海创图网络科技发展有限公司 一种基于dem数据的金字塔切片方法
CN103455624A (zh) * 2013-09-16 2013-12-18 湖北文理学院 一种轻量级全球多维遥感影像网络地图服务实现方法
CN106815807A (zh) * 2017-01-11 2017-06-09 重庆市地理信息中心 一种基于gpu‑cpu协同的无人机影像快速镶嵌方法
US9723269B2 (en) * 2010-07-07 2017-08-01 Pictometry International Corp. Real-time moving platform management system
CN107229742A (zh) * 2017-06-20 2017-10-03 华中科技大学 一种基于遥感大数据确定城市易涝点的方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8775142B2 (en) * 2010-05-14 2014-07-08 Conocophillips Company Stochastic downscaling algorithm and applications to geological model downscaling
CN111243091B (zh) * 2020-04-09 2020-07-24 速度时空信息科技股份有限公司 基于分布式系统的海量dem金字塔切片并行构建方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP1819044B1 (en) * 2001-06-26 2009-04-22 Nokia Corporation A multi-level quantizer with current mode DEM switch matrices and separate DEM decision logic for a multibit sigma delta modulator
US9723269B2 (en) * 2010-07-07 2017-08-01 Pictometry International Corp. Real-time moving platform management system
CN103034999A (zh) * 2012-12-05 2013-04-10 上海创图网络科技发展有限公司 一种基于dem数据的金字塔切片方法
CN103455624A (zh) * 2013-09-16 2013-12-18 湖北文理学院 一种轻量级全球多维遥感影像网络地图服务实现方法
CN106815807A (zh) * 2017-01-11 2017-06-09 重庆市地理信息中心 一种基于gpu‑cpu协同的无人机影像快速镶嵌方法
CN107229742A (zh) * 2017-06-20 2017-10-03 华中科技大学 一种基于遥感大数据确定城市易涝点的方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
An Efficient Tile-Pyramids Building Method for Fast Visualization of Massive Geospatial Raster Datasets;Ning GUO等;《Advances in Electrical and Computer Engineering》;20160430;第16卷(第4期);第3-8页 *
一种支持云计算的遥感影像数据组织模型研究;赖积保等;《计算机科学》;20130731;第40卷(第7期);第80-83页 *
面向Spark的遥感影像金字塔模型的并行构建方法;黄冬梅等;《计算机应用与软件》;20170531;第34卷(第5期);第175-181页 *

Also Published As

Publication number Publication date
WO2021204026A1 (zh) 2021-10-14
CN111243091A (zh) 2020-06-05

Similar Documents

Publication Publication Date Title
CN111243091B (zh) 基于分布式系统的海量dem金字塔切片并行构建方法
CN105118090B (zh) 一种自适应复杂地形结构的点云滤波方法
CN110827398B (zh) 基于深度神经网络的室内三维点云自动语义分割方法
Chen et al. Multiscale grid method for detection and reconstruction of building roofs from airborne LiDAR data
CN113781667B (zh) 三维结构简化重建方法、装置、计算机设备和存储介质
CN102930594B (zh) 基于保持边界的计算机三维地形网格的绘制方法
CN112115534A (zh) 一种将三维房屋模型转换带高度属性的二维矢量面的方法
CN108022307B (zh) 基于增材再制造点云模型的自适应平面分层方法
CN103700088A (zh) 一种基于可变形图结构表示的图像集无监督共分割方法
CN115797592B (zh) 一种基于倾斜摄影三维模型自动生成建筑物体块的方法和装置
US20230274466A1 (en) Point cloud polar coordinate coding method and device
CN102799898A (zh) 大背景下高分辨率遥感影像目标识别高效处理方法
CN105931256A (zh) 基于cuda的大幅面遥感影像快速分割方法
CN114117702A (zh) 基于点云的输电线路自动化逆向建模方法
CN110570445A (zh) 一种图像分割方法、装置、终端及可读介质
CN103679808A (zh) 利用二维核磁图像重建头部三维模型的方法及系统
CN114119902A (zh) 一种基于无人机倾斜三维模型的建筑物提取方法
CN114332291A (zh) 一种倾斜摄影模型建筑物外轮廓规则提取方法
CN108986218A (zh) 一种基于pmvs的建筑物密集点云快速重建方法
CN113379919A (zh) 一种基于无人机rgb相机的植被冠层高度快速提取方法
CN113221895A (zh) 小目标检测方法、装置、设备及介质
CN116721228B (zh) 一种基于低密度点云的建筑物高程提取方法及系统
CN109816781A (zh) 一种基于图像细节和结构增强的多视图立体几何方法
CN108269228B (zh) 基于gpu并行计算的无人机影像拉花区域自动探测方法
CN105653881A (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
CP01 Change in the name or title of a patent holder
CP01 Change in the name or title of a patent holder

Address after: 210042 8 Blocks 699-22 Xuanwu Avenue, Xuanwu District, Nanjing City, Jiangsu Province

Patentee after: Speed Technology Co.,Ltd.

Address before: 210042 8 Blocks 699-22 Xuanwu Avenue, Xuanwu District, Nanjing City, Jiangsu Province

Patentee before: SPEED TIME AND SPACE INFORMATION TECHNOLOGY Co.,Ltd.