Sliding mode corrector for jet UAV control

Abstract

A sliding mode corrector is presented for disturbance rejection in position sensing using relatively accurate velocity measurement. The corrector design is based on a robust second-order sliding mode (2-sliding mode), which makes the fusion of position and velocity on a sliding surface to reject disturbance. Even when the frequency bands of disturbance and actual position signal overlap, or large-magnitude disturbance exists, the corrector can still provide the accurate and smoothed estimate of position. The proposed corrector is applied to a jet UAV navigation and control. In the unmanned aerial vehicle (UAV) system, the disturbances exist in position and attitude measurements, and the uncertainties exist in the system dynamics. For the UAV trajectory tracking control, the system model is constructed in the earth-fixed frame, and the constructed model is fit for observer design to estimate system uncertainties. The control laws are designed according to the correction of position and attitude by the correctors and the estimation of system uncertainties by an existing observer. Finally, the flight experiment demonstrates the effectiveness of the proposed method.

Nomenclature

${p_m}(t)$

position measurement

${p_0}(t)$

actual position

${d_1}(t)$

disturbance in position measurement

${L_1}$

upper-bound of position disturbance

${v_m}(t)$

velocity measurement

${v_0}(t)$

actual velocity

${d_2}(t)$

disturbance in velocity measurement

${L_2}$

upper-bound of velocity disturbance

$\varepsilon $

upper-bound of sensor accuracy ratio

${e_1}$

sliding variable

${e_2}$

sliding variable

${k_1}$

corrector parametre

${k_2}$

corrector parametre

${k_3}$

corrector parametre

${\omega _1}$

position disturbance frequency

$\rho ({\omega _1})$

disturbance rejection ratio

$x$

position in earth-fixed frame x-direction

$y$

position in earth-fixed frame y-direction

$z$

position in earth-fixed frame z-direction

$\phi $

roll angle

$\theta $

pitch angle

$\psi $

yaw angle

$\alpha $

angle-of-attack

$\beta $

sideslip angle

U

linear velocity in body frame axis ${x_b}$

$V$

linear velocity in body frame axis ${y_b}$

$W$

linear velocity in body frame axis ${z_b}$

$\Gamma $

the earth-fixed frame

$\Lambda $

the body frame

$m$

UAV mass

$g$

gravity acceleration

${\Omega _\Lambda }$

angular rate vector

${\Omega _\Gamma }$

Euler angle derivative vector

F

total external force

${F_{jet}}$

thrust by jet engine

${F_w}$

aerodynamic forces on the fixed wing

${F_f}$

aerodynamic forces on the fuselage

${F_r}$

forces created by the rudders

${F_e}$

forces created by the elevators

${F_d}$

uncertainties and external disturbances

$\tau $

total moment

${\tau _w}$

moments created by the fixed wings

${\tau _r}$

moments created by the rudders

${\tau _e}$

moments created by the elevators

${\tau _d}$

moments due to the uncertainties and external disturbances

$\rho $

air density

${S_w}$

area of the half wing

${C_{L0}}$

fixed wing lift coefficient when the angle-of-attack $\alpha $ equals zero

${C_{L\alpha }}$

fixed wing lift coefficient due to the angle-of-attack $\alpha $

${\delta _i}$

fixed wing aileron deflection

${C_{L{\delta _i}}}$

lift coefficient due to the aileron deflection ${\delta _i}$

${C_{D0}}$

fixed wing drag coefficient when $\alpha = {\delta _i} = 0$

${A_w}$

aspect ratio of the fixed wing

${e_w}$

value of the Oswald’s efficiency factor

${\tau _{wa}}$

fixed wing aerodynamic moment

${\tau _{wc}}$

fixed wing control torque

${S_f}$

fuselage equivalent cross-sectional area

${L_f}$

lift force generated by the fuselage

${D_f}$

drag force generated by the fuselage

${C_{lf}}$

fuselage lift coefficient

${C_{df}}$

fuselage drag coefficient

${C_{df0}}$

fuselage constant in the coefficient of drag force

${S_e}$

area of the elevator

${\delta _e}$

elevator deflection

${C_{le\alpha }}$

eleviator lift coefficient due to the angle-of-attack $\alpha $ and the deflection ${\delta _e}$

${C_{de0}}$

drag coefficient when $\alpha + {\delta _e} = 0$

${A_e}$

aspect ratio of the elevator

${e_e}$

the Oswald’s efficiency factor

${\tau _{ea}}$

elevator aerodynamic moment

${\tau _{ec}}$

elevator control torque

${S_r}$

area of the rudders

${C_{lr\beta }}$

rudder lift coefficient due to the sideslip angle $\beta $

${\delta _r}$

rudder deflection

${C_{lr{\delta _r}}}$

rudder lift coefficient due to the deflection ${\delta _r}$

${C_{dr0}}$

rudder drag coefficient when $\beta = {\delta _r} = 0$

${A_r}$

aspect ratio of the rudder

${e_r}$

rudder Oswald’s efficiency factor

${\tau _{ra}}$

rudder aerodynamic moment

${\tau _{rc}}$

rudder control torque

${\Delta _x}$

uncertainty in x-direction

${\Delta _y}$

uncertainty in y-direction

${\Delta _z}$

uncertainty in z-direction

${\Delta _\phi }$

uncertainty in roll

${\Delta _\theta }$

uncertainty in pitch

${\Delta _\psi }$

uncertainty in yaw

${y_{ * 1}}$

measurement for position/angle

${y_{ * 2}}$

measurement for velocity/angular velocity

${x_d}$

reference position in x-direction

${y_d}$

reference position in y-direction

${z_d}$

reference position in z-direction

${\phi _d}$

desired roll angle

${\theta _d}$

desired pitch angle

${\psi _d}$

desired yaw angle

${\lambda _{ * 1}}$

parametre of extended observer

${\lambda _{ * 2}}$

parametre of extended observer

${\alpha _ * }$

parametre of extended observer

${k_{p1}}$

controller parametre in position dynamics

${k_{p2}}$

controller parametre in position dynamics

${k_{a1}}$

controller parametre in attitude dynamics

${k_{a2}}$

controller parametre in attitude dynamics

1.0 Introduction

This paper considers correction of stochastic disturbance in position sensing and application to jet UAV navigation and control. This interest was motivated by the enormous civil and military applications of such fixed wing UAVs. It is one of the most attractive research focuses because the dynamical system of a jet UAV has many prominent features including powerful thrust provision, payload augmentation, high-speed flight and a high manoeuverability [1–3]. UAV large-range flight needs information of global position, attitude and dynamic model, also flying velocity and angular velocity are necessary. However, in many cases, disturbances exist in position and attitude sensing, and uncertainties are inevitable in system modelling. These bring challenge for control.

In flight control systems, rational desired attitude is important for safe flight, and the determination of desired attitude needs the information of actual position and attitude [4]. However, disturbances in position and attitude sensing render the incorrect desired attitude, and the unwanted control command is generated. Constant sensing disturbance can be overcome through initial calibration. Comparing to constant sensing disturbance, time-varying position disturbance is more likely to rend a serious mismatch between desired attitude and actual position, and it causes dangerous flight. Furthermore, the frequency bands of disturbance and actual position signal may overlap, and disturbance cannot be separated from actual position signal using the usual low-pass filters.

GPS (global positioning system) can provide global position information with accuracy of several metres or even tens of metres [5, 6]. Adverse environmental influences may contaminate GPS signals [6], and the position accuracy may become worse. Velocity is also important for UAV navigation and control. GPS can measure device velocity with two different accuracies: (1) large-error velocity by the difference method with accuracy of a metre per second due to position accuracy and noise effect; (2) accurate velocity by Doppler shift measurement with accuracy of a few centimetres, or even the accuracy approaching 5mm/s is possible [7, 8]. Alternatively, accurate velocity of device can be measured by a Doppler radar sensor with accuracy of a few centimetres [9]. Hence, measuring Doppler shift is a preferred way to get velocity. Except for sensing, velocity can be estimated from position using the observers or differentiators [10, 11]. However, relatively accurate measurement of position is required.

INS (inertial navigation system) can estimate position and velocity through integrations from acceleration measurement. However, measurement error or non-zero mean noise in acceleration through integrations cause velocity and position to drift over time. The observer-based INS methods were used to estimate unknown variables in navigation [12, 13]. However, position signals are limited to be local, but not global. For attitude information, an IMU (inertial measurement unit) can determine attitude angle from accurate angular velocity through integration, but angle drift happens. Meanwhile, the outputs of the accelerometres and the magnetometre in IMU can determine the large-error pitch, roll and yaw angles [14].

Uncertainties in UAV flight dynamics include: aerodynamic disturbance, unmodelled dynamics and parametric uncertainties. These uncertainties bring challenges for control system design. The uncertainty in a system can be estimated by an extended state observer [15, 16]. However, accurate position measurement is required as the input of observer. Even velocity can be use for estimation, disturbance in position cannot still be corrected. For a jet UAV flight control system, an integral-uncertainty observer was designed to estimate the attitude angles and attitude dynamic uncertainty, and an augmented observer was used to estimate the flying velocity and position dynamic uncertainty [17]. However, drift may happen for long-time flight due to effect of disturbance and actuator vibrations on the IMU. Meanwhile, the augmented observer can only reduce high-frequency noise, low- or mid-frequency disturbances still exist.

In order to reduce disturbances in position and attitude, the popular methods of GPS/INS based on KF (Kalman filter) or EKF (extended Kalman filter) are used for signal fusion to overcome the limits of individual measurements based on optimisation of a recursive least mean square error [18–21]. Thus, measurement accuracy is improved. For KF or EKF, the relatively accurate system models are needed. Furthermore, the uncertainty in noise statistics limits the performance. In addition, for EKF, system model linearisation may cause filtering divergence, and the derivation of the Jacobian matrices are nontrivial. A finite-time-convergent signal corrector was designed for position correction in a quadrotor UAV control system [22]. The signal corrector is complex, and the finite-time convergence cannot be implemented in engineering practice. Furthermore, the parametres’ selection is sensitive to the estimate performance.

In this paper, a corrector based on robust 2-sliding mode is presented to correct position disturbance using relatively accurate velocity measurement. The 2-sliding mode can reduce the estimate errors of corrector and make the fusion of position and velocity on a linear sliding surface. Position disturbance is reduced further on the sliding surface. Not only the corrector can reject high-frequency noise, but also the low- and mid-frequency disturbances are reduced largely. Therefore, the corrector can reject low/mid/high frequency disturbances, and it is unrelated to the types of actual position signals. Due to the existence of linear sliding surface in the 2-sliding mode, the estimate outputs from the corrector are accurate and smoothed.

The contributions of the proposed corrector include: (1) the corrector can reject position disturbance in low/mid/high frequency bands; (2) due to the continuity of 2-sliding mode, the estimate outputs from the corrector are smoothed and accurate, and they can be used directly for control without any additional filters; (3) due to the robust sliding mode, the corrector parametres are highly inclusive to change of disturbance and signal; (4) because only switch logic and linear functions are used in the corrector, the corrector can be implemented easily in the current hardware of computational environments.

The proposed corrector is applied to navigation and control of a jet UAV. In the UAV flight test, the following adverse conditions are considered: disturbances in the measurements of GPS position and IMU attitude angles, and uncertainties in the UAV flight dynamics. For the UAV trajectory tracking control, the UAV system model is constructed in the earth-fixed frame [23]. Furthermore, the model is fit for observer design to estimate the system uncertainties. The correctors are adopted to correct the disturbances in GPS position and IMU attitude angles. In addition, an existing extended state observer [24] is used to estimate the uncertainties in the UAV flight dynamics. The performance of corrector is compared to the KF-based signal fusion methods [21, 25]. Moreover, based on the correction and estimation, the desired attitude is determined, and the control laws are designed to drive the UAV to achieve the flight mission.

2.0 Problem description

The problem considered in this paper is to reject disturbance in position and attitude sensing.

2.1 Position and velocity sensing

GPS provides position of a device, and accurate velocity can be determined by GPS with Doppler shift measurement or by a Doppler radar sensor.

Define the position measurement: ${p_m}(t) = {p_0}(t) + {d_1}(t)$ , where, ${p_0}(t)$ is the actual position; ${d_1}(t)$ is the disturbance in position measurement, and $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_1}(t)} \right| \le {L_1} \lt \infty $ .

Define the velocity measurement: ${v_m}(t) = {v_0}(t) + {d_2}(t)$ , where, ${v_0}(t)$ is the actual velocity; ${d_2}(t)$ is the disturbance in velocity measurement, and $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_2}(t)} \right| \le {L_2} \lt \infty $ .

Remark 2.1: The accuracy ${L_1}$ of GPS position sensing is usually a metre, a few metres or even tens of metres. Doppler shift measurement enables velocity accuracy ${L_2}$ of a few centimetres per second, even the accuracy approaching 5mm/s (i.e. 0.005m/s) is possible [8]. Therefore, ${L_2} \ll 1$ holds. When using the consistent unit standard, we can get ${L_2} \ll {L_1}$ , i.e. the sensor accuracy ratio $\frac{{{L_2}}}{{{L_1}}} \ll 1$ .

Furthermore, due to the reliability of Doppler measurement, velocity accuracy usually remains unchanged. However, the position accuracy may become worse because of different environmental influences, i.e. ${L_1}$ may increase. Therefore, there exists a small constant $\varepsilon \gt 0$ , i.e. the upper-bound of sensor accuracy ratio, such that the sensor accuracy inequality $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ holds.

2.2 Attitude sensing

The gyroscopes in IMU provide the relatively accurate angular velocities, e.g. their accuracy is about ${L_2}{ = 10^ \circ }$ /hr = 0.003º/s. The accuracy of attitude angles from IMU is relatively large, e.g. about ${L_1}{ = 1.0^ \circ }$ . Therefore, $\frac{{{L_2}}}{{{L_1}}} = \frac{{0.003}}{1} = 0.003 \ll 1$ , and we can select an upper-bound of sensor accuracy ratio $0.003 \le \varepsilon = 0.003 \ll 1$ to satisfy $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ .

2.3 Effect of sensing disturbance on safe flight

Sensing disturbance has a serious impact on safe flight. Disturbance in position sensing may generate incorrect desired attitude angles, and they are mismatched to actual position trajectory. Therefore, the determined attitude is unwanted, and it may be dangerous.

In addition, in position measurement, the frequency bands of actual position ${p_0}(t)$ and sensing disturbance ${d_1}(t)$ may have intersections: the disturbance ${d_1}(t)$ may be in low-, mid- or high-frequency bands. It is impossible for the usual filters to separate ${d_1}(t)$ from the actual position signal ${p_0}(t)$ .

Questions: How to reject disturbance ${d_1}(t)$ in position measurement ${p_m}(t)$ using the relatively accurate velocity ${v_m}(t)$ ? Also, how to reject the disturbance in angle measurement using the relatively accurate angular velocity?

3.0 Preliminary, corollary and notation

3.1 Preliminary

The related background is presented here.

Lemma 3.1 (sliding mode with prescribed convergence law) [26, 27]: The following system is considered:

(1) \begin{align}\dot{e}_1 & = {e_2}\nonumber\\ \dot{e}_2 & = \varphi ( t ) + \gamma ( t )\dot{u} \end{align}

where, ${e_1}$ and ${e_2}$ are the states; u is the controller, and $\dot{u}$ in the system; $\left| {\varphi \left( t \right)} \right| \le \Phi $ , $0 \lt {\Gamma _m} \le \gamma \left( t \right) \le {\Gamma _M}$ , $\Phi \gt 0$ . A 2-sliding control algorithm is as follows:

(2) \begin{equation}\dot{u} = \begin{cases} - u,\ {\rm{if}}\left| u \right| \gt 1;\\ - {V_M}{\rm{sign}}\left[ {{e_2} + g({e_1})} \right]\!,\ {\rm{if}}| u | \le 1\end{cases}\end{equation}

where, $g({e_1})$ is smooth everywhere except on ${e_1} = 0$ , for example, $g({e_1}) = {k_1}{\left| {{e_1}} \right|^\alpha }{\rm{sign}}({e_1})$ , $\alpha \in \left[ {0.5,1} \right)$ ; ${V_M} \gt \frac{\Phi + \sup \left[ {{g^\prime}\left( {e_1} \right)g\left( {e_1} \right)} \right]}{\Gamma_m}$ . Then, we get the finite-time convergence law ${\dot{e}_1} = - g({e_1})$ (i.e. sliding surface ${e_2} + g({e_1}) = 0$ ), and there exists a finite time ${t_s} \gt 0$ , for $t \geqslant {t_s}$ , such that

(3) \begin{equation}{e_1} = 0\ {\rm{and}}\ {e_2} = 0\end{equation}

Remark 3.1: For system (1), when we select $\varphi \left( t \right) = 0$ and $\gamma \left( t \right) = 1$ , it becomes

(4) \begin{align}{\dot{e}_1} & = {e_2}\nonumber\\{\dot{e}_2} & = \dot{u}\end{align}

Then, the 2-sliding control algorithm (2) is expressed by

(5) \begin{equation}\dot{u} = {\dot{e}_2} = \begin{cases} - {e_2},\ {\rm{if}}\left| {{e_2}} \right| \gt 1;\\ - {V_M}{\rm{sign}}\left[ {{e_2} + g({e_1})} \right]\!,\ {\rm{if}}\left| {{e_2}} \right| \le 1\end{cases}\end{equation}

3.2 Corollary on 2-sliding mode system

Combining the system (4) and the 2-sliding control algorithm (5), we get the following corollary on a 2-sliding mode with prescribed convergence law.

Corollary 3.1 (sliding mode with prescribed finite-time convergence law): A 2-sliding mode system is as follows:

(6) \begin{align}{\dot e_1} & = {e_2}\nonumber\\{\dot e_2} & = \left\{ {\begin{array}{*{20}{l}}{ - {e_2},\ {\rm{if}}\left| {{e_2}} \right| \gt 1;}\\{ - {V_M} \cdot {\rm{sign}}\left[ {{e_2} + g({e_1})} \right]\!,\ {\rm{if}}\left| {{e_2}} \right| \le 1}\end{array}} \right.\end{align}

where, ${e_1}$ and ${e_2}$ are the sliding variables; $g({e_1})$ is smooth everywhere except on ${e_1} = 0$ , for example, $g({e_1}) = {k_1}{\left| {{e_1}} \right|^\alpha }{\rm{sign}}({e_1})$ , $\alpha \in \left[ {0.5,1} \right)$ ; ${V_M} \gt \sup \left[ {{g^\prime}\left( {{e_1}} \right)g\left( {{e_1}} \right)} \right]$ . Then, we get the finite-time convergence law ${\dot e_1} = - g({e_1})$ (i.e. sliding surface ${e_2} + g({e_1}) = 0$ ), and there exists a finite time ${t_s} \gt 0$ , for $t \geqslant {t_s}$ , such that

(7) \begin{equation}{e_1} = 0\ {\rm{and}}\ {e_2} = 0\end{equation}

Remark 3.2: For system (6), the parametre selection condition ${V_M} \gt \sup \left[ {{g^\prime}\left( {{e_1}} \right)g\left( {{e_1}} \right)} \right]$ is too strict. In order to relax the parametre selection conditions, we can use a linear convergence law ${\dot e_1} = - {k_1}{e_1}$ for the nonlinear ${\dot e_1} = - g({e_1})$ , and only ${V_M} \gt {k_1} \gt 0$ will be required. In the following section, we will give a theorem on 2-sliding mode system with linear convergence law to be exponentially stable.

3.3 Notation

“ $a( \omega )\;:\;{b_1} \to {b_2}$ as $\omega \;:\; {c_1} \to {c_2}$ ” means that function $a\left( \omega \right)$ varies monotonically increasing or decreasing from ${b_1}$ to ${b_2}$ as $\omega $ increases from ${c_1}$ to ${c_2}$ .

4.0 Robust 2-sliding mode system

Before we present the design of sliding mode corrector, we give a 2-sliding mode system, and a Theorem is presented as follows.

Theorem 4.1 (sliding mode with prescribed linear convergence law): The 2-sliding mode system is as follows:

(8a) \begin{align}{\dot e_1} & = {e_2}\nonumber\\[5pt]{\dot e_2}& = \left\{ {\begin{array}{*{20}{l}}{ - {e_2},\ {\rm{if}}\left| {{e_2}} \right| \gt 1;}\\{ - {k_2}{\rm{sign}}\left( {{e_2} + {k_1}{e_1}} \right),\ {\rm{if}}\left| {{e_2}} \right| \le 1}\end{array}} \right.\end{align}

or

(8b) \begin{align}{\dot e_1} & = {e_2}\nonumber\\[5pt]{\dot e_2} & = \left\{ {\begin{array}{*{20}{l}}{ - {k_3}{\rm{sign}}\left( {{e_2}} \right),\ {\rm{if}}\left| {{e_2}} \right| \gt 1;}\\{ - {k_2}{\rm{sign}}\left( {{e_2} + {k_1}{e_1}} \right),\ {\rm{if}}\left| {{e_2}} \right| \le 1}\end{array}} \right.\end{align}

where, ${e_1}$ and ${e_2}$ are the sliding variables; ${k_2} \gt {k_1} \gt 0$ , and ${k_3} \gt 0$ . Then, we get the linear convergence law ${\dot e_1} = - {k_1}{e_1}$ (i.e. sliding surface ${e_2} + {k_1}{e_1} = 0$ ), and the system (8a) or (8b) is exponentially stable, i.e.

(9) \begin{equation}\mathop {\lim }\limits_{t \to \infty } {e_1} = 0\ {\rm{and}}\ \mathop {\lim }\limits_{t \to \infty } {e_2} = 0\end{equation}

The proof of Theorem 4.1 is presented in Appendix. $\blacksquare $

In fact, for (8b), when $\left| {{e_2}} \right| \gt 1$ , we use ${\dot e_2} = - {k_3}{\rm{sign}}\left( {{e_2}} \right)$ for ${\dot e_2} = - {e_2}$ to speed up ${e_2}$ convergence and overcome disturbance, where, ${k_3} \gt 0$ .

Simulation example (sliding mode with prescribed linear convergence law): For the sliding mode system (8b), we select ${k_1} = 1$ , ${k_2} = 10$ , and ${k_3} = 5$ . Then, we get ${e_1}$ and ${e_2}$ in Fig. 1. Figure 1 illustrates the fast convergence of ${e_1}$ and ${e_2}$ .

Figure 1. Sliding variables ${e_1}$ and ${e_2}$ .

In the following, we give a robust 2-sliding mode system considering existence of multiple disturbances, and a Theorem is presented as follows.

Theorem 4.2 (robust 2-sliding mode): A 2-sliding mode system considering the unknown bounded disturbances is as follows:

(10) \begin{align}{\dot e_1}& = {e_2}\nonumber\\{\dot e_2} & = \left\{ {\begin{array}{*{20}{l}}{ - {k_3}{\rm{sign}}\left[ {{e_2} - {d_2}(t)} \right] - {d_3}(t),\ {\rm{if}}\left| {{e_2} - {d_2}(t)} \right| \gt 1;}\\{ - {k_2}{\rm{sign}}\left[ {{e_2} - {d_2}(t) + {k_1}({e_1} - {d_1}(t))} \right] - {d_4}(t),\ {\rm{if}}\left| {{e_2} - {d_2}(t)} \right| \le 1}\end{array}} \right.\end{align}

where, ${e_1}$ and ${e_2}$ are the sliding variables; the unknown bounded disturbances ${d_i}(t)$ satisfy $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_i}(t)} \right| \le {L_i} \lt \infty $ ( $i = 1,2,3,4$ ); ${L_2} \ll 1$ , and there exists a constant $\varepsilon \gt 0$ such that $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ ; ${k_1} = 2{\varepsilon ^{1 - r}}$ , where $r \in \left( {0,\frac{1}{2}} \right]$ ; ${k_2} \gt {k_1} + {L_4}$ , and ${k_3} \gt {L_3}$ . Then, the effect of disturbances is rejected, and the variables ${e_1}$ and ${e_2}$ of system (10) are in the bounds as follows:

(11) \begin{equation}\mathop {\lim }\limits_{t \to \infty } \left| {{e_1}} \right| \le \rho ({\omega _1}){L_1};\ \mathop {\lim }\limits_{t \to \infty } \left| {{e_2}} \right| \le {L_2}\end{equation}

where, ${\omega _1}$ is the angular frequency of disturbance ${d_1}(t)$ ; and the rejection ratio is expressed by

(12) \begin{equation}\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + \frac{1}{4}{\varepsilon ^{2r - 2}}\omega _1^2} }} + \frac{1}{2}{\varepsilon ^r}\end{equation}

The rejection ratio $\rho ({\omega _1})$ is a monotonically decreasing function of ${\omega _1} \in \left[ {0,\infty } \right)$ , and it satisfies:

1. (i) In $\left[ {{\omega _0},\infty } \right)$ , $\rho ({\omega _1})\;:\;{\varepsilon ^r} \to \frac{1}{2}{\varepsilon ^r}$ as ${\omega _1}\;:\;{\omega _0} \to \infty $ , where, ${\omega _0} = 4{\varepsilon ^{1 - 2r}}\sqrt {1 - \frac{1}{4}{\varepsilon ^{2r}}} $ ;

2. (ii) In $\left( {{\omega _c},{\omega _0}} \right)$ , $\rho ({\omega _1})\;:\;1 \to {\varepsilon ^r}$ as ${\omega _1}\;:\;{\omega _c} \to {\omega _0}$ , where, ${\omega _c} = \frac{{{\varepsilon ^{1 - \frac{1}{2}r}}\sqrt {4 - {\varepsilon ^r}} }}{{1 - \frac{1}{2}{\varepsilon ^r}}} \lt 4{\varepsilon ^{1 - \frac{1}{2}r}} \ll 1$ ;

3. (iii) In $\left[ {0,{\omega _c}} \right]$ , $\rho ({\omega _1})\;:\;1 + \frac{1}{2}{\varepsilon ^r} \to 1$ (i.e. $\rho ({\omega _1}) \approx 1$ due to $0 \lt \varepsilon \ll 1$ ) as ${\omega _1}\;:\;0 \to {\omega _c}$ , and this frequency band is sufficiently small due to ${\omega _c} \ll 1$ .

The proof of Theorem 4.2 is presented in Appendix. $\blacksquare $

Remark 4.1: From $0 \lt \varepsilon \ll 1$ and $r \in \left( {0,\frac{1}{2}} \right]$ , we can get that both the rejection ratio $\rho ({\omega _1})$ in frequency band $\left[ {{\omega _0},\infty } \right)$ and the frequency ${\omega _0}$ are small enough. Therefore, the disturbance ${d_1}\left( t \right)$ is reduced at very small rejection ratio in the large frequency band $\left[ {{\omega _0},\infty } \right)$ . In fact, the disturbance bound ${L_2}$ in velocity sensing through Doppler effect is a few centimetres or a few millimetres, and the disturbance bound ${L_1}$ in position sensing is usually about a few metres or even tens of metres. Thus, the sensing accuracy inequality $\frac{{{L_2}}}{{{L_1}}} \ll 1$ holds, and there exists $\varepsilon \gt 0$ such that $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ . When disturbance ${d_1}\left( t \right)$ becomes worse, i.e. ${L_1}$ increases, the inequality $\frac{{{L_2}}}{{{L_1}}} \le \varepsilon \ll 1$ still holds. Therefore, the frequency band $\left[ {{\omega _0},\infty } \right)$ covers the low/mid/high frequency bands, and the disturbance ${d_1}(t)$ in position sensing in $\left[ {{\omega _0},\infty } \right)$ is rejected sufficiently by the corrector. Furthermore, the disturbance ${d_1}( t )$ can still be rejected largely in the other frequency bands.

In the following, we consider the position disturbance ${d_1}(t)$ is rejected to the maximum extent in a given frequency band, and the disturbance in the other bands can still be rejected largely. We will determine the corrector parametre ${k_1}$ to get the minimum value of the rejection ratio in the given frequency band $\left[ {{\omega _{req}},\infty } \right)$ , and a Theorem is presented as follows.

Theorem 4.3 (sufficient disturbance rejection in given frequency band): The sliding mode system (10) is considered, where, the unknown bounded disturbances ${d_i}(t)$ satisfy $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_i}(t)} \right| \le {L_i} \lt \infty $ ( $i = 1,2,3,4$ ); ${L_2} \ll 1$ , and there exists a constant $\varepsilon \gt 0$ such that $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ . For the a given ${\omega _{req}} \gt 0$ , if we select ${k_3} \gt {L_3}$ , ${k_2} \gt {k_1} + {L_4}$ and ${k_1} = 1/{x_{\min }}$ , where, ${x_{\min }}$ is the unique solution to

(13) \begin{equation}\omega _{req}^2{x_{\min }}{(1 + \omega _{req}^2x_{\min }^2)^{ - \frac{3}{2}}} - \varepsilon = 0\end{equation}

in the range $(\frac{1}{{\sqrt 2 {\omega _{req}}}},\infty )$ , then, the effect of disturbances is rejected, especially it is rejected sufficiently in the frequency band $\left[ {{\omega _{req}},\infty } \right)$ , and the variables ${e_1}$ and ${e_2}$ of system (10) are in the bounds as follows:

(14) \begin{equation}\mathop {\lim }\limits_{t \to \infty } \left| {{e_1}} \right| \le \rho ({\omega _1}){L_1};\ \mathop {\lim }\limits_{t \to \infty } \left| {{e_2}} \right| \le {L_2}\end{equation}

where, the rejection ratio is expressed by

(15) \begin{equation}\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + x_{\min }^2\omega _1^2} }} + \varepsilon \cdot {x_{\min }}\end{equation}

The rejection ratio $\rho ({\omega _1})$ is a monotonically decreasing function of ${\omega _1} \in \left[ {0,\infty } \right)$ , and it satisfies:

1. (i) In $\left[ {{\omega _0},\infty } \right)$ , $\rho ({\omega _1})\;:\;{\rho _{\min }} \to \varepsilon \cdot {x_{\min }}$ as ${\omega _1}\;:\;{\omega _{req}} \to \infty $ .

2. (ii) In $\left( {{\omega _c},{\omega _0}} \right)$ , $\rho ({\omega _1})\;:\;1 \to {\rho _{\min }}$ as ${\omega _1}\;:\;{\omega _c} \to {\omega _{req}}$ , where, ${\omega _c} = \frac{{\sqrt {\frac{{2\varepsilon }}{{{x_{\min }}}}} \sqrt {1 - \frac{1}{2}\varepsilon {x_{\min }}} }}{{1 - \varepsilon \cdot {x_{\min }}}}$ .

3. (iii) In $\left[ {0,{\omega _c}} \right]$ , $\rho ({\omega _1})\;:\;1 + \varepsilon \cdot {x_{\min }} \to 1$ as ${\omega _1}\;:\;0 \to {\omega _c}$ .

where,

(16) \begin{equation}{\rho _{\min }} = \frac{1}{{\sqrt {1 + \omega _{req}^2x_{\min }^2} }} + \varepsilon \cdot {x_{\min }}\end{equation}

The proof of Theorem 4.3 is presented in Appendix. $\blacksquare $

5.0 Design of sliding mode corrector

According to Theorem 4.2, and considering disturbance in position sensing, a sliding mode corrector is designed to reject the disturbance, and a theorem is presented as follows.

Theorem 5.1 (sliding mode corrector): Suppose position measurement is ${p_m}(t) = {p_0}(t) + {d_1}(t)$ , and velocity measurement is ${v_m}(t) = {v_0}(t) + {d_2}(t)$ ; ${p_0}(t)$ is the actual position, and ${v_0}(t)$ is the actual velocity; the unknown bounded disturbances ${d_1}(t)$ and ${d_2}(t)$ satisfy $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_i}(t)} \right| \le {L_i} \lt \infty $ , where, $i = 1,2$ ; ${L_2} \ll 1$ , and there exists a constant $\varepsilon \gt 0$ such that $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ ; $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{{\dot v}_0}(t)} \right| \le {L_3} \lt \infty $ . The sliding mode corrector is designed as follows:

(17) \begin{align}{\dot x_1}& = {x_2}\nonumber\\ {\dot x_2} & = \left\{ {\begin{array}{*{20}{l}}{ - {k_3}{\rm{sign}}\left[ {{x_2} - {v_m}(t)} \right]\!,\ {\rm{if}}\left| {{x_2} - {v_m}(t)} \right| \gt 1;}\\{ - {k_2}{\rm{sign}}\left[ {{x_2} - {v_m}(t) + {k_1}({x_1} - {p_m}(t))} \right]\!,\ {\rm{if}}\left| {{x_2} - {v_m}(t)} \right| \le 1}\end{array}} \right.\end{align}

where, ${x_1}$ and ${x_2}$ are the corrector variables; the corrector parametres satisfy ${k_1} = 2{\varepsilon ^{1 - r}}$ (where $r \in \left( {0,\frac{1}{2}} \right]$ ), ${k_2} \gt {k_1} + {L_3}$ and ${k_3} \gt {L_3}$ . Then, the disturbance ${d_1}(t)$ is rejected, and the corrector estimate outputs satisfy:

(18) \begin{equation}\mathop {\lim }\limits_{t \to \infty } \left| {{x_1} - {p_0}(t)} \right| \le \rho ({\omega _1}){L_1};\ \mathop {\lim }\limits_{t \to \infty } \left| {{x_2} - {v_0}(t)} \right| \le {L_2}\end{equation}

where, ${\omega _1}$ is the angular frequency of disturbance ${d_1}(t)$ , and the rejection ratio is expressed by

(19) \begin{equation}\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + \frac{1}{4}{\varepsilon ^{2r - 2}}\omega _1^2} }} + \frac{1}{2}{\varepsilon ^r}\end{equation}

The rejection ratio $\rho ({\omega _1})$ is a monotonically decreasing function of ${\omega _1} \in \left[ {0,\infty } \right)$ , and it satisfies:

1. (i) In $\left[ {{\omega _0},\infty } \right)$ , $\rho ({\omega _1})\;:\;{\varepsilon ^r} \to \frac{1}{2}{\varepsilon ^r}$ as ${\omega _1}\;:\;{\omega _0} \to \infty $ , where, ${\omega _0} = 4{\varepsilon ^{1 - 2r}}\sqrt {1 - \frac{1}{4}{\varepsilon ^{2r}}} $ ;

2. (ii) In $\left( {{\omega _c},{\omega _0}} \right)$ , $\rho ({\omega _1})\;:\;1 \to {\varepsilon ^r}$ as ${\omega _1}\;:\;{\omega _c} \to {\omega _0}$ , where, ${\omega _c} = \frac{{{\varepsilon ^{1 - \frac{1}{2}r}}\sqrt {4 - {\varepsilon ^r}} }}{{1 - \frac{1}{2}{\varepsilon ^r}}} \lt 4{\varepsilon ^{1 - \frac{1}{2}r}} \ll 1$ ;

3. (iii) In $\left[ {0,{\omega _c}} \right]$ , $\rho ({\omega _1})\;:\;1 + \frac{1}{2}{\varepsilon ^r} \to 1$ (i.e. $\rho ({\omega _1}) \approx 1$ ) as ${\omega _1}\;:\;0 \to {\omega _c}$ . This frequency band is sufficiently small due to ${\omega _c} \ll 1$ .

The proof of Theorem 5.1 is presented in Appendix. $\blacksquare $

Remark 5.1 (analysis of corrector (17)): The sensing accuracy inequality $0 \lt \varepsilon \ll 1$ and $r \in \left( {0,\frac{1}{2}} \right]$ guarantee that the rejection ratio $\rho ({\omega _1})$ in frequency band $\left[ {{\omega _0},\infty } \right)$ is small enough, and $\left[ {{\omega _0},\infty } \right)$ covers the low/mid/high frequency bands. Therefore, the disturbance ${d_1}(t)$ in $\left[ {{\omega _0},\infty } \right)$ can be rejected sufficiently by the corrector. Furthermore, the disturbance ${d_1}(t)$ in the low frequency band $\left( {{\omega _c},{\omega _0}} \right)$ is still reduced largely. Only the approximate constant disturbances in the extreme low-frequency band $\left[ {0,{\omega _c}} \right]$ are not rejected. In fact, signals in the extreme low-frequency band $\left[ {0,{\omega _c}} \right]$ is approximate constant due to sufficiently small ${\omega _c} \ll 1$ . Therefore, the position disturbance ${d_1}(t)$ can be rejected sufficiently by the corrector even when the disturbance frequency covers low/mid/high frequency bands.

From Theorem 4.3, considering sufficient disturbance rejection in a given frequency band, we get the optimal sliding mode corrector, and a theorem is presented as follows.

Theorem 5.2 (optimal sliding mode corrector): The corrector (17) and the measurement signals in Theorem 5.1 are considered, where, the unknown bounded disturbances ${d_i}(t)$ satisfy $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_i}(t)} \right| \le {L_i} \lt \infty $ ( $i = 1,2$ ); ${L_2} \ll 1$ , and there exists a constant $\varepsilon \gt 0$ such that $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ ; $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{{\dot v}_0}(t)} \right| \le {L_3} \lt \infty $ . For the a given ${\omega _{req}} \gt 0$ , if we select ${k_3} \gt {L_3}$ , ${k_2} \gt {k_1} + {L_3}$ and ${k_1} = 1/{x_{\min }}$ , where, ${x_{\min }}$ is the unique solution to

(20) \begin{equation}\omega _{req}^2{x_{\min }}{(1 + \omega _{req}^2x_{\min }^2)^{ - \frac{3}{2}}} - \varepsilon = 0\end{equation}

in the range $(\frac{1}{{\sqrt 2 {\omega _{req}}}},\infty )$ , then, the disturbance ${d_1}(t)$ is rejected, especially it is rejected sufficiently in the frequency band $\left[ {{\omega _{req}},\infty } \right)$ , and the corrector estimate outputs satisfy:

(21) \begin{equation}\mathop {\lim }\limits_{t \to \infty } \left| {{x_1} - {p_0}(t)} \right| \le \rho ({\omega _1}){L_1};\ \mathop {\lim }\limits_{t \to \infty } \left| {{x_2} - {v_0}(t)} \right| \le {L_2}\end{equation}

where, the rejection ratio is expressed by

(22) \begin{equation}\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + x_{\min }^2\omega _1^2} }} + \varepsilon \cdot {x_{\min }}\end{equation}

The rejection ratio $\rho ({\omega _1})$ is a monotonically decreasing function of ${\omega _1} \in \left[ {0,\infty } \right)$ , and it satisfies:

1. (i) In $\left[ {{\omega _{req}},\infty } \right)$ , $\rho ({\omega _1})\;:\;{\rho _{\min }} \to \varepsilon \cdot {x_{\min }}$ as ${\omega _1}\;:\;{\omega _{req}} \to \infty $ .

2. (ii) In $\left( {{\omega _c},{\omega _{req}}} \right)$ , $\rho ({\omega _1})\;:\;1 \to {\rho _{\min }}$ as ${\omega _1}\;:\;{\omega _c} \to {\omega _{req}}$ , where, ${\omega _c} = \frac{{\sqrt {\frac{{2\varepsilon }}{{{x_{\min }}}}} \sqrt {1 - \frac{1}{2}\varepsilon {x_{\min }}} }}{{1 - \varepsilon \cdot {x_{\min }}}}$ .

3. (iii) In $\left[ {0,{\omega _c}} \right]$ , $\rho ({\omega _1})\;:\;1 + \varepsilon \cdot {x_{\min }} \to 1$ as ${\omega _1}\;:\;0 \to {\omega _c}$ .

where,

(23) \begin{equation}{\rho _{\min }} = \frac{1}{{\sqrt {1 + \omega _{req}^2x_{\min }^2} }} + \varepsilon \cdot {x_{\min }}\end{equation}

According to Theorem 4.3 and the system error (123) in the proof of Theorem 5.1, we can get the bounds of the correction errors in (21). This concludes the proof of Theorem 5.2. $\blacksquare $

Remark 5.2 (parametres regulation of corrector (17)):

1. (1) The selection of ${k_3} \gt {L_3}$ , ${k_2} \gt {k_1} + {L_3}$ and ${k_1} = 2{\varepsilon ^{1 - r}}$ (where, $\max \left\{ {\frac{{{L_2}}}{{{L_1}}}} \right\} \le \varepsilon \ll 1$ and $r \in \left( {0,\frac{1}{2}} \right]$ ) makes the corrector stable: ${k_3} \gt {L_3}$ and ${k_2} \gt {k_1} + {L_3}$ make the corrector estimate errors satisfy the convergence law (90); ${k_1} = 2{\varepsilon ^{1 - r}}$ makes the sliding surface stable, and it further reduces the corrector estimate error.

2. (2) The upper-bound of sensor accuracy ratio $0 \lt \varepsilon \ll 1$ and $r \in \left( {0,\frac{1}{2}} \right]$ make the rejection ratio $\rho ({\omega _1})$ in frequency band $\left[ {{\omega _0},\infty } \right)$ and the frequency ${\omega _0}$ sufficiently small, and the sufficiently rejectable frequency ${\omega _1} \in \left[ {{\omega _0},\infty } \right)$ may be in the low/mid/high frequency bands.

3. (3) Furthermore, the selection of parametre r affects the rejection ratio $\rho ({\omega _1})$ and the sufficiently rejectable frequency band $\left[ {{\omega _0},\infty } \right)$ . In fact:

   1. (i) Due to $0 \lt \varepsilon \ll 1$ , the smaller r is, the bigger ${\varepsilon ^r} \in \left( {0,1} \right)$ is, and ${\omega _0}$ decreases. Therefore, the frequency band $\left[ {{\omega _0},\infty } \right)$ becomes relatively wider.

   2. (ii) Conversely, the larger r is, the smaller ${\varepsilon ^r} \in \left( {0,1} \right)$ is, and ${\omega _0}$ increases. Therefore, the frequency band $\left[ {{\omega _0},\infty } \right)$ is reduced.

   3. (iii) Minimum range and minimum value of $\rho ({\omega _1})$ : From $\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + \frac{1}{4}{\varepsilon ^{2r - 2}}\omega _1^2} }} + \frac{1}{2}{\varepsilon ^r}$ , we know that $\rho ({\omega _1}) \in \left( {\frac{1}{2}{\varepsilon ^r},1 + \frac{1}{2}{\varepsilon ^r}} \right]$ when ${\omega _1} \in \left[ {0,\infty } \right)$ . Due to $0 \lt \varepsilon \ll 1$ and $r \in \left( {0,\frac{1}{2}} \right]$ , when we select $r = \frac{1}{2}$ , we can get the minimum value of $\min \left\{ {{\varepsilon ^r}} \right\} = {\varepsilon ^{\frac{1}{2}}}$ . Therefore, the minimum value of $\rho ({\omega _1})$ is $\min \left\{ {\rho ({\omega _1})} \right\} = \frac{1}{2}{\varepsilon ^{\frac{1}{2}}}$ , and the minimum range is $\rho ({\omega _1}) \in \left( {\frac{1}{2}{\varepsilon ^{\frac{1}{2}}},1 + \frac{1}{2}{\varepsilon ^{\frac{1}{2}}}} \right]$ when ${\omega _1} \in \left[ {0,\infty } \right)$ .

6.0 Application to jet UAV navigation and control

An RC-model-based F/A-18 Hornet prototype is used [17], which is shown in Fig. 2, and the forces and torques of UAV are described in Fig. 3.

Figure 2. Jet UAV prototype.

Figure 3. UAV aerodynamic mesh, forces and moments.

6.1 Modelling of jet UAV flight dynamics

For the UAV trajectory tracking control, the modelling is considered in the earth-fixed frame [23]. Furthermore, the constructed model is fit for observer design to estimate system uncertainties.

Let $\Gamma = ({E_x},{E_y},{E_z})$ denote the earth-fixed frame, and $\Lambda = (E_x^b,E_y^b,E_z^b)$ denote the body frame of the UAV. ${\Theta _\Gamma } = [\begin{array}{*{20}{c}}\phi \theta \psi \end{array}{]^T} \in {\Re ^3}$ describes the UAV roll, pitch and yaw angles (Euler angles), and ${\Omega _\Gamma } = [\begin{array}{*{20}{c}}{\dot \phi } {\dot \theta } {\dot \psi }\end{array}{]^T}$ . We use ${s_\theta }$ for $\sin \theta $ and ${c_\theta }$ for $\cos \theta $ . ${R_{\Gamma \Lambda }}$ is the transformation matrix representing the orientation of the body frame $\Lambda $ with respect to the earth-fixed frame $\Gamma $ :

(24) \begin{align}{R_{\Gamma \Lambda }} = \left[ \begin{matrix} {c_\theta }{c_\psi } &&& {{c_\psi }{s_\theta }{s_\phi } - {s_\psi }{c_\phi }} &&& {{c_\psi }{s_\theta }{c_\phi } + {s_\psi }{s_\phi }}\\{c_\theta }{s_\psi } &&& {{s_\psi }{s_\theta }{s_\phi } + {c_\psi }{c_\phi }} &&& {{s_\psi }{s_\theta }{c_\phi } - {c_\psi }{s_\phi }}\\ - {s_\theta } &&& {{s_\phi }{c_\theta }} &&& {{c_\phi }{c_\theta }}\end{matrix} \right]\end{align}

Let $\alpha $ and $\beta $ be the angle-of-attack and the sideslip angle, respectively, we can get

(25) \begin{equation}\alpha = \theta - \mathop {\arctan }\nolimits^{ - 1} (W/U),\beta = \mathop {\arcsin }\nolimits^{ - 1} (V/{V_T})\end{equation}

where, ${\upsilon _\Lambda } = (U,V,W)$ is the linear velocity in body frame $\Lambda $ , and ${V_T} = \sqrt {{U^2} + {V^2} + {W^3}} $ .

Define ${p_\Gamma } = (x,y,z)$ and ${\upsilon _\Gamma } = (\dot x,\dot y,\dot z)$ as the position and velocity vectors of centre of gravity, respectively, relative to the earth-fixed frame $\Gamma $ ; ${\Omega _\Lambda }$ is the angular rate vector of the airframe expressed in the body frame $\Lambda $ ; $m \in \Re $ is the UAV mass, and $J = {\rm{diag}}\{ {J_{xb}}$ , ${J_{yb}}$ , ${J_{zb}}\} \in {\Re ^{3 \times 3}}$ is the UAV inertia matrix. Then, the dynamic equations for the jet UAV subjected to force $F \in {\Re ^3}$ and torque $\tau \in {\Re ^3}$ are given by

(26) \begin{align}{\dot p_\Gamma }& = {\upsilon _\Gamma }\nonumber\\[3pt]m \cdot {\dot \upsilon _\Gamma } & = F + mg{E_z}\nonumber\\[3pt]J \cdot {\dot \Omega _\Lambda } & = - {\Omega _\Lambda } \times \left( {J{\Omega _\Lambda }} \right) + \tau \end{align}

where, ${E_z} = {\left[ {\begin{matrix}0 && 0&& 1\end{matrix}} \right]^T}$ . The relation between the angular rate vector ${\Omega _\Lambda } = {\left[ {\begin{matrix}{{p_\Lambda }}&& {{q_\Lambda }}&& {{r_\Lambda }}\end{matrix}} \right]^T}$ and the Euler angle derivative vector ${\Omega _\Gamma } = [\begin{matrix}{\dot \phi }&& {\dot \theta }&& {\dot \psi }\end{matrix}{]^T}$ is given by

(27) \begin{equation}{\Omega _\Lambda } = \mathbb{Z}{\Omega _\Gamma }\ {\rm{or}}\ {\Omega _\Gamma } = {\mathbb{Z}^{ - 1}}{\Omega _\Lambda }\end{equation}

where,

(28)

The total external force F consists of the thrust ${F_{jet}}$ generated by the jet engine, aerodynamic forces on the fixed wing ${F_w}$ , aerodynamic forces on the fuselage ${F_f}$ , the forces created by the rudders ${F_r}$ , the forces created by the elevators ${F_e}$ , and uncertainties and external disturbances ${F_d}$ . These forces are expressed in body frame $\Lambda $ , and they are transformed by ${R_{\Gamma \Lambda }}$ to be expressed in the earth-fixed frame $\Gamma $ as follows:

(29) \begin{equation}F = {R_{\Gamma \Lambda }}({F_{jet}} + {F_w} + {F_f} + {F_r} + {F_e} + {F_d})\end{equation}

The total moment $\tau $ consists of the moments created by the fixed wings ${\tau _w}$ , the moments created by the rudders ${\tau _r}$ , the moments created by the elevators ${\tau _e}$ , and moments due to the uncertainties and external disturbances ${\tau _d}$ :

(30) \begin{equation}\tau = {\tau _w} + {\tau _r} + {\tau _e} + {\tau _d}\end{equation}

The aerodynamic parametres of the UAV are from [17]. The CFD (computational fluid dynamics) simulation was performed, and the results from the wind tunnel tests were compared. For the boundary conditions, the method of free-stream boundary condition based on Riemann invariants was utilised [28]: no-slip viscous flow condition; the linearised one-dimensional Euler equations; the free-stream values for determination of the value of the Riemann invariants; and the symmetrical boundary condition for the symmetrical UAV. The ANSYS Fluent was used for the CFD simulation, and the design steps included meshing fluid field, fluent simulation and post-processing [29].

1. (1) Thrust by jet engine: The thrust of engine in body frame is expressed by

   (31) \begin{equation}{F_{jet}} = {\left[ {\begin{matrix}{{F_c}} && 0 &&0\end{matrix}} \right]^T}\end{equation}

2. (2) Aerodynamics of fixed wings

Define $q = 0.5\rho ({U^2} + {W^2})$ , where, $\rho $ is the air density. The lift and drag forces generated by the fixed wings are, respectively

(32) \begin{align}{L_i} & = q{S_w}{C_{Li}},{C_{Li}} = {C_{L0}} + {C_{L\alpha }}\alpha + {C_{L{\delta _i}}}{\delta _i}\nonumber\\[3pt]{D_i} & = q{S_w}{C_{Di}},{C_{Di}} = {C_{D0}} + C_{Li}^2/(\pi {A_w}{e_w})\nonumber\\[3pt]{e_w} & = 1.78(1 - 0.045A_w^{0.68}) - 0.46\end{align}

where $i = 1,2$ ; ${S_w}$ is the area of the half wing; ${C_{L0}}$ is the lift coefficient when the angle-of-attack $\alpha $ equals zero; ${C_{L\alpha }}$ is the lift coefficient due to the angle-of-attack $\alpha $ ; ${\delta _i}$ is the aileron deflection, and ${C_{L{\delta _i}}}$ is the lift coefficient due to the aileron deflection ${\delta _i}$ ; ${C_{D0}}$ is the drag coefficient when $\alpha = {\delta _i} = 0$ ; ${A_w}$ is the aspect ratio of the fixed wing; ${e_w}$ is the value of the Oswald’s efficiency factor. The expression of lift and drag coefficients is considered as valid for low angles of attack.

Then the aerodynamic force vector ${F_w}$ on the fixed wings in body frame can be written as

(33) \begin{equation}{F_w} = \left[ {\begin{matrix}{({L_1} + {L_2}){s_\alpha } - ({D_1} + {D_2}){c_\alpha }}\\[3pt]0\\[3pt]{ - ({L_1} + {L_2}){c_\alpha } - ({D_1} + {D_2}){s_\alpha }}\end{matrix}} \right]\end{equation}

The fixed-wing moment ${\tau _w}$ includes the aerodynamic moment ${\tau _{wa}}$ and control torque ${\tau _{wc}}$ around the body axis $E_x^b$ , i.e. ${\tau _w} = {\tau _{wa}} + {\tau _{wc}}$ , where,

(34) \begin{equation}{\tau _{wa}} = \left[ {\begin{matrix}{{l_w}({D_1} - {D_2}){s_\alpha }}\\[3pt]{{l_c}[({L_2} + {L_1}){c_\alpha } + ({D_2} + {D_1}){s_\alpha }]}\\[3pt]{{l_w}[({L_1} - {L_2}){s_\alpha } + ({D_2} - {D_1}){c_\alpha }]}\end{matrix}} \right]\end{equation}

and

(35) \begin{equation}{\tau _{wc}} = \left[ {\begin{array}{*{20}{c}}{{l_w}q{S_w}{C_{L{\delta _{1.2}}}}({\delta _1} - {\delta _2}){c_\alpha }}\\[3pt]0\\[3pt]0\end{array}} \right]\end{equation}

1. (3) Fuselage

The parametres of fuselage lift and drag are described as follows:

(36) \begin{align}{L_f}& = q{S_f}{C_{lf}},{C_{lf}} = {C_{lf\alpha }}\alpha ,\nonumber\\[3pt]{D_f} & = q{S_f}{C_{df}},{C_{df}} = {C_{df0}} + {C_{df\alpha }}\alpha \end{align}

where, ${S_f}$ is the fuselage equivalent cross-sectional area; ${L_f}$ and ${D_f}$ are the lift and drag forces generated by the fuselage, respectively; ${C_{lf}}$ is the lift coefficient; ${C_{df}}$ is the drag coefficient; ${C_{df0}}$ is the constant in the coefficient of drag force. Then the force vector ${F_f}$ on the fuselage in body frame is expressed by

(37) \begin{equation}{F_f} = \left[ {\begin{array}{*{20}{c}}{{L_f}{s_\alpha } - {D_f}{c_\alpha }}\\[3pt]0\\[3pt]{ - {L_f}{c_\alpha } - {D_f}{s_\alpha }}\end{array}} \right]\end{equation}

1. (4) Elevator

The parametres of elevator lift and drag are described as follows:

(38) \begin{align}{L_e} & = q{S_e}{C_{le}},{C_{le}} = {C_{le\alpha }}(\alpha + {\delta _e})\nonumber\\{D_e} & = q{S_e}{C_{de}},{C_{de}} = {C_{de0}} + C_{le}^2/(\pi {A_e}{e_e})\nonumber\\{e_e} & = 1.78(1 - 0.045A_e^{0.68}) - 0.46\end{align}

where ${S_e}$ is the area of the elevator, ${\delta _e}$ is the elevator deflection; ${C_{le\alpha }}$ is the lift coefficient due to the angle-of-attack $\alpha $ and the deflection ${\delta _e}$ ; ${C_{de0}}$ is the drag coefficient when $\alpha + {\delta _e} = 0$ ; ${A_e}$ is the aspect ratio of the elevator; ${e_e}$ is the value of the Oswald’s efficiency factor. Then the force vector ${F_e}$ on the elevator in body frame is expressed by

(39) \begin{equation}{F_e} = \left[ {\begin{array}{*{20}{c}}{{L_e}{s_\alpha } - {D_e}{c_\alpha }}\\0\\{ - {L_e}{c_\alpha } - {D_e}{s_\alpha }}\end{array}} \right]\end{equation}

The elevator moment ${\tau _e}$ includes the aerodynamic moment ${\tau _{ea}}$ and control torque ${\tau _{ec}}$ in the body axis $E_y^b$ , i.e. ${\tau _e} = {\tau _{ea}} + {\tau _{ec}}$ , where,

(40) \begin{equation}{\tau _{ea}} = \left[ {\begin{array}{*{20}{c}}0\\{ - {l_e}{D_e}{s_\alpha }}\\0\end{array}} \right]\end{equation}

and

(41) \begin{equation}{\tau _{ec}} = \left[ {\begin{array}{*{20}{c}}0\\{ - {l_e}q{S_e}{C_{le\alpha }}(\alpha + {\delta _e}){c_\alpha }}\\0\end{array}} \right]\end{equation}

1. (5) Rudders

Define $\bar q = 0.5\rho ({U^2} + {V^2})$ . The lift and drag forces generated by the rudders, respectively

(42) \begin{align}{L_r} & = \bar q{S_r}{C_{lr}},{C_{lr}} = {C_{lr\beta }}\beta + {C_{lr{\delta _r}}}{\delta _r}\nonumber\\{D_r} & = \bar q{S_r}{C_{dv}},{C_{dr}} = {C_{dr0}} + C_{lr}^2/(\pi {A_r}{e_r})\nonumber\\{e_r} & = 1.78(1 - 0.045A_r^{0.68}) - 0.46\end{align}

where, ${S_r}$ is the area of the rudders; ${C_{lr\beta }}$ is the lift coefficient due to the sideslip angle $\beta $ ; ${\delta _r}$ is the rudder deflection, and ${C_{lr{\delta _r}}}$ is the lift coefficient due to the deflection ${\delta _r}$ ; ${C_{dr0}}$ is the drag coefficient when $\beta = {\delta _r} = 0$ ; ${A_r}$ is the aspect ratio of the rudder; ${e_r}$ is the value of the Oswald’s efficiency factor. Then the aerodynamic force vector ${F_r}$ on the rudders in body frame can be expressed by

(43) \begin{equation}{F_r} = \left[ {\begin{array}{*{20}{c}}{{L_r}{s_\beta } - {D_r}{c_\beta }}\\{{L_r}{c_\beta } + {D_r}{s_\beta }}\\0\end{array}} \right]\end{equation}

The rudder moment ${\tau _r}$ includes the aerodynamic moment ${\tau _{ra}}$ and control torque ${\tau _{rc}}$ in body axis $E_z^b$ , i.e. ${\tau _r} = {\tau _{ra}} + {\tau _{rc}}$ , where,

(44) \begin{equation}{\tau _{ra}} = \left[ {\begin{array}{*{20}{c}}0\\0\\{{l_r}{D_r}{s_\beta }}\end{array}} \right]\end{equation}

and

(45) \begin{equation}{\tau _{rc}} = \left[ {\begin{array}{*{20}{c}}0\\0\\{{l_r}\bar q{S_r}({C_{lr\beta }}\beta + {C_{lr{\delta _r}}}{\delta _r}){c_\beta }}\end{array}} \right]\end{equation}

1. (6) UAV motion equations in the earth-fixed frame considering system uncertainties

According to [23], from (27) and (26), we get

(46) \begin{align}{\dot \Omega _\Gamma }& = {\left( {{\mathbb{Z}^{ - 1}}} \right)^{\prime}}{\Omega _\Lambda } + {\mathbb{Z}^{ - 1}}{\dot \Omega _\Lambda }\nonumber\\& = {\left( {{\mathbb{Z}^{ - 1}}} \right)^{\prime}}\mathbb{Z}{\Omega _\Gamma } - {\mathbb{Z}^{ - 1}}{J^{ - 1}}\left[ {\mathbb{Z}{\Omega _\Gamma } \times \left( {J\mathbb{Z}{\Omega _\Gamma }} \right)} \right] + {\mathbb{Z}^{ - 1}}{J^{ - 1}}\tau \end{align}

The total moment $\tau $ can be categorised into the control torque ${\tau _c}$ , uncertain moment ${\tau _d}$ and other moments ${\tau _{other}}$ , i.e.

(47) \begin{equation}\tau = {\tau _c} + {\tau _{other}} + {\tau _d}\end{equation}

From (35), (41) and (45), we get

(48) \begin{equation}{\tau _c} = {\tau _{wc}} + {\tau _{ec}} + {\tau _{rc}}\end{equation}

and from (34), (40) and (44), we get

(49) \begin{equation}{\tau _{other}} = {\tau _{wa}} + {\tau _{ea}} + {\tau _{ra}}\end{equation}

In system (46), considering (47), we define

(50) \begin{equation}{\left( {{\mathbb{Z}^{ - 1}}} \right)^{\prime}}\mathbb{Z}{\Omega _\Gamma } - {\mathbb{Z}^{ - 1}}{J^{ - 1}}\left[ {\mathbb{Z}{\Omega _\Gamma } \times \left( {J\mathbb{Z}{\Omega _\Gamma }} \right)} \right] + {\mathbb{Z}^{ - 1}}{J^{ - 1}}{\tau _{other}}\mathop = \limits^{{\rm{define}}} \left[ {\begin{array}{*{20}{c}}{{a_\phi }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right)}\\[3pt]{{a_\theta }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right)}\\[3pt]{{a_\psi }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right)}\end{array}} \right]\end{equation}

(51) \begin{equation}{\mathbb{Z}^{ - 1}}{J^{ - 1}}{\tau _c}\mathop = \limits^{{\rm{define}}} {\left[ {\begin{matrix}{{\tau _\phi }} &&& {{\tau _\theta }}&&& {{\tau _\psi }}\end{matrix}} \right]^T}\end{equation}

In (29), we define

(52) \begin{equation}{R_{\Gamma \Lambda }}({F_w} + {F_f} + {F_r} + {F_e})/m\mathop = \limits^{{\rm{define}}} \left[ {\begin{array}{*{20}{c}}{{F_{xa}}(\phi ,\theta ,\psi ,\alpha ,\beta )}\\[3pt]{{F_{ya}}(\phi ,\theta ,\psi ,\alpha ,\beta )}\\[3pt]{{F_{za}}(\phi ,\theta ,\psi ,\alpha ,\beta )}\end{array}} \right]\end{equation}

Then, the jet UAV motion equations written in terms of the centre of mass C in the earth-fixed frame are

(53) \begin{align}\ddot x & = {c_\theta }{c_\psi }{F_c}/m + {F_{xa}}(\phi ,\theta ,\psi ,\alpha ,\beta ) + {\Delta _x}\nonumber\\\ddot y & = {c_\theta }{s_\psi }{F_c}/m + {F_{ya}}(\phi ,\theta ,\psi ,\alpha ,\beta ) + {\Delta _y}\nonumber\\\ddot z & = - {s_\theta }{F_c}/m + {F_{za}}(\phi ,\theta ,\psi ,\alpha ,\beta ) + g + {\Delta _z}\end{align}

(54) \begin{align}\ddot \phi& = {a_\phi }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) + {\tau _\phi } + {\Delta _\phi }\nonumber\\\ddot \theta & = {a_\theta }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) + {\tau _\theta } + {\Delta _\theta }\nonumber\\\ddot \psi & = {a_\psi }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) + {\tau _\psi } + {\Delta _\psi }\end{align}

where, m is the mass of the UAV; g is the gravity acceleration; ( ${\Delta _x}$ , ${\Delta _y}$ , ${\Delta _z}$ ) and ( ${\Delta _\phi }$ , ${\Delta _\theta }$ , ${\Delta _\psi }$ ) are the bounded uncertainties in the position and attitude dynamics, respectively; ${\tau _\phi }$ , ${\tau _\theta }$ and ${\tau _\psi }$ are the control torques for roll, pitch and yaw dynamics, respectively, defined in (51).

6.2 Measurements

For the jet UAV, a GPS receiver provides the global position and the velocity. An IMU gives the attitude angle and angular velocity. The measurement outputs are expressed by

(55) \begin{equation}{y_{ * 1}} = * + {d_{ * 1}}(t),{y_{ * 2}} = \dot * + {d_{ * 2}}(t)\end{equation}

where, ${ * _i} = \{ x,y,z,\phi ,\theta ,\psi \} $ , and $\dot * = \{ \dot x,\dot y,\dot z,\dot \phi ,\dot \theta ,\dot \psi \} $ ; ${d_{ * 1}}(t)$ denotes the disturbances in position and angle measurements, and $\mathop {\sup }\nolimits_{t \in [0,\infty )} \left| {{d_{ * 1}}(t)} \right| \le {L_{ * 1}} \lt \infty $ ; ${d_{ * 2}}(t)$ denotes the disturbances in the measurements of flying velocity and angular velocity, and $\mathop {\sup }\nolimits_{t \in [0,\infty )} \left| {{d_{ * 2}}(t)} \right| \le {L_{ * 2}} \lt \infty $ ; ${L_{ * 2}} \ll 1$ , and $\max \left\{ {\frac{{{L_{ * 1}}}}{{{L_{ * 2}}}}} \right\} \le {\varepsilon _ * } \ll 1$ . The corrector (17) is used to reject the measurement disturbances and to estimate the actual (x, y, z, $\phi $ , $\theta $ , $\psi $ ).

6.3 Controller design

In this section, the control laws are derived for UAV trajectory tracking and attitude stabilisation. The position, attitude and system uncertainties are reconstructed by the presented corrector and an existing extended state observer.

The scheme of control system with correction and estimation is shown in Fig. 4: (1) the correctors estimate position and attitude angles, and the disturbances are rejected from the measurements; (2) the extended state observers estimate the uncertainties in position and attitude dynamics, respectively; (3) according to the reference trajectory and estimation from the correctors and observers, the position controller drives the position dynamics; (4) from the estimation of position, attitude, uncertainties and the reference trajectory, the desired attitude is determined; (5) according to the desired attitude and estimation from the correctors and observers, the attitude controller drives the attitude dynamics.

1. (1) Error systems

Figure 4. Scheme of control system.

Suppose the reference trajectory and its finite order derivatives are bounded, and can be generated directly. For the reference trajectory ${X_d} = ({x_d},{y_d},{z_d})$ , let ${e_x} = x - {x_d}$ , ${e_y} = y - {y_d}$ , and ${e_z} = z - {z_d}$ , then the system error for position dynamics (53) is

(56) \begin{equation}{\ddot e_p} = {u_p} + {\Xi _p} + {\Delta _p}\end{equation}

where,

(57) \begin{align}{e_p} & = \left[ {\begin{array}{*{20}{c}}{{e_x}}\\[3pt]{{e_y}}\\[3pt]{{e_z}}\end{array}} \right]\!,\ {u_p} = \left[ {\begin{array}{*{20}{c}}{{u_{px}}}\\[3pt]{{u_{py}}}\\[3pt]{{u_{pz}}}\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}{{c_\theta }{c_\psi }}\\[3pt]{{c_\theta }{s_\psi }}\\[3pt]{ - {s_\theta }}\end{array}} \right]{F_c}/m,\nonumber\\{\Xi _p} & = \left[ {\begin{array}{*{20}{c}}{{F_{xa}}(\phi ,\theta ,\psi ,\alpha ,\beta ) - {{\ddot x}_d}}\\[3pt]{{F_{ya}}(\phi ,\theta ,\psi ,\alpha ,\beta ) - {{\ddot y}_d}}\\[3pt]{{F_{za}}(\phi ,\theta ,\psi ,\alpha ,\beta ) + g - {{\ddot z}_d}}\end{array}} \right]\!,\ {\Delta _p} = \left[ {\begin{array}{*{20}{c}}{{\Delta _x}}\\[3pt]{{\Delta _y}}\\[3pt]{{\Delta _z}}\end{array}} \right]\end{align}

For the desired attitude angle ${\Theta _d} = ({\phi _d},{\theta _d},{\psi _d})$ , let ${e_\phi } = \phi - {\phi _d}$ , ${e_\theta } = \theta - {\theta _d}$ , and ${e_\psi } = \psi - {\psi _d}$ , then the system error for attitude dynamics (54) is

(58) \begin{equation}{\ddot e_a} = {u_a} + {\Xi _a} + {\Delta _a}\end{equation}

where,

(59) \begin{align}{e_a} & = \left[ {\begin{array}{*{20}{c}}{{e_\phi }}\\[3pt]{{e_\theta }}\\[3pt]{{e_\psi }}\end{array}} \right]\!,\ {u_a} = \left[ {\begin{array}{*{20}{c}}{{\tau _\phi }}\\[3pt]{{\tau _\theta }}\\[3pt]{{\tau _\psi }}\end{array}} \right]\!,\ {\Delta _a} = \left[ {\begin{array}{*{20}{c}}{{\Delta _\phi }}\\[3pt]{{\Delta _\theta }}\\[3pt]{{\Delta _\psi }}\end{array}} \right]\!,\nonumber \\{\Xi _a} & = \left[ {\begin{array}{*{20}{c}}{{a_\phi }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) - {{\ddot \phi }_d}}\\[3pt]{{a_\theta }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) - {{\ddot \theta }_d}}\\[3pt]{{a_\psi }\left( {\phi ,\theta ,\psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) - {{\ddot \psi }_d}}\end{array}} \right]\end{align}

1. (2) Extended state observers for the uncertainties ${\Delta _p}$ and ${\Delta _a}$

The extended state observers [24] are used to estimate the uncertainties in position dynamics (53) and attitude dynamics (54). For the observers, the measurements of flying velocity and angular velocity ${y_{ * 2}} = \dot * + {d_{ * 2}}(t)$ (where, $\dot * = \{ \dot x,\dot y,\dot z,\dot \phi ,\dot \theta ,\dot \psi \} $ ) are the input signals. Then, the uncertainty ${\Delta _p} = {\left[ {\begin{array}{*{20}{c}}{{\Delta _x}} {{\Delta _y}} {{\Delta _z}}\end{array}} \right]^T}$ in the position dynamics and the uncertainty ${\Delta _a} = {\left[ {\begin{array}{*{20}{c}}{{\Delta _\phi }} {{\Delta _\theta }} {{\Delta _\psi }}\end{array}} \right]^T}$ in the attitude dynamics are estimated.

The continuous extended state observers can provide smooth and accurate estimations of the uncertainties in the system dynamics, reducing high-frequency vibrations.

The following extended state observers are used [24]:

(60) \begin{align}{\dot x_{ * 1}}& = {x_{ * 2}} - {\lambda _{ * 1}}{\left| {{x_{ * 1}} - {y_{ * 2}}} \right|^{\frac{{1 + {\alpha _ * }}}{2}}}sign({x_{ * 1}} - {y_{ * 2}}) + {H_ * }\nonumber\\{\dot x_{ * 2}}& = - {\lambda _{ * 2}}{\left| {{x_{ * 1}} - {y_{ * 2}}} \right|^{{\alpha _ * }}}sign({x_{ * 1}} - {y_{ * 2}})\end{align}

From Theorem 1 in, [24] there exist a finite time ${t_s} \gt 0$ such that, for $t \geqslant {t_s}$ ,

(61) \begin{equation}{x_{ * 1}} = \dot * ,{x_{ * 2}} = {\Delta _ * }\end{equation}

where, $ * = \left\{ {x,y,z,\psi ,\theta ,\phi } \right\}$ and $\dot * = \{ \dot x,\dot y,\dot z,\dot \phi ,\dot \theta ,\dot \psi \} $ ; ${\lambda _{ * 1}},{\lambda _{ * 2}} \gt 0$ , and ${\alpha _ * } \in (0,1)$ ; The measurement ${y_{ * 2}}$ defined in (55) is the observer input signal;

(62) \begin{align}{H_x} & = {u_{px}} + {F_{xa}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta )\nonumber\\{H_y} & = {u_{py}} + {F_{ya}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta )\nonumber\\{H_z} & = {u_{pz}} + {F_{za}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) - g\nonumber\\{H_\phi } & = {\tau _\phi } + {a_\phi }\left( {\widehat \phi ,\widehat \theta ,\widehat \psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right)\nonumber\\{H_\theta } & = {\tau _\theta } + {a_\theta }\left( {\widehat \phi ,\widehat \theta ,\widehat \psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right)\nonumber\\{H_\psi } & = {\tau _\psi } + {a_\psi }\left( {\widehat \phi ,\widehat \theta ,\widehat \psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right)\end{align}

and ( $\widehat \phi $ , $\widehat \theta $ , $\widehat \psi $ ) are from the outputs of the sliding mode correctors. From (60) and (61), we get

(63) \begin{equation}{\widehat \Delta _p} = {\left[ {\begin{array}{*{20}{c}}{{x_{x2}}} {{x_{y2}}} {{x_{z2}}}\end{array}} \right]^T}\end{equation}

and

(64) \begin{equation}{\widehat \Delta _a} = {\left[ {\begin{array}{*{20}{c}}{{x_{\psi 2}}} {{x_{\theta 2}}} {{x_{\phi 2}}}\end{array}} \right]^T}\end{equation}

1. (3) Controller design for position dynamics

For position dynamics (53), to track the reference trajectory ${X_d} = ({x_d},{y_d},{z_d})$ , when we select the following controller, the position system error (56) will converge asymptotically to zero:

(65) \begin{equation}{u_p} = - {\widehat \Xi _p} - {\widehat \Delta _p} - {k_{p1}}{\widehat e_p} - {k_{p2}}{\widehat {\dot e}_p}\end{equation}

where ${\widehat e_x} = \widehat x - {x_d}$ , ${\widehat {\dot e}_x} = \widehat {\dot x} - {\dot x_d}$ , ${\widehat e_y} = \widehat y - {y_d}$ , ${\widehat {\dot e}_y} = \widehat {\dot y} - {\dot y_d}$ , ${\widehat e_z} = \widehat z - {z_d}$ , ${\widehat {\dot e}_z} = \widehat {\dot z} - {\dot z_d}$ ; ${k_{p1}},{k_{p2}} \gt 0$ ; ( $\widehat x$ , $\widehat y$ , $\widehat z$ , $\widehat {\dot x}$ , $\widehat {\dot y}$ , $\widehat {\dot z}$ , $\widehat \phi $ , $\widehat \theta $ , $\widehat \psi $ ) are from the outputs of the correctors; ${\widehat \Delta _p}$ is from the outputs of the extended state observer; and

(66) \begin{align}{\widehat e_p} & = \left[ {\begin{array}{*{20}{c}}{{{\widehat e}_x}}\\[3pt]{{{\widehat e}_y}}\\[3pt]{{{\widehat e}_z}}\end{array}} \right]\!,\ {\widehat {\dot e}_p} = \left[ {\begin{array}{*{20}{c}}{{{\widehat {\dot e}}_x}}\\[3pt]{{{\widehat {\dot e}}_y}}\\[3pt]{{{\widehat {\dot e}}_z}}\end{array}} \right]\!,\ \nonumber\\[6pt]{\widehat \Xi _p} & = \left[ {\begin{array}{*{20}{c}}{{F_{xa}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) - {{\ddot x}_d}}\\[3pt]{{F_{ya}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) - {{\ddot y}_d}}\\[3pt]{{F_{za}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) - g - {{\ddot z}_d}}\end{array}} \right]\end{align}

From (57), we get the engine thrust

(67) \begin{equation}{F_c} = m{\left\| {{u_p}} \right\|_2} = m\sqrt {u_{px}^2 + u_{py}^2 + u_{pz}^2} \end{equation}

1. (4) Desired attitude angles

From (56) and (65), we get

(68) \begin{equation}{\ddot e_z} = - {s_\theta }{F_c}/m + {F_{za}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) + g - {\ddot z_d} + {\widehat \Delta _z} = - {k_{p1}}{\widehat e_z} - {k_{p2}}{\widehat {\dot e}_z}\end{equation}

Then, the desired pitch angle is

(69) \begin{equation}{\theta _d} = \arcsin \frac{{m({F_{za}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) + g - {{\ddot z}_d} + {{\widehat \Delta }_z} + {k_{p1}}{{\widehat e}_z} + {k_{p2}}{{\widehat {\dot e}}_z})}}{{{F_c}}}\end{equation}

Also, from (56) and (65), we get

(70) \begin{align}{\ddot e_x} & = {c_\theta }{c_\psi }{F_c}/m + {F_{xa}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) + {\widehat \Delta _x} = - {k_{p1}}{\widehat e_x} - {k_{p2}}{\widehat {\dot e}_x}\nonumber\\{\ddot e_y} & = {c_\theta }{s_\psi }{F_c}/m + {F_{ya}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) + {\widehat \Delta _y} = - {k_{p1}}{\widehat e_y} - {k_{p2}}{\widehat {\dot e}_y}\end{align}

Then, the desired yaw angle is

(71) \begin{equation}{\psi _d} = \arctan \frac{{{F_{ya}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) + {{\widehat \Delta }_y} + {k_{p1}}{{\widehat e}_y} + {k_{p2}}{{\widehat {\dot e}}_y}}}{{{F_{xa}}(\widehat \phi ,\widehat \theta ,\widehat \psi ,\alpha ,\beta ) + {{\widehat \Delta }_x} + {k_{p1}}{{\widehat e}_x} + {k_{p2}}{{\widehat {\dot e}}_x}}}\end{equation}

The combination of lift force from the wings, elevator and fuselage can provide the centripetal force. The radius of curvature is

(72) \begin{equation}r = \frac{{{{\dot x}^2} + {{\dot y}^2}}}{{\left| {\dot x\ddot y - \ddot x\dot y} \right|}}\end{equation}

and the centripetal force is

(73) \begin{equation}{f_{centri}} = \frac{{m({{\dot x}^2} + {{\dot y}^2})}}{r} = m\left| {\dot x\ddot y - \ddot x\dot y} \right|\end{equation}

Also, the centripetal force can be expressed by

(74) \begin{equation}{f_{centri}} = {L_{wef}}\cos {\phi _d}\end{equation}

where, ${L_{wef}} = - ({L_1} + {L_2}){c_\alpha } - ({D_1} + {D_2}){s_\alpha } - {L_e}{c_\alpha } - {D_e}{s_\alpha } - {L_f}{c_\alpha } - {D_f}{s_\alpha }$ . Then, we get the desired roll angle

(75) \begin{equation}{\phi _d} = \arccos \frac{{{f_{centri}}}}{{{L_{wef}}}}\end{equation}

1. (5) Controller design for attitude dynamics

For attitude dynamics (54), to track the desired attitude ${\Theta _d} = ({\psi _d},{\theta _d},{\phi _d})$ , when we select the following controller, the attitude system error (58) will converge asymptotically to zero:

(76) \begin{equation}{u_a} = - {\widehat \Xi _a} - {\widehat \Delta _a} - {k_{a1}}{\widehat e_a} - {k_{a2}}{\widehat {\dot e}_a}\end{equation}

where, ${k_{a1}},{k_{a2}}$ ; ${\widehat e_\phi } = \widehat \phi - {\phi _d}$ , ${\widehat {\dot e}_\phi } = \widehat {\dot \phi } - {\dot \phi _d}$ , ${\widehat e_\theta } = \widehat \theta - {\theta _d}$ , ${\widehat {\dot e}_\theta } = \widehat {\dot \theta } - {\dot \theta _d}$ , ${\widehat e_\psi } = \widehat \psi - {\psi _d}$ , ${\widehat {\dot e}_\psi } = \widehat {\dot \psi } - {\dot \psi _d}$ ; ${\widehat e_a} = {\left[ {\begin{matrix}{{{\widehat e}_\phi }}&& {{{\widehat e}_\theta }} &&{{{\widehat e}_\psi }}\end{matrix}} \right]^T}$ ; ${\dot e_a} = {\left[ {\begin{matrix}{{{\widehat {\dot e}}_\phi }}&& {{{\widehat {\dot e}}_\theta }}&& {{{\widehat {\dot e}}_\psi }}\end{matrix}} \right]^T}$ ; ( $\widehat \phi $ , $\widehat \theta $ , $\widehat \psi $ ) are from the outputs of the correctors; ${\widehat \Delta _a}$ is from the outputs of the extended state observers; and

(77) \begin{equation}{\widehat \Xi _a} = \left[ {\begin{array}{*{20}{c}}{{a_\phi }\left( {\widehat \phi ,\widehat \theta ,\widehat \psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) - {{\ddot \phi }_d}}\\[3pt]{{a_\theta }\left( {\widehat \phi ,\widehat \theta ,\widehat \psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) - {{\ddot \theta }_d}}\\[3pt]{{a_\psi }\left( {\widehat \phi ,\widehat \theta ,\widehat \psi ,\dot \phi ,\dot \theta ,\dot \psi ,\alpha ,\beta } \right) - {{\ddot \psi }_d}}\end{array}} \right]\end{equation}

7.0 Simulation examples

We use two examples to illustrate the sliding mode corrector presented in Theorems 5.1 and 5.2.

Example 1 (sliding mode corrector design from Theorem 5.1):

1. (1) Sensor outputs and actual values

Measurement signals for position and velocity: ${p_m}(t) = {p_0}(t) + {d_1}(t)$ , ${v_m}(t) = {v_0}(t) + {d_2}(t)$

where, ${d_1}(t)$ and ${d_2}(t)$ are the disturbances in position and velocity measurements, respectively;

Suppose the actual position: ${p_0}\left( t \right) = 10 + 20\sin \left( t \right)$ ;

and the actual velocity: ${v_0}\left( t \right) = 20\cos \left( t \right)$ .

1. (2) Disturbance in position measurement

Position sensing disturbance ${d_1}(t) = {d_{11}}(t) + {d_{12}}(t)$ includes disturbance ${d_{11}}(t)$ and stochastic noise ${d_{12}}(t)$ . We consider the following three types of disturbance ${d_{11}}(t)$ :

1. (a) ${d_{11}}(t) = 2\sin (4t) + \cos (9t)$ ;

2. (b) ${d_{11}}(t) = 6\sin (4t) + 3\cos (9t)$ , and stochastic noise $3{d_{12}}(t)$ (disturbance magnitude increases, i.e. ${L_1}$ increases)

3. (c) ${d_{11}}(t) = 1.5\cos (0.2t) + 0.5\sin (0.6t) + 2\sin (4t) + \cos (9t)$ (disturbance in low and mid frequency bands)

1. (3) Disturbance in velocity measurement

Velocity sensing disturbance ${d_2}(t) = {d_{21}}(t) + {d_{22}}(t) + {d_{23}}(t)$ includes time-varying disturbance ${d_{21}}(t)$ , constant disturbance ${d_{22}}(t)$ and stochastic noise ${d_{23}}(t)$ . We suppose ${d_{21}}(t) = 0.05\cos (0.3t) + 0.03\sin (0.6t)$ and ${d_{22}}(t) = 0.02$ .

1. (4) Determination of upper-bound of sensor accuracy ratio

From the sensor accuracy, we can get $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_1}(t)} \right| \le {L_1} = 3\,(m)$ , $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_2}(t)} \right| \le {L_2} = 0.1\,(m/s)$ . The upper-bound of sensor accuracy ratio $\varepsilon = 0.034$ can be selected such that the sensor accuracy inequality $\frac{{{L_2}}}{{{L_1}}} = \frac{{0.1}}{3} \le \varepsilon = 0.034 \ll 1$ holds.

1. (5) Corrector parametres’ selection

From the corrector parametre selection in the Remark 5.2 of Theorem 5.1, we select $r = \frac{1}{2}$ to get the minimum value and minimum range of rejection ratio $\rho ({\omega _1})$ : $\min \left\{ {\rho ({\omega _1})} \right\} = \frac{1}{2}{\varepsilon ^{\frac{1}{2}}}$ , and $\rho ({\omega _1}) \in \left( {\frac{1}{2}{\varepsilon ^{\frac{1}{2}}},1 + \frac{1}{2}{\varepsilon ^{\frac{1}{2}}}} \right]$ when ${\omega _1} \in \left[ {0,\infty } \right)$ .

Therefore, from $\varepsilon = 0.034$ and $r = \frac{1}{2}$ , we can determine the corrector parametre ${k_1} = 2{\varepsilon ^{1 - r}} = 2 \times {0.034^{1 - 0.5}} = 0.36$ .

According to $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{{\dot v}_0}(t)} \right| \le {L_3} = 10$ , ${k_2} \gt {k_1} + {L_3}$ and ${k_3} \gt {L_3}$ , we select ${k_2} = 100$ , ${k_3} = 100$ .

1. (6) Rejection ratio and disturbance frequency bands

From $\varepsilon = 0.034$ and $r = \frac{1}{2}$ , the rejection ratio is expressed by

$$\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + \frac{1}{4}{\varepsilon ^{2r - 2}}\omega _1^2} }} + \frac{1}{2}{\varepsilon ^r} = \frac{1}{{\sqrt {1 + 7.35\omega _1^2} }} + 0.09$$

and we get:

\begin{align*}{\omega _0}& = 4{\varepsilon ^{1 - 2r}}\sqrt {1 - \frac{1}{4}{\varepsilon ^{2r}}} = 4\sqrt {1 - \frac{1}{4} \times 0.034} = 3.98\,({\rm{rad}}/{\rm{s}})\\ {\omega _c}& = \frac{{{\varepsilon ^{1 - \frac{1}{2}r}}\sqrt {4 - {\varepsilon ^r}} }}{{1 - \frac{1}{2}{\varepsilon ^r}}} = \frac{{{{0.034}^{1 - \frac{1}{4}}}\sqrt {4 - {{0.034}^{0.5}}} }}{{1 - \frac{1}{2} \times {{0.034}^{0.5}}}} = 0.17\,({\rm{rad}}/{\rm{s}})\end{align*}

Therefore, the rejection ratio in the different frequency bands of $\left[ {0,\infty } \right)$ can be described by:

1. (i) In $\left[ {3.98{\rm{rad/s}},\infty } \right)$ , $\rho ({\omega _1})\;:\;0.18 \to 0.09$ as ${\omega _1}\;:\;3.98 \to \infty \,({\rm{rad}}/{\rm{s}})$ ;

2. (ii) In $\left( {0.17,3.98{\rm{rad/s}}} \right)$ , $\rho ({\omega _1})\;:\;1 \to 0.18$ as ${\omega _1}\;:\;0.17 \to 3.98\,({\rm{rad}}/{\rm{s}})$ ;

3. (iii) In $\left[ {0,0.17{\rm{rad/s}}} \right]$ , $\rho ({\omega _1})\;:\;1.09 \to 1$ (i.e. $\rho ({\omega _1}) \approx 1$ ) as ${\omega _1}\;:\;0 \to 0.17\,({\rm{rad}}/{\rm{s}})$ .

Comparison with signal fusion based on Kalman filter:

We compare the sliding mode corrector with KF-based method. For this example, there are only two measurement signals, and no system model is given, we use the direct KF-based signal fusion method [25].

The position measurement is ${p_m}(t) = {p_0}(t) + {d_1}(t)$ , and the velocity measurement is ${v_m}(t) = {v_0}(t) + {d_2}(t)$ . According to the Taylor’s expansion, position and velocity in discrete system can be expressed approximately by

(78) \begin{align}{p_0}(k) &\approx {p_0}(k - 1) + {v_0}(k - 1) \cdot \Delta T\nonumber\\{v_0}(k) &\approx {v_0}(k - 1)\end{align}

where, $\Delta T$ is the sampling time, and k is the sample step.

Define $X\left( k \right) = {\left[ {\begin{matrix}{{p_0}(k)}&& {{v_0}(k)}\end{matrix}} \right]^T}$ , $A = \left[ {\begin{matrix}1 && {\Delta T}\\0&& 1\end{matrix}} \right]$ , the above relation can be described by a matrix system:

(79) \begin{equation}X\left( k \right) = A \cdot X\left( {k - 1} \right)\end{equation}

For ${p_m}(t) = {p_0}(t) + {d_1}(t)$ and ${v_m}(t) = {v_0}(t) + {d_2}(t)$ , we get

(80) \begin{equation}{p_m}(k) = {p_0}(k) + {d_1}(k);\ {v_m}(k) = {v_0}(k) + {d_2}(k)\end{equation}

Define $H = \left[ {\begin{matrix}1 &&& 0\\0&&& 1\end{matrix}} \right]$ , $V\left( k \right) = {\left[ {\begin{matrix}{{d_1}(k)} &&{{d_2}(k)}\end{matrix}} \right]^T}$ , then the measurement outputs can be expressed by

(81) \begin{equation}Y\left( k \right) = H \cdot X\left( k \right) + V\left( k \right)\end{equation}

Therefore, the Kalman filter for signal integration is designed as follows:

(82) \begin{align}X(k|k - 1) & = A \cdot X(k - 1|k - 1)\nonumber\\P(k|k - 1) & = A \cdot P(k - 1|k - 1){A^T} + Q\nonumber\\K(k) & = \frac{{P(k|k - 1){H^T}}}{{H \cdot P(k|k - 1){H^T} + R}}\nonumber\\X(k|k) & = X(k|k - 1) + K(k)(Y(k) - H \cdot X(k|k - 1))\nonumber\\P(k|k) & = (I - K(k)H)P(k|k - 1)\end{align}

where, $I = \left[ {\begin{matrix}1&& 0\\0&& 1\end{matrix}} \right]$ ; $Q = \left[ {\begin{matrix}{\frac{1}{3}{q_c}\Delta {T^3}}&& {\frac{1}{2}{q_c}\Delta {T^2}}\\{\frac{1}{2}{q_c}\Delta {T^2}}&& {{q_c}\Delta T}\end{matrix}} \right]$ is the process noise covariance matrix, and ${q_c}$ is the power spectral density of the input white noise; and R is the measurement noise covariance matrix. In the simulation, the power spectral density of the input white noise is selected as ${q_c} = 1$ , [25] $R = \left[ {\begin{matrix}{0.8}&& 0\\0&& {0.8}\end{matrix}} \right]$ , $P\left( {0|0} \right) = \left[ {\begin{matrix}{0.1}&& 0\\0 && {0.1}\end{matrix}} \right]$ and the sampling time is $\Delta T = 0.008$ (sec).

Analysis of simulation results:

The disturbance rejections in position sensing are presented in Fig. 5. Figure 5(a) describes the disturbance rejection when the disturbance ${d_1}(t)$ is in the frequency band that is rejected sufficiently. In fact, the frequency band of disturbance ${d_1}(t)$ is ${\omega _1} \in \left[ {4{\rm{rad/s}},\infty } \right)$ , and $\left[ {4{\rm{rad/s}},\infty } \right) \subset \left[ {3.98{\rm{rad/s}}, + \infty } \right)$ . We know that in frequency band $\left[ {3.98{\rm{rad/s}},\infty } \right)$ , the minimum rejection ratio is obtained, i.e. $\rho ({\omega _1}) \in \left( {\frac{1}{2}{\varepsilon ^{\frac{1}{2}}},{\varepsilon ^{\frac{1}{2}}}} \right] = \left( {0.09,0.18} \right]$ . Therefore, the position disturbance is rejected sufficiently.

Figure 5. Example 1 – Simulation on position disturbance rejection. (a) Disturbance rejection when ${L_1} = 3$ . (b) Disturbance rejection when ${L_1} = 9$ . (c) Disturbance rejection when low-frequency disturbance is also included.

Figure 5(b) presents the disturbance rejection when the magnitude of disturbance ${d_1}(t)$ increases. Even the position sensing accuracy becomes worse ( ${L_1}\;:\;3\,{\rm{m}}$ $ \to 9\,{\rm{m}}$ ), the corrector with the original parametres can still reject the disturbance sufficiently.

Figure 5(c) shows the disturbance rejection when the position disturbance ${d_1}(t)$ covers the low/min/high frequency bands. In the disturbance ${d_1}(t)$ , the mid/high frequency part is rejected sufficiently. Even the low-frequency disturbance exists, the effect of disturbance can still be rejected largely. In fact, in the position disturbance ${d_1}(t) = {d_{11}}(t) + {d_{12}}(t)$ (where, ${d_{11}}(t) = 1.5\cos (0.2t) + 0.5\sin (0.6t) + 2\sin (4t) + \cos (9t)$ and the high-frequency stochastic noise ${d_{12}}(t)$ ):

1. (1) the part $2\sin (4t) + \cos (9t)$ is within the frequency band $\left[ {3.98{\rm{rad/s}},\infty } \right)$ , and the rejection ratio is minimum, i.e. $\rho ({\omega _1}) \in \left( {0.09,0.18} \right]$ ; therefore, this part of disturbance is rejected sufficiently; ${d_{12}}(t)$ is also rejected sufficiently due to its high frequency.

2. (2) the part $1.5\cos (0.2t) + 0.5\sin (0.6t)$ is within the frequency band $\left( {0.17,3.98{\rm{rad/s}}} \right]$ , and the rejection ratio is $\rho ({\omega _1}) \in \left[ {0.18,1{\rm{rad/s}}} \right)$ ; therefore, this part of disturbance is still reduced.

From Fig. 5, we can also find that the estimate outputs of corrector are accurate and smoothed even stochastic noise exists in position and velocity measurements. In addition, if initial calibration is done for position sensing, the corrector error keeps small from the beginning; and if there is no initial calibration for position sensing, the corrector error can still converge to the small bound. The corrector performance is compared with the estimate results of the KF-based method. Comparing to the corrector, due to the existence of widely frequency-band disturbance in position sensing, the obviously large estimate errors exist in the outputs of the KF, although it can reduce the effect of disturbance to some extent.

Example 2 (sliding mode corrector design from Theorem 5.2):

This example illustrates the position disturbance ${d_1}(t)$ is rejected to the maximum extent within a given frequency band, and the disturbance in the other bands can still be reduced largely. The corrector parametres are determined to get an optimal value of rejection ratio in the given frequency band.

1. (1) Sensor outputs and actual values

Sensing signals for position and velocity are ${p_m}(t) = {p_0}(t) + {d_1}(t)$ and ${v_m}(t) = {v_0}(t) + {d_2}(t)$ , respectively, where, ${d_1}(t)$ and ${d_2}(t)$ are the disturbances in position and velocity sensing, respectively;

the actual position: ${p_0}\left( t \right) = 10 + 20\sin \left( t \right)$ ;

and the actual velocity: ${v_0}\left( t \right) = 20\cos \left( t \right)$ .

1. (2) Disturbance in position sensing

Position sensing disturbance ${d_1}(t) = {d_{11}}(t) + {d_{12}}(t)$ includes time-varying disturbance ${d_{11}}(t)$ and stochastic noise ${d_{12}}(t)$ , and we suppose ${d_{11}}(t) = 2\sin (1.5) + \cos (3t)$ .

1. (3) Disturbance in velocity sensing

Velocity sensing disturbance ${d_2}(t) = {d_{21}}(t) + {d_{22}}(t) + {d_{23}}(t)$ includes time-varying disturbance ${d_{21}}(t)$ , constant disturbance ${d_{22}}(t)$ , and stochastic noise ${d_{23}}(t)$ . We suppose ${d_{21}}(t) = 0.05\cos (0.3t) + 0.03\sin (0.6t)$ and ${d_{22}}(t) = 0.02$ .

1. (4) Corrector parametres’ selection

In this example, we suppose the disturbance ${d_1}(t)$ is required to be rejected sufficiently in the given frequency band $\left[ {{\omega _{req}},\infty } \right) = \left[ {1{\rm{rad/s}},\infty } \right)$ , and a small rejection ratio is obtained in this frequency band. From ${\omega _{req}} = 1$ rad/s, we can get the unique solution ${x_{\min }}$ to

$${1^2}{x_{\min }}{(1 + {0.5^2}x_{\min }^2)^{ - \frac{3}{2}}} - 0.034 = 0$$

in the range $(\frac{1}{{\sqrt 2 \times 1}},\infty )$ , i.e. ${x_{\min }} = 5.26$ , and ${k_1} = 1/{x_{\min }} = 0.19$ .

According to $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{{\dot v}_0}(t)} \right| \le {L_3} = 20$ , ${k_2} \gt {k_1} + {L_3}$ and ${k_3} \gt {L_3}$ , we select ${k_2} = 100$ , ${k_3} = 100$ .

1. (5) Rejection ratio and disturbance frequency bands

From the sensor accuracy ${L_1} = 3\,({\rm{m}})$ and ${L_2} = 0.1\,({\rm{m}}/{\rm{s}})$ , the upper-bound of sensor accuracy ratio $\varepsilon = 0.034$ is selected such that the sensor accuracy inequality $\frac{{{L_2}}}{{{L_1}}} = \frac{{0.1}}{3} \le \varepsilon = 0.034 \ll 1$ holds.

From $\varepsilon = 0.034$ and ${x_{\min }} = 5.26$ , the rejection ratio can be expressed by

$$\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + x_{\min }^2\omega _1^2} }} + \varepsilon \cdot {x_{\min }} = \frac{1}{{\sqrt {1 + 27.67\omega _1^2} }} + 0.18$$

Therefore, the rejection ratio in the different frequency bands of $\left[ {0,\infty } \right)$ can be described by:

1. (i) In $\left[ {1{\rm{rad/s}},\infty } \right)$ , $\rho ({\omega _1})\;:\;0.37 \to 0.18$ as ${\omega _1}\;:\;1 \to \infty \,({\rm{rad}}/{\rm{s}})$ ;

2. (ii) In $\left( {0.13,1{\rm{rad/s}}} \right)$ , $\rho ({\omega _1})\;:\;1 \to 0.37$ as ${\omega _1}\;:\;0.13 \to 1\,({\rm{rad}}/{\rm{s}})$ , in which, ${\omega _c} = \frac{{\sqrt {\frac{{2\varepsilon }}{{{x_{\min }}}}} \sqrt {1 - \frac{1}{2}\varepsilon {x_{\min }}} }}{{1 - \varepsilon \cdot {x_{\min }}}} = 0.13\,({\rm{rad}}/{\rm{s}})$ ;

3. (iii) In $\left[ {0,0.13{\rm{rad/s}}} \right]$ , $\rho ({\omega _1})\;:\;1.18 \to 1$ as ${\omega _1}\;:\;0 \to 0.13\,({\rm{rad}}/{\rm{s}})$ .

The disturbance rejection for Example 2 is presented in Fig. 6. Due to the relatively small rejection ratio $\rho ({\omega _1}) \in \left( {0.18,0.37} \right]$ when ${\omega _1} \in \left[ {1{\rm{rad/s}},\infty } \right)$ , the disturbance ${d_1}\left( t \right)$ in position sensing is rejected sufficiently in the given frequency band $\left[ {1{\rm{rad/s}},\infty } \right)$ . Also, we can find that the corrector has strong robustness against stochastic noise from the measurements of position and velocity, and the estimate output is smoothed. Thus, the plot performance confirms the results of numerical calculation.

Figure 6. Example 2 – Simulation on position disturbance rejection for the given frequency band.

For example 1, when ${\omega _1} = 1\,{\rm{rad}}/{\rm{s}}$ , we can get $\rho ({\omega _1}) = \frac{1}{{\sqrt {1 + 7.35\omega _1^2} }} + 0.09 = 0.44\,{\rm{rad}}/{\rm{s}}$ . Therefore, in the given frequency band $\left[ {1{\rm{rad/s}},\infty } \right)$ , the rejection ratio $\rho ({\omega _1})$ it satisfies $\rho ({\omega _1})\;:\;0.44 \to 0.09$ as ${\omega _1}\;:\;1 \to \infty \,({\rm{rad}}/{\rm{s}})$ . The rejection performance of Example 1 is a little worse than that of Example 2 near the frequency ${\omega _1}$ due to $0.44 \gt 0.37$ . However, in low frequency band, the method in Example 1 can get better performance because $\rho ({\omega _1})\;:\;1.09 \to 1$ (i.e. $\rho ({\omega _1}) \approx 1$ ) as ${\omega _1}\;:\;0 \to 0.17\,({\rm{rad}}/{\rm{s}})$ ; while, for Example 2, we can get $\rho ({\omega _1})\;:\;1.18 \to 1$ as ${\omega _1}\;:\;0 \to 0.13\,({\rm{rad}}/{\rm{s}})$ . For the whole frequency range $\left[ {0,\infty } \right)$ , the corrector in Example 1 can get $\rho ({\omega _1})\;:\;1.09 \to 0.09$ as ${\omega _1}\;:\;0 \to \infty $ , while the corrector in Example 2 can get $\rho ({\omega _1})\;:\;1.18 \to 0.18$ as ${\omega _1}\;:\;0 \to \infty $ .

8.0 Experiment on jet UAV navigation and control

In this section, an experiment on a jet UAV flight is presented. The jet UAV prototype (an RC-model-based F/A-18 Hornet) shown in Fig. 2 is used for the flight test [17]. A JetCat P200-SX jet engine is adopted to provide the thrust, and the engine starter includes: Jet-tronic ECU for fuel control; electronic valve; electronic starting gas valve; electronic fuel valve; fuel tubing, tubing connector set, filters and cable set; 2-cell, 3,300mA LiPoly battery pack; and starting gas tank. The engine can provide 220N (52 lbs) thrust for 112,000 RPM, and RPM range: $33000 \sim 112000$ RPM. A Gumstix microcomputer is used for data collection and signal processing from sensors. The flight control system implementation on the hardware is shown in Fig. 7. An Arduino Mega 2560 is taken as the driver board, which has multiple PWM output channels. The input voltage is $7 \sim 12\,{\rm{V}}$ . The control update time is 5ms. The FUTABA S3001 servos are adopted to control the deflections of ailerons, elevators and rudders. A 10Hz GPS MediaTek MT3329 is selected as the GPS receiver. A 9Hz VTI SCP1000 altimetre with 10cm resolution is utilised for above the sea level altitude measurements at higher altitudes. A 12Hz SF02-F laser altimetre is used for altitude measurements at lower altitudes with 40m range. A 10kHz Xsens MTI AHRS provides the 3-axial accelerations, the angular rates and the earth’s magnetic field. A 192kHz kpilot 32 digital air speed sensor is utilised to obtain the relative wind speed. A 100kHz 4239-01 AOA sensor is used to measure the angle-of-attack.

Figure 7. Control system hardware.

The jet UAV parametres are described in Table 1.

Table 1. UAV parametres

Measurements of position and velocity: A 10Hz GPS MediaTek MT3329 (without aid) provides position at accuracy of 3m and velocity at accuracy of 0.1m/s.

Measurements of attitude angle and angular velocity: A 10kHz Xsens MTI AHRS provides attitude angles and angular velocity, in which: roll/pitch accuracy: ${0.5^ \circ }$ , yaw accuracy: ${1.0^ \circ }$ ; angular accuracy: 9º/hr = 0.0025º/s.

Desired flight trajectory: The desired flight trajectory consists of takeoff, climb, cruise in a circle with the radius 500m and the height 300 m, and landing back, which is shown in Fig. 8(a).

Figure 8. UAV flight based on correction. (a) 3D navigation trajectories. (b) Position comparison in the three directions.

In the experiment, considering measurement disturbances and the uncertainties in the UAV flight dynamics, the jet UAV is controlled to track the reference trajectory. The position and velocity are obtained from the GPS receiver, and the attitude angle and the angular velocity are measured by the IMU. The corrected positions from the correctors and the system uncertainty estimations from the extended state observers are used for determination of the desired attitude and design of the controllers. The controllers (65) and (76) drive the UAV to track the reference trajectory. The performance of position correction by the correctors is compared with the EKF-based method [21].

8.1 Design of correctors

8.1.1 Determination of upper-bound of sensor accuracy ratio

From the position sensor accuracy, we get $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_{ * 1}}(t)} \right| \le {L_{ * 1}} = 3\,({\rm{m}})$ , $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_{ * 2}}(t)} \right| \le {L_{ * 2}} = 0.1\,({\rm{m}}/{\rm{s}})$ , where, $ * = \left\{ {x,y,z} \right\}$ . The upper-bound of sensor accuracy ratio $\varepsilon = 0.034$ can be selected such that the sensor accuracy inequality $\frac{{{L_2}}}{{{L_1}}} = \frac{{0.1}}{3} \le {\varepsilon _ * } = 0.034 \ll 1$ holds.

From the attitude sensor accuracy, we get $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_{ * 1}}(t)} \right| \le {L_{ * 1}}{ = 0.5^ \circ }$ , $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_{ * 2}}(t)} \right| \le {L_{ * 2}}{ = 0.0025^ \circ }/{\rm{s}}$ , where, $ * = \left\{ {\phi ,\theta } \right\}$ ; and $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_{\psi 1}}(t)} \right| \le {L_{\psi 1}}{ = 1^ \circ }$ , $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{d_{\psi 2}}(t)} \right| \le {L_{\psi 2}}{ = 0.0025^ \circ }/{\rm{s}}$ . For the roll/pitch, the upper-bound of sensor accuracy ratio ${\varepsilon _ * } = 0.005$ can be selected such that the sensor accuracy inequality $\frac{{{L_2}}}{{{L_1}}} = \frac{{0.0025}}{{0.5}} \le {\varepsilon _ * } = 0.005 \ll 1$ holds, where, $ * = \left\{ {\phi ,\theta } \right\}$ . For the yaw, the upper-bound of sensor accuracy ratio ${\varepsilon _\psi } = 0.0025$ can be selected such that the sensor accuracy inequality $\frac{{{L_2}}}{{{L_1}}} = \frac{{0.0025}}{1} \le {\varepsilon _\psi } = 0.0025 \ll 1$ holds.

8.1.2 Corrector parametres selection

For the position, we select $r = \frac{1}{2}$ to get the minimum value and minimum range of rejection ratio $\rho ({\omega _1})$ : $\min \left\{ {\rho ({\omega _{ * 1}})} \right\} = \frac{1}{2}\varepsilon _ * ^{\frac{1}{2}} = \frac{1}{2}{0.034^{0.5}} = 0.09$ , and $\rho ({\omega _{ * 1}}) \in \left( {\frac{1}{2}\varepsilon _ * ^{\frac{1}{2}},1 + \frac{1}{2}\varepsilon _ * ^{\frac{1}{2}}} \right] = \left( {0.09,1.09} \right]$ when ${\omega _{ * 1}} \in \left[ {0, + \infty } \right)$ , where, $ * = \left\{ {x,y,z} \right\}$ .

From ${\varepsilon _ * } = 0.034$ and $r = \frac{1}{2}$ , we can determine the corrector parametre ${k_{ * 1}} = 2\varepsilon _ * ^{1 - r} = 2 \times {0.034^{1 - 0.5}} = 0.36$ .

According to $\mathop {\sup }\nolimits_{t \in \left[ {0,\infty } \right)} \left| {{{\dot v}_{ * 0}}(t)} \right| \le {L_{ * 3}} = 50$ , ${k_{ * 2}} \gt {k_{ * 1}} + {L_{ * 3}}$ and ${k_{ * 3}} \gt {L_{ * 3}}$ , we select ${k_{ * 2}} = 100$ , ${k_{ * 3}} = 100$ , where, $ * = \left\{ {x,y,z} \right\}$ .

For the attitude, we know that ${\varepsilon _\phi } = 0.005$ , ${\varepsilon _\theta } = 0.005$ and ${\varepsilon _\psi } = 0.0025$ , we select $r = \frac{1}{2}$ . Then, we get the minimum rejection ratio $\min \left\{ {\rho ({\omega _{ * 1}})} \right\} = \frac{1}{2}\varepsilon _ * ^{\frac{1}{2}} = \frac{1}{2}{0.005^{0.5}} = 0.035$ , and $\rho ({\omega _{ * 1}}) \in \left( {0.035,1.035} \right]$ when ${\omega _{ * 1}} \in \left[ {0,\infty } \right)$ , where, $ * = \left\{ {\phi ,\theta } \right\}$ ; and $\min \left\{ {\rho ({\omega _{\psi 1}})} \right\} = \frac{1}{2}\varepsilon _\psi ^{\frac{1}{2}} = \frac{1}{2}{0.0025^{0.5}} = 0.025$ , and $\rho ({\omega _{\psi 1}}) \in \left( {0.025,1.025} \right]$ when ${\omega _{\psi 1}} \in \left[ {0,\infty } \right)$ .

Therefore:

From ${\varepsilon _ * } = 0.005$ and $r = \frac{1}{2}$ , we can determine ${k_{ * 1}} = 2{\varepsilon ^{1 - r}} = 2 \times {0.005^{1 - 0.5}} = 0.14$ , where, $ * = \left\{ {\phi ,\theta } \right\}$ .

From ${\varepsilon _\psi } = 0.0025$ and $r = \frac{1}{2}$ , we can determine ${k_{\psi 1}} = 2{\varepsilon ^{1 - r}} = 2 \times {0.0025^{1 - 0.5}} = 0.1$ .

For the other parametres, to overcome the effect of angular accelerations on the correctors, we select the relatively large ${k_{ * 2}} = 10$ , ${k_{ * 3}} = 10$ , where, $ * = \left\{ {\phi ,\theta ,\psi } \right\}$ .

Therefore, we get the corrector parametres.

Correctors for position: ${k_{ * 1}} = 0.36$ , ${k_{ * 2}} = 100$ , ${k_{ * 3}} = 100$ , where, $ * = \left\{ {x,y,z} \right\}$ ;

Correctors for attitude: ${k_{\phi 1}} = 0.14$ , ${k_{\theta 1}} = 0.14$ , ${k_{\psi 1}} = 0.1$ ; ${k_{ * 2}} = 10$ , ${k_{ * 3}} = 10$ , where, $ * = \left\{ {\phi ,\theta ,\psi } \right\}$ .

8.1.3 Rejection ratio and disturbance frequency bands

1. (1) Rejection ratios in the disturbance frequency bands for position

From ${\varepsilon _ * } = 0.034$ and $r = \frac{1}{2}$ , where, $ * = \left\{ {x,y,z} \right\}$ , the rejection ratio is expressed by

$$\rho ({\omega _{ * 1}}) = \frac{1}{{\sqrt {1 + \frac{1}{4}\varepsilon _ * ^{2r - 2}\omega _{ * 1}^2} }} + \frac{1}{2}\varepsilon _ * ^r = \frac{1}{{\sqrt {1 + 7.35\omega _{ * 1}^2} }} + 0.09$$

Therefore, the rejection ratios in the different frequency bands of $\left[ {0, + \infty } \right)$ can be described by:

1. (i) In $\left[ {3.98{\rm{rad/s}},\infty } \right)$ , $\rho ({\omega _{ * 1}})\;:\;0.18 \to 0.09$ as ${\omega _{ * 1}}\;:\;3.98 \to \infty \,({\rm{rad}}/{\rm{s}})$ ;

2. (ii) In $\left( {0.17,3.98{\rm{rad/s}}} \right)$ , $\rho ({\omega _{ * 1}})\;:\;1 \to 0.18$ as ${\omega _{ * 1}}\;:\;0.17 \to 3.98\,({\rm{rad}}/{\rm{s}})$ ;

3. (iii) In $\left[ {0,0.17{\rm{rad/s}}} \right]$ , $\rho ({\omega _{ * 1}})\;:\;1.09 \to 1$ (i.e. $\rho ({\omega _{ * 1}}) \approx 1$ ) as ${\omega _{ * 1}}\;:\;0 \to 0.17\,({\rm{rad}}/{\rm{s}})$ .

1. (2) Rejection ratios in the disturbance frequency bands for attitude ( $\phi $ , $\theta $ )

From ${\varepsilon _ * } = 0.005$ and $r = \frac{1}{2}$ , where, $ * = \left\{ {\phi ,\theta } \right\}$ , the rejection ratio is expressed by

$$\rho ({\omega _{ * 1}}) = \frac{1}{{\sqrt {1 + \frac{1}{4}\varepsilon _ * ^{2r - 2}\omega _1^2} }} + \frac{1}{2}\varepsilon _ * ^r = \frac{1}{{\sqrt {1 + 50\omega _{ * 1}^2} }} + 0.035$$

Therefore, the rejection ratios in the different frequency bands can be described by:

1. (i) In $\left[ {4{\rm{rad/s}},\infty } \right)$ , $\rho ({\omega _{ * 1}})\;:\;0.07 \to 0.035$ as ${\omega _{ * 1}}\;:\;4 \to \infty \,({\rm{rad}}/{\rm{s}})$ ;

2. (ii) In $\left( {0.04,4{\rm{rad/s}}} \right)$ , $\rho ({\omega _{ * 1}})\;:\;1 \to 0.07$ as ${\omega _{ * 1}}\;:\;0.04 \to 4\,({\rm{rad}}/{\rm{s}})$ ;

3. (iii) In $\left[ {0,0.04{\rm{rad/s}}} \right]\rho ({\omega _{ * 1}})\;:\;1.035 \to 1$ (i.e. $\rho ({\omega _{ * 1}}) \approx 1$ ) as ${\omega _{ * 1}}\;:\;0 \to 0.04\,({\rm{rad}}/{\rm{s}})$ .

1. (3) Rejection ratios in the disturbance frequency bands for attitude ( $\psi $ )

From ${\varepsilon _\psi } = 0.0025$ and $r = \frac{1}{2}$ , the rejection ratio is expressed by

$$\rho ({\omega _{\psi 1}}) = \frac{1}{{\sqrt {1 + \frac{1}{4}\varepsilon _\psi ^{2r - 2}\omega _{\psi 1}^2} }} + \frac{1}{2}\varepsilon _\psi ^r = \frac{1}{{\sqrt {1 + 100\omega _{\psi 1}^2} }} + 0.025$$

Therefore, the rejection ratios in the different frequency bands can be described by:

1. (i) In $\left[ {4{\rm{rad/s}},\infty } \right)$ , $\rho ({\omega _{\psi 1}})\;:\;0.05 \to 0.025$ as ${\omega _{\psi 1}}\;:\;4 \to \infty \,({\rm{rad}}/{\rm{s}})$ ;

2. (ii) In $\left( {0.023,4{\rm{rad/s}}} \right)$ , $\rho ({\omega _{\psi 1}})\;:\;1 \to 0.05$ as ${\omega _{\psi 1}}\;:\;0.023 \to 4\,({\rm{rad}}/{\rm{s}})$ ;

3. (iii) In $\left[ {0,0.023{\rm{rad/s}}} \right]$ , $\rho ({\omega _{\psi 1}})\;:\;1.025 \to 1$ (i.e. $\rho ({\omega _{\psi 1}}) \approx 1$ ) as ${\omega _{\psi 1}}\;:\;0 \to 0.023\,({\rm{rad}}/{\rm{s}})$ .

8.2 Parametres of observers and controllers

According to the selection rules of observer parametres [24], we select the extended state observer parametres: ${\lambda _{1 * }} = 4$ , ${\lambda _{ * 2}} = 20$ , ${\alpha _ * } = 0.6$ , where, $ * = \left\{ {x,y,z,\theta ,\phi ,\psi } \right\}$ . According to the properties and tests of engine and digital servos, we select the control law parametres: ${k_{p1}} = 16$ , ${k_{p2}} = 8$ , ${k_{a1}} = 25$ , ${k_{a2}} = 8$ .

8.3 Analysis of UAV navigation and control performance

Figure 8(a) shows the comparison of the flight trajectories, including the measured from GPS, the reference trajectory, and the estimations by the corrector and the EKF-based method. Meanwhile, the trajectory comparisons in the three directions are shown in Fig. 8(b): due to the effect of adverse conditions, e.g. engine vibration and communication, the actual measurement disturbances in position from GPS were about 10m. The estimate errors by the corrector were less than 1m, while the estimate errors by the KF were about 3m. From the estimate errors and the above numerical calculation, we can find that the position disturbances were mainly within the frequency band $\left[ {3.98{\rm{rad/s}}, + \infty } \right)$ . Therefore, the disturbances in position measurements were rejected sufficiently by the correctors, and the correctors provided the relatively accurate and smoothed correction outputs. The attitude angle comparisons in the three directions are shown in Fig. 9. During the UAV flight, the actual measurement disturbances in attitude from the IMU were about ${3^ \circ } \sim {4^ \circ }$ . The corrector estimate errors for attitude angle were less than ${0.2^ \circ }$ . Then, comparing the above calculation, we can find that the attitude disturbances were mainly within the frequency band $\left[ {4{\rm{rad/s}}, + \infty } \right)$ . Therefore, the disturbances in the attitude measurements were rejected sufficiently by the correctors due to the very small rejection ratio in this frequency band. From the flight test, we can find that the correctors also reduced the sensing disturbances from the effect of UAV vibrations, and the jet UAV remained in the safe flight condition throughout the flight.

Figure 9. Attitude correction.

9.0 Conclusions

In this paper, a sliding mode corrector has been presented, which can correct disturbance in position measurement using relatively accurate velocity. The performance of the corrector was demonstrated by two simulation examples and a jet UAV flight test: (i) it succeeded in rejecting the disturbances largely in position and attitude sensing, even though the disturbances are in the low/mid/high frequency bands; (ii) the experimental test verified the validity of the corrector’s providing accurate and smoothed estimate of position and attitude; (iii) the estimate outputs from the correctors can be used directly by the control system without any additional filters. The merits of the corrector include its model free, bounded corrector gains, the accurate and smoothed estimate outputs and strong parametre inclusion to change of disturbance and signal.