Multi-target tracking method in unknown, non-uniform and time-varying clutter environment

Abstract

According to the multi-target tracking method in the unknown, non-uniform and time-varying clutter environment, the clutter density is calculated in real time through adaptive bandwidth kernel density estimation (KDE), and the performance of a probability hypothesis density (PHD) filter is optimized; comprising the steps of initializing filter parameters, predicting a target state, updating the target state, trimming low-weight targets, combining similar targets and the like. Compared with a traditional GM-PHD algorithm, the method has the advantages that the accuracy of clutter density estimation is remarkably improved, the false track rate is reduced, and the target cardinal number estimation and tracking precision is improved. Experimental results show that the method has better tracking performance and robustness in a complex clutter environment, and is suitable for radar, sonar and other sensor application scenes.

Description

Technical Field

[0001] The present invention belongs to the technical field of multi-target tracking, and particularly relates to a multi-target tracking method in an unknown, non-uniform, and time-varying clutter environment. Background Art

[0002] With the rapid development of sensor technologies such as radar, sonar, infrared, and optoelectronics, multi-target tracking (MTT) has been increasingly widely applied in fields such as military, traffic monitoring, unmanned driving, and security. Multi-target tracking aims to perform real-time positioning and state estimation of multiple targets in a dynamic environment. The number and state of targets often change with time, and together with the interference of clutter and noise in the electromagnetic environment, a complex tracking scenario is formed.

[0003] To address the problems of low computational efficiency and poor tracking performance of traditional multi-target tracking algorithms in complex scenarios, in 1994, Mahler proposed a theoretical system of finite set statistics (FISST) for processing multi-sensor multi-target filtering based on random finite set (RFS) theory. This theory provides a comprehensive, unified, and explicit statistical model that can integrate the two steps of target detection and state estimation in multi-target tracking into a Bayesian optimal step, thereby improving computational efficiency and tracking accuracy. Mahler proposed the probability hypothesis density (PHD) algorithm at the beginning of this century, and its core is to use the posterior intensity to approximate the uncomputable multi-target probability density. Vo B-N provides a closed-form implementation of the PHD algorithm in the linear Gaussian (GM) model, namely the Gaussian mixture (PHD Gaussian Mixture PHD, GM-PHD) filter.

[0004] GM-PHD relies on estimating the probability distribution of the target state and assumes that the clutter density is known and uniformly distributed. However, in practical applications, the clutter conditions are often unknown and may exhibit time-varying and non-uniform characteristics due to the complex electromagnetic environment. The limitations of this assumption make GM-PHD perform poorly in processing multi-target tracking in real scenarios, especially in the case of strong clutter, where the performance of the filter significantly degrades.

[0005] Some related literature shows that the difficulty in dealing with clutter environments stems from the complexity of clutter models: in practical applications, clutter can come from multiple sources, and its distribution characteristics are complex and difficult to accurately model. For example, clutter in the ocean environment changes due to factors such as waves and meteorological changes, making fixed clutter models inapplicable. There are also dynamically changing clutter characteristics. In some cases, the characteristics of clutter change over time, causing tracking algorithms based on static assumptions to fail. For example, as the target moves and the environment changes, the intensity and distribution of background clutter may change significantly.

[0006] In view of the above problems, aiming at the limitations of traditional multi-target tracking algorithms in non-uniform and time-varying clutter environments, a new solution is proposed. In this paper, the kernel density estimation method (KDE) with adaptive bandwidth is used to obtain clutter parameters in real time, and the filter is optimized based on these parameters. The simulation experimental results show that this method significantly improves the algorithm performance in dealing with unknown, non-uniform and time-varying clutter environments, demonstrating its application potential in the field of multi-target tracking. Summary of the Invention

[0007] The purpose of the embodiments of the present invention is to provide a multi-target tracking method in an unknown, non-uniform and time-varying clutter environment to improve the tracking performance of multi-target tracking algorithms in unknown, non-uniform and time-varying clutter environments.

[0008] The technical solution adopted by the present invention is a multi-target tracking method in an unknown, non-uniform and time-varying clutter environment, and the steps include:

[0009] Step S1, set the filter-related parameters according to the multi-target scenario, and initialize the state estimation, estimation covariance and weight of the initial frame of the filter;

[0010] Step S2, perform the filter prediction step to obtain the target prediction state, prediction covariance and prediction weight of the k-th frame;

[0011] Step S3, perform the filter update step to obtain the target update state, estimation covariance and weight of the k-th frame;

[0012] Step S4, set the pruning weight threshold;

[0013] Step S5, set the merging threshold;

