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

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

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

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

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

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

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

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

0引言

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

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

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

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

1Kalman濾波方程

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

1.1Kalman濾波的狀態(tài)方程

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

(1)

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

(2)

(3)

則式(1)可寫成

(4)

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

(5)

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

(6)

式(6)即為采用CORS RTK模式進行動態(tài)監(jiān)測時，動態(tài)系統(tǒng)的Kalman濾波狀態(tài)方程。

1.2Kalman濾波的觀測方程

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

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

以CORS RTK測量獲得的監(jiān)測點3維位置為觀測值時，某一監(jiān)測點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為監(jiān)測點i在第k+1歷元的觀測時刻;tk+1為監(jiān)測點i本次測量時所有歷元(如30個歷元)觀測的中心時刻。

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

狀態(tài)方程(6)和觀測方程(7)共同構(gòu)成了監(jiān)測點3維位置序列的Kalman濾波模型

(8)

式中：Φk，k-1為k-1到k歷元的轉(zhuǎn)移矩陣；Γk,k-1為系統(tǒng)噪聲矩陣；Ωk-1為k-1歷元的系統(tǒng)噪聲；Bk為k歷元系統(tǒng)的觀測矩陣；Δk為k歷元系統(tǒng)的觀測噪聲；Xk為k歷元的系統(tǒng)待估狀態(tài)參數(shù)；Lk為k歷元系統(tǒng)的觀測向量矩陣；Xk和Lk均為監(jiān)測點的3維位置和速度向量。

2位置序列的Kalman濾波算法

2.1Kalman濾波遞推公式

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

(9)

式中I為單位矩陣，并有

(10)

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

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

2.2濾波初值的確定

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

1)狀態(tài)參數(shù)初值的確定。

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

(11)

式中以X分量為例，有：

即狀態(tài)參數(shù)初值為

(12)

2)方差陣初值的確定。

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

(13)

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

(14)

式中：

由協(xié)方差陣傳播定律，狀態(tài)參數(shù)中的3維速率的初始方差陣為

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

(15)

式中△t為相鄰歷元間的采用間隔，則狀態(tài)參數(shù)的初始協(xié)方差陣為

(16)

3)動態(tài)噪聲的初始方差陣。

當以測點的3維位置和瞬時速率作為狀態(tài)參數(shù)時，則瞬時加速率為動態(tài)噪聲，根據(jù)協(xié)方差陣傳播定律，其方差陣為

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

(17)

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

3試驗與結(jié)果分析

3.1試驗簡介

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

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

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

3.2數(shù)值結(jié)果與分析

(18)

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

則平均值中誤差為:

有

(19)

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

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

表1　RTK測量平面精度統(tǒng)計信息

表2　RTK測量高程精度統(tǒng)計信息

表3　高程測量差值統(tǒng)計信息

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

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

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

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

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

4結(jié)束語

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

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

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

參考文獻

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

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

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

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

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

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

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

[8]余學祥，呂偉才.GPS監(jiān)測網(wǎng)動態(tài)數(shù)據(jù)處理抗差Kalman濾波模型[J].中國礦業(yè)大學學報，2000，29(6)：553-557.

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

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

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

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

[13]徐彥田，程鵬飛，蔡艷輝，等.單頻RTK動態(tài)解算的卡爾曼濾波算法研究[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]余學祥，呂偉才，柯福陽，等.煤礦開采沉陷自動化監(jiān)測系統(tǒng)[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)；淮南礦業(yè)(集團)有限責任公司資助項目(HNKY-JTJS(2013)- 28)；安徽省國土資源廳科技資助項目(2011-K -22；2011-K-18) 。

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

中圖分類號：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.