袁小翠,陳華偉

(1.南昌工程學院 江西省精密驅(qū)動與控制重點實驗室,南昌 330099; 2.貴州師范大學 機電工程學院,貴陽 550025)

0 概述

在逆向工程應用領(lǐng)域,重建的產(chǎn)品外形比較復雜,由多個連續(xù)性或種類不同的特征面按照不同的拼接條件構(gòu)成,如平面、球面、圓柱面和過渡面。對一些復雜曲面，單純用一種幾何模型去擬合得到的擬合曲面準確性差。在對產(chǎn)品的CAD模型重構(gòu)之前,一般需要將復雜產(chǎn)品的不連續(xù)曲面分割成若干個連續(xù)曲面,準確識別每個特征曲面的類型后對不同的曲面分別建模,再將這些曲面拼接以構(gòu)成完整產(chǎn)品。因此,點云模型特征面分割與識別是曲面準確重建的基礎。

散亂點云數(shù)據(jù)分割主要分為基于三角網(wǎng)格模型和基于散亂點云的數(shù)據(jù)分割方法[1]?；谌蔷W(wǎng)格模型分割需要對點云進行曲面重建,然而三維點云曲面重建比較耗時,基于點云的數(shù)據(jù)分割直接對點云模型進行分割,不需要對曲面進行重建,吸引了越來越多的關(guān)注。目前,點云數(shù)據(jù)分割方法可以歸納為3類,分別是基于邊界檢測、區(qū)域增長和聚類的分割算法[2]。基于邊界的分割方法主要通過查找點云模型的不同特征面的過渡邊界,對邊界線擬合從而將一個復雜曲面劃分為多個獨立的特征面。文獻[3]通過提取點云的幾何屬性值如法向量和曲率等,將點云的法向量映射到高斯球上,不同連續(xù)曲面映射到高斯球上的不同區(qū)域,同一片連續(xù)曲面映射到高斯球上的同一片區(qū)域,再根據(jù)邊界線的區(qū)域分割出如二次曲面、拉伸面和直紋面等特征曲面?；趨^(qū)域生長的模型分割方法的關(guān)鍵在于種子點的選取和生長規(guī)則,通過選定種子點,按照一定的生長規(guī)則將滿足條件點合并在同一區(qū)域中。比如:文獻[4]采用區(qū)域生成算法對LiDAR點云提取建筑物和植被區(qū)域,且取得了滿意的結(jié)果;文獻[5]計算每個采樣點的法向量,然后根據(jù)點云法向量檢測出點云模型中平滑區(qū)域,利用區(qū)域增長法對點云模型進行分割，與基于區(qū)域增長的分割算法相比,聚類方法并不需要初始化種子點和設計生長規(guī)則,其關(guān)鍵在于確定點云聚類模式;文獻[6]采用K-means聚類法對三維點云模型進行分割,但是該方法的缺點在于聚類中心初始化和K值的確定;文獻[7]基于3D活動輪廓模型對點云進行分割,通過提取模型的輪廓,將點云數(shù)據(jù)模型分割成若個干不同的區(qū)域;文獻[8]利用譜聚類算法對點云數(shù)據(jù)分塊,根據(jù)歸一化的非對稱Laplacian矩陣構(gòu)造譜聚類空間,通過特征精簡,在更低維的空間中進行點云分割。

目前雖然已有不少點云數(shù)據(jù)分割的方法,但許多點云處理專業(yè)軟件，如Geomagic和Imageware,都不能實現(xiàn)完全點云自動分割,只能通過手動和算法相結(jié)合的方法實現(xiàn)數(shù)據(jù)分塊,模型分割和特征面識別速度和準確性有待改進。為此,本文在準確估算點云法矢的基礎上采用高次曲面估算點云曲率,并根據(jù)平均曲率和高斯曲率將點云劃分為8種類型,通過引入二維圖像處理的連通區(qū)域標記法對同種類型的點云數(shù)據(jù)進行分割,并根據(jù)概率統(tǒng)計法和點云曲率特性,識別連通區(qū)域特征面點云的所屬特征曲面類型。

1 點云模型特征面

1.1 微分信息估算

1.1.1 法向量估算

