Mars Networks-Based Navigation: Observability and Optimization

Abstract

In order to achieve more scientific returns for Mars, future Mars landers will be required to land at certain landing point with special scientific interest. Therefore, autonomous navigation is indispensable during the Mars approach, entry, and landing phase. However, the number of beacons or the Mars orbiters which can provide the navigation service is so limited and the line-of-sight visibility cannot be guaranteed during the landing period. So the navigation scheme especially the beacon configuration has to be optimized in order to efficiently use the limited navigation information. This chapter aims to analyze the feasibility and optimize the performance of the Mars Networks-based navigation scheme for the Mars pinpoint landing. The observability of navigation system is used as an index describing the navigation capability. Focusing on the relationship between the configuration of radio beacons and observability, the Fisher information matrix is introduced to analytically derive the degree of observability, which gives valuable conclusions for navigation system design. In order to improve the navigation performance, the navigation scheme is optimized by beacon configuration optimization, which gives the best locations of beacons (or the best orbit of navigation orbiters). This is the main approach to improve the navigation capability.

1. Introduction

As the most similar planet to the Earth in the Solar system, Mars is considered as an ideal target for planetary exploration [1, 2]. Since the 1960s, humans have investigated the Mars exploration missions in the near distance. With the development of aerospace science and technology, the manner of Mars exploration has shifted from flyby/orbiting to landing and roving explorations. Considering scientific returns and exploration capabilities, Mars landing exploration is also essential and is one of the most popular tasks of human deep space exploration in the near future. The representative Mars landing missions including NASA’s Viking 1 and 2, Mars Pathfinder (MPF), Mars Exploration Rovers (MER, including the Spirit and Opportunity rovers), Phoenix, Mars Science Laboratory (MSL, including the Curiosity rover), and ESA’s Mars Express/Beagle 2 mission. All of these greatly inspire the development of advanced guidance, navigation, and control (GNC) technologies.

During the past 50 years of Mars exploration, 46 Mars exploration spacecraft have been launched. The overall success rate is only 41.3% though. Furthermore, among the 20 Mars landing attempts, only 7 robotic rovers were successful. The success rate for Mars landing missions is only 35%. Among the failed landing missions, most failures occur during the landing phase. The pinpoint landing has to be based on the precise autonomous navigation technology.

In the entry phase of a Mars landing, the lander is covered by a heat shield which blocks the optical sensor measurement, causing that all landers relied on the Inertial Measurement Unit (IMU) recursion. The initial errors of the lander cannot be corrected by IMU data. Even worse, the recursion errors using IMU are accumulated due to the sensor bias and noise. To overcome the incapability of IMU, the Mars Network-based Mars entry navigation is developed based on high frequency radio communication between the lander with ground or orbiting radio beacons [3, 4, 5]. Involving the radio measurement date into a navigation filter, the position and velocity of the lander can be optimally estimated.

The Mars Network-based Mars entry navigation is faced with two challenging. One is that the geometric configuration of the radio beacons affects the navigation performance. The other is that the available beacons at present are very limited. Considering these two factors, effort should be devoted to optimizing the configuration of radio beacons to maximize the function of the limited beacons. In [7], the navigation accuracy from the Extend Kalman Filter (EKF) by processing the radio measurements is analyzed, and the optimal configuration of ground beacons is selected among potential beacon position. Yu focused on the navigation observability and take it as a performance index to optimize the configuration of radio beacons [8]. The research on ground beacons, to some extent, inspired the future Mars landing navigation. However, the practice application of ground beacon-based navigation is hardly applied in practice. The first concern is that no ground beacon is available. Even if several beacons are distributed on Mars surface, it’s still a tough job to place them exactly at the optimal locations. Moreover, the accurate positions of the beacons are hardly obtained accurately. Considering the immovability of ground beacons, the potential location areas are constrained by the line-of-sight visibility, resulting in an unsatisfactory beacon configuration during the entry phase.

As a substitution of ground radio beacons, the Mars orbiters which can also serve as beacons for Mars Network-Based Navigation are of more practice value. Currently, the operational orbiter around Mars includes 2001 Mars Odyssey and 2005 Mars Reconnaissance Orbiter. With another forthcoming spacecraft Mars Atmosphere and Volatile Evolution (MAVEN) [9], the capability of Mars network can be further increased. Focusing on how to fulfill the function of a Mars network, Ely firstly established the basic principle to design a constellation for navigation [10]. Then, taking the Mean of the Position Accuracy Response Time (MPART) as the performance index, the constellation configuration was optimized [11]. In [12], the number of orbiters and the coverage was considered to design the Martian navigation constellations envisaged in the ESA’s Martian Constellation for Precise Object Location program. The optimization method of the above researches is inherited from the Global Positioning System (GPS). The global navigation performance was emphasized. For the limited amount of Mars orbiters, global coverage is difficult to realized, and local navigation performance should be investigated thoroughly for specific missions. Moreover, the effect of geometric configuration of the Mars network on the navigation performance should be revealed clearly. Inspired by these requirements, Yu et al. optimized the orbits of Mars orbiters in the observability point of view, and tried to explain the relationship between the configuration of beacons and orbiters and the navigation capability [13, 14].

To optimize the configuration of the radio beacons, a performance index should be firstly setup. The observability of the navigation system is selected as the performance index since it reflects the navigation capability directly. A lot of work has investigated the observability of linear and nonlinear dynamic systems [6, 15, 16, 17]. However, the analytic relationship between geometric configuration and observability has never been revealed. According to Cramér-Rao inequality [18], the inverse of the Fisher Information Matrix (FIM) estimates the lower bound of the estimation error. Therefore, FIM can be used to quantify the observability of the navigation system [19, 20, 21]. In this circumstance, some valuable analytic conclusions about the navigation design can thus be obtained.

Based on the requirement of the navigation optimization for Mars pinpoint landing, this chapter discusses the design and optimization of the Mars Networks-based navigation during Mars entry phase. Firstly, the Mars Networks-based navigation scheme is introduced, and the dynamic model and the observation model are given. Based on the navigation system, the observability of the Mars entry navigation analysis, and the analysis methods based on the quadratic approximation and Fisher information matrix are proposed. The relationship between the observability and the beacon configuration is derived, and the theoretically optimal configuration is given. Considering the constraints of Mars entry scenario, the ground beacons and the orbit of Mars orbiters are optimized based on observability based on an entry trajectory. The simulations also indicate the improved navigation performance.

2. Mars networks-based navigation scheme

2.1. Dynamic model of Mars entry phase

In the dynamical model with respect to a stationary atmosphere of a rotating planet, the 6 dimensional states x of the entry vehicle include r (radius from the center of Mars to the vehicle’s center of mass), θ (longitude), $\varphi$ (latitude), V (relative velocity), γ (flight path angle), and Ψ (heading angle, with Ψ = 0 as due east). The motion of the entry vehicle is governed by the following state equations:

$$\begin{array}{l}\dot{r}=V\text{sin}\gamma \\ \dot{\theta }=V\text{cos}\gamma \text{cos}\Psi /\left(r\text{cos}\varphi \right)\\ \dot{\varphi }=V\text{cos}\gamma \text{sin}\Psi /r\\ \dot{V}=-d-g\text{sin}\gamma \\ \dot{\gamma }=\left[l\text{cos}\sigma -\left(g-V^2/r\right)\text{cos}\gamma \right]/V+2\omega \left(\text{tan}\gamma \text{sin}\Psi \text{cos}\varphi -\text{sin}\varphi \right)\\ \dot{\Psi }=-\left(l\text{sin}\sigma +V^2{\text{cos}}^2\gamma \text{cos}\Psi \text{tan}\varphi /r\right)/\left(V\text{cos}\gamma \right)+2\omega \text{cos}\Psi \text{cos}\varphi \end{array}$$E1

In the equation, $\sigma$ is the banking angle, which is fixed at 0 in the following analysis. $\omega$ refers to the rotation rate of Mars. For simplicity, the second order terms of $\omega$ are neglected, which is feasible because the value of $\omega$ is quite small. Then the gravity acceleration g, lift and drag accelerations l and d are given by

$$g=\mu /r^2$$E2

$$l=0.5\rho V^2{C}_{l}S/m$$E3

$$d=0.5\rho V^2{C}_{d}S/m$$E4

where $\mu$ is the Martian gravitational constant. S and m denotes the reference area and mass of the entry vehicle, and $C_l$ and $C_d$ are the lift and drag coefficients respectively. Furthermore, the Mars atmospheric density $\rho$ is approximated by the conventional exponential model

$$\rho ={\rho }_0\text{exp}\left[\left(r_0-r\right)/h_s\right]$$E5

where ${\rho }_0=2\times {10}^{-4}$ kg/m3 is the reference density, $r_0=3437.2$ km is the reference radial position, and $h_s=7500$ m refers to the atmospheric scale height. The dynamical model of the entry vehicle is abbreviated as $\dot{\mathit{x}}=\mathit{f}\left(\mathit{x}\right)$.

2.2. Observation model

The radio ranging and velocity data between the lander and the radio beacon can be measured through radio communication, given by

$$\begin{array}{l}{y}_{R_i}=R_i+{\epsilon }_R^i\\ =\sqrt{{\left(x_B^i-x\right)}^2+{\left(y_B^i-y\right)}^2+{\left(z_B^i-z\right)}^2}+{\epsilon }_R^i\\ \begin{array}{lll}x=r\text{cos}\varphi \text{cos}\theta ,& y=r\text{cos}\varphi \text{sin}\theta ,& z=r\text{sin}\varphi \end{array}\end{array}$$E6

where $R_i$ is the real range between the lander and the i^th beacon, $x_B^i$, $y_B^i$, and$z_B^i$represent respectively the triaxial position components of the beacon, and ${\epsilon }_R^i$ is the radio ranging measurement noise.

The relative velocity model is given by

$$y_{V_i}=V_i+{\epsilon }_V^i={\mathit{dR}}_i/\mathit{dt}+{\epsilon }_V^i$$E7

where $V_i$ is the real line-of-sight relative velocity between the lander and the i^th radio beacon, and ${\epsilon }_V^i$ is the velocity measurement noise.

With different radio beacons come different navigation scenarios. Without losing the generality, the observation model can be summarized as $\mathit{y}=\mathit{h}\left(\mathit{x}\right)$. Obviously, both radio measurements in Eqs. (6) and (7) are nonlinear. Moreover, the navigation performance is closely related to the geometric configuration of radio beacons. Therefore, the beacon configuration needs to be optimized based on the observability analysis.

3. Observability of the navigation system

3.1. Observability analysis based on the quadratic approximation

Consider the following nonlinear system:

$$\mathit{\Sigma }:\left\{\begin{array}{l}\dot{\mathit{x}}=\mathit{f}\left(\mathit{x}\right)\\ \mathit{y}=\mathit{h}\left(\mathit{x}\right)\end{array}\right.$$E8

where $\mathit{x}\in {\mathbb{R}}^n$ is the n-dimensional state vector and $\mathit{y}\in {\mathbb{R}}^m$ is the m-dimensional observation vector. Define $\mathit{h}:{\mathbb{R}}^n\to {\mathbb{R}}^m$ as the nonlinear measurement operator.

The Lie algebra is an efficient tool for observability analysis. For the k^th order Lie derivative of the j^th measurement function, which can be expressed as $L_{\mathit{f}}^k{h}_j$, the k + 1th order Lie derivative $L_{\mathit{f}}^{k+1}{h}_j$ with respect to the state equation $\mathit{f}$ can be computed as:

$$\begin{array}{l}{L}_{\mathit{f}}^{k+1}{h}_j=\sum _{i=1}^n\frac{\partial L_{\mathit{f}}^k{h}_j}{\partial x_i}{f}_i=\nabla L_{\mathit{f}}^k{h}_j\mathit{f},\\ \begin{array}{ll}& \begin{array}{ll}k=0,1,\cdots & j=1,2\cdots ,m\end{array}\end{array}\end{array}$$E9

The differential of $L_{\mathit{f}}^k{h}_j$ is defined as

$$\nabla L_{\mathit{f}}^k{h}_j=\left[\begin{array}{lll}\frac{\partial L_{\mathit{f}}^k{h}_j}{\partial x_1},& \cdots ,& \frac{\partial L_{\mathit{f}}^k{h}_j}{\partial x_n}\end{array}\right]$$E10

Regarding the zero-order Lie derivative of the jth measurement function $h_j$ as $h_j$ itself, the matrix$\nabla L_{\mathit{f}}^k\mathit{h}$ is given as

$$\nabla L_{\mathit{f}}^k\mathit{h}={\left[\begin{array}{lll}{\left(\nabla L_{\mathit{f}}^k{h}_1\right)}^T,& \cdots ,& {\left(\nabla L_{\mathit{f}}^k{h}_m\right)}^T\end{array}\right]}^T$$E11

It is proven that the dynamical system $\mathit{\Sigma }$ at state ${\mathit{x}}_0$ is locally observable if the observability matrix ${\mathit{O}}_{\mathit{\Sigma }}$ given below has the rank of n.

$${\mathit{O}}_{\mathit{\Sigma }}={\left.{\left[\begin{array}{ll}\begin{array}{lll}{\left(\nabla L_{\mathit{f}}^0\mathit{h}\right)}^T,& {\left(\nabla L_{\mathit{f}}^1\mathit{h}\right)}^T,& \cdots ,\end{array}& {\left(\nabla L_{\mathit{f}}^{n-1}\mathit{h}\right)}^T\end{array}\right]}^T\right|}_{\mathit{x}={\mathit{x}}_0}$$E12

It’s a heavy burden to calculate the observability matrix in Eq. (12) due to the existence of high order differential, especially for the 6-dimensional dynamics of Mars entry phase which requires the calculation of 5th order Lie derivatives. Next, a quadratic approximation method is developed to simplify the computation of the observability matrix.

First of all, the quadratic approximation of the kth order Lie derivative $L_{\mathit{f}}^k{h}_j$ is given as

$$L_{\mathit{f}}^k{h}_j\approx L_{\mathit{f}}^k{h}_{j0}+{\mathit{J}}_{\mathit{Lj}}^k\left(\mathit{x}-{\mathit{x}}_0\right)+\frac{1}{2}{\left(\mathit{x}-{\mathit{x}}_0\right)}^T{\mathit{H}}_{\mathit{Lj}}^k\left(\mathit{x}-{\mathit{x}}_0\right)$$E13

where $L_{\mathit{f}}^k{h}_{j0}$ is the value of $L_{\mathit{f}}^k{h}_j$ at ${\mathit{x}}_0$, and ${\mathit{J}}_L^k$and ${\mathit{H}}_L^k$ refer to, respectively, the Jacobian and Hessian matrix of $L_{\mathit{f}}^k{h}_j$ at ${\mathit{x}}_0$. The linearized state equation is given by

$$\mathit{f}\approx {\mathit{f}}_0+{\mathit{J}}_{\mathit{f}}\left(\mathit{x}-{\mathit{x}}_0\right)$$E14

in which ${\mathit{f}}_0$ refers to the value of $\mathit{f}$ at ${\mathit{x}}_0$, and ${\mathit{J}}_{\mathit{f}}$ is the Jacobi matrix of $\mathit{f}$ at ${\mathit{x}}_0$.

According to Eq. (9) and Eq. (13), the relationship between the kth and k + 1th order Lie derivative can be rewritten as

$$\begin{array}{l}{L}_{\mathit{f}}^{k+1}{h}_j=\nabla L_{\mathit{f}}^k{h}_j\cdot \mathit{f}\hfill \\ =\left[{\mathit{J}}_{\mathit{Lj}}^k+\frac{1}{2}{\left(\mathit{x}-{\mathit{x}}_0\right)}^T\left({\mathit{H}}_{\mathit{Lj}}^k+{\left({\mathit{H}}_{\mathit{Lj}}^k\right)}^T\right)\right]\left[{\mathit{f}}_0+{\mathit{J}}_{\mathit{f}}\left(\mathit{x}-{\mathit{x}}_0\right)\right]\begin{array}{lll}& & =\end{array}{\mathit{J}}_{\mathit{Lj}}^k{\mathit{f}}_0\\ +\left[{\mathit{J}}_{\mathit{Lj}}^k{\mathit{J}}_{\mathit{f}}+\frac{1}{2}{\left(\left[{\mathit{H}}_{\mathit{Lj}}^k+{\left({\mathit{H}}_{\mathit{Lj}}^k\right)}^T\right]{\mathit{f}}_0\right)}^T\right]\left(\mathit{x}-{\mathit{x}}_0\right)+\frac{1}{2}{\left(\mathit{x}-{\mathit{x}}_0\right)}^T\left[{\mathit{H}}_{\mathit{Lj}}^k+{\left({\mathit{H}}_{\mathit{Lj}}^k\right)}^T\right]{\mathit{J}}_{\mathit{f}}\left(\mathit{x}-{\mathit{x}}_0\right)\hfill \\ \begin{array}{lll}& & =L_{\mathit{f}}^{k+1}{h}_{j0}+{\mathit{J}}_{\mathit{Lj}}^{k+1}\left(\mathit{x}-{\mathit{x}}_0\right)+{\left(\mathit{x}-{\mathit{x}}_0\right)}^T{\mathit{H}}_{\mathit{Lj}}^{k+1}\left(\mathit{x}-{\mathit{x}}_0\right)\end{array}\hfill \end{array}$$E15

This leads to

$$\begin{array}{l}{L}_{\mathit{f}}^{k+1}{h}_{j0}={\mathit{J}}_{\mathit{Lj}}^k{\mathit{f}}_0\\ {\mathit{J}}_{\mathit{Lj}}^{k+1}={\mathit{J}}_{\mathit{Lj}}^k{\mathit{J}}_{\mathit{f}}+\frac{1}{2}{\left(\left[{\mathit{H}}_{\mathit{Lj}}^k+{\left({\mathit{H}}_{\mathit{Lj}}^k\right)}^T\right]{\mathit{f}}_0\right)}^T\\ {\mathit{H}}_{\mathit{Lj}}^{k+1}=\frac{1}{2}\left[{\mathit{H}}_{\mathit{Lj}}^k+{\left({\mathit{H}}_{\mathit{Lj}}^k\right)}^T\right]{\mathit{J}}_{\mathit{f}}\end{array}$$E16

The observability matrix can be computed as

$$\begin{array}{l}{\mathit{O}}_{\mathit{\Sigma }}={\left.{\left[\begin{array}{ll}\begin{array}{lll}{\left(\nabla L_{\mathit{f}}^0\mathit{h}\right)}^T,& {\left(\nabla L_{\mathit{f}}^1\mathit{h}\right)}^T,& \cdots ,\end{array}& {\left(\nabla L_{\mathit{f}}^{n-1}\mathit{h}\right)}^T\end{array}\right]}^T\right|}_{\mathit{x}={\mathit{x}}_0}\\ \begin{array}{ll}& =\end{array}{\left[\begin{array}{ll}\begin{array}{lll}{\left({\mathit{J}}_L^0\right)}^T,& {\left({\mathit{J}}_L^1\right)}^T,& \cdots ,\end{array}& {\left({\mathit{J}}_L^{n-1}\right)}^T\end{array}\right]}^T\end{array}$$E17

where ${\mathit{J}}_L^k={\left[\begin{array}{lll}{\left({\mathit{J}}_{L1}^k\right)}^T,& \cdots ,& {\left({\mathit{J}}_{\mathit{Lm}}^k\right)}^T\end{array}\right]}^T$.

Obtaining ${\mathit{J}}_{\mathit{Lj}}^0$ and ${\mathit{H}}_{\mathit{Lj}}^0$, the observability matrix can be iteratively calculated. Only 2nd order differential of $\mathit{h}$ is needed here to compute the Jacobian and Hessian matrices, reducing largely the computation cost.

Linearize the dynamical and observation model by first-order approximation

$$\begin{array}{l}\mathit{f}\approx {\mathit{f}}_0+{\mathit{J}}_{\mathit{f}}\left(\mathit{x}-{\mathit{x}}_0\right)\\ \mathit{h}\approx {\mathit{h}}_0+{\mathit{J}}_{\mathit{h}}\left(\mathit{x}-{\mathit{x}}_0\right)\end{array}$$E18

Construct the observability matrix according to the linear system theory

$${\mathit{O}}_{\mathit{\Sigma }}^l={\left[\begin{array}{ll}\begin{array}{lll}{\left({\mathit{J}}_{\mathit{h}}\right)}^T,& {\left({\mathit{J}}_{\mathit{h}}{\mathit{J}}_{\mathit{f}}\right)}^T,& \cdots ,\end{array}& {\left({\mathit{J}}_{\mathit{h}}^{n-1}{\mathit{J}}_{\mathit{f}}\right)}^T\end{array}\right]}^T$$E19

The Hessian matrix is involved in the quadratic approximation, improving the accuracy of observability analysis compared with the linearized observability analysis. However, the higher order terms of $\mathit{x}-{\mathit{x}}_0$ may appear when computing $L_{\mathit{f}}^{k+1}{h}_j$ in Eq.(9) if the state equation is approximated to a higher order. In this case, the predetermined presentation form in Eq.(13) is no longer valid. One way to defeat this case is to increase the approximation order of Lie derivatives. Note that tensor calculus can be involved and the computation complexity is increased. Thus, the trade between accuracy and computation cost is balanced by the quadratic approximation of Lie derivatives and the linearization of state equation.

In the optimization of observability, the condition number of observability matrix is selected as the performance index, given by

$$\mathit{\text{cond}}\left(\mathit{M}\right)=\frac{{\sigma }_{\max }\left(\mathit{M}\right)}{{\sigma }_{\min }\left(\mathit{M}\right)}$$E20

where ${\sigma }_{\max }$ and ${\sigma }_{\min }$ are, respectively, the maximum and minimum singular value of the matrix. The condition number measures the singularity of the matrix. A larger condition number means a more singular matrix. Here we take the inverse of the condition number to quantify the system observability.

$$\delta =\frac{1}{\mathit{\text{cond}}\left({\mathit{O}}_{\mathit{\Sigma }}\right)}=\frac{{\sigma }_{\min }\left({\mathit{O}}_{\mathit{\Sigma }}\right)}{{\sigma }_{\max }\left({\mathit{O}}_{\mathit{\Sigma }}\right)}$$E21

Obviously, the observability degree $\delta$ is in the interval [0, 1]. When $\delta =0$, the observability matrix is rank defect, and the navigation system is locally unobservable. When $\delta >0$, the observability is full rank, indicating an observable navigation system.

3.2. Observability analysis based on the fisher information matrix

Without loss of generality, we will consider the nonlinear observation models

$$\begin{array}{ll}{y}_i=h_i\left(\mathit{x}\right)+{\epsilon }_i,& i=1,\cdots ,N\end{array}$$E22

This equation may describe the measurement of relative range and range-rate according to Eq. (6) and (7). Meanwhile, in order to investigate the impact of different measurement methods on the observability of position and velocity of the entry vehicle separately, the 3-dimensional state x may be r or v of the entry vehicle. The likelihood function of x is defined as the joint probability density function of multiple measurements given by

$$L\left(y_1,\cdots ,y_3|\mathit{x}\right)=\prod _{i=1}^N\frac{1}{\sqrt{2\pi }{\sigma }_i}\text{exp}\left(-\frac{1}{2}{\sigma }_i^{-2}{‖y_i-h_i\left(\mathit{x}\right)‖}^2\right)$$E23

Then, take the negative of the natural log of Eq. (23) and omitting the terms not related to x, and the loss function can be derived as

$$J\left(\mathit{x}\right)=\frac{1}{2}\sum _{i=1}^N{\sigma }_i^{-2}{‖y_i-h_i\left(\mathit{x}\right)‖}^2$$E24

Find a state vector to minimize $J\left(\mathit{x}\right)$ and the state vector is the optimal estimation of the lander’s states. The FIM of the state is given by

$$\mathit{F}=E\left\{\frac{{\partial }^2}{\partial \mathit{x}\partial {\mathit{x}}^T}J\left(\mathit{x}\right)\right\}=\sum _{i=1}^N{\sigma }_i^{-2}\frac{\partial h_i\left(\mathit{x}\right)}{\partial \mathit{x}}{\left(\frac{\partial h_i\left(\mathit{x}\right)}{\partial \mathit{x}}\right)}^T$$E25

The estimate error covariance and FIM satisfy the following equation

$$\mathit{P}\ge {\mathit{F}}^{-1}$$E26

where $\mathit{P}$ is the estimate error covariance, and “$\ge$” means that ($\mathit{P}-{\mathit{F}}^{-1}$) is positive semidefinite. According to Eq. (26), the FIM can be used to evaluate the lower bound of the estimation error covariance, and further the system observability. Give the trace of ${\mathit{F}}^{-1}$ in Eq. (27).

$$\mathrm{tr}\left({\mathit{F}}^{-1}\right)=\sum _{i=1}^3\frac{1}{{\lambda }_i}$$E27

where ${\lambda }_i\left(i=1,2,3\right)$ are the eigenvalues of F. It’s illustrated from Eq. (27) that larger eigenvalues of the FIM leads to smaller trace of estimation error covariance and stronger system observability. Quantify the observability by the determinant of FIM $\text{det}\left(\mathit{F}\right)=\prod _{i=1}^3{\lambda }_i$. The following relationship can be obtained.

$$\mathrm{tr}\left(\mathit{P}\right)\ge \mathrm{tr}\left({\mathit{F}}^{-1}\right)=\sum _{i=1}^3\frac{1}{{\lambda }_i}>\frac{3}{\sum _{i=1}^3{\lambda }_i}=\frac{3}{\mathrm{tr}\left(\mathit{F}\right)}$$E28

Eq. (28) means that the trace of FIM measures the lower bound of estimation errors.

4. Observability analysis of Mars networks-based navigation

4.1. Observability analysis using only range measurements

In this subsection, the system observability using only range measurements between the lander and ground beacons is analyzed. Since no velocity information is included in Eq. (6), only the observability of the position vector is studied. The cases with different amount of beacons are studied.

4.2. One-beacon case

In this case, the FIM is given by

$${\mathit{F}}_1={\sigma }_{R1}^{-2}\frac{\partial R_1\left(\mathit{r}\right)}{\partial \mathit{r}}{\left(\frac{\partial R_1\left(\mathit{r}\right)}{\partial \mathit{r}}\right)}^T={\sigma }_{R1}^{-2}{\mathit{n}}_1{\mathit{n}}_1^T\equiv {\sigma }_{R1}^{-2}{\mathit{N}}_1$$E29

The rank of the matrix ${\mathit{N}}_1$is only one. Solving the following equation

$$\text{det}\left(\lambda {\mathit{I}}_{3\times 3}-{\mathit{N}}_1\right)=0$$E30

Clearly, the eigenvalues of ${\mathit{N}}_1$ are given by twice repeated 0 and $n_{1x}^2+n_{1y}^2+n_{1z}^2=1$. Therefore, the eigenvalues of ${\mathit{F}}_1$ are given by $\begin{array}{ll}{\lambda }_1={\lambda }_2=0,& {\lambda }_3={\sigma }_{R1}^{-2}\end{array}$.

Next, we have the eigenvector corresponding to ${\lambda }_3$

$${\mathit{w}}_3=\frac{1}{n_{1z}}{\left[\begin{array}{lll}{n}_{1x},& n_{1y},& n_{1z}\end{array}\right]}^T=\frac{1}{n_{1z}}{\mathit{n}}_1$$E31

The vector ${\mathit{w}}_3$ corresponds to the observable state combination, and means that only the state component along the vector ${\mathit{n}}_1$ can be observable.

According to Eq. (28), the lower bound of estimation errors can be obtained as

$$\frac{3}{\mathrm{tr}\left({\mathit{F}}_1\right)}=3{\sigma }_{R1}^2$$E32

Eq. (32) means the lower bound of estimation errors is higher than the estimation accuracy. In another word, the estimation accuracy cannot be higher than the measurement accuracy. Note that, even if multiple beacons are involved in the navigation system, the observability is still deteriorated if the beacons are located in similar direction.

4.3. Two-beacon case

Assume two non-collinear beacons, the FIM in Eq. (25) is derived by

$${\mathit{F}}_2=\sum _{i=1}^2{\sigma }_{\mathit{Ri}}^{-2}\frac{\partial R_i\left(\mathit{r}\right)}{\partial \mathit{r}}{\left(\frac{\partial R_i\left(\mathit{r}\right)}{\partial \mathit{r}}\right)}^T=\sum _{i=1}^2{\sigma }_{\mathit{Ri}}^{-2}{\mathit{n}}_i{\mathit{n}}_i^T$$E33

Involving one more measurement, the rank of ${\mathit{F}}_2$ is increased to two. The observable state combinations can be obtained by solving the eigenvalue and eigenvector. In this case, the eigenvalues of ${\mathit{F}}_2$ are given by ${\lambda }_1={\lambda }_2\ne 0$, ${\lambda }_3=0$. The eigenvector corresponding to the zero eigenvalue is obtained as

$${\mathit{w}}_3={\left[\begin{array}{lll}\frac{n_{1y}{n}_{2z}-n_{1z}{n}_{2y}}{n_{1x}{n}_{2y}-n_{1y}{n}_{2x}},& \frac{n_{1z}{n}_{2x}-n_{1x}{n}_{2z}}{n_{1x}{n}_{2y}-n_{1y}{n}_{2x}},& 1\end{array}\right]}^T=\frac{1}{n_{1x}{n}_{2y}-n_{1y}{n}_{2x}}{\mathit{n}}_1\times {\mathit{n}}_2$$E34

The vector ${\mathit{w}}_3$ gives the unobservable state component which is in the direction perpendicular to the plane constructed by ${\mathit{n}}_1$ and ${\mathit{n}}_2$. From an opposite view, all state components in plane are observable.

Since, in this case, the observability matrix is still zero, the navigation system is unobservable. According to Eq.(28), the lower bound of the estimation errors can be obtained as

$$\frac{3}{\mathrm{tr}\left({\mathit{F}}_2\right)}=\frac{3}{\sum _{i=1}^2{\sigma }_{\mathit{Ri}}^{-2}\left(n_{\mathit{ix}}^2+n_{\mathit{iy}}^2+n_{\mathit{iz}}^2\right)}=\frac{3}{\sum _{i=1}^2{\sigma }_{\mathit{Ri}}^{-2}}\ge \frac{3{\sigma }_{R\text{min}}^2}{2}$$E35

where ${\sigma }_{R\text{min}}$ is the smaller standard deviation among ${\sigma }_{R1}$ and ${\sigma }_{R2}$. It’s known by comparing Eqs. (32) and (35) that the estimation accuracy can be improved by using one more radio beacon.

4.4. More-than-two-beacon case

In this case, the FIM is given by

$$\begin{array}{ll}{\mathit{F}}_N=\sum _{i=1}^N{\sigma }_{\mathit{Ri}}^{-2}\frac{\partial R_i\left(\mathit{r}\right)}{\partial \mathit{r}}{\left(\frac{\partial R_i\left(\mathit{r}\right)}{\partial \mathit{r}}\right)}^T=\sum _{i=1}^N{\sigma }_{\mathit{Ri}}^{-2}{\mathit{n}}_i{\mathit{n}}_i^T,& N\ge 3\end{array}$$E36

The matrix ${\mathit{F}}_N$ has a full rank, indicating an observable system. The determinant of ${\mathit{F}}_N$ is given in Eq. (37).

$$\text{det}\left(F_N\right)={\displaystyle \sum _{1\le k_3<k_2<k_1\le N}{\sigma }_{R{k}_1}^{-2}{\sigma }_{R{k}_2}^{-2}{\sigma }_{R{k}_3}^{-2}{\left[n_{k_1}\cdot \left(n_{k_2}\times n_{k_3}\right)\right]}^2}$$E37

The detailed derivation can be found in Ref. [14]. From Eq.(7), we can know that more radio beacons, no matter where they are, increase the determinant of the FIM, thus increase the system observability. To analyze the maximum value of $\text{det}\left({\mathit{F}}_N\right)$, Eq. (37) is reorganized as

$$\text{det}\left({\mathit{F}}_N\right)\le {\sigma }_{R\text{min}}^{-6}\sum _{1\le k_3<k_2<k_1\le N}{\left[{\mathit{n}}_{k_1}\cdot \left({\mathit{n}}_{k_2}\times {\mathit{n}}_{k_3}\right)\right]}^2$$E38

where ${\sigma }_{R\text{min}}$ is the minimum value among ${\sigma }_{\mathit{Ri}}$. The selection of the direction of radio beacons to maximize the observability can be described by the following optimization problem

$$\begin{array}{l}\text{max}\begin{array}{ll}& \begin{array}{ll}& \end{array}\end{array}\sum _{1\le k_3<k_2<k_1\le N}{\left[{\mathit{n}}_{k_1}\cdot \left({\mathit{n}}_{k_2}\times {\mathit{n}}_{k_3}\right)\right]}^2\\ \text{subject to}\begin{array}{ll}& \end{array}\begin{array}{ll}‖{\mathit{n}}_i‖=1,& i=1,\cdots ,N\end{array}\end{array}$$E39

Note that the locations of radio beacons are not constrained. In cases with three beacons, the determinant of ${\mathit{F}}_3$ is maximized if and only if n₁, n₂, and n₃ are orthogonal to each other. However, no analytic results can be obtained when there are more than three beacons. Thus, a Genetic Algorithm is exploited to solve the optimization problem. The maximum determinants are listed in Table 1.

Number of beacons	Maximum determinant of FIM
3	1.000 ${\sigma }_{R\text{min}}^{-6}$
4	2.3704 ${\sigma }_{R\text{min}}^{-6}$
5	4.6296 ${\sigma }_{R\text{min}}^{-6}$
6	8.0000 ${\sigma }_{R\text{min}}^{-6}$
7	12.7037 ${\sigma }_{R\text{min}}^{-6}$
8	18.9630 ${\sigma }_{R\text{min}}^{-6}$

Table 1.

Maximum determinants of FIM related to different number of beacons.

According to the results in Table 1, the relationship between the maximum determinant and the number of beacons can be induced by an exponential formulation, given by

$$\text{det}{\left({\mathit{F}}_N\right)}_{\max }={\left(\frac{N}{3{\sigma }_{R\text{min}}^{-2}}\right)}^3$$E40

The lower bound of estimation errors is derived as

$$\frac{3}{\mathrm{tr}\left({\mathit{F}}_N\right)}=\frac{3}{\sum _{i=1}^N{\sigma }_{\mathit{Ri}}^{-2}\left(n_{\mathit{ix}}^2+n_{\mathit{iy}}^2+n_{\mathit{iz}}^2\right)}\ge \frac{3{\sigma }_{R\text{min}}^2}{N}$$E41

The change of lower bound of estimation errors with number of beacons is shown in Figure 1. It’s shown that with more beacons comes more accurate estimation. However, the increasing rate of accuracy is slowed down, indicating that the navigation accuracy cannot be improved endlessly by only increasing the number of beacons.

Figure 1.

Lower bound of estimation errors with beacon number.

4.5. Observability analysis of the navigation using range-rate measurements

4.5.1. Observability analysis of vehicle’s velocity

The FIM of vehicle’s velocity using range-rate data is given by

$$\begin{array}{ll}{\mathit{F}}_N=\sum _{i=1}^N{\sigma }_{\mathit{Vi}}^{-2}\frac{\partial V_i\left(\mathit{r},\mathit{v}\right)}{\partial \mathit{v}}{\left(\frac{\partial V_i\left(\mathit{r},\mathit{v}\right)}{\partial \mathit{v}}\right)}^T=\sum _{i=1}^N{\sigma }_{\mathit{Vi}}^{-2}{\mathit{n}}_i{\mathit{n}}_i^T,& N\ge 1\end{array}$$E42

Eq. (42) has a similar form with Eq. (36) which describes the FIM of position. The only difference lies in the measurement deviation. Hence the same conclusion of the observability of velocity can be obtained as that in Section 4.1. The detailed analysis is omitted here.

4.6. Observability analysis of vehicle’s position

Using the range-rate measurements, the FIM of the lander’s position is derived as

$$\begin{array}{ll}{\mathit{F}}_N=\sum _{i=1}^N{\sigma }_{\mathit{Vi}}^{-2}\frac{\partial V_i\left(\mathit{r},\mathit{v}\right)}{\partial \mathit{r}}{\left(\frac{\partial V_i\left(\mathit{r},\mathit{v}\right)}{\partial \mathit{r}}\right)}^T=\sum _{i=1}^N{\sigma }_{\mathit{Vi}}^{-2}{\mathit{L}}_i\mathit{v}{\mathit{v}}^T{\mathit{L}}_i^T,& N\ge 1\end{array}$$E43

where ${\mathit{L}}_i$ is given by

$${\mathit{L}}_i=\frac{1}{R_i}\left[\begin{array}{lll}{n}_{\mathit{iy}}^2+n_{\mathit{iz}}^2& -n_{\mathit{ix}}{n}_{\mathit{iy}}& -n_{\mathit{ix}}{n}_{\mathit{iz}}\\ -n_{\mathit{iy}}{n}_{\mathit{ix}}& n_{\mathit{iz}}^2+n_{\mathit{ix}}^2& -n_{\mathit{iy}}{n}_{\mathit{iz}}\\ -n_{\mathit{iz}}{n}_{\mathit{ix}}& -n_{\mathit{iz}}{n}_{\mathit{iy}}& n_{\mathit{ix}}^2+n_{\mathit{iy}}^2\end{array}\right]$$E44

The FIM here is much more complicated than that in Section 4.1 due to the involvement of both range and velocity information in FIM. Define ${\mathit{V}}_i={\mathit{L}}_i\mathit{v}$, the following equation can be obtained.

$$\begin{array}{ll}{\mathit{F}}_N=\sum _{i=1}^N{\sigma }_{\mathit{Vi}}^{-2}{\mathit{V}}_i{{\mathit{V}}_i}^T,& N\ge 1\end{array}$$E45

When there are one or two beacons, the FIM is rank defect, and the navigation system is also unobservable. With three or more beacons comes the full-rank FIM. In this section, only the observable cases are focused on.

It is also concluded that the determinant of FIM will be zero if only one or two beacons is used, which indicates that the position of entry vehicle will be observable if more than two beacons are used. Furthermore, we focus on three-beacon and more-than-three beacon cases. Compare Eq. (46) with Eq. (36), we can find that the determinant of FIM for range-rate measurement cases has a similar format as Eq. (37)

$$\begin{array}{ll}\text{det}\left({\mathit{F}}_N\right)=& \sum _{1\le k_3<k_2<k_1\le N}{\sigma }_{{\mathit{Vk}}_1}^{-2}{\sigma }_{{\mathit{Vk}}_2}^{-2}{\sigma }_{{\mathit{Vk}}_3}^{-2}{\left[{\mathit{V}}_{k_1}\cdot \left({\mathit{V}}_{k_2}\times {\mathit{V}}_{k_3}\right)\right]}^2\le {\sigma }_{V\text{min}}^{-6}\sum _{1\le k_3<k_2<k_1\le N}{\left[{\mathit{V}}_{k_1}\cdot \left({\mathit{V}}_{k_2}\times {\mathit{V}}_{k_3}\right)\right]}^2,N\ge 3\hfill \end{array}$$E46

It’s shown that Eq.(46) has the similar format with Eq. (37). Thus, the change of the observability with the number of radio beacons is similar with the results in Table 1. However, due to involving relative range and velocity information, the optimal geometric configuration is different with the cases using only range measurements.

The lower bound of estimation errors in this case is evaluated by

$$\begin{array}{l}\frac{3}{\mathrm{tr}\left({\mathit{F}}_N\right)}=\frac{3}{\sum _{i=1}^N{\sigma }_{\mathit{Vi}}^{-2}{{\mathit{V}}_i}^T{\mathit{V}}_i}\\ =\frac{3}{\sum _{i=1}^N\frac{{\sigma }_{\mathit{Vi}}^{-2}}{R_i^2}\left\{\begin{array}{l}{\left[v_x\left(n_{\mathit{iy}}^2+n_{\mathit{iz}}^2\right)-v_y{n}_{\mathit{ix}}{n}_{\mathit{iy}}-v_z{n}_{\mathit{ix}}{n}_{\mathit{iz}}\right]}^2+{\left[-v_x{n}_{\mathit{ix}}{n}_{\mathit{iy}}+v_y\left(n_{\mathit{ix}}^2+n_{\mathit{iz}}^2\right)-v_z{n}_{\mathit{iy}}{n}_{\mathit{iz}}\right]}^2+\\ {\left[-v_x{n}_{\mathit{ix}}{n}_{\mathit{iz}}-v_y{n}_{\mathit{iy}}{n}_{\mathit{iz}}+v_z\left(n_{\mathit{ix}}^2+n_{\mathit{iy}}^2\right)\right]}^2\end{array}\right\}}\\ >\frac{3}{\sum _{i=1}^N\frac{{\sigma }_{\mathit{Vi}}^{-2}}{R_i^2}\left\{{\left[v_x\left(1+\frac{1}{2}\left(n_{\mathit{iy}}^2+n_{\mathit{iz}}^2\right)\right)+v_y\left(1+\frac{1}{2}\left(n_{\mathit{ix}}^2+n_{\mathit{iz}}^2\right)\right)+v_z\left(1+\frac{1}{2}\left(n_{\mathit{ix}}^2+n_{\mathit{iy}}^2\right)\right)\right]}^2\right\}}\\ >\frac{3}{\frac{3}{2}\sum _{i=1}^N\frac{{\sigma }_{\mathit{Vi}}^{-2}}{R_i^2}{\left(v_x+v_y+v_z\right)}^2}\end{array}$$E470

$$\begin{array}{l}\ge \frac{2}{3\sum _{i=1}^N{\sigma }_{\mathit{Vi}}^{-2}\frac{v^2}{R_i^2}}\end{array}$$E47

where $v=\sqrt{v_x^2+v_y^2+v_z^2}$ is the lander’s velocity value. Obviously, more radio beacons lead to more accurate estimation. Since the value of relative range is much bigger than relative velocity, the lower bound of estimation errors using range-rate date is larger than that using range data. Besides, it’s concluded that more accurate range-rate measurement, closer relative range, and slower velocity can realizer more accurate position estimation.

5. Orbit optimization based on observability analysis

5.1. Optimization of navigation using ground beacons

The configuration radio beacons is expressed by the following set

$$\mathit{C}=\left\{\left.{\mathit{p}}_B^i\right|i=1,\cdots ,l\right\}$$E48

where ${\mathit{p}}_B^i={\left[\begin{array}{lll}{x}_B^i,& y_B^i,& z_B^i\end{array}\right]}^T$is the position of the i^th beacon. Considering the time-varying observability, the minimum value of the observability in the entry phase is taken as the optimization performance index.

$$D\left(\mathit{C}\right)=\underset{\mathit{x}\in {\mathit{T}}_{\mathit{x}}}{\text{min}}\delta$$E49

To realize the Mars network-based navigation, the visibility of the beacons to the lander should be guaranteed. Define two unit vectors as follows

$$\begin{array}{ll}{\mathit{n}}_{\mathit{Bi}}=\frac{{\left[\begin{array}{lll}{x}_B^i,& y_B^i,& z_B^i\end{array}\right]}^T}{\sqrt{{\left(x_B^i\right)}^2+{\left(y_B^i\right)}^2+{\left(z_B^i\right)}^2}},& {\mathit{n}}_C=\frac{{\left[\tilde{x},,,\tilde{y},,,\tilde{z}\right]}^T}{\sqrt{{\tilde{x}}^2+{\tilde{y}}^2+{\tilde{z}}^2}}\end{array}$$E50

where ${\left[\tilde{x},,,\tilde{y},,,\tilde{z}\right]}^T$is the relative position vector from the lander to the radio beacon, obtained as

$${\left[\tilde{x},,,\tilde{y},,,\tilde{z}\right]}^T={\left[\begin{array}{lll}x,& y,& z\end{array}\right]}^T-{\left[\begin{array}{lll}{x}_B^i,& y_B^i,& z_B^i\end{array}\right]}^T$$E51

To guarantee the visibility, the two vectors in Eq. (50) should satisfy

$$\begin{array}{ll}\mathrm{arc}\text{cos}\left({\mathit{n}}_{\mathit{Bi}}\cdot {\mathit{n}}_C\right)<\frac{\pi }{2},& \mathit{x}\in {\mathit{T}}_{\mathit{x}}\end{array}$$E52

The schematic of visibility is shown in Figure 2.

Figure 2.

Principle of the line-of-sight visibility.

The optimization problem of beacon configuration is given as

$$\begin{array}{l}\text{max}\begin{array}{ll}& D\end{array}\left(\mathit{C}\right)\\ \begin{array}{ll}s.t.\begin{array}{ll}& \end{array}{\mathit{p}}_B^i\in \mathit{\Omega },& i=1,\cdots ,l\end{array}\end{array}$$E53

where $\mathit{\Omega }$ is the set of the areas of radio beacons that satisfy the visibility during the whole entry phase. In this optimization problem, the global optimization algorithm is selected to obtain the optimal beacon configuration.

The initial states of the lander are listed in Table 2.

Initial state	$r\left(0\right)$km	$\theta \left(0\right)$deg	$\phi \left(0\right)$deg	$V\left(0\right)$m/s	$\gamma \left(0\right)$deg	$\psi \left(0\right)$deg
Value	3518.2	−89.872	−28.02	5515	−11.8	5.156

Table 2.

Initial states of the lander.

It is assumed that the Mars entry phase lasts approximately 240 seconds. The entry trajectory and the corresponding visible area are shown in Figure 3.

Figure 3.

Entry trajectory and the visible area.

Three radio ranging measurements at a certain time can geometrically determine the position of the lander. Thus the navigation scenario with three beacons is first analyzed with respect to the observability. The optimal locations of beacons are displayed in Figure 4.

Figure 4.

Optimal configuration for the scenario with three beacons.

The optimal three beacons are located close to the edge of both sides of the visible area. The beacon on the east side is almost along the entry trajectory, while the west two beacons are separated on the north and south side of the entry trajectory. The observability degree in this situation calculated by different methods is illustrated in Figure 5, and the computation time for each method is listed in Table 3.

Figure 5.

Degree of observability for the optimal scenario with three beacons.

Analysis approach	Computation time, s
Method with Lie algebra	>10,000
Linearization method	1.3987
Method based on quadratic approximation	2.1558

Table 3.

Computation time for each approach.

Figure 5 shows a huge undulation in observability degree during the Marts entry phase. The maximum and minimum value are $1.413\times {10}^{‐8}$ and $2.945\times {10}^{‐7}$ respectively. Considering the machine precision, the navigation system is observable only if the observability exceeds $1\times {10}^{‐16}$. The observability degree during the entire entry phase passes through the threshold, and thus, the navigation system is observable. The minimum degree of observability occurs at the beginning of the entry phase when the entry vehicle is at its greatest distance from radio beacons, while the maximum degree of observability occurs when the entry vehicle approaches two beacons on the west side. In order to explain the evolution of the degree of observability. An index related to the geometric configuration of the lander and radio beacons is given in Eq. (54) to explain the evolution of the degree of observability.

$$I=\sum _{1\le i<j<k\le N}{\left[{\mathit{n}}_i\cdot \left({\mathit{n}}_j\times {\mathit{n}}_k\right)\right]}^2$$E54

where ${\mathit{n}}_i$, ${\mathit{n}}_j$, and ${\mathit{n}}_k$ are the unit vectors from the beacon to the lander, N is the number of beacons. The evolution of index I is displayed in Figure 6, showing an identical variation trend with observability degree and backing up the observability analysis conclusion.

Figure 6.

The value of I for the optimal scenario with three beacons.

The observability degree obtained from the three methods is quite close to each other. However, the method based on Lie algebra consumes the most time. The linearization method provides the largest deviations, especially at the peak time, indicating a relatively low accuracy. The proposed quadratic approximation method achieves a performance balance in accuracy and complexity. To analyze the navigation accuracy, the Extended Kalman Filter (EKF) is used to estimate the lander’s states. The range measurement error is assumed to be Gaussian white noise with a standard deviation of 100 m. The initial errors are listed in Table 4. The estimation errors and the 1-sigma uncertainty bounds are depicted in Figure 7.

Initial state	$r\left(0\right)$m	$\theta \left(0\right)$deg	$\phi \left(0\right)$deg	$V\left(0\right)$m/s	$\gamma \left(0\right)$deg	$\psi \left(0\right)$deg
Error	1000	0.2	0.2	10	0.2	0.2

Table 4.

Errors of initial states.

Figure 7.

Navigation results for the optimal scenario with three beacons.

It’s illustrated that θ and $\varphi$ have the most accurate estimation and the fastest convergence. The convergence of the states V, γ, and Ψ is relatively slow at the beginning of the Mars entry phase due to the weak observability. With the increase of the observability degree comes the rapid convergence of the uncertainty bounds and the state estimation errors from about 90 to 115 seconds. The max deceleration of the lander also contributes to the rapid convergence.

5.2. Optimization of navigation using Mars orbiters

Compared with ground beacons, the Mars obiters are constrained by the orbital dynamics, which is considered to be two-body dynamics here. In this subsection, the initial states of the Mars orbiters are considered as the optimized variables. Furthermore, assuming that the Mars orbiters moves in a circular orbit, the variables to be optimized are simplified as inclination i, longitude of ascending node Ω, and the true anomaly f. The initial states of the orbiter can be expressed by the optimized variables, given by

$$\left\{\begin{array}{l}{\mathit{r}}_{i0}^B=\left(R_M+a_i\right)\text{cos}{f}_i{\mathit{P}}_i+\left(R_M+a_i\right)\text{sin}{f}_i{\mathit{Q}}_i\\ {\mathit{v}}_{i0}^B=-\sqrt{\mu /\left(R_M+a_i\right)}\text{sin}{f}_i{\mathit{P}}_i+\sqrt{\mu /\left(R_M+a_i\right)}\text{cos}{f}_i{\mathit{Q}}_i\end{array}\right.$$E55

where R_M is the radius of Mars, a_i is the orbit altitude, and P_i and Q_i are given by

$$\left\{\begin{array}{l}{\mathit{P}}_i={\left[\begin{array}{lll}\text{cos}{\mathrm{\Omega }}_i,& \text{sin}{\mathrm{\Omega }}_i,& 0\end{array}\right]}^T\\ {\mathit{Q}}_i={\left[\begin{array}{lll}-\text{sin}{\mathrm{\Omega }}_i\text{cos}{i}_i,& \text{cos}{\mathrm{\Omega }}_i\text{cos}{i}_i,& \text{sin}{i}_i\end{array}\right]}^T\end{array}\right.$$E56

Given the initial states of the Mars obiter, the subsequent states can be obtained by propagating the two-body dynamics. Likewise, the trajectory of the lander can be also obtained by propagating the entry dynamics. To evaluate the overall performance of the observability of the entry phase, the integration of the observability is taken as the performance index, given by

$$I\left(\mathit{e}\right)=\underset{t=0}{\overset{t_f}{\int }}O\left(t\right)\mathit{dt}$$E57

where t_f is the final time of entry phase, $\mathit{e}=\left\{\begin{array}{lll}{\mathit{e}}_1,& \cdots ,& {\mathit{e}}_n\end{array}\right\}$, ${\mathit{e}}_i={\left[\begin{array}{lll}{\mathrm{\Omega }}_i,& i_i,& f_i\end{array}\right]}^T$denotes the optimization variables. Similar to ground beacon-based navigation, the visibility between the lander and the Mars obiters should be also guaranteed. Define two angles as follows:

$$\left\{\begin{array}{l}{\theta }_0=\text{arccos}\left(\frac{R_M}{‖\mathit{r}‖}\right)\\ {\theta }_{1i}=\text{arccos}\left(\frac{R_M}{‖{\mathit{r}}_i^B‖}\right)\end{array}\right.$$E58

The angle between the position vectors r and ${\mathit{r}}_i^B$ is given by

$${\theta }_i=\text{arccos}\left(\frac{{\mathit{r}}_i^B\cdot \mathit{r}}{‖{\mathit{r}}_i^B‖‖\mathit{r}‖}\right)$$E59

The visibility requires that

$$\begin{array}{ll}{\theta }_i<{\theta }_0+{\theta }_{1i},& t\in \left[\begin{array}{ll}0,& t_f\end{array}\right]\end{array}$$E60

The schematic of the visibility is illustrated in Figure 8. The gray part represents the area in which the Mars orbiter is invisible to the lander.

Figure 8.

The schematic of the visibility.

Then the orbit optimization problem is given by

$$\begin{array}{l}\text{max}\begin{array}{ll}& \begin{array}{ll}& \end{array}I\left(\mathit{e}\right)=\underset{t=0}{\overset{t_f}{\int }}O\left(t\right)\mathit{dt}\end{array}\\ \begin{array}{ll}\text{subject to}\begin{array}{ll}& \end{array}\begin{array}{ll}{\theta }_i<{\theta }_0+{\theta }_{1i},& t\in \left[\begin{array}{ll}0,& t_f\end{array}\right],\end{array}& i=1,\cdots ,n\end{array}\end{array}$$E61

In the optimization problem, the performance index cannot be expressed explicitly by the optimization variables, and the gradient cannot be obtained. Thus, the heuristic global optimization algorithm is chosen to solve the optimization problem. The lander’s initial states are listed in Table 5 with the assumption of a ballistic entry having a banking angle of zero. The duration of entry phase is setup as 240 seconds.

State	Value	Unit
x	−3.92	km
y	−3103.37	km
z	−1665.41	km
vx	5775.31	m/s
vy	1124.27	m/s
vz	1175.48	m/s

Table 5.

Initial states of the lander.

The navigation scenario with three Mars orbiters is analyzed. The nominal orbit altitude of the three orbiters is 725 km. The observability is quantified by

$$O=\text{det}\left({\mathit{N}}_3\right)={\left[{\mathit{n}}_1\cdot \left({\mathit{n}}_2\times {\mathit{n}}_3\right)\right]}^2$$E62

At a certain epoch, the maximum value of O is 1 when and only when three unit vectors n₁, n₂, and n₃ are orthogonal to each other. Considering the overall observability of the entry phase, the orbits of the Mars orbiters are optimized and shown in Figure 9, and the optimal initial elements are listed in Table 6.

Figure 9.

The optimal orbits of three orbiters.

Orbit element	Obiter 1	Orbiter 2	Orbiter 3
Ω (deg)	49.329	16.136	36.562
i (deg)	24.209	35.889	18.901
f (deg)	240.219	256.141	229.294

Table 6.

Initial orbit elements of three orbiters.

It’s shown that the three orbiters keep a relatively stable configuration, and stays orthogonal approximately to each other. The value of maximized performance index is 237.963. The observability almost reaches the maximum value all the time during the Mars entry phase. The comparison of Mars obiters-based navigation and ground beacon-based navigation is performed. The observability degree of these two scenarios is shown in Figure 10.

Figure 10.

Degree of observability in two navigation schemes.

The fixed ground beacons have limited locations due to the visibility constrain and the geometric configuration cannot remain optimal during the entry phase. Thus, the observability is undulated to a large extent. The Mars orbiters overcome this defect with its moving property. To show straightforward the geometric configuration., the observability degree is close to maximum value at each epoch during the Mars entry phase. The angles between the vectors n₁, n₂, and n₃ are depicted in Figure 11.

Figure 11.

Angles between three unit vectors in two navigation schemes.

It’s shown that, using the ground beacons, the angles between the three vectors change dramatically in the entry phase. The optimal configuration can be met only at the epoch of 75 s. However, for the orbiter-based navigation scheme, n₁, n₂, and n₃ are almost orthogonal throughout the entry phase. The advantages of orbiter-based navigation scheme in the configuration and observability performance improve the navigation capability.

Next, 500-time Monte Carlo simulations of navigation systems based on EKF are carried out. The initial position and velocity have standard deviations of 1 km and 0.5 m/s respectively. The measurement error is set to be 50 m, and considered as Gaussian white noise. The simulation results are shown in Figure 12.

Figure 12.

1σ error bounds of states in two navigation schemes.

Since no information of entry vehicle’s velocity is provided from range measurements, the convergence of velocity estimation is not as quick as position estimation. A much better navigation performance can be achieved by the Mars orbiter-based navigation. It can be concluded that the configuration of orbiters is a main contributor to the navigation performance. The Mars orbiter-based navigation, which can achieve a better configuration, is more practical for Mars entry navigation.

6. Conclusions

This chapter introduced the Mars Networks-based navigation for the Mars entry phase. Based on the navigation scheme, the observability of the navigation system was analyzed using the proposed two novel observability analysis methods. Furthermore, the beacon configuration was optimized based on observability considering the line-of-sight constraints were concluded that the beacon configuration is a main contributor to the Mars Networks-based navigation. The observability analysis showed that an improved behavior of observability and more flexibility of beacon configuration determination can be achieved using more beacons. Navigation also demonstrated this conclusion. Meanwhile, compared with the ground beacons, Mars orbiters may be a better choice as Mars Network which gives a more accurate navigation result.