Non-linear robust tracking method for space-based radar aerial target

Abstract

A space-based radar aerial target nonlinear robust tracking method mainly comprises the following steps: constructing a space-based radar target observation geometric model, and establishing a conversion relation between a target state and measurement; a VB approximation method is introduced into a PHD recursive framework, and a Gaussian closed solution of PHD filtering under the condition of unknown measurement noise covariance is deduced, so that the problem that the measurement noise is uncertain under the space-based observation condition is solved; on the basis, volume Kalman filtering is introduced into a VB-PHD filter based on a Gaussian mixture model, PHD recursive density function integral operation under a nonlinear system model is completed by using a third-order spherical radial volume rule, and the tracking precision under a strong nonlinear observation condition is improved. According to the method, the filtering precision can be effectively improved, and the robustness of space-based radar aerial target tracking is improved.

Description

Technical Field

[0001] The invention relates to the technical field of radar target tracking, and in particular to a nonlinear robust tracking method for aerial targets of a space-based radar. Background Art

[0002] Space-based radar uses satellites as payload platforms, which are not restricted by land, territorial waters and climate. It can detect and continuously and stably track low-altitude aerial targets from a long distance, which can effectively improve early warning capabilities. However, due to the influence of factors such as the platform's high-speed movement, long-distance detection, perturbation interference, and nonlinear observation, space-based radar target tracking poses severe challenges, mainly manifested in the following aspects: 1) The space-based radar payload platform and the aerial target are in different inertial systems, resulting in nonlinear radar observation equations; 2) Under long-distance observation conditions, the measurement noise uncertainty increases, making high-precision stable tracking difficult.

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

[0004] Probability Hypothesis Density (PHD) is a first-order moment approximation of the optimal multi-objective Bayesian filter. It has been widely used due to its low computational complexity. The article "The Gaussian mixture probability hypothesis density filter" derived the analytical solution of PHD recursion under the linear Gaussian dynamic model, and extended it to the nonlinear dynamic model using the extended Kalman filter (Extended Kalman Filter, EKF) and unscented Kalman filter (Unscented Kalman Filter, UKF), but it is only applicable to weak nonlinear and known measurement noise scenarios. For actual space-based radar target tracking scenarios, the observation equation is strongly nonlinear, and the measurement noise is unknown and time-varying.

[0005] The variational Bayesian method (VB) can be used to estimate the target state under unknown measurement noise variance. 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 paper "An improved PHD filter based on variationalBayesian method for multi-target tracking" introduces the variational Bayesian approximation method into the PHD recursive filter to achieve the joint estimation of the posterior distribution of multi-target states and 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 is not applicable to nonlinear observation models. Summary of the invention

[0006] In order to solve the problem of strong nonlinear observation and unknown measurement noise variance in space-based radar aerial target tracking, the present invention provides a space-based radar aerial target nonlinear robust tracking method based on VB-PHD filter:

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

[0008] Then, the VB approximation method is introduced into the PHD recursive framework to derive the Gaussian closed-form solution of PHD filtering with unknown measurement noise covariance, so as to solve the problem of measurement noise uncertainty under space-based observation conditions.

[0009] On this basis, the cubic Kalman filter (CKF) is introduced into the VB-PHD filter based on the Gaussian mixture model (GM), and the third-order spherical radial volume rule is used to complete the PHD recursive density function integral operation under the nonlinear system model to improve the tracking accuracy under strong nonlinear observation conditions.

[0010] Compared with the prior art, the beneficial effects of the present invention are: (1) a geometric model of aerial target observation of space-based radar is established, and the VB approximation method is introduced into the PHD filter framework, which solves the problem of measurement noise uncertainty under space-based observation conditions; (2) the volumetric Kalman filter algorithm is introduced into the VB-GM-PHD filter framework, and the density function integral operation under the nonlinear system model in the recursive process of the PHD filter is completed by using the third-order spherical radial volume rule, thereby improving the filtering accuracy under strong nonlinear observation conditions; (3) the introduction of variational processing can effectively improve the filtering accuracy and enhance the robustness of space-based radar aerial target tracking. BRIEF DESCRIPTION OF THE DRAWINGS

