GNSS multipath parameter estimation method based on multi-correlator output

Abstract

The invention discloses a GNSS multipath parameter estimation method based on multi-correlator output, and belongs to the field of satellite signal multipath resistance. The invention discloses a GNSS multipath signal parameter estimation method, and aims to provide a novel GNSS multipath signal parameter estimation method which comprises the following steps: performing correlation operation on a local spreading code and a received signal through a GNSS correlator, and mining internal correlation between a correlator interval and a multipath error based on a correlation output result of a composite signal; the method comprises the following steps: acquiring multipath error data of multiple groups of correlators to construct a relation curve by combining an analytical formula of correlator spacing and multipath errors, reversely deducing multipath parameters by solving coordinates of turning points of the curve, and finally realizing accurate estimation of the multipath parameters. According to the method, the time delay, phase and amplitude parameters of the short multipath can be accurately identified, and the capability bottleneck of a traditional algorithm in the aspects of short multipath suppression and estimation is effectively broken.

Description

Technical Field

[0001] This invention relates to a parameter estimation algorithm, specifically to a GNSS multipath parameter estimation method based on the output of a multicorrelator, belonging to the field of satellite signal anti-multipath. Background Technology

[0002] GNSS technology is a core tool for modern navigation, positioning, and timing, and is widely used in transportation, agriculture, autonomous driving, and many other fields. However, multipath effect, as a key interference, severely restricts its performance—satellite signals propagate through multiple paths after reflecting off objects, and the superposition of these paths leads to decreased receiver positioning accuracy, signal loss, and error accumulation. This effect is particularly significant in complex environments such as urban high-rise buildings and tunnels. Therefore, multipath suppression technology has become a core research direction in the field of GNSS.

[0003] Current research on multipath suppression mainly falls into two categories: nonparametric and parametric algorithms. Nonparametric algorithms suppress multipath interference by optimizing the receiver loop architecture or designing novel phase detector functions. Parametric algorithms, on the other hand, use multicorrelator arrays to obtain the characteristic parameters of multipath signals and combine them with parameter estimation algorithms to separate direct signals from multipath interference.

[0004] Despite some progress in related research, existing technologies still have significant shortcomings: they struggle to effectively suppress short multipath events (such as short-delay multipath events with path differences less than 30m), and some methods suffer from high complexity or even non-convergence. Because multipath events in complex environments often exhibit short-range, short-delay characteristics, the time interval between direct and multipath signals is extremely small, drastically increasing the difficulty of parameter estimation and causing errors to peak. Therefore, exploring high-precision, low-complexity methods to combat short multipath events is crucial for improving the reliability of GNSS in complex scenarios. Summary of the Invention

[0005] To address the aforementioned technical problems, this invention proposes a GNSS multipath parameter estimation method based on the output of multiple correlators.

[0006] A method for estimating GNSS multipath parameters based on the output of multiple correlators, the specific steps of which are as follows:

[0007] Step 1: Collect N sets of "correlator spacing dS curve zero-crossing offset (SCB) ε" SCB "Value" data, forming a set:

[0008] {d p ,ε SCB,p}(p=1,2,3,...,N);

[0009] Step 2: Transfer the second set of data (d) n ,ε SCB,n(n=2) is used as the initial candidate turning point. Curve fitting is performed on the first half and the second half of the data respectively, and the sum of the fitting errors of the two segments is calculated.

[0010] Step 3: Select each subsequent data set (d) in sequence. n+1 ,ε SCB,n+1 (n = 3, 4, ..., N-1) are used as new candidate inflection points. Step 2 is repeated to obtain the fitting error corresponding to each candidate point, forming an error set {e2, e3, ..., e...} N-1}, the candidate inflection point (d) corresponding to the minimum fitting error in the error set. k ,ε SCB,k This is the true turning point;

[0011] Step 4: Set the actual inflection point coordinates (d) k ,ε SCB,k Substitute the values ​​into the multipath delay analytical formula to estimate the delay of the multipath signal;

[0012] Step 5: Change the phase of the local carrier and repeat steps 1 to 4 twice to obtain the inflection points of the local carrier with three different phases. Combine the multipath delay results to estimate the phase and amplitude of the multipath signal.

[0013] The step 1 involves collecting N sets of "correlation spacing dS curve zero-crossing offset ε". SCB Value data, forming a set: {d p ,ε SCB,p}(p=1,2,3,...,N), the process is as follows: ε is derived from the relationship between multipath delay and code phase error. SCB The piecewise function expression for d, thus yielding {d p ,ε SCB,p Data set}(p=1,2,3,...,N).

[0014]

[0015] Step 2 describes transferring the second set of data (d) n ,ε SCB,n (n=2) is used as the initial candidate inflection point. Curve fitting is performed on the first half and the second half of the data respectively, and the sum of the fitting errors of the two segments is calculated. The process is as follows: For any inflection point The bisegmental linear model is constructed as follows:

[0016] and

[0017] The optimal coefficients for each segment are solved using the Iterative Weighted Least Squares (IRLS) method of Huber regression. Initially, the weight matrix is ​​set as the identity matrix, and the initial coefficients are obtained using the ordinary least squares solution. Iterative updates are performed, and the residual r is calculated. (p) =Y-Xβ (p) Update the weights according to the group rules that Huber lost. Using weighted least squares to update coefficient β (t+1) =(X T W (t+1) X) -1 X T W (t+1) The iteration stops when the coefficient changes or the weighted loss decreases to a threshold. Applying the above IRLS to the two datasets respectively yields the Huber regression coefficients for the two datasets. and And calculate the total fitting error.

[0018] Step 3 describes the sequential selection of each subsequent data set (d). n+1 ,ε SCB,n+1 (n = 3, 4, ..., N-1) are used as new candidate inflection points. Step 2 is repeated to obtain the fitting error corresponding to each candidate point, forming an error set {e2, e3, ..., e...} N-1}, the candidate inflection point (d) corresponding to the minimum fitting error in the error set. k ,ε SCB,k The true turning point is reached by repeating step 2 of claim 1 to obtain the total fitting error set of all SCB results. With minimum error k opt =argmin k∈p e (k) The SCB result for p = {2, 3, ..., N-1} is the coordinate of the inflection point.

[0019] Step 4 describes the process of obtaining the actual inflection point coordinates (d). k ,ε SCB,k Substituting the multipath delay analytical formula, the delay of the multipath signal is estimated as follows: The multipath delay can be estimated using the inflection point coordinates of the SCB result set. Since the inflection point coordinates of the SCB result set change with the correlator spacing, to simplify the calculation, the inflection point coordinates are redefined as... in and Numerically equal to d k and ε SCB,k (d k ), different multipath delay scenarios And Δτ1 is expressed as

[0020] Will and Substitute ε SCB The piecewise function expression for d can be expressed as a multipath delay:

[0021]

[0022] Step 5 involves changing the phase of the local carrier. Steps 1 to 4 are repeated twice to obtain three sets of inflection points under local carriers with different phases. Combined with the multipath delay results, the phase and amplitude of the multipath signal are estimated. The process is as follows: using different phases... The inflection point positions of the SCB result set obtained from the local carrier can be obtained with different values. Substituting this into the following formula, the phases of the direct signal and the multipath signal can be estimated.

[0023]

[0024] By establishing Mathematical relationship with amplitude ratio a1 / a0 By substituting the expressions for the phase estimate φ0 of the direct signal and the phase estimate φ1 of the multipath signal, the amplitude of the multipath signal can be estimated.

[0025] Based on the above steps, the GNSS multipath parameter estimation method proposed in this invention, which is based on the output of multiple correlators, can achieve the following effects:

[0026] I. It can accurately identify the time delay, phase and amplitude parameters of short multipath, effectively solving the problem of the lack of short multipath suppression and estimation capabilities of traditional algorithms;

[0027] Second, while ensuring the accuracy of the estimation, the computational complexity is significantly reduced, balancing accuracy and efficiency.

[0028] Third, based on the multi-correlator output architecture, it has good compatibility with new GNSS system signals;

[0029] In summary, compared with existing GNSS multipath parameter estimation algorithms, this invention has stronger versatility and higher robustness, and has lower computational complexity compared with other parameter-based algorithms, thus demonstrating wider applicability and stronger practicality. Attached Figure Description

[0030] Figure 1 This is a flowchart illustrating the implementation process of the method described in this invention;

[0031] Figure 2 This is a structural block diagram of the method described in this invention;

[0032] Figure 3This is a graph showing the relationship between correlator spacing and code ring tracking error under different phase local carriers; Specific implementation methods

[0033] The following description, in conjunction with the accompanying drawings and specific embodiments, provides a further introduction to a GNSS multipath parameter estimation method based on multiple correlator output according to the present invention.

[0034] like Figure 1 As shown, the present invention provides a method for estimating S-curve fitting parameters based on the output of multiple correlators. The specific implementation steps are as follows:

[0035] Step 1: Collect N sets of "correlation spacing dS curve zero-crossing offset ε" SCB Value data, forming a set: {d p ,ε SCB,p}(p=1,2,3,...,N).

