利用“L-R”剖分及三角形聚合的近似等面積菱形格網(wǎng)建模

苑爭一，趙學勝，趙 靜

(1. 中國地震臺網(wǎng)中心，北京 100045； 2. 中國礦業(yè)大學(北京)地球科學與測繪工程學院，北京 100083)

全球離散格網(wǎng)(discrete global grid,DGG)是一種基于球面或橢球面的可以無限細分的地球擬合格網(wǎng)，該模型可以有效克服傳統(tǒng)地圖投影的束縛，從根本上解決傳統(tǒng)平面格網(wǎng)模型在全球多尺度空間數(shù)據(jù)管理上的數(shù)據(jù)斷裂、幾何變形和拓撲不一致等問題[1-2]，為解決越來越多全球一體化背景下的宏觀應(yīng)用問題提供思路[3-4]。球面菱形格網(wǎng)因其幾何結(jié)構(gòu)簡單，而且具有有一致的方向性、徑向?qū)ΨQ性及平移相和性，因而在空間操作和數(shù)據(jù)集成方面有著廣泛的應(yīng)用。文獻[5]基于球面菱形格網(wǎng)模型，發(fā)展了具有固定方向的“塊”層次編碼技術(shù)及鄰近搜索算法；文獻[4]實現(xiàn)了在球面菱形格網(wǎng)上集成和管理中國海潮波系統(tǒng)模型相關(guān)數(shù)據(jù)；文獻[6]提出了一種基于菱形格網(wǎng)的，用于層級間數(shù)據(jù)集中與分散的預(yù)測模型(Kalman過濾模型)。

格網(wǎng)等積性作為提高空間統(tǒng)計分析精度和空間采樣合理性的有效手段[7]，已經(jīng)成為球面格網(wǎng)模型構(gòu)建的研究熱點之一。目前的菱形剖分通常以球面三角剖分為基礎(chǔ)，通過合并相鄰的格網(wǎng)單元生成。傳統(tǒng)的等面積球面三角剖分模型大多基于等面積投影，該類模型需要復(fù)雜的迭代計算，而且投影邊界扭曲嚴重，為了方便球面點的定位，通常利用大圓弧線代替投影邊界，而替代后的模型會引起格網(wǎng)單元的面積變形，如基于Snyder等面積多面體投影的ISEA剖分模型[8]。另外，Song等[9]提出了一種基于小圓弧的等面積剖分模型，該模型克服了投影法邊界扭曲嚴重的缺陷，但其計算方式復(fù)雜，坐標轉(zhuǎn)換困難；Seong[10]采用等長度的經(jīng)線、緯線弧段構(gòu)建了等積球面四邊形格網(wǎng)，該方案破壞了格網(wǎng)單元間簡潔的臨近關(guān)系，增加了臨近搜索的復(fù)雜性；Holhos等[11]針對傳統(tǒng)投影方法計算復(fù)雜的缺陷，推導(dǎo)出了一種新的等面積投影公式，該公式具有對稱性，計算簡便，但是仍然沒有解決傳統(tǒng)投影方法邊界扭曲嚴重的缺陷，而且球面點的定位較為復(fù)雜。文獻[7]在Snyder等積投影的基礎(chǔ)上，沿著面積誤差率較大的3條軸線(由八面體中心到3個頂點)對相鄰的球面三角形格網(wǎng)進行了兩兩合并，構(gòu)建了一種近似等面積的球面菱形格網(wǎng)模型，該方法有效減小了格網(wǎng)單元的面積變形，但由于格網(wǎng)合并的方向不同，增加了模型構(gòu)建的復(fù)雜性。文獻[12]提出了一種計算簡便，同時又具有層次嵌套性、格網(wǎng)統(tǒng)一性和全球無縫連續(xù)性等特點的近似等面積QTM剖分方法——緯線環(huán)法(latitude-ring (L-R) method)。

