An Extended Kalman Filter-Based Technique for On-Line Identification of Unmanned Aerial System Parameters

: The present article deals with the identification, at the same time, of aircraft stability and control parameters taking into account dynamic damping derivatives. Such derivatives, due to the rate of change of the angle of attack, are usually neglected. So the damping characteristics of aircraft dynamics are attributed only on pitch rate derivatives. To cope with the dynamic effects of these derivatives, authors developed devoted procedures to estimate them. In the present paper, a complete model of aerodynamic coefficients has been tuned-up to identify simultaneously the whole set of derivatives. Besides, in spite of the employed reduced order model and/or decoupled dynamics, a six degrees of freedom model has been postulated without decoupling longitudinal and lateral dynamics. A recursive non-linear filtering approach via Extended Kalman Filter is proposed, and the filter tuning is performed by inserting the effects of dynamic derivatives into the mentioned mathematical model of the studied aircraft. The tuned-up procedure allows determining with noticeable precision the stability and control derivatives. In fact, either by activating maneuvers generated by all the control surfaces or by inserting noticeable measurement noise, the identified derivatives show very small values of standard deviation. The present study shows the possibility to identify simultaneously the aircraft derivatives without using devoted procedures and decoupled dynamics. The proposed technique is particularly suited for on-line parametrical identification of Unmanned Aerial Systems. In fact, to estimate both state and aircraft parameters, low power and time are required even using measurement noises typical of low-cost sensors.

ABSTRACT: The present article deals with the identification, at the same time, of aircraft stability and control parameters taking into account dynamic damping derivatives. Such derivatives, due to the rate of change of the angle of attack, are usually neglected. So the damping characteristics of aircraft dynamics are attributed only on pitch rate derivatives. To cope with the dynamic effects of these derivatives, authors developed devoted procedures to estimate them. In the present paper, a complete model of aerodynamic coefficients has been tuned-up to identify simultaneously the whole set of derivatives. Besides, in spite of the employed reduced order model and/or decoupled dynamics, a six degrees of freedom model has been postulated without decoupling longitudinal and lateral dynamics. A recursive non-linear filtering approach via Extended Kalman Filter is proposed, and the filter tuning is performed by inserting the effects of dynamic derivatives into the mentioned mathematical model of the studied aircraft. The tuned-up procedure allows determining with noticeable precision the stability and control derivatives. In fact, either by activating maneuvers generated by all the control surfaces or by inserting noticeable measurement noise, the identified derivatives show very small values of standard deviation. The present study shows the possibility to identify simultaneously the aircraft derivatives without using devoted procedures and decoupled dynamics. The proposed technique is particularly suited for on-line parametrical identification of Unmanned Aerial Systems. In fact, to estimate both state and aircraft parameters, low power and time are required even using measurement noises typical of low-cost sensors.
KEYWORDS: Aircraft dynamic derivatives, Extended Kalman Filter, On-line identification, Unmanned Aerial System.

INTRODUCTION
In spite of their effects on aircraft damping characteristics, the alpha-dot derivatives (C Lά and C mά ) are often not included into the dynamic model of aircraft. Therefore, in developing system identification techniques, neither off-line nor on-line procedures take into account dynamic stability derivatives (Tischler and Remple 2006;Klein and Morelli 2006;Jategaonkar 2006).
Usually, due to the strong correlation between pitch rate and alpha-dot derivatives, only C Lq and C mq are inserted into the mathematical model of the aircraft. In fact, because of such derivatives depend either on the aerodynamic transient or on the down wash lag, to perform their determination, especially devoted tests have to be designed. Besides, the alpha rate is a non-measurable variable, therefore on-line parameter identification procedures have to perform numerical derivation to determine such a quantity.
Besides, in the application of identification theories to Unmanned Aerial Systems (UAS or Remote Piloted Vehicles -RPV), reduced order models of aircraft dynamics are employed.
Thus, longitudinal and lateral aircraft parameters, usually, are identified separately by using devoted flight test procedures. Rimal et al. (2009)   (1) a non-linear global model. They have used model reference control architecture to simulate the non-linear identification and control of the system. Jameson and Cooke (2012) have proposed a post-maneuver parameters estimation with postulated models for the reduced order Short Period and Dutch Roll modes. Equation Error Method in the frequency domain has been applied and validated with Cranfield Jetstream 31 flight data. Kallapur and Anavatti (2006) have made a parametrical estimation in the time domain through Extended Kalman Filter (EKF) but only considering the three moment equations.
Other authors decoupled longitudinal and lateral dynamics. Dorobantu et al. (2013) have made parametrical identification of a 6-DoF small UAV through frequency domain method: they started from flight data and calculated parameters in frequency domain; identified models are then validated with flight data in the time domain. Nicolosi et al. (2010) have estimated aircraft stability derivatives from acquired flight data using the Output Error Method (OEM) technique; in this study, longitudinal and lateral dynamics are decoupled and they have used two of the lateral dynamics to determine the derivatives.
In a previous paper (Grillo and Montano 2014), the authors, by using the Filter Error Method (FEM) approach, have developed a procedure to determine aircraft longitudinal parameters including the angle of attack variation derivatives.
The aim of this paper is to identify, without decoupling longitudinal and lateral dynamics, the whole set of stability and control parameters of a non-linear 6-DoF model of UAV including dynamic derivatives. The original contribution of the present study is: (i) the simultaneously identification of longitudinal and lateral stability and control derivatives; (ii) the identification of dynamic derivatives without using a devoted procedure.
Because of the EKF has been successfully employed for off-line and on-line identification of stability and control derivatives for manned aircraft (Alonge et al. 2015;Shim et al. 2014, Kokolios 1994Speyer and Crues 1987;Kobayashi et al. 2005;Huang et al. 1996;Jategaonkar andPlaetenschke 1988, 1989;García-Velo and Walker 1997), by considering the augmented computation capability of on-board UAS, a procedure based on EKF has been implemented.
Such a procedure could be successfully employed for real-time system identification leading to costs and time reduction in design and development of UAS.

