Cycle slip detection method based on combination of improved Kalman filtering and wavelet transform

Abstract

The invention belongs to the technical field of satellite navigation signal processing, and particularly relates to an improved Kalman filtering and wavelet transform combined cycle slip detection method, which comprises the following steps of: firstly, acquiring carrier phase observation data of a global navigation satellite system by using a single-frequency GNSS (Global Navigation Satellite System) receiver, and performing first-order differential preprocessing on the observation data; improving a Kalman filter by adopting a moving average filtering method, calculating a residual error and a standard deviation of differential data for judging and eliminating a gross error, and smoothing the differential data by using the improved Kalman filter; performing three-layer decomposition and reconstruction on the smoothed differential data by utilizing bior3.1 wavelet transform, extracting and calculating a detailed coefficient sum of the wavelet transform, and performing cycle slip detection according to the change difference of different frequency components of the wavelet coefficient; and finally, the detected cycle slip is replaced by the average value of the values before and after the cycle slip to carry out cycle slip repair, so that the data becomes smoother, and the influence of errors on the result is reduced.

Description

Technical Field

[0001] The invention relates to the technical field of satellite navigation signal processing, and in particular to a cycle slip detection method combining improved Kalman filtering with wavelet transform. Background Art

[0002] In the data processing of high-precision satellite navigation systems, ground receivers need to obtain the satellite's carrier phase, which is the main observable quantity for parameter estimation. However, in the carrier phase measurement, some areas are often blocked by obstacles such as high-rise buildings or tall trees, resulting in gross errors and cycle slips in the satellite navigation system's signals, causing the quality of the observation data to deteriorate, which will bring difficulties to the subsequent resolution of whole-cycle ambiguities and also affect the accuracy and reliability of satellite positioning.

[0003] With the rapid development of carrier phase high-precision positioning technology, many cycle slip detection algorithms have been proposed in recent years. The most commonly used algorithms in single-frequency carrier phase cycle slip detection are Doppler integration method, polynomial fitting method, high-order difference method, Kalman filter method and wavelet transform method. The detection effect of the Doppler integral method is affected by the Doppler observation error and the carrier phase observation value; the polynomial fitting method needs to ensure that there is no cycle slip in the fitting window in order to accurately detect the cycle slip. The cycle slip detection capability is limited by the sampling rate and the observation error of the carrier phase observation value, and is very insensitive to the detection of small cycle slips; the high-order difference method is simple and easy to use, but when using this method, the errors caused by clock jumps and epoch intervals must be considered. While amplifying the cycle slip, the noise will be amplified, and it is only suitable for large cycle slip detection; the Kalman filter method needs to properly set the accuracy of the prior observation noise and process noise to ensure the filtering effect, and is only suitable for high sampling rate data, and has certain defects; the wavelet transform method can accurately detect the location of the cycle slip. When the cycle slip occurs frequently, the wavelet transform method can still detect the cycle slip, but the wavelet transform method can only detect the location of the cycle slip and cannot directly calculate the size of the cycle slip. Many solutions have also been proposed in the dual-frequency carrier phase cycle slip detection, such as the ionospheric residual method, the MW combination method, the TurboEdit method, and the adaptive particle filter method. The ionospheric residual method is greatly affected by the multipath effect and the ionospheric change rate. It can detect small cycle slips, but this method has multi-valued properties and is difficult to detect when the carrier phase observations of L1 and L2 have cycle slips at the same time; the MW combination method is affected by pseudorange noise, and there are problems such as inaccurate detection of 1 to 2 small cycle slips and threshold diffusion after cycle slip detection failure; the TurboEdit method is suitable for baselines of arbitrary length and is suitable for cycle slip detection and repair in precise single-point positioning mode, but it is easily interfered by the original data and cannot detect cycle slips within 6 weeks; the adaptive particle filter method is easy to initialize and converges quickly, but is easily affected by the sampling rate. Summary of the invention

