A multi-target track identity maintenance method and system in a strong clutter background

Abstract

The application relates to the technical field of radar data processing, and discloses a multi-target track identity maintaining method and system under a strong clutter background, wherein the multi-target track identity maintaining method under the strong clutter background comprises the following steps: acquiring current point track and historical track data to form a candidate point track set; calculating a point track quality score; extracting a motion mode feature vector of a track to construct a similarity matrix between tracks; calculating the weight of an edge in a bipartite graph and solving global optimal correlation to obtain a correlation matrix of the track and the point track; starting a track identity correction mechanism when a switching condition is triggered to obtain a corrected correlation matrix; through a global correlation optimization method, global optimal correlation of multi-target tracks is realized, the accuracy of track correlation in a multi-target scene is improved; through a track identity verification and correction mechanism, a track switching problem in multi-target cross maneuvering is solved, and the reliability of multi-target track identity maintaining is improved.

Description

Technical Field

[0001] This application relates to the field of radar data processing technology, and more specifically, to a method and system for maintaining the identity of multiple targets under strong clutter backgrounds. Background Technology

[0002] The simultaneous flight and overlapping maneuvers of multiple drones are becoming increasingly common, making continuous monitoring and accurate identification of such multi-target scenarios a key technology for low-altitude safety management. Radar, as a primary sensor for all-weather, long-range detection, plays a crucial role in drone monitoring. However, the strong clutter in urban environments and the high maneuverability of small drones pose challenges to radar target tracking.

[0003] Chinese patent CN113311398B discloses a method for tracking highly maneuverable, small targets against a complex background with strong clutter. The radar signal processing unit uses three parallel processing methods to detect the target signal from the received echo signal. After angle calculation to form a target point track, it is sent to the radar data processing unit. The track processing software performs secondary processing on the received target point track to form a target track. The target track enters the track fitting module via a sliding window, where a one-dimensional second-order fitting is performed on the target track information in a Cartesian coordinate system. The target fitting module calculates the target maneuvering acceleration variance and sends it to the track filtering module. The module continuously calculates the target motion parameters in real time and estimates the target motion state parameters for the next moment, completing the iterative update of the filter gain coefficients. This adaptively achieves stable tracking of highly maneuverable targets against a complex background with strong clutter.

[0004] However, existing technologies do not take into account the increased uncertainty in prediction when targets intersect. When the distance between targets is reduced to 2-3 times the radar resolution unit (about 30-50 meters), the association gates of multiple tracks overlap, resulting in one-to-many or many-to-one association ambiguity. After the targets intersect, the track of target A may be incorrectly continued as the track of target B, leading to confusion in identification and making it impossible to conduct effective threat assessment and evidence tracing. Summary of the Invention

[0005] The purpose of this application is to address the problem that existing radar multi-target tracking methods are prone to track misassociation and target identity continuation errors due to overlapping association gates and similar target movements in strong clutter backgrounds and multi-target cross-maneuver scenarios, making it difficult to stably maintain the track identity of multiple targets. Therefore, this application provides a method and system for maintaining the track identity of multiple targets in strong clutter backgrounds.

[0006] This application provides a method for maintaining the identity of multiple targets' tracks against a background of strong clutter, including the following steps:

[0007] Acquire the current UAV target's point traces and historical flight path data, perform flight path prediction based on the historical flight path data to obtain the predicted flight path position, calculate the distance between the point traces and the predicted flight path position, and form a candidate point trace set for points whose distance is less than the Mahalanobis distance threshold.

[0008] The quality of the traces in the candidate trace set is evaluated, the quality feature vector of the traces is extracted, and the trace quality score is calculated based on the quality feature vector.

[0009] Motion pattern modeling is performed for each track, and motion pattern feature vectors of the track are extracted based on historical position sequences in historical track data. A similarity matrix between tracks is constructed based on the motion pattern feature vectors.

[0010] Global association optimization is performed based on the candidate point set and point quality score. A bipartite graph including the track set and point set is constructed. The weights of the edges in the bipartite graph are calculated and the global optimal association is solved to obtain the association matrix of track and point.

[0011] The correlation matrix is ​​verified based on the similarity matrix between the tracks. When the exchange condition is triggered, the track identity correction mechanism is activated to obtain the corrected correlation matrix.

[0012] Furthermore, methods for forming a set of candidate traces include:

[0013] Obtain distance, azimuth, and elevation angle information from the point data;

[0014] The three-dimensional position coordinates of the point are determined based on the distance, azimuth, and pitch angle information, and the Mahalanobis distance between the three-dimensional position coordinates of the point and the predicted position of the flight path is calculated.

[0015] Points whose Mahalanobis distance is less than the Mahalanobis distance threshold are included in the candidate point set of the corresponding track.

[0016] Furthermore, methods for calculating dot quality scores include:

[0017] Extract the quality feature vector, which includes signal-to-noise ratio, Doppler consistency, radar cross-section stability, and clutter probability;

[0018] Each component of the quality feature vector is fuzzified. An S-shaped membership function is used for the signal-to-noise ratio. The original values ​​of Doppler consistency and radar cross-section stability are directly used as membership. The reverse membership function is used for clutter probability.

[0019] The quality score of the dots is calculated by the inner product of the preset basic weight vector and the fuzzy quality feature vector.

[0020] Furthermore, methods for extracting motion pattern feature vectors from tracks include:

[0021] The average velocity is calculated based on the three-dimensional velocity components in the historical location sequence and the preset number of historical moments.

[0022] Based on the three-dimensional velocity components, the velocity vector is extracted, and the difference between the velocity vectors at adjacent time points is divided by the sampling period to obtain the acceleration magnitude. Within the adjacent time points corresponding to the preset number of historical time points, the maximum value of the acceleration magnitude is the maximum acceleration.

[0023] The mean value of the sine of the angle between the velocity vectors at each adjacent time point is calculated to obtain the track curvature.

[0024] The straight-line flight index is calculated based on the angle between the velocity vectors of adjacent moments, and the straight-line flight probability is obtained by calculating the mean of the straight-line flight index at each moment within a preset historical time period.

[0025] The motion pattern feature vector includes average velocity, maximum acceleration, trajectory curvature, and probability of straight flight.

[0026] Furthermore, methods for constructing the similarity matrix between tracks include:

[0027] For any two tracks, calculate the Euclidean distance between the motion pattern feature vectors of the two tracks, and calculate the motion pattern similarity by multiplying the Euclidean distance by a natural exponential function with a negative preset attenuation coefficient.

[0028] The similarity matrix between the tracks is constructed using the similarity of the motion patterns as elements.

[0029] Furthermore, methods for calculating the weights of edges in a bipartite graph include:

[0030] Construct a bipartite graph, which includes a vertex set and an edge set. The vertex set includes a track set and a point set, and the edge set includes multiple edges. If point i is in the candidate point set of track j, then there is an edge between track j and point i.

[0031] The normalized distance similarity is calculated based on the distance between the predicted location of the point and the track, and the Mahalanobis distance threshold.

[0032] The consistency of trajectory prediction is calculated based on the three-dimensional position coordinates of the points, the predicted trajectory position, and the preset position deviation threshold.

[0033] The weight of the edge is a weighted sum of normalized distance similarity, point quality score, and trajectory prediction consistency.

[0034] Furthermore, methods for obtaining the correlation matrix between tracks and points include:

[0035] When a point is in the candidate point set of multiple paths at the same time, a competition penalty factor is calculated based on a preset competition penalty coefficient, and the weight of the original edge is multiplied by the competition penalty factor to obtain the adjusted edge weight.

[0036] Set up an association matrix. The association variable in the j-th row and i-th column of the association matrix represents the association relationship between track j and point i. The association relationship corresponds to the edge in the bipartite graph.

[0037] The association weights are calculated based on the adjusted weights and related variables;

[0038] With the goal of maximizing the total association weight, and under the constraint that each point can be associated with at most one track and each track can be associated with at most one point, the association matrix is ​​obtained by using the Hungarian algorithm.

[0039] Furthermore, the methods for determining the triggering of the swap condition include:

[0040] Calculate the velocity change and heading change angle based on the three-dimensional position coordinates of the point trace;

[0041] For each association pair in the association matrix, if the speed change is greater than the product of the preset maximum maneuvering acceleration capability and the sampling period, or the heading change angle is greater than the preset maximum heading change angle and the straight flight probability of the trajectory is greater than the preset straight flight probability threshold, then the exchange condition is triggered; wherein, the association pair is the trajectory j and point i that have an association relationship in the association matrix;

[0042] Calculate the exchange risk index based on the Mahalanobis distance between the point and the track;

[0043] For any two associated tracks in the current association matrix, if the exchange risk index is less than the preset exchange risk threshold, and the similarity of the motion patterns of the two associated tracks is greater than the preset motion pattern similarity threshold, then the exchange condition is triggered.

[0044] Furthermore, methods for obtaining the corrected correlation matrix include:

[0045] For track j and track l that trigger the exchange condition, and their associated points i and p in the current association matrix, determine the original association relationship and the association relationship to be exchanged; the original association relationship is that track j is associated with point i and track l is associated with point p, and the association relationship to be exchanged is that track j is associated with point p and track l is associated with point i.

[0046] Calculate the comprehensive score of the original association and the association to be exchanged. If the comprehensive score of the association to be exchanged is greater than the comprehensive score of the original association, then update the original association to the association to be exchanged in the current association matrix to obtain the corrected association matrix. The comprehensive score is the sum of the edge weights of the two tracks and their respective associated point tracks.

