Abstract

As part of the GPS modernization, the GPS L1 C/A signal will be augmented by the L1C signal. With this improvement, for the first time, several signals from the same constellation will be available at the same frequency. In this paper, an acquisition method is implemented to combine the GPS L1 C/A and L1C signals before correlation. The combined acquisition succeeds to acquire the signal at low C/N0, whereas the acquisition of the L1C data channel alone fails. Concerning the tracking, a method to combine the GPS L1 C/A and L1C signals before the discriminator is developed. This method shows better performance than the traditional tracking using only one signal. Finally, a Kalman filter to combine the signals in the tracking is developed. It shows better performance than the traditional tracking in all conditions. Since the L1C signal will not be broadcast before 2013, these methods are tested using a software signal simulator.

1. Introduction

The current Global Positioning System (GPS) has been widely studied since its creation. Nevertheless, with the recent or upcoming availability of several new signals from GPS, Galileo, GLONASS and Compass, the beginning of a new era in the field of GNSS will soon arrive. With such an improvement in the number of signals available to the user and the geometry offered by the combination of the various constellations and signals, the current capabilities of Global Navigation Satellite System (GNSS) receivers will increase. At the same time, however, the receiver implementation will also become more complex. For instance, concerning GPS L1, many signals will be available in the next decade: L1 C/A (currently available) and L1C (divided into two channels). The focus of this paper is to combine the GPS L1 C/A and L1C signals at the acquisition and tracking level in order to improve the performance relative to the L1 C/A-only case.

In the literature, different methods of acquisition have been developed and their performance determined. These techniques can be divided in three main categories [1]: coherent and noncoherent integration [2, 3], differentially coherent integration [4], or a combination of these. In the case of coherent acquisition, the size of the frequency bins during acquisition needs to be decreased with increasing integration time [5]. Thus, in the case of weak signals (e.g., indoors) where longer coherent integration times are desired, the processing time increases dramatically. Another solution is to use a noncoherent correlation to acquire lower power signals. This way, the size of the frequency bin does not need to be changed but the overall performance is not as good in terms of signal-to-noise ratio (SNR) as with coherent integration. Concerning the combination of a noncoherent or a differential technique with a coherent integration, it permits to acquire signals attenuated by as much as 20 dB [6, 7].

The GPS L1C signal should be less affected by the tracking errors (thermal noise, multipath, interference, user dynamics, etc.) than GPS L1 C/A due to the length of the spreading codes and the modulations [810]. Indeed, in the case of a BOC (Binary Offset Carrier) modulation, the autocorrelation function contains a main peak that is sharper than for a BPSK (Binary Phase Shift Keying) modulation and thus reduces the noise [9]. However, the BOC side peaks can create an ambiguity that can be wrongly tracked, leading to degraded tracking performance, and ultimately biased positioning. Some research has been done on code tracking discriminators to avoid the tracking of these side peaks [11].

Unlike the studies mentioned above, the objective of this paper is to improve signal acquisition and tracking performance of the GPS L1 signal by acquiring and tracking several signal modulations (L1 C/A, L1C data and L1C pilot) simultaneously. The combined method of acquisition takes into account the relative powers of the signals in order to combine them. Two tracking methods are also developed. The first one combines the signals in a traditional scheme using their relative transmitted powers. The second one is a single-signal Kalman filter modified to use all available signals.

The contributions of this paper are threefold. First, it presents new methods of acquiring and tracking all three of the GPS L1 signals. For tracking, a preliminary multipath error assessment is also included. Second, although the focus is on GPS L1 signals, the methods, with appropriate modifications, are applicable to any GNSS frequency containing multiple signals. Third, some of the already existing methods present in the literature can be used in combination with the acquisition and first tracking method introduced here.

The paper begins with an overview of the GPS L1 C/A and L1C signals and their main properties. Then, the method used to test the algorithms is briefly described. Next, the combined acquisition method developed is presented and some results are shown. Afterwards, the combined tracking method developed is described and results in various conditions are shown. Finally the combined Kalman filter is introduced and tests in various conditions are shown.

2. Description of the GPS L1 C/A and L1C Signals

The structure and properties of the L1 C/A signals are not presented here since they can be readily found in the literature. For example, the reader can refer to Van Dierendonck [5] or Kaplan [12]. Unless otherwise specified, all the descriptions below concerning the structure of the GPS L1C signal are from IS-GPS-800 [13].

The GPS L1C signal is transmitted at the same central frequency as GPS L1 C/A (1575.42 MHz), with a nominal minimum received power 1.5 dB above L1 C/A, namely, −157 dBW. It is composed of two channels: a pilot (dataless) channel (denoted L1CP) and a data channel (denoted L1CD) transmitted in phase quadrature. The pilot channel includes a spreading code and an overlay, or secondary code, denoted as L1CO. The overlay code is generated using LFSR (Linear Feedback Shift Register). It is transmitted at 100 bits/s and contains 1800 bits, thus lasting 18 seconds. The overlay code is unique for each PRN (Pseudo-Random Noise). The data channel includes a spreading code and a navigation message. The spreading codes of the pilot and data channels are time synchronized. The spreading codes are broadcast at the same chipping rate as L1 C/A, that is, 1.023 Mchips/s. On both the pilot and data channels, the spreading codes have a period of 10 ms, and therefore contain 10,230 chips. The L1C spreading codes are generated using a modified Weil code (generated from a Legendre sequence). The generation is not reproduced in detail here.