[0004] The purpose of the present invention is to provide a cycle slip detection method for improving Kalman filtering combined with wavelet transform, aiming to solve the cycle slip detection problem in carrier phase high-precision positioning technology and improve the positioning accuracy and reliability of satellite navigation systems.

[0005] To achieve the above object, the present invention provides a cycle slip detection method of improving Kalman filtering combined with wavelet transform, comprising the following steps:

[0006] Step 1: Acquire and preprocess the carrier phase observation data of the global navigation satellite system;

[0007] Step 2: Use the moving average filter method to improve the Kalman filter and smooth the data;

[0008] Step 3: Use bior3.1 wavelet transform to perform three-layer decomposition and reconstruction of the data;

[0009] Step 4: Extract and calculate the sum of wavelet coefficients, and detect and repair cycle slips through the difference in changes in different frequency components of the wavelet coefficients.

[0010] Optionally, during the execution of step 1, a single-frequency GNSS receiver is used to continuously observe satellite signals, save the observation data, and calculate the first-order difference of the observation data.

[0011] Optionally, the improved Kalman filter in step 2 uses system noise and observation noise as statistical characteristics, observations as input, and estimated values ​​as output, eliminates interference in signal observations, and estimates the signal to be processed based on the state equation and measurement equation.

[0012] Optionally, the process of smoothing the data in step 2 specifically involves calculating the residual and standard deviation of the differential data for judging gross errors, setting a threshold to identify gross error observation points, and eliminating gross errors.

[0013] Optionally, in the execution process of step 3, the double difference sequence of carrier phase observation values ​​is first used as a detection sequence for cycle slip detection, the detection sequence is decomposed into three layers of detailed coefficients and one approximation coefficient by wavelet, and then the decomposed signal sequence is reconstructed.

[0014] Optionally, the cycle slip detected in step 4 is repaired by replacing it with the average value of the values ​​before and after the cycle slip point. The expression is as follows:

[0015]

[0016] Among them, x jump Indicates the value of the current week jump point, x jump+1 Indicates the value of the next cycle jump point, x jump-1 The week jump point indicates the value at the previous moment.

[0017] The present invention provides a cycle slip detection method of improving Kalman filtering and combining wavelet transform. Firstly, a single-frequency GNSS receiver is used to obtain carrier phase observation data of a global navigation satellite system, and the observation data is preprocessed by first-order difference. Then, a moving average filtering method is used to improve the Kalman filter, and the residual and standard deviation of the differential data are calculated to judge and eliminate gross errors, and the improved Kalman filter is used to smooth the differential data. Then, bior3.1 wavelet transform is used to perform three-layer decomposition and reconstruction on the smoothed differential data, and the detailed coefficients and of the wavelet transform are extracted and calculated. The cycle slip is detected through the difference in the changes of different frequency components of the wavelet coefficients. Finally, the cycle slip is repaired by replacing the detected cycle slip with the average value of the values ​​before and after the cycle slip, so that the data becomes smoother and the influence of errors on the results is reduced. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0019] Figure 1 It is a schematic diagram of the steps of an improved Kalman filter combined with wavelet transform cycle slip detection method of the present invention.

[0020] Figure 2 It is a specific flow chart of a cycle slip detection method of an improved Kalman filter combined with wavelet transform of the present invention.

[0021] Figure 3 It is a schematic diagram of the curve of the observed value and the Kalman filter estimated value without adding the cycle slip and the gross error in the specific embodiment of the present invention.

[0022] Figure 4 This is a graph of cycle slip detection results of observed values ​​without adding cycle slips and gross errors in a specific embodiment of the present invention.

[0023] Figure 5 It is a schematic diagram of a curve of observation values ​​and Kalman filter estimation values ​​with cycle slips and gross errors added in a specific embodiment of the present invention.

[0024] Figure 6 It is a residual sequence diagram when the sampling rate of the observation value with cycle slips and gross errors added is 1 second in a specific embodiment of the present invention.

[0025] Figure 7 This is a graph of cycle slip detection results after filtering out gross errors when the sampling rate of observation values ​​with cycle slips and gross errors is 1 second in a specific embodiment of the present invention.

[0026] Figure 8 It is a residual sequence diagram when the sampling rate of the observation value with cycle slips and gross errors added is 3 seconds in a specific embodiment of the present invention.

[0027] Fig. 9 This is a graph of cycle slip detection results after filtering out gross errors when the sampling rate of observation values ​​of cycle slips and gross errors is 3 seconds in a specific embodiment of the present invention.

[0028] Fig.10 It is a residual sequence diagram when the sampling rate of the observation value with cycle slips and gross errors added is 5 seconds in a specific embodiment of the present invention.

[0029] Fig.11 It is a graph of cycle slip detection results after filtering out gross errors when the sampling rate of observation values ​​of cycle slips and gross errors is 5 seconds in a specific embodiment of the present invention. DETAILED DESCRIPTION

[0030] Embodiments of the present invention are described in detail below, examples of which are shown in the accompanying drawings, wherein the same or similar reference numerals throughout represent the same or similar elements or elements having the same or similar functions. The embodiments described below with reference to the accompanying drawings are exemplary and are intended to be used to explain the present invention, and should not be construed as limiting the present invention.

[0031] See also Figure 1 The present invention provides a cycle slip detection method of improving Kalman filtering combined with wavelet transform, comprising the following steps:

[0032] S1: Acquire and preprocess the carrier phase observation data of the global navigation satellite system;

[0033] S2: Use the moving average filter method to improve the Kalman filter and smooth the data;

[0034] S3: Three-layer decomposition and reconstruction of data using bior3.1 wavelet transform;

[0035] S4: Extract and calculate the sum of wavelet coefficients, and detect and repair cycle slips through the difference in changes in different frequency components of the wavelet coefficients.

[0036] The following is a further explanation of the specific implementation steps. Figure 2 As shown:

[0037] Step S1: Acquire and pre-process the carrier phase observation data of the global navigation satellite system. Use a single-frequency GNSS receiver to continuously observe satellite signals, save the observation data, and calculate the first-order difference of the observation data.

[0038] Step S2: Use the moving average filtering method to improve the Kalman filter, and calculate the residual and standard deviation of the differential data to judge the gross error, set the threshold to identify the gross error observation point, and eliminate the gross error. The improved Kalman filter uses system noise and observation noise as statistical characteristics, system observations as input, and estimated values ​​as output. It can eliminate interference in signal observations and estimate the signal to be processed based on the state equation and measurement equation. Use the filter to smooth the differential data, and when encountering a gross error, it is determined not to update the state estimate.

[0039] The moving average filter is a simple and effective filter commonly used in digital signal processing. It performs convolution operations on the input signal, calculates the average value in a certain number of sample windows, and uses the calculated result as the output of the current point. When the window moves along the signal, the new average value of the data in the window is calculated to obtain a smoothed signal and reduce the random noise of the data. The mathematical formula of the moving average filter is:

[0040]

[0041] In the above formula, x[n] is the input signal, y[n] is the output result of the sliding average filter, and N is the sliding window size, which determines the filtering effect. The larger the N value, the more obvious the filtering effect, but it will increase the signal delay.

[0042] The Kalman filter is an efficient recursive filter that uses the minimum mean square error as the estimation criterion, adopts the state equation of signal and noise, and uses the estimated value at the previous moment and the observed value at the current moment to update the estimate of the state variable and derive the estimated value at the current moment. It is mainly divided into two steps: prediction and update. It requires inputting the initial value of the system and constructing the state equation and observation equation of the system. Initial estimate of the filter state is the mean of the true value of the state at that time, that is:

[0043]

[0044] In the above formula, is the initial estimate of the system state, x 0 The initial true value of the system state.

[0045] The corresponding initial estimate of the state mean square error P 0 satisfy:

[0046]

[0047] In the above formula, P 0 represents the initial mean square error of the system state, Initial estimate of the system state, x 0 The initial true value of the system state.

[0048] In the prediction phase, the Kalman filter uses the estimated value of the previous state and the state equation of the system to make an estimate of the current state. The mathematical formula is:

[0049]

[0050] In the above formula, It represents the Kalman filter response to the system state x at time k k The prior estimate of It represents the Kalman filter for the system state x at time k-1 k-1 The optimal estimate of A represents the state transfer matrix, B represents the relationship matrix between the system input and the system state, and u k-1 Represents the input quantity of the system at time k-1. represents the covariance matrix of the prior estimation error at time k, P k-1 represents the covariance matrix of the estimation error at time k-1, and Q represents the covariance matrix of the process noise. When predicting the state, since the value of the process noise is an unknown number, the covariance matrix Q of the process noise is added to the calculation formula of the mean square error matrix of the state prior estimation value to reduce the correctness of the state prior estimation value and maintain the reliability of the mean square error matrix. Usually, after the prediction process is completed, the mean square error value of the system state estimation usually increases.

[0051] In the update phase, the filter uses the current actual measurement value to correct the predicted value obtained in the prediction phase. The mathematical formula is:

[0052]

[0053] In the above formula, K k represents the Kalman gain, represents the covariance matrix of the prior estimation error at time k, C represents the relationship matrix between the actual measurement value and the system state, R represents the covariance matrix of the measurement noise, It represents the Kalman filter response to the system state x at time k k The best estimate of It represents the Kalman filter response to the system state x at time k k The prior estimate of y k Indicates the actual measured value, P k represents the covariance matrix of the estimation error at time k, and I represents the unit matrix. When updating the state, since the state estimate is verified by the actual measurement value, its mean square error value decreases and its reliability increases. The Kalman filter integrates and balances the two information of the prior estimate and the actual measurement, so that the state estimate after the measurement update has the smallest mean square error.

[0054] Step S3: Use bior3.1 wavelet transform to perform three-layer decomposition and reconstruction of data. bior3.1 wavelet is a biorthogonal wavelet suitable for multi-scale decomposition, smoothing and denoising of signals. First, the double difference sequence of carrier phase observations is used as the detection sequence for cycle slip detection. The detection sequence is decomposed into three layers of detailed coefficients and one approximation coefficient by wavelet, and then the decomposed signal sequence is reconstructed.

[0055] Wavelet transform definition: Let ψ(t) = L 2 (R), if the conditions are met:

[0056]

[0057] In the above formula, ψ(t) is a basic wavelet. * =R-{0} represents all non-zero real numbers, and ψ(ω) is the Fourier transform of ψ(t).

[0058] Scale and translate ψ(t), let the scale factor be a and the translation factor be τ, then:

[0059]

[0060] In the above formula, ψ a,τ (t) is the wavelet basis function that depends on the parameters a and τ.

[0061] The bior wavelet is a type of Biorthogonal wavelet system. It is symmetrical, compactly supported, and has the characteristics of biorthogonality and linear phase. The general expression is biorNr.Nd, and there is no clear mathematical expression. Nr represents the order of the reconstruction filter, and Nd represents the order of the low-pass decomposition filter. Biorthogonal wavelets are mainly used for signal decomposition and reconstruction. The use of the Biorthogonal wavelet with Nr=3 and Nd=1 combination has high stability and good detection effect for cycle slip detection.

[0062] Step S4: Extract the detailed coefficients of the first, second and third layers of the wavelet and the approximation coefficients of the third layer, calculate the sum of the detailed coefficients by adding the detailed coefficients of different layers, and detect the cycle slips through the difference in the changes of different frequency components of the wavelet coefficients. For the detected cycle slips, they are repaired by replacing them with the average value of the values ​​before and after the cycle slip point to improve the data quality and accuracy. The mathematical formula is:

[0063]

[0064] In the above formula, x jump Indicates the value of the current week jump point, x jump+1 Indicates the value of the next cycle jump point, x jump-1 The week jump point indicates the value at the previous moment.

[0065] Furthermore, the present invention is further explained and supplemented by specific embodiments:

[0066] Implementation Example 1: Analysis of Cycle Slip Detection Effects at Different Sampling Rates

[0067] In order to verify the cycle slip detection effect of the improved Kalman filter combined with wavelet transform method proposed in the present invention, a simulation experiment is conducted using data observed in a certain place. The original observation data collected a total of 3000 epochs of GPS data. The original observation data is firstly processed by first-order difference, and then gross errors and cycle slips are artificially added. Then, the method proposed in the present invention is used to remove gross errors and detect cycle slips when the sampling rate of the observation data is 1 second, 3 seconds, and 5 seconds, respectively, to verify the performance of the proposed cycle slip detection method.

[0068] In order to test the sensitivity of the present invention to cycle slip detection, the original observation data without adding simulated gross errors and cycle slips are selected for detection. Figure 3 As shown in Figure 1, the improved Kalman filter model has a high degree of fit between the original observation value estimation curve without bolding errors and cycle slips and the original observation value curve without bolding errors and cycle slips, except for a certain gap at the beginning of iterative estimation. Figure 4 From the residual sequence shown in , it can be concluded that when no cycle slip occurs, the residual value of the detection model is less than 0.6. It can be seen that the model has a high sensitivity to cycle slip detection and can detect small cycle slips of less than 1 cycle in the differential sequence. In order to test the smoothing effect of the improved Kalman filter model on the observed data, gross errors and cycle slips are first added to the original observed data, and then the model is used to smooth the data, such as Figure 5 As shown in the figure, the observation data with gross errors and cycle slips have some spikes, which will eventually affect the fixation of the whole cycle ambiguity and reduce the positioning accuracy. The improved Kalman filter model smoothes the original observation values ​​with gross errors and cycle slips. Figure 5 It can be seen that the estimated value curve eliminates the spikes in the original observations with gross errors and cycle slips, removes the gross errors, and the smoothing effect of the model is better.

[0069] In order to verify the accuracy of the effectiveness of the present invention in detecting cycle slips, when the sampling rate is 1 second, different sizes of gross errors and cycle slips are artificially added to certain epochs of the first-order difference sequence, as shown in Table 1.

[0070] Table 1. The addition of cycle slips to differential sequences

[0071]

[0072] Depend on Figure 6As shown in FIG. 1 , when the gross errors are not eliminated, there are gross errors at the 600 epoch and the 800 epoch of the data. If the gross errors are not processed, the subsequent cycle slip detection effect will be affected. After the data is smoothed using the model of the present invention, as shown in FIG. Figure 7 As shown, gross errors are completely eliminated. When the sampling rate is 1 second, the present invention can detect cycle jumps of different sizes and special cycle jumps, such as small cycle jumps of 1 cycle or more, continuous epoch cycle jumps, and large cycle jumps.

[0073] In order to verify the cycle slip detection performance of the present invention under low sampling rate conditions, gross errors and cycle slips of different sizes are artificially added to certain epochs of the first-order difference sequence when the sampling rates are 3 seconds and 5 seconds, as shown in Tables 2 and 3, respectively.

[0074] Table 2: Adding cycle slips to differential sequences

[0075]

[0076] Depend on Figure 8 As shown in Figure 1, when the rough errors are not removed, there are gross errors at 150 epochs and 400 epochs in the data. After smoothing the data using the model of the present invention, as shown in Figure 1 Fig. 9 As shown, gross errors are completely eliminated. When the sampling rate is 3 seconds, the present invention can detect cycle jumps of different sizes and special cycle jumps, such as cycle jumps of 2 cycles or more, continuous epoch cycle jumps, and large cycle jumps of more than 2 cycles.

[0077] Table 3. The addition of cycle slips to differential sequences

