Optimized GNSS/INS Integration with Time Synchronization for High-Accuracy Smartphone Positioning

2.1 DGNSS Post-Processing

A conceptual diagram of our DGNSS post-processing methodology is shown in Figure 1. Figure 1 shows the overall processing flow, which follows the standard RTKLIB post-processing pipeline. Before DGNSS processing, raw GNSS log data were converted into observational data in receiver independent exchange (RINEX) format using the “Android GNSS logger to RINEX converter” tool (Everett, 2022), to optimize the number of GNSS measurements fed into DGNSS. The filtering rules proposed by Everett et al. (2022) for removing invalid or low-quality raw measurements were adopted without modification, with the exception of the threshold for the “maximum carrier-phase uncertainty.” Specifically, this uncertainty corresponds to the per-measurement indicator AccumulatedDeltaRangeUncertaintyMeters provided by the Android GNSS Raw Measurement API (Google, 2024), which represents the 1σ Accumulated Delta Range (ADR) uncertainty in meters. Carrier-phase measurements are rejected when this value exceeds a threshold; in this study, the threshold was relaxed from 0.1 to 0.3 m to retain more usable measurements and increase observation availability for DGNSS in the noisier GSDC data sets. Navigation messages (BRDM), containing navigation parameters for all GNSS constellations, were obtained from the data archive of the Crustal Dynamics Data Information System. Base station observations at time intervals of 15 s were obtained from UNAVCO. The selected nearby stations include SLAC, P217, P222, and P225 in San Francisco and TORP, VDCY, and MLFP in the Los Angeles area. Common errors in smartphone GNSS observations, such as satellite orbit error, satellite clock error, atmospheric delay, satellite signal bias, and tidal terms, were mitigated by differentiating the GNSS observations between a smartphone and base station. The EKF computes forward and backward solutions, encompassing position, velocity, acceleration, and float ambiguities.

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 1

Schematic of our carrier-phase DGNSS methodology

Integer ambiguity resolution was not applied in this DGNSS processing; consequently, carrier-phase ambiguities were estimated as float states. Although the feasibility of ambiguity resolution with modern smartphones has been demonstrated under controlled zero- or short-baseline conditions (Li et al., 2022; Odolinski et al., 2024; Yong et al., 2021), ambiguity-fixing performance remains highly sensitive to the antenna environment, multipath, and device placement. Specifically, the integer least-squares fixing success rate is substantially degraded when native smartphone antennas are used in close proximity to surrounding surfaces (Yong et al., 2021). Moreover, in-vehicle placement behind a windshield can significantly reduce the carrier-to-noise density ratio ( C/N0 ) and the fixing rate owing to signal attenuation; for instance, a reduction from approximately 90.5% to 36.2% has been reported in such environments compared with open-sky conditions (Li et al., 2022). Given that the GSDC data sets were collected with smartphones located inside vehicles under dynamic driving conditions, we prioritized robustness and solution availability across heterogeneous devices and environments. Therefore, float solutions were adopted in this study.

Finally, a smoothed solution was generated by weighting the state estimates of the forward and backward filters based on their inverse covariances.

The observation vector comprises a double-differenced carrier-phase vector ϕ, double-differenced pseudorange vector P, and single-differenced doppler vector d, expressed as follows:

y=[ϕ1Tϕ5TP1TP5Td1Td5T]T1

The subscript indices (1 and 5) indicate the corresponding indices of the GNSS signals. The state variables include position r, velocity v, acceleration a, and the rover-base single-differenced ambiguities B1,B5 defined in Equation (3). The residual ionospheric and tropospheric delays were not estimated because the baseline lengths in this data set were relatively short (typically <20km ). For such baselines, the differential atmospheric residuals are generally negligible compared with the dominant smartphone multipath and noise contributions:

x=[rTvTaTB1TB5T]T2

Bf=[Brb,f1Brb,f2…Brb,fm]T3

For a given receiver r, base station b, satellite indices p and q, signal index f, and measurement error ϵ, the measurement equation can be expressed as follows:

ϕrb,fpq=ρrbpq+λf(Brb,fp−Brb,fq)+ϵϕ4

Prb,fpq=ρrbpq+ϵP5

