汪明

(重慶市勘測院，重慶 400020)

1 引言

在解決大規(guī)模地形數(shù)據(jù)三維可視化和實時繪制的問題中，國內外許多學者都進行了廣泛、深入的研究，主要集中于地形和紋理數(shù)據(jù)的分頁管理、紋理數(shù)據(jù)和地形數(shù)據(jù)的LOD控制、可見域的裁剪等多個方面。此外，隨著圖形處理器(Graphic Processing Unit，GPU)性能的提高，出現(xiàn)了基于硬件GPU的粗粒度的離散層次細節(jié)方法［1］。

本文擬研究基于OpenGL的大規(guī)模地形三維地形構建關鍵技術。由于高分辨率的DTM或DEM的數(shù)據(jù)集龐大，若直接使用原始數(shù)據(jù)進行可視化，則需要繪制的三角形或四邊形將達到上百萬個。因此，為了防止顯示滯后，緩解內存和顯卡緩存的壓力，需要考慮以下幾個方面:對復雜的地形做數(shù)據(jù)簡化處理;三維地形數(shù)據(jù)LOD構建等。

2 三維地形構建技術算法及實現(xiàn)

本文基于多分辨率高程數(shù)據(jù)文件的層次細節(jié)構建地形，主要思想是對大范圍的DEM數(shù)據(jù)先做預處理，即裁剪分層(見3.1)。由于文件按行列組織，因此某一塊的索引則為i*m+j(i為行號，j為列號，m為列數(shù))。我們只讀取和顯示在可視范圍內的文件塊，這解決了把DEM數(shù)據(jù)全部讀入內存而實際只需要顯示地形中很小一部分的矛盾。

2.1 高程數(shù)據(jù)組織方法

本文中三維地形的構建，將采用多層次文件的形式存儲與讀取數(shù)據(jù)，而Esri的GIRD數(shù)據(jù)卻在保存上稍顯復雜，為了簡化存儲與后期處理的方便，將采用自定義文件存儲方法，對數(shù)據(jù)進行多層次組織與保存。

如圖1所示，A為低層級數(shù)據(jù)L，它所存儲的當前層級下的數(shù)據(jù)相對于高層級下的數(shù)據(jù)B(級別L－1)來說，其存儲范圍為2×2個B層級數(shù)據(jù)，但是A的分辨率只有B的1/4。A和B均為所在層級下的一個文件，而DEM數(shù)據(jù)由很多個這樣的文件組成，而且不管是同一層級還是不同層級，每個文件的大小一樣(除去邊緣上不足塊文件柵格數(shù)的文件)。本文采用每塊文件柵格數(shù)為長寬均為256，文件大小為“長×寬×浮點數(shù)長度”，浮點數(shù)為4個字節(jié)，即文件存儲大小為256×256×4=262144Byte=256kb。所分層級總數(shù)是按照如下公式計算得出:

每一層級(level)文件中相對于源文件重采樣點的間隔數(shù)XYDiff為:

圖1 分級數(shù)據(jù)示意圖

在保存文件時，采用每一層級存放在單獨文件夾下面，每一層級文件夾下有多個文件塊，分別代表當前級別下256×256個柵格的數(shù)據(jù)。如圖1中A和B中的紅色框。每塊數(shù)據(jù)以二進制格式寫入文件，后綴名為.bk，文件命名方式為:“當前級別_所在行號_所在列號.bk”，如圖1中A級別為L，并在i行j列，因此在級別為L的文件夾下，有文件名為L_i_j.bk。在L同級文件夾中包含一個cfg.cfg的二進制文件，其中包含了DEM的頭文件信息，結構體如下:

如果要對所構建的地形進行補縫，則需要在生成數(shù)據(jù)時在每個塊文件中插入相鄰的最后一行(列)的數(shù)據(jù)，即數(shù)據(jù)需要有重疊部分，如圖2中D包含了B有重疊行，而和C有重疊列。

圖2 補縫數(shù)據(jù)分塊示意圖

2.2 數(shù)據(jù)裁剪及重采樣

數(shù)據(jù)的裁剪程序使用了GDAL函數(shù)庫中的重采樣方法:

該函數(shù)中，nXOff、nYOff表示讀取柵格的起始行列;nBufXSize、nBufYSize表示讀取柵格的行列數(shù)。

通過設定這4個參數(shù)，可以實現(xiàn)柵格數(shù)據(jù)的剪裁和重采樣。

2.3 三維地形構建算法

(1)分塊分級的地形構建算法

在地形數(shù)據(jù)分級存儲的基礎上，確定觀察點與觀察目標點后，可以根據(jù)視點(Object)到觀察(Target)點的距離計算出當前讀入數(shù)據(jù)所需要的文件級別以及當前的視域。為保證三維地形顯示與DEM大小無關，規(guī)定顯示的地形最多只能是觀察點及其周圍8個文件塊。這樣當讀入文件都在同一級別時，讀入內存的數(shù)據(jù)最多為256×9 kb，從而在內存方面保證了實現(xiàn)海量數(shù)據(jù)的實時快速地形構建的可行性。

在地形顯示時，為保證地形顯示的效果，讀取文件時，不一定只讀取一個級別的文件，而更多的是通過四叉樹與級別判定公式讀取當前顯示塊、當前顯示級別下的文件。如圖3中，設L(i，j)為觀察點塊文件，而周圍三個塊文件則是經(jīng)過一定計算會顯示的文件塊。

圖3 地形構建實時分級

首先計算出視點到L(i，j)的空間距離Lij，然后轉換成柵格個數(shù)RLij，通過公式(3)計算得到當前塊的級別(其中blockSize為分塊文件的柵格數(shù))。把L(i，j)分為相同的4塊，以相同的計算方法計算L(i，j)4塊到視點的距離，并計算層級，如果有一個層級比L(i，j)中心點的級別更高，就需要讀取更高級別的文件數(shù)據(jù)，如圖3中L(i，j)繼續(xù)分級后的文件塊為圖4所示。圖5為某一次繪制時，讀入的不同級別的文件塊。

圖4 各級數(shù)據(jù)文件示意圖

圖5 視圖中的地形數(shù)據(jù)

為了保存每個文件塊(block)文件的數(shù)據(jù)，我們設計了如下結構體:

該數(shù)據(jù)結構中，iX、iY分別表示當前塊在當前級別(iLevel)下的列、行號;fData表示保存在內存中的文件塊中的高層數(shù)據(jù);pBlocks表示文件塊劃分的內存子塊，這些子塊將通過自繪完成地形的可視化;Width、Height表示內存中pBlocks所包含的內存子塊的列、行數(shù)。

當需要獲得某點高程數(shù)據(jù)時，只需要知道當前點P在總的DEM柵格中的行列號(x，y)，將行列號減去當前文件塊的起始行列號(ox，oy)，可得到P在當前文件夾中的行列號(m，n)∶(m，n)=(x－ox，y－oy)，再由n*Width為索引得到柵格點的值。

(2)塊級別評價函數(shù)

(3)式定義了數(shù)據(jù)分塊級別判定函數(shù)，RLij表示視點到當前文件塊或內存子塊中心點的柵格個數(shù)，中心點的高程值取整個DEM的平均值，且level均為整數(shù)。我們定義最精細的層級為1，視點距離目標點的柵格數(shù)小于文件分塊的柵格數(shù)的兩倍。(3)式中l(wèi)og2(RLij)求取RLij的對數(shù)，－6表示－7+1，即為－log2(BlockSize)+1，文中的BlockSize為文件子塊所劃分的柵格數(shù)，即256。因此，判別公式實際上是(4)式。

(5)式表示一般的視距柵格數(shù)求取方法，由公式可以看出相對于當前文件來說的柵格視距RL只與層級level與采樣間隔數(shù)XYDiff相關，在裁剪數(shù)據(jù)生成后BlcokSize為常數(shù)。由于當前視圖讀取的文件層級為level，并且當前文件中采樣間隔數(shù)XYDiff由level確定(2)式，因此可確定RL只與BlockSize相關，且為一常數(shù)(6)式。

由于RL為一個常數(shù)BlockSize，就可以通過目標點坐標與RL計算出當前級別下視域范圍(當前文件級別下的柵格矩形)，就可以在計算出目標點所在的文件塊后(圖6中文件塊5)，根據(jù)RL就可以判定，只需要再讀入文件塊5周圍的8個文件塊，就可以覆蓋整個視域，達到只讀取和構建視域范圍內地形的目的。

圖6 視域中的文件塊

(3)四叉樹分塊過程

3 算法性能分析與討論

3.1 數(shù)據(jù)分塊分析

測試中所用GRID數(shù)據(jù)為列12020，行12618，柵格大小 22.5 ×22.5，計算面積為 76782.10725 km2，文件大小578.57 MB，文件裁剪耗時5 min。文件裁剪后分為7個級別的數(shù)據(jù)，即包含7個文件夾，共3166個文件，總的數(shù)據(jù)大小為771 Mb，占用空間774 Mb。數(shù)據(jù)大小增加了33.25%。第一級別的數(shù)據(jù)大小為578 Mb，占用空間580 Mb，每降低一個級別，數(shù)據(jù)大小為上一級的四分之一。

3.2 三維地形顯示性能分析

在內存消耗方面，相對于把所有的地形數(shù)據(jù)讀入內存中，本文中提出的算法只讀入少量數(shù)據(jù)，僅僅是讀取有限的數(shù)十個文件塊?？梢詫λ惴ㄖ袛?shù)據(jù)讀入大小做一定假設，每次讀入的數(shù)據(jù)最多不超過40個小文件塊，即不超過10 Mb?？偟膩碚f，讀入內存的數(shù)據(jù)趨于常數(shù)，不會因為地形數(shù)據(jù)變大而變大。

在實時構建的地形中，同一時間所顯示的頂點數(shù)在256×256=65536到256×256×40=2621400之間。在具體實際中也可以將分開的柵格數(shù)降為128，可以減少場景中繪制的頂點數(shù)。在實時顯示時，場景動畫的幀頻為十幀以上。場景中繪制頂點數(shù)，與讀入的地形數(shù)據(jù)相關，會在一定范圍內變動。

3.3 地形縫隙修補

由于在地形構建中有不同級別的數(shù)據(jù)顯示，在不同級別數(shù)據(jù)的臨邊上，細節(jié)較豐富的塊上必定多出部分點，如圖7中的兩個圓圈。圓圈處不能保證左右兩邊的面在同一個面上，因此就需要對這些點和面做單獨處理。

圖7 地形縫隙示意圖

對于這種處于邊緣，并且相鄰子塊處于不同級別的子塊來說，需要用三角形繪制出采樣間隔造成的不重合的點，如圖8。但是這種情況只有在相鄰級別不超過1的情況下才能滿足條件，因此，需要在前期判斷和處理，使相鄰級別不超過1。

圖8 地形補縫原理示意圖

在每個64×64的子塊自繪時，進行對當前塊與其下部塊與左邊塊進行補縫，那就需要知道每個子塊的級別以及子塊所在的文件塊，以便讀取高程數(shù)據(jù)。我們把要顯示的文件塊用鏈表存儲，并為存儲在文件子塊結構體中的子塊提供相關的信息，比如，子塊所在文件在鏈表中的索引，再通過索引獲取結構體中的信息，這樣就可以獲得所需要的信息。

4 結論

本文提出的基于多級別柵格數(shù)據(jù)的多層次LOD劃分方法，對原有DEM數(shù)據(jù)文件做層級劃分與重采樣，生成具有層次和不同分辨率的小文件塊。在三維地形實時構建時，計算出視圖中要顯示的文件與級別，實現(xiàn)了部分地形數(shù)據(jù)的讀取，使內存占用大大減少。在顯示時，對內存中數(shù)據(jù)進一步做LOD分級，讓顯示的面片更少。

與已有的三維地形構建方法相比，本文中算法實時，由于只對每個塊進行層級判別，因此可視化過程中的計算量少，可視化過程中幾乎不受計算過程的影響，級別判斷時間大概為0.0083 s左右;而且相對于其他算法把全部數(shù)據(jù)讀入內存，本文算法在顯示中只需讀取可視范圍內的DEM數(shù)據(jù)與影像數(shù)據(jù)，在可視化過程中占用很少的內存，大概在3Mb～10 Mb之間。在這種算法基礎上，基于格網(wǎng)表達的DEM模型數(shù)據(jù)管理簡單易于實現(xiàn)交互操作。

由于在內存使用及計算量方面減少，本文提出的算法總的來說達到了預期目標，提出了一種新的三維地形構建方案，并優(yōu)化了現(xiàn)有的三維地形在某些方面的不足。

［1］呂?？?，周小平.實戰(zhàn)OpenGL三維可視化系統(tǒng)開發(fā)與源碼精解［M］.北京:電子工業(yè)出版社，2009:125－239.

［2］孫敏，薛勇，馬藹乃.基于格網(wǎng)劃分的大數(shù)據(jù)集DEM三維可視化［J］.計算機輔助設計與圖形學學報.2002，14(6):566－580.

［3］羅景馨，唐琎.基于改進四叉樹分割和結點存儲的LOD算法［J］.計算機工程.2009，35(20):202 －204.

［4］徐文學，江濤，姚兵.基于OpenGL的DEM三維仿真方法研究［J］.國土資源信息化.2008(1):7－11.

［5］黃超超，凌永順，呂相銀.ROAM動態(tài)地形渲染算法研究［J］.計算機仿真，2005，22(1):216 －219.

［6］周宏偉，楊平，陳琦.實時地形可視化ROAM算法的分塊改進［J］.測繪通報，2003，(8):61 －63.

［7］江流長.ROAM地形渲染算法的實現(xiàn)［J］.農業(yè)網(wǎng)絡信息，2006，(5):42 －44.

［8］涂超.ROAM算法原理及其應用研究［J］.遼寧工程技術大學學報，2003，22(2):176 －179.

［9］解志剛，馬秋禾，趙金萍.基于二叉樹結構的地形分塊實時漫游的研究［J］.海洋測繪，2004，24(5):35－38.

［10］劉修國，張劍波.基于DEM庫的地表模型實時簡化方法［J］.小型微型計算機系統(tǒng)，2004，25(2):280 －282.

［11］http://www.gdal.org［OL］