## 2 Materials and Methods
### Low-Complexity Correlation Approach  
In (6), the direct correlation calculation at the signal level is the cause of the high computational complexity of the DPE method. Theoretical studies have shown that although the direct correlation calculation can effectively suppress the cross-correlation interference, its effect on accuracy improvement is relatively limited in terms of the main source of localization errors [4]. In contrast, an alternative correlation calculation method can reduce the computational complexity while retaining the framework of the DPE and its joint processing in the navigation domain [13]. Specifically, using the autocorrelation function of the differences in time delay and Doppler shift between a candidate point and the satellite has been proposed to reduce complexity while ensuring that the correlation values correspond to the position domain, thus maintaining the dimensional consistency of DPE estimation. Moreover, it has been demonstrated that the correlation values for time delay and Doppler shift can be estimated independently with minimal accuracy loss [28]. By applying (4), the eight-dimensional optimization problem in (6) can be separated into two four-dimensional correlation functions for time delay and Doppler shift, respectively:  
\[R_{i}^{\tau}\left(\Delta\tau_{i}\right)= \begin{cases}1-\dfrac{ \left|\Delta\tau_{i}\right|}{T_{c}},&0\leq\dfrac{\left|\Delta\tau_{i}\right|} {T_{c}}\leq 1\\ 0,&\text{otherwise}.\end{cases}  \tag{7}\] \[R_{i}^{f}\left(\Delta f_{d_{i}}\right)=\operatorname{sinc}(\pi \Delta f_{d_{i}}T_{c}), \tag{8}\]  
where \(\Delta\tau_{i}=\tau_{i}-\tau_{i}\) and \(\Delta f_{d_{i}}=f_{d_{i}}-f_{d_{i}}\) represent the deviations of the candidate's time delay and Doppler shift from the true values for the \(i\)-th satellite. \(T_{c}\) represents the chip period of the PRN code. Notably, in static or low-speed dynamic scenarios where the velocity and clock drift are considered fixed, the MLE of position and time bias can be derived from (7). By leveraging the relationship between the code phase of the GNSS signal and the candidate's position and clock bias, the correlation function can be approximated by a pseudorange-based function to optimize the calculation:  
\[\Lambda_{i}(\hat{\gamma})=1-\left|\dfrac{||\hat{p}-\boldsymbol{p}||-\rho_{i} }{cT_{c}}+\dfrac{\delta\bar{t}-\delta t_{i}}{T_{c}}\right|. \tag{9}\]  
Beyond ensuring dimensional consistency, DPE inherently incorporates the optimal weighting of different satellites in (6), which is crucial for maintaining the validity of the MLE criterion [29]. Consequently, when performing correlation approximation, it is necessary to reassign weights for individual satellites. Since the approximated correlation function is based on pseudoranges, factors influencing pseudorange measurements must be considered. By analyzing the variation of pseudorange errors with satellite elevation and carrier-to-noise ratio (\(C/N_{0}\)), [PERSON] and [PERSON] [30] proposed a weighting model for LOS signal measurements, given by  
\[W_{i}=10^{\frac{CN_{i}}{100}}\times\sin^{2}(e_{i}), \tag{10}\]  
where \(e_{i}\) and \(CN_{i}\) represent the elevation and \(C/N_{0}\) of the \(i\)-th satellite, respectively. To better replicate the implicit weights of the various satellites derived from the signal correlations received by the DPE, these weights are incorporated into (9). Consequently, the low-computation correlation approach accounting for both signal power and elevation is given by  
\[\Lambda(\hat{\gamma})=\sum_{i=1}^{M}W_{i}\Lambda_{i}(\hat{\gamma}). \tag{11}\]  
In this approach, the correlation value of a candidate is directly calculated based on its time delay relative to the visible satellite. Additionally, the two factors that most affect the magnitude of the correlation value, the elevation and power of the signal, are taken into account to ensure the accuracy of the approximate correlation value.