# Wind and Payload Disturbance Rejection Control Based on Adaptive Neural Estimators: Application on Quadrotors

Authors: Jesús Enrique Sierra; Matilde Santos
Journal: Complexity
 (2019)

---

## Abstract

In this work, a new intelligent control strategy based on neural networks is proposed to cope with some external disturbances that can affect quadrotor unmanned aerial vehicles (UAV) dynamics. Specifically, the variation of the system mass during logistic tasks and the influence of the wind are considered. An adaptive neuromass estimator and an adaptive neural disturbance estimator complement the action of a set of PID controllers, stabilizing the UAV and improving the system performance. The control strategy has been extensively tested with different trajectories: linear, helical, circular, and even a lemniscate one. During the experiments, the mass of the UAV is triplicated and winds of 6 and 9 in Beaufort’s scale are introduced. Simulation results show how the online learning of the estimator increases the robustness of the controller, reducing the effects of the changes in the mass and of the wind on the quadrotor.


---

## Body

## 1. Introduction

In recent years, new and valuable applications of unmanned aerial vehicles (UAV) have emerged in different sectors such as defense, security, construction, agriculture, entertainment, and shipping [1–3]. These and other applications demand the design of efficient and robust controllers for those autonomous vehicles. That is why the modelling and control of these complex and unstable systems still motivate the research and interest of the scientific community [4–11].Nevertheless, the modelling and control of quadrotor vehicles are not an easy task. The complexity comes from the randomness of the airstreams and of the exogenous forces, the high nonlinearity dynamics, the coupling between the internal variables, the uncertainty of the measurements, etc. These factors make the techniques based on artificial intelligence a promising approach for the identification and control of these systems [12].Moreover, these intelligent strategies are especially interesting when the model parameters change while the system is working [13]. For example, the total mass will undergo variations when the vehicle is performing logistic tasks, since the mass depends on the loads that are shipped.There are few works that study the effect of the payload variation on the quadrotor dynamics and that take it into account. In [14], an adaptive control is used to mitigate the impact of the parameter variation by estimating them under guaranteed performance. In [15], least square and gradient methods are implemented for adaptive parameter estimation, which are used to update the control output to the current UAV mass and inertia moment. Recently, Wang [16] applied again the same strategy to estimate the variations in the payload and the effect of the wind gusts.We propose to use estimators based on neural networks. They offer three main advantages with respect to other techniques: the fact that knowledge about the internal structure of the system to be estimated is not necessary, unlike other estimation techniques such as PEM (parametric error models), Hammerstein-Weiner, and Volterra; its online learning ability; and the easy parallelization. Even more, using adaptive neural networks allows the online learning of the estimation not to be so biased by the selection of the training data. This is a relevant advantage against other offline estimation methods which are very sensitive to the dataset used during the training.There are other papers more focused on a closely related problem, the UAV stabilization and path tracking when its center of gravity changes, and what happens by the manipulation of suspended loads by cables. In [17], using the flatness property, a trajectory generation method is presented that enables finding nominal trajectories with various constraints that result in minimal load swing. Also, the same authors present a very interesting cooperative control strategy to manipulate suspended loads by several quadrotors at the same time [18]. In [19], a controller based on a set of connected PD regulators is used to tackle this issue. An adaptive tracking controller based on output feedback linearization is used in [20]; this controller compensates for dynamical changes in the center of gravity of the quadrotor. In [21], this issue is again addressed using the dynamic programming approach. In [22], the complexity of the aerial vehicle is incremented by considering the elasticity of the cable in the system equations.Other studies are focused on the rejection of wind disturbances. In [23], Lyapunov-based observers are used to estimate the external force disturbances. In [24], a control strategy based on sliding mode and adaptive control techniques is proposed to deal with slow and fast time-varying wind conditions. In [25], a switching model predictive attitude controller for an unmanned aerial vehicle subject to atmospheric disturbances is presented. In [26], a nonlinear adaptive state feedback controller for thrust and torque actuation is designed, so that it guarantees global convergence of the closed-loop path following in the presence of constant wind disturbances. In [27], a sliding mode control driven by sliding mode disturbance observer (SMC-SMDO) approach is used to design a robust flight controller for a small quadrotor.So far, the studies are mainly focused on payload variations and on wind disturbances rejection. Only few recent papers have been found that address both problems at the same time, such as in [16]. This challenging issue demands further research. Furthermore, although there are some papers where neural networks are applied to model quadrotors [28, 29], and to control them [30–33], these techniques have not been explored to solve this specific research problem.Therefore, in this work, we propose the design of an intelligent control strategy based on neural networks to cope with these external disturbances, payload changes and wind, that can affect quadrotor dynamics. The final goal is to stabilize the UAV and to improve the system performance. The control strategy has been extensively tested by simulation with different trajectories. Indeed, the online learning estimator that has been implemented increases the robustness of the controller, reducing the effects produced for these variations.The paper is organized as follows. In Section2, the dynamic behavior of the system is described. The design of the controller and the adaptive neural estimators are presented in Section 3. Simulation results are discussed in Section 4. The document ends with the conclusions and future works.


## 2. System Model

A quadrotor vehicle is composed by four perpendicular arms, each one with a motor and a propeller (Figure1, left). The four motors drive the lift and direction control.Figure 1
A quadrotor vehicle (left) and UAV’s coordinate system (right).The UAV absolute position is described by three coordinates,(x,y,z), and the attitude is given by the three Euler’s angles (ϕ,θ,ψ), under conditions: (-π≤ψ<π) for the yaw angle, (-π/2≤ϕ<π/2) for the roll, and (-π/2≤θ<π/2) for the pitch, all angles in radians.The system is based on two couples of propellers opposed to each other, (1, 3) and (2, 4) (Figure1, right). To keep the balance of the system, one pair of motors turns clockwise while the other one spins counterclockwise. The increment of the speed of rotor 3 with respect to rotor 1 produces a positive pitch (θ>0), while increasing the speed of rotor 4 regarding rotor 2 produces a positive roll (ϕ>0). The increment of the speeds of rotors 1 + 3 with respect to rotors 2 + 4 produces a positive yaw (ψ>0).Regarding the UAV modelling, in the related literature, there are typically two approaches to obtain the mathematical model of the quadrotor: the Lagrangian method and the one based on the representation of the translational and angular dynamics. As we need to develop a control-oriented model, we have used the latter, that is, the Newton-Euler method, which describes the dynamic systems in terms of force and momentum.The Newton dynamic equation states that the sum of forces applied in a system is equal to the variation of the lineal momentum:(1)∑Fi=ddtmv.While the mass is constant, this equation is equivalent to(2)∑Fi=mv˙.In our case, the forces Fi are the vector of forces T produced by the rotors and by the gravity, considering an Earth reference system; thus, the previous equations give the translational dynamic (3) considering the assumption mv˙≫mv˙, implicitly assumed in other papers [33]:(3)mv˙=RT-mge3,where m is the mass of the quadrotor in Kg, R is the rotation matrix which is dimensionless, g is the gravitational acceleration in m·s−2, T is a vector of forces in N, and e3=[0,0,1]T is a unit vector which describes the rotor orientation.The Euler dynamic equation expresses that the sum of torques is equivalent to the variation of the angular moment:(4)∑μi=ddtJ×w=ddtJ×w+J×ddtw.This sum of torques represents the vector of torques τ (N·m) in the three axes produced by the rotors, that is, the angular dynamic of system (5) [34], considering the assumption J×d/dtw≫d/dtJ×w, also implicitly assumed in other papers [33]:(5)τ=Jω˙+ω×Jω,where J is the inertia tensor in Kg·m2 (6), ω is the angular velocities vector in rad/s, and × represents the vector product: (6)J=Ix000Iy000Iz.The vectors τ (7) and T (8) are a function of the velocities of the propellers:(7)τ=blΩ42-Ω22blΩ32-Ω12dΩ22+Ω42-Ω12-Ω32(8)T=00bΩ12+Ω22+Ω32+Ω42,where b is the thrust coefficient in N·s2, d is the drag coefficient, l is the longitude of each arm in m, and Ω1,…,Ω4 are the velocities in rad/s of the rotors 1 to 4, respectively.To simplify the calculations, instead of using the speed of the rotors, it is possible to define a set of control signalsu1,u2,u3 y u4 as follows (9):(9)u1u2u3u4=11110-101-10101-11-1Ω12Ω22Ω32Ω42.This matrix is invertible, so it is possible to generate speed references for the rotors from a set of control signals.Finally, from (1) to (9), the following system of equations is derived:(10)ϕ¨=θ˙ψ˙Iy-IzIx+lbIxu2(11)θ¨=ϕ˙ψ˙Iz-IxIy+lbIyu3(12)ψ¨=ϕ˙θ˙Ix-IyIz+dIzu4(13)X¨=-sin⁡θcos⁡ϕbmu1(14)Y¨=sin⁡ϕbmu1(15)Z¨=-g+cos⁡θcos⁡ϕbmu1.The constants of (10) to (15) that are used during the simulations are listed in Table 1. The values have been extracted from [35].Table 1
Parameter values of the model.
Parameter Description Value/Units l Longitude of an arm 0.232 m m Mass of the quadrotor 0.52 Kg d Drag coefficient 7.5e−7 N.m.s2 b Thrust coefficient 3.13e−5 N.s2 I x Inertia in X 6.228e−3 Kg.m2 I y Inertia in Y 6.225e−3 Kg.m2 I z Inertia in Z 1.121e−2 Kg.m2 ρ a i r Density of the air 1.2 Kg/m3 A Area in the direction of the wind 0.0186 m2 C d Wind drag coefficient 1


## 3. Control Strategy Design

### 3.1. First Approach: Adaptive Inverse Controller

