# Ionospheric Delay Handling for Relative Navigation by Carrier-Phase Differential GPS

Authors: A. Renga; U. Tancredi; M. Grassi
Journal: International Journal of Aerospace Engineering
 (2015)

---

## Abstract

The paper investigates different solutions for ionospheric delay handling in high accuracy long baseline relative positioning by Carrier-Phase Differential GPS (CDGPS). Standard literature approaches are reviewed and the relevant limitations are discussed. Hence, a completely ionosphere-free approach is proposed, in which the differential ionospheric delays are cancelled out by combination of dual frequency GPS measurements. The performance of this approach is quantified over real-world spaceborne GPS data made available by the Gravity Recovery and Climate Experiment (GRACE) mission and compared to the standard solution.


---

## Body

## 1. Introduction

Carrier-Phase Differential GPS (CDGPS) is a proven technology in several fields of application. CDGPS has been already employed for relative positioning of Low Earth Orbit (LEO) satellites flying in formation [1–4], of aircraft with respect to runways [5] and for cooperative self-separation of general aviation aircraft [6]. The capability to achieve high accuracy by CDGPS is based on the possibility to exploit the integer nature of Double Difference (DD) carrier-phase ambiguities [7]. However, as the separation among the satellites increases, the correlation of ionospheric delays among the receivers decreases [8]. As a result, DD GPS observables are affected by significant errors that complicate the integer resolution task. This paper investigates the effects of different strategies for ionospheric delay compensation on the accuracy in the relative positioning of GPS receivers in LEO over long baselines. Its results can be extended, at least in principle, to other formation flying applications [9], such as those involving Very Light Jets [10] and/or Unmanned Aerial Systems [11–13].Different approaches exist in the literature for dealing with ionospheric delays. In high accuracy, postprocessing applications with dual frequency data [14, 15], the DD ionospheric delays are estimated within a dynamic filter [16], for example, the Extended Kalman Filter (EKF), and are modelled by very simple stochastic models, typically using random walk processes in the filter’s state vector. As an alternative, delays are modelled by Lear’s model [1, 17] which allows relating the n slant ionospheric delays to the Vertical Total Electron Content (VTEC) above the receivers. Even though modelling the ionospheric delays helps to increase their observability, and thus to aid in the ambiguity resolution task, Lear’s model is known to be structurally capable of reproducing actual ionospheric delays only to a limited extent [8].It is well known that linear combinations of GPS dual frequency measurements can be used to delete the ionospheric delays to first order. When those combinations are used as an input to the EFK, it can be expected that the magnitude of the ionospheric delays does not affect the achievable relative positioning accuracy. However, the use of an ionosphere-free approach is known to complicate significantly the integer ambiguities resolution task if compared to approaches that attempt to model the ionosphere [18, 19]. Hence, the choice between model-based methods (e.g., based on Lear’s model) and ionosphere-free approaches is in general not trivial [20], depending on the relative distance between the receivers, on relative dynamics, and on the status of ionospheric activity. In this sense, the ionospheric activity plays a major role in determining the set of GPS measurements and combinations to process to improve the relative positioning accuracy.In this paper, a completely ionosphere-free approach is pursued, in which the ionospheric delays are cancelled out by combination of dual frequency GPS measurements. Several alternative combinations are investigated based on the ionosphere-free combination of pseudorange and carrier-phase observables, but also on GRoup And PHase Ionospheric Correction (GRAPHIC) and Melbourne-Wubbena combinations and thus the best combinations are selected. Based on a relative positioning scheme previously developed by the authors [17], the performance of each approach is quantified over real-world spaceborne GPS data made available by the Gravity Recovery and Climate Experiment (GRACE) mission.The paper is organized as follows. First the conventional approach is presented in which the ionospheric delays are estimated as part of the state vector through Lear’s model. Then ionosphere-free observables are derived from the DD observation models, and two different combinations are selected for being integrated in an ionosphere-free formulation of the EKF. The developed filter is finally tested on GRACE data.


## 2. Ionospheric Delays Estimated by Lear’s Model

