2. FUNDAMENTALS OF SF-PPP AND SF-RTK

In this section, the fundamentals of SF-PPP and SF-RTK are reviewed in Section 2.1 and Section 2.2 respectively.

2.1. SF-PPP

PPP utilises precise satellite products and calibrates the error sources to achieve centimetre-level to decimetre-level positioning accuracy with a single receiver. Since single-frequency low-cost GNSS receivers dominate the market, SF-PPP has been investigated in recent years (de Bakker and Tiberius, Reference de Bakker and Tiberius2017; Nie et al., Reference Nie, Zhou, Liu, Wang and Gao2019; Nie et al., Reference Nie, Liu and Gao2020). The single-frequency PPP positioning model can be given as

(1)$$\begin{array}{l} P_1 =\rho +c\cdot dt_r -c\cdot dt^s+T+I+B_{pr} -B^{ps}+\varepsilon _P \\ L_1 =\rho +c\cdot dt_r -c\cdot dt^s+T-I+\lambda _1 \cdot N_1 +b_{lr} -b^{ls}+\varepsilon _L \\ \end{array} $$

where P ₁ and L ₁ are the code and carrier-phase measurements of L1 (unit: metre); ρ is the geometric distance between the satellite and receiver; c is the light speed; dt _r represents the receiver clock offset; dt ^s is the satellite clock offset; T is the tropospheric delay; I is the ionospheric effect on L1; λ₁ and N ₁ represent the wavelength and the integer ambiguity; B _pr and B ^ps are the receiver and satellite code biases; similarly, b _lr and b _ls are the phase biases at the receiver- and satellite-end; $\varepsilon_{P}$ and $\varepsilon_{L}$ are the corresponding measurement noise.

In Equation (1), the satellite orbit and clock errors can be corrected with precise satellite products. Since most precise satellite products are referred to dual-frequency ionosphere-free code measurements, the code and phase biases on the right side of Equation (1) should be considered carefully. The satellite code bias can be corrected with the application of C1-P1 and P1-P2 differential code biases (DCB) products. Considering that the receiver code bias cannot be separated from the receiver clock offset, they are estimated as one unknown in the estimator. The phase biases are absorbed into ambiguities in the estimator. The tropospheric delay can be corrected by the Saastamoninen model (Saastamoinen, Reference Saastamoinen, Henriksen, Mancini and Chovitz1972) while the ionospheric delay can be relieved with the application of ionospheric products (e.g. IGS Global Ionospheric Map, CNES, VTEC, SSR etc.) (Schaer, Reference Schaer1997; Schaer et al,, Reference Schaer, Gurtner and Feltens1998; CNES, 2019). In addition to the errors listed in Equation (1), other error sources (e.g., the Sagnac effect, special relativistic effect, Shapiro time delay, satellite/receiver antenna phase centre offset and phase centre variation, phase wind-up effect etc.) also have to be calibrated before estimation (Wu et al., 1992; Ashby, Reference Ashby2003; Schmid et al., Reference Schmid, Steigenberger, Gendt, Ge and Rothacher2007). In other words, only the position and clock errors together with ambiguities are estimated. With the methods applied in this work, the positioning accuracy can reach within a few decimetres with a few seconds of convergence time. Due to the existence of ionosphere residuals, it is difficult to improve the performance without external aids (e.g. map matching).

2.2. SF-RTK

RTK applies the double difference between the base and rover observations to remove the common errors. In this work, a short baseline is applied and the atmospheric delays are assumed to be removed after double difference. The short-baseline SF-RTK positioning model between the base (a) and rover (b) for GNSS satellite i and j can be given as (Tian et al., Reference Tian, Ge and Neitzel2015)