[0036] Build as Figure 2 The code phase tracking loop structure of the GNSS shown is represented by the GNSS receiver receiving signal as follows:

[0037] Where c(t) represents the spreading code sequence; a i τ i φ i Let represent the amplitude, delay, and phase of the i-th path, respectively; i = 0 represents a direct signal, and the rest represent multipath signals; ω is the carrier frequency. This refers to the carrier phase.

[0038] Receivers typically use a correlator to perform correlation calculations between the local spreading code and the received signal. The correlation output R of the composite signal... Σ (ε) can be expressed as: If multipath signals exist, the discriminator output is an "S"-shaped curve formed by the superposition of the direct signal and the multipath signals.

[0039]

[0040] Collect the zero-crossing frequency shift values ​​of the S-curve under different correlator spacings to form a set: {d p ,ε SCB,p}(p=1,2,3,...,N).

[0041] Step 2: Transfer the second set of data (d) n ,ε SCB,n (n=2) is used as the initial candidate turning point. Curve fitting is performed on the first half and the second half of the data respectively, and the sum of the fitting errors of the two segments is calculated.

[0042] For any turning point The bisegmental linear model is constructed as follows, and the optimal coefficients of each segment are solved by iterative weighted least squares regression of Huber regression.

[0043] The first model is:

[0044] The second model is:

[0045] For residuals Define the Huber loss function after threshold δ

[0046]

[0047] Its corresponding influence function is

[0048]

[0049] To minimize Iterative weighted least squares can be used:

[0050] Initialize by letting the weight matrix W (0) =I, using the ordinary least squares solution as the initial value.

[0051] Iterate t = 0, 1, 2, ... until convergence:

[0052] Calculate the residuals:

[0053]

[0054] Calculate the weight for each sample:

[0055]

[0056] Update using weighted least squares:

[0057] make Update the optimal coefficients:

[0058] β (t+1) =argmin β (Y-Xβ) T W (t+1) (Y-Xβ)=(X T W (t+1) X) -1 X T W (t+1) Y

[0059] Convergence criterion: When ||β|| (t+1) -β (t) Stop when || = 0 or the weighted loss decreases less than a preset threshold.

[0060] The parameter vector and the XY matrix are defined as follows:

[0061]

[0062] The IRLS algorithm described above is applied to each of the two segments until β1 and β2 converge, ultimately yielding the Huber regression coefficients for the first and second segments:

[0063]

[0064] And calculate the total fitting error.

[0065]

[0066] The composite error is:

[0067]

[0068] Step 3: Select subsequent data for each group (d) n+1 ,ε SCB,n+1 (n = 3, 4, ..., N-1) are used as new candidate turning points. The second step is repeated to obtain the fitting error of the line segments to the left and right of each candidate point, forming an error set {e2, e3, ..., e...}. N-1}

[0069] The candidate inflection point k corresponding to the minimum fitting error in the error set opt =argmin k∈p e (k) The SCB result for p = {2, 3, ..., N-1} is the coordinate of the inflection point.

[0070] Step 4: Set the actual turning point coordinates (d) k ,ε SCB,k Substitute the values ​​into the multipath delay analytical formula to estimate the delay of the multipath signal.

[0071] Since the inflection point coordinates of the SCB result set change with the correlator spacing, the inflection point coordinates are redefined to simplify the calculation. in and Numerically equal to d k and ε SCB,k (d k ), different multipath delay scenarios And Δτ1 are expressed as:

[0072]

[0073] Will and Substitute ε SCB The piecewise function expression for d can be expressed as a multipath delay:

[0074]

[0075] Step 5: Change the phase of the local carrier and repeat steps 1 to 4 to obtain the result.

[0076] Using different phases The inflection point position of the SCB result set obtained from the local carrier is as follows: Figure 3 As shown, by The expressions for Δτ1 and Δτ1 can yield different values. Substituting this into the following formula, the phases of the direct signal and the multipath signal can be estimated.

[0077]

[0078] By establishing The mathematical relationship between the amplitude ratio a1 / a0 and the amplitude ratio:

[0079]

[0080] By substituting the expressions for the phase estimate φ0 of the direct signal and the phase estimate φ1 of the multipath signal, the amplitude of the multipath signal can be estimated.

Claims

1. A method for estimating GNSS multipath parameters based on the output of multiple correlators, characterized in that: The specific steps are as follows: Step 1: Collect N sets of "correlator spacing dS curve zero-crossing offset (SCB) ε" SCB "Value" data, forming a set: {d} p ,he SCB,p }(p=1,2,3,...,N); Step 2: Transfer the second set of data (d) n ,ε SCB,n (n=2) is used as the initial candidate turning point. Curve fitting is performed on the first half and the second half of the data respectively, and the sum of the fitting errors of the two segments is calculated. Step 3: Select each subsequent data set (d) in sequence. n+1 ,ε SCB,n+1 (n = 3, 4, ..., N-1) are used as new candidate inflection points. Step 2 is repeated to obtain the fitting error corresponding to each candidate point, forming an error set {e2, e3, ..., e...} N-1 }, the candidate inflection point (d) corresponding to the minimum fitting error in the error set. k ,ε SCB,k This is the true turning point; Step 4: Set the actual inflection point coordinates (d) k ,ε SCB,k Substitute the values ​​into the multipath delay analytical formula to estimate the delay of the multipath signal; Step 5: Change the phase of the local carrier and repeat steps 1 to 4 twice to obtain the inflection points of the local carrier with three different phases. Combine the multipath delay results to estimate the phase and amplitude of the multipath signal.

2. As described in claim 1, characterized in that: The zero-crossing offset ε of the N sets of correlator spacing dS curves mentioned in step 1 is collected. SCB Value data, forming a set: {d p ,ε SCB,p }(p=1,2,3,...,N), the process is as follows: ε is derived from the relationship between multipath delay and code phase error. SCB The piecewise function expression for d, thus yielding {d p ,ε SCB,p Data set}(p=1,2,3,...,N).

3. As described in claim 1, characterized in that: Step 2 describes transferring the second set of data (d) n ,ε SCB,n (n=2) is used as the initial candidate inflection point. Curve fitting is performed on the first half and the second half of the data respectively, and the sum of the fitting errors of the two segments is calculated. The process is as follows: For any inflection point The bisegmental linear model is constructed as follows: and The optimal coefficients for each segment are solved using the Iterative Weighted Least Squares (IRLS) method of Huber regression. Initially, the weight matrix is ​​set as the identity matrix, and the initial coefficients are obtained using the ordinary least squares solution. Iterative updates are performed, and the residual r is calculated. (p) =Y-Xβ (p) Update the weights according to the group rules that Huber lost. Using weighted least squares to update coefficient β (t+1) =(X T W (t+1) X) -1 X T W (t+1) The iteration stops when the coefficient changes or the weighted loss decreases to a threshold. Applying the above IRLS to the two datasets respectively yields the two Huber regression coefficients. and And calculate the total fitting error.

4. As described in claim 1, characterized in that: Step 3 describes the sequential selection of each subsequent data set (d). n+1 ,ε SCB,n+1 (n = 3, 4, ..., N-1) are used as new candidate inflection points. Step 2 is repeated to obtain the fitting error corresponding to each candidate point, forming an error set {e2, e3, ..., e...} N-1 }, the candidate inflection point (d) corresponding to the minimum fitting error in the error set. k ,ε SCB,k The true turning point is reached by repeating step 2 of claim 1 to obtain the total fitting error set of all SCB results. With minimum error k opt =argmin k∈p e (k) The SCB result for p = {2, 3, ..., N-1} is the coordinate of the inflection point.

5. As described in claim 1, characterized in that: Step 4 describes the process of obtaining the actual inflection point coordinates (d). k ,ε SCB,k Substituting the multipath delay analytical formula, the delay of the multipath signal is estimated as follows: The multipath delay can be estimated using the inflection point coordinates of the SCB result set. Since the inflection point coordinates of the SCB result set change with the correlator spacing, to simplify the calculation, the inflection point coordinates are redefined as... in and Numerically equal to d k and ε SCB,k (d k ), under different multipath delay scenarios And Δτ1 is expressed as Will and Substitute ε SCB The piecewise function expression for d can be expressed as a multipath delay:

6. As described in claim 1, characterized in that: Step 5 involves changing the phase of the local carrier. Steps 1 to 4 are repeated twice to obtain three sets of inflection points under local carriers with different phases. Combined with the multipath delay results, the phase and amplitude of the multipath signal are estimated. The process is as follows: using different phases... The inflection point positions of the SCB result set obtained from the local carrier can be obtained with different values. Substituting this into the following formula, the phases of the direct signal and the multipath signal can be estimated. By establishing Mathematical relationship with amplitude ratio a1 / a0 By substituting the expressions for the phase estimate φ0 of the direct signal and the phase estimate φ1 of the multipath signal, the amplitude of the multipath signal can be estimated.

Images