程志全, 葉永凱, 李 寶
(國防科學技術大學計算機學院,湖南 長沙 410073)
將數據采集設備獲取的實物樣件數據進行處理和三維重構,反求工程在計算機上復現實物樣件的幾何形狀,并在此基礎上進行原樣復制、修改或重設計,在工程實驗、計算機圖形和計算機輔助幾何設計等方面有著廣泛的應用。其中,如何從掃描得到的離散數據(通常為點云模型)中提取原始模型外表的幾何基元(多為二次曲面),是對模型進行各種處理的基礎的、關鍵性的工作。基元提取擬合給定點云模型的最佳二次曲面,用點云與提取出的二次曲面之間的“最小距離”對擬合性能進行評價,以確定擬合曲面是否真正充分表達了整個模型[1]。當前各基元提取算法多針對平面、球、圓柱、圓錐和可擴展曲面等相對簡單的二次曲面,很少涉及橢球等復雜的基元[2-3]。面向幾何基元提取問題,本文提出了一種從點云模型中有效提取復雜橢球的算法。
基元幾何基元的檢測和提取一直備受研究人員的關注。多年來,有大量不同的方法被提出,限于篇幅,本文將不在此對這些方法進行詳細地討論,僅對一些最重要的方法做一個簡要的回顧。
RANSAC(RANdom SAmpling Consensus,隨機采樣一致性)[2]和Hough變換是用于幾何基元提取的最常用的兩種方法。有大量的噪聲和外點的情況下,二者都被證明能成功地檢測出二維和三維空間的形狀,但是低效和高存儲容量要求是它們的主要缺點。為此,COHEN-STEINER等[3]提供了一個通用的變分框架,使用多個平面來逼近表面。WU等[4]將此方法進行了擴展,使用一組更加復雜的基元(平面、球面和可擴展的圓環)來逼近表面。他們的目標不僅僅是檢測幾何基元,而且還試圖尋找給定基元數目時的全局最優表示。RUWEN等[2,5]增加基元的類型(平面、球體、圓柱、圓錐和圓環),并提出了一種高性能的 RANSAC算法,能在提取不同類型的基元形狀的同時,保持基本 RANSAC算法框架的良好特性,如一般性和簡潔性。橢球可以看作廣義圓柱體的一個特例,具有一定的對稱性,相對于球也能更靈活地表現較復雜的幾何信息。本文基于文獻[2]中的算法框架,實現點云中橢球基元的有效檢測與提取。
本文采用RANSAC框架(偽代碼見圖1)完成橢球的提取工作:使用一種局部的采樣策略進行采樣[2],從點云數據(點數為 N)中隨機選取最小點集,即能唯一地確定橢球的最少數目的點的集合,本文設點數為6(2.1節);利用最小二乘的方法計算橢球的參數,即橢球的中心、半軸長和偏轉角度(2.2節);基于距離和偏角參數驗證求得的橢球的正確性(2.3節);將通過驗證的橢球與點云中的所有點進行匹配,以判斷有多少點可以由此候選點基元近似表示,點的數目被稱作分數,完成了分數函數(2.4節)的計算,得到了評估橢球的方法,最后依據分數方程篩選出最好的候選點m。這個最佳的候選點只有在這種情況下才會被接受:由m中點的個數|m|和已經選出的候選點形狀的個數得到的能正確檢測出形狀的置信概率大于用戶定義的閾值pt。進而通過求投影點的方法繪制出檢測到的橢球。如果一個候選點被接受,那么相應的點會從點云中刪除,并這個候選也會從集合 C中刪除。算法在概率大于pt的時候停止,其中τ是用戶定義的最小的形狀包含的點的數目,這表明所包含的點超過τ的形狀已經被檢測出來的概率大于pt。算法用了一種標準的分數方程來計算一個候選點能容許的點的數目,即RANSAC算法所說的一致集。分數方程的輸入包括評估的候選形狀的參數,以及有兩個用戶定義的閾值ε和α:ε表示一個能容許的點到候選形狀的最大距離;α是一個點的法向到候選形狀的法向的最大偏移量。綜上所述,在檢測開始前,用戶需要定義以下的閾值:計算分數方程時的距離閾值ε、角度閾值α、最小的橢球包含的點的個數τ和置信概率閾值pt。
對于置信概率,我們假設任意6個點都可以產生一個候選點,那么一次能夠正確檢測到形狀的概率是

