Unmanned aerial vehicle PPK post-processing method and device integrating filtering and ambiguity processing

Abstract

The invention provides an unmanned aerial vehicle PPK post-processing method and device fusing filtering and ambiguity processing, the positioning data is filtered for multiple times, the multiple filtering results are fused in the subsequent processing step, the result types are fused in the fusion process, the numerical values are fused, and the accuracy of the result types is improved. Compared with the prior art that only result types are fused, the accuracy and anti-interference performance of the fusion result are effectively improved. Meanwhile, composite fixation is carried out in the ambiguity fixation step, and the fixation rate is obviously superior to that of an existing scheme.

Description

Technical Field

[0001] The present invention relates to the field of data processing, in particular to the field of unmanned aerial vehicle positioning, and discloses a method and device for unmanned aerial vehicle PPK post-processing that integrates filtering and ambiguity processing. Background Art

[0002] UAV PPK positioning technology is an advanced positioning method that allows UAVs to obtain high-precision position information during flight. PPK technology, namely Post-Processing Kinematic positioning technology, synchronously receives GPS satellite signals through a base station and at least one mobile station (UAV), and then processes the data in a computer to determine the precise three-dimensional coordinates of the mobile station.

[0003] However, the current solution only has the result fusion of two-way filtering, which cannot really give play to the advantages of post-processing. At the same time, the current solution ambiguity fixation is only a simple partial ambiguity fixation, and there is no composite fixation method, so the fixation rate is low. In addition, the current solution result fusion technology is limited to the comparison of result types, and does not examine more detailed parameters. Summary of the invention

[0004] The present disclosure at least provides a UAV PPK post-processing method and device that integrates filtering and fuzzy processing to solve at least one of the aforementioned technical problems.

[0005] According to one aspect of the present disclosure, a UAV PPK post-processing method for fusion filtering and ambiguity processing is provided, comprising:

[0006] Acquire first positioning data received by the base station, and convert the first positioning data and second positioning data into data in a predetermined data format; wherein the second positioning data is positioning data received by the drone;

[0007] Performing N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format; wherein N is a positive integer greater than or equal to 4;

[0008] Perform cycle slip detection processing, floating point solution processing, integer ambiguity fixation processing, baseline solution, coordinate conversion and adjustment processing on the first positioning data and the second positioning data after each filtering, and obtain the UAV positioning result corresponding to each filtering;

[0009] For the drone positioning results corresponding to each filtering, if the types of N drone positioning results are the same, the types of the N drone positioning results are used as the initial positioning result type; if the types of two drone positioning results are the same, the types of the two drone positioning results are used as the initial positioning result type; if the types of N drone positioning results are different, the drone positioning results are sorted according to the filtering order, and the types of the two middle drone positioning results are used as the initial positioning result type; the values ​​of the N drone positioning results are merged to obtain the target value of the drone positioning result; the types of the N drone positioning results are upgraded and downgraded to obtain the target positioning result type of the drone positioning result;

[0010] Among them, the types of N drone positioning results are upgraded and downgraded, including:

[0011] Subtract the values ​​of the positioning results of two drones in the types of N drone positioning results to obtain multiple difference values. If all the obtained difference values ​​exceed the initial threshold, the type of the initial positioning result is downgraded;

[0012] Alternatively, the values ​​of the positioning results of the two drones are subtracted from each other, and if both differences exceed the initial threshold, the type of the initial positioning result is downgraded;

[0013] Among them, the numerical values ​​of N UAV positioning results are fused, including:

[0014] According to the values ​​of N UAV positioning results, and using the acceptance residual sigma as the inverse proportional weighting, the target value of the UAV positioning result is obtained;

[0015] Alternatively, the target value of the UAV positioning result is obtained by taking the positioning results of two UAVs of the same type and using the posterior residual sigma as an inverse proportional weight.

[0016] In a possible implementation manner, performing N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format includes:

[0017] The first positioning data and the second positioning data converted according to the predetermined data format are subjected to four filtering processes.

[0018] In a possible implementation, the following steps are also included:

[0019] The first filtering is a forward Kalman filter, and the first epoch is the initial state, and the initial state matrix and the initial variance matrix are parameters of the initial state; wherein the epoch is a frame of data.

[0020] In a possible implementation, the following steps are also included:

[0021] The second filtering is Kalman filtering, and the last epoch is the initial state, and the initial state matrix and initial variance matrix are the parameters of the initial state.

[0022] In a possible implementation, the following steps are also included:

[0023] The third filtering is Kalman filtering, in which the first epoch of the third filtering is the last epoch, and the initial state matrix and the initial variance matrix are the results of the last epoch of the first filtering. The third filtering is performed halfway, and the third filtering ends at the middle time point epoch.

[0024] In a possible implementation, the following steps are also included:

[0025] The fourth filtering is Kalman filtering, in which the first epoch of the fourth filtering is the first epoch, and the initial state matrix and the initial variance matrix are the results of the last epoch of the second filtering. The fourth filtering is only performed halfway, and the fourth filtering ends at the middle time point epoch.

[0026] In a possible implementation manner, the integer ambiguity fixing process includes:

[0027] Screening satellite frequencies of all full-frequency navigation satellites, selecting a navigation satellite of the first frequency, and performing ambiguity search and fixing on the corresponding ambiguity floating-point solutions obtained by floating-point solution processing;

[0028] Arrange the parameters of all full-frequency points, assign scores to the parameters according to weights, arrange all scores, perform star-kicking operations according to multiple preset ratios, and perform ambiguity search and fixation on the floating-point solutions of the ambiguities obtained by the floating-point solution processing corresponding to the remaining navigation satellites;

[0029] For the remaining navigation satellites, single satellites are eliminated, and single satellites are eliminated cyclically until the number of eliminated navigation satellites accounts for a first preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of ambiguities obtained by floating-point solution processing corresponding to the remaining navigation satellites; for the navigation satellites that cannot be fixed among the remaining navigation satellites, two satellites are eliminated, and so on, until the number of eliminated navigation satellites accounts for a second preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of ambiguities obtained by floating-point solution processing corresponding to the remaining navigation satellites; wherein the parameters include navigation satellite altitude angle, signal-to-noise ratio, azimuth angle, and number of continuous tracking epochs; the ambiguity search and fixation are used to determine an ambiguity combination; the ambiguity combination includes integer solutions of the floating-point solutions of each navigation satellite obtained by fixing;

[0030] Based on probability theory and through statistical testing, a target integer ambiguity combination is searched out from the ambiguity combinations.

[0031] According to another aspect of the present disclosure, a UAV PPK post-processing device for fusion filtering and ambiguity processing is provided, comprising:

[0032] A data receiving and converting module, used for acquiring first positioning data received by the base station, and converting the first positioning data and the second positioning data into data in a predetermined data format; wherein the second positioning data is positioning data received by the drone;

[0033] A filtering processing module, used to perform N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format; wherein N is a positive integer greater than or equal to 4;

[0034] The fixed positioning solution module is used to perform cycle slip detection processing, floating point solution processing, integer ambiguity fixation processing, baseline solution, coordinate conversion and adjustment processing on the first positioning data and the second positioning data after each filtering, so as to obtain the UAV positioning result corresponding to each filtering;

[0035] The filtering fusion module is used for the UAV positioning results corresponding to each filtering. If the types of N UAV positioning results are the same, the types of the N UAV positioning results are used as the initial positioning result type; if the types of two UAV positioning results are the same, the types of the two UAV positioning results are used as the initial positioning result type; if the types of N UAV positioning results are different, the UAV positioning results are sorted according to the filtering order, and the types of the two middle UAV positioning results are used as the initial positioning result type; the values ​​of the N UAV positioning results are fused to obtain the target value of the UAV positioning result; the types of the N UAV positioning results are upgraded and upgraded to obtain the target positioning result type of the UAV positioning result;

[0036] Among them, the types of N drone positioning results are upgraded and downgraded, including:

[0037] Subtract the values ​​of the positioning results of two drones in the types of N drone positioning results to obtain multiple difference values. If all the obtained difference values ​​exceed the initial threshold, the type of the initial positioning result is downgraded;

[0038] Alternatively, the values ​​of the positioning results of the two drones are subtracted from each other, and if both differences exceed the initial threshold, the type of the initial positioning result is downgraded;

[0039] Among them, the numerical values ​​of N UAV positioning results are fused, including:

[0040] According to the values ​​of N UAV positioning results, and using the acceptance residual sigma as the inverse proportional weighting, the target value of the UAV positioning result is obtained;

[0041] Alternatively, the target value of the UAV positioning result is obtained by taking the positioning results of two UAVs of the same type and using the posterior residual sigma as an inverse proportional weight.

[0042] According to another aspect of the present disclosure, an electronic device is provided, including a memory, a processor, and a computer program stored in the memory, wherein the processor implements any of the above methods when executing the computer program.

