Multipath
Multipath error typically occurs near the receiver and is caused by the refraction or reflection of GNSS signals by various obstacles such as metal structures, glass, trees, and the ground surface. It affects both code and phase measurements, albeit in different ways. The extent of this error depends on the distance between the receiver and the satellite and varies at different frequencies. For instance, in code measurements, multipath error can reach up to 1.5 times the wavelength, while in phase measurements, the effect can be as small as a quarter of the wavelength (Lopez, 2011). The impact of multipath is observed differently during position determination, where the C/A code observations are accurate up to several meters, while the carrier wave phase observations provide accuracy up to several centimeters (Georgiadou and Kleusberg, 1988).
In dense urban areas, multipath manifests in three distinct types (Kubo et al., 2020):
-
Complete signal blockage occurs when the signal path is entirely obstructed and cannot be received by the receiver. This situation reduces the visible satellite geometry for the receiver, sometimes making it impossible to determine the position (Fig. 1a).
-
Indirect signals are received when an obstacle blocks the direct signal path, causing the signal to reach the receiver with a delay or distortion. This type of multipath is known as an indirect signal (Fig. 1b).
-
Both the direct and indirect signals are received by the receiver. This type of multipath can be mitigated through changes in receiver design or signal processing techniques (Fig. 1c).
In this article, the focus is on investigating the multipath effect at the level of raw measurements. Various methods exist to address multipath at different stages of measurement processing, but the specific aim of this study is to examine its impact on raw measurements (Kubo et al., 2020).
Code Minus Carrier-phase observable
To study the multipath effect, the carrier phase minus code (CMCD) observation technique is employed, which involves subtracting the time difference between code and phase observations (both measured in meters). This approach eliminates the need to resolve phase ambiguity for continuously sampled observations without cycle slips (Beitler et al., 2015; Caamano et al., 2020; Pirsiavash et al., 2018).
$$\begin{gathered} CMCD_{r}^{s}({t_r},{t_{r+1}})=dR_{r}^{s}({t_r},{t_{r+1}}) - \lambda d\Phi _{r}^{d}({t_r},{t_{r+1}}) \hfill \\ \hfill \\ \end{gathered}$$
1
Where \(dR_{r}^{s}({t_r},{t_{r+1}})\)is the code observation change rate (meters), \(\lambda\)is wavelength (meters) and \(d\Phi _{r}^{d}({t_r},{t_{r+1}})\) is the phase observation change rate per cycle unit. Simplifying Eq. (1) leads to:
$$CMCD_{r}^{s}({t_r},{t_{r+1}})=2{\dot {\delta }_{ion}}+{\dot {\delta }_{mu{l_R}}} - \lambda \dot {N}_{r}^{s}+{\varepsilon _R}$$
2
Where \({\dot {\delta }_{ion}}\) is ionosphere changes in meters in two consecutive epochs, \({\dot {\delta }_{mu{l_R}}}\) is multipath effect changes in meters in two consecutive epochs, \(\dot {N}_{r}^{s}\) is cycle slip in cycles and \({\varepsilon _R}\) are observations noise in meters. Although the ionosphere changes in the observation sampling intervals can be assumed equal to zero (Beitler et al., 2015, Pirsiavash et al., 2018), but due to higher-order effects of the ionosphere and the purpose of research to investigate the CMCD index solely, this effect is calculated by implementing ionosphere global model. Therefore, if the observations are continuous and without any cycle slips, CMCD will be equal to the multipath changes (Kim et al., 2016).
Detecting and removing cycle slips
In this research, the Doppler integration method has been utilized to detect and eliminate cycle slips. This method takes advantage of the fact that Doppler measurements are typically free from cycle slips (Karaim et al., 2014). Correcting phase observations for cycle slips generally involves two steps: first, detecting the occurrence of the cycle slip in time, and second, calculating and applying the corresponding correction to the observations (Karaim et al., 2014). Doppler measurements, which capture instantaneous frequency changes, serve as a suitable indicator for detecting and correcting cycle slips.
If the Doppler equation is expressed discretely (Kim et al., 2016), below equation estimates the computational phase at moment k + 1:
$${\hat {\Phi }_{k+1}}={\Phi _k}+\frac{{{{\dot {\Phi }}_{k+1}}+{{\dot {\Phi }}_k}}}{2}\Delta t$$
3
Where \({\dot {\Phi }_k}\) and \({\dot {\Phi }_{k+1}}\)are the rate of phase observations in cycles at k and k + 1, \({\Phi _k}\) is the value of the observed phase at moment k and\({\hat {\Phi }_{k+1}}\) is the value of the calculated phase for moment k + 1.
By comparing the calculated and observed phase at k + 1, it is possible to estimate the amount of cycle slip at k + 1 (Fujita et al., 2013, Ji et al., 2018, Karaim et al., 2014, Zhao et al., 2020):
$${\dot {N}_{k+1}}={\hat {\Phi }_{k+1}} - {\Phi _{k+1}}$$
4
If \(\dot {N}\) (cycle slip) exceeds the defined threshold, it is possible to detect the occurrence of a cycle slip equal to \(\dot {N}\); otherwise, the period of measurement under investigation is free of any phase jump (Karaim et al., 2014). The cycle slip detection threshold in this research considered \(\lambda\) (the wavelength of each signal for the respective satellite system).
Multipath Detection
By incorporating ionosphere variations and the magnitude of cycle slip into Eq. 2, the CMCD index is able to estimate the multipath rate and measurement noise within the specific time interval.\(CMCD_{r}^{s}({t_r},{t_{r+1}})={\dot {\delta }_{mu{l_R}}}+{\varepsilon _R}\) (5)
Equation 5 introduces the noise parameter, which is an unknown and requires a threshold to differentiate between multipath and observation noise (Beitler et al., 2015; Pirsiavash et al., 2018). To address this, the method employed by Caamano et al. (2020) involves grouping the estimated CMCD values for satellites based on their elevation angles. Then, Eq. 6 is utilized to calculate the standard deviation within each elevation bin (Caamano et al., 2020).
$${\sigma _{CMC{D^q}}}=\sqrt {\frac{{\sum {{{(CMC{D_i} - \overline {{CMCD}} )}^2}} }}{N}}$$
6
where \({\sigma _{CMC{D^q}}}\) is the standard deviation of the CMCD measurement in the q height bin, \(CMC{D_i}\) is desired CMCD value (meters), \(\overline {{CMCD}}\) is average of CMCDs in the desired height bin (meters), and is number of samples in the desired height category. If Eq. 7 holds true, the value of CMCD is accepted as the multipath rate in that time period.
$$\left| {CMCD_{r}^{s}} \right| \geqslant \kappa {\sigma _{CMC{D^q}}}$$
7
where \(\kappa\) is an experimental coefficient and is chosen based on environmental conditions and measurements (Beitler et al., 2015, Pirsiavash et al., 2018).
Datasets and processing
In this research we tried first to increase the number of satellite observations through multi-system measurements and to improve the satellites' geometry, to provide the possibility of determining the precise position of smartphone even in challenging conditions. In the next step, the main goal was to detect, calculate and remove the multipath effect. The later goal was based on the GNSS signals received by the smartphone receiver and checking the effect of this error on the final position of the receiver.
However, due to the impossibility of solving the phase ambiguity in smartphone GNSS measurements, it was impossible to use the code and phase differential method. therefore, to estimate the multipath effect CMCD algorithm was considered.
On the other hand, in the observations of smartphones, there are many phase interruptions and sudden jumps in the carrier phase tracking that took away the possibility of correct estimation of CMCD; therefore, another goal of this research was to find cycle slips and remove them from the measurements.
Data collecting
In this article, the data utilized for the research was obtained from the Android GNSS-Logger software, which was recorded during the "Google Smartphone Decimeter Challenge" competition hosted by Google and Kaggle. The data was collected in the highways of the San Francisco Bay area and comprised a total of 2175 epochs with a sampling rate of 1 second (Kaggle, 2020). The GNSS-Logger software produces a RINEX 3.0.3 formatted file containing various measurements such as pseudo-ranges, accumulated delta-range, Doppler frequencies, and CN0 values for GPS, GLONASS, Galileo, BeiDou Navigation Satellite System (BDS), and QZSS satellites. These measurements are obtained across different frequencies, including L1, L5, E1B, E1C, and E5A (Google, 2022). The observations were acquired using a Pixel 4 smartphone, which began recording on September 4, 2020, at 17:59:26.4297477. During the data storage process, the smartphone was placed on the car's dashboard, while a geodetic antenna located on the car's roof was used for data validation.
In multipath investigation, it is crucial to assess the magnitude of the error in each measurement epoch. The calculated multipath rate is related to the pseudo-range rate, enabling the correction of pseudo-range values in each epoch by subtracting the CMCD (Carrier Minus Code) from the pseudo-range rate. By comparing the measured pseudo-range with the calculated pseudo-range (after CMCD removal) in each epoch, it becomes possible to determine the amount of multipath error on the code observations at that moment. This comparison allows for an evaluation of the consistency between the observed and calculated pseudo-range values in each epoch. the amount of multipath error on the code observations at each moment is obtained using equations below:
$$\begin{gathered} d\hat {R}_{r}^{s}({t_r})=dR_{r}^{s}({t_r},{t_{r+1}}) - {{\dot {\delta }}_{mu{l_R}}} \hfill \\ =dR_{r}^{s}({t_r},{t_{r+1}}) - CMCD \hfill \\ \end{gathered}$$
8
$$\hat {R}_{r}^{s}({t_{r+1}})=R_{r}^{s}({t_r})+d\hat {R}_{r}^{s}$$
9
$${\delta _{mu{l_R}}}=\hat {R}_{r}^{s}({t_{r+1}}) - R_{r}^{s}({t_{r+1}})$$
10
As mentioned earlier, there are various methods to deal with the multipath effect, from detecting and reducing the weight of multipath observations in the positioning process to calculating and removing its effect from the measurements and then solving the equations by the desired method. In this research, we tried to estimate This effect using the CMCD index and to remove its effects from the measurements. The positioning method used in this research was Single Point Positioning (SPP) using correction files. For this obtained the best possible positioning accuracy by using the corrected code in the GPS, GLONASS, Galileo, and BDS. Since the estimation of the receiver position depends on the position of the satellites, the precise positions of the satellites were calculated using precise orbit products from International GNSS Service (IGS).
Data processing
The raw GNSS measurements obtained from smartphones require filtering and pre-processing to transform them into observations comparable to those from a geodetic receiver. Figures 3 to 6 in the article illustrate the observations of code rate, phase rate, and the removal of cycle slips for four tracked satellites throughout the experiment. These Figures provide a visual representation of the variations in code rate and phase rate over time, as well as the detection and correction of cycle slips in the measurements. By analyzing these Figures, researchers can gain insights into the behavior of the observed satellite signals and the effectiveness of the applied processing techniques.
To remove the jumps (cycle slips) in the phase observations, the calculated values for the cycle slips are added to the phase observation rate. By incorporating these corrected values, the cycle slips are effectively eliminated from the phase observations along the entire path of the receiver's movement. This process helps ensure the accuracy and continuity of the phase observations, allowing for more reliable analysis and positioning calculations.
The carrier-phase observations after removed cycle-slip are displayed in Fig. 7 to Fig. 10, after removing cycle slips, the phase rate observations free of cycle slips and were entered into CMCD algorithm.
The estimated values for the CMCD observables (Fig. 11 to Fig. 14), which were calculated by inserting the code rate and phase rate observations at each instant into Eq. 1, were applied to Eq. 10 to obtain the multipath error rate
Threshold values according to the satellite height category and acceptable values as multipath rate for GNSS satellites based on their elevation angles are shown in Fig.s 11 to 14. In this research, (in Eq. 7) has been selected as 1 according to the remaining observations of the code (\(1\sigma\)). The estimated CMCD of the satellites at each epoch are depicted in Figs. 15 to 18:
Multipath corrected code
In Eq. 10, the multipath error is quantified by comparing the pseudo-range measured by the receiver with the pseudo-range estimated using the accepted CMCD measurements. The difference between these two values represents the impact of multipath error, which has been effectively removed through the CMCD algorithm. By applying this process, a pseudo-range free from multipath error is calculated. This corrected pseudo-range, which is not influenced by multipath effects, is then used in the positioning process. By accounting for and correcting the multipath effect, the accuracy and reliability of the positioning results can be improved.
Positioning results
The investigations conducted on determining the position using only GPS code observations revealed that there were epochs where the receiver received a limited number of observations, resulting in the inability to accurately estimate the position. This issue is evident in the receiver's movement path map shown in Fig. 19. To address this problem, multi-GNSS positioning was employed, which involved utilizing data from all available satellites, including GPS, GLONASS, Galileo, and Beidou. By incorporating observations from multiple GNSS systems, the number of available observations was increased, enhancing the accuracy and reliability of the position estimation process. This approach helps overcome the limitations associated with a single GNSS system and improves the overall positioning performance.
Positioning by correcting all errors except Multipath
By incorporating multi-GNSS measurements and employing the least squares technique, the receiver is able to estimate the positions of epochs that could not be determined with GPS-only positioning. Additionally, by considering satellites with observations having residuals of less than 5 meters, more reliable positioning results are obtained. Table 1 provides a comparison of the estimated position differences between GPS-only positioning and multi-GNSS positioning for epochs 350 to 360 (18:05:15.4297493 to 18:05:25.4297493). This time period corresponds to when the car is moving through underpasses, where direct GNSS signals are not received. The table highlights the improvement in positioning accuracy achieved by utilizing multi-GNSS measurements, as it allows for better estimation of positions in challenging environments where direct signals are obstructed.
Table 1
Position estimation error in different time epochs
Time epoch | GPS estimation error (m) | GNSS estimation error without multipath correction (m) | GNSS estimation error with multipath correction (m) |
| dE | dN | dE | dN | dE | dN |
350 | -1.6369 | -4.2028 | -3.428 | -2.5131 | -3.317 | 1.418 |
351 | -6.4109 | -0.3964 | -4.3004 | 1.71973 | -4.4953 | -0.207 |
352 | -5.3875 | -2.8614 | -4.9326 | 1.06864 | -6.5319 | 0.407 |
353 | -6.9505 | -7.7602 | -5.2078 | 1.43748 | -6.984 | 2.835 |
354 | NaN | NaN | -40.309 | -7.0156 | -10.736 | 4.030 |
355 | NaN | NaN | -48.542 | -16.324 | -16.3 | 3.564 |
356 | NaN | NaN | -64.731 | -31.643 | -21.407 | 1.792 |
357 | NaN | NaN | -47.061 | -7.9282 | -23.179 | 0.419 |
358 | NaN | NaN | -39.167 | 7.80952 | -19.622 | 0.451 |
359 | -4.509 | -4.868 | -5.712 | -0.865 | -5.362 | -0.668 |
360 | 2.325 | 1.5924 | 10.355 | 3.915 | 3.204 | 1.558 |
Positioning by correcting all errors for multi-GNSS satellites signal
The results demonstrate that correcting the multipath error in GNSS signals and employing the weighted least-squares algorithm for position estimation led to significant improvements in the horizontal position differences. Figure 24 provides a visual representation of the distance differences between the estimated positions obtained using three different determination methods (including multipath correction) and the reference position obtained from the GNSS receiver. This Fig. illustrates the improvements in position accuracy during the period from epoch 350 to 360.
Validation of the estimated positions
Figure 21 displays the estimated positions obtained from the smartphone receiver using three different positioning methods: GPS-only, multi-GNSS without multipath correction, and multi-GNSS with multipath correction. These positions are compared with the positions obtained from the reference NovAtel SPAN GNSS receiver mounted on the roof-top. The comparison results are presented in Table 1, which provides a quantitative assessment of the differences between the estimated positions from the smartphone receiver and the reference positions. Additionally, Fig. 21 visually represents the variations in position estimates among the three different positioning methods, highlighting the impact of multipath correction on improving the accuracy of the smartphone receiver's estimated positions.
In Table 2, the RMSE values are presented for the results of positioning without multipath correction and positioning with multipath correction.
Table 2
Comparison of RMSE values in positioning with two with and without multipath correction scenarios
Positioning Mode | RMSE (m) |
| E | N |
GNSS without multipath correction | 13.36576 | 8.455281 |
GNSS with multipath correction | 11.49932988 | 7.975923989 |