在選出了|C|個候選點之后可以正確的檢測出形狀的概率是

在進行采樣時,本文采用局部的采樣策略進行采樣[2],使用一個八叉樹來有效地建立樣本點間的空間近似度。當要選擇樣本點來產生新的候選點時,首先沒有任何約束地從空間中選出一個點,然后從八叉樹中任意一級隨機選擇一個包含該點的單元格,其他的5個點從這個單元格內選出。
橢球的估計是指從上面采樣的點中產生一個與點云中實際隱藏的橢球相近的一個橢球模型。橢球面的參數化方程如下

要解出這個方程中A到I這9個參數,需要在橢球面上的9個點,將這9個點帶入到公式(1)中,問題就變成了一個解線性方程組的問題了。而實際上,因為不能保證所抽取的樣本點都在橢球面上,就算有9個樣本點的數據,構成一個線性的方程組,這個方程組也可能無法解出這9個參數。本文利用橢球的幾何性質,先估計出橢球的中心,然后用最小二乘方法計算出橢球的軸長和偏轉角度。
橢球中心的估算方法的基礎為橢球中心的求解引理《空間解析幾何》[6]:對于橢球上的 4個不共面點和p4,任兩個點的中心和兩個點的橢球切平面相交所得的直線建立平面Pi,所有Pi相交于橢球的中心點。任取4個帶法向的點,建立相應的平面Pi(i=1,…4),平面方程為:。因為這4個點不一定在橢球面上,所以從中任取4個平面不一定交于一點。因此,本文選擇距離這4個平面的距離平方之和最小的點作為我們估計出來的合理的橢球中心,即求取下式的最小值

求出橢球中心后,將坐標系的原點平移到橢球中心得到公式

參照《空間解析幾何》[6]中所述的方法,取6個點分別代入公式,求出橢球的3個半長軸的長度和偏角。
設對稱矩陣S為

如果I, J, K都大于0,則證明所求系數代入上式中得到的二次曲面的確是橢球。設λ1,λ2,λ3為S的特征值,v1,v2,v3為對應的特征向量,那么 3個半軸的長度
由于估計出來的橢球可能存在誤差,在進一步地評估候選的好壞之前,先驗證用于估計橢球的樣本點本身是否能很好地符合估計出來的橢球,若估計時使用的所有樣本點都能很好符合,則驗證通過。為了驗證,必須得到兩個數據:樣本點到橢球的距離d,樣本點的法向與樣本點在橢球面上的投影處的法向之間的偏角θ。當兩個數據滿足所給的閾值條件即d<ε,θ<α時,即認為此樣本點可以很好符合橢球。
橢球的分數就是點云中距離和法向偏角滿足一定的閾值要求,且其投影點屬于橢球的最大連通分量內的點的個數。直接在橢球面上求最大連通分量是困難的。因為計算橢球面上兩點間沿橢球面的最短弧的長度是比較復雜的。在這里,為了簡化,我們先將點都投影到一個平面的參數域上,通過計算平面參數域上的最大連通分量求得橢球面上的最大連通分量。
將橢球面上的點投影到參數平面之后,先要將參數平面位圖化,然后在進行連通分量的檢測。位圖的像素的大小由原始數據的采樣分辨率確定。如果一個像素的區域內存在橢球面上點的投影,那么這個像素將被著色。連通分量的計算過程是這樣的,首先將每個像素的連通分量編號都設為-1,表示不屬于任何分量。逐行掃描位圖的每個像素,在掃描的過程中要記錄每個已檢測出的連通分量包含的像素點,如果掃描到的這個像素點未被著色,不做任何處理,如果這個像素著色了,則考察已經掃描過的與這個像素點相鄰的像素點,將這些相鄰像素點所屬的連通分量合并,再加入當前的這個點得到一個新的連通分量,當掃描完畢之后就找到了所有的連通分量。包含像素點最多的就是最大連通分量。
本文算法在 VC++平臺上用 OPENGL實現。為了測試本文橢球提取算法的有效性,本文在人工合成和掃描儀獲取的點云模型上展開了實驗,并使用運行時間、擬合誤差兩個指標對算法的性能進行了定量描述。
首先,本文使用了一個手工創建的橢球模型并記錄下其各個參數,并對此模型進行隨機采樣得到10萬個點;然后運行本文算法,將得到橢球參數與真實參數進行比對。參數設置為:最大距離 ε = 0 .008*BB (BB為包圍盒的對角線長度),最大偏角α=30°,最小形狀的大小為τ=200,置信概率 pt= 9 9.9999%。圖2中給出了球的檢測結果,算法運行時間為1.1233秒,平均距離誤差(指的是原始數據點到其在檢測出的形狀上的對應點的平均距離)為0.0003。表1則給出了原始橢球的數據與檢測出來的橢球的數據的比較,結果表明,本文的方法能夠正確地檢測出橢球。

