A robust extended object tracking method based on time-series double-ellipse constraint

Abstract

The application discloses a kind of robust extended target tracking method based on timing double-ellipse constraint, belongs to passive perception extended target tracking technical field, it obtains the angle of arrival and propagation distance observation value of direct path and target reflection path;Estimate the position of non-cooperative transmitter, and accordingly solve the spatial position of each reflection point;With parameter vector indicating the quadratic curve of extended target, construct algebraic residual;Single-ellipse estimation is carried out to obtain core parameter estimation value;When non-initial time, adjacent time core parameter and algebraic residual are used to build double-ellipse estimation target function, with the constraint that quadratic curve of each time is ellipse and adjacent time major and minor axis are equal;Adaptive weight is introduced to obtain joint optimization problem;Decomposed into parameter vector subproblem and adaptive weight subproblem, and alternately iterated solution is solved;The center position, size and attitude of extended target are recovered from the current time parameter vector;The application has robustness to measurement bias and wild value, and high tracking stability.

Description

Technical Field

[0001] This invention relates to the field of passive sensing extended target tracking technology, and in particular to a robust extended target tracking method based on temporal double ellipse constraints. Background Technology

[0002] Extended targets typically generate multiple spatially resolved measurement points at a single sampling moment. Therefore, in addition to estimating the center position and trajectory of the extended target, the tracking system also needs to simultaneously acquire spatial extension information such as the target's physical dimensions and orientation. This information is crucial for low-altitude aircraft surveillance, vehicle collision avoidance, ship monitoring, and situational awareness in complex environments. Elliptical models describe the center, major and minor axes (physical dimensions), and attitude angles (orientation) of extended targets with fewer parameters, achieving a good balance between model expressive power and real-time computational overhead.

[0003] In passive sensing scenarios, receivers utilize non-cooperative transmitters in unknown locations (such as broadcast stations, communication base stations, or navigation signal sources in unknown locations) for detection. Typically, they can only obtain one angle of arrival and time of arrival from the direct path, and multiple angles of arrival and times of arrival from multiple target reflection paths, but cannot directly obtain the Cartesian coordinates of each reflection point on the extended target. Because the location of the non-cooperative transmitter is unknown, it is necessary to first estimate its position, then use the geometric relationships of the target reflection paths to calculate the positions of each reflection point, and finally perform a joint estimation of the extended target's center, major and minor axes, and attitude angles.

[0004] Due to occlusion, shadowing effects, and illumination direction, reflection points tend to concentrate on the side of the extended target facing the non-cooperative transmitter, forming an incomplete and significantly asymmetrical measurement distribution, resulting in a systematic bias in the residuals. Furthermore, outliers may be introduced into the measurements. When using only current-time data for least-squares elliptic fitting, the center and semi-major and semi-minor axes of the extended target are prone to shifts. Tracking methods relying on zero-mean Gaussian error assumptions or uniform reflection assumptions experience significant performance degradation under conditions of sparse measurements, unilateral distribution, or the presence of impulse interference. Moreover, while existing methods construct the current-time prior from the previous time-wise posterior result, they do not fully utilize the physical property that the semi-major and semi-minor axes of the extended target remain essentially unchanged between adjacent time-wise moments.

[0005] Therefore, there is a need for an extended target tracking method that can calculate the position of the reflection point from observations of the angle of arrival and time of arrival of the direct path and the target reflection path, maintain robustness against measurement bias and outliers, and improve tracking stability through consistency constraints of the semi-major and semi-minor axes at adjacent time points. Summary of the Invention

[0006] For passive sensing extended target tracking scenarios where the location of the non-cooperative transmitter is unknown, the distribution of reflection points on the extended target is incomplete, and the measurement may contain outliers, this invention provides a robust extended target tracking method based on temporal double ellipse constraints. It can calculate the position of the reflection point from the arrival angle and arrival time observations of the direct path and the target reflection path. It can address the measurement bias caused by the unilateral distribution of reflection points through kernel parameter adaptive estimation, suppress the influence of outliers through adaptive weights, and improve tracking stability through the consistency constraint of the major and minor semi-axes of adjacent time moments.

[0007] The technical solution adopted by this invention to solve the above-mentioned technical problems is: a robust extended target tracking method based on temporal double ellipse constraints, comprising the following steps: S1: The receiver acquires the current observation values ​​of the angle of arrival and propagation distance of the direct path and the reflection paths of each target; S2: Estimate the position of the non-cooperative transmitter at the current moment based on the observed angle of arrival and propagation distance of the direct path; S3: Based on the observed angle of arrival and propagation distance of each target reflection path at the current moment, and the position of the non-cooperative transmitter at the current moment, calculate the spatial position of each reflection point at the current moment; S4: Use a parameter vector to represent the quadratic curve of the extended target at the current moment, construct a data vector using the spatial position of each reflection point at the current moment, and then construct the algebraic residual of each reflection point relative to the extended target at the current moment. S5: Perform single ellipse estimation based on the spatial position of each reflection point at the current moment to obtain the temporary estimate of the parameter vector and the kernel parameter estimate at the current moment; If the current time is the initial time, output the temporary estimate of the parameter vector at the current time; if the current time is not the initial time, use the kernel parameter estimates and algebraic residuals of the previous and current times to construct the bielliptic estimation objective function, with the constraint that the quadratic curve at each time is an ellipse, and that the lengths of the major and minor axes of the previous and current times are equal. S6: Transform the constraints and introduce adaptive weights into the bielliptic estimation objective function to obtain the joint optimization problem; S7: Decompose the joint optimization problem into a parameter vector subproblem and an adaptive weight subproblem; use alternating iteration to solve the problem and obtain the parameter vector estimates of the previous and current time steps; S8: Convert the current parameter vector estimate into the center position, major and minor axis lengths, and attitude angles of the extended target.

