# Autonomous Phase Control Combining EKF and Adaptive Neural Network for Remote Sensing Satellites

Authors: Xin Wang; Guogang Wang; Hongru Zheng
Journal: International Journal of Aerospace Engineering
 (2022)

---

## Abstract

In this paper, a novel, effective, and feasible autonomous phase control algorithm based on the extended Kalman filtering (EKF) and neural network is firstly developed for addressing the problems of configuration maintenance of remote sensing satellite constellation. A balanced moment arm optimization method is employed to design the installation structure layout of the chemical propulsion system in the satellite. On that basis, an autonomous orbit control strategy is presented for controlling the phase of satellites, and the EKF algorithm is utilized to determine the orbit used to calculate the satellite phase. A radial basis function (RBF) neural network-based attitude control method is proposed to solve the attitude disturbance problem in the course of the phase control, and the RBF neural network is utilized to approximate the coupling torque of the orbit control. The simulation results demonstrate the feasibility of the designed automatic phase control strategy of the satellite.


---

## Body

## 1. Introduction

With the rapid development of commercial earth-observation satellites, some application requirements such as rapid revisit and coverage have been improved, while a single satellite could not meet these requirements. Therefore, satellite networks like Jilin-1 constellation, Zhuhai-1 constellation, OneWeb constellation, and Aleph constellation have been established by many commercial satellite companies. At present, the orbit control of satellite networks is entirely dependent on control instructions sent from the ground stations, increasing the burden of daily ground operation and the ground station pressure [1]. The ground station fails can affect the processing capacity of the satellite network. This clarifies the importance of autonomous orbit control of satellite networks. Although continuous tracking of the satellite by the ground monitoring and control network is unnecessary, a few data transceiver stations should be established to complete the uplink and downlink data transmission. Therefore, autonomous orbit control provides various advantages, including human investment saving and the flexibility improvement of space tasks [2].There are lots of satellites in the orbit, such as earth-observation satellite [3, 4], global navigation satellite [5, 6], and communication satellite [7, 8]; the orbit control of these satellites includes orbit maintenance, phase adjusting, phase maintenance, position keeping, and local time maintenance of descending node. The preparatory work of orbit control includes orbit determination, orbit control strategy, orbit maneuver command processing, and operation execution. These frequent maneuvers may become burdensome for daily ground operations. Various studies have been devoted to autonomous orbit control. In the early days, the GPS receiver introduction leads to several theoretical works on autonomous orbit control of LEO satellites [9], in which its feasibility has been further verified by Plam et al. [10] and Gaias and Ardaens [11]. A ground-based autonomous orbit control strategy has been proposed by Christian et al. [12] and successfully applied to TerraSAR-X with a special computing software to perform orbit control tasks. Although the mentioned strategy successfully solved the burden of daily ground operations, there are limited resources in the ground station of an extensive satellite network to meet frequent satellite operations requirements. Thus, an onboard autonomous orbit control strategy has been presented by Bolandi and Abrehdari [13], Florio et al. [14], Garulli et al. [15], De Florio et al. [16], and Zhong and Gurfil [17]. The autonomous orbit control of low-orbit satellite includes orbit determination, orbit, and attitude control.In terms of orbit determination, least squares and Kalman filter methods are the most common methods. Since the Kalman filtering algorithm does not require historical observation data and a large amount of computation, it has been widely utilized in real-time processing systems. In Kalman filtering, both the system state and measurement equations should be linear, and the dynamic and observed noises should be zero-mean Gaussian white noise sequences. However, the orbit motion and GPS observation equations are highly nonlinear equations. Thus, the improved versions of the Kalman filter algorithm have been employed in the literature to determine the real-time orbit. Li and Xu [18] employed the onboard image of the regular-shaped object to design an orbit determination method. The 3D reconstruction method has been adopted to measure the object center’s visual axis vector, while the EKF algorithm has been utilized to estimate the orbit parameters. Zhai et al. [19] adopted Lagrangian multiplier techniques to develop a novel constrained Kalman filter framework for autonomous orbit determination. Jiang et al. [20] employed the reduction dynamic EKF to construct an improved GPS+BeiDou system (BDS) combination model and utilized the orbital dynamics model and EKF algorithm for smoothing the global navigation satellite system (GNSS) positioning results to improve the orbit accuracy. The Newton-Raphson method (NRM), combined with the EKF algorithm, named the NRM/EKF method, has been employed to present an orbit determination method for low-earth orbit satellites. Rough position values calculated in the NRM pre-processing step were provided as measurement inputs to the EKF algorithm to estimate the satellite’s position and velocity vector [21]. Sun et al. [22] employed the random model, measurement equation, and model structure of orbit motion for filtering design and proposed an enhanced EKF to deal with unknown significant measurement deviations. A real-time high-precision orbit determination has been successfully realized using the improved Kalman filter algorithm in the above studies. Thus, the Kalman filter algorithm is appropriate to determine the orbit in real-time.In terms of orbit and attitude control, due to a large deviation of the propulsion system, the satellite attitude will be greatly disturbed in the orbit control process, resulting in flywheel saturation. Therefore, it is necessary to adopt the propulsion system for attitude control. A thrusters-based attitude control method has been proposed for autonomous orbit control by Sato et al. [23], employing a low-thrust, 4-nozzle propulsion system and a pre-calculated thruster selection table for optimizing thrust and torque to achieve effective orbit control and maintain high-precision attitude stability for a long time. Gao et al. [24] combined autonomous orbit determination with autonomous control to present an integrated design method in which the neural network state observer has been designed to estimate the inter-satellite ranging information, and the controller has been constructed autonomously by another neural network based on the estimated state and the predefined correction trajectory. Bock and Tajmar [25] adopted a self-developed electric propulsion system to realize orbit and two-axis attitude control. The propulsion system can provide a large control torque suitable for the control system with a large disturbance torque. Yang and Liu [26] proposed a RBF neural network-based control algorithm to solve the unknown external disturbance problem for the satellite system. Cheng and Shu [27] combined the neural network with a genetic algorithm to estimate the unknown parameters of the satellite system. However, the neural network is a slow algorithm. Sliding mode control has been widely utilized in many fields due to its advantages, such as its simple algorithm, fast response, and strong robustness [28]. Asadi and Bagherzadeh [29] proposed a sliding mode control approach for nonlinear uncertain systems with external disturbances. The sliding mode control has been utilized to solve the spacecraft attitude stabilization problem in the presence of unknown external disturbances [30]. The mentioned problems have been solved through sliding mode control and neural networks. Therefore, the autonomous phase control using EKF and adaptive neural network is an effective method for low-orbit satellite constellation maintenance.Finally, the essential contributions of this paper can be summarized as follows:(1)
In order to get the reasonable torque in the three directions, the structural layout of the chemical propulsion system is optimized by using the balanced force arm optimization method(2)
The EKF algorithm is utilized to determine the orbit used to calculate the satellite phase, and an autonomous orbit control strategy is designed for controlling the satellite phase(3)
A RBF neural network-based attitude control method is presented to solve the attitude disturbance problem in the course of the phase control, and the RBF neural network is employed to approximate the coupling torque of the orbit control. Simulation results show the effectiveness of the proposed algorithm. Compared to the PD control method, the proposed control strategy eliminates the influence of coupling torque and improves the attitude control accuracyThis paper is organized as follows. The propulsion system’s structure layout and the model of the orbit and attitude dynamics are introduced in Section2. This is followed by designing an autonomous phase control algorithm in Section 3, including EKF design, autonomous orbit control algorithm, and attitude control algorithm. The autonomous orbit control strategy associated with a comparison with different control strategies is presented in Section 4. Finally, conclusions are provided in Section 5.


