基于向量旋轉(zhuǎn)的擬合空間圓方法研究

張振虎

上海勘察設(shè)計(jì)研究院（集團(tuán)）有限公司 上海 200082

引言

隨著科技的發(fā)展、信息技術(shù)的進(jìn)步，以及精密工程技術(shù)的發(fā)展，產(chǎn)生了很多新的測(cè)量技術(shù)，三維測(cè)量建模技術(shù)就是時(shí)代的產(chǎn)物。在核電廠、光源儲(chǔ)存環(huán)、離子加速器實(shí)驗(yàn)裝置等項(xiàng)目建造過程中，三維設(shè)計(jì)、設(shè)備制造、核心設(shè)備精密安裝及重要設(shè)備在線監(jiān)測(cè)均進(jìn)行了大量建模，建模技術(shù)的采用促使高精度測(cè)量工作日益轉(zhuǎn)向可視化、智能化及精密化。逆向建模作為一種最常用的建模方法得到了越來越廣泛的應(yīng)用，它主要通過激光掃描、激光跟蹤或攝影測(cè)量等技術(shù)，獲得目標(biāo)物體的三維測(cè)量數(shù)據(jù)，然后根據(jù)測(cè)量數(shù)據(jù)構(gòu)建三維高精度模型。以核電廠為例，其在建造過程中采用了激光跟蹤技術(shù)對(duì)一回路主要設(shè)備，如反應(yīng)堆壓力容器、蒸汽發(fā)生器、主泵、主管道及波動(dòng)管等進(jìn)行了精密建模，這些設(shè)備的接口多為圓形，而這些設(shè)備接口的圓形并不在嚴(yán)格的水平面上，因此，對(duì)這些設(shè)備建模通常需要對(duì)空間圓進(jìn)行擬合。

1 數(shù)學(xué)模型及算法

空間圓無直接公式描述，基于公式擬合的方法很難實(shí)施，故空間圓參數(shù)擬合一般采取間接方式進(jìn)行。根據(jù)參考文獻(xiàn)[1-5]，擬合空間圓主要采用平面與圓球相截法。由圓球性質(zhì)可以得知，任何一個(gè)平面與圓球相截得到的截面均為圓形，故平面與圓球相截法是通過擬合平面與圓球得到空間圓參數(shù)。其中，文獻(xiàn)[1-4]采取擬合平面和一個(gè)圓球中心位于擬合平面上的圓球得到空間圓參數(shù)；文獻(xiàn)[5]也采取擬合平面和圓球方式進(jìn)行，該方法不要求圓球圓心位于擬合平面上，在平面和圓球擬合后將球心投影到擬合平面上以得到圓心，過程較文獻(xiàn)[1-4]復(fù)雜。擬合空間圓需要確定7個(gè)參數(shù)，包括空間圓所在平面的法向量(a，b， c)、空間圓的圓心(x0， y0， z0)和空間圓的半徑R。本文通過擬合平面求取平面法向量，構(gòu)建旋轉(zhuǎn)矩陣，然后將空間圓擬合問題轉(zhuǎn)化為平面圓擬合問題。在工程實(shí)例中，圓形截面的特征數(shù)據(jù)采集通常采用全站儀、激光掃描儀或激光跟蹤儀完成。部分物體圓形截面因受視線阻隔往往難以一次測(cè)量就采集到所有特征數(shù)據(jù)，通常需要轉(zhuǎn)站1～2次完成。轉(zhuǎn)站需要進(jìn)行大角度坐標(biāo)轉(zhuǎn)換，故在擬合空間圓前需要將多站的采集數(shù)據(jù)轉(zhuǎn)換到同一個(gè)坐標(biāo)系中。

綜上，本文算法主要分為以下五步：①大角度坐標(biāo)轉(zhuǎn)換；②擬合空間圓所在平面求取平面法向量；③根據(jù)平面法向量構(gòu)建旋轉(zhuǎn)矩陣；④將擬合坐標(biāo)點(diǎn)旋轉(zhuǎn)，求取平面圓圓心和半徑；⑤平面圓圓心旋轉(zhuǎn)得到空間圓圓心。

1.1 大角度坐標(biāo)轉(zhuǎn)換

坐標(biāo)轉(zhuǎn)換模型可表示為：

其中，（X'、Y'、Z'）為轉(zhuǎn)換前坐標(biāo)，（X、Y、Z）為轉(zhuǎn)換后坐標(biāo)，λ為尺度系數(shù)，R為含三個(gè)旋轉(zhuǎn)角α、β、γ的旋轉(zhuǎn)矩陣，△X、△Y、△Z為平移參數(shù)。R=RxRyRz，Rx是把原坐標(biāo)繞X軸旋轉(zhuǎn)α角得到的旋轉(zhuǎn)矩陣，Ry是把原坐標(biāo)繞Y軸旋轉(zhuǎn)β角得到的旋轉(zhuǎn)矩陣，Rz是把原坐標(biāo)繞Z軸旋轉(zhuǎn)γ角得到的旋轉(zhuǎn)矩陣，α、β、γ均為繞軸逆時(shí)針旋轉(zhuǎn)。

