CN105446352B - A kind of proportional navigation law recognizes filtering method - Google Patents

A kind of proportional navigation law recognizes filtering method Download PDF

Info

Publication number
CN105446352B
CN105446352B CN201510822999.2A CN201510822999A CN105446352B CN 105446352 B CN105446352 B CN 105446352B CN 201510822999 A CN201510822999 A CN 201510822999A CN 105446352 B CN105446352 B CN 105446352B
Authority
CN
China
Prior art keywords
mrow
mtd
msub
msubsup
mtr
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201510822999.2A
Other languages
Chinese (zh)
Other versions
CN105446352A (en
Inventor
邹昕光
周荻
朱蕊蘋
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin Institute of Technology Shenzhen
Original Assignee
Harbin Institute of Technology Shenzhen
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Institute of Technology Shenzhen filed Critical Harbin Institute of Technology Shenzhen
Priority to CN201510822999.2A priority Critical patent/CN105446352B/en
Publication of CN105446352A publication Critical patent/CN105446352A/en
Application granted granted Critical
Publication of CN105446352B publication Critical patent/CN105446352B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G05CONTROLLING; REGULATING
    • G05DSYSTEMS FOR CONTROLLING OR REGULATING NON-ELECTRIC VARIABLES
    • G05D1/00Control of position, course, altitude or attitude of land, water, air or space vehicles, e.g. using automatic pilots
    • G05D1/10Simultaneous control of position or course in three dimensions
    • G05D1/107Simultaneous control of position or course in three dimensions specially adapted for missiles

Landscapes

  • Engineering & Computer Science (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Automation & Control Theory (AREA)
  • Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)
  • Other Investigation Or Analysis Of Materials By Electrical Means (AREA)

Abstract

一种比例导引制导律辨识滤波方法,本发明涉及比例导引制导律辨识滤波方法。本发明是为了解决现有MMAE滤波器假定PN制导律导航常数已知和未考虑pursuer控制器饱和情况而导致模型不准确,估计精度低的问题。本发明首先建立PN制导律和饱和运动模型在俯仰和偏航平面的状态方程;接着计算系统的状态转移矩阵和设计Kalman滤波器方程;然后计算pursuer当前时刻PN制导律运动模型和饱和运动模型的后验概率;最后根据后验概率计算pursuer运动模型切换的时刻,在此时刻前采用饱和运动模型Kalman滤波器方程的估计结果否则采用PN制导律运动模型Kalman滤波器方程的估计结果。本发明应用于航天领域。

A proportional guidance guidance law identification filtering method, the invention relates to a proportional guidance guidance law identification filtering method. The invention aims to solve the problem that the existing MMAE filter assumes that the PN guidance law navigation constant is known and does not consider the saturation of the pursuer controller, resulting in inaccurate models and low estimation precision. The present invention first establishes the state equation of the PN guidance law and the saturation motion model in the pitch and yaw planes; then calculates the state transition matrix of the system and designs the Kalman filter equation; then calculates the PN guidance law motion model and the saturation motion model of the pursuer at the current moment Posterior probability; finally, calculate the moment when the pursuer motion model is switched according to the posterior probability. Before this time, the estimated result of the Kalman filter equation of the saturated motion model is used; otherwise, the estimated result of the Kalman filter equation of the PN guidance law motion model is used. The invention is applied in the aerospace field.

Description

一种比例导引制导律辨识滤波方法A Guidance Law Identification and Filtering Method for Proportional Guidance

技术领域technical field

本发明涉及比例导引制导律辨识滤波方法。The invention relates to a proportional guidance guidance law identification and filtering method.

背景技术Background technique

在大气层外弹道导弹突防应用中,为了提高弹道的突防能力,弹道导弹可以在大气层外释放防御导弹与其同速伴飞。当防御导弹发现拦截导弹时通过制导律控制与其碰撞,从而提高弹道导弹的突防概率。弹道导弹的最优躲避策略和防御导弹的最优制导律的设计均假定拦截导弹的制导律已知。虽然基于滑模控制制导律可以将拦截导弹的机动当做未知的外部扰动来处理,但如果能估计出拦截导弹的加速度,可以大幅提高滑模制导律的性能。因此研究拦截导弹制导律的辨识问题具有实际的意义。In the application of extra-atmospheric ballistic missile penetration, in order to improve the ballistic penetration capability, ballistic missiles can release defensive missiles outside the atmosphere to accompany them at the same speed. When the defensive missile finds the intercepting missile, the collision with it is controlled by the guidance law, thereby increasing the penetration probability of the ballistic missile. The optimal evasion strategy for ballistic missiles and the optimal guidance law for defensive missiles both assume that the guidance law for intercepting missiles is known. Although the maneuvering of the interceptor missile can be treated as an unknown external disturbance based on the sliding mode control guidance law, if the acceleration of the interceptor missile can be estimated, the performance of the sliding mode guidance law can be greatly improved. Therefore, it is of practical significance to study the identification of the guidance law of intercepting missiles.

Shaferman等人在《JOURNAL OF GUIDANCE CONTROL AND DYNAMICS》杂志2010年第6期的文章“Cooperative Multiple-Model Adaptive Guidance for an AircraftDefending Missile”中设计了一种基于MMAE(Multiple-Model Adaptive Estimator)的飞机主动防御协同制导律。假设拦截导弹采用比例导引(PN)制导律,增强比例导引(APN)制导律或者最优制导律(OGL)中的一种,采用多模型自适应滤波(MMAE)对拦截导弹的制导律进行辨识。但该工作存在两个问题。1)该方法假定拦截导弹采用了已知的制导律常数N。这限制了该MMAE的适用范围。在实际场景中,拦截导弹制导律的导航常数是未知的。2)该方法没有考虑拦截导弹执行器饱和的情况。一般情况下,在制导阶段初期,拦截导弹会由于快速补偿初始的对齐误差而暂时进入执行器饱和状态,在该误差被补偿后控制器退出饱和状态。不考虑执行器饱和的情况,将使得MMAE滤波器模型不准确,估计精度低。In the article "Cooperative Multiple-Model Adaptive Guidance for an AircraftDefending Missile" in "JOURNAL OF GUIDANCE CONTROL AND DYNAMICS" magazine No. 6, 2010, Shaferman et al. designed an active aircraft defense based on MMAE (Multiple-Model Adaptive Estimator) Coordinated Guidance Law. Assuming that the intercepting missile adopts one of the proportional guidance (PN) guidance law, the enhanced proportional guidance (APN) guidance law or the optimal guidance law (OGL), the guidance law of the intercepting missile is controlled by the multi-model adaptive filtering (MMAE) to identify. But there are two problems with this work. 1) This method assumes that the interceptor missile uses a known guidance law constant N. This limits the applicability of this MMAE. In practical scenarios, the navigation constants of the interceptor missile guidance law are unknown. 2) This method does not consider the saturation of interceptor missile actuators. Generally, at the beginning of the guidance phase, the interceptor missile will temporarily enter the actuator saturation state due to the rapid compensation of the initial alignment error, and the controller will exit the saturation state after the error is compensated. If the actuator saturation is not considered, the MMAE filter model will be inaccurate and the estimation accuracy will be low.

发明内容Contents of the invention

本发明的目的是为了解决以下问题:1)现有方法假定拦截导弹采用了已知的制导律常数N,限制了MMAE的适用范围。而在实际场景中,拦截导弹制导律的导航常数是未知的;2)拦截导弹执行器在制导初始时刻可能出现饱和情况。现有方法没有考虑以上实际情况,使得MMAE滤波器模型不准确,估计精度低的问题,针对以上问题,本发明提出了一种比例导引制导律辨识滤波方法。The purpose of the present invention is to solve the following problems: 1) the existing method assumes that the intercepting missile has adopted a known guidance law constant N, which limits the scope of application of MMAE. However, in actual scenarios, the navigation constant of the guidance law of the interceptor missile is unknown; 2) The actuator of the interceptor missile may be saturated at the initial moment of guidance. The existing method does not consider the above actual situation, which makes the MMAE filter model inaccurate and the estimation accuracy low. To solve the above problems, the present invention proposes a proportional guidance guidance law identification filtering method.

上述的发明目的是通过以下技术方案实现的:Above-mentioned purpose of the invention is achieved through the following technical solutions:

步骤一:pursuer有两种运动模型:控制器饱和时的饱和运动模型和控制器未饱和时的PN制导律运动模型;建立PN制导律运动模型在俯仰平面的状态方程、PN制导律运动模型在偏航平面的状态方程、控制器饱和运动模型在俯仰平面的状态方程和控制器饱和运动模型在偏航平面的状态方程;evader为弹道导弹;pursuer为拦截导弹;PN制导律为比例导引制导律;Step 1: The pursuer has two motion models: the saturated motion model when the controller is saturated and the PN guidance law motion model when the controller is not saturated; establish the state equation of the PN guidance law motion model in the pitch plane, and the PN guidance law motion model in the The state equation of the yaw plane, the state equation of the controller saturation motion model in the pitch plane, and the state equation of the controller saturation motion model in the yaw plane; evader is a ballistic missile; pursuer is an intercept missile; PN guidance law is proportional guidance guidance law;

步骤二:根据步骤一的PN制导律运动模型在俯仰平面的状态方程、PN制导律运动模型在偏航平面的状态方程、控制器饱和运动模型在俯仰平面的状态方程和控制器饱和运动模型在偏航平面的状态方程,设计PN制导律运动模型在俯仰平面的状态转移矩阵、PN制导律运动模型在偏航平面的状态转移矩阵、控制器饱和运动模型在俯仰平面的状态转移矩阵、控制器饱和运动模型在偏航平面的状态转移矩阵;Step 2: According to the state equation of the PN guidance law motion model in the pitch plane in step 1, the state equation of the PN guidance law motion model in the yaw plane, the state equation of the controller saturation motion model in the pitch plane, and the controller saturation motion model in The state equation of the yaw plane, design the state transition matrix of the PN guidance law motion model in the pitch plane, the state transition matrix of the PN guidance law motion model in the yaw plane, the state transition matrix of the controller saturation motion model in the pitch plane, and the controller The state transition matrix of the saturated motion model in the yaw plane;

步骤三:根据步骤二的PN制导律运动模型在俯仰平面的状态转移矩阵、PN制导律运动模型在偏航平面的状态转移矩阵、控制器饱和运动模型在俯仰平面的状态转移矩阵、控制器饱和运动模型在偏航平面的状态转移矩阵,设计PN制导律运动模型在俯仰平面的Kalman滤波器方程、PN制导律运动模型在偏航平面的Kalman滤波器方程、控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程,Kalman滤波器方程为卡尔曼滤波器方程;Step 3: According to the state transition matrix of the PN guidance law motion model in the pitch plane in step 2, the state transition matrix of the PN guidance law motion model in the yaw plane, the state transition matrix of the controller saturation motion model in the pitch plane, and the controller saturation The state transition matrix of the motion model in the yaw plane, design the Kalman filter equation of the PN guidance law motion model in the pitch plane, the Kalman filter equation of the PN guidance law motion model in the yaw plane, and the controller saturation motion model in the pitch plane The Kalman filter equation and the Kalman filter equation of the controller saturation motion model in the yaw plane, and the Kalman filter equation is the Kalman filter equation;

步骤四:根据步骤三的PN制导律运动模型在俯仰平面的Kalman滤波器方程、PN制导律运动模型在偏航平面的Kalman滤波器方程、控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程,计算pursuer当前时刻PN制导律运动模型的后验概率和控制器饱和运动模型的后验概率;Step 4: According to the Kalman filter equation of the PN guidance law motion model in step 3 in the pitch plane, the Kalman filter equation of the PN guidance law motion model in the yaw plane, the Kalman filter equation of the controller saturation motion model in the pitch plane, and The Kalman filter equation of the controller saturated motion model in the yaw plane, calculate the posterior probability of the PN guidance law motion model and the posterior probability of the controller saturated motion model at the current moment of the pursuer;

步骤五:根据步骤四得到的pursuer当前时刻PN制导律运动模型的后验概率和饱和运动模型的后验概率,计算pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻;在pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻前采用控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程的估计结果;否则采用PN制导律运动模型在俯仰平面的Kalman滤波器方程和PN制导律运动模型在偏航平面的Kalman滤波器方程的估计结果。Step 5: According to the posterior probability of the PN guidance law motion model and the posterior probability of the saturated motion model obtained in step 4, calculate the PN guidance when the pursuer switches from the saturated motion model when the controller is saturated to the unsaturated controller The moment of the law motion model; before the moment when the pursuer switches from the saturated motion model when the controller is saturated to the PN guidance law motion model when the controller is not saturated, use the Kalman filter equation of the controller saturated motion model in the pitch plane and the controller The estimation result of the Kalman filter equation of the saturated motion model in the yaw plane; otherwise, the Kalman filter equation of the PN guidance law motion model in the pitch plane and the estimation result of the Kalman filter equation of the PN guidance law motion model in the yaw plane.

