Study of the contribution of an ionospheric model embedded on a dual frequency receiver

(1)

Study of the contribution of an ionospheric

model embedded on a dual frequency

receiver

Alice Tourtier

Master of Science Thesis in Geodesy No. 3142

TRITA-GIT EX 16-002

School of Architecture and the Built Environment

Royal Institute of Technology (KTH)

Stockholm, Sweden

(2)

Abstract English

GNSS (Global Navigation Satellite System) is a system that allows a receiver to know its position. It includes GPS (Global Positioning System), a system designed by the United States, Galileo (European), Beidou (Chinese) and others.

(3)

Abstract Swedish

GNSS (Global Navigation Satellite System) är system som gör det möjligt för en mottagare att känna sin position. Det täcker över GPS (Global Positioning System), ett system som drivs av USA, Galileo (europeisk), Beidou (kinesiska) och andra.

(4)

Acknowledgments

First, I deeply thank my supervisors at CNES, Marion Aubault and Sebastien Rougerie for their guidance during my thesis work and allowing me to perform this thesis.

I express my gratitude to Anna Jensen, professor at KTH Royal Institute of Technology for approving this master thesis and for providing useful comments and remarks during my master thesis. I also thank Nickolay Ivchenko, my master director, for approving this subject.

My thanks also go to all my colleagues from the Signaux et equipement de Radionavigation service, especially Yoan Gregoire, for answering my questions and their support during the time I worked at CNES.

(5)

Content

1. Introduction ... 10

2. Overview on GNSS ... 11

2.1. Principle of GNSS ... 11

2.1.1. Clock bias ... 11

2.1.2. Computation of the position ... 11

2.2. GNSS Signal description ... 12 2.2.1. PRN Code ... 12 2.2.2. Navigation message ... 13 2.2.3. Signal model ... 13 2.3. Structure of a receiver ... 14 2.3.1. Signal processing ... 15

2.3.2. The Lock Loop ... 18

3. Propagation ... 20

3.1. The ionosphere ... 20

3.1.1. Ionosphere influence on electromagnetic waves. ... 20

3.1.2. Description of the electron density of the ionosphere ... 21

3.1.3. Thin layer assumption ... 22

3.1.4. Examples of ionosphere models ... 23

3.2. Presence of multipath ... 24

3.2.1. Effect on the DLL ... 24

3.2.2. Effect on the PLL ... 26

3.3. The observation model ... 27

3.3.1. The observables ... 27

3.3.2. Propagation effects ... 27

3.3.3. Obsevation model ... 27

3.4. Problem ... 29

4. Techniques to estimate the ionospheric delays ... 30

4.1. Dual and single frequency measurments ... 30

4.1.1. Comparison of the code’s delays ... 30

4.1.2. Comparison of the phase’s delay ... 30

4.1.3. Comparison of both phase’s and code’s delay ... 30

4.1.4. Observable possible ... 31 4.2. Estimation Strategy ... 31 4.2.1. Model 0 ... 31 4.2.2. Model 1 ... 31 4.2.3. Model 2 ... 32 4.2.4. Model 3 ... 33 4.2.5. Model 4 ... 36

4.3. Improvement of the models ... 37

4.3.1. Model 0 ... 37

4.3.2. Model 1 ... 37

(6)

4.3.5. Model 4 ... 41

4.3.6. Remarks ... 43

4.4. Loss of the frequency of a satellite. ... 44

4.4.1. Model 1 ... 44 4.4.2. Model 2 ... 44 4.4.3. Model 3 ... 44 4.4.4. Model 4 ... 44 5. Simulation ... 46 5.1. Tools ... 46 5.1.1. Spring software ... 46

5.1.2. The GNSS receptor simulator ... 47

5.2. Simulation environment ... 47

5.2.1. Environment 1: Fixed point at mid latitudes at midday. ... 47

5.2.2. Environment 2: Fixed point at low latitudes at midday. ... 48

5.2.3. Environment 3: Fixed point at mid latitudes in the morning ... 49

5.2.4. Environment 4: Fixed point at low latitudes in the morning ... 49

5.2.5. Comparison of the different simulations ... 49

5.3. Results ... 49

5.3.1. Parameters of the Kalman filter ... 49

5.3.2. Environment 1 ... 50 5.3.3. Environment 2 ... 57 5.3.4. Environment 3 ... 62 5.3.5. Environment 4 ... 66 5.3.6. Other factors ... 70 6. Final conclusion ... 72 6.1. Conclusion ... 72 6.2. Future work ... 72 7. Bibliography ... 73

8. Annex 1: Calculation of the shape of the signal after correlation. ... 74

9. Annex 2: Kalman filter ... 76

9.1. Theory on Kalman filter ... 76

9.2. Expressions of the inputs parameters in the models used in this study ... 76

(7)

List of figures

Figure 1 : Presentation of multipath ... 10

Figure 2 : Autocorrelation and correlation of PRN 1 and 17 [11]_{... 13}

Figure 3 : Structure of a receiver ... 15

Figure 4 : Tracking of the phase’s delay ... 16

Figure 5 : Tracking of the code's delay ... 17

Figure 6 : Early, Prompt, Late synchronization [2]_{... 18}

Figure 7 : Block diagram of a lock loop [10]_{... 18}

Figure 8 : Electron density in the ionosphere [7]_{... 21}

Figure 9 : Geographic fluctuation [7]_{... 22}

Figure 10 : Thin Layer Model [3]_{... 22}

Figure 11 : Difference between STEC and VTEC ... 23

Figure 12 : Correlation function in presence of multipath [2]_{... 25}

Figure 13: Tracking error envelope of the code's delay [2]_{... 26}

Figure 14 : Maximum and minimum on the error of the PLL [2]_{... 26}

Figure 15: Situation with two satellites ... 33

Figure 16 : Representation of the local ionospheric model ... 34

Figure 17: Process to collect data ... 46

Figure 18 : Example of SPRING processing ... 46

Figure 19 : Approximate position of the receiver of environment 1 ... 48

Figure 20: Location of the warehouse and the receiver ... 48

Figure 21: Approximate position of the receiver of environment 2 ... 49

Figure 22: Legend of the results ... 50

Figure 23: Results for Model 0 in Toulouse at midday ... 51

Figure 28: Comparison of the different results for environment 1 ... 56

Figure 29: Results for Model 0 in Tamanrasset at midday ... 57

Figure 33: Comparison of the different models for the second environment ... 61

Figure 34: Results for Model 0 in Toulouse in the morning ... 62

Figure 38: Comparison of the different models for the third environment ... 65

Figure 39: Results for Model 0 in Tamanrasset in the morning... 66

Figure 43: Comparison of the different models for the fourth environment ... 70

Figure 44: Comparison of the different models ... 72

(8)

List of tables

(9)

Glossary

DLL: Delay Lock Loop

GNSS: Global Navigation Satellite System HF: High Frequency

FLL: Frequency Lock Loop IPP: Ionospheric Pierce Point LS: Least Square

LOS: Line of sight PLL: Phase Lock Loop PSD: Power Spectral Density PRN: Pseudo Random Noise PVT: Position Velocity Time STEC: Slant Total Electron Content TEC: Total Electron Content

(10)

1. Introduction

CNES (Centre National d’Etudes Spatiales) is the French National Space Center. It conducts research in five different fields: design of Ariane (a European launcher); exploration of the space; observation of the Earth; defense and telecommunication.

The study I have carried out for my master thesis, “Study of the contribution of an ionospheric model embedded on a GNSS receiver”, is within the telecommunication field. The work that was carried out is part of an internship.

The signal transmitted between a satellite and a receiver is altered by different factors such as delays when spreading through the troposphere and ionosphere. To estimate these alterations, GNSS multi-frequencies receivers are currently used due to the fact that ionospheric delay is frequency dependent. However, the signal can also be diffracted / reflected / scattered by the nearby environment of the receiver (building, trees, cars …). An example of a wave scattered by a house is presented in Figure 1. As it is well know that all these contributions (called multipath) lead to errors on the receiver position estimation, the impact of multipath on the estimations of ionospheric delays needs to be tackled. The goal of my study is to propose a new multi-frequency receiver architecture which may be robust to multipath in order to improve the ionospheric delay estimation. The main idea is to use an on-board model for the ionosphere to estimate a value for the ionospheric delay close to the real one even if the signal is modified by the presence of multipath.

The CNES leads research activities about GNSS multi-frequency receivers robust to multipath and ionospheric effects. This project is a part of these activities, and compensate the multipath in the ionospheric delay estimation is a first step to improve at the end the estimation of the position. Frequency diversity may also be used to compensate multipath, but this point is out of the scope.

Figure 1 : Presentation of multipath

(11)

2. Overview on GNSS

Global Navigation Satellite System (GNSS) is used to obtain the position of a receiver whether it is fixed on Earth or orbiting around it at an altitude lower than the one of the GNSS satellites. It includes systems such as GPS or Galileo. For each system, a constellation of satellites is orbiting around the Earth. Each satellite sends a number of signals which can be used for positioning.

2.1. PRINCIPLE OF GNSS

Each receiver receives from all the satellites that are in its field of view, a navigation message. In this navigation messages, there is information for the receiver to know the position of the satellites. The receiver is also able to estimate the propagation time of the signal (see section 2.3). This propagation time multiplied by the speed of light is equal to the distance between the satellite and the receiver.

Mathematically, if 𝜏 is the propagation time, the distance 𝑅𝑋 between the satellite X and the receiver is equal to [3]_:

𝑅𝑋= 𝑐 ∗ 𝜏 + 𝑛𝑥

Equation 2.1-1 Where:

- c is the speed of light

- 𝑛𝑥 is the thermal noise in the receiver

