移動最小二乘重采樣在三維重建中的應用

康傳利,時滿星,程 耀,張臨煒,顧峻峰,陳 洋

(桂林理工大學 a.廣西空間信息與測繪重點實驗室;b.測繪地理信息學院,廣西 桂林 541006)

隨著三維激光掃描儀的普及和計算機硬件水平的提高,地面三維激光掃描儀被廣泛應用于數(shù)字化城市[1]、古建筑文物保護[2]、虛擬現(xiàn)實[3]等領(lǐng)域。為了真實重建具有不規(guī)則特征的實體,通常直接構(gòu)建其表面三角網(wǎng),這種曲面重建的方法主要分為兩類,即顯式曲面重建和隱式曲面重建。顯式曲面重建主要是基于Delaunay三角剖分和voronoi圖的計算幾何方法,如Crust算法[4]、tight Cocone算法[5]、3D alpha-shape算法都是在Delaunay三角剖分的基礎(chǔ)上發(fā)展而來的,若點云數(shù)據(jù)非均勻且含有噪聲,使用該種重建方法會丟失曲面重建的精度，并導致幾何模型變形。三維激光掃描儀采集的點云數(shù)據(jù)按點與點之間的排列方式可以分為:線式點云數(shù)據(jù)、陣列式點云數(shù)據(jù)、部分有序點云數(shù)據(jù)和完全散亂點云數(shù)據(jù)。地面三維激光掃描儀的點云數(shù)據(jù)屬于陣列式結(jié)構(gòu),具有明顯的空間拓撲關(guān)系，但為了全方位地獲得實體表面的空間數(shù)據(jù)信息,需進行多站點間的配準。因每兩個測站至少要有20%以上的重疊區(qū)域,經(jīng)配準后的點云數(shù)據(jù)拓撲關(guān)系將會丟失,數(shù)據(jù)冗余度變大;且受配準精度的影響,任意兩站點的數(shù)據(jù)不能完全融合,小尺度噪聲將會增加,不利于曲面重建。Aurich等[6]研究了基于線性濾波器的高斯濾波算法,缺點是無法適用于不存在拓撲關(guān)系的散亂點云數(shù)據(jù)。馬先明等[7]對無人機影像生成的點云數(shù)據(jù)利用雙邊濾波法進行了去噪實驗和分析,通過建立點云鄰域?qū)π⌒驮肼朁c進行了糾正,缺點是若噪聲點影響到鄰域特征的估計,則會減弱該噪聲算法的穩(wěn)定性。宋陽等[8]針對點云數(shù)據(jù)中體外孤點不易辨別難以去噪的問題,提出了一種改進的C均值算法,有效去除了大尺度噪聲,但該方法對點云數(shù)據(jù)中小尺度噪聲的去噪效果不佳。

針對上述散亂點云數(shù)據(jù)噪聲難以用統(tǒng)計分析的方法消除以及數(shù)據(jù)存在冗余從而影響曲面重建質(zhì)量和效率的問題,本文提出了一種基于移動最小二乘(MLS)平滑去噪重采樣的方法,既獲得了隱式曲面高抗噪性的光順效果,又保留了實體表面的細節(jié)特征。

1 MLS算法描述

Lancaster等以最小二乘法為理論基礎(chǔ)最早提出了移動最小二乘法,之后，相關(guān)學者對其進行了改進,提高了MLS的自適應能力,同時也繼承了最小二乘法高精度的數(shù)值計算能力。本文MLS算法的計算流程如圖1所示：首先建立空間數(shù)據(jù)點Pi的鄰域,有基于規(guī)則采樣的歐氏鄰域和基于不規(guī)則采樣的k近鄰鄰域;其次構(gòu)建以Pi鄰域為局部區(qū)域的擬合函數(shù),擬合函數(shù)中基函數(shù)的階數(shù)不同，其擬合精度也不同;之后確定權(quán)函數(shù),擬合函數(shù)能夠很好地繼承權(quán)函數(shù)的連續(xù)性,如果權(quán)函數(shù)是C1階連續(xù)的,則擬合函數(shù)也是C1階連續(xù)的[9],不同的權(quán)函數(shù)平滑效果也不同;通過設(shè)置體素化網(wǎng)格(VoxelGrid)的大小,將點云數(shù)據(jù)所在的空間進行體素柵格化處理,用體素網(wǎng)格內(nèi)所有點的重心近似代替網(wǎng)格中所有的數(shù)據(jù)點,從而簡化冗余數(shù)據(jù)點,最后得到重采樣后的點云。目前大多數(shù)平滑去噪算法在進行去噪時都沒有考慮到點模型的優(yōu)化,點云平滑前后數(shù)量未發(fā)生改變,本文算法將體素化網(wǎng)格下采樣模型引入到MLS平滑去噪中,與單獨使用MLS重采樣算法相比,在點云平滑去噪的同時實現(xiàn)三維點模型的優(yōu)化,能夠有效去除小尺度噪聲并解決點云冗余問題。

