EN中文

地球信息科学学报 >

2023 , Vol. 25 >Issue 2: 239 - 251

DOI: https://doi.org/10.12082/dqxxkx.2023.220725

地球信息科学理论与方法

正二十面体六边形全球离散格网编码运算

  • 赵龙 , 1, 2 ,
  • 李国庆 1, 2 ,
  • 姚晓闯 , 3, * ,
  • 马跃 1, 2
展开
  • 1.中国科学院空天信息创新研究院，北京 100094
  • 2.中国科学院大学电子电气与通信工程学院，北京 100049
  • 3.中国农业大学土地科学与技术学院，北京 100193
* 姚晓闯（1986— ），男，河南郑州人，副教授，主要从事空间大数据技术及农业应用方向研究。 E-mail:

赵 龙（1992— ），男，安徽宿州人，博士生，主要从事空间数据管理方面的研究。E-mail:

收稿日期: 2022-09-25

  修回日期: 2022-11-11

  网络出版日期: 2023-04-19

基金资助

中国科学院战略性先导科技专项（A类）(XDA19020103)

收起

Code Operation Scheme for the Icosahedral Hexagonal Discrete Global Grid System

  • ZHAO Long , 1, 2 ,
  • LI Guoqing 1, 2 ,
  • YAO Xiaochuang , 3, * ,
  • MA Yue 1, 2
Expand
  • 1. Aerospace Information Research Institute, Chinese Academy of Sciences, Beijing 100094, China
  • 2. School of Electronic, Electrical and Communication Engineering, University of Chinese Academy of Sciences, Beijing 100049, China
  • 3. College of Land Science and Technology, China Agricultural University, Beijing 100193, China
* YAO Xiaochuang, E-mail:

Received date: 2022-09-25

  Revised date: 2022-11-11

  Online published: 2023-04-19

Supported by

Strategic Priority Research Program of the Chinese Academy of Sciences, Grant(XDA19020103)

Fold

摘要

全球离散格网系统是指把地球表面按照一定规则离散分割成多分辨率层次结构的格网单元，广泛应用于海量多源空间数据的组织、管理和分析中。六边形全球离散格网具有优良的几何特性，非常适合于空间数据的处理，如何进一步提高六边形全球离散格网编码运算的效率仍是当前研究的重点。本文采用正二十面体施耐德投影四孔径六边形全球离散格网模型，基于六边形三轴坐标与编码的二进制数的对应关系构建四孔六边形的基础编码结构，将二十面体划分为32个基础六边形，并将之分为3种基础六边形剖分瓦片，在每个六边形剖分瓦片采用基础编码结构进行编码，建立了四孔六边形全球离散格网编码，同时设计了并实现了四孔六边形编码与六边形三轴坐标之间的快速转换，基于此构建了一种高效的四孔六边形全球离散格网编码运算方案，包括编码的算数运算、空间拓扑运算和邻域检索运算及跨面运算。与现有的六边形全球离散格网编码运算方案相比，本文的方案进一步提高了编码算数运算、空间拓扑运算和邻域检索运算的效率，编码加法运算是HLQT的2~3倍，邻域检索运算分别是HLQT的3~5倍和H3的2~3倍，且受格网编码层次的影响较小，编码的跨面邻域检索运算时间略高于面内的运算，可以为全球离散格网系统的研究应用提供支撑。

关键词： 全球离散格网系统; 正二十面体; 六边形; 格网编码; 编码运算; 六边形三轴坐标系; 坐标转换; 格网邻域检索

本文引用格式

赵龙 , 李国庆 , 姚晓闯 , 马跃 . 正二十面体六边形全球离散格网编码运算[J]. 地球信息科学学报, 2023 , 25(2) : 239 -251 . DOI: 10.12082/dqxxkx.2023.220725

Abstract

The discrete global grid system refers to the discrete partitioning of the earth's surface into grid cells with multi-resolution hierarchical structure according to certain rules, which is widely used in organization, management, and analysis of massive multi-source spatial data. The hexagonal global discrete grid has excellent geometric properties and is well suited for spatial data processing. However, how to further improve the efficiency of the hexagonal global discrete grid coding operation is still the focus of current research. In this paper, we adopt the model of icosahedral snyder equal-area projection aperture 4 hexagonal discrete global grid system and construct the base coding structure of aperture 4 hexagon based on the correspondence between the hexagonal triaxial coordinates and the coded binary numbers, consisting of 7 base digits in the first layer and 4 base digits int the other layer. We divide the icosahedron into 3 base hexagonal subdivision tiles according to the different subdivision structures and adopt the base coding structure for coding scheme in each hexagonal subdivision tile to establish the aperture 4 hexagonal discrete global grid coding scheme. Besides, this paper designs and implements a fast conversion between aperture 4 hexagonal code and hexagonal triaxial coordinates, based on which an efficient aperture 4 hexagonal discrete global grid encoding operation scheme is constructed, including arithmetic operation of encoding, spatial topology operation, and neighbourhood retrieval operation and cross-plane operation of encoding. Compared with the existing hexagonal discrete global grid coding scheme, the coding scheme proposed in this paper has fewer base code digits, is more concise, and facilitates faster conversion to the hexagonal triaxial coordinates of the grid. Compared with the existing coding operation scheme, the proposed scheme further improves the efficiency of coding arithmetic operation, spatial topology operation, and neighbourhood retrieval operation. The coding addition operation is 2~3 times more efficient than HLQT. The neighbourhood retrieval operation is 3~5 times and 2~3 times more efficient than HLQT and H3, respectively, and is less affected by the coding level of the grid coding. The proposed coding scheme in this paper has the same efficiency of additive operation and subtractive operation, and the efficiency of spatial topology operation is 2 times that of arithmetic operation. The coding cross-plane neighbourhood retrieval operation time is slightly longer than that of the in-plane operation, and the impact on the overall operation time is not significant. This study provides support for the research application of discrete global grid system.

