王劍輝，尹 彤，趙有兵，張 瑩，李旭洋

（1.廣東省地質(zhì)測(cè)繪院，廣州 510800；2.自然資源部 測(cè)繪標(biāo)準(zhǔn)化研究所，西安 710054；3.中鐵二院工程集團(tuán)有限責(zé)任公司，成都 610031；4.自然資源部 第一地形測(cè)量隊(duì)，西安 710054）

0 引言

高速鐵路框架控制網(wǎng)（high-speed railway frame control network,CP0）應(yīng)在初測(cè)前采用全球衛(wèi)星導(dǎo)航系統(tǒng)（global navigation satellite system,GNSS）測(cè)量方法建立空間3 維（3D）控制網(wǎng)，全線一次性布網(wǎng)，統(tǒng)一測(cè)量，整體平差，作為全線（段）的平面坐標(biāo)框架基準(zhǔn)[1]。CP0 外業(yè)測(cè)量時(shí)，應(yīng)聯(lián)測(cè)國(guó)家A 級(jí)、B 級(jí)GNSS 控制網(wǎng)或國(guó)際GNSS 服務(wù)組織（International GNSS Service,IGS）連續(xù)運(yùn)行參考站等高等級(jí)已知點(diǎn)進(jìn)行整網(wǎng)基線解算和平差計(jì)算。

我國(guó)目前使用的坐標(biāo)系統(tǒng)為2000 國(guó)家大地坐標(biāo)系（China geodetic coordinate system 2000,CGCS2000），其坐標(biāo)框架為國(guó)際地球參考框架1997（international terrestrial reference frame,ITRF1997），歷元為2 000.0。其起始?xì)v元至今已跨度二十余年，由于亞歐板塊的絕對(duì)運(yùn)動(dòng)及內(nèi)部子板塊的相對(duì)運(yùn)動(dòng)，CGCS2000的很多框架點(diǎn)位置已經(jīng)發(fā)生較大變化。所以在CP0 與CGCS2000 框架點(diǎn)或 IGS 站聯(lián)測(cè)時(shí)，大跨度的高等級(jí)點(diǎn)的CGCS2000 坐標(biāo)跟其當(dāng)前歷元坐標(biāo)是否具有一致性，是確?；鶞?zhǔn)穩(wěn)定的關(guān)鍵。為了避免將基準(zhǔn)誤差帶入平差結(jié)果中去，必須提前對(duì)高等級(jí)已知點(diǎn)進(jìn)行一致性驗(yàn)證。為此，本文研究赫爾默特（Helmert）模型算法，并使用Helmert 法對(duì)高等級(jí)已知點(diǎn)穩(wěn)定性進(jìn)行一致性檢驗(yàn)，探測(cè)并剔除有問題的高等級(jí)已知點(diǎn)，以避免在CP0的解算中引入基準(zhǔn)誤差，從而得到更為精確的坐標(biāo)。

1 基準(zhǔn)一致性檢驗(yàn)方法

Helmert 法能夠計(jì)算輸入坐標(biāo)集到參考坐標(biāo)集的變換參數(shù)，從而體現(xiàn)出2 套坐標(biāo)集之間在平移、尺度、定向等方面的相關(guān)變化特性[2]，并通過轉(zhuǎn)換后的殘差探測(cè)和定位異常坐標(biāo)點(diǎn)。Helmert 模型主要是通過對(duì)2 套坐標(biāo)集進(jìn)行對(duì)比，通過2 套坐標(biāo)集公用點(diǎn)求取Helmert 模型7 參數(shù)，即3 個(gè)平移參數(shù)、3 個(gè)旋轉(zhuǎn)參數(shù)和1 個(gè)比例因子[3]，并按照7 參數(shù)將檢驗(yàn)坐標(biāo)系坐標(biāo)轉(zhuǎn)化到參考坐標(biāo)下。該模型建立的前提是承認(rèn)這2 套系統(tǒng)的質(zhì)量重心是一致的[4]。采用最小二乘約束平差方法，將2 套坐標(biāo)的差異導(dǎo)入Helmert 模型，然后采用參數(shù)平差模型，最終估計(jì)出Helmert 參數(shù)，同時(shí)求取坐標(biāo)殘差[5]。