[0008] In a preferred embodiment, step S2, the position estimation of the non-cooperative transmitter includes: Calculate the coarse estimate of the position of the non-cooperative transmitter at the current time, and express the coarse estimate of the position of the non-cooperative transmitter at the k-th time as t. ini,k , ,in, and d k These represent the observed angle of arrival and propagation distance of the direct path at time k, respectively. The coarse position estimate is recursively averaged to obtain the fine position estimate of the non-cooperative transmitter at the current moment.

[0009] The position of the non-cooperative transmitter is obtained directly by using the observations of the angle of arrival and propagation distance of the direct path. The position of the non-cooperative transmitter is then refined step by step by running average recursion to suppress the impact of single observation noise on positioning accuracy.

[0010] A preferred embodiment is that, in step S3, the calculation of the spatial position of each reflection point includes: First-order approximations are made for the observed angle of arrival and propagation distance of each target reflection path at the current moment, forming a system of linear equations about the position of each reflection point at the current moment; The linear equations are solved using a coarse estimate of the position of the non-cooperative transmitter at the current moment to obtain the spatial position of each reflection point at the current moment.

[0011] By making a first-order approximation to the nonlinear observation equations for the angle of arrival and propagation distance, the geometric positioning problem is transformed into a system of linear equations. Then, by substituting the coarse estimate of the position of the non-cooperative transmitter into the closed-form solution, the Cartesian coordinates of each reflection point can be quickly calculated.

[0012] In a preferred embodiment, in step S3, the current time is set as the k-th time, and the spatial position of the i-th reflection point at the k-th time is represented by u. i,k u i,k =(A ini,i,k ) -1 b ini,i,k , where A ini,i,k Let b represent the approximate coefficient matrix consisting of the observed angle of arrival and propagation distance of the i-th target reflection path at time k, and the coarse estimate of the position of the non-cooperative transmitter at time k. ini,i,k This represents the approximate coefficient vector consisting of the observed propagation distance of the i-th target reflection path at time k and the coarse estimate of the position of the non-cooperative transmitter at time k.

[0013] In a preferred embodiment, in step S4, the algebraic residual of each reflection point relative to the extended target at the current moment is the product of the transpose of the parameter vector and the data vector.

[0014] The preferred embodiment is that, in step S5, the single ellipse estimation adopts an ellipse fitting method based on the variable center maximum correlation entropy criterion of the Laplace kernel function; the kernel parameter estimates include kernel center estimates and kernel bandwidth estimates; the objective function for the double ellipse estimation is: , where p j The parameter vector represents the quadratic curve of the expanded target at time j, where time k is the current time and time (k-1) is the previous time. Q j v represents the number of reflection points on the extended target at time j. i,j This represents a data vector constructed using the spatial location of the i-th reflection point at time j. and Let δ represent the kernel center estimate and kernel bandwidth estimate at time j, respectively. i,j (p j ) = p j T v i,j This represents the algebraic residual of the i-th reflection point at time j relative to the extended target; constraint one is an ellipticity criterion constraint: , where p j,(1) p j,(2) p j,(3) They represent p respectively j The first element, the second element, and the third element; Constraint condition two is the consistency constraint of the major and minor semi-axes: a k =a k-1 , b k =b k-1 , where a k and b k Let a represent the lengths of the major and minor axes of the expanded target at time k, respectively. k-1 and b k-1 Let represent the lengths of the major and minor semi-axis of the extended target at time k-1, respectively. Single-ellipse estimation employs the variable-center maximum correlation entropy criterion, using adaptive estimation of the kernel center to align the kernel function with the center of the residual distribution, effectively addressing residual bias caused by the unilateral distribution of reflection points. Double-ellipse estimation aims to maximize the sum of kernel similarities between two consecutive time points. Ellipticity discrimination constraints are introduced to ensure the quadratic curve is an ellipse, and semi-major and minor semi-axis consistency constraints are introduced to leverage the physical characteristic that the extended target's size remains essentially unchanged between adjacent time points, improving the temporal stability of shape estimation. No equality constraints are set for the center position and attitude angle to allow for translation and rotation of the extended target.

[0015] The preferred embodiment is that, in step S6, the transformation of the constraint conditions includes: Transform the ellipticity criterion constraint into a second-order cone constraint: ||[p j,(2) ,ε, p j,(1) -p j,(3)] T ||≤p j,(1) + p j,(3) Where ε is a preset non-zero constant; For the consistency constraint of the major and minor semi-axes, the first-order Taylor expansion of the major and minor semi-axes functions at the previous time step is performed using the parameter vector estimate at the previous time step, and the first-order Taylor expansion of the major and minor semi-axes functions at the current time step is performed using the temporary parameter vector estimate at the current time step, so as to obtain the linearized consistency constraint of the major and minor semi-axes. The objective function of the joint optimization problem is: , where w i,j ρ(w) represents the adaptive weight corresponding to the i-th reflection point at time j. i,j ) represents the convex conjugate function corresponding to the Laplace kernel function; The constraints of the joint optimization problem are the second-order cone constraint and the linearized major and minor semi-axes consistency constraint.

