Triple-Frequency Code-Phase Combination Determination: A Comparison with the Hatch-Melbourne-Wübbena Combination Using BDS Signals

: Considering the inﬂuence of the ionosphere, troposphere, and other systematic errors on double-differenced ambiguity resolution (AR), we present an optimal triple-frequency code-phase combination determination method driven by both the model and the real data. The new method makes full use of triple-frequency code measurements (especially the low-noise of the code on the B3 signal) to minimize the total noise level and achieve the largest AR success rate (model-driven) under different ionosphere residual situations (data-driven), thus speeding up the AR by directly rounding. With the triple-frequency Beidou Navigation Satellite System (BDS) data collected at ﬁve stations from a continuously-operating reference station network in Guangdong Province of China, different testing scenarios are deﬁned (a medium baseline, whose distance is between 20 km and 50 km; a medium-long baseline, whose distance is between 50 km and 100 km; and a long baseline, whose distance is larger than 100 km). The efﬁciency of the optimal code-phase combination on the AR success rate was compared with that of the geometry-free and ionosphere-free (GIF) combination and the Hatch-Melbourne-Wübbena (HMW) combination. Results show that the optimal combinations can always achieve better results than the HMW combination with B2 and B3 signals, especially when the satellite elevation angle is larger than 45 ◦ . For the wide-lane AR which aims to obtain decimeter-level kinematic positioning service, the standard deviation (STD) of ambiguity residuals for the suboptimal combination are only about 0.2 cycles, and the AR success rate by directly rounding can be up to 99%. Compared with the HMW combinations using B1 and B2 signals and using B1 and B3 signals, the suboptimal combination achieves the best results in all baselines, with an overall improvement of about 40% and 20%, respectively. Additionally, the STD difference between the optimal and the GIF code-phase combinations decreases as the baseline length increases. This indicates that the GIF combination is more suitable for long baselines. The proposed optimal code-phase combination determination method can be applied to other multi-frequency global navigation satellite systems, such as new-generation BDS, Galileo, and modernized GPS.