整理以后，旋轉(zhuǎn)矩陣R為以下形式：

兩坐標(biāo)系之間轉(zhuǎn)換需要求取7個(gè)轉(zhuǎn)換參數(shù)，一般需要三個(gè)以上的公共點(diǎn)坐標(biāo)來進(jìn)行轉(zhuǎn)換。

為便于數(shù)據(jù)處理，一般需要對(duì)所使用的公共點(diǎn)坐標(biāo)進(jìn)行重心化處理，將兩套坐標(biāo)系的坐標(biāo)均轉(zhuǎn)化為重心化坐標(biāo)。根據(jù)參考文獻(xiàn)[8]，重心化坐標(biāo)計(jì)算公式如下：

將公共點(diǎn)的中心化坐標(biāo)代入式（1），可得：

從上式可以看出，坐標(biāo)轉(zhuǎn)換參數(shù)計(jì)算時(shí)需要先計(jì)算尺度系數(shù)和旋轉(zhuǎn)矩陣，然后代入式（8）就可以得到平移參數(shù)。

根據(jù)參考文獻(xiàn)[9]，對(duì)式（7）兩邊取2-范數(shù)，由于R為正交矩陣，及λ＞0，可得

式（7）變化為以下形式：

求取正交矩陣R，使A-λRB的Frobenius范數(shù)最小時(shí)所得的R即為最佳轉(zhuǎn)換矩陣，即求解以下問題：

將上式進(jìn)行跡函數(shù)計(jì)算并展開，得：

根據(jù)參考文獻(xiàn)[6]，當(dāng)且僅當(dāng)R =VUT時(shí)，式（14）得到極小值。得到旋轉(zhuǎn)矩陣后，將旋轉(zhuǎn)矩陣代入到公式（8），就可以得到平移參數(shù)。

1.2 平面法向量計(jì)算

將所有坐標(biāo)采集數(shù)據(jù)轉(zhuǎn)換至同一坐標(biāo)系后，就可以開始計(jì)算空間圓所在平面的法向量計(jì)算。

平面方程一般可表示為：

為減少參數(shù)的計(jì)算，可以將方程簡(jiǎn)化，兩側(cè)同除以D（D≠0），方程可以簡(jiǎn)化為：

因方程的三個(gè)參數(shù)之間存在著線性關(guān)系，故可以采用求解線性方程組參數(shù)的方法進(jìn)行結(jié)算。假設(shè)空間圓上的n個(gè)觀測(cè)點(diǎn)為( xi， yi， zi)（1≤i≤n），以矩陣形式表示如下：

上式可以簡(jiǎn)化為：

其最小二乘解為：

上述算法為擬合平面最常用算法，但該算法并未考慮系數(shù)矩陣N 中存在的誤差，且造成了固定項(xiàng)l 中數(shù)據(jù)的改動(dòng)，得到的結(jié)果必然不是最優(yōu)的。

假設(shè)系數(shù)矩陣N中存在獨(dú)立同分布誤差，根據(jù)參考文獻(xiàn)[6]，其參數(shù)求解可以采用以下算法。

將( xi， yi， zi)重心化，得到重心坐標(biāo)重心化公式如下：

