Parameter-Independent Event-Triggered Implicit UKF for the Celestial Navigation Using Time Delay Measurement

Abstract

Celestial navigation using time delay measurement is an innovative autonomous navigation method. To calculate the equivalent measurement, the numerical method needs to be applied, which is time-consuming. The event-triggered mechanism intermittently and aperiodically processes measurements by judging if the update error has changed drastically. However, its performance is greatly affected by the constant threshold. To solve this problem, a parameter-independent event-triggered implicit unscented Kalman filter (UKF) is proposed and applied to the celestial navigation using time delay measurement. The innovation at the current moment and the updated estimate covariance at the last moment are compared with the previous value instead of the constant threshold. The event is automatically triggered when the accuracy of the state estimate is low. Simulation results indicate that the proposed parameter-independent event-triggered implicit UKF can reduce the running time by reducing unnecessary measurement updates, whose performance will not be affected by any parameter or window size. In a word, the proposed method substitutes the dynamic threshold for the constant threshold, ensuring that its performance will not be affected by any parameter or window size.

1. Introduction

Deep space exploration is an important indicator of a country’s comprehensive strength and technology level. For the success of a deep space exploration mission, accurate, prompt, and dependable navigation information is essential [1,2]. With the increase in the distance between the spacecraft and the earth, the delay caused by the long roundtrip communication distance becomes an obstacle to the real-time navigation of the ground tracking system [3]. The Sun transit could result in the outage of the spacecraft’s communication links. Furthermore, substantial spacecraft greatly increase the burden and cost of supporting this system.

All these objections can be circumvented if the spacecraft has autonomous navigation capabilities. Celestial navigation is a suitable autonomous navigation method for deep space exploration [4,5]. Commonly used celestial navigation measurements include star angle [6,7], pulsar time of arrival (TOA) [8,9], and Doppler velocity [10,11]. Solar oscillation time delay is an innovative celestial navigation measurement [12,13]. The solar oscillation results in dramatic changes in the intensity and the spectral central wavelength of sunlight [14,15]. Two atomic resonance spectrometers pointing to the Sun and the reflecting celestial body simultaneously detect the spectral central wavelength of sunlight and record the time. The spectral central wavelength of the directly received sunlight can be compared with that of the sunlight reflected by the nearby celestial body to obtain the corresponding time delay. The spacecraft’s position information with respect to the nearby celestial body can be provided by the time delay measurement. The measurement model of the time delay is an implicit function, and thus the implicit unscented Kalman filter (IUKF) is applied to acquire the state estimate [16].

To calculate the equivalent measurement, numerical methods such as the dichotomy method need to be applied to figure out the equation set, which is time-consuming. In fact, real-time performance is as important as accuracy for the autonomous navigation of deep space probes. The probe flies at a speed of tens of kilometers per second during the planetary capture segment, and thus it needs to quickly react based on navigation information. This requires the navigation information to be solved in a short time. Normally, the real-time performance and accuracy of navigation cannot be optimal at the same time. In addition to the above two aspects, factors such as hardware volume, weight, cost, power consumption, etc. should also be considered. It is imperative to find a balance between these factors. In other words, it is necessary to reduce the amount of computation while maintaining high accuracy.

Some measurements can hardly provide valuable information and need not be transmitted and processed. The event-triggered mechanism is intermittent aperiodic sampled data, devoted to a desirable compromise between the real-time performance and accuracy of navigation [17,18,19,20]. It was first introduced into state estimation in Ref. [21] and was shown to outperform periodic sampling at the same sampling rate. For nonlinear systems, the event-triggered extended Kalman filter [22,23], the event-triggered unscented Kalman filter [24], and the event-triggered cubature Kalman filter [25,26] were successively proposed. Ref. [27] proposed a distributed UKF algorithm based on consensus with an event-triggering communication mechanism for multiple unmanned aerial vehicles. Ref. [28] proposed a nonlinear stochastic event-triggered estimator based on UKF for controllable and uncontrollable systems. Ref. [29] designed an event-triggered orbit estimator for a spacecraft with intermittent sensor measurements. To deal with the implicit measurement model, an event-triggered IUKF is presented for celestial navigation using time delay measurement [30]. However, the efficiency of the aforementioned event-triggered mechanism is influenced by the constant threshold. If the threshold is not set properly, a serious decrease in navigation accuracy or less computation load decrease will occur. This forms the incentive for our work.

To sum up, this paper proposes a parameter-independent event-triggered implicit UKF and implements it for celestial navigation using time delay measurement. Compared with the existing works, the main contributions have the following two aspects:

(1)

The dynamic threshold related to previous moments is substituted for the constant threshold. By comparing the innovation at the current moment and the updated estimate covariance at the last moment with the previous value, the event is automatically triggered when the accuracy of the state estimate is low. Different from the traditional event-triggered mechanisms, the performance of the proposed mechanism is not affected by any parameter.

(2)

Considering that large measurement errors will lead to large innovation, we introduce the updated estimate covariance at the previous time in the event-triggered condition. The parameter-independent event-triggered mechanism considering both innovation and updated estimate covariance can get higher navigation accuracy and less running time than the parameter-independent event-triggered mechanism only considering innovation, which will be verified by the simulation given below.

The rest of this paper is organized as follows: Section 2 shows the basic principle of celestial navigation using time delay measurement. The parameter-dependent event-triggered mechanism and the proposed parameter-independent event-triggered IUKF are introduced in Section 3. Section 4 compares the simulation results of the proposed method and other existing methods. The conclusions are given in Section 5.

2. Celestial Navigation Using Time Delay Measurement

It has been known that solar oscillations occur frequently, which causes dramatic changes in the solar spectral wavelength. This can be regarded as a feature in acquiring the difference between the arrival time of direct sunlight and that of reflected sunlight. The time delay of the reflected sunlight is related to the relative position of the spacecraft, the reflecting celestial body, and the Sun. Therefore, it can be adopted as a navigation measurement to provide the distance information of the spacecraft relative to the target celestial body.

2.1. State Model

The state model of the navigation system is constructed based on the orbital dynamics of the Mars probe. When the Mars probe is in the Mars approach stage, its motion can be described as a perturbed three-body model with Mars as the central body. The state model in the Mars-centered inertial coordinate system (J2000.0) can be written as follows:

$$\left\{\begin{array}{l}\dot{\mathit{r}}=\mathit{v}\hfill \\ \dot{\mathit{v}}=-{\mu }_{\mathrm{m}}\frac{\mathit{r}}{r^3}-{\mu }_{\mathrm{s}}\left[\frac{{\mathit{r}}_{\mathrm{ts}}}{r_{\mathrm{ts}}^3}+\frac{{\mathit{r}}_{\mathrm{sm}}}{r_{\mathrm{sm}}^3}\right]+{\mathit{w}}_{\mathrm{v}}\hfill \end{array},\right.$$

(1)

where ${\mu }_{\mathrm{m}}$ and ${\mu }_{\mathrm{s}}$ are the gravitational constants of Mars and the Sun, respectively; $r$ and ${\mathit{r}}_{\mathrm{ts}}$ are the position vectors of the spacecraft with respect to Mars and the Sun, respectively; $r$ and $r_{\mathrm{ts}}$ are the magnitude of $r$ and ${\mathit{r}}_{\mathrm{ts}}$, respectively; $\mathit{v}$ is the velocity vector of the spacecraft with respect to Mars; ${\mathit{r}}_{\mathrm{sm}}$ is the position vector of the Sun with respect to Mars; $r_{\mathrm{sm}}$ is the magnitude of ${\mathit{r}}_{\mathrm{sm}}$; and ${\mathit{w}}_{\mathrm{v}}$ is the process noise that comes from miscellaneous perturbations. Equation (1) can be shown as follows:

$$\dot{\mathit{X}}\left(t\right)=f\left(\mathit{X}\left(t\right),t\right)+\mathit{W}\left(t\right),$$

(2)

where $\mathit{X}={\left[\mathit{r},\mathit{v}\right]}^{\mathrm{T}}$ and $\mathit{W}={\left[0,{\mathit{w}}_{\mathrm{v}}\right]}^{\mathrm{T}}$.

2.2. Measurement Model

Mars is surrounded by an atmosphere, which changes the speed and path of sunlight. Thus, Phobos is adopted as the reflective celestial body for the Mars probe. Two atomic resonance spectrometers are employed simultaneously to record the solar spectral wavelength. One is aimed at the Sun to observe the direct sunlight, and the other is aimed at Phobos to observe the reflected sunlight. By comparing the wavelength features of the sunlight propagating along different paths, the time delay measurement can be acquired:

$$\mathit{Z}=\left[\Delta t\right]=\left[t_2-t_1\right],$$

(3)

where $t_2$ is the feature’s moment of the reflected sunlight, and $t_1$ is the feature’s moment of the direct sunlight.

The basic principle of celestial navigation using time delay measurement is given in Figure 1. Solar photons eject at $t_0$. Some of them transmit along path 1 and are captured by the spectrometer at $t_1$. Some other photons transmit along path 2 and reach the reflective celestial body at $t_r$. The position and velocity of Phobos related to the Sun are ${\mathit{r}}_{\mathrm{psr}}$ and ${\mathit{v}}_{\mathrm{psr}}$, respectively. They are captured by the spectrometer at $t_2$. The position vectors of the Mars probe with respect to the Sun are ${\mathit{r}}_{\mathrm{ts}1}$ and ${\mathit{r}}_{\mathrm{ts}2}$ when the time is $t_1$ and $t_2$. The position vectors of the Mars probe with respect to Mars are ${\mathit{r}}_1$ and ${\mathit{r}}_2$ when the time is $t_1$ and $t_2$. The velocity vectors of the Mars probe with respect to the Sun are $v_{\mathrm{ts}1}$ and $v_{\mathrm{ts}2}$ when the time is $t_1$ and $t_2$. The measurement model of time delay can be established based on the relative position relationship [13]

$$\Delta t=t_2-t_1=\frac{‖{\mathit{r}}_{\mathrm{psr}}‖+‖{\mathit{r}}_{\mathrm{ts}2}-{\mathit{r}}_{\mathrm{psr}}‖-‖{\mathit{r}}_{\mathrm{ts}1}‖}{c},$$

(4)

where $c$ is the speed of light and $‖ \cdot ‖$ represents the norm operator. Based on vector relationships, Equation (4) can be written as follows:

$$\Delta t=\frac{‖{\mathit{r}}_{\mathrm{pmr}}-{\mathit{r}}_{\mathrm{smr}}‖+‖{\mathit{r}}_2-{\mathit{r}}_{\mathrm{sm}2}-{\mathit{r}}_{\mathrm{pmr}}+{\mathit{r}}_{\mathrm{smr}}‖-‖{\mathit{r}}_1-{\mathit{r}}_{\mathrm{sm}1}‖}{c},$$

(5)

where ${\mathit{r}}_{\mathrm{pmr}}$ and ${\mathit{r}}_{\mathrm{smr}}$ are the position vectors of Phobos and the Sun with respect to Mars at $t_r$, respectively; ${\mathit{r}}_{\mathrm{sm}1}$ and ${\mathit{r}}_{\mathrm{sm}2}$ are the position vectors of the Sun with respect to Mars at $t_1$ and $t_2$, respectively.

Because filtering is executed at $t_2$, ${\mathit{r}}_{\mathrm{ts}1}$ and ${\mathit{r}}_{\mathrm{psr}}$ have to be represented by ${\mathit{r}}_{\mathrm{ts}2}$. This process is given in Ref. [13]. The calculations of ${\mathit{r}}_{\mathrm{ts}1}$ and ${\mathit{r}}_{\mathrm{psr}}$ are related to the measurement $\Delta t$, and thus the measurement model is an implicit function. Equation (5) can be written as follows:

$$0=\frac{‖{\mathit{r}}_{\mathrm{pmr}}-{\mathit{r}}_{\mathrm{smr}}‖+‖{\mathit{r}}_2-{\mathit{r}}_{\mathrm{sm}2}-{\mathit{r}}_{\mathrm{pmr}}+{\mathit{r}}_{\mathrm{smr}}‖-‖{\mathit{r}}_1-{\mathit{r}}_{\mathrm{sm}1}‖}{c}-\Delta t,$$

(6)

Equation (6) can also be shown as follows:

$$0=h\left(\mathit{X}\left(t\right),\mathit{Z}\left(t\right)+\mathit{V}\left(t\right)\right),$$

(7)

where $0$ is the equivalent measurement and $\mathit{V}$ is the measurement error. The implicit unscented Kalman filter (IUKF) [16] is employed to achieve the optimal estimate.

3. Parameter-Independent Event-Triggered Implicit Unscented Kalman Filter

The calculation process of the celestial navigation system using time delay measurement is time-consuming: $t_r$ is obtained by solving the equation set with a numerical method; the equivalent measurement noise covariance needs to be worked out through unscented transformation (UT) [31]. Due to the limited computing resources on the Mars probe and high real-time requirements for navigation, it makes sense to reduce the unnecessary computational load.

3.1. Traditional Parameter-Dependent Event-Trigger Method

3.1.1. Parameter-Dependent Event-Triggered Mechanism Based on Measurement

The event-triggered mechanism intends to achieve an advisable compromise between the navigation real-time requirement and navigation accuracy. Commonly, the difference between the last released sensor measurement and the instantaneous sensor measurement is predefined as the update error. If the instantaneous sensor measurement is not much different from the last released sensor measurement, the predicted state estimate and the predicted estimate covariance are directly treated as output. IUKF is run only if the instantaneous sensor measurement is significantly different from the last released sensor measurement. The event-triggered condition is given as follows [24,27]:

$${\gamma }_k=\left\{\begin{array}{c}1, {\left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)}^{\mathrm{T}}\cdot \left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)>\delta \\ 0, {\left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)}^{\mathrm{T}}\cdot \left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)⩽\delta \end{array}\right.,$$

(8)

where $\tilde{\mathit{Z}}$ denotes the last released time delay measurement and $\delta \in {ℝ}^{+}$ is the constant threshold that needs to be set appropriately. $\delta$ can also be replaced with a function dependent on the last released sensor measurements, and the corresponding event-triggered condition can be set as follows [32,33]:

$${\gamma }_k=\left\{\begin{array}{c}1, {\left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)}^{\mathrm{T}}\cdot \left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)>\sigma {\tilde{\mathit{Z}}}^{\mathrm{T}}\cdot \tilde{\mathit{Z}}\\ 0, {\left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)}^{\mathrm{T}}\cdot \left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)⩽\sigma {\tilde{\mathit{Z}}}^{\mathrm{T}}\cdot \tilde{\mathit{Z}}\end{array}\right.,$$

(9)

where $\sigma \in {ℝ}^{+}$ denotes the threshold parameter.

3.1.2. Parameter-Dependent Event-Triggered Mechanism Based on Innovation

In the above cases, an event is triggered once the update error ${\left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)}^{\mathrm{T}}\cdot \left(\tilde{\mathit{Z}}-{\mathit{Z}}_k\right)$ exceeds the constant threshold $\delta$ or the released measurement-dependent threshold $\sigma {\tilde{\mathit{Z}}}^{\mathrm{T}}\cdot \tilde{\mathit{Z}}$, which essentially makes judgments based on changes in measurements. Another event-triggered condition makes judgments based on the accuracy of the state estimate. The innovation in the Kalman filter is the difference between the predicted measurement calculated by the predicted state estimation and the actual measurement, which reflects the information quantity of the measurement. Thus, the event-triggered condition can be designed as follows [25,34,35]:

$${\gamma }_k=\left\{\begin{array}{c}1, {\mathit{\upsilon }}_k^{\mathrm{T}}\cdot {\mathit{\upsilon }}_k>\delta \\ 0, {\mathit{\upsilon }}_k^{\mathrm{T}}\cdot {\mathit{\upsilon }}_k⩽\delta \end{array}\right.,$$

(10)

where ${\mathit{\upsilon }}_k={\mathit{Z}}_k-{\widehat{\mathit{Z}}}_{k\left|k-1\right.}$ is the innovation in the Kalman filter. ${\widehat{\mathit{Z}}}_{k\left|k-1\right.}$ is the predicted measurement obtained by predicted state estimate and measurement model. For the implicit measurement model, the innovation can be written as follows:

$${\mathit{\upsilon }}_k=h\left({\widehat{\mathit{X}}}_{k\left|k-1\right.},{\mathit{Z}}_k\right),$$

(11)

where ${\widehat{\mathit{X}}}_{k\left|k-1\right.}$ is the predicted state estimate. Then, the event-triggered condition for the implicit measurement model can be written as follows [30]:

$${\gamma }_k=\left\{\begin{array}{c}1, h^{\mathrm{T}}\left({\widehat{\mathit{X}}}_{k\left|k-1\right.},{\mathit{Z}}_k\right)\cdot h\left({\widehat{\mathit{X}}}_{k\left|k-1\right.},{\mathit{Z}}_k\right)>\delta \\ 0, h^{\mathrm{T}}\left({\widehat{\mathit{X}}}_{k\left|k-1\right.},{\mathit{Z}}_k\right)\cdot h\left({\widehat{\mathit{X}}}_{k\left|k-1\right.},{\mathit{Z}}_k\right)⩽\delta \end{array}\right.,$$

(12)

3.2. Parameter-Independent Event-Trigger Method

It can be seen from the above section that the threshold is either related to $\delta$ or $\sigma$. The setting of the parameter will affect the efficiency of the event-triggered mechanism. If the parameter is too large, the measurement update runs too few times, resulting in a serious decrease in navigation accuracy. When the parameter is set too small, the measurement update runs too many times. The computational load increased significantly without any significant improvement in navigation accuracy.

The parameter-independent event-triggered mechanism is expected to be set to automatically trigger an event when the accuracy of the state estimate is low. To get rid of the influence of the parameter value, the dynamic threshold is substituted for the constant threshold. The idea of a sliding window is applied. The innovation at the current moment is compared with the maximum value of the innovations at the previous $M$ moments. $M$ denotes the window size. When the innovation at the current moment exceeds the dynamic threshold, the measurement update needs to be run. The event-triggered condition can be set as follows:

$${\gamma }_k=\left\{\begin{array}{l}1, {\mathit{\upsilon }}_k^{\mathrm{T}}\cdot {\mathit{\upsilon }}_k>\max \left\{{\mathit{\upsilon }}_{k-i}^{\mathrm{T}}\cdot {\mathit{\upsilon }}_{k-i}\right\}\hfill \\ 0, \mathrm{Otherwise}\hfill \end{array}\right.,$$

(13)

where $i=1,2,\dots ,M$.

It is worth noting that large measurement errors will also cause innovation to increase, and the event should not be triggered at this time. Updated estimate covariance can reflect the deviation of the state estimate. When the updated estimate covariance at the previous time is small and the innovation at the next time is large, it can be considered that the large innovation is caused by the large measurement noise. Thus, the updated estimate covariance at the previous time can be introduced in the judgment to eliminate the influence of large measurement errors. The parameter-independent event-triggered condition can be set as follows:

$${\gamma }_k=\left\{\begin{array}{l}1, {\mathit{\upsilon }}_k^{\mathrm{T}}\cdot {\mathit{\upsilon }}_k>\max \left\{{\mathit{\upsilon }}_{k-i}^{\mathrm{T}}\cdot {\mathit{\upsilon }}_{k-i}\right\}\&P_{r,k-1}>\max \left\{P_{r,k-1-i}\right\}\hfill \\ 0, \mathrm{Otherwise}\hfill \end{array}\right.,$$

(14)

where $i=1,2,\dots ,M$. $P_{r,k-1}$ denotes the position estimate error calculated from the updated estimate covariance ${\mathit{P}}_{k-1}$:

$$P_{\mathit{r},k-1}=\sqrt{P_{k-1,11}+P_{k-1,22}+P_{k-1,33}}$$

(15)

$P_{k-1,11}$, $P_{k-1,22}$ and $P_{k-1,33}$ represent the first three elements on the main diagonal of ${\mathit{P}}_{k-1}$. When ${\mathit{\upsilon }}_k$ satisfies the condition in Equation (14) and $P_{r,k-1}$ does not, it may be caused by too large measurement errors, and the event is not triggered. When both conditions in Equation (14) are met, it signifies that the large innovation is caused by the large state estimate error, and an event is triggered.

3.3. Filtering Process of the Parameter-Independent Event-Triggered IUKF

Calculate Sigma Points

$$\begin{array}{c}{\chi }_{0,k-1}={\widehat{\mathit{X}}}_{k-1},{\omega }_0=\tau /\left(n+\tau \right)\\ {\chi }_{i,k-1}={\widehat{\mathit{X}}}_{k-1}+\sqrt{n+\tau }{\left(\sqrt{{\mathit{P}}_{k-1}}\right)}_i,{\omega }_i=1/\left[2\left(n+\tau \right)\right]\\ {\chi }_{i+n,k-1}={\widehat{\mathit{X}}}_{k-1}-\sqrt{n+\tau }{\left(\sqrt{{\mathit{P}}_{k-1}}\right)}_i,{\omega }_{i+n}=1/\left[2\left(n+\tau \right)\right]\end{array},$$

(16)

where $\tau \in \mathit{R}$. $n$ is the dimension of the state vector; ${\chi }_{i,k-1}$ is the $2n+1$ sigma points whose mean and covariance are ${\widehat{\mathit{X}}}_{k-1}$ and ${\mathit{P}}_{k-1}$, respectively; and ${\left(\sqrt{{\mathit{P}}_{k-1}}\right)}_i$ is the $i$th row of the matrix square root.

Time Update

$${\chi }_{i,k\left|k-1\right.}=f\left({\chi }_{i,k-1},k-1\right),$$

(17)

$${\widehat{\mathit{X}}}_{k\left|k-1\right.}={\displaystyle \sum _{i=0}^{2n}{\omega }_i{\chi }_{i,k\left|k-1\right.}},$$

(18)

$${\mathit{P}}_{k\left|k-1\right.}={\displaystyle \sum _{i=0}^{2n}{\omega }_i\left({\chi }_{i,k\left|k-1\right.}-{\widehat{\mathit{X}}}_{k\left|k-1\right.}\right){\left({\chi }_{i,k\left|k-1\right.}-{\widehat{\mathit{X}}}_{k\left|k-1\right.}\right)}^{\mathrm{T}}}+{\mathit{Q}}_k,$$

(19)

where ${\chi }_{i,k\left|k-1\right.}$ is the propagated sigma points; ${\mathit{P}}_{k\left|k-1\right.}$ is the predicted estimate covariance; and ${\mathit{Q}}_k$ is the covariance of the process noise.

Measurement Update

The sigma points of the predicted measurement are calculated as follows:

$${\mathit{z}}_{i,k\left|k-1\right.}=h\left({\chi }_{i,k\left|k-1\right.},{\mathit{Z}}_k\right)$$

(20)

The predicted measurement can be obtained as

$${\widehat{\mathit{Z}}}_{k\left|k-1\right.}={\displaystyle \sum _{i=0}^{2n}{\omega }_i{\mathit{z}}_{i,k\left|k-1\right.}}$$

(21)

The predicted error covariance of the measurement and the cross-covariance of the state and measurement are calculated as follows:

$${\mathit{P}}_{zz,k}={\displaystyle \sum _{i=0}^{2n}\begin{array}{c}{\omega }_i\left({\mathit{z}}_{i,k\left|k-1\right.}-{\widehat{\mathit{Z}}}_{k\left|k-1\right.}\right){\left({\mathit{z}}_{i,k\left|k-1\right.}-{\widehat{\mathit{Z}}}_{k\left|k-1\right.}\right)}^{\mathrm{T}}\end{array}}+{\mathit{S}}_k$$

(22)

$${\mathit{P}}_{xz,k}={\displaystyle \sum _{i=0}^{2n}\begin{array}{c}{\omega }_i\left({\chi }_{i,k\left|k-1\right.}-{\widehat{\mathit{X}}}_{k\left|k-1\right.}\right){\left({\mathit{z}}_{i,k\left|k-1\right.}-{\widehat{\mathit{Z}}}_{k\left|k-1\right.}\right)}^{\mathrm{T}}\end{array}}$$

(23)

where ${\mathit{S}}_k$ is the covariance of the equivalent measurement noise, whose detailed computation process is given in Ref. [16].

The updated state estimate and the updated estimate covariance can be acquired as follows:

$${\widehat{\mathit{X}}}_k=\left\{\begin{array}{l}{\widehat{\mathit{X}}}_{k\left|k-1\right.},\hspace{1em}\hspace{1em}\hspace{1em}\hspace{0.17em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{\gamma }_k=0\hfill \\ {\widehat{\mathit{X}}}_{\mathit{k}\left|k-1\right.}-K_k{\widehat{\mathit{Z}}}_{k\left|k-1\right.},\hspace{1em}\hspace{1em}\hspace{1em}{\gamma }_k=1 \parallel k\le M\hfill \end{array}\right.,$$

(24)

$${\mathit{P}}_k=\left\{\begin{array}{l}{\mathit{P}}_{k\left|k-1\right.},\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{0.17em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{\gamma }_k=0\hfill \\ {\mathit{P}}_{k\left|k-1\right.}-{\mathit{K}}_k{\mathit{P}}_{zz,k}{\mathit{K}}_k^{\mathrm{T}},\hspace{1em}\hspace{1em}\hspace{1em}{\gamma }_k=1 \parallel k\le M\hfill \end{array}\right.,$$

(25)

where ${\mathit{K}}_k={\mathit{P}}_{xz,k}{\mathit{P}}_{zz,k}^{-1}$ is the Kalman gain.

Remark

The main improvement of the parameter-independent event-triggered IUKF is the event-triggered condition. Its proof can refer to Ref. [30].

The process of the parameter-independent event-triggered IUKF is given in Figure 2 For the first $M$ filtering cycles, ${\widehat{\mathit{X}}}_k$ and ${\mathit{P}}_k$ are figured out by the time update and the measurement update. After $M$ filtering cycles, $P_{r,k-1}$ is compared with the maximum value at the last $M$ moments. If $P_{r,k-1}$ does not exceed the maximum value, the predicted state estimate and predicted estimate covariance are regarded as the updated value and input to the next filtering cycle. If $P_{r,k-1}$ exceeds the maximum value, ${\mathit{\upsilon }}_k$ is calculated from Equation (11) and compared with the maximum value at the last $M$ moments. If ${\mathit{\upsilon }}_k$ does not exceed the maximum value, the predicted state estimate and predicted estimate covariance are regarded as the updated value and input to the next filtering cycle. If ${\mathit{\upsilon }}_k$ exceeds the maximum value, the measurement update is run, and the updated state estimate and updated estimate covariance are input to the next filtering cycle.

4. Simulation Results

In this section, the navigation performance of the proposed method is compared with some existing methods through simulation to demonstrate its effectiveness and superiority.

4.1. Simulation Conditions

A desktop is employed in the simulation, whose processor is 3.50 GHz Intel Core i7-9700K and RAM is 16 GB. The Systems Tool Kit (STK) astrogator is adopted to generate the spacecraft’s trajectory data.

Table 1. Initial orbital parameters.

Parameter	Value
Launch date	20 July 2020
Arrival date	8 March 2021
C3 energy	$18.2287 {\mathrm{km}}^2/{\mathrm{s}}^2$
Right ascension of outgoing asymptote	6.80876°
Declination of outgoing asymptote	45.9439°

gives the initial orbital parameters.

The orbit data of planets are set up based on JPL DE421 [36], and the orbit data of Phobos are built based on SPICE ephemeris [37]. The standard deviation of the time delay measurement error is set as $1\times {10}^{-7}$ s. Other filter parameters are given in

Table 2. Filter Parameters.

Parameter	Value
Initial state errors	$\begin{array}{l}\delta {\mathit{X}}_0={\left[\delta x_1,\delta x_1,\delta x_1,\delta x_2,\delta x_2,\delta x_2\right]}^{\mathrm{T}}\\ \delta x_1=5 \mathrm{km},\delta x_2=0.1 \mathrm{m}/\mathrm{s}\end{array}$
Initial estimation error covariance	$\begin{array}{l}{\mathit{P}}_0=\mathrm{diag}{\left[p_1,p_1,p_1,p_2,p_2,p_2\right]}^{\mathrm{T}}\\ p_1=25 {\mathrm{km}}^2,p_2={\left(0.1 \mathrm{m}/\mathrm{s}\right)}^2\end{array}$
Covariance matrix of process noise	$\begin{array}{l}\mathit{Q}=\mathrm{diag}{\left[q_1,q_1,q_1,q_2,q_2,q_2\right]}^{\mathrm{T}}\\ q_1={10}^{-3} {\mathrm{m}}^2,q_2={10}^{-7} {\left(\mathrm{m}/\mathrm{s}\right)}^2\end{array}$

.

4.2. Results and Analysis

4.2.1. Results of the Celestial Navigation Using Time Delay

The navigation results of the celestial navigation using time delay measurement with different sampling periods are given in Figure 3 and

Table 3. Navigation results with different sampling periods.

Sampling Period (min)	Mean Position Error (km)	Mean Velocity Error (m/s)	Number of Measurement Update Runs	Running Time (s)
1	1.09	0.02	5760	123
5	1.59	0.02	1152	62
10	1.73	0.02	576	35
100	3.99	0.03	57	26
300	8.01	0.03	19	25

. The data given in Table III and the following parts are the navigation results of the second half of the simulation. It can be seen that the highest navigation accuracy is achieved when the filtering period is 1 min. The highest navigation accuracy also results in the longest run time. The position error and velocity error of the celestial navigation with a period of 10 min is similar to that of the celestial navigation with a period of 1 min, while the running time is reduced by about 72%. When the filtering period exceeds 100 min, the navigation error increases significantly and even diverges while the running speed is not significantly improved. It can be seen from the results that the filtering period affects the navigation performance. Too short a filtering period leads to excessive running time, while too long a filtering period leads to poor navigation accuracy.

4.2.2. Results of the Parameter-Dependent Event-Triggered Mechanism Based on Measurement

Figure 4 and

Table 4. Navigation results with different $\sigma$.

$\mathit{\sigma }$	Mean Position Error (km)	Mean Velocity Error (m/s)	Number of Measurement Update Runs	Running Time (s)
$\sigma =1\times {10}^{-4}$	1.13	0.02	4326	103
$\sigma =1\times {10}^{-3}$	1.38	0.02	949	42
$\sigma =1\times {10}^{-2}$	1.86	0.03	110	28
$\sigma =1\times {10}^{-1}$	19.59	0.08	15	24

show the navigation results of the parameter-dependent event-triggered mechanism based on measurement. Although the highest navigation accuracy can be obtained when $\sigma =1\times {10}^{-4}$, it is also the most time-consuming. When $\sigma =1\times {10}^{-1}$, the calculation takes the least time, but the navigation error diverges due to too few measurement update runs. When $\sigma =1\times {10}^{-2}$ or $\sigma =1\times {10}^{-3}$, the running speed is improved while maintaining high navigation accuracy. Thus, the parameter $\sigma$ has a great influence on the effect of the event-triggered mechanism. It is difficult and crucial to choose an accurate parameter in practice.

4.2.3. Results of the Parameter-Dependent Event-Triggered Mechanism Based on Innovation

The navigation results of the parameter-dependent event-triggered mechanism based on innovation with different $\delta$ are compared in Figure 5.

Table 5. Navigation results with different $\delta$.

$\mathit{\delta }$	Mean Position Error (km)	Mean Velocity Error (m/s)	Number of Measurement Update Runs	Running Time (s)
$\delta =1\times {10}^{-7}$	1.09	0.02	4812	117
$\delta =3\times {10}^{-7}$	1.10	0.02	2854	80
$\delta =5\times {10}^{-7}$	1.23	0.02	1010	50
$\delta =8\times {10}^{-7}$	2.37	0.04	59	33
$\delta =1\times {10}^{-6}$	6.38	0.04	22	31

gives the detailed navigation results and running time. Its navigation results are similar to those of the event-triggered mechanism based on measurement. The parameter $\delta$ has a great impact on the effect of the event-triggered mechanism. When $\delta =1\times {10}^{-7}$, many unnecessary measurement updates are run in the filtering process and the running time is too long. On the contrary, the measurement updates are run too few times and the navigation error is too large when $\delta =1\times {10}^{-6}$. When $\delta =5\times {10}^{-7}$, the running time is reduced while maintaining high navigation accuracy. It is worth noting that the mean position error of the celestial navigation with a period of 100 min is about 1.7 times larger than that of the event-triggered mechanism when $\delta =8\times {10}^{-7}$. This indicates that higher navigation accuracy can be achieved by the event-triggered mechanism when the number of measurement update runs is the same.

4.2.4. Results of the Parameter-Independent Event-Triggered Mechanism Only Considering Innovation

To verify the effect of introducing the updated estimate covariance in the event-triggered condition, we first investigate the navigation performance of the parameter-independent event-triggered mechanism only considering innovation (whose event-triggered condition is Equation (13)). Its navigation results are given in Figure 6 and

Table 6. Navigation results of the parameter-independent event-triggered mechanism only considering innovation, with different M.

M	Mean Position Error (km)	Mean Velocity Error (m/s)	Number of Measurement Update Runs	Running Time (s)
3	1.09	0.02	1438	55
5	1.13	0.02	972	46
10	1.29	0.02	528	40
20	1.50	0.02	301	35
30	2.27	0.02	217	33

. From Figure 6 it can be seen that the position errors and velocity errors with different $M$ are not much different, which demonstrates that the navigation accuracy of the parameter-independent event-triggered mechanism is also robust to the window size. Comparing Table 6 with Table 4 and Table 5, it can be seen that the navigation accuracy and running time can achieve a desirable compromise no matter what $M$ is taken. Different from the constant threshold in the traditional event-triggered mechanisms, the dynamic threshold related to previous moments is applied to automatically trigger an event when the innovation is large. When ${\rm M}=10$, the mean position error of the parameter-independent event-triggered mechanism is about 1.2 times larger than that of the celestial navigation with a period of 1 min, while the running time is reduced by about 70%.

4.2.5. Results of the Parameter-Independent Event-Triggered Mechanism

Then, we investigate the navigation performance of the parameter-independent event-triggered mechanism considering both innovation and updated estimate covariance (whose event-triggered condition is Equation (14)). The navigation results of the parameter-independent event-triggered mechanism with different $M$ are given in Figure 7 and

Table 7. Navigation results of the parameter-independent event-triggered mechanism with different M.

M	Mean Position Error (km)	Mean Velocity Error (m/s)	Number of Measurement Update Runs	Running Time (s)
3	0.88	0.02	1024	50
5	1.11	0.02	632	42
10	1.19	0.03	342	37
20	1.24	0.03	183	35
30	2.21	0.02	148	29

. Comparing Table 6 and Table 7, it can be seen that the parameter-independent event-triggered mechanism considering both innovation and updated estimate covariance can get higher navigation accuracy and shorter running time than the parameter-independent event-triggered mechanism only considering innovation when $M$ is the same. The time delay measurement error is random noise. When the updated estimate covariance at the previous time is small and the innovation at the next time is large, it can be considered that the large innovation is caused by the large measurement noise. The measurement update will not be triggered in this case after introducing the updated estimate covariance in the event-triggered condition. In other words, the measurement update is triggered only when large innovation is caused by large state estimation errors. By eliminating the measurement update trigger caused by large measurement noise, the running time is further reduced while ensuring navigation accuracy. Comparing Table 3 and Table 7, it can be seen that the parameter-independent event-triggered mechanism has a longer running time than the uniform periodic sampling mechanism at the same number of measurement update runs. The reason is that determining whether the event-triggered conditions are met takes a certain amount of running time.

5. Conclusions

This paper proposes a parameter-independent event-triggered implicit UKF for celestial navigation using time delay measurement. Parameter-independent event-triggered implicit UKF automatically triggers an event when the accuracy of the state estimate is low. To get rid of the influence of the parameter value, the innovation at the current moment and the updated estimate covariance at the last moment are compared with the previous value. Simulation results indicate that parameter value has a great effect on the navigation performance of the traditional parameter-dependent event-triggered mechanism, which is difficult and crucial to choose in practice. The proposed parameter-independent event-triggered implicit UKF can reduce the running time while maintaining high navigation accuracy. Meanwhile, its performance will not be affected by any parameter value or window size. It is notable that the following aspects deserve future study: (1) The proposed parameter-independent event-triggered mechanism is based on the system model of the celestial navigation using time delay measurement in this paper. Extending the proposed method to the general class of nonlinear systems with order n or strict-feedback systems requires future research. (2) A stability analysis of the parameter-independent event-triggered implicit UKF needs to be performed in the future. (3) It is worthy of follow-up research to verify the effectiveness of the proposed method through more scenarios.