Key words： discrete global grid systems; icosahedron; hexagon; grid coding; code operation; hexagonal triaxial coordinate system; encoding and coordinate conversion; neighbourhood retrieval operation

1 引言

全球离散格网系统（Discrete Global Grid System，DGGS）是指把地球表面按照一定规则离散分割成多分辨率层次结构的格网单元，具有层次性、全球连续性等特点，可以突破传统的基于投影的平面格网模型的限制[1-2]，并能够实现海量、多源异构、多分辨率的对地观测数据的融合[3]。目前全球离散格网系统已经被广泛应用于海量多源空间数据的组织、管理和分析中，如地理空间数据管理[4⇓-6]、全球邮政建设[7]、全球油棕榈监测[8]、综合环境分析[9]、仇恨犯罪预测[10]等。
编码运算是全球离散格网系统的核心，按照编码原理可以划分为层次编码运算、填充曲线编码运算和整数坐标编码运算[11]。六边形全球离散格网系统主要采用层次编码运算，可以分为三孔径、四孔径和七孔径3种编码结构。在三孔径六边形编码方面，Sahr[12⇓-14]在“一般平衡三进制”（Generalized Balanced Ternary，GBT）六边形结构的基础上构建了“改进的广义平衡三进制”（Modified Generalized Balanced Ternary，MGBT）编码方案，并将其扩展到二十面体上；Peterson[15]提出了PYXIS方案来编码和索引正二十面体三孔径六边形全球离散格网系统；Vince[16]和 Zheng[17]提出了基于PYXIS的编码运算和傅里叶变换算法。在四孔径六边形编码方面，童晓冲[18]等将2种类型的四孔径六边形格网系统组合在一起，构建了六边形四元平衡结构（Hexagonal Quaternary Balanced Structure，HQBS）；根据四孔径六边形格网系统的特点，王蕊等[19-20]设计了六边形格点四叉树结构（Hexagon Lattice Quad Tree, HLQT），编码运算的效率优于PYXIS和HQBS。在七孔径六边形编码方面，目前的研究较少，其中Uber H3将正二十面体划分为122个基本单元（110个六边形和12个五边形），并分配了从0到121的数字，在每个基本单元上按七孔径进行剖分，采用十六进制编码来标识每一个格网[21]。目前的六边形全球离散格网层次编码运算规则大多采用归纳总结并产生编码运算表的形式构建，而整数坐标编码运算的效率高于基于归纳总结的层次编码运算，并有很大差距[22]，如何才能更进一步提高六边形层次编码运算的效率。本文计划将六边形层次编码转换为整数坐标编码来实现编码运算，其中编码与整数坐标的转换是影响运算效率最大的因素，因此如何构建可与整数坐标快速转换的层次编码，并实现编码与坐标的高效转换算法是本文研究的重点。
正二十面体四孔径六边形全球离散格网具有格网变形小[23]、六边形剖分层次的方向一致等特点，不需要考虑奇偶层格网方向变换的问题，便于实现层次编码与整数坐标的快速转换。本文以正二十面体四孔六边形全球离散格网为基础模型，基于六边形三轴坐标与编码的二进制数的对应关系构建四孔六边形的基础编码结构，并将其扩展到正二十面体上，设计了四孔六边形全球离散格网编码方案，实现了编码与六边形三轴坐标之间的快速转换，并基于此构建了编码的算数运算、空间拓扑运算、邻域检索运算及跨面运算，最后通过对比实验证明本文方案的有效性和高效性。

2 四孔六边形全球离散格网编码与编码运算方案设计

本文研究中采取的技术路线如图1所示，分为六边形全球离散格网编码和编码运算2个部分，下面将分为5个小节进行详细阐述。
图1 正二十面体六边形全球离散格网编码运算技术路线

Fig. 1 Technical route of code operation scheme for the icosahedral hexagonal discrete global grid system

2.1 基于六边形三轴坐标系的平面六边形格网编码

为了便于表示和计算全球离散格网编码，我们采用六边形三轴坐标系（ i j k坐标系）来描述和定位每一个六边形格网，六边形三轴坐标系的原点位于初始六边形的中心点，3个坐标轴之间的夹角均为120°，且分别与六边形的右上边、左上边和正下方边垂直，如图2所示。每个象限的格网坐标均由相邻的两个坐标轴决定，另外一个坐标轴的值为0。
图2 不同格网层次上的六边形三轴坐标系

