Generalized schur decomposition time delay estimation method for GNSS array tensor modeling

Abstract

The present invention discloses a generalized Schur decomposition delay estimation method for GNSS array tensor modeling. Specifically, the method comprises the following steps: first, calculating the HOSVD of the signal tensor, compressing the core tensor and singular matrix, performing a generalized Schur decomposition on the two frontal slices of the core tensor, and then upper triangulating the matrix of the frontal slice of the core tensor; using LSKRF to decompose the matrix and solve the factor matrix, and then using the factor matrix to find the corresponding Loss of Sight (LOS) signal components; then, calculating the signal amplitude based on the structure of the signal tensor; interpolating to obtain the post-correlation interpolation function of the correlator and the LOS signal, and estimating the LOS delay based on the interpolation function peak search. The present invention achieves high-precision and low-computational-complexity delay estimation of GNSS signals, with higher robustness and lower computational complexity.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of navigation signals, and particularly relates to a generalized Schur decomposition time delay estimation method for GNSS array tensor modeling. BACKGROUND

[0002] Global Navigation Satellite System (GNSS), such as GPS and Beidou of China, provides high-precision positioning services for the fields of aviation, autonomous driving, precision agriculture, etc., and has achieved wide application. At present, the error factors affecting the positioning accuracy include satellite orbit error, clock error, ionospheric delay, tropospheric delay, and multipath, etc. Among them, the satellite clock error, orbit error, ionospheric and tropospheric delay errors have space-time correlation, and can be eliminated by using differential positioning technology. However, the multipath error belongs to accidental error, and when the receiver antenna is located at different positions, the errors caused by multipath are generally not correlated with each other, and it is not practical to model the environment of each receiver antenna. Since the multipath effect does not have space-time correlation, it cannot be eliminated by using differential technology, which becomes the main bottleneck for high-precision users to improve the positioning accuracy and reliability. In order to eliminate the influence of multipath, the multipath signal suppression algorithm has become a hot spot in the field of satellite navigation. The ranging data quality provided by the GNSS receiver depends largely on the accuracy of the LOS signal propagation time delay estimation for each satellite. However, due to the influence of the propagation environment, the reflection of satellite signals on trees, lamps and buildings will produce multipath components, and at the receiver, the superposition of Line Of Sight (LOS) signals and None Line Of Sight (NLOS) signals will reduce the time delay estimation accuracy, thereby reducing the accuracy of position estimation. Therefore, a GNSS receiver equipped with an antenna array is considered, and a correlator is used to process to reduce the influence of multipath components. In the face of the multi-dimensional structure of the GNSS array receiving system, the time delay estimation algorithm based on tensor modeling shows obvious advantages: high accuracy and low computational complexity.

[0003] Currently, there are several tensor-based high-precision time delay estimation algorithms as follows: (1) FBA+ESPS+HOSVD algorithm. The algorithm adopts forward-backward averaging (FBA) and expanded spatial smoothing (ESPS) to preprocess the tensor data, and then estimates the output data of the correlator of the data tensor by using the HOSVD algorithm. Finally, the number of points is increased by interpolation, and one-dimensional peak search is performed to find the time delay corresponding to the maximum power. (2) DoA / KrF algorithm. The algorithm also uses FBA and ESPS to preprocess the tensor data, and then estimates the DoA of the array by using the estimation of signal parameter via rotational invariance technique (ESPRIT) algorithm to recover the third-dimensional factor matrix (steering matrix) of the tensor. Then, the remaining two-dimensional matrix factors are estimated by using the least squares Khatri-Rao factorization (LSKRF), and finally, the time delay is estimated by interpolation. (3) CPD-GEVD-based time delay estimation algorithm. The algorithm solves the first-dimensional matrix factor by using HOSVD and generalized eigenvalue decomposition (GEVD). Then, the remaining factor matrix is estimated by using LSKRF. Then, the LOS signal with the maximum power is selected for interpolation, and finally, one-dimensional peak search is performed to complete the time delay estimation. (4) SECSI+HOSVD-based time delay estimation algorithm. The algorithm is based on the structural characteristics of the signal tensor, and solves the first-dimensional matrix factor by using SECSI and HOSVD, and then separates the remaining factor matrix by using LSKRF. Then, the number of points is increased by interpolating the correlator output data, and finally, one-dimensional peak search is performed to estimate the signal time delay.

[0004] Among the current tensor-based time delay estimation algorithms, the FBA+ESPS+HOSVD and DoA / KrF algorithms are easily affected by array position defects, the CPD-GEVD algorithm has poor performance when the time delays of the LOS signal and the NLOS signal are similar, and the SECSI+HOSVD algorithm has high computational complexity although it has good performance. SUMMARY

