加權(quán)整體最小二乘坐標(biāo)匹配算法在機(jī)場(chǎng)道面測(cè)量中的應(yīng)用

戴中東孟 良高永攀項(xiàng) 偉

1空軍研究院工程設(shè)計(jì)研究所，北京，100068

2北京特種工程設(shè)計(jì)研究院，北京，100028

機(jī)場(chǎng)道面高程圖（也稱為板角高程圖）測(cè)量一般在機(jī)場(chǎng)改擴(kuò)建設(shè)計(jì)前或者機(jī)場(chǎng)建設(shè)竣工后進(jìn)行。板角高程測(cè)量精度至毫米［1］，一般尾數(shù)為估讀，且板角點(diǎn)相對(duì)臨近水準(zhǔn)點(diǎn)高程中誤差不大于5 mm，極限誤差不大于1 cm，因此可認(rèn)為，板角高程測(cè)量的精度要求優(yōu)于1 cm。傳統(tǒng)的測(cè)量方法將平面數(shù)據(jù)和高程數(shù)據(jù)分開進(jìn)行采集，然后一一對(duì)應(yīng)合并（一個(gè)機(jī)場(chǎng)大概3萬個(gè)數(shù)據(jù)），過程繁瑣，工作效率低下，且很容易出現(xiàn)錯(cuò)誤。地面三維激光掃描儀（terrestrial laser scanner，TLS）發(fā)展已經(jīng)30多年，特別是近10年來，精度高、測(cè)速快的掃描儀已經(jīng)逐漸走進(jìn)各個(gè)應(yīng)用領(lǐng)域［2］，比如工程測(cè)量［3］、文物保護(hù)［4］、工業(yè)測(cè)量［5］、交通道路測(cè)量［6-8］等，極大提高了工作效率。

但根據(jù)機(jī)場(chǎng)道面高程測(cè)量中的精度要求，采用常規(guī)掃描儀流程，主要存在以下3個(gè)主要問題：一是靶標(biāo)控制點(diǎn)的高程精度要足夠高，且測(cè)量過程中需要盡量減少“量高誤差”；二是由于道面測(cè)量對(duì)平面精度和高程精度要求不一致，特別是掃描儀沒有整平裝置時(shí)，平面坐標(biāo)誤差會(huì)影響高程匹配的精度，且配套軟件沒有對(duì)權(quán)值進(jìn)行設(shè)定；三是由于機(jī)場(chǎng)周邊地勢(shì)平坦，且無連續(xù)制高點(diǎn)，距離儀器越遠(yuǎn)，道面測(cè)量點(diǎn)的反射角就越小，大于70 m幾乎就沒有回波數(shù)據(jù)，這樣就使得相鄰測(cè)站的距離不能太大。因此采用常規(guī)的三維掃描儀測(cè)量方法，不僅效率低，而且也不能滿足厘米級(jí)以內(nèi)的高程精度。

本文首先分析了使用加權(quán)整體最小二乘［9-11］進(jìn)行坐標(biāo)匹配轉(zhuǎn)換的算法，該算法公式經(jīng)過推導(dǎo)和向量化后，采用LM（Levenberg-Marquardt）算法求出最優(yōu)轉(zhuǎn)換參數(shù)，并對(duì)目標(biāo)坐標(biāo)和測(cè)站坐標(biāo)的誤差來源進(jìn)行分析，確定了權(quán)值；其次，為了減少靶標(biāo)量高誤差和提高地面回波率，設(shè)計(jì)了一種控制點(diǎn)布設(shè)方案和測(cè)量方法；最后，通過實(shí)例數(shù)據(jù)的分析和比較，得出結(jié)論。

1 坐標(biāo)匹配算法和權(quán)值的確定

1.1 加權(quán)整體最小二乘法坐標(biāo)匹配算法

TLS外業(yè)測(cè)量時(shí)采用的是測(cè)站坐標(biāo)系統(tǒng)，而目標(biāo)坐標(biāo)為測(cè)圖時(shí)采用的坐標(biāo)系統(tǒng)，兩者之間需要轉(zhuǎn)換。由于測(cè)站坐標(biāo)和目標(biāo)坐標(biāo)均含有隨機(jī)誤差，直接采用經(jīng)典最小二乘法求解轉(zhuǎn)換參數(shù)時(shí)，求出的結(jié)果會(huì)有偏差［12］。因此需要使用變量中存在誤差的EIV（errors-in variables）模型，求出最優(yōu)估計(jì)參數(shù)解。從測(cè)站坐標(biāo)到目標(biāo)坐標(biāo)系統(tǒng)EIV模型公式為：