[0043] According to another aspect of the present disclosure, a computer-readable storage medium is provided, wherein a computer program is stored in the computer-readable storage medium, and when the computer program is executed by a processor, any of the above-mentioned methods is implemented.

[0044] The disclosed UAV PPK post-processing method and device for fusion filtering and ambiguity processing performs multiple filtering on the positioning data, and fuses the multiple filtering results in the subsequent processing steps. In the fusion process, the disclosed method not only fuses the result type, but also fuses the value. Compared with the prior art that only fuses the result type, the accuracy and anti-interference of the fusion result are effectively improved. At the same time, the disclosed method performs composite fixation in the ambiguity fixing step, and the fixing rate is significantly better than the existing solution.

[0045] It should be understood that the content described in this section is not intended to identify the key or important features of the embodiments of the present disclosure, nor is it intended to limit the scope of the present disclosure. Other features of the present disclosure will become easily understood through the following description. BRIEF DESCRIPTION OF THE DRAWINGS

[0046] The accompanying drawings are used to better understand the present solution and do not constitute a limitation of the present disclosure.

[0047] Figure 1 is a flow chart of a UAV PPK post-processing method according to the present disclosure integrating filtering and ambiguity processing;

[0048] Figure 2A is a flowchart of cycle slip detection according to an embodiment of the present disclosure;

[0049] Figure 2B is a single-station data preprocessing flow chart according to an embodiment of the present disclosure;

[0050] Figure 2C is a workflow diagram of a polynomial fitting algorithm according to an embodiment of the present disclosure;

[0051] Figure 3 is a schematic diagram of the structure of a UAV PPK post-processing device according to the present disclosure that incorporates filtering and ambiguity processing;

[0052] Figure 4 is a schematic structural diagram of an electronic device according to the present disclosure. DETAILED DESCRIPTION

[0053] The following is a description of exemplary embodiments of the present disclosure in conjunction with the accompanying drawings, including various details of the embodiments of the present disclosure to facilitate understanding, which should be considered as merely exemplary. Therefore, it should be recognized by those of ordinary skill in the art that various changes and modifications may be made to the embodiments described herein without departing from the scope and spirit of the present disclosure. Similarly, for the sake of clarity and conciseness, descriptions of well-known functions and structures are omitted in the following description.

[0054] The present disclosure is due to the fact that the current UAV PPK post-processing solution only has the result fusion of two-way filtering, which cannot really bring into play the advantages of post-processing, and the current solution ambiguity fixation is only a simple partial ambiguity fixation, and there is no composite fixation method, the fixation rate is low, and the current solution result fusion technology is limited to the comparative fusion of result types, and does not examine the defects of more detailed parameters. The present disclosure provides a UAV PPK post-processing method and device with fusion filtering and ambiguity processing. The present disclosure performs multiple filtering on the unmanned positioning data, and fuses the multiple filtering results in the subsequent processing steps. The present disclosure not only fuses the result type in the fusion process, but also fuses the numerical value. Compared with the prior art that only fuses the result type, the accuracy and anti-interference of the fusion result are effectively improved. At the same time, the present disclosure performs composite fixation in the ambiguity fixing step, and the fixation rate is significantly better than the existing solution.

[0055] The technical solution of the present disclosure is described below through specific embodiments.

[0056] like Figure 1 As shown, it is a flowchart of the UAV PPK post-processing method of fusion filtering and fuzzy processing of this embodiment. The execution subject of this embodiment is a computing device or component with data processing capability. The specific method of this embodiment may include the following steps:

[0057] S110. Acquire first positioning data received by the base station, and convert the first positioning data and second positioning data into data in a predetermined data format; wherein the second positioning data is positioning data received by the drone.

[0058] The UAV can be regarded as a mobile station, receiving the UAV's second positioning data sent by the navigation satellite.

[0059] The base station here is a site whose position does not change, and it is used to receive the first positioning data of the drone sent by the navigation satellite.

[0060] The predetermined format here may be the RINEX format.

[0061] S120, performing N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format; wherein N is a positive integer greater than or equal to 4.

[0062] Specifically, the first positioning data and the second positioning data converted according to the predetermined data format can be filtered 4 times. The first filtering is a forward Kalman Kalman filter, and the first epoch is the initial state, and the initial state matrix and the initial variance matrix are the parameters of the initial state. The second filtering is a Kalman Kalman filter, and the last epoch is the initial state, and the initial state matrix and the initial variance matrix are the parameters of the initial state. The third filtering is a Kalman Kalman filter, wherein the first epoch of this filtering is the last epoch of the positioning data, and the initial state matrix and the initial variance matrix are the results of the last epoch of the first filtering, and the third filtering is performed halfway, and the third filtering ends at the middle time point epoch. The fourth filtering is a Kalman Kalman filter, wherein the first epoch of this filtering is the first epoch of the positioning data, and the initial state matrix and the initial variance matrix are the results of the last epoch of the second filtering, and the fourth filtering is only performed halfway, and the fourth filtering ends at the middle time point epoch.

[0063] The parameters will have different effects on the overall filtering accuracy and convergence time, so special processing is required for the number of filtering times and the setting of the initial state matrix and the initial variance matrix. The setting method disclosed in the present invention can effectively improve the filtering effect and is conducive to improving the positioning accuracy of subsequent drones.

[0064] S130, respectively perform cycle slip detection processing, floating point solution processing, integer ambiguity fixation processing, baseline solution, coordinate conversion and adjustment processing on the first positioning data and the second positioning data after each filtering to obtain the UAV positioning result corresponding to each filtering.

[0065] The cycle slip detection process, floating point solution process, baseline solution, coordinate conversion and adjustment process are all the same as the processing steps in the prior art.

[0066] The difference is that the above integer ambiguity fixing process may specifically include the following steps:

[0067] Step 1: Screen the satellite frequencies of all full-frequency navigation satellites, select the navigation satellite of the first frequency, and perform ambiguity search and fixation on the corresponding floating-point solution of the ambiguity obtained by the floating-point solution processing.

[0068] Step 2: Arrange the parameters of all full-frequency points, assign scores to the parameters according to weights, arrange all scores, perform star-kicking operations according to multiple preset ratios, and perform ambiguity search and fixation on the floating-point solutions of the ambiguities obtained by the floating-point solution processing of the remaining navigation satellites.

[0069] Step three, for the remaining navigation satellites, single satellites are eliminated, and single satellites are eliminated cyclically until the number of eliminated navigation satellites accounts for a first preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of the ambiguities obtained by the floating-point solution processing corresponding to the remaining navigation satellites; for the navigation satellites that cannot be fixed among the remaining navigation satellites, two satellites are eliminated, and so on, until the number of eliminated navigation satellites accounts for a second preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of the ambiguities obtained by the floating-point solution processing corresponding to the remaining navigation satellites; wherein the parameters include navigation satellite altitude angle, signal-to-noise ratio, azimuth, and number of continuous tracking epochs; the ambiguity search and fixation are used to determine an ambiguity combination; the ambiguity combination includes integer solutions of the fixed floating-point solutions of each navigation satellite.

[0070] The second preset percentage can be set according to actual conditions, for example, 25%. The first preset percentage is higher than the second preset percentage.

[0071] Step 4: Based on probability theory, through statistical testing, search for the target integer ambiguity combination from the ambiguity combinations.

[0072] S140, fusing the drone positioning results corresponding to each filtering: for the drone positioning results corresponding to each filtering, if the types of the N drone positioning results are the same, then the types of the N drone positioning results are used as the initial positioning result type; if the types of two drone positioning results are the same, then the types of the two drone positioning results are used as the initial positioning result type; if the types of the N drone positioning results are different, then the drone positioning results are sorted according to the filtering order, and the type of the intermediate drone positioning result is used as the initial positioning result type; the values ​​of the N drone positioning results are fused to obtain the target value of the drone positioning result; the types of the N drone positioning results are upgraded or downgraded to obtain the target positioning result type of the drone positioning result;

[0073] Among them, the types of N drone positioning results are upgraded and downgraded, including:

[0074] Subtract the values ​​of the positioning results of two drones in the types of N drone positioning results to obtain multiple difference values. If all the obtained difference values ​​exceed the initial threshold, the type of the initial positioning result is downgraded;

[0075] Alternatively, the values ​​of the positioning results of the two UAVs are subtracted from each other, and if both differences exceed the initial threshold, the type of the initial positioning result is downgraded.

[0076] The result types from high to low are fixed solution, floating point solution, and single point solution.

[0077] Among them, the numerical values ​​of N UAV positioning results are fused, including:

[0078] According to the values ​​of N UAV positioning results, and using the acceptance residual sigma as the inverse proportional weighting, the target value of the UAV positioning result is obtained;

[0079] Alternatively, the target value of the UAV positioning result is obtained by taking the positioning results of two UAVs of the same type and using the posterior residual sigma as an inverse proportional weight.

