基于球面QTM的大范圍有源淹沒算法研究

邢華橋，侯妙樂，王　磊，姜曉軼

(1. 中國礦業(yè)大學(北京)地球科學與測繪工程學院，北京 100083； 2. 北京建筑大學，

北京 100044； 3. 國家海洋信息中心，天津 300171)

Submergence Analysis Algorithm for Large Area Based on Spherical QTM

XING Huaqiao，HOU Miaole，WANG Lei，JIANG Xiaoyi

?

基于球面QTM的大范圍有源淹沒算法研究

邢華橋1，侯妙樂2，王磊1，姜曉軼3

(1. 中國礦業(yè)大學(北京)地球科學與測繪工程學院，北京 100083； 2. 北京建筑大學，

北京 100044； 3. 國家海洋信息中心，天津 300171)

Submergence Analysis Algorithm for Large Area Based on Spherical QTM

XING Huaqiao，HOU Miaole，WANG Lei，JIANG Xiaoyi

摘要：面對海平面上升后沿海區(qū)域大范圍淹沒模擬的應用需求，提出了一種基于球面QTM的大范圍有源淹沒算法，首先構(gòu)建了研究區(qū)域的QTM三角格網(wǎng)集，并利用海陸分離的思想確定了初始淹沒單元；然后基于QTM的12鄰近搜索設計了面向多目標三角格網(wǎng)的淹沒范圍擴張方法，將初始淹沒單元遞歸擴張得到最終的淹沒范圍；最后應用Visual C#語言及Microsoft DirectX三維圖形接口設計開發(fā)了相應的試驗原型系統(tǒng)，并利用天地圖高分辨率影像和SRTM地形數(shù)據(jù)對試驗區(qū)域進行了海平面上升后的大范圍淹沒模擬試驗。結(jié)果表明，該算法具有較好的實用性，可以有效地對海平面上升后的大區(qū)域范圍進行動態(tài)地淹沒模擬，對輔助防災減災具有一定意義。

關鍵詞：球面QTM；海平面上升；海陸分離；12鄰近搜索；有源淹沒

一、引言

隨著全球氣候變暖等極端天氣的頻發(fā)，海平面異常出現(xiàn)的次數(shù)明顯增多，異常值明顯增大，近30年來全球海平面上升速率有加大的趨勢[1]，而同期中國沿海海平面的上升速率高于全球的平均值[2-3]。海平面上升帶來的最為直接和嚴重的影響是沿海地區(qū)的淹沒問題。因此，準確、科學地預測海平面上升后的淹沒范圍，對海洋災害的預測預報、災害過程的動態(tài)模擬，以及促進海洋災害機理的研究、輔助防災減災等都具有重要意義。

近年來，國內(nèi)外研究學者對淹沒分析的研究取得了一定的進展，淹沒分析主要分為無源淹沒和有源淹沒[4-6]，無源淹沒將所有高程值低于給定水位的區(qū)域均記為淹沒區(qū)，有源淹沒除了考慮上述情況外，還考慮到了“流通”的淹沒情形，即只有與淹沒區(qū)域相連通且高程值低于給定水位的區(qū)域才被視為淹沒區(qū)。目前的淹沒分析算法多以研究區(qū)域的DEM數(shù)據(jù)為基礎[7]，利用種子蔓延算法進行淹沒區(qū)域的模擬。如張東華等提出一種基于DEM的洪水有源淹沒算法，在GIS技術的基礎上應用數(shù)字高程模型(DEM)的格網(wǎng)模型進行洪水淹沒分析[8]；丁志雄等提出基于格網(wǎng)模型的洪水淹沒分析算法，在給定洪水水位和洪量兩種條件下，分別基于三角形格網(wǎng)模型和任意多邊形格網(wǎng)模型進行了水淹模擬[9]；聶漢江等提出了基于庫群聯(lián)合調(diào)度和DEM的區(qū)域洪水淹沒范圍模擬方法，以聯(lián)合調(diào)度模型得到的累積成災水量為淹沒水量，以研究區(qū)域的DEM為數(shù)據(jù)基礎，運用體積法對一定體量的洪水淹沒進行了模擬[10]。但這些研究大都是針對小范圍區(qū)域，采用傳統(tǒng)的平面數(shù)據(jù)模型進行的洪水淹沒分析，而海平面上升帶來的淹沒通常范圍較大，傳統(tǒng)的淹沒分析方法并不能滿足這種大范圍淹沒分析的需求。

針對這一問題，本文基于球面四元三角網(wǎng)(quaternary triangular mesh，QTM)的12鄰近搜索，提出一種面向大范圍區(qū)域的有源淹沒算法，以實現(xiàn)海平面上升后研究區(qū)域的淹沒動態(tài)模擬。

