A multi-stage joint enhanced positioning method and system assisted by a low earth orbit satellite constellation

Abstract

The application discloses a kind of low-orbit satellite constellation assisted multi-stage joint enhanced positioning method and system, specifically: receiver is carried out consistency test to pseudo-range, carrier phase and Doppler observation value, and rough error of pseudo-range and carrier phase phase slip are detected by Doppler observation value;Joint GNSS observation data and low-orbit satellite pseudo-range information constructs RAIM algorithm, enhances the failure detection rate of GNSS;GNSS phase epoch difference, kalman filter process noise adaptive adjustment and satellite pseudo-range / pseudo-range rate confidence method are used for ambiguity parameter estimation, LAMBDA method is used for ambiguity fixing, and the ambiguity is verified by Ratio test, to obtain GNSS high-precision positioning result;Pressure altimeter assisted GNSS spoofing interference detection and AHRS assisted GNSS spoofing interference detection are carried out.The application can effectively identify and eliminate satellite disturbance signal, improve satellite navigation positioning speed and accuracy, enhance positioning performance, improve the safety and stability of navigation system.

Description

Technical Field

[0001] The present invention relates to the technical field of multi-source combined navigation, and in particular to a multi-stage joint enhanced positioning method and system assisted by a low-orbit satellite constellation. Background Art

[0002] In recent years, the navigation, positioning, and timing functions of global satellite navigation systems have played a vital role in many fields, including emergency rescue. However, due to the inherent fragility of GNSS signals, interference, particularly through deceptive jamming techniques, can significantly impact their normal operation and cause significant losses. Anti-interference positioning systems based on a single GNSS system often perform poorly in environments with strong interference, high dynamics, and multipath conditions, and may even fail to function properly.

[0003] Although GNSS can provide highly accurate user terminal location information, its output frequency is low. In complex obstructed environments such as canyons, tunnels, and forests, GNSS signals are prone to loss of lock, resulting in capture difficulties or signal loss, and the position precision factor is poor. At this time, if a deceptive jammer transmits a deceptive signal, it is easy to deceive user location information and other information. Summary of the Invention

[0004] The purpose of the present invention is to provide a low-orbit satellite joint positioning method and system with low algorithm complexity, high recognition accuracy, fast recognition speed, and the ability to perform multi-level joint detection and anti-interference on satellite disturbance signals.

[0005] The technical solution for achieving the purpose of the present invention is: a multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation, comprising the following steps:

[0006] Step 1: Receiver cycle slip detection and pseudorange gross error detection: Perform consistency check on pseudorange, carrier phase, and Doppler observations, and detect pseudorange gross errors and carrier phase cycle slips through Doppler observations.

[0007] Step 2: Enhanced GNSS fault detection using low-orbit satellites: Build a RAIM algorithm combining GNSS observation data with pseudo-range information from low-orbit satellites to increase the GNSS fault detection rate and enhance positioning performance.

[0008] Step 3: Estimate ambiguity parameters using the GNSS relative positioning model: GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method are used to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the ratio test is used to verify the ambiguity to obtain the GNSS high-precision positioning result.

[0009] Step 4: Barometric altimeter-assisted GNSS spoofing detection: After obtaining the GNSS positioning results, the barometric altimeter is used to assist in GNSS spoofing detection. An adaptive threshold detection method based on moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within a time window for spoofing detection.

[0010] Step 5: AHRS-assisted GNSS spoofing interference detection: An AHRS system is constructed by combining magnetometer data with inertial measurement unit data. The receiver's GNSS plane acceleration is compared with the plane acceleration consisting of the east and north directions measured by the accelerometer in the AHRS system. The difference between the GNSS plane acceleration and the AHRS plane acceleration is calculated to see if there is any jump data within the time window, thereby detecting whether there is spoofing interference.

[0011] A low-orbit satellite constellation-assisted multi-stage joint enhanced positioning system, which is used to implement the low-orbit satellite constellation-assisted multi-stage joint enhanced positioning method described above, and the system includes:

[0012] Detection module, used for receiver cycle slip detection and pseudorange gross error detection: consistency check of pseudorange, carrier phase and Doppler observation values, and detection of pseudorange gross errors and carrier phase cycle slips through Doppler observation values;

[0013] Fault detection module for enhanced GNSS fault detection on LEO satellites: RAIM algorithm is constructed by combining GNSS observation data with LEO satellite pseudorange information to enhance GNSS fault detection rate and positioning performance;

[0014] The ambiguity parameter estimation module estimates ambiguity parameters through the GNSS relative positioning model. It uses GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the ratio test is used to verify the ambiguity to obtain GNSS high-precision positioning results.

[0015] The first auxiliary interference detection module uses a barometric altimeter to assist in GNSS spoofing interference detection: After obtaining the GNSS positioning results, the barometric altimeter is used to assist in GNSS spoofing detection. The adaptive threshold detection method of moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within a time window for spoofing interference detection;

