Track-before-detect method for radar weak target

Abstract

The invention relates to the technical field of radar target tracking, in particular to a track-before-detect method for a radar weak target. In order to estimate the target state and the measurement noise covariance at the same time, firstly, the joint probability density of the target state and the measurement noise covariance needs to be calculated. And then Gaussian inverse gamma mixed distribution is introduced to model the joint probability density, and due to the fact that the target state and the measurement noise covariance are coupled in a joint likelihood function, the joint posterior probability density is difficult to analyze and solve, so that the separable approximate solution of the joint posterior probability density is approximately solved by adopting a variational Bayesian method. And finally, in a filtering updating stage, information exchange is carried out based on the separation approximate solution, namely, each Bernoulli component is updated by using prediction state information shared by other Bernoulli components. Through simulation verification, in a low signal-to-noise ratio scene, the method can adaptively estimate and measure the noise covariance, and the tracking precision is improved.

Description

Technical Field

[0001] The invention relates to the technical field of radar target tracking, and in particular to a pre-detection tracking method for radar weak targets. Background Art

[0002] In radar systems, the signal-to-noise ratio (SNR) of targets is low due to distance attenuation effects, environmental noise, system performance limitations, and other issues, making it difficult to estimate the motion state and trajectory of targets based on observed data. As a traditional weak target tracking algorithm, the detect-before-track (DBT) algorithm mainly sets a detection threshold to segment the measurement data, and then extracts target information from the segmented data. Compared with the DBT algorithm, track-before-detect (TBD) does not require setting a detection threshold, so the measurement data can be used to the maximum extent for target tracking. The current mainstream TBD algorithm is based on the random finite set (RFS) theory, which can effectively avoid complex data association and simultaneously jointly estimate the number and state of multiple targets.

[0003] Radar systems operate under different environmental conditions, such as urban, forest or marine environments. The noise characteristics of these environments may change with time and location, making it difficult to accurately estimate the measurement noise covariance, and weak target signals are easily masked by background noise, which in turn leads to a decrease in algorithm tracking accuracy. Existing TBD algorithms based on RFS theory all assume that the measurement noise covariance is known, such as the IEMB-PLF-TBD algorithm. Because it uses a fixed measurement noise covariance and uses a common target pruning method, its tracking accuracy will be reduced. Summary of the invention

[0004] The object of the present invention is to provide a method for tracking a weak target before detection by radar, aiming to solve the technical problem of poor accuracy of target tracking methods based on unknown measurement noise covariance in the prior art.

[0005] To achieve the above object, the present invention provides a method for tracking a weak target before detection by radar, comprising the following steps:

[0006] Step 1: Model the target set as a random finite set and use the point spread function to generate measurements to obtain the target state model and measurement model;

[0007] Step 2: In the process of calculating the joint probability density of the target state and the measurement noise covariance, the Gaussian inverse gamma mixture distribution is introduced to model the joint probability density, and the variational Bayesian method is used to approximate the separable approximate solution of the joint posterior probability density;

[0008] Step 3: Based on the target state model and the measurement model, a prediction step is performed using the Gaussian inverse gamma mixture implementation of the variational Bayesian information exchange multi-Bernoulli posterior linearization filter detection-before-tracking algorithm to obtain the target prediction state at time k;

[0009] Step 4: Perform the update step of the Gaussian inverse gamma mixture implementation of the variational Bayesian information exchange multi-Bernoulli posterior linearization filter detection before tracking algorithm, and perform information exchange during the update process, that is, each Bernoulli component uses the predicted state information shared by other Bernoulli components to perform the update, and obtain the target update state at time k;

[0010] Step 5: Prune and merge the targets, and extract the states based on the target existence probability.

[0011] Optionally, in step 1, assume that at time k (k = 1, ..., K) there are M k The target moves in the observation area and is represented by a random finite set: where F(X) is all finite subsets of the target states.

[0012] kth time l(l=1,…,M k The discrete motion equation of the target can be expressed as in is the state vector of target l, and Represent the position, speed and strength of target l respectively. k Represents the state transition matrix of the target. Represents process noise, which has a mean vector of 0 and a covariance matrix of Q k Gaussian distribution.

[0013] The infrared sensor provides a two-dimensional image sequence of the observation area. Each image contains N×M resolution units, and each resolution unit corresponds to a Δ x ×Δ y The observation area corresponding to the (i,j)th resolution unit is (iΔ x ×jΔ y ), i = 1, ..., N, j = 1, ..., M. The measured image is recorded at intervals of T, and the observed intensity of the (i, j)th resolution unit at the kth time It can be expressed as Where C represents the target impact diffusion area. represents the background noise of the (i, j)th resolution unit, which has a mean vector of 0 and a variance of (σ (i,j) ) 2 Gaussian distribution, and are independent between frames and resolution units. It represents the contribution intensity of the lth target at time k to the (i, j)th resolution unit, generally in the form of a point spread function Where ∑ usually represents the known fuzzy coefficient, and the measurement at time k can be expressed as For the sake of subsequent discussion, z is represented as an NM×1 column vector. k ,Right now The corresponding

[0014] Optionally, in step 2, since the inverse gamma distribution is the conjugate prior distribution of the variance of the Gaussian distribution, the inverse gamma distribution is often used to model the variance of the Gaussian distribution. Therefore, in the process of calculating the joint probability density of the target state and the measurement noise covariance, the Gaussian inverse gamma mixture distribution is introduced to model it. Since the target state and the measurement noise covariance are coupled in the joint likelihood function, this will make the joint posterior probability density difficult to solve analytically. Therefore, the variational Bayes method is used to approximate the separable approximate solution of the joint posterior probability density.

[0015] Optionally, in step 3, assume that the joint posterior probability density at time k-1 is expressed using a multi-Bernoulli parameter set as and denote the existence probability and spatial probability density of the i-th Bernoulli component respectively. The predicted joint probability density is expressed as in and The multi-Bernoulli parameter sets representing the survival target and the new target, M k-1 and M Γ,k They represent the number of Bernoulli components of the surviving target and the new target respectively, and the total number of predicted Bernoulli components is M k|k-1 =M k-1 +M Γ,k .

[0016] Optionally, in step 4, assume that the predicted joint probability density at time k is expressed using a multi-Bernoulli parameter set as When the measurement z at time k is given k , then the updated joint posterior probability density function is approximately expressed by the multi-Bernoulli parameter set as During the update process, each Bernoulli component uses the forecast state information shared by other Bernoulli components to calculate the forecast measurement value and the innovation covariance matrix Then according to the measurement z k Perform the update.

[0017] Optionally, in step 5, set the pruning threshold to T, and record is the pruned target set; set the merging threshold to U, record is the merged target set. Set the target state extraction threshold to ε, and extract the target state estimation with a probability greater than ε.

[0018] The present invention provides a pre-detection tracking method for radar weak targets. First, in order to simultaneously estimate the target state and the measurement noise covariance, the joint probability density of the target state and the measurement noise covariance is calculated. Then, the Gaussian inverse gamma mixture distribution is introduced to model the joint probability density, and the variational Bayesian method is used to approximate the separable approximate solution of the joint posterior probability density. Finally, in the filter update stage, information exchange is performed based on the separated approximate solution, that is, each Bernoulli component uses the predicted state information shared by other Bernoulli components to perform updates. It has been verified by simulation that in low signal-to-noise ratio scenarios, the present invention can adaptively estimate the measurement noise covariance, and the tracking accuracy is improved. BRIEF DESCRIPTION OF THE DRAWINGS

[0019] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0020] Figure 1 The present invention is a schematic flow chart of the steps of a method for tracking a weak radar target before detection.

[0021] Figure 2 It is a schematic diagram of the actual target trajectory of a specific embodiment of the present invention.

[0022] Figure 3 It is a comparison diagram of target number estimation of PHD-TBD, IEMB-IPLF-TBD, IEMB-UKF-TBD and the method proposed in the present invention (VB-IEMB-PLF-TBD) under different measurement noise standard deviations in this specific embodiment when the SNR is 3dB.

[0023] Figure 4 It is a comparison diagram of the overall OSPA errors of PHD-TBD, IEMB-IPLF-TBD, IEMB-UKF-TBD and the method proposed in the present invention (VB-IEMB-PLF-TBD) under different measurement noise standard deviations in this specific embodiment when the SNR is 3 dB.

[0024] Figure 5It is a comparison diagram of target number estimation of PHD-TBD, IEMB-IPLF-TBD, IEMB-UKF-TBD and the method proposed in the present invention (VB-IEMB-PLF-TBD) under different measurement noise standard deviations in this specific embodiment when the SNR is 5dB.

[0025] Figure 6 It is a comparison diagram of the overall OSPA errors of PHD-TBD, IEMB-IPLF-TBD, IEMB-UKF-TBD and the method proposed in the present invention (VB-IEMB-PLF-TBD) under different measurement noise standard deviations in this specific embodiment when the SNR is 5 dB. DETAILED DESCRIPTION

[0026] Embodiments of the present invention are described in detail below, examples of which are shown in the accompanying drawings, wherein the same or similar reference numerals throughout represent the same or similar elements or elements having the same or similar functions. The embodiments described below with reference to the accompanying drawings are exemplary and are intended to be used to explain the present invention, and should not be construed as limiting the present invention.

[0027] In the present application, there are the following terms: "SNR" means signal-to-noise ratio, "DBT" means detection before tracking, "TBD" means tracking before detection, "VB" means variational Bayes, "VB-IEMB-PLF-TBD" means information exchange multi-Bernoulli posterior linearization filtering detection and tracking algorithm based on variational Bayes, "PHD-TBD" means probability hypothesis density-based detection and tracking algorithm, "IEMB-UKF-TBD" means information exchange multi-Bernoulli unscented Kalman filtering detection and tracking algorithm, "IEMB-IPLF-TBD" means information exchange multi-Bernoulli iterative posterior linearization filtering detection and tracking algorithm, and "OSPA" means optimal sub-pattern allocation.

[0028] The present invention provides a method for tracking a weak target before detection by radar, comprising the following steps:

[0029] S1: Model the target set as a random finite set and use the point spread function to generate measurements to obtain the target state model and measurement model;

[0030] S2: In the process of calculating the joint probability density of the target state and the measurement noise covariance, the Gaussian inverse gamma mixture distribution is introduced to model the joint probability density, and the variational Bayesian method is used to approximate the separable approximate solution of the joint posterior probability density;

[0031] S3: Based on the target state model and the measurement model, the prediction step of the Gaussian inverse gamma mixture implementation of the variational Bayesian information exchange multi-Bernoulli posterior linearization filter detection and tracking algorithm is performed to obtain the target prediction state at time k;

[0032] S4: Perform the update step of the Gaussian inverse gamma mixture implementation of the variational Bayesian information exchange multi-Bernoulli posterior linearization filter detection tracking algorithm, and perform information exchange during the update process, that is, each Bernoulli component uses the predicted state information shared by other Bernoulli components to perform the update, and obtain the target update state at time k;

[0033] S5: Prune and merge the targets, and extract the states based on the target existence probability.

[0034] Please refer to the execution process and steps Figure 1 The following is a further explanation of the specific steps:

[0035] In step 1, assume that at time k (k = 1, ..., K) there are M k The target moves in the observation area and is represented by a random finite set: where F(X) is all finite subsets of the target states.

[0036] kth time l(l=1,…,M k The discrete motion equation of the target can be expressed as in is the state vector of target l, and Represent the position, speed and strength of target l respectively. k Represents the state transition matrix of the target. Represents process noise, which has a mean vector of 0 and a covariance matrix of Q k Gaussian distribution.

[0037] The infrared sensor provides a two-dimensional image sequence of the observation area. Each image contains N×M resolution units, and each resolution unit corresponds to a Δ x ×Δ y The observation area corresponding to the (i,j)th resolution unit is (iΔ x ×jΔ y ), i = 1, ..., N, j = 1, ..., M. The measured image is recorded at intervals of T, and the observed intensity of the (i, j)th resolution unit at the kth time It can be expressed as Where C represents the target impact diffusion area. represents the background noise of the (i, j)th resolution unit, which has a mean vector of 0 and a variance of (σ (i,j) )2 Gaussian distribution, and are independent between frames and resolution units. It represents the contribution intensity of the lth target at time k to the (i, j)th resolution unit, generally in the form of a point spread function Where ∑ usually represents the known fuzzy coefficient, and the measurement at time k can be expressed as For the sake of subsequent discussion, z is represented as an NM×1 column vector. k ,Right now The corresponding

[0038] In step 2, since the inverse gamma distribution is the conjugate prior distribution of the variance of the Gaussian distribution, the inverse gamma distribution is often used to model the variance of the Gaussian distribution. Therefore, in the process of calculating the joint probability density of the target state and the measurement noise covariance, the Gaussian inverse gamma mixture distribution is introduced to model it. Since the target state and the measurement noise covariance are coupled in the joint likelihood function, this will make the joint posterior probability density difficult to solve analytically. Therefore, the variational Bayes method is used to approximate the separable approximate solution of the joint posterior probability density.

[0039] In step 3, assume that the joint posterior probability density at time k-1 is expressed by a multi-Bernoulli parameter set: and denote the existence probability and spatial probability density of the i-th Bernoulli component respectively. The predicted joint probability density is expressed as in and The multi-Bernoulli parameter sets representing the survival target and the new target, M k-1 and M Γ,k They represent the number of Bernoulli components of the surviving target and the new target respectively, and the total number of predicted Bernoulli components is M k|k-1 =M k-1 +M Γ,k .

[0040] In step 4, assume that the predicted joint probability density at time k is expressed by a multi-Bernoulli parameter set as When the measurement z at time k is given k , then the updated joint posterior probability density function is approximately expressed by the multi-Bernoulli parameter set as During the update process, each Bernoulli component uses the forecast state information shared by other Bernoulli components to calculate the forecast measurement value and the innovation covariance matrix Then according to the measurement z k Perform the update.

[0041] In step 5, set the pruning threshold to T, and record is the pruned target set; set the merging threshold to U, record is the merged target set. Set the target state extraction threshold to ε, and extract the target state estimation with a probability greater than ε.

[0042] Furthermore, the present invention also provides a specific embodiment for auxiliary explanation. Through simulation experiments, VB-IEMB-PLF-TBD (the method of the present invention), IEMB-IPLF-TBD, IEMB-UKF-TBD and PHD-TBD are compared. In order to verify the performance of different algorithms, a multi-target scenario with different measurement noise standard deviations and SNRs of 3dB and 5dB is set, and for VB-IEMB-PLF-TBD, its R=diag((σ 1 ) 2 ,…(σ N×M ) 2 ) is unknown, and the measurement noise standard deviation of the comparison algorithm is set to σ 1 =…=σ N×M =σ=0.8,1,1.2, and analyzed through 100 Monte Carlo simulation experiments. For multi-target scenarios, the target number estimation and Optimal SubPattern Assignment (OSPA) distance are used to evaluate the algorithm, where the parameters of OSPA distance are set as: order p=1, truncation error c=10.

[0043] 1. Simulation conditions: The present invention is simulated on a computer with an Intel(R) Core(TM) i5-14600K@3.50 processor and 32.0GB memory.

[0044] 2. Simulation scene setting: multi-target motion trajectory such as Figure 2 As shown in the target discrete motion model and where q 1 =0.001 represents the target motion noise, q 2 = 0.010 represents the target signal intensity noise.

[0045] The sensor receives 50 frames of images at intervals of T = 1s. The actual measurement noise standard deviation σ = σ 1 =…=σ N×M =1, other parameter values ​​are Δ x =Δ y =1, n=m=20, ∑=0.9. SNR is expressed as Where I is set according to the simulation scenario.

[0046] Assume that there are 4 weak targets in total, and the target parameters are shown in Table 1, where the target signal strength I is set according to the simulation scenario.

[0047] Table 1 Target parameters

[0048]

[0049] The distribution of new targets uses multi-Bernoulli parameter sets Indicates that r Γ =10 -4 , about It uses the observed intensity M Γ = 20 resolution units as potential positions, speed is 0, and intensity is set according to the simulation scenario; the covariance matrix is ​​P Γ =diag([0.5,0.1,0.5,0.1,0.1]) 2 ; In order to ensure the stability of the adaptive estimation of the measurement noise covariance, the parameter of the initial inverse gamma distribution is set to α 0 =β 0 =100, degradation factor ρ = 0.99; target survival probability is set to P s =0.95.

[0050] About the setting of target pruning and merging parameters: pruning threshold T = 10 -3 , merging threshold U=6.

[0051] Experiment 1 considers a multi-target scenario with an SNR of 3dB. Figure 3-4 The target number estimation of different algorithms and the overall error comparison of OSPA are shown respectively under different measurement noise standard deviations σ. Figure 3 It can be seen that the target number estimation result of VB-IEMB-PLF-TBD remains stable and is generally better than the comparison algorithm. When σ = 0.8, because the value of the existence probability is inversely proportional to the measurement noise covariance, σ is smaller than the true value at this time, resulting in a larger existence probability, so the comparison algorithm will overestimate the number of targets. On the contrary, when σ = 1.2, because the measurement noise covariance is large, the existence probability is small, so the comparison algorithm underestimates the number of targets. Figure 4 It can also be seen that at the five moments k = 1, 10, 20, 30, and 40, the error of VB-IEMB-PLF-TBD will suddenly increase compared with other moments. This is because the target state estimation extraction sets an extraction threshold. The target is only considered to exist when the existence probability is greater than the extraction threshold. Therefore, errors will inevitably occur in tracking new targets and removing dead targets, which will lead to larger errors at these five moments.

[0052] Table 2 shows the comparison of the overall OSPA error averaged at all times and the algorithm time under different measurement noise standard deviations. As can be seen from Table 2, because VB-IEMB-PLF-TBD adds an iterative estimation process for parameter update during operation, its running time is longer than that of IEMB-UKF-TBD. It can also be seen that the error of the comparison algorithm fluctuates under different measurement noise standard deviations, while the error of the VB-IEMB-PLF-TBD algorithm is relatively small, because it adaptively estimates the measurement noise covariance and updates the mean vector and covariance matrix of the target state during the iteration process.

[0053] Table 2 Algorithm performance at 3dB SNR

[0054]

[0055]

[0056] Experiment 2 considers a multi-target scenario with an SNR of 5dB. Figure 5-6 The target number estimation of different algorithms and the overall error comparison of OSPA are shown respectively under different measurement noise standard deviations σ. Figure 5 It can be seen that the target number estimation result of VB-IEMB-PLF-TBD remains stable, while the comparison algorithm will over-estimate and underestimate the target number when σ=0.8 and σ=1.2 respectively. Compared with the 3dB scenario, the target number estimation of VB-IEMB-PLF-TBD is smoother and more stable. Figure 6 It can also be seen that in the case of target birth and death, the error of VB-IEMB-PLF-TBD at the five moments k = 1, 10, 20, 30, and 40 is reduced compared with the 3dB scenario, and the error convergence speed becomes faster.

[0057] Table 3 shows the comparison of the overall OSPA error and algorithm time averaged at all times under different measurement noise standard deviations. It can be seen from Table 3 that after the SNR is improved, the operating efficiency of each algorithm is improved, and the overall OSPA error shows a downward trend compared with the 3dB scenario; due to the stronger robustness of VB-IEMB-PLF-TBD, it has better tracking accuracy than the comparison algorithm.

[0058] Table 3 Algorithm performance at SNR 5dB

[0059]

[0060] In summary, compared with the prior art, the present invention has the following advantages:

[0061] 1. Compared with the traditional detection-then-tracking method, the measurement information can be used to the maximum extent for target tracking.

[0062] 2. When the measurement noise covariance is unknown, the measurement noise covariance can be estimated adaptively, which has better robustness.

[0063] 3. Fully consider the mutual impact between goals and effectively improve the utilization of information.

[0064] What is disclosed above is only one or more preferred embodiments of the present invention, which certainly cannot be used to limit the scope of rights of the present invention. Ordinary technicians in this field can understand that all or part of the processes of implementing the above embodiments and making equivalent changes according to the claims of the present invention still fall within the scope of the invention.

Claims

1. A method for tracking a weak target before detection by radar, characterized in that: The following steps are involved: Step 1: Model the target set as a random finite set and use the point spread function to generate measurements to obtain the target state model and measurement model; Step 2: In the process of calculating the joint probability density of the target state and the measurement noise covariance, the Gaussian inverse gamma mixture distribution is introduced to model the joint probability density, and the variational Bayesian method is used to approximate the separable approximate solution of the joint posterior probability density; Step 3: Based on the target state model and the measurement model, a prediction step is performed using the Gaussian inverse gamma mixture implementation of the variational Bayesian information exchange multi-Bernoulli posterior linearization filter detection-before-tracking algorithm to obtain the target prediction state at time k; Step 4: Perform the update step of the Gaussian inverse gamma mixture implementation of the variational Bayesian information exchange multi-Bernoulli posterior linearization filter detection before tracking algorithm, and perform information exchange during the update process, that is, each Bernoulli component uses the predicted state information shared by other Bernoulli components to perform the update, and obtain the target update state at time k; Step 5: Prune and merge the targets, and extract the states based on the target existence probability.

2. The method for tracking a weak radar target before detection as claimed in claim 1, characterized in that: In step 1, assume that at time k (k = 1, ..., K) there are M k The target moves in the observation area and is represented by a random finite set: Where F(X) is all finite subsets of the target state; kth time l(l=1,…,M k The discrete motion equation of the target can be expressed as in is the state vector of target l, and represent the position, velocity and strength of the target l, respectively, and F k represents the state transition matrix of the target, Represents process noise, which has a mean vector of 0 and a covariance matrix of Q k Gaussian distribution of The infrared sensor provides a two-dimensional image sequence of the observation area. Each image contains N×M resolution units, and each resolution unit corresponds to a Δ x ×Δ y The observation area corresponding to the (i,j)th resolution unit is (iΔ x ×jΔ y ), i = 1, ..., N, j = 1, ..., M; the measured image is recorded at intervals of T, and the observed intensity of the (i, j)th resolution unit at the kth time It can be expressed as Where C represents the target impact diffusion area, represents the background noise of the (i, j)th resolution unit, which has a mean vector of 0 and a variance of (σ (i,j) ) 2 Gaussian distribution, and is independent between frames and resolution units. It represents the contribution intensity of the lth target at time k to the (i, j)th resolution unit in the form of a point spread function Where ∑ usually represents the known fuzzy coefficient, and the measurement at time k can be expressed as Represent z in the form of an NM×1 column vector k ,Right now The corresponding 3. The method for tracking a weak radar target before detection as claimed in claim 2, characterized in that: In the process of calculating the joint probability density of the target state and the measurement noise covariance in step 2, a Gaussian inverse gamma mixture distribution is introduced to model it, and a variational Bayesian method is used to approximately solve the separable approximate solution of the joint posterior probability density.

4. The method for tracking a weak radar target before detection as claimed in claim 3, characterized in that: In step 3, assume that the joint posterior probability density at time k-1 is expressed by a multi-Bernoulli parameter set as and They represent the existence probability and spatial probability density of the i-th Bernoulli component respectively; the predicted joint probability density is expressed by the multi-Bernoulli parameter set as in and The multi-Bernoulli parameter sets representing the survival target and the new target, M k-1 and M Γ,k They represent the number of Bernoulli components of the surviving target and the new target respectively, and the total number of predicted Bernoulli components is M k|k-1 =M k-1 +M Γ,k .

5. The method for tracking a weak radar target before detection as claimed in claim 4, characterized in that: In step 4, assume that the predicted joint probability density at time k is expressed by a multi-Bernoulli parameter set as When the measurement z at time k is given k , then the updated joint posterior probability density function is approximately expressed by the multi-Bernoulli parameter set as During the update process, each Bernoulli component uses the forecast state information shared by other Bernoulli components to calculate the forecast measurement value and the innovation covariance matrix Then according to the measurement z k Perform the update.

6. The method for tracking a weak radar target before detection as claimed in claim 5, characterized in that: In step 5, set the pruning threshold to T, record is the pruned target set; set the merging threshold to U, record is the merged target set; set the target state extraction threshold to ε, and extract the target state estimation with a probability greater than ε.

Images