## 2. Dynamics Model

The satellite orbit control is subject not only to various perturbation forces but also to the thrust of the propulsion system. Similarly, the satellite attitude depends on the effect of the coupling torque. Hence, this section models the orbital dynamics, propulsion system layout, and attitude dynamics, respectively.


### 2.1. Orbit Dynamics Model

Since the satellite is in a synchronous solar orbit of 500~1000 km, the atmospheric influence should be considered. The orbital dynamics model [31] is obtained as
(1)M0‐MR¨=Fc+Fa+Fs+Fg.(2)M=Ft1Is,where R=xyzT is the position vector of the satellite, M0 is the initial mass of the satellite, M is the mass of fuel consumption, F is the thrust provided by each thruster, t1 is the working time of the orbit-controlled thruster, Is is vacuum specific impulse, Fa is the atmospheric perturbation force, Fs is the solar pressure, Fg is the gravity gradient, and Fc is the control force of the propulsion system.


### 2.2. Structure Layout of Propulsion System

The chemical propulsion system consists of five thrusters; one thruster is utilized for orbit control, while the others are employed for an attitude adjustment. The body coordinate system is denoted byO−XYZ. The attitude-controlled thruster is installed symmetrically on the X and Y axes of the satellite, while the geometric axis of the orbit-controlled thruster passes through the mass center. In order to make each attitude-controlled thruster of the propulsion system generate torque in three directions and better adjust the attitude of the satellite, an inclination angle is designed when the attitude-controlled thruster is installed. This optimal inclination angle is selected according to the method and the principles in the Appendix, and the inclination angle of the attitude-controlled thruster is -20°, as shown in Figure 1. The installation structure of the chemical propulsion system is successfully applied to an in-orbit satellite. According to the above structure model, the triaxial torque of ith thruster acting on the satellite coordinate system can be described as
(3)Tix=FLacosθTiy=F−Lbcosθ+LcsinθTiz=FLcsinθ,i=1,2,3,4,where La and Lb are, respectively, the distances from each thruster to the X and Y axes, Lc is the projected distance from each thruster to the mass center, F is the chemical propulsion thrust, and θ is the inclination angle of the attitude-controlled thruster.Figure 1
Diagram of propulsion system structure.


### 2.3. Attitude Dynamics Model

Defining attitude quaternionq=q1q2q3q4T. The kinematics and dynamics equations of the satellite can be written as
(4)q̇=12BwIẇ+w×Iw+τd=τ.whereB=q1q2q3q4T, q1=q4q3q2−q1T, q2=−q3q4q1−q2T, q3=q2−q1q4−q3T, q4=q1q2q3q4T, I is the moment of inertia, w is the angular velocity of the satellite in the inertial system, τ is the control torque of the satellite, and τd is the interference torque of the chemical propulsion system. Since it is difficult to obtain an accurate system deviation from the measuring parts, the adaptive controller is designed to compensate for the uncertainties, while the sliding mode control will be utilized to improve the stability and robustness of the control system.


## 2.1. Orbit Dynamics Model