通過分析緯線環(huán)法QTM模型的格網(wǎng)單元面積變形規(guī)律，發(fā)現(xiàn)該剖分模型的格網(wǎng)單元面積變形具有鋸齒分布特征，即相鄰的上下兩個格網(wǎng)單元總是一大一小交替分布，而且變形方向一致。因此，本文以該模型為基礎(chǔ)，通過合并上下兩個相鄰的三角形單元，得到一種新的面積更加接近的球面菱形剖分模型，最后設(shè)計相關(guān)變形試驗，分析該模型的格網(wǎng)單元面積變形特征及分布規(guī)律。

1 基于緯線環(huán)法的等面積QTM模型

文獻[12]提出了一種基于緯線法QTM剖分及變經(jīng)緯度剖分的近似等面積剖分方案——緯線環(huán)法QTM剖分，其基本思路為：通過相鄰兩條緯線所圍成的緯線環(huán)的理想面積(該緯線環(huán)所包含的三角形個數(shù)與理想剖分單元面積的乘積)，構(gòu)建方程依次求解每條分割緯線的位置，進而通過平分經(jīng)度的方法確定格網(wǎng)間的經(jīng)度間隔，完成對整個球面的QTM剖分(具體過程詳見文獻[12—13])。

按照傳統(tǒng)的緯線法和上述緯線環(huán)法分別對球面進行了QTM剖分，同時計算了格網(wǎng)單元的面積變形率。以1/8球面為例，將球面剖分6層，把格網(wǎng)中心點投影到球面內(nèi)接正八面體的一個面上，用Z軸表示格網(wǎng)單元的面積變形率，得到了圖1所示的兩種剖分模型的格網(wǎng)單元面積變形率的對比圖。從中可以看出，“L-R”QTM剖分(圖1①部分)除了極點附近外的其他大部分區(qū)域都非常接近于理想QTM(圖1③部分，其格網(wǎng)單元無面積變形)，表明該剖分有效控制了格網(wǎng)單元的面積變形，能夠達到近似等面積剖分的效果。

圖1 兩種剖分方法的面積變形對比(第6層)

為了進一步研究該模型面積變形的空間分布情況，分別計算每個格網(wǎng)單元的面積，根據(jù)其大小賦予不同的灰度值，借助DirectX繪制在球面上得到了第4—6層的QTM格網(wǎng)單元面積分布灰度圖和面積分布直方圖，如圖2所示。從中可以看出：①格網(wǎng)單元面積變形分布呈現(xiàn)出規(guī)則的鋸齒狀特性，即相鄰兩個格網(wǎng)單元的面積總是一大一小的交替分布，該變形特征便于沿同一方向進行三角形聚合；②盡管由大圓弧線代替投影曲線產(chǎn)生了面積變形，但觀察面積變形分布直方圖可以發(fā)現(xiàn)，格網(wǎng)單元面積分布較為集中，即格網(wǎng)單元的面積近似性較好。綜上所述，緯線環(huán)法QTM剖分不但能夠有效控制剖分單元的面積變形，而且其不需要復(fù)雜的投影及迭代計算，可以大大提高剖分的效率。因此，本研究擬在此剖分方法的基礎(chǔ)上，通過三角形聚合的方法實現(xiàn)對球面的菱形剖分。

2 球面菱形格網(wǎng)剖分及變形分析

2.1 剖分過程

菱形格網(wǎng)生成的詳細過程如圖3所示。首先按上一節(jié)中的緯線環(huán)法對八分之一球面進行近似等面積剖分；進而通過旋轉(zhuǎn)實現(xiàn)對整個球面的QTM剖分；將球面沿0°、90°、180°、270°經(jīng)線展開成4個基礎(chǔ)菱形塊；最后在每個基礎(chǔ)菱形塊中合并相鄰的上下兩個三角形格網(wǎng)單元，從而得到覆蓋整個球面的菱形剖分。

