手機中衛星導航數據的精度及隨機特性

曾樹林，匡翠林，李燕杰

（中南大學 地球科學與信息物理學院，長沙 410083）

0 引言

智能手機作為1 種使用最普遍的全球衛星導航系統（global navigation satellite system, GNSS）設備，其最大的不足是定位精度低。在2016 年以前，用戶僅能從手機中獲取位置、速度和時間，以及衛星高度角、方位角等信息[1]，難以利用手機中的GNSS 原始觀測數據開展研究，這在一定程度上限制了手機端GNSS 定位精度的提升。2016 年5 月，谷歌公司宣布從安卓（Android）7.0 起，對外開放GNSS 的原始觀測數據，用戶可直接獲得手機偽距、載波、多普勒和載噪比等原始觀測信息[2]。目前已有一些有關手機中GNSS 數據質量分析和定位效果評估的研究，但對數據質量的分析大多都停留在信號載噪比、觀測數據殘差和占空比（duty cycle）對周跳的影響等方面[3-4]，對觀測數據隨機誤差特性仍缺乏深入分析。精確分析觀測數據隨機誤差特性，將有利于在高精度定位時建立更合適的隨機模型。

分析GNSS 導航數據隨機誤差特性的方法有極大似然估計[5]、方差分量估計[6]和阿倫（Allan）方差[7]等。Allan 方差因計算簡單，能夠識別多種隨機噪聲并提取噪聲特性參數，已廣泛地應用到鐘穩定度和慣導隨機誤差特性的分析中，目前也逐漸應用到GNSS 定位誤差或觀測值誤差的分析中。文獻[7]利用Allan 方差分析了多種定位模式下GNSS 的定位誤差，證明了Allan 方差用于分析GNSS 隨機誤差的可用性和有效性。文獻[8-9]利用Allan 方差，分別對測量型接收機中全球定位系統（global positioning system, GPS）和北斗衛星導航系統（BeiDou navigation satellite system, BDS）觀測值驗后殘差隨機特性進行了分析。文獻[10]利用Allan 方差驗證了GNSS 導航數據未模型化誤差的隨機特性。為深入分析手機中GNSS 導航數據的誤差特性，本文運用超短基線差分法，計算出手機中雙頻GNSS 觀測數據的單差殘差，以單差殘差作為觀測數據的誤差，分析其精度，并運用Allan 方差分析殘差隨機特性。

1 觀測數據殘差提取和分析方法

本文首先計算出觀測數據單差殘差，再利用Allan 方差分析單差殘差隨機特性。

1.1 雙差殘差恢復單差殘差原理

雙差法是先將2 個測站的GNSS 觀測值通過差分計算，得到站間單差差分值，再將站間單差差分值，在不同衛星間進行差分計算。在零基線或超短基線下，雙差法消除了衛星星歷誤差、大氣延遲和接收機鐘差等誤差，僅剩余接收機內部相關誤差。在雙差法的計算過程中，會引入參考星誤差影響，不便于分析單顆衛星誤差特性[11]，故一般采用站間單差法計算殘差。在利用站間單差法計算殘差時，得到的站間差分值中，還含有接收機相對鐘差和相對模糊度這2 個未知參數，需要采用最小二乘方法估計。而采用雙差殘差恢復單差殘差[12]的方法，可直接將雙差殘差轉換為站間單差殘差，無需進行參數估計，計算更簡單。

若已知1 個歷元站間的單差（single-difference,SD）殘差序列為（i=1, 2,…, n，n 為當前歷元共視衛星數），則以1 號衛星為參考星時，由單差殘差計算雙差（double-difference，DD）殘差的公式為

已知雙差殘差時，由式（1）即可反算出單差殘差，公式推導過程可參見文獻[12]。

1.2 Allan 方差計算及分析方法

Allan 方差[13]開始主要用于振蕩器頻率穩定度的評估，現多用于慣導器件的隨機誤差建模，目前也用于GNSS 數據隨機誤差分析。

1.2.1 Allan 方差計算過程

