https://theses.gla.ac.uk/ Theses Digitisation: https://www.gla.ac.uk/myglasgow/research/enlighten/theses/digitisation/ This is a digitised version of the original print thesis. Copyright and moral rights for this work are retained by the author A copy can be downloaded for personal non-commercial research or study, without prior permission or charge This work cannot be reproduced or quoted extensively from without first obtaining permission in writing from the author The content must not be changed in any way or sold commercially in any format or medium without the formal permission of the author When referring to this work, full bibliographic details including the author, title, awarding institution and date of the thesis must be given Enlighten: Theses https://theses.gla.ac.uk/ research-enlighten@glasgow.ac.uk (c) C.G. Black B.Sc. , 1988A METHODOLOGY FOR THE IDENTIFICATION OF HELICOPTER MATHEMATICAL MODELS FROM FLIGHT DATA BASED ON THE FREQUENCY DOMAIN.Colin Gordon BlackSubmitted for tlie Degree of Ph.D. in the Dept. Aeronautics and Fluid Mechanics University of GlasgowProQuest Number: 10970806All rights reservedINFORMATION TO ALL USERS The quality of this reproduction is dependent upon the quality of the copy submitted.In the unlikely event that the author did not send a com p le te manuscript and there are missing pages, these will be noted. Also, if material had to be removed,a note will indicate the deletion.uestProQuest 10970806Published by ProQuest LLC(2018). Copyright of the Dissertation is held by the Author.All rights reserved.This work is protected against unauthorized copying under Title 17, United States C odeMicroform Edition © ProQuest LLC.ProQuest LLC.789 East Eisenhower Parkway P.O. Box 1346 Ann Arbor, Ml 48106- 1346Gantents(1.1) lotoodoctiori.1.1.1 The Nature of System Identification.1.1.2 The Need for System Identification in Helicopter Applications.(1.2) The Beli.copter Model.1.2.1 The Nature of a Model1.2.2 Adequate Models.1.2.3 Model Structure.1.2.4 The Helicopter Model.1.2.5 A Suitable Model for the Identification of Quasi- Static Derivatives.1.2.6 The Theoretical Model (HELISTAB) .- ii -(1.3) Difficulties Associated with. Helicopter System 25Tdeinif-i f i r a t i o n .1.3.1 General Discussion of Identification Difficulties. 251.3.2 Kinematic Consistency Checking. 311.3.3 The Role of Extended Kalman Filtering Techniques. 33(1.4) Obtaining an Adequate Response for BeMoopter 35System Identification.1.4.1 The Design of Suitable Test Inputs. 351.4.2 The Problem of Identifying a Linear Model from 41Non-Linear Data.2. Sonne Pm1imireiT~ios on the Use of the Frequency 47Donum o _(2.1) Introduction: 472.1.1 An Outline of the Previous Use of the Frequency 47Domain in Aircraft System Identification Work.(2.2) The Use of Frequency Donnaim Data. 502.2.1 The Spectral-Analysis Approach. 502.2.2 An Alternative Use of Fourier Transforms. 51- iii -2.2.3 An Expression for the Discrete Fourier Transform 61of a Time Derivative.2.2.4 A Frequency-Domain Representation of the Model. 672.2.5 The Modelling of Noise in the Frequency Domain. 73(2.3) The Jhelusico of Time Belays in the Model. 772.3.1 Representation of Delays in the Model. 772.3.2 Effect of Sampling Interval on the Estimates of 83the Delay.2.3.3 A Demonstration of the Use of a Time Delay Using 88Real Flight Data.3 Tnltial Stages in the Identification IfrMhodolooy., 98and a COn^i derati m of Sonne Relevant Identification Techniques.(3.1) Mtrodoctioo and Overview of the Methodology. 98(3.2) The ErequeocyHDannaiii Eiquatioo-Error Stage in the 98Tdpjpt-i f i nalt-ii on Methodology.3.2.1 The Basic Method and Some Comments on the 100Frequency-Domain Formulation.3.2.2 Subset Regression. 1153.2.3 The Method of Successive Residuals. 120- iv -3.2.4 Parameter Significance Values.3.2.5 Implemen .tation of the Equation-Error Method Using Singular-Value Decomposition.3.2.6 Practical Application of the Singular-Value- Decomposition Approach to Frequency-Domain System Identification.3.2.7 Some Comments on the Use of Zero Pads in Frequency Domain Indentification.3.2.8 Summary of Section 3.2.(3-3) ^plicatioos of the Singulajr-Value-DeoanrpositianL Ftequency-DanriBijni Ekgpatioo—Error Method to Flight3.3.1 Flight Data Sets and Basic Model Equations Used.3.3.2 Selection of Frequency Range.3.3.3 Length of Time-Domain Record Used.3.3.4 Selection of an Appropriate Number of Orthogonal Components.3.3.5 Implementation of Estimates with Theoretical Values.3.3.6 An Equation-Error Multirun Example.3.3.7 Consideration of the Predictive Qualities of a Model.4. The Fg^qpefMryHDoinrrpiin O utput-K rror S ta g e i n t h e125127138141145146146150153154160166167169- V -(4.1) lolbrodoctLoe. 1694.1.1 Background to the Use of the Output-Error 169Technique.4.1.2 Previous Applications of the Technique. 170(4.2) Development of ttoe Frequency-Domain Output-Error 1724.2.1 Formulation of the Frequency-Domian Cost Function. 1724.2.2 Implementation of the Gauss-Newton Minimisation 180Algorithm.4.2.3 Constrained Optimization. 1974.2.4 Controlling Convergence of the Algorithm. 2034.2.5 Some Further Comments on Neglecting Process-Noise 204Terms.4.2.6 Summary of Section 4.2. 209(4.3) AppMcation of the FfequmcyHDoraoiii Guttput-Erxor 210Method to Flight Bata.4.3.1 Applications Involving Longitudinal Test Inputs. 2104.3.2 Applications Involving Lateral Test Inputs with a 227Consideration of Trends in the Estimates forDifferent Flight Conditions.5. Final .sragp of the Identification Methodology and 241a Sunmp^ Ty of the Complete Identification Approach.- vi -(5.1) Estiinati-CEi of Zero Offsets and 3jnd.ti.al Gooditions. 241‘ 5.1.1 Introduction. 2415.1.2 Minimisation of Time-Domain Cost Function. 2415.1.3 Correlations Between the Sets of Parameters. 245(5.2) Summery of the Complete Identification Methodolgy. 2475.2.1 Stage 1: Frequency-Domain Singular-Value-Decomp- 248osition Equation-Error Identification5.2.2 Stage 2: Frequency-Domain Output-Error Identi- 249fication (with options for the inclusion of constraints)5.2.3 Stage 3: Time-Domain Output-Error Identification 251(with options for the inclusion of constraints).6. Goor-llinisioanis and Recennmneodations for Future Mbrk. 253Tables 260Figures 296TNtomentclature 3477. appendices 354Appendix 1: Ti inward ?sed Egnjaflous of Aircraft Motion. 354- vii -a p p e n d ix 2: fininriii-inns -For U n b iased E s t im a te s from t ± e 359F!rec|pQniCvHDom i^iin OutyufcHRrrnr MeffJhod. a p p e n d i x 3 : "It¥* TVS inf*—f^ m in c A i M.qoritlim. 363J^ppeodix 4: Tfr¥* Kalimim-F-i Algorithm. 367Refereix^s/Biblionirarfiv. 371I should first like to express my sincere thanks to Professor David Murray-Smith, of the Department of Electronics and Electrical Engineering, for the encouragement, help, and perceptive guidance given to me over the course of the research. His sharp intellect and clearness of thought has been a great source of inspiration.I am also greatly indebted to Mr Roy Bradley of the Department of Aeronautics and Fluid Mechanics for the help and advice, often outside the call of duty, given to me. His insightfulness into many aspects of the work has been invaluable.Thanks are also due to my external supervisor Doctor Gareth Padfield, of the Royal Aircraft Establishment, Bedford, who, in effect, set me on my present course many years ago. His enthusiasm and wide-ranging knowledge must be acknowledged as an essential contribution to the research.I should also like to acknowledge Professor Brian Richards, of the Department of Aeronautics and Fluid Mechanics, for his efforts in supporting the research project.Thanks are also due to the computing staff: MrsLinda McCormick and Mrs Carol Carmichael, for their cheerful helpfulness over the years.Finally I wish to thank Mr Frank Lombos and hisassistants for typing this document.This research was supported by the U.K. Ministry of Defence (Procurement Executive) through Agreement No. 2048/028 XR/STR'.This thesis is dedicated to my parents.- X -ABSTRACTThere is considerable need for the application of system identification techniques to helicopters. These include their use in the validation and improvement of existing theoretical flight-mechanics models, and for development flight testing. In both cases, estimates of stability and control parameters are sought. Most applications of system identification techniques to helicopters have involved time-domain methods which use reduced-order mathematical models representingsix-degrees-of-freedom rigid-body motion. In thisdocument, an identification methodology which uses the frequency-domain to obtain estimates of the stability and control parameters is advocated.For applications to helicopters, the use of the frequency domain has some considerable advantages for the identification. For example, in the identification of six-degrees-of-freedom rigid-body models, where anadequate description of the fuselage motion is required,the ability to use a restricted frequency range (i.e. one which excludes the higher-order frequencies associated with faster rotor dynamics) is advantageous. In addition to providing a basis for establishing reduced-ordermodels valid over a defined range of frequencies, the use of the frequency domain also results in a significant data reduction in comparison with time-domain methods, thus resulting in a speedier identification. The fact- xi -that the model is represented in terms of Fourier transformed quantities is advantageous for the development of expressions required for the minimisation algorithm, and also facilitates the incorporation of time delays into the model. A suitable means of symbolically representing the model equations in the frequency-domain is presented. Models valid for small perturbations about a nominal trim state are used in the identification.The identification methodology is based upon the use of three distinct identification stages similar to those advocated by Fu and Marchand in their paper presented at the Ninth European Rotocraft Forum at Stresa, Italy, in 1983. Firstly, a frequency-domain equation-error type identification is used to provide initial parameter estimates. Secondly, a frequency-domain output-error type estimation is used to provide under ideal conditions unbiased estimates of the stability and control derivatives. And finally, a time-domain output-error type estimation is used to obtain initial state conditions and zero offsets, and to provide a time-domain verification of the identified model.For the identification methodology developed in this document, additional features to the basic identification methods are advocated. In broad terms these relate to the form of the model used, to the determination of a suitable model structure, and to the implementation of the identification algorithms. An important and original feature of the identification methodology developed hereis the incorporation and estimation of time delays in some of the controls. This is seen primarily as a means of accounting for the effect of rotor-transient effects which are ignored in the quasi-static formulation of the rigid-body motion. It is known from previous studies that unmodelled rotor effects can have a corrupting influence on some important parameters when a simplified form of model is used in the identification; the improvement in the identification resulting from this new approach is demonstrated using real flight-data from a Puma helicopter.The use of techniques (singular-value decomposition for the equation-error method and rank deficiency for the output-error method) that enable the analyst to explore the effects, and likely improvements, as a result of effectively changing the model structure through the introduction of constrained relationships between parameters within the model, is an approach advocated here. It is believed that these techniques have not previously found an application in helicopter system identification. Results demonstrating the usefulness of these techniques are presented, and are obtained using real flight data. For the frequency-domain output-error method, additional means of introducing constraints into the identification are considered.The identification results are discussed in terms of a comparison between the estimated values and those suggested by theory, as well as through the frequency-- xiii -domain fits obtained and the time-domain verifications of the identified model. Results are presented for real flight-data, corresponding to the application of test-input signals on both longitudinal and lateral controls.Some recommendations that can assist in overcoming some of the practical difficulties traditionally associated with helicopter system identification are made. In addition, some new theoretical results and observations relating to the frequency-domain formulation of the identification are made.-1-CH&PTER 11.0 INTRODUCTION: THE HELICOPTER SYSTEM IDENTIFICATION PROBLEM.1.1 Introduction.1.1.1 The Nature of System Identification.In a large number of branches of science, the requirement for a mathematical formulation of the behaviour of a system based on measurements of the response to various inputs, has led to the development of system identification. System identification is commonly referred to as an inverse problem, because instead of computing the response of a system with known characteristics, an attempt is made to solve the opposite problem: that of obtaining the system characteristics from measured responses. In the reference MAIN002, it is suggested that the inverse problem might be phrased as: "Given the answer, what was the question?"In successfully applying the techniques of system identification, the analyst must be able to draw upon results and analytical tools from a wide range of different fields of mathematics (categorized under such titles as probability, statistics, optimization, numerical analysis, linear algebra, mechanics etc.). In addition, the analyst needs to be familiar with the-2-details of the data collection, such as: thecharacteristics of the instrumentation; the level of noise on the measurements; the limitation of information content as a result of the sampling rate, or record length, etc. Even in the case of familiarity with all these aspects, much more is still required if the analyst is to be successful in applying the techniques of system identification. Insight into the system under consideration is important, and subjective judgements (i.e. engineering judgement) have an important role to play in the work.System identification, it could be said, is as much an art as a science - that is, knowledge acquired through experience and intuitive insight is important, and is complimentary to the more exact aspects of this applied subject. In this context, the following sentence taken from MAIN002, relating to the justification for the application of a particular algorithm (axial iteration, discussed infact in Chapter 4) is relevant: "...Oftenwith little more justification than it seems to work well. This is, of course, the final and most important justification."The term: system identification, is the most general description of the inverse modelling procedure, involving both the determination of an appropriate model structure, and the estimation of the coefficients (i.e. parameters) within that model structure; the latter process is known more commonly as parameter estimation, though the termsare often used to mean the same thing. The term, parameter identification may also be used, but this has come to imply a situation where the form of the system model is assumed to be known.1.1.2 The Meed for System Identification in Helicopter Applications.In the context of helicopters, system identification techniques have considerable potential for the validation, or improvement, of theoretical flight mechanics models; this in turn leads to its usage as an aid to the testing of new designs. In addition, during development flight-testing using clinical flight qualities tests, measurements made can be used to estimate stability and control parameters of interest, in order to demonstrate compliance with requirement specifications. A further application involves the use of an on-line identification procedure in an adaptive flight-control system; this however, has yet to be fully exploited.1.2 The Helicopter Model.1.2.1 The Mature of a Model.Since the purpose of system identification is to obtain models of a particular system (in the current-4-context, the system is a helicopter), it would be usefulto define what is meant by a model. In the reference/KRAM001, the following words are said on the subject of models: "The essence of using models is that a material, or formal image of a system is made which is easier to study than the system itself. This image is then used as a model of the system. The model must then obviously contain information about the system. Hence there must be a certain resemblance [emphasis added by the current author] between the model and the system." The same reference gives the following formal definition of a model: "If a system M, independent of a system S, is used to obtain information about system S, we say that M is a model of S."The last definition begs the question as to what is a system. The same reference gives the definition of a system (actually before defining a model): "A system is a set of interrelated entities, of which no subset is unrelated to any other subset." The system M - that is the model shall in the case of the current investigationbe an abstract model, consisting of constant-ccoefficient, coupled, ordinary differential equations. Although the form of the equations used arise from linearisations of the equations of motion, use of the adjective linear is avoided to allow for time delays in some of the terms, rendering the description as linear, inexact. The models identified, however, are based on the requirement of small-scale perturbations in the states,-5-and so it would be equally inexact to talk of non-linear models; instead I shall henceforth refer to them as state-space models (e.g. GELB001) with delays.Model Simplifications.In the reference MAIN002, it is pointed out that a system is never described exactly by the simplified models used for analysis, and unexplained sources of modelling will always be present, since there is no unique correct model. The models used for the helicopter system identification in this document, for example, are based on the rigid-body motion (with six degrees of freedom: three of body translation, and three of body rotation). The most important aspect of the helicopter apparently excluded from this type of model, is that involving the behaviour of the main and tail rotors, and their interaction or influence on the rigid-body motion. However, by using rigid-body models, in conjunction with a new modelling technique developed in this document (viz the use of time delays in some of the controls), an attempt can be made to model the influence of the rotor(s) on the rigid-body motion. The dynamics of the states associated with the rotor(s) (i.e. flap, lag,torsion etc.) are excluded from the state-space representation of the model. This approach is consistent with another feature of the rigid-body models used in the identification work presented in this document: that of-6-using some state variables as "pseudo-controls" in a reduced-order model. For example, in the identification of the longitudinal subset of the six-degrees-of-freedom model, important lateral states may be included as forcing terms along with the controls and this facilitates the formulation of the identification problem in,a manageable form. Examples showing the use of this form of model reduction are presented in Sections 2.3.3 and 4.3.2.Most models require simplifying assumptions, and the analyst is forced by practical constraints to define the system boundary at some point, that is, where the relations between entities are "less concentrated"(KRAM001); in the present context, we might say: "almostdecoupled."1.2.2 Adequate Models.As a result of the practical constraints imposed on the model (e.g. its size), it can be appreciated that there is no such thing as a "true" model. Consequently, the emphasis should be on the estimation of an "adequate" model (i.e. one that is suitable for the purposerequired) and not on the estimation of the "true” model.The myth of the "true" model is very strong, and inMAIN002, the authors quote the comments made by someone on this matter: "A favourite form of lunacy amongaeronautical engineers produces countless attempts to-7-decide what differential equation governs the motion of some physical object, such as a helicopter rotor... But arguments about which differential equation represents truth, together with their fitting calculations, are wasted time."In the context of the work carried out in this document, some degree of validation of a theoretical model is the primary motivation for the identification work. The theoretical model concerned is HELISTAB (SMIT001). This has six degrees-of-freedom, and assumes that the rotor-disc tilts instantaneously for a given stick input (and similarly for the tail-rotor following a pedal input). Consequently, the most adequate model for the purpose of the identification (and validation) will have six degrees-of freedom (with the same states as the theoretical model, so that direct comparisons can be made easily between parameters). However, since the real system, i.e. the helicopter, has a rotor-disc which does not respond instantaneously to pilot stick inputs, but has a finite transient response (which is, in fact, very quick in comparison to the transient behaviour of the rigid-body), it suggests some further modifications to the basic model are required in order to have an adequate representation of the system for identification. In fact, results have indicated that by ignoring the short term transient effects of the rotor on the rigid-body states, the unmodelled "contaminating" effects will result in physically unrealistic (biased, or wrong, call them whatyou may) estimates of some parameters; the most striking example of which, is in the case of the primary-rate pitch-damping parameter Mg (this is considered, in detail, in Section 2.3.4).Accuracy of EstimatesTo make effective use of parameter estimates we must be able to gauge their accuracy, and in order to define accuracy, it is necessary to have an assumed model structure. Once we have decided upon an adequate model structure to use in the identification, then several measures of accuracy can be found. A commonly used statistical measure is the estimated standard deviations of the parameter estimates, also known as the Cramer-Roo bound (MAIN001). As pointed out in the reference, this measure combines two commonly used indicators of parameter accuracy, that of sensitivity of a parameter, and correlation between any two parameters. The estimated deviation is the measure of parameter accuracy used for the work covered in this document. It should be stressed that the derivation of these error bounds is based on assumptions regarding the noise characteristics on the model; these assumptions will be pointed out in due course for the equation-error and output-error techniques used here.-9-1.2.3 Model Structure.It is said in the reference MAIN002, that the real model-structure determination problem is not to determine some non-existent "correct” model structure, but to determine an adequate model structure. If too detailed a model structure is demanded for a data set containing insufficient information, then the result will be a model which appears to fit the given data set closely, but which may be, in reality, a poor and inadequate model, with poor predictive qualities when extrapolated to data different from that used in the estimation. The inclusion of more and more parameters to get a better fit to the observed data, with the additional parameters not really reflecting any physical effect, is known in identification parlance as "over parameterization."Sometimes it is wiser to use a simple, rather than detailed, model; this is known as the "principle of parsimony.” In the reference KLEI003, the authors speak of a "parsimonious model", and seek to obtain, using the equation-error approach, a model (i.e. an equation) that has good predictive qualities, and is not over-parameterized. They define a criterion for testing the predictive quality of a model (the prediction sum of squares), and use it in conjunction with an automatic model structure determination process known as optimal subset regression.In the optimal subset regression procedure, the-10-"best'’ model is selected from a pool of possible candidate models. On the basis of statistical tests performed on the model fit and individual parameter estimates, variables can be incorporated into, or removed from, a model structure. Hence the "optimal" model will consist of a subset of the available variables. In most implementations of the method, it is usually possible to force variables into the model (sometimes called modified stepwise regression; e.g. KLEI002) regardless of whether or not they would be included by the statistical tests, thus introducing a subjective element into the model selection procedure. The subset regression procedure is discussed in more detail in Section 3.2.2.There are many other methods for determining the adequacy of a model structure (e.g. SODEOOl). One such method is the testing of residuals; this has been used successfully by the current author in studies using simulated data, to show the superiority of one model structure over another. Plots of the residuals and their autocorrelation function (which should be impulsive in appearance for uncorrelated stationary Guassian noise) are used (BLAC003).With the frequency-domain identification methodology developed in this document, one important means for determining an adequate model structure, both for the equation-error and output-error techniques, is based on the "degree of singularity" (i.e. the condition number--11-see Section 3.2.5) of the information matrix. In the case of the equation-error method, use is made of the singular-value-decomposition technique to enable solutions to be found, that correspond to a fit obtained using, in effect, a parameter sub-space; this is discussed in detail in Section 3.2.5. There is a direct analogy with the method of subset-regression mentioned earlier (and discussed in Section 3.2.2) in the sense that the "best" model is arrived at using only a subset (or subspace) of available variables. However, in the case of the singular-value-decomposition approach, orthogonal variables - constructed as linear combinations of the original model variables - are used, and it is combinations of variables rather than individual variables, which in general, are excluded from the model.An interpretation that the exclusion of insignificant orthogonal variables from the model, corresponds to the introduction of relational constraints between the parameters in the original model, is offered by the current author in Section 3.2.5. For the equation-error identification technique, using singular-valuedecomposition, the decision as to what constitutes the "best" model can be based (as with most practical applications of the subset-regression procedure), on both objective and subjective criteria.Both the subset-regression and the singular-value- decomposition approaches, ultimately use the correlations that exist between the responses of the variables-12-considered, in the selection of a suitable model structure, and any attempt to say that one method is more correct than another would be inappropriate. The most striking difference between the two approaches, however, is that estimates are provided for all the parameters stipulated toobe in the model in the case of the singular-value- decomposition approach (regardless of the dimension of the parameter sub-space used in the estimation), whilst in the case of the subset-regression approach estimates are only provided for those parameters incorporated into the model on completion of the algorithm. For example, for the pitching-moment equation we know from the physics of the system that an appropriate equation has the form:0^ — M u d -*■ M w w + + M^> p -*■ (1.2—1)For the singular-value-decomposition approach, we will obtain estimates of Mu, Mw, Mq, Mv, Mp and M/\V5 , regardless of the model structure ( defined by the constraints relating dependent sets of parameters, and determined on the basis of the data used in the identification).Applying the subset-regression procedure to the same set of data, we may find that the proposed adequate model structure is:h = Mod + M w w f\\$ (1.2-2)-13-Thus providing estimates of the parameters Mu, Mw, Mq, and M/\1S) , but not of Mv and Mp (except that is it implied that they can be considered to be zero for the data used in identification).The equation-error method is seen, for the identification methodology developed in this document, as a means of providing initial guesses for the more advanced output-error method. The singular-value- decomposition approach provides a full set of parameter values for use in the output-error method, whereas the subset-regression procedure, in general, does not. It is also of interest to point out that the speed of application of the singular-value decompostion procedure to the problem of helicopter system identification, is greatly facilitated by the use of the frequency domain, because of the small number of data used.The equation-error approach, in addition to providing initial estimates for the output-error method, can also assist in the determination of a suitable model structure through the use of some measures of parameter significance, as defined and described in Section 3.2.4. In addition, some investigation into the use of different model structures for the output-error identification stage can be performed using rank deficiency. This is described in Section 4.2.3, where the method is implemented using a singular-value decomposition of the information matrix. Here, it is found that it may be appropriate to reduce the effective dimension of the-14-parameter space used in the iterative estimation technique. This amounts to introducing relational constraints on the update increments.The preferred techniques for obtaining a suitable model structure for both the equation-error and output- error methods (referred to by the current author thorughout this document as singular-value-decomposition and rank-deficiency respectively, for the purpose of distinction) are demonstrated in Sections 3.3.4 and 4.3.2 using real flight data.1.2.4 the Helicopter Model.The helicopter equations of motion used for the work in this document, correspond to the form of model defined in the reference PADF009. Here the applied forces and moments acting on the vehicle are incorporated into equations, with a right-hand set of body-fixed axes X,Y,Z defining the rigid body, having vehicle translational velocity components u,v,w, and vehicle rotational componenets p,q,r, along and about the X,Y and Z axes respectively. The Euler angles Y , © ,