[0016] The second auxiliary interference detection module uses AHRS to assist in GNSS deception interference detection: by combining magnetometer data with inertial measurement unit data to build an AHRS system, the receiver's GNSS plane acceleration is compared with the plane acceleration composed of east and north directions measured by the accelerometer in the AHRS system, and the difference between the GNSS plane acceleration and the AHRS plane acceleration is calculated to see whether there is jump data within the time window, thereby detecting whether there is deception interference.

[0017] A mobile terminal includes a memory, a processor, and a computer program stored in the memory and runnable on the processor, characterized in that when the processor executes the program, the multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation is implemented.

[0018] A computer-readable storage medium stores a computer program, which, when executed by a processor, implements the steps of the multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation.

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

[0020] (1) The GNSS positioning performance is enhanced through the low-orbit satellite constellation, and multiple sensors such as inertial measurement units, magnetometers, and altimeters are combined to perform multi-level joint detection and anti-interference of satellite disturbance signals, thereby improving the robustness, reliability, and stability of the multi-source positioning and navigation system;

[0021] (2) The algorithm is simple, with high recognition accuracy and fast recognition speed, which improves the positioning speed and accuracy and improves the positioning efficiency in complex environments. BRIEF DESCRIPTION OF THE DRAWINGS

[0022] Figure 1 It is a flow chart of a multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation of the present invention.

[0023] Figure 2 Schematic diagram of the principle of conversion between GNSS altitude and barometric altimeter altitude in the present invention. DETAILED DESCRIPTION

[0024] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0025] The present invention discloses a multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation. By combining multi-sensor data such as an inertial measurement unit, a magnetometer, and an altimeter with low-orbit satellite pseudo-range information, a multi-stage joint detection and anti-interference method for satellite disturbance signals is realized. This method provides a satellite disturbance signal detection and anti-interference algorithm with low algorithm complexity, high recognition accuracy, and fast recognition speed for unmanned platforms in future urban scenarios.

[0026] The overall flow chart of the multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation of the present invention is as follows: Figure 1 As shown, the following steps are included:

[0027] Step 1: Receiver cycle slip detection and pseudorange gross error detection: Perform consistency check on pseudorange, carrier phase, and Doppler observations, and detect pseudorange gross errors and carrier phase cycle slips through Doppler observations.

[0028] Step 2: Enhanced GNSS fault detection using low-orbit satellites: Build a RAIM algorithm combining GNSS observation data with pseudo-range information from low-orbit satellites to increase the GNSS fault detection rate and enhance positioning performance.

[0029] Step 3: Estimate ambiguity parameters using the GNSS relative positioning model: GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method are used to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the ratio test is used to verify the ambiguity to obtain the GNSS high-precision positioning result.

[0030] Step 4: Barometric altimeter-assisted GNSS spoofing detection: After obtaining the GNSS positioning results, the barometric altimeter is used to assist in GNSS spoofing detection. An adaptive threshold detection method based on moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within a time window for spoofing detection.

[0031] Step 5: AHRS-assisted GNSS spoofing interference detection: An AHRS system is constructed by combining magnetometer data with inertial measurement unit data. The receiver's GNSS plane acceleration is compared with the plane acceleration consisting of the east and north directions measured by the accelerometer in the AHRS system. The difference between the GNSS plane acceleration and the AHRS plane acceleration is calculated to see if there is any jump data within the time window, thereby detecting whether there is spoofing interference.

[0032] As a specific example, in step 1, receiver cycle slip detection and pseudorange gross error detection are as follows:

[0033] In complex scenarios, satellite obstruction is severe and navigation signals are fragile. GNSS observations are susceptible to gross errors and cycle slips. Failure to accurately detect and process gross errors in observations in real time will severely impact positioning accuracy and reliability. To address this issue, a consistency check is performed on pseudorange, carrier phase, and Doppler observations. Doppler observations are used to detect pseudorange gross errors and carrier phase cycle slips, as follows:

[0034] When the phase observation value has no cycle slip, the inter-epoch difference can eliminate the ambiguity parameter; when the interval between the previous and subsequent epochs is short, the ionospheric and tropospheric delay errors change little, and the inter-epoch difference can also basically eliminate their influence. Therefore, the pseudorange and phase inter-epoch difference observation equation can be expressed as:

[0035]

[0036] Where, Δ k is the inter-epoch differential operator; dt is the receiver clock error equivalent range error; P and φ are the pseudorange and phase observation values, respectively; ρ is the geometric distance between the satellite and the receiver; ε and ξ are the combined effects of pseudorange and phase noise and unmodeled errors, respectively;

[0037] Since the Doppler observations are less affected by external interference, the consistency test of the Doppler observations between the previous and next epochs will be carried out. The following test quantities are constructed respectively through the differential pseudorange between epochs, the differential phase between epochs, and the average Doppler observation value between two epochs:

[0038]

[0039] Where D P+fd and D φ+fd are the consistency test quantities of pseudorange / Doppler and carrier phase / Doppler respectively; f d is the Doppler frequency shift observation value; λ j is the wavelength; Δt is the time interval between epochs;

[0040] When the interval between the previous and next epochs is short and there are no gross errors or cycle slips in the observations, the consistency check quantity is a stable sequence. The Kalman filter position process noise can be adaptively adjusted in real time based on its size. The carrier phase observation noise adjustment factor is set to:

[0041]

[0042] The pseudorange observation noise adjustment factor is set as:

[0043]

[0044] Where α is the noise adjustment factor; C is the TDCP displacement threshold, and setting C φ 50cm, C P is 5m; σ is the observation noise spectral density, set σ φ is 10cm, σ P is 1m; t is the time interval.

[0045] As a specific example, in step 2, low-orbit satellite enhanced GNSS fault detection is as follows:

[0046] Low-orbit satellites orbit less than 2,000 km above the Earth's surface. Compared to medium- and high-orbit satellites, low-orbit satellites have lower orbital altitudes, higher signal power, and stronger anti-interference capabilities. By combining GNSS observation data with pseudo-range information from low-orbit satellites to build a RAIM algorithm, the GNSS fault detection rate and positioning performance can be enhanced as follows:

[0047] The linearized GNSS observation equation is:

[0048] y=Ax+ε

[0049] Where y is the pseudorange observation; A is the observation matrix; x is the unknown parameter vector; ε is the pseudorange observation error;

[0050] The pseudorange residual vector v is:

[0051]

[0052] Where P is the observation weight matrix, which can be expressed as:

[0053]

[0054] Where, σ L0 is the unit weight standard deviation of the low-orbit satellite, set to 1m; σ G0 is the unit weighted standard deviation of GNSS, set to 3m; and Satellite weight allocation according to satellite elevation angle;

[0055] The cofactor matrix Q of the pseudorange residual vector v for:

[0056] Q v =P -1 -A(A T PA) -1 A T

[0057] Post-test unit weight error for:

[0058]

[0059] Where n is the number of observations; m is the number of satellite systems. Low-orbit satellites are set as system platforms other than BDS to transmit navigation signals and keep consistent with the GPS time and space reference.

[0060] As a specific example, in step 3, the GNSS relative positioning model estimates the ambiguity parameters as follows:

[0061] GNSS positioning is inevitably affected by unmodeled errors such as observation noise, which can lead to abnormal fluctuations in GNSS monitoring sequences. To suppress observation noise and improve monitoring accuracy and reliability, GNSS phase inter-epoch differences, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence methods are used to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the ambiguity is verified through the Ratio test. The GNSS high-precision positioning results are as follows:

[0062] Step 3.1: Use GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method to estimate ambiguity parameters, as follows:

[0063] For the base station r and mobile station u, the reference satellite s1 and non-reference satellite s of the GNSS system p , the double-difference pseudorange and carrier phase observation equations at frequency j are:

[0064]

[0065] Where, represents the inter-station double difference operator; N is the integer ambiguity; ε and ξ are the combined effects of pseudorange, phase noise and unmodeled error, respectively;

[0066] The state equation and observation equation of Kalman filter can be expressed as:

[0067]

[0068] Where, X k and X k-1 are the state vectors of the kth and k-1th epochs, including the mobile station coordinate vector, velocity vector, and double-difference ambiguity vector; Z k is the observation vector of epoch k; H k is the observation matrix; V k is the observation noise; F k,k-1 is the state transfer matrix; W k-1 is the system noise vector;

[0069] The state prediction value is measured and updated according to the observation value of the current epoch to obtain the final estimate of the state vector of the current epoch and its variance-covariance matrix P k ,Right now:

[0070]

[0071] in:

[0072]

[0073] Where I is the unit matrix; K k is the gain matrix, R is the noise adjustment factor α based on the carrier phase observation φ , pseudorange observation noise adjustment factor α P and posterior unit weight error The set measurement noise variance matrix;

[0074] Step 3.2: Use the LAMBDA method to narrow the search space and fix the ambiguity;

[0075] Step 3.3: Due to the inaccuracy of the positioning model and errors in the observation data, especially when the observation noise is large or the satellite geometry is poor, the candidate vectors that meet the requirements in the search space are not necessarily the correct ambiguity solutions, so the results need to be tested and confirmed. Ambiguity confirmation refers to judging whether the searched ambiguity is accepted as the correct ambiguity value based on certain criteria, which is a hypothesis testing problem. Only the ambiguity that has been confirmed to be correct can be used to solve the baseline vector fixed solution. Otherwise, the ambiguity solution fails and only a baseline vector floating-point solution with poor accuracy can be output. The most commonly used test method is based on the ratio test between the second smallest and the minimum a posteriori variance in the fixed solution, that is, the Ratio value test method. The larger the Ratio, the higher the credibility of the search result. In this paper, when the Ratio value is greater than 2, the fixed solution is judged to be the correct ambiguity.

[0076] As a specific example, in step 4, the barometric altimeter assists in GNSS spoofing interference detection as follows:

[0077] After obtaining the GNSS positioning results, the barometric altimeter is used to assist GNSS spoofing detection. The adaptive threshold detection method of moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within the time window to detect spoofing interference. The details are as follows:

[0078] Since the change of air pressure with altitude is affected by many factors, the most common influencing factors are temperature, time, latitude, etc., the international standard atmosphere ISO2533-1975 is used, with mean sea level as the reference altitude, that is, zero altitude. The relationship between the temperature and altitude at a certain point is:

[0079] T r =T0+β(H r -H0)

[0080] Where, H0 is the height of the standard sea level, which is set to 0m; H r is the altitude of the test point; T0 is the temperature at standard sea level, which is 288.15K; β is the temperature gradient, which is -0.0065K / m;

[0081] Use the pressure-height equation to express the relationship between air pressure and altitude:

[0082]

[0083] Where R is the air gas constant, which is 287.05m 2 / (s 2 K); g is the acceleration due to gravity, which is 9.8 m / s 2 ; P0 is the air pressure value at H0, which is 101.325kPa; P r It is the barometric pressure value measured in real time by the barometric altimeter;

[0084] In the measurement process of the barometric altimeter, there will be interference noise. In order to reduce the influence of interference noise in the measurement process of the barometric altimeter, the ARIMA (p, d, q) model noise reduction algorithm after the d-order difference is established to filter the random noise in the non-wide-sense stationary signal. The formula is:

[0085] Φ(L)Δ d x t =δ+Θ(L)ε t

[0086] In the formula, {x t} is a time series that obeys normal distribution and has zero mean; L is the backshift operator, Lx t =x t-1 ;

[0087]

[0088] The GPS / BDS satellite navigation system uses the WGS-84 / CGCS2000 coordinate system. The reference surface in this coordinate system is the reference ellipsoid. The height of the point to be measured is the distance from the point to the ellipsoid along the normal line of the ellipsoid.

[0089] The GNSS altitude and the altitude of the barometric altimeter are used for test verification. Therefore, the altitude conversion is first performed to convert the GNSS geodetic height into altitude, such as Figure 2 As shown, the conversion formula is:

[0090] H p =h ideal +Δh

[0091] Among them, h ideal is the altitude, H p is the geoid height, Δh is the difference between the height of the ellipsoid and the geoid;

[0092] The altitude of the carrier can be obtained by calculating the carrier coordinate position and converting it to the user coordinate system. Therefore, the altitude measurement equation of GNSS can be expressed as:

[0093] h g =h ideal +n s

[0094] Where n s is the white noise generated by measurement error;

[0095] The Euclidean distance D(i) between the GNSS altitude data and the filtered barometric altimeter data is calculated as:

[0096] D(i)=|h g -h q |,i=1,2,…,n

[0097] Where h q Indicates the altitude of the barometric altimeter;

[0098] Calculate the mean η of the moving variance within a time window of width w:

[0099]

[0100] Set the moving variance threshold to:

[0101]

[0102] Where, is the variance of the error sequence when there is no deception signal; T m To set the threshold value when the false alarm rate is 0.135%;

[0103] When the mean of the moving variance is greater than or equal to the moving variance threshold, that is, η≥T h There is a deceptive jamming signal.

[0104] As a specific example, in step 5, AHRS assists GNSS spoofing interference detection as follows:

[0105] An AHRS system is constructed by combining magnetometer data with inertial measurement unit data. The receiver's GNSS acceleration is compared with the acceleration in the east and north directions measured by the accelerometer in the AHRS system. The difference between the GNSS acceleration and the AHRS acceleration is calculated to see if there is a jump in the time window. This allows detection of spoofing interference, as follows:

[0106] The AHRS system uses calibrated output data from the accelerometer, gyroscope, and magnetic sensor as system input, then fuses the data through a Kalman filter to calculate the attitude angle output. The Kalman filter algorithm employed in this invention uses a seven-dimensional vector composed of quaternions and gyroscope drift errors as the state vector. Compared to the Euler angle method, the quaternion approach avoids numerous trigonometric calculations, improves the update rate, ensures real-time performance, and eliminates singular values. Furthermore, introducing gyroscope errors into the state vector allows for real-time prediction of gyroscope drift errors and prompt correction. The outputs of the magnetic sensor and accelerometer are calibrated separately and used as observation vectors.

[0107] Step 5.1: Build the AHRS system by combining the magnetometer data with the inertial measurement unit data as follows:

[0108] The state equation of the AHRS system is:

[0109] x n =F nln-1 x n-1 +v n-1

[0110] Where x n Represents the state vector of the system at time n, including quaternion and gyro drift error; v n-1 is the state noise of the system at time n-1; F nln-1 is the system state transfer matrix;

[0111] When the gyroscope data output rate is high, the angular rate can be approximately considered constant within a short time interval. According to the quaternion differential equation in the formula, the angular rate can be integrated to update the quaternion. The update equation of the attitude quaternion is:

[0112]

[0113] in

[0114]

[0115] Where w b =[ω x ω y ω z ] T is the angular rate output of the gyroscope in the carrier coordinate system; Q n is the quaternion value at time n; Δt is the time interval for gyroscope data update;

[0116] Therefore, the system state transfer matrix is:

[0117]

[0118] Step 5.2: Use the accelerometer in the AHRS system to measure the plane acceleration composed of the east and north directions as follows:

[0119] The observation equation of the AHRS system is expressed as:

[0120] z n =C n x n +w n