[0080] After the receiver captures the navigation satellite signal, as long as the tracking is not interrupted (lost lock), the receiver will automatically give the change in the number of integer cycles of the carrier phase during the tracking period. However, in the actual process, due to the temporary blocking of the navigation satellite signal or the influence of external interference factors, the satellite signal tracking is often temporarily interrupted, resulting in a cycle slip. When a cycle slip occurs, the carrier phase ranging accuracy will be seriously reduced, causing the RTK algorithm to lose the centimeter-level measurement accuracy. This step identifies and repairs the cycle slip, improves the accuracy of the data, and is conducive to improving the positioning accuracy of subsequent drones.

[0081] Among them, the ambiguity search and fixation of the navigation satellite can be achieved by the following steps:

[0082] The WL combination is used to search and fix the ambiguity of the navigation satellite; wherein the WL combination is a wide lane combination with a wavelength of 0.86m and an integer ambiguity.

[0083] The above steps can be implemented using specific modules. Specifically, the conversion of the data format in step S110 can be implemented using a data conversion module, which is mainly used to convert the raw data of the GNSS board into data in RINEX format so that the subsequent algorithm can read and process the data. The filtering in step S120 can be completed using a filtering processing module. The cycle slip detection in step S130 can be completed using a data preprocessing module. In step S130, the common-view satellite needs to be determined in the floating-point solution processing, and this step can be completed using a common-view satellite processing module. In step S130, the pseudo-range observation equation needs to be determined in the floating-point solution processing, and this step can be completed using a double-difference pseudo-range processing module. In step S130, the pseudo-range single-difference observation equation needs to be determined in the floating-point solution processing, and this step can be completed using a double-difference carrier phase processing module. In step S130, the pseudo-range double-difference equation needs to be determined in the floating-point solution processing, and the floating-point solution is performed in combination with the pseudo-range double-difference equation and the double-difference equation. This step can be completed using a floating-point solution processing module. The integer ambiguity fixation process in step S130 can be completed by the ambiguity search fixation module. The ambiguity fixation solution verification module can be used to verify the fixed solution in step S130. Step S140 can be completed by the filtering result fusion module.

[0084] The data preprocessing module can be used to complete the following processing:

[0085] To address the cycle slip problem, the present invention adopts a method that combines single-station non-difference data preprocessing (TurboEdit algorithm), differential data processing (polynomial fitting algorithm) and post-test residual analysis to perform preliminary cycle slip detection in the non-difference and differential observation data processing stages, and mark the detected cycle slips; after the differential carrier phase is solved, a post-test residual analysis is used to determine whether there are undetected cycle slips. For satellites where the ambiguity has been fixed and the residual exceeds the limit, the satellites that may have cycle slips are eliminated.

[0086] According to the signal frequency involved in the solution, the carrier phase differential algorithm can be divided into two categories, one is single-frequency differential positioning, and the other is dual-frequency differential positioning. The main difference lies in the different methods of carrier cycle slip detection and integer ambiguity search strategies. The dual-frequency differential positioning algorithm has many types of original observations, including dual-frequency carrier phase and dual-frequency pseudorange values, which can form a variety of dual-frequency phase and pseudorange combinations, which is conducive to cycle slip detection and repair, and speeds up the ambiguity search speed.

[0087] Several commonly used dual-frequency carrier linear combinations are given below, and their characteristics are explained using the GPS system as an example. Assuming that the phase observation values ​​of the L1 and L2 carriers of the GPS system are φ1 and φ2 respectively, the general form of the linear combination is φ n,m =n·φ1+m·φ2. m and n are the coefficients of the linear combination. The corresponding frequency of the linear combination observation is fn,m , wavelength λ n,m 、Integer ambiguity N n,m and measurement noise The relationship between the L1 and L2 response values ​​is:

[0088]

[0089] f1 and f2 represent the frequencies of the two carriers, N1 and N2 represent the integer ambiguities of the two carriers, c represents the speed of light, σ φ1 , σ φ2 The measured noise of the two carriers respectively.

[0090] The above method can be used to perform a variety of linear combinations of dual-frequency observations. Different combination observations have different application characteristics. The linear combinations of commonly used observations are shown in the following table:

[0091] Linear combination of GPS dual-frequency carrier phase observations

[0092]

[0093] The LI combination is called the ionosphere-free combination, which eliminates the first-order ionosphere effect and can significantly improve the accuracy of the medium and long baseline solution. For medium and long baselines (baseline distance greater than 10km), the residual ionospheric delay error after single difference is still large, which is the main error source affecting positioning accuracy. This combination is used for baseline calculation in this disclosure. For short baselines (baseline distance within 10km), the ionospheric delay error can be basically eliminated after single difference, and observation noise is the main error. The observation noise of the LI combination is large, so single-frequency observations are used for baseline calculation in this disclosure.

[0094] The LG combination is called the ionospheric residual combination. This combination has nothing to do with the geometric distance from the receiver to the navigation satellite, and eliminates the geometric distance from the receiver to the navigation satellite, orbit error, receiver clock error, satellite clock error and tropospheric delay error. It is only related to the ionospheric delay, combined ambiguity and observation value noise, so the observation value should be very smooth in theory. In the present disclosure, the LG combination is used to detect and repair dual-frequency cycle slips.

[0095] The WL combination is called a wide lane combination, with a wavelength of 0.86 m, an integer ambiguity, and a small ionospheric delay. In the present disclosure, the LW combination is used to accelerate the fixation of the ambiguity.

[0096] The LN combination is called LW, which is a narrow lane combination with a wavelength of 0.11m. The relative observation noise is close to that of L1, and the ambiguity is also an integer. The ionospheric effect is equal to that of LN, but the sign is opposite. Although the relative noise of LN is close to that of L1, the ratio of the ionospheric error of LN to the wavelength is too large, so it is only used for short baseline high-precision positioning. This combination is not used in this disclosure.

[0097] In addition, through the combination of dual-frequency carrier and dual-frequency pseudorange, various combinations can be formed, such as MW combination, GF combination, etc., which can be used for the detection of carrier cycle slips.

[0098] In the observation data processing stage, a polynomial fitting algorithm is used for single-frequency data to perform cycle slip detection on the double-difference single-frequency carrier phase data. For dual-frequency data, since there are relatively more observation quantities and a variety of dual-frequency combinations can be used to improve the detection efficiency, the TurboEdit algorithm is used in the design, combined with dual-frequency pseudorange data, to perform cycle slip detection on the non-difference carrier phase data, and then the polynomial fitting method is used to perform cycle slip detection on the single-difference L4 and L5 combinations.

[0099] In the carrier phase positioning algorithm, after carrier cycle slip detection, the cycle slip repair method can be further used for processing. However, there are often risks in carrier cycle slip repair, so another solution is adopted in the present disclosure: once a satellite detects a cycle slip, the ambiguity of the satellite is set to a new unknown. In this case, there is no need to repair the size of the cycle slip value to avoid the impact of the cycle slip value determination error on subsequent calculations. At the same time, this method also has some disadvantages. When the data quality decreases and there are too many satellites with cycle slips, the number of new ambiguity unknowns will also increase, which will affect the positioning results. However, considering that the errors caused by cycle slip repair errors are often even more unacceptable, the strategy of cycle slip detection and repair (marking) of single-station data is adopted here, and in differential data processing, only detection but no repair is performed. This not only enables the algorithm to have a certain cycle slip repair capability, but also reduces the risk of cycle slip repair errors. Figure 2A The figure shows the workflow of cycle slip detection before solution.

[0100] TurboEdit Algorithm

[0101] At present, TurboEdit is generally used in dual-frequency RTK positioning technology to detect and repair cycle slips. Generally, the measurement type of dual-frequency GPS receivers contains the following observations: CA code pseudorange C1; L1 precise code pseudorange P1; L2 precise code pseudorange P2; L1 carrier phase L2 carrier phase Doppler observation values ​​D1, D2. The following will introduce the cycle slip detection workflow under different combinations based on the above observations.

[0102] (1) Melbourne-Wübbena combination

[0103] The properties and characteristics of the Melbourne-Wübbena combined observations: ① It eliminates the influence of errors such as the ionosphere, troposphere, receiver and satellite clock errors, and station geometry, and is only affected by observation noise and multipath effects; ② It has a longer wavelength (about 86 cm) and smaller measurement noise; ③ The calculation results of the combined observations only contain wide-lane ambiguity parameters. Therefore, the MW combination can well complete the detection and repair of wide-lane cycle slips, eliminate gross errors, etc. The MW combination value and its variance are:

[0104]

[0105] Wherein, L7(i) is the MW combined observation value at epoch i, L1(i) and L2(i) are the dual-frequency carrier phase observation values, f1 and f2 are the dual-frequency carrier frequencies, λ7 is the MW combined carrier wavelength, and N7(i) is the MW combined carrier integer ambiguity; They are the measurement variance of MW combined observations, the measurement variance of L1 observations, the measurement variance of L2 observations, the measurement variance of P1 observations, and the measurement variance of P2 observations at epoch i respectively.

