基于GPU并行處理的地形三維重建技術研究

楊秀峰，靳海亮，臧文乾

(1.河南理工大學，河南 焦作454000;2.中國科學院遙感與數字地球研究所，北京100101)

0 引言

隨著計算機科學和計算機圖形學的發展，三維地形的應用范圍不斷擴大，越來越多的涉及虛擬與現實、戰場環境仿真、土地管理與利用、地理信息系統、娛樂與游戲等大規模的三維地形應用中，因此，如何快速的實現大規模的地形重建是一個值得研究的問題［1］。

近年來，很多學者對基于點云的三維重建技術進行了研究。目前基于點云的三維重建技術大致可劃分為基于Delaunay三角剖分的方法、基于隱式曲面的方法和基于區域增長的方法［2］。在眾多基于Delaunay三角剖分的方法中，Amenta［3］等人提出的基于中軸逆變換的Power Crust算法最具代表性。針對點云密度不均勻、帶孔洞以及尖銳特征的任意點云，該方法都可重建出精密網格模型而不需要任何后期再處理，但是該算法需要進行復雜的三維Delaunay三角化處理，時間復雜度相當高，因此難以處理大規模地形的三維重建。基于隱式曲面的方法的核心是如何選取適當的隱曲面模型，其中比較著名的算法有Carr［4］等人提出的基于徑向基函數的方法，Alexa［5］等人提出的基于移動最小二乘的方法和Kazhdan［6］等人提出的基于泊松方程的方法。這類方法能夠較好的處理帶有噪聲的點云數據，但是容易丟失點云模型的原始數據信息和過渡平滑地物細節，并且針對大規模的點云，處理時間長。基于區域增長的方法一般都是從一個初始種子三角形出發，根據一定的判斷原則添加點，不斷向周邊擴張直到所有數據點都被處理。其技術難點在于初始三角形的確定和新加入三角形點的判斷。該方法可以實現大規模點云的三維重建，但是點云的重建質量過分依賴種子三角形的選取，并且在點云密度不均勻的情況下容易產生孔洞。Bemardin［7］提出的基于α-shape理論的滾球法(Ball Pivoting Algorithm，BPA)是比較著名的。目前，除了對上述單一方法的研究外，針對各種方法的融合使用也進行了大量的研究。王先澤［2］等人提出了一種基于Delaunay三角剖分和區域生長相結合的方法。首先對點云模型進行Delaunay三角剖分，然后從區域的邊界邊上添加新的三角面片，直到生成一張完整的三角網格模型。該方法在重建效率上有所提高，并能在一定程度上避免狹長三角形的產生從而提高生成網格的質量，但是針對大規模點云的重建，耗時仍是難以克服的問題。

整體而言，無論是單一的方法，還是多種方法的融合使用都無法有效的解決大規模點云三維重建的耗時問題。因此，本文提出一種基于GPU并行化的局部平面投影的地形三維重建方法，實現了三維地形的快速重建。

1 基于局部平面投影的地形三維重建算法

基于局部平面投影的方法是一種降維的三角剖分計算，在二維平面上實現三維重建。首先計算每個點的鄰域信息，然后將每個點的鄰域點投影到過該點的切平面上進行該點局部的Delaunay三角剖分，最后返回三維空間合并網格以及網格法向歸一化調整。在本文基于點云的地形三維重建中主要包括k近鄰計算、法向量計算、投影、徑向排序和Delaunay近鄰計算。

1.1 k近鄰計算

k近鄰計算一直是計算機圖形學、機器學等領域的難點，特別當點云數據量大時，k近鄰計算更是一個耗時的工作。因此，選擇一個好的k近鄰計算方法對于提高整個算法性能是至關重要。近似最近鄰(approximate nearest neighbor，ANN)算法庫是近年來運用比較廣泛的一種k近鄰計算方法，它是基于kd樹搜索的一種方法，相對與傳統的方法在性能上有75%提高，因此本文采用ANN算法庫實現點k近鄰。

1.2 法向量計算

法向量是點云具有的一個重要的局部特征，是進行點云處理必不可少的信息。本文中法向量的求解采用平面擬合的方法，主要用于鄰域點的投影和排序，以及網格法向歸一化調整。假設一點為p，其k近鄰集合為Q={qi}，i=1，2，3，…，k，令=qi-p，應用最小二乘法，該點的法向量使得式(1)最小:

可得到3×3矩陣C:

容易證明，對于協方差矩陣C，其最小特征值對應的特征向量進行單位化可作為法向量的近似值。

1.3 投影

如圖1所示，假設T為過某點的切平面，Q'= {q'i}，i=1，2，3，…，k為該點的領域點在其切平面上的投影點集合，則q'i可由在與所確定的平面內旋轉得到

其中，q^i是qi在切平面上的正交投影

圖1 旋轉投影

1.4 徑向排序

1.5 Delaunay近鄰計算