[0121] Where z n represents the observation vector of the system at time n, including the three-axis accelerometer velocity and the three-axis geomagnetic field measurement value; w n C is the AHRS system observation noise at time n; n is the AHRS system observation matrix;

[0122] Taking the northeast celestial coordinate system as the geographic coordinate system, the gravity field vector and the geomagnetic field vector in the geographic coordinate system are expressed as:

[0123] a n =[0 0 g] T

[0124] m n =[m nx m ny m nz ] T

[0125] Where g is the acceleration due to gravity; m nx 、m ny and m nz are the measured values ​​of the geomagnetic field in the northeast celestial coordinate system;

[0126] The coordinate transformation matrix from the navigation system to the carrier system is:

[0127]

[0128] Taking the geographic coordinate system as the navigation coordinate system, the gravity field vector in the navigation coordinate system is converted to the carrier coordinate system, and the gravity field vector in the carrier coordinate system is obtained as follows:

[0129]

[0130] Convert the geomagnetic field vector in the navigation coordinate system to the carrier coordinate system, then we have

[0131]

[0132] Therefore, the system observation matrix is:

[0133]

[0134] Based on the Kalman filter algorithm, the accelerometer and magnetometer sensor data are fused to correct the attitude quaternion. The coordinate transformation matrix is ​​constructed based on the attitude quaternion:

[0135]

[0136] The measurement value of the carrier's accelerometer in the northeast direction under the navigation system is:

[0137]

[0138] Step 5.3: Calculate the receiver's GNSS plane acceleration as follows:

[0139] Based on the GNSS high-precision positioning results calculated in step 3, the current east, north, and celestial displacements are calculated as:

[0140]

[0141] Where x k ,y k ,z k and x k-1 ,y k-1 ,z k-1 are the positions of the carrier in the ECEF coordinate system at the current moment and the previous moment respectively; k-1 and L k-1 is the latitude and longitude of the carrier at the previous moment; Δe k ,Δn k ,Δu k is the displacement of the carrier in the east, north and celestial directions;

[0142] The speed of the carrier in the east and north directions is:

[0143]

[0144] Where V e represents the eastward velocity of the carrier, V n Indicates the northbound velocity of the carrier, Δt indicates the time between two positions. The GNSS position calculation frequency is 10 Hz, so Δt is 0.1s;

[0145] The acceleration of the carrier in the east and north directions is further calculated as:

[0146]

[0147] Where A e represents the eastward velocity of the carrier, A n represents the northbound velocity of the carrier;

[0148] Step 5.4: Compare the receiver's GNSS acceleration with the acceleration in the east and north directions measured by the accelerometer in the AHRS system. Calculate whether there is any jump data in the difference between the GNSS acceleration and the AHRS acceleration within the time window to detect spoofing interference. The details are as follows:

[0149] Acceleration error E between GNSS and IMU k for:

[0150] E k =|A GNSS,k -A acc,k |

[0151] Where A GNSS,k A is the vector composed of GNSS east acceleration and north acceleration; acc,k The vector consisting of the east acceleration and north acceleration measured by the IMU in the AHRS;

[0152] Compute the standard deviation of the error series within a time window:

[0153]

[0154] Where, σ k is the standard deviation of the error sequence within the time window; The mean of the error series within the time window;

[0155] The false alarm rate is set to 0.135%, and the judgment threshold is T a =3σ k , when E k ≥T a It is determined that there is a deceptive interference signal at this moment.

[0156] The present invention also provides a low-orbit satellite constellation-assisted multi-stage joint enhanced positioning system, which is used to implement the low-orbit satellite constellation-assisted multi-stage joint enhanced positioning method described above, and the system includes:

[0157] Detection module, used for receiver cycle slip detection and pseudorange gross error detection: consistency check of pseudorange, carrier phase and Doppler observation values, and detection of pseudorange gross errors and carrier phase cycle slips through Doppler observation values;

[0158] Fault detection module for enhanced GNSS fault detection on LEO satellites: RAIM algorithm is constructed by combining GNSS observation data with LEO satellite pseudorange information to enhance GNSS fault detection rate and positioning performance;

[0159] The ambiguity parameter estimation module estimates ambiguity parameters through the GNSS relative positioning model. It uses GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the ratio test is used to verify the ambiguity to obtain GNSS high-precision positioning results.

[0160] The first auxiliary interference detection module uses a barometric altimeter to assist in GNSS spoofing interference detection: After obtaining the GNSS positioning results, the barometric altimeter is used to assist in GNSS spoofing detection. The adaptive threshold detection method of moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within a time window for spoofing interference detection;

[0161] The second auxiliary interference detection module uses AHRS to assist in GNSS deception interference detection: by combining magnetometer data with inertial measurement unit data to build an AHRS system, the receiver's GNSS plane acceleration is compared with the plane acceleration composed of east and north directions measured by the accelerometer in the AHRS system, and the difference between the GNSS plane acceleration and the AHRS plane acceleration is calculated to see whether there is jump data within the time window, thereby detecting whether there is deception interference.

[0162] The present invention also provides a mobile terminal, including a memory, a processor, and a computer program stored in the memory and runnable on the processor, characterized in that when the processor executes the program, the multi-stage joint enhanced positioning method assisted by the low-orbit satellite constellation is implemented.