[0016] The ellipticity discrimination constraint is transformed into a second-order cone form, making the constraint convex. The consistency constraint of the major and minor semi-axis is linearized through a first-order Taylor expansion, and the original nonlinear equation is approximated by a linear equation at local expansion points, so that the constraint can be incorporated into the convex optimization framework. Based on the convex conjugate property of the kernel function, adaptive weights are introduced, transforming the original optimization problem into a joint optimization problem of parameter vectors and adaptive weights. This allows outliers that deviate significantly from the extended target contour to automatically receive low weights during the optimization process, thus achieving automatic suppression of the influence of outliers.

[0017] The preferred solution is that, in step S7, before the alternating iteration, the parameter vector subproblem is transformed into a convex optimization subproblem by introducing slack variables and aiming to minimize the sum of the weighted slack variables. The convex optimization subproblem is a second-order cone programming subproblem, and its constraints include: the absolute value of the difference between the algebraic residual of each reflection point relative to the extended target and the kernel center estimate is not greater than the corresponding slack variable, as well as second-order cone constraints and linearized major and minor semi-axes consistency constraints. The convex optimization subproblem and the adaptive weighting subproblem are solved iteratively, including: By fixing the adaptive weights of the current iteration, solve the second-order cone programming subproblem to obtain the parameter vector estimates of the previous and current time steps of the current iteration; Fix the parameter vector of the current iteration; according to The adaptive weight estimates for the previous and current iterations are calculated, where, This represents the adaptive weight estimate corresponding to the i-th reflection point at time j of the current iteration. This represents the parameter vector estimate at time j of the current iteration; Repeat the above process of solving the second-order cone programming subproblem with fixed adaptive weights and updating the adaptive weights with fixed parameter vectors until the termination condition is met.

[0018] By introducing slack variables, the parameter vector subproblem is transformed into a second-order cone programming subproblem. The objective is to minimize the sum of the weighted slack variables, thus transforming the original non-convex problem into an efficiently solvable convex optimization problem. The adaptive weight subproblem directly updates the weights using the closed-form solution of the convex conjugate function. The two subproblems are solved iteratively, with the adaptive weights fixed to solve the parameter vector and the parameter vector fixed to update the adaptive weights. This transforms the original nonlinear, non-convex joint optimization problem into a series of convex subproblems, balancing robustness and computational efficiency.

[0019] The preferred approach is that the objective function of the second-order cone programming subproblem is: ,in, Let represent the slack variable corresponding to the i-th reflection point at time j.

[0020] In a preferred embodiment, the termination condition for the alternating iterative solution in step S7 is: the change in the parameter vector estimate obtained from two adjacent iterations at the current time is less than a preset threshold, or the number of iterations reaches the preset maximum number of iterations.

[0021] Compared with the prior art, the advantages of the present invention are as follows: Based on the observed angle of arrival and propagation distance of the direct path, the position of the non-cooperative transmitter can be estimated without prior knowledge of the non-cooperative transmitter's position. This overcomes the limitation of traditional passive sensing methods that rely on the prior position of the non-cooperative transmitter and is suitable for passive monitoring of extended targets such as low-altitude aircraft, vehicles, and ships under unknown non-cooperative transmitter conditions.

[0022] Based on the observed angle of arrival and propagation distance of the target reflection path, as well as the position of the non-cooperative transmitter, the spatial position of the reflection point is calculated, realizing the transformation from the observed angle of arrival and propagation distance to the Cartesian coordinates of the reflection point, providing a data foundation for extended target shape estimation.

[0023] By estimating kernel parameters (kernel center and kernel bandwidth), the kernel center can be adaptively adjusted according to the algebraic residual distribution, effectively addressing the measurement bias problem caused by the concentration of reflection points on one side of the extended target.

[0024] By introducing adaptive weights into the objective function of biellipse estimation, outliers that deviate significantly from the extended target contour are assigned smaller adaptive weights, thus automatically suppressing the influence of outliers without the need for pre-detection or removal of outliers.

[0025] By introducing a constraint that the lengths of the major and minor semi-axis are equal at adjacent time points, and taking advantage of the fact that the lengths of the major and minor semi-axis remain basically unchanged at adjacent time points, historical shape information is used for compensation when the number of measurement points is insufficient or poorly distributed, thereby improving tracking stability.

[0026] By decomposing the joint optimization problem into two subproblems and applying alternating iterative solutions, the original non-convex and nonlinear optimization problem is transformed into a convex optimization subproblem that can be solved efficiently, thus balancing robustness and computational efficiency. Attached Figure Description

[0027] Figure 1 This is a flowchart illustrating the overall implementation of the method of the present invention; Figure 2 This is a schematic diagram of a passive sensing extended target tracking scenario to which the method of the present invention is applicable; Figure 3 A schematic diagram of the measurement points and the true contour trajectory of the extended target in a Monte Carlo experiment with outliers. Figure 4 This diagram illustrates the measurement points and estimation results of the method of this invention in a Monte Carlo experiment with outliers. Figure 5 The graph shows the average Gauss-Wasestein distance of the method of the present invention as a function of sampling time in the presence of outliers. Figure 6 The graph shows the root mean square error of the method for estimating the position of a non-cooperative transmitter as a function of sampling time in the presence of outliers. Detailed Implementation

