葛寶臻 ，彭 博 ，田慶國

(1. 天津大學(xué)精密儀器與光電子工程學(xué)院，天津 300072；2. 光電信息技術(shù)教育部重點實驗室，天津 300072)

三維點云數(shù)據(jù)的配準與融合是計算機視覺領(lǐng)域一個基本而重要的問題．其研究成果在醫(yī)學(xué)、逆向工程、虛擬現(xiàn)實、計算機視覺等領(lǐng)域得到了廣泛的應(yīng)用[1]．根據(jù)不同視圖點云數(shù)據(jù)間是否存在畸變，現(xiàn)有的配準算法大致可以分為基于剛性變換的配準和非剛性配準兩類．本文只討論基于剛性變換的配準算法．給定具有重疊部分的三維點云數(shù)據(jù)，其任務(wù)是求取點云數(shù)據(jù)間的最佳剛性變換，使不同視圖的點云數(shù)據(jù)準確地融合在一起[2-4]．

迭代最近點(iterative closest point，ICP)算法[5]是當前最流行的剛性配準算法，該算法以相對位置經(jīng)過粗略配準的兩片點云數(shù)據(jù)為基礎(chǔ)，迭代求取兩片點云數(shù)據(jù)間的對應(yīng)點及相應(yīng)的剛性變換，直到對應(yīng)點之間的距離誤差評價函數(shù)最?。撍惴ǜ鶕?jù)兩點云數(shù)據(jù)點之間的距離來確定對應(yīng)點，所以需要待配準的兩點云數(shù)據(jù)之間基本已經(jīng)配準．實際上待配準的點云數(shù)據(jù)間的相對位置往往未知，怎樣利用粗略配準算法將點云數(shù)據(jù)間的相對位置變換到ICP 算法的收斂域成為問題的關(guān)鍵．

目前普遍的解決方法是通過不依賴于三維旋轉(zhuǎn)和平移變換的點云特征信息來建立點云數(shù)據(jù)間的對應(yīng)關(guān)系．這些特征信息包括法向[6]、曲率[7]、體積分[8]等簡單的一維特征描述函數(shù)，還包括旋轉(zhuǎn)圖像[9]、形狀內(nèi)容[10]等高維特征描述函數(shù)，有的甚至還包括顏色、強度值等輔助信息．但存在的問題是，低維特征描述函數(shù)計算簡單，但所包含的特征信息較少，且對噪聲敏感；高維特征函數(shù)能很好地反映特征信息，但計算復(fù)雜，時間、空間開銷較大，且怎樣合理地進行特征比較是一個難題．

筆者致力于通過新穎的特征信息和匹配算法快速、準確地得到點云數(shù)據(jù)間的對應(yīng)關(guān)系．受Gatzke等[11]的啟發(fā)，在以點幾何表示的三維點云數(shù)據(jù)中引入曲率圖的概念，并將其作為三維點云數(shù)據(jù)的特征描述函數(shù)．該描述函數(shù)的優(yōu)點在于通過選擇不同的環(huán)數(shù)及各環(huán)半徑值可得到不同的曲率圖，很容易實現(xiàn)多比例空間的特征保持分析．有了曲率圖的概念以后，本文的主要算法思想是通過一維曲率圖粗選特征點和潛在對應(yīng)點，然后對特征點和潛在對應(yīng)點進行高維曲率圖計算，并再次通過高維曲率圖進一步精簡潛在對應(yīng)點．這樣既保證了豐富的特征信息以獲得正確的對應(yīng)點，也大大減小了由于對整片點云進行高維曲率圖計算所花費的時間和空間開銷．在特征點的匹配算法上，對已有的Mitra 貪婪算法[12]進行了改進，大大提高了其穩(wěn)定性．最后通過對真實三維點云數(shù)據(jù)進行配準，驗證了該算法的有效性．

1 三維點云數(shù)據(jù)的曲率及曲率圖

在微分幾何中，對于三維歐幾里得空間中可微曲面上的任一點A，給定了兩個主曲率 1k 、2k[13]來衡量該點在不同方向上的彎曲程度．三維點云數(shù)據(jù)是對真實物體表面進行掃描采樣獲得的離散數(shù)據(jù)，點云數(shù)據(jù)中的每個頂點同樣具有主曲率，不但如此，它還是點云數(shù)據(jù)的一項重要特性．經(jīng)常提及和主曲率相關(guān)的幾個曲率概念還有高斯曲率 K = k1k2、均值曲率

