Generic Model of a Satellite Attitude Control System

Abstract

A generic model of a nanosatellite attitude control and stabilization system was developed on the basis of magnetorquers and reaction wheels, which are controlled by PID controllers with selectable gains. This approach allows using the same architectures of control algorithms (and software) for several satellites and adjusting them to a particular mission by parameter variation. The approach is illustrated by controlling a satellite attitude in three modes of operation: detumbling after separation from the launcher, nominal operation when the satellite attitude is subjected to small or moderate disturbances, and momentum unloading after any reaction wheel saturation. The generic control algorithms adjusted to each mode of operation were implemented in a complete attitude control system. The control system model was embedded into a comprehensive simulation model of satellite flight. The simulation results proved the efficiency of the proposed approach.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

The advancements in miniaturization of sensors and mechanical devices and computing hardware and software lead to enhanced performance of small satellites. The novel concepts of cooperating satellites flying in a formation give prospects for a further increase of the possibility of providing services similar to those offered now by large and expensive spacecraft and also of offering new services or functionalities [1]. An important requirement for satellite operation in all flight phases, crucial for collaboration in formation flying, is effective attitude stabilization and control.

Satellite control torques may be generated in a passive or active way. Passive attitude control does not consume power but may be applied only for dumping of unwanted satellite motion resulting from external disturbances [2]. Active control provides higher control authority than passive methods and allows performing rotational maneuvers by using actuators on board a spacecraft. The most commonly used satellite actuators are magnetorquers and reaction wheels [3–5], sometimes also thrusters [6], but for attitude control, they are rarely used due to fuel mass and volume penalty and its irreversible depletion. When low pointing accuracy is sufficient, magnetorquers may be the only actuators used for attitude control [7, 8]. CubeSats attitude control systems are overviewed in [9]. Most CubeSats use reaction wheels and magnetorquers, although passive magnetic stabilization is also applied, especially in small units.

After separation from a launcher, the satellite may tumble, i.e., spin with a high angular rate, which [10–12] may be moderate, like 0.1 rad/s in each axis, but also may be above 95rad/s (600°/s). In nominal operating conditions, a satellite is subjected to various disturbances, which require quick and effective reactions by the attitude control system, usually using reaction wheels as actuators. Uncontrolled reaction wheels saturate, i.e., they reach the limits of maximal angular velocity, due to large external disturbances or when counteracting disturbing torques accumulate satellite angular velocities for a long time. In this research, three modes of ACS operation are considered: detumbling, nominal, and momentum unloading. In the detumbling mode, usually only magnetorquers are used, as the reaction wheels may quickly saturate. However, in [4, 10, 13], reaction wheels were considered also for CubeSat detumbling. In [14], B-dot and bang-bang B-dot control laws were applied for magnetorquers, and bang-bang B-dot provided quicker stabilization but with a higher final angular velocity error than the B-dot one. In [15], satellite angular velocity feedback was implemented for the detumbling mode.

The nominal control laws provide stabilization counteracting external disturbances, which, even when small (like gravity-gradient) but acting for a long time, may drift a satellite away from its required attitude. Another control task in the nominal mode is satellite rotational maneuvers. As primary actuators, reaction wheels are commonly used in the nominal mode, as they may assure a high attitude accuracy and force a relatively quick response (although slower than thrusters). To control reaction wheels, proportional-derivative (PD) controllers are used [16–24], sometimes with modifications like including additional nonlinear terms [16, 23]. The control laws for the nominal mode are often developed using a simplified satellite model (for instance, only rotational motion is considered neglecting disturbance torques) before being applied to nonlinear systems [16].

In such a case, reaction wheels require diminishing their momentum. For small satellites in low Earth orbits (LEO), magnetorquers are used to unload momentum, as cheap, reliable, and effective torque-generating devices [3]. In [24, 25], the crossproduct control law was applied to control three reaction wheels, and in [3], its modification called a revisited-product control law was investigated. A more general approach was presented in [5] to control a set of three or more reaction wheels.

For each mode of operation, the signals from the control law block should be passed as input to particular actuators. It may be achieved in various ways, for instance, using one type of actuators only (reaction wheels in the nominal mode or magnetorquers in the detumbling mode) and then minimizing the control effort [3, 4].

The main objective of the study was to develop an attitude control system to be implemented on board CubeSats for prospective formation flying. The ease of applicability to similar satellite sizes was the main requirement to diminish the cost of control system development and tuning. This explains why the traditional PID approach was implemented with focus on flexibility of implementation to similar satellites and scalability for a number of actuators. To evaluate feasibility of the developed concept, a comprehensive WUT simulation model of a satellite flight was supplemented by a control system model, allowing evaluation of the control quality. The paper will describe the satellite simulation model, and then next, the generic structure of the control system will be presented. The simulation results will confirm the efficiency of the proposed control system and algorithms.

2. Satellite Simulation Model

The satellite simulation software developed at the Warsaw University of Technology (Figure 1) contains satellite dynamic and kinematic equations, models of sensors and actuators, and algorithms of the navigation and control system. The modular structure allows flexibility in model extension and modifications. The model operates in the time domain.

[figure omitted; refer to PDF]

A satellite bus is modeled as a rigid body with six degrees of freedom and time-varying mass and inertia. The dynamic equations of motion are derived in the Body Frame Fixed (BFF) coordinate system (Figure 2), the kinematic equations are derived in the Earth Centered Inertial (ECI) coordinate system, and the control loads are calculated first in local actuator coordinate systems and then transferred to BFF (Figure 2).

[figure omitted; refer to PDF]

The rigid spacecraft dynamic equations of motion are described as $$\begin{array}{}\text{(1)}& \underset{\text{Inertia}\text{loads}}{\underbrace{\mathbf{A}\left(t\right)\dot{\mathbf{x}}+\mathbf{B}\left(t,\mathbf{x}\right)\mathbf{x}}}=\underset{\text{External}\text{loads}}{\underbrace{{\mathbf{f}}_G+{\mathbf{f}}_A+{\mathbf{f}}_{\text{SRP}}+{\mathbf{f}}_M}}+\underset{\text{Control}\text{loads}}{\underbrace{{\mathbf{f}}_R+{\mathbf{f}}_T+{\mathbf{f}}_Q}},\end{array}$$where $t$ is time, $\mathbf{x}={\left[UVWPQR\right]}^T$ is the satellite state vector composed of components of satellite linear ${\mathbf{v}}_F$ and angular ${\omega }_F$ velocity vectors, and ${\mathbf{f}}_{\left(\right)}={\left[\begin{array}{cc}{\mathbf{F}}_{\left(\right)}^T& {\mathbf{M}}_{\left(\right)}^T\end{array}\right]}^T$ is the load vectors composed of force ${\mathbf{F}}_{\left(\right)}^T$ and moment ${\mathbf{M}}_{\left(\right)}^T$ components. Indices refer to the following: $G$: gravity; $A$: aerodynamics; $\text{SRP}$: solar radiation pressure; $M$: magnetic field (acting on the bus); $R$: reaction wheels (or other momentum exchange devices); $T$: propulsion (not used in this study); $Q$: magnetorquers.

