# Navigation Facility for High Accuracy Offline Trajectory and Attitude Estimation in Airborne Applications

Authors: A. Renga; G. Fasano; D. Accardo; M. Grassi; U. Tancredi; G. Rufino; A. Simonetti
Journal: International Journal of Navigation and Observation
 (2013)

---

## Abstract

The paper focuses on a navigation facility, relying on commercial-off-the-shelf (COTS) technology, developed to generate high-accuracy attitude and trajectory measurements in postprocessing. Target performance is cm-level positioning with tenth of degree attitude accuracy. The facility is based on the concept of GPS-aided inertial navigation but comprises carrier-phase differential GPS (CDGPS) processing and attitude estimation based on multiantenna GPS configurations. Expected applications of the system include: (a) performance assessment of integrated navigation systems, developed for general aviation aircraft and medium size unmanned aircraft systems (UAS); (b) generation of reference measurements to evaluate the flight performance of airborne sensors (e.g., radar or laser); and (c) generation of reference trajectory and attitude for improving imaging quality of airborne remote sensing data. The paper describes system architecture, selected algorithms for data processing and integration, and theoretical performance evaluation. Experimental results are also presented confirming the effectiveness of the implemented approach.


---

## Body

## 1. Introduction

Field test refers to the testing of a device or sensor in the conditions under which it will be actually used. Field testing becomes necessary when numerical and indoor testing may fail, that is, in all the cases in which the operative conditions are difficult to be reproduced with high fidelity by software or laboratory simulations.A typical example is the performance assessment of integrated navigation systems for airborne applications. In this case several factors must be taken into account, including but not limited to: (a) the intrinsic dynamical and statistical models of sensor and systems, (b) the selected data fusion strategy, and (c) the typical and the worst case manoeuvres that must be considered to determine a trustworthy dynamical model of the aircraft where the system is installed. Indeed, high performance heading and attitude determination units are needed both in general aviation and in unmanned aircraft systems (UAS) applications to attain an adequate control performance. UAS are more demanding in terms of attitude determination performance than manned aircraft for a series of issues. First, the absence of a human pilot onboard prevents the aircraft from using the human senses such as vision and coordination to integrate onboard systems. Furthermore, UAS are required to perform autonomous flight in case of loss of data link with the ground control station where the remote pilot is located [1, 2]. This condition may have to be managed even during the most critical phases of the missions, such as takeoff and landing. In these phases, navigation performance requirements are typically tighter than during cruise flight.In the latest years, the technological evolution of microelectromechanical system (MEMS) inertial sensors has been characterized by increasing accuracy and reliability, so that they are finding increasing usage in high-accuracy integrated navigation systems which also include global navigation satellite system (GNSS) receivers and magnetic sensors. Their most recent versions can attain accuracy levels that are close to the one of expensive fiber optic gyros (FOGs). Indeed, the concurrent adoption of payloads installed in strap-down configuration, integrated navigation systems based on MEMS sensors and GPS guarantees a significant cost reduction with the same levels of performance for many UAS.As MEMS accuracy improves, it becomes more and more challenging to verify navigation performance of the overall integrated system using commercial-off-the-shelf (COTS) components. Nonlinear dynamical models, non-Gaussian sensor error models, and the wide range of operative conditions cannot be simulated with the required fidelity level. Consequently, field testing is the best option to validate the developed navigation units. In this case, a proper selection of the maneuvers to be executed during the flight missions allows determining a good coverage of typical and worst case conditions.The generation of reference measurements to evaluate the performance of flight sensors, such as radar or laser sensors, is another research area requiring field testing. In this concern, the flying platform is supposed to embark one or more sensors that must be used in specific flight conditions, such as autonomous landing or sense and avoid [3]. The case of performance assessment of planetary landing radars [4] represents an additional example of a similar problem.Based on the above discussion, this paper presents a navigation facility, relying on COTS technology, under development at the Department of Industrial Engineering, University of Naples, aimed at generating offline high-accuracy attitude and trajectory measurements for airborne applications. The navigation facility is based on the concept of GPS-aided inertial navigation and offline processing of sensor data. The system fuses carrier-phase differential GPS (CDGPS) and multiantenna GPS techniques with FOG and MEMS measurements. In line with the expectations of airborne field tests for general aviation and medium-small size UAS applications [5, 6], the navigation facility shall fulfill the requirement of cm-level positioning accuracy and attitude estimation error better than 0.15°. Such a performance can potentially widen the range of applications of the developed navigation facility to sectors different from field testing. For instance, the capability of delivering accurate reference trajectories and attitude histories is required for enhancing the products and the imaging quality [7] of airborne remotely sensed data, for example, synthetic aperture radar or electro-optical data, both in near-real time and in postprocessing applications. Indeed, it is worth noting that an attitude measurement accuracy of the order of less than 10−1 deg is achievable by high-performance inertial systems based on ring laser gyros, but size, weight, and cost budgets are in these cases incompatible with installation onboard light aircraft and medium/close range UAS.The paper is organized as follows. First, system architecture, components, and sensors are presented (Section2). Then the selected algorithms for CDGPS, attitude estimation by multiantenna GPS, and integration of GPS-based measurements with inertial ones are described in detail (Section 3). Section 4 derives the theoretical performance of the navigation facility through simulated scenarios and covariance analysis. Finally, the results of the first experimental tests are discussed in Section 5.


## 2. System Architecture

GPS-aided inertial navigation [8] exploits the features of GPS and inertial sensors to achieve performance exceeding that of the individual sensors. It is well known that inertial sensors suffer from poor long-term accuracy without calibration, whereas the limitations of GPS include integrity or interruption of the signal and sensitivity to interferences from external sources. In GPS-aided inertial navigation, long-term accuracy is guaranteed because GPS limits the navigation error and is also used to calibrate inertial sensors [9]. The developed navigation facility implement a postprocessing version of this principle with important distinguishing characteristics. The integration between GPS and INS is not performed assuming a single receiver but considering differential GPS techniques and exploiting high-accuracy carrier-phase measurements. As a consequence, the derivation of position, velocity, and attitude of the moving vehicle is based on a series of GPS and inertial measurements collected by the sensors mounted on the vehicle, completed by GPS data gathered by one base station. Specifically, the application of CDGPS techniques between the moving vehicle and the base station makes it possible to generate high accuracy, that is, cm-scale, debiased pseudorange measurements for the moving vehicle. In addition, FOGs and MEMS accelerometers are used to providing high-frequency navigation estimates, also during eventual losses of high-accuracy GPS solution, and a multiantenna configuration is adopted to achieve GPS-based attitude estimation, that is, high-accuracy heading measurements, primarily. Finally, an extended Kalman filter (EKF) integrates those measurements to generate the accurate navigation solution.For the sake of completeness, it is important to note that alternative methods exist to implement GPS-aided inertial navigation. These solutions include the following.(1)
Standalone GPS, that is, a single GPS antenna mounted on the moving vehicle.(2)
Code-only differential GPS (DGPS), that is, a single GPS antenna mounted on the moving vehicle and a base station antenna. The differential solution can be provided in real time if a radio link with the base station is implemented.(3)
Real-time kinematic (RTK), that is, single or multiantenna systems mounted on the moving vehicle and connected to the base station through a radio link, with exploitation of carrier-phase differential measurements.Considering an inertial unit of the same class of the one selected for the system presented in the paper, the first and the second solutions can guarantee a pitch and roll accuracy worse than 0.5° and a position accuracy on meter and decimeter scale, respectively. In addition, they cannot accurately measure the heading angle unless significant horizontal accelerations are present. Using magnetic sensors can avoid the need for such acceleration levels by providing magnetic heading information but with a degraded accuracy due to local variations in Earth’s magnetic field, modeling inaccuracies, and external magnetic sources. RTK schemes could be of interest because they can achieve the desired position and attitude accuracy. Anyway, they do not work without a radio link between the vehicle and the base station or when the link is lost. Moreover, RKT operation is not reliable when there is a significant altitude difference, larger than 500 m, between the vehicle and the based station such as in the considered airborne applications. Indeed, in this case the integer resolution can fail due to the differential tropospheric delays (see Section3.1). In conclusion, these methods cannot be considered a valid alternative because they are characterized by unsatisfactory performance for the intended applications.According to these considerations, the system architecture reported in Figure1 has been selected. The basic hardware components of the vehicle are the inertial unit, the GPS receivers and three GPS antennas. The adopted inertial measurement unit (IMU) is the IMU700CA200 manufactured by Crossbow. It is equipped with three MEMS accelerometers that are based on differential capacitance, and three FOGs. The GPS receivers are the models BD960 and BD982 manufactured by Trimble. Both receivers are able to provide dual-frequency carrier-phase and pseudorange observables from the connected antennas. In more detail, BD960 is a single antenna receiver, while the BD982 is a double antenna receiver which can provide raw GPS observables from one antenna and the baseline vector in the North-East-Down reference frame between the two antennas [10]. The GPS antennas mounted on the vehicle are the model AV59, manufactured by Trimble and characterized by an aviation type of design, that is, including rugged package and bulkhead mounting, suitable for vehicle applications. The two receivers and the IMU are connected to a data acquisition computer by means of dedicated RS-232 links. GPS data can be saved at a frequency up to 20 Hz (50 Hz for the BD982 receiver). The IMU output frequency is about 200 Hz. The data acquisition software saves all the data with a time tag based on the CPU clock. This time tag allows the synchronization among GPS and inertial data in the postprocessing phase to be performed.Figure 1
System hardware architecture.Concerning the base station, it is equipped with a BD960 GPS receiver with its antenna, model GA810 by Trimble, and a computer which is used to saving all the GPS observables as it is done on the vehicle. The GPS antenna of the base station is static, as a consequence its position can be considered known with very high accuracy, because either the antenna is placed at a known position or proper base station initialization techniques are applied, as in GPS surveying. Finally, it is worth noting that since carrier-phase differential GPS processing is carried out offline no radio links are foreseen between the vehicle and the base.The nominal performance parameters for all the introduced sensors and components are reported in Table1.Table 1
Performance parameters of the selected sensors and components.
Sensor/component
Performance parameter
Value
BD982 GPS receiver
Carrier-phase noise
<1 mm
AV59 GPS Antenna
Phase center stability
sub-centimeter
BD960 GPS receiver
Carrier-phase noise
<1 mm
GA810 GPS Antenna
Horizontal phase stability
<0.4 cm
Vertical phase stability
<1.4 cm
700CA-200 IMU
Gyro bias stability (constant temperature)
<20°/h
Angle random walk
<0.4°/hr1/2
Acceleration bias
<±8.5 mg
Velocity random walk
<0.1 m/s/hr1/2