[0011] The above and other objects, features and advantages of the present disclosure will become more apparent through a more detailed description of exemplary embodiments of the present disclosure in conjunction with the accompanying drawings, wherein like reference numerals generally represent like components throughout the exemplary embodiments of the present disclosure.

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

[0013] Figure 2 A schematic diagram of a space-based radar target tracking scenario in an ECEF coordinate system in an exemplary embodiment;

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

[0015] Figure 4 is a comparison diagram of multi-target OSPA distances in an exemplary embodiment;

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

[0017] Figure 6 is a flow chart of an exemplary embodiment according to the present disclosure. DETAILED DESCRIPTION

[0018] The preferred embodiments of the present disclosure will be described in more detail below with reference to the accompanying drawings. Although the preferred embodiments of the present disclosure are shown in the accompanying drawings, it should be understood that the present disclosure can be implemented in various forms and should not be limited by the embodiments set forth herein. On the contrary, these embodiments are provided to make the present disclosure more thorough and complete, and to fully convey the scope of the present disclosure to those skilled in the art.

[0019] The present disclosure proposes a nonlinear robust tracking method for space-based radar aerial targets based on VB-PHD filter. The exemplary embodiment process is shown in the attached figure. Figure 6 As shown, the specific steps include:

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

[0021] Assume that the earth is a regular uniform sphere, and the radar carrier satellite platform is in a circular low-Earth orbit. The radar air target tracking system model is usually described by state equations and observation equations. The target state equation is a description of the target motion state, which is established in the Earth-Centered and Earth-Fixed (ECEF) rectangular coordinate system. The coordinate origin O is located at the center of the earth, 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 is a description of the real-time measurement value of the radar on the target, which is usually established in the polar coordinate system. It is assumed that the coordinate origin O' of the space-based radar body coordinate system O'X'Y'Z' is the center position of the radar, the O'X' axis is the extension of the line connecting the center of the earth O and the center of the radar O', the O'Z' axis is in the plane formed by the O'X' axis and the earth's rotation axis and is perpendicular to the O'X' axis, and the O'Y' is determined by the O'Z' axis and the O'X' axis according to the right-hand screw law, such as Figure 1 shown.

[0022] Assume that the position of the satellite in the geocentric rectangular coordinate system is (x s ,y s ,z s ), the aircraft target coordinate position is (x a ,y a ,z a ), in the radar payload satellite body coordinate system, the aircraft target measurement value is Then the conversion relationship from the radar body coordinate system measurement value to the rectangular coordinate system state value is:

[0023]

[0024] Then the conversion relationship between the aircraft target from the geocentric rectangular coordinate system to the radar body coordinate system is:

[0025]

[0026] Based on formulas (1) and (2), it can be concluded that the state conversion relationship is strongly nonlinear. Therefore, the present 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 and derive the Gaussian closed-form solution of PHD filtering under the condition of unknown measurement noise covariance 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:

[0029]

[0030] Where: B kis the new target state set; ∪ is the set union operation.

[0031] The target motion state is expressed as:

[0032]

[0033] Where: x k,i ,y k,i 、z k,i is the spatial position of the target in the ECEF coordinate system; is the target movement speed.