Fig. 2 Hexagonal triaxial coordinate system at different grid levels

每一层格网的六边形三轴坐标系统之间的单位均不相同，第 n层格网坐标系统的 i、 j、 k 3个坐标轴的基本单位 u i n ,   u j n ,   u k n满足如下公式：
u i n - 1 ,   u j n - 1 ,   u k n - 1 = 1 2 u i n ,   1 2 u j n ,   1 2 u k n u i n ,   u j n ,   u k n = 3 2 n L ,   3 2 n L ,   3 2 n L
式中： u i n 、 u j n 、 u k n和 u i n - 1 、 u j n - 1 、 u k n - 1分别是 第 n层和第 n - 1层格网的 i、 j、 k 3个坐标轴的基本坐标单位； L为图2中初始六边形的边长。
由于四孔径六边形格网并不能实现完全的嵌套剖分，为了能构建覆盖初始六边形格网的剖分编码结构，本文将编码分为第1层和第2层—第n层编码，并基于 i j k坐标与编码的二进制数的对应关系构建四孔六边形的基础编码结构，便于实现编码与六边形三轴坐标系间的高效转换。其中第1层格网的 i   j   k坐标从 i轴开始沿逆时针方向对应的二进制数依次为{100，110，010，011，001，101}，对应的基础码元为{4，6，2，3，1，5}，如图3和表1所示，第2层—第n层编码沿 i轴、 j轴和 k轴方向的二进制数依次为{100，010，001}，对应的基础码元为{4，2，1}，如图4和表2所示，最终可构建如 图5所示的四孔六边形编码。
图3 四孔六边形格网第一层编码

Fig. 3 The first code of four-aperture hexagonal grid

表1 第1层格网编码与 i   j   k坐标对应

Tab. 1 Correspondence table between first code and i j k coordinate

( i ,   j ,   k ) ( 0,0 , 0 ) ( 0,0 , 1 ) ( 0,1 , 0 ) ( 0,1 , 1 ) ( 1,0 , 0 ) ( 1,0 , 1 ) ( 1,1 , 0 )
与 i轴夹角 无 4 π 3 2 π 3 π 0 5 π 3 π 3
二进制 000 001 010 011 100 101 110
a 1 0 1 2 3 4 5 6
图4 四孔六边形格网第2层—第n层编码

Fig. 4 The 2 ~ nth code of four-aperture hexagonal grid

表2 第2层—第 n层格网编码与 i j k坐标对应

Tab. 2 Correspondence table between 2 ~ nth code and i j k coordinate

( i , j , k ) ( 0,0 , 0 ) ( 0,0 , 1 ) ( 0,1 , 0 ) ( 1,0 , 0 )
与 i轴夹角 无 4 π 3 2 π 3 0
二进制 000 001 010 100
a 2 ~ a n 0 1 2 4
图5 四孔六边形编码

Fig. 5 Aperture 4 hexagonal encoding

2.2 平面六边形格网编码在正二十面体上的扩展

为了构建全球离散格网编码，需要将平面六边形编码结构在正二十面体上展开。正二十面体有20个正三角形，我们可以将其划分为三角面上完整的六边形面瓦片和由顶点三角形组成的六边形顶点瓦片，如图6所示，其中六边形面瓦片共有20个，编号为0~19，六边形顶点瓦片共有12个，编号为20~31。通过上述方法，将正二十面体划分为32个基础六边形瓦片，然后在每个基础六边形上进行剖分和编码。
图6 基础六边形瓦片结构

Fig. 6 Basic hexagonal tile structure diagram

在剖分的过程中，不同瓦片之间存在一些跨面的子六边形，为了便于编码和计算，我们根据瓦片相接的结构将其划分到特定的瓦片中，如图7所示，最终形成4种不同类型的瓦片结构。1号瓦片结构保留了右上角和正下方的跨面子六边形，2号瓦片保留了右下角的跨面子六边形，3号瓦片保留了所有的跨面子六边形，同时去除正上方三角形内的六边形，4号瓦片则去除了正下方三角形内的六边形。
图7 第3层六边形编码方案在正二十面体上的扩展

Fig. 7 The extension of the third-level hexagonal encoding scheme on the icosahedron

在4种瓦片结构的基础上，构建了图8所示的3种瓦片编码结构。图8（a）为1号瓦片所对应的编码结构，图8（b）为2号瓦片所对应的编码结构，因1号瓦片与2号瓦片所在的三角面方向相反，故图8（a）与图8（b）的编码方向也相反。3号与4号瓦片结构相同，方向相反，故采用同一种编码结构图8（c）。
图8 3种六边形剖分瓦片与编码结构

Fig. 8 Three Hexagonal tile and coding structures

2.3 六边形格网编码与六边形三轴坐标的快速转换