[0047] This application provides a multi-target track identification system under strong clutter background, which is used to implement the aforementioned multi-target track identification method under strong clutter background. The system includes:

[0048] The dynamic threshold module acquires the current UAV target's traces and historical trajectory data, performs trajectory prediction based on the historical trajectory data to obtain the predicted trajectory position, calculates the distance between the traces and the predicted trajectory position, and forms a candidate trace set for traces whose distance is less than the Mahalanobis distance threshold.

[0049] The quality assessment module assesses the quality of the traces in the candidate trace set, extracts the quality feature vectors of the traces, and calculates the trace quality score based on the quality feature vectors.

[0050] The motion modeling module models the motion pattern of each track, extracts the motion pattern feature vector of the track based on the historical position sequence in the historical track data, and constructs a similarity matrix between tracks based on the motion pattern feature vector;

[0051] The association optimization module performs global association optimization based on the candidate point set and point quality score, constructs a bipartite graph including the track set and the point set, calculates the weight of the edges in the bipartite graph and solves the global optimal association, and obtains the association matrix of track and point.

[0052] The verification and correction module verifies the association matrix based on the similarity matrix between the tracks. When the exchange condition is triggered, the track identity correction mechanism is activated to obtain the corrected association matrix.

[0053] The beneficial effects of this application are as follows: This application achieves adaptive tracking of maneuvering targets through a dynamic correlation threshold calculation mechanism. This adaptive mechanism effectively solves the problem that traditional fixed threshold methods are prone to losing tracks or generating false correlations when the target is maneuvering, significantly improving the stability and accuracy of maneuvering target tracking. Furthermore, by using a fuzzy comprehensive evaluation method to quantitatively assess the quality of the tracking points, the assessment of the tracking point quality is made more consistent with the actual motion state of the target, effectively improving the ability to distinguish between real target tracking points and clutter tracking points in a strong clutter background, and reducing the interference of clutter on track correlation.

[0054] This application achieves globally optimal association of multi-objective tracks through a graph-based global association optimization method. The method constructs a bipartite graph comprising a track set and a point set, calculates edge weights that comprehensively consider Mahalanobis distance, point quality, and track prediction consistency, and employs the Hungarian algorithm to globally match and solve the edge weights after competition penalty adjustment. Specifically, when a point is simultaneously in the candidate point sets of multiple tracks, a competition penalty factor is calculated to reduce the edge weights between all competing tracks and that point, thus favoring non-competitive association schemes in the global optimization. This global optimization strategy effectively avoids association conflicts and suboptimal solutions that are prone to occur in local association methods, significantly improving the accuracy of track association in multi-objective scenarios.

[0055] By employing a track authentication and correction mechanism, the problem of track exchange during multi-target cross-maneuvering is effectively resolved. The mechanism verifies the rationality of the association matrix, checks whether the changes in speed and heading of the associated tracks conform to the target's maneuverability and motion pattern feature vectors, and detects track exchange risks by calculating an exchange risk index combined with motion pattern similarity. When a track exchange risk is detected, the mechanism compares the combined scores of the original association scheme and the exchange association scheme, selecting the higher-scoring scheme to automatically correct the track identity. This authentication and correction mechanism ensures the consistency between the track identifier and the target's true identity, effectively avoiding track exchange and merging errors, and significantly improving the reliability of maintaining the track identity of multiple targets. Attached Figure Description

[0056] Figure 1 This is a flowchart illustrating a method for maintaining the identity of multiple targets under strong clutter backgrounds, as proposed in this application.

[0057] Figure 2 This is an example diagram of the correlation matrix generation for a multi-target track identity preservation method under strong clutter background proposed in this application;

[0058] Figure 3 This is a module example diagram of a multi-target track identification system under strong clutter background according to this application. Detailed Implementation

[0059] The subject matter described herein will now be discussed with reference to exemplary embodiments. It should be understood that these embodiments are discussed only to enable those skilled in the art to better understand and implement the subject matter described herein, and changes may be made to the function and arrangement of the elements discussed without departing from the scope of this specification. Various processes or components may be omitted, substituted, or added as needed in the examples. Furthermore, features described in some examples may be combined in other examples.

[0060] Example 1:

[0061] A method for preserving the identity of multiple targets' tracks against a background of strong clutter, such as Figure 1 As shown, it includes the following steps:

[0062] This embodiment applies to a monostation radar system in an urban low-altitude environment. The radar system operates in the X-band and is used to monitor the cross-maneuvering flight of multiple small multi-rotor UAVs in adjacent airspaces. These UAVs fly simultaneously in adjacent airspaces, and the background environment is subject to strong clutter interference such as multipath reflections from buildings and moving vehicle clutter. The radar system has completed the initial track establishment of the target and entered the track maintenance phase.

[0063] Step 100: Obtain the current UAV target's point trace and historical trajectory data; perform trajectory prediction based on the historical trajectory data to obtain the predicted trajectory position; calculate the distance between the point trace and the predicted trajectory position; and form a candidate point trace set by collecting points whose distance is less than the Mahalanobis distance threshold.

[0064] The radar system acquires the current moment's target trace data and historical trajectory data of UAV K. The trace data includes parameters such as range, azimuth, pitch, Doppler velocity, signal-to-noise ratio, and radar cross-section. The historical trajectory data includes the trajectory identifier, historical position sequence, velocity sequence, and historical statistical values ​​of radar cross-section for each trajectory. The selection of these input data is based on the actual needs of multi-target tracking in urban low-altitude environments. Among them, range, azimuth, and pitch provide spatial position information of the target, Doppler velocity reflects the target's motion state, and signal-to-noise ratio and radar cross-section are used to evaluate the reliability of the trace and distinguish real targets from clutter. The introduction of historical trajectory data can utilize the temporal continuity of target motion to provide prior information for current-moment correlation decisions. Especially in the context of strong clutter, historical information helps to maintain the stability and identity consistency of the trajectory.

[0065] For each trajectory j, an adaptive Kalman filter is used for trajectory prediction. The reason for choosing the adaptive Kalman filter algorithm is that small multi-rotor UAVs often perform maneuvering flights in urban low-altitude environments. Traditional fixed-parameter Kalman filters cannot effectively cope with the increased model uncertainty caused by target maneuvering. Adaptive Kalman filters can dynamically adjust the process noise covariance matrix according to the target's maneuvering state, thereby improving the tracking accuracy of maneuvering targets while maintaining filtering stability. The trajectory state vector includes the UAV target's position coordinates x, y, z in three-dimensional space, as well as the UAV target's velocity components vx, vy, vz in the three coordinate axis directions.

[0066] The prediction process is achieved through a state transition matrix, a control input matrix, and a control input vector. The predicted state vector at time k is equal to the product of the state transition matrix and the filtered estimated state vector at time k minus 1, plus the product of the control input matrix and the control input vector at time k minus 1. Here, the predicted state vector refers to the predicted state vector value at time k obtained through the Kalman filter prediction step, and the filtered estimated state vector refers to the estimated state vector value at time k minus 1 obtained through the Kalman filter update step. The predicted state vector includes the predicted position coordinates at time k and the predicted velocity components at time k, where the predicted velocity in the x-direction, y-direction, and z-direction corresponds to the fourth, fifth, and sixth components of the predicted state vector, respectively. These predicted velocity components are obtained from the velocity components in the filtered estimated state vector at time k minus 1 through the state transition matrix and are corrected according to the maneuvering acceleration in the control input vector.

[0067] The state transition matrix adopts a uniform motion model and is a 6x6 matrix. The first row of the matrix contains 1, 0, 0, sampling period, 0, 0; the second row contains 0, 1, 0, 0, sampling period, 0; the third row contains 0, 0, 1, 0, 0, sampling period; the fourth row contains 0, 0, 0, 1, 0, 0; the fifth row contains 0, 0, 0, 0, 1, 0; and the sixth row contains 0, 0, 0, 0, 0, 1. Since the motion of a small multi-rotor UAV within a short time interval can be approximated as uniform motion, and the uniform motion model has a simple structure and high computational efficiency, it can meet the requirements of real-time tracking. The preset sampling period is determined according to the data update rate of the radar system, that is, the preset sampling period is equal to the time interval corresponding to the reciprocal of the radar data update period.

[0068] The control input matrix consists of two parts: the upper part is the square of the preset sampling period divided by 2 and then multiplied by a third-order identity matrix; the lower part is the preset sampling period multiplied by a third-order identity matrix. The control input vector includes the target's maneuvering accelerations ax, ay, and az in the three coordinate axes. When a target maneuver is detected, the control input vector is set to the estimated acceleration vector calculated based on the difference between two adjacent historical filtered estimated velocity vectors; when no target maneuver is detected, the control input vector is set to zero. The purpose of introducing the control input is to compensate for the deviation of the target maneuver from the uniform motion model. By adding the maneuvering acceleration term during the prediction process, the prediction error during maneuvering can be reduced, and the accuracy of trajectory prediction can be improved.

