# Research on the Effectiveness of Different Estimation Algorithm on the Autonomous Orbit Determination of Lagrangian Navigation Constellation

Authors: Youtao Gao; Junkang Chen; Bo Xu; Jianhua Zhou
Journal: International Journal of Aerospace Engineering
 (2016)

---

## Abstract

The accuracy of autonomous orbit determination of Lagrangian navigation constellation will affect the navigation accuracy for the deep space probes. Because of the special dynamical characteristics of Lagrangian navigation satellite, the error caused by different estimation algorithm will cause totally different autonomous orbit determination accuracy. We apply the extended Kalman filter and the fading–memory filter to determinate the orbits of Lagrangian navigation satellites. The autonomous orbit determination errors are compared. The accuracy of autonomous orbit determination using fading-memory filter can improve 50% compared to the autonomous orbit determination accuracy using extended Kalman filter. We proposed an integrated Kalman fading filter to smooth the process of autonomous orbit determination and improve the accuracy of autonomous orbit determination. The square root extended Kalman filter is introduced to deal with the case of inaccurate initial error variance matrix. The simulations proved that the estimation method can affect the accuracy of autonomous orbit determination greatly.


---

## Body

## 1. Introduction

Deep space exploration has become a hot spot of aerospace. Several deep space probes have been launched. The autonomous navigation is important for deep space probes to deal with communication delay as well as reducing the dependency on ground stations. As early as 1968, the sextant had been used for autonomous navigation in “Apollo program” [1]. In 1999, “Deep Space 1” achieved autonomous orbit determination by tracking small celestial bodies with an optical sensor, which is the first successful on-orbit application of the deep space autonomous navigation technology [2]. The comet probe “Deep Impact” which was launched in 2005 also carried out its navigation and control automatically based on an optical navigation system with a high resolution imager [3]. In [4], Downs proposed to use X-ray pulsar radiation signal to navigate spacecraft. The rotation period of X-ray pulsar is extremely stable; therefore, time and the location of spacecraft can be determined by tracking several X-ray pulsars with given and fixed frequency [5]. The satellite navigation constellation can also provide navigation information for deep space probes. GPS can navigate the deep space probes when they are running in low-Earth orbits and medium Earth orbits. For deep space transfer orbits and deep space target orbits, the GPS is not good enough. Several researchers investigated weak GNSS signal navigation for the deep space probes [6–8]. Witternigg et al. introduced how GPS and Galileo could be used for orbit determination in future missions to the Moon [8]. Farquhar introduced a concept of using Earth-Moon libration point satellites for lunar navigation [9, 10]. In 2005, Hill suggested placing navigation constellation on the periodic orbits in the vicinity of libration points of the Earth-Moon system to support deep space navigation [11]. Zhang and Xu analyzed the architecture and navigation performance of the Lagrangian point satellite navigation system [12–14]. The Lagrangian navigation constellation is introduced to navigate the deep space probes autonomously. Hence, the navigation constellation itself should have the ability of autonomous orbit determination (AOD). The methods introduced in [1–5] can be considered as absolute navigation (or absolute autonomous orbit determination) because the estimated orbit refers to an inertial or quasi-inertial frame. Methods introduced in [6–14] can be classified as relative navigation. Relative navigation seeks optimal estimates for the position and velocity of one satellite relative to the other one.Relative navigation usually is applied to satellites in a formation or constellation involved using GPS which restricts the spacecraft formation to near-Earth applications, such as Deep Space Mission 3 [15] and Grace project [16]. Relative navigation is primarily proposed for formation configuration control and formation reconfiguration. However, Lagrangian navigation constellation should provide absolute navigation information to deep space probes to achieve absolute navigation. Therefore, autonomous orbit determination of the Lagrangian navigation constellation is actually using the relative measurement to achieve absolute navigation. The laser interferometer space antenna (LISA) mission is an example which uses the relative range to assist the absolute orbit determination. The LISA mission which consists of three spacecraft separated by 5 million kilometers forming an equilateral triangle is a huge Michelson interferometer in space for gravitational wave detection [17]. The deep space-network provides a raw estimation of the absolute positioning for the three satellites. An accurate relative positioning will be provided with a laser-based ranging measurement in order to obtain an accuracy of positioning of tens of meters [18, 19]. Psiaki [20] and Markley [21] suggested using crosslink range, attitude information, and an optical tracker to determine the orbits autonomously. Yim et al. proposed using optical tracking and attitude information to find the direction vector between the two spacecraft and determine both orbits [22]. But these methods require extensive hardware development. In order to reduce the operational cost, size, and weight of spacecraft for formation missions, crosslink range can be used as the only measurement for orbit determination of a constellation. However, for the Earth navigation satellite constellation, there is a rank deficiency problem when only crosslink range is used to determine the orbit [23, 24]. Hill’s study illustrated that the rank deficiency problem does not exist for the Lagrangian navigation satellites because of the special dynamics near the libration points [11]. Thus, the Lagrangian navigation satellites can autonomously determine their orbits using only crosslink range. In [11], Hill discussed the unique distribution of Lagrangian orbit from the view of dynamics, which theoretically proved the autonomy of the Lagrangian navigation constellation. From the perspective of identifiability of epoch state, Qian et al. verified the feasibility of AOD for satellites in quasiperiodic orbits about the Earth-Moon libration point [25]. Based on circular restricted three-body problem (CR3BP), Du et al. researched the autonomous orbit determination method of satellites in halo orbits, and only crosslink range was used as observation [26]. In [23], Gao et al. discussed the feasibility of autonomous orbit determination using only the crosslink range measurement for a combined Lagrangian navigation constellation and GNSS. The most widely used algorithm for Lagrangian navigation satellite is extended Kalman filter (EKF) method. However, EKF is based on linearization of the system dynamics and the assumption of Gaussian process/measurement noise. These can seriously affect the performance of the state estimation and even lead to divergence. Unscented Kalman filter (UKF) which is based on the unscented transform can achieve higher accuracy than EKF while the added computational cost is not significant. More importantly, UKF is robust with respect to the initial conditions. Therefore, Sun et al. introduced UKF to relative navigation for multiple spacecraft formation flying [27]. Giannitrapani et al. analyze the performance of EKF and UKF for the localization of a spacecraft [28]. In order to improve the robustness and stability accuracy, Wang and Gu applied fault tolerant UKF in autonomous determination of relative orbit for satellite formation flying [29]. Rigatos introduced the technical analysis and implementation cost assessment of sigma-point Kalman filtering and particle filtering in autonomous navigation systems [30]. Reali and Palmerini provided a preliminary comparison of different estimation techniques to be used in formation flying navigation [31].Since the CR3BP is sensitive to the state error and calculation error, the AOD of Lagrangian navigation satellite may refer to the accuracy of the estimation algorithms. One factor which must be considered in the AOD of Lagrangian navigation constellation is to prevent the divergence of AOD error. Therefore, we introduce four estimation methods to achieve the AOD of Lagrangian navigation constellation to analyze the effect on AOD by estimation method.


