Distribution of Zero-effort Miss Distance Estimation Errors in Continuous-time Controlled System With Mode Mismatch
-
摘要: 针对大机动目标拦截中零控脱靶量（Zero-effort miss distance，ZEM）估计误差分布的求解问题，在线性估计器与独立模式辨识器的配置下，提出了一种解析求解方法.在存在固定的模式辨识延迟下，ZEM的估计误差服从有偏的高斯分布.最后，通过与蒙特卡洛仿真的比较，验证了本文方法的正确性.Abstract: For the problem of solving the error distribution of zero-effort miss distance (ZEM) in highly maneuvering target interception, an analytical approach is proposed for the configuration of a linear estimator and a separate mode decision-maker. In the case of a fixed mode decision delay, the error of ZEM follows a biased Gaussian distribution. Finally, the correctness of the proposed method is verified by comparison with Monte Carlo simulation.
-
海湾战争以后, 战术弹道导弹(Tactical ballistic missiles, TBM)逐渐成为防空系统的拦截对象, 而具有潜在机动能力的TBM的出现, 则给现有防空系统提出严峻的挑战.当前, 高速大机动目标拦截问题引起了国内外各研究团队的广泛关注, 其中最具代表的是以色列的Shinar团队[1-4].此类目标拦截的最大挑战在于它要求特别小的脱靶量(Miss distance, MD)甚至是直接碰撞.对于该问题, 由于导弹加速度饱和、非高斯噪声以及系统的非线性, 确定性等价原理[5]已不再适用.然而, 当导弹和目标的相对状态满足可观测性条件时, 部分分离定理仍然成立, 此时估计器仍可独立于导引律进行设计[6].
在估计器给定后, 建立零控脱靶量(Zero-effort miss distance, ZEM)的估计误差模型便具有重要意义.一方面, ZEM是导引律中的一项关键输入参数, 不同导引律间的区别仅在于不同的ZEM计算方式.因此, ZEM的估计精度将直接影响制导性能, 合适的ZEM估计误差模型可以有效指导制导系统的设计.另一方面, 作为评价高度大机动目标拦截性能的重要指标, 脱靶量的解析计算方法中需要用到ZEM的估计误差模型[7-10].
目前, 在给定系统动态模型、估计器、控制策略、以及具体的扰动和噪声模型后, ZEM的估计误差分布主要通过大量蒙特卡洛仿真得到.然而, 这种后验方法在控制系统设计阶段并不十分适用. Moldavskaya等[11]提出了一种求解ZEM估计误差分布的解析方法.采用成形滤波器技术近似目标加速度指令, 近似的指令由白噪声通过一个一阶线性系统后得到, 且与原指令具有相同的自相关函数.假定初始的估计误差为零, 测量噪声和状态噪声为互相独立的零均值高斯白噪声, 则得出ZEM的估计误差服从零均值高斯分布, 其中状态估计误差的方差矩阵满足Riccati方程.
在末制导系统中, 由于系统误差、观测噪声等不确定因素以及不完全的状态观测, 估计器需同时扮演观测器和滤波器的角色, 即要求其同时兼顾未知状态重建及测量噪声滤波这两方面的性能, 但使用固定带宽的成形滤波器难以同时兼顾估计精度和响应速度.相关的研究工作表明, 组合使用独立模式辨识器和一个低带宽高精度估计器更适合高速大机动目标的拦截问题, 可显著提升制导性能[12-14].实际上, 对于高速大机动目标拦截, 雷达和光电导引头能够观测到目标的特征信息, 这些特征为模式辨识器的设计提供了额外的信息且可改善系统对目标机动的响应速度[15-18]. Fan等在引入独立模式辨识器的基础上, 分析了状态估计的误差特性[19].考虑到ZEM对大机动目标拦截的重要意义, 本文在此基础上推导模式失配条件下ZEM的估计误差分布.
1. 问题描述
1.1 系统运动模型
与大多数机动目标拦截文献类似[7-11, 20], 本文仅考虑一弹一目的平面拦截情形.如图 1所示, 用$P$和$E$分别表示导弹(追方)和目标(逃方), 并作如下假设[11, 19-21]:
1) $P$和$E$的控制动态可用一阶转移函数来近似, 相应的时间常数分别为${\tau _p}$和${\tau _e}$;
2) $P$和$E$的速度是恒定的, 分别用${V_p}$和${V_e}$表示;
3) $P$和$E$的横向加速度有界, 其最大横向加速度分别用$a_p^{\max}$和$a_e^{\max}$表示.
在图 1中, $X$轴沿弹目初始视线方向; $Y$垂直于$X$轴; $({x_p}, {y_p})$和$({x_e}, {y_e})$分别为$P$和$E$的当前坐标; ${array}{*{20}{c}}{{\phi _i}, }&{i = p, e}{array}$分别是$P$和$E$的偏航角, 它表示速度矢量和$X$轴的夹角.假设偏航角满足$(\sin {\phi _p} = \phi _p, \sin \phi _e = \pi - {\phi _e}), $则弹目运动轨迹可以沿着初始视线进行线性化.假定接近速度恒定, 起始时刻${t_0} = 0$ s, 给定弹目起始距离${r_0}$后, 拦截的终止时刻满足:
$ \begin{equation} {t_f} \approx \frac{{{r_0}}}{{{V_p}\cos {\phi _{p}(0)} - {V_e}\cos {\phi _{e}(0)}}} \end{equation} $
(1) 剩余飞行时间定义为${t_{go}} = {t_f} - t$.
定义状态矢量${\pmb x} = {[{x_1}, {x_2}, {x_3}, {x_4}]^{\rm T}} = {[y, \dot y, a_y^e, a_y^p]^{\rm T}}$.基于上述三条假定, 整个拦截过程的时间$t \in [0, {t_f}]$, 且具有下述线性动态模型:
$ \begin{equation} \begin{array}{*{20}{l}} {{{\dot x}_1} = {x_2}, }&{{x_1}(0) = 0}\\ {{{\dot x}_2} = {x_3} - {x_4}, }&{{x_2}(0) = {V_e}{\phi _{e}(0)} - {V_p}{\phi _{p}(0)}}\\ {{{\dot x}_3} = \dfrac{{u_e} - {x_3}}{{\tau _e}}, }&{{x_3}(0) = 0}\\[2mm] {{{\dot x}_4} = \dfrac{{u_p} - {x_4}}{{\tau _p}}, }&{{x_4}(0) = 0} \end{array} \end{equation} $
(2) 其中, $x_1 = {y_e} - {y_p}$为$P$和$E$之间沿Y轴的相对距离; $x_2$为相对的横向速度; ${array}{*{20}{c}}{a_y^i, {u_i}, }&{i= p, e}{array}$分别为$i$的横向加速度和加速度指令且满足:
$ \begin{equation}\ \begin{array}{*{20}{c}} {|{u_i}(t)| \le |a_i^{\max}|, }&{i = p, e} \end{array} \end{equation} $
(3) 定义
$ \begin{equation} \gamma = \frac{{a_p^{\max }}}{{a_e^{\max }}} \end{equation} $
(4) $ \begin{equation} \varepsilon = \frac{{{\tau _e}}}{{{\tau _p}}} \end{equation} $
(5) 其中, $\gamma$称为机动过载比, $\varepsilon $通常称为敏捷系数.为了保证较好的拦截精度, 传统最优控制框架下的导引律如OGL (Optimal guidance law)通常需要过载优势满足$\gamma>3$, 比例导引和增广比例导引则需要4到5倍的过载优势, 本文面向TBM拦截应用, 考虑过载比$\gamma < 2.5$的大机动目标且$\varepsilon \approx 1$.
将动态模型写成矢量形式
$ \begin{equation}\ \begin{array}{l} {\bf{\dot {\pmb x}}}(t) = {A}{\pmb x}(t) + {{\pmb B}_1}{u_p}(t) + {{\pmb B}_2}{u_e}(t)\\ {\pmb x}(0) = {(0, {x_2}(0), 0, 0)^{\rm T}} \end{array}\ \end{equation} $
(6) 其中
$ \begin{align} &{{A}} = \left[{\begin{array}{*{20}{c}} 0&1&0&0\\ 0&0&1&{-1}\\ 0&0&{\dfrac{{-1}}{{{\tau _e}}}}&0\\ 0&0&0&{\dfrac{{-1}}{{{\tau _p}}}} \end{array}} \right]\nonumber\\ &{{{\pmb B}_1} = \left[{\begin{array}{*{20}{c}} 0\\ 0\\ 0\\ {\dfrac{1}{{{\tau _p}}}} \end{array}} \right]}, ~{{\pmb B}_2} = \left[{\begin{array}{*{20}{c}} 0\\ 0\\ {\dfrac{1}{{{\tau _e}}}}\\ 0 \end{array}} \right] \end{align} $
(7) 对状态方程(6)采用终端投影变换
$ \begin{equation}\ z(t) = {{\pmb D}^{\rm T}}{\Phi}({t_f}, t){\pmb x}(t) \end{equation} $
(8) 则可将平面拦截问题转换为一个标量问题.系统的新状态变量为零控脱靶量$z(t)$, 而系统的脱靶量则为终止时刻$t_f$的$z(t)$, 即$z(t_f)$.式(8)中, ${\pmb D} = {[1, 0, 0, 0]^{\rm T}}$; ${\Phi }({t_f}, t)$为满足齐次方程${{\dot {\pmb x}}}(t) = {A\pmb x}(t)$的状态转移矩阵, 求解得到:
$ \begin{equation}\ {\Phi }({t_f}, t) = {{\text{e}}^{{A}({t_f} - t)}} \end{equation} $
(9) 因此, 可将零控脱靶量$z(t)$表示为
$ \begin{equation}\ z(t) = {{\pmb g}^{\rm T}}(t){\pmb x}(t) \end{equation} $
(10) 其中, ${{\pmb g}^{\rm T}}(t) = [{g_1}(t), {g_2}(t), {g_3}(t), {g_4}(t)]$且
$ \begin{align} &{g_1}(t) = 1\nonumber\\ & {g_2}(t) = {t_{go}}\nonumber\\ & {g_3}(t) = \tau _e^2\left\{ {\exp \left( {\frac{{{t_f} - t}}{{{\tau _e}}}} \right) + \frac{{{t_f} - t}}{{{\tau _e}}} - 1} \right\}\nonumber\\ & {g_4}(t) = - \tau _p^2\left\{ {\exp \left( {\frac{{{t_f} - t}}{{{\tau _p}}}} \right) + \frac{{{t_f} - t}}{{{\tau _p}}} - 1} \right\} \end{align} $
(11) 1.2 马尔科夫跳变的加速度指令模型
与采用成形滤波器不同, 本文引入一个独立的模式辨识器, 采用马尔科夫跳变模型来描述目标的横向加速度控制指令:
$ \begin{equation}\ {u_e}(t) = m(t) + w(t) \end{equation} $
(12) 如图 2所示, 这里将目标横向加速度指令所在的控制空间量化为一系列离散点构成的模式集.假定: ${u_e}(t)$在这些点之间跳变; $m(t)$为离散化的目标横向加速度指令, 即目标的运动模式; $w(t)$为量化误差, 假定其为零均值的高斯白噪声, 功率谱密度为${s_w}$.有关机动目标跟踪模型集的设计方法可以参见文献[15, 22-24].
不失一般性, 假定$[0, {t_f}]$内目标只发生一次模式切换, 令${t_{sw}}$表示模式切换时刻, $m_1$, $m_2$分别表示模式切换前后目标的运动模式量, 则$m(t)$可表示为
$ \begin{equation}\qquad\quad m(t) = {m_1} + ({m_2} - {m_1})u(t - {t_{sw}}) \end{equation} $
(13) 其中, $u(t)$为阶跃函数, 定义为
$ \begin{equation} u(t) = \left\{ {\begin{array}{*{20}{c}} {1, }&{t \ge 0}\\ {0, }&{t < 0} \end{array}} \right. \end{equation} $
(14) 1.3 系统观测模型
与文献[19]类似, 本文采用的观测模型为
$ \begin{equation}\ {\pmb Y}(t) = {H\pmb x}(t) + \pmb v(t) \end{equation} $
(15) 其中观测矩阵
$ \begin{equation}\ H = \left[{\begin{array}{*{20}{l}} 1&0&0&0\\ 0&0&0&1 \end{array}} \right]\ \end{equation} $
(16) 观测噪声$\pmb v(t)$为零均值的高斯白噪声, 其协方差矩阵为$R(t)$.
2. ZEM估计误差模型的推导
考虑如图 3所示的一种典型的制导系统结构. 图 3中, 独立的模式辨识器为估计器和导引律的设计提供目标的运动模式信息.文献[14]基于该架构给出的导引律为基于逻辑的联合估计导引律.本文在这个框架下推导模式失配条件下ZEM的估计误差分布形式.
在此, 我们对模式辨识器的行为作如下假定:目标运动模式切换后经过$\Delta t$的时间延迟, 模式辨识器可给出正确的模式估计结果.目标模式切换和模式辨识器的模式决策过程如图 4所示.
将式(12)的目标加速度指令代入系统状态方程, 可以得到:
$ \begin{align} {{\pmb x}}(t) =\,&{A\pmb x}(t) + {{\pmb B}_1}{u_p}(t)+ {{\pmb B}_2}{m_1}+ \nonumber\\ & {{\pmb B}_2}({m_2} - {m_1})u(t - {t_{{\rm{sw}}}}) + \pmb w(t) \end{align} $
(17) 其中, ${\pmb w}(t) = {{\pmb B}_2}w(t)$, 其协方差矩阵$Q = {{\pmb B}_2}{\pmb B}_2^{\rm T}{s_w}$.
估计器采用Kalman滤波器, 由图 4可知估计器的动态方程为
$ \begin{align} {\dot {\pmb x}}(t) =\,&{A\pmb x}(t) + {{\pmb B}_1}u_P(t)+ {{\pmb B}_2}{m_1}+\nonumber\\ & {{\pmb B}_2}({m_2} - {m_1})u(t - {t_{{sw}}} - \Delta t) + {\pmb w}(t) \end{align} $
(18) 下面分三种情况讨论:
情形1. ${t_f} < {t_{sw}}$.在这种情形下估计器将一直保持正确的模式直到整个末制导过程结束, 即不存在模式切换.此时, 系统的状态方程为
$ \begin{equation}\ {\dot{\pmb x}}(t) = {A\pmb x}(t) + {{\pmb B}_1}u_P(t) + {{\pmb B}_2}{m_1} + {\pmb w}(t) \end{equation} $
(19) 结合式(15)观测模型可以得到下述滤波方程
$ \begin{align} {{\dot{\hat {\pmb x}}}}(t) =\,&{{A\hat {\pmb x}}}(t) + {{\pmb B}_1}u_P(t) + {{\pmb B}_2}{m_1} +\nonumber\\ & k(t)({\pmb Y}(t) - {{H\hat {\pmb x}}}(t)) \end{align} $
(20) 其中, ${{k}}(t)$为系统的连续Kalman增益矩阵且满足
$ \begin{equation}\ {{k}}(t) = {{P}}(t){{{H}}^{\rm T}}{{{R}}^{ - 1}}(t) \end{equation} $
(21) ${{P}}(t)$为预测误差的协方差矩阵, 满足下述Riccati方程
$ \begin{equation}\ {{\dot P}}(t) = {{AP}}(t) + {{P}}(t){{{A}}^{\rm T}} + {{Q}} - {{P}}(t){{{H}}^{\rm T}}{{{R}}^{ - 1}}(t){{HP}}(t)\ \end{equation} $
(22) 令${\tilde {\pmb x}}(t) = {\hat {\pmb x}}(t) - {\pmb x}(t)$表示状态估计误差, 则可以得到下述状态估计误差传递方程:
$ \begin{equation} \begin{array}{l} {\dot{\tilde {\pmb x}}}(t) = (A - k(t)H){\tilde {\pmb x}}(t) + (k(t)\pmb v(t) - \pmb w(t))\\ {\tilde {\pmb x}}(0) = {{{\tilde {\pmb x}}}_0} \end{array} \end{equation} $
(23) 若${\pmb w}(t)$和${\pmb v}(t)$相互独立且与${{\tilde {\pmb x}}_0}$无关, 并令: $F(t) = A - k(t)H$, ${\pmb \zeta }(t) = k(t){\pmb v}(t) - \pmb w(t)$, ${\Phi _F}({t_2}, {t_1}) = \exp \left\{ {\int_{{t_1}}^{{t_2}} {F(\theta )} {\rm{d}}\theta } \right\}$, 由定积分的性质很容易验证${{\Phi_F}}({t_2}, {t_1}) = {{\Phi_F}}({t_2}, {t_3})\cdot{{\Phi_F }}({t_3}, {t_1})$.求解方程(23)可以得到:
$ \begin{equation} {\tilde {\pmb x}}(t) = {{\Phi_F}}(t, 0){{\tilde {\pmb x}}_0} + \int\limits_0^t {{\Phi_F}(t, s){\pmb \zeta}(s){\rm d}s} \end{equation} $
(24) 其中, ${{\tilde {\pmb x}}_0}$为初始的状态估计误差.
定义${{\pmb \xi }(t) = {\rm E}\{ {\tilde {\pmb x}}(t)\} }$, ${{\Sigma }(t) = {\rm{var}} \{ {\tilde {\pmb x}}(t)\} }$, 求解得到:
$ \begin{equation} {\pmb \xi }(t) = {\Phi_F}(t, 0){\rm E}({{\tilde {\pmb x}}_0}) \end{equation} $
(25) $ \begin{equation} \begin{split} &~~~~~~~~~~~~~~{\Sigma }(t) = {\Phi_F}(t, 0){{\tilde P}_0}{\Phi_F}^{\rm T}{(t, 0)}+ \\ &\int\limits_0^t {{\Phi_F}(t, s)[{Q} + {k}(s){R}(s){{k}^{\rm T}}(s)]{\Phi_F}^{\rm T}{{(t, s)}}{\rm d}s} \end{split} \end{equation} $
(26) 其中, ${{{\tilde P}}_0} = {\rm E}\left\{ {[{{{\tilde {\pmb x}}}_0}-{\rm E}({{{\tilde {\pmb x}}}_0})]{{[{{{\tilde {\pmb x}}}_0}-{\rm E}({{{\tilde {\pmb x}}}_0})]}^{\rm T}}} \right\}$为初始的估计误差协方差矩阵.上述结果的推导过程见附录A.
情形2 ${t_{sw}} \leq {t_f} < {t_{sw}} + \Delta t$.需要分两段讨论:
1) $t \in [0, {t_{sw}})$, 此时估计器一直保持正确的模式, 状态估计误差的均值和方差矩阵分别如式(25)和(26)所示.
2) $t \in [{t_{sw}}, {t_f}]$, 此时存在模式失配, 系统的状态方程为
$ \begin{equation} {\dot {\pmb x}}(t) = {A\pmb x}(t) + {{\pmb B}_1}u_P(t) + {{\pmb B}_2}{m_2} + {\pmb w}(t) \end{equation} $
(27) 估计器的滤波方程仍为式(20).根据定义, 可以得到状态估计误差的时间传递方程:
$ \begin{equation} \begin{array}{l} {\dot {\tilde {\pmb x}}}(t) = F(t){\tilde {\pmb x}}(t) - {{\pmb B}_2}({m_2} - {m_1}) + {\pmb \zeta}(t)\\ {\tilde {\pmb x}}({t_{sw}}) = {{{\tilde {\pmb x}}}_{sw}} \end{array} \end{equation} $
(28) 求解方程(28)得到
$ \begin{align} {\tilde {\pmb x}}(t) =\,&{\Phi_F}(t, {t_{sw}}){\tilde {\pmb x}}({t_{sw}}) + \int\limits_{{t_{sw}}}^t {{\Phi_F}(t, s){\pmb \zeta }(s){\rm d}s}-\nonumber\\ &({m_2} - {m_1})\int\limits_{{t_{sw}}}^t {{\Phi_F}(t, s){\rm d}s} {{\pmb B}_2} \end{align} $
(29) 利用式(24)可得:
$ \begin{equation} {\tilde {\pmb x}}({t_{sw}}) = {\Phi_F}({t_{sw}}, 0){{\tilde {\pmb x}}_0} + \int\limits_0^{{t_{sw}}} {{\Phi_F}({t_{sw}}, s){\pmb \zeta }(s){\rm d}s} \end{equation} $
(30) 将其代入式(29), 则
$ \begin{align} {\tilde {\pmb x}}(t) =\,&{\Phi_F}(t, 0){{\tilde {\pmb x}}_0} + \int\limits_0^t {{\Phi_F}(t, s){\pmb \zeta }(s){\rm d}s} - \nonumber\\ &({m_2} - {m_1})\int\limits_{{t_{sw}}}^t {{\Phi_F}(t, s){\rm d}s{{\pmb B}_2}} \end{align} $
(31) 根据定义, 状态估计误差的均值为
$ \begin{equation} {\pmb \xi }(t) = {\Phi_F}(t, 0){\rm E}\{ {{\tilde {\pmb x}}_0}\} - ({m_2} - {m_1})\int\limits_{{t_{sw}}}^t {{\Phi_F}(t, s){\rm d}s} {{\pmb B}_2} \end{equation} $
(32) 状态估计误差的方差矩阵同式(26).详细的推导过程见附录A.
情形3 ${t_{sw}} + \Delta t \leq {t_f}$.需要分三段进行讨论:
1) $t \in [0, {t_{sw}})$, 此时系统的状态方程和估计器的滤波方程分别如式(19)和(20)所示, 状态估计误差的均值和方差矩阵则分别见式(25)和(26).
2) $t \in [{t_{sw}}, {t_{sw}} + \Delta t)$, 存在模式失配, 此时系统的状态方程和估计器的滤波方程分别如式(27)和(20)所示, 状态估计误差的均值和方差则矩阵分别见式(32)和(26).
3) $t \in [{t_{sw}} + \Delta t, {t_f}]$, 估计器回到正确的目标模式上, 此时的系统状态方程为式(27), 而估计器的滤波方程为
$ \begin{equation} {\dot{ \hat {\pmb x}}}(t) = {A\hat {\pmb x}}(t) + {{B}_1}u_P(t) + {{B}_2}{m_2} + {k}(t)({\pmb Y}(t) - {H\hat {\pmb x}}(t)) \end{equation} $
(33) 根据定义, 可以得到状态估计误差满足方程
$ \begin{equation} \begin{split} &{\dot {\tilde {\pmb x}}}(t) = {F}(t){\tilde {\pmb x}}(t) + {\pmb \zeta}(t) \\ &{\tilde {\pmb x}}({t_{sw}} + \Delta t) = {{{\tilde {\pmb x}}}_{{t_{sw}} + \Delta t}} \end{split} \end{equation} $
(34) 求解得到
$ \begin{align} {\tilde {\pmb x}}(t) =\,&{\Phi_F}(t, {t_{sw}} + \Delta t){{\tilde {\pmb x}}_{{t_{sw}} + \Delta t}} +\nonumber\\ & \int\limits_{{t_{sw}} + \Delta t}^t {{\Phi_F}(t, s){\pmb \zeta }(s){\rm d}s} \end{align} $
(35) 由式(31)可得:
$ \begin{align} {{{\tilde {\pmb x}}}_{{t_{sw}} + \Delta t}} =\,&{\Phi_F}({t_{sw}} + \Delta t, 0){{{\tilde {\pmb x}}}_0}+\nonumber\\ &\int\limits_0^{{t_{sw}} + \Delta t} {{\Phi_F}({t_{sw}} + \Delta t, s){\pmb \zeta }(s){\rm d}s}- \nonumber\\ &({m_2} - {m_1})\int\limits_{{t_{sw}}}^{{t_{sw}} + \Delta t} {{\Phi_F}({t_{sw}} + \Delta t, s){\rm d}s} {{\pmb B}_2} \end{align} $
(36) 代入式(35)后得到
$ \begin{align} {\tilde {\pmb x}}(t) =\,&{\Phi_F}(t, 0){{\tilde {\pmb x}}_0} + \int\limits_0^t {{\Phi_F}(t, s){\pmb \zeta}(s){\rm d}s} -\nonumber\\ &({m_2} - {m_1})\int\limits_{{t_{sw}}}^{{t_{sw}} + \Delta t} {{\Phi_F}(t, s){\rm d}s} {{\pmb B}_2} \end{align} $
(37) 根据定义, 状态估计误差的均值为
$ \begin{align} {\pmb \xi}(t) =\,&{\Phi_F}(t, 0){\rm E}({{\tilde {\pmb x}}_0}) -\nonumber\\ & ({m_2} - {m_1})\int\limits_{{t_{sw}}}^{{t_{sw}} + \Delta t} {{\Phi_F}(t, s){\rm d}s} {{\pmb B}_2} \end{align} $
(38) 状态估计误差的方差矩阵满足式(26), 详细推导见附录A.
令$\tilde z(t) = \hat z(t) - z(t)$表示ZEM的估计误差, 则
$ \begin{align} \tilde z(t) =\,&{{\pmb g}^{\rm T}}(t){\hat {\pmb x}}(t) - {{\pmb x}^{\rm T}}(t){\pmb x}(t) = {{\pmb g}^{\rm T}}(t){\tilde {\pmb x}}(t)=\nonumber\\ &{\tilde x_1}(t) + {\tilde x_2}(t)({t_f} - t) + {\tilde x_3}(t) \tau _e^2\psi\left( {\frac{{{t_f} - t}}{{{\tau _e}}}}\right)-\nonumber\\ &{\tilde x_4}(t)\tau _p^2\psi\left( {\frac{{{t_f} - t}}{{{\tau _p}}}}\right) \end{align} $
(39) 其中, $\psi (\theta ) = {\text{exp}}( - \theta ) + \theta - 1$.
令${\mu (t) = {\rm E}\{ \tilde z(t)\} }$, ${{\sigma ^2}(t) = {\mathop{\rm var}} \{ \tilde z(t)\} }$, 则由式(39)可以得到:
$ \begin{equation} \mu (t) = {{\pmb g}^{\rm T}}(t){\rm E}\{{\tilde {\pmb x}}(t)\} ={{\pmb g}^{\rm T}}(t){\pmb \xi }(t) \end{equation} $
(40) $ \begin{align} {\sigma ^2}(t) =\, &{\mathop{\rm var}} \{ \tilde z(t)\} = {\rm E}\{ {\{ \tilde z(t) - {\rm E}[\tilde z(t)]\} ^2}\}=\nonumber\\ &{\rm E}\{ {\{ {{\pmb g}^{\rm T}}(t){\tilde {\pmb x}}(t) - {{\pmb g}^{\rm T}}(t){\rm E}[{\tilde {\pmb x}}(t)]\} ^2}\}=\nonumber\\ &{\rm E}\{ {\{ {{\pmb g}^{\rm T}}(t)\{ {\tilde {\pmb x}}(t) - {\rm E}[{\tilde {\pmb x}}(t)]\} \} ^2}\}=\nonumber\\ &{{\pmb g}^{\rm T}}(t){\rm E}\{ \{ {\tilde {\pmb x}}(t) -\nonumber\\ &{\rm E}[{\tilde {\pmb x}}(t)]\} \cdot {\{ {\tilde {\pmb x}}(t) - {\rm E}[{\tilde {\pmb x}}(t)]\} ^{\rm T}}\} {\pmb g}(t)=\nonumber\\ &{{\pmb g}^{\rm T}}(t){\mathop{\rm{var}}} \{ {\tilde {\pmb x}}(t)\}{\pmb g}(t)=\nonumber\\ &{{\pmb g}^{\rm T}}(t){\Sigma }(t){\pmb g}(t) \end{align} $
(41) 因此, 存在模式失配时每一时刻ZEM的估计误差均服从有偏的高斯分布, 其均值和方差分别为$\mu (t)$和${\sigma ^2}(t)$.从式(40)和(41)的结果来看:
1) 制导系统对ZEM估计误差的的影响主要体现在导引律、弹目时间常数、剩余飞行时间${t_{go}}$以及观测精度4个方面; 目标的影响则主要体现由机动导致的模式失配上.
2) 模式失配只影响ZEM估计误差的均值$\mu (t)$, 对${\sigma ^2}(t)$没影响.
3) $\mu (t)$和${\sigma ^2}(t)$与系统所用的导引律无关, 在后面的仿真验证中, 不失一般性, 导引律选用DGL/1.
4) ${t_{go}}$的估计精度及弹目时间常数${\tau _e}$和${\tau _p}$通过投影向量${{\pmb g}}(t)$影响脱靶量.在具体的拦截问题中, 弹目时间常数通常可假定为确定已知的, 而雷达导引头可直接获得高精度的$t_{go}$测量, 因此本文分析中不考虑它们对ZEM估计误差的影响.
5) 估计器的观测精度将直接影响到Kalman增益系数$k(t)$, 见式(21), 进而影响ZEM估计误差的均值和方差.
3. 仿真实验
本节通过一个典型的TBM拦截场景验证前面理论推导的正确性, 仿真参数设置如表 1, 蒙特卡洛仿真次数设置为1 000.
表 1 仿真参数Table 1 Simulation parameters参数类型 参数名称 单位 值(范围) 弹目参数 ${V_p}$ m/s 2 300 ${V_e}$ m/s 2 700 $a_p^{\max}$ g 30 $a_e^{\max}$ g 12, 15 ${\tau _p}$ s 0.2 ${\tau _e}$ s 0.2 观测参数 T s 0.01 ${\sigma _\theta }$ mrad 5 ${\sigma _a }$ $\rm m/{\rm s}^2$ 1 场景参数 $r_0$ m 15 000 ${\phi_p}(0)$ rad 均匀分布($-\pi /18, \pi /18$) ${\phi_e}(0)$ rad $ > \pi /2$且满足碰撞三角形 ${u _p}(0)$ g 0 ${a_p}(0)$ g 0 ${a_e}(0)$ g $a_e^{\max}$ 目标机动方式 - 随机乒乓 估计器参数 $s_w$ $\rm g^2\rm{/Hz}$ 1 $t_{sw}$ s 2 $\Delta t$ s 0.1 初始状态 - ${{\hat {\pmb x}}_0} = {[0, 0, 0, 0]^{\rm T}}$ 初始估计误差 - ${{\tilde {\pmb x}}_0} = {[0, 0, a_e^{\max}, 0]^{\rm T}}$ 初始估计协方差 - ${{{\tilde P}}_0} = \left[{{array}{*{20}{c}} 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0\\ 0 & 0 &{{{(a_e^{\max})}^2}}& 0\\ 0 & 0 & 0 & 0 {array}} \right]$ 图 5给出了两种不同过载比$\gamma=2$和$\gamma=2.5$下ZEM估计误差的均值变化曲线.从实验结果可以看出, 本文推导的理论结果与蒙特卡洛仿真的曲线基本吻合.由图 5还可以看出, 当目标的运动模式改变时($t = 2$ s), ZEM的估计误差会迅速增大, 当运动模式被正确识别后($t = 2.1$ s), ZEM的估计误差将逐渐减小. 图 6给出了这两种情形下ZEM估计误差的方差分布, 可以看出本文理论推导结果与蒙特卡洛仿真的曲线同样也是吻合的. 图 5和图 6的仿真结果充分说明了本文理论推导的正确性.
图 7给出了各状态分量的估计误差.由该图可见, 当目标运动模式改变后, 各状态分量的估计误差都迅速增大, 当模式匹配后, 估计误差逐渐减小; 导弹自身加速度估计误差分量不受模式失配的影响, 这与导弹自身的加速度模型是完美可知且可精确测量的假设相一致.
4. 结论
本文针对高速大机动目标拦截问题, 推导了模式失配条件下ZEM估计误差的分布形式.在过程噪声和测量噪声均为零均值高斯白噪声且相互独立的假定下, 每一时刻ZEM的估计误差服从有偏的高斯分布.本文得到了各时刻ZEM估计误差均值和方差的解析表达式, 并与蒙特卡洛仿真实验进行了对比, 验证了理论推导的正确性.
将本文的ZEM估计误差模型应用于模式失配条件下脱靶量模型的推导, 以及研究基于特征辅助的目标运动模式辨识算法将是下一步工作的方向.
附录A
情形1中$\pmb \xi(t)$和$\Sigma(t)$形式证明:
因为
$ \begin{align} {\rm E}\{ {\pmb \zeta}(t)\} = \, &{\rm E}\{ {k}(t){\pmb v}(t) - {\pmb w}(t)\} =\nonumber\\ &{k}(t){\rm E}\{ {\pmb v}(t)\} - {\rm E}\{ {\pmb w}(t)\} = {\pmb 0} \end{align} $
(A1) 所以
$ \begin{equation} \begin{split} {\rm E}\{ {\tilde {\pmb x}}(t)\} &= {\Phi_F}(t, 0){\rm E}\{ {{{\tilde {\pmb x}}}_0}\} + \int\limits_0^t {{\Phi_F}(t, s){\rm E}\{ {\pmb \zeta}(s)\} {\rm d}s} = \\ &{\Phi_F}(t, 0){\rm E} ({{\tilde {\pmb x}}}_0) \end{split} \end{equation} $
(A2) 同理, 很容易推出式(32)和(38), 这里不再赘述.
由${\pmb w}(t)$与${\pmb v}(t)$相互独立, 所以
$ \begin{align} {\rm E}\{ {\pmb \zeta }(t){{\pmb \zeta}^{\rm T}}(t)\} = \, &{\rm E}\{ [{k}(t){\pmb v}(t)- \nonumber\\ & {\pmb w}(t)]{[{k}(t){\pmb v}(t)-{\pmb w}(t)]^{\rm T}}\} =\nonumber\\ &{k}(t){R}(t){{k}^{\rm T}}(t) + {Q} \end{align} $
(A3) 根据定义且由假设条件${\pmb w}(t)$、${\pmb v}(t)$和${{\tilde {\pmb x}}_0}$之间互相独立
$ \begin{align} &{\Sigma }(t) = {\rm E}\{ \tilde {\pmb x}(t) - {\rm E}\{ \tilde {\pmb x}(t)\} \} ^2={\rm E}\Bigg\{ {\Phi_F}(t, 0){{{\tilde {\pmb x}}}_0} +\nonumber\\ & \int\limits_0^t {{\Phi_F}(t, s){\pmb \zeta}(s){\rm d}s} - {\Phi_F}(t, 0){\rm E}\{ {{{\tilde {\pmb x}}}_0}\} \Bigg\}^2=\nonumber\\ & {\rm E}{\left\{ {{\Phi_F}(t, 0)[{{{\tilde {\pmb x}}}_0}- {\rm E}({{{\tilde {\pmb x}}}_0})] + \int\limits_0^t {{\Phi_F}(t, s){\pmb \zeta }(s){\rm d}s} } \right\}^2}=\nonumber\\ & {\Phi_F}(t, 0){\rm E}\left\{ {[{{{\tilde {\pmb x}}}_0}-{\rm E}({{{\tilde {\pmb x}}}_0})]{{[{{{\tilde {\pmb x}}}_0}- {\rm E}({{{\tilde {\pmb x}}}_0})]}^{\rm T}}} \right\}{\Phi_F}^{\rm T}{(t, 0)}+\nonumber\\ &\int\limits_0^t {{\Phi_F}(t, s){\rm E}\{ {\pmb \zeta } (s){{\pmb \zeta }^{\rm T}}(s)\} {\Phi_F}^{\rm T}{{(t, s)}}{\rm d}s}=\nonumber\\ & {\Phi_F}(t, 0){{{\tilde P}}_0}{\Phi_F}^{\rm T}{(t, 0)}+\nonumber\\ &\int\limits_0^t {{\Phi_F}(t, s)[{Q} + {k}(s){R}(s){{k}^{\rm T}}(s)]{\Phi_F}^{\rm T}{{(t, s)}}{\rm d}s} \end{align} $
(A4) 同理, 很容易得出在情形2和情形3下, $\Sigma(t)$的表达式与式(26)相同, 这里不再赘述.
附录B
符号说明 $P$ 导弹 $E$ 目标 $\tau_p$, $\tau_e$ 导弹和目标控制系统的时间常数 $a_p^{\max}, a_e^{\max}$ 导弹和目标最大横向加速度 ${V_p}, {V_e}$ 导弹和目标的飞行速度 ${u_p}, {u_e}$ 导弹和目标的横向加速度指令 $r$ 弹目相对距离 ${t_{sw}}$ 目标模式切换时刻 $t$ 仿真时间 ${t_f}$ 终止时刻 g 重力加速度, $9.8\rm m/{\rm{s}^2}$ $m$ 目标的运动模式 ${m_1}, {m_2}$ 目标在模式切换时刻前后的运动模式 $\Delta m$ 目标运动模式改变量, $\Delta m = {m_2} - {m_1}$ $T$ 离散采样时间间隔 ${\sigma _\theta }$ 测角精度 ${\sigma _a}$ 导弹加速度测量精度 ${s_w}$ 目标指令加速度误差的功率谱密度 $\Delta t$ 目标运动模式辨识延迟 ${\tilde {\pmb x}}$ 状态估计误差 ${\pmb \xi }, {\Sigma}$ 状态估计误差的均值和方差 $\eta (t)$ ZEM估计误差 $\mu, {\sigma ^2}$ ZEM估计误差的均值和方差
-
表 1 仿真参数
Table 1 Simulation parameters
参数类型 参数名称 单位 值(范围) 弹目参数 ${V_p}$ m/s 2 300 ${V_e}$ m/s 2 700 $a_p^{\max}$ g 30 $a_e^{\max}$ g 12, 15 ${\tau _p}$ s 0.2 ${\tau _e}$ s 0.2 观测参数 T s 0.01 ${\sigma _\theta }$ mrad 5 ${\sigma _a }$ $\rm m/{\rm s}^2$ 1 场景参数 $r_0$ m 15 000 ${\phi_p}(0)$ rad 均匀分布($-\pi /18, \pi /18$) ${\phi_e}(0)$ rad $ > \pi /2$且满足碰撞三角形 ${u _p}(0)$ g 0 ${a_p}(0)$ g 0 ${a_e}(0)$ g $a_e^{\max}$ 目标机动方式 - 随机乒乓 估计器参数 $s_w$ $\rm g^2\rm{/Hz}$ 1 $t_{sw}$ s 2 $\Delta t$ s 0.1 初始状态 - ${{\hat {\pmb x}}_0} = {[0, 0, 0, 0]^{\rm T}}$ 初始估计误差 - ${{\tilde {\pmb x}}_0} = {[0, 0, a_e^{\max}, 0]^{\rm T}}$ 初始估计协方差 - ${{{\tilde P}}_0} = \left[{\begin{array}{*{20}{c}} 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0\\ 0 & 0 &{{{(a_e^{\max})}^2}}& 0\\ 0 & 0 & 0 & 0 \end{array}} \right]$ 
下载: 导出CSV
符号说明 $P$ 导弹 $E$ 目标 $\tau_p$, $\tau_e$ 导弹和目标控制系统的时间常数 $a_p^{\max}, a_e^{\max}$ 导弹和目标最大横向加速度 ${V_p}, {V_e}$ 导弹和目标的飞行速度 ${u_p}, {u_e}$ 导弹和目标的横向加速度指令 $r$ 弹目相对距离 ${t_{sw}}$ 目标模式切换时刻 $t$ 仿真时间 ${t_f}$ 终止时刻 g 重力加速度, $9.8\rm m/{\rm{s}^2}$ $m$ 目标的运动模式 ${m_1}, {m_2}$ 目标在模式切换时刻前后的运动模式 $\Delta m$ 目标运动模式改变量, $\Delta m = {m_2} - {m_1}$ $T$ 离散采样时间间隔 ${\sigma _\theta }$ 测角精度 ${\sigma _a}$ 导弹加速度测量精度 ${s_w}$ 目标指令加速度误差的功率谱密度 $\Delta t$ 目标运动模式辨识延迟 ${\tilde {\pmb x}}$ 状态估计误差 ${\pmb \xi }, {\Sigma}$ 状态估计误差的均值和方差 $\eta (t)$ ZEM估计误差 $\mu, {\sigma ^2}$ ZEM估计误差的均值和方差
下载: 导出CSV
-
[1] Shinar J, Turetsky V. Three-dimensional validation of an integrated estimation/guidance algorithm against randomly maneuvering targets. Journal of Guidance, Control, and Dynamics, 2014, 32(3):1034-1039 doi: 10.2514/1.40542 [2] Shinar J, Shima T. Nonorthodox guidance law development approach for intercepting maneuvering targets. Journal of Guidance, Control, and Dynamics, 2002, 25(4):658-666 doi: 10.2514/2.4960 [3] Shinar J, Turetsky V. Meeting the challenges of modern interceptor guidance by non-conventional approaches. In: Proceedings of the 17th Mediterranean Conference on Control and Automation. Thessaloniki, Greece: IEEE, 2009. 1563-1568 http://ieeexplore.ieee.org/document/5164770/ [4] Shinar J, Oshman Y, Turetsky V. On the need for integrated estimation/guidance design for hit-to-kill accuracy. In: Proceedings of the 2003 American Control Conference. Denver, CO, USA: IEEE, 2003, 1: 402-407 [5] Witsenhausen H S. Separation of estimation and control for discrete time system. Proceedings of the IEEE, 1971, 59(11):1557-1566 doi: 10.1109/PROC.1971.8488 [6] Alexandre P. Multiple Model Estimation and Detection for Adaptive Guidance of Hybrid Systems [Master dissertation], McGill University, Canada, 2004. [7] Glizer V Y, Turetsky V, Shinar J. Terminal cost distribution in discrete-time controlled system with disturbance and noise-corrupted state information. International Journal of Applied Mathematics, 2012, 42(1):52-59 http://d.old.wanfangdata.com.cn/OAPaper/oai_doaj-articles_9728187a2dae195516565edae35386af [8] Shinar J, Glizer V Y, Turetsky V. Distribution of terminal cost functional in continuous-time controlled system with noise-corrupted state information. In: Proceedings of the 27th Convention of Electrical & Electronics Engineers in Israel. Eilat, Israel: IEEE, 2011. 1-5 http://ieeexplore.ieee.org/xpls/icp.jsp?arnumber=6377008 [9] Shinar J, Glizer V Y, Turetsky V. Terminal state distribution of continuous-time system with random disturbance and noise-corrupted information. International Journal of Applied Mathematics, 2015, 45(2):77-84 [10] Glizer V Y, Turetsky V, Shinar J. Distribution of terminal cost functional in discrete-time controlled system with noise-corrupted state information. In: Proceedings of the 2012 IEEE 27th Convention of Electrical & Electronics Engineers in Israel (IEEEI). London, UK: IEEE, 2012. [11] Moldavskaya E, Shinar J. Distribution of the zero-effort miss distance estimation error in interception problems. Applied Mathematics and Computation, 2015, 269:217-231 doi: 10.1016/j.amc.2015.07.073 [12] Shinar J, Turetsky V, Oshman Y. New logic-based estimation/guidance algorithm for improved homing against randomly maneuvering targets. In: Proceedings of the 2014 AIAA Guidance, Navigation and Control Conference. Providence, Rhode Island: AIAA, 2004. 1-17 doi: 10.2514/6.2004-4886 [13] Dionne D, Michalska H, Shinar J, et al. Decision-directed adaptive estimation and guidance for an interception endgame. Journal of Guidance, Control, and Dynamics, 2015, 29(4):970-980 doi: 10.2514/1.16050 [14] Shinar J, Turetsky V, Oshman Y. Integrated estimation/guidance design approach for improved homing against randomly maneuvering targets. Journal of Guidance, Control, and Dynamics, 2007, 30(1):154-161 doi: 10.2514/1.22916 [15] 范红旗.主动寻的制导中机动目标运动模式辨识技术博士学位论文, 国防科学技术大学, 中国, 2008.Fan Hong-Qi. Technology on Maneuvering Target Motion Mode Identification in Active Homing Guidance [Ph. D. dissertation], National University of Defense Technology, China, 2008. [16] Author(s) Zhu Y L, Fan H Q, Fan J P, Lu Z Q, Fu Q. Target turning maneuver detection using high resolution doppler profile. IEEE Transactions on Aerospace and Electronic Systems, 2012, 48(1): 762-779 http://ieeexplore.ieee.org/xpl/articleDetails.jsp?arnumber=6129669 [17] Bizup D F, Brown D E. Maneuver detection using the radar range rate measurement. IEEE Transactions on Aerospace and Electronic Systems, 2004, 40(1):330-336 doi: 10.1109/TAES.2004.1292169 [18] Hughes E J, Leyland M. Target manoeuvre detection using radar glint. Electronics Letters, 1998, 34(17):1695-1696 doi: 10.1049/el:19981188 [19] Fan H Q, Zhu Y L, Fu Q. Impact of mode decision delay on estimation error for maneuvering target interception. IEEE Transactions on Aerospace and Electronic Systems, 2011, 47(1):702-711 doi: 10.1109/TAES.2011.5705700 [20] Turetsky V, Shinar J. Missile guidance laws based on pursuit-evasion game formulations. Automatica, 2003, 39(4):607-618 doi: 10.1016/S0005-1098(02)00273-X [21] 李运迁, 齐乃明.基于零控脱靶量的大气层内拦截弹制导律.宇航学报, 2010, 31(7):1768-1774 doi: 10.3873/j.issn.1000-1328.2010.07.011Li Yun-Qian, Qi Nai-Ming. A zero-effort miss distance-based guidance law for endoatmoshperic interceptor. Journal of Astronautics, 2010, 31(7):1768-1774 doi: 10.3873/j.issn.1000-1328.2010.07.011 [22] Blair W D, Watson G A, Kirubarajan T, Bar-Shalom Y. Benchmark for radar allocation and tracking in ECM. IEEE Transactions on Aerospace and Electronic Systems, 1998, 34(4):1097-1114 doi: 10.1109/7.722694 [23] Li X R, He C. Model-set design, choice, and comparison for multiple-model estimation. In: Proceedings of SPIE 3809 International Symposium on Optical Science, Engineering, and Instrumentation. Denver, CO, USA: SPIE, 1999. 501-513 http://spie.org/x648.xml?product_id=364047 [24] Baram Y, Sandell N. An information theoretic approach to dynamical systems modeling and identification. IEEE Transactions on Automatic Control, 1978, 23(1):61-66 doi: 10.1109/TAC.1978.1101690 期刊类型引用(1)
1. Shengwen Xiang，Hongqi Fan，Qiang Fu. Distribution of Miss Distance for Pursuit-Evasion Problem. IEEE/CAA Journal of Automatica Sinica. 2020(04): 1161-1168 . 
必应学术其他类型引用(1)
-
下载:
计量
- 文章访问数: 2168
- HTML全文浏览量: 250
- PDF下载量: 836
- 被引次数: 2
