基于新息速率抗差估計的INS/GNSS組合導航系統(tǒng)欺騙檢測算法

張 超，呂志偉,2，張倫東，高揚駿,2

（1. 信息工程大學，鄭州 450000；2. 地理信息工程國家重點實驗室，西安 710054）

慣性導航系統(tǒng)（Inertial Navigation System, INS）和全球衛(wèi)星導航系統(tǒng)（Global Navigation Satellite System, GNSS）的組合導航系統(tǒng)在軍用和民用領域得到了廣泛而深入的應用。但由于GNSS 信號微弱且民碼結構公開，使其容易受到欺騙干擾的影響。欺騙干擾是指干擾源播發(fā)虛假衛(wèi)星信號，誘導目標接收機進行捕獲與跟蹤，進而輸出錯誤信息的技術[1]。GNSS欺騙技術的分類方式有很多，按攻擊階段分類可分為捕獲階段欺騙干擾和跟蹤階段欺騙干擾，而跟蹤階段欺騙干擾按欺騙信號的碼相位和載波頻率與真實信號是否保持同步，又可分為同步式欺騙干擾和異步式欺騙干擾[2]。在欺騙干擾環(huán)境下，組合導航系統(tǒng)的GNSS 模塊鎖定欺騙信號并輸出錯誤的偽距和偽距率，進而影響測量更新階段的狀態(tài)誤差估計值并輸出錯誤的導航結果。因此，欺騙檢測對于提高組合導航系統(tǒng)的完好性和可用性顯得十分重要。

INS/GNSS 組合導航系統(tǒng)的欺騙檢測算法可分為非耦合和耦合兩類[3]。非耦合欺騙檢測是指GNSS 與INS 獨立工作，將GNSS 輸出的位置、速度轉換為加速度和角速度，與INS 輸出的加速度和角速度進行一致性檢驗。耦合欺騙檢測是指GNSS 與INS 輸出的信息通過Kalman 濾波等方式融合，通常使用殘差或新息序列作為檢驗統(tǒng)計量，檢測方式采用假設檢驗。這類方法的缺點是在欺騙干擾影響下組合導航系統(tǒng)的反饋校正機制會使殘差或新息序列偏離正常值導致誤警率和漏警率升高[4]。耦合欺騙檢測又可分為“快照法”和“連續(xù)法”。“快照法”是以當前時刻的新息序列構成檢驗統(tǒng)計量，適合檢測GNSS 測量值中階躍類型的誤差，對緩變類型的誤差檢測時間較長，例如新息欺騙檢測、殘差卡方檢測屬于該類算法，在單路衛(wèi)星信號受到0.1 m/s 緩變類型的干擾時，新息欺騙檢測算法的檢測時間為360 s[5]；“連續(xù)法”是將檢測窗口內的新息序列構成檢驗統(tǒng)計量，適合檢測緩變類型的誤差，但檢測窗口長度的設置需要在告警頻率和檢測成功率之間權衡，例如自主完好性檢測外推法屬于該類算法，在單顆衛(wèi)星受到0.1 m/s 緩變類型的干擾時，檢測窗口設置為150 s 的自主完好性檢測外推法的檢測時間為162 s[5]。

針對GNSS 測量值中存在緩變類型的干擾，Bhatti提出了一種新息速率檢測算法[6]，將新息的變化率作為檢測量，在單顆衛(wèi)星受到0.1 m/s 緩變類型的干擾時檢測時間為110 s，但該方法需要與自主完好性檢測外推法結合使用，無法針對單個GNSS 測量值檢測。針對組合導航系統(tǒng)反饋校正對耦合欺騙檢測的影響，Shizhuang Wang 將抗差估計與基于新息的檢測算法結合，單顆衛(wèi)星受到0.1 m/s 緩變類型的干擾時檢測時間縮短至65 s[7]；Chuang Zhang 提出的檢測算法在單顆衛(wèi)星受到0.5 m/s 緩變類型的干擾時漏警率為26%[8]，但這兩種算法需要設置檢測窗口；許睿等提出了MEDLL 輔助的GNSS/INS 系統(tǒng)欺騙信號檢測方法[9]，成功檢測出牽引速率為2 m/s 的欺騙信號。近幾年有學者將人工智能（例如神經網絡[10]、支持向量機[4]）、芯片級原子鐘[11]等應用到欺騙檢測算法中，但對計算能力要求較高，與現(xiàn)有系統(tǒng)兼容性不強。本文主要針對同步式欺騙干擾，其造成的影響可視為緩變類型的誤差[12]。通過將基于新息速率的檢測算法與抗差估計結合，在無需設置檢測窗口的情況下有效提高了檢測時間，降低了漏警率和誤警率。