對(duì)于比較大型的區(qū)域或全球網(wǎng)絡(luò)的坐標(biāo)，則應(yīng)采用空間直角坐標(biāo)（X，Y，Z）形式。在這種情況下，估計(jì)的Helmert 轉(zhuǎn)換參數(shù)是指地心笛卡爾坐標(biāo)系。該坐標(biāo)系下的轉(zhuǎn)換模型為

式中：Xlarge-0為大區(qū)域的參考坐標(biāo)矩陣；Xlarge為大區(qū)域的輸入坐標(biāo)矩陣；X0為平移向量矩陣；Ri(ri)為是繞i軸的旋轉(zhuǎn)角度ri的旋轉(zhuǎn)矩陣，i=x、y、z；k為縮放參數(shù)。

對(duì)于地方網(wǎng)絡(luò)和小區(qū)域網(wǎng)絡(luò)，適合采用測(cè)站空間直角坐標(biāo)（N，E，U）形式，則估計(jì)的變換參數(shù)將引用以輸入坐標(biāo)集重心（質(zhì)量中心）為原點(diǎn)的局部地平系統(tǒng)。輸入坐標(biāo)集的重心計(jì)算公式為

借助于重心X的橢球坐標(biāo)和，將這2 個(gè)坐標(biāo)集轉(zhuǎn)換為局部地平系統(tǒng)，即

式中：Xlocal為局部地平坐標(biāo)系坐標(biāo)矩陣；Ry和Rz為y軸或z軸的選擇矩陣；和為重心X的橢球坐標(biāo)；X為須轉(zhuǎn)換的坐標(biāo)向量。

按式（4）變換輸入坐標(biāo)集與參考坐標(biāo)集之間差的平方和，通過最小二乘法來估計(jì)方程中的Helmert 變換參數(shù)。

式中：Xlocal-0為小區(qū)域的參考坐標(biāo)集局部地平系坐標(biāo)；Xlocal-1為小區(qū)域輸入坐標(biāo)集局部地平系坐標(biāo)。

無論大區(qū)域還是地方小區(qū)域，在確定了3 個(gè)平移參數(shù)、3 個(gè)旋轉(zhuǎn)參數(shù)、1 個(gè)比例因子之后，可以根據(jù)式（5）將輸入坐標(biāo)轉(zhuǎn)換系至參考坐標(biāo)系，得到轉(zhuǎn)換后的參考坐標(biāo)系坐標(biāo)，即

式中：Xout為轉(zhuǎn)換后參考坐標(biāo)矩陣；Xin為輸入坐標(biāo)矩陣。

然后與原參考坐標(biāo)系坐標(biāo)進(jìn)行比較，得到轉(zhuǎn)換后的殘差，為：

式中：ΔXlarge為大區(qū)域轉(zhuǎn)換后殘差矩陣；ΔXlacal為小區(qū)域轉(zhuǎn)換后殘差矩陣。

計(jì)算轉(zhuǎn)換得到的參考系坐標(biāo)集與對(duì)應(yīng)的原參考坐標(biāo)集的差異值，即N（北）、E（東）、U（天頂）3 個(gè)方向的模型殘差分量[6]，獲得Helmert 轉(zhuǎn)換殘差，當(dāng)N、E、U 3 個(gè)方向的坐標(biāo)模型殘差大于該方向上綜合中誤差的3 倍的時(shí)候，即可判斷該坐標(biāo)不符合2 個(gè)坐標(biāo)集的整體一致性。

在高速鐵路框架控制網(wǎng)平差基準(zhǔn)一致性檢驗(yàn)中，可采用以上方法準(zhǔn)確判斷和定位有問題的高等級(jí)已知點(diǎn)。

2 實(shí)驗(yàn)與結(jié)果分析

本文采用2017 年年積日第214—222 天，某鐵路CP0的數(shù)據(jù)進(jìn)行實(shí)驗(yàn)分析。該數(shù)據(jù)分2 輪進(jìn)行外業(yè)觀測(cè)，每輪同步觀測(cè)5 個(gè)時(shí)段，白天和夜間觀測(cè)均勻分布，時(shí)間不低于5.5 h，輪間搭接2 個(gè)CP0 控制點(diǎn)，每個(gè)時(shí)段數(shù)據(jù)合格。求解14 個(gè)CP0控制點(diǎn)，基線解算和平差計(jì)算的起算點(diǎn)選用周邊10 個(gè) IGS 站，即北京房山站（BJFS）、長(zhǎng)春站（CHAN）、拉薩站（LHAZ）、波利根站（POL2）、上海天文臺(tái)站（SHAO）、韓國(guó)站（SUWN）、中國(guó)臺(tái)灣桃園站（TWTF）、蒙古站（ULAB）、烏魯木齊站（URUM）、武漢站（WUHN），在2017 年年積日第214—222 天期間內(nèi)，這些IGS 站的數(shù)據(jù)均比較完整，數(shù)據(jù)質(zhì)量良好，P1 多路徑影響mp1和P2 多路徑影響mp2 值均小于0.5 m。CP0 及IGS基準(zhǔn)站網(wǎng)型如圖1 所示。

圖1 CP0 及IGS基準(zhǔn)站分布網(wǎng)型

2.1 數(shù)據(jù)處理

高精度基線解算采用麻省理工學(xué)院（Massachusetts Institute of Technology,MIT）和美國(guó)斯克里普斯海洋研究所（Scripps Institution of Oceanography,SIO）聯(lián)合開發(fā)的加米特（GAMIT）/格洛布克（GLOBK）軟件10.70 版本[7]。GAMIT/GLOBK 軟件主要功能模塊包括軌道積分程序、觀測(cè)方程構(gòu)建程序、單差自動(dòng)修復(fù)周跳程序、雙差自動(dòng)修復(fù)周跳程序、交互式修復(fù)周跳程序、最小二乘求解程序等[8]。軟件基于多種改正模型，采用優(yōu)化的策略，可針對(duì)上千千米的超長(zhǎng)基線，實(shí)現(xiàn)毫米級(jí)的高精度基線解算[9]。

在基線解算中，固定站的N、E、U 3 個(gè)方向分別進(jìn)行松弛約束，約束分別為5、5、10 cm，分別得到年積日第214—222 天共9 個(gè)時(shí)段的10 個(gè)IGS 站的單時(shí)段松弛解，這10 個(gè)IGS 站為BJFS、CHAN、LHAZ、POL2、SHAO、SUWN、TWTF、ULAB、URUM、WUHN。

2.2 Helmert 轉(zhuǎn)換坐標(biāo)殘差指標(biāo)定位異常點(diǎn)

在得到單天松弛約束解后，分別提取BJFS、CHAN、LHAZ、POL2、SHAO、SUWN、TWTF、ULAB、URUM、WUHN 等10 個(gè)站的每個(gè)時(shí)段的松弛約束解的空間直角坐標(biāo)，作為第1 個(gè)坐標(biāo)集；收集該10 個(gè)站在國(guó)際ITRF2014 框架內(nèi)標(biāo)準(zhǔn)歷元坐標(biāo)作為參考數(shù)據(jù)集；對(duì)第1 坐標(biāo)集和參考坐標(biāo)集進(jìn)行Helmert 轉(zhuǎn)換；根據(jù)Helmert 轉(zhuǎn)換參數(shù)計(jì)算坐標(biāo)轉(zhuǎn)換殘差。根據(jù)Helmert的坐標(biāo)殘差指標(biāo)可以定位異常點(diǎn)。

本次采用Helmert 模型計(jì)算10 個(gè)IGS 站單個(gè)時(shí)段松弛解坐標(biāo)與 ITRF2014 標(biāo)準(zhǔn)歷元坐標(biāo)的Helmert 轉(zhuǎn)換坐標(biāo)殘差。殘差如表1 所示。

表1 各時(shí)段松弛解Helmert 轉(zhuǎn)換坐標(biāo)殘差 mm

由表可知：2017 年年積日第216 天SHAO 站在E 方向的殘差為43.4 mm，大于2 倍中誤差，U方向殘差為64.3 mm，大于3 倍中誤差；2017 年年積日第221 天，CHAN 站在U 方向的殘差為88.9 mm，大于3 倍中誤差。

由此可以判斷出2017 年年積日第216 天的SHAO 站、2017 年年積日第221 天的CHAN 站存在數(shù)據(jù)質(zhì)量問題，導(dǎo)致了解算的坐標(biāo)精度降低，需要將其剔除掉重新計(jì)算，才能保證基準(zhǔn)平差不引入基準(zhǔn)誤差。

2.3 單天解卡方檢驗(yàn)

卡方檢驗(yàn)是一種非連續(xù)性資料的假設(shè)檢驗(yàn)方法。基于卡方檢驗(yàn)可實(shí)現(xiàn)3 種檢驗(yàn)：①擬合度檢驗(yàn)，比較2 個(gè)及2 個(gè)以上樣本率；②相關(guān)性分析，2 個(gè)分類變量之間有無關(guān)聯(lián)性；③統(tǒng)一性檢驗(yàn)，檢驗(yàn)2 個(gè)或2 個(gè)以上總體的某一特性分布[10]。

在GLOBK 平差中，可以通過卡方檢驗(yàn)來宏觀了解平差基準(zhǔn)誤差的大小。平差后的卡方檢驗(yàn)值越小越好，如果卡方檢驗(yàn)值大于1，說明平差基準(zhǔn)存在誤差。本實(shí)驗(yàn)將通過卡方檢驗(yàn)方法來驗(yàn)證前期Helmert的分析結(jié)果。

實(shí)驗(yàn)對(duì)9 個(gè)時(shí)段分別進(jìn)行了單天平差計(jì)算，約束10 個(gè)IGS 站，采用GAMIT/GLOBK 軟件對(duì)每天的基線文件進(jìn)行平差計(jì)算，在獲取單天解的同時(shí)也可以得到單天解的卡方檢驗(yàn)值。該數(shù)據(jù)平差后的單天解卡方檢驗(yàn)值如表2 所示。

表2 單時(shí)段平差的卡方檢驗(yàn)

由表可以看出：2017 年年積日第216 天和第221 天這2 個(gè)時(shí)段的卡方檢驗(yàn)大于1，說明這2 個(gè)時(shí)段約束平差的固定基準(zhǔn)存在不一致的情況，但不能定位是哪些站出了問題。根據(jù)前面Helmert 殘差分析的情況，可以更加準(zhǔn)確地定位到出現(xiàn)問題的約束站，說明利用Helmert 殘差可以確定異?；鶞?zhǔn)點(diǎn)。

后面的平差過程需要針對(duì)有問題的時(shí)段，將有問題的站松弛開，約束其他沒有問題的基準(zhǔn)站重新進(jìn)行單天約束平差。如果卡方檢驗(yàn)通過，則須在單天數(shù)據(jù)中剔除掉有問題的基準(zhǔn)站，或者替換成其他基準(zhǔn)站，重新計(jì)算基線和平差。

3 結(jié)束語(yǔ)

本文通過與約束平差的卡方檢驗(yàn)指標(biāo)進(jìn)行對(duì)比，證明利用Helmert 模型計(jì)算殘差來發(fā)現(xiàn)和定位有問題的基準(zhǔn)站是可行的。在CP0 數(shù)據(jù)處理中，尤其是對(duì)大區(qū)域或者含IGS 站在內(nèi)的超大區(qū)域的CP0 數(shù)據(jù)進(jìn)行處理時(shí)，必須嚴(yán)格驗(yàn)證所選國(guó)家A級(jí)、B 級(jí)GNSS 控制網(wǎng)或IGS 站等高等級(jí)已知點(diǎn)的坐標(biāo)一致性。通過Helmert 模型能夠分析出高等級(jí)已知點(diǎn)異常坐標(biāo)，準(zhǔn)確定位異常數(shù)據(jù)，從而避免將誤差帶入控制網(wǎng)中，最終提高控制網(wǎng)的計(jì)算精度。