二、QTM理論知識

QTM是由內(nèi)接于球體的正八面體經(jīng)多次遞歸剖分而成的球面擬合三角格網(wǎng)集[11]，具有符合球面的特性和本質(zhì)是柵格的特點，可以用QTM模擬地球表面；另外，QTM形狀和面積的近似相等能保證格網(wǎng)之間具有簡潔的鄰近關系，便于格網(wǎng)間的鄰近查詢、分析等處理[12]。

1. QTM剖分與編碼

在QTM初始化剖分時，內(nèi)接于球體的正八面體頂點占據(jù)球面主要點(包括兩極)，而邊的投影與赤道、主子午線，以及90°、180°、270°子午線重合，如圖1所示。在對初始的QTM進一步細化剖分時，用大弧平分法找到三角形邊的中點，通過連接它們將球面三角形分成4個小的三角形，以此類推進行遞歸，對整個球面進行近似均勻的剖分[13]，形成了全球的QTM離散格網(wǎng)集。

圖1　QTM初始剖分單元

QTM編碼不僅暗含著格網(wǎng)的空間位置，而且表示了格網(wǎng)的剖分層次，因此QTM編碼是進行鄰近搜索及連通性分析的基礎。本文采用的編碼方案為固定方向編碼[12]，與其他編碼方案相比，該方案的所有格網(wǎng)都具有固定方向，更利于鄰近搜索。如圖2所示，編碼的第1位為八分體號，從第2位開始，0表示中間三角格網(wǎng)，1表示頂(底)三角格網(wǎng)，2表示左三角格網(wǎng)，3表示右三角格網(wǎng)，即可以用一個四進制的Morton碼來標示。因此，QTM的編碼M由八分碼M0加上Morton碼組成，即M=M0q1q2…qi。

圖2　固定方向編碼方案

2. QTM的鄰近關系

鄰近關系是GIS空間分析中必不可少的一類空間關系，是球面實體擴張和空間索引的基礎[14-15]，在實際應用中(如水淹分析、最短路徑分析等)具有重要的意義。如圖3所示，QTM三角格網(wǎng)的鄰近三角格網(wǎng)分為具有公共邊的邊鄰近(edge-adjacent)三角格網(wǎng)(編號為E)和具有公共點的角鄰近(vertex-adjacent)三角格網(wǎng)(編號為V)。將球面按QTM剖分后，每個球面三角格網(wǎng)都有3個邊鄰近三角格網(wǎng)，位于南、北極的8個頂角三角形和赤道上經(jīng)度為0、90°、180°(-180°)、-90°的16個頂角三角形均具有7個角鄰近三角形(如圖3(a)所示)，其他位置的非頂角三角形均具有9個角鄰近三角形(如圖3(b)所示)。根據(jù)QTM的鄰近關系，其鄰近搜索分為三鄰近搜索和12鄰近搜索，由于海水淹沒不具有方向性，即海水可能向所有方向蔓延，因此，利用12鄰近搜索更符合實際情況。本文主要考慮QTM格網(wǎng)的12鄰近搜索[16]。

圖3　QTM的鄰近三角格網(wǎng)

三、基于QTM 12鄰近搜索的有源淹沒算法

1. 總體思路

首先構(gòu)建研究區(qū)域的一定層次的QTM格網(wǎng)集，基于海陸分層的思想確定初始淹沒單元，并以此為基礎進行12鄰近搜索獲取其鄰近格網(wǎng)，如果鄰近格網(wǎng)的高程低于給定的上升水位，即HQTM

圖4　基于QTM 12鄰近搜索的有源淹沒算法流程

2. 基于海陸分離思想的初始淹沒單元判定

確定初始淹沒區(qū)域是海平面上升分析與模擬的首要問題，海平面上升后的淹沒過程是從陸地與海洋相接處開始的。因此，本文利用海陸分離的思想，將研究區(qū)域的QTM按照3個頂點的高程劃分為陸地和海洋，并確定陸地邊緣格網(wǎng)。具體方法為：對于每個格網(wǎng)，如果其3個頂點均位于海洋(頂點的高程值均為-9999 m，如圖5(a)所示)則為海洋格網(wǎng)；如果其3個頂點均位于陸地(頂點的高程值均大于-9999 m，如圖5(b)所示)則為陸地格網(wǎng)；如果其3個頂點并不完全位于陸地或海洋，那么將該格網(wǎng)視為陸地邊緣格網(wǎng)(如圖5(c)所示)，所有陸地邊緣格網(wǎng)構(gòu)成的集合即為初始淹沒單元。