本文首先根據(jù)同步式欺騙原理建立了基于GNSS測量值的欺騙仿真模型，并介紹了INS/GNSS 緊耦合的新息欺騙檢測算法；然后針對新息欺騙檢測算法存在的缺陷，提出了一種基于新息速率抗差估計的INS/GNSS 組合導航系統(tǒng)欺騙檢測算法；最后通過仿真驗證了新算法的有效性。

1 GNSS 欺騙仿真模型

對于INS/GNSS 組合導航系統(tǒng)的欺騙仿真分析，基于GNSS 測量值的欺騙模型可以滿足仿真需求，避免了復雜的欺騙信號仿真，容易實現(xiàn)。本文模擬同步式欺騙信號。

第i 顆衛(wèi)星在t 時刻的真實偽距測量值可表示為：

其中τ(i)為信號傳播延遲，c 為光速， δtu和δt(i)為接收機鐘差和衛(wèi)星鐘差。

異步式欺騙信號的碼相位和載波頻率與真實信號無法保持同步，所以其偽距測量值是任意的。同步式欺騙信號初期與真實信號保持同步，當接收機跟蹤環(huán)路鎖定欺騙信號后，為防止接收機自主完好性檢測，采取逐步誘導的方式使偽距測量值逐漸偏離真實值。

綜上所述，基于GNSS 測量值的欺騙模型可表示為：

tInit為欺騙信號鎖定接收機跟蹤環(huán)路的時刻，a為誘導速率，b 為真實偽距與欺騙偽距的偏差，且a≠ 0。當b= 0時， s ( t )代表同步式欺騙；當b≠0時， s ( t )代表異步式欺騙。

2 新息欺騙檢測算法

新息序列是實際觀測向量與一步預測狀態(tài)向量計算的觀測向量的差值，由系統(tǒng)模型誤差和觀測量誤差共同決定，當系統(tǒng)模型誤差較小時，新息序列可以反映觀測量誤差。

設觀測矩陣為Hk，觀測向量為 Zk，一步預測狀態(tài)向量協(xié)方差矩陣為Pk|k-1，觀測噪聲協(xié)方差矩陣為Rk，則新息序列 rk與其協(xié)方差矩陣 Vk的表達式如下：

歸一化后的新息可以表示為：

其中rki為k 時刻新息序列的第i 個值（i=1…n，n為可見衛(wèi)星個數(shù)），Vkii為k 時刻rki的方差。 ωi代表k時刻歸一化后的第i 個新息值，可以反映第i 個GNSS測量值誤差。

如果接收到的第i 路衛(wèi)星信號是真實信號，在無其他故障的情況下，第i 個測量值不存在粗差， ωi服從標準正態(tài)分布：

如果接收到的第i 路衛(wèi)星信號是欺騙信號，則第i個測量值存在粗差， ωi服從非中心正態(tài)分布：

根據(jù)實際應用對于導航系統(tǒng)完好性的需求，設定誤警率 Pfa，則對應的第i 個測量值的誤警率為[9]：

構造假設H0:ωi～N(0,1)代表第i 路衛(wèi)星信號狀態(tài)正常，H1: ωi～N(δ ,1)代表第i 路衛(wèi)星信號狀態(tài)異常，根據(jù)標準正態(tài)分布概率密度函數(shù)，由 α0可計算出門限值 Td，公式如下：

根據(jù)假設檢驗原理，可由下式判斷第i 路衛(wèi)星信號是否存在異常：

3 基于新息速率抗差估計的INS/GNSS 組合導航系統(tǒng)欺騙檢測算法

3.1 基于新息速率的檢測算法

新息速率欺騙檢測算法通過歸一化新息 ωi的變化速率來判斷GNSS 測量值是否存在異常?？紤]到測量噪聲的影響，通常采用Kalman 濾波實時估計新息ωi的速率。

Kalman 濾波估計的狀態(tài)向量定義為：

系統(tǒng)模型定義為：

觀測模型定義為：

其中， ωi作為測量值輸入， ni為觀測噪聲。

在無其他干擾情況下，新息速率 vi服從高斯分布，給定誤警率Pfd，由式(9)計算出單個測量值對應的誤警率 α0，則 vi對應的檢測門限為：

其中，Q-1是高斯分布的逆，Pvi是狀態(tài)協(xié)方差陣中 vi的方差。

同步欺騙導致GNSS 測量值誤差隨時間積累，如式(3)所示，a 越大，誤差積累越快。新息欺騙檢測算法存在最小可檢測誤差（Minimum Detectable Blunders, MDB）[8]，只有當誤差積累到MDB 大小時才能檢測出同步欺騙。由于誤差累積會導致新息呈增長或減少趨勢，新息速率欺騙檢測算法通過檢測新息增長或減少的速率判定是否存在欺騙干擾，無需等誤差累積到一定程度才能被檢測。所以，新息速率欺騙檢測算法比新息欺騙檢測算法的檢測時間更短。但計算新息速率的濾波器在收斂之前無法進行檢測。

