新非平稳信号分析方法-自适应径向墨西哥草帽核时频分布

第一作者李志农男，博士，教授，1966年8月生

新非平稳信号分析方法—自适应径向墨西哥草帽核时频分布

李志农，朱明，龙盛蓉(南昌航空大学无损检测技术教育部重点实验室, 南昌330063)

摘要：据自适应最优核设计准则，构造新的自适应径向墨西哥草帽型核函数时频分析方法。该方法能据分析信号分布自适应调节核函数的扩展方向及宽度，使在信号自项方向尽可能延伸，而不在信号自项方向的互项尽可能抑制，能克服传统时频分布中固定核函数缺乏对信号自适应能力的不足。给出自适应径向墨西哥草帽型核函数时频分布定义及算法，并与传统的短时傅里叶变换、Wigner-Ville分布及小波变换方法对比分析。仿真结果表明，所提方法优于传统的时频分布方法，能更有效处理非平稳信号，可获得较传统时频分布分辨率及抗噪更高性能。通过用于转子裂纹故障诊断实验结果表明，该方法能有效反映裂纹故障的严重程度。

关键词：径向墨西哥草帽核；自适应时频分布；自适应优化；故障诊断；裂纹转子

基金项目：国家自然科学基金(51075372，50775208， 51265039);江西省教育厅科技计划项目(GJJ12405);湖南科技大学机械设备健康维护湖南省重点实验室开放基金(201204) ;江西省研究生创新基金项目(YC2013-S214)

收稿日期：2014-02-08修改稿收到日期：2014-05-23

中图分类号:TN911;TH165.3;TH17;TP391文献标志码：A

A new method for non-stationary signal analysis-adaptive radial mexican-hat kernel time-frequency distribution

LIZhi-nong,ZHUMing,LONGSheng-rong(MOE Key Laboratory of Nondestructive Testing, Nanchang Aeronautical University, Nanchang 330063, China)

Abstract:A new adaptive radial Mexican-hat kernel (RMK) time-frequency distribution method was proposed according to the design criteria of adaptive optimal kernel. The characteristic of the proposed method was that RMK could self-adaptively adjust the expansion direction and width of the kernel function according to the distribution of the analyzed signal. The kernel function expanded as far as possible in auto-term direction, and was suppressed as far as possible in cross-term direction to overcome the deficiency of lack of self-adaptability for the fixed kernel function in the traditional time-frequency distribution. Here, the definition and the algorithm of RMK were presented, and the proposed method was compared with the traditional time-frequency distribution, such as, short Fourier transformation, Wigner-Ville distribution and wavelet transformation. The simulation results showed that the proposed method is superior to the traditional time-frequency distribution, it can more effectively be used to process non-stationary signals, and obtain the higher time-frequency resolution and anti-interference performance. Finally, the proposed method was applied in the fault diagnosis of rotor crack. The test results showed that the proposed method is very effective and can reflect the severity of rotor crack fault.

Key words:radial mexican-hat kernel; adaptive time-frequency distribution; adaptive optimization; fault diagnosis; cracked rotor

工程中所得动态信号往往具有非平稳性。由于非平稳动态信号的统计特性与时间有关，故对其须同时进行时、频分析[1]。时频分析能提供频谱内容随时间变化的信息，为非平稳信号分析的有力工具[2]。然而，传统的短时傅里叶变换、小波变换、Cohen类双线性时频分布等分析方法均采用固定的核函数，对信号缺乏自适应性，导致分析结果往往失去实际物理意义而不能正确提取信号本质特征[3-4]。Baraniuk等[5-6]针对双线性时频分布中核函数固定不变、对信号缺乏自适应缺点，提出自适应优化核函数的设计准则，大大增强对非平稳信号的处理能力，可据不同类型信号进行自适应处理[7-11]。本文针对固定核缺陷，据自适应优化核设计准则，提出新的自适应墨西哥草帽核时频分布。该方法能随信号模糊函数不同而自适应调节核的方向，使之与模糊函数的自项相匹配且最大限度扩展自项方向核函数，对其它方向进行最大抑制，自项获得最大保留，交叉项得到较好抑制，时频分辨率高，亦能较好抑制噪声。将该方法与传统时频分析方法在时频分辨率及抑制噪声能力方面进行对比分析，并用于裂纹转子故障特征分析，实验结果验证该方法的有效性。