[0014] Step S6, according to the weights of the updated targets after merging, output and store the target states, count the number of targets, enter the next frame, and repeat steps S1 to S5 until all frames are traversed to find all targets.

[0015] Further, the specific steps of S2 are as follows:

[0016] S21. Using the Kalman filter prediction formula, obtain the predicted state, predicted state estimation covariance, and weight of the surviving targets in the (k-1)-th frame. The specific formulas are as follows:

[0017]

[0018] P survival.g,k = Q + F·P g,k-1 ·F T

[0019] w survival.g,k = P s ·w g,k-1

[0020] where x survival.g,k represents the predicted state of the surviving target, P survival.g,k represents the predicted state estimation covariance of the surviving target, represents the predicted state of the g-th surviving target in the (k-1)-th frame, P g,k-1 represents the predicted state estimation covariance of the g-th surviving target in the (k-1)-th frame, Q represents the process noise covariance, F represents the state transition matrix, w survival.g,k represents the predicted state weight of the surviving target, P s represents the target survival probability, w g,k-1 represents the weight of the g-th surviving target in the (k-1)-th frame, T represents the transpose;

[0021] S22. Combine the surviving targets and the newly born targets as the predicted targets;

[0022] S23. Calculate the Mahalanobis distance between the measurement in the k-th frame and the target predicted state. According to the noise data of the target measurement and the threshold value d γ set by the association probability, compare the threshold value d γ with the Mahalanobis distance, filter out the measurements that are easily associated with the predicted target state, and delete the measurements with low association possibility. The calculation formula of the Mahalanobis distance d mahal.n,i,k is as follows:

[0023]

[0024] where, represents the predicted state of the i-th predicted target in the k-th frame, z m,k represents the m-th measurement in the k-th frame, H represents the measurement matrix, R represents the measurement noise, represents the predicted covariance of the i-th predicted target among them, T represents the transpose, represents to z m,k the square of the Mahalanobis distance. When d mahal.n,i,k > d γ , filter out the corresponding measurement zm,k When d mahal.n,i,k ≤ d γ and this value is the minimum value corresponding to the measured z m,k give the measured z m,k a virtual weight is the predicted weight. After filtering all the measurements, the next update step is carried out.

[0025] Furthermore, the set of predicted states, predicted weights, and predicted covariances of the predicted targets described in S22 is specifically as follows:

[0026]

[0027] Among them, represents the set of predicted states of the k-th frame target, represents the set of predicted weights of the k-th frame, represents the set of predicted covariances of the k-th frame, x birth.j,k represents the predicted state of the j-th newly born target, t represents the total number of newly born targets, x survival.g,k represents the predicted state of the g-th surviving target, v represents the total number of surviving targets, w birth.j,k represents the predicted weight of the j-th newly born target, w survival.g,k represents the predicted state of the g-th surviving target, P birth.j,k represents the predicted covariance of the j-th newly born target, P survival.g,k represents the predicted covariance of the g-th surviving target.

[0028] Furthermore, the updated targets in S3 include the targets that have been successfully detected and the targets that have not been detected. The specific steps of S3 are as follows:

[0029] S31, Obtain the states of the missed detection targets, the state estimation covariances, and the weights. The specific formulas are as follows:

[0030]

[0031] Among them, represents the state of the missed detection target, P qd.i,k represents the state estimation covariance of the missed detection target, w qd.i,k represents the weight of the missed detection target, P D represents the detection probability, represents the predicted state of the i-th predicted target in the k-th frame, represents the state estimation covariance of the i-th predicted target in the k-th frame, represents the predicted weight of the i-th predicted target in the k-th frame;

[0032] S32. Use the Kalman filter update formula to obtain the detection target state estimate, the detection target estimation covariance, and the detection target estimation likelihood;

[0033] S33. Calculate the detection target weight according to the clutter density;

[0034] S34. Combine the missed detection target set and the detection target set as the updated target set.

[0035] Furthermore, in the above S32, the specific calculation formulas for the detection target state estimate, the detection target estimation covariance, and the detection target estimation likelihood are as follows:

[0036]

[0037] Among them, represents the detection target state estimate, P pd.n,i,k the detection target estimation covariance, represents the detection target estimation likelihood, I is the unit diagonal matrix, H is the measurement matrix, R is the measurement noise, represents the target prediction covariance, S pd.i,k represents the prediction variance, K pd.i,k represents the Kalman gain, represents the target prediction state, z m,k represents the m-th measurement data in the k-th frame, T represents the transpose, and e represents the base of the natural logarithm.