The inertia matrix $\mathbf{A}\left(t\right)$ has the form $$\begin{array}{}\text{(2)}& \mathbf{A}\left(t\right)=\left[\begin{array}{cccccc}m& 0& 0& 0& S_z& -S_y\\ 0& m& 0& -S_z& 0& S_x\\ 0& 0& m& S_y& -S_x& 0\\ 0& -S_z& S_y& I_{xx}& -I_{xy}& -I_{xz}\\ S_z& 0& -S_x& -I_{xy}& I_{yy}& -I_{yz}\\ -S_y& S_x& 0& -I_{xz}& -I_{yz}& I_{zz}\end{array}\right],\end{array}$$where $m$, $S_x$, $S_y$, $S_z$, $I_{xx}$, $I_{yy}$, $I_{zz}$, $I_{xy}$, $I_{yz}$, and $I_{xz}$ are the satellite mass, static mass moments, and moments and products of inertia, respectively.

The $\mathbf{B}\left(t,\mathbf{x}\right)$ matrix is calculated as $$\begin{array}{}\text{(3)}& \mathbf{B}\left(t,\mathbf{x}\right)=\mathbf{O}\left(\mathbf{x}\right)\mathbf{A}\left(t\right),\end{array}$$where $$\begin{array}{}\text{(4)}& \mathbf{O}\left(\mathbf{x}\right)=\left[\begin{array}{cccccc}0& -R& Q& 0& 0& 0\\ R& 0& -P& 0& 0& 0\\ -Q& P& 0& 0& 0& 0\\ 0& -W& V& 0& -R& Q\\ W& 0& -U& R& 0& -P\\ -V& U& 0& -Q& P& 0\end{array}\right].\end{array}$$

The satellite kinematic equations for position and attitude have the form $$\begin{array}{c}\text{(5)}& {\dot{\mathbf{r}}}_I={\mathbf{C}}_F^I{\mathbf{v}}_F,{\dot{\mathbf{q}}}_I=\frac{1}{2}{\mathbf{Q}\mathbf{q}}_I,\end{array}$$where ${\mathbf{r}}_{FI}={\left[\begin{array}{ccc}{x}_{FI}& y_{FI}& z_{FI}\end{array}\right]}^T$ is the satellite position in the inertial coordinate system (ECI), ${\mathbf{v}}_F={\left[\begin{array}{ccc}U& V& W\end{array}\right]}^T$ is the satellite linear velocity in the satellite body frame, ${\omega }_F={\left[\begin{array}{ccc}P& Q& R\end{array}\right]}^T$ is the satellite angular velocity in the satellite body frame, ${\mathbf{q}}_I={\left[\begin{array}{cccc}{q}_0& q_1& q_2& q_3\end{array}\right]}^T$ is the quaternion of a satellite attitude in the inertial frame, and ${\mathbf{C}}_F^I$ is the matrix of transformation from the satellite body to inertial frame in the form $$\begin{array}{}\text{(6)}& {\mathbf{C}}_F^I={\mathbf{C}}_F^I\left(\mathbf{q}\right)=\left[\begin{array}{ccc}{q}_0^2+q_1^2-q_2^2-q_3^2& 2\left(q_1{q}_2-q_3{q}_0\right)& 2\left(q_1{q}_3+q_2{q}_0\right)\\ 2\left(q_1{q}_2+q_3{q}_0\right)& q_0^2-q_1^2+q_2^2-q_3^2& 2\left(q_2{q}_3-q_1{q}_0\right)\\ 2\left(q_1{q}_3-q_2{q}_0\right)& 2\left(q_2{q}_3+q_1{q}_0\right)& q_0^2-q_1^2-q_2^2+q_3^2\end{array}\right].\end{array}$$

Matrix $\mathbf{Q}$ is composed of satellite angular velocities and has the form $$\begin{array}{}\text{(7)}& \mathbf{Q}=\left[\begin{array}{cccc}0& -P& -Q& -R\\ P& 0& R& -Q\\ Q& -R& 0& P\\ R& Q& -P& 0\end{array}\right].\end{array}$$

It is assumed that on a satellite, $N$ control actuators are mounted, producing torques ${\mathbf{M}}_{Ci}$ only, so the total control torque is $$\begin{array}{}\text{(8)}& {\mathbf{f}}_C=\sum _{i=1}^N{\mathbf{f}}_{Ci}-={\left[\begin{array}{cc}\mathbf{0}& \sum _{i=1}^N{\mathbf{M}}_{Ci}\end{array}\right]}^T\hspace{1em}\left(i=1,\cdots ,N\right),\end{array}$$where ${\mathbf{M}}_{Ci}$ is the torque generated by the $i$-th actuator (reaction wheel or magnetorquer) in the BFF coordinate system.

A reaction wheel is modelled as a flywheel fixed to the satellite body spinning with angular velocity ${\Omega }_R$ and moment of inertia $I_{Rx}$. Its angular momentum in the actuator coordinate system is $$\begin{array}{}\text{(9)}& {\mathbf{H}}_R={\left[\begin{array}{ccc}{I}_{Rx}\left({\Omega }_R+{\omega }_{Rx}\right)& 0& 0\end{array}\right]}^T,\end{array}$$where ${\mathbf{H}}_R$ is the angular momentum of the reaction wheel, $I_{Rx}$ is the moment of inertia of the reaction wheel about its spin axis, ${\Omega }_R$ is the spin velocity of the reaction wheel, and ${\omega }_R={\left[\begin{array}{ccc}{\omega }_{Rx}& {\omega }_{Ry}& {\omega }_{Rz}\end{array}\right]}^T$ is the satellite angular velocity in the actuator coordinate system.

A torque acting on the satellite due to a single reaction wheel (in the BFF frame) is calculated as $$\begin{array}{}\text{(10)}& {\mathbf{M}}_R={\left[\begin{array}{ccc}-I_{Rx}\left({\dot{\Omega }}_R+{\dot{\omega }}_{Rx}\right)& -I_{Rx}{\omega }_{Rz}\left({\Omega }_R+{\omega }_{Rx}\right)& I_{Rx}{\omega }_{Ry}\left({\Omega }_R+{\omega }_{Rx}\right)\end{array}\right]}^T.\end{array}$$

Angular velocity ${\Omega }_R$ of a reaction wheel is the control variable for the generated torque.

The reaction wheel is driven by an electric motor, the dynamics of which in general is described by combined electromechanical equations. As control signals, COTS (commercial off-the-shelf) reaction wheels use current, torque, or angular velocity, and their control electronics converts them to the input voltage of the motor. In this study, it is assumed that the armature current $i_a$ is a control variable, and the reaction wheel model is described by mechanical state equations in the form $$\begin{array}{}\text{(11)}& I_{Rx}\left({\dot{\Omega }}_R+{\dot{\omega }}_{Rx}\right)=k_R{i}_R-b_R{\Omega }_R,\end{array}$$where $i_R$ is the armature current, $k_R$ is the motor torque constant, and $b_R$ is the viscous friction coefficient of the motor shaft.

The torque from a magnetorquer is a crossproduct of the magnetic dipole generated by a coil and the Earth magnetic field vector: $$\begin{array}{}\text{(12)}& {\mathbf{M}}_Q={\mathbf{D}}_Q\times {\mathbf{B}}_Q={\left[\begin{array}{ccc}0& -D_Q{B}_{Qz}& D_Q{B}_{Qy}\end{array}\right]}^T,\end{array}$$where ${\mathbf{D}}_Q={\left[\begin{array}{ccc}{D}_Q& 0& 0\end{array}\right]}^T$ is the magnetic dipole generated by the coil and ${\mathbf{B}}_Q={\left[\begin{array}{ccc}{B}_{Qx}& B_{Qy}& B_{Qz}\end{array}\right]}^T$ is a vector of the Earth magnetic field in the magnetorquer coordinate system.

A magnetic dipole is calculated as $$\begin{array}{}\text{(13)}& D_Q=D_0{u}_Q,\end{array}$$where $D_0$ is the nominal dipole moment of the magnetorquer at a given nominal voltage and $u_Q$ is the control variable of a magnetorquer.

3. Attitude Control System

The general architecture of a satellite attitude control system investigated in this study is shown in Figure 3. The blocks “actuators” and “satellite” (Figure 3) are the models described in the previous chapter used for simulation of bus and actuators dynamics. “Sensors” are not considered in this study. The block “current state” data contains state variables of a satellite, reaction wheels and magnetorquers, quaternion ${\mathbf{q}}_I={\left[\begin{array}{cccc}{q}_0& q_1& q_2& q_3\end{array}\right]}^T$ representing satellite attitude in the inertial frame, vector ${\omega }_F={\left[\begin{array}{ccc}P& Q& R\end{array}\right]}^T$ of a satellite angular velocity expressed in the satellite body coordinate system, vector ${\Omega }_{RW}$ of reaction wheels angular velocities, and vector ${\mathbf{u}}_M$ of the magnetorquer control signal. The required values of these quantities, denoted with index $\text{cmd}$, are taken from the block named “flight management.”

[figure omitted; refer to PDF]

The differences between the required and current states are inputs to the satellite attitude control system (ACS), composed of two blocks denominated in Figure 3 as “satellite attitude control” and “reaction wheel desaturation.” The task of the satellite attitude control block is to bring a satellite to the required attitude state and/or to keep it at that state. The task of the reaction wheel desaturation block is to decelerate them to angular velocities determined by vector ${\Omega }_{\text{cmd}}$. The desaturation process should not influence significantly the satellite rotational dynamics, so it is included as a separate line of control.

For both control blocks, full PID control is envisaged and implemented in the algorithm. $$\begin{array}{}\text{(14)}& {\mathbf{u}}_{\left(\ast \right)}={\mathbf{K}}_{P\left(\ast \right)}{\sigma }_{\left(\ast \right)}+{\mathbf{K}}_{D\left(\ast \right)}{\dot{\sigma }}_{\left(\ast \right)}+{\mathbf{K}}_{I\left(\ast \right)}\int {\sigma }_{\left(\ast \right)}dt,\end{array}$$where $\sigma$ is the vector of the differences between the required and current values; ${\mathbf{K}}_{P\left(\ast \right)}$, ${\mathbf{K}}_{D\left(\ast \right)}$, and ${\mathbf{K}}_{I\left(\ast \right)}$ are the proportional, derivative, and integral gain matrices, respectively; and $\left(\ast \right)$ is the index of the actuator used.

The applied PID control is illustrated in Figure 4 showing the selection of actuators and control gain matrix for three system operating modes: detumbling, nominal, and momentum unloading.

[figure omitted; refer to PDF]

3.1. Detumbling Mode

The detumbling mode is activated after separation of the satellite from the launch vehicle or at any time when the satellite’s angular velocity is too large. The satellite angular velocity is reduced by the attitude controller, and the momentum management controller is inactive. In the detumbling mode, only magnetorquers are used for spacecraft stabilization, and PID is reduced to differential control related to the actual value of satellite angular velocity $$\begin{array}{}\text{(15)}& {\mathbf{u}}_Q={\mathbf{K}}_{DQ}{\omega }_e={\mathbf{K}}_{DQ}\omega ,\end{array}$$as the required satellite velocity reduction ${\omega }_{\text{cmd}}=0$.

The gain matrix ${\mathbf{K}}_{DQ}$ depends on the magnetorquer arrangement, magnetic field vector, satellite inertia, and orbital period (similar to [15]) and is calculated as $$\begin{array}{}\text{(16)}& {\mathbf{K}}_{DQ}=\left(\frac{4\pi }{T_{\text{orb}}}\mathbf{I}\right)\frac{{\mathbf{D}}_0^{-1}{\mathbf{G}}_{\text{MTQ}}^{+}{\mathbf{B}}_F^{\times }}{{\left|{\mathbf{B}}_F\right|}^2},\end{array}$$where $N_Q$ is the number of magnetorquers, ${\mathbf{D}}_0=\mathrm{diag}\left(\begin{array}{ccc}{D}_{01}& \cdots & D_{0N_Q}\end{array}\right)$ is the diagonal matrix composed of nominal dipole moments $D_{0i}$ of magnetorquers, ${\mathbf{G}}_{\text{MTQ}}$ is the matrix ($3\times N_Q$) representing the configuration of magnetorquers in the satellite body frame (the “+” in the superscript indicates a pseudoinverse or an inverse of matrix ${\mathbf{G}}_{\text{MTQ}}$ if it exists), ${\mathbf{B}}_F$ is the magnetic field vector ($3\times 1$) in BFF (the “×” in the superscript indicates a skew symmetric matrix ($3\times 3$) satisfying ${\mathbf{B}}_F^{\times }\omega ={\mathbf{B}}_F\times \omega$), $T_{\text{orb}}$ is the orbital period of the satellite, and $\mathbf{I}$ is the matrix of the satellite moment of inertia matrix about its center of mass.