从上文中六边形三轴坐标系的定义可知，全球离散格网编码与六边形三轴坐标之间是一一对应关系，且四叉树全球离散格网编码本身就隐含了格网之间的空间位置关系，故可以根据编码的原理实现编码 a 1 a 2   .   .   .   a n与 i   j   k坐标之间的快速转换。全球离散格网编码转换为六边形三轴坐标的运算如算法1所示：首先确定第1层码元 a 1与 i j k坐标的对应关系，根据表3实现第1层码元的转换，然后确定第2层—第 n层码元 a 2   .   .   .   a n与 i   j   k坐标的对应关系，根据表4实现第 2 层 — 第 n层码元的转换并累加，最后将转换结果修正为一个坐标值为0的 i ,   j ,   k坐标，如算法2所示。
算法1：全球离散格网编码转换为六边形三轴坐标
1. 输入：格网编码 a 1 a 2   .   .   .   a n
2. 根据表3实现第1层码元的转换
i 1 , j 1 , k 1 = 2 n - 1 ∙ φ 1 a 1 （2）
式中： n表示编码的长度， i 1 , j 1 , k 1表示第一层格网所对应的六边形三轴坐标， a 1表示首层码元 φ 1 ( a 1 )如表3所示。
3. 根据表4实现第 2 ~ n层码元的转换并累加；
i , j , k = i 1 , j 1 , k 1 + ∑ r = 2 n 2 n - r + 1 ∙ φ 2 a r （3）
式中： n表示编码的长度， i , j , k表示第 n层格网所对应的六边形三轴坐标， a r表示第 r个码元， φ 2 ( a 2 )在表4中给出。
4. 修正转换结果
m = m i n i , j , k ;   i = i - m ;   j = j - m ;   k = k - m （4）
5. 输出：六边形三轴坐标 ( i , j , k )
表3 第1层码元对应的 i j k坐标 φ 1

Tab. 3 The i   j   k coordinates φ 1 associated with first symbol

编码 0 1 2 3 4 5 6
φ 1 ( 0,0 , 0 ) ( 0,0 , 1 ) ( 0,1 , 0 ) ( 0,1 , 1 ) ( 1,0 , 0 ) ( 1,0 , 1 ) ( 1,1 , 0 )
表4 第2层—第n层码元对应的 i   j   k坐标 φ 2

Tab. 4 The i j k coordinates φ 2 associated with 2 — nth symbols

编码 0 1 2 4
φ 2 ( 0,0 , 0 ) ( 0,0 , 1 ) ( 0,1 , 0 ) ( 0,1 , 1 )
算法2：六边形三轴坐标转为全球离散格网编码
1. 输入：六边形三轴坐标 ( i , j , k )；
2. f o r   r = n   t o   2
3. ( v i , v j , v k ) = ( i 1 , j 1 , k 1 )
4. a r = ω v i , v j , v k
5. i f   Ɵ ( v i , v j , v k )
6. i = i ≫ 1 ;   j = j ≫ 1 ;   k = k ≫ 1
7. e l s e
8. i = ( i + 1 ) ≫ 1 ;   j = ( j + 1 ) ≫ 1 ;   k = ( k + 1 ) ≫ 1
9. a 1 = ω 0 v i , v j , v k
10. 输出：格网编码 a 1 a 2   .   .   .   a n
在六边形三轴坐标转为全球离散格网编码的运算中，首先计算第2层—第n层中每一层码元对应的 i j k坐标，共有7种坐标，而编码只有4种，故需要在原有4种对照关系基础上加入3种新的对照关系，如表5所示，根据表5确定相应的码元，并在出现这3种坐标时对六边形三轴坐标 i , j , k进行修正，然后再根据表6确定第1层码元，最后将所有码元合并为格网编码 a 1 a 2   .   .   .   a n。
表5 第2层—第n层中 i   j   k坐标对应的码元 ω

Tab. 5 The 2 — nth symbol ω associated with the i   j   k coordinate

( i , j , k ) ( 0,0 , 0 ) ( 0,0 , 1 ) ( 0,1 , 0 ) ( 0,1 , 1 ) ( 1,0 , 0 ) ( 1,0 , 1 ) ( 1,1 , 0 )
ω 0 1 2 4 4 2 1
Ɵ true true true false true false false
表6 第1层中 i j k坐标对应的码元 ω 0

Tab. 6 The first symbol ω 0 associated with the i j k coordinate

( i , j , k ) ( 0,0 , 0 ) ( 0,0 , 1 ) ( 0,1 , 0 ) ( 0,1 , 1 ) ( 1,0 , 0 ) ( 1,0 , 1 ) ( 1,1 , 0 )
ω 0 0 1 2 3 4 5 6
其中， n表示编码的长度， i , j , k表示第 n层格网所对应的六边形三轴坐标， v i , v j , v k表示第 n个码元的二进制坐标， a r表示第 r个码元， ω v i , v j , v k和 Ɵ ( v i , v j , v k )如表5所示， a 1表示第 1个码元， ω 0 v i , v j , v k如表6所示。

2.4 六边形全球离散格网编码运算