The most common approach for handling the differential ionospheric delays in high accuracy long baseline applications is to estimate the DD ionospheric delays within the EKF modelling them using random walk processes in the state vector [14, 15]. This scheme usually integrates an EKF with an Integer Least Squares estimator based on the LAMDBA method. This approach can be not sufficiently accurate in real-time on-board implementations when accurate stochastic and dynamic models cannot be used. In such conditions, a reliable way to proceed is to model the ionospheric delays through the VTEC above the receivers by Lear’s mapping function [1, 8]. With specific reference to a formation of two satellites, this approach leads to the following state and measurement vectors:(1)x∈R8+2n,x=bb˙VTECNWLN1T;y∈R4n,y=P1P2L1L2T,where b and b˙ are the baseline and the relative velocity vectors expressed in the Earth-Centred Earth-Fixed (ECEF) reference frame and N1 stack together the cycle ambiguities on the L1 frequency for the n DD couples:(2)N1∈Rn,N1=N1AB01⋯N1AB0nT,where the pivot satellite, that is, the reference GPS satellite selected for calculating DD observables, is indicated with 0 for simplicity. NWL, instead, stacks together wide-lane cycle ambiguities. In (1) P indicates DD pseudorange observables and L represents DD carrier-phase ones. The subscripts 1 and 2 stand for L1 and L2 GPS frequencies, respectively, whereas A and B subscripts refer to the receivers onboard the satellites. Finally, VTEC is the vector including the two vertical total electron contents for the two receivers. According to Lear’s model the zero difference ionospheric delay in meters for the GPS satellite i and the receiver A can be written as [1, 8](3)IAi=40.3m3/s2f122.037sin⁡EAi+sin2⁡EAi+0.076VTECA,where f1 is L1 carrier frequency in Hz, VTECA is the (scalar) vertical total electron content for modeling the ionospheric delay of the receiver A and is expressed in number of electrons per square meter, and EAi is the elevation of the GPS satellite i with respect to the receiver A.The main advantage of Lear’s model is the capability to predict zero difference ionospheric delays, relevant to different tracked GPS satellites, as a function of a single unknown parameter (i.e., VTECA), which is a desirable property for navigation filters [1].This reference model of (1)–(3) has shown satisfactory observability features [17], since it is capable of delivering good estimates of the Integer Ambiguities (IA) even in case of only 3 DD observations. However, this VTEC-based EKF has some inherent accuracy limitations, due to its inability of rejecting deviations of the true ionosphere from Lear’s model, which appear as additional error terms in the baseline estimate [8]. In what follows, the possibility is discussed of overcoming the limitations of the VTEC-based EKF by deleting the ionospheric delays through measurement combinations.


## 3. Ionosphere-Free Observables

Dual frequency DD carrier-phase and pseudorange observables can be combined in different ways to generate ionosphere-free measurements. The complete DD observation model is thus presented before deriving the relevant ionosphere-free observables.


### 3.1. DD Observation Model

The DD measurements have the following expressions [21]:(4)P1ABjk=ρABjk+IABjk+η1ABjk,P2ABjk=ρABjk+γ-2IABjk+η2ABjk,L1ABjk=ρABjk-IABjk+λ1N1ABjk+β1ABjk,L2ABjk=ρABjk-γ-2IABjk+λ2N2ABjk+β2ABjk,where(i)
γ=f2/f1 is ratio between the L2 and L1 frequencies, and the j and k superscripts refer to the GPS satellites radiating the navigation signal;(ii)
λ1, λ2 are the wavelengths of L1 and L2 signals, respectively;(iii)
ρABjk=ρBk-ρAk-ρBj+ρAj is the DD geometric term between the two receivers;(iv)
IABjk=IBk-IAk-IBj+IAj is the DD ionospheric delay on the L1 frequency, which is denoted simply by ionospheric delay in the following;(v)
N1ABjk=N1Bk-N1Ak-N1Bj+N1Aj is the integer ambiguity on the L1 frequency (on L2 which is analogous);(vi)
η1ABjk=η1Bk-η1Ak-η1Bj+η1Aj is the P1 noise term on the L1 frequency (on L2 which is analogous);(vii)
β1ABjk=β1Bk-β1Ak-β1Bj+β1Aj is the L1 noise term on the L1 frequency (on L2 which is analogous).Each of the four observables in (4) is assumed to be independent from the other ones. However, the n DD measurements of the same kind at a certain time epoch are mutually correlated due to the presence of the pivot satellite in all the measurements. More precisely, denoting generically by X the observation type, that is, X=P1,…,L2, and by σX its standard deviation, we have(5)cov⁡XAB01⋯XAB0n=σX2D,D≔42⋯2242⋮⋱⋮2⋯⋯4.


### 3.2. Ionospheric Free Combination

The most common combination for eliminating the ionospheric delay, referred to as Ionosphere-Free (IF) [21], is concerned with combining observations of the same type on the two carrier frequencies, exploiting the frequency dependence of the first-order ionospheric delay effect. More precisely, the IF combinations are obtained as(6)XIF=f22X2-f12X1f22-f12=γ2X2-X1γ2-1.Thus, two IF observables per each DD couple jk out of the four measurements in (4) are(7)PIFABjk=11-γ2P1ABjk-γ2P2ABjk,(8)LIFABjk=11-γ2L1ABjk-γ2L2ABjk.The IF observation model can be obtained combining (4), (6), and (7) which yields(9)PIFABjk=ρABjk+ηIFABjk,LIFABjk=ρABjk+11-γ2λ1N1ABjk-γ2λ2N2ABjk+βIFABjk.Since the observation types are assumed to be independent, assuming white Gaussian measurements noises, the IF combinations are affected by Gaussian white noise with variance(10)Var⁡XIFABjk=4σX12+γ4σX221-γ22,X=P,L.Hence, the noise is increased compared to the original uncombined observations.