[0069] The process of calculating the prediction covariance matrix is ​​as follows: the prediction covariance matrix at time k is equal to the product of the state transition matrix and the filter covariance matrix at time k minus 1, multiplied by the transpose of the state transition matrix, and finally added to the process noise covariance matrix. The process noise covariance matrix adopts a block diagonal matrix form, with the two submatrices on the diagonal being the position noise covariance submatrices and the velocity noise covariance submatrices, respectively, and the off-diagonal elements being all zero. These two submatrices are adaptively adjusted according to the target maneuver characteristics. Since the process noises of position and velocity are physically independent, the block diagonal matrix form not only conforms to the actual physical process but also simplifies the computational complexity. At the same time, by adaptively adjusting the submatrix parameters, the response characteristics of the filter can be flexibly adjusted for different maneuver states.

[0070] Based on the predicted covariance matrix and the target maneuver state, a Mahalanobis distance threshold is determined, and the process noise covariance matrix is ​​adaptively adjusted. Mahalanobis distance is used as the correlation threshold because it considers the covariance information of the prediction error and can automatically adjust the shape and size of the correlation gate according to the magnitude of the prediction uncertainty. Compared with Euclidean distance, it is more suitable for handling anisotropic error distributions. First, the estimated acceleration amplitude is calculated. This amplitude is based on the difference between the velocity components in the filtered estimated state vector at time k-1 and time k-2. Specifically, it is the square root of the sum of the squares of the differences between the filtered estimated velocity components at the two times, divided by a preset sampling period. The filtered estimated velocity components include the filtered estimated velocity in the x-direction and the y-direction... The velocity is estimated by filtering in the z-direction and in the z-direction. The reason for using filtered velocity estimates instead of predicted velocity to calculate acceleration is that filtered velocity estimates are corrected by observation information and are more accurate and reliable than pure predicted velocity estimates. This avoids the positive feedback bias problem that may occur when using the difference between predicted and actual velocity estimates to estimate acceleration, especially when there is strong clutter or correlation errors. Acceleration judgment based on filtered estimates is more robust. The estimated acceleration amplitude is used to determine whether the target is in a maneuvering state. When the estimated acceleration amplitude is greater than a preset maneuvering judgment acceleration threshold, the target is considered to be in a maneuvering state. This preset maneuvering judgment acceleration threshold is determined based on the target's normal flight acceleration level, so that it can effectively distinguish between normal flight and maneuvering flight.

[0071] When the target is in a maneuvering state, the process noise covariance matrix is ​​adaptively amplified. The amplified process noise covariance matrix is ​​equal to the original process noise covariance matrix multiplied by a maneuver amplification factor. The maneuver amplification factor is equal to 1 plus a preset maneuver noise gain coefficient multiplied by the ratio of the estimated acceleration amplitude to the preset maneuver judgment acceleration threshold. This preset maneuver noise gain coefficient is determined based on the degree of increase in uncertainty during maneuvering; the larger the coefficient, the greater the increase in process noise during maneuvering, thus increasing the predicted covariance matrix. By amplifying the process noise covariance matrix, the predicted covariance matrix naturally increases during maneuvering, thereby increasing the denominator in the Mahalanobis distance calculation, resulting in a relative decrease in the Mahalanobis distance value. This allows more candidate points to pass through the gating and enter the candidate set, avoiding maneuver-induced... To prevent false associations from being missed by real target points, although increasing prediction uncertainty would effectively loosen the association gate, this method suppresses false associations introduced by the loosened gate through point quality assessment in step 200 and global association optimization in step 400, thereby controlling the false alarm rate while ensuring that no maneuvering targets are missed. The Mahalanobis distance threshold is set as the threshold of the q-dimensional chi-square distribution that the square of the Mahalanobis distance follows under a preset confidence level, where q equals 3 corresponding to the three-dimensional position space. This preset confidence level is determined according to the reliability requirements of target detection. The higher the confidence level, the larger the association threshold, thereby reducing the probability of missed detection but potentially increasing false alarms. The threshold of the chi-square distribution of the corresponding dimension under the preset confidence level is obtained by looking up a table, and the square root of the threshold is then taken to obtain the Mahalanobis distance threshold.

[0072] For each point i, calculate the Mahalanobis distance between its observation vector and the predicted position of the track j. The original radar system measurements are in polar coordinates, including range, azimuth, and elevation, which need to be converted to three-dimensional position coordinates in a rectangular coordinate system. The coordinate transformation method is as follows: the x-coordinate of the point is equal to the range multiplied by the azimuth cosine and then by the elevation cosine; the y-coordinate is equal to the range multiplied by the azimuth sine and then by the elevation cosine; and the z-coordinate is equal to the range multiplied by the elevation sine. Since the coordinate transformation is nonlinear, the observation noise covariance matrix needs to be calculated for error propagation using the Jacobian matrix. Specifically, this is the product of the Jacobian matrix and the measurement noise covariance matrix in the polar coordinate system, multiplied by the transpose of the Jacobian matrix. The measurement noise covariance matrix in the polar coordinate system is a diagonal matrix, with the diagonal elements being the squares of the standard deviations of the range, azimuth, and elevation measurements, respectively. The Jacobian matrix is ​​a matrix composed of the partial derivatives of the coordinate transformation function with respect to the range, azimuth, and elevation.

[0073] The standard deviations of range measurement, azimuth measurement, and elevation measurement are all determined based on the measurement accuracy indicators of the radar system. The observation vector of track i includes the three-dimensional position coordinates of track i after coordinate transformation: px, py, and pz. The predicted position of track j is obtained by multiplying the observation matrix by the predicted state vector of track j. The observation matrix is ​​a 3x6 matrix, with the first three columns being a 3x3 identity matrix and the last three columns being a 3x3 zero matrix, used to extract the position component from the predicted state vector. The Mahalanobis distance is calculated by first calculating the difference between the observation vector and the predicted position, then multiplying the transpose of this difference by the inverse of the predicted observation covariance matrix, multiplying it by the difference itself, and finally taking the square root of the result. The predicted observation covariance matrix is ​​equal to the product of the observation matrix and the predicted covariance matrix multiplied by the transpose of the observation matrix, and finally adding the observation noise covariance matrix. The observation noise covariance matrix is ​​calculated using the Jacobian matrix error propagation method described above.

[0074] Select points whose Mahalanobis distance is less than the Mahalanobis distance threshold to form a candidate point set for track j.

[0075] Step 200: Evaluate the quality of the traces in the candidate trace set, extract the quality feature vectors of the traces, and calculate the trace quality score based on the quality feature vectors.

[0076] The quality of the traces in the candidate trace set for each track obtained in step 100 is evaluated. The quality feature vector of trace i is extracted, which includes signal-to-noise ratio, Doppler consistency, radar cross-section stability, and clutter probability. The reason for introducing trace quality evaluation is that in a strong clutter background, the candidate trace set contains both real target traces and clutter traces. By quantitatively evaluating the trace quality, a more reliable basis can be provided for subsequent global association optimization, reducing the interference of clutter on track association and improving the accuracy of track identity maintenance. The quality feature vector reflects the reliability of the trace from different perspectives. The signal-to-noise ratio reflects the detection intensity of the trace, the Doppler consistency reflects the degree of matching between the trace motion state and the track prediction, the radar cross-section stability reflects the consistency of the trace target characteristics, and the clutter probability combines environmental prior information and trace features to determine the possibility that the trace is clutter.

[0077] The signal-to-noise ratio is obtained directly from the radar spot data; Doppler consistency represents the degree of consistency between the Doppler velocity of the spot and the predicted velocity of the track. The calculation method is as follows: first, extract the predicted velocity in the x-direction, y-direction, and z-direction from the predicted state vector of track j; then, calculate the radar radial unit vector based on the radar position and the predicted target position; project the predicted velocity vector onto this radial unit vector to obtain the predicted velocity component of track j in the radar radial direction; next, calculate the absolute value of the difference between the observed Doppler velocity of the spot and the predicted velocity component; then, calculate the square of this difference and divide it by the square of twice the preset velocity measurement standard deviation; take the negative sign of the result and use it as the exponent of the natural exponential function for calculation; the preset velocity measurement standard deviation is determined according to the Doppler velocity measurement accuracy of the radar system; the higher the measurement accuracy, the smaller the standard deviation.

[0078] Radar cross-section (RCS) stability indicates the degree of closeness between the RCS value of a point track and the historical RCS mean of the track. It is calculated as follows: first, the absolute value of the difference between the observed RCS value of the point track and the historical RCS mean of track j is calculated; then, the square of this difference is divided by the square of twice the preset RCS measurement standard deviation. The result is then negative and used as the exponent of the natural exponential function. The preset RCS measurement standard deviation is determined based on the radar system's accuracy in measuring the RCS of small targets and is set with reference to the fluctuation range of the target's average RCS to reflect the fluctuation characteristics of the small target's RCS. It should be noted that due to multipath effects and rapid changes in target attitude in urban low-altitude environments, the RCS of small multi-rotor UAVs can fluctuate significantly, and the measurement noise of the RCS of small targets by a single-station radar is relatively high. Therefore, RCS stability is only used as an auxiliary feature for point track quality assessment, with a relatively low weight. Signal-to-noise ratio (SNR) and Doppler consistency are the primary criteria. RCS stability is used in conjunction with clutter probability to improve the reliability of point track quality assessment.