對于三維點云數(shù)據(jù)，估算其頂點主曲率 1k 、2k 的方法有很多，具代表性的有多項式擬合(polynomial fitting)法[15]、主成分分析(principal component analysis)法[16]以及移動最小二乘(moving least-squares)法[17]等．其中多項式擬合法對噪聲敏感程度最小，穩(wěn)定性高，但計算相對復(fù)雜；主成分分析法求解過程相對簡單，但只適用于采樣點稠密的大型點云；移動最小二乘法具有較好的法向估算精度，但需要求解非線性方程．本文根據(jù)不同的點云數(shù)據(jù)采取與之相適合的主曲率估算方法．但僅利用單一點的均方根曲率C 并不能很好地反映該點周圍點的特征，并且對噪聲比較敏感．曲率圖利用該點及其周圍點的曲率信息分布，能夠準確地反映特征信息的同時有效地抑制噪聲．下面分別介紹以點幾何表示的三維點云數(shù)據(jù)中一維曲率圖和高維曲率圖的構(gòu)造方式．這里以降采樣后的點云數(shù)據(jù)Stanford bunny00[18]為例，對其表面上的任意一頂點A，一維曲率圖的構(gòu)造比較簡單，如圖1(a)所示，以目標點A 為球心，對半徑為R 的球體內(nèi)所有頂點的均方根曲率 Ci進行平均，得到的平均值即為該點的一維曲率圖．

圖1 曲率圖的概念示意Fig.1 Sketch of curvature map

式中n 為距離A 點的距離小于R 的頂點個數(shù)．對于高維曲率圖，如圖1(b)所示，同樣以目標點A 為球心，利用半徑為 rj(0 ＜ j ≤m)的同心球?qū)⑧徲蚩臻g劃分為m 個球殼層，然后對每個球殼內(nèi)的所有頂點的均方根曲率進行平均得到 f(rj)(0＜ j≤m)，高維曲率圖 kmap即為 f(rj)隨 rj的變化關(guān)系．用向量表示為

2 基于一維曲率圖提取特征點

有了曲率圖的概念以后，對三維點云數(shù)據(jù)的每個頂點進行一維曲率圖運算．根據(jù)所有頂點的一維曲率圖分布情況進行特征點的提?。@里采用的思想是對所有頂點的一維曲率圖進行統(tǒng)計分析，求取一維曲率圖的平均值μ和方差σ，選取一維曲率圖 f(R)處于μ± kσ之外的點作為特征點．一片點云的數(shù)據(jù)量很大，使得一維曲率圖往往近似服從正太分布，實驗中發(fā)現(xiàn) k = 1時，獲得的特征點數(shù)量大致為點云總數(shù)的10%，通過調(diào)整系數(shù)k 可選擇不同程度的特征點．為了使選取的特征點能更好地反映點云的特征信息，抑制噪聲，對每個頂點選取不同的鄰域范圍做多比例空間的特征分析．這里通過改變一維曲率圖的構(gòu)造半徑R ，例如對于半徑為 Ri和 Ri+1時的鄰域范圍分別做特征分析得到特征點集 Pi和 Pi+1，選取 Pi和Pi+1的交集作為特征點．最后的特征點集為

式中n 為做多比例空間特征分析選取的半徑個數(shù)．對Stanford bunny00 進行了基于一維曲率圖的特征點提取及特征保持分析．算法流程為：采用空間包圍盒法將三維點云數(shù)據(jù)結(jié)構(gòu)化；在此基礎(chǔ)上按歐式距離進行點云的K 鄰域搜索；采用主成分分析法估算主曲率及均方根曲率(K＝20)；對每個半徑R 分別計算每個頂點的一維曲率圖 f(R)(R ＝0.000 5、0.001、0.002、0.005,mm)；對每個半徑R 分別進行曲率圖統(tǒng)計分析，選取 f(R)處于μ± kσ之外的點作為該曲率圖半徑下的特征點集(k ＝1)；對上述4 個曲率圖半徑R進行多比例空間的特征保持分析，得到特征保持點集；對特征保持點集按一定間隔距離值(0.005,mm)進行采樣．

圖2所示為該算法的特征點提取過程．圖2(a)為R ＝0.000 5,mm 時提取到的特征點，可以看到，通過一維曲率圖提取到的特征點基本上都是處在兔子的耳朵、脖子、腿等輪廓特征比較明顯的部位，驗證了一維曲率圖提取特征點的有效性；圖2(b)顯示了通過多比例空間的特征保持分析后得到的特征保持點集，對比圖2(a)可知，特征保持分析能有效去除一些非特征點或特征不夠顯著的點；圖 2(c)為按0.005,mm 間隔對特征保持點集采樣的結(jié)果，對特征保持點集按一定距離采樣是為了進一步減少特征點數(shù)量，且讓特征點間有一定的距離，為后面的配準算法打下了基礎(chǔ)．

