高超聲速滑翔飛行器多模型魯棒跟蹤方法

王瑞鵬，王小剛，張豪杰

（1.哈爾濱工業(yè)大學(xué)航天學(xué)院，哈爾濱 150001；2.北京電子工程總體研究所，北京 100854）

0 引言

隨著跟蹤、預(yù)報(bào)以及攔截等反導(dǎo)技術(shù)的快速發(fā)展，彈道形式相對固定的傳統(tǒng)彈道導(dǎo)彈已經(jīng)難以突破現(xiàn)有的導(dǎo)彈防御體系，近年來，高超聲速滑翔飛行器（Hypersonic gliding vehicle，HGV）憑借其大升阻比的外形，可在臨近空間環(huán)境內(nèi)以超過5 的馬赫數(shù)實(shí)現(xiàn)遠(yuǎn)距離滑翔及大范圍機(jī)動(dòng)，成為工程研究的重點(diǎn)，現(xiàn)已開展多次HGV 飛行試驗(yàn)［1-2］。HGV 對現(xiàn)有的導(dǎo)彈防御體系造成了極大的壓力，對HGV 的跟蹤難點(diǎn)在于其彈道形式多樣，在滑翔段可以以平衡滑翔或跳躍滑翔方式飛行［3］，因此難以對其進(jìn)行準(zhǔn)確建模。此外，HGV 在臨近空間內(nèi)高速飛行產(chǎn)生時(shí)變的等離子鞘套，對地基雷達(dá)的探測過程也產(chǎn)生了較大影響，使得探測噪聲呈現(xiàn)出非高斯分布的特性［2］。因此，在非高斯噪聲條件下研究臨近空間HGV 軌跡跟蹤方法就顯得格外重要。

為解決HGV 的軌跡跟蹤問題，一個(gè)重要的方面就是如何建立符合HGV 運(yùn)動(dòng)特點(diǎn)的目標(biāo)運(yùn)動(dòng)模型，目前的研究方向主要分為兩個(gè)方向，分別是基于動(dòng)力學(xué)的跟蹤方法以及基于運(yùn)動(dòng)學(xué)的跟蹤方法。基于動(dòng)力學(xué)的跟蹤方法其主要思想是對未知的氣動(dòng)參數(shù)進(jìn)行辨識，即將氣動(dòng)參數(shù)增廣到狀態(tài)方程中，基于氣動(dòng)力模型對HGV 進(jìn)行跟蹤［4-5］。魏世君等［6］根據(jù)天基紅外衛(wèi)星的測量信息設(shè)計(jì)了機(jī)動(dòng)觀測器，對未知的氣動(dòng)加速度進(jìn)行估計(jì)，并分析了觀測器的誤差界。文獻(xiàn)［7］將氣動(dòng)加速度建模為Singer模型，同時(shí)利用正交原理判斷建模誤差，自適應(yīng)修改模型的機(jī)動(dòng)頻率。但此類方法的跟蹤精度嚴(yán)重依賴于氣動(dòng)參數(shù)的辨識結(jié)果，并往往需要目標(biāo)更多的先驗(yàn)信息，當(dāng)氣動(dòng)力變化較快或先驗(yàn)信息與實(shí)際相差較大時(shí)難以保證估計(jì)精度，因此在工程中通常難以直接應(yīng)用［8］。在運(yùn)動(dòng)學(xué)方面，常用的假設(shè)加速度為常值的常加速度（CA）模型、Singer模型、對加速度進(jìn)行自適應(yīng)的“當(dāng)前”統(tǒng)計(jì)模型（Current statistical model，CSM）以及對加加速度進(jìn)行估計(jì)的Jerk 模型均難以準(zhǔn)確描述HGV 的運(yùn)動(dòng)特點(diǎn)。針對HGV 的周期性的跳躍滑翔式運(yùn)動(dòng)，王國宏等［9］將加速度建模為具有正弦波自相關(guān)的零均值過程，進(jìn)而提出了SW（Sine wave）模型，該模型合理的描述了HGV 運(yùn)動(dòng)的周期性特征。文獻(xiàn)［10-11］將加速度的均值假設(shè)為正弦自相關(guān)隨機(jī)過程，并利用一階或二階馬爾可夫過程對機(jī)動(dòng)進(jìn)行補(bǔ)償，進(jìn)而提出了自適應(yīng)非零均值衰減振蕩（ANMDO）模型，該模型同時(shí)響應(yīng)了HGV 跳躍滑翔階段的周期性和衰減性的特點(diǎn)，是目前描述HGV 跳躍滑翔運(yùn)動(dòng)較為準(zhǔn)確的運(yùn)動(dòng)學(xué)模型。