[0078]

[0079] Depend on Fig.10 As shown in FIG. 1 , when the gross errors are not eliminated, there are gross errors at the 400th and 550th epochs of the data. After smoothing the data using the model of the present invention, as shown in FIG. Fig.11 As shown, gross errors are completely eliminated. When the sampling rate is 5 seconds, the present invention can detect cycle jumps of different sizes and special cycle jumps such as cycle jumps of 3 cycles or more, continuous epoch cycle jumps, and large cycle jumps of more than 3 cycles.

[0080] In summary, the present invention provides a cycle slip detection method of an improved Kalman filter combined with a wavelet transform. First, a single-frequency GNSS receiver is used to obtain the carrier phase observation data of the global navigation satellite system, and the observation data is preprocessed by a first-order difference. Then, a moving average filtering method is used to improve the Kalman filter, and the residual and standard deviation of the differential data are calculated to judge and eliminate gross errors, and the filter is used to smooth the differential data. Then, the bior3.1 wavelet transform is used to perform a three-layer decomposition and reconstruction on the smoothed differential data, and the detailed coefficients and of the wavelet transform are extracted and calculated. The cycle slip is detected by the change difference of different frequency components of the wavelet coefficients. Finally, the cycle slip is repaired by replacing the detected cycle slip with the average value of the values ​​before and after the cycle slip, so that the data becomes smoother and the influence of errors on the results is reduced.

[0081] What is disclosed above is only a preferred embodiment of the present invention, and it certainly cannot be used to limit the scope of rights of the present invention. Ordinary technicians in this field can understand that all or part of the processes of the above embodiment and equivalent changes made according to the claims of the present invention still fall within the scope of the invention.

Claims

1. A cycle slip detection method based on improved Kalman filtering combined with wavelet transform, characterized in that: The following steps are involved: Step 1: Acquire and preprocess the carrier phase observation data of the global navigation satellite system; Step 2: Use the moving average filter method to improve the Kalman filter and smooth the data; Step 3: Use bior3.1 wavelet transform to perform three-layer decomposition and reconstruction of the data; Step 4: Extract and calculate the sum of wavelet coefficients, and detect and repair cycle slips through the difference in changes in different frequency components of the wavelet coefficients.

2. The cycle slip detection method of improved Kalman filtering combined with wavelet transform as claimed in claim 1, characterized in that: During the execution of step 1, a single-frequency GNSS receiver is used to continuously observe satellite signals, save the observation data, and calculate the first-order difference of the observation data.

3. The cycle slip detection method of improved Kalman filtering combined with wavelet transform as claimed in claim 2, characterized in that: The improved Kalman filter in step 2 uses system noise and observation noise as statistical characteristics, observations as input, and estimates as output. It eliminates interference in signal observations and estimates the signal to be processed based on the state equation and measurement equation.

4. The cycle slip detection method of improved Kalman filtering combined with wavelet transform as claimed in claim 3, characterized in that: The process of smoothing the data in step 2 is to calculate the residual and standard deviation of the differential data to judge the gross error, set the threshold to identify the gross error observation points, and eliminate the gross error.

5. The cycle slip detection method of improved Kalman filtering combined with wavelet transform as claimed in claim 4, characterized in that: In the execution process of step 3, the double difference sequence of the carrier phase observation value is first used as the detection sequence for cycle slip detection, the detection sequence is decomposed into three layers of detailed coefficients and one approximation coefficient by wavelet, and then the decomposed signal sequence is reconstructed.

6. The cycle slip detection method of improved Kalman filtering combined with wavelet transform as claimed in claim 5, characterized in that: In step 4, the detected cycle slip is repaired by replacing it with the average value of the values ​​before and after the cycle slip point. The expression is as follows: Among them, x jump Indicates the value of the current week jump point, x jump+1 Indicates the value of the next cycle jump point, x jump-1 The week jump point indicates the value at the previous moment.

Images