對切平面上的投影點進行徑向排序之后，利用二維平面三角剖分原理，從投影點集合中確定每個點的Delaunay鄰近點，從而實現每個點的局部重建。如圖2所示，l2是線段pq'i的垂直平分線，M是其中點，I是 pq'i-1的垂直平分線 l1與pq'i+1的垂直平分線l3的交點。如果I點不在p點與M點之間，則q'i是p點的有效Delaunay鄰近點，然后繼續下一點的Delaunay近鄰判斷;否則刪除該點，并回溯判斷〈q'i-2，q'i-1，q'i+1〉的有效性從而再次判斷p'i-1是否為有效的Delaunay鄰近點。直到沒有點標記為無效的Delaunay鄰近點，Delaunay鄰近計算完畢。

圖2 Delaunay近鄰點計算

2 基于GPU并行處理的地形三維生成

近年來，圖形處理器(Graphic Processing U-nit，GPU)在通用計算構架方面的飛速發展，使其計算性能不斷提高，計算能力已遠遠超過了通用處理器CPU，為加速三維重建提供了新的技術手段。2007年Nvidia公司正式發布的CUDA編程模型引發了GPU通用計算的革命，其將CPU作為主機，而GPU作為協處理器或者設備(device)。在模型中，CPU和GPU協同工作，其中CPU主要負責處理邏輯性強的事物和串行計算; GPU主要負責執行高度線程化的并行處理，實現大數據量的密集運算。

2.1 算法框架

通過對局部平面投影方法的分析，基于CUDA框架實現地形的三維重建中，鄰域點的計算、網格合并和法向歸一化在CPU完成，而針對每一個點的局部重構則在GPU中實現，其流程圖如圖3所示。

圖3 基于CUDA的地形三維重建流程

首先，在CPU上讀取點云數據，并調用ANN算法庫進行k近鄰計算。然后，將數據傳入GPU中，接下來正式在GPU上進行局部平面投影算法。利用CPU控制內存分配和啟動內核函數，而針對每一個點的局部平面投影算法的每個步驟都在GPU上執行，從而實現了并行處理。最后，在CPU上完成網格合并和法向歸一化調整。

2.2 算法實現及結果分析

所采用實驗設備為:Intel(R)Core(TM) i7CPU，主頻1.6 GHz，Windows7操作系統，顯卡為NVIDIA Quadro FX 880M。分別采用CPU模式下的局部平面投影方法和GPU模式下的局部平面投影方法對文家溝(圖4所示)進行了三維地形重建，其結果如圖5所示，并將2種重建方法的時間進行了比對，如表1所示。從統計的結果可以看出，利用基于GPU的局部平面投影算法重建地形的三維網格比基于CPU的方法有顯著提高，特別在法向量計算和角度計算2個步驟中，分別有30多倍和40多倍速度的提升。

圖4 文家溝三維點云模型

圖5 文家溝三維重建網格模型

3 總結與展望

面對大規模地形三維重建耗時的問題，通過GPU并行處理的方式，使大規模網格生成的時間大幅降低。描述了局部平面投影算法的前提下，提出了一種基于GPU的局部平面投影算法，并基于CUDA框架實現了汶川文家溝三維地形的快速重建。實驗結果表明，利用基于GPU的局部平面投影算法比傳統的CPU模式下的局部投影算法有顯著的提高，可見，在大規模地形三維重建中，GPU相對于CPU都有著相當大的優勢。然而，由于實際點云數據的多樣性和復雜性，總是會有點云稀疏或不均勻的情況，因此，重建過程或多或少都會導致洞的產生，孔洞修補成為必不可少的步驟。下一步工作，如何在GPU上實現孔洞的修補，從而快速、高質量的實現大規模地形的三維重建。

表1 運行時間結果對比

［1］ 徐 超.三維地形快速生成算法［D］邯鄲:河北工程大學，2011.

［2］ 王先澤，李忠科，馬亞奇，等.基于Delaunay三角剖分和區域生長的散亂點云重構［J］.四川工兵學報，2012，33(5):108-111.

［3］ Amenta N，Choi S，Kolluri R K.The power crust，unions of balls，and the medial axis transform［J］.Computing Geometry，2001，19(2):127-153.

［4］ Carr J C，Beatson R K，Cherrie J B，et al.Reconstruction and representation of 3D objects with radial basis functions［C］.Proceedings of the 28thannual Conference on Computer Graphics and Interactive Techniques，2001:67-76.

［5］ Alexa M，Behr J，Cohen-Or D，et al.Computing and rendering point set surfaces［J］.Visualization and Computer Graphics，IEEE Transactions on，2003，9 (1):3-15.

［6］ Kazhdan M，Bolitho M，Hoppe H.Poisson surface reconstruction［C］.Eurographics Symposium on Geometry Processing，2006.

［7］ Bernardini F，Mittleman J，Rushmeier H.The ball-Pivoting algorithm for surface reconstruction［J］.IEEE Transactions on Visualization and Computer Graphics，1999，5(4):349-359.