圖2 特征點的提取過程Fig.2 Process of extracting feature points

3 高維曲率圖精簡潛在對應(yīng)點

對于要配準的兩片點云數(shù)據(jù)，通常稱之為源點集和目標點集．理論上，只要在源點集和目標點集之間至少找出3 對對應(yīng)點，就能夠計算出這兩片點云數(shù)據(jù)間的剛性變換關(guān)系，從而將這兩片點云數(shù)據(jù)配準到一起．因此，準確找到點集間的對應(yīng)點成為配準的關(guān)鍵．

本文的思想是對源點集進行基于一維曲率圖的特征點提取，然后對目標點集進行相應(yīng)的一維曲率圖運算，一維曲率圖半徑為R ．對源點集中的每個特征點根據(jù)半徑為R 時一維曲率圖的相似性在目標點集中求取相應(yīng)的對應(yīng)點，即選取 fp(R) ? fq(R )≤Δ的點作為該特征點的潛在對應(yīng)點，Δ為設(shè)定閾值，設(shè)定原則是保證理論上具有對應(yīng)點的特征點都能在潛在對應(yīng)點中找到正確的對應(yīng)點．通過實驗發(fā)現(xiàn)，該值設(shè)為所有源點集中特征點的一維曲率圖平均值的10%時，基本能夠滿足上述要求．此時源點集中的一點 pi往往對應(yīng)目標點集中的多個點 qi1，qi2，…，qik．為進一步減少潛在對應(yīng)點的數(shù)量，有效去除錯誤對應(yīng)點，對pi及其潛在對應(yīng)點 qi1，qi2，…，qik進行高維曲率圖運算，高維曲率圖的半徑大小及環(huán)數(shù)根據(jù)點云采樣密度決定．然后將 pi的高維曲率圖分別與 qi1，qi2，…，qik的高維曲率圖做比較，比較的方式直接采用高維空間的歐式距離，即E2中選擇與 ie 沒有交集的2 點對 je ，使得 ie 、je 組成的4 點對的距離均方差dRMS 最小．每構(gòu)成一個新的4 點對，便從E2中刪除包含這4 點對中任意一點的2 點對，遍歷直到E2為空．最后得到4 點對集E4，同樣按照距離均方差的值對E4升序排列．這里距離均方差dRMS 的定義為

(3) 粗略匹配對應(yīng)點．同理可以對E4進行組合，得到8 點對集E8，乃至16 點對集E16等．在對應(yīng)點對不是太多的情況下一般得到E4、E8即可．選擇E4或E8中距離均方差最小的4 點對或8 點對作為粗略的對應(yīng)點．

(4) 精確匹配對應(yīng)點．在目標點集中以粗略對應(yīng)點為球心，半徑為 CR 做球，處于球內(nèi)的點作為潛在對應(yīng)點，遍歷所有的對應(yīng)點組合，選擇使dRMS 最小的點對組合作為最終的對應(yīng)點．CR 的選取值越大得到的結(jié)果更精確，但耗時也越大．

對Stanford bunny00 和Stanford bunny45 進行上述匹配算法，點對組合至E4，CR 取0.005,mm，得到4對對應(yīng)點，如圖3 所示．從圖中可以看到，兩視圖間的對應(yīng)點對基本上處于兔子的同一個位置，證實了該匹配算法的有效性．

式中：m 為高維曲率圖的維數(shù)；fp(rj)、fq(ri)分別為點 pi、qi的高維曲率圖分量．對 pi的所有潛在對應(yīng)點按d 的大小進行升序排列，選取靠前的M 個對應(yīng)的 qi作為 pi最終的潛在對應(yīng)點．M 的選取原則是盡量保證正確的對應(yīng)點包括在其中，而且值還不能過大，不然會給后面的匹配算法帶來困難．經(jīng)過實驗驗證M 值選取為10，即一個點具有10 個對應(yīng)點時，點云重疊部分的特征點基本都能找到正確的對應(yīng)點，而且對后續(xù)的匹配算法不會增加太多負擔(dān)．

4 精確匹配對應(yīng)點

目前主要匹配算法思想是利用形狀的空間變換不變性來有效匹配對應(yīng)點．在眾多的算法中，Mitra的貪婪算法雖然不能保證找到最佳的對應(yīng)點，但它能高效快速地找到接近最佳匹配的對應(yīng)點對．為了較為精確地匹配對應(yīng)點，筆者提出下面的兩步匹配算法：首先，通過貪婪算法粗略匹配對應(yīng)點，然后同樣利用形狀空間變換不變性，遍歷這些點某個鄰域范圍內(nèi)的所有點對組合，找出最佳對應(yīng)點．具體算法有如下4 個步驟．