[0163] The present invention also provides a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the steps in the multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation.

[0164] The above are only preferred embodiments of the present invention. It should be pointed out that for ordinary technicians in this technical field, several improvements and modifications can be made without departing from the principles of the present invention. These improvements and modifications should also be regarded as the scope of protection of the present invention.

Claims

1. A multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation, characterized in that: The following steps are involved: Step 1: Receiver cycle slip detection and pseudorange gross error detection: Perform consistency check on pseudorange, carrier phase, and Doppler observations, and detect pseudorange gross errors and carrier phase cycle slips through Doppler observations. Step 2: Enhanced GNSS fault detection using low-orbit satellites: Build a RAIM algorithm combining GNSS observation data with pseudo-range information from low-orbit satellites to increase the GNSS fault detection rate and enhance positioning performance. Step 3: Estimate ambiguity parameters using the GNSS relative positioning model: GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method are used to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the ratio test is used to verify the ambiguity to obtain the GNSS high-precision positioning result. Step 4: Barometric altimeter-assisted GNSS spoofing detection: After obtaining the GNSS positioning results, the barometric altimeter is used to assist in GNSS spoofing detection. An adaptive threshold detection method based on moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within a time window for spoofing detection. Step 5: AHRS-assisted GNSS spoofing interference detection: An AHRS system is constructed by combining magnetometer data with inertial measurement unit data. The receiver's GNSS plane acceleration is compared with the plane acceleration consisting of the east and north directions measured by the accelerometer in the AHRS system. The difference between the GNSS plane acceleration and the AHRS plane acceleration is calculated to see if there is any jump data within the time window, thereby detecting whether there is spoofing interference.

2. The multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation according to claim 1 is characterized in that: In step 1, the consistency check of pseudorange, carrier phase and Doppler observations is performed, and the Doppler observations are used to detect pseudorange errors and carrier phase cycle slips, as follows: When the phase observation value has no cycle slip, the inter-epoch difference can eliminate the ambiguity parameter; when the interval between the previous and next epochs is shorter than the set value, the change in the ionospheric and tropospheric delay errors is less than the threshold, and the inter-epoch difference eliminates the influence. Therefore, the pseudorange and phase inter-epoch difference observation equation is expressed as: Where, Δ k is the inter-epoch differential operator; dt is the receiver clock error equivalent range error; P and φ are the pseudorange and phase observation values, respectively; ρ is the geometric distance between the satellite and the receiver; ε and ξ are the combined effects of pseudorange and phase noise and unmodeled errors, respectively; Ignoring the external interference of Doppler observations, the consistency test of Doppler observations between the previous and next epochs is performed. The following test quantities are constructed respectively through the differential pseudorange between epochs, the differential phase between epochs, and the average Doppler observation value between two epochs: Where, and are the consistency test quantities of pseudorange / Doppler and carrier phase / Doppler respectively; f d is the Doppler frequency shift observation value; λ j is the wavelength; Δt is the time interval between epochs; When the interval between the previous and next epochs is shorter than the set value and there are no gross errors or cycle slips in the observations, the consistency check quantity is a stable sequence. The position process noise of the Kalman filter is adaptively adjusted in real time based on the size of the consistency check quantity. The carrier phase observation noise adjustment factor is set to: The pseudorange observation noise adjustment factor is set as: Where α is the noise adjustment factor; C is the TDCP displacement threshold, and setting C φ 50cm, C P is 5m; σ is the observation noise spectral density, set σ φ is 10cm, σ P is 1m; t is the time interval.

3. The multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation according to claim 2 is characterized in that: In step 2, the RAIM algorithm is constructed by combining GNSS observation data with pseudorange information of low-orbit satellites to enhance the GNSS fault detection rate and positioning performance, as follows: The linearized GNSS observation equation is: y=Ax+ε Where y is the pseudorange observation; A is the observation matrix; x is the unknown parameter vector; ε is the pseudorange observation error; The pseudorange residual vector v is: Where P is the observation weight matrix, expressed as: Where, σ L0 is the unit weight standard deviation of the low-orbit satellite, set to 1m; σ G0 is the unit weighted standard deviation of GNSS, set to 3m; and Satellite weight allocation according to satellite elevation angle; The cofactor matrix Q of the pseudorange residual vector v for: Q v =P -1 -A(A T PA) -1 A T Post-test unit weight error for: Where n is the number of observations; m is the number of satellite systems. Low-orbit satellites are set as system platforms other than BDS to transmit navigation signals and keep consistent with the GPS time and space reference.