## 3. Algorithms Selection

According to the system architecture presented in the previous section, the navigation software is based on three main steps:(i)
offline derivation of high accuracy debiased pseudorange measurements for the moving platform by CDGPS;(ii)
offline derivation of attitude measurements for the moving platform by multiantenna CDGPS;(iii)
integration of GPS data and FOG/MEMS measurements by EKF.Such a multistep scheme works properly only if the derived GPS-based debiased pseudoranges and attitude estimates do not have significant time correlation, and they are sufficiently uncorrelated with one another. As shown in the following subsections, the desired uncorrelation properties result from the kinematic point-wise scheme developed to generate both pseudorange observables and attitude measurements for the moving vehicle by differential techniques. Alternatively, a different, “single-step”, scheme could be designed in which the Kalman filter processes raw pseudorange and carrier observables of both base station and moving platform together with FOG/MEMS data. This approach has the advantage of taking into account all the correlations but has the drawback of requiring very complex handling of carrier-phase ambiguities and biases, involving ambiguity validation [11] along with residual bias computation [12]. The selected multistep approach is simpler to implement, relying on more assessed and documented algorithms, but also robust because it takes full advantage of postprocessing and offline ambiguities and bias estimation.


### 3.1. CDGPS

The capability to generate accurate, cm-scale, debiased pseudorange measurements for the moving platform is based on the possibility to exploit the integer nature of carrier-phase double difference (DD) ambiguities. To this end, dual frequency carrier-phase and pseudorange measurements are collected by both the base station and the moving vehicle. Those indifferent observables are then combined to calculate single and DD observables. Specifically, one double difference is formed by subtracting two single difference equations of the same type and frequency, taken by the same two GPS receivers at the same time from two different GPS satellites, one of which, named pivot, is taken as a reference. DD observables represent the main input to carrier-based differential GPS processing (see Figure2). The first step is a recursive least-square estimation based on the following nonlinear DD observation model [13]:
(1)(P1)SMjk=ρSMjk+TSMjk+ISMjk+(σP1)SMjk,(P2)SMjk=ρSMjk+TSMjk+η-2ISMjk+(σP2)SMjk,(L1)SMjk=ρSMjk+TSMjk-ISMjk+λ1(N1)SMjk+(σL1)SMjk,(L2)SMjk=ρSMjk+TSMjk-η-2ISMjk+λ2(N1)SMjk-λ2(NWL)SMjk+(σL2)SMjk,
where P, L indicate pseudorange and carrier-phase observables, σ indicates measurement noise, the subscripts S and M indicate the base station and the moving platform, the subscripts 1, 2, and WL indicate GPS L1 and L2 signal carrier frequencies and wide lane combination, respectively, the superscripts j and k are the pivot satellite and a generic GPS satellite, ρ is the DD geometric term depending on the positions of GPS satellites, r, of the base station, rb, and the baseline vector B from the base station to the moving platform,
(2)ρMSjk=∥rk-(rb+B)∥-∥rk-rb∥-∥rj-(rb+B)∥+∥rj-rb∥,
in addition, T is the tropospheric delay, I is the ionospheric delay, η is the ratio between L1 and L2 wavelengths, λ is the wavelength, and N the integer ambiguity.Figure 2
Block diagram of the implemented CDGPS scheme.Apart from the measurement noise, the main error sources affecting DD observables include GPS ephemeris errors, ionospheric errors, and tropospheric errors. According to the analysis reported in [13], differential ephemeris errors can be neglected if IGS product [14] is used. In addition, differential ionospheric delays are of the order of a few centimeters for baseline up to 50 km. Since the goal of the recursive least-square algorithm is to generate a float estimation of the integer cycle ambiguities, DD ionospheric errors can be also neglected. Any residual differential ionospheric delays can be easily cancelled by measurement combination after integer ambiguity resolution [15].Differential tropospheric delays, instead, deserve more attention. Actually when the moving vehicle is an airplane important vertical separations, up to a few kilometers, can be established between the base station and the moving platform, thus, increasing DD tropospheric delays to decimeter scale [13]. In this case tropospheric errors must be accounted for to correctly estimate the integer ambiguities, and different strategies can be used. In detail, when meteorological measurements are collected simultaneously to GPS ones, that is, both the base station and the moving vehicle are equipped with proper weather stations, tropospheric delays can be measured and the relevant GPS observable compensated. For the sake of reducing the hardware complexity of the moving vehicle, the utilization of meteorological data is not foreseen in the framework of the present activity. Nonetheless, tropospheric models exist that are able to estimate the delays without meteorological sensors. Based on the accuracy and the reliability, UNB3 troposphere model [16] was selected. The model exploits temperature and humidity profiles contained in the US standard atmosphere to derive five meteorological parameters relevant to the current user position and the day of year. Then, Saastamoinen’s model [17] can be used to calculate wet and dry zenith tropospheric delays (ZTD) that can be related to the slant tropospheric delays through Niell’s mapping function [18]. UNB3 model completed with Niell’s mapping function can guarantee a few centimeter accuracy in the prediction of the tropospheric delay [19].According to this consideration, the DD observation model, after tropospheric delay compensation, can be simplified to(3)(P1)SMjk=ρSMjk+(σP1)SMjk,(P2)SMjk=ρSMjk+(σP2)SMjk,(L1)SMjk=ρSMjk+λ1(N1)SMjk+(σL1)SMjk,(L2)SMjk=ρSMjk+λ2(N1)SMjk-λ2(NWL)SMjk+(σL2)SMjk.Thus, the recursive least-square estimation is applied considering the following state, measurement vectors:(4)xGPS∈ℝ(3+2p)×1,xGPS=[B(N1)SMjk(NWL)SMjk],yGPS∈ℝ(4p)×1,yGPS=[(P1)SMjk(P2)SMjk(L1)SMjk(L2)SMjk],
and modeling [20]
(5)KtGPS=Pt-1GPS(HtGPS)T[HtGPSPt-1GPS(HtGPS)T+RtGPS]-1,xtGPS=xt-1GPS+KtGPS(ytGPS-HtGPSxt-1GPS),PtGPS=(I-KtGPSHtGPS)Pt-1GPS,
where t indicates a generic time instant, I is the identity matrix, HGPS is the Jacobian of the observation function reported in (3), RGPS is the covariance matrix of yGPS, and PGPS is the covariance matrix of xGPS.Recursive application of (5) allows data accumulation over multiple time epochs to be performed also considering correlation between different integer ambiguities. It is worth noting that (5) must be implemented taking into account changes in the observed GPS satellites and therefore, reorganizing and reordering both xGPS and PGPS when new GPS satellites are acquired [21] or pivot satellite changes [22]. The estimated float ambiguities and the accompanying covariance matrix are fed to LAMBDA [23] that derives an integer estimate of the ambiguities.Assuming that these ambiguities are correct, carrier-phase measurements on both frequencies can be debiased, resulting in unambiguous pseudorange measurements, in which the largest remaining error source, when significant vertical separations are present, is due to residual tropospheric delays. In [24] it is shown that this effect can be modeled as a residual ZTD, and the latter can be calculated together with the baseline between the base station and the moving vehicle as the solution of a point-wise kinematic least-square algorithm processing the debiased carrier-phase observables. This approach is implemented in this work.The main problem with the selected algorithm is that the integer estimates derived by LAMBDA are not deterministic but stochastic variables, so there is in general a nonnull probability to derive wrong integer estimates. For increasing the chances to avoid using an erroneous estimate, integer validation tests can be applied to discriminate estimates that are likely to be correct from others that are probably wrong. The prominent integer validation tests in the open literature, described in great detail in [11], can be termed as vector tests, because they screen the whole vector of ambiguities. In case only a few ambiguities in the vector are erroneous, they discard the whole vector. These tests are based on the assumption of debiased Gaussian float ambiguity distribution, and they test the closeness of the float solution to its nearest integer vector. Under this assumption, they provide a theoretical framework for deciding which integer estimates are more likely to be correct and to filter out erroneous ones. The most common vector test is the so-called ratio test [11] that has been recently revised [25, 26] to achieve constant and fixed failure rate. The design failure rate value is chosen as a tradeoff between the requirement of not fixing erroneous integer estimates and the need for fixing the highest possible number of correct ones. Even though the acceptable failure rates depend on the application at hand, typical values are in the order 0.1–1%. This means that, independently from the actual LAMBDA estimates, the test will necessarily discard an amount of integer estimates. In other words, in all the conditions in which the failure rate that can be achieved by the combination of the selected GPS models and the LAMBDA method is lower than that of the validation test, the latter becomes of limited practical usefulness and can be avoided with no performance losses.GPS models can be classified as strong or weak based on the failure rate they achieve. The selected dual-frequency measurement model of (1)–(3), fed with measurements collected by high performance carrier-phase GPS receivers and subcentimeter phase center stability antennas, falls in the class of the strong models when baselines up to several tens of km are considered. The experimental results presented in [25, 26] clearly demonstrate that the ratio test can be discarded for strong models even in the extreme case of epoch-by-epoch (EBE) ambiguity estimation. In EBE the integer ambiguities are calculated every epoch anew. This approach is less conservative and robust than that selected one in this work in which the ambiguities are accumulated over many time epoch, and correlations among the ambiguities and the other parameters of the state vector are exploited; hence, the accumulation is able to cut down the percentage of wrong estimates with respect to EBE, notably [25–27]. In conclusion, the accumulation achieved by recursive least-square estimation and the strong GPS model are expected to reduce the failure rates to values that make standard vector tests useless. In such a condition, the detection of the residual wrong integers can be based on checking individual measurement residuals after residual ZTD computation. When measurement residuals show significant biases, it is likely that one or more integer ambiguities are wrong. As a consequence, the computation of the residual ZTD is performed again considering only a subset of the integer estimates. The procedure is repeated until a sufficient number of integer ambiguities can be assumed valid and correct.The result is an accurate estimation of the baseline vector between the base station and the moving vector that can be used to calculate debiased high-accuracy pseudorange measurements for the moving vehicle as(6)PPRi=∥rk-(rb+B)∥,k=1⋯s,
where s is the number of validated integer ambiguities.In conclusion, it is important to point out that the generated pseudorange measurements can be considered noontime correlated. Indeed, information relevant to multiple time instants is used only to estimate the integer ambiguities. However, once fixed ambiguities are validated they can be considered as known deterministic, and hence uncorrelated, parameters. Moreover, the derived baseline vector and residual ZTD are also time uncorrelated because they are the result of kinematic point-wise least-square estimation. Actually, the sole source of time correlation is due to receiver code and carrier tracking loops that is typically negligible for present generation high-performance GPS receivers.


### 3.2. Attitude Estimation by Multiantenna CDGPS

Multiple GPS antennas must be properly mounted on the moving vehicle to perform GPS-based attitude determination. In the present system assembly, GPS measurements from one antenna of the triplet are also those used to calculate the position of the moving vehicle with respect to the base station. This antenna is called “master” in what follows, while the other two antennas are named “slave”. GPS-based attitude determination strongly relies on the principle of CDGPS. Actually, the baseline vectors between the GPS antennas embarked on the moving vehicle are accurately known in a local body reference frame (BRF), but they can be derived with very high accuracy by CDGPS in Earth-centred reference frame (ECEF) or equivalently in North-East-Down (NED) reference frame with origin at the center of mass of the moving platform. The rotation matrix from BRF to NED can be then calculated, thus, leading to a measure of the attitude angles.CDGPS processing for attitude determination is conceptually similar to that presented in the previous section to estimate the baseline vector between the moving vehicle and the base station. More in detail, DD observables can be derived by combining GPS observables collected by couples of GPS antennas embarked on the moving vehicle, and the baseline vectors between those couples can be then estimated. When addressing the attitude determination problem the distance between GPS antennas is limited to very short values, several meters at most. In such a case, the atmospheric effects, that is, ionospheric and tropospheric delays, are completely cancelled when forming DD observables. Thus, the observation model represented by (3) can be applied. carrier-phase integer ambiguities have to be calculated for each couple of GPS antennas. As detailed in the previous subsection, the integer estimation scheme must include (1) derivation of a float solution, (2) estimation of the corresponding integer solution, and (3) validation of the integer estimate. The validation phase is again based on checking the measurement residuals, but it also takes advantage from the fact that the geometry of the antenna array is rigid and known. Indeed the rigidity constraint reduces notably the number of potential valid integer ambiguities for each couple of measurements [28]. Alternatively, the rigidity constraint could be included in the integer estimation step. In this regard, a multivariate constrained version of the LAMBDA method has been proposed [29] as a reliable solution for singleepoch, real time, and attitude estimation. With specific reference to the test system presented in this paper, however, since no real time operation is foreseen, rigidity constraint is only considered in the validation phase of integer ambiguities.The resulting relative positioning accuracy between the antennas, after ambiguity validation and compensation, is dominated by carrier-phase noise and GPS satellite geometry. It is, therefore, reasonable to assume a horizontal baseline error ranging from 2.5 to 5 mm and a vertical baseline accuracy between 5 mm and 1 cm.Once the antenna baseline vectors have been computed in the NED reference frame, different methods can be used to estimate vehicle attitude. For instance, the same methods adopted for spacecraft attitude determination, such as TRIAD or QUEST, can be used [30, 31]. An interesting solution is represented by the implicit least-squares approach [32], that is, an iterative procedure where both the baseline vectors in NED and in BRF are treated as observations with proper covariance matrices, so that baseline uncertainties in BRF and in NED can be both considered. Starting from a first tentative solution, the method works by successive corrections until a convergence criterion is satisfied. In addition, this approach also computes the covariance matrix for the least-squares solution, which can be straightforwardly included in the navigation filter, and therefore, it has been selected herein for implementation.The first tentative solution can be computed by a suboptimal approach such as the direct computation method [32]. This is briefly discussed in the following because it allows deriving a basic error budget for attitude measurements and understanding the antenna configuration impact on the accuracy of attitude determination. Let us consider the antenna geometry depicted in Figure 3. The master antenna is located in the origin of BRF. In this way, we can neglect the distance of the master antenna from the vehicle center of mass. The baseline from the master to the first slave antenna is along the vehicle longitudinal axis (body axis 1), while the baseline from the master antenna to the second slave antenna is along the body axis 2. Thus, the antenna geometry is such that the two baselines are orthogonal and of equal length L. This is an optimal installation geometry as discussed in [33].Figure 3
Antennas configuration geometry.Given this geometry, the baseline components in BRF and in NED are as follows:(7)b_1=[L00]b_2=[0L0]u_1=[n1e1d1]u_2=[n2e2d2].
By definition, the attitude matrix M converts vectors from NED to BRF.Consider(8)b_1=Mu_1,b_2=Mu_2.Considering the rotation sequence 3-2-1 (headingγ-pitch β-roll α) and using (8), it is possible to derive the following expressions for the attitude angles:
(9)γ=atan(e1n1)+k1π,β=-atan(d1n12+e12),α=asin(d2Lcosβ)+k2π,
where in general k1 and k2 can equal 0 or 1 depending on the sign of NED baseline components.These equations can be used to calculate a first-order estimate of the uncertainty on the attitude angles, under the assumption of uncorrelated errors. For example, heading uncertainty can be estimated as follows:(10)σγ2≅(∂γ∂n1)2σn12+(∂γ∂e1)2σe12.Explicit calculation of the derivatives shows that the propagation of CDGPS errors into the attitude solution depends on the actual vehicle attitude. As an example, Figure4 depicts heading and pitch uncertainty as a function of baseline length L, considering γ = 45°, β = 5°, and α = 1°. Different values for the CDGPS errors are considered which are based on the previous discussion. It can be seen that availability of mm-level accuracy on the NED components makes it possible to achieve about 0.1° accuracy even with relatively short baselines. Pitch (and roll) error is larger than heading due to the effect of vertical accuracy, antenna geometry, and vehicle attitude.Heading (a) and pitch (b) uncertainty as a function of antenna baseline.
(a)
(b)Finally, it is worth noting that in order to ensure accurate GPS-based attitude determination, after system installation, a calibration procedure has to be implemented to calculate with high accuracy the BRF components of baselines among the antennas.


### 3.3. EKF

The integration among the accurate pseudoranges generated by CDGPS algorithm, the attitude measurements delivered by the multiantenna GPS configuration and the inertial measurements are performed by extended Kalman filter. The state vector of the EKF include 9 components which represent estimated error on vehicle absolute position, that is, latitudeϕ, longitude ψ, and altitude h, with respect to the WGS-84 ellipsoid, ground velocity, that is, North vn, East ve, and Down vd components, and attitude with respect to the so-called geographic or navigation frame [20].Consider(11)δx=[δpδvρ],δp=[δϕδψδh]T,δv=[δvnδveδvd]T,ρ=[εnεeεd]T.The origin of this frame is the projection of vehicle position onto the WGS84 ellipsoid, the geographicz-axis points toward the interior of the ellipsoid, and x-axis points toward north. The y-axis completes the orthogonal right-handed reference frame pointing east. The attitude vector error ρ comprises the three small-angle rotations from the true geographic frame to that calculated from the computed state vector [20].The state vector can be propagated through the following linearized model that can be derived from INS mechanization equations:(12)[δp˙δv˙ρ˙]=F[δpδvρ]+[00-Rbn00-Rbn][δfδω],
where 0 is the null matrix, Rbn is the rotation matrix from the body to the navigation reference frame, f is the vector of specific forces measured by the accelerometer, and ω is the vector of relative angular rate measurements. The expression of the matrix F is not reported herein for brevity, but it can be found in [20] to which we refer the interested reader.The observation model of the EKF accounts for the contributions of CDGPS and multiantenna GPS. Specifically, the measurement vector, representing the difference between the true measurements and the measurements estimated from the propagated state vector can be written as(13)δy∈ℝ(s+3)×1,δy=[δPPR1⋮δPPRsδε^nδε^eδε^d],
where ε^ is the attitude rotation measured by multiantenna GPS. As a consequence the linearized observation model is
(14)δy=HEKFδx,HEKF∈ℝ(s+3)×9,HEKF=[Hs0000000I3],Hs=[h1RneD⋮hsRneD],
where Ip is p-dimensional identity matrix, hk is the line-of-sight vector in the ECEF reference frame from the kth GPS satellite to the vehicle, Rne is the rotation matrix from the navigation frame to ECEF, and D is the matrix converting polar coordinates (i.e., ϕ, ψ, and h) to linear coordinates within the navigation frame.Starting from the initial condition, standard nonlinear INS mechanization equations are used to propagate vehicle position, velocity, and attitude, whereas the linear model of (12) propagates the navigation error. When CDGPS and multiantenna GPS measurement are available, EKF estimation theory is applied considering the dynamic model of (12) and the measurement model of (14) leading to an estimate of the measurement error, that is, then subtracted to vehicle position, velocity, and attitude to provide the final GPS-aided accurate estimate (see Figure 5). Finally, concerning filter tuning, process noise matrix can be calculated on the basis of IMU specifications, while measurement covariance matrix is built on the basis of the residuals estimated by GPS processing and the attitude angles uncertainties estimated in the least-squares solution. In this way, less accurate pseudorange and/or attitude measurements have a smaller impact on the accuracy of the final navigation solution. The filter solution is provided at the IMU output frequency, thus, enabling a significant increase with respect to GPS data availability.Figure 5
Block diagram of the aiding process.


## 3.1. CDGPS

The capability to generate accurate, cm-scale, debiased pseudorange measurements for the moving platform is based on the possibility to exploit the integer nature of carrier-phase double difference (DD) ambiguities. To this end, dual frequency carrier-phase and pseudorange measurements are collected by both the base station and the moving vehicle. Those indifferent observables are then combined to calculate single and DD observables. Specifically, one double difference is formed by subtracting two single difference equations of the same type and frequency, taken by the same two GPS receivers at the same time from two different GPS satellites, one of which, named pivot, is taken as a reference. DD observables represent the main input to carrier-based differential GPS processing (see Figure2). The first step is a recursive least-square estimation based on the following nonlinear DD observation model [13]:
(1)(P1)SMjk=ρSMjk+TSMjk+ISMjk+(σP1)SMjk,(P2)SMjk=ρSMjk+TSMjk+η-2ISMjk+(σP2)SMjk,(L1)SMjk=ρSMjk+TSMjk-ISMjk+λ1(N1)SMjk+(σL1)SMjk,(L2)SMjk=ρSMjk+TSMjk-η-2ISMjk+λ2(N1)SMjk-λ2(NWL)SMjk+(σL2)SMjk,
where P, L indicate pseudorange and carrier-phase observables, σ indicates measurement noise, the subscripts S and M indicate the base station and the moving platform, the subscripts 1, 2, and WL indicate GPS L1 and L2 signal carrier frequencies and wide lane combination, respectively, the superscripts j and k are the pivot satellite and a generic GPS satellite, ρ is the DD geometric term depending on the positions of GPS satellites, r, of the base station, rb, and the baseline vector B from the base station to the moving platform,
(2)ρMSjk=∥rk-(rb+B)∥-∥rk-rb∥-∥rj-(rb+B)∥+∥rj-rb∥,
in addition, T is the tropospheric delay, I is the ionospheric delay, η is the ratio between L1 and L2 wavelengths, λ is the wavelength, and N the integer ambiguity.Figure 2
Block diagram of the implemented CDGPS scheme.Apart from the measurement noise, the main error sources affecting DD observables include GPS ephemeris errors, ionospheric errors, and tropospheric errors. According to the analysis reported in [13], differential ephemeris errors can be neglected if IGS product [14] is used. In addition, differential ionospheric delays are of the order of a few centimeters for baseline up to 50 km. Since the goal of the recursive least-square algorithm is to generate a float estimation of the integer cycle ambiguities, DD ionospheric errors can be also neglected. Any residual differential ionospheric delays can be easily cancelled by measurement combination after integer ambiguity resolution [15].Differential tropospheric delays, instead, deserve more attention. Actually when the moving vehicle is an airplane important vertical separations, up to a few kilometers, can be established between the base station and the moving platform, thus, increasing DD tropospheric delays to decimeter scale [13]. In this case tropospheric errors must be accounted for to correctly estimate the integer ambiguities, and different strategies can be used. In detail, when meteorological measurements are collected simultaneously to GPS ones, that is, both the base station and the moving vehicle are equipped with proper weather stations, tropospheric delays can be measured and the relevant GPS observable compensated. For the sake of reducing the hardware complexity of the moving vehicle, the utilization of meteorological data is not foreseen in the framework of the present activity. Nonetheless, tropospheric models exist that are able to estimate the delays without meteorological sensors. Based on the accuracy and the reliability, UNB3 troposphere model [16] was selected. The model exploits temperature and humidity profiles contained in the US standard atmosphere to derive five meteorological parameters relevant to the current user position and the day of year. Then, Saastamoinen’s model [17] can be used to calculate wet and dry zenith tropospheric delays (ZTD) that can be related to the slant tropospheric delays through Niell’s mapping function [18]. UNB3 model completed with Niell’s mapping function can guarantee a few centimeter accuracy in the prediction of the tropospheric delay [19].According to this consideration, the DD observation model, after tropospheric delay compensation, can be simplified to(3)(P1)SMjk=ρSMjk+(σP1)SMjk,(P2)SMjk=ρSMjk+(σP2)SMjk,(L1)SMjk=ρSMjk+λ1(N1)SMjk+(σL1)SMjk,(L2)SMjk=ρSMjk+λ2(N1)SMjk-λ2(NWL)SMjk+(σL2)SMjk.Thus, the recursive least-square estimation is applied considering the following state, measurement vectors:(4)xGPS∈ℝ(3+2p)×1,xGPS=[B(N1)SMjk(NWL)SMjk],yGPS∈ℝ(4p)×1,yGPS=[(P1)SMjk(P2)SMjk(L1)SMjk(L2)SMjk],
and modeling [20]
(5)KtGPS=Pt-1GPS(HtGPS)T[HtGPSPt-1GPS(HtGPS)T+RtGPS]-1,xtGPS=xt-1GPS+KtGPS(ytGPS-HtGPSxt-1GPS),PtGPS=(I-KtGPSHtGPS)Pt-1GPS,
where t indicates a generic time instant, I is the identity matrix, HGPS is the Jacobian of the observation function reported in (3), RGPS is the covariance matrix of yGPS, and PGPS is the covariance matrix of xGPS.Recursive application of (5) allows data accumulation over multiple time epochs to be performed also considering correlation between different integer ambiguities. It is worth noting that (5) must be implemented taking into account changes in the observed GPS satellites and therefore, reorganizing and reordering both xGPS and PGPS when new GPS satellites are acquired [21] or pivot satellite changes [22]. The estimated float ambiguities and the accompanying covariance matrix are fed to LAMBDA [23] that derives an integer estimate of the ambiguities.Assuming that these ambiguities are correct, carrier-phase measurements on both frequencies can be debiased, resulting in unambiguous pseudorange measurements, in which the largest remaining error source, when significant vertical separations are present, is due to residual tropospheric delays. In [24] it is shown that this effect can be modeled as a residual ZTD, and the latter can be calculated together with the baseline between the base station and the moving vehicle as the solution of a point-wise kinematic least-square algorithm processing the debiased carrier-phase observables. This approach is implemented in this work.The main problem with the selected algorithm is that the integer estimates derived by LAMBDA are not deterministic but stochastic variables, so there is in general a nonnull probability to derive wrong integer estimates. For increasing the chances to avoid using an erroneous estimate, integer validation tests can be applied to discriminate estimates that are likely to be correct from others that are probably wrong. The prominent integer validation tests in the open literature, described in great detail in [11], can be termed as vector tests, because they screen the whole vector of ambiguities. In case only a few ambiguities in the vector are erroneous, they discard the whole vector. These tests are based on the assumption of debiased Gaussian float ambiguity distribution, and they test the closeness of the float solution to its nearest integer vector. Under this assumption, they provide a theoretical framework for deciding which integer estimates are more likely to be correct and to filter out erroneous ones. The most common vector test is the so-called ratio test [11] that has been recently revised [25, 26] to achieve constant and fixed failure rate. The design failure rate value is chosen as a tradeoff between the requirement of not fixing erroneous integer estimates and the need for fixing the highest possible number of correct ones. Even though the acceptable failure rates depend on the application at hand, typical values are in the order 0.1–1%. This means that, independently from the actual LAMBDA estimates, the test will necessarily discard an amount of integer estimates. In other words, in all the conditions in which the failure rate that can be achieved by the combination of the selected GPS models and the LAMBDA method is lower than that of the validation test, the latter becomes of limited practical usefulness and can be avoided with no performance losses.GPS models can be classified as strong or weak based on the failure rate they achieve. The selected dual-frequency measurement model of (1)–(3), fed with measurements collected by high performance carrier-phase GPS receivers and subcentimeter phase center stability antennas, falls in the class of the strong models when baselines up to several tens of km are considered. The experimental results presented in [25, 26] clearly demonstrate that the ratio test can be discarded for strong models even in the extreme case of epoch-by-epoch (EBE) ambiguity estimation. In EBE the integer ambiguities are calculated every epoch anew. This approach is less conservative and robust than that selected one in this work in which the ambiguities are accumulated over many time epoch, and correlations among the ambiguities and the other parameters of the state vector are exploited; hence, the accumulation is able to cut down the percentage of wrong estimates with respect to EBE, notably [25–27]. In conclusion, the accumulation achieved by recursive least-square estimation and the strong GPS model are expected to reduce the failure rates to values that make standard vector tests useless. In such a condition, the detection of the residual wrong integers can be based on checking individual measurement residuals after residual ZTD computation. When measurement residuals show significant biases, it is likely that one or more integer ambiguities are wrong. As a consequence, the computation of the residual ZTD is performed again considering only a subset of the integer estimates. The procedure is repeated until a sufficient number of integer ambiguities can be assumed valid and correct.The result is an accurate estimation of the baseline vector between the base station and the moving vector that can be used to calculate debiased high-accuracy pseudorange measurements for the moving vehicle as(6)PPRi=∥rk-(rb+B)∥,k=1⋯s,
where s is the number of validated integer ambiguities.In conclusion, it is important to point out that the generated pseudorange measurements can be considered noontime correlated. Indeed, information relevant to multiple time instants is used only to estimate the integer ambiguities. However, once fixed ambiguities are validated they can be considered as known deterministic, and hence uncorrelated, parameters. Moreover, the derived baseline vector and residual ZTD are also time uncorrelated because they are the result of kinematic point-wise least-square estimation. Actually, the sole source of time correlation is due to receiver code and carrier tracking loops that is typically negligible for present generation high-performance GPS receivers.


## 3.2. Attitude Estimation by Multiantenna CDGPS

Multiple GPS antennas must be properly mounted on the moving vehicle to perform GPS-based attitude determination. In the present system assembly, GPS measurements from one antenna of the triplet are also those used to calculate the position of the moving vehicle with respect to the base station. This antenna is called “master” in what follows, while the other two antennas are named “slave”. GPS-based attitude determination strongly relies on the principle of CDGPS. Actually, the baseline vectors between the GPS antennas embarked on the moving vehicle are accurately known in a local body reference frame (BRF), but they can be derived with very high accuracy by CDGPS in Earth-centred reference frame (ECEF) or equivalently in North-East-Down (NED) reference frame with origin at the center of mass of the moving platform. The rotation matrix from BRF to NED can be then calculated, thus, leading to a measure of the attitude angles.CDGPS processing for attitude determination is conceptually similar to that presented in the previous section to estimate the baseline vector between the moving vehicle and the base station. More in detail, DD observables can be derived by combining GPS observables collected by couples of GPS antennas embarked on the moving vehicle, and the baseline vectors between those couples can be then estimated. When addressing the attitude determination problem the distance between GPS antennas is limited to very short values, several meters at most. In such a case, the atmospheric effects, that is, ionospheric and tropospheric delays, are completely cancelled when forming DD observables. Thus, the observation model represented by (3) can be applied. carrier-phase integer ambiguities have to be calculated for each couple of GPS antennas. As detailed in the previous subsection, the integer estimation scheme must include (1) derivation of a float solution, (2) estimation of the corresponding integer solution, and (3) validation of the integer estimate. The validation phase is again based on checking the measurement residuals, but it also takes advantage from the fact that the geometry of the antenna array is rigid and known. Indeed the rigidity constraint reduces notably the number of potential valid integer ambiguities for each couple of measurements [28]. Alternatively, the rigidity constraint could be included in the integer estimation step. In this regard, a multivariate constrained version of the LAMBDA method has been proposed [29] as a reliable solution for singleepoch, real time, and attitude estimation. With specific reference to the test system presented in this paper, however, since no real time operation is foreseen, rigidity constraint is only considered in the validation phase of integer ambiguities.The resulting relative positioning accuracy between the antennas, after ambiguity validation and compensation, is dominated by carrier-phase noise and GPS satellite geometry. It is, therefore, reasonable to assume a horizontal baseline error ranging from 2.5 to 5 mm and a vertical baseline accuracy between 5 mm and 1 cm.Once the antenna baseline vectors have been computed in the NED reference frame, different methods can be used to estimate vehicle attitude. For instance, the same methods adopted for spacecraft attitude determination, such as TRIAD or QUEST, can be used [30, 31]. An interesting solution is represented by the implicit least-squares approach [32], that is, an iterative procedure where both the baseline vectors in NED and in BRF are treated as observations with proper covariance matrices, so that baseline uncertainties in BRF and in NED can be both considered. Starting from a first tentative solution, the method works by successive corrections until a convergence criterion is satisfied. In addition, this approach also computes the covariance matrix for the least-squares solution, which can be straightforwardly included in the navigation filter, and therefore, it has been selected herein for implementation.The first tentative solution can be computed by a suboptimal approach such as the direct computation method [32]. This is briefly discussed in the following because it allows deriving a basic error budget for attitude measurements and understanding the antenna configuration impact on the accuracy of attitude determination. Let us consider the antenna geometry depicted in Figure 3. The master antenna is located in the origin of BRF. In this way, we can neglect the distance of the master antenna from the vehicle center of mass. The baseline from the master to the first slave antenna is along the vehicle longitudinal axis (body axis 1), while the baseline from the master antenna to the second slave antenna is along the body axis 2. Thus, the antenna geometry is such that the two baselines are orthogonal and of equal length L. This is an optimal installation geometry as discussed in [33].Figure 3
Antennas configuration geometry.Given this geometry, the baseline components in BRF and in NED are as follows:(7)b_1=[L00]b_2=[0L0]u_1=[n1e1d1]u_2=[n2e2d2].
By definition, the attitude matrix M converts vectors from NED to BRF.Consider(8)b_1=Mu_1,b_2=Mu_2.Considering the rotation sequence 3-2-1 (headingγ-pitch β-roll α) and using (8), it is possible to derive the following expressions for the attitude angles:
(9)γ=atan(e1n1)+k1π,β=-atan(d1n12+e12),α=asin(d2Lcosβ)+k2π,
where in general k1 and k2 can equal 0 or 1 depending on the sign of NED baseline components.These equations can be used to calculate a first-order estimate of the uncertainty on the attitude angles, under the assumption of uncorrelated errors. For example, heading uncertainty can be estimated as follows:(10)σγ2≅(∂γ∂n1)2σn12+(∂γ∂e1)2σe12.Explicit calculation of the derivatives shows that the propagation of CDGPS errors into the attitude solution depends on the actual vehicle attitude. As an example, Figure4 depicts heading and pitch uncertainty as a function of baseline length L, considering γ = 45°, β = 5°, and α = 1°. Different values for the CDGPS errors are considered which are based on the previous discussion. It can be seen that availability of mm-level accuracy on the NED components makes it possible to achieve about 0.1° accuracy even with relatively short baselines. Pitch (and roll) error is larger than heading due to the effect of vertical accuracy, antenna geometry, and vehicle attitude.Heading (a) and pitch (b) uncertainty as a function of antenna baseline.
(a)
(b)Finally, it is worth noting that in order to ensure accurate GPS-based attitude determination, after system installation, a calibration procedure has to be implemented to calculate with high accuracy the BRF components of baselines among the antennas.


## 3.3. EKF

The integration among the accurate pseudoranges generated by CDGPS algorithm, the attitude measurements delivered by the multiantenna GPS configuration and the inertial measurements are performed by extended Kalman filter. The state vector of the EKF include 9 components which represent estimated error on vehicle absolute position, that is, latitudeϕ, longitude ψ, and altitude h, with respect to the WGS-84 ellipsoid, ground velocity, that is, North vn, East ve, and Down vd components, and attitude with respect to the so-called geographic or navigation frame [20].Consider(11)δx=[δpδvρ],δp=[δϕδψδh]T,δv=[δvnδveδvd]T,ρ=[εnεeεd]T.The origin of this frame is the projection of vehicle position onto the WGS84 ellipsoid, the geographicz-axis points toward the interior of the ellipsoid, and x-axis points toward north. The y-axis completes the orthogonal right-handed reference frame pointing east. The attitude vector error ρ comprises the three small-angle rotations from the true geographic frame to that calculated from the computed state vector [20].The state vector can be propagated through the following linearized model that can be derived from INS mechanization equations:(12)[δp˙δv˙ρ˙]=F[δpδvρ]+[00-Rbn00-Rbn][δfδω],
where 0 is the null matrix, Rbn is the rotation matrix from the body to the navigation reference frame, f is the vector of specific forces measured by the accelerometer, and ω is the vector of relative angular rate measurements. The expression of the matrix F is not reported herein for brevity, but it can be found in [20] to which we refer the interested reader.The observation model of the EKF accounts for the contributions of CDGPS and multiantenna GPS. Specifically, the measurement vector, representing the difference between the true measurements and the measurements estimated from the propagated state vector can be written as(13)δy∈ℝ(s+3)×1,δy=[δPPR1⋮δPPRsδε^nδε^eδε^d],
where ε^ is the attitude rotation measured by multiantenna GPS. As a consequence the linearized observation model is
(14)δy=HEKFδx,HEKF∈ℝ(s+3)×9,HEKF=[Hs0000000I3],Hs=[h1RneD⋮hsRneD],
where Ip is p-dimensional identity matrix, hk is the line-of-sight vector in the ECEF reference frame from the kth GPS satellite to the vehicle, Rne is the rotation matrix from the navigation frame to ECEF, and D is the matrix converting polar coordinates (i.e., ϕ, ψ, and h) to linear coordinates within the navigation frame.Starting from the initial condition, standard nonlinear INS mechanization equations are used to propagate vehicle position, velocity, and attitude, whereas the linear model of (12) propagates the navigation error. When CDGPS and multiantenna GPS measurement are available, EKF estimation theory is applied considering the dynamic model of (12) and the measurement model of (14) leading to an estimate of the measurement error, that is, then subtracted to vehicle position, velocity, and attitude to provide the final GPS-aided accurate estimate (see Figure 5). Finally, concerning filter tuning, process noise matrix can be calculated on the basis of IMU specifications, while measurement covariance matrix is built on the basis of the residuals estimated by GPS processing and the attitude angles uncertainties estimated in the least-squares solution. In this way, less accurate pseudorange and/or attitude measurements have a smaller impact on the accuracy of the final navigation solution. The filter solution is provided at the IMU output frequency, thus, enabling a significant increase with respect to GPS data availability.Figure 5
Block diagram of the aiding process.


## 4. Performance Assessment

