Network RTK user request fast matching method based on row and column number mechanism

Abstract

The application discloses a network RTK user request fast matching method based on a row-column number mechanism, comprising the following steps: acquiring satellite ephemeris, reference station original observation values and reference station coordinate data; a grid division module divides regional grid points according to the minimum value of the left lower corner, the maximum value of the right upper corner and the division interval of a reference station group; Delaunay triangular network elements are divided according to the accurate coordinates of the reference station group, and the gravity center coordinates of each triangle are calculated; the ambiguity between reference stations in each triangular network element is calculated according to the division result; the space correlation error values between reference stations are calculated and generated; the pseudo-range observation values and the carrier phase observation values of virtual reference stations are interpolated; all difference correction products are sent to a broadcasting platform; when a user requests product service, the broadcasting platform firstly verifies the user, refuses service if the verification fails, and performs product fast matching if the verification passes; after the matching is completed, the product is broadcast to the user to complete high-precision positioning.

Description

Technical Field

[0001] The present invention belongs to the field of satellite navigation and positioning technology, in particular to regional satellite navigation differential enhancement technology, and specifically relates to a network RTK user request fast matching method based on a row and column number mechanism. Background Art

[0002] In recent years, my country's Beidou satellite navigation system has seen rapid development. It has played a vital role in areas of national economy and livelihood, such as engineering surveying, geological exploration, deformation monitoring, vehicle navigation, and smart logistics. With the development of human society, people's demand for location-based services is increasing. Currently, satellite navigation positioning can only provide a horizontal positioning accuracy of approximately 3 meters. Differential positioning technology can significantly improve positioning accuracy. Real-time carrier phase differential positioning (RTK) is currently the most widely used high-precision positioning technology. It uses double-difference between a base station and a rover to eliminate common-view errors, providing high-precision positioning results for rover users. However, as the distance between the base station and the rover increases, the spatial correlation between the two weakens, causing conventional RTK differential positioning solutions to lose accuracy due to distance limitations.

[0003] Because conventional RTK positioning technology is limited by the distance between the reference station and the rover, network RTK positioning technology has been proposed. Network RTK positioning technology divides all reference stations in a region into a grid. By calculating the comprehensive errors of the ionosphere, troposphere, and ephemeris of each reference station in the grid, it interpolates the error corrections for each virtual point and calculates the pseudorange and carrier phase observation values ​​based on this. The virtual observation values ​​of each virtual point are sent to the broadcast platform after standard encoding. The rover users in the region receive the virtual point differential corrections sent by the broadcast platform and calculate the terminal's high-precision position coordinates through differential calculation. Network RTK technically makes up for the shortcomings of conventional RTK technology, enabling cm-level real-time dynamic positioning solutions within the region, expanding the service scope of high-precision positioning. Because virtual reference station technology (VRS) has great advantages in terms of operating distance, accuracy, reliability, and real-time performance, it has become the most widely used and successful solution in the field of network RTK. Therefore, the present invention also adopts VRS technology to implement network RTK differential service functions.

[0004] Network RTK differential services provide differential services using virtual observations calculated from virtual points within a grid. When a user requests differential data products from the service center, the service center locates the grid point based on the approximate coordinates provided by the user, finds the grid point closest to the user, and sends the coordinates of that grid point and the differential data product to the user. The user then receives the product for high-precision positioning. For large areas with small grid intervals, the number of grid points within the area is very large. Traversing the grid point database based on the user's location is not only inefficient but also introduces delays to regional services. Summary of the Invention

[0005] The object of the present invention is to provide

[0006] The technical solution for achieving the purpose of the present invention is:

[0007] A method for quickly matching network RTK user requests based on a row and column number mechanism, comprising:

[0008] Step 1: Obtain satellite ephemeris, base station original observation values ​​and base station coordinate data;

[0009] Step 2: Delineate a rectangular area based on the base station coordinates and calculate the minimum latitude and longitude value BL of the lower left corner of the rectangular area. min and the maximum latitude and longitude value BL in the upper right corner max , and divide the grid points in sequence according to the interval inval;

[0010] Step 3: Divide the reference stations in the area into triangulated elements in the form of Delaunay triangulation elements;

[0011] Step 4: Calculate the inter-reference station ambiguity within one of the triangulated network elements, first perform wide lane ambiguity Fixed, then use the ionosphere-free combined model to fix the wide lane and L1 ambiguity Restore its integer properties and finally determine the L2 ambiguity using the fixed wide lane and L1 ambiguity size;

