一種高超聲速飛行器在線反饋濾波算法*

高長(zhǎng)生，王越欣，荊武興，胡玉東

(哈爾濱工業(yè)大學(xué) 航天學(xué)院 自主空間系統(tǒng)實(shí)驗(yàn)室，黑龍江 哈爾濱 150001)

0 引言

臨近空間是指距離地面高度為20～100 km的空域，該空域處在航空和航天之間，比衛(wèi)星飛行高度低，比飛機(jī)飛行高度高，探測(cè)信息相對(duì)空缺，掌握該空域?qū)θ〉弥瓶諜?quán)有很大的作用，臨近空間特殊的戰(zhàn)略地位對(duì)未來(lái)戰(zhàn)爭(zhēng)體系的影響不可估量[1-3],近年來(lái)得到了各國(guó)軍方的高度重視。以美國(guó)為代表，如HTV2，CAV等高超聲速飛行器，其高速、遠(yuǎn)程、精確的特點(diǎn)可以實(shí)現(xiàn)全球戰(zhàn)略打擊目的，給其他國(guó)家的國(guó)防帶來(lái)了嚴(yán)肅的考驗(yàn)[4-5]。在此背景下，針對(duì)該類飛行器的攔截研究已迫在眉睫。

在解決臨近空間目標(biāo)跟蹤問(wèn)題時(shí)，一個(gè)重要的方面就是構(gòu)建高超聲速飛行器的運(yùn)動(dòng)模型，經(jīng)過(guò)近幾十年各國(guó)科學(xué)家和學(xué)者的深入研究，其相關(guān)理論和模型日趨完善。相關(guān)的運(yùn)動(dòng)學(xué)模型主要有：勻速模型(CV)，狀態(tài)量只有位置和速度2方面，適合跟蹤目標(biāo)作簡(jiǎn)單的勻速直線運(yùn)動(dòng)，無(wú)任何機(jī)動(dòng)情況下的運(yùn)動(dòng)形式[6]。勻加速模型(CA)，相較于勻速模型，把目標(biāo)的加速度加入到狀態(tài)量中，適合跟蹤目標(biāo)作勻加速運(yùn)動(dòng)對(duì)應(yīng)的跟蹤模型[7]。圓周轉(zhuǎn)彎模型(CT),指物體的速度和加速度大小均不變化，只有方向改變，主要用于跟蹤目標(biāo)轉(zhuǎn)彎，常用于描述二維運(yùn)動(dòng)狀態(tài)[8]。上述模型均是將加速度假設(shè)為高斯白噪聲系統(tǒng)，這在實(shí)際應(yīng)用中會(huì)出現(xiàn)很大的誤差，此外一階時(shí)間相關(guān)模型(Singer)，將目標(biāo)的加速度作為具有零均值指數(shù)衰減自相關(guān)函數(shù)的隨機(jī)過(guò)程，且具有指數(shù)自相關(guān)的特性，彌補(bǔ)了這一假設(shè)誤差[9]。在此基礎(chǔ)上，國(guó)內(nèi)的周宏仁教授對(duì)Singer模型進(jìn)行了改進(jìn)，認(rèn)為下一時(shí)刻的加速度不能隨意取值，只能在當(dāng)前加速度的領(lǐng)域內(nèi)取值，目標(biāo)加速度的統(tǒng)計(jì)特性用修正的瑞利分布來(lái)表示。在目標(biāo)運(yùn)動(dòng)過(guò)程中實(shí)現(xiàn)了加速度均值自適應(yīng)變化，并且加速度概率分布與均值相關(guān)，方差則由加速度均值決定。目標(biāo)的加速度只能在前一時(shí)刻加速度的鄰域范圍內(nèi)變化，建模為非零均值的自相關(guān)指數(shù)衰減過(guò)程[10]。

在高超聲速飛行器跟蹤系統(tǒng)中，工程上應(yīng)用最為廣泛的是卡爾曼濾波。應(yīng)用卡爾曼濾波算法完成系統(tǒng)狀態(tài)估計(jì)，需要依賴構(gòu)建的系統(tǒng)模型和實(shí)時(shí)的測(cè)量數(shù)據(jù)。然而，在跟蹤未知臨近空間目標(biāo)時(shí)，飛行器的運(yùn)動(dòng)規(guī)律及參數(shù)無(wú)法得知，必然會(huì)引起跟蹤模型的不匹配問(wèn)題，同時(shí)目標(biāo)運(yùn)動(dòng)過(guò)程中隨時(shí)可能存在機(jī)動(dòng)模式切換情況[11]。目標(biāo)運(yùn)動(dòng)的不確定性嚴(yán)重影響了算法的精度和穩(wěn)定性。目標(biāo)機(jī)動(dòng)運(yùn)動(dòng)實(shí)質(zhì)上是一種非線性現(xiàn)象，可以利用神經(jīng)網(wǎng)絡(luò)的強(qiáng)非線性映射能力對(duì)目標(biāo)運(yùn)動(dòng)的不確定性建模，對(duì)目標(biāo)機(jī)動(dòng)在線辨識(shí)，實(shí)時(shí)修正濾波估計(jì)值[12-14]。BP神經(jīng)網(wǎng)絡(luò)是一種按照誤差逆向傳播算法訓(xùn)練的多層前饋神經(jīng)網(wǎng)絡(luò)，是應(yīng)用最廣泛的神經(jīng)網(wǎng)絡(luò)，可以利用BP神經(jīng)網(wǎng)絡(luò)對(duì)濾波算法的濾波誤差進(jìn)行預(yù)測(cè)，并補(bǔ)償給濾波估計(jì)值，從而提高濾波精度[15]。

綜上所述，本文針對(duì)高超聲速飛行器難以進(jìn)行高精度跟蹤這一問(wèn)題，提出了一種高超聲速飛行器在線反饋濾波算法。在當(dāng)前統(tǒng)計(jì)模型的基礎(chǔ)上，利用BP神經(jīng)網(wǎng)絡(luò)與卡爾曼濾波相結(jié)合進(jìn)行濾波器設(shè)計(jì)，實(shí)現(xiàn)對(duì)高超聲速飛行器高精度跟蹤。最后采用數(shù)學(xué)仿真手段對(duì)所提出的算法進(jìn)行了驗(yàn)證，并進(jìn)行相關(guān)分析。

1 運(yùn)動(dòng)學(xué)跟蹤模型與量測(cè)模型

1.1 運(yùn)動(dòng)學(xué)跟蹤模型

目標(biāo)跟蹤模型是跟蹤算法的基礎(chǔ)，是系統(tǒng)的狀態(tài)方程部分，即當(dāng)前時(shí)刻目標(biāo)的狀態(tài)量由前一時(shí)刻狀態(tài)量的表示方法，目標(biāo)的運(yùn)動(dòng)特性由位置、速度、加速度等狀態(tài)量表征。

(1)

(2)

假設(shè)紅外采樣間隔為T，對(duì)上述方程離散化，則離散時(shí)間狀態(tài)方程為

(3)

式中：X(k)為k時(shí)刻目標(biāo)運(yùn)動(dòng)的狀態(tài)量;F(k)為狀態(tài)轉(zhuǎn)移矩陣;G(k)為輸入控制矩陣；W(k)為噪聲。

狀態(tài)轉(zhuǎn)移矩陣F(k)的表達(dá)式為

(4)

輸入控制矩陣G(k)的表達(dá)式為

(5)

當(dāng)前時(shí)刻的加速度均值一般情況下難以確定，用濾波遞推過(guò)程中的加速度一步預(yù)測(cè)代替，有

(6)

目標(biāo)當(dāng)前加速度方差估計(jì)：

(7)

式中：amax和a-max分別為目標(biāo)加速度的最大值和最小值。

噪聲W(k)的協(xié)方差為

(8)

Q為對(duì)稱矩陣，而且有

(9)

綜上，“當(dāng)前”統(tǒng)計(jì)模型將上一個(gè)時(shí)刻的加速度輸入到當(dāng)前時(shí)刻來(lái)估算加速度方差，進(jìn)而自適應(yīng)地調(diào)整過(guò)程噪聲。

1.2 量測(cè)模型

首先定義探測(cè)坐標(biāo)系OxPyPzP，簡(jiǎn)記為P。OyP軸在探測(cè)系原點(diǎn)水平面內(nèi)，指向北極方向，OzP軸垂直于水平面指向上方，OxP構(gòu)成右手坐標(biāo)系，坐標(biāo)系固連在地球表面。如圖1所示。

設(shè)目標(biāo)在探測(cè)系下的位置矢量：

r=(x,y,z).

(10)

紅外探測(cè)器基點(diǎn)在探測(cè)系下的位置矢量：

Si=(xi,yi,zi).

(11)

目標(biāo)在本體系下的矢量為

(12)

所以有

(13)

如圖2所示。

圖2 定位原理Fig.2 Positioning principle

定位：

(14)

整理為矩陣形式：

