An improved direct position estimation method for satellite navigation

Abstract

The present invention discloses an improved satellite navigation direct position estimation method, which relates to the field of satellite navigation technology. The method comprises the following steps: step S1: constructing a multi-snapshot model of satellite signals at a user receiving end in combination with prior information; step S2: obtaining a multivariate optimization parameter estimation problem based on the multi-snapshot model; step S3: performing a coarse search based on a spatial alternating generalized expectation maximization algorithm with a variable step size; and step S4: performing a fine search based on a spatial alternating generalized expectation maximization algorithm with a fixed step size. The present invention obtains more accurate received signal information through prior information, and then performs estimation through the spatial alternating generalized expectation maximization algorithm, thereby achieving accurate positioning in a multipath environment.

Description

Technical Field

[0001] The present invention relates to the technical field of satellite navigation, and in particular to an improved satellite navigation direct position estimation method. Background Art

[0002] The statements in this section merely provide background information related to the present disclosure and may not constitute prior art.

[0003] Traditional Global Navigation Satellite System (GNSS) receivers are typically designed to operate under open sky conditions. Their positioning speed and reliability in simple signal environments meet people's needs. However, when there is no line of sight between mobile devices and satellites, the positioning performance of traditional receivers will be severely degraded. Direct Position Estimation (DPE) is a new satellite navigation positioning method that can meet this challenge.

[0004] The general principle of the direct position estimation method is to first use the known satellite ephemeris information under the current sky state as the prior information required for model construction, and then combine the information of multiple satellites, position and velocity information, pseudo code and Doppler relationship to construct the signal form closest to the received signal, and estimate the optimal position and velocity state information through the maximum likelihood method. However, it is not suitable for accurate positioning of satellite navigation in harsh multipath environments. Summary of the Invention

[0005] The purpose of the present invention is to address the problem that current direct position estimation methods are not suitable for accurate positioning of satellite navigation in harsh multipath environments. An improved direct position estimation method for satellite navigation is provided. The method utilizes known satellite ephemeris prior information to reconstruct satellite signals at the receiving end and, in combination with coarse and fine search, estimates the user's accurate position. Direct position estimation is achieved through an improved space alternating generalized expectation maximization (Sage) algorithm, achieving accurate estimation of the user's position in harsh environments such as multipath, thereby solving the above-mentioned problem.

[0006] The technical solutions of the present invention are as follows:

[0007] An improved satellite navigation direct position estimation method, comprising:

[0008] Step S1: Construct a multi-snapshot model of satellite signals at the user receiving end by combining prior information;

[0009] Step S2: Based on the multi-snapshot model, a multivariate optimization parameter estimation problem is obtained;

[0010] Step S3: performing a coarse search based on a spatial alternating generalized expectation maximization algorithm with a variable step size;

[0011] Step S4: Perform a fine search based on a fixed-step-size spatial alternating generalized expectation-maximization algorithm.

[0012] Furthermore, the step S1 includes:

[0013]

[0014] Where: M is the number of satellites, the first term in the formula is the direct signal, the second term is the noise signal, a i is the signal amplitude of the i-th satellite, s i is a pseudo code signal, n(t) is a complex Gaussian white noise signal, γ is the user information to be estimated, τ i is the time delay between the i-th satellite and the user, Δf i is the Doppler shift between the i-th satellite and the user.

[0015] Furthermore, the i-th satellite signal amplitude a i =α i exp{jφ i,0};

[0016] The complex white Gaussian noise signal

[0017] The user information to be estimated

[0018] Where: x, y, z correspond to the user's position coordinates, δt is the time difference between the user's receiver clock and the GPS system, Corresponding to the user's speed, is the user receiver clock drift;

[0019] The time delay between the i-th satellite and the user

[0020] Among them, r i =[xx i yy i zz i ] · is the user's position relative to the satellite, δt i is the time difference between the satellite clock and the GPS system, and c is the speed of light;

[0021] The Doppler frequency shift between the i-th satellite and the user

[0022] Where: f c is the carrier frequency of the satellite signal. is the user's velocity relative to the satellite, is the clock drift of the user receiver clock relative to the GPS system time, It is the clock drift of the satellite clock relative to the GPS system time.

[0023] Furthermore, the step S1 further includes:

[0024] x=aD(γ)+n

[0025] Where: x is the observed signal of K snapshots, a is the composite signal amplitude matrix, D(γ) is the delay-Doppler matrix, and n represents the zero mean and variance within the observation interval. k snapshots of complex Gaussian white noise signal.

[0026] Furthermore, the observation signals of the K snapshots

[0027] The composite signal amplitude matrix a=[a1 a2...a M ];

[0028] The delay-Doppler matrix D(γ)=[d(t0)…d(t K-1 )];

[0029] The observation interval has zero mean and variance k snapshots of complex Gaussian white noise signal

[0030] Furthermore, the step S2 includes:

[0031] Maximizing the likelihood function x = aD(γ) + n is equivalent to minimizing the nonlinear least squares problem Λ(a,γ) = ‖x-aD(γ)‖ 2 , we get a multivariate optimization parameter estimation problem where Λ(γ)=xx H -xD H (γ)(D(γ)D H (γ)) -1 D(γ)x H .

[0032] Furthermore, the step S3 includes:

[0033] Set the variables to be estimated in the model to x, y, and z, and the other four velocity variables to 0. To simplify the algorithm search, the clock error is also set to 0;

[0034] In the variable step size part, first determine the initial search range, the x, y, and z ranges are all [-L, L], and the search step size is 10 a , the search range after each update will also change with the corresponding search factor, and the specific value is {x,y,z}∈[p0(i)-bound,p0(i)+bound], where p0=[x,y,z] is the last estimated value, bound is the search boundary, and bound=10 a+2 ;

[0035] In the Sage algorithm part, each step size corresponds to a Sage algorithm search, and the number of Sage iterations is set to N0. After all searches are completed, a rough estimate of the position is obtained.

[0036] Furthermore, the step S4 includes:

[0037] First, the initial parameter estimate is given according to the result of step S3, and i=1 is set to calculate the cost function

[0038] Second, minimize the cost function to update the current parameters

[0039] Then let i=i+1 and determine whether all current variables have been searched. If i≤L, return to search the next variable. Otherwise, proceed to the next step.

[0040] Finally, the number of iterations k = k + 1, if all iterations are completed, the final estimation result is obtained

[0041] Compared with the existing technology, the beneficial effects of the present invention are:

[0042] An improved satellite navigation direct position estimation method uses prior information to obtain more accurate received signal information, and then uses the spatial alternating generalized expectation maximization algorithm to perform estimation, achieving accurate positioning in multipath environments. BRIEF DESCRIPTION OF THE DRAWINGS

[0043] Figure 1 This is a rough search flow chart;

[0044] Figure 2 This is the Sage algorithm flow chart;

[0045] Figure 3 This is the positioning error diagram under different carrier-to-noise ratios. DETAILED DESCRIPTION

[0046] It should be noted that relational terms such as "first" and "second" are used only to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any actual relationship or order between these entities or operations. Moreover, the terms "comprises," "comprising," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus comprising a series of elements includes not only those elements, but also other elements not explicitly listed, or elements inherent to such process, method, article, or apparatus. In the absence of further limitations, an element defined by the phrase "comprising a ..." does not exclude the presence of additional identical elements in the process, method, article, or apparatus comprising the element.

[0047] The features and performance of the present invention are further described in detail below with reference to the embodiments.

[0048] Example 1

[0049] See also Figure 1-3 , an improved satellite navigation direct position estimation method, specifically comprising the following steps:

[0050] Step S1: Construct a multi-snapshot model of satellite signals at the user receiving end by combining prior information;

[0051] Step S2: Based on the multi-snapshot model, a multivariate optimization parameter estimation problem is obtained;

[0052] Step S3: performing a coarse search based on a spatial alternating generalized expectation maximization algorithm with a variable step size;

[0053] Step S4: Perform a fine search based on a fixed-step-size spatial alternating generalized expectation-maximization algorithm.

[0054] In this embodiment, specifically, step S1 includes:

[0055]

[0056] in:

[0057] M is the number of satellites;