AIRCRAFT MATHEMATICAL MODEL
As previous stated, the aim of this paper is the identification of unknown aerodynamic parameters for an UAS.
The take-off weight of the analyzed vehicle is 24.5 N. It is equipped with: In this study, the rigid body dynamics of a 6-DoF aircraft is of interest. The aircraft motion can be described by the following equations in body axes (Etkin 1972): where: where: u, v, and w are velocity components in body reference frame; p, q, r are roll, pitch, and yaw angular speeds [rad/s]; φ, ϑ, ψ are Euler angles; m is the aircraft mass; I x , I y , I z are aircraft moments of inertia; V is the air speed [m/s]; C L , C D , C m , C T are the lift coefficient, drag coefficient, pitch moment coefficient, and thrust coefficient; C Lα , C mα , C Lq , C mq , C Lἀ , C mἀ , C Yβ , C lβ , C nβ , C Yp , C lp , C np , C Yr , C nr , C TV , are the stability derivatives; C l , C n , C Y are coefficients of roll moment, yaw moment, and side-force; C Lδe , C mδe , C Yδr , C lδa , C lδr , C nδa , C nδr , C Tδth , are control derivatives.
Notice that the over line individuates non-dimensional quantities; besides, in Eq. 10, altitude variation is considered.
Finally, instead of modeling the drag coefficient by means of stability derivatives, an analytical expression of the drag polar is employed. Such equation has been obtained by Computational Fluid Dynamics (CFD) analysis of the whole aircraft.
Based on these assumptions, the state of the system is given by x = [u v w p q e φ ϑ ψ h] T , while the set of inputs in = [δ e , δ th , δ a , δ r ] T is made up of the control surfaces and throttle deflections, where δ e is the elevator deflection; δ a is the aileron deflection; δ r is the rudder deflection and δ th is the throttle displacement.
Equations 1 -10 represent the aircraft state equations. Because of the on-board instrumentation is located very close to the center of gravity of the aircraft, no correction has to be made to the set of observation equations. Besides, to avoid errors due to poor quality of the pitot probe, instead of using anemometric data, the airspeed is obtained by integrating the IMU data. Finally, the heading ψ and the altitude h have been measured by using the magnetometer and the altimeter, respectively.
The observation equations are: where: a x , a y , and a z are the acceleration of aircraft mass center along the body-fixed x-axis, y-axis, and z-axis.
Based on these assumptions, the output vector of the system is given by y = [a x a y a z V p q r ψ h] T .