[0106] Therefore, the integer ambiguity and variance of the combined observation value are expressed as:

[0107] N7(i)=L7(i) / λ7 (4)

[0108]

[0109] In actual calculation, the recursive method can be used to calculate the ambiguity prediction value of each epoch integer <N7> i and its variance value:

[0110]

[0111] For the observation data of the i-th epoch, if |N7(i)-<N7> i-1 |<4σ i-1 , it is assumed that no cycle slip occurs in epoch i, and the mean and variance of epoch i+1 are then recursively deduced.

[0112] If |N7(i)-<N7> i-1 |≥4σ i-1 , then the i epoch may have a cycle slip or be an outlier point. If a cycle slip occurs, by calculating the difference ΔN7 between the integer ambiguities of the multiple epochs before and after the cycle slip, ΔN7 is the cycle slip value between epochs, and the cycle slip ΔN7 has the following relationship with the cycle slips ΔN1 and ΔN2 of the L1 and L2 observation values:

[0113] ΔN7=ΔN1-ΔN2(8)

[0114] It should be noted that if the cycle slips occurring on L1 and L2 are equal in size, that is, ΔN1=ΔN2, then ΔN7=0, the cycle slips in the observed values ​​cannot be detected using the MW combination. Therefore, it is necessary to combine other dual-frequency combination methods for cycle slip detection.

[0115] (2) GF combination without geometric structure

[0116] The ionospheric residual combination L4 of the dual-frequency phase is combined with the ionospheric residual combination P4 of the dual-frequency code measurement to form the geometric structure-free combined observation value LP4 of the dual-frequency phase and pseudorange, and its formula is as follows:

[0117] LP4=L4+P4=(L1-L2)+(P1-P2)=λ2N2-λ1N1 (9)

[0118] This combination eliminates the influence of errors such as ionosphere, troposphere, receiver and satellite clock errors, and station satellite geometry, and the result is the difference between the ambiguities of L1 and L2. Since the whole cycle ambiguity remains unchanged when there is no cycle slip, this combined observation value is suitable for data processing such as gross error elimination, cycle slip detection and repair, and when the cycle slips on L1 and L2 are equal, that is, ΔN1 = ΔN2, this combination can also detect cycle slips.

[0119] By eliminating the gross error observations and calculating the cycle slip value ΔN4 of L4 at the cycle slip based on the multi-epoch observation data before and after the cycle slip, the relationship between ΔN4 and the cycle slips ΔN1 and ΔN2 of the L1 and L2 observations is as follows:

[0120] ΔN4=λ2ΔN2-λ1ΔN1(10)

[0121] According to ΔN7 and ΔN4, ΔN1 and ΔN2 can be calculated as follows:

[0122]

[0123] (3) Ionospheric-free IF combination

[0124] The ionosphere elimination combination L3 of the dual-frequency phase and the ionosphere elimination combination P3 of the dual-frequency code pseudorange are combined again to form the observation value LP3:

[0125]

[0126] This combined observation eliminates the effects of geometric distance, troposphere, ionosphere, etc., and only includes the effects of noise. Its disadvantage is that the noise is significantly amplified (three times the P1 code noise). However, LP3 can be used to detect gross errors caused by the receiver's own system errors and eliminate poor quality observations that are not deleted in the MW combination.

[0127] (4) Based on the above three traditional carrier and pseudorange combination algorithms, the present disclosure adopts a simplified TurboEdit algorithm, that is, the MW combination and the simplified ionospheric residual combination LG combination are used to perform cycle slip detection respectively. Section 3.2.1 gives the LG combination method, namely:

[0128] L4=L1-L2=λ2N2-λ1N1+ΔION L4 (13)

[0129] where ΔION L4 is the dual-frequency ionospheric residual, whose value is the L1 and L2 dual-frequency ionospheric error ΔION L1 , ΔION L2 The difference, that is:

[0130]

[0131] The observations of P1 and P2 are not used in the LG combination, and the results are no longer dependent on the measurement accuracy of the P code, which improves the measurement accuracy from decimeter level to millimeter level. However, the ionospheric residual value is introduced at the same time. However, the dual-frequency ionospheric residual value changes very little in the previous and next epochs. The following is an estimate of the change rate of the ionospheric residual value:

[0132] The rate of change of the ionospheric residual error between the previous and next epochs Δε I The single-frequency L1 ionospheric error change rate Δε IL1 The relationship is:

[0133]

[0134] The MEO satellite period is 11h56min. During this period, the satellite rotates 360° around the earth. The change in the observed altitude angle of the GPS satellite within 1 second is about 0.5′ (the corresponding radian value is ). Taking GPS L1 as an example, when the satellite moves Δφ, the corresponding ionospheric delay change rate Δε IL1 Approximately:

[0135]

[0136] Since Δφ is a very small quantity, we can get:

[0137]

[0138] Single frequency L1 ionospheric error ΔION L1 The error value is affected by factors such as satellite elevation angle and electron concentration, and the error value is about 5m-100m; α is the satellite elevation angle, and the elevation angle is greater than 15°. Therefore, for the extreme case, the ionospheric error is 100m and the elevation angle is 15°, we can get:

[0139]

[0140] Δε I ≈0.7Δε IL1 =0.0384m (19)

[0141] It can be seen that the change rate of the ionospheric residual of the LG combination observation value within 1 second does not exceed 4cm, so the error of the LG combination observation value is also at the centimeter level, which is easier to detect cycle slips than the GF combination. In this disclosure, the threshold for cycle slip detection of the LG combination is set to 5cm. The LG combination value can detect most cycle slips, except for some L1 and L2 that meet the special relationship

[0142] λ1ΔN1-λ2ΔN2|<5cm cycle slip.

[0143] When encountering a cycle slip that cannot be detected by the GF combination, the MW combination can basically achieve effective detection of single-station non-difference data.

[0144] (5) The above-mentioned combined observation values ​​can largely detect the cycle slips on the L1 and L2 carriers and their sizes ΔN1 and ΔN2. After completing the cycle slip detection, the processing of the carrier observation values ​​at the cycle slip location mainly includes two methods: the first is to correct the carrier phase observation values ​​after the cycle slip location based on the cycle slip ΔN1 and ΔN2; the second is to reset the ambiguity parameters at each cycle slip location (many GNSS post-processing software currently use this method, such as Bernese 5.0 software).

[0145] In this project, when the receiver observation values ​​have dual-frequency P1, P2, L1, and L2, the single-station observation data is preprocessed based on the simplified TurboEdit algorithm. When a cycle slip occurs in the detection, the cycle slip values ​​on L1 and L2 can be solved according to the corresponding formula, the carrier phase observation values ​​can be corrected, and the corresponding satellite can be marked.

[0146] After simplifying the TurboEdit data preprocessing, the observation data entering the main RTK algorithm is relatively "clean" and carries the corresponding "cycle slip repair mark". At this time, the main RTK algorithm can process these "marked" satellites according to the corresponding specific needs, using single-difference-based The polynomial fitting method of the observations detects, repairs and marks the cycle slips, thereby reconfirming the cycle slip repair of the TurboEdit algorithm. For single-frequency receivers, the single-station preprocessing step will be skipped because there are no corresponding P2, L2, etc. observations.

[0147] Taking GPS processing as an example, the single-station data preprocessing process based on the TurboEdit algorithm is as follows: Figure 2B shown.

[0148] Polynomial Fitting Algorithm

[0149] Based on the observation data of single-station dual-frequency (including P1 and P2), simplified TurboEdit data preprocessing can be performed, and relatively "clean" carrier observation data can be obtained. However, due to the limitations of P code observation noise and simplified TurboEdit algorithm, the effect of preprocessing cannot be fully guaranteed. To this end, after performing inter-station data differentiation, a cycle slip detection and marking module based on a polynomial fitting algorithm was designed and implemented, and ambiguity parameters were added at the cycle slip.

[0150] The principles and methods of polynomial fitting are as follows:

[0151] Substitute the m cycle slip-free difference observations into the following formula:

[0152]

[0153] Where n is the polynomial fitting order; t0 is the time base; t i is the time variation; t i The fitted carrier phase difference observation value corresponding to the time. Where i = 1, 2, 3, ..., m (m> n + 1).

[0154] Step 1: Use the least squares adjustment principle

[0155]

[0156] Obtain the polynomial coefficients a0, a1, a2, a3 in the formula, and calculate the residual v after fitting i Calculation error

[0157]

[0158] Step 2: Use the obtained polynomial coefficients to extrapolate the double-difference carrier phase value of the next epoch and compare it with the actual observation value. When the difference between the two is less than 3δ, it is considered that the double difference observation has no cycle slip, the earliest observation is removed, the actual double difference observation is added, and then the first step is returned to continue the polynomial fitting to extrapolate the next epoch.

[0159] When the difference between the extrapolated value and the actual double-difference observation value is greater than 3δ, the actual double-difference observation value is considered to contain a cycle slip. At this time, the whole cycle of the extrapolated value should be used to replace the whole cycle number of the observation value with a cycle slip error. Then remove the earliest observation value and add the corrected Go back to the first step to continue to calculate the polynomial coefficients and extrapolate the differential value of the next epoch. Usually, the change of the fitting value above the third order is very small, and the contribution to the fitting accuracy is limited. Therefore, the current algorithm uses a second-order polynomial fitting based on 5 epoch observations for cycle slip detection. Figure 2C The figure shows the workflow of the polynomial fitting algorithm.

[0160] When the observation data is single-frequency, the single difference between stations cannot eliminate the receiver clock error, so the double difference method is used to further eliminate this error and improve the accuracy of cycle slip detection; when the observation data is dual-frequency, the selected dual-frequency combination Single difference can eliminate all kinds of major errors including receiver clock error, so polynomial fitting cycle slip detection can be performed in the single difference stage. In this stage, the strategy of only detection without repair is adopted, and the differential observation data with possible cycle slips are "marked" and new ambiguity parameters are added. Subsequently, the double difference ambiguity corresponding to the "marked" satellite will be reset to the parameters to be estimated, and the ambiguity parameters will be re-solved.

[0161] The double-difference pseudorange processing module is used to complete the following processing:

[0162] Pseudorange difference uses pseudorange observations as the basic input. It is assumed that the base station and the mobile station synchronously observe a group of navigation satellites to obtain the pseudorange observations of n common-view satellites synchronously observed by the base station (station A) and the mobile station (station B). and i=1,2,3......n, select the navigation satellite r with the highest elevation angle among the n navigation satellites as the reference satellite to form a relationship between any navigation satellite j (j=1,2,3......n and j≠r) and the reference satellite r at any t i The pseudo-range observation equation at time is as follows:

[0163]

[0164] The meanings of the parameters in the above formula are as follows:

[0165] c: speed of light (m / s);

[0166] λ: Satellite navigation signal carrier wavelength (m);

[0167] f: Satellite navigation signal carrier frequency (Hz);

[0168] t i Pseudorange observation value of navigation satellite r observed by the reference station at the moment (m);

[0169] t i Pseudorange observation value of navigation satellite j observed by the reference station at time (m);

[0170] t i The geometric distance between the time reference station and the navigation satellite r (m);

[0171] t i The geometric distance between the time reference station and navigation satellite j (m);

[0172] δt A (t i ): t i The receiver clock error of the reference station at the moment (s);

[0173] δt r (t i ): t i Satellite clock error of navigation satellite r at time (s);

[0174] δt j (t i ): t i Satellite clock error of navigation satellite j at time (s);

[0175] t i Earth rotation error between the time reference station and the navigation satellite r (m);

[0176] t i Earth rotation error between the time reference station and navigation satellite j (m);

[0177] t i Ionospheric delay error between the time reference station and navigation satellite r (m);

[0178] t i Ionospheric delay error between the time reference station and navigation satellite j (m);

[0179] t i Tropospheric delay error between the time reference station and navigation satellite r (m);

[0180] t i Tropospheric delay error between the time reference station and navigation satellite j (m);

[0181] t i Thermal noise of pseudo-range measurement of the reference station receiver about navigation satellite r at the moment (m);

[0182] t iThermal noise of pseudo-range measurement of navigation satellite j by the base station receiver at time (m);

[0183] t i Pseudorange observation value of navigation satellite r observed by the mobile station at the moment (m);

[0184] t i Pseudorange observation value of navigation satellite j observed by the mobile station at time (m);

[0185] t i The geometric distance between the mobile station and the navigation satellite r at any moment (m);

[0186] t i The geometric distance between the mobile station and navigation satellite j at the moment (m);

[0187] δt B (t i ): t i The mobile station receiver clock error at the moment (s);

[0188] t i Earth rotation error between the mobile station and the navigation satellite r at the moment (m);

[0189] t i Earth rotation error between the mobile station and navigation satellite j at the moment (m);

[0190] t i Ionospheric delay error between the mobile station and navigation satellite r at the moment (m);

[0191] t i Ionospheric delay error between the mobile station and navigation satellite j at time (m);

[0192] t i Tropospheric delay error between the mobile station and navigation satellite r at time (m);

[0193] t i Tropospheric delay error between the mobile station and navigation satellite j at time (m);

[0194] t i Thermal noise of the pseudorange measurement of the mobile station receiver about the navigation satellite r at time (m);

[0195] ti Thermal noise of pseudorange measurement of rover receiver about navigation satellite j at time (m);

[0196] Linearization:

[0197] In the above formula,

[0198]

[0199] Where (X r ,Y r ,Z r ) is the satellite position, (X A ,Y A ,Z A ) is the receiver position

[0200] Let (X A0 ,Y A0 ,Z A0 ), (δX A (t),δY A (t),δZ A (t)) are the approximate values ​​and corrections of the observation station coordinates, respectively. A0 ,Y A0 ,Z A0 ) can be linearized by using Taylor expansion as follows:

[0201]

[0202] In the formula,

[0203]

[0204] Then the formula in (23) can be transformed into:

[0205]

[0206] The double-difference carrier phase processing module is used to complete the following processing:

[0207] Due to the existence of multiple errors such as user receiver clock error, navigation satellite clock error, earth rotation error, ionospheric delay and tropospheric delay, the pseudorange observation equation and carrier phase observation equation have many unknown parameters, complex solutions, and poor solution accuracy. Since the common errors such as navigation satellite clock error, earth rotation error, ionospheric delay and tropospheric delay between the base station and the mobile station are highly correlated, the influence of the common errors can be eliminated by the difference between the observation equations to improve the solution accuracy. By subtracting the two pseudorange observation equations corresponding to the same navigation satellite of the base station and the mobile station, the corresponding pseudorange single difference observation equation can be obtained as follows:

[0208]

[0209] The meanings of the parameters in the above formula are as follows:

[0210] t i The single difference of pseudorange between the reference station and the mobile station with respect to the navigation satellite r at the moment (m);

[0211] t i The single difference of pseudorange between the reference station and the mobile station about navigation satellite j at the moment (m);

[0212] δt AB (t i ): t i The difference between the base station receiver clock error and the mobile station receiver clock error at the moment (s);

[0213] t i The single difference of pseudorange noise between the reference station and the mobile station about the navigation satellite r at the moment (m); t i The single difference of pseudorange noise between the reference station and the mobile station about navigation satellite j at the moment (m);

[0214] The single-difference observation equation eliminates the influence of common errors such as navigation satellite clock error, ionospheric delay error, and tropospheric delay error. However, due to the existence of the user receiver clock error in the single-difference observation equation, the complexity of the solution is increased, affecting the solution accuracy. By re-difference of the single-difference observation equation, the influence of the user receiver clock error can be eliminated, further improving the accuracy of the relative position solution.

[0215] Subtract the pseudo-range single-difference observation equation of navigation satellite j observed by the base station and the mobile station from the pseudo-range single-difference observation equation of the reference satellite r, and obtain n-1 pseudo-range double-difference equations as follows:

[0216]

[0217] The meanings of the parameters in the above formula are as follows:

[0218] t i Pseudorange measurement residual between the reference station and the mobile station at the moment (m);

[0219] t i The double difference measurement coefficient matrix of the mobile station about navigation satellites r and j at the moment;

[0220] t i Position correction vector between mobile stations at the moment (m);

[0221] The constant term (m) of the pseudorange measurements of the base station and the rover station with respect to the navigation satellites r and j;

[0222] The double-difference pseudorange observation equation can be simplified as:

[0223] L ρ =AX+V ρ (30)

[0224] The meanings of the parameters in the above formula are as follows:

[0225] X: position correction vector of the rover;

[0226] L ρ : Pseudorange double difference measurement constant vector between the base station and the rover;

[0227] A: relative position vector coefficient matrix between the base station and the mobile station;

[0228] V ρ : Pseudorange double difference measurement residual between the base station and the rover;

[0229] The floating point solution processing module is used to complete the following processing:

[0230] For n GPS solution satellites, the carrier phase double difference observation equation contains 3 unknown relative position parameters and n-1 unknown integer ambiguity parameters. The number of unknown parameters is greater than the number of equations. In order to solve the above carrier phase double difference observation equation, it is necessary to comprehensively utilize the two observation quantities of pseudorange and carrier. The pseudorange carrier phase joint solution process is given below.

[0231] Combining the double-difference pseudorange equation and the carrier phase double-difference equation, we can get the following set of equations:

[0232]

[0233] The meanings of the parameters in the above formula are as follows:

[0234] a: relative position correction vector between the base station and the rover;

[0235] b: integer ambiguity vector between the base station and the rover;

[0236] L ρ : Pseudorange double difference measurement constant vector between the base station and the rover;