U = ^ P)T(The states X in the model, and the measured quantities Z are in general, related in a non-linear manner; this is discussed in Section 2.2.4).A more general representation of the helicopter, will also include equations with states representing the behaviour of the rotor, and its coupling with the rigid- body motion of the helicopter. It is convenient to represent the rotor by means of the attitude and shape of the tip-path plane, with respect to the shaft, instead of considering individual blade motions (whose description requires a system of coordinates rotating with the rotor). This is achieved using a multi-blade coordinate transformation (e.g. PADF009) and corresponds to the representation of the rotor as a continuous and periodic function of azimuth position. For a four-bladed rotor, the three terms of a Fourier series representation, corresponding to the average, and first harmonic sine and cosine terms, and a fourth term corresponding to the average of the second harmonic sine and cosine terms, can be obtained. In the case of blade flap, these correspond to j8o (coning angle) ; fl\c. (tip-path plane tilt in pitch) ;-16-P,* (tip-path plane tilt in roll) ; and /3d (differential coning). They are obtained using (PADF009):A - A | aS is. = 2A £ JSi t»nYt(1.2-4)s.c = zA i & cos n' LI tjB*. -- A i / J i t - i ifor a four-bladed rotor, with individual blade flapping measurements jBi., and rotor azimuth position YJ, .Identical expressions can be defined for the laggingmotion of the rotors.If we define the rigid-body states to be , andthe rotor states to be Xft / then we can write thefollowing linear constant coefficent model (KALE001):XR = A rf X f + A*R X* + Br U (1.2-6)Where the matrices A ff and A rR represent theuncoupled systems for the fuselage and rotor. Thematrices A fr and A rf represent the coupling between thefuselage and the rotor; Bp and B R are the fuselage androtor control dispersion matrices.It is pointed out in the reference KALE001 for a model with 14 degrees of freedom (i.e. 6 rigid-body, 4= A FF Xp + A ff, XR + BfV)-17-rotor flap, and 4 rotor lag), with the rotor modelled as a second-order system, that there are 32 6 state matrix coefficients (of which, - about 300 have to be identified) . And this is after two approximations have been made: namely the small perturbation assumption, and the averaging of periodic terms in (1.2-5) and (1.2-6) to give constant coefficients. This is obviously an unrealistic number of parameters to even consider attempting to estimate. The comments made in Section1.2.1 about the practical constraints which force the analyst to make necessary simplifications, are relevant in this context. The simplifications which can be carriedout here include neglecting lead-lag states (on the basis/that the frequencies associated with them are well separated from those associated with the rigid-body states), thus representing the rotor dynamics in terms of flapping motion only. In turn, the rotor flap motion can undergo a series of simplifications, whereby the second- order system can be simplified to a first-order system, and ultimately to a zeroth order system. The final simplification (i.e. zeroth order system) is a form of quasi-static model. It assumes an instantaneous response of the rotor following a pilot control input, and is a form of model which has found applications in real time helicopter simulations (e.g. PADF009). Again, this assumption is justified on the basis that the frequency seperation between the rotor modes and the fuselage rigid-body modes is high: for the time scales typical of-18-the rigid-body motion, the behaviour of the rotor states appears to be instantaneous.The Quasi-Static Assumption.Using the assumption that the rotor disc tilts instantaneously following the application of a pilot control input, a form of the model given in (1.2-5) and (1.2-6) can be written which considers fuselage statesonly. Setting X & in (1.2-6) to zero gives (KALE001) :Xfl = - A r*'*A*f X f - (1.2-7)And substituting for in (1.2-5) we have:^ ( Apf - Apfv A rr - A rf) Xp +■( B f " A rr - U (1.2-8)- A*F X f + B f U (1.2-9)The contribution to the rotorcraft motion made by therotor is lumped into the fuselage coefficients. The point can be made that the derivatives, identified for a model of the form given in (1.2.9), but with data which is not quasi-static, will in general be different to thetheoretical quasi-static values.The effects of rotor/fuselage coupling on derivative estimation was demonstrated in MOLU003, by attempting an-19-identification of a six-degrees-of-freedom model from simulation data that included rotor flapping modes. Substantial differences existed between the predictions and the quasi-static values. However by including the multiblade flapping coordinates f*>or /3IC, and /31S , in the assumed model, and performing the identification, and then reducing to a six-degrees-of-freedom form (as in1.2-8), the values obtained agreed very well with the calculated quasi-static values. The conclusion was that the quasi-static mode was quite adequate in representing the long-term helicopter motion, but that profound difficulties could be expected in trying to estimate quasi-static model parameters using fuselage states only, in the estimation. The most serious discrepancy was the underestimation of primary-rate damping derivatives such as Mq.The inclusion of rotor states in the estimation is one possibility suggested by the above conclusion. However, this would require a fairly large increase in the number of parameters to be estimated. The technique developed by the current author in the following section, requires only the inclusion of one extra parameter, namely a time delay, for an improved estimation of quasistatic derivatives.-20-1.2.5 A Suitable Model for the Identification of Qnasi- Static Derivatives-The quasi-static model containing fuselage states only is given in (1.2-8), and is written compactly in(1.2-9). To be exact, an additional term is required in(1.2-9) to represent the transient rotor effects. We may rearrange (1.2-6), to obtain (c.f. 1.2-7):A rr ‘ A rr X F - A r r - B ^ U + A r r X r (1.2-10)Substituting for ^ in (1.2-5), and rearranging as in (1.2-8), we have:= A FF Xp * B F u * ( 1 . 2 - n )Hence it is the unmodelled input - A Ff\* A kr X r representing the rotor transient effects, which is responsible for the poor estimates of some quasi-static parameters in the matrices A FF and Bp . It can be appreciated that the unmodelled input A pr' A rk X k will become significant for a small time following a movement of the stick (or pedal) by the pilot: the time span of this transient term will be characterized by the time-constants of the corresponding important main rotor (or tail rotor) modes. There is clearly a requirement for aremodelling of the effective control input. We requireA ^the estimates of the quasi-state derivatives: A ff anci-21-PU p for direct comparison with theory, and we shall seek a model for the identification with the following form:X F = A f f X f + b f* y eP (1.2-12)where Ue-p is the effective control input to the model. A comparison of (1.2-11) and (1.2-12) implies that:Bpu A*cA** X* = B p ( i . 2-i3)When - 0 (or near enough) then we have the effectivecontrol input Uep equals the measured pilot input U; this will be the case for a step input after a time-span ofthe order of the time constants of the important main rotor (or tail rotor) modes has elapsed. When X r ^ 0, theonly meaningful interpretation of (1.2-13) is aleast-squares solution for The dimensions of Bp aren x m (for n fuselage states and m controls). Since n > m, we have an overdetermined system of equations which may be solved for U<=f using the following pseudo inverse (e.g GELB001);(1.2-14)Hence using (1.2-14) in (1.2-13), we have:U cf - U + ( Bp B f Bp A fr A rr Xf r\FR h Aft (1.2-15)-22-An improvement on the assumption of an instantaneouscontrol input, is to use a first-order model. If weHence (1.2-15) can be written as a scalar equation.The effective input is composed of the pilot input and an additional transient term related to the rotor states.An approximate approach adopted here, incorporating transient and steady-state behaviour, involved modelling the effective controlinput Uef as a first-order lag, driven by the applied or measured pilot input U. The time constant would be expected to have a value characteristic of the important main rotor (or tail rotor) modes. We can thus write:Transforming (1.2-17) into the frequency domain (see Section 2.2.2) and for a periodic measurement window of length T (U(o) = U(T) = 0 ) we haveresponse of the rotor states X«, following a pilotignore all but the most important element in X r (assuming Xft represents rotor flapping states), this will leave, for example, /Sic for a longitudinal cyclic input /"|is(1.2-16)Uef * Uef = / x U (1.2-17)(1.2-18)where Uep(0) = O’, since U(0) = 0 and using (1.2-16).Now for |o*c| << 1, (i.e. for small frequencies and small time constants) we can write (1.2-18) as:-23-U epCco) + O Cco^ -c1)) U Cod) (1.2-19)Since for |00TT I « 1 it is true that, r*. r *.-r-3Ac. ■£* \ - joot + O C ^ ; (1.2-20)we have:A - J C O X AU e f M — e. UCo) (1.2-21)In the reference ISER001, and in the context of transfer-function representations, the replacement of small time constants in the denominator by delays in the numerator, is a suggested simplification; this corresponds to the approximation given by (1.2-20). Converting (1.2-21) back into the time domain we have.Uef CO = U ( t - r ) (1.2-22)Hence a more suitable model for the identification of quasi-static derivatives ( A ff and B f ) using fuselage states and the measured or applied pilot control input, is given by:%p = Apr Xf + BpUCt-t) (1.2-23)where T is an unknown parameter to be estimated alongside-24-elements of App and Bp . The estimated time delay T may also contain a "transport" component in the control.1.2.6 The Theoretical Model (HELISTAB).The helicopter flight-mechanics software package- HELISTAB, developed at the Royal Aircraft Establishment, Bedford (SMIT001), is used to provide the theoretical quasi-static parameter values for comparison with estimates obtained from flight data. In addition, thepackage can be used to generate simulated non-linear orlinear time responses, which is useful both for the validation of estimation techniques and software. It has - options for a range of different degrees of freedom (DOF) in the model: 6D0F (quasi-static model); 9 DOF (withfirst order modelling of main rotor flap); 12 DOF (with second order modelling of rotor flap).The theoretical non-linear model (PADF00 9) - for a specified helicopter configuration - is linearised for a given flight condition in order to provide the theoretical quasi-static linear model values. A trim routine calculates the control settings for the givenflight condition, and then the force and momentderivatives are calculated numerically using forward and backward differencing about the trim point. These are then combined with the linearised kinematic and gravitational effects to produce the complete linear model.-25-1.3 Difficulties Associated with Helicopter System Identification.1.3.1 General Discussion of Identification Difficulties.The application of system-identification techniques, especially in the time domain, has a relatively long history of success in the context of fixed-wing aircraft (e.g. ROSSOOl, KLEI007, FOST002, FOST003). The pioneering work of Molusis in helicopter system identification (MOLUOOl, MOLU002, MOLU003), and work published by subsequent authors (e.g. PADF002, PADF005, BLAC007, DUVA001), highlighted some of the difficulties of applying identification techniques to helicopters. These difficulties, including those encountered for the frequency-domain methodology developed in this document, will be considered now:1. System Complexity: As was discussed in Section 1.2.4, the helicopter represents a complex, and highly coupled system. Strong coupling exists between the longitudinal and lateral fuselage dynamics. Consequently a large number of parameters are required for a model approaching something of a global representation of the helicopter. It is pointed out in KALE001 that such "global" models may be very inefficient for use in problems such as control system design, stability and control analysis,-26-and direct correlation with flight test data for validation purposes. There is therefore a need for simplifications in the model representation.The justification for the assumed form of model used in the identification work covered in this document, has been given in the previous sections. It was also mentioned that a further reduction in the model complexity can be obtained by excluding lateral dynamics, but incorporating lateral states in an extended control vector for a longitudinal input (and vice versa). This reduces the identification problem to a manageable size for the output-error method, and is an approach advocated in the helicopter system identification methodology developed in this document.Identification problems resulting from fuselage- rotor interaction were discussed in Section 1.2.4, and a new feature of the assumed model structure, namely the use of a time delay, was presented as a means of overcoming to some extent these problems; this circumvents the need for introducing the additional complexity in the estimation model that would occur with the inclusion of rotor states.2. Hiah-Vibration Environment: Helicopter vibration is caused mainly by the higher harmonic components of the blade flapping and lagging motion, in addition to blade bending modes. The high-vibration environment reduces the signal-to-noise ratio. A high degree of uncertainty-27-or noise in the measurements results in a greater degree of uncertainty in the parameter estimates finally obtained.3. Instabilities: The inherent instability of the helicopter restricts the length of the data recordsavailable for estimation. Data records which are short in comparison to long-period modes, may result in identification problems; the convergence of parameterestimates with record length is discussed and demonstrated in Section 3.3.3. The smallest frequency which can be used in the identification is determined solely by the length of the data record; an artificial means of reducing the spacing between frequencies, namely by the use of zero pads, is shown in Section 3.2.7 to be worthless for identification purposes.Because of the instab .ility of the helicopter, difficulties can be experienced in the practical application by the pilot, of test input signals which may have desirable spectral properties for systemidentification work, but which may result in an unacceptable response of the aircraft.The use of stability augmentation systems can enable longer data records to be obtained. However a decreased accuracy in the estimates of the parameters can be expected because of the likely existence of strong correlations between the states and the controls, when feedback loops are used. The identification of a highly-28-augment ed airplane is considered in the reference BATT001, where time-domain stepwise regression and maximum-likelihood procedures are used; it was concluded for the identified model, that some degree of correlation between the parameters must be accepted. The use of rank-deficiency as a means of establishing a suitable model structure when strong correlations or nearly perfectly defined relationships exist between groups of parameters, is one aspect of the identification methodology used by the current author. The use of rank- deficiency, together with the longer data records that would be possible as a result of using a stability augmentation system, may well offer some new possibilities for improved helicopter systemidentification.The use of a stability and control augmentation system in obtaining response data for system identification work, is also mentioned in TISC001, where it is stated as a key requirement that the total surface deflection (i.e. the control input signal) must contain a significant component from the pilot which is uncorrelated with the output.4. Non Linearities: One of the drawbacks of the frequency-domain identification methodology is that it is restricted to linearised state-space models with time delays. As a result there is a need to avoid large excursions of the helicopter from the trim, in order not-29-to violate the small-perturbation assumptions on which the model is based. Consequently there is a requirement at the data-collection stage to use control input signals which will induce a favourable type of response in the helicopter. The effect of control input shape and the effect of non-linearities is considered in more detail in Sections 1.4.1 and 1.4.2.5. Errors Associated with Measurement Devices: Some measurement devices (e.g vanes measuring sideslip and incidence angles, and the speed probe) are susceptible to rotor wake and fuselage flow-field effects. Dynamic lags may also be present in some measurement devices, such as accelerometers, and as such, should be taken into consideration when the measurement system is modelled. If measurement channels with lags were to be used in an identification, then for the frequency-domain output- error program: OUTMOD, developed by the current author, the corresponding measurement channel could be modelled with a time delay, for both small frequencies and small accelerometer time constants (See Equations 1.2-17 to1.2-22 in Section 1.2.5 for some justification of this).It may be the case that scale factor or bias errors are present in the data; these are considered in more detail in Section 2.2.4. For the use of frequency-domain data, constant biases have no effect on the estimation, when the CO =0 frequency is excluded from the stipulated range. In addition, it is possible to estimate-30-measurement scale factors occurring in the linearised measurement transition matrix (relating measured responses to states in the model) using the frequency-domain output-error program: OUTMOD.The locations of measurement devices, relative to the centre of gravity, have to be taken into consideration when modelling the measurement system, although these can also be estimated as unknown parameters. Offsets relative to the centre of gravity are discussed in Section 2.2.4.6. Correlations in the Data: A specific difficulty resulting from highly correlated lateral response variables, associated with a "Dutch-Roll" type mode, is considered in Section 4.3.2. This is thought to be largely responsible for the physically unrealistic estimates obtained for some lateral parameters in previous studies (PADF002). The use of rank-defLciency in the output-error method is demonstrated in Section 4.3.2, as a technique, which through the establishment of a more appropriate model structure, leads to improved estimates.Some Further Remarks: A list of further, and possibly insidious, sources of error in the identification would be endless. However, some additional problems encountered by the current author include: errors in the calibration (introduced at the data preparation stage); the use of incorrect measurement units (i.e. a scale- factor type-31-error) the use of an inappropriate data window for the identification (e.g. the inclusion of a portion of data, usually at the end of a measurement data file, corresponding to the pilot recovery period, and which is not adequately modelled using only one control input in the identification model).1.3.2 Kinematic Consistency Checking.Errors in the flight due to instrumentation failure, or errors introduced at the pre-estimation processing stage (i.e. conversion into standard measurement units using calibration values) can cause a lot of wasted hours at the estimation stage. However, the presence of such errors can be revealed by performing a Kinematic consistency check (REID001, FEIK002) on the data. A Kinematic-consistency check program was developed by the current author (KINECON) for use in such cases.The Kinematic equations of translational motion and associated Euler equations of rotational motion, shown below, from the basis of the computer program - KINECON.uvw (1.3-1)p + o 0 h w )<9 + r c o s j z f h x n ©<3Cos0 - rsm^J i3 S e c © + ^ C o s $ S e c ©-32-The above set of equations may be integrated (using Euler's method or a Runge-Kutta method) to obtain time histories of the velocity components: U,V,W, and Eulerangles: 0 , © , Y • T^e integration is driven by measurements of the accelerations: ax, ay, az, and rotational velocities: p, q, r. A comparison is made with estimates of U, V, W, obtained independently of ax, ay, and az, by solving iteratively the following non-linear set of coupled equations (shown in a form that includes the offsets relative to the centre of gravity in Section 2.2.5, and implemented in the program in such a form).(u* * v' + v j1)'/z* ( 1 . 3 - 2 )yS -- tan ( v /u )CK * tan' ( w /u )Measurements of speed, flank angle and incidence angle: are obtained for use in the aboveprocedure.A graphical comparison of "measured1 time histories of: U, V, W, , © , V5p , jB - (X and the correspondingquantities obtained by integration of (1.3-1) can then be made. This can highlight bias errors which show up as drifts in the computed time history (see Figure 1-1) .A procedure for estimation of accelerometer biases is also a feature of the program KINECON. The biases areestimated by an iterative regressive procedure whichattempts at each stage to fit a least squares regressionline through the residuals obtained for the U, V and W comparisons. Using the slopes obtained, estimates for the biases are found, and the corresponding accelerometer measurements are corrected. Integration of Equations (1. 3-1) is then repeated, and a new comparison of the observed and estimated quantities U, V and W is made. The procedure is repeated until a satisfactory match between the observations and the estimates is obtained. The use of the technique is illustrated in Figure 1-1, and has proven its worth in application to a large number of flight test records. A more general facility for handling bias and scale factor errors is given by the identification methodology developed in this document.1.3.3 The Role of Extended Kalman Filtering Techniques.A problem encountered in the application of system identification techniques to helicopters, is the presence on some measurement channels of large amounts of noise. State estimation, or Kalman filtering techniques seek to reduce the amount of uncertainty associated with a given signal. A computer program specifically for rotorcraft state estimation was developed by NASA, called DEFKIS (HALLO02). The program represents the implementation of an Extended Kalman filter and smoother (e.g. GELB001), where the non-linear model is based on the Kinematic equations of motion, similar to those given in (1.3-1), but supplemented by further equations incorporating-34-additional states. The program can also be seen to play the role of a Kinematic consistency checker, and in addition, allows for the reconstruction of unmeasured states and state time derivatives. The DEFKIS program is an integral part of the time-domain helicopter system identification package PEP, currently implemented at the Royal Aircraft Establishment Bedford, and outlined in the reference PADF005.A significant disadvantage in the use of DEFKIS is the sheer complexity of the input data required for successful operation of the program. The user is required to specify certain noise statistics (measurement and process) which may not be known apriori with confidence. Execution time of the program may also be considerable.The frequency-domain methodology developed here is capable of using the untreated measured flight data, where it is transformed directly into the frequency domain. For the frequency-domain equation-error method, it is shown in Section 3.2.1 that measurements of time derivatives need not be available for the estimation to proceed; in addition, the Fourier transforms of states to be included in the model, are easily obtained from the transformed measurements.The possible need for the use of an Extended Kalman Filter/Smoother Program, such as DEFKIS, within the frequency-domain identification methodology, is discussed in Section 4.2.5 in the context of noisily observed-35-inputs to the output-error method.1.4 Obtaining an Adequate Response for Helicopter System ^Identification.1.4.1 The Design of Suitable Test Inputs.Mode ParticipationIn order to be able to perform a satisfactory identification, it is necessary for the measured responses to contain sufficient information about the system. For a linear system of the form given in 1.2-3, with n states, the following expression can be written for the free response (PADF003):y (0 - ^ c v L • x(o))c d.4-1)where* Xc o ) is the value of the state vector following the application of the test-input signal; U i. are the eigenvectors of the system matrix A; VL is the ith row of U t where U is the n x n matrix formed by closing together in ordered form the n eigenvectors of A.The terms in (1.4-1), in parenthesis, are scalars known as mode participation factors. If X ( o ) is proportional to one of the Ul, then only that mode will respond. For the free response, the relative degree to which each mode participates in the motion is determined-36-by the participation factor.It is because of the fact that the application of a test input signal, on a given control, may not adequately excite all the required rigid body modes, that the combination of several data sets for the identification is advocated, and even stated as necessary, by some authors (e.g. KALE001). This is not, however, a universally held belief: a methodology based on single run equation-error evaluation in the time domain is described in the reference (DUVA001), and is considered in Section 3.2.3. The use of multirun data for the equation-error method is demonstrated in Section 3.3.6 using real flight data.Multi Step InputsThe design of suitable control inputs for aircraft system identification studies is an area which has received considerable attention (e.g. GUPT002, STEP001, AGAR001). The object of these test inputs is to excite the system so that the modes of interest are participating significantly in the response. In the context of helicopter system identification, the use of pseudo-stochastic multistep inputs has been proposed (e.g. MURR002). Originally developed for fixed-wing applications, these signals can be tailored so that the frequency components of the signals span the passband of the system under test. However, it should also be pointed-37-out that sharp-edged inputs might excessively excite the higher-order rotor modes, which is an undesirable feature for the identification of rigid-body models (DUVA001).A commonly used signal of this type is the n3211"; this is shown in Figure 1-2, together with its auto spectrum, for clock periods of 1.0 and 0.5 seconds. The "3211" input is often favoured as a test input because it has an auto spectrum which closely approximates that of a band-limited white noise source, which has an equal power distribution over a certain frequency range, and so excites all modes within the band equally.For the flight data sets analysed by the current author, the "3211" inputs occurred with a 1.0 second clock period. It can be appreciated that as the clock period of a multi-step is increased, there will be an increase in the overall power content of the signal; this is shown in Figure 1-2 for a "3211" input. The problem with applying a signal which has a relatively high power level throughout the range of frequencies associated with the rigid-body modes (0-3 rads/s), is that it will induce large excursions of the aircraft, invalidating the small- perturbation assumptions of the model. This is likely to be a significant reason why better results were obtained for doublet inputs (shown in Figure 1-3, together with auto spectrum), of shorter duration, by the current author, for the flight data sets considered.The above findings have been supported by some additional research, subsequently carried out at Glasgow-38-University (LEIT001), in which it was shown that the "3211" input, with a 1.0 second clock period had a lower coherence function (BEND001) than the 1.0 second clock period doublet, over the range of frequencies associated with the rigid-body motion. A lowering of the coherency at particular frequencies can be caused by insufficient excitiation and non linearities; both of which are undesirable features from the point of view of identification.An identification carried out by the current author (BLAC003) using simulated data generated by the HELISTAB program, sought to compare the 0.5 second "3211” input and a doublet of similar duration. For the regression fits obtained, it was concluded by examination of the residuals (and their auto-correlation functions, where for uncorrelated Guassian residuals it is an impulse function) that the 0.5 second "3211" was the better of the two. This conclusion was also backed up by consideration of the coherency functions of the two inputs, shown in (LEIT001).Another factor which has to be taken into account, for a given multistep shape, is the amplitude of the signal; increasing the amplitude, for a given clock period, has the effect of scaling the autospectrum. A large-amplitude input will thus also induce an undesirably large response of the helicopter. Some results of an investigation carried out by the current author (BLAC003) into the effect of amplitude and clock-39-period on estimates of lateral derivatives, using non-linear simulation data generated by HELISTAB, are shown in Figure 1-4. The deteriorating effect of both increasing clock-period and increasing amplitude for "3211" input is demonstrated, where normalised lateral derivative estimates (normalised with respect to the corresponding linear model values) are showrv. In reality, of course, very small amplitude inputs will be unrealistic in that they need to be large enough to overcome the dead zone of the control; in addition, a very small response would be overwhelmed by measurement noise.It is clear from the aforegoing discussion, that in the design of suitable test inputs for helicopter system identification, more consideration than has been given in the past is needed of the factors likely to induce severe non-linearities in the response. Previously, frequency characteristics based on the general shape of the input have been given primary consideration. Recently, however, a procedure for the design of multi- step inputs that are suitable for helicopter system identification, has been developed, which takes into consideration such factors, and some new input types have been designed (LEIT001)Other Test Inputs UsedIn addition to multistep inputs, some other types of test input signals have been used in the past. In-40-(TISC001) the authors use frequency sweeps (or swept sine waves) in the identification of single-input singleoutput transfer functions for the XV-15 tilt rotor aircraft; these inputs are cited as ideal because they result in bounded and reasonable excursions of the airgraft, and suitably excite the important rigid-body modes over the entire frequency range of interest.However it is believed that fatigue problems may arise with frequency-sweep inputs when applied manually through the pilot's controls. With manual inputs it is difficult to limit the magnitude of the perturbation at the resonant frequencies of the aircraft. With an appropriate control input device, these problems can be overcome by limiting the amplitude of the test signal at critical parts of the frequency range.In the context of time-domain helicopter system identification, the authors in the reference. DUVA001 consider a M3211", a doublet, and a low frequency sinusoid as test-input signals. Simulated data representing the Rotor Systems Research Aircraft, with rotor degrees of freedom was used. The simulations included cases with a stability augmentation system (SAS) switched on. It was found that the use of a "3211” and doublet, both of which have significant high frequency components capable of excessive excitation of rotor modes, resulted in considerable SAS activity. (The detrimental effect of the use of stablility augmentation systems on an identification was discussed in Section-41-1.3.1). In contrast, the sinusoid input provided only low frequency excitation, reducing the SAS activity, which in turn lead to better regression estimates of longitudinal parameters. However, it was also concluded in DUVA001, that in the event that the SAS can be left off for lateral excitation, it may prove more effective to use inputs with higher frequency components (i.e. a ”3211").1.4.2 The Problem of Identifying a Linear Model from Mon-Linear Data.The extent to which the response of a helicopter in flight-test conditions is non-linear is determined, to a large extent, by the test input signal used; this was considered in the last section. A disadvantage of the frequency-domain methodology advocated here, is that it is limited to the identification of small-perturbation models. However, for the state-space formulation assumed for the identification, it is possible in theory to include specific non-linearterms in the model by including non-linear'states in an extended control input vector. The state-space model would have the form:rucolX C t ) = A X ( t ) + [ B S 6*1 i r c t ) ( i . 4 - 2 )where S Ct) is the non-linear input vector (corresponding to higher-order terms in the Taylor Expansion of the non-linear model); and B is thecorrresponding partition of the enlarged control dispersion matrix. The use of higher-order aerodynamic terms in time-domain fixed-wing system identification is discussed in KLEI002 using stepwise regression and in ROSSOOl for the output-error approach where a formulation of the model identical to that given in Equation (1.4-2) is considered. The inclusion of such non-linear terms in the model will mean an increase in the number of parameters to be estimated. In the case of helicopter system identification, an increase in the number of parameters should be avoided at all costs, and it is thus better to reduce non-linear effects at the outset than to introduce additional complexity into the estimation.Non linearities corresponding to the known kinematic terms (see Equations 1.3-1) and gravitational effects (known components of the accelerations ax,ay and az given in Equations 1.3-1) could also be included as non-linear inputs (assumed to be noise-free). This would mean that linearisations of the aerodynamic force and moment terms only are considered in the identification.An investigation was carried out by the current author using non-linear simulation data from the HELISTAB package, in order to gain an appreciation of any difficulties that can be expected in the estimation of a linearised version of a non-linear system. The results are presented in detail in the reference BLAC003. The general conclusions were although the linearised reponse (when added to the trim values) may be very-43-similar in appearance to the non-linear response, for a given control input, the estimates obtained, using the equation-error approach, can be quite different for some of the weakly-defined parameters. Differences were observed in both sign and order of magnitude.The time-domain residuals from model-fits when there was a substantial amount of non-linearity were quite deterministic in form. The presence of deterministic unmodelled non-linear effects is effectively the same as having a correlated equation-error term, and as shown in Section 3.2.1, will lead to biased parameter estimates. It was demonstrated in the simulation study, that small non-linear effects have a much more detrimental effect on the estimates of parameters, than a 'larger' amount of white measurement noise on the variables. When estimates obtained either from non-linear simulation studies or from actual flight data are compared with the theoretical HELISTAB linear model values, it should be noted that the HELISTAB values are obtained by numerical differencing of the non-linear model and as such, are dependent on the prescribed perturbation values used in the differencing. The perturbation values used in the linearisation can be varied.One observation which is encouraging for the identification of small-perturbation models, was that for the longitudinal derivatives considered for the investigation, important derivatives such as the pitching moment derivatives: Mu, Mw, Mq, Mv, Mp, M: l^s , and the-44-normal force derivative Zw, gave estimates which were close to the calculated linear model values. A reduction in the amplitude of the control input used, brought about an even closer match between the two sets of values. Estimates of the aforementioned parameters obtained from flight data, using the frequency-domain estimation techniques developed in this document, are presented in Chapters 2 and 3.Sumnmniarv of Chapter 1./The general concept of system identification was discussed for the specific application considered here, namely the identification of quasi-static helicopter models involving fuselage states only.-. Systemidentification was presented as a means of verifying, and updating,q. theoretical helicopter simulation model.A form of model was proposed for use in the identification which was considered to be more appropriate than those used hitherto for the identification of quasi-static derivatives from real flight-data. For the model, which involves fuselage states only, justification was given for theincorporation and estimation of a time delay in some of the controls as a means of accounting for contaminating rotor-transient effects.The role of the helicopter flight-mechanics package HELISTAB in the identification methodology as a source of-45-theoretical values for comparison with those obtained from flight-data, and its facility for the generation of simulated responses corresponding to known models which can be then used to test identification techniques and their software implementations, was pointed out.Some of the difficulties which are associated with the application of system identification techniques to helicopter flight data were enumerated, and some means of alleviating them were suggested. For example, it is possible to derive very comprehensive and high-order linear models for a helicopter; however, their complexity precludes the use of such models for system identification purposes. It was stated that modifications, such as the use of a subsystem in the identification (e.g. longitudinal motion states only) driven by, in addition to the pilot control, pseudo control terms (e.g. lateral states), could allow the identification problem to be cast in a manageable form for an output-error identification. The removal of rotor states from explicit representation in the model, and the use of time delays to model the influence of rotor- transient effects, represents a simplification to the model which is suggested here.Some means of checking the flight data, prior to identification, were outlined. The role of an Extended Kalman Filter Procedure in the context of helicopter system identification was discussed, particularly as a means of constructing unmeasured states and removing-46-noise levels prior to a time-domain identification. It was stated that the frequency-domain approach advocated here, also has the facility to effectively construct unmeasured states, and that the ability to successfully apply the identification techniques using raw flight-data, thus avoiding the complexities of an ExtendedKalman Filter, is advantageous.Finally, the problem of obtaining a response from the helicopter which would facilitate a successful identification, were discussed. These included the need for adequate model excitation in addition to the requirement that the small-perturbation assumptions, underlying the form if model used, were not to beviolated. Some commonly used test-input signals weredescribed, and it was argued that more attention should be paid to the undesirable effects on the motion of the helicopter of large amplitude and large clock periods for multistep inputs, such as the commonly used '3211'.-47-CHAP.TER 22-0 Some. Pce\\mtfio.o«.s on VVxe. o s e o? Wve. Tre^o«2.n.c,y Docneun2.1 Introduction2.1.1 An Outline of the Previous Use of the Frequency Domain in Aircraft System Identification Work.Most published accounts of applications of system identification techniques to helicopters have been concerned with time-domain methods using a reduced-order mathematical model representing six degrees-of-freedom rigid-body motion (e.g. MOLUOOl, HALLO01, DUVA001, PADF005). However, interest in the use of the frequency domain for aircraft and helicopter system identification has increased recently as evidenced by published work (e.g. KLEI005, FUKH001, TISC001, BLAC007).In the reference TISC001, transfer function models are used for the modelling of XV-15 tilt-rotor aircraft dynamics, and the identification in the frequency domain is approached using spectral-analysis techniques which are considered briefly in Section 2.2.1. For the identification work carried out in this document, the problem was formulated in terms of a state-space model with stability and control derivatives explicitly represented; it can be shown that an equivalent state- space representation can also be found for a transfer-function model.A frequency-domain maximum-likelihood method of-48-system identification, based on a state-space model representation is developed in KLEI005, with a view to application to aircraft system identification problems. The frequency-domain output-error method developed in this document in Chapter 4 is similar, in many respects to the method developed in KLEI005, in that the same form of cost-function is used, whose frequency-domain residual terms are calculated using a state-space type model. Also, a quasi-Newton method is used in the minimisation of the cost-function in both cases, though some additional useful features for introducing various types of constraints into the identification are available for the software implementation written by the current author. The main differences between the method developed in KLEI005 and the method developed in this document, are asociated with both the form of the frequency-domain state-space model used and the way in which it is represented symbolically.In FUKH001, the authors apply the frequency-domain output-error method, developed in KLEI005, to some helicopter flight data. A modification term for nonperiodic measurement windows, to be included in the state equation, is presented as a means of improving the basic method (this modification term is discussed in Section 2.2.3). The results for the application reported in FUKH001, are presented only in terms of the final frequency-domain fits obtained, without any discussion of the parameter values found. Nevertheless, the work-49-presented by the authors in FUKH001 (and in a similar paper: MARCO02) was very important for the development of the frequency-domain methodology presented in this document. To the knowledge of the current author, no applications of a frequency-domain output-error method to helicopter flight data have been reported, with an accompanying assessment of the estimated stability and control derivatives, other than those given in this document (and by the current author in the referencs BLAC005 and BLAC007).The authors in the paper FUKH001 advocate an identification methdology which uses both an equation- error and a frequency-domain output-error method, with the former used to provide initial estimates to the latter. A final stage involving a time-domain output- error identification of biases and zero-offsets leading on to a time-domain verification of the identified model is also part of the methodology.The three identification stages given above also form the basis of the identification methodology developed in this document. However, it is believed that the additional features described in this document relating to model representation, model-structure determination, computer-implementation and the software facilities developed, constitute a significantly improved methodology for helicopter system identification.-50-2.2 The Use of Freqoencv-Domain Data.2.2.1 The Spectral—Analysis Approach.Some of the early approaches to the extraction of airplane stablility and control parameters were carried out using frequency-domain data, (a brief outline of the earlier work is provided in the reference KLEI005).Frequency-response functions, formulated as rationalpolynomials of jcu, were estimated by least-square minimisation of the equation-errors, (e.g. MARCO01). The frequency-response values for the measured data wereobtained through the use of auto and cross spectral densities. The equation-error term can be written as:FesVCci) - Few) ~ x^yCcS)! SxxCco) “* £ D K (jco)* (2.2"1) HrsV (go) 55 5 xyCcS)/ S x* Cco) (2.2-2)|\Jt and Dk are the coefficients in the numerator and denominator polynomials, and Sxx (co) and Sxy (co) are the calculated auto and cross spectral density functions. Windowing techniques can be used to help reduce the errors introduced as a result of using finite-length records - i.e. leakage of power from different frequencies, or smearing (loss of resolution) (e.g. PAPOOOl, BEND001).-51-The frequency-response errors are:A F Cco) = F^y Cco) - F Cco) (2.2-3)A F C cS) - R* C, AFCco^l + J A.ryj F. A (2.2-4)The least-squares criterion is defined as the minimisation of:Frequency-dependent weights can also be incorporated into the above equation.The majority of texts, when defining frequency- domain estimation, use approaches like the one outlined above (e.g. BENN001), where there is a requirement for the accurate estimation of power-spectral-density functions. The estimated models are in the form of frequency-response functions, characterised by the identified poles and zeros.2.2.2 An Alternative Use of Fourier Transforms.There is,an alternative way of using the Fourier Transform, and frequency-domain data, to that given in Section 2.2.1. Fourier's thereom (e.g. PAPOOOl) states that any periodic function f(t), with period T, consists(2.2-5)-52-of a sequence of impulses, separated from each other in the frequency domain by C0o = Z n / j . The Fourier-series expansion of f(t) is given by:COPCt) = ^Q-n exp C iaoo01) (2.2-6)ns-w 1 Uwhere the coefficients CL are given by: r/zA n - '/t S FCh) expC-jnoot") At (2.2-7)This is a common example of a more general result for expansions using orthogonal functions (e.g. WYLI001), where the length of the interval is T seconds, and defines the fundamental frequency F0 = = l/T.For a periodic function f(t), it can easily be shown that for At = T/N (e.g. PAPOOOl):w - \ _A * e x p (j 2X*mn/ru) (2 .2- 8)where Q n are known as the aliased coefficients, and are related to the CL in (2.2-6), which is defined in (2.2-7), byoOOln = ^ (2.2-9)r= -&0The finite Fourier-series expansion is a name often given to Equation (2.2-8) . From (2.2-8) we see that determining the Fourier series is reduced to the solution of system-53-of N equations; this corresponds to the well known result (PAPOOOl):Taken together, Equations (2.2-8) and (2.2-10) correspond to the well-known finite Fourier-transform pair.It is perhaps more enlightening to rewrite (2.2-8) in a form which explicitly includes the frequencies con , defined below.Folding FrequencyFrom Shannon's sampling thereom, the highest frequency which can be uniquely distinguished is given by half the sampling frequency 1/At (e.g. OTNEOOl); this corresponds to the frequency defined in (2.2-12),and frequencies above this are aliased into a lower value. It can also be shown that:- Vru Z f(maO e.kp (-j /ro) (2 .2-10)(2 .2-11)6Ja = nzTr/rOAt (2 .2-12)Re U a rl - ReCOU-rl (2.2-13)-54-Averaae ValueThe co = o (i.e. C0o) value represents the estimated average, or D.C. level, of the data. From (2.2-10) we have:(M- 1CLo = '/rv) Z PCmat) (2.2-14)rn- oFor the current application, these represent approximate values of the trim levels of the measured variables transformed into the frequency domain; that is, the points about which the small perturbation model is being considered. In addition, constant biases or zero offsets in the measurements, will manifest themselves in this term. The co = o value is excluded from the range of frequencies used in the identification; this reduces the number of parameters that have to be identified using frequency-domain data.Application of Fourier Transforms to State-Soace Model For a state-space model given by:J(ct) .= A X O O + B y c O (2.2-15)we can apply the discrete Fourier transformation to obtain a frequency-domain representation of the model:-55-^ ^ ^Xccon) = AXCcj^) + Bycco^ (2.2-16)where the notation ^ signifies the transform of a time derivative. The model given in (2.2-16) is defined for discrete values of frequency given by (2.2-12).The frequency-domain estimation techniquessubsequently developed, and used in the work covered in this document, use complex-valued frequency-domain quantities, such as those given in (2.1-16) in theestimation algorithms.Errors Due to Finite Measurement Window Length and Non- PeriodicitvWithin the limits of machine accuracy, a set ofmeasured data transformed into the frequency domain using Equation (2.2-10), can be obtained once again in itsoriginal time-domain form using Equation (2.2-8) (i.e.inverse Fourier transformation). It can be seen, however, using Equation (2.2-8) and the expressions in (2.2-13), that we could perform some sort of filtering operation by setting to zero those calculated Fourier coefficients ^which correspond to a particular frequency range. For the removal of frequencies in the range 63^ t0 we havefor the "filtered” signal fq:CrnAb)(2.2-17)OJ-Z j l r * £ x p (j2 rtn m /l\j) - Z dn exp ( j X ^ r^ /ru )-56-By performing an estimation with a restricted frequency range, the effective time-domain response, with which the model is being fitted to in the identification, is given by (2.2-17).In representing a function in the form of a Fourier- series expansion such as Equation (2.2-6) with an infinite number of terms, or as Equation (2.2-8) for a finite number of terms (with aliased coefficients) it was assumed that the function was periodic. However, in the present context we have to note that helicopter responses are transient, and can be considered to be periodic functions with an infinite period. If the transient response of a state has decayed to a value near zero(assuming the response to be that of a stable system) atthe end of the record, then we can consider it to beapproximately a periodic function within the measurementwindow. For the test-input signal, we can ensure that it is zero (i.e. its perturbed value above the trim setting is zero) prior to the application of the input, and at the end of the record; thus it is periodic within the measurement window.For a finite measurement window, the Fourier series approximation assumes that the function duplicates itself at regular intervals, for all times outside the measurement window (i.e. assumes a periodic function). Inorder to overcome this problem, we would require recordsJof infinite length for the responses and the test-input signal, corresponding to the limiting case of functions-57-with infinite period. In reality, we can only use finite-length data records. This is equivalent to multiplying the "true infinite" signal by a box-car function, amounting to a convolution of the true Fourier coefficients with the transform of the box-car function (e.g. BEND001). Moreover, the use of a series representation is an approximation, since an infinite- period function can only be described using a continuous range of frequencies (i.e. an integral expression).The effect of using a finite window - and by implication the assumption of periodicity in the analysis - is apparent when "filtering" operations, as defined by Equation (2.2-17) are performed. This is shown in Figure 2-1 using real flight-data, where the original signal and the "filtered" signal are superimposed for a test-input signal and some measured responses. Uncharacteristicoscillations are present at the beginning and the end of the "filtered" signals in cases where the original unfiltered signals are not periodic within themeasurement window. This corresponds to the discreteequivalent of Gibb's phenomenon (e.g. PAPOOOl) which occurs in analytical work using continuous Fourierintegrals, when a restricted low frequency range is used in describing a piecewise continuous function in the vicinity of a discontinuity: the function begins to oscillate rapidly as the discontinuity is approached from both directions, and the concentration of the ripple increases as the frequency range used increases.-58-The above effects represent the second maindisadvantage of the indicated use of the frequency domain in helicopter system identification: namely, theeffective approximation of the response using a Fourier- series representation, whose coefficients are calculated using a finite-length record. (The other disadvantage was the restriction of the technique to smallperturbation systems represented by a state-space model with delays).Hence, for the identification, it is desirable to have - for a stable system - long data records, andperiodicity in the measurement window, resulting in a representation more characteristic of the real system. However, as was pointed out in Section 1.3.1, obtaining sufficiently-long records,'without the use of a stability augmentation system, is a major problem for helicopter system identification.Least-Squares Sense of Orthoaonal-Function ApproximationIt is a property of the discrete Fourier-series expansion obtained for a given data set, that it represents uniquely the "best” - in the least-squares sense - trigonometric (and hence periodic) approximation for a specified order (ISAA001). The maximum order is determined by the number of data points available, and all other orders of approximation, which are also best in the least-squares sense, correspond to Equation (2.2-17).-59-This is a consequence of the well known, and more general result for orthogonal series approximations, whereby increasing the order of the approximation does not require the re-calculation of any coefficients, only the inclusion of an additional term. Hence, for the use of a given data set in an identification, carried out over a stipulated frequency-range, the effective time-domain response, to which the model is being fitted, is the "best1 periodic approximation in the least-squares sense.(NAG) Library Routines Used to Obtain Fast Fourier TransformsThe identification software developed by the currentauthor, P^es Hve f\)/\Cr (FFT) routine C06FAF(NAG001) to obtain the discrete transforms. The Fourier-transform pair so defined are similar to (2.2-8) and (2.2-10) except for multiplicative factor of l/\/N; they are:Equivalence of Time and Freauency-Domain ApproachesQn - '//n Z^PCmatO^xp rnn/l\)) (2.2-19)Using all the available frequency-domain points in an identification is equivalent to using the original-60-time- domain sequence. The cost functions which are minimised for the identification techniques (i.e. equation error or output error) can be represented as functions involving sums of squares. Parseval's theorem (PAP0001) states that the following relationship exists between the squared frequency-domain and squared time-domain quantities (it is written in terms of the fourier transform pairs defined in 2.2-8 and 2.2-10):The equivalence between the time and frequency domain approaches no longer holds if a restricted frequency range is used. This is, however, one advantage in using the frequency-domain: the frequency range over which the cost function is minimised can be varied as required, within the estimation algorithm, and there is only the requirement for the storage of the original data set. Time-domain estimation for a variety of different frequency ranges would require the production and storage of at least one additional data set.Some Further Advantages of Using the Frequency DomainThe ease with which the frequency range used in the identification can be restricted was cited above as one(2 .2-20)-61-advantage in using the frequency-domain; this is particularly advantageous for helicopter system identification when models for rigid-body dynamics, valid for low frequencies, are to be obtained. The exclusion of high frequencies, representing unmodelled effects, is also beneficial from the point of view of data reduction.Another advantage is that transformation of the model equations into the frequency-domain has some important implications for the estimation algorithms themselves: operations in the time-domain such asdifferentiation are replaced by multiplication; and equations, which in the time domain require numerical integrations over the entire time record for each iteration (such as the equations for the partial derivatives of the cost function with respect to the unknown parameters), are replaced by expressions which can be solved by algebraic manipulation. This is brought out in Section 4.2.2, where the expressions are developed.In addition, by formulating the estimation problem in the frequency-domain, it is easy to include time delays in the set of parameters to be estimated; this is also brought out in Section 4.2.2.2.2.3 An Expression for the Discrete Fourier Transform of a Time DerivativeFrom the state-space model given in Equation-62-(2.2-15), we can write the equation shown in (2.2-21) involving the finite Fourier integral, for an interval of T seconds. The reason for displacing the integration limits from their usual values by half a sampling increment, will be indicated shortlyt- T-&PexpC-joJtf) d\r = A 3 X Ct) dXp dt +(2 .2-21)B 5 UCOexpC-^oOcit■^cThe integral expressions in Equation (2.2-21) are known as finite Fourier integrals (BEND001), and can be defined for a continuous range of frequency values c o .By integrating the left-hand side of Equation (2.2-21) by parts we can obtain an expression for thefinite Fourier integral of a quantity differentiated withcirespect to time, in terms of the integral of the original quantity. For periodic functions, the result is well known. In the references FUKH001 and MARCO02, the authors present an expression for the discrete Fourier transform of a quantity that is differentiated with respect to time which involves both the discrete Fourier transform of the original quantity and some additional points pot used in the transformation. An important feature of the result is that it does not assume the measurement data window to be periodic (which is the usual case for helicopter flight data).The justification for the non-periodic measurement --63-window result given in the references FUKH001 and MARC002 relies on using an approximation between a discrete Fourier transform and a Fourier integral defined between the limits (-At/2, .T-At/2) as its starting point. No mention is made of the numerical integration technique implied in making the approximation (given in 2.2-24), or why the integration limits are so defined. A justification for the approximation, which is the starting point of the derivation, will be provided by the - current author in due course. Firstly, however, let us consider the derivation of the result itself. For thediscrete Fourier transform defined by (2.2-10) we have:£ ru - \XCccO - '/r\) X J( CmaO exp ( ( 2.2-22)we will write:T- AtA ^X(oO = ' / r ) X c O e* pC - \c ont ) d t (2.2-23)-AtXand similarly we can write:f- AtA 3*X(co„') * V t ) X ( t ) e.xp (-'\cor>t') A t (2.2-24)-At 1 JZIntegrating the right-hand side of (2.2-24) by parts, and using the fact that we are considering only discretevalues of frequency: , (defined in 2.2-12) we thushave:-64-\V • A/ ✓ j cQq A tA C ly > U are the offsets from the centre of gravity of the speed, flank and incidence measurement devices respectively, and where Ue , Ve. ,-69-We represent the trim values of the velocity components in the body-fixed axes along the X, Y, Z directions respectivelyThe (time-domain) linearised measurement equations, obtained using a first-order Taylor expansion about the trim, are as follows for a typical set of flight measurements, where measurement bias terms:( bv, bp, btf, bp be) bp )b^ > br and measurement noise terms ( Hv , ri/3, 0*0 r\> He, np> rty > Dr )/ are included.V.AmOiftil m0mn) , etc. we have:AX(cOn') - COn XCcOn') + C^On G CcOrD +■ Cj Coin) (2.3~5)Awhere X(^n) is the discrete Fourier transform of theAsecond time derivative of X(t), and (r (co^ ) is analogous toA(b(cc>(\) given by (2.2-8) but evaluated at the corresponding time derivatives. Hence we may write the right-hand side of (2.3-4) as:X'Cwj + 'c(ico0 XCoO + &CaS) + -cV2i +A -A .J 6Jr\ Gr CcJn') + G C W n ) ) + ••• =( l + tjOJn - ^ T ^ c O n / * ! + .•OX'CcoO + ( t + T j© n /Z l * . .0 Orton) *CrV^l ...) g-Ccoo + o c t 3)(2.3-6)Hence we have thatn>-\»/ro Z X ( U t + r ) e J“"Lat = er“n- »/w X XCiat) e iu"Ut*■" o a 5. -©+ G^CoO(2.3-7)Awhere Gt(con) is a complex (both in form and in value) modification term. In the application of the time delay to a control, it is reasonable to assume that the (perturbed) control state is zero (corresponding to the trim position) prior to the application of a test input-80-signal, and is zero at the end of a time recordcorresponding to the free response, in both cases for a time span of length at least IT seconds: hence the signalAis periodic and Cxr(con) is zero since all the higherderivatives are zero. For time delays in the controls (and for the test input signals used in practice) and for delays in other quantities periodic in the measurement window, we can use the result given by (2.3-7), with theAmodification term Gc(coO set equal to zero.Changes to the Model EquationsBy formulating the estimation problem in thefrequency domain, the inclusion of time delays (in the controls U and measurements Z ) results in a relatively simple modification to the frequency-domain stateequations given by (2.2-37) and (2.2-38). Additional diagonal matrices of the form:— d = C os ( Ti. coO - j j l - v/Ue -15/Ue. O WCOo 0 1 O l it ')e M(0 o o o 1 .©CO4*vorf4(of F (JcFF _6©fp(2.3-12)The terms in the above equation were explained in Section 2.2.4.A number of the parameters in these equations are known to be insignificant (measures of parameter significance as explained in Section 3.2.4 were used) and were excluded from the estimation process. Initial parameter estimates ' for the output-error method were-90-obtained from equation-error fits for each of the rows. Description of Results ObtainedThree different cases for the estimation were considered:(1) estimation without a time delay in the longitudinal- cyclic control (i.e. the delay in the above model is effectively fixed at a zero value), using the flight-data sampling time of 0.015625 seconds;(2) estimation with a time delay in the longitudinal-cyclic control using the flight-data sampling time of 0.015625 seconds; and(3) estimation without a time delay, but using anincreased sampling time of (6 x 0.015625 =) 0.09375seconds. The final estimates obtained on convergence of the estimation algorithm, are given in Table 2-2, along with theoretical values from the HELISTAB program.Examination of the results for cases (1) and (2) reveals that the inclusion of the time delay resulted in a smaller cost-function value (the frequency-domain output- error cost function is derived in Section 4.2.1) at convergence, for the sampling interval of 0.015625 seconds: that is to say a better fit with the observedflight data was obtained for case (2). The estimate ofthe important pitch-rate-damping parameter without the-91-delay is very much underestimated in comparison to theory, but with time delay estimated in the model, the agreement is excellent. The underestimation of this parameter in helicopter parameter identification studies is reported by some authors (e.g. PADF002), where the assumption of quasi-static rotor dynamics in the simplified six-degrees-of-freedom models is cited as the most probable cause. The inclusion, and estimation, of the time delay in the longitudinal-cyclic control appears to have made a considerable improvement to the situation. For the case with the time delay, excellent estimates are obtained for the cross-coupling parameters: and Mp..The control sensitivity M^1S agrees very well with theory for the inclusion of the delay. The time delay itself is estimated with confidence; its estimated value is of the same order of magnitude as the time constant for longitudinal and lateral cyclic flapping. For the derivative Mw, it is seen that the theoretical HELISTAB prediction represents a more stable aircraft that is suggested by the results; the change in sign obtained for this parameter, when the time delay is included in the model, indicates that there may be some strong correlation between and Mw. The normal-forcederivative: Zw, is estimated with a small error-bound,and is unaffected by the time delay.Looking now at case (3) , it can be seen that by increasing the sampling interval to 0.0 9375 seconds, the most significant features of the results obtained in case-92-(2) are repeated: namely, the significant improvements in the estimates of and M^1S over those obtained for case (1); and the reversal in sign for the estimate- of M w. Improvements in the quasi-static estimates of someparameters, such as pitch-rate damping and controlsensitivity, obtained by increasing the sampling interval, were demonstrated in the reference PADF002 using simulated data and time-domain equation-error estimation. These observations have now been confirmed using real flight-data, and a more advanced frequency-domain output-error technique. This improvement is because, as was indicated in Section 2.3.3, thehigher-order rotor transient effects are effectivelyfiltered out for large sampling-intervals. It was also pointed out, however, that this is not the approach advocated here: the- use of a time delay with a smallsampling-interval is preferred - the superiority of this approach will be made apparent in the rest of this section.The cost-function value at convergence for case (3)is smaller than case (2) ; however, because the sampling interval in case (3) is six times that in case (2), thetwo values are not directly commensurable. Examination ofTable 2-2 reveals that the cross-coupling derivatives: Mp and Mp, have a less satisfactory agreement with theory than case (2). In fact, plots of the frequency-domain fits for the three cases, shown in Figures 2-5a), b) and c), reveal that case (2) (with the time delay in the-93-control and the small sampling interval) gives the best frequency-domain fit; this is especially true for the pitch rate.Using the parameters identified in the frequency- domain, a time-domain output-error procedure can then be used to estimate initial state conditions, and measurement zero offsets - quantities not identifiable in the frequency- domain because of the exclusion of the zero frequency to obtain the time-domain fits of the estimated model. This was done using the program: OFBIT, developed by the current author; the theoretical basis of the method, and the computer implementation, is explained in detail in Chapter 5 and in the reference BLAC006. The time-domain fits obtained for three cases above, are shown in Figures 2-6 a), b) and c). Case (2) is again the most satisfactory, having the best time-domain reconstruction. The inclusion of the time delay results in a much tighter fit between the observed and predicted time-domain responses, especially for the pitch rate: ^ . A comparison of the pitch-rate fits in Figures 2-6 a) and b), shows clearly the improvement resulting from the inclusion, and estimation of, the time delay in the longitudinal-cyclic control. The lateral variables used as 'pseudo controls' in the extended control vector, along with the longitudinal- cyclic control are shown in Figure 2-6 d).-94-Sum m ary o f C h a p t e r 2The use of the the discrete Fourier transform as a means of representing a linearised state-space model of a helicopter in the frequency-domain, was presented. The model was extended to include time delays in themeasurements and controls; and the ease with which the time delays can be incorporated into the frequency-domain model was emphasized.In addition to the state-space equation,consideration was given to the frequency-domainrepresentation of the measurement system: it wasexplained that a linearised representation was necessary because of the practical difficulties encountered in applying a linear transformation, such as the discrete Fourier transformation, to non-linear equations in general. It was pointed out that a drawback of thefrequency-domain approach was its limitation tosmall-perturbation models (non-linear only in the sense that time delays are permissible). Another drawback of the frequency-domain approach was explained: namely that errors are introduced into the identification as a result of representing the- response as a Fourier series whose coefficients are calculated using a finite-length record.However, the potential advantages for helicopter system identification, gained by formulating the problem in the frequency domain are manifold. Some examples of the advantages were given in both this and the previous-95-chapter, and are brought out in detail, or demonstrated, in the remaining chapters of this document. These include:1) Noisy flight data can be used directly in the system identification. There is no pressing requirement for pre-estimation techniques, such as Extended-Kalman-Filter state estimation, in order to remove high-frequency noise, or to construct unmeasured quantities, such as the time derivatives of states; thus cutting down considerably on the time required to perform an identification.2) The ability to easily introduce selectivity in the frequencies used in the identification is useful for models valid for a specified frequency range; however, the resultant reduction in the number of data is also computationally beneficial. For the helicopter flight data used by the current author, there was found to be typically a 1:50 reduction in the amount of data actually used in the frequency-domain identification algorithm, in comparison to the original time-domain data. This kind of reduction in the amount of data augurs well for identification techniques based on cost functions involving summations. It should also be pointed out that obtaining the Fourier transforms of the measured time- domain data can be done very quickly using standard FFT routines, such as the NAG routine C06FAF, mentioned in Section 2.2.2. The availability of frequency-domain-96-records provides a useful indication of the degree of excitation of the system at frequencies of interest, and there is no requirement for the creation of new data sets each time the frequency range used in the identification changes.3) By excluding the zero frequency value from the range used in the identification, there is a reduction in the number of parameters that have to be estimated using the frequency domain. Zero-offset terms are decoupled from the rest of the model. These can be estimated at a later stage, along with initial state conditions, using the time domain.4) Expressing the identification problem in the frequency-domain has some important practical benefits for the estimation algorithms, and their computer implementation. The problem of obtaining the estimated model output, and other quantities, becomes algebraic; in the time-domain, numerical integration is required for the corresponding operations.Finally in this chapter, the usefulness of a new approach (i.e. the incorporation and identification of a time delay in the control used) in overcoming some of the problems associated wiVh quasi-static derivatives was demonstrated using both simulated data and real flight-data. For real flight-data, the parameter estimates obtained compared favourably, in most cases, with the theoretical values provided by the HELISTAB-97-program, when a time delay was identified in the longitudinal cyclic input. The improvement was also apparent in the time-domain reconstructions of the models identified in the frequency-domain.Results from simulated data indicated that in -some cases, there would be no requirement for the inclusion of a time delay in the collective input channel because of its association with high-frequency coning effects, and because of the sampling interval used in the flight data; results obtained by the current author using real flight- data, have confirmed-this. The HELISTAB program does not include the modelling of tail-rotor dynamics, however, results for pedal inputs using real flight-data, and presented in Section 4.3.2, show the inclusion of a time delay in the pedal control, to be an important modelling requirement.-98- CHAPTER 33.0 Initial Stages In The Identification Methodology, And A Consideration Of Some Relevant Identification Techniques.3.1 Introduction And Overview Of The Methodology.The three basic elements of the identification methodology are developed and discussed separately in each of the remaining chapters of this document, culminating in a summary of the salient points for each stage in the final Sections of Chapter 5. These three elements are: a frequency-domain equation-error method;a frequency-domain output-error method; and a time-domain output-error method (see Figure 5-1). The ordering of the chapters represents the natural sequence of application for each of these stages, where the results from one stage are used to initiate the next stage.Initial parameter estimates are obtained at the equation-error stage; these estimates, which are biased in the presence of noise on the independent variables used in the equations, may then be used as initial guesses for the iterative frequency-domain output-error method, which is capable of producing unbiased estimates in the absence of process (or model) noise. Both the equation-error and output-error identifications are carried out over a restricted frequency range appropriate-99-to the assumed form of model, namely one which describes the rigid-body motions of a helicopter. Convergence of the frequency-domain output-error method is assisted by having good initial estimates of the unknown parameters; the singular-value-decomposition implementation of the equation-error method, presented in this Chapter, is seen as a means of improving the estimates obtained at the equation-error stage. For the frequency-domain output- error method implemented here, options exist for investigating the effects of different model structures on the estimates. Results from the equation-error stage can also assist in isolating weak, or insignificant, parameters to be fixed, or excluded, from the state-space model used in the output-error identification.The final stage of the identification methodology involves the estimation of initial conditions and zero- offsets which are uncoupled from the frequency-domain identification, but which are necessary for a final time- domain verification of the model. Stability and control derivatives, and time delays, estimated at the frequency- domain output-error stage (together with those elements of the model fixed for that identification) are fixed for the time-domain output-error identification.There are some issues which relate to both the frequency-domain equation-error and output-error methods, such as the selection of an appropriate frequency range, and the selection of an appropriate length of. time-domain record for transformation. These are discussed in this-100-Chapter since it is at the equation-error stage that they are first encountered. In addition, the comments made about attempts to apparently increase the frequency resolution through the use of zero pads in this Chapter, also apply to their use at the output-error stage.The equation-error method is a widely used approach to system identification. This is so, mainly because of its simplicity and easy application to linear or nonlinear models. In the reference KLEI001, the author pays tribute to the equation-error method: "Despite all these degradations in the accuracy of the estimates resulting from real flight-data [because of violations of the assumed noise properties], the equation-error method is often used, sometimes with very consistent results in comparison with more sophisticated techniques". Some previously reported applications of the equation-error technique are discussed in appropriate Sections of this Chapter.3.2 The Freoruency-Domain Equation-Error Stage In The Identification Methodology.3.2.1 The Basic Method And Some Comments On The Frequency-Domain Formulation.The equation-error technique is based on the least- squares principle, and is well-known for its application to curve fitting or regression analysis. With thisapproach, the solution aims to minimise the sum of squares of deviations between measured data points and corresponding points obtained from the solution, where the measured data is to be represented by a functional relationship, or smooth curve.For the estimation of parameters of a dynamical system with state-space form (as given by Equation 2.2- 15), it is necessary to consider one row of the matrix- equation at a time for the equation-error method. Consider an element of the state time-derivative vector as the dependent variable, with unknown parameters CUj corresponding to the ith row of state matrix A , and h.* corresponding to the Lth row of the control dispersion matrix B . It is assumed that the state X and control vector U are measured without error, whereas the measured values Y of the dependent variable Xu are corrupted with zero-mean uncorrelated Gaussian measurement noise such thatVCO - XuCr) * £cr)(3.2-1)— Q . X . C r ) + OLi.2 , X2 .C1O Xrj CfOb L| U,Cr) + b^mU^Cr) + £.cr) (3.2-2)for an nxn state matrix A and nxm control dispersion matrix B. The error term B.Cr) in (3.2-1) can represent process (i.e. modelling) noise on the model, in addition to measurement noise on the dependent variable Xt_-102-assuming, of course, that a zero-mean uncorrelatedGaussian representation is appropriate. The least-squares criterion in defined as the minimisation of thecost function: uT ( © ) = ' /Z. Z C'/Cr) - OUiXiCr) - C 3m y c ^ )]") (3.2-21)for a frequency range, given by r=l, 2..q, used in theidentification. By similarly constructing vectors forA/ A/ A the independent variables: Xi (r) , X;?(r) 0m(r) i-n(3.2-16) and (3.2-17), we can perform the identificationusing both the real and imaginary parts, and are thusminimising the cost function:The cost function defined in (3.2-22) can only have zero contribution at a given frequency r, when (3.2-20) holds.Instrumental Variable Technique.The requirement for a zero-mean uncorrelated Gaussian error-term in the equation-error method (in addition to error-free measurements of the independent variables) is pointed out in references to the equation- error method (e.g. KLEI001); this was expressed by Equations (3.2-9) and (3.2-10) as a necessary condition for obtaining unbiased estimates (i.e. estimates without systematic errors). In the context of time-domain estimation, a technique which attempts to overcome the difficulties of correlated noise is. through the use of an(3.2-22)-109-instrumental-variable matrix (BOOMOOl, KLEI001) ; theinstrumental-variable matrix Zxv' premultiplies all the original quantities in the equation-error estimation problem given by Equation (3.2-4) to give, by analogy with (3.2-9):E t © , 2 1 = © + EC.(Czlv X)T(ZI ,X ))"(2 I ,X')T(Z x v t ) l(3.2-23)where the properties of the instrumental-variable matrix Zxvr are such that:l(2rvX)TCZivX)\ > O (3.2-24)E t Z i v t ) = O (3.2-25)For a positive-definite Hermitian matrix (e.g. STOEOOl) as defined by (3.2-24), and for an uncorrelated effective noise source, as given by (3.2-25), we have that the second term on the right-hand side of Equation (3.2-23) has an expected value of zero, resulting in unbiased estimates of 0.An important property of equation-error estimates (in the case of uncorrelated noise) and of instrumental- variable estimates is that of consistency. Basically, a consistent estimator in the time-domain is one whose covariance matrix (i.e. the parameter error bounds) tends to the null matrix (i.e. zero error bounds) as the number of points used in the identification increases. This-110-will be demonstrated for parameter estimates obtained using real flight data in Section 3.3.3. The property of consistency is defined by:In the reference KLEI001, the consistency of equation- error estimates (in the case of uncorrelated and stationary noise) is proven; the consistency of instrumental-variable estimates is also stated, but not proven, though it can be easily shown that the proof of consistency follows in a similar fashion.With regard to the requirements given by (3.2-24) and (3.2-25), it is suggested in the reference KLEI001, that the instrumental variables (i.e. the rows of "Zxv)/ could be obtained as the output from a Kalman filter with known input variables, and approximate values for the stability and control derivatives.Use Of A Restricted Frequency Range.It was observed by the current author that the solution to the frequency-domain equation-error estimation problem, carried out over a restricted frequency range, could be written in a form analogous to the instrumental-variable formulation for time-domain estimation, given by Equation (3.2-23). Unbiasedestimates are obtained for requirements similar to thoseO (3.2-26).-111-given by (3.2-24) and (3.2-25). In addition, the estimates are consistent for an increase in the number of frequency-domain points within the frequency range valid for the model; this corresponds to (3.2-26) .To show that the above is true, consider first that the discrete Fourier transformation of a time-domain quantity which is represented as a column vector, isachieved through pre-multiplying it by a complex-valued_?/ / matrix - /p^ , whose elements - Zp* are:- '/V w exp(.-j2tf (.L-OCvt-O/ro) (3.2-27)LjK . ... NNow the method devised by the current author for using what are essentially .complex-valued quantities in the estimation, was to alternate the real and imaginary parts of a column vector, as shown in (3.2-21). If we use a real-valued matrix Z f*>, where alternative rows correspond respectively to the real and imaginary parts of a discrete Fourier transform, defined by:= ' / / N Cos 2 * Ck - O / n(3.2-28)= "'/v/iv S m 2 K O . - 0 ( . K - 0 / wL, K = 1 , 2 ^ ••• IV)then the frequency-domain equation-error problem can be written in terms of the matrix Z fk and the original time- domain quantities. Pre-multiplying Equation (3.2-4) by Z fr,/ we have:-112-= Z f^ X S + ZfPv £ (3.2-29)By excluding the rows in (3.2-28) corresponding to 1=1, we will exclude the zero frequency from the identification; in addition, if we define Zpft such that it only has rows corresponding to Z ^ i ^ 4 , then weare performing the identification on a restricted range of lower frequencies.By analogy with Equation (3.2-23), we have that:assuming the frequency range implied by the rows of Z f ^ is the appropriate one for the identification of the model.Using the definition of Z fk given by (3.2-28), constructed for a reduced frequency range, we can see that provided the sums of the squared magnitudes of the transforms (i.e. the power within the stipulated frequency range - see Parseval's theorem discussed in Section 2.2.2) are 'sufficiently' large for each of the independent variables in Equation (3.2-29), then no terms having a small magnitude will occur in the leading diagonal of the real-valued matrix product, whose diagonal terms are defined below:- © - E E.cc zF(kx)T(zfr,x))'(Zfrx)(zFR£)](3.2-30)(3.2-31)-113-/»where are elements of the matrix defined by:- X (3.2-32)In addition, provided the frequency-domain responses of the independent variables are 'sufficiently'uncorrelated, we will then have a condition analogous to (3.2-24) :( Z f* X) (ZfpvX)I > 0 (3.2-33)The above discussion relating to the expression given in (3.2-33) is in fact centred upon the condition number of the frequency-domain information matrix (i.e. the matrix whose determinant is considered in 3.2-33), which is affected by both the sensitivities of, and correlations between, the independent variables; these concepts are considered in more detail in Section 3.2.5.In order to finally complete the analogy between the time-domain instrumental variable technique and the frequency-domain equation-error method (using a restricted frequency range), consider the matrix For helicopter flight' test data, the original time-domain noise source E, for a six-degrees-of-freedom model that includes fuselage states only, will not, in general, be adequately modelled as a zero-mean uncorrelated Gaussian disturbance. Consequently, biased parameter estimates-114-and poor fits will be obtained when the original time- domain data, or equivalently the entire frequency range available, is used in the identification; this is shown graphically in Section 3.3.2 using real flight data. However, if the restricted frequency range used in the identification of a six-degrees-of-freedom model is such that the frequency-domain noise source is adequatelymodelled as a zero-mean uncorrelated noise source, then we will have a property for the real-valued matrix ,identical to that for the time-domain instrumental- variable matrix 7 . ^ , given by (3.2-25), that is:E E Z FR£ l = O (3.2-34)Hence we have, for an appropriate frequency range implied by the rows of ZFR , and using (3.2-33) and (3.2-34), that:E C (C ZF* x)T( zFR x))"(zFR X) a * t)] = EE((zFRx)T(zFRx))"UFRX)] • ELzFRt) - O (3-2_35)From Equation (3.2-30) we can see that this results inAunbiased parameter estimates (i.e.E u O FJ = § ) .Despite the violation of the assumption of no noise for the identification op real flight data, it will be demonstrated in later Sections of this Chapter using real flight data, that useful results can still be obtained using the frequency-domain equation-error approach over a restricted frequency range.-115-Direct Transformation Of Flight Data.When the states to be used in the identification and the flight-test measurements are related linearly, we have:R e L f C c o y } = H_ ~ tl a (3.2-36)I m U Z C c o V ] - H X m U XCco)lThe real and imaginary parts of the transformed quantities can be arranged in the same way as the vector defined in (3.2-21), and a single matrix-linear Equation solved for the relatively few frequency-domain points used in helicopter system identification, in order toAobtain X (w) . The same approach is used to obtain the transforms of time derivatives.3.2.2 Subset Regression.The equation-error approach to helicopter and fixed-wing system identification problems is mostly discussed in the available literature in the context of a subset (or stepwise) regression procedure, using time-domain data (e.g. PADF005 [current author is a co-author], KLEI002, KLEI003, DUVA001).The subset regression procedure allows for some-116-investigation in determining the 'best' model structure from a pool of possible model structures; that is, the 'best' subset of independent variables from a set of variables specified by the analyst. The procedure applies least-squares fits in a sequence of steps, where at each step a variable is added to, or removed from the model, on the basis of objective statistical tests. The process terminates when, according to some criteria specified by the analyst, the 'best' structure is obtained, and no further improvement can be made either by the addition, or the removal, of a variable from the Equation. The subset-regression procedure is described in the reference KLEI002, and is outlined below.From a pool of possible independent variables to be included in the Equation, the variable chosen at each stage of the procedure, for entry into the Equation, will have the highest partial correlation coefficient defined by:is subtracted from the dependent variable. At the start of the algorithm this will correspond to the dependent variable to which the model is being fitted - i.e. p=o in (3.2-28). For stage p+1 we have:rJV, = f t x * c l > - x j U i y V o - (3.2-37)\/* .where / is the residual left when the current model fity c o — XxCc} £>

> , fix... ftp , finally using the residual defined as:Xj = — y8p Vp £r)(3.2-39)At every stage of the procedure, variables which are already incorporated into the model, and the variable selected for possible inclusion in the model (on the basis of its partial correlation) are subjected to the partial-F test. The partial-F statistic is defined in (3.2-40):Fp = t t l / sVtfO (3.2-40)where (XL i s a parameter estimate and ^ ( C i c ) is theestimated variance of the estimate O i t. This statistic isevaluated for each of the variables under consideration,and is compared with a pre-defined value set by theanalyst. If it is less than the pre-defined value, then the variableis removed from, or prevented from entering into, the model. This amounts to performing a hypothesis test using the appropriate F-distribution (e.g. MOODOOl), where the null hypothesis postulates that the true value of the parameter is zero (i.e. the variable is not-118-present in the model). The procedure continues until no more variables can be admitted to, or removed from, the model Equation.The partial-F statistic is an indicator of the relative accuracy of a parameter estimate, since it is defined as the square of the ratio of estimated parameter value to estimated standard error.The squared multiple-correlation coefficient - K provides a measure of the accuracy of the fit (the fractional portion of the dependent variable which can be attributed to linear effects), and is defined by:Rl 4 i < y c o - y u ) f (3.2-41)fU2 C YCC) - y c o fuAwhere /CC) are the model predictions and Yc'C) the measured responses.The total F-ratio provides a measure of the confidence that can be attributed to the fit, and is defined by:Ft™ - -R / ( P - 0 (3.2-42)( l - R ' ) / ( W - P )where N is the number of data points used in the identification and p is the number of parameters in the linear model.Whilst the tests for deciding the 'best' model are objective, a subjective component of the process of model-structure determination exists; this is because the-119-selection of the 'best' model is made from a finite set of possibilities specified by the analyst. In addition, most applications of the method allow the analyst to force terms into the model, regardless of the fact that they would not be selected by the automatic process. This feature is present in the optimal-subset-regression program OSR, developed at NASA Ames specifically for helicopter system identification (HALLO02). In the reference KLEI002, a procedure termed "modified stepwise regression" is described; here the aerodynamic force and moment coefficients are described using multivariable polynomials in response and control variables (i.e. a Taylor-series expansion about the initial steady-state flight condition, including terms higher than first order), and in applying the stepwise regression procedure, the significant linear terms are forced into and kept in the Equation, with the testing for entry or removal carried out only on the higher-order terms. It should be pointed out that some analysts prefer not to use the statistical approach for the selection of significant non-linear terms; in the reference ROSSOOl the author points out that some flight dynamicists have reservations in 'leaving out the physics' from the identification.Results obtained in the application of the subset regression procedure to real flight-data, using the time domain, are presented in the reference PADF005 [current author is a co-author]. In that paper, attention is-120-focussed on the rolling and yawing moment Equations; it was observed that the 'best' parameter estimates - in terms of their physical realism, and accuracy indicated by individual partial-F statistics - were obtained for the primary stiffness (N/, Lv), primary damping (Nr, L p), and control derivatives (NAp , L^p) when each of the Equations contained only the three indicated terms. When further terms were brought into the Equation by the subset regression procedure, a deterioration in the physical realism, and a collapse in the value of the partial-F statistics was observed. As would be expected for an increasing number of parameters, the squared correlation coefficient continued to increase for the higher-order, but less satisfactory models.The reduced-order models corresponding to the first three steps of the subset regression procedure were judged to be the 'best' representative models. The concept of obtaining the best representative model based on the choice of one particular model structure from a set of possible structures, will be discussed in Section 3.2.5, using an approach known as singular-value decomposition.3.2.3 The Method Of Successive Residuals.In addition to the subset-regression method, another approach to the use of the equation-error technique in helicopter system identification, called "the method ofsuccessive residuals'1, has been suggested in the references - DUVA001 and WANG001. This method involves partitioning the fuselage six-degrees-of-freedom model into longitudinal, lateral, and cross-coupling subsets, shown here:fx UotvJOlYAn A,a\ ^21 Asa!v \A UorvJOX'B., B,z\+uf\r / B*, BU ^X -^ococri y uf\r /(3.2-43)It is proposed that the partition matrices: Au and Bn ,are best identified by applying the equation-error technique to data generated by a longitudinal input, and neglecting the effects of the lateral states and controls. Once these matrices have been identified, a constraint Equation is constructed as:¥ r Xu> N& A»\ XLotvJO- - B u U LOrOCr (3.2-44)The new dependent variable Y is generated from a lateral- input manoeuvre, and is the residual due to the unidentified cross-coupling terms ( Aw and Bw were identified at the previous longitudinal-input stage). The cross-coupling terms are then estimated by applying the equation-error technique to the data-set generatedfrom the lateral input,using the Equation:y = a, xf} LftT Bia U (3.2-45)The cross-coupling matrices A 12 and are thus-122-estimated in such a way that they are consistent with the An and Bn matrices obtained from the longitudinal manoeuvre. A similar approach is used in identifying the lateral force and moment matrices.It is pointed out in the reference DUVA001, that in theory, the method could be used to identify rotor derivatives corresponding to fast modes and fuselage derivatives corresponding to slow modes independently,, whilst combining them- in a consistent fashion. The frequency content of each manoeuvre would be separated into a high and low frequency portion, and a process similar to that outlined above for the fuselage derivatives alone, could be carried out.In the reference DUVA001, some results are presented where "the method of successive residuals" is applied to time-domain data generated for a linear simulation. A comparison was made of the estimates of the crosscoupling derivates obtained from this approach, with those obtained by identifying uncoupled and cross-coupled terms simultaneously on a data set containing longitudinal and lateral manoeuvres. This data set involved one set of time-history records placed directly after another, and the stacking of data in this way is permissible b e causet he regression looks only for correlations, and is unaffected by discontinuities. The estimates of the cross-coupling terms were found to be superior when "the method of successive residuals" was used.The use of combined data sets from several manoeuvres has been a commonly used, and advocated, approach for helicopter system identification (e.g. KALE001, MOLUOOl); however, as is pointed out in DUVA001: the problem with this approach is that some states are poorly excited in one manoeuvre, and well excited in another, so that combining the data tends to degrade the identification, in comparison to that obtained using only the well-excited data."The method of successive residuals" was investigated by the current author using simulated linear data obtained from the HELISTAB program. Noise was added to the dependent variables (i.e. the state time derivatives) to make the simulation more realistic. A variation on the originally proposed method was also tried, where for example with a longitudinal input, the lateral states were not assumed to be negligible, and all parameters known to non-zero were included in the initial stage of the identification procedure. The crosscoupling estimates obtained at the initial stage were then effectively disgarded; the intention being to estimate these at the second stage, using a lateral test input.Results obtained from this investigation are presented in the reference BLAC003. It was also concluded by the current author, that "the method of successive residuals" was better in estimating the crosscoupling derivatives than the use of stacked data.-124-However, in some cases the variation on the original method (as described before) gave the best estimates, whilst in others the originally proposed method gave the best estimates. When there was a substantial difference in accuracy (i.e. when compared to the true values) this was seen to be reflected in the relative values of the partial-F statistics.The main results obtained from this study, which were carried forward into the identification methodology developed in this document - remembering that the equation-error technique is seen more in terms of providing good initial estimates for use in the output- error method, than in being an end in itself - were as follows:Firstly, single-run identification, where only one control at a time is used in applying the test input signal, is not an unreasonable approach to the problem of helicopter system identification. Data sets at exactly the same flight condition may not be available to apply the multirun approach, and also as was pointed out earlier, the use of such data often leads to the deterioration of some estimates in comparison to the single-run case.Secondly, provision for the effects of lateral states for models identified using longitudinal inputs (and vice versa) should be made. This is accomplished in the implementation of the equation-error technique described in Section 3.2.5, and demonstrated using real-125-flight-data in the sub-Sections corresponding to 3.3. In the application of the frequency-domain output-error method, where for example a longitudinal test input is used, the lateral states are included in an extended control vector; this was shown for the example given in Section 2.3.3, and the corresponding case for the application of a lateral test input is shown in Section 4.3.2, where the cross-coupling derivatives are estimated as if they were elements of the control dispersion matrix.And finally, the fact that some derivatives were best estimated with one model structure, and some with another (e.g. the inclusion or exclusion of crosscoupling terms) emphasized the importance of investigating the effect of model structure, and of determining a suitable one. The method for the determination of a suitable model structure, proposed in this document for the equation-error technique, has similarities with some of the features associated with the two applications of the equation-error method described in this Section (method of successive residuals), and in the previous Section (subset regression).3-2.4 Parameter Significance Values.There is a requirement to determine which parameters are contributing significantly towards the helicopter-126-motion, and which parameters are not, and can therefore be excluded from the identification, where the identification technique in mind is the output-error method. We would expect that parameters which only have a minor influence on the vehicle motion would not be estimated with great confidence. An indication of whether a parameter is weak or strong can be provided by the equation-error stage.In the reference KALE001, the author defines measures of significance for individual parameters, based on the absolute integral of the time response over the time span used in the identification multiplied by the estimated parameter value, and normalised by the absolute integral of the dependent variable. Any other 'vector norm' which gives some indication of the size of the response, such as the root-mean square, could equally well be used instead. Taking the example of the pitching-moment Equation in the time domain, this leads to measures such as:1 H A *) \ uco \ I 5 \ ^CO\ 6 . Y(3.2-46)I Mwl 5 I W(b") u t / ) i^ooidLtThe same concept was extended into the frequency domainby the current author. Here, the integration is carried out over the range of frequencies used in the identification, and the magnitudes of the Fourier transform are used instead of the magnitudes of the time responses. These frequency-domain measures of parameter-127-significance take the form:1 M u \ *) I Occo \^ cico / S ^ | c c j ) \ dco(3.2-47)1 M w I 5 ^ ^ / $" I ^ ^^Parameter significance values are demonstrated in Section3.3.5 using real flight-data.3.2.5 An Implementation Of The Equation—Error MethodUsing Simgular-Value Decomposition.Ill Conditioning Of The Information Matrix.We can see from Equation (3.2-7) that the parameterAestimates 0 , for the equation-error method, are obtained formally by solving a linear matrix Equation; that is by inverting the square matrix ( X V ) (i.e. theequation-error information matrix - though strictly speaking t i^e maximum-likelihood informationmatrix, whose inverse gives the parameter covariance matrix). From the expression for the parameter covariance matrix given in (3.2-12), we can see that in order to have uncorrelated parameter estimates, we would require the responses of the independent variables, written as columns of the matrix X/ t0 k>e orthogonal; this would result in a diagonal parameter covariance matrix. The precision of the parameter estimates (i.e. the error bounds) are given by the diagonal elements of the matrix" •-128-Correlations between some of the response variables. . Twill result in some columns of a X being almost linearly dependent; in addition, the presence of insensitive variables in the model will result in some of the columns of X X being all zero: in both cases the result is the same - the matrix XX will be 'nearly singular'.If we write the matrix linear Equation which represents the solution to the equation-error problem in the formA§ =■ b (3.2-48)we have from Equation (3.2-7)A i ) T X (3.2-49)M X Y (3.2-50)AA measure of the sensitivity of the relative error in © to changes in the right-hand-side vector b , and to changes in the matrix A / is given by the condition number of A (STOEOOl). There are several ways of defining the condition number. It can be shown that the ratio of the largest and smallest eigenvalues of the matrix A can be used. An inequality relating the relative change in the solution, to the condition number of the matrix A, and to a relative change in the right-hand-side vector b is derived in the reference STOEOOl; in the context of the-129-equation-err or solution considered here, it can be written as:It AO ll/ll 011 4 II All II All Ubll/llfeU =Cond (A') \\ L b W I Hb\\ (3.2-51)where 1| 11 represents a norm, and the condition number: G o o d (A) = 11 A W H All* It is also shown in STOEOOl that the condition number can be used to measure the sensitivity of the solution to changes in the matrix A .If the information matrix (i.e. ) “ representedin the foregoing discussion as /{ - is 'nearly singular', it will have at least one relatively small eigenvalue, and consequently has a large condition number.From (3.2-51) it can be seen that it is desirable from the point of view of the accuracy of the parameter estimates, to have a small condition number for theinformation matrix. This is true of the equation-error technique discussed here, and the output-error technique discussed in Chapter 4 (the conditioning of the output- error information matrix is discussed in Section 4.2.3). It is observed by the current author that the lowest upper-bound of the parameter-estimate errors - defined in terms of the norm lj ]j - is proportional both to thecondition number of the information matrix and therelative-error in the vector b (again defined in terms of the norm 11 11 ) . Furthermore, we can use the definition of b, given in (3.2-50), to associate errors in b with the-130-noise on the dependent variable Y . These observations are consistent with the expression for the parameter- estimates error given in Equation (3.2-12), which has two independent components: a scalar componentestimated as 5* given in Equation (3.2-14) - and whichcan be associated with llAbll/llbW in (3.2-51); and a matrix component (XX) which can be associated with Cond(A) in (3.2-51) .The implementation of the equation-error method used by the current author, and which is described in detail below, allows information matrices, which are of reduced rank, but which have a lower condition number, to be used in the identification. The ill-conditioning of the higher-rank problem is seen as an indication of the requirement to restructure the model; and it is observed by the current author that this is equivalent to the introduction of relational constraints between the parameter values. The derivation of the solution of the least-squares problem using singular-value decomposition, is in essence, the same as is presented in the reference NASH001; some additional references to singular-value decomposition are JACOOOl and KERR001.Sinaular-Value Decomposition.Consider the NxP matrix X given in Equation (3.2-4) containing the values for each of the independent variables arranged in columns. We can find a PxP-131-orthogonal matrix V which transforms X into another NxPmatrix B, whose columns are orthogonal. We will have:B = XV = ( b, , b*; . . .bp' ) (3.2-52)where bJ*bj = O LS? represents the squared magnitude of each of the NX1 column vectors b;_; the positive squared roots of these are referred to as the singular values of the matrix X • For non-zero singular values we may obtain unit- orthogonal vectors from:UL = bu/s,. (3.2-53)5c ^ OWe may separate the matrix of unit-orthogonal vectors, and the corresponding magnitudes of the original orthogonal vectors, into two matrices.B - U S (3.2-54)5 is a diagonal matrix with the singular values arranged in descending order down the leading diagonal. Equating (3.2-52) and (3.2-54), and using the orthogonal property of y, we have the singular-value decomposition of the matrix X--132-^ ~ (3.2-55)If we write the least-squares problem given in (3.2.4) as * X® (3.2-56)and use (3.2-54) and (3.2-55) in the above, we havey * U S V © = U S § = B & (3.2-57)V © = ^ (3.2-58)The matrix 3 contains in its columns the principal components of the matrix X/ these correspond to new sets of orthogonal independent-variable responses that relate to the new set of parameters 0 . Estimates of the original model parameters © can be obtained from estimates of 0, remembering that the matrix V is orthogonal.From (3.2-57), it can be seen that the normal Equations corresponding to the least-squares solution of y~ BjZjare:B Ty = B rB$£ 3 5'$ (3.2-59)cx .Since o is a square diagonal matrix, least-squares solutions for 0 are easily obtained as:-133-A -Z T ,S B / (3.2-60)And using (3.2-54) we have:$ -- s ' u V (3.2-61)where it is assumed that all the singular values are non zero. In the case of real data this will almost certainly be the case. However, the problem of zero singular values can still be overcome in the quest for a solution, by effectively setting the corresponding element of 0 in (3.2-59) to zero. In fact, for relatively small non-zero singular values, it may be desirable to assume that the corresponding element of 0is zero. The solution can still be written in a formwhere r is the number of non-zero singular values retained in the diagonal matrix S and arranged in descending order of magnitude down the leading diagonal. Solutions corresponding to the inclusion of the r most significant singular values in the solution are known as principal component solutions. This may be summarised by the following.equivalent to Equation (3.2-61), where S a is replaced by Siy which is defined as:-134-$ = 5 U / (3.2-63)^ s 0 I >r (3.2-64)= C --- ^r,0j0;0) (3.2-65)a.Using (3.2-58), we can relate the parameter estimates 0for the original least-squares problem given in (3.2-56), in terms of the estimates of the transformed variables0*(r) given by (3.2-63 - 3.2-65):.T A© - V g C O (3.2-66)It is noted by the current author that by using 0(r), asAdefined in (3.2-65), in order to estimate © , P-r relational constraints are introduced into the estimation; this means that the estimation is carried out in a subspace of the original parameter space.Covariance Matrix.Expressions for the covariance matrix of the parameter estimates corresponding to principal component solutions are not given in the reference NASH001, but are easily derived as follows. Using (3.2-54) and (3.2-57), we have from (3.2-12) thatCoV($-$) = cr'S* (3.2-67)-135-And using (3.2-58) we have thatC o v C © - © ) - V C o v C ^ ” ^ ) V T (3.2-68)In the case of principal component solutions indicated by (3.2-66) S is replaced by 5 defined in (3.2-62).Reduction In The Expected Value Of The Residual Sum Of Squares.For estimates 0, and measurements of the dependent variable y, we can write the residual sum of squares as:Ze'co ( y -X © )T( ' / -X0 ) = ( / - B g j ( y - B g l )(3.2-69)This can be easily rearranged to give:Xdco =■ yTy - y T(3$ - <$Tb ( .y -b$) 0.2-70)In NASH001, the expression:Ze'co - y y - / 6 0 <3.2-71)is written, without comment, as a strict equality;however, this assumes that the term ^6(^*6$) is exactlyzero, which is not true. Certainly, if the residuals: are zero-mean and uncorrelated, the expected value of the residual sum of squares is given by the right-hand-136-side of (3.2-71), that is:E H i e r c c Q = y Ty - y Ts § (3.2-72)Using (3.2-54) and (3.2-63) we can rewrite (3.2-72) as:where r is the number of orthogonal components used. Hence as r increases, it is the expected value of the residual sum of squares which must decrease rather than the calculated residual sum of squares, for an increase in the number of orthogonal components used. As can be seen from (3.2-70), the predicted values of the dependent variable (i.e.B0) would have to differ substantially from the measured values / over a sizeable portion of the data for the calculated residual sum of squares to increase, for the inclusion of an extra orthogonal component in the solution. This is possible for helicopter flight-data represented in the frequencydomain even when good fits, with large values of the. , . 2 ,squared correlation coefficient R , are obtained. Ifthere is a region of the frequency range used in the identification, where the discrete Fourier transforms have a relatively small magnitude in comparison to the important modal effects which are modelled by the fit, a small decrease in the calculated residual sum of squares could occur for the inclusion of an additional orthogonal(3.2-73)component that does not contribute to any real improvement in the fit; this would indicate that further orthogonal components should not be included in the model.Some Comments On The Condition Number And On Ill-Conditioning Due To Bad Scaling.Equation-error solutions for the use of the r most significant orthogonal components correspond to the use of an information matrix which is better conditioned than that for the usual equation-error solution where all the orthogonal components are present. This is because the condition number of the information matrix can be defined as the ratio of its largest and smallest eigenvalues. It is shown below that the use of only a subset of the most significant orthogonal components corresponds to the removal of the smallest eigenvalues from the information matrix; whether this is a valid modelling step or not is a separate issue.TFor the equation-error information matrix X X , we have, using the singular-value decomposition of the matrix X given in Equation (3.2-55) that:X X - ( U S V T) ( U S V T) - V-TLVT (3.2-74)where _TLtL = Su. • This is the familiar eigenvalue- eigenvector factorisation result often given in texts on-138-linear algebra (e.g. MURD001). The P eigenvalues of the information matrix xx are equal to the squared singular values of the matrix X.If the dependent variable Y is best approximated by an orthogonal variable corresponding to a column of 0 associated with a small singular value, then the problem is ill-conditioned due to bad scaling of the data. For applications of regression, it is often recommended that if the independent variables differ substantially in order of magnitude that a correlation transformation be applied, where the variables are described as perturbations about their estimated mean values, normalised by the respective estimated standard deviations.3.2.6 Practical Application Of The Singnlar-Value- Decomposition Approach To Frequency—Domain System Identification.The singular-value-decomposition approach was implemented in a FORTRAN-77 program - SINGVAL; details of the program are given in BLAC006. The program was written especially to perform frequency-domain identification, and incorporates all the features outlined in this document that are special to the frequency-domain implementation of the equation-error method.The program enables the user to investigate with ease, the effect on the identification of varying the-139-number of orthogonal components. That is, it provides information concerning the effect of changing the model structure on the estimated parameter values, on the indicated accuracy of individual parameters (using partial-F statistics), on the indicated significance of individual parameters (using parameter significance values), and on indicators of the accuracy and goodness of the overall fit (using the total F statistic and the squared correlation coefficient).For the practical application of the method tohelicopter flight-data using the frequency domain, themost convenient approach was found to involve building upthe model from a small number of orthogonal componentspossibly 1 - up to the case where all the orthogonalcomponents are used, corresponding in form to the usualleast-squares solution obtained using the matrix pseudoinverse. It is the value of the squared-correlation2 , ,coefficient R which is the most immediate and most useful indicator of hpw additional orthogonal components are affecting the solution. A progressive increase in the value of R z as the number of orthogonal componentsincreases up to a certain number, followed by negligibleincreases for additional components, may indicate that the orthogonal components corresponding to the smaller singular values should be excluded from the solution, since these components will consist essentially of noise. In addition, factors including the physical realism of the estimates themselves, and their indicated accuracies,-140-can be used by the analyst in choosing the 'optimum' number of orthogonal components for the solution. In the reference NASH001, the use of tolerances relating to the magnitudes of the singular values, is presented as a means of deciding if a particular orthogonal component should be excluded from the least-squares solution. For the relatively few points used in helicopter frequency- domain identification it is possible to consider the effect of all possible combinations of orthogonal components on the solution, very quickly, thus rendering the specification of tolerances on the singular values unnecessary.The ideas given above on the selection of an appropriate number of orthogonal components are demonstrated using real flight-data in Section 3.3.4.(NAG) Library Routines Used To Carry Out The Sinaular- Value Decomposition And Computational Details.The program SINGVAL uses the NAG routine F02WAF (NAG004) to carry out the singular-value decomposition of the frequency-domain independent-variable matrix X .As a result of the reduction in the amount of data used in the identification, when the frequency domain is used, the dimensions of the vector and matrix quantities involved in performing the singular-value decomposition will be relatively small when compared to those involving the use of time-domain data. For example, the number of-141-data points in the flight-data sets used for the work carried out in this document, is typically 1500; when transformed into the frequency domain, the number of data points (i.e. real and imaginary parts) used in the identification is about 20. For N data values and P parameters, the timing of the singular-value decomposition of the NxP matrix for the NAG routine F02WAF is quoted as being approximately proportional to P2 (N+6P) (NAG004); this means that for the Figures quotedbefore, the time to perform a singular-value decomposition using frequency-domain data is only 1.3% of the time required to carry out the corresponding calculation in the time domain. For equation-error identification, this is a strong point in favour of the use of frequency-domain data.3.2.7 Some Comments On The Use Of Zero Pads In Frequency-Domain Identification.For the discrete Fourier transform, the frequency resolution is determined solely by the length of the time-domain record used. This is unfortunate in the context of helicopter system identification in view of the difficulty of obtaining long records. In the reference PADF002, the use of zero-pads is suggested as a possible means of reducing the spacing between the frequency-domain points; the motivation being a more accurate estimation of the low-frequency modes which are-142-difficult to identify from a data record of limited duration. Zero-padding has been used in spectral analysis to smoothen the appearance of calculated power spectra (e.g. BUCK002). However, the use of zero-padding is, in the opinion of the current author, not an appropriate technique for system identification. A proof is offered below which shows that when zero-padding is used, the modelled noise term is in fact correlated with itself over a range of frequencies. It follows from the discussion given in Section 3.2.1, and in particular from Equation (3.2-30), that this will lead to biased estimates.The argument proposed by the current author first shows that zero-padding the time-domain data is equivalent to using a complex-valued, unit magnitude, set of weights on the original time-domain data in order to obtain an expression for the discrete Fourier transform at a frequency value lying between two integer multiples of the fundamental. Consider a data record of N sampled values padded with an additional (L-l) x N zeros, where L is a positive integer. We have for the discrete Fourier transform of a quantity V(i.At) :where CO* = 2flK/tr0At, K = 0,1,2 ... LN-1; V(lAt) = 0 for L > N-l; and 00F - 27v/rUAt, We can write:LW - \(3.2-75)-143-GJK . H kA 3 cjf + ^ / L cof (3.2-76)Hence (3.2-75) can be rewritten as:-j cck/l1) asdescribed in Section 2.2.3.Lateral Equations.For the pedal-doublet data set, the two lateral Equations considered (rolling moment and yawing moment) have the following frequency-domain representation:p CooT) - Lu OccS) + Lw w CcS) + L(j_ c^CcS) + Lv \/c.oS) * Lp pccj) (3.3-3)Lr r Cci) + L Ap f\p CcS) rCco) = Nu uccS) *• IMw wcci) + IXJ^ C^co) + lVJv + IXJppcci) +(3.3-4)r\jr fccod * i\JopAp6ci)The portion of data used was about 2 6.6 seconds in length, giving a spacing between points in the frequency domain of about 0.037 6 Hz.-149-Incorporation Of Time Delays Into The Equations.The incorporation and estimation of a time delay in the control channel used by the pilot is an important feature of the identification methodology developed in this document. It is primarily seen as a feature most appropriately and most easily accomplished at the frequency-domain output-error stage of the identification methodology. However, using a first-order approximation- \GO"Cto the complex exponential C (valid for small delays *C and small frequency values go, as explained in Section1.2.5) the following means for the estimation of delays using the equation-error method was developed by the current author. For a time-delayed quantity C(t-t), where X is the delay, and whose coefficient in an Equation (such as those given in 3.3-1 - 3.3-4) is wehave for the frequency-domain representation that:Kc C v- CicS) — K c e J ^co^o (3.3-5)Taking the first-order form and splitting it into its real and imaginary components, we have:Kc. 0 - jco“0 CCtt) = Kc C fU I. C Coil + 63X I«nZCCo}3) + j Kc ( Im CCCco)3 - COX Ret CCoXl) ( 3 .3 -6 )The right-hand-side of Equation (3.3-6) can be rewritten with terms having the same coefficient:-150-Kc ( l - J U - O CCco) = Kc R e t C C c o ^ + j K c I m C C C c O ] *K5 Re. C S CcS\\ + j Ks loa t bCo)l(3.3-7)whereRe X SCco)3 ~ <*3 Xm X CCcj)3Im C S Cco)3 = ~ CJ Rc X CCcj)3 (3.3-8)KsWe may thus construct another independent variable in theimplementation of the frequency-domain equation-error method described in detail in Section 3.2.1, we have from the definitions given in (3.3-8) that the estimate of the delay can be obtained by taking the ratio of the estimates of the two coefficients Ks and Kc:For the incorporation of a time delay in the control input, an additional independent variable, constructed as explained above, is included in Equations (3.3-1) to (3.3-4).frequency domain 5 (co) , with coefficient Ks. Using thek s / K c a x (3.3-9)3.3.2 Selection Of Frequency Range.The magnitudes of the Fourier-transform quantities-151-can be used as an indication of a suitable frequency range for use in the identification of rigid-bodyparameters. For example, the magnitudes of transformsused in the pitching-moment Equation, are shown in Figure 3-1 (for a magnitude log-frequency scale). It can be seen that above about 0.54Hz. (i.e. log10 ^ -0.25) themagnitudes of the Fourier transforms are small incomparison.The number of data points used in frequency-domain helicopter system identification is relatively small, and it is an easy matter to investigate the effect of frequency range on the estimates. For the pitching- moment Equation, the effect of using a progressively higher upper-limit of frequency in the identification (the lower limit of frequency is fixed at the smallest non-zero value of 0.0752 Hz) is shown in Figures 3-2 and 3-3, where some parameter estimates and their partial-F statistics are shown in terms of both the frequency range and the number of orthogonal components used. The squared-correlation coefficient for the pitching-moment Equation is shown similarly in Figure 3-4. Focussing on the results corresponding to 5 orthogonal components, it can be seen in Figure 3-3 that peaks in the partial-F statistics (indicating highest confidence in the corresponding parameter estimates) are observed for an uppermost frequency of about 0.75 Hz.. The parameter estimates shown in Figure 3-2 have stabilised over this region for 5 orthogonal components.-152-Hence all the indicators - the magnitudes of the transforms, and the partial-F statistics - point to the use of an upper frequency value in the range 0.5 to 0.75 Hz., with the parameter estimates themselves showing no significant changes over this range. In Figure 3-5, the squared-correlation coefficient R2' is presented as in Figure 3-4 except that the entire frequency range available (with an upper limit corresponding to the Nyquist folding frequency of 32 Hz.) is used. The value of Rz can be seen to fall away from the high levels associated with a low frequency range as higher-order effects, not accounted for in the model, are brought into play. It is interesting to note that the final low value of the squared-correlation coefficient corresponds to carrying out the identification on the original raw flight-data in the time domain.Equivalent Time-Domain Reconstruction For A Selected Frequency Range.In Section 2.2.2, an interpretation of the use of a restricted frequency range in terms of an equivalent time-domain response was given. For one quantity (forward velocity) the effect of progressively increasing the number of frequency-domain points used (i.e. terms used in the Fourier series approximation) is shown in Figure 3-6, up to the point where all available points are used (corresponding to the original measured-153-response) .3.3.3 Length Of Time-Oomain Record Used.Figure 3-7 shows, for the pitching-moment Equation, the effect of the record length used on the parameter estimates obtained; the frequency range used in the identification is approximately up to 0.5 Hz., with 5 orthogonal components used for the solutions (this choice is discussed in the next Section). The estimates are seen to reach almost constant values as the record length approaches the time span used for the results presented here; however, it is clear in some cases, such as for Mp, that a longer record than that which was available would have been desirable for the identification. The error bounds on the estimates are clearly decreasing with increasing record length: this is the consistencyproperty of the estimator given by (3.2-26) in Section 3.2.1.It is of interest to note that the parameters estimated with the greatest confidence, namely M u and / approach their final estimated values for much shorter record lengths than some of the other parameters.If any credence is to be given to the identifiedwmodel, it is important to ensure that the data records 'are of a duration which has allowed the parameter estimates to converge.-154-3.3.4 Selection Of to appropriate Member Of OrthogonalPitching Moment EquationConsider first the pitching-moment Equation, where the frequency-domain form of the Equation is as given by (3.3-1). The model incorporating an approximation to a time-delay in the control input will be considered secondly, and the form of Equation (3.3-1) corresponds to the results shown in Figures 3-2 - 3-4; tabulatedparameter estimates for different numbers of- orthogonal components are presented in Table 3-1 a).For this example, the squared-correlation coefficient (R2 ) is seen to increase. This is as would be expected in the statistical sense and is consistent with the discussion in Section 3.2.5 on the expected reduction of the residual sum of squares as more orthogonal components (associated with singular values of a progressively smaller magnitude) are used in the solution. The size of the increase in R2 in going from one structure to another, can be seen from Table 3-2 a), to be associated with the 'significance' of the newly introduced orthogonal component. A convenient measure of 'parameter significance' (see Section 3.2.4) was seen by the current author to be given by the absolute product of the estimated coefficient of the orthogonal component which remains unchanged for the addition of further-155-orthogonal components - and the associated singular value (which, as was explained in Section 3.2.5, represents the squared magnitude of the orthogonal component). In Table 3-2 b) partial-F statistics of the orthogonal variables are presented for the principal component solutions. It is interesting to note that for the solution corresponding to six orthogonal components (i.e. the conventional least squares solution), the least significant orthogonal component is seen to have a partial-F statistic of less than unity, implying that the estimate is an approximation to zero: the sixth component may well be associated mostly with noise. The sixth orthogonal component can also be seen, from Table 3-2 a), to have a relatively small measure of 'parameter significance', defined by the absolute product of the orthogonal-variable estimate and the corresponding singular value (given in the Table as ABS. PRODUCT a) x b) ) .The relatively small difference between the high values of F? given by 5 and 6 orthogonal components is one indicator that the solution corresponding to 5 orthogonal components may be selected. In deciding on the number of orthogonal components to use, the physical realism of the model estimates can also be taken into account, along with indicators of parameter accuracy (partial-F statistics and parameter significance). For example, the solution corresponding to 5 orthogonal components is seen to give the highest values of-156-partial-F statistics for the important parameters, M u, M^, (see Figure 3-3) ; five orthogonal componentsalso give the most physically realistic estimates - over all the other solutions shown - of pitch-damping Mpand control sensitivity .The parameters and Mp do not play a significant part in the solutions - and consequently, are estimated poorly - for the inclusion of up to 4 orthogonal components. For 5_ orthogonal components, relatively large entries in the 3rd and 5th columns of the 5th row of the matrix V"*" (relate to Equation 3.2-58 in Section3.2.5), which is shown in Table 3 c), corresponds to the effective use of pitch rate and roll rate in the estimation. In Figure 3-8 'parameter significance' values are given for the model variables, for all the principal component solutions. They show the importance of M0, M w and M/j in the first few principal component solutions. However, the solution corresponding to 5 orthogonal components shows for the stability derivatives, 'parameter significance' values which are all of similar magnitude; these are included in Table 3 - l b ) .The incorporation of a time delay in the control has been advocated throughout this document as a more appropriate means of producing parameter estimates that can be compared with theoretical values obtained from a quasi-static model. An approximate method for the incorporation, and estimation, of time delays using the frequency-domain equation-error method was offered iniK ©-157-Section 3.3.1. Parameter estimates are given in Table 3-3 a) and these are discussed in detail in the nextSection. Using the same approach as for the case withoutthe delay, it can be seen that this time, all theevidence points to the selection of the model with 7orthogonal components. There are no principal-componentnzsolutions, with high values of K , which appear to have (associated with the smaller singular values) orthogonal components that are relatively insignificant (see Table 3-4 a), or that have very small partial-F statistics (see Table 3-4 b). 'Parameter significance' values and partial-F statistics of the model parameters, for the selected number of orthogonal variables are given by Table 3-3 b).Tables 3-5 a) and 3-5 b) contain parameter estimates for the pitching-moment Equation without and with a time delay in the control respectively, for a wider frequency range of 0.075 to 0.753 Hz.; it can be seen that there is no significant difference from the estimates obtained for 0.075 to 0.527 Hz., for the selected number of orthogonal components.Normal Force Equation.The parameters in the normal force Equation, where the frequency-domain form is given by Equation (3.3-2), were estimated using the same frequency range as the pitching-moment Equation. The parameter estimates are-158-presented in Table 3-6 a); it can be seen that the squared-correlation coefficient levels out in value once5 orthogonal components are used. The use of 7 orthogonal components would correspond to the conventional least- squares solution. However, from Table 3-6 c), it can beseen that there are 2 relatively insignificant orthogonal components associated with the 2 smallest singularvalues. It should be stressed that it is not possible to make the final choice of the number of orthogonalcomponents to be used in the solution completelyobjective, but on the basis of the information given in Table 3-6 a) and 3-6 c), 5 orthogonal components is ajustifiable choice. For the selected number of orthogonal components, the 'parameter significance' values and partial-F statistics are presented in Table 3-6 b). It can be seen here that Zw i-s the mostsignificant parameter by far in the Equation; the control sensitivity is the next most significant parameter,and has a fairly large partial-F statistic. The'parameter significance' values are shown in Figure 3-9.It is useful to remember that the equation-error approach is seen in the identification methodology developed in this document both as a means of obtaining initial guesses for the output-error method and inisolating insignificant parameters through the use of 'parameter significance' values. If the normal-force Equation, with parameters identified using the equation- error approach, was itself intended for use, then further-159-model simplifications would be required in order to give the model satisfactory predictive qualities.Insignificant parameters, estimated with great error, would be removed to give the simplest model adequate for the intended application.For the normal-force Equation with a time delay estimated in the control, the parameter estimates are presented in Table 3-7 a). The evidence points to the selection of 6 orthogonal components: the levelling-off of the squared-correlation coefficient; the presence of 2 relatively insignificant orthogonal components associated with the 2 smallest singular values (see Table 3-7 c)); and the emergence of an estimate of the time delay in the control which is very close to the value also given from the corresponding pitching-moment results (c.f Table 3-3a) ) .For the selected number of orthogonal components, the 'parameter significance' values and partial-F statistics are presented in Table 3-7 b): the mostimportant and most confidently estimated stability and control derivatives are once again 2^ / and .Yawinq-Moment And Rollina-Moment Equations.Tables 3-8 a), b), c) relate to the yawing-moment Equation (frequency-domain form given by Equation (3.3— 4)). The corresponding results for the model with the delay in the control are given in Tables 3-9 a) and b).-160-The selected number of orthogonal components are 6 and 7 respectively: for the selected models - particularly for the model with the delay - the important yaw-damping parameter estimate Nr is very close to theory. Results for the rolling-moment Equation, when a time delay is estimated in the control, are given in Tables 3-10 a) andb). The selected model corresponding to the use of 7 orthogonal components gives the estimate of the important (as evidenced by Table 3-10 b)) derivative, Lv r that is closest to theory. It should be pointed out that the rolling-moment Equation has exactly the same independent variables as the yawing-moment Equation: consequently the associated singular values are identical in both cases (compare Tables 3-10 c) and 3-9 c)); for the two Equations, the smallest singular value was judged to be an approximation to zero.3.3.5 Comparison Of Estimates With Theoretical Values.Since actual experimental flight data have been analysed in the application of the equation-error method demonstrated here, there is no set of 'true' parameter values with which the estimates can be compared. The helicopter flight mechanics package HELISTAB provides theoretical parameter values which may be considered alongside the estimates obtained from flight data.For the pitching-moment results presented in Table 3-1 a), it can be seen that the selected model,-161-corresponding to the use of 5 orthogonal components, results in estimates of the parameters M^, Mp and I%IS which are closer to theory than the conventional least- squares solution (given by 6 orthogonal components). Mu compares very well with theory, whilst Mw anc* Mv are estimated with relatively small error bounds and are of the same order of magnitude as the theoretical values.At this stage in the discussion of pitching-moment results it is of interest to consider the results given in Table 3-11 a) obtained from flight data for a Puma helicopter at a trim speed of 80 knots with a longitudinal cyclic '3211' input applied. For the selection of 5 orthogonal components, the results for the 80 knots case display some similar characteristics to Y those obtained for the 100 knots case. The emergency^ of the most physically realistic estimate of issignificant and an estimate of Mp which compares very well with theory, and which has a small error bound is another important feature of the results which contrasts with the estimate for 6 orthogonal components. In addition, the estimate of compares very well withtheory, whilst My, is estimated with a relatively small error bound and the magnitude-trend for Mv (i.e. the change in going from 100 to 80 knots) is consistent with theory.The most significant discrepancy between the flight- data estimates and theoretical quasi-static values in the case of the results presented in Tables 3-1 a) and 3-11-162-a) was the underestimation of the parameter . From Table 3-3 a) it can be seen that when an approximation to a time delay is incorporated and estimated in the longitudinal cyclic control, the pitch-rate damping parameter is significantly bigger in magnitude, and compares very well with theory. The control sensitivity also compares much better with theory than the case without the delay. It is reassuring to note that the approximate time delay estimate is almost identical to the value obtained for the application of the output- error method to the same data set, described in Section 2.3.3; in fact, a great many similarities can be found when the two sets of results are compared. It is worthwhile pointing out that in the case of a time delay in the control, the iterative procedure for the output- error results of Section 2.3.3 was initiated using equation-error results corresponding to an estimation without a time delay in the control; an initial guess of zero was used for the time delay.For the normal-force-Equation results presented in Tables 3-6 a) and 3-7 a) - without and with a time delay in the control respectively - it is only worth discussing estimates of the dominant parameters: 7w and Z^1S . For the selected models, the time delay is shown not to affect the estimates of Z w t this observation was also made for the output-error identification. The selected Zw estimates have relatively small error bounds, and are seen to be of a slightly greater magnitude than the-163-theoretical value. The selected Z i\IS estimates compare favourably, in both cases, with the theoretical value, and have relatively small error bounds. The time delay is estimated with a value that is almost identical to that obtained for the pitching-moment Equation, and ipso facto it is in close agreement with the value obtained for the application of the output-error method in Section 2.3.3. There is a large difference between the estimated delay for the selected model (6 orthogonal components) and that which is offered for the conventional least- squares solution (8 orthogonal components); this observation constitutes a strong argument in favour of the selected model in preference to the conventional least-squares solution.For the yawing-moment results presented in Table 3- 8a), it can be seen that the selected model, corresponding to the use of 6 orthogonal components, results in an estimate of the yaw-rate damping parameter f\jr that is both substantially different from the value offered for the conventional least-squares solution (7 orthogonal components), and is closer to the theoretical value of the parameter. The estimate of 1\JV is slightly greater than theory, and this parameter is estimated with a relatively small.error bound; conversely, the estimate of the control sensitivity is slightly smaller inmagnitude than theory, and is also estimated with a relatively small error bound. The estimate of R)p has a relatively large accompanying error bound, and differs-164-somewhat from theory. On the basis of the 'parameter significance' values presented in Table 3-8 b) the lateral/longitudinal coupling parameters: HJU/ KJw and i\)^ play a relatively insignificant role in the Equations; this contrasts with the pitching-moment Equation where longitudinal/lateral coupling intensity was shown to be significant for the flight condition considered.With a time delay in the control there is a noticeable reduction in the indicated error bounds of the important yawing-moment derivatives: !\)v, f\)r an<^ !\)r\pfor the selected model over' the case without a delay, as shown in Table 3-9 a). In addition, the estimate of the yaw-rate damping LJr is in very close agreement with theory as a result of the time delay being used in the control.Estimates of the rolling-moment Equation are presented in Table 3-10 a), and it can be seen from Table 3-10 b), that like the yawing-moment results, lateral/longitudinal coupling effects are relatively insignificant. For the selected model (7 orthogonal components) it can be seen that the estimate of Lv matches very well with theory, and has a relatively small error bound. The estimate of Lp is smaller in magnitude than its theoretical counterpart, whilst conversely the estimate of Lr is larger than theory suggests. The tendency for the derivatives Lp and Lr to have this reciprocal relationship is considered again later in Section 4.3.2, where the estimation of yawing-moment and-165-rolling-moment derivatives is examined in detail using the frequency-domain output-error method with rank- deficient solutions; it is shown that more satisfactory estimates of lateral parameters, particularly Lo and L r , can be obtained, than by considering the yawing-moment and rolling-moment Equations in isolation, as is the case with the equation-error method.Measured And Predicted Frecruency-Domain Values.In addition to consideration of the individual parameter estimates themselves, further evidence of the overall validity of an identified model can be obtained by comparing measured frequency-domain values with corresponding values predicted by the model. A convenient way of displaying the fits is to use the magnitudes of the complex-valued quantities, remembering that the identification is performed by fitting both the real and imaginary parts of both sides of an Equation. In Figures 3-10 a), b) the frequency-domain fits for the pitching-moment Equation, without and with a time delay estimated in the control, respectively, are presented for solutions obtained for different numbers of orthogonal components, including the selected model. similarly, Figures 3-11 a), b); 3-12 a), b); and 3-13 a),b) areprovided for the normal-force, yawing-moment, and rolling-moment fits respectively.-166-3.3.6 Am Equation—Error Multirun Example.An investigation into the use of four combined manoeuvres for frequency-domain equation-erroridentification was carried out. The manoeuvres involve the use of all four controls, and consisted of: acollective doublet; a longitudinal cyclic '3211'; a pedal doublet; and a lateral cyclic step input. The responses were stacked to produce a record length of just over 60 seconds as shown in Figure 3-14. Frequency-domain fits for the pitching-moment, normal-force, yawing-moment and rolling-moment Equations are presented in Figure 3-15. The frequency range used extended up to about 0.5 Hz. with 35 different frequencies used in the identification. It was observed that the estimates of some longitudinal parameters were closer to theoretical values than for the single-input longitudinal-cyclic case described earlier, for the selected model, namely M.w and - estimated as -0.0040 [0.0009] and -0.1965 [0.07] respectively (wherethe square brackets indicate the estimated error bounds). The corresponding theoretical values are given in Table 3-1 a). On the other hand, some important longitudinal parameters, namely and - estimated as 0.0015[0.0003] and -0.0593 [0.15] respectively, had estimates which were worse than the single-inputlongitudinal-cyclic case in terms of physical realism.tThese findings appear to support the comments made in the reference DUVA001 - mentioned in Section 3.2.3 - about-167-the shortcomings of using several stacked data sets in the identification. It should be noted however that the data used for the investigation do not constitute an ideal set. The fact that a step input, rather than a (small amplitude) multistep input which involve a return to the trim setting, was used, is one notable deficiency. This choice of lateral input was dictated by the available test records for the flight condition. The results for the multirun identification are presented in more detail in the reference BLAC005.The presence of discontinuities presents no problems for the use of time-domain data in equation-error identification, however, for the frequency-domain approach exemplified in this Section, the Gibb's oscillation effect - described in Section 2.2.2 - will be a source of error in the identification.3.3.7 Consideration Of The Predictive Qualities Of A Model.If there are sufficient flight-data sets available to the analyst, the identified model should be tested on a completely independent data set. Figures 3-16 a), b) provide examples of this type of assessment, where frequency-domain values are shown for the response to a longitudinal cyclic '3211', together with predictions based on the identified model using the longitudinal cyclic doublet (discussed in Section 3.3.5 and-168-corresponding to Table 3-1 a)). In addition, the response to the longitudinal-cyclic doublet input, and predictions based on the multirun model, considered in the last Section, are also shown. The overall agreement between the measured and predicted responses is reasonable in both cases.- 169 -CHAPTER 44.0 THE FREQUENCY-DOMAIN QUTPUT-ERRQR STAGE IN THE IPEUTIFICATIQU METHODOLOGY.4.1 Introduction.4.1.1 Background To The Use Of The Output—Error Technique.The frequency-domain output-error method is presented here as a special case of the more general maximum-likelihood method, where there is the assumption of no process (i.e. modelling) noise. The original features of the method presented here are given inSections 4.2.2 to 4.2.4, and are associated with modelrepresentation and computer implementation; the cost function used as the basis for the method is the same as that derived in KLEI004.A significant distinction between the output-error and equation-error techniques is that the former requires start-up values, or initial guesses, to perform the estimation. For the implementation used in this work the minimisation of the cost-function is performediteratively using the Gauss-Newton method (with anadditional line-search modification). Good initial guesses are required for rapid convergence because of the representation of the cost-function in the vicinity of- 170 -the minimum by a quadratic. Whereas the equation-error method can only be used to identify one equation at a time, the output-error method identifies state-space type models. The frequency-domain output-error method is presented here, as the second, and most important stage in the helicopter system identification methodology that is developed.4.1.2 Previous Applications Of The Technique.Maximum-likelihood/output-error methods havefeatured often in previous helicopter and aircraft system identification work. Most reported applications have used the time domain (e.g. FOSTOOl, FOST002, ROSSOOl, KLEI007, STEP001, FEIK002). Applications to the problem of helicopter system identification include: HALL001,MOLUOOl, MOLU002, MOLU003. Real progress with these techniques for helicopter system identification, however, have been disappointing, with few reported successes in the literature. The large computation times required, coupled with the dificulties of obtaining an adequate model structure in the time-domain, (without introducing a formidable amount of complexity into the problem through, for example the inclusion of rotor states), have severely limited the use of this method for helicopter system identification.The use of the frequency domain, for the output- error method, is advantageous for several reasons,- 171 -including the ease with which innovative modelling features can be incorporated (e.g. the use of delays in the controls, and the use of a specific frequency range for the estimation). In addition, use of the frequency domain can lead to a reduction in the amount of data (and hence computation time) required. In the reference FUKH 001, the authors report some success in the application to helicopter flight data of a frequency-domain identification technique based on the same cost-function as that developed in section 4.2.1. In MARC002, the authors present some results, again using the frequency- domain technique, to data generated from an in-flight- simulator. Defined time delays existed between the commanded control inputs (pitch stick, roll stick) and the corresponding control surfaces (elevator, rudder), and these were successfully identified leading to improved fits. In Chapter 2, the idea of using time- delays to account for higher-order rotor effects, in addition to any purely transport effects, was advanced as a new and original tactic in helicopter system identification. Some results using real flight data were used to demonstrate the resulting improvement. Further examples of the use of a time delay with this interpretation, are presented later in this chapter (and in the references BLAC005, BLAC007).- 172 -A Further Comment.Although not used to estimate the stability and control derivatives, time-domain estimation does have a place in the estimation methodology presented in this document. It is used to estimate zero-offset and initial state conditions (see chapter 5), following a frequency domain identification, of the main model parameters.4.2 Development Of The Freouencv-Domain Output-Error Algorithm.4.2.1 Formulation Of The Frequency—Domain CostFunction.Set Of Transformed Measurements.Consider the following equations representing the state-space model and the measurement system (delays in the controls and measurements have been omitted from the notation here for clarity), and using the same notation as in section 2.2.4:= A X C O + Butt) + WCt)/v /V A.2 CO = H X C t ) + K * b VCt")^ A *WCt) is a process (or modelling) noise term.(4.2-1)(4.2-2)- 173 -Transforming the above into the frequency domain, we have:VV A AX C gV ) = A X Coj*') + B U C o O * ( 4 . 2 - 3 )^ /% r*- 'A a AZ . ( O k') = H X C < 0 , 0 + V C c^ k') ( 4 . 2 - 4 )where coK^o , and Cok = KX7t/l\iAt for N time-domain .-samples of interval A t . Let be the set of q frequency-domain measurement values. If W t o and Vet) are independent Gaussian random variables, then so are the transformedA Aquantities W a n d VCcjO (SCHOOOl) . We can define theA- Afollowing set:Z , A A A£ = I ZCcoO , • • • Z C o ^ j ( 4 . 2 - 5 )where o < Maximum-Likelihood Principle.The maximum-likelihood principle can be defined as follows: Given an experiment whose outcome depends on a set of unknown parameters Y the maximum-likelihood estimates of the unknown parameters are those values ofY, for which the observed values in the experiment, aremost likely to occur. This involves maximization of the conditional probability-density function of the observations given V ; and is called the likelihood- 174 -function. For the set of measured values this can bestated as:Km* -- max P (4.2-6)PtZj/sri is the probability density of ' h conditioned on Y , and ^ are the maximum-likelihood estimates of Y • The probability of obtaining the set *Z<^ can be obtained by concatenating conditional probabilities (using Baye's theorem):pczt/n = ptzc^/zr ,,n • pEz*-,/n(4.2-7)Repeated application results in:PEzt/n = a Ptzcoo/z«,r: <4.2-8)To find the probability distribution of 2 (coj given 21u_, and Y , the mean value and variance are determined first. Let us define:r\E L Z t o o / Z . . , X I - Zcoj./l-0 (4.2-9)E tC Z c c o o - £ c « . / l - o ) ( . Z c < j j > - Z c w . / L - o f ] * 4 - 2 - 1 0 )= 5 Cart- 175 -In order to obtain Z (£){./ l-\ ) , consider first the discrete steady-state Kalman Filter representation of the system given by (4.2-1) and (4.2-2) (e.g. KLEI004; K and b addedby current author)X (n+t = d) X(nAt') + D U C n a O * Kv^cnAt) (4.2-11)/ .(4 . 2-14)A. *where the covariance matrix S i s independent of frequency for white time-domain Gaussian noise sources V(t) andA*W ( t ) .Output Errors Instead Of Kalman-Filter Innovations.If we assume that there is no process noise term W(t) (i.e. no modelling errors) in (4.2-1), then the steady-state Kalman-gain matrix K is zero (APPENDIX 4). The innovations \T{coK) become output-errors V(<^ *) • HenceA n .we have for (4.2-13):ZCco*") = H XGp*D * $ ( < & * ) ) <=->* £ O (4.2-15)AFor 2(C0i./>--\ ) defined in (4.2-9), we can see from(4.2-15) that it is dependent solely on the model output at that particular frequency, and is consequently onlyconditioned on the model parameters- 177 -Log Likelihood Function.Returning to (4.2-8), we may now write the conditional probabilities:P I T . o / n = <4.2-16)Using the Gaussian probability density function, the following argument for the derivation of the likelihood function is similar to'that presented in the references KLEI004 and KLEI005, except that output- errors have been used instead of innovations, following the assumption of no process noise given earlier. The probability densityAfunction for the Gaussian variable Z(coJ conditioned on Y is given by:L ( 2 7 0 l (b\3 * e x p E r H C V o ajj*S' VCao)3 (4.2-17)Including terms of the form given in (4.2-17) into the product given in (4.2- 16), and taking the naturallogarithm of both sides, we have:loge P L Z , / n = i l o 3 e P C f C o O / t l (4.2-18)- 178 -And using (4.2-17) where:PC'ZCcdi)I^ CJ - \o^\S\ “ Vz V Coif 5 VCoO- V z ( 2 * ) (4.2-19)we have:C^-looje P tZ l / r3 = '/a I ( 02m£)S*\/cmo + \o3 * \ s 0-V C orvsV a 'aV(4.2-20)Ignoring the constant term, which is of no consequence for the optimization, and noting that the log function is monotonic increasing, a suitable cost function is given by:T ( ^ , S') - /^z 2 C VCcu) 5 Vccj) + \oc^ls\)(4.2-21)This is the cost function for the frequency-domain output-error method used in this document. The change of sign in Equation (4.2-21), as compared with Equation (4.2-20), should be noted. It can be seen that the cost function is a function of the model parameters ^ (i.e.elements of the matrices A., B and H given in (4.2-1) and (4.2-2), in addition to any time delays included in the model), and elements of the frequency-domain error-- 179 -covariance matrix S (defined in (4.2-10), but independent of frequency for white measurement noise). The cost- function summation. is carried out over the desired frequency range for the estimation.Some Further Comments.It is appreciated by the current author that the above justification for the cost function, based on a probabilistic argument, may be considered unsatisfactory for a number of reasons. Firstly, the assumption that there are no modelling errors might seem to be patently unjustified. In the context of helicopters, obvious modelling errors which spring to mind, are ones arising from non-linearities (i.e. large excursions from the trim), and the coupling of rigid-body dynamics with rotor dynamics. In the case of non-linearities, the use of more appropriate control inputs (discussed in section 1.4.2) can reduce this problem. This is a situation where prevention is better than a 'cure', because non- linearities are deterministic effects which are obviously not appropriately modelled as a Gaussian random effect. For the effect of rotor dynamics on rigid-body states, the use of time delays to re-model the effective control input has already been discussed, and it has been emphasized that this is only an approximate method. There is also the problem of having noise on the controls, and this is commented on in section 4.2.5.- 180 -The cost-function given in (4.2-21) is of the familiar least-squares form when the error-covariance matrix S is known, and this is the sort of ad-hoc cost function which one might define as the basis for a model fitting technique. However, when Sis unknown, and has to be estimated in addition to the model parameters Y , then the logarithmic component must be included. Finally, it should be mentioned that the cost-function given in (4.2-21) can be modified, when a statistical model is available for the parameters as well as the measurements, by the inclusion of an additional term. This corresponds to a Bayesian estimation approach (e.g. GELB001), and for the estimates (4.2-6) is replaced by:= r o a * P t L Z • P m (4.2-22)As is pointed out in MAIN002, the probability distribution of parameters in a model, such as doesnot represent any randomness in the value of Y , but reflects knowledge or information about the value of Y • For the work covered in this report, however, there is no assumption of an a priori distribution of the parameter values.4.2.2 Implementation Of Gauis s-Newt on MinimisationThe cost-function given in (4.2-21) is to be- 181 -minimised with respect to the unknown parameters in the model. The values of the parameters at the minimum of the cost-function are taken as the estimates. The minimization scheme used is a modified Newton-Raphson technique; this is also known as a Gauss-Newton method (e.g. MAIN002, WALS001). An additional scalar line- search algorithm (FOST002 and Appendix 3) is also used.A Convenient Symbolic Representation Of The Model.A convenient symbolic representation of the frequency-domain model, which facilitates the formulation of required terms for the minimisation, is that developed in Section 2.2.4. This involves the stacking of the real components of a complex-valued vector on top of the imaginary components; and the development of the corresponding matrix terms through equating separately the real and imaginary parts of both sides of an expression. Let us define the following shorthand notation:V. &Re t ttoii AccJ) * - A - « XI m L X Cco)] col - A .U CoS) = R e C ( W/*• i i B 0Imltocoti- /V « 0 B.G CcS) * _ ,- W t(4 .2-23)A X C o s 5^ K =1X o ---1AX / o H.- 182 -For the frequency-domain model (i.e. system and measurement equations) we can thus write:JK-CcS) XCcS) - U CcS) -v & C o )' A /S. (4.2-24)ZeeS) - y t )(CcS) -+■ Vgo)A. '(4.2-25)If we include time delays in the model, then the VJ(co) and 2(co) vectors are pre-multiplied by matrices whose elements are frequency-dependent sinusoidal and cosinusoidal functions of the time delays (See Equation 2.3-8 in Section 2.3-1). Representing these matrices by \ ( Q ) and 3~l(<^>) f then in place of the above equations we have:A CcS) X CcS) = Ha.C6oO.ccs) + G-CcS) (4.2-26)H zCco)ZCco') = H. X’Cco) + ){CcS) (4.2-27)For the case of r controls, the definition of 1^®°) must be:Cos COtiC o i W t l QCos cj X .f— Sir\ojr»— S»rv oj TZ\o*bm wr^ qS\n cJTTrCoScJTTi Cos w C z _\ O \O \C o S c jf r(4.2-28)- 183 -The matrix is defined in an identical fashion toA Convenient Symbolic Representation Of The Cost FunctionUsing (4.2-15) and (4.2-21), and the definitions:(4.2-29)Cco) 4 R e C E t l c ® - ) ! !XnfEtiCcoXfl (4.2-30)where, for white measurement noise - £[V(w)] = O r we can re-express the frequency-domain cost function as:t ( r ,s ') * ' A L ( c t e a ) - c A03 ~ ~ ~(4.2-31)For white Gaussian measurement noise, the error-covariance matrix is independent of frequency, and so for [\L. frequency values used in the estimation, we have:T A , S ) = '/z 2. ((ZccS) - Kccsfi$~X % - %&>)')I f CO ~ ~(4.2-32)- 184 -Minimisation Of The Cost Function.Since the objective is to minimise the cost function with respect to the unknown parameters, which are (the set of model parameters) and 5 (the set of elements of the error covariance matrix 5), then we must have:( ... ^ 7C1f,S) - Q (4.2-34)Assuming that the matrix 5 (and hence ^ ) is completely unknown, except for the restriction that it be po'sitive semi-definite, then for (4.2-33) and (4.2-34) respectively (using matrix differentiation results for (e.g. MAINO07, GELB001) , we havea V a ^ = o - Z(Clcco-)-^cco-))T $ ~ ' V s L ' { i l f c < s ) - % c c S ) S )S ^ A. ' ' * A*(4.2-35)= o = - '/a 6= O (4.2-33)+ /2 l\)cj 6 (4.2-36)- 185 -Using (4.2-36), givesP>esY r '/rOa) Il C Z C cS) ~ ^ Cco))(ZCcd) -}iCcS)) (4.2-37)C O ^ ^The result given in (4.2-37) is consistent with the definition of 5 (with the dependence on frequency omitted) given in (4.2-10) because of the equi- distribution of power between the real and imaginaryparts of the transform of a Gaussian random variable,discussed in section 2.2.5.Separate Algorithms For The Estimation Of X And S .A distinction has already been made between the model parameters: X , and elements of the error-covariance matrix: S • In terms of the minimisation of the cost- function, the conditions for the minimum, in terms ofeach of these sets of parameters, were given in (4.2-35) and (4.2-36). It would be possible to apply the sameoptimization technique to both sets of parameters simultaneously, and for a gradient minimisation scheme the right-hand sides of (4.2-35) and (4.2-36) would beused to supply the required first-order partial derivatives. An alternative approach, and one which is used here because of its ease of application, uses separate algorithms for the estimation of X and S . This- 186 -approach is called axial iteration, and is one which has been widely, and successfully, used in optimization problems analogous to the present one (e.g. FEIK002, FOST002, MAIN002). Convergence in both algorithms is required for the estimation procedure to terminate. The Gauss-Newton method is used in the estimation of ^ , and is written as:ie.Y = - M % (4.2-38)where the super-script refers to the iteration number. M K is a linearised approximation to the Hessian matrix (this will be explained in more detail shortly), at the Kih iteration, such that:*(4.2-39)where the second-order partial derivatives are evaluated for the parameter estimates at the Kth iteration. 9 is the gradient vector of the cost function at the Kth iteration (similarly evaluated for the current parameter estimates), such that:- 187 -g K _ C aT/asr,, 3T/ a ; . . . ^ / a r - G (4.2-40)Elements of ^ are estimated at each iteration step using(4.2-37), where the model output terms' 7 C (w) are calculated using the current estimates of {^. By using two separate algorithms, the dimension of the unknown parameter vector that would have to be used in (4.2-38), is reduced.Summary Of The Axial Iteration Algorithm.The procedure may be summarised as follows:1) Estimate 5 using (4.2-37) and the current estimateof Y2) Perform one iteration of the Gauss-Newton algorithm (using 4.2-38) to revise the estimate of Y3) Repeat steps 1 and 2 until some convergencecriterion is satisfiedOnly one iteration of the Gauss-Newton algorithm is carried out at stage 2, because it is pointless improving on the estimates of Y for a given 5 any further, if 5 itself is to be changed at the next pass through the overall procedure.The situation could be further simplified if S was- 188 -known and fixed throughout the estimation. This wouldcorrespond to a weighted least-squares type estimation(FOSTOOl), with the logarithmic component of the costfunction becoming a constant, irrelevant to theoptimization. The elements of S are fixed by the analyst- \in this type of estimation (in fact 5 is the weighting matrix used), and can be chosen to reflect the importance attached to the fitting of each measured quantity. For the applications of the method presented in this document, this will not be the mode of operation used: 5 and X will both be estimated using the axial iterationprocedure. However, the ability to fix elements of thematrix S at user-defined - values for the first fewiterations of the procedure, is a very useful feature.This is standard procedure for the application of minimisation schemes similar to the one described here. Experience has shown that if the initial parameter guesses (necessary to start the iterative procedure) are far from their true values, then the corresponding estimate of S (and hence J~>) given by (4.2-37) will be very poor because of the significant modelling errors. This, consequently, could result in a divergence of the algorithm.Convergence Of The Non-Loaarithmic Component Of The Cost Function.It can be seen by substituting for the converged- 189 -estimate of {d into the non-logarithmic component of the cost function in (4.2-32), that for convergence of the estimation algorithm, the non-logarithmic component of the cost function should tend towards the integer value - rUx x IUj, where ID* is the number of states in the state- space equation, and f\Jw is the number of frequency points used in the estimation.Zz - W x ' W w (4.2-41)GOThis was found by the current author to be a useful check that the algorithm was converging correctly.Expressions For The Sensitivities (Or PartialDerivatives).Expressions for the elements of the matrix rl and the vector 9 given in (4.2-39) and (4.2-40) respectively are obtained as follows. Using (4.2-35) for a point away from the minimum (thus dropping the equality with zero), we have for an element of 9 , corresponding to theunknown parameter K. :^ ( ( fecS) - ^ CcSj) $ IL^Cco) - K CcS)](4 .2-42)By partially differentiating (4.2-42) with respect to asecond parameter - ^ , we have:- 190 -■a. ^^ " X (C ZCw) -Xto)) $> 3^*Cu3lk'm C 7 c■d ' C . ~ ^-V H - a X C o ' ) (4.2-48)The steps required to solve for B_ [?(( (f)j^ > ^ 3 ) • • • j O i O j O ) (4.2-58)For the symmetrical matrix M we can formally write the following expression for the inverse of a rank- deficient M:M = X Xc (VcV ur) (4.2-59)U-ft 1where r is the number of eigenvalues removed from the information matrix M. Whilst the above expression is not used in calculating a rank-deficient 1^., it highlights the fact that the smallest eigenvalues are associated with directions in the parameter space that have the most uncertainty associated with them.An example of the application of rank-deficiency is given in section 4.3.2. Using (4.2-58) in (4.2-55) it can be seen that the removal of insignificant eigenvalues corresponds to the introduction of relational constraints on the update increments A^. As with the case where there are known relationships between sets of parameters, there is a reduction in the dimension of the parameter space effectively used in the estimation. The advantage of the use of rank-deficiency is that it may not be possible for the analyst to easily define the relationships (or correlations) that exist between a large set of- 203 -parameters. However, with the . rank-deficient approach, the relationships (associated with near-zero eigenvalues) are effectively constructed automatically from the measured data used in the identification; more will besaid of this in the context of the practical examplesgiven in section 4.3.2.4.2.4 Controlling Convergence Of The Algorithm.The Gauss-Newton iterative technique used in the minimisation of the cost function J requires the specification of criteria to indicate when the estimates are acceptable. There are several possible ways of doing this. These include the following, suggested in thereference FOST002: (1) inspect the size of 3T, and enda tthe search when the elements of this vector become small; (2) end the search when the change in J is small; and (3) end the search for small changes in each of the elements of T. As pointed out in FOST002, the problem with thefirst approach is that it is difficult to define what constitutes a small gradient vector in a way that is readily appreciable. It is also pointed out that the problem with the second approach is that, for parameters weakly defined by the observation data, the algorithm can terminate at values which are highly dependent on initial guesses. The third approach is given as the most desirable, where what constitutes a small change is specified individually for each parameter. This is- 204 -better than using percentage changes for all the parameters, because it avoids difficulties that could be experienced for parameters whose true value is near zero; in addition, it indicates acceptable levels of accuracy for each of the parameters, which as the author points out in FOST002, is no bad thing.A problem that may be experienced is that the algorithm may converge not to the desired global minimum, but to a local minimum. It should be apparent when this occurs from the un-realistic (or non-physical) parameter estimates obtained on completion; in addition, the cost value will be higher than might be expected, based on experience from similar runs. As was stressed in section4.2.2, the performance of the identification algorithm is very dependent on the initial guesses - hence the requirement for good initial guesses to ensure rapid convergence.4.2.5 Some Further Comments On Neglecting Process- Noise Terms.The Problem Of Noisily Observed Inputs.In the development of the cost function for the frequency-domain output-error method presented in section 4.2.1, it was assumed that there was no process noise. This assumption was made to avoid the considerable complexity that would be incurred as a result of- 205 -formulating the model output in terms of a Kalman filter. This decision was then rationalized by saying that the effect of some unmodelled inputs could be lessened: forexample, avoiding the use of large amplitude and/or large clock-period test inputs to lessen problems due to non- linearities. In addition, in the case of unmodelled rotor effects, the incorporation and estimation of a time delay in some controls, could improve the modelling of the effective control input. However, the problem of having noise on the control inputs is not one which is easily dismissed. This is particularly so for the approach presented in this document: that of using an extended control vector to include in addition to the applied pilot control (possibly with a time delay), additional (noisy) terms to drive the state-space model. This approach was convenient because it enabled the helicopter identification problem to be cast in a manageable size.The problem of having noisily observed inputs is equivalent to having a process noise term in the model. As in the case of measurement noise on the independent variables in the equation-error approach, the presence of process noise in the output-error approach will lead, in general, to biased estimates. One possible approach, which if not a solution, is at least a palliative, is to use a state estimator, such as the Extended Kalman Filtering program DEFKIS, mentioned in section 1.3.3. This technique may be able to remove some of the noise (or uncertainty) associated with those quantities to be- 206 -included in the extended control vector, used to drive the reduced-order model. It will, however, add considerably to the time and cost required to carry out an identification.Estimation Of Kalman Gains.For the state-space model given by equations (4.2-1) and (4.2-2) the corresponding continuous-time Kalman- Filter equations (which enable the expected-value quantities which are defined in (4.2-9) and (4.2-10), and are required for the maximum-likelihood formulation of the problem, to be obtained) are given by (e.g. GELB001, and APPENDIX 4):j? CO = A X CO + Byct) ■+' Kv'co ( 4 . 2 - 6 O )Z C O - H X C O + \fCt7 ( 4 . 2 - 6 DThe zero-offset terms are dropped from the aboveformulation; since they disappear anyway for transformation into the frequency domain, with thecondition that co^ 0. In addition, the steady-stateKalman gains K are used (transients associated withuncertainties about the initial state conditions havedecayed, and the noise statistics are assumed stationary).The steady-state Kalman-Filter representation of the- 207 -system given by (4.2-60) and (4.2-61) can be rearranged to give:yet") = A* yet) + b z (4.2-62)where: A Z= A - K H > B Z= (B,K) ; U,(t) = (U(t) , 7 (t) fBy transforming equations (4.2-61) and (4.2-62) into the frequency domain, the frequency-domain innovations vector V^ co) can be obtained for use in the probabilityz'*/vdensity function given in (4.2-17), instead of V (co) . In addition, the covariance matrix of the innovations 5/ can also be obtained from the Kalman-filter algorithm. Strictly speaking, the Kalman-gain matrix K, implicit in (4.2-62), is not an independent entity which, could be estimated as if it was part of a control dispersion matrix but is subject to the constraints of both itsdefinition and the matrix Riccati equation (governing the propagation of the error-covariance matrix of the state estimates through time) (e.g. GELB001).The above discussion of the inclusion of Kalman gains (and by implication, the inclusion of process noise terms) in the estimation, is not presented in detail, because no results were obtained applying the method. It was included to show the possibility of future work in this direction, though it would also be important to assess whether the considerable increase in complexity incurred as a result, would justify its use for helicopter system identification.- 208 -In the references KLEI004 and KLEI005, an algorithm for the frequency-domain identification of a linear system using maximum likelihood estimation, was developed. The algorithm is based on the transformation of discrete-time system equations into the frequency domain. However, the elements of the Kalman-gain matrix are estimated as if they were free parameters and are not subject to the constraints of the discrete-time Kalman filter algorithm such as the definition of the steady- state Kalman gain matrix, the error-covariance extrapolation and error-covariance update equations (see Appendix 4). (As an aside, it is also assumed in the derivation that the measurement window is periodic: theoutput-error method developed in this document incorporates a non-periodic measurement window correction term at the outset). Before the implications of neglecting the constraints and treating the Kalman-gain matrix elements as if they were ordinary parameters are discussed, it is worth pointing out that such an estimation could be carried out using the frequency- domain output-error method developed in this document. By formulating the estimation problem in terms of equation (4.2-62), and using the control vector extended to include the measured responses, and making use of the ability to include any relationships between different parts of the model, an identification could be carried out. There will be a substantial increase in the number of parameters to be estimated because of the need to- 209 -estimate the gain matrix K as part of the control dispersion matrix. Because this approach requires no modification to existing software, it is a possible area for future research.It is pointed out in the reference FOST002, that by using the simple method, whereby the Kalman gains are treated as if they were independent parameters, the parameter search will take place in a space with many more dimensions than necessary. For a time-domain maximum likelihood method, the author in FOST002 reports some success with the simple method, but recommends a more rigourous implementation J of the Kalman-filter algorithm, and develops some analytical expressions for the implementation of the algorithm.4.2.6 Summary Of Section 4.2.A frequency-domain output-error method has been developed and implemented in a FORTRAN-77 program: OUTMOD. The Gauss-Newton minimisation method, with an additional line-search facility, was used. A separate estimation algorithm was used for the elements of the error-covariance matrix.The output-error method was implemented with many unique and potentially useful features. The frequency- domain representation of the model equations, developed by the current author, facilitated the development and implementation of both the required expressions for the- 210 -estimation algorithm, and the useful features of theprogram OUTMOD. These include: the ability to definerelations between different parts of the model structure (i.e. a relational constraint); the ability to restrict the allowable range of individual parameter values; and the facility for obtaining rank-deficient solutions (linear relational constraints) using a singular-value decomposition of the information matrix.In addition, time delays in the observed control inputs and measured responses can be estimated, as well as elements of the model matrices: A ., 15 , H . The correction term for a non-periodic measurement window, discussed in section 2.2.4, is also incorporated in thealgorithm. The program OUTMOD, is to the best of thecurrent author's knowledge, the best and most versatile frequency-domain state-space identification software available to date.4.3 Applications Of The Fregpuency—Domain Output—Error Method To Flight Data.4.3.1 Applications Involving Longitudinal Test Inputs.The frequency-domain output-error method has been applied to a large number of data sets, although, only a representative set of results is included in this document. These examples serve to illustrate important features of the method and the computer implementation.- 211 -In Section 2.3.3, results were presented using real flight data for a longitudinal-cyclic input. It was shown that the output-error method produced estimates, which in the case with a time delay incorporated, compared favourably with theoretical estimates; the plots showingthe frequency-domain fits and the time-domainverification also indicated a successful identification. Figure 4-1 shows the estimates of the parameters, and the total cost-function value and its non-logarithmic component, as a function of iteration number. The estimated diagonal elements of the inverse of the error- covariance matrix as a function of iteration number are shown in Figure 4-2; these elements were fixed at specified values (obtained by estimating the relative noise levels using the original flight-data records) for the first iteration, to prevent possible initialdivergence of the algorithm due to the use of anestimated error-covariance matrix which could be greatly in error due to the initial parameter estimates being far from their true values.Twelve model parameters and four diagonal elementsof the error-covariance matrix (estimated using the axial iteration algorithm explained in Section 4.2.2) were estimated. It can be seen from Figure 4-1 thatconvergence is reached after about 16 iterations althoughthe exact number will vary according to the tolerancesspecified on the update increments for each of the parameters individually. This is of the order of one- 212 -iteration per unknown parameter. The total C.P.U. time for the run was about 4 minutes on a VAX 11/750 computer with V.M.S. operating system. The cost-function value shown in Figure 4-2 is progressively reduced after each iteration until it levels out to a value of about -131.4; and the non-logarithmic component (as explained in section 4.2.2) converges to the value of 64, corresponding to the product of the number of states in the model and the number of frequency points used in the identification. It should be noted that an additional factor of 2 is present because the cost-function minimised in the implementation here is not premuliplied by 1/2; this is, however, taken into consideration when the error-bounds are calculated. The non-logarithmic component should always reach the value indicated above, and is a useful additional check that the identification algorithm has converged correctly. The results shown in Figure 4-2, for the error-covariance matrix, also confirm that its elements have converged.A Transfer Function Example.An alternative to the formulation of the identification in terms of the linearised equations of motion in state-space form, is provided by a single-input single-output transfer function representation. The use of single-input single-output transfer functions valid over a defined frequency range for chosen flight- 213 -conditions, is an approach which has been used by some researchers (e.g. TISC001) using the chirp Ztransform. The transfer-function response for the measured data is obtained from cross and input auto spectra (as outlined in section 2.2.1). The approach adopted by the current author was to use the frequency- domain output-error method developed in this document: this is possible because the transfer-function models can be written in an equivalent state-space form.The classical pitch rate and normal acceleration responses to a longitudinal-cyclic input for the short period mode are given by:— Ste__ (s + yVe) e+ 2]CA<, pS + co?p (4.3-1)^ r -v , - 5 ITaxCLzCS) €.S* + 2]CJspb + GJ1? (4.3-2)where q(s)/f\ys (s) is the Laplace-transformed pitch-rate response to a longitudinal-cyclic input, and 5 is the Laplace operator. The parameter M^is is the longitudinal- cyclic pitch sensitivity, and T© is the effective time delay on the input for pitch rate. The parameter ~jf@ is given by (PADF006):T 9 - M A i* (4.3-3)- 214 -where ] and GOsp are the equivalent short-period mode damping and natural frequency respectively. In equation (4.3-2), the term Q z(s) //\vs (s) is the Laplace-transformed normal acceleration response to a longitudinal-cyclic input, whilst is the longitudinal-cyclic normal-force sensitivity. The effective time delay on the input, I ai, was assumed to be negligible for this investigation. The denominator parameters are identical to those in the pitch-rate transfer function.Equations (4.3-1) and (4.3-2) may be written in time-domain state-space form as:ito O 1 0 O iCt)i w -co?p -Z\cOsp O O i c oQ*(t) 0 O O 1 cuct)tact) O a,ct)o o o dis Ct-t®)MWTe Mfe O n.5Ct-te)o o o C b - Xaz). O ° *1* •(4.3-4)The identification problem is thus formulated now in a way that allows the frequency-domain output-error method to be used. Use can be made of the facility in the software to relate parameters within the model structure, because , 2?coSp and X© all occur twice. By specifying the equalities which exist between the elements in the second and fourth rows of the state- 215 -matrix in equation (4.3-4), we are effectively imposing equality in the denominator coefficients for the transfer functions shown in equations (4.3-1) and (4.3-2).Using flight-data for a longitudinal-cyclic doublet input (input and response time histories shown in Figure2-6), estimates were obtained, using the transfer-function formulation, for C3sp, , etc .. Resultsare shown in Table 4-1 for three cases, representing :(1) an identification using only the pitch-rate terms that is the subsystem of equation (4.3-4) corresponding to a state-space formulation of the transfer function given by (4.3-1); (2) an identification using (4.3-4,with the equalities existing between elements in the model structure specified; and finally (3), anidentification also using (4.3-4), but with the acceleration measurements modelled such that the actual acceleration is lagged relative to the measuredacceleration, with the acceleration time-constant estimated as an unknown parameter. Before the results are considered, some additional remarks concerning the third case would be of use. The acceleration time- constant is modelled as a delay in the corresponding measurement channels; this is a satisfactory approximation when \coXm \ « \ (see equations (1.2-17) (1.2-22) in section 1.2.5) .Comparing the results for cases (1) and (2) in Table 4-1, it can be seen that by using both pitch-rate and acceleration measurements to estimate the transfer-- 216 -function parameters (using equation 4.3-4), a substantial improvement is obtained in the parameter estimates, both in terms of their comparison with theoretical values and their estimated error bounds, when compared to those obtained for the single pitch-rate transfer-function. Finally, it can be seen that the case (3) estimates are in essence very similar to those obtained for case (2). Both sets of estimates agree within the indicated 1 crerror bounds, however, it is remarkable that 2 1 andM^is compare more favourably with theory for case (3) ,and have smaller error bounds. The acceleration time- constant is estimated with great uncertainty in this example, - in fact, a zero value for T™ , correspondingeffectively to case (2), lies within the indicated \crerror bound. However, the example shown here serves toillustrate a potentially useful idea which could be used in future applications for lagged (or phase-shiftedmeasurements.A noticeable feature of all three cases shown in Table 4-1, is that the time delay T© in the longitudinal- cyclic control is estimated consistently to be around 0.2 seconds, and has a relatively small error bound for cases(2) and (3) .It can be seen from equation (4.3-4) that implicit in the formulation is the assumption that measurements of <^ (t) and c\z(t) are available. In fact, for the example presented, they were not, and the problem of having to construct time-derivative signals was tackled by- 217 -differentiating a signal constructed using only the first few Fourier coefficients (corresponding to a frequency range similar to that used in the identification). The facility for differentiating a filtered signal over a specified frequency range is available in the frequency- domain output-error program OUTMOD. One drawback with this approach, is that, as pointed out in section 2.2.2, the ' filtered' signal might suffer from Gibb's- oscillation effects, resulting in some errors in the differentiated signal if the beginning and end of the response are very different. There is no requirement to differentiate a control-signal time history in the same way, because when the model equations are transformed into the frequency domain, the relationship between the Fourier transform of a given quantity and the Fourier transform of its time derivative can be used, as explained in section 2.2.3.The normal acceleration signal used in the identification, was effectively obtained from measurements of the incidence angle, using:W + W e ~ De.Wm*The frequency-domain fits (corresponding to case (3) in Table 4-1) are shown in Figure 4-4; and apart from the lowest frequency value, (possibly associated with non-linearities) the correspondence between measured and predicted value is good for and W.A Run With A Defined Non-Linear Relationship (Using A Combination Of Data For Two Test-Input Signals Applied ToThe Lonaitudinal-Cvclic Control).An identification was carried out using a combination of two longitudinal-cyclic inputs - a '3211' and a doublet. Both sets of data represent the same nominal flight condition: a speed of 100 knots instraight and level flight. The form of model used in the identification is the same as that given in section2.3.3. The data corresponding to each input is combined using the linear principle of superposition. This is possible in the present context, even when a non-linear feature, namely a time delay, is incorporated and estimated in the longitudinal-cyclic control, because we- 219 -can assume that the time delay is the same for both test- input signals. Some evidence for the validity of this assumption can be provided by the investigation carried out in section 2.3.2 using simulated data from the HELISTAB program, where, as shown in Figure 2-3, the estimated delays for both a doublet and a 3211, are identical within the limits of the indicated accuracy.The combined time-domain (averaged) responses are shown as part of Figure 4-5. For a linear system, whether the data is transformed into the frequency domain and then combined, or the perturbation data is combined in the time domain and then transformed into the frequency domain, should be immaterial; both approaches are attempted and are referred to as cases (1) and (2) respectively.A feature of the identification carried out here for the combined data set, was the specification of a nonlinear relationship between an element of the system matrix f\ and an element of the measurement transition matrix VI- The system matrix and measurement transition matrix are defined in section 2.3.3 by equations (2.3-11) and (2.3-12) respectively. The non-linear relationship defined was:A „ = ' / R » (4.3-7)where Aij and VVq are the Lj Yh elements of the matrices. It was assumed that A 23 was equal to the trim value of- 220 -the forward velocity component Ue , and that thederivative Z^ was negligible in comparison - from which, the relationship specified by equation (4.3-7) follows.The identification was carried out over a frequency range up to 0.60 Hz, with the transformed record length being about 13.3 seconds in length (see Figure 4-5). Thirteen independent model parameters, in addition to the parameter Ars defined by equation (4.3.7), were estimated. A third case was attempted where '/ Ue (and\Je) were fixed at the values indicated for the flightcondition.The parameter estimates for all three cases aregiven in Table 4-2. It can be seen that the results, on the basis of individual parameter estimates, are not as good as those obtained for the use of the doublet-input data alone (results shown in Table 2-2, case (2), in section 2.3.3). This is apparent for several parameters which previously had estimates which compared favourably with theory, and which had relatively small error bounds, namely: M^, M p, M^, Mq|4 and Z<\Xi . The parameter '/lieestimated for the relationship defined by equation 4.3.7 is, however, for case (1), estimated with a relatively high degree of confidence, and is in close agreement with what would be expected for a nominal trim speed of 168.9 ft/s (100 knots) . For case (2), the estimate of '/l)e is smaller in value; the reason for the difference is not clear, but may be associated with the strong non- linearities (i.e. large excursions from the trim)- 221 -apparent for the speed time history shown as part of Figure 4-7. In fact a comparison of the responses shown in Figures 4-7 and 2-6, shows that the combined run has much larger excursions in the measurements, which is not good for the identification of a model based on theassumption of small perturbations.The estimates of Zw , shown in Table 4-2, areconsistent for all three cases and compare very favourably with theory. The values of the time delay obtained are in agreement with that obtained for thesingle-input doublet run, within the limits of indicated accuracy.The frequency-domain fits obtained at convergence are shown in Figure 4-6, and, unfortunately, provide little useful information. It is apparent that thecombined run has resulted in responses which are poorly excited above 0.15Hz.. The poor excitation above 0.15 Hz., together with the large excursions in the measurements, are major factors in causing the estimates to be inferior to the single doublet-input case.When the one-second clock-period '3211' input (i.e. a total of 7 seconds in length) alone was used, theidentification algorithm diverged. It was pointed out in section 1.4.1 that the use of this input is far from ideal for the identification of small perturbation models. The effective input for the combined run, shown in Figure 4-5, represents an improvement, but is demonstrably inferior to the two-second clock-period- 222 -doublet input. The estimated auto-spectra for the two cases are shown in Figures 4-8 a), b).The time-domain predictions of the model, following a time-domain output-error estimation of initial conditions and zero offsets are shown as part of Figure 4-7.An Qutput-Error Multirun Example.In the reference MARC0 02, the combination of data sets in the frequency domain, corresponding to the application of test-inputs on different controls, is an approach which is advocated for multirun identification. This requires all the data sets - representing the same flight condition - to be of equal length. Time histories of the data used in a multirun identification presented in the reference MARC002 reveal that the test-input signals, applied on all four controls, are all of the same general shape. Although no parameter values are presented and final fits between measured and estimated time histories only are given, it is the opinion of the current author that this is likely to lead to identification difficulties, particularly for elements of the true control dispersion matrix. The reason for this can be appreciated by examination of the state-equation, where it can be seen that if two controls have similarly shaped test-input signals and are dispersed onto a given state, then the two corresponding row entries in the- 223 -control dispersion matrix are indistinguishable for identification purposes.A multirun identification was carried using two data sets; one corresponding to the application of a longitudinal-cyclic doublet (as used in section 2.3.3), and the other corresponding to a collective doublet. The test-input signals are shown in Figure 4-9, along withthe time-histories of the effective lateral-state responses incorporated as pseudo-controls in an extended control vector. The form of the model used is as given by equations (2.3-11) and (2.3-12) in section 2.3.3, except that an additional collective control term is present.If we consider Figure 4-11, and for the time being concentrate only on the measured responses of the longitudinal variables, it can be seen that theexcursions from the initial trim values are small in comparison to those given by Figure 4-7. This indicates that the small-perturbation modelling assumption is likely to be more satisfactory in this case. It may,however, have been more satisfactory to use different test-input signals rather than two two-second clock- period doublets, because of the reasons discussed earlier. The choice of data sets used was dictated byavailability.Results for the identification, which was carried out over a frequency range up to 0.60 Hz., are given in Table 4-3. Also shown are results corresponding to the use of an information matrix of rank 13 (14 model- 224 -parameters in total were estimated). The full-rank and rank-deficient solutions are referred to in the table as Cases 1 and 2 respectively.Examination of the results given for Case 1 in Table 4-3 for the multirun data, and a comparison of the results given in Table 2-2 (Case 2) for the single-run data, shows that for this example, the use of multirun data in the identification has tended to degrade the estimates in comparison to the single-run case. Estimates of parameters which had previously compared well with theory, and were estimated with a relatively high degree of confidence, no longer correspond with theory, in the majority of cases.The use of a rank-deficient information matrix resulted in an improvement in the estimates obtained in a few cases, especially for . The value of the time delay estimated in the longitudinal-cyclic control is in agreement with the value obtained for the single-run case, and is estimated with a relatively high degree of confidence. It is also interesting to note that a smaller value, which is less significant on the basis of the indicated error bounds, is obtained for the delay on the collective control channel. The investigation carried out in section 2.3.2 using simulated data, indicated that a time delay in the collective was likely to be of less importance than a time delay in the longitudinal-cyclic control, for the identification of a quasi-static model using flight-data contaminated with- 225 -rotor transient effects. A further reduction in the rank of the information matrix did not give, in general, appreciably better results: more physically realisticestimates of Mp and emerged, but this was found to coincide with less satisfactory estimates of other parameters, particularly for the time delay in thelongitudinal-cyclic control.The frequency-domain fits obtained for the rank-13 solution are given in Figure 4-10. The leastsatisfactory frequency-domain fit is obtained for the speed channel. A time-domain verification of the model, following a time-domain output-error estimation ofinitial conditions and zero-offsets, is shown as part of Figure 4-11.It can be concluded that if multirun data sets are to be used as a means for identifying more robust models than single-run data sets, then more consideration of the test-input signals used on each of the control channels is required. The use of multirun data sets for theattempted identification of a full 8-state model,involving the application of uncorrelated test-input signals on the longitudinal and lateral control channels is a task to which the identification software could be applied.- 226 -Other Applications To The Identification Of Longitudinal Models.Some other attempted applications of the frequency- domain output-error method to the identification of longitudinal models involved mostly the use of one-second clock-period 3211 test inputs. The unsatisfactory nature of this test-input signal for the identification of models relating to the Puma helicopter has already been emphasized. Nevertheless, despite obtaining in most cases poorer parameter estimates, some of the features of the more successful applications discussed here for longitudinal models, and in the next section for lateral models, could be seen in the results obtained. In particular, an improvement was found in the fits, and in some of the parameter estimates, for the incorporation and estimation of a time delay in the longitudinal-cyclic control. Results are presented in Table 4-4 for twoflight-data sets corresponding to a Puma helicopter,trimmed at a speed of 80 knots in straight and levelflight.- 227 -4.3.2 Applications Involving Lateral Test Inputs, With A Consideration Of Trends In The Estimates For Different Flight Conditions.Previously Reported Problems.The tendency for primary-rate damping derivatives, such as M<^ and Lp to be underestimated, is a recurring theme in previously published studies on helicopter parameter estimation studies (e.g. PADF002). Theincorporation and estimation of a time delay in the applied control for the output-error method is an approach advocated in this document for theidentification of quasi-static derivatives when using a six-degrees-of-freedom model; this approach wasdemonstrated in section 2.3.3 for a longitudinal model. In addition to the problem of the underestimation of roll-damping Lp, the results obtained in PADF002 using the equation-error approach, show that the accompanying estimate of the derivative Lr was very high compared to theory; this anomaly has also been observed by the current author, for a range of flight conditions between 60 and 100 knots in level flight, using the equation- error approach.If we consider the response (shown for the lateral variables as part of Figures 4-17, 4-18, and 4-19) to a pedal doublet input (Figure 4-12), for a Puma helicopter flying at a nominal trim sped of 100 knots, in straight- 228 -and level flight, altitude 6000 ft., it can be seen that there is a strong correlation between the roll and yaw rate responses. The damped sinusoidal roll and yaw rate responses are almost 7\ radians out of phase with each other, and are associated with a 'Dutch-Roll' type mode. It can therefore be expected that difficulties inestimating some of the parameters associated with these variables will occur due to the ill-conditioning of theinformation matrix. Ill-conditioning of the information matrix was considered in section 3.2.5 in the context of the equation-error method; and in section 4.2.3 thefacility for producing a rank-deficient - and therefore better conditioned - information matrix for the output- error method was explained. The presence of one, or more, relatively small eigenvalues in the information matrix is an indication that the parameter set used does notconstitute a truly independent set of parameters, andthat a better identification may be carried out using a rank-deficient information matrix; indeed, convergence problems can be experienced in the use of an ill-conditioned - nearly singular - information matrix.Before results obtained in applying the rank- deficient frequency-domain output-error technique for the identification of lateral derivatives are presented, it is worth noting that, in addition to the problems resulting from linear dependence, the information matrix will also be nearly singular, if one or more of theparameters to be estimated are weakly defined. In terms- 229 -of the cost- function surface defined in the parameter space, this means that in the vicinity of the minimum, the surface is relatively flat in at least one direction; a relatively insignificant change in the cost-function value would occur for a relatively large change in the parameter value. For the iterative output-error estimation scheme where a rank-deficient information matrix was used, the final estimate of a weakly-defined parameter will depend very much on the initial guess. Consequently, weak or insignificant parameters should be excluded from the identification at the outset; 'parameter significance' values obtained at the equation-error stage can be used as a basis for judgement.Pedal Input - 100 Knots.Consider the pedal-doublet run described earlier. The full-rank, and rank-deficient results, obtained using the frequency-domain output-error estimation technique are shown in Tables 4-5 a) and b). The results are also represented graphically for the important lateral derivatives, the delay, and the cost-function value in figure 4-13. Full-rank solutions with, and without, a delay in the control are also presented for comparison.In total, 12 stability and control parameters were estimated for each of the rank-deficient solutions. The number of time-domain points input to the estimation- 230 -program, and transformed into the frequency domain was 1700, sampled at 64 Hz., making a record of length 26.6 seconds. The frequency range used in the estimation was 0.03765 - 0.4894 Hz., corresponding to 13 complex-valued frequency-domain points, and was chosen on the basis of magnitude plots of the Fourier transforms, obtained at the equation-error stage. Initial guesses for the parameters were also obtained at the equation-error stage, except for the delay which had an initial guess of zero.The incorporation and estimation of a delay in the control, results in a substantially lower cost value at convergence. Figures 4-14 and 4-15 show the frequency- domain fits obtained for these two cases. Theimprovement obtained as a result of the delay is particularly visible for frequencies on either side of the peak at about 0.22 Hz.. In general, the agreement between measured and predicted frequency-domain responses is very good, especially the rolling and yawing moment fits. The delay itself is estimated to be about 0.2 seconds, and has a relatively small error bound.In the case of the rolling-moment parameters, the inclusion of a delay results in estimates that are in much better agreement with theory, than in the case without the delay. As the rank, of the solution is decreased to 9, there is a noticeable change in the estimates of Lv and Lr; L v agrees very well with theory, whilst the Lr estimate is much closer to theory than the- 231 -higher-rank cases. The roll damping Lp is lower than the theoretical prediction, but it is larger than corresponding estimates obtained from the equation-error approach (values of -0.9 and 0.8 6 for Lp and Lr respectively). The incorporation of a delay in the control has thus increased the estimate of Lp. The combination of the use of the output-error technique, the incorporation and estimation of a delay in the control, and the use of rank-deficiency in the information matrix - in particular, for an information matrix of rank 9 has led to rolling moment parameter estimates that are in generally good agreement with theory.Consider now the yawing-moment derivatives: therank-9 estimate of Nr is in excellent agreement with theory. Np differs somewhat from theory, but is estimated with a relatively small error bound. Nv is estimated to be larger than the theoretical prediction, but is still of comparable magnitude. The pedal control sensitivity to yaw N^p is smaller than theory suggests; however, the estimate obtained from the rank-9 solution is the closest to theory. The frequency-domain fits obtained at convergence for the rank-9 solution are shown in Figure 4-16. It can be seen that they are very similar to those for the full-rank case with delay presented in Figure 4-15.Following a frequency-domain estimation of the stability and control derivatives, the next stage in the identification scheme is to perform a time-domain output-- 232 -error estimation to obtain estimates of the zero-offsets, and initial state conditions, with a view to obtaining a time-domain verification of the model identified in the frequency domain. This stage in the identification methodology is explained in Chapter 5. This was done for the estimated model obtained in the following threecases: 1) full rank with no delay in the control 2) fullrank with delay in the control and 3) rank-9 solutionwith delay in the control. The time-domain verification results following from the time-domain estimation are shown in Figures 4-17 to 4-19. First comparing Figures 4-17 and 4-18 for the full-rank solutions: it can be seen that for the roll rate channel in particular, theinclusion of the delay leads to a much tighter fit over the first few seconds of data, when the control input is applied; the rank-9 solution also shows this, and in comparison to the full rank case in Figure 4-18, the time-domain fit is only slightly degraded towards the end of the time record.The preference for the rank-9 solution was based on comparisons of the predicted theoretical values with corresponding estimates. It is accepted that all the parameter estimates obtained from flight data need not equal the theoretical values, since the purpose of system identification in the current context is both to confirm some aspects of the theoretical model, and to update others. However, important primary effects should be able to be predicted by relatively simple theory, and so the- 233 -estimates of parameters strongly influencing these effects may be used as an indicator of how good the model is, along with the time-domain reconstructions and predictions of the model.In going cfrom rank 9 to rank 8 there is asubstantial degradation in the estimates of most of the important parameters, such as: Lp, Lr and Nr . Figure 4-20 shows for one iteration of the frequency-domain output-error method a typical reduction i n cost-function value (normalised to the £ull-rank case) that wouldobtained as the rank of the information matrix isincreased from 1 up to the full-rank case; this isobtained using the total differential:A T =* a? A'C, + ZJ A.Yi + M A'A. + ... 3 7 •a-e, a'Cz b 'Cj ■a'Cp(4.3-8)There is clearly a distinction between the rank-9, andrank-8 and lower-rank cost-reductions, whereas the rank-9 reduction is of comparable magnitude to the higher-rank cost-reductions. This observation seems to reflect the degradation observed in the parameter estimates for ranks lower than 9.Longitudinal measurements were included in anextended deterministic control vector for the four-state lateral model used in the estimation, with the significant cross coupling terms identified as elements of the control dispersion matrix B. Table 4-5 b) shows- 234 -the values of the estimated lateral-longitudinal crosscoupling terms. In regard to the agreement with theory for these parameters, it can be seen that whilst there is not even any approximate matching, there is some evidence that rank-deficiency produces estimates which are at least of the correct order of magnitude. It can be appreciated that the degree of coupling-intensity between lateral and longitudinal states, in the 6 degrees-of-freedom (6D0F) model, will determine the ease with which satisfactory estimates of these parameters can be obtained. As pointed out earlier, the use of rank-deficiency in situations where the information matrix is 'near singular' because of the presence of one or more weakly-defined parameters, can result in situations where the estimates of the weak parameters are dependent on initial guesses. Indeed, for the current data set, there was found to be some evidence that for the cross-coupling terms, the rank-deficient solutions were dependent on the initial guesses, although estimates of the expected order of magnitude were still found. A contributory factor to the disagreement between estimates and theory for the cross-coupling terms could also be due to the fact that in directly incorporating longitudinal measurements into the extended control vector, any offsets relative to the centre of gravity of the corresponding measurement devices were not taken into consideration. In the case of the measurements that relate to the states in the linear model, offsets are accounted for through the measurement transition matrix- 235 -Y\. For the extended control vector, this would requireextra measurements to be included in the control vector, and extra parameters to be estimated in the 3 matrix.It should also be noted that the use of noisymeasurements - such as C/(t) and C^(t) - as deterministic pseudo-controls is a possible source of error in theestimates, where in the estimation algorithm there is the inherent assumption that these are noise free. It is assumed that there is no process noise on the model. Larger models incorporating all longitudinal and lateral states would avoid the problem of noise on the controls, but would mean the estimation of a larger number of parameters, and would require the use data sets generated from control inputs that excite both the longitudinal and lateral modes.Lateral Cyclic Input - 60 Knots.Consider now the results obtained for a lateral- cyclic doublet input, for a Puma helicopter flying at a nominal trim level of 60 knots, in straight and level flight, altitude 1000 ft.. The input is shown in Figure 4-21, and the lateral response variables are shown as part of Figure 4-24. The length of record available for estimation is much shorter than in the previous case, with 800 points transformed into the frequency domain, sampled at 64 Hz., making a record of length 12.5 seconds. The frequency range used in the estimation was- 236 -0.08 to 0.56 Hz., corresponding to 7 complex-valued frequency-domain points,.The important lateral stability derivatives obtained from the frequency-domain output-error estimation are shown in Figure 4-22 (a), together with the estimates obtained from the previous 100-knots case, in order to clearly visualise any trends that may be apparent in the estimated values. Theoretical HELISTAB values are also shown for comparison. Error bounds are shown only for the rank-9 case in order to avoid the figure becoming too cluttered; the error bounds for the other cases are of a similar magnitude.Concentrating first on the rolling-moment parameters, it can be observed that in the case of L r the substantial improvement in the estimate that was obtained in the 100-knots case, as a result of using rank- deficiency, is not repeated for the 60-knots case; the lower-rank solutions are, however, smaller than the full- rank case. Once again, the estimate of Lv is improved by rank-deficiency, with the rank-10 and rank-9 estimates in excellent agreement with theory. The rank-10 'and rank-9 estimates of Lp, whilst as with the 100-knots case are lower than theory, are consistent with the magnitude order predicted.by theory.For the yawing-moment derivatives, there is close agreement with theory for estimates of Nr , with the rank- 10 and rank-9 solutions. A value higher than the theoretical prediction is obtained once again for Np ,- 237 -though the value obtained for all the ranks is on the whole larger than the 100-knots case, and this trend is predicted by theory. The Ny estimate is in excellent agreement with theory for the rank-10 case.In Figure 4-22(b) are shown estimates of the lateral- cyclic control sensitivity with respect to roll rate; the estimated delay in the control; and the final cost- function value at convergence. The rank-10 and rank-9 estimates of the control sensitivity are identical within the range of error. The estimated delay is not as large for the lateral-cyclic input as it was for the pedal input, and is not estimated with the same degree of confidence, although it is of a magnitude comparable to the time constant of the main-rotor longitudinal and lateral cyclic flapping modes.It appears that for the lateral-cyclic case at 60 knots, the rank-10 solution gives the most satisfactory agreement with theory. The frequency-domain fits obtained at convergence for the rank-10 solution are shown in Figure 4-23. In Figure 4-24, the time-domain reconstruction is shown for the rank-10 estimates, following a time-domain output-error estimation of the zero offsets, and initial state conditions.Tail-Rotor - 60 Knots.Data obtained for a pedal-doublet input, for a Puma helicopter flying at a nominal trim level of 60 knots, in- 238 -straight and level flight, altitude 1000 ft., was also analysed. The input is shown in Figure 4-25. As with the 60-knots lateral-cyclic case, the length of record available for use was 12.5 seconds; this is considerably shorter than the 26.5 seconds of data available for the 100-knots tail-rotor case. The magnitude of the doublet input, however, in this case is larger than the 100-knots case: the result being that the corresponding excursions from the nominal trim levels are large. This is not good for the estimation of a small-perturbation model, and highlights the important point concerning control-input design made in Section 1.4.1: attention should beaddressed not only to the shape or frequency content of any applied input signal, but also to its amplitude and the magnitudes of the excursions likely to be produced.The frequency range used in the estimation was the same as the lateral-cyclic case - 0.08 to 0.56 Hz.,corresponding to 7 complex-valued frequency-domain points. Estimates of the important lateral stability derivatives obtained from the frequency-domain output- error estimation are shown in Figure 4-26. The full-rank case failed to converge, but by turning to rank-deficient solutions, convergence was obtained for the output-error method.The results shown in Figure 4-26 also include the estimates obtained for a full-rank solution when the Lr parameter is considered to be linearly related to Lp; the Lp parameter is estimated freely and the Lr estimate is- 239 -constrained using the theoretical HELISTAB ratio of the two parameters:Lr = L p (4-3-9)The result is that there are 11 free stability and control parameters to be estimated, and one additional related parameter which is updated at each iteration; the sensitivities are calculated within the output-error algorithm taking the defined relation into consideration. The ability to define relations between sets of parameters, is one of the features of the estimation program OUTMOD.It can be seen that the cost-value obtained for the case with the relation between Lp and Lr is almost identical to the rank-11 cost-value obtained at convergence. In addition, th^ parameter estimates shown are identical, within the bounds of accuracy. For the rolling-moment parameters, there is good agreement with theory for Lv , Lp and Lr, in both cases. In reality, the correlation is not usually between pairs of parameters, but may involve a large number of unknown parameters, and so the technique of fixing relationships between parameters is not a practical solution to the problem of correlations between parameters. The example shown, however, does perhaps reinforce earlier statements about likely problems in the estimation caused by strong correlations between the roll and yaw responses. The fact- 240 -that for the results presented in the two previous cases, solutions of rank-9 and rank-10 respectively, gave the best estimates, where the full-rank case was of rank 12, does indicate that the existing correlations were indeed between more than simply Lp and Lr.If we consider the yawing-moment parameters, it is seen that they are not estimated very well for the rank- 11, and full-rank case with the defined relationship (i.e. rank 11), but require lower-rank solutions in order to approach the indicated theoretical values. On the whole, the estimated delay value is very similar to that obtained for the 100-knots tail-rotor case.Comparing the cost-function values given in Figures4-26 and 4-22 b), it can be seen that for the 60-knots cases, they are, for all the rank-deficient solutions presented, much greater for the pedal input than for the lateral-cyclic input. This is reflected in the better parameter estimates obtained in the latter case. It should be noted that the larger the cost value obtained at convergence, the poorer is the fit, and that comparisons of this type are valid here because the same number of frequency-domain points were used in the estimation process in both cases. It is felt that the large amplitude of the pedal input, and as a result, the less satisfactory adherence to the small-perturbation modelling assumption of the responses is to blame for this.-241-CHAPTER 55.0 FINAL STAGE OF THE IDENTIFICATION METHODOLOGY AMD ASUMMARY OF THE COMPLETE IDEMTIFICATIOM APPROACH.Estimation of Zero Offsets and Initial Conditions.5.1.1 Introduction.By excluding the zero frequency value from those values used in the identification, the estimation of the stability and control derivatives, and time delays, is carried out independently of the zero offsets and initial state conditions. However, for a time-domain verification of the model identified in the frequency- domain, it is necessary to estimate these quantities. This is seen as the final stage of the identification methodology, and is carried out using a time-domain output-error method. Those parameters estimated using the frequency domain are considered fixed for the purpose of the time-domain identification.5-1.2 Minimisation of Time-Bomain Cost Function.The Time-Domain Cost Function and Model.The time-domain maximum-likelihood cost function is well known (eg. FEIK002). For the implementation of the-242-method by the current author, Kalman-filter innovations are replaced by output-errors; this is analogous to the simplifications made in the frequency-domain application of the method, discussed in Section 4.2.1.T t = ' A 2 ( Vet) R':Vet) + loje.lRl) (5.1-1)V(t) are the output-errors and R is the time-domain error-covariance matrix. The state equation and the output-errors are given by:i o o - AXct) + B C u c t , x * ) - k h (5.i-2)Vcf) - HXct) + k m - (s.i-3)^ / s r—where T c and are the set of time delays in thecontrols and measurements respectively; Kc , are the zero-offset vectors in the controls and measurements, and together with the initial state conditions: Xcd) , theyrepresent the unknowns in the estimation.It is assumed that elements of the state matrix - A , the control dispersion matrix - B / the measurement transition matrix - H , and the sets of time delays T , T are known, or have been estimated using thecfrequency domain; these are fixed throughout the estimation.The discrete representation of (5.1-2) and (5.1-3)may be written as:-243-)(Cn+\At) = ^ XCnz^t) + D ( L K n&t v Tc ) ~ V O(5.1-4)Vtn^t") = HXCoak) + - ZCn^ t, )(5.1-5)where the following are defined (e.g . GELB001):T /t . A ACtz-tOQKti.W) = e (5.1-6)X Aq) = e (5.1-7)n+i*tD UCn^t, Tc ) = 5 5) (n+\ At,©) B U ( e 3 T c I d ©nat (5.1-8)For small At we may use the following approximation:.D - Bat (5.1-9)3k A AtAnd for (t) = Tm ) ^ H X Co) Km (5.1-17)V ( n & O ^ Q (5.1-18)As an example, consider the measurement system for longitudinal quantities, described in Section 2.3.3, evaluated at time t=o, and assuming Vcn&tO- O :Vm Co) 1 o i ; 0 UCcD VoFF(Y„C0) — 0 ‘/Ue “15/ Ue o W Co) + 0(OFF 4-J0 ou. 'O 4-> 0 CO J2 *->H co03oo03 CO CDin LO t-H o o-o CO o 1-—1 coo CO o 02 o« • • • •o o o o oCOo- 03 CO LO _ _ „— .co CO t" CO -— . 03) CO -— - CO 'Ci-i LO03 O CO CO o LO 03 CO o o CO 03o O o o CO O' o o LO CO CO o 03o o c o T-* 1 o o ’— I o o 03 LO# # # • . . • • • • • • • •o o o1 o o1 oo1 o o1 o o1 oo COLO03 LO .—. -— - CO .— . -—-oo CO 03 CO o oo o -—* LO O- T^l03 o 03 o o co 03 CO 03 03 03 CD O'o o o o CO CO o o 03 O- CO G CD COo o o o CO o o o 03 G o G 03 03« . • . • • • • • • • • •o o o o o o o o c G G G o LO" ^ 1 1 1 1 1 ICO ._ „ .—« ,—* 'cr1 03LO 03 o CO LO 03 -—- G .—i -—~03 03 CD 03 CO o LO T—1 G O LO O- o03 o 03 o o o CO CO G G CO co G COo o o o o o o o G O 03 G *—< 03o o o o o o o G G O G O 03« • • • • • • • • • •o o o o o o o o G G G G G CD1 -- "—" 1 -■ I — I ,-HCOo03oCOm0303LOoCOCOCOCOCOoCOCOo% I03-263-Table 3-1 b) Parameter S i g n i f i c a n c e Values and P a r t i a l - F S t a t i s t i c s for the S e l e c t i o n o f 5 Orthogonal Components.Parameter S i g n i f . Value Part i a l - FMi 0 .2 4 2 0 .0 0Mv 0.21 14 .07Mq 0 .2 5 14.07Mv 0 .1 5 56 .21Mp 0 .1 9 14 .10Mqls 0 .8 5 147.80Table 3-2a) Equation-Error Results Pitching-IVfament Equation Puma, 100KN, Run R0201A, Long-Cyclic Input. Estimates of Orthogonal Variables with Singular-264-CO 0 +-> b£> CO CIX (1I—CO O—< c4-J 0s- 3co cr CL 05-hX! U-, 04-> *,CO COco bO■3JSCM CMt o COc o CO— t o CM CM o• *d o oin CO 1fH COi n CMo d doCO o(M CO 05‘S C *odCOoiCM dcoCM CO COO t - oo.c o o • oo .diCO oCO o 05o c o ^_ cm o • in0 s CM .diCO oCOino •- ,“H o CMC k • CO »o1oCO CO II0 0 -— -+-j 0 X3 co 3CO X.-H CO COS-l 4-> 0 > -— -CO CO 3 CO> w — * 0o •> c XCO CO oc u o S-4o CO •*■H CLbO <4-1O OS •f-H .X bO c CO4-> c bO -Q-— - S_ •!—1 •-4CO o X l CO COTable 3-2b)-265-cn UOCM CO CO CD CDr - CO CD uo O CDCD • • • • • •CM CM i CD Oi— i CD CM 4— 1OCO CO CD ot"- c o CD UO 4-HID • • • • •CM CM CD*—i CD CMCOCDO CD CD l>-ID CD OO« • • •UO CD oCO UOCD "^4 CDCO CD TTco CD• • •o CD ot>o ▼—(CM o CO• •t^ -c oCDCO COCD o1— 4 •f-HJD 4->CO CO• r—1 •1—1J-, ■l_>CO CO i— <> 4-» COC/3 C1— 1 OCO U hibO COc O CDo —i -C — IbO CO +j -OO 1—4 5-4 COd: 4—> o — 44-> }-. uU< CO *4-1 coo CL o >C O uo C D-266-Table 3 - 2 c ) M atr ix0 .6016 0 .7959 0 .0 0 4 7 - 0 .0 5 6 7 - 0 .0 0 3 2 0 .0 3 8 5- 0 .7 9 5 6 0 .5963 0 .0 0 0 3 0 .0022 - 0 .0 0 2 8 0 .1 0 7 3- 0 .0 2 0 7 - 0 .0 6 4 4 - 0 .0 0 4 6 - 0 .9 7 1 9 0.0061 0 .2 2 5 30 .0691 - 0 .0 8 2 7 - 0 .0 1 2 7 0 .2 2 8 4 0 .0 0 0 4 0 .9 6 7 50 .0 0 1 5 0 .0008 - 0 .7 4 3 2 -0 .0 0 2 9 - 0 .6 6 9 0 - 0 . 0 0 8 9- 0 .0 0 1 0 - 0 .0 0 7 0 0 .6 6 8 9 -0 .0051 - 0 .7 4 3 2 0 .0 0 9 8Table 3-3a) Equation-Error Results. Pi tching-Mment Equation (with Delay in Control). Puma, 100 KN, Run R0201A, Long Cyclic Input.Parameter EstimatesFrequency Range: (0.075-0.527 Hz.)- 267 -iCOCO 4—>0) CO£55oCMooCM CO CDin in o o-o CO o COo 00 o CM o• • • • •o o o o oCMCMCO rr - -■r— 1 oo O CM CO 03 t> CO 03D o O CM oo in -—« 00 t^- COr*- o o O o CO o o CO t"- CO o in oo coo o o o oo c o o o CM ▼-H o o —'H 03 COa # a • . • • • « • • • • . •o o o o o o o o o o o o o O CO' 1 1 1 1 COt" CM —H CO CD —-^ '—-co 1-H CM in *—t CO oCM o o o CD CM 03 -— - o CM o in ooco o o o o o m o o CO CO o t"- CMo o o o CM o o o in »— 1 o o o 03 r-a • a . . • • • • . • • • • •o o o o o o o o o o o o o o ini 1 1 1 1CO-h~. .— . in oo 03 COCO t^ - CM in CM ▼-H - - 1 CM05 05 CO o CM CM o o oo oo CO 03CM o CM O o o CO c o o o 00 CO oo coin o o O O o o o o o o CM o in » H -"tfo o o O o o o o o o o o o 03 in• a a a • • • • • • • . • • •o o o o o o o o o o o o o o CMi 1 1 1t*- in1 c o• • • • • • • ooo oo mT—t COc o • • • • • • m ■rP 00O ooC M ooC M• • • • • • • ooO oCM03o oo o• •o oxoo*£coCMPC 4“ -If- °Estimated lcrerror bound.See Section 3.3.1 for explanation of how approx."^ value is obtained Selected number of orthogonal components (=7)-268-T ab le 3 - 3 b) Parameter S i g n i f i c a n c e Values and P a r t i a l - F S t a t i s t i c s for the S e l e c t i o n o f 7 Orthogonal Components.Parameter S i g n i f . Value Part i a l - FMi 0 . 3 4 22 .8 6My 0 .0039 0 .00038Mq 0 .6 5 6 .6 6Mv 0 .1 7 7.71Mp 0 . 1 7 2 .3 3“ 1* 1.01 68 .9 2Table 3-4a) Equation-Error Results. Pitching-Mcment Equation (with Delay in Control), 100 KN, Run R0201A, Long-Cyclic Input.-269-3CO — 'coc- oco bOCD O3 X— 4->CO> ou <+-■CO o3bO co3 s-0C O X4->- cCO 0 -—0 S-. .— 0 NXI *+* JLCO <+*• «1—1u. Q C MCO in> .oo—< c-l—I iCOc inO CONbOOOO —' •X 4J O4-> CO■— -Si •-!o «CO 0■m bOi+i c/3 3O COVCCO [j_,>,0 i o— iceg co 03—H +J cr+J l-l0CO CO SiWCLUhCOCD o*—1• •o oCDO©CD +—1in +—i COin CM• • •o o oCOco oo XCM in Ciin O • CD'CL •o CM OCMCM CM CDo O CDo . Oo .o1 CO OCio COo CO^ CO o • oo .©1oCOCDo CD co_ CM o . in©1CMC3•ot'- CD 05o . CMo CM©1CO ©CO CO I0 0•M 0 XX eg 3 X!3 3 i—i• H CO 3u 4-> 0 > -—-3 CO 3 3> w 3 0O_ > 3 TO3 3 O3 u O S—iO 3 •—< O hbO ■— i 3 bO -9-—- i-, <3 O X CO COTable 3-4b)-270-t-Ht-" CO CO CM 1—1t'- CM c- ’-H o• « • • • • •05 C o *■“4 CM»-H CO 05 CM1t> 05CO 1 "■ < LOCO uo CM LO CO• • • • • •05 o i H 05 uor—1 05 COcoLO coCO 'cT' *-H 05 oLO CM "CT CO• • • • •LO co o o COCO LOCM 05CO CM COco o o• • • •o CO o oCOCM l"-CO oo co05 CO o• • •o CO oCMCOCOo CM• •1— 1COoco COCO CO•f-4 •1—1s_ 4-1CO CO —•> 4-1 COin cC+-I « oO CO [l, M coc l O 0co O •—i jo —u b£) CO 4J X5CD O —' J-1 CO-Q 4-> O —H 4-1 S-. Si5 S-i CO ^ COz o O, o >’-I CM CO UO COTable 3-5a) Equation-Error Results. Pitching-Mment Equation, Puma, 100 KN, Run-R0201A, Long-Cyclic-271-cxEin%oE0EO'0S-HUhcoS-, •o CO<4-1 •l_>E0CO E0 O3 &CO£> U_ _ ,S-4 30 E■u> O05 bOE O3 JEU 4—>3 i-,CLOF tot71.345CMCM»—1o o o od oo'CM t>CO Od dCO 00CMCM OO OO Od o1in cmoo inc m oo o o o• •o o.X)E30 0 &4-> 4 J cSE E£3 •E oS-H 4-> i->3 CO S-hX w u- 272 -Ou-.ou4-»C3CO&x:NXC Oint^-oiint^-o.oinicoo£CO 054->o 054-> • •U-, 05Oc—oi ,od 05•oinXo. .D. o£— ■CO£ CO COfo o ooI oO COCOo o05 —Tf COin *-hJ > o oo o• •o o^ — ■t'- COft>— CMo o05uo —»in t>-o *- <2mE H Uro .5 ou ■u> uco C/5 l_Oh w wSee Section 3.3.1 for explanation of how approx. X value is obtained.Table 3-6a) Equation-Error Results. Nanral Force Equation, Plim, 100 I 00CO • in• CM •o 05CO^aooCOCM .CD 00. CO •O I COooCO CO I 1-1 .— .CO CO G CO> w 0CO o •> c X3CO CO Oc Si o Sio CO CL,bO ■—1 <41O a •1-1 .x: bO G CO4-> c bO -9Si 1—1 •pH • i-S 0 0 0 0 » - < 0 0 0 0 c 0 . —' CO 0 5 CMo c d o o o o o o o o o c d o o o o q1 1 1 1 CM1_ [>. -— - CO -— - * - 1 -—' -—- t -~ V - - -— - CO CO - H c c OD Q COC— 0 5 ’—i CO '—i CO -— ' 1 ^ CM O CM CO 0 5 — <8 1 P 3 8 8 2 3 8 8 8 8 S 8 8 f e £o o c o o o o o o o o o o o r p o c >CO0.973831.918CM *0.972229.9310.82444.024Number ofOrthogonal Components:c/5Vi 0 Xl£ CO 33 Ch — i X• <-H •1—1 CO 3l-l ■M 0 > -— -ca C/5 3 3> u 0ca O •> c XIca 3 oca i_ O s_o ca QhbO •—H <+Ho 3 —H •x: bO 3 CO■i_> C bO -Q-- * Ih s tf-H i-iCO3bO355-o0)4-1COooc/3c/DCOcobOO43+->S-HoCM o oCM CM Oc o 03t— l> CMS^L o d oCOCM t'- 03CM o oto T—Idi d dCO o 0303 o o !>•ID' a di CM oCMC"- CMCM CM oo • 03'Q . oiI>CO doCO oo o o too • CM'Ql to •o1o o ot"-COoo CO in^-*cm o • of"o co di y~“1co03 co 03O • l—HO to' a t •diCO COCO CO 1!0 0 .— .f—1 4-> 0 JOX> co 3CO 3 X•i—< COJ-H 4-> 0 .— -CO CO 3 3> DJ 0CO 0 •> 3 *0CO 3 oc S-H O uo CO CLbO <4-,O 3 •— I •43 bo 3 CO4-> c bO -9S—i ___ • f-H ing £ s 2 g g g g 8 8 Q g 3 8 $ S w 8 8 8 8 8 8 8 8 8 8 8 8 8 8 3 S 8O 3 O 3 O 3 O 3 O 3 O 3 O 3 O O CM0.916428.196co0.916228.114CM0.910826.250-0.910626.192Mrrber ofOthogonal Cmponents:CO 4—>I Ha . 1 12 2 f 2 2: 2 £ l i t 4- aEstimated lcrerror bound.Selected nuiber of orthogonal ccnponents (=7).-283-T a b le 3 - 9 b) Parameter S i g n i f i c a n c e Values and P a r t i a l - F S t a t i s t i c s for the S e l e c t i o n o f 7 Orthogonal Components.Parameter S i g n i f . Value P a r t i a l - FNu 0 .0 1 9 0 .4 0 0Nw 0 . 1 0 11.53Nq 0 .0 4 9 31 .3 7Nv 0 .7 6 185.6Np 0 .2 3 7 .187Nr 0 .3 7 35 .55Nrjp 0 .3 8 30 2 .9-284-CO in oCO »-H CMoo COoo CM r-H'S d o•ot>- CM CMt^- inGOo- in CMdi ©Ooo COCM OOCM CO »—1CM CO in'a O CM ocoGbJOGOOX3CD4-1COooViinCOcobOOJZ4-45-HoCM COCOCM CM oin O • 05^ a oir-CO dCOCO COo COo CM ino • o^ a d oo oooo inCOo oo coCO o * CMdiCO00 dooCO CMo CMo CO inCM o • o^ a oit>-CO ot '-co05 in i HO • CMo COoit"-CO COVi C/3 iiCD 0 -—-4—> 0 X!X £9 GCO G i—1 X—h •-H CO COs_ 4-> 0 > H—HCO CO G G> U 0CO o •—H > c TDCO G Oc S-H O i io CO DubO —1 <4-1O G —4 •XI bJD G CO4—> G bO -9._ S-H __ .*—« c/d C/DTable 3-10a) Equation-Error Results. Fbl 1 ing-Mrrent Equation (withEtelay inGontrol). Flrra, 100 ICOcoco oCM CD ot CM 00in o • Mf"Q. dt-co dCDo oo00 CD oo CM COo • CD' S L 00 *d o o05co COCD 00o 00 05CO o • in73. oCD00 oCMCMCM CO ooCM o • 05TQt d1t"-co CMCOo05 in CMo . 00o 'Cf COot>CO COco C/5 IICD 0 -— -4-1 0 JOCO 3CO s —- X■—< •—1 CO 3S-. 4-> 0 >CO C/5 3 3> w 0CO O .— > 3 "OCO 3 oc S-. O l-lo CO •*H D-,bJO c+1O 3 ..-1 .x: bO 3 CO■i_> 3 bO -9- i-i .-- •*-« % CM CD o co COs- CO in 05 o o CMO o o CO o o in o * •K-0 . • • . . . .So1CO1o o o1o1o14 ^ CO CO 4—4CO o »—4 -—- 00 O CO "CT 4—4 — ^ 4-—oo ^ o CO 1-H in 05 in r - o CO — —H t " fH rH CO COcs o —< -sF CD CO co o o o in ■'cf t '- 1 CO o CM CO .o o CO o 05 O o o o o CO o CD O o o CM O CDC • • . . • • . . • • • • . . • • • •cq o o -H O o o o o o o o o o o O O o o& 1 s— 1 '’ w ^11 '' .1CD 4—4 iH iH 4----oo o .—4 -—- 'sSH co ^ CD 4-— ----- CD CD 4-—05 co o — CO o - in CO O CO CD in co CM — O 05 —HCO o 05 — o O O co in cm in CO O — CM •o o 'cf O o o o O O CM o in o O O CM O CMC • • » • • • • • • • • • . . . . . . ino o -H O o o o o O O o o o o o o O o —iC-I 1 ■" 1 ' 1 s—^ 1 '44'1i4"— .—. CO t> 4—4o CD CO 4-—-V in co coco 1-H ----- 4—4 CO CD ----4* CO *—i 'M1'—- CO 05 05 in CO o co ^ LO t"~ CM — CD — i"~CO o CO 05 CO 05 — o o o CM C- co in CO O CM CO .o o LO — in o o o o o • o in o o o CM O COG • • • • • • . . • . o • . . • • • • ino o .-4 O o o o o o o i o o o o o o o »"H& 1 1 ^ ^ 11i4—4 4—4 \r~ co 4—41—H CO CD 05 CD CD CD 4.—V ----- 05 ----41 00 ^ OO.— . CO — CO 'si’1 CO o Tjl CO co t"- —H H—4 TCis — 05^ o o CD CO o o o in cd — CD CO O O CO .o o Xf* »-H C- »—• o o o o CM O in o o o CM OC . . . . . . • • • • • • . . . . . . CDo o o o o o o o o o o o o o o o o —1 ' 1 44—4 . — 1 4 4 1 4 4 1c -—- 4«—*4 in o -----CO oo x f in CD CD 4—4 4—4 CM rcf ----45- CO 1-H 0 5 — 4 [4-4-4 CO ^ CO O C- CD O CM —H —H "M* — 05*—1 o CO ^ 05 CO CO o o o CM CD O CD CO o O CO •— CO o o in — O • o o o o CM O in o o o CM O CO1—• 1—1 • • • • . . . . • • . . . . . . • •o o —• o o o o o o o o o o o o o o —HL . d 1 "4—' 1 1' ' 1 1 4—4 1-f-C -—4 LO CD 4—4CO CO CO CO in cd in 4—4 ----- inU —i r^ — — in o CO c m CD ’- i CD0 1 o — CO CO ^ o o o in cd 05 CD CM O * .— TJ o o CO — CO — o o o o CO o ^ c o O O or—H • • • • . . • • • • . . . . • • CDo o — o -H O o o o o o o o o o oIj-i 2 1 s—' 1 ' 1 1 4441l5-404->O-h£ “0 4-> 35h a 0-4 £ C/5 —COOh LvLp LrcrNv§ Nr £ P Q g H—Estimated lcrerror bound-295-JQinij=i-0.0460.339-0.328o O 03 o c^o in o CMcm o CO O o» • • • • •o o o o o oi 1CM m CO 00 03CM o CM oo in03 O O o« • • • . •o o o o o o1 103 o COCO © o CO o oCO CO W' 1A • • • • •o o o O o o1 ’ 1t'- oin t- 1— 1 CO T— l CO03 t—H iH 1-H CM. • • • • •o t-H o o o1 1o CMo 00 CM o' o CMt" CM co CO CO t• • . • • •CO o CM o o o1 100COin CO CM CM03 CO CO o ▼-H. • . • • •CO o CM o o o1 1*nJ LqNqKEY TQ PLQTS.F.EE.QHEKCY,~P.QM Ifi... ,E US.(Magnitudes of FFT’s are used to display fits) Measured : —Predicted : XTIME-DQMAIM FITS Measured :Predictedurns.Angular measurements are in radians. Speed is measured in Knots.State-space angular quantities are in radians. State-space translational velocities are in ft/s.- - - - - - - Measurements+ + + +-+• + Estimates without bias correction Estimates with bias correctionX18W(n /s)l—011V( f t / s )0w10 1u(«/•)17T T T - n " [ T i *T19f f X ‘ii-1iffy A 1 _ - “ l c~r'Trr7T7T_T"T 1— r t— 7-7 r— n^TT— 1— i— 1— 1— t— 1— 1— r— 1— 1— r t— 1—(rad) 0 j -j— — |--i— 1— 1— 1 - r ~p^f*T t t20X 10<#>( r a d )0I-1-4:f™— |~“f—F~r T T H n - n - r --- r~~r— j— t— r iferr'T2yjf-p-j— r_T~p=S--r6 3'4T 107*1- r16 13 20FIGURE 1-1 DETECTION AND REMOVAL OF BIASES IN ACCELEROMETERS USING THE PROGRAM KINECON1 SEC 3211 1/2 SEC 3211A !■»-0.4TIME (SEC)0.8.40.00•0■08TIME (SEC)2. 1.30.30.50FREQ (RADS/SEC)987G532G 90 3FREQ (RADS/SEC)FIGURE 1-2 TIME HISTORIES AND AUTO-SPECTRA OF 1/2 SECOND AND 1 SECOND 3211 INPUTS1 SEC DOUBLET 1/2 SEC08.400 0■0.4•0,8TIME (SEC)0 80.40,0■0■0 8Q_COoTIME (SEC).5000.30.20G 9FREQ (RADS/SEC)2.0.6.2080.4G 93FREQ (RADS/SEC)FIGURE 1-3 TIME HISTORIES AND AUTO-SPECTRA OF 1/2 SECOND AND 1 SECOND DOUBLET INPUTS'boOI2 4 r1141040-94AMPLITUDE OF INPUT X I0 3(RADS)0-0 25 SO 7-5 100 12-5 150 17-5 2 0 0 AMPLITUDE OF INPUT X I0 3 (RAD)N n 3 00 0 25 50 7-5 100 12-5 150 17-5 200 AMPLITUDE OF INPUT X 103 (RADS)3 02-520N,0-50 0OO 25 50 7-5 IOO 12-5 150 175 20 0AMPLITUDE OF INPUT X 10 3 (RADS)FIGURE 1-4 THE EFFECT OF INPUT SIZE AND CLOCK PERIOD ONESTIMATED YAWING-MOMENT DERIVATIVES. (ESTIMATES NORMALISED WITH RESPECT TO LINEAR MODEL VALUES: HELISTAB SIMULATED DATA; PUMA 80 KNOTS,LATERAL CYCLIC 3211 INPUT62i*1s (%fwd;measured /— reconstructedreconstructedmeasured-o.oa0.12rreconstructedmeasuredmeasured /0-0* ■/ reconstructedmeasuredreconstructedFIGURE 2-1 MEASURED TIME HISTORIES TOGETHER WITH LOW-FREQUENCY RECONSTRUCTIONS:?o 'a2At 4-Afc $6b 6At 7At SAt <=\AtAtT = \OAfcFIGURE 2-2 THE FOURIER TRANSFORM AS AN APPROXIMATION TO THE FOURIER INTEGRAL USING MIDPOINT INTEGRATIONX 1 0 ~0 . 9 0 .10 . 7 5 J. i. 11oCDO0 . 4 5 i0 . 3 0 . = l D0 . 1 51o . o o i- . 1 5 i ox- . 3 0 i0 . 0 0* L o n g . c y c l i c d o u blet o C o l l e c t i v e 3211 o L o n g . c y c l i c 3211 n;L a t . eye I ic 321 1 xTa.il r o t o r 321 1i i i i i i i i i0 . 3 0 0 . 6 0 0 . 9 0INTERVAL (S)i r1 .20X10“ 1T.50FIGURE 2-3 ESTIMATED TIME DELAY VERSUS SAMPLING INTERVALUSING 9 DOF SIMULATED DATA AND A 6 OOF ESTIMATION MODEL (HELISTAB)^0 310'ft/stfs)rad0-0 t(s)0ft/s00-0 3-60*0 6rad/stCs)0-0FIGURE 2-4 a) TIME-DOMAIN FITS. 9 DOF DATA AND G DOF MODEL WITHOUT TIME DELAY. HELISTAB, 80 KNOTS, LONGITUDINAL-CYCLIC DOUBLET9 DOF R E S P O N S E :-- --6 DOF PREDICTION WITHOUT DELAY:ft/s0-0rad7-U0-0ft/s0-0radfe0-0FIGURE 2-4 b) TIME-DOMAIN FITS. 9 DOF DATA AND G DOF MODEL WITH TIME DELAY. HELISTAB, 80 KNOTS, LONGITUDINAL-CYCLIC DOUBLET9 DOF RESP O N S E :S DOF PREDICTION WITH DELAY:FIGURE 2-5 a) FREQUENCY-DOMAIN FITS (OUTPUT-ERROR METHOD) FIGURE 2-5 b) FREQUENCY-DONAIN FITS (OUTPUT-ERROR METHOD) CASE 1: SAMPLING INTERVAL 0.015625 SECONDS; CASE 2: SAMPLING INTERVAL 0.015625 SECONDS;NO TIME DELAY IN MODEL. TIME DELAY IN MODEL.PUMA, 100 KNOTS, RUN R0201A ' ; PUMA, 1QQ KNOTS, RUN RQ2Q1A00 JJ.00_____0.?0____0*. 13_____0,PC Or90<>J.00 0.10PJT,S■>RUTO.SPECP-i;l=1dl ©>£FIGURE 2-5 c) FREQUENCY-DONAIN FITS (OUTPUT-ERROR METHOD) CASE 3: SAMPLING INTERVAL 0.09375 SECONDS; NO TIME DELAY IN MODEL.PUMA, 100 KNOTS, RUN R0201AP4"i2p.|(ros*nC7>cr70*o Z CO — io >CO z> —1 m mizz. —* am ■■ oo =so o COm >r~ zz zz “OCD —< *<— 1 mCO z 70z CD*n✓a ZZcr o Z oz: a — i >*m m — i70 r* 70o ohO >* zO n*> o oczo — 1"Ocn c:cn — iro icn m70CO 70m oo 70oz zza mCO -H31rse> &r \cc .—I C~) —I >—< > 1—1 =t cn rs m m mo 3* cn >-<>>•—i7^ 2 2z: >— “o o z r <- O CT) i—icu nto m i—ii—icn to m o o 70oJB.m 3?.08*O -rooPO .—I< o > z r*~o o crC5 — 1 CD “D CO d■s^ -HCO PUm 70 o oO ?D< /—Nr*. O'rc>STATE1.09 -0.03 -0.01 0.01 0.03^ da zs>• 3: .— 1 •—» m co 70> o r~ su■o m co co mcz o o -nCD/-srr.ecoFIGURE 3-1 LOG-MAGNITUDE LOG-FREQUENCY PLOTS OF VARIABLES USED IN THE PITCHING-T10HENT EQUATION si.08 -2.00 O.CQ 2.x__ U.00<>&1*10 xJO'S'J'0 - 0 7 50-2>7fcHz.0-00-075 \ * 5\Hz.0 - 0 7 5FIGURE 3-2 THE DEPENDENCE OF PARAMETER ESTIMATES ON THE FREQUENCY RANGE AND THE NUMBER OF ORTHOGONAL COMPONENTS USED. RUN R0201A~?\oParV\a\-fW w0-0750-5760-0Hz0-0751- 510-0Hz.f aHz.FIGURE 3-3 THE DEPENDENCE OF PARTIAL-F RATIOS ON THE FREQUENCY RANGE ANO THE NUMBER OF ORTHOGONAL CONPONENTS USED.RUN R0201A0 075 - 0*%7€>pEFIGURE 3-4 THE DEPENDENCE OF THE SQUARED-CORRELATION COEFFICIENT ON THE FREQUENCY RANGE AND THE NUMBER OF ORTHOGONAL COMPONENTS USEO FOR THE PITCHING-MOMENT EQUATION.RUN R0201AFREQUENCY(Hz)NUMBER OF ORTHOGONAL COMPONENTSFIGURE 3-5 THE DEPENDENCE OF THE SQUARED-CORRELATION COEFFICIENT ON THE FREQUENCY RANGE AND THE NUMBER OF ORTHOGONAL COMPONENTS USED FOR THE PITCHING-flOMENT EQUATION (FOR FULL RANGE OF FREQUENCIES. RUN R0201A)a s s ®\0.016cr error boond0 000 13X) F(S)5 0 7.0 11.09.00.000-0.008, 13.0 t(s)5.0 7.0 11.09.01.20.45.0 13.0 t (s)7.0 9.0 11.00.0050.015 13.0 t(s)9.0 11.07.05.00.090.397.05.0 9.0 13.0 t(s)11.00.020\ __0-036 J13.0 t(s)5.0 7.0 11.09.0FIGURE 3-7 VARIATION OF ESTIMATES AND ERROR BOUNDS (BROKEN LINES) WITH LENGTH OF RECORD USED. RUN R0201Abi -0-oFIGURE 3-11 b) FREQUENCY-DOMAIN EQUATION-ERROR FITSFOR DIFFERENT NOS. OF ORTHOGONAL COMPONENTS (SINGULAR-VALUE DECOMPOSITION).MODEL WITH TIME DELAY. NORMAL-FORCE EQUATION. RUN RQ2Q1A_______-o—o1^0.00•"“0FIGURE 3-11 b) FREQUENCY-DOMAIN EQUATION-ERROR FITSFOR DIFFERENT NOS. OF ORTHOGONAL COMPONENTS (SINGULAR-VALUE DECOMPOSITION).MODEL WITH TIME DELAY.l?0.00 180.00—o120.001^ 0.00 160.00FIGURE 3-12 a) FREQUENCY-DOMAIN EQUATION-ERROR FITSFOR DIFFERENT NOS. OF ORTHOGONAL COMPONENTS (SINGULAR-VALUE DECOMPOSITION).MODEL WITHOUT TIME DELAY.YAWING-MOMENT EQUATION. RUN R0201AFIGURE 3-12 b) FREQUENCY-DOMAIN EQUATION-ERROR FITSFOR DIFFERENT NOS. OF ORTHOGONAL COMPONENTS (SINGULAR-VALUE DECOMPOSITION).MODEL WITH TIME DELAY.YAWING-MOMENT EQUATION. RUN R0201A-o o"OFIGURE 3-13 a) FREQUENCY-DOMAIN EQUATION-ERROR FITSFOR DIFFERENT NQS. OF ORTHOGONAL COMPONENTS (SINGULAR-VALUE DECOMPOSITION).MODEL WITHOUT TIME DELAY. ROLLINS-MOMENT EQUATION. RUN R02Q1A_____FIGURE 3-13 b) FREQUENCY-DOMAIN EQUATION-ERROR FITSFOR DIFFERENT NOS. OF ORTHOGONAL COMPONENTS (SINGULAR-VALUE DECOMPOSITION).MODEL WITH TIME DELAY.ROLLING-MOHENT EQUATION. RUN R0201A-oFIGURE 3-14 MULTIRUN TIME-DOMAIN RECORDS TOGETHER WITH LOW-FREQUENCY RECONSTRUCTIONSoNJro rooo310co(bClCO •+•v-/0) d>orofo.CO "h £0 -*•0 00.0 0.80.2 0.4 0.60.8 f (Hz)0-00.80.0 0.4 0.60.2f (Hz)0.0 0.80.60.40.20.0f (Hz)1800.80.4 0.60.0 0.2f (Hz)MEASURED -----PREDICTED xx xxFIGURE 3-15 FREQUENCY-DOMAIN EQUATION-ERROR FITS FOR MULTIRUN DATA SET3 ^1.20.0_0.05 0.80.60.2 0.4f (Hz)3200.05 0.2 0.6 0.80.4f (Hz)MEASURED -------PREDICTED * *FIGURE 3-16 q ) FREQUENCY-DOMAIN VERIFICATION.MEASURED VALUES OBTAINED FROM RESPONSE TO LONGITUDINAL-CYCLIC 3211 INPUT. PREDICTED VALUES BASED ON MODEL IDENTIFIED FROM RESPONSE TO LONGITUDINAL-CYCLIC DOUBLET0.80.0 _ 0.05 0.2 0.60.4 0.8f (Hz)160,05 0.6 0.80.2 0.4f (Hz)MEASURED -------PREDICTED * *FIGURE 3-16 b) FREQUENCY-DOMAIN VERIFICATION.MEASURED VALUES OBTAINED FROM RESPONSE TO LONGITUDINAL-CYCLIC DOUBLET. PREDICTED VALUES BASED ON MODEL IDENTIFIED FROM MULTIRUN DATAFIGURE 4-! PARAMETER ESTIMATES AND COST-FUNCTION VALUE AS A FUNCTION OF ITERATION NUMBER. RUN R0201ArHMHnciu; -0.08 -£PRRRHETER.ESTHS 0.20 o.M Q.PARAMETER.ESTHSO.SO 1,00 I.PRRRHETER.ESTHS -o .e o -o .s o -oPRRRMETER.ESTHS *10'-O.OJ 0.00 0.01FIGURE 4-1 CONTDPRRRMETER.ESTHS •132.00 ~,I20.0Q —.12*1.00 -,120.00 -118.00PRRRHETER.ESTHS20.00 HO.00 60lpo.oo 2UD.00 SpP.OO Jt. 00 2p.0Q *UO.OQ 60.00 M.OO J.00 Up. 00 BO.00 IpD.OO IpO.OOsM+s- \s- I2.Z.FIGURE 4-2 ESTIMATES OF THE DIAGONAL ELEMENTS OF INVERSE OF ERROR-COVARIANCE MATRIX AS A FUNCTION OF ITERATION NUMBER. RUN R020IASECONDSV OTUoSECONDSw c o'o .o oSECONOSwet)jr \ wet)FIGURE 4-3 MEASURED SIGNALS FOR TRANSFER-FUNCTION EXAMPLE. RUN R0201A0.520.H8 0.56 0.600.280.12 0.16). 82 0.86 CFREQUENCT.(HZ)0.0.52 0.580.32 0.88 CFREQUENCT. (HZ)0.280.08 0.12FIGURE 4-4 FREQUENCY-DOMAIN OUTPUT-ERROR FITS.TRANSFER FUNCTION EXAMPLES. RUN R0201A5EC0N0SSECONOS'o .o o o '. flo iTe o 2'. 110 a'. 20 u'. oo y eo s .e o e '.uo >'.20 e ' . c o a ’. eo a ' . e o i b . n o i ' i . 2 o 1 2 .0 0 ib '.a o iVT 6 o u .n oFIGURE 4-5 TIME HISTORIES OF THE LONGITUDINAL-CYLIC CONTROL AND LATERAL PSEUDO CONTROL TERMS.COMBINED DATA RUN: RQ201A, R050IAFIGURE 4-6 FREQUENCY-DOMAIN OUTPUT-ERROR FITS. • FIGURE 4-7 TIME-DOMAIN FITS.COMBINED DATA RUN: R0201A, R050IA COMBINED DATA RUN: R0201A, R0501A60 eP,. 00___0.00 1. CC____2.10IG>FIGURE 4-8 a) MAGNITUDE OF FFT's FOR DOUBLET INPUTFIGURE 4-8 b) MAGNITUDE OF FFT's FOR THE EFFECTIVEINPUT. COMBINED DATA RUN p (t)tyfrSECONDS12.00 12.80 13.60 l^.MOFIGURE 4-9 H U E HISTORIES OF CONTROLS AND LATERAL PSEUDO CONTROLS. I1ULTIRUN EXAMPLE: R0201A, R1501UFIGURE 4-10 FREQUENCY-D0MA1N OUTPUT-ERROR FITS (RANK 13). FIGURE 4-11 TIME-DOMAIN FITS (FOR RANK 13 ESTIMATES). MULTIRUN EXAMPLE: R0201A, R1501U ; MULTIRUN EXAMPLE: R0201A, R1501U<>£9?. 00FIGURE 4-12 PEDAL DOUBLET INPUT. PUNA, 100 KNOTS, RUN R1201L.-*-Delay. -Npedal o No delay. -NpedalXI0-23.83:2.50:RANK-*-Delay. Nv o No delay. NvX10-20.90.0.80.0.700.60 XHELISTAB10 912RANK $ERRGR 80UN0-*-Delay. -Lv o No delay. -Lv2.40XHELISTAB1 .87,1 .330.80,RANK-^-Estimate oF DelayXI 0_11 .83:RANK ERROR-*-Delay. -Np o No delay. -Np0.40,0.270 .00, -I XHELISTAB91012RANK $ERR0R "U™-*-Delay. -Lp o No delay. -Lp2.10,XHELISTAB.73:.00error boundRANK-*-Delay. Cost o No delay. Cost-150,-160:-170.-180RANK-*-Delay. -Nr o No delay. -NrX10“16.20,5.60:4.40;ERROR BOUNDRANK-*-Delay. Lr o No delay. Lr1 .70,0.63:ERROR BOUNDRANKFIGURE 4-13 ESTIMATE VERSUS RANK OF INFORMATION MATRIX USED. PUMA, 100 KNOTS, RUN R1201L5FREQUENCY. (HzfRo'tro.ou 0.08 0.12 0.16 0.20sFIGURE 4 - H OUTPUT-ERROR FREQUENCY-DOMAIN FITS. NO DELAY IN ESTIMATION MODEL.PUMA, 100 KNOTS, RUN R1201LFIGURE 4-15 OUTPUT-ERROR FREQUENCY-DOMAIN FITS. FIGURE 4-1G OUTPUT-ERROR FREQUENCY-DOMAIN FITS.DELAY IN ESTIMATION MODEL. DELAY IN ESTIMATION MODEL, RANK-9 SOLUTION.P UMA, 100 KNOTS, RUN R1201L PUMA, 100 KNOTS, RUN R1201L£—T>~o>FIGURE 4-17 TIME—DOMAIN VERIFICATION. FIGURE 4-18 TIME-DOMAIN-VERIFICATION.NO DELAY IN ESTIMATION MODEL. ' DELAY IN ESTIMATION MODEL. PUMA, 100 KNOTSj RUN R1201L PUr1A’ 100 KN0TSs RUN R1201L“D/ NIT.*-■' w r iSECONDSgSECONDSsP(t)2b.oo sU .0 0ill.00SECONDSgill.00FIGURE 4-19 TIME-DONAIN VERIFICATION.DELAY IN ESTIMATION MODEL,RANK-9 SOLUTION. PUMA, 100 KNOTS, RUN R1201L1 .00_0.80:0 = 60:0 = 40:0 = 20:o = oo:RANKFIGURE 4-20 TYPICAL REDUCTION IN COST VALUE FOR USE OF INCREASING RANK OF INFORMATION MATRIX. NORMALISED WITH RESPECT TO FULL-RANK SOLUTION,FIGURE 4-21 LATERAL-CYCLIC DOUBLET INPUT. PUMA, GO KNOTS, RUN R2501R100 KN Nv -b -60 KN NvX10-20 .90,,0.80; □ helibtas0.70;0 .GO;12 10 91 1-*-100 KN -Np -b -60 KN -Np-*-100 KN -Nr -b -60 KN -NrXI 0-1XHELISTAB0.40,0.27;1012 9□HELISTAB ! XHELISTAB6.205.60;912 10RANK RANK RANK-B-100 KN -Lv 60 KN -LvX10-2.87; □ HELISTAB1 .33.12 10 9100 KN -Lp -b -60 KN -LpRANK21 .731 .00.10RANK9XHELISTAB □ HELISTAB100 KN Lr b -60 KN Lr70S B - -0.63.ERRORRANKFIGURE 4-22 a) ESTIMATE VERSUS RANK OF INFORMATION MATRIX USED. PUMA, 60 KNOTS (LATERAL-CYCLIC RUN R2501R) & 100 KNOTS (PEDAL RUN R1201L) COMPARISONSX10-260 KN -L Lateral cyclic5 .00,.3.53ERRORRANK—*— 60 KN Estimate oF Delay-*-60 KN CostXI O' X102.00., -9.50,.67233; -10.63J.00;12 1 1 910RANKERROR BOUNDFIGURE 4-22 b) ESTIMATE VERSUS RANK OF INFORMATION MATRIX USED. PUMA, 60 KNOTS, RUN R2501RUS I —H >- >• ~o“ -< d— r~ oT)>&~ >■ 1 >-H?C —h “T1hO Ocn o o m-0.09SECONDSFIGURE 4-25 PEDAL DOUBLET INPUT.PUMA, 60 KNOTS, RUN R340IR-*-60 KN -Npedal o Lp - Lr relationXI O'2.50,.83:1 .50:ERROR 80UN0RANK-*-Estimate oF Delay o Lp - Lr relation-*-60 KN Cost o Lp - Lr relationXI O'3.00,.33:.67:10 912RANKX10-7 .00,-7.33-7.67:-8 .00:RANKBOUND-*-60 KN Nv o Lp - Lr relationX10-2 1 .30.1 .070.830.6012 1 1 10RANK-*-60 KN -Np o Lp - Lr relation[ERROR BOUND0.50,0.23,XHELISTAB-.03,-.30,ERRORRANK60 KN -Nr o Lp - Lr relation0.62,0.24:0.05:.91012 1 1RANK [ERROR BOUND60 KN -Lv o Lp - Lr relation.87: XHELISTAB1 .33:0.80:RANK60 KN -Lp o Lp - Lr relation21 .73,1 .37,1 .00,RANK60 KN Lr o Lp - Lr relationXHELISTABERROR BOUND0.80,0.570.33,ERROR BOUNORANKFIGURE 4-26 ESTIMATE VERSUS RANK OF INFORMATION MATRIX USED. PUMA, SO KNOTS, RUN R3401RjrFLIGHT DATA^ Time - d o m a i n veriPication oF modelOne equation at a time estimatedEQUATION-ERROR ESTIMATION (FREQUENCY DOMAIN):Estimates oF measurement zero-oFFsets & initial conditionsOUTPUT-ERROR ESTIMATION (TIME DOMAIN):State-space model estimatedOUTPUT-ERROR ESTIMATION (FREQUENCY DOMAIN):FIGURE 5-! BASIC STAGES OF THE IDENTIFICATION METHODOLOGY- 347-KQMElitGLAlIiifoA State matrix.A,b Matrices used to represent a linear system of equations.A ff, A t r ... Bp Partition matrices of the 9 DOF model relating fuselageand rotor effects.A n >Au> ••• B 22 Partitions of the six-degrees-of-freedom state andcontrol matrices,Az,Gh,U 2 State matrix, control dispersion matrix, and controlvector for steady-state Kalman-filter representation of the system.True values and estimates of the state-space matrix (equation error),Q^ n Coefficients and aliased coefficients of the Fourier-series expansion.CL/ Components of acceleration._A.Cco”) Frequency-dependent matrix associated with the statematrix A , and used in the frequency-domain model y representation developed for quantities of the form C 3 . Line-search scalar parameter,...^i .J3p Estimates of coefficients associated with various stages in the subset-regression procedure.,(3 Control dispersion matrix,0 Matrix with orthogonal columns used in singular-valuedecomposition.B w (Approximate) Bandwidth of measurement noise.Frequency-dependent matrix relating to delays in the controls.b Vector of biases in measurements.Columns of matrix B used in singular-value decomposition.bij, True values and estimates of elements of the controldispersion matrix, (equation error),by,b/3 br Constant biases in indicated measured quantities."fr Frequency-independent matrix associated with the controldispersion matrix B * and used with vector quantities of the form Cv3*B.6w ) Frequency-dependent matrix associated with the matrix^Bjto) * and used with vector quantities of the formC 3./3o , /3iCj Pa, fit. Coning, longitudinal and lateral-cyclic flapping, and differential coning multiblade flapping angles./sC.CcS)j DCco) Frequency-domain equation-error control terms used inv the approximate estimation of time delays.aCcS) Expected model-output vector (output-error method) inthe form £ v3 •C o n 4 C 3 Condition number of indicated matrix.Cov C ) Covariance matrix of indicated vector quantity,- 348 -^ Input vector containing higher-order terms of Taylor~ expansion.Diaq ( ") Diagonal matrix with indicated quantities in the leadingdiagonal.e C O , e C r ) Estimated equation-error residuals (time and frequencydomain).£ Vector of equation errors.a a f\ a Collective, longitudinal-cyclic, lateral-cyclic, and c) ,c> p pedal control positions.E C 3 Expected value of indicated quantity,^Xp E 1 Complex exponential.P'TOT A W, Veco) Frequency-response values (estimated and predicted)Partial-F statistic (or ratio)Total-F statistic (or ratio)?(t)J General function of time; discrete function of time andfiltered form.AFcco) Frequency-response errors.Ft! Fourier transform operator (integral defined betweeninfinite limits).&CcoJ} G^COn) Correction terms for a non-periodic measurement window ~K ~ (frequency-dependent).9 Gradient vector of the output-error cost function at the~ iteration.9 Acceleration due to gravity.*C Vector of true output-error model parameter values.Update increment vector for model parameter values (output error).Update increment vector required to take estimates to ^ K the true model parameter values (output error),Vector of output-error model parameter estimates and the vector of update increments at the iteration.V Set of output-error model parameters.H Measurement transition matrix (relating measuredquantities to states in the model).H aCci) Frequency-dependent matrix relating to delays in themeasurements.H Frequency-independent matrix associated with themeasurement transition matrix H • and used with vector quantities of the form C l .H xCco) Frequency-dependent matrix associated with the matrixHjto)» and used with vector quantities of the form C 1 .Identity matrix.- 349 -T(e) Equation-error cost function (time and frequencydomain),TCY.S) Frequency-domain output-error (maximum-likelihaod) costfunction.T„c©) Equation-error cost function (Hartley transform).T t Time-domain output-error (maximum-likelihood) costfunction,\ Complex number such that J* = -1.A 7 Change in cost-function value.K Matrix of Kalman gains.K Vector of trim constants for measurements.w \s Vectors of control and measurement zero offsets. ^h r, Coefficients of the equation-error control terms C (w)l\c J K* , £ , Nand (cu) .K c Set of control zero offsets.K m Set of measurement zero offsets.^ Eigenvalues of indicated matrix.loae Natural logarithm.i;jL" Offsets (in body-fixed axes) of indicated measurementdevices relative to the centre of gravity,JVf Output-error information matrix at the Kth iteration,Mho . Output-error information matrix evaluated at the trueparameter values.I V •• Pitching-moment derivatives with respect to subscripted quantities./^, Mean and variance of measurement noise.J\| Number of sampled points in time-domain record.|\J Number of data points used in equation-erroridentification.IUi .Dk Numerator and denominator transfer-functioncoefficients,fyjgj Number of frequency-domain points used in output-erroridentification.f\Jx Number of states in the model for output-erroridentification,fW j 0/5 , Hr Measurement noise terms in indicated quantities.Vcnnt') Vca>*) Vec'tors discrete Kalman-filter innovations; time * ^ domain (real-valued) and frequency domain (complexvalued) .CD 3_con Angular frequency and discrete values of angularfrequency. cof Fundamental frequency.cj5p Short-period mode natural frequency.-350-«J\. Diagonal matrix containing eigenvalues of the equation-error information matrix.O o Quantities of the indicated order.AP Number of unknown parameters to be estimated.Y, Rotor azimuth position.0 Vectors of true orthogonal parameter values and^ estimated orthogonal parameter values used in singular-value decomposition. AVector whose first r elements equal those of 0 , theremaining elements being zero, ^The it'h element of $ .Discrete versions of the state matrix A and control dispersion matrix 0 .Orthogonal-update increment vector related to 0 ,Vector whose first r elements equal those of A ® , theremaining elements being zero. ^Conditional probability (a conditioned on b),P Co./ bJAngular rates of pitch, roll and yaw. f\ C*) Laplace transforms of indicated quantities.Time-domain error-covariance matrix,Squared-correlation coefficient.RaV Theoretical ratio of parameter estimates.Partial-correlation coefficient used in subset regression procedure.Real and imaginary parts of the complex-valued quantity [ ].Sj ‘DO-jO Frequency-domain error-covariance matrix (frequency- independent and frequency-dependent versions).S Diagonal square matrix with singular values in theleading diagonal.5 Matrix obtained from the inverse of 5 used in singular-value decomposition.Variance of estimated coefficients obtained for the subset-regression procedure.SXyCco)j Sx*Cw) Cross and auto-spectral density functions, c! Frequency-domain error-covariance matrix used withvector quantities of the form[[v3 •S Laplace operator.S\. Singular values.£ Set of unknown parameters of error covariance matrixcr 3 cT* Standard deviation and variance of parameter estimate.(1).0"% Equation-error variance of residuals (true and estimatedvalues). (2).- 351 -T Time length of measurement window."f© Transfer-function parameter.t Time.At Sampling interval.§*/> ©FF Vector of parameter estimates relating to the~ instrumental variable technique and its frequency-domainequivalent.<9 Vector of true parameter values (equation error).0 Vector of parameter estimates (equation error),Euler pitch, roll, and yaw angles (and perturbations from trim).T Time delay; time constant."C©,Taz Effective transfer-function time delays.© Set of equation-error parameters.Xc Set of time delays in the controls.T m Set of time delays in the measurements.U Matrix formed by closing together the eigenvectors ofthe system matrix. (1).u Matrix with orthogonal columns; used in singular-valuedecomposition. (2),U Control vector.0 11 _ Portions of the control vector relating to theX uow&, u«vr °longitudinal and lateral controls respectively.\JL Columns of the matrix U used in singular-valuedecomposition.Ul Eigenvectors of system matrix.\J,V,W Aircraft translational velocity components (totalvalues).Uef Effective control input.w Aircraft translational velocity components(perturbations from trim).V Orthogonal matrix used in singular-value decomposition.Vet"), VCCij VCcS), Continuous and discrete measurement noise (time domain* 'Cfclf’ ~ and frequency domain).V5p, Speed, flank angle and incidence angle.Vi. The itM row of VW e t ), W C o K) Process noise vector (time and frequency domain).X , ft Matrix of independent-variable frequency (or time)response data arranged in columns (for equation-error method). The form of the matrix with a hat is used to indicate frequency-domain quantities when a distinction with time-domain quantities is required.X State vector.X(b), UCtlj ZCt) State, control and output vectors (time domain).XCcS), GccJ);^CcS) State, control and output vectors (frequency domain),vJ" ' CT v t/%CcS) dccoy.^ CcS) State, control and output vectors in the form X*vlA. * «. t- 352 -Xco) Vector of initial state conditions.X u >no X LnT Partitions of the state vector relating to longitudinaland lateral states respectively.Rotor state vector and fuselage state vector.A X , A X Nft<> Vector associated with correction term for a non-~ periodic measurement window, and version for NAGdefinition of discrete Fourier transform.X,'/,Z A Body-fixed axes components.XtCr), UjCO, Xc(0 Independent variables relating to states and controlsGj(r) ^ (equation error); time and frequency domain.X,.Cr) j X«.C.r) Dependent variable relating to time derivative of state(equation error); time and frequency domain.XjCr) Subset-regression independent variables.X)*<0 Subset-regression variables associated with X, (D .Xu, Xw •• Forward force derivatives.X0 Set of initial state conditions.) Short-period mode dampingV Equation-error dependent variable (vector form).ycr) , 9CO Equation-error dependent variable and model prediction(scalar form).y(r) Subset-regression dependent variable constructed fromresiduals.y,,,'/* •« Lateral force derivatives.Z fr Complex-valued matrix used to illustrate the frequency-domain equivalence to the instrumental-variable technique.Z fp. Real-valued matrix associated with 2V* .~ZXv Instrumental-variable matrix.Z Model output vector.2ccS) Model output vector (output-error method) in the form* 1: v i'~ZCcou/t-0 Expected value of Z (“0 conditioned on the sets andY ♦Vectors associated with the correction term for a nonperiodic measurement window (relate to A X and A X WfV&.).ZujZvv*. Normal force derivatives.2L Set of measured frequency-domain values corresponding tothe first q frequency values.Z 3 Derivative of I 3 with respect to time.3 ^ ^ C 3 Partial derivative of C 3 with respect to indicatedquantity.H 3 t Least-squares pseudo inverse of C 3.£ 1 Inverse of [ 3.£ 3 T Transpose of C 3.E 3* Transpose of complex conjugate of C 3.\ 1 Magnitude of indicated quantity.) | Determinant of indicated quantity.£*3 Fourier transform of the time derivative of C 3.Measured, value of C 3.Offset term for C 3 combining trim constants and measurement biases.Fourier transform of [ 3 - unless otherwise defined separately.Expected (or mean) value of C 3.Matrix and vector norm of indicated quantity.Definition.Integer value of indicated quantity (i.e. truncated not rounded).Frequency-domain vector quantity relating to [ 3 and obtained by stacking the real on top of the imaginary part of the otherwise complex-valued vector quantity. Subscript e signifies trim values of indicated scalar quantities, e.g. 0®.-354-Appendix 1 Linearised Equations of Aircraft MotionThe force, moment and kinematic equations of motion are given below:Translational MotionmU = - m C W ^ - V r ) + Xcx3 ~ mcj S\n© (ai.dmV = - m C U r - W p ) + Yc x) + rqg Cos© Svn$zS (ai.2)mW = ~m(Vp - Ucp + Z(x) + rngCose Cos0 (A1..3)Rotational MotionIxxp = ( Iyy - I 2ZV + Ixz (r + + Lex) (Al.4)Iyy^ ' ( Izz - IxxV p + Ixz ( r*" p1) + MCx) (A1.5)Izzf" = Clxx -Iy/)p<£ + Ixz(p-^r) + IMoO ■ \n 0s) + /ex)/m(A1.14)W = - ( Ve p - Ue. - q (© Cos 0 1 ■sin ©<. + Sin 0 C cos ©e) + Zcx) / m (A1.15)Xxx p = Ixz r + Lw) (A1.16)lyy^ Mcx) (A1.17)Iz i r = 1 x2 p + l\)cx) (A1.18)II p + a sin 0s. t"an Qs. + r c o s ^ c. \'an <9e. (A1.19)© * a cos 0s. - r sin (A1.20)Cos 0c Y - <^ Sin 0^ + r cos 0^ . (A1.21)The linearised equations of motion can be written instate -space form as:X - AX + Bu/v* r* (A1.22)whereTX = Cu» v5 Pj 0,r)^ 1 A»c >-357-The system matrix A and control dispersion matrix B follow from the linearised equations of motion given by Equations (A1.13) to (A1.21,) and are as follows:A =X „ X w X<£ We. Xe+C, X v X p x ^ X r + V eZ o z . Z ^ + U e Z * * G z . Z F- V e Z ^ * C a Z rM „ M w M e M w M P Mr0 0 Cos 0e 0 0 O O -" S m ^ cy „ Yw Ye +C X , Yp + W e Yr-UcC C L i l ;U Lp L^s L r0 0 C-5 0 0 1-0 0 c 6. n ; l \ CW l i u lro; [ \ i ; i d ;c , = - 0 Cos ©€ c * -• q COS 0e ©eC;5 ~ “ 2 *010 CoS ©e ©4. = - 0 ‘bvn ©e. ‘binC5 - Sin (^e ^o.n ©c. C g = Cos ta n ©cB -X 0 X'l.i V X/\Pz,„ Z>„ Z/),,. Z)el;. C . C U/)p0 0 0 0w;. V * w ; f-358-where, for a general subscript c , Lc = l_c + • Ixz/Izz/ .and I\jc = [\jc + fc . Ixz/Ixx. The derivatives relating to the force equations (i.e. for a general subscript c , Xc , V Q , and Zc) are divided by the aircraft mass m. The pitching-moment derivatives are divided by Iyy, and therolling-moment and yawing-moment derivatives are dividedz aby Izz/(lzzlxx - Ixz ) and Ixx/(Izzlxx - Ixz )respectively.-359-Aippeodix 2 Conditions for Unbiased Estimates from theFreqnency-Domain Output—Error MethodThe proof given below to show, in the absence of process (or model) noise, and for uncorrelated Graussian measurement noise, that unbiased parameter estimates are obtained using the output-error method, is in essence the same as that given in KLEI001; however, the proof is presented here for frequency-domain quantities and some additional explanation and development is provided. The same notation is used as in Chapter 4.Consider Equation (4.2-44) giving the linearised approximatioon for elements of the Hessian matrix at the Kth iteration (defined by Equation (4.2-39). We can write:M = Z A*Cco) p> K c . c S ) (A2.1)cowhere A«.,(co) = (A2.2)And similarly for the gradient vector 3 (defined by Equation (4.2-40)) we can write, using Equation (4.2-42) and the definition given by (A2.2):9 - 2 A CcS) - ZCcS l \ (A2.3)~ CO ^ ~For the Gauss-Newton method, given in Equation (4.2-38), and using (A2.1) and (A2.3), we have:-360-rei: - - a c ; » - w y -( Z A* G») 5) K ccS)) ( H P \ C g 5 ) $ L 2CcS) - ^ lCco)D (A2 . 4)GO W ~The value of the output vector *X(co) at the Kth iteration is dependent upon the vector of parameter estimates K. obtained, and indicating this fact in the notation, we can expand the true model output vector in a Taylor expansion about the output vector obtained at the kth iteration:#c. is the parameter increment required to move from the current parameter estimates to their true values. The measurements and the true model output vector are related (using Equation (4.2-25) and the definition (4.2-30)) by:ZCci) = A! (co/() + VCcn) (A2 . 6)Hence, by using the first-order terms only in (A2.5), we can substitute into (A2.6) to obtain:v1 ^ZCci) - - A c d ) A X t + V CcJ") (A2.7)Now using the first-order approximation for the left-hand-361-side of (A2.7), we can substitute it into (A2.4), and obtain:A 0 for a minimum.3. Re-order (X,, CXz , and CX3 so that J()4. Fit a quadratic A 0 for a minimum.5. A series of tests are now carried out to see if oT is to be accepted as the value of the line-search parameter. These are:-365-a) I c C - (X,I 4 TESTS (3)b) Absolute value of the gradient of the line from(0(i / J(

t

Abstract

Extracted data

t

Abstract

Extracted data

Related items

Related items