Allan 方差的基本思想是對完整數據序列以一定采樣率進行采樣，計算該采樣率下的Allan 方差值，再改變采樣率依次計算各采樣率下的Allan 方差值，以不同采樣率表征序列隨機特性。Allan 方差的采樣方式有重疊采樣和非重疊采樣2 種，其中重疊采樣的計算流程為：

1）設數據序列總時長為T，時間間隔為0τ ，包含數據點個數為N。以m 個數據為采樣間隔，將數據序列按重疊采樣方式劃分為M = N －m +1個數據子集，對應采樣時間為 τ= m·τ0。

4）計算全部數據子集平方值的平均值，即可得出對應采樣率下的Allan 方差計算值2σ 。Allan方差的計算公式為

5）改變采樣間隔，重復步驟1）～步驟4），計算不同采樣間隔下的Allan 方差值。

1.2.2 Allan 方差隨機誤差識別方法

隨機誤差的類型通過Allan 方差平方根值的雙對數（log-log）圖曲線斜率進行識別。常見的隨機誤差有白噪聲（white noise，WN）、閃爍噪聲（ficker noise，FN）、隨機游走（random walk，RW）和1 階高斯馬爾科夫過程（Gauss-Markov，GM）等。根據隨機誤差的理論功率譜密度（power spectral density，PSD），Allan 方差理論值的計算式[14]為

式中 SΩ( f )為隨機誤差的PSD。根據理論Allan 方差即可得出log-log 圖理論斜率，并提取噪聲特性參數，如表1 所示。

表1 4 種典型隨機誤差特性表

2 數據采集及殘差計算

2.1 實驗數據采集

實驗采用2 臺型號相同的雙頻手機小米8，實驗場地為中南大學某辦公樓天臺，實驗場地視野開闊。由于手機中內置GNSS 天線的準確位置未知，假定天線位于手機正面上邊緣中心。將2 臺小米8 手機并排緊貼放置，基線長約為7 cm，可視天線位于同一位置，實驗布置如圖1 所示。GNSS數據記錄軟件為Geo++ RINEX Logger(v2.1.3)[15]，該軟件可直接將手機中GNSS 數據記錄為與接收機無關的交換格式（receiver independent exchange format，RINEX），數據采樣率為1 Hz，采集時長6 h。實驗前，采用2 臺天寶R9 接收機，以相對定位技術測得手機觀測位置參考坐標，用于后續固定坐標計算觀測數據殘差；因手機準確天線相位中心未知，不能將測量型天線與手機天線相位中心位置精確對中，故不考慮測量型天線對中和天線高誤差。

圖1 手機GNSS 數據采集實驗

2.2 觀測數據單差殘差計算

手機導航觀測數據的單差殘差值通過雙差恢復單差法計算得出。首先采用爾特克利布（RTKLIB）軟件固定（Fixed）解算模式計算雙差殘差，在Fixed 模式中，將基準站和流動站坐標均固定為已知值，只對其他測量誤差參數進行估計，輸出結果為觀測值殘差，可專門用于殘差分析。其中固定坐標由測量型接收機采用相對定位方式測量得到。衛星截止高度角設為5°。計算出雙差殘差后，根據雙差恢復單差原理，從GPS、BDS、伽利略衛星導航系統（Galileo navigation satellite system,Galileo）和格洛納斯衛星導航系統（global navigation satellite system, GLONASS）這4 個系統中選取每個系統高度角最大的衛星作為參考星，然后將雙差殘差轉換為單差殘差。

2.2.1 L1 頻率偽距和相位殘差

手機中L1 頻率（為表述方便，Galileo E1、BDS B1 頻率均稱為L1 頻率）導航觀測數據的偽距和相位單差殘差結果如圖2 所示。從圖2 可以看出，不同系統導航數據的的偽距和相位殘差基本分布在零參考線附近，均無明顯的系統偏差或線性趨勢項。從偽距殘差的數值上看，GPS 數據的殘差大小在10 m 左右，Galileo 和BDS 在5 m 左右，GLONASS 中約5%的殘差超過了15 m。對于相位殘差，GPS 中97%的殘差集中在20 mm 以內，但部分殘差的數值波動較大，其波動幅值達到了50 mm，Galileo 的相位殘差在10 毫米級，GLONASS 和BDS 在15 mm 左右，但BDS 中存在較多粗差。