[0028] The present invention will be further described in detail below with reference to the accompanying drawings and embodiments.

[0029] This invention proposes a robust extended target tracking method based on temporal bielliptic constraints, applicable to passive sensing extended target tracking scenarios. These scenarios include a receiver with a known location, a non-cooperative transmitter with an unknown location that remains stationary during observation, and an extended target with multiple reflection points. Without loss of generality, the example is a two-dimensional scenario where the receiver is positioned at the origin, and the receiver, non-cooperative transmitter, and any reflection point are not collinear. However, this invention is not limited to this. For multi-receiver systems, the reflection point localization results obtained by each receiver can be fused before performing the same bielliptic estimation, which still falls within the scope of this invention.

[0030] See Figure 1 The robust extended target tracking method of the present invention includes the following steps: S1: The receiver receives the signals of the direct path and the reflection paths of each target, and obtains the observation values ​​of the angle of arrival and propagation distance of the direct path and the reflection paths of each target at the current time; wherein, each target reflection path corresponds to a reflection point on the extended target.

[0031] It should be noted that the physical quantity directly acquired by the receiver is the time of arrival observation value. Multiplying the time of arrival observation value by the signal propagation speed (usually the speed of light) can be used to calculate the propagation distance observation value.

[0032] Specifically, if we set the current time as the k-th time, then the observed angle of arrival and propagation distance of the direct path at the k-th time are respectively... and d k The observed angle of arrival and propagation distance of the i-th target reflection path at the k-th time are θ, respectively. i,k and r i,k , , , , Where k=1,2,……,i=1,2,……,Q k Q k This represents the number of reflection points on the extended target at time k. This represents the true angle of arrival of the direct path at time k. express The observation noise included Let n represent the true propagation distance of the direct path at time k. d,k d k The observation noise included This represents the true angle of arrival of the i-th target reflection path at time k. Represents θ i,k The observation noise included This represents the true propagation distance of the i-th target reflection path at time k. Indicates r i,k The observation noise included is so that all true values ​​are unknown. o This represents the actual location of the non-cooperative transmitter. The superscript T indicates the transpose operation. Let x and y represent the true values ​​of the x-axis and y-axis of the non-cooperative transmitter, respectively. This represents the true position of the i-th reflection point on the extended target at time k. , Let x and y represent the true values ​​of the x-coordinate and y-coordinate of the i-th reflection point on the extended target at the k-th time, respectively. Let ||·|| represent the L2 norm operation, and arctan2(·,·) be the arctangent function.

[0033] Figure 2 A schematic diagram of a passive sensing extended target tracking scenario to which the method of the present invention is applicable is shown, where 's' represents the receiver.

[0034] S2: Estimate the position of the non-cooperative transmitter at the current moment based on the observed angle of arrival and propagation distance of the direct path.

[0035] Further specifying, the current position estimate of the non-cooperative transmitter includes: Calculate the coarse estimate of the position of the non-cooperative transmitter at the current time. Let the current time be the k-th time, and express the coarse estimate of the position of the non-cooperative transmitter at the k-th time as t. ini,k , , where the initial value of k is 1, that is, the initial time is when k=1.