There are different control strategies based on neural networks [29–34]. In a first approach [32], a variant of the generalized learning algorithm (GLA) was used to control this system (Figure 2). The procedure was as follows.Figure 2
Offline training to identify the plant inverse dynamic.The first step is the application of the GLA algorithm to offline train the neural network in order to identify the inverse dynamic of the plant (Figure2). Once the network has been offline trained, it is placed in cascade connection with the plant and a PID controller. In Figure 3, this control strategy is shown for the altitude variable of the UAV. Then, the configuration of the network is online refined. In order to do this, during each control interval, two processes are sequentially applied to the network (first the simulation, later the online learning):(1)
Simulation: The output of the PID, PIDOUTti, feeds one of the inputs of the artificial neural network; the rest of the inputs are past values of the plant output, PLANTOUTti-j∗Ts. The network generates the control input, u1 (16), which is the input of the plant (Figure 3, switch in the upper position): (16)u1ti=fNETPIDOUTti,PLANTOUTti-j∗Ts,parNETti-Tsj=1…Ninputs-1,whereTs is the sampling time in seconds, parNET denotes the configuration parameters of the network, and Ninputs means the number of inputs of the neural network.(2)
Online learning: The neural network is trained again with the current and previous outputs of the plant, in order to generate the control output, u1, obtaining the new configuration parameters parNET (17). The network input dataset is made up of the past values of the plant output, PLANTOUTti-j∗Ts. The output dataset is the current value of the plant input u1(ti) (Figure 3, switch in the lower position):(17)parNETti=fPLANTOUTti-j∗Ts,u1ti,parNETti-Ts,j=0…Ninputs-1.In order to test the validity of this first approach, we firstly focused only on the altitude control. UAVs are normally provided with accelerometers, so it is assumed that the acceleration in the z-axis (Z¨) is available. The network must be able to simulate the control signal u1 by using acceleration measurements.Figure 3
Adaptive altitude neurocontrol strategy.In this example,PIDOUT (16) is the reference of the acceleration in the z-axis, Z¨PID; PLANTOUT is the acceleration in the z-axis in m·s−2, and u1 is the control signal. Thanks to the artificial neural network, the PID does not need to include the plant gain. The network is able to learn the plant gain and work with it. In other words, with this approach, it is not necessary to know the system parameters to control it [32].


### 3.2. New Advanced Strategy: Controller with Adaptive Neural Estimators

The generalized learning algorithm proposed in the previous section is especially useful when there is not any knowledge about the real dynamic of the plant (black-box system). However, when some knowledge about the system dynamics is available, even if it is not complete or accurate, it is positive to include it in some way in the controller. In this section, the system equations ((10) to (15)) are introduced in the controller, and the neural networks are focused on the uncertain terms: the mass, the wind disturbance, and the nonmodelled dynamics.The control system has been designed to track trajectories defined by tuples of three coordinates(Xref,Yref,Zref).The model of the whole system with the controllers is shown in Figure4. For the sake of simplicity, several ports that route the signals haven been introduced to reduce the number of lines in the diagram.Figure 4
Control system with mass and disturbance estimators.The UAV model control inputs are the four control signals,ui, i=1,..,4, that represent the power of the rotors. Four main controllers are defined to obtain these model inputs, the controllers of the X, Y, and Z coordinates, and another for the yaw angle. This is because the pitch and roll angles are used to track the Xref and Yref reference coordinates. The control of Z is carried out by the control signal u1. The Y controller generates the roll reference, and the tracking of the roll angle is performed by the control signal u2. In the same way, the controller of X generates the pitch angle reference value, and control signal u3 will be in charge of getting this value. The control signal u4 is used to stabilize the yaw angle around zero.The aim of the PIDs controllers that appear in Figure4 is to generate the acceleration references to make the attitude error (ϕr-ϕ,θr-θ,ψr-ψ) and the tracking error converge to zero. The rest of the controllers are used to compensate the nonlinearities of the system. The PID tuning parameters have been set by trial and error and the values of the gains are shown in Table 2.Table 2
PID parameters for each variable.
Signal Kp Kd Ki X 8 8 0 Y 8 8 0 Z 2 2 0.9 Roll 8 8 0 Pitch 8 8 0 Yaw 8 8 0A mass estimator and a wind disturbance estimator have been added (Figure4). They are implemented by neural networks. These estimators are used to compensate the variations of the system’s mass and the influence of the wind disturbances. These adaptive neural estimators feed the inputs of the X, Y, and Z controllers.The controllers are defined as follows.The performance of theZ controller is given by(18)u1ti=rZti+gmestti-1bcos⁡θi-1cos⁡ϕi-1=0rZti+gmestti-1b1cos⁡θi-1cos⁡ϕi-1cos⁡θi-1cos⁡ϕi-1≠0,where rZ is the output of the Z PID controller, mest is the estimation of the mass in Kg, (ϕi-1,θi-1) represents the roll and pitch signals at time ti-1, and the rest of the parameters and variables have been previously defined. In (18), it is possible to see how there is a discontinuity at cos⁡θcos⁡ϕ=0 that has been taken into account.TheY controller is defined by(19)ϕRti=ϕRti-1u1ti=0asin⁡rYti+distYestti-1mestti-1bu1tiu1ti≠0,where ϕR is the reference of the roll PID controller in rad, rY is the output of the Y PID controller in m·s−2, distYest is the estimation of the wind disturbance in the y-axis in m·s−2, and the rest of the parameters and variables have been already defined. The discontinuity at u1(ti)=0 has been considered.The control of theX coordinate is given by the expression:(20)θRti=θRti-1u1ticos⁡ϕi-1=0-asin⁡rXti+distXestti-1mestti-1bu1ticos⁡ϕi-1u1ticos⁡ϕi-1≠0,where θR is the reference of the pitch PID controller, rX is the output of the X PID controller in m·s−2, distXest is the estimation of the wind disturbance in the X coordinate in m·s−2, and the rest of the parameters and variables have been already cited. Again, there is a discontinuity at u1ticos⁡ϕi-1=0.


### 3.3. Adaptive Neural Estimator for Disturbances

The control scheme proposed in Figure4 of this paper uses three different estimations: the estimation of the total mass of the system, the estimation of the disturbances in the X coordinate, and the estimation of the disturbances in the Y coordinate. These approximations, as may be observed in (18)-(20), feed the inputs of the different nonlinear controllers in order to reject the effect of the changes in the mass (mest) and the wind external disturbances (distXest and distYest).These parameters that affect the dynamic of the UAV must be estimated if there are no sensors that could measure these disturbances, as it is the case. Although the estimation of each one of these variables (mest, distXest, and distYest) has been implemented in a different function block, the inner structure of all of them is the same. They differ in some parameters of the configuration and the input and output signals. Thus, they can be jointly explained.Theestimator is based on an artificial neural network with online learning. Thus, there is one neural network to, let say, model each of the three parameters or disturbances considered. The following figure represents its generic structure (Figure 5). It is based on the diagram shown in [28], but in this case a new parameter model has been included. Furthermore, in this case, the inputs of the neural network are the outputs of the UAV, and its output are the output of the parameter model (when in [28], the output of the neural network is fed by the outputs of the UAV).Figure 5
Neural estimator structure and configuration.The parameter model receives the inputs and outputs of the UAV. The output of the parameter model is used as target output of the neural network during the training. The output of the proposed estimator is always the output of the neural network. The parameter model is needed because we are using supervised artificial networks and the inputs, and their corresponding outputs must be known.We will use one parameter model for the mass, one parameter model for the wind disturbance in the x-axis, and another for the y-axis. The complete process is further explained below.Each element of the dataset used to offline train the network is composed of the following:(i)
Target: Parameter (disturbance) value calculated with the model at instantti, param(ti)(ii)
Network inputs: for each output signali, the previous Ni values to ti are collected and structured as an array. A theoretical example would be the following: At ti=10 and for M = 3 outputs with a configuration Nout1=3, Nout2=2, Nout3=1, the network inputs are [out1(7), out1(8), out1(9), out2(8), out2(9), out3(9)] and the target is: param(10).If the parameter cannot be calculated (division by zero, squared root of negative numbers, or any other singularity), that element is not included into the training dataset.In this offline learning, the training dataset has as many elements as previous instants of time are considered. Longer time will normally produce better accuracy but will require more computational effort. The selection of the data for training is a delicate task. Indeed, the accuracy of the model depends on the data used to train the network. Another disadvantage of exclusively using offline approaches is that they do not capture the dynamics when it is changing over time [28]. For these reasons, in our proposal, we use adaptive learning for the mass and wind disturbance estimators’ calculation.Once the offline learning has finished, at each instant of time, a new training element is added{target,networkinputs}, as has been previously explained. The target is obtained as the output of the model and the network inputs the output signals of the UAV. This new element is used to teach the networks how to adapt its parameters according to the new input. That is, the function to estimate the parameter is continuously changing over time.Pseudocode1 details the algorithm which relates to offline learning, the simulation, the online learning, and how the parameters are updated.Pseudocode 1
Net= configureNet () Dataset= ∅ For i=1 to toff    #Generating the dataset for off-line learning inputNet =[out1(i-1),… out1(i-n1), out2(i-1),…, outM (i-1)…. outM (i- nM )] ParM (i)=model(out1(i-1),….., outM (i-1)]) Element =  parM (i), inputNet Dataset=Dataset ∪ Element endFor Net= offlineTraining (Net, Dataset) For i= toff to tend   #on-line learning inputNet =[out1(i-1),… out1(i-n1), out2(i-1),…, outM (i-1)…. outM (i- nM )] ParS (i)=simulate(Net, inputNet ) [out1(i),….., outM (i)]= executeUAV ( ParS (i))     #Execute UAV + Controller ParM (i)=model(out1(i-1),….., outM (i-1)]) If ∃ ParM(i) then Element=  parM (i), inputNet Net= onlineTraining (Net, Element) endIf endForAs it has been commented before, when it is not possible to measure or to calculate the parameter with the model, the input values are not defined, or there is a singularity in the calculation, the artificial neural network of our approach is very useful to estimate the disturbance. In these situations, we could say roughly speaking that the network generates new knowledge.Once the inner structure of the estimator has been introduced, we explain the specific configurations for the mass and the wind disturbances.Mass Estimation (i)
Target output: The mass is approximated with the model:(21)m~ti=Z¨ti-1+gu1ti-1cos⁡θi-1cos⁡ϕi-1. (ii)
Network inputs:
(a)
Acceleration in the z-axisZ¨ti-1 (b)
Roll and pitch cosine angles multiplied (cos⁡θi-1cos⁡ϕi-1).Disturbance in X-Axis Estimation (i)
Target output: The disturbance is approximated with the model:(22)distX~ti=-sin⁡θi-1cos⁡ϕi-1bu1ti-1mestti-1-X¨ti-1. (ii)
Network inputs:
(a)
Acceleration in the x-axisX¨ti-1 (b)
Roll angleθi-1 (c)
Pitch angleϕi-1.Disturbance in Y-Axis Estimation (i)
Target output: The disturbance is approximated with the model:(23)distY~ti=sin⁡ϕi-1bu1ti-1mestti-1-Y¨ti-1. (ii)
Network inputs:
(a)
Acceleration in the y-axisY¨ti-1 (b)
Roll angleθi-1 (c)
Pitch angleϕi-1.For the three different estimators, the artificial network implemented is a multilayer perceptron (MLP) with a hidden layer. The number of neurons of the hidden layer has been set to 20. The Levenberg-Marquardt algorithm withμ=0.001 has been used for the training. The network is offline trained for the first 2 seconds and then the online learning is applied for the remaining 13 seconds.


