羅德里格矩陣在盾構(gòu)姿態(tài)解算中的應(yīng)用

魏云豹 劉恒杰 周小利

摘 要：在大量實(shí)際隧道工程中，采用羅德里格矩陣進(jìn)行盾構(gòu)機(jī)軸線精確解算，其精度可達(dá)到毫米級(jí)，完全可滿足日常隧道施工需求。由此，本文介紹了利用羅德里格矩陣模型和反對(duì)稱矩陣，把傳統(tǒng)七參數(shù)模型中的三個(gè)旋轉(zhuǎn)參數(shù)用反對(duì)稱矩陣中三個(gè)獨(dú)立參數(shù)代替，使用Excel表格工具解算出盾構(gòu)機(jī)軸線坐標(biāo)。

關(guān)鍵詞：盾構(gòu)機(jī)軸線坐標(biāo);羅德里格矩陣;三維坐標(biāo)轉(zhuǎn)換

中圖分類號(hào)：U452.13 文獻(xiàn)標(biāo)識(shí)碼：A 文章編號(hào)：1003-5168（2019）28-0076-04

Application of Rodrigue's Matrix in Attitude Calculation of Shield Machine

WEI Yunbao LIU Hengjie ZHOU Xiaoli

（China Railway Engineering Equipment Co.， Ltd.，Zhengzhou Henan 450000）

Abstract： In a large number of practical tunnel projects， Rodrigo's matrix is used to accurately calculate the axis of shield machine， and its accuracy can reach millimeter level， which can fully meet the needs of daily tunnel construction. Therefore， this paper introduced using Rodrigo's matrix model and antisymmetric matrix， replacing three rotation parameters in the traditional seven parameter model with three independent parameters in antisymmetric matrix， and using Excel table tool to solve the axis coordinates of shield machine.

Keywords： shield machine axis coordinate;Rodrigo's matrix;3-D coordinate conversion

盾構(gòu)機(jī)是一種被廣泛應(yīng)用于城市地鐵建設(shè)的大型工程機(jī)械。精確計(jì)算盾構(gòu)機(jī)軸線坐標(biāo)精度是保證其正常掘進(jìn)的關(guān)鍵。本文主要以羅德里格矩陣計(jì)算盾構(gòu)機(jī)軸線精度。

1 數(shù)學(xué)模型

零位坐標(biāo)系中，盾構(gòu)機(jī)的軸線坐標(biāo)和特征點(diǎn)坐標(biāo)是已知的。在施工過程中，當(dāng)?shù)刈鴺?biāo)系下的特征點(diǎn)坐標(biāo)是可以通過全站儀實(shí)測(cè)得到的，這樣就可以通過公共特征點(diǎn)建立2個(gè)坐標(biāo)系的空間關(guān)系式，即求解7參數(shù)，包括3個(gè)平移參數(shù)、3個(gè)旋轉(zhuǎn)參數(shù)、1個(gè)尺度參數(shù)（[Δx]、[Δy]、[Δz]、[θ1]、[θ2]、[θ3]、[K]），進(jìn)而解算出盾構(gòu)機(jī)的軸線坐標(biāo)[1-4]。

根據(jù)坐標(biāo)轉(zhuǎn)換的物理過程，可得到如式（1）所示的數(shù)學(xué)模型：

[xyzi=kΔxΔyΔz+kRxyzI]? ? ? ? ? ? ? ? （1）

根據(jù)羅德里格矩陣原理，具有3個(gè)自由度的正交旋轉(zhuǎn)矩陣[R]可以看成是由反對(duì)稱矩陣[S]構(gòu)成的羅德里格矩陣：

[R=I+SI-S-1]? ? ? ? ? ? ? ? ? ? ?（2）

[S=0-c-bc0-aba0]? ? ? ? ? ? ? ? ? ? ? ? ? ? （3）

式中，[a]、[b]、[c]是3個(gè)未知參數(shù)，用于替代3個(gè)旋轉(zhuǎn)參數(shù)，簡化計(jì)算過程，由式（1）可知2個(gè)公共點(diǎn)可以列6個(gè)方程式。

在7個(gè)參數(shù)模型中，[R]為3個(gè)旋轉(zhuǎn)矩陣的乘積：

[R=RxαRyβRzγ=1000cosαsinα0-sinαcosαcosβ0-sinβ010sinβ0cosβcosγsinγ0-sinγ-cosγ0001]? ? ? ? ? ?（4）

空間中有兩個(gè)不同的坐標(biāo)系，即零位坐標(biāo)系[Oi-xyzxiyizi]和當(dāng)?shù)刈鴺?biāo)系[OI-xyzxIyIzI]，同時(shí)有[n]個(gè)公共點(diǎn)[i，I=1，2，3，…]。

由于在實(shí)際盾構(gòu)機(jī)軸線坐標(biāo)解算中尺度變換為零，故[k=1]，這樣只需要解算出旋轉(zhuǎn)參數(shù)和平移參數(shù)即可。

旋轉(zhuǎn)參數(shù)解算：將公共點(diǎn)[ni]與[nI]在[Oi-xyzxiyizi]和[OI-xyzxIyIzI]中的坐標(biāo)值帶入（1）式得：

[xi'yi'zi'=kΔxΔyΔz+kRxiyizi]? ? ? ? ? ? ? ? ? ? （5）

[xI'yI'zI'=kΔxΔyΔz+kRxIyIzI]? ? ? ? ? ? ? ? ? ? ?（6）

將式（5）和式（6）相減，消去平移參數(shù)，整理后可得：

[xI'-xi'yI'-yi'zI'-zi'=kRxI-xiyI-yizI-zi]? ? ? ? ? ? ? ? ? ? ? ? ?（7）

將式（2）與式（7）結(jié)合，整理后為：

[0-zIi-zIi'-yIi-yIi'-zIi-zIi'0xIi+xIi'yIi+yIi'xIi+xIi'0abc=xIi'-xIiyIi'-yIizIi'-zIi]? ? ? ? ? ? （8）

由于式（8）的系數(shù)陣為反對(duì)稱矩陣，行列式為0，秩等于2，式（8）中只有兩個(gè)獨(dú)立方程，無法解出3個(gè)參數(shù)[a]、[b]、[c]，因此需要將3個(gè)公共點(diǎn)[n1]、[n2]、[n3]代入式（1），按照同樣原理與式（8）聯(lián)合解出未知參數(shù)[a]、[b]、[c]，方程組如下：

[0-z21-z21'-y21-y21'-z21-z21'0z21+z21'y21+y21'x21+x21'00-z31-z31'-y31-y31'-z31-z31'0z31+z31'y31+y31'x31+x31'0abc=x21'-x21y21'-y21z21'-z21x31'-x31y31'-y31z31'-z31]? ? ? ? （9）

根據(jù)式（9）即可解出未知參數(shù)[a]、[b]、[c]，再結(jié)合式（2）可解出旋轉(zhuǎn)矩陣[R]。

平移參數(shù)計(jì)算：根據(jù)上述公式解出旋轉(zhuǎn)矩陣[R]，代入式（1），即可解出平移參數(shù)[Δx]、[Δy]、[Δz]。同理，求出每個(gè)點(diǎn)的平移參數(shù)，再計(jì)算平均值。

2 計(jì)算步驟

第一步：把零位坐標(biāo)系和當(dāng)?shù)刈鴺?biāo)系下的公共點(diǎn)（特征點(diǎn)）輸入表1中。

第二步：根據(jù)式（9）列出矩陣列表，如表2所示。

第三步：用Excel中的TRANSPOSE命令算出[B]矩陣的轉(zhuǎn)置矩陣[BT]。首先在Excel中確定矩陣[B]的范圍（列數(shù)和行數(shù)），然后在空白區(qū)域選擇矩陣[B]的行數(shù)和列數(shù)，然后在選定范圍后輸入轉(zhuǎn)置矩陣函數(shù)，例如：選中范圍后輸入函數(shù)=TANSPOSE（B2：D19），然后按Ctrl+Shift+Enter鍵，即可求出轉(zhuǎn)置矩陣[BT]。

第四步：用Excel中的MMULT命令算出[BTB]。首先在Excel表確定乘積矩陣的行數(shù)和列數(shù)范圍。然后在輸入欄中輸入矩陣乘積函數(shù)，例如，=Mmult（F4：W6，B2：D19），然后按Ctrl+Shift+Enter鍵，得出[BTB]。

第五步：用Excel中的MMULT命令算出[BTL]。

第六步：用Excel中的MINVERSE命令算出[（BTB）-1]。

第七步：用Excel中的MMULT命令求轉(zhuǎn)換參數(shù)[a]、[b]、[c]。

[abc=（BTB）（BTB）-1]? ? ? ? ? ? ? ? ? ? ? ? ? ?（10）

第八步：根據(jù)式（1）算出每個(gè)點(diǎn)的平移參數(shù)[Dx]、[Δy]、[Δz]后取平均值。

運(yùn)行界面見圖1。

3 工程應(yīng)用驗(yàn)算

本算例中，盾構(gòu)機(jī)軸線坐標(biāo)和特征點(diǎn)坐標(biāo)為現(xiàn)場(chǎng)實(shí)際坐標(biāo);零位坐標(biāo)系為任意設(shè)站轉(zhuǎn)換為機(jī)器坐標(biāo)系（盾構(gòu)機(jī)盾首坐標(biāo)[x=0]、[y=0]、[z=0]，盾尾坐標(biāo)[x=-5.4]、[y=0]、[z=0]）;施工坐標(biāo)系為當(dāng)?shù)刈鴺?biāo)系。

3.1 模型建立

零位坐標(biāo)系下的坐標(biāo)，按CAD法整體轉(zhuǎn)換到當(dāng)?shù)刈鴺?biāo)系中，然后提取各特征點(diǎn)坐標(biāo)和盾構(gòu)機(jī)軸線坐標(biāo)，特征點(diǎn)和盾構(gòu)機(jī)軸線坐標(biāo)相對(duì)空間關(guān)系固定有唯一解。

3.2 驗(yàn)算過程

隨機(jī)在當(dāng)?shù)刈鴺?biāo)系（[x]，[y]，[z]）下加入-0.003～0.003m的誤差（實(shí)際工作中，誤差不能超過0.003m限值）;圖2為羅德里格矩陣解與真值比較之差（盾尾），圖3為羅德里格矩陣解與真值比較之差（盾首），數(shù)據(jù)采集為69組。通過分析可知在限差范圍-0.003～0.003m內(nèi)，盾尾各坐標(biāo)（[x]，[y]，[z]）波動(dòng)都在-0.004 7～0.004 8m，盾首各坐標(biāo)（[x]，[y]，[z]）波動(dòng)都在-0.004 7～0.004 8m，可滿足日常施工要求。

4 結(jié)語

在大量實(shí)際隧道工程中，采用羅德里格矩陣進(jìn)行盾構(gòu)機(jī)軸線精確解算，其精度可達(dá)到毫米級(jí)，完全可滿足日常隧道施工需求（隧道施工測(cè)量要求為公分級(jí)）。同時(shí)，把羅德里格矩陣解算過程通過Excel表格形式逐步展示出來，可以為其他研究人員提供良好的參考價(jià)值。

參考文獻(xiàn)：

[1]姚吉利，韓保民，楊元喜.羅德里格矩陣在三維坐標(biāo)轉(zhuǎn)換嚴(yán)密解算中的應(yīng)用[J].武漢大學(xué)學(xué)報(bào)（信息科學(xué)版），2006（12）：63-65，88.

[2]潘國榮，周躍寅.兩種坐標(biāo)系轉(zhuǎn)換計(jì)算方法的比較[J].大地測(cè)量與地球動(dòng)力學(xué)，2011（3）：62-66.

[3]潘國榮，汪大超，周躍寅.兩種大轉(zhuǎn)角空間坐標(biāo)轉(zhuǎn)換模型研究[J].山東科技大學(xué)學(xué)報(bào)（自然科學(xué)版），2015（1）：61-67.

[4]原玉磊，蔣理興，劉靈杰.羅德里格矩陣在坐標(biāo)系轉(zhuǎn)換中的應(yīng)用[J].測(cè)繪科學(xué)，2010（2）：178-179.