### 3.3. GRAPHIC Combinations

GRAPHIC combinations exploit the asymmetry of the ionospheric effect on group and phase propagation [22]. In practice, they combine pseudorange and carrier-phase measurements on each frequency, as follows:(11)G1ABjk=0.5P1ABjk+L1ABjk,G2ABjk=0.5P2ABjk+L2ABjk.From the above equations, the GRAPHIC observation model and variance read(12)G1ABjk=ρABjk+λ12N1ABjk+ϕ1ABjk,G2ABjk=ρABjk+λ22N2ABjk+ϕ2ABjk,Var⁡GiABjk=σPi2+σLi2,i=1,2,where ϕ1, ϕ2 indicate the noise terms of G1 and G2 combinations, respectively.


### 3.4. Melbourne-Wubbena Combinations

Melbourne-Wubbena (MW) combinations combine all four observable types for cancelling out the ionospheric delay [23, 24]. They build upon the definition of the wide and narrow lane (NL) wavelengths:(13)λWL=cf1-f2;λNL=cf1+f2,where c is the speed of light in vacuum. The MW combinations are obtained as(14)MW′ABjk=λWLL1ABjkλ1-L2ABjkλ2-λNLP1ABjkλ1+P2ABjkλ2and have the following observation model and variance:(15)MWABjk=λWLNWLABjk+ωABjk;NWLABjk=N1ABjk-N2ABjk,Var⁡MWABjk=4λNL2σP12λ12+σP22λ22+4λWL2σL12λ12+σL22λ22,where ω represents the noise term of the combination.


## 3.1. DD Observation Model

The DD measurements have the following expressions [21]:(4)P1ABjk=ρABjk+IABjk+η1ABjk,P2ABjk=ρABjk+γ-2IABjk+η2ABjk,L1ABjk=ρABjk-IABjk+λ1N1ABjk+β1ABjk,L2ABjk=ρABjk-γ-2IABjk+λ2N2ABjk+β2ABjk,where(i)
γ=f2/f1 is ratio between the L2 and L1 frequencies, and the j and k superscripts refer to the GPS satellites radiating the navigation signal;(ii)
λ1, λ2 are the wavelengths of L1 and L2 signals, respectively;(iii)
ρABjk=ρBk-ρAk-ρBj+ρAj is the DD geometric term between the two receivers;(iv)
IABjk=IBk-IAk-IBj+IAj is the DD ionospheric delay on the L1 frequency, which is denoted simply by ionospheric delay in the following;(v)
N1ABjk=N1Bk-N1Ak-N1Bj+N1Aj is the integer ambiguity on the L1 frequency (on L2 which is analogous);(vi)
η1ABjk=η1Bk-η1Ak-η1Bj+η1Aj is the P1 noise term on the L1 frequency (on L2 which is analogous);(vii)
β1ABjk=β1Bk-β1Ak-β1Bj+β1Aj is the L1 noise term on the L1 frequency (on L2 which is analogous).Each of the four observables in (4) is assumed to be independent from the other ones. However, the n DD measurements of the same kind at a certain time epoch are mutually correlated due to the presence of the pivot satellite in all the measurements. More precisely, denoting generically by X the observation type, that is, X=P1,…,L2, and by σX its standard deviation, we have(5)cov⁡XAB01⋯XAB0n=σX2D,D≔42⋯2242⋮⋱⋮2⋯⋯4.


## 3.2. Ionospheric Free Combination

The most common combination for eliminating the ionospheric delay, referred to as Ionosphere-Free (IF) [21], is concerned with combining observations of the same type on the two carrier frequencies, exploiting the frequency dependence of the first-order ionospheric delay effect. More precisely, the IF combinations are obtained as(6)XIF=f22X2-f12X1f22-f12=γ2X2-X1γ2-1.Thus, two IF observables per each DD couple jk out of the four measurements in (4) are(7)PIFABjk=11-γ2P1ABjk-γ2P2ABjk,(8)LIFABjk=11-γ2L1ABjk-γ2L2ABjk.The IF observation model can be obtained combining (4), (6), and (7) which yields(9)PIFABjk=ρABjk+ηIFABjk,LIFABjk=ρABjk+11-γ2λ1N1ABjk-γ2λ2N2ABjk+βIFABjk.Since the observation types are assumed to be independent, assuming white Gaussian measurements noises, the IF combinations are affected by Gaussian white noise with variance(10)Var⁡XIFABjk=4σX12+γ4σX221-γ22,X=P,L.Hence, the noise is increased compared to the original uncombined observations.


## 3.3. GRAPHIC Combinations

GRAPHIC combinations exploit the asymmetry of the ionospheric effect on group and phase propagation [22]. In practice, they combine pseudorange and carrier-phase measurements on each frequency, as follows:(11)G1ABjk=0.5P1ABjk+L1ABjk,G2ABjk=0.5P2ABjk+L2ABjk.From the above equations, the GRAPHIC observation model and variance read(12)G1ABjk=ρABjk+λ12N1ABjk+ϕ1ABjk,G2ABjk=ρABjk+λ22N2ABjk+ϕ2ABjk,Var⁡GiABjk=σPi2+σLi2,i=1,2,where ϕ1, ϕ2 indicate the noise terms of G1 and G2 combinations, respectively.


## 3.4. Melbourne-Wubbena Combinations

Melbourne-Wubbena (MW) combinations combine all four observable types for cancelling out the ionospheric delay [23, 24]. They build upon the definition of the wide and narrow lane (NL) wavelengths:(13)λWL=cf1-f2;λNL=cf1+f2,where c is the speed of light in vacuum. The MW combinations are obtained as(14)MW′ABjk=λWLL1ABjkλ1-L2ABjkλ2-λNLP1ABjkλ1+P2ABjkλ2and have the following observation model and variance:(15)MWABjk=λWLNWLABjk+ωABjk;NWLABjk=N1ABjk-N2ABjk,Var⁡MWABjk=4λNL2σP12λ12+σP22λ22+4λWL2σL12λ12+σL22λ22,where ω represents the noise term of the combination.


## 4. Ionosphere-Free Relative Positioning

This section is concerned with establishing which of the measurement combinations presented in Section3 is suitable for computing the relative position of the two receivers in long baseline applications. The measurement combinations available with no ionospheric effects are the nPIF combinations, the nLIF combinations, the n GRAPHIC combinations on L1 and L2, and the n MW combinations, for a total of 5n measurements. However, the 5n measurements are not linearly independent. In particular, each group of 5 observables per each of the n DD couples can be seen as a linear transformation of the 4 uncombined measurements of (4). More precisely, the cancellation of the ionosphere from (4) can be seen as a linear projection of the 4-dimensional measurement vector P1,P2,L1,L2T onto a 3-dimensional hyperplane. This implies that no more than three linear combinations of the four dimensional vector exist being linearly independent, and thus suitable for use as measurement vector in an EKF. In addition, whatever set of three linearly independent vectors lying on the hyperplane can be used as a basis for describing all the vectors belonging to the hyperplane. Even though any of such bases would yield theoretically equivalent positioning solutions, significant differences exist in practical applications due to unmodeled systematic errors. In order to find such “best” basis, its performance in accurately estimating both the geometric term and the cycle ambiguity is analyzed. A natural indicator to evaluate the accuracy achieved by a specific measurement combination is the standard deviation (STD) of the noise affecting the measurements.Hence, a specific analysis has been performed. Starting from different candidate noise levels for the uncombined measurements, that is,P1, P2, L1, and L2 (see the first three columns in Table 1), reflecting typical performance of spaceborne GPS receivers, the expected STD for the estimation of the geometric term ρjk is derived as listed in Table 1. The rationale of this analysis is that different combinations of GPS measurements may modify, that is either amply or reduce, the uncertainty of the geometric term. Combinations able to reduce this uncertainty should be preferred for implementation in an ionosphere-free EKF. The same analysis is then repeated for the integer ambiguity on L1. Again, the smaller the uncertainty on this ambiguity, the larger the probability to correctly fix it. The results are shown in Table 2. In this case, a scale factor is introduced to represent the STD of combined measurements as a fraction of L1 cycles. In addition to this theoretical analysis, the dispersion of actual GRACE data is computed and the relevant STD is calculated (see the last row in Tables 1 and 2) for gaining further insight into the true-world accuracy of the various combinations.Table 1
Expected and estimated standard deviations of the various combinations forρjk estimation (NA: not available).
P
1
P
2
L
1, L2
L
IF
G
1
G
2
P
IF
MW
Formula
σ
P
1
2
σ
P
2
2
σ
L
1
2
=
σ
L
2
2
2
σ
L
1
2
+
γ
4
σ
L
2
2
1
-
γ
2
2
σ
P
1
2
+
σ
L
1
2
σ
P
2
2
+
σ
L
2
2
2
σ
P
1
2
+
γ
4
σ
P
2
2
1
-
γ
2
2
NA
Unit
m
M
m
m
m
m
m
m
Example1
0.150
0.200
0.010
0.060
0.150
0.200
0.983
NA
Example2
0.100
0.150
0.005
0.030
0.100
0.150
0.689
NA
Example3
0.100
0.150
0.001
0.006
0.100
0.150
0.689
NA
GRACE data
0.005
0.082
0.149
0.708
NA
(Jan. 09)Table 2
Expected and estimated standard deviations of the various combinations forN1jk estimation (NA: not available).
P
1
P
2
L
1, L2
L
IF
G
1
G
2, PIF
MW
Formula
σ
P
1
2
σ
P
2
2
σ
L
1
2
=
σ
L
2
2
2
σ
L
1
2
+
γ
4
σ
L
2
2
λ
1
1
-
γ
2
2
λ
1
σ
P
1
2
+
σ
L
1
2
NA
2
λ
NL
2
λ
WL
σ
P
1
2
λ
1
2
+
σ
P
2
2
λ
2
2
+
λ
WL
σ
L
1
2
λ
1
2
+
σ
L
2
2
λ
2
2
Unit
m
m
m
Cycles
Cycles
Cycles
Cycles
Example1
0.150
0.200
0.010
0.123
1.580
NA
0.312
Example2
0.100
0.150
0.005
0.061
1.052
NA
0.211
Example3
0.100
0.150
0.001
0.012
1.051
NA
0.201
GRACE data
0.010
0.862
NA
0.194
(Jan. 09)According to the results presented in Tables1 and 2, the set of measurements including LIF, G1, and MW combinations is the one allowing a higher potential accuracy on both the DD geometric term and the cycle ambiguities. This set of measurements is called set “3n” in what follows. A potential limitation of set 3n is that cycle ambiguities are present in all the three observation models. Hence, when the number of fixed ambiguities is smaller than three, the observation model provides no useful information to the EKF, which is thus forced to rely only on its propagation model. In this particular situation, the EKF is expected to perform poorly. To cope with this problem an additional set is considered, called set “4n” and made of LIF, G1, and MW combinations, plus the uncombined P1 measurement. The addition of P1 observations allows the filter to rely on relative positioning information even when less than three cycle ambiguities are fixed.Because the various combinations have different observation models selecting different sets of measurements affects the observability of the EKF state vector differently. The observability of the system at hand is extremely involved to be analyzed. Indeed, the dynamical system is nonlinear and evolves on a time-varying trajectory. Its linearization is thus time-varying and depends on the state estimate. The observation model is nonlinear, due to the DD geometry terms. Its linearization depends on the state trajectory, and thus it is time-varying as well. Both the dynamical and the observation models depend on an exogenous vector parameter (e.g., the absolute position of one of the receivers) [17]. Finally, both the EKF state vector and the measurements depend on the number of DD couples n, which is time-varying as well. In such a condition, instead of developing a complicated observability analysis a different approach is pursued, which is based on numerically evaluating the EKF performances using both set 3n and set 4n as shown in the next section.


### 4.1. No-Iono EKF Dynamic Model

The no-iono EKF state vector is defined as(16)x∈R6+2n,x=bb˙NWLN1T.Concerning the process model, the baseline dynamics is described by the nonlinear relative Keplerian model with J2 effects included. This model is equal to the one used in the VTEC-based EKF [17] and is denoted at the time epoch tk as follows:(17)b¨=fb,b˙,rA,bb˙k+1=Φb,kbb˙k+wb;EwbwbT=Qb,where f is the nonlinear baseline propagation function, rA denotes the receiver A position vector in ECEF, Φ is the time-discrete baseline propagation matrix, wb is the baseline process noise vector, and E indicates the mathematical expectation.Cycle ambiguities are instead assumed to be constant. However, a small process noisewN might be injected, letting them be modeled as a random constant plus random walk. The model is equal to the one used in the VTEC-based EKF [17] and is denoted by (1n is the n-by-n identity matrix)(18)N˙WLN˙1=0,NWLN1k+1=1nNWLN1k+wN;EwNwNT=QN=σN21n.


### 4.2. No-Iono EKF Observation Model

