[發明專利]基于指數加權衰減記憶濾波的彈道軌跡形成方法

無效

• 申請號: 201310723456.6
• 申請日: 2013-12-24
• 公開（公告）號: CN103744058A
• 公開（公告）日: 2014-04-23
• 發明（設計）人: 胡居榮;周晶;歐陽廣帥;韓亞紅;曹寧
• 申請（專利權）人: 河海大學
• 主分類號: G01S7/02
• 分類號: G01S7/02
• 代理公司: 南京經緯專利商標代理有限公司 32200
• 代理人: 朱小兵

1.一種基于指數加權衰減記憶濾波的彈道軌跡形成方法，其特征在于，包括步驟如下：

步驟1、根據彈道模型建立基于狀態變量的彈道方程和觀測方程，將彈道方程和觀測方程離散化、線性化，從而建立彈道方程和觀測方程的近似線性方程；

步驟2、將狀態變量與標稱狀態作差得到狀態偏差，建立基于狀態偏差的卡爾曼濾波遞推方程，對雷達觀測的一段飛行彈道參數進行濾波，在濾波過程中減少彈道方程在線性化過程中產生的誤差；

步驟3、采用指數加權衰減記憶濾波方法以抑制步驟2中對飛行彈道參數濾波過程中產生的濾波發散，為彈道外推提供可靠的濾波數據；

步驟4、根據濾波后的飛行彈道參數，利用彈道模型外推彈道落點，形成彈道軌跡。

2.根據權利要求1所述的一種基于指數加權衰減記憶濾波的彈道軌跡形成方法，其特征在于，步驟1中所述建立基于狀態偏差的近似線性方程具體如下：

步驟101，建立彈道方程和觀測方程；

在地面坐標系中建立彈道方程：

dXdt=dxdtdydtdzdtdvxdtdvydtdvzdt=vxvyvz-KDSCx(Ma)2mρvr(vx-wx)+KDSCx(Ma)2mρvrwx′-KDSCx(Ma)2mρvrvy-g-KDSCx(Ma)2mρvr(vz-wz)+KLBZ+KDSCx(Ma)2mρvrwx′]]>

式中：X為t時刻的狀態變量；表示對時間t求導數；x，y，z分別為彈丸的距離、高度、側偏；v_x,v_y,v_z分別為彈丸速度的三軸分量；v_r為含風速在內的彈丸實際速度；w_x,w_z分別為縱風和橫風風速，在標準條件下均取為0；w′_x，w′_z分別為縱風和橫風的隨機量；S為彈丸最大截斷面積，d為彈徑；m為彈體質量；Ma為馬赫數，C_s為聲速，其中k為比熱比，對于空氣k＝1.404；R_d為干空氣氣體常數；τ為虛溫，是把濕空氣折合成干空氣時對氣溫的修正；C_x(Ma)為空氣阻力系數，它是馬赫數Ma的函數；ρ為空氣密度，其中，p為氣體壓強；g為海拔y處的重力加速度，其中g₀為海平面重力加速度，r₀為地球半徑；K_D為阻力符合系數；B_Z為側向升力加速度；K_L為由動力平衡角引起的側向加速度時彌補誤差的修正系數，又稱靜力矩符合系數；

令狀態變量X＝[x,y,z,v_x,v_y,v_z,c,ac]^T＝[x₁,x₂,x₃,x₄,x₅,x₆,x₇,x₈]^T，式中c為彈道系數，ac為修正的彈道系數，x₁～x₈為狀態變量X的8個狀態，分別與x,y,z,v_x,v_y,v_z,c,ac對應；

令x_4r＝x₄-w_x表示彈丸在x方向受風速影響后的實際速度；x_6r＝x₆-w_z，表示彈，丸在z方向的實際飛行速度，則

取作為中間變量；同時將聲速帶入，馬赫數表示為Ma=vrCs=vr20.0468τ,]]>將彈道方程改寫為：

dXdt=f(X)+ΓW=x4x5x6-x7Nx4r-x7Nx5-g-x7Nx6r+x8Bz00+03×3I3×302×3x7Nwx′0wz′]]>

式中，f(X)為狀態變量X的函數；Γ為噪聲驅動陣，Γ＝[0_3×3??I_3×3??0_2×3]^Tx₇N；W為系統噪聲，W＝[w′_x??0??w′_z]^T；0_3×3表示3×3的零矩陣；I_3×3表示3×3的單位陣；0_2×3表示2×3的零矩陣；

建立觀測方程：

Z=h(X)+V=x2+y2+z2arctanzxarctanyx2+z2+V=x12+x22+x32arctanx3x1arctanx2x12+x32+V]]>

式中，h(X)=x2+y2+z2arctanzxarctanyx2+z2T;]]>V為雷達量測噪聲，且V服從均值為零、方差為R的正態分布：

R=σr2σB2σe2]]>

式中，σ_γ、σ_β、σ_ε分別為炮位雷達在天線坐標系中的隨機測量誤差；

步驟102，對彈道方程和觀測方程離散化、線性化，得到彈道方程和觀測方程的近似線性方程；

在采樣區間[t_k,t_k+1]上將f(X(t))近似展開并將彈道方程離散化，得到：

Xk+1=Xk+f(Xk)T+A(Xk)f(Xk)T22+Γkqk]]>

式中，X_k為k時刻的狀態變量；X_k+1為k+1時刻的狀態變量；A(X)是f(X)關于X的偏導，即模型噪聲T為時間間隔，T＝t_k+1-t_k；

令F(Xk)=Xk+f(Xk)T+A(Xk)f(Xk)T22,]]>則

X_k+1＝F(X_k)+Γ_kq_k

若不考慮系統噪聲，則系統k時刻的標稱狀態為

Xk*=F(Xk-1*)]]>

圍繞標稱狀態將F(X_k)泰勒展開，并取一階近似，即得到彈道方程的近似線性化方程：

Xk+1=F(Xk*)+∂F(X)∂X|X=Xk*(Xk-Xk*)+Γkqk;]]>

由于F(Xk*)=Xk+1*,]]>并令φk=∂F(X)∂X|X=Xk*=I+A(Xk*)T+d[A(X)f(X)]dX|X=Xk*T22,]]>則彈道方程的近似線性化方程改寫為：

Xk+1=Xk+1*+φk(Xk-Xk*)+Γkqk;]]>

將觀測方程中的非線性函數h(X)在標稱狀態處展開成泰勒級數，并取一階近似，即得到觀測方程的近似線性化方程：

Zk=Zk*+∂h(X)∂X|X=Xk*[Xk-Xk*]+Vk]]>

令觀測矩陣則觀測方程的近似線性化方程改寫為：

Zk=Zk*+Hk(Xk-Xk*)+Vk]]>

式中，Z_k為k時刻的觀測變量，為k時刻的觀測變量標稱值，V_k為k時刻的觀測噪聲。