If 𝑃 = ( 𝑋 𝑌 𝑍) and 𝑃𝑆𝑋= ( 𝑋𝑆_𝑋 𝑌𝑆_𝑋

𝑍𝑆_𝑋) are respectively the position of the receiver and the satellite X, the distance between the satellite and the receiver is, from a geometric point of view, equal to [3]_:

𝑅𝑋,𝑔𝑒𝑜= √(𝑋 − 𝑋𝑆𝑋) 2 + (𝑌 − 𝑌𝑆𝑋) 2 + (𝑍 − 𝑍𝑆𝑋) 2 Equation 2.1-2

2.1.1. CLOCK BIAS

A specific reference time for GPS was introduced in 1980 and is used to synchronize all the clocks in GPS satellites and receivers. However, there is always a bias between time of such a clock and the reference time. The two distances computed in the previous section are not exactly the same but differ from a factor equal to the sum of the satellite and receiver clock bias multiplied by the speed of light.

𝑅𝑋= 𝑅𝑋,𝑔𝑒𝑜+ 𝑐 ∗ (𝛿𝑇𝑆_𝑋+ 𝛿𝑡) + 𝑛𝑋 Equation 2.1-3 𝑅𝑋= √(𝑋 − 𝑋𝑆_𝑋)2+ (𝑌 − 𝑌𝑆_𝑋) 2 + (𝑍 − 𝑍𝑆_𝑋)2+ 𝑐 ∗ (𝛿𝑇𝑆_𝑋+ 𝛿𝑡) + 𝑛𝑋 Equation 2.1-4 With, 𝛿𝑇𝑆_𝑥,, 𝛿𝑡, respectively the satellite and receiver clock bias.

Thanks to the data contained in the navigation message, an estimation of the satellite clock bias can be found. However, the receiver clock bias is considered as an unknown parameter. The pseudo-range between the receiver and the satellite 𝑑𝑋 without the error of the satellite clock is equal to:

𝑑𝑋= √(𝑋 − 𝑋𝑆𝑋) 2 + (𝑌 − 𝑌𝑆𝑋) 2 + (𝑍 − 𝑍𝑆𝑋) 2 + 𝑐 ∗ 𝛿𝑡 + 𝑛𝑋 Equation 2.1-5

2.1.2. COMPUTATION OF THE POSITION

(12)

𝑃 = ( 𝑋 𝑌 𝑍 𝛿𝑡 )

At least four satellites are needed to compute the position of the receiver. If more than four satellites are available, algorithms such as the least square‘s algorithm estimation are used to compute the position while reducing the errors.

2.2. GNSS SIGNAL DESCRIPTION

In this section, it is explained the shape of the signal and how it propagates in the atmosphere. GPS satellites emit waves at three different frequencies:

- 𝑓1= 1575.42 𝑀𝐻𝑧, called the L1 band, - 𝑓2= 1227.60 𝑀𝐻𝑧, called the L2 band

- 𝑓5= 1176.45 𝑀𝐻𝑍, called the L5 band, launched since 2010.

In this study, only data coming from signals at frequencies 𝑓1 and 𝑓2 are used. Here, we will describe the Coarse Acquisition (C/A) signal broadcast on the L1 band. This signal is based on two components: the PRN code (section 2.2.1) and the navigation message (section 2.2.2).

2.2.1. PRN CODE

A family of Pseudo Random Noise (PRN) codes is used to transmit the signal from the satellite to the receiver. The length of the code used in GPS satellites is of 1023 chips. Its frequency is 𝑓 = 1.023 𝑀ℎ𝑧. Thus, the length of one chip is approximately 𝑇𝑐ℎ𝑖𝑝≈ 1 µ𝑠 and it takes about 𝑇𝑐𝑜𝑑𝑒= 1 𝑚𝑠 for the code to repeat itself. It has the two following properties:

- For each different code, the autocorrelation of the code is close to zero everywhere except in zero where the value of the autocorrelation is 1

- The cross-correlation of two different codes from the same family is close to zero everywhere. For PRN codes used in GPS satellite, the results of the autocorrelation and the correlation of two codes 𝐶𝑖 and 𝐶𝑗 are [2]: 𝐾𝐶𝑖,𝐶𝑗 (𝑚)= { −1/1023 𝑜𝑟 −65/1023 𝑜𝑟 63/1023 Equation 2.2-1

Unless 𝑖 = 𝑗 (autocorrelation of the same code), and 𝑚 = 0, 𝐾𝐶𝑖,𝐶𝑖(0) = 1

(13)

Figure 2 : Autocorrelation and correlation of PRN 1 and 17 [11]

In practice, the information (see section 2.2.2) emitted by a satellite is multiplied by a different code of the same family. The Power Spectral Density of the GNSS signal is therefore spread on a larger band and gets a value lower than the thermal noise. Thus, GNSS signals do not disturb any other telecommunications systems, but can be difficult to track in noisy environment

However the receiver is able to recover the signal emitted by a satellite by multiplying the signal it receives by the same PRN code as the one used by the satellite (see section 2.3.1).

2.2.2. NAVIGATION MESSAGE

Each satellite sends messages with information regarding the satellite. For example, the GPS navigation message contains the information the receiver must know about the satellite to compute its position: the satellite’s coordinates and the satellite clock bias. The length of the navigation message is 750 𝑠 and is composed of 37500 bits. Thus, there is one bit every 𝑇𝑏𝑖𝑡 = 20 𝑚𝑠 [2].

As seen in section 2.2.1, the length of the PRN code is 𝑇𝑐𝑜𝑑𝑒 = 1𝑚𝑠, the PRN repeats itself 20 times between two navigation bits.

2.2.3. SIGNAL MODEL

Signal delivered by the satellite

The signal delivered by the satellite has the following shape [2]_:

𝑠(𝑡) = 𝐴 ∗ 𝑐(𝑡) ∗ 𝑑(𝑡) ∗ cos (2πf0𝑡 + ϕ0)

Equation 2.2-2

With:

- A: the amplitude of the signal - c: the PRN code

- d: the navigation message - 𝜙0: the initial random phase shift - 𝑓0: the carrier frequency

Propagation of the signal in the atmosphere

(14)

from the receiver to the satellite. However, the satellite is not fixed in space and has a certain velocity relative to the receiver. The frequency seen by the receiver is the sum of the carrier frequency 𝑓0 and the frequency of the Doppler Effect 𝜈.

Thus, the signal reaching the receiver coming from the satellite 𝑆𝑥 is [2]_:

sreceived,Sx(t) = A ∗ c(t − τ) ∗ d(t − τ) ∗ cos (2π(f0+ ν)(t − τ) + ϕ0) + 𝑛(𝑡)

Equation 2.2-3 𝑠𝑟𝑒𝑐𝑒𝑖𝑣𝑒𝑑,𝑆𝑥(𝑡) = 𝐴 ∗ 𝑐(𝑡 − 𝜏) ∗ 𝑑(𝑡 − 𝜏) ∗ cos (2π(f0+ 𝜈)𝑡 + ϕ) + 𝑛(𝑡)

Equation 2.2-4 With:

- 𝐴, the amplitude of the signal

- 𝜏, the propagation time between satellite 𝑆𝑥 and the receiver.

- 𝜈, the Doppler effect due to the motion of satellite 𝑆𝑥 relative to the receiver - 𝜙 = 𝜙0− 2𝜋(f0+ 𝜈)𝜏, the phase delay

- 𝑛, the thermal noise.

In this signal, there is information that the receiver must figure out in order to compute its position or velocity: the propagation time 𝜏, the phase delay 𝜙, the Doppler Effect 𝜈 and the information signal. The propagation time 𝜏 can also be called the code’s delay. In section 2.3, the structure of a receiver and the method used to recover the data are presented.

Signal received by the receiver

In this section, it is assumed that there is no multipath on any of the satellites. Section 3.2 will describe the modification induced by the presence of multipath. The signal received by the receiver is the sum of all the signals sent by the different visible satellites of the GNSS constellation plus a background noise due both to the presence of other electromagnetic waves in the atmosphere and the electronic in the antenna.

The signal received by the receiver is:

𝑠𝑟𝑒𝑐𝑒𝑖𝑣𝑒𝑑(𝑡) = ∑𝑁𝑠𝑎𝑡(sreceived,S_i(t) 𝑖=1 + 𝑛𝑖(𝑡)) Equation 2.2-5 𝑆𝑟𝑒𝑐𝑒𝑖𝑣𝑒𝑑(𝑡) = [∑ 𝐴i∗ 𝑐𝑖(𝑡 − 𝜏𝑖) ∗ 𝑑𝑖(𝑡 − 𝜏𝑖) ∗ cos (2π(f0+ 𝜈𝑖)𝑡 + ϕi) 𝑁_𝑠𝑎𝑡 𝑖=1 ] + 𝑛(𝑡) Equation 2.2-6 Where: 𝑛(𝑡) = ∑ 𝑛𝑖(𝑡) 𝑁_𝑠𝑎𝑡 𝑖=1 With:

- 𝑁𝑠𝑎𝑡, the number of visible satellite,

- 𝐴𝑖, the amplitude of the signal coming from satellite 𝑆𝑖

- 𝜏𝑖, the propagation time between satellite 𝑆𝑖 and the receiver,

- 𝜈𝑖, the Doppler effect due to the motion of satellite 𝑆𝑖 relative to the receiver, - 𝜙𝑖, the phase delay between satellite 𝑆𝑖 and the receiver,

- 𝑛, the thermal noise coming from all the satellite, - 𝑐𝑖, the PRN code of satellite 𝑆𝑖

- 𝑑𝑖, the information regarding the navigation message sent by satellite 𝑆𝑖

2.3. STRUCTURE OF A RECEIVER

The goal of a GNSS receiver is to compute its position thanks to the signal received by its antenna. A GNSS receiver can be separated in three main blocks plus an antenna as presented in Figure 3: the high frequency (HF) block, the signal processing block and the data processing block

The HF block:

(15)

The signal processing block:

The goal of this block is to figure out, for each of the satellite contributing to the signal received by the receiver, the values of the code’s and phase’s delays. The contribution of each satellite is separated from the others by correlating the signal received by the receiver with a local version of the same PRN code as the one used by the satellite to emit the signal. The contribution of each satellite is processed by a different channel and for each of them; lock loops for the phase’s and for the code’s delays are used to track these values. A short description of the theory on look loops is presented in section 2.3.2.

The data processing block:

The last block of a GNSS receiver is the data processing block. For each satellite, the receiver does two tasks in parallel: demodulation of the navigation message and calculation of the pseudo-range between the receiver and the satellite from the values of the code’s and phase’s delays (see section 3.3). Part of the treatment to demodulate the signal is done in the signal processing block but it is in this block that the values of data bits are figured out. Then it gathers all this information to compute the position or the velocity (PVT- Position, Velocity, Time) of the receiver by using the method described in section 2.1. In Figure 3, the capital D’s are the distance calculated by the receptor between itself and a satellite.

Figure 3 : Structure of a receiver

2.3.1. SIGNAL PROCESSING

The goal of this block is to track the values of the code’s and phase’s delays for each satellite. The receiver has a certain number of channels and analyses the contribution from each satellite in parallel. In this section, it is described what happens within one channel and the signal it is mentioned is only the part of the received signal coming from the satellite that is analyzed (see Equation 2.2-4).

Within the channel, there are two modes possible: the acquisition mode and the tracking mode. Indeed, before the tracking of the code’s and phase’s delays can be possible, a rough estimation of the Doppler frequency and the code’s delay must be found in the acquisition mode. Every time the signal is lost, the receiver must redo the acquisition mode before the tracking can restart.

2.3.1.1. ACQUISITION MODE

The goal of the acquisition is to detect a satellite and provide a first estimation of the code’s delay and the Doppler frequency as the tracking mode needs to be initialized. The acquisition is done by evaluating the value of an output parameter for all the values possible for the couple (𝜏, 𝜈). The value of the output parameter must be high when the right couple (𝜏, 𝜈) is found and low otherwise so that the software can know when it has tested the right values.

2.3.1.2. TRACKING MODE

(16)

of the Doppler frequency so if one parameter is known, the other one is also known (possibly with an ambiguity). Frequency Lock Loop (FLL) is used to track the Doppler frequency. Given the situation, the receiver tracks the phase’s delay or the Doppler frequency, the FLL is less precise but more robust than the Phase Lock Loop (PLL).

The block diagrams (inspired from [2]) for the tracking are: Tracking of the phase’s delay

Figure 4 : Tracking of the phase’s delay

The treatment done on the signal is the same as in the acquisition mode at the beginning. The expressions of 𝐼 and 𝑄 are the same as in Equation 2.2-3 and Equation 2.2-4.

The parameter 𝜖𝜙_𝑖 is the error (also called the discriminator) of the phase’s delay for satellite 𝑆𝑖. This parameter is the input parameter of the PLL (see section 2.3.2 to see the use of this parameter). The expression of this parameter is a function of 𝐼 and 𝑄. Depending on the receiver, there are different ways to express it but in 𝐼 and 𝑄, only the cosine and the sine of this error is present so no matter the way to express 𝜖𝜙_𝑖, an ambiguity of an integer multiple of 2𝜋 appears.

In the receiver used in this study and presented in section 5.1.2, the discriminator 𝐷𝑃𝐿𝐿 used for the PLL is: 𝐷𝑃𝐿𝐿 = atan (𝑄

𝐼)

Equation 2.3-1 Assuming that there is no noise in the expression of 𝐼 and 𝑄 given in equations Equation 2.2-3 and Equation 2.2-4, 𝐷𝑃𝐿𝐿= atan (𝑠𝑖𝑛 (𝜙𝑖(𝑘) − 𝜙̂ (𝑘))𝑖 𝑐𝑜𝑠 (𝜙𝑖(𝑘) − 𝜙̂ (𝑘))𝑖 ) Equation 2.3-2 𝐷𝑃𝐿𝐿= 𝜙𝑖(𝑘) − 𝜙̂ (𝑘) = 𝜖𝜙𝑖 _𝑖 Equation 2.3-3 The discriminator is indeed equal to the error of the phase’s delays.

(17)

Figure 5 : Tracking of the code's delay

The tracking of the code’s delay requires more calculations as the one for the phase’s delay. The error or discriminator of the code’s delay 𝜖𝜏_𝑖for satellite 𝑆𝑖 needs more input parameters. It is based on the comparison of the values of the autocorrelation function of the PRN code at two different times: one in advance of 𝑑₂ compare to the error of the code’s delay 𝜖𝜏_𝑖 and one delayed of 𝑑₂. 𝑑 is called the spacing chip.

There are three cases possible (illustration is given in Figure 6): - 𝐾𝑐_𝑖(𝜖𝜏_𝑖−𝑑

2) = 𝐾𝑐𝑖(𝜖𝜏𝑖+ 𝑑

2), the error 𝜖𝜏𝑖 is null, - 𝐾𝑐_𝑖(𝜖𝜏_𝑖−𝑑

2) ≥ 𝐾𝑐𝑖(𝜖𝜏𝑖+ 𝑑

2), the error 𝜖𝜏𝑖 is positive and the synchronization is in advance, - 𝐾𝑐_𝑖(𝜖𝜏_𝑖−𝑑₂) ≤ 𝐾𝑐_𝑖(𝜖𝜏_𝑖+𝑑₂), the error 𝜖𝜏_𝑖 is negative and the synchronization is late.

The values of 𝐾𝑐_𝑖(𝜖𝜏_𝑖−𝑑

2) and 𝐾𝑐𝑖(𝜖𝜏𝑖+ 𝑑

2) are known thanks to the parameters 𝐼𝐸, 𝐼𝑃, 𝐼𝐿, 𝑄𝐸, 𝑄𝑃 and 𝑄𝐿 (E stands for Early, P for Prompt and L for Late).

As for the PLL, the discriminator of the code’s delay 𝜖𝜏_𝑖 is the input of the DLL used to track the code’s delay. In the receiver used in this study and presented in 5.1.2, the discriminator 𝐷𝐷𝐿𝐿 used for the DLL is:

𝐷𝐷𝐿𝐿 =2 − 𝑑 ∗ 𝑘𝑑𝑖𝑠𝑐𝑟𝑖 2 ∗ 𝑘𝑑𝑖𝑠𝑐𝑟𝑖 ∗ ( √𝐼𝐸2+ 𝑄𝐸2− √𝐼𝐿2+ 𝑄𝐿2 √𝐼𝐸2+ 𝑄𝐸2+ √𝐼𝐿2+ 𝑄𝐿2) ) Equation 2.3-4 Assuming there is no noise,

𝐷𝐷𝐿𝐿=2 − 𝑑 ∗ 𝑘𝑑𝑖𝑠𝑐𝑟𝑖 2 ∗ 𝑘𝑑𝑖𝑠𝑐𝑟𝑖 ∗ ( |𝐾𝑐_𝑖(𝜖𝜏_𝑖−𝑑_{2) | − |𝐾}𝑐_𝑖(𝜖𝜏_𝑖+𝑑_{2) |} |𝐾𝑐𝑖(𝜖𝜏𝑖− 𝑑 2)| + |𝐾𝑐𝑖(𝜖𝜏𝑖+ 𝑑 2) | ) Where:

(18)

Several other discriminator formulations exist for the DLL, but will not be listed in this report.

Figure 6 : Early, Prompt, Late synchronization [2]

2.3.2. THE LOCK LOOP

Lock loops are used to track the value of a parameter Φ (here code, phase or Doppler) over time by using a filtered measure of the error (thanks to the discriminator) of this parameter 𝜖Φ̂ .This measure is an input parameter whereas the filter is implemented inside the loop. The output parameter of the loop is the estimation of the value Φ.

The algorithm for the computation of Φ̂ is [1]_:

Φ̂(𝑛 + 1) = Φ̂(𝑛) + 𝑉(𝑛) Where 𝑉(𝑛) is the result of the filtering of 𝜖Φ̂.

Assuming that the error is altered by a noise 𝑁, the error 𝜖Φ̂ is equal to [1]_: 𝜖̂ (𝑛) = Φ(𝑛) − ΦΦ ̂ (𝑛) + 𝑁(𝑛)

The scheme of the lock loop is:

(19)

(20)

3. Propagation

3.1. THE IONOSPHERE

The ionosphere is the part of the atmosphere located between 80 and 1000 km of altitude. It is a neutral media composed of ionized gasses as well as free negatively charged electrons and positively charged ions. The presence of electrons in the atmosphere has an influence on electromagnetic waves. For high frequencies, such as the ones used in GNSS, the presence of electrons in the ionosphere leads to a modification of the group and phase velocities of the waves. This phenomenon adds a delay on the propagation time of the wave and an advance on the phase of the wave. More information is available in [7].

3.1.1. IONOSPHERE INFLUENCE ON ELECTROMAGNETIC WAVES.

The electrons present in the ionosphere have an influence on the time it takes for a wave to propagate through the ionosphere.

The time it takes for a wave to cover a distance 𝑑𝑠 at the speed 𝑣 is: 𝑑𝑡 =𝑑𝑠

𝑣 Along a path s, the propagation time 𝑡 is:

𝑡 = ∫1 𝑣𝑑𝑠

However, 𝑣 = 𝑐/𝑛, with 𝑐, the speed of light and 𝑛 the refraction index. Thus,

𝑡 =1

𝑐∫ 𝑛 ∗ 𝑑𝑠

The delay of propagation 𝜏 comparing to a wave traveling at the speed of light is: 𝜏 = ∫1 𝑣∗ 𝑑𝑠 − ∫ 1 𝑐∗ 𝑑𝑠 𝜏 =1 𝑐∗ ∫(𝑛 − 1) ∗ 𝑑𝑠

The first order of the Appleton-Hartree equation states that, the refractive indexes 𝑛𝑔 and 𝑛𝜙 respectively for the group and the phase of a signal are [6]_:

𝑛𝑔= 1 + 40.3

𝑓2 ∗ 𝑁𝑒 𝑛𝜙= 1 −

40.3 𝑓2 ∗ 𝑁𝑒 where 𝑓 is the frequency of the wave and 𝑁𝑒 is the electron density

Thus, the group delay is:

𝜏𝑔=1 𝑐∗

40.3

𝑓2 ∗ ∫ 𝑁𝑒∗ 𝑑𝑠 (𝑠) The ionospheric delay in meter for the code is:

𝐼𝑃= 40.3

𝑓2 ∗ ∫ 𝑁𝑒∗ 𝑑𝑠 (𝑚) The phase delay is:

𝜏𝜙= − 1 𝑐∗

40.3

𝑓2 ∗ ∫ 𝑁𝑒∗ 𝑑𝑠 (𝑠) The ionospheric delay in meter for the phase is:

𝐼𝜙= −40.3

𝑓2 ∗ ∫ 𝑁𝑒∗ 𝑑𝑠 (𝑚)

The absolutes values of the group and phase delays are the same but the signs are opposite: the group delay leads to an actual delay whereas the phase delay leads to an advance.

The value ∫ 𝑁𝑒∗ 𝑑𝑠 is called Total Electron Content (TEC), it corresponds to the number of electrons present over the path of a signal for a cross section of one square meter.

(21)

3.1.2. DESCRIPTION OF THE ELECTRON DENSITY OF THE IONOSPHERE

The electron density within the atmosphere is not constant but depends mostly on the altitude and the solar radiations. The ionosphere can be described by using three layers (see Figure 8): the D-layer with a low ionization of 109_{𝑒𝑙/𝑚}3_{, the E-layers with medium ionization closed to}₁₀11_{𝑒𝑙/𝑚}3_{and the F-layer with high} ionization around 1012_{𝑒𝑙/𝑚}3_{. Those values are for daytime. During nighttime, the values are lower as there is} no solar radiation.

Figure 8 : Electron density in the ionosphere [7]

Other parameters can influence the electron density:

- The time of the day : peaks at 2 pm in the afternoon - The season : peaks appear at spring and fall equinoxes, - The solar activities : the 11 year solar cycle and solar storms,

- The latitude: the electron density is higher at low latitudes as the solar radiations increase.

(22)

Figure 9 : Geographic fluctuation [7]

The dimension of the TEC parameter is electron per square meter but as its values are always above 1016_{𝑒𝑙/𝑚}2_{, a new unit called TEC unit (TECU) was defined:}

1 𝑇𝐸𝐶𝑈 = 1016_{𝑒𝑙/𝑚}2

Because of the many unpredictable parameters from which the TEC parameter depends on, it is not easy to forecast or model its values. Thus, it is difficult to estimate the ionospheric delays in advance. Nevertheless, in the section 3.1.4, two models of the ionosphere delay, the Klobuchar and NeQuick models are explained.

3.1.3. THIN LAYER ASSUMPTION

The Thin Layer Assumption is a way to model the ionosphere. This assumption is used in many ionosphere models such as the Klobuchar model (section 3.1.4.1). The idea is to represent the ionosphere as a thin shell layer located at an altitude of 350 km as presented in Figure 10. The height was chosen to be equal to the height of the maximum electron density in the ionosphere as seen in Figure 8.

The point where the signal crosses the single layer of ionosphere is called the Ionospheric Pierce Point (IPP). From this model, another parameter describing the signal going from the satellite to the receiver can be defined: the Vertical Total Electron Content (VTEC).

(23)

3.1.3.1. THE VTEC PARAMETER

The VTEC parameter of a signal is equal to the value of the TEC parameter of a signal having the same IPP but an elevation of 90°. For this reason, the parameter TEC can also be called STEC for Slant Total Electron Content. From now on, in order not to have any confusion, it will be written STEC and not TEC. Figure 11 presents the difference between the STEC and VTEC parameters.

Figure 11 : Difference between STEC and VTEC

3.1.3.2. THE OBLIQUITY FACTOR

The STEC value is dependent on the elevation of the satellite: the lower the elevation angle is, the longer the time it takes for the signal to go through the atmosphere. For example, in Figure 11, the signal sent from satellite 1 (in green) spends less time in the ionosphere as the one sent by satellite 2 (in orange). Thus, if the elevation decreases, the ionospheric delay increases. The obliquity factor is a factor depending only on the elevation that states the relation between the STEC and VTEC parameters.

The obliquity factor OF is equal to [6]_:

𝑂𝐹(𝜉) = (1 − (𝑅𝐸∗ sin(𝜉) 𝑅𝐸+ ℎ𝐼 ) 2 ) −1₂ Equation 3.1-3 Where:

- 𝑅𝐸= 6371 𝑘𝑚, the Earth radius, - ℎ𝐼= 350 𝑘𝑚, the mean isotropic height,

- 𝜉 = 90 − 𝑒, with 𝑒 the elevation angle of the satellite. Thus,

𝑆𝑇𝐸𝐶 = 𝑂𝐹(𝜉) ∗ 𝑉𝑇𝐸𝐶

Equation 3.1-4

3.1.4. EXAMPLES OF IONOSPHERE MODELS

3.1.4.1. KLOBUCHAR MODEL

The Klobuchar model [3]_{is based on the thin layer assumption. At a given point, the VTEC parameter is} modeled like a sinusoid function in time for which the amplitude A and period P depend on the position of the point and the time. The parameters A and P are modeled as polynomials for which the coefficients are broadcast in the navigation message of the signal.

In general, the Klobuchar model corrects only half of the ionospheric delay. It is usually more efficient for mid-latitudes. More information in [3]

3.1.4.2. NEQUICK MODEL

(24)

3.2. PRESENCE OF MULTIPATH

In case of multipath more than one signal coming from the same satellite reaches the receiver. Each of them has different parameters.

The model for the signal received by the receiver is [2]_: 𝑆𝑟𝑒𝑐𝑒𝑖𝑣𝑒𝑑(𝑡) = ∑ ∑ 𝐴_𝑗𝑖_𝑐_{𝑗(𝑡 − 𝜏} 𝑗𝑖) ∗ 𝑑𝑗(𝑡 − 𝜏𝑗𝑖) ∗ cos (2𝜋(𝑓0+ 𝜈𝑗𝑖)𝑡 + 𝜙𝑗𝑖) 𝑁_𝑗 𝑖=1 𝑁𝑠𝑎𝑡 𝑗=1 + 𝑛(𝑡) Equation 3.2-1 With:

- 𝑁𝑠𝑎𝑡, the number of visible satellite,

- 𝑁𝑗, the number of paths (LOS + multipath) for satellite 𝑆𝑗,

- 𝐴𝑗𝑖, the amplitude of the multipath 𝑖 signal coming from satellite 𝑆𝑗

- 𝜏𝑗𝑖, the propagation time of the multipath 𝑖 between satellite 𝑆𝑗 and the receiver,

- 𝜈𝑗𝑖, the Doppler effect of the multipath 𝑖 due to the motion of satellite 𝑆𝑗 relative to the receiver, - 𝜙𝑗𝑖, the phase delay of the multipath 𝑖 between satellite 𝑆𝑗 and the receiver,

- 𝑐𝑖, the PRN code of satellite 𝑆𝑖,

- 𝑑𝑖, the information regarding the navigation message sent by satellite 𝑆𝑖, - 𝑛, the total noise,

- 𝑓0, the carrier phase frequency.

The presence of multipath has effects on the DLL and PLL. Each multipath has an effect on the lock loops and leads to errors on the determination of the code’s and phase’s delays. To explain the effect of multipath, only one satellite will be taken into account for simplification 𝑆𝑗 (𝑁𝑠𝑎𝑡 = 1 in Equation 3.2-1) with one disturbed signal (𝑁1= 2 in Equation 3.2-1).The index 1 is for the direct path and the index 2 is for the other path.

𝑆𝑟𝑒𝑐(𝑡) = 𝐴_𝑗1_𝑐_{𝑗(𝑡 − 𝜏}

𝑗1)𝑑𝑗(𝑡 − 𝜏𝑗1)cos (2𝜋(𝑓0+ 𝜈𝑗1)𝑡 + 𝜙𝑗1) + 𝐴𝑗2𝑐𝑗(𝑡 − 𝜏𝑗2)𝑑𝑗(𝑡 − 𝜏𝑗2)cos (2𝜋(𝑓0+ 𝜈𝑗2)𝑡 + 𝜙𝑗2) + 𝑛(𝑡) Equation 3.2-2

3.2.1. EFFECT ON THE DLL

In presence of multipath, the formula of the parameter 𝐼 of the Figure 4 is the sum of the formulas of 𝐼 for each of the above component (see section 8 – Annex 1 for the details of the calculations).

𝐼(𝑘) = 𝐴𝑗 1_𝑑 𝑗(𝑘) 2 𝐾𝑐𝑖(𝜖𝜏_𝑗1) cos (𝜖𝜙_𝑗1) sinc (𝜋 (𝜖𝜈_𝑗1) 𝑇𝑖𝑛𝑡) + 𝐴𝑗2𝑑𝑗(𝑘) 2 𝐾𝑐𝑖(𝜖𝜏_𝑗2) cos (𝜖𝜙2_𝑗) sinc (𝜋 (𝜖𝜈_𝑗2) ∗ 𝑇𝑖𝑛𝑡) + 𝑛𝐼(𝑘) Equation 3.2-3 With: - 𝜖_𝜏 𝑗 𝑖 = 𝜏𝑗𝑖− 𝜏̂ 𝑗 - 𝜖_ϕ 𝑗 𝑖 = ϕ_𝑗𝑖 − 𝜙𝑗̂ - 𝜖_𝜈 𝑗 𝑖 = 𝜈_𝑗𝑖− 𝜈𝑗̂

In the discriminator of code, the presence of the autocorrelation of the PRN code will appear twice, at two different moments: 𝜖_𝜏 𝑗 1 and 𝜖_𝜏 𝑗 2. 𝜖_𝜈 𝑗 2 is equal to 𝜖_𝜈 𝑗

1 plus a positive delay 𝑑𝑀𝑃 due to the fact that the multipath

arrives after the direct path. There are two possibilities depending on the chip period 𝑇𝑐ℎ𝑖𝑝: - 𝑑𝑀𝑃> 𝑇𝑐ℎ𝑖𝑝, 𝐾𝑐_𝑖(𝜖_𝜈

𝑗2) ≈ 0 and the estimation of 𝜏𝑗̂(𝑘) is not biased by the presence of multipath,

- 0 < 𝑑𝑀𝑃< 𝑇𝑐ℎ𝑖𝑝, the estimation of 𝜏̂(𝑘) is biased by the presence of multipath. 𝑗

(25)

Figure 12 : Correlation function in presence of multipath [2]

Figure 12 is an example of the shape of the correlation function for two paths that are in phase. On the above picture, the synchronization is assumed perfect (i.e 𝜖_𝜏

𝑗1= 0). The goal of the DLL is to balance the Early and Late correlation output. Here, although the estimation is perfect (the one of the direct signal), the discriminator will measure 𝐾𝑐_𝑖(𝜖𝜏_𝑖−𝑑₂) ≥ 𝐾𝑐_𝑖(𝜖𝜏_𝑖+𝑑₂), which mean the synchronization is in advance.. Thus, the DLL will react to this wrong information, which will lead to a bias on the measurement of the code’s delay due to multipath.

(26)

Figure 13: Tracking error envelope of the code's delay [2]

3.2.2. EFFECT ON THE PLL

For a given difference of propagation time between the direct signal and the multipath 𝑑𝑀𝑃, the error induced on the PLL varies between a maximum and a minimum depending on the difference 𝜖_𝜙

𝑗2− 𝜖𝜙𝑗1.

The maximum is achieved when the two signals are in phase, i.e. 𝜖_𝜙 𝑗

2− 𝜖_𝜙

𝑗1= 0 [2𝜋] and the minimum is achieved when the two signals are in opposition of phase, i.e. 𝜖_𝜙

𝑗2− 𝜖𝜙𝑗1 = 𝜋 [2𝜋].

It is possible to plot the two functions maximum and minimum of the error induced on the PLL versus the multipath delay. Figure 14 is a plot of those functions for a multipath signal whose power is the same as the one of the direct signal. The maximum of error is equal to one-quarter of the wavelength, which corresponds to 5cm at L1.

(27)

3.3. THE OBSERVATION MODEL

In section 2.3, it is written that a receiver is composed of three main blocks. The main goal of the two first ones is to track the values of the code’s and phase’s delays. The last one must compute the PVT of the receiver thanks to these two values and the information message carried by the signal. As written in section 2.1.2, the receiver needs to have access to the distances between itself and the visible satellites, it must figure them out from the values of the delays.

3.3.1. THE OBSERVABLES

Both code’s and phase’s delays are due to the time the signal needs to propagate from the satellite to the receiver and thus depend on the distance between the receiver and the satellite. The units of the code’s and phase’s delays are respectively in time and radians. In order to compare them, they must be expressed in the same unit. The easiest way is to transform the delays in a distance expressed in meter:

- The code’s delay is multiplied by the speed of light 𝑐 to be changed from time to distance.

- The phase’s delay is multiplied by _2𝜋𝜆 (if it is expressed in radians, otherwise must first be conversed from the unit to radian) where 𝜆 is the wavelength of the wave in the vacuum.

Eventually, two observables 𝑃 and Φ, expressed in meter 𝑚, are obtained called code and phase. Numerically, 𝑃 = 𝑐 ∗ 𝜏 Equation 3.3-1 Φ = 𝜆 2𝜋∗ 𝜙 Equation 3.3-2 However, these parameters are not equal to the distance between the satellite and the receiver as they are biased by some effects due to for instance the propagation of the signal in the atmosphere.

3.3.2. PROPAGATION EFFECTS

The signal does not propagate in the vacuum but in the Earth’s atmosphere. Thus, the speed changes and it adds a delay. This delay is decomposed into several parts:

- The ionospheric delay presented in section 3.1.1

- The tropospheric delay: unlike the ionospheric delay, this one does not depend on the frequency of the wave.

- The presence of multipath due to the nearby environment.

The clock bias presented in section 2.1.1 has also an effect on the observables.

3.3.3. OBSEVATION MODEL

A model for each observable must be defined with the different phenomena disturbing the measure of the delays.

The observable model for one satellite 𝑆𝑥 at the time k and frequency 𝑦 is [5]:

𝑃𝑆𝑥 ,𝑦(𝑘) = 𝜌𝑆𝑥(𝑘) + 𝑐 (𝛿𝑇𝑆𝑥(𝑘) − 𝛿𝑡(𝑘)) + 𝑇𝑆𝑥(𝑘) + 𝐼𝑆𝑥 ,𝑦(𝑘) + 𝑀𝑃𝑃 ,𝑆𝑥 ,𝑦 (𝑘) + 𝑛𝑃 ,𝑆𝑥 ,𝑦(𝑘) + 𝑏𝑃 ,𝑆𝑥 ,𝑦(𝑘) Equation 3.3-3 Φ𝑆_𝑥 ,𝑦(𝑘) = 𝜌𝑆_𝑥(𝑘) + 𝑐 (𝛿𝑇𝑆_𝑥(𝑘) − 𝛿𝑡(𝑘)) + 𝑇𝑆_𝑥(𝑘) − 𝐼𝑆_𝑥 ,𝑦(𝑘) + 𝑀𝑃Φ ,𝑆_𝑥 ,𝑦 (𝑘) + 𝑛Φ ,𝑆_𝑥 ,𝑦(𝑘) + 𝑏Φ ,𝑆_𝑥 ,𝑦(𝑘) + 𝜆𝑦𝐴𝑆_𝑥 ,𝑦 Equation 3.3-4 With,

(28)

- Φ𝑆_𝑥 ,𝑦, the phase’s delay in meter

- 𝜌𝑆𝑥, the geometric distance between the receiver and the satellite Sx

- 𝛿𝑇𝑆_𝑥,, 𝛿𝑡𝑆_𝑥, respectively the satellite and receiver clock error (see section 2.1.1) - 𝑇𝑆𝑥, the tropospheric delay

- 𝐼𝑆_𝑥 ,𝑦, the ionospheric delay (see section 3.1)

- 𝑀𝑃𝑃, 𝑀𝑃Φ, respectively the error induced by multipath on the code and on the phase - 𝑛𝑃, 𝑛Φ, respectively the thermic noise on the code and on the phase

- 𝑏𝑃, 𝑏Φ, respectively both the satellite and the receiver clock bias at the frequency Y - 𝐴𝑦, ambiguity factor (see section 2.3.1.2)

- 𝜆𝑦, wave length of the signal at the frequency Y

The difference between the satellite and receiver clock errors 𝛿𝑇𝑆𝑥,, 𝛿𝑡𝑆𝑥 and the satellite and receiver clock biases 𝑏𝑃, 𝑏Φ is that the first ones do not depend on the frequency whereas the second ones depend on it. It is possible to gather the terms that are independent of the frequency:

𝐷𝑆_𝑥(𝑘) = 𝜌𝑆𝑥(𝑘) + 𝑐 (𝛿𝑇𝑆𝑥(𝑘) − 𝛿𝑡𝑆𝑥(𝑘)) + 𝑇𝑆𝑥(𝑘) Equation 3.3-5 Thus, 𝑃𝑆_𝑥 ,𝑦(𝑘) = 𝐷𝑆_𝑥(𝑘) + 𝐼𝑆_𝑥 ,𝑦(𝑘) + 𝑀𝑃𝑃 ,𝑆_𝑥 ,𝑦 (𝑘) + 𝑛𝑃 ,𝑆_𝑥 ,𝑦(𝑘) + 𝑏𝑃 ,𝑆_𝑥 ,𝑦(𝑘) Equation 3.3-6 Φ𝑆_𝑥 ,𝑦(𝑘) = 𝐷𝑆_𝑥(𝑘) − 𝐼𝑆_𝑥 ,𝑦(𝑘) + 𝑀𝑃Φ ,𝑆_𝑥 ,𝑦 (𝑘) + 𝑛Φ ,𝑆_𝑥 ,𝑦(𝑘) + 𝑏Φ ,𝑆_𝑥 ,𝑦(𝑘) + 𝜆𝑦𝐴𝑆_𝑥 ,𝑦(𝑘) Equation 3.3-7

As the goal of the study is to estimate the ionosphere, for each equation, there are three sources of errors, the multipath, the noise and the clocks bias. It is possible to gather them under a term called error (the index 𝑖 represents the code or the phase), 𝐸:

𝐸𝑖,𝑆_𝑥 ,𝑦(𝑘) = 𝑀𝑃𝑖 ,𝑆_𝑥 ,𝑦 (𝑘) + 𝑛𝑖 ,𝑆_𝑥 ,𝑦(𝑘) + 𝑏𝑖 ,𝑆_𝑥 ,𝑦(𝑘) Equation 3.3-8 Thus, 𝑃𝑆𝑥 ,𝑦(𝑘) = 𝐷𝑆𝑥(𝑘) + 𝐼𝑆𝑥 ,𝑦(𝑘) + 𝐸𝑃,𝑆𝑥 ,𝑦(𝑘) Equation 3.3-9 Φ𝑆_𝑥 ,𝑦(𝑘) = 𝐷𝑆𝑥(𝑘) − 𝐼𝑆𝑥 ,𝑦(𝑘) + 𝜆𝑦𝐴𝑆𝑥 ,𝑦(𝑘) + 𝐸𝜙,𝑆𝑥 ,𝑦(𝑘) Equation 3.3-10 Standard deviations

When we assume the observable are perturbed only by noises, the expressions of the standard deviations of the code’s and phase’s delays are available (needed later in this study). Those values depend on the parameters of the tracking mode of the receiver and of the signal power indicator. The formulas are [5]_:

𝜎𝑃,𝑛𝑜𝑖𝑠𝑒= 𝑐 𝑓𝑃𝑅𝑁 √( 𝐵𝐿𝐷𝐿𝐿𝑑 2𝛼_𝑁0𝐶 (1 + 1 𝐶 𝑁0𝑇𝑖𝑛𝑡 ) ) (𝑚) Equation 3.3-11 𝜎Φ,noise,Y=𝜆𝑌 2𝜋 √ 𝐵𝐿𝑃𝐿𝐿 𝐶 𝑁0 (𝑚) Equation 3.3-12 Where:

(29)

- 𝑇𝑖𝑛𝑡 is the integration time in s

- 𝑑 is the DLL correlator spacing in chip - c is the speed of light

- 𝑓𝑃𝑅𝑁 cadency of the PRN code (1.023 MHz for the code’s delay)

- 𝛼 is the slope of the normalized spreading sequence autocorrelation function in 𝑑/2, - 𝐶

𝑁0 is the carrier to noise PSD ratio in dB-Hz

- 𝜆𝑌 is the wavelength of the signal in m (𝜆1= 0.1903 𝑚 for 𝐿1 and 𝜆2= 0.1898 𝑚 for 𝐿2) The order of magnitude of the carrier to noise _𝑁𝐶

0 for the transmission of a signal between a receiver and a satellite in normal conditions (no degradations of the signal in the atmosphere) is 50 𝑑𝐵𝐻𝑧. Thus, the orders of magnitude of the standard deviations for the code’s and phase’s delays for the receiver used in this study (presented in section 5.1.2) are:

𝜎𝑃,𝑛𝑜𝑖𝑠𝑒≈ 1 𝑚 Equation 3.3-13 𝜎Φ,noise,L₁≈ 0.3 𝑚𝑚 Equation 3.3-14 𝜎Φ,noise,L2≈ 0.4 𝑚𝑚 Equation 3.3-15 The PLL is indeed more precise than the DLL.

3.4. PROBLEM

As stated in the introduction, the signal transmitted between a satellite and a receiver is altered by different factors such as delays when spreading through the troposphere and ionosphere. As seen in section 3.1.1, the ionospheric delays are frequency dependent. If one combines Equation 3.3-9 and Equation 3.3-10 at two different frequencies, it is possible to have access to an estimation of the ionospheric delay (more details in section 4.1). However, in presence of multipath, the estimations of the ionospheric delays are biased. My study focuses on establishing multi-frequency algorithms in order to estimate the ionospheric delays and testing them in different environment to see if they are robust to multipath. To do so, for each observation model and for each environment, two sets of simulations are done: one with multipath and one without. Then, the errors on the ionospheric delays are estimated in both cases and compares to see if the model is robust to multipath.

The next section presents the different estimation strategy used in this study. Each solution is composed of one observation vector 𝑌 with the observables used for this model, one state vector 𝑋 with the unknown parameters, an observation matrix 𝐻 and a noise 𝑁.

The equation of the model is:

𝑌 = 𝐻 ∗ 𝑋 + 𝑁 The Least Square (LS) solution of this equation is:

𝑋 = (𝐻𝑇_{∗ 𝐻)}−1_𝐻𝑇_{∗ 𝑌}

(30)

4. Techniques to estimate the ionospheric delays

More information about this section is available in [5]. All the equations are taken from this paper.

4.1. DUAL AND SINGLE FREQUENCY MEASURMENTS

As seen in section 3.1.1, the absolute value of the ionospheric delays is: 𝐼 =40.3

𝑓2 ∗ 𝑆𝑇𝐸𝐶

)

Equation 4.1-2

The difference of the phase’s delays multiplied by a relevant factor is equal to the STEC value biased by the difference of the errors on both frequencies and the difference of ambiguities.

4.1.3. COMPARISON OF BOTH PHASE’S AND CODE’S DELAY

If we compare both the code’s and phase’s delays at the same frequency Y,

𝑘𝑌( 𝑃𝑆_𝑥 ,𝑌(𝑘)

−

Φ𝑆_𝑥 ,𝑌(𝑘)) = 𝑆𝑇𝐸𝐶𝑆_𝑥(𝑘) − 𝑘𝑌(𝐸𝜙,𝑆_𝑥 ,𝑌(𝑘) − 𝐸𝑃,𝑆_𝑥 ,𝑌(𝑘) + 𝜆𝑦𝐴𝑆_𝑥 ,𝑦)

(31)

4.1.4. OBSERVABLE POSSIBLE

According to the three previous sections, there are four ways to estimate the values of the STEC parameter from the values of the code’s and the phase’s delays both at L1 and at L2 : difference of the code’s delays, difference of the phase’s delays, difference of the code’s delay and the phase’s delay at L1 and at L2. Those four differences are possible observables, the notations used to represent them are:

Δ𝑃(𝑘) =

𝑘

12[𝑃𝑆𝑥 ,1(𝑘) − 𝑃𝑆𝑥 ,2(𝑘)], difference of code’s delays ΔΦ(k) =

𝑘

12[Φ𝑆𝑥 ,2(𝑘) − Φ𝑆𝑥 ,1(𝑘)], difference of phase’s delays 𝐶𝑀𝐶1(𝑘) = 𝑘1[ 𝑃𝑆_𝑥 ,1(𝑘)

−

Φ𝑆_𝑥 ,1(𝑘)], difference of code’s and phase’s at L1 𝐶𝑀𝐶2(𝑘) = 𝑘2[ 𝑃𝑆𝑥 ,2(𝑘)

−

Φ𝑆𝑥 ,2(𝑘)], difference of code’s and phase’s at L2

In the observation vector, the differences are multiplied by the relevant factor k in order to simplify the expression of the observation matrix.

These four differences are not exactly equal to the ionospheric delays, the error terms (see Equation 3.3-8) due to multipath and clock biases are still present.

Table 1 presents the different advantages and disadvantages of each observable.

Strengths Weaknesses

Difference of code’s delays Without ambiguity Need of two frequencies Less precise (order of m) Difference of phase’s delays More precise (order of mm) With ambiguity

Need of two frequencies Difference of code’s and phase’s Only one frequency needed With ambiguity

Less precise (order of m) Table 1 : Strengths and weaknesses of the different observables

4.2. ESTIMATION STRATEGY

The first part of my study was to choose the observation model I will use to compute the ionospheric delays and test if they are robust to multipath or not. This section presents the models I have defined as well as the reasons why. Except for the local ionospheric model presented in section 4.2.4.1, the models presented in this section are my own derivation and I have implemented them myself in order to run the simulations

4.2.1. MODEL 0

Model 0 is the simplest model that allows estimating the STEC parameter. It only uses one satellite at a time and one observable, the difference of code’s delays.

For this model, the observation vector is:

𝑌(𝑘) = Δ𝑃(𝑘)

Equation 4.2-1 The state vector is:

𝑋(𝑘) = 𝑆𝑇𝐸𝐶(𝑘)

Equation 4.2-2 The observation matrix is:

𝐻 = 1

Equation 4.2-3 The main advantage of this model is that it is really easy to implement and allows a fast estimation of the ionospheric delays. This model is used as a reference model.

4.2.2. MODEL 1

(32)

𝑌(𝑘) = ( Δ𝑃(𝑘) ΔΦ(k) 𝐶𝑀𝐶1(𝑘) 𝐶𝑀𝐶2(𝑘) ) Equation 4.2-4

The state vector is:

𝑋(𝑘) = ( 𝑆𝑇𝐸𝐶(𝑘) 𝑏𝑖𝑎𝑠Φ(𝑘) 𝑏𝑖𝑎𝑠𝐶𝑀𝐶1(𝑘) 𝑏𝑖𝑎𝑠𝐶𝑀𝐶₂(𝑘)) Equation 4.2-5

The observation matrix is:

𝐻(𝑘) = ( 1 0 0 0 1 𝑘12 0 0 1 0 𝑘1 0 1 0 0 𝑘2 ) Equation 4.2-6 This model uses the four observables possible so it uses all the information available for each satellite. In case of the loss of one frequency, there is still one observable possible and thus, information regarding the changes of the 𝑆𝑇𝐸𝐶 parameter. However, in this model, all the satellites are analyzed separately so it is not possible to compute the value of the 𝑆𝑇𝐸𝐶 parameter of a satellite with respect to another one.

4.2.3. MODEL 2

Model 2 uses one satellite at a time and only two observables for each satellite: the code’s and phase’s difference delays. As presented in section 4.1, there are two unknown parameters: the STEC parameter and a bias due to the ambiguities on the measurement of the phase’s delay.

Thus, the observation vector is:

𝑌(𝑘) = (_ΔΦ(k))Δ𝑃(𝑘)

Equation 4.2-7

The state vector is:

𝑋(𝑘) = (𝑆𝑇𝐸𝐶(𝑘)_{𝑏𝑖𝑎𝑠 (𝑘) )}

Equation 4.2-8

The observation matrix is:

𝐻(𝑘) = (1_{1 −𝑘}0 12)

Equation 4.2-9

Model 2 uses only two observables so it is expected to be more precise than model 1 as there is less incertitude. As for model 1, all the satellites are analyzed separately. In case of the loss of frequency of a satellite, this model does not allow to keep track of the changes of the 𝑆𝑇𝐸𝐶 parameter.

(33)

4.2.4. MODEL 3

In the three previous models, the satellites were analyzed one by one. However, all the satellites seen by a receiver are not necessarily impacted by multipath. For example, as presented in Figure 15, satellite 1 is not impacted by multipath whereas satellite 2 is impacted by them. The idea of model 3 is to figure out a model for the ionosphere where not impacted satellites can help impacted ones to estimate the ionospheric delays.

Figure 15: Situation with two satellites

4.2.4.1. LOCAL IONOSPHERIC MODEL

The model chosen for the ionosphere is a local model of the VTEC parameter stating that the VTEC factor evolves linearly in the North, South, East and West directions. Thus, the VTEC parameter of a point close to the user is equal to the VTEC of the user plus a factor depending on the relative position of the point from the user in the North, South, East and West directions. This model was chosen for this study after reading [5]. Mathematically, the value of the VTEC parameter at a point P (see Figure 16 for the position of the point P and the receiver) is equal to:

𝑉𝑇𝐸𝐶𝑝(𝑘) = ( 𝑉𝑇𝐸𝐶𝑢(𝑘) + max(𝑙𝑎𝑡𝑝𝑆_{(𝑘) − 𝑙𝑎𝑡} 𝑢(𝑘), 0) ∗ 𝑔𝑁(𝑘) + min(𝑙𝑎𝑡𝑝𝑆_{(𝑘) − 𝑙𝑎𝑡} 𝑢(𝑘), 0) ∗ 𝑔𝑠(𝑘) + max(𝑙𝑜𝑛𝑔𝑝𝑆(𝑘) − 𝑙𝑜𝑛𝑔𝑢(𝑘), 0) ∗ 𝑔𝐸(𝑘) + min(𝑙𝑜𝑛𝑔𝑝𝑆(𝑘) − 𝑙𝑜𝑛𝑔𝑢(𝑘), 0) ∗ 𝑔𝑊(𝑘)) Where:

- 𝑉𝑇𝐸𝐶𝑢 is the VTEC value of the receiver,

- 𝑙𝑎𝑡𝑝𝑆(𝑘) 𝑎𝑛𝑑 𝑙𝑎𝑡𝑢(𝑘) are respectively the point and receiver latitudes, - 𝑙𝑜𝑛𝑔𝑝𝑆(𝑘) 𝑎𝑛𝑑 𝑙𝑜𝑛𝑔𝑢(𝑘) are respectively the point and receiver longitudes.

- 𝑔𝑁, 𝑔𝑆, 𝑔𝐸, 𝑔𝑊 are the gradient respectively in the North, South, East and West directions. The following notations are used from now on in this study:

- Δ𝑙𝑎𝑡𝑆,𝑢,𝑝= 𝑙𝑎𝑡𝑝𝑆− 𝑙𝑎𝑡𝑢

- Δ𝑙𝑜𝑛𝑔𝑆,𝑢,𝑝= 𝑙𝑜𝑛𝑔𝑝𝑆(𝑘) − 𝑙𝑜𝑛𝑔𝑢(𝑘)

(34)

only five unknowns to compute all the values of the ionospheric delays. The number of unknowns to compute the ionospheric delays is diminished if there are at least six satellites seen by the receiver.

Figure 16 : Representation of the local ionospheric model

4.2.4.2. MATRICES OF MODEL 3

Model 3 uses all the satellites available at a time and only two observables for each satellite: the difference of code’s and of phase’s. It is based on the local ionospheric model.

(35)

𝑋(𝑘) = ( 𝑉𝑇𝐸𝐶𝑢(𝑘) 𝑔𝑁(𝑘) 𝑔𝑆(𝑘) 𝑔𝐸(𝑘) 𝑔𝑊(𝑘) 𝑏𝑖𝑎𝑠𝑁₁(𝑘) : 𝑏𝑖𝑎𝑠𝑁_𝑖(𝑘) : 𝑏𝑖𝑎𝑠𝑁𝑛(𝑘)) Equation 4.2-11 The observation matrix is:

𝐻(𝑘) = (𝑂𝐹𝑛∗1(𝑘) Δ𝑙𝑎𝑡𝑛∗2(𝑘) Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) 0𝑛∗𝑛 𝑂𝐹𝑛∗1(𝑘) Δ𝑙𝑎𝑡𝑛∗2(𝑘) Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) 𝐼𝑛 )

Equation 4.2-12 With,

- 𝑂𝐹𝑛∗1(𝑘), the matrix with the obliquity factor of all the visible satellites.

𝑂𝐹𝑛∗1(𝑘) = ( 𝑂𝐹𝑁1(𝑘) : 𝑂𝐹𝑁_𝑖(𝑘) : 𝑂𝐹𝑁𝑛(𝑘)) Equation 4.2-13 - Δ𝑙𝑎𝑡𝑛∗2 (𝑘) the matrix of the difference of latitude between the receiver and all the satellite

Δ𝑙𝑎𝑡𝑛∗2(𝑘) = 𝑂𝐹𝑁1∗ max (Δ𝑙𝑎𝑡𝑁1,𝑢,𝑝1, 0) 𝑂𝐹𝑁1∗ min (Δ𝑙𝑎𝑡𝑁1,𝑢,𝑝1, 0) : : 𝑂𝐹𝑁_𝑖∗ max (Δ𝑙𝑎𝑡𝑁_𝑖,𝑢,𝑝_𝑖, 0) 𝑂𝐹𝑁_𝑖∗ min (Δ𝑙𝑎𝑡𝑁_𝑖,𝑢,𝑝_𝑖, 0) : : 𝑂𝐹𝑁_𝑛∗ max (Δ𝑙𝑎𝑡𝑁_𝑛,𝑢,𝑝𝑛, 0) 𝑂𝐹𝑁𝑛∗ min (Δ𝑙𝑎𝑡𝑁𝑛,𝑢,𝑝𝑛, 0) Equation 4.2-14

- Δl𝑜𝑛𝑔𝑛∗2(𝑘) the matrix of the difference of longitude between the receiver and all the satellite

Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) = ( 𝑂𝐹𝑁1∗ max (Δ𝑙𝑜𝑛𝑔𝑁1,𝑢,𝑝1, 0) 𝑂𝐹𝑁1∗ min (Δ𝑙𝑜𝑛𝑔𝑁1,𝑢,𝑝1, 0) : : 𝑂𝐹𝑁_𝑖∗ max (Δ𝑙𝑜𝑛𝑔𝑁_𝑖,𝑢,𝑝_𝑖, 0) 𝑂𝐹𝑁_𝑖∗ min (Δ𝑙𝑜𝑛𝑔𝑁_𝑖,𝑢,𝑝_𝑖, 0) : : 𝑂𝐹𝑁𝑛∗ max (Δ𝑙𝑜𝑛𝑔𝑁𝑛,𝑢,𝑝𝑛, 0) 𝑂𝐹𝑁𝑛∗ min (Δ𝑙𝑜𝑛𝑔𝑁𝑛,𝑢,𝑝𝑛, 0)) Equation 4.2-15

The obliquity factors of all the visible satellites are known by using the approximate position of the receiver given in the observation file and the positions of all the satellites given in the navigation message.

The number of unknowns in this model is 𝑛 + 5 while the number of observable is 2𝑛. To be able to inverse the system, the number of observables must be superior to the number of unknowns. Thus

2𝑛 ≥ 𝑛 + 5

Equation 4.2-16 𝑛 ≥ 5

(36)

4.2.5. MODEL 4

Model 4 uses all the satellites available at a time, the same local ionospheric model as in Model 3, and the four observables available for each satellite. The unknowns are the five unknowns of the local ionospheric model and, for each visible satellite, three biases due to the ambiguities on the measurement of the phase’s delay. If 𝑛 is the number of satellite available and 𝑁1, … , 𝑁𝑖, … . , 𝑁𝑛 the numbers of these satellites at a time 𝑘; the state and observation vectors are:

𝑌(𝑘) = ( Δ𝑃𝑁1(𝑘) : Δ𝑃𝑁_𝑖(𝑘) : Δ𝑃𝑁_𝑛(𝑘) ΔΦ𝑁1(k) : ΔΦ𝑁_𝑖(k) : ΔΦ𝑁_𝑛(k) 𝐶𝑀𝐶𝑌1,𝑁1(𝑘) : 𝐶𝑀𝐶𝑌₁,𝑁_𝑖(𝑘) : 𝐶𝑀𝐶𝑌_1,𝑁𝑛(𝑘) 𝐶𝑀𝐶𝑌2,𝑁1(𝑘) : 𝐶𝑀𝐶𝑌2,𝑁𝑖(𝑘) : 𝐶𝑀𝐶𝑌₂,𝑁_𝑛(𝑘)) Equation 4.2-18 And 𝑋(𝑘) = ( 𝑉𝑇𝐸𝐶𝑢(𝑘) 𝑔𝑁(𝑘) 𝑔𝑆(𝑘) 𝑔𝐸(𝑘) 𝑔𝑊(𝑘) 𝑏𝑖𝑎𝑠Φ,𝑁1(𝑘) : 𝑏𝑖𝑎𝑠Φ,𝑁_𝑛(𝑘) 𝑏𝑖𝑎𝑠𝐶𝑀𝐶1,𝑁1(𝑘) : 𝑏𝑖𝑎𝑠𝐶𝑀𝐶1,𝑁𝑛(𝑘) 𝑏𝑖𝑎𝑠𝐶𝑀𝐶2,𝑁1(𝑘) : 𝑏𝑖𝑎𝑠𝐶𝑀𝐶2,𝑁𝑛(𝑘)) Equation 4.2-19

The observation matrix is: 𝐻(𝑘) = ( 𝑂𝐹𝑛∗1(𝑘) Δ𝑙𝑎𝑡𝑛∗2(𝑘) Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) 0𝑛∗𝑛 0𝑛∗𝑛 0𝑛∗𝑛 𝑂𝐹𝑛∗1(𝑘) Δ𝑙𝑎𝑡𝑛∗2(𝑘) Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) 𝐼𝑛 0𝑛∗𝑛 0𝑛∗𝑛 𝑂𝐹𝑛∗1(𝑘) Δ𝑙𝑎𝑡𝑛∗2(𝑘) Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) 0𝑛∗𝑛 𝐼𝑛 0𝑛∗𝑛 𝑂𝐹𝑛∗1(𝑘) Δ𝑙𝑎𝑡𝑛∗2(𝑘) Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) 0𝑛∗𝑛 0𝑛∗𝑛 𝐼𝑛 ) Equation 4.2-20 With,

(37)

𝑂𝐹𝑛∗1(𝑘) = ( 𝑂𝐹𝑁1(𝑘) : 𝑂𝐹𝑁_𝑖(𝑘) : 𝑂𝐹𝑁𝑛(𝑘)) Equation 4.2-21 - Δ𝑙𝑎𝑡𝑛∗2 (𝑘) the matrix of the difference of latitude between the receiver and all the satellite

Δ𝑙𝑎𝑡𝑛∗2(𝑘) = ( 𝑂𝐹𝑁1∗ max (Δ𝑙𝑎𝑡𝑁1,𝑢,𝑝1, 0) 𝑂𝐹𝑁1∗ min (Δ𝑙𝑎𝑡𝑁1,𝑢,𝑝1, 0) : : 𝑂𝐹𝑁_𝑖∗ max (Δ𝑙𝑎𝑡𝑁_𝑖,𝑢,𝑝_𝑖, 0) 𝑂𝐹𝑁_𝑖∗ min (Δ𝑙𝑎𝑡𝑁_𝑖,𝑢,𝑝_𝑖, 0) : : 𝑂𝐹𝑁_𝑛∗ max (Δ𝑙𝑎𝑡𝑁_𝑛,𝑢,𝑝_𝑛, 0) 𝑂𝐹𝑁_𝑛∗ min (Δ𝑙𝑎𝑡𝑁_𝑛,𝑢,𝑝𝑛, 0)) Equation 4.2-22

- Δl𝑜𝑛𝑔𝑛∗2(𝑘) the matrix of the difference of longitude between the receiver and all the satellite

Δ𝑙𝑜𝑛𝑔𝑛∗2(𝑘) = ( 𝑂𝐹𝑁1∗ max (Δ𝑙𝑜𝑛𝑔𝑁1,𝑢,𝑝1, 0) 𝑂𝐹𝑁1∗ min (Δ𝑙𝑜𝑛𝑔𝑁1,𝑢,𝑝1, 0) : : 𝑂𝐹𝑁_𝑖∗ max (Δ𝑙𝑜𝑛𝑔𝑁_𝑖,𝑢,𝑝_𝑖, 0) 𝑂𝐹𝑁_𝑖∗ min (Δ𝑙𝑜𝑛𝑔𝑁_𝑖,𝑢,𝑝_𝑖, 0) : : 𝑂𝐹𝑁_𝑛∗ max (Δ𝑙𝑜𝑛𝑔𝑁_𝑛,𝑢,𝑝_𝑛, 0) 𝑂𝐹𝑁_𝑛∗ min (Δ𝑙𝑜𝑛𝑔𝑁_𝑛,𝑢,𝑝𝑛, 0)) Equation 4.2-23 As for model 3, this model can only be used if the number of visible satellite is above 5.

4.3. IMPROVEMENT OF THE MODELS

4.3.1. MODEL 0

As there is only one observation parameter and one state parameter, it is not relevant to use a Kalman filter on this model. The impact of multipath on the observation parameter can be seen as a noise that sometimes is positive and sometimes is negative. By averaging the values of the observation parameter, the positive and negatives impacts of the multipath cancel themselves and the values of the ionospheric delays are more accurate. Thus, the filter chosen to improve this model is simply a mean filter that averages the value of the observation parameter at a time k with M the previous ones.

The value chosen for M results of a compromise between the improvement of the accuracy added by the filter and its costs in terms of memory and speed. In this study, it has been chosen to average the 𝑆𝑇𝐸𝐶 values over 50 𝑠. As these values are computed every second the value, 𝑀 = 50.

𝑆𝑇𝐸𝐶𝑓𝑖𝑙𝑡𝑒𝑟(𝑘) = 1 𝑀∑ Δ𝑃(𝑘 − 𝑖) 𝑀−1 𝑖=0 Equation 4.3-1

4.3.2. MODEL 1

The accuracy of the values of the state vector of model 1 can be improved by using a Kalman filter. It is presented in Annex 2: Kalman filter. As presented in this annex, in order to use a Kalman filter, new parameters must be defined: the transition matrix and the covariance matrices of the observation and system errors. Also, the filter needs to be initialized.

4.3.2.1. EXPRESSIONS OF THE MATRICES

For a given path, the values of the STEC parameters change slowly so it is assumed that the STEC parameter is constant between two close consecutive times.

The ambiguity factor on the phase’s delay is constant over a period, then changes suddenly and is again constant over a period. Thus, it is supposed in the system model that the values of the biases factors are constant between two consecutive times.

(38)

𝐴 = ( 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 ) Equation 4.3-2

The covariance matrix of the observation error at any time 𝑡𝑘 depends on the values of the standard deviations of the code’s and phase’s delays given in equations Equation 3.3-11 and Equation 3.3-12. Its expression is demonstrated in section 9.2. 𝑅𝑘 = ( 𝑘2_(𝜎 𝑐𝑜𝑑𝑒,𝐿1 2 _{+ 𝜎} 𝑐𝑜𝑑𝑒,𝐿2 2 ₎ ₀ _{𝑘𝑘1𝜎} 𝑐𝑜𝑑𝑒,𝐿1 2 _{𝑘𝑘2𝜎} 𝑐𝑜𝑑𝑒,𝐿2 2 0 𝑘2_(𝜎 𝑝ℎ𝑎𝑠𝑒,𝐿2 ₁+ 𝜎𝑝ℎ𝑎𝑠𝑒,𝐿2 ₂) 𝑘1𝑘2𝜎𝑝ℎ𝑎𝑠𝑒,𝐿2 ₁ 𝑘1𝑘2𝜎𝑝ℎ𝑎𝑠𝑒,𝐿2 ₂ 𝑘𝑘1𝜎𝑐𝑜𝑑𝑒,𝐿2 1 𝑘1𝑘2𝜎𝑝ℎ𝑎𝑠𝑒,𝐿2 1 𝑘12(𝜎𝑐𝑜𝑑𝑒,𝐿2 1+ 𝜎𝑝ℎ𝑎𝑠𝑒,𝐿2 1) 0 𝑘𝑘1𝜎𝑐𝑜𝑑𝑒,𝐿2 ₂ _{𝑘1𝑘2𝜎𝑝ℎ𝑎𝑠𝑒,𝐿} 2 2 ₀ _𝑘22_{(𝜎𝑐𝑜𝑑𝑒,𝐿} 2 2 _{+ 𝜎𝑝ℎ𝑎𝑠𝑒,𝐿} 2 2 ₎₎ Equation 4.3-3 Section 3.3.3 gives orders of magnitude for the standard deviations of the code’s and phase’s delays . Using these values, the order of magnitude for the matrix 𝑅𝑘 (with the values introduced in section 3.3.3) is:

𝑅𝑘 = ( 150 0 25 15 0 2 ∗ 10−5 _{5 ∗ 10}−7 _{9 ∗ 10}−7 25 5 ∗ 10−7 ₈ ₀ 15 9 ∗ 10−7 ₀ ₃ ) (𝑚2₎ Equation 4.3-4 The values of the covariance matrix of the system error 𝑄𝑘 at any time 𝑡𝑘 cannot be calculated or estimated before running the filter. Here we assume that the three parameters of the state vector are independent so 𝑄𝑘 is diagonal. It is also assumed in this study that no matter if the ambiguities come from 𝐿1 or 𝐿2, the variances of these parameters are the same.

Finally, the expression of 𝑄𝑘 is:

𝑄𝑘= ( 𝑞𝑆𝑇𝐸𝐶 0 0 0 0 𝑞𝑝ℎ𝑎𝑠𝑒 0 0 0 0 𝑞𝑏 0 0 0 0 𝑞𝑏 ) Equation 4.3-5 With, 𝑞𝑆𝑇𝐸𝐶 and 𝑞𝑏 being two parameters to tune.

4.3.2.2. INITIALIZATION OF THE FILTER.

As presented in Annex 2: Kalman filter, the values of the state vector and its covariance matrix must be initialized. To do so, it has been chosen in this study to use the values of the state vector obtained thanks to the observation model presented in section Model 1.

𝑋𝑘𝑎𝑙𝑚𝑎𝑛𝑓𝑖𝑙𝑡𝑒𝑟 (𝑘 = 0) = 𝑋(𝑘 = 0) = ( 𝐻𝑡 _{∗ 𝐻)}−1 𝑡_𝐻_{∗ 𝑌(𝑘 = 0)} Its covariance is initialized to a diagonal matrix,

𝑃(𝑘 = 0) = 𝑝 (1 0 00 1 0 0 0 1

) With, 𝑝 being a parameter to tune.

4.3.3. MODEL 2

4.3.3.1. EXPRESSIONS OF THE MATRICES

Model 2 is closed to Model 1. For the same reasons as the ones described section 4.3.2, the transition matrix is equal to the identity matrix.

𝐴 = (1 0 0 1)