[0005] In order to overcome the defects of the existing tensor-based time delay estimation algorithm, the application provides a generalized Schur decomposition time delay estimation method for GNSS array tensor modeling.

[0006] The application discloses a generalized Schur decomposition time delay estimation method for GNSS array tensor modeling.

[0007] The complex baseband signal output response received at M antennas is expressed as a tensor model:

[0008]

[0009] In the formula, is a complex amplitude of K time periods, is a pseudo-random binary sequence (PRBS) sample of L signals, wherein, corresponds to a satellite number d, a time delay is a pseudo-random binary sequence sample of L1, L2,..., L D D respectively are the number of signals generated by the 1st, 2nd,..., Dth star, L d includes 1 LOS signal and L d -1 NLOS signal, is a steering matrix, is a steering vector corresponding to an antenna array with a bearing angle l is a steering vector corresponding to an antenna array with a bearing angle is a Gaussian white noise tensor.

[0010] The tensor model of the data has three dimensions, wherein the first dimension and the second dimension are time dimensions, and the third dimension is a space dimension; the first dimension is a time period of sampling, and the size is K; the second dimension is N sampling points of each time period; and the third dimension corresponds to M antennas of the array in space.

[0011] The signals emitted by the satellites adopt code division multiple access technology, each satellite has a specific and fixed pseudo-random binary sequence. The receiver is provided with a correlator group to separate the satellite signals that are desired to be received from different satellite signals, and to separate the LOS signals and NLOS multipath components; assuming that the correlator tap number is R, the correlator is expressed as:

[0012]

[0013] In the formula, indicates a pseudo-random binary sequence sample of a satellite number d and a time delay .

[0014] Since directly introducing the correlator R d will change the noise into colored noise, R d singular value decomposition is performed to obtain instead of R d ; the tensor data model after correlation is expressed as:

[0015]

[0016] where, is the L d complex amplitudes of the d is the array steering matrix of the d d satellite in the L is the correlated Gaussian white noise, is the multipath interference tensor of other satellite signals, after correlation processing, the interference of other satellite signals can be approximately ignored. The correlated data tensor also has three dimensions: the first dimension is the time period of sampling, the size is K, the second dimension is the output of the sampled signal after correlation with the correlator group, the size is R, and the third dimension corresponds to the M antennas of the array in space.

[0017] The HOSVD low-rank approximation of the signal model is represented as:

[0018]

[0019] where, is the truncated core tensor, and are the factor matrices of the tensor decomposition truncated, considering the modulo 1 expansion of the data model as Combining the tensor data model (3) and equation (4), the following equation is established:

[0020]

[0021] Comparing equation (5) and equation (6), it can be seen that the column vectors of U (1) and Γ T span the same subspace, so there is an invertible matrix such that the following equation is established:

[0022] Γ T = U (1) T1 (7)

[0023] Similarly, the data model is expanded modulo 2 and modulo 3 respectively, and through similar analysis, it can be obtained that there is an invertible matrix and such that the following equation is established:

[0024]

[0025] A = U (3) T3 (9)

[0026] From equations (5), (6) and equations (7), (8), (9), the core tensor is represented as:

[0027]

[0028] For the core tensor there exist unitary matrices and such that:

[0029]

[0030] are upper triangular matrices, where matrices Q and Z are the unitary matrices of the matrix pencil decomposition.

[0031] Define the data tensor:

[0032]

[0033] Now and are upper triangular, and the frontal slice is expressed as:

[0034]

[0035] where is a matrix with the rth row element of T3 being 1 and other elements being 0, r = 1,..., L d ; obviously, (Z T T2) -T can be regarded as the generalized eigenvector of the matrix pencil , by sorting the column vectors appropriately, (Z T T2) -T is also an upper triangular matrix, and the corresponding Z T T2is a lower triangular matrix.

[0036] From equation (13), it can be deduced that QT1is also an upper triangular matrix, so the following equation is valid:

[0037]

[0038] Substituting the matrix obtained from equation (14) into equation (9) gives:

[0039]

[0040] Next, the data tensor is unfolded according to modulo 3, and multiplied by to obtain:

[0041]

[0042] The matrix F is the matrix and Γ TThe Khatri-Rao product result of the matrix is decomposed by using the LSKRF to solve the factor matrix and Next, the corresponding LOS signal component is found by using the factor matrix. Since the NLOS is obtained by reflection, refraction and other ways of the original signal through obstacles, the power thereof is usually lower than that of the LOS, so the maximum signal power method is used to select the LOS signal. The maximum signal power method is used to select the LOS signal, and all factor matrices are unitized:

[0043]

[0044]

[0045]

[0046] The matrix X is defined as:

[0047]

[0048] Wherein According to the structure of the signal tensor, the amplitude of the signal is calculated:

[0049]

[0050] Finally, the LOS component is selected:

[0051]

[0052] v :,ld is l d The corresponding amplitude.

[0053] Because the data model is introduced into the correlator, the correlator is preprocessed by SVD, so it is necessary to multiply Sigma V H The matrix restores the original data signal:

[0054]

[0055] Wherein, q represents the correlation value size of each tap of the correlator and the LOS signal; the vector q is interpolated to obtain the interpolation function F(τ) of the correlator and the LOS signal after correlation, and finally the F(τ) is used to estimate the LOS delay:

[0056]

[0057] At this point, the time delay estimation is completed.

[0058] Compared with the prior art, the present application has the beneficial technical effects that:

[0059] The application realizes high-precision and low-complexity time delay estimation of GNSS signals, and has higher robustness and lower calculation complexity. BRIEF DESCRIPTION OF DRAWINGS

[0060] Figure 1 A generalized Schur decomposition time delay estimation method framework for GNSS array tensor modeling of the application.

[0061] Figure 2 A simulation diagram for the case of 2 signals, different time delay differences between LOS signals and NLOS signals.

[0062] Figure 3 A simulation diagram for the case of 3 signals, different time delay differences between LOS signals and NLOS signals.

[0063] Figure 4 A simulation diagram for the case of 2 signals, different array position disturbances.

[0064] Figure 5 A simulation diagram for the case of 3 signals, different array position disturbances.

[0065] Figure 6 A simulation diagram for the case of 2 signals, different signal incidence angle differences.

[0066] Figure 7 A simulation diagram for the case of 3 signals, different signal incidence angle differences. DETAILED DESCRIPTION

[0067] The application will be further described in detail below in combination with the drawings and specific embodiments.

[0068] A generalized Schur decomposition time delay estimation method for GNSS array tensor modeling of the application is shown in Figure 1 The specific steps are as follows:

[0069] 1: Calculate the HOSVD of the signal tensor

[0070]

[0071] 2: Compress the core tensor and the singular vector matrix:

[0072]

[0073] 3: Generalized Schur decomposition of the two frontal slices of the core tensor to obtain Q and Z matrices:

[0074]

[0075]

[0076] ​where Q, Z are unitary matrices.

[0077] 4: Triangularize the core tensor frontal slice matrix:

[0078]

[0079] 5:

[0080] 6: Compute steering matrix

[0081]

[0082] 7:

[0083] 8: For l d = 1 to L d do

[0084] 9: Compute SVD of rank 1 matrix

[0085]

[0086] 10: Adjust using previously computed singular values so that u ld and are the corresponding singular vectors: ld

[0087]

[0088] 11: End for

[0089] 12: Concatenate previously computed eigenvectors to estimate and

[0090]

[0091]

[0092] 13: For l d = 1 to L d do

[0093] 14: Unitize T , and

[0094]

[0095]

[0096] ​​

[0097] 15: End for

[0098] 16: where, l d = 1,...,L d .

[0099] 17: According to the structure of the signal tensor, the amplitude of the signal is calculated:

[0100]

[0101] where, is the l d corresponding amplitude.

[0102] 18: Finally, the LOS component is selected:

[0103]

[0104] 19: where q represents the correlation value of each tap of the correlator and the LOS signal.

[0105] 20: The vector q is interpolated to obtain the interpolation function F(τ) after correlation of the correlator and the LOS signal.

[0106] 21: Finally, the LOS delay is estimated using F(τ):

[0107]

[0108] The embodiment is simulated by using MATLAB tool, and numerical simulation verifies the correctness and feasibility of the application. Figures 2 to 7 The simulation is the root mean square error of time delay estimation under different algorithms and different scenes in Monte Carlo simulation 2000 times. Among them, the Proposed HOSVD+QZ curve is the root mean square error of time delay estimation using HOSVD+QZ algorithm, the HOSVD+SECSI curve is the root mean square error of time delay estimation using HOSVD+SECSI algorithm, the CPD-GEVD curve is the root mean square error of time delay estimation using CPD-GEVD algorithm, the FBA+ESPS+HOSVD curve is the root mean square error of time delay estimation using FBA+ESPS+HOSVD algorithm, the DoA / KrF curve is the root mean square error of time delay estimation using DoA / KrF algorithm, the Known A and Γ curve is the root mean square error of time delay estimation with noise when A and Γ are known, and the Noiseless case curve is the root mean square error of time delay estimation without noise when A and Γ are known. Figure 2 It is a simulation diagram of different time delay differences between the LOS signal and the NLOS signal when the number of signals is 2. Figure 3Figure 3 is a simulation diagram of the time delay difference when the number of signals is 3 and the LOS signal and the NLOS signal are different. Figure 4 Figure 4 is a simulation diagram when the number of signals is 2 and the array position is different. Figure 5 Figure 5 is a simulation diagram when the number of signals is 3 and the array position is different. Figure 6 Figure 6 is a simulation diagram when the number of signals is 2 and the signal incidence angle difference is different. Figure 7 Figure 7 is a simulation diagram when the number of signals is 3 and the signal incidence angle difference is different.