[0036] By performing a running average recursion on the coarse position estimate, we obtain the fine position estimate of the non-cooperative transmitter at the current time. The fine position estimate of the non-cooperative transmitter at the k-th time is denoted as t. fin,k , t fin,k =t fin,k-1 +(1 / k)(t ini,k -t fin,k-1 ), where t fin,k-1 This represents the precise position estimate of the non-cooperative transmitter at time k-1, where 1 / k represents the update coefficient of the running average recursion. When k=1, t fin,k = t fin,k-1 +(1 / k)( t ini,k -t fin,k-1 This can be simplified to t fin,k =t ini,k That is, the precise position estimate of the non-cooperative transmitter at the first moment is its coarse position estimate.

[0037] S3: Based on the observed angle of arrival and propagation distance of each target reflection path at the current moment, as well as the position of the non-cooperative transmitter, calculate the spatial position of each reflection point at the current moment.

[0038] Further specifying, the calculation of the spatial position of each reflection point at the current moment includes: First-order approximations are made to the expressions for the angle of arrival and propagation distance observations of each target reflection path at the current moment, forming a system of linear equations about the position of each reflection point at the current moment.

[0039] By using a coarse estimate of the position of the non-cooperative transmitter at the current moment, a system of linear equations is solved to obtain the spatial position of each reflection point at the current moment.

[0040] Specifically, if we set the current time as the k-th time, then the observed angle of arrival θ of the i-th target reflection path at the k-th time is... i,k and the observed propagation distance r i,k First-order approximations are performed separately to obtain corresponding approximate expressions. These two approximate expressions are then combined to form a system of linear equations concerning the position of the i-th reflection point at time k. The coarse estimate t of the non-cooperative transmitter's position at time k is then used. ini,k Substituting the equations into the linear system and solving them, we obtain the spatial location of the i-th reflection point at the k-th time.

[0041] Specifically, let the current time be the k-th time, and let the spatial position of the i-th reflection point at the k-th time be denoted as u. i,k u i,k =(A ini,i,k ) -1 b ini,i,k , where A ini,i,k This represents the angle of arrival θ of the i-th target reflection path observed at time k. i,k and the observed propagation distance r i,k And the coarse estimate t of the position of the non-cooperative transmitter at time k. ini,k The approximate coefficient matrix formed, b ini,i,k This represents the observed propagation distance r of the i-th target reflection path at time k. i,k And the coarse estimate t of the position of the non-cooperative transmitter at time k. ini,k The approximate coefficient vector formed , , t ini,k =[x t,ini,k ,y t,ini,k ] T x t,ini,k and y t,ini,k These represent the coarse estimates of the x-coordinate and y-coordinate of the non-cooperative transmitter at time k, respectively. In actual implementation, these are the observed angle of arrival θ of the i-th target reflection path at time k. i,k and the observed propagation distance r i,k and the true position t of the non-cooperative transmitter o Constituent coefficient matrix A i,k The propagation distance r of the i-th target reflection path at the k-th time is given by the observed value. i,k and the true position t of the non-cooperative transmitter o Constituent coefficient vector b i,k , using t ini,k Replace A i,k and b i,k t in o That is, we get A. ini,i,k and b ini,i,kSince the receiver, the non-cooperative transmitter, and any reflection point are not collinear, therefore A ini,i,k Reversible.

[0042] S4: Use a parameter vector to represent the quadratic curve of the extended target at the current moment, construct the corresponding data vector using the spatial position of each reflection point at the current moment, and thus construct the algebraic residual of each reflection point relative to the extended target at the current moment.

[0043] Specifically, let the current time be the k-th time, and let the parameter vector (six-dimensional vector) of the quadratic curve of the extended target at the k-th time be represented as p. k p k =[A k B k C k D k E k ,F k ] T , where A k B k C k D k E k ,F k The six coefficients of the quadratic curve of the extended target at time k, the center position of the extended target, the lengths of its major and minor axes, and the attitude angles can be recovered from the six-dimensional parameter vector of the quadratic curve; the data vector constructed using the spatial position of the i-th reflection point at time k is represented as v. i,k , , where x i,k and y i,k Let u represent the x-coordinate and y-coordinate of the i-th reflection point at the k-th time, respectively. i,k =[x i,k , y i,k ] T The data vector contains x-coordinates, y-coordinates, quadratic terms, cross terms, and constant terms; the algebraic residual of each reflection point relative to the extended target at the current time is the product of the transpose of the parameter vector and the data vector, and the algebraic residual of the i-th reflection point at the k-th time relative to the extended target is expressed as δ. i,k (p k ), δ i,k (p k )=p k T v i,k .

[0044] S5: Based on the spatial positions of each reflection point at the current moment, perform single ellipse estimation to obtain temporary estimates of the parameter vector at the current moment, as well as estimates of the kernel parameters at the current moment. Here, the single ellipse estimation adopts an ellipse fitting method based on the variable center maximum correlation entropy criterion of the Laplace kernel function; the kernel parameter estimates include kernel center estimates and kernel bandwidth estimates. The kernel center estimates are used to reflect the center of the algebraic residual distribution, and the kernel bandwidth estimates are used to reflect the dispersion of the algebraic residual distribution.

[0045] If the current time is the initial time, the temporary estimate of the parameter vector at the current time is directly output; if the current time is not the initial time, i.e., starting from the second time, the kernel parameter estimates of the previous and current times and the algebraic residuals are used to construct a bielliptic estimation objective function for the parameter vector by maximizing the sum of the kernel similarities of two consecutive times. The constraints of the bielliptic estimation objective function are that the quadratic curve at each time is an ellipse, and the lengths of the major and minor axes of the previous and current times are equal. The center position and attitude angle of the extended target are not subject to the constraint of being equal at the previous and current times, so as to allow the extended target to translate and rotate, i.e., allow its center position and attitude angle to change with time.

[0046] In this embodiment, the ellipse fitting method based on the variable center maximum correlation entropy criterion of the Laplace kernel function aims to minimize the weighted sum of the algebraic residuals of each reflection point relative to the extended target at the current time, and uses the quadratic curve at the current time as the constraint of the ellipse to solve for the temporary estimate of the parameter vector at the current time and the estimate of the kernel parameters at the current time, including the kernel center and kernel bandwidth estimates.

[0047] In this embodiment, the objective function for estimating the double ellipse is: Where max is the maximum value operation, p j The parameter vector represents the quadratic curve of the expanded target at time j, where time k is the current time and time (k-1) is the previous time. Q j Let v represent the number of reflection points on the extended target at time j, where exp is an exponential function with base e, and v i,j This represents a data vector constructed using the spatial location of the i-th reflection point at time j. and Let δ represent the kernel center estimate and kernel bandwidth estimate at time j, respectively. i,j (p j ) = p j T v i,j This represents the algebraic residual of the i-th reflection point at time j relative to the extended target; constraint one is an ellipticity criterion constraint: , where p j,(1) p j,(2)p j,(3) They represent p respectively j The first element, the second element, and the third element; Constraint condition two is the consistency constraint of the major and minor semi-axes: a k =a k-1 , b k =b k-1 , where a k and b k Let a represent the lengths of the major and minor axes of the expanded target at time k, respectively. k-1 and b k-1 Let represent the lengths of the major and minor axes of the extended target at time k-1, respectively. Constraint 2 utilizes the physical property that the major and minor axes of the extended target remain basically unchanged between adjacent time points.

[0048] S6: Transform the constraints of the bielliptic estimation objective function and combine them with the bielliptic estimation objective function to form a nonlinear nonconvex optimization problem; through the properties of the kernel function, introduce adaptive weights into the objective function of the nonlinear nonconvex optimization problem, i.e., the bielliptic estimation objective function, according to the degree of deviation between the algebraic residual and the kernel center estimate, to obtain a joint optimization problem about the parameter vector and the adaptive weights.

[0049] In this embodiment, the transformation of constraints includes: Transform the ellipticity criterion constraint into a second-order cone constraint: ||[p j,(2) ,ε, p j,(1) -p j,(3) ] T ||≤p j,(1) + p j,(3) Where ε is a preset nonzero constant used to ensure strict ellipticity, this constraint is related to The ellipse discrimination condition corresponds to this.

[0050] For the consistency constraint of the major and minor semi-axes, the first-order Taylor expansion of the major and minor semi-axe functions at the previous time step is performed using the parameter vector estimate at the previous time step, and the first-order Taylor expansion of the major and minor semi-axe functions at the current time step is performed using the temporary parameter vector estimate at the current time step, resulting in the linearized consistency constraint of the major and minor semi-axes, described as follows: , Here, the major axis function and the minor axis function represent the mapping relationship between the parameter vector of the conic section and the corresponding ellipse's major and minor axes, respectively. Specifically, these relationships can be obtained through the coefficient correspondence between the general equation and the standard equation of the ellipse. This represents the parameter vector estimate at the previous time step, when k=2. , This represents the temporary estimate of the parameter vector at the previous time step. This represents the temporary estimate of the parameter vector at the current moment. and Let represent the gradients of the major and minor semi-axis functions at the current time step, respectively, with respect to the parameter vector at the current time step. and Let represent the gradients of the major and minor semi-axis functions at the previous time step, respectively, with respect to the parameter vector at the previous time step. The condition to the right of the vertical line represents the local expansion point used in the first-order Taylor expansion.

[0051] In this embodiment, the objective function of the joint optimization problem is: , where w i,j ρ(w) represents the adaptive weight corresponding to the i-th reflection point at time j. i,j ) represents the convex conjugate function corresponding to the Laplace kernel function.