圖5　QTM三角格網(wǎng)陸地與海洋的判別

3. 面向多目標三角格網(wǎng)的淹沒范圍擴張

由淹沒單元通過鄰近搜索擴張得到其鄰近格網(wǎng)時，淹沒單元往往不是單獨的一個QTM格網(wǎng)，而是多個連續(xù)的格網(wǎng)，這就對鄰近搜索提出了特殊的需求，即需要對多個三角格網(wǎng)同時進行擴張，而面向單目標三角格網(wǎng)的擴張方法會造成格網(wǎng)的重復擴張，從而增加了系統(tǒng)的計算量。本文在淹沒范圍的擴展過程中，剔除擴張得到的鄰近格網(wǎng)中的重復格網(wǎng)，同時為進一步減小工作量，每一次擴張只將上一次擴張得到的被淹沒格網(wǎng)向外擴張，而不再處理已經(jīng)判斷過的格網(wǎng)。面向多目標三角格網(wǎng)的淹沒范圍擴張過程如圖6所示。該擴張方法的過程具體如下：

1) 利用12鄰近搜索算法依次計算初始淹沒單元S中所有格網(wǎng)的鄰近格網(wǎng)，并從中剔除重復格網(wǎng)，將不重復的鄰近格網(wǎng)集合記作G。

2) 利用雙線性插值算法計算格網(wǎng)集G中每個格網(wǎng)中心點的高程。

3) 將格網(wǎng)集G中每個格網(wǎng)的高程與當前水位數(shù)據(jù)進行對比，高程低于當前水位的格網(wǎng)為被淹沒的格網(wǎng)，記作G′。

4) 以G′代替步驟1)中的初始淹沒單元S，執(zhí)行步驟1)—步驟3)，直至G中所有格網(wǎng)的高程均高于當前水位(即G′為空)或整個試驗區(qū)的格網(wǎng)搜索完畢。

圖6　面向多目標三角格網(wǎng)的淹沒范圍擴張過程

四、試驗與分析

為了對算法進行驗證，本文以C# 3.0為開發(fā)語言，DirectX為三維開發(fā)包，從底層開發(fā)了試驗原型系統(tǒng)。試驗硬件環(huán)境為Pentium(R) Dual-core CPU T4200 @2.00 GHz，2.00 GB內(nèi)存，采用SRTM 地形數(shù)據(jù)作為DEM數(shù)據(jù)源，天地圖高分辨率影像為影像數(shù)據(jù)源，在給定淹沒水位的情形下對海平面上升后研究區(qū)域進行動態(tài)淹沒模擬。試驗中首先構(gòu)建了一定層次的QTM三角格網(wǎng)集作為研究區(qū)域的模擬圖層，然后利用海陸分離思想提取試驗區(qū)域的初始淹沒單元，如圖7所示。

圖7　試驗區(qū)域的初始淹沒單元

進而利用本文提出的算法在給定淹沒水位下獲取所有被淹沒的格網(wǎng)存入淹沒格網(wǎng)集，并將淹沒格網(wǎng)集中所有格網(wǎng)進行著色和渲染，然后與無源淹沒算法的淹沒效果進行了對比。圖8為淹沒水位為5 m的試驗區(qū)域淹沒范圍的模擬結(jié)果，圖9為利用無緣淹沒算法得到的試驗區(qū)域淹沒范圍。通過淹沒模擬對比圖可以看到，利用本文算法得到的淹沒區(qū)域具有連通性，更加符合淹沒的實際情況。

為了進一步比較本文算法與無源淹沒算法的區(qū)別，并驗證本文算法在大區(qū)域范圍進行淹沒模擬的有效性，在相同的試驗硬件環(huán)境下，應用不同的淹沒水位分別在不同面積的大范圍試驗區(qū)域進行有源淹沒和無源淹沒模擬，并比較了不同淹沒算法的淹沒面積，見表1。試驗結(jié)果表明，隨著淹沒水位的升高，淹沒范圍逐漸增大；在利用相同淹沒水位對同一的試驗區(qū)域進行模擬時，本文算法得到的淹沒范圍相對于無源淹沒算法較小，這是由于本文算法考慮到了淹沒區(qū)域的連通性，只有高程低于給定水位且與初始淹沒單元相連通的區(qū)域才被視為最終的淹沒區(qū)域。此外，與傳統(tǒng)的有源淹沒算法相比，本文算法所應用的試驗區(qū)域達到了數(shù)萬平方千米的級別，可以有效地應用于更大范圍的淹沒分析中， 并且能將淹沒結(jié)果動態(tài)地在三維球面上進行模擬和展示。