[0109] Table 1 is a comparison of the calculation complexity of different time delay estimation algorithms

[0110]

[0111] Table 1 is the calculation complexity of different time delay estimation algorithms, wherein I K , I R , I M , I Ms , I Ls , I E , I F represents the power iteration number of the randomized SVD, J is the iteration number of the convergence of the simultaneous joint diagonalization matrix, and N is the number of slices. Numerical simulation shows that the present application has the characteristics of low calculation complexity and high precision.

Claims

1. A generalized Schur decomposition delay estimation method for GNSS array tensor modeling, characterized by: Specifically: The complex baseband signal output response received at the M antennas is represented as a tensor model: Where, is the complex amplitude of K time periods, is a pseudo-random binary sequence sampling of L signals, where Corresponding to the satellite number d, the time delay is Pseudo-random binary sequence sampling, L1, L2, ..., L D are the number of signals generated by the 1st, 2nd, ..., Dth stars, L d Includes 1 LOS signal and L d -1 NLOS signal, is the orientation matrix, The azimuth angle is φ l Corresponding to the steering vector of the antenna array, is the Gaussian white noise tensor; The tensor model of data has three dimensions, of which the first and second dimensions are time dimensions, and the third dimension is space dimension; the first dimension is the sampling time period, size K; the second dimension is N sampling points in each time period; the third dimension corresponds to the M antennas arrayed in space; The receiver is equipped with a correlator group to separate the desired satellite signal from different satellite signals, and to separate the LOS signal and NLOS multipath components. Assuming the number of correlator taps is R, the correlator is expressed as: Among them, c d [κ r ] indicates that the satellite number is d and the delay is κ r Pseudo-random binary sequence sampling of ; R d Singular value decomposition is obtained Instead of R d ; The relevant post-tensor data model is expressed as: in, yes L of the dth satellite d The complex amplitude of a signal, yes The corresponding satellite L d An array steering matrix of signals; is the correlated Gaussian white noise, is the multipath interference tensor of other satellite signals; the correlated data tensor also has three dimensions: the first dimension is the sampling time period, whose size is K, the second dimension is the output after the sampled signal is correlated with the correlator group, whose size is R, and the third dimension corresponds to the M antennas in the spatial array; Perform HOSVD low-rank approximation on the signal model, and the signal model is expressed as: in, is the truncated core tensor, and is the factor matrix of tensor decomposition truncation; combined with the tensor data model (3) and formula (4), the following formula holds: There exists an invertible matrix Make the following equation true; Γ T =U (1) T1(7) Expand the data model modulo 2 and modulo 3 respectively, and we can get the following: there is a reversible matrix and So that the following formula holds: A=U (3) T3(9) Combining Equations (5), (6) and (7), (8), and (9), the core tensor S can be expressed as: For the core tensor S, there exists a unitary matrix and So that: At the same time, it is an upper triangular matrix, where the matrices Q and Z are the matrix bundles Decomposed unitary matrix; define the data tensor: at this time and is upper triangulated, and the frontal slice is represented as: in, Represents a matrix whose diagonal elements are the r-th row elements of T3 and whose other elements are 0, r=1,...,L d ;(Z T T2) -T Matrix bundle The generalized eigenvectors of , sorted by column vectors, (Z T T2) -T It is also an upper triangular matrix, and the corresponding Z T T2 is a lower triangular matrix. From equation (13), we can infer that QT1 is also an upper triangular matrix. Therefore, the following equation holds: The matrix obtained by formula (14) Substituting into formula (9) we get: Next, the data tensor is expanded modulo 3 and compared with Multiplying them together gives: The matrix F is the matrix and Γ T The Khatri-Rao product result is used to decompose the matrix using LSKRF to solve the factor matrix and Next, we use the factor matrix to find the corresponding LOS signal components; we use the maximum signal power method to select the LOS signal and normalize all factor matrices: Define the matrix X: in According to the structure of the signal tensor, calculate the amplitude of the signal: Finally, select the LOS components: v :,ld Yes d The corresponding amplitude; Multiply ΣV H The matrix restores the original data signal: Where q represents the correlation value between each tap of the correlator and the LOS signal. The vector q is interpolated to obtain the correlation interpolation function F(τ) between the correlator and the LOS signal. Finally, F(τ) is used to estimate the LOS delay: At this point, the delay estimation is complete.

Images