基于網格特征臨界點的三維工程模型檢索算法
摘要：隨著計算機輔助設計(Computer Aided Design，CAD)技術和三維圖形硬件的不斷發展，專業化CAD軟件在工業中得到了廣泛使用。三維工程模型已成為工程分析和生產制造的基礎，是現代工程企業產品數據事實上的標準，為工程信息的構建、分析和重用提供了新的手段，大大提高了設計和制造的效率。由于產品結構越來越復雜，產品類型不斷增加，需要設計的模型越來越多，造成工程三維模型爆炸式的增長。統計顯示，在產品設計中只有20%的零部件是需要全新設計的，40%可以從現有設計中直接得到，剩下的40%可以從現有設計中修改得到;75%的新設計都需要參考已有的設計和知識口]。
現今許多企業正在建立企業內部的三維工程模型數據庫，方便了產品開發人員及時有效地獲得所需的三維模型，加快了產品開發的步伐。在客戶需求多樣化的今天，有效檢索并重用已有的三維模型及相關設計知識已成為實現產品快速研發、提高企業競爭力的重要手段。
傳統的檢索方式是將CAD模型中附帶的文件名、零部件數量或內容等信息作為關鍵詞進行檢索，這種方法相對簡單易行，但已不能滿足日益增長的檢索需求 [z]。許多學者采用基于圖(graph)的方法對模型進行檢索[3q]，并將其應用于基于實例的產品設計中。他們將零件本身的結構特征(如幾何、加工精度特征等)、工藝特征(如外圓、內孔、平面、槽等)及其相互間的關系提取出來用有向圖表示，進而通過子圖同構來檢索需要的模型。這種方法有效地利用了零件自身的信息，與領域知識關聯緊密。但前提是必須對模型進行特征識別，才能準確提取出模型的特征信息。由于不同商業CAD系統內部三維模型表示方法以及建模方式不同，阻礙了CAD系統問的產品數據交換和模型共享。目前的通用加工特征識別算法不穩定，特征識別只能針對某種CAD系統單獨進行二次開發，工作量大，且缺乏通用性和一般性。況且子圖同構算法是NP難問題，一旦零件復雜，對應的有向圖急劇膨脹，檢索效率將大大降低。
為此，本文提出一種與CAD系統無關的基于網格特征臨界點的三維工程模型檢索算法。該算法以三維模型的網格化表示作為檢索輸入，通過對網格模型的分析，找出表征網格形狀的關鍵點，即特征臨界點，以這些點為基礎計算三維模型的形狀度量，通過相似性比較，從模型數據庫中檢索出與輸入模型相似的模型。
1、原理與方法
1.1 Morse理論和網格特征臨界點
1934年，美國數學家M.Morse提出用分析方法研究空間拓撲性質，即Morse理論[5]，成為微分拓撲學的一個重要分支。空間是幾何研究的對象，而函數是分析研究的對象。因此，定義于空間上的函數與空間的形狀有著緊密的聯系。Morse理論正是研究空間上的函數與空間關系的，它特別關注的是函數的臨界點，并從臨界點的信息中獲取空間形狀的信息，即研究流形上光滑函數的臨界點與流形本身拓撲結構之間的關系。本文以Morse理論作為理論基礎，將網格模型映射為Morse臨界點的集合。
設MCR“為行維空間下的緊致流形，函數，：
M—R為作用在M上的一個光滑實值函數，在點P(越)∈M處建立局部坐標系，“=(U。，?，U。)，則：①如果函數廠在點P的梯度為零，即 3f/aU。一0，i∈[1，咒]，則P是函數，的臨界點。②如果P是函數，的臨界點，且，在P點處的Hessian矩陣日(，)(P;32，越)5(蠢女)非奇異，則P是Morse。$i界-A。
本文只考慮M是三角網格、函數廠是分段線性實值函數、且作用在網格M的頂點上的情況，三角面片內部的函數值可由頂點處的函數值插值得到。
假設對任意的邊<夕，Pz>∈M，均有廠(P，)≠廠(Pz)，可推知在任意三角面片內部，梯度都是恒定的非零實數，且I臨界點只出現在網格的頂點上[5]。
對于二維流形，Morse臨界點分為極大值點、極小值點和鞍點三類。若P為極小值點，則廠在P的/J,lI臺i域的各個方向均上升;若P為鞍點，則P 為廠在P的小臨域內上升下降轉換的過渡點;若P為極大值點，則廠在P的小臨域的各個方向均下降。圖1表示了網格頂點P的分類情況，其中圓圈表示與P相連的頂點可i，0表示f(v。)>廠(p)，e表示廠(讓)<廠(夕)。
1.2網格特征區域
Chang Ha Lee等人于2005年提出了網格特征區域(mesh saliency)的概念[6]。網格特征區域是對網格模型局部區域重要性的一種度量，即該網格模型局部區域能否體現網格形狀特點，其基本思想是函數，將模型表示為采樣點處函數值的概率分布。
計算網格頂點的平均曲率，并對其進行高斯加權平該方法易于理解，實現簡單，計算效率高，對網格退均(Gaussian-weighted average)，然后計算頂點的化的情況魯棒性好，但存在以下不足：
特征因子，以此對頂點進行濾波，找出能體現網格模(1)由于該方法中采樣點的數量是一定的，不能型外形特點的區域。根據模型的形狀與復雜程度自適應地產生采樣點，本文結合網格特征區域和Morse理論，找出網特別針對工程模型，一些很重要的局部幾何特征(如格上表征網格局部形狀的特征臨界點，并以基于這槽、凸臺)很難被充分采樣，而對于平面特征往往被些臨界點之間的距離和角度值的統計規律作為形狀過采樣，因此采樣點不能很好地體現工程模型的形描述子，對網格進行形狀比較。狀特點。
1.3基于形狀分布的圖形檢索算法
普林斯頓大學Robert Osada等人提出的基于形狀分布的方法(D2)E73是圖形檢索領域中十分重要的算法，其基本思想是：首先對模型進行采樣，生成一系列的采樣點，然后通過作用在模型上的形狀(2)此算法采用單一形狀分布函數，即任意兩采樣點間的歐氏距離，雖然計算簡單，但沒有充分利用模型的其他幾何信息，如三角面片和網格頂點處的法矢和曲率信息，所以單一的形狀分布函數很難充分反映出模型的外形特點，檢索結果差強人意。圖2所示為三個模型的比較結果，可以看出，雖然三個模型的外形差異很大，但通過此算法得到的形狀分布曲線是相似的�？梢娫摲椒ú⒉荒芎芎玫貐^分這三種模型，此算法有待改進和增強。
C.Y.IP等人根據任意兩點的連線與模型的位置關系，將距離分為In distances，Out distances和Mixed distances，改進了D2算法，并用于比較CAD模型[8]，改進方法較傳統的形狀分布算法有效，但計算過程復雜，計算量過大，實時性有待加強，且不能避免傳統形狀分布算法不足對檢索結果的影響。
2、基于特征臨界點的形狀檢索算法
為了克服傳統基于形狀分布檢索算法的不足，本文提出了基于特征臨界點的檢索算法。其基本思想是以三維工程網格模型作為輸入，通過分析網格頂點的離散法矢和平均曲率，找出表征其結構的特征臨界點，然后計算臨界點間的形狀函數值，并對其進行概率統計，給出形狀分布圖，通過計算形狀分布相對于傳統算法，該算法的優勢在于：
(1)用網格模型的特征臨界點代替傳統的采樣點進行形狀函數計算，避免了采樣不足對算法造成的影響。根據Morse理論，特征臨界點反映了空間模型的拓撲幾何信息，且其對模型中有重要工程意義的結構敏感，因此比隨機生成的采樣點更能體現模型的結構特點，更具針對性。
(2)取網格益面上兩點間近似測地線距離與兩點處法矢夾角的余弦值作為聯合形狀函數，對函數值進行概率統計，給出二維形狀函數分布圖。較之單一的形狀函數D2，近似測地距離更能表征模型表面形狀的變化，法矢夾角余弦反映了模型局部區域內形狀的分布，因此聯合形狀函數在兩個方面綜合考慮了模型的結構，更全面地反映出模型的外形特點。
(3)不計算任意兩點間的聯合形狀函數值，而是根據Morse理論，將II缶界點分為極大值點、極小值點和鞍點，針對每類臨界點計算聯合形狀函數值，模型的整體描述由三者加權平均得到。這樣提高了計算效率，使統計結果更具可比性，能夠更合理地描述模型結構。
2.1 網格微分量的計算
三角網格模型是一種離散曲面表示形式。法矢和曲率作為重要的微分幾何特性。描述了三角網格模型的局部幾何特征。由于離散形式的曲面是一種分片線性曲面，沒有連續的法矢和曲率，通常只考慮頂點處的法矢和曲率，其余各點的幾何特性可通過對頂點進行線性插值獲得。近幾年來，眾多學者提出了一系列估算離散曲面微分量的新方法。浙江大學方惠蘭等人對國際上提出的三角網格曲面上估算平均曲率的七種方法和估算高斯曲率的四種方法進行了總結和比較，指出2002年 Meyer等人提出的Voronoi方法對高斯曲率和平均曲率估算最優砷]，因此本文將該方法擴展到頂點法矢的計算，并用其計算網格頂點的離散曲率。
Voronoi方法的基本思想是把光滑曲面看作是一族網格的極限或者線性逼近，把三角網格每個頂點的微分量看作是此三角網格在此頂點一個小鄰域的平均度量。與一般算法等同看待小鄰域內所有三角形不同，該算法根據小鄰域內三角形的不同類型(銳角、直角和鈍角三角形)，分別采用不同的面積計算方法，更好地體現了面片上的微分量對頂點處微分量的影響。具體過程可參考文獻[10]。
圖4表示按照Voronoi方法計算得到的三個工程模型頂點處的離散曲率分布，顏色越深代表此處的曲率越大。2.2網格特征臨界點的計算根據 Morse理論，筆者首先取作用在網格模型上的實值函數廠為網格頂點p處的離散平均曲率H(p)與該點影響因子的代數和，則函數，的Morse臨界點體現了網格模型的空間拓撲幾何結構。
Morse臨界點的計算如下‘11。：
步驟1對每個頂點戶∈M，采用Voronoi方法計算其平均曲率H(p)。
步驟2計算頂點P對局部形狀的影響因子。
根據頂點P的鄰域對H(戶)進行加權平均，以此將頂點處的曲率與頂點的鄰域聯系起來，通過計算鄰域頂點間對局部形狀影響的差異，體現該頂點對局部形狀的重要性。采用雙向平滑算子(bilateralsmoothing operation)對H(戶)進行鄰域相關的加權平均。影響因子的計算如下：
B(H(戶)，口)一[Σ(H(z)一H(夕))J∈F‘p·2a)W。(fI z—p J)w。(I H(z)一H(夕)1)]/[Σz∈F(p·2a)Wc(1l z—P 11)Ws(1 H(z)一H(戶)})]。(1)其中，F(p，盯)為距離頂點戶在盯范圍內的頂點的集合;Wc(z)=exp[-一z2/(2口冬)]為標準高斯濾波函數，W。(z)=exp[一，/(2盯§)]為特征保持加權函數。
步驟3根據B(H(p)，盯)更新函數廠的值。
步驟4控制計算的迭代次數，控制最終得到的特征臨界點的數量。迭代次數越多，生成的特征臨界點的數量越少。
步驟5根據Morse臨界點的分類，將得到的臨界點分為極大值點、極小值點和鞍點。
上述過程的偽代碼如下：Procedure SalientCriticalPoints(M，d，iteration)輸入：M為三角網格模型，d為頂點鄰域半徑，iteration為迭代次數。
輸出：Sc為特征臨界點集合。
局部變量：N為M中頂點的數量，Pl為M中的第i個頂點，f(pi)為頂點pl處的實值函數值，a，為頂點Pl的影響因子。
Begin集合{P1)=M的所有頂點;N一集合{Pl}的勢;for(j=1 to iteration)計算頂點P，處的平均曲率H(pO，令f(pi)=H(pi)lfor(i一1 to N)計算頂點pl的影響因子ai(M，f(pi)，d);更新f(Pi)+一8i;endend將得到的實值函數值{f(P1)}，利用Morse臨界點的定義判斷Pi的類型，return網格特征臨界點集合Sc，end圖5是對3個零件網格模型進行臨界點計算的結果，其中第一列是模型所有的臨界點，第二列是極大值點，第三列是極小值點，第四列是鞍點，圖片下方的數字表示此類臨界點的數量。
2.3形狀函數
形狀函數定義了能夠表征模型形狀的幾何特征，通過計算臨界點處的函數值來反映模型的結構特點，因此形狀函數的選擇對提高算法的性能至關重要。傳統的 D2距離雖然計算簡單，但并不足以反映模型的外形特點，因此檢索結果不理想。本文采用兩點間近似測地距離和兩點處法矢夾角余弦作為聯合形狀函數，近似測地距離反映了整個模型表面形狀的變化，法矢夾角余弦體現了模型局部區域形狀的分布。
(1)近似測地距離函數測地線在微分幾何中有著嚴格的定義，即曲面上的一條曲線，如果它的每一點處的測地曲率為零，則稱為測地線，在工程中可將其直觀理解為曲面上連接兩點的最短路徑。由于精確計算測地線算法復雜、時間復雜度較高，本文采用近似測地線，網格上組成近似測地線的所有小線段的歐式距離和即為近似測地距離。
本文采用Takashi Kanai等人提出的算法，計算網格上兩頂點間的近似測地路徑，其思想是首先將網格的頂點看作圖的節點，網格頂點間的連接看作圖的邊，采用計算圖最短路徑的Dijkstra算法獲得初始最短路徑，然后對此路徑所在的局部鄰域進行網格細分，再對細分后的區域采用Dijkstra算法計算最短路徑，以此迭代，直至近似最短路徑的長度變化小于預先設定的閾值。具體實現細節可參考文獻[12]。為了更快速地計算，本文采用最短路徑快速算法(Shortest Path Faster Algorithm，SPFA)代替Dijkstra算法[131。
由近似測地距離的定義可知，此形狀函數是平移、旋轉不變的，但因距離是對模型大小的一種度量，故是縮放可變的。圖6表示了從模型中某一臨界點到其余同類臨界點間的近似測地路徑。
極大值點極小值點鞍點圖6網格模型l司類臨界點同的近似測地路徑(2)頂點法矢夾角余弦定義NA為法矢夾角余弦，設臨界點p;，p，處的單位法矢分別為露，。和n，。，則由向量內積可得n以·n聲，==I刀以I×I n戶;I×cosL(力以，n聲，)=》NA=cosL(肛p，np)一np·np。(2)因為在區間[o，丌]內余弦值是單調下降的，所以可以用來唯一度量兩矢量的夾角。由定義可知此形狀函數是平移、旋轉和縮放不變的。
2.4形狀分布矩陣
根據特征臨界點的類型，分別計算每類臨界點(即極大值點、極小值點和鞍點)中兩兩頂點之間的聯合分布函數。因為特征點的類型不同，其所表征的網格上局部區域的凸凹性不同，所以分類型計算聯合分布函數使得到的形狀分布矩陣更具可比性，并且把計算聯合分布函數局限到每一類中，比傳統的計算任意兩點間的形狀函數的計算量大大減少，效率更高。這樣做雖然忽略了兩組l臨界點間的拓撲關系，但因為本檢索算法本質上是基于概率統計的，并且每類臨界點的分布并沒有局限到某一區域，而是在模型整體上分布的，所以只要有足夠的臨界點，足夠的臨界點間形狀函數值，該模型的形狀特點就能通過概率統計反映出來。絕大多數情況下不同組內的函數值加權得到的統計規律已能很好地表現出模型的形狀分布規律。
計算每類臨界點兩兩之間的測地線距離和法矢夾角余弦值。由于近似測地線距離函數是縮放可變的，要對其進行歸一化處理，常用的方法有最大值歸一和平均值歸一�？紤]到最大值歸一易受噪聲數據影響，故采用平均值歸一化：設每一類中計算得到近似測地距離的最大值為D……，平均值為D……，需要將距離 [o，D……]平均劃分為L個單元，則將[o，n。]和[D……，D……]分別等分為L。/2個單元，即得到平均值歸一化的分布。NA值域為[一 1，1]，設需要將其平分為L，個單元。對每對臨界點的聯合形狀函數值進行統計，可生成形狀分布矩陣MD=(z州)，×L√z。為近似測地距離￡屬于單元 i且NA屬于單元�！旱呐R界點對總數占所有臨界點對總數的比例。將P，NA，l：值分別對應z，Y，z軸，可得到聯合形狀函數的二維概率分布圖。
分別對三類臨界點進行上述計算得到三個矩陣，則一個網格模型即可抽象為三個聯合形狀分布矩陣的線性組合，即M=謝DIn。;+凸 MDmi。+7，MDs“。(3)其中：M為網格模型的形狀分布矩陣;MD加剃MDm"M腑a分別表示與極大值點、極小值點和鞍點對應的形狀分布矩陣;a，J9，y分別表示三類矩陣的系數，且滿足歸一化條件口+口+y=1。設N。，N。t。，N州分別表示模型中極大值點、極小值點和鞍點的數量，本文取口5瓦i可吒而N’?盧一—Nma --1-N—minN-。}I_n-N,.a’
(4)圖7表示按照本算法得到的歸一化后的模型形狀分布。其中三維柱形圖分別表示了一個模型總體在極大值點集、極小值點集和鞍點集層次上的形狀分布。從圖中可以看出，分屬不同類別的兩個模型，其形狀分布具有較大差異，因此本算法具有較好的區分能力。
圖7也比較了統一計算任意兩臨界點的形狀函數值和按組計算形狀函數值得到的模型形狀分布的差異。前者因用于統計的函數值大大增加，統計結果較分組計算的結果分布稍均勻一些，但最后得出的統計規律與分組計算形狀函數值并無本質區別。平面著色圖反映了在一個柱形圖中單元數值的分布情況，單元的顏色越亮，數值越大，可見同一模型的三個柱形圖內部數值分布規律并不相同，說明了將特征臨界點分為三類分別進行計算的合理性。在工程模型中存在大量的相互之間成直角的面，形狀分布矩陣中的非零值會集中在表示180。，90。，o。的三行中(Y軸度量的是夾角的余弦值，分別對應Y軸的0，10，20)，針對該情況筆者進行了專門的實驗。圖8說明了分布在上述三行中的存在大量互成直角面的模型的形狀分布特點，它與存在其他分布特點的模型并無本質不同。實驗證明，此類模型的形狀函數值雖然主要分布在這三行中，但每一行的具體分布規律因模型不同而有所不同，因此本算法也能較好地區分其間的差別。
圖8大量存在互成直角面的模型及其形狀分布2.5形狀比較經上述處理，模型的比較轉化為形狀分布矩陣的比較。本文仍然采用L。范數，即歐式距離計算兩矩陣間的距離。令MA=(m。)I，xL，B=(m6)f。×L。v M分別為模型A，B對應的形狀分布矩陣，則模型A，B間的距離為D(A，B)=D(MA，MB)一
3、算法復雜度分析
由算法過程可知，提取網格特征臨界點和計算近似測地距離是時間復雜度較高的兩個步驟。
3.1提取網格特征臨界點的時間復雜度
(1)網格微分量計算需計算每個三角面片的面積，最壞情況下，即每個三角面片為銳角三角形時，按照Voroini方法，一個三角面片的面積需計算三次(三次面積各不相同)。由于面積計算僅是單獨的乘法，與其他步驟相比，可認為此步驟的時間復雜度為常數。
(2)頂點影響因子計算最耗時部分為搜索任意頂點在某鄰域內的所有頂點。本文將頂點數據存入k-d樹來加速此過程。設網格頂點數為N，計算臨界點時的循環次數為艴，文獻[-14]中已證明，在一棵肛d樹中插入和搜索節點的平均時間復雜度為0(109：N)，故此步驟最終的時間復雜度為 o(nN·l092N)。
3.2計算近似測地距離的時間復雜度
本文采用SPFA算法計算從任意頂點到其余頂點的近似最短路徑。設網格模型的邊數為E，求得的所有臨界點數量為Nc，文獻E13]已證明，計算從任意點出發到所有頂點的近似最短距離的時間復雜度為0(E)，故此步驟最終的時間復雜度為o(Nc·E)。
4、算法驗證與討論
以Visual C++6.0為集成開發環境，結合Matlab 6.5實現了本文提出的算法，并在普渡大學建立的工程標準模型庫ESB(engineering shapebenchmark)‘153上進行了驗證。實驗采用PC機，AMD 2500+CPU，512 M內存。
ESB中包含866個STL格式的工程模型，筆者以其中任意一個模型作為輸入，欲檢索出模型庫中相似的模型，并與傳統圖形分布算法(D2 shape distribution)和增強圖形分布算法[16](D-IA shape distributionwith simplification)進行了對比。本文描述的算法以模型的STL格式作為檢索輸入，實驗中的參數為L，=64，L，=20;根據實驗結果，取盯=0.24%￡，￡為網格模型包圍盒對角線的長度，并取UC一如一口，iteration一10。圖9列出了前五個檢索結果，其中的數字表示檢索到的模型與輸入模型的距離。從圖中可以看出，基于特征臨界點的算法能夠比另外兩種方法檢索出更多有效的模型，這是因為特征臨界點本身就是模型的頂點，由 Morse理論可知，臨界點表征了模型的空間拓撲結構和幾何結構，因此比基于采樣點的方法更能反映模型自身的特點;且本算法采用的近似測地距離代替傳統的 D2距離，測地線本身就是對模型表面形狀的反映，其長度亦體現了模型局部形狀的變化，用于做形狀函數效果較好，但其計算比D2距離要復雜。因此，鑒于對時間復雜度的考慮，近似測地距離是一種較好的折衷。
現對不同檢索輸入的檢索時間和有效性進行統計，可得到算法的平均性能。針對ESB庫中的三大類模型，即薄壁件、棱柱類零件和旋轉件分別統計，將其綜合得到算法在整體ESB庫上的平均性能，做出查準率一查全率(Precision—Recall，PR)曲線，并與其他形狀檢索算法(表面積與體積描述 [151、角度分布直方圖[17‘、3D球諧描述子[18]、2D形狀直方圖[1?、光場描述子[20])進行了比較，如圖10所示。
從PR曲線可以看出，本方法雖然目前還達不到光場描述子的效果，但它有效提高了傳統D2形狀分布算法的檢索性能，在計算復雜度和檢索效果之間找到了一種平衡。
對圖10的實驗數據進行綜合，表1定量地比較了D2形狀分布、帶網格簡化的D-IA形狀分布和本文算法的檢索精度和效率，FT，ST，NN的含義見參考文獻E16]。從中可以看出，基于特征臨界點的方法雖然總的檢索時間有所增加(主要是計算近似測地距離需時較多)，但仍在可接受的范圍內，且其檢索結果要明顯優于另外兩種方法。
5、結束語
針對傳統基于形狀分布檢索算法的不足，提出了基于特征臨界點的檢索算法。該算法以三維工程網格模型作為輸入，以Morse臨界點理論為依據，用網格模型的特征臨界點代替傳統的采樣點進行形狀函數計算，避免了采樣不足對算法造成的影響;取網格曲面上兩臨界點問近似測地線距離與兩點處法矢夾角的余弦值作為聯合形狀函數，更好地表征了模型整體表面和局部區域的變化;分別針對極大值點、極小值點和鞍點計算聯合形狀函數值，模型的整體描述由三者加權平均得到。實驗結果表明，本算法客觀反映了工程模型的相似程度，明顯提高了傳統圖形分布方法檢索的有效性。
以下幾項工作在未來研究中應該著重加強：①本文算法的時間效率雖然可以接受，但還有很大的提升余地來找到效率更高的近似測地距離算法。②采用歐拉距離計算兩矩陣間的距離存在缺點。形狀矩陣中的每一列(行)代表一形狀度量向量，不同度量向量對于區分形狀有著不同的重要性，而歐拉距離將矩陣不同列(行)之問的差別等同看待，因此最終的距離只反映了平均綜合的結果，不能表示出某一形狀度最向量對總體差異的影響。這樣導致了兩矩陣即使距離較近，但其代表的圖形可能并不相似。
因此可以嘗試其他的度量方法，從而不僅反映矩陣間的距離，而且能體現矩陣自身的分布。③由于特征臨界點體現了模型局部的形狀信息，可以考慮將其用于局部形狀檢索。
【基于網格特征臨界點的三維工程模型檢索算法】相關文章：
基于QBASIC環境下的數學算法教學11-14
基于滾動計劃的動態企業資源優化模型03-29
探析基于VaR模型的證券投資組合風險12-05
基于勝任力模型的組織生涯管理探析12-12
《基于導納的圖像加密算法的研究》開題報告12-03
基于教學知識點的模型框架與結構分析12-05
淺談基于知識的網格技術應用研究11-20
- 相關推薦