4. The multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation according to claim 3 is characterized in that: In step 3, GNSS phase epoch difference, Kalman filter process noise adaptive adjustment and satellite pseudorange / pseudorange rate confidence method are used to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the Ratio test is used to verify the ambiguity to obtain the GNSS high-precision positioning result, as follows: Step 3.1: Use GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method to estimate ambiguity parameters, as follows: For the base station r and mobile station u, the reference satellite s1 and non-reference satellite s of the GNSS system p , the double-difference pseudorange and carrier phase observation equations at frequency j are: Where Δ▽ represents the inter-station double difference operator; N is the integer ambiguity; ε and ξ are the combined effects of pseudorange and phase noise and unmodeled errors, respectively; The state equation and observation equation of Kalman filter are expressed as: Where, X k and X k-1 are the state vectors of the kth and k-1th epochs, including the mobile station coordinate vector, velocity vector, and double-difference ambiguity vector; Z k is the observation vector of epoch k; H k is the observation matrix; V k is the observation noise; F k,k-1 is the state transfer matrix; W k-1 is the system noise vector; The state prediction value is measured and updated according to the observation value of the current epoch to obtain the final estimate of the state vector of the current epoch and its variance-covariance matrix P k ,Right now: in: Where I is the unit matrix; K k is the gain matrix, R is the noise adjustment factor α based on the carrier phase observation φ , pseudorange observation noise adjustment factor α P and posterior unit weight error The set measurement noise variance matrix; Step 3.2: Use the LAMBDA method to narrow the search space and fix the ambiguity; Step 3.3: Verify the ambiguity through the Ratio test. When the Ratio value is greater than 2, the fixed solution is determined to be the correct ambiguity.

5. The multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation according to claim 4 is characterized in that: In step 4, after obtaining the GNSS positioning results, the barometric altimeter is used to assist the GNSS in spoofing detection. The adaptive threshold detection method of moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within the time window to perform spoofing interference detection. The details are as follows: Using the International Standard Atmosphere ISO2533-1975, with mean sea level as the reference height, i.e. zero altitude, the relationship between the temperature and altitude at a certain point is: T r =T0+β(H r -H0) Where, H0 is the height of the standard sea level, which is set to 0m; H r is the altitude of the test point; T0 is the temperature at standard sea level, which is 288.15K; β is the temperature gradient, which is -0.0065K / m; Use the pressure-height equation to express the relationship between air pressure and altitude: Where R is the air gas constant, which is 287.05m 2 / (s 2 K); g is the acceleration due to gravity, which is 9.8 m / s 2 ; P0 is the air pressure value at H0, which is 101.325kPa; P r It is the barometric pressure value measured in real time by the barometric altimeter; In order to reduce the influence of interference noise in the measurement process of the barometric altimeter, the ARIMA (p, d, q) model noise reduction algorithm after d-order difference is established to filter the random noise in the non-wide sense stationary signal. The formula is: F(L)D d x t =δ+Θ(L)ε t In the formula, {x t } is a time series that obeys normal distribution and has zero mean; L is the backshift operator, Lx t =x t-1 ; The GPS / BDS satellite navigation system uses the WGS-84 / CGCS2000 coordinate system. The reference surface in this coordinate system is the reference ellipsoid. The height of the point to be measured is the distance from the point to the ellipsoid along the normal line of the ellipsoid. The GNSS altitude and the altitude of the barometric altimeter are used for test verification. Therefore, the altitude conversion is first performed to convert the GNSS geodetic height to the altitude. The conversion formula is: H p =h ideal +Δh Among them, h ideal is the altitude, H p is the geoid height, Δh is the difference between the height of the ellipsoid and the geoid; The altitude of the carrier is obtained by calculating the carrier coordinate position and converting it to the user coordinate system. Therefore, the altitude measurement equation of GNSS is expressed as: h g =h ideal +n s Where n s is the white noise generated by measurement error; The Euclidean distance D(i) between the GNSS altitude data and the filtered barometric altimeter data is calculated as: D(i)=|h g -h q |,i=1,2,…,n Where h q Indicates the altitude of the barometric altimeter; Calculate the mean η of the moving variance within a time window of width w: Set the moving variance threshold to: Where, is the variance of the error sequence when there is no deception signal; T m To set the threshold value when the false alarm rate is 0.135%; When the mean of the moving variance is greater than or equal to the moving variance threshold, that is, η≥T h There is a deceptive jamming signal.