发明效果Invention effect

本发明提出的比例导引制导律辨识滤波器采用了多模型滤波方法。设计了分别针对PN制导律运动模型和饱和运动模型的扩展卡尔曼滤波器。其中PN制导律运动模型将PN制导律的导航常数作为待估计的状态,从而解决了PN制导律导航常数未知的实际情况。上述两类扩展卡尔曼滤波器并行运行,通过贝叶斯推理框架计算两个模型的概率,从而辨识出拦截导弹采用的运动模型,进而估计出拦截导弹的加速度和PN制导律导航常数。这解决了拦截导弹在制导初期其控制器可能饱和的问题。本发明提出的滤波方法相对于传统的只采用PN制导律模型的扩展卡尔曼滤波方法其估计精度提高了40-50%左右。The proportional guidance guidance law identification filter proposed by the invention adopts a multi-model filtering method. Extended Kalman filters are designed for PN guidance law motion model and saturated motion model respectively. The motion model of the PN guidance law takes the navigation constant of the PN guidance law as the state to be estimated, thus solving the actual situation that the navigation constant of the PN guidance law is unknown. The above two types of extended Kalman filters run in parallel, and calculate the probability of the two models through the Bayesian inference framework, thereby identifying the motion model adopted by the interceptor missile, and then estimating the acceleration and PN guidance law navigation constant of the interceptor missile. This solves the problem that the controller of the interceptor missile may saturate during the initial period of guidance. Compared with the traditional extended Kalman filtering method using only PN guidance law model, the filtering method proposed by the invention improves the estimation accuracy by about 40-50%.

考虑了在制导过程中pursuer控制器在初期产生饱和的可能性,使用了多模型滤波方法对pursuer制导律进行了辨识,对PN制导律下的导航常数和在俯仰平面和偏航平面加速度进行了估计。能更准确的对制导律进行辨识。在本发明的实施例中,设MMAE滤波器对pursuer加速度估计产生的均方误差为VarMMAE;全程使用PN制导律对pursuer加速度估计产生的均方误差为VarPN。设α=VarMMAE/VarPN。当0<α<1,MMAE滤波器对pursuer加速度估计较说明全程使用PN制导律对pursuer加速度估计更精确;α越小,说明相对于全程使用PN制导律对pursuer加速度估计来说,MMAE对pursuer加速度估计越精确。Considering the possibility of saturation of the pursuer controller at the initial stage in the guidance process, the multi-model filtering method is used to identify the pursuer guidance law, and the navigation constants and accelerations in the pitch plane and yaw plane under the PN guidance law are calculated. estimate. The guidance law can be identified more accurately. In an embodiment of the present invention, it is assumed that the mean square error generated by the MMAE filter for the estimation of the pursuer acceleration is Var MMAE ; the mean square error generated by the PN guidance law for the estimation of the pursuer acceleration throughout the whole process is Var PN . Let α=Var MMAE /Var PN . When 0<α<1, the MMAE filter is more accurate in estimating the pursuit acceleration than that using the PN guidance law in the whole process; the smaller the α, it means that compared with using the PN guidance law in the whole process to estimate the pursuer acceleration, the MMAE is more accurate in the pursuit acceleration estimation. Acceleration estimates are more accurate.

当evader无机动时,在pursuer饱和阶段,对pursuer加速度在Evader惯性坐标系y轴的分量ay的估计时α=0.59;对pursuer加速度在Evader惯性坐标系z轴的分量az的估计时α=0.56。When the evader is not maneuvering, in the saturation stage of the pursuer, the estimation of the component a y of the pursuer acceleration on the y-axis of the Evader inertial coordinate system is α=0.59; the estimation time of the component a z of the pursuer acceleration on the z-axis of the Evader inertial coordinate system is α = 0.56.

当evader进行常值机动时,在pursuer饱和阶段,对pursuer加速度ay的估计时α=0.55;对pursuer加速度az的估计时α=0.62。When the evader performs a constant maneuver, in the saturation stage of the pursuer, the estimation of the pursuit acceleration a y is α=0.55; the estimation of the pursuit acceleration a z is α=0.62.

当evader进行正弦机动时,在pursuer饱和阶段,对pursuer加速度ay的估计时α=0.56;对pursuer加速度az的估计时α=0.60。When the evader performs a sinusoidal maneuver, in the saturation stage of the pursuer, the estimation of the pursuit acceleration a y is α=0.56; the estimation of the pursuit acceleration a z is α=0.60.

附图说明Description of drawings

图1为pursuer和evader相对运动关系示意图,E为evader在evader惯性系下的坐标,P为pursuer在evader惯性系下的坐标,q为视线俯仰角,q为视线偏航角,r为evader到pursuer的相对距离,LOS0为evader到pursuer的初始视线;Figure 1 is a schematic diagram of the relative motion relationship between the pursuer and the evader, E is the coordinate of the evader in the evader inertial system, P is the coordinate of the pursuer in the evader inertial system, q is the pitch angle of the line of sight, q is the yaw angle of the line of sight, and r is The relative distance from evader to pursuer, LOS 0 is the initial line of sight from evader to pursuer;

图2为Evader常值机动情况下对Pursuer俯仰平面导航常数Nε估计结果图;Figure 2 is the estimation result of Pursuer pitch plane navigation constant N ε under the condition of Evader constant value maneuver;

图3为Evader常值机动情况下对Pursuer偏航平面导航常数Nβ估计结果图;Fig. 3 is the estimation result diagram of Pursuer yaw plane navigation constant N β in the case of Evader constant value maneuver;

图4为Evader无机动情况下使用PN制导律Kalman辨识滤波器对Pursuer加速度在Evader惯性坐标系y轴的分量ay估计结果图,g为重力加速度;Fig. 4 is the estimation result diagram of the component a y of the Pursuer acceleration on the y-axis of the Evader inertial coordinate system using the PN guidance law Kalman identification filter in the case of Evader without maneuvering, and g is the gravitational acceleration;

图5为Evader无机动情况下使用PN制导律Kalman辨识滤波器对Pursuer加速度在Evader惯性坐标系z轴的分量az估计结果图,g为重力加速度;Fig. 5 is the estimation result diagram of the component a z of the Pursuer acceleration on the z axis of the Evader inertial coordinate system using the PN guidance law Kalman identification filter in the case of Evader without maneuvering, and g is the acceleration of gravity;

图6为Evader无机动情况下使用MMAE滤波器对Pursuer加速度ay估计结果图;Fig. 6 is the estimation result diagram of the Pursuer acceleration a y using the MMAE filter in the case of Evader without maneuvering;

图7为Evader无机动情况下使用MMAE滤波器对Pursuer加速度az估计结果图;Fig. 7 is the estimation result diagram of Pursuer's acceleration az using MMAE filter under Evader's non-maneuvering situation;

图8为Evader常值机动情况下使用PN制导律Kalman辨识滤波器对Pursuer加速度ay估计结果图;Fig. 8 is a diagram of the estimation results of the Pursuer acceleration a y using the PN guidance law Kalman identification filter in the case of Evader constant maneuvering;

图9为Evader常值机动情况下使用PN制导律Kalman辨识滤波器对Pursuer加速度az估计结果图;Fig. 9 is a graph of the estimation results of the Pursuer acceleration az using the PN guidance law Kalman identification filter under the condition of Evader constant value maneuver;

图10为Evader常值机动情况下使用MMAE滤波器对Pursuer加速度ay估计结果图;Fig. 10 is the estimation result diagram of Pursuer acceleration a y using MMAE filter under Evader constant maneuvering situation;

图11为Evader常值机动情况下使用MMAE滤波器对Pursuer加速度az估计结果图;Fig. 11 is the estimation result diagram of Pursuer acceleration az using MMAE filter under Evader constant value maneuvering situation;

图12为Evader正弦机动情况下使用PN制导律Kalman辨识滤波器对Pursuer加速度ay估计结果图;Fig. 12 is a graph of the estimation results of the Pursuer acceleration a y using the PN guidance law Kalman identification filter in the case of Evader sinusoidal maneuver;

图13为Evader正弦机动情况下使用PN制导律Kalman辨识滤波器对Pursuer加速度az估计结果图;Fig. 13 is a diagram of the estimation results of the Pursuer acceleration a z using the PN guidance law Kalman identification filter in the case of the Evader sinusoidal maneuver;

图14为Evader正弦机动情况下使用MMAE滤波器对Pursuer加速度ay估计结果图;Fig. 14 is the estimation result diagram of Pursuer acceleration a y using MMAE filter in the case of Evader sinusoidal maneuver;

图15为Evader正弦机动情况下使用MMAE滤波器对Pursuer加速度az估计结果图。Fig. 15 is a diagram showing the results of estimation of the Pursuer acceleration a z using the MMAE filter in the case of Evader sinusoidal maneuvering.

具体实施方式Detailed ways

具体实施方式一:本实施方式的一种比例导引制导律辨识滤波方法,具体是按照以下步骤制备的:Specific implementation mode 1: A proportional guidance guidance law identification and filtering method in this implementation mode is specifically prepared according to the following steps:

步骤一:pursuer有两种运动模型:控制器饱和时的饱和运动模型和控制器未饱和时的PN制导律运动模型;建立PN制导律运动模型在俯仰平面的状态方程、PN制导律运动模型在偏航平面的状态方程、控制器饱和运动模型在俯仰平面的状态方程和控制器饱和运动模型在偏航平面的状态方程;evader为弹道导弹;pursuer为拦截导弹;PN制导律为比例导引制导律;Step 1: The pursuer has two motion models: the saturated motion model when the controller is saturated and the PN guidance law motion model when the controller is not saturated; establish the state equation of the PN guidance law motion model in the pitch plane, and the PN guidance law motion model in the The state equation of the yaw plane, the state equation of the controller saturation motion model in the pitch plane, and the state equation of the controller saturation motion model in the yaw plane; evader is a ballistic missile; pursuer is an intercept missile; PN guidance law is proportional guidance guidance law;

步骤二:根据步骤一的PN制导律运动模型在俯仰平面的状态方程、PN制导律运动模型在偏航平面的状态方程、控制器饱和运动模型在俯仰平面的状态方程和控制器饱和运动模型在偏航平面的状态方程,计算PN制导律运动模型在俯仰平面的状态转移矩阵、PN制导律运动模型在偏航平面的状态转移矩阵、控制器饱和运动模型在俯仰平面的状态转移矩阵、控制器饱和运动模型在偏航平面的状态转移矩阵;Step 2: According to the state equation of the PN guidance law motion model in the pitch plane in step 1, the state equation of the PN guidance law motion model in the yaw plane, the state equation of the controller saturation motion model in the pitch plane, and the controller saturation motion model in The state equation of the yaw plane, calculate the state transition matrix of the PN guidance law motion model in the pitch plane, the state transition matrix of the PN guidance law motion model in the yaw plane, the state transition matrix of the controller saturation motion model in the pitch plane, and the controller The state transition matrix of the saturated motion model in the yaw plane;

步骤三:根据步骤二的PN制导律运动模型在俯仰平面的状态转移矩阵、PN制导律运动模型在偏航平面的状态转移矩阵、控制器饱和运动模型在俯仰平面的状态转移矩阵、控制器饱和运动模型在偏航平面的状态转移矩阵,设计PN制导律运动模型在俯仰平面的Kalman滤波器方程、PN制导律运动模型在偏航平面的Kalman滤波器方程、控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程,Kalman滤波器方程为卡尔曼滤波器方程;Step 3: According to the state transition matrix of the PN guidance law motion model in the pitch plane in step 2, the state transition matrix of the PN guidance law motion model in the yaw plane, the state transition matrix of the controller saturation motion model in the pitch plane, and the controller saturation The state transition matrix of the motion model in the yaw plane, design the Kalman filter equation of the PN guidance law motion model in the pitch plane, the Kalman filter equation of the PN guidance law motion model in the yaw plane, and the controller saturation motion model in the pitch plane The Kalman filter equation and the Kalman filter equation of the controller saturation motion model in the yaw plane, and the Kalman filter equation is the Kalman filter equation;

步骤四:根据步骤三的PN制导律运动模型在俯仰平面的Kalman滤波器方程、PN制导律运动模型在偏航平面的Kalman滤波器方程、控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程,计算pursuer当前时刻PN制导律运动模型的后验概率和控制器饱和运动模型的后验概率;Step 4: According to the Kalman filter equation of the PN guidance law motion model in step 3 in the pitch plane, the Kalman filter equation of the PN guidance law motion model in the yaw plane, the Kalman filter equation of the controller saturation motion model in the pitch plane, and The Kalman filter equation of the controller saturated motion model in the yaw plane, calculate the posterior probability of the PN guidance law motion model and the posterior probability of the controller saturated motion model at the current moment of the pursuer;