[0052] In this embodiment, the constraints of the joint optimization problem are the second-order cone constraint and the linearized major and minor semi-axes consistency constraint.

[0053] S7: The joint optimization problem is still a nonlinear and nonconvex optimization problem, which is difficult to solve. In order to solve this problem, the joint optimization problem is decomposed into a parameter vector subproblem and an adaptive weight subproblem. The alternating iterative solution method is used to obtain the parameter vector estimates of the previous time step and the current time step.

[0054] In this embodiment, before the alternating iterative solution, the parameter vector subproblem is transformed into a convex optimization subproblem by introducing slack variables and aiming to minimize the sum of the weighted slack variables. The convex optimization subproblem is a second-order cone programming subproblem with the objective function: ,in, This represents the slack variable corresponding to the i-th reflection point at time j; its constraints include: the absolute value of the difference between the algebraic residual of each reflection point relative to the extended target and the kernel center estimate is not greater than the corresponding slack variable (expressed as: ), as well as second-order cone constraints and linearized major and minor semi-axes consistency constraints.

[0055] Further specifying, the convex optimization subproblem and the adaptive weighting subproblem are solved iteratively, including: By fixing the adaptive weights of the current iteration, solve the second-order cone programming subproblem to obtain the parameter vector estimates of the previous and current time steps of the current iteration.

[0056] Fix the parameter vector of the current iteration; according to The adaptive weight estimates for the previous and current iterations are calculated, where, This represents the adaptive weight estimate corresponding to the i-th reflection point at time j of the current iteration. This represents the parameter vector estimate at time j of the current iteration. As can be seen from the formula for calculating the adaptive weight estimate, the adaptive weights are based on the centered residual. The adaptive weighting mechanism automatically adjusts its size. When a reflection point deviates significantly from the extended target contour, the absolute value of its adaptive weight automatically decreases, thus reducing its contribution to the parameter vector estimation of the extended target. Conversely, reflection points closer to the extended target contour have larger absolute values ​​of their adaptive weights, playing a dominant role in parameter vector estimation. This mechanism simultaneously performs outlier suppression during parameter vector estimation, eliminating the need for pre-detection or removal of outliers and avoiding the difficulties of threshold selection. Furthermore, the consistency constraint between the major and minor semi-axles utilizes historical shape information to constrain the current estimation, preventing outliers from causing drastic deviations in the estimation results. This, in conjunction with the adaptive weighting, ensures the robustness of the parameter vector estimation.

[0057] Specifically, in solving the second-order cone programming subproblem, the initial value of the adaptive weights corresponding to all reflection points at time j is -1 / Q. j (Adaptive weights are negative), meaning that the adaptive weights for all reflection points at time j in the first iteration are -1 / Q. j Starting from the second iteration, the fixed adaptive weights are the adaptive weight estimates of all reflection points at time j of the previous iteration, and the introduced slack variables are non-negative slack variables. In the process of solving the adaptive weight subproblem, the fixed parameter vector of the current iteration is the parameter vector estimate of the previous time and the current time obtained by solving the second-order cone programming subproblem.

[0058] The process of solving the second-order cone programming subproblem with fixed adaptive weights and updating the adaptive weights with fixed parameter vectors is repeated until a termination condition is met. The termination condition is: the change in the parameter vector estimate obtained between two consecutive iterations is less than a preset threshold, or the number of iterations reaches a preset maximum number of iterations. In this embodiment, the change in the parameter vector estimate obtained between two consecutive iterations refers to the L2 norm of the difference between the parameter vector estimate obtained in the current iteration and the parameter vector estimate obtained in the previous iteration. The preset threshold is preferably 10. -5 The maximum number of iterations is preferably set to 50.