[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 is k It can be represented by the union of the measurements and clutter generated by the target:

[0035]

[0036] In the formula, K k is the false alarm clutter measurement set.

[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 multi-target state and measurement noise variance is v k-1 (x, R), and the two are independent of each other, then the joint prediction strength can be obtained as:

[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] Among them, γ k It is the intensity of simultaneous birth.

[0042] ·renew

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

[0044]

[0045] Among them, 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 target likelihood function, and R k Unknown and time-varying, difficult to get v k According to 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 ) is Gaussian distribution, Q R (R k ) is the inverse gamma Gaussian distribution, expressed as:

[0047]

[0048] Here, 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.

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

[0051] The VB-PHD non-Gaussian recursive implementation process 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 measurement noise variance 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 strength v k|k-1 (x,R) can be expressed as:

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

[0057]

[0058] Among them, ρ l ∈(0,1] is the forgetting factor, which is used to adjust the variance, γ k (x,R) are the simultaneous birth intensities represented by a mixture of Gaussian and inverse gamma distributions.

[0059] For each mean The covariance is The Gaussian component of , whose predicted mean and covariance are calculated as follows:

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

[0061] ii) According to the dimension n of the target state vector x x , the volume point number 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 volumetric points via a single target nonlinear motion model:

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

[0066] iv) The mean and covariance of the predictions can be calculated as:

[0067]

[0068] ·renew

[0069] Assume that the joint prediction strength of multi-target states and measurement noise variance is v k|k-1 (x, R), and represented by an inverse gamma and Gaussian mixture distribution, then the joint posterior v k (x,R) can be obtained:

[0070]

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

[0072]

[0073] Given the predicted mean and the covariance matrix Under these conditions, the measurement value update calculation is as follows:

[0074] i) The predicted value of u volume points is expressed as:

[0075]

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

[0077]

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

[0079]

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

[0081]

[0082] v) Estimation of the innovation covariance matrix:

[0083]

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

[0085]

[0086] pass The Cholesky decomposition of The inverse gamma Gaussian distribution parameter β is updated as:

[0087]

[0088] until The iteration ends and sets and

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

[0090]

[0091] Application Examples

[0092] The above method is explained and analyzed based on the scenario of multi-target tracking in the sky by space-based radar. The ADS-B (Automatic Dependent Surveillance-Broadcast) data of real aircraft targets and low-orbit space-based radar data are imported into the satellite tool software STK (Satellite Tool Kit). STK is used to generate the trajectory data of aircraft targets in the satellite observation area, and this is used as the true value of the trajectory. In the space-based radar observation coordinate system, distance, azimuth, and pitch dimension measurement noise are added to generate target measurement points. At the same time, clutter measurements are added in 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 λ c represents the average number of clutter per unit volume. In this embodiment, V is set to 2×10 10 ,λ c =2.5×10 -9 In the ECEF coordinate system, the trajectory of the satellite and the target is as follows: Figure 2 As shown in Figure 2, the simulation scene after adding clutter measurement is as follows: Figure 3 shown.

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

[0096] from Figure 4 From (a) to (d), it can be seen that under the condition of certain azimuth and elevation errors, as the distance error increases, the tracking accuracy of both deteriorates. However, compared with PHD-CKF, the OSPA distance of PHD-CKF-VB is smaller and the improvement gradually increases. Therefore, it is more robust to changes in distance noise.

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

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

[0099] The above technical scheme is only an exemplary embodiment of the present invention. For those skilled in the art, it is easy to make various types of improvements or modifications based on the application methods and principles disclosed in the present invention, and it is not limited to the method described in the above specific embodiment of the present invention. Therefore, the method described above is only preferred and does not have a restrictive meaning.

Claims

1. A nonlinear robust tracking method for aerial targets using a space-based radar, comprising the following steps: S1, construct the geometric model of space-based radar target observation and establish the conversion relationship between target state and measurement; S2, the VB approximation method is introduced into the PHD recursive framework to derive the Gaussian closed-form solution of PHD filtering with unknown measurement noise covariance; S3, introduces the volumetric Kalman filter into the VB-PHD filter based on the Gaussian mixture model, and uses the third-order spherical radial volume rule to complete the PHD recursive density function integral operation under the nonlinear system model.