Since the satellite is in a synchronous solar orbit of 500~1000 km, the atmospheric influence should be considered. The orbital dynamics model [31] is obtained as
(1)M0‐MR¨=Fc+Fa+Fs+Fg.(2)M=Ft1Is,where R=xyzT is the position vector of the satellite, M0 is the initial mass of the satellite, M is the mass of fuel consumption, F is the thrust provided by each thruster, t1 is the working time of the orbit-controlled thruster, Is is vacuum specific impulse, Fa is the atmospheric perturbation force, Fs is the solar pressure, Fg is the gravity gradient, and Fc is the control force of the propulsion system.


## 2.2. Structure Layout of Propulsion System

The chemical propulsion system consists of five thrusters; one thruster is utilized for orbit control, while the others are employed for an attitude adjustment. The body coordinate system is denoted byO−XYZ. The attitude-controlled thruster is installed symmetrically on the X and Y axes of the satellite, while the geometric axis of the orbit-controlled thruster passes through the mass center. In order to make each attitude-controlled thruster of the propulsion system generate torque in three directions and better adjust the attitude of the satellite, an inclination angle is designed when the attitude-controlled thruster is installed. This optimal inclination angle is selected according to the method and the principles in the Appendix, and the inclination angle of the attitude-controlled thruster is -20°, as shown in Figure 1. The installation structure of the chemical propulsion system is successfully applied to an in-orbit satellite. According to the above structure model, the triaxial torque of ith thruster acting on the satellite coordinate system can be described as
(3)Tix=FLacosθTiy=F−Lbcosθ+LcsinθTiz=FLcsinθ,i=1,2,3,4,where La and Lb are, respectively, the distances from each thruster to the X and Y axes, Lc is the projected distance from each thruster to the mass center, F is the chemical propulsion thrust, and θ is the inclination angle of the attitude-controlled thruster.Figure 1
Diagram of propulsion system structure.


## 2.3. Attitude Dynamics Model

Defining attitude quaternionq=q1q2q3q4T. The kinematics and dynamics equations of the satellite can be written as
(4)q̇=12BwIẇ+w×Iw+τd=τ.whereB=q1q2q3q4T, q1=q4q3q2−q1T, q2=−q3q4q1−q2T, q3=q2−q1q4−q3T, q4=q1q2q3q4T, I is the moment of inertia, w is the angular velocity of the satellite in the inertial system, τ is the control torque of the satellite, and τd is the interference torque of the chemical propulsion system. Since it is difficult to obtain an accurate system deviation from the measuring parts, the adaptive controller is designed to compensate for the uncertainties, while the sliding mode control will be utilized to improve the stability and robustness of the control system.


## 3. Autonomous Phase Control Algorithm Design

Autonomous phase control includes filter design, autonomous orbit control algorithm design, and attitude control design in the process of orbit control, as shown in Figure2.Figure 2
The flow chart of autonomous phase control.


### 3.1. EKF Design

Considering the central gravitation term, and higher order termJ2, J3, and J4, the gravity potential function of the earth is
(5)U=μR1−J2Re22R23sin2φ−1−J3Re32R35sin3φ−3sinφ−J4Re48R435sin4φ−30sin2φ+3+Wg,where Re is the earth radius, μ is the gravitational constant, φ is the geocentric latitude, φ=arcsin)z/R, R=R is the orbital position module, and Wg is the higher order uncertainties.The motion of the satellite can be formulated as(6)d2xdt2=−μxR31+3J2Re22R21−5z2R2+5J3Re32R33zR−7z3R3−5J4Re48R43−42z2R2+63z4R4+Wx,(7)d2ydt2=∂U∂y=d2xdt2yx+Wy,(8)d2zdt2=−μzR31+3J2Re22R23−5z2R2+5J3Re32R36zR−7z3R3−3R5z−5J4Re48R415−70z2R2+63z4R4+Wz,where WxWyWzT=∂Wg/∂x∂Wg/∂y∂Wg/∂zT is model uncertainties. Consider that the satellite velocity vector is denoted by V=vxvyvzT, while the state variable of the determined orbit state filter is represented by X=xyzvxvyvzT. Now, the state and observation equations of the satellite orbit can be obtained as
(9)Ẋ=FX,t+WtZ=HX,t+Vt,where W=000WxWyWzT, vt is the measurement noise, and HX,t and FX,t are the vector function, FX,t=F1F2F3F4F5F6T, which can be described as
(10)F1=vx,F2=vy,F3=vzF4=d2xdt2,F5=d2ydt2,F6=d2zdt2.After linearizing and discretizing Equation (9), we have
(11)Xk=AXk−1+WkZk=HkXk+Vk,where Wk and Vk are independent Gaussian white noise with zero mean, E=WkWjT=Qkδkj, and E=VkVjT=Rkδkj, in which Qk and Rk are the covariance of Gaussian white noises. Here, the matrix A is given by
(12)A=∂FX,t∂X=000100000010000001∂F4∂x∂F4∂y∂F4∂z000∂F5∂x∂F5∂y∂F5∂z000∂F6∂x∂F6∂y∂F6∂z000.The state transition matrix can be calculated as(13)Φk,k−1=I6×6+AT+A2T2/2+⋯+AkTk/k!+⋯where I6×6 is the identity matrix, and T is the sampling period.The observation matrixHk is described as
(14)Hk=∂HX,t∂X=∂HX,t∂R∂HX,t∂V.Then, the EKF equation can be expressed as(15)X^k/k−1=X^k−1+FX^k−1,tk−1TPk/k−1=Φk/k−1Pk−1Φk/k−1T+QkKk=Pk/k−1HkTHkPk/k−1HkT+Rk−1δX^k=KkZk−HX^k/k−1,kPk=I−KkHkPk/k−1I−KkHkT+KkRkKkTX^k=X^k/k−1+δX^k,where Pk is the state estimate covariance. Initialize X^0,P0,T,Qk,Rk, the observed data Zk is the system input. Recursion is carried out through formula (15) until the end. The algorithm description of EKF is shown as Figure 3. Equation (15) can be employed to determine the satellite position and velocity in real-time.Figure 3
The structure of EKF.