Many efforts have been made to discuss how to select useful linear combinations from an early stage.Following a systematic search for possible wide-lane (WL) combinations using dual-frequency signals by Cocard and Geiger [11], systematic investigations of optimal carrier phase combinations for the triple-frequency case were conducted.Richert and El-Sheimy [12] discussed the optimal combinations of triple-frequency GNSS signals, and gave the mathematical theory involved in creating linear combinations and developed three categories of combinations from the perspective of the ionosphere effects, the effects of thermal noise and multipath, and the troposphere effects.Cocard et al. [13] presented the concepts of lane-number to describe the wavelength of a particular combination, and analyzed the noise and ionosphere amplification factors with respect to the resulting wavelength for GPS signals.Urquhart [14] analyzed the multi-frequency carrier phase linear combination for GPS and Galileo and pointed out three factors to determine the optimal combinations, including the atmospheric errors, the observation noise, and the corresponding wavelength.Similar work for BDS signals has been done by Li et al. [15] and Zhang and He [16].For the application of AR using triple-frequency combinations, Han and Rizos [17] first extended the linear combination theory of carrier phase measurements, and proposed some typical carrier phase combinations with the free ionosphere and some with longer wavelengths.Feng [18] identified three combinations suffered from minimal or relatively-low ionosphere effects under different observational conditions.Geng and Bock [10] determined a new triple-frequency ionosphere-free (IF) carrier phase combination using the wide-lane carrier phase combinations with L1/L2 and L2/L5 signal pairs, helping implement rapid triple-frequency PPP AR.However, the above proposed combinations are all geometry-based, thus, an integer least-squares estimation is needed to resolve the combination ambiguities.In the meantime, if the geometry-free (GF) carrier phase combinations are used, it is difficult to guarantee the integer nature of the combination ambiguities, or it is impossible to determine the combination ambiguities by directly rounding because of the very large total noise level.
In the above discussions only the carrier phase measurements are adopted to form the linear combinations.If the code measurements are involved, a difference can be made.The best-known code-phase combination should be the Hatch-Melbourne-Wübbena (HMW) combination [19][20][21] for dual-frequency signals, which is expressed as: where P and L represent the code and carrier phase in meters, f is the frequency, the subscript numbers denote the frequency numbers, N W L = N 1 − N 2 is the WL ambiguity, and λ W L is the corresponding wavelength.From Equation (1) we can see that the HMW combination is free from geometry, the ionosphere, and the troposphere, and is only affected by the noise of the code and carrier phase measurements.Moreover, the WL ambiguity can be easily determined by rounding the average after several epochs, whose efficiency in PPP-AR was discussed in detail by Geng et al. [22].This combination was extended to the triple-frequency case and generated the "extra-wide-lane" (EWL) combination with L2 and L5 for GPS [10], as well as B2 and B3 for BDS, whose wavelengths can reach up to 5.861 m and 4.884 m, respectively.They were used in the first step of triple-frequency carrier phase ambiguity resolution (TCAR) [15].
Except for the HMW combination, many researchers tried to find other code-phase combinations with similar characteristics.Henkel and Günther [23] recommended a group of code-phase combinations with a wavelength of several meters and a noise level of a few centimeters, to maximize the ambiguity discrimination and ensure the reliability.Wang and Rothacher [24] presented a simplified method to minimize the noise level for the triple-frequency GF and ionosphere-free (GIF) linear combinations with both code and carrier phase measurements.However, no validation was made with real BDS data in the above two papers.Deo et al. [25] gives the triple-frequency ionosphere-free combinations for code and phase observations separately, but only the dual-frequency signals were used in the mixed code-phase combinations.Zhao et al. [26] introduced different criteria to select the optimal combined signals for GF and geometry-based (GB) TCAR, respectively.For GF TCAR the sum of the ionosphere scale factor (ISF), the measurement noise, the wavelength to total noise ratio, and the success rate are all considered, and the coefficients of the code combination in code-EWL and code-WL signal pairs were restricted to be integer numbers but, in fact, this is needless.
Since the second and third steps of general TCAR method is less reliable, Xu et al. [8] proposed a new GIF combination with P1, P2, and the ambiguity-resolved EWL observation to determine the WL ambiguity, and Zhao et al. [27] employed all three original code measurements and the ambiguity-resolved carrier phase combinations in previous steps to round off the WL ambiguity and the original ambiguity respectively.Inevitably, both the above proposed code-phase combinations rely on the HMW combination.Besides, for cycle slip detection, researchers paid more attention on the coefficients of the carrier phase part, and only P1 or a simple average of P1, P2 and P3 was used in the code part of the GF code-phase combination [3,28].
Table 1 gives an overview of some important works related to the optimal combinations, along with their characteristics and limitations.No employed parameters are listed since in some contributions there are a series of combinations to select for a specific system.Generally, all researchers proposed the EWL combination ϕ (0,−1,1) as the optimal one thus more attention should be paid on the code coefficients in the mixed code-phase combinations.A mixed code-phase combination may be preferred because the geometry-free model can be established.In this situation the ambiguity can be properly determined by directly rounding, rather than being searched in the geometry-based model.On the other hand, the HMW combination on B2 and B3 signals can reliably fix the EWL ambiguity within a few epochs, and in other proposed code-phase combinations the EWL combination should be used.From this perspective, for a mixed code-phase combination, a comparison with the HMW combination is necessary and sufficient to evaluate the performance of this new combination.
Table 1.Overview of several contributions to the optimal combinations for reliable AR.

Reference
Characteristics Limitation [12][13][14][15][16] (1) Theoretically discussed the factors that effects the performance of combinations (2) Gave the general selection procedure and specific useful combinations (1) Only carrier phase is used (2) The ambiguity should be searched due to the used geometry-based model [10,17,18] select combinations that reduce or eliminate the effect of ionosphere delay [23,24] (1) Mixed code-phase combinations (2) Minimize the noise level of the combinations and fix ambiguity by rounding with high success rate No BDS data are validated [8,26,27] (1) Mixed code-phase combinations (2) Fix ambiguity by rounding WL ambiguity fixing relies on the resolved EWL ambiguity [3,28] (1) Mixed code-phase combinations (2) Fix ambiguity by rounding no specific discussion of the code coefficients in the combination From the above discussion we can see that, currently, the HMW combination is still an efficient approach for AR with both dual-frequency and triple-frequency cases, and it can be applied in cycle-slip detection and AR with high reliability.However, for the BDS system, attention should be paid to two special problems: (1) The noise levels of three original code measurements are no longer identical [29], thus, the calculated total noise should be modified; (2) There are satellite-induced code variations in BDS system which can up to several meters [30].
These systematic biases have a significant impact on AR using HMW combinations, as thoroughly discussed by Zhang et al. [31] and a model correction is required before data processing [32].
For these reasons we propose an optimal code-phase combination determination method, considering all the factors that may affect the AR reliability for different situations.These factors include the combined ionosphere delay, the wavelength, and the noise level of the code and carrier phase measurements.For easy implementation, a sophisticated traversal searching approach is employed, in which the range of the phase coefficients and the combined ISF is limited, and three conditions are restricted with given a priori values of the ionosphere delay on B1 and noise levels on code and carrier phase measurements.
This paper is organized as follows: we first briefly present the general functional model of triple-frequency code-phase combinations and then give the GF code-phase combination.The two-step method to determine the optimal triple-frequency GF code-phase combinations is introduced afterwards.The coefficients of the determined optimal combinations are presented in the first subsection of the Results section.Then the collected real triple-frequency BDS data are described and the ambiguity residuals and their standard deviation (STD) for the determined optimal, suboptimal combinations and the HMW combination are calculated and compared.Followed by this, we discuss the results sufficiently, and then draw some conclusions in the last section, including better performance of the determined optimal and suboptimal combinations than the HMW combinations.

