*3.2.3.3.1 Satellite attitude and angular velocity estimation algorithms*

This group of algorithms was presented above in 3.2.1.2 and can be used here.

For example, let us consider single-axis stabilized satellite that should keep one axis (e.g., *Z*) permanently pointed to the Sun as in **Figure 28**. Only two angles of the satellite deviation from this direction and their angular velocities are required to know to point and keep it in this direction. The satellite has two-axis Sun sensor that

**Figure 28.** *Satellite pointed by the Z-axis to the sun.*

can measure two angles *αx*, *α<sup>y</sup>* of satellite deviation from the sun direction. Its axes coincide with the satellite axes *xyz*. The axis *z* is the sensitivity axis that nominally should point into the Sun's direction (center of brightness), and *xy* is the focal plane. The Sun vector is referenced in the Sun frame as *Sr* <sup>¼</sup> ½ � 00S *<sup>T</sup>* and is measured in the Sun sensor measured frame as *Sm* ¼ *Sx Sy Sz* � �*<sup>T</sup>* .

In **Figure 28**, *xyz* is the satellite body frame, *xsys zs* is the Sun reference frame, SS is the two-axis Sun sensor and *α<sup>x</sup>* and *α<sup>y</sup>* are the turn angles of satellite *x* and *y* axis accordingly.

The following formula represents the mathematical transformation of the Sun vector from the reference into the body frame:

$$\mathbf{S\_m} = \mathbf{C\_{bs}} \mathbf{S\_r} \tag{41}$$

where **Sr** <sup>¼</sup> ½ � 00S *<sup>T</sup>*, **Sm** <sup>¼</sup> *Sx Sy Sz* � �*<sup>T</sup>* , and **Cbs** is the DCM between the reference (Sun) and measured (satellite) frames. Let us consider that the order of rotation from the Sun to the satellite frame is 3-2-1 (*αz*, *αy*, *αx*); the DCM matrix **Cbs** is as follows [9]:

$$\mathbf{C\_{bb}} = \begin{bmatrix} \cos a\_{\mathcal{V}} \cos a\_{\mathcal{z}} & \cos a\_{\mathcal{V}} \sin a\_{\mathcal{z}} & -\sin a\_{\mathcal{V}} \\ -\cos a\_{\mathcal{x}} \sin a\_{\mathcal{z}} + \sin a\_{\mathcal{x}} \sin a\_{\mathcal{y}} \cos a\_{\mathcal{z}} & \cos a\_{\mathcal{x}} \cos a\_{\mathcal{z}} + \sin a\_{\mathcal{x}} \sin a\_{\mathcal{y}} \sin a\_{\mathcal{z}} & \sin a\_{\mathcal{x}} \cos a\_{\mathcal{y}} \\ \sin a\_{\mathcal{x}} \sin a\_{\mathcal{z}} + \cos a\_{\mathcal{x}} \sin a\_{\mathcal{y}} \cos a\_{\mathcal{z}} & -\cos a\_{\mathcal{x}} \cos a\_{\mathcal{z}} + \cos a\_{\mathcal{x}} \sin a\_{\mathcal{y}} \sin a\_{\mathcal{z}} & \cos a\_{\mathcal{x}} \cos a\_{\mathcal{y}} \end{bmatrix} \tag{42}$$

Then from (41), (42) can derive the following formulas:

$$\begin{cases} \mathbb{S}\_{xm} = -\mathbb{S}\sin a\_{\mathcal{Y}} \\ \mathbb{S}\_{ym} = \mathbb{S}\sin a\_{\mathcal{X}}\cos a\_{\mathcal{Y}} \\ \mathbb{S}\_{xm} = \mathbb{S}\cos a\_{\mathcal{X}}\cos a\_{\mathcal{Y}} \end{cases} \tag{43}$$

From (43), desired angles and can be derived that can be used for satellite attitude control.

$$\begin{cases} a\_{\rm x} = \tan^{-1} \frac{\mathcal{S}\_{\rm ym}}{\mathcal{S}\_{\rm xm}} \\\\ a\_{\rm y} = -\sin^{-1} \frac{\mathcal{S}\_{\rm xm}}{\sqrt{\mathcal{S}\_{\rm ym}^2 + \mathcal{S}\_{\rm xm}^2}} \end{cases} \tag{44}$$

#### *3.2.3.3.2 Angular velocities (body rates)*

Let us also assume that the satellite does not have angular velocity sensors RS and its angular velocities should be derived from the measured angles *α<sup>x</sup>* and *αy*. Simple low-frequency first-order differentiating fitters can be applied for this purpose. Laplace operator s-form (transfer functions) of these filters are presented below:'

$$\begin{cases} \hat{a}\_{\mathbf{x}} = \frac{s}{T\_{f\mathbf{x}}s + \mathbf{1}} a\_{\mathbf{x}} \\ \hat{a}\_{\mathbf{y}} = \frac{s}{T\_{f\mathbf{y}}s + \mathbf{1}} a\_{\mathbf{y}} \end{cases} \tag{45}$$

where *Tfx*, *Tfy* are filter time constants, typically, *T <sup>f</sup>* <ð Þ 3 � 10 *s*.

*Satellite Control System: Part I - Architecture and Main Components DOI: http://dx.doi.org/10.5772/intechopen.92575*

#### *3.2.3.3.3 Satellite control algorithms*

If not the optimization criterion to characterize the control quality [13] is required, then conventional negative feedback closed control loop with linear PID (proportional, integral, and damping) control law [9] that provides a good performance for many practical satellite control applications can be used to satisfy the requirements. They are typical for any automatic control system requirements: such as transfer process decay time and overshooting, residual static error caused by the permanent external disturbance, etc. Today, attitude control system performance can be verified mainly on ground with simulation. If we try to evaluate it in flight, then only onboard attitude sensors TLM data can be used for postprocessing, and it should be taken into account that mainly sensors that detect high-frequency noise (perceived errors) will be observable and low-frequency components (sensor biases) are compensated in the closed control attitude stabilization loop. Simple example of single-axis satellite attitude stabilization control loop is presented below. It is a simplified linear model; however, it presents the stabilization principle and essential features. Let us assume that a simple, positional, and damping control law is used to stabilize satellite axis *Z* in Sun direction *S* as in the example C1 for attitude and angular rate determination above. Let us assume that only the one-axis control channel *X* is considered; small angle *α<sup>x</sup>* ð Þ tan *α<sup>x</sup>* ≈*α<sup>x</sup>* is measured with SS (44), and its rate is derived with linear differentiating filter (45). Then requested PD control torque can be presented by the following formula:

$$T\_{\rm cx} = -k\_{\rm px} a\_{\rm xm} - k\_{\rm dx} \frac{s}{T\_{\rm fx} + \mathbf{1}} a\_{\rm xm} \tag{46}$$

where *kpx* is the position control coefficient, *kdx* is the damping control coefficient, *<sup>s</sup>* <sup>¼</sup> *<sup>d</sup> dt* is the Laplace operator, *Tfx* is the differentiating filter time constant, *αxm* ¼ *α<sup>x</sup>* þ Δ*α<sup>x</sup>* is the measured angle *αx*,*α<sup>x</sup>* is the true value and Δ*α<sup>x</sup>* is the measured error. Let us assume that this torque is generated by only one MTR *Y*. Eq. (27) is determined as:

$$T\_{c\mathbf{x}} = B\_{\mathbf{z}} \mathbf{M}\_{\mathbf{y}} = B\_{\mathbf{z}} k\_{\text{TR}\mathbf{y}} \mathbf{u}\_{\mathbf{y}} \tag{47}$$

where *Bz* ¼ *const* is the local *Z* component of Earth magnetic induction vector, *My* ¼ *kTRyuy* is the *Y* MTR magnetic moment, *kTRy* is the *Y* MTR control gate and *iy* is the control voltage applied to *Y* MTR winding. Then as it follows from (46), (47) requested from the AODCS OBC control voltage to the winding of *Y* MTR is:

$$
\mu\_{\rm y} = -K\_{\rm px} a\_{\rm xm} - K\_{\rm dx} \frac{s}{T\_{f\rm x} + \mathbf{1}} a\_{\rm xm} \tag{48}
$$

where *Kpx* <sup>¼</sup> *kpx BzkTRy* and *Kdx* <sup>¼</sup> *kdx BzkTRy* are position and damping magnetic control coefficients.

Let us take a ball-shaped satellite with the inertia matrix as follows:

*J* ¼ *Jx* 0 0 0 *Jy* 0 0 0 *Jz* 2 6 4 3 7 5 where *Jx* ¼ *Jy* ¼ *JzJx* ¼ *Jy* ¼ *Jz*. Then in inertial spatial, its linear

angular dynamical equations for the axis *X* can be approximately written as follows:

$$J\_{\mathbf{x}}\mathbf{s}^{2}a\_{\mathbf{x}} = T\_{c\_{\mathbf{x}}} + T\_{d\_{\mathbf{x}}} \tag{49}$$

where *<sup>s</sup>* <sup>¼</sup> *<sup>d</sup> dt* is the Laplace operator, *Tcx* is the control torque, and *Tdx* is the disturbing external torque (satellite residual and induction magnetism torque, atmosphere drug torque, solar pressure torque). Then substituting in Eq. (49) Eqs. (47) and (48), we can rewrite it as follows:

$$J\_{\mathbf{x}}\mathbf{s}^{2}a\_{\mathbf{x}} = -k\_{\mathbf{p}\mathbf{x}}a\_{\mathbf{x}m} - k\_{d\mathbf{x}}\frac{\mathbf{s}}{T\_{f\mathbf{x}}\mathbf{s}+\mathbf{1}}a\_{\mathbf{x}m} + T\_{d\mathbf{x}}\tag{50}$$

Let us divide all terms in Eq. (50) by the coefficient *kp* and substitute *αxm* value, then it can be represented in the following form:

$$\left(T\_{\rm x}^{2}\mathbf{s}^{2} + 2d\_{\rm x}T\_{\rm x}\frac{\mathbf{s}}{T\_{\rm f\mathbf{x}}\mathbf{s} + \mathbf{1}} + \mathbf{1}\right)a\_{\rm x} = \frac{\mathbf{1}}{k\_{\rm px}}T\_{\rm dx} - \Delta a\_{\rm x} \tag{51}$$

where *Tx* ¼ ffiffiffiffiffi *Jx kpx* <sup>q</sup> is the *<sup>X</sup>* control channel time constant and *dx* <sup>¼</sup> *kdx* <sup>2</sup> ffiffiffiffiffiffiffi *kpxJx* p is the

*X*control channel-specific damping coefficient. Eq. (51) is a third-order linear differential equation and could be analytically analyzed. In particular its stability can be analyzed with algebraic Hurwitz criterion [13]. However, more simple and general results can be obtained with the following approximate consideration. If *Tfx* < < *Tx*, then by the filter time constant, *Tfx* can be neglected, and (51) can approximately be considered as a standard second-order control unit, presented by the second-order linear time invariant (LTI) differential equation and rewritten as follows:

$$(T\_\mathbf{x}^2 \mathbf{s}^2 + 2d\_\mathbf{x} T\_\mathbf{x} \mathbf{s} + \mathbf{1})a\_\mathbf{x} \simeq \frac{\mathbf{1}}{k\_{\text{px}}} T\_{d\mathbf{x}} - \Delta a\_\mathbf{x} \tag{52}$$

As it follows from Eq. (52), steady-state error in attitude stabilization can be calculated with the formula:

$$a\left(T\_{\mathbf{x}}^2\mathbf{s}^2 + 2d\_{\mathbf{x}}T\_{\mathbf{x}}\mathbf{s} + \mathbf{1}\right)a\_{\mathbf{x}} \simeq \frac{1}{k\_{\text{px}}}T\_{d\mathbf{x}} - \Delta a\_{\mathbf{x}} \tag{53}$$

where *Tdx*<sup>0</sup> ¼ *const* and Δ*α<sup>x</sup>*<sup>0</sup> ¼ *const*Δ*α<sup>x</sup>*<sup>0</sup> ¼ *const* (ATT sensor bias). For Eq. (52), *the optimal* damping coefficient is *dx* <sup>¼</sup> ffiffi 2 p <sup>2</sup> ¼ 0*:*707 [13]. **Numerical example**

Let us evaluate satellite time constant *Tx*. Let us assume that for the LEO satellite magnetic field induction vector *B* has the following value of the projection on the Sun direction.

and MTR has the following parameters: maximal magnetic moment *My* max ¼ 35*Am*<sup>2</sup> , maximal control current *Iy* max ¼ 100 mA, winding resistance *R* ¼ 280 Ohm; and maximal control voltage *<sup>u</sup>*ymax <sup>¼</sup> *<sup>R</sup>* � *Iy* max <sup>¼</sup> 28 V. MTR gate is *kTRy* <sup>¼</sup> *My* max *uy* max ¼ 35*Am*<sup>2</sup> <sup>28</sup>*<sup>V</sup>* <sup>¼</sup> <sup>1</sup>*:*25*Am*<sup>2</sup> *<sup>=</sup>V*. Then maximal available magnetic torque is *Tcx* max <sup>¼</sup> <sup>10</sup>�<sup>3</sup> *Nm*. If maximal linear zone for this control channel is *<sup>α</sup><sup>y</sup>* max <sup>¼</sup> *<sup>π</sup>* <sup>2</sup> rad, then the position control coefficient *kpx* can be calculated with the following formula:

$$k\_{p\text{x}} = \frac{T\_{c\text{x}\text{max}}}{a\_{\text{max}}} \tag{54}$$

For the data above, it has the value of *kpx* <sup>¼</sup> <sup>0</sup>*:*<sup>001</sup> <sup>3</sup>*:*14�0*:*<sup>5</sup> <sup>¼</sup> <sup>6</sup>*:*<sup>366</sup> � <sup>10</sup>�<sup>4</sup> Nm*=*rad. Let us consider example of the first Soviet satellite "Sputnik" (SS-1) that had the *Satellite Control System: Part I - Architecture and Main Components DOI: http://dx.doi.org/10.5772/intechopen.92575*

mass *m* ¼ 83*:*6 kg and the radius*R* ¼ 0*:*29 m and take the assumption that its mass was homogeneously distributed within its spherical volume of *<sup>V</sup>* <sup>¼</sup> <sup>4</sup> <sup>3</sup> *<sup>π</sup>R*<sup>3</sup> <sup>¼</sup> <sup>0</sup>*:*1 m3, and then its inertia ( *<sup>J</sup>* <sup>¼</sup> *Jx* <sup>¼</sup> *Jy* <sup>¼</sup> *Jz*<sup>Þ</sup> can be calculated with the following formula:

$$J = \frac{2}{5} mR^2 \tag{55}$$

as for a homogeneous sphere. Substituting into Eq. (55) the data above, we can calculate that for SS-1 *<sup>J</sup>* <sup>¼</sup> <sup>2</sup>*:*82 kgm<sup>2</sup> , then its time constant with the considered MTR control might be:

$$T\_x = \sqrt{\frac{I\_x}{k\_{px}}} = \sqrt{\frac{2.82 \text{kg} m^2}{6.366 \cdot 10^{-4} \text{kg} m^2 s^{-2}}} = 66.56 \text{ s.}$$

Now damping coefficient can be calculated with the following formula:

$$k\_{d\mathbf{x}} = 2d\_{\mathbf{x}}T\_{\mathbf{x}}k\_{p\mathbf{x}} \tag{56}$$

It has the following value: *kdx* <sup>¼</sup> <sup>2</sup> � <sup>0</sup>*:*<sup>707</sup> � <sup>66</sup>*:*<sup>56</sup> � <sup>6</sup>*:*<sup>366</sup> � <sup>10</sup>�<sup>4</sup> <sup>s</sup> � Nm*=*rad <sup>¼</sup> 0*:*06 Nms*=*rd**.**

Finally, *Kpx* and *Kdx* can be calculated. They are as follows: *Kpx* <sup>¼</sup> *kpx BzkTRy* ¼ <sup>6</sup>*:*366�10�4*Nm=rad* <sup>2</sup>*:*86�10�5�1*:*25*T*�*A*�*m*2*=<sup>V</sup>* <sup>¼</sup> <sup>17</sup>*:*<sup>8</sup> *<sup>V</sup>=rad* and *Kdx* <sup>¼</sup> *kdx BzkTRy* <sup>¼</sup> <sup>0</sup>*:*06*Nms=rad* <sup>2</sup>*:*86�10�5�1*:*25*T*�*A*�*m*2*=<sup>V</sup>* <sup>¼</sup> <sup>1</sup>*:*<sup>678</sup> � 10<sup>3</sup> *sV=rad*.

When *Tx* is determined, then the time constant for the differencing filter *Tfx* can be chosen from the condition that *Tfx* < <*Tx*. In our example *Tx* ¼ 66*:*56 s, let us take that *Tfx* ¼ 5 s.

#### **Simulation**

Eqs. (51) and (52) were simulated using Simulink (see **Figure 29**).

**Figure 29.** *Satellite single-axis attitude control Simulink block scheme.*

Blocks in the pink color present the satellite model, the dark green color is for control law blocks, the cyan blocks are registration oscilloscopes, and the display and the orange color are the disturbances. The red manual switch allows to implement the differentiating filter, transforming the scheme from the approximation (52) to the accurate presentation (51). Disturbing external torque *Md* is constant; attitude sensor error is represented by the constant value *ALP*<sup>0</sup> and limited range white noise *V* that has spectral density *SV* <sup>¼</sup> <sup>2</sup>*σ*<sup>2</sup> *VTv* (*σ<sup>v</sup>* is the standard deviation (SD), *Tv* is the correlation time). The results of the simulation (**Figure 29**) with and without the differentiating filter (with the assumption that **α**\_ *<sup>x</sup>α*\_ *<sup>x</sup>* is directly measured without any errors) are presented below in **Figures 30**–**34** (left A, (52), without the filter; right B, (51), with the filter). The numerical data for the simulation are as follows:

$$J = 2.82 \text{ kg} \text{m}^2, T\_x = \sqrt{\frac{f\_x}{k\_{px}}} = 66.56 \text{ s}, d\_x = \frac{k\_{dx}}{2T\_x k\_{px}} = 0.707,$$

$$k\_{px} = 6.366 \cdot 10^{-4} \text{ Nm/rad}, k\_{dx} = 0.06 \text{ Nm} \cdot \text{s/rad}, T\_{fx} = 5 \text{ s},$$

$$M\_d = 10^{-5} \text{Nm}, \Delta a\_{x0} = dltALP\_0 = 0.1^\circ, \Delta a\_{x0} = dltALP\_0 = 0.1^\circ$$

Simulation of ACS (**Figure 29**) is presented in **Figures 30**–**34**. Units: vertical axis (deg), horizontal axis: (s).

**Figure 30.**

*Response to initial deviation angle <sup>α</sup><sup>x</sup>*<sup>0</sup> <sup>¼</sup> 10*. (a) without dif. filter and (b) with dif. filter.*

**Figure 31.** *Response to initial angular velocity α*\_ *<sup>x</sup>*<sup>0</sup> ¼ 0*:*01deg*=s.*

*.*

**Figure 32.** *Response to external disturbance torque Td* <sup>¼</sup> <sup>10</sup>�5Nm*. Static error <sup>α</sup>*<sup>∗</sup> *<sup>x</sup>* <sup>¼</sup> <sup>0</sup>*:*898<sup>∘</sup>

#### *Satellite Control System: Part I - Architecture and Main Components DOI: http://dx.doi.org/10.5772/intechopen.92575*

**Figure 33.**

*Response to attitude sensor bias* <sup>Δ</sup>*<sup>x</sup>*<sup>0</sup> <sup>¼</sup> <sup>1</sup><sup>∘</sup> *plus white noise <sup>σ</sup><sup>V</sup>* <sup>¼</sup> <sup>0</sup>*:*1<sup>∘</sup> , *TV* ¼ 1 s*. Satellite attitude stabilization errors, ALP.*

**Figure 34.** *Satellite attitude measured errors ALPm.*

Decay time: *τ* ¼ 195*s*, *αx*ð Þ¼ 3*τ* 5%*αx*0(deg), (s) *τ* ¼ 275*s***.**

As it can be seen, measured noise is filtered effectively in the control loop, and stabilization error is equal to the sensor bias with opposite sign.

In **Figure** 34, we can see that the measured (perceived) errors that TLM data provide to ground after the decay time do not present sensor bias and present only measured noise. It is because satellite stabilization error with opposite sign compensates the bias. In general, it can also be seen that the simulation of the approximate second-order model (52) is very close to the accurate model (51). Hence, at least for the analytical representation, (52) can be successfully used.

## **4. Conclusion**

Part I of this chapter presents an overview of practical satellite control system, satellite guidance, navigation and control equipment. The work presented here is based on the author's point of view of integration of this GN&C equipment in the integrated AODCS system (satellite GN&C Spacetronics System). Main work principles, architecture, and components of the satellite control system were briefly highlighted.

The chapter can serve to a wide pool of space system specialists as an introduction to satellite control system development.