[0058] The first term in the formula is the direct signal;

[0059] The second term is the noise signal;

[0060] a i is the signal amplitude of the i-th satellite;

[0061] s i It is a pseudo code signal;

[0062] n(t) is a complex Gaussian white noise signal;

[0063] γ is the user information to be estimated;

[0064] τ i is the time delay between the i-th satellite and the user;

[0065] Δf i is the Doppler shift between the i-th satellite and the user.

[0066] In this embodiment, specifically, the i-th satellite signal amplitude a i =α i exp{jφ i,0};

[0067] The complex white Gaussian noise signal

[0068] The user information to be estimated

[0069] Where: x, y, z correspond to the user's position coordinates, δt is the time difference between the user's receiver clock and the GPS system, Corresponding to the user's speed, is the user receiver clock drift;

[0070] The time delay between the i-th satellite and the user

[0071] Among them, r i =[xx i yy i zz i ] · is the user's position relative to the satellite, δt i is the time difference between the satellite clock and the GPS system, and c is the speed of light;

[0072] The Doppler frequency shift between the i-th satellite and the user

[0073] Where: f c is the carrier frequency of the satellite signal. is the user's velocity relative to the satellite, is the clock drift of the user receiver clock relative to the GPS system time, It is the clock drift of the satellite clock relative to the GPS system time.

[0074] In this embodiment, specifically, step S1 further includes:

[0075] x=aD(γ)+n

[0076] Where: x is the observed signal of K snapshots, a is the composite signal amplitude matrix, D(γ) is the delay-Doppler matrix, and n represents the zero mean and variance within the observation interval. k snapshots of complex Gaussian white noise signal.

[0077] In this embodiment, specifically, the K snapshot observation signals

[0078] The composite signal amplitude matrix a=[a1 a2...a M ];

[0079] The delay-Doppler matrix D(γ)=[d(t0)…d(t K-1 )];

[0080] The observation interval has zero mean and variance k snapshots of complex Gaussian white noise signal

[0081] In this embodiment, specifically, step S2 includes:

[0082] Maximizing the likelihood function x = aD(γ) + n is equivalent to minimizing the nonlinear least squares problem Λ(a,γ) = ‖x-aD(γ)‖ 2 , we get a multivariate optimization parameter estimation problem where Λ(γ)=xx H -xD H (γ)(D(γ)D H (γ)) -1 D(γ)x H .

[0083] In this embodiment, specifically, step S3 includes:

[0084] Set the variables to be estimated in the model to x, y, and z, and the other four velocity variables to 0 (i.e., the user is a fixed point). To simplify the algorithm search, the clock error is also set to 0;

[0085] In the variable step size part, first determine the initial search range, the x, y, and z ranges are all [-L, L], and the search step size is 10 a (a is called the search factor, and it is decremented by a=a-2 each time through the loop, and the loop repeats three times.) The search range after each update will also change accordingly with the search factor. The specific value is {x,y,z}∈[p0(i)-bound,p0(i)+bound], where p0=[x,y,z] is the last estimated value, bound is the search boundary, and bound=10 a+2 ;

[0086] In the Sage algorithm part, each step size corresponds to a Sage algorithm search, and the number of Sage iterations is set to N0. After all searches are completed, a rough estimate of the position is obtained.

[0087] In this embodiment, specifically, step S4 includes:

[0088] First, the initial parameter estimate is given according to the result of step S3, and i=1 is set to calculate the cost function

[0089] Second, minimize the cost function to update the current parameters

[0090] Then let i=i+1 and determine whether all current variables have been searched. If i≤L, return to search the next variable. Otherwise, proceed to the next step.

[0091] Finally, the number of iterations k = k + 1, if all iterations are completed (k>N t ), and obtain the final estimation result

[0092] Example 2

[0093] The second embodiment is a specific application of an improved satellite navigation direct position estimation method proposed in the first embodiment.

[0094] Assume that the current scenario is fixed-point user positioning. There are 6 satellites in the sky environment. There are reflected signals from 3 satellites reaching the user receiving end. The ephemeris information of 6 satellites is known. To obtain the user location information in the current situation, the specific implementation method is as follows.

[0095] Step 1: Combine the known satellite position and velocity information to build a multi-snapshot model of the satellite signal at the user receiving end:

[0096]

[0097] Among them, the number of satellites is 6, a i is the signal amplitude of the i-th satellite, s i is a pseudo-code signal, γ=[xyz 0 0 0 00] T is the user location information to be estimated.

[0098] Step 2: The K snapshot received signals are:

[0099] x=aD(γ)+n

[0100] Among them, the signal amplitude matrix The value of and the noise n are determined according to the specific carrier-to-noise ratio;

[0101] is the delay-Doppler matrix;

[0102] D(γ)=[d(t0)…d(t K-1 )]

[0103] in, d i =s i (t-τ i (γ))exp{j2πΔf i (γ)t}.

[0104] Step 3: Minimize the nonlinear least squares problem Λ(a,γ)=‖x-aD(γ)‖ 2 , we get a multivariate optimization parameter estimation problem where Λ(γ)=xx H -xD H (γ)(D(γ)D H (γ)) -1 D(γ)x H .

[0105] Step 4: Perform a rough search based on the variable step-size space alternating generalized expectation maximization algorithm. The algorithm flowchart is as follows: Figure 1 As shown in the figure, here the variables to be estimated in the model are set to 3, which are x, y, and z coordinates respectively. The other four speed variables are set to 0 (i.e. the user is a fixed point). In order to simplify the algorithm search, the clock error is also set to 0. For the variable step size, first determine that the initial search range x, y, and z are all in the range of [-6386137, 6386137], and the search step size is 10 a (where a is called the search factor, the initial value of a is 5, and it is decremented by a=a-2 each time through the loop, and the loop repeats three times.) The corresponding step sizes are 100km, 1km, and 10m, respectively. The search range after each update will also change with the search factor. The specific value is {x,y,z}∈[p0(i)-bound,p0(i)+bound], where p0=[x,y,z] is the last estimated value, bound is the search boundary, and bound=10 a+2 In the Sage algorithm part, each step size corresponds to a Sage algorithm search. The number of Sage iterations is set to 10. The first step size searches according to the initial given range. The number of searches corresponding to a 100km step size is The second and third steps are searched according to the updated search range, and the number of searches is and After all searches are completed, a rough estimate of the location is obtained.

[0106] Step 5: Fine search, perform fine search based on fixed-step spatial alternating generalized expectation maximization algorithm, Sage algorithm such as Figure 2 First, according to the result of step 4, the initial parameter estimate is given, let i = 1, and calculate the cost function Second, minimize the cost function to update the current parameters Then let i = i + 1, then determine whether all current variables have completed the search. If i ≤ 3, return to search the next variable, otherwise, proceed to the next step; finally, the number of iterations k = k + 1, if all iterations are completed (k > 11), then the final estimation result is obtained

[0107] Description of simulation conditions.

[0108] The simulation experiment was conducted using MATLAB R2021a on a Windows 10 operating system. The simulation setup included six satellites visible in the sky, all with direct signals. Three of these satellites had reflected signals, with the amplitude of the reflected signal being half that of the direct signal. The amplitudes of the three reflected signals were consistent. The multipath delay was one chip, and the carrier-to-noise ratio varied incrementally from 30 to 55. The reference coordinates were x = -1723586.255, y = 5005657.9520, and z = 3547005.7570.

[0109] The positioning error under different carrier-to-noise ratios was simulated and analyzed. Figure 3 The following are positioning error curves for different carrier-to-noise ratios. The horizontal axis represents the carrier-to-noise ratio, and the vertical axis represents the absolute value of the difference between the estimated coordinates and the true position coordinates. As can be seen, the higher the carrier-to-noise ratio, the smaller the positioning error. When the carrier-to-noise ratio reaches 35dB-Hz, the x and y coordinate positioning accuracy is already relatively low, within 5 meters. At this point, the height error is relatively large, approximately 20 meters. When the carrier-to-noise ratio reaches 45dB-Hz, the x and y coordinate positioning accuracy is very low, around 1 meter. The height error is also smaller, at approximately 3 meters. Therefore, the DPE method has relatively good position estimation performance in multipath environments and also provides relatively good position estimation at low carrier-to-noise ratios.

[0110] The above-described embodiments merely represent specific implementation methods of the present application. While the descriptions are relatively specific and detailed, they should not be construed as limiting the scope of protection of the present application. It should be noted that a person skilled in the art would be able to make numerous variations and improvements without departing from the technical concept of the present application, and all such variations and improvements fall within the scope of protection of the present application.

[0111] This background section is provided to generally present the context of the invention, and the work of the presently named inventors, the work to the extent described in this background section, and aspects of the description in this section that did not constitute prior art at the time of filing are neither explicitly nor implicitly admitted to be prior art to the present invention.

Claims

1. An improved satellite navigation direct position estimation method, characterized in that: include: Step S1: Construct a multi-snapshot model of satellite signals at the user receiving end by combining prior information; Step S2: Based on the satellite signal multi-snapshot model, a multivariate optimization parameter estimation problem is obtained; Step S3: performing a coarse search based on a spatial alternating generalized expectation maximization algorithm with a variable step size; Step S4: performing a fine search based on a fixed-step-size spatial alternating generalized expectation-maximization algorithm; The satellite signal multi-snapshot model in step S1 is expressed as follows: Where: M is the number of satellites, It's a direct signal. It is i Satellite signal amplitude, is a pseudo-code signal, is a complex white Gaussian noise signal, is the user information to be estimated, It is i The time delay between the satellite and the user, It is i Doppler shift between the satellite and the user; The said i Satellite signal amplitude ; The complex white Gaussian noise signal ; The user information to be estimated ; Among them: x, y, z correspond to the user's position coordinates, is the time difference between the user receiver clock and the GPS system clock, 、 、 Corresponding to the user's speed, is the user receiver clock drift; The said i The delay between the satellite and the user ; in, is the user's position relative to the satellite, is the time difference between the satellite clock and the GPS system, and c is the speed of light; The said i Doppler shift between satellites and user ; in: is the carrier frequency of the satellite signal. is the user's velocity relative to the satellite, is the clock drift of the user receiver clock relative to the GPS system time, It is the clock drift of the satellite clock relative to the GPS system time.

2. The improved satellite navigation direct position estimation method according to claim 1, characterized in that: The step S1 further includes: the observation signals of the K snapshots are expressed as follows: in: is the observation signal of K snapshots, is the composite signal amplitude matrix, is the delay-Doppler matrix, represents a zero mean within the observation interval and a variance of k snapshots of complex Gaussian white noise signal.

3. The improved satellite navigation direct position estimation method according to claim 2, characterized in that: The K snapshot observation signals ; The composite signal amplitude matrix , ; The delay-Doppler matrix , ; The observation interval has zero mean and variance k snapshots of complex Gaussian white noise signal , , , .

4. The improved satellite navigation direct position estimation method according to claim 3, characterized in that: The step S2 includes: By maximizing The likelihood function is equivalent to minimizing the nonlinear least squares problem , we get a multivariate optimization parameter estimation problem ,in .

5. The improved satellite navigation direct position estimation method according to claim 4, characterized in that: The step S3 includes: The variables to be estimated in the satellite signal multi-snapshot model are set as , the other four speed variables are set to 0. In order to simplify the algorithm search, the clock error is also set to 0; In the variable step size part, first determine the initial search range. The range is [-L, L], and the search step is After each update, the search range will also change accordingly with the search factor. The specific value is ,in, is the previous estimate, is the search boundary, ; In the Sage algorithm part, each step size corresponds to a Sage algorithm search, and the number of Sage iterations is set to N0. After all searches are completed, a rough estimate of the position is obtained.

6. The improved satellite navigation direct position estimation method according to claim 5, characterized in that: The step S4 comprises: First, the initial parameter estimate is given according to the result of step S3, and , calculate the cost function ; Second, minimize the cost function to update the current parameters ; Again , and judge whether all current variables have completed the search, if , it returns to search for the next variable, otherwise, it goes to the next step; Finally, the number of iterations , if all iterations are completed, the final estimated result is obtained .

Images