[0038] Furthermore, in the above S33, the clutter density at the location of each target is calculated in real time through kernel density estimation with an adaptive bandwidth. The adaptive bandwidth calculation method is as follows:

[0039]

[0040] Among them, h represents the fixed bandwidth, h m represents the adaptive bandwidth, α represents the sensitivity factor, 0 ≤ α ≤ 1. When α takes 0, the adaptive bandwidth will become the fixed bandwidth; M represents the number of measurements, is the virtual weight corresponding to the measurement z m m represents the m-th measurement, f(z m ) represents the kernel density estimate value at the measurement z m at the fixed bandwidth;

[0041] The adaptive kernel density estimate after using the adaptive bandwidth is:

[0042]

[0043] The above result is used as the probability density of the clutter, that is:

[0044]

[0045] Substitute this value into the formula for calculating the detection and reporting target weight \(w\) qd.n,i,k to calculate the weight value of the detection and reporting target. The specific formula is as follows:

[0046]

[0047] where \(w\) pd.m,i,k represents the detection and reporting target weight, \(P\) D represents the detection probability, represents the estimated likelihood of the detection and reporting target, \(N\) is the number of predicted targets, represents the prediction weight, \(M\) represents the number of samples, represents the adaptive kernel density estimation function at the multi-dimensional coordinate \(z\), \(h\) m represents the adaptive bandwidth, \(n\) represents the measurement dimension, \(z\) m is the multi-dimensional coordinate of the \(m\)-th measurement, represents the probability density of clutter, \(H\) is the measurement matrix, represents the target prediction state, \(K\) h (·) is the multi-dimensional kernel function.

[0048] Furthermore, in S34, update the target update state set of the target set Update the weight set \(W\) k and the updated covariance set \(P\) k , and the specific sets are as follows:

[0049]

[0050]

[0051] where represents the update state set of the \(k\)-th frame, \(W\) k represents the update weight set of the \(k\)-th frame, \(P\) k represents the updated covariance set of the \(k\)-th frame, represents the update state of the \(i\)-th predicted target, \(N\) represents the number of predicted targets, represents the update state of the target obtained by associating the \(m\)-th measurement with the \(i\)-th predicted target, \(M\) represents the number of measurements, \(w\) qd.i,k represents the update weight of the \(i\)-th predicted target, \(w\) pd.m,i,k represents the update weight of the target obtained by associating the \(m\)-th measurement with the \(i\)-th predicted target, \(P\) qd.i,k represents the updated covariance of the \(i\)-th predicted target, \(P\) pd.m,i,k represents the updated covariance of the target obtained by associating the \(m\)-th measurement with the \(i\)-th predicted target.

[0052] Furthermore, in S4, set the pruning weight threshold \(\tau\) prune, when w i,k <τ prune , delete the corresponding target data.

[0053] Furthermore, in S5, a merging threshold τ is set. merge , the Mahalanobis distance between the a-th target and the b-th target When , the two targets are merged, and the weighted state mean and weighted covariance of the two targets after merging, as well as the sum of the weights, are calculated as the new merged components. The specific formula is as follows:

[0054] w c,k =w a,k +w b,k

[0055]

[0056] Among them, w a,k 、w b,k , P a,k , P b,k Represent the updated weight, updated target state, and updated covariance of the a-th target and the b-th target respectively; w c,k , P c,k They represent the merged updated weights, updated target states, and updated covariance respectively.

[0057] The beneficial effects of the present invention are:

[0058] 1. In the process of multi-target tracking, the present invention ensures the tracking effect in a complex clutter environment by estimating the clutter density in real time.

[0059] 2. The present invention uses an adaptive kernel density estimation method to estimate clutter density, taking into account factors such as clutter characteristics and target weight, as well as local and overall characteristics of clutter density, avoiding the estimation bias problem of fixed bandwidth kernel density estimation and improving the accuracy of clutter distribution estimation.

[0060] 3. The present invention is more accurate in target cardinality estimation, and the false track rate is lower than that of the traditional GM-PHD algorithm. When the clutter characteristics are complex, the present invention can also obtain the clutter density in time, especially in high clutter density areas, the false track rate is much lower than that of the traditional algorithm.

[0061] 4. In terms of tracking accuracy, the present invention avoids the influence of clutter as much as possible by accurately estimating the clutter distribution. Compared with the traditional GM-PHD algorithm, the tracking accuracy is higher and the distance error is significantly reduced. BRIEF DESCRIPTION OF THE DRAWINGS

[0062] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the following will briefly introduce the accompanying drawings required for the description of the embodiments or the prior art. Obviously, the accompanying drawings in the following description are only some embodiments of the present invention. For those of ordinary skill in the art, without creative efforts, other accompanying drawings can be obtained based on these drawings.

[0063] Figure 1 It is the trajectory diagram of the real target simulated in the simulation experiment.

[0064] Figure 2 It is the measurement diagram simulated in the simulation experiment.

[0065] Figure 3 It is the tracking result diagram of the traditional GM-PHD method.

[0066] Figure 4 It is the tracking result diagram of the method proposed in this paper.

[0067] Figure 5 It is the comparison diagram of the optimal sub-pattern assignment distance between the method proposed in this paper and the traditional GM-PHD method.

[0068] Figure 6 It is the comparison diagram of the cardinality estimation error between the method proposed in this paper and the traditional GM-PHD method.

[0069] Figure 7 It is the comparison diagram of the distance estimation error between the method proposed in this paper and the traditional GM-PHD method.

[0070] Figure 8 It is the target quantity estimation diagram between the method proposed in this paper and the traditional GM-PHD method.

[0071] Among them, Figures 3 - 4 in, the black trajectory is the tracking result, and the gray dots are the measurements. Specific implementation manners

[0072] Embodiment

[0073] As Figures 1 - 8 shown, the embodiment of the present invention provides a multi-target tracking method in an unknown, non-uniform, and time-varying clutter environment. The steps include:

[0074] Step S1, set the filter-related parameters according to the multi-target scenario, and initialize the state estimation of the initial frame of the filter Estimate the covariance P i,0 and the weight w i,0 .

[0075] Step S2, perform the filter prediction step to obtain the predicted state, predicted covariance, and prediction weight of the k-th frame. The predicted targets include the surviving targets from the (k-1)-th frame to the k-th frame and the newly born targets in the k-th frame. Then, filter the measurement data of the k-th frame using an association threshold to remove measurements with low association probabilities, reducing the computational complexity of the association calculation. The specific steps are as follows:

[0076] S21, using the Kalman filter prediction formula, obtain the predicted state, predicted state estimation covariance, and weight of the surviving targets in the (k-1)-th frame. The specific formulas are as follows:

[0077]

[0078] P survival.g,k = Q + F·P g,k-1 ·F T

[0079] w survival.g,k = P s ·w g,k-1

[0080] where, x survival.g,k represents the predicted state of the surviving target, P survival.g,k represents the predicted state estimation covariance of the surviving target, represents the predicted state of the g-th surviving target in the (k-1)-th frame, P g,k-1 represents the predicted state estimation covariance of the g-th surviving target in the (k-1)-th frame, Q represents the process noise covariance, F represents the state transition matrix, w survival.g,k represents the prediction state weight of the surviving target, P s represents the target survival probability, w g,k-1 represents the weight of the g-th surviving target in the (k-1)-th frame, T represents the transpose.

[0081] S22, combine the surviving targets and the newly born targets as the predicted targets. The sets of the predicted state, prediction weight, and prediction covariance of the predicted targets are as follows:

[0082]

[0083] where, represents the set of predicted states of the k-th frame targets, represents the set of prediction weights of the k-th frame, represents the set of prediction covariances of the k-th frame, x birth.j,k represents the predicted state of the j-th newly born target, t represents the total number of newly born targets, x survival.g,k represents the predicted state of the g-th surviving target, v represents the total number of surviving targets, w birth.j,k represents the prediction weight of the j-th newly born target, w survival.g,krepresents the predicted state of the g-th survival target, P birth.j,k represents the predicted covariance of the jth new target, P survival.g,k represents the predicted covariance of the g-th survival target.

[0084] S23, calculating the Mahalanobis distance between the measurement of the kth frame and the target prediction state. The threshold value d is set according to the noise data of the target measurement and the associated probability. γ , threshold value d γ By comparing with the Mahalanobis distance, we can filter out the measurements that are easily associated with the predicted target state and delete the measurements with low association probability, so as to reduce the computational burden in the data association process and improve the efficiency of the multi-target tracking algorithm. mahal.n,i,k The calculation formula is as follows:

[0085]

[0086] in, represents the predicted state of the i-th predicted target in the k-th frame, z m,k represents the mth measurement of the kth frame, H represents the measurement matrix, R represents the measurement noise, represents the predicted covariance of the i-th predicted target, T represents the transposition, express to z m,k The square of the Mahalanobis distance, when d mahal.n,i,k >d γ When the corresponding measurement z is filtered out m,k , when d mahal.n,i,k ≤d γ And this value is the measurement z m,k When the corresponding minimum value is given, the measurement z m,k Assign a virtual weight To predict the weights, all measurements are filtered out before proceeding to the next update step.