[0237] The carrier phase double difference measurement constant vector between the base station and the rover;

[0238] A: relative position vector coefficient matrix between the base station and the mobile station;

[0239] B: integer ambiguity vector coefficient matrix between the base station and the mobile station;

[0240] V ρ : Pseudorange double difference measurement residual between the base station and the rover;

[0241] Carrier phase double-difference measurement residuals between the base station and the rover;

[0242] P ρ : Pseudorange double difference measurement weight matrix between the base station and the mobile station;

[0243] The carrier phase double difference measurement weight matrix between the base station and the rover;

[0244] P: basic weight coefficient matrix, n is the number of stars to be solved, lin and col represent the rows and columns of the matrix P, is the measurement noise variance;

[0245] q: The weight coefficient factor of the pseudorange double difference measurement between the base station and the rover. Its size depends on the pseudorange measurement accuracy and error level, and its typical value is 10-4 to 10-6;

[0246] The above equations can be expressed in the following matrix form:

[0247]

[0248] Adding the coefficient matrix, we can get the following using the least squares method:

[0249]

[0250] The above formula can be simplified as:

[0251]

[0252] Let N aa =A T PA, N ab =A T PB, N ba =B T PA, N bb =B T PB, The above formula can be simplified as:

[0253]

[0254] The relative position and integer ambiguity floating point solution can be obtained by using the least squares method:

[0255]

[0256]

[0257] The fuzzy search fixed module can complete the following data processing:

[0258] A three-step method is used to fix the ambiguity: 1. Ambiguity search is performed only on the first frequency point of all ambiguities; 2. Satellite elevation angle, signal-to-noise ratio, azimuth angle and other information are used to control and screen satellites to fix some ambiguities; 3. Finally, the partial satellite elimination method is used to perform cyclic satellite removal and then fix the ambiguity.

[0259] The first step is to screen the satellite frequencies of all satellites in the system and select only the satellites with the first frequency for ambiguity search. This can reduce the fixed search time for ambiguity while also taking into account the positioning accuracy of the fixed solution.

[0260] The second step is to screen the satellite frequencies based on information such as satellite altitude angle, signal-to-noise ratio, and number of continuous tracking epochs. First, the satellite altitude angle, signal-to-noise ratio, and number of continuous tracking epochs of the entire system and all frequencies are arranged, and these parameters are scored according to their weights. The overall scores are arranged, and satellites are kicked out as a whole according to the three ratios of 90%, 85%, and 80%, and the search for the remaining satellites is fixed.

[0261] The third step is to use the method of partially eliminating satellites to perform ambiguity search and fixation. The first step is to eliminate a single satellite, and then eliminate a single satellite in a loop, and then perform ambiguity search and fixation. If this method cannot fix it, then eliminate two satellites, and so on, until the number of eliminated satellites reaches 25% of all satellites, and no more satellite elimination is performed.

[0262] The ambiguity fixed solution verification module is used to complete the following processing:

[0263] Fixed solution verification method

[0264] Ambiguity confirmation is the last step of integer ambiguity resolution. Its main function is to select the best from the set of possible ambiguity combinations obtained during the ambiguity search process and determine its correctness. Generally speaking, all ambiguity resolution methods select the best from the set of ambiguity combinations through some kind of test, so the choice of test method becomes a crucial issue. If the conditions are too loose, the entire process will be prolonged, which is not conducive to the rapid resolution of ambiguity; if the conditions are too strict, the correct ambiguity is likely to be excluded.

[0265] The ambiguity confirmation algorithm used here is based on probability theory. Through three statistical tests, it is finally confirmed that the searched integer ambiguity combination is the correct integer ambiguity combination to be found. After substituting this group of integer ambiguities into the normal equation, the solution obtained is the correct fixed solution.

[0266] (1) Check the consistency of the baseline vector obtained by the integer solution and the initial solution.

[0267] Let the baseline vector obtained by the integer solution be X, and the baseline vector obtained by the initial solution be The corresponding cofactor array If the following equation holds, then From the perspective of statistical testing, it is consistent and compatible with X:

[0268]

[0269] where β = ξ F (u,f,1-α) is the one-tail quantile of the Fisher distribution with confidence level 1-α and degrees of freedom f and u; u is the number of unknown parameters; f is the degrees of freedom in parameter estimation; σ0 is the unit weighted mean error in the initial solution.

[0270] Since the ambiguity parameters are already within the corresponding confidence intervals, as long as the baseline vector is also consistent, it means that the integer solution and the initial solution vector are consistent.

[0271] (2)Consistency check of the unit weighted mean error between the integer solution and the initial solution.

[0272] Assume the unit weighted mean error of the integer solution is σ A , the unit weighted mean error of the initial solution is σ0, if the following formula holds, it means that the two statistical tests are consistent:

[0273]

[0274] The above test is also called the χ test of the variance factor. 2 The symbols in the above formula have the same meanings as above.

[0275] (3) Ratio test: the smallest unit weighted mean error σ and the second smallest unit weighted mean error σ in the integer solution 次 Significance test between .

[0276] Ratio test is the most commonly used test method. Its basic idea is to compare the minimum residual sum of squares with the second minimum residual sum of squares. Since the double difference residual calculated with the correct ambiguity group is significantly smaller than other incorrect ones, a threshold can be set according to various error factors such as measurement error and multipath error level. If it meets:

[0277]

[0278] That is, if it is greater than the preset threshold value, it is considered to have passed the test and the fixed solution is correct, otherwise it is considered that the solution still has deviations. In this design, the ratio value is set to ≥3.

[0279] If any of the above three tests fails, the search fails and re-testing or other measures (such as expanding the confidence interval) are required.

[0280] The filtering result fusion module is used to complete the following processing:

[0281] Since the present invention has four filtering results to be fused, the result positioning types are arranged first. Each epoch will have four positioning results. The result type can be preliminarily determined according to the following three schemes: first, the result type is defined first. If the four result types are the same, this positioning result type is temporarily determined; second, if two are the same, these two same types are temporarily used as the positioning result type; if the four result types are not the same, the result in the middle of the order shall prevail.

[0282] Fusion of result values: The first method is to use the acceptance residual sigma as inverse proportional weighting according to the existing four results to finally obtain the positioning result; the second method is to use the post-test residual sigma as inverse proportional weighting according to two positioning results of the same type to finally obtain the positioning result.

[0283] Upgrade or downgrade the result type: The first type is to perform mutual subtraction of the two results of the four positioning results to obtain 12 difference values. If all 12 difference values ​​exceed the initial threshold, the result is downgraded. The second type is to perform mutual subtraction of two positioning results. If both difference values ​​exceed the initial threshold, the result is downgraded. The third type is not to process the result type of the data.

[0284] At this point, the target positioning result type and target value of the drone are obtained, that is, the positioning result of the drone is obtained.

[0285] In the present disclosure, in the UAV PPK positioning, the mobile station, i.e. the UAV, obtains the relevant observation data of the base station through post-processing, and then combines the local original observation data to form a double difference equation between the two, eliminating most of the correlation errors, so that the carrier phase double difference ambiguity can be fixed quickly to obtain a high-precision positioning result. In order to ensure the positioning accuracy and reliability, according to the actual application requirements of RTK positioning of GNSS measurement receivers. The algorithm mainly includes the analysis and acquisition of dual-station original observation data (the data includes GPS: L1L1L5 / BDS: B1I B2I B3I B1C B2A / QZSS: L1 L2 L5 / GAL: E1 E5aE5b / GLN: LG1 G2), and the original data includes pseudorange, carrier phase, Doppler, signal-to-noise ratio information, etc., single-station data preprocessing based on simplified TurboEdit, base station and mobile station single-point positioning and single-point velocity measurement, pseudorange differential positioning, cycle slip detection, repair and marking based on polynomial fitting algorithm, establishment of carrier phase double difference observation and adjustment, solution of floating point ambiguity and covariance matrix, LAMBDA ambiguity search and fixation, solution based on fixed ambiguity baseline, post-verification residual estimation analysis and result output, etc. The present invention uses four-time multi-directional filtering results in the PPK solution process and fuses the results. Different from the existing two-way filtering method, the present invention uses four-time multi-directional result fusion in the PPK solution to comprehensively consider the baseline length, ratio, nb, rms and other parameters. At the same time, the present invention uses a composite fixing method in the existing traditional fixing method in the PPK solution, and the fixing rate is significantly better than the existing method.

[0286] Based on the same inventive concept, the present disclosure provides a UAV PPK post-processing device with fusion filtering and fuzzy processing. The steps performed by the components of the device are the same or similar to those of the above method, so the similar parts are not repeated. Figure 3 As shown, the UAV PPK post-processing device for fusion filtering and fuzzy processing in this embodiment includes:

[0287] The data receiving and converting module 310 is used to obtain the first positioning data received by the base station, and convert the first positioning data and the second positioning data into data in a predetermined data format; wherein the second positioning data is the positioning data received by the drone.