[0079] Clutter probability is determined by combining prior environmental information and clutter characteristics of the location of the track. For tracks located in known clutter areas, such as building reflection zones, the clutter probability is set to a higher value; for tracks with stable Doppler velocity and radar cross-section characteristics, the clutter probability is set to a lower value. The mechanism for acquiring and updating clutter areas is as follows: During the initialization phase, the system establishes an initial clutter area map based on building location and height information from the urban geographic information database. This clutter area map is stored in a three-dimensional grid format, with each grid recording the prior clutter probability at that location. During system operation, the clutter area map is updated using a sliding time window method by statistically analyzing the spatial distribution of tracks that have not been successfully associated with flight paths in historical track data. Specifically, spatial clustering is performed on unrelated points within a certain number of time periods of a preset clutter statistical time. Regions with point density exceeding a preset clutter density threshold are marked as clutter regions. The number of preset clutter statistical time periods is determined based on the time-varying characteristics of the clutter environment. The more time periods, the more stable the clutter region map is statistically, but the slower the response to environmental changes. The preset clutter density threshold is determined based on the false alarm rate of the radar system and the intensity of environmental clutter. Through this combination of offline calibration and online updating, the clutter region map can adapt to the non-stationary characteristics of clutter in urban environments. The clutter probability is calculated by first calculating the clutter indication value, which is equal to the weighted sum of the position indication term, the Doppler consistency term, and the radar cross-section stability term.

[0080] The position indication term is equal to a preset first clutter weighting coefficient multiplied by a position indication function. The input of the position indication function is the three-dimensional position coordinates of the point. It takes a value of 1 when the point is located in a clutter region, and a value of 0 otherwise. The Doppler consistency term is equal to 1 minus the Doppler consistency coefficient multiplied by a preset second clutter weighting coefficient. The radar cross-section stability term is equal to 1 minus the radar cross-section stability coefficient multiplied by a preset third clutter weighting coefficient. The three preset clutter weighting coefficients sum to 1. These three preset clutter weighting coefficients are determined based on the environmental clutter characteristics and the contribution of each point's characteristics to clutter discrimination. In the case of building multipath reflection and vehicle... In urban low-altitude environments with strong mobile clutter, prior position information contributes significantly to clutter discrimination. Therefore, the first clutter weighting coefficient should be set to the largest of the three. The contributions of Doppler consistency and radar cross-section stability are relatively low, so the second and third clutter weighting coefficients should be set to smaller values ​​accordingly. Since the position indication function takes the value of 0 or 1, and the values ​​of Doppler consistency and radar cross-section stability are both between 0 and 1, and the sum of the three weighting coefficients is equal to 1, the calculated clutter indication value naturally ranges between 0 and 1 and can be directly used as the clutter probability.

[0081] The fuzzy comprehensive evaluation method is used to calculate the quality score of the trace. Since the various features of trace quality have different dimensions and value ranges, and the influence of each feature on quality is fuzzy and uncertain, the fuzzy comprehensive evaluation method maps different features to a unified evaluation space through fuzzy membership functions. This effectively handles the uncertainty problem in multi-feature fusion and yields a more reasonable comprehensive evaluation result. First, the components of the quality feature vector are fuzzified, and a fuzzy membership function is defined. For the signal-to-noise ratio (SNR), an S-shaped membership function is used. This function is calculated by dividing 1 by the exponent term, where the exponent term is 1 plus the natural exponent function. The base of the natural exponent function is the natural constant. A preset steepness parameter with a negative exponent is multiplied by a difference term, which is the SNR minus a preset center value. The preset steepness parameter is determined based on the sensitivity of the membership function to changes in the SNR. A larger parameter results in a steeper change in the membership function, while a smaller parameter results in a smoother transition.

[0082] The preset center value is determined based on the typical signal-to-noise ratio (SNR) level of the radar system under strong clutter background. This preset center value corresponds to the inflection point where the SNR membership degree is 0.5, that is, when the SNR equals the preset center value, the membership degree is exactly 0.5. For Doppler consistency and radar cross-section stability, since their value range is already in the 0 to 1 interval, the original values ​​are directly used as the membership degree. For clutter probability, the inverse membership degree is adopted, that is, 1 minus the clutter probability, so that the lower the clutter probability, the higher the quality score.

[0083] The method for calculating the trace quality score is as follows: the inner product of the transpose of the preset basic weight vector and the fuzzy quality feature vector, after expansion, equals the preset first quality weight multiplied by the fuzzy membership degree of the signal-to-noise ratio, plus the preset second quality weight multiplied by the Doppler consistency, plus the preset third quality weight multiplied by the radar cross-section stability, plus 1 minus the difference of clutter probability, and then multiplied by the preset fourth quality weight. The preset basic weight vector includes the preset first quality weight, preset second quality weight, preset third quality weight, and preset fourth quality weight, and the sum of the four weights equals 1. These four preset quality weights are determined according to the contribution of each quality feature to the trace quality assessment. In target tracking under strong clutter background, the signal-to-noise ratio and Doppler consistency contribute more to distinguishing real targets from clutter, so the preset first quality weight and preset second quality weight should be set to larger values, while the contribution of radar cross-section stability and clutter probability is relatively low, so the preset third quality weight and preset fourth quality weight should be set to smaller values.

[0084] The weight vector is adjusted based on the motion mode feature vector of track j; based on the maneuvering state and straight flight tendency of the track, the weight of Doppler consistency is increased when the track is in a maneuvering state, and the weight of radar cross-section stability is increased when the track is in a straight flight state. The reason for introducing the adaptive weight adjustment mechanism is that the importance of each quality feature is different under different motion modes. Under maneuvering state, the target speed changes rapidly, and Doppler consistency can better reflect the matching relationship between the point track and the track. Under straight flight state, the target attitude is relatively stable, and radar cross-section stability becomes a more reliable discrimination criterion. By dynamically adjusting the weight according to the motion mode, the accuracy and adaptability of the point track quality assessment can be improved.

[0085] The adjusted weight vector equals the preset base weight vector plus a preset weight adjustment coefficient multiplied by the weight adjustment vector. This preset weight adjustment coefficient is determined based on the influence of the motion mode feature vector on the weight adjustment; a larger coefficient indicates a more significant adjustment effect of the motion mode feature vector on the weights, while a smaller coefficient indicates a more conservative weight adjustment. The weight adjustment vector is calculated based on the motion mode feature vector of track j. Specifically, it is calculated by first determining whether track j is in a maneuvering state. When the maximum acceleration of track j exceeds a preset maneuvering acceleration threshold, the track is considered to be in a maneuvering state. At this point, the second weight adjustment vector... The component corresponding to Doppler consistency is set to a positive value to increase its weight. When the probability of straight flight of track j is greater than the preset straight flight probability threshold, the track is considered to be in a straight flight state. At this time, the third component of the weight adjustment vector, which corresponds to the radar cross-section stability, is set to a positive value to increase its weight. When the weight of a certain quality feature is increased, the increase is deducted from the remaining components according to the basic weight ratio of the remaining quality features, and the weight vector after deduction is normalized so that each component is not less than zero and the weight sum is 1. This step uses the motion mode feature vector of the previous moment to adjust the weight.

[0086] The quality score of each candidate dot is obtained, with a value ranging from 0 to 1. The higher the value, the higher the quality of the dot.

[0087] Step 300: Model the motion pattern of each track, extract the motion pattern feature vector of the track based on the historical position sequence in the historical track data, and construct the track similarity matrix based on the motion pattern feature vector.

[0088] Motion pattern modeling is performed for each track, and motion pattern feature vectors are extracted based on the historical position sequence of the track. The motion pattern feature vectors include average velocity, maximum acceleration, track curvature, and straight flight probability. The motion pattern similarity between tracks is calculated based on the motion pattern feature vectors, and a track similarity matrix is ​​constructed.

[0089] The reason for introducing motion pattern modeling is that in multi-target cross-maneuver scenarios, the motion patterns of different targets often differ. By extracting and quantifying the motion pattern features of each track, important evidence can be provided for track identity verification. When the motion patterns of two tracks are highly similar, they are more prone to identity swapping errors, which need to be given special attention in association decision-making. The reason for choosing the four features of average speed, maximum acceleration, track curvature, and straight flight probability is that they can comprehensively describe the motion characteristics of the target. Average speed reflects the overall speed of the target's motion, maximum acceleration reflects the target's maneuverability, track curvature reflects the target's turning characteristics, and straight flight probability reflects the target's motion regularity. The combination of these features can effectively distinguish the motion patterns of different targets.

[0090] The historical position sequence of track j over several preset historical times includes all filtered estimated state vectors from k minus the preset historical time number to k minus 1. Each filtered estimated state vector includes three-dimensional position coordinates and three-dimensional velocity components. The preset historical time number is determined according to the time span of the target motion pattern. The more time numbers there are, the longer the motion pattern can be captured, but the computational load is also greater. It should be selected according to the typical duration of the target maneuvering flight to ensure that the motion pattern characteristics can be fully captured.

[0091] Calculate the average velocity of track j. The average velocity is the average velocity of track j over a number of preset historical times. The calculation method is as follows: for each time t from k minus the preset number of historical times to k minus 1, extract the three-dimensional velocity components from the filtered estimated state vector at that time, calculate the square root of the sum of squares of the three-dimensional velocity components at that time to obtain the velocity amplitude, then sum the velocity amplitudes over all times, and finally divide by the preset number of historical times. The three-dimensional velocity components include the x-direction velocity, y-direction velocity, and z-direction velocity of track j at time t.