(1) 構(gòu)造點對．對于源特征點集中的每個點對組合 pi和 pj，從其潛在對應(yīng)點中選擇 qik和 qjl，使得

圖3 精確匹配對應(yīng)點Fig.3 Accurately matching of corresponding points

5 結(jié)果與討論

(2) 組合點對．順序從E2中取出2 點對 ei，再從

針對Stanford bunny00 和Stanford bunny45，配準前如圖4(a)所示，利用匹配算法得到對應(yīng)點以后，運用奇異值分解(SVD)法[19]直接求取剛性變換，將Stanford bunny45 粗略配準到Stanford bunny00 上，得到圖4(b)所示結(jié)果，可以看到兩片點云數(shù)據(jù)基本上被配準在一起，但仍出現(xiàn)大片單一視圖的點云數(shù)據(jù)，需要進一步的精確配準．

為此，將粗配準后的兩片點云數(shù)據(jù)作為輸入，利用經(jīng)典ICP 算法進一步精配準，ICP 算法的具體算法流程為：

(1) 通過隨機采樣的方式從源點集中選擇一部分點集P；

中日醫(yī)院呼吸中心煙草病學(xué)與戒煙中心主任肖丹教授指出，無論吸什么樣的煙，只要吸煙就是有害健康的。大量研究表明，煙草煙霧中有超過7000種化學(xué)物質(zhì)，包括69種致癌物，依然含有尼古丁、丙二醇、丙三醇、重金屬等有害物質(zhì)，會增加患慢性阻塞性肺疾病、肺癌、心血管疾病的風(fēng)險，可以引起人體中多個器官的損害。

(2) 對點集P 中的任何一點，在目標點集中找尋與之歐式距離最近的點作為對應(yīng)點，得到點集P的對應(yīng)點集Q；

(3) 從對應(yīng)點集P 和Q 中刪除距離值大于平均距離的點，剩下的點作為最后的對應(yīng)點；

(4) 通過對應(yīng)點求取坐標變換矩陣(四元素法或SVD)；

(5) 構(gòu)造并計算誤差評價函數(shù)，其實就是第4 節(jié)提到的距離均方差；

(6) 判斷誤差評價函數(shù)的大小是否滿足精度要求，如果滿足則算法終止，返回結(jié)果；如不滿足，迭代進入步驟(1)，此時的目標點集為原始目標點集經(jīng)過步驟(4)獲得的變換后的結(jié)果．

通過ICP 算法的精配準得到最終的配準結(jié)果如圖4(c)所示，可以看到兩片點云數(shù)據(jù)交錯有致，達到了較好的配準效果．

圖4 兔子點云的配準過程Fig.4 Process of registration for Stanford bunny

為了驗證本文配準算法的通用性和適應(yīng)性，利用本實驗室的掃描系統(tǒng)，對圖5 所示的人頭模型進行掃描實驗．但這次掃描條件比較特殊，只利用其中的一個攝像機，分別從不同角度進行掃描，得到一系列具有重疊部分的三維點云視圖．

圖5 人頭模型Fig.5 Model of head

圖 6 所示為以任意兩角度拍攝的兩點云視圖．運用本文的配準算法模塊，對這兩點云視圖進行配準，得到如圖7 所示的結(jié)果．圖7(a)為配準前的兩點云視圖在同一坐標系中的表示，可以看出兩點云的視圖旋轉(zhuǎn)了一定的角度，且存在少量的平移；圖7(b)所示為本文算法粗配準后的點云視圖，這時兩點云數(shù)據(jù)已經(jīng)相當接近配準狀態(tài)；圖7(c)為經(jīng)過ICP 精配準的結(jié)果，可以看到視圖間的點云交錯有致，達到了良好的配準效果，驗證了本文算法的有效性．

為說明該算法的抗噪性能，對Stanford bunny00和Stanford bunny45 加入不同程度的白高斯噪聲．噪聲方差δ分別設(shè)為原始點云相鄰頂點間平均距離的0.1、0.3 和1.0 倍．為了和同類算法進行比較，在1.0倍噪聲時，分別利用旋轉(zhuǎn)圖像算法[6]和本文的算法進行粗配準．圖8(a)、8(b)分別表示在0.1 倍和0.3 倍噪聲下本文算法的配準結(jié)果，圖8(c)、8(d)分別表示在1.0 倍噪聲下利用本文算法和旋轉(zhuǎn)圖像算法的粗配準結(jié)果．從圖中可以看到，本文的配準算法具有很好的抗噪性，即使在1.0 倍噪聲下也能達到不錯的配準效果，遠遠優(yōu)越于旋轉(zhuǎn)圖像算法．

圖6 不同角度的兩點云視圖Fig.6 Two views of head with different angles

圖7 人頭點云的配準過程Fig.7 Process of registration for the head

圖8 不同高斯噪聲下的點云配準結(jié)果Fig.8 Results of registration for Stanford bunny with different noise levels

為了進一步說明該算法在尋找對應(yīng)點上的時間復(fù)雜度，以體積分(integral volume descriptors)算法[4]作對比．針對點云數(shù)據(jù)Bunny、Chicken、Dragon[15]，以Bunny000、Chicken_view1 和Dragonstandright_0作為源點集，Bunny045、Chicken_view2 和Dragonstandright_24 作為目標點集．在相同硬件條件下，在源點集中選取相同的100 個特征點，然后分別利用曲率圖和體積分算子在對應(yīng)的目標點集中尋找對應(yīng)點，每個特征點的對應(yīng)點數(shù)量設(shè)定為20，實驗結(jié)果如表1所示．由于本文的匹配算法是通過兩步法的思想尋找對應(yīng)點，只對很少的部分點集進行了高維曲率圖計算，所以耗費的時間也大大少于體積分算子．

表1 曲率圖和體積分算子找到對應(yīng)點的耗時對比Tab.1 Comparison of the point matching time for some models

6 結(jié) 語

本文在點幾何表示的三維點云數(shù)據(jù)中，構(gòu)造了一維、高維曲率圖．并利用曲率圖分兩步確立了具有重疊部分的三維點云數(shù)據(jù)間的對應(yīng)點．先利用一維曲率圖粗略找對應(yīng)點，然后對找出的點進行高維曲率圖運算，剔除大量錯誤對應(yīng)點．在此基礎(chǔ)上，通過改進的貪婪算法實現(xiàn)了對應(yīng)點之間的精確匹配．找到匹配點對后，通過SVD 求解剛性變換，將兩片點云數(shù)據(jù)粗略配準到ICP 算法的收斂域，最后借助經(jīng)典ICP算法實現(xiàn)了點云數(shù)據(jù)的精確配準．通過各種配準實驗，體現(xiàn)了該算法不但具有很強的抗噪性，在計算時間上也大大優(yōu)于其他基于特征點的配準算法．

致 謝：

特別感謝卡內(nèi)基美隆大學(xué)機器人研究所的Daniel Huber 教授提供了利用旋轉(zhuǎn)圖像進行三維點云數(shù)據(jù)配準的可執(zhí)行代碼．同時感謝斯坦福大學(xué)的計算機圖形學(xué)實驗室提供的三維點云數(shù)據(jù)．

2006.

［9］Johnson A，Hebert M. Using spin images for efficient object recognition in cluttered 3D scenes [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence，1999，21(5)：674-686.

［10］Frome A，Huber D，Kolluri R，et al. Recognizing objects in range data using regional point descriptors[J].Computer Vision-ECCV，2004，3023：224-237.

［11］Gatzke T，Grimm C，Garland M，et al. Curvature maps for local shape comparison[C]//IEEE International Conference on Shape Modeling and Applications. 2005：244-253.

［12］Mitra N J，Gelfand N，Pottmann H，et al. Guibas registration of point cloud data from a geometric optimization perspective [C]//Proc Symposium on Geometry Processing. 2004：23-32.

［13］Gray A. Modern Differential Geometry of Curves and Surfaces[M]. Boca Raton，F(xiàn)L，USA：CRC Press，1993.

［14］Huy Tho Ho，Gibbins Danny. Multi-scale feature extraction for 3D models using local surface curvature [C]//Digital Image Computing：Techniques and Applications.2008：16-23.

［15］Cazals F，Pouget M. Estimating differential quantities using polynomial fitting of osculating jets[J]. Computer Aided Geometric Design，2005，22(2)：121-146.

［16］Lange C，Polthier K. Anisotropic smoothing of point sets[J]. Compute Aided Geometric Design ， 2005 ，22(7)：680-692.

［17］Alexa M，Behr J，Cohen-Or D，et al. Point set surfaces[C]// Proceedings of IEEE Visualization. 2001：21-28.

［18］Stanford University. Stanford 3D Scanning Repository[EB/OL]. http：//www-graphics. stanford. edu/data/3Dscanrep/.

［19］Arun K S，Huang T S，Blostein S D. Least-squares fitting of two 3-D point sets[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence ， 1987 ，9(5)：698-700.