圖1 MLS算法流程圖Fig.1 Algorithm flow chart of MLS

2 MLS算法基本原理

2.1 擬合函數(shù)的建立

在擬合區(qū)域的一個局部子域U上,其擬合函數(shù)f(x)表示為

(1)

式中:P∈U;α(x)=[α1(x),α2(x),…,αm(x)]T為待求系數(shù),它是一個關(guān)于坐標x的函數(shù);P(x)=[P1(x),P2(x),…,Pm(x)]T稱為基函數(shù),是一個k階完備多項式;m是基函數(shù)的項數(shù)。

對于二維區(qū)域,通常基函數(shù)的形式為

線性基P=[1,x,y]T,m=3;

二次基P=[1,x,y,x2,xy,y2]T,m=6。

在實際運用中,二次基函數(shù)使用較多,所以式(1)也可以描述為

f(x)=α0(x)+α1(x)x+α2(x)y+α3(x)x2+

α4(x)xy+α5(x)y2。

(2)

(3)

式中:n為受影響區(qū)域內(nèi)節(jié)點的數(shù)目;f(x)為擬合函數(shù);yi=y(xi)是x=xi處的節(jié)點值;w(x-xi)是節(jié)點xi的權(quán)函數(shù)。

α=(BWBT)-1BWy,

(4)

W為n×n的對角矩陣;y=[y(x1),y(x2),…,y(xn)]T。

由式(4)能夠求出式(1)中的系數(shù), 再將待定點代入式(1), 可得到該點的擬合函數(shù)值。

2.2 確定權(quán)函數(shù)

移動最小二乘法、最小二乘和加權(quán)最小二乘的主要區(qū)別是在權(quán)函數(shù)的運用上,如引入緊支(compact support)的概念。認為點x處的值僅受x子域中相鄰節(jié)點的影響,該子域稱為點x的影響區(qū)域,影響區(qū)域一般是半徑為h的圓形,權(quán)函數(shù)在x的影響區(qū)域內(nèi)不為0,在影響區(qū)域外全為0，并且權(quán)函數(shù)w(x-xi)應該是非負的, 權(quán)值隨著‖x-xi‖2的增加而單調(diào)遞減。常用的權(quán)函數(shù)主要有:

①三次樣條權(quán)函數(shù)

(5)

②高斯權(quán)函數(shù)

(6)

式中:di(x)表示第i個原始采樣點Pi與重采樣點P之間的距離；u表示影響區(qū)域中的特征與重采樣間的影響因子。

在空間點云數(shù)據(jù)重采樣的過程中,重采樣點的位置相對于原始采樣點會發(fā)生變化,因此式(4)中權(quán)因子w(x-xi)應隨著待求點位置的不同而進行移動,本文選用高斯權(quán)函數(shù),能夠?qū)κ?4)進行重新計算,得到擬合函數(shù)新的系數(shù)矩陣

α(x)=(BW(x)BT)-1BW(x)y,

(7)

式中:W(x)是由式(6)計算得到的對角矩陣。

2.3 體素化(VoxelGrid)網(wǎng)格下采樣模型

體素化網(wǎng)格下采樣的基本原理是[10-11]:在點云所在坐標系下建立該輸入點云的最小包圍盒,包圍盒是與各坐標平面平行且完全包含點云數(shù)據(jù)的最小六面體,根據(jù)用戶輸入采樣距離的大小,將六面體等間距均勻分解成m×n×l個小網(wǎng)格立方體,所有的點云數(shù)據(jù)被再次劃分到小網(wǎng)格立方體中,這些小網(wǎng)格立方體可看作是子包圍盒,刪除子包圍盒之外的小網(wǎng)格立方體,再分別計算子包圍盒中所有點的重心,以該重心來近似代替這些點。相對于其他點云簡化方法,該種方法簡單易實現(xiàn),特別是對于散亂點云不需要再次建立點云間的拓撲關(guān)系,減小了計算的復雜度。該采樣模型的實現(xiàn)步驟如下:

①最小包圍盒的確定。找出輸入點云坐標系下各個坐標軸方向、點云坐標的最大值和最小值,分別記為Xmax、Xmin、Ymax、Ymin、Zmax、Zmin。

②根據(jù)用戶輸入的采樣距離L,確定最小包圍盒在各坐標軸方向被均勻分割出的立方體的個數(shù),即m、n、l的值。

(8)

③建立子包圍盒的索引。點云數(shù)據(jù)中的任一點Pi(xi,yi,zi)都被點云數(shù)據(jù)最小包圍盒的子包圍盒所包圍, 所以能夠根據(jù)這個點的空間位置確定其子包圍盒的空間索引(mindex,nindex,lindex)。

(9)

④子包圍盒內(nèi)近似代表點的確定。建立好子包圍盒的索引后,根據(jù)該子包圍盒內(nèi)所有點的個數(shù),求出它們的重心,并用這個重心來代表這些點。若某個子包圍盒內(nèi)有n個點,則代表點的坐標為

(10)

3 實驗結(jié)果與分析

3.1 實驗環(huán)境

在VS2010環(huán)境下采用C++語言和PCL實現(xiàn)了本文的算法,即完成空間點云數(shù)據(jù)的平滑去噪和重采樣。實驗的硬件配置為Intel(R) Xeon(R) CPU,2.40 GHz,4.00 GB內(nèi)存,Windows 7操作系統(tǒng),對采樣后的點分別在Geomagic Studio 2013和Matlab 2014中進行了顯式曲面重建。

3.2 實驗數(shù)據(jù)

主要測試數(shù)據(jù)來源于地面三維激光掃描儀Riegl vz1000對桂林理工大學雁山校區(qū)東門前景觀石的掃描, 為了全方位的獲取景觀石表面的空間數(shù)據(jù)信息, 此次掃描共架設(shè)5個站點, 掃描角分辨率為0.13°, 掃描時長30 min, 在掃描儀自帶的聯(lián)機控制軟件Pro riscan進行了多站點的ICP(iterative closest point)配準, 整體配準精度為0.003 2 m。圖2a為景觀石前方表面單個站點的掃描數(shù)據(jù),點云排列整齊、凹凸有致;圖2b為配準后的掃描數(shù)據(jù),點云的空間關(guān)系被完全打亂,且受配準精度的影響,具有相同位置的點云并不能完全重合,點云的冗余度較大。

圖2 配準前后的點云Fig.2 Point cloud before and after registration

3.3 基于Delaunay三角剖分的Crust三維重建

設(shè)置體素化網(wǎng)格小立方體的邊長為0、2.0、2.5、3.0 cm,分別與圖3和表1中下采樣1(僅MLS重采樣)、下采樣2、下采樣3、下采樣4相對應,原始點云是指景觀石的原始掃描數(shù)據(jù)。

點云名稱點云個數(shù)三角形個數(shù)點數(shù)簡化率/%Delaunay三角化耗時/s曲面重建總耗時/s原始點云455092905936085.84185.95下采樣1455092910172095.88194.44下采樣211033822066875.822.2846.99下采樣37403314806083.715.0531.61下采樣45301310602088.410.6522.43

從圖3b—d可以看出,經(jīng)體素化網(wǎng)格下采樣后的點云在空間結(jié)構(gòu)上具有一定的拓撲關(guān)系,且下采樣點云密度在具有表面特征的區(qū)域較大,如景觀石表面刻字部分的點個數(shù)比同等范圍平緩區(qū)域要多,即曲率變化大的區(qū)域點云分布比較集中,避免了細節(jié)特征信息的丟失。從圖4a和表1可以看出,存在冗余數(shù)據(jù)的點云會嚴重影響曲面重建的時間開銷，且三維重建效果不佳,正常分辨率下因重建網(wǎng)格過于稠密無法顯示出三維立體效果;從圖4d中可以看出景觀石表面“字跡”較為模糊,主要因采樣距離設(shè)置過大,劃分進入子包圍盒的點數(shù)較多,相應刪除了較多的點,從而造成局部特征不明顯。從表1可以看出,經(jīng)MLS重采樣的點云,其中滿足Delaunay三角剖分空外接圓和最小內(nèi)角最大準則條件的點數(shù)增加,故下采樣1重建三角形的個數(shù)要高于原始點云,曲面重建時間有所增加，通過對Delaunay三角剖分和曲面重建總耗時的比較,采樣間距越大點云簡化率越高，其表面重建時間越短。