式中，[XY Z]T為目標(biāo)坐標(biāo)向量；[Vx Vy Vz]T為目標(biāo)坐標(biāo)誤差向量；[x y z]T表示測(cè)站坐標(biāo)向量；[vx vy vz]T為測(cè)站坐標(biāo) 誤差向 量，[ΔxΔyΔz]T為平移參數(shù)向量；λ代表尺度縮放參數(shù)，TLS的坐標(biāo)轉(zhuǎn)換中，一般認(rèn)為λ等于1；M為3個(gè)旋轉(zhuǎn)角組成的旋轉(zhuǎn)矩陣。

根據(jù)式（1）、式（2）可以得出：

由式（3）可以看出，在儀器整平的情況下，εx、εy為小角度，sinεx、sinεy也為極小值，cosεx、cosεy接近1，這時(shí)，目標(biāo)坐標(biāo)高程主要與測(cè)站坐標(biāo)的z方向坐標(biāo)分量、誤差、平移參數(shù)相關(guān)。但是，TLS一般沒有整平裝置，這時(shí)，目標(biāo)坐標(biāo)高程與測(cè)站3個(gè)坐標(biāo)分量以及誤差、平移參數(shù)均相關(guān)。

由式（1）、式（2）可以看出旋轉(zhuǎn)參數(shù)和系數(shù)矩陣呈非線性關(guān)系，無法直接使用整體最小二乘算法求解，一般將式（1）轉(zhuǎn)換成線性化的模型求解，但該方法僅僅適合于兩套坐標(biāo)系間的旋轉(zhuǎn)角較小時(shí)的轉(zhuǎn)換，若旋轉(zhuǎn)角較大時(shí)，模型誤差會(huì)很大，甚至完全得不到正確的結(jié)果［13，14］。

將式（1）和式（4）聯(lián)立，坐標(biāo)轉(zhuǎn)換模型轉(zhuǎn)化為非線性GM模型：

式中，l為式（1）和式（4）中目標(biāo)坐標(biāo)和測(cè)站坐標(biāo)的向量化；V為誤差的向量化。假設(shè)A為測(cè)站坐標(biāo)系統(tǒng)的n個(gè)坐標(biāo)，βc為6個(gè)坐標(biāo)轉(zhuǎn)換參數(shù)，vec表示逐行向量化，則變量參數(shù)β可表示為：

其中，f(β)為式（1）和式（4）聯(lián)立后的方程的右半部分。

經(jīng)過向量化變化后，式（5）轉(zhuǎn)換為非線性方程的最優(yōu)化問題，加權(quán)整體最小二乘的最優(yōu)估計(jì)的函數(shù)可表示為：

式中，P表示權(quán)值。

對(duì)于非線性方程的最優(yōu)化問題，常見的方法有梯度法和高斯牛頓法。梯度法主要從負(fù)梯度方向進(jìn)行迭代搜索，并最終找到極小值，該方法開始時(shí)下降較快，但靠近極小值時(shí)收斂速度很慢。高斯牛頓法中由雅克比（Jacobian）矩陣近似計(jì)算海森（Hessian）矩陣，避免了直接計(jì)算海森矩陣的復(fù)雜性，收斂速度快，但是該算法對(duì)初始值較敏感，容易陷入局部極小值。LM算法引入阻尼因子［15］，每次迭代時(shí)調(diào)節(jié)阻尼因子μ，使算法介于梯度法和高斯牛頓法之間，能有效克服常見缺點(diǎn)。對(duì)于本文的加權(quán)非線性方程的求解中，在梯度向量G和海森矩陣H更新時(shí)需要加入權(quán)值項(xiàng)。

每一次迭代時(shí)，h表示β值的更新值：

當(dāng)梯度的無窮范數(shù)小于預(yù)定的一個(gè)極小值ε1，或者通過更新值h和變量參數(shù)的2-范數(shù)進(jìn)行比較，當(dāng)小于變量參數(shù)的一個(gè)極小比例ε2時(shí)，程序找到非線性最小二乘的最優(yōu)解β，由式（6）可知，解的前6個(gè)值就是旋轉(zhuǎn)角和平移量參數(shù)。本文中ε1、ε2均設(shè)置為1×10-10。

1.2 權(quán)值的確定

機(jī)場(chǎng)道面很平整，橫坡和縱坡最大也不超過2%［16］，即平面偏差5 cm，對(duì)高程影響也就1 mm。為了保證高精度的道面高程測(cè)量成果，同時(shí)又提高效率。靶標(biāo)控制點(diǎn)的平面坐標(biāo)使用GPS-RTK采集，高程則使用二等水準(zhǔn)進(jìn)行測(cè)量。道面測(cè)量采用RIEGL VZ-1000三維激光掃描儀（角度分辨率1.8″，小于100 m距離的測(cè)距精度5 mm），是一款高密度、高精度、高效率的掃描設(shè)備。