dr,fpq=ρ˙rpq+ϵd6

The measurement model vector h can be expressed as follows:

h(x)=[hϕ,1Thϕ,5ThP,1ThP,5Thd,1Thd,5T]T7

hϕ,f=[ρrb12+λi(Brb,f1−Brb,f2)ρrb13+λi(Brb,f1−Brb,f3)⋮ρrb1m+λi(Brb,f1−Brb,fm)],hP,f=[ρrb12ρrb13⋮ρrb1m],hd,f=[ρ˙r12ρ˙r13⋮ρ˙r1m]8

where ρ is the geometric range and λ and i are the wavelength and corresponding index of the GNSS L1 and L5 signals, respectively. The partial derivative matrix Hk and measurement noise covariance matrix Rk for epoch k can be expressed as follows:

Hk=∂h(x)∂x=[−DEOOλ1DO−DEOOOλ5D−DEOOOO−DEOOOOO−DEOOOO−DEOOO]9

E=[er1er2⋯erm]T10

Rk=diag[DRϕDT,DRPDT,DRdDT]11

where e is the unit line-of-sight vector from the receiver to the satellite. D represents the between-satellite differencing matrix with respect to the pivot (reference) satellite. Note that B is not fully observable under double-differenced (DD) modeling because only DB appears in the residuals (one common-mode/gauge component is unobservable). We fix this gauge following RTKLIB by applying a common offset to all active ambiguity states in the phase-bias update that enforces phase-code coherence.

Rϕ,RP, and RD are the measurement noise matrices for the carrier phase, pseudorange, and Doppler, respectively. The time update of the Kalman filter is expressed as follows:

x―k+1=Fx^k12

P¯k+1=FP^kFT+QτGNSS13

F=[I3×3I3×3τGNSSI3×3τGNSS22O3×mO3×mO3×3I3×3I3×3τGNSSO3×mO3×mO3×3O3×3I3×3O3×mO3×mOm×3Om×3Om×3Im×mOm×mOm×3Om×3Om×3Om×mIm×m]14

Q=diag[O3×3,O3×3,Qa,QB1,QB5]15

The discrete-time state transition matrix F is derived from a constantacceleration kinematic model, and the parameter τGNSS denotes the discretization time step corresponding to the GNSS observation interval (1 s). Process noise is applied to the acceleration and ambiguity states, which are both modeled as random-walk processes.

The components of both ambiguity process noise matrices, QB1 and QB5, are set to 0.00012m2/s, consistent with the default configurations of RTKLIB. For the acceleration states, the process noise matrix Qa is defined in the Earth-centered, Earth-fixed (ECEF) frame for state propagation. This matrix is derived by transforming a locally defined covariance matrix into the ECEF frame to reflect typical vehicle dynamics. Specifically, the acceleration process noise is first defined anisotropically in a local-level frame, with the horizontal components set to 0.22(m/s2)2/s and the vertical component set to 0.12(m/s2)2/s. The resulting local-level covariance is then rotated into the ECEF frame to form Qa.

The measurement update is as follows:

Sk=HkP¯kHkT+Rk16

vk=yk−h(x¯k)17

Kk=P¯kHkTSk−118

x^k=x―k+Kkvk19

P^k=(I−KkHk)P―k20

where vk represents the innovation vector, Sk denotes the innovation covariance matrix, and ∧ and −denote the a posteriori and a priori estimates, respectively.

A smoothed solution is obtained by combining the forward and backward estimates:

x^kfb=(P^kf−1+P^kb−1)−1(P^kf−1x^kf+P^kb−1x^kb)21

where (x^kf,P^kf) and (x^kb,P^kb) denote the forward and backward filter estimates and their associated error covariances, respectively.

To enhance robustness against outliers and accommodate low-quality smartphone GNSS measurements, several enhancements were incorporated into our processing engine. The following minor changes from RTKLIB 2.4.3 were implemented:

• Doppler measurements were included in the observation vector, as expressed in Equation (1), to impose additional range-rate constraints, thereby directly refining the estimation of receiver velocity. It is well established that enhancing velocity observability improves instantaneous positioning performance, particularly under conditions of weak measurement geometry or limited satellite availability (Kubo, 2009). In the specific context of smartphone navigation, velocity-aided single-frequency positioning leveraging Dopplerderived range rates has been demonstrated to enhance accuracy in real-world driving scenarios (Han et al., 2024).

• Cycle-slip detection was enhanced for single-frequency observations. A detection process was implemented in which a cycle slip is identified if the discrepancy between the time-differenced carrier phase and the Dopplerderived range rate exceeds a threshold L (Everett et al., 2022; Suzuki, 2023):

|λi(ϕk+1−ϕk)/Δt−dk+1|>L22

where λi is the carrier wavelength, ϕ is the carrier phase in cycles, Δt is the measurement interval, and d denotes the Doppler-derived range rate expressed in meters per second. Accordingly, the threshold L is defined in units of meters per second. In this study, the threshold was set to L=1m/s, empirically determined based on the Doppler noise level observed in representative GSDC data sets. Furthermore, this threshold value aligns with those successfully employed in prior GSDC studies conducted for similar signal characteristics and driving environments (Suzuki, 2023).

• Innovation testing was refined to detect and eliminate outliers. In the RTKLIB code, the threshold Tmax for innovation vk was set at a constant value as follows:

|vk|<Tmax23

This study implements a normalized innovation test to account for the varying quality of measurements. Specifically, we set σs,k=[Sk]ii as the standard deviation of the innovation vk, derived from the i-th diagonal element of the innovation covariance matrix Sk (Equation (16)). The term vk represents the innovation scalar for a given pseudorange, carrier-phase, or Doppler measurement. Outliers are detected and rejected using a 4σ criterion, defined by the normalized threshold γ=4:

|vkσs,k|<γ24

The processing blocks in Figure 1 follow the standard RTKLIB workflow; the differences from the baseline RTKLIB processing are implemented as localized code modifications inside the corresponding blocks. The primary modifications in RTKLIB 2.4.3 b34 are summarized in Table 1.

2.2 Tightly Coupled GNSS/INS Integration

The DGNSS post-processing method was expanded to include tightly coupled GNSS/INS integrated navigation systems by incorporating inertial navigation processing and IMU bias estimation/correction. Figure 2 highlights the extension from the DGNSS processing to the tightly coupled GNSS/INS scheme; the blocks newly introduced for GNSS/INS integration are emphasized with thick borders. While the EKF structure used in DGNSS (time update and measurement update) is retained, the major modification is that the state time propagation is performed by inertial navigation mechanization using IMU measurements, instead of the constant-acceleration kinematic model. The GNSS/INS integration filter receives pseudorange, Doppler, and carrier-phase measurements from the GNSS observation correction. These measurements are used to calculate corrections to attitude, velocity, and position estimates, along with accelerometer and gyroscope biases. The estimated accelerometer and gyroscope biases are then employed in the IMU bias correction to calculate the true specific force f and true angular rate ω. The timestamp in the IMU data set is pre-adjusted using the offset calculated according to the method described in Section 2.3. During the inertial navigation processing, the position is maintained by integrating velocity, which is derived by integrating acceleration measurements obtained from the IMU. Attitude is maintained by integrating angular rate measurements. Furthermore, forward and backward solutions are computed in the GNSS/INS integration process.

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 2

Tightly coupled integration of the GNSS/INS architecture

The observation vector comprises the attitude θ, velocity v, position r, accelerometer bias ba, gyroscope bias bg, and ambiguity B. The attitude vector θ comprises roll, pitch, and yaw angles:

x=[θtbTvtbbTrtbtTbaTbgTB1TB5T]T25

Here, the Euler angle vector θβα indicates the orientation of frame α with respect to the reference frame β, and vβαγ represents the velocity of frame α with respect to frame β, expressed in frame γ⋅rβαγ represents the position of frame α with respect to frame β, as expressed in frame γ. Frames b,t,e, and i are defined as follows:

• b : The IMU frame is centered at and aligned with the IMU accelerometer triad.

• t : The local tangent plane frame has a fixed origin relative to the starting point of the drive. Its z-axis aligns vertically downward, and its x- and y-axes align northward and eastward, respectively.

• e: ECEF coordinate frame.

• i: Earth-centered inertial coordinate frame.

