基于法向量區域聚類分割的點云特征線提取

史紅霞 王建民

太原理工大學礦業工程學院，太原，030024

0 引言

隨著三維掃描技術的飛速發展，點云模型已成為幾何模型的通用表示形式，并廣泛應用于計算機可視化、逆向工程、三維重建、文物保護等領域[1-5]。特征是反映點云模型幾何外觀的最小基元，可精準地表達三維模型表面特性和幾何形狀[1]，因此，特征提取已成為點云數據處理中的一個基礎環節，如何有效準確地提取模型特征是當前點云數據處理的關鍵。

為解決點云模型特征提取問題，部分學者通過點云局部鄰域特性及幾何特征提取模型特征。幾何特征主要指點間距離特征、角度特征及點云法矢、曲率等信息。文獻[6]采用空間距離和角度特性提取車輛行駛道路和礦井道路邊界，具有一定的準確性，但未考慮道路障礙物的影響。文獻[7-8]利用點云法向量信息識別特征點，但均未考慮法向量估算的準確度。文獻[9]通過對各點局部鄰域協方差矩陣特征值的區域增長生成分割帶，在各分割區域內根據主方向提取特征點。文獻[10]引入多尺度分類算子分析局部鄰域點，計算不同尺度下采樣點成為特征點的概率，但多尺度使其計算效率下降。文獻[11]根據曲率突變點隱含特征點，設置突變點比率選取備選特征點，但由于比率具有不確定性，故可能會造成細節特征的缺失。

基于點云局部鄰域特性及幾何特征的方法大多依賴全局閾值進行特征點識別，全局固定閾值的設置易造成特征點的誤判或漏判，為此，有些學者運用聚類分割或局部曲面網格重建的方法識別特征點。例如，文獻[12]對各點局部鄰域三角區域的單位法向量進行高斯映射，對映射點進行聚類分析獲取特征點；文獻[13] 通過對遺跡碎片表面法向量的聚類實現碎片分割，結合法向量和碎片表面粗糙度提取邊界線；文獻[14]對已標記的潛在谷脊點構建局部幾何特征信息三角網格，然后再進行谷脊點識別，該方法計算復雜度較高。

上述兩類方法均不易準確完整地提取所有模型特征點，尤其是模型過渡特征。為此，有些學者引入數學形態學[15]、投影映射[16]、切片技術[17]、神經網絡[18]、張量投票理論[19]等方法，使某些模型特征提取效果達到最優。

本文針對現有模型過渡線及細節特征線提取不完整問題，基于區域分割思想進行特征線提取，在綜合利用主成分分析(principal component analysis,PCA)法、螢火蟲算法[20](firefly algorithm,FA)和模糊C均值(fuzzy C-mean, FCM)聚類算法[21]對點云模型進行聚類分割的基礎上，將各分割塊邊界點作為備選點集識別特征點，采用三次B樣條法擬合生成特征線[22]。

1 FCM聚類算法

已知n個l維空間數據集S={s1,s2,…,sj,…,sn}?Rl，將數據集S劃分為c(2≤c≤n)個子集，即聚類數為c，各子集聚類中心C={C1,C2,…,Ci,…,Cc}。最佳聚類準則是極小化目標函數值，使具有相似特征的數據點聚為一類。目標函數Jm的一般定義如下：

(1)

dij=‖sj-Ci‖Uij∈[0,1]

式中，dij為空間數據集S中第j個樣本點sj與第i個聚類中心Ci的歐氏距離；m為模糊加權指數；Uij為第j個樣本點屬于第i個聚類中心的隸屬度。

目標函數Jm取極小值的約束條件是

最小化Jm過程中，可運用下式計算隸屬度和更新聚類中心：

(2)

(3)

式中，dgj為空間數據集S中第j個樣本數據點sj與第g個聚類中心Cg的歐氏距離。

FCM聚類算法的基本思路就是通過不斷迭代更新隸屬度Uij和聚類中心Ci，使目標函數Jm最小，然后依據Jm最小時的隸屬度矩陣劃分數據集S。