The two different observation models share some common structure, which is discussed first. Then, details on each of the two models will be given. The no-iono EKF measurement vectory is made of a series of measurement combinations, depending on the chosen set. Each measurement combination stacks together the combinations for each of the n DD couples:(19)X∈Rn,X=XAB01⋯XAB0nT.The covariance matrix of y has a nonsparse structure, due to the presence of the various uncombined measurements in more combinations, and because of the full covariance matrix of each of the uncombined measurements X. In order to have the desired covariance matrix, cov(y), it is possible to introduce the vector z that stacks the uncombined DD observations P1,…,L2:(20)z=P1AB01⋯P1AB0nP2AB01⋯P2AB0nL1AB01⋯L1AB0nL2AB01⋯L2AB0nT,covz=σP12D0σP22DσL2D0σL2D.The chosen combined measurements set can be expressed as a linear function of the uncombined measurement vector z, by means of a matrix T depending on the chosen set. The matrix cov(y) can be obtained by cov(z) applying the linear transformation:(21)y=Tz;R=cov⁡y=Tcov⁡zTT.The observation model relating the EKF state to the observables is nonlinear due to the geometric term. The geometric term is independent from the chosen measurement set and thus is discussed herein. The n DD geometric terms are related to the baseline by the following equations (Ri denotes the ith GPS satellite position vector in ECEF):(22)ρABjk=Rk-rA+b-Rk-rA-Rj-rA+b+Rj-rA.The nonlinear observation model and the linear delta observation model, needed for the computation of the Kalman gain in the EKF, are obtained by linearization of the propagated state vector, and it can be expressed generically as(23)y∈R3n,y=hx,Δy∈R3n,Δy=HΔx.Thus, at each time epoch, from the observation equations, the measurement vector y is computed combining uncombined measurements. The choice of the measurement set to consider thus changes only the following: (1)T matrix for obtaining the EKF combined measurement vector y from the uncombined one z, (2)h(x) nonlinear observation function, and (3)H matrix of the linearized observation model. These items are discussed in the following subsections for each measurement set.


#### 4.2.1. Set3n

TheT matrix, the h(x) function, and the H matrix take the following form:(24)
T
=
0
n
0
n
1
1
-
γ
2
1
n
-
γ
2
1
-
γ
2
1
n
1
2
1
n
0
n
1
2
1
n
0
n
-
λ
NL
λ
1
1
n
-
λ
NL
λ
2
1
n
λ
WL
λ
1
1
n
-
λ
WL
λ
2
1
n
,
h
x
=
ρ
A
B
j
b
ρ
A
B
j
b
0
n
×
1
+
λ
1
γ
1
-
γ
2
1
n
λ
1
1
+
γ
1
n
0
n
λ
1
2
1
n
λ
WL
1
n
0
n
N
WL
N
1
,
H
=
∇
b
ρ
A
B
j
0
3
λ
1
γ
1
-
γ
2
1
n
λ
1
1
+
γ
1
n
∇
b
ρ
A
B
j
0
3
0
n
λ
1
2
1
n
0
3
0
3
λ
WL
1
n
0
n
,where ρABj vector stacks together n DD geometric terms, ∇b is the gradient operator with respect to the baseline vector, and 0n is the n-by-n null matrix.


#### 4.2.2. Set4n

TheT matrix, the h(x) function, and the H matrix take the following form:(25)
T
=
1
n
0
n
0
n
0
n
T
3
n
,
h
x
=
ρ
A
B
j
b
h
3
n
x
,
H
=
∇
b
ρ
A
B
j
0
3
0
n
0
n
H
3
n
,where the subscript “3n” is used to indicate matrix and functions relevant to set 3n. It is worth noting that, since the ionospheric delay is cancelled and not included in the state vector, the implemented observation model for P1 measurement is (26)P1ABjk=ρABjk+η1ABjk.Comparing (4) to (26), it is evident that set 4n implicitly assumes the DD ionospheric delay modelled as an additional Gaussian white noise affecting the P1 measurements. In other words, considering that the DD ionospheric delay is only a small fraction of the total ionospheric delays, the DD P1 observable noise is modelled as if the ionospheric delay was not time correlated but included within the white noise error. This approach resembles that one used in [7, 14] where ionospheric delays are modelled as random walk processes characterized by large process noises and correlation times in the order of the filter time step.


## 4.1. No-Iono EKF Dynamic Model

The no-iono EKF state vector is defined as(16)x∈R6+2n,x=bb˙NWLN1T.Concerning the process model, the baseline dynamics is described by the nonlinear relative Keplerian model with J2 effects included. This model is equal to the one used in the VTEC-based EKF [17] and is denoted at the time epoch tk as follows:(17)b¨=fb,b˙,rA,bb˙k+1=Φb,kbb˙k+wb;EwbwbT=Qb,where f is the nonlinear baseline propagation function, rA denotes the receiver A position vector in ECEF, Φ is the time-discrete baseline propagation matrix, wb is the baseline process noise vector, and E indicates the mathematical expectation.Cycle ambiguities are instead assumed to be constant. However, a small process noisewN might be injected, letting them be modeled as a random constant plus random walk. The model is equal to the one used in the VTEC-based EKF [17] and is denoted by (1n is the n-by-n identity matrix)(18)N˙WLN˙1=0,NWLN1k+1=1nNWLN1k+wN;EwNwNT=QN=σN21n.


## 4.2. No-Iono EKF Observation Model

The two different observation models share some common structure, which is discussed first. Then, details on each of the two models will be given. The no-iono EKF measurement vectory is made of a series of measurement combinations, depending on the chosen set. Each measurement combination stacks together the combinations for each of the n DD couples:(19)X∈Rn,X=XAB01⋯XAB0nT.The covariance matrix of y has a nonsparse structure, due to the presence of the various uncombined measurements in more combinations, and because of the full covariance matrix of each of the uncombined measurements X. In order to have the desired covariance matrix, cov(y), it is possible to introduce the vector z that stacks the uncombined DD observations P1,…,L2:(20)z=P1AB01⋯P1AB0nP2AB01⋯P2AB0nL1AB01⋯L1AB0nL2AB01⋯L2AB0nT,covz=σP12D0σP22DσL2D0σL2D.The chosen combined measurements set can be expressed as a linear function of the uncombined measurement vector z, by means of a matrix T depending on the chosen set. The matrix cov(y) can be obtained by cov(z) applying the linear transformation:(21)y=Tz;R=cov⁡y=Tcov⁡zTT.The observation model relating the EKF state to the observables is nonlinear due to the geometric term. The geometric term is independent from the chosen measurement set and thus is discussed herein. The n DD geometric terms are related to the baseline by the following equations (Ri denotes the ith GPS satellite position vector in ECEF):(22)ρABjk=Rk-rA+b-Rk-rA-Rj-rA+b+Rj-rA.The nonlinear observation model and the linear delta observation model, needed for the computation of the Kalman gain in the EKF, are obtained by linearization of the propagated state vector, and it can be expressed generically as(23)y∈R3n,y=hx,Δy∈R3n,Δy=HΔx.Thus, at each time epoch, from the observation equations, the measurement vector y is computed combining uncombined measurements. The choice of the measurement set to consider thus changes only the following: (1)T matrix for obtaining the EKF combined measurement vector y from the uncombined one z, (2)h(x) nonlinear observation function, and (3)H matrix of the linearized observation model. These items are discussed in the following subsections for each measurement set.


### 4.2.1. Set3n

TheT matrix, the h(x) function, and the H matrix take the following form:(24)
T
=
0
n
0
n
1
1
-
γ
2
1
n
-
γ
2
1
-
γ
2
1
n
1
2
1
n
0
n
1
2
1
n
0
n
-
λ
NL
λ
1
1
n
-
λ
NL
λ
2
1
n
λ
WL
λ
1
1
n
-
λ
WL
λ
2
1
n
,
h
x
=
ρ
A
B
j
b
ρ
A
B
j
b
0
n
×
1
+
λ
1
γ
1
-
γ
2
1
n
λ
1
1
+
γ
1
n
0
n
λ
1
2
1
n
λ
WL
1
n
0
n
N
WL
N
1
,
H
=
∇
b
ρ
A
B
j
0
3
λ
1
γ
1
-
γ
2
1
n
λ
1
1
+
γ
1
n
∇
b
ρ
A
B
j
0
3
0
n
λ
1
2
1
n
0
3
0
3
λ
WL
1
n
0
n
,where ρABj vector stacks together n DD geometric terms, ∇b is the gradient operator with respect to the baseline vector, and 0n is the n-by-n null matrix.


### 4.2.2. Set4n

TheT matrix, the h(x) function, and the H matrix take the following form:(25)
T
=
1
n
0
n
0
n
0
n
T
3
n
,
h
x
=
ρ
A
B
j
b
h
3
n
x
,
H
=
∇
b
ρ
A
B
j
0
3
0
n
0
n
H
3
n
,where the subscript “3n” is used to indicate matrix and functions relevant to set 3n. It is worth noting that, since the ionospheric delay is cancelled and not included in the state vector, the implemented observation model for P1 measurement is (26)P1ABjk=ρABjk+η1ABjk.Comparing (4) to (26), it is evident that set 4n implicitly assumes the DD ionospheric delay modelled as an additional Gaussian white noise affecting the P1 measurements. In other words, considering that the DD ionospheric delay is only a small fraction of the total ionospheric delays, the DD P1 observable noise is modelled as if the ionospheric delay was not time correlated but included within the white noise error. This approach resembles that one used in [7, 14] where ionospheric delays are modelled as random walk processes characterized by large process noises and correlation times in the order of the filter time step.


## 4.2.1. Set3n

TheT matrix, the h(x) function, and the H matrix take the following form:(24)
T
=
0
n
0
n
1
1
-
γ
2
1
n
-
γ
2
1
-
γ
2
1
n
1
2
1
n
0
n
1
2
1
n
0
n
-
λ
NL
λ
1
1
n
-
λ
NL
λ
2
1
n
λ
WL
λ
1
1
n
-
λ
WL
λ
2
1
n
,
h
x
=
ρ
A
B
j
b
ρ
A
B
j
b
0
n
×
1
+
λ
1
γ
1
-
γ
2
1
n
λ
1
1
+
γ
1
n
0
n
λ
1
2
1
n
λ
WL
1
n
0
n
N
WL
N
1
,
H
=
∇
b
ρ
A
B
j
0
3
λ
1
γ
1
-
γ
2
1
n
λ
1
1
+
γ
1
n
∇
b
ρ
A
B
j
0
3
0
n
λ
1
2
1
n
0
3
0
3
λ
WL
1
n
0
n
,where ρABj vector stacks together n DD geometric terms, ∇b is the gradient operator with respect to the baseline vector, and 0n is the n-by-n null matrix.


## 4.2.2. Set4n

TheT matrix, the h(x) function, and the H matrix take the following form:(25)
T
=
1
n
0
n
0
n
0
n
T
3
n
,
h
x
=
ρ
A
B
j
b
h
3
n
x
,
H
=
∇
b
ρ
A
B
j
0
3
0
n
0
n
H
3
n
,where the subscript “3n” is used to indicate matrix and functions relevant to set 3n. It is worth noting that, since the ionospheric delay is cancelled and not included in the state vector, the implemented observation model for P1 measurement is (26)P1ABjk=ρABjk+η1ABjk.Comparing (4) to (26), it is evident that set 4n implicitly assumes the DD ionospheric delay modelled as an additional Gaussian white noise affecting the P1 measurements. In other words, considering that the DD ionospheric delay is only a small fraction of the total ionospheric delays, the DD P1 observable noise is modelled as if the ionospheric delay was not time correlated but included within the white noise error. This approach resembles that one used in [7, 14] where ionospheric delays are modelled as random walk processes characterized by large process noises and correlation times in the order of the filter time step.


## 5. Results

The reference EKF and the no-iono solutions have been run on GRACE data. Comparisons among the solutions are performed looking at the baseline magnitude estimation error and at the capability to obtain correct estimates of the integer ambiguities. Both of them are estimated thanks to the high accuracy reference solutions, which can be computed processing GRACE data [8]. As expected, application to mild and moderate ionospheric conditions, not shown herein for brevity, suggest that the filter with Lear’s model is capable to achieve a baseline estimation error which is the half of the no-iono approach (both set 3n and set 4n). This behaviour changes under strong ionospheric activity. In order to show this, DOY286 (October 13th) 2011 has been selected for the analysis, in which the two GRACE spacecrafts are separated of more than 200 km, and ionospheric activity is extremely intense. Figures 1–3 report the baseline magnitude estimation error and the IA estimation performance for the reference solution and the no-iono solutions, respectively. Specifically, only L1 estimation performance is shown because wide lane fail rates are significantly lower than those on L1.Figure 1
Baseline magnitude estimation error and L1 ambiguities estimation performance for the reference solution.Figure 2
Baseline magnitude estimation error and L1 ambiguities estimation performance for set3n.Figure 3
Baseline magnitude estimation error and L1 ambiguities estimation performance for set4n.Figure1 shows that, under strong ionospheric activity, Lear’s model is not able to support the integer estimation, so a significant number of ambiguities is wrong and baseline magnitude estimation errors larger than 1 m occur. IA estimation performance increases notably when set 3n is considered (see Figure 2), and the baseline estimation error reduces to decimetre scale in this case. Finally, with reference to set 4n, it is evident that the addition of the P1 measurements is not a successful strategy when the ionospheric activity is intense. Indeed, DD ionospheric delays dominate over P1 measurement noise and so they cannot be modelled as additional Gaussian white noise.


## 6. Concluding Remarks

This paper has compared various strategies for compensating ionospheric delays in long baseline relative positioning applications for LEO spacecraft using dual-frequency CDGPS techniques. Based on a relative positioning technique previously developed by the authors, problem settings that do not require ionospheric delays estimation have been investigated. The ionosphere-free alternatives have been also compared with a reference solution making use of a standard ionospheric delay model for spaceborne applications, originally developed by Lear. Results suggest that, among the possible alternatives and under strong ionospheric activity, best performance is achieved using the observation model made of ionosphere-free combinations of carrier-phase measurements, GRAPHIC combinations on the L1 frequency, and Melbourne-Wubbena combinations. Future work will deal with the extension of these results to a more complete data set, including different levels of the ionospheric activity. The ultimate goal is the introduction of quantitative tools to guide the selection between model-based and ionosphere-free approaches depending on the ionospheric activity and the available ionospheric models.



---

*Source: 570382-2015-07-16.xml*