## 2. Dynamical Model of Lagrangian Navigation Satellite

For satellites in Lagrangian point orbits, the equation of CR3BP should be an appropriate model to describe the satellites’ dynamical characteristics. Consider two massive bodiesm1 and m2 moving under the action of just their mutual gravitation, and let their orbit around each other be a circle of radius r12. As shown in Figure 1, a noninertial, comoving frame of referenceo-xyz is defined. The origin of frameo-xyz lies at the center of mass of the two-body system. The positive x direction goes from m1 to m2. The positive y-axis is parallel to the velocity vector. The z-axis is perpendicular to the orbital plane. Now the third body of mass m which is vanishingly small compared to the primary masses m1 and m2 is introduced. We assume that the mass m is so small that it has no effect on the motion of the primary bodies. This is called the restricted three-body problem.Figure 1
Diagram of the circular restricted three-body problem ino-xyz frame.The nondimensional equations of motion for the CR3BP are shown as follows [32]:(1)x¨-2y˙=x-1-μx+μr13-μx+μ-1r23,y¨+2x˙=1-1-μr13-μr23y,z¨=μ-1r13-μr23z,where m1 and m2 are the masses of the two primaries. r1 and r2 are the distances of the spacecraft from the massive and secondary primaries, respectively. r1, r2, and μ are given by (2)r1=x+μ2+y2+z2,r2=x+μ-12+y2+z2,μ=m2m1+m2.The nondimensional scales are defined as(3)M=m1+m2,L=r12,T=r123Gm1+m21/2,where r12 is the distance between the two primaries. G is the Universal Constant of Gravitation.Equations (1) have five equilibrium points which are called Lagrange points (or libration point), as shown in Figure 2.Figure 2
Schematic location of the five Lagrange points in the CR3BP.Three Lagrange points on thex-axis are unstable, and two Lagrange points which form an equilateral triangle with the two primary bodies in thex-y plane are stable. Many interesting periodic orbits exist in the vicinity of the Lagrange points no matter whether they are stable. The Lagrange navigation satellite constellation discussed in this paper is distributed in these periodic orbits.


## 3. AOD of Lagrangian Navigation Constellation with Different Estimation Methods

EKF method is a commonly used method for AOD of satellite. EKF is summarized as follows [33]:(4)X^k=Φk,k-1X^k-1+Kkρk-HkΦk,k-1X^k-1,Kk=PkHkTRk-1,Pk=I-KkHkPk,k-1,Pk,k-1=Φk,k-1Pk-1Φk,k-1T+Γk,k-1Qk-1Γk,k-1T,where X^k-1 is the estimated state at tk-1. Φtk,tk-1 is the transition matrix. Pk is the covariance matrix. Hk is the mapping matrix that relates the observation deviation vector to the state deviation vector. Kk is the Kalman gain matrix. In this paper, only crosslink range between two satellites is used as observation. The crosslink range is described as(5)ρi,j=xi-xj2+yi-yj2+zi-zj2.So (6)Hk=∂ρxk,k∂xkTxk,k-1.Two periodic orbits around L1 and L2, respectively, are chosen to achieve the simulation. The initial states of the two orbits are shown in Table 1. The values are nondimensional. 10 m initial state errors and 10 m measurement noise are added. In Figures 3 and 4, the orbits based on different initial state is shown. As we can see, the trajectory will diverge within a short time if the initial state error is 10 m. Therefore, the initial state error and the measurement noise are chosen based on the navigation requirement, not the ability of the existing technology.Table 1
Initial states of the periodic orbits around Lagrange pointsL1 and L2.
Satellites L 1 L 2 x ( L ) 0.8440212401521470 1.1726789595745000 y ( L ) 0.0000000000000000 0.0000000000000000 z ( L ) 0.0592695845762629 0.0834298981288340 x ˙ ( L / T ) 0.00000000000000000 0.0000000000000000 y ˙ ( L / T ) - 0.0053009560909076 - 0.1864489128036080 z ˙ ( L / T ) 0.0000000000000000 0.0000000000000000Figure 3
The periodic orbit based on accurate initial state.Figure 4
The orbit based on the initial state with 10 m error.The initial covariance matrix is defined as(7)Pt0=σx02000000σy02000000σz02000000σx˙02000000σy˙02000000σz˙02,where σx0, σy0, and σz0 are the initial position error in three coordinate axes and σx˙0, σy˙0, and σz˙0 are the initial velocity error in three coordinate axes.Figures5 and 6 show the results of AOD of Lagrangian navigation satellite using EKF. As we can see, the maximum errors in three axes of L1 are 541 m, 254 m and 168 m. The maximum errors in three axes of L2 are 555 m, 195 m, and 263 m. And the estimation process is divergent.Figure 5
The AOD position error ofL1 using EKF.Figure 6
The AOD position error ofL2 using EKF.Since CR3BP is sensitive to calculation error, the state estimation error covariance matrix of EKF may lose its properties of positive definite when the computations are carried out with finite digit arithmetic of the computer. Nonpositive definite covariance matrix will cause the gain matrix to lose its function gradually and lead to divergence of AOD error. In order to suppress the divergence of estimation process, the fading-memory filter (FMF) is introduced to improve the AOD characteristics. The FMF is essentially the old estimate plus a gain times a residual (difference between current measurement and previous estimate). The difference between the EKF and FMF is that the FMF is recursive and weights new measurements more heavily than older measurements [34].The basic equations of FMF are listed as follows [33]:(8)X^k=Φk,k-1X^k-1+Kkρk-HkΦk,k-1X^k-1,Kk=PkHkTRk-1,Pk=I-KkHkPk,k-1,Pk,k-1=Φk,k-1Pk-1Φk,k-1TeCk-1+Γk,k-1Qk-1Γk,k-1T.FMF raises the role of new observation data and reduces the influence of old observation to the filter in order to suppress the filtering divergence. Figures 7 and 8 show the AOD position errors of L1 and L2 using FMF. It can be seen that the divergence of errors is suppressed by using FMF method. The accuracy of the AOD is improved as well.Figure 7
The AOD position error ofL1 using FMF.Figure 8
The AOD position error ofL2 using FMF.One disadvantage of FMF is that the oscillation amplitude of filtering error at early stage is larger. Comparing Figure5 with Figure 7, we can see that estimation process of EKF is relatively ideal during the first 100 s. Therefore, we can use EKF to determine the orbit of Lagrangian navigation satellite at the early stage of the AOD. When the AOD error is bigger than a threshold value, the estimation method will be switching to FMF. We call this method integrated Kalman fading filter (IKFF) which is illustrated as (9)X^k=Φk,k-1X^k-1+Kkρk-HkΦk,k-1X^k-1,Kk=PkHkTRk-1,Pk=I-KkHkPk,k-1,Pk,k-1=Φk,k-1Pk-1Φk,k-1T+Γk,k-1Qk-1Γk,k-1T,ife≥es,Φk,k-1Pk-1Φk,k-1TeCk-1+Γk,k-1Qk-1Γk,k-1T,ife<es,where es is the threshold value of the AOD error. Since the estimation error cannot be got in practical applications, the residual can be used to approximate the threshold value of estimation error. The linearized relationship between the residual and the AOD error can be expressed as(10)y=Ce,where y is the residual and C is denoted as(11)C=x^-x^2ρy^-y^2ρz^-z^2ρ000.x^y^z^T and x^2y^2z^2T are the estimated position variable of Lagrangian satellite 1 and Lagrangian satellite 2, respectively. ρ is the crosslink range between these two satellites.The AOD results of IKFF are shown in Figures9 and 10. We list the maximum error of the above three methods in Tables 2 and 3. As shown in Table 2, the maximum error of L1 satellite in three axes using IKFF reduced to 204 m, 104 m, and 62 m. Meanwhile, oscillation characteristic at the early stage of AOD is improved in Figure 9. The similar improvement of AOD can be seen about L2 satellite. Therefore, as long as the threshold value is selected appropriately, IKFF will produce better AOD results than FMF or EKF.Table 2
The maximum AOD error ofL1 using different methods.
Maximum error inx-axis/m Maximum error iny-axis/m Maximum error inz-axis/m EKF 541 253 168 FMF 273 129 92 IKFF 204 104 62Table 3
The maximum AOD error ofL2 using different methods.
Maximum error inx-axis/m Maximum error iny-axis/m Maximum error inz-axis/m EKF 554 195 262 FMF 296 226 163 IKFF 215 81 128Figure 9
The AOD position error ofL1 using IKFF.Figure 10
The AOD position error ofL2 using IKFF.To overcome the divergence of AOD error caused by the state estimation error covariance matrix losing its properties of positive definite symmetric, we also introduce the square root extended Kalman filter (SR-EKF) to estimate the state of the Lagrangian satellite using only crosslink range. The SR-EKF is described as follows [33]:(12)X^k=Φk,k-1X^k-1+Kkρk-HkΦk,k-1X^k-1,X^0=EX0,Kk=akSk,k-1Fk,Fk=Sk,k-1THkT,ak=FkTFk+Rk-1,Sk,k-1=Φk,k-1Sk-1,P0=S0S0T,Sk=Sk,k-1I-akrkFkFkT,rk=11±akRk.The state estimation error covariance matrix P is replaced by SST where S is the square root of P [35]. Here we define the order of σx02, σy02, and σz02 as P0. This method can guarantee the P is symmetric nonnegative definite matrix at any time. In Tables 4 and 5, we give the maximum error of AOD of L1 satellite and L2 satellite with different initial state estimation error covariance matrix. As we can see, with a small initial value of error variance matrix, SR-EKF algorithm cannot estimate the state very well. The estimation error using SR-EKF is bigger than the EKF results. However, the accuracy of AOD is improved with the increase of initial error covariance matrix. Thus, the effectiveness of SR-EKF method is affected by the initial value of state estimation error covariance matrix P. SR-EKF is useful when the initial error variance matrix is not accurate. However, the AOD error will increase when the initial error variance matrix is bigger than a boundary. In our simulation, the best AOD result coincides with the P0 = 10−10. Figures 11 and 12 exhibit the AOD error using SR-EKF under the same initial condition as the above three methods with P0 = 10−10. Even when the initial error variance matrix is not accurate, the precision of AOD of two Lagrangian satellites is better than EKF.Table 4
The maximum AOD error ofL1 using SR-EKF with different P0.
Maximum error inx-axis/m Maximum error iny-axis/m Maximum error inz-axis/m P 0 = 1 0 - 16 17936 8507 4378 P 0 = 1 0 - 14 853 393 184 P 0 = 1 0 - 12 404 186 88 P 0 = 1 0 - 10 377 178 104Table 5
The maximum AOD error ofL2 using SR-EKF with different P0.
Maximum error inx-axis/m Maximum error iny-axis/m Maximum error inz-axis/m P 0 = 1 0 - 16 17290 5216 6532 P 0 = 1 0 - 14 790 283 301 P 0 = 1 0 - 12 374 134 154 P 0 = 1 0 - 10 365 107 134Figure 11
The AOD position error ofL1 using SR-EKF with P0=10-10.Figure 12
The AOD position error ofL2 using SR-EKF with P0=10-10.


## 4. Conclusion

Since the special dynamical characteristics of CR3BP, the AOD accuracy of the Lagrangian navigation constellation is sensitive to all the errors introduced by the AOD process. Different estimation methods will cause great different AOD accuracy.Four estimation methods are used to determinate the orbit of Lagrangian navigation constellation. EKF is a commonly used method, but it is not the best choice for the AOD of Lagrangian navigation constellation. The accuracy of AOD using FMF can improve 50% more than using EKF. And the divergence of the AOD error is suppressed.We proposed a new method, IKFF, to estimate the state of Lagrangian navigation satellite. The accuracy of AOD using IKFF can improve 20% more than using FMF. Furthermore, IKFF can smooth the process of the AOD.For the case of not accurate initial error variance matrix, SR-EKF can be chosen as the AOD estimation method. The accuracy of SR-EKF is close to the FMF.In order to improve the AOD accuracy, all the factors which will cause errors should be analyzed. As the estimation methods, we will find other advanced technologies to get a much more accurate AOD of Lagrangian navigation satellite.



---

*Source: 8392148-2016-10-25.xml*