李 治

1珠海市測繪院，廣東 珠海，519015

在工業(yè)、工程檢測和逆向工程中，曲面擬合是經常遇到的問題。一般實物在空間的模型都有著嚴格的數學公式表達，通過掃描獲取實物表面點的三維坐標，可以擬合出實物表面在空間坐標系中的幾何方程［1，2］，這對于發(fā)現實物的整體形變以及生成模型設計CAD圖紙是至關重要的。谷川等［3，4］提出了采用遺傳算法進行曲面擬合，并且通過工程實例說明了該算法的有效性，但是對于柱面這一類的空間曲面進行坐標平移旋轉以及標準曲面表達參數的求取存在一定的困難。張益澤等［5］利用圓柱的幾何原理列誤差方程，并對其做平方變換，相較于谷川的方法較為簡單，但要加兩個限制條件，而且要迭代兩次。田曉等［6］利用坐標轉換的方法達到減少了參數的目的，但也僅僅是對圓柱面進行了擬合。王解先［7，8］針對工業(yè)測量中的二次曲面擬合提出了按特征值、旋轉角和平移量為參數擬合二次曲面，其對于旋轉矩陣是分別求3個旋轉角的3個旋轉矩陣的乘積，不能直接得到旋轉矩陣。據此，本文提出了基于豪斯霍爾德變換的三維空間柱面擬合法，豪斯霍爾德變換可以直接得到旋轉矩陣，將提高擬合的效率，通過豪斯霍爾德變換將柱面的軸向轉換到與測量坐標系Z軸平行，點云投影到二維平面，再擬合投影的曲線，最后拓展到三維空間達到柱面擬合，為后面點云三維建模中一般柱面的建模提供一種新的方法，并結合實例說明該方法可行且精度可靠。

1 三維空間柱面擬合的數學模型

常見的三維空間柱面包括圓柱面、橢圓柱面、拋物柱面和雙曲柱面［9］。這些柱面在一個測站的數據在二維平面都可以用式（1）來表達。

一般地，在Rn中，將向量x映射為關于“與單位向量u正交的n-1維子空間”對稱的向量y的鏡像變換定義如下：設單位向量ω∈Rn，稱R=E-2ωωT為豪斯霍爾德矩陣（初等反射矩陣），由豪斯霍爾德矩陣確定的線性變換成為豪斯霍爾德變換（初等反射變換、鏡像變換或反射變換）［10］?？臻g中的三維柱面的姿態(tài)都是隨機的，將柱面的軸向利用豪斯霍爾德變換至與測量坐標系中的Z軸平行，就可以將柱面投影至XOY平面，利用點的X、Y坐標擬合曲線。最終表達一個柱面模型的參數分別是（α，β，a0，a1，a2…an），其中α、β分別是柱面軸向與X、Y軸 的 夾 角，a0、a1、a2…an是 二 維 曲 線 表 達 式 的系數。

2 三維空間柱面擬合

2.1 計算法向量

對于求解法向量的方法有很多，本文采用主成分分析法（principal component analysis，PCA），該方法利用搜索查詢點的k（k≥3）鄰域構造協(xié)方差矩陣，協(xié)方差矩陣反映了該鄰域局部曲面的變化規(guī)律，由此可以求解該矩陣特征值和特征向量來估計查詢點的法向量。

查詢點pi的協(xié)方差矩陣C為：

利用雅克比方法求協(xié)方差矩陣C，其特征值為λ1、λ2、λ3(λ1≤λ2≤λ3)，其對應的特征向量為n1、n2、n3。3個向量中最小特征值對應的特征向量n1就是查詢點的法向量。為了方便，對于所有的法向量都標準化處理。

2.2 計算軸向向量

利用柱面軸向與柱面點的法向量都垂直的特點求取柱面軸向向量l=(n，m，t)，利用克拉默法則求取l，對 于 兩 個 不 共 線 上 向 量a=(x1，y1，z1)，b=(x2，y2，z2)。

2.3 軸向轉換

式中，R為旋轉矩陣，由豪斯霍爾德變換求得R.

2.4 曲線擬合

在得到旋轉矩陣R后，將點云通過旋轉調整至與x、y軸平面垂直的姿態(tài)。此時柱面物體在x、y軸平面的投影就是一條曲線。由點云的x，y坐標通過最小二乘求多項式系數就能把柱狀物體的大小確定下來了。

3 實例應用與分析

3.1 規(guī)則柱面擬合

本文采用Riegl VZ-1000對某建筑物傾斜圓柱體進行掃描，Riegl VZ-1000是脈沖式掃描，標稱作業(yè)距離為1 400 m，單次掃描精度為5 mm/100 m，本次數據采集的距離小于100 m，故點云的精度在5 mm以內。在掃描得到的點云數據中裁剪出研究對象點云，利用統(tǒng)計分析技術剔除點云數據中的離群點，計算每個點到它所有鄰域內點的平均距離，假設得到的結果是一個高斯分布，其形狀由均值和標準差定義，平均距離在標準范圍（有全局距離平均值和方差定義）之外的點，可被定義為離群點并剔除，最終采用體素格網化法進行下采樣后點云如圖1所示，計算點云的法向量之后得到軸向量為（0.418，0.114，0.485），經過豪斯霍爾德轉換得到旋轉矩陣R。

圖1 原始圓柱點云Fig.1 Original Cylindrical Point Cloud

通過坐標轉換將傾斜圓柱點云數據轉換值與xoy平面垂直如圖2所示，將點云數據投影至xoy平面，利用最小二乘算法求取圓的圓心和半徑，得到圓的表達式為：(X+456.828)2+(Y-214.149)2=29.8152，最終圓柱面擬合如圖3所示。

圖2 旋轉后圓柱點云Fig.2 Cylindrical Point Cloud After Rotation

圖3 擬合圓柱面Fig.3 Fitting a Cylindrical Surface

對擬合后的圓柱面利用旋轉矩陣進行姿態(tài)調整，旋轉回初始姿態(tài)，利用初始點云數據與擬合后旋轉柱面參數進行誤差分析，計算每個點至圓柱面的距離，如圖4所示，除少部分噪聲點影響，大部分的誤差控制在0.05 cm以內。

圖4 擬合殘差Fig.4 Fitting Residuals

3.2 一般柱面擬合

本文采用Riegl VZ-1000對某大學行政樓掃描的部分數據進行驗證，如圖5所示，部分點云數據如表1所示（已經過重心化處理）。計算點云的法向量之后得到軸向量為（-0.443，-0.110，-0.512），經過豪斯霍爾德轉換得到旋轉矩陣R。

圖5 行政樓影像圖Fig.5 Image of the Administrative Building

表1 部分點云數據Tab.1 Some Point Cloud Data

由 此 得 到α=arccos0.985=9.93°，β=arccos0.761=40.45°。通過坐標轉換以后利用X、Y坐標擬合曲線，采用5階至10階多項式對坐標轉換后投影至xoy平面點云數據進行擬合，經過實驗采用10階多項式擬合效果最好，各項系數（見表2）如下：a0=-1.034×10-6，a1=3.18×10-7，a2=5.98×10-5，a3=-0.000 1，a4=-0.000 9，a5=0.002，a6=0.007，a7=-0.03，a8=0.046，a9=0.01，a10=1.582，旋轉后點云如圖6所示，其中擬合的曲線如圖7所示。

圖6 旋轉后點云Fig.6 Point Cloud After Rotation

圖7 擬合曲線Fig.7 Fitting Curve

表2 多項式系數Tab.2 Polynomial Coefficients

最終擬合的曲面如圖8所示。最后根據擬合點云到擬合曲面的距離評定擬合精度，空間距離投影到xoy平面上就是點到曲線的距離，最后得到的距離殘差圖如圖9所示，絕大部分點云偏離曲面都在0.03 cm以內。擬合的精度滿足一般的工程曲面擬合精度。

圖8 擬合曲面Fig.8 Fitting Surface

圖9 擬合殘差Fig.9 Fitting Residuals

4 結束語

針對空間中任意姿態(tài)一般的柱面，提出了先計算其軸向向量，通過豪斯霍爾德變換進行投影，降維處理來達到擬合一般曲面的曲線公式，在從二維擴展到三維空間。并以實例進行驗證，取得了較好的效果。在曲線擬合的時候對于不同的曲面需要綜合考慮曲線的階數，以達到最佳效果。