(2)$$\begin{array}{l} P_{ab}^{ij} =\rho _{ab}^{ij} +\delta _{ab}^{ij} +\xi _{ab}^{ij} \\ \lambda ^i \varphi _{ab}^i -\lambda ^j \varphi _{ab}^j =\rho _{ab}^{ij} +\lambda ^iN_{ab}^i -\lambda ^jN_{ab}^j +\gamma _{ab}^{ij} +\varepsilon_{ab}^{ij} \\ \end{array}$$

where $P_{ij}^{ab}$ (unit: metre) and $\phi_{ij}^{ab}$ (unit: cycle) are the double-difference code and phase measurements; $\rho_{ij}^{ab}$ is the double-difference geometric distance; $\varphi_{i}^{ab}$ and $\varphi_{j}^{ab}$ are the single-difference between-station phase measurements of satellite i and j respectively; λⁱ and λ^j are the wavelengths of satellite i and j respectively; $N_{i}^{ab}$ and $N_{j}^{ab}$ are the single-difference ambiguities of corresponding satellites; $\delta_{ij}^{ab}$ and $\gamma_{ij}^{ab}$ are the inter-frequency biases (IFB) in code and phase measurements; $\xi_{ij}^{ab}$ and $\varepsilon_{ij}^{ab}$ are the noise of double-difference measurements.

In this paper, GPS measurements are applied. Therefore, $\delta_{ij}^{ab}$ and $\gamma_{ij}^{ab}$ are eliminated after double difference. Moreover, the single-difference ambiguities can be added together as double-difference ambiguities due to the same wavelength. After removing all the error sources in RTK, the position unknowns and the double-difference ambiguities are estimated in the estimator. The integer ambiguities in Equation (2) can be fixed and the positioning accuracy can be greatly improved. For dual-frequency RTK (DF-RTK), it is much easier to fix the ambiguities due to redundant measurements. However, it is challenging for SF-RTK to achieve the instantaneous fixed solution in some circumstances especially with low-cost GNSS receivers.

3. ROAD LANE DETECTION AND MAP MATCHING

3.1. Road lane detection

Road lane detection and recognition is an important part of the lane departure warning system in a modern automobile including an autonomous car. It also becomes increasingly significant in precise positioning, especially with a precise digital map. The most commonly applied sensor for lane detection is a camera, and the methods based on images can be categorised as mark-feature-based method and machine-learning-based method (Neven et al., Reference Neven2018; Pan et al., Reference Pan, Shi, Luo, Wang and Tang2018; Zou et al., Reference Zou, Jiang, Dai, Yue, Chen and Wang2019). Narote et al. Reference Narote, Bhujbal, Narote and Dhane2018 summarised and compared the lane detection methods based on the mark features. Generally, edge detection, colour segment, Hough transform or Kalman filter are applied to detect the lane boundaries. The main drawback of these methods is that they are affected by road scene variations (e.g. lane numbers). Although machine-learning-based methods usually need to train the complicated model with a big dataset, the detection can be more robust and accurate. LaneNet method is applied for lane detection in this paper. The model is trained with TuSimple datasets (TuSimple, 2017). The general procedures are summarised in the following paragraph (Neven et al., Reference Neven2018; MaybeShewill-CV, 2019).

The convolutional network consists of two branches, namely binary segmentation and instance segmentation, both of which are based on encoder-decoder network Enet (Paszke et al., Reference Paszke, Chaurasia, Kim and Culurciello2016). The network architecture is shown in Figure 1. The binary segmentation branch is to distinguish the lanes from the background while the instance segmentation branch can distinguish the lanes from each other. The two branches share most of the encoder stages except for the last one, which makes the last stage in encoder and decoder the backbone of each branch. The loss function for the binary segmentation branch is cross-entropy with bounded inverse class weighting (Paszke et al., Reference Paszke, Chaurasia, Kim and Culurciello2016) to avoid problems caused by high imbalance between the lanes and backgrounds. For the instance segmentation branch, distance metric learning-based discriminative loss functions are applied to disentangle the lane pixels identified by the segmentation branch (De Brabandere et al., Reference De Brabandere, Neven and Van Gool2017). Clustering is applied to separate the pixel embedding of different road lanes. Specifically, the loss function consists of three terms, namely variance term, distance term and regularisation term, given as