6. The multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation according to claim 5, characterized in that: In step 5, the AHRS system is constructed by combining the magnetometer data with the inertial measurement unit data as follows: Step 5.1: Build the AHRS system by combining the magnetometer data with the inertial measurement unit data as follows: The state equation of the AHRS system is: x n =F nln-1 x n-1 +v n-1 Where x n Represents the state vector of the system at time n, including quaternion and gyro drift error; v n-1 is the state noise of the system at time n-1; F nln-1 is the system state transfer matrix; When the gyroscope data output rate is higher than the set value, the angular rate is considered unchanged within the set time interval. According to the quaternion differential equation in the formula, the angular rate is integrated to update the quaternion. The update equation of the attitude quaternion is: in Where w b =[ω x ω y ω z ] T is the angular rate output of the gyroscope in the carrier coordinate system; Q n is the quaternion value at time n; Δt is the time interval for gyroscope data update; Therefore, the system state transfer matrix is: Step 5.2: Use the accelerometer in the AHRS system to measure the plane acceleration composed of the east and north directions as follows: The observation equation of the AHRS system is expressed as: z n =C n x n +w n Where z n represents the observation vector of the system at time n, including the three-axis accelerometer velocity and the three-axis geomagnetic field measurement value; w n C is the AHRS system observation noise at time n; n is the AHRS system observation matrix; Taking the northeast celestial coordinate system as the geographic coordinate system, the gravity field vector and the geomagnetic field vector in the geographic coordinate system are expressed as: to n =[00g] T m n =[m nx m ny m nz ] T Where g is the acceleration due to gravity; m nx 、m ny and m nz are the measured values ​​of the geomagnetic field in the northeast celestial coordinate system; The coordinate transformation matrix from the navigation system to the carrier system is: Taking the geographic coordinate system as the navigation coordinate system, the gravity field vector in the navigation coordinate system is converted to the carrier coordinate system, and the gravity field vector in the carrier coordinate system is obtained as follows: Convert the geomagnetic field vector in the navigation coordinate system to the carrier coordinate system, then we have Therefore, the system observation matrix is: Based on the Kalman filter algorithm, the accelerometer and magnetometer sensor data are fused to correct the attitude quaternion. The coordinate transformation matrix is ​​constructed based on the attitude quaternion: The measurement value of the carrier's accelerometer in the northeast direction under the navigation system is: Step 5.3: Calculate the receiver's GNSS plane acceleration as follows: Based on the GNSS high-precision positioning results calculated in step 3, the current east, north, and celestial displacements are calculated as: Where x k ,y k ,z k and x k-1 ,y k-1 ,z k-1 are the positions of the carrier in the ECEF coordinate system at the current moment and the previous moment respectively; k-1 and L k-1 is the latitude and longitude of the carrier at the previous moment; Δe k ,Δn k ,Δu k is the displacement of the carrier in the east, north and celestial directions; The speed of the carrier in the east and north directions is: Where V e represents the eastward velocity of the carrier, V n Indicates the northbound velocity of the carrier, Δt indicates the time between two positions. The GNSS position calculation frequency is 10 Hz, so Δt is 0.1s; The acceleration of the carrier in the east and north directions is further calculated as: Where A e represents the eastward velocity of the carrier, A n Indicates the northbound velocity of the carrier.

7. The multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation according to claim 6, characterized in that: In step 5, the receiver's GNSS acceleration is compared with the acceleration in the east and north directions measured by the accelerometer in the AHRS system. The difference between the GNSS acceleration and the AHRS acceleration is calculated to see if there is any jump data within the time window, thereby detecting the presence of spoofing interference. The details are as follows: Acceleration error E between GNSS and IMU k for: AND k =|A GNSS,k -THE acc,k | Where A GNSS,k A is the vector composed of GNSS east acceleration and north acceleration; acc,k The vector consisting of the east acceleration and north acceleration measured by the IMU in the AHRS; Compute the standard deviation of the error series within a time window: Where, σ k is the standard deviation of the error sequence within the time window; The mean of the error series within the time window; The false alarm rate is set to 0.135%, and the judgment threshold is T a =3σ k , when E k ≥T a It is determined that there is a deceptive interference signal at this moment.

8. A multi-stage joint enhanced positioning system assisted by a low-orbit satellite constellation, characterized in that: The system is used to implement the multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation according to any one of claims 1 to 7, and the system includes: Detection module, used for receiver cycle slip detection and pseudorange gross error detection: consistency check of pseudorange, carrier phase and Doppler observation values, and detection of pseudorange gross errors and carrier phase cycle slips through Doppler observation values; Fault detection module for enhanced GNSS fault detection on LEO satellites: RAIM algorithm is constructed by combining GNSS observation data with LEO satellite pseudorange information to enhance GNSS fault detection rate and positioning performance; The ambiguity parameter estimation module estimates ambiguity parameters through the GNSS relative positioning model. It uses GNSS phase epoch difference, Kalman filter process noise adaptive adjustment, and satellite pseudorange / pseudorange rate confidence method to estimate ambiguity parameters. The LAMBDA method is used to fix the ambiguity, and the ratio test is used to verify the ambiguity to obtain GNSS high-precision positioning results. The first auxiliary interference detection module uses a barometric altimeter to assist in GNSS spoofing interference detection: After obtaining the GNSS positioning results, the barometric altimeter is used to assist in GNSS spoofing detection. The adaptive threshold detection method of moving variance is used to compare the moving variance value of the distance difference between the barometric altimeter and the GNSS altitude within a time window for spoofing interference detection; The second auxiliary interference detection module uses AHRS to assist in GNSS deception interference detection: by combining magnetometer data with inertial measurement unit data to build an AHRS system, the receiver's GNSS plane acceleration is compared with the plane acceleration composed of east and north directions measured by the accelerometer in the AHRS system, and the difference between the GNSS plane acceleration and the AHRS plane acceleration is calculated to see whether there is jump data within the time window, thereby detecting whether there is deception interference.

9. A mobile terminal comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein: When the processor executes the program, the multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation as described in any one of claims 1 to 7 is implemented.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the program is executed by a processor, the steps of the multi-stage joint enhanced positioning method assisted by a low-orbit satellite constellation are implemented as described in any one of claims 1 to 7.

Images