第六讲参数化时频分析理论方法

1、1896192019872006第五讲参数化变换General Parameterized Time-frequency Transform彭志科: z.peng上海交通大学机械系统与振动研究背景及意义雷达信号旋转机械振动信号非平稳信号波 导波生物电信号瞬时电流、语音信号非平稳信号特点：非平稳信号的频率常随时间变化例1：(0 t 15) f (t) = 10 + 2.5t + t2 / 3- t3 / 30 (Hz) 0.080.060.040.02Fourier变换采样频率120Hz1)能反映信号的频率范围不能反映频率随时间变化的变化规律2)0204060Freq /

2、Hzs(t) = sin æ 2p æ10t + 5 t2 + 1 t3 -1t4 ööçç49120÷÷èèøø非平稳信号例2：1.51.5Fourier变换110.50.50000102030102030Freq / Hz准确反映信号所含频率分量Freq / Hz1)2)不能反映频率分量的时间段x2x1x (t) = ì2sin(6p t) + sin(12p t)0 £ t < 2s2íîsin(12p t) + 2sin(

3、18p t)2 £ t £ 4sx1(t) = sin(6pt) + sin(12pt) + sin(18pt)0 £ t < 2s分析短时变换连续小波变换(STFT)(CWT)D. Gabor 1946J. Morlet 1984Nobel Prize inPhysics 1971非平稳信号Wigner-Ville分布线调频小波变换(WVD)(Chirplet)E.P Wigner 1932S. Mann,J. Ville 1948S. Haykin, 1991短时变换0.35a=1定义：0.150.1本质：加窗变换0.050适用对象：

4、分段平稳信号-50t/sec53例2- x1例2- x2210-1-2-3t/sect/sec012t/sec34Freq/Hz窗函数g (t) =1e-t2 / 4aa2 pa¥- jwtSTFTx (t ,w) = ò-¥ x(t)g(t -t )edt连续小波变换定义：0 = 510.50-0.510.50-0.5-1-1-50t/sec5-50t/sec5本质：变分辨率带通滤波适用对象：局部奇异性信号示例转子碰摩故障典型特征0.40.2000010rad/seca = 0.52010rad/seca = 120母波函数y

5、 (t) = p -1/4e- jw te-t2 /20CWT (a, b;y ) =1¥ x(t)y æ t - b ö dtxa ò-¥ça÷èøWigner-Ville分布定义：本质：瞬时相关函数的适用对象：单分量信号变换示例交叉项多分量信号x单分量信号x2单分量信号x1WVD (t,w) = 1 ò+¥ x* æ t - 1 t öx æ t + 1 t öe- jtwdt x2p-¥ç2÷ ç2&

6、#247;èø èø线调频小波变换定义：调频窗函数本质：加调频窗的变换适用对象：线性调频信号STFT示例脉冲反射波Chirplet反射波信号简化模型s(t) = Ae jwc (t-l )+ jg (t-l )2- j a (t -t )2Y(t,a ) (t) = ga (t -t )e2+¥- jwtxò-¥(t,a )CT (t,w,a) =x(t)Y*(t )edt分析的不足例1信号真10STFTCWT实曲线00不足：1)不能正确反映分量的瞬时幅值集中性较差Chirplet( = 3)2)WVD非线性调频分量特点：频

7、率是时间的非线性函数(A)(B)(C)Opt Exp, 19(2011) 26174(D)J. Ac. Soc. Am. 107(2000), Pt.1App. Ac. 71 (2010) 10701080(E)J S. Vib. 330 (2011)12251243J. G. Con. DY 21(1998)375-382(A) 激光脉冲信号(B) Lamb波信号(C) 鲸鱼声波(D) 水轮机停机振动信号(E) 战斗机机翼测试信号参数化分析-原理线调频信号正弦信号STFTSTFT参数化分析-原理Chirplet 定义新表达Chirplet 工作原理t0TimeFrequencyìz

8、 (t ) = s(t )FR (t )FM (t , t)ïaa2ïFR (t ) = e- jat /2旋转算子íaïFM (t , t) = e jatt平移算子aï A(t) = e- jt2a /2 î旋转平移CTs (t,w,a )+¥- jwt= A(t)ò-¥ z (t )g(t - t)edt参数化分析-原理例3：Chirp Rate = 5 / s = 0 = 5 = 2.5 = -2.5x(t) = sin(2p (10 + 2.5t)t) +sin(2p (12 + 2.5t)t)

9、 (0 £ t £ 15s)参数化分析-原理信号模型:参数化瞬率: IF(t) = w0旋转算子+ f (t) 分析原理图解平移算子FM (t ,t) = e- jy (t )tjætöFR (t) = e ç ò y (t )dt ÷è 0ø- jæw t +t f (t )dt öx(t) = eç0ò0÷èø参数化分析-原理性质定义线性可加 ax + by ® aPTFT + bPTFT xy时移不变频移不变x (t )

10、e jw0t ® PTFT (t,w + w ,y )x0x (t - t0 ) ® PTFTx (t - t0 ,w,y )ìïz (t) = FRFM x(t)ïFM (t , t)= e- jy (t )tíïætöïîFR (t) = e jç ò0y (t )dt ÷èøPTFTx (t,w,y )+¥- jwt= ò-¥ z (t )g(t - t)edt参数化分析-(PCT)多项式调频小波变换

11、旋转算子思想：闭区间上的连续函数可用多项式函数一致逼近平移算子FM (t ,t) = e- jy (t )t示例STFTPCT50the original data45the polynomial approximation 4035302520150246810n+1y (t) = åatk-1k -1k =2jætöFR (t) = e ç ò y (t )dt ÷è 0ø参数化分析-PCT示例-STFTCWTPCTChirplet ( = 3)WVD参数化分析-样条调频小波变换 (SCT)旋转算子平移算子FM