## 3.1. First Approach: Adaptive Inverse Controller

There are different control strategies based on neural networks [29–34]. In a first approach [32], a variant of the generalized learning algorithm (GLA) was used to control this system (Figure 2). The procedure was as follows.Figure 2
Offline training to identify the plant inverse dynamic.The first step is the application of the GLA algorithm to offline train the neural network in order to identify the inverse dynamic of the plant (Figure2). Once the network has been offline trained, it is placed in cascade connection with the plant and a PID controller. In Figure 3, this control strategy is shown for the altitude variable of the UAV. Then, the configuration of the network is online refined. In order to do this, during each control interval, two processes are sequentially applied to the network (first the simulation, later the online learning):(1)
Simulation: The output of the PID, PIDOUTti, feeds one of the inputs of the artificial neural network; the rest of the inputs are past values of the plant output, PLANTOUTti-j∗Ts. The network generates the control input, u1 (16), which is the input of the plant (Figure 3, switch in the upper position): (16)u1ti=fNETPIDOUTti,PLANTOUTti-j∗Ts,parNETti-Tsj=1…Ninputs-1,whereTs is the sampling time in seconds, parNET denotes the configuration parameters of the network, and Ninputs means the number of inputs of the neural network.(2)
Online learning: The neural network is trained again with the current and previous outputs of the plant, in order to generate the control output, u1, obtaining the new configuration parameters parNET (17). The network input dataset is made up of the past values of the plant output, PLANTOUTti-j∗Ts. The output dataset is the current value of the plant input u1(ti) (Figure 3, switch in the lower position):(17)parNETti=fPLANTOUTti-j∗Ts,u1ti,parNETti-Ts,j=0…Ninputs-1.In order to test the validity of this first approach, we firstly focused only on the altitude control. UAVs are normally provided with accelerometers, so it is assumed that the acceleration in the z-axis (Z¨) is available. The network must be able to simulate the control signal u1 by using acceleration measurements.Figure 3
Adaptive altitude neurocontrol strategy.In this example,PIDOUT (16) is the reference of the acceleration in the z-axis, Z¨PID; PLANTOUT is the acceleration in the z-axis in m·s−2, and u1 is the control signal. Thanks to the artificial neural network, the PID does not need to include the plant gain. The network is able to learn the plant gain and work with it. In other words, with this approach, it is not necessary to know the system parameters to control it [32].


## 3.2. New Advanced Strategy: Controller with Adaptive Neural Estimators

The generalized learning algorithm proposed in the previous section is especially useful when there is not any knowledge about the real dynamic of the plant (black-box system). However, when some knowledge about the system dynamics is available, even if it is not complete or accurate, it is positive to include it in some way in the controller. In this section, the system equations ((10) to (15)) are introduced in the controller, and the neural networks are focused on the uncertain terms: the mass, the wind disturbance, and the nonmodelled dynamics.The control system has been designed to track trajectories defined by tuples of three coordinates(Xref,Yref,Zref).The model of the whole system with the controllers is shown in Figure4. For the sake of simplicity, several ports that route the signals haven been introduced to reduce the number of lines in the diagram.Figure 4
Control system with mass and disturbance estimators.The UAV model control inputs are the four control signals,ui, i=1,..,4, that represent the power of the rotors. Four main controllers are defined to obtain these model inputs, the controllers of the X, Y, and Z coordinates, and another for the yaw angle. This is because the pitch and roll angles are used to track the Xref and Yref reference coordinates. The control of Z is carried out by the control signal u1. The Y controller generates the roll reference, and the tracking of the roll angle is performed by the control signal u2. In the same way, the controller of X generates the pitch angle reference value, and control signal u3 will be in charge of getting this value. The control signal u4 is used to stabilize the yaw angle around zero.The aim of the PIDs controllers that appear in Figure4 is to generate the acceleration references to make the attitude error (ϕr-ϕ,θr-θ,ψr-ψ) and the tracking error converge to zero. The rest of the controllers are used to compensate the nonlinearities of the system. The PID tuning parameters have been set by trial and error and the values of the gains are shown in Table 2.Table 2
PID parameters for each variable.
Signal Kp Kd Ki X 8 8 0 Y 8 8 0 Z 2 2 0.9 Roll 8 8 0 Pitch 8 8 0 Yaw 8 8 0A mass estimator and a wind disturbance estimator have been added (Figure4). They are implemented by neural networks. These estimators are used to compensate the variations of the system’s mass and the influence of the wind disturbances. These adaptive neural estimators feed the inputs of the X, Y, and Z controllers.The controllers are defined as follows.The performance of theZ controller is given by(18)u1ti=rZti+gmestti-1bcos⁡θi-1cos⁡ϕi-1=0rZti+gmestti-1b1cos⁡θi-1cos⁡ϕi-1cos⁡θi-1cos⁡ϕi-1≠0,where rZ is the output of the Z PID controller, mest is the estimation of the mass in Kg, (ϕi-1,θi-1) represents the roll and pitch signals at time ti-1, and the rest of the parameters and variables have been previously defined. In (18), it is possible to see how there is a discontinuity at cos⁡θcos⁡ϕ=0 that has been taken into account.TheY controller is defined by(19)ϕRti=ϕRti-1u1ti=0asin⁡rYti+distYestti-1mestti-1bu1tiu1ti≠0,where ϕR is the reference of the roll PID controller in rad, rY is the output of the Y PID controller in m·s−2, distYest is the estimation of the wind disturbance in the y-axis in m·s−2, and the rest of the parameters and variables have been already defined. The discontinuity at u1(ti)=0 has been considered.The control of theX coordinate is given by the expression:(20)θRti=θRti-1u1ticos⁡ϕi-1=0-asin⁡rXti+distXestti-1mestti-1bu1ticos⁡ϕi-1u1ticos⁡ϕi-1≠0,where θR is the reference of the pitch PID controller, rX is the output of the X PID controller in m·s−2, distXest is the estimation of the wind disturbance in the X coordinate in m·s−2, and the rest of the parameters and variables have been already cited. Again, there is a discontinuity at u1ticos⁡ϕi-1=0.


## 3.3. Adaptive Neural Estimator for Disturbances

The control scheme proposed in Figure4 of this paper uses three different estimations: the estimation of the total mass of the system, the estimation of the disturbances in the X coordinate, and the estimation of the disturbances in the Y coordinate. These approximations, as may be observed in (18)-(20), feed the inputs of the different nonlinear controllers in order to reject the effect of the changes in the mass (mest) and the wind external disturbances (distXest and distYest).These parameters that affect the dynamic of the UAV must be estimated if there are no sensors that could measure these disturbances, as it is the case. Although the estimation of each one of these variables (mest, distXest, and distYest) has been implemented in a different function block, the inner structure of all of them is the same. They differ in some parameters of the configuration and the input and output signals. Thus, they can be jointly explained.Theestimator is based on an artificial neural network with online learning. Thus, there is one neural network to, let say, model each of the three parameters or disturbances considered. The following figure represents its generic structure (Figure 5). It is based on the diagram shown in [28], but in this case a new parameter model has been included. Furthermore, in this case, the inputs of the neural network are the outputs of the UAV, and its output are the output of the parameter model (when in [28], the output of the neural network is fed by the outputs of the UAV).Figure 5
Neural estimator structure and configuration.The parameter model receives the inputs and outputs of the UAV. The output of the parameter model is used as target output of the neural network during the training. The output of the proposed estimator is always the output of the neural network. The parameter model is needed because we are using supervised artificial networks and the inputs, and their corresponding outputs must be known.We will use one parameter model for the mass, one parameter model for the wind disturbance in the x-axis, and another for the y-axis. The complete process is further explained below.Each element of the dataset used to offline train the network is composed of the following:(i)
Target: Parameter (disturbance) value calculated with the model at instantti, param(ti)(ii)
Network inputs: for each output signali, the previous Ni values to ti are collected and structured as an array. A theoretical example would be the following: At ti=10 and for M = 3 outputs with a configuration Nout1=3, Nout2=2, Nout3=1, the network inputs are [out1(7), out1(8), out1(9), out2(8), out2(9), out3(9)] and the target is: param(10).If the parameter cannot be calculated (division by zero, squared root of negative numbers, or any other singularity), that element is not included into the training dataset.In this offline learning, the training dataset has as many elements as previous instants of time are considered. Longer time will normally produce better accuracy but will require more computational effort. The selection of the data for training is a delicate task. Indeed, the accuracy of the model depends on the data used to train the network. Another disadvantage of exclusively using offline approaches is that they do not capture the dynamics when it is changing over time [28]. For these reasons, in our proposal, we use adaptive learning for the mass and wind disturbance estimators’ calculation.Once the offline learning has finished, at each instant of time, a new training element is added{target,networkinputs}, as has been previously explained. The target is obtained as the output of the model and the network inputs the output signals of the UAV. This new element is used to teach the networks how to adapt its parameters according to the new input. That is, the function to estimate the parameter is continuously changing over time.Pseudocode1 details the algorithm which relates to offline learning, the simulation, the online learning, and how the parameters are updated.Pseudocode 1
Net= configureNet () Dataset= ∅ For i=1 to toff    #Generating the dataset for off-line learning inputNet =[out1(i-1),… out1(i-n1), out2(i-1),…, outM (i-1)…. outM (i- nM )] ParM (i)=model(out1(i-1),….., outM (i-1)]) Element =  parM (i), inputNet Dataset=Dataset ∪ Element endFor Net= offlineTraining (Net, Dataset) For i= toff to tend   #on-line learning inputNet =[out1(i-1),… out1(i-n1), out2(i-1),…, outM (i-1)…. outM (i- nM )] ParS (i)=simulate(Net, inputNet ) [out1(i),….., outM (i)]= executeUAV ( ParS (i))     #Execute UAV + Controller ParM (i)=model(out1(i-1),….., outM (i-1)]) If ∃ ParM(i) then Element=  parM (i), inputNet Net= onlineTraining (Net, Element) endIf endForAs it has been commented before, when it is not possible to measure or to calculate the parameter with the model, the input values are not defined, or there is a singularity in the calculation, the artificial neural network of our approach is very useful to estimate the disturbance. In these situations, we could say roughly speaking that the network generates new knowledge.Once the inner structure of the estimator has been introduced, we explain the specific configurations for the mass and the wind disturbances.Mass Estimation (i)
Target output: The mass is approximated with the model:(21)m~ti=Z¨ti-1+gu1ti-1cos⁡θi-1cos⁡ϕi-1. (ii)
Network inputs:
(a)
Acceleration in the z-axisZ¨ti-1 (b)
Roll and pitch cosine angles multiplied (cos⁡θi-1cos⁡ϕi-1).Disturbance in X-Axis Estimation (i)
Target output: The disturbance is approximated with the model:(22)distX~ti=-sin⁡θi-1cos⁡ϕi-1bu1ti-1mestti-1-X¨ti-1. (ii)
Network inputs:
(a)
Acceleration in the x-axisX¨ti-1 (b)
Roll angleθi-1 (c)
Pitch angleϕi-1.Disturbance in Y-Axis Estimation (i)
Target output: The disturbance is approximated with the model:(23)distY~ti=sin⁡ϕi-1bu1ti-1mestti-1-Y¨ti-1. (ii)
Network inputs:
(a)
Acceleration in the y-axisY¨ti-1 (b)
Roll angleθi-1 (c)
Pitch angleϕi-1.For the three different estimators, the artificial network implemented is a multilayer perceptron (MLP) with a hidden layer. The number of neurons of the hidden layer has been set to 20. The Levenberg-Marquardt algorithm withμ=0.001 has been used for the training. The network is offline trained for the first 2 seconds and then the online learning is applied for the remaining 13 seconds.


## 4. Results and Discussion

Simulation results have been obtained with Matlab/Simulink software. The duration of each simulation is 15 s. The controller is offline trained during the first 2 s. Then, the online learning algorithms are applied for the remaining 13 s. The sample timeTs is set to 10 ms.In the experiments, in order to simplify the system, the yaw angle is set to 0, but it could be set to any other value by the user.


### 4.1. Trajectory Tracking without Disturbances

The control system has been first tested and validated with several trajectories without considering any disturbances. For each trajectory, the path followed by the UAV with adaptive neural estimators (red line), without them (yellow line), the reference (blue line), and the tracking error, are shown (Figures6–11). The tracking error is calculated by the following equation: (24)tERRti=Xi-Xrefi2+Yi-Yrefi2+Zi-Zrefi2,where Xrefi,Yrefi,Zrefi are the references for X, Y, and Z in ti, respectively.Figure 6
UAV tracking of a linear trajectory (left) and its tracking error (right).Figure 7
UAV tracking of a circular trajectory (left) and its tracking error (right).Figure 8
UAV tracking of a helical trajectory (left) and its tracking error (right).Figure 9
UAV tracking of a cyclic helical trajectory (left) and its tracking error (right).Figure 10
UAV tracking of a lemniscate trajectory (left) and its tracking error (right).Figure 11
UAV tracking of a helical lemniscate trajectory (left) and its tracking error (right).The trajectories used to validate the controller are described by the following equations, with the corresponding parameters:Linear:(25)XYZ=0.2t+10.2t+15∗stept+Z0.Circular:(26)XYZ=cos⁡0.4t+1sin⁡0.4t+15∗stept+Z0.Helical:(27)XYZ=cos⁡0.8t⁡sin⁡0.8t0.4t+Z0.Cyclic helical:(28)XYZ=cos⁡1.6t⁡sin⁡1.6tsin⁡0.4t⁡+Z0.Lemniscate:(29)XYZ=cos⁡0.4t⁡1+sin⁡0.4t⁡2cos⁡0.4tsin⁡0.4t1+sin⁡0.4t2Z0.Helical lemniscate:(30)XYZ=cos⁡0.8t1+sin⁡0.8t2cos⁡0.8tsin⁡0.8t1+sin⁡0.8t20.1t+Z0.In Figures 6–11, it is possible to observe how the tracking error decreases along the time. In t = 0, the error is high because the value of the reference starts in a very different value than the initial position of the system, that is, (X0,Y0,Z0)=(0,0,20). Then, the controller starts to work, and the tracking error is reduced.For example, in the linear trajectory (Figure6), the reference in t0 is (1, 1, 25); thus, the initial tracking error at t=0 is (1+1+25)=5.196.For the same reason, the initial error of the helical and lemniscate trajectories is much lower. The initial value reference is (1, 0, 20), and thus the initial error is 1.In these figures, there is not any significant visual difference between the results with and without the neural network due to the fact that if the PIDs are well tuned and there are no changes either in the mass or in the external disturbances, as it is the case in this section, the neural networks do not provide relevant advantages. Nevertheless, the networks in the control scheme of Figure4 show their full potential, when they must tackle disturbances, as it can be seen in the following sections. Consequently, the red and yellow lines seem to be overlapped.Other appreciable result is that the stationary tracking error is higher when the trajectory has helical component. Predictably, the circular and helical trajectories show a cyclic nature.As a conclusion, the controller performs a good tracking for a wide range of different trajectories. Thus, it can be said that the control strategy has been validated.


### 4.2. Control Robustness with Mass Variations

#### 4.2.1. Mass Disturbance Model

Adding a payload in the quadrotor typically has three effects: the total mass of the system is increased, and the centre of the gravity can be modified and therefore also the inertia. In this work, we assume that the payload is an isotropic symmetric rigid solid attached to the UAV, not suspended, with dimensions much smaller than the dimensions of the quadrotor. The distance between the centre of gravity of the UAV and the centre of gravity of the load depends on the shape and the weight distribution of the manipulator and the shape and weight distribution of the load. In our experiment, we assume that this distance is zero. Under these circumstances, the effect of the inertia and centre of gravity variation can be neglected. Therefore, our paper is only addressing the effect of the mass variation.When the UAV is performing a logistic task, there are two possible stable states regarding the mass: one is before the load is in contact with the quadrotor, and then only the mass of the UAV is considered; the second one is when the payload mass is part of the system, and the sum of both masses is then considered as an only system. Between these two states, there may be several profiles of mass variation depending on the grasping technology and the properties of the surface of the load.We assume that the grasping and the load are nondeformable; thus, one step profile may be applied regarding the mass disturbance, as in other papers [16].The mass variation is simulated by adding a new termdistm (34) to (13) to (15), resulting in (31) to (34). The modelling of the mass disturbance is a step function. The total mass then is triplicated at t = 4s (34) meaning that a payload has been attached to the UAV. After 4 sec, the total mass of the system is m+mL=3m, with m being the mass of the quadrotor and the mass of the load mL=2m. Even if it can be considered a simple approach of dealing with this disturbance, the final effect is well represented:(31)X¨=-sin⁡θcos⁡ϕbm+distmu1(32)Y¨=sin⁡ϕbm+distmu1(33)Z¨=-g+cos⁡θcos⁡ϕbm+distmu1(34)distm=mL·stept-4.This experiment represents a possible situation while performing a logistic task, where the payload is heavier than the quadrotor itself.


#### 4.2.2. Simulation Results with Mass Variation

Now, we test the control proposal with mass disturbances for the helical lemniscate trajectory because it is the most challenging one.Figure12, left, shows the reference trajectory in blue, the trajectory obtained by the controller with the neuromass estimator in red, and without the mass estimation in yellow (only PID). It is possible to observe how the trajectory obtained without the adaptive neural estimators (yellow line) moves away the reference even if it later comes closer to it. At the right side of Figure 12, the tracking error is presented (red line with mass estimation and blue line without mass estimation). In both figures, it is possible to observe how the performance with the mass estimator is much better. The tracking error is mostly the same until the mass changes at t = 4 s, but then the error without mass estimation increases significantly.Figure 12
UAV tracking of a helical lemniscate trajectory with mass variation (left) and its tracking error (right).The trajectory tracking of each coordinate has been represented in a separated figure to study the performance of the controller. Figure13, top left, shows the tracking of the X coordinate, at the top, right, the tracking of the Y coordinate, and at the bottom the tracking of the Z coordinate, with reference (magenta line), the trajectory with mass estimator (red line), and the trajectory without mass estimator (blue line). For every coordinate, it is possible to see how the trajectory obtained by the controller with the mass estimator fits better the reference. Even though, the deviation mainly comes from the Z coordinate due to the fact that the acceleration in the z-axis is more sensible to changes in the mass.Figure 13
UAV tracking trajectory of X, Y, and Z coordinates, respectively, with mass variation.Figure14 shows the mass obtained by the neuromass estimator in this experiment (real mass in blue and the estimated mass in red). The estimate of the mass is very similar to the real one but slightly noisy. Clearly, it can be appreciated that the mass is triplicated following a step variation profile.Figure 14
Mass obtained by the neural estimator.The controller is robust even with other types of variations in the mass like a sine function, for instance, if the mass changes according to(35)mt=mmodelstept-stept-4+stept-42·mmodel+mmodel·sin⁡2π8t.The results of considering this mass variation profile are shown in Figure 15. The tracking error of each component is like the previous ones, but the system experiences a stronger change with the sinusoidal reference profile. In this case, the tracking error of the system without the adaptive neural estimators tends to increase; meanwhile, in the case of the step change, the tracking error decreased over the time.Figure 15
UAV tracking of a helical lemniscate trajectory with sinusoidal mass variation (left) and its tracking error (right).Therefore, it can be seen how the controller with the neural estimator also works well with different mass change profiles.Figure16 shows the mass obtained by the neuromass estimator (red line) in this case and the corresponding real mass (blue line). The estimate of the mass is again very similar to the real one but slightly noisy but now the sinusoidal nature of the profile can be clearly observed.Figure 16
Mass obtained by the neural estimator with a sinusoidal mass change profile.The controller has been extensively tested for different trajectories and the MSE numerical results are summarized in Table3 for the different trajectories and for each coordinate, with the neural estimator (neuro) and without it (PID). The last column shows the absolute error. The best result for each component and trajectory is boldfaced.Table 3
Comparison of the MSE of the tracking error for different trajectories with mass variation.
Trajectory M S E X M S E Y M S E Z M S E T Neuro PID Neuro PID Neuro PID Neuro PID Linear 0.0707 0.0707 0.0708 0.0708 1.6099 21.9021 0.6635 3.8267 Circular 0.2766 0.2775 0.0720 0.0729 1.6101 21.9026 0.7076 3.8657 Helical 0.0934 0.1053 0.0018 0.0148 0.0029 18.2737 0.1899 3.2131 Cyclic helical 0.1067 0.2145 0.0126 0.1212 0.0046 17.5445 0.2715 3.1897 Lemniscate 0.0902 0.0916 0.003 0.0020 0.0010 18.2297 0.1589 3.1981 Helical lemniscate 0.0923 0.1060 0.0018 0.0166 0.0034 18.2335 0.1956 3.2087Table4 compares the values of the maximum error in this experiment for the same cases with before.Table 4
Comparison of the MAX of the tracking error for different trajectories with mass variation.
Trajectory M A X X M A X Y M A X Z M A X T Neuro PID Neuro PID Neuro PID Neuro PID Linear 1.0270 1.0270 1.0270 1.0270 5.0010 9.7835 5.1986 9.7835 Circular 2 2 1.0562 1.0562 5.0010 9.7836 5.4796 9.7837 Helical 1 1 0.1207 0.2074 0.1582 9.2655 1.0033 9.2675 Cyclic helical 1 1 0.2060 0.6163 0.1591 9.1763 1.0006 9.1937 Lemniscate 1 1 0.0614 0.0962 0.0913 9.2541 1 9.2544 Helical lemniscate 1 1 0.1128 0.2703 0.1322 9.2549 1 9.2553In Tables3 and 4, it is possible to observe that the controller with adaptive neural estimators provides less or equal tracking error for every tested trajectory. A general trend in the “PID column” (without mass disturbance estimation) is that the worst tracking error is obtained for the Z component. One of the reasons of this may be the fact that the maximum amplitude of the reference signal is in the z-axis.


### 4.3. Control Robustness with Wind Disturbances

#### 4.3.1. Wind Disturbance Model

The effect of the wind in an UAV in flight can be considered as a drag force and torque depending on the aerodynamics properties of the object and the characteristics of the air flow.The wind speed is a vector field; that is, its value may be different in each coordinate (x, y, z) of the space. For small UAVs like the one we are using in the experiments, we can assume a planar air flow, so the wind is the same in the region of the space where the UAV is flying. Under these circumstances, the wind influence in the torque can be neglected.On the other hand, the following equation represents the drag force [36]:(36)FD-=A-·Cd·ρair·V-22,where Cd is the drag coefficient, which is usually determined experimentally and collects the complex dependency, A is the area exposed to the wind, ρair is the density of the air, and V is the relative velocity of the flying object with respect to the wind. Since the air flow is planar, it can be assumed that the drag force is fully exerted in the centre of the gravity of the vehicle.In our work, this equation is divided by the mass. This way the expression is transformed to an acceleration ((40) to (43)) that can be easily introduced in the equations of the translational dynamics ((37) to (39)), where wind disturbances have been represented as distwX, distwY, and distwZ.As the wind speed is a vectorial magnitude, in our experiment, we suppose the wind components X and Y are twice the Z one. Due to this assumption, there appears a factor of 2 multiplying the wind speed in (41) and (42):(37)X¨=-distwX-sin⁡θcos⁡ϕbmu1(38)Y¨=-distwY+sin⁡ϕbmu1(39)Z¨=-distwZ-g+cos⁡θcos⁡ϕbmu1.In the free atmosphere, the wind is a balance between the Coriolis, centrifugal, and pressure forces acting on the air mass. But in the boundary layer, the mean wind velocity is also a function of the height [37], and it can be expressed as a logarithmic function. This fact is shown in the following:(40)vwZ=vwZ=20·log⁡Z/Clog⁡20/C(41)distwX=sgn⁡vwZ·ρair·Ax·Cd·X˙-2∗vwZ22m(42)distwY=sgn⁡vwZ·ρair·Ay·Cd·Y˙-2∗vwZ22m(43)distwZ=sgn⁡vwZ·ρair·Az·Cd·Z˙-vwZ22m,where vw(Z=20) is the wind speed at 20 m of altitude in m/s, vw(Z) is the wind speed at Z altitude in m/s, C is a constant related to the flight (in this experiment, the value is set to 1.5) which is dimensionless, ρair is the air density in Kg·m3, Ax to Az are the effective area of the quadrotor exposed to each component of the wind in m2, Cd is the drag coefficient with respect to the wind which is dimensionless, X˙, Y˙, and Z˙ are the velocities in the x-axis, y-axis, and z-axis in m/s, and sgn denotes the sign function.The mean wind velocity can be considered constant during the experiments, but not its instant value. The most common assumption is to consider the wind turbulence as a stationary Gaussian random process [37]. Considering this fact, in our experiment, the wind speed is simulated by a step with Gaussian noise at t = 4 s. The SNR between the average wind and the noise is 10 dB. The average wind speed is 12 m/s in the z-axis and 24 m/s in the x-axis and y-axis. These values match numbers 6 and 9 in Beaufort’s scale (strong breeze and strong gale) [38].


#### 4.3.2. Simulation Results with Wind Variation

Figure17, left, shows the tracking of the trajectory, the reference in blue, the trajectory with adaptive neural estimators in red, and without estimators in yellow. At the right side of Figure 17, the tracking error is shown with estimators (red line) and without estimators (blue line). In both figures, it is possible to see how the performance of the control strategy with the neural estimators is much better. Like in the case of the mass variation, the tracking error is the same until t = 4 s because before there is no wind disturbance, but from that moment on, the tracking error of the controller without the estimators increases significantly.Figure 17
UAV tracking of a helical lemniscate trajectory with wind disturbance (left) and its tracking error (right).If the tracking according to each coordinate, X, Y, and Z, is represented (Figure18), the contribution to the tracking error seems to be more balanced (reference, magenta; trajectory with adaptive neural estimators, red; trajectory without neural estimators, blue). For every coordinate, it is possible to see how the trajectory obtained by the controller with the neural estimators better fits the reference. Figure 18 also shows how the disturbance produces an important deviation around t=7s due to the big peak of disturbance in the y-axis which cannot be compensated without the adaptive neural estimator and it becomes a stationary error. This deviation can be also easily observed in Figure 17 in the PID line.Figure 18
UAV tracking of X, Y, and Z coordinates, respectively, with wind disturbances.Another interesting result which can be drawn from the previous figures is that we are not using a wind disturbance estimator for the Z coordinate (see Figure4) but the tracking is still good enough. The reason is that adaptive neuromass estimator senses the wind disturbances in the z-axis as a virtual mass variation and it can compensate it. This is shown in Figure 19, where the estimate of the wind disturbance in the x-, y-, and z-axis is interpreted by the mass estimator as a mass variation. Indeed, the neural estimator (red line) fits reasonably well the real disturbances (blue line).Figure 19
Wind disturbance estimation in the x- and y-axis, top, and mass estimation (bottom).The controller has been extensively tested for different trajectories and the MSE numerical results are summarized in Table5 for the different trajectories and for each coordinate, with the neural estimator (neuro) and without it (PID). The last column shows the absolute error. The best result for each component and trajectory is boldfaced.Table 5
Comparison of the MSE of the tracking error for different trajectories with wind disturbances.
Trajectory M S E X M S E Y M S E Z M S E T Neuro PID Neuro PID Neuro PID Neuro PID Linear 0.0707 0.0707 0.0708 0.0708 1.6163 2.1214 0.6730 1.0130 Circular 0.2805 0.8612 0.0739 0.2596 1.6186 2.2711 0.7591 1.6014 Helical 0.0951 0.4850 0.0042 0.2009 0.0048 0.2847 0.2181 0.8666 Cyclic helical 0.1053 0.3435 0.0132 0.2000 0.0069 0.2596 0.2746 0.8084 Lemniscate 0.0911 0.1928 0.0011 0.0335 0.0005 0.3036 0.1748 0.6008 Helical lemniscate 0.0940 0.2614 0.0032 0.0695 0.0010 0.2460 0.2015 0.6903Table6 compares the values of the maximum error in this experiment for the same cases with before.Table 6
Comparison of the MAX of the tracking error for different trajectories with wind disturbances.
Trajectory M A X X M A X Y M A X Z M A X T Neuro PID Neuro PID Neuro PID Neuro PID Linear 1.0270 1.0270 1.0270 1.0270 5.0010 5.0010 5.1986 5.1986 Circular 2 2 1.0562 1.1864 5.0010 5.0010 5.4796 5.4796 Helical 1 1.4086 0.1834 0.9557 0.1582 1.4093 1.0033 1.6068 Cyclic helical 1 1.0948 0.2060 1.0506 0.1591 1.5356 1.0006 1.6874 Lemniscate 1 1 0.1035 0.5375 0.0710 1.6160 1 1.6166 Helical lemniscate 1 1.0582 0.1228 0.7613 0.0771 1.3190 1 1.3371According to Tables5 and 6, the controller with adaptive neural estimators provides less or equal tracking error for every tested trajectory. The minimal tracking error is achieved for the lemniscate trajectory in the z-axis, with a very small error value of 0.0005, 600 times less than the error value without the neural estimator. The maximum tracking error is obtained with the circular trajectory due to the high error of the initial conditions.The controller with the adaptive neural estimators works also well when the wind follows other different profiles, such as a sinusoidal function. In the next experiment, the wind average speed at Z = 20 m has been set to(44)vwZ=20=12+6∗sin⁡2π4t.Figure 20 shows the results. It can be seen how the performance of the controller without estimators is much worse than with the step wind profile. Nevertheless, the performance of the neurocontroller is like the previous ones. It is also possible to see the sinusoidal shape of the tracking error in Figure 20, right.Figure 20
UAV tracking of a helical lemniscate trajectory with sinusoidal wind profile disturbance (left) and its tracking error (right).


## 4.1. Trajectory Tracking without Disturbances

The control system has been first tested and validated with several trajectories without considering any disturbances. For each trajectory, the path followed by the UAV with adaptive neural estimators (red line), without them (yellow line), the reference (blue line), and the tracking error, are shown (Figures6–11). The tracking error is calculated by the following equation: (24)tERRti=Xi-Xrefi2+Yi-Yrefi2+Zi-Zrefi2,where Xrefi,Yrefi,Zrefi are the references for X, Y, and Z in ti, respectively.Figure 6
UAV tracking of a linear trajectory (left) and its tracking error (right).Figure 7
UAV tracking of a circular trajectory (left) and its tracking error (right).Figure 8
UAV tracking of a helical trajectory (left) and its tracking error (right).Figure 9
UAV tracking of a cyclic helical trajectory (left) and its tracking error (right).Figure 10
UAV tracking of a lemniscate trajectory (left) and its tracking error (right).Figure 11
UAV tracking of a helical lemniscate trajectory (left) and its tracking error (right).The trajectories used to validate the controller are described by the following equations, with the corresponding parameters:Linear:(25)XYZ=0.2t+10.2t+15∗stept+Z0.Circular:(26)XYZ=cos⁡0.4t+1sin⁡0.4t+15∗stept+Z0.Helical:(27)XYZ=cos⁡0.8t⁡sin⁡0.8t0.4t+Z0.Cyclic helical:(28)XYZ=cos⁡1.6t⁡sin⁡1.6tsin⁡0.4t⁡+Z0.Lemniscate:(29)XYZ=cos⁡0.4t⁡1+sin⁡0.4t⁡2cos⁡0.4tsin⁡0.4t1+sin⁡0.4t2Z0.Helical lemniscate:(30)XYZ=cos⁡0.8t1+sin⁡0.8t2cos⁡0.8tsin⁡0.8t1+sin⁡0.8t20.1t+Z0.In Figures 6–11, it is possible to observe how the tracking error decreases along the time. In t = 0, the error is high because the value of the reference starts in a very different value than the initial position of the system, that is, (X0,Y0,Z0)=(0,0,20). Then, the controller starts to work, and the tracking error is reduced.For example, in the linear trajectory (Figure6), the reference in t0 is (1, 1, 25); thus, the initial tracking error at t=0 is (1+1+25)=5.196.For the same reason, the initial error of the helical and lemniscate trajectories is much lower. The initial value reference is (1, 0, 20), and thus the initial error is 1.In these figures, there is not any significant visual difference between the results with and without the neural network due to the fact that if the PIDs are well tuned and there are no changes either in the mass or in the external disturbances, as it is the case in this section, the neural networks do not provide relevant advantages. Nevertheless, the networks in the control scheme of Figure4 show their full potential, when they must tackle disturbances, as it can be seen in the following sections. Consequently, the red and yellow lines seem to be overlapped.Other appreciable result is that the stationary tracking error is higher when the trajectory has helical component. Predictably, the circular and helical trajectories show a cyclic nature.As a conclusion, the controller performs a good tracking for a wide range of different trajectories. Thus, it can be said that the control strategy has been validated.


## 4.2. Control Robustness with Mass Variations

### 4.2.1. Mass Disturbance Model

Adding a payload in the quadrotor typically has three effects: the total mass of the system is increased, and the centre of the gravity can be modified and therefore also the inertia. In this work, we assume that the payload is an isotropic symmetric rigid solid attached to the UAV, not suspended, with dimensions much smaller than the dimensions of the quadrotor. The distance between the centre of gravity of the UAV and the centre of gravity of the load depends on the shape and the weight distribution of the manipulator and the shape and weight distribution of the load. In our experiment, we assume that this distance is zero. Under these circumstances, the effect of the inertia and centre of gravity variation can be neglected. Therefore, our paper is only addressing the effect of the mass variation.When the UAV is performing a logistic task, there are two possible stable states regarding the mass: one is before the load is in contact with the quadrotor, and then only the mass of the UAV is considered; the second one is when the payload mass is part of the system, and the sum of both masses is then considered as an only system. Between these two states, there may be several profiles of mass variation depending on the grasping technology and the properties of the surface of the load.We assume that the grasping and the load are nondeformable; thus, one step profile may be applied regarding the mass disturbance, as in other papers [16].The mass variation is simulated by adding a new termdistm (34) to (13) to (15), resulting in (31) to (34). The modelling of the mass disturbance is a step function. The total mass then is triplicated at t = 4s (34) meaning that a payload has been attached to the UAV. After 4 sec, the total mass of the system is m+mL=3m, with m being the mass of the quadrotor and the mass of the load mL=2m. Even if it can be considered a simple approach of dealing with this disturbance, the final effect is well represented:(31)X¨=-sin⁡θcos⁡ϕbm+distmu1(32)Y¨=sin⁡ϕbm+distmu1(33)Z¨=-g+cos⁡θcos⁡ϕbm+distmu1(34)distm=mL·stept-4.This experiment represents a possible situation while performing a logistic task, where the payload is heavier than the quadrotor itself.


### 4.2.2. Simulation Results with Mass Variation

Now, we test the control proposal with mass disturbances for the helical lemniscate trajectory because it is the most challenging one.Figure12, left, shows the reference trajectory in blue, the trajectory obtained by the controller with the neuromass estimator in red, and without the mass estimation in yellow (only PID). It is possible to observe how the trajectory obtained without the adaptive neural estimators (yellow line) moves away the reference even if it later comes closer to it. At the right side of Figure 12, the tracking error is presented (red line with mass estimation and blue line without mass estimation). In both figures, it is possible to observe how the performance with the mass estimator is much better. The tracking error is mostly the same until the mass changes at t = 4 s, but then the error without mass estimation increases significantly.Figure 12
UAV tracking of a helical lemniscate trajectory with mass variation (left) and its tracking error (right).The trajectory tracking of each coordinate has been represented in a separated figure to study the performance of the controller. Figure13, top left, shows the tracking of the X coordinate, at the top, right, the tracking of the Y coordinate, and at the bottom the tracking of the Z coordinate, with reference (magenta line), the trajectory with mass estimator (red line), and the trajectory without mass estimator (blue line). For every coordinate, it is possible to see how the trajectory obtained by the controller with the mass estimator fits better the reference. Even though, the deviation mainly comes from the Z coordinate due to the fact that the acceleration in the z-axis is more sensible to changes in the mass.Figure 13
UAV tracking trajectory of X, Y, and Z coordinates, respectively, with mass variation.Figure14 shows the mass obtained by the neuromass estimator in this experiment (real mass in blue and the estimated mass in red). The estimate of the mass is very similar to the real one but slightly noisy. Clearly, it can be appreciated that the mass is triplicated following a step variation profile.Figure 14
Mass obtained by the neural estimator.The controller is robust even with other types of variations in the mass like a sine function, for instance, if the mass changes according to(35)mt=mmodelstept-stept-4+stept-42·mmodel+mmodel·sin⁡2π8t.The results of considering this mass variation profile are shown in Figure 15. The tracking error of each component is like the previous ones, but the system experiences a stronger change with the sinusoidal reference profile. In this case, the tracking error of the system without the adaptive neural estimators tends to increase; meanwhile, in the case of the step change, the tracking error decreased over the time.Figure 15
UAV tracking of a helical lemniscate trajectory with sinusoidal mass variation (left) and its tracking error (right).Therefore, it can be seen how the controller with the neural estimator also works well with different mass change profiles.Figure16 shows the mass obtained by the neuromass estimator (red line) in this case and the corresponding real mass (blue line). The estimate of the mass is again very similar to the real one but slightly noisy but now the sinusoidal nature of the profile can be clearly observed.Figure 16
Mass obtained by the neural estimator with a sinusoidal mass change profile.The controller has been extensively tested for different trajectories and the MSE numerical results are summarized in Table3 for the different trajectories and for each coordinate, with the neural estimator (neuro) and without it (PID). The last column shows the absolute error. The best result for each component and trajectory is boldfaced.Table 3
Comparison of the MSE of the tracking error for different trajectories with mass variation.
Trajectory M S E X M S E Y M S E Z M S E T Neuro PID Neuro PID Neuro PID Neuro PID Linear 0.0707 0.0707 0.0708 0.0708 1.6099 21.9021 0.6635 3.8267 Circular 0.2766 0.2775 0.0720 0.0729 1.6101 21.9026 0.7076 3.8657 Helical 0.0934 0.1053 0.0018 0.0148 0.0029 18.2737 0.1899 3.2131 Cyclic helical 0.1067 0.2145 0.0126 0.1212 0.0046 17.5445 0.2715 3.1897 Lemniscate 0.0902 0.0916 0.003 0.0020 0.0010 18.2297 0.1589 3.1981 Helical lemniscate 0.0923 0.1060 0.0018 0.0166 0.0034 18.2335 0.1956 3.2087Table4 compares the values of the maximum error in this experiment for the same cases with before.Table 4
Comparison of the MAX of the tracking error for different trajectories with mass variation.
Trajectory M A X X M A X Y M A X Z M A X T Neuro PID Neuro PID Neuro PID Neuro PID Linear 1.0270 1.0270 1.0270 1.0270 5.0010 9.7835 5.1986 9.7835 Circular 2 2 1.0562 1.0562 5.0010 9.7836 5.4796 9.7837 Helical 1 1 0.1207 0.2074 0.1582 9.2655 1.0033 9.2675 Cyclic helical 1 1 0.2060 0.6163 0.1591 9.1763 1.0006 9.1937 Lemniscate 1 1 0.0614 0.0962 0.0913 9.2541 1 9.2544 Helical lemniscate 1 1 0.1128 0.2703 0.1322 9.2549 1 9.2553In Tables3 and 4, it is possible to observe that the controller with adaptive neural estimators provides less or equal tracking error for every tested trajectory. A general trend in the “PID column” (without mass disturbance estimation) is that the worst tracking error is obtained for the Z component. One of the reasons of this may be the fact that the maximum amplitude of the reference signal is in the z-axis.


## 4.2.1. Mass Disturbance Model

Adding a payload in the quadrotor typically has three effects: the total mass of the system is increased, and the centre of the gravity can be modified and therefore also the inertia. In this work, we assume that the payload is an isotropic symmetric rigid solid attached to the UAV, not suspended, with dimensions much smaller than the dimensions of the quadrotor. The distance between the centre of gravity of the UAV and the centre of gravity of the load depends on the shape and the weight distribution of the manipulator and the shape and weight distribution of the load. In our experiment, we assume that this distance is zero. Under these circumstances, the effect of the inertia and centre of gravity variation can be neglected. Therefore, our paper is only addressing the effect of the mass variation.When the UAV is performing a logistic task, there are two possible stable states regarding the mass: one is before the load is in contact with the quadrotor, and then only the mass of the UAV is considered; the second one is when the payload mass is part of the system, and the sum of both masses is then considered as an only system. Between these two states, there may be several profiles of mass variation depending on the grasping technology and the properties of the surface of the load.We assume that the grasping and the load are nondeformable; thus, one step profile may be applied regarding the mass disturbance, as in other papers [16].The mass variation is simulated by adding a new termdistm (34) to (13) to (15), resulting in (31) to (34). The modelling of the mass disturbance is a step function. The total mass then is triplicated at t = 4s (34) meaning that a payload has been attached to the UAV. After 4 sec, the total mass of the system is m+mL=3m, with m being the mass of the quadrotor and the mass of the load mL=2m. Even if it can be considered a simple approach of dealing with this disturbance, the final effect is well represented:(31)X¨=-sin⁡θcos⁡ϕbm+distmu1(32)Y¨=sin⁡ϕbm+distmu1(33)Z¨=-g+cos⁡θcos⁡ϕbm+distmu1(34)distm=mL·stept-4.This experiment represents a possible situation while performing a logistic task, where the payload is heavier than the quadrotor itself.


## 4.2.2. Simulation Results with Mass Variation

Now, we test the control proposal with mass disturbances for the helical lemniscate trajectory because it is the most challenging one.Figure12, left, shows the reference trajectory in blue, the trajectory obtained by the controller with the neuromass estimator in red, and without the mass estimation in yellow (only PID). It is possible to observe how the trajectory obtained without the adaptive neural estimators (yellow line) moves away the reference even if it later comes closer to it. At the right side of Figure 12, the tracking error is presented (red line with mass estimation and blue line without mass estimation). In both figures, it is possible to observe how the performance with the mass estimator is much better. The tracking error is mostly the same until the mass changes at t = 4 s, but then the error without mass estimation increases significantly.Figure 12
UAV tracking of a helical lemniscate trajectory with mass variation (left) and its tracking error (right).The trajectory tracking of each coordinate has been represented in a separated figure to study the performance of the controller. Figure13, top left, shows the tracking of the X coordinate, at the top, right, the tracking of the Y coordinate, and at the bottom the tracking of the Z coordinate, with reference (magenta line), the trajectory with mass estimator (red line), and the trajectory without mass estimator (blue line). For every coordinate, it is possible to see how the trajectory obtained by the controller with the mass estimator fits better the reference. Even though, the deviation mainly comes from the Z coordinate due to the fact that the acceleration in the z-axis is more sensible to changes in the mass.Figure 13
UAV tracking trajectory of X, Y, and Z coordinates, respectively, with mass variation.Figure14 shows the mass obtained by the neuromass estimator in this experiment (real mass in blue and the estimated mass in red). The estimate of the mass is very similar to the real one but slightly noisy. Clearly, it can be appreciated that the mass is triplicated following a step variation profile.Figure 14
Mass obtained by the neural estimator.The controller is robust even with other types of variations in the mass like a sine function, for instance, if the mass changes according to(35)mt=mmodelstept-stept-4+stept-42·mmodel+mmodel·sin⁡2π8t.The results of considering this mass variation profile are shown in Figure 15. The tracking error of each component is like the previous ones, but the system experiences a stronger change with the sinusoidal reference profile. In this case, the tracking error of the system without the adaptive neural estimators tends to increase; meanwhile, in the case of the step change, the tracking error decreased over the time.Figure 15
UAV tracking of a helical lemniscate trajectory with sinusoidal mass variation (left) and its tracking error (right).Therefore, it can be seen how the controller with the neural estimator also works well with different mass change profiles.Figure16 shows the mass obtained by the neuromass estimator (red line) in this case and the corresponding real mass (blue line). The estimate of the mass is again very similar to the real one but slightly noisy but now the sinusoidal nature of the profile can be clearly observed.Figure 16
Mass obtained by the neural estimator with a sinusoidal mass change profile.The controller has been extensively tested for different trajectories and the MSE numerical results are summarized in Table3 for the different trajectories and for each coordinate, with the neural estimator (neuro) and without it (PID). The last column shows the absolute error. The best result for each component and trajectory is boldfaced.Table 3
Comparison of the MSE of the tracking error for different trajectories with mass variation.
Trajectory M S E X M S E Y M S E Z M S E T Neuro PID Neuro PID Neuro PID Neuro PID Linear 0.0707 0.0707 0.0708 0.0708 1.6099 21.9021 0.6635 3.8267 Circular 0.2766 0.2775 0.0720 0.0729 1.6101 21.9026 0.7076 3.8657 Helical 0.0934 0.1053 0.0018 0.0148 0.0029 18.2737 0.1899 3.2131 Cyclic helical 0.1067 0.2145 0.0126 0.1212 0.0046 17.5445 0.2715 3.1897 Lemniscate 0.0902 0.0916 0.003 0.0020 0.0010 18.2297 0.1589 3.1981 Helical lemniscate 0.0923 0.1060 0.0018 0.0166 0.0034 18.2335 0.1956 3.2087Table4 compares the values of the maximum error in this experiment for the same cases with before.Table 4
Comparison of the MAX of the tracking error for different trajectories with mass variation.
Trajectory M A X X M A X Y M A X Z M A X T Neuro PID Neuro PID Neuro PID Neuro PID Linear 1.0270 1.0270 1.0270 1.0270 5.0010 9.7835 5.1986 9.7835 Circular 2 2 1.0562 1.0562 5.0010 9.7836 5.4796 9.7837 Helical 1 1 0.1207 0.2074 0.1582 9.2655 1.0033 9.2675 Cyclic helical 1 1 0.2060 0.6163 0.1591 9.1763 1.0006 9.1937 Lemniscate 1 1 0.0614 0.0962 0.0913 9.2541 1 9.2544 Helical lemniscate 1 1 0.1128 0.2703 0.1322 9.2549 1 9.2553In Tables3 and 4, it is possible to observe that the controller with adaptive neural estimators provides less or equal tracking error for every tested trajectory. A general trend in the “PID column” (without mass disturbance estimation) is that the worst tracking error is obtained for the Z component. One of the reasons of this may be the fact that the maximum amplitude of the reference signal is in the z-axis.


## 4.3. Control Robustness with Wind Disturbances

### 4.3.1. Wind Disturbance Model

The effect of the wind in an UAV in flight can be considered as a drag force and torque depending on the aerodynamics properties of the object and the characteristics of the air flow.The wind speed is a vector field; that is, its value may be different in each coordinate (x, y, z) of the space. For small UAVs like the one we are using in the experiments, we can assume a planar air flow, so the wind is the same in the region of the space where the UAV is flying. Under these circumstances, the wind influence in the torque can be neglected.On the other hand, the following equation represents the drag force [36]:(36)FD-=A-·Cd·ρair·V-22,where Cd is the drag coefficient, which is usually determined experimentally and collects the complex dependency, A is the area exposed to the wind, ρair is the density of the air, and V is the relative velocity of the flying object with respect to the wind. Since the air flow is planar, it can be assumed that the drag force is fully exerted in the centre of the gravity of the vehicle.In our work, this equation is divided by the mass. This way the expression is transformed to an acceleration ((40) to (43)) that can be easily introduced in the equations of the translational dynamics ((37) to (39)), where wind disturbances have been represented as distwX, distwY, and distwZ.As the wind speed is a vectorial magnitude, in our experiment, we suppose the wind components X and Y are twice the Z one. Due to this assumption, there appears a factor of 2 multiplying the wind speed in (41) and (42):(37)X¨=-distwX-sin⁡θcos⁡ϕbmu1(38)Y¨=-distwY+sin⁡ϕbmu1(39)Z¨=-distwZ-g+cos⁡θcos⁡ϕbmu1.In the free atmosphere, the wind is a balance between the Coriolis, centrifugal, and pressure forces acting on the air mass. But in the boundary layer, the mean wind velocity is also a function of the height [37], and it can be expressed as a logarithmic function. This fact is shown in the following:(40)vwZ=vwZ=20·log⁡Z/Clog⁡20/C(41)distwX=sgn⁡vwZ·ρair·Ax·Cd·X˙-2∗vwZ22m(42)distwY=sgn⁡vwZ·ρair·Ay·Cd·Y˙-2∗vwZ22m(43)distwZ=sgn⁡vwZ·ρair·Az·Cd·Z˙-vwZ22m,where vw(Z=20) is the wind speed at 20 m of altitude in m/s, vw(Z) is the wind speed at Z altitude in m/s, C is a constant related to the flight (in this experiment, the value is set to 1.5) which is dimensionless, ρair is the air density in Kg·m3, Ax to Az are the effective area of the quadrotor exposed to each component of the wind in m2, Cd is the drag coefficient with respect to the wind which is dimensionless, X˙, Y˙, and Z˙ are the velocities in the x-axis, y-axis, and z-axis in m/s, and sgn denotes the sign function.The mean wind velocity can be considered constant during the experiments, but not its instant value. The most common assumption is to consider the wind turbulence as a stationary Gaussian random process [37]. Considering this fact, in our experiment, the wind speed is simulated by a step with Gaussian noise at t = 4 s. The SNR between the average wind and the noise is 10 dB. The average wind speed is 12 m/s in the z-axis and 24 m/s in the x-axis and y-axis. These values match numbers 6 and 9 in Beaufort’s scale (strong breeze and strong gale) [38].


### 4.3.2. Simulation Results with Wind Variation

Figure17, left, shows the tracking of the trajectory, the reference in blue, the trajectory with adaptive neural estimators in red, and without estimators in yellow. At the right side of Figure 17, the tracking error is shown with estimators (red line) and without estimators (blue line). In both figures, it is possible to see how the performance of the control strategy with the neural estimators is much better. Like in the case of the mass variation, the tracking error is the same until t = 4 s because before there is no wind disturbance, but from that moment on, the tracking error of the controller without the estimators increases significantly.Figure 17
UAV tracking of a helical lemniscate trajectory with wind disturbance (left) and its tracking error (right).If the tracking according to each coordinate, X, Y, and Z, is represented (Figure18), the contribution to the tracking error seems to be more balanced (reference, magenta; trajectory with adaptive neural estimators, red; trajectory without neural estimators, blue). For every coordinate, it is possible to see how the trajectory obtained by the controller with the neural estimators better fits the reference. Figure 18 also shows how the disturbance produces an important deviation around t=7s due to the big peak of disturbance in the y-axis which cannot be compensated without the adaptive neural estimator and it becomes a stationary error. This deviation can be also easily observed in Figure 17 in the PID line.Figure 18
UAV tracking of X, Y, and Z coordinates, respectively, with wind disturbances.Another interesting result which can be drawn from the previous figures is that we are not using a wind disturbance estimator for the Z coordinate (see Figure4) but the tracking is still good enough. The reason is that adaptive neuromass estimator senses the wind disturbances in the z-axis as a virtual mass variation and it can compensate it. This is shown in Figure 19, where the estimate of the wind disturbance in the x-, y-, and z-axis is interpreted by the mass estimator as a mass variation. Indeed, the neural estimator (red line) fits reasonably well the real disturbances (blue line).Figure 19
Wind disturbance estimation in the x- and y-axis, top, and mass estimation (bottom).The controller has been extensively tested for different trajectories and the MSE numerical results are summarized in Table5 for the different trajectories and for each coordinate, with the neural estimator (neuro) and without it (PID). The last column shows the absolute error. The best result for each component and trajectory is boldfaced.Table 5
Comparison of the MSE of the tracking error for different trajectories with wind disturbances.
Trajectory M S E X M S E Y M S E Z M S E T Neuro PID Neuro PID Neuro PID Neuro PID Linear 0.0707 0.0707 0.0708 0.0708 1.6163 2.1214 0.6730 1.0130 Circular 0.2805 0.8612 0.0739 0.2596 1.6186 2.2711 0.7591 1.6014 Helical 0.0951 0.4850 0.0042 0.2009 0.0048 0.2847 0.2181 0.8666 Cyclic helical 0.1053 0.3435 0.0132 0.2000 0.0069 0.2596 0.2746 0.8084 Lemniscate 0.0911 0.1928 0.0011 0.0335 0.0005 0.3036 0.1748 0.6008 Helical lemniscate 0.0940 0.2614 0.0032 0.0695 0.0010 0.2460 0.2015 0.6903Table6 compares the values of the maximum error in this experiment for the same cases with before.Table 6
Comparison of the MAX of the tracking error for different trajectories with wind disturbances.
Trajectory M A X X M A X Y M A X Z M A X T Neuro PID Neuro PID Neuro PID Neuro PID Linear 1.0270 1.0270 1.0270 1.0270 5.0010 5.0010 5.1986 5.1986 Circular 2 2 1.0562 1.1864 5.0010 5.0010 5.4796 5.4796 Helical 1 1.4086 0.1834 0.9557 0.1582 1.4093 1.0033 1.6068 Cyclic helical 1 1.0948 0.2060 1.0506 0.1591 1.5356 1.0006 1.6874 Lemniscate 1 1 0.1035 0.5375 0.0710 1.6160 1 1.6166 Helical lemniscate 1 1.0582 0.1228 0.7613 0.0771 1.3190 1 1.3371According to Tables5 and 6, the controller with adaptive neural estimators provides less or equal tracking error for every tested trajectory. The minimal tracking error is achieved for the lemniscate trajectory in the z-axis, with a very small error value of 0.0005, 600 times less than the error value without the neural estimator. The maximum tracking error is obtained with the circular trajectory due to the high error of the initial conditions.The controller with the adaptive neural estimators works also well when the wind follows other different profiles, such as a sinusoidal function. In the next experiment, the wind average speed at Z = 20 m has been set to(44)vwZ=20=12+6∗sin⁡2π4t.Figure 20 shows the results. It can be seen how the performance of the controller without estimators is much worse than with the step wind profile. Nevertheless, the performance of the neurocontroller is like the previous ones. It is also possible to see the sinusoidal shape of the tracking error in Figure 20, right.Figure 20
UAV tracking of a helical lemniscate trajectory with sinusoidal wind profile disturbance (left) and its tracking error (right).