圖2 手機中L1 頻率導航數據的偽距和相位單差殘差結果

2.2.2 L5 頻率偽距和相位殘差

小米8 手機目前僅支持接收GPS 和Galileo衛星的雙頻信號，第2 個頻率為L5/E5a 波段，為表述方便，以下均稱為L5 頻率。但手機中雙頻數據可用率還較低，本實驗中僅接收到6 顆GPS和3 顆Galileo 衛星的L5 頻率信號，其原因是目前能夠播發L5 頻率信號的衛星還較少，如GPS中僅Block IIF 衛星播發L5 頻率信號。圖3 為L5 頻率導航觀測數據的偽距和相位單差殘差，同樣地，L5 頻率觀測數據的殘差序列也無明顯的系統偏差和線性趨勢項。GPS 數據的偽距殘差大小在5 m 以內，相位殘差在20 mm 左右。實驗中接收到Galileo 衛星L5 頻率信號的時間比L1頻率更短，可用L5 頻率觀測數據很少。在數值上，L5 頻率偽距殘差在3 m 以內，相位殘差在10 mm 左右。

圖3 手機中L5 頻率導航數據的偽距和相位單差殘差結果

2 臺手機型號和觀測環境完全相同，可認為觀測精度相同，根據誤差傳播定律，將單差殘差除以即可得到單臺手機觀測精度[16]。表2 為手機中所有衛星不同頻率觀測數據的標準方差平均值。對于偽距觀測數據，L5 頻率數據的精度均高于L1頻率，L1 頻率數據中除GLONASS 外，其他系統數據的精度在2 m 左右，L5 頻率數據的精度優于1 m。而相位觀測數據中雙頻數據的精度無明顯差別，數值均在3～6 mm 范圍內。但Galileo L5 頻率相位數據的結果比L1 頻率略差，這可能是由于Galileo 觀測數據過少，導致統計精度出現了偏差。

表2 所有衛星不同頻率偽距及相位精度統計結果

此外，由于接收到L1 頻率數據的衛星數比接收到L5 頻率數據的衛星數多，且部分衛星中L1頻率數據比L5 頻率數據的觀測時段更長，這可能會對2 個頻率數據的精度統計產生影響。對此，需選取同時包含L1 和L5 2 個頻率的衛星導航數據，在相同衛星數和相同時段長度的條件下進行精度統計。因Galileo 衛星觀測時段過短，不利于統計分析，則僅對GPS 的2 個頻率數據在相同條件下進行統計，結果如表3 所示。

表3 相同條件下GPS 雙頻偽距及相位精度統計結果

對比表3 與表2，對于相同頻率數據的統計結果，2 個表中偽距精度結果僅有厘米級差別，相位精度僅亞毫米級差別。在表3 中，GPS L5 頻率數據的偽距精度仍為L1 頻率的3 倍左右，L5 頻率的相位精度比L1 頻率約高1 mm，該結果與表2結果基本一致。

綜上所述，手機中Galileo 導航數據的偽距和相位精度均較高，但不足的是本次實驗中可用的Galileo 數據過少，GPS 和BDS 數據的精度次之，GLONASS 精度最低。L5 頻率數據中偽距和相位精度整體上均高于L1 頻率數據。

不斷強化政策制度保障。結合修訂制定國防動員、民兵工作相關法律法規，從國家層面進一步明確經濟功能區和企業在國防建設方面的法律責任，強化頂層設計，提高制度權威。各地區、各行業健全完善地方性、行業性法規制度，對后備力量建設相關內容作出具體規定，強化剛性措施，規范經濟功能區武裝工作秩序。探索企業稅收減免、發展資金扶持等經濟利益與履行國防義務掛鉤的措施辦法，引導企業自覺履行國防義務、主動支持民兵建設。堅持嚴格執法，適時通報批評落實武裝工作不力的黨政機關和個人，追究相關領導責任，對拒絕履行國防義務的企業實施嚴厲的行政和經濟處罰，形成震懾效應，督促各級黨政機關和企業切實抓好民兵建設。