對圖3、 圖4和表1在點云分布情況、 曲面重建立體化效果和曲面重建時間開銷等進行綜合性分析可以得出, 體素化網(wǎng)格下采樣3(2.5 cm×2.5 cm×

圖4 Crust三維重建部分截圖Fig.4 Partial screenshot of the Crust 3D reconstruction

2.5 cm)能夠在保證重建質(zhì)量下同時提高重建效率, 點數(shù)簡化率為83.7%, 重建時間縮短了83%。 因此，針對點云數(shù)據(jù)冗余的情況, 能夠根據(jù)點云的平均密度和局部細節(jié)特征選取合適的采樣間距進行簡化處理,提高曲面重建效率。

3.4 基于Geomagic Studio的三維重建

為了驗證本文算法的適用性及在實際建模中的魯棒性,對MLS重采樣的點云在逆向工程軟件Geomagic Studio中進行了真三維重建,并截取了景觀石表面帶有明顯細節(jié)特征(“理”的字跡)的部分點云作進一步驗證和分析(圖2b)。圖5a—c是原始點云分別在無平滑去噪、平滑去噪強度低、平滑去噪強度高的條件下封裝成的三角網(wǎng)格表面；圖5d為無去噪三角網(wǎng)格圖5a表面經(jīng)Laplace網(wǎng)格平滑去噪后的效果圖；圖5e是經(jīng)MLS重采樣后的點云在無平滑去噪條件下封裝成的網(wǎng)格表面。

在無去噪的情況下對配準后的點云數(shù)據(jù)進行網(wǎng)格化重建,其重建表面布滿了鱗片狀的褶皺,字跡處也出現(xiàn)了漏洞(圖5a),說明因配準生成的小尺度噪聲會嚴重影響曲面重建的質(zhì)量;比較圖5b、c表面的重建效果,隨著去噪強度的提高,重建表面也愈趨近于光滑,但圖5c表面的字跡因過于光滑已經(jīng)變得模糊,造成細節(jié)特征信息的丟失;比較圖5d和e,因Laplace是采用多次迭代的方法將數(shù)據(jù)點向該點鄰域中心移動基于網(wǎng)格的平滑去噪算法,故無法修復圖5a中由噪聲形成的漏洞；但從圖5e中可以看出，經(jīng)本文算法重采樣后的點云其表面重建效果較好,不僅整體表面較為平滑,同時字跡特征也能清晰地顯現(xiàn)出來且無漏洞,視覺效果優(yōu)于Laplace平滑去噪算法。

基于本文算法,將景觀石重采樣的點云數(shù)據(jù)導入到Geomagic Studio中直接封裝成網(wǎng)格模型,如圖6a所示,可以看出表面的整體光滑度高,局部細節(jié)特征清晰;為了給模型增加真實感的視覺效果,對圖6a進行紋理貼圖,生成景觀石的真三維實體模型如圖6b所示。

圖6 景觀石實體模型Fig.6 Solid model of landscape stone

4 結(jié) 論

點云數(shù)據(jù)的平滑去噪和簡化,對模型三維立體重建的質(zhì)量和快束存儲、處理、顯示等具有重要的意義。將移動最小二乘重采樣運用到地面三維激光掃描儀對大型構(gòu)筑物的三維重建上,能夠很好地解決因配準誤差不可避免而造成數(shù)據(jù)冗余和小尺度噪聲增加的問題。經(jīng)實例分析表明:

(1)選擇合適的體素化網(wǎng)格小立方體的邊長,能夠有效簡化海量點云數(shù)據(jù),且重采樣后的點云在曲率變化大的地方分布多，平緩區(qū)域分布少,在提高曲面重建效率的同時也保證了曲面重建的質(zhì)量。

(2)將移動最小二乘重采樣后的點云直接在逆向工程軟件中封裝成網(wǎng)格模型,其表面細節(jié)特征清晰、光滑度好,跳過了后期打磨、松弛和去除特征等修理環(huán)節(jié),再對網(wǎng)格模型進行紋理貼圖生成真三維實體模型,逼真程度高、視覺效果好,優(yōu)于Laplace網(wǎng)格平滑去噪算法。

圖5 平滑去噪效果Fig.5 Results of smoothing and denoising

因此,移動最小二乘重采樣在地面三維激光掃描儀三維重建中具有很高的實用價值。