法向量是點云的重要屬性之一,點云法向量的有效估計是點云分割的基礎。點云法向估計方法可分為局部鄰域擬合法和Voronoi/Dalaunay方法兩類[9]。局部鄰域擬合法在點云法向量估算中應用比較廣泛,由文獻[10]提出,稱之為主成分分析法,該方法能快速有效地估算光滑曲面的法向量,但當曲面包含尖銳特征時估算的法向量不準確。一般來說產(chǎn)品的CAD模型是由多個基面組成的復雜模型,在多個基面交界的區(qū)域(特征區(qū)域)主成分分析方法估算的法矢被平滑,為了更準確地估算點云法向量,研究人員提出了許多改進方法[11-13]。其中,文獻[11]所提方法能快速、準確地估算包含尖銳特征復雜模型的點云法矢,且參數(shù)具有自適應性,無需人工設置參數(shù)。本文采用文獻[11]的方法進行法矢估計。

給定點集X={x1,x2,…,xN},其中,N為點云總數(shù),任意一點xi的最近k鄰域表示為Nb(xi),文獻[11]對點xi的k鄰域擬合的平面表示為:

(1)

wr(ri)=exp(-(ri/σr)2),wn(ni)=exp(-‖ni-nj‖/σn2)和wn(xi)=exp(-‖xi-xj‖/σd2)分別為鄰域點對當前點的殘差、法向偏差和距離高斯權(quán)重,σd、σr和σn分別是距離、殘差和法矢偏差帶寬,用來控制各鄰域點對當前點作用的大小。式(1)可以轉(zhuǎn)化為對式(2)中半正定協(xié)方差矩陣C進行特征值分解。

(2)

C可以分解為3個特征向量v2、v1和v0,3個特征向量對應的特征值分別是λ2,λ1,λ0,其中,λ2≥λ1≥λ0。最小特征值對應的特征向量為平面的法向量,即點xi的法向量ni=v0。

通過分解式(2)得到的法向量方向不一致,為了后續(xù)計算需要,把各點的法向量方向調(diào)整為一致方向。本文采用最小成本路徑傳播方法進行法矢方向調(diào)整[10]。圖1是2種不同模型法矢估算和方向調(diào)整的結(jié)果,參數(shù)t=5,圖中估算的法矢與相應局部面近似垂直。

圖1 2種不同模型散亂點云法向量估算結(jié)果

1.1.2 離散點云曲率估算

通過法矢計算和方向調(diào)整,獲得了比較準確的法矢。以法矢為鄰域曲面的局部支撐,即可通過鄰域構(gòu)造離散點的局部參數(shù)化曲面,從而對離散點估算其曲面特性參數(shù)。

設待擬合鄰域曲面為z(u,v),通過Tylor級數(shù)展開,去掉高階無窮小,即可獲得近似曲面表達式:

z(u,v)=TA,d(u,v)+O(d+1)

(3)

本文構(gòu)造4次曲面估算點云曲率,曲面Tylor表達式為:

z(u,v) ≈TA,4(u,v) =

A0,0+ (A1,0u+A0,1v) +

(A2,0u2+2A1,1uv+A0,2v2)/2 +

(A3,0u3+3A2,1u2v+3A1,2uv2+A0,3v3)/6+

(A4,0u4+4A3,1u3v+6A2,2u2v2+

4A1,3uv3+A0,4v4)/24

(4)

引入新的矢量表達式:

P= (1,u,v,u2,uv,v2,u3,u2v,uv2,v3,u4,u3v,u2v2,uv3,v4)

則4次參數(shù)曲面可表示為如式(5)所示的矩陣形式。

PTQ=Z

(5)

其中,矩陣P中的u,v值為局部坐標系下的(u,v)參數(shù),Q為待定系數(shù)Ak-i,i構(gòu)成的矩陣,Z為局部坐標系下的z坐標值。將點及其鄰域點的坐標(ui,vi,zi)代入矩陣方程即可求出所有待定系數(shù)Ak-i,i。

求解曲面基本參數(shù)E、F、G和L、M、N,并構(gòu)造Weingarten曲率矩陣W[14]:

(6)

Weingarten矩陣的特征值就是主曲率k1、k2,特征向量就是主曲率方向,由k1、k2直接求高斯曲率K和平均曲率H。

曲率推導過程中,平均曲率H與法矢向量正相關(guān),因而其符號與法矢方向正相關(guān)。法矢實際是由曲率主方向叉乘而來,如果曲率在法矢方向調(diào)整之前就已計算,則在法矢方向調(diào)整后,曲率方向和符號也要進行相應調(diào)整。如圖2所示,法矢n反向后,對應最小主曲率k1和最大主曲率k2應反向互換。由H=(k1+k2)/2關(guān)系式可知,平均曲率H也要變號,由關(guān)系式K=k1×k2知,在k1和k2均變號的情況下,高斯曲率K不變號。

