三維地形數據塔式結構模型與可視化方法研究

李鶴元，陳 剛，柯希林

（1.信息工程大學 地理空間信息學院，河南 鄭州 450052；2.地理信息工程國家重點實驗室，陜西 西安 710054）

三維地形數據塔式結構模型與可視化方法研究

李鶴元1,2，陳 剛1，柯希林2

（1.信息工程大學 地理空間信息學院，河南 鄭州 450052；2.地理信息工程國家重點實驗室，陜西 西安 710054）

針對在地形繪制中三維地形數據量巨大和繪制速度慢等問題，提出了一種基于GPU加速的大區域三維地形數據的塔式結構模型及實時可視化新方法。首先構建細節層次模型來減少大區域地形實時繪制的地形數據量，然后建立DEM數據和DOM數據的金字塔式數據模型來實現大容量三維數據的實時裝載和動態卸載，并采用CUDA編程實現了大區域三維地形可視化的GPU加速。

地形可視化；細節層次模型；計算機視覺

大區域三維地形數據的交互可視化中一直存在著一對矛盾，即地形數據大容量和可視化速度無延遲要求之間的矛盾。描述一個大區域地形通常需要10 G以上的數據，如果對每一塊數據進行詳細處理顯示會是一個很大的工程，而在地形可視化過程中，視域內的地形范圍只是整個地形的一部分，為提高地形繪制效率，將三維地形數據分塊存儲，在地形繪制時只處理視域內的部分，一直是解決矛盾的突破口。除此之外，由于各地形塊距視點的遠近不同，繪制時繪制分辨率也不必要求一致。因此在地形可視化過程中，根據視域的地形裁剪以及視相關的地形簡化能提高地形繪制速度，同時地形繪制精度也能得到保證。

海量三維地形數據可視化算法的要求為：①便于生成和管理海量三維數據；②便于根據場景窗口查找可視的地形區域；③便于快速裝載、卸載相關地形數據；④便于在實時地形繪制中快速計算地形塊的相應細節層次（LOD）等級；⑤便于對地形進行交互式高精度繪制。為達到上述要求，就需對三維地形數據進行預處理。通過地形數據的金字塔式管理，從原始數據派生得到多分辨率數據，由此可生成連續LOD四叉樹數據結構，從而簡單地實現三維地形的快速動態繪制。

1996年，Lindstrom等首次提出了實時連續LOD地形繪制方法[1-6]。該方法基于自下向上的地形分塊策略，第一次根據地形塊在屏幕上的投影精度誤差來確定地形塊繪制時的層次等級，在地形繪制中，增加那些誤差精度大于一定閾值的頂點，并通過在相鄰塊間采用不同的繪制等級來消除塊間的地形裂縫，在相鄰地形塊的接邊處共用相同的頂點。該算法通過預處理減少幀之間的頂點數從而實現幀之間的連貫性；通過地形的四叉樹分割來自動生成地形格網。1997年，Duchaineauy等使用ROAM方法繪制地形[7-14]。該方法繼承了Lindstrom算法中的一些思想，地形繪制速度更快。ROAM主要是基于三角形樹結構進行地形三角形的分割或合并，這一點和Lindstrom的算法中頂點二元化處理相似，都是采用自下向上的方式進行處理。ROAM算法使用了分割隊列和合并隊列2個棱形隊列的地形三角形網格，三角形網格根據精度誤差或其他誤差數值確定優先級，再進行排序。該算法不但嚴格控制了誤差精度，而且能控制場景繪制過程中三角形的個數，從而實現場景繪制中幀之間的連貫性和連續的幀速。此外還有其他一些關于實時地形三維可視化的算法[15-19]，Hoppe在地形繪制中運用連續視相關LOD控制技術[20]，Rottger在1998年提出高度場連續LOD模型實時生成算法，Jonathan在2000年提出高精度LOD模型的地形繪制。但這些方法依然存在繪制大數據量高精度地形圖時易出現時延較大、精度不夠高等問題，且有些是相似的，如場景的視景體裁剪、LOD方法的使用和地形裂縫的消除等。針對這些問題，本文結合LOD模型、多分辨率金字塔結構模型和基于GPU的加速處理，提出了一種三維地形數據的塔式管理和實時可視化新方法。

1 金字塔結構模型

1.1 多分辨率結構數據模型

金字塔數據模型是在地形數據和影像數據預處理基礎上建立的[21-23]，其策略是建立DEM和DOM數據的多分辨率數據模型。最頂層設為0級，分辨率最低；級數越高，數據分辨率越高；最高級數據的分辨率與原始數據的分辨率一致，如圖1。

圖1 金字塔結構數據模型

假定地形水平方向的范圍是（MinCol，MaxCol）、垂直方向范圍是（MinRow，MaxRow），地形的根節點是M×N塊，每塊大小是（2n+1）×（2n+1）， 通常為17×17 或33×33，相鄰地形節點間有條公共邊。對于DOM數據，紋理數據的大小通常強制設為2k×2k。如果當前地形層次等級是t，那么該層地形水平方向和垂直方向的地形分辨率為：

為了使地形節點大小固定，DEM和DOM數據的分辨率比率通常為1∶16或1∶32。為確保地形繪制效果，地形水平和垂直方向的分辨率比率采用同一個值。如圖2所示，在金字塔結構數據模型構建時，下一級地形塊的塊數是上一級地形塊的4倍，依次類推建立金字塔數據模型。

數據管理采用四叉樹結構來進行，高一級地形節點數是低一級的4倍，如圖3所示。

圖4是地形四叉樹結構層中的數據，根節點為0級，節點格網大小為n=2k+1。t級節點分辨率為：

圖2 金字塔模型地形分層分塊策略

圖3 地形數據管理中的四叉樹結構

其頂點位于上一級t-1頂點中心位置（圖4中的白點）。

圖4 金字塔模型中的四叉樹結構

1.2 多層紋理壓縮技術

在海量地形數據可視化時，DEM和DOM數據量很大，如果將原始數據一次性載入顯示列表，數據量肯定大大超過計算機內存和性能，從而降低了地形可視化的速度，影響交互可視化效果。而且用高分辨率影像數據可視化時，原始高分辨率與視覺習慣不相適應，因為視距的不同，有不同的視覺效果。事實上，離視點近的場景看起來細節詳細些，而離視點遠的場景細節幾乎是不體現的。因此對于紋理圖像數據，LOD技術非常有用。

1976年，Clark[24]提出了LOD模型的構建方法。當地形覆蓋的屏幕區域較小時，采用較粗略的模型描述要進行繪制的地形，并構建了可見面判定算法的幾何層次模型來對復雜地形進行快速繪制。在不影響地形觀察者的視覺效果的前提下，LOD技術通過逐次簡化地形場景的表面細節來降低地形場景的復雜性，從而提高繪制算法的效率。對地形原始數據建立幾個不同逐次逼近真實地形精度的幾何模型，與原始地形數據相比，每個LOD等級模型均保留了一定層次的細節。在繪制時，根據不同的選擇標準確定適當的層次地形LOD等級來繪制地形場景。地形場景LOD等級確定方法通常有3類：基于距離和物體尺寸標準的方法，基于偏心率、視野深度、運動速度等標準的方法和考慮保持恒定幀率的方法。本文忽略視距較遠的地形細節，因此選擇距離標準計算對應的分辨率。

首先對紋理數據進行分塊分級，建立金字塔結構模型。根據地形區域范圍和位置對紋理數據進行分塊，再對每塊數據進行不同分辨率的分級處理，分辨率分為256×256、128×128、64×64、32×32和16×16五 級。紋理數據塊的分辨率取決于塊到視點的距離，距視點越遠分辨率越低；反之越高。

即使使用多層紋理，每塊地形的紋理數據量約為256 KB，若繪制地形區域為60塊，則地形紋理高達15 MB。而實際地形塊一般都多于60塊，所以在地形繪制過程中，存儲、裝載和可視化這么大的紋理數據的效率是很低的，往往會消耗計算機性能和內存。為很好解決這個問題，需要對紋理數據進行壓縮，本文采用DXT1格式對紋理數據進行壓縮。DXT是一種圖片壓縮格式[25]，使用了由S3 公司發明的一種有損圖像壓縮算法，現在已為絕大多數 3D 顯卡硬件所支持。DXT1格式主要適用于不具透明度的貼圖或Alpha貼圖，其中Alpha貼圖是一個8位的灰度圖片，用256 級灰度來記錄圖像中的透明度，定義透明、不透明和半透明區域，分別用白、黑和灰色表示。DXT1將每4×4個像素塊視為一個壓縮單位，壓縮后占用64位，其中有2個16位的RGB顏色和16個2 位索引，格式描繪如圖5所示。