目標(biāo)坐標(biāo)系統(tǒng)下靶標(biāo)點(diǎn)的測(cè)量誤差分為平面和高程兩個(gè)部分。由于平面坐標(biāo)由GPS-RTK采集，有固定誤差，對(duì)中桿手握對(duì)中以及其他因素引起的偶然誤差，要準(zhǔn)確計(jì)算每一個(gè)點(diǎn)的平面誤差比較困難。由于機(jī)場(chǎng)環(huán)境空曠，高度角大于15°的衛(wèi)星一般都能接收到信號(hào)，本文計(jì)算中統(tǒng)一設(shè)置為2.1 cm中誤差（兩個(gè)坐標(biāo)分量都是1.5 cm的中誤差）。由于電子水準(zhǔn)儀的普遍應(yīng)用，水準(zhǔn)測(cè)量的效率和精度都得到了很大提高，有的儀器甚至已經(jīng)達(dá)到了0.2 mm/km的高程中誤差，本案例的控制點(diǎn)采用二等水準(zhǔn)測(cè)量，平差后高程中誤差為0.6 mm，加上放置靶標(biāo)的誤差后，統(tǒng)一設(shè)置為2 mm。

測(cè)站坐標(biāo)系統(tǒng)是掃描儀外業(yè)測(cè)量時(shí)采用的坐標(biāo)系統(tǒng)，靶標(biāo)位置在該系統(tǒng)的誤差包括兩部分：一是儀器的測(cè)量誤差，另一個(gè)是靶標(biāo)識(shí)別的誤差。其中在測(cè)站坐標(biāo)系統(tǒng)下，如圖1所示，儀器的觀測(cè)量為：水平角α、垂直角θ、以及斜距S。則儀器的測(cè)量中誤差在3個(gè)方向上分別為［17］：

圖1 測(cè)站坐標(biāo)示意圖Fig.1 Schematic Diagram of Station Coordinate

靶標(biāo)識(shí)別誤差主要和靶標(biāo)上的點(diǎn)云密度有關(guān)，掃描儀道面掃描時(shí)使用的是6 cm×6 cm的平面靶標(biāo)，放置在10 m位置時(shí)約有6 000個(gè)點(diǎn)，點(diǎn)間距約為0.8 mm；放置在120 m位置時(shí)就只有200個(gè)點(diǎn)，點(diǎn)間距約為4.2 mm。雖然掃描儀沒有對(duì)中和整平裝置，但是一般情況下，掃描儀在測(cè)量時(shí)基本保持水平，本文中將靶標(biāo)識(shí)別誤差設(shè)定為點(diǎn)間距σxy≈σz2≈±d，其中x、y代表水平方向，即使用隨機(jī)軟件RiSCAN中的點(diǎn)云個(gè)數(shù)求出各靶標(biāo)的點(diǎn)間距后，將水平方向的誤差、垂直方向的誤差設(shè)定為點(diǎn)間距。

如圖2所示，x、y、z方向的識(shí)別誤差分別為：

圖2 靶標(biāo)誤差平面投影Fig.2 Plane Projection of Target Error

測(cè)站坐標(biāo)的3個(gè)方向總誤差分別為：

目標(biāo)坐標(biāo)和測(cè)站坐標(biāo)的誤差都計(jì)算出來后，就可以對(duì)整體最小二乘的權(quán)值進(jìn)行確定，pi=σ02/σi2，其中σ02可以選定為σi2的最大值。權(quán)矩陣可表示為：

式中，diag代表生成對(duì)角化矩陣的符號(hào)，即將向量的成員作為矩陣對(duì)角線上的元素。這樣，權(quán)值即可確定。

2 實(shí)例數(shù)據(jù)采集和分析

2.1 數(shù)據(jù)的獲取

機(jī)場(chǎng)道面周邊無連續(xù)制高點(diǎn)，數(shù)據(jù)采集時(shí)使用專用腳架提高儀器的架設(shè)高度。普通腳架只能把儀器高度架設(shè)到1.5～1.7 m，使用工業(yè)級(jí)的專用腳架，把儀器高度架設(shè)到2.6～3.0 m。提高儀器高度，不僅可以測(cè)量得更遠(yuǎn)，提高效率；同時(shí)，由于反射角的增大，反射率增加，也提高了測(cè)量精度。