圖8　本文算法得到的淹沒模擬圖

圖9　無源淹沒算法得到的淹沒模擬圖

km2

五、結(jié)論與展望

本文針對海平面上升后的沿海區(qū)域淹沒模擬問題，提出了基于球面QTM 12鄰近搜索的大范圍有源淹沒算法，并設計開發(fā)了相應的原型試驗系統(tǒng)。試驗證明該方法可以用于球面大范圍海平面上升后的水淹動態(tài)模擬，能夠為防災減災提供決策依據(jù)。但該算法也有一定的局限性，海平面上升是個漸變的過程，其年變化范圍在毫米量級，而本文假定的海平面上升水位實際是疊加了平均潮位和極值潮位的最終上升水位。如何根據(jù)真實的海平面上升高度模擬淹沒范圍，并對海平面上升后給淹沒區(qū)域帶來的社會、經(jīng)濟影響等進行分析將是下一步的研究重點。

參考文獻：

[1]李響，劉克修， 董軍興，等. 中國沿海海平面上升風險區(qū)劃方法研究[J].海洋預報， 2014，31(2)：41-49.

[2]周子鑫.我國海平面上升研究進展與前瞻[J]. 海洋地質(zhì)動態(tài)，2008， 24(10)：14-18.

[3]王慧，劉克修，范文靜，等. 2012年中國沿海海平面上升顯著成因分析[J].海洋學報， 2014，36(5)：8-17.

[4]劉仁義，劉南. 基于GIS的復雜地形洪水淹沒區(qū)計算方法[J].地理學報，2001，56(1)：1-6.

[5]劉小生，陳英俊，黃玉生.基于GIS技術的洪水淹沒區(qū)確定[J]. 測繪科學，2007，32(5)：136-137.

[6]王建軍，吳北平，朱萬里，等.基于ArcGIS的風暴潮災害疏散系統(tǒng)的研究與開發(fā)[J].測繪通報，2013(9)：91-94.[15]國家測繪局．GB/T13990—921∶5000、1∶10000地形圖航空攝影測量內(nèi)業(yè)規(guī)范[S].北京：中國標準出版社，1992.

[7]韓麗蓉，張濤.基于RS/GIS的大南川水庫潰壩災害模擬方法研究[J].測繪通報，2014(1)：97-99.

[8]張東華，劉榮，張詠新，等.一種基于DEM的洪水有源淹沒算法的設計與實現(xiàn)[J].東華理工大學學報:自然科學版，2009，32(3)：181-184.

[9]丁志雄，李紀人，李琳.基于GIS 格網(wǎng)模型的洪水淹沒分析方法[J].水利學報，2004(6)：56-60.

[10]聶漢江，劉彬，何立新，等. 基于庫群聯(lián)合調(diào)度和DEM的區(qū)域洪水淹沒范圍模擬[J].水電能源科學，2014，32(5)：43-46.

[11]DUTTON G. Encoding and Handling Geospatial Data with Hierarchical Triangular Meshes [C]∥Proceeding of 7th International Symposium on Spatial Data Handling. Netherlands: Delft University of Technology, 1996. 34-43.

[12]SAHR K, WHITE D，KIMERLING A. Geodesic Discrete Global Grid Systems[J]. Cartography & Geographical Information Science, 2003,30(2): 121-134.

[13]LEE M，SAMET H. Navigating through Triangle Meshes Implemented as Linear Quadtree[J]. ACM Transactions on Graphics, 2000,19(2): 79-121.

[14]陳軍，趙仁亮.GIS 空間關系的基本問題與研究進展[J]. 測繪學報，1999，28(2)：96-102.

[15]侯妙樂，邢華橋，趙學勝，等.球面四元三角網(wǎng)的復雜拓撲關系計算[J].武漢大學學報:信息科學版，2012，37(4)：468-471.

[16]趙學勝，侯妙樂，白建軍.全球離散格網(wǎng)的空間數(shù)字建模[M].北京：測繪出版社，2007：40-46.引文格式： 顏金彪，胡最，禹信. CORS與靜態(tài)相對定位技術下的坐標轉(zhuǎn)換研究[J].測繪通報，2015(12)：54-56.DOI:10.13474/j.cnki.11-2246.2015.377

作者簡介：邢華橋(1988—)，男，博士生，研究方向為全球離散格網(wǎng)、服務計算。E-mail：xinghuaqiao@126.com

基金項目：國家自然科學基金(41171304；40701152；41171306)

收稿日期：2015-02-26； 修回日期： 2015-10-27

中圖分類號：P208

文獻標識碼：B

文章編號：0494-0911(2015)12-0046-04