Estimation of Pico-Satellite Attitude Dynamics and External Torques via Unscented Kalman Filter

ABSTRACT: In this study, an Unscented Kalman Filter (UKF) algorithm is designed for estimating the attitude of a picosatellite and the in-orbit external disturbance torques. The estimation vector is formed by the satellite’s attitude, angular rates, and the unknown constant components of the external disturbance torques acting on the satellite. The gravity gradient torque, residual magnetic moment, sun radiation pressure and aerodynamic drag are all included in the estimated external disturbance torque vector. The satellite has magnetometers and gyros onboard as the attitude sensors. Because of the inherent nonlinear dynamics and the nonlinear measurement model, the UKF, which is a nonlinear version of the Kalman Filter, is selected as the filter algorithm. Performance of the proposed algorithm is demonstrated via simulations for a cube pico-satellite and the results are analyzed for different scenarios.


INTRODUCTION
Although there are numerous researches on cubesats, and this number is increasing day by day, the investigations are still far from being concluded.The Cubesat is basically a cubic picosatellite which has a volume of 1 liter and mass of no more than 1.3 kg.These types of satellites are the outcomes of a search for lighter, smaller and cheaper spacecrafts, and recently, they have mostly been considered as a part of research projects of organizations like universities (Toorian et al., 2008).
The biggest difficulty of cubesat applications is the limitations on the size and mass of the satellite.All the subsystems must be designed regarding these limitations and, unquestionably, a trade-off between the performance and applicability might be required.Specifically for the attitude determination and control system (ADCS), we cannot use highly accurate sensors and actuators onboard a cubesat, since they are usually heavy and large.The attitude of the satellite must be determined and controlled using miniaturized sensors and actuators, which are usually coarser and less accurate (Candini et al., 2012).In this case, one technique to increase the system performance is to properly select the onboard running ADCS algorithm.
The Extended Kalman Filter (EKF) is a nonlinear Kalman Filter (KF), which has been widely used for satellite attitude estimation (Psiaki et al., 1990;Lefferts et al., 1982).On the other hand, the EKF has some disadvantages, especially for the highly nonlinear systems.Generally, this is caused by the mandatory linearization phase of the EKF procedure and so the Jacobians derived with that purpose.For most of the applications, the generation of the Jacobians is time consuming, difficult and prone to human errors (Julier and Uhlmann, 2004;Sekhavat et al., 2007).Nonetheless, the linearization brings about an unstable filter performance when the time steps for the update is not sufficiently small and so the estimation or identification procedure fails as the filter diverges (Julier et al., 1995).Per contra, small time steps increases the computational burden because of the larger number of Jacobian calculations.As a result of these facts, the EKF may be efficient only if the system is almost linear on the timescale of update intervals (Julier and Uhlmann, 2004).
A relatively new Kalman filtering technique, which does not have the shortcomings of the EKF for the nonlinear systems, is the Unscented Kalman Filter (UKF).The UKF generalizes the Kalman filter for both the linear and nonlinear systems and, in case of nonlinear dynamics, it provides relatively more accurate estimation results than other known observer design methodologies such as the EKF.The essence of the UKF is the fact that the approximation of a nonlinear distribution is easier than the approximation of a nonlinear function or transformation (Julier et al., 2000).The UKF introduces sigma points for catching higher order statistics of the system.It satisfies both better estimation accuracy and convergence characteristics by securing higher order information (Sekhavat et al., 2007).Furthermore, the UKF is more robust against the initial estimation errors, so even in case of inaccurate a priori knowledge about the initial condition of the states, it performs well.The only disadvantage of the UKF is the increase in the computational burden when compared to the EKF.Yet, the computational load depends on the type of application and there are several researches aiming at making the UKF computationally more efficient (Li et al., 2013).
There are many documented studies for the UKF as an estimation algorithm in astronautics.In (Crassidis and Markley, 2003) it is used as a state estimator, while both the states and the parameters of the satellite are estimated by the UKF in Dyke et al. (2004), Sekhavat et al. (2007) and Sekhavat et al. (2009).Moreover, in Vinther et al. (2011) the UKF is preferred as a part of the inexpensive cubesat attitude estimation method, where also the magnetometer biases are estimated and in Inamori et al. (2009) it is used for in-orbit magnetic disturbance estimation and compensation.In Soken and Hajiyev (2011) and Soken and Sakai (2011), both the magnetometer and the gyro biases are estimated as well as the attitude of the satellite by using the UKF.In Soken and Sakai (2013), two different UKFs are run onboard a nanosatellite for attitude estimation, sensor calibration and residual magnetic moment compensation.In Oliveria et al.
(2014), the UKF is used for attitude estimation of a nanosatellite and the results are compared with the EKF.In Li et al. (2013), an adaptive version of the UKF is implemented as a part of the ADCS of a nanosatellite.Lastly, in Cornejo et al. (2010), the UKF is used for state estimation and fault detection purposes within the nanosatellite attitude control system.
In this study, an UKF algorithm is designed for estimating a pico-satellite' s attitude and the in-orbit external disturbance torques acting on the pico-satellite.Gravity-gradient torque, sun pressure, aerodynamic drag and residual magnetic moment are included in the estimated disturbance torque vector.The satellite, for which the algorithm is proposed, has magnetometers and gyros onboard as attitude sensors.The main contribution of this study is to show that it is possible to estimate the unknown external torques as well as the attitude dynamics parameters of the satellite via a simple UKF based algorithm.In addition to the investigations in (Soken and Hajiyev, 2009), a more detailed discussion is given for the torque estimation with an extended simulation scenario.Moreover the case for uncertainty in the satellite's moments of inertia is covered.
The paper proceeds as follows: the mathematical model of the pico-satellite attitude dynamics is presented in the following section.In "The sensors measurement models" the sensors' measurement models are given."Unscented Kalman filter (UKF) for the estimation of attitude dynamics and external torques" contains information about the UKF algorithm for the attitude state and unknown external torque estimation.In "Simulations", the performance of the proposed algorithm is demonstrated via simulations.And in the last section of this study, there is a brief summary of the obtained results and the conclusion.