[0012] Step 5: Calculate the spatial correlation error between reference stations. When the ambiguity between reference stations is fixed, generate the double-difference phase observation value Δ▽φ, double-difference pseudorange observation value Δ▽ρ, double-difference ionospheric delay Δ▽I, double-difference tropospheric delay Δ▽T and orbit comprehensive error Δ▽O between each reference station in real time.

[0013] Step 6: Use linear interpolation to interpolate the double-difference ionospheric residual Δ▽I between the main reference station and the virtual reference station at the grid point in each triangulation element. i,VRS , double-difference tropospheric residual Δ▽T i,VRS and double difference orbit residual Δ▽O i,VRS On this basis, the pseudorange observation value P(VRS) and carrier phase observation value of the virtual reference station are calculated through the pseudorange observation value and carrier phase observation value of the main reference station.

[0014] Step 7, repeating step 6, sequentially calculating the pseudorange observation values ​​and carrier phase observation values ​​of all grid point virtual reference stations in the triangulation element, and mapping the pseudorange observation values ​​and carrier phase observation values ​​of the grid point virtual reference station to the grid point row and column numbers one by one;

[0015] Step 8: Repeat steps 4-7 to calculate the pseudorange observation values ​​and carrier phase observation values ​​of the virtual reference station at all grid points in sequence;

[0016] Step 9: When a user makes a service request, first verify their identity. If the verification fails, the process ends immediately. If the verification passes, the user will be assigned a service request based on the latitude and longitude information (B i ,L i ), and the minimum longitude and latitude of the region BL min (B min , L min ) and divide it by the interval inval to calculate the grid point position closest to the user position (H i , R i );

[0017] Step 10, grid points (H i , R i )The differential correction product mapped by the GNSS is sent to the user for positioning.

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

[0019] When a user requests differential data from the service center, the service center calculates the row and column numbers based on the user's approximate coordinates and the latitude and longitude of the lower left corner of the region, dividing the difference by the partition interval inval. The resulting values ​​are then sent directly to the user along with the grid point coordinates and differential correction products corresponding to these values. This method maps row and column numbers to correction products during differential correction generation, efficiently and effectively finding the grid point closest to the user without traversing the grid point database. This improves the service efficiency of regional RTK software, especially in large regions with small partition intervals, where the efficiency of matching user requests is even more pronounced.

[0020] The present invention is further described in detail below with reference to the accompanying drawings. BRIEF DESCRIPTION OF THE DRAWINGS

[0021] Figure 1 Schematic diagram of a fast matching method for service requests based on row and column numbers.

[0022] Figure 2 This is a schematic diagram of the system flow adopted by the present invention. DETAILED DESCRIPTION

[0023] Combine Figure 1 、 Figure 2 The system flow of a network RTK user request fast matching method based on a row and column number mechanism in this embodiment mainly includes steps such as data access, data solution and data broadcast service;

[0024] Step 1: Obtain satellite ephemeris, original observation values ​​of the reference station, and precise coordinate data of the reference station, send these data to the software data receiving module, decode the ephemeris and original observation data, and store the precise coordinates of the reference station;

[0025] All base station data in the region are connected to the local CASTER server via the Ntrip protocol, and the precise coordinates of each base station are obtained after calculation. The data stream is accessed in real time, the ephemeris and raw observation data are decoded, and the precise coordinates of the base station are stored;

[0026] Step 2: Calculate the minimum latitude and longitude value BL of the lower left corner of the rectangular area based on the regional reference station group min and the maximum latitude and longitude value BL in the upper right corner max , and divide the grid points GW in sequence according to a certain interval inval (i,j) , the maximum row and column numbers are (H max , R max );

[0027] According to the precise coordinates of the base station, the rectangular area is delineated and the minimum longitude and latitude value BL of the lower left corner of the rectangular area is calculated. min and the maximum latitude and longitude value BL in the upper right corner max , divide the regional grid points GW in sequence according to a certain interval inval (i,j) , the maximum row and column numbers are (H max , R max ), according to formula (1) and (2), we can obtain BL min , BL max , according to formula (3) and (4), we can get the maximum row and column number (H max , R max ), and the grid point longitude and latitude GW is obtained according to formula (5) (i,j) .