[0092] Calculate the maximum acceleration of trajectory j. The maximum acceleration represents the maximum acceleration amplitude of trajectory j over a preset number of historical moments. The calculation method is as follows: for each moment t from k minus the preset number of historical moments plus 1 to k minus 1, extract the velocity vector from the filtered estimated state vector at moment t and t minus 1. Specifically, the filtered estimated state vector includes position and velocity components, and the velocity vector is composed of the three-dimensional velocity components in the filtered estimated state vector. When extracting the velocity vector, read the velocity components in the three coordinate directions according to the preset arrangement order of the components in the state vector, and combine the three read velocity components into the velocity vector at the corresponding moment. Calculate the difference between the velocity vector at moment t and the velocity vector at t minus 1, divide it by the preset sampling period to obtain the acceleration vector, then calculate the magnitude of the acceleration vector, and finally take the maximum value of the acceleration magnitude among all moments.

[0093] Calculate the curvature of trajectory j. The curvature represents the degree of turning of trajectory j. It is calculated by summing the curvature values ​​at each time t from k minus a preset historical time number plus 1 to k minus 1, and then dividing by the time difference term. The curvature value at time t is obtained by calculating the sine of the angle between the velocity vectors at adjacent times. The magnitude of the cross product of the velocity vector at time t and the velocity vector at time t minus 1 is calculated and then divided by the product of the magnitudes of the velocity vectors at time t and t minus 1. The time difference term is the preset historical time number minus 1. This curvature value reflects the degree of change in velocity direction between adjacent times. A larger cross product magnitude indicates a more drastic change in velocity direction. Normalization by dividing by the product of the magnitudes of the two velocity vectors ensures that the curvature is not affected by the magnitude of the velocity. Its value ranges from 0 to 1; a larger value indicates a more drastic turn.

[0094] Calculate the probability of straight flight for trajectory j. The probability of straight flight represents the tendency of trajectory j to maintain straight flight. It is calculated by summing the straight flight indices for each time t from k minus a preset historical time number plus one to k minus one, and then dividing by the time difference term. The straight flight indices for time t are calculated as follows: Calculate the angle between the velocity vectors at time t and t minus one, then calculate the square of this angle and divide it by an angle threshold term. The result is then negative and used as the exponent of the natural exponential function to calculate the straight flight indices for that time. The angle threshold term is twice the square of a preset angle threshold. The preset angle threshold is determined based on the heading stability of the target during straight flight and the accuracy of radar angle measurement. The smaller the threshold, the more sensitive it is to heading changes, and the larger the threshold, the more tolerant it is to heading changes. The selection of the preset angle threshold needs to take into account the influence of radar angle measurement error. When the standard deviation of radar azimuth and pitch angle measurements is large, the preset angle threshold should be appropriately increased to avoid misjudging heading jitter caused by measurement noise as maneuvering. Specifically, it can be adjusted according to the velocity direction estimation error obtained by converting the radar angle measurement accuracy, to ensure that the preset angle threshold is several times greater than the velocity direction estimation error.

[0095] The average speed, maximum acceleration, trajectory curvature, and straight-line flight probability obtained above are combined to form the motion pattern feature vector of trajectory j; this motion pattern feature vector is used for subsequent trajectory similarity calculation and association decision.

[0096] The similarity of motion patterns between tracks is calculated. For tracks j and l, the method for calculating their motion pattern similarity is as follows: First, the Euclidean distance between the motion pattern feature vectors of the two tracks is calculated. Specifically, after normalizing each component of the motion pattern feature vectors of the two tracks, the Euclidean distance between the normalized motion pattern feature vectors is calculated. Then, a negative preset decay coefficient is calculated and multiplied by the Euclidean distance, and the result is used as the exponent of the natural exponential function. The reason for using an exponential decay function based on Euclidean distance to calculate similarity is that Euclidean distance can comprehensively reflect the degree of difference between two feature vectors in a multidimensional feature space, while the exponential decay function can map distance to similarity. The smaller the distance, the higher the similarity. Moreover, the exponential function has good mathematical properties, which facilitates subsequent optimization calculations. The preset decay coefficient is determined according to the influence of the Euclidean distance of the motion pattern feature vectors on the similarity. The larger the coefficient, the faster the similarity decays with increasing distance, and the smaller the coefficient, the slower the similarity decays. The value range of similarity is greater than 0 and less than or equal to 1. The larger the value, the more similar the motion patterns of the two tracks are.

[0097] Construct a similarity matrix between tracks, where the element in the j-th row and l-th column of the matrix represents the similarity of motion patterns between track j and track l; this matrix is ​​used in step 500 to determine the similarity of motion patterns between candidate exchange tracks.

[0098] Step 400: Based on the candidate point set and point quality scores, perform global association optimization, construct a bipartite graph including the track set and point set, calculate the weights of the edges in the bipartite graph, and solve for the globally optimal association to obtain the association matrix between the track and point. Specifically, as follows... Figure 2 As shown.

[0099] Global association optimization is performed using the candidate point set obtained in step 100, the point quality score obtained in step 200, the motion pattern feature vector and the similarity matrix between tracks obtained in step 300. The reason for adopting a graph theory-based global association optimization method is that multi-target track association is essentially a combinatorial optimization problem, which requires finding the global optimal solution among all possible association schemes. The graph theory-based method can transform the association problem into a bipartite graph matching problem and solve it efficiently using mature graph theory algorithms such as the Hungarian algorithm. Compared with local optimization methods such as greedy association, global optimization can avoid local optimum traps and achieve better association performance in complex scenarios such as multi-target cross maneuvers.

[0100] A bipartite graph is constructed, comprising a vertex set and an edge set. The vertex set includes a track set and a point set, where the track set contains m tracks and the point set contains n points. The edge set contains multiple edges. If point i is in the candidate point set of track j, then there is an edge between track j and point i, indicating a candidate association between track j and point i. The reason for using a bipartite graph structure is that the association between tracks and points naturally constitutes a bipartite graph, with tracks and points located on opposite sides of the graph, and the association corresponding to the edges in the graph. This modeling method can clearly express the constraints of the association problem, namely, each track is associated with at most one point, and each point is associated with at most one track. When the number of tracks is not equal to the number of points, or when unassociated tracks are allowed, virtual points or virtual tracks are added to the cost matrix, enabling the Hungarian algorithm to output matching results that include unassociated states.

[0101] Calculate the weight of each edge, which comprehensively considers Mahalanobis distance and spot quality. The weight of the edge between track j and spot i is equal to the product of the preset first edge weight coefficient and the normalized distance similarity, plus the product of the preset second edge weight coefficient and the spot quality score, plus the product of the preset third edge weight coefficient and the track prediction consistency. The three preset edge weight coefficients add up to 1. These three preset edge weight coefficients are determined based on the importance of Mahalanobis distance, spot quality, and track prediction consistency to the association decision. In multi-target tracking against strong clutter backgrounds, Mahalanobis distance contributes the most to the association decision, so the preset first edge weight coefficient should be set to the largest of the three. The contribution of spot quality is second, so the preset second edge weight coefficient should be set to the middle value. The contribution of track prediction consistency is relatively low, so the preset third edge weight coefficient should be set to the smallest of the three. The normalized distance similarity is equal to 1 minus the distance term, which is the Mahalanobis distance divided by the Mahalanobis distance threshold. The smaller the Mahalanobis distance, the larger this term. The spot quality score is the result calculated in step 200.

[0102] The trajectory prediction consistency is the similarity between trajectory j and its historical trajectory. It is calculated by multiplying a first exponential term by a second exponential term. The base of the first exponential term is a natural constant. The square of the Euclidean distance between the three-dimensional position coordinates of the point with a negative exponent and the predicted trajectory position is divided by a position deviation term, which is the square of a preset position deviation threshold multiplied by 2. The base of the second exponential term is a natural constant. The square of the absolute value of the difference between the Doppler velocity of the point with a negative exponent and the radial predicted velocity component of the trajectory is divided by a velocity measurement term, which is the square of a preset velocity measurement standard deviation multiplied by 2. The radial predicted velocity component of the trajectory is obtained using the same method as in step 200 when calculating the Doppler consistency, i.e., by extracting the predicted velocity component from the predicted state vector and projecting it onto the radar radial unit vector. The preset position deviation threshold is determined based on the distance between the target in adjacent airspace and the radar measurement accuracy. This term integrates the historical consistency of position and velocity.

[0103] A Hungarian matching method based on competition penalty and edge weight adjustment is used to solve for the globally optimal association. Specifically, this embodiment does not change the row and column reduction, zero-coverage, and matching search rules of the Hungarian algorithm. Instead, before solving the Hungarian algorithm, the candidate edge weights are preprocessed according to the situation where multiple tracks compete for the same point, so that the cost matrix input to the Hungarian algorithm can reflect the risk of competitive association. Then, the Hungarian algorithm is used to solve the cost matrix after competition penalty adjustment to obtain the globally optimal association matrix that satisfies the one-to-one constraint. The reason for choosing the Hungarian algorithm is that it can solve the maximum weight matching problem of bipartite graphs in polynomial time, with a computational complexity of cubic level, which can meet the real-time requirements for the number of targets in practical applications. The reason for introducing the edge weight adjustment mechanism is that in multi-target close-range intersection scenarios, multiple tracks may simultaneously include a certain point in the candidate set. If they are associated indiscriminately, it is easy to cause track swapping errors. By penalizing the edge weights of competing points, the optimization algorithm can be guided to choose non-competitive association schemes, thereby reducing the risk of track swapping.