全球离散格网编码记录了六边形之间的相对位置关系，每一个格网编码都隐含了它相对于中心点的方位和距离，可以被视为一个向量并定义相关的运算规则。本文将编码运算分为3类：编码的算数运算、编码的空间拓扑运算以及编码的邻域检索运算。编码的算数运算主要包括加法和减法运算，与向量运算类似，编码的算数运算满足平行四边形法则，如算法3所示。编码的空间拓扑运算是计算2个格网单元之间的最小连通距离，用于描述格网之间的连通性和邻接性，如算法4所示。编码的邻域检索运算是计算格网单元的6个相邻格网的编码，如算法5所示。本文利用全球离散格网编码与六边形三轴坐标之间的转换算法来构建编码的算数运算、空间拓扑运算以及邻域检索运算（图9）。
算法3：编码算数运算
1. 输入：格网编码 a 1 a 2   .   .   .   a n和 b 1 b 2   .   .   .   b n
2. 利用算法1，分别计算 a 1 a 2   .   .   .   a n和 b 1 b 2   .   .   .   b n所在格网
的六边形三轴坐标 i 1 , j 1 , k 1和 i 2 , j 2 , k 2
3. i 3 , j 3 , k 3 = i 1 + i 2 , j 1 + j 2 , k 1 + k 2或者
i 3 , j 3 , k 3 = i 1 - i 2 , j 1 - j 2 , k 1 - k 2
4. m = m i n i 3 , j 3 , k 3
5. i 3 , j 3 , k 3 = i 3 - m , j 3 - m , k 3 - m
6. 利用算法2，将六边形三轴坐标 i 3 , j 3 , k 3转换为格网
编码 c 1 c 2   .   .   .   c n
7. 输出：格网编码
算法4：编码空间拓扑运算
1. 输入：格网编码 a 1 a 2   .   .   .   a n和 b 1 b 2   .   .   .   b n
2. 利用算法1，分别计算 a 1 a 2   .   .   .   a n和
b 1 b 2 … b n 所在格网的六边形三轴坐标
i 1 , j 1 , k 1和 i 2 , j 2 , k 2；
3. i 3 , j 3 , k 3 = i 1 - i 2 , j 1 - j 2 , k 1 - k 2
4. m = m i n i 3 , j 3 , k 3；
5. i 3 , j 3 , k 3 = i 3 - m , j 3 - m , k 3 - m
6. L = m a x i 3 , j 3 , k 3 ;
7. 输出： 2个格网之间的最短连通距离
算法5：编码的邻域检索运算
1. 输入：格网编码 a 1 a 2   .   .   .   a n
2. 利用算法1，计算 a 1 a 2   .   .   .   a n所在格网的六边形三轴
坐标 i , j , k
3. 分别计算格网的6个相邻格网的六边形三轴坐标
i 1 , j 1 , k 1 = i + 1 ,   j ,   k               i 2 , j 2 , k 2 = i + 1 ,   j + 1 ,   k i 3 , j 3 , k 3 = i ,   j + 1 ,   k             i 4 , j 4 , k 4 = i ,   j + 1 ,   k + 1 i 5 , j 5 , k 5 = i ,   j ,   k + 1               i 6 , j 6 , k 6 = i + 1 ,   j ,   k + 1 （5）
4. f o r   r = 1   t o   6
5. m = m i n i r , j r , k r
6. i r , j r , k r = i r - m , j r - m , k r - m
7. 利用算法2，将六边形三轴坐标 i r , j r , k r转换为格网编码 a 1 r a 2 r … a n r
8. 输出：相邻的6个格网编码，，
，，，
图9 全球离散格网编码邻域检索运算

Fig. 9 Discrete global grid code heighborhood retrieval operation

2.5 六边形全球离散格网编码跨面运算

上文中的全球离散格网编码运算均是在每一个单独的瓦片上构建的，而不同瓦片上的编码运算是不可或缺的，本文利用不同瓦片上的六边形三轴坐标系之间的转换来实现编码的跨面运算，如算法6所示。首先将待运算的编码转换到同一个六边形三轴坐标系下，然后利用 i   j   k坐标进行相应的运算，最后将运算的结果转换到原六边形三轴坐标系下，并将坐标转换为格网编码，本文以编码的跨面邻域检索运算为例。对于一个六边形瓦片来说，可以根据六边形边界和坐标轴将其划分为8个边界区域，如图10所示。8个边界区域上的格网可以划分为临边、沿边和跨边3种类型，每一种类型的格网的 i 、 j 、 k坐标均可以利用三轴坐标方程来进行归纳和表达，如图11和表7所示，同时不同瓦片的坐标系之间也可以实现转换。在此基础上，可以构建全球离散格网编码的跨面邻域检索运算算法。
算法6：编码的跨面邻域检索运算
1. 输入：格网编码 a 1 a 2   .   .   .   a n
2. 利用算法1，计算 a 1 a 2   .   .   .   a n所在格网的六边形三轴坐标         i , j , k
3. 根据不同的瓦片类型及表2中格网的六边形三轴坐标所 属位置，计算在同一瓦片内的相邻格网的六边形三轴坐 标，并利用算法2将其转换为格网编码
4. 将 i , j , k坐标转换为以相邻瓦片为中心的坐标系内的六 边形三轴坐标 i 0 , j 0 , k 0；
5. 计算 i 0 , j 0 , k 0的相邻格网的六边形三轴坐标，并利用算 法2将其转换为格网编码
6. 输出：相邻的6个格网编码 a 11 a 21   .   .   .   a n 1， a 12 a 22   .   .   .   a n 2，
a 13 a 23   .   .   .   a n 3， a 14 a 24   .   .   .   a n 4， a 15 a 25   .   .   .   a n 5， a 16 a 26   .   .   .   a n 6
图10 跨面检索运算的六边形瓦片分区