[0087] Step S3, perform filter update step, obtain the target update state, estimated covariance and weight of the kth frame, and the updated targets include successfully detected targets (i.e., reported targets) and undetected targets (i.e., missed targets). When calculating the reported target weight, adaptive kernel density estimation is used to estimate the clutter density to improve the calculation accuracy. The specific steps are as follows:

[0088] S31, obtain the missed target state, state estimation covariance and weight, the specific formula is as follows:

[0089]

[0090] in, Indicates the missed target state, Pqd.i,k Indicates the covariance of the missed detection target state estimate, w qd.i,k Indicates the weight of the missed detection target, P D Indicates the detection probability Indicates the predicted state of the i-th predicted target in the k-th frame Indicates the covariance of the state estimate of the i-th predicted target in the k-th frame Indicates the predicted weight of the i-th predicted target in the k-th frame

[0091] S32. Using the Kalman filter update formula, obtain the detected target state estimate, the detected target estimation covariance, and the detected target estimation likelihood. The specific calculation formulas are as follows:

[0092]

[0093]

[0094] Among them, Indicates the detected target state estimate, P qd.n,i,k The detected target estimation covariance Indicates the detected target estimation likelihood, I is the unit diagonal matrix, H is the measurement matrix, R is the measurement noise Indicates the target prediction covariance, S pd.i,k Indicates the prediction variance, K pd.i,k Indicates the Kalman gain Indicates the target predicted state, z m,k Indicates the m-th measurement data in the k-th frame, T represents the transpose, and e represents the base of the natural logarithm

[0095] S33. Calculate the detected target weight according to the clutter density. Different from the uniform clutter prior of the traditional algorithm, the present invention calculates the clutter density at the location of each target in real time through kernel density estimation, and uses Indicates in the measurement space The estimated density clutter at the location. Kernel Density Estimation is a non-parametric estimation method commonly used in statistics to infer its probability density function through a finite sample. On this basis, the kernel density estimation with an adaptive bandwidth further exploits the advantages of this method. It fully excavates the internal characteristics of the sample data without making any prior assumptions about the data distribution, thereby significantly improving the accuracy and robustness of the probability density function estimation. The specific calculation formula of the kernel density estimation is as follows:

[0096]

[0097] If the Gaussian kernel function is used for kernel density estimation and the normal distribution approximation (Gaussian approximation, Silverman's rule of thumb) is used, the optimal choice of the multi-dimensional bandwidth h (i.e., the bandwidth that minimizes the average integrated squared error) is:

[0098]

[0099] The Gaussian kernel function is:

[0100]

[0101] where K h (·) is the multi-dimensional kernel function, h represents the fixed bandwidth, M represents the number of measurements, represents the probability density at the multi-dimensional coordinate z, and z m is the multi-dimensional coordinate of the m-th measurement, represents the standard deviation of the measurement, n represents the dimension of the measurement, M represents the number of measurements, h represents the fixed bandwidth, x = z - z s , and x represents the independent variable of K h (·).

[0102] When the density distribution of this empirical formula is relatively complex, there may be serious estimation biases. The kernel density estimation method with an adaptive bandwidth is obtained by modifying the bandwidth parameter on the basis of the kernel density function with a fixed bandwidth, which can avoid this problem to a certain extent. The calculation method of the adaptive bandwidth is as follows:

[0103]

[0104] where h represents the fixed bandwidth, h m represents the adaptive bandwidth, α represents the sensitivity factor, 0 ≤ α ≤ 1, usually taking 0.5. When α takes 0, the adaptive bandwidth will become the fixed bandwidth; M represents the number of measurements, is the virtual weight corresponding to the measurement z m , m represents the m-th measurement, and f(z m ) represents the kernel density estimation value at the measurement z m under the fixed bandwidth.

[0105] When the virtual weight m corresponding to the measurement z is higher, the local clutter density is higher, and more attention needs to be paid to the specific clutter characteristics. Therefore, the adaptive bandwidth is narrower. On the contrary, when the virtual weight m corresponding to the measurement z is lower and the local clutter density is lower, a larger adaptive bandwidth is required to obtain the clutter characteristics of a larger area. The adaptive kernel density estimation after using the adaptive bandwidth is:

[0106]

[0107] This result can be approximately used as the probability density of clutter, that is:

[0108]

[0109] Substitute this value into the formula for calculating the detection and reporting target weight w qd.n,i,k in the formula, and the weight value of the detection and reporting target can be calculated. The specific formula is as follows:

[0110]

[0111] where w pd.m,i,k represents the detection and reporting target weight, P D represents the detection probability, represents the estimated likelihood of the detection and reporting target, N is the number of predicted targets, represents the prediction weight, M represents the number of samples, represents the adaptive kernel density estimation function at the multi-dimensional coordinate z, h m represents the adaptive bandwidth, n represents the measurement dimension, z m is the multi-dimensional coordinate of the m-th measurement, represents the probability density of clutter, H is the measurement matrix, represents the target prediction state, K h (·) is the multi-dimensional kernel function.

[0112] S34, merge the missed detection target set and the detection and reporting target set as the updated target set, and the target update state set of the updated target set Update the weight set W k and the updated covariance set P k , and the specific sets are as follows:

[0113]

[0114]

[0115] where, represents the updated state set of the k-th frame, W k represents the updated weight set of the k-th frame, P k represents the updated covariance set of the k-th frame, represents the updated state of the i-th predicted target, N represents the number of predicted targets, represents the updated state of the target obtained by associating the m-th measurement with the i-th predicted target, M represents the number of measurements, w qd.i,k represents the updated weight of the i-th predicted target, w pd.m,i,k represents the updated weight of the target obtained by associating the m-th measurement with the i-th predicted target, P qd.i,k represents the updated covariance of the i-th predicted target, Ppd.m,i,k Represents the updated covariance of the target obtained by associating the m-th measurement with the i-th prediction target.

[0116] Step S4, remove the estimated targets with smaller weights, remove the components that may be considered unimportant or noise and have less impact on the model, thereby reducing the computational amount in the subsequent process, and set the pruning weight threshold τ prune , when w i,k <τ prune , delete the corresponding target data.

[0117] Step S5, set the merging threshold τ merge , the Mahalanobis distance between the a-th target and the b-th target , merge these two targets, calculate the weighted state mean, weighted covariance, and the sum of weights after the merger of the two targets as the new merged component. The specific formulas are as follows:

[0118] w c,k = w a,k + w b,k

[0119]

[0120] Among them, w a,k , w b,k , P a,k , P b,k respectively represent the updated weights, updated target states, and updated covariances of the a-th target and the b-th target; w c,k , P c,k respectively represent the updated weights, updated target states, and updated covariances after the merger.

[0121] Step S6, output and store the target states with updated weights greater than 0.5 after the merger, count the number of targets, enter the next frame, and repeat Steps S1 to S5 for prediction, update, pruning, and merging until all frames are traversed to find all targets, thus completing the multi-target tracking process.

[0122] Through the above steps of updating, predicting, pruning, and merging of the probability hypothesis filter, using the adaptive kernel density algorithm to estimate the clutter density in real time, the accurate estimation of the target state is realized, and thus the accuracy of target tracking is improved.

[0123] Experimental verification

[0124] Conduct a simulation experiment on the true target trajectory, and the results are as Figure 1As shown, the target true trajectory is the black solid line, the starting point is marked with a black "○", and the ending point is marked with a black "△". In the observation space, there are clutters with the following clutter characteristics:

[0125] From frame 0 to 100, there are clutters whose quantity follows a Poisson distribution with an expectation of 100 and is uniformly distributed in the space where x ∈ [-1000, 1000] and y ∈ [-1000, 1000];

[0126] From frame 20 to 100, there are clutters whose quantity follows a Poisson distribution with an expectation of 50 and is uniformly distributed in the space where x ∈ [-500, 500] and y ∈ [-500, 500];

[0127] From frame 40 to 100, there are clutters whose quantity follows a Poisson distribution with an expectation of 25 and is uniformly distributed in the space where x ∈ [-250, 250] and y ∈ [-250, 250];

[0128] The obtained measurements are as Figure 2 shown. The z-axis is time, and the x-axis and y-axis are the x coordinate and y coordinate of the measurement respectively.

[0129] It can be seen from the figure and the clutter characteristics that the clutter in the entire observation space conforms to the characteristics of time-varying and non-uniform. The traditional GM-PHD algorithm and the KDE-GM-PHD algorithm are respectively used to track the target, and the tracking results are respectively as Figure 3 and Figure 4 shown. By comparing the tracking results of the two figures, it can be seen that the tracking error of the KDE-GM-PHD algorithm is smaller and the false tracks are fewer. Multiple Monte Carlo experiments are carried out to compare the specific tracking performances of the two algorithms.

[0130] Figure 5 It is the comparison of the optimal sub-pattern assignment (OSPA) distances of the two algorithms, Figure 6 It is the comparison of the cardinality estimation errors of the two algorithms, Figure 7 It is the comparison of the distance estimation errors of the two algorithms. It can be seen from the above three figures that the distance error of the KDE-GM-PHD algorithm is significantly lower than that of the traditional GM-PHD algorithm, Figure 8 It is the comparison of the cardinality estimation contrasts of the two algorithms. It can be proved that the cardinality estimation value of the KDE-GM-PHD algorithm has a smaller gap with the true value compared with the traditional GM-PHD algorithm. And the peak value of the cardinality estimation error of the KDE-GM-PHD algorithm is higher than that of the traditional GM-PHD algorithm. This is because in the high clutter density area, the target starting requires more frames to exclude false targets, and the prior of the PHD algorithm is uniform clutter, and there will be many false targets in the high clutter density area as Figure 3As shown, the accuracy of cardinality estimation is significantly better than that of the GM-PHD algorithm in the case of no new target start, which reflects that the KDE-GM-PHD algorithm improves the performance of multi-target tracking.

[0131] Each embodiment in this specification is described in a related manner. For the same or similar parts among the embodiments, reference can be made to each other. Each embodiment focuses on the differences from other embodiments. In particular, for the system embodiment, since it is basically similar to the method embodiment, the description is relatively simple, and the relevant parts can be referred to the partial description of the method embodiment.

[0132] The above description is only a preferred embodiment of the present invention and is not intended to limit the protection scope of the present invention. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention are included in the protection scope of the present invention.

Claims

1. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment, characterized in that the steps Including: Step S1: Set the filter-related parameters according to the multi-objective scenario, and initialize the state estimation, estimation covariance, and weight of the initial frame of the filter; Step S2: Perform the filter prediction step to obtain the target prediction state, prediction covariance, and prediction weight of the k-th frame; Step S3: Perform the filter update step to obtain the target update state, estimation covariance, and weight of the k-th frame; Step S4: Set the pruning weight threshold; Step S5: Set the merging threshold; Step S6: According to the weight of the updated target after merging, output and store the target state, count the number of targets, enter the next frame, and repeat steps S1 to S5 until all frames are traversed to find all targets.

2. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 1, characterized in that, The specific steps of S2 are as follows: S21: Use the Kalman filter prediction formula to obtain the prediction state, prediction state estimation covariance, and weight of the surviving targets in the (k - 1)-th frame. The specific formula is as follows: P survival.g,k = Q + F·P g,k-1 ·F T w survival.g,k = P s · w g,k-1 where x survival.g,k represents the predicted state of the survival target, P survival.g,k represents the covariance estimate of the predicted state of the survival target, represents the predicted state of the g-th survival target in the (k-1)-th frame, P g,k-1 represents the covariance estimate of the predicted state of the g-th survival target in the (k-1)-th frame, Q represents the process noise covariance, F represents the state transition matrix, w survival.g,k represents the weight of the predicted state of the survival target, P s represents the target survival probability, w g,k-1 represents the weight of the g-th survival target in the (k-1)-th frame, T represents the transpose; S22: Merge the surviving targets and the newly born targets as the prediction targets; S23. Calculate the Mahalanobis distance between the measurement of the k-th frame and the target predicted state, and set the threshold value d according to the noise data and association probability of the target measurement γ , the threshold value d γ is compared with the Mahalanobis distance to screen out the measurements that are easily associated with the predicted target state, and delete the measurements with low association possibility. The Mahalanobis distance d mahal.n,i,k is calculated as follows: Among them, represents the prediction state of the i-th predicted target in the k-th frame, z m,k represents the m-th measurement in the k-th frame, H represents the measurement matrix, R represents the measurement noise, represents the prediction covariance of the i-th predicted target among them, T represents the transpose, represents to z m,k the square of the Mahalanobis distance. When d mahal.n,i,k > d γ , filter out the corresponding measurement z m,k . When d mahal.n,i,k ≤ d γ and this value is the minimum value corresponding to the measurement z m,k , assign a virtual weight to the measurement z m,k . is the prediction weight. After filtering out all measurements, proceed to the next update step.

3. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 2, characterized in that The sets of the prediction state, prediction weight, and prediction covariance of the prediction targets described in S22 are specifically as follows: Among them, represents the set of target prediction states for the k-th frame, represents the set of prediction weights for the k-th frame, represents the set of prediction covariances for the k-th frame, x birth.j,k represents the predicted state of the j-th newly born target, t represents the total number of newly born targets, x survival.g,k represents the predicted state of the g-th surviving target, v represents the total number of surviving targets, w birth.j,k represents the prediction weight of the j-th newly born target, w survival.g,k represents the predicted state of the g-th surviving target, P birth.j,k represents the predicted covariance of the j-th newly born target, P survival.g,k represents the predicted covariance of the g-th surviving target.

4. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 1, characterized in that, The updated targets in S3 include the targets that have been successfully detected and the targets that have not been detected. The specific steps of S3 are as follows: S31: Obtain the missed detection target state, state estimation covariance, and weight. The specific formula is as follows: Among them, represents the missed detection target state, P qd.i,k represents the covariance of the missed detection target state estimation, w qd.i,k represents the weight of the missed detection target, P D represents the detection probability, represents the predicted state of the i-th predicted target in the k-th frame, represents the covariance of the state estimation of the i-th predicted target in the k-th frame, represents the predicted weight of the i-th predicted target in the k-th frame; S32: Use the Kalman filter update formula to obtain the detected target state estimation, detected target estimation covariance, and detected target estimation likelihood; S33: Calculate the detected target weight according to the clutter density; S34: Merge the missed detection target set and the detected target set as the updated target set.

5. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 4, characterized in that In S32, the detected target state estimation, detected target estimation covariance, and detected target estimation likelihood are specifically calculated as follows: Among them, represents the detection and reporting target state estimation, P pd.n,i,k is the covariance of the detection and reporting target estimation, represents the likelihood of the detection and reporting target estimation, I is the unit diagonal matrix, H is the measurement matrix, and R is the measurement noise, represents the covariance of the target prediction, S pd.i,k represents the prediction variance, K pd.i,k represents the Kalman gain, represents the predicted state of the target, z m,k represents the m-th measurement data in the k-th frame, T represents the transpose, and e represents the base of the natural logarithm.

6. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 4, characterized in that, In S33, the clutter density at the location of each target is calculated in real time through kernel density estimation with an adaptive bandwidth. The adaptive bandwidth calculation method is as follows: where h represents the fixed bandwidth, and h m represents the adaptive bandwidth, α represents the sensitivity factor, 0 ≤ α ≤ 1. When α is 0, the adaptive bandwidth becomes the fixed bandwidth; M represents the number of measurements, is the virtual weight corresponding to the measurement z m , m represents the m-th measurement, and f(z m ) represents the kernel density estimate value of the measurement z m at the fixed bandwidth. The adaptive kernel density estimation after using the adaptive bandwidth is: The above results are used as the probability density of the clutter, that is: Substitute this value into the formula for calculating the weight \(w\) of the detection target qd.n,i,k to calculate the weight value of the detection target. The specific formula is as follows: Among them, w pd.m,i,k represents the detection target weight, P D represents the detection probability, represents the detection target estimated likelihood, N is the number of predicted targets, represents the prediction weight, M represents the number of samples, represents the adaptive kernel density estimation function at the multi-dimensional coordinate z, h m represents the adaptive bandwidth, n represents the measurement dimension, z m is the multi-dimensional coordinate of the m-th measurement, represents the probability density of clutter, H is the measurement matrix, represents the target prediction state, K h (·) is the multi-dimensional kernel function.

7. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 4, characterized in that, In S34, update the target update status set of the target set Update the weight set W k And update the covariance set P k , and the specific sets are as follows: Among them, represents the update status set of the k-th frame, W k represents the update weight set of the k-th frame, P k represents the update covariance set of the k-th frame, represents the update status of the i-th predicted target, N represents the number of predicted targets, represents the update status of the target obtained by associating the m-th measurement with the i-th predicted target, M represents the number of measurements, w qd.i,k represents the update weight of the i-th predicted target, w pd.m,i,k represents the update weight of the target obtained by associating the m-th measurement with the i-th predicted target, P qd.i,k represents the update covariance of the i-th predicted target, P pd.m,i,k represents the update covariance of the target obtained by associating the m-th measurement with the i-th predicted target.

8. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 1, characterized in that In S4, a pruning weight threshold τ is set. prune , when w i,k < τ prune , the corresponding target data is deleted.

9. A multi-target tracking method in an unknown, non-uniform, time-varying clutter environment according to claim 1, characterized in that, The S5 sets a merging threshold τ merge , the Mahalanobis distance between the a-th target and the b-th target When, merge these two targets, and calculate the weighted state mean, weighted covariance, and the sum of weights after the merger of the two targets as the new merged component. The specific formulas are as follows: w c,k = w a,k + w b,k Among them, w a,k , w b,k , P a,k , P b,k respectively represent the update weights, update target states, and update covariances of the a-th target and the b-th target; w c,k , P c,k respectively represent the combined update weights, update target states, and update covariances.

Images