3 殘差隨機特性分析

3.1 殘差自相關性檢驗

首先利用簡單的自相關性分析檢驗殘差隨機特性，采用矩陣實驗室（MATLAB）自相關函數計算偽距和相位殘差自相關系數。以G32 號衛星為例，圖4 為不同頻率的偽距和相位殘差在200 s 相關時間內的自相關系數圖。

圖4 G32 不同頻率導航數據的偽距和相位殘差自相關系數

從圖4 可以看出，對于L1 頻率，偽距殘差和相位殘差自相關系數迅速減小到0.2 以內，無明顯自相關性，但偽距比相位相關性波動更大。L5 頻率導航數據的偽距殘差在20 s 內自相關系數均超過了0.4，并在100 s 后才逐漸減小至0.3 以下，存在較強的時間相關性，但相位殘差基本無自相關性。這些結果表明，手機中GNSS 導航數據并不只包含白噪聲，還含有有色噪聲，需要采取新方法進一步對有色噪聲類別進行識別。本文以Allan 方差方法分析觀測數據殘差隨機特性。

3.2 殘差Allan 方差分析

利用Allan 方差分析觀測數據殘差隨機特性。在圖2 和圖3 所示的整個觀測過程中，單顆衛星殘差序列會出現不連續，但文獻[17]指出序列的不連續性不會對Allan 方差分析結果造成影響，故可對單顆衛星計算Allan 方差。另外，粗差會引起Allan 方差圖像變形，進而影響GM 這類噪聲的識別和噪聲參數的提取[17]，在計算Allan 方差前以3σ 原則剔除粗差。

以G32 衛星為例，其偽距殘差Allan 方差log-log 曲線如圖5 所示，根據表1 方法，識別噪聲類型并提取特性參數值。對于 L1 頻率導航數據的偽距殘差，可以看出，在采樣時間小于5 s 的初始階段，曲線斜率與斜率為的WN 參考線一致，殘差表現為白噪聲。白噪聲特性參數是在采樣時間τ=1 s 處取得的，對應的WN 參數值為W=2.18 m。在采樣時間大于6 s 以后，曲線斜率趨于0，表現出閃爍噪聲FN，FN 的噪聲參數是在曲線斜率為0 的任意采樣時間下取得的，圖5 中斜率最早趨近于0 的時刻為τ=6 s，此時對應的FN 噪聲參數B=1.30 m。但在采樣時間達到100 s 后，噪聲類型不再表現為FN。對于L5 頻率數據，在采樣時間小于10 s 時，Allan 方差曲線幾乎與WN 參考線完全重合，表現出較純凈的WN，對應的噪聲參數值為W=0.63 m。在采樣時間大于10 s 以后，曲線斜率為0，此后一直表現為FN，取斜率最早趨近于0的時刻，即τ=10 s，此時對應的噪聲參數B=0.22 m。

圖5 G32 不同頻率數據的偽距殘差Allan 方差曲線

FN 是1 種低頻有色噪聲，FN 的存在表明殘差序列存在時間相關性，而在殘差自相關性分析中，G32 衛星L5 頻率數據的偽距中表現出的自相關性也表明殘差序列中存在時間相關性，這與Allan 方差分析結果一致，證明了Allan 方差結果的正確性。另一方面，在自相關分析中，G32 L1 偽距殘差自相關性不明顯，但在Allan 方差分析中表明其存在有色噪聲，即存在時間相關性，說明Allan方差對于識別相關誤差具有更強的敏感性。