PICO-SATELLITE ATTITUDE DYNAMICS MATHEMATICAL MODEL
Although the Euler angles may have singularity in some certain cases as an attitude kinematics representation technique, they are physically more significant and easier to interpret than the quaternions.Moreover, when the quaternions are used for kinematics representation, implementing the UKF becomes more difficult because of the quaternion norm constraint (Crassidis and Markley, 2003).Hence, in this study, we preferred to use the Euler angles for simplicity and presumably to reduce the computational load.
When we use the Euler angles the mathematical model of the satellite can be expressed with a 9 dimensional system vector.

Estimation of Pico-Satellite Attitude Dynamics and External Torques via Unscented Kalman Filter
The Euler angles representing the attitude of the satellite (φ is the roll angle about x axis; θ is the pitch angle about y axis; ψ is the yaw angle about z axis), the body angular rates with respect to the inertial axis frame and the constant components of the external disturbance torques form the state vector: (1) where (2) (3) ω is the angular velocity vector of the body frame with respect to the inertial frame and N d is the external disturbance torques vector.Since we only estimate the constant components of the torques in this study: (4) The dynamic equations of the satellite can be derived by the use of the angular momentum conservation law (Wertz, 1998): (5) where J is the inertia matrix, consisting of the main moments of inertia as J=diag(J x ,J y ,J z ) and N c is the applied control torque.In this study it is assumed that there is no attitude controller so N c =0 Kinematic equations of the pico-satellite in terms of Euler angles can be given as: . (6) Here, c (.), s (.) and t (.) are the cosine, sine and tangent functions, respectively.Besides, p, q and r are the components of the ω BR vector, which indicates the angular velocity of the body frame with respect to the reference frame (orbit frame).ω BI and ω BR can be related via: (7) where ω 0 denotes the angular velocity of the orbit with respect to the inertial frame, found as .A represents the direction cosine matrix (Wertz, 1998):

THE SENSORS' MEASUREMENT MODELS
The investigated cubic pico-satellite has two types of sensors onboard: the magnetometers, which measure the strength of the Earth's magnetic field; and the gyros, which provide the angular rates of the satellite with respect to the inertial frame.
In this section, measurement models of these two sensors are presented.

THE MAGNETOMETER MEASUREMENT MODEL
The Earth magnetic field vector components can be modeled in the orbit frame as a function of time (Sekhavat et al., 2007): Here, M e is the magnetic dipole moment of the Earth (M e = 7.943 x 10 15 Wb.m); μ is the the Earth Gravitational constant (μ = 3.98601 x 10 14 m 3 /s 2 ); i is the orbit inclination; ε is the magnetic dipole tilt (ε = 11.7 0 ); ω e is the spin rate of the Earth (ω e = 7.29 x 10 -5 rad/s); r 0 is the distance between the centre of mass of the satellite and the Earth.
Three onboard magnetometers of the pico-satellite measures the components of the magnetic field vector in the body frame.Therefore, for the measurement model, which characterizes the measurements in the body frame, the magnetic field vector terms must be transformed by the use of the direction cosine matrix A. The overall measurement model may be given as: (12) where, H 1 (t), H 2 (t) and H 3 (t) represent Earth's magnetic field vector components in the orbit frame as a function of time, and H x (ϕ, θ, ψ, t) H y (ϕ, θ, ψ, t) and H x (ϕ, θ, ψ, t) show the measured Earth magnetic field vector components in the body frame as a function of time and varying attitude.Furthermore, η 1 is the zero average Gaussian white noise with the characteristic of ( 13) Here I 3x3 is the identity matrix with the dimension of 3 x 3, σ m is the standard deviation of each magnetometer error and δ kj is the Kronecker delta.
In this study, we assume that the magnetometers are calibrated using either an on-ground or in-orbit technique so no bias terms are included in Eq.( 12).

THE GYRO MEASUREMENT MODEL
As aforementioned, the other attitude sensor onboard the satellite is the gyro.Widely used model for the gyro measurements is as follows: where,ω BI,meas is the measured angular rates of the satellite, and η 2 is the zero mean Gaussian white noise with the characteristic of: Here, σ v is the standard deviation of each rate gyro random error.Same as the magnetometers, the gyros are assumed to be calibrated.

UNSCENTED KALMAN FILTER (UKF) FOR THE ESTIMATION OF ATTITUDE DYNAMICS AND EXTERNAL TORQUES
The essence of the UKF is the unscented transform, a deterministic sampling technique, that we use for obtaining a minimal set of sample points (or sigma points) from the a priori average and covariance of the states.Then, these sigma points go through nonlinear transformation.The posterior average and the covariance are determined using the transformed sigma points (Julier et al., 1995;Crassidis and Markley, 2003).
The UKF is derived for discrete-time nonlinear equations, so the system model is given by; (16) (17) Here x k is the state vector and is the measurement vector.Moreover w k and v k are the process and measurement error noises, which are assumed to be Gaussian white noise processes with the covariances of Q (k) and R (k) respectively.
The UKF is based on the determination of 2n+1 sigma points with an average of x (k|k) and a covariance of P (k|k).For an n dimensional state vector, these sigma points are obtained by: where, x 0 (k|k), x 1 (k|k) and x 1+n (k|k) are sigma points, n is the state number and k is the scaling parameter which is used for fine tuning, and the heuristic for choosing this parameter is n+k=3 (Julier and Uhlmann, 2004).corresponds to the l th column of the indicated matrix and l is given as l=1... n.
The next step of the UKF process is transforming each sigma point by the use of system dynamics: Here, is the predicted average and P (k+1|k) is the predicted covariance.Furthermore, the predicted observation vector is: (24) where ( 25) After that, the observation covariance matrix is determined as: where the innovation covariance is ( 27) Here R (k+1) is the measurement noise covariance matrix.On the other hand, the cross correlation matrix can be obtained as: (28) The following part is the update phase of the UKF algorithm.At that phase, first by using measurements, y (k+1), the innovation sequence is found as: Here, (k+1|k+1) is the estimated state vector, and P (k+1|k+1) is the estimated covariance matrix.