In this study, the filter estimates the smartphone attitude (roll, pitch, and yaw) rather than the vehicle attitude. Because vehiclemotion constraints, such as non-holonomic constraints or zero-velocity updating (ZUPT), are not incorporated, explicit alignment between the smartphone and vehicle frames is not required in the present formulation. Consequently, the initial attitude is treated as unknown and initialized to zero, allowing the angles to converge gradually through the Kalman filtering and smoothing process. Although utilizing the gravity vector or magnetometer measurements for initialization could potentially accelerate convergence, this approach was not prioritized, as a forward-backward smoother was applied to ensure trajectory consistency.

The error model for an accelerometer triad is simplified as follows:

f∼ibb=fibb+ba+na26

where f∼ibb denotes the IMU-output specific force vector, fibb represents the true specific force, ba represents the accelerometer bias vector, and na represents the random noise vector of the accelerometer. The accelerometer scale factor and cross-coupling error were not considered. Similarly, for a gyroscope, we have the following:

ω∼ibb=ωibb+bg+ng27

where ω∼ibb denotes the IMU-output angular rate vector, ωibb denotes the true angular rate, bg denotes the gyroscope bias vector, and ng denotes the gyroscope random noise vector. ωβαγ represents the orientation rates of the α-frame axes with respect to the β-frame axes, expressed in the γ-frame axes. The gyroscope scale factor and cross-coupling error were not considered in this study.

The three-dimensional navigation equations are expressed as follows. These equations assume a strapdown configuration, where the IMU is fixed relative to the smartphone:

r˙tbt=Ctbvtbb28

v˙tbb=fibb+Cebgeb(reb)−(Ωibb+Ωibe)vtbb29

q˙ib=12qibωibb30

Here, geb represents the acceleration due to gravity and is the sum of the gravitational and centrifugal accelerations. Ctb and Ceb represent the transformation matrices from frames b to t and b to e, respectively. qib denotes the quaternion attitude of frame b with respect to frame i.Ωβαγ represents the skew-symmetric matrix of the angular rate vector ωβαγ :

Ωβαγ=[0−ωβα,zγωβα,yγωβα,zγ0−ωβα,xγ−ωβα,yγωβα,xγ0]31

The error covariance matrix can be expressed as follows:

P―k+1=ΦkP^kΦkT+GQGT32

where Φ represents the state transition matrix, G represents the system noise distribution matrix, and the time interval of the IMU data is 0.01 s. We have the following:

Q=diag[SaI3×3,SgI3×3,SbaI3×3,SbgI3×3,SBIm×m,SBIm×m]33

where Sa,Sg,Sba,Sbg, and SB denote the squared standard deviations that define the discrete-time process-noise terms for the accelerometer random noise, gyroscope random noise, accelerometer-bias variation, gyroscope-bias variation, and ambiguity variation, respectively. In this study, these parameters were set to Sa=0.012(m/s2)2,Sg=0.0012(rad/s)2,Sba=0.00012(m/s2)2,Sbg=0.000012(rad/s)2, and SB=0.0000012m2. The noise variances Sa and Sg were empirically determined from representative stationary IMU segments by removing the sample mean for each axis and calculating the variance of the resulting time series. The ambiguity process noise SB was set to be consistent with the default configuration in RTKLIB. Finally, the bias process noise parameters, Sba and Sbg, were selected via a parameter sweep over representative GSDC data sets to optimize positioning performance.

The observation matrix is expressed as follows:

Hk=∂h(x)∂x=[OO−DECetOOλ1DOOO−DECetOOOλ5DOO−DECetOOOOOO−DECetOOOOCebΛ−DECebOOOOOCebΛ−DECebOOOOO]34

Zero-sideslip constraints and ZUPT were not incorporated, although both can be effective for low-cost IMUs. Non-holonomic (zero-sideslip) constraints were omitted because the smartphone IMU lacked a prescribed alignment to the vehicle frame and its rigidity relative to the vehicle throughout the drive could not be guaranteed. Applying such constraints without a reliable mounting model may introduce biased pseudo-measurements, potentially degrading the overall estimation accuracy:

Λ=[0−vtb,zbvtb,ybvtb,zb0−vtb,xb−vtb,ybvtb,xb0]35

