球面菱形離散格網正二十面體剖分法

周良辰,盛業華,林冰仙,閭國年,趙志鵬

1.南京師范大學虛擬地理環境教育部重點實驗室,江蘇南京 210023;2.江蘇省地理信息資源開發與利用協同創新中心,江蘇南京 210023

球面菱形離散格網正二十面體剖分法

周良辰1,2,盛業華1,2,林冰仙1,2,閭國年1,2,趙志鵬1,2

1.南京師范大學虛擬地理環境教育部重點實驗室,江蘇南京 210023;2.江蘇省地理信息資源開發與利用協同創新中心,江蘇南京 210023

基于菱形單元的球面離散格網具有結構簡單和方向性一致等優點,使得其適用于全球多尺度空間數據建模與分析。本文提出一種球面菱形離散格網正二十面體剖分法。首先根據地球的主要地理特征確定正二十面體各個頂點在球面上的位置;其次針對南北兩極處的極點奇異性,建立了球面經緯度坐標與三維直角坐標之間的轉換關系,在此基礎上提出球面菱形離散格網層次遞歸剖分算法,并從菱形格網單元面積、長短軸比兩方面對剖分結果進行幾何形變分析;最后驗證了本文方法的可行性,試驗結果表明,該方法可有效地構建基于正二十面體的球面菱形離散格網,其總體幾何特性明顯優于基于正八面體構建的球面菱形離散格網。

全球離散格網;菱形;正二十面體;層次遞歸剖分;幾何形變

1 引 言

全球離散格網(discrete global gird,DGG)主要研究如何將地球表面剖分成一個等面積、等形狀、具有多分辨率的空間層次格網結構[1],采用每個格網對應的地址碼代替地理坐標在球面上進行操作。由于地址碼既表示位置,又表達了比例尺和精度,因此已被廣泛應用于全球尺度空間數據集成與模式計算中[2-3]。最常見的DGG構建方法是以經緯弧線為基礎,使得基于經緯坐標系的數據及算法無需進行坐標變換即可使用,但該格網存在單元面積從赤道向南北兩極逐漸收斂造成格網粒度不均勻以及極點收斂等問題[4]。另一類常見的DGG構建方法是基于5種基本柏拉圖立體,目前研究較多的是基于正八面體或正二十面體的球面三角形離散格網[5-6]、球面菱形離散格網[7-8]與球面六邊形離散格網[9-10]。該類DGG因具有粒度近似多層次多分辨率等特性,目前成為研究的重點。

球面菱形離散格網具有類似于平面正方形格網的幾何結構,在剖分方法與編碼機制上比三角形或六邊形格網簡單,且具有方向性一致、徑向對稱和平移相和性,易于實現鄰近搜索和層次關聯等空間操作,適應于全球空間數據管理與模式計算[11]。正八面體6個頂點可以定位在兩個極點和4個等分的赤道點上,任意的經緯度坐標能很容易地定位到其中一個面上[12]。相關學者對球面菱形離散格網正八面體剖分方法、編碼模型、空間數據建模以及空間分析等方面均進行了較為深入的研究[5,7-8]。然而,正八面體的形狀與球體相差較大,因此格網系統的單元形狀會出現較大變形(圖1)。相比之下,正二十面體是最接近球面的柏拉圖立體[12],基于正二十面體所構建的DGG將具有更均勻的幾何性質,有利于全球空間數據建模與模式模擬[13],因此逐步引起相關研究者的關注[14-16]。

圖1 正八面體和正二十面體對比Fig.1 Comparison between octahedron and icosahedron

剖分方法是構建DGG的基礎,以菱形作為基本格網單元的DGG,通常采用四叉樹剖分策略[19],但正二十面體的幾何性質、拓撲性質比正八面體更為復雜,無法直接采用現有的球面菱形離散格網正八面體剖分方法[5,7-8]。

2 正二十面體與球面的對應關系

球面菱形離散格網正二十面體剖分首先要確定正二十面體各個頂點在球面上的位置。文獻[17]分別在兩極各放置一個頂點,其中一條通過北極的邊線和0°經線重合(圖2(a))。圖2(b)、(c)將大部分多面體頂點放在海洋上,減少頂點對陸地數據的影響[16-18]。地理現象客觀存在著水平地帶性特征(即經度地帶性和緯度地帶性的合稱),同時地球是以南北兩極所確定的地軸為中心進行自轉。后一種放置方法只能針對陸地應用,且沒能顧及地球的自轉特性。因此本文采用文獻[17]提出的方法建立正二十面體與球面的對應關系。將正二十面體中頂點位于極點處的三角形與其極點對邊相鄰的三角形進行合并,可以將正二十面體劃為10個初始菱形格網單元。這些格網單元沿赤道對稱,符合地球南北兩極的方向性、對稱性和水平地帶性特征,便于后續的剖分、編碼與地理坐標轉換。

圖2 正二十面體的3種不同放置方法Fig.2 Icosahedron placed in three different ways

3 菱形格網單元層次遞歸剖分

初始菱形格網單元確定后,即可進行格網單元的層次遞歸剖分。球面菱形離散格網單元形狀、大小近似規則、均勻,拓撲連接關系類似于正方形格網,可以按照四叉樹或其他方式進行層次遞歸剖分[19]。正二十面體結構的復雜性以及南北兩極處的極點奇異性給格網單元與球面經緯度的轉換帶來問題。分析發現,球面菱形離散格網單元的幾何位置既可以用經緯坐標表達,也可以由以球心為原點的三維直角坐標表達,且在三維直角坐標系中不存在極點收斂,有利于實現球面菱形離散格網的層次遞歸剖分。

3.1 球面經緯度坐標與三維直角坐標的相互轉換

3.1.1 球面經緯度坐標向三維直角坐標轉換

設地球半徑為R,格網點(i,j)的經緯度坐標分別為G(i,j)＝(λ(i,j),φ(i,j)),對應的空間直角坐標為C(i,j)＝(X(i,j),Y(i,j),Z(i,j)),借助空間幾何關系即可建立從球面經緯度坐標向三維直角坐標的坐標轉換關系如下

3.1.2 三維直角坐標向球面經緯度坐標轉換

球面格網點的三維直角坐標向球面經緯度坐標轉換,是將一個不收斂空間轉換到一個收斂空間。對于非極點處(即Z(i,j)＝?{R,－R}),可基于式(2)建立轉換關系。地球南北兩極處(當Z(i,j)＝R或－R時),極點奇異性帶來一對多的映射沖突,但與極點同屬于一條格邊的另一格點并不存在映射沖突的問題,且這兩個格點處在同一條經線上,即λ(i,j＋1)＝λ(i,j),可基于式(3)建立轉換關系

3.2 格網單元層次遞歸剖分算法

基于球面經緯度坐標與三維直角坐標的相互轉換關系,進行球面菱形離散格網層次遞歸剖分,算法主要步驟包括(圖3):①根據正二十面體與球面的對應關系,確定初始菱形格網的球面經緯坐標;②將初始菱形格網從球面經緯坐標映射為三維直角坐標,根據剖分層次,采用大圓弧平分方法進行層次遞歸剖分,建立剖分層級與對應比例尺的關系(表1);③將剖分結果從三維直角坐標轉換至球面經緯度坐標。

圖3 基于正二十面體的球面菱形離散格網層次遞歸剖分算法流程Fig.3 The recursive subdivision algorithms flow of spherical diamond discrete grids based on icosahedron

表1 格網剖分層級對應的地圖比例尺Tab.1 The relationship between subdivision level of grid and map scale

4 菱形格網單元幾何形變分析

從球面幾何可知,沒有一種剖分方法能使球面格網在每個層次上具有像平面柵格那樣完全相同的幾何特征(如:面積、長度、角度),只能達到近似相等[12,20]。目前對于球面離散格網的幾何特性研究通常以對格網單元幾何特性的統計分析為主[20]。下面分別從球面菱形離散格網單元面積比和標準差、長/短軸比(角度)和標準差為度量指標對本文提出的球面菱形離散格網正二十面體剖分的結果進行分析(表2)。

4.1 最大/小面積比和面積標準差

本文采用球面大圓弧剖分方法,菱形單元的每條邊均為球面大圓弧,因此菱形單元的面積等同于兩個球面三角形面積之和,即SABCD＝SABCD＋SBCD。根據式(4)即可求得球面三角形的面積和格網單元面積的標準差

式中,∠A、∠B、∠C、∠D為菱形單元的4個頂點的內角;σ表示格網面積的標準差;Si表示每個菱形單元的面積;S表示該剖分層次下所有菱形格網

的平均面積;N表示該層次下菱形單元的個數。

表2 正二十面體球面菱形格網統計指標Tab.2 The statistical indicators of spherical diamond grid based on icosahedron

從表2可以看出,隨著剖分層次的增加,面積比首先增大然后逐漸收斂到1.3左右,標準差先增加后減,逐漸收斂到0(圖4)。說明隨著剖分層次的增加,格網單元的大小能夠保持近似均勻。同時,菱形單元的面積分布和面積標準差分布在初始菱形中由內到外呈現出類似于五邊形的分布規律,五邊形外部的菱形單元面積標準差大于內部(圖5)。

圖4 最大/最小面積比和面積標準差Fig.4 Max/min area ratio and area standard deviation of diamond cell

圖5 菱形單元面積和面積標準差全球分布圖Fig.5 The global distribution of diamond element area and the area＇s standard deviation

4.2 長/短軸比和軸長標準差

菱形單元的長短軸之比和標準差指標的組合可以用來衡量格網單元角度變形的情況[26]。對于球面菱形的長短軸可以作如下的定義,長軸方向即為連接豎直方向上兩個頂點的大圓弧,短軸為水平方向上兩個頂點的大圓弧,長軸、短軸球面長度可定義為

式中,OA、OC、OB、OD為球心至菱形格網單元4個頂點的距離。

同樣,從表2可以看出,最長最短軸長比和軸長標準差表現出了與面積特性極為相似的變化的特征,即隨著剖分層次的增加逐漸收斂穩定(圖6),這一特性表明隨著剖分層級的增加,正二十面體格網能夠保持穩定的形狀變形特征。格網的形狀以及面積越穩定、均勻,越有利于格網編碼的實現以及數據的集成。從軸長比的分布圖和軸長的分布圖,可以發現:軸長比和軸長標準差的分布呈現菱形分布的特點,菱形水平方向兩頂點處的變形較大;豎直方向兩頂點處的變形較小(圖7)。

圖6 長短軸比和軸長標準差Fig.6 The axial ratio and axial standard deviation of diamond cell

4.3 與其他球面菱形離散格網剖分方法對比

為了進一步分析基于正二十面體的球面菱形離散格網(正二十大圓弧)的幾何形變性質,將其與正八面體大圓弧平分法[19](正八大圓弧法)、經緯線平分法[7](正八經緯法)以及大圓弧和經緯線混合平分法[8](正八混合法)所得到的3種球面菱形離散格網進行幾何性質對比(圖8)。對比分析可以看出,基于正二十面體的球面菱形格網隨著剖分層次的遞增,角度變形和面積變形均保持為最小,進而表明格網的幾何性質更好,更有利于空間數據集成和提高空間分析的精度。這顯然是由于正二十面體相較于正八面體更加接近球面。

圖7 菱形單元長短軸比和長短軸比標準差全球分布圖Fig.7 The global distribution of diamond element’s ratio of major/minor axis and the ratio’s standard deviation

圖8 不同剖分方法的幾何形變對比Fig.8 Comparisons geometry distortion between different subdivision methods

5 實例驗證與分析

本文在VC＋＋2013環境下,基于OpenGL開發了試驗系統。以90 m分辨率數字高程模型SRTM數據的球面空間建模為例,將本文所提出的剖分模型與前文介紹的3種基于正八面體的球面菱形離散格網剖分模型進行效果對比(圖9)。

圖9 4種球面菱形離散格網建模效果對比Fig.9 Comparisions of modeling effect between four spherical diamond discrete grids

從圖9中可以看出:正八面體大圓弧剖分方法在同一基礎菱形區域內格網單元分布較均勻,但在不同區域間的格網單元面積有較明顯的變化,格網整體分布不均勻(圖9(a)、(b));正八面體經緯網剖分方法在不同區域間變化非常明顯,相對于正八面體大圓弧剖分方法,面積形變更大(圖9(c)、(d));正八面體混合剖分方法在同一基礎菱形區域內的格網單元分布比較均勻,在不同區域間的格網單元面積有較明顯的變化(圖9(e)、(f)),并且在某些區域內存在突變的情況,如圖9(f)中紅色箭頭所指的地方;本文剖分方法得到的格網整體分布均勻,格網單元面積形變小、排列緊致,不存在明顯的形狀突變(圖9(g)、(h))。對比以上4種格網的可視化效果,可以看出本文所提出的基于正二十面體的球面菱形格網在空間數據建模與集成方面具有明顯的優勢。

6 結論與展望

本文提出了一種球面菱形離散格網正二十面體剖分法。從地球的主要地理特征出發確定正二十面體與球體的對應關系;針對南北兩極處的極點奇異性,設計了收斂的球面經緯度坐標與不收斂的三維直角坐標的轉換關系,實現了球面菱形離散格網的層次遞歸剖分與地理坐標確定;基于菱形格網單元的最大/最小面積比和面積標準差、菱形長/短軸比和軸長標準差進行了剖分單元幾何變形性質的分析。分析結果表明,與現有各類球面菱形離散格網剖分方法相比,本文方法所得到的格網在不同剖分層次,角度變形和面積變形均保持為最小,具有良好的幾何性質。通過數據集成實例驗證,本文方法所得到的格網整體分布均勻,格網單元面積形變小、排列緊致,不存在明顯的形狀突變,在空間數據建模與集成方面具有明顯的優勢。下一步的工作是基于正二十面體菱形格網剖分模型探索菱形格網編碼模型以及空間數據的表達和管理方法。

[1] GOODCHILD F.Discrete Global Grids:Retrospect and Prospect[J].Geography and GeoInformation Science, 2012,28(1):1-6.

[2] VINCE A.Indexing the Aperture 3 Hexagonal Discrete Global Grid[J].Journal of Visual Communication and Image Representation,2006,17(6):1227-1236.

[3] LAMBRECH TS J,COMBLEN R,LEGAT V,et al.Multiscale Mesh Generation on the Sphere[J].Ocean Dynamics,2008,58(3):461-473.

[4] STANIFORTH A,THUBURN J.Horizontal Grids for Global Weather and Climate Prediction Models:A Review [J].Quarterly Journal of the Royal Meteorological Society, 2012,138:1-26.

[5] ZH AO Xuesheng,BAI Jianjun,WANG Zhipeng.An Adaptive Visualized Model of the Global Terrain Based onQTM[J].Acta Geodaetica et Cartographica Sinica,2007,36 (3):316-320.(趙學勝,白建軍,王志鵬.基于QTM的全球地形自適應可視化模型[J].測繪學報,2007,36(3): 316-320.)

[6] BAI Jianjun,SUN Wenbin,ZHAO Xuesheng.Character Analysis and Hierarchical Partition of WGS-84 Ellipsoidal Facet Based on QTM[J].Acta Geodaetica et Cartographica Sinica,2011,40(2):243-248.(白建軍,孫文彬,趙學勝.基于QTM的WGS-84橢球面層次剖分及其特點分析[J].測繪學報,2011,40(2):243-248.)

[7] ZHANG Yumei,CHEN Weihua,NIE Hongshan,et al.Study on Sphere Rhombus Grid Recursive Subdivision[J].Geography and Geo-Information Science,2010,26(6):34-37.(張玉梅,陳維華,聶洪山.球面菱形網格遞歸剖分方法研究[J].地理與地理信息科學,2010,26(6):34-37.)

[8] ZHAO Xuesheng,BAI Jianjun.Hierarchical Model of Global Discrete Grids Based on Diamonds[J].Journal of China University of Mining and Technology,2007,36(3):397-401.(趙學勝,白建軍.基于菱形塊的全球離散格網層次建模[J].中國礦業大學學報,2007,36(3):397-401.)

[9] BEN Jin,TONG Xiaochong,YUAN Chaopeng.Indexing Schema of the Aperture 4 Hexagonal Discrete Global Grid System[J].Acta Geodaetica et Cartographica Sinica, 2011,40(6):785-789.(賁進,童曉沖,元朝鵬.孔徑為4的全球六邊形格網系統索引方法[J].測繪學報,2011,40 (6):785-789.)

[10] BAI Jianjun.Location Coding and Indexing Aperture 4 Hexagonal Discrete Global Grid Based on Octahedron[J].Journal of Remote Sensing,2011,15(6):1131-1146.(白建軍.基于正八面體的四孔六邊形球面格網編碼及索引[J].遙感學報,2011,15(6):1131-1146.)

[11] TODD R,PETERSEN M,ROBERT L.A Multi-resolution Approach to Global Ocean Modeling[J].Ocean Modelling, 2013,69:211-232.

[12] WHITE D.Comparing Area and Shape Distortion on Polyhedral-based Recursive Partitions of the Sphere[J].International Journal of Geographical Information Science, 1998,12(8):805-807.

[13] SATOH M,MATSUNO T,TOMITA H,et al.Nonhydrostatic Icosahedral Atmospheric Model(NICAM)for Global Cloud Resolving Simulations[J].Journal of Computational Physics,2008,227:3486-3514.

