周明志，胡志剛，李鵬博

(武漢大學(xué)衛(wèi)星導(dǎo)航定位技術(shù)研究中心，湖北 武漢 430079)

?

BDS/GPS實時變形監(jiān)測技術(shù)的軟件實現(xiàn)

周明志，胡志剛，李鵬博

(武漢大學(xué)衛(wèi)星導(dǎo)航定位技術(shù)研究中心，湖北 武漢 430079)

Realization of BDS/GNSS Real-time Deformation Monitoring Software

ZHOU Mingzhi，HU Zhigang，LI Pengbo

摘要：闡述了BDS/GPS變形監(jiān)測實時算法，并利用C#和C語言研發(fā)了BDS/GPS實時變形監(jiān)測軟件。該軟件支持多種格式的數(shù)據(jù)采集、存儲、解算和圖形顯示。精密導(dǎo)軌試驗結(jié)果表明，5 min的BDS/GPS實時解在E、N、U 3個方向的重復(fù)性精度分別為1.0～1.3 mm、1.4～1.7 mm和2.7～3.5 mm，探測精度在平面和高程方向分別為0.2～0.4 mm和1.1 mm。

關(guān)鍵詞：BDS；實時監(jiān)測；Kalman濾波

以GPS為代表的衛(wèi)星導(dǎo)航定位技術(shù)由于具有高精度、全天候、準(zhǔn)實時的優(yōu)點，在橋梁、尾礦、大壩、高層建筑等方面得到了廣泛應(yīng)用[1-4]。1998年清江隔河巖大壩建立了我國第一個GPS變形監(jiān)測系統(tǒng)[5]，并在當(dāng)年的抗洪搶險決策中發(fā)揮了巨大作用。我國北斗衛(wèi)星導(dǎo)航系統(tǒng)于2012年12月27日完成區(qū)域階段部署，目前有5顆GEO、5顆IGSO和4顆MEO衛(wèi)星在軌運行，可為亞太大部分地區(qū)提供公開服務(wù)[6]。近年來圍繞北斗系統(tǒng)的高精度導(dǎo)航定位應(yīng)用研究正在逐步開展。已有研究表明，BDS/GPS聯(lián)合定位精度比GPS單系統(tǒng)的精度提高20%～30%，尤其在動態(tài)定位模式其精度提升明顯[7]。在GNSS變形監(jiān)測方面，現(xiàn)代化的超高層建筑、大型橋梁等不但要求毫米級的定位精度，還需要盡量縮短觀測周期來監(jiān)測其受大風(fēng)等外部策動力影響時建(構(gòu))筑物的擺動幅度和振動周期信息。這就對GNSS變形監(jiān)測的精度、穩(wěn)定性和實時性提出了較高的要求。

大壩、水庫、滑坡體變形具有連續(xù)、緩慢、量級小的特點，常規(guī)監(jiān)測中一般利用一段時間數(shù)據(jù)(數(shù)小時)作為一組，采用靜態(tài)相對定位方式處理，但此方法比較滯后，實時性不強。為了提高數(shù)據(jù)處理的實時性，常采用動態(tài)單歷元處理模式(OTF)[8]。OTF需要初始化幾十分鐘并且保證模糊度連續(xù)跟蹤。李博峰、韓保民等從提高法方程病態(tài)性角度求得單歷元穩(wěn)定解[9]。余學(xué)祥從觀測方程入手，提出直接求得變形量的方法[10]。上述兩種方法雖然避免了載波周跳探測，但都需要外界提供相對準(zhǔn)確的未知參數(shù)先驗信息。

本文首先探討了BDS/GPS多模數(shù)據(jù)融合應(yīng)用于高精度GNSS形變監(jiān)測的關(guān)鍵技術(shù)，包括BDS/GPS時空基準(zhǔn)統(tǒng)一、實時探測周跳、相位中心改正、Kalman濾波參數(shù)估計；然后運用C#和C語言開發(fā)了BDS/GPS實時變形監(jiān)測軟件；最后基于精密導(dǎo)軌試驗平臺采用人工增加變形量的方法驗證了算法的正確性和軟件的解算精度。

一、BDS/GPS多模數(shù)據(jù)融合關(guān)鍵算法

1. BDS/GPS時空基準(zhǔn)統(tǒng)一

由于廣播星歷擁有特定的坐標(biāo)基準(zhǔn)和時間基準(zhǔn)，進行不同系統(tǒng)組合定位時，往往需要考慮系統(tǒng)間的時間基準(zhǔn)和坐標(biāo)基準(zhǔn)差異。GPST和BDT都采用原子時作為時間基準(zhǔn)，但是其起算點不同：BDT起始?xì)v元分別為UTC時2006年1月1日0時0分0秒，而GPST為UTC時1980年1月6日0時0分0秒。由于UTC閏秒的存在，BDT與GPST之間有14 s的差異。除此之外，兩者之間還存在微小的同步誤差τ，在相對定位中可以忽略不計。兩者之間的轉(zhuǎn)換關(guān)系為

GPST=BDT+14 s+τ

在多模接收機采集的數(shù)據(jù)中，需要考慮上式以換算正確的衛(wèi)星發(fā)射時刻來計算正確的衛(wèi)星位置。

另外，BDS廣播星歷采用CGCS2000框架，GPS廣播星歷采用WGS-84框架，兩者之間存在一定差異。文獻[11]指出，在坐標(biāo)系的實現(xiàn)精度范圍內(nèi)，無須考慮坐標(biāo)系差異引起的結(jié)果偏差。實際上，由于大壩等構(gòu)筑物變形監(jiān)測局限于幾十千米范圍以內(nèi)，坐標(biāo)系之間的差異對基線解的影響不超過1 mm，故本文未考慮系統(tǒng)間的坐標(biāo)轉(zhuǎn)換。

時空基準(zhǔn)統(tǒng)一后的BDS/GPS雙差數(shù)學(xué)模型為

(1)

線性化后函數(shù)模型為

(2)

BDS與GPS是兩個獨立系統(tǒng)，觀測精度和空間分布不同導(dǎo)致對待估參數(shù)貢獻不同，在組合定位中需要進行非等權(quán)處理。如表1所示，Odolinski通過大量的數(shù)據(jù)統(tǒng)計出BDS和GPS非差偽距噪聲σp和非差載波噪聲σΦ經(jīng)驗值[13]。因此，BDS/GPS組合定位隨機模型為

(3)

表1　非差偽距載波噪聲

2. 接收機天線相位中心改正

接收機天線相位中心一般表述為相對于天線參考點的天線相位中心偏差(antenna phase center offser，PCO)和相對于PCO的天線相位中心變化(phase center variations，PCVs)。研究表明，不同的頻率其天線相位中心存在幾毫米甚至幾厘米的差異。不考慮相位中心改正項將會造成毫米級甚至厘米級的位置偏差[14-15]，特別是高程方向。然而，目前還沒有官方發(fā)布北斗接收機天線相位中心精確的改正模型。在這種情況下建議所有參與變形監(jiān)測的接收機天線都進行指北，從而在進行相對定位時可以消除天線相位中心變化(PCV)對基線解的影響，但PCO仍然無法消除。此時得到的基線是相對于北斗PCO和GPS PCO相應(yīng)的頻率組合PCO值而不是監(jiān)測站標(biāo)石。對于關(guān)心監(jiān)測站相對位置的變形監(jiān)測來說影響不大，但對于精密單點定位來說，如果要求計算測站標(biāo)位置時，接收機天線相位中心改正是強制修改項。

3. 實時周跳處理

周跳的產(chǎn)生中斷了高精度載波觀測值跟蹤弧段的連續(xù)性。1周的周跳將對載波產(chǎn)生約20 cm的偏差?？煽康闹芴綔y算法，尤其是小周跳探測在目前仍然具有一定的挑戰(zhàn)性。準(zhǔn)確可靠和能夠?qū)崟r周跳探測是高精度實時變形監(jiān)測解決模糊度的關(guān)鍵之一。目前學(xué)者提出了許多方法，如高次差法、小波分析法等，但都是基于事后數(shù)據(jù)，不能應(yīng)用于實時[16]。軟件采用劉志趙提出的TECR(total electron contents rate)方法[17]：采用電離層變化率并結(jié)合寬巷組合可實現(xiàn)實時周跳探測。寬巷Nwl表達式為

(4)

(5)

進一步推導(dǎo)k歷元TECR為

(6)

將電離層變化率推導(dǎo)k歷元電離層總量(TEC)變化值與載波觀測值推導(dǎo)出的電離層變化值作比較，可得

(7)

由式(5)—式(7)推導(dǎo)TECR探測周跳噪聲為

(8)

(9)

4.Kalman濾波參數(shù)估計

卡爾曼濾波算法是利用觀測向量來估計隨時間不斷變化的隨機向量，是一種對動態(tài)監(jiān)測數(shù)據(jù)進行處理的有效方法，具有提高數(shù)據(jù)處理效率、減少噪聲、提高預(yù)測精度等優(yōu)點，在變形監(jiān)測中常用來實時估計監(jiān)測站點位。其公式如下。

狀態(tài)更新

(10)

測量更新

(11)

(12)

初始方差陣P0為對角矩陣，其中坐標(biāo)參數(shù)方差由偽距單點定位求得，單差模糊度方差由式(12)通過誤差傳播定律確定。當(dāng)Pk陣設(shè)置不正確將會導(dǎo)致濾波發(fā)散，對于上一步探測出有周跳的衛(wèi)星，先驗Pk陣不能反映其模糊度精度，需要對其相應(yīng)的單差模糊度參數(shù)方差重新初始化。

(13)

再根據(jù)MLAMBDA[18]方法求雙差模糊度整數(shù)解N，根據(jù)模糊度殘差平和的次小值與最小值的比率ratio判斷模糊度固定是否正確[19]。檢驗通過后固定解為

(14)

值得注意的是，當(dāng)有新衛(wèi)星進入視野時，ratio值偏低，致使檢驗無法通過。這是因為新參與解算的模糊度精度較低，若和其他模糊度一起在整數(shù)最小二乘估計下(IIS)估計會影響其他模糊度固定。Verhagen提出部分模糊度固定方法[20]：從精度最高的模糊度開始，依次尋找子集模糊度，使它們理論固定成功率不低于給定的數(shù)值p0，通常是99%。只固定這部分模糊度，其余模糊度保持不變。這種方法需要反復(fù)計算模糊度成功率，直至找到最佳子模糊度組合。此處參考Verhagen提出的模糊度精度越高，整數(shù)解越可靠的思想，設(shè)定只有衛(wèi)星高度角大于25°其模糊才參與固定，25°以下模糊度只參與實數(shù)解。經(jīng)過15°到25°的濾波收斂，其精度不斷提高，理論成功固定率也提高，這時再參與LAMBDA解算[21]。

二、軟件實現(xiàn)

軟件采用客戶端/服務(wù)端(C/S)模式編寫，客戶端和服務(wù)端通過TCP/IP協(xié)議進行數(shù)據(jù)交換?？蛻舳素?fù)責(zé)結(jié)果圖表顯示和用戶命令交換，采用C#編寫。第1節(jié)所闡述的關(guān)鍵算法，以及測站管理、數(shù)據(jù)解碼等核心功能通過服務(wù)端完成，采用C語言編寫。服務(wù)端在后臺連續(xù)運行，只要輸入正確的賬號和密碼，客戶端就可以獲得服務(wù)端發(fā)送過來的結(jié)果數(shù)據(jù)，并進行圖表顯示，軟件界面如圖1所示。軟件流程如圖2所示，主要包括測站管理、數(shù)據(jù)解碼、BDS/GPS實時定位算法、結(jié)果顯示等功能。其中，測站管理是對測站進行添加刪除等操作，數(shù)據(jù)解碼是對實時二進制數(shù)據(jù)解碼。雖然RTCM 3.x修訂版開始支持BDS數(shù)據(jù)，message type為1121—1127。但現(xiàn)在市場上接收機大部分還是RTCM 2.0版本或不是最新版本的RTCM 3.0，并不包含BDS數(shù)據(jù)。監(jiān)測中采用BINEX格式，它是瑞士unavco機構(gòu)協(xié)調(diào)眾多接收機廠商制定的標(biāo)準(zhǔn)二進制數(shù)據(jù)格式，支持BDS/GNSS實時數(shù)據(jù)編碼，北斗觀測值和導(dǎo)航星歷Record ID分別是0x01-05和0x7f-05，可對解算的結(jié)果進行實時顯示，包括星空圖、信噪比、監(jiān)測站位移等。

圖1　軟件界面

圖2　軟件流程

三、精密導(dǎo)軌試驗驗證

為了驗證軟件的實時監(jiān)測精度，在國家光電中心樓頂對本文研發(fā)的變形監(jiān)測軟件進行了測試。測試樓頂有2個觀測墩，都安置TRM59800天線，其中一個監(jiān)測墩安置導(dǎo)軌，導(dǎo)軌每次移動大小可以通過千分尺準(zhǔn)確量取，如圖3所示，這樣可以通過移動導(dǎo)軌來模擬變形量。觀測時段模擬變形量見表2。在測試過程中，設(shè)Kalman濾波時間為5 min。當(dāng)時刻為整5 min時，重新初始化濾波器；最后統(tǒng)計5 min濾波結(jié)果精度。這樣既可以保證精度又可以保證探測靈敏度。

圖3　試驗環(huán)境

序號觀測時段狀態(tài)時段12015-04-0916:40:57—2015-04-1015:25:00 初始狀態(tài)時段22015-04-1015:30:00—2015-04-1116:00:20 E方向移動10mm時段32015-04-1116:40:10—2015-04-1217:25:00 W方向移動10mmU方向挪動30mm時段42015-04-1217:25:25—2015-04-1317:20:25 E方向移動10mm時段52015-04-1317:20:25—2015-04-1415:20:20 W方向移動5mm

系統(tǒng)連接接收機IP后，自動對流動站和監(jiān)測站傳過來的BINEX格式BDS/GPS實時數(shù)據(jù)解碼、同步、組基線、Kalman濾波、LAMBDA解算、圖形顯示并將結(jié)果保存到數(shù)據(jù)庫中。從數(shù)據(jù)庫中提取濾波5 min結(jié)果，其時間序列圖和統(tǒng)計表分別見圖4和表3。

圖4　模擬形變E/N/U三維方向時間序列圖

時段平均值/m均方根RMS/mmENUENU解個數(shù)1-3.925-0.479-0.2061.01.52.82752-3.915-0.479-0.2061.11.52.92863-3.925-0.479-0.2371.21.73.52964-3.915-0.478-0.2381.31.42.72885-3.920-0.479-0.2381.01.52.7233

注：每個時段起始截止觀測時間不一致，網(wǎng)絡(luò)傳輸可能存在丟包情況，因而每個時段解的個數(shù)不一樣。

表4　各時段模擬形變探測精度　mm

從圖4可以看出，當(dāng)沒有發(fā)生變形時，實時解穩(wěn)定，當(dāng)導(dǎo)軌移動，實時解有一個突變，可以很明顯看出一個形變過程。從表3中可以看出，BDS/GPS實時解3個方向解的穩(wěn)定性分別在1.0～1.3 mm、1.4～1.7 mm、2.7～3.5 mm之間，E方向解的穩(wěn)定性高于N方向，這是BDS衛(wèi)星幾何結(jié)構(gòu)決定的。從表4可以得到軟件實時探測精度在平面和高程方向分別是0.2～0.4 mm、1.1 mm。

四、結(jié)束語

高精度實時變形監(jiān)測軟件對實時性要求很高，傳統(tǒng)的最小二乘往往是基于數(shù)小時的數(shù)據(jù)，結(jié)果滯后嚴(yán)重，不能應(yīng)用于實時。單歷元解不穩(wěn)定，需要額外提供相對準(zhǔn)確的先驗信息。Kalman濾波能夠基于狀態(tài)域和觀測域進行狀態(tài)參數(shù)實時更新，非常適合應(yīng)用于實時變形監(jiān)測。本文介紹了基于Kalman濾波的BDS/GPS定位的關(guān)鍵算法，并在實現(xiàn)中用單差模糊度替代雙差模糊度，用部分模糊度替代全部模糊度進行固定?；谠撍惴?，采用C和C#語言，筆者開發(fā)了BDS/GPS實時變形監(jiān)測軟件。通過精密導(dǎo)軌試驗證明該軟件BDS/GPS實時解穩(wěn)定性為1.0～1.3 mm、1.4～1.7 mm和2.7～3.5 mm之間，探測精度在平面和高程方向分別為0.2～0.4 mm和1.1 mm，滿足變形監(jiān)測的需求。在高精度形變監(jiān)測中，多路徑是影響定位結(jié)果一個重要因素，在后續(xù)工作，需對其影響作進一步研究。

參考文獻：

[1]劉超. 礦區(qū) GPS 變形監(jiān)測及其偽衛(wèi)星增強技術(shù) [D]. 徐州：中國礦業(yè)大學(xué), 2011.

[2]李偉. GNSS 算法研究及其在變形監(jiān)測中的應(yīng)用 [D]. 成都：西南交通大學(xué), 2013.

[3]任常. 超高層建筑 GPS 變形監(jiān)測單歷元算法研究[D]. 北京：清華大學(xué), 2013.

[4]丁培江. 水庫潛在滑坡區(qū) GNSS 變形監(jiān)測管理信息系統(tǒng)研究與開發(fā)[D]. 昆明：昆明理工大學(xué), 2014.

[5]姜衛(wèi)平, 劉經(jīng)南. GPS技術(shù)在隔河巖大壩監(jiān)測中的應(yīng)用研究[J]. 武漢測繪科技大學(xué)學(xué)報, 1998,23(S0):20-22.

[6]Chinese Satellite Navigation System Management Office. Beidou Satellite Navigation System Development Report(2.1)[R]. Beijing：[s.n.], 2012.

[7]SHI C, ZHAO Q, HU Z, et al. Precise Relative Positioning Using Real Tracking Data from COMPASS GEO and IGSO Satellites[J]. GPS solutions, 2013, 17(1): 103-119.

[8]陳永奇,JAMES L. 單歷元GPS變形監(jiān)測數(shù)據(jù)處理方法的研究[J]. 武漢測繪科技大學(xué)學(xué)報, 1998，23(4): 45-49,84.

[9]李博峰, 沈云中. 附有約束條件的GPS模糊度快速解算[J]. 武漢大學(xué)學(xué)報(信息科學(xué)版), 2009,34(1): 117-121.

[10]余學(xué)祥, 徐紹銓, 呂偉才. GPS變形監(jiān)測信息的單歷元解算方法研究[J]. 測繪學(xué)報, 2002,31(2):123-128.

[11]魏子卿. 2000中國大地坐標(biāo)系及其與WGS84的比較[J].大地測量與地球動力學(xué), 2008, 28(5):1-5.

[12]盧祥. GPS/BDS短基線實時快速相對定位研究[D]. 北京：中國科學(xué)院大學(xué)，2014.

[13]ODOLINSKI R, TEUNISSEN P J G, ODIJK D. An Analysis of Combined COMPASS/BeiDou-2 and GPS Single-and Multiple-frequency RTK Positioning[C]∥Proceedings of the ION 2013 Pacific PNT Meeting. Honolulu, Hawaii：[s.n.], 2013: 69-90.

[14]ROTHACHER M, SCHAER S, MERVART L, et al. Determination of Antenna Phase Center Variations Using GPS Data[C]∥Proc.. IGS Workshop on Special Topics and New Directions. [S.l.]: IGS,1995: 205-220.

[15]MADER G L. GPS Antenna Calibration at the National Geodetic Survey[J]. GPS Solutions, 1999, 3(1): 50-58.

[16]方榮新, 施闖, 魏娜,等. GPS數(shù)據(jù)質(zhì)量控制中實時周跳探測研究[J]. 武漢大學(xué)學(xué)報(信息科學(xué)版), 2009,34(9): 1094-1097.

[17]LIU Z. A New Automated Cycle Slip Detection and Repair Method for a Single Dual-frequency GPS Receiver[J]. Journal of Geodesy, 2011, 85(3): 171-183.

[18]CHANG X W, YANG X, ZHOU T. MLAMBDA: a Modified LAMBDA Method for Integer Least-squares Estimation[J]. Journal of Geodesy, 2005, 79(9): 552-565.

[19]TEUNISSEN P J G, VERHAGEN S. The GNSS Ambiguity Ratio-test Revisited: A Better Way of Using It[J]. Survey Review, 2009, 41(312): 138-151.

[20]VERHAGEN A A, TEUNISSEN P J G, VAN DER MAREL H, et al. GNSS Ambiguity Resolution: Which Subset to Fix?[C]∥IGNSS Symposium 2011, International Global Navigation Satellite Systems Society. Sydney, Australia: University of New South Wales, 2011.

[21]TAKASU T, YASUDA A. Kalman-filter-based Integer Ambiguity Resolution Strategy for Long-baseline RTK with Ionosphere and Troposphere Estimation[C]∥ION NTM. [S.l.]: ION, 2010: 161-171.

中圖分類號：P228.4

文獻標(biāo)識碼：B

文章編號：0494-0911(2016)02-0026-06

作者簡介：周明志(1990—), 男，碩士，主要從事形變監(jiān)測研究。E-mail: nanwang_luojia@foxmail.com

基金項目：國家863計劃(2013AA122501)；國家自然科學(xué)基金(2014AA123101)；中央高校基本科研業(yè)務(wù)費專項資金(2042014kf0085)

收稿日期：2015-06-30

引文格式： 周明志，胡志剛，李鵬博. BDS/GPS實時變形監(jiān)測技術(shù)的軟件實現(xiàn)[J].測繪通報，2016(2)：26-31.DOI:10.13474/j.cnki.11-2246.2016.0042.