2.2 格網(wǎng)單元面積變形分析

按上述方法構(gòu)建了球面菱形剖分模型，計算出了每個菱形格網(wǎng)面積變形率(相對于理想剖分)。為了便于比較聚合前后格網(wǎng)單元面積變形情況的差異，同時計算了基礎(chǔ)QTM剖分的格網(wǎng)單元面積變形率。本節(jié)將從格網(wǎng)單元面積變形分布區(qū)段及格網(wǎng)單元面積變形位置分布兩個方面出發(fā)，對比分析剖分模型聚合前后的面積變形情況。

2.2.1 剖分單元面積變形分布區(qū)段研究

根據(jù)式(1)計算格網(wǎng)單元相對于理想剖分單元的面積變形率(Area_dis_rate)，式中：area為剖分單元面積；areanorm為理想剖分單元面積。

圖2 格網(wǎng)單元面積分布灰度圖、面積分布直方圖

圖3 菱形格網(wǎng)生成過程

(1)

基礎(chǔ)QTM剖分的格網(wǎng)單元面積變形率分布區(qū)段見表1，從中可以看出：基礎(chǔ)三角剖分模型在前5層中面積變形率的分布無明顯的規(guī)律，第6層以后格網(wǎng)單元面積變形率分布近似呈正態(tài)分布(如圖4所示)。隨著剖分層次的增加，模型的面積變形情況繼續(xù)改善，第8層以后剖分單元面積變形率被控制在±1%以內(nèi)的格網(wǎng)單元達到99%以上。

表1 基礎(chǔ)三角剖分單元面積變形率分布區(qū)段 (%)

圖4 緯線環(huán)法QTM剖分面積變形的統(tǒng)計特征

聚合菱形剖分單元面積變形分布區(qū)段見表2，可以看出：聚合菱形剖分模型在前3層中面積變形率的分布無明顯規(guī)律，第4層以后90%以上的格網(wǎng)單元面積變形率分布在±1%以內(nèi)；將面積變形區(qū)段(-1%，1%]進一步細分，并統(tǒng)計了剖分單元的面積變形率在第6～10層的分布情況，從中可以看出：格網(wǎng)單元面積變形情況繼續(xù)得到改善，第8層以后面積變形率被控制在±0.1%內(nèi)的格網(wǎng)達到99%以上。

表2 聚合菱形剖分單元面積變形率分布區(qū)段 (%)

2.2.2 剖分單元面積變形位置分布研究

以八分之一球面和對應(yīng)的聚合后的一個菱形塊為例，將球面剖分4、5、6層，給每個面積變形率分布區(qū)段對應(yīng)一個特定的灰度值，借助DirectX繪制了剖分模型面積變形率的位置分布灰度圖(如圖5、圖6所示)。圖中淺色區(qū)域的格網(wǎng)單元面積變形率較小，深色區(qū)域較大。

圖5 基礎(chǔ)QTM剖分格網(wǎng)單元面積變形率的位置分布

圖6 聚合菱形剖分單元面積變形率的位置分布

從圖6及其層次變化可以看出：①隨著剖分層次的增加，三角剖分模型的面積變形情況逐漸得到改善，變形邊界(變形率=1%)逐漸向極點移動；②隨著剖分層次的增加，聚合菱形剖分模型的面積變形情況也不斷改善，第6層以后，中低緯度絕大部分區(qū)域的面積變形被控制在1%以內(nèi)；③對于相同層次的兩種剖分模型，聚合菱形剖分模型的面積變形率更小，分布也更加均勻。

3 結(jié) 語

本文利用緯線環(huán)法QTM剖分方法，構(gòu)建了近似等面積的球面四元三角剖分模型，并分析了該剖分模型的格網(wǎng)單元面積變形。在此基礎(chǔ)上通過三角形聚合，構(gòu)建了一種面積近似效果很好的球面菱形剖分模型。通過對比聚合前后的格網(wǎng)單元面積變形情況，得出了以下結(jié)論：