1自适应径向墨西哥草帽核时频分布构造

由于直角坐标系中传统二次型时频分布的核函数Φ(θ,τ)固定不变，不能随信号不同而改变，对信号缺乏自适应性。为此，本文据自适应最优核设计准则，构造出二维径向非增低通的墨西哥草帽(Mexican-hat)核函数，其极坐标形式为

(1)

该核函数为径向非增低通函数，随远离原点而光滑下降，引入参数ξ使核函数通带与阻带之前的宽度可调；亦可据模糊域中信号的自项方向，自适应调节核方向与自项方向匹配，尽可能使核在自项方向扩展,其它方向尽量收缩，使信号自项获得最大的保留，交叉项得到最大抑制；且能较好抑制噪声。

由式(1)可知，二维径向函数是一维扩展函数的参数化表示，可将二维函数的优化求解问题转换成一维函数的优化求解问题，优化目的简化为寻找最优的一维扩展函数σ(ψ)。径向Mexican-hat核函数的优化设计准则在极坐标形式下可表示为

(2)

约束条件为

(3)

(4)

式中：AF(r,ψ)为信号模糊函数极坐标形式。

式(3)将核函数限制为径向Mexican-hat函数形式，式(4)中a限制核函数的能量体积，且直接影响自项的聚集程度及交叉项的抑制程度。因模糊函数关于原点对称，式(3)中σ(ψ)积分范围可变为0≤ψ≤π，故对式(4)化简为

(5)

通过符号积分可得

(6)

式中：

η=[2+ξ2-2ξ-2exp(-ξ)]ξ-3

(7)

所以

(8)

对自适应最优Mexican-hat核函数的求解转变为优化求解问题，即

(9)

约束为式(4)、式(8)。由式(9)可知，要实现最大化，模糊函数AF(r,ψ)越大，核函数Φ(θ,τ)越大。然而受a限制，其核函数能量体积有限。为实现极大化优化准则，核函数Φ(θ,τ)将尽可能对靠近模糊域原点的自分量方向作远的扩展而对远离模糊域原点的互分量方向尽可能收缩。参数a通过控制最优核体积决定交叉项抑制程度及自分量聚集程度。文献[3-4]给出a的选择原则，即通常情况下取1≤a≤5。a越小对互分量抑制越厉害，但会滤掉一些自分量，使分辨率降低；a越大时频分辨率越高，但互分量则得不到较好抑制。因此，利用该方法需折中选择。经多次试验，本文认为选a=2最合适。

在所提出算法中需考虑调节系数ω，ξ的选择。由式(3)知，ω及ξ值对Φ有直接影响，但整个优化选择中a控制着核函数的能量体积，优化时σ2(ψ)初值为式(8)取等号情况。a不变时σ2(ψ)与ω成正比，故其结果对Φ影响不大。而随ξ增大核函数过渡带变小，形状亦会发生一定变化，对时频分辨率有一定影响。经大量试验对比证明，ω对时频分布影响不大，本文取ω=0.1。对ξ选择，取ξ=1时能获得较好结果。

自适应Mexican-hat核函数参数确定后，须进行约束条件下优化问题求解。工程应用中时频分布计算总建立在离散化基础上，因此须对优化目标函数及约束条件离散化。通常，模糊函数计算及从模糊域到时频域转换均在直角坐标系中进行，可利用快速傅里叶变换实现。而优化问题在极坐标中计算则更简便。因此计算过程中涉及两个坐标系的尺度变换。自适应径向Mexican-hat核时频分布计算步骤如下：

(1)对直角坐标系下信号的模糊函数AF(θ,τ)进行矩形采样。设采样长度为L，采样间隔为T，则离散模糊函数为

AFd(m,n)=AF(θ;τ)|θ=mΔθ,τ=mΔτ

(10)

离散模糊函数表达式为

AFd(m,n)=

(11)

(2)直角坐标系下不便于信号的模糊函数优化，可通过插值转化到极坐标下。为便于变换，设L×L矩形网格宽高比相等，对两坐标系进行尺度转换，令Δτ=Δθ，则

(12)

(13)

r=pΔr,(p=0,1,…,P-1)

ψ=qΔψ,(q=0,1,…,Q-1)

式中：p，q分别为极半径、极角的采样点。

获得离散化模糊函数极坐标形式AFp(p,q)，此时核函数可表示为

(14)

(3)由上两步得待优化目标函数及约束离散形式

(15)

(16)

由于最优核函数优化过程为非线性过程，直接求解计算量非常大，故采用梯度上升迭代算法求其最优解σq，令

(17)

迭代式为

σ(k+1)=σ(k)+μ(k)f(k)

(18)

(19)

(20)

为使算法收敛，用归一化方法使约束条件成立，即

(4)将优化的最优核函数σq用插值法转换到直角坐标系中。

(5)在直角坐标系中对模糊函数及最优核函数的乘积进行一次FFT及IFFT，即得最优时频分布。

2性能研究

为验证自适应径向墨西哥草帽核时频分布的有效性，构造仿真信号为

(22)

该信号由两个复线性调频信号之和构成，选采样频率fs=100 Hz，采样点数N=256。其时域波形见图1(a)，自适应径向墨西哥草帽核时频分布见图1(b)。为对比分析，给出传统的时频分析方法所得Wigner-Ville分布、小波变换及短时傅里叶变换的时频分布，见图1(c)、(d)、(e)。从时频分辨率及抗噪能力两方面考察自适应径向墨西哥草帽核时频分布处理非平稳信号能力。

2.1时频分辨率分析

对比图1各时频分布结果可知，对线性调频信号，自适应墨西哥草帽核时频分辨率明显高于小波变换及短时傅里叶变换。Wigner-Ville分布虽分辨率最高，但存在严重的交叉干扰项，而自适应墨西哥草帽核时频分布、小波变换及短时傅里叶变换均不存在交叉项。因此由时频分辨率及交叉项抑制来看，自适应径向墨西哥草帽核时频分布明显优于其它时频分析方法。此可从理论上对该结果给出充分解释。

短时傅里叶变换由于受Heisenberg测不准原的约束，时间、频率分辨率不可能同时达到最佳，且其时频窗口大小固定不变，缺乏自适应性。Wigner-Ville 时频分布在模糊域中全通，时频聚集性好，能较好反映信号的时频特征，但存在远大于自项的交叉干扰项，且其时频分布仍采用固定核函数。即适合于某种特定类型信号分析的固定核函数不一定适合其它类型信号分析。况且Wigner-Ville 时频分布缺乏对信号的自适应能力。小波虽具有可调时频聚窗口，但存在一定局限，如小波基函数选择对非平稳信号处理影响较大。故选合适的小波基满足不同类型信号较困难[12-13]，往往只适用某些类型信号，而对其它类型信号不一定适用，对信号缺乏广泛地自适应性。且小波变换的时移、频移是固定变化的，只是对时频平面进行机械的格型分割，因此在小波变换中小波基一经选择，在整个分解、重构过程中均无法更改；小波基可能在全局上是最佳的，但对某个局部可能最差，故小波变换不具备自适应的信号分解特性。

与传统短时傅里叶变换、Wigner-Ville分布、小波变换相比，自适应墨西哥草帽核时频分布能据不同信号模糊函数自适应调节核的扩展方向及宽度，对在信号自项方向尽可能延伸，而对不在信号自项方向的互项尽可能抑制。因此该时频分布自适应性较强。

图1(f)为信号x(t)对应的自适应Mexican-hat核等高图，可见Mexican-hat核能自动随信号方向扩展，而在其它方向收缩，从而能得到最佳的时频分辨率。

(a) 信号x(t)的时域波形　　　　　　(b)自适应径向墨西哥草帽核分布　　　　　　　(c) Wigner-Ville分布

(d) 小波变换　　　　　　　　　　(e) 短时傅里叶变换　　　　　　　　(f) 信号x(t)对应的墨西哥草帽核 图1　信号的时频分布 Fig.1 The time-frequency distribution of signal

2.2抗噪性分析

自适应径向墨西哥草帽核在二维平面上为低通核，而低通核有助于减少噪声敏感度。为验证该草帽核抑制噪声的能力，对信号x(t)分别施加不同信噪比(SNR) 的高斯白噪声进行测试。

分别加入信噪比5 dB、-1 dB的加性高斯白噪声时自适应径向墨西哥草帽核时频分布见图2、图3。为对比分析，亦给出相应的不同噪声环境下传统Wigner-Ville分布、小波变换及短时傅里叶变换时频分布。由图2、图3可知，Wigner-Ville分布存在严重的噪声污染，信号特征完全被淹没，此因信号在模糊域全通；小波变换及短时傅里叶变换时频图均存在不同程度的噪声污染，不能较好分辨信号的时频特征，且随噪声增强，反映信号的时频特征越模糊。原因为两种变换采用的固定小波函数、窗函数所致，缺乏自适应性。而自适应墨西哥草帽核时频分布(图2(b)、3(b))可清晰看出频率随时间线性变化特性，且绝大部分噪声被抑制。明显优于其它时频分布，尤其信噪比越低该优势越明显。其原因为该自适应核在二维平面上为低通核，且能随信号特征朝信号自项方向自动作尽可能延伸，在其它方向尽可能收缩，故噪声能被尽可能抑制。

从时频分辨率及抗噪能力两方面均明显反映出自适应墨西哥草帽核时频分布在处理非平稳信号方面优于Wigner-Ville分布、小波变换、短时傅里叶变换。该草帽核时频分布为优秀的自适应时频分析方法。

3实验研究

为进一步验证方法的有效性，利用自适应径向墨西哥草帽核时频分布对转子裂纹故障进行分析。实验在Bently转子试验台上进行，为分析不同裂纹的自适应径向墨西哥草帽核时频分布特性，加工两根裂纹深度不同的轴，一根裂纹深度为直径的25%，另一根为直

(a) 信号x(t)染噪后的时域波形,信噪比为5dB(b) 自适应径向墨西哥草帽核分布(c) Wigner-Ville分布

(d) 小波变换(e) 短时傅里叶变换图2 各种时频分布Fig.2Thetime-frequencydistribution

(a) 信号x(t)染噪后的时域波形,信噪比为-1dB(b) 自适应径向墨西哥草帽核分布(c) Wigner-Ville分布

(d) 连续小波变换分布(e) 短时傅里叶变换分布图3 各种时频分布Fig.3Thetime-frequencydistribution

(a) 深裂纹转子振动信号时域及频域(b) 自适应径向墨西哥草帽核分布(c) Wigner-Ville分布

(d) 小波变换(e) 短时傅里叶变换图4 浅裂纹的时频分布Fig.4Thetime-frequencyofshallowcrack

(a)深裂纹转子振动信号时域及频域(b)自适应径向墨西哥草帽核分布(c)Wigner-Ville分布

(d)小波变换(e)短时傅里叶变换图5 深裂纹的时频分布Fig.5Thetime-frequencyofdeepcrack

径的50%。裂纹轴加工方法为：用线切割机在轴上指定位置切割出给定深度槽，切槽宽度约0.12 mm；在切槽中嵌入厚度0.10 mm金属片，形成近似开闭裂纹。为防止薄金属片实验中脱落，嵌入切槽前在薄金属片某面涂502胶水与切槽粘牢[14-15]。两条不同程度裂纹转子振动信号的时域图、频谱图、自适应径向墨西哥草帽核时频分布及传统的Wigner-Ville分布、小波变换、短时傅里叶变换时频分布分别见图4、图5。其中图4(a)、图5(a)为不同裂纹深度时的时域图、频域图，可知，频谱中能看出倍频信息，但却看不出信号频率随时间变化特征。由图4(b)、图5(b)看出，浅裂纹时转子系统自适应径向墨西哥草帽核时频分布的一、二倍频占主要成分且持续存在，而随裂纹深度增加，裂纹转子系统除一、二倍频等主要成分仍持续存在外，四、六、八倍等高倍频成分相继出现，较一、二倍频弱，且只在一定时区内存在，与理论分析结果吻合。由图4(c)、图5(c)知，Wigner-Ville分布虽时频分辨率高但存在严重的交叉干扰项，难以看清时频特征。由图4(d)、5(d),4(e)、5(e)知，小波变换及短时傅里叶变换的时频分辨率均不及本文时频分布方法。因此，自适应径向墨西哥草帽核时频分布能有效反映裂纹故障的严重程度。

4结论

针对传统时频分布中固定核函数缺乏对信号的自适应能力的不足，据自适应最优核设计准则，提出新的自适应墨西哥草帽核时频分布。给出该核函数的定义、优化方法及算法，考察该方法对非平稳信号处理能力，并与传统的Wigner-Ville分布、小波变换、短时傅里叶变换进行对比分析,结论如下：

(1)所提时频分布无论对时频分辨率或噪声抑制均优于传统的时频分析方法，因核函数为二维的径向非增低通核函数，且能随信号特征朝信号自项方向自动尽可能延伸，在其它方向尽可能收缩，交叉项及噪声被尽可能抑制。

(2)通过进行转子裂纹故障诊断实验，结果进一步验证所提方法的有效性。

参考文献

[1]何正嘉,訾艳阳,孟庆丰,等. 机械设备非平稳信号的故障诊断原理及应用[M].北京: 高等教育出版社, 2001.

[2]向玲,唐贵基,胡爱军. 旋转机械非平稳振动信号的时频分析比较[J]. 振动与冲击, 2010,29(2): 42-45.

XIANG Ling,TANG Gui-ji, HU Ai-jun. Vibration signal’s time-frequency analysis and comparison for a rotating machinery[J]. Journal of Vibration and Shock, 2010,29(2): 42-45.玲,唐贵基,胡爱军. 旋转机械非平稳振动信号的时频分析比较[J]. 振动与冲击, 2010,29(2): 42-45.

XIANG Ling,TANG Gui-ji, HU Ai-jun. Vibration signal’s time-frequency analysis and comparison for a rotating machinery[J]. Journal of Vibration and Shock, 2010,29(2): 42-45.

[3]Sheng Chen-xing, Li Zhi-xiong, Li Qin, et al. Recent progress on mechanical condition monitoring and fault diagnosis[J]. Procedia Engineering, 2011,15: 142-146.

[4]Zhao X M, Patel T H, Zuo M J. Multivariate EMD and full spectrum based condition monitoring for rotating machinery[J]. Mechanical Systems and Signal Processing,2012,27:712-728.

[5]Baraniuk R G,Jones D L. A signal dependent time-frequency representation: optimal kernel design[J]. IEEE Trans on Signal Processing, 1993, 41(4) :1959-1602.

[6]Baraniuk R G, Jones D L. Signal dependent time-frequency analysis using a radially Gaussian kernel[J].Signal Processing, 1993, 32(3) : 263-284.

[7]Alp Y K, Arikan O. Time-frequency analysis of signals using support adaptive Hermite-Gaussian expansions[J]. Digital Signal Processing, 2012,22:1010-1023.

[8]Jones D L, Baraniuk R G. An adaptive optimal-kernel time-frequency representation[J].IEEE Trans. on Signal Processing,1995,43(10):2361-2371.

[9]Zhao Y, Altlas L E ,Marks R J. The use of cone-shaped kernels for generalized time-frequency representations of nonstationary signals[J].IEEE Transactions on Acoustics, Speech, and Signal Processing,1990,38(7):1084-1091.

[10]Wu D, Morris J M.Time-frequency representations using a radial butterworth kernel[C]. IEEE-SP Int’l. Symp. on Time-Freq. and Time-Scale Anal, Philadelphia, PA. 1994: 60-63.

[11]王晓凯, 高静怀, 何洋洋. 基于时频自适应最优核的时频分析方法[J]. 系统工程与电子技术, 2010, 32(1): 22-26.

WANG Xiao-kai, GAO Jing-huai, He Yang-yang, Time-frequency analysis based on time-frequency-adaptive optimal-kernel[J]. Systems Engineering and Electronics, 2010, 32(1): 22-26.

[12]Yan R, Gao R X. Base wavelet selection for bearing vibration signal analysis[J]. International Journal of Wavelets, Multiresolution and Information Processing, 2009, 7(4): 411-426.

[13]赵学智,陈统坚,叶邦彦,等. 基于参数方程的小波基自适应选择[J]. 机械工程学报, 2005, 40(11): 123-128.

ZHAO Xue-zhi, CHENG Tong-jian ,YE Bang-yan,et al. Adaptive selection of wavelet based on parametr equation [J]. Chinese Journal of Mechanical Engineering, 2005, 40(11): 123-128.

[14]曾复. 裂纹转子非线性振动机理及实验研究[D]. 杭州:浙江大学, 2001.

[15]李志农. 时序模型盲辨识在旋转机械故障诊断中的应用研究[D]. 杭州：浙江大学, 2002.