2.3 Synchronizing IMU Data with GNSS Observations

Accurate synchronization between IMU data and GNSS observations is crucial to achieve high-precision positioning. GNSS measurements were timestamped using a GNSS clock, and IMU measurements were timestamped using the system clock (Sharma et al., 2021). In practice, the system clock is not always precisely synchronized with the GNSS clock. Although synchronization information between the system and GNSS clocks, such as ChipsetElapsedRealtimeNanos in GNSS raw data and ElapsedRealtimeNanos in IMU data, is available for certain driving data sets, participants typically ensured temporal consistency between GNSS measurements and IMU data using their own methods for most data sets. To address this issue, we proposed a novel timestamp adjustment approach. In particular, this method estimates the time offset (δt) between the system clock and GNSS times and adjusts the IMU timestamps accordingly, as expressed in Equation (36):

Adjusted_utcTimeMills=utcTimeMills +δt36

where utcTimeMills represents the original timestamp included in the IMU data and Adjusted_utcTimeMills denotes the timestamp adjusted using the offset estimated by δt. Figure 3 illustrates the concept of the proposed time-offset estimation method. In this figure, the vehicle moves from the bottom to the top of the slide. The red and black lines represent the true trajectory and GNSS/ INS-integrated solution, respectively. The orange dots indicate epochs at which GNSS measurements are available. Between these epochs, the vehicle position is propagated by integrating IMU measurements. At each GNSS epoch, the propagated position is corrected by applying a GNSS measurement update, resulting in a position correction vector Δr. The position correction can be expressed as follows:

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 3

Time synchronization method

Δrk(δt)=Er(δt)Kk(δt)vk(δt)37

where Er=[I30] extracts the position components from the state vector. If the system clock is not accurately synchronized with the GNSS clock, the magnitude of Δr increases. Taking advantage of this characteristic, the proposed method estimates the time offset δt by minimizing the evaluation function, which is defined as the RMS of the position corrections, as expressed in Equation (38):

Jk(δt)=1N∑k=1NΔrk(δt)238

δt^=arg⁡minδtJk(δt)39

where N represents the total number of epochs for which GNSS measurements are available in each driving data set. In this study, we assume that the time offset remained constant throughout each individual drive. The time offset δt is not assumed to be constant across all drives. Instead, it is determined on a per-drive basis, as the optimal offset exhibits significant variability among individual drives and across different smartphone models. The observed ranges of δt for all training drives, categorized by smartphone type, are reported in the final paragraph of Section 3.1.

Algorithm 1 summarizes the proposed procedure for estimating the time offset between GNSS and IMU data. Starting from an initial guess δt0, the method explores the candidate offsets within a predefined search range M using a step size Δt. The IMU timestamps are adjusted accordingly for each candidate δti, and the innovation-based position correction Δrk is evaluated across all GNSS epochs. Then, the evaluation function J(δti) is computed. The offset δt^ that minimizes J(δti) is selected as the final estimate, which represents the optimal synchronization between the GNSS and IMU clocks. The initial offset was initialized to δt0=0s. The grid search employed a step size of Δt=20ms spanning the range of [-1.00, +1.00]s, yielding a set of candidate offsets with 20-ms spacing.

The relationship between the evaluation function value and the corresponding positioning score with respect to the timestamp adjustment δt is illustrated in Figure 4. As shown, a decrease in position correction results in a decrease in positioning score, indicating a positive relationship between the two variables. In this example, both values reach their minima at a timestamp adjustment of 0.2 s. Specifically, the positioning score without timestamp adjustment (i.e., a time offset of 0.0 s ) is 1.60 m , whereas applying a time offset of 0.20 s to the IMU utcTimeMillis timestamps reduces the score to 0.59 m. This example confirms the validity of the proposed approach, in which the time offset is estimated based on the magnitude of position corrections.

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 4

(Top) Relationship between timestamp adjustment (δt) and evaluation function value J(δt); (bottom) relationship between time timestamp adjustment (δt) and score, with means of 50% and 95%

Both relationships were calculated using the entire “2020-12-10-22-52-us-ca-sjc-c/mi8” driving data set.

2.4 Adaptive Measurement Weighting

