连续小波变换应用于密集模态参数识别

我的第一篇文章：
摘要：
        用 连续小波变换方法进行模态参数识别时,减小小波函数的带宽可以提高小波变换的频率分辨能力,但同时也会增大边沿效应,影响密集模态参数识别的精度.本文提 出了极值位置方法,在模态耦合严重的小波脊附近,找到模态解耦最好的尺度,在此尺度上进行频率和阻尼的识别,可以获得较好的识别精度.使用Morlet小 波函数,对GARTEUR飞机模型中的三阶密集模态进行识别,验证了该方法的有效性。

全文链接：
柳小勤，岳林，朱如鹏. 连续小波变换应用于密集模态参数识别. 南京航空航天大学学报，2007，39（4），496-500.
并没有很好地解决问题，惭愧。参数的拟合区间始终是个大问题，文中没有提到，因为我不知道该怎么比较好地解决。

欢迎到我的讨论组发表意见：http://groups.google.com/group/liuxq-works?hl=zh-CN

连续小波变换应用于密集模态参数识别

柳小勤，岳林，朱如鹏

关键词：模态分析，小波变换，密集模态， Morlet

中图分类号：032           TB112         文献标识码：A

Applying continuous wavelet transform to close modal parameter estimation

LIU Xiao-qin，YUE Lin，ZHU Ru-peng

(College of Mechanical and Electrical Engineering, Nanjing University of Aeronautics and Astronautics,     Nanjing 210016, China)

Abstract: Continuous wavelet transform based on the Morlet function is used to identify the frequencies and damping ratios of close modes. By minifying the bandwidth of wavelet function, close modes will be easier to extract. But on the other hand, the edge-effect will be larger, and the precision of damping ratios will be of influence. It is very difficult to decompose very compact modes into separated wavelet ridges. Numerical simulation shows that there are two scales beside the composed wavelet ridge, where modes are separated better. A method was developed to find those extreme positions, and then modal parameters are obtained from the wavelet coefficients at those positions. The accuracy of this technique is confirmed by applied to simulated GARTEUR plane model.

Keywords: modal analysis, continuous wavelet transform (CWT), close mode, Morlet

引  言

通过连续小波变换（CWT），可以将多自由度系统的自由响应分解为各个单自由度系统的自由响应，从每个单自由度对应的小波系数中可获得各阶模态的频率和阻尼^[1]，这种方法已经应用于环境激励下建筑、桥梁等的模态参数识别^[2][3]，以及导弹弹翼的颤振实验分析^[4]，具有较好的识别精度。

颤振试验是飞机试验中一个重要环节，试验过程中各种干扰和测量传输噪声对信号造成的污染，使得对密集模态频率和阻尼的精确识别非常困难^[5]。本文将连续小波变换方法应用于GARTEUR飞机模型中密集模态的参数识别，计算精度明显优于正交多项式拟合方法。

1．基本原理

如果函数φ(t) ÎL²(R)，并且满足容许条件

则φ(t)是一个基本小波。对任意的函数f(t)ÎL²(R)的连续小波变换定义为^[6]：

                                                                        

式中：*表示共轭，a和b分别为尺度因子和平移变量。

在模态参数识别中使用的基本小波有Morlet小波，Cauchy小波，和谐波小波等，它们的计算结果比较接近 ^[7]。本文选择最常用的Morlet小波，其时域表达式

（N>0, ω₀ ≥5        ）

频域    对于单自由度系统的自由振动

式中：A为幅值，ζ为阻尼比，ω_n为固有频率，ω_d为阻尼固有频率。其Morlet小波变换可以表示为^[2]  

                     

当a_i = ω₀ / ω_d时，小波系数在整个时间范围内，幅值都达到极值，(a_i, b)形成一个小波脊，W_f (a_i, b)称为小波骨架^[1]。这时 

                                  （1）

                                    （2）

对单自由度系统，通过求极值的方法找出小波脊，－ζω_n对应小波骨架包络线的自然对数的斜率， ω_d对应小波骨架相位的斜率。由ω_d =ω_n sqrt(1-ζ²)，就可以求出系统的固有频率ω_n和阻尼比ζ。

对于具有M个自由度的线性系统，其自由响应

的小波变换可以表示为

                                                （3）

当a_i = ω₀ / ω_di时，在相邻频率间隔比较大的情况下，第i阶模态f_i对小波系数的贡献最大，其它模态对应的小波系数的幅值很小，可以忽略不计，得到

                                                           

这时，在a_i处，小波系数取得局部极值，形成小波脊。这样每一阶模态都对应一个小波脊，利用（1）（2）式，就可以计算出每一阶的频率和阻尼。