### 3.2. Autonomous Orbit Control Algorithm

For commercial Earth-observation constellation, the phase distribution maintenance is particularly essential for the revisit and coverage of satellites. If the altitude difference of two satellites in the same orbital plane is large, the phase will exceed a threshold within a few days. Therefore, two satellites should be kept in one orbital altitude as far as possible. The following assumptions are considered:Assumption 1.
A satellite can acquire orbit data from another satellite.Assumption 2.
The GPS navigation receiver data is usually available during orbit control.Assumption 3.
The deviation quaternion and angular velocity are continuous and bounded.The orbital dynamics equation can be converted to a perturbation equation with six elements. The average semi-major axis of the satelliteNA and the satellite NB are a1 and a2, respectively, with phase σ distributed in an orbital plane. Their phases should be kept within σ±δ, while the orbital period difference ΔT is given by
(16)ΔT2π=a13/μ−a23/μ,ifa1>a2a23/μ−a13/μ,ifa1<a2.where 2/R−V2/μ=1/a, R=sqrtx2+y2+z2, V=sqrtvx2+vy2+vz2, where X=xyzvxvyvzT are obtained by formula (15). The phase shift of the two satellites is determined as
(17)φ=ΔTT12×360×T,where T1=2πa13/μ.Autonomous logical judgment strategy is given by(18)Δa=−λa1−a2,a1>a2&φ>δ+λa1−a2,a1≤a2&φ>δ0,φ≤δ,where λ is the gain coefficient and δ is the phase limit of two satellites, The symbols + and – indicate the rise and fall of the orbit altitude of satellite NA.The semi-major axis change caused by orbital perturbation is obtained as(19)Δa1=2n1−e2eΔVSsinf+ΔVTcosf+ΔVTΔa2=2n1−e2ΔVT,where ΔVS and ΔVT are, respectively, the radial and lateral velocity increments provided by the orbit-controlled thruster, ΔVT is the lateral velocity increment provided by the attitude-controlled thruster, n is the orbital angular velocity, ee≈0 is the eccentricity, and f is the true anomaly.Then, the controller of the chemical propulsion system is given by(20)Δa=Δa1+Δa2,where Δa2=KaΔa1. Ka is the proportion coefficient of the attitude-controlled thruster, which should be calibrated after the satellite enters into the orbit.The impulse control model is applied to the chemical propulsion system. Thus, the control time is allocated to calculate the velocity increment, and the lateral velocity increment provided by the thruster is obtained as(21)ΔVT=Ft1M0ΔVt=Ft2M0,where t2=∑i=14t2i is the total working time of the attitude-controlled thruster acting on the lateral direction, and t2i is the working time of each thruster acting on the lateral direction. During the working time t2 of the propulsion system, there will be coupling torque τd in the satellite attitude.


### 3.3. Attitude Control Algorithm during Orbit Control