[0288] The filtering processing module 320 is used to perform N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format; wherein N is a positive integer greater than or equal to 4.

[0289] The fixed positioning solution module 330 is used to perform cycle slip detection processing, floating point solution processing, integer ambiguity fixation processing, baseline solution, coordinate conversion and adjustment processing on the first positioning data and the second positioning data after each filtering, so as to obtain the UAV positioning result corresponding to each filtering;

[0290] The filtering fusion module 340 is used for the drone positioning results corresponding to each filtering. If the types of the N drone positioning results are the same, the types of the N drone positioning results are used as the initial positioning result type; if the types of two drone positioning results are the same, the types of the two drone positioning results are used as the initial positioning result type; if the types of the N drone positioning results are different, the drone positioning results are sorted according to the filtering order, and the types of the two middle drone positioning results are used as the initial positioning result type; the values ​​of the N drone positioning results are fused to obtain the target value of the drone positioning result; the types of the N drone positioning results are upgraded or downgraded to obtain the target positioning result type of the drone positioning result;

[0291] Among them, the types of N drone positioning results are upgraded and downgraded, including:

[0292] Subtract the values ​​of the positioning results of two drones in the types of N drone positioning results to obtain multiple difference values. If all the obtained difference values ​​exceed the initial threshold, the type of the initial positioning result is downgraded;

[0293] Alternatively, the values ​​of the positioning results of the two drones are subtracted from each other, and if both differences exceed the initial threshold, the type of the initial positioning result is downgraded;

[0294] Among them, the numerical values ​​of N UAV positioning results are fused, including:

[0295] According to the values ​​of N UAV positioning results, and using the acceptance residual sigma as the inverse proportional weighting, the target value of the UAV positioning result is obtained;

[0296] Alternatively, the target value of the UAV positioning result is obtained by taking the positioning results of two UAVs of the same type and using the posterior residual sigma as an inverse proportional weight.

[0297] In some embodiments, when the filtering processing module 320 performs N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format, it is specifically used to:

[0298] The first positioning data and the second positioning data converted according to the predetermined data format are subjected to four filtering processes.

[0299] In some embodiments, the first filtering is a forward Kalman filter, and the first epoch is an initial state, and the initial state matrix and the initial variance matrix are parameters of the initial state; wherein the epoch is a frame of data.

[0300] In some embodiments, the second filtering is a Kalman filter, and the last epoch is the initial state, and the initial state matrix and the initial variance matrix are parameters of the initial state.

[0301] In some embodiments, the third filtering is a Kalman filter, wherein the first epoch of the third filtering is the last epoch, and the initial state matrix and the initial variance matrix are the results of the last epoch of the first filtering, and the third filtering is performed halfway, and ends at the middle time point epoch.

[0302] In some embodiments, the fourth filtering is a Kalman filter, wherein the first epoch of the fourth filtering is the first epoch, and the initial state matrix and the initial variance matrix are the results of the last epoch of the second filtering, and the fourth filtering is only performed halfway, and the fourth filtering ends at the middle time point epoch.

[0303] In some embodiments, the fixed positioning module 330 is specifically used to:

[0304] Screening satellite frequencies of all full-frequency navigation satellites, selecting a navigation satellite of the first frequency, and performing ambiguity search and fixing on the corresponding ambiguity floating-point solutions obtained by floating-point solution processing;

[0305] Arrange the parameters of all full-frequency points, assign scores to the parameters according to weights, arrange all scores, perform star-kicking operations according to multiple preset ratios, and perform ambiguity search and fixation on the floating-point solutions of the ambiguities obtained by the floating-point solution processing corresponding to the remaining navigation satellites;

[0306] For the remaining navigation satellites, single satellites are eliminated, and single satellites are eliminated cyclically until the number of eliminated navigation satellites accounts for a first preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of ambiguities obtained by floating-point solution processing corresponding to the remaining navigation satellites; for the navigation satellites that cannot be fixed among the remaining navigation satellites, two satellites are eliminated, and so on, until the number of eliminated navigation satellites accounts for a second preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of ambiguities obtained by floating-point solution processing corresponding to the remaining navigation satellites; wherein the parameters include navigation satellite altitude angle, signal-to-noise ratio, azimuth angle, and number of continuous tracking epochs; the ambiguity search and fixation are used to determine an ambiguity combination; the ambiguity combination includes integer solutions of the floating-point solutions of each navigation satellite obtained by fixing;

[0307] Based on probability theory and through statistical testing, a target integer ambiguity combination is searched out from the ambiguity combinations.

[0308] According to an embodiment of the present disclosure, the present disclosure also provides an electronic device and a computer-readable storage medium.

[0309] Figure 4 A schematic block diagram of an example electronic device 400 that can be used to implement an embodiment of the present disclosure is shown. The electronic device is intended to represent various forms of digital computers, such as laptop computers, desktop computers, workstations, personal digital assistants, servers, blade servers, mainframe computers, and other suitable computers. The electronic device can also represent various forms of mobile devices, such as personal digital processing, cellular phones, smart phones, wearable devices, and other similar computing devices. The components shown herein, their connections and relationships, and their functions are merely examples and are not intended to limit the implementation of the present disclosure described and / or required herein.

[0310] like Figure 4 As shown, the device 400 includes a computing unit 410, which can perform various appropriate actions and processes according to a computer program stored in a read-only memory (ROM) 420 or a computer program loaded from a storage unit 480 into a random access memory (RAM) 430. In the RAM 430, various programs and data required for the operation of the device 400 can also be stored. The computing unit 410, the ROM 420, and the RAM 430 are connected to each other via a bus 440. An input / output (I / O) interface 450 is also connected to the bus 440.

[0311] A number of components in the device 400 are connected to the I / O interface 450, including: an input unit 460, such as a keyboard, a mouse, etc.; an output unit 470, such as various types of displays, speakers, etc.; a storage unit 480, such as a disk, an optical disk, etc.; and a communication unit 490, such as a network card, a modem, a wireless communication transceiver, etc. The communication unit 490 allows the device 400 to exchange information / data with other devices through a computer network such as the Internet and / or various telecommunication networks.

[0312] The computing unit 410 may be a variety of general and / or special processing components with processing and computing capabilities. Some examples of the computing unit 410 include, but are not limited to, a central processing unit (CPU), a graphics processing unit (GPU), various dedicated artificial intelligence (AI) computing chips, various computing units running machine learning model algorithms, digital signal processors (DSPs), and any appropriate processors, controllers, microcontrollers, etc. The computing unit 410 performs the various methods and processes described above. For example, in some embodiments, any of the above methods may be implemented as a computer software program, which is tangibly contained in a machine-readable medium, such as a storage unit 480. In some embodiments, part or all of the computer program may be loaded and / or installed on the device 400 via ROM 420 and / or a communication unit 490. When the computer program is loaded into RAM 430 and executed by the computing unit 410, one or more steps of any of the methods described above may be performed. Alternatively, in other embodiments, the computing unit 410 may be configured to perform any of the methods described above in any other appropriate manner (e.g., by means of firmware).

[0313] Various implementations of the systems and techniques described above herein can be implemented in digital electronic circuit systems, integrated circuit systems, field programmable gate arrays (FPGAs), application specific integrated circuits (ASICs), application specific standard products (ASSPs), systems on chips (SOCs), load programmable logic devices (CPLDs), computer hardware, firmware, software, and / or combinations thereof. These various implementations can include: being implemented in one or more computer programs that can be executed and / or interpreted on a programmable system including at least one programmable processor, which can be a special purpose or general purpose programmable processor that can receive data and instructions from a storage system, at least one input device, and at least one output device, and transmit data and instructions to the storage system, the at least one input device, and the at least one output device.

[0314] The program code for implementing the method of the present disclosure may be written in any combination of one or more programming languages. These program codes may be provided to a processor or controller of a general-purpose computer, a special-purpose computer, or other programmable data processing device, so that the program code, when executed by the processor or controller, enables the functions / operations specified in the flow chart and / or block diagram to be implemented. The program code may be executed entirely on the machine, partially on the machine, partially on the machine and partially on a remote machine as a stand-alone software package, or entirely on a remote machine or server.

[0315] In the context of the present disclosure, a machine-readable medium may be a tangible medium that may contain or store a program for use by or in conjunction with an instruction execution system, device, or equipment. A machine-readable medium may be a machine-readable signal medium or a machine-readable storage medium. A machine-readable medium may include, but is not limited to, an electronic, magnetic, optical, electromagnetic, infrared, or semiconductor system, device, or equipment, or any suitable combination of the foregoing. A more specific example of a machine-readable storage medium may include an electrical connection based on one or more lines, a portable computer disk, a hard disk, a random access memory (RAM), a read-only memory (ROM), an erasable programmable read-only memory (EPROM or flash memory), an optical fiber, a portable compact disk read-only memory (CD-ROM), an optical storage device, a magnetic storage device, or any suitable combination of the foregoing.