EXTENDED KALMAN FILTER DESIGN
As it is well known, the EKF has been successfully employed as recursive non-linear filtering algorithm for aircraft parameter estimation (Speyer and Crues 1987;Kobayashi et al. 2005;Jategaonkar and Plaetenschke 1988).
As previous stated, in the present study, an EKF has been designed to estimate both the state and stability and control derivatives for the UAS described by the model outlined in "Aircraft Mathematical Model" section.
To obtain joint estimation of state and parameters, the state vector has been augmented by defining the unknown parameters as additional state variables. Such a process leads to the following state vector: where: w(t) and w p (t) are, respectively, a 10 x 1 and a 22 x 1 random process with unknown statistics and t is time. These have been inserted to represent uncertainties in the system model. For the augmented model, the output equation is: where: x = [u v w p q e φ ϑ ψ h] T and p is the set of unknown stability and control derivatives that, according to Eqs. 1 -9, is: Discrete measurement equations have been used: with: where: ν(t) is the measurement noise vector and T is the sampling interval of the instrumentation.
According to the theory (Lewis et al. 2008), the structure of the continuous discrete Kalman filter consists of two subsystems: a subsystem (predictor) which gives a time update of the estimation and a subsystem (corrector) which adjusts the prediction by means of measurements.
A first estimate of the state is obtained by means of the known dynamics of the system: This estimate is the corrected using: Besides, we have postulated the constant system parameter p, and the following equations have been inserted into the UAS model: The following augmented dynamic model of the system has been obtained by joining Eqs. 1 -10 and Eq. 20: where: X pred are predicted states, X est are esteemed states, z is the measurement vector, and y pred is predicted output, with the initial condition: where: K k is the gain matrix at step k, which is computed to minimize the covariance of the estimation error: In the present paper, Eq. 24 has been integrated by using a second-order Runge-Kutta algorithm; the Jacobean matrices A and H have been evaluated by using a numerical implementation of the central difference formula; the state transition matrix has been obtained by Taylor series expansion: The estimation error covariance matrix P pred and P est are given by: with the initial condition P estk = P(0) . where R is the measurement noise covariance matrix, Q is the process noise covariance matrix and I is the identity matrix.
In the previous equations, ϕ k is the discrete time state transition matrix: and, according to Equation (21), Jacobeans A k and H k are defined as next: where: f is a non-linear relation between states, inputs, disturbances and time (see Eq. 21).
Equations 24 -28 describe the EKF for the system given by Eqs. 21 and 22. In particular, Eqs. 21 and 27 are the prediction equations and Eqs. 25, 26 and 28 are the correction equations.
Obviously, to perform the estimation process of the nonlinear system given by Eqs. 21 and 22, numerical integration algorithm has to be used for integration of Eq. 24.
Besides, because of Eqs. 29 and 30 denote the linearization process for the state matrix (Eq. 29) and the measurement matrix (Eq. 30), numerical linearization has to be performed.
Finally, the state transition matrix ϕ k must be approximated.
The core of the tuning of the EKF is the calibration of both measurement and process noise covariance matrix and the selection of P(0), which represents the confidence in the initial state estimates.
Physical airframe size and cost restrictions limit the availability and quality of UAS on-board sensors. Besides, such aircraft is often significantly affected by wind disturbance, both due to their relatively small dimensions and to their mission profiles. Finally, small wind components represent relevant non-measured inputs. Therefore, measured flight data contain considerable amount of noise.
To take into account either the poor quality of the on-board instruments or the relevant effects of the wind disturbance on the flight path characteristics, in the present paper, very high values have been used for the measurement noise covariance matrix R. In particular, instead of inserting into the measurement noise covariance matrix the characteristics of the on-board instruments, the elements of R have been selected by determining errors due to realistic values of the wind components.
Because of the UAS mathematical models are usually very inaccurate to take into account both high measurement errors in the aircraft outputs and poor knowledge of the parameter values, high values have been assumed for P(0). In this way, it is possible to test the EKF robustness in parameters estimation. Therefore, the measurement noise has been postulated as Gaussian band limited white noise. To simulate the influence of the atmospheric turbulence on the UAS flight, the amplitude of such a noise is comparable to realistic values of atmospheric wind components.
Finally, a trial and error procedure has been employed to determine the process noise covariance matrix Q by imposing satisfaction of requirements assigned in terms of mean value of estimation error for the non-augmented state of the aircraft.

FILTER TUNING
The identification procedure has been tuned-up by means of simulation of the non-linear model discussed in "Aircraft Mathematical Model" section. A first set of stability and control derivatives of the studied UAS has been calculated by means of analytical relations (Etkin 1972). Such a set represents the initial condition of the parameters into the augmented state.
To perform the parametrical identification, the non-linear model postulated in Eqs. 1 -10 has been employed to take into account every correlation between variables. In order to excite both dynamics, elevator, ailerons and rudder, deflections have been simultaneously applied. The comparison has been made between simulation results and estimates.
Simulation has been performed by selecting a trim condition of rectilinear horizontal flight at a reference speed V = 24.63 m/s, corresponding to an angle of attack α = 0.0471 rad and an elevator deflection δ e = -0.2610 rad.
The so-called 3-2-1 input has been selected (alternating pulses with time and width in the ratio 3-2-1).
To cope with transient characteristics of aircraft dynamic modes, elevator deflection has been chosen constant (Δδ e = ± 0.1 rad) and scaled on time (3-2-1 seconds). Such an input, as it is well known, is the best way to excite the short period mode. In fact, the short period and phugoid frequencies are: 2.4 x 10 -3 Hz and 4.44 x 10 -4 , respectively. Because of the dutch roll frequency is 7.53 x 10 -4 Hz, the rolling convergence and spiral mode time constants are 0.054 s and 2.04 s, respectively; aileron and rudder deflections have been chosen in 3-2-1 ratio amplitudes (δ a , δ r = ± 0.12, ± 0.08 and ± 0.04 rad; duration: 1s). The selected inputs are shown in Figs. 2 and 3.
The initial condition of augmented state has been formed by the trim value of state variables (x) and the analytical values The measurement vector (z) has been obtained by simulation, integrating state equations 1 -10. To simulate real measurements, a Gaussian band limited white noise has been added to the obtained results.
Process and measurement noises have been postulated as uncorrelated; in this way, either the process noise covariance matrix or the measurement noise covariance matrix are diagonal.
The selected constrains are presented in Table 1. At the end of the tuning procedure, the following values of the process noise covariance matrix have been obtained:  Q = diag [1; 1; 10; 1; 1; 1; 1; 1; 1; 1; 0.1; 0.01; 10; 10; 1; 1; 10; 10; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 10; 1] x 10 -2 The non-augmented state variables estimated with the tuned-up filter parameters are shown in  In these figures, the values of the state variables calculated by simulation (drawn in green line) and those obtained by performing state estimation (drawn in red line) are compared. Figures 4 -13 show very good estimation accuracy for all states. Mean errors are largely below the maximum values admitted by the imposed specification. Besides, the maximum estimation error is very small.      The obtained values of maximum and mean error are shown in Table 2. Grillo

RESULTS
The tuned-up EKF has been employed to determine the unknown aircraft parameters. To verify robustness properties of the tuned-up procedure, a Gaussian white noise has been added to the true outputs in order to better simulate real measurement.
The standard deviations of the measurement errors have been set as follows: In spite of the very poor qualities of the postulated measures, for the non-augmented state, the obtained mean errors are smaller than those shown in Table 2. The obtained results confirm the suitability of the tuning procedure. In fact, as it is known, small values of the process noise covariance matrix correspond to a better filtering of high-frequency disturbances.
The time evolution of the estimated aircraft damping parameters is shown in Figs. 14 -19.    Notice that these parameters are shown because of their relationship with the time constant of aircraft dynamic modes. These figures clearly show that the parameters go to their asymptotic value after few seconds.
To demonstrate the suitability of the proposed procedure, Figs. 20 and 21 show the fundamental aircraft parameters C Lα and C mα . Table 3 shows the identified values of the parameters compared with their analytical values, as well as the standard deviations.   The parametrical identification has been carried out by using low computational power and considering sensor characterized by high measurement noise.
The tuned-up procedure affords:

•
The simultaneous identification of longitudinal and lateral stability and control derivatives.

•
The identification of dynamic derivatives without using a devoted procedure.

•
The reduction of costs and time in performing parametrical identification. The second item represents a very interesting result. In fact, as it is well known, in off-line identification procedures (based on filter error methods or output error method, for example), the parametrical identification is performed by minimizing the selected cost function. Such a function is related to the measurement noise. Obviously, is not measurable in a simple way; therefore, devoted flight test campaign or devoted instruments are necessary to obtain valuable values of such a quantity.
Similarly, in on-line identification procedures to determine dynamic derivatives, it is necessary to identify first the angle of attack change rate.
By using the proposed procedure, it is sufficient to evaluate numerical values of by using measured values of α. The obtained results show that the implemented algorithm affords to determine aircraft parameters with noticeable precision.
Very similar values of the estimated aircraft parameters have been identified by varying the measurement noise characteristics. In particular, very small standard deviations have been obtained even in presence of noticeable measurement noise.
The tuned-up procedure is particularly suited for UAS because a low computational power is required; besides, lowcost sensors may be employed.
Because a non-linear 6-DoF mathematical model of the aircraft has been used, the tuned-up maneuvers afford to determine the whole set of aircraft stability and control derivatives, leading to reduction of costs and time in performing parametrical identification.
Moreover, the implemented algorithm has good robustness properties. In fact, by adding a Gaussian white noise with noticeable standard deviation to the system output, the obtained results show very small errors.
Finally, the ἀ derivatives (C Lά and C mά ), which often are not included into the dynamic model of aircraft, have been The obtained results show that it is possible to identify at the same time the whole set of the aircraft derivatives; besides, it is not necessary to tune-up dedicated procedure to estimate the dynamic damping derivatives. In fact, the obtained values of and , that are shown in Table 3, have very small values of standard deviation.
It is noticeable the negative value of the parameter , which is usually positive for conventional aircraft. This means that it will be necessary, in performing flight test campaign, to study devoted manoeuvres in order to determine the real sign of the mentioned derivative.

CONCLUSIONS
The identification process of a 6-DoF UAS model has been performed developing a procedure based on the Extended Kalman Filter.