圖2 曲率方向調(diào)整示意圖

采用二次曲面法計算各點曲率,曲面坐標系在點的k鄰域局部構(gòu)建,各點之間的坐標構(gòu)建相對獨立,因而計算結(jié)果中的曲率值標準并不統(tǒng)一。根據(jù)點云分割需要,需對各曲率參數(shù)H、K、k1、k2進行歸一化。同一物體的點云分割應采用同一判斷標準,采用曲率法對點所屬曲面類型進行判斷,還涉及曲率零值和符號的確定問題,因而必須事先對所有點的曲率做全局歸一化處理,并保證曲率在零值附近的分布規(guī)律不變。對此,本文采用零值對齊的歸一化方法,從而保證了曲率數(shù)值符號的不變性。圖3所示為4種歸一化情況,其中M和m分別為曲率最大值和最小值,y為任一點曲率,不同情況的曲率歸一化方法如下:

1)如果M>m>0,則置最大曲率為1,最小曲率為m/M,其他點的曲率按y=y/M計算。

2)如果M>0,m<0,M>abs(m),處理方法同第1種情況。

3)如果M>0,m<0,M

4)如果m

圖3 不同情況曲率歸一化示意圖

1.2 特征面的分割與識別

1.2.1 點云初始聚類

基于離散點的高斯曲率K和平均曲率H,可以將點附近的曲面形狀劃分為峰、阱、脊、谷、鞍形脊、鞍形谷、平面和極小曲面8種類型,如表1所示。據(jù)此可將點云劃分為8種類型,實現(xiàn)點云的初始聚類。

表1 曲面分類

高斯曲率由最小和最大主曲率乘積而得,會引入高階誤差,一般采用平均曲率作為上述劃分的依據(jù)。在曲率歸一化處理的前提下,根據(jù)表1的判斷條件即可確定點所屬面域特征。此時,曲率零值的設定將對曲面類型劃分結(jié)果產(chǎn)生較大影響,本文處理方法是預設極小值e作為零值。經(jīng)實驗測定,推薦平均曲率H的零值閾值為[0.001,0.005],一般可預設eH=0.005,高斯曲率K的零值可按eK=eH/4設定,并將這2個閾值參數(shù)作為用戶輸入?yún)?shù),便于調(diào)整。

1.2.2 點云連通區(qū)域分割

點云的初始聚類只是將所有離散點歸屬至8種曲率類型,但是同種類型的點云可能分屬于不同的模型表面,如圖4所示,黑色點集和白色點集都屬于平面類型點,但是分別位于不同的面片,2個平面點集之間并不連通。因而,需要對同一種類型的點判斷其連通性,將不連通的區(qū)域分割出來,實現(xiàn)點云的特征面分割。

圖4 同種類型點云不連通情況示意圖

從初始點云的角度看,2個分離區(qū)域必然存在鄰域不連續(xù)的點,同一連通區(qū)域(局部相鄰)內(nèi)的點必定是全局相鄰的。本文將二維圖像處理中連通區(qū)域標記理論引入至三維點云區(qū)域分割[15],對三維點云連通性進行判斷。

連通區(qū)域識別步驟如下:

1)當前點xi加入當前連通區(qū)域Ccidx_cur。

2)采用上述判斷條件,對xi的局部鄰域點xj進行全局連通性判斷。

3)如果兩者全局相鄰,則判斷xj是否已追加至其他連通區(qū)域Ccidx_merge,如果是則合并Ccidx_cur和Ccidx_merge。

4)否則,將xj加入當前連通區(qū)域Ccidx_cur。

經(jīng)過連通區(qū)域分割后同種類型的面被分割成不同的塊,需對各塊點云進行特征面識別。

1.2.3 分塊點云的曲面類型識別

CAD模型一般由規(guī)則型面組成,例如平面、圓柱面等規(guī)則特征曲面用于安裝面、孔、倒角等加工特征,逆向工程中的特征面識別也主要是針對此類規(guī)則面。經(jīng)過連通區(qū)域分析,已經(jīng)獲得分片點云,可以進一步展開典型造型面(平面、球面、圓柱面等)的類型判斷。從統(tǒng)計角度上看,經(jīng)分塊后的點集同屬一個曲面類型。因此,本文采用數(shù)學統(tǒng)計法對分塊點云所屬曲面類型進行判斷,對區(qū)域點集中的點逐一判斷其所屬曲面類型,統(tǒng)計該區(qū)域中屬于不同曲面類型的點的數(shù)目,將最大數(shù)目對應的曲面類型作為整個區(qū)域的曲面類型。

