Flutter Test Data Processing Based on Improved Hilbert-Huang Transform

Flutter tests are conducted primarily for the purpose of modal parameter estimation and flutter boundary prediction, the accuracy of which is severely affected by the acquired data quality, structural modal density, and nonstationary conditions. An improved Hilbert-Huang Transform (HHT) algorithm is presented in this paper which mitigates the typical mode mixing effect via modulation. The algorithm is validated by theory, by numerical simulation, and per actual flight flutter test data. The results show that the proposed method could extract the fluttermodel parameters and predict the flutter speedmore accurately,which is feasible for the current flutter test data processing.


Introduction
"Flutter" is the self-excited vibration of an elastic structure under the coupling of aerodynamic force, elastic force, and inertial force; it is often accompanied by catastrophic structural damage [1].Flutter analysis is a crucial aspect in the design of new or modified aircraft.Wind tunnel tests, aeroelastic models, and flight flutter tests are important components of flutter analysis.Such tests are risky, timeconsumptive, and costly to conduct.To this effect, accurate and effective flutter test data processing techniques are in high demand.One of the most popular flutter boundary prediction (FBP) methods currently is the damping-based method [2,3], which serves to extrapolate a curve fitted by the estimated structural modal damping factor against the airspeed to the abscissa while zero damping is considered the critical flutter criterion.The key to this type of FBP method is the accurate identification of modal damping factors.
Many algorithms have been previously developed for identifying flutter test modal parameters, including fast Fourier transform-(FFT-) based methods [4,5], Random Decrement Technique (RDT) [6], natural excitation technique combined with the eigensystem realization algorithm (NExT-ERA) [7], time series analysis based on the Autoregressive (AR) model [8,9], and Stochastic Subspace Identification (SSI) [10,11].Although some of these algorithms are effective, none is ideal; for instance, nonstationary measured data and low signal to noise ratio (SNR) affect the Fourier-based methods.The AR model requires an appropriate model order and algorithm to function properly; RDT functionality is limited by the number of the main modals contained in the structure [12].SSI and NExT-ERA are problematic in terms of their spurious mode [13,14].The Hilbert-Huang Transform (HHT) proposed in 1998 is an adaptive scheme well-suited to nonlinear, nonstationary time series analysis; however, the mode mixing effect [15] which emerges when dealing with signals over multifrequencies in each frequency band severely limits its application to flutter test data processing due to the inherent density of modal problems.This paper proposes an improved HHT which applies to flutter test modal parameter identification.

Hilbert-Huang Transform Theory
HHT is a relatively new nonstationary and nonlinear signal processing approach which is not limited under linear and stationary spectral analysis theory based on the Fourier transform.HHT has two basic steps.First, the original signal is decomposed into a series of Intrinsic Mode Function (IMF) components via Empirical Mode Decomposition (EMD); second, a Hilbert transform is performed on each IMF 2 Mathematical Problems in Engineering component to obtain instantaneous amplitude and phase information.
According to classical Fourier theory, the local frequency only can be defined when there is at least one complete sine or cosine oscillation.This strict definition is not suitable for nonstationary signals which have ever-changing frequencies.A Hilbert transformation for the original signal () is performed as follows to reveal the instantaneous frequency of the signal at any moment: where  denotes the Cauchy principal value; the transformation is available for all   classes.According to this definition, when a complex conjugate is formed by () and (), the analytical signal () is where () = √ 2 () +  2 () denotes the instantaneous amplitude and () =  tan(()/()) denotes the instantaneous phase.
In this context, the instantaneous frequency can be defined as follows: In a sense, the local properties of () are emphasized because the Hilbert transform is defined as a convolution between () and 1/ under (1).The local characteristics are further expressed as a polar coordinate in (2), which is the optimal approximation of the trigonometric functions with variable amplitude and phase.
Although the definition of the instantaneous frequency and phase is given, some extra conditions are required to ensure the instantaneous frequencies have physical significance.The signals must be symmetrical with a local mean of zero and must have the same zero crossings and extreme points.A general signal can be described by the Hilbert transform in the frequency domain only after EMD.
The EMD decomposition is based on the assumption that any signal is comprised of a series of different IMF components.The envelope of each IMF is defined by the local maximum and the minimum is symmetrical about the abscissa.The number of extreme points and zero crossing points should be equal or no more than one in the entire data sequence.Ideally, each IMF component contains a single modal and the number of IMFs equals the modal number of the original signal.
The EMD steps to decompose any signal () are as follows.The EMD procedure is also shown in Figure 1.The IMF, per its theoretical definition, meets the conditions of the Hilbert transform as the instantaneous frequencies have physical significance.Thus, the HHT can be used for complex nonstationary signal analysis.

Improved HHT and Modal Parameter Identification
3.1.Adjacent Mode Separation.Mode mixing problems are inherent to the EMD method.When the signal contains multiple modals with similar frequencies, the energy between the IMF components is widely divergent and the instantaneous frequencies of the IMFs may intermittently fluctuate, so the IMF components will not be decomposed correctly.
As an example, consider the signal () =  1 cos(2 1 ) +  2 cos(2 2 ).If the two frequencies are adjacent, the lowfrequency component will be mixed into the high-frequency IMF component during the EMD process to create a mode mixing effect.After several simulation experiments, when 0.5 <  1 / 2 < 1, mode mixing is inevitable during the decomposition process and grows increasingly severely as the ratio approaches 1; this is a serious problem in regard to flutter test data analysis due to the dense modal and modal coupling.Without loss of generality, given  1 <  2 , for any positive number  0 less than  1 , ( 1 − 0 )/( 2 − 0 ) <  1 / 2 .The mode mixing effect can be reduced when the two close frequency components shift leftward simultaneously.Therefore, any signal () = ∑  =1   () can be processed via the following steps.(b) Modulation.Multiply the constructed complex analytical signal with the complex exponential signal () =  −2 0 ⋅ .
Equation (5) shows that the modulated signal z() is a complex signal.Its real and imaginary parts are a pair of Hilbert transform pairs.Therefore, if the original signal () is processed as above, the corresponding modulated signal X() can be obtained as follows: and ( 5) can be rewritten as By definition, the instantaneous frequency of each IMF component of the modulated signal is The instantaneous frequency is reduced by the above process.As long as the proper selection is made, the alias mode can be completely separated.

Modal Parameter Identification.
The aeroelastic system is a typical multiple degree-of-freedom (DOF) dynamic system.The impulse response function can be expressed as follows: where   denotes modal natural frequencies,   denotes modal damping ratios, and   denotes natural damping frequencies with   =   √ 1 −   2 .
The proposed "improved" HHT method was designed to overcome the mode mixing problem caused by the dense modal and modal coupling for effective modal parameter identification in flutter tests.It works in the following stepwise process.
(a) Select an appropriate complex exponential signal to modulate the flutter test signal to reduce the denseness of the adjacent modal by (5).
(b) Obtain the IMF components by using EMD: (c) Apply the Hilbert transform to each IMF component to calculate the modal parameters.The instantaneous amplitude   () and instantaneous phase   () can be expressed as Logarithms and differential operators are performed for ( 11) and ( 12), respectively: Then the natural frequency and damping ratio are identified by straight line fitting with the instantaneous frequency   () and the logarithm of instantaneous amplitude ln   ().
(d) The IMF component is the decomposition of X(), so the actual natural frequencies and damping ratios can be calculated from ( 15) and ( 16) by demodulation:

Numerical Simulation and Application
4.1.Modal Parameter Identification.The effects of flutter boundary prediction are dependent on accurate modal parameter identification.The following simulation was primarily run to test the performance of HHT in terms of modal parameter estimation.The response of a 3DOF dynamic system was constructed according to (9) with the modal parameters shown in Table 1.The temporal history of the response waveform excited by the impulse signal is shown in Figure 2, where the duration is 1 s and sampling rate is 1000 Hz.Since the first and second modal frequencies are adjacent,  1 / 2 = 0.7, the traditional HHT method inevitably results in modal aliasing as shown in Figure 3, where each IMF component and its frequency spectrum are listed.Figure 3 also shows that the frequencies in the first and third IMFs  are pure of any other frequency components, but the second contains both 70 Hz and 50 Hz frequency components.
In the modulation frequency  0 = 30Hz, all the frequencies contained are shifted by 30 Hz to the left; this reduces the ratio of  1 / 2 to 0.5 and prevents mode mixing.The improved HHT was performed subsequently to obtain the decomposition results shown in Figure 4.Each IMF component was separated clearly with an individual modal.After the right decomposition, the instantaneous frequency and instantaneous amplitude curve of each IMF were calculated according to ( 13) and ( 14) as shown in Figure 5. Linear fitting was then conducted separately to obtain each modal natural frequency and damping ratio of the modulated signal.The final modal parameters of the original signal were obtained by demodulation according to (15)     in Figure 6.The data have strong noise in the lower frequency rage.The third mode was not effectively activated due to limited excited energy.The first and second modal parameters including the frequency and damping at each test stage were estimated by the improved HHT method as shown in Table 4; the second mode is the primary flutter mode with smaller damping and higher energy.We conducted extrapolation and curve fitting based on the second modal damping factors; the FBP result is 824 Km/h, as shown in Figure 7, which is consistent with the design value.

Conclusion
An improved HHT method was proposed in this paper for identifying flight flutter modal parameters and predicting flutter boundaries.Typical mode mixing effect issues were resolved by applying a modulation algorithm that solely requires a structural response in the form of measured acceleration.Simulation and flight flutter test data results altogether indicated that this approach is simple, effective, and feasible for flutter test data processing.
(d) Determine whether () meets the two IMF conditions given above.If so, () is an IMF, as   .Otherwise, let () = () and repeat steps (a)-(d) until the conditions are met.(e) Let the residual  +1 = () −   be the new signal to be decomposed, and repeat steps (a)-(d) for other IMFs.Repeat the steps above until the residual is a monotonic signal or falls below a certain threshold; then the decomposition is complete.

Figure 6 :
Figure 6: Spectrums of vibration signals at every airspeed.
and (16), as shown in Table2.Compared with the true values, the results in Table2indicate that the proposed method yields very accurate modal parameters.There is slight (tolerable) error We next applied the proposed method to actual flutter flight test data (China Flight Test Establishment) to determine its practical applicability.Four main acceleration sensors were installed on the auxiliary tank and wingtip of the aircraft symmetrically, and the structural responses were excited by a small rocket as pulse excitations.
improved HHT method.The flutter boundary designed for the aircraft is 850 Km/h.The lowest three natural frequencies measured by the vibration test and the FEM analysis are given in Table3.The signals from the acceleration sensor located at the left wingtip were selected and the power spectral density of every airspeed was determined as shown 6 Mathematical Problems in Engineering

Table 4 :
Changes in modal parameters with wind speed.