A New Method for Multipath Filtering in GPS Static High-Precision Positioning
Next Article in Journal
A Low-Cost Chamber Prototype for Automatic Thermal Analysis of MEMS IMU Sensors in Tilt Measurements Perspective
Previous Article in Journal
Achieving Tiered Model Quality in 3D Structure from Motion Models Using a Multi-Scale View-Planning Algorithm for Automated Targeted Inspection
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A New Method for Multipath Filtering in GPS Static High-Precision Positioning

School of Electronic Engineering, Beijing University of Posts and Telecommunications, No. 10 Xitucheng Road, Beijing 100876, China
*
Author to whom correspondence should be addressed.
Sensors 2019, 19(12), 2704; https://doi.org/10.3390/s19122704
Submission received: 23 April 2019 / Revised: 12 June 2019 / Accepted: 13 June 2019 / Published: 16 June 2019
(This article belongs to the Section Remote Sensors)

Abstract

:
It is well known that multipath is one of the main sources of errors in GPS static high precision positioning of short baselines. Most algorithms for reducing multipath manipulate the GPS double difference (DD) observation residuals as input signal in GPS signal processing. In the traditional multipath mitigation methods, applying the wavelet transform (WT) to decompose the GPS DD observation residuals for identifying the multipath disturbance cannot effectively filter out the white noise of the high frequency part of the signal, and it is prone to edge effect. In this paper, for extracting multipath, a wavelet packet algorithm based on two-dimensional moving weighted average processing (WP-TD) is proposed. This algorithm can not only effectively filter out the white noise of the high frequency part of the signal, but also weaken the influence of the edge effect. Furthermore, considering the repeatability of multipath error in static positioning, we propose a method for determining the level of wavelet packet decomposition layers which make multipath extraction more effectively. The experimental results show that the corrected positioning accuracy is 14.14% higher than that of the traditional wavelet transform when applying the obtained multipath to DD coordinate sequences for position correction.

1. Introduction