SIMULATIONS
The proposed UKF based estimation method is tested in this section.The simulations are performed for 40,000 s, which is the period for almost 7 orbits of the satellite.The orbit of the cubesat is a circular orbit with an altitude of r=550km and an inclination of i=97º.The inertia matrix of the satellite is J=diag(2.1x10 - ,2.0x10 -3 ,1.9x10 -3 )kg.m 2 , which corresponds to a 10cm cubic satellite with an approximate mass of 1.2kg, as mentioned.
The magnetometer sensor noise is characterized by zero mean Gaussian white noise, with a standard deviation of σ m =300 nT.The rate gyro random error is taken as .As the filter parameter for the UKF, ĸ is selected as ĸ=-3, which is different than the suggested heuristic.
The initial attitude errors for all the simulations are set to few degrees.The initial estimates for the angular rates and the external disturbance torques are all zero.Besides, the initial value of the covariance matrix is taken as P 0 =10 -10 for all the estimated states, while the process noise covariance matrix is selected as Q=10 -29 for the attitude and Q=10 -20 for the rest of the states.The R matrix is composed of the sensor noise covariances for the magnetometers and gyros, which means 9x10 -12 T 2 as 3 diagonal components corresponding to magnetometer measurements and 64 x 10 -10 rad 2 /s 2 as 3 diagonal elements corresponding to the gyro measurements.The constant torque term used for modeling purposes in the simulations (the values which will be estimated) is Nd=[5 -3 4] T x 10 -1 μNm.
The algorithm was tested in two cases.We are running a dynamics based model for the estimations (in other words, since we are using the satellite dynamics in addition to the kinematics) so there might be uncertainties, mainly caused by the mismatch between the real dynamics and the model used in the filter.The only source for such uncertainty in our simulations might be the inertia terms, which are not exactly known.Hence the two scenarios for the simulations are the cases with and without the uncertainty in the inertia terms.

l=0... 2n
In the first scenario, there is no uncertainty in inertia terms.In Fig. 1, an Euler angle estimation example is given.The top plot compares the estimation with the actual value, and the lower one presents the estimation error.As seen, the UKF accurately estimates the attitude, and after almost 1.5 orbits (≈11000s.)the estimation error falls below ±1deg.If the attitude determination accuracy requirement for such pico-satellite without attitude control is considered (±1deg of determination accuracy may be accepted as sufficient), then we may state that the UKF algorithm works well and provides sufficient accuracy for pico-satellite missions.Besides, as presented in Fig. 2, the UKF also has a good estimation characteristic for the angular rates.
As to the main purpose of the proposed algorithm, the Fig. 3 presents an example for the estimation of the external disturbance torques.As clearly seen, the UKF converges to the real values after 10,000th s, and gives sufficiently accurate estimations for the torque terms.Therefore, we can say that the given UKF algorithm estimates all of the states accurately.The estimation results for the rest of the parameters are similar and can be seen in Figs.4-6.The key contribution of the method is proving that such estimation can be achieved by a simple algorithm that is computationally nondemanding and efficient enough for cubesat applications.
The second simulation is performed for a scenario where we do not know the inertia of the satellite accurately, and there are 5% and 10% uncertainty in inertia terms, respectively.For understanding sake, the results are given in table with comparison with the results obtained for the 0.5 0 1 1.5 2 2.5 3 3.5 4 x 10 4 x 10 4 first scenario, beforehand.In Table 1, root mean square errors are tabulated for 10,000 s. between 20,000 th and 30,000 th seconds such that; As it can be seen, even for a small uncertainty in inertia terms of the satellite, the estimation error increases for all the estimated states.Yet, for the simulation case with 5% uncertainty, estimation error is still within the acceptable bounds for cubesat applications, regarding that the attitude estimation accuracy is better than ±1deg.However, if the level of uncertainty is 10%,    then it is not possible to satisfy these requirements and the attitude estimation error becomes more than ±1deg.Further examinations show that if there is even more uncertainty in the inertia terms of the satellite, the filter may diverge in long term.

CONCLUSION
In this study, an Unscented Kalman Filter (UKF) algorithm is designed for estimating a pico-satellite's attitude and the in-orbit external disturbance torques acting on the pico-satellite.The gravitygradient torque, sun pressure, aerodynamic drag and residual magnetic moment are included in the estimated disturbance torque vector.The satellite, for which the algorithm is proposed, has magnetometers and gyros onboard as attitude sensors.The main contribution of this study is to show that it is possible to estimate the unknown external torques as well as the attitude dynamics parameters of the satellite via a simple UKF based algorithm.The results show that by using the presented algorithm it is possible to estimate the attitude of the cubesat with an accuracy better than ±1deg, which is accurate enough for such satellites.Moreover the case for uncertainty in the satellite' s moments of inertia is addressed with additional simulations.It is shown that if the uncertainty is less than 5%, the requirements for the attitude estimation performance are still satisfied although the estimation error increases.

Figure 4 .
Figure 4. Estimation of the pitch and yaw angles.

Figure 5 .
Figure 5. Estimation of the angular velocities about "y" and "z" axes.

Figure 6 .
Figure6.Estimation of the constant external torques about "x" and "z" axes.

Table 1 .
Root mean square error for the cases with and without uncertainty in the inertia terms of the satellite.