A space-based radar air target nonlinear robust tracking method

Abstract

A space-based radar air target nonlinear robust tracking method, the main steps include: constructing a space-based radar target observation geometry model, establishing the conversion relationship between the target state and the measurement; VB approximation method is introduced into the PHD recursive framework, the Gaussian closed solution of PHD filter under unknown measurement noise covariance is derived to solve the problem of measurement noise uncertainty under space-based observation conditions; On this basis, the cubature Kalman filter is introduced into the VB-PHD filter based on Gaussian mixture model, and the third-order spherical radial cubature rule is used to complete the PHD recursive density function integration operation under the nonlinear system model, and the tracking accuracy under strong nonlinear observation conditions is improved. The method can effectively improve the filtering accuracy and improve the robustness of space-based radar air target tracking.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of radar target tracking, and particularly relates to a space-based radar air target nonlinear robust tracking method. BACKGROUND

[0002] The space-based radar takes a satellite as a load platform, can be free from the restriction of territory, territorial sea and climate, can discover and continuously and stably track a low-altitude air target which is on the attack, and can effectively improve early warning capability. However, due to the influence of factors such as high-speed movement of the platform, long-distance detection, perturbation interference and nonlinear observation, the target tracking of the space-based radar is brought severe challenges, mainly manifested as follows: 1) the space-based radar load platform and the air target are in different inertial systems, causing the radar observation equation to be nonlinear; 2) under the condition of long-distance observation, the measurement noise uncertainty increases, resulting in difficulty in high-precision stable tracking.

[0003] The target tracking algorithm based on a random finite set (RFS) utilizes the finite set statistics (FISST) theory to construct an optimal multi-target Bayesian filter under a Bayesian framework, provides a unified and complete theoretical basis for a multi-target tracking problem, and is suitable for solving a time-varying multi-target filtering problem in a complex environment. Among them:

[0004] The probability hypothesis density filter (PHD) is a first-order moment approximation of the optimal multi-target Bayesian filter, has been widely applied due to the advantage of low computational complexity. In the document “The Gaussian mixture probability hypothesis density filter”, an analytic solution of PHD recursion under a linear Gaussian dynamic model is derived, and the extended Kalman filter (EKF) and unscented Kalman filter (UKF) two methods are used to extend it to a nonlinear dynamic model, but it is only applicable to a weakly nonlinear and known measurement noise scene. For an actual space-based radar target tracking scene, the observation equation presents strong nonlinearity, and the measurement noise is unknown and time-varying.

[0005] The variational Bayesian (VB) method can be used to estimate the target state under unknown measurement noise variance, and the main idea is to parameterize the joint posterior distribution of the measurement noise covariance and the target state, and give its recursive form. Based on this, the document "An improved PHD filter based on variational Bayesian method for multi-target tracking" introduces the variational Bayesian approximation method into the PHD recursive filter to realize the joint estimation of the posterior distribution of the multi-target state and the measurement noise variance, and derives the closed-form solution of the recursive equation based on the inverse gamma and Gaussian mixture distribution under the linear Gaussian model, but it cannot be applied to the nonlinear observation model. SUMMARY

[0006] For the problem of strong nonlinear observation and unknown measurement noise variance in space-based radar air target tracking, the present disclosure provides a space-based radar air target nonlinear robust tracking method based on a VB-PHD filter:

[0007] First, a space-based radar target observation geometric model is constructed to establish the conversion relationship between the target state and the measurement;

[0008] Then, the VB approximation method is introduced into the PHD recursive framework, and the Gaussian closed solution of the PHD filter under unknown measurement noise covariance is derived to solve the problem of measurement noise uncertainty under space-based observation conditions;

[0009] On this basis, the Cubature Kalman Filter (CKF) is introduced into the VB-PHD filter based on the Gaussian Mixture (GM), and the third-order spherical radial cubature rule is used to complete the integration operation of the PHD recursive density function under the nonlinear system model, thereby improving the tracking accuracy under strong nonlinear observation conditions.