In GPS static high-precision positioning of short baselines, errors such as ionospheric delay, tropospheric delay, satellite orbit error, receiver clock error, and satellite clock error can be eliminated or attenuated by differential techniques while multipath disturbance are not correlated at both ends of the baseline which make it impossible to be eliminated by differential techniques and become the main source of error affecting the positioning accuracy [1,2].
GPS multipath disturbance occurs when GPS signals travel from a satellite to a receiver via several paths due to reflection or diffraction of signals by nearby obstacles [3]. Multipath can distort signal modulation and carrier phase, thereby reducing the accuracy of GPS static high-precision positioning. In addition, this effect may hinder the fixation of ambiguities and even lead to erroneous solutions.
In order to mitigate and eliminate multipath, a variety of different strategies have been proposed, including appropriate site selection, hardware-based methods and software-based methods.
Appropriate site selection is the simplest and the most effective way to solve the multipath effect. The construction of the base station should be chosen in an open area, avoiding various obstacles causing excessive reflection and obtaining as many visible satellite signals as possible.
PozoRuz et al. (1998) proposed a method calculating the positions in the plane using the three highest satellites based on a new satellite selection criterion which can reduce multipath [4].
The hardware-based method refers to improvements of antenna and receiver design.
Scire Scappuzzo et al. (2009) designed a low-multipath wideband GPS antenna with cutoff or non-cutoff corrugated ground plane which can operate uniformly between 1.15 GHz and 1.60 GHz while maintaining the required low-multipath performance in the whole bandwidth [5]. Ray et al. (1999) developed a system for reducing the effect of carrier-phase multipath on static GPS applications using multiple closely spaced antennas [6]. Bryan R. Townsend et al. (1994) introduced a new tracking loop which takes full advantage of the Narrow Correlator spacing design, but, in addition, is much more resistant to multipath effects on the correlation function and thereby reduces the multipath bias on the pseudorange measurements [7]. R.D.J. van Nee et al. (1992) designed a specific receiver structure which simultaneously estimates the parameters of line-of-sight plus multipath signals for reducing GPS code and carrier multipath errors [8].
However, the hardware is inaccessible for most users and it is expensive. Despite this, we cannot deny that hardware methods have better performance in other aspects, such as stability and efficiency. In contrast, software methods are more convenient and cheaper for most people.
The software-based method mainly focuses on the post-processing of signal by using algorithms and filters.
Axelrad et al. (1996) proposed an improved technique which can adaptively estimate the spectral parameters (frequency, amplitude, phase offset) of multipath in the associated signal-to-noise ratio (SNR) [9]. Linlin Ge et al. (2000) proposed an adaptive finite-duration impulse response filter, based on a least-mean-square algorithm, that has been developed to derive a relatively noise-free time series from the CGPS results [10]. There are some other filters, such as Kalman filters (Nce and Sahin 2000) [11], FIR filters (Han and Rizos (1997)) [12], and Vondrak filters (Vondrak (1977) [13], have also been developed to reduce GPS multipath effects.
In addition, wavelet transform has a good localization property within both the time and frequency domains. Particularly, after Stephane Mallat proposed the multi-resolution algorithm, also named Mallat algorithm [14], wavelet transform is widely applied signal processing, including extracting multipath.
Huang et al. (2003) applied wavelet transform in dynamic deformation monitoring for high-rise buildings [15]. Zhang and Bartone (2004) applied the wavelet technique to mitigate errors for satellite based navigation systems which can mitigate multipath error in a real-time conditions [16]. E. M. Souza, J. F. G. Monico (2004) apply the wavelet transform to decompose the pseudorange and carrier phase DD signals to separate the high frequencies, which are due to multipath from long delays, and the low frequencies effects, associated with multipath from short delays [17]. P. Zhong et al. (2008) proposed a method based on the technique of cross-validation for automatically identifying wavelet signal layers is developed and used for separating noise from signals in data series, and applied to mitigate GPS multipath effects [3]. M. R. Azarbad and M. R. Mosavi (2014) proposed a new multipath mitigation method based on the wavelet transform (WT). The method uses the stationary wavelet transform (SWT) to decompose the double difference residuals [18]. Lawrence Lau (2017) proposed a generic and robust three-level wavelet packets based denoising method for repeat-time-based carrier phase multipath filtering in relative positioning; the results show that the wavelet packets based method is better than the DWT-based method in the repeat time-based multipath filtering [19]. Souza et al. (2017) proposed a new approach for structure monitoring from GPS multipath effect and wavelet spectrum, and the experience investigated the feasibility of using wavelet spectra analysis of the multipath signal to monitor structure movement [20].
As we know, based on the Mallat algorithm, a signal can be decomposed into scale space and wavelet space. With the increase of the decomposition levels, the scale space would be decomposed continuously, whereas the wavelet space cannot be decomposed, which results in a low resolution in high frequency regions of wavelet transform [21]. This means that the multipath in the high frequency part of the signal cannot be effectively extracted by traditional wavelet transform.
Hence, wavelet packet transform (WPT) was pioneered by Coifman et al. based on wavelet transform [22]. WPT keeps the property of wavelet transform naturally, good localization property within both the time and frequency domains. This property greatly improves the accuracy of WPT in analyzing non-stationary and non-periodic signals. In addition, WPT effectively mitigates the defects in wavelet transform of low resolution in high frequency regions by decomposing both the scale space and wavelet space [21]. This means we can further subdivide the high frequency portion of signal and extract multipath more effectively.
The multipath extraction efficiency partially depends on the chosen mother wavelet for WPT. In this paper, we won’t attempt to compare the performance of mother wavelets on retrieving multipath errors from noisy coordinate residual sequences. E. M. Souza et al. (2004) had verified that SYM12, SYM8 and DAUB8 presented the better performance for reconstruction of the GPS DD signal [17]. In this paper, we will choose DAUB8 as the mother wavelet for WPT.
The other challenging issue related to the multipath extraction efficiency is the level of decomposition which depends on the sub-band frequency in which multipath lies [18]. Thence, we also propose a method for determining the number of wavelet packet decomposition levels considering the repeatability of multipath error in static positioning.
In addition, edge effects in signal processing using wavelet packet transform is another un-neglected problem. P. Zhong et al. (2008) choose about 70% of the data in the middle of the observational series for cross-validation to prevent edge effects due to poorer filtering results at the ends of a data series [3]. This method solves the edge effect problem in a certain condition, but, when the coordinate sequence residual is short, it will become invalid. Therefore, referring to the idea of moving weighted average method and the principle of bilateral filters in image processing [23,24], we propose a new algorithm named the two-dimensional moving weighted average algorithm. The double-difference coordinate residual sequences are smoothed by the proposed algorithm to achieve the purpose of edge-preserving and pre-denoising.
In general, we propose a new method named wavelet packet algorithm based on two-dimensional moving weighted average processing (TDMWA), compared with the traditional wavelet algorithm, which can not only more effectively mitigate the multipath of the DD coordinate residual sequences, but also effectively weaken the influence of the edge effect.
This study is organized as follows. The principle of position correction based on multipath periodicity is introduced in the Section 2. The theory of the two-dimensional moving weighted average algorithm and the preprocessing of smoothing the double-difference coordinate sequence residual obtained by static observation for three consecutive days using TDWMA algorithm are introduced in Section 3. The theory of wavelet packet transform, the theory of correlation analysis of smoothed double-difference coordinate residual sequences, the specific steps for determining the number of wavelet packet decomposition layers and the principle of wavelet packet threshold denoising will be introduced in Section 4. The overall program of mitigating the multipath in GPS static high-precision positioning will be introduced in Section 5. The experimental verification will be introduced in Section 6 and Section 7. The data set used in Section 6 is the simulation data for verifying the feasibility of the proposed method. Section 7 uses the GPS data acquisition system built to collect the GPS data to further verify the feasibility of the actual environment and the denoising performance of the Vondrak algorithm, traditional WT algorithm and the WP-TD algorithm is compared from the correlation and root mean square error. Conclusions are given in Section 8.

2. The Principle of Position Correction Based on Multipath Periodicity

In this paper, we take the double difference carrier phase observation residuals obtained from the reference day (Day 1) as the correction to the observations on the adjacent days (Day 2 and Day 3) considering that multipath is repetitive in static observations of the fixed point. Methods for position correction based on multipath repeatability can be found in Khelifa et al. (2011) [25], Ye et al. (2013) [26], and Lawrence Lau et al. (2017) [19].
The double difference carrier phase observations [27] indicated as Equation (1):
Φ r b n m = λ 1 ρ r b n m + N r b n m + ε ϕ , r b n m + M ϕ , r b n m ,
where Φ r b n m denotes the double difference carrier phase observation between satellites n and m, and stations r and b. ρ r b n m denotes the geometric distance from the satellite center to the antenna phase center. N r b n m denotes the double difference integer ambiguity. ε ϕ , r b n m denotes the double difference carrier phase measurement noise. M ϕ , r b n m denotes the double difference carrier phase multipath. λ denotes the carrier wavelength.
In GPS high-precision relative positioning of short baseline (less than 3 km). Since the DD carrier phase multipath error is always less than a quarter of the carrier wavelength, the observation residuals obtained on the reference day do not take the double difference integer ambiguity into consideration [21].
Based on the double difference carrier phase observations, the positioning solution of each epoch on the reference day indicated as Equation (2):
[ X ˙ ] R e = [ X ¯ ] + [ M ] R e + [ ε ] R e ,
where X ˙ denotes the best estimated positioning solution of the fixed point in static observation. X ¯ denotes the known positioning solution of the fixed point in static observation. M denotes the multipath, and ε denotes the observation noise (white noise). Subscript R e indicates the reference day.
Moving the known item of Equation (2) to the left side to get Equation (3):
[ X ˙ ] R e [ X ¯ ] = [ M ] R e + [ ε ] R e .
The right side of the Equation (3) is the positioning residuals that consist of multipath and observation noise (white noise), which will be regarded as the DD coordinate residual sequences in the next section.
Because the repeat-time-based multipath increases the observation noise level in the positioning solutions of the fixed point on the adjacent days, it is necessary to eliminate the noise of the DD coordinate residual sequences and obtain the noiseless multipath as the correction of the positioning solutions of the adjacent days.
The corrected position of each epoch on the adjacent day indicated as Equation (4):
[ X ^ ] A d + [ ε ] A d = [ X ˙ ] A d [ M ] R e ,
where X ^ denotes the multipath-corrected positioning solutions. Subscript A d indicates the adjacent days.

3. Preprocessing of the DD Coordinate Residual Sequences by the TDMWA Algorithm

3.1. The Basic Principle of the TDMWA Algorithm

In GPS static high-precision positioning, the multipath tends to be stable between the adjacent epochs, but the DD coordinate residual sequences fluctuation intensified between adjacent epochs due to the existence of white noise. Even worse, there exist poorer filtering results at the end of the DD coordinate residual sequences because of the edge effects.
In order to weaken the influence of edge effects and weaken the white noise, we propose a new method named two-dimensional moving weighting average algorithm (TDMWA) for smoothing the DD coordinate residual sequences. Considering that the white noise is independent and obeys the normal distribution, the algorithm performs moving weighted average processing on each of the observation epoch in the time domain and the value domain, respectively, which can be described as the following formulation equation:
{ X r T = α [ T 2 ] x j [ T 2 ] + α [ T 2 ] 1 x j [ T 2 ] + 1 + α 0 x j + α 1 x j + 1 + + α [ T 2 ] x j + [ T 2 ] T X r N r = k = 1 , x k N r β k x k N r ,
where X r T is the time domain component of the original DD coordinate residual sequence of the r-th observation epoch, and T is the moving average period of the time domain component, taken as an odd number, and T I ( I is the total number of observation epochs). X r N r is the value domain component of the original DD coordinate residual sequence of the r-th observation epoch, and N r is the moving average period of the value domain component, and N r I , r represents the current epoch, and r I .
[] denotes an integer, x j is the decentralized epoch coordinate of the time domain component which determines the smoothed DD coordinate residual of each observation epoch in the time domain, x k is the decentralized epoch coordinate of the value domain component which determines the smoothed DD coordinate residual of each observation epoch in the value domain, and is the set consisting of the decentralized epoch coordinates.
α j is the weighting coefficient of the j-th delay epoch in the moving average period T , and β k is the weighting coefficient of the k-th decentralized epoch in the moving average period N r , then
{ α j = α 0 , α 1 , , α [ T 2 ] , α 0 + 2 ( α 1 + α [ T 2 ] ) = 1 β k = β 1 , ,   β [ N r ] , β 1 + + β [ N r ] = 1 .

3.2. The Process of Smoothing the DD Coordinate Sequence Residual by the TDMWA Algorithm

3.2.1. Determine the Moving Average Periods T and the Decentralized Epoch Coordinate x j of Time Domain Component

Wherein the time domain components have the same moving average period T , and T is consistent with the influence period of the multipath effect.
The current epoch coordinate x r is taken as the median value of the corresponding decentralized epoch coordinate of the r-th observation epoch and the number of the corresponding decentralized epochs is T.

3.2.2. Determine the Moving Average Periods N r and the Decentralized Epoch Coordinate x k of Value Domain Component

Based on the idea of bilateral filters in image processing, the moving average period N r of each value domain component and the corresponding decentralized epoch coordinate x k is determined according to the discrete distribution of the observation epochs.
Firstly, we need calculate the average value of all the observation epoch coordinate x ¯ I
x ¯ I = i = 1 I x i I ,
where x i represents the all observation epoch coordinates, and i I .
Then, calculating the standard deviation of the observation epoch coordinate, and using the standard deviation as the threshold θ I for selecting decentralized epoch coordinate constituting :
{ σ I = i = 1 I ( x i x ¯ I ) 2 I θ I = σ I .
Then, the value of the current observation epoch coordinate x r becomes the new mean of all the observation epochs,
x ¯ I = x r , r I .
Next, calculating the standard deviation σ i of all the observation epochs coordinate x i according to the new mean obtained in Equation (9),
σ i = i = 1 I ( x i x r ) 2 I .
The observation epoch coordinate whose standard deviation σ i is not greater than the threshold θ I is used as the decentralized epoch coordinate x k , which will be selected to determine the smoothed value of the r-th observation epoch in the value domain, and the number of the decentralized epoch coordinate x k is the moving average period N r , and N r I .

3.2.3. Determine the Weighting Coefficient α j and β k

In the time domain, calculating the epoch delay of each decentralized epoch coordinate x j and the current observation epoch coordinate x r , and weighting the corresponding decentralized epoch coordinate x j according to the size of the delayed epoch, the smaller the epoch delay, the greater the influence of the proximity effect, and the greater the weight of the corresponding decentralized epoch coordinate x j .
Let τ ¯ be the average epoch delay of each decentralized epoch coordinate x j relative to the current observation epoch x r , Δ τ j is the epoch delay of each decentralized epoch coordinate x j relative to the current observation epoch coordinate x r , and ϑ j is the variation coefficient of each decentralized epoch coordinate x j relative to the current observation epoch coordinate x r , and j T , j = t [ T 2 ] , , 0 , , t + [ T 2 ] , the weighting coefficient of the decentralized epoch coordinate x j in the current epoch r can be calculated by the following formula:
{ ϑ j = ( Δ τ j τ ¯ ) 2 τ ¯ τ ¯ = i = t [ T 2 ] t + [ T 2 ] Δ τ j T α j = ϑ j j = t [ T 2 ] t + [ T 2 ] ϑ j .
In the value domain, the difference between each decentralized epoch coordinate x k and the current observation epoch coordinate x r is calculated, and the weighting coefficient are assigned to the respective decentralized epoch coordinate x k according to the magnitude of the difference, and, the smaller the difference, the greater the influence of similarity effect, the greater the weight of the corresponding decentralized epoch coordinate x k . Let μ ¯ be the average of the magnitude of the difference of the decentralized epoch coordinate x k , Δ μ k is the difference of the each decentralized epoch coordinate x k with respect to the current observation epoch coordinate x r , and ϑ k is the variation coefficient of the each decentralized epoch coordinate x k . Then, the weighting coefficient of the decentralized epoch coordinate x k in the current epoch r can be calculated by the following formula:
{ ϑ k = ( Δ μ k μ ¯ ) 2 μ ¯ μ ¯ = k = 1 , N r Δ μ k N r β k = k k = 1 , N r ϑ k .

3.2.4. Calculating the Time Domain Component X r T and the Value Domain Component X r N r

Then, we substitute T , N r , x j , x k ,   α j , β k into Equation (5), the current observation epoch coordinate x r is subjected to moving weighted averaging processing in the time domain and in the value domain, respectively, and the time domain component X r T and the value domain component X r N r of the current epoch coordinate are obtained.

3.2.5. Calculating the Smoothed Observation Epoch X r ^

Next, the time domain component X r T and the value domain component X r N r of the current observation epoch r are comprehensively weighted:
X r ^ = γ T X r T + γ N r X r N r ,
where γ T is the weighting coefficient of the time domain component, γ N r is the weighting coefficient of the value domain component, and X r ^ is the new value smoothed in the current observation epoch r .
In the current epoch r, the difference of the time domain component and the current observation epoch coordinate x r , and the difference of the value domain component and the current observation epoch coordinate x r are calculated, respectively, and the values of the corresponding components are weighted according to the magnitude of the difference.
Assuming that the average value of the current epoch coordinate x r in two-dimensional space is Υ ¯ , the variation coefficient of the time domain component is ϑ T , and the variation coefficient of the value domain is ϑ N r . Then, the weighting coefficients of the time domain component and the value domain component of the current observation epoch r are calculated by the following formula:
{ ϑ T = ( Δ X T Υ ¯ ) 2 Υ ¯ ϑ N r = ( Δ X N r Υ ¯ ) 2 Υ ¯ γ T = ϑ T ϑ T + ϑ N r γ N r = ϑ N r ϑ T + ϑ N r ,
where Υ ¯ = X r T + X r N r 2 , Δ X T = X r T Υ ¯ , Δ X N r = X r N r Υ ¯ . Finally, we substitute the obtained γ T and γ N r into Equation (13) to obtain the new value X r ^ of the current observation epoch r , which is the smoothed observation epoch coordinate.
Following the steps as above, this method is expected to solve the edge effect and achieve the purpose of edge preservation and denoising, but using this method will retain too much high frequency information, so that the white noise in high frequency will not be completely eliminated, so this paper will add wavelet packet transform based on this method to further eliminate the white noise in high frequency.

4. Time-Frequency Conversion by WPT

4.1. The Theory of WPT

From the introduction, we have known that wavelet packet transform provides a more sophisticated analysis of the signal. Wavelet packet transform divides the time-frequency space into more detail, and it has higher resolution for the high-frequency part of the signal than the binary wavelet transform. Moreover, based on the theory of wavelet analysis, it introduces the concept of optimal basis selection. That is, after the frequency band is divided into multiple levels, according to the characteristics of the analyzed signal, the optimal basis function is adaptively selected to match the signal to improve the signal analysis capability.
In order to further eliminate the white noise of the high frequency part of the smoothed DD coordinate residual sequences, which wavelet packet transform is applied to decompose, the wavelet packet decomposition algorithm [21] can be indicated as Equation (15):
{ d M 2 ω 1 = k Z H ( k 2 i ) d M 1 ω d M 2 ω = k Z G ( k 2 i ) d M 1 ω ,
where d M 2 ω 1 denotes the 2ω − 1 wavelet packet coefficient of the M-th layer, d M 1 ω represents the ω-th wavelet packet coefficient on the M-1 layer, and p M 2 ω represents the 2ω-th wavelet packet coefficients on the M-th layer, G is the decomposition filter relation to scale function, H is the decomposition filter of wavelet function, and M 1 ,   i = 1 ~ I , ω = 1 ~ 2 M .
The wavelet packet reconstruction algorithm [21] can be indicated as Equation (16):
d M ω = 2 [ k Z h ( i 2 k ) d M + 1 2 ω 1 + k Z g ( i 2 k ) d M + 1 2 ω ] ,
where d M + 1 2 ω 1 represents the 2 ω 1 wavelet packet coefficient on the M + 1 layer, d M + 1 2 ω indicates the 2 ω wavelet packet coefficients on the J + 1 layer, g is the reconstruction filter related to scaling function, and h is the reconstruction filter related to wavelet function.
From the perspective of time-frequency conversion, the DD coordinate residual sequences after the first WPT decomposition could be converted into two equal frequency bands by Equation (15), i.e., ( 0 ~ N f 2 ) and ( N f 2 ~ N f ) , where N f is Nyquist frequency. The WPT coefficients d M 1 ω are divided into d M 2 ω 1 and d M 2 ω as well. Moreover, both ( 0 ~ N f 2 ) and ( N f 2 ~ N f ) can be divided continuously.
An example of three level decomposition of WPT is shown as Figure 1.
In Figure 1, A represents a low frequency, D represents a high frequency, and the serial number at the end represents the number of layers of wavelet decomposition.
The decomposition relationship can be expressed as Equation (17):
S = AAA 3 + DAA 3 + ADA 3 + DDA 3 + AADD 3 + DAD 3 + ADD 3 + DDD 3 .
After the DD coordinate residual sequences are decomposed according to the structure shown in Figure 1, using the thresholding denoising mentioned in the Section 4.4, the wavelet packet coefficients that satisfy the condition will be used to reconstruct the compressed DD coordinate residual sequences while others will be set to zero.

4.2. The Theory of Correlation Analysis

In view of the strong correlation of multipath in the same time between the adjacent days in static positioning, though the smoothed DD coordinate residual sequences still contain residual white noise, but the multipath error is dominant, which contributes to strong correlation between the smoothed DD coordinate residual sequences.
We have the decomposed DD coordinate residual sequences in Section 3.1, which can be applied to correlation analysis for determining the best decomposition level of the WPT.
In the three consecutive days, regard one day as the reference day and the rest as the adjacent days.
A is the wavelet packet coefficient of the low-frequency part of the decomposed DD coordinate sequences of the reference day, and B is the wavelet packet coefficient of the corresponding low-frequency part of the decomposed DD coordinate sequence of the adjacent day, and ρ A B is the cross-correlation coefficient.
The formula for calculating the correlation coefficient can be expressed as Equation (18):
ρ AB = c o v ( A , B ) D ( A ) D ( B ) ,
where c o v ( A , B ) = ω = 1 2 M ( a ω u a ) ( b ω u b ) is the covariance of variables A and B , D ( A ) = ω = 1 2 M ( a ω u a ) 2 is the variance of the variable A , D ( B ) = ω = 1 2 M ( b ω u b ) 2 is the variance of the variable B , u a = 1 2 M ω = 1 2 M a ω is the average value of the variable A , u b = 1 2 M ω = 1 2 M b ω is the average value of the variable B .
The magnitude of ρ AB reflects the correlation degree between the analyzed sequences. The closer the correlation coefficient | ρ AB | is to 1, the closer A and B are to linear correlation, and the greater the correlation of A and B .
However, due to the influence of residual white noise, satellite visible conditions and the distance between the antenna of the monitoring station and the reflecting surface, the correlation coefficient of the smoothed DD coordinate residual sequences between adjacent days becomes smaller.

4.3. The Method for Determining the Best Decomposition Level of the WPT

Determining the best decomposition level of the WPT that is related to the validity of multipath extraction, if there are too many decomposition layers, the useful information in the smoothed DD coordinate residual sequences will be lost. If the number of decomposition layers is too small, the white noise in the smoothed DD coordinate residual sequences cannot be eliminated absolutely.
The method based on the correlation analysis for determining the best decomposition level of the WPT can be described as follows:
Step1: Calculating the cross-correlation coefficient of the smoothed DD coordinate residual sequences between the reference day and other adjacent days, for any reference day, and selecting the smoothed DD coordinate residual sequences that have the largest cross correlation coefficients.
Step2: Performing M-layers decomposition on the set of the smoothed DD coordinate residual sequences obtained in step 1 by using the WPT algorithm, and calculating the cross-correlation coefficient of the wavelet packet coefficients of the approximate part (low-frequency part) of the each decomposition layer between the smoothed DD coordinate residual sequences obtained in step 1.
Step3: Determining whether the target layer J is smaller than M, where M is the preset number of decomposition layers, and the target layer J is the number of layers in which the maximum value of the cross-correlation coefficient is located.
If yes, determining the target layer J be the best decomposition level of the WPT;
If not, the value of M is incremented by one, and the process returns to step 3 until the target layer J is smaller than M. By this way, we can figure out the best decomposition level of the WPT that can be applied to analyze the smoothed DD coordinate residual sequences in the optimal state.

4.4. The Principle of Wavelet Packet Threshold Denoising

We have the decomposed DD coordinate residual sequences with the best decomposition level of the WPT and the wavelet packet coefficients obtained should be denoised by thresholding for reconstructing the DD coordinate residual sequences.
The process of denoising by thresholding is carried out by comparing the magnitude of the wavelet packet coefficients d i with a threshold λ . The process can be described in two aspect: one is the choice of the threshold function F and the other is the choice of the threshold parameter λ .
There are three choices of threshold function F including hard thresholding F λ H ( d i ) , soft thresholding F λ S ( d i ) and quantitative thresholding F λ Q ( d i ) [28]:
F λ H ( d i ) = { 0 , | d i | < λ d i , | d i | λ ,
F λ S ( d i ) = sign ( d i ) ( | d i | λ ) ,     | d i | < λ ,
F λ Q ( d i ) = { 0 , | d i | < Ρ d i , | d i | Ρ ,
where d i denotes the wavelet coefficients, and Ρ is the value for which a certain percentage of coefficients d i is eliminated [17].
In this paper, the hard thresholding function is chosen for denoising the wavelet packet coefficients, which had been proved be the most suitable function for GPS applications [17].
As for the choice of the threshold parameter λ , Donoho (1995) had proposed the universal threshold [29]:
λ = 2 σ 2 log ( n ) ,
where n denotes the number of samples in white noise vector N and σ is the white noise level and needs to be estimated in each DD coordinate residual sequence. For determining the value of σ , Donoho and Johnstone (1994) proposed the following estimator [30]:
σ = Median [ | N | ] 0.6745 .
The noise vector N in the first level of decomposition layers is used for determining the value of σ .
Applying the thresholding method mentioned as above, all of the wavelet packet coefficients are denoised, which contribute to the reconstruction of the DD coordinate residual sequences and obtain the noiseless multipath.

5. The Overall Program

The overall program of mitigating the multipath in GPS static high-precision positioning is shown in Figure 2.

6. Verify the Feasibility of the Algorithm under Simulated Conditions

The DD coordinate residual sequences is simulated with the following model:
S ( t ) = cos ( 2 π t 1200 ) + sin ( 2 π t 900 ) + sin ( 2 π t 300 ) + e ( t ) ,
where S ( t ) consists of a cosine signal with a period of 1200 s and two sinusoidal signals with periods of 900 s and 300 s which represents typical GPS multipath wavelengths, and e ( t ) denotes a Gaussian white noise sequence.
Take the noise level of e ( t ) as N (0, 1.52). The data sampling rate is 1 s, and the sample size is 6000.
The simulated signal is denoised by WT transform, the TDMWA algorithm and WP-TD algorithm, respectively. DAUB8 is taken as the wavelet function, and we assume that the best decomposition level is 5.
The denoising result is shown in Figure 3, and the relative distribution of the signal denoised by the three different algorithms with the original signal is shown in Figure 4.
The residual signal generated by the three different algorithms is shown in Figure 5.
From Figure 4 and Figure 5, we can visually see that the signal denoised by the WP-TD algorithm better restores the original signal. The denoised result of the signal processed using the TDMWA algorithm is also slightly better than the traditional WT algorithm. In addition, the signal denoised using the WT algorithm has an edge effect, but the signal denoised using the TDMWA algorithm does not appear.
The signals denoised by the above three algorithms and the residual signals generated are analyzed quantitatively. The analysis indicators include the root mean square (RMS) value of the noise part of denoised signal N R M S , the RMS value of the signal part of the denoised signal S R M S and the correlation coefficient ρ between the denoised signal and the original signal. In order to scientifically evaluate the performance of the three algorithms, simulation analysis was carried out in four different Gaussian white noise simulation environments. The results are shown in Table 1, Table 2 and Table 3.
On one hand, it can be seen from Table 1, Table 2 and Table 3 that the correlation coefficients of the three algorithms mostly exceed 0.95 at different noise levels, indicating that all of them can restore the information of the original very well. However, the value of S R M S of the WP-TD algorithm is relatively optimal at different noise levels, and its correlation coefficient is greater than 0.98.
On the other hand, the value of N R M S of the three algorithms are basically equivalent while there is some difference in the value of S R M S . According to the tables, as the noise level increases, the value of S R M S of the WT algorithm and the TDMWA algorithm increase significantly compared with WP-TD algorithm. However, at high noise levels, the TDMWA algorithm increases slowly while the WT algorithm is relatively fast, and its S R M S is still superior to the WT algorithm.
In summary, we can draw such a conclusion that the WP-TD algorithm can denoise the noise signal more effectively, and can effectively weaken the influence of the edge effect in signal filtering.

7. Analyze the Performance of the Algorithm in a Measured Environment

In order to further verify the ability of the WP-TD algorithm to eliminate noise in the actual environment, we established a GPS data acquisition system in the laboratory and conducted experiments in a multipath environment.

7.1. Set up GPS Data Acquisition System

The block diagram of the overall scheme and key modules of the positioning data acquisition system are shown in Figure 6 and Figure 7, which shows the scene map and physical device diagram.
In this system, the base station and mobile station are all equipped with a DC-RTK-00 Beidou/GPS/Galileo three-mode single-frequency RTK module, a power model (5 V), a microcontroller (Raspberry3 b+) and a data transmission module (Ethernet card/4G card). The base station sends the differential data to the server over the wired network and is cached by the lab server. The mobile station accesses the server through the 4G network to obtain differential data, and, at the same time, the mobile station equipped with the 4G network card transmits the positioning result back to the lab server.

7.2. Data Collection and Processing

The fixed point is located at the top of the research building of Beijing University of Posts and Telecommunications. The settings of related parameters of data acquisition are shown in Table 4.
The DD coordinate residual sequences to be processed are part of the elevation data sequence measured of the fixed point for three consecutive days in a multipath environment.
Figure 8 shows the DD coordinate residual sequences of the elevation direction for three consecutive days. As can be seen from Figure 8, under the influence of the observation noise, the DD coordinate residual sequences for three consecutive days still have significant repeatability.
According to the overall program of mitigating the multipath in GPS static high-precision positioning established in this paper, we need to smooth the DD coordinate residual sequences first using the TDMWA algorithm.
Secondly, we need to analyze the correlation of the DD coordinate sequences residuals for three consecutive days and determine the best decomposition level of WPT based on the method proposed in this paper.
The process of determining the best decomposition level of WPT is shown in Figure 9.
From right to left, and from top to bottom are the DD coordinate residual sequence of the reference day preprocessed by the TDMWA algorithm, the seven-layer decomposition tree structure of the DD coordinate residual sequence of the reference day, the cross-correlation colored coefficients of the corresponding nodes and the value of the first wavelet packet coefficient of the decomposed DD coordinate residual sequence of the reference day located on the sixth layer. It can be seen from the distribution of the colored coefficients in the figure that the larger values of the correlation coefficient are mostly distributed in the sixth layer and the fifth layer, and the maximum value appears in the sixth layer.
According to the judgment conditions in the method proposed in this paper, we can determine that the best denoise result can be obtained when the decomposition level of the wavelet packet is six layers.
Next, we perform a six-layer wavelet packet decomposition on the smoothed DD coordinate residual sequences and reconstruct the decomposed DD coordinate residual sequences based on the wavelet packet coefficients obtained by threshold denoising.
Finally, we can obtain almost noiseless multipath, which can be applied to correct the position sequences of the adjacent days.
The DD coordinate residual sequences denoised by the WP-TD algorithm, Vondrak algorithm and WT algorithm for three consecutive days are shown in Figure 10, and Figure 11 shows the difference of denoised DD coordinate residual sequences for the second day and the third day compared to the first day, respectively, and Figure 12 shows the noise residuals denoised by the three algorithms on the first day.
It can be seen from the Figure 10 and Figure 11 that the difference of DD coordinate residual sequence processed by the WP-TD algorithm is much smaller than the WT algorithm and the Vondrak algorithm, which shows that applying the DD coordinate residual sequence denoised by the WP-TD algorithm to correct the position sequence of the adjacent day will be more precise.
Moreover, the DD coordinate residual sequence has a larger extremum at the end of sequence, which was denoised by the WT algorithm without the preprocessing of TDMWA algorithm, while the DD coordinate residual sequence presents a smooth curve at the end of the sequence which is preprocessed by the TDMWA algorithm before being denoised by the WPT algorithm.

7.3. Data Analysis

In order to further compare the performance for extracting multipath of the three algorithms in the actual environment, the RMS value of the DD coordinate residual sequences before and after filtering, the RMS value of the DD coordinate residual sequences of the adjacent days corrected by the almost noiseless multipath obtained in the reference day and the correlation coefficient of the DD coordinate residual sequence of the three days are used as the evaluation indicator.
The comparison results of the three algorithms are shown in Table 5, Table 6 and Table 7, respectively.
First of all, it can be seen from the results in Table 5 that the RMS value of the DD coordinate residual sequence after filtering by the WP-TD algorithm is basically equivalent to the RMS value of the DD coordinate residual sequence after filtering by the Vondrak algorithm, but both are slightly better than the RMS value of the DD coordinate residual sequence after filtering by the WT algorithm. In addition, the noise residuals denoised by the three algorithms on the first day are not much different as shown in Figure 12, indicating that the observation noise (white noise) is only a small part of the DD coordinate residuals sequence, and multipath dominates.
Secondly, the results in Table 6 show that the GPS multipath of the fixed point is highly repetitive, and the correlation coefficient increases after filtering observation noise. Moreover, the performance of filtering of the WP-TD algorithm is better than WT algorithms and Vondrak algorithms more or less. Through calculation, we can know that the correlation of the DD coordinate residual sequence of the three days filtered by WP-TD is 3.02% and 1.78% higher than that of WT algorithms and Vondrak algorithms, respectively.
Last but not least, according to the data in Table 7, after correcting by the denoised multipath of the reference day (Day 1) which was extracted by the WT algorithms, the Vondrak algorithms, and the WP-TD algorithm, we can figure out that the accuracy of the position sequence of the adjacent days (Day 2 and Day 3) is increased by 57.47%, 65.98%, and 71.61%, respectively, compared with before correction, which further verified that the WP-TD algorithm can obtain an more accurate multipath correction model.

8. Conclusions

The WP-TD algorithm is proposed, in which the method of determining the optimal decomposition level of WPT is given, which can avoid the redundancy and loss of information after signal reconstruction to the greatest extent.
The WP-TD algorithm is proposed, which makes use of the TDMWA algorithm, and is proven to have good performance that weakens the influence of the edge effect through simulation analysis and measured data analysis.
The WP-TD algorithm is proposed, which uses the WPT to divide the signal frequency band into multiple layers, further decomposes the high frequency part without subdivision in the traditional WT algorithm, and eliminates the observation noise in the high frequency part.
The WP-TD algorithm is proposed, which applies the denoised multipath of the reference day to correct the position sequence of the adjacent days, improving positioning accuracy by 14.4%, compared to the traditional WT algorithm.
In summary, the WP-TD algorithm is proposed, which can perform well for eliminating the multipath in GPS static high-precision positioning.

Author Contributions

Conceptualization, K.H. and C.T.; Methodology, K.H. and C.T.; Data collection and analysis, C.T., K.H. and Z.D.; Investigation, C.T.; Writing-Original Draft Preparation, C.T.; Writing-Review & Editing, C.T., K.H. and Z.D.; Supervision, Z.D. and K.H.; Project Administration, Z.D. and K.H.

Funding

This research was funded by [The authors acknowledge the National Key Research and Development Program of China] grant number [2016YFB0502001].

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Ji, X.; Lu, C. Comparative study on GPS multipath error filtering method. Bull. Surv. Mapp. 2015, 4, 10–13. [Google Scholar]
  2. Zhou, D. Research on Multipath Error Based on Wavelet Transform; Jiangxi University of Technology: Nanchang, China, 2010. [Google Scholar]
  3. Zhong, P.; Ding, X.L.; Zheng, D.W.; Chen, W.; Huang, D.F. Adaptive wavelet transform based on cross-validation method and its application to GPS multipath mitigation. GPS Solut. 2008, 12, 109–117. [Google Scholar] [CrossRef]
  4. Pozoruz, A.; Martinez, J.L.; Garciacerezo, A. A new satellite selection criterion for DGPS using two low-cost receivers. In Proceedings of the IEEE International Conference on Robotics & Automation, Leuven, Belgium, 20 May 1998. [Google Scholar]
  5. Scire-Scappuzzo, F.; Makarov, S.N. A low-multipath wideband GPS antenna with cutoff or non-cutoff corrugated ground plane. IEEE Trans. Antennas Propag. 2009, 57, 33–46. [Google Scholar] [CrossRef]
  6. Ray, J.K.; Cannon, M.E.; Fenton, P. Mitigation of static carrierphase multipath effects using multiple closely spaced antennas. Navig 1999, 46, 193–202. [Google Scholar] [CrossRef]
  7. Townsend, B.R.; Fenton, P.C. A practical approach to the reduction of pseudorange multipath errors in a Ll GPS receiver. In Proceedings of the ION GPS-94, U.S. Inst. of Navigation, Salt Lake City, UT, USA, 20–23 September 1994; pp. 143–148. [Google Scholar]
  8. Nee, R.D.J.V. The multipath estimating delay lock loop. In Proceedings of the IEEE Second International Symposium on Spread Spectrum Techniques and Applications, ISSTA 92, Yokohama, Japan, 29 November–2 December 1992; IEEE: Piscataway, NJ, USA, 1992. [Google Scholar]
  9. Comp, C.; Axelrad, P. Adaptive SNR-based carrier phase multipath mitigation technique. IEEE Trans. Aerosp. Electr. Syst. 1998, 34, 264–276. [Google Scholar] [CrossRef]
  10. Ge, L.; Han, S.; Rizos, C. Multipath mitigation of continuous GPS measurements using an adaptive filter. GPS Solut. 2000, 4, 19–30. [Google Scholar] [CrossRef]
  11. Nce, C.D.; Sahin, M. Real-time deformation monitoring with GPS & Kalman filter. Earth Planets Space 2000, 52, 837–840. [Google Scholar]
  12. Han, S.; Rizos, C. Multipath effects on GPS in mine environments. In Proceedings of the Xth International Congress of the International Society for Mine Surveying, Fremantle, Australia, 2–6 November 1997; pp. 447–457. [Google Scholar]
  13. Vondrák, J. Problem of smoothing observational data II. Bull. Astron. Inst. Czechoslov. 1977, 28, 84. [Google Scholar]
  14. Mallat, S.G. A theory for multiresolution signal decomposition: The wavelet representation. IEEE Trans. Pattern Anal. Mach. Intell. 1989, 11, 674–693. [Google Scholar] [CrossRef]
  15. Huang, S.; Liu, J.; Liu, X. Deformation analysis based on wavelet and its applicationin dynamic monitoring for high-rise buildings. Acta Geod. Cartogr. Sin. 2003, 32, 153–157. [Google Scholar]
  16. Zhang, Y.; Bartone, C. Real-time multipath mitigation with WaveSmoothTM technique using wavelets. In Proceedings of the ION GNSS 2004, Long Beach, CA, USA, 21–24 September 2004; pp. 1181–1194. [Google Scholar]
  17. Souza, E.M.; Monico, J.F.G. Wavelet Shrinkage: High frequency multipath reduction from GPS relative positioning. GPS Solut. 2004, 8, 152–159. [Google Scholar] [CrossRef]
  18. Azarbad, M.R.; Mosavi, M.R. A new method to mitigate multipath error in single-frequency GPS receiver with wavelet transform. GPS Solut. 2014, 18, 189–198. [Google Scholar] [CrossRef]
  19. Lau, L. Wavelet packets based denoising method for measurement domain repeat-time multipath filtering in GPS static high-precision positioning. GPS Solut. 2017, 21, 461. [Google Scholar] [CrossRef]
  20. Souza, E.M.; Negri, T.T. First prospects in a new approach for structure monitoring from GPS multipath effect and wavelet spectrum. Adv. Space Res. 2017, 59, 2536–2547. [Google Scholar] [CrossRef]
  21. Chen, G.; Li, Q.-Y.; Li, D.-Q.; Wu, Z.-Y.; Liu, Y. Main frequency band of blast vibration signal based on wavelet packet transform. Appl. Math. Model. 2019, 74, 569–585. [Google Scholar] [CrossRef]
  22. Coifman, R.R.; Wickerhauser, M.V. Entropy-based algorithms for best basis selection. IEEE Trans. Inf. Theory 1992, 38, 713–718. [Google Scholar] [CrossRef]
  23. Alessio, E.; Carbone, A.; Castelli, G.; Frappietro, V. Second-order moving average and scaling of stochastic time series. Eur. Phys. J. B Condens. Matter 2002, 27, 197–200. [Google Scholar] [CrossRef]
  24. Khelifa, S.; Kahlouche, S.; Belbachir, M.F. Application of wavelet analysis to GPS stations coordinate time series. In Proceedings of the FIG Working Week 2011, Bridging the Gap Between Cultures, Marrakech, Morocco, 18–22 May 2011. [Google Scholar]
  25. Ye, X.; Hu, Y.; Ao, M. GPS multipath mitigation based on wavelet packet decomposition. Comput. Inf. Syst. 2013, 9, 477–484. [Google Scholar]
  26. Strang, G.; Borre, K. Linear Algebra, Geodesy, and GPS; Wellesley-Cambridge Press: Wellesley, MA, USA, 1997; ISBN 0961408863. [Google Scholar]
  27. Elad, M. On the origin of the bilateral filter and ways to improve it. IEEE Trans. Image Process. 2002, 11, 1141–1151. [Google Scholar] [CrossRef]
  28. Percival, D.B.; Walden, A.T. Wavelet Method for Time Series Analysis; Cambridge University Press: Cambridge, UK, 2000; pp. 398–440. [Google Scholar]
  29. Donoho, D.L. De-noising by soft-thresholding. IEEE Trans. Inf. Theory 1995, 41, 613–627. [Google Scholar] [CrossRef]
  30. Donoho, D.L.; Johnstone, J.M. Ideal spatial adaptation by wavelet shrinkage. Biometrika 1994, 81, 425–455. [Google Scholar] [CrossRef]
Figure 1. Wavelet packet decomposition tree.
Figure 1. Wavelet packet decomposition tree.
Sensors 19 02704 g001
Figure 2. The overall program of mitigating the multipath in GPS static high-precision positioning.
Figure 2. The overall program of mitigating the multipath in GPS static high-precision positioning.
Sensors 19 02704 g002
Figure 3. The signals denoised by the three different algorithms.
Figure 3. The signals denoised by the three different algorithms.
Sensors 19 02704 g003
Figure 4. The relative distribution of the signal denoised by the three different algorithms with the original signal.
Figure 4. The relative distribution of the signal denoised by the three different algorithms with the original signal.
Sensors 19 02704 g004
Figure 5. The residual signals generated by the three different algorithms.
Figure 5. The residual signals generated by the three different algorithms.
Sensors 19 02704 g005
Figure 6. Data acquisition system.
Figure 6. Data acquisition system.
Sensors 19 02704 g006aSensors 19 02704 g006b
Figure 7. Scene map and physical device diagram.
Figure 7. Scene map and physical device diagram.
Sensors 19 02704 g007
Figure 8. Original DD coordinate residual sequences of the elevation direction for the three consecutive days.
Figure 8. Original DD coordinate residual sequences of the elevation direction for the three consecutive days.
Sensors 19 02704 g008
Figure 9. The process of determining the best decomposition level of WPT.
Figure 9. The process of determining the best decomposition level of WPT.
Sensors 19 02704 g009
Figure 10. The DD coordinate residual sequences denoised by the three algorithms for three consecutive days.
Figure 10. The DD coordinate residual sequences denoised by the three algorithms for three consecutive days.
Sensors 19 02704 g010
Figure 11. The difference of DD coordinate residual sequences denoised by the three algorithms.
Figure 11. The difference of DD coordinate residual sequences denoised by the three algorithms.
Sensors 19 02704 g011
Figure 12. Noise residuals denoised by the three algorithms on the first day.
Figure 12. Noise residuals denoised by the three algorithms on the first day.
Sensors 19 02704 g012
Table 1. Performance indicators of the WT algorithm in different Gaussian white noise simulation environments.
Table 1. Performance indicators of the WT algorithm in different Gaussian white noise simulation environments.
Noise Level0.51.52.53.5
S R M S 0.13820.26480.32320.3631
N R M S 0.49121.45142.46232.9827
ρ 0.98860.96390.96240.9593
Table 2. Performance indicators of the TDMWA algorithm in different Gaussian white noise simulation environments.
Table 2. Performance indicators of the TDMWA algorithm in different Gaussian white noise simulation environments.
Noise Level0.51.52.53.5
S R M S 0.12680.22310.26380.2985
N R M S 0.48891.44852.468842.9816
ρ 0.98920.98020.97420.9698
Table 3. Indicators of the WP-TD algorithm in different Gaussian white noise simulation environments.
Table 3. Indicators of the WP-TD algorithm in different Gaussian white noise simulation environments.
Noise Level0.51.52.53.5
S R M S 0.08870.14010.20220.2436
N R M S 0.47221.43532.468412.9707
ρ 0.99980.98820.98320.9801
Table 4. The setting of data acquisition parameters.
Table 4. The setting of data acquisition parameters.
Main ParametersValue
Receive signal typeGPS L1
Sampling frequency (HZ)1
Baud rate (bps)115,200
Elevation mask angle (deg)20.000
SNR mask (dB)30.000
Baseline length (km)<3
Sampling time9:00 p.m.–10:30 p.m. (1–3 November)
Total epoch5400
Table 5. RMS value of the DD coordinate residual sequences before and after filtering by the three algorithms (mm).
Table 5. RMS value of the DD coordinate residual sequences before and after filtering by the three algorithms (mm).
TimeBefore FilteringWTVondrakWP-TD
D116.442216.322416.008815.8628
D217.389517.334716.765516.6834
D316.987216.782516.445616.2218
Table 6. Correlation coefficient of the DD coordinate residual sequence of the three days after being denoised by the three algorithms.
Table 6. Correlation coefficient of the DD coordinate residual sequence of the three days after being denoised by the three algorithms.
TimeBefore FilteringWTVondrakWP-TD
D1-D20.93880.95860.97220.9873
D2-D30.92360.95210.96360.9822
D3-D10.92480.95490.96470.9826
Table 7. The RMS value of the DD coordinate residual sequences of the adjacent days corrected by the almost noiseless multipath obtained in the reference day (mm).
Table 7. The RMS value of the DD coordinate residual sequences of the adjacent days corrected by the almost noiseless multipath obtained in the reference day (mm).
TimeBefore FilteringWTVondrakWP-TD
D215.34766.54525.45434.5343
D316.82487.13475.46684.5839

Share and Cite

MDPI and ACS Style

Han, K.; Tang, C.; Deng, Z. A New Method for Multipath Filtering in GPS Static High-Precision Positioning. Sensors 2019, 19, 2704. https://doi.org/10.3390/s19122704

AMA Style

Han K, Tang C, Deng Z. A New Method for Multipath Filtering in GPS Static High-Precision Positioning. Sensors. 2019; 19(12):2704. https://doi.org/10.3390/s19122704

Chicago/Turabian Style

Han, Ke, Canyang Tang, and Zhongliang Deng. 2019. "A New Method for Multipath Filtering in GPS Static High-Precision Positioning" Sensors 19, no. 12: 2704. https://doi.org/10.3390/s19122704

APA Style

Han, K., Tang, C., & Deng, Z. (2019). A New Method for Multipath Filtering in GPS Static High-Precision Positioning. Sensors, 19(12), 2704. https://doi.org/10.3390/s19122704

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop