改善移動終端測量精度的卡爾曼濾波算法

呂偉才，蔣法文，杭玉付，朱亞洲，楊　旭，陳小軼，汪　洋，汪　澤

(1.安徽理工大學　測繪學院，安徽 淮南　232001；2.淮南礦業(集團)有限責任公司，安徽 淮南　232001)

改善移動終端測量精度的卡爾曼濾波算法

呂偉才1，蔣法文2，杭玉付2，朱亞洲1，楊旭1，陳小軼1，汪洋1，汪澤1

(1.安徽理工大學測繪學院，安徽 淮南232001；2.淮南礦業(集團)有限責任公司，安徽 淮南232001)

摘要：針對開采沉陷監測中常規RTK測量精度偏低的問題，結合開采沉陷監測的特點，以煤礦開采沉陷自動化監測系統的實時數據采集終端系統2 min內各歷元采集的坐標及精度為基礎，構建卡爾曼濾波模型，以進一步改善移動終端測量精度。利用地表移動觀測站的實測數據，通過與常規RTK測量、水準測量結果的比較分析，結果表明濾波后的測量精度得到了較好的改善。從內部符合精度來看，濾波RTK測量的平面位置平均測量精度約為±0.3 cm，高程方向平均測量精度約為±0.5 cm，完全滿足開采沉陷監測的精度要求；與水準測量成果相比，高程方向平均外部測量精度約為±0.8 cm，基本滿足開采沉陷監測的精度要求。

關鍵詞：開采沉陷監測；移動終端；RTK測量；卡爾曼濾波；精度分析

0引言

《煤礦測量規程》要求，礦區的第一個回采工作面必須建立開采沉陷監測地表移動觀測站(以下簡稱地表移動觀測站)，通過重復性的數據采集、處理和分析來掌握開采引起的地表移動變形特征與規律，為滿足煤礦安全生產、恢復與重建礦區生態環境及研究重復采動時的地表移動變形規律提供基礎資料。在建立地表移動觀測站的過程中，涉及到連接測量、全面觀測、巡視測量、日常觀測等階段的數據采集工作[1]，特別是全面觀測和日常觀測階段的數據采集工作更為頻繁。在監測點較多的情況下(如，顧橋煤礦11171首采面地表移動觀測站含18個控制點和174個監測點)[2-4]，采用常規方法就需要組織大量的人力物力來實施，建立觀測站的成本居高不下。

對于礦山開采沉陷監測，要求相鄰兩期間平面點位相對中誤差≤±2 cm，最弱點高程中誤差≤±1 cm，以保證解算的開采沉陷的關鍵參數(如邊界角、移動角等)的精度。這就要求一次測量的平面點位相對中誤差≤±1.4 cm，最弱點高程中誤差≤±0.7 cm。雖然從外業數據采集速度上來說，目前有多種技術，如連續運行參考站(continuously operating reference stations，CORS)系統、測量機器人、3維激光掃描等可以達到這一要求，但對于開采沉陷監測而言，這些手段還難以推廣應用。如與常規實時動態(real-time kinematic，RTK)定位相比，CORS RTK技術具有數據采集速度快、精度高的特點，在20 km范圍內公認的平面位置精度約±2 cm，高程精度約±3 cm[5-6]。這種測量精度(特別是高程方向的精度)還難以滿足開采沉陷監測的要求，同時對目前的礦區來說，還沒有全面建立網絡CORS系統的條件。

為解決當前開采沉陷數據采集技術手段落后、效率低、信息化程度低等問題，淮南礦業(集團)有限責任公司與安徽理工大學、南京科博空間信息科技有限公司聯合開展“地表移動自動化監測系統研究”項目，其研究內容之一是結合煤礦開采沉陷監測的特點，研究集成CORS和移動地理信息系統(geographic information system，GIS)的高精度快速數據采集終端系統。地表移動自動化監測系統共布設1個基準站(PYDCDP)、9個連續實時監測站(CORS1～CORS9)和約60個非連續實時監測站(地表移動觀測站中的監測點)[7]。其中:PYDCDP布設于礦區辦公樓頂，提供監測基準及差分數據；CORS1～CORS9分別布設于觀測線關鍵部位，采用連續實時監測、定時實時監測模式；非連續實時監測站采用CORS RTK監測模式，也可采用全站儀、數字水準儀的測量模式。

在“地表移動自動化監測系統研究”項目中，為提高移動終端的CORS RTK測量的精度，特別是高程方向的測量精度，采用了事后卡爾曼(Kalman)濾波算法。目前Kalman濾波被廣泛地應用于各種動態測量系統中，特別是在變形監測[8-9]、動態數據處理[10-11]、全球定位系統(global positioning system，GPS)定位定軌等方面的應用更為多見。基于Kalman濾波的RTK 測量算法討論的較多，但一般是以虛擬參考站的虛擬觀測值和移動站的實測觀測值為基礎，通過快速解算模糊度來實現的[12-14]。本文在介紹對常規CORS RTK測量獲得的測點位置序列進行濾波的基本原理基礎上，通過與常規CORS RTK測量、水準測量結果的比較分析，獲得濾波CORS RTK測量的精度信息，為采用Kalman濾波算法的移動終端進行開采沉陷監測提供參考依據。

1Kalman濾波方程

采用地表移動自動化監測系統的移動終端進行CORS RTK測量時，各測點每一歷元的主要測量信息包括點名、參考坐標系下的平面坐標(x，y)、擬合高程h、WGS-84大地坐標(B，L，H)、平面坐標中誤差和擬合高程中誤差。對某一測點進行了多次測量(如30次)時，Kalman濾波算法的主要思路為：首先將各歷元在WGS-84坐標系下的大地坐標轉換為高斯平面坐標，并利用Kalman濾波模型進行濾波，獲得濾波后的高斯平面坐標和大地高；然后利用連接測量建立的坐標系統轉換模型和高程系統轉換模型，獲得各歷元在參考坐標系下(如1954北京坐標系)下的坐標和擬合高程(如1985國家高程基準下的高程)；最后取各歷元轉換后的平面坐標和高程的平均值并進行精度評定。

1.1Kalman濾波的狀態方程

采用CORS RTK測量獲得的測點3維坐標(WGS-84坐標系下高斯平面直角坐標和大地高)序列中，將t時刻測點i的3維坐標ξi(t)和3維速率λi(t)作為狀態向量，將其瞬時加速率Ωi(t)看作隨機干擾，則有以下微分關系式[14]

(1)

記i點的狀態向量為Xi(t)(在本文中，用X、Y、Z表示測點3維坐標)，即：

(2)

(3)

則式(1)可寫成

(4)

式中0和E分別為3階零矩陣和3階單位陣。對式(4)進行拉普拉斯變換并將其離散化，可得Kalman濾波的狀態方程為

(5)

式中：Δtk=tk+1-tk；tk和tk+1分別為第k歷元和第k+1歷元的觀測時刻。將式(5)記為

(6)

式(6)即為采用CORS RTK模式進行動態監測時，動態系統的Kalman濾波狀態方程。

1.2Kalman濾波的觀測方程

利用CORS RTK測量監測點的3維位置時，采用的是虛擬參考站(virtual reference stations，VRS)技術。移動站開機啟動后，先向控制中心發送一個概略位置信息，控制中心接收到該位置信息后，利用基準站(PYDCDP)的觀測數據在移動站附近模擬出一個虛擬參考站，整體地改正全球衛星導航系統(global navigation satellite system，GNSS)的軌道誤差以及大氣折射引起的誤差;在此基礎上摸擬出高精度的虛擬觀測信息并發給移動站。虛擬參考站與移動站構成短基線，從而高精度、快速地解算出移動站的3維位置及相關精度信息[15]。

本文所討論的Kalman濾波的觀測方程，是以VRS技術解算的流動站的3維坐標為觀測值，而不是以虛擬參考站的虛擬觀測值和移動站的實測觀測值為依據。這是因為:一方面，在進行CORS RTK測量時，一般不保留移動站的實測觀測值；另一方面，以3維坐標為觀測值時有利于Kalman濾波的快速實施。

以CORS RTK測量獲得的監測點3維位置為觀測值時，某一監測點i在第k+1個歷元的觀測方程為

Li/k+1=ξi/k+1+Δti/k+1λj/k+1+Δi/k+1。

(7)

式中：Δti/k+1=ti/k+1-tk+1;ti/k+1為監測點i在第k+1歷元的觀測時刻;tk+1為監測點i本次測量時所有歷元(如30個歷元)觀測的中心時刻。

采用CORS RTK技術進行開采沉陷地表移動變形監測的目的之一，是在保證測量精度的前提下盡量縮短作業時間，如一個測點的測量時間不超過2 min。此時Δti/k+1不超過1 min，對于開采沉陷監測而言，在如此短的時間內，地表下沉速度是可以忽略不計的。

狀態方程(6)和觀測方程(7)共同構成了監測點3維位置序列的Kalman濾波模型

(8)

式中：Φk，k-1為k-1到k歷元的轉移矩陣；Γk,k-1為系統噪聲矩陣；Ωk-1為k-1歷元的系統噪聲；Bk為k歷元系統的觀測矩陣；Δk為k歷元系統的觀測噪聲；Xk為k歷元的系統待估狀態參數；Lk為k歷元系統的觀測向量矩陣；Xk和Lk均為監測點的3維位置和速度向量。

2位置序列的Kalman濾波算法

2.1Kalman濾波遞推公式

與式(8)相應的Kalman濾波遞推公式為

(9)

式中I為單位矩陣，并有

(10)

式中：X(k/k-1)為一步預測值；DX(k/k-1)為一步預測方差陣；Jk為狀態增益矩陣；Ek為預測殘差。

式(9)和式(10)即為監測點空間位置序列的Kalman濾波遞推計算公式。在按一定方法確定了系統的初始狀態后，利用新的觀測值(下1個歷元測點的空間位置)即可求得新的狀態參數濾波值。

2.2濾波初值的確定

從Kalman濾波方程可以看出，要確定動態系統在tk時刻的狀態，首先必須知道系統的初始狀態，即應賦給系統初值。

1)狀態參數初值的確定。

對于CORS RTK測量獲得測點空間位置的動態系統而言，狀態參數為測點的3維坐標和速度。由于此處的Kalman濾波是采用事后處理模式，因此：3維坐標的初值可取為測點各歷元3維位置的平均值；3維速度的初值可取為相鄰歷元間3維位置變化的平均值。若某一測點觀測了n個歷元，則有

(11)

式中以X分量為例，有：

即狀態參數初值為

(12)

2)方差陣初值的確定。

某一歷元觀測噪聲的方差陣DΔ， 即CORSRTK測得的測點該歷元3維位置的方差陣，由該歷元3維位置的中誤差mX、mY、mZ確定，即

(13)

狀態參數中的3維位置的初始方差陣，可由CORS RTK測得的測點各歷元3維位置的中誤差平均值確定，即

(14)

式中：

由協方差陣傳播定律，狀態參數中的3維速率的初始方差陣為

Dλ0=2Δt-2Dξ0。

(15)

式中△t為相鄰歷元間的采用間隔，則狀態參數的初始協方差陣為

(16)

3)動態噪聲的初始方差陣。

當以測點的3維位置和瞬時速率作為狀態參數時，則瞬時加速率為動態噪聲，根據協方差陣傳播定律，其方差陣為

DΩ(0)=4Δt-4Dξ0。

(17)

這樣，利用CORSRTK測量獲得測點3維位置序列即可獲得系統的初始狀態；然后利用式(9)和式(10)對此3維位置序列進行濾波。

3試驗與結果分析

3.1試驗簡介

為驗證本文提出的濾波算法的正確性，于2014-11-17在研究區的潘一東區1252(1)首采工作面觀測站和12421(1)綜采面觀測站的部分監測點(共45點)上進行試驗。試驗時，分別采用1臺LeciaSprinter150數字水準儀和1臺套CORSRTK移動終端(采用兼容BDS/GPS的接收機)進行外業數據采集。采用CORSRTK測量時，每點觀測30次以上。

對于RTK測量，坐標系統轉換和高程系統轉換參數采用連接測量階段的成果。觀測時，要求對中桿對點準確、立桿穩定，兩人配合，一人對點，一人采集數據；每點觀測2min，采樣率為2s，當觀測條件較差時(如樹底下)，適當延長時間。

對于RTK測量，分別采用常規數據處理方法和本文提出的濾波處理方法。采用常規數據處理方法時，取各測點所有觀測歷元的平均值并進行精度評價；采用濾波方法，數據處理方法參見第一節所述“Kalman濾波算法的主要思路”。

3.2數值結果與分析

(18)

對于濾波RTK測量，若某一測點i觀測了n個歷元，取各歷元濾波后平面坐標(分量)和高程的平均值為該測點的最終成果，即:

則平均值中誤差為:

有

(19)

從表1和表2中可以看出：對于RTK測量本身的精度而言，與常規RTK測量相比，采用本文提出的Kalman濾波RTK測量，其平面位置精度和高程位置精度得到了大幅度的提升(約提高了4倍)；同時精度分布更為均勻，即濾波后各點的測量成果更為穩定和可靠。這就為獲得穩定的正常高提供了基礎。

表3中列出了這45個測點上濾波RTK高程測量、常規RTK高程測量與水準測量高程測量之間結果差異的統計信息，圖1繪出了該差異圖形。

表1　RTK測量平面精度統計信息

表2　RTK測量高程精度統計信息

表3　高程測量差值統計信息

圖1　RTK高程測量與水準測量的差異

從表3和圖1可以看出，對于本次試驗而言：

1)對于常規RTK測量，其轉換后的正常高與水準測量成果相比，差異很大，且這種差異存在系統性,最大相差約10 cm，平均差值約8 cm，差值中誤差約1 cm，難以滿足開采沉陷監測對高程測量的精度要求。

2)對于濾波RTK測量，與常規RTK測量相比，除其內部精度得到大幅度提升外(參見表1和表2)，其轉換后的正常高的精度也得到了大幅度的提高，且消除了常規RTK高程測量的系統性誤差。其轉換后的正常高與水準測量成果相比，差異較小，最大相差約1.8 cm，平均差值約0.3 cm，差值中誤差約0.9 cm。從中誤差來看，略超過要求的±0.7 cm，基本滿足開采沉陷監測對高程測量的精度要求。

3)對于濾波RTK測量，與水準測量高程相比，差值在±0.7 cm以內的點有21個，占46.7%，差值在±1.0 cm以內的約占70%以上，差值在±1.5 cm以內的點有44個，占97.7%，極少出現差值超過±1.5 cm的情況。這說明濾波RTK高程測量，具有較好的穩定性和可靠性。

4結束語

在“地表移動自動化監測系統研究”項目中，為提高實時數據采集終端系統采用CORS RTK測量模式獲得的各測點空間位置(特別是高程)信息的精度和可靠性，結合開采沉陷監測的特點，本文建立了Kalman濾波算法，對RTK測量獲得的測點空間位置序列進行進一步處理。研究區域的實測數據表明:采用Kalman濾波RTK測量，其平面位置精度和高程位置精度得到了大幅度的提升(約為常規RTK測量的4倍);同時精度分布更為均勻，即濾波后各點的測量成果更為穩定和可靠。這就為獲得穩定的正常高提供了基礎。

對于開采沉陷監測，為確定移動變形邊界，要求平面點位相對中誤差不超過±1.4 cm，邊界部分最弱點高程中誤差不超過±0.71 cm。對于本試驗的濾波RTK測量而言，平面位置精度平均不超過±2.5 mm，最大不超過±5.0 mm，完全滿足開采沉陷監測對平面位置測量的精度要求；對于濾波RTK高程測量，與水準測量高程相比，高程差值的中誤差為±0.84 cm，與±0.71 cm相接近。若為確定正確的移動變形邊界，該種算法還需要進一步改進；但對于沉陷區的非邊界區域測量，濾波RTK高程測量仍可以獲得良好的監測成果，這就為高效、快速采集移動變形信息提供了技術保障。

影響濾波RTK高程測量精度的因素主要有RTK測量本身的精度(平面和高程)、天線相位中心偏差、測區高程(似大地水準面)擬合模型的精度、觀測條件等，特別是前2種因素，這是后續改進中需要考慮的主要問題。

參考文獻

[1]何國清，楊倫，凌賡娣，等.礦山開采沉陷學[M].徐州：中國礦業大學出版社，1991：53-56.

[2]余學祥，秦永洋，孫興平，等.相鄰工作面綜合地表移動觀測站的設計與連接測量[J].大連大學學報，2008，29(6)：74-79.

[3]余學祥，秦永洋，孫興平，等.顧橋煤礦11-2煤綜采面地表移動變形基本特征分析[J].礦山測量，2009(12)：8-12.

[4]余學祥，鄧蓉蓉，張美微，等.基于似單差法的井筒沉降監測試驗與結果分析[J].合肥工業大學學報：自然科學版，2012，35(6)：804-808.

