基于空氣動力學模型含雷達散射截面觀測的聯合目標跟蹤識別

郭鵬，白亮，武夢潔，蔣宏

(1.北京航空航天大學飛行器控制一體化技術重點實驗室，北京100191;2.洛陽電光設備研究所 火力控制技術國防科技重點實驗室，河南 洛陽471009)

0 引言

聯合目標跟蹤識別(JTC)技術是信息融合領域新興的一個研究方向［1］，大致可分為以文獻［2 -4］為代表的基于質點運動模型的JTC 技術和以文獻［5 -8］為代表的基于剛體運動模型的JTC 技術兩類?；谫|點運動模型的JTC 技術主要是依賴運動學信息來進行目標分類的，通常只能實現廣泛意義上的分類，如“戰斗機/客機”;而要識別目標飛機型號，如F-15 戰斗機/VFY-218 戰斗機，則需要采用基于剛體運動模型的JTC 技術，并在觀測中至少包含一個目標特征，該特征有可能區分特定的目標型號。例如，文獻［5］建議使用高分辨距離成像特征，文獻［6 -7］選擇了雷達散射截面(RCS)特征，因為這些分類特征不僅依賴運動學信息，還包含目標的屬性信息，能夠使分類結果更準確，同時，目標飛機型號的正確識別也有助于提高目標跟蹤精度。

為了正確識別出目標的型號，向指揮和制導系統提供更精確的信息。本文以地基被動警戒雷達為背景，利用地面上多個商業廣播所發射的無線電調頻(FM)信號，提出一種有效的基于空氣動力學模型含RCS 觀測的聯合目標跟蹤識別技術。

1 建立目標飛機的狀態模型

在使用真實的空氣動力學模型來描述目標的平動和轉動時，需要考慮特定飛機的具體參數，因此目標狀態和目標類型之間是完全耦合的。

1.1 對飛機目標的平動建模

有3 種類型的力作用于飛機上:重力，推力和氣動力。重力表示在慣性坐標系下，向下為正，W =mg，其中m 為飛機質量。推力的大小由飛行員控制，對于一個給定類型的飛機，推力有上界和下界，上界由引擎最大的輸出決定，下界由維持足夠升力的最低飛行速度決定。推力用T 表示，并假定它沿著機體坐標系的x 軸。氣動力表示在氣流坐標系下:

式中:D、C、L 分別為阻力、側力、升力;CD、CC、CL分別為無量綱的阻力系數、側力系數、升力系數;ρ 為飛機飛行高度上的空氣密度;v 為飛機飛行速度;S機翼面積。只要側滑角β 足夠小，可假定CC≈0.對亞音速和超音速飛行，在一定的α(迎角)范圍內，可以用以下近似:

式中:CLα為升力線斜率;CDmin為零升阻力系數;KD為展弦比，它們主要是飛機形狀和速度的函數。在有限的速度變化范圍內，可忽略這些參數對速度的依賴?？紤]飛機所受的重力、推力和氣動力，將飛機平動的總力F 表示在慣性坐標系下:

式中:Oib為機體坐標系到慣性坐標系變換矩陣;Obw為氣流坐標系到機體坐標系變換矩陣。由于有CD、CL，(4)式不僅是一個剛體運動模型，還是一個設計明確的飛行。由(4)式可知，要計算F 需要兩套參數。第一套是目標類型參數Ac={CLα，CDmin，KD，S，m}，下標c 表示這些參數所代表的目標類型。目標身份類型一旦被確定，這套參數就完全被確定了;第二套是目標運動學參數{pz，v，ψ，θ，φ，T}，其中ψ、θ、φ 為3 個歐拉角，z 方向的高度pz和空氣密度ρ 相關。

1.2 對飛機目標的轉動建模

使飛機產生滾動，俯仰和偏航的力矩都是在飛行員的直接控制下。因此，可以通過對飛行員的意圖建模來代替對力矩建模。飛行員的意圖對轉動的影響可以定義為一組轉動模態，每個轉動模態對應一個不同類型的飛行，如轉彎或爬升。這樣，飛機的轉動只不過是飛行員命令的一系列未知機動{γ1，γ2，…，γk}的結果，其中γk為飛機從kΔt 到(k +1)Δt時刻的操作模態。