運用常規的FCM聚類算法對點云模型進行分割進而提取特征線時，需要考慮兩個問題。首先，若直接利用各點空間坐標信息進行聚類劃分，則所得模型分割線將與模型特征線無直接聯系，不易于后續特征線提取，為此，引入自適應鄰域的主成分分析(PCA)法[23]求取模型各點法向，以各點法向量方向的接近程度為聚類特征實現點云模型的有效分割。其次，由于初始聚類中心選擇的隨機性以及梯度下降式的最優解搜索法，致使算法容易陷入局部最優，為解決這一問題，在FCM聚類算法中引入FA算法，利用具有全局搜索能力的FA算法為FCM聚類算法選取最優初始聚類中心，有效解決其對初始值敏感及過早收斂的問題。

2 基于法向量估計的FA-FCM組合聚類算法

法向量是散亂點云模型的一個重要屬性，可用來描述模型特征信息。點云模型中任意采樣點局部區域較為平坦時，鄰域內各點法向量變化平緩；相反，采樣點局部區域起伏變化時，各點法向量變化較大，因此，可估算模型法向量，進而提取模型特征點。

2.1 法向量估計

目前，在常見的點云法向量估計算法中，PCA方法因其算法簡單高效、穩定性強而被廣泛使用。PCA算法實質是對點云數據局部鄰域的主成分分析，然而，固定鄰域尺寸的選擇致使其估算結果不精確。以二維情況(三維類似)為例，說明因鄰域值選擇不合理而造成的法向量估算偏差問題，如圖1所示。

圖1 局部鄰域值大小對法向量估算的影響Fig.1 Effect of neighborhood size on normal estimation

圖1中圓形區域代表鄰域范圍，虛直線代表虛線圓形鄰域內數據擬合形成的直線，虛線箭頭表示根據虛線圓形鄰域范圍計算所得的法向量，實直線代表實線圓形鄰域內數據擬合形成的直線，實線箭頭表示根據實線圓形鄰域范圍計算所得的法向量。圖中點s1為尖銳特征點，應選擇較小的鄰域值擬合估算該點法向量，故圖中s1點處虛線法向量因鄰域值過大而造成偏差，而s2點所在區域較為平坦，只有局部略微凸起，應選擇較大的鄰域值擬合估算該點法向量，故圖中s2點處實線法向量因鄰域值過小而造成偏差。為解決這一問題，采用點云局部鄰域幾何維度特征最優鄰域選擇法[23]查找各點最優鄰域值，然后結合PCA算法估算各點法向。

對點云樣本數據集而言，空間維度l的值為3，此時，點云樣本數據集S={s1,s2,…,sj,…,sn}?R3，j=1,2，…,n,n為樣本總點數，取數據集S中任意樣本點s，其k近鄰域點集Ns={s1,s2,…,sb, …,sk}，b=1,2，…,k，點s的協方差矩陣可表示為

(4)

由于協方差矩陣M為對稱正定陣，因此其特征值λ0、λ1、λ2具有非負性，同時，可用三個特征值描述樣本點局部鄰域維度特性，其特性如下：

(1)當λ0?λ1,λ2>0時，局部鄰域呈線性特征分布。

(2)當λ0,λ1?λ2>0時，局部鄰域呈平面特征分布，此時最小特征值λ2所對應的特征向量即為鄰域點法向量。

(3)當λ0≈λ1≈λ2時，局部鄰域呈三維曲面特征分布。

根據各樣本點局部鄰域維度特性可定義其特征模型：

(5)

式中,a1、a2、a3分別表示當前點局部鄰域為線性(一維)、平面(二維)和曲面(三維)特征。

基于信息熵理論，依據上述維度特征模型定義點s的近鄰點數為k時信息熵函數：

(6)

結合上述基于信息熵理論的自適應鄰域選擇法，采用PCA方法進行法向量估算，步驟如下：

(1)對于點云樣本數據集S中任意采樣點s，設置鄰域搜索上下限值分別為kmax和kmin(kmax,kmin均表示鄰域點數)，搜索步長為1，初始化鄰域值k←kmin。

(3)設置搜索鄰域值k←k+1，當k≤kmax時，返回步驟(2)；否則執行下一步。

2.2 FA-FCM組合聚類算法

為解決傳統FCM聚類算法因隨機性的初始聚類中心選擇法造成的算法過早局部收斂問題，引入FA，將每只螢火蟲的位置作為聚類中心，通過比較螢火蟲亮度值進行算法迭代循環，從而為FCM聚類算法尋找最優初始聚類中心。在已求解模型法向的基礎上，介紹引入FA算法后FCM聚類算法的初始聚類中心選擇方法。