步骤五:根据步骤四得到的pursuer当前时刻PN制导律运动模型的后验概率和饱和运动模型的后验概率,计算pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻;在pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻前采用控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程的估计结果;否则采用PN制导律运动模型在俯仰平面的Kalman滤波器方程和PN制导律运动模型在偏航平面的Kalman滤波器方程的估计结果。Step 5: According to the posterior probability of the PN guidance law motion model and the posterior probability of the saturated motion model obtained in step 4, calculate the PN guidance when the pursuer switches from the saturated motion model when the controller is saturated to the unsaturated controller The moment of the law motion model; before the moment when the pursuer switches from the saturated motion model when the controller is saturated to the PN guidance law motion model when the controller is not saturated, use the Kalman filter equation of the controller saturated motion model in the pitch plane and the controller The estimation result of the Kalman filter equation of the saturated motion model in the yaw plane; otherwise, the Kalman filter equation of the PN guidance law motion model in the pitch plane and the estimation result of the Kalman filter equation of the PN guidance law motion model in the yaw plane.

具体实施方式二:本实施方式与具体实施方式一不同的是:所述步骤一中pursuer有两种运动模型:控制器饱和时的饱和运动模型和控制器未饱和时的PN制导律运动模型;建立PN制导律运动模型在俯仰平面的状态方程、PN制导律运动模型在偏航平面的状态方程、控制器饱和运动模型在俯仰平面的状态方程和控制器饱和运动模型在偏航平面的状态方程;具体过程为:Embodiment 2: The difference between this embodiment and Embodiment 1 is that the pursuer in the step 1 has two motion models: a saturated motion model when the controller is saturated and a PN guidance law motion model when the controller is not saturated; Establish the state equation of the PN guidance law motion model in the pitch plane, the state equation of the PN guidance law motion model in the yaw plane, the state equation of the controller saturation motion model in the pitch plane, and the state equation of the controller saturation motion model in the yaw plane ; The specific process is:

在导弹突防场景中,突防的弹道导弹要逃避对方导弹防御系统发射的拦截导弹的拦截,是逃逸方,称为evader,和它相关的所有物理量都含有下标e;拦截导弹是追捕方,称为pursuer,和它相关的所有物理量都含有下标p;假设拦截导弹的制导控制可解耦为纵向平面和侧向平面。图1描绘了pursuer和evader在俯仰平面的相对运动关系。In the missile penetration scenario, if the penetrating ballistic missile wants to evade the interception of the intercepting missile launched by the opponent’s missile defense system, it is the evader, called evader, and all physical quantities related to it contain the subscript e; the intercepting missile is the hunting party , called the pursuer, and all physical quantities related to it contain the subscript p; it is assumed that the guidance and control of the interceptor missile can be decoupled into the longitudinal plane and the lateral plane. Figure 1 depicts the relative motion relationship between the pursuer and the evader in the pitch plane.

Evader视线坐标系定义;视线坐标系o′x4y4z4原点o′位于导引头回转中心,o′x4轴与目标—导弹视线一致,由导引头回转中心指向目标为正,o′y4轴位于包含o′x4轴的铅锤面内,与o′x4轴垂直,指向上方为正,o′z4轴按右手定则确定;而evader的末制导惯性坐标系定义为制导初始时刻的evader视线坐标系o0x0y0z0Evader line-of-sight coordinate system definition; line-of-sight coordinate system o′x 4 y 4 z 4 origin o′ is located at the center of rotation of the seeker, and the axis o′x 4 is consistent with the target-missile line of sight, and it is positive from the center of rotation of the seeker to the target. The o′y 4 -axis is located in the plumb plane including the o′x 4 -axis, perpendicular to the o′x 4 -axis, and pointing upward is positive, and the o′z 4 -axis is determined by the right-hand rule; while the evader’s terminal guidance inertial coordinate system Defined as the evader line of sight coordinate system o 0 x 0 y 0 z 0 at the initial moment of guidance;

pursuer在制导初期需要消除对齐误差,这有可能导致pursuer控制器进入饱和状态,待误差消除后,pursuer控制器会退出饱和状态;因此pursuer有两种运动模型:控制器饱和时的饱和运动模型和控制器未饱和时的PN制导律运动模型;下面给出pursuer两种运动在俯仰平面和偏航平面的模型;Pursuer needs to eliminate the alignment error in the early stage of guidance, which may cause the pursuer controller to enter a saturated state. After the error is eliminated, the pursuer controller will exit the saturated state; therefore, the pursuer has two motion models: the saturation motion model when the controller is saturated and The motion model of the PN guidance law when the controller is not saturated; the models of the two motions of the pursuer in the pitch plane and the yaw plane are given below;

pursuer的PN制导律运动方程建模,rx为evader和pursuer相对位置在evader惯性系x轴下的分量;ry为evader和pursuer相对位置在evader惯性系y轴下的分量;rz为evader和pursuer相对位置在evader惯性系z轴下的分量;evader为弹道导弹,pursuer为拦截导弹;The PN guidance law motion equation of the pursuer is modeled, r x is the component of the relative position of the evader and the pursuer under the x-axis of the evader inertial system; r y is the component of the relative position of the evader and the pursuer under the y-axis of the evader inertial system; r z is the component of the evader The component relative to the position of the pursuer under the z-axis of the evader inertial system; the evader is a ballistic missile, and the pursuer is an interceptor missile;

rx=xp-xe ry=yp-ye rz=zp-ze (1)r x =x p -x e r y =y p -y e r z =z p -z e (1)

式中,[xp,yp,zp]T和[xe,ye,ze]T分别是pursuer和evader在evader惯性系下的位置,T为转置,xp为pursuer在evader惯性系x轴下的位置,yp为pursuer在evader惯性系y轴下的位置,zp为pursuer在evader惯性系z轴下的位置,xe为evader在evader惯性系x轴下的位置,ye为evader在evader惯性系y轴下的位置,ze为evader在evader惯性系z轴下的位置;In the formula, [x p , y p , z p ] T and [x e , y e , z e ] T are the positions of the pursuer and evader in the evader inertial system respectively, T is the transpose, and x p is the position of the pursuer in the evader The position under the x-axis of the inertial system, y p is the position of the pursuer under the y-axis of the evader inertial system, z p is the position of the pursuer under the z-axis of the evader inertial system, x e is the position of the evader under the x-axis of the evader inertial system, y e is the position of evader under the y-axis of the evader inertial system, and z e is the position of the evader under the z-axis of the evader inertial system;

PN制导律运动模型在俯仰平面的状态方程为The state equation of the PN guidance law motion model in the pitch plane is

PN制导律运动模型在偏航平面的状态方程为The state equation of the PN guidance law motion model in the yaw plane is

式中,vx为evader和pursuer相对速度在evader惯性系x轴下的分量,vy为evader和pursuer相对速度在evader惯性系y轴下的分量,vz为evader和pursuer相对速度在evader惯性系z轴下的分量;apx为pursuer的加速度在evader惯性系x轴下的分量,apy为pursuer的加速度在evader惯性系y轴下的分量,apz为pursuer的加速度在evader惯性系z轴下的分量;aex为evader的加速度在evader惯性系x轴下的分量,aey为evader的加速度在evader惯性系y轴下的分量,aez为evader的加速度在evader惯性系z轴下的分量;Nε是pursuer在俯仰平面的导航常数;Nβ是pursuer在偏航平面的导航常数;τ是时间常数,为rx的一阶导数,为ry的一阶导数,为rz的一阶导数,为vx的一阶导数,为vy的一阶导数,为vz的一阶导数,为apx的一阶导数,为apy的一阶导数,为apz的一阶导数,为Nε的一阶导数,为Nβ的一阶导数;In the formula, v x is the component of the relative velocity of evader and pursuer under the x-axis of the evader inertial system, v y is the component of the relative velocity of evader and pursuer under the y-axis of the evader inertial system, and v z is the relative velocity of evader and pursuer under the evader inertial system The component under the z-axis of the system; a px is the component of the acceleration of the pursuer under the x-axis of the evader inertial system, a py is the component of the acceleration of the pursuer under the y-axis of the evader inertial system, a pz is the acceleration of the pursuer under the z-axis of the evader inertial system The component under the axis; a ex is the component of the acceleration of the evader under the x-axis of the evader inertial system, a ey is the component of the acceleration of the evader under the y-axis of the evader inertial system, a ez is the acceleration of the evader under the z-axis of the evader inertial system Component; N ε is the navigation constant of the pursuer in the pitch plane; N β is the navigation constant of the pursuer in the yaw plane; τ is the time constant, is the first derivative of r x , is the first derivative of r y , is the first derivative of r z , is the first derivative of v x , is the first derivative of v y , is the first derivative of v z , is the first derivative of a px , is the first derivative of a py , is the first derivative of a pz , is the first derivative of N ε , is the first derivative of N β ;

假设evader能获得pursuer在惯性系下的位置,进而计算出相对自身(即evader)的位置;计算出相对自身的位置的测量方程为Assume that evader can obtain the position of the pursuer in the inertial system, and then calculate the position relative to itself (ie evader); the measurement equation for calculating the position relative to itself is

h=[rx,ry,rz]T (4)h=[r x ,r y ,r z ] T (4)

控制器饱和运动模型在俯仰平面的状态方程为The state equation of the controller saturation motion model in the pitch plane is

控制器饱和运动模型在偏航平面的状态方程为The state equation of the controller saturation motion model in the yaw plane is

式中,d1和d2是控制器饱和运动模型的填充状态;目的是为了和PN制导律运动模型的维数保持一致。In the formula, d 1 and d 2 are the filling states of the saturated motion model of the controller; the purpose is to keep consistent with the dimension of the PN guidance law motion model.

其它步骤及参数与具体实施方式一相同。Other steps and parameters are the same as those in Embodiment 1.

具体实施方式三:本实施方式与具体实施方式一或二不同的是:所述步骤二中根据步骤一的PN制导律运动模型在俯仰平面的状态方程、PN制导律运动模型在偏航平面的状态方程、控制器饱和运动模型在俯仰平面的状态方程和控制器饱和运动模型在偏航平面的状态方程,计算PN制导律运动模型在俯仰平面的状态转移矩阵、PN制导律运动模型在偏航平面的状态转移矩阵、控制器饱和运动模型在俯仰平面的状态转移矩阵、控制器饱和运动模型在偏航平面的状态转移矩阵;具体过程为:Specific embodiment three: the difference between this embodiment and specific embodiment one or two is: in described step 2, according to the state equation of the PN guidance law motion model of step one in the pitch plane, the PN guidance law motion model in the yaw plane The state equation, the state equation of the controller saturation motion model in the pitch plane and the state equation of the controller saturation motion model in the yaw plane, calculate the state transition matrix of the PN guidance law motion model in the pitch plane, and the PN guidance law motion model in the yaw plane The state transition matrix of the plane, the state transition matrix of the controller saturation motion model in the pitch plane, and the state transition matrix of the controller saturation motion model in the yaw plane; the specific process is:

PN制导律运动模型在俯仰平面的状态转移矩阵为The state transition matrix of the PN guidance law motion model in the pitch plane is

式中,T为测量周期,单位为秒;In the formula, T is the measurement cycle in seconds;

PN制导律运动模型在偏航平面的状态转移矩阵为The state transition matrix of the PN guidance law motion model in the yaw plane is

控制器饱和运动模型在俯仰平面的状态转移矩阵为The state transition matrix of the controller saturation motion model in the pitch plane is

控制器饱和运动模型在偏航平面的状态转移矩阵为The state transition matrix of the controller saturation motion model in the yaw plane is

其它步骤及参数与具体实施方式一或二相同。Other steps and parameters are the same as those in Embodiment 1 or Embodiment 2.

具体实施方式四:本实施方式与具体实施方式一至三之一不同的是:所述步骤三中根据步骤二的PN制导律运动模型在俯仰平面的状态转移矩阵、PN制导律运动模型在偏航平面的状态转移矩阵、控制器饱和运动模型在俯仰平面的状态转移矩阵、控制器饱和运动模型在偏航平面的状态转移矩阵,设计PN制导律运动模型在俯仰平面的Kalman滤波器方程、PN制导律运动模型在偏航平面的Kalman滤波器方程、控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程,Kalman滤波器方程为卡尔曼滤波器方程;具体过程为:Specific embodiment four: this embodiment is different from one of specific embodiments one to three in that: in the step 3, according to the state transition matrix of the PN guidance law motion model of step 2 in the pitch plane, the PN guidance law motion model is in the yaw The state transition matrix of the plane, the state transition matrix of the controller saturated motion model in the pitch plane, the state transition matrix of the controller saturated motion model in the yaw plane, design the Kalman filter equation of the PN guidance law motion model in the pitch plane, and the PN guidance The Kalman filter equation of the law motion model in the yaw plane, the Kalman filter equation of the controller saturation motion model in the pitch plane, and the Kalman filter equation of the controller saturation motion model in the yaw plane, and the Kalman filter equation is Kalman filter The device equation; the specific process is:

PN制导律辨识滤波器和饱和辨识滤波器均采用EKF,EKF为扩展卡尔曼滤波器,根据步骤二计算得到的状态转移矩阵,Both the PN guidance law identification filter and the saturation identification filter use EKF, and EKF is an extended Kalman filter. According to the state transition matrix calculated in step 2,

PN制导律运动模型在俯仰平面的Kalman滤波器方程为The Kalman filter equation of the PN guidance law motion model in the pitch plane is

PN制导律运动模型在偏航平面的Kalman滤波器方程为The Kalman filter equation of the PN guidance law motion model in the yaw plane is

控制器饱和运动模型在俯仰平面的Kalman滤波器方程为The Kalman filter equation of the controller saturation motion model in the pitch plane is

控制器饱和运动模型在偏航平面的Kalman滤波器方程为The Kalman filter equation of the controller saturation motion model in the yaw plane is

式中,是对k+1时刻状态的估计值;是对k+1时刻状态的预报值;分别是k时刻和k+1时刻的PN制导律运动模型在俯仰平面的滤波器增益;分别是k时刻和k+1时刻的PN制导律运动模型在偏航平面的滤波器增益;In the formula, is the estimated value of the state at time k+1; is the predicted value of the state at k+1 time; and are the filter gains of the PN guidance law motion model in the pitch plane at time k and k+1, respectively; and are the filter gains of the PN guidance law motion model in the yaw plane at time k and k+1, respectively;

分别是k时刻和k+1时刻的控制器饱和运动模型在俯仰平面的滤波器增益; and are the filter gains of the controller saturation motion model in the pitch plane at time k and k+1, respectively;

分别是k时刻和k+1时刻的控制器饱和运动模型在偏航平面的滤波器增益; and are the filter gains in the yaw plane of the controller saturated motion model at time k and k+1, respectively;

为k+1时刻PN制导律运动模型在俯仰平面的状态预报的协方差矩阵; is the covariance matrix of the state prediction of the PN guidance law motion model in the pitch plane at time k+1;

为k+1时刻PN制导律运动模型在偏航平面的状态预报的协方差矩阵; is the covariance matrix of the state prediction of the PN guidance law motion model in the yaw plane at time k+1;

为k+1时刻控制器饱和运动模型在俯仰平面的状态预报的协方差矩阵; is the covariance matrix of the state prediction of the controller saturated motion model in the pitch plane at time k+1;

为k+1时刻控制器饱和运动模型在偏航平面的状态预报的协方差矩阵; is the covariance matrix of the state prediction of the controller saturated motion model in the yaw plane at time k+1;

R(1)为PN制导律运动模型在俯仰平面的测量噪声协方差矩阵;R (1) is the measurement noise covariance matrix of the PN guidance law motion model in the pitch plane;

R(2)为PN制导律运动模型在偏航平面的测量噪声协方差矩阵;R (2) is the measurement noise covariance matrix of the PN guidance law motion model in the yaw plane;

R(3)为控制器饱和运动模型在俯仰平面的测量噪声协方差矩阵;R (3) is the measurement noise covariance matrix of the controller saturation motion model in the pitch plane;

R(4)为控制器饱和运动模型在偏航平面的测量噪声协方差矩阵;R (4) is the measurement noise covariance matrix of the controller saturation motion model in the yaw plane;

为PN制导律运动模型在俯仰平面的k时刻到k+1时刻的状态转移矩阵; is the state transition matrix of the PN guidance law motion model from time k to time k+1 in the pitch plane;

为PN制导律运动模型在偏航平面的k时刻到k+1时刻的状态转移矩阵; is the state transition matrix of the PN guidance law motion model from time k to k+1 in the yaw plane;

为控制器饱和运动模型在俯仰平面的k时刻到k+1时刻的状态转移矩阵; is the state transition matrix of the controller saturation motion model from time k to time k+1 in the pitch plane;

控制器饱和运动模型在偏航平面的k时刻到k+1时刻的状态转移矩阵; The state transition matrix of the controller saturation motion model from time k to k+1 in the yaw plane;

Q(1)为PN制导律运动模型在俯仰平面的过程噪声协方差矩阵;Q (1) is the process noise covariance matrix of the PN guidance law motion model in the pitch plane;

Q(2)为PN制导律运动模型在偏航平面的过程噪声协方差矩阵;Q (2) is the process noise covariance matrix of the PN guidance law motion model in the yaw plane;

Q(3)为控制器饱和运动模型在俯仰平面的过程噪声协方差矩阵;Q (3) is the process noise covariance matrix of the controller saturated motion model in the pitch plane;

Q(4)为控制器饱和运动模型在偏航平面的过程噪声协方差矩阵;Q (4) is the process noise covariance matrix of the controller saturation motion model in the yaw plane;

为PN制导律运动模型在俯仰平面的k时刻状态估计协方差矩阵; Estimating the covariance matrix for the state of the PN guidance law motion model at time k in the pitch plane;

为PN制导律运动模型在偏航平面的k时刻状态估计协方差矩阵; Estimate the covariance matrix for the state of the PN guidance law motion model at time k in the yaw plane;

为控制器饱和运动模型在俯仰平面的k时刻状态估计协方差矩阵; Estimating the covariance matrix for the state of the controller's saturated motion model at time k in the pitch plane;

为控制器饱和运动模型在偏航平面的k时刻状态估计协方差矩阵; Estimating the covariance matrix for the state of the controller's saturated motion model at time k in the yaw plane;

为PN制导律运动模型在俯仰平面的k+1时刻状态估计协方差矩阵; Estimate the covariance matrix for the state of the PN guidance law motion model at time k+1 in the pitch plane;

为PN制导律运动模型在偏航平面的k+1时刻状态估计协方差矩阵; Estimate the covariance matrix for the state of the PN guidance law motion model at time k+1 of the yaw plane;

为控制器饱和运动模型在俯仰平面的k+1时刻状态估计协方差矩阵; Estimating the covariance matrix for the state of the controller saturated motion model at time k+1 in the pitch plane;

为控制器饱和运动模型在偏航平面的k+1时刻状态估计协方差矩阵; Estimate the covariance matrix for the state of the controller saturated motion model at time k+1 of the yaw plane;

为PN制导律运动模型在俯仰平面的k时刻到k+1时刻的状态转移矩阵的转置; is the transposition of the state transition matrix of the PN guidance law motion model from time k to time k+1 in the pitch plane;

为PN制导律运动模型在偏航平面的k时刻到k+1时刻的状态转移矩阵的转置; is the transposition of the state transition matrix of the PN guidance law motion model from time k to k+1 in the yaw plane;

为控制器饱和运动模型在俯仰平面的k时刻到k+1时刻的状态转移矩阵的转置; is the transposition of the state transition matrix of the controller saturation motion model from time k to k+1 in the pitch plane;

为控制器饱和运动模型在偏航平面的k时刻到k+1时刻的状态转移矩阵的转置; is the transposition of the state transition matrix of the controller saturation motion model from time k to k+1 in the yaw plane;

Zk+1为k+1时刻的测量值;I为7×7单位矩阵;H为测量矩阵,H=[1 1 0 0 0 0 0],HT为H的转置,k为采样时刻。Z k+1 is the measured value at time k+1; I is the 7×7 unit matrix; H is the measurement matrix, H=[1 1 0 0 0 0 0], H T is the transpose of H, and k is the sampling time .

其它步骤及参数与具体实施方式一至三之一相同。Other steps and parameters are the same as those in Embodiments 1 to 3.

具体实施方式五:本实施方式与具体实施方式一至四之一不同的是:所述步骤四中根据步骤三的PN制导律运动模型在俯仰平面的Kalman滤波器方程、PN制导律运动模型在偏航平面的Kalman滤波器方程、控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程(2个俯仰一组为一个MMAE,2个偏航一组为一个MMAE),计算pursuer当前时刻PN制导律运动模型的后验概率和饱和控制器运动模型的后验概率;具体过程为:Specific embodiment five: this embodiment is different from one of specific embodiments one to four in that: in the step 4, according to the Kalman filter equation of the PN guidance law motion model of step 3 in the pitch plane, the PN guidance law motion model is in the deflection The Kalman filter equation of the flight plane, the Kalman filter equation of the controller saturation motion model in the pitch plane, and the Kalman filter equation of the controller saturation motion model in the yaw plane (two pitches are one MMAE, two yaw One group is one MMAE), calculate the posterior probability of the motion model of the PN guidance law of the pursuer and the posterior probability of the motion model of the saturated controller at the current moment; the specific process is:

设模型集合为M={Mj|j=1,…,r},式中,r为模型的个数,r=2,M1表示PN制导律运动模型,M2表示控制器饱和运动模型,在k时刻模型j的先验概率为P{Mj|Zk-1},在k时刻模型j的后验概率为P{Mj|Zk},第k时刻的测量值为z(k),前k个时刻的测量值序列定义为Let the model collection be M={M j |j=1,...,r}, where r is the number of models, r=2, M 1 represents the PN guidance law motion model, M 2 represents the controller saturation motion model , the prior probability of model j at time k is P{M j |Z k-1 }, the posterior probability of model j at time k is P{M j |Z k }, and the measured value at time k is z( k), the sequence of measured values at the first k moments is defined as

根据贝叶斯定理,According to Bayes' theorem,

式中,P{Mi|Zk-1}为模型i的先验概率,p(z(k)|Zk-1,Mi)为模型i在k时刻的似然函数,P{Mj|z(k),Zk-1}为k时刻模型Mj的后验概率,p(z(k)|Zk-1,Mj)为模型j在k时刻的似然函数;因此模型j的后验概率可通过所有模型的先验概率和似然函数计算得到,模型i的似然函数是该模型对应滤波器新息的分布函数;假定新息服从均值为零的正态分布,设其k时刻的方差为Si(k),即In the formula, P{M i |Z k-1 } is the prior probability of model i, p(z(k)|Z k-1 ,M i ) is the likelihood function of model i at time k, and P{M j |z(k),Z k-1 } is the posterior probability of model M j at time k, p(z(k)|Z k-1 ,M j ) is the likelihood function of model j at time k; therefore The posterior probability of model j can be calculated by the prior probability and likelihood function of all models, and the likelihood function of model i is the distribution function of the model corresponding to the filter innovation; it is assumed that the innovation follows a normal distribution with a mean value of zero , let the variance at time k be S i (k), that is

p(z(k)|Zk-1,Mi)=N(0,Sj(k)) (39)p(z(k)|Z k-1 ,M i )=N(0,S j (k)) (39)

式中,N(0,Sj(k))代表均值为0,方差为Sj(k)正态分布的概率密度函数;Sj(k)是模型j新息的方差,即In the formula, N(0,S j (k)) represents the mean value is 0, and the variance is the probability density function of the normal distribution of S j (k); S j (k) is the variance of model j innovation, namely

Sj(k)=E[v(k)vT(k)] (41)S j (k)=E[v(k)v T (k)] (41)

式中,v(k)为新息,vT(k)是v(k)的转置,E[·]是数学期望函数;In the formula, v(k) is the new information, v T (k) is the transposition of v(k), E[ ] is the mathematical expectation function;

设模型i对应的Kalman滤波器的状态估计值为则MMAE滤波器第k时刻的状态估计值Let the state estimate of the Kalman filter corresponding to model i be Then the estimated value of the state of the MMAE filter at the kth moment for

i取值为1或2。i takes the value of 1 or 2.

其它步骤及参数与具体实施方式一至四之一相同。Other steps and parameters are the same as in one of the specific embodiments 1 to 4.

具体实施方式六:本实施方式与具体实施方式一至五之一不同的是:所述步骤五中根据步骤四得到的pursuer当前时刻PN制导律运动模型的后验概率和饱和运动模型的后验概率,计算pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻;在pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻前采用控制器饱和运动模型在俯仰平面的Kalman滤波器方程和控制器饱和运动模型在偏航平面的Kalman滤波器方程的估计结果;否则采用PN制导律运动模型在俯仰平面的Kalman滤波器方程和PN制导律运动模型在偏航平面的Kalman滤波器方程的估计结果;具体过程为:Embodiment 6: This embodiment differs from one of Embodiments 1 to 5 in that: the posterior probability of the PN guidance law motion model and the posterior probability of the saturated motion model obtained in step 5 according to the current moment of the pursuer , calculate the moment when the pursuer switches from the saturated motion model when the controller is saturated to the PN guidance law when the controller is not saturated; when the pursuer switches from the saturated motion model when the controller is saturated to the PN guidance law when the controller is not saturated The Kalman filter equation of the controller saturation motion model in the pitch plane and the estimation result of the Kalman filter equation of the controller saturation motion model in the yaw plane are used before the moment of the motion model; otherwise, the Kalman filter equation of the PN guidance law motion model in the pitch plane is used. Estimation results of the filter equation and the Kalman filter equation of the PN guidance law motion model in the yaw plane; the specific process is:

考虑到在制导的第一阶段,pursuer外在行为由于饱和的影响不再是PN制导律的运动模型,之后解除饱和进入PN制导律的运动模型。我们采用两个MMAE滤波器分别对两个阶段pursuer的状态进行估计;将这两个MMAE滤波器分别称为MMAE_A和MMAE_B滤波器,其中MMAE_A滤波器估计第一阶段,即pursuer处于控制器饱和时的状态;而MMAE_B滤波器估计第二阶段,即pursuer处于控制器未饱和时的状态;另外需要估计出pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻,从而确定在某时刻t采用哪一个MMAE滤波器的估计值来代表当前pursuer的状态;Considering that in the first stage of guidance, the external behavior of the pursuer is no longer the motion model of the PN guidance law due to the influence of saturation, and then the saturation is released to enter the motion model of the PN guidance law. We use two MMAE filters to estimate the state of the two-stage pursuer respectively; these two MMAE filters are called MMAE_A and MMAE_B filters respectively, where the MMAE_A filter estimates the first stage, that is, when the pursuer is in controller saturation state; while the MMAE_B filter estimates the second stage, that is, the state of the pursuer when the controller is not saturated; in addition, it is necessary to estimate that the pursuer switches from the saturated motion model when the controller is saturated to the PN guidance law motion model when the controller is not saturated , so as to determine which estimated value of the MMAE filter is used at a certain time t to represent the state of the current pursuer;

若pursuer未进入控制器饱和时的状态,则pursuer全程的运动模型保持一致,MMAE_A和MMAE_B估计出的结果一致;If the pursuer does not enter the state when the controller is saturated, the motion model of the pursuer is consistent throughout the whole process, and the estimated results of MMAE_A and MMAE_B are consistent;

pursuer从控制器饱和时的饱和运动模型切换到控制器未饱和时的PN制导律运动模型的时刻的估计方法具体过程为:The specific process of estimating the moment when the pursuer switches from the saturated motion model when the controller is saturated to the PN guidance law motion model when the controller is not saturated is as follows:

每个MMAE滤波器中均含有两个模型:控制器未饱和时的PN制导律运动模型和控制器饱和时的饱和运动模型;根据多模型滤波方法计算出当前滤波器两个模型的后验概率,然后根据当前滤波器两个模型的后验概率,确定pursuer目前是否处于控制器饱和时的饱和阶段;MMAE_B滤波器采用扫描的方式来检查饱和退出点,即每间隔check_period长的时间,根据式(48)判断当前pursuer的运动模型,假设MMAE滤波器中控制器未饱和时的PN制导律运动模型和控制器饱和时的饱和运动模型的后验概率分别为P{MPNG}和P{MSat},则判断过程为:Each MMAE filter contains two models: the PN guidance law motion model when the controller is not saturated and the saturated motion model when the controller is saturated; the posterior probability of the two models of the current filter is calculated according to the multi-model filtering method , and then according to the posterior probability of the two models of the current filter, it is determined whether the pursuer is currently in the saturation stage when the controller is saturated; the MMAE_B filter uses a scanning method to check the saturation exit point, that is, the interval of check_period is long, according to the formula (48) Judging the motion model of the current pursuer, assuming that the posterior probabilities of the PN guidance law motion model when the controller is not saturated and the saturated motion model when the controller is saturated in the MMAE filter are P{M PNG } and P{M Sat }, the judgment process is:

式中,ε为门限值。当P{MPNG}-P{MSat}≥ε时,则当前pursuer采用的模型为PN制导律模型,判定pursuer退出了饱和状态;否则重新启动和初始化MMAE_B滤波器。每隔check_period时间进行一次判断,直到判定pursuer退出了饱和状态位置。check_period为检验周期,一般check_period设置为0.5秒。In the formula, ε is the threshold value. When P{M PNG }-P{M Sat }≥ε, the model adopted by the current pursuer is the PN guidance law model, and it is determined that the pursuer has exited the saturation state; otherwise, restart and initialize the MMAE_B filter. A judgment is made every check_period time until it is determined that the pursuer has exited the saturation state. check_period is the inspection period, and generally check_period is set to 0.5 seconds.

当ε取较大的值时,对pursuer当前模型的判断准确率更高;反之当ε取较小值时,对pursuer当前模型的判断准确率更低。但这并不意味着ε的取值越大越好,后面会针对此进行分析。When ε takes a larger value, the judgment accuracy of the current model of the pursuer is higher; on the contrary, when ε takes a smaller value, the judgment accuracy of the current model of the pursuer is lower. But this does not mean that the larger the value of ε, the better, which will be analyzed later.

check_period和参数ε有关当ε取很大值时,需要的check_period也较大,这就会导致pursuer退出饱和的时间点估计不准确;当ε取很小值时,check_period可以取更小的值。表面上看这可以提高pursuer退出饱和时间点估计的精度,但是当ε取很小值时对pursuer当前模型的判断准确度会降低。因此需要选择合适的ε和check_period值才能兼顾pursuer当前模型判断准确度和退出饱和时刻估计的精确度。check_period is related to the parameter ε. When ε takes a large value, the required check_period is also large, which will lead to inaccurate estimation of the time point when the pursuer exits saturation; when ε takes a small value, check_period can take a smaller value. On the surface, this can improve the accuracy of the estimation of the time point when the pursuer exits saturation, but when ε takes a small value, the accuracy of the judgment of the current model of the pursuer will decrease. Therefore, it is necessary to select the appropriate ε and check_period values to take into account the accuracy of the current model judgment of the pursuer and the accuracy of the estimation of the exit saturation time.

下面是MMAE_B对pursuer退出饱和时刻估计功能的伪代码,算法计算出pursuer退出饱和时刻的估计值quit_saturation_time之后退出。The following is the pseudocode of MMAE_B’s function of estimating the exit saturation time of the pursuer. The algorithm calculates the estimated value quit_saturation_time of the exit saturation time of the pursuer and then exits.

其它步骤及参数与具体实施方式一至五之一相同。Other steps and parameters are the same as one of the specific embodiments 1 to 5.

采用以下实施例验证本发明的有益效果:Adopt the following examples to verify the beneficial effects of the present invention:

实施例一:Embodiment one:

本实施例一种比例导引制导律辨识滤波方法具体是按照以下步骤制备的:In this embodiment, a proportional guidance guidance law identification and filtering method is specifically prepared according to the following steps:

本专利提出的多模型PN制导律辨识滤波器在一个导弹突防场景中进行了仿真验证。为了验证该制导律辨识滤波器的性能,将其与PN制导律Kalman辨识滤波器的估计性能进行了对比。该滤波器未考虑控制器饱和情况,认为pursuer全程在PN制导律约束下运动。The multi-model PN guidance law identification filter proposed in this patent has been simulated and verified in a missile penetration scenario. In order to verify the performance of the guidance law identification filter, it is compared with the estimated performance of the PN guidance law Kalman identification filter. This filter does not consider the saturation of the controller, and considers that the pursuer moves under the constraint of the PN guidance law throughout the whole process.

在该导弹突防场景中,pursuer拦截evader。考虑到PN制导律是常见的制导律,假定evader事先知道pursuer采用PN制导律,但是不知道其在俯仰和偏航通道上的导航常数。evader可以获得pursuer的位置,使用多模型PN制导律辨识滤波器估计pursuer的PN制导律导航常数和在俯仰平面与偏航平面上的加速度。这些信息会传递给保护evader的防御导弹,防御导弹利用这些信息对pursuer进行拦截,以保护evader。pursuer的仿真配置在表1中描述,evader的仿真配置在表2中描述。In this missile penetration scenario, the pursuer intercepts the evader. Considering that the PN guidance law is a common guidance law, it is assumed that the evader knows in advance that the pursuer adopts the PN guidance law, but does not know its navigation constants on the pitch and yaw channels. The evader can obtain the position of the pursuer, and use the multi-model PN guidance law identification filter to estimate the pursuer's PN guidance law navigation constants and accelerations on the pitch plane and yaw plane. This information will be passed to the defensive missiles protecting the evader, and the defensive missiles will use this information to intercept the pursuer to protect the evader. The simulation configuration of the pursuer is described in Table 1, and the simulation configuration of the evader is described in Table 2.

表1:pursuer仿真配置信息Table 1: Pursuer simulation configuration information

表2:evader仿真配置信息Table 2: evader simulation configuration information

PN制导律Kalman辨识滤波器和饱和Kalman辨识滤波器的过程噪声矩阵和测量噪声矩阵设置为The process noise matrix and measurement noise matrix of PN guidance law Kalman identification filter and saturated Kalman identification filter are set as

QPNG=QSat=diag(10,10,10,20,20,20,2)Q PNG =Q Sat =diag(10,10,10,20,20,20,2)

RPNG=RSat=diag(10,10,10)R PNG =R Sat =diag(10,10,10)

图2-15是在该仿真场景下的仿真结果。图2和图3是evader常值机动情况下,滤波器MMAE_A对pursuer制导律的导航常数的估计值。MMAE_A和MMAE_B对导航常数的估计值都在pursuer退出饱和后很短时间内收敛到了真值。Pursuer在俯仰通道和偏航通道控制器退出饱和的时间分别是2.7s和4.6s。从图2和图3中可以看到在pursuer控制器处于饱和状态时,PN制导律Kalman辨识滤波器由于大的模型误差导致对导航常数的估计值不收敛到真值。在pursuer退出饱和之后很快对导航常数的估计值收敛到了真值,分别是俯仰通道导航常数为4,偏航通道导航常数为5。Figure 2-15 is the simulation result in this simulation scenario. Figure 2 and Figure 3 are the estimated value of the navigation constant of the pursuer guidance law by the filter MMAE_A in the case of evader constant maneuvering. Both MMAE_A and MMAE_B's estimates of navigation constants converge to the true value within a short time after the pursuer exits saturation. The time for Pursuer to exit saturation in pitch channel and yaw channel controller is 2.7s and 4.6s respectively. It can be seen from Figure 2 and Figure 3 that when the pursuer controller is in a saturated state, the PN guidance law Kalman identification filter does not converge to the true value of the estimated value of the navigation constant due to a large model error. After the pursuer exits saturation, the estimated value of the navigation constant converges to the true value, respectively, the pitch channel navigation constant is 4, and the yaw channel navigation constant is 5.

上述仿真过程考虑了在制导过程中pursuer控制器在初期产生饱和的可能性,使用了多模型滤波方法对pursuer制导律进行了辨识,对PN制导律下的导航常数和在俯仰平面和偏航平面加速度进行了估计。能更准确的对制导律进行辨识。在本发明的实施例中,设MMAE滤波器对pursuer加速度估计产生的均方误差为VarMMAE;全程使用PN制导律对pursuer加速度估计产生的均方误差为VarPN。设α=VarMMAEVarPN。当0<α<1,MMAE滤波器对pursuer加速度估计较说明全程使用PN制导律对pursuer加速度估计更精确;α越小,说明相对于全程使用PN制导律对pursuer加速度估计来说,MMAE对pursuer加速度估计越精确。The above simulation process takes into account the possibility of saturation of the pursuer controller at the initial stage in the guidance process, and uses the multi-model filtering method to identify the pursuer guidance law. Acceleration was estimated. The guidance law can be identified more accurately. In an embodiment of the present invention, it is assumed that the mean square error generated by the MMAE filter for the estimation of the pursuer acceleration is Var MMAE ; the mean square error generated by the PN guidance law for the estimation of the pursuer acceleration throughout the whole process is Var PN . Let α=Var MMAE Var PN . When 0<α<1, the MMAE filter is more accurate in estimating the pursuit acceleration than that using the PN guidance law in the whole process. Acceleration estimates are more accurate.

当evader无机动时,在pursuer饱和阶段,对pursuer加速度ay的估计时α=0.59;对pursuer加速度az的估计时α=0.56。When the evader is not maneuvering, in the saturation stage of the pursuer, the estimation of the pursuit acceleration a y is α=0.59; the estimation of the pursuit acceleration a z is α=0.56.

当evader进行常值机动时,在pursuer饱和阶段,对pursuer加速度ay的估计时α=0.55;对pursuer加速度az的估计时α=0.62。When the evader performs a constant maneuver, in the saturation stage of the pursuer, the estimation of the pursuit acceleration a y is α=0.55; the estimation of the pursuit acceleration a z is α=0.62.

当evader进行正弦机动时,在pursuer饱和阶段,对pursuer加速度ay的估计时α=0.56;对pursuer加速度az的估计时α=0.60。When the evader performs a sinusoidal maneuver, in the saturation stage of the pursuer, the estimation of the pursuit acceleration a y is α=0.56; the estimation of the pursuit acceleration a z is α=0.60.

图4、图5、图6、图7是在evader无机动情况下,整个制导过程都使用PN制导律Kalman辨识滤波器对pursuer进行辨识时对pursuer在俯仰通道和偏航通道过载的估计。在pursuer控制器处于饱和状态时,PN制导律Kalman辨识滤波器在俯仰和偏航通道上的过载估计值均与真值有较大的差异。而一旦pursuer控制器退出饱和,PN制导律Kalman辨识滤波器的估计值很快就收敛到了真值。本文提出的多模型估计方法考虑了饱和对加速度辨识的影响。综合考虑了MMAE_A和全程PN制导律辨识卡尔曼滤波器的估计值如图4所示。根据计算出的pursuer退出饱和的时间,在pursuer的控制器退出饱和之前使用MMAE_A滤波器的估计值,而在pursuer的控制器退出饱和之后使用全程PN制导律辨识卡尔曼滤波器的估计结果。Figure 4, Figure 5, Figure 6, and Figure 7 are the estimation of the pursuer's overload in the pitch channel and yaw channel when the evader is not maneuvering and the PN guidance law Kalman identification filter is used to identify the pursuer during the entire guidance process. When the pursuer controller is in a saturated state, the overload estimation values of the PN guidance law Kalman identification filter on the pitch and yaw channels are quite different from the true values. And once the pursuer controller exits saturation, the estimated value of the Kalman identification filter of the PN guidance law converges to the true value very quickly. The multi-model estimation method proposed in this paper takes into account the influence of saturation on acceleration identification. Considering MMAE_A and full-range PN guidance law, the estimated value of Kalman filter identification is shown in Fig. 4. According to the calculated time when the pursuer exits saturation, the estimated value of the MMAE_A filter is used before the pursuer controller exits saturation, and the estimated result of the Kalman filter is identified using the global PN guidance law after the pursuer controller exits saturation.

对比图4、图5、图6、图7看到在饱和阶段本文提出的滤波器对pursuer俯仰通道和偏航通道的过载估计精度均高于PN制导律Kalman辨识滤波器的估计精度。而在pursuer控制器退出饱和后,两种滤波器的估计精度相同。当evader执行常值机动或者正弦机动时,结果和evader不机动的结果一致。如图8、图9、图10、图11、图12、图13、图14、图15所示。Comparing Fig. 4, Fig. 5, Fig. 6, and Fig. 7, it can be seen that in the saturation stage, the overload estimation accuracy of the filter proposed in this paper for the pursuer pitch channel and yaw channel is higher than that of the PN guidance law Kalman identification filter. After the pursuit controller exits saturation, the estimation accuracy of the two filters is the same. When the evader performs a constant maneuver or a sine maneuver, the result is consistent with the result of the evader not maneuvering. As shown in Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14, Figure 15.

本发明还可有其它多种实施例,在不背离本发明精神及其实质的情况下,本领域技术人员当可根据本发明作出各种相应的改变和变形,但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。The present invention can also have other various embodiments, without departing from the spirit and essence of the present invention, those skilled in the art can make various corresponding changes and deformations according to the present invention, but these corresponding changes and deformations are all Should belong to the scope of protection of the appended claims of the present invention.

Claims (6)