如圖3所示，將靶標(biāo)貼在平整長(zhǎng)方體一面，直接將靶標(biāo)的投影中心對(duì)準(zhǔn)控制點(diǎn)平面中心，放置在控制點(diǎn)的最高處，靶標(biāo)面朝向儀器中心，優(yōu)點(diǎn)在于避免了過大的對(duì)中誤差和量高誤差。需要注意的是，匹配該靶標(biāo)點(diǎn)的高程應(yīng)該為控制點(diǎn)高程加上3 cm（平面靶標(biāo)的一半尺寸），平面坐標(biāo)不變。

圖3 架設(shè)儀器和放置靶標(biāo)Fig.3 Instrument Erection and Target Placement

靶標(biāo)布設(shè)時(shí)要分布均勻，且避免靶標(biāo)點(diǎn)共線和共面的情況，防止坐標(biāo)匹配時(shí)無窮多解的情況。待匹配的靶標(biāo)點(diǎn)布設(shè)在道面兩邊線外的助航燈基座螺母上，航燈的間距約為50 m，跑道寬度一般和機(jī)場(chǎng)等級(jí)有關(guān)，多數(shù)為45～60 m，如圖4所示。一個(gè)測(cè)站可以掃描12個(gè)靶標(biāo)，其中8個(gè)在本測(cè)站掃描范圍內(nèi)，另外4個(gè)在范圍外附近。機(jī)場(chǎng)道面測(cè)量時(shí)，將掃描儀架設(shè)在跑道中心線上，為便于相鄰測(cè)站的結(jié)果進(jìn)行檢核，使后一測(cè)站的點(diǎn)云和前一測(cè)站有重疊，使測(cè)站重疊區(qū)域離相鄰測(cè)站的最大距離在100 m左右。

圖4 靶標(biāo)點(diǎn)布設(shè)示意圖Fig.4 Schematic Diagram of Target Placement

2.2 分析和比較

實(shí)驗(yàn)1：采用8種坐標(biāo)匹配方案，比較各種方案的匹配精度，具體匹配方案和匹配精度如表1、表2所示。

表1 坐標(biāo)匹配方案Tab.1 Coordinate Matching Scheme

表2 坐標(biāo)匹配精度統(tǒng)計(jì)比較/mmTab.2 Statistical Comparison of Coordinate Matching Precision/mm

在表2中，方案1只有必要的3個(gè)點(diǎn)參與參數(shù)計(jì)算，沒有多余點(diǎn)，無法平差，雖然內(nèi)符合精度很高，但從檢核點(diǎn)的精度統(tǒng)計(jì)可以看出高程的差值范圍很大，標(biāo)準(zhǔn)差最大，證實(shí)了采用本文算法的必要性。

由于EIV模型把誤差作為變量進(jìn)行平差，考慮了目標(biāo)坐標(biāo)系統(tǒng)中存在誤差，所以在參與計(jì)算點(diǎn)的內(nèi)符合精度統(tǒng)計(jì)中，方案3和4的高程差值標(biāo)準(zhǔn)差比較大，但是檢核點(diǎn)的高程標(biāo)準(zhǔn)差最小。方案3的檢核點(diǎn)高程精度最高，平面精度符合實(shí)際RTK測(cè)量的情況，因此8個(gè)臨近的點(diǎn)的參數(shù)估計(jì)最穩(wěn)定。方案4的整體精度和方案3相近，也符合式（3）的εx、εy為小角度的情況。

方案5～方案8為匹配點(diǎn)在一側(cè)的情況，另一側(cè)的檢核點(diǎn)高程出現(xiàn)較大的偏差，這是由于靶標(biāo)分布不均勻造成的，這個(gè)結(jié)論也符合文獻(xiàn)［18］的結(jié)論。方案6、方案8使用了兩個(gè)近130 m的靶標(biāo)，靶標(biāo)的分辨率過低引起這兩個(gè)方案的高程整體精度不如方案5和方案7。

實(shí)驗(yàn)2：TLS多數(shù)沒有整平和對(duì)中系統(tǒng)，儀器的坐標(biāo)平面和目標(biāo)坐標(biāo)平面產(chǎn)生夾角在所難免，主要體現(xiàn)在以坐標(biāo)x軸和y軸發(fā)生旋轉(zhuǎn)，本文假設(shè)產(chǎn)生0～10°的旋轉(zhuǎn)，使用測(cè)站周邊的8個(gè)點(diǎn)作為參數(shù)計(jì)算點(diǎn)，剩下的4個(gè)點(diǎn)作為檢核點(diǎn)，對(duì)不加權(quán)和加權(quán)進(jìn)行比較。