圖5 DXT1格式壓縮示意圖

圖6為采用DXT1格式的多分辨率紋理數據，數據分辨率逐漸成倍減小，數據塊的大小也隨之減小，地形越來越模糊。

圖6 DXT1格式的多分辨率紋理數據

紋理數據經過DXT1格式壓縮后，存儲和實時裝載速度大大提高，內存消耗少，極大提高了大場景地形的繪制效率。

1.3 GPU的加速處理

近年來，GPU廣泛應用于圖像處理和科學計算領域，GPU的大規模并行硬件使得數以萬計的線程并行成為可能。同時，采用CUDA和OpenGL圖形渲染混合編程方法是可視化技術的重要應用方向，也是提高性能和編程效率的理想選擇。

壓縮后的紋理數據塊提高了裝載和存儲速度，但對每塊數據的處理仍然是串行的。對于大區域海量數據而言，隨著紋理數據塊數量級的增加，程序的延時導致了實時性差和響應慢，因此采用GPU加速技術，對每塊紋理數據并行處理，提高數據處理的實時性和效率。

1.4 大區域三維地形的動態裝載和繪制

首先確定處于可視場景中的地形節點，再根據地形節點到視點的距離采用不同的細節層次進行繪制。紋理數據也同樣使用不同的分辨率來對地形塊進行映射，如圖7所示。

圖7 同一地形節點在不同場景中的分辨率

計算地形節點是否處于可視場景內，首先要計算場景可視矩陣，再判斷地形節點的4個角點是否在可視矩陣內。只要4個角點中有1個在可視矩陣內，則該地形節點被劃分到可視場景內。

確定需要進行繪制的地形節點后，根據地形節點距視點的距離和地形的精度誤差來計算地形節點的繪制等級，以實現地形快速繪制。為了適用于任意范圍的地形區域，本算法利用動態裝載和卸載相關地形塊數據的方法優化了計算機性能和內存。當某塊地形離開場景時，從內存中卸載掉該塊地形節點數據；而當某塊地形正從場景外進入場景時，則實時動態裝載該塊數據。這種優化方法使得算法中用來調度數據的內存始終只占固定的內存容量。

1.5 地形相鄰節點間裂縫的消除

根據地形節點距視點距離以及自身精度誤差的差異，繪制時使用的LOD等級層次也不同。當地形離視點近時，繪制的分辨率較高；反之，則采用較低的分辨率進行繪制。如圖8所示，顏色越深表示分辨率越高。本文采用一種適合GPU加速的裂縫消除算法。該算法要求地形模型符合2個要求：邊界點個數滿足n2+1（n≥3）；如果是RSG/TIN 混合結構，要求TIN 邊界點與RSG邊界點一致。該算法的基本思想是在低級地塊的邊界上等間隔插入備用點，地形場景實時繪制時，通過低等級LOD地形塊中備用點來參與地形場景繪制，從而消除相鄰地形塊繪制時的裂縫；最大的特點是既可用于規則格網RSG，也可用于TIN，且適用于RSG/TIN 混合結構。

圖8 地形可視化中多分辨率繪制

本算法中有關頂點和三角形的分類，如圖9。

邊界點：所有位于地形塊凸殼上的點為邊界點，如點E1、E2和E3。

圖9 頂點及三角形分類

內部點：除邊界點外所有參與地形場景構網的地形節點為內部點，如點I1、I2和I3。如果分別以E、I表示邊界點集合和內部點集合，那么構網點集合Q 就等于E 和I 的并集，即Q=E∪I，且E∩I=Φ。

備用點：每條邊界上兩相鄰邊界點的中點為備用點，如點S1、S2和S3。

邊界三角形：構成地形格網三角形中的3點至少有2點是邊界點的地形場景三角形格網為邊界三角形，如△E1E2I1和△E2E3I1。

內部三角形：構成地形格網三角形中的3點至少有2點是內部點的地形場景三角形格網為內部三角形，如△I1I2I3和△I1I3E3。