The modulation on L1CD (data channel) is a BOC(1,1) modulation (Binary Offset Carrier). A BOC signal is created by modulating a sine wave carrier with the product of the PRN spreading code and the sign of the sinus of the subcarrier. The first parameter of the BOC is the ratio between the subcarrier frequency and the reference frequency, and the second parameter is the ratio between the spreading code frequency and the reference frequency. Due to this modulation, the spectrum of the signal is split in two symmetrical components with no remaining power at the central frequency. The two main lobes are centered at +/− 1.023 MHz relative to the central frequency, thus a bandwidth of 4.092 MHz is needed to keep most of the power. The modulation on L1CP (pilot channel) is a TMBOC (Time Multiplexed BOC), which consists of 29/33 chips modulated using a BOC(1,1) and others using a BOC(6,1). When BOC(6,1) is used, the number of chips (and correspondingly the chipping rate) is increased by a factor of six compared to the BOC(1,1) for a time interval equivalent to one chip of the original ranging code. This increased chipping rate has the effect of requiring an increased sampling frequency to extract all the information (because the minimum necessary bandwidth is now 14.322 MHz).

The power spectrum of the MBOC and BOC(1,1) modulations is shown in Figure 1.

The L1C navigation message is transmitted on the data channel at 100 symbols/s. The message is broadcast at a higher chipping rate than on L1 C/A (50 bits/s). Indeed, the navigation message of L1 C/A was not broadcast at a higher rate in order to have enough redundancy to guarantee a certain level of integrity of the extracted data bits. However, the L1C navigation message is transmitted using Cyclic Redundancy Check (CRC) and Forward Error Correction (FEC) techniques which improve the integrity even at the higher bit rate.

3. Description of the Simulations

Since the GPS L1C signal is not yet broadcast, a software signal simulator was developed in order to test the acquisition and tracking algorithms described here. The software signal simulator is a C language-based program and is an enhanced version of the IFG2Sim software (Intermediate Frequency GPS and Galileo Simulator) developed by Macchi et al. [14] in the PLAN Group at the University of Calgary. This software is able to generate samples of the L1 frequency at a specifiable intermediate frequency (IF). The software generates IF samples for static or mobile users using one bit quantization. The software simulates the signal at the output of a front-end taking into account only white noise. Possible interference sources and satellite clock errors, multipath, orbital errors, receiver clock errors, and atmospheric errors are not simulated. An overview of the signal generation portion of the software is presented in Figure 2. First, an initialization phase is realized (read configuration parameters), then the positions of the satellites are computed; the ones visible are determined and their position is approximated using cubic spline coefficients over an interval. If the satellite just appeared in view, the navigation message is generated. All the signals from visible satellites are then generated and added together. Finally, Gaussian noise is generated and added to the signal.

4. Combined Acquisition Method

In this section, a new method of acquisition is developed to combine the GPS L1 C/A and L1C (data and pilot channels) signals for the acquisition. Note that a similar method can be applied to Galileo L1-B and L1-C channels as shown in Mattos [15]. First, the method is explained in detail. Then, some results are shown using a signal generated by the software signal simulator described in the previous section.

4.1. Presentation of the Method

The combined acquisition algorithm proposed herein is shown in Figure 3.

Most of the scheme in Figure 3 is similar to a classic parallel code phase search acquisition scheme. However, there is an important difference in the generation of the local replica of the ranging code. Specifically, instead of generating only a single ranging code, a composite code consisting of the L1 C/A, L1CD and L1CP ranging codes is generated. The acquisition is realized as follows.(1)A carrier wave is generated at a fixed Doppler frequency 𝑓D using a local oscillator.(2)The incoming signal is multiplied by the carrier wave generated in step (1) to generate the in phase component 𝐼.(3)The incoming signal is multiplied by the carrier wave generated in step (1) shifted by 90° to generate the quadrature component 𝑄.(4)The Fourier transform of the combination of the 𝐼 and 𝑄 components is computed: 𝐹(𝐼+𝑗𝑄) (where 𝐹 represents the Fourier transform).(5)The local replica of the three ranging codes (GPS L1 C/A, L1C pilot, and L1C data channels) are generated.(6)The three ranging codes generated in step (5) are multiplied, respectively, by 𝐴, 𝑗𝐵, and 𝐶, and then added together. The values of the scalars 𝐴, 𝐵, and 𝐶 are computed according to the relative power of the GPS L1 C/A and L1C signals; details are shown later in this section. The L1C pilot channel is multiplied by 𝑗 because the pilot channel is transmitted with a 90° shift compared to the L1 other channels.(7)The Fourier transform of the combination created in step (6) is computed and its complex conjugate is multiplied by the Fourier transform computed in step (4). (8)The inverse Fourier transform of the vector computed in step (7) is computed and the results are squared.(9)Steps (1) through (8) are repeated for all Doppler shifts until the signal is found. Note that the whole process can be repeated for all Doppler shifts, or the frequency shift property of the FFT can be used to reduce the processing time. The result would be identical.

Herein, a frequency step of 50 Hz has been chosen due to the length of L1C (10 ms). Indeed, the minimum size of the frequency step is the inverse of the double of the integration time.