When the specified low angular velocity is achieved, the attitude control system is switched from the detumbling mode to the more precise nominal mode.

3.2. Nominal Mode

The nominal mode is used to maintain or change satellite attitude as required with adequate accuracy, despite external disturbances. The control is performed by reaction wheels in a full PID form $$\begin{array}{}\text{(17)}& {\mathbf{u}}_R={\mathbf{K}}_{PR}{\sigma }_R+{\mathbf{K}}_{DR}{\omega }_e+{\mathbf{K}}_{IR}\int {\sigma }_{R}dt.\end{array}$$

In (17), the attitude error vector $\sigma$ is calculated as a rotation from the actual attitude ${\mathbf{q}}_I={\left[\begin{array}{cccc}{q}_0& q_1& q_2& q_3\end{array}\right]}^T$ to the commanded attitude ${\mathbf{q}}_{I\text{cmd}}={\left[\begin{array}{cccc}{q}_{0c}& q_{1c}& q_{2c}& q_{3c}\end{array}\right]}^T$ in a quaternion notation as [24] $$\begin{array}{c}\text{(18)}& \left[\begin{array}{c}{q}_{0e}\\ q_{1e}\\ q_{2e}\\ q_{3e}\end{array}\right]=\left[\begin{array}{cccc}{q}_{0c}& q_{1c}& q_{2c}& q_{3c}\\ -q_{1c}& q_{0c}& q_{3c}& -q_{2c}\\ -q_{2c}& -q_{3c}& q_{0c}& q_{1c}\\ -q_{3c}& q_{2c}& -q_{1c}& q_{0c}\end{array}\right]\left[\begin{array}{c}{q}_0\\ q_1\\ q_2\\ q_3\end{array}\right],{\sigma }_R=2q_{0e}{\left[\begin{array}{ccc}{q}_{1e}& q_{2e}& q_{3e}\end{array}\right]}^T.\end{array}$$

The derivative of attitude error ${\dot{\sigma }}_R$ is the difference of satellite angular velocities: $$\begin{array}{}\text{(19)}& {\omega }_{Fe}={\omega }_{F\text{cmd}}-{\omega }_F.\end{array}$$

The gain matrices ${\mathbf{K}}_{PR}$, ${\mathbf{K}}_{DR}$, and ${\mathbf{K}}_{IR}$ are calculated as $$\begin{array}{}\text{(20)}& {\mathbf{K}}_{PR}=\left({\omega }_n^2+\frac{2{\zeta \omega }_n}{T}\right){\mathbf{R}}^{-1}{\mathbf{G}}_{RW}^{+}\mathbf{I},\text{(21)}& {\mathbf{K}}_{DR}=\left(2{\zeta \omega }_n+\frac{1}{T}\right){\mathbf{R}}^{-1}{\mathbf{G}}_{RW}^{+}\mathbf{I},\text{(22)}& {\mathbf{K}}_{IR}=\frac{{\omega }_n^2}{T}{\mathbf{R}}^{-1}{\mathbf{G}}_{RW}^{+}\mathbf{I},\end{array}$$where $N_R$ is the number of reaction wheels; $\mathbf{R}=\mathrm{diag}\left(\begin{array}{ccc}{k}_{R1}& \cdots & k_{{RN}_R}\end{array}\right)$ is the diagonal matrix ($N_R\times N_R$) containing reaction wheel motor constants; ${\mathbf{G}}_{RW}$ is the matrix ($3\times N_R$) representing the configuration of reaction wheels, relative to satellite body axes; and ${\omega }_n$, $\zeta$, and $T$ are the required linear control bandwidth, dumping ratio, and time constant of the integral controller, respectively.

3.3. Momentum Unloading Mode

The momentum unloading mode is used to prevent reaction wheel saturation by reducing their angular velocities. In the momentum unloading mode, reaction wheels are used for satellite attitude control, and simultaneously, magnetorquers are switched on for reaction wheel desaturation. The mode is activated when the velocity of at least one reaction wheel exceeds its design limit. However, for better mission performance, the momentum unloading system may be switched on from time to time to reduce accumulating wheel rotating velocity.

In the momentum unloading mode, the control of reaction wheels should assure satellite attitude control and diminish the spin of reaction wheels, so the total reaction wheel control ${\mathbf{u}}_R$ in this mode is done as $$\begin{array}{}\text{(23)}& {\mathbf{u}}_R={\mathbf{u}}_{RS}+{\mathbf{u}}_{RR},\end{array}$$where ${\mathbf{u}}_{RS}$ is applied to perform the first task and ${\mathbf{u}}_{RR}$ the second one.

As the integral term in PID might slow down the control reaction, for satellite attitude control ${\mathbf{u}}_{RS}$ in the momentum unloading mode, the PD controller for reaction wheels is used: $$\begin{array}{}\text{(24)}& {\mathbf{u}}_{RS}={\mathbf{K}}_{PRS}{\sigma }_R+{\mathbf{K}}_{DRS}{\omega }_e,\end{array}$$with modified $P$ and $D$ gain matrices (eq. (20) and (21)) calculated here as $$\begin{array}{}\text{(25)}& {\mathbf{K}}_{PRS}={\omega }_n^2{\mathbf{R}}^{-1}{\mathbf{G}}_{RW}^{+}\mathbf{I},\text{(26)}& {\mathbf{K}}_{DRS}=2{\zeta \omega }_n{\mathbf{R}}^{-1}{\mathbf{G}}_{RW}^{+}\mathbf{I}.\end{array}$$

To decelerate angular velocities of reaction wheels to a specified value ${\Omega }_{\text{cmd}}$, ${\mathbf{u}}_{RR}$ is calculated as $$\begin{array}{}\text{(27)}& {\mathbf{u}}_{RR}={\mathbf{K}}_{DRR}{\Omega }_e,\end{array}$$where ${\Omega }_e={\Omega }_{\text{cmd}}-{\Omega }_{RW}$ is the reaction wheel angular velocity error vector, and matrix ${\mathbf{K}}_{RR}$ is calculated as $$\begin{array}{}\text{(28)}& {\mathbf{K}}_{RR}=k_R{\mathbf{R}}^{-1}{\mathbf{I}}_{RW},\end{array}$$where $k_R$ is the scalar control gain.

The magnetorquers are also controlled by the D controller: $$\begin{array}{}\text{(29)}& {\mathbf{u}}_Q={\mathbf{K}}_{DQ}{\Omega }_e,\end{array}$$and gain ${\mathbf{K}}_{DQ}$ is calculated as [3, 24, 25] $$\begin{array}{}\text{(30)}& {\mathbf{K}}_{DQ}=\frac{k_Q}{{\left|{\mathbf{B}}_F\right|}^2}{\mathbf{D}}_0^{-1}{\mathbf{G}}_{MTQ}^{+}{\mathbf{B}}_F^{\times }{\mathbf{G}}_{MTQ}{\mathbf{I}}_{RW},\end{array}$$where $k_Q$ is the scalar control gain.

Proper selection of control gains $k_Q$ and $k_R$ in the momentum unloading mode results in lower disturbances of spacecraft attitude and shorter duration of the reaction wheel desaturation process.

4. Simulation Study

To investigate the efficiency of the developed attitude control system, a simulation study was done for the generic satellite mass data given in Table 1.

Table 1

Satellite mass properties.

Description	Value	Unit
Mass	6.2	[kg]
Moments of inertia: $I_{xx}$, $I_{yy}$, $I_{zz}$	0.0756, 0.0763, 0.0209	[kg·m2]
Products of inertia: $I_{xy}$, $I_{yz}$, $I_{xz}$	-0.0002, 0.0020, -0.0019	[kg·m2]
Location of center of mass relative to geometrical center ($x$, $y$, $z$)	0.0034, -0.0036, -0.0007	[m]
Dimensions	$$0.7\times 0.7\times 0.4$$	[m]

The initial satellite orbit elements are provided in Table 2.

Table 2

Initial satellite orbital parameters.

Description	Value	Unit
Semimajor axis	6371 000	[m]
Eccentricity	0	[-]
Longitude of ascending node	331.7739	[°]
Orbit inclination	51.6393	[°]
Argument of perigee	0	[°]
True anomaly	20.5285	[°]
Orbital period	5500	s

Three configurations of reaction wheels: orthogonal, pyramidal, and tetrahedral (Figure 5), were investigated accompanied by a set of three orthogonal magnetorquers.

[figures omitted; refer to PDF]

The reaction wheel parameters used in simulations are based on Sinclair RW-0.03 data [26] and are given in Table 3.

Table 3

RW data.

Description	Value	Unit
Moment of inertia about spin axis	51.16	[kg·mm2]
Armature current	0.35	[A]
Viscous friction	3.837	[kg·mm2/s]
Motor torque constant	51.6393	[°]

The major parts of simulations were performed for the orthogonal configuration of reaction wheels as the most commonly used one [9]. The other two configurations were investigated mainly for validating the scalability of the developed algorithms.

The matrices ${\mathbf{G}}_{RW}$ describing reaction wheel configuration, used in equations (20), (21), (22) and (25), (26), (27) are

Orthogonal configuration: identity matrix dimension ($3\times 3$)

Pyramid configuration: $$\begin{array}{}\text{(31)}& {\mathbf{G}}_{RW}=\left[\begin{array}{cccc}0.5774& -0.5774& -0.5774& 0.5774\\ 0.5774& 0.5774& -0.5774& -0.5774\\ 0.5774& 0.5774& 0.5774& 0.5774\end{array}\right]\end{array}$$

Tetrahedron configuration: $$\begin{array}{}\text{(32)}& {\mathbf{G}}_{RW}=\left[\begin{array}{cccc}0.9428& -0.4714& -0.4714& 0\\ 0& 0.8165& -0.8165& 0\\ -0.3333& -0.3333& -0.3333& 1\end{array}\right].\end{array}$$

The set of magnetorquers contains three devices in an orthogonal configuration with nominal dipoles as in [27]: in the $O_F{x}_F$ and $O_F{y}_F$ axes, $D_{01}=D_{02}=0.2\text{A}{\text{m}}^2$, and in the $O_F{z}_F$ axis, $D_{03}=0.24\text{A}{\text{m}}^2$.

For the detumbling mode, the constant part of the gain matrix was calculated purely based on the satellite inertia and orbital parameters as per equation (16). The selective parameters for the nominal (equations (20), (21), (22)) and momentum unloading (equations (28) and (30)) modes were chosen using the trial and error method. The attitude controller parameters were set to ${\omega }_n=0.8$ rad/s, $\zeta =1$, and $T=12.5$ s. For the nominal mode, the aim was to allow relatively quick maneuvers without significant overshoot of the attitude and sufficiently high pointing accuracy for the initial wheel rates of 0 rpm/s. The gains in the momentum unloading mode, $k_Q=k_R=0.05$, were chosen to have relatively small values in order to limit the deviation from the required attitude during that phase.

Initial attitudes for all simulated cases of yaw, pitch, and roll angle values were zero, i.e., ${\Psi }_0={\Theta }_0={\Phi }_0=0$.

The simulation rationale for all cases was to check the applicability of the algorithms developed. In particular, for the detumbling mode, the ability to suppress the satellite rotation was verified. For the normal operating mode, the ability to maintain or change the satellite attitude was investigated, for various angular velocities and configurations of reaction wheels. For the momentum unloading case, the ability to desaturate high rates of reaction wheels was verified.

4.1. Detumbling Mode

The detumbling mode of the ACS operation was simulated to verify the effectiveness of reducing the satellite’s high value of 10°/s per axis of initial angular velocity, which may appear after separation from the launcher. The magnetorquers were actuators used in this mode. The time profiles of the satellite angular velocity module and magnetorquer dipole moments module are presented in Figures 6(a) and 6(b).

[figures omitted; refer to PDF]

The satellite angular velocity was reduced below 0.5°/s after 3200 s and below 0.2°/s after 6000 s, i.e., in time about 9% longer than the orbital period which was 5500 s. The result is satisfying, in terms, for instance, of the requirement formulated in [10] that a satellite should detumble in about one orbital period from the initial angular velocity of 5.73°/s per axis.

As indicated by other simulations, using the detumbling mode controller did not reduce the satellite rate below 0.2°/s. To cross this threshold, the ACS nominal operation mode should be used.

4.2. Nominal Mode

In the simulations of the nominal mode, only reaction wheels were active. The quantities monitored for this mode are satellite Euler angles: $\Psi$ (yaw), $\Theta$ (pitch), and $\Phi$ (roll), angular rates of reaction wheels, and optionally satellite angular velocity. The list of scenarios simulated for the nominal mode of ACS operation is given in Table 4.

Table 4

Nominal mode simulation cases.

Case	Initial satellite rates $P=Q=R$	Commanded attitude(Initial)(After 100 s)	Initial RW spin	RW configuration
1	0.5°/s	0°0°	0	Orthogonal 3 wheels
2	0	0°, 0°, 0°120°, 60°, 150	0, 2000, 4000 rpm	Orthogonal 3 wheels
3	0	0°, 0°, 0°Nadir pointing	0 rpm	Orthogonal 3 wheels
4	0	120°, 60°, 150°	0 rpm	Pyramid 4 wheels
				Tetrahedron 4 wheels

The simulation of case 1 illustrates the transition from detumbling to nominal mode. The initial angular velocity vector is assumed to be 0.5°/s per axis, which is the residual of the detumbling mode. The requirement was to keep the zero Euler angles, i.e., $\Psi =0^{°}$, $\Theta =0^{°}$, and $\Phi =0^{°}$.

At the beginning of the control, there was an increase in the reaction wheel angular rate, but after 10 seconds (Figure 7), the satellite angular rate was reduced from 0.5°/s to almost 0°/s in each axis, which was the result of reaction wheel activity. At that moment, the velocities of the first and the second reaction wheels were stabilized at 120 rpm. The third wheel stabilized at 30 rpm, because the satellite moment of inertia about the respective axis was lower than those on $X$ and $Y$ axes which are approximately equal.

[figures omitted; refer to PDF]

The simulation case 2 was performed to investigate attitude control system performance for different initial reaction wheel rates. In this case, the initial attitude was $\Psi =0^{°}$, $\Theta =0^{°}$, and $\Phi =0^{°}$ and initial angular velocity vector was 0°/s per axis. The satellite was required to keep this status for 100 s and then change attitude to $\Psi ={120}^{°}$, $\Theta ={60}^{°}$, and $\Phi ={150}^{°}$. The simulation was performed for a set of 3 orthogonal reaction wheels, for 3 different initial velocities of reaction wheels.

In case 2 for initial RW rates 0 and 2000 rpm, the performance is acceptable. The high initial reaction wheel rate such as 4000 rpm does not allow achieving a large commanded angle of rotation, which may be the result of reaction wheel saturation, as the second reaction wheel reached the maximal angular velocity of 6200 rpm (Figure 8). The remedy is that for the required large satellite attitude variation, the maneuver should be divided into a few smaller attitude changes to avoid high angular velocity of reaction wheels. Another solution could be to use different values of bandwidth and dumping in the gain matrices for this controller.

[figure omitted; refer to PDF]

The case 3 simulation was performed to check the ACS performance for the cases similar for a time-varying nadir commanded attitude. In the first 100 seconds, the commanded attitude was $\Psi =0^{°}$, $\Theta =0^{°}$, and $\Phi =0^{°}$. After that, the commanded satellite attitude was varying for nadir pointing. The simulation was performed for 3 different sets of initial reaction wheel velocities (0, 2000, and 4000 rpm). The results are presented in Figures 9(a) and 9(b).

[figures omitted; refer to PDF]

In this case, the control system response is similar to that in case 2. The system has a better performance for lower initial reaction wheel rates, as the commanded attitude was reached more quickly and with better accuracy. For initial rates of 4000 rpm, the commanded attitude is not achieved and a significant pointing error occurs. It may be a symptom of reaction wheel saturation, as a result of summing up the angular velocities of reaction wheels and satellite, which makes the controller ineffective.

The simulation case 4 was performed to verify whether the developed control laws may be also used for reaction wheel configurations other than the orthogonal one. Initial satellite angular velocity vector was 0°/s per axis. The pyramid and tetrahedron configurations of four reaction wheels were investigated. In control gain calculations, the configurations differ by matrices ${\mathbf{G}}_{RW}$. In these cases, the satellite was expected to change attitude from initial $\Psi =0^{°}$, $\Theta =0^{°}$, and $\Phi =0^{°}$ to the commanded one: $\Psi ={120}^{°}$, $\Theta ={60}^{°}$, and $\Phi ={150}^{°}$. The results are presented in Figures 10(a) and 10(b).

[figures omitted; refer to PDF]

For both configurations, the maneuver was completed after about 15 seconds. The satellite motion was different for each configuration, but the final results were similar. The overshoot for pyramid configuration was larger on the pitch axis than for the tetrahedron configuration, but on the roll axis, an opposite tendency occurred.

The systems using 4 reaction wheels reach the commanded attitude significantly earlier with smaller overshoot than the one with 3 reaction wheels (case no. 3), as using more reaction wheels allows to generate a higher torque on the satellite. The performance of reaction wheel arrays is different; however, in both cases, the maximal speed achieved by any reaction wheel was similar (above 2000 rpm for the 1^st wheel of the pyramid configuration and for the 2^nd wheel of the tetrahedron configuration).

The results presented in Figure 10(a) illustrate that for different configurations, the performance of the control system may be similar. This poses a difficult task for the designer which configuration should be chosen if 4 reaction wheels are used. A wider set of scenarios would indicate more differences between those 2 configurations. However, a detailed comparison between different reaction wheel configurations is not covered in this study.

4.3. Momentum Unloading Mode

In this mode of ACS operation, both reaction wheels and magnetorquers are active. To verify the effectiveness of diminishing reaction wheel angular velocities, a desaturation process was performed for initial velocities of reaction wheels equal to 6000 rpm which is near their maximal value of operation—6200 rpm. Although it is a rather unrealistic value in practice, this case was studied to check the efficiency of the control. Simulations were done for orthogonal and pyramid configurations of reaction wheels.

For the orthogonal configuration of three reaction wheels (Figures 11(a) and 11(b)), their spin velocities were reduced to near zero values after 6400 s, i.e., 115% of the orbital period, which is a satisfying result. During the control, the rates of reaction wheels were not decreasing monotonically. It may be attributed to the fact that magnetorquers were operating close to their saturation region for most of the time. The satellite pointing error was decreasing with decreasing rates of reaction wheels. It may be explained by viscous dumping of reaction wheel rotation (providing residual torque) which was not compensated for in the momentum unloading mode.

[figures omitted; refer to PDF]

To verify whether it is possible to apply control laws developed during this study to more than three reaction wheels, the simulation for a set of 4 reaction wheels in the pyramid configuration was performed.

The desaturation process was successful for this case too (Figures 12(a) and 12(b)), as angular rates of all four reaction wheels were reduced to near zero. It took more time compared to the operation of 3 orthogonal reaction wheels, most likely due to the fact that a system of 4 reaction wheels has a higher total angular momentum to be reduced. This might also result in a higher initial pointing error.

[figures omitted; refer to PDF]

4.4. Simulation Result Summary

The objective of the simulation study was to verify whether the attitude control system developed on the basis of linear control may be applied to the nonlinear satellite model. Three modes of operations were simulated, and the results of the control are satisfactory.

In the detumbling mode, the satellite angular velocity was reduced even for very high initial angular velocities, such as 10°/s per axis. The other simulation results indicated that using only magnetorquers for high accuracy pointing in the presence of external torque acting on the satellite may not be effective, which suggests that reaction wheels should be used in the final tuning of satellite attitude.

The simulations for the nominal mode of operation indicated that the control algorithm, applied in this mode, provides sufficient attitude control. The control algorithm developed in this study may be applied to various configurations of reaction wheels, which was illustrated by simulation results for orthogonal, pyramid, and tetrahedron reaction wheel configurations. The study revealed limits of using reaction wheels at high angular rates. As varying external loads acting on a satellite results in increase of reaction wheel angular velocities, the momentum unloading mode is necessary to maintain the required performance of reaction wheel actuation during a mission.

The simulations of the momentum unloading mode confirmed that the developed control algorithm diminishes angular velocities of reaction wheels in their various configurations. The desaturation process has a significant impact on satellite attitude control and indicates that in the momentum unloading mode, a satellite should not be commanded to perform tasks which require high accuracy pointing.

5. Conclusions

An attitude control system was developed which fulfils two main tasks: provides a satellite with sufficient attitude control capabilities in the detumbling and normal modes of operation and ensures adequate performance of control actuators by the momentum unloading control process. The actuators of the control system are magnetorquers and reaction wheels. The software architecture allows scalability of the system to various numbers and configurations of actuators as well as changing PID regulator parameters, adjusting to a particular case. The model of the attitude control system was integrated with the satellite simulation software developed at the Warsaw University of Technology. The results of simulations of the control system operation are satisfying. The model may be a useful tool for development of an attitude control system to be operating in real time and embedded into nanosatellite hardware. Also, the simulation results presented in the paper may be useful for similar system development. In further research, the software is intended to be used for other missions and other control strategies.

Acknowledgments

This research was supported by the European Commission’s Horizon 2020 Program [grant number 687490].

Glossary

Nomenclature

Vectors and matrices:Bold

$\mathbf{A}\left(t\right)$, $\mathbf{B}\left(t,\mathbf{x}\right)$:Satellite inertia matrices

$b_R$:Viscous friction coefficient of the RW motor shaft

${\mathbf{B}}_F$:Magnetic field vector

${\mathbf{B}}_Q$:The Earth magnetic field in a magnetorquer coordinate system ${\mathbf{B}}_Q={\left[\begin{array}{ccc}{B}_{xQ}& B_{yQ}& B_{zQ}\end{array}\right]}^T$

${\mathbf{C}}_A^B$:Transformation matrix from coordinate system $A$ to system $B$

$D_0$:Nominal magnitude of magnetic dipole generated by one magnetorquer

${\mathbf{D}}_Q$:Magnetic dipole generated by the magnetorquer coil, ${\mathbf{D}}_Q={\left[\begin{array}{ccc}{D}_Q& 0& 0\end{array}\right]}^T$

${\mathbf{D}}_0$:Diagonal matrix composed of magnetorquers nominal dipole moments ${\mathbf{D}}_0=\mathrm{diag}\left(\begin{array}{ccc}{D}_{01}& \cdots & D_{0N_Q}\end{array}\right)$

${\mathbf{F}}_{\left(\ast \right)}$:Force generated by ($\ast$)

${\mathbf{f}}_{\left(\ast \right)}$:Satellite loads ${\mathbf{f}}_{\left(\ast \right)}={\left[\begin{array}{cc}{\mathbf{F}}_{\left(\ast \right)}& {\mathbf{M}}_{\left(\ast \right)}\end{array}\right]}^T$

${\mathbf{G}}_{\left(\ast \right)}$:Matrix representing actuator configuration

${\mathbf{H}}_R$:Angular momentum of the reaction wheel

$i_R$:Armature current of the reaction wheel motor

$I_{xx},I_{yy},I_{zz}$:Satellite moments of inertia

$I_{xy},I_{yz},I_{xz}$:Satellite products of inertia

$I_{Rx}$:Moment of inertia of the reaction wheel about its spin axis

$\mathbf{I}$:Satellite moments of inertia matrix

${\mathbf{I}}_{\text{RW}}$:Diagonal matrix containing the moments of inertia of reaction wheels

$k_R$:Torque constant of the RW motor

${\mathbf{K}}_{P\left(\ast \right)}{\mathbf{K}}_{D\left(\ast \right)}{\mathbf{K}}_{I\left(\ast \right)}$:Gain matrices: proportional, derivative, and integral

${\mathbf{M}}_{\left(\ast \right)}$:Torque

$m$:Satellite mass

$N$:Total number of actuators

$N_{\left(\ast \right)}$:Number of devices

${\mathbf{q}}_I$:Quaternion of a satellite attitude in the inertial frame ${\mathbf{q}}_I={\left[\begin{array}{cccc}{q}_0& q_1& q_2& q_3\end{array}\right]}^T$

$\mathbf{Q}$:Matrix of angular velocities in attitude equation

${\mathbf{r}}_I$:Satellite position in the inertial coordinate system (ECI), ${\mathbf{r}}_I={\left[\begin{array}{ccc}{x}_I& y_I& z_I\end{array}\right]}^T$

$R$:Satellite yaw rate

$\mathbf{R}$:Diagonal matrix containing reaction wheel motor constants, $\mathbf{R}=\mathrm{diag}\left(\begin{array}{ccc}{k}_{R1}& \cdots & k_{{RN}_R}\end{array}\right)$

$S_x$, $S_y$, $S_z$:Static mass moments of the satellite

$t$:Time

$T$:Time constant of the integral controller

$T_{\text{orb}}$:Satellite orbital period

${\mathbf{u}}_{\left(\ast \right)}$:Actuator control signals

${\mathbf{v}}_F$:Satellite linear velocity in the satellite body frame, ${\mathbf{v}}_F={\left[\begin{array}{ccc}U& V& W\end{array}\right]}^T$

$\mathbf{x}$:Satellite state vector

$\zeta$:Specified control dumping ratio

$\sigma$:Attitude error vector

${\omega }_n$:Specified linear control bandwidth

${\omega }_F$:Satellite angular velocity in the satellite body frame, ${\omega }_F={\left[\begin{array}{ccc}P& Q& R\end{array}\right]}^T$

${\omega }_R$:Satellite angular velocity in the actuator coordinate system ${\omega }_R={\left[\begin{array}{ccc}{\omega }_{Rx}& {\omega }_{Ry}& {\omega }_{Rz}\end{array}\right]}^T$

${\Omega }_R$:Spin velocity of a single reaction wheel

$\mathbf{O}$:Satellite linear and angular velocity matrix

${\Omega }_{\left(\ast \right)}$:Reaction wheel angular velocities.

($\ast$):Reflects the device or coordinate system

$A$:Aerodynamics

$C$:Control

$\text{cmd}$:Commanded

$e$:Difference between the commanded and current values

$F$:In the satellite body coordinate system

$G$:Gravity

$M$:Magnetic field (acting on the bus)

$Q$:Magnetorquers

$R$:Reaction wheels

$\text{SRP}$:Solar radiation pressure

$T$:Propulsion (not used in this study).

+:Matrix pseudoinverse

$\times$:Matrix from vectors of cross-product

$T$:Vector, matrix transpose

dot:Differentiation in time.

References

[1] J. Narkiewicz, S. Topczewski, M. Sochacki, "Nanosatellite formation flying for maritime applications," 7th International Science and Technology Conference on Naval Technologies for Defence and Security NATCON, pp. 345-351, .

[2] D. T. Gerhardt, S. E. Palo, "Passive magnetic attitude control for CubeSat spacecraft," 24th annual AIAA/USU Conference on Small Satellites, .

[3] J. F. Trégouët, D. Arzelier, D. Peaucelle, C. Pittet, L. Zaccarian, "Reaction wheels desaturation using magnetorquers and static input allocation," IEEE Transactions on Control Systems Technology, vol. 23 no. 2, pp. 525-539, DOI: 10.1109/TCST.2014.2326037, 2015.

[4] E. A. Hogan, H. Schaub, "Three-axis attitude control using redundant reaction wheels with continuous momentum dumping," Journal of Guidance, Control, and Dynamics, vol. 38 no. 10, pp. 1865-1871, DOI: 10.2514/1.G000812, 2015.

[5] F. Giulietti, A. A. Quarta, P. Tortora, "Optimal control laws for momentum-wheel desaturation using magnetorquers," Journal of Guidance, Control, and Dynamics, vol. 29 no. 6, pp. 1464-1468, DOI: 10.2514/1.23396, 2006.

[6] X. Chen, W. H. Steyn, S. Hodgart, Y. Hashida, "Optimal combined reaction-wheel momentum management for Earth-pointing satellites," Journal of Guidance, Control, and Dynamics, vol. 22 no. 4, pp. 543-550, DOI: 10.2514/2.4431, 1999.

[7] R. Wisniewski, "Linear time-varying approach to satellite attitude control using only electromagnetic actuation," Journal of Guidance, Control, and Dynamics, vol. 23 no. 4, pp. 640-647, DOI: 10.2514/2.4609, 2000.

[8] M. L. Psiaki, "Magnetic torquer attitude control via asymptotic periodic linear quadratic regulation," Journal of Guidance, Control, and Dynamics, vol. 24 no. 2, pp. 386-394, DOI: 10.2514/2.4723, 2001.

[9] H. C. Polat, J. Virgili-llop, M. Romano, "Survey, statistical analysis and classification of launched CubeSat missions with emphasis on the attitude control method," Journal of Small Satellites, vol. 5, pp. 513-530, 2016.

[10] Y. K. Chang, B. H. Lee, J. W. Choi, M. Y. Yun, S. J. Kang, "Rapid initial detumbling strategy for micro/nanosatellites using a pitch bias momentum system," Transactions of the Japan Society for Aeronautical and Space Sciences, vol. 50 no. 167, pp. 63-69, DOI: 10.2322/tjsass.50.63, 2007.

[11] W. H. Steyn, Y. Hashida, V. Lappas, "An attitude control system and commissioning results of the SNAP-1 nanosatellite," 14th annual AIAA/USU Conference on Small Satellites, .

[12] W. Steyn, Y. Hashida, "An attitude control system for a low-cost earth observation satellite with orbit maintenance capability," 13th annual AIAA/USU Conference on Small Satellites, .

[13] Y. Y. Lin, C. T. Wang, "Detumbling of a rigid spacecraft via torque wheel assisted gyroscopic motion," Acta Astronautica, vol. 93,DOI: 10.1016/j.actaastro.2013.06.021, 2014.

[14] W. Flatley, W. Morgenstern, A. Reth, F. Bauer, "A B-dot acquisition controller for the RADARSAT spacecraft," NASA, pp. 79-90, 1996.

[15] G. Avanzini, F. Giulietti, "Magnetic detumbling of a rigid spacecraft," Journal of Guidance, Control, and Dynamics, vol. 35 no. 4, pp. 1326-1334, DOI: 10.2514/1.53074, 2012.

[16] B. Wie, Space Vehicle Dynamics and Control, Second Edition,DOI: 10.2514/4.860119, 2008.

[17] Z. Ismail, R. Varatharajoo, "A study of reaction wheel configurations for a 3-axis satellite attitude control," Advances in Space Research, vol. 45 no. 6, pp. 750-759, DOI: 10.1016/j.asr.2009.11.004, 2010.

[18] C.-H. Won, "Comparative study of various control methods for attitude control of a LEO satellite," Aerospace Science and Technology, vol. 3 no. 5, pp. 323-333, DOI: 10.1016/S1270-9638(00)86968-0, 1999.

[19] J. Li, M. Post, T. Wright, R. Lee, "Design of attitude control systems for CubeSat-class nanosatellite," Journal of Control Science and Engineering, vol. 2013,DOI: 10.1155/2013/657182, 2013.

[20] M. Fakhari Mehrjardi, H. Sanusi, M. A. M. Ali, M. A. Taher, "PD controller for three-axis satellite attitude control using discrete Kalman filter," 2014 International Conference on Computer, Communications, and Control Technology (I4CT), pp. 83-85, DOI: 10.1109/I4CT.2014.6914151, .

[21] S. R. Omar, "An inverese dynamics satellite attitude determination and control system with autonomous calibration," 29th annual AIAA/USU Conference on Small Satellites, .

[22] M. Moradi, "Satellite neuro-PD three-axis stabilization based on three reaction wheels," Journal of Aerospace Engineering, vol. 27 no. 1, pp. 177-184, DOI: 10.1061/(ASCE)AS.1943-5525.0000231, 2014.

[23] B. Wie, J. Lu, "Feedback control logic for spacecraft eigenaxis rotations under slew rate and control constraints," Journal of Guidance, Control, and Dynamics, vol. 18 no. 6, pp. 1372-1379, DOI: 10.2514/3.21555, 1995.

[24] M. J. Sidi, Spacecraft Dynamics and Control,DOI: 10.1017/CBO9780511815652, 2014.

[25] P. J. Camillo, F. L. Markley, "Orbit-averaged behavior of magnetic control laws for momentum unloading," Journal of Guidance, Control, and Dynamics, vol. 3 no. 6, pp. 563-568, DOI: 10.2514/3.19725, 1980.

[26] "Sinclair Interplanetary, Nanosatellite reaction wheels (RW-0.03)," . January 2018, http://www.sinclairinterplanetary.com/reactionwheels/30mNm-secwheel2016a.pdf?attredirects=0

[27] "Innovative Solutions In Space, Magnetorquer board (iMTQ)," , . January 2018, https://www.cubesatshop.com/wp-content/uploads/2016/06/iMTQ-Brochure-v1.pdf