對于狀態模型，假設目標的轉動可以用兩個縱向模態和一個橫向模態描述。這3 個模態為:模態1(CV)，常速固定高度飛行;模態2(PT)，縱向爬升或下降;模態3(CT)，恒速率協調轉彎。

狀態模型需要隨時間變化的歐拉角，因此可以繞過力矩建模問題，直接對每個轉動模態進行更有意義的歐拉角建模。

表1總結了每個轉動模態下的歐拉角模型。其中，帶“* ”變量表示驅動飛機轉動的未知飛行員輸入，ω*表示轉彎角速率，Δθ*表示俯仰角的變化，αss表示以規定速度維持水平飛行所需的穩態迎角。最后，假設模態序列{γk}按照模型轉移概率矩陣為P的一階馬爾可夫鏈演變，矩陣P 中的元素被定義為P(i，j)=Pr (γk=j|γk-1=i).

表1 歐拉角模型Tab.1 Euler angle model

1.3 完整的飛行模型

圖1給出一個由xk-1到xk的偽碼描述。每次掃描，首先判斷該模型是否發生了模態改變:如果沒有改變，上一時刻的駕駛員輸入被添加少量噪聲而作為當前時刻的駕駛員輸入;如果有改變，需要獲取新的飛行員輸入。圖中Tss為以規定速度維持水平飛行所需的穩態推力，分別為縱向爬升的俯仰角變化、縱向爬升的推力以及恒速率協調轉彎的轉彎角速率，都服從均勻分布，且新的飛行員輸入從這些均勻分布中取值。通過讓這些分布依賴于目標類型，每個飛機實際的運動學限制可以被使用，從而有利于目標識別。

圖1中偽碼的最后一行調用一個數值積分子程序，它利用k-1 時刻的狀態和當前飛行員的輸入，完成了飛機的平動和轉動，并將狀態從k-1 時刻推進到k 時刻。圖2給出了這個數值積分子程序的偽碼。

圖1 飛機狀態模型偽碼Fig.1 Pseudo-code of state model for flight

到此為止，建立了一個基于飛機實際的平移和轉動的狀態模型。

2 建立被動雷達的觀測模型

被動雷達可以觀測到時延、多普勒頻移及RCS.在三維空間要唯一地定位目標，至少需要3 個源于不同發射器的延遲測量;同樣，要唯一地確定目標的三維速度，也至少需要3 個源于不同發射器的多普勒頻移測量。因此，為簡單起見，采用3 個發射器，單被動雷達接收器。這樣，k 時刻被動雷達的觀測為Zk={(τm，k，dm，k，σm，k)}，m =1，2，3，其中(τm，k，dm，k，σm，k)分別為k 時刻源于第m 個廣播電臺的時延、多普勒及RCS 觀測量。

圖2 狀態數值積分子程序偽碼Fig.2 Pseudo-code of state integration subroutine

2.1 對時延和多普勒頻移建模

式中:τm，k為時延;pk為k 時刻目標的位置;pr為被動雷達接收站的位置;c 為光速。

式中:dm，k為多普勒頻移;vk為k 時刻目標的速度矢量;nm為位置的單位向量;λm為電磁波的波長。

式中:ptm為第m 個廣播電臺的位置;pr為被動雷達接收站的位置，如果讓慣性坐標系的原點與被動雷達接收器重合，則pr=0.時延和多普勒頻移的測量噪聲可假設為零均值高斯白噪聲，可得時延、多普勒觀測量的似然函數為

式中:h(τ)m，k、h(d)m，k由(6)式和(7)式定義。

2.2 對RCS 建模

