三維點云法向量的模糊估值算法*

劉艷菊 張永德 楊波

(1.哈爾濱理工大學 智能機械研究所，黑龍江 哈爾濱 150080;2.齊齊哈爾大學 計算中心，黑龍江 齊齊哈爾 161006)

應用三維點云的曲面重建技術有很多，但大多數技術在三維重建過程中都要求點云數據含有法向量.法向量是點云數據的一個基本屬性，用于標識曲面的內/外側.雖然法向量可以通過三維激光掃描儀獲得的深度圖像而得到，但掃描過程中不可避免的噪聲、振動等因素使得到的法向量不準確甚至有缺失，很難滿足曲面重建的要求，有必要對法向量重新進行估值.

目前法向量估值算法總體分為兩類:局部估值和全局估值.在局部估值算法中，文獻［1］的原則組元分析(PCA)法通過計算兩點間歐式距離得到法向量所在直線，但該方法只適用于規則采樣表面，要求點云必須均勻，實際上很難達到這樣的要求.文獻［2］在此方法基礎上提出局部權值優化投影因子(WLOP)算法，從原始數據中篩選提取能夠滿足PCA 要求的點云數據，然后進行法向量估值，雖然效果較好但程序復雜度增加一倍.文獻［3］用Mahlanobis 距離替代歐式距離的算法，可以同時求出法向量所在直線及其方向，但對點云數據要求也非常高，同樣難以實現.全局估值法向量的算法主要是由輸入點云數據構建Voronoi 圖［4］，將Voronoi 圖中的每個極點看做是數據點的法向量;這類算法可以簡單高效地估算出法向量，但如果點云數據中含有噪聲點估算出來的法向量就會有較大的偏差.文獻［5］將局部算法和全局算法相結合，具有二者的優勢但效率較低.

近期又有一些研究新思路:文獻［6］利用魯棒統計算法進行法向量估值，該算法適用于含有尖銳特征的模型，而且需要人工調整參數.文獻［7］利用圖的優化理論對法向量的方向進行推導，但此算法不適用于薄片類模型.文獻［8］采用變分算法既可以估算法向量所在直線又可以確定其方向，該算法對較薄模型推導效果較好，但對于含有尖銳特征的曲面容易產生錯誤結果.總體上，每種算法都在某些模型上能夠實現較好的法向量估值而忽略其他情況，不能適用于任意形狀模型的法向量估值.

由于不同模型差別很大，很難用一個統一的數學模型表達出來.筆者在文獻［9］的基礎上提出基于模糊推理的法向量估值算法，對獲取的點云數據進行單元劃分，并且求出每個單元中數據點的k 鄰近距離和曲率，將這兩項數據輸入到模糊系統中，根據人類經驗給出模糊規則，通過模糊推理確定該部分數據點適合的估值算法，然后按照相應的法向量估值算法進行計算.

1 法向量模糊估值系統及相應算法

1.1 模糊推理系統

模糊推理系統［9］是由IF-THEN 規則構建而成的，可以通過系統的方法將專家的語言信息結合到系統中.該系統包括模糊化、模糊規則集合、模糊推理機制和反模糊化，如圖1 所示.模糊化是將非模糊的輸入空間映射到模糊集合;反模糊化是執行相反映射從模糊集合映射到非模糊的輸出空間;模糊集合的特征是由隸屬函數μ(·)表示;模糊推理機制是執行基于模糊規則的一個從輸入的模糊集合空間U∈Rn向輸出模糊集合空間V∈Rm的映射.設U =U1U2…Un，Ui∈R，i =1，2，…，n;V =V1V2…Vm，Vk∈R，k=1，2，…，m.下面詳細介紹模糊推理系統在法向量估值上的具體應用.

圖1 模糊推理系統組成Fig.1 Composition of fuzzy inference system

對于點云集合中某個數據點pk，其k 個相鄰的數據點集合為數據點的k 鄰近距離和曲率的組合能表示出模型所在區域的特征，所以將這兩個變量作為模糊規則的條件給出一組IF_THEN規則:

其中:x 是模糊推理系統的輸入向量，x =［xknn，xcur］T，xknn是點pk的k 個鄰近點的歐式距離之和，，距離之和大說明點云分布比較松散，距離之和小說明點云分布密集;xcur是點pk處的曲率;y 是模糊推理系統的輸出向量，(i =1，2)，Bk(k =1，2，…，m)分別是模糊集合Ui和Vk的語言變量，由隸屬函數(xi)和(yk)(j=1，2，…，M)表示.其中:xi是xknn或xcur中的任意一個，根據i 的取值來確定;yi是y 的某個取值;M 是模糊規則的總數量.使用推理機制、單值模糊化和中心平均反模糊化策略的模糊推理系統的輸出表示為

ξ(x)是模糊基函數向量，ξ(x)=［ξ1(x)，ξ2(x)，…，ξM(x)］T，其中ξj(x)可以被定義為

為了系統的計算效率，xknn和xcur的隸屬函數采用三角型隸屬函數.

1.2 法向量的模糊估計算法

經過對物體形狀的分析，總體上將物體形狀分成3類:局部平滑部分、薄片突出部分和不規則尖銳特征部分.一個物體可能含有其中1 種情況，也可能具有上述3 種情況.只需要將物體的點云數據中對應的k 鄰近距離及其曲率輸入到模糊系統中，按照模糊規則就可以將此物體每一部分屬于哪種類型區分出來，根據不同類型給出不同的法向量估值算法.

對于平滑部分的法向量估值可以直接應用PCA法計算法向量;對于薄片部分文中增加了一個檢測器，專門對此部分的法向量進行特定估值;不規則尖銳部分采用增加附加點的方式進行估值，具體方式為

其中，［1 α γ］是由模糊推理系統的輸出y 確定的.公式(4)中第1 項NPCA是每個數據點都要進行的估算，所以其系數是1，對于平滑曲面只采用此估算就可以了.NPCA采用常用的原則組件分析方法，參見文獻［10］.Nchk和Natt是專門分別處理薄片部分和尖銳特征部分的算法.

1.2.1 檢測器估值算法

在曲率變化很大，而k 鄰近距離也很大時，說明該區域是屬于薄片特征的部分，可以增加檢測器［11］對法向量進一步估值.檢測器沿兩個相鄰數據點的局部切線方向，按照k 鄰近搜索方式簡單地逼近薄片關鍵點.算法步驟如下:

(1)通過數據點pi的k 鄰近距離投影到其切平面上，來確定當前pi的法向量的方向.

(2)如果pi的k 鄰近距離投影是在凸出部分的外側，那么pi被認為是在薄片位置的數據點.

(3)取pi的k 鄰近點pj，(ni，nj)分別是這兩點的法向量所在直線，按照公式(5)推導:

式中:Dij是pi、pj兩點的投影距離;fbs是線段的中點是沿著各自法向量所在直線距離pi和pj的單位長度;obs是fbs投影到線段的點.

(5)在薄片位置的數據點pj被重新定向，即ninj=-1，但它的方向不允許傳播給其他點.

(6)按照這種方式，在法向量方向迭代檢測過程中，執行檢測器可以將相鄰的曲面方向翻轉.

1.2.2 附加點估值算法

在曲率變化大且不規則、k 鄰近距離之和小的區域，點云數據分布密集，且含有尖銳特征.文中采用在曲面的內側添加附加點［12］的方式幫助識別曲面的方向.附加點是從給定的點P 沿著法向量N 的方向投影產生的點，附加點與原始點距離是β，附加點集合={pn+1，pn+2，…，p2n}，按如下公式產生:

其中，N 是由PCA 方法計算每個數據點得到的無方向的法向量直線，N = {n1，n2，…，nn};通過設置β(0 ＜β≤0.001)的取值，附加點被定位在期望曲面的內側，與原始數據非常接近;公式(6)中設pi是原始點云.將所有的數據點連接生成四面體模型，如圖2 所示.

圖2 加入附加點后的四面體Fig.2 Tetrahedral with attach points

在圖2 中，四面體可以分為3類，分別是:①四面體的各頂點只含有原始點云集合P 中的點;②四面體的各頂點同時含有集合P 和集合的點;③四面體的各頂點只含有集合中的點.在這3類四面體中，①和③類在法向量估計中是冗余的，可以刪掉.刪除①和③類四面體后，就能夠通過含有P 和的三角網格推導出恰當的方向了.具體算法步驟如下:

(1)按照式(6)，求出附加點pn+i，構建附加點集合.

(2)應用合并集合中的數據點構建四面體模型，具體方法見參考文獻［12-13］.

(3)在點云集合中，設一個數據點pk為起始點，構建pk周圍的所有三角片集合{}，g 是環繞在點pk周圍的三角片的數量;構建以pk為頂點的各角集合{}，如果角度之和為360°，pk點的法向量按公式(7)估算:

(4)構建pk周圍的所有三角片邊長集合}，其中是第i 個三角片的3 條邊中的最長邊.因為冗余三角片與其他三角片相比總是有較長的邊，所以對三角片集合中各最長邊集合按升序排序.排序后求的和，j 是使等式成立的值.同時認為是冗余三角片，將其刪除.接下來按式(7)求pk點的法向量值，然后執行步驟(6).如果不存在j，按步驟(5)處理.

(5)采用加權求平均的方式計算法向量:

式中:hik為兩個數據點之間的距離，w(h)為權值函數、w(h)=1-6J2+8J3-4J4，J=h/R，R 為點pk的支撐半徑;pk的k 鄰近點相應的法向量集合是已經按照上面步驟(3)、(4)估算出的法向量的值;是所包含的鄰近點數量.

(6)循環檢測此區域的數據點，即可估算出相應的法向量.

2 算法驗證

牙齒模型是一個典型的不規則模型(如圖3 所示)，牙齒的咬合面部分含有不規則的尖銳特征，切齒、犬齒含有薄片特征，而牙冠和牙齦部分的曲面較平滑.

圖3 不規則模型—牙頜Fig.3 Irregular model—denture

2.1 牙齒模型的點云數據分析

文中得到牙頜模型的點云數據后，計算每個數據的k 鄰近距離之和(i =1，2，…，n)，按升序排序并添加到優先級隊列，從而確定點云的疏密程度.計算數據點的曲率xcur，在點云集中的區域沿著曲率變化最小的方向調整法向量的方向可以避免最大的可能誤差以及在具體實施過程中的“死鎖”問題，例如，含有尖銳特性的牙尖或牙窩.從xknn和xcur的取值可以知道點云數據密度和尖銳特征的分布情況，如圖4 所示.

圖4 點云的分布和曲率Fig.4 Distribution and curvature changes of point clouds

2.2 牙齒模型對應的模糊規則

經過對牙頜模型和實體牙齒的點云數據測量分析，給出相應的模糊規則:

模糊規則中xknn和xcur的隸屬函數取值如圖5所示.圖5(a)中，“大”、“小”對應的取值范圍分別是(0.05，+∞)、(-∞，0.05］;圖5(b)中“負無窮大”、“正無窮大”對應的取值范圍分別是(-∞，-25］，［25，+∞)，“小”的范圍是(-25，25).

圖5 三角型隸屬函數Fig.5 Triangle membership functions

2.3 牙齒模型的法向量估值

文中采用21 號標準牙齒和牙頜模型的數據，按照前面給出的模糊推理系統，分別由切齒、犬齒、臼齒和牙頜模型的點云數據計算出相應的k 鄰近距離和曲率，并輸入到模糊系統中，按照模糊推理規則將點云數據進行分類估值法向量，其中計算含有尖銳特征時增加附加點使用的參數β =0.000 75.同時采用文獻［2］的WLOP 算法對上述數據進行法向量估算.表1 給出文中算法與WLOP 算法執行時間的比較結果.可以看出，數據量越大估值時間越長，當模型中出現多種特征時，文中算法明顯優于WLOP 算法的法向量估值.

表1 文中算法與WLOP 算法的執行時間比較Table 1 Running time comparison between WLOP algorithm and the proposed one

兩種方法估算的法向量在進行隱式曲面擬合［14-15］時，設擬合精度δ =5 ×10-4，圖6 分別給出上述4 種模型利用兩種法向量估值進行曲面擬合的時間比較.從圖中可以看出，由于WLOP 算法法向量估算精度低，所以曲面擬合時間相對較長，尤其是對于數據量大的牙頜模型.曲面擬合的結果如圖7所示.

圖6 兩種算法的曲面擬合時間比較Fig.6 Comparison of surface fitting time between two algorithms

圖7 切牙、犬齒、臼齒和牙頜模型擬合結果Fig.7 Models of cutting，canine teeth，molar and denture

3 結語

文中提出了可對任意形狀物體的點云數據進行法向量模糊估值的算法，利用模糊推理系統將點云數據進行分類，將物體可能含有薄片特征和尖銳特征的點云區域區分處理:對于平滑的區域直接應用PCA 方法估值其法向量，而含有薄片特征的區域由檢測器對關鍵點的法向量進行推導，不規則尖銳特征的區域由增添附加點的算法對法向量估值.這樣使檢測器和附加點算法只應用在局部區域，減少了計算量，而且針對不同區域應用不同算法的法向量估值準確率和曲面擬合效率均較高.

［1］Hoppe H，Tony D R，Tom D.Surface reconstruction from unorganized points［C］∥Proceedings of the 19th Annual Conference on Computer Graphics and Interactive Techniques.New York:ACM SIGGRAPH Computer Graphics，1992:71-78.

［2］Huang H，Li D，Zhang H，et al.Consolidation of unorganized point clouds for surface reconstruction ［J］.ACM Transactions on Graphics，2009，28(5):1-7.

［3］Jaakko L，Matthias Z，Emmanuel T，et al.A meshless hierarchical representation for light transport ［J］.ACM Transactions on Graphics，2008，27(3):1-37.

［4］Amenta N，Bern M.Surface reconstruction by Voronoi filtering［J］.Discrete ＆ Computational Geometry，1998，22(3):481-504.

［5］Alliez P，David C S，Tong Y Y，et al.Voronoi-based variational reconstruction of unoriented point sets ［C］∥Proceedings of the Fifth Eurographics Symposium on Geometry Processing.Barcelona:Eurographics Association，2007:39-48.

［6］Li B，Schnabel R W，Klein R，et al.Robust normal estimation for point clouds with sharp features［J］.Computer ＆Graphics，2010，34(2):94-106.

［7］Chen Y L，Lai S H.An orientation inference framework for surface reconstruction from unorganized point clouds［J］.IEEE Transactions on Image Processing，2011，20(3):762-775.

［8］Wang J，Yang Z W，Chen Fa-lai.A variational model for normal computation of point clouds［J］.The Visual Computer，2012，28(2):163-174.

［9］Liu Y J，Zhang Y D.A novel normal estimation based on fuzzy inference for mass point clouds of denture in the three-dimensional reconstruction ［J］.ICIC Express Letters Part B:Applications，2012，3(4):931-938.

［10］?？?，達飛鵬，蔡亮.基于特征融合的三維人臉識別［J］.東南大學學報:自然科學版，2011，41(1):47-51.Chang Jun-yan，Da Fei-peng，Cai Liang.3D face recognition based on feature fusion ［J］.Journal of Southeast University:Natural Science Edition，2012，41(1):47-51.

［11］Yukie N，Yutaka O，Hiromasa S.Smoothing of partition of unity implicit surfaces for noise robust surface reconstruction［J］.Computer Graphics Forum，2009，28(5):1339-1348.

［12］Taku I.Surface reconstruction with an interactive modification of point normals［C］∥Proceedings of the 2nd International Conference on Computational Modeling of Objects Represented in Images.Buffalo:CESER Publications，2010:263-274.

［13］Taku I，Yoshifumi K.Surface reconstruction from 3D scattered points with normals using both delaunay tetrahedralization and implicit function［C］∥Proceedings of the 7th Asia Simulation Conference on System Simulation and Scientific Computing.Beijing:Institute of Electrical and Electronics Engineers，2010:780-785.

［14］Stelldinger P，Latecki J L，Siqueira M.Topological equivalence between a 3D object and the reconstruction of its digital image［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2007，29(1):126-140.

［15］Lange C，Polthier K.Anisotropic smoothing of point sets［J］.Computer Aided Geometric Dessign，2005，22(7):680-692.