Figure 1. LaneNet architecture. The encoding stage is denoted by blue arrows while the decoding stage is represented by dash orange arrows.

(3)$$\begin{array}{l} L_{var} =\frac{1}{C}\sum\limits_{c=1}^C {\frac{1}{N_c }} \sum\limits_{i=1}^{N_C } {[ {\vert \vert \mu_c -x_i \vert \vert -\delta v} ] } _+^2 \\ L_{dist} =\frac{1}{C(C-1)}\sum\limits_{c_A =1}^C \sum\limits_{c_B =1}^C [ {2\delta_d -\vert \vert \mu_{c_A } -\mu_{c_B } \vert \vert } ] _+^2 \\ L_{reg} =\frac{1}{C}\sum\limits_{c=1}^C {\vert \vert \mu _c \vert \vert } \\ L=L_{var} +L_{dist} +0.001L_{reg} \\ \end{array}$$

where L is the discriminative loss function; L _var, L _dist and L _reg are the variance term, distance term and regularisation term respectively; C is the number of clusters in the ground truth; N _C is the number of elements in cluster c; x _i is the pixel embedding of the lane; μ_c is the mean of embedding, namely the centre of the cluster; μ_cA and μ_cB represent different centres; [x]₊ = max(0, x); δ v and δ d are the margins of variance and distance loss.

The variance term imposes an intra-cluster pull force on embedding pixels towards the centre while the distance term pushes the clusters away from each other. The regularisation term is a small pull that keeps the activations bounded. With such a network model, the road lanes can be detected for further applications.

3.2. Map matching

Map matching aims to determine the road link on the digital map on which the land vehicle is moving. Based on the land vehicle state (e.g., position, velocity, azimuth etc.) and digital map information, the land vehicle position can be projected on the road link of the digital map. Generally, map matching algorithms can be classified into three categories, namely geometric map matching, topological map matching and advanced map matching (e.g. fuzzy logic map matching). Geometric map matching determines the correct road link based on the relationship between the vehicle localisation and digital map entities. For topological map matching, the topology of the road network is also taken into consideration. Fuzzy logic map matching utilises fuzzy logic to find out the best match of road link, which has been proved to be the one with least mismatching (Quddus, Reference Quddus2006). With an accessible precise digital map, map matching provides external position information in addition to GNSS positioning (Liu, Reference Liu2017; Balazadegan Sarvrood et al., Reference Balazadegan Sarvrood, Liu and Gao2018). In the past, the accuracy of the GNSS solution was low, therefore, it was necessary to adopt more complicated algorithms to avoid mismatching. However, with the development of GNSS positioning technology, the horizontal accuracy could reach lane level (sub-metre) even with low-cost receivers, which enables correct map matching using simply geometric map-matching algorithms. In this research, geometric point-to-point matching is applied. This method finds the closest node of the road network based on the position of the land vehicle, given that the SF-PPP and SF-RTK float solutions can reach decimetre accuracy. Usually, there are multiple lanes within one road, which makes map matching useless without lane information in precise GNSS positioning. In this paper, the precise lane-level digital map generated based on RTK fixed solution is applied. In addition, the lane detection method introduced before is adopted to find out the lane that the vehicle is moving on. The procedure is summarised in Figure 2.

Figure 2. SF-GNSS with aid of lane-level digital map matching.

In Figure 2, GNSS first provides the initial solution, which is projected on the closest road link with two lanes. The digital map used contains the information for both lanes. Therefore, with the detected lane number obtained by lane recognition, the GNSS solution can be projected on the correct lane link.

4. INTEGRATION OF SF-GNSS AND MAP MATCHING

The projected point on the lane link, of which the accuracy can reach a few decimetres on the cross-road direction, can be further applied to reduce the SF-PPP positioning error or used to improve the SF-RTK fixing rate. To be more specific, the two-dimensional (2D) coordinates of the projected points can be applied as additional observations other than the GNSS measurements in the Kalman filter. The processing procedure of the integration is summarised in Figure 3.

Figure 3. Processing procedures for SF-GNSS with aid of lane-level map matching.

In Figure 3, the inputs for map matching are the same as in Figure 2, namely the initial GNSS solution and detected lane number. For map matching, a modified OpenStreetMap (OSM) shapefile is loaded in an open-source relational database management system (RDBMS) named PostgreSQL, which contains several types of the lane link information (e.g. the link ID, link geometry, link type, lane number, link azimuth, start and end node IDs and coordinates in Universal Transverse Mercator etc.). In map matching, the PostgreSQL is connected to manage the spatial queries. One thing that needs to be mentioned is that the accuracy of the original OSM shapefile is in metres, which cannot reach the requirement of this work. Besides, the lane links in each road are unavailable. Therefore, to access lane-level digital maps, the modified shapefile based on RTK fixed solution on each lane is made by QGIS, which is a free geographic information system application that allows creating, editing, visualising and analysing geospatial data. In this shapefile, lane information on several roads is provided. In Figure 3 a feedback filter is designed to exclude those possible mismatched points in map matching, which is important since faulty projected points lead to completely incorrect results. In this work, there are three principles of the feedback filter which are listed as follows. The first one is that if the difference between the heading of the land vehicle and the lane link is larger than a threshold, map matching stops. The second one requires the distance from the land vehicle to a junction less than a threshold (e.g. 5 m) to have valid map matching projected points. The last one stops map matching when lane detection fails.

The additional positioning information from map matching is further applied to improve the SF-GNSS positioning results. Specifically, the Kalman filter is applied to integrate the SF-GNSS and map matching. The projected points on the road lanes from map matching are applied as pseudo measurements in the Kalman filter functioning as constraints of the positioning results. The general Kalman filter calculation procedure is defined in Equation (4) as (Welch and Bishop, Reference Welch and Bishop1995)

(4)$$\begin{array}{l} X_{k+1}^- =\Phi _{k+1,k} X_k +w_k \\ P_{k+1}^- =\Phi _{k+1,k} P_k \Phi _{k+1,k}^T +Q_k \\ K_{k+1} =P_{k+1}^- H_{k+1}^T (H_{k+1} P_{k+1}^- H_{k+1}^T +R_{k+1} )^{-1} \\ X_{k+1} =X_{k+1}^- +K_{k+1} (Z_{k+1} -H_{k+1} X_k^- ) \\ P_{k+1} =(I-K_{k+1} H_{k+1} )P_{k+1}^- \\ \end{array} $$

where X _k+1⁻ is the predicted state vector; w _k is the process noise vector; P _k+1⁻ is the covariance of the predicted state vector; Φ_{k+1, k} is the transition matrix between two consecutive epochs; Q _k is the process noise covariance; K _k+1 is the optimal gain matrix; H _k+1 is the corresponding design matrix; R _k+1 is the covariance matrix of the measurements; Z _k+1 is the measurement vector; I is the identity matrix; X _k+1 is the updated state vector after applying the measurement; P _k+1 is the updated state vector covariance.

In this paper, the position errors, ambiguities and receiver clock offset (for SF-PPP) are estimated as the state vector in the Kalman filter. In addition to the GNSS measurements in Equations (1) and (2) the difference between the GNSS solution and the map matched point is also applied in the measurement update in the Kalman filter. By projecting the GNSS solution on the lane link, the map matching measurement, namely the coordinate difference in the cross-lane direction, can be given as

(5)$$Z_{mm} =\delta c$$

where Z _mm is the additional measurement formed by map matching, namely the distance between the GNSS solution and the map matched point by projection. The corresponding design matrix can be formed as

(6)$$H_{mm} =\left[\begin{array}{@{}ccccc@{}} {r_{11} } & {r_{12} } & {r_{13} } & 0 & \cdots \end{array} \right]$$

where r _1i represents the elements for the position corrections in the state vector while the other elements are 0. The r _1i is the first row of the rotation matrix R from ECEF (Earth-centred, Earth-fixed frame) to the lane frame that is defined as the x-axis pointing to the cross-lane direction, y-axis pointing to the along-lane direction while z-axis completes the right-handed coordinate system. The lane frame is defined by the rotation of the local-level frame with a yaw angle y, which is shown in Figure 4.

Figure 4. Definition of lane frame.

The rotation matrix R from ECEF to the lane frame can be given as

(7)$$ R=R_l^m R_e^l =\left[\begin{array}{@{}ccc@{}} {r_{11} } & {r_{12} } & {r_{13} } \\ {r_{21} } & {r_{22} } & {r_{23} } \\ {r_{31} } & {r_{32} } & {r_{33} } \end{array} \right]$$

where R is the rotation matrix from ECEF to the lane frame; r ₁₁, r ₁₂ and r ₁₃ are the elements in the design matrix in Equation (6); $R_{l}^{m}$ and $R_{e}^{l}$ are the rotation matrices from ECEF to local-level frame and from local-level frame to the lane frame respectively, given as

(8)$$ R_l^m =\left[\begin{array}{@{}ccc@{}} {\cos y} & {\sin y} & 0 \\ {-\sin y} & {\cos y} & 0 \\ 0 & 0 & 1 \\ \end{array} \right]R_e^l =\left[\begin{array}{@{}ccc@{}} {-\sin \lambda } & {\cos \lambda } & 0 \\ {-\sin \varphi \cos \lambda } & {-\sin \varphi \sin \lambda } & {\cos \varphi } \\ {\cos \varphi \cos \lambda } & {\cos \varphi \sin \lambda } & {\sin \varphi } \end{array} \right]$$

where y is the yaw angle; φ and λ are the latitude and longitude.

5. FIELD TESTS AND RESULTS

To verify the effectiveness of the proposed integration of map matching and SF-GNSS, a field kinematic test was conducted in Calgary. The test setup is shown in Figure 5. A u-blox M8U receiver and a NovAtel receiver were connected to a geodetic antenna on the roof of a land vehicle through a signal splitter. The u-blox receiver is used to conduct the SF-PPP and SF-RTK while the NovAtel receiver is used to conduct dual-frequency RTK with IGS (International GNSS Service) station UCAL as the base. The dual-frequency fixed RTK solution is applied as the reference to evaluate the performance of different methods. For the SF-RTK, the base station is also UCAL. In SF-GNSS processing, GPS code and phase measurements were applied. A PIXELINK camera is used to take images at 1 FPS (frame per second) for lane detection. The u-blox receiver is used to trigger the camera to take images with GPS time tag. All the datasets are logged by a laptop and the integration processing is post-mission processed using the collected datasets. The current implementation of SF-GNSS and map matching cannot achieve real-time solutions.

Figure 5. Kinematic test setup.

During the test, the land vehicle was moving on 32 Avenue NW, Shaganappi Trail NW, and 33 Street NW. The test environment is suburban, with trees and buildings along the road. For the land vehicle trajectory, 32 Avenue NW and Shaganappi Trail NW are two-lane roads with the more accurate lane information in the digital map while the others are just one-lane roads. When the land vehicle was approaching intersections or ramps, map matching was not applied. The trajectory and lane information are indicated in Figure 6. The roads with two lanes, one lane, and no map matching are marked in magenta, cyan and yellow respectively.

Figure 6. Kinematic test trajectory.