1. a kind of proportional navigation law recognizes filtering method, it is characterised in that a kind of proportional navigation law recognizes filtering method Specifically follow the steps below:
Step 1:Pursuer has two kinds of motion models:When saturation motion model and controller during controller saturation are unsaturated PN Guidance Law motion models;State equation, the PN Guidance Law motion models that PN Guidance Laws motion model is established in pitch plane exist The state equation of the plane of yaw, controller saturation motion model pitch plane state equation and controller saturation motion model In the state equation of the plane of yaw;Evader is ballistic missile;Pursuer is interception guided missile;PN Guidance Laws are proportional guidance system Lead rule;
Step 2:According to state equation of the PN Guidance Laws motion model of step 1 in pitch plane, PN Guidance Law motion models State equation, controller saturation motion model in the plane of yaw move mould in the state equation and controller saturation of pitch plane Type is in the state equation of the plane of yaw, state-transition matrix, PN Guidance Law of the calculating PN Guidance Laws motion model in pitch plane Motion model the state-transition matrix of the plane of yaw, controller saturation motion model pitch plane state-transition matrix, State-transition matrix of the controller saturation motion model in the plane of yaw;
Step 3:Moved according to state-transition matrix, PN Guidance Law of the PN Guidance Laws motion model of step 2 in pitch plane Model is in the state-transition matrix of the state-transition matrix of the plane of yaw, controller saturation motion model in pitch plane, control Device saturation motion model is in the state-transition matrix of the plane of yaw, Kalman of the design PN Guidance Laws motion model in pitch plane Filter equation, PN Guidance Laws motion model exist in Kalman filter equation, the controller saturation motion model of the plane of yaw The Kalman filter equation and controller saturation motion model of pitch plane the plane of yaw Kalman filter equation, Kalman filter equation is Kalman filter equation;
Step 4:According to Kalman filter equation of the PN Guidance Laws motion model of step 3 in pitch plane, PN Guidance Laws Motion model the Kalman filter equation of the plane of yaw, controller saturation motion model pitch plane Kalman filter Device equation and controller saturation motion model are at the Kalman filter equation of the plane of yaw, calculating pursuer current times PN The posterior probability of Guidance Law motion model and the posterior probability of controller saturation motion model;
Step 5:Posterior probability and the saturation fortune of the pursuer current time PN Guidance Law motion models obtained according to step 4 The posterior probability of movable model, calculate pursuer from controller saturation when saturation motion model be switched to controller unsaturation when PN Guidance Law motion models at the time of;Saturation motion model when pursuer is from controller saturation is switched to controller not Before at the time of PN Guidance Law motion models during saturation using controller saturation motion model pitch plane Kalman filter The estimated result of device equation and controller saturation motion model in the Kalman filter equation of the plane of yaw;Otherwise PN systems are used Rule motion model is led in the Kalman filter equation and PN Guidance Law motion models of pitch plane in the plane of yaw
The estimated result of Kalman filter equation.
A kind of 2. proportional navigation law identification filtering method according to claim 1, it is characterised in that:In the step 1 Pursuer has two kinds of motion models:PN Guidance Laws fortune when saturation motion model and controller during controller saturation are unsaturated Movable model;PN Guidance Laws motion model is established in the state equation of pitch plane, PN Guidance Law motion models in the plane of yaw State equation, controller saturation motion model pitch plane state equation and controller saturation motion model in the plane of yaw State equation, pursuer is interception guided missile;Detailed process is:
In missile breakthrough scene, anti-ballistic missile of dashing forward will escape blocking for the interception guided missile of other side's missile defense systems transmitting Cut, be escape side, be known as evader, interception guided missile is the side of chasing, and is known as pursuer, and pursuer there are two kinds of motion models:Control PN Guidance Law motion models when saturation motion model and controller during device saturation processed are unsaturated;
The PN Guidance Laws equation of motion modeling of pursuer, rxIt is evader and pursuer relative positions in evader inertial system x-axis Under component;ryFor component of the evader and pursuer relative positions under evader inertial system y-axis;rzFor evader and Component of the pursuer relative positions under evader inertial system z-axis;Evader is ballistic missile, and pursuer is interception guided missile;
rx=xp-xe ry=yp-ye rz=zp-ze (1)
In formula, [xp,yp,zp]T[xe,ye,ze]TIt is the position of pursuer and evader under evader inertial systems respectively, T is Transposition, xpFor positions of the pursuer under evader inertial system x-axis, ypFor positions of the pursuer under evader inertial system y-axis Put, zpFor positions of the pursuer under evader inertial system z-axis, xeFor positions of the evader under evader inertial system x-axis, ye For positions of the evader under evader inertial system y-axis, zeFor positions of the evader under evader inertial system z-axis;
PN Guidance Laws motion model is in the state equation of pitch plane
<mrow> <mtable> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>x</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>y</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>y</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>x</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>y</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>y</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>y</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mi>&amp;tau;</mi> </mfrac> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>y</mi> </mrow> </msub> <mo>=</mo> <msub> <mi>N</mi> <mi>&amp;epsiv;</mi> </msub> <mfrac> <mrow> <mo>(</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>+</mo> <msub> <mi>r</mi> <mi>y</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>)</mo> <mo>(</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>-</mo> <msub> <mi>r</mi> <mi>y</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>)</mo> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mfrac> <mn>3</mn> <mn>2</mn> </mfrac> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mfrac> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>y</mi> </mrow> </msub> <mi>&amp;tau;</mi> </mfrac> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>N</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>&amp;epsiv;</mi> </msub> <mo>=</mo> <mn>0</mn> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow>
PN Guidance Laws motion model is in the state equation of the plane of yaw
<mrow> <mtable> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>x</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>z</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>z</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>x</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>z</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>z</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>z</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mi>&amp;tau;</mi> </mfrac> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>z</mi> </mrow> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <msub> <mi>N</mi> <mi>&amp;beta;</mi> </msub> <mi>&amp;tau;</mi> </mfrac> <mfrac> <mrow> <mo>(</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>+</mo> <msub> <mi>r</mi> <mi>z</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>)</mo> <mo>(</mo> <msub> <mi>r</mi> <mi>z</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>-</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>)</mo> </mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mfrac> <mn>3</mn> <mn>2</mn> </mfrac> </msup> </mfrac> <mo>-</mo> <mfrac> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>z</mi> </mrow> </msub> <mi>&amp;tau;</mi> </mfrac> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>N</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>&amp;beta;</mi> </msub> <mo>=</mo> <mn>0</mn> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </mrow>
In formula, vxFor component of the evader and pursuer relative velocities under evader inertial system x-axis, vyFor evader and Component of the pursuer relative velocities under evader inertial system y-axis, vzIt is evader and pursuer relative velocities in evader Component under inertial system z-axis;apxFor component of the acceleration under evader inertial system x-axis of pursuer, apyFor pursuer's Component of the acceleration under evader inertial system y-axis, apzFor component of the acceleration under evader inertial system z-axis of pursuer; aex, aeyAnd aezIt is component of the acceleration of evader under evader inertial systems, NεIt is navigation of the pursuer in pitch plane Constant;NβIt is navigation constants of the pursuer in the plane of yaw;τ is time constant,For rxFirst derivative,For rySingle order Derivative,For rzFirst derivative,For vxFirst derivative,For vyFirst derivative,For vzFirst derivative,For apx First derivative,For apyFirst derivative,For apzFirst derivative,For NεFirst derivative,For NβSingle order Derivative;
Assuming that evader obtains positions of the pursuer under inertial system, and then calculates the opposite position of itself, calculate opposite The measurement equation of the position of itself is
H=[rx,ry,rz]T (4)
Controller saturation motion model is in the state equation of pitch plane
<mrow> <mtable> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>x</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>y</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>y</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>x</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>y</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>y</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>y</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mi>&amp;tau;</mi> </mfrac> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>y</mi> </mrow> </msub> <mo>=</mo> <mn>0</mn> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>d</mi> <mo>&amp;CenterDot;</mo> </mover> <mn>1</mn> </msub> <mo>=</mo> <mn>0</mn> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>5</mn> <mo>)</mo> </mrow> </mrow>
Controller saturation motion model is in the state equation of the plane of yaw
<mrow> <mtable> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>x</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>r</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>z</mi> </msub> <mo>=</mo> <msub> <mi>v</mi> <mi>z</mi> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>x</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>x</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>v</mi> <mo>&amp;CenterDot;</mo> </mover> <mi>z</mi> </msub> <mo>=</mo> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>z</mi> </mrow> </msub> <mo>-</mo> <msub> <mi>a</mi> <mrow> <mi>e</mi> <mi>z</mi> </mrow> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <msub> <mi>a</mi> <mrow> <mi>p</mi> <mi>x</mi> </mrow> </msub> <mi>&amp;tau;</mi> </mfrac> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>a</mi> <mo>&amp;CenterDot;</mo> </mover> <mrow> <mi>p</mi> <mi>z</mi> </mrow> </msub> <mo>=</mo> <mn>0</mn> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mover> <mi>d</mi> <mo>&amp;CenterDot;</mo> </mover> <mn>2</mn> </msub> <mo>=</mo> <mn>0</mn> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>6</mn> <mo>)</mo> </mrow> </mrow>
In formula, d1And d2It is the occupied state of controller saturation motion model.
A kind of 3. proportional navigation law identification filtering method according to claim 2, it is characterised in that:In the step 2 According to the PN Guidance Laws motion model of step 1 in the state equation of pitch plane, PN Guidance Law motion models in the plane of yaw State equation, controller saturation motion model pitch plane state equation and controller saturation motion model in the plane of yaw State equation, calculate PN Guidance Laws motion model in the state-transition matrix of pitch plane, PN Guidance Law motion models inclined The navigate state-transition matrix of plane, controller saturation motion model is transported in the state-transition matrix of pitch plane, controller saturation State-transition matrix of the movable model in the plane of yaw;Detailed process is:
PN Guidance Laws motion model is in the state-transition matrix of pitch plane
<mrow> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mn>1</mn> </msubsup> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mrow> <mn>1</mn> <mo>-</mo> <mi>T</mi> <mo>/</mo> <mi>&amp;tau;</mi> </mrow> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <msub> <mi>a</mi> <mn>1</mn> </msub> </mtd> <mtd> <msub> <mi>a</mi> <mn>2</mn> </msub> </mtd> <mtd> <msub> <mi>a</mi> <mn>3</mn> </msub> </mtd> <mtd> <msub> <mi>a</mi> <mn>4</mn> </msub> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mrow> <mn>1</mn> <mo>-</mo> <mi>T</mi> <mo>/</mo> <mi>&amp;tau;</mi> </mrow> </mtd> <mtd> <msub> <mi>a</mi> <mn>5</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>1</mn> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;epsiv;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>3</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>-</mo> <mn>2</mn> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>r</mi> <mi>y</mi> </msub> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <mn>2</mn> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>r</mi> <mi>y</mi> </msub> <msubsup> <mi>v</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>-</mo> <mn>5</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>+</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>-</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>5</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>8</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>2</mn> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;epsiv;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>-</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>+</mo> <mn>5</mn> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>r</mi> <mi>y</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>-</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <msubsup> <mi>v</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>-</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>3</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>5</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>9</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>3</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;epsiv;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>-</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>r</mi> <mi>y</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>-</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>3</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>10</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>4</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;epsiv;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>+</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>r</mi> <mi>y</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>-</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>3</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>11</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>5</mn> </msub> <mo>=</mo> <mfrac> <mrow> <mi>T</mi> <mrow> <mo>(</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>+</mo> <msub> <mi>r</mi> <mi>y</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>)</mo> </mrow> <mrow> <mo>(</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>y</mi> </msub> <mo>-</mo> <msub> <mi>r</mi> <mi>y</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>y</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>3</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>12</mn> <mo>)</mo> </mrow> </mrow>
In formula, T is measurement period, and unit is the second;
PN Guidance Laws motion model is in the state-transition matrix of the plane of yaw
<mrow> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mn>2</mn> </msubsup> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mrow> <mn>1</mn> <mo>-</mo> <mi>T</mi> <mo>/</mo> <mi>&amp;tau;</mi> </mrow> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <msub> <mi>a</mi> <mn>6</mn> </msub> </mtd> <mtd> <msub> <mi>a</mi> <mn>7</mn> </msub> </mtd> <mtd> <msub> <mi>a</mi> <mn>8</mn> </msub> </mtd> <mtd> <msub> <mi>a</mi> <mn>9</mn> </msub> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mrow> <mn>1</mn> <mo>-</mo> <mi>T</mi> <mo>/</mo> <mi>&amp;tau;</mi> </mrow> </mtd> <mtd> <msub> <mi>a</mi> <mn>10</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>13</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>6</mn> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;beta;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>3</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>-</mo> <mn>2</mn> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>r</mi> <mi>z</mi> </msub> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <mn>2</mn> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>r</mi> <mi>z</mi> </msub> <msubsup> <mi>v</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>-</mo> <mn>5</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>+</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>-</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>5</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>14</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>7</mn> </msub> <mo>=</mo> <mo>-</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;beta;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>-</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>3</mn> </msubsup> <msubsup> <mi>v</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>+</mo> <mn>5</mn> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>r</mi> <mi>z</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>-</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <msubsup> <mi>v</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <msubsup> <mi>v</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>-</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>3</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>5</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>15</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>8</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;beta;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>-</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>r</mi> <mi>z</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>-</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>3</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>16</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>9</mn> </msub> <mo>=</mo> <mfrac> <mrow> <msub> <mi>TN</mi> <mi>&amp;beta;</mi> </msub> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>+</mo> <mn>2</mn> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>r</mi> <mi>z</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>-</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>3</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>17</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>a</mi> <mn>10</mn> </msub> <mo>=</mo> <mfrac> <mrow> <mi>T</mi> <mrow> <mo>(</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>+</mo> <msub> <mi>r</mi> <mi>z</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>)</mo> </mrow> <mrow> <mo>(</mo> <msub> <mi>r</mi> <mi>x</mi> </msub> <msub> <mi>v</mi> <mi>z</mi> </msub> <mo>-</mo> <msub> <mi>r</mi> <mi>z</mi> </msub> <msub> <mi>v</mi> <mi>x</mi> </msub> <mo>)</mo> </mrow> </mrow> <mrow> <msup> <mrow> <mo>(</mo> <msubsup> <mi>r</mi> <mi>x</mi> <mn>2</mn> </msubsup> <mo>+</mo> <msubsup> <mi>r</mi> <mi>z</mi> <mn>2</mn> </msubsup> <mo>)</mo> </mrow> <mrow> <mn>3</mn> <mo>/</mo> <mn>2</mn> </mrow> </msup> <mi>&amp;tau;</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>18</mn> <mo>)</mo> </mrow> </mrow>
Controller saturation motion model is in the state-transition matrix of pitch plane
<mrow> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mn>3</mn> </msubsup> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mrow> <mn>1</mn> <mo>-</mo> <mi>T</mi> <mo>/</mo> <mi>&amp;tau;</mi> </mrow> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>19</mn> <mo>)</mo> </mrow> </mrow>
Controller saturation motion model is in the state-transition matrix of the plane of yaw
<mrow> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mn>4</mn> </msubsup> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mi>T</mi> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mrow> <mn>1</mn> <mo>-</mo> <mi>T</mi> <mo>/</mo> <mi>&amp;tau;</mi> </mrow> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>1</mn> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>20</mn> <mo>)</mo> </mrow> <mo>.</mo> </mrow>
A kind of 4. proportional navigation law identification filtering method according to claim 3, it is characterised in that:In the step 3 It is flat in yaw in state-transition matrix, the PN Guidance Laws motion model of pitch plane according to the PN Guidance Laws motion model of step 2 State-transition matrix, controller saturation movement mould of the state-transition matrix, controller saturation motion model in face in pitch plane Type the plane of yaw state-transition matrix, design PN Guidance Laws motion model pitch plane Kalman filter equation, PN Guidance Laws motion model is in the Kalman filter equation of the plane of yaw, controller saturation motion model in pitch plane Kalman filter equation and controller saturation motion model are in the Kalman filter equation of the plane of yaw, Kalman filter Equation is Kalman filter equation;Detailed process is:
PN Guidance Laws recognize wave filter and saturation identification wave filter uses EKF, and EKF is extended Kalman filter, according to step Rapid two state-transition matrixes being calculated,
PN Guidance Laws motion model is in the Kalman filter equation of pitch plane
<mrow> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>=</mo> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>+</mo> <msubsup> <mi>K</mi> <mi>k</mi> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>-</mo> <mi>H</mi> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>21</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>HP</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>R</mi> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msup> <mo>)</mo> </mrow> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>22</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>Q</mi> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>23</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <mrow> <mo>(</mo> <mi>I</mi> <mo>-</mo> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mi>H</mi> <mo>)</mo> </mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>24</mn> <mo>)</mo> </mrow> </mrow>
PN Guidance Laws motion model is in the Kalman filter equation of the plane of yaw
<mrow> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>=</mo> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>+</mo> <msubsup> <mi>K</mi> <mi>k</mi> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>-</mo> <mi>H</mi> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>25</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>HP</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>R</mi> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msup> <mo>)</mo> </mrow> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>26</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>Q</mi> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>27</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <mrow> <mo>(</mo> <mi>I</mi> <mo>-</mo> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <mi>H</mi> <mo>)</mo> </mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </msubsup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>28</mn> <mo>)</mo> </mrow> </mrow>
Controller saturation motion model is in the Kalman filter equation of pitch plane
<mrow> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>=</mo> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>+</mo> <msubsup> <mi>K</mi> <mi>k</mi> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>-</mo> <mi>H</mi> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>29</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>HP</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>R</mi> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msup> <mo>)</mo> </mrow> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>30</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>Q</mi> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>31</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <mrow> <mo>(</mo> <mi>I</mi> <mo>-</mo> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <mi>H</mi> <mo>)</mo> </mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </msubsup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>32</mn> <mo>)</mo> </mrow> </mrow>
Controller saturation motion model is in the Kalman filter equation of the plane of yaw
<mrow> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>=</mo> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>+</mo> <msubsup> <mi>K</mi> <mi>k</mi> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <mrow> <mo>(</mo> <msub> <mi>Z</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>-</mo> <mi>H</mi> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>33</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>HP</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mi>H</mi> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>R</mi> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msup> <mo>)</mo> </mrow> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>34</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <msup> <mrow> <mo>(</mo> <msubsup> <mi>&amp;Phi;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>+</mo> <msup> <mi>Q</mi> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>35</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <mrow> <mo>(</mo> <mi>I</mi> <mo>-</mo> <msubsup> <mi>K</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <mi>H</mi> <mo>)</mo> </mrow> <msubsup> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </msubsup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>36</mn> <mo>)</mo> </mrow> </mrow>
In formula,It is the estimate to k+1 moment states;It is the predicted value to k+1 moment states;
WithIt is the filter gain of the PN Guidance Laws motion model in pitch plane at k moment and k+1 moment respectively;
WithIt is the filter gain of the PN Guidance Laws motion model in the plane of yaw at k moment and k+1 moment respectively;
WithIt is the filter gain of the controller saturation motion model in pitch plane at k moment and k+1 moment respectively;
WithIt is the filter gain of the controller saturation motion model in the plane of yaw at k moment and k+1 moment respectively;
For k+1 moment PN Guidance Law motion model the state forecast of pitch plane covariance matrix;
For k+1 moment PN Guidance Law motion model the state forecast of the plane of yaw covariance matrix;
For k+1 moment controller saturation motion models the state forecast of pitch plane covariance matrix;
For k+1 moment controller saturation motion models the state forecast of the plane of yaw covariance matrix;
R(1)For PN Guidance Laws motion model pitch plane measurement noise matrix;
R(2)For PN Guidance Laws motion model the plane of yaw measurement noise matrix;
R(3)Measurement noise matrix of the device saturation motion model in pitch plane in order to control;
R(4)Measurement noise matrix of the device saturation motion model in the plane of yaw in order to control;
For PN Guidance Laws motion model k moment to the k+1 moment of pitch plane state-transition matrix;
For PN Guidance Laws motion model k moment to the k+1 moment of the plane of yaw state-transition matrix;
State-transition matrix of the device saturation motion model at k moment to the k+1 moment of pitch plane in order to control;
State-transition matrix of the device saturation motion model at k moment to the k+1 moment of the plane of yaw in order to control;
Q(1)For PN Guidance Laws motion model pitch plane process noise matrix;
Q(2)For PN Guidance Laws motion model the plane of yaw process noise matrix;
Q(3)Process noise matrix of the device saturation motion model in pitch plane in order to control;
Q(4)Process noise matrix of the device saturation motion model in the plane of yaw in order to control;
For PN Guidance Laws motion model pitch plane k moment state estimation covariance matrixes;
For PN Guidance Laws motion model the plane of yaw k moment state estimation covariance matrixes;
K moment state estimation covariance matrix of the device saturation motion model in pitch plane in order to control;
K moment state estimation covariance matrix of the device saturation motion model in the plane of yaw in order to control;
For PN Guidance Laws motion model pitch plane k+1 moment state estimation covariance matrixes;
For PN Guidance Laws motion model the plane of yaw k+1 moment state estimation covariance matrixes;
K+1 moment state estimation covariance matrix of the device saturation motion model in pitch plane in order to control;
K+1 moment state estimation covariance matrix of the device saturation motion model in the plane of yaw in order to control;
For PN Guidance Laws motion model pitch plane the k moment to the state-transition matrix at k+1 moment transposition;
For PN Guidance Laws motion model the plane of yaw the k moment to the state-transition matrix at k+1 moment transposition;
Device saturation motion model turns at the k moment of pitch plane to the state-transition matrix at k+1 moment in order to control Put;
Device saturation motion model turns at the k moment of the plane of yaw to the state-transition matrix at k+1 moment in order to control Put;
Zk+1For the measured value at k+1 moment;I is 7 × 7 unit matrixs;H is calculation matrix, H=[1 10000 0], HTFor H Transposition, k is sampling instant.
A kind of 5. proportional navigation law identification filtering method according to claim 4, it is characterised in that:In the step 4 According to the PN Guidance Laws motion model of step 3 in the Kalman filter equation of pitch plane, PN Guidance Law motion models inclined The Kalman filter equation of plane, controller saturation motion model navigate in the Kalman filter equation of pitch plane and control Device saturation motion model is in the Kalman filter equation of the plane of yaw, calculating pursuer current times PN Guidance Law movement mould The posterior probability of type and the posterior probability of saturated controller motion model;Detailed process is:
If model set is M={ Mj| j=1 ..., r }, in formula, r is the number of model, r=2, M1Represent PN Guidance Laws movement mould Type, M2Represent controller saturation motion model, be P { M in the prior probability of k moment models jj|Zk-1, after k moment models j It is P { M to test probabilityj|Zk, the measured value at kth moment is z (k), and the measured value sequence definition at preceding k moment is
<mrow> <msup> <mi>Z</mi> <mi>k</mi> </msup> <mover> <mo>=</mo> <mi>&amp;Delta;</mi> </mover> <mo>{</mo> <mi>z</mi> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>|</mo> <mi>k</mi> <mo>=</mo> <mn>1</mn> <mo>,</mo> <mo>...</mo> <mo>,</mo> <mi>k</mi> <mo>}</mo> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>37</mn> <mo>)</mo> </mrow> </mrow>
According to Bayes' theorem,
<mrow> <mi>P</mi> <mo>{</mo> <msub> <mi>M</mi> <mi>j</mi> </msub> <mo>|</mo> <msup> <mi>Z</mi> <mi>k</mi> </msup> <mo>}</mo> <mo>=</mo> <mi>P</mi> <mo>{</mo> <msub> <mi>M</mi> <mi>j</mi> </msub> <mo>|</mo> <mi>z</mi> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>,</mo> <msup> <mi>Z</mi> <mrow> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>}</mo> <mo>=</mo> <mfrac> <mrow> <mi>p</mi> <mrow> <mo>(</mo> <mi>z</mi> <mo>(</mo> <mi>k</mi> <mo>)</mo> <mo>|</mo> <msup> <mi>Z</mi> <mrow> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>,</mo> <msub> <mi>M</mi> <mi>j</mi> </msub> <mo>)</mo> </mrow> <mi>P</mi> <mo>{</mo> <msub> <mi>M</mi> <mi>j</mi> </msub> <mo>|</mo> <msup> <mi>Z</mi> <mrow> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>}</mo> </mrow> <mrow> <munderover> <mi>&amp;Sigma;</mi> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>r</mi> </munderover> <mi>p</mi> <mrow> <mo>(</mo> <mi>z</mi> <mo>(</mo> <mi>k</mi> <mo>)</mo> <mo>|</mo> <msup> <mi>Z</mi> <mrow> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>,</mo> <msub> <mi>M</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mi>P</mi> <mo>{</mo> <msub> <mi>M</mi> <mi>i</mi> </msub> <mo>|</mo> <msup> <mi>Z</mi> <mrow> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mo>}</mo> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>38</mn> <mo>)</mo> </mrow> </mrow>
In formula, P { Mi|Zk-1Be model i prior probability, p (z (k) | Zk-1,Mi) for model i in the likelihood function at k moment, P {Mj|z(k),Zk-1It is k moment model MsjPosterior probability, p (z (k) | Zk-1,Mj) it is likelihood functions of the model j at the k moment;Cause The posterior probability of this model j can be calculated by the prior probability and likelihood function of all models, and the likelihood function of model i is The distribution function that the model respective filter newly ceases;It is assumed that new breath obeys the normal distribution that average is zero, if the variance at its k moment For Si(k), i.e.,
p(z(k)|Zk-1,Mi)=N (0, Sj(k)) (39)
In formula, N (0, Sj(k)) average is represented as 0, variance Sj(k) probability density function of normal distribution;Sj(k) it is model j The variance newly ceased, i.e.,
<mrow> <mi>v</mi> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>z</mi> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>-</mo> <mi>H</mi> <msub> <mover> <mi>X</mi> <mo>^</mo> </mover> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>,</mo> <mi>k</mi> </mrow> </msub> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>40</mn> <mo>)</mo> </mrow> </mrow>
Sj(k)=E [v (k) vT(k)] (41)
In formula, v (k) ceases to be new, vT(k) be v (k) transposition, E [] is mathematic expectaion function.
A kind of 6. proportional navigation law identification filtering method according to claim 5, it is characterised in that:In the step 5 The posterior probability of pursuer current time PN Guidance Law motion models and the posteriority of saturation motion model obtained according to step 4 Probability, calculate pursuer from controller saturation when saturation motion model be switched to controller unsaturation when PN Guidance Laws fortune At the time of movable model;Saturation motion model when pursuer is from controller saturation is switched to PN systems during controller unsaturation Controller saturation motion model is used before at the time of leading rule motion model in the Kalman filter equation of pitch plane and control Estimated result of the device saturation motion model in the Kalman filter equation of the plane of yaw;Otherwise PN Guidance Law motion models are used Pitch plane Kalman filter equation and PN Guidance Laws motion model the plane of yaw Kalman filter equation Estimated result.Detailed process is:
The state of two stage pursuer is estimated respectively using two MMAE wave filters;By the two MMAE wave filters MMAE_A and MMAE_B wave filters are referred to as, wherein MMAE_A wave filters estimation first stage, i.e. pursuer is in controller State during saturation;And MMAE_B wave filters estimation second stage, i.e. pursuer are in state during controller unsaturation;Separately It is outer need to estimate pursuer from controller saturation when saturation motion model be switched to controller unsaturation when PN Guidance Laws At the time of motion model, so that it is determined that representing current pursuer using the estimate of which MMAE wave filter in certain moment t State;
If pursuer is introduced into state during controller saturation, the motion model of pursuer whole process is consistent, MMAE_A It is consistent with the result that MMAE_B is estimated;
Pursuer from controller saturation when saturation motion model be switched to controller unsaturation when PN Guidance Law motion models At the time of method of estimation detailed process be:
Contain two models in each MMAE wave filters:PN Guidance Laws motion model and controller when controller is unsaturated are satisfied With when saturation motion model;The posterior probability of two models of present filter is calculated according to multiple model filtering method, then According to the posterior probability of two models of present filter, saturation rank when whether pursuer is in controller saturation at present is determined Section;MMAE_B wave filters check saturation exit point, the i.e. time at interval of check_period long by the way of scanning, Check_period is round of visits, the motion model of current pursuer is judged according to formula (48), it is assumed that controlled in MMAE wave filters The posterior probability of saturation motion model during PN Guidance Laws motion model and controller saturation when device processed is unsaturated is respectively P {MPNGAnd P { MSat, then deterministic process is:
<mrow> <mtable> <mtr> <mtd> <mrow> <mi>mod</mi> <mi>e</mi> <mo>=</mo> <mi>P</mi> <mi>N</mi> <mi>G</mi> </mrow> </mtd> <mtd> <mrow> <mi>i</mi> <mi>f</mi> <mi> </mi> <mi>P</mi> <mo>{</mo> <msub> <mi>M</mi> <mrow> <mi>P</mi> <mi>N</mi> <mi>G</mi> </mrow> </msub> <mo>}</mo> <mo>-</mo> <mi>P</mi> <mo>{</mo> <msub> <mi>M</mi> <mrow> <mi>S</mi> <mi>a</mi> <mi>t</mi> </mrow> </msub> <mo>}</mo> <mo>&amp;GreaterEqual;</mo> <mi>&amp;epsiv;</mi> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mi>mod</mi> <mi>e</mi> <mo>=</mo> <mi>S</mi> <mi>a</mi> <mi>t</mi> <mi>u</mi> <mi>r</mi> <mi>a</mi> <mi>t</mi> <mi>i</mi> <mi>o</mi> <mi>n</mi> </mrow> </mtd> <mtd> <mrow> <mi>e</mi> <mi>l</mi> <mi>s</mi> <mi>e</mi> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>48</mn> <mo>)</mo> </mrow> </mrow>
In formula, ε is threshold value;As P { MPNG}-P{MSatDuring } >=ε, then the model that current pursuer is used is PN Guidance Law mould Type, judges that pursuer has exited saturation state;Otherwise restart and initialize MMAE_B wave filters, check_period is set It is set to 0.5 second.
CN201510822999.2A 2015-11-23 2015-11-23 A kind of proportional navigation law recognizes filtering method Active CN105446352B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201510822999.2A CN105446352B (en) 2015-11-23 2015-11-23 A kind of proportional navigation law recognizes filtering method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201510822999.2A CN105446352B (en) 2015-11-23 2015-11-23 A kind of proportional navigation law recognizes filtering method

Publications (2)

Publication Number Publication Date
CN105446352A CN105446352A (en) 2016-03-30
CN105446352B true CN105446352B (en) 2018-04-24

Family

ID=55556672

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201510822999.2A Active CN105446352B (en) 2015-11-23 2015-11-23 A kind of proportional navigation law recognizes filtering method

Country Status (1)

Country Link
CN (1) CN105446352B (en)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108009358B (en) * 2017-12-01 2021-06-15 哈尔滨工业大学 Three-dimensional guidance law identification and filtering method based on IMM_UKF
CN108052112B (en) * 2017-12-01 2020-10-02 哈尔滨工业大学 Multi-aircraft threat degree obtaining method based on PN guidance law identification
CN110187640B (en) * 2019-06-29 2022-04-29 东南大学 Multi-missile cooperative combat guidance law design method for maneuvering target and allowable communication time lag
CN110686564B (en) * 2019-10-15 2020-08-11 北京航空航天大学 Infrared semi-strapdown seeker guidance method and system
CN112731965B (en) * 2020-12-17 2022-04-08 哈尔滨工业大学 A guidance method based on target maneuver identification

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
USH1980H1 (en) * 1996-11-29 2001-08-07 The United States Of America As Represented By The Secretary Of The Air Force Adaptive matched augmented proportional navigation
CN104266546A (en) * 2014-09-22 2015-01-07 哈尔滨工业大学 Sight line based finite time convergence active defense guidance control method

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
USH1980H1 (en) * 1996-11-29 2001-08-07 The United States Of America As Represented By The Secretary Of The Air Force Adaptive matched augmented proportional navigation
CN104266546A (en) * 2014-09-22 2015-01-07 哈尔滨工业大学 Sight line based finite time convergence active defense guidance control method

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
cooperative multiple-model adaptive guidance for an aircraft defending missile;Vitaly Shafeman etc,;《journal of guidance,control and dynamics》;20101231;第1801-1813页 *
Pursuer Identification and Time-To-Go Estimation using Passive Measurements from an Evader;L. LIN etc,;《IEEE TRANSACTIONS ON AEROSPACE AND ELECTRONIC SYSTEMS》;20050131;第190-201页 *
一种导弹导引律及参数的非线性MMAE辨识方法;王小平等;《飞行力学》;20151031;第430-434页 *
寻的导弹制导中距离测量可能错误时的多模型滤波;周荻等;《宇航学报》;20051031;第66-69页 *
比例导引下的机动目标自适应滤波算法研究;邹冈等;《计算机仿真》;20070930;第26-29页 *

Also Published As

Publication number Publication date
CN105446352A (en) 2016-03-30

Similar Documents

Publication Publication Date Title
CN105446352B (en) A kind of proportional navigation law recognizes filtering method
Lee et al. Polynomial guidance laws considering terminal impact angle and acceleration constraints
Kim et al. Time-to-go polynomial guidance with trajectory modulation for observability enhancement
CN109472418B (en) Prediction and optimization method of maneuvering target state based on Kalman filter
CN107941087A (en) A kind of superb steady gliding reentry guidance method of high lift-drag ratio based on resistance profiles
CN105202972B (en) Multi-missile cooperative engagement guidance method based on model predictive control technique
CN112819303B (en) Aircraft tracking efficiency evaluation method and system based on PCE proxy model
CN108267731B (en) Construction method and application of UAV target tracking system
CN107255924A (en) Method for extracting guidance information of strapdown seeker through volume Kalman filtering based on dimension expansion model
CN107618678A (en) Attitude control information consolidation method of estimation under attitude of satellite angular deviation
CN109001699B (en) A Tracking Method Based on Constraints of Destination Information with Noisy
KR101658464B1 (en) Real-time prediction method of impact point of guided missile
CN103528449B (en) Missile formation control method based on disturbance observer and finite time control
CN112099348B (en) Collision angle control guidance method based on observer and global sliding mode
CN115047769A (en) Unmanned combat platform obstacle avoidance-arrival control method based on constraint following
CN111121770A (en) Interactive multi-missile multi-model flight path fusion method
CN108052112A (en) Multi-aircraft Threat acquisition methods based on the identification of PN Guidance Laws
CN109780933B (en) Dynamic target prediction guidance method for individual-soldier guided rocket
CN114740760B (en) Semi-physical simulation method and system for strapdown guided missile
Kumar et al. Target tracking using adaptive Kalman Filter
Kumar et al. Adaptive extended kalman filter for ballistic missile tracking
Zakharin et al. Concept of navigation system design of UAV
CN117889868B (en) Missile position accurate estimation method integrating infrared seeker information
CN108009358A (en) Three-dimensional guidance rule identification filtering method based on IMM_UKF
ATES Lyapunov based nonlinear impact angle guidance law for stationary targets

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant