Positioning method for uniformly moving targets with asynchronous transmitters and unknown positions

Abstract

The present invention discloses a method for positioning a uniformly moving target under conditions where transmitters are asynchronous and positions are unknown. Aiming at multi-base positioning scenarios, the method utilizes time delay information at consecutive moments to replace a traditional hybrid method of time delay and Doppler information, thereby reducing system hardware requirements. Aiming at a distance measurement model, the method utilizes numerical characteristics between measurement values ​​to eliminate some coupling variables, then separates unknown variables, and transforms the positioning problem into a constrained weighted least squares problem. By performing a first-order Taylor expansion on the objective function of the problem, a constrained weighted least squares problem with bias reduction capability is obtained. The problem is then transformed into a convex semi-definite programming problem. Finally, optimal estimates of the target's initial coordinates and moving speed are obtained. The method has the advantages of utilizing distance measurement values ​​observed at consecutive time points to estimate the target's position and speed when transmitters are asynchronous and positions are unknown, and has low equipment cost, low computational complexity, and a significant bias reduction capability.

Description

Technical Field

[0001] The present invention relates to a target positioning technology, in particular to a method for positioning a uniformly moving target with an asynchronous transmitter and an unknown position. Background Art

[0002] In the field of target positioning, locating moving targets is more challenging than locating stationary targets, especially in applications such as target tracking, sensor networks, and search and rescue. Traditional positioning techniques primarily rely on measurements such as time delay (TD), angle of arrival (AOA), and received signal strength (RSS). To more accurately estimate the velocity of moving targets, modern positioning systems must also utilize Doppler frequency shift (DFS) information.

[0003] In recent years, mobile target localization methods that combine time and frequency information have been widely studied and developed. These methods mainly include the combined application of Time Difference of Arrival (TDOA) and Frequency Difference of Arrival (FDOA).

[0004] Based on these studies, multi-base station positioning models suitable for mobile targets have been proposed. These models typically utilize passive coherent positioning technology to estimate the position and velocity of a moving target by analyzing the time delay and Doppler shift information between signals transmitted from base stations such as TV and radio towers and the echo signals reflected from the target. Because transmitters such as TV and radio towers are not specifically designed for positioning, the system is generally unable to directly determine their exact location. However, the timing function of signal transmission can be used to determine the corresponding signal transmission time, thereby indirectly obtaining time delay information. The main advantage of multi-base station positioning technology is that it can utilize existing infrastructure, eliminating the need for additional dedicated transmitting and receiving equipment, thereby reducing system costs.

[0005] Multi-base positioning technology not only excels in common applications but also holds great potential in scenarios requiring high-precision, wide-coverage positioning. Multi-base positioning systems utilize multiple transmitters and receivers, leveraging various measurement methods such as time difference of arrival (TDOA) and Doppler shift (DFS), to accurately locate moving targets.

[0006] Although the inclusion of Doppler frequency shift information can significantly improve the reliability of position and velocity estimates for moving targets, this approach is primarily applicable to narrowband signals. For wideband signals, the accuracy of frequency information is often insufficient for precise positioning. To overcome this limitation, some studies have proposed models that assume a moving target moves at a constant speed over a short period of time. By collecting corresponding distance or angle measurements at multiple consecutive time points, the initial position and velocity of the moving target are determined, thereby delineating the target's trajectory. However, the inclusion of multiple measurements increases the complexity of existing positioning methods.

[0007] While there has been some research on multi-base mobile target positioning, there has been no research on estimating the position and velocity of a mobile target using distance measurements from consecutive time points based on a multi-base positioning model. This new positioning problem poses a challenge to existing positioning technologies, but also provides an important opportunity for the development of new adaptive positioning algorithms, and is of great research significance. Summary of the Invention

[0008] The technical problem to be solved by the present invention is to provide a method for locating a uniformly moving target when the transmitter is asynchronous and the position is unknown. When the transmitter is asynchronous and the position is unknown, the method uses the distance measurement values ​​observed at consecutive time points to estimate the target position and speed, solving the problem that the target speed estimation cannot be achieved due to equipment cost limitations. The method has low computational complexity and the ability to reduce deviations.

[0009] The technical solution adopted by the present invention to solve the above technical problems is: a method for locating a uniformly moving target with an asynchronous transmitter and an unknown position, characterized by comprising the following steps:

[0010] Step 1: In a multi-base positioning scenario, assume there are several sensors and a target. Some sensors act as receivers with known coordinates and positions. Different receivers share the same reference clock, meaning they are synchronized. The remaining sensors act as transmitters with unknown coordinates and positions. The transmitters and receivers are asynchronous, meaning they have different reference clocks and a clock synchronization error exists between them. The target's initial coordinates are unknown and it moves at an unknown constant speed.

[0011] Step 2: At different time points, the signal transmitted by the transmitter reaches the receiver via a direct path or an indirect path reflected by the target. The receiver can obtain the direct path distance measurement values ​​and the indirect path distance measurement values ​​observed at the consecutive time points. Based on these values, a distance measurement model is constructed, including the direct path distance measurement part and the indirect path distance measurement part. The signal transmitted by the transmitter carries the clock information at the time of transmission.

[0012] Step 3: For the indirect path distance measurement portion of the distance measurement model, the indirect path distance measurement values ​​obtained by the remaining receivers except the first receiver are subtracted from the indirect path distance measurement value obtained by the first receiver, and the unknown variables are separated through mathematical processing. For the direct path distance measurement portion of the distance measurement model, the unknown variables are separated through mathematical processing. Then, using all the unknown variables, the distance measurement model after the unknown variables are separated is transformed into the original optimization problem, which is a weighted least squares problem with constraints.

[0013] Step 4: Perform a first-order Taylor expansion on the objective function in the original optimization problem; then, obtain a constraint on the bias term based on the first-order Taylor expansion; then, combine the original optimization problem and the constraint on the bias term to obtain a constrained weighted least squares problem with bias reduction capability;

[0014] Step 5: Use the semi-definite relaxation technique to relax the constrained weighted least squares problem with bias reduction capability into a convex semi-definite programming problem that is easy to solve.

[0015] Step 6: Directly solve the convex semi-positive programming problem to obtain the optimal estimate of the target's initial coordinate position, the optimal estimate of the target's moving speed, the optimal estimate of the coordinate position of each transmitter, and the optimal estimate of the clock synchronization error of each transmitter.

[0016] In step 1, the multi-base positioning scenario is two-dimensional or three-dimensional, the number of receivers is N, and the coordinate position of the nth receiver is s n , n=1,2,…,N, the number of transmitters is M, and the coordinate position of the mth transmitter is The clock synchronization error of the mth transmitter is The initial coordinate position of the target is p o , the target's moving speed is v o .

[0017] In step 2, the distance measurement model constructed is described as: Among them, K represents the total number of time points, r mn,k The indirect path distance measurement value obtained when the signal transmitted by the mth transmitter at time k reaches the nth receiver through the indirect path reflected by the target, d mn,k represents the direct path distance measurement value obtained by the signal transmitted by the mth transmitter at time k to the nth receiver through the direct path, ε rmn,k Represents r mn,k The corresponding measurement noise, ε dmn,k Indicates d mn,k The corresponding measurement noise, p o +kvo represents the coordinate position of the target at time k, and “||||” is the Euclidean norm operator.

[0018] In step 3, the specific process of separating the unknown variables in the indirect path distance measurement part of the distance measurement model is as follows:

[0019] Step 3.1a: Add the indirect path distance measurement part of the distance measurement model

[0020]

[0021] Split into two formulas, namely

[0022] and

[0023]

[0024] Step 3.2a: Subtract the two equations obtained in step 3.1a to obtain

[0025]

[0026] Step 3.3a: For r in step 3.1a m1,k In the expression Perform mathematical processing of transposing the terms, and then perform mathematical processing of square expansion on the result after transposing the terms, ignoring the second-order noise term, and obtain: The unknown variables are separated, where The superscript "T" indicates the transpose operation of a matrix or vector; for r in step 3.2a mn,k -r m1,k ||p in the expression o +kv o -s1|| performs mathematical processing of transposing the term, and then performs mathematical processing of square expansion on the result after transposing the term, and we get: The unknown variables are separated, where

[0027] In step 3, the specific process of separating the unknown variables in the direct path distance measurement part of the distance measurement model is as follows:

[0028] Direct path distance measurement part of the distance measurement model in Perform mathematical processing of transposing the terms, and then perform mathematical processing of square expansion on the result after transposing the terms, ignoring the second-order noise term, and obtain: The unknown variables are separated, where

[0029] In step 3, the process of obtaining the original optimization problem is as follows:

[0030] Step 3.1b: Combine Equations 1, 2, and 3 to obtain the distance measurement model after separating the unknown variables, which can be described as: Then all unknown variables are expressed as vectors, denoted as x o , in, All are coupled variables;

[0031] Step 3.2b: Exploit x o The distance measurement model after separating the unknown variables is transformed into the original optimization problem. The original optimization problem is a weighted least squares problem with constraints, which is described as: Among them, x is the optimization variable in the original optimization problem, x and x o Corresponding to, (g-Hx) T W(g-Hx) is the objective function in the original optimization problem, L represents the dimension of the multi-base positioning scenario, L is 2 or 3, and x is defined i Represents the i-th element of x, define x (i:j) represents the column vector consisting of the i-th element to the j-th element of x, B=blkdiag(B1,...,B K ),B k =blkdiag(B r1,k ,B r2,k ,B d ), Definition I i Represents the identity matrix of dimension i×i, Definition 1 i Represents a vector of all 1s of length i, and defines 0 i×j Represents an all-zero matrix of dimension i×j, Diag() represents the element diagonalization operation function, that is, the vectors are arranged diagonally to form a matrix, blkdiag() represents the matrix diagonalization operation function, that is, each matrix is ​​arranged diagonally to form a matrix, Q k Represents ε k The corresponding covariance matrix, ε k represents the composite noise vector at time k, ε k Subject to mean 0 and covariance matrix Q k Gaussian distribution, r n,k =[r 1n,k ,...,r Mn,k ] T , d n,k =[d 1n,k ,...,d Mn,k ] T , Represented by the matrix The submatrix consisting of the M+1th row to the MNth row and all columns, is the Kronecker product operator.

[0032] The specific process of step 4 is as follows:

[0033] Step 4.1: Define optimization variables in, α is an intermediate variable; then the objective function (g-Hx) in the original optimization problem is T W(g-Hx) is transformed into y T R T WRy, where R represents the coefficient matrix carrying noise, R = [-H, g]; then the objective function y T R T WRy performs a first-order Taylor expansion to obtain a first-order Taylor expansion, which is described as: in, Indicates the mathematical expectation, R o Represents the coefficient matrix without noise, r mn,k The corresponding indirect path distance measurement without noise d mn,k The corresponding direct path distance measurement without noise Replace r in R=[-H,g] mn,k d mn,k Get R o , “⊙” is the Hadamard product operator, ΔR=[-ΔH,Δg], A k,r1 =[0 M×(M+2)L ,I M ,0 M×(k-1) ,1 M ,0 M×(K-k) ,0 M×M(K+1) ], A k,r2 =[0 M(N-1)×((M+2)L+(M+k-1)) ,-1 M(N-1) ,0 M(N-1)×(K-k) ,0 M(N-1)×M(K+1) ], Represents the composite noise vector after deformation at K moments, represents the composite noise vector after deformation at time k,

[0034] Step 4.2: According to the first-order Taylor expansion, and introduce Form a constraint on the deviation term, described as: in,

[0035] Step 4.3: Combining the original optimization problem and the constraints on the bias term, we obtain a constrained weighted least squares problem with bias reduction capability, which can be described as: Among them, y end Represents the last element of y, define y (i:j) Represents the column vector consisting of the i-th element to the j-th element of y, and defines y i Represents the i-th element of y.

[0036] In step 5, the convex semi-positive programming problem is described as: Where Y represents the newly introduced optimization variable, and Y=yy T , tr{R T WRY} is the objective function in the convex semi-positive programming problem, tr{·} represents taking the matrix diagonal elements and summing them, Y is the optimization variable of the convex semi-positive programming problem, and is defined as Indicates that from row i to row j, column l to column j of Y The submatrix consisting of the elements corresponding to the columns is defined as Y (i,j) Represents the element corresponding to the i-th row and j-th column of Y, Y (end,end) is the index of the last row or column of Y. Y ≥ 0 means that the matrix Y is positive semidefinite.

[0037] In step 6, the initial coordinate position p of the target is o The optimal estimate of * , Set the target's moving speed v o The optimal estimate of * , The coordinate position of the mth transmitter The best estimate of The clock synchronization error of the mth transmitter is The best estimate of in, Y * represents the optimal estimate of Y, and defines Represents the submatrix consisting of the elements corresponding to the i-th row to the j-th row and the last column of Y, and defines Represents the element corresponding to the i-th row and j-th column of Y.

[0038] Compared with the prior art, the advantages of the present invention are:

[0039] 1) The method of the present invention considers a completely new positioning problem in which distance measurements at consecutive time points are used instead of the traditional method of using mixed distance and Doppler measurements to simultaneously estimate the target position and moving speed. This avoids the inability to estimate the moving target speed due to performance limitations of the signal acquisition equipment, and at the same time reduces the requirements for signal parameters compared to the mixed distance and Doppler method.

[0040] 2) The method of the present invention utilizes the numerical relationship between distance measurement values ​​to eliminate some coupled variables, thereby constructing a weighted least squares problem with constraints and a smaller number of coupled unknown variables; the reduction in the dimension of unknown variables can effectively reduce the search dimension of the problem, thereby reducing the computational complexity of the problem and facilitating the tightening of the optimization problem.

[0041] 3) The MSE analysis performed by the present invention shows that the constrained weighted least squares problem with bias reduction capability before relaxation can estimate unknown variables close to the performance of CRLB under Gaussian noise. This indicates that the method of the present invention can effectively estimate the target position and target velocity in this positioning scenario. At the same time, the elimination of coupling variables and reduction of estimation bias applied in forming the optimization problem optimize the results without reducing the accuracy of the estimation results. In addition, a theoretical bias expression is derived to evaluate the magnitude of the bias. BRIEF DESCRIPTION OF THE DRAWINGS

[0042] Figure 1 This is a block diagram of the overall implementation of the method of the present invention;

[0043] Figure 2 Schematic diagram of the system model of the present invention;

[0044] Figure 3 Comparison diagram of the logarithmic mean square error (MSE) of target position and target velocity estimation using the method of the present invention (SDP, BR-SDP) and the Cramer-Rao lower bound (CRLB) as the measurement noise power changes;

[0045] Figure 4 The figure is a comparison of the logarithmic bias (Bias) of the target position and target velocity estimation using the method of the present invention (SDP, BR-SDP) and the theoretical bias (Theoretical) as the measurement noise power changes. DETAILED DESCRIPTION

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

[0047] The present invention proposes a method for locating a uniformly moving target with an asynchronous transmitter and unknown position. The overall implementation block diagram is as follows: Figure 1 As shown, it includes the following steps:

[0048] Step 1: If Figure 2 As shown in the figure, in a multi-base positioning scenario, there are several sensors and one target. Some sensors act as receivers and their coordinate positions are known. Different receivers have the same reference clock, that is, the clocks of different receivers are synchronized. The remaining sensors act as transmitters and their coordinate positions are unknown. The time between the transmitter and the receiver is asynchronous, that is, the reference clocks are different, and there is a clock synchronization error between the transmitter and the receiver. The initial coordinate position of the target is unknown, and it moves at an unknown constant speed.

[0049] Specifically, in step 1, the multi-base positioning scenario is two-dimensional or three-dimensional, the number of receivers is N, and in this embodiment, N=6, and the coordinate position of the nth receiver is s n , n=1,2,…,N, the number of transmitters is M, in this embodiment, M=2, the coordinate position of the mth transmitter is The clock synchronization error of the mth transmitter is The initial coordinate position of the target is p o , the target's moving speed is v o .

[0050] Step 2: At different consecutive time points, the signal transmitted by the transmitter reaches the receiver via a direct path or an indirect path reflected by the target. The receiver can obtain the direct path distance measurement values ​​and indirect path distance measurement values ​​observed at consecutive time points. Based on this, a distance measurement model is constructed, including the direct path distance measurement part and the indirect path distance measurement part. Among them, the signal transmitted by the transmitter carries the clock information when the signal was transmitted.

[0051] It is worth noting that this positioning scenario is different from the existing positioning scenario. When the target is moving, it considers collecting only the time delay measurement value. The present invention considers a new positioning scenario for the problem of mobile target positioning. Specifically: in traditional positioning methods, if it is necessary to achieve simultaneous estimation of the position and speed of the moving target, it is usually necessary to use time delay information and Doppler information at the same time. The use of Doppler information requires that the signal acquisition equipment can achieve high-precision time-frequency measurement, and at the same time have the ability to sample and process signals at high speed. The method of the present invention innovatively uses only the time delay information (i.e., the distance information measurement value r mn,k , d mn,k), replacing the traditional method of mixing time delay and Doppler information. This significantly reduces the requirements and complexity of system hardware, reduces dependence on precision equipment, and improves the feasibility of the system.

[0052] Furthermore, the constructed distance measurement model is described as: Among them, K represents the total number of time points, r mn,k The indirect path distance measurement value obtained when the signal transmitted by the mth transmitter at time k reaches the nth receiver through the indirect path reflected by the target, d mn,k It represents the direct path distance measurement value obtained when the signal transmitted by the m-th transmitter at time k reaches the n-th receiver via the direct path. Specifically, after the signal transmitted by the m-th transmitter at time k reaches the n-th receiver via the direct path or the indirect path reflected by the target, the n-th receiver determines the signal flight time by comparing its own clock information at the time of signal arrival with the clock information carried by the signal, thereby obtaining the direct path distance measurement value or the indirect path distance measurement value. Represents r mn,k The corresponding measurement noise is, Indicates d mn,k The corresponding measurement noise, p o +kv o represents the coordinate position of the target at time k, "||||" is the Euclidean norm operator, used to represent the distance, ||p o +kv o -s n || represents the distance between target k and the nth receiver at time instant, represents the distance between target k and the mth transmitter at time instant, Indicates the distance between the mth transmitter and the nth receiver.

[0053] Step 3: For the indirect path distance measurement part of the distance measurement model, the indirect path distance measurement values ​​obtained by the remaining receivers except the first receiver are subtracted from the indirect path distance measurement value obtained by the first receiver, and the unknown variables are separated through mathematical processing; for the direct path distance measurement part of the distance measurement model, the unknown variables are separated through mathematical processing; then, using all the unknown variables, the distance measurement model after the separation of the unknown variables is transformed into the original optimization problem, which is a weighted least squares problem with constraints.

[0054] Specifically, in step 3, the specific process of separating the unknown variables in the indirect path distance measurement part of the distance measurement model is as follows:

[0055] Step 3.1a: Add the indirect path distance measurement part of the distance measurement model

[0056]

[0057] Split into two formulas, namely

[0058] and

[0059]

[0060] Step 3.2a: Subtract the two equations obtained in step 3.1a to obtain This difference step effectively reduces the number of coupled variables, achieving the goal of tightening the final optimization problem. The present invention cleverly utilizes the numerical relationship between the target reflection path distance information measurements between different receivers, and eliminates the common norm term before forming the optimization problem. This process has the following advantages: (1) By eliminating the common norm term, the number of coupled variables in the final optimization problem is effectively reduced, thereby reducing the complexity of the calculation. The reduction in the number of variables means that when solving the optimization problem, the calculation time and resource consumption can be shortened, and the solution efficiency can be improved. (2) Removing the common norm term can make the optimization problem more compact and the constraints more stringent. This has the effect of "tightening" the optimization problem to a certain extent, making the final solution more accurate and stable. Such a tightening effect helps to improve the convergence of the optimization problem, reduce the feasible domain of the solution, and thus enhance the robustness and accuracy of the algorithm.

[0061] Step 3.3a: For r in step 3.1a m1,k In the expression Perform mathematical processing of transposition, and then perform mathematical processing of square expansion on the result after transposition, ignoring the second-order noise term get: The unknown variables are separated, where The superscript "T" indicates the transpose operation of a matrix or vector, s1, r m1,k and Respectively represent s n 、r mn,k and The result when n=1; for r in step 3.2a mn,k -r m1,k ||p in the expression o +kv o -s1|| performs mathematical processing of transposing the term, and then performs mathematical processing of square expansion on the result after transposing the term, and we get: The unknown variables are separated, where

[0062] Specifically, in step 3, the specific process of separating the unknown variables in the direct path distance measurement part of the distance measurement model is as follows:

[0063] Direct path distance measurement part of the distance measurement model in Perform mathematical processing of transposing the terms, and then perform mathematical processing of square expansion on the result after transposing the terms, ignoring the second-order noise term, and obtain: The unknown variables are separated, where

[0064] Specifically, in step 3, the process of obtaining the original optimization problem is as follows:

[0065] Step 3.1b: Combine Equations 1, 2, and 3 to obtain the distance measurement model after separating the unknown variables, which can be described as: Then all unknown variables are expressed as vectors, denoted as x o , in, are coupled variables.

[0066] Step 3.2b: Exploit x o The distance measurement model after separating the unknown variables is transformed into the original optimization problem. The original optimization problem is a weighted least squares problem with constraints, which is described as: Among them, min means to minimize the objective function, st means "subject to the constraints...", considering x o The constraints generated by the numerical relationship between the unknown variables in the original optimization problem are x, x and x. o correspond, p represents the target position variable, v represents the target velocity variable, t1 represents the first transmitter position variable, t M represents the position variable of the Mth transmitter, τ1 represents the clock synchronization error variable of the first transmitter, τ M represents the clock synchronization error variable of the Mth transmitter, ξ 1,1 ,...,ξ 1,K , They represent the coupled variables formed by different combinations of the previous independent variables, (g-Hx) T W(g-Hx) is the objective function in the original optimization problem, L represents the dimension of the multi-base positioning scenario, L is 2 or 3, and x is defined i represents the i-th element of x, x ((M+2)L+M+k) Represents x iThe result when i=(M+2)L+M+k, x ((M+2)L+K+kM+m) Represents x i The result when i=(M+2)L+K+kM+m, x ((M+2)L+m) Represents x i The result when i=(M+2)L+m, x ((M+2)L+(M+1)K+M+m) Represents x i The result when i=(M+2)L+(M+1)K+M+m, define x (i:j) represents the column vector consisting of the i-th element to the j-th element of x, x (1:L) Represents x (i:j) The result when i=1, j=L, x (L+1：2L) Represents x (i:j) The result when i=L+1,j=2L,x ((m+1)L+1:(m+2)L) Represents x (i:j) The result when i=(m+1)L+1,j=(m+2)L is: g, g r1,k 、g r2,k and g d,k are introduced coefficient vectors, s2 represents s n The result when n=2, s N Indicates s n The result when n=N, r 11,k Represents r mn,k The result when n=1,m=1, r M1,k Represents r mn,k The result when n=1,m=M, r 12,k Represents r mn,k The result when n=2,m=1, r M2,k Represents r mn,k The result when n=2,m=M, r 1N,k Represents r mn,k The result when n=N,m=1, r MN,k Represents r mn,k The result when n=N,m=M,d 11,k Indicates d mn,k The result when n=1,m=1,d M1,k Indicates d mn,k The result when n=1,m=M,d 1N,k Indicates d mn,k The result when n=N,m=1, d MN,k Indicates d mn,k The result when n=N,m=M is: B=blkdiag(B1,...,B K ),Bk =blkdiag(B r1,k ,B r2,k ,B d ), B.B k 、B r1,k 、B r2,k 、B d are the introduced coefficient matrices, and I is defined as i Represents the identity matrix of dimension i×i, I MN Indicates I i The result when i=MN, I M Indicates I i The result when i=M, Definition 1 i represents a vector of all 1s of length i, 1 N-1 Indicates 1 i The result when i=N-1, define 0 i×j Represents a matrix of all zeros with dimension i×j, 0 M×MN Represents 0 i×j The result when i=M,j=MN is 0 M(N-1)×M(N-1) Represents 0 i×j The result when i=M(N-1), j=M(N-1), Diag() represents the element diagonalization operation function, that is, the vectors are arranged diagonally to form a matrix, blkdiag() represents the matrix diagonalization operation function, that is, each matrix is ​​arranged diagonally to form a matrix, represents the distance between the target and the first transmitter at time k, represents the distance between the target k and the Mth transmitter at time instant, Indicates the distance between the first transmitter and the first receiver, represents the distance between the Mth transmitter and the first receiver, Indicates the distance between the 1st transmitter and the Nth receiver, represents the distance between the Mth transmitter and the Nth receiver, represents the distance between the target and the second receiver at time k, represents the distance between target k and the Nth receiver at time k, Q k Represents ε k The corresponding covariance matrix, ε k represents the composite noise vector at time k, Respectively represent r 11,k 、r M1,k 、r 1N,k 、r MN,k d 11,k dM1,k d 1N,k d MN,k The corresponding measurement noise, ε k Subject to mean 0 and covariance matrix Q k Gaussian distribution, H, H k 、H r1,k 、H r1,k1 、H r1,k2 、H r1,k3 、H r2,k 、H r2,k1 、H r2,k2 、H d,k 、H d,k1 、H d,k2 、H d,k3 are the introduced coefficient matrices, 1 M Indicates 1 i The result when i=M is ​​1 N Indicates 1 i The result when i=N is 0 M×ML , 0 M×M , 0 M×(k-1) , 0 M×(K-k) , 0 M×(k-1)M , 0 M×(K-k)M , 0 M(N-1)×M(L+1) , 0 M(N-1)×M(K+1) , 0 M(N-1)×(k-1) , 0 M(N-1)×(K-k) , 0 MN×2L and 0 MN×K(M+1) Represents 0 respectively i×j When i=M,j=ML,i=M,j=M,i=M,j=k-1,i=M,j=Kk,i=M,j=M(k-1),i=M,j=M(Kk),i=M(N-1),j=M(L+ 1), the result when i=M(N-1),j=M(K+1),i=M(N-1),j=k-1,i=M(N-1),j=Kk,i=MN,j=2L,i=MN,j=K(M+1), r n,k =[r 1n,k ,...,r Mn,k ] T , r 1n,k Represents r mn,k The result when m=1, r Mn,k Represents r mn,k The result when m=M, d n,k =[d 1n,k ,...,d Mn,k ] T , d 1n,k Indicates d mn,k The result when m=1, dMn,k Indicates d mn,k The result when m=M is, Represented by the matrix The submatrix consisting of the M+1th row to the MNth row and all columns, is the Kronecker product operator.

[0067] Step 4: Perform a first-order Taylor expansion on the objective function in the original optimization problem; then, obtain a constraint on the bias term based on the first-order Taylor expansion; and then combine the original optimization problem and the constraint on the bias term to obtain a constrained weighted least squares problem with bias reduction capability.

[0068] Specifically, the specific process of step 4 is:

[0069] Step 4.1: Define optimization variables in, α is an intermediate variable; then the objective function (g-Hx) in the original optimization problem is T W(g-Hx) is transformed into y T R T WRy, where R represents the coefficient matrix carrying noise, R = [-H, g]; considering that the measurement values ​​in the coefficient matrix R carry measurement errors, the objective function y T R T WRy performs a first-order Taylor expansion to obtain a first-order Taylor expansion, which is described as: in, Indicates the mathematical expectation, R o Represents the coefficient matrix without noise, r mn,k The corresponding indirect path distance measurement without noise d mn,k The corresponding direct path distance measurement without noise Replace r in R=[-H,g] mn,k d mn,k Get R o , “⊙” is the Hadamard product operator, ΔR=[-ΔH,Δg], A k,r1 =[0 M×(M+2)L ,I M ,0 M×(k-1) ,1 M ,0 M×(K-k) ,0 M×M(K+1) ], A k,r2 =[0 M(N-1)×((M+2)L+(M+k-1)) ,-1 M(N-1) ,0 M(N-1)×(K-k) ,0 M(N-1)×M(K+1) ], Represents the composite noise vector after deformation at K moments, represents the composite noise vector after deformation at time k, express The result when k=1 is: express The result when k=K is, express The result when m=1 is: express The result when m=M is, express The result when m=1 is: express The result when m=M.

[0070] Step 4.2: According to the first-order Taylor expansion, and introduce Form a constraint on the deviation term, described as: in,

[0071] Step 4.3: Combining the original optimization problem and the constraints on the bias term, we obtain a constrained weighted least squares problem with bias reduction capability, which can be described as: Among them, y end Represents the last element of y, define y (i:j) Represents the column vector consisting of the i-th element to the j-th element of y, y (1:L) represents y (i:j) The result when i=1,j=L, y ((L+1):2L ) represents y (i:j) The result when i=L+1,j=2L, y (((m+1)L+1):(m+2)L ) represents y (i:j) The result when i=(m+1)L+1,j=(m+2)L is defined as y i Represents the i-th element of y, y ((M+2)L+M+k) represents y i The result when i=(M+2)L+M+k is y ((M+2)L+K+kM+m) represents y i The result when i=(M+2)L+K+kM+m, y ((M+2)L+m) represents y i The result when i=(M+2)L+m, y ((M+2)L+(M+1)K+M+m) represents y i The result when i=(M+2)L+(M+1)K+M+m.

[0072] Step 5: Use the semi-definite relaxation technique to relax the constrained weighted least squares problem with bias reduction capability into a convex semi-definite programming problem that is easy to solve.

[0073] Specifically, in step 5, the convex semi-positive programming problem is described as: Where Y represents the newly introduced optimization variable, and Y=yy T , tr{R T WRY} is the objective function in the convex semi-positive programming problem, tr{·} represents taking the matrix diagonal elements and summing them, Y is the optimization variable of the convex semi-positive programming problem, and is defined as Indicates that from row i to row j, column l to column j of Y The submatrix Y consists of the elements corresponding to the columns. (1:L,1:L) express When i=1, j=L, l=1, The result when Y ((L+1):2L,(L+1):2L) express When i=L+1,j=2L,l=L+1, The result when Y (1:L,(L+1):2L) express When i=1, j=L, l=L+1, The result when Y (1:L,end) express When i=1, j=L, The result when Y ((L+1):2L,end) express In i=L+1,j=2L, The result when Y (((m+1)L+1):(m+2)L,((m+1)L+1):(m+2)L )express In i=(m+1)L+1, j=(m+2)L, l=(m+1)L+1, The result when Y (1:L,((m+1)L+1):(m+2)L) express When i=1, j=L, l=(m+1)L+1, The result when Y ((L+1):2L,((m+1)L+1):(m+2)L) express In i=L+1,j=2L,l=(m+1)L+1, The result when Y (i,j) Represents the element corresponding to the i-th row and j-th column of Y, Y ((M+2)L+m,(M+2)L+m) Indicates Y (i,j) The result when i=(M+2)L+m,j=(M+2)L+m,Y ((M+2)L+M+k,(M+2)L+M+k) Indicates Y (i,j) The result when i=(M+2)L+M+k,j=(M+2)L+M+k is Y (end,end) Indicates Y (i,j) The result when i=j=end, Y ((M+2)L+K+kM+m,end) Indicates Y (i,j)The result when i=(M+2)L+K+kM+m,j=end, Y ((M+2)L+m,(M+2)L+M+k) Indicates Y (i,j) The result when i=(M+2)L+m,j=(M+2)L+M+k is Y ((M+2)L+(M+1)K+M+m,end) Indicates Y (i,j) The result when i=(M+2)L+(M+1)K+M+m,j=end, Y (end,end) is the index of the last row or column of Y. Y ≥ 0 means that the matrix Y is positive semidefinite.

[0074] Step 6: Use existing toolboxes such as the interior point method to directly solve the convex semi-positive programming problem to obtain the optimal estimate of the target's initial coordinate position, the optimal estimate of the target's moving speed, the optimal estimate of the coordinate position of each transmitter, and the optimal estimate of the clock synchronization error of each transmitter.

[0075] Specifically, in step 6, the initial coordinate position of the target p o The optimal estimate of * , Set the target's moving speed v o The optimal estimate of * , The coordinate position of the mth transmitter The best estimate of The clock synchronization error of the mth transmitter is The best estimate of in, Y * represents the optimal estimate of Y, and defines represents the submatrix consisting of the elements corresponding to the i-th row to the j-th row and the last column of Y, express The result when i=1,j=L is: express The result when i=L+1,j=2L is: express The result when i=(m+1)L+1,j=(m+2)L is defined as represents the element corresponding to the i-th row and j-th column of Y, express The result when i=(M+2)L+m,j=end, express The result when i=j=end.

[0076] In order to verify the feasibility and effectiveness of the method of the present invention, a simulation test was carried out on the method of the present invention.

[0077] Assume that there are two transmitters (i.e., M = 2) and six receivers (i.e., N = 6) in a three-dimensional space, i.e., a three-dimensional coordinate system (i.e., L = 3). Their true Cartesian coordinate positions in the three-dimensional space are randomly generated. The true Cartesian coordinate position (initial coordinate position) and true movement speed of the target are also randomly generated. The range is shown in Table 1:

[0078] Table 1 Scenario range

[0079]

[0080] At the same time, the scene and measurement values ​​are scaled by 10 times (to avoid calculation problems in the toolbox due to excessively large measurement values ​​during the problem solving process). The measurement noise is assumed to have a mean of 0 and a covariance matrix Q at time k. k =σ 2 I 2MN Gaussian noise, where σ 2 represents the variance of the measurement noise, I 2MN Represents the identity matrix of dimension 2MN×2MN.

[0081] Based on the above parameter settings, the simulation randomly generated 10 test scenarios. In each test scenario, the proposed method was tested for the logarithmic mean square error (MSE) and logarithmic bias (Bias) between the optimal estimate of the target's initial coordinate position and the optimal estimate of the target's movement speed and the true value. To verify the stability of the method, 1000 Monte Carlo runs were performed on each test scenario, and the final result was the average of 1000 Monte Carlo runs and 10 test scenarios. To verify the performance of the method, the Cramer-Rao bound (CRLB) and theoretical bias were introduced for comparison.

[0082] Figure 3 Comparison graphs showing the variation of the log mean square error (MSE) of target position and target velocity estimates with the proposed method (SDP, BR-SDP) and the Cramer-Rao lower bound (CRLB) as the measurement noise power changes are shown. The top graph shows the variation of the log mean square error (MSE) of target position estimation with the measurement noise power, and the bottom graph shows the variation of the log mean square error (MSE) of target velocity estimation with the measurement noise power. Figure 4 Comparison graphs showing how the logarithmic bias (Bias) of target position and target velocity estimates varies with measurement noise power using the method of the present invention (SDP, BR-SDP) and the theoretical bias (Theoretical) are used. The top graph is a schematic diagram showing how the logarithmic bias (Bias) of target position estimation varies with measurement noise power, and the bottom graph is a schematic diagram showing how the logarithmic bias (Bias) of target velocity estimation varies with measurement noise power. Figure 3 and Figure 4 In the above, SDP refers to the semi-definite relaxation method, that is, the method obtained by the present invention without the process of reducing the deviation ability, and BR-SDP refers to the semi-definite relaxation method with the ability to reduce the deviation. Figure 3 As can be seen, when the measured noise power is 10log 10 σ 2 When θ changes, the method of the present invention can approach the Cramer-Rao lower bound under the entire noise condition, indicating that the method of the present invention can accurately estimate the target position and target velocity using only the time delay information. At the same time, the elimination of the coupling variables does not affect the performance of the final estimation result and also reduces the computational complexity of the problem. Figure 4 The results show that the proposed method can effectively reduce the logarithmic deviation of the final estimation result by reducing the deviation, which has an advantage of about 15dB compared with the result without bias reduction, and can be close to the theoretical deviation in the case of medium and below medium noise.

[0083] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A method for locating a uniformly moving target with an asynchronous transmitter and unknown position, characterized in that The following steps are involved: Step 1: In a multi-base positioning scenario, assume there are several sensors and a target. Some sensors act as receivers with known coordinates and positions. Different receivers share the same reference clock, meaning they are synchronized. The remaining sensors act as transmitters with unknown coordinates and positions. The transmitters and receivers are asynchronous, meaning they have different reference clocks and a clock synchronization error exists between them. The target's initial coordinates are unknown and it moves at an unknown constant speed. Step 2: At different time points, the signal transmitted by the transmitter reaches the receiver via a direct path or an indirect path reflected by the target. The receiver can obtain the direct path distance measurement values ​​and the indirect path distance measurement values ​​observed at the consecutive time points. Based on these values, a distance measurement model is constructed, including the direct path distance measurement part and the indirect path distance measurement part. The signal transmitted by the transmitter carries the clock information at the time of transmission. Step 3: For the indirect path distance measurement portion of the distance measurement model, the indirect path distance measurement values ​​obtained by the remaining receivers except the first receiver are subtracted from the indirect path distance measurement value obtained by the first receiver, and the unknown variables are separated through mathematical processing. For the direct path distance measurement portion of the distance measurement model, the unknown variables are separated through mathematical processing. Then, using all the unknown variables, the distance measurement model after the unknown variables are separated is transformed into the original optimization problem, which is a weighted least squares problem with constraints. Step 4: Perform a first-order Taylor expansion on the objective function in the original optimization problem; then, obtain a constraint on the bias term based on the first-order Taylor expansion; then, combine the original optimization problem and the constraint on the bias term to obtain a constrained weighted least squares problem with bias reduction capability; Step 5: Use the semi-definite relaxation technique to relax the constrained weighted least squares problem with bias reduction capability into a convex semi-definite programming problem that is easy to solve. Step 6: Directly solve the convex semi-positive programming problem to obtain the optimal estimate of the target's initial coordinate position, the optimal estimate of the target's moving speed, the optimal estimate of the coordinate position of each transmitter, and the optimal estimate of the clock synchronization error of each transmitter.

2. The method for locating a uniformly moving target with an asynchronous transmitter and unknown position according to claim 1 is characterized in that In step 1, the multi-base positioning scenario is two-dimensional or three-dimensional, the number of receivers is N, and the coordinate position of the nth receiver is s n , n=1,2,…,N, the number of transmitters is M, and the coordinate position of the mth transmitter is The clock synchronization error of the mth transmitter is The initial coordinate position of the target is p o , the target's moving speed is v o .

3. The method for locating a uniformly moving target with an asynchronous transmitter and unknown position according to claim 2 is characterized in that In step 2, the distance measurement model constructed is described as: Among them, K represents the total number of time points, r mn,k The indirect path distance measurement value obtained when the signal transmitted by the mth transmitter at time k reaches the nth receiver through the indirect path reflected by the target, d mn,k represents the direct path distance measurement value obtained when the signal transmitted by the m-th transmitter at time k reaches the n-th receiver through the direct path, Represents r mn,k The corresponding measurement noise is, Indicates d mn,k The corresponding measurement noise, p o +kv o represents the coordinate position of the target at time k, and "||| |" is the Euclidean norm operator.

