Simulation research on high precision multimode GNSS positioning algorithm

Abstract

The traditional Gauss-Newton iterative method is highly dependent on the initial value when locating the multimode GNSS receiver. If the difference between the initial value and the true value is higher, the algorithm has the problem of increasing number of iterations, and the algorithm lacks the self main monitoring process of the GNSS receiver, which leads to a great reduction in the positioning accuracy. A high precision multi-mode GNSS positioning algorithm is proposed. It is based on the composition and working principle of multimode GNSS multimode receiver, and the pseudo distance positioning distance is obtained by using GNSS multi constellation combined location algorithm. It uses a direct algorithm without initial value and iteration through the new algorithm of high precision positioning. After linearizing the pseudo range location distance equation, the user’s general position is calculated. After the pseudo range location distance equation is carried out in the general position of the user, the user’s position correction is calculated by weighted least squares, and the exact location of the user is obtained. The receiver autonomous integrity monitoring (RAIM) algorithm based on the least square residual method of GNSS receiver is used to realize the self-improvement monitoring of GNSS, and to further improve the precision of the multi-mode GNSS positioning algorithm. Experimental results show that the proposed location algorithm has high location accuracy and stability.

Keywords

High accuracy multimode GNSS positioning algorithm receiver pseudo range fusion RAIM algorithm perfection

1 Introduction

Since the steady progress of the construction of BDS and Galileoleo systems, the application of the multimode GNSS data processing method has become the research hot spots in all fields, especially in the satellite navigation system [1]. Its high precision, global, and real-time characteristics have been paid attention to by the satellite navigation system, which are widely used in the world, and opening a new era of navigation. The current research results showed that the integrated multimode GNSS system can effectively improve the reliability and stability of GNSS positioning, especially in areas such as mountain areas, urban canyons, and so on [2]. The traditional Gauss-Newton iterative method is highly dependent on the initial value in the positioning process of GNSS receiver. If the difference between initial value and true value is higher, the algorithm has a problem of increasing number of iterations, and the algorithm lacks the process of independent monitoring of the GNSS receiver, which leads to the problem of the multimode GNSS positioning precision greatly reduced. To solve this problem, a high-precision multimode GNSS location algorithm is proposed in this paper to achieve high-accuracy positioning of multi-mode GNSS receiver.

2 Simulation of High precision multimode GNSS Positioning algorithm

2.1 The composition and working principle of multimode GNSS receiver

The specific process of multimode GUSS receiver [4] is acquisition, tracking, bit synchronization, bit extraction, frame synchronization, extraction of satellite parameters, satellite position speed calculation, and user position speed calculation. First, the receiver makes up for the satellite signal. After tracking the satellite signal, the multimode GUSS receiver can be measured to settle according to the positioning settlement scheme [5], which can calculate the longitude, dimension, height, speed, time and other information of the user’s location. The composition and working principle of multimode GUSS receiver, receiver function module, hardware and software architecture are the necessary parameter indicators of multimode GUSS receiver.

The multimode GUSS receiver [6] receives satellite signals through the antenna unit. After the signal is amplified by low noise, the intermediate frequency signal is output through two frequency conversion, low pass filter, medium frequency amplification, and automatic gain control (AGC). The intermediate frequency signal is quantized by A/D, the signal is processed by A/D, pseudo distance measurement is obtained by solution expansion and conciliation [7]. Based on this, the location, speed, clock error and clock error rate of users were calculated, and the navigation business under high dynamic conditions was implemented. The block diagram of its composition is as shown in Fig. 1.

Fig.1

Structure diagram of multimode GNSS receiver.

The antenna element is composed of active antenna, passive antenna, acoustic meter filter and low noise amplifier (LNA) [8]. The satellite signals of different frequency points need an independent antenna. The function of the antenna is to transform the very weak electromagnetic wave energy of the satellite signal into the corresponding current. The sound meter filter is to filter out the external interference signal and suppress the frequency of the cross harmonic mirror, while the LNA is to amplify the signal and improve the signal to noise ratio.

The receiving and processing module is composed of a radio frequency receiving circuit and a baseband acquisition and processing circuit. The radio frequency receiving circuit is composed of filter, amplifier, down converter and wood vibration source. It is responsible for amplifying, filtering, down conversion to intermediate frequency signals received from the antenna, to ensure stable high gain and make the power of the intermediate frequency signal in the range required by the baseband acquisition. The baseband acquisition and processing circuit is composed of A/D, FPGA, memory, DSP and so on. It uses 8 bit A/D for intermediate frequency sampling, and the signal channel uses FPGA. Its function is to retrieve, pull and track the satellite, and provide the microprocessors with correlation values, pseudo range measurement values, carrier phase measurement values and Doppler shift measurement values. Memory mainly stores correlation value, pseudo range observation value, carrier phase observation value and Doppler shift. The microprocessor uses DSP chip, which is the control part of the receiver. The function is to self-detect, detect and search the satellite signal for the receiver, coordinate the channel to track the satellite signal, the bit synchronization, the frame synchronization, the despreading of the navigation message data, the demodulation of the navigation message, and the related operation [9]. (For example, the three-dimensional coordinates of the station, the navigation parameters, etc.). Receiving user input information, such as measuring station, station number, observer name, antenna, etc., and data communication at the same time. The control circuit is used for the interaction between the receiving processing module and the human interface. It displays some data and state through the liquid crystal display and the LED light emitting diode, receives the control signal of the chassis panel at the same time, sets up the receiving processing module and regulates the state [10]. The power module provides the power needed for the receiver based on the power interface requirements of the system.

2.2 GNSS multi constellation positioning algorithm

GNSS multi constellation combined positioning solution mode is to select satellites in three satellite navigation systems, including GPS, GLONASS and BD-2, and to locate satellites and locate [11]. GNSS multi constellation combined localization algorithm is divided into two ways: location result fusion and pseudo range fusion. In this study, the multimode GNSS receiver is solved by pseudo range fusion. The pseudo range fusion positioning solution is to select satellites in three satellite navigation systems of GPS, GLONASS and BD-2. The corresponding pseudorange of each satellite is unified spatio-temporal and with simultaneous equations.

The multi constellation combined location and solution mode signal is combined to solve [12]. As shown in Fig. 2, at this time, each receiving channel gives Doppler information, satellite launch time and navigation message information of the tracking satellite. The program is based on the satellite launching time and the ephemeris information to get the pseudo distance and the satellite position. After the implementation of the above parameters, the simultaneous pseudo range formulas can obtain the user’s position calculation results.

Fig.2

Pseudo range fusion flow chart.

The location principle of GNSS multi constellation combined location algorithm is basically the same with the positioning principle of single constellation location solution system (such as GPS system) [13]. Through the single constellation positioning algorithm, at least 4 satellite information are needed to calculate. When the two GNSS systems use multi constellation combination location solving model, the two modes of satellite system use different system time, so the unknown quantity increases the time difference between two time system and at least 5 satellite information are needed for location calculation; in the same way, in three GNSS system uses multi constellation combined location solution mode, at least 6 satellites are needed to solve the location problem. More visible satellites can improve positioning accuracy. Although the method of pseudo range fusion is relatively complex, the effect of information fusion is better. It can make better use of GDOP and RAIM technology for satellite selection and integrity detection [14].

2.2.1 Localization of gnss pseudo range fusion method

It is assumed that at some time, the multimode GNSS receiver receives a total of n visible satellites in GPS, GLONASS and BD-2 systems. By obtaining the pseudo range and navigation message of all the visible satellites from each receiving channel, the simultaneous formulas are obtained:

$\begin{matrix} ρ_{1} & = & \sqrt{(x_{1} - x_{u})^{2} + (y_{1} - y_{u})^{2} + (z_{1} - z_{u})^{2}} \\ + b_{u} + k_{1} δ_{b} + l_{1} δ_{g} \\ ρ_{2} & = & \sqrt{(x_{2} - x_{u})^{2} + (y_{2} - y_{u})^{2} + (z_{2} - z_{u})^{2}} \\ + b_{u} + k_{2} δ_{b} + l_{2} δ_{g} \\ ρ_{3} & = & \sqrt{(x_{3} - x_{u})^{2} + (y_{3} - y_{u})^{2} + (z_{3} - z_{u})^{2}} \\ + b_{u} + k_{3} δ_{b} + l_{3} δ_{g} \\ ρ_{n} & = & \sqrt{(x_{n} - x_{u})^{2} + (y_{n} - y_{u})^{2} + (z_{n} - z_{u})^{2}} \\ + b_{u} + k_{n} δ_{b} + l_{n} δ_{g} \end{matrix}$ (1)

Among them, ρ_i represents the pseudo (x_i, y_i, z_i) range of the satellite i, (x_i, y_i, z_i) represents the position of the satellite i, (x_u, y_u, z_u) represents the user coordinates, and b_u represents the distance deviation between the user receiver and the clock difference of the GPS system, b_u = ct_u; δ_b and δ_g represent the distance deviation caused by the clock difference between the BD-2 system and the GLONASS system, respectively. In the BD-2 system satellite; k_i and l_i represent that when the satellite i is GPS system satellite, it will take zero. When the satellite i belongs to the BD-2 system satellite, then k_i = 1, l_i = 0. When the satellite i belongs to the GL k_i = 0 ONASS system satellite, then k_i = 0, l_i = 1.

The formula (1) is written in a matrix form: $Δ ρ = H Δ x$ (2)

According to the least square estimation, it is obtained that: $Δ x = {[H^{T} H]}^{- 1} H^{T} Δ ρ$ (3)

The result of weighted least square estimation is: $Δ x = {[H^{T} H]}^{- 1} H^{T} P Δ ρ$ (4)

2.2.2 Pseudorange location formula

The pseudo range observations of the satellite are subject to multiple disturbances. The pseudo range of the GPS observed from the receiver with the satellite i can be described as follows:

$\begin{matrix} {\tilde{ρ}}_{i} & = & \sqrt{(x_{i} - x_{u})^{2} + (y_{i} - y_{u})^{2} + (z_{i} - z_{u})^{2}} \\ + Δ D_{i} + Δ b_{i} + Δ v_{i} + Δ I_{i} + Δ T_{i} + b_{u} + e_{i} \\ (i = 1, 2, \dots m, m \geq 4) \end{matrix}$ (5)

(x, y, z) ^T and (x_i, y_i, z_i) ^T are the WGS-84 coordinates of the receiver and the satellite i, respectively, ΔD_i is the influence of the satellite position error to the distance, Δb_i is the distance error caused by the satellite clock difference, Δv_i is the distance error caused by the relativistic time correction, ΔI_i and Δt_i are the ionosphere and troposphere delay error respectively, b_u is the distance difference caused by the receiver clock difference, e_i is the receiver measurement noise error. ΔD_i, Δb_i, Δv_i, ΔI_i and Δt_i can be corrected through ephemeris data. In the single frequency receiver, it can be corrected mainly according to ephemeris. It can improve the ΔI_i and Δt_i correction in dual frequency receiver through frequency delay characteristics. Assuming that the corrected positive term is completely corrected, the pseudo range distance formula is:

${\begin{matrix} ρ_{1} = \sqrt{(x_{1} - x_{u})^{2} + (y_{1} - y_{u})^{2} + (z_{1} - z_{u})^{2}} + b_{u} + e_{1} \\ . . . \\ ρ_{m} = \sqrt{(x_{m} - x_{u})^{2} + (y_{m} - y_{u})^{2} + (z_{1} - z_{u})^{2}} + b_{u} + e_{m} \end{matrix}$ (6)

The range error e_i caused by the receiver noise mainly depends on the tracking error of the code loop, and the tracking error of the code tracking loop obeys the Gauss distribution that the mean value is 0 and the variance is $σ_{i}^{2}$ , that is, the corresponding error $e_{i} \sim N (0, σ_{i}^{2})$ of the distance $σ_{i}^{2}$ of the satellite i, different satellites are different in distance and signal shielding. The difference of ranging error variance caused by signal arrival and receiving is different from that of $σ_{i}^{2}$ . The least square solution has the same variance and the least variance of the variance only in the task component. Because the error variance $σ_{i}^{2}$ of the satellite distance error is different, it is impossible to obtain the least solution of the location error by using the least square method in the positioning. The weighted least two multipliers should be adopted. It can be proved that the weighted least square solution can be proved when the weight of the component is the standard deviation of the corresponding error. It is optimal in the sense of minimum variance.

2.3 New algorithm for high precision positioning

The Gauss-Newton iterative method has local convergence [15], it has high dependence on initial value, and the convergence of Gauss-Newton is very fast when the initial value is close to true value. But when the difference between initial value and true value is high, the number of iterative times increases dramatically. The direct solution location method does not need initial value and iteration, the operation amount is small but the positioning precision is low. The new algorithm of high precision positioning combines the advantages of the direct solution algorithm and the high precision of the Gauss-Newton iterative method. Based on the GNSS multi constellation positioning algorithm, the initial positioning distance is obtained, and the precise location of the multimode GNSS is realized.

2.3.1 Initial location of direct solution algorithm

For the pseudo range location formula group (6), we get the square of the transfer item:

${\begin{matrix} (ρ_{1} - b_{u})^{2} \\ = (x_{1} - x_{u})^{2} + (y_{1} - y_{u})^{2} + (z_{1} - z_{u})^{2} + δ_{1} \\ . . . \\ (ρ_{m} - b_{u})^{2} \\ = (x_{m} - x_{u})^{2} + (y_{m} - y_{u})^{2} + (z_{m} - z_{u})^{2} + δ_{m} \end{matrix}$ (7)

Among them, $δ_{i} = 2 (ρ_{i} - b_{u}) e_{i} - e_{i}^{2} (i = 1, 2 \dots m)$ is the corresponding random error of the I equation.

The first formula in the formula (7) is subtracted from the two term of the unknown number, thus the purpose of linearization of the formula group is reached. $r_{i}^{2} = x_{i}^{2} + y_{i}^{2} + z_{i}^{2} (i = 1, 2 \dots m)$ , the formula after subtraction is collated:

${\begin{matrix} 2 (x_{2} - x_{1}) x_{u} + 2 (y_{2} - y_{1}) y_{u} + 2 (z_{2} - z_{1}) z_{u} - 2 (ρ_{2} - ρ_{1}) b_{u} \\ = (r_{2}^{2} - r_{1}^{2}) - (ρ_{2}^{2} - ρ_{1}^{2}) + (δ_{2} - δ_{1}) \\ ⋮ \\ 2 (x_{m} - x_{1}) x_{u} + 2 (y_{m} - y_{1}) y_{u} + 2 (z_{m} - z_{1}) z_{u} - 2 (ρ_{m} - ρ_{1}) b_{u} \\ = (r_{m}^{2} - r_{1}^{2}) - (ρ_{m}^{2} - ρ_{1}^{2}) + (δ_{m} - δ_{1}) \end{matrix}$ (8)

Put formula (8) into matrix form: ${Bx}_{0} = c + e_{1}$ (9)

Where: $\begin{matrix} B & = & [\begin{matrix} 2 (x_{1} - x_{2}) 2 (y_{1} - y_{2}) 2 (z_{1} - z_{2}) 2 (ρ_{1} - ρ_{2}) \\ ⋮ ⋮ ⋮ ⋮ \\ 2 (x_{1} - x_{m}) 2 (y_{1} - y_{m}) 2 (z_{1} - z_{m}) 2 (ρ_{1} - ρ_{m}) \end{matrix}] \\ = & [b_{1} b_{2} b_{3} b_{4}] \end{matrix}$ (10) $x_{0} = [\begin{matrix} x_{u} \\ y_{u} \\ z_{u} \\ b_{u} \end{matrix}], c = [\begin{matrix} (ρ_{2}^{2} - ρ_{1}^{2}) + (r_{1}^{2} - r_{2}^{2}) \\ ⋮ \\ (ρ_{m}^{2} - ρ_{1}^{2}) + (r_{1}^{2} - r_{m}^{2}) \end{matrix}], e_{1} = [\begin{matrix} δ_{2} - δ_{1} \\ ⋮ \\ δ_{m} - δ_{1} \end{matrix}]$

If there is m satellites, m - 1 formulas can be obtained from the formula (9). It is obvious that the position of the receiver can be solved directly when m > 4, and there are only three formulas at m = 4, so the formula (9) is divided into two cases:

When m > 4, the formula (9) coefficient matrix B is full rank, and the location of the receiver is to solve the linear formulas with e₁ perturbation. The variance of the e₁ component (δ_i+1 - δ₁) is different because of the different ranging errors of each satellite, and the linear formulas are weighted least squares, and both sides of formula (9) are multiplied by a weight matrix W₁: $W_{1} {Bx}_{0} = W_{1} c + W_{1} e_{1}$ (11)

The weighted least squares solution is as follows: $x_{0} = {[{(W_{1} B)}^{T} W_{1} B]}^{- 1} {(W_{1} B)}^{T} W_{1} c$ (12)

Where:

$W_{1} = diag [\frac{1}{\sqrt{D (δ_{2} - δ_{1})}} \frac{1}{\sqrt{D (δ_{3} - δ_{1})}} \dots \frac{1}{\sqrt{D (δ_{m} - δ_{1})}}]$ (13)

D (·) represents the variance. From the formula (13), we can see that the weight of the first i (i = 1, 2 ⋯ m - 1) formula is the reciprocal (δ_i+1 - δ₁) of the corresponding standard error of ${(\sqrt{D (δ_{i + 1} - δ_{1})})}^{- 1}$ . So the key to the solution of the minimum variance solution is to seek ${(\sqrt{D (δ_{i + 1} - δ_{1})})}^{- 1}$ . An estimation method for ${(\sqrt{D (δ_{i + 1} - δ_{1})})}^{- 1}$ is given below.

Because δ_i is independent of each other and b_u is small relative to ρ_i, so:

$\begin{matrix} D (δ_{i + 1} - δ_{1}) = D (δ_{i + 1}) + D (δ_{i}) \\ = D [2 (ρ_{i} - b_{u}) e_{i} - e_{i}^{2}] \\ + D [2 (ρ_{1} - b_{u}) e_{1} - e_{1}^{2}] \\ \approx D (2 ρ_{i} e_{i} - e_{i}^{2}) + D (2 ρ_{1} e_{1} - e_{1}^{2}) \end{matrix}$ (14)

Then:

$\begin{matrix} D (2 ρ_{i} e_{i} - e_{i}^{2}) \\ = E {(2 ρ_{i} e_{i} - e_{i}^{2})}^{2} - {[E (2 ρ_{i} e_{i} - e_{i}^{2})]}^{2} \\ = E (4 ρ_{i}^{2} e_{i}^{2} + e_{i}^{4}) - {[E (e_{i}^{2})]}^{2} = 2 δ_{i}^{4} + 4 ρ_{i}^{2} δ_{i}^{2} \end{matrix}$ (15)

So:

$\begin{matrix} {(\sqrt{D (δ_{i + 1} - δ_{1})})}^{- 1} \\ = {(\sqrt{(2 δ_{i + 1}^{4} + 4 ρ_{i + 1}^{2} δ_{i + 1}^{2}) + (2 δ_{1}^{4} + 4 ρ_{1}^{2} δ_{1}^{2})})}^{- 1} \\ (i = 1, 2 \dots m - 1) \end{matrix}$ (16)

Handle formula (16) right into the matrix (13), then the weighting matrix into formula (12) to obtain the initial value of the position of the receiver.

When m = 4, the first formula of the simultaneous (7) and formula (7) (no redundant information is not weighted at this time, so the writing formula does not consider the random error term of the formula): ${\begin{matrix} x_{u}^{2} + y_{u}^{2} + z_{u}^{2} - b_{u}^{2} - 2 x_{1} x_{u} - 2 y_{1} y_{u} - 2 z_{1} z_{u} \\ + 2 ρ_{1} b_{u} = ρ_{1}^{2} - r_{1}^{2} \\ {Bx}_{0} = c \end{matrix}$ (17)

First, in the formula group Bx₀ = c, the b_u is known as the known number, x_u, y_u and z_u are solved, and then x_u, y_u and z_u are substituted into formula (17) to solve b_u. Finally, the form of the solution can be obtained: ${\begin{matrix} b_{u} = (- b \pm \sqrt{b^{2} - 4 ac}) / - 2 a \\ x_{u} = a_{1} + b_{1} b_{u} \\ y_{u} = a_{2} + b_{2} b_{u} \\ z_{u} = a_{3} + b_{3} b_{u} \end{matrix}$ (18)

Where: ${\begin{matrix} a = b_{1}^{2} + b_{2}^{2} + b_{3}^{2} - 1 \\ b = 2 a_{1} b_{1} + 2 a_{2} b_{2} + 2 a_{3} b_{3} - 2 x_{1} b_{1} \\ - 2 y_{1} b_{2} - 2 z_{1} b_{3} + 2 ρ_{1} \\ c = a_{1}^{2} + a_{2}^{2} + a_{3}^{2} - 2 x_{1} a_{1} - 2 y_{1} a_{2} - 2 z_{1} a_{3} \\ - ρ_{1}^{2} + r_{1}^{2} \end{matrix}$

$\begin{matrix} {[a_{1} a_{2} a_{3}]}^{T} = {[b_{1} b_{2} b_{3}]}^{- 1} c \\ {[b_{1} b_{2} b_{3}]}^{T} = {[b_{1} b_{2} b_{3}]}^{- 1} b_{4} \end{matrix}$ (19)

In the formula (19), we can see that the time error b_u of the receiver has two solutions, resulting in the user’s position also having two sets of solutions. It is known from the geometry knowledge that the two sets of solutions are about the plane symmetry of the three satellites. Obviously, the group near the earth is the location of the receiver on the ground, that is, the smaller group of $\sqrt{x_{u}^{2} + z_{u}^{2} + y_{u}^{2}}$ solutions.

2.3.2 Precise positioning

The first order Taylor expansion of the distance formulas of formula (6) is carried out at the x₀ of formula (11) or formula (18). $A δ x = δ ρ + e_{2}$ (20)

Where, $A = [\begin{matrix} \frac{x_{1} - x_{u}}{ρ_{1} - ρ_{u}} \frac{y_{1} - y_{u}}{ρ_{1} - b_{u}} \frac{z_{1} - z_{u}}{ρ_{1} - b_{u}} 1 \\ ⋮ ⋮ ⋮ ⋮ \\ \frac{x_{m} - x_{u}}{ρ_{m} - ρ_{u}} \frac{y_{m} - y_{u}}{ρ_{m} - b_{u}} \frac{z_{m} - z_{u}}{ρ_{m} - b_{u}} 1 \end{matrix}]$ (21) $δ ρ = {[δ ρ_{1} δ ρ_{2} \dots δ ρ_{n}]}^{'}$ (22) $δ x = {[δ x_{u} δ y_{u} δ z_{u} δ b_{u}]}^{'}$ (23) $e_{2} = {[e_{1} e_{2} \dots e_{m}]}^{'}$ (24)

The pseudo moment ρ and δρ calculated from x₀ and satellite positions are calculated from the difference between observed pseudo moment and ρ. Because the distance error of each satellite is different, the variance of each component e₂ of the error term e_i of formula (20) is different (i = 1, 2 ⋯ m), and the same is solved by the weighted least squares, and the formula (20) is multiplied by the weighted matrix W₂: $W_{2} A δ x = W_{2} δ ρ + W_{2} e_{2}$ (25)

The weighted least squares solution is obtained as follows: $δ x = {[{(W_{2} A)}^{T} W_{2} A]}^{- 1} {(W_{2} A)}^{T} W_{2} δ ρ$ (26)

Where, $W_{2} = diag [\frac{1}{δ_{1}} \frac{1}{δ_{2}} \dots \frac{1}{δ_{m}}]$ (27)

The weight of the first i (i = 1, 2… m) equation is the reciprocal $σ_{i}^{- 1}$ corresponding to the standard deviation of the error e_i. The smaller the ranging error is, the more reliable the satellite is, the greater the weight is, so as to reduce the location error. The user location correction volume δx is obtained from formula (26): $x = x_{0} + δ x$ (28)

The process of obtaining the new algorithm is to use the direct solution algorithm to calculate the user’s general position x₀, to carry out Taylor expansion at the distance equation group at x₀, to calculate the user position correction δx by weighted least squares, and to obtain the user location x = x₀ + δx.

3 GNSS receiver autonomous integrity monitoring (RAIM) algorithm based on least squares residual method

Perfect monitoring is one of the most important technical problems in the development of GNSS receiver, and its function is to provide alarm information for users who cannot be used in navigation service time [16]. Specifically, the RAIM algorithm is used to monitor the normal positioning results of the GNSS receiver, and to find and eliminate the data of the fault satellite when it is not working properly. It is mainly divided into two steps, the first step is the error detection, if the measurement data exceeds a certain threshold, the error of the task data is stored and the obstacle stars exist; the second step is to determine and eliminate the fault satellite, if the fault star is found, it is found according to certain rules. At present, the GNSS receiver’s self-improvement monitoring is usually based on the least square residual [17]. This study analyzes and studies the receiver autonomous integrity monitoring algorithm of GNSS receiver based on the least squares residual method, and explains the problems in the application of the algorithm to the GPS/GLONASS/BD-2 multimode GNSS receiver.

3.1 The basic principle of the least squares residual method

The observation equations of the GNSS system are as follows: $y = G \cdot X + ɛ$ (29)

y ∈ Rⁿ⁺¹ indicates the difference between the measured pseudo range and the estimated value; X ∈ R⁴⁺¹ indicates the true position of the user receiver and the clock difference; ɛ ∈ Rⁿ⁺¹ represents the pseudo distance measurement error of the satellites; and G ∈ Rⁿ⁺⁴ is a linear matrix of the direction cosine vector process from each satellite to the user receiver.

Tracking least squares algorithm, the least square solution vector X is: ${\overset{⌢}{X}}_{LS} = (G^{T} G)^{- 1} \cdot G^{T} \cdot y \leq$ (30)

Then: $\overset{⌢}{y} = G \cdot {\overset{⌢}{X}}_{LS} = G \cdot (G^{T} \cdot G)^{- 1} \cdot G^{T} \cdot y$ (31)

With the description of the residual vector, there are:

$\begin{matrix} w & = & y - \overset{⌢}{y} = y - G \cdot {\overset{⌢}{X}}_{LS} \\ = & [I - G \cdot (G^{T} \cdot G)^{- 1} G^{T}] \cdot y \\ [I - G \cdot (G^{T} \cdot G)^{- 1} G^{T}] \cdot [G \cdot X + ɛ] \\ = & [G - G \cdot (G^{T} \cdot G)^{- 1} G^{T} G] \cdot X \\ + [I - G \cdot (G^{T} \cdot G)^{- 1} G^{T}] \cdot ɛ \\ = & [I - G \cdot (G^{T} \cdot G)^{- 1} G^{T}] \cdot ɛ = S \cdot ɛ \end{matrix}$ (32)

The analysis formula (32) is available. Vector w contains satellite ranging error information, which can be used as a basis for judging whether there is a fault star. $SSE = w^{T} w = u_{1}^{2} + u_{2}^{2} + . . . + u_{u + 4}^{2}$ is defined as the sum of squares of the satellite pseudo range residuals, which is a χ² distribution of a degree of freedom (n-4). In practical applications, it usually uses $T_{s} = \sqrt{SSE / (n - 4)}$ as a statistical test, and compares it with the threshold value T₀ at each location. If T_s ≥ T₀, there is a fault star. If T_s < T₀, there is no fault star. When the statistical detection volume exceeds the threshold value, it indicates the presence of the failure star. At this time, we can find the most likely satellite with the greatest likelihood of failure based on the maximum likelihood estimation method [18]. Because w = S · y, if the $w_{i}^{2} / s_{n}$ value of a satellite reaches the maximum, the mission of the satellite is most likely to fail.

3.2 Operation flow of raim algorithm based on least squares residuals

Based on the least squares residuals of the RAIM algorithm, the operation process is described in Fig. 3.

Fig.3

RAIM algorithm.

The description and analysis of the GNSS receiver RAIM algorithm showed that the algorithm is characterized by self-improvement monitoring through the current time data; the proposed algorithm performs error detection and error correction on the assumption that only one satellite fails at the same time. If a navigation system has more than one fault star at the same time, then RAIM can only give warning and cannot put forward the fault star. Under the single constellation positioning solution mode, 5 satellites can detect the fault star, and more than 6 satellites can be observed to put forward the fault star, and the more satellite number, the better error correction performance is [19].

The RAIM algorithm of the GPS/GLONASS/BD-2 multimode GNSS receiver is basically the same as the single constellation GNSS system, but due to the equivalent distance error of the GLONASS satellite and the slight difference between the GPS system, in the operation detection threshold, it should take the system of equivalent distance error to carry out the operation [20], which can reduce the false alarm probability of the system. At this time, there are 6 unknowns (X, Y, Z, Δt₁, Δt₂, Δt₃) that 6 satellites are needed to locate. In this way, the degree of freedom of the system is changed from n-4 to n-6 under a single system. Therefore, in the RAIM algorithm of GPS/GLONASS/BD-2 combination location, it needs 7 satellites to detect, and 8 satellites to correct errors.

4 Experimental analysis

4.1 Comparison of three GNSS location algorithms

The Gauss-Newton iterative method needs to be iterated many times and each iteration also revise the position of the satellite in addition to the update of the receiver position. The tropospheric ionosphere delay also needs to be reestimated according to the latest pitching angle. It involves a large number of matrix multiplications, matrix inversion and the calculation of the mathematical formula. The number of iterations has great influence on the computation of Gauss-Newton iteration. The Gauss-Newton iteration location algorithm, the two step location algorithm and the computation amount of the proposed algorithm are statistically analyzed. The results are described in Table 1. m is the number of satellites, and k is the number of iterations.

Table 1
Computational complexity of three GNSS positioning algorithms

Algorithm Multiplicative times Addition times Number of division Open root number trigonometric function

Gauss Newton iterative algorithm (50m + 60) + (346m + 111) (k – 1) (39m + 16) + (137m + 76) (k – 1) m + (39m + 15) (k-) m + (5m + 7) (k-1) (8m + 3) (k-1)

Two step location algorithm 57m + 16 44m-26 m-1 0 0

Algorithm in this paper 403m + 127 181 + 56 40m + 14 5m + 7 8m + 3

Algorithm	Multiplicative times	Addition times	Number of division	Open root number	trigonometric function
Gauss Newton iterative algorithm	(50m + 60) + (346m + 111) (k – 1)	(39m + 16) + (137m + 76) (k – 1)	m + (39m + 15) (k-)	m + (5m + 7) (k-1)	(8m + 3) (k-1)
Two step location algorithm	57m + 16	44m-26	m-1	0	0
Algorithm in this paper	403m + 127	181 + 56	40m + 14	5m + 7	8m + 3

4.2 Simulation results and analysis

In the simulation experiment, the initial value of the Gauss-Newton iterative algorithm uses the method of literature [7]: the initial value is set as (0, 0, 0, 0), the location of the user is 40 degrees, the longitude is 120 degrees, the height is 100 meters, the satellite is set up of 7, and the GDOP is 2.109. Table 2 gives the positioning error of the three location algorithms, and the number of trigonometric functions represents the number of times for all kinds of trigonometric and inverse trigonometry. From Table 2, we can see that the computation amount of the algorithm is about 1/20 of the Gauss-Newton iterative algorithm. The two step localization algorithm is about 1/3 of the Gauss-Newton iterative method. It can be seen that the computation of this algorithm is greatly reduced and the location efficiency is high.

Table 2
Computational complexity of three GNSS positioning algorithms

Algorithm Multiplicative times Addition times Number of division Root number Trigonometric function times

Gauss Newton iterative algorithm 8009 3394 871 133 177

Two step location algorithm 415 282 6 0 0

Algorithm in this paper 2948 1323 294 42 59

Algorithm	Multiplicative times	Addition times	Number of division	Root number	Trigonometric function times
Gauss Newton iterative algorithm	8009	3394	871	133	177
Two step location algorithm	415	282	6	0	0
Algorithm in this paper	2948	1323	294	42	59

The Gauss-Newton iterative algorithm, the two step location algorithm and the proposed algorithm are used to locate the data above 10 000 times. The statistical results of positioning accuracy are shown in Table 3. From the statistical data of Table 3, the positioning accuracy of the two step positioning algorithm is about 15 m, the positioning precision of Gauss-Newton is 2.44 m, and the positioning accuracy of the two step positioning algorithm is nearly 7 times worse than the Gauss-Newton iterative method. The location accuracy of the proposed algorithm is 1.65 m, which is higher than the other two algorithms, which indicates that the proposed algorithm has high location accuracy.

Table 3

Positioning precision of three positioning algorithms

Positioning algorithm/Three dimensional	2DRMS
Gauss Newton iterative algorithm	2.44 rice
Two step location algorithm	15.84 rice
Algorithm in this paper	1.65 rice

4.3 Accuracy evaluation of positioning result for high dynamic gnss receiver

In this study, the accuracy of the algorithm for a high dynamic GNSS receiver is tested, and the stability and accuracy of the algorithm are verified. The detailed results are as follows:

4.3.1 Evaluation of level classification accuracy

According to the fitting coefficient of the horizontal component and the horizontal coordinates of each ten meter pile through each ten meter pile, the horizontal component precision in the high dynamic environment can be compared with the known horizontal coordinates, as shown in Fig. 4.

Fig.4

Schematic diagram of the deviation between the fitting values of X and Y components and the known values.

Most of the deviations of horizontal components X and Y are within 1 cm. The RMS of X components is 4 mm, the RMS of Y components is 5 mm, and the RMS of horizontal components is 7 mm.

4.3.2 Accuracy evaluation of elevation component

The experiment selects the average value of the height observation value of the last two epochs of each ten meter pile, and compares the accuracy of the GNSS horizontal component in the high dynamic environment by comparing with the known elevation value.

As shown in Fig. 5, the result of the proposed algorithm to locate GNSS in high dynamic environment is that the deviation of the H component is within 5 cm and RMS is 2.9 cm, the accuracy of the proposed algorithm for high dynamic GNSS positioning is higher, and the higher stability of the proposed algorithm is also explained in the other aspect.

Fig.5

Schematic of the deviation between the fitting value of H component and the known values.

4.4 Comparison between test results and traditional commercial software calculation results

The final solution results of this experiment (Fig. 5) are compared with the results of the mature commercial software GrafNav. The concrete method is that the time system of the whole test is first unified into GPS time, then the data of this test is calculated by GrafNav, and the results of the experiment are analyzed in two ways:

Comparison and analysis of the results of location calculation with the algorithm in this paper, the contrast result is shown in Fig. 6.

The analysis of the standard error (standsrd deviation, SD) of GrafNav software and the “true value difference” (the difference between the results of the proposed algorithm and the true value) were analyzed.

Fig.6

CrafNav software solution value and self research software.

As shown in Fig. 6, the results showed that the results of the commercial software CrafNav calculation are in centimeter level and the difference between the estimation error and the actual calculation error is within the centimeter level, the difference between X direction and Y direction is within 3 cm, and the H direction is basically within 5 cm, only a few epochs have gross error. Through comparison, we can see that the high precision and multi-mode GNSS localization algorithm proposed in this paper has high positioning accuracy.

5 Conclusion

In this study, a high precision multi-mode GNSS location algorithm is proposed, which solves the drawbacks of the traditional Gauss Newton iterative algorithm for the high dependence of the initial value and the independent monitoring of the GNSS in the GNSS receiver positioning process, and improves the positioning accuracy of the multimode GNSS receiver. The proposed algorithm obtains the pseudo range location distance through the GNSS multi constellation positioning algorithm, and obtains the accurate location of the user through the new algorithm based on the distance. It adopts the GNSS receiver autonomous integrity monitoring (RAIM) algorithm based on the least square residual method to realize the autonomous and perfect monitoring of GNSS, which greatly enhances the positioning accuracy of the multimode GNSS receiver.

References

Garcia

J.G.

, Roncagliolo

P.A.

and Muravchik

C.H.

, A bayesian technique for real and integer parameters estimation in linear models and its application to GNSS high precision positioning, IEEE Transactions on Signal Processing64(4) (2016), 923–933.

Dong

, Wang

, Chen

et al., Mitigation of multipath effect in GNSS short baseline positioning by the multipath hemispherical map, Journal of Geodesy90(3) (2016), 255–262.

Chen

, Zhang

, Wang

et al., A simplified and unified model of multi-GNSS precise point positioning, Advances in Space Research55(1) (2015), 125–134.

Shen

, Li

, Xu

et al., Spatiotemporal filtering of regional GNSS network’s position time series with missing data using principle component analysis, Journal of Geodesy88(1) (2014), 1–12.

Huang

, Yao

, Cui

et al., Dilution of precision analysis for GNSS collaborative positioning, IEEE Transactions on Vehicular Technology65(5) (2016), 3401–3415.

, Li

X.C.

and Zhang

R.F.

, Research on algorithm of time delay estimation in satellite positioning systems, Computer Simulation32(4) (2015), 61–65.

, Ge

, Dai

et al., Accuracy and reliability of multi-GNSS real-time precise positioning: GPS, GLONASS, BeiDou, and Galileo, Journal of Geodesy89(6) (2015), 607–635.

Paziewski

and Wielgosz

, Investigation of some selected strategies for multi-GNSS instantaneous RTK positioning, Advances in Space Research59(1) (2017), 12–23.

Shu

, Shi

, Xu

et al., Error analysis of high-rate GNSS precise point positioning for seismic wave measurement, Advances in Space Research59(11) (2017), 2691–2713.

10.

Liu

, Xu

, Shu

et al., A new method to improve the performance of multi-GNSS pseudorange positioning in signal-degraded environment, Advances in Space Research58(4) (2016), 577–586.

11.

, Ding

, Dai

et al., Systematic error mitigation in multi-GNSS positioning based on semiparametric estimation, Journal of Geodesy91(12) (2017), 1491–1502.

12.

Liang

, Walter

, Enge

et al., GNSS multipath and jamming mitigation using high-mask-angle antennas and multiple constellations, IEEE Transactions on Intelligent Transportation Systems16(2) (2015), 741–750.

13.

Kashcheyev

, Nava

and Radicella

S.M.

, Estimation of higher-order ionospheric errors in GNSS positioning using a realistic 3-D electron density model, Radio Science47(4) (2016), 1–7.

14.

Chen

, Huang

, Liu

et al., Feasibility analysis of GNSS multi-constellation positioning for lunar spacecraft, Advances in Space Research32(6) (2015), 1155–1159.

15.

Gorni

and Zampieri

, Nonlocal and nonvariational extensions of killing-type equations, Discrete and Continuous Dynamical Systems-Series S11(4SI) (2018), 675–689.

16.

Cheng

, He

and Hu

, Quantitative risk analysis method of information security-Combining fuzzy comprehensive analysis with information entropy, Journal of Discrete Mathematical Sciences and Cryptography20(1) (2017), 149–165.

17.

Xue

, Jiang

and Sun

, Modeling of safety evaluation of super-large deep-water group pile foundation, Journal of Interdisciplinary Mathematics20(1) (2017), 101–111.

18.

Butt

, Ahmad

S.S.

, Shabbir

and Erum

, GIS based surveillance of road traffic accidents (RTA) risk for Rawalpindi city: A geo-statistical approach, Kuwait Journal of Science44(4) (2017), 129–134.

19.

Gao

, Farahani

M.R.

, Aslam

and Hosamani

, Distance learning techniques for ontology similarity measuring and ontology mapping, Cluster Computing-The Journal of Networks Software Tools and Applications20(2SI) (2017), 959–968.

20.

, Liu

, Furuta

and Peng

, Adsorption characteristics of sulfur solution by acticarbon against drinking-water toxicosis, Saudi Journal of Biological Sciences24(6) (2017), 1355–1360.