Pseudorange residuals differ significantly, depending on the smartphone model and signal frequency. As shown in Figure 5, the double-differenced pseudorange residuals for GPS L1 observed with Pixel 5 are approximately twice those from the Xiaomi Mi 8. However, the differences between the smartphone models were minimal for GPS L5. Regarding the driving environment, Figure 6 shows that the pseudorange residuals recorded by the s908b under urban street driving conditions are approximately 1.5 times those recorded under highway (open-sky) conditions.

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 5

Comparison of double-differenced pseudorange residuals of each satellite signal and smartphone in the “2021-03-16-18-59-us-ca-mtv-a” trip with respect to the ground truth Panels (a) and (c) represent GPS L1C/A and L5 residuals for Pixel 5, respectively, and panels (b) and (d) represent GPS L1C/A and L5 residuals for Mi 8, respectively.

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 6

Comparison of double-differenced pseudorange residuals of Samsung sm-s908b GPS L1 in (a) the “2023-09-06-18-04-us-ca” open-sky trip and (b) the “2023-05-19-20-10-us-ca-mtv-ie2” urban street trip with respect to the ground truth

Therefore, achieving high-precision positioning requires careful adjustment of measurement noise weighting based on the smartphone model and driving environment. To address the discrepancies in pseudorange residuals, an adaptive measurement-weighting scheme was introduced. This study proposed an algorithm that adjusts pseudorange, carrier-phase, and Doppler measurement weighting using an IAE method (Akhlaghi et al., 2017). Similarly, innovation-driven covariance adaptation has also been investigated for smartphone GNSS measurements, where standardized innovations are assessed via hypothesis testing (e.g., a t-test) and the measurement covariance is scaled accordingly (Raghuvanshi & Bisnath, 2025). In contrast, this study adopts an IAE-based covariance-matching update to continuously track the time-varying measurement quality. The algorithm can be expressed as follows:

Rk+1=αRk+(1−α)(vkvkT−HkP―kHkT)40

where α is the forgetting factor ( 0≤α≤1 ). In this study, α=0.5 was employed for pseudorange and Doppler measurements to facilitate rapid tracking of environmental noise fluctuations encountered during driving. Conversely, a larger value of α=0.9 was selected for the carrier phase to maintain stability in the covariance updates, reflecting the relatively consistent noise characteristics of smartphone carrier-phase observations. Notably, the updated covariance matrix Rk+1 must remain positive-definite. However, Equation (40) does not inherently guarantee positive-definiteness, because Rk+1 is estimated by subtracting two positive-definite matrices. Therefore, to ensure that the resulting covariance matrix is positive-definite, only the diagonal elements of Rk+1 are updated, based on Equation (40), and the upper and lower limits are applied to these elements. Noise covariances were determined for the pivot satellite using an elevation-dependent weighting scheme implemented in RTKLIB.

3.1 Evaluation Using the Training Data Set

To assess the impact of the proposed method on horizontal positioning accuracy, we analyzed the positioning solutions computed by each method using the complete GSDC 2023-2024 training data set, with respect to the provided ground truth. Consistent with the GSDC scoring scheme, we evaluated horizontal positioning error only. Before evaluation, data sets corresponding to smartphone models not included in the test data set were excluded from the analysis. Five distinct positioning solutions were compared:

1. The PPK solution computed by RTKLIB demo5 (Everett, 2025), denoted as “RTKLIB demo5"

2. The DGNSS solution, as discussed in Section 2.1, denoted as “DGNSS"

3. The GNSS/INS integration method without timestamp adjustment or adaptive measurement weighting, denoted as “DGNSS/INS"

4. The GNSS/INS integration method with timestamp adjustment, denoted as “DGNSS/INS with time-adj."

5. The GNSS/INS integration method with timestamp adjustment and adaptive measurement weighting, denoted as “GNSS/INS with time-adj. and AEKF"

The RTKLIB demo5, which is a derivative of RTKLIB2.4.3 managed by Tim Everett, supports the latest raw receiver data and features enhanced algorithms for calculating high-precision positioning solutions for low-cost receivers. The availability of the source code and calculation modules for GSDC 2022 further confirms RTKLIB demo5 as a crucial benchmark and baseline solution. Regarding measurement weighting, the RTKLIB demo5 solution employs a hybrid scheme that combines elevation-dependent weighting with a C/N0-dependent term. In contrast, the proposed DGNSS solution utilizes solely an elevation-dependent weighting scheme to determine measurement variances.

The overall positioning scores, 50th percentiles, and 95th percentiles across all rides are provided in Figure 7, Table 2, and Table 3, respectively. The score in Figure 7 is defined as the mean of the 50th and 95th percentiles of horizontal error over the evaluation set:

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 7

Comparison of positioning errors among different smartphones and methods for 136 training data sets Filter parameters were not individually optimized for each phone type; however, optimized common parameters were applied to all smartphone models.

Score=12(P50+P95)41

where Pq denotes the q-th percentile of horizontal error. The evaluated Samsung smartphones included S-series smartphones (e.g., sm-g988b and sm-s908b). The Samsung S series includes models, such as the sm-g988b and sm-s908b, equipped with dualfrequency GNSS chips. The Samsung A series comprises models, such as sm-a325f, samsunga32, samsunga325g, sm-a205u, and sm-a505g, equipped with single-frequency GNSS chips. For the Samsung A50 (sm-a505u) and A20 (sm-a205u) models, a reversed accumulated delta range sign issue was confirmed; therefore, their GNSS positioning was processed after the sign had been reversed. Although the filter parameters were not optimized individually for each smartphone model, a set of common parameters optimized for general performance was applied to all smartphones. GNSS/INS integration with timestamp adjustment and adaptive measurement weighting achieved a score of 0.92 m for all smartphone models. The score decreased from 1.42 to 1.06 m following the introduction of IMU integration and further decreased to 0.92 m with adaptive measurement weighting. The timestamp adjustment in the GNSS/INS integration proposed in this study also contributed significantly to improving the score, further reducing the score by 0.11 m.

Comparative analyses of the positioning solutions of DGNSS and DGNSS/INS with time-adj. and AEKF are illustrated in Figures 8 and 9. Although the horizontal positioning error in the DGNSS solution increased near a tall building, DGNSS/ INS with time-adj. and AEKF significantly mitigated this error increase. Because not all phones were utilized during every test ride, the score differences among the smartphone models could not be easily analyzed. However, the positioning accuracies of Pixel 6 Pro and Samsung A series (single-frequency) were inferior to those of the other smartphones. This disparity can be attributed to the absence of a carrier phase in Pixel 6 Pro and the limitation of single-frequency measurements in the Samsung A series.

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 8

Comparison of vehicle trajectories calculated using (a) DGNSS post-processing and (b) GNSS/INS integration in the “2020-12-10-22-52-us-ca-sjc-c/pixel5” data set The yellow and red lines represent the solution and the ground truth, respectively.

• Download figure
• Open in new tab
• Download PowerPoint

FIGURE 9

Comparison of horizontal positioning errors of DGNSS (DGNSS) and GNSS/ INS integration (DGNSS/INS with time-adj. and AEKF) Errors were calculated using the entire “2020-12-10-22-52-us-ca-sjc-c/pixel5” driving data set. “Baseline (WLS)” represents the baseline solution, which is included in the provided data sets, derived from weighted least squares.

Finally, for the Pixel 4/4XL/5 models, synchronizing the system clock with GPS time resulted in minimal score differences between GNSS/INS integration solutions with and without timestamp adjustment. However, the optimal timestamp adjustments required to minimize the positioning errors varied significantly among smartphone models: the ranges were -0.04 to 0.10 s for the Pixel 4/4XL/5, -0.04 to 0.84 s for the Pixel 6 Pro, 0.00 to 0.28 s for the Pixel 7 Pro, 0.16 to 0.32 s for the Mi 8,0.08 to 0.20 s for the Samsung S series, and 0.10 to 0.16 s for the Samsung A series.

3.2 Evaluation Using a Private Data Set

Position solutions were calculated using a tightly coupled GNSS/INS integration method, with common parameters optimized using the training data set. The final scores achieved were 0.905 and 0.883 m for the public and private leaderboards, respectively, rendering a significant milestone, as the first time that accuracy within 1 m was achieved, securing first place at GSDC 2023-2024.