一種雙星故障條件下RAIM可用性評估的改進方法

馬下平 余科根 賀小星 嚴 麗 趙立都

1 西安科技大學測繪科學與技術學院，西安市雁塔路58號，710054 2 中國礦業大學環境與測繪學院，江蘇省徐州市大學路1號，221116 3 江西理工大學土木與測繪工程學院，江西省贛州市紅旗大道86號，341000 4 江西省防震減災與工程地質災害探測工程研究中心，南昌市廣蘭大道418號，330013 5 重慶交通大學土木工程學院，重慶市學府大道66號，400074

BDS是我國自主建設、獨立運行的全球衛星導航系統[1-4]，可為用戶提供更多高定位精度和可靠性的導航信息，同時也可為民用航空接收機自主完好性監測(receiver autonomous integrity monitoring，RAIM)提供有力支撐。RAIM對衛星進行故障監測會受到衛星個數和衛星幾何分布的影響。當可見衛星數過少或幾何分布不理想時，由衛星幾何結構引起的定位誤差會遮掩衛星故障引起的定位誤差，造成完好性監測結果不可信，漏檢率增加。因此在進行故障監測之前必須進行可用性判斷，以保證不會影響故障監測性能[5]。RAIM可用性計算中涉及到漏檢率和誤警率參數，王爾申等[6]研究發現，RAIM可用性隨漏檢率和誤警率的降低而降低。

目前，基于單星故障的RAIM可用性評估方法主要有水平保護法(horizontal protection level method，HPLM)、最大水平精度因子法(δHDOPmax)和近似徑向誤差保護法(approximate radial-error protection method，ARPM)，理論上這3種方法具有等價性[7-8]；考慮到2顆衛星同時發生故障的可能性，倪育德等[9]從1顆衛星故障的可用性出發，推導出2顆衛星故障的可用性；針對多故障衛星的可用性，陳金平等[10]提出基于漏檢概率和圓概率誤差的可用性分析方法。隨著GNSS的發展與應用以及飛機更高階段精密進近的需求，具備垂直導航能力的高級接收機自主完好性監測(advanced RAIM，ARAIM)被引入。ARAIM旨在全球范圍內提供垂直導向-200 ft定標性能(LPV-200)級別的航空導航服務，眾多學者對此進行了相關算法改進和定位性能方面研究[11-17]。

現有的雙星故障條件下RAIM可用性評估大多沿用已有的矩陣最大特征值方法(matrix maximum eigenvalue method，MMEM)[18]，該方法可解決2個粗差比值的計算，但該值還可采用高等數學中求極值方法得到。為此，本文提出針對RAIM可用性評估的極大值方法(maxima method，MM)，并基于完好性風險參數和IGMAS中8個GNSS站的實測數據對這2種方法的RAIM可用性進行驗證，得到最新的中國境內BDS可用性性能。

1 RAIM方法原理

1.1 RAIM觀測方程建立

RAIM主要是基于GNSS中的偽距觀測量。若在歷元t時接收機r與衛星s進行同步觀測，并觀測m顆衛星，則觀測方程可表示為：

y=Hx+ε

(1)

式中，y為m維偽距觀測值與近似衛地距的差向量；H為m行4列設計矩陣，表示各衛星對用戶的投影向量；x為測站坐標和接收機鐘差構成的四維向量，ε為m維觀測偽距的誤差向量。

給定觀測量權陣為P，則式(1)的最小二乘解為：

=x+(HTPH)-1HTPε

(2)

=(I-H(HTPH)-1HTP)ε

(3)

令

(4)

(5)

1.2 檢驗統計量構建

(6)

(7)

一般取檢驗統計量(T)為：

(8)

在導航過程中，將實時解算的T值與監測門限(TD)進行比較，TD計算可參考文獻[19]。若T

2 RAIM可用性分析方法

在進行RAIM之前，首先根據性能指標對當前可見星的幾何分布進行判斷，分析是否適合進行完好性監測，即判斷RAIM算法是否可用。

2.1 監測統計量與定位誤差的關系

RAIM可用性判斷時需要構建T與定位誤差之間的數學關系。由于RAIM不具備垂直導航能力，因此本文只討論水平保護級別(horizontal protection level，HPL)。

圖1中橫軸表示T，縱軸表示水平徑向誤差(horizontal radial error，HRE)，將這2個變量取比值后得到圖中斜線斜率(SLOPE)，即

(9)

圖1 T與HPL之間的關系Fig.1 The relationship between T and HRE

從圖1可以看出，T與HRE呈線性關系，沿此傾斜線，當偽距觀測量中存在偏差時，T與HRE線性增加。偏差較小時，HRE小于HPL且T小于TD，系統處于正常狀態(圖1左下角)，PMD較低；偏差較大時，HRE大于HPL且T大于TD，系統能正確監測偏差(圖1右上角)，PMD也較低；當偏差處于兩者之間時，PMD可能達到最大(圖1左上角)。對應相同T值時，SLOPE越大，衛星發生故障時MD概率越高。對應衛星的SLOPEs為：

(10)

在m個觀測量中，HPL可表示為：

HPL=max(SLOPEs)·T=SLOPEmax·T

(11)

RAIM可用性分析就是利用SLOPEs在各歷元的最大值SLOPEmax進行判定，依據SLOPEmax求出HPL，并與水平告警限值(horizontal alarm limit，HAL)進行比較，如果HPL

2.2 雙星故障時MM求解HPL

在忽略觀測偽距隨機誤差的情況下，當觀測誤差ε在第i個和第j個觀測值中存在粗差bi和bj時，則有：

(12)

根據式(5)可以得到定位解的誤差向量為：

(13)

則有：

(14)

將式(12)代入式(5)可得：

(15)

(16)

以往對k值主要采用矩陣最大特征值方法(MMEM)[9,18,20]進行求解，但其也可采用高等數學中求極值的方法獲取，具體如下：

為求式(16)的極大值，兩邊取對數，并對式中k求偏導，可得一元二次方程式：

ak2+bk+c=0

(17)

設k1和k2分別為式(17)兩根，將k1、k2分別代入式(16)，則當前歷元的(SLOPEij)max可表示為：

(SLOPEij)max=

max((SLOPEij)k1,(SLOPEij)k2)

(18)

即可得到雙星故障時的HPL為[19]：

(19)

2.3 雙星故障時MMEM求解HPL

目前求解HPL主要采用MMEM方法，該方法主要是構建關于粗差向量b的二次型矩陣，根據矩陣的最大值為矩陣最大特征值這一定律求解HPL，具體計算過程可參考文獻[18]。

3 數據處理

3.1 實驗數據和處理策略

實驗數據采用國際GNSS監測評估系統(international GNSS monitoring and assessment system，IGMAS)提供的中國境內BJF1、CHU1、GUA1、KUN1、LHA1、SHA1、WUH1、XIA1共8個連續跟蹤站2020-09-06的觀測數據和星歷文件。通過自編程序統計這些測站在截止高度角5°以上的BDS衛星可見性、DOP值和基于2種方法的雙星故障HPL值，比較并分析中國境內BDS可用性和2種雙星故障RAIM可用性的差異。具體數據情況及處理策略為：1)BJF1和CHU1站觀測時段為00:00:00～09:59:30；GUA1、KUN1、LHA1、SHA1站觀測時段為00:00:00～15:59:30；WUH1和XIA1站觀測時段為00:00:00～16:59:30；2)考慮到每個測站不完全包含雙頻或三頻觀測值，因此采用BDS B1I頻率的單頻觀測值；3)全球星歷文件中每隔1 h給出用于計算每個歷元的BDS衛星軌道參數、衛星鐘差參數和衛星狀態信息，其1 d的全球星歷文件可從MGEX網站下載；4)電離層延遲采用Klobuchar模型進行改正，對流層延遲采用Saastamoinen模型進行改正；5)式(2)中P采用已知的高度角定權法[21]，定權時考慮3種異構衛星之間的差異，將GEO和IGSO+MEO衛星的權比關系設置為1∶5[22]。

3.2 衛星可見數和DOP值

利用8個IGMAS站的BDS觀測數據和全球星歷文件，計算出每個測站每個歷元對應衛星的高度角，并統計可見衛星數、平均衛星數、GDOP值及其平均值。結果表明：1)中國境內8個IGMAS站的BDS可見衛星數最小值為17顆，完全能夠滿足RAIM可用性評估；2)GDOP最小值、平均值和最大值的均值分別為1.117、1.437、2.186，說明BDS衛星在中國境內空間分布情況很好。

3.3 2種方法計算雙星故障時HPL比較

目前BDS在中國區域的定位精度為5m左右[4]。保守起見，將σ0設置為8 m，PMD設置為10-3，PFA設置為10-5[19]。利用2種方法分別計算BDS雙星故障時HPL值，圖2～5分別為8個IGMAS站的HPL坐標時間序列，表1(單位m)

圖2 雙星故障條件下BJF1和CHU1站HPL時間序列Fig.2 HPL time series of BJF1 and CHU1 stations under double-satellite faults conditions

圖3 雙星故障條件下GUA1和KUN1站HPL時間序列Fig.3 HPL time series of GUA1 and KUN1 stations under double-satellite faults conditions

圖4 雙星故障條件下LHA1和SHA1站 HPL時間序列Fig.4 HPL time series of LHA1 and SHA1 stations under double-satellite faults conditions

圖5 雙星故障條件下WUH1和XIA1站HPL時間序列Fig.5 HPL time series of WUH1 and XIA1 stations under double-satellite faults conditions

表1 2種方法計算的HPL值統計

為8個IGMAS站的HPL最小值、最大值和平均值。通過分析可知：

1)從整體上看，MMEM計算的HPL小于MM計算的HPL，表明MMEM得到的RAIM可用性高。利用MM求解的8個IGMAS站的HPL最小值、平均值、最大值的均值分別為124.921 m、233.141m、3 020.005 m；MMEM求解的HPL最小值、平均值和最大值的均值分別為219.093 m、515.110 m、1 302.518 m。

2)2種方法得到的HPL時間序列中，大部分歷元的HPL值一致性很好，HPL平均值的差值約為281.969 m，主要原因可能為這2種方法的計算公式中所用參數個數和精度存在差異。

3)在中國區域，雙星故障條件下HPL值的空間分布存在差異，表明RAIM可用性與地理位置相關。

3.4 計算時間比較

表2(單位s)為2種方法的計算時間，從表中可以看出，MM和MMEM的平均計算時間分別為11.733 s和21.387 s，MM比MMEM計算時間短。

3.5 2種方法雙星故障時可用性統計

飛機在飛行階段包含非精密進近(non-precision approach，NPA)、終端、本土航路和遠洋航路4個飛行階段[23]，4個階段的HAL值分別為555.6 m、1 852 m、5 556 m和7 408 m。基于2種方法得到雙星故障條件下HPL值，根據每個測站計算出觀測歷元中各階段的可用性百分比，再計算8個測站的可用性比值均值，得到每個階段的可用性結果(表3)。從表中可以看出：1)MM計算的可用性在航路和遠洋階段達到100%，而MMEM計算的RAIM可用性除在NPA階段外，其他階段的可用性均達到100%；2)結合所計算的HPL值，雖然MMEM在NPA階段的可用性不及MM，但整體上MMEM可用性優于MM。

表2 2種方法的計算時間比較

表3 單星和雙星故障條件下RAIM可用性統計

4 結 語

本文從雙星存在故障的前提條件出發，在總結已有的最大特征值計算RAIM可用性方法的基礎上，提出RAIM可用性評估的極大值方法，并以2020-09-06中國區域8個IGMAS連續跟蹤站的BDS觀測數據為例，比較和分析中國境內HPL變化及RAIM可用性情況，得到以下結論：

1)中國境內可見衛星數和GDOP完全能夠滿足可用性計算需求，中國境內BDS目前至少可觀測到17顆可見衛星，GDOP值最大為2.186，可為RAIM可用性計算和故障監測提供豐富的數據源。

2)雖然MM計算的RAIM可用性在NPA階段高于MMEM，且耗時較少，但總體來說，MMEM計算HPL的時間序列小于MM，即對應的RAIM可用性更高。

3)2種方法計算的HPL存在差異，主要原因可能為各方法所用參數個數及本身精度存在差異，且RAIM可用性計算中涉及到完好性風險參數，不同的參數初值會對RAIM可用性評估產生不同影響，因此參數取值應引起重視。