備用三角形：構成地形格網三角形中的3點至少有1點是備用點的地形場景三角形格網為備用三角形，如△E1I1S1、△E2I1S1、△E2I1S2和△E3I1S2。

邊角三角形：構成地形格網三角形中的3點全為邊界點的地形場景三角形格網為邊角三角形，如△E4E5E6。

該算法消除相鄰地形節點間的地形裂縫的步驟為：

1）備用點插入。在每條邊界的任意兩個邊界點中點處按同樣的順序插入等間隔備用點，備用點插入后所有邊界點仍然保持從左到右、從下到上、從右到左、從上到下的順序，如圖10所示。

圖10 接邊算法之插入備用點

2）邊角三角形分割。首先判斷地形塊中是否存在邊角三角形。若存在，則在邊角三角形的斜邊上強行插入圖11b中的點I1，對邊角三角形進行分割，分割后的三角形與原邊角三角形共斜邊。

圖11 邊角三角形分割

3）邊界三角形分割。使用備用點將所有邊界三角形一分為二，如圖12所示，并將分割后的三角形追加到地形塊結尾。

經過備用點插入和邊角三角形分割后，只需判斷是否接邊就可實現不同分辨率地形模型接邊處的裂縫消除。如圖13所示，當左右兩邊地形LOD等級不等時，在低等級的地形LOD塊使用備用點可很好地解決地形繪制接邊裂縫問題。

圖12 邊界三角形分割

圖13 經過接邊處理后的地形塊拼接

圖14是采用本文算法實現的地形場景，可以看出，靠近視點的區域使用高級LOD地形塊，遠離視點的區域使用低級地塊，同時在接邊處使用了其邊界備用點。

2 實驗與分析

實驗是在HP Z600 Workstation、9 GB RAM、NVIDIA GeForce GTX 750顯卡計算機上，用Windows Visual C++、 OpenGL和CUDA編程環境實現的。實驗設備的CPU、GPU配置分別為：CPU 8核Intel(R) Xeon(R)處理器，內存容量9 GB；GPU GeForce GTX 750，DDR5，1G顯存，80．2 GB/s。

實驗中，地形數據包括DEM和DOM，約為16 GB，對實驗數據建立了7層金字塔結構數據模型。紋理數據經過DXT1格式壓縮后，存儲數據量不到700 MB，僅使用CPU處理，數據預處理時間約為30 min，而使用GPU加速后，數據預處理時間為17 828 ms（約18 s）。地形繪制中的速度約為60 fps，每幀顯示三角形個數約為20 000個，地形繪制速度可滿足地形實時交互式瀏覽的要求。

圖14 地形繪制效果圖

3 結 語

本文針對三維地形數據量大，難以快速繪制的問題，把地形數據劃分為可視場景與非可視場景。將可視場景內部的數據動態載入內存，非可視場景的數據從內存卸載，減少了內存容量。首先，通過LOD技術將可視場景區域劃分為不同細節層次，根據視距劃分的不同場景確定相應的分辨率。然后，根據不同的分辨率對三維數據進行塔式結構模型表示，壓縮了數據量，并充分利用了顯卡的GPU加速，提高了大區域三維地形場景數據處理和快速渲染速度，較好地滿足了大區域三維地形場景實時交互式瀏覽的要求。參考文獻

[1] 游雄．基于虛擬現實技術的戰場環境仿真[J]．測繪學報, 2002,31(1)∶7-11

[2] Koller D,Lindstrom P,Ribarsky W,et al．Virtual GIS∶ a Real-time 3D Geographic Information System[C]．IEEE Computer Society Proceedings of the 6th Conference on Visualization'95, Atlanta,1995

[3] Lindstrom P,Koller D,Hodges L F,et al．Level-of-detail Management for Real-time Rendering of Phototextured Terrain[J]．GVU TR,1995(2)∶15-20

[4] Lindstrom P,Koller D,Ribarsky W,et al．An Integrated Global GIS and Visual Simulation System[J]．GVU TR,2007(2)∶38-44

[5] Lindstrom P,Koller D,Ribarsky W,et al．Real-time, Continuous Level of Detail Rendering of Height Fields[C]．Proceedings of the 23rd Annual Conference on Computer Graphics and Interactive Techniques, 2006

[6] Lindstrom P,Pascucci V．Visualization of Large Terrains Made Easy[C]．IEEE Visualization Conference,San Diego,2001

[7] Lindstrom P, Pascucci V． Terrain Simplification Simplified∶a General Framework for View-dependent Out-of-core Visualization[J]．IEEE Transactions on Visualization and Computer Graphics,2002,8(3)∶239-254

[8] Asano T,Ranjan D,Roos T,et al．Space-filling Curves and Their Use in the Design of Geometric Data Structures[J]． Theoretical Computer Science,1997,181(1)∶3-15

[9] Balmelli L,Ayer S,Vetterli M．Efficient Algorithms for Embedded Rendering of Terrain Models[C]．IEEE International Conference on Image Processing,Chicago,1998

[10] Balmelli L,Liebling T,Vetterli M．Computational Analysis of 4-8 Meshes with Application to Surface Simplification Using Global Error[C]．Proc． of the 13th Canadian Conference on Computational Geometry (CCCG),Waterloo,2001

[11] Cignoni P, Montani C, Scopigno R． A Comparison of Mesh Simplification Algorithms[J]．Computers & Graphics, 1998,22(1)∶37-54

[12] De Floriani L, Marzano P, Puppo E． Multiresolution Models for Topographic Surface Description[J]．The Visual Computer,1996,12(7)∶317-345

[13] Duchaineau M,Wolinsky M,Sigeti D E,et al． Roaming Terrain∶ Teal-time Optimally Adapting Meshes[C]．Proceedings of the 8th Conference on Visualization'97,AZ,1997

[14] De Floriani L, Kobbelt L, Puppo E． A Survey on Data Structures for Level-of-detail Models[M]． Springer Berlin Heidelberg∶Advances in Multiresolution for Geometric Modeling,2005

[15] Gerstner T,Pajarola R．Topology Preserving and Controlled Topology Simplifying Multiresolution Isosurface Extraction[C]．Proceedings of the Conference on Visualization, Salt Lake City,2000

[16] Pajarola R． Large Scale Terrain Visualization Using the Restricted Quadtree Triangulation[C]． IEEE Visualization'98 Conference,North Carolina,1998

[17] Hwa L M, Duchaineau M A, Joy K I． Real-time Optimal Adaptation for Planetary Geometry and Texture∶ 4-8 Tile Hierarchies[J]． IEEE Transactions on Visualization and Computer Graphics,2005,11(4)∶355-368

[18] Pajarola R, Gobbetti E．Survey of Semi-regular Multiresolution Models for Interactive Terrain Rendering[J]． The Visual Computer,2007,23(8)∶583-605

[19] Levenberg J．Fast View-dependent Level-of-detail Rendering Using Cached Geometry[C]． IEEE Visualization Conference,Boston,2002

[20] Losasso F,Hoppe H． Geometry Clipmaps∶ Terrain Rendering Using Nested Regular Grids[J]．ACM Transactions on Graphics (TOG),2004,23(3)∶769-776

[21] Hoppe H．Smooth View-dependent Level-of-detail Control and Its Application to Terrain Rendering[C]． IEEE Visualization'98,North Carolina,1998

[22] Pajarola R,Antonijuan M,Lario R． Quadtin∶ Quadtree Based Triangulated Irregular Networks[C]．Proceedings of the Conference on Visualization'02, Boston,2002

[23] Pajarola R,Ohler T,Stucki P,et al． The Alps at Your Fingertips∶Virtual Reality and Geoinformation Systems[C]． IEEE 14th International Conference on Data Engineering,1998

[24] Clark J A．Greenland's Rapid Postglacial Emergence∶ a Result of Ice-water Gravitational Attraction[J]．Geology,1976,4(5)∶310-312

[25] JIANG Yifei,HUAN Dandan． Improved Texture Compression for S3TC[C]．Picture Coding Symposium (PCS),Poznan,2010

P208

B

1672-4623（2016）09-0012-05

10.3969/j.issn.1672-4623.2016.09.004

李鶴元，博士，高級工程師，主要從事三維地形仿真以及網絡地理信息服務研究。

2015-06-23。

項目來源：國家科技支撐計劃課題資助項目（2012BAI32B07）； 國家自然科學基金青年基金資助項目（41001314）。