如圖5所示，等權(quán)的情況下，參與計(jì)算點(diǎn)的內(nèi)符合高程轉(zhuǎn)換差值的標(biāo)準(zhǔn)差在1.3～1.8 mm之間，檢核點(diǎn)的高程轉(zhuǎn)換差值的標(biāo)準(zhǔn)差在3.5～5.5 mm之間。不論是內(nèi)符合還是檢核點(diǎn)的精度，都隨著旋轉(zhuǎn)角的加大而增加，最大的情況出現(xiàn)在x軸和y軸同時(shí)旋轉(zhuǎn)10°的情況?？梢钥闯?，在2°旋轉(zhuǎn)角的范圍內(nèi)，高程的轉(zhuǎn)換誤差還算穩(wěn)定，也就是旋轉(zhuǎn)角度較小的情況。

圖5 等權(quán)情況下的精度統(tǒng)計(jì)圖Fig.5 Precision Statistical Diagram Under Equal Weight

如圖6所示，加權(quán)的情況下，參與計(jì)算點(diǎn)的高程差值的標(biāo)準(zhǔn)差非常平穩(wěn)，保持在1.3 mm，檢核點(diǎn)的轉(zhuǎn)換差值標(biāo)準(zhǔn)差也較穩(wěn)定，保持在3.5～3.7 mm。對(duì)比圖5、圖6可以看出，符合式（3）的推導(dǎo)，即沒有整平裝置時(shí)，目標(biāo)坐標(biāo)高程和測(cè)站3個(gè)坐標(biāo)分量以及誤差、平移參數(shù)都相關(guān)的結(jié)論。加權(quán)的情況下保持著等權(quán)差值的最小值，并非常穩(wěn)定，說明了該算法具有通用性和強(qiáng)壯性。

圖6 加權(quán)情況下的精度統(tǒng)計(jì)圖Fig.6 Precision Statistical Diagram Under Weighted

由于兩站重疊部分的點(diǎn)云距離兩測(cè)站都最遠(yuǎn)，該處的掃描精度最低，外業(yè)時(shí)使用常規(guī)的方法（水準(zhǔn)儀）抽樣測(cè)量了幾個(gè)測(cè)站重疊部分的高程，共計(jì)119個(gè)點(diǎn)，與TLS點(diǎn)云過濾分類后的結(jié)果進(jìn)行比較得出，差值范圍在-12.1～8.7 mm之間，大于1 cm差值的點(diǎn)只有3個(gè)，且都位于道面邊緣處，是源于該處有覆土或者植被的影響，所有差值的均值為2.2 mm，差值標(biāo)準(zhǔn)差為±4.6 mm，也就是精度優(yōu)于1 cm，可以認(rèn)為本文的方法完全能滿足機(jī)場(chǎng)道面測(cè)量的要求。

3 結(jié)束語

本文分析了加權(quán)整體最小二乘算法的理論和實(shí)現(xiàn)方法，并且對(duì)TLS的目標(biāo)坐標(biāo)系統(tǒng)和測(cè)站坐標(biāo)系統(tǒng)誤差進(jìn)行了分析確定，并針對(duì)某機(jī)場(chǎng)實(shí)測(cè)的數(shù)據(jù)進(jìn)行計(jì)算比較，得出結(jié)論：由于測(cè)量效率的需要，匹配點(diǎn)的目標(biāo)坐標(biāo)采用不同的設(shè)備采集，這樣就需要使用加權(quán)整體最小二乘法進(jìn)行轉(zhuǎn)換參數(shù)的求解，這樣才能得到最優(yōu)的參數(shù)估計(jì)，并有效減小三維激光掃描儀不能整平帶來的誤差。同時(shí)在高精度的道面高程測(cè)量的應(yīng)用中，采用此方法可以有效的提高高程精度，最終達(dá)到優(yōu)于1 cm的高程精度要求。

傳統(tǒng)方式測(cè)量一個(gè)機(jī)場(chǎng)的板角高程需要4人（一名觀測(cè)員、一名記錄員、兩名跑尺員），15 d左右才能完成工作。采用本文方法只需要3人（兩名儀器操作員，一名靶標(biāo)擺設(shè)人員），3 d就可完成同樣工作。無論從人力需求上還是工作效率上都明顯得到了提高，而且避免了傳統(tǒng)方式平面和高程分開采集，最后融合易產(chǎn)生的錯(cuò)誤。本文方法對(duì)于其他一些對(duì)單一方向有高精度要求的測(cè)量工程，如道路、特殊構(gòu)件的平整度測(cè)量等，都有很好的借鑒作用。