The 𝐴, 𝐵, and 𝐶 coefficients above are computed using the relative powers of the transmitted signals as follows. If the correlations of the three signals were computed separately, before squaring the results (step (8) above), the distributions of the three correlator outputs in the presence of signal would be [16](i)𝑋C/A𝑁(𝜇C/A,𝜎2C/A) for the distribution of L1C/A𝑋C/A,(ii)𝑋L1CD𝑁(𝜇L1CD,𝜎2L1CD) for the distribution of L1C data channel 𝑋L1CD,(iii)𝑋L1CP𝑁(𝜇L1CP,𝜎2L1CP) for the distribution of L1D pilot channel 𝑋L1CP,

where𝑁(𝜇,𝜎2) represents a Gaussian distribution with mean 𝜇 and variance 𝜎2, and 𝜇C/A, 𝜇L1CD, and 𝜇L1CP are the amplitudes of the correlation peak of the L1 C/A, and L1C data and pilot channels, respectively.

𝜎2C/A, 𝜎2L1CD, and 𝜎2L1CP are the variances of the noise associated with the L1 C/A, and L1C data, and pilot channels, respectively. They are computed using the relative broadcast power of these signals.

Therefore, using the relationship of noise/variance with power, the optimal combination—in a minimum variance sense—for the three channels is 𝑋L1=𝜎2C/A𝜎2C/A+𝜎2L1C×𝜎2L1CD𝜎2L1CD+𝜎2L1CP𝑋L1CP+𝜎12L1CD𝜎2L1CD+𝜎2L1CP𝑋L1CD+𝜎12C/A𝜎2C/A+𝜎2L1C𝑋C/A,𝜎𝐴=12C/A𝜎2C/A+𝜎2L1C𝜎=0.4142,𝐵=2C/A𝜎2C/A+𝜎2L1C𝜎2L1CD𝜎2L1CD+𝜎2L1CP𝜎=0.4388,𝐶=2C/A𝜎2C/A+𝜎2L1C𝜎12L1CD𝜎2L1CD+𝜎2L1CP=0.1470.(1)

Note that the detail of this derivation can be found in Macchi [17].

Finally, the acquisition presented in Figure 3 needs to be performed in parallel four times to account for all possible relative signs between the L1 C/A and L1C pilot and data channels. Specifically, due to the different navigation messages on the L1 C/A and L1C data channels, and the secondary code on the L1C pilot channel, it is not possible to know (without aiding information) the relative sign between these signals during acquisition. Therefore, a sign is applied in front of the 𝐴, 𝐵, and 𝐶 scalars to compute all the possible relative signs. The following possibilities are tested:(i)𝐴+𝑗𝐵+𝐶 called Comb1, (ii)𝐴𝑗𝐵𝐶 called Comb2, (iii)𝐴+𝑗𝐵𝐶 called Comb3, (iv)𝐴𝑗𝐵+𝐶 called Comb4.

These combinations can also be computed at the correlator outputs level after correlating the three codes separately. This last method would require less processing time but use more memory. The choice of the method therefore depends on the platform and whichever parameter is the most critical: processing time or memory.

The combination can be realized on the correlator outputs and in this case just requires 4 more operations. Since the results are squared at the end of the acquisition, the absolute sign is not important, only the relative sign. The maximum of the four combinations is selected.

4.2. Results Using the Software Signal Simulator

In order to test the method presented in the previous section, GPS L1 C/A and L1C data are simulated using the software signal simulator described previously. The data was processed with a software receiver that used the above-mentioned algorithms.

Two different sets of data are simulated with the parameters shown in Table 1. Due to a possible change in the sign of the data bit, or the secondary bit each ranging code period in the case of L1C, the correlation is done using 20 ms of incoming data using the zero padding technique [18].

The correlator outputs from the combined acquisition using the values of the first simulation are shown in Figure 4. Note that the curves have been slightly shifted (of 10 chips for each combination) for a better view of the various peaks. During the first 10 ms, 10 peaks are visible for all combinations. These 10 peaks correspond to the correlation of GPS C/A each millisecond. The 4th peak (around the 3rd ms) is much higher than the other peaks for two of the combinations. One of these combinations corresponds to the acquisition of the three channels with the right sign. This is determined by the highest peak (blue curve). The other combination (green curve) corresponds to the combination with the correct sign of L1 C/A with the pilot channel of L1C but a destructive combination (due to the wrong sign) with the data channel of L1C. Indeed, the broadcast power of the data channel of L1C is much smaller than for L1 C/A and the pilot channel of L1C. Therefore, a high peak is visible, but still smaller than when the three channels are combined with the correct sign. Concerning the second 10 ms of data, there is again a peak every 10 ms for all the combinations corresponding to the correlation of L1 C/A. One of the peaks is higher for one of the combinations (black curve) and this peak is located 10 ms after the highest peak mentioned above. This peak is smaller than the main peak of the first 10 ms due to a destructive combination due to a change in the data bit sign of L1C.

The corresponding result of the combined acquisition using the values of the second simulation is shown in Figure 5. As with the first simulation, during the first 10 ms, 10 peaks corresponding to the correlation of GPS C/A are visible. However, the ratio between the peaks and the noise level is much smaller than in the previous case due to the lower C/N0. The 4th peak (around the 3rd ms) is again much higher than the other peaks for two of the combinations. Overall, the same conclusions as above can be drawn, only the ratio between the peaks and the noise level changes. For comparison, Figure 6 shows the acquisition of the L1C data channel alone using the same data (second simulation). The L1C data channel is chosen for this comparison because it is the signal transmitted with the lowest power. This acquisition is done using a zero padding technique with 20 ms of data. It is clearly visible that the acquisition of the L1C data channel alone fails whereas it succeeds in the case of the combined acquisition. However, the secondary peaks in the combined acquisition (appearing C/A ranging code period) can become an issue at low C/N0. Particularly, the probability of false alarm increases due to these peaks.

5. Combined Tracking Method

In this section, a new method of tracking is developed to combine the GPS L1 C/A and L1C (data and pilot channels) signals. Note that a similar method can be applied to Galileo L1-B and L1-C channels. First, the method is explained in detail; then some results are shown using a signal generated by the software signal simulator described earlier.

5.1. Method of Tracking

The combined tracking method developed is shown in Figure 7 and consists of the following steps.(1)Generate the local carrier (corresponding to the intermediate frequency plus Doppler).(2)Remove the carrier from the incoming signal. Note that a shift of 𝜋/2 is applied to the local carrier for the L1C pilot channel because the pilot channel is transmitted in quadrature.(3)Separately correlate the outputs from step (2) with the L1 C/A, L1C data and pilot channels ranging codes. This correlation is computed for the prompt, early, and late correlators for the inphase 𝐼 and in quadrature 𝑄 components. These correlations are computed over 10 ms for each ranging code.(4)The correlator outputs are then multiplied by the scalars 𝐴, 𝐵, and 𝐶. These scalars have been computed in the section on acquisition (1).(5)Compute the following combinations (same as in the combined acquisition method except the “𝑗” is removed because the local carrier is already shifted by 𝜋/2 for the pilot channel):(a)𝐴+𝐵+𝐶 called Comb1, (b)𝐴𝐵𝐶 called Comb2, (c)𝐴+𝐵𝐶 called Comb3, (d)𝐴𝐵+𝐶 called Comb4. (6)Take the maximum of the four combinations in step (5).(7)Use the usual discriminators and filters for the DLL (Delay Lock Loop) and PLL (Phase Lock Loop).

The combination proposed—based on the minimum variance—is optimal. A similar combination can be done at the output of the discriminators. However, the output of the discriminators is no longer Gaussian, especially at low C/N0. Therefore, the combination will no longer be optimal.

5.2. Results Using the Software Signal Simulator

In order to test the method presented in the previous section, three different sets of data are simulated with the parameters shown in Table 2 and processed with the software receiver. The tracking parameters used in the receiver are shown in Table 3. Concerning the third simulation, the motion is as follows:(i)static for 5 seconds,(ii)acceleration of 10 m/s2 for 10 seconds,(iii)constant speed for 5 seconds,(iv)deceleration of 10 m/s2 for 5 seconds,(v)acceleration of 10 m/s2 for 5 seconds.

The simulated signals are tracked using the combined tracking method presented in the previous section as well as using a traditional tracking method able to track the GPS L1C data channel alone (lowest power of the three channels). One metric used to assess the performance of the proposed method relative to the traditional tracking method is the Phase Lock Indicator (PLI). The PLI used in this work is given in Van Dierendonck [5] and ranges from −1 (worst case) to 1 (perfect phase tracking).

The results of the static tests (first and second simulations) are shown in Figures 8 and 9, respectively. Note that the results of the first simulation for the combined method are reproduced in both figures for comparison purposes. As shown by the Doppler measurements and the PLI in these figures, the combined tracking is able to track the signals at 45 dB-Hz. Tracking of only the L1C data channel is still possible at 45 dB-Hz but its performance is not as good as that of the combined tracking. Specifically, the Doppler measurements are noisier and the PLI is not as close to unity. These results are expected given that only a single signal is used. Nevertheless, it serves to illustrate one of the benefits of the combined tracking algorithm, namely, lower measurement noise.

Concerning the test at 32 dB-Hz, the Doppler measurements are noisier and the PLI is not as close to unity, but the combined tracking algorithm is still able to track the signal. However, using only the L1C data channel, tracking is not possible. This shows that the combined approach offers increased sensitivity relative to tracking signal individually. This was expected since the powers of the three signals are effectively added when they are combined with the proper sign.

Concerning the test with motion (third simulation), the Doppler and PLI measurements for the combined tracking and the tracking using the L1C data channel only are shown in Figure 10. As for the static test, the Doppler measurements from the combined method are less noisy than the ones of the L1C data channel only. Concerning the PLI, again the L1C data channel-only tracking is not as close to unity as the combined tracking and is similarly affected by the user motion.

From these tests, one can deduce that the combined tracking performs better than the traditional tracking but does not show a high improvement in case of dynamics. Indeed, since the same PLL bandwidth is used, the combined tracking is not able to perform better under dynamics. However, the increase in the number of observations increases the tracking accuracy and reliability.

6. Combined Kalman Filter for Tracking

The Kalman filters have been used for signal tracking in several works to replace the discriminators and filters of the PLL and DLL [1922]. In all these cases, only one signal at a time was tracked, but results showed improved accuracy and precision in the carrier phase measurements (relative to the traditional approaches) which directly implies a better solution in position. Therefore, a Kalman filter is developed here which uses all the GPS signals in the L1 band in order to again improve the carrier phase measurements.

6.1. Method

This section describes a Kalman filter to simultaneously track the L1 C/A signal and the L1C (pilot and data channels) signal. In the following, the TMBOC modulation is considered as a BOC(1,1) modulation (i.e., the BOC(6,1) is ignored) since the loss created by this approximation is inconsequential but results in a receiver that is easier to implement [14].

One possible method to track all signals is to augment the state vector to include the L1C signal tracking errors. However, this approach implicitly considers the different signals as independent. Since the signals come from the same satellite and are time synchronized, one effectively loses information by doing so. Therefore, this approach was not selected.

Another approach is to augment the measurement model using the outputs of the 18 correlators (6 correlators per signal/channel). The method developed here uses this approach. However, for the sake of efficiency, the Kalman filter is implemented using a sequential measurement processing method [20]. The method is described below.

The measurements to the filter are the outputs of the 18 correlators, six for L1 C/A and six for each the pilot and data channels of L1C. The output of these correlators depends on the amplitude of the signal, the autocorrelation function, the error in code delay and the error in phase [23], namely, 𝐼=𝐴𝑅(𝛿𝜏Δ)cos(𝛿𝜙)+𝑛𝐼,𝑄=𝐴𝑅(𝛿𝜏Δ)sin(𝛿𝜙)+𝑛𝑄(2) with 𝐴=𝐴sin(𝜋𝛿𝑓𝑇)/𝜋𝛿𝑓𝑇 and 𝛿𝜙=𝛿𝜙0+(𝑇/2)𝛿𝑓0+(𝑇2/6)𝛿𝑎0 using the following notations.(i)𝐴 is the signal amplitude.(ii)𝛼D and 𝛼P are the scale factors to account for the difference in transmitted power of the L1C data and pilot channels relative to the L1 C/A signal (known). These coefficients are used to normalize the outputs of the correlators.(iii)𝑅 is the autocorrelation function of the corresponding signal (taking into account the number of code periods considered for each coherent integration).(iv)𝛿𝜏 is error in the local code delay.(v)Δ is spacing between the prompt and the early/late correlators.(vi)𝛿𝜙 is average error in the local carrier phase over the integration interval.(vii)𝑇 is coherent integration time interval.(viii)𝛿𝜙0 is error in the local carrier phase at the beginning of the integration interval.(ix)𝛿𝑓0 is error in the local carrier frequency at the beginning of the integration interval.(x)𝛿𝑎0 is at phase acceleration error (frequency rate error) the beginning of the integration interval.(xi)𝑛𝐼 are the noise values in the in phase correlator (Gaussian).(xii)𝑛𝑄 are the noise values in the in quadrature correlator (Gaussian).

The correlators outputs 𝐼 and 𝑄 are then normalized (for L1C) using the coefficients 𝛼D and 𝛼P.

The difference from the Kalman filter for one signal only (presented in [19, 21, 24]) is that in this case, there are three measurement vectors, one for each channel. The state vector is the same as that of the Kalman filter with one signal only since all the errors are the same between channels. The state vector is given by 𝑋𝑘=𝐴𝛿𝜏𝛿𝜙0𝛿𝑓0𝛿𝑎0𝑘.(3)

Only the amplitude of the signal changes between signals due to the difference in the transmitted power. This is compensated by using a scale factor in the amplitude term (see later in this section for more detail). For each channel, all the outputs of the correlators are normalized by the square root of the sum of the square of the in-phase and in-quadrature prompt correlators in order to match them with the proposed model of the autocorrelation function. The different signal powers are accounted for in the Kalman filter observation covariance matrix.

Concerning the autocorrelation of L1 C/A, the derivative is estimated first, and then is integrated to obtain the autocorrelation function. The derivative is generated as piecewise-continuous as follows:(i)𝑑𝑅(𝑡)/𝑑𝑡=1 for 𝑡=[10.01],(ii)𝑑𝑅(𝑡)/𝑑𝑡=100𝑡 for 𝑡=[0.010.01],(iii)𝑑𝑅(𝑡)/𝑑𝑡=1 for 𝑡=[0.011].

Integrating the above using appropriate initial conditions gives(i)𝑅(𝑡)=𝑡+1 for 𝑡=[10.01],(ii)𝑅(𝑡)=50𝑡2+0.995 for 𝑡=[0.010.01],(iii)𝑅(𝑡)=𝑡+1 for 𝑡=[0.011].

The result of the above approach deviates slightly from the ideal correlation function within the range of ±0.01 chips. This is necessary because the Kalman filter implementation (details to follow) requires that the autocorrelation function to be differentiable everywhere.

Concerning L1C, a similar approach for developing the autocorrelation function model as for L1 C/A was developed. However, the model becomes complicated due to the inflexion points resulting from the side peaks (around 0.5 chips from the main peak) and since it does not have an impact on the performance, these areas are approximated by a constant over very small intervals and the derivative is modeled without discontinuities. The model of the autocorrelation function and its derivative is presented in Figure 11.

With the exceptions of the small approximations in the autocorrelation functions mentioned above, the front-end bandwidth was essentially assumed to be infinite. In practice, however, the use of a finite bandwidth front-end will cause the autocorrelation function peaks to be even more rounded. This additional rounding of the correlation function will have a larger impact at low C/N0 but should be small overall. Specifically, the rounding will only affect the prompt correlator and will primarily affect the estimate of the signal amplitude. In contrast, the early and late correlators should be virtually unaffected because they are located at locations of constant slope and thus will not affect tracking. Overall, this approximation is not expected to introduce significant adverse effects.

To distinguish between channels, the channel observation vector is denoted as 𝑍channel,𝑘, where “channel” is either L1C/A, L1CD or L1CP: 𝑍channel,𝑘=𝐼_𝑃𝐼_𝐸𝐼_𝐿𝑄_𝑃𝑄_𝐸𝑄_𝐿channel,𝑘,(4) where(i)𝐼_𝑃, 𝐼_𝐸, and 𝐼_𝐿 are, respectively, the in phase prompt, early and late correlators outputs (for one channel);(ii)𝑄_𝑃, 𝑄_𝐸, and 𝑄_𝐿 are, respectively, the quadrature prompt, early and late correlators outputs (for one channel).

As for one channel only case, the measurement model has to be linearized as follows:𝛿𝑧channel,𝑘=𝜕channel𝑥𝑘,𝑘𝜕𝑥𝑘||||𝑥𝑘=̂𝑥𝑘𝛿𝑥𝑘+𝑣𝑘=𝐻channel,𝑘𝛿𝑥𝑘+𝑣𝑘,(5) where(i)𝑧𝑘 and 𝑥𝑘 are the real values of the measurements and the state vector;(ii)̂𝑧𝑘 and ̂𝑥𝑘 are the estimated values of the measurements and the state vector;(iii)𝛿𝑧𝑘 are the measurement misclosures;(iv)𝛿𝑥𝑘 is the perturbed state vector;(v)𝐻𝑘 is the design matrix and is detailed in the following.(vi)𝑣𝑘 is the noise associated with the measurement model (assumed Gaussian with covariance matrix 𝐶𝑣𝑘 detailed in the following).

Therefore, the measurement model can be expressed as 𝛿𝑧channel,𝑘=𝐼_𝑃channel𝐼_𝐸channel𝐼_𝐿channel𝑄_𝑃channel𝑄_𝐸channel𝑄_𝐿channel𝐴𝑅channel𝛿𝜏1cos𝛿𝜙1𝑅channel𝛿𝜏1Δcos𝛿𝜙1𝑅channel𝛿𝜏1+Δcos𝛿𝜙1𝑅channel𝛿𝜏1sin𝛿𝜙1𝑅channel𝛿𝜏1Δsin𝛿𝜙1𝑅channel𝛿𝜏1+Δsin𝛿𝜙1.(6) The covariance matrix 𝐶𝑣𝑘 associated to the measurement model can be expressed as follows:𝐶𝑣𝑘=𝜎2𝑛𝑅(0)𝜎2𝑛𝑅(Δ)𝜎2𝑛𝜎𝑅(Δ)0002𝑛𝑅(Δ)𝜎2𝑛𝑅(0)𝜎2𝑛𝜎𝑅(2Δ)0002𝑛𝑅(Δ)𝜎2𝑛𝑅(2Δ)𝜎2𝑛𝑅(0)000000𝜎2𝑛𝑅(0)𝜎2𝑛𝑅(Δ)𝜎2𝑛𝑅(Δ)000𝜎2𝑛𝑅(Δ)𝜎2𝑛𝑅(0)𝜎2𝑛𝑅(2Δ)000𝜎2𝑛𝑅(Δ)𝜎2𝑛𝑅(2Δ)𝜎2𝑛;𝑅(0)(7) where 𝜎2𝑛 is the variance of the noise and 𝑅 is the known autocorrelation function.

Because the measurement vector is the same for all the signals, we can use the same model as for L1 C/A, and apply it three times (in succession) per integration interval, using the measurements of the three channels. This is possible because the noise values associated with the measurement model are sufficiently uncorrelated between the three channels due to the orthogonality properties of the codes. For the receiver implemented herein, the integration time chosen is 10 ms due to the ranging code period of the L1C codes.

The signs of the outputs of the prompt correlators are used to determine the signs of the data bits on L1 C/A and L1C.

By using the sequential processing method, the sizes of the corresponding update matrices (i.e., 𝐻, 𝑅, and 𝐾) are not increased compared to the traditional Kalman filter used in previous research, thus providing an efficiency gain. Moreover, using this method (with appropriate modifications), one can update the different channels asynchronously and, for example, integrate longer on the pilot channel (no data bits) than on the data channels to track weaker signals.

Concerning the process noise associated with the dynamic model, the values used in static mode are presented in Table 4. The effect of dynamics is taken into account by increasing the variance of the phase acceleration process noise (higher values means higher dynamics).

6.2. Results

In order to analyze the performance of the combined Kalman filter, results are compared with those of a C/A code-only Kalman filter and to traditional tracking. The values of the process noise of the C/A code-only Kalman filter are the same as for the combined Kalman filter and the parameters of the traditional tracking are the same as in the previous section. Results of several tests in static and dynamic mode for each method are shown below.

It would have been nice to compare the performance of the combined tracking above with the Kalman filter at the same C/N0. However, we made the choice to have to show the performance at high C/N0 for all methods (45 dB-Hz) and at the limit of each method. Moreover, each test with the Kalman filter is compared to a traditional tracking.

Note that through the test, we choose to not simulate multipath in this simulator because a realistic model of multipath is very complex and is a topic by itself.

6.2.1. Tests in Static Mode

To analyze the performance of the combined Kalman filter and compare it with the performance of the Kalman filter with C/A only and the traditional tracking, two tests were simulated. The simulation parameters for the static tests are shown in Table 5.

For each test, the performance is assessed using the PLI and the Doppler values. Moreover, since all the results are generated using the software simulator, the true (simulated) Doppler and code delay are known. Therefore, these values are compared with the ones output by the different tracking strategies.

The Doppler error of the first test is shown in Figure 12. Both Kalman filters are less noisy than the traditional tracking. Indeed, the noise associated with the traditional tracking is about 1 Hz around the true Doppler, whereas the noise associated with both Kalman filters is similar and almost null. Figure 13 shows the code tracking error for the first test for the three tracking methods. All three methods converge to the true code delay. Secondly, the traditional tracking is much noisier than both Kalman filters: the error in code delay varies from about −2.5 m to +2.5 m in the case of traditional tracking and between −50 cm and +50 cm for both Kalman filters. The variance of code delay error of the combined Kalman is smaller than that for the Kalman filter with only L1 C/A.

Finally, Figure 14 shows the PLI for the first test for the three tracking methods. All three PLI converge to approximately unity and the differences between the three approaches are insignificant. Therefore, with a high C/N0 and without any motion or error, the combined Kalman filter is able to track slightly better than the two other methods. Moreover, the three signals are tracked simultaneously using the combined Kalman filter. The performance does not change significantly compared to the Kalman filter with L1 C/A alone. However, this is not surprising because the results were already very good using only L1 C/A because of the strong signals considered, and it must be borne in mind that the simulator does not simulate any errors except white noise.

Concerning the second static test, the results are shown in Figure 15. Here, only the results for the combined Kalman filter and the Kalman filter for C/A only are shown. The traditional filter was not able to track the signal in this case. As shown on the top figures (Doppler and Doppler error), both filters are able to track at low C/N0 (30 dB-Hz), and the Kalman filter for C/A only is slightly noisier than the combined Kalman filter. The same conclusion can be drawn from the bottom left figure, showing the code delay error. The code delay error stays larger for the Kalman filter for C/A only than the combined Kalman filter. This is again confirmed by the PLI which experiences marginally lower values for the Kalman filter for C/A only than the combined Kalman filter, even if the PLI is still very close to one. Therefore, despite the increased noise, both Kalman filters are able to track at low C/N0 values, whereas the traditional tracking method cannot. Moreover, the combined Kalman filter shows better performance than the Kalman filter for C/A only for all the parameters considered. Nevertheless, one has to keep in mind that these results have been generated using simulated data, thus only white noise is simulated. With real data, the difference is expected to be larger since other effects will be present. In particular, weak signals will likely also contain multipath, which can degrade the performance of the filter.

The improvement brought by the Kalman filter relative to the traditional tracking method can be explained by the selection of the process noise parameters and by the fact that the code, frequency, and phase are all tracked together, instead of in separate loop filters (linked only through the loop aiding). Note that in the case of varying C/N0, the Kalman filter can weight measurements according to the C/N0 (whereas a regular loop gives equal weights to all measurements) and therefore should exhibit better performance. Moreover, the use of the three signals with the combined Kalman filter increases the number of observations threefold and therefore improves the tracking performance and is able to track better at high and low C/N0.

6.2.2. Tests in Dynamic Mode

The performance of the combined Kalman filter is now compared with the performance of the L1 C/A code Kalman filter and the traditional Kalman filter in terms of tracking accuracy and the maximum level of dynamics supported. For all these tests, the variances of the noise associated with the dynamic model for the single signal Kalman filter and the combined Kalman filter are set to the same values in order to compare their performance with the same level of process noise. However, the value of the process noise corresponding to phase acceleration is increased depending on the dynamics as explained above. Toward this goal, three simulations are realized with the parameters shown in Table 6. The motion parameters are summarized in Table 7. The values of the phase acceleration process noise are the maximum acceleration values squared.

The velocity in function of time for these two tests is shown in Figure 16.

The Doppler, Doppler error, code delay error, and PLI over time for the first dynamic test are shown in Figure 17. In the case of the satellite used in this test, about one third of the user dynamics are projected on the line of sight between the satellite and the receiver. As seen in this figure, the three tracking methods are able to track under the simulated dynamics at high C/N0. The variance of the Doppler error is much smaller in the case of the Kalman filters than for the traditional filter. Moreover, at each change in dynamics (every 5 seconds), the traditional tracking is clearly affected (as demonstrated by the PLI) whereas the two other methods are only slightly affected. Furthermore, in the figures showing the Doppler error and the PLI over time, we can see that the combined Kalman filter is less affected by the change in dynamics than the C/A code Kalman filter. Concerning the code delay error, the noise of the combined Kalman filter is slightly smaller than for the C/A code Kalman filter which is already much smaller than the traditional tracking. Finally, on the plots of the PLI and the Doppler error over time, one can see that the biggest impact of the dynamics is around 20 seconds, which corresponds to a change of dynamics from −10 m/s2 to +10 m/s2. Therefore, in presence of reasonable dynamics and at high C/N0, the combined Kalman filter performs best.

The Doppler, Doppler error, code delay error and PLI over time for the second dynamic test for the both Kalman filters are shown in Figure 18. The traditional method is not able to track at low C/N0. The combined Kalman filter is able to follow the medium dynamics at low C/N0. However, as it can be seen with the Doppler error, the code delay error, and the PLI plots, the noise is higher at lower C/N0. Nevertheless, the combined Kalman filter is able to track medium dynamics at low C/N0 better than the C/A code Kalman filter (which is more affected by the dynamic changes every 5 seconds as shown by the PLI).

The Doppler, Doppler error, code delay error, and PLI over time for the third dynamic test for the three methods are shown in Figure 19. Both Kalman filters are able to track the high dynamics for the length of the simulation with medium Doppler and code delay errors. Both Kalman filters experience a jump at every change of dynamics (every 5 seconds) as visible on the plots of the Doppler error and the code delay error. However, both Kalman filters recover much faster than the traditional tracking (as evidenced by the Doppler error for the highest dynamic change). Moreover, as demonstrated by the PLI, the traditional tracking experiences a lot of difficulties at each dynamic change. In addition, the Doppler noise and code delay error are smaller for the combined Kalman filter than for the C/A code Kalman filter. Therefore, the combined Kalman filter yields the best tracking and error performance in presence of high dynamics at high C/N0.

6.3. Initial Multipath Performance

An important source of error in any GNSS receiver is multipath. Although a comprehensive multipath analysis is beyond the scope of the paper, a basic simulation was performed to illustrate the performance of the different tracking approaches in the presence of multipath. Specifically, only a single reflection, constant delay model was simulated with the software signal simulator described earlier. In addition to the direct signal, a multipath signal with half the amplitude of the direct signal (6 dB less power) and a delay of exactly 0.3 C/A code chips (~87.9 m) was generated on each of the three signals. Since the delay is an integer number of L1 cycles, the multipath signal is perfectly in phase with the direct signal. All other simulation parameters were the same as before. The multipath corrupted signal was then processed using the same receiver parameters as before.

Figure 20 shows the code delay error for each of the traditional C/A code tracking algorithm (blue line), the Kalman filter using only the C/A code (red line) and the combined Kalman filter (black line). As can be seen, the C/A code-only results obtained using either traditional or Kalman filter tracking are approximately the same after the convergence period and have a steady state error of approximately 10 m. The main differences are the convergence period (longer for the traditional tracking) and the level of noise (smaller in the case of the Kalman filter). This is expected in this case since the Kalman filter does not include a specific model to mitigate (estimate) multipath.

Compared to the C/A code-only error, the combined Kalman filter error is much smaller with a steady state error of only about 3 m. This improvement is explained by the shape of the autocorrelation function of the L1C signal. More specifically, these results are consistent with the weighted average of the multipath error envelopes found in Van Dierendonck [5] for L1 C/A and Macchi et al. [14] for the approximation an MBOC by a BOC(1, 1). Since all three signals have the same delay and phase, no significant improvement is expected beyond simply averaging the multipath error of each signal.

7. Conclusions and Future Work

This paper first presented a method to combine the GPS L1 C/A and L1C (data and pilot channels) signals during the acquisition process. This method was tested using a software signal simulator at low and high C/N0 values. At low C/N0 (32 dB-Hz), the combined acquisition was able to acquire signals that could not be acquired using the L1C data channel alone using a coherent integration time of 10 ms in both cases. Moreover, the three channels (L1 C/A and L1C data and pilot channels) were acquired simultaneously.

For tracking, a new method was developed to combine the GPS L1 C/A and L1C (data and pilot channels) signals. This algorithm was tested using a signal software simulator at high (45 dB-Hz) and low C/N0 (30 dB-Hz) and in static and dynamic modes (about 10 g). The combined tracking showed better performance than the traditional tracking of the L1C data channel. Indeed, at high C/N0 (45 dB-Hz), the Doppler measurements were less noisy and the PLI was closer to unity. At low C/N0 (30 dB-Hz), only the combined tracking was able to track the signals. However, as expected, the combined tracking did not show a significant improvement compared to the traditional tracking in case of dynamics.

Finally, using a single signal Kalman filter as a starting point, a Kalman filter combining the GPS L1 C/A and L1C (pilot and data channels) channels was developed to track all three signals simultaneously. The combined Kalman filter was compared with the single-signal Kalman filter and the traditional tracking under various conditions. Both Kalman filters showed better performance than the traditional tracking (more accurate tracking and far less noise) in all cases. Furthermore, the Kalman filters were able to track at low C/N0 (30 dB-Hz) whereas the traditional tracking could not. Moreover, the combined Kalman filter showed better performance than the single signal Kalman filter for one signal in static mode at low and high C/N0. The combined Kalman filter also showed better performance in presence of reasonable dynamics as well as high dynamics. Finally, a simple multipath assessment showed that the combined tracking approach produces errors consistent with the weighted average of the multipath errors of each signal.

A summary of the performance for the dynamic and static tests and for each method is shown in Tables 8 and 9.

Acknowledgments

The first author would like to thank Cyrille Gernot, Ph.D. student in geomatics engineering, for his help and support. The authors would like to thank the Informatics Circle Of Research Excellence and the GEOIDE Networks of Centres of Excellence for their financial support.