GPS-based attitude determination

(1)

Institutionen för systemteknik

Department of Electrical Engineering

Examensarbete

GPS-based attitude determination

Examensarbete utfört i Reglerteknik vid Tekniska högskolan i Linköping

av

Johan Bejeryd

LiTH-ISY-EX--07/4127--SE

Linköping 2007

Department of Electrical Engineering Linköpings tekniska högskola

Linköpings universitet Linköpings universitet

(2)

(3)

GPS-based attitude determination

Examensarbete utfört i Reglerteknik

vid Tekniska högskolan i Linköping

av

Johan Bejeryd

LiTH-ISY-EX--07/4127--SE

Handledare: Björn Gabrielsson / Fredrik Neregård

Saab Bofors Dynamics AB

David Törnqvist

isy, Linköpings universitet

Examinator: Fredrik Gustafsson

isy, Linköpings universitet

(4)

(5)

Avdelning, Institution

Division, Department

Division of Automatic Control Department of Electrical Engineering Linköpings universitet

SE-581 83 Linköping, Sweden

Datum Date 2007-12-09 Språk Language Svenska/Swedish Engelska/English Rapporttyp Report category Licentiatavhandling Examensarbete C-uppsats D-uppsats Övrig rapport

URL för elektronisk version

http://www.control.isy.liu.se http://urn.kb.se/resolve?urn=urn:nbn:se:liu:diva-4127 ISBN — ISRN LiTH-ISY-EX--07/4127--SE

Serietitel och serienummer

Title of series, numbering

ISSN

—

Titel

Title

Attitydmätning med GPS GPS-based attitude determination

Författare

Author

Johan Bejeryd

Sammanfattning

Abstract

Inertial sensors and magnetometers are often used for attitude determination of moving platforms. This thesis treats an alternative method; GPS-based at-titude determination. By using several GPS-antennas, and with carrier phase measurements determining the relative distance between them, the attitude can be calculated.

Algorithms have been implemented in Matlab and tested on real data. Two commercial GPS-based attitude determination systems have also been tested on a mobile platform and compared to a navigation grade Inertial Navigation System (INS). The results from the tests show that GPS-based attitude determination works well in open areas, but would require support from additional sensors in urban and forest environments.

Nyckelord

(6)

(7)

Abstract

Inertial sensors and magnetometers are often used for attitude determination of moving platforms. This thesis treats an alternative method; GPS-based at-titude determination. By using several GPS-antennas, and with carrier phase measurements determining the relative distance between them, the attitude can be calculated.

Algorithms have been implemented in Matlab and tested on real data. Two commercial GPS-based attitude determination systems have also been tested on a mobile platform and compared to a navigation grade Inertial Navigation System (INS). The results from the tests show that GPS-based attitude determination works well in open areas, but would require support from additional sensors in urban and forest environments.

(8)

(9)

Acknowledgments

First of all I would like to thank the involved people at Saab Bofors Dynamics for giving me the oppurtunity of working with this thesis. The work has been both challenging and rewarding. I would also like to thank Mårten Lindgren at FMV for providing a GPS receiver that has been well used. Additionally I would like to thank my examiner Fredrik Gustafsson and my supervisor David Törnqvist at Linköpings Universitet for helpful comments on the thesis.

Finally I would especially like to thank my supervisors at Saab Bofors Dynam-ics; Fredrik Neregård for starting up the work, and Björn Gabrielsson for support and guidance troughout the work.

(10)

(11)

Introduction

1.1 Background

For attitude determination of moving platforms it is common to use inertial sensors and magnetometers. Basically two distinctions can be made; high grade gyros and accelerometers are combined to form an INS which can operate autonomously and find the initial attitude, lower grade gyros are combined with accelerometers and use magnetometers for initial attitude estimates as well as for long term support. Inertial navigation systems are often very expensive but on the other hand deliver high precision measurements, while magnetometers are sensitive to anomalies and disturbances from for example technical equipment. An interesting alternative is GPS-based attitude determination which could have a lower cost, be free from drifts in the attitude estimates and deliver good precision.

1.2 Purpose

The purpose of this work has been to examine the principles of GPS-based attitude determination. Saab Bofors Dynamics (SBD) expressed a need to understand if GPS-based attitude determination could be an alternative in different applications.

Two main goals have been defined:

• To understand the principles for GPS-based attitude determination,

includ-ing basic theory and description of the algorithms that are used.

• To gain practical experience of what performance that can be expected.

1.3 Solution method

The work started with a literature study to gain knowledge of GPS-based attitude determination. As practical experience of the expected performance was a goal, a data acquisition system was assembled. Data was then recorded using a com-mercial GPS attitude system which also logged raw data that could be used for

(16)

2 Introduction

algorithm development and evaluation. To evaluate the performance of the GPS attitude solution a navigation grade INS was used as reference. Algorithms were later implemented and the results evaluated against the commercial GPS attitude system, as well as the INS. A simple simulation model was also implemented to facilitate the algorithm development.

1.4 Simplifications

The focus of this master thesis has been on the principles of GPS-based attitude de-termination. Therefore, the implementation of the algorithms has been simplified and not intended for real use. The simplifications concern real time performance and the efficiency of handling events where the composition of tracked satellites changes. Another simplification is that only two GPS-antennas have been used, which implies that only two of the attitude angles have been computed. This thesis has not considered using signals from the Russian GLONASS satellite navigation system.

1.5 Related work

GPS-based attitude determination has been a subject for research for many years and still is. The topic is closely related to Carrier Phase Differential GPS (CPDGPS), where carrier phase measurements from several observers are used for high preci-sion positioning. The GPS-based attitude determination problem is a special case of CPDGPS with the extra constraint of the known baseline length.

The major challenge is to find the unknown initial number of whole carrier phase cycles between observers. Much of the research has treated this problem and also the problem of validating the result. Examples of algorithms are the LAMBDA method developed at the University of Delft and the LSAST method with derivates developed at the University of Calgary [20] [21].

Several commercial GPS-based attitude determination systems are available, but the number is limited. Presently most of the systems are based upon high quality receivers with higher costs. Included here are the JNSGyro receivers from Javad Navigation Systems and the different PolaRx2e receivers from Septen-trio Satellite Navigation [18] [11]. There are also receivers available from Hemi-sphere GPS that come at a lower price. Older receivers include the TANS Vector, Ashtech/Thales/Magellan ADU5 and Surrey Satellite Technology Limited SGR-20 [2].

GPS-based attitude determination is used in a variety of applications, for ex-ample in automatic control of agricultural machines. Other exex-amples are marine navigation and attitude determination for small unmanned aerial vehicles [7]. Even space craft are viable platforms [2].

(17)

1.6 Outline 3

1.6 Outline

The subsequent part of the report consists of the following chapters.

Chapter 2 is an introduction to the Global Positioning System, including a brief

description of signals and measurements.

Chapter 3 treats the attitude determination problem using GPS measurements.

Chapter 4 describes methods for integer ambiguity resolution.

Chapter 5 treats the algorithms that have been implemented during the

the-sis work.

Chapter 6 describes the equipment used for the practical tests.

Chapter 7 describes the achieved results from practical tests of commercial

GPS-based attitude determination systems and the implemented algorithms.

(18)

(19)

Chapter 2

Introduction to GPS

The Global Positioning System (GPS) is a satellite based navigation system ini-tially developed by the U.S. Department of Defense (DoD) and operated by the U.S. Air Force. The system provides accurate three-dimensional positioning and velocity for an unlimited number of users. GPS also offers very precise timing for its users since timing is at the heart of the technology that makes it work. The system is continuously online with satellite coverage of the whole globe [15]

Although GPS is a military system, there are now many civilian applications, such as automobile navigation and geodesy. GPS receivers have recently also started to appear in cell phones. The development of the system began in the seventies and it was declared operational in 1995, but the development has not halted; instead major upgrades are planned [12].

2.1 The satellite constellation

The space segment of GPS is made up of the satellites orbiting the earth at a Medium Earth Orbit (MEO) with a height above the earth surface of approxi-mately 20200 km. The orbit time is 11h 58m, thus the satellite constellation for a given time is repeated four minutes earlier the next day. The satellites travel at a speed of approximately 4 km/s; however the satellites positions can be estimated with errors less than a few meters when up to two day’s old predictions are used [15].

There are at least 24 satellites orbiting the earth and they are divided in six orbital planes, where each plane has at least four satellites. At the present time of writing 31 satellites are orbiting the earth [1]. The satellites are not evenly spaced in the planes; instead the constellation is designed to be robust if there would be a satellite failure. The constellation is also designed so there are continuously four visible satellites over the horizon at the places where the system is intended to work.

(20)

6 Introduction to GPS

2.2 The control segment

Spread around the world is a network of monitor stations. The main control station is located in Colorado at the Schriever Air Force Base. Data from the monitor stations is sent there and processed to update the satellites’ navigation status. Through uplinks the satellites are time synchronized and their navigation messages are uploaded [15].

2.3 The user segment

A usual application is made up of a receiving antenna, a signal processing unit, and an external interface for presenting information. There are a variety of GPS receivers available on the market. For civilian purposes there is everything from small OEM receiver modules, at rates from less than $10 and above, to high quality geodetic grade receivers at prices ranging from thousands of dollars and above.

2.4 GPS signals and measurements

Presently the GPS satellites emit two carrier signals, the L1 with the frequency 1575.42 MHz and the L2 with the frequency 1228.6 MHz. In general the L1 signal is used by both military and civil applications, and the L2 signal only for military signals. However some civil dual frequency receivers can partly use the L2 frequency.

Modulated on the L1 carrier is the Coarse / Acquisition (C/A) code which is freely available, but also the Precise code which is usually encrypted and dedicated for military applications only. The C/A code is a chip pseudo-random code (PRN) with a chip rate of 1.023 million chips/sec. Each satellite has its own unique PRN code that can be distinguished among the others, which also is a requirement as all satellites transmit on the same frequencies. The receivers can generate the PRN codes and by measuring the difference between the on board generated signal and the received signal the transmission time can be calculated. A measure of distance is attained by multiplying the signal propagation time with the speed of light. This distance is commonly named as the pseudorange as it is an approximation.

Also modulated on the carrier is the satellite navigation message, which con-tains information of satellite status and also ephemeris data which is used to calculate the satellite’s positions. The ephemerides contains lists of parameters which are used for predicting the satellite positions.

The width of the C/A code is 1 ms which corresponds to roughly 300 meters of length. Receivers can usually measure the C/A code within a couple of meters. The carrier phase which also can be measured has a wavelength of approximately 19 cm for L1, and provides a more accurate distance measurement. In fact the receivers measure the phase within a centimeter or even down to tenths of millime-ters. However the phase measurement is ambiguous as only the fractional part of phase is measured. But the receivers accumulate the phase since start of tracking and only the initial integer number of cycles is unknown.

(21)

2.4 GPS signals and measurements 7

Figure 2.1. A simplified block diagram of a carrier tracking loop.

The receiver uses carrier tracking loops to keep track of the phase for the sig-nals. There are two main types of carrier loop discriminators, the phase lock loop (PLL), and the frequency lock loop (FLL). Costas PLL is a PLL that can operate with data modulated on the processed signal, and it is named after its original inventor. The PLLs are the most accurate discriminators with long integration time, while the FLL is more tolerant to dynamic stress; however the FLL only tries to track the frequency and not the phase. A good receiver design should have a combination of both types to give good accuracy and tolerance for dynamic stress. Initially the FLL should be used and when the replica carrier is aligned enough the loop should transform into a PLL, that can track the phase. [15]

The GPS receiver processes the signals in several steps, a thorough description can be found in [12]. A simplified block diagram is shown in Figure 2.1. When the signal has been converted to an intermediary frequency and the C/A code has been stripped off then the carrier is removed by using an onboard generated replica of the carrier including the Doppler. It is this replica driven by the onboard oscillator that is adjusted in the loop to attain the phase measurement. By multiplication with the sine the in-phase data (I) is attained and by using the cosine the quadraphase data (Q) is attained. The I and Q data still contain mostly noise and are further processed before used to measure the phase. Suppose that the phase loop is in lock, then the in-phase data should contain the data, and the quadraphase data should ideally be zero. In reality both will also contain noise. If the loop is not in lock then the Q data will contain a part of the signal. As the I and Q data are 90◦ apart, the phase angle from the I-axis can be calculated as tan−1Q_I. This value is used as feedback to adjust the replica carrier and minimize the error, if a PLL is used. When the PLL is in lock the off-phase Q data will only contain noise and the in-phase I data will contain the signal with noise.

When the phase lock loop is designed the bandwidth must be taken into con-sideration. High bandwidth gives short rise times and enables the receiver to keep track in high dynamics; however the noise level is then increased. Usually a second order system is used, but a second order system has a steady state error when the acceleration is constant. This is a dilemma for the receiver designer. For very high dynamic environments one solution could be to tightly couple the receiver to inertial sensors. The PLL could also be extended to a third order system instead of the common second order system. A third order system removes the steady state

(22)

8 Introduction to GPS

Figure 2.2. Pseudo-trilateration, a 2D-example

error for constant acceleration, but might introduce other difficulties. In spite of the challenges of obtaining a good phase measurement, the high grade commercial receivers today succeed well. [15]

2.5 Principle of GPS positioning

In general four visible satellites are required to calculate the receiver position. The unknowns are the three coordinates in space and the receiver clock error. The receiver clock error needs to be resolved as the receiver must be synchronized to be able to generate the PRN codes and measure the pseudoranges with good accuracy.

The positioning is a matter of trilateration, or actually pseudo-trilateration as the pseudorange is not the real distance. It is affected by several errors, some that can be modeled and corrected for and some that are left unknown. If only one satellite is used the receiver position would be known to lie on the surface of a sphere centered around the satellite, and with a radius of the pseudorange. As the pseudorange is not a precise measurement an uncertainty is added and the receiver position would instead be known to lie in a thin shell. If three satellites were to be used then the position would be known to lie in a small volume where all three shells intersect, (assuming the clock offset is known). A two dimensional example is plotted in Figure 2.2 where the black area is the intersection where the receiver is known to be within.

The equation for the pseudorange can be written as

p = ks − uk2+ cto+ e0+ e1 (2.1)

(23)

2.5 Principle of GPS positioning 9

p := pseudorange measurement

u := receiver position

s := satellite position

to := receiver to satellite time offset

c := speed of light in vacuum

e0 := unmodeled errors

e1 := modeled errors.

Expanding the equation yields

pi=

p

(xi− xu)2+ (yi− yu)2+ (zi− zu)2+ cto+ e0+ e1. (2.2)

Where xi, yi, zi denotes the position of satellite i. The four unknowns in the

equation are the receiver coordinates xu, yu, zu and the clock offset to. The

mod-elled errors can be approximated and partly corrected for; however this is not described here. By combining measurements from at least four satellites the equa-tion system can be solved and the posiequa-tion calculated. However note that equaequa-tion (2.2) is non-linear and must be linearized, or treated in other ways, before solved.

(24)

(25)

Chapter 3

Attitude determination

3.1 Coordinate systems

GPS-systems usually uses an Earth Centered Earth Fixed (ECEF) coordinate system named WGS-84. When determining attitude a local coordinate system is needed. This normally coincides with the north, east and up or down directions. Examples of such systems are the right handed East North Up (ENU) and North East Down (NED) coordinate systems. In this work the NED system was chosen, with the coordinates defined as

x := north,

y := east,

z := down.

To transform coordinates from the earth fixed ECEF system to the local NED system a rotation matrix is used.

CECEFN ED =





− sin(φ) cos(λ) − sin(φ) sin(λ) cos(φ)

− sin(λ) cos(λ) 0

− cos(φ) cos(λ) − cos(φ) sin(λ) − sin(φ)



, (3.1) where φ denotes the latitude and λ denotes the longitude, both expressed in WGS-84.

3.2 Determining attitude

To determine the attitude the idea is to calculate the relative positions between two or more antennas. If two relative positions are known, the vector between them contains information of heading and elevation angle. But with only two relative positions the roll angle can not be determined, as this is the rotation around the own axis. To solve for the roll angle another relative position must be

(26)

12 Attitude determination

Figure 3.1. The phase difference between two antennas

known. This position must also form a non-collinear baseline to give the needed extra information.

The relative position between the antennas is estimated using carrier phase measurements, as the C/A code measurements are not precise enough. The phase difference between the antennas is used to determine the distance between the antennas in the satellite direction, see Figure 3.1. By combining measurements from several satellites the relative positions can be calculated. The satellites are so far away from the receiving antennas that the incoming wave fronts can be approximated to be planar. However as mentioned earlier the initial number of integer cycles is ambiguous. This will be described more in the following sections.

3.3 Measurement model

The fractional part of the phase is measured by the carrier tracking loops, but the receiver also counts the number of whole cycles that passes by; however the initial whole number of cycles between the receiver and satellite antennas is unknown. To attain a precise distance measurement these integer ambiguities must be resolved. [12]

The phase measurements are affected by several error sources. When propa-gating from the satellites to the observers GPS signals are affected by interference from the atmosphere. The higher refractive index of the atmosphere makes the signals refract when they enter from space, and the passage through the tropo-sphere and ionotropo-sphere introduces other errors such as for example phase delay. Timing is another error source, the satellite and receiver clocks are not perfectly synchronized.

The following equation summarizes the observation of the carrier phase, [2]:

φi_A= ri_A+ βA− λNAi +  i₊i

(27)

3.3 Measurement model 13

where

φi

A := integrated Doppler (meters),

ri_A := distance from observer A to satellite i (meters),

βA := clock bias for observer A (meters),

λ := carrier wavelength (meters),

N_Ai := integer ambiguity (cycles),

i := measurement error specific for satellite i (meters),

i

A := measurement error specific for satellite i and observer A (meters).

For attitude determination the relative distance between several observers is interesting. As mentioned before the attitude can be calculated if several relative positions are known. The vector between two such positions is called a baseline and used to calculate the attitude. Motivated by this the single difference of phase measurements is formed by subtracting the phase measurement of observer A for satellite i with the phase measurement of observer B for satellite i.

∆φi_AB = φi_A− φi B = r_Ai + βA− λNAi +  i₊i A− (r i B+ βB− λNBi +  i₊i B) = rAi − rBi − λ∆NABi + ∆βAB+ ∆iAB (3.3) where λ∆NABi := λNAi − λNBi ∆βAB := βA− βB ∆iAB := i A−  i B. (3.4)

The atmospheric disturbances on GPS signals can be assumed to be equal for different observers if they are close. Thus the measurement error specific for the observers and satellite i is assumed to be eliminated in the single differences. However the clock bias is not eliminated and this is a good reason for forming the double differences, which are the difference of single differences for two different satellites. More generally also a linear combination of single differences could be used. The ∇∆ denotes the double differences of a measurement, ∇ is not the gradient. ∇∆φij_AB = ∆φi_AB− ∆φj_AB = riA− r i B− r j A+ r j B− λ∇∆N ij AB+ ∇∆ ij AB (3.5) where λ∇∆N_ABij := λ∆NABi − λ∆N j AB (3.6) ∇∆ij_AB := ∆i_AB− ∆j_AB. (3.7)

(28)

The clock bias between the receivers or observers is now eliminated but at the cost of increased measurement noise. But note that it is a linearized model. The clock bias is not allowed to be too large or the equation is not longer correct which means that the clock bias can not be expected to be eliminated.

3.3.1 Error sources

Noise

Thermal noise depends on the power of the incoming signal, the antenna signal amplifier and the receiver amplifier. However it is usually a small error source of millimeter or submillimeter size.

Multipath

Multipath is caused by reflections from surrounding objects. The reflected signals have often lower power, but usually causes errors from a millimeter up to two centimeters, or in the worst case a maximum error of a quarter of a wave length. Diffuse multipath from rough surfaces has reflections with random phase and low power. But specular multipath from smooth surfaces generates coherent waves, whose systematic effect can make the error be mistaken for attitude motion [2].

Multipath is not eliminated for closely spaced antennas. Multipath errors can be dampened by using choke ring antennas, or better receivers. If the antenna site is static then the effects from it can be surveyed, and later corrected. This also works for flying aircraft or satellites where multipath is caused by reflections from the own body. As multipath is a dominant error source it is important to choose a good antenna placement.

Oscillator stability

The oscillator in the receiver is not perfectly stable. Instabilities in the oscillator generates errors in the phase tracking loop. The stability can be adjusted but is weighted against the bandwidth in the tracking loops. These errors are of millimeter size [2].

Satellite position

The calculated satellite positions are not perfectly correct, based upon predictions of up to two days old data the position can be estimated within a couple of meters. However the antennas are usually closely spaced and the error should be the same for both. Thus the satellite position error is removed when forming the single differences.

Antenna phase center variations

Antennas may not have a stable phase center; instead it can vary between different incident directions. The phase center for different frequencies (L1 & L2) does not necessarily match either. If the same type of antennas is used and they are aligned

(29)

3.3 Measurement model 15

then the errors should cancel out. But if the phase center variations differ between individual antennas of the same type, then the error will of course not cancel out. The error could be a few millimeters large, or up to 20 centimeters if the antennas are misplaced [2].

3.3.2 Variance and covariance of phase measurements

The accuracy of phase measurements varies between different receivers and man-ufacturers. Cheaper low end receivers measure with 1 cm or better precision [16], and for high quality receivers figures as low as 0.1 mm [11] are mentioned although 1 mm is more common [17]. Those numbers only apply to perfect conditions and the type of error is not always specified. A more dominant error source is multi-path that usually causes measurement errors of millimetre size but in the worst case scenario can cause errors of a quarter of a wavelength which for the L1 band is close to 5 cm.

The phase measurements of different satellites with the same receiver are partly correlated as the same antenna and signal amplifier is used. But for simplicity they are assumed to be uncorrelated. This is probably a good approximation as receiver noise is a minor error source, compared to multipath.

Example 3.1: Covariance matrix for phase measurements

Assume that four common satellites are tracked by two receivers. The single differences are combined from the phase measurements.

∆Φ = ΦA− ΦB. (3.8)

Assume also for simplicity that the receivers operate synchronously with zero line bias. The error then becomes

∆ = A− B. (3.9)

The measurement errors are assumed independent with zero mean. The co-variance is

R∆ = E{(A− B)(A− B)T}

