A doppler velocity driven PPP-B2b robust adaptive precise positioning method and system

Abstract

The application discloses a kind of PPP-B2b anti-difference adaptive precision positioning method and system driven by Doppler velocity measurement, for the problem that prediction model deviation is large and is susceptible to sudden observation value interference in dynamic scene, the application utilizes Doppler frequency shift estimation receiver three-dimensional instantaneous velocity;Utilize three-dimensional instantaneous velocity to calculate pseudo-range prediction variation, and cross consistency test is carried out with pseudo-range actual variation comparison to remove abnormal observation value;In extended Kalman filter, three-dimensional instantaneous velocity is used as speed control input to drive state prediction, and process noise is determined adaptively by velocity measurement accuracy;According to the calculation of standard residual error, the equivalent weight factor is used to adjust the observation noise covariance matrix, and the filtering solution is completed.The application makes the prediction model accurately fit the real motion trajectory of receiver, effectively suppresses multipath interference and sudden noise, significantly shortens the filtering convergence time, and improves the dynamic positioning accuracy and robustness in complex observation environment.

Description

Technical Field

[0001] This invention relates to the field of satellite navigation and positioning technology, specifically to a Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method and system. Background Technology

[0002] Precise point positioning technology utilizes a single satellite navigation receiver in conjunction with precise orbit and clock bias products to achieve high-precision coordinate acquisition. The PPP-B2b service launched by the BeiDou-3 system broadcasts precise ephemeris correction information via space signals, providing a core data source for real-time precise point positioning without reference station support.

[0003] In conventional real-time precision positioning data processing, the extended Kalman filter is a fundamental mathematical tool for achieving optimal parameter estimation. Its computation is highly dependent on the construction of the state prediction model. For user terminals in motion, existing technologies typically employ random walk or uniform velocity models to extrapolate states between epochs, relying on pre-defined empirical process noise parameters to cover unknown physical displacement changes.

[0004] However, existing filtered state prediction models are disconnected from the actual physical motion of the receiver, resulting in deficiencies in their positioning performance in dynamic scenarios. Because traditional prediction models do not incorporate real-time velocity information as a physical constraint, they passively adapt to position changes by relying solely on empirically amplified process noise. This makes the state variables output in the prediction step highly prone to deviating from the actual motion trajectory. This prediction bias leads to longer convergence times for the filter during the solution process and causes over-constraint in parameter estimation, ultimately resulting in a significant degradation in the accuracy of dynamic continuous positioning. Summary of the Invention

[0005] To address the shortcomings of existing technologies, this invention provides a Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method and system, which solves the problems of large prediction deviations and susceptibility to abnormal and sudden changes in observation values ​​caused by the static prediction model and fixed observation weights in dynamic scenes and complex observation environments in existing PPP-B2b real-time precision single-point positioning.

[0006] To achieve the above objectives, the present invention provides the following technical solution: The first aspect of the present invention provides a Doppler velocity measurement driven PPP-B2b robust adaptive precision positioning method.

[0007] The method includes: performing multi-source GNSS observation data acquisition and PPP-B2b precision correction decoding; performing Doppler instantaneous velocimetry and receiver velocity field construction; and using weighted least squares method to estimate the receiver's three-dimensional instantaneous velocity and receiver clock drift.

[0008] The velocity-pseudorange cross-consistency test and outlier detection are performed. The pseudorange prediction change is calculated using the receiver's three-dimensional instantaneous velocity, and the difference between the actual pseudorange change and the pseudorange prediction change is used as the consistency test statistic to eliminate outlier observations.

[0009] The adaptive dynamic model driven by Doppler velocity was constructed. The receiver's three-dimensional instantaneous velocity was used as the velocity control input in the state prediction equation of the extended Kalman filter, and the process noise was adaptively determined by the Doppler velocity measurement accuracy.

[0010] Perform robust adaptive Kalman filtering based on equivalent weights, calculate the equivalent weight factors according to the standardized residuals, and apply the equivalent weight factors to the observation noise covariance matrix to complete the Kalman filtering solution; then perform convergence judgment and accuracy evaluation of the positioning results.

[0011] Furthermore, the process of performing Doppler instantaneous velocimetry and constructing the receiver velocity field includes: converting the Doppler frequency shift observations of the visible satellite into satellite-to-ground radial velocities; simultaneously establishing the radial velocity observation equations of the visible satellite, wherein the radial velocity observation equations include the satellite direction vector, the satellite velocity vector, the receiver's three-dimensional instantaneous velocity, the receiver clock drift, and the satellite clock drift; and solving the radial velocity observation equations using the weighted least squares method to estimate the receiver's three-dimensional instantaneous velocity and the receiver clock drift.

[0012] Furthermore, the execution of velocity-pseudorange cross-consistency test and abnormal observation detection includes: using the receiver's three-dimensional instantaneous velocity to predict the displacement increment between epochs, and combining it with the receiver clock error change based on the receiver clock drift prediction, to calculate the pseudorange prediction change of each satellite channel.

[0013] The difference between the pseudorange observation value of the current epoch and the pseudorange observation value of the previous epoch is calculated to obtain the actual change in pseudorange. When the consistency test statistic exceeds the anomaly detection threshold, the corresponding pseudorange observation value is determined to be an anomaly observation value, and the anomaly observation value is removed from the localization solution of the current epoch. The anomaly detection threshold is calculated as follows:

[0014] Calculate the square of the product of the Doppler velocity noise standard deviation and the sampling interval, add it to the square of twice the pseudorange noise standard deviation, take the square root of the sum, and multiply it by the detection sensitivity factor to obtain the anomaly detection threshold.

[0015] Furthermore, the construction of the adaptive dynamic model driven by Doppler velocity includes: constructing a state vector containing receiver position, receiver clock error, tropospheric wet delay, and carrier phase ambiguity; in the state prediction equation of the extended Kalman filter, setting the position component of the control input matrix as the sampling interval multiplied by the identity matrix, and substituting the receiver's three-dimensional instantaneous velocity into the control input matrix to drive state prediction; in the covariance propagation of state prediction, adaptively determining the process noise of the position component in the state vector by calculating the square of the product of the Doppler velocity noise standard deviation and the sampling interval.

[0016] Furthermore, in the process of performing robust adaptive Kalman filtering based on equivalent weights, the calculation steps of the equivalent weight factors include: establishing a joint observation equation that includes pseudorange observations and carrier phase observations; calculating the ratio of the filtered prediction residual to the standard deviation of the prediction residual to obtain the standardized residual;

[0017] When the absolute value of the standardized residual is less than or equal to the normal observation boundary, the equivalent weight factor is set to one; when the absolute value of the standardized residual is greater than the normal observation boundary and less than or equal to the complete elimination boundary, the equivalent weight factor is smoothly calculated based on the absolute value of the standardized residual, the normal observation boundary, and the complete elimination boundary; when the absolute value of the standardized residual is greater than the complete elimination boundary, the equivalent weight factor is set to zero.

[0018] The Kalman filter solution steps include: dividing the diagonal elements corresponding to the original observation noise covariance matrix by the equivalent weight factor to obtain a robust equivalent weight matrix; substituting the robust equivalent weight matrix into the Kalman filter update equation to calculate the Kalman gain matrix, thereby completing the posterior estimation of the state vector and the state covariance matrix.

[0019] Furthermore, the process of performing multi-source GNSS observation data acquisition and PPP-B2b precise correction decoding includes: receiving pseudorange observations, carrier phase observations, Doppler shift observations, and signal-to-noise ratios acquired by a multi-frequency GNSS receiver; synchronously receiving and decoding PPP-B2b signals to obtain precise orbit corrections, precise clock corrections, and ionospheric corrections; and using the precise orbit corrections and precise clock corrections to correct the broadcast ephemeris and obtain precise satellite positions and precise satellite clock errors.

[0020] Furthermore, the method also includes a multi-constellation joint processing step: adding inter-system deviation parameters of the Global Positioning System and the Galileo system relative to the BeiDou system to the state vector; adding the corresponding inter-system deviation parameter terms to the observation equations of the Global Positioning System and the Galileo system satellites respectively, thereby extending the positioning processing to the multi-constellation joint processing solution of the BeiDou-3, Global Positioning System and Galileo system.

[0021] A second aspect of the present invention provides a Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning system for implementing the method described in any one of the first aspects, the positioning system comprising:

[0022] The data acquisition module is used to perform multi-source GNSS observation data acquisition and PPP-B2b precision correction data decoding.

[0023] The Doppler velocity measurement module is used to perform Doppler instantaneous velocity measurement and receiver velocity field construction, and uses the weighted least squares method to estimate the receiver's three-dimensional instantaneous velocity and receiver clock drift.

[0024] The anomaly detection module is used to perform velocity-pseudorange cross-consistency tests and anomaly detection. It calculates the pseudorange prediction change using the receiver's three-dimensional instantaneous velocity and uses the difference between the actual pseudorange change and the pseudorange prediction change as a consistency test statistic to eliminate anomalies.

[0025] The dynamic model building module is used to complete the construction of the adaptive dynamic model driven by Doppler velocity. In the state prediction equation of the extended Kalman filter, the receiver's three-dimensional instantaneous velocity is used as the velocity control input, and the process noise is adaptively determined by the Doppler velocity measurement accuracy.

[0026] The robust filtering module is used to perform robust adaptive Kalman filtering based on equivalent weights. It calculates the equivalent weight factors based on the standardized residuals and applies the equivalent weight factors to the observation noise covariance matrix to complete the Kalman filtering solution.

[0027] The accuracy assessment module is used to determine the convergence of positioning results and output the accuracy assessment.

[0028] This invention provides a Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method and system. It has the following beneficial effects:

[0029] 1. This invention uses the receiver's three-dimensional instantaneous velocity as the velocity control input for extended Kalman filter state prediction, and adaptively determines the process noise of the position component by Doppler velocity measurement accuracy. This enables the state prediction model to accurately fit the receiver's actual motion trajectory, achieving rapid convergence of the filter in dynamic scenarios. It avoids the prediction deviation and over-constraint problems caused by empirical settings in traditional random walk models, significantly shortening the convergence time and improving the dynamic continuous positioning accuracy.

[0030] 2. This invention utilizes the receiver's three-dimensional instantaneous velocity to calculate the pseudorange prediction change to construct a consistency test statistic to eliminate outlier observations. It also combines an equivalent weighting factor based on standardized residual calculation to adaptively adjust the observation noise covariance matrix, thereby achieving real-time detection of pseudorange mutations and automatic weight reduction of observations affected by multipath interference. This overcomes the shortcomings of fixed-weight models that cannot adapt to time-varying noise and improves the robustness of positioning results in complex observation environments.

[0031] 3. This invention extends single-system positioning processing to the joint processing and calculation of multiple constellations of BeiDou-3, GPS, and Galileo by adding inter-system deviation parameters of GPS and Galileo relative to BeiDou system to the state vector. This effectively increases the number of available satellites and improves the satellite spatial geometry configuration, solves the problem of insufficient visible satellites in a single system under obstructed environment, and ensures the continuity of positioning while reducing the position accuracy factor. Attached Figure Description

[0032] Figure 1 This is an overall flowchart of a Doppler velocimetry-driven PPP-B2b robust adaptive precision positioning method according to an embodiment of the present invention;

[0033] Figure 2 This is a graph showing the verification results of the Doppler velocity measurement accuracy and inter-epoch displacement prediction accuracy of each station in the embodiments of the present invention;

[0034] Figure 3 This is a graph showing the anomaly detection performance based on the velocity-pseudorange cross-consistency test in an embodiment of the present invention.

[0035] Figure 4 This is a bar chart comparing the positioning accuracy and convergence time of the velocity-driven dynamic model and the traditional random walk model in this embodiment of the invention.

[0036] Figure 5 This is a comparison chart showing the performance improvement of robust adaptive estimation under different levels of multipath interference in the embodiments of the present invention;

[0037] Figure 6 This is a comparison chart of the number of visible satellites and positioning accuracy of BeiDou-3 single system and multi-constellation joint processing in an embodiment of the present invention. Detailed Implementation

[0038] The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0039] Please see the appendix Figure 1 , Figure 1 This is a flowchart illustrating a Doppler velocity measurement-driven PPP-B2b robust adaptive precise positioning method according to an embodiment of the present invention. The overall workflow of this Doppler velocity measurement-driven PPP-B2b robust adaptive precise positioning method consists of the following processing stages:

[0040] Step S1: Perform multi-source GNSS observation data acquisition and PPP-B2b precise correction decoding. Receive the acquired pseudorange observations, carrier phase observations, Doppler shift observations, and signal-to-noise ratio. Simultaneously receive and decode the PPP-B2b signal to obtain precise orbit corrections, precise clock corrections, and ionospheric corrections. Use the obtained precise orbit corrections and precise clock corrections to correct the broadcast ephemeris and obtain precise satellite positions and precise satellite clock errors.

[0041] Step S2 involves performing Doppler instantaneous velocimetry and constructing the receiver velocity field. The Doppler frequency shift observations of the visible satellite are converted into satellite-to-ground radial velocities. The radial velocity observation equations of the visible satellite are then combined, and the weighted least squares method is used to estimate the receiver's three-dimensional instantaneous velocity and receiver clock drift.

[0042] Step S3: Perform velocity-pseudorange cross-consistency check and outlier detection. Calculate the predicted pseudorange change using the receiver's three-dimensional instantaneous velocity to predict the displacement increment between epochs. The difference between the actual pseudorange change and the predicted pseudorange change is used as the consistency check statistic. When the consistency check statistic exceeds the outlier detection threshold, the corresponding pseudorange observation is determined to be an outlier and removed from the positioning solution for the current epoch.

[0043] Step S4 completes the construction of the Doppler velocity-driven adaptive dynamic model. A state vector is constructed, including receiver position, receiver clock error, tropospheric wet delay, and carrier phase ambiguity. In the extended Kalman filter's state prediction equation, the receiver's three-dimensional instantaneous velocity is used as the velocity control input, and the process noise of the position component is adaptively determined by the Doppler velocimetry accuracy.

[0044] Step S5: Perform robust adaptive Kalman filtering based on equivalent weights. Establish a joint observation equation including pseudorange and carrier phase observations. In the filter update, calculate the equivalent weight factor based on the standardized residual determined by the standard deviation of the prediction residual. Apply the equivalent weight factor to the observation noise covariance matrix to achieve adaptive adjustment of the observation weights and complete the Kalman filtering solution.

[0045] Step S6: Perform convergence determination and accuracy assessment of the positioning results. The filtered output position sequence is converged; convergence is achieved when the position change for a consecutive preset number of epochs is less than the convergence threshold. The positioning accuracy in the east, north, and zenith directions is calculated, and the positioning results, including three-dimensional position coordinates, receiver velocity vector, convergence status flag, estimated positioning accuracy, number of satellites used, and position accuracy factor (PDOP) value, are output.

[0046] Step S7: Perform multi-constellation joint processing and inter-system bias elimination. Add inter-system bias parameters between the GPS and Galileo systems and the BeiDou system to the state vector, extending the positioning processing to the multi-constellation joint processing solution of BeiDou-3, GPS, and Galileo.

[0047] The specific steps of the Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method are as follows:

[0048] Step S1: Perform multi-source GNSS observation data acquisition and PPP-B2b precision correction data decoding.

[0049] The S101 employs a multi-frequency GNSS receiver that supports BeiDou-3 B2b signals, with a sampling interval Collect satellite observation data. Sampling interval. The value ranges from 1 to 30 seconds. The collected observation data includes pseudorange observations. Carrier phase observations Doppler frequency shift observations and signal-to-noise ratio Pseudorange observations This represents the pseudorange measurement from the receiver to the satellite, in meters; carrier phase observations. This represents the carrier phase measurement, in cycles; Doppler frequency shift observation. The signal-to-noise ratio (SNR) represents the carrier Doppler frequency shift, measured in Hertz. This indicates the carrier-to-noise ratio, measured in decibels (dB / Hertz).

[0050] S102 synchronously receives and decodes PPP-B2b signals to obtain precise track correction data. Precision clock error correction and ionospheric correction Precision track correction numbers These are corrections for the radial, tangential, and normal components of the satellite orbit, in meters; precision clock error corrections. Satellite clock correction, in nanoseconds; ionospheric correction. This is the grid correction value for the total vertical electron content of the ionosphere, expressed in TECU.

[0051] S103, utilizing precision track correction numbers and precision clock error correction Correct the broadcast ephemeris to obtain precise satellite positions. and precision satellite clock bias :

[0052]

[0053]

[0054] In the formula, Satellite positions calculated for broadcast ephemeris. For broadcast ephemeris clock errors.

[0055] Step S2: Perform Doppler instantaneous velocities measurement and construct the receiver velocity field.

[0056] S201, the first The satellite in the Doppler shift of epochs Converted to satellite-to-ground radial velocity :

[0057]

[0058] In the formula, The speed of light is 2.99792458 × 10^8 m / s; The carrier frequency is expressed in Hertz (Hz).

[0059] S202, the radial velocity observation equation is:

[0060]

[0061] In the formula, For the receiver to point to the first The unit direction vector of a satellite; This is the satellite velocity vector, obtained from the precise ephemeris time derivative, in m / s; This is the receiver velocity vector, in m / s; This refers to receiver clock drift, measured in seconds (s / s). This refers to satellite clock drift, measured in seconds (s / s). This is for Doppler observation noise.

[0062] S203, solve the simultaneous radial velocity equations for all visible satellites, and use the weighted least squares method to estimate the receiver's three-dimensional velocity. And Zhong Piao :

[0063]

[0064] In the formula, This is the design matrix for the Doppler observation equations, where each row of the design matrix consists of satellite direction vectors. Composed of clock drift coefficient; This is the weight matrix of the Doppler observations; To observe the residual vector.

[0065] Step S3: Perform velocity-pseudorange cross-consistency test and outlier detection.

[0066] S301, utilizing receiver speed Predict the displacement increment of the current epoch relative to the previous epoch. :

[0067]

[0068] S302, Calculate the first displacement increment based on the predicted displacement increment. The predicted change in pseudorange of each satellite at the current epoch. :

[0069]

[0070] In the formula, This represents the change in receiver clock bias based on clock drift prediction.

[0071] S303, Calculate the actual change in pseudorange The difference between the predicted change and the pseudorange is used to define the consistency test statistic. :

[0072]

[0073] S304, when the consistency test statistic Exceeding the anomaly detection threshold At that time, the judgment of the first The pseudorange observations of these satellites at the current epoch are outliers. Anomaly detection threshold. Determined jointly by Doppler velocity measurement accuracy and pseudorange noise level:

[0074]

[0075] In the formula, This represents the standard deviation of Doppler velocity measurement noise, ranging from 0.001 to 0.01 m / s, with a default value of 0.003 m / s. This represents the standard deviation of pseudorange noise, ranging from 0.3 to 3.0 meters, with a default value of 1.0 meter. The detection sensitivity factor has a value range of 3.0 to 5.0, with a default value of 3.5. Detected abnormal satellites are removed from the positioning calculation in the current epoch and are not included in the filter update.

[0076] Step S4: Complete the construction of the Doppler velocity-driven adaptive dynamic model.

[0077] S401, Introducing Doppler velocity into the extended Kalman filter state prediction model, defining the state vector. :

[0078]

[0079] In the formula, These are the receiver's position coordinates, in meters. This refers to the receiver clock error, expressed in seconds. This is an estimate of the tropospheric wet delay, in meters; This is the carrier phase ambiguity parameter vector, in cycles.

[0080] S402, the state prediction equation for speed-driven operation is:

[0081]

[0082] In the formula, The state transition matrix is ​​defined as follows: the position component of the state transition matrix is ​​the identity matrix, the clock error component is 1, and the tropospheric and ambiguity components are 1. To control the input matrix, the position components of the input matrix correspond to... Multiply by the identity matrix, and the remaining components are 0; The speed control input is obtained from Doppler velocimetry.

[0083] S403, the covariance propagation of state prediction is:

[0084]

[0085] In the formula, the process noise matrix process noise of position components Determined by the accuracy of Doppler velocity measurement:

[0086]

[0087] Clock difference component process noise ; Tropospheric wet delay component process noise to The ambiguity parameter and process noise are set to 0.

[0088] Step S5: Perform robust adaptive Kalman filtering based on equivalent weights.

[0089] S501, construct the joint observation equation of pseudorange and carrier phase. For the first... For each satellite, the pseudorange observation equation is:

[0090]

[0091] The carrier phase observation equation is:

[0092]

[0093] In the formula, The geometric distance between the satellite and the receiver, in meters; For the first The clock bias of each satellite; Tropospheric delay, in meters; Ionospheric delay, measured in meters; The carrier wavelength is measured in meters. For integer ambiguity; This is pseudorange observation noise. This is carrier phase observation noise.

[0094] S502, in the Kalman filter update step, introduces an equivalent weight function. For the ... Observation values, based on standardized residuals Determine the equivalent weight factor :

[0095]

[0096] In the formula, For the filtered prediction residual, This represents the standard deviation of the corresponding predicted residuals. The equivalent weighting factor adopts the IGGIII scheme:

[0097] when hour:

[0098]

[0099] when hour:

[0100]

[0101] when hour:

[0102]

[0103] In the formula, This represents the boundary of normal observations, with a value range of 1.0 to 2.0 and a default value of 1.5. To completely eliminate boundaries, the value ranges from 2.5 to 4.5, with a default value of 3.0.

[0104] S503 applies equivalent weighting factors to the observation noise covariance matrix. Obtain the robust equivalent weight matrix :

[0105]

[0106] The filter update equation is:

[0107]

[0108]

[0109]

[0110] In the formula, The Kalman gain matrix; The observation matrix; For observation vectors; It is an identity matrix.

[0111] Step S6: Perform a convergence determination and accuracy assessment of the positioning results.

[0112] S601 performs a convergence check on the position sequence of the filtered output. The convergence condition is defined as continuity. The positional changes of each epoch are all less than the convergence threshold. :

[0113]

[0114] In the formula, To determine the convergence window, the value range is 10 to 60 epochs, with a default value of 20 epochs; The convergence threshold ranges from 0.1 to 0.5 meters, with a default value of 0.3 meters.

[0115] S602, calculate the positioning accuracy in the three directions of east, north, and zenith:

[0116]

[0117] In the formula, , , These are the diagonal elements of the covariance matrix corresponding to the components in the local coordinate system.

[0118] S603 outputs positioning results, including: three-dimensional position coordinates, receiver velocity vector, convergence status flag, estimated positioning accuracy in the three directions of East / North / Zenith, number of satellites used, and Position Accuracy Factor (PDOP) value.

[0119] Step S7: Perform multi-constellation joint processing and eliminate inter-system biases.

[0120] S701 extends positioning processing to joint processing of multiple constellations including BeiDou-3, Global Positioning System, and Galileo. An inter-system bias parameter is added to the state vector, resulting in the extended state vector:

[0121]

[0122] In the formula, This refers to the receiver clock bias with reference to the BeiDou system. The inter-system deviation of the Global Positioning System relative to the BeiDou system is expressed in meters. The system deviation between the Galileo system and the BeiDou system is expressed in meters.

[0123] S702, for GPS and Galileo system satellites, respectively, adds corresponding... or The inter-system deviation parameter is used as the parameter to be estimated in the Kalman filter, with an initial value set to 0, and the process noise set to... The Doppler velocities in step S2, the anomaly detection in step S3, and the robust estimation in step S5 are simultaneously applicable to both GPS and Galileo satellite systems.

[0124] Based on the same inventive concept as the aforementioned Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method, this invention provides a Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning system. The positioning system includes a data acquisition module, a Doppler velocity measurement module, an anomaly detection module, a dynamic model construction module, a robust filtering module, an accuracy evaluation module, and a multi-constellation processing module.

[0125] The data acquisition module is used to receive multi-frequency GNSS observation data and PPP-B2b precise corrections. Specifically, the data acquisition module receives pseudorange observations, carrier phase observations, Doppler shift observations, and signal-to-noise ratios acquired by the multi-frequency GNSS receiver. It simultaneously receives and decodes PPP-B2b signals to obtain precise orbit corrections, precise clock corrections, and ionospheric corrections. It then uses the precise orbit corrections and precise clock corrections to correct the broadcast ephemeris and obtain precise satellite positions and precise satellite clock errors.

[0126] The Doppler velocimetry module is used to estimate the receiver's three-dimensional instantaneous velocity and clock drift based on Doppler frequency shift. Specifically, the Doppler velocimetry module converts the Doppler frequency shift observations of visible satellites into satellite-to-ground radial velocities, combines the radial velocity observation equations of the visible satellites, and uses the weighted least squares method to estimate the receiver's three-dimensional instantaneous velocity and receiver clock drift.

[0127] The anomaly detection module is used to perform a cross-consistency check between the predicted pseudorange change and the actual pseudorange change to identify and remove abnormal observations. Specifically, the anomaly detection module uses the displacement increment between epochs of the receiver's three-dimensional instantaneous velocity prediction to calculate the predicted pseudorange change, and uses the difference between the actual pseudorange change and the predicted pseudorange change as a consistency check statistic. When the consistency check statistic exceeds the anomaly detection threshold, the corresponding pseudorange observation is determined to be an abnormal observation, and the abnormal observation is removed from the positioning solution of the current epoch.

[0128] The dynamic model building module is used to drive the extended Kalman filter state prediction using Doppler velocity as a control input, and to determine process noise based on the velocity measurement accuracy. Specifically, the dynamic model building module constructs a state vector including receiver position, receiver clock error, tropospheric wet delay, and carrier phase ambiguity, and uses the receiver's three-dimensional instantaneous velocity as the velocity control input in the extended Kalman filter state prediction equation. Simultaneously, the process noise of the position component is adaptively determined by the Doppler velocity measurement accuracy.

[0129] The robust filtering module is used to adaptively weight the observations using an equivalent weight function to achieve robust Kalman filtering. Specifically, the robust filtering module establishes a joint observation equation including pseudorange and carrier phase observations, calculates equivalent weight factors based on the standardized residuals, applies these equivalent weight factors to the observation noise covariance matrix to obtain a robust equivalent weight matrix, and uses this matrix to complete the Kalman filter update solution.

[0130] The accuracy assessment module is used for convergence determination and accuracy index calculation. Specifically, the accuracy assessment module performs convergence determination on the position sequence of the filtered output, calculates the positioning accuracy in the three directions of east, north, and zenith, and outputs the positioning result including three-dimensional position coordinates, receiver velocity vector, convergence status flag, estimated positioning accuracy, number of satellites used, and position accuracy factor (PDOP) value.

[0131] The multi-constellation processing module is used for joint calculation of multiple systems and elimination of inter-system biases. Specifically, the multi-constellation processing module adds inter-system bias parameters between the Global Positioning System (GPS) and Galileo system and the BeiDou system to the state vector, extending positioning processing to the joint calculation of multiple constellations of BeiDou-3, GPS, and Galileo systems.

[0132] To provide hardware support for the above method, this embodiment of the invention also provides a GNSS positioning device. The GNSS positioning device includes: a multi-frequency GNSS receiver, a PPP-B2b signal decoding unit, an embedded processor, a storage unit, and a communication module.

[0133] Multi-frequency GNSS receivers are used to acquire satellite observation data and Doppler shift.

[0134] The PPP-B2b signal decoding unit is used to receive and decode precision correction data.

[0135] An embedded processor is used to run the aforementioned Doppler velocimetry-driven PPP-B2b robust adaptive precision positioning method.

[0136] Storage unit, used to store observation data and positioning results.

[0137] The communication module is used to output positioning results and speed information.

[0138] This invention also provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the aforementioned Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method.

[0139] See attached document Figure 2 , Figure 2 This figure shows the verification results of Doppler velocity measurement accuracy and displacement prediction accuracy according to an embodiment of the present invention. Five GNSS monitoring stations distributed at different latitudes were selected, and BeiDou-3 observation data were collected continuously for 72 hours, with a sampling interval of 1 second. The accuracy of Doppler velocity measurement and the accuracy of velocity-based inter-epoch displacement prediction were evaluated using the precise post-processing results as the reference true value. The velocity RMSE at station A was 2.1 mm / s, and the displacement prediction RMSE was 1.9 mm; the velocity RMSE at station B was 2.3 mm / s, and the displacement prediction RMSE was 2.1 mm; the velocity RMSE at station C was 1.8 mm / s, and the displacement prediction RMSE was 1.6 mm; the velocity RMSE at station D was 3.2 mm / s, and the displacement prediction RMSE was 2.9 mm; and the velocity RMSE at station E was 2.5 mm / s, and the displacement prediction RMSE was 2.3 mm. The average RMSE of the three-dimensional velocity measured by Doppler at each station was 2.4 mm / s. The average RMSE of displacement prediction within one epoch based on velocity measurement results is 2.2 mm. Specific velocity measurement and displacement prediction accuracy data are shown in Table 1. The experimental data above confirm that adaptively determining process noise based on Doppler velocity measurement accuracy can ensure a strict match between process noise and the actual motion state of the receiver. When the receiver is stationary, the filter fully utilizes historical information; when in motion, it effectively avoids the over-constraint problem caused by traditional empirical settings.

[0140] Table 1: Accuracy of Doppler velocity and displacement prediction at each station

[0141]

[0142] See attached document Figure 3 , Figure 3This is a performance chart of velocity-pseudorange cross-consistency check anomaly detection according to an embodiment of the present invention. Pseudorange anomalies of different amplitudes were artificially injected into the same station data. For pseudorange anomalies with an amplitude greater than or equal to 3 meters, the detection rate was 98.7% with a false positive rate of 0.3%; for pseudorange anomalies with an amplitude greater than or equal to 1 meter, the detection rate was 92.1% with a false positive rate of 0.8%; and for anomalies with an amplitude greater than or equal to 0.5 meters, the detection rate was 78.5% with a false positive rate of 1.2%. Specific detection performance statistics for different anomaly amplitudes are shown in Table 2. Experimental results show that by using Doppler velocity measurements to predict pseudorange changes and comparing them with actual pseudorange changes, real-time detection of pseudorange anomaly observations on a satellite-by-satellite and epoch-by-epoch basis can be achieved. The cross-consistency check mechanism effectively eliminates the dependence on continuous carrier phase locking and can maintain extremely high detection accuracy even in complex obstruction environments where satellite signals frequently lose lock.

[0143] Table 2: Detection performance under different anomaly amplitudes

[0144]

[0145] See attached document Figure 4 , Figure 4 This is a comparison chart of the positioning accuracy of a velocity-driven dynamic model and a traditional model according to an embodiment of the present invention. The positioning performance of the random walk model scheme in traditional PPP-B2b positioning and the velocity-driven combined robust estimation scheme of the present invention are compared. Data from five consecutive days were processed using a single BeiDou-3 system configuration. The average positioning accuracy of traditional PPP-B2b positioning in the east, north, and zenith directions was 14.3 cm, 9.3 cm, and 21.3 cm, respectively, with an average convergence time of 16.3 minutes. The positioning method of the present invention achieved average positioning accuracies of 10.8 cm, 6.5 cm, and 16.2 cm in the east, north, and zenith directions, respectively, representing improvements of 24.5%, 30.1%, and 23.9% compared to the traditional scheme, with an average convergence time of 5.3 minutes, a reduction of 67.5%. Detailed comparison data of the positioning accuracy and convergence time between the velocity-driven dynamic model and the traditional random walk model are shown in Table 3. Comparative data confirms that introducing the Doppler instantaneous velocity into the prediction step of the extended Kalman filter, replacing the traditional random walk model, significantly accelerates the filter's convergence process. In dynamic scenarios, the velocity-driven prediction model can accurately track the receiver's actual motion trajectory, eliminating the accuracy degradation caused by excessive prediction bias in traditional models.

[0146] Table 3: Comparison of Positioning Accuracy and Convergence Time

[0147]

[0148] See attached document Figure 5 , Figure 5This is a graph illustrating the performance of robust adaptive estimation under multipath conditions according to an embodiment of the present invention. The localization performance of fixed-weight extended Kalman filtering and the robust adaptive extended Kalman filtering of the present invention is compared under normal conditions and different levels of multipath interference. Under normal conditions, the performance of the two methods is similar. Under mild multipath conditions, the RMSE of fixed-weight extended Kalman filtering degrades to 32.6 cm in the zenith direction, while the method of the present invention achieves 18.3 cm, an improvement of 43.9%. Under moderate multipath conditions, the RMSE of fixed-weight extended Kalman filtering degrades to 74.5 cm in the zenith direction, while the method of the present invention achieves 25.6 cm, an improvement of 65.6%. Under severe multipath conditions, the RMSE of fixed-weight extended Kalman filtering degrades to 125.8 cm in the zenith direction, while the method of the present invention achieves 29.7 cm, an improvement of 76.4%. Detailed comparison data of RMSE in the zenith direction (U direction) under different multipath conditions are shown in Table 4. The above comparison demonstrates the performance advantage of robust adaptive Kalman filtering based on IGGIII equivalent weights under degraded environments. The observation weights are automatically adjusted in real time based on the standardized residuals. Without the need for a priori noise model, the system can automatically reduce or eliminate abrupt noise caused by multipath interference, ionospheric scintillation, etc., which significantly enhances the robustness of the positioning results.

[0149] Table 4: Comparison of RMSE in the U-direction under different multipath environments (unit: cm)

[0150]

[0151] See attached document Figure 6 , Figure 6 This is a performance comparison chart of multi-constellation joint positioning according to an embodiment of the present invention. In the multi-constellation joint positioning performance verification experiment, the performance of the BeiDou-3 single system was compared with that of the BeiDou-3, Global Positioning System (GPS), and Galileo system joint positioning. The average number of visible satellites in the three-system joint positioning increased from 11.8 in the single system to 27.3, and the PDOP value decreased from 2.21 to 1.38. The positioning accuracy of the three-system joint positioning in the east, north, and zenith directions was 8.3 cm, 4.8 cm, and 12.5 cm, respectively, representing improvements of 23.1%, 26.2%, and 22.8% compared to the single system. The convergence time was further reduced from 5.3 minutes to 3.6 minutes. The performance indicators of the single system and multi-constellation joint positioning are compared in Table 5. The test results show that by expanding and estimating the inter-system bias parameters in the state vector, the multi-constellation joint processing increases the number of available satellites, optimizes the satellite spatial geometry, and significantly reduces the PDOP value, ensuring that the equipment still has a sufficient number of visible satellites and continuous positioning accuracy even in obstructed environments.

[0152] Table 5: Performance Comparison of Single-System and Multi-Constellation Joint Positioning

[0153]

[0154] Furthermore, in practical engineering deployment, the method provided by this invention has moderate computational complexity, the Doppler velocity measurement module outputs independent least squares solutions, and the robust estimation stage only adds weight update operation logic to the original covariance matrix. The memory and computing power consumption are both within a controllable range, meeting the requirements for real-time high-precision calculation on embedded hardware platforms.

[0155] For the various test data acquisition interface formats, root mean square error and detection rate general statistical calculation formulas, and conventional GNSS monitoring reference station setup specifications involved in the embodiments of this invention, those skilled in the art can refer to existing mature surveying and mapping engineering specifications and GNSS data statistical evaluation criteria for operation. The relevant operation basis and statistical methods are well-known technologies in the field and will not be elaborated here.

[0156] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method, characterized in that, Includes the following steps: Perform multi-source GNSS observation data acquisition and PPP-B2b precision correction data decoding; Doppler instantaneous velocimetry and receiver velocity field construction were performed, and the three-dimensional instantaneous velocity of the receiver and receiver clock drift were estimated using the weighted least squares method. The velocity-pseudorange cross-consistency test and outlier detection are performed. The pseudorange prediction change is calculated using the receiver's three-dimensional instantaneous velocity, and the difference between the actual pseudorange change and the predicted pseudorange change is used as the consistency test statistic to eliminate outlier observations. The adaptive dynamic model driven by Doppler velocity was constructed. The three-dimensional instantaneous velocity of the receiver was used as the velocity control input in the state prediction equation of the extended Kalman filter, and the process noise was adaptively determined by the Doppler velocity measurement accuracy. Perform robust adaptive Kalman filtering based on equivalent weights, calculate the equivalent weight factors according to the standardized residuals, and apply the equivalent weight factors to the observation noise covariance matrix to complete the Kalman filtering solution; Perform convergence determination and accuracy evaluation of positioning results and output.

2. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 1, characterized in that, The process of performing Doppler instantaneous velocities and constructing the receiver velocity field includes: Convert Doppler shift observations of visible satellites into satellite-to-ground radial velocities; Combine the radial velocity observation equations of the visible satellite, wherein the radial velocity observation equations include the satellite direction vector, the satellite velocity vector, the receiver's three-dimensional instantaneous velocity, the receiver clock drift, and the satellite clock drift; The radial velocity observation equation is solved using the weighted least squares method to estimate the receiver's three-dimensional instantaneous velocity and the receiver's clock drift.

3. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 1, characterized in that, The execution speed-pseudorange cross-validation and outlier detection include: The interepoch displacement increment is predicted using the receiver's three-dimensional instantaneous velocity, and combined with the receiver clock error change based on the receiver clock drift prediction, the pseudorange prediction change for each satellite channel is calculated. The difference between the pseudorange observation value of the current epoch and the pseudorange observation value of the previous epoch is calculated to obtain the actual change in pseudorange. When the consistency test statistic exceeds the anomaly detection threshold, the corresponding pseudorange observation is determined to be an anomaly observation, and the anomaly observation is removed from the localization solution in the current epoch.

4. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 3, characterized in that, The method for calculating the anomaly detection threshold is as follows: Calculate the square of the product of the Doppler velocity noise standard deviation and the sampling interval, add it to the square of twice the pseudorange noise standard deviation, take the square root of the sum, and multiply it by the detection sensitivity factor to obtain the anomaly detection threshold.

5. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 1, characterized in that, The construction of the adaptive dynamic model driven by Doppler velocity includes: Construct a state vector that includes receiver position, receiver clock error, tropospheric wet delay, and carrier phase ambiguity; In the state prediction equation of the extended Kalman filter, the position component of the control input matrix is ​​set as the sampling interval multiplied by the identity matrix, and the three-dimensional instantaneous velocity of the receiver is substituted into the control input matrix to drive the state prediction. In the covariance propagation of state prediction, the process noise of the position component in the state vector is adaptively determined by calculating the square of the product of the standard deviation of Doppler velocity noise and the sampling interval.

6. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 1, characterized in that, In the execution of the robust adaptive Kalman filter solution based on equivalent weights, the calculation steps of the equivalent weight factors include: Establish a joint observation equation that includes pseudorange and carrier phase observations; The standardized residual is obtained by calculating the ratio of the filtered prediction residual to the standard deviation of the prediction residual. When the absolute value of the standardized residual is less than or equal to the boundary of the normal observation, the equivalent weight factor takes the value of one. When the absolute value of the standardized residual is greater than the normal observation boundary and less than or equal to the complete elimination boundary, the equivalent weighting factor is smoothly calculated based on the absolute value of the standardized residual, the normal observation boundary, and the complete elimination boundary. When the absolute value of the standardized residual is greater than the complete elimination boundary, the equivalent weight factor is zero.

7. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 6, characterized in that, In the execution of the robust adaptive Kalman filter solution based on equivalent weights, the Kalman filter solution steps include: Divide the diagonal elements of the original observation noise covariance matrix by the equivalent weight factor to obtain the robust equivalent weight matrix. The robust equivalent weight matrix is ​​substituted into the Kalman filter update equation to calculate the Kalman gain matrix, thus completing the posterior estimation of the state vector and the state covariance matrix.

8. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 1, characterized in that, The process of acquiring multi-source GNSS observation data and decoding PPP-B2b precise correction data includes: Receive pseudorange observations, carrier phase observations, Doppler frequency shift observations, and signal-to-noise ratios collected by a multi-frequency GNSS receiver; Simultaneously receive and decode PPP-B2b signals to obtain precise orbital corrections, precise clock error corrections, and ionospheric corrections; The broadcast ephemeris is corrected using the precise orbital correction and the precise clock error correction to obtain precise satellite position and precise satellite clock error.

9. The Doppler velocity measurement-driven PPP-B2b robust adaptive precision positioning method according to claim 1, characterized in that, The method also includes a multi-constellation joint processing step: Add inter-system deviation parameters between the Global Positioning System and Galileo system and the BeiDou system to the state vector; The observation equations for GPS and Galileo satellites are supplemented with corresponding inter-system bias parameters, extending the positioning process to the multi-constellation joint processing and solution of BeiDou-3, GPS, and Galileo systems.

10. A Doppler velocimeter-driven PPP-B2b robust adaptive precision positioning system, used to implement the Doppler velocimeter-driven PPP-B2b robust adaptive precision positioning method according to any one of claims 1-9, characterized in that, include The data acquisition module is used to perform multi-source GNSS observation data acquisition and PPP-B2b precision correction data decoding; The Doppler velocimetry module is used for Doppler instantaneous velocimetry and receiver velocity field construction, and uses the weighted least squares method to estimate the receiver's three-dimensional instantaneous velocity and receiver clock drift. Anomaly detection module is used to perform velocity-pseudorange cross-consistency test and anomaly detection. It calculates the pseudorange prediction change using the receiver's three-dimensional instantaneous velocity and uses the difference between the actual pseudorange change and the predicted pseudorange change as a consistency test statistic to eliminate anomalies. The dynamic model building module is used to complete the adaptive dynamic model building of Doppler velocity drive. In the state prediction equation of the extended Kalman filter, the three-dimensional instantaneous velocity of the receiver is used as the velocity control input, and the process noise is adaptively determined by the Doppler velocity measurement accuracy. The robust filtering module is used to perform robust adaptive Kalman filtering based on equivalent weights. It calculates the equivalent weight factors based on the standardized residuals and applies the equivalent weight factors to the observation noise covariance matrix to complete the Kalman filtering solution. The accuracy assessment module is used to determine the convergence of positioning results and output the accuracy assessment.

Images