但是，当相邻的两阶模态频率f_i，f_i+1的间距f_∆比较小时， f_i+1和f_i的小波系数互相耦合，对f_i

               （4）

上式中的第二项不能忽略，令

( z₁ÎR)，当z₁足够大时，（4）式中的第二项近似为0。得到

                                                                                                       （5）

从上式得到可以解耦的最小频率间隔满足 

在尺度为a_i=ω₀/ω_di时，Morlet小波φ(t)的带宽为^{[7 ]}         

                   （6）

从（5）、（6）式可以看出，通过减小小波函数的带宽（增大ω₀或N），使得z₁取得较大值，可以减小（4）式中的第二项，减小f_i+1对f_i的耦合，提高小波变换对密集模态的解耦能力。

对有限长度的数据进行小波变换时，边沿效应总是存在^[6]。根据（1）（2）式，要对小波系数进行直线拟合来求斜率，边沿效应内不准确的数据显然要排除在拟合区间之外，因此边沿效应应当尽量小。边沿效应的时间长度与持续时间有关，持续时间越长，边沿效应的时间也越长^[8]。

在尺度为a_i=ω₀/ω_di时，Morlet小波的持续时间为^[7]        

       （7）

从 上式可以看出，减小小波函数的带宽，会造成持续时间的延长，边沿效应的增大，影响频率和阻尼识别的精度。可见，在频率过于接近的情况下，通过减小小波带宽 来提高模态解耦能力的方法受到边沿效应的限制，往往不足以将各阶频率分解为单独的小波脊。这时的小波脊含有两阶或更多的模态信息（图1），难以识别出准确的结果，尤其是阻尼的识别误差更大（图2）。

 

      图1 两阶模态耦合的小波系数                                 图 2 含二阶密集频率的小波脊包络线

2. 极值位置方法

在含有耦合模态的小波脊附近，小波系数的包络线呈现较大的波动（图2）。用ln|W(a_i, b)|偏离其拟合直线L(a_i, b)的标准差σ_e来衡量波动幅度的大小

σ_e (a_i)=std(ln|W(a_i, b)|- L(a_i, b))                                                                            （8）

式中：std表示标准偏差。在尺度a_i处，f_i+1对f_i的耦合越严重，ln|W(a_i, b)|的波动幅度越大，偏离拟合直线的幅度也越大，σ_e也越大；相反，f_i，f_i+1的解耦越好，则ln|W(a_i, b)|的波动幅度越小，σ_e也越小。

两个单自由度系统f₁，f₂，分别进行小波变换得到W_f1(a, b)和W_f2(a, b)。将这两个系统叠加，得到一个两自由度系统，其小波变换：W(a, b)= W_f1(a, b)+W_f2(a, b)。

对W_f1(a, bi)和W_f2(a, bi)求比值      R_W(a)=| W_f2(a, bi)/ W_f1(a, bi)|  

通过数值计算，可以得到如图3所示的R_W分布曲线（与平移变量b_i无关）。该曲线在频率分别为f′和f″时达到最小值和最大值。在f′处，R_W最小，W_f1在W中所占的比重最大，W_f1的影响最小，可以忽略，对此位置的小波系数W(f′, b)进行f₁的参数估计，计算精度更高。同样，可以在f″处进行f₂的参数识别。即在f′和f″处，两阶模态的解耦最好。

 

图3 两个单自由度系统的小波系数及其比值                                      图4 σ_e和R_W的极值位置一致

在f′ 处，其它模态的干扰很小，可以认为W(f′, b)只含有f₁的模态信息，ln|W(f′, b)|与其拟合直线的符合较好，即σ_e在f′处取得极小值。同理，σ_e在f″处也取得极小值（图4）。因此通过求σ_e的极小值可以间接地得到f′和f″，然后在这两个位置进行模态参数估计，从而获得更精确的结果。

3.数值模拟

为了验证上述识别方法的有效性，构造一个两自由度密集模态系统，第一阶频率f₁=39Hz，ζ₁=0.015；第二阶频率f₂=40Hz，阻尼比ζ₂=0.01。对其脉冲响应进行分析，采样频率512Hz，采样时间2s。Morlet小波中取N=4, ω₀=12π。小波系数|W(a, b_i)|如图1所示，两阶频率的小波脊耦合在一起。对f₁和f₂分别进行小波变换得到小波系数W_f1(a, b)和W_f2(a, b)，两者幅值的比值R_W如图3实线所示。两个极值点分别对应频率38.1Hz、40.9Hz。

图5极值位置38.1Hz的小波包络线及其拟合直线

对小波脊两侧的小波系数，根据式（8）， 计算ln|W(a_i, b)| 到其拟合直线的标准差σ_e，其分布曲线中两个极小值点对应的频率分别为38.1Hz、40.8Hz，同R_W的极值点位置基本一致（图4）。根据（1）（2）式，对这两个位置的小波系数进行参数识别，拟合直线与小波系数的符合较好（图5）。识别结果见表1。

表1                      两自由度密集频率识别结果

频率

阻尼

理论值（Hz）

计算值（Hz）

相对误差（%）

理论值（%）

计算值（%）

相对误差（%）

39

39.036

0.09

1.5

1.54

2.67

40

40.004

0.01

1.0

0. 97

3.0

4. GARTEUR模型仿真

GARTEUR飞机模型是欧洲航空科技研究集团于20世纪九十年代设计的，具有真实飞机模态频率密集且低的特点。该模型真实地模拟了实际飞机的动态特性^[9]，第3～5阶为密集模态，频率依次为33.01，33.66和35.14Hz，阻尼比均为1％。

根据模态叠加原理，构造出具有三阶密集模态的系统。对该系统进行线性扫频激励，同时加20%白噪声和20%有色噪声，采样频率f_s=160Hz。由频响函数（图6）的IFFT变换得到脉冲响应函数，Morlet小波取N=4，ω₀=15π。

从图7看出，第三阶的小波系数区分得很明显，而前两阶则耦合在一起。根据极值位置方法对33.5Hz附近的小波系数包络线进行分析，σ_e分布如图8，两个极小值分别对应频率32.3Hz和33.95Hz。由这两个尺度处的小波系数（图9）得到两阶模态的频率和阻尼。三阶模态的识别结果见表2，与正交多项式拟合方法相比，小波变换方法的识别精度更高。

  

图6 加噪声后频响函数幅频曲线                 图7 GARTEUR模型中密集模态的小波系数幅值

   

图8 σ_e分布曲线中的两个极值位置                                       图9 极值位置32.3Hz的小波包络线及其拟合直线

表2          GARTEUR模型的识别结果

连续小波变换

正交多项式曲线拟合

频率（Hz）

相对误差（%）

阻尼（%）

相对误差（%）

频率（Hz）

相对误差（%）

阻尼（%）

相对误差（%）

32.95

0.18

0.91

9.0

33.33

0.97

0.87

13

33.55

0.33

1.07

7.0

34.47

2.41

0.68

32

35.13

0.03

0.98

2.0

35.10

0.11

0.86

14

 

5. 结论

        本文探讨了连续小波变换方法进行密集模态参数识别的一些问题，得到以下结论：

（1）通过减小小波带宽来提高模态解耦能力的方法受到边沿效应的限制，往往不足以将密集模态分解为单独的小波脊。这时的小波脊含有两阶或更多的模态信息，难以直接识别出准确的结果。

（2）通过极值位置方法，寻找相邻频率小波系数比值的极值位置，得到模态解耦对应的最佳尺度，在此尺度上进行模态参数估计，可以提高频率和阻尼的识别精度。

 

参考文献

[1]      Staszewski W J. Identification of damping in MDOF systems using time-scale decomposition[J]. Journal of Sound and Vibration, 1997, 203(2): 283-305.

[2]      Lardies J, Ta M N, Berthillier M. Modal parameter estimation based on the wavelet transform of output data[J]. Archive of Applied Mechanics, 2004, 73(9-10): 718-733.

[3]      Piombo B A D, Fasana A, Marchesiello S, et al. Modeling and identification of the dynamic response of a supported bridge[J]. Mechanical Systems and Signal Processing, 2000, 14(1): 75-89.

[4]      Yua Kaiping, Yea Jiyuan, Zou Jingxiang, et al. Missile flutter experiment and data analysis using wavelet transform[J]. Journal of Sound and Vibration[J]. 2004, 269(3-5): 899–912.

[5]      Kehoe M A. Historical overview of flight flutter testing[R]. NASA TM 4720, 1995.

[6]      杨福生. 小波变换的工程分析与应用[M]. 北京:科学出版社,1999：1-14.

[7]      Le T P, Argoul P. Continuous wavelet transform for modal identification using free decay response[J]. Journal of Sound and Vibration, 2004, 277(1-2): 73-100.

[8]      Kijewski T, Kareem A. Wavelet transforms for system identification in civil engineering[J]. Computer-Aided Civil and Infrastructure Engineering, 2003, 18(5):339-355.

[9]      Balmes E. GARTEUR group on ground vibration testing-Results from the test of a single structure by 12 laboratories in Europe[C]//Proceedings of 15th IMAC, California, USA: SEM，1997: 1346 – 1352.