4. The method for locating a uniformly moving target with an asynchronous transmitter and unknown position according to claim 3 is characterized in that In step 3, the specific process of separating the unknown variables in the indirect path distance measurement part of the distance measurement model is as follows: Step 3.1a: Add the indirect path distance measurement part of the distance measurement model Split into two formulas, namely and Step 3.2a: Subtract the two equations obtained in step 3.1a to obtain Step 3.3a: For r in step 3.1a m1,k In the expression Perform mathematical processing of transposing the terms, and then perform mathematical processing of square expansion on the result after transposing the terms, ignoring the second-order noise term, and obtain: (Formula 1) realizes the separation of unknown variables, where The superscript "T" indicates the transpose operation of a matrix or vector; for r in step 3.2a mn,k -r m1,k ||p in the expression o +kv o -s1|| performs mathematical processing of transposing the term, and then performs mathematical processing of square expansion on the result after transposing the term, and we get: (Equation 2) realizes the separation of unknown variables, where In step 3, the specific process of separating the unknown variables in the direct path distance measurement part of the distance measurement model is as follows: Direct path distance measurement part of the distance measurement model in Perform mathematical processing of transposing the terms, and then perform mathematical processing of square expansion on the result after transposing the terms, ignoring the second-order noise term, to obtain: (Equation 3) realizes the separation of unknown variables, where 5. The method for locating a uniformly moving target with an asynchronous transmitter and unknown position according to claim 4 is characterized in that In step 3, the process of obtaining the original optimization problem is as follows: Step 3.1b: Combine Equations 1, 2, and 3 to obtain the distance measurement model after separating the unknown variables, which can be described as: Then all unknown variables are expressed as vectors, denoted as x o , in, All are coupled variables; Step 3.2b: Exploit x o The distance measurement model after separating the unknown variables is transformed into the original optimization problem. The original optimization problem is a weighted least squares problem with constraints, which is described as: Among them, x is the optimization variable in the original optimization problem, x and x o Corresponding to, (g-Hx) T W(g-Hx) is the objective function in the original optimization problem, L represents the dimension of the multi-base positioning scenario, L is 2 or 3, and x is defined i Represents the i-th element of x, define x (i:j) represents the column vector consisting of the i-th element to the j-th element of x, B=blkdiag(B1,...,B K ),B k =blkdiag(B r1,k ,B r2,k ,B d ), Definition I i Represents the identity matrix of dimension i×i, Definition 1 i Represents a vector of all 1s of length i, and defines 0 i×j Represents an all-zero matrix of dimension i×j, Diag( ) represents the element diagonalization operation function, that is, the vectors are arranged diagonally to form a matrix, blkdiag( ) represents the matrix diagonalization operation function, that is, each matrix is ​​arranged diagonally to form a matrix, Q k Represents ε k The corresponding covariance matrix, ε k represents the composite noise vector at time k, ε k Subject to mean 0 and covariance matrix Q k Gaussian distribution, r n,k =[r 1n,k ,...,r Mn,k ] T , d n,k =[d 1n,k ,...,d Mn,k ] T , Indicated by the moment Formation The submatrix consisting of the M+1th row to the MNth row and all columns, is the Kronecker product operator.

6. The method for locating a uniformly moving target with an asynchronous transmitter and unknown position according to claim 5 is characterized in that The specific process of step 4 is as follows: Step 4.1: Define optimization variables in, α is an intermediate variable; then the objective function (g-Hx) in the original optimization problem is T W(g-Hx) is transformed into y T R T WRy, where R represents the coefficient matrix carrying noise, R = [-H, g]; then the objective function y T R T WRy performs a first-order Taylor expansion to obtain a first-order Taylor expansion, which is described as: in, Indicates the mathematical expectation, R o Represents the coefficient matrix without noise, r mn,k The corresponding indirect path distance measurement without noise d mn,k The corresponding direct path distance measurement without noise Replace r in R=[-H,g] mn,k d mn,k Get R o , "⊙" is the Hadamard product operator, ΔR=[-ΔH,Δg], A k,r1 =[0 M×(M+2)L ,I M ,0 M×(k-1) ,1 M ,0 M×(K-k) ,0 M×M(K+1) ], A k,r2 ＝[0 M(N-1)×((M+2)L+(M+k-1)) ,-1 M(N-1) ,0 M(N-1)×(K-k) ,0 M(N-1)×M(K+1) ]， Represents the composite noise vector after deformation at K moments, represents the composite noise vector after deformation at time k, Step 4.2: According to the first-order Taylor expansion, and introduce Form a constraint on the deviation term, described as: in, Step 4.3: Combining the original optimization problem and the constraints on the bias term, we obtain a constrained weighted least squares problem with bias reduction capability, which can be described as: Among them, y end Represents the last element of y, define y (i:j) Represents the column vector consisting of the i-th element to the j-th element of y, and defines y i Represents the i-th element of y.

7. The method for locating a uniformly moving target with an asynchronous transmitter and unknown position according to claim 6 is characterized in that In step 5, the convex semi-positive programming problem is described as: AND ((M+2)L+K+kM+m,end) =-tr{Y (((m+1)L+1):(m+2)L,((m+1)L+1):(m+2)L) }+2tr{And (1:L,((m+1)L+1):(m+2)L) }+2ktr{And ((L+1):2L,((m+1)L+1):(m+2)L) }+Y ((M+2)L+m,(M+2)L+m) +2Y ((M+2)L+m,(M+2)L+M+k) AND ((M+2)L+(M+1)K+M+m,end) =And ((M+2)L+m,(M+2)L+m) -tr{And (((m+1)L+1):(m+2)L,((m+1)L+1):(m+2)L) } Where Y represents the newly introduced optimization variable, and Y=yy T , tr{R T WRY} is the objective function in the convex semi-positive programming problem, tr{·} represents taking the matrix diagonal elements and summing them, Y is the optimization variable of the convex semi-positive programming problem, and is defined as Indicates that from row i to row j, column l to column j of Y The submatrix consisting of the elements corresponding to the columns is defined as Y (i,j) Represents the element corresponding to the i-th row and j-th column of Y, Y (end,end) is the index of the last row or column of Y. Y ≥ 0 means that the matrix Y is positive semidefinite.

8. The method for locating a uniformly moving target with an asynchronous transmitter and unknown position according to claim 7 is characterized in that In step 6, the initial coordinate position p of the target is o The optimal estimate of * , Set the target's moving speed v o The optimal estimate of * , The coordinate position of the mth transmitter The best estimate of The clock synchronization error of the mth transmitter is The best estimate of in, Y * represents the optimal estimate of Y, and defines Represents the submatrix consisting of the elements corresponding to the i-th row to the j-th row and the last column of Y, and defines Represents the element corresponding to the i-th row and j-th column of Y.

Images