2. The method according to claim 1, characterized in that The step S1 specifically includes: Assume that the earth is a regular uniform sphere, the radar carrier satellite platform is in a circular low-Earth orbit, and the radar air target tracking system model is described by the state equation and observation equation. The target state equation is established in the geocentric rectangular coordinate system, and the observation equation is established in the polar coordinate system. Assume that the position of the satellite in the geocentric rectangular coordinate system is (x s ,y s ,z s ), the aircraft target coordinate position is (x a ,y a ,z a ), in the radar payload satellite body coordinate system, the aircraft target measurement value is Then the conversion relationship from the radar body coordinate system measurement value to the rectangular coordinate system state value is: Then the conversion relationship between the aircraft target from the geocentric rectangular coordinate system to the radar body coordinate system is: The multi-objective state at time k is described by RFS: Where: B k is the new target state set; ∪ is the set union operation; The target motion state is expressed as: Where: x k,i ,y k,i 、z k,i is the spatial position of the target in the ECEF coordinate system; is the target movement rate; 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 is k The union of the measurements and clutter generated by the target is expressed as: In the formula, K k is the false alarm clutter measurement set.

3. The method according to claim 2, characterized in that In step S2, the VB-PHD filtering recursive process is as follows: ·predict Assume that at time k-1, the joint posterior density of the multi-target state and measurement noise variance is v k-1 (x, R), and the two are independent of each other, then the joint prediction strength is: v k|k-1 (x,R)=∫p S,k (g)f k|k-1 (x|ζ)f k|k-1 (R|Σ)×v k-1 (ζ,Σ)dζdΣ+γ k (x,R) (6) Among them, γ k For simultaneous birth intensity; ·renew Assume that at time k-1, the joint prediction density of multi-target states and measurement noise variance is v k|k-1 (x,R), then the joint posterior strength is expressed as: Among them, p D,k (x) represents the detection probability, κ k (z) represents the Poisson clutter density; In the above formula, due to g k (z|x,R) represents the target likelihood function, and R k Unknown and time-varying, difficult to get v k (x, R), so according to the VB approximation principle, the joint posterior strength is approximately expressed as v k (x,R)≈Q x (x k )Q R (R k ), where Q x (x k ) is Gaussian distribution, Q R (R k ) is the inverse gamma Gaussian distribution, expressed as: Here, IG(·; α, β) represents the inverse gamma distribution.

4. The method according to claim 3, characterized in that In step S2, the optimal approximation of the posterior strength is obtained by minimizing the KL divergence between the true posterior value and the approximation.

5. The method according to claim 3 or 4, characterized in that: The recursive process in step S3 includes: ·predict At time k-1, the joint prediction strength v k|k-1 (x, R) is expressed as a Gaussian mixture, based on the third-order spherical radial volume rule, for each mean The covariance is The Gaussian component of , whose predicted mean and covariance are calculated as follows: i) Given Based on Cholesky decomposition, we get ii) According to the dimension n of the target state vector x x , the volume point number is determined as u = 2n x ; Then, the volume point χ p,k-1 , p=1,2,...,u is generated based on the following formula: in, iii) Propagating volumetric points via a single target nonlinear motion model: x p,k|k-1 =f k|k-1 (x p,k-1 ),p=1,2,...,u (18) iv) The mean and covariance of the predictions can be calculated as: ·renew Joint posterior v k (x, R) represents the product of Gaussian distribution and inverse gamma distribution. To obtain the parameters of inverse gamma and Gaussian mixture, the following calculation is performed: First, set Then, iterate the following steps: Given the predicted mean and the covariance matrix Under these conditions, the measurement value update calculation is as follows: i) The predicted value of u volume points is expressed as: ii) Calculate the propagation volume points based on the observation model: iii) The estimated value of the forecast measurement is: iv) The estimated value of the cross matrix is: v) Estimation of the innovation covariance matrix: Furthermore, the state update value is calculated as follows: pass The Cholesky decomposition of The inverse gamma Gaussian distribution parameter β is updated as: until The iteration ends and sets and

Images