[0059] S8: The parameter vector estimate at the current moment is converted into the center position, semi-major and semi-minor axis lengths and attitude angles of the extended target, and then the center trajectory and shape sequence of the extended target in time order can be obtained.

[0060] To verify the feasibility and stability of the method of the present invention under conditions containing outliers, simulation experiments were conducted.

[0061] The receiver is positioned at the origin of a two-dimensional coordinate system; the actual position of the non-cooperative transmitter is [30m, 10m].T The extended target was modeled as an ellipse with a major semi-axis length of 6m and a minor semi-axis length of 3m, and its attitude angle changed with the motion process. The total observation time was 32s, the sampling interval was 0.2s, and the number of reflection points on the extended target at each moment followed a Poisson distribution with a mean of 20.

[0062] At each moment, 70% of the reflection points were set on the side of the extended target facing the non-cooperative transmitter, 20% were set on the side facing away from the non-cooperative transmitter, and the remaining 10% were set as outliers. Outliers were generated by superimposing pulse errors uniformly distributed within an 8m × 8m square area centered on the center of the extended target onto the selected reflection points. Unless otherwise specified, all statistical curves were obtained by averaging 1000 Monte Carlo experiments.

[0063] Figure 3 The measurement point (i.e., reflection point) and the extended true contour trajectory of the target are given in a typical Monte Carlo experiment. Figure 4 The trajectory of the measurement point (i.e., the reflection point) and the estimation result of the method of this invention in a typical Monte Carlo experiment is given. From Figure 3 and Figure 4 As can be seen, although some measurement points deviate significantly from the extended target boundary, the adaptive weights effectively reduce the impact of outliers on ellipse fitting, and the consistency constraint of major and minor semi-axes compensates for the lack of information caused by unilateral sparse measurements by utilizing historical shape information, and the estimated trajectory maintains a high degree of consistency with the true trajectory.

[0064] The mean Gauss-Wasestein distance is used to comprehensively measure the estimation errors of the center position, semi-major and semi-minor axis lengths, and attitude angles of the extended target. Figure 5 The curves showing the variation of the average Gauss-Wasestein distance with sampling time are presented. It can be seen that the distance remains stable throughout the tracking process and does not diverge significantly due to outliers. This indicates that the method of the present invention can continuously obtain reliable estimates of the center position, semi-major axis length, and attitude angle of the extended target.

[0065] Figure 6 The root mean square error (RMSE) of the non-cooperative transmitter position estimation is presented as a function of sampling time. It can be seen that the RMS error decreases rapidly and stabilizes as direct path observations accumulate. Since the non-cooperative transmitter position is mainly estimated based on direct path observations, outliers in the reflection point do not directly affect the estimation process, thus maintaining good positioning accuracy.

[0066] The simulation parameters described above are only used to illustrate one embodiment of the present invention and do not constitute a limitation on the extended target size, sampling period, number of reflection points, outlier ratio, and application scenario.

Claims

1. A robust extended target tracking method based on temporal double elliptic constraints, characterized in that, Includes the following steps: S1: The receiver acquires the current observation values ​​of the angle of arrival and propagation distance of the direct path and the reflection paths of each target; S2: Estimate the position of the non-cooperative transmitter at the current moment based on the observed angle of arrival and propagation distance of the direct path; S3: Based on the observed angle of arrival and propagation distance of each target reflection path at the current moment, and the position of the non-cooperative transmitter at the current moment, calculate the spatial position of each reflection point at the current moment; S4: Use a parameter vector to represent the quadratic curve of the extended target at the current moment, construct a data vector using the spatial position of each reflection point at the current moment, and then construct the algebraic residual of each reflection point relative to the extended target at the current moment. S5: Perform single ellipse estimation based on the spatial position of each reflection point at the current moment to obtain the temporary estimate of the parameter vector and the kernel parameter estimate at the current moment; If the current time is the initial time, output the temporary estimate of the parameter vector at the current time; If the current time is not the initial time, the kernel parameter estimates and algebraic residuals of the previous and current times are used to construct the double ellipse estimation objective function, with the constraints being that the quadratic curve at each time is an ellipse, and the lengths of the major and minor axes of the previous and current times are equal. S6: Transform the constraints and introduce adaptive weights into the bielliptic estimation objective function to obtain the joint optimization problem; S7: Decompose the joint optimization problem into a parameter vector subproblem and an adaptive weight subproblem; use alternating iteration to solve the problem and obtain the parameter vector estimates of the previous and current time steps; S8: Convert the current parameter vector estimate into the center position, major and minor axis lengths, and attitude angles of the extended target.

2. The robust extended target tracking method based on temporal double elliptic constraints according to claim 1, characterized in that, In step S2, the position estimation of the non-cooperative transmitter includes: Calculate the coarse estimate of the position of the non-cooperative transmitter at the current time, and express the coarse estimate of the position of the non-cooperative transmitter at the k-th time as t. ini,k , ,in, and d k These represent the observed angle of arrival and propagation distance of the direct path at time k, respectively. The coarse position estimate is recursively averaged to obtain the fine position estimate of the non-cooperative transmitter at the current moment.

3. The robust extended target tracking method based on temporal double elliptic constraints according to claim 2, characterized in that, In step S3, the calculation of the spatial position of each reflection point includes: First-order approximations are made for the observed angle of arrival and propagation distance of each target reflection path at the current moment, forming a system of linear equations about the position of each reflection point at the current moment; The linear equations are solved using a coarse estimate of the position of the non-cooperative transmitter at the current moment to obtain the spatial position of each reflection point at the current moment.

4. The robust extended target tracking method based on temporal double elliptic constraints according to claim 3, characterized in that, In step S3, the current time is set as the k-th time, and the spatial position of the i-th reflection point at the k-th time is represented by u. i,k u i,k =(A ini,i,k ) -1 b ini,i,k , where A ini,i,k Let b represent the approximate coefficient matrix consisting of the observed angle of arrival and propagation distance of the i-th target reflection path at time k, and the coarse estimate of the position of the non-cooperative transmitter at time k. ini,i,k This represents the approximate coefficient vector consisting of the observed propagation distance of the i-th target reflection path at time k and the coarse estimate of the position of the non-cooperative transmitter at time k.

5. The robust extended target tracking method based on temporal double elliptic constraints according to claim 1, characterized in that, In step S4, the algebraic residual of each reflection point relative to the extended target at the current moment is the product of the transpose of the parameter vector and the data vector.

6. A robust extended target tracking method based on temporal double elliptic constraints according to claim 1 or 5, characterized in that, In step S5, the single ellipse estimation adopts an ellipse fitting method based on the variable center maximum correlation entropy criterion of the Laplace kernel function. The kernel parameter estimates include the kernel center estimate and the kernel bandwidth estimate; the objective function for the bielliptic estimation is: , where p j The parameter vector represents the quadratic curve of the expanded target at time j, where time k is the current time and time (k-1) is the previous time. Q j v represents the number of reflection points on the extended target at time j. i,j This represents a data vector constructed using the spatial location of the i-th reflection point at time j. and Let δ represent the kernel center estimate and kernel bandwidth estimate at time j, respectively. i,j (p j ) = p j T v i,j This represents the algebraic residual of the i-th reflection point at time j relative to the extended target; constraint one is an ellipticity criterion constraint: , where p j,(1) p j,(2) p j,(3) They represent p respectively j The first element, the second element, and the third element; Constraint condition two is the consistency constraint of the major and minor semi-axes: a k =a k-1 , b k =b k-1 , where a k and b k Let a represent the lengths of the major and minor axes of the expanded target at time k, respectively. k-1 and b k-1 These represent the lengths of the major and minor axes of the extended target at time k-1, respectively.

7. A robust extended target tracking method based on temporal double elliptic constraints according to claim 6, characterized in that, In step S6, the transformation of the constraint conditions includes: Transform the ellipticity criterion constraint into a second-order cone constraint: ||[p j,(2) ,ε, p j,(1) -p j,(3) ] T ||≤p j,(1) + p j,(3) Where ε is a preset non-zero constant; For the consistency constraint of the major and minor semi-axes, the first-order Taylor expansion of the major and minor semi-axes functions at the previous time step is performed using the parameter vector estimate at the previous time step, and the first-order Taylor expansion of the major and minor semi-axes functions at the current time step is performed using the temporary parameter vector estimate at the current time step, so as to obtain the linearized consistency constraint of the major and minor semi-axes. The objective function of the joint optimization problem is: , where w i,j ρ(w) represents the adaptive weight corresponding to the i-th reflection point at time j. i,j ) represents the convex conjugate function corresponding to the Laplace kernel function; The constraints of the joint optimization problem are the second-order cone constraint and the linearized major and minor semi-axes consistency constraint.

8. A robust extended target tracking method based on temporal double elliptic constraints according to claim 7, characterized in that, In step S7, before the alternating iteration solution, the parameter vector subproblem is transformed into a convex optimization subproblem by introducing slack variables and aiming to minimize the sum of weighted slack variables. The convex optimization subproblem is a second-order cone programming subproblem, and its constraints include: the absolute value of the difference between the algebraic residual of each reflection point relative to the extended target and the kernel center estimate is not greater than the corresponding slack variable, as well as second-order cone constraints and linearized major and minor semi-axes consistency constraints. The convex optimization subproblem and the adaptive weighting subproblem are solved iteratively, including: By fixing the adaptive weights of the current iteration, solve the second-order cone programming subproblem to obtain the parameter vector estimates of the previous and current time steps of the current iteration; Fix the parameter vector of the current iteration; according to The adaptive weight estimates for the previous and current iterations are calculated, where, This represents the adaptive weight estimate corresponding to the i-th reflection point at time j of the current iteration. This represents the parameter vector estimate at time j of the current iteration; Repeat the above process of solving the second-order cone programming subproblem with fixed adaptive weights and updating the adaptive weights with fixed parameter vectors until the termination condition is met.

9. A robust extended target tracking method based on temporal double elliptic constraints according to claim 8, characterized in that, The objective function of the second-order cone programming subproblem is: ,in, Let represent the slack variable corresponding to the i-th reflection point at time j.

10. A robust extended target tracking method based on temporal double elliptic constraints according to claim 1, characterized in that, In step S7, the termination condition for alternating iterative solution is: the change in the parameter vector estimate obtained in two adjacent iterations at the current time is less than a preset threshold, or the number of iterations reaches the preset maximum number of iterations.

Images