The performance of the developed navigation facility is preliminarily evaluated by analyzing the behaviour of the covariance of the EKF. Covariance propagation can be conducted considering the following EKF equations:(15)Pt|t-1EKF=Φt-1Pt-1|t-1EKFΦt-1T+Qt-1EKF,KtEKF=Pt|t-1EKF(HtEKF)T[HtEKFPt|t-1EKF(HtEKF)T+RtEKF]-1,Pt|tEKF=(I-KtEKFHtEKF)Pt|t-1EKF,
where the notation Πi|i-1 is used for the variable Π at the time instant i conditioned on all the measurements up to the time instant i-1. In (15) PEKF is the covariance matrix of δx, QEKF is the process noise covariance matrix of the aiding measurements, and Φ is the discrete time state propagation matrix depending on the matrix F and the time step Δt.Consider(16)Φ∈ℝ9×9,Φ=eFΔt.The process noise covariance matrix can be evaluated on the basis of the process noise spectral density matrixQc according to Von Loan’s formula. In the present case, a diagonal process noise spectral density matrix is assumed, with entries depending on the estimated variance of the accelerometers, σa2, and gyroscopes, σg2 estimated based on inertial sensors datasheet the following:
(17)Qc∈ℝ9×9,Qc=[0000σa2I3000σg2I3].REKF is the covariance of the measurement noise and depends on the variance of the accurate pseudoranges, σPPR2, that are assumed uncorrelated, and the variance of attitude estimations. According to the considerations reported in Section 3.2, the latter can be represented as two different contributions, the tilt error determined by pitch and roll accuracy and the yaw error depending on heading accuracy.Consider(18)REKF∈ℝ(s+3)×(s+3),REKF=[σPPR2Is000σTilt2I2000σYaw2].Inertial measurement rate and GPS measurement rate are assumed to be of 100 Hz and 1 Hz, respectively. The values reported in Table2 are assumed for the parameters defining the matrix REKF. With specific reference to the selected value for GPS-based attitude estimation accuracy, they correspond to GPS antennas mounted with about 1 m baseline length and therefore compatible with the installation on a general aviation aircraft or medium size UAVs.Table 2
GPS-based measurement parameters.
Parameter
Value
σ
PPR
3 cm
σ
Tilt
0.3°
σ
Yaw
0.15°As expected, covariance propagation is not affected by the true measurements but only by their statistical characterization, and the resulting diagonal terms of the state covariance matrix represent a direct measurement of the accuracy the EKF is able to attain. Nevertheless, the practical computation of the evolution of the state covariance matrix requires reference scenarios to be set in terms of time variation of the matrixF, or equivalently position, velocity, and attitude of the moving vehicle. GPS satellite geometry must be also simulated for evaluating the observation matrix HEKF.Two different scenarios are considered in this section. The first scenario simulates a small aircraft or an unmanned aerial vehicle (UAV) flying a linear flight path at constant velocity and altitude, whereas the second scenario considers the air vehicle loitering over a specific area of interest, for example, to perform observation or surveillance of that area. Specifically, the second scenario is based on a coordinated turn with the vehicle flying a circular path in about 80 seconds, with constant altitude and roll angle, 500 m ground radius, and heading angle describing the full circle. The detailed list of the parameters relevant to the selected scenarios is reported in Table3. In both scenarios a typical GPS geometry is assumed characterized by 5 available accurate pseudoranges and yielding values of 4 for the position dilution of precision (DOP) and 3.8 for the vertical DOP.Table 3
Aiding parameters.
Scenario 1straight path
Scenario 1loitering
Velocity norm
40 m/s
Altitude
1000 m
Starting latitude
40.8°
Starting longitude
14.8°
Starting heading
45°(constant)
0°
Heading rate
0°/s
~4.5°/s
Roll
0°(constant)
17.75°(constant)
Pitch
2°(constant)
0°(constant)Figure6 shows the resulting height and attitude errors for the first scenario as derived from the diagonal terms of the propagated PEKF. It is immediate to note that the filter readily converges in a few seconds. The role of GPS aiding can be also clearly recognized from the saw-tooth shape of the diagrams, meaning that the accumulating error due to the inertial measurements is drastically reduced when GPS-based correction is available. In addition, the added value of combining GPS-based and inertial measurements is also verified because both the height and the attitude errors are lower than those achievable by individual measurements. The same considerations hold for the remaining components of the state vector that are not reported for brevity.Figure 6
Height and attitude errors for scenario 1 (straight path).As far as the second scenario is concerned, Figure7 shows the attitude errors. Convergence is still achieved in few seconds, and the difference between roll, pitch, and heading accuracy is essentially originated by the actual value of the attitude angles needed to perform the coordinated turn. The main effect of the loitering maneuver is that the errors slightly oscillate depending on the variation of the specific forces while the air vehicle is turning. Due to the coupling between the different terms of the EKF state vector, this effect is present for all the variables even if, as expected, it is more pronounced for north and east velocity components as shown in Figure 8.Figure 7
Attitude errors scenario 1 for scenario 2 (loitering).Figure 8
Oscillation of the estimation errors in scenario 2 (loitering) generated by the variation of the specific forces.In conclusion, the conducted covariance analysis confirms all the expected results; namely, the proposed solution for GPS-aided inertial navigation, implemented through the described system architecture and the selected algorithms, is theoretically able to match the desired position, velocity, and attitude requirements.


## 5. Experimental Results

First, experimental data were collected during a van test conducted to validate the system concept and also to evaluate the performance of AX1-GNS unit manufactured by Axitude [5, 6]. AX1-GNS is a compact combined navigation unit based on three separated components: attitude sensors relying on MEMS technology, remote magnetic sensing unit, and single frequency GPS receiver with satellite-based augmentation system (SBAS) capabilities.In this test, a partial version of the developed navigation facility was integrated which comprises only the IMU and the BD960 receiver with a single onboard antenna. Having a single antenna configuration on the vehicle hinders attitude calculation based on multiantenna GPS from being performed, with negative effects on heading accuracy, especially, in absence of significant accelerations. However, the setup was used for a first assessment of the system potential and reliability with particular concern to position and velocity measurements.The test was conducted on April 28th 2011 in the area of the Axitude factory (Giugliano in Campania, to the North of Naples, Italy), at about 40.9° latitude and 14.1° longitude. The acquisitions started in the morning at the time 10 UTC. The maximum distance between the vehicle and the base station during the test was about 11 km, while the maximum velocity was of the order of 36 m/s. Raw inertial and GPS data from the test system were acquired together with the standalone position fix of the BD960 receiver and the outputs from the AX1-GNS system. GPS-carrier-phase and position data were acquired at 5 Hz frequency in this test. A stream of about 45 min data has been processed by the designed algorithms.


### 5.1. Navigation System Performance

GPS satellites availability was satisfactory during raw data acquisition. The number of tracked satellites ranged from 6 to 10 for both the base station and the moving receiver. The relatively short distances between the base station and the rover, with not meaningful vertical separation, allows tropospheric and ionospheric delays to be neglected. For this reason a simplified version of the CDGPS algorithm presented in Figure2 was implemented in which residual ZTD is not computed. As detailed in Section 3.1, crucial steps in CDGPS schemes are integer ambiguity estimation and validation. Figure 9 shows an example of single-epoch float estimates for both wide-lane and L1 ambiguities. The zero level corresponds to the validated ambiguity. As far as wide-lane ambiguity is concerned the mean of the residuals is very close to one integer value, with a difference lower than a few hundredth of WL cycle; moreover, the float estimate are well centered about the mean value with more than 80% of the residuals within half a wide lane cycle. L1 residuals are even closer to the integer solution thanks to the very low noise of L1 carrier-phase observables. As a matter of fact all the samples falls within half a L1 cycle. In these cases the accumulation process performs well and correct integer estimates can be derived.Example of (a) wide lane and (b)L1 float residuals. The zero-level in red refers to the validated ambiguity (graphs in cycles, i.e., with L1 cycles being more than 4 times shorter than wide lane one).
(a)
(b)As noted above, the implemented CDGPS scheme is based on the assumption that both ionospheric and tropospheric DD delays are negligible in the considered test. The validity of this hypothesis can be verified after integer ambiguity validation. Actually, once wide lane andL1 ambiguities are resolved, also L2 ambiguities can be derived. With dual-frequency debiased carrier-phase observables, ionospheric-free combinations can be formed [20] and the relevant residuals computed. Figure 10 compares L1 and ionospheric free residuals derived for the analyzed van test data. It is immediate to note that L1 residuals are smaller than ionospheric free ones. Both the residuals are dominated by measurements noise that is about a factor 3 higher for ionospheric free combination. This result actually suggests that DD ionospheric delays are not significant in the considered test. In addition, since L1 residuals are of the order of the expected DD L1 carrier-phase noise, also the possibility to neglect DD tropospheric delays can be considered verified.Comparison betweenL1 (a) and ionospheric free (b) measurement residuals.
(a)
(b)For the sake of assessing, the quality of the results produced by the developed Kalman filter, the estimated vehicle positions can be compared with the internal stand-alone solution computed by the BD960 receiver. During all the test, a good agreement (difference of the order of a few meters at most) was observed between the two solutions. It is clear however that BD960 stand-alone solution cannot achieve the same accuracy of the EKF, because it cannot exploit either differential corrections or INS measurements. Figure11 reports the two solutions (North components) computed in the last part of the acquisition when the rover was standstill. Very stable solution, with variations limited to cm level, is attained by the developed EKF. The same type of comparison can be carried out for the velocity estimates. Again, a good agreement between receiver and filter output was observed in the entire test duration, but the CDGPS/INS solution exhibits a lower noise, some cm/s, as shown in Figure 12 which depicts the North velocity components, again during the final, static part of the van test.Figure 11
Comparison between North-components estimated by the developed EKF and the BD960 stand-alone fixing (standstill vehicle).Figure 12
North velocity component as estimated by the developed EKF and by the stand-alone BD960 receiver (standstill vehicle).