(15)

則有

(16)

綜上，使用最小二乘算法定位，則有

M·X=Y?X=(MT·M)-1·MT·Y.

(17)

2 基于BP神經(jīng)網(wǎng)絡(luò)的濾波器設(shè)計(jì)

2.1 BP神經(jīng)網(wǎng)絡(luò)

BP神經(jīng)網(wǎng)絡(luò)由3部分組成，分別為輸入層、隱含層和輸出層，其中輸入層接受外界信息和數(shù)據(jù)輸入，輸出層則輸出網(wǎng)絡(luò)對(duì)輸入信息和數(shù)據(jù)的處理結(jié)果，隱含層作為輸入層和輸出層的媒介，完成信息和數(shù)據(jù)的分析、處理過(guò)程。如圖3所示。

圖3 BP神經(jīng)網(wǎng)絡(luò)基本結(jié)構(gòu)Fig.3 Basic structure of BP neural network

神經(jīng)元作為人工神經(jīng)網(wǎng)絡(luò)中信息、數(shù)據(jù)處理的獨(dú)立單元，主要完成3個(gè)功能：加權(quán)、求和以及激活，其模型如圖4所示。

圖4 神經(jīng)元模型Fig.4 Neuron model

(1) 加權(quán)：接受前一層傳來(lái)的數(shù)據(jù)，并且每一個(gè)輸入信號(hào)都對(duì)應(yīng)一個(gè)不同的權(quán)重；

(2) 求和：確定全部輸入在不同權(quán)重下的組合效果；

(3) 激活：將數(shù)據(jù)的加權(quán)和輸入到激活函數(shù)，并將函數(shù)輸出結(jié)果傳遞給下一層。

激活函數(shù)采用tansig函數(shù)，可以將一個(gè)實(shí)數(shù)映射到[-1,1]區(qū)間內(nèi)，表達(dá)式為

y=2/[1+exp(-2x)]-1.

(18)

如圖5所示。

圖5 tansig函數(shù)Fig.5 Tansig function

2.2 混合濾波器設(shè)計(jì)

將BP神經(jīng)網(wǎng)絡(luò)的學(xué)習(xí)和適應(yīng)能力應(yīng)用到卡爾曼濾波中，使濾波器兼?zhèn)鋵W(xué)習(xí)能力和估計(jì)性能，從而提高系統(tǒng)的跟蹤性能。

輸出濾波器的狀態(tài)估計(jì)方程為

(19)

將估計(jì)方程做等價(jià)變形，有

(20)

圖6 算法原理Fig.6 Algorithm principle

2.3 算法流程

(1) 網(wǎng)絡(luò)結(jié)構(gòu)與參數(shù)配置

X=[xyzvxvyvzaxayaz].

(21)

狀態(tài)估計(jì)9維，殘差3維，增益為9×3=27維，則輸入節(jié)點(diǎn)數(shù)為39個(gè)，輸出節(jié)點(diǎn)數(shù)為9個(gè)。本文設(shè)置隱含層為2層，節(jié)點(diǎn)數(shù)分別為45，9。學(xué)習(xí)率采用Levenberg-Marquardt算法。

(2) 采集訓(xùn)練數(shù)據(jù)集

本文采用有監(jiān)督學(xué)習(xí)的方式離線訓(xùn)練網(wǎng)絡(luò)，典型彈道及一些常見(jiàn)的機(jī)動(dòng)彈道等彈道的跟蹤結(jié)果進(jìn)行數(shù)據(jù)采集，包含15 000個(gè)輸入輸出時(shí)間序列數(shù)據(jù)。

(3) 數(shù)據(jù)預(yù)處理

BP神經(jīng)網(wǎng)絡(luò)的激活函數(shù)將數(shù)據(jù)映射到(-1,1)，需要對(duì)輸出數(shù)據(jù)進(jìn)行歸一化處理，考慮到神經(jīng)網(wǎng)絡(luò)的輸入數(shù)據(jù)之間存在很大的量級(jí)差別，為了避免因量級(jí)差異引起網(wǎng)絡(luò)誤差預(yù)測(cè)較大，需要對(duì)輸入數(shù)據(jù)進(jìn)行預(yù)處理。采用最大最小值方法進(jìn)行數(shù)據(jù)歸一化。用神經(jīng)網(wǎng)絡(luò)得到的預(yù)測(cè)結(jié)果還需要進(jìn)行數(shù)據(jù)的反歸一化處理，反歸一化是歸一化的逆過(guò)程。

(4) 網(wǎng)絡(luò)訓(xùn)練

設(shè)置網(wǎng)絡(luò)最大訓(xùn)練次數(shù)為3 000，訓(xùn)練精度為1×10-5，是輸入輸出數(shù)據(jù)歸一化后的精度要求。設(shè)置隱含層的激活函數(shù)為tansig函數(shù)，輸出層的激活函數(shù)為purelin函數(shù)。

(5)卡爾曼BP神經(jīng)網(wǎng)絡(luò)融合輸出

將訓(xùn)練好的BP神經(jīng)網(wǎng)絡(luò)嵌入到卡爾曼濾波器中，形成混合濾波器。在對(duì)機(jī)動(dòng)目標(biāo)狀態(tài)估計(jì)時(shí)，經(jīng)過(guò)卡爾曼濾波得到網(wǎng)絡(luò)的輸入數(shù)據(jù)，需要一步數(shù)據(jù)歸一化處理。而后將歸一化數(shù)據(jù)輸入網(wǎng)絡(luò)得到輸出數(shù)據(jù)，需要輸出數(shù)據(jù)反歸一化處理，得到網(wǎng)絡(luò)預(yù)測(cè)的估計(jì)誤差。

同時(shí)需要對(duì)網(wǎng)絡(luò)的估計(jì)誤差設(shè)置閾值檢測(cè)，若估計(jì)誤差超過(guò)閾值限制，則認(rèn)為此次的網(wǎng)絡(luò)預(yù)測(cè)值錯(cuò)誤，對(duì)濾波估計(jì)值不進(jìn)行修正。即：

(22)

3 仿真校驗(yàn)

3.1 仿真場(chǎng)景

為了驗(yàn)證本文所提算法的有效性，采用數(shù)學(xué)仿真方法進(jìn)行驗(yàn)證，具體仿真場(chǎng)景如下：起點(diǎn)高度h0為45 km，起點(diǎn)經(jīng)度λ0為152°E，起點(diǎn)緯度φ0為25°N，初速度v0為6 km/s，初始航跡角γ0為0°，初始航向角ψ0為260°，初始攻角α0為5°，初始傾側(cè)角σ0為5°，仿真場(chǎng)景如圖7所示。

圖7 仿真場(chǎng)景Fig.7 Scene simulation

3.2 仿真結(jié)果

設(shè)置系統(tǒng)測(cè)量噪聲服從干擾高斯分布：

vk～(1-α)N(0,Rk)+αN(0,λRk),

(23)

是2種高斯噪聲的疊加。對(duì)目標(biāo)的運(yùn)動(dòng)軌跡進(jìn)行跟蹤測(cè)試。本文設(shè)置兩臺(tái)紅外探測(cè)器(S1,S2)，高度為30 km，分別布置在(120.0 E,42.9 N)，(116.5 E,40.1 N)完成對(duì)目標(biāo)的探測(cè)定位。

選擇指標(biāo)為均方根誤差(RMSE):

(24)

式中:N為蒙特卡羅仿真次數(shù)；下標(biāo)i代表狀態(tài)向量的第i個(gè)分量；j代表濾波第幾步；k為第k次蒙特卡羅仿真。

跟蹤時(shí)間為150 s，進(jìn)行50次蒙特卡羅仿真結(jié)果如圖8，9所示。

圖8 位置RMSE對(duì)比圖Fig.8 Speed RMSE comparison

圖9 速度RMSE對(duì)比圖Fig.9 Location RMSE comparison

由上述仿真結(jié)果可知，經(jīng)過(guò)BP神經(jīng)網(wǎng)絡(luò)補(bǔ)償后的混合濾波器明顯收斂速度更快，而且跟蹤誤差要小很多。跟蹤精度有大幅度的提高。軌跡更加光滑。說(shuō)明輸入預(yù)測(cè)估計(jì)值與濾波估計(jì)值的差值、濾波增益和新息，通過(guò)神經(jīng)網(wǎng)絡(luò)可以很準(zhǔn)確地預(yù)測(cè)估計(jì)誤差。

4 結(jié)束語(yǔ)

本文針對(duì)高超聲速飛行器運(yùn)動(dòng)軌跡跟蹤，在當(dāng)前統(tǒng)計(jì)模型的基礎(chǔ)上，提出了一種高超聲速飛行器在線反饋濾波算法。仿真結(jié)果表明，此在線反饋濾波算法應(yīng)對(duì)高超聲速飛行器具有良好的跟蹤能力，利用此算法可以獲得更高的濾波精度和穩(wěn)定性。