如果假設被動雷達系統是相干的(通過使用相當有限的調頻無線電信號帶寬)，理論上可以還原k時刻源自第m 個發射臺的復數散射系數sm，k.然而，用于前端濾波器組的積分時間不能用散射系數的相位來分類，故使用|sm，k|，m =1，…，M 來分類，因為RCS 被定義為σm，k= |sm，k|2，這相當于使用RCS.RCS 是一個非負量，假設其測量噪聲為高斯分布不恰當，因此可以對(相干)匹配濾波器的同步和正交分量的熱噪聲建模。對RCS 測量噪聲，采用的具體模型為

式中:上標代表sm，k的實部和虛部;w1，w2為獨立同分布、零均值、方差為的高斯隨機變量。假定熱噪聲的方差獨立于收到的FM 信號，可去掉的下標m.最終生成的似然函數可以表示為

式中:I0是第一種零階修正貝塞爾函數;σm，k是k 時刻對源自第m 個發射臺的RCS 觀測量;σ2s 是影響散射系數sm，k測量的熱噪聲方差，m=1，…，M.

2.3 完整的觀測模型

假設各觀測量τm，k、dm，k、σm，k的噪聲是獨立的，粒子濾波所需要的觀測似然函數將是各觀測量似然函數的乘積，即

式中:p(σm，k|xk)由(12)式決定。

3 Matlab 實時調用FEKO 軟件

FEKO 是一種用于RCS 理論計算的電磁仿真軟件?？梢杂赡繕说膶崟r位置、姿態角、電磁波發射站和雷達接收站的位置得到FEKO 軟件的輸入信息，即入射電磁波的方位角、高低角，散射電磁波的方位角、高低角。然后經Matlab 調用FEKO，可以得到不同型號飛機的實時RCS.調用過程如下:

1)運行PREFEKO，對飛機模型進行預處理，飛機模型如圖3所示。

2)運行EDITFEKO，對“.pre”文件中的控制量、激勵和遠場等信息進行配置。配置好后的參數如圖4所示。

3)運行FEKO 軟件，對配置好的模型求解RCS值。

圖4 EDITFEKO 界面Fig.4 The interface of EDITFEKO

4)運行Matlab 的m 函數實現調用，獲取RCS.

4 粒子濾波執行

由于RCS 觀測量，目標姿態和目標位置等狀態量之間沒有解析關系，并且飛機實際平動和轉動的狀態模型沒有簡單的顯式表達式，因此不能使用擴展卡爾曼濾波，然而粒子濾波提供了解決方案。

正則化粒子濾波［9］(RPF)是一種基于重采樣的改進算法，它與傳統粒子濾波的區別在于:傳統粒子濾波從離散近似分布中重采樣，而正則化粒子濾波從連續近似分布中重采樣。正則化粒子濾波可以有效緩解重采樣造成的粒子多樣性喪失問題，能夠獲得較高的精度。其算法如圖5所示，S 表示統計協方差陣，K 表示Epanechnikov Kernel，一般情況下取近似值，用高斯核代替，hopt表示最優帶寬。

目標有3 種飛行模態，本文將采用交互式多模型正則化粒子濾波算法來實現，并用交互式多模型粒子濾波算法與其比較。

5 仿真實驗

5.1 仿真場景及參數設置

在仿真過程中，對戰斗機VFY-218 進行跟蹤和識別，無人機UAV 充當目標庫中的參考目標。采用Monte Carlo 仿真來驗證算法的有效性，仿真的次數為NMC=50，每次仿真的采樣次數為Nk=75.仿真軟件為Matlab 和FEKO，計算機的硬件為AMD Athlon(tm)II X2 250 Processor 3.01 GHz，1.75 GHz 內存。

圖5 正則化粒子濾波Fig.5 Regularized particle filter