通常情況下，模糊加權指數m為2，對點云模型法向量進行聚類劃分時，距離d(vj,Ci)和隸屬度Uij的表達如下：

(7)

(8)

則改進后的模糊聚類目標函數為

(9)

FA算法是尋找亮度最亮的螢火蟲個體，其中亮度和吸引度是該算法的兩個重要概念[24]，而FCM聚類算法旨在尋找目標函數最小時的聚類中心，因此，螢火蟲亮度值可定義為

(10)

螢火蟲吸引度為

β=β0exp(-γr2)

(11)

式中，r為兩螢火蟲間距離;β0為最大吸引度，即r=0時的吸引度；γ為介質對光的吸收因子。

兩螢火蟲A、B間距離為

(12)

i=1,2，…,c

式中，CAi、CBi分別為螢火蟲A和螢火蟲B的第i個聚類中心。

螢火蟲A受亮度較高的螢火蟲B吸引而向B移動，A移動后的位置由下式決定：

CAi=CAi+β(CBi-CAi)+α(rand-(0.5))

(13)

式中，α為步長因子，α∈[0,1]；rand是由3個隨機因子構成的向量，在區間[0,1]上屬于均勻分布。

FA-FCM組合聚類算法實現步驟如下：

(1)初始化FA算法參數β0,γ,α,c,T1，N(T1為FA算法最大迭代次數，N為螢火蟲數)及螢火蟲位置，設置FCM算法收斂閾值ε(ε>0)。

(2)將初始化的螢火蟲位置作為初始聚類中心，聯合式(7)、式(8)計算隸屬度矩陣U(N×c)×(n)和距離d(vj,Ci)，聯合式(9)、(10)計算N只螢火蟲亮度值矩陣Is(N×1)。

(3)比較各螢火蟲亮度，聯合式(11)、式(12)、式(13)更新各螢火蟲位置。

(4)判斷算法迭代次數是否達到T1，若是，則按式(7)～式(10)計算各螢火蟲亮度，取最亮螢火蟲位置為最優初始聚類中心；否則，返回步驟(2)，繼續迭代。

(5)確定初始聚類中心后，使用FCM聚類算法進行聚類分割，令FCM算法迭代次數t=0。

(7)按式(3)更新聚類中心位置，此時，t←t+1。

2.3 FA-FCM聚類有效性檢驗

為驗證聚類方法的有效性,本文選取的評價指標為類內緊致度Var和類間離散度Sep[25]。

根據點云模型隸屬度矩陣及聚類中心，類內緊致度為

式中，U為隸屬度值Uij構成的隸屬度矩陣，i=1,2,…,c,j=1,2,…,n；n(i)為第i個聚類中心樣本數。

Var表示樣本類內集中程度，Var越小，類內越緊致。

類間離散度為

0≤Sep(c,U)≤1

式中，FX、FY為兩個模糊集合；vj表示樣本點向量，j=1,2,…,n;min(UXj，UYj)表示該樣本點分別屬于FX和FY隸屬度的最小值；S(FX,FY)表示兩集合間相似度；當FX=FY時，Sep(c,U)=0。

Sep表示兩聚類集群離散度，Sep越大，各類間越分散。

為綜合評價聚類結果，將Var與Sep結合，定義綜合評價指標為

VS越小，聚類效果越好。

3 特征點提取

3.1 邊界點剔除與合并

應用FA-FCM聚類算法對點云模型進行分割后，需要對各分割塊邊界點進行特征點提取，在特征點提取過程中，所得分割塊的邊界點有兩種情況需要處理。首先，過分割問題的存在容易致使模型平坦區域出現分割線，顯然該分割線為非特征線，該分割線處的點為明顯非特征點，表示該分割線的兩組邊界點相鄰且共面，應從邊界點集中剔除；其次，若分割線為特征線，則表示該分割線的兩組邊界點相鄰但不共面，應對其進行合并處理。

本文在利用鄰近投影點相鄰向量夾角的方法[26]提取各分割塊邊界點的基礎上，對邊界點進行優化處理，剔除相鄰且共面的邊界點，合并相鄰但不共面的邊界點。

