黃勝利 卜 昆 程云勇 周麗敏

西北工業(yè)大學(xué)現(xiàn)代設(shè)計與集成制造技術(shù)教育部重點實驗室,西安,710072

0 引言

在渦輪葉片模具型腔反變形優(yōu)化設(shè)計中,計算分析葉片測量數(shù)據(jù)相對于CAD模型的三維偏差至關(guān)重要。葉片測量數(shù)據(jù)的主要獲取方法為三坐標(biāo)測量和光學(xué)掃描測量。近年來,隨著激光技術(shù)的發(fā)展,光學(xué)掃描儀獲取數(shù)據(jù)的速度和精度都有顯著的提高,應(yīng)用越來越廣泛,但其獲得的海量數(shù)據(jù)給后續(xù)的數(shù)據(jù)處理、分析、計算帶來很大的不便。由于測量數(shù)據(jù)是在測量設(shè)備坐標(biāo)系下獲得的,故在對測量數(shù)據(jù)進行偏差分析前,需把測量數(shù)據(jù)與CAD模型的坐標(biāo)系統(tǒng)一起來,也就是實現(xiàn)測量數(shù)據(jù)與CAD模型的配準(zhǔn)定位。

目前應(yīng)用比較廣泛的配準(zhǔn)算法是Besl等[1]提出的迭代最近點(iterative closest point,ICP)算法。該方法通過最小二乘法構(gòu)造配準(zhǔn)的目標(biāo)函數(shù),具有收斂速度快、配準(zhǔn)精度高的優(yōu)點?；贗CP的改進算法[2]在含有自由曲面的工件配準(zhǔn)中得到了廣泛的應(yīng)用。ICP算法要求待配準(zhǔn)兩模型的初始位置要足夠接近,否則不能得到好的收斂結(jié)果,因此,使用ICP算法進行精確配準(zhǔn)前,需要對兩待配準(zhǔn)模型進行預(yù)配準(zhǔn)。武殿梁等[3]采用遺傳算法實現(xiàn)了曲面之間的預(yù)配準(zhǔn),此方法比較耗時。吳鋒等[4]基于輪廓的力矩主軸實現(xiàn)了醫(yī)學(xué)圖像的預(yù)配準(zhǔn),該方法的缺點是對數(shù)據(jù)的缺失比較敏感,需要整個物體的全部信息。劉晶[5]通過建立模型的包圍盒實現(xiàn)三維測量模型與CAD模型之間的預(yù)配準(zhǔn),該方法對測量模型的完整性要求較高,且包圍盒的建立有一定的難度。

光學(xué)測量的密集點云數(shù)據(jù)的配準(zhǔn)計算非常耗時,為了縮短配準(zhǔn)計算的時間,本文提出一種基于簡化模型曲率計算的配準(zhǔn)方法,首先對密集點云模型進行簡化,估算簡化模型的曲率,根據(jù)曲率提取對應(yīng)特征面;然后尋找匹配點對從而計算出變換矩陣,把求出的變換矩陣應(yīng)用于簡化前的原始點云模型以實現(xiàn)模型的預(yù)配準(zhǔn);再通過奇異值分解和最近點迭代相結(jié)合(singular value decompo-sition and iterative closest point,SVD-ICP)算法[6-7]進行精確配準(zhǔn);最后以某型渦輪葉片的蠟?zāi)槔?驗證方法的有效性和實用性。

1 點云模型簡化

由于光學(xué)測量獲取的數(shù)據(jù)量巨大,為了提高計算的速度需對點云模型進行適當(dāng)?shù)暮喕?。圖1所示分別為某型葉片及其蠟?zāi)５墓鈱W(xué)測量模型。目前,常用的點云簡化方法有包圍盒法、均勻網(wǎng)格法、聚類法、隨機采樣法、曲率采樣法、均勻采樣法等[8],它們各有優(yōu)缺點和不同的適用范圍。本文對點云模型進行簡化的主要目的是提高計算的速度,因此采用均勻采樣法。

圖1 光學(xué)測量密集點云模型

2 曲率估算

曲率是曲面的重要幾何特征,具有平移、旋轉(zhuǎn)和縮放不變性,可以作為曲面特征識別的重要依據(jù)。本文采用鄰域內(nèi)擬合二次參數(shù)曲面的方法估算散亂點云的曲率[9]。

2.1 散亂數(shù)據(jù)參數(shù)化

對原始點云進行空間劃分,計算點云中每個點pi的K 鄰域,K鄰域記為Knbhd(pi)。每個Knbhd(pi)可以擬合成一曲面且此曲面可以被一個局部基面逼近,該局部基面記為LBS(Knbhd(pi)),本文選擇點pi處的切平面為局部基面。切平面可通過各點的鄰近關(guān)系和最小二乘原理[10]求得。對于不同的模型,合理選取K的值(K一般取10～30),以保證鄰域內(nèi)凹凸一致,這樣得到的切平面才能更好地逼近原始點集,使投影點在局部基面的參數(shù)化更好地反映點云的參數(shù)化。

把Knbhd(pi)中的各點在相應(yīng)局部基面上的投影點記為PRO(Knbhd(pi)),投影點一般是隨機分布的。為了參數(shù)化PRO(Knbhd(pi)),首先需要確定 u、v兩個方向。在 PRO(Knbhd(pi))中 ,p′i是pi在局部基面上的投影點,求出距離p′i最遠的點p′,將連接p′i、p′兩點的直線方向作為u 方向,垂直于u方向的直線方向作為v方向。這里v方向可由u矢量和局部基面的法矢作向量積求得。然后將PRO(Knbhd(pi))中的每一點p′j(j=1,2,…,K)都與p′i相連,得到K 條有向線段,把這些有向線段分別與u矢量作點積,記為Dj,把得到的點積進行排序,最大值為Dmax,最小值為Dmin。由式(1)就可求得PRO(Knbhd(pi))中各點的u參數(shù)值:

其中,j=2,3,…,K+1且U1=0。同理可以求得各點的v參數(shù)值。把所求的PRO(Knbhd(pi))中各點的參數(shù)值作為Knbhd(pi)中對應(yīng)點的參數(shù),就實現(xiàn)了局部散亂數(shù)據(jù)的參數(shù)化。

2.2 擬合二次參數(shù)曲面

由2.1節(jié)中方法得到每個點pi的Knbhd(pi)的參數(shù)u、v,就可以建立二次參數(shù)曲面:

其中,A為3×3的系數(shù)矩陣。運用Yang[11]提出的二次參數(shù)曲面逼近法求得散亂數(shù)據(jù)點的二次參數(shù)曲面。首先引入新的矢量表達式:

則二次參數(shù)曲面可表示為

其中,系數(shù)a、b、c的下標(biāo)為式(2)中A的行號和列號。采用最小二乘法使Knbhd(pi)上各點到二次參數(shù)曲面的歐幾里德距離之和最小,可得

式中,X為鄰域Knbhd(pi)的矩陣。

2.3 曲率計算

由于計算得到的曲面法矢方向一般是不協(xié)調(diào)的,即并不是所有的法矢都指向曲面的同一側(cè),因此有必要對求得的法矢進行歸一化調(diào)整,具體方法參見文獻[9]。

根據(jù)曲面的曲率性質(zhì)可以得出高斯曲率:

平均曲率:

上述物理量可由曲面的第一、第二基本形式推出。

3 提取特征面

渦輪葉片主要由自由曲面和一些平面組成,而自由曲面的曲率特征沒有一定的規(guī)律性,因此,曲面的提取有一定的難度。本文中,將測量模型和CAD模型上某一組對應(yīng)平面作為尋找匹配點對的依據(jù)。高斯曲率I和平均曲率H均分別反映了曲面的凸凹形狀,I和H的不同組合表現(xiàn)的是不同的幾何形狀,表1給出了幾種常見的曲面形狀與I、H組合的關(guān)系。

表1 點的局部曲面類型

由表1可以看出,當(dāng)I、H 都為0時,曲面就是一平面,根據(jù)此條件可提取模型的某一平面為特征面,為尋找對應(yīng)匹配點作好準(zhǔn)備。

4 模型配準(zhǔn)

4.1 預(yù)配準(zhǔn)

