Global Navigation Satellite Systems (GNSS) constitute the main information supplier for outdoor positioning and timing. Unfortunately, unintended and malicious radio-frequency interference constitutes a major concern for GNSS. Within the maritime domain, radio attacks deprive the skippers from navigation and lead to the failure of diferent interfaces on a vessel's bridge. Ranging-Mode (R-Mode) is a positioning system which uses signals-of-opportunity and allows to backup GNSS. This paper describes the overall R-Mode architecture and details the types of ranging and velocity observations derived from medium and very high frequency radio signals. This paper contributes with the derivation of recursive estimators, in the form of Cubature Kalman Filters, and the discussion of the particulars required to accommodate R-Mode observations. The theoretical analysis and performance comparison of the proposed filters is addressed via Monte Carlo simulation following a realistic maritime scenario.
2008 [
R-Mode uses two types of existing signals on diferent frequency bands for deriving ranging
observations: 1) diferential GNSS corrections transmitted by maritime radio beacons over
medium frequency (MF); 2) base stations receiving and broadcasting AIS messages over very
high frequency (VHF). Fig.1 illustrates the potential availability of R-Mode compatible stations
over the Baltic Sea [
The rest of the paper is organized as follows: Section 2 provides a brief overview on the R-Mode architecture, its signal models and related ranging and velocity observations. Section 3 discusses on recursive estimators using R-Mode observations and analyses the specifics for a Cubature Kalman Filter solution. Then, Section 4 studies the positioning performance of the previously derived estimators for a simulated R-Mode scenario. Finally, Section 5 presents an outlook of the work and discusses the future lines of work.
The R(anging)-Mode concept emerged in 2008 in order to ofer an alternative PNT terrestrial
system to the maritime users [
The potential of R-Mode has been presented in three feasibility studies published in 2014
[
The frequency of the signals in use is fundamental in order to understand how the
radiowave propagates. In the VHF band the signals propagate as line-of-sight (LOS), hence the
transmitter-receiver distance can be described with the classical formulation of the Euclidean
distance, analogous to the GNSS case. VHF measurements are obtained by using correlation
techniques. In particular, range and velocity observation are gathered from the stations by
evaluating the time of arrival (TOA) and the Doppler frequency of the received signal
respectively [
ρ = τ c 0
vr = −
∆ f
f0
c0
with ∆ f Doppler frequency and f0 carrier reference frequency. In (2) a positive Doppler
means that the distance between the receiver and the transmitter is decreasing, therefore the
resulting radial velocity is negative. In general, velocity measurements present noise of one
order of magnitude lesser than the ranges and, therefore, they are of particular relevance to
achieve good positioning performance. The interested readers can refer to [
Diferently from VHF signals, the MF ones propagate as ground waves which implies that the
signals follow the curvature of the Earth. The main advantage of MF signals is that they are
not limited to LOS and can, therefore, be received at greater distances. However, such curved
propagation implies that the Euclidean distance does not correctly model the range between
receiver and transmitter. Instead, MF-derived ranges are typically modelled as geodesics on
the Earth ellipsoid and described with the Vincenty’s distance model [
gj,k = Nj,kλ j,k + ϕ j,kλ j,k
where j is the station index whereas k = 1, 2 is the CW index. N represents an integer number
of phase cycle, often referred to as ambiguity in the GNSS literature, ϕ is the fractional phase
cycle and λ represents the wavelength of the CW signal. Further information about the signal
characteristics and phase estimation technique for MF R-Mode can be found in [
After this brief description of the signals in use, one clearly identifies that one of the main
challenges for the position solver is to consider and combine the diferent observations with
their pertinent models, especially given the MF ground-wave signal propagation. In [
Let us consider a general discrete state-space model (SSM) whose time evolution is described by a process model f (· ) and the relationship between the state x and a set of measurements y by an observation model h(· ), as xk = f (xk− 1, wk− 1) , yk = h (xk) + nk, with the vectors of process and observation noise assumed to follow Gaussian distributions of known parametrization (i.e., typically zero-mean and known covariance matrix), such that wk− 1 ∼ N (0, Qk− 1), nk ∼ N (0, Σ k).
At a particular time instance k, the state estimate is given by
⎡ pk ⎤ xk = ⎢⎢⎣cvδt kk⎥⎦⎥ , cδ ̇ tk with pk, vk ∈ R3, cδt k, cδ ̇ tk ∈ R1, xk ∈ Rp, (6) where pk, vk denote the three-dimensional vectors of position and velocity, expressed in a global reference frame (such as the Earth-centered Earth-fixed (ECEF)), while cδt and cδ ̇ t express the clock ofset and clock ofset rate, respectively. Then, p is the state dimensionality. Note that, in analogy to the GNSS world, the clock ofset delays at the receiver end shall be estimated while clock deficiencies at the transmitter end are expected to be covered by the data transmission or a prospective R-Mode augmentation system.
As introduced during Section 2, the reception of VHF and MF R-Mode signals leads to two types of ranging observations and a radial velocity measurement, with details on the corresponding observation models provided next.
Propagating in a straight line, VHF signals lead to ranges in an Euclidean space and the following observation model
ρ iVHF = ⃦⃦ pi − p⃦⃦ + cδt + nρi , for i = 1, . . . , nVHF ρ ̇ iVHF = − ︁( vi − v︁) ⊤
i p − p ∥pi − p∥ + cδ ̇ t + nρi̇ , for i = 1, . . . , nVHF, (8) where the time index has been neglected for simplicity, the superscript i refers to the each of the nVHF stations and ∥ · ∥ is an Euclidean norm. Note that the small time deviations among stations due to the TMDA modulation is negligible and, therefore, disregarded in this work.
The application of Doppler measurements for the R-Mode VHF implementation leads to
significant accuracy improvements [
As previously discussed, MF measurements propagate following the Earth’s curvature, leading to a complicated geolocation-dependent observation model, such that
ρ jMF = ⃦⃦ pj − p⃦⃦ V + cδt + nρj , for j = 1, . . . , nMF, where ∥ · ∥ V is the Vincenty’s distance, or the geodesic between two arbitrary points on the surface of the Earth. Presenting high fidelity for the Earth ellipsoid, Vincenty’s inverse problem leads to an iterative procedure based on the two latitude and longitude coordinates. Thus, the use of the celebrated Extended Kalman Filter (EKF) is jeopardized on the derivation of complex Jacobians and sensitive to the high nonlinearity of ellipsoidal distances, especially prone to errors for long distances.
To exploit the series of observations collected over time, one generally makes use of the
Recursive Bayesian Estimation (RBE) framework. Among the families of RBE algorithms, the
Kalman Filter (KF) and its nonlinear extensions –e.g., the EKF, UKF or Sigma-Point Gaussian
Filters (SPGF)– have become the norm for navigation and tracking applications [
given by Considering a constant velocity time evolution of the positioning, the process model in (4) is xˆk = Fk− 1xˆk− 1, Fk− 1 = ⎢⎢⎣ 00
0 ⎡I3 ∆ tI3 0
I3 0 0 0 0 1 ∆ t 0 ⎤ 0 ⎥⎥ ,
⎦ 1 linear KF form, as
with ∆ t the time elapsed between two consecutive discrete time instances. With the prediction step being a linear function, the prediction for the covariance matrix can be realized using the
Pk|k− 1 = Fk− 1Pk− 1|k− 1Fk⊤− 1 + Qk− 1.
Then, the update equations for the state and covariance matrix result from the propagation of the cubature points.
To do so, one first propagates these points factorization of the covariance matrix and their evaluation on the observation model, as Xk|k− 1 based on the with ε the generators for the 2p points. From these, we may reconstruct the mean of our vector of observation models with where n is the total number of observations (i.e., n = 2 · nVHF + nMF to account for the Doppler and ranging measurements of VHF and the ranging observations from MF) and [a]i indicates the ith position on a generic vector a. Then, the update of the state and covariance estimates follow the well-known CKF algorithm [17]:
Xk|k− 1 = xk|k− 1 + P1/2 k|k− 1
ε Yk = h (︁ Xk|k− 1
︁) ˆ yˆk = n i=1 n 1 ∑︂[Yˆ k]i Pyy,k|k− 1 = Pxy,k|k− 1 = n n n i=1 n i=1 1 ∑︂ Yˆ kYˆ k⊤ − yˆkyˆk⊤ + Σ k 1 ∑︂ Xˆ k|k− 1Yˆ k⊤ − xˆk|k− 1yˆk⊤ Kk = Pxy,k|k− 1P−yy1,k|k− 1 xˆk|k = xˆk|k− 1 + Kk (yk − yˆk) Pk|k = Pk|k− 1 −
KkPyy,k|k− 1Kk⊤.
(10) (11) (12) (13) (14) (15) (16) (17) (18) 9 11
12 14 Longitude [Deg] 15 17 18
The algorithm proposed has been tested in a simulated scenario. Fig. 2 shows the map with the VHF (blue triangles) and MF (red circles) base stations used to generate the synthetic observations (i.e. pseudoranges and Doppler measurements). The initial receiver position is also visible in the map represented as a black square. It is important to mention that the location of the MF transmitters used for the simulation are all enabled to send R-Mode signals, whereas the VHF ones are only potential candidates for the future implementation of the system.
To assess the performance of the filter under study, a Monte-Carlo simulation with 100 realizations was performed. The values of the process noise covariance matrix describing the model mismatch of the filter are used to generate the random walk of the reference trajectory in each realization. The trajectory is initialized randomly around the starting point indicated in Fig. 2 with an uncertainty as indicated by the first element of the first row of Tab. 1. The values of the first row of Tab. 1 indicated as initial uncertainties are also used as initial process noise covariance matrix of the filter. In the second row of Tab. 1 the standard deviation used for the process noise covariance matrix are given, for both velocity and clock drift. Last but not least, the third row of Tab. 1 contains the standard deviation used to generate the noise on the simulated measurements (i.e MF range, VHF range, VHF Doppler Velocity) and the measurement noise covariance matrix Σ .
In Fig. 3a, the time evolution for the root mean squared error (RMSE) for the position estimates is depicted. Specifically, the orange line (as indicated by the legend) illustrates the error in the ECEF coordinate system, while the blue line indicates the horizontal 2D positioning error in the East-North coordinate system. The diference between the two represents therefore Clock rate: 0.1 [s/s/s] Velocity (East-North-Up): [0.1, 0.1, 10− 4] [m/s/s] (a) RMSE of position in 3D ECEF and 2D ENU (b) Mean DOP in horizontal and vertical domain. the contribution of the error in the vertical domain (Up), on which we expect poor performance due to the low vertical dilution of precision (VDOP). Nevertheless, for maritime application the interest primarily lays on the horizontal domain, verifying that there are no instabilities in the vertical one. Indeed, for the three-dimensional RMSE in Fig. 3a, it seems apparent the presence of a transient efect of the filter, which is a typical and expected behavior. The convergence to a steady-state appears to be slow but it can be explained by the fact that the values in the process noise covariance matrix are very small as indicated in Tab. 1. Also, the long convergence time observed might be induced by the poor observability over the vertical component. In the future, a deeper analysis on the convergence will be done and techniques to speed-up this process will be investigated.
In terms of performance, the filter should provide close-to-optimal filtering capabilities, conditioned on the stochastic modeling is known. As afore-explained, the same process noise covariance matrix Q is used to create the reference and the measurement noise is assumed to be perfectly known. In order to mitigate the errors in the vertical domain, the associated process noise standard deviation is reduced to 10− 4 m in the process covariance matrix. Following the small acceleration capabilities of commercial vessels, a standard deviation of 0.1 m/s/s is simulated and thus used for filter design.
Even if the three-dimensional position estimate appears to not fully converge to steady-state, an accuracy of 62 m is achieved at the end of the simulation. Nevertheless, such accuracy is surprising given the large VDOP as visible in Fig. 3b. For better assessment, Fig. 4 depicts only the RMSE in the horizontal direction. Following the small HDOP values visible in Fig. 3b, the filter shows excellent performance. Despite the small overshooting efect visible until epoch 1800, the filter is able to fully converge and estimate the user position with the depicted RMSE of approximately 2.5 m. As stressed earlier, the positioning performance, especially in the horizontal domain, is of particular interest for the maritime application. The accuracy achieved by this approach outperforms the predicted accuracy within the testbed. However, it is important to mention that this is the result of a simulation for filter validation. Knowing all simulated uncertainties stated in Tab. 1, especially for the process noise, is not possible in a real-world environment. Still, this performance validates the implemented filter approach in accordance to the low HDOP values.
Both DOP values shown in Fig. 3b show constant performance throughout the whole simulation. The orange graph in depicts the HDOP in the associated East-North coordinate system. The extremely low value of 0.65 shows good performance, which is expected due to the geometrically favourable conditions of the initial position (black square in Fig. 2). As the geometrical conditions experience no significant changes throughout the whole simulation, only small variations can be observed in the HDOP. The vertical DOP (blue graph of Fig. 3b) amounts to 17.25 and demonstrates considerably bad geometrical conditions in the vertical domain. Also this value shows only small changes, which can be explained by the mostly steady geometrical conditions. Using exclusively transmitters placed on the earth’s surface, a poor VDOP is expected. The approach used to limit the positioning uncertainties in the vertical direction addresses this issue. We consider this method suitable due to the harsh geometrical conditions of the vertical positioning but also to incorporate the fact that our receiver will stay on the earth’s surface, hence the approach is application specific and cannot be adopted for other scenarios where the receiver can experience large vertical variation (e.g.in aviation).
Furthermore, the validity of the approach in discussion can be stressed by the low values (a) 3D RMSE of speed (b) RMSE of clock ofset (of 0.3 m/s) of the three-dimensional RMSE on the speed estimation depicted in Fig. 5a. The iflter is suficient dynamic to achieve convergence within approximately 2000 s, which is fast, compared to the convergence time of the RMSE on the position depicted in Fig. 3a.
The convergence time of the RMSE of the receiver clock ofset depicted in Fig. 5b is even smaller. Even though the low ofset is centered around a relatively low value of 1.5 m, its behaviour is comparably noisy. This efect could be mitigated by an increased number of Monte Carlo realizations, which would raise the computation time and therefore was not considered in the current implementation.
In this paper, we described the overall architecture of R-Mode, a GNSS backup system for the maritime domain. The positioning accuracy of the system is predicted to be less than 10 m, 95 % of the time. From a positioning perspective, the main challenges are related to the integration of the two subsystem (MF, VHF) and the low observability on the vertical direction caused by the transmitters geometrical setup.
Due to the high non linearity of the MF range measurement and the complexity of its Jacobian term derivation, we excluded the adoption of the EKF. Therefore, a CKF was implemented as PNT algorithm. Its theory was briefly presented along with the definition of the diferent observation models.
To assess the performance of the filter, a simulation environment was created in order to perform Monte Carlo runs. The realizations were obtained by using real VHF and MF transmitter locations and we showed that the filter is able to accurately estimate all the relevant states by using an adapted process noise covariance matrix. In particular, the usage of very low variance in the vertical speed acts as a bound for the variation in the same direction and reduces the impact of unfavorable VDOP. The results of the horizontal positioning performance are extremely promising: After convergence, the RMSE in the horizontal domain amounts to 2.5 m, which is surprising, even in a simulated environment. Nevertheless, the long convergence time of the filter in the three-dimensional domain shows the poor observability in the vertical direction. However, for the application presented, the practical purpose of a vertical positioning estimate is limited.
In the future, several research directions are foreseen to be considered. First of all, perform a test with the proposed algorithm on the field by using real measurements. Secondly, the comparison of alternative positioning algorithms and finally the integration of additional sensors which can reduce the efects associated to the poor geometry in the vertical domain. cubature rules applied to sensor data fusion for positioning, in: 2010 IEEE international conference on communications, IEEE, 2010, pp. 1–5. [15] D. Medina, A. Heßelbarth, R. Bu¨scher, R. Ziebold, J. Garıac´ , On the Kalman filtering formulation for RTK joint positioning and attitude quaternion determination, in: 2018 IEEE/ION Position, Location and Navigation Symposium (PLANS), IEEE, 2018, pp. 597–604. [16] A. Heßelbarth, D. Medina, R. Ziebold, M. Sandler, M. Hoppe, M. Uhlemann, Enabling Assistance Functions for the Safe Navigation of Inland Waterways, IEEE Intelligent Transportation Systems Magazine 12 (2020) 123–135. [17] I. A. S. Haykin, Cubature Kalman Filters, IEEE Transactions on automatic control (2009).