The digital map is made based on the RTK solution with attributes containing the link number, lane order and total number of lanes on the road. When the test was conducted, the digital map was ready. During the test, there were two lanes at most on any road. The right lane and left lane are marked as lane 1 and lane 2 respectively in the digital map. Based on the GPS initial positioning result, the road can be selected on the digital map. With the lane detection result from the images taken by the camera, the lane order can be determined and the GPS positioning result can be projected on the lane link. Figure 7 shows lane detection based on the image taken by the camera. The lane in the left image is detected and marked as white pixels in the right binary image. Based on the location of the detected lane in the image, the moving lane order can be determined.

Figure 7. Lane detection.

The satellite number, PDOP (position dilution of precision) and HDOP (horizontal dilution of precision) during the test can be found in Figure 8. During the test, more than eight satellites can be observed for most of the time. The drop of satellite number around 370 s is due to signal blockage by a bridge.

Figure 8. Satellite number and DOP(dilution of precision) values during the test.

With the external position information, based on Equations (4)–(7), SF-RTK ambiguity fixing rate can be further improved and SF-PPP positioning accuracy can be slightly enhanced with the removal of large errors. Figure 9 shows the SF-RTK positioning results in each direction. The ambiguity ratio in this test is set as 3. In Figure 9, the periods with different numbers of lanes are marked in different colours as in Figure 6. The red, green and blue lines represent positioning errors in east, north and up directions, respectively. State 1 means ambiguities are fixed while state 2 represents unsuccessful ambiguity resolution. The epoch-wise ambiguity fixing method is applied here in which the ambiguities and position are reinitialised every epoch. In this way, if the ambiguities are not fixed, the performance is similar to DGNSS since the Kalman filter is reinitialised. It can be seen that the ambiguity fixing rate is improved, especially during the period with double lanes where the digital map is more accurate. The number of epochs with ambiguities fixed was improved from 252 to 491. With the total number of epochs as 936, the ambiguity fixing rate increased from 26.9% to 52.5%. Due to the improved ambiguity fixing rate, the overall performance of SF-RTK is improved as well. The positioning errors in the horizontal plane and vertical direction were reduced from 1.21 m to 0.80 m and from 1.11 m to 0.91 m respectively.

Figure 9. Positioning errors of epoch-wise SF-RTK and SF-RTK/MM.

For the 463 epochs when the land vehicle was moving on the road with two lanes, the number of epochs with successful ambiguity resolution increased from 122 to 333 with the aid of map matching. The improvement during these periods is more obvious, compared with single-lane periods, because of higher digital map accuracy. For the 140 epochs of single-lane periods, the number of fixed solutions was improved from 35 to 58. The improvement for the single-lane period is not as good as the double-lane period because of degraded digital map accuracy.

SF-PPP using the u-blox receiver is also processed with the aid of map matching. The processing strategy of PPP is summarised as follows. Precise products of satellite orbit and clock are applied to correct the satellite-end errors. Since the precise products are referred to ionospheric-free (P1, P2) code, C1-P1 and P1-P2 DCB products are needed to correct the group delay for C1 code measurements of the u-blox receiver. Ionospheric delay is corrected by CODE ionospheric global map while the tropospheric delay is corrected using the Saastamoninen model. The position unknowns and receiver clocks are estimated in the Kalman filter. The residuals of ionospheric delays and tropospheric delays are not included in estimation to simplify the estimation model and reduce the convergence time in practical application. With such a strategy, we can achieve fast solutions with an accuracy of about half a metre. The PPP positioning error plot is shown in Figure 10.

Figure 10. Positioning errors of SF-PPP and SF-PPP/MM.

It can be seen that the horizontal accuracy of SF-PPP is improved from 0.45 m to 0.31 m with the aid of map matching. Some large horizontal positioning errors between 450 s and 550 s are reduced with the aid of map matching. Generally, during the 463 epochs of double-lane periods, the horizontal positioning RMS was reduced from 0.52 m to 0.30 m. The improvement in height is less, compared with the horizontal plane, since the 2D digital map is applied without height information. The error in the up direction is reduced from 0.41 m to 0.34 m due to the correlation of position unknowns in different directions.