[0010] Compared with the prior art, the present disclosure has the following beneficial effects: (1) a space-based radar air target observation geometric model is established, and the VB approximation method is introduced into the PHD filter framework, thereby solving the problem of measurement noise uncertainty under space-based observation conditions; (2) the Cubature Kalman Filter algorithm is introduced into the VB-GM-PHD filter framework, and the third-order spherical radial cubature rule is used to complete the integration operation of the density function under the nonlinear system model in the recursive process of the PHD filter, thereby improving the filtering accuracy under strong nonlinear observation conditions; (3) the introduction of the variational processing can effectively improve the filtering accuracy and enhance the robustness of the space-based radar air target tracking. BRIEF DESCRIPTION OF DRAWINGS

[0011] The above and other objects, features and advantages of this disclosure will become more apparent from the more detailed description of exemplary embodiments of this disclosure taken in conjunction with the accompanying drawings, in which the same reference numerals generally represent the same components.

[0012] Figure 1 This is a schematic diagram of the geometric model for space-based radar observation.

[0013] Figure 2 This is a schematic diagram of a space-based radar target tracking scenario in the ECEF coordinate system in an exemplary embodiment.

[0014] Figure 3 The distribution of measurement points in the ECEF coordinate system in an exemplary embodiment;

[0015] Figure 4 This is a comparison chart of OSPA distances for multiple targets in an exemplary embodiment;

[0016] Figure 5 A comparison chart of RMSE tracking for target 1 in an exemplary embodiment;

[0017] Figure 6 This is a flowchart illustrating an exemplary embodiment of the present disclosure. Detailed Implementation

[0018] Preferred embodiments of the present disclosure will now be described in more detail with reference to the accompanying drawings. While preferred embodiments of the present disclosure are shown in the drawings, it should be understood that the present disclosure may be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that the present disclosure will be thorough and complete, and will fully convey the scope of the present disclosure to those skilled in the art.

[0019] This disclosure proposes a nonlinear robust tracking method for airborne targets using space-based radar based on a VB-PHD filter. An exemplary embodiment is shown in the attached flowchart. Figure 6 As shown, the specific steps include:

[0020] Step 1: First, model the geometric model of space-based radar observation and establish the conversion relationship between target state and measurement.

[0021] Assuming the Earth is a regular, uniform sphere, and the radar carrier satellite platform is in a circular near-Earth orbit, radar airborne target tracking system models are typically described using state equations and observation equations. The target state equation describes the target's motion state and is established in an Earth-Centered and Earth-Fixed (ECEF) Cartesian coordinate system. The origin O is located at the Earth's center, the XOY plane is located at the equatorial plane, the X-axis points to the Greenwich Meridian, and the Z-axis points to the North Pole. The observation equation describes the radar's real-time measurements of the target and is usually established in a polar coordinate system. Assuming the origin O′ of the space-based radar's body coordinate system O′X′Y′Z′ is the radar's center, the O′X′ axis is the extension of the line connecting the Earth's center O and the radar center O′, the O′Z′ axis lies in the plane formed by the O′X′ axis and the Earth's rotation axis and is perpendicular to the O′X′ axis, and O′Y′ is determined by the O′Z′ axis and the O′X′ axis according to the right-hand screw law, as shown below. Figure 1 As shown.

[0022] Let the position of the satellite in the geocentric fixed rectangular coordinate system be (x s ,y s ,z s The target coordinates of the aircraft are (x...). a ,y a ,z a In the radar payload satellite body coordinate system, the measured value of the aircraft target is... The transformation relationship between the radar body coordinate system measurement values ​​and the rectangular coordinate system state values ​​is as follows:

[0023]

[0024] The transformation relationship between the aircraft target and the radar body coordinate system from the geocentric fixed rectangular coordinate system is as follows:

[0025]

[0026] Based on formulas (1) and (2), it can be concluded that the state transition relationship exhibits strong nonlinearity. Therefore, this invention introduces the spherical radial volume rule into the VB-PHD filter to improve the nonlinear tracking accuracy.

[0027] Step 2: Introduce the VB approximation method into the PHD recursive framework to derive the Gaussian closed-form solution of the PHD filter under the condition of unknown measurement noise covariance, so as to solve the problem of measurement noise uncertainty under space-based observation conditions.

[0028] The multi-objective state at time k can be described using RFS as follows:

[0029]

[0030] In the formula: B kis the set of target states for new students; ∪ is the union operation of sets.

[0031] The target's motion state is represented as:

[0032]

[0033] In the formula: x k,i y k,i z k,i The spatial position of the target in the ECEF coordinate system; The target speed of motion.

[0034] Assuming that the target observation process and the clutter generation process are independent of each other, the multi-target measurement Z received by the sensor... k The union of measurements and clutter generated by the target can be represented as follows:

[0035]

[0036] In the formula, K k This is a set of false alarm clutter measurements.

[0037] The VB-PHD filtering recursive process is as follows:

[0038] ●Prediction

[0039] Assume that at time k-1, the joint posterior density of the variances of the multi-target states and measurement noise is v. k-1 Given (x, R) and that the two are independent of each other, the joint prediction strength can be derived as follows:

[0040] v k|k-1 (x,R)=∫p S,k (ζ)f k|k-1 (x|ζ)f k|k-1 (R|Σ)×v k-1 (ζ,Σ)dζdΣ+γ k (x,R) (6)

[0041] Where, γ k It refers to the intensity of births occurring simultaneously.

[0042] ·renew

[0043] Assume that at time k-1, the joint prediction density of the multi-target state and measurement noise variance is v. k|k-1 If (x,R), then the joint posterior strength can be expressed as:

[0044]

[0045] Where, p D,k (x) represents the detection probability, κ k (z) represents the Poisson clutter density.

[0046] In the above formula, g k (z|x,R) represents the objective likelihood function, and R k Unknown and time-varying, making it difficult to obtain v k The analytical expression for (x,R). Based on the VB approximation principle, the joint posterior strength is approximated as v. k (x,R)≈Q x (x k )Q R (R k ), where Q x (x k Q follows a Gaussian distribution. R (R k The distribution is an inverse gamma-Gaussian distribution, expressed as:

[0047]

[0048] Where IG(·;α,β) represents the inverse gamma distribution.

[0049] Furthermore, the optimal approximation of the posterior strength can be obtained by minimizing the KL (Kullback-Leibler) divergence between the true posterior value and the approximation value.

[0050] Step 3: Introduce the Cubature Kalman Filter (CKF) into the VB-PHD filter based on the Gaussian Mixture (GM) model, and use the third-order spherical radial volume rule to complete the integral operation of the PHD recursive density function under the nonlinear system model.

[0051] The non-Gaussian recursive implementation process of VB-PHD based on the spherical radial volume rule is as follows:

[0052] ·predict

[0053] At time k-1, the joint posterior v of the target state and the variance of the measurement noise is... k-1 (x,R) can be approximated by the product of a Gaussian mixture and an inverse gamma distribution, expressed as:

[0054]

[0055] Furthermore, the joint prediction intensity v k|k-1 (x,R) can be represented as:

[0056] v k|k-1 (x,R)=v S,k|k-1 (x,R)+γ k (x,R) (11)

[0057]

[0058] Where, ρ l ∈(0,1] is the forgetting factor, used to adjust variance, γ k (x,R) represents the simultaneous birth intensity, expressed by a Gaussian mixture and inverse gamma distribution.

[0059] For each mean, covariance is The Gaussian components, their predicted mean and covariance are calculated as follows:

[0060] i) Given Based on Cholesky decomposition, we can obtain

[0061] ii) Based on the dimension n of the target state vector x x The number of volume points can be determined as u = 2n x Then, the volume point χ p,k-1 p = 1, 2, ..., u can be generated based on the following formula:

[0062]

[0063] in,

[0064] iii) Propagating volume points through a single-objective nonlinear motion model:

[0065] χ p,k|k-1 =f k|k-1 (χ p,k-1 ),p=1,2,...,u (18)

[0066] iv) The predicted mean and covariance can be calculated as follows:

[0067]

[0068] ·renew

[0069] Assume the joint prediction strength of the multi-target state and measurement noise variance is v. k|k-1 Given (x,R) and a mixture of inverse gamma and Gaussian distributions, the joint posterior v k From (x,R), we can conclude that:

[0070]

[0071] To obtain the parameters of the inverse gamma and Gaussian mixture, the following calculations are performed. First, set... Then, iterate through the following steps:

[0072]