為直觀表達，以圖2所示圓柱體與長方體的組合模型為例進行說明。圖2a所示為該組合模型經FA-FCM算法聚類分割后所得結果，由該圖可知，為有效分割長方體四個側面，圓柱體側面也會被劃分為四組，存在過分割現象，導致圓柱體側面邊界點集l1、m1為非特征點集，長方體相鄰兩側面邊界點集l2、m2均為兩側面交線處特征點，如圖2b所示。為此，以各分割塊為處理單元，通過分析比較兩分割塊邊界點間距離和法向夾角，利用設置閾值的方法剔除相鄰且共面的邊界點集(如l1、m1)，將相鄰但不共面的邊界點集(如l2、m2)合并，方法如下。

(a)區域聚類分割結果 (b)邊界點提取結果圖2 邊界點優化準則Fig.2 Optimization rule of boundary points

模型分割后，定義分割塊邊界點集為

F={fi}i=1,2,…,c

式中，fi為第i個分割塊邊界點集。

給定任意兩組邊界點集fw1,fw2(w1,w2=1,2，…,c，w2≠w2)，依據點云密度設置距離閾值D，以及依據邊界點集fw1和fw2所對應分割塊的聚類中心間夾角設置法向夾角閾值T，分別取fw1和fw2中任一點su1(u1=1,2，…,n1)和su2(u2=1,2，…,n2)，n1、n2為各點集總點數,求su1、su2兩點間坐標距離d、法向夾角φ。判定準則如下：

(1)若d≥D，則點su1、su2不做任何處理。

(2)若d

(3)若d

3.2 候選特征點平滑

各分割塊邊界點集經優化處理后所得點集為無序點，稱該無序狀態點集為候選特征點集，為便于曲線擬合，需先對候選特征點進行排序，然后再擬合成線。然而，由于點云數據在獲取過程中不可避免地受外界環境、掃描設備等因素影響，導致候選特征點集中存在噪點，若直接對排序后的候選特征點進行曲線擬合，則噪點的存在會嚴重影響擬合精度，因此，以各點局部鄰域主軸方向[27]為基準對排序后的候選特征點進行平滑處理。

取排序后候選特征點集中的一點s0，確定點s0的鄰域值k0，對于鄰域內任意點su(u=1,2,…,k0)，記s0與su兩點連線方向s0su與平滑參考方向正向夾角為θu(u=1,2，…,k0)，如圖3b所示。平滑參考方向構造法如下：

(1)候選特征點s0為初始點時，無參考方向，取局部鄰域點PCA主軸方向D1為平滑參考方向，如圖3a所示。

(a)采樣點局部鄰域點集及主軸方向

(b)相鄰點方向夾角圖3 候選特征點平滑準則Fig.3 Smoothing rule of candidate feature points

(2)候選特征點s0為非初始點時，取s0的前一點s00與s0的連線方向s00s0為平滑參考方向。

確定參考方向后，平滑準則如下：

(1)距離準則。min(|s0su|)(|s0su|表示點s0與點su的距離，u=1,2，…,k0)，即查找k0鄰域內距s0最近的點。

(2)夾角準則。min(θu)(u=1,2，…,k0)，即查找s0和k0鄰域內任意點su(u=1,2，…,k0)的連線方向s0su與平滑參考方向正向間夾角最小的點。

如圖3b所示，取鄰域值k0為5，搜索s0的下一點，min(|s0su|)=|s0s5|，min(θu)=θ3。僅使用距離準則時，s5被確定為s0的下一點，此時s5與平滑參考方向角度偏差太大；僅使用夾角準則時，s3被確定為s0的下一點，此時距離跨度太長，直線拐點處特征點不易保留，容易出現過度平滑現象。因此，在對候選特征點進行平滑處理的過程中，應組合使用距離準則和夾角準則，設置夾角閾值T2，對滿足條件θu

根據候選特征點分布特性設置夾角閾值T2，多次調整閾值以達到最佳平滑效果，在實驗過程中，對直線和曲線類特征線設置不同的角度閾值T2，分別為T2=5°～10°和T2=20°～30°。

4 實驗結果與分析

4.1 主要參數設置

本文將基于法向量的FA-FCM聚類算法應用于點云模型特征線提取問題中，選取AIM@SHAHE Shape Repository公開數據庫中的Smooth、Fandisk、Vase模型以及Princeton ModelelNet數據庫中的m932-table模型作為實驗對象，驗證所提算法的可行性和優越性。

目前螢火蟲算法的參數設置還沒有明確合理的參考值，為此，首先通過分析比較文獻中的參數值來設置合理的參數范圍，然后測試不同的參數組合，確定該特征線提取算法的初始化參數值，設置見表1。

表1 參數列表Tab.1 The list of parameters

4.2 特征線提取結果的有效性及精度分析

為驗證本文算法的可行性，選取Smooth、Fandisk、Table和Vase模型進行特征線提取。首先估算各模型數據點法向，由FA-FCM混合算法對其進行聚類分割，各模型聚類有效性評價指標見表2。由表2可知，4個模型的Sep(0≤Sep≤1)值均接近于1，表明各模型整體劃分效果較好，但分析比較Var值與VS值可知，Vase模型所對應的Var和VS值均為最大，表明Vase模型的聚類效果較差，還需進一步優化。

表2 各模型聚類有效性指標Tab.2 Cluster validity indexes of each model

4個模型經FA-FCM算法聚類分割后，基于各分割塊邊界點提取特征點，對應該相交的特征線端點進行合并處理，然后擬合生成特征線，各模型特征線提取不同階段效果顯示如圖4～圖7所示。

(a)Smooth模型原始點云 (b)區域聚類結果

(c)特征點提取結果 (d)特征線提取結果圖4 Smooth模型特征線提取Fig.4 Feature line extraction from Smooth model

(a)Fandisk模型原始點云 (b)區域聚類結果

(c)特征線提取結果圖5 Fandisk模型特征線提取Fig.5 Feature line extraction from Fandisk model

Smooth模型結構簡單，主要由平面和曲面構成。圖4b是模型經FA-FCM區域聚類分割所得結果，通過相鄰兩分割塊交線提取該模型特征線，特征線以直線為主，特征信息較為突出，即為尖銳特征，圖4c為Smooth模型特征點識別結果，表明模型尖銳特征點能夠得到準確提取。

(a)Table模型原始點云 (b)區域聚類結果

(c)特征線提取結果圖6 Table模型特征線提取Fig.6 Feature line extraction from Table model

(a)Vase模型原始點云 (b)區域聚類結果

(c)特征線提取結果圖7 Vase模型特征線提取Fig.7 Feature line extraction from the Vase model

圖5為Fandisk模型特征線提取不同階段效果圖。Fandisk模型為機械輪盤模型，相比于Smooth模型結構較為復雜，既包含明顯棱線(即尖銳特征線)，又包含過渡特征線(如在模型底部相鄰曲面相交處的特征點曲率變化較小，法向量方向相近，特征信息不明顯，這類特征即為過渡特征)。圖5c為Fandisk模型特征點識別結果，表明Fandisk模型的尖銳特征得到了完整提取，同時也有效保留了過渡特征。

圖6為Table模型特征線提取不同階段效果圖。Table模型為餐桌模型，結構相對復雜，主要由平面、曲面以及球面構成，特征信息豐富，既包含尖銳特征，也包含細節特征，如模型桌腳及兩球面相交處均包含有細節特征。圖6c所示為Table模型特征點識別結果，表明本文所提方法在完整提取Table模型尖銳特征的同時也能有效識別模型細節特征。

圖7為Vase模型特征線提取不同階段效果圖。Vase模型為構造復雜的花瓶模型，由內外兩層曲面構成，特征線呈“之”字形，主要包括尖銳特征和過渡特征。與Fandisk模型相比，該模型曲面彎曲度較大，過渡特征的特征強度較弱且分布不均。圖7c所示為Vase模型特征點識別結果，表明該算法對Vase模型中強度較弱的過渡特征提取效果不佳，部分特征線斷裂。這與由表2得出的Vase模型聚類效果較差的結論相對應。

為對本文算法的特征線提取結果進行定量評價，分別計算各模型特征線的數量提取率和長度提取率。

模型數量提取率通過計算精確率Rpc和召回率Rrc[28]獲得，精確率表示被正確提取的特征線數目占被提取出特征線總數的比值，召回率表示被正確提取的特征線數目占模型特征線總數的比值。精確率、召回率的計算公式如下：