[0028] BL min =MIN(B i , L i ) (1)

[0029] BL max =MAX(B i , L i ) (2)

[0030]

[0031]

[0032] GW (i,j) =(MIN(B i )+i*inval,MIN(L j)+j*inval)(0≤i≤H max ,0≤j≤R max )(5)

[0033] Among them, MIN() is the minimum function, MAX() is the maximum function, int[] is the integer function, B i Indicates latitude, L i represents longitude, i represents row number, j represents column number, Hmax represents the maximum row number, and Rmax represents the maximum column number.

[0034] Step 3: Divide the reference stations in the area into triangulated elements in the form of Delaunay triangulation elements;

[0035] Step 4: Calculate the inter-reference station ambiguity within one of the triangulated network elements, first perform wide lane ambiguity Fixed, then use the ionosphere-free combined model to fix the wide lane and L1 ambiguity Restore its integer properties and finally determine the L2 ambiguity using the fixed wide lane and L1 ambiguity size;

[0036] To calculate the ambiguity between reference stations within a triangulation network element, the "three-step method" is generally used.

[0037] The first step is to use the wide lane ambiguity N as the wavelength is longer and less affected by the atmospheric delay error. i wide It is relatively easy to fix the wide lane ambiguity using MW combination. The observation value of MW combination is shown in formula (6):

[0038]

[0039] Where, Φ MW is the MW combined carrier observation value, f1 and f2 are the L1 and L2 frequencies, λ1 and λ2 are the L1 and L2 wavelengths, φ1 and φ2 are the L1 and L2 carrier observation values, P1 and P2 are the L1 and L2 pseudorange observation values, λ wide is the wide-lane wavelength, N wide is the wide lane ambiguity, ε P is the pseudorange noise.

[0040] Then we can get the MW combined double difference observation equation as follows:

[0041] Δ▽Φ MW =λ wide Δ▽N wide +Δ▽ε P (7)

[0042] Where Δ▽ is the double difference symbol. It can be seen from formula (3) that when using the MW combined double difference observation value to solve the wide lane ambiguity, it is only affected by the pseudorange observation noise. Therefore, in the data processing process, in order to calculate the correct and reliable wide lane ambiguity, the average value of multiple epochs is generally used to solve the wide lane ambiguity.

[0043] In the second step, the L1 ambiguity is estimated by adaptive Kalman filtering using the ionosphere-free combined model. Floating point solution. The L1 band ambiguity calculation formula is:

[0044]

[0045] Where Δ▽N1 is the L1 double difference ambiguity, Δ▽N W is the wide lane ambiguity, λ n ,λ w are the narrow lane and wide lane wavelengths, respectively; Δ▽ρ, Δ▽O, Δ▽T, and Δ▽M are the double-difference satellite-to-ground distance, double-difference orbit error, double-difference tropospheric delay, and double-difference multipath effect, respectively.

[0046] The third step is to search and fix the L1 ambiguity through the LAMBDA algorithm. After the L1 ambiguity is fixed, the L2 ambiguity is obtained using the known wide lane ambiguity. Fixed solution;

[0047] Step 5. Calculate the spatial correlation error between reference stations. When the ambiguity between reference stations is fixed, generate the double-difference phase observation value Δ▽φ, double-difference pseudorange observation value Δ▽ρ, double-difference ionospheric delay Δ▽I, double-difference tropospheric delay Δ▽T and orbit comprehensive error Δ▽O between each reference station in real time. Generally, the tropospheric delay and orbit comprehensive error are calculated together.

[0048] After the ambiguity between the reference stations is fixed, the geometry-free observation equation is constructed according to the carrier observation equation. The spatial error values ​​between the reference stations, such as the double-difference ionospheric residual Δ▽I, the double-difference tropospheric residual Δ▽T, and the double-difference orbit residual Δ▽O, are generated by real-time calculation. Δ▽I is calculated by equation (9). When generating the double-difference tropospheric delay and orbit error between the reference stations in real-time calculation, they are generally calculated together and calculated by equation (10).

[0049]

[0050]

[0051] in are the carrier phase observations at L1 and L2, N1 and N2 are the ambiguities at L1 and L2, c is the speed of light, and λ1 and λ2 are the wavelengths at L1 and L2.

[0052] Step 6: Use linear interpolation to interpolate the double-difference ionospheric residual Δ▽I between the main reference station and the virtual reference station at the grid point in each triangulation element. i,VRS , double-difference tropospheric residual Δ▽T i,VRS and double difference orbit residual Δ▽O i,VRS On this basis, the pseudorange observation value P(VRS) and carrier phase observation value of the virtual reference station are calculated through the pseudorange observation value and carrier phase observation value of the main reference station.

[0053] According to the grid point coordinates, the double-difference ionospheric residual Δ▽I between the main reference station in each triangulated network element and the virtual reference station at the grid point is interpolated by linear interpolation. i,VRS , double-difference tropospheric residual Δ▽T i,VRS and double difference orbit residual Δ▽O i,VRS On this basis, the pseudorange observation values ​​and carrier phase observation values ​​of the virtual reference station are calculated through the pseudorange observation values ​​and carrier phase observation values ​​of the main reference station;

[0054] The plane interpolation model is the inverse distance weighted interpolation method, as shown in formula (6), (X M ,Y M ,Z M ) is the coordinate value of the virtual reference station, (X i ,Y i ,Z i ) is the coordinate value of the network element base station.

[0055]

[0056] Among them, d i is the distance between the virtual reference station and the network element base station, b i is the coefficient, b is the sum of the coefficients, a i is the coefficient weight.

[0057]

[0058] Where, is the carrier phase observation value of the grid point virtual reference station, is the carrier phase observation value of L1 and L2 corresponding to the main reference station in the triangulation network element (i = 1, 2, corresponding to L1 and L2 respectively), Δρ i,VRS is the difference between the two stations, λ i are the wavelengths of the L1 and L2 frequency points (i=1, 2, corresponding to L1 and L2 respectively). is the multipath error and measurement noise.

[0059]

[0060] Where P(VRS) is the pseudorange observation value of the grid point virtual reference station, Pi is the L1 and L2 pseudorange observation values ​​corresponding to the main reference station in the triangulation network element (i = 1, 2, corresponding to L1 and L2 respectively), Δρ i,VRS is the difference between the two stations, Δ▽M P,i,VRS , ε P,i,VRS is the multipath error and measurement noise.

[0061] Step 7, repeating step 6, sequentially calculating the pseudorange observation values ​​and carrier phase observation values ​​of all grid point virtual reference stations in the triangulation element, and mapping the pseudorange observation values ​​and carrier phase observation values ​​of the grid point virtual reference station to the grid point row and column numbers one by one;

[0062] Step 8: Repeat steps 4-7 to calculate the pseudorange observation values ​​and carrier phase observation values ​​of the virtual reference station at all grid points in sequence;

[0063] According to the above steps, the pseudorange observation values ​​and carrier phase observation values ​​of the virtual reference station at all grid points are calculated in sequence, and the grid point row and column numbers are mapped one by one with the pseudorange observation values ​​and carrier phase observation value products of the virtual reference station at the grid point to form an array.

[0064] Step 9: When a user makes a service request, first verify their identity. If the verification fails, the process ends immediately. If the verification passes, the user will be assigned a service request based on the latitude and longitude information (B i ,L i ), and the minimum longitude and latitude of the region BL min (B min , L min ) and divide it by the interval inval to calculate the grid point position closest to the user position (H i , R i );

[0065] When a user makes a service request, the user's identity is verified first. If the verification fails, the service ends immediately. If the verification passes, the user uploads the approximate coordinates of the latitude and longitude information (B i ,L i ), and the minimum longitude and latitude of the region BL min (B min , L min ) and divide it by the interval inval, and calculate the grid point position closest to the user position (H i , R i );

[0066]

[0067] Step 10: Set the grid point position (H i , R i ) and the corresponding differential correction product is sent to the user for high-precision positioning.

Claims

1. A method for fast matching of network RTK user requests based on row and column number mechanism, characterized in that: include: Step 1: Obtain satellite ephemeris, base station original observation values ​​and base station coordinate data; Step 2: Delineate a rectangular area based on the base station coordinates and calculate the minimum latitude and longitude value BL of the lower left corner of the rectangular area. min and the maximum latitude and longitude value BL in the upper right corner max , and divide the grid points in sequence according to the interval inval; Step 3: Divide the reference stations in the area into triangulated elements in the form of Delaunay triangulation elements; Step 4: Calculate the inter-reference station ambiguity within one of the triangulated network elements, first perform wide lane ambiguity Fixed, then use the ionosphere-free combined model to fix the wide lane and L1 ambiguity Restore its integer properties and finally determine the L2 ambiguity using the fixed wide lane and L1 ambiguity size; Step 5: Calculate the spatial correlation error between the reference stations. When the ambiguity between the reference stations is fixed, generate the double-difference phase observation value between each reference station in real time. Double-difference pseudorange observations Double-difference ionospheric delay Double-difference tropospheric delay and orbital errors Step 6: Use linear interpolation to interpolate the double-difference ionospheric residuals between the main reference station and the virtual reference station at the grid point in each triangulation element. Double-difference tropospheric residual and double-difference orbital residuals On this basis, the pseudorange observation value P(VRS) and carrier phase observation value of the virtual reference station are calculated through the pseudorange observation value and carrier phase observation value of the main reference station. Step 7, repeating step 6, sequentially calculating the pseudorange observation values ​​and carrier phase observation values ​​of all grid point virtual reference stations in the triangulation element, and mapping the pseudorange observation values ​​and carrier phase observation values ​​of the grid point virtual reference station to the grid point row and column numbers one by one; Step 8: Repeat steps 4-7 to calculate the pseudorange observation values ​​and carrier phase observation values ​​of the virtual reference station at all grid points in sequence; Step 9: When a user makes a service request, first verify their identity. If the verification fails, the process ends immediately. If the verification passes, the user will be assigned a service request based on the latitude and longitude information (B i ,L i ), and the minimum longitude and latitude of the region BL min (B min , L min ) and divide it by the interval inval to calculate the grid point position closest to the user position (H i , R i ); Step 10, grid points (H i , R i )The differential correction product mapped by the GNSS is sent to the user for positioning.

2. The method for fast matching of network RTK user requests based on row and column number mechanism according to claim 1, characterized in that: Maximum row and column number (H max , R max ) is calculated as: Among them, MIN() is the minimum function, MAX() is the maximum function, int[] is the integer function, B i Indicates latitude, L i Indicates longitude.

3. The method for fast matching of network RTK user requests based on row and column number mechanism according to claim 2, characterized in that: Grid point longitude and latitude GW (i,j) , calculated by the following formula: GW (i,j) =(MIN(B i )+i*selection,MIN(L j )+j*inval)(0≤i≤H max ,0≤j≤R max ) Where i represents the row number and j represents the column number.

4. The method for fast matching of network RTK user requests based on row and column number mechanism according to claim 1, characterized in that: Step 4 specifically includes the following steps: In the first step, the wide lane ambiguity is fixed using the MW combination to obtain the MW combination double-difference observation equation; In the second step, the L1 ambiguity is estimated by adaptive Kalman filtering using the ionosphere-free combined model. floating point solution; The third step is to search and fix the L1 ambiguity through the LAMBDA algorithm. After the L1 ambiguity is fixed, the L2 ambiguity is obtained using the known wide lane ambiguity. Fixed solution.

5. The method for fast matching of network RTK user requests based on row and column number mechanism according to claim 4, characterized in that: The calculation formula for L1 band ambiguity is: In the formula is the L1 double-difference ambiguity, is the wide lane ambiguity, λ n ,λ w are the wavelengths of narrow lane and wide lane respectively, They are double-difference satellite-to-ground distance, double-difference orbit error, double-difference tropospheric delay, and double-difference multipath effect; f1 and f2 are L1 and L2 frequencies, and φ1 and φ2 are L1 and L2 carrier observation values.

6. The method for fast matching of network RTK user requests based on row and column number mechanism according to claim 5, characterized in that: Calculated by formula (9), the double-difference tropospheric delay and orbit error between reference stations are calculated in real time and combined together to calculate by formula (10): N1 and N2 are L1 and L2 ambiguities respectively, c is the speed of light, are the carrier phase observation values ​​of L1 and L2 frequency points respectively, and λ1 and λ2 are the wavelengths of L1 and L2 respectively.

7. The method for fast matching of network RTK user requests based on row and column number mechanism according to claim 5, characterized in that: The grid point location closest to the user's location (H i , R i )for: Among them (B i ,L i ) is the latitude and longitude information uploaded by the user, (B min , L min ) is the minimum longitude and latitude of the area, and int[] is the rounding function.

Images