目標初始位置p0=［-95 000，50 000，4 000］T，初始速度為v0=［110，160，0］T，初始俯仰角、滾轉角和偏航角分別為θ=4°，φ=0°和ψ=55.49°.3 個發射器和被動雷達的位置分別為pt1=［48 000，-18 000，0］T，pt2=［- 23 000，38 000，0］T，pt3=［42 000，23 000，0］T，pr=［0，0，0］T.采樣周期Δt=0.4 s，積分步長Δi =0.2 s，3 個發射器交替向被動雷達發送信息以獲得觀測值，每輪處理3 個發射器信息的時間為1.2 s.粒子濾波的粒子數Np=600，模型轉移概率矩陣P =［0.98，0.01，0.01;0.01，0.98，0.01;0.01;0.01，0.98］，初 始 模 型 概 率μcv(0)=0.8，μct(0)= 0.1，μpt(0)=0.1;傳輸時延的觀測噪聲σ2τ=0.000 1 ms2，多普勒頻移的觀測噪聲σ2d=1 Hz2，RCS 的觀測噪聲σ2s=64 m4.

目標的運動模態如表2所示，其真實軌跡和跟蹤軌跡如圖6所示。

5.2 跟蹤結果及分析

圖7、圖8和圖9分別給出了姿態角跟蹤圖、速度跟蹤圖和位置均方根誤差(RMSE)，圖中的豎直虛線表示該時刻目標運動模態發生變化。

表2 目標的運動模態及參數Tab.2 The motion model and parameters of target

圖6 目標運動軌跡Fig.6 Target motion trajectory

1 ～20 采樣區間在固定高度做勻速直線運動，飛機受力平衡，速度和姿態角基本不變，此時交互式多模型粒子濾波(IMMPF)和交互式多模型正則化粒子濾波(IMMRPF)在x、y 和z 方向的RMSE 相當，都較小。

21 ～35 采樣區間在豎直平面內做爬升運動，飛行員給飛機施加的推力T 大于空氣阻力與重力的分力，俯仰角變化Δθ*=6° ＞0°，故飛機向上爬升。由于飛機在縱向發生機動，故此時的俯仰角變化和速度變化存在跟蹤滯后現象。與勻速運動相比，爬升運動的IMMPF 和IMMRPF 在x、y 和z 方向的均方根誤差都增大，且增大的幅度相當。

36 ～50 采樣區間在固定高度做勻速直線運動，由于受到爬升運動俯仰角和速度滯后的影響，z 方向速度不能減小到0.在x 方向上，IMMRPF 均方根誤差要小于IMMPF，而在y 和z 方向上均方根誤差相當。與爬升運動相比，勻速運動的IMMPF 和IMMRPF 在x、y 和z 方向的均方根都減小，表明目標在非機動情況下的跟蹤性能要優于機動情況下。

51 ～75 采樣區間做恒速率協調轉彎運動，在水平方向飛行員施加的推力和空氣阻力平衡，飛機機身傾斜使升力在豎直方向的分力平衡重力，在水平方向的分力提供向心力。此時x 方向的速度不斷增大，y 方向的速度不斷減小，滾轉角從0°變到-75°.IMMRPF 的偏航角在協調轉彎過程中都要優于IMMPF，滾轉角在51 ～60 采樣區間也優于IMMPF.在x 和y 方向，IMMRPF 的速度相對于IMMPF 更接近真實值，因此IMMRPF 的均方根誤差明顯小于IMMPF.

圖7 姿態角跟蹤圖Fig.7 Euler angle tracking

由圖7和圖8可以看出，目標發生機動之后，對于目標狀態的跟蹤存在滯后現象。這種滯后反映在模態概率上為機動檢測時刻總是滯后于機動發生時刻。如圖10所示，IMMPF 的模型概率在51 ～60 采樣區間要略優于IMMPF，其他時候相當。

在第n 次Monte Carlo 仿真中，第k 次采樣時刻的位置誤差定義為

圖8 速度跟蹤圖Fig.8 Velocity tracking

式中:NL為跟蹤失敗的次數。

表3 跟蹤性能Tab.3 Tracking performance

從表3可以看出，IMMRPF 的綜合跟蹤性能要優于IMMPF.

圖9 位置均方根誤差Fig.9 The position RMSEs

圖10 3 個模態的概率Fig.10 The probabilities of three models

5.3 分類結果及分析