[0316] To provide interaction with a user, the systems and techniques described herein can be implemented on a computer having: a display device (e.g., a CRT (cathode ray tube) or LCD (liquid crystal display) monitor) for displaying information to the user; and a keyboard and pointing device (e.g., a mouse or trackball) through which the user can provide input to the computer. Other types of devices can also be used to provide interaction with the user; for example, the feedback provided to the user can be any form of sensory feedback (e.g., visual feedback, auditory feedback, or tactile feedback); and input from the user can be received in any form (including acoustic input, voice input, or tactile input).

[0317] The systems and techniques described herein may be implemented in a computing system that includes back-end components (e.g., as a data server), or a computing system that includes middleware components (e.g., an application server), or a computing system that includes front-end components (e.g., a user computer with a graphical user interface or a web browser through which a user can interact with implementations of the systems and techniques described herein), or a computing system that includes any combination of such back-end components, middleware components, or front-end components. The components of the system may be interconnected by any form or medium of digital data communication (e.g., a communication network). Examples of communication networks include: a local area network (LAN), a wide area network (WAN), and the Internet.

[0318] A computer system may include a client and a server. The client and the server are generally remote from each other and usually interact through a communication network. The relationship of client and server is generated by computer programs running on respective computers and having a client-server relationship with each other. The server may be a cloud server, a server of a distributed system, or a server combined with a blockchain.

[0319] It should be understood that the various forms of processes shown above can be used to reorder, add or delete steps. For example, the steps recorded in this disclosure can be executed in parallel, sequentially or in different orders, as long as the desired results of the technical solutions disclosed in this disclosure can be achieved, and this document does not limit this.

[0320] The above specific implementations do not constitute a limitation on the protection scope of the present disclosure. It should be understood by those skilled in the art that various modifications, combinations, sub-combinations and substitutions can be made according to design requirements and other factors. Any modification, equivalent substitution and improvement made within the spirit and principle of the present disclosure shall be included in the protection scope of the present disclosure.

Claims

1. A UAV PPK post-processing method integrating filtering and fuzzy processing, characterized in that: include: Acquire first positioning data received by the base station, and convert the first positioning data and second positioning data into data in a predetermined data format; wherein the second positioning data is positioning data received by the drone; Performing N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format; wherein N is a positive integer greater than or equal to 4; Perform cycle slip detection processing, floating point solution processing, integer ambiguity fixation processing, baseline solution, coordinate conversion and adjustment processing on the first positioning data and the second positioning data after each filtering, and obtain the UAV positioning result corresponding to each filtering; For the drone positioning results corresponding to each filtering, if the types of N drone positioning results are the same, the types of the N drone positioning results are used as the initial positioning result type; if the types of two drone positioning results are the same, the types of the two drone positioning results are used as the initial positioning result type; if the types of N drone positioning results are different, the drone positioning results are sorted according to the filtering order, and the types of the two middle drone positioning results are used as the initial positioning result type; the values ​​of the N drone positioning results are merged to obtain the target value of the drone positioning result; the types of the N drone positioning results are upgraded and downgraded to obtain the target positioning result type of the drone positioning result; Among them, the types of N drone positioning results are upgraded and downgraded, including: Subtract the values ​​of the positioning results of two drones in the types of N drone positioning results to obtain multiple difference values. If all the obtained difference values ​​exceed the initial threshold, the type of the initial positioning result is downgraded; Alternatively, the values ​​of the positioning results of the two drones are subtracted from each other, and if both differences exceed the initial threshold, the type of the initial positioning result is downgraded; Among them, the numerical values ​​of N UAV positioning results are fused, including: According to the values ​​of the positioning results of N drones, and using the acceptance residual sigma as the inverse proportional weighting, the target value of the drone positioning result is obtained; Alternatively, the target value of the UAV positioning result is obtained by taking the positioning results of two UAVs of the same type and using the posterior residual sigma as an inverse proportional weight.

2. The method according to claim 1, characterized in that The performing N times of filtering on the first positioning data and the second positioning data converted according to the predetermined data format comprises: The first positioning data and the second positioning data converted according to the predetermined data format are subjected to four filtering processes.

3. The method according to claim 2, characterized in that Also includes: The first filtering is a forward Kalman filter, and the first epoch is the initial state, and the initial state matrix and the initial variance matrix are parameters of the initial state; wherein the epoch is a frame of data.

4. The method according to claim 3, characterized in that Also includes: The second filtering is Kalman filtering, and the last epoch is the initial state, and the initial state matrix and initial variance matrix are the parameters of the initial state.

5. The method according to claim 4, characterized in that Also includes: The third filtering is Kalman filtering, in which the first epoch of the third filtering is the last epoch, and the initial state matrix and the initial variance matrix are the results of the last epoch of the first filtering. The third filtering is performed halfway, and the third filtering ends at the middle time point epoch.

6. The method according to claim 5, characterized in that Also includes: The fourth filtering is Kalman filtering, in which the first epoch of the fourth filtering is the first epoch, and the initial state matrix and the initial variance matrix are the results of the last epoch of the second filtering. The fourth filtering is only performed halfway, and the fourth filtering ends at the middle time point epoch.

7. The method according to claim 1, characterized in that The integer ambiguity fixing process comprises: Screening satellite frequencies of all full-frequency navigation satellites, selecting a navigation satellite of the first frequency, and performing ambiguity search and fixing on the corresponding ambiguity floating-point solutions obtained by floating-point solution processing; Arrange the parameters of all full-frequency points, assign scores to the parameters according to weights, arrange all scores, perform star-kicking operations according to multiple preset ratios, and perform ambiguity search and fixation on the floating-point solutions of the ambiguities obtained by the floating-point solution processing corresponding to the remaining navigation satellites; For the remaining navigation satellites, single satellites are eliminated, and single satellites are eliminated cyclically until the number of eliminated navigation satellites accounts for a first preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of ambiguities obtained by floating-point solution processing corresponding to the remaining navigation satellites; for the navigation satellites that cannot be fixed among the remaining navigation satellites, two satellites are eliminated, and so on, until the number of eliminated navigation satellites accounts for a second preset percentage of the number of all navigation satellites, and then ambiguity search and fixation are performed on the floating-point solutions of ambiguities obtained by floating-point solution processing corresponding to the remaining navigation satellites; wherein the parameters include navigation satellite altitude angle, signal-to-noise ratio, azimuth angle, and number of continuous tracking epochs; the ambiguity search and fixation are used to determine an ambiguity combination; the ambiguity combination includes integer solutions of the floating-point solutions of each navigation satellite obtained by fixing; Based on probability theory and through statistical testing, a target integer ambiguity combination is searched out from the ambiguity combinations.

8. A UAV PPK post-processing device integrating filtering and fuzzy processing, characterized in that: include: A data receiving and converting module, used for acquiring first positioning data received by the base station, and converting the first positioning data and the second positioning data into data in a predetermined data format; wherein the second positioning data is positioning data received by the drone; A filtering processing module, used to perform N filtering processes on the first positioning data and the second positioning data converted according to the predetermined data format; wherein N is a positive integer greater than or equal to 4; The fixed positioning solution module is used to perform cycle slip detection processing, floating point solution processing, integer ambiguity fixation processing, baseline solution, coordinate conversion and adjustment processing on the first positioning data and the second positioning data after each filtering, so as to obtain the UAV positioning result corresponding to each filtering; The filtering fusion module is used for the UAV positioning results corresponding to each filtering. If the types of N UAV positioning results are the same, the types of the N UAV positioning results are used as the initial positioning result type; if the types of two UAV positioning results are the same, the types of the two UAV positioning results are used as the initial positioning result type; if the types of N UAV positioning results are different, the UAV positioning results are sorted according to the filtering order, and the types of the two middle UAV positioning results are used as the initial positioning result type; the values ​​of the N UAV positioning results are fused to obtain the target value of the UAV positioning result; the types of the N UAV positioning results are upgraded and upgraded to obtain the target positioning result type of the UAV positioning result; Among them, the types of N drone positioning results are upgraded and downgraded, including: Subtract the values ​​of the positioning results of two drones in the types of N drone positioning results to obtain multiple difference values. If all the obtained difference values ​​exceed the initial threshold, the type of the initial positioning result is downgraded; Alternatively, the values ​​of the positioning results of the two drones are subtracted from each other, and if both differences exceed the initial threshold, the type of the initial positioning result is downgraded; Among them, the numerical values ​​of N UAV positioning results are fused, including: According to the values ​​of the positioning results of N drones, and using the acceptance residual sigma as the inverse proportional weighting, the target value of the drone positioning result is obtained; Alternatively, the target value of the UAV positioning result is obtained by taking the positioning results of two UAVs of the same type and using the posterior residual sigma as an inverse proportional weight.

Images