運用上述方法分析所有衛星的偽距和相位殘差，其Allan 方差曲線如圖6 所示。從圖6 可以看出，對于同一類觀測數據，各顆衛星觀測數據殘差的Allan 方差曲線走勢基本一致，說明手機中不同衛星同類觀測數據誤差特性相似，但在數值大小上有所差別，主要是由于衛星高度角不一致或其他因素引起的觀測數據精度不同。在曲線初始階段，曲線斜率約為，L1 與L5 偽距和相位殘差均主要表現為WN。隨著采樣間隔的增大，2 個頻率偽距殘差斜率趨近于0，開始出現FN，尤其是L5 頻率數據的偽距殘差中表現出了更明顯的FN。而對于相位殘差，L1 頻率的結果中部分衛星的殘差表現出較弱的FN，L5 頻率的結果中所有衛星的殘差均一直表現為較純凈的WN，未表現出明顯的有色噪聲。

圖6 不同頻率數據的偽距和相位殘差Allan 方差曲線（圖中虛線為WN 參考線）

不同頻率偽距和相位噪聲特性參數統計值如表4 所示。從偽距中WN 的噪聲參數值來看，Galileo L1 頻率數據中的WN 值最小，GPS 和BDS次之，GLONASS 最大，這與偽距殘差標準方差分析結果一致。同樣，L5 頻率數據的WN 值也是Galileo 最小。另外，除GLONASS 外，其余各系統的偽距殘差還表現出FN，其噪聲參數在數值上有分米級差別。而雙頻相位觀測數據中的噪聲類型均主要為WN，數值上仍是Galileo 最小，但不同系統之間的數值僅相差1～2 mm。

表4 不同頻率偽距和相位噪聲特性參數統計表

以上分析結果中，手機中偽距和相位2 類觀測數據包含的噪聲類型不一致，其中偽距包含WN 和FN 2 類噪聲，而相位則主要包含WN。偽距中存在明顯有色噪聲可能源于以下2 個原因：①手機中使用的Geo++ RINEX Logger 數據采集軟件已經對偽距進行了平滑處理，即采用相位平滑偽距技術減小了偽距誤差，部分測量型接收機中也具有這種功能，如天寶R9；②多路徑效應對測量結果的影響，理由是通過超短基線差分不能完全消除觀測值數據中的多路徑誤差，并且手機中采用的線性極化GNSS 天線，導致GNSS 信號中更易產生多路徑誤差，即使是在弱多路徑環境下，偽距數據也受到了較強的多路徑影響。除多路徑效應外，文獻[7]指出，殘差中殘余的對流層延遲濕分量誤差也可引起有色噪聲，而本實驗中2 臺手機的基線長僅幾厘米，可認為大氣延遲誤差完全被消除，故本實驗中偽距所包含的有色噪聲不是來源于大氣殘余誤差。對手機中GNSS 導航數據包含的有色噪聲的來源有待進一步驗證。

4 結束語

本文分析了手機中雙頻多系統GNSS 原始觀測數據的精度，并采用Allan 方差分析了其隨機特性，得出如下主要結論：

1）從精度統計結果來看，就不同頻率而言，L5 頻率觀測數據的偽距和相位精度整體上高于L1頻率。對于偽距觀測數據，L1 頻率數據的精度在2 m 左右，L5 頻率數據精度優于1 m。而相位觀測數據中，雙頻精度均為3～6 mm。對于手機中不同衛星系統的觀測數據，Galileo 的偽距和相位精度均較高，GPS 和BDS 次之，GLONASS 精度最低。

2）Allan 方差分析結果表明，手機中GNSS 觀測數據并不只含白噪聲。在弱多路徑環境下，L1 和L5 頻率數據的偽距中均包含WN 及FN。L1 頻率相位中主要包含WN，但存在較弱的FN，而L5 頻率相位中僅包含較純凈的WN。

手機中內置GNSS 天線為線性極化天線，這類天線對多路徑的抵抗能力較差，這可能是觀測數據產生有色噪聲的主要原因，對于手機中GNSS 信號的多路徑效應有待于深入研究。另外，傳統的卡爾曼濾波要求過程噪聲和觀測噪聲均為白噪聲，而本文結果表明，手機接收的GNSS 數據中存在有色噪聲，因此僅考慮白噪聲的手機端GNSS 定位算法是不夠嚴密的。在建立適用于手機的GNSS 坐標參數估計隨機模型時，需要考慮有色噪聲。