Jtam-A4.dvi JOURNAL OF THEORETICAL AND APPLIED MECHANICS 54, 1, pp. 87-98, Warsaw 2016 DOI: 10.15632/jtam-pl.54.1.87 INCLUSION OF THE D-OPTIMALITY IN MULTISINE MANOEUVRE DESIGN FOR AIRCRAFT PARAMETER ESTIMATION Piotr Lichota Warsaw University of Technology, Institute of Theoretical and Applied Mechanics, Warsaw, Poland e-mail: plichota@meil.pw.edu.pl This paper is concernedwith the designing of simultaneous flight control deflections for air- craft system identification. The elevator, ailerons and rudder are excited with harmonically relatedmultisine signals. The optimal deflections are designedwhen there is no information about the stability and control derivatives and when this information is available. The in- clusion of the system dynamics in the inputs design phase is done with the D-optimality criterion. Both sets of optimal flight surface deflections are used as excitations of a nonli- near aircraftmodel which is identified through themaximum likelihood estimationmethod. Parameters accuracy for those maneuvers (designed with and without a-priori knowledge) is presented and compared. Keywords: inputs design, system identification, flight dynamics 1. Introduction To obtain precise information about aircraft dynamics, aerodynamic stability and control de- rivatives have to be determined. This can be done by various methods like wind tunnel tests (Hoe et al., 2012), Computational Fluid Dynamics – CFD (Mader andMartins, 2011) and Sys- tem Identification – Sys-ID (Jameson and Cooke, 2012; Jategaonkar, 2006; Lichota and Ohme, 2014). Among this group, the gathering of aerodynamic databases from flight test data is the most reliable as it is based on experimentsmade on real objects. On the other hand, asmultiple flight tests are performed in the Sys-ID, it is of high cost and time consuming.CFDmethods are the least expensive in this group, however final results of aerodynamic derivatives (regardless of the object) have to be compared with experiments (Rogowski andMaroński, 2015). The Sys-ID approach consists of four main steps that are: designing and performingmano- euvres, measuring the data, modelling the object and estimating unknown parameters. Those steps are strictly connected, and if the registered data will be inaccurate or there will be not enough information about the aircraft dynamics in the output signals, reliable estimation will not be possible. As modern sensors allow one to obtain very high accuracy of the data andme- asurement techniques are well developed a strong emphasis should be put on designing inputs used for exciting the optimal aircraft response. The research concerning excitations that were to maximize the information content in the measured data was extensively investigated in the seventies and eighties. Those studies showed that sine-sweep (linear or logarithmical) or typicalmulti-step inputs (pulse, doublet, 3-2-1-1) can be used to develop a good mathematical model of an aircraft if only one flight control surface is deflected through a manoeuvre. It has been shown recently that simultaneous multi-step excitations could be used for this purpose as well (Lichota and Ohme, 2014), but this requires a priori knowledge of the aerodynamic derivatives. If the flight controls can be deflected at the same time, the use of harmonically related multisine signals is possible aswell. This approach does not require the initial information about 88 P. Lichota the system and allows one to obtain estimates with the same quality (in terms of parameters accuracy) as the design with simultaneous multi-step inputs. If the inclusion of the a priori knowledge in themultisine inputsdesign phasewould increase the quality of estimation, itwould mean that multisine signals are more adequate for aircraft Sys-ID purposes than simultaneous multi-step excitations. The incorporation of initial knowledge of aerodynamic parameters can be done by introdu- cing theD-optimality criterion instead of theRelative Peak Factor in themultisine optimization phase. In the present study, a linear aircraft model is used for this purpose as it is less com- putationally demanding than the nonlinear representation. The designed set of simultaneous multisine excitations is used as inputs for a nonlinear aircraft model that has been created in Matlab. On the basis of recorded signals the unknown parameters of the object are estimated by applying the maximum likelihood principle. A similar procedure is applied to multisine signals that are designed without a priori know- ledge of the system dynamics. Evaluated ailerons, elevator and rudder deflections are used as excitations for the nonlinear aircraftmodel. The response of the object is recorded and then the unknown parameters of the system are identified. The results are compared with the estimates obtained from themanoeuvre inwhich initial values of the aerodynamic derivatives are available in the inputs design phase. 2. Multisine input signals The multisine input is an excitation that is composed of summed harmonic sinusoids with individual amplitudesAk and phase shifts φk (Morelli, 2012) δ= M∑ k=1 Ak sin(2πfkt+φk) (2.1) where k = 1, . . . ,M stands for the number of the harmonic, fk is the frequency of the k-th component and δ is the input signal. An important feature of the multisine signals is that they can be designed as mutually orthogonal in the time and frequency domain (Morelli, 2003). This means that simultaneous ailerons, elevator and rudder deflections can be independent. In order to achieve this aim, it is required to assign different harmonics to each flight control (e.g. 2f0, 5f0, 8f0,. . .to ailerons, 3f0, 6f0, 9f0,. . . to elevator, 4f0, 7f0, 10f0,. . . to rudder).This assignment provides orthogonality in the frequency domain as distinct spectral lines form the frequency content of each input. In the time domain, the mutual orthogonality is achieved due to orthogonality properties of the sine function (Morelli, 2012). It is practicable to omit the first harmonic in this assignment in order to optimize the cost function effectively. Frequencies of the consecutive components in the multisine signals are evenly spaced and based on the excitation time T : fk = k/T . This also limits the minimum available frequency which must satisfy the condition: fmin ­ 2/T . The maximum available frequency is limited by the frequency range of interest in the investigated case. Multisine input signals have wide-frequency band and the amplitudes of the different har- monicsAk are chosen to achieve desired power spectrum. If there is no need to put an emphasis on specific frequencies, a uniform power spectrum should be used. For simultaneous aileron, elevator and rudder deflections and flat power spectrum, the amplitudesAk related to the j-th control are given by Aj,k = Aj√ Mj (2.2) Inclusion of the D-optimality in multisine manoeuvre design... 89 where Aj is the amplitude of the j-th flight control and Mj is the number of the harmonic components assigned to the j-th flight surface. For each flight control, input energymaximization is done through proper selection of conse- cutiveharmonicsphaseanglesφj,k.Thephaseangles shouldmaximize theexcitation effectiveness without unnecessary increase in each signal value as this can cause that aircraft will go too far from the trim conditions andwill hinder the Sys-IDprocess. For this purpose, theRelative Peak Factor (RPF) can be used as it expresses the input amplitude range divided by the excitation energy. TheRPF is also scaled so that for a single sinusoid it is equal to 1. For the j-th flight surface deflection δj, theRPF is given by RPF(δj)= maxδj −minδj 2 √ 2rms(δj) (2.3) RPF minimization is equivalent to input effectiveness maximization and can be achieved with a simplex optimization algorithm. In our study, the Schroeder phase angles are used as initial values of the phase angles φk for the j-th flight control (Schroeder, 1970) φj,k =φj,k−1+2π(fj,k−1−fj,k)tj,k−1 tj,k =T k∑ l=1 pj,l (2.4) where fj,k is thek-th harmonic frequency assigned to the j-th inputand tk is time epoch.For the flat power spectrumand the j-th flight control, the power of the k-th component is pj,k =1/Mj. As deflections should start and end with zero amplitude, it is necessary to find a constant time offset for the components of each excitation. This is equivalent to sliding the inputs along the time axis until the zero crossing occurs at t=0. 3. D-optimal criterion Thedesign process ofmultisine input signals does not require knowledge of the systemdynamics in terms of stability and control derivatives. However, it is considered that inclusion of some information about aerodynamic parameters in the design phase could increase the quality of the aircraft response in terms of the Sys-ID. This requires introducing estimator based on the minimummean squared error (Kay, 1993) MSE(Θ̂)=E[(Θ̂−Θ)2] = cov(Θ̂)+bias2(Θ,Θ̂) (3.1) whereΘ are the parameters of the model and Θ̂ are their estimates. When the estimator is designed to be unbiased, the minimum mean squared error consists only of the covariance term. This part can be determined on the basis of Cramér-Rao inequality which states that the covariance of the unbiased estimator is at least as high as the inverse of the Fisher InformationMatrixF cov(Θ̂)­F−1 (3.2) TheFisher InformationMatrix which is ameasure of the amount of information that observable variables carry about the unknown system parameters is defined as F =E {[ ∂ lnL(Θ|z) ∂Θ ][ ∂ lnL(Θ|z) ∂Θ ]T} (3.3) where the likelihood function L(Θ|z) is equal to p(z|Θ) – conditional probability that the measurement vector z is observed for the model parameters Θ. 90 P. Lichota If the likelihood function is twice differentiable with respect to themodel parameters, it can be shown that F =−E [ ∂2 lnL(Θ|z) ∂ΘΘT ] (3.4) Multivariate normal distribution is usually chosen to evaluate the probability density function at a certain time point tk. Due to variables, independence the probability can be defined for all time points in the manoeuvre (Jategaonkar, 2006) p(z|Θ) = 1 √ (2π)n|R|N exp ( − 1 2 N∑ k=1 [z(tk)−y(tk)]TR−1[z(tk)−y(tk)] ) (3.5) where n is the number of model outputs y, k=1, . . . ,N is the index of time samples and R is the measurement noise covariance matrix. The described approach with neglecting small terms leads strictly to the Fisher Information Matrix sensitivity form F ≈ N∑ k=1 [ ∂y(tk) ∂Θ ]T R −1 [ ∂y(tk) ∂Θ ] (3.6) where the output signals gradients ∂y(tk)/∂Θ can be obtained through forward difference for- mula. Introducing central differences does not increase the accuracy significantly and raises the computational time, therefore, it is not used in the evaluations. The columns of the Fisher Information Matrix represent contributions of the model para- meters, so if the columns are independent, the determinant reaches its maximum value. On the contrary, if the columns are linearly dependent, the determinant will accept the minimum value. Therefore, the optimality criterion (D-optimality) is expressed by maximization of the Fisher Information Matrix determinant what means that the estimation error ellipsoid volume is minimized. Inclusion of the D-optimality criterion in the multisine design phase can be used in order to introduce the initial information about stability and control derivatives. This task can be achieved by incorporating the linear aircraft model andminimization of the Fisher Information Matrix inverse (Parameter Error CovarianceMatrixP) determinant. The linear representation has been selected because it well describes the object, and the computational time for single evaluation is relatively short. 4. Linear model The aircraft dynamic equations of motion are derived in a vehicle carried coordinate system Oxyz and linearised in accordance with Etkin (1972). The origin of theOxyz reference frame is located at the center of gravity. TheOx axis coincides with the longitudinal axis of the airplane. TheOy axis is normal to the aircraft symmetryplane and is pointing in the direction of the right wing. TheOz axis is oriented downward, so it completes the right-handed coordinate system. TheOxyz system is related to the vehicle carriedOxgygzg reference frame through rotation angles:Φ (roll angle),Θ (pitch angle),Ψ (yawangle)whichareused todescribe the orientation of the object. TheOxgygzg coordinate system remaines parallel to the earth fixed reference frame O1x1y1z1, whose origin is located at an arbitrary point of the Earth with theOx1 axis pointed north,Oy1 axis pointed east and theOz1 axis is directed toward center of the Earth. Relations between the described coordinate systems are shown in Fig. 1, and the transformations of linear Inclusion of the D-optimality in multisine manoeuvre design... 91 Fig. 1. Coordinate systems andmotion variables and angular quantities from the gravitational reference frame Oxgygzg to the body coordinate system Oxyz are given by the matrices ΛV =    1 0 0 0 cosΦ sinΦ 0 −sinΦ cosΦ       cosΘ 0 −sinΘ 0 1 0 sinΘ 0 cosΘ       cosΨ sinΨ 0 −sinΨ cosΨ 0 0 0 1    ΛΩ =    1 0 −sinΘ 0 cosΦ sinΦcosΘ 0 −sinΦ cosΦcosΘ    (4.1) Dynamic equations of motion are derived from Newton’s second law of motion in Oxyz and kinematic relationships, which leads to m(V̇O+Ω×VO)=F IΩ̇+Ω× (IΩ)=MO Φ̇=Λ−1Ω Ω (4.2) where VO = [U V W ] T is the velocity of the origin, Ω= [P Q R]T is the angular velocity, Φ̇ = [Φ Θ Ψ]T describes the aircraft orientation, m stands for mass, I for inertia matrix, F= [X Y Z]T andMO = [L M N] T are forces andmoments acting on the object. The dot symbol is used to denote derivatives with respect to time. When the level flight is the equilibrium state, it is possible to linearize the equations of motion and describe the system as follows u̇=Xuu+Xαα+(Xq −u0α0)q+gθcosΘ0+XδEδE β̇=Yββ+(Yp+α0)p+(Yr−1)r+ g |V0| cosΘ0+YδAδA+YδRδR α̇=Zuu+Zαα+(Zq+1)q− g U0 sinΘ0+ZδEδE ṗ=Lββ+Lpp+Lrr+LδAδA+LδRδR q̇=Muu+Mαα+Mqq+MδEδE ṙ=Nββ+Npp+Nrr+NδAδA+NδRδR φ̇= p+rtanΘ0 θ̇= q (4.3) where small letters are used to describe the perturbations of the flight state from the trim condition which is denoted by subscript 0 (e.g. U =u0+u). Symbol g is used for gravitational 92 P. Lichota acceleration, α is the angle of attack, β is the angle of sideslip and δA, δE, δR are ailerons, elevator and rudder deflections respectively. The signs of flight control deflections are defined with accordance toHopkin (1970) – a positive flight surface deflection causes a negative aircraft response. This means that: left aileron up, elevator down and rudder to left are considered as positive. The symbols that are not defined yet are known as dimensional stability and control derivatives, e.g. Mα is the pitching moment derivative with respect to the angle of attack. 5. D-optimal multisine inputs Application of the D-optimality criterion for phase angles selection in multisine inputs design requires minimization of the Fisher InformationMatrix inverse determinant (cost function) δ=min φ F −1 (5.1) The cost function hasmultiple local minima, therefore a genetic algorithm is selected for finding the optimal solution.Thephase angles relatedwith j-th flight surface andk-th harmonicφj,k are coded in a binary string that consists of fixed-length strings s. Base-2 floating point representa- tion is used to determine phase angles for different aileron, elevator and rudder deflections φj,k =2π ∑ s bs2 −s (5.2) where s=1, . . . ,12 is a part of the string that contains bit values bs ∈{0,1} used for coding a specific phase angle φj,k. Genetic algorithms mimics the process of natural selection and require multiple executions of four steps that are: population generation, selection mating, exchange of the information and mutation (Mitchel, 1999). The initial population containsP possible solutions that are randomly drawn. The solutions are decoded into input signals and sorted in descending order in accordance with the cost function. Then different selection probability is assigned for each solution. The binary string that represents simultaneous excitations and produces the highest cost function has the smallest drawing probability pmin = 1/ ∑P i=1 i. The binary string that represents excitations with the smallest cost function has the highest drawing probability pmax =P/ ∑P i=1 i. The probabilities are assigned with the linear scale. In the next step, solutions are drawn in accordance to their probabilities and combined in pairs. Exchanging information between the solutions in eachmatedpair (parents) is done by apply- ing uniform crossover which results in creating two new solutions (offspring). Applying uniform crossover requires drawing a mask of flags which has the same length as coded solutions and is filled with zeros and ones. The first offspring is formed from fields of the first parent when the flag of the mask is 0 and from fields of the second parent when the flag of the mask is 1. The second offspring is created in a similar way – from fields of the first parent when the flag of the mask is 1 and fromfields of the second parentwhen the flag of themask is 0. The uniform cross- over idea is shown in Fig. 2. Fields of the first parent have white background and of the second parent – gray. The described procedure is carried for all mated pairs and results in creating a set of new solutions whose size is equal to the initial population. To introduce more diversity in the newly created population, a mutation operator is used. This action requires drawing a number for each field of each solution. When this number is higher than the threshold (0.95), the value of the field is changed from 0 to 1 or from 1 to 0. After performing mutation for all solutions, the initial population is replaced by the new one. In order to preserve the best simultaneous inputs, the parent that has the lowest cost function is also included in this population. Inclusion of the D-optimality in multisine manoeuvre design... 93 Fig. 2. Uniform crossover The described steps: selection, crossing, mutation and population generation are performed until the global maximum is found. 6. Nonlinear model Multisine input signals designed with and without a priori knowledge of stability and control derivatives are used as excitations of the nonlinear aircraft model (Lichota and Ohme, 2014; Raab, 2006). The dynamic equations of motion of the object are given by (4.2). The forces and moments that act on the aircraft consist of aerodynamic, gravitational and propulsion terms F=Fa+Fg+Ft MO =Ma+ 4∑ i=1 ri×Ft (6.1) where r is the thrust force arm and the indices a, g and t indices denote aerodynamic, gravita- tional and propulsion forces andmoments. The aerodynamic forcesFa = [Xa Ya Za] andmomentsMaO = [La Ma Na] are obtained with the use of dimensionless force andmoment coefficients (about quarter-chord point) CD =CD0+kC 2 L CY =CY0+CYββ+(CYp+CYpαα)p ∗+(CYr+CYrαα)r ∗+CYδRδR CL =CL0+CLWBαα+ SH S (CLHαHαH +CLHδEδE) Cl25 =Cl0+Clββ+(Clp+Clpαα)p ∗+(Clr+Clrαα)r ∗+ClδAδA+ClδRδR Cm25 =Cm0+CmWBqq ∗− xLH c SH S (CLHαHαH +CLHδEδE) Cn25 =Cn0+Cnββ+(Cnp+Cnpαα)p ∗+(Cnr+Cnrαα)r ∗+CnδAδA+CnδRδR (6.2) where p∗ = pb/(2VO), q ∗ = qc/(VO) and r ∗ = rb/(2VO) are normalised angular rates, S is the wing area, SH the horizontal tail area, xLH is the horizontal tail arm and k is the drag polar coefficient. The indices WB and H in the equations for the lift force and pitching moment coefficients denoteWing-Body and Horizontal tail, respectively. The angle of attack at the horizontal tail is given by αH =α− ∂ε ∂α α ( t− xLH VO ) + iH +αdyn (6.3) where ε is the downwash angle, iH is the horizontal tail setting angle and αdyn = qxLH/VO is the dynamic angle of attack. The force andmomentcoefficients aredefined in theaerodynamic coordinate systemOxeyeze. TheOxe axis has the direction of the flow, theOye axis is directed towards the right wing and 94 P. Lichota the Oze axis is directed upwards. As the equations of motion are derived in the body fixed coordinate system Oxyz, the following transformation is required Cx =−CD cosα+CL sinα Cy =CY Cz =−CD sinα−CLcosα Cl =Cl25cosα−Cn25 sinα Cm =Cm25 Cn =Cl25 sinα+Cn25cosα (6.4) This allows one to obtain components of the aerodynamic forces andmoments as Xa =CxqS Ya =CyqS Za =CzqS La =ClqSb Ma =CmqSc Na =CnqSb (6.5) The thrust force is a function of the throttle position, altitude and speed: T(δth,h,VO). It is stored in a 3D array and interpolated for specific flight conditions. The gravitational force Fg = mg is defined in the Oxgygzg coordinate system, therefore transformation to Oxyz is required. 7. System Identification The classical and most widely used definition of System Identification was given by Zadeh: “Identification is the determination, on the basis on input and output, of a system within a specified class of systems, to which the system under test is equivalent” (Zadeh, 1962). If the response of the dynamical system and input signals are recorded and the structure of themodel is known then the Sys-ID is equal to adjusting the model parameters Θ so that the model outputsymatch themeasured aircraft response zwhen the inputsu are the same. The process can be done either offline or online in real time (Hendzel and Trojnacki, 2014). Thestructureof themodelweused in theSys-IDwasnonlinearandthe sameas forgenerating the aircraft response what allowed us to eliminate modelling errors in the study. The adjusting of model parameters was done in the time domain when all data was collected (offline). To find the set of unknown model parameters we used the Output Error Method whose aim was to minimize the error between measured system outputs and estimated response. When the Maximum Likelihood Estimation principle is selected for error minimization the task is equivalent to finding the parameters vector that maximize the conditional probability p(z|Θ) Θ̂=argmax Θ p(z|Θ) (7.1) where the hat symbol denotes the estimates, and the conditional probability p(z|Θ) is given by (3.5). As it is frequently done, we replace the probability maximization task by an easier action – minimization of a negative log likelihood function (Jategaonkar, 2006) L(Θ|z) = 1 2 N∑ k=1 [z(tk)−y(tk)]TR−1[z(tk)−y(tk)]+ nN 2 ln(2π)+ N 2 ln(det(R)) (7.2) The unknownmeasurement covariance matrix R is estimated as R̂= 1 N N∑ k=1 [z(tk)−y(tk)][z(tk)−y(tk)]T (7.3) Substitution of the measurement covariance matrix estimate R̂ into the negative log likelihood function and rejection of the constant terms allows one to obtain the cost function J(Θ) =det(R) (7.4) Inclusion of the D-optimality in multisine manoeuvre design... 95 TheGauss-Newton algorithm has been used for cost function minimization. The accuracy of each estimated parameter is determined on the basis of theParameter Error Covariance Matrix P σ(Θi)= √ Pi,i (7.5) where σ is the standard deviation and Pi,i is the element of the Parameter Error Covariance Matrix P =F−1. 8. Results In the study, we investigate two designs with simultaneous multisine inputs: when there is no a priori knowledge of the systemdynamics andwhen this information is available. The frequency range of interest upper bound is set to fmax =2Hz and the lower is determined on the basis of the excitation length T =20s. A total of 39 harmonics is evaluated and assigned to the elevator (f2,f5, . . .), aileron (f3,f6, . . .) and rudder deflections (f4,f7, . . .). The maximum values of all inputs are set to 1deg. In the multisine design without a priori knowledge, the phase angles φk are optimized through Relative Peak Factor minimization. When additional information about the system is available, the D-optimality criterion is used for signal optimization. In this case, the cost function extreme is found by a genetic algorithm. Both designed sets of excitations are used as input signals for a nonlinear aircraft model. In both cases, the response of the model is recorded. The gathered data is used to perform a Sys-ID for a model with unknown aerodynamic parameters. The results of the estimation are shown in Fig. 3, Fig. 4 and in Table 1. Fig. 3. Measured and estimated signals for multisine inputs 96 P. Lichota Fig. 4. Measured and estimated signals for the D-optimal multisine inputs Table 1. Standard deviations Parameter Multisine D-optimal Parameter Multisine D-optimal multisine multisine k 2.5440E-004 3.5897E-004 Cl0 2.2346E-009 7.2915E-010 CL0 3.6488E-005 3.9050E-005 Clβ 2.6671E-006 9.9505E-007 CLWBα 2.6817E-004 3.1839E-004 Clp 1.3783E-004 2.1207E-004 CLHαH 1.1588E-003 1.3114E-003 Clpα 1.6153E-003 2.6471E-003 CLHδE 6.8113E-005 8.6226E-005 Clr 2.3736E-004 2.6195E-004 ∂ε/∂α 1.7339E-004 1.9543E-004 Clrα 3.0391E-003 3.2459E-003 CD0 1.0919E-004 1.6787E-004 ClδA 1.2369E-006 4.8513E-007 Cm0 2.1883E-004 1.8874E-004 ClδR 9.6532E-007 4.0260E-007 CmWBq 1.0397E-002 1.1892E-002 Cn0 2.1906E-009 7.2876E-010 CY0 2.2318E-008 8.4724E-009 Cnβ 2.5441E-006 1.6348E-006 CYβ 1.3871E-005 5.2384E-006 Cnp 1.1095E-004 2.1430E-004 CYp 1.0030E-003 1.6169E-003 Cnpα 1.3218E-003 2.6642E-003 CYpα 1.2198E-002 2.0359E-002 Cnr 1.9973E-004 2.3221E-004 CYr 1.6164E-003 2.0590E-003 Cnrα 2.5561E-003 2.8803E-003 CYrα 2.0288E-002 2.5716E-002 CnδA 1.5451E-006 8.3159E-007 CYδR 4.8513E-006 2.6505E-006 CnδR 1.2249E-006 7.7315E-007 In Figs. 3 and 4, solid lines are used to denote the estimated aircraft response and the input signals. The cross points denote the recorded values of flight parameters. For presenting the measurements, a 10-point data reduction is used. From Figs. 3 and 4 it can be seen that the obtained set of estimates allows one to obtain a good visualmatch for flight parameters whether Inclusion of the D-optimality in multisine manoeuvre design... 97 the D-optimality criterion is used or not. However, as more than one set of model parameters can produce a response that visually fits well the data, another quality of estimation indicator should be used – standard deviations. Table 1 indicates that for both analysed designs the evaluated standard deviations are small (of theorderof 1.0E-002or lower) forall parameterswhatmeans that theestimation is successful. The standard deviations obtained for model parameters when the D-optimality is used in the designphaseareof the sameorderaswhenthe criterion isnotused in theoptimization.Therefore we consider that the inclusion of the D-optimality criterion for phase angles optimization in the simultaneous multisine inputs design does not increase significantly the quality of the Sys-ID manoeuvre. Moreover, as the Parameter Error Covariance Matrix determinant minimization is more computationally demanding thatRelativePeakFactor optimization,we found the inclusion of the D-optimality criterion impractical. Acknowledgments This work has been supported by the European Union in the framework of European Social Fund through the “Didactic Development Program of the Faculty of Power and Aeronautical Engineering of theWarsawUniversity of Technology”. The author would like to thank German Aerospace Center (DLR) for the support and guidance he received during his stay in the Institute of Flight Systems (FT) in Braunschweig. References 1. Etkin B., 1972,Dynamics of Atmospheric Flight, JohnWiley and Sons, NewYork, USA 2. Hendzel Z., Trojnacki M., 2014, Neural network identifier of a four-wheeled mobile robot subject to heel slip, Journal of Automation, Mobile Robotics& Intelligent Systems, 8, 4, 24-30 3. Hoe G., Owensand D.B., Denham C., 2012, Forced oscillationwind tunnel testing for FASER flight researchaircraft,AIAAAtmospheric FlightMechanics Conference,Minneapolis,USA,AIAA 2012-4645 4. Hopkin H.R., 1970, A Scheme of Notation and Nomenclature for Aircraft Dynamics and Asso- ciated Aerodynamics, Aeronautical Research Council, Report 3562 5. Jameson P., Cooke A., 2012, Developing real-time system identification for UAVs, UKACC International Conference on Control 2012, Cardiff, UK, 958-963 6. Jategaonkar R.V., 2006, Flight Vehicle System Identification: A Time Domain Methodology, AIAA, Reston, USA 7. Kay S.M., 1993, Fundamentals of Statistical Signal Process – Estimation Theory, Prentice Hall, Upper Saddle River 8. LichotaP.,OhmeP., 2014,Design andAnalysis of newMultiAxis InputManoeuvres forAircraft Sys-ID, Deutsches Zentrum für Luft- und Raumfahrt (DLR), Report DLR-IB 111-2014/46 9. Mader C.A., Martins J.R.R.A., 2011, Computation of aircraft stability derivatives using an automatic differentiation adjoint approach,AIAA Journal, 49, 12, 2737-2750 10. Mitchel M., 1999,An Introduction to Genetic Algorithms, MIT Press, Cambridge, USA 11. Morelli E.A., 2003, Multiple input design for real-time parameter estimation in the frequency domain, 13th IFAC Symposium on System Identification, Rotterdam, Netherlands, REG-360 12. MorelliE.A., 2012,Efficientglobalaerodynamicmodeling fromflightdata,50thAIAAAerospace Sciences Meeting, Nashville, USA, AIAA-2012-1050 13. Raab C., 2006, A Versatile and Modular Architecture for Aircraft System Simulation and Test, Deutsches Zentrum für Luft- und Raumfahrt (DLR), Report DLR-IB 111-2006/22 98 P. Lichota 14. Rogowski K., Maroński R., 2015, CFD computation of the Savonius rotor, Journal of Theore- tical and Applied Mechanics, 53, 1, 37-45 15. Schroeder M., 1970, Synthesis of low-peak-factor signals and binary sequences with low auto- correlation, IEEE Transactions on Information Theory, 16, 1, 85-89 16. Zadeh L.A., 1962, From circuit theory to control theory,Proceeding of the IRE, 50, 5, 856-865 Manuscript received February 25, 2015; accepted for print July 8, 2015