[14] YUAN Wen,MA Ainai,GUAN Xiaojing.A New Projection for Spherical Triangle:Equal Angle Ratio Projection (EARP)[J].Acta Geodaetica et Cartographica Sinica, 2005,34(1):78-84.(袁文,馬藹乃,管曉靜.一種新的球面三角投影:等角比投影(EARP)[J].測繪學報,2005,34 (1):78-84.)

[15] BEN Jin,TONG Xiaochong,ZHANG Yongsheng,et al.Research on Generating Algorithm and Software Model of Discrete Global Grid Systems[J].Acta Geodaetica et Cartographica Sinica,2006,36(2):187-191.(賁進,童曉沖,張永生,等.球面等積網格系統生成算法與軟件模型研究[J].測繪學報,2007,36(2):187-191.)

[16] TONG Xiaochong,BEN Jin,WANG Yin.A New Effective Hexagonal Discrete Global Grid System:Hexagonal Quad Balanced Structure[C]∥Proceedings of 2010 18th International Conference on Geoinformatics.Beijin:[s.n.],2010:1-6.

[17] FEKETE G,TREINISH L.Sphere Quadtrees:A New Data Structure to Support the Visualization of Spherically Distributed Data[C]∥Proceedings of SPIE on Extracting Meaning from Complex Data:Processing,Display, Interaction.[S.l.]:SPIE,1990:242-253.

[18] SAHR K,WHITE D,KIMERLING J.Geodesic Discrete Global Grid Systems[J].Cartography and Geographic Information Science,2003,30(2):121-134.

[19] WHITE D.Global Grids from Recursive Diamond Subdivisions of the Surface of an Octahedron or Icosahedron[J].Environmental Monitoring and Assessment,2000,64 (1):93-103.

[20] MING Tao,YUAN Wen,PENG Guagnxiong,et al.The Study on Error Analysis of Discretization Area in Discrete Global Grid System[J].Procedia Environmental Science, 2011,10(B):1122-1128.

(責任編輯:宋啟凡)

Diamond Discrete Grid Subdivision Method for Spherical Surface with lcosahedron

ZHOU Liangchen1,2,SHENG Yehua1,2,LlN Bingxian1,2,LüGuonian1,2,ZHAO Zhipeng1,2
1.Key Laboratory of Virtual Geographic Environment for the Ministry of Education,Nanjing Normal University, Nanjing 210023,China;2.Jiangsu Center for Collaborative lnnovation in Geographical lnformation Resource Development and Application,Nanjing 210023,China

Diamond discrete grids of spherical surface can meet the requirements of multi-scale spatial data modeling and analysis due to the advantages of simple structure and uniform directionality.This paper presents a new diamond discrete grid subdivision method for spherical surface with icosahedron.Firstly,the location of the icosahedron’s vertexes on the spherical surface are determined according to the earth’s major geographical features.Secondly,to solve the problem of grid convergence in poles, transformation relation between the latitude and longitude coordinates and 3D Cartesian coordinates are created.Then,a recursive subdivision approach of spherical diamond discrete grids is provided,and the geometric deformation of subdivision result is analyzed in terms of cell area and cell axial ratio.Finally, the feasibility of the method proposed is verified.The experiment results show that our subdivision method can construct spherical diamond discrete grids with icosahedron effectively,of which the general geometric characteristics are much better than that of the diamond discrete grids with octahedron.

discrete global gird;diamond;icosahedron;hierarchy recursively subdivision;geometric

ZHOU Liangchen(1979—),male,PhD, majors in discrete global grid systems.

SHEN Yehua

P208

A

1001-1595(2014)12-1293-07

國家自然科學基金(41301415;41001224;41271383);江蘇高校優勢學科建設工程資助項目

2014-01-03

周良辰(1979—),男,博士,主要研究方向為全球離散格網系統。

E-mail:zhoulch＠gmail.com

盛業華

E-mail:shengyehua＠njnu.edu.cn

ZHOU Liangchen,SHENG Yehua,LIN Bingxian,et al.Diamond Discrete Grid Subdivision Method for Spherical Surface with Icosahedron[J].Acta Geodaetica et Cartographica Sinica,2014,43(12):1293-1299.(周良辰,盛業華,林冰仙,等.球面菱形離散格網正二十面體剖分法[J].測繪學報,2014,43(12):1293-1299.)

10.13485/j.cnki.11-2089.2014.0192

distortion

修回日期:2014-09-03