Fig. 10 Hexagonal tile partition map for cross-face retrieval operations

图11 跨面格网单元与 i j k方程

Fig. 11 Cross-face grid cells and the i j k equation

表7 跨面格网单元所对应的方程

Tab. 7 Table of equations corresponding to cross-face grid cells

分区 临边 沿边 跨边
① 2×i-j-2n+2=0 2×i-j-2n+1=0 2×i-j-2n=0
② i+j+2-2n=0 i+j+1-2n=0 i+j+1-2n=0
③ 2×j-i-2n+2=0 2×j-i-2n+1=0 2×j-i-2n=0
④ 2×j-k-2n+2=0 2×j-k-2n+1=0 2×j-k-2n=0
⑤ j+k-2n+2=0 j+k-2n+1=0 j+k-2n=0
⑥ 2×k-j-2n+2=0 2×k-j-2n+1=0 2×k-j-2n=0
⑦ 2×k-i-2n+2=0 2×k-i-2n+1=0 2×k-i-2n=0
⑧ i+k-2n+2=0 i+k-2n+1=0 i+k-2n=0
⑨ 2×i-k-2n+2=0 2×i-k-2n+1=0 2×i-k-2n=0

注：表中n表示全球离散格网系统第n层，i、j、k分别表示第n层格网所在六边形三轴坐标系中的坐标，ijk方程与图10中内容相对应。

3 实验与对比分析

目前，研究与应用中最常用的六边形全球离散格网编码方案有PYXIS、HQBS、HLQT和Uber H3[16,18⇓⇓-21]，根据之前的研究，HLQT的代码加法效率优于HQBS和PYXIS，HLQT的邻域检索效率优于HQBS[19-20,24]。为了验证本文编码运算的效率和优势，本文将与HLQT之间进行加法效率的对比，并与HLQT、Uber H3进行邻域检索效率的对比，分别设计了全球离散格网编码方案的对比、编码算数运算和空间拓扑运算的效率对比和编码邻域检索运算的效率对比。本文所有试验的计算机配置如下，操作系统：Windows 11 x64，硬件配置：11th Gen Intel(R) Core(TM) i5-11300H @3.10 GHz，16 GB RAM， SKHynix 512GB SSD。

3.1 全球离散格网编码方案的对比

本文将编码分为第1层和第2层—第n层编码， 并基于 i j k坐标与编码的二进制数的对应关系构建四孔六边形的基础编码结构，第1层由{0，1，2，3，4，5，6} 7个基础码元组成，第2层—第n层由{0， 1， 2，4}4个基础码元组成，HLQT为了实现编码运算在正二十面体上的扩 展， 引入了首位码元的概念，首位码元的集合E = {0, 1, 2, 3, 4, 5, 6, 10, 20, 30, 40, 50, 60, 100, 200, 300, 400, 500, 600, 106, 601, 201, 102, 302, 203, 403, 304, 504, 405, 605, 506, 1000, 2000, 3000, 4000, 5000, 6000}，通过对比可发现本文编码更为简洁，且无需利用逗号来隔开首位编码（图12）。
图12 本文编码方案与HLQT的对比

Fig. 12 Comparison of the encoding scheme in this article with HLQT

3.2 编码的算数运算和空间拓扑运算的效率对比

在编码的算数运算和空间拓扑运算实验中，从每个格网层次中随机选择1 000 000对编码进行HLQT加法运算、本文加法运算、减法运算和空间拓扑运算，统计其在第5层—第12层的运算时间，实验结果见表8，编码算数运算和空间拓扑运算的效率对比如图13所示。
表8 编码算数运算和空间拓扑运算的实验结果

Tab. 8 Experimental results of arithmetic operations and spatial topological operations (ms)

层次 HLQT加法
运算		 本文加法
运算		 本文减法
运算		 本文拓扑
运算
5 77 34 37 15
6 99 39 41 15
7 116 45 42 19
8 141 47 52 21
9 168 56 58 24
10 195 63 64 29
11 220 63 67 28
12 244 69 72 30
图13 编码算数运算和空间拓扑运算的效率对比

Fig. 13 Efficiency comparison of arithmetic operations and spatial topology operations