最小奇異值所對(duì)應(yīng)的V向量[ a 'b ' c ']T即為平面方程中參數(shù)A、B、C的值。平面方程中常數(shù)項(xiàng)的值可按下式計(jì)算：

根據(jù)文獻(xiàn)[7]研究，基于奇異值分解的擬合平面方法與基于正交距離擬合平面的方法是等價(jià)的，其求得的平面參數(shù)可以使點(diǎn)到平面的距離平方和最小，其比式(18)—(21)采用的擬合平面常用算法更為精確。

1.3 求取圓心和半徑

（1）構(gòu)建旋轉(zhuǎn)矩陣

假設(shè)空間圓的圓心坐標(biāo)為( x0， y0， z0)，其應(yīng)位于擬合的平面上，故有：

應(yīng)用羅德里格旋轉(zhuǎn)公式矩陣形式，將平面法向量旋轉(zhuǎn)至與Z軸平行，其構(gòu)成的旋轉(zhuǎn)矩陣如下：

其旋轉(zhuǎn)公式表達(dá)為：

（2）擬合平面圓

平面圓方程可以用下式表示：

采用最小二乘法可以得到E、F、G的值，平面圓圓心及半徑值計(jì)算如下：

（3）解算空間圓圓心

1.4 擬合空間圓效果評(píng)價(jià)

計(jì)算擬合點(diǎn)到空間圓的法向偏差：

計(jì)算擬合點(diǎn)到擬合圓的徑向偏差：

計(jì)算擬合點(diǎn)到擬合圓的空間偏差：

計(jì)算RMSE（均方根誤差）：

2 工程算例

某設(shè)備中部的法蘭需要與其他設(shè)備法蘭進(jìn)行精密對(duì)接，為檢測(cè)設(shè)備的匹配性，需要在安裝前進(jìn)行模擬對(duì)接。為精確模擬，采用激光跟蹤儀對(duì)設(shè)備的法蘭圓面進(jìn)行了測(cè)量。因視線阻擋，難以在一站就完成所有數(shù)據(jù)的采集，故法蘭圓面的測(cè)量采取了轉(zhuǎn)站測(cè)量的方式，兩站之間采用公共點(diǎn)進(jìn)行連接。采集的原始數(shù)據(jù)請(qǐng)見下表1。

表1 某設(shè)備中部法蘭圓第一站測(cè)量數(shù)據(jù)（單位：mm）

表2 某設(shè)備中部法蘭圓第二站測(cè)量數(shù)據(jù)（單位：mm）

表3 公共點(diǎn)測(cè)量數(shù)據(jù)（單位：mm）

表4 坐標(biāo)轉(zhuǎn)換參數(shù)

采用本文坐標(biāo)轉(zhuǎn)換算法，計(jì)算第二站轉(zhuǎn)入第一站的坐標(biāo)轉(zhuǎn)換參數(shù)。因激光跟蹤儀的精度很高，且測(cè)量距離較短，其坐標(biāo)轉(zhuǎn)換尺度系數(shù)可視為1。計(jì)算出來的旋轉(zhuǎn)矩陣及平移參數(shù)請(qǐng)見表4。

利用表4的坐標(biāo)轉(zhuǎn)換參數(shù)將第二站的測(cè)量數(shù)據(jù)轉(zhuǎn)入到第一站測(cè)量坐標(biāo)系中，并計(jì)算擬合圓的參數(shù)。根據(jù)本文算法，對(duì)空間圓的數(shù)據(jù)進(jìn)行了擬合，擬合后的平面方程單位法向量為（0.93393557，-0.357440912，0.000587606），平面方程常數(shù)項(xiàng)為-3541.36949，空間圓圓心坐標(biāo)為（2605.2329，-3100.094608，253.8655206），半徑為589.6816919mm。根據(jù)擬合結(jié)果計(jì)算出的擬合后各點(diǎn)的法向偏差、徑向偏差及空間偏差見表5。

表5 某設(shè)備端面圓擬合后偏差統(tǒng)計(jì)（單位：mm）

采用式(40)，根據(jù)表5中所示的各點(diǎn)的法向偏差、徑向偏差及空間偏差計(jì)算均方根誤差（RMSE）。經(jīng)計(jì)算，法向偏差RMSE為0.042167mm，徑向偏差RMSE為0.409195mm，空間偏差RMSE為0.411362mm，說明擬合圓的各項(xiàng)偏差均很小，空間圓參數(shù)解算較準(zhǔn)確。

為驗(yàn)證本文中算法，使用文獻(xiàn)[2]中的算法對(duì)本文數(shù)據(jù)進(jìn)行擬合。為便于擬合后參數(shù)的對(duì)比，本文將采用文獻(xiàn)[2]中算法解算出的平面法向量進(jìn)行了單位化。文獻(xiàn)[2]中算法和本文算法解算的結(jié)果如表6所示。

表6 文獻(xiàn)[2]算法和本文算法解算結(jié)果（單位：mm）

通過表6可發(fā)現(xiàn)，本文所述算法擬合空間圓后法向偏差RMSE和文獻(xiàn)[2]中使用的算法結(jié)果基本一致，但徑向偏差RMSE及空間偏差RMSE均略優(yōu)于文獻(xiàn)[2]中使用的算法，說明本文算法可行，精度較高。

3 結(jié)束語

本文通過奇異值分解法擬合平面，構(gòu)建旋轉(zhuǎn)矩陣化空間圓擬合為平面圓擬合，原理清晰，不用進(jìn)行迭代，且有利于編程計(jì)算，并通過實(shí)踐證實(shí)，該方法實(shí)用且可靠。適合用于核電廠、光源存儲(chǔ)環(huán)實(shí)驗(yàn)裝置、直線加速器實(shí)驗(yàn)裝置、飛機(jī)、造船、汽車、航空航天等行業(yè)大型設(shè)備模型構(gòu)建過程中圓形截面的擬合。