設xj為連通點云數(shù)據(jù)塊X2∈X1中的任一點,Ti(i=0,1,2)對應平面、球面、圓柱面3類基本造型面,根據(jù)曲率條件,將點xj存入相應的類型數(shù)組Ti。3種特征面的判斷條件如下:

1)平面:xj的最小曲率和最大曲率均等于0,即k1=k2=0,亦即K=H=0。如果過濾掉曲率突變點,則法矢夾角法更為直觀,可優(yōu)先識別平面。

2)圓柱面:xj的最小曲率等于0,最大曲率為常數(shù),即k1=0,k2=C,亦即K=0,H=C。

3)球面:xj的最小曲率和最大曲率相等,即k1=k2=C=1/R,亦即K=1/R2,H=1/R(R為球半徑)。

考慮到計算和浮點誤差,在實際程序處理中,需對零值和兩值是否相等的判斷做特殊處理,預設極小值e,對零值和兩值是否相等分別使用x

基于概率統(tǒng)計的曲面類型判斷算法步驟如下:

1)對數(shù)據(jù)塊X2中各點曲率做統(tǒng)計分析,獲得曲率參數(shù)H、K、k1、k2的平均值和方差。

2)遍歷X2中的所有點,按照點云所屬曲面類型判斷方法,將點劃分至Ti。

3)統(tǒng)計Ti中的點數(shù)Ni,比值Ni/N即為當前區(qū)域?qū)儆诿骖愋蚷的概率pi。

4)按標準正態(tài)分布3σ原則預設概率閾值pT(pT=0.68)。如果pi>pT,則認定當前聚類屬于曲面類型i,返回曲面類型參數(shù)i。

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

使用VS2008和OpenGL開發(fā)平臺,在Windows 10下分步實現(xiàn)了上述算法,并使用多種模型驗證算法的有效性。實驗所用的計算機的配置為CPU Intel Core i3 3.4 GHz,內(nèi)存8 GB。

實驗結(jié)果如表2、圖5～圖9所示。其中,在表2中,N為模型點云數(shù)量,面特征列第一個數(shù)字表示特征面類型(1-平面,2-圓柱面,3-球面),括號內(nèi)的數(shù)字為該類型面的數(shù)目,r、eH分別表示曲率突變比例、平均曲率零值閾值。表中識別率表示算法對模型特征面所屬曲面類型識別的正確率,是當前模型中正確識別的面片與總面片數(shù)的比值。圖5～圖9分別是模型1～模型5幾種不同模型的分割和特征面識別結(jié)果。在計算參數(shù)中,模型2選用了較低的r值,是為了避免將部分圓柱面點納入突變點,模型4選用較高的r值,以將前端溝槽處的棱線納入突變點。此外,各模型的eH取值范圍在[0.001,0.005]都能取得較好的分割效果。在接下來的點云識別計算中,規(guī)則模型1～模型3無過渡面,型面曲面類型全部正確識別,識別率達到100%,模型4除了過渡面之外全部正確識別,總體識別率達到95%,模型5除了過渡面未識別之外,其他面域識別均是正確的,主要型面的曲面類型識別率大于90%。面面相交實驗顯示,模型1～模型4的交線全部正確提取,模型5提取了主要平面與圓柱,以及圓柱與圓柱交線(相貫線)。此外,從計算效率來看,點云數(shù)量在50 000以下的模型各節(jié)點耗時都在0.65 s以內(nèi)。

表2 參數(shù)設置及實驗結(jié)果

圖5 模型1點云分割與特征面識別結(jié)果

圖6 模型2點云分割與特征面識別結(jié)果

圖7 模型3點云分割與特征面識別結(jié)果

圖8 模型4點云分割與特征面識別結(jié)果

圖9 模型5點云分割與特征面識別結(jié)果

3 結(jié)束語

本文提出基于連通區(qū)域標記和統(tǒng)計法的散亂點云特征面分割與識別方法,在準確估算點云法矢的基礎上,以估算的法矢作為4次曲面曲率局部坐標系的坐標軸。根據(jù)點云的曲率特性分析曲面的型面并進行點云預分割,引入二維圖像處理中的連通區(qū)域標記法對點云連通區(qū)域分割,并采用統(tǒng)計法對點云所屬特征類型面進行判斷。實驗結(jié)果表明，該方法能夠取得較好的點云分割和識別效果。