12、 (t ,t) = e- jy (t )t思想：分段多项式逼近，避免Runge现象0.6示例PCTSCT0.40.20-0.2t-0.4y =-0.61 + t4-0.8-505jætöFR (t) = e ç ò y (t )dt ÷è 0øn+1y (t) = åli (t - t ) k-1kik =2if ti £ t £ ti +1,i = 1, l -1参数化分析-SCT示例-SCTCWTSTFTPCTWVD参数化分析-广义Warblet变换(GWT)旋转算子思想：任何可积函数都可用

13、三角函数展开平移算子FM (t ,t) = e- jy (t )t示例STFTGWTNy (t) = åY(w )e jwitii=1jætöFR (t) = e ç ò y (t )dt ÷è 0ø参数化分析-GWT示例-STFTCWTWVDGWT参数化分析-如何确定变换核参数多项式调频小波变换样条调频小波变换广义Warblet变换Y1, Yn él1 , l1 ùê 1n úêúêël l , l l úû1na

14、1,an 参数化分析-变换核参数估计-PCT40+最小二乘法拟合35+30+25+20STFT151005101540+35+30+25+PCT20最小二乘法拟合1510051015éëa1,a1 ùû1néëa 0 ,a 0 ùû1n参数化分析-变换核参数估计-GWTFourier变换STFTGWTFourier变换Freq/HzFreq/HzéëY1, Y1 ùû1n35302520151050102030405060Time/SecéëY0, Y0

15、ùû1n35302520151050102030405060Time/Sec参数化分析-如何分析多频率分量信号？示例-多频率分量信号IF1IF2 IF2(t) = 0.7t +12.5IF1(t) = 30 + p cos(p t 6)jp é60t+12sinæ p t öùs(t) = eêç 6 ÷úëèøû+ e j(0.7p t +25p t ) + n2参数化分析-多频率分量信号分析-融合法1IF1IF22参数化分析-多频率分量信号分析-融合法

16、高通滤波高通滤波æ 0.11110.11110.1111öH = ç 0.11110.88890.1111÷ç÷ç 0.11110.11110.1111÷èø参数化分析-融合法示例蝙蝠回波定位信号STFTPWVDPCTChirplet50-0.05-0.1-0.15-0.2-0.2501234567Time / secx 10-3参数化分析-多频率分量信号分析-信号分解法信号模型:问题:Nx(t) » å xk (t)k =1- jæwö

17、;tòt +y(t ;a k ,a k )dt +fçkk ÷k1lx (t) = eèø0k()(t )yt ;a ,a» fkkk1lkN- jæw t +t f (t )dt +f öx(t) = åeçkò0kk ÷ + nèøk =1参数化分析-多频率分量信号分析-信号分解法+ e j(25p t+0.7p t2 ) + n示例信号 x(t) = e jp éë60t+12sin(p t/6)ùûFR (t

18、) = e- j12p sin(p t /6)FR(t) = 1y0y20.250.81FF (Fs)A- F (s)(F s)RRB-C-y1y0.050.200000050Freq/Hz10050Freq/Hz10050Freq/Hz100Ind:(A)-0.1007 ; (B)-0.4650; (C)-0.6060Ind = ò F (FRs) 4 dwkFR (t) = e- j(0.7p t2 )y1参数化分析-多频率分量信号分析-信号分解法参数估计:初始值设定:342422322030181628142601251

19、015051015Time /sTime /s粒子群优化算法或遗传算法估计算法:Freq /HzFreq /Hz(a 2 ,a 2 )1l20(a1,a1 )1l10(a k ,.,a k ) = arg max (ò F (FRs) 4 dw)1l(a k ,.,a k )k1l参数化分析-多频率分量信号分析-信号分解法相位估计: z(t) = FR x(t) » e j(wkt +fk ) + Dykjw t yk (t) = e k相关分析:( (t ) )t maxt max= arg max rz, ytk® fkkòrz, y(t ) =yk

20、 (t)z(t -t )dt(w ,a k ,a k ) k1l参数化分析-多频率分量信号分析-信号分解法+ e j(25p t+0.7p t2 +0.7) s(t) = e jéëp (60t+12sin(p t/6)+0.5ùû分量1分量2参数估计值12.50464481779405910.69676176552245720.000257614908471参数估计值33.1489358178738701-0.0408034080162722-0.3773289802169103-0.02781192253533440.0165482590746405

21、-0.0006590658501026-0.00010817446963370.000009353945751参数化分析-如何刻画局部频率时延？例:频散曲线瞬率曲线局部频率时延曲线4030352530202515201015510000510155101520Time / secFreq / HzFreq / HzTime / sec参数化分析-参数化频率时延分析定义工作原理图解ìïz (w) = FRFM X (w)ïFM (W,w)= e- jy (W)wíïæwöïîFR (w) = e j

22、1; ò0 y (W)d W÷èøPTFTx (w, t,y )+¥= òz (w)g(W - w)e jWtdW-¥参数化分析-示例STFTWVD多项式频率时延变换ChirpletCWT参数化分析-应用应用一：转子瞬时转速估计实验装置起停机过程信号电机DAiscpclac/memmen/st/2mm度传转子感器0.40.20-0.2-0.405101520Time/Sec参数化分析-应用应用一：转子瞬时转速估计STFTWVD瞬时转速曲线GWTCWTFreq/Hz302826242205101520Time/SecTFR Peak Estimated Speed参数化分析-应用应用二：水轮机振动信号精细特征提取水轮机简图123456789上导轴承发电机转子推力轴承 主轴轴承蜗壳导叶水轮机叶片尾水管传感器测量位置V1 上导轴承V2 推力轴承V3 轴承参数化分析-应用应用二：水轮机