(1) 面積變形?。壕酆锨暗娜切纹史帜Ｐ?，其相鄰的上下兩個格網(wǎng)單元總是一大一小交替分布，通過將其聚合成菱形格網(wǎng)后，有效地減小了這一變形。相比同層次的三角剖分，聚合菱形剖分模型的格網(wǎng)單元面積變形率更小。

(2) 位置分布明確：不論是基礎(chǔ)三角剖分，還是聚合后的菱形剖分，其面積變形都會隨剖分層次的增加而不斷改善，中低緯度區(qū)域不斷接近等面積剖分。相比之下，聚合菱形剖分的改善效果更好，第6層以后，除了極點處的極少數(shù)格網(wǎng)單元，其他區(qū)域基本達到了近似等面積剖分的效果，能夠滿足絕大部分統(tǒng)計分析的應(yīng)用需求。

另外，由于該剖分方案不涉及復(fù)雜的投影變換和迭代計算，顯著加快了計算速度，從而有效提高了剖分的效率。

[1] 趙學勝，王磊，王洪斌，等.全球離散格網(wǎng)的建模方法及基本問題[J].地理與地理信息科學,2012, 28(1)：29-34.

[2] 趙學勝，侯妙樂，白建軍.全球離散格網(wǎng)的空間數(shù)字建模[M].北京：測繪出版社,2007.

[3] 邢華橋，侯妙樂，王磊，等.基于QTM的大范圍有源淹沒算法研究[J].測繪通報,2015(12)：46-49.

[4] 劉光鑫.基于菱形離散格網(wǎng)的中國海潮波系統(tǒng)數(shù)據(jù)集成方法[D].贛州：江西理工大學,2013.

[5] 趙學勝，白建軍.基于菱形塊的全球離散格網(wǎng)層次建模[J].中國礦業(yè)大學學報,2007, 36(3)：397-401.

[6] HUANG H, GRESSIE N. Multiscale Spatial Modeling[C]∥Proceeding of the Section on Statistics and Environment.[S.l.]：IEEE,1997.

[7] 孫文彬，周長江.一種近似等面積球面菱形格網(wǎng)的構(gòu)建方法[J].武漢大學學報(信息科學版),2016, 41(8)：1040-1045.

[8] SNYDER J P. An Equal-area Map Projection for Polyhedral Globes[J]. Cartographica: The International Journal for Geographic Information and Geovisualization,1992, 29(1)：10-21.

[9] SONG L, KIMERLING A J, SAHR K. Developing an Equal Area Global Grid by Small Circle Subdivision[C]∥Second International Conference on Discrete Global Grids. Santa Barbara: National Center for Geographic Information & Analysis, 2002.

[10] SEONG J C. Implementation of an Equal-area Gridding Method for Global-scale Image Archiving[J]. Photogram-metric Engineering& Remote Sensing,2003, 71(5)：623-627.

[12] 趙學勝，苑爭一，趙龍飛，等.一種改進的近似等面積QTM剖分模型[J].測繪學報,2016, 45(1)：112-118.

[13] 趙學勝，孫文彬，陳軍.基于QTM的全球離散格網(wǎng)變形分布及收斂分析[J].中國礦業(yè)大學學報,2005, 34(4)：438-442.

[14] GOODCHILD M F. Criteria for Evaluation of Global Grid Models for Environmental Monitoring and Analysis[J]. Handout from NCGIA Initiative,1994, 15：9-47.

[15] ZHOU L C, LIN B X, LV G N, et al. Distortion Distribution and Convergence Analysis of Spherical Diamond Discrete Grids[J]. ISPRS-International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences,2013, 1(2)：89-93.

[16] 孫文彬，趙學勝，高彥麗，等.球面似均勻格網(wǎng)的剖分方法及特征分析[J].地理與地理信息科學,2009, 25(1)：53-56.