# Microsatellite Attitude Determination and Control Subsystem Design and Implementation: Software-in-the-Loop Approach

Authors: Ho-Nien Shou
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

The paper describes the development of a microsatellite attitude determination and control subsystem (ADCS) and verification of its functionality by software-in-the-loop (SIL) method. The role of ADCS is to provide attitude control functions, including the de-tumbling and stabilizing the satellite angular velocity, and as well as estimating the orbit and attitude information during the satellite operation. In Taiwan, Air Force Institute of Technology (AFIT), dedicating for students to design experimental low earth orbit micro-satellite, called AFITsat. For AFITsat, the operation of the ADCS consists of three modes which are initialization mode, detumbling mode, and normal mode, respectively. During the initialization mode, ADCS collects the early orbit measurement data from various sensors so that the data can be downlinked to the ground station for further analysis. As particularly emphasized in this paper, during the detumbling mode, ADCS implements the thrusters in plus-wide modulation control method to decrease the satellite angular velocity. ADCS provides the attitude determination function for the estimation of the satellite state, during normal mode. The three modes of microsatellite adopted Kalman filter algorithm estimate microsatellite attitude. This paper will discuss using the SIL validation ADCS function and verify its feasibility.


---

## Body

## 1. Introduction

The purpose of the satellite project is to establish indigenous capability for spacecraft development, including the design, analysis, manufacturing, assembly, integration, tests, and in-space validation. The paper is relevant to the design and implement of the attitude determination and control subsystem (ADCS) for a microsatellite simulation in software. The simulation of the system in software, settle SIL, is the aim of the paper. In a satellite mission, the ADCS plays an important role as it provides the attitude information of the satellite and stabilizes the satellite. The design of the ADCS is typically divided into several phases, including simulation of the system in software, implementation into hardware, realization and implementation of coding in hardware, and the verification test of the system. SIL for research development easy and do not have to spend too much hardware funding; it is suitable for students in the lab as a topic of special purpose.Micro-satellite Attitude Determination and Control Subsystem Design and Implementation: Software-in-the-loop Approach.To support the satellite mission, the ADCS is required to provide the in-orbit attitude control and determination functions. The attitude control function includes the stabilization and detumbling control of the satellite to the prerequisite value, while the attitude determination function is to facilitate the estimated attitude information. ADCS will also operate and record the pertinent data after the separation of satellite from the launcher.SIL is on a mathematical model, in theory, to complete the preliminary microsatellite ADCS design. The SIL platform contains a dynamic simulator and a real-time microcontroller, as well as some interfacing circuitry. The dynamic simulator is capable of performing simulation of the space environment, orbit dynamics, attitude dynamics, and sensor/actuator models. The controller/estimator of ADCS is a realization of the embedded microcontroller for attitude determination and control. In this paper, the requirements of AFITsat ADCS will be briefly reviewed. The design and implementation of AFITsat ADCS will also be described [1–3]. Furthermore, the design of the ADCS design can then be verified by the SIL test and the performance of the ADCS can be observed from the simulation results.In designing the ADCS, the satellite environment, dynamics, control, and estimation functions are simulated through a SIL simulation endeavor. SIL simulation is a virtual platform which can provide the simulation of the satellite working environment, satellite dynamics, and ADCS instruments including attitude sensors and actuators. Besides offering a flexible working platform for the early design stage of ADCS, SIL also eases the modification of the design by modulized all the models. The attitude determination filter and attitude control law are verified through this SIL. After the simulation of ADCS in software, the ADCS algorithm has to be implemented in the ADCS on board microcontroller. To fully develop a verifiable ADCS algorithm, a processor-in-the-loop (PIL) test is developed. PIL test is a way to bridge the gap between the simulation results and the final system construction. It increases the realism of the simulation and provides access to hardware features that are currently not available in the software-only simulation. PIL enables the testing of the actual control and estimation codes running on the microcontroller in the verification environment.This paper is organized into six sections, which are introduction, ADCS design, ADCS verification, verification results, and conclusions, respectively. In Section2, the satellite model and space environment are briefly reviewed. The dynamics and kinematics models of the satellite are derived. Moreover, the design of attitude controller and filter are described. The verification methods and platform setups are addressed in Section 3. Section 4 reviews nonlinear H
∞ robust control law. Section 5 is brief introduction of the algorithm of unscented Kalman filter (UKF) and attitude rate estimator is described. After that, the results from ADCS simulation and verification test are provided and discussed in Section 4. At last, the paper is concluded with the summary and some remarks for the future research are depicted.


## 2. ADCS Design

### 2.1. Satellite Model

The dynamics of the spacecraft in inertial space are governed by Euler’s equations of motion [4–6]
(1)
J
ω
˙
+
ω
×
J
ω
=
τ
c
+
τ
d
,
where J is the moment of inertia which is symmetric:
(2)
J
=
J
0
+
J
Δ
=
[
J
x
x
0
0
0
J
y
y
0
0
0
J
z
z
]
+
[
Δ
J
x
x
0
0
0
Δ
J
y
y
0
0
0
Δ
J
z
z
]
,
with J
i
i being the moment of inertia along the i-axis and moment of inertia uncertain elements are expanded around the nominal point in ith direction with uncertainty, J
i
i
+
Δ
J
i
i, Δ
J
i
i the inertia moment decrease progressively in submodes transition, ω is the vector of angular rates of the satellite in body frame, and ω
i is the satellite angular rates in the ith axis, the moment acting on satellite which consists of control torque applied to the satellite τ
c and disturbance torques τ
d, τ
d including gravity gradient force, air resistance force, and magnetic disturbed force various sources.ω
×
=
[
0
-
ω
z
ω
y
ω
z
0
-
ω
x
-
ω
y
ω
x
0
] is a corresponding ω opposing skew matrix.The kinematic update of the satellite is realized by using the quaternion representation. The following vector set of differential equations is used:(3)
q
˙
=
[
-
ω
×
ω
-
ω
T
0
]
[
q
q
4
]
,
where q is the attitude quaternion vector q
=
[
q
q
4
]
T, q
=
[
q
1
q
2
q
3
]
T, for convenient. The derivative of quaternion in (3) can also be arranged as
(4)
q
˙
=
1
2
Q
ω
with
(5)
Q
=
[
q
4
-
q
3
q
2
q
3
q
4
-
q
1
-
q
2
q
1
q
4
-
q
1
-
q
2
-
q
3
]
.


### 2.2.  Environment Model

The space environment must be considered in modeling and analyzing the motion of spacecraft. The space environment models including magnetic field model, Earth’s upper atmosphere model, and sun position model are described in the following. In the magnetic field model, the International Geomagnetic Reference Field (IGRF) is an attempt by the International Association of Geomagnetism and Aeronomy (IAGA) to provide a geomagnetic field mathematical model [7]. In the SIL design, the latest IGRF model, IGRF 2010, is adopted.Models of the upper atmosphere are usually based on either empirical or theoretical work. The main concern is about the total atmospheric density in upper atmospheric, which will affect the satellite attitude. The atmospheric density is determined by a modified analytical expression of the Jacchia-Roberts theory [8]. For microsatellite which is planned to operate at the altitude of 550 km, the atmospheric density is about 1.454
×
10
-
13 kgm−3.By knowing the Sun’s orbit and position with respect to the satellite, the eclipse/shadowing status of satellite can be determined. The microsatellite ADCS is concerned with the satellite solar eclipse because the sun sensor reading is affected. To this end, the sun position model is also incorporated.Satellites are subject to many types of disturbances which will affect the attitude of the satellite in space. The dominant sources of attitude disturbance torques are the gravity-gradient torque, aerodynamic torque, and magnetic disturbance. The gravity gradient torque is one of the largest sources of disturbance that would affect the low earth orbit satellite. In orbit, any object which is nonsymmetrical with finite dimensions will be subjected to a gravitational torque due to the variation in the Earth’s gravitational force over the object. The gravity-gradient torque act on the satellite can be expressed as follows:(6)
τ
g
g
=
3
μ
∥
r
s
∥
3
r
s
×
(
I
r
s
)
,
where μ is the earth’s gravity constant and r
s is the location of the satellite mass center relative to the earth center in body frame [4, 8]. For the microsatellite operation, the gravity-gradient torque is computed. The magnitude is of the order 10−7 Nm.The interaction of the upper atmosphere with the surface of satellite will produce a torque about the center of mass. The aerodynamic torque is given by the following equation:(7)
τ
ad
=
∑
k
=
1
6
r
k
×
f
ad
,
k
with
(8)
f
ad
,
k
=
{
1
2
C
d
ρ
V
2
A
k
(
n
k
·
v
→
)
v
→
,
n
k
·
v
→
>
0
0
,
n
k
·
v
→
≤
0
.In the above, as the microsatellite is of rectangular shape, there are six surfaces that need to be considered. The location of thekth plate is denoted by r
k. The corresponding force due to the impact of molecules on the surface is denoted by f
ad
,
k. This force component is a function of the drag coefficient C
d, atmospheric density ρ, magnitude of the translational velocity of the satellite V, area of the plate A
k, unit normal vector of the surface n
k, and the unit vector of the velocity vector v
→. The magnitude of the aerodynamic torque acting on the microsatellite is computed to be of the order 10−8 Nm.Finally, magnetic disturbance torques result from the interaction between the satellite’s residual magnetic field and the earth geomagnetic field. The dominant source of the residual magnetic field is the spacecraft magnetic moments. This instantaneous magnetic disturbance torque due to the spacecraft effective magnetic moment and the geocentric magnetic field can be described as follows [9–13]:
(9)
τ
mag
=
m
×
b
body
,
where m is the sum of the internal residual magnetic dipole moment and b
body is the geomagnetic flux intensity in body frame.The magnetic field of microsatellites in orbit coordinates frame is as shown in Figure1.Figure 1
The magnetic field of microsatellites in orbit coordinates frame.


## 2.1. Satellite Model

The dynamics of the spacecraft in inertial space are governed by Euler’s equations of motion [4–6]
(1)
J
ω
˙
+
ω
×
J
ω
=
τ
c
+
τ
d
,
where J is the moment of inertia which is symmetric:
(2)
J
=
J
0
+
J
Δ
=
[
J
x
x
0
0
0
J
y
y
0
0
0
J
z
z
]
+
[
Δ
J
x
x
0
0
0
Δ
J
y
y
0
0
0
Δ
J
z
z
]
,
with J
i
i being the moment of inertia along the i-axis and moment of inertia uncertain elements are expanded around the nominal point in ith direction with uncertainty, J
i
i
+
Δ
J
i
i, Δ
J
i
i the inertia moment decrease progressively in submodes transition, ω is the vector of angular rates of the satellite in body frame, and ω
i is the satellite angular rates in the ith axis, the moment acting on satellite which consists of control torque applied to the satellite τ
c and disturbance torques τ
d, τ
d including gravity gradient force, air resistance force, and magnetic disturbed force various sources.ω
×
=
[
0
-
ω
z
ω
y
ω
z
0
-
ω
x
-
ω
y
ω
x
0
] is a corresponding ω opposing skew matrix.The kinematic update of the satellite is realized by using the quaternion representation. The following vector set of differential equations is used:(3)
q
˙
=
[
-
ω
×
ω
-
ω
T
0
]
[
q
q
4
]
,
where q is the attitude quaternion vector q
=
[
q
q
4
]
T, q
=
[
q
1
q
2
q
3
]
T, for convenient. The derivative of quaternion in (3) can also be arranged as
(4)
q
˙
=
1
2
Q
ω
with
(5)
Q
=
[
q
4
-
q
3
q
2
q
3
q
4
-
q
1
-
q
2
q
1
q
4
-
q
1
-
q
2
-
q
3
]
.


## 2.2.  Environment Model

The space environment must be considered in modeling and analyzing the motion of spacecraft. The space environment models including magnetic field model, Earth’s upper atmosphere model, and sun position model are described in the following. In the magnetic field model, the International Geomagnetic Reference Field (IGRF) is an attempt by the International Association of Geomagnetism and Aeronomy (IAGA) to provide a geomagnetic field mathematical model [7]. In the SIL design, the latest IGRF model, IGRF 2010, is adopted.Models of the upper atmosphere are usually based on either empirical or theoretical work. The main concern is about the total atmospheric density in upper atmospheric, which will affect the satellite attitude. The atmospheric density is determined by a modified analytical expression of the Jacchia-Roberts theory [8]. For microsatellite which is planned to operate at the altitude of 550 km, the atmospheric density is about 1.454
×
10
-
13 kgm−3.By knowing the Sun’s orbit and position with respect to the satellite, the eclipse/shadowing status of satellite can be determined. The microsatellite ADCS is concerned with the satellite solar eclipse because the sun sensor reading is affected. To this end, the sun position model is also incorporated.Satellites are subject to many types of disturbances which will affect the attitude of the satellite in space. The dominant sources of attitude disturbance torques are the gravity-gradient torque, aerodynamic torque, and magnetic disturbance. The gravity gradient torque is one of the largest sources of disturbance that would affect the low earth orbit satellite. In orbit, any object which is nonsymmetrical with finite dimensions will be subjected to a gravitational torque due to the variation in the Earth’s gravitational force over the object. The gravity-gradient torque act on the satellite can be expressed as follows:(6)
τ
g
g
=
3
μ
∥
r
s
∥
3
r
s
×
(
I
r
s
)
,
where μ is the earth’s gravity constant and r
s is the location of the satellite mass center relative to the earth center in body frame [4, 8]. For the microsatellite operation, the gravity-gradient torque is computed. The magnitude is of the order 10−7 Nm.The interaction of the upper atmosphere with the surface of satellite will produce a torque about the center of mass. The aerodynamic torque is given by the following equation:(7)
τ
ad
=
∑
k
=
1
6
r
k
×
f
ad
,
k
with
(8)
f
ad
,
k
=
{
1
2
C
d
ρ
V
2
A
k
(
n
k
·
v
→
)
v
→
,
n
k
·
v
→
>
0
0
,
n
k
·
v
→
≤
0
.In the above, as the microsatellite is of rectangular shape, there are six surfaces that need to be considered. The location of thekth plate is denoted by r
k. The corresponding force due to the impact of molecules on the surface is denoted by f
ad
,
k. This force component is a function of the drag coefficient C
d, atmospheric density ρ, magnitude of the translational velocity of the satellite V, area of the plate A
k, unit normal vector of the surface n
k, and the unit vector of the velocity vector v
→. The magnitude of the aerodynamic torque acting on the microsatellite is computed to be of the order 10−8 Nm.Finally, magnetic disturbance torques result from the interaction between the satellite’s residual magnetic field and the earth geomagnetic field. The dominant source of the residual magnetic field is the spacecraft magnetic moments. This instantaneous magnetic disturbance torque due to the spacecraft effective magnetic moment and the geocentric magnetic field can be described as follows [9–13]:
(9)
τ
mag
=
m
×
b
body
,
where m is the sum of the internal residual magnetic dipole moment and b
body is the geomagnetic flux intensity in body frame.The magnetic field of microsatellites in orbit coordinates frame is as shown in Figure1.Figure 1
The magnetic field of microsatellites in orbit coordinates frame.


## 3. Attitude Control

Attitude control is the process of changing the orientation of satellite. This process is vital because a satellite will always experience internal and external disturbance torques. The attitude control in microsatellite is mainly concerned about the attitude stabilization. Approaches of attitude control have been discussed in this.


### 3.1. Production of Satellite Thrust [14–16]

The satellite exercise frame is as in Figure2. Four thrusters are put in x
-
y plane, including angle, the jet moves toward -
z direction, and lines to be a square distance d. Fuel is put the jet x
-
y plane l distance, and +
z direction is toward satellite center, in Figure 3.Figure 2
Microsatellite geometric appearance.Figure 3
The geometric arrangement of the thrusters.TheΔ
λ, Δ
d, and Δ
l, respectively, represent the jet including angle, jet relative distance, and the uncertainty of the fuel gravity position.The following is productionτ
u from the thruster:
(10)
τ
u
=
[
τ
x
τ
y
τ
z
]
=
[
r
→
1
×
F
→
1
r
→
2
×
F
→
2
r
→
3
×
F
→
3
r
→
4
×
F
→
4
]
[
u
1
u
2
u
3
u
4
]
=
[
-
β
1
Δ
-
β
1
Δ
β
1
Δ
β
1
Δ
-
β
2
Δ
β
2
Δ
-
β
2
Δ
-
β
2
Δ
β
3
Δ
-
β
3
Δ
-
β
3
Δ
β
3
Δ
]
[
u
1
u
2
u
3
u
4
]
=
β
u
,
where
(11)
1
2
3
4
r
→
=
[
-
(
d
+
Δ
d
)
-
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
-
(
d
+
Δ
d
)
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
(
d
+
Δ
d
)
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
(
d
+
Δ
d
)
-
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
]
x
y
z
F
→
1
=
[
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
F
→
2
=
[
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
F
→
3
=
[
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
F
→
4
=
[
ρ
sin
⁡
(
λ
+
Δ
λ
)
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
α
=
1
2
,
β
1
=
d
+
α
λ
l
,
β
2
=
d
-
α
λ
l
,
β
3
=
λ
d
α
,
Δ
β
1
=
Δ
d
-
λ
Δ
λ
(
d
+
Δ
d
)
+
α
(
λ
Δ
l
+
l
Δ
λ
+
Δ
λ
Δ
l
)
,
Δ
β
2
=
Δ
d
-
λ
Δ
λ
(
d
+
Δ
d
)
-
α
(
λ
Δ
l
+
l
Δ
λ
+
Δ
λ
Δ
l
)
,
Δ
β
3
=
1
α
(
λ
Δ
d
+
d
Δ
d
+
Δ
λ
Δ
d
)
,
β
i
Δ
=
β
i
+
Δ
β
i
,
i
=
1,2
,
3
.


### 3.2. Pulse-Width Modulation Rules [17–22]

The controller requires a continuum of positive and negative thrust values, which are physically impossible with the four fixed thrusters. This requirement can be met by pulse-width modulating the thrusters. The four thrusters span the three-dimensional control space when generating a positive thrust along the satellitez-axis. The thruster commands are all biased upwards so that each thruster receives a positive thrust command. Each thruster can be commanded on for a variable time during each control system sample period, with a resolution of 50 msec. The variable thrust level from the linear control signal is realized by turning on each thruster for less than the full sample period, such that the impulse of the thruster matches the commanded impulse.The procedure to convert the controller to pulse-width modulation is as follows.(1)
The thruster commandsu.
(2)
Convertu
(
n
) from variable thrusts for the fixed total sample period to fixed thrusts for variable times: t
on
⁡ = (sample time/thrust magnitude) ×  u.
(3)
Limit on-times to ± half the sample period.
(4)
Compute a bias term; add the bias to all on-times such that at least one thruster is on for the complete sample period.According to the PWM duty rule, the duty period can be divided into two modulations: off-modulation and on-modulation. In this project, an on-modulation is adopted. Microsatellite attitude control pulse-width modulation block diagram is as shown in Figure4. Figure 5 shows PWM processing; from the control flow chart of pulse-width on/off modulation algorithm, there is at least one thruster that works in the off-modulation, while there is at least one thruster that does not work in the on-modulation. In Figures 4 and 5, M
thruster, T
s represents the maximum amplitude of the thruster, sampling time.Figure 4
Microsatellite attitude control pulse-width modulation block diagram.Figure 5
Pulse-width on/off modulation control flow chart.


## 3.1. Production of Satellite Thrust [14–16]

The satellite exercise frame is as in Figure2. Four thrusters are put in x
-
y plane, including angle, the jet moves toward -
z direction, and lines to be a square distance d. Fuel is put the jet x
-
y plane l distance, and +
z direction is toward satellite center, in Figure 3.Figure 2
Microsatellite geometric appearance.Figure 3
The geometric arrangement of the thrusters.TheΔ
λ, Δ
d, and Δ
l, respectively, represent the jet including angle, jet relative distance, and the uncertainty of the fuel gravity position.The following is productionτ
u from the thruster:
(10)
τ
u
=
[
τ
x
τ
y
τ
z
]
=
[
r
→
1
×
F
→
1
r
→
2
×
F
→
2
r
→
3
×
F
→
3
r
→
4
×
F
→
4
]
[
u
1
u
2
u
3
u
4
]
=
[
-
β
1
Δ
-
β
1
Δ
β
1
Δ
β
1
Δ
-
β
2
Δ
β
2
Δ
-
β
2
Δ
-
β
2
Δ
β
3
Δ
-
β
3
Δ
-
β
3
Δ
β
3
Δ
]
[
u
1
u
2
u
3
u
4
]
=
β
u
,
where
(11)
1
2
3
4
r
→
=
[
-
(
d
+
Δ
d
)
-
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
-
(
d
+
Δ
d
)
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
(
d
+
Δ
d
)
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
(
d
+
Δ
d
)
-
(
d
+
Δ
d
)
-
(
l
+
Δ
l
)
]
x
y
z
F
→
1
=
[
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
F
→
2
=
[
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
F
→
3
=
[
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
F
→
4
=
[
ρ
sin
⁡
(
λ
+
Δ
λ
)
ρ
sin
⁡
(
λ
+
Δ
λ
)
-
cos
⁡
⁡
(
λ
+
Δ
λ
)
]
,
α
=
1
2
,
β
1
=
d
+
α
λ
l
,
β
2
=
d
-
α
λ
l
,
β
3
=
λ
d
α
,
Δ
β
1
=
Δ
d
-
λ
Δ
λ
(
d
+
Δ
d
)
+
α
(
λ
Δ
l
+
l
Δ
λ
+
Δ
λ
Δ
l
)
,
Δ
β
2
=
Δ
d
-
λ
Δ
λ
(
d
+
Δ
d
)
-
α
(
λ
Δ
l
+
l
Δ
λ
+
Δ
λ
Δ
l
)
,
Δ
β
3
=
1
α
(
λ
Δ
d
+
d
Δ
d
+
Δ
λ
Δ
d
)
,
β
i
Δ
=
β
i
+
Δ
β
i
,
i
=
1,2
,
3
.


## 3.2. Pulse-Width Modulation Rules [17–22]

The controller requires a continuum of positive and negative thrust values, which are physically impossible with the four fixed thrusters. This requirement can be met by pulse-width modulating the thrusters. The four thrusters span the three-dimensional control space when generating a positive thrust along the satellitez-axis. The thruster commands are all biased upwards so that each thruster receives a positive thrust command. Each thruster can be commanded on for a variable time during each control system sample period, with a resolution of 50 msec. The variable thrust level from the linear control signal is realized by turning on each thruster for less than the full sample period, such that the impulse of the thruster matches the commanded impulse.The procedure to convert the controller to pulse-width modulation is as follows.(1)
The thruster commandsu.
(2)
Convertu
(
n
) from variable thrusts for the fixed total sample period to fixed thrusts for variable times: t
on
⁡ = (sample time/thrust magnitude) ×  u.
(3)
Limit on-times to ± half the sample period.
(4)
Compute a bias term; add the bias to all on-times such that at least one thruster is on for the complete sample period.According to the PWM duty rule, the duty period can be divided into two modulations: off-modulation and on-modulation. In this project, an on-modulation is adopted. Microsatellite attitude control pulse-width modulation block diagram is as shown in Figure4. Figure 5 shows PWM processing; from the control flow chart of pulse-width on/off modulation algorithm, there is at least one thruster that works in the off-modulation, while there is at least one thruster that does not work in the on-modulation. In Figures 4 and 5, M
thruster, T
s represents the maximum amplitude of the thruster, sampling time.Figure 4
Microsatellite attitude control pulse-width modulation block diagram.Figure 5
Pulse-width on/off modulation control flow chart.


## 4. NonlinearH
∞ Robust Control Law

The microsatellite kinematic equation is described by(12)
q
˙
=
1
2
[
-
ω
×
ω
-
ω
T
0
]
[
q
q
0
]
=
1
2
[
q
×
+
q
0
I
-
q
T
]
ω
=
1
2
Q
(
q
)
ω
,
where
(13)
Q
(
q
)
≜
1
2
[
q
×
+
q
0
I
-
q
T
]
.In the attitude control design, it is desired that the system is stable and the excursions of the angular rate and control input are minimized. This motivates the use of the following penalty function to be minimized:(14)
z
=
[
ρ
1
ω
T
J
Δ
ω
+
ρ
2
η
2
(
q
0
)
ρ
3
u
]
,
where ρ
1, ρ
2, and ρ
3 are weighting coefficients introduced for the trade-off between performance and control effort. The function η
(
q
0
) is defined as η
(
q
0
)
=
2
cos
⁡
-
1
|
q
0
|.The attitude control problem can then be formulated as(15)
x
˙
=
(
f
(
x
)
+
Δ
f
(
x
)
)
+
(
g
1
(
x
)
+
Δ
g
1
(
x
)
)
w
+
(
g
2
(
x
)
+
Δ
g
2
(
x
)
)
u
(16)
z
=
[
h
(
x
)
u
]
=
h
1
(
x
)
,
where x
=
[
ω
T
q
T
q
0
]
Tis the state x and w is the disturbance. The Δ
f, Δ
g
1, and Δ
g
2 are uncertain smooth vector fields and mappings. Denote δ
J
=
J
-
1
Δ
J
J
-
1, for small perturbations; f, g
1, g
2, h
1, Δ
f, Δ
g
1, and Δ
g
2 can be expressed as
(17)
f
(
x
)
=
[
-
J
-
1
q
×
J
q
Q
(
q
)
ω
]
,
g
1
(
x
)
=
[
J
-
1
0
]
,
g
2
(
x
)
=
[
J
-
1
B
0
]
,
h
1
(
x
)
=
ρ
1
q
T
J
q
+
ρ
2
η
2
(
q
0
)
,
Δ
f
(
x
)
=
[
-
J
-
1
q
×
J
q
+
J
-
1
δ
J
q
×
J
q
+
J
-
1
δ
J
q
×
Δ
J
q
0
]
,
Δ
g
1
(
x
)
=
[
-
δ
J
0
]
,
Δ
g
2
(
x
)
=
[
J
-
1
Δ
B
-
δ
J
(
B
+
Δ
B
)
0
]
.


### 4.1. NonlinearH
∞ Robust Control Theory

Consider a nonlinear state-space systemG
0
(18a)
x
_
˙
=
f
(
x
_
)
+
g
(
x
_
)
w
_
,
(18b)
z
_
=
h
(
x
_
)
,wherex
_ is the state vector, w
_ is the exogenous disturbance to be rejected, and z
_ is the penalized output signal. We assume that f
(
x
_
), g
(
x
_
), and h
(
x
_
) are C
∞ functions and x
_
=
0 is the equilibrium point of the system; that is, f
(
0
)
=
h
(
0
)
=
0.Given a positive numberγ
>
0, system (18a) and (18b) has a finite L
2—gain less than γ if, for all w
_
∈
L
2
[
0
,
T
] with 0
≤
T
<
∞,
(19)
∫
0
T
∥
z
_
(
t
)
∥
2
d
t
≤
γ
2
∫
0
T
∥
w
_
(
t
)
∥
2
d
t
.Lemma 1.
Assume thatΔ
f
(
x
_
)
=
0, Δ
g
1
(
x
_
)
=
0, and Δ
g
2
(
x
_
)
=
0 in (15). If (
f
(
x
_
)
,
h
1
(
x
_
)
) is zero-detectable and there exists a positive-definite function V
(
x
_
) such that the following Hamilton-Jacobi partial differential inequality holds [23–25]
(20)
H
γ
=
V
x
_
T
f
+
1
2
V
x
_
T
(
1
γ
2
g
1
g
1
T
-
g
2
g
2
T
)
V
x
_
+
1
2
h
1
T
h
1
<
0
then the system has a finite L
2 gain less than γ. Once V
(
x
_
) is found, the control signal u
_ can be synthesized as
(21)
u
_
=
-
g
2
T
V
x
_
,
where V
x
_
=
[
∂
V
/
∂
x
1
∂
V
/
∂
x
2
⋯
∂
V
/
∂
x
n
]
T.Theorem 2.
Consider that the nonlinear uncertain system (17) and the uncertainty of system parameters satisfy the following equation [26, 27]:
(22)
Δ
f
(
x
_
)
=
Φ
(
x
_
)
Ψ
T
(
x
_
)
,
Δ
g
1
(
x
_
)
Δ
g
1
T
(
x
_
)
≤
G
1
(
x
_
)
G
1
T
(
x
_
)
,
Δ
g
2
(
x
_
)
Δ
g
2
T
(
x
_
)
≤
G
2
(
x
_
)
G
2
T
(
x
_
)
,
where Φ
(
x
_
)
∈
ℝ
n
×
m, Ψ
(
x
_
)
∈
ℝ
1
×
m, and G
i
(
x
_
)
∈
ℝ
n
×
p
i, i
=
1,2, are known vector, for which we assume that there exist scale functions Θ
(
x
_
)
≥
0 with Ψ
(
x
_
) described by ∥
Ψ
(
x
_
)
∥
≤
Θ
(
x
_
) with Θ
(
0
)
=
0. Under this assumption, if system (15) satisfies the following conditions(i)
(
f
(
x
_
)
+
Δ
f
(
x
_
)
,
h
1
(
x
_
)
) is zero-state detectable;
(ii)
there exists a smooth functionσ
=
σ
(
x
_
)
>
0 such that the following Hamilton-Jacobi partial differential inequality
(23)
V
x
_
T
f
+
1
2
V
x
_
T
(
1
γ
2
(
1
+
σ
)
g
1
g
1
T
-
(
1
-
σ
)
g
2
g
2
T
)
V
x
_
+
1
2
V
x
_
T
(
1
γ
2
1
+
σ
σ
G
1
G
1
T
+
1
-
σ
σ
G
2
G
2
T
+
σ
Φ
Φ
T
)
V
x
_
+
1
2
σ
Ψ
T
Ψ
+
1
2
h
1
T
h
1
≤
0has a smooth positive solution V
(
x
_
) where x
_
∈
X and X is a neighborhood of origin in ℝ
n, then by the following state feedback law
(24)
u
(
x
_
)
=
-
g
2
T
V
x
_
the closed-loop system of (15) has a robust H
∞-performance.


## 4.1. NonlinearH
∞ Robust Control Theory

Consider a nonlinear state-space systemG
0
(18a)
x
_
˙
=
f
(
x
_
)
+
g
(
x
_
)
w
_
,
(18b)
z
_
=
h
(
x
_
)
,wherex
_ is the state vector, w
_ is the exogenous disturbance to be rejected, and z
_ is the penalized output signal. We assume that f
(
x
_
), g
(
x
_
), and h
(
x
_
) are C
∞ functions and x
_
=
0 is the equilibrium point of the system; that is, f
(
0
)
=
h
(
0
)
=
0.Given a positive numberγ
>
0, system (18a) and (18b) has a finite L
2—gain less than γ if, for all w
_
∈
L
2
[
0
,
T
] with 0
≤
T
<
∞,
(19)
∫
0
T
∥
z
_
(
t
)
∥
2
d
t
≤
γ
2
∫
0
T
∥
w
_
(
t
)
∥
2
d
t
.Lemma 1.
Assume thatΔ
f
(
x
_
)
=
0, Δ
g
1
(
x
_
)
=
0, and Δ
g
2
(
x
_
)
=
0 in (15). If (
f
(
x
_
)
,
h
1
(
x
_
)
) is zero-detectable and there exists a positive-definite function V
(
x
_
) such that the following Hamilton-Jacobi partial differential inequality holds [23–25]
(20)
H
γ
=
V
x
_
T
f
+
1
2
V
x
_
T
(
1
γ
2
g
1
g
1
T
-
g
2
g
2
T
)
V
x
_
+
1
2
h
1
T
h
1
<
0
then the system has a finite L
2 gain less than γ. Once V
(
x
_
) is found, the control signal u
_ can be synthesized as
(21)
u
_
=
-
g
2
T
V
x
_
,
where V
x
_
=
[
∂
V
/
∂
x
1
∂
V
/
∂
x
2
⋯
∂
V
/
∂
x
n
]
T.Theorem 2.
Consider that the nonlinear uncertain system (17) and the uncertainty of system parameters satisfy the following equation [26, 27]:
(22)
Δ
f
(
x
_
)
=
Φ
(
x
_
)
Ψ
T
(
x
_
)
,
Δ
g
1
(
x
_
)
Δ
g
1
T
(
x
_
)
≤
G
1
(
x
_
)
G
1
T
(
x
_
)
,
Δ
g
2
(
x
_
)
Δ
g
2
T
(
x
_
)
≤
G
2
(
x
_
)
G
2
T
(
x
_
)
,
where Φ
(
x
_
)
∈
ℝ
n
×
m, Ψ
(
x
_
)
∈
ℝ
1
×
m, and G
i
(
x
_
)
∈
ℝ
n
×
p
i, i
=
1,2, are known vector, for which we assume that there exist scale functions Θ
(
x
_
)
≥
0 with Ψ
(
x
_
) described by ∥
Ψ
(
x
_
)
∥
≤
Θ
(
x
_
) with Θ
(
0
)
=
0. Under this assumption, if system (15) satisfies the following conditions(i)
(
f
(
x
_
)
+
Δ
f
(
x
_
)
,
h
1
(
x
_
)
) is zero-state detectable;
(ii)
there exists a smooth functionσ
=
σ
(
x
_
)
>
0 such that the following Hamilton-Jacobi partial differential inequality
(23)
V
x
_
T
f
+
1
2
V
x
_
T
(
1
γ
2
(
1
+
σ
)
g
1
g
1
T
-
(
1
-
σ
)
g
2
g
2
T
)
V
x
_
+
1
2
V
x
_
T
(
1
γ
2
1
+
σ
σ
G
1
G
1
T
+
1
-
σ
σ
G
2
G
2
T
+
σ
Φ
Φ
T
)
V
x
_
+
1
2
σ
Ψ
T
Ψ
+
1
2
h
1
T
h
1
≤
0has a smooth positive solution V
(
x
_
) where x
_
∈
X and X is a neighborhood of origin in ℝ
n, then by the following state feedback law
(24)
u
(
x
_
)
=
-
g
2
T
V
x
_
the closed-loop system of (15) has a robust H
∞-performance.


## 5. Kalman Filter

Attitude determination involves a two-part process: estimation of a vehicle’s orientation from body measurements and known reference observations and filtering of noisy measurements. The Kalman filter (KF) is essentially a set of mathematical equations that implement a predictor-corrector type estimator so that it can minimize the estimated error covariance. However Kalman filter is only applicable for linear systems and not suitable for the nonlinear satellite dynamic system. The UKF (unscented KF) is a standard state estimation scheme for nonlinear system. The UKF extends the scope of Kalman filter to nonlinear optimal filtering problems by unscented condition estimation which means it is based on the input unscented system and the output measuring data to estimate the unpredictable interior condition. Unscented condition estimation is important in scientific research and constructional applied field. The attitude determination algorithm implemented in microsatellite is UKF. The UKF starts with a given initial estimate for the state and its error covariance matrix. In order to obtain the optimal estimates of the real system, UKF will predict the system state recursively by measurement update and propagation update. The propagation update will be propagated as the prior state and covariance for the next time step. In the design of the microsatellite attitude estimation scheme, several different filter formulations have been assessed. In its simplest form, the states of the UKF include the satellite quaternion and angular rates (both in satellite body frame) [28, 29]:
(25)
x
=
[
q
1
q
2
q
3
q
4
ω
x
ω
y
ω
z
]
T
.Based on the discussions in Section2, the following process and measurement models can be established:(26a)
x
k
+
1
=
f
(
x
k
)
+
w
k
,
(26b)
z
k
=
h
(
x
k
)
+
v
k
.The variables in the above equation will be described.x
k
is a system condition vector in thek  moment.
f
(
·
) is unscented system model.
w
k is a dynamic mixed signal in the k moment.
z
k is a measuring dynamic vector in the k moment.
h
(
·
) is a unscented system measuring model.
v
k is dynamic measuring mixed signal in the k moment.The measurement vector is denoted byz
k. The measurements that are used for attitude determination in the microsatellite ADCS include sensor readings of three-axis gyroscope, sun sensor, and three-axis magnetometer. The process noise w
k and the measurement noise v
k are assumed to be zero-mean white noise. The process noise covariance matrix and the measurement noise covariance are given by Q
j and R
j, respectively.The system error covariance matrixQ
j is as follows:
(27)
E
{
w
i
w
j
T
}
=
{
Q
j
,
i
=
j
0
,
i
≠
j
.The measuring error covariance matrixR
j is as follows:
(28)
E
{
v
i
v
j
T
}
=
{
R
j
,
i
=
j
0
,
i
≠
j
.Here, and are independent and unrelated:(29)
E
{
w
i
v
j
T
}
=
0
,
∀
i
,
j
.


### 5.1. Unscented Kalman Filter Processing [30, 31]

(1)
Producing standard points and calculation value:(30)
χ
0
=
x
^
k
-
,
χ
i
=
x
^
k
-
+
(
(
n
+
ε
)
P
k
)
i
T
i
=
1,2
,
…
,
n
,
χ
i
+
n
=
x
^
k
-
-
(
(
n
+
ε
)
P
k
)
i
T
i
=
1,2
,
…
,
n
,
W
i
m
=
ε
(
n
+
ε
)
,
W
0
c
=
W
0
m
+
(
1
-
ε
1
2
+
ε
2
)
,
W
i
m
=
W
i
c
=
1
2
(
n
+
ε
)
.
The parameter is a scaling parameter defined as(31)
ε
2
=
ε
1
2
(
n
+
ε
3
)
-
n
.
The positive constantsε
i, i
=
1,2
,
3, are used as parameter of the method; a priori and a posteriori estimates of the state are denoted by x
^
k
- and x
^
k.
(2)
Time updating:
condition predicted value:(32)
(
ς
k
-
)
i
=
f
(
(
χ
k
)
i
)
,
condition predicted average:(33)
x
^
k
-
=
∑
i
=
0
2
n
W
i
m
(
ς
k
-
)
i
,
covariance matrix:(34)
P
k
-
=
∑
i
=
0
2
n
W
i
c
[
(
ς
k
-
)
i
-
x
^
k
-
]
[
(
ς
k
-
)
i
-
x
^
k
-
]
T
+
Q
k
.
(3)
Observation updating:
observation measurement predicted value:(35)
(
z
k
-
)
i
=
h
(
(
ς
k
-
)
i
)
,
observation measurement predicted average:(36)
z
^
k
-
=
∑
i
=
0
2
n
W
i
m
(
z
k
-
)
i
,
P
x
y and P
y
y update:
(37)
P
x
y
=
∑
i
=
0
2
n
W
i
c
[
(
ς
k
-
)
i
-
x
^
k
-
]
[
(
ς
k
-
)
i
-
z
^
k
-
]
T
,
P
y
y
=
∑
i
=
0
2
n
W
i
c
[
(
z
k
-
)
i
-
x
^
k
-
]
[
(
z
k
-
)
i
-
z
^
k
-
]
T
+
R
k
.
(4)
Calculating Kalman gain value:(38)
K
k
=
P
x
y
P
y
y
-
1
.
(5)
Updating estimated value to measurement value:(39)
x
^
k
=
x
^
k
-
+
K
k
(
z
k
-
z
^
k
-
)
.
(6)
Updating condition error covariance matrix:(40)
P
k
+
1
=
P
k
-
-
K
k
P
y
y
K
k
T
.


### 5.2. Unscented Kalman Filter Flow Chart

The flow chart of the UKF is as shown in Figure6.Figure 6
The algorithm flow chart of UKF.


## 5.1. Unscented Kalman Filter Processing [30, 31]

(1)
Producing standard points and calculation value:(30)
χ
0
=
x
^
k
-
,
χ
i
=
x
^
k
-
+
(
(
n
+
ε
)
P
k
)
i
T
i
=
1,2
,
…
,
n
,
χ
i
+
n
=
x
^
k
-
-
(
(
n
+
ε
)
P
k
)
i
T
i
=
1,2
,
…
,
n
,
W
i
m
=
ε
(
n
+
ε
)
,
W
0
c
=
W
0
m
+
(
1
-
ε
1
2
+
ε
2
)
,
W
i
m
=
W
i
c
=
1
2
(
n
+
ε
)
.
The parameter is a scaling parameter defined as(31)
ε
2
=
ε
1
2
(
n
+
ε
3
)
-
n
.
The positive constantsε
i, i
=
1,2
,
3, are used as parameter of the method; a priori and a posteriori estimates of the state are denoted by x
^
k
- and x
^
k.
(2)
Time updating:
condition predicted value:(32)
(
ς
k
-
)
i
=
f
(
(
χ
k
)
i
)
,
condition predicted average:(33)
x
^
k
-
=
∑
i
=
0
2
n
W
i
m
(
ς
k
-
)
i
,
covariance matrix:(34)
P
k
-
=
∑
i
=
0
2
n
W
i
c
[
(
ς
k
-
)
i
-
x
^
k
-
]
[
(
ς
k
-
)
i
-
x
^
k
-
]
T
+
Q
k
.
(3)
Observation updating:
observation measurement predicted value:(35)
(
z
k
-
)
i
=
h
(
(
ς
k
-
)
i
)
,
observation measurement predicted average:(36)
z
^
k
-
=
∑
i
=
0
2
n
W
i
m
(
z
k
-
)
i
,
P
x
y and P
y
y update:
(37)
P
x
y
=
∑
i
=
0
2
n
W
i
c
[
(
ς
k
-
)
i
-
x
^
k
-
]
[
(
ς
k
-
)
i
-
z
^
k
-
]
T
,
P
y
y
=
∑
i
=
0
2
n
W
i
c
[
(
z
k
-
)
i
-
x
^
k
-
]
[
(
z
k
-
)
i
-
z
^
k
-
]
T
+
R
k
.
(4)
Calculating Kalman gain value:(38)
K
k
=
P
x
y
P
y
y
-
1
.
(5)
Updating estimated value to measurement value:(39)
x
^
k
=
x
^
k
-
+
K
k
(
z
k
-
z
^
k
-
)
.
(6)
Updating condition error covariance matrix:(40)
P
k
+
1
=
P
k
-
-
K
k
P
y
y
K
k
T
.


## 5.2. Unscented Kalman Filter Flow Chart

The flow chart of the UKF is as shown in Figure6.Figure 6
The algorithm flow chart of UKF.


## 6. The SIL Test of Microsatellite ADCS

The SIL test of microsatellite ADCS is conducted to verify the function of each module. To this end, the Labview [32] with is used to facilitate the simulation.SIL system is divided into two parts: one is for the host computer and the other is the on-board computer (real-time controller NI cRIO-9014, National Instruments) [33], Ethernet connection between the Host computer and on board computer with a group and as for satellite attitude simulation tests. Host computer simulates microsatellite dynamic equations of motion, in orbit on the environment (environmental disturbance torque, 10 orders of the Earth’s magnetic field) and attitude measurement sensing element mathematical model. The on-board computer not only is a real-time hardware system to replace the satellite attitude controller but also simulated satellite attitude control element mathematical model (on/off modulator PWM, thruster). Host computer with on-board computer through the Ethernet connection into a loop, they linked for the satellite attitude (angle, angular velocity) and a control command. The SIL test can be conducted in a subsystem level or closed-loop system level. In SIL testing, the system under consideration is composed of six functional blocks as depicted in Figure 7. The function of each block is described as follows [28–31].(1)
Satellite Attitude and Orbit Information. It provides three main functions including the computation of satellite attitude and angular rates, satellite position and velocity, and coordinate frame transformation.
(2)
Sensor Model. Three types of sensors are modelled and simulated in this block, including magnetometer, sun sensor, and gyro.
(3)
Disturbance Torques. Three types of disturbance torques are considered, including the aerodynamics torque, gravity gradient torque, and magnetic disturbance torque.
(4)
Attitude Determination. This block provides the function of attitude estimation by using UKF.
(5)
Attitude Control. This block executes the attitude control function.
(6)
Attitude Determination Error. This block calculates the state error estimated by the UKF.(a) Block diagram of SIL (b) The photos of SIL experimental.
(a)
(b)The attitude determination and control functions are controlled by switches and can be enabled or disabled anytime during the SIL test.


### 6.1. Parameters of Experimental

The parameters of microsatellite are shown in Table1 and initial conditions of microsatellite are shown in Table 2.Table 1
Parameters of microsatellite.
Item
Uncertainty bound
Value
Unit
λ
Δ
λ
/
λ
=
±
0.1
5
deg.
l
Δ
l
/
l
=
±
0.1
0.24
m
d
Δ
d
/
d
=
±
0.1
0.25
m
J
x
x
(
Δ
J
x
x
)
/
(
J
x
x
=
±
0.1
)
5.5384
kg-m2
J
y
y
(
Δ
J
y
y
)
/
(
J
y
y
=
±
0.1
)
5.6001
kg-m2
J
z
z
(
Δ
J
z
z
)
/
(
J
z
z
=
±
0.1
)
4.2382
kg-m2Table 2
Initial conditions of microsatellite.
Item
Description
Value
Unit
[
φ
,
θ
,
ψ
]
Initial Euler angle
[−120, 30, 60]
deg.
[
q
1, q
2, q
3, q
4
]
Initial quaternion
[0.0474, 0.5303, −0.6597, 0.5303]
[
ω
x, ω
y, ω
z
]
Initial angular rate
[3.5, −3.5, 5.23]
rad/s
d
t
Step time
1.0
sec
T
sim
Simulation time
1200
secIn the next step, we will use the satellite parameters of Tables1 and 2 of the initial conditions to simulate the satellite attitude response.


### 6.2. Experimental Results

In this section, we will implement the SIL simulation using Table1 parameters of microsatellite and initial conditions of microsatellite of Table 2. Figures 8(a) to 8(b), using the Matlab/Simulink simulation results. Then, the SIL simulation results are as shown in Figures 9 to Figure 10.(a) Time response of quaternion. (b) Time response of thruster command.
(a)
(b)(a) Time response ofφ
˙. (b) Time response of θ
˙. (c) Time response of ψ
˙.
(a)
(b)
(c)(a) Time response of thruster 1. (b) Time response of thruster 2. (c) Time response of thruster 3. (d) Time response of thruster 4.
(a)
(b)
(c)
(d)


## 6.1. Parameters of Experimental

The parameters of microsatellite are shown in Table1 and initial conditions of microsatellite are shown in Table 2.Table 1
Parameters of microsatellite.
Item
Uncertainty bound
Value
Unit
λ
Δ
λ
/
λ
=
±
0.1
5
deg.
l
Δ
l
/
l
=
±
0.1
0.24
m
d
Δ
d
/
d
=
±
0.1
0.25
m
J
x
x
(
Δ
J
x
x
)
/
(
J
x
x
=
±
0.1
)
5.5384
kg-m2
J
y
y
(
Δ
J
y
y
)
/
(
J
y
y
=
±
0.1
)
5.6001
kg-m2
J
z
z
(
Δ
J
z
z
)
/
(
J
z
z
=
±
0.1
)
4.2382
kg-m2Table 2
Initial conditions of microsatellite.
Item
Description
Value
Unit
[
φ
,
θ
,
ψ
]
Initial Euler angle
[−120, 30, 60]
deg.
[
q
1, q
2, q
3, q
4
]
Initial quaternion
[0.0474, 0.5303, −0.6597, 0.5303]
[
ω
x, ω
y, ω
z
]
Initial angular rate
[3.5, −3.5, 5.23]
rad/s
d
t
Step time
1.0
sec
T
sim
Simulation time
1200
secIn the next step, we will use the satellite parameters of Tables1 and 2 of the initial conditions to simulate the satellite attitude response.


## 6.2. Experimental Results

In this section, we will implement the SIL simulation using Table1 parameters of microsatellite and initial conditions of microsatellite of Table 2. Figures 8(a) to 8(b), using the Matlab/Simulink simulation results. Then, the SIL simulation results are as shown in Figures 9 to Figure 10.(a) Time response of quaternion. (b) Time response of thruster command.
(a)
(b)(a) Time response ofφ
˙. (b) Time response of θ
˙. (c) Time response of ψ
˙.
(a)
(b)
(c)(a) Time response of thruster 1. (b) Time response of thruster 2. (c) Time response of thruster 3. (d) Time response of thruster 4.
(a)
(b)
(c)
(d)


## 7. Conclusions

In this paper, the satellite attitude determination and control techniques are studied and the verification results are presented. The ADCS is designed with the aid of software and being tested in SIL simulation for the verification of the ADCS model. The results of SIL show that the attitude estimation errors are within 3°and the UKF error can converge after 240 seconds even though the initial attitude error is large. The designed ADCS can fulfill the mission requirement by achieving attitude determination error of less than 5°.In the future, some additional subjects will be investigated. The sensors, magnetometer, and sun sensor should undergo calibration tests in order to increase the precision of the simulated sensor model. The classic Kalman filter (UKF) was superseded by unscented Kalman filter that can be used as the algorithm for attitude estimation.



---

*Source: 904708-2014-05-06.xml*