在第k 次掃描時，目標類型的概率可以表示為式中:Nc表示類型的數目;c 表示從離散類型標簽中獲取信息的判斷準則;z1:k表示1 到k 時刻的觀測序列。而Pr(c(xk)=j|z1:k)的表達式可以利用粒子濾波的權值近似表示為

從表4可以看出，在聯合目標跟蹤識別中，IMMPF 和IMMRPF 正確識別目標的概率都較高，因而能夠較準確識別出目標的類型，并且IMMRPF 的類型平均識別率優于IMMPF.IMMRPF 利用了正則化粒子濾波的優勢，從后驗密度的連續近似采樣得到新的粒子，因此避免了大權值粒子被多次復制，保持了粒子的多樣性，故類型平均概率要優于IMMPF.

表4 類型平均識別率Tab.4 The average identification probability

在正確識別出戰斗機的類型后，分類器將戰斗機的RCS 值和類型參數Ac={CLα，CDmin，KD，S，m}輸入到跟蹤器中，其中CLα= 2.75 rad-1，CDmin=0.013，KD=2.68，S=42.7 m2，m=13 397 kg.

6 結論

本文對飛機目標建立了平動和轉動的空氣動力學模型，利用各飛行模態的特點將全狀態變量簡化為關注目標類型、飛行模態和飛行員輸入;通過Matlab 調用FEKO 電磁仿真軟件，將跟蹤器的輸出與各傳感器的位置信息綜合換算得到FEKO 電磁仿真軟件的輸入，有效縮短了求解RCS 的時間，然后將分類器得到的實時RCS 值和類型參數輸入到跟蹤器當中，完成目標的聯合跟蹤識別。

仿真實驗結果表明，粒子濾波能夠較好地處理高度非線性的系統模型，IMMRPF 在位置均方根誤差、跟蹤成功的次數和類型平均概率等方面都要優于IMMPF，具有較好的綜合跟蹤性能。

References)

［1］ 單甘霖，梅衛，王春平.聯合目標跟蹤與分類技術的進展及存在問題［J］.兵工學報，2007，28(6):733 -738.SHAN Gan-lin，MEI Wei，WANG Chun-ping.The development and problem of joint target tracking and classification technology［J］.Acta Armamentarii，2007，28(6):733 -738.(in Chinese)

［2］ Yang W，Fu Y W，Long J Q，et al.Joint detection，tracking and classification of multiple targets in clutter using the PHD filter［J］.IEEE Transactions on Aerospace and Electronics Systems，2012，48(4):3594 -3609.

［3］ Melzi M，Ouldali A.Joint multiple target tracking and classification using the unscented Kalman particle PHD filter［C］∥2011 IEEE 9th International New Circuits and System Conference.Bordeaux，FR:IEEE，2011:534 -537.

［4］ MEI W，SHAN G L，LI X R.Simultaneous tracking and classification:a modularized scheme［J］.IEEE Transactions on Aerospace and Electronic Systems，2007，43(2):487 -506.

［5］ Miller M I，Srivastava A，Grenander U.Conditional-mean estimation via jump-diffusion process in multiple target tracking/recognition［J］.IEEE Transactions on Signal Processing，1995，43(11):2678 -2690.

［6］ Guo Y F，Peng D L，Chen H J，et al.A kernel particle filter algorithm for joint tracking and classification［C］∥2012 15th International Conference on Information Fusion.Singapore:IEEE，2012:2386 -2391.

［7］ Herman S，Moulin P.A particle filtering approach to FM-band passive radar tracking and automatic target recognition［C］∥2002 IEEE Aerospace Conference Proceedings.Big Sky，Montana:IEEE，2002:1789 -1808.

［8］ Herman S M.A particle filter approach to joint passive radar tracking and target classification［D］.US:Graduate College of the University of Illinois，2002.

［9］ 羅飛騰.目標跟蹤的粒子濾波技術研究［D］.合肥:中國科學技術大學，2010:16 -17.LUO Fei-teng.Research on particle filtering in target tracking［D］.Hefei:University of Science and Technology of China，2010:16 -17.(in Chinese)