通过以上实验可知：① 本文编码方案在第5层—第12层的加法运算效率是HLQT的2~3倍，且随着全球离散格网编码层次的增加，本文方案运算效率所受到的影响较小。主要是因为HLQT每一位码元的计算都可能会在左侧多个加法位上产生进位码元，且首位码元的设计较为复杂，增加了运算的次数，而本文方案每一位码元都只运算一次，计算复杂度为 O ( N )。② 本文方案的加法运算效率与减法运算效率相同，空间拓扑运算效率是算数运算的2倍。加法与减法运算的步骤基本相同，运算时间也相当，而空间拓扑运算与算数运算相比，不需要将计算的六边形三轴坐标转换为格网编码。

3.3 编码的邻域检索运算的效率对比

为了验证本文编码的邻域检索运算及跨面运算的效率，本文设计了编码的邻域检索运算实验，从每个格网层次中随机选择1 000 000个编码进行HLQT、H3和本文面内及跨面的编码邻域检索运算，统计其在第5层—第12层的运算时间，实验结果如表9所示。图14为编码邻域检索运算的效率对比结果。
表9 编码邻域检索运算的实验结果

Tab. 9 Experimental results of neighborhood retrieval operations (ms)

层次 HLQT H3 本文方案 本文方案跨面
5 458 404 152 186
6 584 403 162 203
7 697 398 185 213
8 836 423 204 234
9 970 403 224 247
10 1130 418 251 270
11 1261 403 258 281
12 1424 396 278 296
图14 编码邻域检索运算的效率对比

Fig. 14 Efficiency comparison of coded neighborhood retrieval operations

通过以上实验可知:① 本文编码方案在第5 层—第12层的邻域检索运算效率是HLQT的3~5倍，H3的2~3倍。随着全球离散格网编码层次的增加，HLQT的编码运算时间快速增加，H3基本保持不变，本文方案则缓慢增加。HLQT和本文的邻域检索运算均是由6次加法运算组成，本文方案简化了编码转换为六边形三轴坐标的过程，H3在编码运算的过程中需要考虑到初始的122个基本单元和七孔径六边形剖分中的方向变化，增加了运算的时间，同时直接采用二进制数进行运算减少了层次的影响。② 编码的跨面邻域检索运算时间略高于面内的运算时间，主要是因为编码的跨面运算加入了编码跨面位置的确定与不同瓦片坐标系的转换，相较于编码与六边形三轴坐标之间的转换，其对运算时间的影响并不大。

4 结论

本文基于六边形三轴坐标与编码的二进制数的对应关系构建四孔六边形的基础编码结构，将二十面体划分为3种瓦片结构，构建正二十面体四孔径六边形全球离散格网编码，并利用六边形三轴坐标系实现编码与六边形三轴坐标的快速转换，基于此构建编码的算数运算、空间拓扑运算、邻域检索运算和跨面运算。通过与现有方案对比，本文方案有以下优点：
（1）本文基于六边形三轴坐标系构建四孔六边形编码结构，第一层由7个基础码元组成，第2层—第n层由4个基础码元组成，并把二十面体划分为 3种类型的瓦片，实现编码在正二十面体上的扩展。与现有编码方案相比，本文编码方案的基础码元更少，编码更为简洁，并便于与格网所在的六边形三轴坐标的快速转换。
（2）本文编码方案在第5—第12层的加法运算效率是HLQT的2~3倍，且随着全球离散格网编码层次的增加，本文方案运算效率所受到的影响较小，加法运算效率与减法运算效率相同，空间拓扑运算效率是算数运算的2倍。
（3）本文编码方案在第5—第12层的邻域检索运算效率是HLQT的3~5倍，H3的2~3倍，编码的跨面邻域检索运算时间略高于面内的运算时间，相较于编码与六边形三轴坐标之间的转换，其对运算时间的影响并不大，可以为全球离散格网上的应用提供技术支撑。

参考文献

原文顺序 | 文献年度倒序 | 文中引用次数倒序
[1]
Bauer-Marschallinger B, Sabel D, Wagner W. Optimisation of global grids for high-resolution remote sensing data[J]. Computers & Geosciences, 2014, 72:84-93. DOI:10.1016/j.cageo.2014.07.005

DOI

[2]
Lin B X, Zhou L C, Xu D P, et al. A discrete global grid system for earth system modeling[J]. International Journal of Geographical Information Science, 2018, 32(4):711-737. DOI:10.1080/13658816.2017.1391389

DOI

[3]
Yao X C, Li G Q, Xia J S, et al. Enabling the big earth observation data via cloud computing and DGGS: Opportunities and challenges[J]. Remote Sensing, 2019, 12(1):62. DOI:10.3390/rs12010062

DOI

[4]
Li M, Stefanakis E, Mcgrath H. National Terrain Data Management on Discrete Global Grids in Canada[C]. AutoCarto2020, 2020.

[5]
李爽, 程承旗, 童晓冲, 等. 基于多级信息网格的海量遥感数据存储管理研究[J]. 测绘学报, 2016, 45(S1):106-114.

[ Cheng C Q, Tong X C, et al. A study on data storage and management for massive remote sensing data based on multi-level grid model[J]. Acta Geodaetica et Cartographica Sinica, 2016, 45(S1):106-114. ] DOI:10.11 947/j.AGCS.2016.F013