In the orbit control process, due to the thrust direction deviation, thruster installation deviation, mass center deviation, and other factors in the propulsion system, the orbit-controlled thrust will generate the attitude coupling torque to the satellite, leading to the satellite’s attitude instability. Now, an orbit control algorithm is designed for the satellite system (3) and orbit dynamics model (1)~(2) to ensure attitude stability.Consider thatqd is the desired attitude angle, and wd=w1dw2dw3dT is the desired angular velocity. Now, the deviation quaternion e1 and the deviation angular velocity e2 can be calculated as
(22)e1=q⊗qde2=w−wd.The following sliding mode surface is appropriately defined to tend the deviation quaternion and deviation angular velocity to zero(23)S=e2+K1e1=w−wr,where wr=wd−K1e1, when S⟶0, the deviation quaternion and deviation angular velocity tend to zero along the sliding mode surface.The derivative of the sliding surfaceS is calculated as
(24)Ṡ=ẇ−ẇr.In the orbit control process, the uncertain disturbance part of the dynamic model is approximated using a neural network. Defining(25)tdq,w,W=WTφ+ε,where W is the neural network weight value, φ is the neural network basis function, and ε is the approximation error. The neural network weight is defined as
(26)W=argminW^∈Usupq,w∈Dt^dq,w,W^−tdq,w,W,where U and D are the ranges of q, w, and W^, and t^dq,w,W^ is the estimated value of tdq,w,W. As shown in Figure 4, a neural network is employed to approximate the coupling torque τd of the satellite. It can be expressed as
(27)t^dq,w,W=W∧Tφ^+ε^.Figure 4
Diagram of the attitude control algorithm.The estimated error is given by(28)t~d=td−t^d=WTφ+ε−W∧Tφ^−ε^=WTφ−WTφ^+WTφ^−W∧Tφ^+ε−ε^=WTφ~+W~Tφ^+ε~,where W~=W−W^, and φ~=φ−φ^.The control law is defined as(29)τ=Iẇr+w×Iw+t^d−KdS−KssatS,where the saturation functionsatS can be expressed as
(30)satS=1Si>Δi−1Si<−ΔiSi/ΔiSi≤Δi,where S=S1S2S3S4S5S6T, and Δii=1,⋯6 are arbitrary known constants. The saturation function is employed to suppress the chattering phenomenon caused by repeated switching on the sliding mode surface during the convergence process.Equation (29) can be rewritten with the following form, which is similar to Equation (4)
(31)Iẇ−ẇr+td−t^d+KdS+KssatS=0.According to relations (24), (28), and (Proof), we have
(32)IṠ=−t~d−KdS−KssatSTheorem.
For satellite system (4), consider that the system interference and the satellite’s desired orbit and attitude are bounded. Now, if an adaptive attitude controller can be designed using Equations (29) and (30), then the tracking error is bounded, while the adaptation update law can be described as
(33)W^̇=−Γφ^STwhere Γ=diagΓiii=1,2⋯7, and Γ>0. When t⟶∞, and S⟶0, then q⟶qd and w⟶wd. This means that the tracking error is bounded.Proof.
Consider the following Lyapunov candidate function(34)V=12STIS+trW~TΓ−1W~.
The time derivative of Equation (34) is calculated as
(35)V̇=STIṠ+trW~TΓ−1W~̇=ST−t~d−KdS−KssatS+trW~TΓ−1W~̇=−STKdS−SKssatS−STWTφ~+W~Tφ^+ε~+tr−WTΓ−1W^̇=−STKdS−SKssatS−STWTφ~+ε~+tr−W~TΓ−1W^̇+φ^ST=−STKdS−SKssatS−STWTφ~+ε~.If matrixN6×6 is positive definite, WTφ~+ε~≥N6×6S. Thus, we have
(36)ST⋅WTφ~+ε~≥ST⋅N6×6S.The following inequalities are valid(37)ST⋅WTφ~+ε~≥STWTφ~+ε~ST⋅N6×6S≥STN6×6S.Then whenSTWTφ~+ε~≥ST⋅N6×6S, the following inequalities should be fulfilled
(38)STWTφ~+ε~≥STN6×6S.Thus, we have(39)V̇=−STKdS−SKssatS−STWTφ~+ε~≤−STKdS−SKssatS−STN6×6S≤0.Case 1.
IfS=0, then V̇=0;Case 2.
IfS≠0, then V̇<0. According to the Lyapunov stability theory, the satellite control system of the satellite is bounded and stable.SinceS is uniformly continuous, integrable, and bounded as t⟶∞, it is obtained by
(40)limt⟶∞s=0.Sincee1⟶0 and e2⟶0, according to Barbalat Lemma, the system is asymptotically convergent, besides being stable.


## 3.1. EKF Design

Considering the central gravitation term, and higher order termJ2, J3, and J4, the gravity potential function of the earth is
(5)U=μR1−J2Re22R23sin2φ−1−J3Re32R35sin3φ−3sinφ−J4Re48R435sin4φ−30sin2φ+3+Wg,where Re is the earth radius, μ is the gravitational constant, φ is the geocentric latitude, φ=arcsin)z/R, R=R is the orbital position module, and Wg is the higher order uncertainties.The motion of the satellite can be formulated as(6)d2xdt2=−μxR31+3J2Re22R21−5z2R2+5J3Re32R33zR−7z3R3−5J4Re48R43−42z2R2+63z4R4+Wx,(7)d2ydt2=∂U∂y=d2xdt2yx+Wy,(8)d2zdt2=−μzR31+3J2Re22R23−5z2R2+5J3Re32R36zR−7z3R3−3R5z−5J4Re48R415−70z2R2+63z4R4+Wz,where WxWyWzT=∂Wg/∂x∂Wg/∂y∂Wg/∂zT is model uncertainties. Consider that the satellite velocity vector is denoted by V=vxvyvzT, while the state variable of the determined orbit state filter is represented by X=xyzvxvyvzT. Now, the state and observation equations of the satellite orbit can be obtained as
(9)Ẋ=FX,t+WtZ=HX,t+Vt,where W=000WxWyWzT, vt is the measurement noise, and HX,t and FX,t are the vector function, FX,t=F1F2F3F4F5F6T, which can be described as
(10)F1=vx,F2=vy,F3=vzF4=d2xdt2,F5=d2ydt2,F6=d2zdt2.After linearizing and discretizing Equation (9), we have
(11)Xk=AXk−1+WkZk=HkXk+Vk,where Wk and Vk are independent Gaussian white noise with zero mean, E=WkWjT=Qkδkj, and E=VkVjT=Rkδkj, in which Qk and Rk are the covariance of Gaussian white noises. Here, the matrix A is given by
(12)A=∂FX,t∂X=000100000010000001∂F4∂x∂F4∂y∂F4∂z000∂F5∂x∂F5∂y∂F5∂z000∂F6∂x∂F6∂y∂F6∂z000.The state transition matrix can be calculated as(13)Φk,k−1=I6×6+AT+A2T2/2+⋯+AkTk/k!+⋯where I6×6 is the identity matrix, and T is the sampling period.The observation matrixHk is described as
(14)Hk=∂HX,t∂X=∂HX,t∂R∂HX,t∂V.Then, the EKF equation can be expressed as(15)X^k/k−1=X^k−1+FX^k−1,tk−1TPk/k−1=Φk/k−1Pk−1Φk/k−1T+QkKk=Pk/k−1HkTHkPk/k−1HkT+Rk−1δX^k=KkZk−HX^k/k−1,kPk=I−KkHkPk/k−1I−KkHkT+KkRkKkTX^k=X^k/k−1+δX^k,where Pk is the state estimate covariance. Initialize X^0,P0,T,Qk,Rk, the observed data Zk is the system input. Recursion is carried out through formula (15) until the end. The algorithm description of EKF is shown as Figure 3. Equation (15) can be employed to determine the satellite position and velocity in real-time.Figure 3
The structure of EKF.