3.2 基于新息速率抗差估計的檢測算法

組合導航系統(tǒng)的反饋校正機制會使前一時刻GNSS 測量值誤差傳遞至下一時刻系統(tǒng)更新階段的狀態(tài)向量估計值中，進而使下一時刻原本狀態(tài)正常的新息偏離正常值，對基于新息序列的檢測算法造成干擾。欺騙干擾對新息序列的影響可表達為[4]：

其中 Φk+1,k為狀態(tài)轉移矩陣，Kk為濾波增益矩陣，Δk=（Δ1…Δi…Δn）為 k 時刻欺騙干擾對GNSS 測量值的影響，Δr～k+1為k+1 時刻新息序列的偏移量。從式(16)可以看出，k 時刻欺騙干擾的影響會使k+1 時刻的新息序列發(fā)生偏移。

油菜葉片在不同生育時期呈現(xiàn)不同的形態(tài)，主要包括長柄葉、短柄葉和無柄葉3種類型。選用3種類型油菜葉片在室內和戶外環(huán)境下的圖像，如圖2(a)4幅圖像所示。利用本文方法所提取的結果(為便于觀察，將葉脈圖像作求補運算)如圖2(b)4幅圖像所示。

通過引入抗差估計可以很好地削弱這種影響。選取IGG-3 等價權函數(shù)[13]，使用新息速率計算等價權，具體公式如下：

vi為第i 個GNSS 測量值對應的新息速率， wi為第i 個測量值的等價權。將 k1設定為新息速率的門限值 vD， k0設為0.5 倍的 k1。

等價權矩陣可表示為：

文獻[7]對增益矩陣Kk進行調整，具體公式如下：

將KR代替式(16)中的Kk，可得：

從式(20)可以看出，Δk中的元素Δi增大或減少，導致k 時刻新息 ωi增大或減少，進而新息速率vi會增大，等價權 wi降低。因此無論 Δi增大或減小其對應的等價權均會降低，從而降低k 時刻的欺騙干擾對k+1 時刻新息序列的影響。

基于新息速率抗差估計的INS/GNSS 組合導航系統(tǒng)欺騙檢測算法流程如圖1 所示。組合導航濾波器輸出各路歸一化新息，計算對應的新息速率和等價權，等價權組成等價權矩陣用于計算KR，將KR替代原增益矩陣Kk輸入到組合導航濾波器。

圖1 檢測算法流程圖Fig.1 Flow chart of detection algorithm

3.3 多路衛(wèi)星信號受欺騙時的算法改進

同步欺騙影響下的GNSS 測量值誤差隨時間累積的快慢與誘導速率有關。誘導速率越小，誤差累積的越慢，對應的新息速率 vi變化的就越慢， vi到達門限值 vD的時間就越久，對應的等價權 wi減至0 的時間變久，最終導致抗差估計的抑制效果降低。尤其多路衛(wèi)星信號受到低誘導速率的同步欺騙時，抗差估計的抑制效果會顯著降低。為此，改進等價權函數(shù)中參數(shù)k1，令 k1=κvD，0＜ κ≤1 。通過調節(jié)參數(shù)κ ，可以達到縮短等價權 wi減至0 的時間，從而提高抗差估計的抑制效果。

4 仿真試驗

為簡化表達，M1表示新息欺騙檢測算法，M2表示新息速率欺騙檢測算法，M3 表示基于新息速率抗差估計的INS/GNSS 組合導航系統(tǒng)欺騙檢測算法。

針對INS/GNSS 緊耦合模式進行仿真實驗，設計了3 個算例。1）在1 顆衛(wèi)星受到同步欺騙時比較M1和M2 的檢測能力；2）在2 顆衛(wèi)星受到同步欺騙時比較M1、M2 和M3 的檢測能力；3）在3 顆衛(wèi)星受到同步欺騙時比較等價權函數(shù)中參數(shù) k1的取值對M3 檢測能力的影響。

4.1 仿真條件

設30 顆GNSS 衛(wèi)星均勻分布在6 個等圓軌道上。GNSS 模塊采用雙頻數(shù)據(jù)消除電離層影響。接收機誤差設為高斯白噪聲序列。誤警率設置為10-5，M1 的門限值為 4.95 m，M2 的門限值為 0.00276 m/s。GNSS 模塊參數(shù)見表1；IMU 模塊參數(shù)配置見表2；Kalman 濾波器配置參數(shù)見表3。

表1 GNSS 模塊參數(shù)Tab.1 Parameters of GNSS model

表2 IMU 模塊參數(shù)Tab.2 Parameters of IMU model