DOI

[6]
Li S, Pu G L, Cheng C Q, et al. Method for managing and querying geo-spatial data using a grid-code-array spatial index[J]. Earth Science Informatics, 2019, 12(2):173-181. DOI:10.1007/s12145-018-0362-6

DOI

[7]
Adams B. Wāhi, a discrete global grid gazetteer built using linked open data[J]. International journal of digital earth, 2017, 10(5):490-503. DOI:10.1080/17538947.2016. 1229819

DOI

[8]
Cheng Y, Yu L, Zhao Y, et al. Towards a global oil palm sample database: Design and implications[J]. International Journal of Remote Sensing, 2017, 38(14):4022-32. DOI:10.1080/01431161.2017.1312622.

DOI

[9]
Robertson C, Chaudhuri C, Hojati M, et al. An Integrated Environmental Analytics System (IDEAS) based on a DGGS[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2020, 162:214-228. DOI:10.1016/j.isprsjpr s.2020.02.009

DOI

[10]
Jendryke M, Mcclure S C. Spatial prediction of sparse events using a discrete global grid system; a case study of hate crimes in the USA[J]. International Journal of Digital Earth, 2021:1-17. DOI:10.1080/17538947.2021.1886 356

DOI

[11]
赵学胜, 贲进, 孙文彬, 等. 地球剖分格网研究进展综述[J]. 测绘学报, 2016, 45(S1):1-14.

[ Zhao X S, BEN J, SUN W B, et al. Overview of the Research Progress in the Earth Tessellation Grid[J]. Acta Geodaetica et Cartographica Sinica, 2016, 45(S1):1-14. ] DOI: 10.11947/j.AGCS.2016.F001

DOI

[12]
Sahr K. Location coding on icosahedral aperture 3 hexagon discrete global grids[J]. Computers, Environment and Urban Systems, 2008, 32(3):174-187. DOI:10.1016/j.compenvurbsys.2007.11.005

DOI

[13]
Sahr K. Hexagonal discrete global grid systems for geospatial computing[J]. Archiwum Fotogrametrii, Kartografii i Teledetekcji, 2011, 22:363-76.

[14]
Sahr K M. Discrete global grid systems: A new class of geospatial data structures[M].  Oregon: University of Oregon ProQuest Dissertations Publishing, 2005.

[15]
Peterson P.Close-packed uniformly adjacent, multiresolutional overlapping spatial data ordering[P]. US Patent No 8018458,2011-09-13.

[16]
Vince A, Zheng X. Arithmetic and Fourier transform for the PYXIS multi-resolution digital Earth model[J]. International Journal of Digital Earth, 2009, 2(1):59-79. DOI:10.1080/17538940802657694

DOI

[17]
Zheng X. Efficient Fourier transforms on hexagonal arrays[D]. Florida: University of Florida, 2007.

[18]
Tong X C, Ben J, Wang Y, et al. Efficient encoding and spatial operation scheme for aperture 4 hexagonal discrete global grid system[J]. International Journal of Geographical Information Science, 2013, 27(5):898-921. DOI:10.1080/13658816.2012.725474

DOI

[19]
王蕊, 贲进, 杜灵瑀, 等. 正二十面体四孔六边形格网系统编码运算[J]. 武汉大学学报·信息科学版, 2020, 45(1):89-96.

[ Wang R, Ben J, Du L Y, et al. Code operation scheme for the icosahedral aperture 4 hexagonal grid system[J]. Geomatics and Information Science of Wuhan University, 2020, 45(1):89-96. ] DOI:10.13203/j.whugis20180191

DOI

[20]
王蕊, 贲进, 杜灵瑀, 等. 平面四孔六边形格网系统编码运算[J]. 测绘学报, 2018, 47(7):1018-1025.

[ Wang R, Ben J, Du L Y, et al. Encoding and operation for the planar aperture 4 hexagon grid system[J]. Acta Geodaetica et Cartographica Sinica, 2018, 47(7):1018-1025. ] DOI:10.11947/j.AGCS.2018.20170374

DOI

[21]
Brodsky I. H3: Uber's hexagonal hierarchical spatial index[EB/OL]. [2018-06-27].

[22]
Zhao L, Li G Q, Yao X C, et al. Comparison and analysis of hexagonal discrete global grid coding[C]. International Conference on Spatial Data and Intelligence. Springer, Cham, 2021:127-133. DOI:10.1007/978-3-030-85462-1_11

DOI

[23]
White D, Kimerling A J, Sahr K, et al. Comparing area and shape distortion on polyhedral-based recursive partitions of the sphere[J]. International Journal of Geographical Information Science, 1998, 12(8):805-827. DOI:10.10 80/136588198241518

DOI

[24]
Wang R, Ben J, Zhou J B, et al. A generic encoding and operation scheme for mixed aperture three and four hexagonal discrete global grid systems[J]. International Journal of Geographical Information Science, 2021, 35(3):513-555. DOI:10.1080/13658816.2020.1763363

DOI

Options
文章导航

/

〈 〉