## 3.2. Autonomous Orbit Control Algorithm

For commercial Earth-observation constellation, the phase distribution maintenance is particularly essential for the revisit and coverage of satellites. If the altitude difference of two satellites in the same orbital plane is large, the phase will exceed a threshold within a few days. Therefore, two satellites should be kept in one orbital altitude as far as possible. The following assumptions are considered:Assumption 1.
A satellite can acquire orbit data from another satellite.Assumption 2.
The GPS navigation receiver data is usually available during orbit control.Assumption 3.
The deviation quaternion and angular velocity are continuous and bounded.The orbital dynamics equation can be converted to a perturbation equation with six elements. The average semi-major axis of the satelliteNA and the satellite NB are a1 and a2, respectively, with phase σ distributed in an orbital plane. Their phases should be kept within σ±δ, while the orbital period difference ΔT is given by
(16)ΔT2π=a13/μ−a23/μ,ifa1>a2a23/μ−a13/μ,ifa1<a2.where 2/R−V2/μ=1/a, R=sqrtx2+y2+z2, V=sqrtvx2+vy2+vz2, where X=xyzvxvyvzT are obtained by formula (15). The phase shift of the two satellites is determined as
(17)φ=ΔTT12×360×T,where T1=2πa13/μ.Autonomous logical judgment strategy is given by(18)Δa=−λa1−a2,a1>a2&φ>δ+λa1−a2,a1≤a2&φ>δ0,φ≤δ,where λ is the gain coefficient and δ is the phase limit of two satellites, The symbols + and – indicate the rise and fall of the orbit altitude of satellite NA.The semi-major axis change caused by orbital perturbation is obtained as(19)Δa1=2n1−e2eΔVSsinf+ΔVTcosf+ΔVTΔa2=2n1−e2ΔVT,where ΔVS and ΔVT are, respectively, the radial and lateral velocity increments provided by the orbit-controlled thruster, ΔVT is the lateral velocity increment provided by the attitude-controlled thruster, n is the orbital angular velocity, ee≈0 is the eccentricity, and f is the true anomaly.Then, the controller of the chemical propulsion system is given by(20)Δa=Δa1+Δa2,where Δa2=KaΔa1. Ka is the proportion coefficient of the attitude-controlled thruster, which should be calibrated after the satellite enters into the orbit.The impulse control model is applied to the chemical propulsion system. Thus, the control time is allocated to calculate the velocity increment, and the lateral velocity increment provided by the thruster is obtained as(21)ΔVT=Ft1M0ΔVt=Ft2M0,where t2=∑i=14t2i is the total working time of the attitude-controlled thruster acting on the lateral direction, and t2i is the working time of each thruster acting on the lateral direction. During the working time t2 of the propulsion system, there will be coupling torque τd in the satellite attitude.


## 3.3. Attitude Control Algorithm during Orbit Control