表3 Kalman 濾波器參數(shù)Tab.3 Parameters of Kalman filter

模擬飛機飛行軌跡，時長746 s，并據(jù)此產生GNSS 和IMU 測量值。飛行軌跡如圖2 所示。

圖2 飛行軌跡Fig.2 Flight trajectory

欺騙場景均采用同步式欺騙干擾，設置如表4 所示。

表4 欺騙場景設置Tab.4 Spoofing scene setup

4.2 仿真結果與分析

1）實驗1

通過實驗1 驗證M2 比M1 對同步欺騙的檢測時間短。M1 和M2 的仿真結果如圖3(a)和3(b)所示。

圖3 M1 與M2 仿真結果比較Fig.3 Comparison of simulation results between M1 and M2

由圖3(a)可知，M1 在300～500 s 內檢測失敗，將欺騙時間延長至 600 s，得到 M1 的檢測時間為210 s；由圖3(b)可知，M2 的檢測時間為80 s。所以，在1 路衛(wèi)星信號受到同步欺騙時，M2 比M1的檢測時間縮短了61.9%。

2）實驗2

通過實驗2 驗證抗差估計抑制新息序列偏離正常值的效果。M1、M2和M3的仿真結果如圖4(a)、4(b)和4(c)所示。由圖4(a)可知，M1 在300～500 s 內未檢測出1、2 號衛(wèi)星所受欺騙，若將欺騙時間延長至600 s，得到M1 的檢測時間為225 s；由圖4(b)可知，M2 的檢測時間為80 s，成功檢測出了1、2 號衛(wèi)星，但狀態(tài)正常的3、4 號衛(wèi)星的新息速率卻超過了門限；由圖4(c)可知，M3的檢測時間為76 s，其余衛(wèi)星的新息速率正常。

圖4 M1、M2 和M3 仿真結果比較Fig.4 Comparison of simulation results among M1, M2 and M3

為更好說明抗差估計的抑制效果，對實驗2 進行100 循環(huán)的蒙特卡洛仿真，得到三種算法的漏警率和誤警率如表5 所示。

表5 實驗2 蒙特卡洛仿真結果Tab.5 Monte Carlo simulation results of Exp.2

實驗2 表明：在2 路衛(wèi)星信號受到同步欺騙時，M2 相比M1 的檢測時間縮短了64.4%，漏警率降低至0%，但誤警率卻顯著上升。M3 不僅漏警率保持0%，誤警率也降低至0%，同時檢測時間在M2 的基礎上又縮短了5%。

3）實驗3

為驗證參數(shù) k1的取值對M3 檢測能力的影響，通過公式 k1=κvD調節(jié)參數(shù) k1的大小，κ =1 時的檢測結果如圖5(a)所示，κ =0.6 的檢測結果如圖5(b)所示。由圖5(a)可知，1 號、2 號和3 號衛(wèi)星的檢測時間分別為101 s、82 s 和108 s，但狀態(tài)正常的4 號衛(wèi)星卻超過了門限值；由圖5(b)可知，1 號、2 號和3 號衛(wèi)星的檢測時間分別為88 s、75 s 和88 s，其余狀態(tài)正常的衛(wèi)星新息速率均正常。

圖5 κ =1 與κ =0.6 的仿真結果比較Fig.5 Comparison of simulation results between κ =1 and κ =0.6

為更好說明參數(shù) k1的取值對M3 檢測能力的影響，對實驗3 進行100 循環(huán)的蒙特卡洛仿真，得到κ =1 和κ =0.6 的漏警率和誤警率如表6 所示。

表6 實驗3 蒙特卡洛仿真結果Tab.6 Monte Carlo simulation results of Exp. 3

實驗3 表明，κ =0.6 的檢測時間在κ =1 的基礎上進一步縮短了12%左右，同時誤警率降低了81%。但是漏警率有所升高。

5 結 論

本文針對基于新息的欺騙檢測算法對同步式欺騙干擾的檢測時間較長、漏警率和誤警率較高的問題，提出了一種基于抗差估計的INS/GNSS 導航系統(tǒng)欺騙檢測算法。該算法將新息的變化速率作為檢測量，并采用抗差估計的方法削弱了欺騙干擾對新息序列產生的影響，通過調參數(shù) k1，進一步提高了檢測能力。仿真實驗驗證了新算法的可行性,相比新息欺騙檢測算法，新算法對同步欺騙的檢測時間縮短了60%以上，漏警率和誤警率維持在4%以內。

未來研究工作展望：1）新算法的檢測時間、誤警率和漏警率隨欺騙施加時長、誘導速率的變化規(guī)律有待進一步研究；2）新算法在不同精度等級的慣導系統(tǒng)下檢測能力的變化有待進一步研究。