式中，Tp為被正確提取的特征線數目；Fp為被錯誤提取的特征線數目；Tz為模型特征線總數。

模型長度提取率通過計算模型中各特征線段長度提取率的平均值獲得，模型中某一特征線的長度提取率為已提取出的該特征線長度與其準確長度的比值。因為三維點云模型中特征曲線的準確長度不易求取，故本文利用“化曲為直”的思想，近似估算曲線長度。模型整體特征線長度提取率的計算公式如下：

式中，hq為第q條特征線的長度提取率(q=1,2,…,Tz)；lq為已提取的該特征線長度；Lq為該特征線準確長度；H為模型整體特征線長度提取率。

對4種模型的精確率、召回率和模型長度提取率H進行統計，結果如表3所示。分析表3可知，Smooth模型和Fandisk模型的特征線提取結果更加完整和準確，所提取特征線不存在誤判和漏判現象，由精確率、召回率及長度提取率H可判定Smooth和Fandisk的特征線基本得到了準確完整的提取，經細微調整后可投入實際應用中。Table模型特征線提取準確性較好，不存在誤判現象，但由其Rrc值可判定該模型特征線未得到完整提取。此外，4種模型相比之下，Vase模型特征線提取的完整性和準確性較低，由其Rpc值和Rrc值可判定該模型的特征線提取既存在誤判現象又存在漏判現象。但Table模型和Vase模型的召回率Rrc均達到90%，兩模型長度提取率H均可達85%，表明本文方法在進行相對復雜模型的特征線提取時也可取得較好成效，可進行局部二次特征提取，從而有效識別強度信息不明顯的特征點，提高模型特征線提取的完整性和準確性。

表3 各模型特征線數量提取率和長度提取率Tab.3 The extraction rate of the number and the length of feature lines of each model %

4.3 特征線提取方法對比

為進一步驗證算法的可行性及優越性，選擇Fandisk模型特征線提取結果與文獻[29]提取結果進行對比，結果如圖8所示。選擇Vase模型特征線提取結果與文獻[29]、[30]和[31]的提取結果進行對比，結果如圖9所示。

(a)本文方法特征線提取(b)文獻[29]特征線提取圖8 Fandisk模型特征線提取方法對比Fig.8 Comparison of the feature extraction methods used in Fandisk model

在圖8中，對比長方形框標定區域可以看出，文獻[29]的方法不易識別特征強度逐漸減弱的平滑過渡區特征，從而導致特征線提取不完整，而本文方法所提取特征線的完整度以及線條光滑度相對較好。

由圖9可知，本文算法能較好地識別內外兩層特征線，兩層線間不存在相互扭曲現象，曲線線條性更好，但不能完整識別特征強度較弱的特征點，致使曲線出現斷裂和不完整現象。

(a)本文方法特征線提取 (b)文獻[29]方法特征線提取

(c)文獻[30]方法特征線提取 (d)文獻[31]方法特征線提取

5 結論

本文基于自適應鄰域的FA-FCM算法實現點云模型的有效分割，通過構造邊界點集剔除與合并準則以及候選特征點集平滑準則從各分割塊邊界點集中識別特征點。結合實驗分析得出以下結論：

(1)相對于固定閾值的主成分分析(PCA)法，自適應鄰域的PCA法可較為精準地估計模型尖銳特征點處法向，再通過熒火蟲算法優化后的FCM算法可使模型分割效果達到最佳，從而提高各分割區邊界點中特征點比率。

(2)邊界點集剔除與合并準則可有效去除邊界點集中多數非特征點，從而提高算法效率；同時，結合候選特征點集平滑準則的應用可使所提算法能夠識別更多準確的模型特征，提高所擬合特征線的光滑度。

(3)通過4種實驗模型的特征線提取精度分析可知，Smooth模型和Fandisk模型的特征線可得到準確完整的提取，結構相對復雜的Table模型和Vase模型特征線的數量提取率可達90%，長度提取率也可達85%，其中，主要特征線均可得到完整提取。結果表明，本文算法簡單有效且具有一定的通用性，在準確提取點云模型顯著特征和細節特征的同時能盡可能多地保留模型過渡特征。