為了分析鑄件相對于其CAD設(shè)計模型的彎扭變形,首先需要實現(xiàn)測量模型與CAD模型的配準(zhǔn)。運用SVD-ICP算法進行精確配準(zhǔn),要想獲得較好的效果,首先需要實現(xiàn)測量模型與CAD模型的預(yù)配準(zhǔn)。本文通過曲率提取特征面的方法實現(xiàn)預(yù)配準(zhǔn),具體方法如下:

(1)運用均勻采樣法對原始測量點云{P}進行簡化得到{P′}。

(2)計算簡化后點云{P′}和CAD模型的重心(分別記為O1、O2),平移兩模型,使兩模型的重心與坐標(biāo)原點O重合。

(3)運用第2節(jié)中所述方法估算簡化后測量模型和CAD模型的高斯曲率K與平均曲率H。

(4)根據(jù)所求得的K和H提取兩模型中相對應(yīng)的特征面,分別記為L1、L2,并求取兩特征面的

形心 P1、Q1。

(5)把P1、Q1沿所在特征面的法線方向移動一定距離 d(d值可適當(dāng)選取,本文實驗中 d=10mm),得點P2、Q2的計算公式:

式中,P1、P2、Q1、Q2為點 P1、P2、Q1、Q2的位置矢量;N1、N2分別為特征面的單位法矢。

(6)由三組對應(yīng)點{O,P1,P2}和{O,Q1,Q2}求取旋轉(zhuǎn)矩陣R和平移矩陣T,方法如下:以O(shè)點為坐標(biāo)系原點,由單位矢量組 e1 、e2 、e3 與e′1 、e′2 、e′3按右手法則構(gòu)建兩個局部坐標(biāo)系;設(shè)存在旋轉(zhuǎn)變換矩陣R與平移變換矩陣T使兩個局部坐標(biāo)系重合,于是有

由式(13)可求得旋轉(zhuǎn)變換矩陣:

待配準(zhǔn)模型中的任一點Pi變換到對應(yīng)點Qi的關(guān)系式為

把P1、Q1代入式(15)可得到平移變換公式:

(7)把求得的旋轉(zhuǎn)矩陣和平移矩陣應(yīng)用到簡化前的原始測量點云模型上,得到變換后的測量模型,再按向量OO2進行整體平移就實現(xiàn)了測量模型與CAD模型的預(yù)配準(zhǔn)。

4.2 精確配準(zhǔn)

實現(xiàn)了測量模型與CAD模型的預(yù)配準(zhǔn)后,就可以運用SVD-ICP算法進行精確配準(zhǔn)。記測量點云數(shù)據(jù)為{pi},與其對應(yīng)的CAD模型上的對應(yīng)點集為{p′i}。配準(zhǔn)目標(biāo)函數(shù)表示為

對配準(zhǔn)點集{pi}和{p′i}運用文獻[7]中奇異值分解法求取變換矩陣R和T,并將變換矩陣作用到測量點云數(shù)據(jù)上,迭代地進行這個操作,直到使目標(biāo)函數(shù)值滿足設(shè)定的閾值為止。這樣就實現(xiàn)了測量點云模型與CAD模型坐標(biāo)系的統(tǒng)一,即實現(xiàn)了測量模型到CAD模型的精確配準(zhǔn)。

5 實例驗證

為了驗證本文方法的有效性,以VC++6.0為平臺,以VTK可視化類庫為輔助工具編程實現(xiàn)了上述方法,并以某型渦輪葉片的蠟?zāi)槔M行了驗證。試驗中所用的測量點云模型是由ATOS流動式光學(xué)掃描儀獲取的,其數(shù)據(jù)點數(shù)為102 908。試驗所用的電腦配置為Intel(R)Core(TM)2 Duo CPU E6550@2.33GHz;內(nèi)存為2.00GB。圖2所示為蠟?zāi)y量模型與CAD模型的初始位置,圖3所示為經(jīng)上述方法計算得到的匹配點對。圖4和圖5分別為預(yù)配準(zhǔn)和精配準(zhǔn)的結(jié)果圖。

圖2 蠟?zāi)Ｅ錅?zhǔn)前初始位置

圖3 對應(yīng)匹配點對

圖4 預(yù)配準(zhǔn)結(jié)果

圖5 精配準(zhǔn)結(jié)果

運用本文方法對原始蠟?zāi)y量模型進行不同程度的簡化,實現(xiàn)預(yù)配準(zhǔn)和精確配準(zhǔn),并輸出配準(zhǔn)后的數(shù)據(jù),通過Geomagic Studio 7.0對配準(zhǔn)后的測量數(shù)據(jù)和CAD模型進行三維偏差計算,結(jié)果如圖6所示。表2記錄了不同簡化程度下,配準(zhǔn)計算所消耗的時間;表3記錄了不同簡化程度下對應(yīng)的三維偏差結(jié)果。

圖6 三維偏差結(jié)果

表2 不同簡化程度下配準(zhǔn)所耗時間對照表

表3 不同簡化程度時三維偏差對照表

由表2可以看出,對原始測量點云數(shù)據(jù)進行不同程度的簡化可以減少預(yù)配準(zhǔn)所消耗的時間,精確配準(zhǔn)時由于仍然是對原始測量點云數(shù)據(jù)進行計算,所以,此階段所耗時間是相同的。由表3可知,與由原始點云計算的偏差相比,簡化到不同點數(shù)時得到的三維偏差結(jié)果相差很小,以簡化到12 864個點時的結(jié)果來看,與以原始點云作為控制點集計算得出的結(jié)果相比,其配準(zhǔn)精度并沒有受到影響,而其預(yù)配準(zhǔn)階段的時間從30s降到了4s,可見速度有了很大的提高。綜合表2和表3的結(jié)果可知,對原始測量點云進行不同程度的簡化,可以提高配準(zhǔn)計算的速度,并且保證了配準(zhǔn)的精度。

6 結(jié)語

本文提出一種基于簡化模型的曲率計算配準(zhǔn)方法。實例表明,該方法既保證了配準(zhǔn)的精度又提高了配準(zhǔn)的速度。但該方法只提高了預(yù)配準(zhǔn)階段的速度,并沒有縮短精確配準(zhǔn)階段所耗的時間。

[1]Besl P J,MaKay N D.A Method of Registration of 3D Shapes[J].IEEE Transaction on Pattern Analysis and Machine Intelligence,1992,14(2):239-255.

[2]Li Y D,Gu P H.Free-form Surface Inspection Techniques State of the Artreview[J].Computer-Aided Design,2004,36(3):1395-1417.

[3]武殿梁,黃海量,丁玉成,等.基于遺傳算法和最小二乘法的曲面匹配[J].航空學(xué)報,2002,23(3):285-288.

[4]吳鋒,錢宗才,杭洽時,等.基于輪廓的力矩主軸法在醫(yī)學(xué)圖像配準(zhǔn)中的應(yīng)用[J].第四軍醫(yī)大學(xué)學(xué)報,2001,22(6):567-569.

[5]劉晶.葉片數(shù)字化檢測中的模型配準(zhǔn)技術(shù)及應(yīng)用研究[D].西安:西北工業(yè)大學(xué),2006.

[6]戴靜蘭,陳志楊,葉修梓.ICP算法在點云配準(zhǔn)中的應(yīng)用[J].中國圖像圖形學(xué)報,2007,12(3):517-521.

[7]張二虎,卞正中,張燕,等.基于ICP和 SVD的視網(wǎng)膜圖像特征點配準(zhǔn)算法[J].小型微型計算機系統(tǒng),2004,25(10):1810-1813.

[8]戴靜蘭.海量點云預(yù)處理算法研究[D].杭州:浙江大學(xué),2006.

[9]賀美芳,周來水,神會存.散亂點云數(shù)據(jù)的曲率估算及應(yīng)用[J].南京航空航天大學(xué)學(xué)報,2005,37(4):515-519.

[10]張麗艷,周儒榮,蔡煒斌,等.海量測量數(shù)據(jù)簡化技術(shù)研究[J].計算機輔助設(shè)計與圖形學(xué)學(xué)報,2001,13(11):1019-1023.

[11]Yang M,Lee E.Segmentation of Measured Datausing a Parametric Quadric Surface Approximation[J].Computer-Aided Design,1999,31(7):449-457.