General Functional Model
Given the original code and carrier phase measurements in meters with frequency i, i.e., P i and L i , the corresponding observation equations can be expressed as: where the superscript s and the subscript r represent the satellite and the receiver, and ρ denotes the distance between them.dT, dI, and dO are the troposphere delay, the ionosphere delay, and the orbital error, respectively, c is the speed of light in vacuum.dt r and dt s represent the clock errors of the receiver and the satellite, dH and dh are the hardware delays of the code and the carrier phase.N and λ are the integer ambiguity and the corresponding wavelength, and ε is the measurement noise.Generally, the linear combination from the code and carrier phase measurements of triple-frequency signals can be formulated as follows: where the subscripts of the measurements denote the frequency, a i and b i (i = 1, 2, 3) are the corresponding coefficients.Derived from Equations ( 2) and (3), the observation equation of the code-phase combination can be expressed as: where G is the overall geometric part and written as: dH (a 1 ,a 2 ,a 3 ) and dh (b 1 ,b 2 ,b 3 ) are the combination of hardware delay for code and carrier phase respectively, dI 1 is the ionosphere delay on B1, and β is the ISF defined by Li et al. [7] and expressed as: The geometric part in the first term on the right side of Equation ( 5) can be scaled by its coefficients.A GF combination can be obtained if the Equation ( 8) holds true, and a geometry-preserving one if Equation ( 9) is favorable:

GF Code-Phase Combination
To simplify the observation of Equation ( 5) and ignore the effect of the hardware delays, the double-differenced measurements are used in this paper.Furthermore, since the optimal code-phase combination is used for fast AR, the condition Equation ( 8) for GF combinations are satisfied.Therefore, we can get the observation equation of the double-differenced GF code-phase combination as follows: where ∆ is the double-differencing operator.
From Equation (4) we can see that the triple-frequency GF code-phase combination can be divided into two parts: the code combination part and the phase combination part.In order to simplify the discussion, in each part the geometry-based mode is adopted, i.e., the code and phase coefficients should satisfy the following condition: Then, according to Feng [18], the geometry-based phase part can be reformulated and Equation ( 4) is thereby rewritten as: where i, j, and k are the integer coefficients.The double-differenced observation Equation ( 10) is, therefore, modified to: where is the combined ISF of the code-phase combination and is the ISF of the phase part.
If the combined ionosphere residual, together with the measurement noise of the code and carrier phase, is regarded as the uncertainty of the integer ambiguity, the GF code-phase combination makes it possible to determine the integer ambiguity by directly rounding with the following formula: where [*] is the rounding operator.Provided the difference between the calculated float ambiguity and its true value follows the normal distribution centered at zero, the AR success rate P s equals the probability of rounding the float ambiguity to its true value.This theoretical value is computed in Appendix A.

A Two-Step Method
The optimal GF code-phase combination could achieve the largest AR success rate.The wavelength, the combined ISF and phase noise factor (PNF) of code and carrier phase all depend on the combination coefficients, thus from Equations (A1) and (A2) we can see that the AR success rate depends only on the combination coefficients, the measurement noise, and the ionosphere residual on B1.The last two factors could be set as a prior values, so the vital issue is to determine the real-valued coefficients (a 1 , a 2 , a 3 ) and the integer coefficients (i, j, k).In this section, a two-step method is carried out, different combination coefficients and their corresponding success rates are calculated, and the optimal combinations are thereby screened out with the largest success rate.
Considering the coefficients (i, j, k) are integers and the combined ISF β 0 should be smaller than 1 to reduce the effect of ionosphere delays, a traversal searching method is adopted in the first step.The ranges for i, j, k are all limited from −50 to 50, and the value of β 0 is taken over the range [−1.0, 1.0] with an interval of 0.01.The real-valued coefficients (a 1 , a 2 , a 3 ) can be determined for a certain set of values (i, j, k, β 0 ) using the following three conditions: where the first condition means the code part of the combination is geometry-based, the second condition defines the combined ionosphere residual, and the last condition denotes the noise level of the code part should be minimized.Absolutely, when ∆dI 1 = 0 the second condition is redundant.In this situation the coefficients (a 1 , a 2 , a 3 ) are determined as an extreme value of the binary function problem prior to the combined ISF, and values of 1+2n 2 , n 2 1+2n 2 , 1 1+2n 2 are obtained.Additionally, when β 0 = 0 the code-phase combination turns ionosphere-free.
In the second step, the AR success rate for each triple-frequency code-phase combination determined in the previous step is calculated in Appendix A under the given ionosphere residual on B1 and the measurement noise of the code and carrier phases.By sorting all the combinations whose phase coefficients are linear-independent, we define the combinations with the largest and second largest AR success rate as the optimal and suboptimal code-phase combinations, respectively.

Optimal and Suboptimal Combinations
Taking the current BDS as an example and providing the measurement noise of the original code and carrier phase as 0.3 m and 0.003 m, respectively, we give the selected optimal and suboptimal code-phase combinations under different ionosphere residuals on B1, as shown in Tables 2 and 3.The assumptions that ∆dI 1 = 0, 0.1, 0.5, and 1.0 m represent four typical cases, i.e., zero or short baselines when the ionosphere residual can be ignored, medium baselines when station's distance is between 20 km and 50 km, medium-long baselines when station's distance is between 50 km and 100 km, and long baselines when the station's distance is larger than 100 km.Along with the combination coefficients are the corresponding standard deviations and AR success rates.Note that the effect of code variations for BDS has not been considered in the proposed combinations, thus a previous correction of BDS code variations is necessary.Meanwhile, the proposed method can be used not only for the current BDS system, but also for other systems with multi-frequency signals, such as the new-generation BDS with B1I, B1C (1575.42MHz), and B2a (1176.45MHz) signals [33,34], modernized GPS with L1, L2, and L5 signals, or Galileo with E1, E5 and E6 signals.
From Table 2 we can see that the phase coefficients of the optimal code-phase combination remains (0, −1, 1) for all cases because of its long wavelength of 4.8842 m.The theoretical success rates of rounding can always be up to 100.0%, though the code coefficients may be changed in different cases.On the other hand, there is a trend that the combined ISF of the optimal combination gradually decrease to zero along with the increase of the ionosphere residual on B1.We can infer that the optimal code-phase combination will become the GIF combination for long baselines with very large ionosphere residuals on B1.In fact, the GIF code-phase combination that (i, j, k) = (0, −1, 1) and (a 1 , a 2 , a 3 ) = (−0.1217,0.0912, 1.0305) can also achieve nearly a 100.0%success rate no matter how large the ionosphere residual on B1 is.
Different from the optimal code-phase combination, the phase coefficients of the suboptimal combination (i, j, k) may be changed when the ionosphere residual on B1 becomes larger to 1.0 m.The relative shorter wavelength can lead to smaller measurement noise of phase combination, and the ionosphere residual of the phase part can be reduced by the code counterpart.The theoretical success rate of rounding when the ionosphere residual on B1 is smaller than 1.0 m can remain over 95%.
It should be noted that the values of (a 1 , a 2 , a 3 ) in Tables 2 and 3 cannot be directly used since they are still not precise enough and will cause considerable computational errors.Therefore, more precise values of code coefficients calculated from Equation (15) with integer coefficients (i, j, k) and combined ISF β 0 in Tables 2 and 3 are necessary for practical data processing.

Data Description
In order to validate the proposed algorithm and demonstrate the performance of the determined combinations, triple-frequency BDS data with 1 Hz sampling rate on day of year (DOY) 339 of 2015 is used, and the satellite cut-off elevation angle is set to 15 • .The data were collected from ten stations in Guangdong Continuous Operational Reference System (GDCORS), whose distributions are plotted in Figure 1.The UR370 receivers manufactured by Unicore Communications, Inc. (Beijing, China) are used at stations ZHGT and QIAO, and the PD318 receivers manufactured by the Panda Space and Time Technology Co., Ltd.(Wuhan, China) are used at another eight stations.Several baselines can be formed: the baseline from ZHGT to QIAO (20 km), the baseline from GMGT to ZSGT (66 km), and the baseline from YFGT to ZSGT (147 km) are selected as examples of medium, medium-long, and long baselines for specific discussion.The values of ∆dI 1 for the three baselines are assumed to be 0.1 m, 0.5 m, and 1.0 m, respectively, thus, the optimal and suboptimal combinations are selected as shown in Table 4. for the three baselines are assumed to be 0.1 m, 0.5 m, and 1.0 m, respectively, thus, the optimal and suboptimal combinations are selected as shown in Table 4.

Optimal Combination
The ambiguity residuals of the optimal combination, the GIF code-phase combination, and the HMW combination with B2 and B3 signals for three baselines are shown from Figures 2-4.In each figure the blue dot represents the ambiguity residual, and the orange one represents the satellite elevation angle.The optimal combination (marked as optimal), the GIF code-phase combination (marked as GIF), and the HMW combination with B2 and B3 signals (marked as HMW23) are shown in the left panels, the middle panels, and the right panels, respectively.The corresponding STD values for each combination are also shown in each panel.

Optimal Combination
The ambiguity residuals of the optimal combination, the GIF code-phase combination, and the HMW combination with B2 and B3 signals for three baselines are shown from Figures 2-4.In each figure the blue dot represents the ambiguity residual, and the orange one represents the satellite elevation angle.The optimal combination (marked as optimal), the GIF code-phase combination (marked as GIF), and the HMW combination with B2 and B3 signals (marked as HMW23) are shown in the left panels, the middle panels, and the right panels, respectively.The corresponding STD values for each combination are also shown in each panel.From the three figures we can see that for the medium baseline the optimal combination can achieve better results than those by the HMW23 combination, especially when the satellite elevation angle is larger than 45°.For the medium-long and long baselines the performance of the optimal combination becomes remarkably better than that of the HMW23 combination.Specifically, an improvement of 24.5% and 15% are achieved for satellite C10.
Furthermore, the STD values for the three combinations of all IGSO and MEO satellites, the improvement percentages of the optimal combination compared with the HMW23 combination ( ), and the STD differences between the optimal and the GIF code-phase combinations ( ), were calculated and are shown in Table 5.Generally, the results in the table coincide with the results for satellites C10 and C14.An overall improvement of optimal combination compared with the HMW23 combination is achieved, with almost 100% AR success rate.In fact, from Figures 2-4 we can see that the ambiguity residuals of the optimal combination at high satellite elevation angles are much smaller than those of the HMW23 combination.Additionally, as the baseline length increases, the STD difference becomes smaller and smaller.This indicates that the GIF combination is more suitable for long baselines.

Suboptimal Combination
Since the EWL ambiguities can be reliably fixed to their right values, we can use them to verify the correctness of the determined integer ambiguity of the suboptimal combination, the HMW combination with B1 and B2 signals, and with B1 and B3 signals.The corresponding ambiguity residuals of each combination for three baselines are shown from Figures 5-7, together with the corresponding STD values.In these three figures the blue dots and orange dots also represent the ambiguity residual and the satellite elevation angle, respectively.The suboptimal combination From the three figures we can see that for the medium baseline the optimal combination can achieve better results than those by the HMW23 combination, especially when the satellite elevation angle is larger than 45 • .For the medium-long and long baselines the performance of the optimal combination becomes remarkably better than that of the HMW23 combination.Specifically, an improvement of 24.5% and 15% are achieved for satellite C10.
Furthermore, the STD values for the three combinations of all IGSO and MEO satellites, the improvement percentages of the optimal combination compared with the HMW23 combination (P ip ), and the STD differences between the optimal and the GIF code-phase combinations (δ), were calculated and are shown in Table 5.Generally, the results in the table coincide with the results for satellites C10 and C14.An overall improvement of optimal combination compared with the HMW23 combination is achieved, with almost 100% AR success rate.In fact, from Figures 2-4 we can see that the ambiguity residuals of the optimal combination at high satellite elevation angles are much smaller than those of the HMW23 combination.Additionally, as the baseline length increases, the STD difference δ becomes smaller and smaller.This indicates that the GIF combination is more suitable for long baselines.

Suboptimal Combination
Since the EWL ambiguities can be reliably fixed to their right values, we can use them to verify the correctness of the determined integer ambiguity of the suboptimal combination, the HMW combination with B1 and B2 signals, and with B1 and B3 signals.The corresponding ambiguity residuals of each combination for three baselines are shown from Figures 5-7, together with the corresponding STD values.In these three figures the blue dots and orange dots also represent the ambiguity residual and the satellite elevation angle, respectively.The suboptimal combination (marked as suboptimal), the HMW combination with B1 and B2 signals (marked as HMW12), and the HMW combination with B1 and B3 signals (marked as HMW13) are shown in the left panels, the middle panels, and the right panels, respectively.All the selected suboptimal combinations in the three baselines are WL combinations whose wavelength is between 0.75 and 2.93 m [7], so the comparison between the suboptimal combination and the HMW12 or the HMW13 combination is reasonable.From the above three figures we can see that the suboptimal combination achieves the best results in all baselines with a dramatic improvement.The STD of ambiguity residuals for the suboptimal combination of C10 and C14 are only about 0.2 cycles, thus, the AR success rate by directly rounding can be up to 99%.Additionally, no apparent difference is observed between the IGSO and MEO satellites, except for the shortest baseline where the residual for C14 is amplified much larger than the residual for C10.This is because, for short baselines, the effect of the code measurement noise, rather than the ionosphere delay, plays a major role in determining the selected combinations.The code measurement noise suffers more at low satellite elevation angles, and there is a greater percentage of the observation data when C14 at a low satellite elevation angle than that of C10.As the baseline distance increases, the effect of the code measurement noise reduces, and the effect of the ionosphere delay enhances.
On the contrary, the HMW12 combination shows the worst performance, especially for satellite C10 at low satellite elevation angles.This indicates that the selected suboptimal combination has the lowest total noise level, and we can reliably fix the ambiguity of the suboptimal combination by rounding with a short time span, even in a single epoch.However, relatively more time is required to reliably fix the HMW12 or HMW13 ambiguity through a smoothing technique.
In addition, the STD values of the suboptimal, HMW12, and HMW13 combinations for all IGSO and MEO satellites, the AR success rate of the suboptimal combination Ps, and the improvement percentages of the suboptimal combination compared with the HMW12 combination and the HMW13 combination, which are marked as ∆ and ∆ , are calculated and summarized in Table 6.Generally speaking, the STD values of all three combinations for all IGSO and MEO satellites coincide with the results of C10 and C14 shown from Figures 5-7.An improvement of approximately 50% has been achieved compared with the HMW12 combination, and about 20% compared with the HMW13 combination.
What is interesting is that for the longest baseline the phase part of the selected suboptimal combination is identical to that of the HMW13 combination.However, the STD value of the suboptimal combination still outperforms the HMW13 combination by 22.078%.This further demonstrates the better performance of the selected combination using our method compared with the HMW combinations.All the selected suboptimal combinations in the three baselines are WL combinations whose wavelength is between 0.75 and 2.93 m [7], so the comparison between the suboptimal combination and the HMW12 or the HMW13 combination is reasonable.From the above three figures we can see that the suboptimal combination achieves the best results in all baselines with a dramatic improvement.The STD of ambiguity residuals for the suboptimal combination of C10 and C14 are only about 0.2 cycles, thus, the AR success rate by directly rounding can be up to 99%.Additionally, no apparent difference is observed between the IGSO and MEO satellites, except for the shortest baseline where the residual for C14 is amplified much larger than the residual for C10.This is because, for short baselines, the effect of the code measurement noise, rather than the ionosphere delay, plays a major role in determining the selected combinations.The code measurement noise suffers more at low satellite elevation angles, and there is a greater percentage of the observation data when C14 at a low satellite elevation angle than that of C10.As the baseline distance increases, the effect of the code measurement noise reduces, and the effect of the ionosphere delay enhances.
On the contrary, the HMW12 combination shows the worst performance, especially for satellite C10 at low satellite elevation angles.This indicates that the selected suboptimal combination has the lowest total noise level, and we can reliably fix the ambiguity of the suboptimal combination by rounding with a short time span, even in a single epoch.However, relatively more time is required to reliably fix the HMW12 or HMW13 ambiguity through a smoothing technique.
In addition, the STD values of the suboptimal, HMW12, and HMW13 combinations for all IGSO and MEO satellites, the AR success rate of the suboptimal combination Ps, and the improvement percentages of the suboptimal combination compared with the HMW12 combination and the HMW13 combination, which are marked as P ∆12 and P ∆13 , are calculated and summarized in Table 6.Generally speaking, the STD values of all three combinations for all IGSO and MEO satellites coincide with the results of C10 and C14 shown from Figures 5-7.An improvement of approximately 50% has been achieved compared with the HMW12 combination, and about 20% compared with the HMW13 combination.
What is interesting is that for the longest baseline the phase part of the selected suboptimal combination is identical to that of the HMW13 combination.However, the STD value of the suboptimal combination still outperforms the HMW13 combination by 22.078%.This further demonstrates the better performance of the selected combination using our method compared with the HMW combinations.

Discussion
Rapid and reliable AR is the prerequisite to provide precise positioning services, and the systematic biases, including orbit error, troposphere delay, and ionosphere delay make it difficult to fix the ambiguity if these biases are not removed or reduced.Since these biases can be divided into two classes, i.e., the frequency-independent geometric term (orbit error and troposphere delay) and the frequency-related term (ionosphere delay), the GF and IF combinations are used accordingly to cancel each term, thus, the GIF combination is preferred.In this paper, we proposed a two-step method, using both the triple-frequency code and carrier phase observations to select the optimal code-phase combination which is also geometry-and ionosphere-free.The selected combinations can instantaneously and reliably fix the ambiguity by simply rounding.From the results of the above section, the phase part of the optimal combination is the EWL, and this agrees with a large number of research results, both in theory [12,13,18] and in practice [27,28].The HMW combination with the B2 and B3 signals shows an excellent capability in speed and reliability of AR.However, from Table 5, we see the improvement using the optimal combination, though not much, in the way of further reducing the total noise level.
On the other hand, the reliability of WL AR using the suboptimal combination has a remarkable improvement.The STD of the ambiguity residual after rounding using the suboptimal combination is much smaller than those using HMW12 or using HMW13 combinations, with more than 20%.No matter the ambiguity N(1, 1, −2) for short or medium baselines, or the ambiguity N(1, 0, −1) for long baselines, once fixed, they can be used together with the fixed EWL ambiguity N(0, −1, 1) to determine an arbitrary WL ambiguity since they are linear-dependent [12,13].In addition, the ionosphere delay on B1 can be estimated using two WL combinations [35].In this way, decimeter-level kinematic positioning services can be provided, which are comparable to that of carrier phase-smoothed code differential positioning [36][37][38].
All the above test results, analysis and discussions focus on the performance of the proposed method on AR success rates, which are indicated by the STD values of the ambiguity residuals in Tables 5 and 6.If the STD values are much smaller, the corresponding combined ambiguity can be reliably fixed after a few seconds' smoothing, or even in a single epoch.Furthermore, the resolved ambiguity can be immediately used for positioning at the decimeter level (with resolved EWL or WL ambiguity) or at the centimeter level (with resolved original ambiguity).We emphasize how to select optimal combinations in order to achieve fast and reliable AR, since it is the prerequisite for precise positioning.How to achieve precise positioning with the resolved ambiguity is not an issue discussed in this paper.
There are three different admissible integer estimators in use, i.e., integer rounding, integer bootstrapping, and integer least-squares [39].Theoretically speaking, the integer least-squares method is the most rigorous and the least-squares ambiguity decorrelation adjustment (LAMBDA) method [40] is widely used.However, the search efficiency will dramatically decrease for high-dimensional integer ambiguity resolution [41], this problem may occur in the multi-frequency case.Even if, currently, the so-called "partial AR" strategy [42][43][44] arises to fix part of the ambiguities, how to select the subset of the ambiguity is another introduced problem.Additionally, for integer least-squares estimation, the procedure includes an integer ambiguity search and ambiguity validation, and the latter is a troublesome problem to be solved [45].If sufficiently reliable, the integer rounding could be an alternative approach, and the selected optimal and suboptimal combinations in this paper make this possible.
Compared with the current code-phase combinations in Table 1 and the HMW combination using signals on two arbitrary frequencies, in our method not only are the selected combinations free from geometry and insusceptible to the ionosphere delay (which is the same as those in most of the current combinations and the HMW combination), but the noise level is also minimized and the largest AR success rate of the selected combination by rounding is guaranteed.In different scenarios the coefficients of the optimal combinations may be different, but the AR performance remains the best, and the suboptimal combination remains independent of the optimal combination.From this point of view, our proposed method is a collection of both model-and data-driven selection strategies.The model-driven strategy is used to determine the optimal combination coefficients and achieve the largest AR success rate under a certain level of measurement noise and ionosphere delay, and the data-driven strategy is used to determine the proper a priori information of measurement noise and ionosphere delay for baselines with various lengths.The mathematical model and the used data make equivalent contributions to the optimal combination determination.For the BDS system, specifically, different noise levels on different signal frequencies are taken fully into account; thus, an improvement of around 10% can be achieved.
Finally, the overall performances of our method for several baselines with various lengths using data at all ten stations (plotted in Figure 1) are shown in Table 7. From the table we can find that in all baselines the optimal combination outperforms the HMW23 combination.Additionally, the suboptimal combination is also effective, whose STD values are much smaller than those of the HMW12 and HMW13 combinations, and the success rates by directly rounding are nearly 99%.

Conclusions
In this paper we proposed a method to determine the optimal triple-frequency code-phase combinations, in which the factors that affect the AR success rate were taken into account.A sophisticated traversal search was employed to calculate the coefficients of the combinations and select the optimal combinations under different ionosphere residuals.In order to validate the algorithm and demonstrate the performance of the determined combinations, triple-frequency BDS data were used and the standard deviations of ambiguity estimation for the optimal, the suboptimal, the geometry-ionosphere-free, and the HMW combination were calculated and compared.Results show that: 1.
Although the optimal combination and the HMW23 combination can both achieve almost 100% AR success rate by rounding, the STD values of the ambiguity residual for the optimal combination are slightly better than those of the HMW23 combination, especially when the satellite elevation angle is higher than 45 • .Specifically, an improvement of nearly 25% can be achieved for C10 satellite in the medium-long baseline.

2.
The suboptimal combination can achieve an AR success rate around 99% by rounding.Compared with the HMW12 and HMW13 combinations, the suboptimal combination achieves best results in all baselines, with a dramatic improvement of about 50% and 20%, respectively.3.
Using the correctly fixed optimal and suboptimal ambiguities, the WL ambiguity can then be instantaneously determined and the ionosphere delay on B1 can be estimated.Therefore, decimeter-level kinematic positioning service can be provided.
As discussed in this paper, the factors of combinations that affect the AR success rate include the combination coefficients, the measurement noise, and the ionosphere residual on B1.To achieve the largest AR success rate, the former is determined by the model and the last two are up to the real data.In the present method we built an explicit relationship between the combination coefficients and the AR success rate, but intensive research on the relationship between the data-driven factors and the AR success rate is insufficient.For example, the measurement noise is affected by the satellite elevation angle and the ionosphere residual is related to the baseline length.We only give a different prior ionosphere residual according to different baselines, and the satellite elevation angle is also an important factor, as we can see from Figures 2-7.Followed this discussion, further study may be conducted to obtain the relationship between precise data-driven factors, such as satellite elevation angles and the AR success rate, and to derive a more rigorous optimal combination.
A small modification should be made for the code measurements of BDS system since the noise of P3 is much smaller than that of P1 and P2 [29].As a result the PNF of the code part for BDS system is redefined as: (a 1 ,a 2 ,a 3 ) = a 2 1 + a 2 2 + (na 3 ) 2 (A4) where n is a factor smaller than 1, and a practical value of 0.2 is used [29].Additionally, the code variations of IGSO and MEO satellites are corrected using the models developed by Zou et al. [47], while those of GEO satellites are not since no effective correction models are available.
Remote Sens. 2018, 10, x FOR PEER REVIEW 8 of 17 China) are used at stations ZHGT and QIAO, and the PD318 receivers manufactured by the Panda Space and Time Technology Co., Ltd (Wuhan, China) are used at another eight stations.Several baselines can be formed: the baseline from ZHGT to QIAO (20 km), the baseline from GMGT to ZSGT (66 km), and the baseline from YFGT to ZSGT (147 km) are selected as examples of medium, mediumlong, and long baselines for specific discussion.The values of Δ

Figure 1 .
Figure 1.Distribution of the used stations in GDCORS and the selected baselines.

Figure 1 .
Figure 1.Distribution of the used stations in GDCORS and the selected baselines.

Table 4 .
The selected optimal and suboptimal combinations for three baselines in this experiment.2493, −0.0354, 0.7861)

Figure 2 .
Figure 2. Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the medium baseline (20 km).

Figure 3 .
Figure 3.Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the medium-long baseline (66 km).

Figure 2 .
Figure 2. Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the medium baseline (20 km).

Figure 2 .
Figure 2. Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the medium baseline (20 km).

Figure 3 .
Figure 3.Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the medium-long baseline (66 km).

Figure 3 .
Figure 3.Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the medium-long baseline (66 km).

Figure 4 .
Figure 4. Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the long baseline (147 km).

Figure 4 .
Figure 4. Comparison of ambiguity residuals and the corresponding standard deviations of optimal, GIF, and HMW23 combinations for the long baseline (147 km).
Remote Sens. 2018, 10, x FOR PEER REVIEW 11 of 17 (marked as suboptimal), the HMW combination with B1 and B2 signals (marked as HMW12), and the HMW combination with B1 and B3 signals (marked as HMW13) are shown in the left panels, the middle panels, and the right panels, respectively.

Figure 5 .
Figure 5.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (20 km).

Figure 6 .
Figure 6.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (66 km).

Figure 5 .
Figure 5.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (20 km).

Figure 5 .
Figure 5.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (20 km).

Figure 6 .
Figure 6.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (66 km).

Figure 6 .
Figure 6.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (66 km).

Figure 7 .
Figure 7.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (147 km).

Figure 7 .
Figure 7.Comparison of ambiguity residuals and the corresponding standard deviations of suboptimal, HMW12, and HMW13 combinations for the medium baseline (147 km).

Table 2 .
Optimal BDS triple-frequency code-phase combinations and the corresponding AR success rate for different code noises and ionosphere residuals.

Table 3 .
Suboptimal BDS triple-frequency code-phase combinations and the corresponding AR success rate for different code noises and ionosphere residuals.

Table 4 .
The selected optimal and suboptimal combinations for three baselines in this experiment.

Table 5 .
Performance comparison for all IGSO and MEO satellites among the optimal, the GIF, and the HMW23 combinations.

Table 5 .
Performance comparison for all IGSO and MEO satellites among the optimal, the GIF, and the HMW23 combinations.

Table 6 .
Performance comparison for all IGSO and MEO satellites among the suboptimal, the HMW12, and the HMW13 combinations.

Table 7 .
Performance comparison for baselines with various lengths among the optimal, the HMW23, the suboptimal, the HMW12, and the HMW13 combinations.