Acta Polytechnica doi:10.14311/AP.2018.58.0077 Acta Polytechnica 58(2):77–91, 2018 © Czech Technical University in Prague, 2018 available online at http://ojs.cvut.cz/ojs/index.php/ap ESTIMATION OF THE LATERAL AERODYNAMIC COEFFICIENTS FOR SKYWALKER X8 FLYING WING FROM REAL FLIGHT-TEST DATA Rahman Mohammadi Farhadi∗, Vyacheslav Kortunov, Andrii Molchanov, Tatiana Solianyk National Aerospace University “KhAI”, 17, Chkalova st., Kharkiv, Ukraine ∗ corresponding author: rmfarhadi.ua@gmail.com Abstract. Stability and control derivatives of Skywalker X8 flying wing from flight-test data are estimated by using the combination of the output error and least square methods in the presence of the wind. Data is collected from closed loop flight tests with a proportional-integral-derivative (PID) controller that caused data co-linearity problems for the identification of the unmanned aerial vehicle (UAV) dynamic system. The data co-linearity problem is solved with a biased estimation via priori information, parameter fixing and constrained optimization, which uses analytical values of aerodynamic parameters, the level of the identifiability and sensitivity of the measurement vector to the parameters. Estimated aerodynamic parameters are compared with the theoretically calculated coefficients of the UAV, moreover, the dynamic model is validated with additional flight-test data and small covariances of the estimated parameters. Keywords: flight system identification; aerodynamic coefficient estimation; data co-linearity; least square method; output error method. 1. Introduction A system identification can be used to characterize forces and moments acting on a UAV. The UAV dy- namic system identification is applied to find the func- tional dependence of aerodynamic forces and moments for the UAV motion and control variables based on flight-test data. The aerodynamic characteristics of a UAV can be obtained by these methods: (1.) computational methods, (2.) wind-tunnel testing (3.) UAV dynamic system identification from flight test data. There are several advantages to perform the UAV dynamic system identification from flight-test data [1– 11]: • To verify and interpret theoretical, numerical and experimental (wind-tunnel test) predictions and results of the UAV dynamic characteristics; • To find more accurate and realistic uncertain math- ematical models of the UAV dynamics, used in the flight control systems design; • To develop flight simulators, which require accu- rate representation of the UAV dynamic model in different flight regimes. Many UAV manoeuvres and flight conditions cannot be conducted simply in wind tunnel nor by analytical or numerical com- putations with sufficient accuracy or efficiency; • To expand the flight envelope for a new UAV, which can include the stability quantification and control the impact of UAV modifications, configuration changes, or special flight conditions; • To verify the UAV specification complied with the flight conditions. There are four important principles known as Quad- M that must be considered in the identification of the flight vehicle (Figure 1): good Manoeuvres or a persistent excitation, sufficient Measurements, a selection of suitable Model and Method [5, 12]. In general, an estimation of the aerodynamic coefficients can be performed in two ways: (1.) Estimation after modelling, some methods in- clude: equation error method, least square method (LSM), output error method (OEM), optimization methods (such as genetic algorithm, simulated an- nealing, pattern search), Kalman filter, filter error method or maximum likelihood method (Kalman filter+OEM), model error estimation; (2.) Estimation before modelling: singular system approach for unknown inputs estimation, Gauss- Markov stochastic models for unknown inputs, com- bination of the Kalman filter and neural network, combination of stochastic models and an optimiza- tion methods [3–7, 9, 11–15]. Particular application of the dynamic model deter- mines its accuracy requirements and estimation or identification methods of the UAV dynamic model [1, 16]. A high accuracy level is required to process the pa- rameter monitoring and failure detection. A medium to high level of accuracy is required for the verification of the theoretical model (off-line). A low to medium 77 http://dx.doi.org/10.14311/AP.2018.58.0077 http://ojs.cvut.cz/ojs/index.php/ap R. Mohammadi Farhadi, V. Kortunov, A. Molchanov, T. Solianyk Acta Polytechnica Figure 1. Quad-M principles for flight system identi- fication. level of accuracy is required for the tuning of the control system parameters (off-line), adaptive control (on-line) and digital control system design (on-line and off-line) [1]. Usually, optimal persistent exciting inputs (ma- noeuvres) are designed so that an assumed aerody- namic model structure (usually linear) would be ad- equate to characterize the data for estimating the aerodynamic coefficients of the flight vehicles [3, 5– 7]. Flight tests are performed for the open loop conditions in the calm weather [3]. This elimi- nates the need to include the process noise in the model so that a simpler output-error method can be used in the data analysis and modelling. How- ever, it is not always possible to fulfil all these conditions. Usually, wind is present in the flight tests and persistent exciting manoeuvres cannot be designed. In addition, in some models, feedback presence can lead to completely incorrect results of system identification [5, 6, 13]. Therefore, pa- rameters estimates diverge from actual values or are multivalued. However, sometimes it is neces- sary to perform the identification in the presence of a feedback due to the instability of the object it- self, the economy, the presence of inherent feedback, safety, production considerations, etc. Therefore, it is necessary to estimate the unknown parameters for weakly excited inputs in wind presence. Persistent excited parts of the data must be chosen in order to identify the unknown parameters based on the optimization criterion such as the sensitivity of the system output to the unknown parameters [6, 11]. In this paper, the combination of the OEM, LSM, and constrained optimization will be used in the presence of the wind and weakly exciting manoeu- vres. It must be considered that the time history match is always a necessity, but not a sufficient condition for the system identification [5]. An unrealistic value or sign for some parameters can occur due to an inade- quate model or existence of some un-excited modes of the system. Inclusion of physical signs and bounds on parameters, fixing some parameters to priori values and interval estimation due to uncertainties must be taken into consideration in the practical application of parameter estimation methods based on the plant dy- namics and the allowable range of parameters [5, 13]. These problems lead to a correlation or a lack of identifiability, singularity of information matrix, un- realizable parameters or large error covariances. It means that all parameters cannot be estimated sepa- rately while some of their linear combinations can be estimated. Identifiability problems can be solved by fixing parameters, a priori weighting, constraint opti- mization, rank deficient solutions, input design and model structure selection [5, 7]. Simulated Annealing optimization algorithm (SAOA) (which imitates the annealing process of metals) is used for the constrained optimization algorithm [17]. By using the SAOA, the probability of getting trapped in local optimums is avoided by adding randomness to the acceptance of a better direction in optimization procedure and initial acceptance of direction with worse value of the cost function. The optimization procedure can be restarted from the obtained minimum with a high probability of accepting worse values of the cost function. Modelling the UAV’s dynamics of motion can be performed on the basis of benchmark models, i.e. the model, that has passed multilevel checks and refine- ment of the physical experiments results in wind tun- nels and flight tests. Since conducting a flight experi- ment is costly, for the development of identification algorithms, it is necessary to perform the simulation of the motion and take into account the noise and dis- placements of all sensors [3]. Then, the identification algorithms that have been tested for such models are recognized as reliable and suitable for a practical ap- plication. Therefore, linear models are generated for the nonlinear simulation of the Aerosonde, Skywalker X8 and a Zagi UAVs [18, 19]. The suggested methods are implemented and validated for the linear stan- dard models to estimate the aerodynamic parameters. Next, the combination of these methods are applied for estimating the aerodynamic coefficients of Sky- walker X8 flying wing (SX8FW) from real flight-test data. Aerodynamics coefficients of the SX8FW were cal- culated theoretically in [19]. In this paper, its aerody- namic coefficients are estimated based on flight-test data and compared with theoretical coefficients of the 78 vol. 58 no. 2/2018 Estimation of Lateral Aerodynamic Coefficients for Skywalker X8 Parameter ZFW SX8FW Mass m kg 1.56 4.5 Moment of inertia Jx kg m2 0.1147 0.45 Moment of inertia Jy kg m2 0.0576 0.325 Moment of inertia Jz kg m2 0.1712 0.75 Moment of inertia Jxz kg m2 0.0015 0.06 Wing area S m2 0.2589 0.75 Wing span b m 1.4224 2.12 Mean aerodynamic chord c m 0.3302 0.3571 Propeller area Sprop m2 0.0314 0.1018 Air density % kg/m3 1.2682 1.2682 Motor constant Kmotor – 20 40 Propeller aerodynamic coefficient Cprop – 1 0.5 Table 1. Airframe parameters for the ZFW and SX8FW. Figure 2. SX8FW (with the mini autopilot developed in at the National Aerospace University (KhAI)). SX8FW and the Zagi flying wing (ZFW). A special case of a parametric identification is a so-called stiff system identification [3]. The presence of the fast and slow modes in such systems creates this feature, which leads to a divergence of identification algorithms. In [3], aerodynamic parameters for a longitudinal mo- tion of the DHC-2 “Beaver” aircraft were estimated with the LSM and Kalman filter in the computer simulation environment. Practically, it is not pos- sible to use the LSM due to the measurement and process noise. In addition, it is difficult to tune the Kalman filter. Therefore, this study proposes a prac- tical approach to apply deterministic equations for the FW motion and OEM for estimating aerodynamic parameters from real flight-test data. The problem is defined in the second section. OEM and identifiability problems of aerodynamic coeffi- cients are discussed in the third section. Fourth sec- tion presents data co-linearity and its solutions for a lateral acceleration regression. All aerodynamic coeffi- cients for the lateral dynamics are estimated with five different approaches in the fifth section. The conclu- sion and some suggestions for future work are given in the last section. 2. Problem statement It is necessary to develop an identification procedure to estimate aerodynamic parameters of the SX8FW. This procedure needs to rely on the results of the flight experiments in the presence of wind and measurement noises. The measurements of Angular rates and lin- Figure 3. The ZFW UAV [18]. ear acceleration are measured with MEMS sensors, while the longitudinal and lateral ground velocities are measured by a GPS receiver. In addition, the ground velocity and Euler angles are calculated with an inertial measurement unit of a mini autopilot. Air- speed is measured with the airspeed sensor at the pivot location, which is installed in front of the UAV. The SX8FW and the ZFW are similar flying wings (FW). Therefore, the ZFW is chosen for comparison. Corresponding airframe parameters are shown in Fig- ures 2 and 3, Table 1 respectively [18, 19]. These FW belong to mini UAVs [20]. The SX8FW uses the mini autopilot AP-AVIA, which was developed at the National Aerospace University named after N. Zhukovsky (KhAI) [21, 22]. The errors of the AP-AVIA autopilot sensors are presented in work [22]. The frequency of the data recording can be 20 Hz (saved to an SD card by the autopilot) or 5 Hz (saved to a log file by the ground control station). The study of the influence of the sensors accuracy on the UAV dynamic model identification is not car- ried out in detail since the accuracy of angles, angu- lar velocities, and linear velocities seems to be suffi- cient for solving this problem. However, the standard deviations of these parameters were higher in flight due to the vibrations of the UAV engine and struc- ture. Body-axis and wind-axis reference frames for the 79 R. Mohammadi Farhadi, V. Kortunov, A. Molchanov, T. Solianyk Acta Polytechnica Figure 4. Body-axis and wind-axis reference frames of the UAV. UAV are shown in Figure 4. The motion of the FW is controlled by elevons, which are moved by the au- tomatic control system. A differentially angular de- flection of the right (δer) and left (δel) elevons has the same effect as aileron (δA) and a commonly angu- lar deflection of them has the same effect as elevator (δE) [18]: [ δE δA ] = [ 1 1 1 −1 ] [ δer δel ] . (1) Nonlinear dynamical equations of a motion in the presence of the wind can be written as [23]: [V ′A]B = fB–aero(. . .) + fB–thrust(δT ) m(t) +HB I gI − ωB × [VI ]B , dωB dt = I−1 B × ( mB–aero(. . .) +mB–thrust(δT ) − ωB × IbωB ) , dΘ dt = LEBωB , (2) where fB–aero = fB–aero(|VA|, α−αW , β−βW , ωB , δA, δE , δR), fB–thrust, mB–aero = mB–aero(|V |A, α − αW , β − βW , ωB , δA, δE , δR), mB–thrust are aerodynamic and thrust forces and moments in the UAV body frame, α, β are angles of attack and side slip, αW , βW are angles of attack and side slip due to wind in the inertial frame, VA, VI are airspeed and ground speed, ωB is angular velocity vector; δA, δE , δR, δT are control signals for aileron, elevator, rudder and thrust side slip angle due to wind, gI is the gravity vector in the inertial frame; IB is matrix of inertial moments, Θ is vector of the Euler angles, HB I and LEB are the corresponding rotation matrices. It means that aerodynamic forces and moments depend on airspeed where VA = VI −W = v1 v2 v3  I − W1 W2 W3  I , [VA]B = uAvA wA  = HB I VA −HB I W, |VA|βA αA  =  √ u2 A + v2 A + w2 A sin−1 vA VA tan−1 wA uA  , (3) where W is the wind vector in the inertial frame. The linearized lateral equations of motion, including the effect of the wind gusts, are: ∆xLatB = [ ∆v ∆p ∆r ∆φ ]T , ∆uLat = ∆δA, ∆wLat = ∆vW , ∆x′LatB = ALat∆xLatB +BLat∆uLat + ELat∆wLat, ∆yLatB = CLat∆xLatB +DLat∆uLat + FLat∆wLat, (4) where matrices ALat, BLat ELat are ALat =  Yv Lv Nv + YvNv′ 0 Lp + w0 Lp Np + (Yp + w0)Nv′ 1 Yr − u0 Lr Nr + (Yr − u0)Nv′ − sin θ0 gI cos θ0 Lφ 0 0  T , BLat =  YδA LδA NδA + YδANv′ 0  , ELat =  −ALat(1, 1) −ALat(2, 1) −ALat(3, 1) 0  =  −Yv −Lv −(Nv + YvNv′) 0  (5) 80 vol. 58 no. 2/2018 Estimation of Lateral Aerodynamic Coefficients for Skywalker X8 and matrices CLat, DLat FLat are CLat =  1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1  , DLat = ELat = [ 0 0 0 0 ]T . (6) It must be taken into account that the velocity com- ponents u, v, w in the state equations are computed at the centre of gravity. Therefore, all quantities in the state as well as observation equations should be defined with respect to the centre of gravity. Although the aircraft rates and the roll and pitch attitudes are not affected by the location of the centre of gravity, the measurements of linear accelerations and velocity components are influenced by the distance between the centre of gravity and the sensor position. The airspeed is measured at the pivot location which is in- stalled in front of the UAV. It must be considered that the airspeed is measured almost along the x-axis, so the normal and lateral components of the wind are not observable and these components only affect the mea- sured normal and lateral accelerations. However, the amplitude and angle of the wind in inertial frame are not changed rapidly; therefore, the UAV flies in four directions at the beginning of the flight to estimate the amplitude and angle of the wind. The amplitude and angle of the wind are practically used during the flight to compensate the effect of the wind. In other words, the lateral wind is estimated in advance. The airspeed components along the body frame axis at the sensor location are given by [5] us = u− rys + qzs, vs = v − pzs + rxs, ws = w − qxs + pys, (7) where [ xs ys zs ]T is the sensor location in the body frame. In addition, sometimes, linear accelerometers are not mounted exactly at the centre of gravity. The accelerations at the centre of gravity can be derived from the measured accelerations at the sensor location using the following equations: ax = axs + (q2 + r2)xs − (pq − r′)ys − (pr + q′)zs, ay = ays − (pq + r′)xs + (r2 + p2)ys − (qr − p′)zs, az = azs − (pr − q′)xs − (qr + p′)ys + (p2 + q2)zs, (8) where ax, ay, az are the accelerations at the cen- tre of gravity, axs, ays, azs the accelerations at the accelerometer location and [ xs ys ys ]T the ac- celerometer location in the body frame. Therefore ax = Xv(u− uW ) +Xpp+Xrr +XδEδE +XT δT , ay = Yv(v − vW ) + Ypp+ Yrr + YδAδA, az = Zv(w − wW ) + Zpp+ Zrr + ZδEδE + ZT δT . (9) It must be noted that there is no rudder control surface in the FW and in most cases, Yp amd Yr are assumed to be zero. 3. OEM and identifiability problem For a practical use in the aircraft parameter es- timation, no process noise is considered in the maximum-likelihood algorithm. Therefore, the practi- cal maximum-likelihood estimation algorithm can be developed for the state space model with a discrete measurement model in a discrete form. In case of no process noise, the negative log-likelihood function is formulated as follows [7]: J(Θ) = 1 2 N∑ k=1 eT (k)R−1(k,Θ)e(k) + N 2 ln|R(k,Θ)|, R(k,Θ) = E{e(k)eT (k)}, (10) where e(k) = y(k) − ŷ(k, Θ̂) is the innovation error between measurements and model outputs. In the expression for J(Θ), it is assumed that the unknowns are the aircraft aerodynamic parameters. Elements of the R matrix and the initial conditions could also be included among the unknowns, but it is preferred to estimate them a priori as part of the data compatibility checking in order to minimize the number of unknown parameters [24]. If the measurement noise is assumed to be Gaussian with zero mean, then E{e(k)} = 0, E{e(k)eT (l)} = Rδkl (11) The unknown R can be estimated by minimizing the likelihood function with respect to R: R̂ = 1 N N∑ k=1 e(k)eT (k). (12) After substituting the estimated R, the cost function for the output error is obtained in the following way: J(Θ) = N∑ k=1 eT (k)R−1e(k). (13) Using the approximation for the second-order gradient of the cost function for the i-th iteration, the estimate vector Θ̂i+1 is obtained from Θ̂i+1 = Θ̂i + ∆Θ̂i+1, ∆Θ̂ = M−1 ∑ k = 1NHT (k)R̂−1e(k), (14) where H(k) = ∂y(k) ∂Θ(k) , Hij(k) = ∂yi(k) ∂Θj(k) , i = 1, 2, . . . ,m, j = 1, 2, . . . , p, (15) 81 R. Mohammadi Farhadi, V. Kortunov, A. Molchanov, T. Solianyk Acta Polytechnica Figure 5. Block diagram for the output-error param- eter estimator. where m is the number of output signals and p is the number of unknown parameters, and matrix M is the approximation of the Fisher information matrix: M = ∑ k = 1NHT (k)R̂−1H(k). (16) In order to have stable estimates with respect to a small deviation in the measurement data, it is nec- essary for the sensitivity coefficients ∂Θ̂ ∂y to be nil or very small. Therefore, one criterion for the stability of the estimations can be written as [6]: J = N∑ k=1 (∂Θ̂ ∂y )−1 ∂Θ̂ ∂y . (17) The elements of the sensitivity matrix H can be com- puted by solving sensitivity equations d dt ∂x ∂Θj = ∂f ∂x ∂x ∂Θj + ∂f ∂Θj , ∂y ∂Θj = ∂h ∂x ∂x ∂Θj + ∂h ∂Θj . (18) Another way is to compute the sensitivities by a numerical method. The simplest one is to use a finite difference approximation. The accuracy of the param- eters is given by the inverse of the information matrix. As it was pointed out, the diagonal elements of the matrix M form the Cramer-Rao lower bound of the parameter variance. Cramer-Rao bound sometimes does not accurately reflect the true parameter variance. Incorrect assump- tions about the measurement and process noise, mod- elling errors and non-linearity can cause a the lower bound and the actual parameter variance to differ [7]. This shows the importance of the accuracy and ade- quacy of the mathematical dynamic model. The maximum likelihood parameter estimation method can be simplified when applied to a determin- istic linear dynamic system. In this case, there is no process noise. The block diagram for the output-error Figure 6. Block diagram for the output-error param- eter estimator with wind estimation. Figure 7. Identification of the UAV and actuator dynamic model. parameter estimator is shown in Figure 5. Close-to- constant wind disturbance and high frequency pro- peller noise practically acts on the UAV in the flight. It is possible to estimate the amplitude and direction of the near-to-constant wind in the inertial frame with the output error estimation algorithm (Figure 6). It is worth noting that in this case, the actuator angular deviation is not measured. Therefore, the dynamic model of the UAV along with the actuator are iden- tified together. Due to the rapidity of the actuator dynamics relative to the UAV, this approximation is acceptable (Figure 7). The mathematical model in the state space, ẋ(t) = A(Θ)x(t) +B(Θ)u(t), y(t) = C(Θ)x(t) +D(Θ)u(t), (19) is identifiable if and only if the rank of the Jacobian 82 vol. 58 no. 2/2018 Estimation of Lateral Aerodynamic Coefficients for Skywalker X8 ∂Q(Θ) ∂Θ equals p [6], where: QT (Θ) =  DT (Θ) [C(Θ)B(Θ)]T [C(Θ)A(Θ)B(Θ)]T ... [C(Θ)A2n−1(Θ)B(Θ)]T  , (20) where n is the length of the state vector and p is the length of the parameter vector Θ. Based on (19) and by using the OEM, all aero- dynamic coefficients of the lateral dynamics for the FW are identifiable and the identifiability matrix is a 50× 12 full rank matrix. It is possible to determine the quantity of the identifiability for these parameters based on the analysis of the singular values of the identifiability matrix. Singular values of the identifi- ability matrix for the lateral dynamics of the ZFW have a very wide range from 0.14 to 1.7 · 1010. The parameters, which correspond to the small singular values of the identifiability matrix, are called poorly identifiable parameters: Clp ↔ 1.7 · 1010, Clβ ↔ 48.4, Cnp ↔ 7.7 · 109, CnδA ↔ 16.9, CYp ↔ 3.8 · 108, Cnr ↔ 12.9, ClδA ↔ 1.6 · 107, CYδA ↔ 1.06, Clnβ ↔ 3.1 · 104, Clr ↔ 0.54, ClYβ ↔ 1.5 · 103, CYr ↔ 0.14. (21) The set of singular values is divided in two subsets of small and large values. It can be said that the parame- ters Clβ , CnδA , Cnr , CYδA , Clr and CYr are the poorly identifiable ones. Therefore, it is difficult to guarantee the convergence of the parameter identification algo- rithm even in the presence of the persistent exciting input. However, there are usually co-linearity and poorly exciting inputs in the closed loop flight tests that further complicate the identification algorithm. 4. Co-linearity analysis for lateral acceleration regression In order to demonstrate the detection of a co-linearity and the biased estimation techniques to solve it, the lateral data of the SX8FW were used. Part of the UAV lateral response is illustrated in Figure 8. Time histories of the control surface δA, and four response variables v−vw, p, r and −aY after removing their mean values are shown. They show that aY depends mainly on δA and r terms, whereas the contribution of the remaining parameters is very small. In addition to the co-linearity, a low sensitivity of several terms in the aerodynamic model equation could also aggravate the estimation procedure and the accuracy of the estimates. 1000 1100 1200 1300 1400 1500 -5 0 5 v- v W [ m /s ec .] 1000 1100 1200 1300 1400 1500 -0.5 0 0.5 p [r ad /s ec .] 1000 1100 1200 1300 1400 1500 1600 -0.2 0 0.2 0.4 r [r ad /s ec .] 1000 1100 1200 1300 1400 1500 -0.05 0 0.05  A [ ra d] 1000 1100 1200 1300 1400 1500 -0.5 0 0.5 -a Y [ m /s 2 ] Time [sec.] Figure 8. Time histories of lateral response variables. Because of small changes in the input and output variables, the equation for the lateral force is formu- lated as ay = Y0 + Yv(v − vw) + Ypp+ Yrr + YδAδA, Yv = %S 2m √ u2 0 + w2 0CYβ , Yp = %Va0Sb 4m CYp , Yr = %Va0Sb 4m CYr , YδA = %V 2 a0 S 2m CYδA (22) where u0 and w0 are longitudinal and vertical veloc- ities in the body frame for trim conditions. Aero- dynamic parameters Y0, Yv, Yp, Yr and YδA can be estimated by using the LSM instead of the OEM. The singular values of the information matrix and the correlation matrix are investigated to detect and assist the co-linearity. The singular values, condi- tion indexes and correlation matrix of the measured variables (without normalization) are presented in Table 2. The singular values vary between 3.22 and 604.95 and there is one condition index greater than 50 that belongs to CYδA . But from the correlation matrix, it is possible to obtain three linear depen- dencies involving CYβ , CYp and CYr because elements of the correlation matrix for these parameters are greater than 0.7. Without a normalization, the con- dition number, with respect to the inversion for the information matrix (the ratio of the largest singular value of the matrix to the smallest), equals to 187.8. Large condition numbers indicate a nearly singular matrix [5, 7, 25]. For further discussion and analysis, it will be more convenient and necessary to deal with regressor vari- ables, which are centred and scaled to a unit length. Therefore, it is suitable to calculate singular values, condition indexes and a correlation matrix for normal- ized measured variables (Table 3). It is clear that the condition number of the matrix for normalized data is reduced from 187.8 to 9.05 and, perhaps, there is a strongly linear dependency only between CYp and CYδA . In most cases, Yp and Yr are assumed to be zero. Therefore, they are deleted systematically from the re- 83 R. Mohammadi Farhadi, V. Kortunov, A. Molchanov, T. Solianyk Acta Polytechnica Singular value σj Condition index σmax σj Correlation Matrix Y0 CYβ CYp CYr CYδA Y0 244.87 2.47 1 -0.2511 0.1055 -0.3117 -0.5426 CYβ 604.95 1 -0.2511 1 -0.9217 0.972 0.5694 CYp 21.98 27.53 0.1055 -0.9217 1 -0.9425 -0.213 CYr 23.89 25.32 -0.3117 0.972 -0.9425 1 0.4613 CYδA 3.22 187.8 -0.5426 0.5694 -0.213 0.4613 1 Table 2. Detection and assessment of data co-linearity (without normalization). Singular value σj Condition index σmax σj Correlation Matrix Y0 CYβ CYp CYr CYδA Y0 62 3.95 1 -0.2778 -0.3073 -0.2967 -0.3862 CYβ 244.87 1 -0.2778 1 -0.5366 0.1108 -0.3966 CYp 46 5.32 -0.3073 -0.5366 1 -0.5821 0.7492 CYr 52.54 4.66 -0.2967 0.1108 -0.5821 1 -0.4039 CYδA 27 9.05 -0.3862 -0.3966 0.7492 -0.4039 1 Table 3. Detection and assessment of data co-linearity (with normalization). gression of the lateral acceleration and, again, singular values, condition indexes and correlation matrix are calculated for remained parameters. In Table 4, aero- dynamic coefficients are estimated for the following six cases [5, 7, 25, 26]: (1.) Principle components estimates (PCE): Using only one parameter in the regression (2.) PCE: deleting the aerodynamic parameters CYp and CYr from the regression (3.) PCE: deleting the aerodynamic parameter CYp from the regression (4.) Ordinary LSM (OLSM): all of the aerodynamic parameters are used in the regression (5.) Mixed estimates (ME): Using a priori values for CYp and CYr (6.) PCE: deleting the aerodynamic parameter CYδA from the regression Based on the condition index analysis, the relative quality of estimation for these aerodynamic parame- ters are as follows: CYβ ↔ 1, Y0 ↔ 2, CYr ↔ 3, CYp ↔ 4, CYδA ↔ 5. (23) Equations (21) and (23) show that the relative quality of estimation for ,CYβ , CYp and CYr is quite different with the OEM and LSM based on the condition index analysis. On the one hand, according to the LSM in Case 1, there is no linear dependency, so estimated parameters, in this case, can be interpreted as values close to reality. On the other hand, if the estimated parameters in the other cases are close or equal to these values, they can be accepted as good estimates. In the fourth case, with the ordinary LSM, the es- timated parameters CYβ , CYp and CYr are different from the first case with 57% , 20% and 7.7% respec- tively. Parameter CYp is completely estimated with a different sign in these cases. In the third case, with removing CYp from the fourth regression, estimated parameters CYβ , CYr and CYδA are different from the first case with 49%, 14% and 1% respectively. And in the second case, with removing CYp and CYr from the fourth regression, estimated parameters CYβ and CYδA are different from the first case with 15% and 0.2% respectively. In the fifth case, the parameters CYp and CYr are fixed with the values from separate regressions and are estimated from the other param- eters using the mixed estimate method. Finally, the parameter CYβ is removed from the regression and the other parameters in the sixth case are estimated. CYδA has the worst condition index, but there is a strong correlation between the lateral acceleration and this parameter. Therefore, the fit error is increased more than twice with removing CYδA from the regres- sion. The parameter estimates for the lateral-force coef- ficient from an ordinary linear regression, principal components regression and mixed estimation are sum- marized in Table 4. These estimates are compared with the theoretically calculated aerodynamic param- eters of the ZFW and SX8FW in the first and second rows of Table 4. The use of the ordinary LSM can result in non-physical values for parameter CYp and a small value for CYβ . The principal components estimates of these parameters are close to the sepa- rately estimated parameters for CYβ and CYδA . The mixed estimation with a priori values for CYp and CYr also improves the parameter values. The esti- mate of CYδA has approximately the same value for all three techniques, which is the result of the very high sensitivity of this parameter despite of its small 84 vol. 58 no. 2/2018 Estimation of Lateral Aerodynamic Coefficients for Skywalker X8 Y0 (3σ) CYβ (3σ) CYp (3σ) CYr (3σ) CYδA (3σ) Criteria SX8FW 0 −0.1949 −0.1172 0.0959 −0.0696 ZFW 0 −0.07359 0 0 0 Case 1 PCE 0 −0.0033 (0.0002) −0.2423 (0.0046) −0.1914 (0.0045) −0.2299 (0.0017) Case 2 PCE 0 (0.0008) −0.0028 (0.0002) 0 0 −0.2294 (0.0017) 0.0661 Case 3 PCE 0 (0.0008) −0.0017 (0.0002) 0 −0.1645 (0.0046) −0.2275 (0.0017) 0.0592 Case 4 OLSM 0 (0.0008) −0.0014 (0.0002) +0.1004 (0.0056) −0.1547 (0.0046) −0.2476 (0.0020) 0.0574 Case 5 ME 0 (0.0008) −0.0024 (0.0002) −0.2423 −0.1914 −0.1795 (0.0017) 0.0766 Case 6 PCE 0 (0.0008) −0.0028 (0.0002) −0.2676 (0.0047) −0.2051 (0.0046) 0 0.1202 Table 4. Estimated aerodynamic coefficients and their 3σ deviations. 1050 1100 1150 1200 1250 1300 1350 1400 -0.6 -0.4 -0.2 0 0.2 0.4 Time [sec.] aY [ m /s 2 ] aY Yv*(v-vW)+Yda*da Figure 9. Part of the time history match of the lateral acceleration and its regressor prediction. singular value. However, there is a high correlation be- tween the lateral acceleration and the aileron control signal. As indicated by the last column in Table 4, the root mean square fit error, given as the estimation error of the lateral acceleration, is not significantly degraded by the principal components estimates and mixed estimation methods compared with the ordinary LSM. It is concluded that it is possible to reduce the co- linearity by removing the dependent parameters. In order to obtain some feel for the importance of parame- ters in (22), the estimated parameters in this equation for the SX8FW are compared with the data for the aerodynamic parameters of the ZFW and SX8FW in Table 4. Parameters estimation errors are three times larger than the standard deviation (3σ) shown in Ta- ble 4. It is suggested, with ignoring CYp and CYr , to use the results of Case 2 as good and accepted results for the next uses. It can be noticed that estimated values for the parameters CYβ and CYδAare different from their theoretically calculated values. The time history match of the lateral acceleration and prediction of its regressor, Yv(v − vw) + YδAδA is shown in Figure 9. 5. Estimation of aerodynamic coefficients of the lateral dynamics Theoretical and Numerical values for the aerodynamic parameters of the ZFW and SX8FW were given in [18, 19].These values are used as initial values for the parameters in the optimization algorithms. With an expert, it is possible to determine the accuracy and range of the uncertainty of the theoretical value for each parameter. Aerodynamic parameters Clβ and Cnβ , which determine the static stability of the UAV, are more accurate than the other parameters can be theoretically calculated. However, there is a strong correlation between the roll angular rate p and aileron control signal δA that makes the co-linearity problem. Therefore, it is preferable to cancel the co-linearity problem because it negatively affects the estimates of other parameters. These problems must be solved: • To find physically acceptable estimates for param- eters (a constrained optimization based on the ac- curacy of theoretical values of the aerodynamic parameters from expert). • To solve the co-linearity problem by removing some dependent variables with the principle component estimate and mixed estimate methods in the con- strained optimization algorithm. 85 R. Mohammadi Farhadi, V. Kortunov, A. Molchanov, T. Solianyk Acta Polytechnica Parameters OEM-IAV (3σ) (SX8FW) LSM-OEM-AC (3σ) Theoretical values ZFW SX8FW CYβ −0.0972 (11%) −0.0008 (0.0029) −0.0736 −0.1949 CYδA −0.0160 (23%) −0.2477 (0.0070) 0 −0.0696 Clβ −0.0993 (2%) −0.0068 (0.0094) −0.0285 −0.0765 Clp −0.4669 (21%) −0.4715 (0.0386) −0.3209 −0.4018 Clr 0.1650 (57%) 0.0250 (0.0038) 0.3066 0.025 ClδA 0.1410 (2%) 0.0355 (0.0022) 0.1682 0.2987 Cnβ 0.0583 (2%) 0.0020 (0.0034) 0.0004 0.0403 Cnp −0.1721 (56%) −0.0246 (0.0246) −0.01297 −0.0247 Cnr −0.1216 (4%) −0.0228 (0.0018) −0.0043 −0.1252 CnδA −0.0635 (6%) −0.0013 (0.0139) −0.0032 −0.0076 Criteria 117.86 116.6849 σv 2.4826 2.5438 σp 0.0869 0.0918 σr 0.0838 0.0638 σϕ 0.106 0.1054 σaY 0.0597 0.0592 Table 5. Aerodynamic parameters estimates and their deviations for SX8FW. These approaches are used to apply the optimization methods for estimating the lateral aerodynamic pa- rameters: • OEM with only Modified Newton-Raphson or Gauss-Newton method with Initial Analytical Val- ues from SX8FW (OEM-IAV): Analytical values of the aerodynamic parameters are used as initial conditions and the Gauss-Newton method is used to minimize the cost function for the unconstrained optimization. There are two shortcomings of this approach: 1. It is not able to find the global opti- mum with the Gauss-Newton method. 2. It leads to a large and unacceptable linear and angular ac- celerations. This means that absolute values for the time history of the right side of the dynamic equa- tions are larger than the ones of the left side (for example, absolute values of Lpp and Nrr greater than the absolute values of p′ and r′). In this case, covariances for all estimated parameters are lower than 20%. • The LSM for the lateral acceleration and the OEM for lateral dynamic equations of motion with Ana- lytic Constraints (LSM-OEM-AC): Analytical val- ues of the aerodynamic parameters and their a priori analytical uncertainties based on an expert evaluation are used for the SAOA in constrained optimization. For removing the co-linearity, Clr and CnδA are assigned zero and for the order reduc- tion, estimated values for CYβ and CYδA in lateral acceleration regression are used in the constrained optimization stage. Then Gauss-Newton method is used to finalize the optimum values and covariances of the estimated parameters. • The OEM with the SAOA and Gauss-Newton Methods with Measured Constraints (OEM-MC): Putting upper or lower limits on the aerodynamic parameters based on the linear and angular acceler- ation for the constraint optimization. It is possible to compute the mean square deviation for the vari- ables v′, p′, r′ from the measurements v, p and r (it means cov(Lpp) ≤ cov(p′) and cov(Nrr) ≤ cov(r′)). The range of the changes of linear and angular accel- 86 vol. 58 no. 2/2018 Estimation of Lateral Aerodynamic Coefficients for Skywalker X8 Parameters OEM-MC (3σ) LMS-OEM-MC (3σ) OEM-IAV (3σ) (ZFW) Theoretical values ZFW SX8FW CYβ −0.0334 (0.0033) −0.0022 (0.0035) −0.0736 (1%) −0.0736 −0.1949 CYδA −0.2089 (0.0100) −0.2432 (0.0066) −0.0 (0.007) 0 −0.0696 Clβ −0.0084 (0.0025) −0.0054 (0.0081) −0.0526 (1%) −0.0285 −0.0765 Clp −0.2376 (0.0304) −0.1586 (0.0167) −0.3209 (1%) −0.3209 −0.4018 Clr 0.0266 (0.0093) 0.0313 (0.0089) 0.1658 (1%) 0.3066 0.025 ClδA 0.0261 (0.0053) 0.0178 (0.0035) 0.1682 (1%) 0.1682 0.2987 Cnβ 0.0037 (0.0004) 0.0028 (0.0120) 0.0006 (33%) 0.0004 0.0403 Cnp −0.0881 (0.0363) −0.0899 (0.0191) −0.0131 (5%) −0.01297 −0.0247 Cnr −0.0301 (0.0024) −0.0329 (0.0044) −0.0041 (80%) −0.0043 −0.1252 CnδA −0.0009 (0.0012) 0.0004 (0.0092) −0.0000 (0.0005) −0.0032 −0.0076 Criteria 117.447 117.1416 126.98 σv 2.4837 2.4759 2.6556 σp 0.0921 0.0918 0.0924 σr 0.0651 0.065 0.0819 σϕ 0.1055 0.1052 0.1095 σaY 0.0755 0.0756 0.0769 Table 6. Aerodynamic parameters estimates and their deviations for SX8FW. erations is an effective rule to validate the estimated parameters. It is possible to find an upper limit for absolute values of the parameters. Then, the SAOA for constrained optimization is used with bounds of this approach. And finally, Gauss-Newton method is applied to find final optimum parameters and their covariances. • The LSM for the lateral acceleration and the OEM for the lateral dynamic equations of motion with Measured Constraints (LMS+OEM+MC): In this case, early estimated values for Crβ and CrδA by the LSM from the lateral acceleration are used in the OEM at the second stage. In this case, for the optimization, at first, the parameters CnδA , Cnp and Clp are also fixed to remove the data co-linearity. Then, the parameters Cnp and Clp are made free in the optimization algorithm. And finally, Gauss- Newton method is applied to find final optimum parameters and their covariances. It can be seen that in this case, covariances of the parameters CnδA , CYβ , Clβ and Cnβ are very large and they are not estimated accurately. • The OEMwith the only Gauss-Newton method with Initial Analytical Values from the ZFW (OEM-IAV). From (21), for weakly identifiable parameters, CYδA and CnδA are put to zero. These parameters have small contributions on outputs and it is possible to use zero value for them. The weakly identifiable parameters Clv and Clr are set to values between analytical ones for the ZFW and SX8FW. IAV-OEM is used with these conditions and initial conditions based on analytical values for other parameters of the ZFW. In Tables 5 and 6, estimation parameters with differ- ent suggested approaches are given. Output signals for mathematical models with estimated parameters of all of the approaches are given in Figures 10–14. Identification quality indicators and properties of es- timated aerodynamic parameters with five different approaches are discussed below: • The largest value for the cost function is obtained 87 R. Mohammadi Farhadi, V. Kortunov, A. Molchanov, T. Solianyk Acta Polytechnica 900 1000 1100 1200 1300 1400 1500 1600 1700 -20 0 20 v [m /s ec ] Flight-Test Data Mathematical Model 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1 p [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -0.5 0 0.5 r [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -0.5 0 0.5  [r ad ] 900 1000 1100 1200 1300 1400 1500 1600 1700 -2 0 2 a Y [ m /s 2 ] Time [sec.] Figure 10. Flight data and responses of the model with the OEM-IAV parameters with initial condition from analytical aerodynamic parameters of SX8FW. 900 1000 1100 1200 1300 1400 1500 1600 1700 -20 0 20 v [m /s ec ] Flight-Test Data Mathematical Model 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1 p [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -0.5 0 0.5 r [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1  [r ad ] 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1 a Y [ m /s 2 ] Time [sec.] Figure 11. Flight data and responses of the model with the LSM-OEM-AC parameters. with the OEM-IAV (ZFW) approach. For this ap- proach, the time history match, especially for the yaw angular rate and roll angle, is reduced. The esti- mated aerodynamic coefficients with this approach are the most similar to the ZFW analytical ones. • There is a tendency to increase the parameter Clp in the optimization process due to the correlation between the roll angular rate, p(t) and the aileron control signals, δA(t). This tendency can be seen es- pecially in the LSM-OEM-AC approach. Therefore, this tendency can be controlled with the constrained optimization on the basis of the dispersion of the measured linear and approximately calculated angu- lar accelerations about their mean values (OEM-MC and LMS+OEM+MC). However, theoretical val- ues of the aerodynamic parameters, especially for Clp and Cnr , are large. Although the lowest value for the cost function belongs to this case (LSM- OEM-AC), the estimated values for the coefficients separately are not consistent with the bounds of the measured-signals-based LSM for angular and 900 1000 1100 1200 1300 1400 1500 1600 1700 -20 0 20 v [m /s ec ] Flight-Test Data Mathematical Model 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1 p [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -0.5 0 0.5 r [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1  [r ad ] 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1 a Y [ m /s 2 ] Time [sec.] Figure 12. Flight data and responses of the model with the OEM-MC parameters. 900 1000 1100 1200 1300 1400 1500 1600 1700 -20 0 20 v [m /s ec ] Flight-Test Data Mathematical Model 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1 p [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -0.5 0 0.5 r [r ad /s ec .] 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1  [r ad ] 900 1000 1100 1200 1300 1400 1500 1600 1700 -1 0 1 a Y [ m /s 2 ] Time [sec.] Figure 13. Flight data and responses of the model with the LSM-OEM-MC parameters. linear accelerations. Results of the LSM-OEM-AC and OEM-IAV cases are shown in Table 5. • Due to the slight change of the amplitude and di- rection of the wind in the inertial frame and distur- bances due to the actuator (because the actuator worked without the position feedback), there is not a complete history matching between mathematical models’ outputs to control signals and real flight- test measurements. • A long time history of flight-test data can be used for an estimation of the parameters, which produce the DC gains of the different transfer functions. • It must be noticed that due to the lack of the iden- tifiability and the linear dependency between some flight variables, it is not possible to estimate all of the aerodynamic parameters. Depending on the degree of identifiability, the linear dependency be- tween the variables is eliminated in the constrained optimization step and then the other aerodynamic parameters are estimated. • It is suggested to use analytical values for some 88 vol. 58 no. 2/2018 Estimation of Lateral Aerodynamic Coefficients for Skywalker X8 0 200 400 600 800 1000 1200 1400 1600 1800 -20 0 20 v [m /s ec ] Flight-Test Data Mathematical Model 0 200 400 600 800 1000 1200 1400 1600 1800 -1 0 1 p [r ad /s ec .] 0 200 400 600 800 1000 1200 1400 1600 1800 -0.5 0 0.5 r [r ad /s ec .] 0 200 400 600 800 1000 1200 1400 1600 1800 -1 0 1  [r ad ] 0 200 400 600 800 1000 1200 1400 1600 1800 -2 0 2 a Y [ m /s 2 ] Time [sec.] Figure 14. Flight data and responses of the model with the OEM-IAV parameters with initial condi- tion from analytical aerodynamic parameters of the ZFW. Figure 15. Typical model prediction and fitting properties as a function of its number of parameters. parameters, which make the DC gains, to remove the co-linearity problem and estimate the remaining parameters. Because, theoretically, these static stability derivatives are calculated more accurately than the damping and primary control derivatives. However, there are significant differences between estimated and theoretical ones regarding the LSM and OEM. • It is visible that the roll rate is very noisy due to vibrations, especially from the propeller thrust. • It is worth mentioning that the linear model predicts the dynamic behaviour of the SX8FW acceptably. In order to determine the structure and the number of parameters for the mathematical model, it must be considered that usually fitting error is decreased monotonically with the number of the parameters, but prediction quality of the model is not increased monotonically (Figure 15) [27]. • For the purposes of the prediction, the co-linearity 0 200 400 600 800 1000 1200 1400 -20 0 20 v [ m /s e c ] Flight-Test Data Mathematical Model 0 200 400 600 800 1000 1200 1400 -2 0 2 p [ ra d /s e c .] 0 200 400 600 800 1000 1200 1400 -1 0 1 r [r a d /s e c .] 0 200 400 600 800 1000 1200 1400 -2 0 2  [ ra d ] 0 200 400 600 800 1000 1200 1400 -2 0 2 a Y [ m /s 2 ] Time [sec.] Figure 16. Prediction capability of the model with estimated parameters for the OEM-AIV (with initial condition from analytical aerodynamic parameters of SX8FW case) and for additional flight-test data. does not cause a serious problem. It must be men- tioned that mathematical models with different groups of estimated parameters have a good predic- tion capability. This shows the accuracy and ade- quacy of the mathematical dynamic model. How- ever, there are differences between estimated pa- rameters and analytical ones for the LSM-OEM-AC and LSM-OEM-MC and covariances of estimated parameters are large in comparison with the OEM- IAV approaches. However, it is possible to achieve physically meaningful estimations and low estima- tion covariances by allowing a tolerance to some increase of the cost function. • Physical meaning, small estimation covariances and the dynamic model prediction capability are used for the model validation as effective indicators. In Figure 16, the accepted prediction capability of the model using the estimated parameters with the OEM-IAV case (initial condition from analytical aerodynamic parameters of SX8FW) is shown. • The OEM is a nonlinear optimization method. It is not possible to prove the achievement of the global minimum for it especially when there are many pa- rameters. Therefore, it is helpful to use the LSM in the first stage for the lateral acceleration to decrease the number of parameters in the optimization stage (use of the OEM). However, LSM did not give esti- mates close to analytical parameters due to noisy measurements. • In the OEM-IAV approach (with the initial con- dition from analytical aerodynamic parameters of SX8FW) estimated values of CYβ , CYδA and CnδA are close to the ZFW aerodynamics and there are significant differences between estimated val- ues from the LSM for the lateral acceleration. The OEM gave better results than the LSM, because it uses the multi-dimensional equations of motion 89 R. Mohammadi Farhadi, V. Kortunov, A. Molchanov, T. Solianyk Acta Polytechnica and more effectively estimates parameters in the presence of the measurements noise. Aerodynamic coefficients, Clβ , Clp , Cnβ and Cnr are estimated close to the SX8FW analytical ones. Estimated values for Cnp and CnδA are not matched with ana- lytical ones for the SX8FW and ZFW. Estimated parameter, Clr has the value between the SX8FW and ZFW one. • In the OEM-IAV approach (with the initial con- dition from analytical aerodynamic parameters of the ZFW), estimated values except for Clβ and Clr are very close to the ZFW aerodynamic parameters with small covariances. For this approach, there are also significant differences between estimated values from the LSM for the lateral acceleration .Based on a singular value analysis, it is expected that parameters CYδA , Clβ , CnδA and Clr are badly estimated. • The LSM-OEM-AC approach could find the small- est value for the cost function, because for this approach, large bound limits are used for the con- straint optimization. As it is previously stated, for the purposes of predic- tion, the co-linearity does not cause a serious problem. However, for the purposes of the estimation of the parameters, biased estimation methods can be used. The increase in the number of flights alone does not lead to a significant increase in accuracy, since flights are similar in manoeuvre and, accordingly, in the na- ture of the input signal. The number of measurements in the experiments depends on the frequency of the data recording (20Hz) and flight time (10–40 min- utes). According to statistical properties, the amount of data is sufficient for an identification, but the spec- tral properties of the input signals depend little on the flight time. 6. Conclusion The Aerodynamic coefficients of the lateral dynamics for the SX8FW were estimated from flight-test data with a co-linearity problem based on the combination of constrained optimization, LSM, and OEM. And several approaches, such as parameter fixing, princi- ple components estimation and a mixed estimation with the LSM and the OEM in constrained optimiza- tion, were applied to solve the co-linearity problem for estimating aerodynamic parameters. The data co- linearity leads to a not repeated and unstable estima- tion of the parameters with large covariances. There- fore, it is necessary to accept some biased estimations and remove the correlation between some variables by fixing or omitting the corresponding parameters to solve the co-linearity problem and estimate the other aerodynamic parameters correctly. As a quality indicator of the estimation to validate the estimated parameters, the physical meanings for aerodynamics parameters, theoretical values of the parameters for the FW being studied and a similar one, small covari- ances for the estimates, repeated and stable estimates for parameters, and prediction capabilities of the esti- mated dynamic model for other flight-test data were used. With additional and persistent flight-test data; it is possible to remove the data co-linearity and near- singularity problems. In addition, the time-varying model for the wind in an inertial frame can be used to improve the model adequacy. The approaches consid- ered for the SX8FW have shown satisfactory results and in the future, we propose to use them for other types of aircrafts, including vertical take off and land- ing vehicles. The work will continue in the direction of developing a real-time identification and control system tuning methods for such UAVs and designing a safe and optimal persistent exciting manoeuvre to remove the co-linearity problem in the UAV system identification. The developed models will later be used to synthesize a robust control and testing in real flights. The modification of the criterion function in the identification process based on the a priori infor- mation is suggested to solve the linear dependency problem. References [1] R. Isermann, “Practical aspects of process identification,” Automatica, vol. 16, no. 5, pp. 575-587, 1980. doi:10.1016/0005-1098(80)90079-5 [2] M. B. Tischler, “System identification methods for aircraft flight control development and validation,” Advances in Aircraft Flight Control, pp. 35-69, 1996. [3] A. A. Tunik, and A. N. Klipa, “Identification of mathematical model of longitudinal motion of the aircraft in the presence of noise measurement and displacement sensors biases (in Russian: Identifikacija matematicheskoj modeli prodol’nogo dvizhenija letatel’nogo apparata pri nalichii shumov izmerenij i smeshhenij datchikov),” Electronic modeling (Jelektronnoe modelirovanie), vol. 6, pp. 3-18, 2001 2001. [4] E. A. Morelli, “System identification programs for aircraft (SIDPAC),” 2002. doi:10.2514/6.2002-4704 [5] J. R. Raol, G. Girija, and J. Singh, Modelling and parameter estimation of dynamic systems. Iet, 2004. [6] V. A. Kas’yanov, Flight simulation (in Russian: Modelirovanie poleta). Moscow: Moscow, 2004, p. 400. [7] V. Klein, and E. A. Morelli, Aircraft system identification: theory and practice. American Institute of Aeronautics and Astronautics Reston, Va, USA, 2006. doi:10.2514/4.861505 [8] B. Mettler, Identification modeling and characteristics of miniature rotorcraft. Springer Science & Business Media, 2013. [9] R. M. Farhadi, V. I. Kortunov, and A. Mohammadi, “UAV motion model and estimation of its uncertainties with flight-test data,” in 22rd Saint Petersburg International Conference on Integrated Navigation Systems, 2015, pp. 131-133. 90 http://dx.doi.org/10.1016/0005-1098(80)90079-5 http://dx.doi.org/10.2514/6.2002-4704 http://dx.doi.org/10.2514/4.861505 vol. 58 no. 2/2018 Estimation of Lateral Aerodynamic Coefficients for Skywalker X8 [10] A. S. Holtsov, R. M. Farhadi, V. I. Kortunov, and A. Mohammadi, “Comparison of the UAV adaptive control with the robust control based on mu-synthesis,” in Methods and Systems of Navigation and Motion Control (MSNMC), 2016 4th International Conference on, 2016, pp. 18-21: IEEE. doi:10.1109/MSNMC.2016.7783096 [11] R. M. Farhadi, and V. I. Kortunov, “Identification of roll channel dynamics for the unmanned aerial vehicle under weakly excited input signal (in Russian: Identifikaciya dinamiki kanala krena bespilotnogo letatel’nogo apparata pri slabo informativnom vxodnom signale),” Radio-electronic and Computer Systems (Radioelektronni i komp’yuterni sistemi), vol. 1, no. 1, pp. 99-106, 2017. [12] R. M. Farhadi, S. K. Yadavar Nikravesh, and A. Mohammadi, “Aerodynamic Model Identification Using Flight Test Data with EKF and NN,” in Proceedings of the 4th Iranian Aerospace Society Conference, Tehran, Iran, 2003, vol. Papers Presented in English, pp. 280-286: Iranian Aerospace Society. [13] M. Verhaegen, and V. Verdult, Filtering and system identification: a least squares approach. Cambridge university press, 2007. doi:10.1017/CBO9780511618888 [14] J. L. Crassidis, and J. L. Junkins, Optimal estimation of dynamic systems. CRC press, 2011. [15] N. Von Hoffer, System identification of a small low-cost unmanned aerial vehicle using flight data from low-cost sensors. Utah State University, 2015. [16] T. Söderström, and P. Stoica, System identification. Prentice-Hall, Inc., 1988. [17] P. J. Van Laarhoven, and E. H. Aarts, “Simulated annealing,” in Simulated annealing: Theory and applications: Springer, 1987, pp. 7-15. doi:10.1007/978-94-015-7744-1_2 [18] R. W. Beard, and T. W. McLain, Small unmanned aircraft: Theory and practice. Princeton university press, 2012. doi:10.1515/9781400840601 [19] K. Gryte, “High Angle of Attack Landing of an Unmanned Aerial Vehicle,” NTNU, 2015. [20] R. Austin, Unmanned aircraft systems: UAVS design, development and deployment. John Wiley & Sons, 2011. [21] V. I. Kortunov, I. Y. Dybska, and G. A. Proskura. (2009, 03 February 2009) Integrated mini INS based on MEMS sensors for UAV control. IEEE Aerospace and Electronic Systems Magazine. 41-43. doi:10.1109/MAES.2009.4772754 [22] V. I. Kortunov, O. V. Mazurenko, A. V. Gorbenko, W. Mohammed, and A. Hussein, “Review and comparative analysis of mini- and micro-UAV autopilots,” in Actual Problems of Unmanned Aerial Vehicles Developments (APUAVD), 2015 IEEE International Conference, Kiev, Ukraine, 2015, pp. 284-289: IEEE. doi:10.1109/APUAVD.2015.7346622 [23] R. F. Stengel, Flight dynamics. Princeton University Press, 2015. [24] E. Tohme, “Initialization of output error identification algorithms,” Université de Poitiers, 2008. [25] N. R. Draper, and H. Smith, Applied regression analysis. John Wiley & Sons, 2014. [26] J. O. Rawlings, S. G. Pantula, and D. A. Dickey, Applied regression analysis: a research tool. Springer Science & Business Media, 2001. [27] D. E. Stepner, and R. K. Mehra, “Maximum likelihood identification and optimal input design for identifying aircraft stability and control derivatives,” 1973. 91 http://dx.doi.org/10.1109/MSNMC.2016.7783096 http://dx.doi.org/10.1017/CBO9780511618888 http://dx.doi.org/10.1007/978-94-015-7744-1_2 http://dx.doi.org/10.1515/9781400840601 http://dx.doi.org/10.1109/MAES.2009.4772754 http://dx.doi.org/10.1109/APUAVD.2015.7346622 Acta Polytechnica 58(2):77–91, 2018 1 Introduction 2 Problem statement 3 OEM and identifiability problem 4 Co-linearity analysis for lateral acceleration regression 5 Estimation of aerodynamic coefficients of the lateral dynamics 6 Conclusion References