由于單模型方法難以匹配HGV 的多種運(yùn)動(dòng)情況，目前的研究重點(diǎn)在如何建立合適的目標(biāo)運(yùn)動(dòng)模型集以及研究基于多模型濾波結(jié)構(gòu)的跟蹤方法上［12］，文獻(xiàn)［13］基于Singer、CA 和CT 模型建立了交互式多模型（Interacting multiple model，IMM）跟蹤濾波器以實(shí)現(xiàn)對HGV 的跟蹤。文獻(xiàn)［14］利用多個(gè)機(jī)動(dòng)角速率不同的Sine 模型對HGV 的跳躍滑翔運(yùn)動(dòng)進(jìn)行匹配。文獻(xiàn)［15］則在模型集中同時(shí)加入了動(dòng)力學(xué)模型和運(yùn)動(dòng)學(xué)模型，利用多模型方法融合了動(dòng)力學(xué)模型建模匹配程度高和運(yùn)動(dòng)學(xué)適應(yīng)性強(qiáng)的優(yōu)點(diǎn)?？紤]到固定結(jié)構(gòu)多模型算法模型集固定且各模型之間存在競爭導(dǎo)致精度下降的問題，文獻(xiàn)［16］利用自適應(yīng)網(wǎng)格變結(jié)構(gòu)多模型對結(jié)構(gòu)和模型參數(shù)進(jìn)行自適應(yīng)調(diào)整。文獻(xiàn)［17］在變結(jié)構(gòu)多模型的基礎(chǔ)上利用有向圖切換的方法對SW 模型的周期機(jī)動(dòng)角速度進(jìn)行自適應(yīng)切換，以實(shí)現(xiàn)對HGV 滑翔段的有效跟蹤。

有關(guān)HGV 軌跡跟蹤，另一個(gè)重要的研究方面是非線性濾波方法。在HGV 的跟蹤系統(tǒng)中狀態(tài)方程和量測方程均具有較強(qiáng)的非線性，因此在工程應(yīng)用中常常使用擴(kuò)展卡爾曼濾波（EKF）對其進(jìn)行跟蹤，但由于EKF 存在較大的截?cái)嗾`差，導(dǎo)致其跟蹤精度有限。因此學(xué)者們逐漸采用基于確定性采樣的非線性濾波方法進(jìn)行跟蹤濾波器的設(shè)計(jì)，包括基于徑向容積規(guī)則的容積卡爾曼濾波（CKF）和基于無跡規(guī)則的無跡卡爾曼濾波（UKF）［18-19］。但是，上述方法均是在噪聲符合高斯分布的假設(shè)下進(jìn)行推導(dǎo)的，當(dāng)探測噪聲呈現(xiàn)非高斯特性時(shí)，這些算法的精度會(huì)嚴(yán)重下降甚至發(fā)散。在處理非高斯噪聲問題方面，最常用的魯棒濾波方法為基于廣義極大似然的Huber濾波，該方法通過混合l1和l2范數(shù)的組合以獲得對量測異常值的有效抑制［20］，其對閃爍噪聲具有良好的漸進(jìn)最優(yōu)魯棒性能。文獻(xiàn)［21］將Huber濾波與五階CKF 進(jìn)行了融合，提出了HCHF 用于助推段軌跡的跟蹤。文獻(xiàn)［22］在Huber濾波處理量測噪聲模型存在誤差問題的基礎(chǔ)上，將UKF 與強(qiáng)跟蹤方法相結(jié)合，用于跟蹤存在傾側(cè)角反轉(zhuǎn)的HGV 目標(biāo)。近年來，也有學(xué)者利用學(xué)生t 分布［23］或高斯混合模型［24］對重尾分布的噪聲進(jìn)行建模，并利用變分貝葉斯方法對噪聲模型的參數(shù)進(jìn)行實(shí)時(shí)估計(jì)以適應(yīng)先驗(yàn)信息的不準(zhǔn)確性。雖然這種方法可以得到較高的估計(jì)精度，但需要利用不動(dòng)點(diǎn)迭代的方法進(jìn)行多次運(yùn)算［25］，往往伴隨著更大的計(jì)算量。

綜上所述，對臨近空間HGV 軌跡跟蹤的難點(diǎn)在于目標(biāo)機(jī)動(dòng)形式復(fù)雜多變且機(jī)動(dòng)時(shí)機(jī)不確定造成的目標(biāo)運(yùn)動(dòng)模型建模困難，以及非高斯噪聲條件下非線性濾波跟蹤方法精度難以保證。因此本文首先設(shè)計(jì)了基于有向圖切換的多模型濾波結(jié)構(gòu)，在傳統(tǒng)交互式多模型算法的基礎(chǔ)上，利用有向圖切換的方法實(shí)現(xiàn)模型集結(jié)構(gòu)及參數(shù)的自適應(yīng)調(diào)整；其次，將Huber 濾波方法與CKF 相結(jié)合，利用Huber 濾波實(shí)現(xiàn)對非高斯噪聲的魯棒性并應(yīng)用CKF 處理跟蹤系統(tǒng)的強(qiáng)非線性；然后，采用不同衰減頻率的ANMDO 模型對HGV 跳躍滑翔階段進(jìn)行近似，同時(shí)利用不同機(jī)動(dòng)頻率的“當(dāng)前”統(tǒng)計(jì)模型描述平衡滑翔運(yùn)動(dòng)并建立了目標(biāo)運(yùn)動(dòng)模型集，之后基于地基雷達(dá)的量測模型完成HGV 跟蹤方法的設(shè)計(jì)；最后利用數(shù)學(xué)仿真的方法對所提算法進(jìn)行了驗(yàn)證。

1 基于有向圖切換的變結(jié)構(gòu)多模型結(jié)構(gòu)

HGV 在滑翔段具有跳躍滑翔和平衡滑翔兩種常用的飛行模式，傳統(tǒng)的單模型方法無法適應(yīng)這種機(jī)動(dòng)樣式復(fù)雜多變的情況，因此需要利用多模型方法以實(shí)現(xiàn)有效跟蹤。最常見的多模型方法即交互式多模型方法，以其結(jié)構(gòu)簡單、計(jì)算效率高等優(yōu)點(diǎn)在工程中得到了廣泛的應(yīng)用。但對于HGV 來說，覆蓋所有可能的機(jī)動(dòng)樣式會(huì)導(dǎo)致模型集過于龐大，不僅導(dǎo)致計(jì)算成本急劇增加，而且由于模型集之間的模型競爭又會(huì)導(dǎo)致跟蹤精度下降［26］，因此需要采用變結(jié)構(gòu)多模型對其進(jìn)行改進(jìn)。本文基于有向圖的思想，在IMM 的基礎(chǔ)上利用有向圖切換的方式對實(shí)時(shí)模型集進(jìn)行自適應(yīng)調(diào)整，在保證對HGV 機(jī)動(dòng)樣式覆蓋性的同時(shí)縮短了算法運(yùn)行時(shí)間。

1.1 交互式多模型結(jié)構(gòu)

考慮如下的非線性多模型狀態(tài)空間系統(tǒng)：

式中：k為離散的時(shí)間序列；xk∈Rn為n維的狀態(tài)向量；yk∈Rm為m維的量測向量；f(·)和h(·)分別代表非線性的狀態(tài)方程和量測方程；wk和vk分別代表過程噪聲和量測噪聲；mk∈{1，…，M}為模型標(biāo)識；(·)代表此時(shí)第mk個(gè)狀態(tài)模型正在生效；mk符合馬爾可夫過程，即：

式中：pij為模型i到模型j的轉(zhuǎn)移概率，可由經(jīng)驗(yàn)確定或由自適應(yīng)算法進(jìn)行調(diào)整。

作為多模型算法中最經(jīng)典的算法，IMM 通過設(shè)計(jì)一個(gè)覆蓋目標(biāo)所有可能運(yùn)動(dòng)情況的固定模型集，每個(gè)模型獨(dú)立運(yùn)行在各自的子濾波器中，并且各子模型在輸入時(shí)進(jìn)行數(shù)據(jù)的相互作用完成交互，當(dāng)子模型得到各自的估計(jì)結(jié)果后，根據(jù)似然函數(shù)計(jì)算相應(yīng)權(quán)值并將所有子濾波器的狀態(tài)估計(jì)結(jié)果進(jìn)行加權(quán)融合，作為多模型的濾波結(jié)果實(shí)時(shí)輸出，其具體步驟如下。

交互式多模型算法首先根據(jù)上一時(shí)刻各子濾波器的估計(jì)結(jié)果進(jìn)行交互計(jì)算，得到當(dāng)前時(shí)刻子濾波器的輸入，包括各子模型的先驗(yàn)?zāi)Ｐ透怕?、狀態(tài)和誤差協(xié)方差。對于第j個(gè)子模型的先驗(yàn)?zāi)Ｐ透怕?/p>

式中：uk-1[i]=p(mk-1=i|y1：k-1)為上一時(shí)刻第i個(gè)模型的后驗(yàn)?zāi)Ｐ透怕省?/p>

對于第j個(gè)子模型的狀態(tài)先驗(yàn)概率p(xk-1|mk=j，y1：k-1)，有：

假設(shè)各子模型的后驗(yàn)概率均符合高斯分布，則第j個(gè)子模型輸入的初始狀態(tài)均值和協(xié)方差分別為

隨后各子濾波器利用當(dāng)前時(shí)刻的量測信息獨(dú)立并行進(jìn)行濾波更新，獲得各自的狀態(tài)估計(jì)、誤差協(xié)方差、新息以及新息協(xié)方差矩陣S(j)。計(jì)算似然函數(shù)對模型概率uk[j]進(jìn)行更新：

式中：N(μ；0，σ)表示高斯分布。

最后IMM 通過融合各子濾波器的估計(jì)結(jié)果獲得最終的狀態(tài)估計(jì)，假設(shè)第j個(gè)子模型的估計(jì)結(jié)果為

1.2 有向圖切換方法

變結(jié)構(gòu)多模型在IMM 的基礎(chǔ)上，提出了模型總集和實(shí)時(shí)模型集的概念。假設(shè)模型總集為ΘM，實(shí)時(shí)模型集為Θk∈ΘM，變結(jié)構(gòu)多模型是通過某種的模型集自適應(yīng)策略，在總模型ΘM中找到可描述目標(biāo)當(dāng)前運(yùn)動(dòng)狀態(tài)的最小模型集Θk，并將Θk作為實(shí)時(shí)模型集運(yùn)行IMM 算法得到機(jī)動(dòng)目標(biāo)的狀態(tài)估計(jì)。基于有向圖切換的變結(jié)構(gòu)多模型方法為通過可能模型集的設(shè)計(jì)，完成對總模型集ΘM的完全覆蓋，并使得相鄰的兩個(gè)可能模型集中存在重復(fù)的模型，同時(shí)保證各可能模型集不完全相同，通過判斷邊緣模型的實(shí)時(shí)模型概率是否超過閾值，從而確定是否對實(shí)時(shí)模型集進(jìn)行切換。

如圖1 所示，首先預(yù)制數(shù)個(gè)具有相同模型數(shù)量的可能模型集Θ={Θ1，Θ2，…，Θm}，可能模型集中包含位于中心位置的中心模型和邊緣模型，可見各子模型集的中心模型同時(shí)為相鄰模型集中的邊緣模型，如果當(dāng)前模型集Θk中的某一邊緣模型mj的后驗(yàn)?zāi)Ｐ透怕蚀笥谒O(shè)定的閾值λ，則實(shí)時(shí)模型集向以模型mj為中心模型的可能模型集移動(dòng)，從而完成模型集的切換。

2 基于Huber更新的魯棒容積濾波方法

當(dāng)量測噪聲不再符合高斯分布時(shí)，若仍用EKF等濾波方法會(huì)導(dǎo)致跟蹤濾波器的估計(jì)精度下降甚至發(fā)散，而Huber 濾波作為一種混合H1/H2的濾波方法，對非高斯噪聲具有良好的魯棒性?？紤]到基于確定性采樣的CKF 具有比EKF 更高的近似精度和良好的數(shù)值穩(wěn)定性，本文將Huber 濾波與CKF 融合進(jìn)行濾波器的設(shè)計(jì)。

2.1 廣義極大似然估計(jì)

首先考慮如下的線性回歸問題：

式中：H為量測矩陣；v為方差陣為R的量測噪聲。利用量測噪聲方差陣進(jìn)行如下變換：

Huber 利用廣義極大似然估計(jì)來求解（15）所示的線性回歸問題，即計(jì)算如下指標(biāo)函數(shù)的極小值：

對式（16）求導(dǎo)并令其等于0，可得：

通過定義權(quán)重函數(shù)ψ(δi)=及權(quán)重矩陣Ψ=diag(ψ(δi))，則其解的隱含數(shù)方程可表示為

Huber 通過選用如下混合l1和l2的ρ函數(shù)以獲得對非高斯噪聲的魯棒性：

式中：γ為可調(diào)參數(shù)。Huber 濾波通過選取迭代初值：

并按下式迭代至收斂獲取估計(jì)的結(jié)果。

2.2 容積Huber-based濾波

容積卡爾曼濾波作為一種確定性采樣的非線性濾波方法，通過徑向容積規(guī)則進(jìn)行容積點(diǎn)的選取，并利用非線性方程對容積點(diǎn)進(jìn)行傳播，以此獲得非線性系統(tǒng)的后驗(yàn)估計(jì)，其不僅具有比擴(kuò)展卡爾曼更高的近似精度，而且相對于無跡卡爾曼來說其數(shù)值穩(wěn)定性更強(qiáng)?？紤]HGV 跟蹤系統(tǒng)為非線性系統(tǒng)，因此使用容積卡爾曼濾波與Huber算法相結(jié)合，處理非線性系統(tǒng)中的非高斯噪聲問題，考慮如下狀態(tài)方程（即mk）已確定的非線性子模型系統(tǒng)，將(xk)簡寫為f(xk)，則式（1）和式（2）可記為

則基于Huber的容積濾波算法步驟如下：

1) 時(shí)間更新

式中：Sk-1為協(xié)方差矩陣Pk-1|k-1的Cholesky 分解。ξj為如下點(diǎn)集的第j列：

將容積點(diǎn)通過非線性狀態(tài)方程進(jìn)行傳播，有：

2) 量測更新

生成用于量測更新的容積點(diǎn)：

定義狀態(tài)預(yù)測誤差為ek=xk-，并將量測方程線性化：

可構(gòu)建如下線性回歸問題：

取迭代初值為：

則迭代求解的狀態(tài)估計(jì)值以及誤差協(xié)方差為：

3 高超聲速滑翔目標(biāo)高精度跟蹤方法

本節(jié)將建立HGV 運(yùn)動(dòng)模型集以及地基雷達(dá)的量測模型，隨后基于有向圖變結(jié)構(gòu)多模型結(jié)構(gòu)并利用容積Huber-based 濾波方法建立HGV 高精度跟蹤濾波器。

3.1 高超聲速滑翔目標(biāo)運(yùn)動(dòng)模型集

HGV 在滑翔段常見的飛行方式包括平衡滑翔和跳躍滑翔，由于這兩種飛行方式的彈道特性有明顯不同，因此需要分別建立目標(biāo)運(yùn)動(dòng)模型。

對于平衡滑翔運(yùn)動(dòng)階段，其彈道變化相對平緩，可用“當(dāng)前”統(tǒng)計(jì)模型對其進(jìn)行描述，該模型假設(shè)目標(biāo)當(dāng)前的加速度在上一時(shí)刻的加速度鄰域范圍內(nèi)。取目標(biāo)的位置、速度以及加速度作為狀態(tài)量，即：

則一維的“當(dāng)前”統(tǒng)計(jì)模型離散化后的狀態(tài)方程為

式中：α為機(jī)動(dòng)頻率參數(shù)，T為采樣步長。噪聲wk的方差為：

σ2代表“當(dāng)前”加速度方差，滿足：

其中，amax和a-max分別為最大正負(fù)加速度。

對于跳躍滑翔運(yùn)動(dòng)階段，考慮高超聲速滑翔目標(biāo)在該階段的加速度存在周期振蕩的特點(diǎn)，可用正弦波模型對其進(jìn)行描述。同時(shí)，由于氣動(dòng)力的存在，目標(biāo)的加速度均值不為零且隨時(shí)間變化，并且在阻力作用下飛行器的能量逐漸耗散，導(dǎo)致周期性逐漸減弱，因此，本文選取文獻(xiàn)［10］所述的ANMDO模型對跳躍滑翔階段進(jìn)行建模。取目標(biāo)的位置、速度、加速度和加加速度作為狀態(tài)量，即：

則一維的ANMDO模型離散化后的狀態(tài)方程為

式中：Fk為狀態(tài)轉(zhuǎn)移矩陣，滿足

其中，

β代表權(quán)重系數(shù)；μ代表衰減系數(shù)；U為加速度補(bǔ)償項(xiàng)，滿足

式中：ω為振蕩頻率。當(dāng)μ為零時(shí)，模型中的加速度僅表現(xiàn)周期性，無衰減；當(dāng)ω為零時(shí)，加速度僅表現(xiàn)衰減性，無周期性。因此通過選擇合適的衰減系數(shù)和振蕩頻率參數(shù)，可以保證角加速度同時(shí)具有周期性和衰減性。

機(jī)動(dòng)頻率α是“當(dāng)前”統(tǒng)計(jì)模型的主要參數(shù)，為保證對平衡滑翔段的覆蓋性，完備模型集中包含兩個(gè)機(jī)動(dòng)頻率不同的“當(dāng)前”統(tǒng)計(jì)模型，其模型參數(shù)分別為α1<α2。在ANMDO 模型中，振蕩頻率ω是其主要的模型參數(shù)，考慮跳躍滑翔階段目標(biāo)可能具有的不同振蕩周期，完備模型集中包含3 個(gè)振蕩頻率不同的ANMDO 模型，對應(yīng)的模型參數(shù)在目標(biāo)可能的振蕩頻率變化范圍內(nèi)選取，即[ωmin<ω1<ω2<ω3<ωmax]，則完備模型集中模型組合的有向圖網(wǎng)絡(luò)為{α2，α1，ω1，ω2，ω3}，取實(shí)時(shí)模型集數(shù)量為3，則共完備模型集可分為3 個(gè)可能模型集，分別為M1={α2，α1，ω1}，M2={α1，ω1，ω2}，M3={ω1，ω2，ω3}，3 個(gè)可能模型集間的切換準(zhǔn)則為：

1) 假設(shè)實(shí)時(shí)模型集為M1，若ω1對應(yīng)的ANMDO模型在多模型算法中的模型概率最大且大于設(shè)定閾值λ，則下一時(shí)刻的實(shí)時(shí)模型集切換為可能模型集中的M2；否則下一時(shí)刻的實(shí)時(shí)模型集仍保持為M1；

2) 假設(shè)實(shí)時(shí)模型集為M2，若α1對應(yīng)的“當(dāng)前”統(tǒng)計(jì)模型在多模型算法中的模型概率最大且大于設(shè)定的閾值λ，則下一時(shí)刻的實(shí)時(shí)模型集切換為M1；若ω2對應(yīng)的ANMDO模型的模型概率最大且大于設(shè)定閾值λ，則下一時(shí)刻的實(shí)時(shí)模型集切換為M3；否則，下一時(shí)刻的模型集保持不變；

3) 假設(shè)實(shí)時(shí)模型集為M3，若ω1對應(yīng)的ANMDO的模型概率最大且大于閾值λ，則下一時(shí)刻的實(shí)時(shí)模型集切換為M2；否則，下一時(shí)刻的模型集不變。

3.2 地基雷達(dá)量測模型

地基雷達(dá)通過主動(dòng)發(fā)射電磁波并接收回波對目標(biāo)進(jìn)行探測，具有探測精度高，作用距離遠(yuǎn)等優(yōu)勢，常用來對HGV 的滑翔段進(jìn)行探測。地基雷達(dá)的量測信息包括目標(biāo)相對于雷達(dá)站的斜距r、高低角γ及方位角η，即：

對應(yīng)的量測方程為：

式中：Δx，Δy和Δz為HGV在雷達(dá)探測坐標(biāo)系中的位置分量。

3.3 跟蹤濾波器算法流程

將建立的HGV 目標(biāo)運(yùn)動(dòng)模型作為狀態(tài)方程，地基雷達(dá)量測模型作為量測方程，則基于有向圖多模型結(jié)構(gòu)和容積Huber-based 濾波算法的HGV 高精度跟蹤算法流程可總結(jié)如下：

1) 數(shù)據(jù)輸入交互：將上一時(shí)刻各子模型的按照式（4）至式（7）進(jìn)行融合得到各子模型的模型概率和當(dāng)前時(shí)刻各子模型的狀態(tài)初值。

2) 狀態(tài)濾波：各子模型按照容積Huber-based方法進(jìn)行濾波，按照式（24）至式（42）得到各子模型的狀態(tài)估計(jì)，估計(jì)誤差協(xié)方差，以及用于模型更新的新息及其協(xié)方差。

3) 模型概率更新：利用獲得的新息及其協(xié)方差，按照式（9）至式（10）更新當(dāng)前時(shí)刻的各子模型的后驗(yàn)?zāi)Ｐ透怕省?/p>

4) 融合輸出：利用更新后的模型概率和各子模型狀態(tài)估計(jì)及協(xié)方差，按式（12）至式（13）進(jìn)行加權(quán)融合輸出。

5) 模型集切換：根據(jù)后驗(yàn)?zāi)Ｐ透怕世糜邢驁D多模型的切換準(zhǔn)則對實(shí)時(shí)模型集進(jìn)行切換，確定新的當(dāng)前實(shí)時(shí)模型集。

4 仿真校驗(yàn)

4.1 仿真場景

為驗(yàn)證本文所提的CHF-DSVSMM 算法在HGV跟蹤中的性能，現(xiàn)以CAV-H 飛行器為例構(gòu)建如下的仿真場景。設(shè)置飛行中最大過載3 倍重力加速度、最大熱流密度1 000 kW/m2、最大動(dòng)壓300 kPa，仿真算例初值條件如表1所示。

表1 飛行軌跡優(yōu)化算例的初值Table 1 Initial value of flight trajectory optimization example

采用高斯偽譜方法進(jìn)行軌跡優(yōu)化［27］，HGV 飛行軌跡如圖2 所示，目標(biāo)飛行速度變化曲線如圖3 所示。地基雷達(dá)部署位置經(jīng)度、緯度、高度為[222.5°，40°，0 m ]，假設(shè)雷達(dá)最大可測量距離為750 km，將可探測弧段同時(shí)繪制于圖2 中，采樣頻率10 Hz，測距精度σr=100 m，測角精度σρ=ση=0.1°。

圖2 HGV飛行軌跡及跟蹤場景Fig.2 HGV flight trajectory and simulation scene

圖3 速度變化曲線Fig.3 Target velocity

4.2 仿真一：高斯噪聲條件下軌跡跟蹤

在所建立的仿真條件下，利用CKF 并結(jié)合CSM單模型、文獻(xiàn)［10］的ANMDO 模型、文獻(xiàn)［13］的交互式多模型、文獻(xiàn)［7］所提的自適應(yīng)動(dòng)力學(xué)模型（Adaptive dynamical model）以及本文所提的CHFDSVSMM 算法進(jìn)行HGV 的軌跡跟蹤，單模型ANMDO 模型中參數(shù)設(shè)置為權(quán)重系數(shù)β=0.7，衰減頻率μ=11 000，周期頻率ω=0.05；IMM 算法由Singer、CA、CT 三個(gè)模型構(gòu)成，Singer 模型的機(jī)動(dòng)頻率為α=140；自適應(yīng)動(dòng)力學(xué)模型中具體參數(shù)見文獻(xiàn)［7］；CHF-DSVSMM 的模型集為“當(dāng)前”統(tǒng)計(jì)模型α1=160，α2=140 和ANMDO 模型ω1=0.07，ω2=0.05，ω3=0.03，切換閾值為λ=0.8。首先在高斯噪聲條件下對驗(yàn)證所提算法對HGV 的跟蹤精度，進(jìn)行500 次蒙特卡洛仿真試驗(yàn)，仿真結(jié)果如圖4至圖6所示。

圖4 位置估計(jì)誤差Fig.4 Position estimation error of different filters

圖5 速度估計(jì)誤差Fig.5 Velocity estimation error of different filters

圖6 加速度估計(jì)誤差Fig.6 Acceleration estimation error of different filters

根據(jù)仿真結(jié)果可知，由于CSM 單模型和文獻(xiàn)［13］中IMM 所包含的模型均難以適應(yīng)HGV 在跳躍滑翔時(shí)的強(qiáng)機(jī)動(dòng)特性，因此跟蹤精度最差。ANMDO單模型同時(shí)具有周期性和衰減性，因此可對跳躍滑翔運(yùn)動(dòng)實(shí)現(xiàn)較為準(zhǔn)確的跟蹤，但單模型方法下其模型參數(shù)固定，不能實(shí)時(shí)準(zhǔn)確地適應(yīng)HGV 真實(shí)的機(jī)動(dòng)情況，導(dǎo)致其跟蹤精度難以長時(shí)間保證。自適應(yīng)動(dòng)力學(xué)模型需要對加速度進(jìn)行自適應(yīng)估計(jì)，因此需要得到一定數(shù)據(jù)量后才能對加速度變化進(jìn)行響應(yīng)。而本文通過建立包含多種參數(shù)的CSM 和ANMDO的模型集，實(shí)現(xiàn)了對HGV 跳躍滑翔運(yùn)動(dòng)較為完整地覆蓋，并通過有向圖切換的方法在模型集中自適應(yīng)地快速調(diào)整實(shí)時(shí)模型集，在保證了計(jì)算量的同時(shí)實(shí)現(xiàn)了對HGV 跳躍滑翔階段持續(xù)穩(wěn)定地高精度跟蹤。

4.3 仿真二：非高斯噪聲條件下軌跡跟蹤

為進(jìn)一步驗(yàn)證所提算法對非高斯噪聲的魯棒性，假設(shè)地基雷達(dá)受到等離子鞘套的干擾，表現(xiàn)出的非高斯探測噪聲滿足閃爍噪聲形式：

進(jìn)行500 次獨(dú)立的蒙特卡洛打靶試驗(yàn)，仿真結(jié)果如圖7至圖9所示。

圖7 位置估計(jì)誤差Fig.7 Position estimation error of different filters

圖8 速度估計(jì)誤差Fig.8 Velocity estimation error of different filters

圖9 加速度估計(jì)誤差Fig.9 Acceleration estimation error of different filters

定義平均均方根誤差（ARMSE）作為性能指標(biāo)，即：

式中：N為蒙特卡洛打靶次數(shù)；M為仿真步長；Xk和分別為狀態(tài)真值與估計(jì)值。

將兩種仿真情況下不同算法的ARMSE 總結(jié)如表2。

表2 不同算法估計(jì)結(jié)果的ARMSETable 2 ARMSE of different filters

根據(jù)仿真結(jié)果可知，在非高斯噪聲條件下，單模型和IMM 中采用的CKF 估計(jì)精度明顯下降，而由于文獻(xiàn)［7］和本文均采用了HCF 處理量測過程中出現(xiàn)的非高斯噪聲，其跟蹤精度在非高斯噪聲條件下精度相當(dāng)。因此所提HCF-DSVSMM 方法不僅可以在高斯噪聲條件下取得穩(wěn)定的跟蹤效果，而且對非高斯噪聲具有良好的魯棒性，非高斯噪聲對所提算法影響最小。

5 結(jié)論

本文研究了非高斯噪聲條件下的臨近空間HGV 軌跡跟蹤方法，通過將容積Huber-based 濾波方法與有向圖變結(jié)構(gòu)多模型方法相結(jié)合，提出了CHF-DSVSMM 算法，并設(shè)計(jì)了HGV 目標(biāo)運(yùn)動(dòng)模型集。與其他方法相比，本文所建立的模型集可以更為準(zhǔn)確地描述HGV 跳躍滑翔階段的運(yùn)動(dòng)，并且利用有向圖切換的方法在保證計(jì)算量的同時(shí)實(shí)現(xiàn)了對HGV 的穩(wěn)定高精度跟蹤，同時(shí)，算法中的HCF 還提供了對非高斯噪聲的魯棒性。仿真結(jié)果驗(yàn)證了該算法的有效性，仿真結(jié)果表明，所提算法可以實(shí)現(xiàn)對具有多種運(yùn)動(dòng)情況的HGV 持續(xù)穩(wěn)定的高精度跟蹤，算法通過模型集切換的方式自適應(yīng)的改變目標(biāo)運(yùn)動(dòng)模型集，達(dá)到了對HGV 實(shí)時(shí)運(yùn)動(dòng)特性的準(zhǔn)確描述，且當(dāng)存在非高斯量測噪聲時(shí)，算法仍保持了良好的魯棒性和跟蹤精度。