圖2 人工合成的橢球的提取結果

表1 提取出的橢球參數與真實數據的比較
本文提出的算法應用到由激光掃描儀得到的原始點云模型,如圖3所示。可以看到,在此點云模型中,每個人物的頭部可以大致由橢球來進行擬合。我們為算法設置的參數為:最大距離ε= 0 .009*BB ,最大偏角α=40°,最小形狀的大小為τ=200,置信概率 pt= 9 9.9999%。實驗結果表明,在 RANSAC算法框架基礎上,本文提出的算法能夠在點云模型中檢測并提取出來可以由橢球擬合的部分(例如人的頭部),并得到很好的擬合結果。算法的運行時間為 58.5s,平均距離誤差為0.0407,檢測出51個基元。
基于RANSAC的形狀檢測已經能從點云中檢測出平面,球,圓錐,圓柱,圓環等形狀,為了能得到更為準確和復雜的基元表示,本文提出了一種在點云模型中提取橢球的算法。該算法使用了RANSAC框架;使用一種局部的采樣策略進行采樣;利用最小二乘的方法計算橢球的參數(中心點、半軸長和偏轉角度);利用距離和偏角進行橢球的初步驗證,并通過一種低扭曲的參數化方法得到參數化平面計算出最大連通分量,完成了分數方程的計算,得到評估橢球的分數;最后依據分數篩選出最好的候選,通過求投影點的方法繪制出提取的橢球。在人工合成和掃描儀獲取的點云數據上,本文驗證了所提算法的有效性和魯棒性。

圖3 真實掃描數據(上)的橢球(下)的提取結果
[1]劉金平, 程志全, 黨 崗, 等. 三角網格中二次曲面的提取方法綜述[C]//第 15屆計算機輔助設計與圖形學學術會議, 大連: 中國計算機學會, 2008: 134-139.
[2]Schnabel R, Wahl R, Klein R. Efficient RANSAC for point-cloud shape detection [J]. Computer Graphic Forum, 2007, 26(2): 214-226.
[3]Cohen S D, Alliez P, Desbrun M. Variational shape approximation [J]. ACM Transactions on Graphics,2004, 23(3): 905-914.
[4]Wu J, Kobbelt L. Structure recovery via hybrid variational surface approximation [J]. Computer Graphics Forum, 2005, 24(3): 277-284.
[5]Li B, Schnabel R, Jin S, et al. Variational surface approximation and model selection [J]. Computer Graphics Forum, 2009, 28(7): 1985-1994.
[6]楊文茂, 李全英. 空間解析幾何[M]. 武漢: 武漢大學出版社, 1998: 128-135.