GNSS Kinematic Position and Velocity
Determination for Airborne Gravimetry
Zur GNSS-basierten Bestimmung von Position
und Geschwindigkeit in der Fluggravimetrie
vorgelegt von
M. Eng.
Kaifei He
geb. in Xianyang, Shaanxi, China
von der Fakultät VI – Planen Bauen Umwelt
der Technischen Universität Berlin
zur Erlangung des akademischen Grades
Doktor der Ingenieurwissenschaften
– Dr.-Ing. –
genehmigte Dissertation
Promotionsausschuss:
Vorsitzender: Prof. Dr.-Ing. Dr. h.c. Harald Schuh
Gutachter: Prof. Dr.-Ing. Frank Neitzel
Gutachter: Prof. Dr.-Ing. Frank Flechtner
Gutachter: Prof. Dr. Phil. nat. Markus Rothacher (ETH Zürich)
Tag der wissenschaftlichen Aussprache: 03. Dezember 2014
Berlin 2015
- i -
Abstract
The Global Navigation Satellite System (GNSS) plays a significant role in the fields of air-
borne gravimetry. The objective of this thesis is to develop reliable GNSS algorithms and
software for kinematic highly precise GNSS data analysis in airborne gravimetry. Based on
the requirements for practical applications in airborne gravimetry and shipborne gravimetry
projects, the core research and the contributions of this thesis are summarized as follows:
Estimation Algorithm: Based on the accuracy requirements for GNSS precise positioning
in airborne gravimetry, the estimation algorithms of least squares including the elimination
of nuisance parameters as well as a two-way Kalman filter are applied to the kinematic
GNSS data post-processing. The goal of these adjustment methods is to calculate non-epoch
parameters (such as system error estimates or carrier phase ambiguity parameters) using all
data in the first step, followed by the calculation of epoch parameters (such as position and
velocity parameters of the kinematic platform) at every epoch. These methods are highly
efficient when dealing with massive amounts of data, and give the highly precise results for
the GNSS data analyzed.
Accuracy Evaluation and Reliability Analysis: The accuracy evaluation and reliability
analysis of the results from precise kinematic GNSS positioning is studied. A special accura-
cy evaluation method in GNSS kinematic positioning is proposed, where the known
distances among multiple antennas of GNSS receivers are taken as an accuracy evaluation
index. The effect of the GNSS receiver clock error in the accuracy evaluation for GNSS kin-
ematic positioning results of a high-speed motion platform is studied and a solution is
proposed.
Kinematic Positioning Based on Multiple Reference Stations Algorithms: In order to
overcome the problem of decreasing accuracy in GNSS relative kinematic positioning for
long baselines, a new relative kinematic positioning method based on a priori constraints for
multiple reference stations is proposed. This algorithm increases the accuracy and reliability
of kinematic positioning results for large regions resp. long baselines.
GNSS Precise Positioning Based on Robust Estimation: In order to solve the problem of
outliers occurring in positioning results which are caused by the presence of gross errors in
the GNSS observations, a robust estimation algorithm is applied to eliminate the effects of
gross errors in the results of GNSS kinematic precise positioning.
- ii -
Kinematic Positioning Based on Multiple Kinematic Stations: In airborne gravimetry,
multiple antennas of GNSS receivers are usually mounted on the kinematic platform. Firstly,
a GNSS kinematic positioning method based on multiple kinematic stations is proposed. Us-
ing the known constant distances among the multiple GNSS antennas, a kinematic
positioning method based on a priori distance constraints is proposed to improve the reliabil-
ity of the system. Secondly, such an approach is also used for the estimation of a common
atmospheric wet delay parameter among the multiple GNSS antennas mounted on the plat-
form. This method does not only reduce the amount of estimated parameters, but also
decreases the correlation among the atmospheric parameters.
Kinematic Positioning Based on GNSS Integration (GPS and GLONASS): To improve
the reliability and accuracy of kinematic positioning, a kinematic positioning method using
multiple GNSS systems integration is addressed. Furthermore, a GNSS integration algorithm
based on Helmert’s variance components estimation is proposed to adjust the weights in a
reasonable way. This improves the results when combining data of the different GNSS sys-
tems.
Velocity Determination Using GNSS Doppler Data: Airborne gravimetry requires instan-
taneous velocity results, thus raw Doppler observations are used to determine the kinematic
instantaneous velocity in high-dynamic environments. Furthermore, carrier phase derived
Doppler observations are used to obtain precise velocity estimates in low-dynamic environ-
ments. Then a method of Doppler velocity determination based on GNSS integration with
Helmert’s variance components estimation and robust estimation is studied.
Software Development and Application: In order to fulfill the actual requirements of air-
borne as well as shipborne gravimetry on GNSS precise positioning, a software system
(HALO_GNSS) for precise kinematic GNSS trajectory and velocity determination for kine-
matic platforms has been developed. In this software, the algorithms as proposed in this
thesis were adopted and applied. In order to evaluate the effectiveness of the proposed algo-
rithm and the HALO_GNSS software, this software is applied in airborne as well as
shipborne gravimetry projects of GFZ Potsdam. All results are compared and examined, and
it is shown that the applied approaches can effectively improve the reliability and accuracy
of the kinematic position and velocity determination. It allows the kinematic positioning with
an accuracy of 1-2 cm and the velocity determination with an accuracy of approximately 1
cm/s using raw and approximately 1 mm/s using carrier phase derived Doppler observations.
Keywords: GNSS, Kinematic Precise Positioning, Airborne Gravimetry, Least Squares,
Nuisance Parameter Elimination, Two-way Kalman filter, Receiver Clock Error, Multiple
Reference Stations, Robust Estimation, Multiple Kinematic Stations, A Priori Distance Con-
- iii -
straint, Common Troposphere Parameterization, GNSS Integration, Helmert’s Variance
Components Estimation, Doppler Velocity Determination, GEOHALO, HALO, HA-
LO_GNSS.
- v -
Zusammenfassung
Das weltumspannende Satelliten-Navigationssystem GNSS spielt eine wichtige Rolle für die
Fluggravimetrie. Gegenstand dieser Arbeit ist die Entwicklung zuverlässiger GNSS-
Algorithmen und Software für die hochgenaue GNSS-Datenanalyse in der Fluggravimetrie.
Ausgehend von den Anforderungen für praktische Anwendungen der Fluggravimetrie lassen
sich die Beiträge und Schwerpunkte dieser Dissertation wie folgt zusammenfassen:
Ausgleichs- bzw. Schätzungs-Algorithmen: Ausgehend von den Genauigkeitsanforderun-
gen an die GNSS-basierte Positionsbestimmung in der Fluggravimetrie werden in einer
kinematischen GNSS-Daten-Auswertung eine Schätzung nach kleinsten Quadraten
einschließlich der Eliminierung von Störparametern sowie ein Zwei-Wege-Kalman-Filter
angewendet. Das Ziel der beiden Ausgleichsverfahren ist es, an jedem Messzeitpunkt
zunächst globale Parameter (wie System-Fehler und Trägerwellen-Ambiguities) und
anschließend lokale Parameter (wie Position und Geschwindigkeit der bewegten
Messplattform) zu bestimmen. Die angewandten Methoden sind sehr effizient und ergeben
hochpräzise Resultate für die GNSS-Datenanalyse.
Analyse von Genauigkeit und Zuverlässigkeit: Die Genauigkeit und Zuverlässigkeit der
Resultate der präzisen kinematischen GNSS-Positionsbestimmung werden untersucht. Dabei
wird eine besondere Methode zur Bewertung der Genauigkeit der kinematischen GNSS-
Positionsbestimmung vorgeschlagen, wo bekannte Entfernungen zwischen mehreren GNSS-
Antennen als Genauigkeits-Maßstab genommen werden. Weiterhin wird der Einfluss der
Uhrenfehler der GNSS-Empfänger auf die Genauigkeit der kinematischen
Positionsbestimmung für die Hochgeschwindigkeits-Plattform untersucht. Für dabei
auftretende Probleme wird eine Lösung vorgeschlagen.
Algorithmen der kinematischen Positionsbestimmung die auf mehreren Referenzstationen
beruhen: Um das Problem der im Falle langer Basislinien abnehmenden Genauigkeit in der
relativen kinematischen GNSS-Positionsbestimmung zu bewältigen, wird ein neuer
Algorithmus vorgeschlagen. Er beruht auf der apriori Einführung von Exzentrizitäts-
Bedingungen für mehrere Referenzstationen. Dieser Algorithmus erhöht die Genauigkeit und
Zuverlässigkeit der Ergebnise in der kinematischen Positionsbestimmung für große
Regionen resp. lange Basislinien.
Präzise GNSS-Positionsbestimmung, beruhend auf robuster Schätzung: Das Vorhanden-
sein von groben Fehlern in den GNSS-Beobachtungen verursacht das Auftreten von
- vi -
Ausreißern in den Ergebnissen der Positionsbestimmung. Um dieses Problem zu überwinden,
wird ein robuster Ausgleichungs-Algorithmus angewendet, der die Auswirkungen von
groben Fehlern in den Ergebnissen der kinematischen GNSS-Positionsbestimmung beseitigt.
Kinematische Positionierung auf der Basis mehrerer bewegter Stationen: In der
Fluggravimetrie werden in der Regel mehrere GNSS-Antennen auf einer bewegten Plattform
installiert. In diesem Zusammenhang wird deshalb erstens ein kinematisches GNSS-
Positionsbestimmungsverfahren vorgeschlagen, das auf mehreren gleichzeitig bewegten
GNSS-Stationen basiert. Aus den bekannten, konstanten Distanzen zwischen den GNSS-
Antennen werden dabei apriori Exzentrizitäts-Bedingungen abgeleitet und in die
Positionsschätzung eingeführt. Dies verbessert die Zuverlässigkeit des Messsystems.
Zweitens wird solch ein Ansatz auch zur Bestimmung eines gemeinsamen
Refraktionsparameters aller GNSS-Antennen der Plattform für den feuchten Teil der
Atmosphäre verwendet. Dieses Verfahren reduziert nicht nur die Menge der geschätzten
Parameter, sondern verringert auch die Korrelation zwischen den atmosphärischen
Parametern.
Kinematische Positionierung basierend auf der Kombination verschiedener GNSS-
Systeme: Um die Zuverlässigkeit und Genauigkeit der kinematischen Positionsbestimmung
zu verbessern, werden die Signale mehrerer GNSS-Systeme (d.h. GPS und GLONASS)
gemeinsam registriert und ausgewertet (sog. GNSS-Integration). Zur Optimierung des
relativen Gewichts zwischen den Daten der verschiedenen GNSS-Systeme wird die
Helmertsche Varianz- Komponenten-Schätzung angewandt. Der auf dieser Basis entwickelte
Kombinationsalgorithmus ermöglicht die Verbesserung der Beiträge von mehreren GNSS-
Systemen.
Geschwindigkeitsbestimmung mit GNSS-Doppler-Daten: Die Auswertung der Schwere-
Messdaten in der Fluggravimetrie verlangt die hochgenaue Bestimmung des
Geschwindigkeitsvektors der bewegten Plattform. Deshalb werden rohe GNSS-Doppler-
Beobachtungen verwendet, um die Geschwindigkeit der bewegten Plattform im Falle hoch-
dynamischer Flugbedingungen kinematisch zu bestimmen. Darüberhinaus werden aus der
Trägerphase abgeleitete Doppler-Beobachtungen verwendet, um präzise
Geschwindigkeitsschätzungen im Falle weniger dynamischer Flugbedingungen zu erhalten.
Die Kombination verschiedener GNSS-Systeme wird auch bei der Doppler-
Geschwindigkeitsbestimmung angewandt. Hierzu wird die Anwendung der Helmertschen
Varianzkomponenten-Schätzung und einer robuste Schätzung untersucht.
Software Entwicklung und Anwendung: Um die aktuellen Anforderungen der GNSS-
basierten Positionsbestimmung in der Flug- sowie Schiffsgravimetrie zu erfüllen, wurde ein
- vii -
Software-System (HALO_GNSS) für die präzise kinematische GNSS-Flugbahn- und
Geschwindigkeitsberechnung kinematischer Plattformen entwickelt. Die in dieser Arbeit
vorgeschlagenen Algorithmen wurden in diese Software integriert. Um die Effizienz der
vorgeschlagenen Algorithmen und der HALO_GNSS Software zu prüfen, wurde diese
Software sowohl in Flug- aus auch in Schiffsgravimetrie-Projekten des GFZ Potsdam
angewandt. Alle Ergebnisse werden verglichen und geprüft und es wird gezeigt, dass die
angewandten Methoden die Zuverlässigkeit und Genauigkeit der kinematischen Positions-
und Geschwindigkeitsbestimmung effektiv verbessern. Die Verwendung der Software
HALO_GNSS ermöglicht kinematische Positionsbestimmung mit einer Genauigkeit von 1-2
cm sowie Geschwindigkeitsbestimmung mit einer Genauigkeit von ca. 1 cm/s mit Roh- und
etwa 1 mm/s mit aus der Trägerphase abgeleiteten Doppler-Beobachtungen.
Stichworte: GNSS, Kinematische präzise Positionsbestimmung, Fluggravimetrie,
Ausgleichung nach kleinsten Quadraten, Paremeter-Eliminierung, Zwei-Wege-Kalman-
Filter, GNSS-Uhrfehler, Mehrfache Referenzstationen, Robuste Schätzung, Mehrfache
bewegte GNSS-Stationen, A Priori Exzentrizitäts-Bedingungen, Gemeinsame Troposphären-
Parametrierung, GNSS-Integration, Helmertsche Varianzkomponenten-Schätzung, Doppler-
basierte Geschwindigkeits-Bestimmung, GEOHALO, HALO, HALO_GNSS.
- ix -
Table of Contents
Abstract ..................................................................................................................................... i
Zusammenfassung ................................................................................................................... v
Table of Contents .................................................................................................................... ix
List of Figures ....................................................................................................................... xiii
List of Tables ......................................................................................................................... xix
List of Abbreviations ............................................................................................................ xxi
1 Introduction ........................................................................................................................ 1
1.1 Motivation .................................................................................................................... 1
1.2 Background .................................................................................................................. 2
1.2.1 The Role of GNSS Application in Airborne Gravimetry .................................... 2
1.2.2 The New Airborne Gravimetric Project .............................................................. 4
1.2.3 The Development Status of GNSS ..................................................................... 5
1.3 Related Research of GNSS Data Processing .............................................................. 12
1.4 Challenges and Research Objectives .......................................................................... 14
1.5 Contributions of this Research ................................................................................... 16
1.6 Overview of Dissertation ........................................................................................... 17
2 GNSS Observations and Error Sources .......................................................................... 19
2.1 Introduction ................................................................................................................ 19
2.2 GNSS Observations .................................................................................................... 19
2.2.1 Pseudorange ...................................................................................................... 19
2.2.2 Carrier Phase .................................................................................................... 20
2.2.3 Doppler ............................................................................................................. 21
2.3 Linear Combinations of GNSS Observations ............................................................ 21
2.3.1 Ionosphere-free Linear Combination (LC) ....................................................... 22
2.3.2 Widelane Combination (WL) ........................................................................... 22
2.3.3 Extra-widelane Linear Combination (EX_WL) ............................................... 23
2.3.4 MW Widelane Linear Combination (MW_WL) .............................................. 23
2.4 Error Sources .............................................................................................................. 24
2.4.1 Satellite Orbit and Clock Errors ....................................................................... 24
2.4.2 Antenna Phase Centre Offset and Variation ..................................................... 26
2.4.3 Satellite and Receiver Hardware Biases ........................................................... 27
2.4.4 Tropospheric Effects ......................................................................................... 27
2.4.5 Ionospheric Effects ........................................................................................... 29
- x -
2.4.6 Site Displacement Effects ................................................................................. 30
2.4.7 Ambiguities and Cycle Slips ............................................................................ 32
2.5 Conclusions ................................................................................................................ 32
3 Algorithms Developing and Quality Analysis for GNSS Kinematic Positioning ........ 35
3.1 Introduction ................................................................................................................ 35
3.2 Nuisance Parameter Elimination in Least Squares ..................................................... 35
3.2.1 Classic Least Squares Adjustment .................................................................... 36
3.2.2 Nuisance Parameter Elimination in Least Squares ........................................... 37
3.3 Two-way Kalman Filter.............................................................................................. 40
3.3.1 Classic Kalman Filter ....................................................................................... 41
3.3.2 Two-way Kalman Filter .................................................................................... 42
3.4 Precision and Reliability Evaluation .......................................................................... 43
3.4.1 Static Experiment ............................................................................................. 43
3.4.2 Kinematic Experiment ...................................................................................... 44
3.4.3 Reliability ......................................................................................................... 47
3.5 Influence of the GNSS Receiver Clock in Accuracy Evaluation ............................... 49
3.5.1 Receiver Clock Jump ........................................................................................ 49
3.5.2 Influence of Receiver Clock in Highly Dynamic Positioning .......................... 51
3.6 Conclusions ................................................................................................................ 53
4 Kinematic Positioning Based on Multiple Reference Stations ...................................... 55
4.1 Introduction ................................................................................................................ 55
4.2 Double Difference Positioning ................................................................................... 55
4.2.1 Principle of Classic Double Difference Positioning ......................................... 55
4.2.2 Experiment and Analysis .................................................................................. 58
4.3 Multiple Reference Stations Kinematic Positioning Based on an A priori Constraint 61
4.3.1 Principle of Multiple Reference Stations Kinematic Positioning Based on an A
priori Constraint ............................................................................................... 62
4.3.2 Experiment and Analysis .................................................................................. 64
4.4 Kinematic Positioning Based on Robust Estimation .................................................. 68
4.4.1 Principle of Robust Kalman Filter .................................................................... 69
4.4.2 Experiment and Analysis .................................................................................. 71
4.5 Conclusions ................................................................................................................ 78
5 Kinematic Positioning Based on Multiple Kinematic Stations ..................................... 81
5.1 Introduction ................................................................................................................ 81
5.2 Kinematic Positioning Based on Multiple Kinematic Stations with Multiple
Reference Stations ...................................................................................................... 82
5.3 Kinematic Positioning Based on an A priori Distance Constraint .............................. 83
- xi -
5.3.1 The A priori Distance Constraints .................................................................... 83
5.4 Kinematic Positioning Based on a Common Tropospheric Parameter....................... 84
5.4.1 Tropospheric Delay Estimation ........................................................................ 85
5.4.2 A Common Tropospheric Wet Delay Parameter ............................................... 86
5.5 Experiment and Analysis ............................................................................................ 87
5.6 Conclusions ................................................................................................................ 94
6 Kinematic Positioning Based on GNSS Integration ...................................................... 95
6.1 Introduction ................................................................................................................ 95
6.2 GPS and GLONASS Integration ................................................................................ 96
6.3 GNSS Integration Based on Helmert’s Variance Components Estimation ................ 98
6.4 Experiment and Analysis .......................................................................................... 100
6.4.1 Kinematic Experiment .................................................................................... 100
6.4.2 Static Experiment ........................................................................................... 105
6.5 Conclusions .............................................................................................................. 110
7 GNSS Velocity Determination Based on the Doppler Effect ........................................ 111
7.1 Introduction ............................................................................................................... 111
7.2 Principle of Doppler-Based Velocity Determination ................................................. 111
7.2.1 Velocity Determination Based on GNSS Raw Doppler Observations ........... 112
7.2.2 Velocity Determination Based on Doppler Observations Derived from GNSS
Carrier Phase Measurements .......................................................................... 114
7.3 Velocity Determination Based on GNSS Integration and Robust Estimation .......... 115
7.4 Experiment and Analysis .......................................................................................... 116
7.4.1 Static Experiment ........................................................................................... 116
7.4.2 Kinematic Experiment .................................................................................... 120
7.4.3 Importance of Precise Velocity Determination for Airborne Gravimetry ....... 127
7.5 Conclusions .............................................................................................................. 130
8 Software Development and Application ....................................................................... 133
8.1 Introduction .............................................................................................................. 133
8.2 HALO_GNSS Software Development ..................................................................... 133
8.2.1 Characteristics of the HALO_GNSS Software............................................... 133
8.2.2 Structure Design of Software ......................................................................... 134
8.3 HALO_GNSS Software Applications ...................................................................... 137
8.3.1 Lake Müritz Shipborne Gravimetry ............................................................... 137
8.3.2 The GEOHALO Italy Mission ....................................................................... 138
8.3.3 Lake Constance Shipborne Gravimetry .......................................................... 139
8.3.4 Baltic Sea Shipborne Gravimetry ................................................................... 140
8.4 Conclusions .............................................................................................................. 141
- xii -
9 Summary and Future Work .......................................................................................... 143
9.1 Summary .................................................................................................................. 143
9.2 Future Work .............................................................................................................. 144
References ............................................................................................................................. 147
Acknowledgments ................................................................................................................ 157
- xiii -
List of Figures
Figure 1.1: Basic principle of the airborne gravimetry (GFZ, 2013)........................................ 3
Figure 1.2: The HALO aircraft (DLR, 2013b) ......................................................................... 4
Figure 1.3: GPS control segment (NOAA, 2014a) ................................................................... 8
Figure 1.4: BDS space constellation (BDS-OS-PS, 2013) ..................................................... 11
Figure 2.1: The principle of troposphere delay (Wang et al., 2009b) ..................................... 28
Figure 3.1: The IGS tracking network (IGS, 2014) ................................................................ 44
Figure 3.2: Equipment for the horizontal movement experiment of a GNSS antenna ........... 45
Figure 3.3: Equipment for the vertical movement experiment of a GNSS antenna ............... 46
Figure 3.4: The continuous steering of GNSS receiver 1 (NOVATEL OEM4) ...................... 50
Figure 3.5: The millisecond jumping of GNSS receiver 2 (JAVAD DELTA G3T) ................ 50
Figure 3.6: The millisecond jumping of GNSS receiver 3 (JAVAD DELTA G3T) ................ 51
Figure 3.7: Variation of the estimated distance between two mechanically fixed antennas,
computed at the same time tags ........................................................................... 52
Figure 3.8: Variation of the estimated distance between two mechanically fixed antennas,
computed at interpolated time epochs .................................................................. 53
Figure 4.1: An illustration for GNSS DD positioning (Stombaugh and Clement, 1999) ....... 56
Figure 4.2: The results of the kinematic positioning for the short baseline POTS-KIN1
during the vertical antenna movement experiment ............................................. 59
Figure 4.3: The long baseline (MATE-KIN1, 1330 km) (IGS, 2014) .................................... 60
Figure 4.4: The results of the kinematic positioning for the long baseline MATE-KIN1
during the vertical antenna movement experiment ............................................. 61
Figure 4.5: The positions of the selected GNSS antennas on the HALO aircraft ................... 65
Figure 4.6: HALO aircraft flight trajectory and the position of reference stations ................ 66
Figure 4.7: The distance between the two antennas AIR5 and AIR6, the trajectories for
these antennas were estimated by using only one reference station (REF6) ........ 67
Figure 4.8: The distance between the two antennas AIR5 and AIR6, the trajectories for
these antennas were estimated including the three reference stations REF6,
RENO and ASCE ................................................................................................. 68
Figure 4.9: The weight factors of Kalman filter robust estimation ......................................... 70
Figure 4.10: The time series for the baseline length between AIR5 and AIR6 based on
pseudorange positioning results without robust estimation ................................. 72
Figure 4.11: The time series for the baseline length between AIR5 and AIR6 based on
pseudorange positioning results with robust estimation....................................... 72
- xiv -
Figure 4.12: The probability distribution of standard residual for the GNSS pseudorange
observations.......................................................................................................... 73
Figure 4.13: The time series for the baseline length between AIR5 and AIR6 based on
carrier phase positioning results without robust estimation ................................. 76
Figure 4.14: The time series for the baseline length between AIR5 and AIR6 based on
carrier phase positioning results with robust estimation ...................................... 77
Figure 4.15: The probability distribution of the standard residuals for the GNSS Carrier
phase observations ............................................................................................... 78
Figure 5.1: Multiple GNSS antennas mounted on the HALO aircraft (DLR, 2013a) ............ 81
Figure 5.2: The relative positions of the kinematic GNSS antennas on the ship .................... 87
Figure 5.3: The track of the ship (green curve) and the positions of the reference stations
in a small region (yellow triangles) during the Baltic Sea shipborne
gravimetric campaign on June 18. 2013 ............................................................... 88
Figure 5.4: The track of the ship (green curve) and the positions of reference stations
(yellow triangles) for the Baltic Sea shipborne gravimetric experiment on
June 18. 2013 ....................................................................................................... 89
Figure 5.5: The distance between two antennas, KIN1 and KIN3; the trajectories of
these antennas were estimated by using two reference stations, 0801 and
0775 ...................................................................................................................... 90
Figure 5.6: The distance between two antennas, KIN1 and KIN3; the tracks of these
antennas were estimated by using two reference stations, WARN and POTS ..... 91
Figure 5.7: The results for KIN1 for the large region with two independent tropospheric
delay parameters compared with scheme 1 (differences) .................................... 91
Figure 5.8: The results for scheme 3 for KIN1 for the large region with two independent
tropospheric delay parameters and distance constraints compared with
scheme 1 (differences) ......................................................................................... 92
Figure 5.9: The results of scheme 4 for KIN1 with reference stations in the large region
with a common tropospheric delay parameter for KIN1 and KIN3 compared
with scheme 1 (differences) ................................................................................. 93
Figure 5.10: The results of scheme 5 for KIN1 with reference stations in the large region
with a common tropospheric delay parameter and distance constraints
compared with scheme 1 (differences) ................................................................. 93
Figure 6.1: GNSS integration (OpeNaviGate, 2012) .............................................................. 95
Figure 6.2: The HALO aircraft flight trajectory on June 6, 2012 and the location of the
selected reference stations .................................................................................. 101
- xv -
Figure 6.3: Number of the selected satellites of GPS (blue line), GLONASS (green line)
and GPS+GLONASS (red line) for the kinematic experiment (GEOHALO
flight on June 6, 2012) ....................................................................................... 102
Figure 6.4: The time series for the apparent baseline length variation between AIR5 and
AIR6 based on GPS as a single GNSS system ................................................... 102
Figure 6.5: The time series for the apparent baseline length variation between AIR5 and
AIR6 based on GLONASS as a single GNSS system........................................ 103
Figure 6.6: The time series for the apparent baseline length variation between AIR5 and
AIR6 based on GPS and GLONASS as an integrated system with 1:1
weights ............................................................................................................... 103
Figure 6.7: The ratios between the weights for the GPS and GLONASS observations ....... 104
Figure 6.8: The time series for the apparent baseline length variation between AIR5 and
AIR6 based on GPS and GLONASS as an integrated system with Helmert’s
VCE .................................................................................................................... 104
Figure 6.9: The number of selected satellites of GPS (blue line), GLONASS (green line)
and GPS+GLONASS (red line) for the static experiment (IGS station TITZ
and FFMJ on January 1, 2013) ........................................................................... 106
Figure 6.10: The differences between the IGS result and the GPS kinematic positioning
results for the baseline TITZ – FFMJ (scheme 1) .............................................. 107
Figure 6.11: The differences between the IGS result and the GLONASS kinematic
positioning results for the baseline TITZ – FFMJ (scheme 2) ........................... 107
Figure 6.12: The differences between the IGS result and the GPS+GLONASS (with 1:1
weights) kinematic positioning results for the baseline TITZ – FFMJ
(scheme 3) .......................................................................................................... 108
Figure 6.13: The differences between the IGS result and the GPS+GLONASS (with
Helmert weights) kinematic positioning results for the baseline TITZ –
FFMJ (scheme 4)................................................................................................ 108
Figure 6.14: The ratios between the weights for the GPS and GLONASS observations
(scheme 4) .......................................................................................................... 109
Figure 7.1: The velocity from raw GPS Doppler observations at the static station RENO
on June 6, 2012 (scheme 1) ................................................................................ 117
Figure 7.2: The velocity from raw GLONASS Doppler observations at the static station
RENO on June 6, 2012 (scheme 2) .................................................................... 117
Figure 7.3: The velocity from the integrated raw GPS and GLONASS Doppler
observations at the static station RENO on June 6, 2012 (scheme 3) ................ 118
- xvi -
Figure 7.4: The velocity from the integrated GPS and GLONASS carrier phase derived
Doppler observations on June 6, 2012 for the static station RENO, without
the usage of the robust estimation (scheme 4) ................................................... 118
Figure 7.5: The velocity from GPS and GLONASS integrated carrier phase derived
Doppler observations on June 6, 2012 for the static station RENO, with the
usage of the robust estimation (scheme 5) ......................................................... 119
Figure 7.6: The components for the flight trajectory of the HALO aircraft on June 6,
2012 .................................................................................................................... 121
Figure 7.7: The velocity components for the flight of the HALO aircraft on June 6, 2012 . 121
Figure 7.8: The acceleration components for the flight of the HALO aircraft on June 6,
2012 .................................................................................................................... 122
Figure 7.9: The number of selected satellites of GPS (blue line), GLONASS (green line)
and GPS+GLONASS (red line) ......................................................................... 122
Figure 7.10: The differences of the velocity estimates between AIR5 and AIR6 ................. 123
Figure 7.11: The difference of the GPS-only velocity results between AIR5 and AIR6 in
the period of a straight smooth flight (scheme 1) ............................................... 124
Figure 7.12: The difference of the GLONASS-only velocity results between AIR5 and
AIR6 in the period of a straight smooth flight (scheme 2) ................................. 124
Figure 7.13: The difference of the velocity results using the integrated GPS and
GLONASS between AIR5 and AIR6 in the period of a straight smooth flight
(scheme 3) .......................................................................................................... 125
Figure 7.14: The differences between the estimated velocities for the raw Doppler
observations and for the carrier phase derived Doppler at station AIR5
without the usage of the robust estimation (integrated GPS and GLONASS
data) .................................................................................................................... 126
Figure 7.15: The differences between the estimated velocities for the raw Doppler
observations and for the carrier phase derived Doppler at station AIR5 by
applying the robust estimation (integrated GPS and GLONASS data) .............. 127
Figure 7.16: Vertical kinematic acceleration of the HALO aircraft as calculated directly
from GNSS observations .................................................................................... 128
Figure 7.17: The black curve is the GNSS-based disturbing vertical kinematic
accelerations of the HALO aircraft after the application of a low-pass filter
to the noisy result as displayed in Figure 7.16. The red curve is a reference
gravity signal computed from the global gravity field mode EIGEN-6C4
minus normal-gravity. The green curve is the CHEKAN measure vertical
acceleration with all corrections except of that from GNSS. ............................. 129
- xvii -
Figure 7.18: The final gravity signal (blue curve) compared with the reference gravity
field model EIGEN-6C4 (red curve) ................................................................. 130
Figure 8.1: The GNSS data processing flowchart of HALO_GNSS software ..................... 135
Figure 8.2: The Lake Müritz location and the tracks of the ship .......................................... 137
Figure 8.3: The ship used and the position of GNSS receiving equipments ........................ 138
Figure 8.4: The GEOHALO mission tracks (red lines) and GNSS reference stations (red
and green dots) in Italy ....................................................................................... 139
Figure 8.5: All tracks of the ship in the Constance Lake in October 2012 ........................... 140
Figure 8.6: All tracks of the ship in the Baltic Sea gravimetric campaign in June 2013 ...... 141
- xix -
List of Tables
Table 1.1: GPS satellite block characteristics in April 2014 (NOAA, 2014b).......................... 6
Table 1.2: GLONASS satellite characteristics in April 2014 (IAC, 2014) ............................... 9
Table 2.1: The relation between admissible orbit error and baseline length for an 1 cm
baseline error for differential GNSS positioning ................................................. 25
Table 2.2: IGS combined orbit and clock products (IGS, 2013) ............................................ 25
Table 4.1: The statistical results of the difference between the long and the short
baseline (Unit: m) ................................................................................................. 61
Table 4.2: Hardware list of the selected stations from the GEOHALO project...................... 65
Table 4.3: The Statistical results for the distance between the two antennas AIR5 and
AIR6 (Unit: m) ..................................................................................................... 67
Table 4.4: The statistical results of the baseline length between AIR5 and AIR6 from
GNSS pseudorange positioning (Unit: m) ........................................................... 71
Table 4.5: The statistical results for the reliability analysis without robust estimation
(for the GNSS DD observations at time epoch 12h:10m:09s) ............................. 74
Table 4.6: The statistical results for the reliability analysis when applying robust
estimation (for the GNSS DD observations at time epoch 12:10:09) ................. 75
Table 4.7: The statistical results of GNSS carrier phase positioning (Unit: m) ...................... 77
Table 5.1: Hardware list of the selected stations from Baltic Sea gravimetric experiment .... 89
Table 5.2: The statistical results for the baseline length between KIN1 and KIN3 (Unit:
m) ......................................................................................................................... 91
Table 5.3: Statistical results of each scheme compared with scheme 1 (Unit: mm) ............... 94
Table 6.1: The statistical results of schemata 1-4 ................................................................. 105
Table 6.2: Hardware list of the selected IGS sites ................................................................ 106
Table 6.3: Statistical results for the IGS station TITZ (Unit: mm) ....................................... 109
Table 7.1: The statistical results of the velocity determination (static experiment) for the
station RENO in the schemata 1-5 (Unit: mm/s) ............................................... 119
Table 7.2: Statistical results for the velocity determination in the kinematic experiment
for schemata 1-3 (Unit: mm/s) ........................................................................... 125
Table 7.3: The statistical results of scheme 4 and 5 as plotted in figures 7.14 and 7.15
(Unit: mm/s) ....................................................................................................... 127
Table 7.4: The statistic results of the comparison between the calculated airborne gravity
signal and the reference gravity field mode EIGEN-6C4 (unit: mGal) ............. 130
- xxi -
List of Abbreviations
AC Analysis Center
AdV Arbeitsgemeinschaft der Vermessungsverwaltungen der Länder der
Bundesrepublik Deutschland (i.e. the Working Committee of the
Surveying Agencies of the States of the Federal Republic of Germany)
ARP Antenna Reference Point
BDS the Chinese BeiDou navigation Satellite System
BKG Bundesamt für Kartographie und Geodäsie (i.e. Federal Agency for
Cartography and Geodesy, Germany)
BMBF Bundesministerium für Bildung und Forschung (i.e. German Federal
Ministry for Education and Research)
C/A Coarse/Acquisition code
CDDIS Crustal Dynamics Data Information System
CDMA Code Division Multiple Access
CODE Centre of Orbit Determination in Europe
CORS Continuously Operating Reference Stations
CS Commercial Service
DCB Differential Code Biases
DD Double Difference
DF Dual-Frequency
DFG Deutsche Forschungsgemeinchaft (i.e. German Research Foundation)
DLR Deutsches Zentrum für Luft- und Raumfahrt (i.e. German Aerospace
Centre)
DOD the U.S. Department of Defense
ECMWF European Centre for Medium-Range Weather Forecasts
ESA European Space Agency
EU European Union
EX_WL Extra-widelane Linear Combination
FDMA Frequency Division Multiple Access
Galileo the European Union satellite navigation system
GEO Satellites in Geostationary Orbit
GEOHALO GEOscience High Altitude and LOng range research project
GFZ Helmholtz-Zentrum Potsdam Deutsches-GeoForschungsZentrum (i.e.
Helmholtz-Centre Potsdam-GFZ German Research Centre for
- xxii -
Geosciences)
GIM Global Ionospheric Map
GLONASS the Russian GLObal NAvigation Satellite System
GLOT GLONASS Time
GMF Global Mapping Function
GNSS Global Navigation Satellite System
GPS Global Positioning System
GPST GPS Time
HALO High Altitude and LOng range airborne gravimetric project
HGF Helmholtz-Gemeinschaft forscht (i.e. Helmholtz Association)
IAR Integer Ambiguity Resolution
IGR IGS Rapid orbit
IGS International GNSS Service
IGSO Inclined Geosynchronous Orbit
IGU IGS Ultra-rapid orbit
INS Inertial Navigation System
IOV In-Orbit Validation
IRNSS Indian Regional Navigation Satellite System
ITRF International Terrestrial Reference Frame
KIT Karlsruher Institut für Technologie (i.e. Karlsruhe Institute for
Technology)
KNITs Coordination Scientific Information Centre of the Russia
LC Ionosphere-free Linear Combination
LEO Low Earth Orbit
MDB Minimal Detectable Bias
MEO Medium altitude Earth Orbit
MET Meteorology
MPG Max-Planck-Gesellschaft (i.e. Max Planck Society)
MW_WL MW Widelane Linear Combination
NASA National Aeronautics and Space Administration
NOAA National Oceanic and Atmospheric Administration
NRTK Network-based Real Time Kinematic positioning
OCS Operational Control Segment
OS Galileo are Open Service
PCO Phase Centre Offsets
PCV Phase Centre Variation
PLL Phase Lock Loop
- xxiii -
PNT Positioning, Navigation and Timing
PPP Precise Point Positioning
PPP-RA Precise Point Positioning Regional Augmentation
PPS Precise Positioning Service
PRS Public Regulated Service
PRN Pseudo-Random Noise
PZ-90 Parametry Zemli 1990
QZSS Quasi-Zenith Satellite System of Japan
RMS Root Mean Square
RTK Real-Time Kinematics
SAPOS Satellite Positioning Service of the German State Survey
SAR Search and Rescue
SD Single Difference
SDev Standard Deviation
SoL Safety-of-Life Service
SP3 IGS Standard Product 3
SPS Standard Positioning Service
STD Slant Total Delay
UPD Un-calibrated Phase Delays
USNO US Naval Observatory
USSR Union of Soviet Socialist Republics
UTC Coordinated Universal Time
VCE Variance Components Estimation
WL Widelane Combination
WGS-84 World Geodetic System 1984
ZHD Zenith Hydrostatic Delay
ZTD Zenith Total Delay
ZWD Zenith Wet Delay
- xxiv -
- 1 -
1 Introduction
1.1 Motivation
Precise kinematic position and velocity determination based on the Global Navigation Satellite
System (GNSS) is widely used in many scientific research and engineering fields, such as ge-
odesy, spatial science, geophysics and meteorology (Chen, 1998; Leick, 2004, p. 2; Hofmann-
Wellenhof et al., 2008, pp. 1-12; Xu, 2007, p. 1). Meanwhile, measuring the Earth’s gravity
field is one of the most important activities in scientific and economic applications, such as
geodesy, geophysics, exploration purposes, geoid determination and satellite orbits prediction
(Kwon and Jekeli, 2001; Kreye and Hein, 2003). In this context, airborne gravimetry plays a
very important role in recovering the Earth’s gravity field in the range of medium to high fre-
quencies, which then closes the gap between the terrestrial gravity filed measurements on the
ground and the global gravity models based on the satellite gravimetry at wavelengths be-
tween 1 and 100-200 km (Hein, 1995; Kwon and Jekeli, 2001).
In the fields of airborne gravimetry, GNSS position and velocity determination plays a sig-
nificant role as well (Forsberg and Olesen, 2010), since the state information of a kinematic
platform with an airborne gravimeter can be obtained independently from the GNSS observa-
tions. The trajectory and attitude of such a kinematic platform are indispensable information
for analysing its airborne gravimetry data. The acceleration information as derived from the
position and/or velocity information for such a kinematic platform can be used to separate the
disturbing kinematic accelerations affecting the airborne platform from the gravitational signal.
Therefore, the estimation of accurate state information for such a kinematic platform by pre-
cise GNSS position and velocity determination is the key factor for any successful
implementation of airborne gravimetry (Torge, 1989; Hannah, 2001; Kreye and Hein, 2003).
More details about the principle of airborne gravimetry are given in Section 1.2.1.
However, there are still many remaining challenges in the application of GNSS position
and velocity determination in airborne gravimetry. Firstly, the novel airborne gravimetric pro-
ject GEOHALO on the German High Altitude and Long Rang Research Aircraft HALO (cf.
Section 1.2.2), which is characterized by its high-altitude, long-range and high dynamics, re-
quires more accurate and reliable GNSS-based state information for the kinematic platform
than before in “traditional” airborne gravimetry. Secondly, the rapid development of multiple
GNSS systems (GPS, GLONASS, etc.) provides more information and new opportunities (cf.
Section 1.2.3) for improving the accuracy and reliability of the GNSS results for kinematic
- 2 -
platforms. In this study, the algorithms and methods for the multiple GNSS systems are ap-
plied and analyzed. Finally, most of the currently available GNSS software tools incl. the
therein adopted methods do not meet the requirements for HALO. For instance, the usage of
multiple GNSS systems is not possible in some software systems like GAMIT (Herring et al.,
2010) or HALO_GPS (Wang et al., 2010). Furthermore, some of the GNSS software systems
like the RTKLIB (Takasu, 2013) cannot exploit information of multiple reference stations for
large regions. And most of the GNSS software cannot provide velocity information of a kine-
matic platform. Therefore, new algorithms and software have to be developed to fulfill the
requirements of the novel airborne gravimetry on HALO. Thus, the objective of this thesis is
the development, application and evaluation of specific precise and reliable GNSS algorithms
and software of kinematic GNSS data analysis for airborne gravimetry on the high-altitude and
long-range aircrafts like HALO.
1.2 Background
1.2.1 The Role of GNSS Application in Airborne Gravimetry
Airborne gravimetry is an emerging technology that provides a powerful tool for both geodesy
and geophysical exploration. Using airborne techniques, detailed local and regional gravity
field information can be collected in a rapid and cost-effective manner (Bruton et al., 1999).
For instance, it is used for geoid determination and for mapping of gravity anomalies in the
contest of geophysical exploration (Forsberg and Olesen, 2010).
The idea of using an airborne platform for gravity measurements is not completely new. It
has been attempted since the 1960s (Lacoste, 1967; Schwarz and Wei, 1995), and recognized
that if an appropriate level of accuracy could be achieved airborne gravimetry would be vastly
superior in economy and efficiency to point-wise terrestrial gravimetry (Alberts, 2009). In the
beginning, airborne gravimetry did not become a major tool for gravity filed mapping, alt-
hough the first experiments gave promising results (Thompson and LaCoste, 1960). The
reason is the disturbing kinematic accelerations of the aircraft which are very difficult to de-
termine due to the lack of sufficient navigation technologies at that time (Hannah, 2001).
Nevertheless, in the late 1980s and early 1990s this situation rapidly improved with the avail-
ability of the Global Positioning System (GPS) and the development of kinematic positioning
using its carrier phase signals. Meanwhile, the accuracy of GNSS-based trajectory and veloci-
ty determination for airborne gravimetry has reached a useful level (Torge, 1989; Hannah,
2001).
In principle, airborne gravimetry requires both a device determining the sum of all gravity
and kinematic accelerations affecting the airborne platform, plus a positioning system (such as
- 3 -
radar, laser altimeter, or GNSS) for the determination of the kinematic accelerations alone
(Hannah, 2001). The gravity vector is determined by differencing between them. The basic
principle is shown in Figure 1.1.
Figure 1.1: Basic principle of the airborne gravimetry (GFZ, 2013)
In Figure 1.1, the yellow arrow means the superposition of all vertical accelerations (
a
)
measured by the gravimeter inside of the airplane. The green arrow indicates the non-
gravitational kinematic acceleration (
x
) to be estimated from the GNSS-recorded trajectory
and/or velocity. The difference between the measured and the kinematic accelerations is the
gravity acceleration (
g
) at the trajectory (GFZ, 2013) which can be expressed in the follow-
ing basic formula as
g a x
. (1.1)
Airborne gravimetry techniques face a number of difficulties that typically do not occur in
the classical field of gravimetry on the ground or on slowly moving vessels like ships. Some
of these difficulties are given by Torge (1989, p. 287) and Hannah (2001):
The high platform velocities require short averaging times and high position and par-
ticularly velocity accuracies.
Elevation changes (uncertainties) directly influence the measured value of gravity.
- 4 -
The power of the gravity field (especially the short-wavelength part), attenuates
steadily with altitude thus causing the signal-to-noise ratio to decrease at a higher alti-
tude. In general, an increase of the altitude of the aircraft decreases the resolution
capability of the gravimetric system at short wavelengths.
As mentioned above, GNSS is used for the reduction of the disturbing non-gravitational
accelerations of the measuring kinematic platform. GNSS allows for a highly accurate estima-
tion of position and velocity. Hence, precise GNSS position and velocity determination is
widely used, and plays a significant role in the field of airborne gravimetry (Kreye and Hein,
2003; Petrovic et al., 2013).
1.2.2 The New Airborne Gravimetric Project
The High Altitude and LOng Range (HALO) project is a joint new project of several German
institutes for atmospheric research and Earth observation (DLR, 2013b). The main project
partners are the German Federal Ministry for Education and Research (BMBF), the German
Centre for Aeronautics and Space Research (DLR), the German Research Foundation (DFG),
the Max Planck Society (MPG) and the Helmholtz Association (HGF) incl. its members Re-
search Centre Juelich (Forschungszentrum Jülich), Karlsruhe Institute for Technology (KIT)
and German Research Centre for Geosciences (GFZ).
The new German HALO aircraft is a business jet G550 produced by the Gulfstream Aero-
space Corporation, shown in Figure 1.2. The aircraft has a spacious cabin for installing various
atmospheric and remote sensing equipment and instruments, which makes the HALO aircraft
an extraordinary platform for atmospheric and geophysical research.
Figure 1.2: The HALO aircraft (DLR, 2013b)
The goal of the GEOHALO project is an aircraft-based geodedic-geophysical survey over
the whole territory of Italy and adjacent regions of the Mediterranean using a variety of meas-
urement systems and remote sensing equipment mounted on the HALO aircraft. The main
- 5 -
fields of the GEOHALO project include: GNSS precise position and velocity determination,
GNSS reflectometry (Semmling et al., 2013; Semmling et al., 2014), airborne gravimetry
(Petrovic et al., 2013), magnetometry and laser altimetry (Scheinert et al., 2013).
From June 6 to 12, 2012, the GEOHALO mission conducted 4 flights, which started from
the German Aerospace Centre (DLR) in Oberpfaffenhofen Germany, scanning around over
Italy and the adjacent regions of the Mediterranean. Every flight lasted up to 10 hours and the
flown trajectories had a total length of 16150 kilometres. The speed of the aircraft on the sur-
vey tracks was approximately 450 km/h. The flight altitude was about 3500 m (Scheinert et al.,
2013; Heyde et al., 2013). GNSS precise position and velocity determination plays an im-
portant role in the whole project, since all remote sensing instruments and measuring devices
are mounted on this kinematic platform. Therefore, precise information about the state of the
platform is a prerequisite and a key factor for the successful accomplishment of the entire re-
search project. In particular, the high speed of HALO is a challenge for the airborne
gravimetry on this aircraft since “classical” airborne gravimetry is usually done on a smaller
and slower aircraft.
1.2.3 The Development Status of GNSS
The advent of satellite navigation technologies has greatly changed the way of human's life
and production. With the development during a few decades, GNSS has been applied widely
in many areas, such as air, sea and land Positioning, Navigation and Timing (PNT), low earth
orbit (LEO) satellite orbit determination, static and kinematic positioning, flight-state monitor-
ing, sea surface monitoring (Galas et al., 2013; Semmling et al., 2014), atmospheric profiling
(Haase et al., 2014), as well as surveying, etc. GNSS has become a necessity in daily life, in-
dustry, research and education (Xu, 2007, p. 1). The application of GNSS is limited by the
human's imagination only.
GNSS includes several different satellite navigation systems, namely the Global Position-
ing System (GPS) developed by the United States, the GLObal NAvigation Satellite System
(GLONASS) developed by Russia, Galileo developed by the European Union, the BeiDou
navigation Satellite System (BDS) developed by China, the Quasi-Zenith Satellite System
(QZSS) developed by Japan and the Indian Regional Navigation Satellite System (IRNSS)
developed by India. As of April 2014, only GPS and GLONASS were fully globally opera-
tional. Galileo and BDS are in the process of growth. QZSS and IRNSS are regional satellite
navigation system. All kinds of GNSS consist of three segments: the space, control and user
segment (Hofmann-Wellenhof et al., 2008, pp. 6, 7; Yang et al., 2011a).
- 6 -
1.2.3.1 The Development Status of GPS
The Global Positioning System (GPS) is the earliest satellite navigation system. It was de-
signed and built, and is operated and maintained by the U.S. Department of Defence (DOD).
GPS permits land, sea, airborne and space users to determine their three-dimensional position,
velocity and time 24 hours a day, in all weather, anywhere in the world with a precision and
accuracy far better than other radio navigation systems available today or in the foreseeable
future (NAVCEN, 2014).
The first GPS satellite was launched in 1978 and the system became fully operational in the
mid-1990s. At the moment, the GPS constellation consists of about 30 satellites in six orbital
planes with about five satellites in each plane. The GPS satellites fly in Medium Earth Orbits
(MEO) at an altitude of approximately 20200 km. The orbital planes are inclined 55 degrees
with respect to the equator. Each GPS satellite is in a nearly circular orbit with a semi-major
axis of 26578 km and a period of about twelve hours, and circles the Earth twice a day (Xu,
2007, p. 2; Hofmann-Wellenhof et al., 2008, p. 322; NOAA, 2014b).
In April 2014, the GPS constellation consisted of 31 satellites, divided into four groups: 8
Block IIA, 12 Block IIR, 7 Block IIR(M) and 4 Block IIF. In the future, GPS satellite Block
III will be launched, due to the GPS modernization. The major characteristics of each Block
are shown in Table 1.1 (NOAA, 2014b).
Table 1.1: GPS satellite block characteristics in April 2014 (NOAA, 2014b)
Block
IIA
IIR
IIR(M)
IIF
III
Operational
Satellite
Number
8
12
7
4
0
Launch
Period
1990-1997
1997-2004
2005-2009
Since 2010
Planning 2016
Design
Life
7.5-year
7.5-year
7.5-year
12-year
15-year
CDMA
Signals
L1, L2
P1, P2
C/A
L1, L2
P1, P2
C/A
L1, L2
P1, P2
C/A
M1, M2
L2C
L1,L 2, L5
P1, P2
C/A
M1, M2
L2C, L5C
L1, L2, L5
P1, P2
C/A
M1, M2
L1C, L2C, L5C
- 7 -
The GPS satellites listed in Table 1.1 use the Code Division Multiple Access (CDMA)
technique to transmit data on the L-band frequencies, L1 (1575.42 MHz), L2 (1227.60 MHz)
and L5 (1176.45 MHz), depending on the block type of the particular satellite. The L1, L2 and
L5 carrier frequencies are generated by multiplying the fundamental frequency (10.23 MHz)
by 154, 120 and 115, respectively. Pseudorandom noise (PRN) codes, along with satellite eph-
emerides, ionospheric models and satellite clock corrections are superimposed onto the carrier
frequencies L1, L2 and L5. The measured transmitting times of the signals that travel from the
satellites to the receivers are used to compute the pseudoranges. The 1st civil signal Coarse-
Acquisition (C/A) code, sometimes called the Standard Positioning Service (SPS), is a pseu-
dorandom noise code that is modulated onto the L1 carrier. The precision (P) code, sometimes
called the Precise Positioning Service (PPS), is modulated onto the L1, L2 and L5 carriers al-
lowing for the removal of the effects of the ionosphere (Xu, 2007, p. 3). The transmission of
the military M-code started with the launch of the Block IIR(M) satellites. The main character-
istics of this military code are a higher resistance against jamming, increased navigation
performance, higher security based on the new cryptography algorithms, and the possibility of
higher transmission power, denoted as flex-power. The M-code is modulated onto the L1 and
L2 carrier frequencies. The new civil signals L1C, L2C and L5C are modulated onto L1, L2
and L5, respectively. The 2nd civil signal L2C has been designed to meet in particular commer-
cial needs. The 3rd civil signal L5C has been designed to especially meet the requirements of
safety-of-life (SoL) applications. The 4th civil signal L1C has been designed to enable interop-
erability between GPS and the other international satellite navigation systems, such as Galileo,
BDS, QZSS and IRNSS (Hofmann-Wellenhof et al., 2008, p. 336; NOAA, 2014b).
The Operational Control Segment (OCS) started to operate in 1985 and consists of a global
network of ground facilities that track the GPS satellites, monitor their transmissions, perform
analyses, and send commands and data to the constellation. The current OCS includes a master
control station in Colorado Springs, an alternate master control station, 12 command and con-
trol antennas, and 16 monitoring sites. The locations of these facilities are shown in Figure 1.3
(NOAA, 2014a). The monitoring stations track all GPS satellites in view and collect ranging
information from the satellite broadcasts. The monitoring stations send the information they
collect from each of the satellites back to the master control station which computes precise
satellite orbits. The information is then formatted into updated navigation messages for each
satellite. The updated information is transmitted to each satellite via the ground antennas,
which also transmit and receive satellite control and monitoring signals (Hofmann-Wellenhof
et al., 2008, pp. 324-327; Xu, 2007, p. 3; Marreiros, 2012; NOAA, 2014a).
- 8 -
Figure 1.3: GPS control segment (NOAA, 2014a)
1.2.3.2 The Development Status of GLONASS
The GLObal Navigation Satellite System (GLONASS) is a satellite navigation system devel-
oped by the former Union of Soviet Socialist Republics (USSR) for military use. Now,
GLONASS is managed by the Russian Space Forces and the system is operated by the Coor-
dination Scientific Information Centre (KNITs) of the Ministry of Defence of the Russian
Federation (Xu, 2007, p. 3).
The first GLONASS satellite was launched into orbit in 1982. The system consists of 24
satellites in three orbital planes, with three in-orbit spares. The ascending nodes of three or-
bital planes are separated by 120 degrees, and the satellites within the same orbit plane are
equally spaced by 45 degrees. The arguments of latitude of satellites in equivalent slots in two
different orbital planes differ by 15 degrees. Each satellite operates in a nearly circular orbit
with a semi-major axis of 25510 km. Each orbital plane has an inclination angle of 64.8 de-
grees, and each satellite completes an orbit in approximately 11 hours 16 minutes (Xu, 2007, p.
3).
As of April 2014, 29 satellites were in orbit with 24 GLONASS-M satellites operational
(IAC, 2014). In the future, GLONASS satellites GLONASS-K1, GLONASS-K2 and
GLONASS-KM will be launched, due to the ongoing GLONASS modernization. The major
characteristics of each type are shown in Table 1.2 (Revnivykh, 2012; IAC, 2014).
- 9 -
Table 1.2: GLONASS satellite characteristics in April 2014 (IAC, 2014)
Type
M
K1
K2
KM
Operational
Satellite Number
24
0
0
0
Launch Period
2003-2016
2011-2014
2015-2024
Planning 2025
Design Life
7-year
10-year
10-year
12-year
FDMA Signals
G1F, G2F
G1OF, G2OF
G1SF, G2SF
G1F, G2F
G1OF, G2OF
G1SF, G2SF
G1F, G2F
G1OF, G2OF
G1SF, G2SF
G1F, G2F
G1OF, G2OF
G1SF, G2SF
CDMA Signals
G3C
G3OC
(from 2014)
G3C
G3OC
G1C, G2C, G3C
G1OC, G2OC, G3OC
G1SC, G2SC
G1C, G2C, G3C
G1OC, G2OC, G3OC
G1SC, G2SC, G3SC
Interoperability
CDMA Signals
G1CM, G3CM, G5CM
G1OCM, G3OCM, G5OCM
Note
F: denotes FDMA. C: denotes CDMA. CM: denotes Interoperability CDMA.
O: denotes standard-accuracy signal, likes C/A code of GPS.
S: denotes high-accuracy signal, likes P code of GPS.
As seen from Table 1.2, GLONASS implements the Frequency Division Multiple Access
(FDMA) technique to differentiate among the signals of different satellites. In this way, the
GLONASS satellites transmit coded signals in two frequencies located on two frequency
bands, which can be computed by the simple formula (Leick et al., 1998; Habrich, 2000; Dach
et al., 2007)
0
1 1 1
k
f f k f
(1.2)
and
0
2 2 2
k
f f k f
, (1.3)
where
k
denote frequency channel number,
7, , 6k
(since 2005), it is part of the navi-
gation message.
0
11602f
MHz,
0
21246f
MHz,
10.5625f
MHz,
20.4375f
MHz. From 2014 onward, the new GLONASS satellites will feature a full suite of modernized
- 10 -
CDMA signals in the existing G1 (1600.995 MHz), G2 (1248.04 MHz), G3 (1202.25 MHz)
bands, which includes standard-accuracy signal and high-accuracy signal, respectively.
GLONASS-KM satellites will increase the open signals (Interoperability CDMA Signals) to
enable interoperability between GLONASS and the other international satellite navigation sys-
tems. The open signal G1OCM is modulated onto the G1 (1575.42 MHz), similar to
modernized GPS signal L1C and Galileo/BDS signal E1. The open signal G3OCM is modu-
lated onto the G3 (1207.14 MHz), similar to Galileo/BDS signal E5b. The open signal
G5OCM is modulated onto the G5 (1176.45 MHz), similar to the GPS safety of life signal
L5C and Galileo/BDS signal E5a (InsideGNSS, 2010).
The ground control stations of the GLONASS are maintained only on the territory of the
former Soviet Union due to the historical reasons. This lack of global coverage is not optimal
for the monitoring of a global navigation satellite system (Xu, 2007, p. 4).
1.2.3.3 The Development Status of BDS
The BeiDou navigation Satellite system (BDS) is a Chinese global navigation satellite system.
Beidou means the constellation of the Big Dipper (or Great Bear) in Chinese. Following the
deployment and operation of the BeiDou Satellite Navigation Experimental System (BeiDou-
1), China is in the progress of deploying BDS (formerly known as COMPASS and BeiDou-2).
The build-up of BDS is divided into two steps. The first step is establishing regional naviga-
tion and positioning services, consisting of 5 satellites in Geostationary Orbit (GEO), 5 in
Inclined Geosynchronous Orbit (IGSO), and 4 in Medium altitude Earth Orbit (MEO). This
step has been completed by the end of October 2012. The second step will complete the con-
stellation, which comprises 5 GEO, 3 IGSO, 27 MEO satellites by the end of 2020. Then, the
BDS system will provide global navigation services similar to GPS, GLONASS and Galileo
(Zhao et al., 2013; Shi et al., 2013; Montenbruck et al., 2013; BDS, 2014).
In April 2014, the BDS constellation consisted of 14 satellites, divided into three types, 5
GEO, 5 IGSO, and 4 MEO. The GEO satellites are operating in orbits with an altitude of
35786 km and are positioned at 58.75 E, 80 E, 110.5 E, 140 E and 160 E, respectively.
The IGSO satellites are operating in orbits with an altitude of 35786 km and an inclination of
55 w.r.t. to the equatorial plane. The phase difference of the right ascensions of the ascending
nodes of the orbital planes is 120 . The MEO satellites are operating in orbits with an altitude
of 21528 km and an inclination of 55 to the equatorial plane. The satellite recursion period is
13 rotations within 7 days. The respective positions of satellites are shown in Figure 1.4 (BDS-
ICD, 2013; BDS-OS-PS, 2013).
- 11 -
The positioning services transmit signals using B1, B2 and B3. The nominal frequency is
1561.098 MHz, 1207.140 MHz and 1268.520 MHz for B1, B2, B3, respectively (BDS-ICD,
2013).
Figure 1.4: BDS space constellation (BDS-OS-PS, 2013)
1.2.3.4 The Development Status of Galileo
The Galileo satellite navigation system is the first civilian GNSS. It has been created by the
European Union (EU) and the European Space Agency (ESA) to provide a highly accurate
navigation and positioning service designed specifically for civilian purposes and run solely
by civilians.
The fully functional Galileo constellation will consist of 27 satellites and 3 spares, posi-
tioned in three orbital planes with nine equally spaced operational satellites in each plane plus
one inactive spare satellite. The ascending nodes of the orbital planes are equally spaced by
120 degrees. The orbital planes are inclined at 56 degrees. Each Galileo satellite is in a nearly
circular orbit with a semi-major axis of 29600 km (Xu, 2007, p. 4; Galileo-ICD, 2010).
The first Galileo test satellite was launched in December 2005. The first four satellites were
launched into orbit in 2011 and 2012. Four satellites is the minimum number needed to obtain
a unique solution. The full completion of the 30-satellite Galileo system is expected by 2019
(ESA, 2014).
The Galileo navigation signals are transmitted in four frequency bands, E5a (1176.450
MHz), E5b (1207.140 MHz), E6 (1278.570 MHz) and E1 (1575.420 MHz). They provide a
wide bandwidth for the transmission of the Galileo Signals (Galileo-ICD, 2010; ESA, 2014).
- 12 -
Mission and services of Galileo are Open Service (OS), Safety-of-Life Service (SoL),
Commercial Service (CS), Public Regulated Service (PRS) and Search and Rescue Service
(SAR) (ESA, 2014).
It should be noted, however, that the GNSS data used for this thesis is based on the current
performance of GNSS.
1.3 Related Research of GNSS Data Processing
With regards to the emergence and development of GNSS, theories and methods of GNSS da-
ta processing have been developed. Scientists and technologists made outstanding
contributions which improve the accuracy of GNSS precise positioning results. Differential (or
relative) and undifferential (or absolute) positioning are the main data processing methods in
GNSS precise positioning. They have specific advantages and disadvantages, differences and
common characteristics.
Double Differencing (DD) is a classical and early GNSS data processing method, which
can reduce common errors (such as receiver clock error, satellite clock error, satellite orbits
error, tropospheric delay, ionospheric effects and so on) by differential observations in a small
region. The DD method can provide a positioning service with an accuracy of centimetre or
even millimetre, which depends on the baseline length. This method has been widely used be-
cause it does not consider the complex model of error correction. In contrast, it has less
unknown parameters, higher accuracy of positioning results and it keeps the integer character-
istics of the DD integer ambiguities (Blewitt, 1989; Dong and Bock, 1989; Wei and Ge, 1998;
Chen, 1998; Hofmann-Wellenhof et al., 2008, p. 175; Wang et al., 2010; Shi et al., 2013). Real
Time Kinematic (RTK) positioning, Network-based Real Time Kinematic (NRTK) positioning
and Continuously Operating Reference Stations (CORS) can provide real-time kinematic posi-
tioning servers with an accuracy of centimetre or millimetre, based on one or more reference
stations or a reference network (Bock and Shimada, 1990; Wübbena et al., 1996; Han, 1997;
Chen et al., 2000; Chen et al., 2001; Grejner-Brzezinska et al., 2007; Snay and Soler, 2008).
In applications for short distances (about 10 km), a common practice is to neglect the
common error effects. However, in the case of longer distances (hundreds or even thousands
of km), differential common error residuals increase and may hamper the differential process,
or may decrease the accuracy (Genrich and Bock, 1992; Wielgosz et al., 2005; Li et al., 2010;
Li et al., 2013; Chu and Yang, 2013; Stephenson et al., 2011). Therefore, the reduction of dif-
ferential error effects is one of the most important steps to improve the relative kinematic
positioning. Hence, the medium and long-range precise kinematic positioning is studied here
in this thesis.
- 13 -
In the most recent decades, an innovative Precise Point Positioning (PPP) technology based
on undifferenced data processing strategies has been extensively studied. Based on the precise
satellite clock and orbit products of the International GNSS Service (IGS), the standard PPP
can provide positioning accuracy of decimetre to a few centimetre levels in post-processing
mode in a global reference frame (Zumberge et al., 1997; Kouba and Héroux, 2001; Dixon,
2006; Geng et al., 2010; Zhang et al., 2011b). Following RTK and NRTK technological revo-
lution, PPP is an innovative and flexible technology, without requiring users to set up their
own ground reference stations, providing unlimited operating distance at low cost. It can be
widely used in LEO satellite precise orbit determination, GNSS meteorology, investigation of
ionospheric effects, precise timing, GNSS seismology, earth plate movement and dynamics
research, engineering applications areas and many other earth science disciplines, with im-
portant application value (Chen et al., 2004; Zhang and Andersen, 2006; Bisnath and Gao,
2008).
In order to improve the accuracy of the standard PPP technology, the PPP Regional Aug-
mentation (PPP-RA) technology (very similar to the PPP-RTK technology) has been
developed to generate undifferential corrections in the observation domain from a regional
reference network which can be disseminated by reference stations and applied to user obser-
vations for an instantaneous ambiguity-fixing. In such a service, instantaneous ambiguity
resolution is accessible for regions with these observation corrections as regional augmenta-
tion information. At the user’s side, the data processing is similar to standard PPP and the
difference is whether the corrections based on the regional network are applied or not. PPP-RA
can provide services of positioning with a few centimetres accuracy after a convergence time
of about 30 minutes or after successful ambiguity fixing, and instantaneous positioning at a
few centimetres level in the regions where augmentation information is available (Collins et
al., 2008; Mervart et al., 2008; Teunissen et al., 2010; Zhang et al., 2011a; Geng et al., 2011;
Ge et al., 2012; Li et al., 2014).
In summary, to reach a higher precision in GNSS precise positioning for large-scale areas
or long baselines more information is needed, which comes from multiple reference stations or
a regional reference network. In this situation, undifferential and differential GNSS data pro-
cessing are essentially different implementations for the usage of the network data and are
integrated into a unified service (Li et al., 2013). In 2002 it was algebraically demonstrated for
the first time that the undifferential and differential algorithms of GNSS are equivalent (Xu,
2002; Xu, 2007, pp. 111, 112; Xu et al., 2010). However, this depends on the way how the sat-
ellite and receiver clocks are analysed. In undifferenced algorithm, the clocks can be estimated,
e.g. precise times corrections. Furthermore, most of these previous studies used the results of
- 14 -
DD processing taken as “true value”, to evaluate their results of the proposed approach. There-
fore, the method of DD processing is deeply studied in this research.
1.4 Challenges and Research Objectives
In the GEOHALO airborne gravimetric project, the new German scientific aircraft HALO was
used. Its characteristics are ultra-high-altitude, ultra-long-range, and high dynamics. The role
of GNSS precise positioning in airborne gravimtery and development of GNSS technology
imposed the following challenges and research objectives.
GNSS precise positioning is required to provide high-precision state information of the
kinematic platform in airborne gravimetry. Since here not only the position but also velocity
information for the platform are needed, the existing estimation method needs to be improved.
Furthermore, since the GNSS data from airborne gravimetric projects are of high frequency,
high dynamic, long range, and contain large amounts of data, the traditional estimation method
is not efficient. Therefore, a GNSS estimation method of high-precision and efficiency needs
to be developed. In this thesis, in order to obtain high-precision results for GNSS kinematic
precise positioning, the estimation methods of least squares including the elimination of nui-
sance parameters and a two-way Kalman filter are applied. Furthermore, accuracy evaluation
and reliability analysis of GNSS kinematic precise positioning is studied in-depth.
From the background of GNSS data processing, we know that both undifferential (PPP-
RTK and PPP-RA) and differential data processing are essentially different implementations
(utilizations) for the data of reference stations or a reference network. In order to obtain the
high-precision results of GNSS precise positioning, the information of far away reference sta-
tions or a large scale reference network should be considered. Since the differential method
shows high-precision results for GNSS precise positioning, by using only the existing IGS
tracking stations or establishment of several temporary reference stations in the survey area,
differential methods are thus widely used in airborne gravimetry. However, the survey area of
possible future airborne gravimetric projects might stretch over thousands of kilometres (like
closing the Antarctic polar gap in satellite-only gravity field models), so that the traditional
single baseline model is facing new challenges. If the length of a single baseline for the rela-
tive positioning mode is too long, the tropospheric and other effects cannot be completely
eliminated by the methods of model correction by applying DD between the kinematic and the
reference station. In addition, the amount of visible common satellites will be reduced with the
increasing length of the single baseline. Therefore, in order to solve this problem, a new meth-
od of GNSS precise positioning based on multiple reference stations is developed.
- 15 -
In airborne gravimetry, multiple antennas of GNSS receivers are usually mounted on the
kinematic platform. The geometric relations between the multiple GNSS antennas and the sim-
ilar characteristic of atmospheric effects within a small area above the kinematic platform
should be taken into account for the improvement of the accuracy and reliability in the posi-
tioning. In this thesis, a GNSS kinematic positioning method based on multiple kinematic
stations is proposed. Also, using the constant distance among multiple GNSS antennas, a kin-
ematic positioning method based on a priori distance constraints is proposed to improve the
reliability of the system. Furthermore, due to the similar characteristic of atmospheric effects
in a small area above the kinematic platform, a common tropospheric wet delay parameter for
the multiple GNSS antennas mounted on the platform is applied.
The trend of the development of GNSS is based on multiple systems and interoperability,
with an increased amount of visible satellites. This makes studying multiple GNSS systems
very essential. In addition, the reliability of a single GNSS system is lower. When there are not
enough observations or no observations, the kinematic positioning results are of a poor accu-
racy or fail. Therefore, a method of kinematic positioning using GNSS integration is employed
to improve the system reliability and accuracy. For this purpose, a GNSS integration algorithm
based on Helmert’s variance components estimation is proposed to reasonably adjust the
weights and the contributions of different GNSS systems.
The velocity information for the kinematic platform is an important state parameter in air-
borne gravimetry. By the methods of position differencing and carrier phase differencing, we
can only get the average velocity between observation epochs. However, a high dynamic and
high speed aircraft can be used in an airborne gravity project such as GEOHALO. The average
velocity (computed by position or carrier phase differencing) cannot meet the requirements of
airborne gravimetry in such a case. Therefore, the instantaneous velocity of the kinematic plat-
form determined by GNSS is much more important in airborne gravimetry on the HALO
aircraft. Thus, raw Doppler observations are used to determine the kinematic instantaneous
velocity in high dynamic environments. Furthermore, the carrier phase derived Doppler obser-
vations are used to obtain precise velocity results in low dynamic environments. Then, a new
method of Doppler velocity determination based on GNSS integration with Helmert’s variance
components estimation and robust estimation is developed to obtain more precise velocity in-
formation.
Finally, based on the algorithm proposed in this thesis, combined with the actual require-
ments of airborne gravimetry and shipborne gravimetry, a high-precision kinematic GNSS
position and velocity determination software (HALO_GNSS) has been developed for airborne
gravimetry.
- 16 -
1.5 Contributions of this Research
The main outcomes of the research done in this thesis are:
The parameter estimation algorithms of least squares adjustment including the elimi-
nation of nuisance parameter as well as the two-way Kalman filter were applied in
kinematic GNSS data post-processing to fulfill the high accuracy requirements of air-
borne gravimetry.
The specific accuracy evaluation methods for GNSS kinematic precise position deter-
mination have been developed and applied. In this context, the receiver clock errors
including clock jumps were analyzed for high dynamics GNSS precise positioning.
A new approach for the kinematic positioning using multiple reference stations based
on an a priori constraint is addressed to increase the accuracy and reliability of the
GNSS precise kinematic positioning in a large region. In this context, a robust estima-
tion theory is used to suppress the impact of observation outliers on the trajectory
estimates for airborne gravimetry.
A new method of GNSS kinematic positioning based on multiple kinematic stations
with an a priori distance constraint and a common tropospheric delay parameter is de-
veloped, which can enhance the accuracy and reliability of the state estimates as well.
A kinematic positioning method using multiple GNSS systems integration based on
Helmert’s variance component estimation is addressed to improve the reliability and
accuracy of GNSS kinematic positioning, and to adjust the weights in a reasonable
way whilst to enhance the contributions among multiple GNSS systems.
Velocity estimates were determined by the use of Doppler observations which were
obtained either from the receiver-generated raw Doppler or from carrier phase derived
Doppler data. A new approach for velocity determination based on the GNSS integra-
tion with Helmert’s variance component estimation and the robust estimation has been
developed.
A reliable and precise kinematic GNSS position and velocity determination software
system (HALO_GNSS) has been developed based on the algorithms proposed in this
thesis. This software can be applied in many fields and allows the GNSS kinematic
positioning with an accuracy of 1-2 cm and GNSS velocity determination with an ac-
curacy of about 1 cm/s using raw and about 1 mm/s using carrier phase derived
Doppler observations.
- 17 -
1.6 Overview of Dissertation
This thesis documents the approach in the development of algorithms, software, and analysis
methods for precise GNSS kinematic position and velocity determination for airborne gravim-
etry. It includes nine chapters and has the following outline.
Firstly, Chapter 1 presents the introduction, motivation, background, challenges and re-
search objectives of this thesis and specifies the contributions of this research.
Chapter 2 introduces the fundamentals of GNSS data processing and explains several dif-
ferent combinations of GNSS observations as used in this thesis for different purposes. The
main error sources of GNSS observations are given in this chapter as well.
Chapter 3 addresses the nuisance parameter elimination in the least squares adjustment as
well as the two-way Kalman filter for kinematic positioning. The methods of accuracy evalua-
tion for static and kinematic positioning results are described. The influence of receiver clock
errors in kinematic positioning is analyzed as well.
Chapter 4 explains the development of a method for multiple reference stations relative po-
sitioning based on a priori constraints. The application of robust estimation in the algorithm is
addressed as well.
Chapter 5 shows the improvements of the method of kinematic positioning based on multi-
ple kinematic stations as well as of the method based on a priori distance constraints among
multiple kinematic stations. The algorithm of using a common tropospheric parameter among
multiple kinematic stations is addressed too.
Chapter 6 highlights the precise positioning based on the integration of different GNSS sys-
tems as well as the Helmert’s variance component estimation used in this algorithm.
Chapter 7 gives the details of the velocity determination from the GNSS Doppler observa-
tions. The used Doppler data were from receiver-generated raw Doppler observations as well
as the carrier phase derived Doppler data. Furthermore, the Doppler velocity determination
based on GNSS integration and robust estimation is addressed.
Chapter 8 describes the characteristics and structure of the HALO_GNSS software based
on the algorithm of this thesis. The application of the kinematic GNSS algorithms and soft-
ware in the kinematic data analysis in an airborne gravimetric project and in shipborne
gravimetry is demonstrated.
Finally, Chapter 9 summarizes the main results as obtained in the previous chapters, pre-
sents the final conclusions and suggests recommendations for future work.
- 19 -
2 GNSS Observations and Error Sources
2.1 Introduction
The basic GNSS observables are the code pseudoranges and carrier phases as well as Doppler
measurements. The principle of the GNSS measurements and their mathematical expressions
are described below. In particular, several different linear combinations of GNSS observations
are used in this thesis for different purposes.
The key point for the GNSS precise positioning is an ability to mitigate all potential error
sources and disturbances in the system. All errors in GNSS observations caused by the space
segment, by the signal propagation, by the environment around the receiver, and by the hard-
ware of the receiver itself need to be mitigated. The mitigation can be carried out by modelling,
estimating, and forming single station observation combinations as well as by using differen-
tial techniques. These observation errors will be discussed in detail in this chapter.
2.2 GNSS Observations
The fundamental measurements recorded by a GNSS receiver are the differences in time or
phase between the signals transmitted by the GNSS satellites and reference signals generated
inside the receiver. The GNSS signals are transmitted at different frequencies (Chen, 1998).
GNSS receivers generate various observables from the GNSS signals. The observation types
of code pseudorange, carrier phase and Doppler are discussed in this section.
2.2.1 Pseudorange
Unlike the terrestrial electronic distance measurements, GNSS satellite navigation and posi-
tioning uses the “one-way concept” where satellite and receiver clocks are involved. Thus the
ranges are biased by satellite and receiver clock errors. Consequently, they are denoted as
pseudoranges (Hofmann-Wellenhof et al., 2008, p. 105). The GPS receiver generates a copy of
the pseudorandom code and compares it to that arriving from the satellite. A time offset is
computed by an autocorrelation function between the received pseudorandom code from the
satellite and that generated by the receiver. This offset contains the signal travel time and the
mis-synchronization of the satellite and receiver clocks. The pseudorange measurements typi-
cally have a precision on the order of 1-10 meters. The equation for the pseudorange
observable is given in (Xu, 2002; Seeber, 2003, p. 255; Deng, 2012)
- 20 -
, , , ,
()
s s s s s s s
r j r r r j r r j j r j
P c dt dt I T b b
, (2.1)
where the superscript
s
refers to a given satellite, subscript
r
refers to the receiver, subscript
j
identifies the frequency of the frequency-dependent terms,
,
s
rj
P
is the pseudorange between
the receiver
r
and the satellite
s
at the frequency
j
,
s
r
is the geometic distance from the
receiver
r
to the satellite
s
, including relativistic corrections as well as phase centre offset
and variations.
c
denotes the speed of light in vacuum,
r
dt
and
s
dt
are the offsets of the re-
ceiver and satellite clock with respect to system time,
,
s
rj
I
is the ionosphere delay on the
signal path at frequency
j
,
s
r
T
is the troposphere delay on the signal path,
,rj
b
and
s
j
b
are the
code bias of receiver
r
and satellite
s
at frequency
j
, respectively,
,
s
rj
represents the effect
of observation noise and all non-modeled error sources, such as errors in the satellite clock and
orbit prediction, inaccuracies in ionospheric and tropospheric modeling, multipath. Unit meter
(m) is used for all terms. In this thesis, the pseudorange measurement is mainly used to obtain
an approximate or initial position and to construct the ionosphere-free, geometry-free meas-
urements.
2.2.2 Carrier Phase
Carrier phase observations are obtained by comparing the phases between a signal transmitted
by a satellite and a similar signal generated by a receiver. The receiver records the fractional
phase of the GNSS satellite and keeps track of the changes of the received carrier phase. The
initial phase, called ambiguity, is unknown. In order to use phase observations the so-called
phase ambiguity must be resolved. The phase observations have a noise of a few millimetres
and are much more accurate than pseudoranges. The equation of the carrier phase measure-
ments can be written as (Chen, 1998; Leick, 2004, p. 256; Hofmann-Wellenhof et al., 2008, p.
106)
, , , , ,
()
s s s s s s s s s s
j r j r r r j r j r j r j j r j
c dt dt I T N d d e
, (2.2)
where
s
j
is the carrier wavelength of satellite
s
at frequency
j
,
,
s
rj
is the carrier phase
observation in cycles between the receiver
r
and the satellite
s
at the frequency
j
,
s
r
is
again the geometric distance from the receiver
r
to the satellite
s
including relativistic cor-
rections as well as phase centre offset and variations,
,
s
rj
I
is the ionosphere delay for carrier
phase observation on the path at frequency
j
scaled to units of length, it has the same magni-
tude as for pseudorange measurements, but is of the opposite sign,
,
s
rj
N
is the integer
- 21 -
ambiguity for a particular receiver-satellite pair at frequency
j
,
,rj
d
and
s
j
d
are the carrier
phase biases of receiver
r
and satellite
s
at frequency
j
, respectively,
,
s
rj
e
represents un-
modeled effects, modeling errors and measurement errors for carrier phase observations, it is
three or four orders of magnitude lower than that of code measurements.
2.2.3 Doppler
The Doppler effect is a phenomenon of frequency shift of the electromagnetic signal caused
by the relative motion of the emitter w.r.t. the receiver. In a first approximation, the Doppler
shift is given as (Hofmann-Wellenhof et al., 2008, p. 108; Xu, 2007, p. 41)
,,
ss
rr
s s s s
r j j r j j s
j
VV
D f f f
c
, (2.3)
where
,
s
rj
D
is the Doppler shift between satellite
s
and receiver
r
at frequency
j
,
s
j
f
de-
notes the emitted frequency
j
of satellite
s
,
,
s
rj
f
the received frequency
j
from satellite
s
,
s
r
V
is the relative velocity along the distance line between satellite
s
and receiver
r
,
s
j
is
carrier wavelength of satellite
s
at frequency
j
. The Eq. (2.3) for the observed Doppler shift
scaled to range rate is given by
,()
s
r
s s s s
j r j r r
V D c dt dt
, (2.4)
where the derivatives with respect to time are indicated by a dot,
is the measurement error.
The Doppler frequency shift is a by-product of the carrier phase measurements, an inde-
pendent observable and a measure of the instantaneous range rate. When a satellite is moving
toward the GNSS receiver, the Doppler shift is positive; so one gets more Doppler counts
when the range is diminishing.
In this thesis, the GNSS Doppler observations are used to estimate the velocity of the kin-
ematic platform and will be discussed in detail in Chapter 7.
2.3 Linear Combinations of GNSS Observations
Several linear combinations of the original GNSS carrier phase and code measurements are
used during data analyses to eliminate or reduce certain components of the observation equa-
tion. For instance, a linear combination can be formed to remove the effect of the ionosphere.
These combinations are listed below, followed by a brief description.
j
represents the phase
- 22 -
observations in cycles at frequency
j
,
j
P
denotes the code observations in meters at frequen-
cy
j
.
2.3.1 Ionosphere-free Linear Combination (LC)
The ionosphere delay caused by the refraction of the electromagnetic GNSS signal when
propagating through the ionospheric layer in the atmosphere ranges from 6 to 150 m. The
normal approach to eliminate the ionospheric delay is forming a dedicated linear combination
of GNSS observations. This combination is called the “ionosphere-free measurement LC”. For
the carrier phase observation Eq. (2.2) and the code observation Eq. (2.1), the ionosphere-free
combination can be written as (cf. Hofmann-Wellenhof et al., 2008, p. 127; Xu, 2007, p. 97)
2
1 1 2
12
2 2 2 2
1 2 1 2
LC
f f f
f f f f
(2.5)
and
22
12
12
2 2 2 2
1 2 1 2
LC
ff
P P P
f f f f
. (2.6)
The ionospheric delay is frequency-dependent. Eq. (2.5) and Eq. (2.6) eliminate the effect
of the first-order ionospheric delay on the observables which is widely used in GNSS data
processing. The disadvantage of this linear combination is that the noise from the
1
and
2
carrier phase measurements is increased almost by a factor of three (Seeber, 2003, p. 263), and
that the ambiguities cannot directly be solved as integer numbers.
For receivers that have dual-frequency capability the LC combination is usually the pre-
ferred method in geodetic and atmospheric applications for the estimation of coordinates,
tropospheric delay values and receiver clock biases.
2.3.2 Widelane Combination (WL)
The widelane observation is a popular linear combination mainly used for ambiguity and cycle
slip fixing which can be described as
12WL
. (2.7)
The combined wavelength
WL
of the L1 and L2 carrier phase measurements amounts to
86 cm. This long wavelength simplifies the ambiguity solution. It is commonly used in the
analysis for GNSS stations which are distant from each other by more than a few tens of km.
In the data pre-processing it can also be applied for detecting cycle slips (Blewitt, 1990).
- 23 -
2.3.3 Extra-widelane Linear Combination (EX_WL)
Since
1
and
2
carry the same geometric information, a position-independent quantity can
be constructed by subtracting the
2
carrier phase observation multiplied by the frequency
ratio from the
1
carrier phase observation as given by
1
_ 1 2
2
EX WL
f
f
. (2.8)
This extra-widelane linear combination eliminates the geometric (orbits, station coordi-
nates), tropospheric, and clock synchronization components of the carrier phase equation.
Thus it is often called geometry-free linear combination or extra-widelane. Since the combina-
tion of the initial phase ambiguities remains, EX_WL can only represent the complete
variation of the ionospheric delay during a continuous tracking (Chen, 1998).
For the GPS data pre-processing it is possible to construct a polynomial fit for EX_WL,
and to identify discontinuities such as cycle slips or outliers. But under high ionospheric activ-
ity conditions it is difficult to detect cycle slips with the ionosphere linear combination
(Blewitt, 1990).
2.3.4 MW Widelane Linear Combination (MW_WL)
The widelane observation in Eq. (2.7) still contains the position information. But the iono-
spheric effects and position information from the widelane observation can be eliminated, due
to the fact that the ionosphere effects on the code and phase measurements are equal but have
opposite signs in Eq. (2.1) and Eq. (2.2). When both, code and phase information are all avail-
able on two frequencies, the position-free and ionosphere-free observation can be given as
(Chen, 1998)
1 2 1 2
_ 1 2
1 2 1 2
MW WL
f f P P
ff
. (2.9)
This quantity is called the Melbourne-Wübbena combination (Melbourne, 1985; Wübbena,
1985). It combines the phase and code observations to eliminate the ionospheric, geometric
and clock effects and will be used for ambiguity initialization in GNSS data processing. Only
the error of multipath and the pseudorange noise are still contained in MW_WL, but they can
be reduced or eliminated by averaging multiple epoch (Blewitt et al., 1988).
- 24 -
2.4 Error Sources
The prerequisite of an accurate GNSS precise positioning is to reduce the errors resp. disturb-
ances contained in GNSS measurements, which include ionospheric effects, tropospheric
effects, relativistic effects, Earth tide and ocean loading tide effects, clock errors, antennas
phase centre corrections, multipath effects, antenna phase wind up as well as hardware biases.
The main errors will be discussed in this section. Further information can be found in Seeber
(2003), Leick (2004), Xu (2007) and Hofmann-Wellenhof et al. (2008).
2.4.1 Satellite Orbit and Clock Errors
The impact of satellite orbit and clock errors depends on the type of the applied GNSS pro-
cessing technique. Based on a simple rule of thumb (Beser and Parkinson, 1982; Seeber, 2003,
p. 302) the impact of the satellite orbit error on single point positioning can be given as
PDOPdp dr
, (2.10)
where
dp
is the point positioning error, PDOP is the Positional Dilution Of Precision
(Langley, 1999) and
dr
is the error of the satellite orbit in meter. For
PDOP 2
, a 2 m posi-
tion error will be caused by a 1 m orbit error.
For differential positioning, the impact of orbit errors on the solutions will grow larger with
the increase of the distance between the reference station and the kinematic GNSS stations
(Abdel-salam, 2005). It has been demonstrated that the relative accuracy of a baseline obtained
from DD GNSS observations can be related to the satellite position error by a rule of thumb
b
db dr
r
(2.11)
as given by Beser and Parkinson (1982), Parrot (1989), Dach et al. (2007) and Seeber (2003, p.
304), where
db
is the baseline error,
b
is the baseline length,
dr
is the orbit error,
r
denotes
the distance between the receiver and satellite (about 23000 km), According to this relation, to
achieve a baseline error below 1 cm, the admissible orbit errors as given in Table 2.1 are re-
quired for specified baseline lengths.
Table 2.1 clearly shows that, for differential positioning over short distances, the required
orbit accuracy is not a critical factor. However, the requirement for 1 cm accuracy over very
large distances, for example, in geodynamic applications over 500 km and more, implies an
orbit accuracy of better than 1 m. Nevertheless, the impact of satellite orbit errors on position-
ing is much larger for single point positioning than for differential positioning.
- 25 -
Table 2.1: The relation between admissible orbit error and baseline length for an 1 cm baseline error for
differential GNSS positioning
Baseline length
Required admissible orbit error
10 km
23 m
50 km
4.6 m
100 km
2.3 m
500 km
0.46 m
1000 km
0.23 m
1500 km
0.15 m
2000 km
0.11 m
In contrast to the orbit error, there is no baseline-dependent component for the satellite
clock correction (Parkinson, 1996). By using global GNSS tracking networks (shown as Fig-
ure 3.1) GNSS satellite orbit and clock parameter can be estimated with high accuracy. The
International GNSS Service (IGS), formerly the International GPS Service, provides orbits
and clocks in different latencies and accuracies as shown in Table 2.2.
Table 2.2: IGS combined orbit and clock products (IGS, 2013)
GNSS Products Type
Accuracy
Latency
Updates
Interval
GPS
Broadcast
orbits
~100 cm
real time
--
daily
Sat. clocks
~5 ns RMS
~2.5 ns SDev
Ultra-Rapid
(predicted half)
orbits
~5 cm
real time
at 03, 09, 15,
21 UTC
15 min
Sat. clocks
~3 ns RMS
~1.5 ns SDev
Ultra-Rapid
(observed half)
orbits
~3 cm
3 - 9 hours
at 03, 09, 15,
21 UTC
15 min
Sat. clocks
~150 ps RMS
~50 ps SDev
Rapid
orbits
~2.5 cm
17 - 41 hours
at 17 UTC daily
15 min
Sat. & Stn.
clocks
~75 ps RMS
~25 ps SDev
5 min
Final
orbits
~2.5 cm
12 - 18 days
every Thursday
15 min
Sat. & Stn.
clocks
~75 ps RMS
~20 ps SDev
Sat.: 30 s
Stn.: 5 min
GLONASS
Final
~3 cm
12 - 18 days
every Thursday
15 min
- 26 -
The IGS has presently 12 Analysis Centres (ACs). The IGS ACs use the collected data
from more than 360 stations worldwide to generate and provide precise GNSS ephemerides
and adjusted clock parameters. GFZ operates one of the IGS Analysis Centres since the very
beginning of the IGS activities (Deng, 2012). Compared to the Broadcast ephemeris, the
GNSS orbit and clock uncertainties can be significantly reduced by using the IGS ultra-rapid
(IGU), rapid (IGR) and final products.
2.4.2 Antenna Phase Centre Offset and Variation
The GNSS observation refers to the distance between the satellite and receiver antenna phase
centres, which varies with the frequency and the receiver-satellite orientation (Mader, 1999).
The need for satellite-related corrections is caused by the difference between the GNSS satel-
lite centre of mass and the phase centre of its antenna. Since the force models used for satellite
orbit modelling refer to the satellite centre of mass, the IGS GNSS precise satellite coordinates
and clock products also refer to the satellite centre of mass. In contrast the orbit ephemerides
in the GNSS broadcast navigation message refer to the satellite antenna phase centre. However,
all GNSS measurements are made to the antenna phase centre, thus one must know the satel-
lite phase centre offsets and has to monitor the orientation of the offset vector in space as the
satellite orbits the earth (Kouba, 2009a). Ignoring these phase centre variations can lead to se-
rious (up to 10 cm) vertical errors (Mader, 1999). Due to the differentiation between the
satellite’s centre of mass and the antenna phase centre position, the corrections for satellite
antenna Phase Centre Offsets (PCO) and Phase Centre Variations (PCV) must be considered in
GNSS data analysis for highly precise applications in particular for long baselines.
Analogously to the satellite, the phase centre of the GNSS receiver antenna is different
from its Antenna Reference Point (ARP) and should be taken into account (Rothacher, 2001).
The mean position of the electrical phase centre is determined for each frequency and the loca-
tion with respect to the ARP in a local reference system is denoted as PCO. However, the
actual electrical phase centre may depend on satellite-receiver direction. The deviation from
the mean phase centre is known as PCV (Moreno Monge, 2012).
Besides individual corrections for a specific antenna, there are also type-specific antenna
corrections (Schmid et al., 2007). These are mean values from calibrations of several antennas
of the same antenna type, which can then be used for all antennas of this type. The IGS cali-
brations of the GNSS satellite and station antennas can be downloaded from the IGS ftp site
(ftp://igscb.jpl.nasa.gov/pub/station/general/).
- 27 -
2.4.3 Satellite and Receiver Hardware Biases
The satellite as well as the receiver hardware systems cause time delays in the measured rang-
es, which differ from each other for every frequency and every modulating code. The
manufacturers spend a lot of work to avoid or at least to calibrate these biases. Yet, these bias-
es cannot be calibrated completely, since they depend on several influences (Roßbach, 2001).
Thus, they must be carefully estimated in precise GNSS applications. These biases are known
as hardware or electronic biases, including satellite code biases
s
j
b
and receiver code biases
,rj
b
in Eq. (2.1), satellite phase biases
s
j
d
and receiver phase biases
,rj
d
in Eq. (2.2)
(Banville et al., 2008).
These biases impose a different clock datum on each observable, so that the true clock er-
rors
r
dt
and
s
dt
cannot really be recovered (Collins et al., 2005). And it must be taken into
account that IGS precise satellite clock corrections are estimated from the ionosphere-free
combination of phase and code observations in P1 and P2, therefore they contain the iono-
sphere-free linear combination of the satellite hardware biases (Dach et al., 2007).
The receiver hardware delays are different for GPS and GLONASS observations. Whereas
for carrier phase observations the hardware delays are absorbed by the initial ambiguities, such
errors can reach several ns for the code observations. Nevertheless, when DD are formed, sat-
ellite and receiver hardware biases can be cancelled out. But in point positioning it’s difficult
to correct these errors caused by receiver code biases since they occur during the analysis of
doubly differential data from more than one receiver type, where their effect is not cancelled
out (Moreno Monge, 2012).
The hardware biases cannot be determined in an absolute sense, but they can be mostly
eliminated from the code observations by using the Differential Code Biases (DCB) correc-
tions, which are provided by the IGS (Collins et al., 2005; Banville et al., 2008). The DCB
corrections can be downloaded from the ftp site of the Centre for Orbit Determination in Eu-
rope (CODE) (ftp://ftp.unibe.ch/aiub/CODE/).
2.4.4 Tropospheric Effects
The troposphere is the lowest part of the atmosphere over the Earth’s surface. In contrast to the
ionosphere, the troposphere is a non-dispersive medium for the GPS carrier frequencies. In
other words, the tropospheric effects on the GPS signal transmission are independent from the
working frequency and the electromagnetic signals are affected by the neutral atoms and mol-
ecules in the troposphere. The effects are called tropospheric delay or tropospheric refraction.
(Xu, 2007, p. 55). The principle of the tropospheric delay is shown in Figure 2.1.
- 28 -
The troposphere causes a transmission delay of the GNSS signals both through signal path
bending and the alteration of the wave propagation velocity. The amount of tropospheric delay
in the zenith direction is about 2 m and it increases with the decrease of the elevation angle of
the sight line from the receiver to the satellite. For satellite elevation angles lower than 10 de-
grees, the tropospheric delay of the GNSS signal can reach up to about 20 meters (Xu, 2007, p.
55; Wang et al., 2009b). Therefore, the tropospheric effect is an important error source in pre-
cise GNSS applications.
The tropospheric delay can be split into the hydrostatic (i.e. dry) and the wet component.
The hydrostatic component is caused by the dry atmospheric gases and induces about 90% of
the total tropospheric delay. The hydrostatic component can be easily modelled (Schüler, 2001;
Abdel-salam, 2005). The wet component is mainly caused by the water vapour in the lower
part of the troposphere up to 11 km above sea level. The wet component causes 10% of the
total tropospheric delay (Abdel-salam, 2005; Hofmann-Wellenhof et al., 2008, p. 129). Due to
the variation of the water vapour density with position and time, modelling of the wet compo-
nent is difficult.
Figure 2.1: The principle of troposphere delay (Wang et al., 2009b)
The excess path length is the troposphere delay and is called Slant Total Delay (STD). The
hydrostatic and wet tropospheric delay are usually modelled as delays in zenith direction and
then mapped down to the satellite elevation using a mapping function. The STD between satel-
lite and receiver at an elevation angle
e
can thus be written in the form given below (Davis et
al., 1985; Collins and Langley, 1997)
STD ( ) ZHD ( ) ZWD
hw
M e M e
, (2.12)
- 29 -
where
STD
is the slant total tropospheric delay between the satellite and receiver,
ZHD
and
ZWD
are the zenith hydrostatic delay and the zenith wet delay, respectively.
()
h
Me
and
()
w
Me
are the hydrostatic and the wet mapping function, respectively, which depend on the
elevation angle
e
. The ZHD can be modelled with sub-millimetre accuracy using surface
pressure measurements under the assumption that the atmosphere is in hydrostatic equilibrium
(Janes et al., 1991). The ZHD can be retrieved with an accuracy of 1-3 mm using global Sea
Level Pressure (SLP) grids of atmospheric models as provided by the European Centre for
Medium-Range Weather Forecasts (ECMWF) together with an adequate model for the height
dependence of the atmospheric pressure (Fernandes et al., 2013; Singh et al., 2014). However,
ZWD cannot be modelled with an accuracy better than 2cm (Mendes and Langley, 1995;
Singh et al., 2014).
Here in the studies for this thesis different tropospheric correction models and mapping
functions based on theoretical approaches as well as on measured atmospheric data are used.
For the tropospheric path delay in zenith direction the models of Saastamoinen (Saastamoinen,
1972), Hopfield (Hopfield, 1969), Black-Eisner (Black and Eisner, 1984), Ifadis (Ifadis, 1986)
and Askne and Mordius (Askne and Nordius, 1987) are applied while the applied mapping
functions comprise Davis, Chao, Marini (Mendes and Langley, 1994), Niell mapping function
(Niell, 1996), Isobaric mapping function (Niell, 2001), the Vienna mapping function (Boehm
and Schuh, 2004) and a Global Mapping Function (GMF) (Bohm et al., 2006) etc.
Basically, three methods are used here for modelling the tropospheric delay. The first
method uses DD processing to eliminate the tropospheric error, but this works only for relative
positioning on short baselines. The second method uses a tropospheric delay model to revise
the error whereas the third method is an adjustment of troposphere parameters to estimate a
residual atmospheric delay. The pros and cons are discussed in Section 5.4.
2.4.5 Ionospheric Effects
The ionosphere is the higher stratum of the atmosphere with an extension from about 50-2000
km. The ionospheric refraction causes a group delay and a phase advance in the propagation of
the GNSS signals. The magnitude of the ionospheric delay or advance of the GNSS signal can
vary from a few meters to more than hundred meters within one day which is larger than the
delay caused by the troposphere (Xu, 2007, p. 48). The advance of the carrier phase and the
delay in the code observation are of the same magnitude but have opposite signs. However the
ionospheric refraction can be eliminated to a first order approximation by forming the men-
tioned ionospheric-free linear combination of dual-frequency (DF) GNSS observations. This
- 30 -
first order ionosphere-free combination can be formed as Eq. (2.5) for the carrier phase obser-
vation and Eq. (2.6) for the pseudorange observation.
Single frequency GNSS users can reduce the ionospheric delay by applying a model, such
as the Klobuchar model (Klobuchar, 1987) or the NeQuick model (Radicella and Leitinger,
2001), or use the values computed by the IGS Analysis Centre for either a global or a regional
ionosphere, such as Global Ionospheric Map (GIM) models (Jee et al., 2010).
2.4.6 Site Displacement Effects
In a global sense, a station undergoes periodic movements (real or apparent) reaching a few
dm that are not included in the International Terrestrial Reference Frame (ITRF) “regularized”
positions, from which “high-frequency” displacements have been removed using models.
Since most of the periodic station movements are nearly the same over broad areas of the
Earth, they almost cancel each other in differential positioning over short (<100km) baselines
and thus they do not need to be considered. However, if one is to obtain a precise station coor-
dinate solution consistent with the current ITRF conventions in PPP, using undifferential
approaches, or in relative positioning over long baselines (> 500 km), the above station
movements must be modelled as per recommendation in the IERS Conventions (Kouba and
Héroux, 2001; Kouba, 2009a).
2.4.6.1 Earth Tides
The Earth responds as an elastic body to external forces caused by the Sun and the Moon.
They cause periodic deformations of the Earth’s crust and lead to vertical and horizontal site
displacement, which can be computed quite accurately from simple Earth models, while the
ocean tides are strongly influenced by the coastal outlines and the shape of the near-coastal
ocean floor. The magnitude of the Earth tides is dependent on station latitude, tide frequency,
and sidereal time. The effect of the tidal variation is larger for the vertical component than for
the horizontal direction and can reach up to 30 cm. The horizontal movement can reach 5 cm
(Kouba and Héroux, 2001; Deng, 2012).
The site displacement effects of the Earth tide due to degree 2 of the tidal potential is
(McCarthy, 1996; Kouba and Héroux, 2001; Petit and Luzum, 2010)
4
32
22
22
3
2
33
22
0.025 sin cos sin
J
j j j
jJ
g
GM hh
r
r l R r R l R r r
GM R
mr
, (2.13)
- 31 -
where
r
is the site displacement vector in Cartesian coordinates
( , , )r x y z
.
GM
and
j
GM
are the gravitational parameters of the Earth, the Moon (
2j
) and the Sun
(
3j
);
,j
rR
are the geocentric state vectors of the station, the Moon and the Sun with the
corresponding unit vectors
r
and
j
R
, respectively;
2
l
and
2
h
are the nominal second degree
Love and Shida dimensionless numbers (about 0.608, 0.085);
,
are the site latitude and
longitude and
g
is Greenwich Mean Sidereal Time.
In the relative airborne kinematic GNSS positioning, the airborne antennas are of course
not affected by the Earth tide effects. However, the static reference antennas mounted on the
Earth surface are not free from tidal effects. In this case, the tidal displacements are usually
independent from the size of the applied area or lengths of the baselines and have to be taken
into account (Xu, 2007, p. 67).
2.4.6.2 Ocean Tide Loading
The ocean tide loading displacements affect mostly only the GNSS stations near the coast. The
displacements at most of the continental stations are less than 1cm. Loading correction is not
commonly considered in GNSS data processing because the computation is more complicated,
and its modelling is less accurate. However, for precise applications, loading effects have to be
taken into account (Xu, 2007, p. 72). When stations are located not too far away from the
nearest coast line (<1000 km), the ocean loading effects should be taken into account for ap-
plications such as point positioning with centimetre to millimetre accuracy or GPS
meteorology. Otherwise, this effect would be mapped into the ZTD and station clock correc-
tion (Kouba and Héroux, 2001; Kouba, 2009b).
The mode of ocean tide loading can be written as (McCarthy, 1996; Petit and Luzum, 2010)
cos( )
j j cj j j j cj
o f A t u
, (2.14)
where
o
is the displacement due to ocean loading,
j
represents the 11 tidal waves ( known
as
2 2 2 2 1 1 1 1
, , , , , , , , ,
fm
M S N K K O P Q M M
,
sa
S
), the
j
f
and
j
u
depend on the longitude of
the lunar node (at 1-3 mm precision
1
j
f
and
0
j
u
),
j
and
j
are the angular velocity
and the astronomical arguments at time
0t
h, corresponding to the tidal wave component
j
.
cj
A
is the station specific amplitude,
cj
is the station specific phase.
- 32 -
2.4.7 Ambiguities and Cycle Slips
As described in Section 2.2.2, the carrier phase measurement is made by shifting the receiver-
generated phase to track the received phase of the satellite signal. The absolute number of full
carrier wave oscillations between the receiver and the satellite cannot be counted at the initial
signal acquisition. Therefore, measuring the carrier phase means to measure a fractional phase
and to keep track of changes in the cycles. Due to the measurement principle in the receiver,
the carrier phase observable is an accumulated carrier phase observation (Xu, 2007, p. 39).
A full carrier wave oscillation is called a cycle. The ambiguous integer number of cycles in
the carrier phase measurement is called ambiguity. It is necessary to correct the initially meas-
ured fractional phase. This is done by setting up an arbitrary integer number at the start epoch
in the observations. Such an arbitrary initial setting will then be adjusted to the correct value
by modelling the ambiguity parameters (Xu, 2007, p. 39). Without an exact determination of
the integer ambiguities, no precise positioning at centimetre level based on GNSS phase ob-
servations can be achieved (Chen, 1998).
If the receiver loses the phase lock of the satellite signal, a cycle slip will occur. The rea-
sons for cycle slips may be obstructions, signal noise, low satellite elevation, weak signals,
antenna inclination in kinematic application (airplane, ship) and are caused by the signal pro-
cessing (Seeber, 2003, p. 277). Cycle slips have either to be removed from the data at the pre-
processing level, or a new ambiguity has to be determined (Seeber, 2003, p. 277).
The ambiguity estimation as well as the cycle slip detection and repair are not in the focus
of this study. More details can be read e.g. in Seeber (2003), Leick (2004), Hofmann-
Wellenhof et al. (2008) and Xu (2007).
2.5 Conclusions
In this chapter, the GNSS observations and their combination have been discussed. The main
error sources and disturbances of the GNSS observations are described as well. That compris-
es the fundamental theories and conceptions for the GNSS signal modelling in this study.
In comparison to pseudoranges the carrier phases are much more precise with a precision
of a few millimetres. For this reason the carrier phases are the most important observation type
for precise GNSS applications like GNSS positioning. The pseudorange observations are only
used to obtain an approximate or initial position and to construct the ionosphere-free, geome-
try-free measurements. The Doppler observations are used to estimate the velocity of the
kinematic platform.
- 33 -
The GNSS data analysis algorithms rely on signal combinations and differential position-
ing. In general, disturbances of the GNSS observations like tropospheric impacts, satellite
orbit and clock errors and receiver clock errors can be reduced or eliminated by differential
positioning processing. For long baselines, some of the errors, especially the tropospheric sig-
nal delay has to be taken into account in precise kinematic positioning.
The following chapters deal with the treatment of these observables and the mathematical
models to obtain a user position and velocity information from these measurements.
- 35 -
3 Algorithms Developing and Quality Analysis for GNSS
Kinematic Positioning
3.1 Introduction
In this chapter, the mostly used adjustment and filtering algorithms for static and kinematic
GNSS data processing are outlined. The adjustment algorithms as discussed here are the least
squares adjustment and an equivalent algorithm which is based on what is known as eliminat-
ed observation equation system. The filtering algorithms discussed here concern the classic
Kalman filter as well as the two-way Kalman filter.
The method used in this thesis for the evaluation of the precision and the reliability of
GNSS kinematic precise positioning is discussed. The influence of the GNSS receiver clock
error on the accuracy evaluation of the high-dynamics GNSS positioning is analyzed as well.
Finally, an advice about the GNSS receiver clock error estimation is given.
3.2 Nuisance Parameter Elimination in Least Squares
In the least squares adjustment, the unknown parameters can be divided into two groups. In
practice, sometimes only one group of unknowns is of interest. Therefore, it is often useful to
“eliminate” the other group of unknowns (called nuisance parameters) because of its size, for
example. In this case, using what is known as the equivalently eliminated observation equation
system could be very beneficial. That means the nuisance parameters can be eliminated direct-
ly from the observation equations instead of from the normal equations (Zhou, 1985; Xu, 2002;
Xu, 2007, p. 146; Shen and Xu, 2008).
In the least squares adjustment for GNSS precise positioning for moving platforms, the un-
known parameters can be divided into two groups as well. One group is the state parameters of
the kinematic platform, such as the time series of the position and velocity parameters. The
other group consists of global parameters, which are connected to all or some part of the data,
such as ambiguity parameters. Normally we are interested in the state parameters of the kine-
matic platform; however, the global parameters such as ambiguity parameters have a huge
impact on the state parameters. That means, if the global parameters are inaccurate, it is hard
to obtain precise solutions of the state parameters.
In order to obtain highly precise state results for kinematic platforms, the principle of the
elimination of nuisance parameters in least squares is used: Firstly, the state parameters are
- 36 -
eliminated for every epoch. Secondly, the global parameters are calculated using all data. Fi-
nally, by backward substitution of the estimated global parameters into the original equations,
the state parameters are calculated at every epoch. This method is based on the classic least
squares adjustment as follows.
3.2.1 Classic Least Squares Adjustment
The classic least squares principle can be summarised as follows (Gelb, 1974, p. 23; Mikhail
and Ackermann, 1976, pp.103-105; Koch, 1999, p. 188; Ghilani, 2010, pp. 178-182):
The GNSS observational model at epoch
i
is given as
i i i i
L A X e
, (3.1)
where
i
L
is an
1n
observation vector at epoch
i
,
i
A
is
nm
design matrix,
i
X
denotes a
1m
unknown parameter vector,
i
e
is the theoretical measurement error vector with zero
mean and covariance matrix
i
, with
21
0ii
P
, where
2
0
is the theoretical variance of
unit weight and
i
P
denotes the weight matrix of the observations.
If the
1n
vector
ˆi
L
denotes the estimator of the expected values
E( )
i
L
of the observa-
tions, then the estimator
i
V
(which is an
1n
vector) of the vector
i
e
of the errors in Eq. (3.1)
is given by
ˆ
i i i
V L L
, (3.2)
with
ˆ ˆ
i i i
L A X
, (3.3)
where
ˆi
X
is the estimator of the vector
i
X
of unknown parameters at epoch
i
. The vector
i
V
is called vector of the residuals, which plays an important role after the adjustment process. It
is basically useful to analyze the elements of
i
V
in order to test the adequacy of the model
(Mikhail and Ackermann, 1976, p. 104).
Following the well-known basic principle of the least squares adjustment (Mikhail and
Ackermann, 1976, p. 104; Koch, 1999, p. 187; Ghilani, 2010, p. 179),
minimum
i i i
V PV
, (3.4)
the solution of Eq. (3.2) and Eq. (3.3) is given by,
1
ˆ()
i i i i i i i
X A PA A PL
, (3.5)
- 37 -
with the covariance matrix
ˆi
X
of the estimated parameter vector
ˆi
X
as,
2
ˆ ˆ
0
ˆ
ii
XX
Q
, (3.6)
with
1
ˆ()
ii i i
X
Q A PA
, (3.7)
where
ˆi
X
Q
is the cofactor matrix. The empirical a posteriori variance of the unit weight
2
0
ˆ
can be computed from the residual vector
i
V
as follows
2
0
ˆ
ˆi i i i i i i i i i
V PV L PL L PA X
n m n m
, (3.8)
where
nm
.
3.2.2 Nuisance Parameter Elimination in Least Squares
Following the principles of the classic least squares adjustment (Gelb, 1974, p. 23; Mikhail
and Ackermann, 1976, pp.103-105; Koch, 1999, p. 188; Ghilani, 2010, pp. 178-182), the un-
known parameters can be divided into two groups, one with the state parameters (i.e. local
parameters), and the other one with the global parameters. The error equation of the GNSS
observational model at epoch
i
is given by
ˆ ˆ
i i i i i
V A X BY L
, (3.9)
where
i
V
is the residual vector of dimension
1n
, the vector
ˆi
X
denotes the
1m
estimated
state vector and
i
A
is the
nm
design matrix of the state vector at the epoch
i
, the vector
ˆ
Y
denotes the
1k
estimated vector of global parameters for all observation periods and
i
B
is
the
nk
design matrix of the global parameter vector at epoch
i
,
i
L
is the
1n
observation
vector with zero mean and the covariance matrix
i
, with
21
0ii
P
,
2
0
is the theoretical
variance of unit weight and
i
P
denotes the weight matrix of the observations.
The normal equation for the least squares adjustment derived from Eq. (3.9) is
ˆ
ˆ
T T T
i
i i i i i i i i i
T T T
i i i i i i i i i
X
A PA A PB A PL
B PA B PB B PL
Y
. (3.10)
- 38 -
Introducing
,,
,,
TT
XX i XY i i i i i i i
TT
YX i YY i i i i i i i
NN A PA A PB
NN B PA B PB
, and
,
,
T
Xi i i i
T
Yi i i i
UA PL
UB PL
, the Eq. (3.10)
can be written as
, , ,
, , ,
ˆ
ˆ
XX i XY i X i
i
YX i YY i Y i
N N U
X
N N U
Y
. (3.11)
Based on the theory of the equivalent elimination of nuisance parameters (Zhou, 1985;
Schaffrin and Grafarend, 1986; Xu, 2002; Schaffrin, 2004; Xu, 2007, p. 146; Shen and Xu,
2008), the estimated state parameter vector
ˆi
X
can be equivalently eliminated at epoch
i
, and
then, the equivalent normal equation of the global parameter vector
ˆ
Y
can be given as
11
, , , , , , , ,
ˆ
()
YY i YX i XX i XY i Y i YX i XX i X i
N N N N Y U N N U
. (3.12)
Introducing
iXYiXXiYXiYYiYY NNNNN ,
1
,,,,
and
iXiXXiYXiYiY UNNUU ,
1
,,,,
, the Eq.
(3.12) can be written as
,,
ˆ
YY i Y i
N Y U
. (3.13)
According to Xu (2002) and Xu (2007, p. 124), let’s define the auxiliary matrix
J
by
1
()
i i i i i i
J A A PA A P
. (3.14)
Hereby,
,( ) ( )( ) ( ) ( )
YY i i i i i i i i i i i i i i i
N B P I J B B P I J I J B B I J P I J B
(3.15)
and
,( ) ( )
Y i i i i i i i i i
U B P I J L B I J PL
. (3.16)
where
I
is the identity matrix. By introducing
()
i i i
D I J B
, Eq. (3.13) can be rewritten as
ˆ
i i i i i i
D PDY D PL
, (3.17)
then, Eq. (3.17) can be seen as the normal equation of the transformed error equation of
ˆ
Y
,
ˆ
i i i
W DY L
(3.18)
where,
i
W
is the respective residual vector having the meaning as
i
V
in Eq. (3.9).
i
L
is the
original observation vector with its associated weight matrix
i
P
.
Following the well-known basic principle of the classic least squares adjustment (cf. Sec-
tion 3.2.1), the solution of Eq. (3.13) is given by,
- 39 -
1
,,
ˆ()
YY i Y i
Y N U
, (3.19)
with the covariance matrix
ˆ
Y
of the estimated parameter vector
ˆ
Y
as,
21
ˆ,
ˆ()
Y YY i
YN
, (3.20)
Following Eq. (3.8) and considering Eq. (3.18), the empirical a posteriori variance of the unit
weight
2
ˆY
at epoch
i
can be computed from the residual vector
i
W
as follows
2ˆ
ˆi i i i i i i i i
Y
W PW L PL L PDY
n k n k
, (3.21)
where
nk
.
Here, the data of all observation periods contribute to the calculation of the global parame-
ter vector
ˆ
Y
. Eq. (3.13) can be summed from the first epoch to the last epoch
z
. The
equivalent normal equation for the global parameter vector
ˆ
Y
is given as
,,
11
ˆ
zz
YY i Y i
ii
N Y U
. (3.22)
The solution of the global parameter vector
ˆ
Y
and its covariance matrix
ˆ
Y
can be written
as
1
,,
11
ˆzz
YY i Y i
ii
Y N U
(3.23)
and
1
2
ˆ,
1
ˆ
z
Y YY i
Y
i
N
. (3.24)
Following Eq. (3.21), the empirical a posteriori variance of the unit weight
2
ˆY
of all period
can be computed as follows
21 1 1
11
ˆ
ˆ
z z z
i i i i i i i i i
i i i
Yzz
ii
ii
W PW L PL L PDY
n k n k
. (3.25)
After the global parameters have been calculated, taking them back to the original normal
equation Eq. (3.11), the state vectors
ˆi
X
( 1,2, , )iz
and their covariance matrices
i
X
ˆ
can be calculated as
- 40 -
1
, , ,
ˆ ˆ
()
i XX i X i XY i
X N U N Y
(3.26)
and
1 1 1
ˆ ˆ
, , , , ,
( ) ( )
iXX i XX i XY i XX i XY i
XY
N N N N N
. (3.27)
The formulas Eq. (3.9) to Eq. (3.27) are the complete formulas for the algorithm for the
elimination of nuisance parameters in the least squares adjustment.
The basic idea of this method is, all of the data are used to calculate global parameters first-
ly, then the local parameters are calculated based on the precise global parameters. The
solution is identical to that of solving Eq. (3.11) and Eq. (3.12) together (Xu, 2002). This elim-
ination process is the same as the Gauss-Jordan algorithm, which is often used for inverting
the normal matrix (or solving the linear equation system) (Xu, 2002).
The elimination of nuisance parameters in least squares adjustment is very useful for GNSS
kinematic data processing while the Integer Ambiguity Resolution (IAR) is a key approach in
high-precision GNSS applications. Without an exact determination of the integer ambiguities,
any precise positioning at centimetre level based on GNSS phase observations cannot be
achieved (Chen, 1998).
Using the described method for the elimination of nuisance parameters in the least squares
adjustment, the ambiguity parameters of the carrier phase observations can be calculated accu-
rately in a first step, where the ambiguity values are obtained as float solutions. These values
are then used to fix the integer ambiguities. Then the fixed ambiguities are taken to transform
the ambiguous carrier phases into ultra-precise pseudoranges which are eventually used for the
precise positioning. Several different methods have been developed to improve the ambiguity
estimation (Blewitt, 1989; Teunissen et al., 1996; Teunissen et al., 1999; Wang et al., 2001; Ge
et al., 2008; Teunissen et al., 2010). One of the most popular algorithms is the LAMBDAD
method (Teunissen et al., 1996) since this approach is numerically efficient and statistically
optimum. This method maximizes the probability of a correct integer ambiguity estimation
and its integer search is efficient.
3.3 Two-way Kalman Filter
The Kalman filter has been widely applied in many fields such as geodetic measurements
(Herring et al., 1990), GNSS kinematic precise positioning (Chen, 1998; Yang, 2010), satellite
orbit determination (He and Xu, 2011; Xu et al., 2012), GPS/INS integration (Mohamed and
Schwarz, 1999; Mercado et al., 2013) etc. Application of the Kalman filter has become stand-
ard for the estimation of linear state space models.
- 41 -
The classic Kalman filter will be described, and then the two-way Kalman filter will be
discussed.
3.3.1 Classic Kalman Filter
The system state equation and observation equation of GNSS kinematic positioning are gener-
ally expressed as (Schwarz et al., 1989; Koch and Yang, 1998; Yang, 2010)
. 1 1i i i i i
X X W
(3.28)
and
i i i i
L A X e
, (3.29)
where the subscript
i
is the time
i
t
,
i
X
and
1i
X
are
1m
state vectors at epoch
i
t
and
1i
t
,
respectively,
.1ii
is an
mm
transition matrix from state
1i
X
to
i
X
,
i
W
is an
1m
error
vector of system state model with zero mean and covariance matrix
i
W
,
i
L
is an
1n
meas-
urement vector at epoch
i
t
,
i
A
is an
mn
design matrix and
i
e
is a measurement error
vector with zero mean and covariance matrix
i
, here
21
0ii
P
,
2
0
is the variance of unit
weight and
i
P
denotes the weight matrix of observations.
The predicted state vector
i
X
and its covariance matrix
i
X
are denoted by
, 1 1
ˆ
i i i i
XX
(3.30)
and
1
ˆ
, 1 , 1 i
ii
T
i i i i W
XX
. (3.31)
The transition matrix
.1ii
of state vector
i
X
and the covariance matrix
i
W
of error vector
i
W
of Kalman filter system state model are given in the literature (Schwarz et al., 1989; Chen,
1998; Abdel-salam, 2005; Yang, 2010).
The error equations for the predicted state vector and the measurement vector are
ˆ
iii
X
V X X
(3.32)
and
ˆ
i i i i
V A X L
, (3.33)
where
i
X
V
and
i
V
denote the estimator of the vector
i
W
and
i
e
, respectively.
- 42 -
According to the classical Kalman filter theory (Kalman, 1960; Kalman and Bucy, 1961;
Huep, 1985), the state estimate at epoch
i
can be expressed as
ˆi i i i i i
X X K L A X
, (3.34)
with its covariance matrix
ˆi
X
,
ˆi
iii X
XI K A
(3.35)
and
1
ii
i i i i i
XX
K A A A
, (3.36)
where
I
is the identity matrix and
i
K
is the so-called Kalman gain matrix.
3.3.2 Two-way Kalman Filter
There are three types of applications of the Kalman filter for GNSS data analysis: smoothing,
filtering and predicting representing the process of providing solutions for past, current, and
future epochs, respectively (Chen, 1998). In this thesis, to improve the accuracy of the state
vector of the kinematic platform, the two-way Kalman filter is applied, where the filter is per-
formed in both forward and backward directions (He and Xu, 2009; He and Xu, 2011; Xu et
al., 2012). A smoothed solution for the two-way Kalman filter and its covariance matrix can be
given by
, , , ,
1 1 1 1 1
ˆ ˆ ˆ ˆ
,two , ,
ˆ ˆ ˆ
( ) ( )
i f i b i f i b
i i f i b
X X X X
X X X
(3.37)
and
, wo , ,
1 1 1
ˆ ˆ ˆ
()
i t i f i b
X X X
, (3.38)
where
,
ˆi two
X
and
,two
ˆi
X
are the smoothed solution and its covariance matrix,
,
ˆif
X
and
,
ˆif
X
are the forward solution and its covariance matrix,
,
ˆib
X
and
,
ˆib
X
are the backward solution
and its covariance matrix.
The two-way Kalman filter can be seen as a weighted average of the estimated states from
the forward and backward runs. By combining the forward and backward results, measure-
ments before and after a given epoch can contribute to the corresponding state estimate, which
improves the accuracy of the state vector of the kinematic platform. This method is also useful
for the estimation of accurate integer ambiguity solutions. It should be pointed out that the
proposed methods are suitable for post-processing of GNSS precise positioning.
- 43 -
3.4 Precision and Reliability Evaluation
The results obtained in the GNSS precise positioning always contain uncertainties, which in-
clude random errors, systematic errors, and outliers as well as coloured noise. There are many
different measures for evaluating the quality of the positioning results, such as precision, accu-
racy, reliability and uncertainty. Precision as expressed by the a posteriori obtained covariance
matrix of the coordinates, means the survey’s characteristics by propagating random errors
(Teunissen and Kleusberg, 1998). Accuracy is the closeness of the obtained positioning results
to the true values. Typically, an independent or alternative measurement technique is required
to assess accuracy. Reliability describes the survey’s ability to check for the presence of mod-
elling errors, refers to the ability to detect blunders and to estimate the effects that undetected
blunders may cause on a solution (Baarda, 1968; Leick, 2004, p. 161). Uncertainty can refer to
vagueness, ambiguity, or anything that is undetermined (Shi, 2010). Although uncertainty
comprises many different factors, this thesis will focus attention to a narrower notion of uncer-
tainty, namely that a (GNSS) survey is considered to be qualified for its purpose when it is
delivered with sufficient precision and accuracy. The terms precision and accuracy are often
used to describe how good the estimated position for a GNSS receiver is (Dixon, 1991).
A distinction must be made between accuracy and precision. Accuracy is the degree of
closeness of an estimate to its truth. Generally, it can be measured using root mean square of
errors (RMS). Precision is the degree of closeness of observations to their means. It can be
measured with standard deviation (SDev).
To evaluate the results of GNSS kinematic precise positioning, based on kinematic pro-
cessing, here in this study a static and a kinematic experiment were carried out as described in
the following sections.
3.4.1 Static Experiment
The GNSS data of a static experiment can be processed kinematically. When the true value of
the static station coordinate is unknown, we can compare the kinematic GNSS results with its
mean value. The precision is described by
2
1
1ˆ
SDev ( )
z
i
i
X
z
， (3.39)
where
SDev
is the precision of the GNSS results,
z
is the number of epochs of the GNSS
observations,
ˆi
X
is position and/or velocity of the static station at epoch
i
,
is the mean
value of
ˆi
X
, (i.e.
1
1ˆ
z
i
i
X
z
). But since precision is an internal statistical evaluation of the
- 44 -
kinematic GNSS results, it cannot show the systematic errors. Therefore an additional external
evaluation is needed.
If the true value of static station parameters is known as
true
X
, we can compare the GNSS
kinematic results with this true value. The more objective accuracy is described as
2
true
1
1ˆ
RMS ( )
z
i
i
XX
z
. (3.40)
Many GNSS tracking stations like those of the International GNSS Service (IGS) are wide-
ly distributed all over the world, shown in Figure 3.1. The highly accurate coordinate estimates
and the GNSS observation data of the IGS stations are available on the IGS website
(http://www.igs.org/). That means IGS tracking stations can be used for static experiments to
check the algorithms and software as applied here in this thesis.
Figure 3.1: The IGS tracking network (IGS, 2014)
3.4.2 Kinematic Experiment
For the kinematic experiment, an evaluation of the results is more difficult than for the static
experiment, since the true parameters of the kinematic station at every time epoch are difficult
to know. In order to check the algorithm and software, the following evaluation methods were
used.
Firstly, a special testing equipment or testing place has been constructed, where a GNSS
antenna can be moved over certain horizontal or vertical distances which are measured with
mm accuracy by rulers. In this thesis, the antenna movement experiment has been carried out
- 45 -
on the roof of building A17 on the Telegrafenberg campus in Potsdam, to check the algorithms
and the software. The experiment equipment is shown in Figure 3.2 for horizontal and in Fig-
ure 3.3 for vertical movement of an antenna. Both equipments have a fixed ruler which can be
used to record the relative position of the moved antenna. The relative positions measured with
the ruler were compared with the GNSS results to get an independent picture of the real accu-
racy of the GNSS positioning. In this thesis only the vertically moving test antenna was
applied (cf. Section 4.2.2).
Figure 3.2: Equipment for the horizontal movement experiment of a GNSS antenna
Secondly, the GNSS results of a short kinematic baseline are used to check the results of a
long kinematic baseline: The kinematic GNSS antenna in the aforementioned experiment was
used to form a short as well as a long baseline. Since the results of the DD kinematic GNSS
positioning for a short baseline of about 2 m are of millimetre accuracy they can be considered
as true position values of the kinematic antenna w.r.t. to those for a long baseline of about
1000 km. The objective, independent accuracy can be obtained by comparison of the GNSS
results for the short and long baselines. The principle can be mathematically written as
2
,long ,short
1
1ˆ ˆ
RMS ( )
z
ii
i
XX
z
, (3.41)
where
,long
ˆi
X
and
,short
ˆi
X
are the state vectors of the same kinematic station obtained for the
long and the short baseline at the same time, respectively.
- 46 -
Figure 3.3: Equipment for the vertical movement experiment of a GNSS antenna
Finally, the methods as described above are effective in the case of kinematic stations for
low dynamic and in a small area, since the position of the kinematic station can be estimated
accurately by other approaches such as short kinematic baseline. However, for a highly dy-
namic and long kinematic baseline situation such as on the HALO aircraft (Figure 1.2), it is
difficult to know the true value of the exact position of the moving station at every epoch.
Therefore, in order to investigate the accuracy of the kinematic results, two relative testing
methods are used.
One relative testing method is based on the fact that the distances between two antennas on
a moving kinematic platform are always constant. This fact is used to check the results of the
highly dynamic and long baseline kinematic station. Firstly, the positions of the kinematic sta-
tions are calculated. Then the length of the baseline between the two kinematic stations was
calculated. The accuracy of the kinematic results was then checked by inspection of the tem-
poral variations of the computed length of the baseline between the two kinematic antennas.
This principle can be mathematically expressed by
2
true
1
1
RMS ( )
z
i
i
dd
z
(3.42)
and
1 2 2 1 2 2 1 2 2
( ) ( ) ( )
i i i i i i i
d x x y y z z
, (3.43)
- 47 -
where
i
d
denotes the distance between the two kinematic antennas at epoch
i
,
true
d
is the true
distance between the two kinematic antennas, which can be measured independently before
the aircraft starts.
( , , )
i i i
x y z
are the coordinates of kinematic station at epoch
i
. If the true
distance is unknown, the mean value of distance is used for this evaluation. The principle can
be expressed as
2
mean
1
1
SDev ( )
z
i
i
dd
z
, (3.44)
where
mean
d
is the mean value over the whole measurement period of the distance between
two kinematic antennas.
Another relative testing method is based on the fact that the difference of positioning re-
sults coming from two GNSS receivers connected with the same GNSS antenna is
theoretically zero. This can also be used to check the results of the highly dynamic and long
baseline kinematic station. The principle can be expressed mathematically as
rec_1 rec_2 2
1
1ˆ ˆ
RMS ( )
z
ii
i
XX
z
, (3.45)
where
rec_1
ˆi
X
and
rec_2
ˆi
X
are the position vectors of same antenna coming from two receivers.
3.4.3 Reliability
For use in geodetic networks, Baarda (1968) developed his testing procedure that ultimately
led to a theory of “reliability” with a number of applications (Baarda, 1976; Förstner, 1983; Li,
1985; Förstner, 1987; Li, 1989; Biacs et al., 1990; Schaffrin, 1997; Hewitson and Wang, 2006;
Knight et al., 2010). The reliability of GNSS is essentially dependent on the redundancy and
geometry of the measurement system as well. Reliability refers to the consistency of the re-
sults provided by a system, dictating the extent to which they can be trusted, or relied upon.
More specifically, in terms of GNSS kinematic positioning, reliability comprises the ability of
the system to detect outliers, referred to as internal reliability, and a measure of the influence
of undetectable outliers on the parameter estimations, referred to as external reliability (Baarda,
1968).
The measure of internal reliability is quantified as the minimal detectable bias (MDB) and
is indicated by the lower bound for detectable outliers. The MDB is the magnitude of the
smallest bias that can be detected for a specific level of confidence and is determined, for cor-
related measurements (Baarda, 1968; Förstner, 1987) by
- 48 -
0
0i
il
i
lr
, (3.46)
where
0i
l
is the lower bound for detectable outliers, the standard deviation
i
l
is the
th
i
(
1,2, ,in
) component of precision of the GNSS observations. The lower bound
0
is the
non-centrality parameter which depends on the given false alarm rate
0
and the delectability
0
. If the significance level is
0
1 99.9%
, then the gross errors
i
l
smaller than
4.13 i
l
are not detectable with a probability larger than
080%
(Förstner, 1987). The
redundancy number
i
r
of the
th
i
component of the GNSS observations is defined (Baarda,
1968; Förstner, 1983; Li, 1985) by
1
( ( ) )
i ii
r I A A PA A P
, (3.47)
where I is the identity matrix, the definitions of matrices A and P can be found in Section 3.2.1.
Hence, the smaller the redundancy number
i
r
of the observation, the larger a gross error has to
be in order to be detectable. This is reasonable, because the residuals are smaller in this case
(Förstner, 1987).
A comparison of the boundary values for heterogeneous observations is not directly possi-
ble. Thus, the controllability of the observation was written as (Li, 1985, p. 20; Förstner, 1987)
00
0,
i
i
i
li
l
r
, (3.48)
Therefore,
0,i
is a factor for the standard deviation
i
l
of the observation needed to obtain
the boundary value
0i
l
. With the controllability factor
0,i
we finally get
0 0, i
i i l
l
. (3.49)
The controllability value
0,i
is thus the factor by which a gross error at least has to be at least
larger than the standard deviation of the observation to be detectable with a probability of at
least
0
(Förstner, 1987).
External reliability of the system is characterized by the extent to which an MDB affects
the estimated parameters (Baarda, 1968). Gross errors smaller than the boundary values
0i
l
may stay undetected and contaminate the result. The maximum influence of undetectable gross
errors on to the estimates can be determined by sensitivity or robustness factor
0,i
(Li, 1985,
p. 21; Förstner, 1987; Leick, 2004, p. 169),
- 49 -
0, 0 0
1
ii
i
ii
ur
rr
, (3.50)
where,
1
1 ( ( ) )
i i ii
u r A A PA A P
, the meaning of matrices A and P is same as in Eq.
(3.47). These sensitivity or robustness factors
0,i
measure the external reliability according to
Baarda (1968). There is one such value for each observation. If the
0,i
are same order of
magnitude, the network is homogeneous with respect to external reliability. If
i
r
is small, the
external reliability factor becomes large and the global falsification caused by a blunder can be
significant. It follows that small redundancy numbers are not desirable (Leick, 2004, p. 169).
3.5 Influence of the GNSS Receiver Clock in Accuracy Evaluation
3.5.1 Receiver Clock Jump
In the application of GNSS precise positioning, most of the GNSS receivers attempt to keep
their internal clocks synchronized to the system time, such as GPS time. This is done by peri-
odically adjusting the clock by applying time jumps. The actual mechanisms of receiver clock
jumps are typically proprietary (Kim and Langley, 2001).
In practice, the clock jumps are largely divided into two categories. The first method is a
more or less continuous steering, where the clock drift is tuned approximately to zero, and the
offset is kept constant within the level of the noise and the tracking jitter. The other method is
an application of millisecond jumps which keeps the receiver clock time synchronized within
1 ms w.r.t. the GNSS system time (Kim and Lee, 2012; Zhang et al., 2013).
In order to illustrate the performance of the clock jumps, three GNSS receivers have been
tested. Figure 3.4 shows for the receiver 1 (NOVATEL OEM4) how the clock of this receiver
is continuous steered. The clocks of the other two receivers 2 and 3 (JAVAD DELTA G3T) are
millisecond jumping, which is shown in the Figure 3.5 and Figure 3.6, respectively.
- 50 -
Figure 3.4: The continuous steering of GNSS receiver 1 (NOVATEL OEM4)
Figure 3.5: The millisecond jumping of GNSS receiver 2 (JAVAD DELTA G3T)
- 51 -
Figure 3.6: The millisecond jumping of GNSS receiver 3 (JAVAD DELTA G3T)
3.5.2 Influence of Receiver Clock in Highly Dynamic Positioning
Due to the influence of GNSS receiver clock errors, which include millisecond jumps, the fi-
nal result time of the GNSS results at epoch
i
is different from the time tag which indicates
the estimated time instant when the measurements are sampled by the receiver. The receiver
clock error must be removed from the time tag of epoch
i
by
result tag
t t t
, (3.51)
where
result
t
is the final result time,
tag
t
is the time tag,
t
is the receiver clock error including
millisecond jumps.
But in some commercial software, the difference between the final results for the time and
the time tag are not considered. Here, usually the time tag is used in the final results. That is
not a problem for static and low dynamics positioning when the velocity is low or nought. But
in the case of the GEOHALO airborne gravimetric project, the positions of the GNSS anten-
nas at the final time epochs
result
t
are different from the positions at the time tags
tag
t
, since the
HALO aircraft has a high velocity of about 450 km/h (corresponding to 125 m/s).
In order to illustrate this issue, highly dynamic GNSS data, coming from the measurement
day June 06, 2012 of the GEOHALO project are selected for a test. Two mechanically fixed
antennas mounted on the HALO aircraft and connected with JAVAD DELTA G3T receivers 2
and 3 were taken and the distance between both antennas was estimated. In this context two
- 52 -
calculation schemata were performed and compared (both receivers were the same as those
used for the test as described in Figure 3.5 and Figure 3.6).
Scheme 1: The final results for the distance between both antennas were calculated at the
same time tags, without considering the clock errors, like in some commercial software. The
result is shown in Figure 3.7.
Figure 3.7: Variation of the estimated distance between two mechanically fixed antennas, computed at
the same time tags
Scheme 2: The final results for the distance between both antennas were calculated at the
epochs of the final results for the receiver time. If the final times of both receivers were differ-
ent, the results have been interpolated to the same time epoch. The result is shown in Figure
3.8.
The above results show that the time tags should not be used for the timing of the final
GNSS results. The described test of distance estimation between two mechanically fixed an-
tennas to evaluate the accuracy of GNSS kinematic positioning shows that the antenna
positions must be interpolated to the same time epochs.
- 53 -
Figure 3.8: Variation of the estimated distance between two mechanically fixed antennas, computed at
interpolated time epochs
3.6 Conclusions
In this chapter 3, the parameter estimation algorithms as used for this thesis and the accuracy
evaluation method were described. Furthermore, GNSS receiver clock errors and their impact
were exemplified.
Firstly, it has been outlined that the algorithm for the elimination of nuisance parameters in
least squares adjustment shall be used for the GNSS kinematic precise positioning. This ad-
justment approach is based on the classic least squares adjustment. The idea of this method is
to use the whole measurement data for the calculation of the global parameters in the first step.
Then the local parameters are calculated based on the already estimated global parameters.
Secondly, based on the classic Kalman filter, the two-way Kalman filter was proposed for
GNSS kinematic precise positioning. The two-way Kalman filter can be regarded as a
weighted average of the estimated states from the forward and backward runs. By combining
the forward and backward results, measurements before and after a given epoch contribute to
the corresponding state estimate, which improves the accuracy of the state vector of the kine-
matic platform.
It should be pointed out that the outlined methods, i.e. least squares including the elimina-
tion of nuisance parameters and the two-way Kalman filter, are suitable to fulfill the high
accuracy requirements for airborne gravimetry in the post-processing of GNSS precise posi-
tioning. That also means that the outlined methods are capable for an accurate integer
ambiguity solution. After the global parameter estimation and an application of the forward
- 54 -
Kalman filter, a precise ambiguity float solution can be obtained, which is the basis for the
final integer ambiguity estimation. Without a precise integer ambiguity determination, precise
positioning on the centimetre level using GNSS phase observations cannot be achieved.
The explained accuracy evaluation method for GNSS kinematic precise positioning has
been applied to data sets recorded in static and kinematic modes, for short as well as long
baselines, and at low and high dynamic states. In this thesis, one or more evaluation methods
have been used for each example depending on the characteristics of the recorded data.
Finally, the receiver clock errors including the clock jumping were analyzed for highly dy-
namic GNSS precise positioning. It has been shown that the receiver clock error including
clock jumps must be removed from the time tag. If a kinematic estimation for the distance be-
tween two mechanically fixed kinematic antennas is used to evaluate the accuracy of GNSS
kinematic precise positioning, the antenna positions must be computed for the same time
epochs.
- 55 -
4 Kinematic Positioning Based on Multiple Reference
Stations
4.1 Introduction
This chapter deals primarily with a GNSS Double Difference (DD) positioning approach
based on multiple reference stations, which are used to determine the baseline vectors between
simultaneously observing receivers.
This chapter begins with the discussion of GNSS DD positioning. The observation equa-
tions for Single Difference (SD) and DD are given in detail. The challenge of DD positioning
in a wide region will be presented by means of ultra-short and ultra-long baseline experiments.
Furthermore, an approach of multiple reference stations kinematic positioning based on an
a priori constraint is addressed for GNSS precise kinematic positioning in a wide region.
Finally, the robust Kalman filter theory is used to suppress the impact of observation outli-
ers on the trajectory estimates for airborne gravimetry.
4.2 Double Difference Positioning
DD positioning is a classic GNSS data processing method (Blewitt et al., 1988; Blewitt, 1989;
Dong and Bock, 1989; Chen, 1998; Leick, 2004, p. 261; Hofmann-Wellenhof et al., 2008, pp.
169-191; Xu, 2007, p. 107). It is regularly used to eliminate or reduce errors between the satel-
lite and the receiver (Bossler et al., 1980; Chu and Yang, 2013). Since the carrier phase
observations to GNSS satellites can be used to determine user positions even more precisely
than by means of the pseudorange measurements (Roßbach, 2001) the principle of the carrier
phase DD positioning is given in this section.
4.2.1 Principle of Classic Double Difference Positioning
A simple example for the DD positioning is illustrated in Figure 4.1.
p
and
q
denote the
GNSS satellites,
r
denotes the reference station and
k
denotes the kinematic station. Accord-
ing to the carrier phase observation Eq. (2.2), the carrier phase observations of the reference
station
r
and the kinematic station
k
to a common satellite
p
can be written as
, , , , ,
()
p p p p p p p p p p
j r j r r r j r j r j r j j r j
c dt dt I T N d d e
(4.1)
- 56 -
and
, , , , ,
()
p p p p p p p p p p
j k j k k k j k j k j k j j k j
c dt dt I T N d d e
. (4.2)
Figure 4.1: An illustration for GNSS DD positioning (Stombaugh and Clement, 1999)
Forming a single difference (Roßbach, 2001; Hofmann-Wellenhof et al., 2008, p. 174), i.e.
subtracting the measurement at the reference station from that at the kinematic station, the sat-
ellite clock
p
dt
and the satellite phase bias
p
j
d
(which are common errors to both observers)
are cancelled out by
, , , , ,
p p p p p p p p
j k r j k r k r k r j k r j k r j k r j k r j
c dt I T N d e
,(4.3)
where the operator
indicates the single difference, for example
,
p
k r j
denotes
,,
pp
k j r j
.
If both antennas are closely located to each other, for example less than 10 km on ground,
the tropospheric delays are approximately the same, so that the difference of the tropospheric
delays
p
kr
T
are small and considerably reduced in Eq. (4.3). But for medium and long base-
lines, the signal travel path through the troposphere is different for each station, especially if
the receivers are placed at different altitudes. This may be the case, e.g., in mountainous re-
gions or when an aircraft is approaching an airport (Roßbach, 2001). Thus the tropospheric
path delay is not cancelled out in single difference positioning for medium and long baselines.
Supposed the kinematic station is sufficiently close to the reference station, the path of the
GNSS satellite signal through the ionosphere will be almost identical for the reference station
and the kinematic station. Hence the ionospheric delay can be considerably reduced. This as-
p
q
r
k
- 57 -
sumption can be made for distances of up to approximately 1000 km (Roßbach, 2001; Dai,
2013). For long baselines, the ionospheric delay can be eliminated to a first order approxima-
tion by forming the ionospheric-free linear combination of dual-frequency observations (cf.
Section 2.3.1). Thus the differential ionospheric delay
,
p
k r j
I
can be neglected. In this case,
Eq. (4.3) can be further simplified to
, , , ,
p p p p p p p
j k r j k r k r k r j k r j k r j k r j
c dt T N d e
. (4.4)
Such a differential processing is called “single difference positioning”, namely differencing
the measurements of two receivers with respect to a common satellite. The single difference
positioning is applied to the satellite
q
with the signal frequency
i
, which is given as
, , , ,
q q q q q q q
i k r i k r k r k r i k r i k r i k r i
c dt T N d e
. (4.5)
Assuming the same frequencies
ji
(i.e.
pq
ji
) for the satellite signals and using the
denotation
p q p q
k r k r k r
, differencing two single difference observations, i.e. Eq.
(4.4) and Eq. (4.5), at the same sites
r
and
k
to two different GNSS satellites
p
and
q
gives
the DD observation as
, , ,
p p q p q p q p p q p q
j k r j k r k r j k r j k r j
T N e
. (4.6)
The receiver clock error
kr
c dt
and the carrier phase bias
,k r j
d
are common error
terms in Eq. (4.4) and Eq. (4.5) and hence can be reduced by DD positioning. That is the main
reason why DD are preferably used.
However, if the signal frequencies between the GNSS satellites are different (i.e.
pq
ji
),
as in the case of the GLONASS system, the carrier phase bias
,k r j
d
and
,k r i
d
cannot be
cancelled out from Eq. (4.6). Introducing a shorthand notation, symbolically
, , ,
p q p q
k r ji k r j k r i
, the DD equation for two different frequency satellites can be
given as
, , , ,
,,
p p q q p q p q p p q q
j k r j i k r i k r k r j k r j i k r i
pq
k r ji k r ji
T N N
de
. (4.7)
It is difficult to separate the DD carrier phase bias
,k r ji
d
from the ambiguity parameter,
thus the float solution for the GLONASS ambiguities is applied.
The procedure described above is called “DD positioning”, where the satellite orbit error,
the satellite clock error and the receiver clock error have been reduced. In case of short base-
- 58 -
lines, the remaining ionospheric and tropospheric delay can be neglected. For medium and
long baselines, the ionosphere-free combination (LC) is used to eliminate the first order iono-
sphere path delays. The remaining tropospheric delay is reduced by the estimation of the
zenith wet delay parameters.
Depending on the principle of the DD processing and the two-way Kalman filter, an exam-
ple is given in the following.
4.2.2 Experiment and Analysis
For testing the capability of the DD processing for short and long baselines, the mentioned
vertical antenna movement experiment (cf. Section 3.4.2) was carried out on the roof of build-
ing A17 at GFZ on March 27, 2012. In this experiment, the IGS station POTS was taken as the
reference station for the short baseline. Its antenna and receiver type are
JAV_RINGANT_G3T NONE and JAVAD TRE_G3TH DELTA, respectively. The experi-
mental station named KIN1 is treated as the moving kinematic station, which is located close
to the IGS reference station POTS with a distance of about 2 m. The station KIN1 is shown in
Figure 3.3. Antenna and receiver type of KIN1 are the same as those of the IGS station POTS.
The experiment was done in the following way: the antenna of KIN1 was kept at the start-
ing point for 18 minutes. Then the antenna was uplifted by 46 cm and kept at this higher
position for 10 minutes. Then the antenna was lowered down to the position of its starting
point. After this, the antenna was uplifted again step by step to heights of 10 cm, 20 cm and 40
cm. From this position the antenna was lowered down to the starting point in two steps of 20
cm each. At the end, the antenna was quickly uplifted and lowered between the largest height
of 46 cm and the starting position. The north, east and up components of the estimated trajec-
tory for the antenna KIN1 w.r.t. POTS (i.e. the short baseline) are shown in Figure 4.2.
The sampling interval for the GNSS observation data was 1 second. The precision of the
pseudorange and carrier phase observations is given to be 1 m and 3 mm, respectively. An ele-
vation cut-off angle of 10 degrees was applied in this test. About 1 hour of GPS and
GLONASS data were processed by the HALO_GNSS software which was developed by the
author for the GEOHALO project at GFZ. The details about the GNSS integration and the
HALO_GNSS software are discussed in Chapter 6 and Chapter 8, respectively.
The results of the kinematic positioning for KIN1 are compared with its initial position on
the starting point. This is shown in Figure 4.2 for the north, east and up directions w.r.t. POTS
(i.e. the short baseline). The RMS of the north and east directions is 2-3 mm, since the antenna
was not moved in the horizontal direction (quadratic mean of north and east). The movement
results of the up direction agree well with the true movement log at 3 mm accuracy. Such a
- 59 -
high accuracy results is achieved since the common errors between the satellites and receivers
are clearly eliminated or reduced by the DD processing in the case of the short baseline. These
results can be treated as the true value of the antenna movement.
Figure 4.2: The results of the kinematic positioning for the short baseline POTS-KIN1 during the
vertical antenna movement experiment
In order to show the performance of the DD positioning in the case of a long baseline, an
IGS tracking station far away from Potsdam named MATE (Matera, Italy) was chosen and
taken as reference station to form the long baseline together with KIN1 (cf. Figure 4.3). The
antenna and receiver type of the IGS station MATE are LEIAT504GG and LEICA
GRX1200GGPRO, respectively. The length of this long baseline MATE-KIN1 is 1330 km.
The GNSS observation data of the IGS station MATE on March 27, 2012, were download-
ed from the IGS ftp-server (ftp://cddis.gsfc.nasa.gov/gps/data/highrate/) with 1 s sampling rate.
This ftp-server is operated by the Crustal Dynamics Data Information System (CDDIS) and
provides also GNSS ephemeris data, orbit information and meteorological parameters of the
observation stations. The precision of the pseudorange and carrier phase observations is given
to be 1 m and 3 mm respectively. An elevation cut-off angle of 10 degrees was applied in this
test. The IGS rapid or final products of GNSS satellite orbits and clocks were chosen to apply
corrections for satellite orbit and clock errors. The HALO_GNSS software was used for the
GNSS data processing. The dual-frequency carrier phase observations were used to form the
ionosphere-free linear combination. The parameters of the wet tropospheric zenith path delay
- 60 -
for MATE and POTS were employed as well and they are assumed to be known with an un-
certainty of 10 cm and a spectral density of
92
10 m /s
. The two-way Kalman filter was used
for the parameter estimation.
Figure 4.3: The long baseline (MATE-KIN1, 1330 km) (IGS, 2014)
The results of the kinematic DD positioning for KIN1 on the long baseline MATE-KIN1
are shown in Figure 4.4 with its north, east and up components. It shows that there are biases
in the trajectory of the kinematic station KIN1.
The results of KIN1 for the long baseline MATE-KIN1 are compared with the “true values”,
of the short baseline POTS-KIN1. The statistics of the comparison given in Table 4.1 shows
that the achieved accuracy is 0.012 m, 0.022 m and 0.068 m for the north, east and up direc-
tions, respectively.
The reason for the lower positioning accuracy of the long baseline compared to the short
baseline is the fact that errors, such as the tropospheric delay, cannot be eliminated in long
baselines. Thus, to improve the accuracy and reliability of the GNSS kinematic positioning in
large regions, a multiple reference stations approach should be taken into account.
KIN1
- 61 -
Figure 4.4: The results of the kinematic positioning for the long baseline MATE-KIN1 during the
vertical antenna movement experiment
Table 4.1: The statistical results of the difference between the long and the short baseline (Unit: m)
Direction
Min
Max
Mean
RMS
North
-0.026
0.004
-0.012
0.012
East
0.009
0.034
0.022
0.022
Up
0.011
0.130
0.062
0.068
4.3 Multiple Reference Stations Kinematic Positioning Based on an A pri-
ori Constraint
The use of a single baseline for kinematic GNSS positioning can provide high positioning ac-
curacies in small areas. But for the ultra-long-range kinematic positioning, the use of multiple
reference stations should be considered (Dai et al., 2001; Lachapelle and Alves, 2009; Wang et
al., 2011; Al-Shaery et al., 2011; Firuzabadì and King, 2012).
Compared to the single baseline users, one of most important advantages of multiple refer-
ence stations is the increase in reliability and availability. If one or two reference stations fail
at the same time, their missed contributions can be substituted by the remaining reference sta-
tions, thus the availability of the service is retained (Fotopoulos and Cannon, 2001). Another
- 62 -
quite important aspect is that it allows to model distance-dependent or spatially correlated er-
rors, such as satellite orbit, tropospheric and ionospheric effects. This effect of distance-
dependent errors can be reduced by combining observations from the multiple reference sta-
tions (Fotopoulos and Cannon, 2001).
4.3.1 Principle of Multiple Reference Stations Kinematic Positioning Based on an A
priori Constraint
Airborne gravimetric projects, such as the GEOHALO project, are usually carried out over
large areas, thus a dedicated approach of the Kalman filtering in kinematic GNSS positioning
based on an a priori constraint for multiple reference stations was developed as described in
this section to improve the accuracy and reliability of the kinematic GNSS positioning. The
idea is the following: When the DD observation equations are formed, there is just only one
station which is used as a formal reference station. The other reference stations are processed
formally as kinematic stations together with the real kinematic stations. But in contrast to the
real kinematic stations an a priori constraint based on the known station information is applied
to the other reference stations. If they do not have those constraints, they would be real kine-
matic stations.
In the Kalman filter DD positioning of a single baseline between the kinematic station
k
and the first reference station
1
r
, the system state equation Eq. (3.28) and observation equation
Eq. (3.29) can be re-written as
, 1 1
k k k k
i i i i i
X X W
(4.8)
and
k k k k
i i i i
L A X e
, (4.9)
where the subscript
i
refers to the time epoch
i
t
. The superscript
k
denotes the kinematic
station.
i
k
X
and
1
k
i
X
are
1m
unknown parameter vectors at time
i
t
and
1i
t
, respectively,
where the state vector of the kinematic station
,,
i i i
x y z
or is included
, , , , ,
i i i i i i
x y z x y z
.
T
is the tropospheric zenith wet delay.
1
12
[ , , , , , , , , ]
rk
i i i i i i n
X x y z T T N N N
where
12
, , , n
N N N
are the DD ambiguity parameters. The subscript
n
is the
number of DD observations.
,1
k
ii
is the
mm
transition matrix from the time
1i
t
to
i
t
,
k
i
W
is the
1m
error vector of the system state model with zero mean and the covariance
matrix
k
i
W
.
k
i
L
is the
1n
DD measurement vector at the time epoch
i
t
.
k
i
A
is the
mn
- 63 -
design matrix and
k
i
e
is the DD measurement error vector with zero mean and covariance ma-
trix
k
i
.
When multiple kinematic stations are adopted, the system state equation and the observa-
tion equation of the Kalman filter DD positioning can be written as
1
1 1 1
2
2 2 2
,1 1
,1 1
,1 1
j
j j j
k
k k k
ii
i i i
k
k k k
ii
i i i
k
k k k
ii
i i i
X X W
X X W
X X W
(4.10)
and
1 1 1 1
2 2 2 2
j j j j
k k k k
i i i i
k k k k
i i i i
k k k k
i i i i
L A X e
L A X e
L A X e
, (4.11)
where the subscript
j
is the number of kinematic stations. These DD observations are formed
w.r.t. the first reference station
1
r
for each of the kinematic stations
12
, , , j
k k k
respectively.
The non-diagonal block matrixes of the design matrix
i
A
of Eq. (4.11) are zero.
When the multiple reference station method including the a priori constraint is applied, the
additional reference stations are processed as kinematic station. This is similar to Eq. (4.10)
and Eq. (4.11). The system state equations and the observation equations of the Kalman filter
DD positioning can be written as
1
2 2 2
,1 1
,1 1
,1 1
s
s s s
r
r r r
ii
i i i
r
r r r
ii
i i i
k
k k k
ii
i i i
X X W
X X W
X X W
(4.12)
and
2
2 2 2
22
2
2
,,
,,
,
,
s
s s s
s s s
s
rr
r r r
r r k
i i i
i i i
r r k
r r r r k
i i i
i i i
kr
kr
k k k
k
i i i
i i i
L X e
A A A
L X e
A A A
L X e
A A A
, (4.13)
where
s
is the number of the multiple reference stations,
23
, , , s
r r r
are the additional refer-
ence stations and
k
is the kinematic station. If these DD observations are formed w.r.t. the
- 64 -
first reference station
1
r
with each other station including the additional reference stations and
the kinematic stations
23
, , , ,
s
r r r k
, respectively, the non-diagonal block matrixes of the de-
sign matrix
i
A
of Eq. (4.13) are zeros. If the DD observations are formed between the
additional reference stations
23
, , , s
r r r
and the kinematic station
k
, the non-diagonal block
matrixes of the design matrix
i
A
of Eq. (4.13) is the design matrix of the corresponding ob-
servations.
The calculation for the additional reference stations
23
, , , s
r r r
is similar as for the kine-
matic station
k
. But a kind of strong a priori constraint for the state vector is applied for the
additional reference stations according to their known state information, such as
2 2 2
0 0 0
( , , )
r r r
x y z
or
2 2 2 2 2 2
0 0 0 0 0 0
( , , , , , )
r r r r r r
x y z x y z
. Thus, a kind of highly precise a priori accuracy for the state vec-
tor and a small enough spectral densities for the Kalman filter are applied to the state vector
during the Kalman filter processing. For example, the elements of a priori variance and the
spectral densities of the position vector are set to
62
1 10 m
and
92
1 10 m /s
, respectively.
This condition guarantees that the state vectors of the additional reference stations are not
modified during the Kalman filter processing.
And then, the estimates for the kinematic station can be obtained according to the theory of
the two-way Kalman filtering discussed in Section 3.3. This new approach is more flexible
when running the software programs since one only needs to fix those kinematic stations
which shall serve as reference stations.
4.3.2 Experiment and Analysis
In order to investigate the capability of the new approach based on multiple reference stations,
highly dynamic and ultra-long-range GNSS observation data were taken from the GEOHALO
project. The GNSS data were collected by GFZ on June 08. 2012 and comprise GPS and
GLONASS observations. The GNSS antennas of the kinematic stations AIR5 and AIR6 were
selected. These antennas were mounted on the HALO aircraft in its front and middle part, re-
spectively (cf. Figure 4.5). The reference station REF6 was set up by GFZ beside the runway
of the DLR airport Oberpfaffenhofen in Germany. The reference stations RENO and ASCE
were located in central Italy and have been chosen as multiple reference stations. The types of
all these receivers and antennas are given in Table 4.2.
- 65 -
Figure 4.5: The positions of the selected GNSS antennas on the HALO aircraft
Table 4.2: Hardware list of the selected stations from the GEOHALO project
Station Name
Receiver type
Antenna type
(with dome)
AIR5
JAVAD TRE_G3TH DELTA
ACCG5ANT_42AT1
NONE
AIR6
JAVAD TRE_G3TH DELTA
ACCG5ANT_42AT1
NONE
REF6
JAVAD TRE_G3TH DELTA
JAV_RINGANT_G3T
NONE
RENO
TPS ODYSSEY_E
TPSCR3_GGD
CONE
ASCE
LEICA GRX1200+
LEIAR25
NONE
The trajectory of the HALO aircraft on June 08, 2012 and the position of the selected mul-
tiple reference stations REF6, RENO and ASCE are shown in Figure 4.6.
In order to investigate the accuracy and reliability of the algorithm as described above,
firstly, the trajectories of the kinematic stations AIR5 and AIR6 were calculated separately.
Then a time series of the baseline length between the two antennas of the kinematic stations
AIR5 and AIR6 was calculated from the estimated trajectories. The precision and reliability of
the kinematic results were analyzed by checking the temporal variations of the baseline length.
For this reason, the following two computation schemes were performed and compared.
AIR5
AIR6
- 66 -
Figure 4.6: HALO aircraft flight trajectory and the position of reference stations
Scheme 1 (single baseline experiment): The GPS and GLONASS measurement data were
used for the calculation of the trajectories for the kinematic stations AIR5 and AIR6 with
REF6 as the only reference station. Then the distance between AIR5 and AIR6 was calculated.
Scheme 2 (multiple reference stations experiment): The GPS and GLONASS data were
used for the calculation of the trajectories for the kinematic stations AIR5 and AIR6 by the
new approach where REF6, RENO and ASCE served as multiple reference stations. Then the
distance between AIR5 and AIR6 was calculated.
In the data processing for the mentioned calculations, the DD processing of the multiple
reference stations kinematic positioning based on an a priori constraint was carried out using
the HALO_GNSS software, which was developed by the author of this thesis for the GEO-
HALO project. Here, the dual-frequency carrier phase observations were used to form the
ionosphere-free linear combination. The tropospheric zenith wet delay was estimated as well
REF6
RENO
ASCE
- 67 -
and an uncertainty of 10 cm and a spectral density of
82
10 m /s
for the kinematic stations
were assumed. The two-way Kalman filter was used for the parameter estimation.
The results of scheme 1 and scheme 2 are given in Figure 4.7 and Figure 4.8, respectively.
The statistical results of the two schemes are given in Table 4.3.
Table 4.3: The Statistical results for the distance between the two antennas AIR5 and AIR6 (Unit: m)
Reference stations
Min
Max
Mean
SDev
REF6
6.892
7.179
7.046
0.039
REF6, RENO, ASCE
6.984
7.141
7.042
0.019
Figure 4.7: The distance between the two antennas AIR5 and AIR6, the trajectories for these antennas
were estimated by using only one reference station (REF6)
- 68 -
Figure 4.8: The distance between the two antennas AIR5 and AIR6, the trajectories for these antennas
were estimated including the three reference stations REF6, RENO and ASCE
Figure 4.7, Figure 4.8 and Table 4.3 show clearly that the GNSS kinematic positioning is
significantly improved for a large region when using the approach of multiple reference sta-
tions.
The reason for this improvement is the reduction of distance-dependent errors, such as sat-
ellite orbit, tropospheric and ionospheric effects when applying the principle of multiple
reference stations.
4.4 Kinematic Positioning Based on Robust Estimation
In GNSS kinematic precise positioning, observation outliers can affect the accuracy and the
reliability of results. In airborne gravimetry, the requirements on the accuracy and the reliabil-
ity of the results are higher than in “usual” geometric GNSS applications. Thus, outliers
should be eliminated. One method for outlier elimination is Baarda’s data snooping (Baarda,
1968). It is based on testing the occurrence of individual residuals under the assumption that
only lone outliers are present in the relevant set of observations (Leick, 2004, p.173). However,
such an approach is complicated to apply, if multitudinous outliers are expected. In GNSS kin-
ematic precise positioning, it is often the case that not only individual outliers occur in the
observations, and is not a good idea to eliminate each observation featuring a small outlier
since this could cause too few remaining GNSS observations in certain epochs. In this context,
to minimize the impact of observation outliers on the estimates, a robust estimation approach
is applied by using all observations including the outliers but reducing the weights of the ob-
servations individually according to the outlier magnitude. This approach keeps the estimation
- 69 -
algorithm stable and can reduce the impact of the outliers when more than one outlier per
GNSS observations epoch exists.
4.4.1 Principle of Robust Kalman Filter
Robust estimation has been widely studied and applied in geodesy (Koch and Yang, 1998;
Yang et al., 1999; Yang et al., 2001; Yang et al., 2010; Koch, 2013). For example, a Danish
method (Kubik, 1982), and an IGG-I scheme (Zhou, 1989) for different weighting of observa-
tions were established. An IGG-III scheme (Yang, 1994) was proposed for the estimation of
dependent observations.
Based on the classical Kalman filter (cf. Section 3.3), if searching for outliers in the GNSS
observations is included, the robust Kalman filter estimation could be described as (Yang et al.,
2001)
ˆi i i i i i
X X K L A X
(4.14)
and
1
ii
i i i i i
XX
K A A A
, (4.15)
where
i
K
denotes the Kalman gain matrix based on the equivalent weight matrix of observa-
tions
i
L
at epoch
i
.
i
is the covariance matrix of the equivalent weight matrix of the
observations
i
L
with
21
0ii
P
, (4.16)
where
2
0
is the variance of unit weight.
i
P
denotes the equivalent weight matrix of observa-
tions
i
L
. In the independent case,
i
P
is a diagonal matrix with the elements
k
i
p
with
k k k
i i i
p c p
, (4.17)
where
1,2, ,kn
,
n
is the number of DD observations,
k
i
p
is the element of the original
weight matrix
i
P
at the observation at epoch
i
and
k
i
c
denotes the weight factor of the robust
estimation (Yang, 1991; Yang et al., 1999). Based on the calculation experiments, the
k
i
c
can
be constructed for three different cases as follows
- 70 -
1
2
1
12
21
2
1
0
k
i
k
i
kk
ii
k
i
k
i
vt
tv
t
c t v t
tt
v
vt
, (4.18)
where
1
t
and
2
t
are two constants with the empirically found values
11.0 1.5t
and
23.5 6.5t
, respectively.
k
i
v
is the standard residual with
k
ki
i
i
v
v
, (4.19)
where
k
i
v
denotes the element of the residual matrix
i
V
(
i
V
can be found in Eq. (3.33)),
2
i
are the robust variances of the observations (Yang et al., 1999).
The weight factor
k
i
c
of Eq. (4.18) depending on the standard residual
k
i
v
is given in Fig-
ure 4.9. In this figure, the middle segment in the weight factor function given by Eq. (4.18) is
the same as that of the LS case. In practice, most residuals of the GNSS observations will sat-
isfy the condition
1
k
i
vt
, and then
1
k
i
c
, which is just the result of a LS estimation. Both
side parts of the weight factor function are down-weighting parts, which smooth the weight
factor from the LS case to a rejection point. The contribution of the corresponding measure-
ments to the state parameter estimations will be reduced, which should make the Kalman filter
more robust.
Figure 4.9: The weight factors of Kalman filter robust estimation
-t1
t1
t2
-t2
- 71 -
The robust estimation theory can be used in airborne gravimetry to control the impact of
GNSS observation outliers on the estimates for the kinematic platform. In the work for this
thesis the experiment as described in the following section has been carried out.
4.4.2 Experiment and Analysis
Here an experiment study has been performed using the same GNSS data as in the experiment
described in Section 4.3.2. The multiple reference stations REF6, RENO and ASCE, as shown
in Figure 4.6, were used to calculate the trajectories of the kinematic stations AIR5 and AIR6.
The hardware types of the GNSS receivers and antennas have been listed in Table 4.2.
Two kinematic experiments have been carried out in order to show the performance of the
robust estimation. To investigate the accuracy and reliability of the algorithm in this experi-
ment, the trajectories of the kinematic stations AIR5 and AIR6 were estimated separately.
Afterwards, a time series of the baseline length between these two antennas was calculated.
The accuracy and reliability of the kinematic results were analyzed by checking the temporal
variations of the baseline length.
4.4.2.1 GNSS Pseudorange observations Experiment
The GNSS pseudorange observations were used to calculate an initial trajectory for the kine-
matic stations AIR5 and AIR6. The results for the baseline length between AIR5 and AIR6
based on the calculation without using the robust estimation are given in Figure 4.10. This fig-
ure shows large spikes in the time series of the baseline length. This is caused by the fact that
outliers in the GNSS observations were included when applying the non-robust estimation.
The results of the baseline length from with robust estimation are given in Figure 4.11. It
shows that now the large spikes have been removed by the robust estimation.
The statistical results of the GNSS pseudorange positioning without and with robust esti-
mation are listed in Table 4.4.
Table 4.4: The statistical results of the baseline length between AIR5 and AIR6 from GNSS pseudor-
ange positioning (Unit: m)
Estimation type
Min
Max
Mean
SDev
Non-robust estimation
3.506
84.400
7.289
1.154
Robust estimation
2.827
12.770
7.298
0.960
- 72 -
Figure 4.10: The time series for the baseline length between AIR5 and AIR6 based on pseudorange po-
sitioning results without robust estimation
Figure 4.11: The time series for the baseline length between AIR5 and AIR6 based on pseudorange po-
sitioning results with robust estimation
Prior to the application of the robust estimation, the probability distribution of the standard
residuals for the GNSS pseudorange observations was calculated as shown in Figure 4.12.
When the robust estimation was applied by using the constants
11.5t
and
26.0t
(cf. Eq.
(4.18)), 87.13% of the pseudorange observations were processed with the same weight as for
the LS case, 12.50% of observations were processed by down-weighting w.r.t. the LS case ac-
- 73 -
cording to Figure 4.9, and 0.37% of the observations (i.e. outliers) were processed with zero
weight. That means outliers can be removed by the application of the robust estimation.
Figure 4.12: The probability distribution of standard residual for the GNSS pseudorange observations
The improvement due to the robust estimation in this experiment can be analyzed by using
Baarda’s reliability theory (Baarda, 1968; Baarda, 1976) (cf. Section 3.4.3). In order to typi-
cally illustrate exemplarily the typical improvement when applying the robust estimation, the
GNSS observation at the time epoch 12h:10m:09s are taken to consider the extremely largest
spike (84.40 m) as visible in Figure 4.10. For this epoch, the redundancy number
i
r
, the con-
trollability factors
0,i
for the internal reliability, the sensitivity factors
0,i
for the external
reliability and the residuals of the GNSS observations were calculated yielding the values giv-
en in Table 4.5 and Table 4.6 obtained without and with applying the robust estimation,
respectively. Here, for the application of Baarda’s theory the significance number
00.1%
and the power
080%
were used. Therefore, the lower bound features
04.13
(Baarda,
1976; Förstner, 1987).
The residuals of the GNSS pseudorange observations in Table 4.6, exhibit two significant
outliers with the 9th observation of GPS and the 15th observation of GLONASS. These outliers
affect the result of the GNSS positioning (cf. Figure 4.10) and the distribution of the residuals
of the GNSS observations (cf. Table 4.5). The reason of this impact is plausibly explained by
Baarda’s reliability theory (Baarda, 1968):
- 74 -
Without robust estimation, the outlier of the 9th GPS observation is of about 5.193 times its
standard deviation. That means, it cannot be detected by Baarda’s data snooping (Baarda, 1968)
with a probability greater than
080%
if a significance level of
0
1 99.9%
is accept-
ed for the test with statistic standard residual. Meanwhile, the effect of the undetected outlier
on the estimated shift is less than 3.148 times its standard deviation.
Table 4.5: The statistical results for the reliability analysis without robust estimation
(for the GNSS DD observations at time epoch 12h:10m:09s)
GNSS DD
observations
The redundancy
number
i
r
The controlla-
bility factors
0,i
The sensitivity
factors
0,i
The residuals of
GNSS obs.
ˆ
V
(Unit: m)
GPS
1
0.896
4.363
1.408
36.687
2
0.859
4.457
1.676
22.055
3
0.808
4.595
2.015
4.636
4
0.783
4.666
2.172
41.432
5
0.957
4.222
0.877
17.105
6
0.899
4.356
1.384
-35.324
7
0.538
5.632
3.830
-15.863
8
0.851
4.478
1.730
-9.445
9
0.632
5.193
3.148
126.008
GLONASS
10
0.940
4.260
1.043
-34.554
11
0.806
4.599
2.024
-6.285
12
0.913
4.432
1.276
-32.567
13
0.875
4.416
1.563
12.765
14
0.757
4.746
2.339
-19.413
15
0.763
4.729
2.303
-25.594
16
0.815
4.573
1.964
-16.676
17
0.799
4.621
2.072
-6.355
When the robust estimation is applied, the design matrix A (cf. Eq. (3.47)) is not changed
since the observation environment or condition are not modified, but the weights of the obser-
vations P are adjusted. The outlier of 9th GPS observation, which is 4.130 times larger than its
standard deviation can be detected by Baarda’s data snooping, with a probability higher than
080%
if a significance level is
0
1 99.9%
for the test with test statistic standard
residual is used. The important improvement is that the effect of the outlier is removed due to
- 75 -
the sensitivity factor
0,i
being zero. In general, very small sensitivity values
0,i
suggest
checking whether the observation is really necessary, or whether it is superfluous (Förstner,
1987). Therefore, the zero weight is given by robust estimation to the 9th GPS observation.
Table 4.6: The statistical results for the reliability analysis when applying robust estimation
(for the GNSS DD observations at time epoch 12:10:09)
GNSS DD
observations
The redundancy
number
i
r
The controllabil-
ity factors
0,i
The sensitiv-
ity factors
0,i
The residuals of
GNSS obs.
ˆ
V
(Unit: m)
GPS
1
0.853
4.472
1.714
2.630
2
0.841
4.503
1.795
1.013
3
0.799
4.620
2.070
-0.427
4
0.702
4.930
2.693
0.522
5
0.963
4.208
0.807
0.548
6
0.902
4.349
1.362
-5.737
7
0.480
5.964
4.302
1.340
8
0.831
4.531
1.864
-2.221
9
1.000
4.130
0.000
187.438
GLONASS
10
0.880
4.402
1.523
-4.640
11
0.798
4.623
2.077
-6.059
12
0.851
4.478
1.731
-1.924
13
0.794
4.635
2.104
-3.000
14
0.706
4.916
2.666
-2.195
15
0.889
4.380
1.460
-10.555
16
0.808
4.595
2.015
-2.125
17
0.795
4.632
2.098
-0.339
The improvement of the reliability can be seen for the 15th observation of GLONASS as
well. The controllability factor of the internal reliability
0,i
is improved from 4.729 to 4.380.
This means that the gross error can be detected from 4.729 times improved to 4.380 times its
standard deviation. The sensitivity factors of the external reliability
0,i
are improved from
2.303 to 1.460. This means that the effect of undetected errors in the 15th observation onto the
result can be reduced from 2.303 times down to 1.460 times its standard deviation. Others re-
- 76 -
sidual of Table 4.6 have the same noise level as the GNSS pseudorange observations, thus they
do not need to be analyzed.
Hence, the robust estimation can improve the reliability in GNSS kinematic positioning as
shown here exemplary for the selected time epoch.
4.4.2.2 GNSS Carrier Phase Observations Experiment
The GNSS carrier phase observations were taken to estimate final precise trajectories for the
kinematic stations AIR5 and AIR6. The results for the time series of the baseline length be-
tween both antennas without and with robust estimation are given in Figure 4.13 and Figure
4.14, respectively. The corresponding statistical results are listed in Table 4.7. The numbers in
this table as well as the plots in both figures show very clearly that the accuracy of the results
from the carrier phase precise positioning based on the classic Kalman filter based on multiple
reference stations is improved by robust estimation.
Figure 4.13: The time series for the baseline length between AIR5 and AIR6 based on carrier phase
positioning results without robust estimation
- 77 -
Figure 4.14: The time series for the baseline length between AIR5 and AIR6 based on carrier phase
positioning results with robust estimation
Table 4.7: The statistical results of GNSS carrier phase positioning (Unit: m)
Estimation type
Min
Max
Mean
SDev
Non-robust estimation
6.984
7.142
7.043
0.019
Robust estimation
6.959
7.101
7.041
0.016
Prior to the application of the robust estimation, the probability distribution of the standard
residuals for the GNSS carrier phase observations was calculated as shown in Figure 4.15.
When the robust estimation was applied by using the constants
11.5t
and
26.0t
in Eq.
(4.18), 86.50% of carrier phase observations were processed with the same weight as for the
LS case, 12.88% of the carrier phase observations were down-weighted w.r.t. the LS case ac-
cording to Figure 4.9, and 0.62% of observations (i.e. outliers) were processed with zero
weight. That means outliers can be removed by the application of the robust estimation.
- 78 -
Figure 4.15: The probability distribution of the standard residuals for the GNSS
Carrier phase observations
The aforementioned results show the following: If the GNSS observations have larger out-
liers, their impact on the errors of the non-robust Kalman filter results (cf. Figure 4.10 and
Figure 4.13) is more significant for the pseudorange positioning than for the carrier phase po-
sitioning. This has the following reason: cycle slips of the carrier phase data will be introduced
whenever outliers occur in the GNSS carrier phase observations. Then, the new carrier phase
ambiguity is applied for the subsequent GNSS carrier phase observations. The other reason is
that the initial position calculated from the pseudoranges is used to form the design matrix of
the carrier phase equation, such as Eq. (4.9) and the final precise results were calculated from
the carrier phase observations.
But larger outliers still have an impact on the results of the carrier phase processing as
shown in Figure 4.13. This should be taken into account in highly precise positioning. Thus
the robust estimation is applied for airborne gravimetry. The impact of outliers in the GNSS
measurements is suppressed by the equivalent weights in the robust estimation. This can be
seen from the comparison with Figure 4.14 and from the statistical results in Table 4.7. That
means the robust estimation should be the preferred estimation approach for precise GNSS
positioning.
4.5 Conclusions
In this Chapter 4, the observation equations of single difference (SD) and DD are described. It
is shown in the studies of this thesis that the common errors between the satellites and the re-
- 79 -
ceivers (such as tropospheric delays, ionospheric delays, GNSS satellite orbit errors and clock
errors etc) can be eliminated or reduced by DD positioning for small areas. But for large re-
gions, these errors cannot be eliminated satisfactorily by DD positioning.
For large regions a significant improvement of the accuracy and reliability for GNSS DD
positioning can be achieved by using a certain algorithm of multiple reference stations which
is based on the application of an a priori constraint. This improvement has been demonstrated
during the experiments as described in this chapter.
In addition, a robust Kalman filter approach is used to suppress the impact of observation
outliers on the trajectory estimates for airborne gravimetry. It is helpful to apply such a robust
estimation in GNSS precise positioning.
- 81 -
5 Kinematic Positioning Based on Multiple Kinematic
Stations
5.1 Introduction
In the application of GNSS precise positioning in the airborne gravimetric project GEOHALO
or any other kinematic positioning, it is often the case that multiple GNSS receiving equip-
ments are mounted on the kinematic platform, thus known functional or theoretical relations
exist among the unknown parameters of the multiple GNSS receiving equipments. For exam-
ple, the distances among the multiple GNSS antennas are known and fixed and the
characteristic of the tropospheric delay in a small area above the HALO aircraft is similar (see
Figure 5.1). This known relations and characteristic can be and should be made use of to im-
prove the accuracy and reliability of the state estimates and to control the divergences.
Figure 5.1: Multiple GNSS antennas mounted on the HALO aircraft (DLR, 2013a)
For this purpose, in this chapter, a method of GNSS kinematic positioning based on multi-
ple kinematic stations with multiple reference stations is developed firstly. And then the
distances among the multiple GNSS antennas are used as constraints within the state parame-
ters. Finally, a common tropospheric zenith delay parameter is used for the multiple kinematic
stations, which depends on the similar characteristic of tropospheric effects.
Antenna 1
Antenna 2
Antenna 3
- 82 -
5.2 Kinematic Positioning Based on Multiple Kinematic Stations with
Multiple Reference Stations
Following the kinematic positioning based on multiple reference stations (cf. Section 4.3), the
Kalman filter system state equation and observation equation of the kinematic positioning
based on multiple kinematic stations with multiple reference stations can be given as
1
2 2 2
1 1 1
,1 1
,1 1
,1 1
,1 1
s
s s s
j
j j j
r
r r r
ii
i i i
r
r r r
ii
i i i
k
k k k
ii
i i i
k
k k k
ii
i i i
X X W
X X W
X X W
X X W
(5.1)
and
2
2
2 2 1
22
21
11
1
1
1 2 1
21
,
,,
,
,,
,
,
,
, , ,
j
s
sj
ss
s s s
j
s
jj
j j s j j
rk
rr
r r k
rr
i i i i
ii
rk
rr
r r r r k
ii
i i i i
kk
kk
kr
k r k
ii
i i i i
kk
k r k r k k k
ii
i i i i
A A A A
LX
LX
A A A A
LX
A A A A
LX
A A A A
2
1
s
j
r
i
r
i
k
i
k
i
e
e
e
e
, (5.2)
where the subscript
i
is the epoch time
i
t
. The superscript
r
and
k
denote the reference sta-
tions and kinematic stations, respectively.
s
and
j
denote the number of reference stations
and kinematic stations.
i
X
and
1i
X
are unknown parameter vectors at epoch
i
and
1i
,
respectively, in which including state vector
,,
i i i
x y z
or
, , , , ,
i i i i i i
x y z x y z
, tropospheric
zenith wet delay
i
T
, which will be discussed in following section, and DD ambiguity pa-
rameters
12
, , , n
N N N
, here subscript
n
is the number of corresponding DD
observations, hereby
12
[ , , , , , , , ]
i i i i i n
X x y z T N N N
.
,1ii
is the transition
matrix from epoch
1i
to
i
and
i
W
is the error vector of the system state model with zero
mean and covariance matrix
i
W
.
i
L
is the DD measurement vector of corresponding observa-
tions at epoch
i
.
i
e
is the DD measurement error vector with zero mean and covariance
matrix
i
.
i
A
is the design matrix, if the DD observations are formed by the first reference
station
1
r
with each station including additional reference stations and kinematic stations
- 83 -
2 3 1 2
, , , , , , ,
sj
r r r k k k
respectively, the non-diagonal block matrixes of the design matrix of
Eq. (5.2) are zeroes. If the DD observations are formed between the additional reference sta-
tions
23
, , , s
r r r
and the kinematic stations
12
, , , j
k k k
, the non-diagonal block matrixes of
design matrix of Eq. (5.2) are the design sub-matrix of the corresponding DD observations.
According to the estimation theory of two-way Kalman filtering, which has been discussed
in Section 3.3, the state estimate of multiple kinematic stations can be obtained.
5.3 Kinematic Positioning Based on an A priori Distance Constraint
When opportunities for the application of certain constraints in the GNSS kinematic state pa-
rameters of multiple GNSS sensors system are given, they should be taken into account for the
improvement of the positioning accuracy and reliability. A typical example for such a con-
straint option is the known baseline length between two GNSS antennas mounted on the
kinematic platform. This distance can be given as
1 2 1 2 1 2 1 2
2 2 2
,k k k k k k k k
i i i i i i i
d x x y y z z
, (5.3)
where
12
,kk
i
d
denotes the distance between the two antennas of the kinematic stations
1
k
and
2
k
,
( , , )
i i i
x y z
is the position vector of the respective kinematic station at epoch
i
.
There exist several approaches to deal with a priori constraints in Kalman filter applica-
tions, which includes linearized constraints, nonlinear constraints, time-varying constraints
and so on (Tahk and Speyer, 1990; Doran, 1992; Simon and Chia, 2002; Julier and LaViola,
2007; Ko and Bitmead, 2007; Yang et al., 2010; Yang et al., 2011b). In this section, a dedicated
kind of a priori distance constraints is applied for the GNSS kinematic positioning.
5.3.1 The A priori Distance Constraints
Due to the presence of the GNSS antenna phase centre offset and its variation (cf. Section
2.4.2) the precise distances between GNSS antennas are difficult to measure in a simple way
(for instance by using a ruler). Therefore in this study, the DD processing on ultra-short base-
lines is applied to measure the distances among GNSS antennas with millimetre accuracy in a
small area. These precise distances are used as a priori constraints with realistic accuracies.
The way to deal with these constraints is the so-called “pseudo-observation method”, that is to
express the constraints as perfect observation equations which strengthen the structure of the
measurement space directly (Yang et al., 2010; Yang et al., 2011b).
When multiple kinematic stations are mounted on a kinematic platform, the linearized con-
straint equations can be written as
- 84 -
ii
D B X
, (5.4)
where
D
denotes the
1u
distance constraint vector at every epoch.
i
B
is an
um
design
matrix at epoch
i
.
i
X
is the
1m
unknown parameter vector at epoch
i
and
is a distance
constraint uncertainty vector with zero mean and covariance matrix
d
.
To express the distance constraints among the kinematic antennas as perfect observation
equations, the error equation of Eq. (5.4) can be combined with the GNSS observation error
equation Eq. (3.33) as
ˆ
ii i
i
d
ii
VA L
X
VB D
. (5.5)
The covariance matrix of Eq. (5.5) can be written as
0
0
i
d
.
According to the classical Kalman filter theory (Kalman, 1960; Kalman and Bucy, 1961;
Huep, 1985), the parameter estimate depending on a priori distance constraints at epoch
i
can
be expressed as (Tahk and Speyer, 1990; Yang et al., 2010)
1
0
ˆ
0
ii
i
ii
i i i i
XXii i i
i i i i
X
dii
i i i i
XX
A A A B L A X
X X A B D B X
B A B B
(5.6)
and the a posteriori covariance matrix
ˆi
X
can be expressed by
1
ˆ
0
0
ii
i i i
i
ii
i i i i
XXii
ii
X X X
X
di
i i i i
XX
A A A B A
AB B
B A B B
. (5.7)
In order to illustrate the performance of the distance constraints in GNSS precise position-
ing for multiple kinematic stations, a realistic shipborne gravimetric experiment is given in
Section 5.5.
5.4 Kinematic Positioning Based on a Common Tropospheric Parameter
The tropospheric delay is an important error source in precise GNSS applications, which has
been introduced briefly in Section 2.4.4. Generally, in GNSS data analysis three methods are
used to estimate the tropospheric delay.
The first method uses the tropospheric delay model to correct the tropospheric error. This
tropospheric model correction method works well for the hydrostatic or dry tropospheric delay
- 85 -
correction (Janes et al., 1991) but not for the wet tropospheric delay, since this is difficult to be
corrected completely (Mendes and Langley, 1998; Wei and Ge, 1998; Singh et al., 2014).
The second method is DD processing. If the reference station and the kinematic station are
closely located to each other, for example less than 10 km on ground, the impact of the tropo-
spheric delay can be considerably reduced, since both tropospheric delays are approximately
the same. However, when applying the DD processing in large regions, where the kinematic
station is too far away from the reference station (more than a few hundred kilometres) (Shi
and Cannon, 1995) or the difference of their elevations is too large (from a few hundred to
thousands of metres) (Tiemeyer et al., 1994; Shi and Cannon, 1995), the difference between
the atmospheric delays on the kinematic station and on the reference station will cause an in-
crease of this error impact. Thus, in the latter cases it is difficult to eliminate the tropospheric
delay completely by the method of DD processing.
Thirdly, after the application of the atmospheric correction model and the DD processing,
the remaining residual tropospheric delay can be corrected by estimating a zenith wet tropo-
spheric delay parameter (Singh et al., 2014).
In order to obtain a high positioning accuracy by DD processing in large regions, the resid-
ual tropospheric delay error should be eliminated by estimating a wet tropospheric zenith
delay parameter.
5.4.1 Tropospheric Delay Estimation
The Slant Total Delay (STD) between satellite and receiver can be written as (Davis et al.,
1985; Chen and Herring, 1997)
STD ( ) ZHD ( ) ZWD
dw
M e M e
, (5.8)
where
ZHD
and
ZWD
are the zenith hydrostatic and the wet delay, respectively.
()
d
Me
and
()
w
Me
are the associated mapping functions for the hydrostatic and the wet atmospheric
delay, which depend on the elevation angle
e
.
Though the hydrostatic and the wet mapping function can be applied to map the tropo-
spheric delay correction, the wet tropospheric delay component is usually hard to model.
Therefore, the state vector element corresponding to the tropospheric delay correction can be
considered as a correction for the wet component as well as for the total tropospheric delay
when using a wet mapping function (Abdel-salam, 2005).
When the tropospheric wet delay parameter is considered, the slant total delay between the
satellite and the antenna
k
at epoch
i
can be given as
- 86 -
STD ( ) ZHD ( ) ZWD
k k k k
i d i w i i
M e M e T
, (5.9)
where the
k
T
is the zenith wet delay parameter to be estimated at station
k
. For this parame-
ter an empirical a priori constraints should be added according to the variations of the local
atmosphere or the model for regional atmospheric corrections.
The nature of the tropospheric wet delay allows for modelling its components as a stochas-
tic process. Therefore, the transition matrix and the noise matrix of the Kalman filter state
equation will consist of elements as given by (Chen, 1998; Abdel-salam, 2005)
trop
,1 1
ii
(5.10)
and
trop
trop
i
wqt
, (5.11)
where the value of the spectral density of the tropospheric delay parameter
trop
q
is the com-
mon value used by many research works to model the stochastic properties of the tropospheric
delay (Zumberge et al., 1997; Niell, 2001). The spectral density of the troposphere was chosen
to be
92
10 m /s
which was observed to perform well with many datasets (Abdel-salam, 2005).
5.4.2 A Common Tropospheric Wet Delay Parameter
It is often the case in airborne gravimetry that multiple GNSS receiving equipments are
mounted on the kinematic platform. Since the tropospheric delay is uniform over such a small
area like a kinematic platform, a common tropospheric delay parameter can be introduced as
described below.
Assuming that all kinematic stations (
12
, , , j
k k k
) in Eq. (5.1) and Eq. (5.2) are mounted
on the kinematic platform, the vector of the unknown parameters vector can be written as
1 1 1 1 1 1 1 1
12
12
[ , , , , , , , ]
[ , , , , , , , ]
j j j j j j j j
k k k k k k k k
i i i i i n
k k k k k k k k
i i i i i n
X x y z T N N N
X x y z T N N N
. (5.12)
The tropospheric zenith wet delay parameters out of Eq. (5.12) are
12
, , , , ,
i
k
kk
i i i
T T T
, (5.13)
where the elements of tropospheric zenith delay parameters are similar in the small area above
the kinematic platform.
- 87 -
Thus, a common tropospheric zenith wet delay parameter is introduced for the multiple
kinematic stations. The common tropospheric zenith wet delay parameter in the unknown pa-
rameter vector
k
X
will be changed as
,,
i
T
. (5.14)
This method not only reduces the amount of estimated parameters, but also harmonizes the
actual conditions.
Generally, for solving the atmospheric delay parameters, there exist several usual methods
such as single parameter estimation, multi-parameter estimation, piecewise linear estimation,
stochastic estimation. In this study, a stochastic estimation is applied (Chen, 1998; Wei and Ge,
1998; Abdel-salam, 2005), where the tropospheric wet delay is computed on the basis of a
first-order Gaussian-Markov model.
5.5 Experiment and Analysis
Since the distance constraints and a common tropospheric delay parameter were used for the
study in this chapter, the accuracy evaluation method which uses the changes of the baseline
length between the kinematic antennas cannot demonstrate the capability of the specific meth-
od as developed in this chapter. Thus, GNSS DD positioning results for a small region are
used as “true value” for the comparison with results of the developed method applied for a
large region. For this purpose, GNSS data of a shipborne gravimetry campaign on the Baltic
Sea were used here.
From June 18 to 27, 2013, ten days of GNSS and gravimetric data were collected by GFZ
and the Federal Agency for Cartography and Geodesy, Germany (BKG) in the Baltic Sea (Ost-
see in German) near Greifswald, Germany. During this campaign, three GNSS antennas were
mounted on the used ship. The positions of the GNSS antennas relative to each other are
shown in Figure 5.2.
Figure 5.2: The relative positions of the kinematic GNSS antennas on the ship
KIN1
KIN2
KIN3
- 88 -
In order to investigate the capability of the strategy, in this chapter the GNSS data of the
first day of the Baltic Sea gravimetric experiment (June 18, 2013) were selected for testing.
The GNSS stations KIN1 and KIN3 were selected as multiple kinematic stations. The baseline
length of 26.342 m was used as the distance constraint in the following experiment. The sta-
tions 0801 and 0775 of the Satellite Positioning Service (SAPOS, www.sapos.de), which is
organized by the Working Committee of the Surveying Agencies of the States of the Federal
Republic of Germany (AdV), were taken as small regional reference stations (i.e. distances <
30 km) (cf. Figure 5.3). The IGS stations WARN and POTS were chosen as large regional ref-
erence stations (distance ~50-200 km). The tracks of the ship and the positions of the reference
stations are shown in Figure 5.4. The hardware types of all GNSS receivers and antennas are
given in Table 5.1.
Figure 5.3: The track of the ship (green curve) and the positions of the reference stations in a small re-
gion (yellow triangles) during the Baltic Sea shipborne gravimetric campaign on June 18. 2013
For the GNSS data processing, the HALO_GNSS software was used based on the methods
described in this chapter. Here, the dual-frequency carrier phase observations were used to
form the ionosphere-free linear combination. The tropospheric zenith wet delay was estimated
as a random walk process, the initial uncertainty of which was assumed to be 10 cm and its
spectral density
12 2
10 m /s
. The reasons for the selection of these values are the slow motion
of the ship and the small changes in the height profile (Abdel-salam, 2005). The two-way
Kalman filter was used for the parameter estimation and the GPS and GLONASS observation
- 89 -
data were used for the calculation of the trajectories for the multiple kinematic stations KIN1
and KIN3.
Figure 5.4: The track of the ship (green curve) and the positions of reference stations (yellow triangles)
for the Baltic Sea shipborne gravimetric experiment on June 18. 2013
Table 5.1: Hardware list of the selected stations from Baltic Sea gravimetric experiment
Station Name
Receiver type
Antenna type
( with dome)
KIN1
JAVAD TRE_G3TH DELTA
LEIAS10
NONE
KIN3
JAVAD TRE_G3TH DELTA
ACCG5ANT_42AT1
NONE
0801
TPS NET-G3A
TPSCR.G3
TPSH
0775
TPS NET-G3A
TPSCR.G3
TPSH
WARN
JPS LEGACY
LEIAR25.R3
LEIT
POTS
JAVAD TRE_G3TH DELTA
JAV_RINGANT_G3T
NONE
- 90 -
Experimental studies were performed with five computational schemes:
Scheme 1 (small region experiment): The trajectories of the multiple kinematic stations
KIN1 and KIN3 were calculated by a single adjustment for the same time with a common
tropospheric delay parameter, where 0801 and 0775 served as multiple reference stations in
small region.
The baseline length between the two kinematic stations KIN1 and KIN3 was calculated for
each epoch. The apparent changes of this baseline length are shown in Figure 5.5 and the sta-
tistical results are given in Table 5.2. These results can be treated as the “true value” (with the
standard deviations of 1.5 cm in the related distance) of the antenna movement. Since the re-
sults of KIN1 and KIN3 are nearly the same, only the KIN1 is used for comparisons with the
following schemes.
Figure 5.5: The distance between two antennas, KIN1 and KIN3; the trajectories of these antennas were
estimated by using two reference stations, 0801 and 0775
Scheme 2 (large region experiment): The trajectories of the multiple kinematic stations
KIN1 and KIN3 were calculated by a single adjustment for the same time with two independ-
ent tropospheric delay parameters, where WARN and POTS served as multiple reference
stations located in the large region.
The baseline length between two kinematic stations KIN1 and KIN3 was calculated for
each epoch. The apparent changes of this baseline length are shown in Figure 5.6 and the sta-
tistical results are given in Table 5.2. The results showed that the precision of scheme 2 is
lower than that of scheme 1. The tracks of KIN1 for the large region were compared with the
“true value” of scheme 1. The comparison results are shown in Figure 5.7 and the statistics of
the comparison as given in Table 5.3 show that the achieved accuracy is 5.8 mm, 4.4 mm and
36.7 mm for the north, east and up directions, respectively.
- 91 -
Figure 5.6: The distance between two antennas, KIN1 and KIN3; the tracks of these antennas were es-
timated by using two reference stations, WARN and POTS
Table 5.2: The statistical results for the baseline length between KIN1 and KIN3 (Unit: m)
Scheme
Region type
Min
Max
Mean
SDev
1
Small region
26.282
26.406
26.342
0.015
2
Large region
26.261
26.428
26.338
0.022
Figure 5.7: The results for KIN1 for the large region with two independent tropospheric delay parame-
ters compared with scheme 1 (differences)
- 92 -
Scheme 3: Based on the scheme 2, the a priori distance constraints were used to calculate
the tracks of KIN1 and KIN3.
Scheme 4: Based on the scheme 2, a common tropospheric delay parameter between the
kinematic stations was used to calculate the tracks of KIN1 and KIN3.
Scheme 5: Based on the scheme 2, a common tropospheric delay parameter and the dis-
tance constraints were used to calculate the tracks of KIN1 and KIN3.
The trajectories of KIN1 for the schemata 3, 4 and 5 were compared with the “true value”
of scheme 1, respectively. The comparisons of the results are shown in Figure 5.8, Figure 5.9
and Figure 5.10, respectively. And the statistics results of the comparisons are given in Table
5.3.
Figure 5.8: The results for scheme 3 for KIN1 for the large region with two independent tropospheric
delay parameters and distance constraints compared with scheme 1 (differences)
From the above computational results, the following conclusions can be drawn:
(1) The precision of the kinematic state estimates provided by the DD processing using
reference stations in the large region is lower than those provided in the small region.
(2) The accuracies of the estimated kinematic state parameters can be improved by the dis-
tance constraints, because of the relationships or known information between states
parameters are used in GNSS precise kinematic positioning.
(3) The method which calculates multiple kinematic stations together and uses a common
tropospheric delay parameter for these kinematic stations, can improve the kinematic
- 93 -
positioning accuracy as well, since it can reduce the amount of unknown parameters
and represents the reality better.
(4) By integration of distance constraints and using a common tropospheric delay parame-
ter in GNSS precise kinematic positioning, the accuracies for the Up direction (cf. the
corresponding mean and RMS values of Table 5.3 ) can be further improved.
Figure 5.9: The results of scheme 4 for KIN1 with reference stations in the large region with a common
tropospheric delay parameter for KIN1 and KIN3 compared with scheme 1 (differences)
Figure 5.10: The results of scheme 5 for KIN1 with reference stations in the large region with a com-
mon tropospheric delay parameter and distance constraints compared with scheme 1 (differences)
- 94 -
Table 5.3: Statistical results of each scheme compared with scheme 1 (Unit: mm)
Scheme
Direction
Min
Max
Mean
RMS
2
Two independent
tropospheric pa-
rameters
North
-27.0
23.0
-6.4
5.8
East
-23.5
14.8
-6.4
4.5
Up
-165.7
165.9
3.5
36.8
3
A priori Distance
constraints
North
-27.7
15.8
-6.4
5.8
East
-21.6
16.4
-6.6
4.2
Up
-161.6
145.4
0.0
33.1
4
A common tropo-
spheric parameter
North
-28.5
14.2
-7.2
5.2
East
-21.5
11.5
-6.4
4.1
Up
-151.4
142.1
0.6
27.8
5
A priori Distance
constraints and a
common tropo-
spheric parameter
North
-28.6
15.9
-6.8
5.4
East
-21.4
11.6
-6.6
4.1
Up
-112.5
108.4
0.0
27.3
5.6 Conclusions
In this chapter a method of GNSS kinematic positioning based on multiple kinematic stations
with multiple reference stations was described. It can be used to provide a basis for the re-
search when multiple kinematic stations exist.
It is often the case that multiple GNSS receiving equipments are mounted on the kinematic
platform, so that known functionals or theoretical relations exist among the unknown parame-
ters of the individual GNSS receiving equipments. The distance among the multiple GNSS
antennas is known and fixed, which has been used as a priori distance constraint to improve
the RMS accuracy of the state estimates. The characteristics of the tropospheric delay in a
small area are similar. Therefore, a common tropospheric delay parameter was used for the
multiple kinematic stations, which can enhance the accuracy and reliability of the state esti-
mates as well, because it can reduce the amount of unknown parameters and represents the
reality better.
Finally, a shipborne gravimetric experiment is used to test above methods since its results
with reference stations in the small region can be used as “true value” to compare with those
based on the reference stations in the large region. Through integration of distance constraints
and a common tropospheric delay parameter accuracies can be further improved.
- 95 -
6 Kinematic Positioning Based on GNSS Integration
6.1 Introduction
The integration of multiple GNSS satellite systems may be considered as a major milestone in
GNSS precise positioning, because it can dramatically improve the reliability and productivity
of GNSS positioning (Wang et al., 2001).
It is well known that, for GNSS based precise positioning systems, the accuracy, availabil-
ity and reliability of the positioning results are strongly dependent on the number of satellites
being tracked, especially in kinematic positioning. However, in some situations, GNSS signals
are blocked or essentially degraded by natural or artificial obstacles, such as in urban canyons
(shown in Figure 6.1), mountainous areas and in deep open-cut mines, so that the number of
visible satellites may not be sufficient for the positioning function (Wang et al., 2001;
Angrisano et al., 2013). In case of a single satellite system, when there are not enough obser-
vations or even no observation, kinematic positioning results are of poor accuracy, unreliable,
or even fail.
Figure 6.1: GNSS integration (OpeNaviGate, 2012)
- 96 -
In applications of long baseline GNSS differential positioning (such as the GEOHALO
project), the number of simultaneously visible satellites will be reduced with the increase of
the baseline length. A possible way to fill this gap is using a GNSS multi-constellation receiver,
considering the combined use of GPS with the signals of other GNSS systems such as
GLONASS, Galileo, the BeiDou navigation Satellite System (BDS) etc. In general, the usage
of multiple GNSS systems will increase the number of simultaneously visible satellites and
improve the reliability.
Furthermore, in this chapter, a GNSS integration algorithm based on Helmert’s variance
components estimation is used to reasonably adjust the weights and the contributions among
multiple GNSS systems.
Currently Galileo has only four satellites in orbit, in the In-Orbit Validation (IOV) phase
(ESA, 2014), while BDS is in the development phase and only the regional phase (Asia-
Pacific area) has been completed at the end of 2012 (Zhao et al., 2013; BDS, 2014). The en-
hancement of the Russian space program has made GLONASS an ideal candidate to form a
multi-constellation with GPS (Angrisano et al., 2013). In the following section, GPS and
GLONASS are considered because they are the only GNSS systems which are declared to be
fully operational.
6.2 GPS and GLONASS Integration
There is a strong interest in including GLONASS satellites in any GPS positioning solution. It
is well known that additional satellites strengthen the solution and increase its reliability
(Leick 1998). Unlike GPS, GLONASS satellites transmit signals at different frequencies,
which causes a significant complexity in terms of the modelling and ambiguity resolution for
an integrated GPS and GLONASS data processing system (Wang et al., 2001).
Referring to the coordinates on the Earth’s surface, the terrestrial reference system of GPS
is the World Geodetic System 1984 (WGS-84) and for GLONASS it is the Parametry Zemli
1990 (PZ-90). GPS time (GPST) is connected to the Coordinated Universal Time (UTC)
which is maintained by the US Naval Observatory (USNO). The GLONASS time (GLOT)
scale is connected to TUC (RU) which is a kind of UTC maintained by Russia. The transfor-
mation formula between GPST and GLOT is
GPST GLOT r u g
, (6.1)
where
r
is the time offset between GLOT and UTC(RU).
u
is the time offset between
UTC(RU) and UTC(USNO).
g
denotes the time offset between GPST and UTC(USNO)
(Hofmann-Wellenhof et al., 2008, p. 315; Angrisano et al., 2013). In this study, WGS84 and
- 97 -
GPST were used as reference systems. More details about the reference systems and their
transformations into one another can be found in Roßbach (2001), Habrich (2000), Xu (2007)
and Hofmann-Wellenhof et al. (2008).
As well known, carrier phase measurements are much more precise than pseudorange ob-
servations. Hence, they are the primary measurements for precise positioning. The DD carrier
phase observation equations for GPS and GLONASS (Roßbach, 2001; Wang et al., 2001) at
epoch
i
are given as
GPS GPS GPS GPS
GPS p q p q GPS p q p q GPS p q
k r k r k r k r k r
N T e
(6.2)
and
ICB
GLO GLO GLO GLO GLO GLO GLO GLO GLO
GLO GLO GLO GLO GLO
p p q q p q p p q q
k r k r k r k r k r
p q p q p p q q
k r k r k r k r
NN
T e e
, (6.3)
where the superscripts
p
and
q
indicate the reference satellite and the other satellite, respec-
tively. The subscripts
r
and
k
mean the reference station and the kinematic station,
respectively.
denotes the single difference operator between the reference station and the
kinematic station, for example
p
kr
denotes
pp
kr
.
is the DD operator among two
satellites and two stations, for example
p q p q
k r k r k r
.
is the wavelength of the
carrier phase of the satellite signal,
is the carrier phase observation between a satellite and a
receiver.
is the geometrical distance between a satellite and a receiver with
2 2 2
p p p p
r r r r
x x y y z z
. (6.4)
N
is the carrier phase ambiguity and
T
is the atmospheric delay along the path.
e
is the
measurement noise of the carrier phase and
ICB
is the inter-channel bias for measurements
on different carrier frequencies for different satellites. It is difficult to separate the DD carrier
phase bias
ICB GLO
pq
kr
from the ambiguity parameter. Thus, the float solution for the
GLONASS DD ambiguity is used in this study.
From Eq. (6.2) and (6.3), the error equations of GPS and GLONASS integration can be
written as
1
1 1 1 1
2 2 2 2
2
0
0
xB
V A N L
I
yT
V A N L
BI
z
, (6.5)
- 98 -
where
i
V
(
1,2i
) denotes the
1
i
n
residual vector.
i
A
is the
3
i
n
design matrix of state
parameters
,,x y z
of the kinematic station.
i
B
is the
1
i
n
coefficient matrix of the tropo-
spheric delay parameter
T
. The state parameters
,,x y z
and the tropospheric delay
parameter
T
are common parameters of the kinematic station between two GNSS systems
at the same epoch. The DD ambiguity parameters
12
NN
between the two systems
are uncorrelated. The observation vectors
1
L
and
2
L
are independent of each other.
The estimates for the kinematic station can be calculated using the estimation theory (cf.
Chapter 3). The accuracy of the multiple GNSS measurement vectors is usually different. In
order to balance the contributions of the various GNSS measurements, the final result is ob-
tained by using the principle of Helmert’s variance component estimation.
6.3 GNSS Integration Based on Helmert’s Variance Components Estima-
tion
The variance component estimation (VCE) was developed by Friedrich Robert Helmert in
1907 (Helmert, 1907) as cited by (Xu et al., 2007) and it has been adopted for a variety of ap-
proaches (Sahin et al., 1992; Crocetto et al., 2000; Koch and Kusche, 2002; Kusche, 2003;
Yang et al., 2005; Vennebusch et al., 2007; Xu et al., 2007; Teunissen and Amiri-Simkooei,
2008; Wang et al., 2009a). The purpose of VCE is to find the realistic and reliable variance
components for the observations in order to construct correctly the a priori covariance matrix
of the observations (Sahin et al., 1992).
In the GNSS kinematic positioning, let us consider
k
types of GNSS measurements as
1
i
n
vectors
i
L
for
1,2, ,ik
. The GNSS observation equations can be written as
1 1 1
2 2 2
k k k
L A e
L A e
X
L A e
. (6.6)
The matrix
i
A
is the
i
nm
design matrix.
X
denotes a
1m
vector of the unknown param-
eters.
i
e
is the theoretical measurement error vector with zero mean and the covariance matrix
is
i
with
21
0,i i i
P
. The value
2
0,i
is the variance of unit weight and
i
P
denotes the
weight matrix of the observations.
- 99 -
The observation vectors may be internally correlated (the noise may be coloured) but we
assume them to be uncorrelated between groups. Moreover, we assume that the covariance
matrices of these observation groups are known a priori up to some scaling factors, the
k
noise levels or variance components
2
0,i
,
21
10,1 1
21
20,2 2
21
0,
00 00
00
00
00 00
kkk
P
P
P
. (6.7)
Eq (6.7) indicates that
12
, , , k
L L L
do not have the same variance of unit weight
2
0
. For-
tunately, the variances for individual measurements or the variance components for the
grouped measurements can a posteriory be estimated based on the measurement residuals from
Least squares solution for the parameters.
The purpose of VCE is to estimate the variance factors
2
0,i
( 1,2, , )ik
and iteratively
or adaptively adjust the measurement weights or variances using the available measurements
(Wang et al., 2009a).
Based on Helmert’s most popular VCE method (Helmert, 1907; Förstner, 1979; Koch,
1986; Kusche, 2003; Bähr et al., 2007), the variance components can be estimated from the
following equation:
2
11 12 1 0,1 1 1 1
2
21 22 2 0,2 2 2 2
2
12 0,
ˆ
ˆ
ˆ
k
k
k k kk kk k k
U U U V PV
U U U V PV
U U U V PV
, (6.8)
where
1 1 2
2 ( ) ( )
ii i i i
U n tr N N tr N N
, (6.9)
and
11
()
ij ji i j
U U tr N N N N
, (6.10)
with
, 1,2, ,i j k
and
ij
,
1
k
i
i
NN
, (6.11)
i i i i
N A PA
(6.12)
- 100 -
and
ˆ
i i i
V A X L
, (6.13)
where
i
V
is the estimated error vector and
ˆ
X
is the estimated parameter vector.
The inverse of the coefficient matrix in Eq. (6.8) makes the method inconvenient. It is the
most time-consuming computation step. Correspondingly, a number of simplified formulas
have been developed (Förstner, 1979; Grafarend et al., 1980; Grodecki, 1997; Bähr et al.,
2007). Starting from the Helmert estimator, Förstner (1979) derived a simplified estimator
with some considerable advantages. The name “Förstner method” has been given by Persson
(1980), and the simplified estimator is obtained by
2
0,
ˆi i i
i
i
V PV
r
(6.14)
with
1
()
i i i
r n tr N N
. (6.15)
Obviously, Eq. (6.14) can be applied easily because it does not need any extra calculations.
This advantage is caused by the fact that the redundant contribution of the measurements is
always required in the data processing for the purpose of outlier detection or reliability analy-
sis (Wang et al., 2009a). It should be pointed out that outliers must be removed prior to the
application of VCE.
6.4 Experiment and Analysis
In order to investigate the capability of the strategy as described in this chapter, the kinematic
GNSS data of the GEOHALO mission and the static GNSS data of the IGS network were used.
6.4.1 Kinematic Experiment
The kinematic GNSS data of the GEOHALO flight on June 6, 2012 were selected for the test.
The selected kinematic stations are AIR5 and AIR6 in the front and back part of the HALO
aircraft, respectively (cf. Figure 4.5). The station REF6, installed by GFZ next to the runway
of the DLR airport in Oberpfaffenhofen, and the station ASCE located in the central region of
Italy were selected as reference stations. The hardware types of all receivers and antennas are
the same as in the experiment described in Section 4.3.2 (cf. Table 4.2). The trajectory of the
HALO aircraft on June 6, 2012 and the positions of the selected reference stations are shown
in Figure 6.2. The selected data contain GPS and GLONASS observations with a sampling
rate of 20 Hz.
- 101 -
For the data processing, the HALO_GNSS software was used. This software for precise
kinematic positioning and velocity determination was developed at GFZ in the frame of this
thesis. The IGS precise ephemeris (Kouba, 2009b) were chosen and the GPS and GLONASS
satellite elevation angle was limited to
10
. The number of GPS, GLONASS and GPS +
GLONASS satellites is shown in Figure 6.3 as blue, green and red line, respectively. As obser-
vations we used the ionospheric-free linear combination (Hofmann-Wellenhof et al., 2008, p.
111; Xu, 2007, p. 97) of pseudorange and carrier phase observations, respectively. The two-
way Kalman filtering method was used for the parameter estimation.
Figure 6.2: The HALO aircraft flight trajectory on June 6, 2012 and the location of the selected refer-
ence stations
REF6
ASCE
- 102 -
Figure 6.3: Number of the selected satellites of GPS (blue line), GLONASS (green line) and
GPS+GLONASS (red line) for the kinematic experiment (GEOHALO flight on June 6, 2012)
In order to investigate the accuracy and reliability of the developed algorithm and software,
firstly the positions of the kinematic stations AIR5 and AIR6 were calculated with REF6 and
ASCE as multiple reference stations, and secondly the baseline length between two kinematic
stations AIR5 and AIR6 was calculated for every epoch. The accuracy and reliability of the
kinematic results were then validated based on the apparent variations of the baseline length.
The following four schemes were compared.
Scheme 1 (GPS as a single GNSS system): GPS as a single system was used and the posi-
tions of the kinematic stations AIR5 and AIR6 were calculated together with REF6 and ASCE.
The result for the baseline is shown in Figure 6.4 as a time series and its statistical characteris-
tics are listed in Table 6.1.
Figure 6.4: The time series for the apparent baseline length variation between AIR5 and AIR6 based on
GPS as a single GNSS system
- 103 -
Scheme 2 (GLONASS as a single GNSS system): GLONASS as a single system was used
and the positions of the kinematic stations AIR5 and AIR6 were calculated together with REF6
and ASCE. The result for the baseline is shown in Figure 6.5 as a time series and its statistical
characteristics are listed in Table 6.1.
Scheme 3 (GPS and GLONASS as an integrated GNSS system with 1:1 weights): GPS and
GLONASS were used together and the positions of the kinematic stations AIR5 and AIR6
were calculated together with REF6 and ASCE. The ratio of the a priori observation accuracy
between GPS and GLONASS systems is set to 1:1. The result for the baseline is shown in
Figure 6.6 as a time series and its statistical characteristics are listed in Table 6.1.
Figure 6.5: The time series for the apparent baseline length variation between AIR5 and AIR6 based on
GLONASS as a single GNSS system
Figure 6.6: The time series for the apparent baseline length variation between AIR5 and AIR6 based on
GPS and GLONASS as an integrated system with 1:1 weights
- 104 -
Scheme 4 (GPS and GLONASS as an integrated system with Helmert’s VCE): GPS and
GLONASS were used together and the positions of the kinematic stations AIR5 and AIR6
were calculated together with REF6 and ASCE. The relation of the observations accuracies
between GPS and GLONASS was determined by Helmert’s VCE. Outliers must be removed
prior to the application of VCE. The ratios between the weights for the GPS and GLONASS
observations are shown in Figure 6.7. The apparent changes of the baseline length are shown
in Figure 6.8 as a time series and their statistical features are listed in Table 6.1.
Figure 6.7: The ratios between the weights for the GPS and GLONASS observations
Figure 6.8: The time series for the apparent baseline length variation between AIR5 and AIR6 based on
GPS and GLONASS as an integrated system with Helmert’s VCE
- 105 -
Table 6.1: The statistical results of schemata 1-4
Scheme
Apparent distance between AIR5 and AIR6 (Unit: m)
Min
Max
Mean
SDev
1
G (GPS)
6.960
7.114
7.034
0.019
2
R (GLONASS)
6.958
7.121
7.034
0.026
3
G+R (1:1 weights)
6.989
7.097
7.038
0.014
4
G+R (Helmert weights)
6.997
7.078
7.036
0.009
These results show the following:
(1) By the comparison of scheme 1 and 2, it can be concluded that the kinematic position-
ing accuracy based on GLONASS as a single system is a little bit lower than the GPS-
only accuracy.
(2) The results for scheme 3 show that the method of GPS and GLONASS integration can
improve the accuracy of GNSS kinematic positioning. This is caused by an increased
number of observations as well as by more available GNSS satellites and improves the
stability and reliability (as an increase in redundancy numbers) of the parameter esti-
mation.
(3) By the comparison of schemes 3 and 4 it can be found that the method which is based
on of Helmert’s VCE can further improve the accuracy of the integrated GNSS kine-
matic positioning, because it can reasonably adjust the weights and the contributions of
the different GNSS systems.
6.4.2 Static Experiment
In order to test the accuracy and reliability of the algorithm and software in a static data exper-
iment, the static GNSS observation data of the IGS sites TITZ and FFMJ (both in Germany)
on January 1, 2013 were randomly selected for testing. The length of this baseline is about 190
km. The hardware types of the receivers and antennas are given in Table 6.2. The number of
GPS, GLONASS and GPS+GLONASS satellites is shown in Figure 6.9 as blue, green and red
line, respectively. FFMJ was selected as the reference station. The station TITZ was processed
by the HALO_GNSS software using kinematic processing. The results were compared with
the coordinates taken from the IGS website (http://www.igs.org), which were regarded as
“truth”. The following four schemes were compared.
- 106 -
Table 6.2: Hardware list of the selected IGS sites
Station Name
Receiver type
Antenna type ( with dome)
TITZ
JPS LEGACY
LEIAR25.R4 LEIT
FFMJ
JPS LEGACY
LEIAR25.R4 LEIT
Figure 6.9: The number of selected satellites of GPS (blue line), GLONASS (green line) and
GPS+GLONASS (red line) for the static experiment (IGS station TITZ and FFMJ on January 1, 2013)
Scheme 1: The GPS system was used alone to calculate the position of TITZ.
Scheme 2: The GLONASS system alone was used to calculate the position of TITZ.
Scheme 3: The GPS and GLONASS integrated system was used to calculate the position
of TITZ with 1:1 weights.
Scheme 4: The GPS and GLONASS integrated system was used to calculate the position
of TITZ with Helmert’s VCE.
The differences between the “true values” and the results of the schemes 1-4 are shown in
Figures 6.10, 6.11, 6.12 and 6.13 as a time series, respectively. The statistics of comparisons
are given in Table 6.3.
- 107 -
Figure 6.10: The differences between the IGS result and the GPS kinematic positioning results for the
baseline TITZ – FFMJ (scheme 1)
Figure 6.11: The differences between the IGS result and the GLONASS kinematic positioning results
for the baseline TITZ – FFMJ (scheme 2)
- 108 -
Figure 6.12: The differences between the IGS result and the GPS+GLONASS (with 1:1 weights) kine-
matic positioning results for the baseline TITZ – FFMJ (scheme 3)
Figure 6.13: The differences between the IGS result and the GPS+GLONASS (with Helmert weights)
kinematic positioning results for the baseline TITZ – FFMJ (scheme 4)
These experiments based on kinematic processing of static data illustrate the effectiveness,
repeatability and stability of the proposed method in static data. The accuracy of GNSS kine-
matic positioning based on GLONASS alone is a little bit lower than the GPS-only accuracy.
The combination of both GPS and GLONASS is better than their use as single systems, and
improves both the accuracy and reliability. The method based on Helmert’s VCE can further
improve the accuracy of the integrated GNSS kinematic positioning. The ratios between the
weights for the GPS and GLONASS observations are shown in Figure 6.14.
- 109 -
Figure 6.14: The ratios between the weights for the GPS and GLONASS observations (scheme 4)
Table 6.3: Statistical results for the IGS station TITZ (Unit: mm)
Scheme
Directions
Min
Max
Mean
RMS
1
G(GPS)
North
-29.8
34.8
6.6
10.0
East
-32.3
17.0
-9.8
7.1
Up
-97.4
71.7
-18.3
21.6
2
R(GLONASS)
North
-19.6
92.0
9.4
10.1
East
-41.8
47.6
9.7
13.3
Up
-111.5
82.2
-19.2
26.8
3
G+R
(1:1 weights)
North
-20.1
55.8
8.0
7.7
East
-39.0
11.0
9.8
6.1
Up
-73.5
87.3
-21.6
17.9
4
G+R
(Helmert weighs)
North
-25.9
23.0
0.0
7.2
East
-20.2
21.9
0.1
6.0
Up
-51.5
62.8
0.9
17.6
- 110 -
6.5 Conclusions
In this chapter it has been shown that the application of the GNSS integration can improve the
accuracy and reliability of GNSS positioning. Since the Galileo system and the BeiDou Navi-
gation Satellite System (BDS) are currently still under construction, only the GPS and
GLONASS systems are used in this study for a multi-system combination. The exercise of
mathematical models and processing methodologies for integrated systems is a valuable topic
as it identifies crucial issues concerning the combination of two or more GNSS positioning
systems. Hence these are experiences that can be applied to the other GNSS systems that
might integrate GPS with Galileo, GLONASS, BDS, or all four. Here it has been shown that
the accuracy of GNSS kinematic positioning based only on GLONASS as a single system is a
little bit worse than for GPS-only. But the combination of GPS and GLONASS is better than
using them as single systems and improves the accuracy and reliability. Finally it has been
shown that Helmert’s VCE method can be applied to estimate the weights of the multiple
GNSS observation data and can further improve the accuracy of integrated GNSS kinematic
positioning.
- 111 -
7 GNSS Velocity Determination Based on the Doppler
Effect
7.1 Introduction
Accurate estimates of the velocity of a platform are often needed in high dynamics positioning,
airborne gravimetry, GNSS/Inertial navigation System (INS) integration and geophysical ex-
ploration, etc (Szarmes et al., 1997; Bruton et al., 1999). GNSS is a cost-effective means to
obtain reliable and highly accurate velocities by exploiting the receiver raw and carrier phase
derived Doppler measurements (Szarmes et al., 1997). GNSS velocity determination can be
based on the Doppler effect. Its principle was first proposed by the Austrian mathematician
and physicist Christian Andreas Doppler (1842).
In the present chapter, firstly the principle of GNSS raw Doppler velocity determination is
discussed, and subsequently the velocity determination by carrier phase derived Doppler ob-
servations is described as well. In order to improve the accuracy and the reliability of the
velocity determination, the algorithms of GNSS integration and robust estimation are applied
as well. Finally, the actually measured GNSS data of the kinematic and static experiments are
applied to validate the presented methods and the importance of precise velocity determination
for airborne gravimetry is discussed.
7.2 Principle of Doppler-Based Velocity Determination
When a wave source is moving relative to an observer the perceived wave frequency is differ-
ent from the emitted frequency (Jonkman, 1980). This effect is named Doppler effect and has
been widely used in velocity determination.
The estimation of the velocity from discrete time signals in GNSS is based on the receiver
generated raw Doppler measurements and Doppler observations derived from the carrier phase
measurements (Cannon et al., 1997; Bruton et al., 1999). The receiver generated raw Doppler
measurement is a measure of near instantaneous velocity, whereas the carrier phase derived
Doppler observation is a measure of the mean velocity between observation epochs. The raw
measurements are usually noisier than the carrier phase derived Doppler measurement, be-
cause they are determined over a very small time interval (Serrano et al., 2004a). Since the
carrier phase derived Doppler shift is computed over a longer time span than the raw one, the
random noise is averaged and suppressed. Therefore, the very smooth velocity estimates in
- 112 -
low dynamic environments can be obtained from carrier phase derived Doppler measurements
if there are no undetected cycle slips (Serrano et al., 2004b).
In this section, the principle of GNSS Doppler velocity determination will be described
based on raw Doppler measurements generated by the GNSS receiver and Doppler measure-
ments derived from the GNSS carrier phase measurements.
7.2.1 Velocity Determination Based on GNSS Raw Doppler Observations
The GNSS raw Doppler shift, which is caused by the relative motion between the receiver and
the satellite, is the measurement of the phase rate (Cannon et al., 1997) directly estimated from
the Phase Lock Loop (PLL) output (Szarmes et al., 1997; Borio et al., 2009). In a first approx-
imation, the Doppler shift
D
between the GNSS satellite
s
and receiver
r
at the frequency
channel
j
can be written as (cf. Eq. (2.3) )
,,
ss
rr
s s s s
r j j r j j s
j
VV
D f f f
c
, (7.1)
where
f
denotes the frequency of the GNSS carrier phase observation.
s
r
V
is the radial ve-
locity of the range between the satellite
s
and the receiver
r
.
c
denotes the speed of light in
vacuum and
is the wavelength.
,
s
rj
D
has a positive sign when the receiver and the transmit-
ter approach each other and a negative sign when they move away from each other. Eq. (7.1)
for the observed Doppler shift scaled to range rate is given by
,()
s
r
s s s s
j r j r r
V D c dt dt
, (7.2)
where the derivatives with respect to time are indicated by a dot.
stands for the geometric
range rate between the satellite and receiver.
r
dt
and
s
dt
denote the receiver clock drift and
satellite clock drift, respectively.
is the effect of the observational noise and all non-
modeled error sources, such as errors in multipath.
Most reference books suggest that the effects with low frequency properties such as iono-
sphere, troposphere, tide and multipath effects are so small that they can be neglected.
However, according to Simsky and Boon (2003) and Zhang (2007), these errors should be tak-
en into account in the precise velocity determination. It is well known that the GNSS
observables have errors and biases that are spatially and temporally correlated. Differential
GNSS techniques provide an efficient solution to account for the correlated errors and biases,
since many un-modelled errors may be significantly reduced by the differencing process
(Zhang et al., 2006). Therefore, similarly to the Double Difference (DD) GNSS positioning (cf.
- 113 -
Section 4.2), the DD GNSS velocity determination is applied in this study for precise velocity
estimation.
Assuming the radial velocity equation of the range of the reference station
r
and the kine-
matic station
k
to a common satellite
q
can be written as
,()
q q q q q
r r j r j r r
e V D c dt dt
(7.3)
and
,()
q q q q q
k k j k j k k
e V D c dt dt
, (7.4)
where
e
is the direction cosine of the satellite and receiver.
q
r
V
denotes the radial velocity of
the range between the satellite
q
and the reference station
r
with
2 2 2
( ) ( ) ( )
q q q q q
r r r r r
V V V x x y y z z
, (7.5)
where
,,
q q q
x y z
denotes the velocity vector of the satellite
q
and
,,
r r r
x y z
denotes the
velocity vector of the reference station
r
.
Forming a single difference, the satellite clock drift
q
dt
can be cancelled out by
, , ,
( ) ( )
q q q q q q q
k k r r j k j r j k r k r
e V e V D D c dt dt
. (7.6)
The single difference for another satellite
p
at the signal frequency
i
is
, , ,
( ) ( )
p p p p p p p
k k r r i k i r i k r k r
e V e V D D c dt dt
. (7.7)
The DD equation for two single difference equations of the satellites
p
and
q
can be writ-
ten as
,
, , , , ,
( ) ( )
p p p p q q q q p p q q p q
k k r r k k r r i k i r i j k j r j k r
e V e V e V e V D D D D
, (7.8)
where the differential receiver clock drift cancels. Since
pp
kk
V V V
,
pp
rr
V V V
,
qq
kk
V V V
and
qq
rr
V V V
, Eq. (7.8) can be rewritten as done by Wang and Xu (2011)
by
, , , ,
, , , ,
p q p q p p q q p q p q
k k r r k r k r k r k r
e V e V e V e V D
, (7.9)
where
,
, , , , ,
( ) ( )
p q p p q q
k r i k i r i j k j r j
D D D D D
and
,p q p q
k k k
e e e
etc. The velocity of the
reference station
r
V
is zero, the velocity of the satellites
p
V
and
q
V
can be obtained from the
GNSS precise ephemerides or broadcast ephemerides with an analytical differentiation of the
position parameter equations. In order to achieve a solution at the mm/s level, satellite posi-
- 114 -
tions have to be known better than 10 m (Serrano et al., 2004a). Therefore, the velocity DD
observation equation for the kinematic station
k
can be given as
, , ,
, , , ,
p q p q p p q q p q
k k k r k r k r k r
e V D e V e V
. (7.10)
The precise velocity estimation of the kinematic station can be directly obtained by the
classic LS adjustment if Doppler observations from more than four GNSS satellites have been
measured (cf. Section 3.2.1).
7.2.2 Velocity Determination Based on Doppler Observations Derived from GNSS
Carrier Phase Measurements
Every GNSS receiver measures Doppler shifts. However, this is primarily used as an interme-
diate process to obtain accurate carrier phase measurements. Thus, some receivers output raw
Doppler shift measurements and some do not output any at all. In the absence of raw Doppler
shift, the carrier phase derived Doppler shift can be used (Zhang et al., 2005). It can be ob-
tained either by differencing carrier phase observations in the time domain, normalizing them
with the time interval of the differenced observations or by fitting a curve to successive phase
measurements, using polynomials of various orders (Serrano et al., 2004a).
At present, the first order central difference is one of the most popular methods for obtain-
ing the virtual Doppler observation (Cannon et al., 1997; Bruton et al., 1999; Serrano et al.,
2004a; Wang and Xu, 2011). Based on the fundamental GNSS carrier phase observation equa-
tion, cf. Eq. (2.2), the carrier phase derived Doppler observation is obtained by
1
22
t t t t t t t t t t
tt t t
, (7.11)
where
t
is the observation time and
t
is the data sample interval.
t
is the variation rate of
the original carrier phase observation
t
, namely the carrier phase derived Doppler observa-
tion.
The resulting observation equation for the velocity determination is
, , .
()
q q q q q q q
j r j r r r r r j r j
e V c dt dt T I
, (7.12)
where
T
and
I
denote the tropospheric and the ionospheric delay rate, respectively.
repre-
sents un-modeled effects, modeling error and measurement error rate for carrier phase
observations. Since the DD velocity determination is used, the error rates of
T
,
I
and
e
can
be summarised as
,.
q q q q
r r r j r j
u T I
(Serrano et al., 2004b).
- 115 -
The raw Doppler velocity DD observation equation Eq. (7.10) can be replaced by the carri-
er phase derived Doppler observation equation
, , ,
, , , ,
p q p q p p q q p q
k k k r k r k r k r
e V e V e V u
. (7.13)
The precise velocity estimation of the kinematic station can be directly obtained by the
classic LS adjustment if the carrier phase derived Doppler observations from more than four
GNSS satellites are available (cf. Section 3.2.1).
The carrier phase derived Doppler observation is an average variation rate of the carrier
phase observation
during a time interval of
2t
, rather than the instantaneous variation
rate of
at the epoch time
t
. However, the functional model for the raw Doppler observation
is stricter than that of the carrier phase derived Doppler measurement. Therefore, the velocity
results of the raw Doppler observations are more reliable than those of carrier phase derived
Doppler observations in the high dynamic environments.
7.3 Velocity Determination Based on GNSS Integration and Robust Esti-
mation
In order to improve the accuracy and the reliability of the velocity determination, the algo-
rithms of GNSS integration based on Helmert’s variance component estimation (VCE) and
robust estimation are applied in GNSS Doppler velocity determination as well.
The accuracy, availability and reliability of GNSS Doppler velocity determination are very
dependent on the number of GNSS Doppler observations, especially in high dynamic envi-
ronments. A possible way to increase the number of observations is to use GNSS multi-
constellation observations. Therefore the GNSS integration is applied in GNSS Doppler veloc-
ity determination as well. The principle of GNSS integration based on Helmert’s VCE has
been discussed in Chapter 6. The performance and the comparison of the application of this
principle in GNSS Doppler velocity determination will be given in detail here.
The carrier phase derived Doppler observation is a measurement of the mean velocity be-
tween observation epochs which means that their random noise is averaged and lowered.
However, the performance of such Doppler data is strongly dependent on the quality of the
carrier phase observation. Cycle slips and gaps are unavoidable in GNSS processing. Usually,
cycle slips can be detected and edited, but sometimes the undetected cycle slips still exist. Cy-
cle slip detection in kinematic applications is difficult because of the increased dynamics
induced by the platform motion of the remote antenna (Cannon et al., 1997). In particular, the
occurrences of carrier cycle slips and gaps cause a significant degradation of the accuracy of
the GNSS application. When undetected cycle slips exist in GNSS observations, outliers will
- 116 -
appear in carrier phase derived Doppler observations, resulting in incorrect velocity determi-
nation. Therefore, the outliers must be eliminated. To minimize the impact of observation
outliers on the estimates, the robust estimation approach is applied here for the GNSS velocity
determination. The principle of robust estimation is discussed in Section 4.4. The performance
and the comparison of the application of robust estimation in GNSS carrier phase derived
Doppler velocity determination will be given in detail.
7.4 Experiment and Analysis
After the implementation of the algorithms for GNSS Doppler velocity determination, some
practical tests were given carried out for a static as well as a kinematic case. Thereby, actually
measured GNSS data of the GEOHALO project were used to validate the developed method
in this study and the importance of a precise velocity determination for airborne gravimetry is
discussed.
7.4.1 Static Experiment
In order to test the accuracy and reliability of the algorithms and software in a static data ex-
periment, the static GNSS observation data of the reference stations REF6 and RENO from
the GEOHALO flight on June 6, 2012 were randomly selected for testing. The distance of this
baseline is about 605 km. The hardware types of all receivers and antennas are the same as in
the experiment of Section 4.3.2 (cf. Table 4.2). The selected data contain GPS and GLONASS
observations with a sampling rate of 1 Hz. REF6 was selected as the reference station. The
station RENO was processed by the HALO_GNSS software as in a low dynamic environment.
Since this test was carried out for a static station, the velocity results were compared with the
expected zero velocity as the truth. The following five schemes were compared.
Scheme 1: The raw Doppler observations of GPS single system were used.
Scheme 2: The raw Doppler observations of GLONASS single system were used.
Scheme 3: The raw Doppler observations of GPS and GLONASS integration were used.
The robust estimation is applied to remove outliers.
Scheme 4: The carrier phase derived Doppler observations of GPS and GLONASS integra-
tion were used.
Scheme 5: The carrier phase derived Doppler observations of the integrated GPS and
GLONASS were used. The robust estimation is applied to remove outliers.
The comparative results are given in Figures 1-5 and Table 7.1.
- 117 -
Figure 7.1: The velocity from raw GPS Doppler observations at the static station RENO on June 6, 2012
(scheme 1)
Figure 7.2: The velocity from raw GLONASS Doppler observations at the static station RENO on June
6, 2012 (scheme 2)
- 118 -
Figure 7.3: The velocity from the integrated raw GPS and GLONASS Doppler observations at the static
station RENO on June 6, 2012 (scheme 3)
Figure 7.4: The velocity from the integrated GPS and GLONASS carrier phase derived Doppler
observations on June 6, 2012 for the static station RENO, without the usage of the robust estimation
(scheme 4)
- 119 -
Figure 7.5: The velocity from GPS and GLONASS integrated carrier phase derived Doppler observa-
tions on June 6, 2012 for the static station RENO, with the usage of the robust estimation (scheme 5)
Table 7.1: The statistical results of the velocity determination (static experiment) for the station RENO
in the schemata 1-5 (Unit: mm/s)
Scheme
Directions
Min
Max
Mean
RMS
1
Raw Doppler
G (GPS)
North
-37.3
30.3
0.3
7.1
East
-22.7
32.8
0.6
5.2
Up
-74.1
72.3
0.7
12.7
2
Raw Doppler
R (GLONASS)
North
-39.1
40.0
-0.1
6.6
East
-33.1
33.8
1.6
6.7
Up
-86.6
78.3
0.3
14.1
3
Raw Doppler
Robust estimation
G+R (Helmert weights)
North
-19.7
21.0
0.3
4.6
East
-18.6
21.7
-0.6
3.7
Up
-52.1
50.0
0.8
9.1
4
Carrier Phase
G+R (1:1 weights)
North
-276.3
171.0
-0.3
3.7
East
-83.6
360.5
0.1
5.0
Up
-485.0
328.5
0.0
6.7
5
Carrier Phase
Robust estimation
G+R (Helmert weights)
North
-3.1
2.7
0.2
0.6
East
-2.7
2.5
0.0
0.5
Up
-7.6
6.3
0.0
1.5
- 120 -
The described static experiment illustrates the effectiveness, repeatability and stability of
the proposed method using static data. The obtained results show the following:
(1) In the case of GPS-only, the accuracy of the velocity determination using the raw Dop-
pler data can achieve about 1 cm/s in the low dynamic environment.
(2) The accuracy of the velocity determination based on raw GLONASS-only Doppler da-
ta is a little bit lower than the GPS-only accuracy. The combination of GPS and
GLONASS based on Helmert’s VCE is better than their respective use as single sys-
tems and improves the accuracy and reliability.
(3) In the case of the velocity determination using the carrier phase derived Doppler data,
the undetected cycle slips influence considerably the velocity results. After the applica-
tion of the robust estimation, the velocity outliers caused by the undetected cycle slips
can be controlled as well.
(4) The highly accurate and very smooth velocity results in low dynamic environments can
be obtained by the carrier phase derived Doppler measurements if there are no unde-
tected cycle slips. The accuracy of the velocity determination for the carrier phase
derived Doppler data can achieve about 1 mm/s in low dynamic environments.
7.4.2 Kinematic Experiment
The kinematic GNSS data of the GEOHALO flight on June 6, 2012 were selected for this test.
These are the same data as used in Section 6.4. The same kinematic stations are selected. AIR5
and AIR6 are mounted on the front and back parts of the HALO aircraft, respectively (cf. Fig-
ure 4.5). The reference station REF6 was installed by GFZ next to the runway of the DLR
airport in Oberpfaffenhofen. The station RENO located in the Italy region was selected as the
second reference station because it includes raw Doppler observations. The hardware types of
all receivers and antennas are the same as in the experiment in Section 4.3.2 (cf. Table 4.2).
The selected data contain GPS and GLONASS observations, and the sampling rate is 1 Hz.
Figure 7.6 shows the HALO aircraft trajectory as latitude, longitude and height compo-
nents on June 6, 2012. The two significant humps on the height curve correspond to crossing
the Alps.
- 121 -
Figure 7.6: The components for the flight trajectory of the HALO aircraft on June 6, 2012
The velocity and the acceleration of this HALO aircraft flight are shown in Figure 7.7 and
Figure 7.8, respectively.
Figure 7.7: The velocity components for the flight of the HALO aircraft on June 6, 2012
- 122 -
Figure 7.8: The acceleration components for the flight of the HALO aircraft on June 6, 2012
The GNSS Doppler velocity determination was done using the HALO_GNSS software
which was developed at GFZ in the frame of this thesis. The IGS precise ephemeris (Kouba,
2009b) were chosen and the GPS and GLONASS satellite elevation angle was limited to
10
.
The number of GPS, GLONASS and GPS+GLONASS satellites is shown in Figure 7.9 as
blue, green and red line, respectively.
Figure 7.9: The number of selected satellites of GPS (blue line), GLONASS (green line) and
GPS+GLONASS (red line)
- 123 -
In order to investigate the capability of the elaborated approach based on the Doppler ve-
locity determination, the velocity results of the kinematic stations AIR5 and AIR6 were
obtained by the HALO_GNSS software using the multiple reference stations REF6 and RENO.
Therefore, five experimental schemes are designed.
Scheme 1: The raw Doppler observations of GPS single system were used.
Scheme 2: The raw Doppler observations of GLONASS single system were used.
Scheme 3: The raw Doppler observations of GPS and GLONASS integration were used.
The robust estimation is applied to remove outliers.
Scheme 4: The carrier phase derived Doppler observations of GPS and GLONASS integra-
tion were used.
Scheme 5: The carrier phase derived Doppler observations of GPS and GLONASS integra-
tion were used. The robust estimation is applied to remove outliers.
Since it is very difficult to obtain an absolutely true velocity value of the moving highly
dynamic platform, an internal comparison and a cross-validation are performed. The differ-
ences of the velocity estimates between two kinematic stations AIR5 and AIR6 are shown in
Figure 7.10. The jumps in Figure 7.10 are caused by the velocity and acceleration changes (cf.
Figure 7.7 and Figure 7.8) when a sudden change of the aircraft’s status occurred, such as
takeoff, landing, braking, turning and accelerating, etc.
Figure 7.10: The differences of the velocity estimates between AIR5 and AIR6
- 124 -
However, when the aircraft is flying straight and smoothly, the difference of the velocity
results between AIR5 and AIR6 should be theoretically zero. This assumption can be used to
check the internal precision of the velocity determination. Therefore, the velocity results for
the second straight and smoothly period of this flight from 6:40:00 am to 7:24:00 am (cf. Fig-
ure 7.7 and Figure 7.8) were selected for this evaluation. The comparison of the results
between AIR5 and AIR6 for the GPS-only scheme 1, the GLONASS-only scheme 2 and the
scheme 3 of the GPS and GLONASS integration are shown in Figure 7.11, Figure 7.12 and
Figure 7.13, respectively. The statistical results of these cases are given in Table 7.2.
Figure 7.11: The difference of the GPS-only velocity results between AIR5 and AIR6 in the period of a
straight smooth flight (scheme 1)
Figure 7.12: The difference of the GLONASS-only velocity results between AIR5 and AIR6 in the pe-
riod of a straight smooth flight (scheme 2)
- 125 -
Figure 7.13: The difference of the velocity results using the integrated GPS and GLONASS between
AIR5 and AIR6 in the period of a straight smooth flight (scheme 3)
Table 7.2: Statistical results for the velocity determination in the kinematic experiment for schemata 1-3
(Unit: mm/s)
Scheme
Directions
Min
Max
Mean
RMS
1
Raw Doppler
G(GPS)
North
-61.4
74.2
0. 9
12.2
East
-44.3
44.7
-1.1
9.8
Up
-113.2
87.6
1.6
23.5
2
Raw Doppler
R(GLONASS)
North
-50.3
49.8
0.3
13.5
East
-81.1
62.6
-0. 3
15.9
Up
-136.0
110.3
2.9
29.4
3
Raw Doppler
Robust estamition
G+R (Helmert weights)
North
-26.6
25.9
0. 3
7.6
East
-27.2
28.0
0.4
6.9
Up
-55.6
54.7
1.8
15.3
The aforementioned results show the following: The raw GNSS Doppler observations can
be used to estimate the precise velocity of a highly dynamic aircraft and the conclusion is the
same as in the recent literature. The precision of the velocity determination based only on the
GLONASS is a little bit worse than for GPS-only. The precision of the velocity determination
using the integrated GPS and GLONASS data based on Helmert’s VCE is better than taking a
single GNSS system.
- 126 -
Theoretically, the difference between different velocity solutions for the same kinematic
station should be zero as well. This assumption can be used to check the internal precision of
the velocity determination. For the kinematic station AIR5, the velocity results for the carrier
phase derived Doppler observations for scheme 4 and 5 were compared with its velocity re-
sults based on the raw Doppler data of scheme 3 (integrated GPS and GLONASS data). The
comparative results are shown in Figure 7.14 and Figure 7.15, respectively. The related statis-
tic results are given in Table 7.3.
Figure 7.14: The differences between the estimated velocities for the raw Doppler observations and for
the carrier phase derived Doppler at station AIR5 without the usage of the robust estimation (integrated
GPS and GLONASS data)
The comparison results show clearly the inability of the carrier phase derived Doppler ve-
locity determination in a high dynamic environment: There are many remaining jumps in
Figure 7.14 which might be caused by undetected cycle slips and high dynamic environments.
After the application of the robust estimation, the number of outliers due to undetected cycle
slips can be obviously reduced (cf. Figure 7.15). But the impact of the highly dynamic motion
still exists in the velocity determination of the carrier phase derived Doppler observations, be-
cause these data are average variation rates of the carrier phase observations
during a time
interval of
2t
which are more disturbed than the instantaneous variation rate of the raw
Doppler data of
at the epoch time
t
(cf. Section 7.2.2). The large spikes in Figure 7.15 cor-
respond to the sudden changes of the state of the aircraft (cf. Figure 7.7 and Figure 7.8).
Therefore, the velocity results for a high dynamic environment based on the raw Doppler
observations are more suitable than those based on carrier phase derived Doppler observations.
- 127 -
Figure 7.15: The differences between the estimated velocities for the raw Doppler observations and for
the carrier phase derived Doppler at station AIR5 by applying the robust estimation (integrated GPS and
GLONASS data)
Table 7.3: The statistical results of scheme 4 and 5 as plotted in figures 7.14 and 7.15 (Unit: mm/s)
Scheme
Directions
Min
Max
Mean
RMS
4
Carrier Phase
G+R (1:1 weights)
North
-212.1
367.2
0. 4
13.2
East
-381.4
330.9
-0. 6
14.9
Up
-606.9
438.2
0. 1
26.6
5
Carrier Phase
Robust estimation
G+R (Helmert weights)
North
-367.0
156.4
-0. 4
12.7
East
-330.7
216.0
0. 4
12.8
Up
-439.3
605.9
-0. 1
26.2
7.4.3 Importance of Precise Velocity Determination for Airborne Gravimetry
The precise velocity information is one of the most critical issues in the application of GNSS
results in airborne gravimetry. It is necessary for deriving the vertical component of the kine-
matic acceleration of the platform. This acceleration is used to separate the gravitational signal
from the disturbing kinematic vertical accelerations affecting the airborne platform (cf. Sec-
tion 1.2.1).
As seen from Table 7.2 and Table 7.3 the accuracy of the determined vertical velocity is
about 2 cm/s. In order to investigate whether this accuracy is sufficient for the above purpose,
- 128 -
the vertical accelerations should be computed from these velocity values. In order to illustrate
the problem, one track flown on the first day (June 6, 2012) of the GEOHALO mission was
chosen and the vertical acceleration was computed by numerical differentiation. The velocity
was calculated from raw Doppler observations (cf. Section 7.4.2). A part of the used velocity
data (second track of the named flight day) is given in the bottom plot in Figure 7.8. The re-
sulting kinematic vertical acceleration of the HALO aircraft is shown in Figure 7.16. The
result is given in units as usual in gravimetry
52
1mGal 10 m/s
, and seems to have a strong
noise component. Taking into account 1) that the resolution of the airborne gravity measure-
ments lies in the order of magnitude of 1 mGal, and 2) that the expected gravity signal along
such a track of 300-400 km length might vary maximally only in the range of several hundred
mGal, the oscillations between -15000 and +15000 mGal as seen in Figure 7.16 seem to imply
that such a result for the kinematic accelerations cannot be used at all. However, this conclu-
sion is overhasty. The accelerations measured by an airborne gravimeter itself show also such
strong oscillations. But when using a low-pass filter, these oscillations can be decomposed into
a feasible signal and a disturbing high-frequency noise.
Figure 7.16: Vertical kinematic acceleration of the HALO aircraft as calculated directly from GNSS
observations
Therefore, in order to obtain a feasible signal for the disturbing kinematic vertical accelera-
tion (Brozena et al., 1989; Bruton et al., 1999; Kreye and Hein, 2003), the same low-pass filter
should be applied to the vertical accelerations derived from the GNSS based vertical velocities.
For this a low-pass filter with a cut-off wavelength of 300 seconds has been applied in the fre-
quency domain (using Fast Fourier Transform). This cut-off wavelength corresponds to what
- 129 -
is typically used in aerogravimetry. The resulting filtered kinematic accelerations for the re-
spective GEOHALO track (produced by Dr. Franz Barthelmes, GFZ) are shown in Figure 7.17
as the black curve. The values are of an expected realistic order of magnitude (several mGal)
and are usable to separate the huge disturbing kinematic accelerations affecting the airborne
platform.
Figure 7.17: The black curve is the GNSS-based disturbing vertical kinematic accelerations of the HA-
LO aircraft after the application of a low-pass filter to the noisy result as displayed in Figure 7.16. The
red curve is a reference gravity signal computed from the global gravity field mode EIGEN-6C4 minus
normal-gravity. The green curve is the CHEKAN measure vertical acceleration with all corrections ex-
cept of that from GNSS.
In Figure 7.17, the black curve is the GNSS-based signal of the disturbing vertical kinemat-
ic accelerations of the HALO aircraft after the application of a low-pass filter to the noisy
result displayed in Figure 7.16. The green curve is the gravity signal of the CHEKAN gravim-
eter including all corrections but without that from GNSS. The red curve denotes the gravity
signal of the global gravity field model EIGN-6C4 minus normal-gravity, which is here used
to check the gravity results. The statistical results of the comparison between the CHEKAN
signal (green curve) and that of the reference gravity field model (red curve) are given in Table
7.4.
The final gravity signal at the HALO trajectory can be obtained by separating the GNSS-
based disturbing kinematic accelerations of the HALO aircraft from the CHEKAN measured
vertical acceleration. This is shown in Figure 7.18 as blue curve. The statistical results of the
comparison between final gravity signal and the reference gravity field model is given in Table
7.3. The comparison of Figure 7.17, Figure 7.18 and Table 7.4 shows the need and the im-
portance of GNSS for airborne gravimetry very significantly. This verification has been done
by the colleagues from GFZ who are dealing with airborne gravitmetry since the processing of
- 130 -
the CHEKAN airborne gravity data itself and the computation of the final gravity signal by
separation of the GNSS-based disturbing kinematic accelerations are not subject of this thesis.
Figure 7.18: The final gravity signal (blue curve) compared with the reference gravity field model
EIGEN-6C4 (red curve)
Table 7.4: The statistic results of the comparison between the calculated airborne gravity signal and the
reference gravity field mode EIGEN-6C4 (unit: mGal)
Gravity signal type
Min
Max
Mean
RMS
CHEKAN measurements gravity
(without GNSS modify)
-15.20
12.54
0.07
6.35
Final gravity (with GNSS modify)
-6.88
4.11
0.08
1.93
7.5 Conclusions
This chapter describes the principles of a series of options for the differential GNSS high accu-
racy velocity determination. Velocity was determined using Doppler observations obtained
either by using the receiver generated raw Doppler shifts or by taking the carrier phase derived
Doppler observations.
The best velocity results for low dynamic environments were obtained by applying the car-
rier phase derived Doppler observations. The velocity results for the raw Doppler observations
are more reliable than those of the carrier phase derived Doppler measurements in high dy-
- 131 -
namic environments, since the carrier phase derived Doppler observation is a measure of the
mean velocity between observation epochs. But they are usually noisier than the velocity re-
sults of the carrier phase derived Doppler observations because they are determined over a
very short time interval.
The accuracy, availability and reliability of GNSS Doppler velocity determination are
strongly dependent on the number of GNSS Doppler observations as well, especially in high
dynamic environments. Therefore, the GNSS integration based on Helmert’s VCE is applied in
GNSS Doppler velocity determination. The accuracies obtained from both the kinematic and
the static experiments have been improved in this way.
Cycle slips and gaps are unavoidable in any application of GNSS. Undetected cycle slips
usually exist in GNSS observations and cause outliers in the carrier phase derived Doppler
observations. This gives incorrect results for the GNSS velocity determination. Applying the
robust estimation, the outliers can be removed.
- 133 -
8 Software Development and Application
8.1 Introduction
In order to fulfill the actual requirements of airborne as well as shipborne gravimetry on
GNSS precise position and velocity determination, a software system (HALO_GNSS) for pre-
cise kinematic GNSS trajectory and velocity determination for kinematic platforms has been
developed. In this software, the algorithms as proposed in this thesis were implemented. In
order to evaluate the effectiveness of the proposed algorithms and the HALO_GNSS software,
this software was applied both in airborne and shipborne gravimetry projects of GFZ Potsdam.
In this chapter, the characteristics, the overall design and the structure of the HALO_GNSS
software are outlined. The applications of the kinematic GNSS algorithm and software are in-
troduced briefly.
8.2 HALO_GNSS Software Development
The HALO_GNSS software is a high-precision GNSS kinematic position and velocity deter-
mination software for scientific research, developed by the author at GFZ during his study on
this thesis. The purpose of this software is to test and evaluate the algorithms as introduced in
the previous chapters and to process the GNSS data for airborne and shipborne gravimetry
projects. All results presented in this thesis were obtained by using this HALO_GNSS soft-
ware.
The HALO_GNSS software has been developed in the programming language FORTRAN
under a Linux operating system by applying the algorithms as described in this thesis and by
learning a lot from the previous HALO_GPS software, which is a GPS single system position-
ing software for single baselines, developed by Wang et al. (2010) at GFZ too. Furthermore,
the development of HALO_GNSS was inspired by the famous software MFGsoft (Xu, 2004)
and some subroutines were adopted from the GAMIT software (Herring et al., 2010).
8.2.1 Characteristics of the HALO_GNSS Software
The most important characteristics of the HALO_GNSS software are:
(1) Multiple GNSS data can be processed. Currently, the GPS-only, GLONASS-only and
GPS/GLONASS integration are available.
- 134 -
(2) Precise position information of a kinematic platform can be obtained. It allows the kin-
ematic positioning with an accuracy of 1-2 cm.
(3) Precise velocity information of a kinematic platform can be obtained as well. It allows
the velocity determination with the accuracy of about 1 cm/s using raw and about 1
mm/s using carrier phase derived Doppler observations. This accuracy seems to be far
too low for the direct application in airborne gravimetry. However, the major error
component is the high-frequency noise, so that a useful signal can be recovered by low-
pass filtering.
(4) Multiple reference stations and multiple kinematic stations can be processed together.
(5) Robust estimation is applied to control the influence of outliers in GNSS position and
velocity determination.
(6) An a priori distance constraint among multiple kinematic stations on the same platform
is used to improve the accuracy and reliability of GNSS kinematic positioning.
(7) A common troposphere zenith delay parameter for multiple kinematic stations mounted
on the same platform is employed to improve the accuracy and reliability of GNSS
kinematic positioning as well.
(8) Helmert’s variance components estimation is used to reasonably adjust the weights and
the contributions among multiple GNSS systems.
(9) This software can be applied in many fields, such as static and kinematic environments,
airborne and shipborne platforms, small and large regions, and in particular in ultra-
high-altitude, ultra-long-range and high dynamic environments.
8.2.2 Structure Design of Software
Most of GPS software systems consist generally of three basic components (e.g. Xu (2004);
Herring et al. (2010); Wang et al. (2010)): Firstly, the functional library provides all possibly
needed physical models, algorithms and tools for use. Secondly, the data platform prepares all
possibly needed data for use. Finally, the data processing core forms the observation equations,
accumulates them and solves the problem (Xu, 2007, p. 219).
The HALO_GNSS software developed as part of this study was designed in accordance
with this general structure. The specific structure design of the HALO_GNSS software is giv-
en in Figure 8.1.
The main steps outlined in Figure 8.1 are:
(1) Start: Running the command “./HALO_GNSS”.
(2) Reading user command file: The most important parameters for controlling the run-
ning of the software are set.
- 135 -
Figure 8.1: The GNSS data processing flowchart of HALO_GNSS software
- 136 -
(3) Reading GNSS satellite orbit and clock files, e.g. for the GPS and GLONASS (cf.
Section 1.2.3).
(4) Read sites GNSS observation files, which include the GNSS data of reference stations
and kinematic stations, GPS and GLONASS, code, carrier phase and Doppler observa-
tions, etc. (cf. Section 2.2).
(5) Read other files if necessary, such as the Differential Code Biases (DCB) corrections,
antenna Phase Centre Offsets (PCO) and Phase Centre Variations (PCV), etc. (cf. Sec-
tion 2.4).
(6) Observation error corrected, reduce the errors resp. disturbances contained in the
GNSS measurements if possible, e.g. antenna offset, atmospheric delay and DCB cor-
rections, etc. (cf. Section2.4).
(7) Code positioning: Data pre-processing for the elevation angles of the satellites, detec-
tion of bad data and initial position of sites (cf. Section 4.2).
(8) Cycle slips detection and Ambiguity initialization which is necessary only for the
carrier phase data processing.
(9) Choice processing mode: Depended on the user command file, the processing mode
will be chosen, such as position or velocity determination, initialization of matrix and
all unknown parameters (cf. Chapter 3-7).
(10) Form the observation mode and/or system state mode, e.g. the DD observation
equations of carrier phase, pseudorange or Doppler observations; system state mode
(needed in Kalman filter only) (cf. Section 4.2, 4.3, 5.2, 5.3, 5.4, 6.2, 7.2 and 7.3).
(11) Form stochastic mode. e.g. the covariance matrix of observations, the covariance ma-
trix of Kalman filter system noise, and the covariance matrix of predicted state vector
(cf. Section 4.4, 6.3 and 7.3).
(12) Parameter estimation (cf. Section 3.2 and 3.3).
(13) Robust estimation to control the outliers if necessary (cf. Section 4.4 and 7.3).
(14) Helmert’s variance components estimation for multiple GNSS systems if necessary
(cf. Section 6.3 and 7.3).
(15) Residual edit. If necessary, set New cycle slips for the carrier phase observation de-
pend on the value of observation residuals.
(16) Output results and summary files.
(17) End of program.
- 137 -
The Operation instructions and the detailed description of the major functional modules of
the HALO_GNSS software will be given in a separate technical report “HALO_GNSS soft-
ware user manual”, which is not included in this thesis.
8.3 HALO_GNSS Software Applications
The HALO_GNSS software can be applied in different fields, such as static and kinematic
environments, airborne and shipborne platforms, small and large regions, in particular in ultra-
high-altitude, ultra-long-range and high dynamic environments. The applications for the
GNSS kinematic position and velocity determination using the HALO_GNSS software are
introduced briefly in this section; some of them have been already presented in the former
chapters.
8.3.1 Lake Müritz Shipborne Gravimetry
This experiment was performed by GFZ section 1.2 from October 25 to 27, 2011 at Lake
Müritz in Mecklenburg-Vorpommern, northern Germany (cf. Figure 8.2). One of the main
purposes of this survey was testing the newly purchased gravimeter CHEKAN. The other pur-
pose was to study the precise GNSS kinematic positioning for gravimetry.
Figure 8.2: The Lake Müritz location and the tracks of the ship
- 138 -
During the experiment, three days of measurements were performed on this lake. Seven
dual-frequency GNSS receivers were used, two of them were set up on the shore as reference
stations and another five were mounted on the ship, see Figure 8.3. The HALO_GNSS soft-
ware has been applied for GNSS data processing. The tracks of the ship on Lake Müritz are
shown in Figure 8.2 as green lines.
Figure 8.3: The ship used and the position of GNSS receiving equipments
8.3.2 The GEOHALO Italy Mission
On June 6, 8, 11 and 12, 2012, the GEOHALO mission was carried out over the Italian area
and adjacent regions of the Mediterranean, starting and ending at the DLR in Oberpfaffenho-
fen Germany (cf. Figure 8.4).
More details about the GEOHALO aircraft flight mission are given in Section 1.2.2. The
main goals of the GEOHALO mission were: GNSS reflectometry (Semmling et al., 2013;
Semmling et al., 2014), airborne gravimetry (Petrovic et al., 2013), magnetometry and laser
altimetry (Scheinert et al., 2013). GNSS precise position and velocity determination plays an
important role in the whole mission, since all remote sensing instruments and measuring de-
vices are fixed on the same kinematic platform. Therefore, the accurate state information of
the kinematic platform is a key factor of the successful implementation of the entire mission.
The HALO_GNSS software was used to provide the precise position and velocity infor-
mation of the HALO aircraft (cf. Figure 5.1) for the whole mission. The flown HALO tracks
(red lines) are shown in Figure 8.4. Some of GNSS data of the GEOHALO-flight have been
compared and analyzed in detail in sections 4.3.2, 4.4.2, 6.4.1 and 7.4.1 in order to demon-
strate the approach proposed in this thesis.
AIR2
AIR3
AIR1
AIR4 and AIR5
- 139 -
Figure 8.4: The GEOHALO mission tracks (red lines) and GNSS reference stations (red and green dots)
in Italy
8.3.3 Lake Constance Shipborne Gravimetry
The Lake Constance (Bodensee in German) campaign was performed by GFZ and BKG from
October 23 to 26, 2012. This lake is located in Germany, Switzerland and Austria near the
Alps (Figure 8.5). One of the main purposes of this survey was to obtain gravity information
of the Lake Constance area. The other purpose was to study the precise GNSS kinematic posi-
tioning for gravimetry.
During this campaign, four days of measurements were accomplished on the lake and the
HALO_GNSS software has been applied for the GNSS data processing. The tracks of the ship
run on Lake Constance are shown in Figure 8.5 as green lines.
- 140 -
Figure 8.5: All tracks of the ship in the Constance Lake in October 2012
8.3.4 Baltic Sea Shipborne Gravimetry
From June 18 to 27, 2013, ten days of GNSS and gravimetric data of a shipborne gravimetric
campaign were collected by GFZ and BKG in Baltic Sea (Ostsee in German) near Greifswald,
Germany (cf. Figure 8.6). The main purpose of this survey was to obtain gravity information
of the Baltic Sea near Greifswald, Germany and to study the precise GNSS kinematic posi-
tioning for gravimetry.
The HALO_GNSS software was used to provide the precise position information of the
ship for the whole campaign. The tracks of the ship during this campaign are shown in Figure
8.6 (green lines). The GNSS data of the first day of the Baltic Sea gravimetric campaign have
been compared and analyzed in detail in Section 5.5 to demonstrate the approach proposed in
chapter 5.
- 141 -
Figure 8.6: All tracks of the ship in the Baltic Sea gravimetric campaign in June 2013
8.4 Conclusions
A high-precision kinematic GNSS position and velocity determination software HALO_GNSS
has been developed based on algorithms proposed in this thesis, to meet the actual require-
ments of airborne as well as shipborne gravimetry. The characteristics, overall design and the
structure of the HALO_GNSS software were outlined firstly. In order to evaluate the effec-
tiveness of the proposed algorithms and the HALO_GNSS software, this software has been
applied in several airborne and shipborne gravimetry projects of GFZ Potsdam. The applica-
tions of the kinematic GNSS algorithm and software are introduced briefly.
- 143 -
9 Summary and Future Work
9.1 Summary
The GNSS kinematic position and velocity determination plays an important role in airborne
gravimetry. Development, analysis and investigation of a reliable GNSS algorithms and soft-
ware for kinematic highly precise GNSS data analysis for airborne gravimetry are the major
goals of this research. It includes investigation of the estimation algorithms, accuracy evalua-
tion analysis, the development of the new methods of integrated GNSS kinematic position and
velocity determination based on multiple reference stations and multiple kinematic stations,
robust estimation, a priori distance constraints, a common tropospheric delay parameter,
Helmert’s variance components estimation, and software development and application. The
core results and the specific contributions of this research can be summarized as follows.
The parameter estimation algorithms of least squares including the elimination of nuisance
parameters as well as the two-way Kalman filter were applied in kinematic GNSS data post-
processing to fulfill the high accuracy requirements of airborne gravimetry.
A specific accuracy and the reliability evaluation method for GNSS kinematic precise posi-
tion determination has been applied in static and kinematic modes, for short as well as long
baselines, and for low and high dynamics states. The receiver clock errors including the clock
jumps were analyzed for highly dynamic GNSS precise positioning. It has been shown that the
receiver clock error including clock jumps must be corrected for each measurement epoch.
A new approach using multiple reference stations based on an a priori constraint is ad-
dressed to increase the accuracy and reliability of GNSS precise kinematic positioning for the
widely arranged GNSS stations. A robust estimation algorithm is used to suppress the impact
of observation outliers on the trajectory estimates.
A new method of GNSS kinematic positioning based on multiple kinematic stations on the
same platform with a priori distance constraints and the estimation of a common tropospheric
delay parameter is developed. In this new method, the distances among the multiple GNSS
antennas are known and used as a priori distance constraints to improve the accuracy and the
reliability of the state estimates. Since the characteristics of the tropospheric delays in a small
area are similar, a common tropospheric delay parameter has been set up for such multiple
kinematic stations, which can enhance the accuracy and reliability of the state estimates.
- 144 -
A kinematic positioning method using the integration of multiple GNSS systems based on
Helmert’s variance components estimation is addressed to adjust the weights in a reasonable
way in order to balance the contributions of multiple GNSS systems, and to improve the relia-
bility and accuracy of GNSS kinematic positioning.
Velocities were determined using Doppler observations. For this purpose two kinds of
Doppler data were used: 1) Receiver generated raw Doppler shift measurements or 2) Doppler
data derived from the carrier phase measurements. The best velocity results for low dynamics
were obtained by applying carrier phase derived Doppler observations. For high dynamic en-
vironments, the velocity results based on the raw Doppler shift measurements are more
reliable than those derived from carrier phase data. A new approach for velocity determination
based on the integration of multiple GNSS systems using Helmert’s variance components es-
timation and the robust adjustment has been developed, which can significantly improve the
accuracy of the results and control the influence of the outliers.
In order to fulfill the actual requirements of airborne and shipborne gravimetry on kinemat-
ic GNSS data analysis, a reliable and precise kinematic GNSS position and velocity
determination software system (HALO_GNSS) has been developed based on the algorithms
proposed in this thesis. This software can be applied for different purposes, such as static and
kinematic processing, airborne and shipborne platforms, small and large regions, in particular
in ultra-high-altitude, ultra-long-range and high dynamic environments. It allows the GNSS
kinematic positioning with an accuracy of 1-2 cm and GNSS velocity determination with the
accuracy of about 1 cm/s using raw and about 1 mm/s using carrier phase derived Doppler ob-
servations. The accuracy of these results for velocity is high enough in GNSS field, but it
cannot be used directly in airborne gravimetry. Because of the high-frequency characters of
the noise, a sufficiently accurate vertical acceleration results can be obtained by applying a
low-pass filter.
These contributions have been supported by numerical results and tested under different
conditions. All results presented in this thesis were obtained by using the HALO_GNSS soft-
ware.
9.2 Future Work
As a continuation of this research, a variety of points can be suggested. Among them are those
described in the following.
The integration of the raw Doppler velocity into the Kalman filter for an improved position
determination should be investigated. The precise instantaneous velocity of a kinematic plat-
form can be obtained by raw Doppler observation. This velocity information should be used to
- 145 -
improve the classic application of the Kalman filter. The observation mode, system mode and
a priori covariance of classic Kalman filter will be improved and enhanced, and will be more
reliable and accurate.
The methodology presented in this study can be in principle applied to any GNSS system,
such as GPS, GLONASS, BDS and Galileo, etc. However, here in this study resp. in the de-
veloped HALO_GNSS software only the two GNSS systems GPS, GLONASS and their
integration/combination were adopted. The exercise of mathematical models and processing
methodologies for integrated systems is valuable since it is capable to identify crucial issues
concerning the combination of any two or more GNSS positioning systems. Hence these expe-
riences can be applied to other GNSS systems in order to integrate GPS with Galileo,
GLONASS, BDS, or all the four together. This interesting topic requires further investigation
and will probably improve the performance of GNSS position and velocity determination.
Since real-time methods play a more and more important role in GNSS applications, the
precise GNSS kinematic position and velocity determination methods for large regions as
studied in this thesis should be adopted in the future for real-time applications.
The attitude determination of a kinematic platform is an important issue for many scientific
applications, like remote sensing and airborne gravimetry. The method of precise GNSS kine-
matic positioning based on multiple stations on the same platform as used in this thesis can be
further developed and used for attitude determination.
In the algorithm of multiple reference stations, the weights of the observations from each
reference station can be introduced depending on the lengths of the distance between the kin-
ematic aircraft and the reference stations.
For the future, it is planned to use the HALO aircraft for airborne gravimetry in central
Antarctica to close the “polar gap” of dedicated satellite gravity field missions (due to their
non-polar orbit inclinations). But for the polar region Antarctica it’s expected that the specific
constellation of the GNSS satellites and the sparse density of reference ground stations will
cause significant erroneous biases in the altitude determination. Therefore, the development of
new related algorithms is urgently needed to fulfill the requirements of HALO in the expected
situation. Furthermore, the HALO data of laser altimetry could be used to compare and com-
bine with the GNSS data.
- 147 -
References
Abdel-salam, M. A. 2005. Precise point positioning using un-differenced code and carrier
phase observations. Ph.D Thesis, University of Calgary.
Al-Shaery, A., Lim, S. & Rizos, C. 2011. Investigation of different interpolation models used
in Network-RTK for the virtual reference station technique. Journal of Global
Positioning Syetems, 10, 136-148
Alberts, B. 2009. Regional gravity field modeling using airborne gravimetry data, Delft, NCG,
Nederlandse Commissie voor Geodesie.
Angrisano, A., Gioia, C., Gaglione, S. & del Core, G. 2013. GNSS Reliability Testing in
Signal-Degraded Scenario. International Journal of Navigation and Observation,
2013, 1-12
Askne, J. & Nordius, H. 1987. Estimation of tropospheric delay for microwaves from surface
weather data. Radio Science, 22, 379-386
Bähr, H., Altamimi, Z. & Heck, B. 2007. Variance component estimation for combination of
terrestrial reference frames, Univ.-Verlag Karlsruhe.
Baarda, W. 1968. A testing procedure for use in geodetic networks. Publications on geodesy,
New Series, Vol. 2, No. 5, Delft
Baarda, W. 1976. Reliability and precision of networks. In: Eichhorn, G. (ed.) VIIth
International Course for Engineering Surveys of high precision. Darmstadt, Germany.
Banville, S., Santerre, R., Cocard, M. & Langley, R. B. 2008. Satellite and receiver phase bias
calibration for undifferenced ambiguity resolution. Proc. ION NTM.
BDS-ICD 2013. BeiDou Navigation Satellite System Signal In Space Interface Control
Docoument. Open Service Signal (Version 2.0). Beijing: China Satellite Navigation
Office.
BDS-OS-PS 2013. BeiDou Navigation Satellite System Open Service Performance Standard.
Version 1.0. Beijing: China Satellite Navigation Office.
BDS. 2014. BeiDou Navigation Satellite System [Online]. Available: http://en.beidou.gov.cn/
[Accessed 03.02 2014].
Beser, J. & Parkinson, B. W. 1982. The application of NAVSTAR differential GPS in the
civilian community. Navigation, 29, 107-136
Biacs, Z., Krakiwsky, E. & Lapucha, D. 1990. Reliability analysis of phase observations in
GPS baseline estimation. Journal of Surveying Engineering, 116, 204-224
Bisnath, S. & Gao, Y. 2008. Current state of precise point positioning and future prospects and
limitations. Observing our changing earth. Springer.
Black, H. & Eisner, A. 1984. Correcting satellite Doppler data for tropospheric effects.
Journal of Geophysical Research: Atmospheres (1984–2012), 89, 2616-2626
Blewitt, G. 1989. Carrier phase ambiguity resolution for the Global Positioning System applied
to geodetic baselines up to 2000 km. Journal of Geophysical Research: Solid Earth
(1978–2012), 94, 10187-10203
Blewitt, G. 1990. An automatic editing algorithm for GPS data. Geophysical Research Letters,
17, 199-202
Blewitt, G., Melbourne, W., Bertiger, W., Dixon, T., Kroger, P., Lichten, S., Meehan, T., Neilan,
R., Skrumeda, L. & Thornton, C. 1988. GPS geodesy with centimeter accuracy. GPS-
Techniques Applied to Geodesy and Surveying. Springer.
Bock, Y. & Shimada, S. 1990. Continuously monitoring GPS networks for deformation
measurements. Global Positioning System: An Overview. Springer.
Boehm, J. & Schuh, H. 2004. Vienna mapping functions in VLBI analyses. Geophysical
Research Letters, 31
- 148 -
Bohm, J., Niell, A., Tregoning, P. & Schuh, H. 2006. Global Mapping Function (GMF): A new
empirical mapping function based on numerical weather model data. Geophysical
Research Letters, 33
Borio, D., Sokolova, N. & Lachapelle, G. 2009. Doppler measurements and velocity
estimation: A theoretical framework with software receiver implementation.
ION/GNSS. Savannah, Georgia.
Bossler, J. D., Goad, C. C. & Bender, P. L. 1980. Using the Global Positioning System (GPS)
for geodetic positioning. Bulletin géodesique, 54, 553-563
Brozena, J., Mader, G. & Peters, M. 1989. Interferometric Global Positioning System: Three‐
dimensional positioning source for airborne gravimetry. Journal of Geophysical
Research: Solid Earth (1978–2012), 94, 12153-12162
Bruton, A. M., Glennie, C. L. & Schwarz, K. P. 1999. Differentiation for high-precision GPS
velocity and acceleration determination. GPS solutions, 2, 7-21
Cannon, M. E., Lachapelle, G., Szarmes, M. C., Hebert, J. M., Keith, J. & Jokerst, S. 1997.
DGPS kinematic carrier phase signal simulation analysis for precise velocity and
position determination. Navigation, 44, 231-245
Chen, G. 1998. GPS kinematic positioning for the airborne laser altimetry at Long Valley,
California. Ph.D Thesis, Massachusetts Institute of Technology.
Chen, G. & Herring, T. A. 1997. Effects of atmospheric azimuthal asymmetry on the analysis
of space geodetic data. Journal of Geophysical Research: Solid Earth (1978–2012),
102, 20489-20502
Chen, W., Hu, C., Chen, Y., Ding, X. & Kwok, S. 2001. Rapid static and kinematic positioning
with Hong Kong GPS active network. ION GPS.
Chen, W., Hu, C., Gao, S., Chen, Y., Ding, X. & Simon, C.-w. K. 2004. Absolute ionospheric
delay estimation based on GPS PPP and GPS active network. 2004 International
Symposium on GNSS/GPS. Sydney, Australia.
Chen, X., Han, S., Rizos, C. & Goh, P. C. 2000. Improving real time positioning efficiency
using the Singapore integrated multiple reference station network (SIMRSN).
Proceedings of the 13th International Technical Meeting of the Satellite Division of
The Institute of Navigation (ION GPS 2000).
Chu, F.-Y. & Yang, M. 2013. GPS/Galileo long baseline computation: method and
performance analyses. GPS Solutions, 18, 1-10
Collins, J. P. & Langley, R. B. 1997. Estimating the residual tropospheric delay for airborne
differential GPS positioning. PROC ION GPS, INST OF NAVIGATION,
ALEXANDRIA, VA,(USA), 1997, 2, 1197-1206
Collins, P., Gao, Y., Lahaye, F., Héroux, P., MacLeod, K. & Chen, K. 2005. Accessing and
processing realtime GPS corrections for precise point positioning—some user
considerations. Proceedings of ION GNSS 18th international technical meeting of the
satellite division. Long Beach, CA, US.
Collins, P., Lahaye, F., Heroux, P. & Bisnath, S. 2008. Precise point positioning with
ambiguity resolution using the decoupled clock model. Fairfax, USA: ION GNSS
Crocetto, N., Gatti, M. & Russo, P. 2000. Simplified formulae for the BIQUE estimation of
variance components in disjunctive observation groups. Journal of Geodesy, 74, 447-
457
Dach, R., Hugentobler, U., Fridez, P. & Meindl, M. 2007. Bernese GPS software version 5.0.
Astronomical Institute, University of Bern.
Dai, L., Han, S., Wang, J. & Rizos, C. 2001. A study on GPS/GLONASS multiple reference
station techniques for precise real-time carrier phase-based positioning. ION GPS
2001. Salt Lake City, USA.
Dai, Z. 2013. On GPS based Attitude Determination. Ph. D, Universität Siegen.
Davis, J. L., Herring, T. A., Shapiro, I. I., Rogers, A. E. E. & Elgered, G. 1985. Geodesy by
radio interferometry: Effects of atmospheric modeling errors on estimates of baseline
length. Radio science, 20, 1593-1607
Deng, Z. 2012. GPS Meteorology with Single Frequency Receivers. Ph.D Thesis, University of
Hannover.
- 149 -
Dixon, K. 2006. StarFire(TM): A Global SBAS for Sub-Decimeter Precise Point Positioning.
ION GNSS 19th International Technical Meeting of the Satellite Divesion. Fort Worth,
TX, USA: Proceedings of ION GNSS 2006.
Dixon, T. H. 1991. An introduction to the Global Positioning System and some geological
applications. Reviews of Geophysics, 29, 249-276
DLR. 2013a. HALO aircraft [Online]. Available:
http://www.dlr.de/dlr/en/desktopdefault.aspx/tabid-10658#gallery/6163 [Accessed
03.03 2014].
DLR. 2013b. HALO Home [Online]. Available: http://www.halo.dlr.de/ [Accessed 01.12 2013].
Dong, D. N. & Bock, Y. 1989. Global Positioning System network analysis with phase
ambiguity resolution applied to crustal deformation studies in California. Journal of
Geophysical Research: Solid Earth (1978–2012), 94, 3949-3966
Doppler, C. 1842. Ueber das farbige Licht der Doppelsterne und einiger anderer Gestirne des
Himmels, Prag, Abh Kgl Bohm Ges Wissench.
Doran, H. E. 1992. Constraining Kalman filter and smoothing estimates to satisfy time-
varying restrictions. The Review of Economics and Statistics, 568-572
ESA. 2014. Galileo General Introduction [Online]. Available:
http://www.navipedia.net/index.php/GALILEO_General_Introduction [Accessed
03.04 2014].
Förstner, W. 1979. Ein verfahren zur schätzung von varianz-und kovarianzkomponenten.
Allgemeine Vermessungsnachrichten, 86, 446-453
Förstner, W. 1983. Reliability and discernability of extended Gauss-Markov models. Deutsche
Geodätische Kommission (DGK). München, Report A 98:79-104
Förstner, W. 1987. Reliability analysis of parameter estimation in linear models with
applications to mensuration problems in computer vision. Computer Vision, Graphics,
and Image Processing, 40, 273-310
Fernandes, M. J., Pires, N., Lazaro, C. & Nunes, A. L. 2013. Tropospheric delays from GNSS
for application in coastal altimetry. Advances in Space Research, 51, 1352-1368
Firuzabadì, D. & King, R. W. 2012. GPS precision as a function of session duration and
reference frame using multi-point software. GPS solutions, 16, 191-196
Forsberg, R. & Olesen, A. V. 2010. Airborne gravity field determination. In: Xu, G. (ed.)
Sciences of Geodesy-I. Heidelberg, London, New York: Springer.
Fotopoulos, G. & Cannon, M. E. 2001. An overview of multi-reference station methods for
cm-level positioning. GPS solutions, 4, 1-10
Galas, R., Schöne, T., Cokrlic, M. & Semmling, M. 2013. On precise GNSS-based sea surface
monitoring systems. ELMAR, 2013 55th International Symposium. IEEE.
Galileo-ICD 2010. European GNSS (Galileo) Open Service Signal in Space Interface Control
Document. 2010.
Ge, M., Douša, J., Li, X., Ramatschi, M., Nischan, T. & Wickert, J. 2012. A novel real-time
precise positioning service system: Global precise point positioning with regional
augmentation. J. Glob. Position. Syst, 11, 2-10
Ge, M., Gendt, G., Rothacher, M., Shi, C. & Liu, J. 2008. Resolution of GPS carrier-phase
ambiguities in precise point positioning (PPP) with daily observations. Journal of
Geodesy, 82, 389-399
Gelb, A. 1974. Applied optimal estimation, Cambridge, Massachusetts, and London, The MIT
press.
Geng, J., Teferle, F. N., Meng, X. & Dodson, A. 2010. Kinematic precise point positioning at
remote marine platforms. GPS solutions, 14, 343-350
Geng, J., Teferle, F. N., Meng, X. & Dodson, A. 2011. Towards PPP-RTK: Ambiguity
resolution in real-time precise point positioning. Advances in space research, 47,
1664-1673
Genrich, J. F. & Bock, Y. 1992. Rapid resolution of crustal motion at short ranges with the
Global Positioning System. Journal of Geophysical Research: Solid Earth (1978–
2012), 97, 3261-3269
- 150 -
GFZ. 2013. Terrestrial and Airborne Gravimetry [Online]. Available: www.gfz-
potsdam.de/en/research/organizational-units/departments/department-1/global-
geomonitoring-and-gravity-field/topics/terrestrial-and-airborne-gravimetry/ [Accessed
31.01 2014].
Ghilani, C. D. 2010. Adjustment computations: spatial data analysis, John Wiley & Sons.
Grafarend, E., Kleusberg, A. & Schaffrin, B. 1980. An introduction to the variance-covariance
component estimation of Helmert type. Zeitschrift für Vermessungswesen, 105, 161-
180
Grejner-Brzezinska, D. A., Kashani, I., Wielgosz, P., Smith, D. A., Spencer, P. S., Robertson, D.
S. & Mader, G. L. 2007. Efficiency and reliability of ambiguity resolution in network-
based real-time kinematic GPS. Journal of Surveying Engineering, 133, 56-65
Grodecki, J. 1997. Estimation of variance-covariance components for geodetic observations
and implications on deformation trend analysis. Geodesy and Geomatics Engineering
Technical Reports; 186,
Haase, J., Murphy, B., Muradyan, P., Nievinski, F., Larson, K., Garrison, J. & Wang, K. N.
2014. First results from an airborne GPS radio occultation system for atmospheric
profiling. Geophysical Research Letters,
Habrich, H. 2000. Geodetic applications of the global navigation satellite system (GLONASS)
and of GLONASS/GPS combinations, Verlag des Bundesamtes für Kartographie und
Geodäsie.
Han, S. 1997. Quality-control issues relating to instantaneous ambiguity resolution for real-
time GPS kinematic positioning. Journal of Geodesy, 71, 351-361
Hannah, J. 2001. Airborne gravimetry: a status report. Dunedin: Department of Surveying,
University of Otago.
He, K. & Xu, T. 2009. Kinematic Orbit Determination of GEO Satellite Based on Kalman
filtering. Journal of Geodesy and Geodynamics, 29, 109-112 (In Chinese)
He, K. & Xu, T. 2011. Orbit determination of GEO satellite based on forward and backward
Kalman filtering with elimination parameter. Science of Surveying and Mapping, 36,
53-55 (In Chinese)
Hein, G. W. 1995. Progress in airborne gravimetry: Solved, open and critical problems.
Proceedings of the IAG Symposium on Airborne Gravity Field Determination, IUGG
XXI General assembly. Boulder, Colorado, July.
Helmert, F. R. 1907. Die Ausgleichungsrechnung nach der Methode der kleinsten Quadrate,
Leipzig, Berlin, BG Teubner.
Herring, T., King, R. & McClusky, S. 2010. GAMIT Reference Manual. GPS Analysis at MIT.
Release 10.4. Massachussetts Institute Technology.
Herring, T. A., Davis, J. L. & Shapiro, I. I. 1990. Geodesy by radio interferometry: The
application of Kalman filtering to the analysis of very long baseline interferometry
data. Journal of Geophysical Research: Solid Earth (1978–2012), 95, 12561-12581
Hewitson, S. & Wang, J. 2006. GNSS receiver autonomous integrity monitoring (RAIM)
performance analysis. GPS Solutions, 10, 155-170
Heyde, I., Barthelmes, F. & Scheinert, M. 2013. First Results of the aerogravity measurements
during the geoscientific flight mission GEOHALO over Italy and the adjacent
Mediterranean EGU General Assembly 2013. Vienna, Austria.
Hofmann-Wellenhof, B., Lichtenegger, H. & Wasle, E. 2008. GNSS–global navigation
satellite systems: GPS, GLONASS, Galileo, and more, Wien, NewYork, Springer.
Hopfield, H. 1969. Two‐quartic tropospheric refractivity profile for correcting satellite data.
Journal of Geophysical research, 74, 4487-4499
Huep, W. 1985. Zur Positionsschätzung im gestörten KALMAN-Filter am Beispiel eines
manövrierenden Wasserfahrzeuges. Ph. D, Fachrichtung Vermessungswesen der Univ.
Hannover.
IAC. 2014. GLONASS Constellation Status [Online]. Available: http://glonass-
iac.ru/en/GLONASS/ [Accessed 02.04 2014].
- 151 -
Ifadis, I. 1986. The atmospheric delay of radio waves: Modeling the elevation dependence on
a global scale. Technical Report. Gothenburg, Sweden: Chalmers University of
Technology.
IGS. 2013. IGS Products [Online]. Available: http://igscb.jpl.nasa.gov/components/prods.html
[Accessed 30 April 2014].
IGS. 2014. The IGS Tracking Network [Online]. Available:
http://www.igs.org/network/complete.html [Accessed 12.02 2014].
InsideGNSS. 2010. Russia to Put 8 CDMA Signals on 4 GLONASS Frequencies [Online].
Available: http://www.insidegnss.com/node/1997 [Accessed 02.04 2014].
Janes, H., Langley, R. & Newby, S. 1991. Analysis of tropospheric delay prediction models:
comparisons with ray-tracing and implications for GPS relative positioning. Bulletin
géodésique, 65, 151-161
Jee, G., Lee, H. B., Kim, Y., Chung, J. K. & Cho, J. 2010. Assessment of GPS global
ionosphere maps (GIM) by comparison between CODE GIM and TOPEX/Jason TEC
data: Ionospheric perspective. Journal of Geophysical Research: Space Physics
(1978–2012), 115
Jonkman, E. 1980. Doppler research in the nineteenth century. Ultrasound in medicine &
biology, 6, 1-5
Julier, S. J. & LaViola, J. J. 2007. On Kalman filtering with nonlinear equality constraints.
Signal Processing, IEEE Transactions on, 55, 2774-2784
Kalman, R. E. 1960. A new approach to linear filtering and prediction problems. Journal of
Fluids Engineering, 82, 35-45
Kalman, R. E. & Bucy, R. S. 1961. New results in linear filtering and prediction theory.
Journal of Fluids Engineering, 83, 95-108
Kim, D. & Langley, R. B. 2001. Instantaneous real time cycle-slip correction of dual-
frequency GPS data. Proceedings of the international symposium on kinematic
systems in geodesy, geomatics and navigation. Banff, Alberta, Canada.
Kim, H.-S. & Lee, H.-K. 2012. Elimination of Clock Jump Effects in Low-Quality Differential
GPS Measurements. Journal of Electrical Engineering & Technology, 7, 626-635
Klobuchar, J. A. 1987. Ionospheric time-delay algorithm for single-frequency GPS users.
Aerospace and Electronic Systems, IEEE Transactions on, 325-331
Knight, N. L., Wang, J. & Rizos, C. 2010. Generalised measures of reliability for multiple
outliers. Journal of Geodesy, 84, 625-635
Ko, S. & Bitmead, R. R. 2007. State estimation for linear systems with state equality
constraints. Automatica, 43, 1363-1368
Koch, K. R. 1986. Maximum likelihood estimate of variance components. Bulletin
Gæodésique, 60, 329-338
Koch, K. R. 1999. Parameter estimation and hypothesis testing in linear models, Springer.
Koch, K. R. 2013. Robust estimation by expectation maximization algorithm. Journal of
Geodesy, 87, 107-116
Koch, K. R. & Kusche, J. 2002. Regularization of geopotential determination from satellite
data by variance components. Journal of Geodesy, 76, 259-268
Koch, K. R. & Yang, Y. 1998. Robust Kalman filter for rank deficient observation models.
Journal of Geodesy, 72, 436-441
Kouba, J. 2009a. A guide to using International GNSS Service (IGS) products. International
GNSS. Ottawa: Geodetic Survey Division of Natural Resources Canada.
Kouba, J. 2009b. Testing of global pressure/temperature (GPT) model and global mapping
function (GMF) in GPS analyses. Journal of Geodesy, 83, 199-208
Kouba, J. & Héroux, P. 2001. Precise point positioning using IGS orbit and clock products.
GPS solutions, 5, 12-28
Kreye, C. & Hein, G. W. 2003. GNSS Based Kinematic Acceleration Determination for
Airborne Vector Gravimetry-Methods and Results. Proceedings of the ION GPS/GNSS
2003 Meeting. Portland, Oregon.
Kubik, K. 1982. An error theory for the Danish method. Symposium of Comm. III of ISP.
Helsinki.
- 152 -
Kusche, J. 2003. Noise variance estimation and optimal weight determination for GOCE
gravity recovery. Advances in Geosciences, 1, 81-85
Kwon, J. H. & Jekeli, C. 2001. A new approach for airborne vector gravimetry using GPS/INS.
Journal of Geodesy, 74, 690-700
Lachapelle, G. & Alves, P. 2009. Multiple reference station approach: overview and current
research. Journal of Global Positioning Syetems, 1, 133-136
Lacoste, L. J. 1967. Measurement of gravity at sea and in the air. Reviews of Geophysics, 5,
477-526
Langley, R. B. 1999. Dilution of precision. GPS world, 10, 52-59
Leick, A. 2004. GPS satellite surveying, John Wiley & Sons.
Leick, A., Beser, J. & Rosenboom, P. 1998. Aspects of GLONASS carrier-phase differencing.
GPs Solutions, 2, 36-41
Li, B., Feng, Y., Shen, Y. & Wang, C. 2010. Geometry‐specified troposphere decorrelation
for subcentimeter real‐time kinematic solutions over long baselines. Journal of
Geophysical Research: Solid Earth (1978–2012), 115
Li, B., Shen, Y., Feng, Y., Gao, W. & Yang, L. 2013. GNSS ambiguity resolution with
controllable failure rate for long baseline network RTK. Journal of Geodesy, 88, 1-14
Li, D. 1985. Theorie und Untersuchung der Trennbarkeit von groben Paßpunktfehlern und
systematischen Bildfehlern bei der photogrammetrischen Punktbestimmung. PhD
thesis, Universität Stuttgart.
Li, D. 1989. A Thought of Optimization and Design of Geodetic Networks in Consideration of
Accuracy and Reliability. Festschrift Friedrich Ackermann zum 60. Geburtstag,
Li, X., Ge, M., Douša, J. & Wickert, J. 2014. Real-time precise point positioning regional
augmentation for large GPS reference networks. GPS solutions, 18, 61-71
Mader, G. L. 1999. GPS antenna calibration at the National Geodetic Survey. GPS solutions, 3,
50-58
Marreiros, J. P. R. 2012. Kinematic GNSS precise point positioning. Ph.D Thesis, University
of Porto.
McCarthy, D. D. 1996. IERS Conventions (1996). IERS Technical Note.
Melbourne, W. G. 1985. The case for ranging in GPS-based geodetic systems. In: Goad, C.
(ed.) Proceedings of first International Symposium on Precise Positioning with the
Global Positioning Systme. Rockville, Maryland: U.S. Department of Commerce.
Mendes, V. B. & Langley, R. B. 1994. A comprehensive analysis of mapping functions used in
modeling tropospheric propagation delay in space geodetic data. Proceeding of the
International Symposium on Kinematic Systems in Geodesy. Banff, Alberta.
Mendes, V. B. & Langley, R. B. 1995. Zenith Wet Tropospheric Delay Determination Using
Prediction Model: Accuracy Analysis. Cartograpia e Cadastro, 2, 41-47
Mendes, V. B. & Langley, R. B. 1998. Tropospheric zenith delay prediction accuracy for
airborne GPS high-precision positioning. Proceedings of The Institute of Navigation
54th Annual Meeting.
Mercado, D. A., Flores, G., Castillo, P., Escareno, J. & Lozano, R. 2013. GPS/INS/optic flow
data fusion for position and Velocity estimation. Unmanned Aircraft Systems (ICUAS),
2013 International Conference on. IEEE.
Mervart, L., Lukes, Z., Rocken, C. & Iwabuchi, T. 2008. Precise Point Positioning with
ambiguity resolution in real-time. Proceedings of ION GNSS.
Mikhail, E. M. & Ackermann, F. E. 1976. Observations and least squares.
Mohamed, A. H. & Schwarz, K. P. 1999. Adaptive Kalman filtering for INS/GPS. Journal of
Geodesy, 73, 193-203
Montenbruck, O., Hauschild, A., Steigenberger, P., Hugentobler, U., Teunissen, P. &
Nakamura, S. 2013. Initial assessment of the COMPASS/BeiDou-2 regional
navigation satellite system. GPS solutions, 17, 211-222
Moreno Monge, B. 2012. Development of algorithms for the GNSS data processing: their
application to the modernized GPS and Galileo scenarios. Ph. D, Universidad
Complutense de Madrid.
- 153 -
NAVCEN. 2014. General Informaiton on GPS [Online]. Available:
http://www.navcen.uscg.gov/?pageName=GPSmain [Accessed 01.04 2014].
Niell, A. E. 1996. Global mapping functions for the atmosphere delay at radio wavelengths.
Journal of Geophysical Research, 101, 3227-3246
Niell, A. E. 2001. Preliminary evaluation of atmospheric mapping functions based on
numerical weather models. Physics and Chemistry of the Earth, 26, 475-480
NOAA. 2014a. The GPS Control Segment [Online]. Available:
http://www.gps.gov/systems/gps/control/ [Accessed 01.04 2014].
NOAA. 2014b. The GPS Space Segment [Online]. Available:
http://www.gps.gov/systems/gps/space/#generations [Accessed 01.04 2014].
OpeNaviGate. 2012. GLONASS / GPS Integration [Online]. Available:
http://openavigate.blogspot.de/2012/01/glonass-gps-integration.html [Accessed 08.05
2014].
Parkinson, B. 1996. Introduction and Heritage of NAVSTAR, the Global Positioning System.
Global Positioning System: Theory and Applications, Volume I. American Institute of
Aeronautics and Astronautics. Inc., Cambridge, Massachusetts,
Parrot, D. 1989. Short-arc orbit improvement for GPS satellites. Geodesy and Geomatics
Engineering Technical Reports; 143,
Persson, C. G. 1980. MINQUE and related estimators for variance components in linear
models, Stockholm, Royal Institute of Technology.
Petit, G. & Luzum, B. 2010. IERS conventions (2010). DTIC Document.
Petrovic, S., Barthelmes, F. & Pflug, H. 2013. Airborne and shipborne gravimetry at GFZ with
emphasis on the GEOHALO project. IAG Scientific Assembly 2013. Potsdam,
Germany.
Radicella, S. & Leitinger, R. 2001. The evolution of the DGR approach to model electron
density profiles. Advances in Space Research, 27, 35-40
Revnivykh, S. 2012. GLONASS Status and Modernization. International GNSS Committee
IGC-7. Beijing.
Roßbach, U. 2001. Positioning and navigation using the Russian satellite system GLONASS.
Ph.D Thesis, Univ. d. Bundeswehr München.
Rothacher, M. 2001. Comparison of absolute and relative antenna phase center variations.
GPS solutions, 4, 55-60
Saastamoinen, J. 1972. Contributions to the theory of atmospheric refraction. Bulletin
Geodesique, 105, 279-298
Sahin, M., Cross, P. A. & Sellers, P. C. 1992. Variance component estimation applied to
satellite laser ranging. Bulletin Geodesique, 66, 284-295
Schüler, T. 2001. On ground-based GPS tropospheric delay estimation. Ph. D, Univ. der
Bundeswehr München.
Schaffrin, B. 1997. Reliability measures for correlated observations. Journal of Surveying
Engineering, 123, 126-137
Schaffrin, B. Some generalized equivalence theorems for least-squares adjustment. In: Sanso,
F., ed. V Hotine-Marussi Symposium on Mathematical Geodesy, 2004. Springer, 107-
112.
Schaffrin, B. & Grafarend, E. 1986. Generating classes of equivalent linear models by
nuisance parameter elimination. Manuscripta geodaetica, 11, 262-271
Scheinert, M., Petrovic, S., Heyde, I., Barthelmes, F. & GEOHALO, G. 2013. The geodetic -
geophysical flight mission GEOHALO: Results of airborne gravimetry and further
geodetic products. IAG Scientific Assembly 2013. Potsdam, Germany.
Schmid, R., Steigenberger, P., Gendt, G., Ge, M. & Rothacher, M. 2007. Generation of a
consistent absolute phase-center correction model for GPS receiver and satellite
antennas. Journal of Geodesy, 81, 781-798
Schwarz, K. P., Cannon, M. E. & Wong, R. V. C. 1989. A comparison of GPS kinematic
models for the determination of position and velocity along a trajectory. Manuscripta
geodaetica, 14, 345-353
- 154 -
Schwarz, K. P. & Wei, M. 1995. Some unsolved problems in airborne gravimetry. In: Sünkel,
H. (ed.) Gravity and Geoid. Berlin, Heidelberg: Springer.
Seeber, G. 2003. Satellite geodesy: foundations, methods, and applications, Walter de Gruyter.
Semmling, A. M., Beckheinrich, J., Wickert, J., Beyerle, G., Schön, S., Fabra, F., Pflug, H., He,
K., Schwabe, J. & Scheinert, M. 2014. Sea surface topography retrieved from GNSS
reflectometry phase data of the GEOHALO flight mission. Geophysical Research
Letters,
Semmling, A. M., Beyerle, G., Schön, S., Beckheinrich, J., Wickert, J. & Scheinert, M. 2013.
Airborne GNSS Reflectometry for Sea Surface Height Estimation as part of the
GEOHALO Mission. IAG Scientific Assembly 2013. Potsdam, Germany.
Serrano, L., Kim, D. & Langley, R. B. 2004a. A single GPS receiver as a real-time, accurate
velocity and acceleration sensor. Proceedings of the 17th International Technical
Meeting of the Satellite Division of the Institute of navigation.
Serrano, L., Kim, D., Langley, R. B., Itani, K. & Ueno, M. 2004b. A gps velocity sensor: how
accurate can it be?–a first look. ION NTM. San Diego, CA.
Shen, Y. & Xu, G. 2008. Simplified equivalent representation of GPS observation equations.
GPS Solutions, 12, 99-108
Shi, C., Zhao, Q., Hu, Z. & Liu, J. 2013. Precise relative positioning using real tracking data
from COMPASS GEO and IGSO satellites. GPS solutions, 17, 103-119
Shi, J. & Cannon, M. 1995. Critical error effects and analysis in carrier phase-based airborne
GPS positioning over large areas. Bulletin géodésique, 69, 261-273
Shi, W. 2010. Principles of modeling uncertainties in spatial data and spatial analyses, CRC
Press.
Simon, D. & Chia, T. L. 2002. Kalman filtering with state equality constraints. Aerospace and
Electronic Systems, IEEE Transactions on, 38, 128-136
Simsky, A. & Boon, F. 2003. Carrier Phase and Doppler-Based Algorithms for Real-Time
Standalone Positioning. Proceedings of the GNSS 2003. Graz, Austria.
Singh, D., Ghosh, J. K. & Kashyap, D. 2014. Precipitable water vapor estimation in India from
GPS-derived zenith delays using radiosonde data. Meteorology and Atmospheric
Physics, 123, 209-220
Snay, R. A. & Soler, T. 2008. Continuously operating reference station (CORS): history,
applications, and future enhancements. Journal of Surveying Engineering, 134, 95-104
Stephenson, S., Meng, X., Moore, T., Baxendale, A. & Edwards, T. 2011. Precision of
Network Real Time Kinematic Positioning for Intelligent Transport Systems.
European Navigation Conference.
Stombaugh, T. S. & Clement, B. R. 1999. Unraveling the GPS Mystery [Online]. Ohio State
University. Available: http://ohioline.osu.edu/aex-fact/0560.html [Accessed 02. 11
2013].
Szarmes, M., Ryan, S., Lachapelle, G. & Fenton, P. 1997. DGPS high accuracy aircraft
velocity determination using Doppler measurements. Proceedings of the International
Symposium on Kinematic Systems (KIS). Banff, AB, Canada.
Tahk, M. & Speyer, J. L. 1990. Target tracking problems subject to kinematic constraints.
Automatic Control, IEEE Transactions on, 35, 324-326
Takasu, T. 2013. RTKLIB ver 2.4.2 [Online]. Available:
http://www.rtklib.com/prog/manual_2.4.2.pdf [Accessed 01.05 2014].
Teunissen, P. J. G. & Amiri-Simkooei, A. R. 2008. Least-squares variance component
estimation. Journal of Geodesy, 82, 65-82
Teunissen, P. J. G., De Jonge, P. J. & Tiberius, C. C. J. M. 1996. The Volume of the GPS
Ambiguity Search Space and its Relevance for Integer Ambiguty Resolution.
PROCEEDINGS OF ION GPS. Citeseer.
Teunissen, P. J. G., Joosten, P. & Odijk, D. 1999. The reliability of GPS ambiguity resolution.
GPS Solutions, 2, 63-69
Teunissen, P. J. G. & Kleusberg, A. 1998. GPS for Geodesy, Berlin, Heidelberg, Springer.
Teunissen, P. J. G., Odijk, D. & Zhang, B. 2010. PPP-RTK: Results of CORS network-based
PPP with integer ambiguity resolution. J Aeronaut Astronaut Aviat Ser A, 42, 223-230
- 155 -
Thompson, L. G. D. & LaCoste, L. J. B. 1960. Aerial gravity measurements. Journal of
Geophysical Research, 65, 305-322
Tiemeyer, B., Cannon, M., Lachapelle, G. & Schanzer, G. 1994. Satellite navigation for high
precision aircraft navigation with emphasis on atmospheric effects. Position Location
and Navigation Symposium. IEEE.
Torge, W. 1989. Gravimetry, Berlin; New York, Walter de Gruyter & Co.
Vennebusch, M., Böckmann, S. & Nothnagel, A. 2007. The contribution of very long baseline
interferometry to ITRF2005. Journal of Geodesy, 81, 553-564
Wübbena, G. 1985. Software developments for geodetic positioning with GPS using TI-4100
code and carrier measurements. In: Goad, C. (ed.) Proceedings of the First
International Symposium on Precise Positioning with the Global Positioning System.
Rockville, Maryland: U.S. Department of Commerce.
Wübbena, G., Bagge, A., Seeber, G., Boeder, V. & Hankemeier, P. 1996. Reducing distance
dependent errors for real-time precise DGPS applications by establishing reference
station networks. In Proceedings of the 9th International Technical Meeting of the
Satellite Division of the Instiute of Navigation. Kansas City, Missouri: Institute of
Navigation.
Wang, J., Gopaul, N. & Scherzinger, B. 2009a. Simplified algorithms of variance component
estimation for static and kinematic GPS single point positioning. Journal of Global
Positioning Systems, 8, 43-52
Wang, J., Rizos, C., Stewart, M. P. & Leick, A. 2001. GPS and GLONASS integration:
modeling and ambiguity resolution issues. GPS solutions, 5, 55-64
Wang, Q. & Xu, T. 2011. Combining GPS carrier phase and Doppler observations for precise
velocity determination. Science China Physics, Mechanics and Astronomy, 54, 1022-
1028
Wang, Q., Xu, T. & Xu, G. 2010. HALO_GPS (High Altitude and Long Range Airborne GPS
Positioning Software)-Software User Manual. Potsdam, Germany:
GeoForschungsZentrum.
Wang, Q., Xu, T. & Xu, G. 2011. Adaptively changing reference station algorithm and its
application in GPS long range airborne kinematic relative positioning. Acta
Geodaetica et Cartographica Sinica, 4, 429-434 (In Chinese)
Wang, X., Ji, J. & Li, Y. 2009b. The applicability analysis of troposphere delay error model in
GPS positioning. Aircraft Engineering and Aerospace Technology, 81, 445-451
Wei, Z. & Ge, M. 1998. Mathematical model for GPS relative positioning, Beijing, China,
Surveying and Mapping Press ( In Chinese).
Wielgosz, P., Kashani, I. & Grejner-Brzezinska, D. 2005. Analysis of long-range network RTK
during a severe ionospheric storm. Journal of Geodesy, 79, 524-531
Xu, G. 2002. GPS data processing with equivalent observation equations. GPS solutions, 6,
28-33
Xu, G. 2004. MFGsoft Software User Manual. Version of 2004. Potsdam, Germany:
Geoforschungszentrum.
Xu, G. 2007. GPS: Theory,Algorithms and Applications, Berlin, Heideberg, New York,
Springer.
Xu, G., Shen, Y., Yang, Y., Sun, H., Zhang, Q., Guo, J. & Yeh, T.-K. 2010. Equivalence of GPS
Algorithms and Its Inference. In: Xu, G. (ed.) Sciences of Geodesy-I Advances and
Future Direcitons. Heidelberg, New York: Springer.
Xu, P., Liu, Y., Shen, Y. & Fukuda, Y. 2007. Estimability analysis of variance and covariance
components. Journal of Geodesy, 81, 593-602
Xu, T., He, K. & Xu, G. 2012. Orbit determination and thrust force modeling for a maneuvered
GEO satellite using two-way adaptive Kalman filtering. Science China Physics,
Mechanics and Astronomy, 55, 738-743
Yang, Y. 1991. Robust bayesian estimation. Bulletin géodésique, 65, 145-150
Yang, Y. 1994. Robust estimation for dependent observations. Manuscripta geodaetica, 19,
10-17
- 156 -
Yang, Y. 2010. Adaptively robust kalman filters with applications in navigation. In: Xu, G. (ed.)
Sciences of Geodesy-I Advances and Future Directions. Heidelberg, Dordrecht,
London, New York: Springer.
Yang, Y., Cheng, M. K., Shum, C. K. & Tapley, B. D. 1999. Robust estimation of systematic
errors of satellite laser range. Journal of Geodesy, 73, 345-349
Yang, Y., Gao, W. & Zhang, X. 2010. Robust Kalman filtering with constraints: a case study
for integrated navigation. Journal of Geodesy, 84, 373-381
Yang, Y., He, H. & Xu, G. 2001. Adaptively robust filtering for kinematic geodetic positioning.
Journal of Geodesy, 75, 109-116
Yang, Y., Li, J., Xu, J., Tang, J., Guo, H. & He, H. 2011a. Contribution of the compass satellite
navigation system to global PNT users. Chinese Science Bulletin, 56, 2813-2819
Yang, Y., Xu, T. & Song, L. 2005. Robust estimation of variance components with application
in global positioning system network adjustment. Journal of surveying engineering,
131, 107-112
Yang, Y., Zhang, X. & Xu, J. 2011b. Adaptively constrained Kalman filtering for navigation
applications. Survey Review, 43, 370-381
Zhang, B., Teunissen, P. J. & Odijk, D. 2011a. A novel un-differenced PPP-RTK concept. J.
Navig, 64, S180-S191
Zhang, J. 2007. Precise velocity and acceleration determination using a standalone GPS
receiver in real time. Ph.D, Royal Melbourne Institute of Technology.
Zhang, J., Zhang, K., Grenfell, R. & Deakin, R. 2006. Short note: On the relativistic Doppler
effect for precise velocity determination using GPS. Journal of Geodesy, 80, 104-110
Zhang, J., Zhang, K., Grenfell, R., Li, Y. & Deakin, R. 2005. Real-time Doppler/Doppler rate
derivation for dynamic applications. Journal of Global Positioning Syetems, 4, 95-105
Zhang, X. & Andersen, O. B. 2006. Surface ice flow velocity and tide retrieval of the Amery
ice shelf using precise point positioning. Journal of Geodesy, 80, 171-176
Zhang, X., Guo, B., Guo, F. & Du, C. 2013. Influence of clock jump on the velocity and
acceleration estimation with a single GPS receiver based on carrier-phase-derived
Doppler. GPS solutions, 17, 549-559
Zhang, X., Li, X. & Guo, F. 2011b. Satellite clock estimation at 1 Hz for realtime kinematic
PPP applications. GPS solutions, 15, 315-324
Zhao, Q., Guo, J., Li, M., Qu, L., Hu, Z., Shi, C. & Liu, J. 2013. Initial results of precise orbit
and clock determination for COMPASS navigation satellite system. Journal of
Geodesy, 87, 475-486
Zhou, J. 1985. On the Jie factor. Acta Geodaetica et Geophysica, 5, (In Chinese)
Zhou, J. 1989. Classical theory of errors and robust estimation. Acta Geodetica Et
Cartographica Sinica, 18, 115-120, (In Chinese)
Zumberge, J., Heflin, M., Jefferson, D., Watkins, M. & Webb, F. 1997. Precise point
positioning for the efficient and robust analysis of GPS data from large networks.
Journal of Geophysical Research, 102, 5005-5017
- 157 -
Acknowledgments
The completion of this thesis would not have been possible without the support of many peo-
ple and organizations.
First of all I would like to express my profound respect and gratitude to my supervisors,
Prof. Dr.-Ing. Frank Neitzel of Berlin University of Technology (TU Berlin) and Prof. h.c. Dr.-
Ing. Guochang Xu (formerly German Research Centre for Geosciences (GFZ), now Shandong
University Weihai Branch (SDUW) in China), for supervising this thesis, for their continuous
support, guidance and encouragement of my studies and research. They were always willing to
share their insights with me and were patient with every question.
I am grateful to Prof. Dr. Frank Flechtner, head of GFZ Section 1.2 for his kindly help, or-
ganization of my work and improvements of my publication. Particularly I would like to thank
my GFZ colleague Dr. Christoph Förste for his critical comments and helpful suggestions to
improve the English in this thesis. Furthermore, the German abstract was translated by Dr.
Christoph Förste since I am not very familiar with the German language.
My special words of thanks should go to Prof. Dr.-Ing. Dr. h.c. Harald Schuh and Prof. Dr.-
Ing. Frank Flechtner of GFZ and TU Berlin, Prof. Dr.-Ing. Frank Neitzel of TU Berlin, and
Prof. Dr. Phil. nat. Markus Rothacher of ETH Zürich for reviewing the thesis and attending
dissertation defense.
I owe a special thank to my former GFZ colleague and advisor Dr. Tianhe Xu, whom I
thank so much for his continuous and extensive support throughout the last years. His selfless
help, detailed guidance and discussion are the source of my forward momentum. I also wish to
thank my former GFZ colleagues Dr. Guanwen Huang and Dr. Qianxin Wang for their discus-
sions, which helped to put ideas together.
I would also like to thank my GFZ colleagues Dr. Franz Barthelmes, Mr. Hartmut Pflug,
Ms. Anglika Svarovsky, M. Sc. Yan Xu and Dr. Maiko Ritschel. I appreciate their kindly help,
cooperation and discussions. I have enjoyed working with all of my colleagues.
Further, many thanks to my GFZ colleagues Dr. Maorong Ge, Dr. Zhiguo Deng, Dr. Jens
Wickert, Dr. Maximilian Semmling, Dr. Ming Shangguan, Dr. Xinging Li, and many PhD stu-
dents of TU Berlin and GFZ, Yumiao Tian, Nan Jiang, LiWah Wong, Rui Tu, Kejie Chen, Hua
Chen and others, for their kindly help and discussions.
I am indebted to my host institute GFZ where I was doing my research from 2010 to 2014.
TU Berlin is also thanked for my study during this period. The organizers and participants of
- 158 -
the GEOHALO campaign are also thanked for kindly providing GNSS data of HALO aircraft.
The partners of the International GNSS Service (IGS) are thanked for providing GNSS data
and products.
Gratefully acknowledged are the China Scholarship Council (CSC) and the Education Sec-
tion of the Embassy of the People’s Republic of China in the Federal Republic of Germany,
which financially supported the work on my PhD research for four years. Also many thanks to
GFZ for financially supporting me in the last four months.
Finally, I would like to express my gratitude to my parents for their endless support and
understanding. Also thanks to my sister-in-law and my elder brother for encouraging me to
study and helping me to take care of our parents. The person I wish to thank most is my girl-
friend, Dr. Lu Zhang, who fills my heart with confidence, courage, and love all the time.
