A method for solving the original observation equation of GNSS observation station network

Abstract

The invention discloses a method for solving the original observation equation of a GNSS observation station network, and belongs to the field of satellite navigation technology. The method is based on the multi-frequency pseudorange and phase observation data of the global or regional GNSS tracking station network, and constructs a full-rank non-difference non-combined full parameter estimation model by adding constraints, thereby realizing the unified estimation of positioning enhancement parameters such as receiver clock error, satellite clock error, tropospheric delay, ionospheric delay and base frequency satellite phase deviation; by giving weights to the parameters to be estimated, the benchmark for parameter estimation can be flexibly selected and adjusted, and the problem of selecting a single benchmark for phase deviation estimation when solving the original non-difference non-combined observation equation is effectively solved. The parameters estimated by this method can be used as PPP-RTK high-precision positioning enhancement services, and can also be used for the service needs of traditional positioning modes such as PPP and RTK.

Description

Technical Field

[0001] The invention relates to a method for solving original observation equations of a GNSS observation station network, belongs to the technical field of satellite navigation, and can be used for solving original (non-difference non-combination) observation equations of a GNSS observation station network. Background Art

[0002] The Global Navigation Satellite System (GNSS) was established to provide outdoor positioning, navigation and timing services. Today, real-time and high-precision satellite navigation and positioning technologies include: real-time differential kinematic positioning (RTK), real-time precise point positioning (Real-TimePPP and Real-Time PPP-AR), real-time kinematic positioning based on a reference station network (network RTK), and PPP-RTK technology. The realization of real-time and high-precision positioning technology not only depends on the corresponding algorithms, but also requires the support of a variety of enhanced information products, such as satellite clock error, satellite phase deviation, and atmospheric delay information. The enhanced information mainly comes from the solution of the GNSS observation station network.

[0003] However, when solving these enhanced information using the original observation equations, existing methods mainly rely on selecting the phase deviation corresponding to a single receiver or a single satellite as the benchmark for all estimations, and fail to use all phase deviation parameters as the benchmark for phase deviation estimation. Summary of the invention

[0004] In view of this, the present invention provides a method for solving the original observation equation of a GNSS observation station network, which can solve the GNSS original observation equation to provide the satellite clock error, satellite phase deviation and ionosphere, troposphere and other enhanced information required for precise navigation and positioning services such as PPP-AR and PPP-RTK.

[0005] In order to achieve the above object, the technical solution adopted by the present invention is:

[0006] A method for solving the original observation equation of a GNSS observation station network comprises the following steps:

[0007] (1) Using a single GNSS dual-frequency or multi-frequency raw observation data, construct the non-differenced and non-combined observation equation:

[0008]

[0009] Where, the superscript s represents the satellite number, the subscript r represents the receiver number, the subscript f represents the frequency band number, e represents the epoch number, P represents the pseudorange observation value in meters, and L represents the carrier phase observation value in meters. is the geometric distance from the receiver r to the satellite s, which is a known value in meters. is the model value of the tropospheric projection function from receiver r to satellite s, Tr is the tropospheric delay in the zenith direction, in meters, dt r is the clock error of the receiver r, in meters, dt s is the clock error of satellite s, in meters, b r,f is the pseudorange hardware delay of receiver r corresponding to frequency f, in meters, is the pseudorange hardware delay of satellite s corresponding to frequency f, in meters, φ r,f is the phase hardware delay of the receiver r corresponding to the frequency f, in meters, is the phase hardware delay of satellite s corresponding to frequency f, in meters, is the ionospheric delay of the first frequency of the receiver r and satellite s, in meters, is the ratio of the ionospheric delay of other frequencies to the first frequency, λ f is the wavelength corresponding to the frequency f, in meters, is the integer ambiguity corresponding to frequency f between receiver r and satellite s, in meters;

[0010] (2) Add constraint type 1: Where Δdt i is the receiver clock correction value of the i-th station, Δdt i is the correction value of the i-th satellite clock error, Pdt i is the weight of the receiver clock correction value of the i-th station, Pdt i is the weight of the i-th satellite clock correction value, nr is the number of stations, and ns is the number of satellites;

[0011] (3) Add constraint type 2: Select the first and second frequencies and add constraint b to each station and satellite. r，f =0 and Where f takes two selected frequencies;

[0012] (4) Add constraint type 3: Add constraints for other frequencies other than the frequency selected in constraint 2. Where f is a frequency other than the frequency selected in constraint 2, Pb i，f is the weight of the receiver pseudorange hardware delay of the i-th station, is the weight of the hardware delay of the ith satellite pseudorange, nr is the number of stations, and ns is the number of satellites;

[0013] (5) Add constraint type 4: Add conditions for each frequency where Pφ i，f is the weight of the receiver phase hardware delay of the i-th station, is the weight of the hardware delay of the i-th satellite phase, nr is the number of stations, and ns is the number of satellites;

[0014] (6) Add constraint type 5: Add conditions for each frequency at each station. Where ns is the number of satellites selected as constraints;

[0015] (7) Add constraint type 6: Add conditions for each frequency of each satellite Where nr is the number of stations selected as constraints;

[0016] (8) Combine the observation equation and the six types of constraints to solve the parameter T r , The specific form of the parameters is as follows:

[0017]

[0018] in, α=1 / (μ2-1), β=μ2 / (μ2-1), γ f =α·(μ f +1), κ f =α·(μ f +μ2),

[0019] Compared with the prior art, the present invention has the following advantages:

[0020] (1) Epoch parameter T r , and They are independent of each other and can be solved in a single epoch;

[0021] (2) Clock difference standard dt base It can be selected as a combination of multiple stations and multiple satellites, and the weight of each station and satellite can be set by setting the Pdt of each station or satellite. r and Pdt s Enable flexible selection of benchmarks;

[0022] (3) Phase deviation reference φ base，f It can be selected as a combination of multiple stations and multiple satellites, and the weight of each station and satellite can be adjusted by setting the Pφ of each station or satellite. r，f and Enable flexible selection of benchmarks;

[0023] (4) First-class ambiguity benchmark and It can be selected as a combination of multiple stations and multiple satellites, and the weight Pφ can be adjusted r，f and As well as the participating stations and satellites in constraint types 5 and 6, flexible benchmark selection is achieved;

[0024] (5) The second type of ambiguity criterion It can be selected as a combination of multiple stations and multiple satellites, and the weight Pφ can be adjusted r，f and As well as the stations and satellites participating in constraint types 5 and 6, flexible selection of benchmarks is achieved.

[0025] In summary, the present invention provides a method for parameter solution by applying constraints using observation data from GNSS global or regional tracking stations, and realizes synchronous estimation of all parameters to be estimated using a non-differenced and non-combined observation model. The benchmark for parameter estimation can be flexibly selected by adjusting the weights between the parameters to be estimated, ensuring that the estimated parameters have good consistency, self-consistency and stability, while providing data support for accurate modeling of the global and regional ionosphere and troposphere. DETAILED DESCRIPTION

[0026] A method for solving the original observation equation of a GNSS observation station network comprises the following steps:

[0027] (1) Pseudorange and phase observation data from a single GNSS system global or regional tracking station network are selected, including BeiDou-2 triple-frequency, BeiDou-3 triple-frequency, GPS dual-frequency or triple-frequency, GLONASS dual-frequency, and Galileo dual-frequency, triple-frequency, or quad-frequency observation data.

[0028] (2) Using the pseudorange and phase observations of all station systems and frequencies in the tracking station network, the non-differenced and non-combined original observation equations are constructed. The known satellite positions and station coordinates are introduced into the observation equations, that is, the geometric distance from the satellite to the receiver is known, and the parameters to be estimated include: tropospheric delay, receiver clock error, satellite clock error, receiver pseudorange hardware delay, satellite pseudorange hardware delay, receiver phase hardware delay, satellite phase hardware delay, ionospheric delay, and integer ambiguity.

[0029] (3) Construct the original observation equation of single system with multiple frequencies:

[0030]

[0031] Where, the superscript s represents the satellite number, the subscript r represents the receiver number, the subscript f represents the frequency band number, e represents the epoch number, P represents the pseudorange observation value in meters, and L represents the carrier phase observation value in meters. is the geometric distance from the receiver r to the satellite s, which is a known value in meters. is the model value of the tropospheric projection function from receiver r to satellite s, T r is the tropospheric delay in the zenith direction, in meters, dt r is the clock error of the receiver r, in meters, dt s is the clock error of satellite s, in meters, b r，f is the pseudorange hardware delay of receiver r corresponding to frequency f, in meters, is the pseudorange hardware delay of satellite s corresponding to frequency f, in meters, φ r，f is the phase hardware delay of the receiver r corresponding to the frequency f, in meters, is the phase hardware delay of satellite s corresponding to frequency f, in meters, is the ionospheric delay of the first frequency of satellite s at station r, in meters, is the ratio of the ionospheric delay of other frequencies to the first frequency, λ f is the wavelength corresponding to the frequency f, in meters, is the integer ambiguity corresponding to frequency f between receiver r and satellite s, in meters.

[0032] (4) Add constraint type 1: Where Δdt i is the receiver clock correction value of the i-th station, Δdt i is the correction value of the i-th satellite clock error, Pdt i is the weight of the receiver clock correction value of the i-th station, Pdt i is the weight of the i-th satellite clock correction value, nr is the number of stations, and ns is the number of satellites.

[0033] (5) Add constraint type 2: Select the first and second frequencies and add constraint b to each station and satellite. r，f =0 and Where f takes two selected frequencies.

[0034] (6) Add constraint type 3: Add constraints for other frequencies other than the frequency selected in constraint 2. Where f is a frequency other than the frequency selected in constraint 2, Pb i，f is the weight of the receiver pseudorange hardware delay of the i-th station, is the weight of the hardware delay of the ith satellite pseudorange, nr is the number of stations, and ns is the number of satellites.

[0035] (7) Add constraint type 4: Add conditions for each frequency where Pφ i，f is the weight of the receiver phase hardware delay of the i-th station, is the weight of the i-th satellite phase hardware delay, nr is the number of stations, and ns is the number of satellites.

[0036] (8) Add constraint type 5: Add conditions for each frequency of each station. Where ns is the number of satellites selected as constraints.

[0037] (9) Add constraint type 6: Add conditions for each frequency of each satellite Where nr is the number of stations selected as constraints.

[0038] (10) Combining the observation equations and the six types of constraints, the column rank deficiency problem in the equations is eliminated, except that the zenith tropospheric delay is still T r In addition, the specific forms of other parameters obtained are as follows:

[0039]

[0040] in, α=1 / (μ2-1), β=μ2 / (μ2-1), γ f =α·(μ f +1), κ f =α·(μ f +μ2),

[0041] In summary, the present invention aims to solve the problem of a single phase deviation estimation benchmark when solving the multi-frequency original observation equations of the GNSS observation station network to provide enhanced information. By adding constraints on the clock error, pseudorange hardware delay and phase deviation, the solution of the GNSS original (non-difference non-combined) observation equations is achieved.

[0042] The present invention is based on the pseudorange and phase observation data of multi-frequency points of the global or regional GNSS tracking station network, and constructs a full-rank non-difference non-combined full parameter estimation model by adding constraints, eliminating the model rank deficiency of the observation equation, maximizing the utilization rate of the observation data, and simultaneously realizing the real-time optimal estimation of multiple state domain parameters including receiver clock error, satellite clock error, phase deviation, pseudorange deviation, atmospheric ionosphere and troposphere; the weights between the parameters to be estimated can be adjusted to flexibly select the benchmark for parameter estimation, effectively solving the problem of selecting a single phase deviation estimation benchmark when solving the original non-difference non-combined observation equation. It can also solve the observation equations of dual frequencies and above, generate pseudorange deviations and phase deviations of the base frequency, and generate satellite clock errors and receiver clock errors corresponding to the ionosphere elimination combination. The method solves the problems of multiple state domain parameters such as satellite phase deviations, ionosphere and troposphere required by PPP-AR and PPP-RTK technologies, and the estimated parameters can be used as PPP-RTK high-precision positioning enhancement services, and can also be used for the service needs of traditional positioning modes such as PPP and RTK.

Claims

1. A method for solving the original observation equation of a GNSS observation station network, characterized in that: The following steps are involved: (1) Using a single GNSS dual-frequency or multi-frequency raw observation data, construct the non-differenced and non-combined observation equation: Where, the superscript s represents the satellite number, the subscript r represents the receiver number, the subscript f represents the frequency band number, e represents the epoch number, P represents the pseudorange observation value in meters, and L represents the carrier phase observation value in meters. is the geometric distance from the receiver r to the satellite s, which is a known value in meters. is the model value of the tropospheric projection function from receiver r to satellite s, T r is the tropospheric delay in the zenith direction, in meters, dt r is the clock error of the receiver r, in meters, dt s is the clock error of satellite s, in meters, b r,f is the pseudorange hardware delay of receiver r corresponding to frequency f, in meters, is the pseudorange hardware delay of satellite s corresponding to frequency f, in meters, f r,f is the phase hardware delay of the receiver r corresponding to the frequency f, in meters, is the phase hardware delay of satellite s corresponding to frequency f, in meters, is the ionospheric delay of the first frequency of the receiver r and satellite s, in meters, is the ratio of the ionospheric delay of other frequencies to the first frequency, λ f is the wavelength corresponding to the frequency f, in meters, is the integer ambiguity corresponding to frequency f between receiver r and satellite s, in meters; (2) Add constraint type 1: Where Δdt i is the receiver clock correction value of the i-th station, Δdt i is the correction value of the i-th satellite clock error, Pdt i is the weight of the receiver clock correction value of the i-th station, Pdt i is the weight of the i-th satellite clock correction value, nr is the number of stations, and ns is the number of satellites; (3) Add constraint type 2: Select the first and second frequencies and add constraint b to each station and satellite. r,f =0 and Where f takes two selected frequencies; (4) Add constraint type 3: Add constraints for other frequencies other than the frequency selected in constraint 2. Where f is a frequency other than the frequency selected in constraint 2, Pb i,f is the weight of the receiver pseudorange hardware delay of the i-th station, is the weight of the hardware delay of the ith satellite pseudorange, nr is the number of stations, and ns is the number of satellites; (5) Add constraint type 4: Add conditions for each frequency Where Pf i,f is the weight of the receiver phase hardware delay of the i-th station, is the weight of the hardware delay of the i-th satellite phase, nr is the number of stations, and ns is the number of satellites; (6) Add constraint type 5: Add conditions for each frequency at each station. Where ns is the number of satellites selected as constraints; (7) Add constraint type 6: Add conditions for each frequency of each satellite Where nr is the number of stations selected as constraints; (8) Combine the observation equation and the six types of constraints to solve the parameter T r , Except for the zenith tropospheric delay, which is still T r In addition, the specific forms of other parameters are as follows: Among them, B r =b r,2 -b r,1 , a=1 / (μ2-1), β=m2 / (μ2-1), γ f =a×(μ f +1),k f =a×(μ f +m2),

Images