[0073] Given the predicted mean Covariance Matrix Under these conditions, the measurement value update calculation is as follows:

[0074] i) The predicted values ​​for u volume points are expressed as follows:

[0075]

[0076] ii) Calculate the propagation volume point based on the observation model:

[0077]

[0078] iii) The estimated value of the predicted measurement is:

[0079]

[0080] iv) The estimated value of the cross matrix is:

[0081]

[0082] v) Estimated value of the new information covariance matrix:

[0083]

[0084] Furthermore, the state update value is calculated as follows:

[0085]

[0086] pass Obtained by Cholesky decomposition, i.e. The inverse gamma-Gaussian distribution parameter β is updated as follows:

[0087]

[0088] until The iteration ends, and settings are set. and

[0089] Finally, the Gaussian component weights can be calculated as follows:

[0090]

[0091] Application Examples

[0092] Based on a space-based radar multi-target airborne tracking scenario, the above method is explained and analyzed. ADS-B (Automatic Dependent Surveillance-Broadcast) data of real aircraft targets and low-Earth orbit space-based radar data are imported into the satellite tool software STK (Satellite Tool Kit). STK is used to generate aircraft target trajectory data within the satellite observation area, and this is used as the true trajectory value. Range, azimuth, and pitch dimension measurement noise is added in the space-based radar observation coordinate system to generate target measurement points. Simultaneously, clutter measurements are added within the observation area and modeled as a Poisson random set K. k Its density is expressed as:

[0093]

[0094] in, λ represents the uniform density covering the entire observation area, V is the volume of the observation area, and λ is the density of the uniform density covering the entire observation area. c This represents the average number of clutter particles per unit volume. In this embodiment, V is set to 2 × 10⁻⁶. 10 , λ c =2.5×10 -9 In the ECEF coordinate system, the trajectories of the satellite and target are as follows: Figure 2 As shown, the simulation scenario after adding clutter measurements is as follows: Figure 3 As shown.

[0095] To comprehensively evaluate the tracking performance of the method, this embodiment designed simulation scenarios under different measurement noise variances. The number of Monte Carlo simulations was 500 for each parameter. The tracking accuracy was evaluated from both multi-target and single-target perspectives. Specifically, OSPA (Optimal Subpattern Assignment) distance was used to evaluate the multi-target tracking performance, and RMSE (Root Mean Square Error) was used to compare the root mean square error values ​​of a single target in terms of range, azimuth, pitch angle, and three-dimensional coordinate system. Figure 4 The figure shows a comparative analysis of tracking accuracy after introducing VB processing under different measurement noise parameters.

[0096] from Figure 4 As can be seen from (a) to (d), under the condition that the azimuth and pitch angle errors are constant, the tracking accuracy of both deteriorates as the distance error increases. However, compared with PHD-CKF, PHD-CKF-VB has a smaller OSPA distance and the improvement gradually increases. Therefore, it has stronger robustness to changes in distance noise.

[0097] Figure 5(a)-(c) show a comparative analysis of the root mean square error of target 1 in the scene for slant range, azimuth, and elevation. It can be seen that after 20s, the tracking filter gradually converges, but the improvement in filtering accuracy for azimuth and elevation is more significant after introducing variational processing.

[0098] also, Figure 5 (d) represents the three-dimensional filtering result in the ECEF coordinate system. It can be seen that as the target slant range increases, the measurement error of the XYZ axes gradually increases. At this time, the introduction of variational processing can effectively improve the filtering accuracy and enhance the robustness of airborne target tracking by space-based radar.

[0099] The above technical solutions are merely exemplary embodiments of the present invention. For those skilled in the art, based on the application methods and principles disclosed in the present invention, it is easy to make various types of improvements or modifications, and not limited to the methods described in the specific embodiments of the present invention. Therefore, the methods described above are merely preferred and not restrictive.

Claims

1. A space-based radar air target nonlinear robust tracking method, comprising the following steps: S1, constructing a space-based radar target observation geometry model, establishing the conversion relationship between the target state and the measurement; S2, introducing the variational Bayesian (VB) approximation method into the probability hypothesis density (PHD) recursive framework, and deriving the Gaussian closed solution of the PHD filter under the condition of unknown measurement noise covariance; S3, introducing the cubature Kalman filter into the VB-PHD filter based on the Gaussian mixture model, and completing the PHD recursive density function integral operation under the nonlinear system model by using the third-order spherical radial cubature rule; The step S1 specifically comprises: The earth is a regular uniform sphere, the radar carrier satellite platform is in a circular near-earth orbit, and the radar air target tracking system model is described by using a state equation and an observation equation, wherein the target state equation is established in a geocentric fixed rectangular coordinate system, and the observation equation is established in a polar coordinate system; The position of the satellite in the fixed right-angle coordinate system is set as (x s ,y s ,z s ), the target coordinate position of the airplane is (x a ,y a ,z a ), and the measured value of the airplane target in the radar load satellite body coordinate system is The conversion relationship from the measured value in the radar body coordinate system to the state value in the right-angle coordinate system is: The conversion relationship between the aircraft target and the radar body coordinate system is: The multi-target state at the k time is described by using the RFS as: In the formula, B k is the new target state set; and ∪ is the set union operation. The target motion state is represented as: wherein: x k,i , y k,i , z k,i are the spatial positions of the target in the ECEF coordinate system; is the target motion velocity; Assuming that the target observation process and the clutter generation process are independent, the multi-target measurement Z received by the sensor can be represented as k The union of the measurements generated by the targets and the clutter can be represented as In the formula, K k is the false alarm clutter measurement set; In the step S2, the VB-PHD filter recursive process is as follows: Prediction Assume that the joint posterior density of the multi-target states and the measurement noise variance at time k - 1 is v k-1 (x, R), and both are independent of each other, then the joint predictive intensity is: v k|k-1 (x, R) = ∫p S,k (ζ)f k|k-1 (x|ζ)f k|k-1 (R|Σ) x v k-1 (ζ, Σ) dζ dΣ + γ k (x, R) (6) where γ k is the simultaneous birth intensity; Update Assume that the joint predictive density of the multi-target states and the measurement noise variance at time k - 1 is v k|k-1 (x, R), then the joint posterior intensity is represented as: where p D,k (x) represents the detection probability, κ k (z) represents the Poisson clutter density; In the above equation, g k (z|x,R) denotes the target likelihood function, and R k is unknown and time-varying, it is difficult to obtain the analytical expression of v k (x,R). Therefore, according to the VB approximation principle, the joint posterior intensity is approximated as v k (x,R)≈Q x (x k )Q R (R k ), where Q x (x k ) is a Gaussian distribution, and Q R (R k ) is an inverse-Gaussian distribution, which is expressed as: Wherein, IG(·; a, b) represents the inverse gamma distribution; The recursive process in the step S3 comprises: Prediction At time k - 1, the joint prediction intensity v k|k-1 (x, R) is represented in a Gaussian mixture form, based on the third-order spherical radial volume rule, for the jth Gaussian component as with the prediction mean and covariance computed as ​ i) Given Based on Cholesky decomposition, we get ii) the dimension n of the target state vector x x The number of volume points is determined as Then, the volume points χ p,k-1 p = 1, 2,..., u are generated based on the following equation: wherein iii) Propagating the cubature points through the single-target nonlinear motion model: χ p,k|k-1 = f k|k-1 (χ p,k-1 ), p = 1, 2,..., u (18) iv) The predicted mean and covariance can be calculated as: Update Joint posterior v k (x, R) denotes the product of a Gaussian distribution and an inverse Gamma distribution, and to obtain the parameters of the inverse Gamma and Gaussian mixture, the following calculations are made: First, set Then, the following steps are iterated: Given the prediction mean and covariance matrix Under the condition, the measurement update is calculated as follows: i) The predicted value of the u cubature points is represented as: ii) The propagated cubature points are calculated based on the observation model: iii) The predicted measurement estimate is: iv) The cross-matrix estimate is: v) The innovation covariance matrix estimate is: Further, the state update value is calculated as: By Cholesky decomposition, i.e. The inverse gamma Gaussian distribution parameter b is updated as: Until The iteration ends and the setting And 2. The method of claim 1, wherein, In the step S2, the optimal approximation of the posterior intensity is obtained by minimizing the KL divergence between the true posterior value and the approximate value.

Images