黃勝利 卜 昆 程云勇 周麗敏
西北工業(yè)大學現(xiàn)代設計與集成制造技術教育部重點實驗室,西安,710072
在渦輪葉片模具型腔反變形優(yōu)化設計中,計算分析葉片測量數(shù)據(jù)相對于CAD模型的三維偏差至關重要。葉片測量數(shù)據(jù)的主要獲取方法為三坐標測量和光學掃描測量。近年來,隨著激光技術的發(fā)展,光學掃描儀獲取數(shù)據(jù)的速度和精度都有顯著的提高,應用越來越廣泛,但其獲得的海量數(shù)據(jù)給后續(xù)的數(shù)據(jù)處理、分析、計算帶來很大的不便。由于測量數(shù)據(jù)是在測量設備坐標系下獲得的,故在對測量數(shù)據(jù)進行偏差分析前,需把測量數(shù)據(jù)與CAD模型的坐標系統(tǒng)一起來,也就是實現(xiàn)測量數(shù)據(jù)與CAD模型的配準定位。
目前應用比較廣泛的配準算法是Besl等[1]提出的迭代最近點(iterative closest point,ICP)算法。該方法通過最小二乘法構造配準的目標函數(shù),具有收斂速度快、配準精度高的優(yōu)點?;贗CP的改進算法[2]在含有自由曲面的工件配準中得到了廣泛的應用。ICP算法要求待配準兩模型的初始位置要足夠接近,否則不能得到好的收斂結果,因此,使用ICP算法進行精確配準前,需要對兩待配準模型進行預配準。武殿梁等[3]采用遺傳算法實現(xiàn)了曲面之間的預配準,此方法比較耗時。吳鋒等[4]基于輪廓的力矩主軸實現(xiàn)了醫(yī)學圖像的預配準,該方法的缺點是對數(shù)據(jù)的缺失比較敏感,需要整個物體的全部信息。劉晶[5]通過建立模型的包圍盒實現(xiàn)三維測量模型與CAD模型之間的預配準,該方法對測量模型的完整性要求較高,且包圍盒的建立有一定的難度。
光學測量的密集點云數(shù)據(jù)的配準計算非常耗時,為了縮短配準計算的時間,本文提出一種基于簡化模型曲率計算的配準方法,首先對密集點云模型進行簡化,估算簡化模型的曲率,根據(jù)曲率提取對應特征面;然后尋找匹配點對從而計算出變換矩陣,把求出的變換矩陣應用于簡化前的原始點云模型以實現(xiàn)模型的預配準;再通過奇異值分解和最近點迭代相結合(singular value decompo-sition and iterative closest point,SVD-ICP)算法[6-7]進行精確配準;最后以某型渦輪葉片的蠟模為例,驗證方法的有效性和實用性。
由于光學測量獲取的數(shù)據(jù)量巨大,為了提高計算的速度需對點云模型進行適當?shù)暮喕D1所示分別為某型葉片及其蠟模的光學測量模型。目前,常用的點云簡化方法有包圍盒法、均勻網(wǎng)格法、聚類法、隨機采樣法、曲率采樣法、均勻采樣法等[8],它們各有優(yōu)缺點和不同的適用范圍。本文對點云模型進行簡化的主要目的是提高計算的速度,因此采用均勻采樣法。
圖1 光學測量密集點云模型
曲率是曲面的重要幾何特征,具有平移、旋轉和縮放不變性,可以作為曲面特征識別的重要依據(jù)。本文采用鄰域內(nèi)擬合二次參數(shù)曲面的方法估算散亂點云的曲率[9]。
對原始點云進行空間劃分,計算點云中每個點pi的K 鄰域,K鄰域記為Knbhd(pi)。每個Knbhd(pi)可以擬合成一曲面且此曲面可以被一個局部基面逼近,該局部基面記為LBS(Knbhd(pi)),本文選擇點pi處的切平面為局部基面。切平面可通過各點的鄰近關系和最小二乘原理[10]求得。對于不同的模型,合理選取K的值(K一般取10~30),以保證鄰域內(nèi)凹凸一致,這樣得到的切平面才能更好地逼近原始點集,使投影點在局部基面的參數(shù)化更好地反映點云的參數(shù)化。
把Knbhd(pi)中的各點在相應局部基面上的投影點記為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)中對應點的參數(shù),就實現(xiàn)了局部散亂數(shù)據(jù)的參數(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的下標為式(2)中A的行號和列號。采用最小二乘法使Knbhd(pi)上各點到二次參數(shù)曲面的歐幾里德距離之和最小,可得
式中,X為鄰域Knbhd(pi)的矩陣。
由于計算得到的曲面法矢方向一般是不協(xié)調(diào)的,即并不是所有的法矢都指向曲面的同一側,因此有必要對求得的法矢進行歸一化調(diào)整,具體方法參見文獻[9]。
根據(jù)曲面的曲率性質(zhì)可以得出高斯曲率:
平均曲率:
上述物理量可由曲面的第一、第二基本形式推出。
渦輪葉片主要由自由曲面和一些平面組成,而自由曲面的曲率特征沒有一定的規(guī)律性,因此,曲面的提取有一定的難度。本文中,將測量模型和CAD模型上某一組對應平面作為尋找匹配點對的依據(jù)。高斯曲率I和平均曲率H均分別反映了曲面的凸凹形狀,I和H的不同組合表現(xiàn)的是不同的幾何形狀,表1給出了幾種常見的曲面形狀與I、H組合的關系。
表1 點的局部曲面類型
由表1可以看出,當I、H 都為0時,曲面就是一平面,根據(jù)此條件可提取模型的某一平面為特征面,為尋找對應匹配點作好準備。
為了分析鑄件相對于其CAD設計模型的彎扭變形,首先需要實現(xiàn)測量模型與CAD模型的配準。運用SVD-ICP算法進行精確配準,要想獲得較好的效果,首先需要實現(xiàn)測量模型與CAD模型的預配準。本文通過曲率提取特征面的方法實現(xiàn)預配準,具體方法如下:
(1)運用均勻采樣法對原始測量點云{P}進行簡化得到{P′}。
(2)計算簡化后點云{P′}和CAD模型的重心(分別記為O1、O2),平移兩模型,使兩模型的重心與坐標原點O重合。
(3)運用第2節(jié)中所述方法估算簡化后測量模型和CAD模型的高斯曲率K與平均曲率H。
(4)根據(jù)所求得的K和H提取兩模型中相對應的特征面,分別記為L1、L2,并求取兩特征面的
形心 P1、Q1。
(5)把P1、Q1沿所在特征面的法線方向移動一定距離 d(d值可適當選取,本文實驗中 d=10mm),得點P2、Q2的計算公式:
式中,P1、P2、Q1、Q2為點 P1、P2、Q1、Q2的位置矢量;N1、N2分別為特征面的單位法矢。
(6)由三組對應點{O,P1,P2}和{O,Q1,Q2}求取旋轉矩陣R和平移矩陣T,方法如下:以O點為坐標系原點,由單位矢量組 e1 、e2 、e3 與e′1 、e′2 、e′3按右手法則構建兩個局部坐標系;設存在旋轉變換矩陣R與平移變換矩陣T使兩個局部坐標系重合,于是有
由式(13)可求得旋轉變換矩陣:
待配準模型中的任一點Pi變換到對應點Qi的關系式為
把P1、Q1代入式(15)可得到平移變換公式:
(7)把求得的旋轉矩陣和平移矩陣應用到簡化前的原始測量點云模型上,得到變換后的測量模型,再按向量OO2進行整體平移就實現(xiàn)了測量模型與CAD模型的預配準。
實現(xiàn)了測量模型與CAD模型的預配準后,就可以運用SVD-ICP算法進行精確配準。記測量點云數(shù)據(jù)為{pi},與其對應的CAD模型上的對應點集為{p′i}。配準目標函數(shù)表示為
對配準點集{pi}和{p′i}運用文獻[7]中奇異值分解法求取變換矩陣R和T,并將變換矩陣作用到測量點云數(shù)據(jù)上,迭代地進行這個操作,直到使目標函數(shù)值滿足設定的閾值為止。這樣就實現(xiàn)了測量點云模型與CAD模型坐標系的統(tǒng)一,即實現(xiàn)了測量模型到CAD模型的精確配準。
為了驗證本文方法的有效性,以VC++6.0為平臺,以VTK可視化類庫為輔助工具編程實現(xiàn)了上述方法,并以某型渦輪葉片的蠟模為例進行了驗證。試驗中所用的測量點云模型是由ATOS流動式光學掃描儀獲取的,其數(shù)據(jù)點數(shù)為102 908。試驗所用的電腦配置為Intel(R)Core(TM)2 Duo CPU E6550@2.33GHz;內(nèi)存為2.00GB。圖2所示為蠟模測量模型與CAD模型的初始位置,圖3所示為經(jīng)上述方法計算得到的匹配點對。圖4和圖5分別為預配準和精配準的結果圖。
圖2 蠟模配準前初始位置
圖3 對應匹配點對
圖4 預配準結果
圖5 精配準結果
運用本文方法對原始蠟模測量模型進行不同程度的簡化,實現(xiàn)預配準和精確配準,并輸出配準后的數(shù)據(jù),通過Geomagic Studio 7.0對配準后的測量數(shù)據(jù)和CAD模型進行三維偏差計算,結果如圖6所示。表2記錄了不同簡化程度下,配準計算所消耗的時間;表3記錄了不同簡化程度下對應的三維偏差結果。
圖6 三維偏差結果
表2 不同簡化程度下配準所耗時間對照表
表3 不同簡化程度時三維偏差對照表
由表2可以看出,對原始測量點云數(shù)據(jù)進行不同程度的簡化可以減少預配準所消耗的時間,精確配準時由于仍然是對原始測量點云數(shù)據(jù)進行計算,所以,此階段所耗時間是相同的。由表3可知,與由原始點云計算的偏差相比,簡化到不同點數(shù)時得到的三維偏差結果相差很小,以簡化到12 864個點時的結果來看,與以原始點云作為控制點集計算得出的結果相比,其配準精度并沒有受到影響,而其預配準階段的時間從30s降到了4s,可見速度有了很大的提高。綜合表2和表3的結果可知,對原始測量點云進行不同程度的簡化,可以提高配準計算的速度,并且保證了配準的精度。
本文提出一種基于簡化模型的曲率計算配準方法。實例表明,該方法既保證了配準的精度又提高了配準的速度。但該方法只提高了預配準階段的速度,并沒有縮短精確配準階段所耗的時間。
[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].航空學報,2002,23(3):285-288.
[4]吳鋒,錢宗才,杭洽時,等.基于輪廓的力矩主軸法在醫(yī)學圖像配準中的應用[J].第四軍醫(yī)大學學報,2001,22(6):567-569.
[5]劉晶.葉片數(shù)字化檢測中的模型配準技術及應用研究[D].西安:西北工業(yè)大學,2006.
[6]戴靜蘭,陳志楊,葉修梓.ICP算法在點云配準中的應用[J].中國圖像圖形學報,2007,12(3):517-521.
[7]張二虎,卞正中,張燕,等.基于ICP和 SVD的視網(wǎng)膜圖像特征點配準算法[J].小型微型計算機系統(tǒng),2004,25(10):1810-1813.
[8]戴靜蘭.海量點云預處理算法研究[D].杭州:浙江大學,2006.
[9]賀美芳,周來水,神會存.散亂點云數(shù)據(jù)的曲率估算及應用[J].南京航空航天大學學報,2005,37(4):515-519.
[10]張麗艷,周儒榮,蔡煒斌,等.海量測量數(shù)據(jù)簡化技術研究[J].計算機輔助設計與圖形學學報,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.