[0104] Set up an association matrix, which is an m x n matrix, where m is the number of tracks and n is the number of points. The association variable in the j-th row and i-th column of this matrix represents the association relationship between track j and point i, with a value of 1 indicating association and a value of 0 indicating no association. The objective function is to maximize the total association weight, which is the sum of the adjusted weights and the association variables for all tracks j from 1 to m and all points i from 1 to n.

[0105] The method for calculating the adjusted edge weights is as follows: First, determine whether point i is in the candidate point set of multiple paths. If point i is only in the candidate point set of path j, the adjusted edge weight is equal to the original edge weight. If point i is in the candidate point set of multiple paths, calculate a competition penalty factor for each competing path. The competition penalty factor is equal to 1 minus the competition term. The competition term is the difference between the number of competing paths and 1, multiplied by a preset competition penalty coefficient. Then, multiply the original edge weight by the competition penalty factor to obtain the adjusted edge weight. The number of competing paths is the number of paths corresponding to the candidate point set where point i is located. The preset competition penalty coefficient is used to avoid multiple paths competing for the same point. The intensity of demand is determined, and the larger the coefficient, the heavier the penalty for competition. To avoid negative or excessively small values ​​of the competition penalty factor leading to abnormal edge weights, a lower limit is set for the competition penalty factor. When the calculated competition penalty factor is less than the preset lower limit, the competition penalty factor is truncated to the preset lower limit. This preset lower limit is determined based on the requirement to ensure the effectiveness of edge weights. It should ensure that the truncated edge weights can still participate in global optimization and will not lose their distinguishability due to being too small. In this way, when multiple tracks compete for the same point, the edge weights of all competing tracks and that point will be reduced, thus favoring the selection of non-competitive association schemes in global optimization.

[0106] The constraints include three terms: The first term is that each point can be associated with at most one track, that is, the sum of the associated variables of all tracks j from 1 to m is less than or equal to 1, and this constraint holds for all points i; The second term is that each track can be associated with at most one point, that is, the sum of the associated variables of all points i from 1 to n is less than or equal to 1, and this constraint holds for all tracks j; The third term is that the associated variables are binary variables, which can only take the values ​​0 or 1.

[0107] Based on the objective function, the correlation matrix, and the constraints, a cost matrix is ​​constructed for the Hungarian algorithm to solve. For track-point pairs with edge connections, the corresponding elements in the cost matrix are set to the negative of the adjusted edge weights to transform the problem of maximizing the total correlation weights into the problem of minimizing the total cost. For track-point pairs without edge connections, the corresponding elements are set to prohibited costs, which are preset positive numbers greater than the absolute values ​​of all valid costs, used to prevent non-candidate correlations from entering the optimal matching result. When the number of tracks is not equal to the number of points, or when tracks are allowed to be unassociated with points at the current time, virtual points or virtual tracks are added to the cost matrix to make the expanded cost matrix a square matrix. Virtual points or virtual tracks correspond to empty correlations, and their costs are set according to the requirements for missing detection retention, so that real valid candidate correlations take precedence over empty correlations, and non-candidate correlations are inferior to empty correlations. The expanded cost matrix is ​​solved using the Hungarian algorithm to obtain the expanded matching result. After removing the virtual node matching relationships in the expanded matching result, the correlation matrix between tracks and points is obtained.

[0108] The Hungarian algorithm outputs an optimal correlation matrix; this optimal correlation matrix is ​​an m x n matrix, where the element in the j-th row and i-th column is 1 if track j is associated with point i, and 0 if not. Specifically, during the solution process, the Hungarian algorithm determines the optimal matching point number for each track. For track j, if the algorithm determines that its optimal matching point number is i, then the element in the j-th row and i-th column of the optimal correlation matrix is ​​set to 1, and all other elements in that row are set to 0; if the algorithm determines that track j does not match any point, then the optimal correlation matrix is ​​set to... All elements in row j are set to 0; the optimal correlation matrix clearly expresses the correlation between each track and point track, providing a basis for subsequent track authentication and status updates; tracks not associated with any point tracks are marked as missed detections; for point tracks not associated with any track, it is determined whether they are new targets or clutter; the rule for determining the start of a new target is that if an unassociated point track has unassociated point tracks with similar positions and consistent movement trends within a certain number of consecutive preset new track confirmation frames, then the point track sequence is considered to correspond to... For a new target, a new track is established and a new track identifier is assigned. The preset number of new track confirmation frames is determined based on the reliability requirements of new target confirmation. More frames result in more reliable new track confirmation but a slower response, while fewer frames result in a faster response but may introduce false tracks. The criterion for judging proximity is that the Euclidean distance between unrelated points at adjacent times is less than a preset new track position threshold. This preset new track position threshold is determined based on the target's maximum displacement between adjacent times and the radar measurement accuracy. The criterion for judging consistent motion trends is that the three-dimensional coordinates of adjacent points... The angle between the velocity vector calculated by the standard deviation and the velocity vector at the previous moment is less than the preset new trajectory heading threshold. The preset new trajectory heading threshold is determined according to the heading continuity requirements of the target motion. The clutter judgment rule is that if the quality score of an unassociated point is lower than the preset clutter quality threshold, or if the point is located in the clutter region and cannot form a continuous point sequence, it is judged as clutter and discarded. The preset clutter quality threshold is determined according to the difference in quality score distribution between clutter points and real target points, and should be set in the boundary region between the quality score distributions of the two types of points.

[0109] Step 500: Verify the association matrix based on the similarity matrix between the tracks. When the exchange condition is triggered, start the track identity correction mechanism to obtain the corrected association matrix.

[0110] The optimal correlation matrix obtained in step 400 is verified and corrected to avoid track swapping and merging errors. The reason for introducing the track identity verification and correction mechanism is that although global correlation optimization can obtain the mathematically optimal solution, in strong clutter backgrounds and multi-target cross-maneuver scenarios, track identity errors may still occur if correlation is based solely on Mahalanobis distance and track quality. In particular, when two tracks are very close and have similar motion patterns, track swapping is likely to occur. By verifying the physical rationality of the correlation results in the correlation matrix and checking the consistency of motion patterns, potential identity errors can be detected and corrected in a timely manner, improving the reliability of track identity maintenance.

[0111] The exchange conditions are judged, including verifying the rationality of each association pair in the association matrix; for each association pair of track j and point i corresponding to each element with a value of 1 in the optimal association matrix, it is checked whether the change in track speed after association exceeds the target's maximum maneuverability; the reason for introducing speed change check is that the target's maneuverability is limited by physical constraints. If the speed change caused by association exceeds the target's maximum maneuverability and acceleration, it indicates that the association is physically unreasonable and is likely an incorrect association; the speed change after association is calculated, which is equal to the difference between the three-dimensional position coordinates of the point and the position component extracted from the filtered estimated state vector of the track at time k minus 1, and then divided by the preset sampling period to obtain the speed vector after association. Then, the speed component is extracted from the filtered estimated state vector of track j at time k minus 1 as the speed vector at time k minus 1, and the magnitude of the difference between the speed vector after association and the speed vector at time k minus 1 is calculated to obtain the speed change.

[0112] If the speed change is greater than the preset maximum maneuvering acceleration capability multiplied by the preset sampling period, the association is considered unreasonable, triggering the exchange condition and initiating the track identity correction mechanism. The preset maximum maneuvering acceleration capability is determined based on the target's maneuverability and should be set with reference to the target type's physical maneuverability limits.

[0113] The system checks whether the heading change of the associated track conforms to the motion mode feature vector. It calculates the heading change angle by first calculating the associated velocity vector as the three-dimensional position coordinates of the track minus the position component extracted from the filtered estimated state vector of the track at time k-1, then dividing by the preset sampling period. Next, it extracts the velocity component from the filtered estimated state vector of track j at time k-1 as the velocity vector at time k-1. Then, it calculates the inner product of the associated velocity vector and the velocity vector at time k-1, divides it by the product of the magnitudes of the two velocity vectors, and takes the inverse cosine function of the result to obtain the heading change angle. If the heading change angle is greater than the preset maximum heading change angle, and the probability of straight flight of track j is greater than the preset straight flight probability threshold, the association is considered unreasonable, triggering the exchange condition and initiating the track identity correction mechanism. The preset maximum heading change angle is determined based on the target's heading stability in straight flight; the smaller the angle, the stricter the constraint on heading change. The preset straight flight probability threshold is determined based on the confidence requirement for judging whether the target is in straight flight; this threshold should be set at a level that reliably distinguishes between straight flight and maneuvering flight.

[0114] Detecting track swapping risk: For any two linked tracks j and l in the optimal association matrix, if track j is associated with track i and track l is associated with track p, calculate the swapping risk index. The reason for introducing swapping risk detection is that in multi-target close-range intersection scenarios, two tracks may simultaneously approach each other's real target tracks. If the association decision is inappropriate, track identity swapping is likely to occur. By calculating the ratio of the total Mahalanobis distance after the swapped association to the current total Mahalanobis distance, the magnitude of the swapping risk can be quantitatively assessed. The reason for using Mahalanobis distance instead of Euclidean distance to calculate the swapping risk index is that Mahalanobis distance considers prediction covariance information and can normalize the distance according to the prediction uncertainty in each direction. Under conditions of clutter background and large single-station radar altimetry error, Mahalanobis distance is more robust to measurement noise than Euclidean distance and can more accurately reflect the statistical distance between point tracks and tracks. The method for calculating the exchange risk index is as follows: First, calculate the Mahalanobis distance between point p and track j. The calculation method for this Mahalanobis distance is the same as in step 100, that is, calculate the difference between the observation vector of point p and the predicted position of track j, then calculate the product of the transpose of this difference and the inverse of the predicted observation covariance matrix of track j, multiply it by the difference itself, and finally take the square root of the result. The Mahalanobis distance between point i and track l, the Mahalanobis distance between point i and track j, and the Mahalanobis distance between point p and track l are calculated in the same way.

[0115] Before calculating the exchange risk index, boundary condition processing of the denominator is required to avoid numerical instability. When the Mahalanobis distance between point i and track j plus the Mahalanobis distance between point p and track l is less than the preset minimum Mahalanobis distance and threshold, it indicates that the predictions of both tracks in the current association scheme are very accurate and the points fall exactly near the predicted positions. At this point, the current association scheme is considered optimal, and there is no track exchange risk. The exchange risk index is directly set to a value greater than the preset exchange risk threshold to skip subsequent exchange risk judgments. The preset minimum Mahalanobis distance and threshold are determined according to the requirements of numerical calculation stability and should be set to positive numbers close to zero to avoid division by zero errors. When the denominator... When the exchange risk index is greater than or equal to the preset minimum Mahalanobis distance and threshold, the exchange risk index is calculated normally. The exchange risk index is equal to the Mahalanobis distance between point p and track j plus the Mahalanobis distance between point i and track l, and then divided by the Mahalanobis distance between point i and track j plus the Mahalanobis distance between point p and track l. This index represents the ratio of the total Mahalanobis distance after the exchange association to the current total Mahalanobis distance. If the exchange risk index is less than the preset exchange risk threshold, it means that the Mahalanobis distance after the exchange association is smaller, and there is a risk of track exchange. The preset exchange risk threshold is determined based on the sensitivity of judging track exchange risk. The larger the threshold, the more sensitive the detection of exchange risk, and the smaller the threshold, the more conservative the detection.

[0116] Simultaneously, the similarity of motion patterns is used for judgment. If the similarity of motion patterns between track j and track l is greater than the preset motion pattern similarity threshold, it indicates that the motion patterns of the two tracks are very similar, and the risk of exchange is further increased. The preset motion pattern similarity threshold is determined based on the confidence requirement for judging whether the motion patterns of the two tracks are similar. The higher the threshold, the stricter the similarity judgment. The comprehensive judgment condition is that the exchange risk index is less than the preset exchange risk threshold and the motion pattern similarity is greater than the preset motion pattern similarity threshold. When the above exchange conditions are triggered, the track identity correction mechanism is activated.

[0117] The correction process of the track identity correction mechanism is as follows: Calculate the comprehensive score after the exchange association; for the original association scheme, the comprehensive score equals the edge weight between track j and point i plus the edge weight between track l and point p; for the exchange association scheme, the comprehensive score equals the edge weight between track j and point p plus the edge weight between track l and point i; if the comprehensive score of the exchange association scheme is greater than the comprehensive score of the original association scheme, then the association exchange is performed, that is, in the optimal association matrix, the elements originally associated with track j to point i are set to 0, and the elements originally associated with track j to point i are changed to 0. The element corresponding to point p is set to 1, while the element originally associated with point p in track l is set to 0, and the element corresponding to point i in track l is set to 1, thereby realizing the exchange of track j associated with point p and track l associated with point i; for the track and point combination in the relationship to be exchanged, if the combination already exists in the edge set of the bipartite graph, the corresponding edge weight is used; if the combination does not exist in the edge set of the bipartite graph, its edge weight is recalculated, and its comprehensive score is set to invalid value when its Mahalanobis distance does not meet the Mahalanobis distance threshold.

[0118] After the above verification and correction, the corrected correlation matrix is ​​obtained. The corrected correlation matrix has the same structure as the optimal correlation matrix, which is still an m-row n-column matrix. The element in the j-th row and i-th column takes a value of 1 to indicate that track j is associated with point i, and takes a value of 0 to indicate that they are not associated. The difference is that some correlation pairs in the corrected correlation matrix may have been swapped to ensure the correctness of track identity.

[0119] Based on the corrected correlation matrix, the state of each track is updated. For track j with a value of 1 in the j-th row and i-th column of the corrected correlation matrix, it indicates that track j has successfully correlated with track i, and an adaptive Kalman filter is used for state update. The reason for using Kalman filtering for state update is that Kalman filtering can fuse prediction and observation information, and provide the optimal state estimate in the sense of minimum mean square error. By automatically balancing the weights of prediction and observation through Kalman gain, it relies more on observation when prediction uncertainty is high and more on prediction when observation noise is high, thus obtaining a stable and reliable state estimate. First, the Kalman gain is calculated. The Kalman gain is equal to the product of the prediction covariance matrix and the transpose of the observation matrix, and then multiplied by the observation term, which is the sum of the observation matrix and the prediction covariance matrix. The product of the two matrices is multiplied by the transpose of the observation matrix and then added to the inverse of the observation noise covariance matrix; the prediction covariance matrix is ​​the result calculated in step 100; the observation matrix is ​​a 3x6 matrix, with the first three columns being a 3x3 identity matrix and the last three columns being a 3x3 zero matrix; the observation noise covariance matrix is ​​set according to the radar measurement accuracy and is a 3x3 diagonal matrix, with the diagonal elements being the squares of the preset x-direction position measurement standard deviation, the preset y-direction position measurement standard deviation, and the preset z-direction position measurement standard deviation, respectively, and the off-diagonal elements being all zero. The preset x-direction position measurement standard deviation and the preset y-direction position measurement standard deviation are determined according to the radar system's measurement accuracy in the horizontal direction, and the preset z-direction position measurement standard deviation is determined according to the radar system's measurement accuracy in the vertical direction.

[0120] The track state is updated. The updated track state vector is equal to the predicted state vector plus the Kalman gain multiplied by the observation residual. The observation residual is equal to the point observation vector minus the product of the observation matrix and the predicted state vector, representing the difference between the actual observation and the predicted observation.

[0121] Update the covariance matrix. The updated covariance matrix is ​​equal to the identity matrix minus the difference between the Kalman gain and the observation matrix, and then multiplied by the prediction covariance matrix.

[0122] For track j in the corrected correlation matrix where all elements in row j are 0, it indicates that track j is not associated with any point and is in a missed detection state. The predicted state remains unchanged, i.e., the updated state vector is equal to the predicted state vector, but the covariance matrix is ​​increased to reflect the increase in uncertainty. The reason for increasing the covariance matrix for missed detection tracks is that when a track is not associated with any observation point, there is a lack of correction based on observation information, and the uncertainty of the state estimation will increase over time. By increasing the covariance matrix, the correlation threshold can be expanded in subsequent time steps, increasing the probability of re-acquiring the target, and also providing a basis for track quality assessment. The updated covariance matrix is ​​equal to the predicted covariance matrix plus the missed detection penalty matrix. The missed detection penalty matrix is ​​equal to the process noise covariance matrix multiple of the preset missed detection penalty multiple. This preset missed detection penalty multiple is determined according to the degree of increase in uncertainty in the missed detection state. The larger the multiple, the greater the increase in uncertainty after the missed detection, and the smaller the multiple, the more conservative the increase in uncertainty.

[0123] The motion pattern feature vector of the track is updated. The filtered estimated state vector at the current moment is added to the historical state sequence. Based on the updated historical position sequence, the average velocity, maximum acceleration, track curvature, and straight flight probability are recalculated to obtain the updated motion pattern feature vector, which is used for the association decision at the next moment. The reason for continuously updating the motion pattern feature vector is that the target's motion pattern may change over time. By updating the motion pattern feature vector in real time through a sliding window, the changes in the target's motion pattern can be captured in a timely manner, providing accurate motion pattern information for the subsequent adjustment of the point quality assessment weights and track identity verification.

[0124] Maintain track identity information; each track retains a unique track identifier, which is assigned during track initialization and remains unchanged throughout the tracking process; through the authentication and correction mechanism in this step, the consistency between the track identifier and the target's true identity is ensured, avoiding track exchange and merging errors; the reason for maintaining a unique track identifier is that the core objective of track identity maintenance is to ensure that each track always tracks the same real target. As the unique identifier of the track identity, the track identifier must remain unchanged throughout the tracking process. Through the synergistic effect of the aforementioned steps, the correspondence between the track identifier and the target's true identity is ensured not to be disrupted.

[0125] Update the quality assessment parameters of the track; calculate the number of consecutive associations and the number of missed detections of the track; if the number of consecutive associations exceeds the preset consecutive association threshold, increase the reliability of the track. The preset consecutive association threshold is determined based on the requirements for judging track stability. The larger the threshold, the higher the requirement for track stability. If the number of missed detections exceeds the preset missed detection threshold, decrease the reliability of the track and consider deleting the track. The preset missed detection threshold is determined based on the fault tolerance of track maintenance. The smaller the threshold, the lower the tolerance for missed detections. It should be set according to the typical duration of target obstruction or missed detection in clutter environments.

[0126] Output updated track status data, including track identifier, position, velocity, acceleration, covariance matrix, motion pattern feature vector, etc.; at the same time, record track processing logs, including correlation decision process, authentication results, status update parameters, etc., for subsequent analysis and system optimization.

[0127] Example 2:

[0128] See Figure 3 As shown, a multi-target track identification system under strong clutter background is provided, which is used to execute the aforementioned multi-target track identification method under strong clutter background. The system includes:

[0129] The dynamic threshold module 101 acquires the current UAV target's traces and historical trajectory data, performs trajectory prediction based on the historical trajectory data to obtain the predicted trajectory position, calculates the distance between the traces and the predicted trajectory position, and forms a candidate trace set for traces whose distance is less than the Mahalanobis distance threshold.

[0130] The quality assessment module 102 assesses the quality of the traces in the candidate trace set, extracts the quality feature vectors of the traces, and calculates the trace quality score based on the quality feature vectors.

[0131] The motion modeling module 103 performs motion pattern modeling for each track, extracts motion pattern feature vectors of the track based on historical position sequences in historical track data, and constructs a similarity matrix between tracks based on the motion pattern feature vectors.

[0132] The association optimization module 104 performs global association optimization based on the candidate point set and point quality score, constructs a bipartite graph including the track set and the point set, calculates the weight of the edges in the bipartite graph and solves the global optimal association, and obtains the association matrix of the track and the point.

[0133] The verification and correction module 105 verifies the association matrix based on the similarity matrix between the tracks. When the exchange condition is triggered, the track identity correction mechanism is activated to obtain the corrected association matrix.

[0134] The embodiments of this application have been described above, but these embodiments are not limited to the specific implementation methods described above. The specific implementation methods described above are merely illustrative and not restrictive. Those skilled in the art can make more equivalent embodiments based on the guidance of these embodiments, and all of them are within the protection scope of these embodiments.

Claims

1. A method for preserving the identity of multiple targets' tracks against a background of strong clutter, characterized in that, Includes the following steps: Acquire the current UAV target's point traces and historical flight path data, perform flight path prediction based on the historical flight path data to obtain the predicted flight path position, calculate the distance between the point traces and the predicted flight path position, and form a candidate point trace set for points whose distance is less than the Mahalanobis distance threshold. The quality of the traces in the candidate trace set is evaluated, the quality feature vector of the traces is extracted, and the trace quality score is calculated based on the quality feature vector. Motion pattern modeling is performed for each track, and motion pattern feature vectors of the track are extracted based on historical position sequences in historical track data. A similarity matrix between tracks is constructed based on the motion pattern feature vectors. Global association optimization is performed based on the candidate point set and point quality score. A bipartite graph including the track set and point set is constructed. The weights of the edges in the bipartite graph are calculated and the global optimal association is solved to obtain the association matrix of track and point. The correlation matrix is ​​verified based on the similarity matrix between the tracks. When the exchange condition is triggered, the track identity correction mechanism is activated to obtain the corrected correlation matrix.

2. The method for maintaining the identity of multiple targets under strong clutter background according to claim 1, characterized in that, Methods for forming a set of candidate points include: Obtain distance, azimuth, and elevation angle information from the point data; The three-dimensional position coordinates of the point are determined based on the distance, azimuth, and pitch angle information, and the Mahalanobis distance between the three-dimensional position coordinates of the point and the predicted position of the flight path is calculated. Points whose Mahalanobis distance is less than the Mahalanobis distance threshold are included in the candidate point set of the corresponding track.

3. The method for maintaining the identity of multiple targets under strong clutter background according to claim 1, characterized in that, Methods for calculating dot quality scores include: Extract the quality feature vector, which includes signal-to-noise ratio, Doppler consistency, radar cross-section stability, and clutter probability; Each component of the quality feature vector is fuzzified. An S-shaped membership function is used for the signal-to-noise ratio. The original values ​​of Doppler consistency and radar cross-section stability are directly used as membership. The reverse membership function is used for clutter probability. The quality score of the dots is calculated by the inner product of the preset basic weight vector and the fuzzy quality feature vector.

4. The method for maintaining the identity of multiple targets under strong clutter background according to claim 1, characterized in that, Methods for extracting motion pattern feature vectors from tracks include: The average velocity is calculated based on the three-dimensional velocity components in the historical location sequence and the preset number of historical moments. Based on the three-dimensional velocity components, the velocity vector is extracted, and the difference between the velocity vectors at adjacent time points is divided by the sampling period to obtain the acceleration magnitude. Within the adjacent time points corresponding to the preset number of historical time points, the maximum value of the acceleration magnitude is the maximum acceleration. The mean value of the sine of the angle between the velocity vectors at each adjacent time point is calculated to obtain the track curvature. The straight-line flight index is calculated based on the angle between the velocity vectors of adjacent moments, and the straight-line flight probability is obtained by calculating the mean of the straight-line flight index at each moment within a preset historical time period. The motion pattern feature vector includes average velocity, maximum acceleration, trajectory curvature, and probability of straight flight.

5. A method for maintaining the identity of multiple targets under strong clutter background according to claim 4, characterized in that, Methods for constructing a similarity matrix between tracks include: For any two tracks, calculate the Euclidean distance between the motion pattern feature vectors of the two tracks, and calculate the motion pattern similarity by multiplying the Euclidean distance by a natural exponential function with a negative preset attenuation coefficient. The similarity matrix between the tracks is constructed using the similarity of the motion patterns as elements.

6. The method for maintaining the identity of multiple targets under strong clutter background according to claim 1, characterized in that, Methods for calculating the weights of edges in a bipartite graph include: Construct a bipartite graph, which includes a vertex set and an edge set. The vertex set includes a track set and a point set, and the edge set includes multiple edges. If point i is in the candidate point set of track j, then there is an edge between track j and point i. The normalized distance similarity is calculated based on the distance between the predicted location of the point and the track, and the Mahalanobis distance threshold. The consistency of trajectory prediction is calculated based on the three-dimensional position coordinates of the points, the predicted trajectory position, and the preset position deviation threshold. The weight of the edge is a weighted sum of normalized distance similarity, point quality score, and trajectory prediction consistency.

7. A method for maintaining the identity of multiple targets under strong clutter background as described in claim 6, characterized in that, Methods for obtaining the correlation matrix between tracks and points include: When a point is in the candidate point set of multiple paths at the same time, a competition penalty factor is calculated based on a preset competition penalty coefficient, and the weight of the original edge is multiplied by the competition penalty factor to obtain the adjusted edge weight. Set up an association matrix. The association variable in the j-th row and i-th column of the association matrix represents the association relationship between track j and point i. The association relationship corresponds to the edge in the bipartite graph. The association weights are calculated based on the adjusted weights and related variables; With the goal of maximizing the total association weight, and under the constraint that each point can be associated with at most one track and each track can be associated with at most one point, the association matrix is ​​obtained by using the Hungarian algorithm.

8. A method for maintaining the identity of multiple targets under strong clutter background according to claim 1, characterized in that, The methods for determining whether a swap condition is triggered include: Calculate the velocity change and heading change angle based on the three-dimensional position coordinates of the point trace; For each association pair in the association matrix, if the speed change is greater than the product of the preset maximum maneuvering acceleration capability and the sampling period, or the heading change angle is greater than the preset maximum heading change angle and the straight flight probability of the trajectory is greater than the preset straight flight probability threshold, then the exchange condition is triggered; wherein, the association pair is the trajectory j and point i that have an association relationship in the association matrix; Calculate the exchange risk index based on the Mahalanobis distance between the point and the track; For any two associated tracks in the current association matrix, if the exchange risk index is less than the preset exchange risk threshold, and the similarity of the motion patterns of the two associated tracks is greater than the preset motion pattern similarity threshold, then the exchange condition is triggered.

9. A method for maintaining the identity of multiple targets under strong clutter background as described in claim 8, characterized in that, Methods for obtaining the corrected correlation matrix include: For track j and track l that trigger the exchange condition, and their associated points i and p in the current association matrix, determine the original association relationship and the association relationship to be exchanged; the original association relationship is that track j is associated with point i and track l is associated with point p, and the association relationship to be exchanged is that track j is associated with point p and track l is associated with point i. Calculate the comprehensive score of the original association and the association to be exchanged. If the comprehensive score of the association to be exchanged is greater than the comprehensive score of the original association, then update the original association to the association to be exchanged in the current association matrix to obtain the corrected association matrix. The comprehensive score is the sum of the edge weights of the two tracks and their respective associated point tracks.

10. A multi-target track identification system under strong clutter background, characterized in that, It is used to implement the multi-target track identification method under strong clutter background as described in any one of claims 1-9, the system comprising: The dynamic threshold module acquires the current UAV target's traces and historical trajectory data, performs trajectory prediction based on the historical trajectory data to obtain the predicted trajectory position, calculates the distance between the traces and the predicted trajectory position, and forms a candidate trace set for traces whose distance is less than the Mahalanobis distance threshold. The quality assessment module assesses the quality of the traces in the candidate trace set, extracts the quality feature vectors of the traces, and calculates the trace quality score based on the quality feature vectors. The motion modeling module models the motion pattern of each track, extracts the motion pattern feature vector of the track based on the historical position sequence in the historical track data, and constructs a similarity matrix between tracks based on the motion pattern feature vector; The association optimization module performs global association optimization based on the candidate point set and point quality score, constructs a bipartite graph including the track set and the point set, calculates the weight of the edges in the bipartite graph and solves the global optimal association, and obtains the association matrix of track and point. The verification and correction module verifies the association matrix based on the similarity matrix between the tracks. When the exchange condition is triggered, the track identity correction mechanism is activated to obtain the corrected association matrix.

Images