In the orbit control process, due to the thrust direction deviation, thruster installation deviation, mass center deviation, and other factors in the propulsion system, the orbit-controlled thrust will generate the attitude coupling torque to the satellite, leading to the satellite’s attitude instability. Now, an orbit control algorithm is designed for the satellite system (3) and orbit dynamics model (1)~(2) to ensure attitude stability.Consider thatqd is the desired attitude angle, and wd=w1dw2dw3dT is the desired angular velocity. Now, the deviation quaternion e1 and the deviation angular velocity e2 can be calculated as
(22)e1=q⊗qde2=w−wd.The following sliding mode surface is appropriately defined to tend the deviation quaternion and deviation angular velocity to zero(23)S=e2+K1e1=w−wr,where wr=wd−K1e1, when S⟶0, the deviation quaternion and deviation angular velocity tend to zero along the sliding mode surface.The derivative of the sliding surfaceS is calculated as
(24)Ṡ=ẇ−ẇr.In the orbit control process, the uncertain disturbance part of the dynamic model is approximated using a neural network. Defining(25)tdq,w,W=WTφ+ε,where W is the neural network weight value, φ is the neural network basis function, and ε is the approximation error. The neural network weight is defined as
(26)W=argminW^∈Usupq,w∈Dt^dq,w,W^−tdq,w,W,where U and D are the ranges of q, w, and W^, and t^dq,w,W^ is the estimated value of tdq,w,W. As shown in Figure 4, a neural network is employed to approximate the coupling torque τd of the satellite. It can be expressed as
(27)t^dq,w,W=W∧Tφ^+ε^.Figure 4
Diagram of the attitude control algorithm.The estimated error is given by(28)t~d=td−t^d=WTφ+ε−W∧Tφ^−ε^=WTφ−WTφ^+WTφ^−W∧Tφ^+ε−ε^=WTφ~+W~Tφ^+ε~,where W~=W−W^, and φ~=φ−φ^.The control law is defined as(29)τ=Iẇr+w×Iw+t^d−KdS−KssatS,where the saturation functionsatS can be expressed as
(30)satS=1Si>Δi−1Si<−ΔiSi/ΔiSi≤Δi,where S=S1S2S3S4S5S6T, and Δii=1,⋯6 are arbitrary known constants. The saturation function is employed to suppress the chattering phenomenon caused by repeated switching on the sliding mode surface during the convergence process.Equation (29) can be rewritten with the following form, which is similar to Equation (4)
(31)Iẇ−ẇr+td−t^d+KdS+KssatS=0.According to relations (24), (28), and (Proof), we have
(32)IṠ=−t~d−KdS−KssatSTheorem.
For satellite system (4), consider that the system interference and the satellite’s desired orbit and attitude are bounded. Now, if an adaptive attitude controller can be designed using Equations (29) and (30), then the tracking error is bounded, while the adaptation update law can be described as
(33)W^̇=−Γφ^STwhere Γ=diagΓiii=1,2⋯7, and Γ>0. When t⟶∞, and S⟶0, then q⟶qd and w⟶wd. This means that the tracking error is bounded.Proof.
Consider the following Lyapunov candidate function(34)V=12STIS+trW~TΓ−1W~.
The time derivative of Equation (34) is calculated as
(35)V̇=STIṠ+trW~TΓ−1W~̇=ST−t~d−KdS−KssatS+trW~TΓ−1W~̇=−STKdS−SKssatS−STWTφ~+W~Tφ^+ε~+tr−WTΓ−1W^̇=−STKdS−SKssatS−STWTφ~+ε~+tr−W~TΓ−1W^̇+φ^ST=−STKdS−SKssatS−STWTφ~+ε~.If matrixN6×6 is positive definite, WTφ~+ε~≥N6×6S. Thus, we have
(36)ST⋅WTφ~+ε~≥ST⋅N6×6S.The following inequalities are valid(37)ST⋅WTφ~+ε~≥STWTφ~+ε~ST⋅N6×6S≥STN6×6S.Then whenSTWTφ~+ε~≥ST⋅N6×6S, the following inequalities should be fulfilled
(38)STWTφ~+ε~≥STN6×6S.Thus, we have(39)V̇=−STKdS−SKssatS−STWTφ~+ε~≤−STKdS−SKssatS−STN6×6S≤0.Case 1.
IfS=0, then V̇=0;Case 2.
IfS≠0, then V̇<0. According to the Lyapunov stability theory, the satellite control system of the satellite is bounded and stable.SinceS is uniformly continuous, integrable, and bounded as t⟶∞, it is obtained by
(40)limt⟶∞s=0.Sincee1⟶0 and e2⟶0, according to Barbalat Lemma, the system is asymptotically convergent, besides being stable.


## 4. Simulation Results and Discussion

In this section, a satellite with the chemical propulsion system of five thrusters is adopted to prove the validity of the proposed control method. The maximum control torque available is 0.5 Nm, where Table1 shows the external input deviation of the propulsion model. The initial orbit and the system parameters as inputs are shown in Table 2 and Table 3, respectively. The EKF parameters are presented in Table 4.Table 1
Deviation parameters of the propulsion system.
ParametersValuesGeometric installation deviationθa=0.5;0;0;0;0.5°Geometric installation deviation azimuthθac=90;90;90;90;90°Thrust deflectionθb=1;1;1;1;1°Thrust deflection azimuthθbc=90;90;90;90;90°Center of massImc=0.01;0.01;0.01mmTable 2
Initial orbit parameters.
ParametersInitial orbit parametersPosition ofx-1897.1 kmPosition ofy-6647.74 kmPosition ofz3.13363 kmVelocity ofvx-0.94481 km/sVelocity ofvy0.273174 km/sVelocity ofvz7.52935 km/sTable 3
Initial parameters of the system.
ParametersValueThrustF=1.3NProportion coefficientKa=1.3Satellite massM0=200 kgVacuum specific impulseIS=2000 N⋅s/kgMoment of inertiaI=diag555030Initial attitude quaternionq0=0.2;0.3;0.4Initial attitude angular velocityw0=‐0.2;0.1;0.3 (°/s)Desired attitude quaternionqd=0;0.8509;0Desired attitude angular velocityw0=0;0;0 (°/s)Table 4
EKF parameters.
ParametersValuePk0diag2003004000.30.50.3Qkdiag151515888×0.0001Rkdiag6006006000.060.060.06×10T0.125sTwo attitude control algorithms, including the proportional-derivative (PD) control (which has been implemented) and the neural networks, are designed for the propulsion system to improve the attitude control accuracy. The PD controller parameters are chosen asKp=diag27.52515 and Kd=diag137.512575, and the RBF neural network parameters are listed in Table 5.Table 5
RBF neural network control parameters.
ParametersValueK1diag0.630.630.63Kddiag88.778.878.8×0.23Ksdiag0.0060.0060.006Γdiag4444444×100A three-layer neural network is employed to approximate the coupling torque of orbit control, while the hidden layer contains seven neurons. The radial basis function of the neural network is chosen as(41)φj=e‐X‐cj/2bj2,j=1,2,⋯7,where cj and bj are the center and width of the RBF neural network.The inboard satellite orbit interfaces with actually random noises, making a large orbital instantaneous fluctuation. Thus, orbit filtering is required for the controller design. The random noises of orbit position and velocity in the system are in the range of [-30 30] and [-0.05 0.05], respectively. The designed EKF algorithm is employed to perform the orbit filtering of the satellite. The EKF effect is shown in Figure5.Figure 5
EKF algorithm effect. (a) The position error. (b) The velocity error.
(a)(b)The following conclusions can be derived from the simulation results:(1)
Figures5(a) and 5(b) show the error between the filtering position and the real position and the error between the filtering velocity and the real velocity, respectively. The solid, the colon, and the dotted lines in Figure 5(a) represent the position filtering errors of the X, Y, and Z axes, respectively. The solid, the colon, and the dotted lines in Figure 5(b) show the velocity filtering errors of the X, Y, and Z axes, respectively. Since the position data error of GPS on the satellite is about 30 m, which is not conducive to the strategy formulation of autonomous orbit control, the EKF algorithm is employed to filter the GPS data. The position noise is attenuated from the peak value of 30 m to within 0.6 m, while the velocity noise is attenuated from the peak value of 0.05 m/s to within 0.016 m/s. It can be seen that the designed EKF algorithm, which is successfully prepared for the formulation of an autonomous control strategy, can significantly improve the orbit determination accuracy(2)
Figures6(a) and 6(b) show the phase change and the average orbit semi-major of the satellite, respectively. After the orbit is determined, the orbit control strategy will be calculated. The phase threshold is adjusted in the range of [117.5°-122.5°], and when the satellite phase exceeds this threshold, the satellite’s propulsion system starts to work automatically. The orbit determination and phase tracking are necessary for autonomous phase maintenance. The average orbit semi-major is 6903.7 ± 0.2 km, the semi-major after applying the EKF algorithm is 6903.7 ± 0.02 km, and the orbital semi-major of the other satellite is 6902.7 ± 0.02 km. Therefore, the phase between the two satellites approaches at a certain speed (see Figure 6(a)), the phase between two satellites exceeds 117.5° at 705 s, the satellite autonomously develops the control strategy, the control time and duration are 705 s, and 62 s, respectively, and the control quantity of the semi-major is 1.04 km. In practice, more than 20 m should be in orbit to ensure that the satellite’s orbit altitude is lower than that of another satellite. Since the orbit after applying the EKF algorithm has an instantaneous fluctuation of ±20 m, the control quantity is increased by 40 m. It can be seen from Figure 6(b) that the satellite starts to control at 705 s, and stops control at 767 s, the orbit semi-major is reduced from 6903.7 ± 0.02 km to 6902.65 ± 0.02 km, and the semi-major control deviation is within 10 m, and then the phase will drift back. It can be concluded that the orbit control precision is very high, which is mainly related to the proportional coefficient of the attitude-controlled thruster. Inaccurate calibration of this coefficient will affect the orbit control accuracy. It proves the validity and feasibility of the whole process from satellite phase tracking to the autonomous determination of orbit control strategy, and then to orbit control. Moreover, autonomous phase maintenance is realized(3)
Figures7 and 8 show the attitude angle and quaternion deviation curves under two control strategies. Attitude adjustment and stability are successfully realized through both attitude control algorithms. As shown in Figure 5, the attitude adjustment is completed in 200 s, the phase exceeds the threshold at 705 s, and the chemical propulsion system begins to work to ensure the attitude stability. In Figures 8(a) and 8(b), the quaternion deviation under neural network sliding mode control is within -2 × 10-3 and 2 × 10-3 in the course of orbit control, while the quaternion deviation under PD control is within -4 × 10-3 and 4 × 10-3. Accordingly, the attitude control precision is significantly improved(4)
Figures9(a) and 9(b) show the control torques in the orbit control process. The control torques of the neural network and PD algorithms are stabilized within 0.15 Nm and 0.5 Nm, respectively. Output torques generated by the propulsion system are shown in Figures 10(a) and 10(b); the maximum torques of the neural network and PD control are 0.9 Nm, and 1.1 Nm, respectively. The designed control torque is smaller than the PD control torque. This demonstrates the superiority of the neural network control algorithmFigure 6
State change of the satellite. (a) The phase. (b) The orbit semi-major.
(a)(b)Figure 7
Attitude angle curve. (a) PD control. (b) RBF Neural network control.
(a)(b)Figure 8
Quaternion deviation curve. (a) PD control. (b) Neural network sliding mode control.
(a)(b)Figure 9
Control torque. (a) PD control. (b) Neural network sliding mode control.
(a)(b)Figure 10
Output torque of propulsion system. (a) PD control. (b) Neural network sliding mode control.
(a)(b)


## 5. Conclusions and Future Work

In this paper, a neural network-based autonomous phase control method is designed to solve the problems of configuration maintenance of remote sensing satellite constellation. Firstly, the balanced moment arm optimization method is utilized to develop the installation structure of the propulsion system to apply in the satellite. Secondly, the EKF algorithm is employed to determine the orbit used to calculate the satellite phase. Finally, an RBFNN-based attitude control method is proposed to solve the attitude disturbance problem in the course of the phase control. The effectiveness and feasibility of the proposed automatic phase control strategy of the satellite are verified through mathematical simulations.The proposed method in this work yields the theoretical command limiting the engineering application of orbit data interaction of small satellites. If the orbit data error occurs, the orbit determination accuracy will be significantly reduced, affecting the phase maintenance. Therefore, an adaptive autonomous orbit control method is necessary to solve this problem. Although the attitude instability problem in the orbit control process should be considered for the successful realization of the orbit control, there is still a long way to achieve autonomous orbit control in orbit. To attain this goal, the reliability improvement of the orbit data in orbit and adaptive autonomous orbit control methods should be studied in the future.



---

*Source: 7153667-2022-04-12.xml*