= E{(ATA+ BTB}

= 2σ2I, (3.10)

where I is the identity matrix.

The single differences are uncorrelated. The double differences are formed by taking the difference between a chosen reference satellite and the other used satellites. This is equivalent to multiply the single difference vector with the matrix

(30)

16 Attitude determination B =   1 −1 0 0 1 0 −1 0 1 0 0 −1  . (3.11) R∇∆= E{B∆∆TBT} = B2σ2IBT = 2σ2BBT (3.12) R∇∆ = 2σ 2   2 1 1 1 2 1 1 1 2   (3.13)

The double differences are indeed correlated.

3.3.3 A simulation model

To initially test the implemented algorithms a simple model that generated sim-ulated phase measurements was created in Matlab. In the model the satellite position and motion is calculated from real ephemerides1 _{taken from the recorded}

data. These positions are used as the truth, and from them the range from observer

k to satellite i, ri_k, and line of sight vector to satellite i, ˆei, can be calculated. At the start of the simulations the exact time, t0, is decided as well as the

exact position of the observer, uk.

The exact distance between satellite and antenna, ri_k, is first determined. Thereafter the integer number of cycles is removed, thus a perfect measurement of the fractional part of the phase is created.

φi_k(to) = rki(t0)%λ (3.14)

The initial number of cycles

N_ki(t0) =

ri

k(t0) − φik(t0)

λ (3.15)

is also stored.

The following measurement epochs2 _{the time is incremented and new satellite}

positions are calculated. The new distances are subtracted by the initial number of cycles, thus any changes of integer size are not lost. The attained measurements are perfect but to make them more realistic white noise and a bias β is added (the bias is different for different observers, otherwise it could be eliminated),

ˆ

φi_k(t) = ri_k(t) − λN_ki(t0) + + β. (3.16)

The simulated phase measurement ˆφ(t)i

k is calculated for each observer and

satel-lite. Together with the line of sight vectors the phase measurements are sent to the algorithm.

1_{Ephemeris: a table of parameters used to describe satellite orbits and clock corrections.} 2_{Epoch: A specific instance in time when measurements are made.}

(31)

3.4 The relation between carrier phase and baselines 17

This simulation model was used for basic debugging of the code, but for eval-uation only real data was used.

3.4 The relation between carrier phase and

base-lines

For simplicity assume that the single differences are used and there is no clock error. The ∆βAB term is then set to zero. If the phase difference including integer

ambiguities is known, then also the relative distance between observers is known with good precision. Equation 3.3 yields,

ri_A− ri B= ∆φ i AB+ λ∆N i AB− ∆AB. (3.17)

The observers are assumed to be close to each other (< 20 meters) and thus the direction to each satellite is approximately the same for both observers. Since the satellites are such far away from the observers this is a reasonable approximation. From Figure 3.1 it is seen that the projection of the baseline onto the line of sight vector equals the phase difference between the observers,

bT_ABˆei= rAi − r i B= ∆φ i AB+ λ∆N i AB− ∆, (3.18)

bAB := baseline from receiver B to A,

ˆ

ei _{:= line of sight vector to satellite i.}

Unknowns for the equations are the three coordinates for the baseline, which can be solved if three or more independent equations are available. This implies that three or more common satellites need to be tracked by both receivers. How-ever in reality four satellites are needed to calculate the master receiver’s position including clock offset, needed to determine the local NED coordinate system and the line of sight vectors. Note that the line of sight vectors have estimation errors as they are derived from predicted satellite positions and measured receiver posi-tions; however the large distance to the satellites motivates the approximation to use them as if they where known.

The equations for several satellites are combined to form an equation system on matrix form.       ˆ e1T ˆ e2T .. . ˆ enT       | {z } :=H b =      ∆φ1 AB ∆φ2 AB .. . ∆φn_AB      | {z } :=∆Φ +λ      ∆N1 AB ∆N2 AB .. . ∆N_ABn      | {z } :=∆N −∆AB Hb = ∆Φ + λ∆N − ∆AB (3.19)

The equivalent for double differences is true, though then the H matrix consist of differenced line of sight vectors,

(32)

Figure 3.2. Extracting attitude from the baseline.

Hb = ∇∆Φ + λ∇∆N − ∆AB. (3.20)

When the integer ambiguities are resolved and four or more equations are available an over-determined equation system is attained, and the baseline can be solved using a regular least squares approach,

HTHb = HT(∇∆Φ + λ∇∆N )

b = (HTH)−1HT(∇∆Φ + λ∇∆N ) (3.21) A weighted least squares estimate can also be used to account for correlation between measurements, which is suitable if the strongly correlated double differ-ences are used. Then the inverse of the covariance matrix is used as the weight matrix,

HTR−1_∇∆Hb = HTR−1_∇∆(∇∆Φ + λ∇∆N )

b = (HTR−1_∇∆H)−1HTR−1_∇∆T (∇∆Φ + λ∇∆N ) (3.22)

.

3.5 Attitude derived from the baseline

The baseline is calculated in a local frame and is used to determine the attitude, but the baseline does not necessarily coincide with the body frame. A rotation matrix from the body fixed frame to the local frame might be needed to extract the attitude. However determining the attitude for a single baseline that coincides with the body fixed frame is simple. The angles are defined in Figure 3.2.

ψ = tan−1(y

x) (3.23)

θ = tan−1(_p z

(33)

3.6 A dynamic state space model 19

ψ := heading

θ := elevation (roll or pitch)

From this it can be seen how errors in the baseline estimate transfers to the attitude angles. Let m denote the maximum measurement error for the relative

positions of the antennas. The displacement m in position gives a maximum

angular error of sin−1 m

L , which is close to m

L for small angles, where L is the

baseline length. Thus the error in attitude is inversely proportional to the baseline length.

3.6 A dynamic state space model

For use with filtering a dynamic model of the problem is suitable,

xn+1 = fn(xn) + Gn(xn)wn, (3.25)

yn = hn(xn) + Jnzn+ vn, (3.26)

where xn ∈ Rn denotes the states which could contain at least the baseline

coordinates or attitude angles, yn ∈ Rn denotes the measurement vector, zn ∈ Zm

denotes the integer ambiguities (note the subscript n which indicates that the integers can vary with time, for example due to cycle slips), and wn∈ Rq together

with vn ∈ Rn denote noise processes. The function fn includes dynamcis and

describes how the current state affects the next, while Gn describes how process

noise affect the states. The function hn describes how the measurements are

related to the states. For example it could be a projection on the line of sight vectors if baseline coordinates were the states (compare to Equation 3.19). Jn

maps the ambiguities to the measurements. For the single difference case using L1 carrier phase measurements the identity matrix would be used, while for the double differences the Jn-matrix could describe the combinations of satellites that

form the double differences.

The model could be used for filtering if the vehicle dynamics is known, but an extension could also be used when additional sensors are available. For example sensor fusion with one or several gyros using an Extended Kalman Filter (EKF) [4]. Other applications are filter banks or particle filters, which could be used for integer ambiguity resolution. An advantage for the particle filter is that the ambiguities can be dealt with as integers directly without any middle stage such as the so called floating ambiguities. However, this is not further elaborated in this thesis, instead the reader is encouraged to view for example the references [4], [5] and [8] for further information.

(34)

(35)

Chapter 4

Ambiguity resolution

The number of whole wavelengths between two antennas is unknown at start of tracking. To determine the total phase difference, which is needed for the attitude solution, the initial integer ambiguities must be found. They can thereafter be held fixed and used for attitude determination until the satellites disappear behind the horizon or a signal is either blocked or a cycle slip occurs. A cycle slip occurs when the receiver fails to register one or several integer cycles, for example due to large measurement errors caused by multipath. Cycle slips can be hard to detect immediately and could potentially corrupt the attitude solution for a period of time, if sufficient error checks are not conducted.

To find the integers different search techniques are used as there is no analytical solution for the problem. The strategy of many methods is to search a limited space of possible candidates and evaluate those using distinguishing tests. There are methods that use a predefined grid of guessed baselines and there are also methods that divide the ambiguities in different sets, where only some of the ambiguities are used in the actual search. Another approach is to collect data for several epochs and solve for the ambiguities directly. Initially a float solution, where the ambiguities are allowed to be real numbers, is attained and thereafter a search is performed around this initial estimate. The search volume is usually based upon the calculated covariance for the float solution. The Least Squares AMbiguity Decorrelation Adjustment method, or just the LAMBDA method, works in this way but in an extra step decorrelates the ambiguities prior to searching around the float solution. The decorrelation step is performed to change the search space from an elliptical one to a more sphere like one, which is easier and faster to search around. This spherical like search space is common when the strongly correlated double difference phase measurements are used.

The choice of method depends partly on what hardware is used. If a dedicated attitude receiver is used the phase measurements could be nearly simultaneous and methods using single differences would be an alternative. Otherwise methods with the double difference of phase measurements are preferred, as they become more general and can be used on many different hardware setups.

(36)

22 Ambiguity resolution

Figure 4.1. 2D example of maximum number of integer cycles that fit

Example 4.1: Brute force search

Assuming that the baseline length is fixed and known, there are only a maximum number of integer cycles that can fit between the antennas. The maximum number of cycles is reached when the baseline is rotated parallel to the specific satellite’s line of sight vector. The maximum number of integers is calculated as the baseline length divided by the carrier wavelength rounded downwards.

Nmax= b

kbk2

λ c (4.1)

The search space for a single ambiguity can be visualized as the two dimensional example in Figure 4.1 where antenna A is held fixed and antenna B, at a distance of the known baseline length, is allowed to rotate.

The baseline can be rotated 180 degrees with the integer being negative, and the baseline can also be rotated in such way that the integer could be zero. The number of candidates for a specific ambiguity thus comes down to 2Nmax+ 1. If

a one meter baseline is used the number of candidates per satellite is 11. Not all combinations of integers are possible when combining several satellites, but if a brute force search was to be performed the total number of candidates for the whole integer set would be (2Nmax+1)pwhere p is the number of ambiguities. Again with

a one meter baseline assume that eight satellites are tracked. The entire search space would then consist of 118 _{= 214358881 candidates, which would require}

considerable computing resources to process. If better accuracy for the attitude solution is wanted a longer baseline would be chosen. If the baseline length is just increased to 2 meters the search space would be increased to 218_{= 37822859361}

candidates. It is evident that a brute force search is not desirable in real time applications, and other techniques that drastically decrease the search space are used instead.

(37)

4.1 Guessed baseline search 23

Figure 4.2. Baseline search space.

4.1 Guessed baseline search

Instead of searching for all combinations of integers the search can be based upon guessed baselines. This technique is mentioned in [4] among other places. Position domain constraints can easily be used as only possible baselines are searched. A car mounted system can be used as an example. As a regular car seldom has large pitch or roll angles it is convenient to only search baselines with small pitch and roll. The search space can be limited to only one dimension of possible heading, or yaw, angles.

The guessed baselines should be equally spaced and not farther away from each other than a single integer difference in phase, but also not too close to each other as that would worsen performance. The angle θ between two baselines is the angle that gives at maximum a whole wavelength of phase difference from one baseline to the other for any given satellite direction. In Figure 4.2 baselines 1 and 2 are plotted with antennas at locations a1, b1 and a2, b2. The requirement that the

phase difference is not larger than λ can be translated to the requirement that the length of the vector b1− b2 is less than or equal to half a wavelength.

kb1− b2k2≤

λ

2 (4.2)

If the known baseline length is d then an isosceles triangle with the angle θ and lengths d/2 and kb1− b2k is formed. The triangle is composed of the vectors

b1, b2and b1− b2. From a geometrical perspective it is given that

d 2sin θ 2 = kb1− b2k2 2 (4.3)

(38)

24 Ambiguity resolution θ = 2 sin−1 kb1− b2k2 2 5 2 sin−1 λ_2d (4.4) Now an upper bound for the angle step between the guessed baselines is at-tained, however it is safer to choose a slightly smaller step if the measurements are noisy. Using half the angle could be a suitable compromise between efficiency and robustness.

Example 4.2

Assume that the baseline length is 0.50 m. The total number of guessed baselines would be

θ = bsin−1λc = 10◦

Ncand=

360

θ = 36

The number of candidates to evaluate is only 36 for a baseline of 0.50 m with the constraint that the pitch and roll angles are small. Observe also that the number of candidates is independent of how many satellites that are tracked and used in the equations. If the pitch and roll angles are not small but known within a few degrees then the search space can be tilted without being extended.

When the search space is created the integers for each baseline candidate can be computed. Remember equation (3.19) where the relationship between the base-line and phase measurements is described. The ambiguities can be calculated by reordering and rounding. The measurement error is not known and therefore set to zero in the equation (the deltas are skipped for convenience).

N = d1

λ(Hbguessed− Φ)e (4.5)

Using the integer set N the baseline can be recalculated and the result will differ from the guessed baseline as the fractional part is now included and also due to the rounding. A test baseline is attained as

btest = H†(λN + Φ) (4.6)

Now a number of test baselines have been calculated and what is left is to figure out which test baseline with corresponding set of integers that is correct. To do this several distinguishing tests are used. These tests are described later in section 4.4.

4.2 Primary and secondary sets

To determine the three baseline coordinates only three ambiguities are required, and the others are for measurement redundancy. According to reference [14] the

(39)

4.3 The LAMBDA method 25

idea to only search the integers for three primary sets was introduced by Ronald Hatch. For the primary set all possible integers could be searched as in the brute force search, but the search space can be made smaller by using the baseline length constraint. The HEADRT +T M software developed at the University of Calgary uses this division of the integers in primary and secondary sets [3].

The secondary integers (subscript s) can be derived from the baseline created from the primary integers (subscript p),

bp = Hp†(λNp+ Φp) (4.7)

Ns = d

1

λ(Hbp− Φs)e, (4.8)

whereby the baseline estimate using all satellites can be made,

N = Np Ns , Φ = Φp Φs , b = H†(Φ + λN ) . (4.9)

4.3 The LAMBDA method

The LAMBDA method developed by P.J.G. Teunissen formulates the problem of finding the baseline as a mixed integer least squares problem [20] [19]. The observation equation is formulated as

y = Aa + Bb + , Cov(y) = Qy (4.10)

and the mixed integer least squares problem is written as min a,b ky − Aa − Bbk 2 Q−1y , a ∈ Zp, b ∈ R3, y ∈ Rp (4.11) where

y := phase measurements, either single or double differences,

a := integer ambiguities,

b := baseline,

A := design matrix for ambiguities,

B := matrix with combinations of line of sight vectors,

:= measurement error,

Qy := variance matrix for y,

p := number of ambiguities.

More generally the measurements can include the pseudoranges as well as the carrier phase measurements.

The function to be minimized is orthogonally decomposed and rewritten as a sum of three squares. The decomposition is attained by inserting the least squares solution of the problem where the integer constraint has been relaxed.

(40)

26 Ambiguity resolution ky − Aa − Bbk2 Q−1y = kêk 2 Qy + kâ − ak 2 Qâ+ kˆb(a) − bk 2 Qˆb(a) (4.12) ˆ

e := the least squares estimate residual, y − Aˆa − Bˆb,

ˆ

a := float ambiguities, ˆ_b _{:= float baseline,}

ˆ_b(a) _{:= conditional least squares solution of b, assuming a is known.}

The residual term cannot be affected and is therefore removed from the mini-mization. The problem is now formulated as

min a∈Zp kˆa − ak2Qˆa+ min_b∈R3kˆb(a) − bk 2 Qˆb(a) . (4.13)

The last term can be made zero for any a as b can be chosen freely and set equal to ˆb(a). Thus when a is found the solved baseline is simply ˆb(a). Finally the

problem is formulated as ˇ a = arg min a∈Zkâ − ak 2 Qâ, (4.14) ˇ_b ₌ ˆ_{b − Q}_ˆ bâQ −1 ˆ a (ˆa − ˇa), (4.15)

where the last equation simply could be rewritten as the least squares estimate for b with ˇa as the fixed ambiguities.

Note that it is the norm weighted by the correlation matrix, Qˆa, that is to be

minimized, the result can therefore differ from simple rounding. The computation of the ambiguities, ˇa, cannot be done analytically; instead a search must be

per-formed, as in the other described methods. However here the search is centered around the float solution and shaped by the covariance of the float solution. If dou-ble difference carrier phase measurements are used, the ambiguities are strongly correlated as mentioned earlier. Then the search space becomes elongated and time consuming to search. The efficiency of the LAMBDA method is gained by first decorrelating the ambiguities. A transformation matrix is calculated, which decorrelates the measurements as much as possible without removing their integer nature.

When the ambiguities are transformed the search is started. It has been proved that for the Gaussian case the integer least squares estimator is optimal, with respect to maximizing the probability of chosing the correct integer set [19].

For the LAMBDA method to work a float solution and its variance must first be calculated. Many measurement epochs of data might be needed for the float solution to be accurate enough that the integers are found. Recently P.J.G. Teu-nissen has proposed a modification of the LAMBDA method that is more suitable for the attitude determination problem. The baseline length constraint is added to the float solution, which improves the reliability. It is stated that this modified version works in the single epoch case when only L1 phase measurements are used, which for the regular method yields too many integer candidates [20].

(41)

4.4 Distinguishing tests 27

4.4 Distinguishing tests

4.4.1 Baseline length check

The baseline length is usually known with the precision of a couple of millimeters up to a couple of centimeters. Candidates with very large length errors can di-rectly be invalidated, but it is a matter of tuning to decide how large the smallest acceptable error should be. With multipath effects the correct solution can have an error of a quarter of a wavelength or so, thus the acceptance bound should be around there, [2].

The length check might not remove all of the correct integer sets immediately, but as the satellite geometry changes it is reasonable to assume that the faulty candidates eventually will show up in the length test. However it is not evident that this will happen as fast as one wants the ambiguity resolution to finish.

4.4.2 Phase residuals

The phase residuals, v, can be defined as the difference between the measured phase and the predicted phase derived from the baseline and satellite directions,

V = Hb − λN − Φ, (4.16)

where V is a vector of size p that contains the residuals for each of the p ambiguities, (each satellite if single differences are used)

V =      v1 v2 .. . vp      .

The measurement errors can be assumed to have zero mean. Therefore the correct integer set should have residuals with zero mean. In reality multipath errors do not have zero mean but they might be averaged out if several epochs of data are used. The residuals should anyway in average be lower for the correct integer set and this can be used to distinguish the correct candidate from the others.

χ2 _test

To combine the residuals for different satellites the sum of squared residuals, VTV ,

is used.

More precisely the residuals are assumed to be statistically normal distributed with zero mean and standard deviation σ. For independent stochastic variables

Xi|i = 1, ..., n their sum of squares is χ2 distributed [6],

1 σ2 n X 1 Xi2∈ χ2(n). (4.17)

(42)

28 Ambiguity resolution

But the double differences are correlated as shown earlier, and a modified version of the expression must be used to account for this, [9],

n X i=1 n X j=1 XiXj σ2 ij (4.18)

which expressed in matrix form is formulated as

VTR−1_∇∆V ∈ χ2(n). (4.19)

A confidence test can be used since the distribution is known. In [22] it is proposed that a confidence level of 75% is used and that the degrees of freedom are set to the number of satellites minus four for double differences, which would mean three for single differences.

For this test to be accurate the standard deviation must be chosen well as a too large value would let too many candidates pass, and a too small value could reject the correct candidate.

The hypotheses are defined as:

H0: VTR−1_∇∆V ∈ χ2(n) ⇒ a probable candidate.

H1: VTR−1∇∆V /∈ χ2(n).

The null hypothesis is rejected if

VTR−1_∇∆V > x. (4.20)

A confidence level of 75% yields

P (VTR−1_∇∆V > x) = 0.25 (4.21)

from which x can be decided using the χ2 distribution.

The phase residuals could be tested in the two ways described above. The correct candidate should pass the χ2 test given that the standard deviation is correct, and the correct candidate should also in average have the lowest residuals.

Ratio test

To speedup the ambiguity resolution another test called the ratio test is popular, [14]. If the residuals are significantly lower for the best set compared to the second best set one can assume that the current best set is the correct one. The test can be defined as VT 1stR −1 ∇∆V1sr VT 2ndR −1 ∇∆V2nd < α. (4.22)

However α must be chosen correctly and it is hard to find a strict mathematical derivation for what value to use. Depending on noise level and multipath the ratio between the best set and the second best set can vary considerably, and it is not even necessarily so that the best set is the correct one. To reduce this difficulty the

(43)

4.4 Distinguishing tests 29

residuals can be accumulated for several epochs and averaged. But α still needs to be chosen, which is a matter of tuning.

The hypotheses are defined as:

H0: The lowest sum of squared residuals is significantly lower than

the others, which makes it probable that it belongs to the correct set

H1: The lowest sum of squared residuals is not significantly lower

than the others.

and the null hypothesis is rejected if equation 4.22 is not fulfilled.

Threshold test

Another test is suggested in [4] where simply the residuals below a chosen constant threshold are assumed to belong to the correct integer set. But then there is the problem to decide what threshold to use. In high noise and multipath environments the residuals can vary considerably.

(44)

(45)

Chapter 5

Implemented algorithms

Based upon the literature study a couple of algorithms were implemented in Mat-lab, and tested using both collected raw data and data from a simple simulation model.

5.1 Logistics

Besides the problem of finding the correct integer set, it is a challenge to find an efficient way of handling events where the number of available tracked satellites changes. It could be events where some of the of satellites temporarily disappear due to signal blockage. It could also be events where new satellites are tracked, or even a combination of new and lost satellites. The simplest, but not the most elegant way to handle this, would be to simply start a new integer search whenever the composition of available satellites changes. This is not a good idea if the search takes long time, and if the environment is such that temporary signal blockage is common, then the integers might never get locked.

However, if a simple and fast algorithm is used to complete an integer search each epoch this is really not a problem. The available measurements for each epoch is used and integers are only searched for those, disregarding what satellites that were tracked the last epoch. But if another method is used it could be a problem, especially if the double differences are used with one satellite as reference, and that particular satellite would be lost from tracking, then all integers would become invalid, and a new search has to be started.

The event of a newly tracked satellite could be handled by working out the new integer directly from the calculated baseline, assuming the baseline is known and the other integers are locked. If it appears in the middle of an integer search then different things can be done. When the search space is based upon guessed baselines the new measurements could be included directly, if it was the guessed baselines that were saved. Even if the baseline is not saved it can be calculated using the measurements and the other integer candidates and the new integer can be estimated thereafter, however normally a rounding operation is performed here and if the result is not perfectly close to an integer before the rounding,

(46)

32 Implemented algorithms

then measurement errors may lead to the selection of the wrong integer. If this is a concern then both the possible integers can be used, and the search space is doubled.

Another way of including the new integer is to expand the search space with all possible entries for the new integer, i.e. for each candidate integer set create new sets with all possible combinations of the added integer. This will imminently make the search space much bigger but the baseline test and residual tests would probably reject many of the new candidates in short time.

A lost satellite could be removed from the equations without affecting the rest of them, unless it is used in combination with other satellites, where also those need to be taken care of.

To handle the events it is necessary to identify which changes that have oc-curred to the satellite composition, and then call the correct routine to deal with them. The code for doing this could easily outgrow the actual algorithm but it depends on how many cases that should be supported before a search is started.

5.2 Cycle slip detection

All phase measurements are checked for cycle slips before sending them to any of the algorithms. The Doppler frequency for the previous epoch together with the current Doppler frequency is used to predict what the phase should be, and if the measurement differs too much from the prediction a cycle slip is assumed found. The integer for the satellite is then reset, and searched for. The predicted phase is calculated as ˆ φ(kT +T )= φ(kT )+ T ˙ φ(kT +T )+ ˙φ(kT ) 2 (5.1)

and the test is

kφ − ˆφk

T < τ (5.2)

where T is the sample time, and τ is a value that must be chosen.

5.3 Guessed baseline search

The guessed baseline search was implemented using both single and double dif-ferences. Based upon the baseline length a search space was created. With short baselines the search space is small enough that a search can be executed at each epoch, however an option of locking to an integer set was added.

The single differences of the carrier phase was used in two ways, either the clock bias including the line bias due to different cable lengths etc. was assumed constant and known, or the bias was assumed unknown and solved for during each epoch. To solve for the bias it was added as an extra unknown together with the baseline. The needed modifications were to add a column of ones in the H-matrix, and an initial guessed, or saved, estimate of the bias in the guessed baseline.

(47)

5.4 Primary and secondary sets 33 H =      e1 x e1y e1z 1 e2x e2y e2z 1 .. . ... ... ... en_x en_y en_z 1      (5.3) .

Now the b vector contains the x, y, z coordinates and the bias

b =     x y z β     (5.4)

and the b vector is solved as a least squares estimate

b = H†(∆Φ + λ∆N ). (5.5)

While the candidates are derived from the guessed baselines a length check is performed and only those candidates that have a small enough error are stored together with their residuals. Later when all candidates are examined they are sorted. If the candidate with the smallest residual passes a ratio test then it is considered the right one. The resolved integers can then be used directly to calculate the attitude, or as an option they can be monitored for several epochs to reduce the possibility of locking to wrong integers. If the ratio test fails, then the problem is unsolved and the search will have to restart when new data arrives.

5.4 Primary and secondary sets

5.4.1 Forming the sets

The baseline can be calculated from three resolved ambiguities, but which ambi-guities and what combination of satellites should be used? To get a good starting estimate the three chosen ambiguities should belong to the best measurements, but it is hard to know exactly which measurements that are the best. In general, however, signals from high elevation satellites travel a shorter distance trough the atmosphere, and are also less probable to be blocked or be subject to multipath. Therefore it would be a good starting point to choose from the satellites with the highest elevation. Beside measurement errors the geometry of the satellites affects the results. The geometric strength can be expressed with the Position Dilution of Precision (PDOP). To calculate the baseline the inverse of the matrix HTH is

multiplied with the measurements. The diagonal elements of HT_{H quantify how}

measurement errors translate to the components of the solution. Let D1, D2, D3

denote the first three diagonal elements, corresponding to the three coordinates in space, then the PDOP is defined as

(48)

34 Implemented algorithms

Using this PDOP measurement the primary set is formed. First the satellites are sorted after elevation and afterwards the combination of the highest elevation satellites that yield the lowest PDOP are chosen.

5.4.2 Forming the search space

The primary integers are searched for all combinations. This is for simplicity, if a speedup is wanted a number of constraints can be used in the search space, for example using Cholesky transformation of the squared line of sight vectors [14]. During the search candidates that have a too large length error are rejected directly. If a candidate passes, the secondary integers are derived from the baseline and a new baseline is calculated. The baseline is subject to a new length test, and if passed the candidate is stored. All of this is performed in the loop where all candidates for the primary set is searched. The procedure follows mostly what is decribed in [14]. Once the search space has been calculated it is stored and reused until either empty or until a candidate is locked.

5.4.3 Distinguishing the candidates

For each epoch the phase residuals for all candidates are calculated and accumu-lated. Also a length check is performed, and if the length error is too large the candidate is removed. Thereafter the size of the search space is checked. If a single candidate remains a χ2 test is performed before it is considered the right one, but it is not locked directly. Instead the phase residuals are observed for another couple of measurement epochs, and then locked if it has passed the tests all times. If the search space is empty a new search will be conducted the next epoch. Otherwise the ratio test is performed and if a candidate passes a χ2 test is also performed, before it is considered the right one. But as with the single candidate case the phase residuals are observed for a couple of epochs before the integers are locked.

GPS-based attitude determination

Institutionen för systemteknik

Department of Electrical Engineering

Examensarbete

GPS-based attitude determination

GPS-based attitude determination

Examensarbete utfört i Reglerteknik

vid Tekniska högskolan i Linköping

av

Abstract

Acknowledgments

Contents

Acronyms

Symbols

Operators and functions

Chapter 1

Introduction

1.1

Background

1.2

Purpose

1.3

Solution method

1.4

Simplifications

1.5

Related work

1.6

Outline

Chapter 2

Introduction to GPS

2.1

The satellite constellation

2.2

The control segment

2.3

The user segment

2.4

GPS signals and measurements

2.5

Principle of GPS positioning

Chapter 3

Attitude determination

3.1

Coordinate systems

3.2

Determining attitude

3.3

Measurement model

3.3.1

Error sources

3.3.2

Variance and covariance of phase measurements

3.3.3

A simulation model

3.4

The relation between carrier phase and

base-lines

3.5

Attitude derived from the baseline

3.6

A dynamic state space model

Chapter 4

Ambiguity resolution

4.1

Guessed baseline search

4.2

Primary and secondary sets

4.3

The LAMBDA method

4.4

Distinguishing tests

4.4.1

Baseline length check

4.4.2

Phase residuals

Chapter 5

Implemented algorithms

5.1

Logistics