[5]馬洪濱，孫軍，周海壯，等.鞍山市連續運行參考站系統建設與定位精度分析[J].導航定位學報，2013，1(1)：89-95.

[6]吳文壇，田摯，李哲，等.河北CORS在大比例尺傾斜航空攝影測量中的應用[J].測繪通報，2013(2)：63-66.

[7]蔣法文，劉可勝，杭玉付，等.煤礦開采沉陷自動化監測實時數據采集終端系統研究[J].礦山測量，2014(5)：93-97.

[8]余學祥，呂偉才.GPS監測網動態數據處理抗差Kalman濾波模型[J].中國礦業大學學報，2000，29(6)：553-557.

[9]余學祥，呂偉才.抗差卡爾曼濾波模型及其在GPS監測網中的應用[J].測繪學報，2001，30(1)：27-31.

[10]楊元喜.動態Kalman濾波模型誤差的影響[J].測繪科學，2006，31(1)：17-18.

[11]尹航，朱紀洪，周池軍，等.基于Kalman預報觀測器的增量動態逆控制[J].清華大學學報：自然科學版，2014，54(12)：1534-1538，1543.

[12]柯福陽，王慶，潘樹國.網絡RTK長基線模糊度解算方法研究[J].大地測量與地球動力學，2012，32(5)：72-77.

[13]徐彥田，程鵬飛，蔡艷輝，等.單頻RTK動態解算的卡爾曼濾波算法研究[J].測繪科學，2012，37(4)：43-44.

[14]LV Weicai，XU Shaoquan.Kalman filtering algorithm research for the deformation information series of the similar single-difference model[J].Journal of China University of Mining and Technology，2004，14(2)：189-194.

[15]余學祥，呂偉才，柯福陽，等.煤礦開采沉陷自動化監測系統[M].北京：測繪出版社，2014：192-217.

Kalman filter algorithm for improving mobile terminal measurement accuracy

LV Weicai1，JIANG Fawen2，HANG Yufu2，ZHU Yazhou1，YANG Xu1，CHEN Xiaoyi1，WANG Yang1，WANG Ze1

(1.School of Geomatics，Anhui University of Science and Technology，Anhui Huainan 232001，China；2.Huainan Mining Group Co.，Ltd.，Anhui Huainan 232001，China)

Abstract：For conventional RTK measurement accuracy is lower in mining subsidence monitoring，according to the characteristics of the coal mining subsidence monitoring，and basing on the coordinate and it’s accuracy collected by the real-time data acquisition terminal system of the coal mining subsidence monitoring automation system，in order to further improve the mobile terminal measurement accuracy，the kalman filtering model was built.Using the measured data of surface movement observation station，and comparing with the results of the conventional RTK measurement and levelling measurement，the filtered measuring accuracy got better improvement.From the internal accuracy of filtering RTK measurement，the plane measurement average accuracy could be up to±0.3 cm，and the height measurement average accuracy was about ±0.5 cm.This internal accuracy could fully meet the accuracy requirement of mining subsidence monitoring.Comparing with the results of the filtering RTK measurement and levelling measurement，the height measurement external accuracy was about ±0.8 cm，which basically meet the accuracy requirement of mining subsidence height monitoring.

Keywords：mining subsidence monitoring；terminal system；RTK surveying；Kalman filter；accuracy analysis

收稿日期：2015-08-12

基金項目：國家自然科學基金資助項目(41474026)；淮南礦業(集團)有限責任公司資助項目(HNKY-JTJS(2013)- 28)；安徽省國土資源廳科技資助項目(2011-K -22；2011-K-18) 。

第一作者簡介：呂偉才(1965—)，女，山東青島人，副教授，碩士生導師，研究方向為GNSS測量與數據處理。

中圖分類號：P228.1

文獻標志碼：A

文章編號：2095-4999(2016)02-0047-06

引文格式：呂偉才，蔣法文，杭玉付，等.改善移動終端測量精度的卡爾曼濾波算法[J].導航定位學報，2016，4(2)：47-52.(LV Weicai，JIANG Fawen，HANG Yufu，et al.Kalman filter algorithm for improving mobile terminal measurement accuracy[J].Journal of Navigation and Positioning，2016，4(2)：47-52.)DOI：10.16547/j.cnki.10-1096.20160210.