### 5.2. Utilization of the Navigation System as Testing Facility

EKF results have been taken as a reference to test the performance of the AX1-GNS system. Thus, in what follows the term “error” will be used to indicate the difference between the AX1-GNS value and the corresponding EKF one. Table4 demonstrates the level of accuracy of the AX1-GNS system measurements. While the maximum value of the positioning error can reach up to 10 m; the values of the error mean and standard deviation (std) are always smaller than 2 m. Concerning the velocity, the error standard deviation is smaller than 10 cm/s. Finally, the performance of the AX1-GNS system while measuring the attitude angles shows about 0.3° error standard deviation for both pitch and roll. Heading results, instead, are not meaningful for the implemented partial version of the test system and therefore are not reported.Table 4
AX1-GNS error analysis for position, velocity, and attitude.
Max error
Mean error
Error Std
x
ECEF
10 m
-
7
·
10
-
2 m
1 m
y
ECEF
−5 m
3
·
10
-
2 m
1 m
z
ECEF
−6 m
3
·
10
-
1 m
2 m
v
n
−1 m/s
<10−2 m/s
8
·
10
-
2 m/s
v
e
2 m/s
<10−2 m/s
10−1 m/s
v
d
6
·
10
-
1 m/s
<10−2 m/s
4
·
10
-
2 m/s
β
4°
<0.2°
0.4°
α
2°
<0.2°
0.3°


## 5.1. Navigation System Performance

GPS satellites availability was satisfactory during raw data acquisition. The number of tracked satellites ranged from 6 to 10 for both the base station and the moving receiver. The relatively short distances between the base station and the rover, with not meaningful vertical separation, allows tropospheric and ionospheric delays to be neglected. For this reason a simplified version of the CDGPS algorithm presented in Figure2 was implemented in which residual ZTD is not computed. As detailed in Section 3.1, crucial steps in CDGPS schemes are integer ambiguity estimation and validation. Figure 9 shows an example of single-epoch float estimates for both wide-lane and L1 ambiguities. The zero level corresponds to the validated ambiguity. As far as wide-lane ambiguity is concerned the mean of the residuals is very close to one integer value, with a difference lower than a few hundredth of WL cycle; moreover, the float estimate are well centered about the mean value with more than 80% of the residuals within half a wide lane cycle. L1 residuals are even closer to the integer solution thanks to the very low noise of L1 carrier-phase observables. As a matter of fact all the samples falls within half a L1 cycle. In these cases the accumulation process performs well and correct integer estimates can be derived.Example of (a) wide lane and (b)L1 float residuals. The zero-level in red refers to the validated ambiguity (graphs in cycles, i.e., with L1 cycles being more than 4 times shorter than wide lane one).
(a)
(b)As noted above, the implemented CDGPS scheme is based on the assumption that both ionospheric and tropospheric DD delays are negligible in the considered test. The validity of this hypothesis can be verified after integer ambiguity validation. Actually, once wide lane andL1 ambiguities are resolved, also L2 ambiguities can be derived. With dual-frequency debiased carrier-phase observables, ionospheric-free combinations can be formed [20] and the relevant residuals computed. Figure 10 compares L1 and ionospheric free residuals derived for the analyzed van test data. It is immediate to note that L1 residuals are smaller than ionospheric free ones. Both the residuals are dominated by measurements noise that is about a factor 3 higher for ionospheric free combination. This result actually suggests that DD ionospheric delays are not significant in the considered test. In addition, since L1 residuals are of the order of the expected DD L1 carrier-phase noise, also the possibility to neglect DD tropospheric delays can be considered verified.Comparison betweenL1 (a) and ionospheric free (b) measurement residuals.
(a)
(b)For the sake of assessing, the quality of the results produced by the developed Kalman filter, the estimated vehicle positions can be compared with the internal stand-alone solution computed by the BD960 receiver. During all the test, a good agreement (difference of the order of a few meters at most) was observed between the two solutions. It is clear however that BD960 stand-alone solution cannot achieve the same accuracy of the EKF, because it cannot exploit either differential corrections or INS measurements. Figure11 reports the two solutions (North components) computed in the last part of the acquisition when the rover was standstill. Very stable solution, with variations limited to cm level, is attained by the developed EKF. The same type of comparison can be carried out for the velocity estimates. Again, a good agreement between receiver and filter output was observed in the entire test duration, but the CDGPS/INS solution exhibits a lower noise, some cm/s, as shown in Figure 12 which depicts the North velocity components, again during the final, static part of the van test.Figure 11
Comparison between North-components estimated by the developed EKF and the BD960 stand-alone fixing (standstill vehicle).Figure 12
North velocity component as estimated by the developed EKF and by the stand-alone BD960 receiver (standstill vehicle).


## 5.2. Utilization of the Navigation System as Testing Facility

EKF results have been taken as a reference to test the performance of the AX1-GNS system. Thus, in what follows the term “error” will be used to indicate the difference between the AX1-GNS value and the corresponding EKF one. Table4 demonstrates the level of accuracy of the AX1-GNS system measurements. While the maximum value of the positioning error can reach up to 10 m; the values of the error mean and standard deviation (std) are always smaller than 2 m. Concerning the velocity, the error standard deviation is smaller than 10 cm/s. Finally, the performance of the AX1-GNS system while measuring the attitude angles shows about 0.3° error standard deviation for both pitch and roll. Heading results, instead, are not meaningful for the implemented partial version of the test system and therefore are not reported.Table 4
AX1-GNS error analysis for position, velocity, and attitude.
Max error
Mean error
Error Std
x
ECEF
10 m
-
7
·
10
-
2 m
1 m
y
ECEF
−5 m
3
·
10
-
2 m
1 m
z
ECEF
−6 m
3
·
10
-
1 m
2 m
v
n
−1 m/s
<10−2 m/s
8
·
10
-
2 m/s
v
e
2 m/s
<10−2 m/s
10−1 m/s
v
d
6
·
10
-
1 m/s
<10−2 m/s
4
·
10
-
2 m/s
β
4°
<0.2°
0.4°
α
2°
<0.2°
0.3°


## 6. Conclusion

A navigation facility generating offline high-accuracy attitude and trajectory measurements for aerial applications was presented. The system achieves high accuracy in position and attitude estimation exploiting carrier-phase differential GPS techniques, multiantenna GPS configurations onboard the vehicle, and integrating GPS-based estimations with inertial data by Kalman filtering. The selected CDGPS algorithm includes a specific procedure to deal with residual differential tropospheric delays resulting from the significant vertical separation that can be established between the base station and the air vehicle depending on the operative scenario. The theoretical performance of the developed system has been evaluated by covariance analysis showing that the navigation facility is potentially able to achieve cm-level navigation accuracy and less than 0.15° errors for the attitude angles, heading angle included.First, experimental results were collected with a partial version of the facility embarked on a van, including a single GPS antenna installed onboard the vehicle. In spite of the limitations of the partial setup, results are promising. Offline processing of carrier-phase differential observables demonstrated high accuracy and very good availability. The latter is a key factor for system effectiveness in flight tests. At the same time, integration of inertial sensors was found to be important for increasing the measurement rate and reducing the noise in the navigation solution. During the static phases of the tests, with standstill vehicle, it has been observed that the estimated position and attitude of the vehicle was very stable featuring less than 2 cm position variations and less than 0.05° attitude variations.Integration of the complete multiantenna test system, and installation onboard an experimental aircraft are currently in progress. An augmented version of the navigation filter, including estimates of accelerometers and gyroscopes biases in the state vector, is also under evaluation to further improve filter performance.



---

*Source: 397686-2013-09-05.xml*