## 4.3.1. Wind Disturbance Model

The effect of the wind in an UAV in flight can be considered as a drag force and torque depending on the aerodynamics properties of the object and the characteristics of the air flow.The wind speed is a vector field; that is, its value may be different in each coordinate (x, y, z) of the space. For small UAVs like the one we are using in the experiments, we can assume a planar air flow, so the wind is the same in the region of the space where the UAV is flying. Under these circumstances, the wind influence in the torque can be neglected.On the other hand, the following equation represents the drag force [36]:(36)FD-=A-·Cd·ρair·V-22,where Cd is the drag coefficient, which is usually determined experimentally and collects the complex dependency, A is the area exposed to the wind, ρair is the density of the air, and V is the relative velocity of the flying object with respect to the wind. Since the air flow is planar, it can be assumed that the drag force is fully exerted in the centre of the gravity of the vehicle.In our work, this equation is divided by the mass. This way the expression is transformed to an acceleration ((40) to (43)) that can be easily introduced in the equations of the translational dynamics ((37) to (39)), where wind disturbances have been represented as distwX, distwY, and distwZ.As the wind speed is a vectorial magnitude, in our experiment, we suppose the wind components X and Y are twice the Z one. Due to this assumption, there appears a factor of 2 multiplying the wind speed in (41) and (42):(37)X¨=-distwX-sin⁡θcos⁡ϕbmu1(38)Y¨=-distwY+sin⁡ϕbmu1(39)Z¨=-distwZ-g+cos⁡θcos⁡ϕbmu1.In the free atmosphere, the wind is a balance between the Coriolis, centrifugal, and pressure forces acting on the air mass. But in the boundary layer, the mean wind velocity is also a function of the height [37], and it can be expressed as a logarithmic function. This fact is shown in the following:(40)vwZ=vwZ=20·log⁡Z/Clog⁡20/C(41)distwX=sgn⁡vwZ·ρair·Ax·Cd·X˙-2∗vwZ22m(42)distwY=sgn⁡vwZ·ρair·Ay·Cd·Y˙-2∗vwZ22m(43)distwZ=sgn⁡vwZ·ρair·Az·Cd·Z˙-vwZ22m,where vw(Z=20) is the wind speed at 20 m of altitude in m/s, vw(Z) is the wind speed at Z altitude in m/s, C is a constant related to the flight (in this experiment, the value is set to 1.5) which is dimensionless, ρair is the air density in Kg·m3, Ax to Az are the effective area of the quadrotor exposed to each component of the wind in m2, Cd is the drag coefficient with respect to the wind which is dimensionless, X˙, Y˙, and Z˙ are the velocities in the x-axis, y-axis, and z-axis in m/s, and sgn denotes the sign function.The mean wind velocity can be considered constant during the experiments, but not its instant value. The most common assumption is to consider the wind turbulence as a stationary Gaussian random process [37]. Considering this fact, in our experiment, the wind speed is simulated by a step with Gaussian noise at t = 4 s. The SNR between the average wind and the noise is 10 dB. The average wind speed is 12 m/s in the z-axis and 24 m/s in the x-axis and y-axis. These values match numbers 6 and 9 in Beaufort’s scale (strong breeze and strong gale) [38].


## 4.3.2. Simulation Results with Wind Variation

Figure17, left, shows the tracking of the trajectory, the reference in blue, the trajectory with adaptive neural estimators in red, and without estimators in yellow. At the right side of Figure 17, the tracking error is shown with estimators (red line) and without estimators (blue line). In both figures, it is possible to see how the performance of the control strategy with the neural estimators is much better. Like in the case of the mass variation, the tracking error is the same until t = 4 s because before there is no wind disturbance, but from that moment on, the tracking error of the controller without the estimators increases significantly.Figure 17
UAV tracking of a helical lemniscate trajectory with wind disturbance (left) and its tracking error (right).If the tracking according to each coordinate, X, Y, and Z, is represented (Figure18), the contribution to the tracking error seems to be more balanced (reference, magenta; trajectory with adaptive neural estimators, red; trajectory without neural estimators, blue). For every coordinate, it is possible to see how the trajectory obtained by the controller with the neural estimators better fits the reference. Figure 18 also shows how the disturbance produces an important deviation around t=7s due to the big peak of disturbance in the y-axis which cannot be compensated without the adaptive neural estimator and it becomes a stationary error. This deviation can be also easily observed in Figure 17 in the PID line.Figure 18
UAV tracking of X, Y, and Z coordinates, respectively, with wind disturbances.Another interesting result which can be drawn from the previous figures is that we are not using a wind disturbance estimator for the Z coordinate (see Figure4) but the tracking is still good enough. The reason is that adaptive neuromass estimator senses the wind disturbances in the z-axis as a virtual mass variation and it can compensate it. This is shown in Figure 19, where the estimate of the wind disturbance in the x-, y-, and z-axis is interpreted by the mass estimator as a mass variation. Indeed, the neural estimator (red line) fits reasonably well the real disturbances (blue line).Figure 19
Wind disturbance estimation in the x- and y-axis, top, and mass estimation (bottom).The controller has been extensively tested for different trajectories and the MSE numerical results are summarized in Table5 for the different trajectories and for each coordinate, with the neural estimator (neuro) and without it (PID). The last column shows the absolute error. The best result for each component and trajectory is boldfaced.Table 5
Comparison of the MSE of the tracking error for different trajectories with wind disturbances.
Trajectory M S E X M S E Y M S E Z M S E T Neuro PID Neuro PID Neuro PID Neuro PID Linear 0.0707 0.0707 0.0708 0.0708 1.6163 2.1214 0.6730 1.0130 Circular 0.2805 0.8612 0.0739 0.2596 1.6186 2.2711 0.7591 1.6014 Helical 0.0951 0.4850 0.0042 0.2009 0.0048 0.2847 0.2181 0.8666 Cyclic helical 0.1053 0.3435 0.0132 0.2000 0.0069 0.2596 0.2746 0.8084 Lemniscate 0.0911 0.1928 0.0011 0.0335 0.0005 0.3036 0.1748 0.6008 Helical lemniscate 0.0940 0.2614 0.0032 0.0695 0.0010 0.2460 0.2015 0.6903Table6 compares the values of the maximum error in this experiment for the same cases with before.Table 6
Comparison of the MAX of the tracking error for different trajectories with wind disturbances.
Trajectory M A X X M A X Y M A X Z M A X T Neuro PID Neuro PID Neuro PID Neuro PID Linear 1.0270 1.0270 1.0270 1.0270 5.0010 5.0010 5.1986 5.1986 Circular 2 2 1.0562 1.1864 5.0010 5.0010 5.4796 5.4796 Helical 1 1.4086 0.1834 0.9557 0.1582 1.4093 1.0033 1.6068 Cyclic helical 1 1.0948 0.2060 1.0506 0.1591 1.5356 1.0006 1.6874 Lemniscate 1 1 0.1035 0.5375 0.0710 1.6160 1 1.6166 Helical lemniscate 1 1.0582 0.1228 0.7613 0.0771 1.3190 1 1.3371According to Tables5 and 6, the controller with adaptive neural estimators provides less or equal tracking error for every tested trajectory. The minimal tracking error is achieved for the lemniscate trajectory in the z-axis, with a very small error value of 0.0005, 600 times less than the error value without the neural estimator. The maximum tracking error is obtained with the circular trajectory due to the high error of the initial conditions.The controller with the adaptive neural estimators works also well when the wind follows other different profiles, such as a sinusoidal function. In the next experiment, the wind average speed at Z = 20 m has been set to(44)vwZ=20=12+6∗sin⁡2π4t.Figure 20 shows the results. It can be seen how the performance of the controller without estimators is much worse than with the step wind profile. Nevertheless, the performance of the neurocontroller is like the previous ones. It is also possible to see the sinusoidal shape of the tracking error in Figure 20, right.Figure 20
UAV tracking of a helical lemniscate trajectory with sinusoidal wind profile disturbance (left) and its tracking error (right).


## 5. Conclusions

Intelligent control strategies are especially useful when the parameters change while the system is performing some tasks, and when the external disturbance are relevant, due to the ability to learn and adapt to the changing conditions. This is a common situation in many UAV applications.In this work, a new intelligent control strategy based on neural networks has been proposed. It includes the design of neural networks that estimate the system parameter variations. They allow the UAV to follow different trajectories with small tracking error when disturbances due to mass changes and wind are included.Simulation results show how the online learning of the neural estimators increases the robustness of the controllers, reducing the effects of the mass variation and the wind on the UAV.Among other possible future works, we may highlight the study of the influence of other disturbances such as the ones generated by the engines. In addition, the analysis of the parallelization of this approach for real-time application could be another interesting research line and help to deal with the high computational demand of these systems.



---

*Source: 6460156-2019-01-03.xml*