自适应分解层数的小波域中值滤波振动信号降噪法

自适应分解层数的小波域中值滤波振动信号降噪法（完整版）实用资料(可以直接使用，可编辑完整版实用资料，欢迎下载）

振动与冲击自适应分解层数的小波域中值滤波振动信号降噪法（完整版）实用资料(可以直接使用，可编辑完整版实用资料，欢迎下载）第28卷第12期JOURNALOFVIBRATIONANDSHOCKVol.28No.122021自适应分解层数的小波域中值滤波振动信号降噪法基金项目:国家自然科学基金资助项目(50678052收稿日期:2021-12-25修改稿收到日期:2021-02-05第一作者滕军男,教授,博士生导师,1962年12月生通讯作者滕军滕军1,朱焰煌1,周峰2,李惠2,欧进萍2(1.哈尔滨工业大学深圳研究生院,深圳518055;2.哈尔滨工业大学土木工程学院,哈尔滨150001摘要:为降低结构健康监测加速度信号中常见的白噪声以及脉冲噪声,提出了中值滤波与小波阈值降噪相结合的方法。分解层数对降噪效果有着重要影响,为取得更好的降噪效果,提出一种分解层数自适应确定法,并给出了各层阈值的取值方法。数值模拟以及国家游泳中心健康监测系统实测数据分析结果表明,冲噪声拟制效果不理想的缺陷,能同时有效降低白噪声以及脉冲噪声,关键词:降噪;小波域中值滤波;分解层数;模极大值;中图分类号:O3271文献标识码:A题。,所获取的振动信号中不可避免地存在一些噪声,其中白噪声和脉冲噪声就是常见的两种噪声。如果不能对信号进行有效的降噪,将会影响到后续模态参数识别等的精度。小波变换具有时频局部化、多分辨率、解相关等特点,因而被广泛应用在信号降噪等领域。基于小波变换的降噪方法主要有模极大值法、空域相关法和阈值法,其中小波阈值降噪方法由于实现简单、计算量小等特点,应用最为广泛[1]。相关学者对该方法中阈值的确定、阈值函数的选取等做了一些改进[2-5],在一定程度上提高了信号的信噪比。基于小波阈值降噪的方法能有效降低信号中白噪声,但对脉冲噪声的抑制效果通常不明显[6],而脉冲噪声也是结构健康监测振动信号中常见的一种噪声,因此有必要研究一种在降低白噪声的同时也能有效抑制脉冲噪声的降噪方法,本文提出了小波阈值降噪与中值滤波相结合的降噪方法。针对分解层数对降噪效果有着重要影响的问题,本文根据噪声与信号在不同分解层上的模极大值特性不同,提出一种分解层数自适应确定的方法,并给出了各层阈值的取值方法。数值模拟以及国家游泳中心健康监测现场实测数据分析表明,所提出的方法能有效降低振动信号中的白噪声以及脉冲噪声,提高信号的信噪比。1降噪方法及其改进111噪声来源在信号传输的所有中间环节,任何电流的存在都会产生热噪声以及其他噪声,产生信号的传感器、处理信号的仪器设备都能将这些噪声引入到信号中。另外传感器常处于恶劣的外部环境中,干扰源类型多而复杂,如交流电源及其谐波干扰,雷电、无线电等的射频干扰、电台电视台的信号干扰等。通常情况下噪声可以分为白噪声、有色噪声和脉冲噪声等。其中白噪声和脉冲噪声是电子仪器仪表中常见的噪声。一维带噪信号通常将噪声以线性叠加的方式作用于真实信号,噪声信号可表示为:x(t=s(t+n(t(1其中s(t为纯信号,n(t为噪声。112中值滤波中值滤波是非线性滤波的典型代表,具有良好的边缘保持特性和抑制脉冲噪声的能力。该方法的本质是一种滑动窗口滤波器,滤波操作是使滑动窗口中心位置的信号抽样值取代当前窗口内的所有抽样的中位值。对于窗口长度n=2k+1(k为正整数的中值滤波器,设第i时刻输入信号序列在窗口内的样本为(xi-k,…,xi,…,xi+k,则此时中值滤波器的输出为:yi=med(xi-k,…,xi,…,xi+k(2式中med(表示取窗口内所有数的中位数。当窗口宽度n=2k+1时,信号序列中脉冲宽度小于或等于k的脉冲信号,经过中值滤波后将被去除;当脉冲宽度大于或等于k+1时,经过中值滤波后脉冲则被保留。这就是中值滤波能去除脉冲噪声而保持信号细节的原因[7]。中值滤波对脉冲宽度满足一定条件下的脉冲噪声有较好的抑制作用,但对随机噪声的抑制能力不够理想。1.3小波阈值降噪设x(k为噪声信号x在k点的离散采样值,令c0,k=x(k,则根据Mallat的小波分解塔式算法,x的正交小波变换的分解公式为:cj,k=∑ncj-1,n・hn-2kdj,k=∑ndj-1,n・gn-2k(3式中,k=0,1,2,…,N-1;N为离散采样数据长度;cj,k为尺度系数;dj,k为小波系数;h、g是一对正交镜像滤波器;j为分解尺度。小波变换是线性的,含噪信号x(t经过小波变换后的小波系数模型可以表示为[8]:y(i=θ(i+z(i(4式中y=wx,θ=ws,z=wn,分别为含噪信号、纯信号以及噪声的小波系数。Donoho[9]是:系数予以保留,,得到小波系数的估计,,对这些估计系数进行重构从而达到降噪的目的。对于脉冲噪声,其能量主要集中在小区域的范围上,其小波系数并不满足远远小于信号的小波系数的条件,因此Donoho的阈值处理方法对脉冲噪声来讲是无效的[6]。1.4改进的降噪法为了同时降低振动信号中的白噪声以及脉冲噪声,结合小波阈值降噪和中值滤波的特点,本文提出了一种小波域中值滤波的降噪方法,具体步骤如下:(1对含噪信号进行适当层数的小波分解;(2由于噪声主要集中在小波细节系数中,对每层的细节系数进行中值滤波,尺度系数不做处理;(3选取合适的阈值,对中值滤波后的细节系数进行阈值降噪处理;(4对经过步骤(3处理后的细节系数和尺度系数进行小波重构,重构后的信号即为降噪后的信号。2分解层数及阈值的自适应确定211分解层数的自适应确定信号和噪声的模极大值在小波域有着不同的性态表现。噪声模极大值的幅值及稠密度随着分解层数的增大而减小,而信号模极大值的幅值随着分解层数的增大而增大,且其个数基本不变。因此可以根据模极大值的幅值来判断小波分解后各层序列是否是噪声控制为主,进而自适应地确定合理的分解层次。对信号f进行m层分解(m一般取3~5,每层细节系数的模极大值序列记为wpeak(j,j=1,…,m。如果满足下式,则可认为信号f是有用信号为主,否则可认为信号f是噪声为主。wpeak(12Φ2wpeak(22Φ…Φwpeak(j2(5分解层数的自适应确定按以下步骤实现:(1对含噪信号进行1层小波分解;(2保留步骤(1中得到的尺度系数c1,k,对细节系数d1,k进行重构,得到新信号sig2new,若sig2new的模极大值满足式(5,则跳到步骤(4,否则对c1,k进一层分解;(3重复步骤(2,即每分解一层对该层的细节系数进行重构得到新信号,检验该新信号的模极大值是否满足式(5,,dj,k重构信号的模极(5;,假设最后一次分解的层数为j,则最终确定的分解层数为j-1。212阈值的选取及处理Donoho等[10]给出了通用阈值的选取,但这种方法在实际应用中效果并不理想,将产生过扼杀现象。本文采用Birge2Massart的惩罚算法来进行阈值选取,该惩罚标准crit(t定义如下:crit(t=-sum[d(k2]+2σ×t×[alpha+log(n/t](6式中d(k表示小波系数,kΦt,该系数按照绝对值递减的顺序存储,n代表系数的个数,alpha为调整参数,为大于1的实数,典型值为2,σ代表估计的噪声的标准值,σ的取值如下:σ=ni=1dji0.6745n(7由上可知,式(6是关于t的函数,设t为tmin时crit(t取最小值,则阈值λ可以表示为:λ=d(tmin(8对阈值的处理一般分为硬阈值函数法、软阈值函数法等。由于硬阈值函数是不连续函数,会产生一些间断点,故在实际应用时常采用软阈值函数法。本文采用软阈值函数法,软阈值函数的表达为:y(x=sgn(d(d-λ,dΕλ0,d<λ(93仿真试验311抑制白噪声以及脉冲噪能力的对比为了验证本文所提方法的优良性,以Matlab软件中带噪声的信号进行仿真分析,并与小波阈值降噪以及中值滤波方法对比。选用‘sym8’小波对四种信噪比情况下的Heavy2sine信号进行降噪,原始信号、加噪信号以及各种方法降噪后的信号时域图如图1所示,以信噪比以及均方误差作为降噪效果的衡量标准,各种方法的降噪效果如表1所示。95第12期滕军等:自适应分解层数的小波域中值滤波振动信号降噪法图1各种方法降噪效果比较表1各种方法降噪效果对比原信号信噪比中值滤波法信噪比/均方差小波阈值法信噪比/均方差本文方法信噪比/均方差5.161710.8569/0.286516.1415/0.155919.7336/0.10319.90915.4830/0.168221.4339/0.084823.5381/0.066515.267220.7913/0.091325.4003/0.053727.6866/0.041319.444725.2312/0.054827.3698/0.042829.2785/0.0344从图1中可以看出,中值滤波法可以有效抑制脉冲噪声,但对白噪声的抑制效果不明显;小波阈值降噪法能有效抑制白噪声,但对脉冲噪声的抑制效果不理想;本文所提出的方法,能同时有效抑制白噪声和脉冲噪声。由表1可知,本文所提出的方法在信噪比和均方误差两个方面均表现出了良好的优越性。图2网壳结构模型图3.2分解层数对降噪效果的影响为验证本文所提出的分解层数自适应法的可行性,现对模拟环境激励(白噪声激励下的网壳模型振动信号进行分析。图2为施威德勒型单层球面网壳模型,跨度为40m,矢高为5m,矢跨比为1/8,均布质量为500kg/m2。杆件均采用钢管,其中径向杆采用<133×4mm,环向杆和斜杆均采用<127×3mm。弹性模量E=2.06×1011N/m2,泊松比μ=0.3,密度ρ=7800kg/m3。提取结构振动响应数据,在响应数据中人工添加白噪声和脉冲噪声,模拟实际信号中的噪声,信号如图3所示。图3原始信号及加噪信号图表2不同分解层数降噪效果对比分解层数中值滤波法信噪比/均方差小波阈值法信噪比/均方差本文方法信噪比/均方差12.1698/0.7790-2.6767/1.3609-2.4688/1.328722.1698/0.7790-0.1139/1.01320.4115/0.953732.1698/0.77902.0399/0.79073.0846/0.701142.1698/0.77903.0174/0.70654.4317/0.600452.1698/0.7790-1.0016/1.1222-0.3981/1.046962.1698/0.7790-0.9240/1.1122-0.3009/1.035206振动与冲击2021年第28卷图3中加噪信号的信噪比为-5.3871dB,根据本文提出的分解层数自适应确定法,得到最优分解层数为4。不同分解层数上的降噪结果如表2所示。由表2可知,分解层数对小波降噪法的降噪效果有着重要的影响。分解层数过多,对每层小波系数都进行阈值处理,会造成有用信号信息的丢失,信噪比反而降低;当分解层数过少,就不能有效地去除噪声,信噪比不能得到很大的提高。当分解层数为4时,取得最优降噪效果,与本文提出的自适应层数确定法确定的分解层数一致,说明本文所提的方法是可行的。由表2还可以看出,在最优分解层数下,本文方法的降噪效果要优于中值滤波法以及小波阈值法,降噪后的信号如图4所示。从图4中可以看出,白噪声和脉冲噪声,。图4利用本文方法降噪后的信号图4工程应用国家游泳中心(“水立方”钢结构工程采用国内外首创仿水泡组合的全新结构体系。根据监测方案,在“水立方”屋盖下弦球节点上安装了加速度传感器。加速度传感器的现场安装如图5所示。图5双轴加速度传感器现场安装图现场采样频率为100Hz,采样时间为1000s,环境激励下加速度传感器所采集到的竖向振动信号如图6所示。经理论计算,国家游泳中心钢结构的前100阶自振频率均小于5Hz,但由图6(b可知,原始信号中有大于5Hz的成分,定性分析认为原信号中可能含有噪声,从图6(a可以看出,原始信号中可能含有脉冲噪声、白噪声等成分。利用本文提出的分解层数自适应确定法,得到最优分解层数为4。前4层上信号的小波系数模极大值均随着分解层数的增大而减小,符合噪声特性,第5层上信号的小波系数模极大值随着分解层数的增大而增大,符合有用信号特性。以上分析表明,原始信号中确实含有噪声和有用信号。利用小波域中值滤波法,对原始信号进行降噪,降噪后的信号及其功率谱如图7所示。从降噪前后的信号图及其功率谱图可知,降噪后的信号保留了有用信号成分,而滤除了杂波干扰,可为健康监测中模态参数识别等提供较好的分析数据。5结论脉冲噪声和白噪声是电子仪器仪表中常见的噪声,小波域中值滤波降噪法结合了中值滤波和小波阈16第12期滕军等:自适应分解层数的小波域中值滤波振动信号降噪法值降噪各自的特点,能同时有效拟制白噪声以及脉冲噪声。根据信号和噪声的模极大值在不同分解层数上传播的特性不同,可以自适应地确定出阈值降噪时的分解层数。数值仿真和国家游泳中心健康监测系统实测数据分析表明,本文所提出的方法能有效地提高信噪比,同时较好地保留了原有信号的主要特征,进一步提高了中值滤波和小波阈值降噪方法的性能。参考文献[1]潘泉,张磊,孟晋丽,等.小波滤波方法及应用[M].北京:清华大学出版社,2005.[2]赵晓燕,李宏男.一种改进的小波分析消噪方法及其在健康监测中的应用[J].振动与冲击,2007,26(10:137-139.[3]陶红艳,秦华峰,余成波.噪算法的研究[J].,(195.[4]刘杰,朱启兵,李允公,等.基于新阈值函数的二进小波变换信号去噪研究[J].东北大学学报(自然科学版,2006,27(5:536-539.[5]张维强,宋国乡.基于一种新的阈值函数的小波域信号去噪[J].西安电子科技大学学报(自然科学版,2004,31(2:296-299.[6]彭玉华.一种改进的小波变换阈值去噪方法[J].通信学报,2004,25(8:119-123.[7]朱凤岗,朱雪田.自适应中值滤波用于色谱信号去噪的研究[J].理化检验———化学分册,2004,40(4:200-203.[8]陈刚,张朝霞,,2006,39(7:797-800.[]DL.bySoft2Thresholding[J].IEEE1995,41(3:613-627.[]果.基于小波分析的转子故障信号自适应降噪技术研究[J].航空动力学报,2021,23(1:9-15.(上接第37页图5三种方法对调幅信号的降噪效果表3三种方法对含有高斯白噪声调幅信号降噪效果降噪前信噪比降噪后信噪比EMD阈值EMD2SG改进的EMD方法1.083.653.643.695结论本文提出了一种改进的EMD降噪方法,基于EMD阈值方法与Savitzky2Golay滤波方法相结合的EMD降噪方法,对频率较高的内蕴模态函数采用阈值降噪的方法,对频率较低的内蕴模态函数采用Savitzky2Golay滤波降噪的方法,这样既能较好地保持信号的高频部分,又能保持信号低频部分的光滑特性。对数值仿真信号的实验表明,改进的EMD降噪方法在信噪比方面要优于单独使用EMD阈值方法或Savitzky2Golay滤波方法。参考文献[1]张贤达.现代信号处理[M].北京:清华大学出版社,2002.[2]DonohoDL.Adaptingtounknownsmoothnessviawavelet2shrinkage[J].J.Amer.Statist.Assoc.,1995,90:1200-1224.[3]HuangNE,ShenZ,LongSR,etal.TheempiricalmodedecompositionandtheHibertSpectrumfornonlinearandnon2stationarytimeseriesanalysis[C].Proc.R.Soc.Lond.A,1998,454:903-995.[4]BoudraaAO,CexusJC,SaidiZ.EMD2BasedSignalNoiseReduction[J].InternationalJournalofSignalProcessing,2004,1(1:33-37.[5]赵雯雯,曾兴雯.一种新的EMD去噪方法[J].电子科技,2021,21(5:30-32.[6]HuangNE,ShenZ,LongSR.Anewviewofnonlinearwa2terwaves:thehibertspectrum[J].Annu.Rev.FluidMech.,1999,31:417-457.[7]PatrickFlandrin,GabrielRilling,PauloGoncalves.Empiri2calModeDecompositionasaFilterBank[J].IEEESignalProcessingLetters,2004,11(2:112-114.[8]DonohoDL,JohnstoneIM.Idealspatialadaptationbywave2letshrinkage[J],Biometrika,1994.81(3:425-455.[9]DonohoDL.De2noisingbysoft2thresholding[J].IEEETrans.Inform.Theory,1995,41(3:613-627.[10]BorgGG,HarrisJH.ApplicationofPlasmaColumnstoRa2diofrequencyAntennas[J].AppliedPhysicsLetters,1999,74(22:3272-3274.26振动与冲击2021年第28卷202JOURNALOFVIBRATIONANDSHOCKVol.28No.122021cussedindetailbasedonthestructuraldamagealarmingmethodforlong2spansuspensionbridges.Thesystemadoptsthedamagealarmingmethodsbasedonenvironmental2condition2normalizedmodalfrequencyandwaveletpacketenergyspec2trum.Theonlineandofflinedamagealarminganalyseswererealizedusingthetwo2stagesdesignmethod.Thesystemcom2biningthemodesofClient/ServerandextendedBrowser/Serverwassuggestedtoconnectthedamagealarmingsystemandothersubsystemsinhealthmonitoringsystemefficientlyandeffectively.ThedamagealarmingsystemsoftwareofRunyangsuspensionbridgewasdevelopedusingtheLabVIEWsoftware,whichincludesthreeparts,i.e.,onlineidentificationofmodalfrequencyandwaveletpacketenergyspectrum,correlationanalysisandstatisticalmodelingofdamagealarmingpa2rametersandenvironmental2conditions,detectionandalarmingofstructuralabnormalconditions.Keywords:suspensionbridge;healthmonitoringsystem;damagealarmingsystem;modalfrequency;wavelet(pp:53-57packetenergyspectrumTENGJun,ZHUYan2huang,ZHOUFeng,LIHui,OUJin2ping1122Vibrationsignaldenoisingmethodbasedonmedianfilterinwaveletdomainwithself2adaptiveleveldecomposition2(1.HarbinInstituteofTechnologyShenzhenGraduateSchool,Shenzhen518055,China;2.HarbinInstituteofTechnologyCivilEngineeringSchool,Harbin150001,ChinaAbstract:Inordertoreducethewhitenoiseandpulsenoiseintheaccelerationsignalofcivilstructuralhealthmoni2toring,acombinativemethodofmedialfilterandwaveletthresholdwasproposed.Becauseofthattheleveldecompositionhasasignificantimpactontheresultofdenoising,inordertogetbetterresult,aself2adaptivemethodtoselectthedecom2positionlevelwasproposed.Themethodtodeterminethevalueofthresholdateachlevelwasalsogivenout.Theanalysisofnumericalsimulationresultandtheon2sitemeasureddatabyuseofthestructuralhealthmonitoringsystemofNationalA2quaticsCenterinBeijing2021Olympicvenuesshowsthat,themethodproposedinthepapercanovercomethedefectofthewavethresholdmethodwhichdoesnπthaveobviouseffecttoreducethepulsenoise.Itcanreducethewhitenoiseandpulsenoiseatthesametimeandprovidedataofhighersignaltonoiseratioformodalparameteridentificationandsoon.Keywords:denoising;medianfilterinwaveletdomain;decompositionlevel;modulusmaximum;structuralhealth(pp:58-62monitoringResearchonconstruction,principle,andperformancesofanoveltwo2degreesoffreedomultrasonicmotorwithanobeliskstatorJINJia2mei,ZHANGJian2hui,ZHAOChun2sheng(PrecisionDriveLaboratory,NanjingUniversityofAeronauticsandAstronautics,Nanjing210016,ChinaAbstract:Anoveltwo2degreesoffreedomultrasonicmotorwasproposed.Themotoriscomposedofanobelisksta2torandabasewithaconcavesphericalsurface.Thestatorispressedtobeincontactwiththeconcavesphericalsurfacebymeansofanopportunepre2pressuresystem.Itiseasytominiaturizethemotor,becausethelocalflexuralvibrationsofthestatorareutilizedastheoperatingmodesofthemotorinsteadofthelongitudinalvibrations.TheellipticalmotionofthetopofthestatorobservedinsimulationwithFEM(finiteelementmethodcanmakethestatormovewhenthetopispressedontheconcavesphericalsurfaceofthebase.Theperformancecharacteristicsachievedbytheprototypeare:thespeedwithoutload145mm/sinx2direction,thespeedwithoutload180mm/siny2direction,thestallthrust4.3Ninx2di2rection,thestallthrust5.0Niny2direction.(pp:63-67Keywords:ultrasonicmotor;two2degreesoffreedom;obeliskstatorPrototypeexperimentandworkingmechanismanalysisofsquareultrasonicmotorbasedonin2planevibrationmodeGUOHui,LIXun2bo,LIHua2feng112(1.SchoolofMechatronicsEngineering,UniversityofElectronicScienceandTechnologyofChina,Chengdu610054,China;2.ResearchCenterofUltrasonicMotors,NanjingUniversityofAeronautics&Astronotics,Nanjing210016,ChinaAbstract:Usingfiniteelementmethod,thevibrationmodeshapesofthesquareultrasonicmotorstatorwerepres2©1994-2021ChinaAcademicJournalElectronicPublishingHouse.Allrightsreserved.://www第四章小波变换降噪分析小波变换是一种崭新的时域(频域信号分析工具。它的发展和思想都来自于傅里叶分析,且在保留了傅里叶分析优点的基础上,较好的解决了时间和频率分辨率的矛盾,在频域与空间域中能够同时具有良好的局部化特性,可进行局部分析。小波去噪的基本原理是根据原始信号和噪声的小波系数在不同尺度上所具有的不同性质,构造相应的规则,在小波域采用其他数学方法对含噪信号的小波系数进行处理。4.1小波变换理论的研究连续小波变换设2((tLRψ∈(2(LR表示平方可积的的空间,即能量有限的信号空间,其傅立叶变换为(ψω。当(ψω满足允许条件(AdmissibleCondition:2(Cφωωω+∞-∞=<∞⎰(4.1时,我们称(tψ为一个基本小波或母小波(MotherWavelet。将母小波函数(tψ经伸缩和平移后,就可以得到一个小波序列。对于连续情况,小波序列为:,((abtbtaψ-=,abR∈0a≠(4.2其中,a——伸缩因子;b——平移因子;——能量归一化因子。这样对于任一信号2011((,(ftbftabdadbCaaφωψ∞∞-∞-=⎰⎰,连续小波变换定义为:,,(,(,(((ababCWTabfttfttdtψ∞-∞==⎰(4.3其逆变换为:2011((,(ftbftabdadbCaaφωψ∞∞-∞-=⎰⎰(4.4离散小波变换实际应用中,尤其是在计算机上实现,如在信号处理领域,必须对连续小波加以离散化。需要强调的是,这一离散化都是针对连续的尺度参数a和连续平移参数b的,而不是针对时间变量t的,这与其它形式的离散化不同。在连续小波中,考虑函数(4.5:,((abtbtaψ-=(4.5这里,,abR∈;0a≠且ψ是容许的,为方便起见,在离散化中限制a取正值,则容许条件变为:2(Cφωωω+∞=<∞⎰(4.6通常,连续小波变换中的尺度因子和平移因子的离散化公式为:000jjaabkab⎧=⎨=⎩(4.7这里,jZ∈,扩展步长01a≠是固定值,且假定01a>。*,,,(,jkjkjkCftdtfψψ∞-∞==⎰(4.8其重构公式为:,,((jkjkjkftCCtψ∞∞=-∞=-∞=∑∑(4.9其中,C是一个与信号无关的常数。然而,怎样选择0a和0b才能够保证重构信号的精度是非常重要的,显然,网格点尽可能密(即0a和0b尽可能小,因为如果网格点越稀疏,使用的小波函数,(jktψ和离散小波系数,jkC就越少,信号重构精确度也就会越低。4.2小波去噪的基本原理的研究原理:根据原始信号和噪声的小波系数在不同尺度上所具有的不同性质,构造相应的规则,在小波域采用其他数学方法对含噪信号的小波系数进行处理。其本质是减小甚至完全剔除由噪声产生的系数,同时最大限度地保留真实信号的系数。小波变换的基本步骤如图4.1所示:图4.1小波去噪的基本步骤4.3小波分析对染噪矩形信号处理处理结果如下图所示:参考信号染噪信号图4.2染噪的矩形波形信号进行降噪处理结果可见,经小波变换降噪之后,噪声水平明显下降,信噪比得到提高。4.4某检测信号降噪读取MATLAB中专用检测信号,对其加入噪声然后进行小波变换降噪,所得结果如下图:原始信号降噪后的信号图4.3染噪的某检测信号进行降噪处理结果可见降噪效果良好。4.5正弦信号降噪生成一段正弦信号,对其加入高斯噪声,进行小波降噪,所得结果如下图所示:原始信号染噪信号消噪信号图4.4染噪的某正弦信号进行降噪处理结果可见,降噪之后信号平滑,信噪比得到改善。4.6带突变信号的小波降噪读取MATLAB中带突变的一段信号,对其加入噪声,然后进行降噪处理,所得结果如图:原始信号噪声信号-信噪比为3降噪信号-SURE降噪信号-Fixedform阈值De-noisedsignal-Minimax图4.5染噪的某正弦信号进行降噪处理结果可见,降噪之后,信号的突变也消失了,从而产生了失真,因此,该种信号不适合用小波变换进行降噪处理。4.7小波包降噪20040060080010001200原始信号20040060080010001200降噪后信号图4.6染噪的某正弦信号进行降噪处理结果可见,降噪之后,信号平滑性增强,有利于信号检测。4.8高斯噪声语音信号的小波降噪语音信号处理的一个重要应用就是进行语音信号的降噪处理,噪声中最常见的是高斯白噪声,故首先研究对高斯白噪声的降噪。读取一段语音信号,在信号中加入高斯白噪声,所得信号的波形及频谱如下图所示:x104图4.7带有高斯噪声的信号波形及频谱可见,叠加高斯白噪声之后信号的变化比较快,而且无规则,在噪声比较小的情况下,虽然可见信号包络,但是有严重的背景噪声,当噪声比较大时,则可能完全听不到信号。宽带噪声与语音信号在时域和频域上基本上重叠,只有在无话期间,噪声分量才单独存在,因此消除这种噪声比较困难。下面采用小波变换法进行滤波降噪。x104图4.8滤波之后波形及频谱图4.8中左图为小波变换后的信号的时域波形图,与处理前的信号相比可以看出,大量的噪声已经去掉,信号已经基本还原为没有噪声干扰的波形了,效果比较理想。右图为处理后信号的频域图,可以看出高频部分大量的噪声谱已经被滤除了,而低频部分的原信号频谱基本全部保留了下来,因此效果比较理想。4.9单频噪声语音信号的小波降噪在语音降噪中,另外一种比较常见的噪声就是单频哨声,这种哨声可能是由于外界某机器产生,也有可能是录音设备稳定性不好造成的。比如,设备中的放大电路中如果有反馈,很肯就会产生正弦振荡产生单频哨声。取一段信号,叠加单频噪声之后的信号波形如下图所示:x104图4.9带有单频噪声的信号波形及频谱可见,信号变化比较快,从频谱中可见,信号中有较强的单频噪声。试听该带噪音乐会发现其中有异常尖锐的单频哨声。频噪声与语音信号在时域上基本上重叠,但在频域是分离的,因此如果知道单频噪声的频率的范围,在频域将其分离比高斯噪声更简单。下面采用小波变换法进行滤波降噪。x104图4.10滤波之后波形及频谱如上图所示,为对含有单频噪声的语音信号进行处理之后的信号的波形和频谱图,与处理前的信号(如图4.10所示相比,可以看出,左图中,信号时域波形中大量的高频跳动噪声信号已经被去掉,信号已经基本还原为没有噪声干扰的波形了,效果比较理想。右图为处理后信号的频域图,可以看出高频部分强度很大的单频噪声谱已经被滤除了,而低频部分的原信号频谱基本全部保留了下来,因此效果比较理想。上述实验展示了小波变换从噪声数据中恢复信号的问题,从图中可以看出,消噪后的信号大体上显示了原始信号的形状,并明显地除去了噪声所引起的干扰。但是,所恢复的信号与原始信号相比有明显的改变。这主要是因为在进行消噪处理的过程中所使用的分析小波和细节参数阈值不恰当所致。4.9小结小波分析是语音信号去噪的新途径,它提供了一种自适应的时域和频域同时局部化的分析方法,无论分析低频或高频局部信号,它都能自动调节时频窗口,以适应实际分析的需要,可以对信号的任意局部细节加以分析。在过去,我们一般使用短时傅里叶变换在频域内对语音信号进行分析去噪,但是对于白噪声,这种方法的效果往往不尽人意。而小波变换在去除语音信号高斯白噪声时,显现出了优势,它的多尺度变换特点,可以由粗到精地逐步观察信号;适当的选择基本小波,便可以使其在时频两域都具有表征信号局部特征的能力,因此非常有利于信号的分析,基于这种特性,它能有效的区分信号中的突变部分和噪声,从而实现信号的去噪。通过以上实验可以看出,在去除语音信号中的高斯白噪声时,小波变换方法明显优于其它方法。附录小波去噪matlab程序******************************************clearclc%在噪声环境下语音信号的增强%语音信号为读入的声音文件%噪声为正态随机噪声sound=wavread('c12345.wav';count1=length(sound;noise=0.05*randn(1,count1;fori=1:count1signal(i=sound(i;endfori=1:count1y(i=signal(i+noise(i;end%在小波基'db3'下进行一维离散小波变换[coefs1,coefs2]=dwt(y,'db3';%[低频高频]count2=length(coefs1;count3=length(coefs2;energy1=sum((abs(coefs1.^2;energy2=sum((abs(coefs2.^2;energy3=energy1+energy2;fori=1:count2recoefs1(i=coefs1(i/energy3;endfori=1:count3recoefs2(i=coefs2(i/energy3;end%低频系数进行语音信号清浊音的判别zhen=160;count4=fix(count2/zhen;fori=1:count4n=160*(i-1+1:160+160*(i-1;s=sound(n;w=hamming(160;sw=s.*w;a=aryule(sw,10;sw=filter(a,1,sw;sw=sw/sum(sw;r=xcorr(sw,'biased';corr=max(r;%为清音(unvoice时,输出为1;为浊音(voice时,输出为0ifcorr>=0.8output1(i=0;elseifcorr<=0.1output1(i=1;endendfori=1:count4n=160*(i-1+1:160+160*(i-1;ifoutput1(i==1switchabs(recoefs1(icaseabs(recoefs1(i<=0.002recoefs1(i=0;caseabs(recoefs1(i>0.002&abs(recoefs1(i<=0.003recoefs1(i=sgn(recoefs1(i*(0.003*abs(recoefs1(i-0.000003/0.002;otherwiserecoefs1(i=recoefs1(i;endelseifoutput1(i==0recoefs1(i=recoefs1(i;endend%对高频系数进行语音信号清浊音的判别count5=fix(count3/zhen;fori=1:count5n=160*(i-1+1:160+160*(i-1;s=sound(n;w=hamming(160;sw=s.*w;a=aryule(sw,10;sw=filter(a,1,sw;sw=sw/sum(sw;r=xcorr(sw,'biased';corr=max(r;%为清音(unvoice时,输出为1;为浊音(voice时,输出为0ifcorr>=0.8output2(i=0;elseifcorr<=0.1output2(i=1;endendfori=1:count5n=160*(i-1+1:160+160*(i-1;ifoutput2(i==1switchabs(recoefs2(icaseabs(recoefs2(i<=0.002recoefs2(i=0;caseabs(recoefs2(i>0.002&abs(recoefs2(i<=0.003recoefs2(i=sgn(recoefs2(i*(0.003*abs(recoefs2(i-0.000003/0.002;otherwiserecoefs2(i=recoefs2(i;endelseifoutput2(i==0recoefs2(i=recoefs2(i;endend%在小波基'db3'下进行一维离散小波反变换output3=idwt(recoefs1,recoefs2,'db3';%对输出信号抽样点值进行归一化处理maxdata=max(output3;output4=output3/maxdata;%读出带噪语音信号,存为'101.wav'wavwrite(y,5500,16,'c101';%读出处理后语音信号,存为'102.wav'wavwrite(output4,5500,16,'c102';第五章小波图像压缩技术•二维正交多分辨分析•图像压缩简介•EZW编码编•SPIHT编码•EBCOT编码与JPEG2000介绍清华大学计算机系孙延奎2021本讲的目的•通过典型的小波图像压缩算法,了解小波图像编码的基本想与关键技术为小波图像编码的研究与应用奠定基础思想与关键技术,为小波图像编码的研究与应用奠定基础。清华大学计算机系孙延奎2021图像压缩简介•图像压缩的目的•两种压缩类型有损压缩与无损压缩•常用压缩方法RLE(Run-LengthEncoding行程长度编码,(ggJEPG,JEPG-LS,JEPEG2000分形编码清华大学计算机系孙延奎2021小波图像压缩图像压缩简介小波图像压缩的般模型•小波图像压缩的一般模型典型的方法:EZW,SPIHT,EBCOT:•小波图像分解的(一种常用方法:可分离二维小波变换滤波器的选取、边界延拓清华大学计算机系孙延奎2021滤波器的选取•主要考虑具有线性相位的、正则性的、完全重构的、有限长度的双正交小波滤波器。滤波器的对称性◆线性相位:滤波器的对称性◆正则性:保证图像的重构质量。在小波分解中,正则性使信获得较的小波系数表在小波构中可以使信号获得较好的小波系数表示;而在小波重构中,正则性导致平滑扰动。◆完全重构、有限长度:一般希望滤波器的分解重构精度越高越好,而要求滤波器对量化误差的反应越不敏感越好,最重要的是需要恢复信号的主、客观质量越高越好。“一般地,自然图像是高度非平稳信源,其功率谱是非对称的,般地,自然图像是高度非平稳信源,其功率谱是非对称的,要求分解滤波器的长度要长,而高通分解滤波器的长度要短。此外,滤波器长度的非对称性可以使能量较好地集中,特别对于高频子图,而当低通合成滤波器的长度较短时,可以减少重构误差”。OlivierEgger,WeiLi.SubbandcodingofimageusingtilfiltbkIEEETIP19954(4asymmetricalfilterbanks.IEEETrans.ImageProc.,1995,4(4:478-485清华大学计算机系孙延奎2021边界处理图像小波变换:边界处理•由于图像小波分解是将图像数据与滤波器系数进行线性卷积运算,使滤波器的输出样点数大于输入的样点数;但从子带压缩编码的效率看,是希望各子图信号的样点数总数不大于原始图像信号的样点数,且子图信号应能完全重构其原始信号。常用的方法是对信号信号应能完全重构其原始信号常用的方法是对信号进行边界延拓。清华大学计算机系孙延奎2021简单的小波图像压缩方案1方案1:只保留低频部分.2:.方案全局阈值法方案3:保留绝对值较大的若干小波系数:利用Matlab编程即可实现。讨论:这些方案利用了小波变换的什么性质?去相关性!小波子带分解中,小波系数还有什么性质?如何开发更清华大学计算机系孙延奎2021先进的小波压缩算法?小波图像编码中的关键问题如何组织小波系数及其位置信息?如何处理小波系数及其位置信息?清华大学计算机系孙延奎2021不同分辨率子带之间数据的相关性.从下图中可观察到不同分解级小波系数之间具有从下图中可观察到,不同分解级小波系数之间具有相似性。,1,2,3jjjjcddd⎡⎤⎢⎥⎣⎦清华大学计算机系孙延奎2021小波系数的四叉树结构表示小波系数的四叉树结构清华大学计算机系孙延奎2021小波系数的零树特性对于给定的阈值T,若小波系数x满足|x|<=T,则称x关于T是不重要的系数。实验统计表明,若一个小波系数关于T是不重要的,则它的所有子孙以上系数关于T是不重要系数的概率在98%以上。[该结论仅供参考]清华大学计算机系孙延奎2021小波零树编码1992,LewisandKnowles利用小波树中各级子带系数之间的相似性,最早提出了小波零树编码算法。该算法将量化后系数为零的系数的子孙系数都指定为零算法将量化后系数为零的系数的子孙系数都指定为零。特点:利用了不同分辨率子带之间数据的相关性.?问题:?清华大学计算机系孙延奎2021EZW嵌入式零树小波图像压缩技术----EZW编码在1993年,Shapiro提出了小波零树编码算法的改进算法,称为嵌入零树小波编码算法,简称为EZW算法。这是一种简单有效的小波编码算法。受该算法的启发,人们后来开发了更为有效的SPIHT及EBCOT算法。☆基本概念☆编码过程清华大学计算机系孙延奎2021----:小波系数的处理扫描顺序为了出现尽可能多的零树根以压缩代码,同时保证先编码重要信息,扫描从最低的精度开始逐渐向高精度级进行。我们教材中采用的是Morton扫描顺序.清华大学计算机系孙延奎2021清华大学计算机系孙延奎2021小波系数位置的处理----由扫描顺序及小波零树共同完成其中,H表示小波系数的绝对值大于或等于当前的阈值;44×图象L表示小波系数的绝对值小于当前的阈值.清华大学计算机系孙延奎2021多遍扫描编码小波系数图像,其中每一遍扫描包含以下的处理步骤:EZW编码过程:1.选择阈值2.主扫描33.辅扫描54.重新排序5.输出编码信号注解:由于”重新排序”需要记录有关信息,且实践表明,这种操作对图像重构的质量影响不大,所以实践中通常忽略这步操作.为此,以下删除教材中EZW编码算法中的第4步”重新排序”.清华大学计算机系孙延奎20211.选择阈值TTT−对于L级小波变换,EZW算法应用一系列的阈值,,…,011L2/1−=iiTTii来确定小波系数的重要性,其中,为扫描次数,=1,2,….,L-1。初始阈值的选择方法如下:log{||}0,22MaxcijT⎢⎥⎢⎥⎣⎦=320=T清华大学计算机系孙延奎20215.5.输出编码信息编码器输出两类信息:一类是给解码器的信息,包括阈值、主扫描表和辅扫描表;的信息包括阈值及获得的重要系数序列320=T1D1S,:PNZTPTTTTZTTZZZZZPZZ;:1010第二类是用于下次扫描的信息,包括阈值及获得的重要系数序列。小波图像数据63320=T,,小波图像数据。{}PPNP−−−−47,49,34,未排序清华大学计算机系孙延奎2021第二次编码输出结果第二次编码输出结果:a为解码器提供的信息161=T2D2S,:NPTTTTTTTTTTTZZZZ;:100110b为下一次扫描的信息161=T,,小波图像数据。{}63,34,49,47,31,23PNPPNP−−−−−−表5.3二次编码的输出结果T3211/SD0PNZTPTTTTZTTZZZZZPZZ/1010NPTTTTTTTTTTTZZZZ/10011022/SD清华大学计算机系孙延奎2021EZW编码与位平面编码之间的关系----小波系数的表示位平面编码=⋅+⋅+⋅+⋅+⋅+⋅63的二进制表示为111111符号位平面及前三个高位的位平面清华大学计算机系孙延奎2021EZW编码与位平面编码之间的关系通过不同阈值如32,16,8,…对小波系数进行扫描,可以发现最重要的小波系数.这等价于对位平面依次进行小波零树编码.由此实现累进编码,也称为嵌入式编码.具体的第一次主扫描本质上是对第一个位平面编码第一次辅扫描具体的,第次主扫描本质上是对第个位平面编码。第次辅扫描本质上是发现”重要系数”在下一个位平面上的改进位.对应当前阈值下4个重要系数在次高位平面的位值,是必然,是巧合?必然.如543210=⋅+⋅+⋅+⋅+⋅+⋅543210=⋅+⋅+⋅+⋅+⋅+⋅第二次辅扫描对应当前阈值6个重要系数在第三高位平面的位值清华大学计算机系孙延奎2021EZW解码EZW算法具有显著编码性能的原因:(1离散小波变换:去相关性与塔式结构(2零树编码:利用不同分辨率同方向子带之间的相关性(3累进逼近:位平面编码(4自适应算法编码:无损编码方式EZW编解码算法的实现,包括C语言与Matlab程序(不考虑排序步骤:清华大学计算机系孙延奎2021分集规则{}(,,max||2(nijijXncSX∈⎧≥⎪=⎨1若(1nSX=X是重要的⎪⎩0其他(0nSX=X是不重要的,(jiO:节点(i,j的坐标集;点(,j所有孩子标集,(jiD:节点(i,j所有子孙的坐标集;H:所有树根的坐标集。,(,(,(jiOjiDjiL−=,(jiL:节点(i,j所有非直系子孙的坐标集;一般地,{}12,12(,2,12(,12,2(,2,2(,(++++=jijijijijiO清华大学计算机系孙延奎2021分集规则1最初坐标集由{}(,|(,ijijH∈和(,|(,DijijH∈且具有非零子孙组成;{}2若,(jiD是重要的,则,(jiD分成iL及4个单节点,,iOlk∈,(j((j3若是重要的,则分成,(jiL,(jiL4个集,(lkD,(,(,jiOlk∈有序表LIP——不重要系数表;LSP——重要系数表;在LIS中,坐标(,ij代表,(jiD或者,(jiL(ij(ijLLIS——不重要子集表。每一个表项都使用坐标(,ij标识分别用(i,jD和(i,jL表示清华大学计算机系孙延奎2021(2排序扫描目的:编码当前位平面的重要系数由以下两个大的步骤构成:对于当前阈值•如果是重要的系数,则输出“1”及其符号位,其中正、负小波系数的符号位分别采用“1”和“0”表示,然后将该系数从LIP中1顺次检查LIP中的所有小波系数(,ij,对于当前阈值,确定其是否重要删除,并添加到有序表LSP的尾部。•如果是不重要的系数,则输出“0”。清华大学计算机系孙延奎2021例题:排序扫描1输出T=32%对LIP中的每个表项顺次进行处理Is(0,0significant?yes:1,gy1(符号位/将(0,0从LIP中删除,添加到LSP的尾部/LSP={(0,0}LIP={(0,1,(1,0,(1,1}Is(0,1significant?yes:10(符号位LSP={(0,0,(0,1}LIP={(1,0,(1,1}Is100Is(1,0significant?no:Is(1,1significant?no:0清华大学计算机系孙延奎20212CheckthesignificanceofalltreesintheLISaccordingtothetypeoftreetype:2对LIS中的每个表项顺次处理,并对D型表项和L型表项分别采用不同的处理方法,具体算法如下:ggypyp√ForatreeoftypeD:°Ifitissignificant,output1,andcodeitschildren:•Ifachildissignificant,output1,thenasignbitandaddittotheLSP•Ifachildisinsignificant,output0andaddittotheendofLIP.•Ifthechildrenhavedescendants,movethetreetotheendofLISastypeL,otherwiseremoveitfromLIS.°Ifitisinsignificant,output0.√ForatreeoftypeL:°Ifitissignificant,output1,addeachofthechildrentotheendofLISasanentryoftypeDandremovetheparenttreefromtheLISofLISasanentryoftypeDandremovetheparenttreefromtheLIS.°Ifitisinsignificant,output0.清华大学计算机系孙延奎2021例5.6第一次SPIHT编码后输出的信息第一次编码过程完成后,编码器输出两类信息:11100011100010000001010110000、精细扫描位流及三个有序表的初始化信息,即LIP,LIS和LSPLIP00011011={(0,0,(0,1,(1,0,(1,1}LIS={(0,1D,(1,0D,(1,1D}LSP={}2用于下次扫描的信息,包括域值52.、三个有序表LIP,LSP,LIS的当前状态信息清华大学计算机系孙延奎2021SPIHT与EZW编码的简单比较共同点:基于(广义小波零树的位平面编码,形成嵌入式位流.其中,在不同分辨率子带之间引入向下的相关性.嵌入式位流的吸引力在于,希望的压缩等级可以在信源已经压缩之后确定.嵌入的一个重要结果是,信源内任意给定的空间区域的信息内容一般必然是分散(dispersed贯穿在压缩表示中.不同点:重要系数排序信息的显式与隐式处理.系数序信式式清华大学计算机系孙延奎2021电子资源•••清华大学计算机系孙延奎2021EBCOT编码最佳截断嵌入码块编码(embeddedblockcodingwithoptimizedtruncation,EBCOT算法是一种非常有效的小波系数编码方法。清华大学计算机系孙延奎2021EBCOT编码质量可伸缩性SNR清华大学计算机系孙延奎2021或可伸缩性或率可伸缩性EBCOT编码EZW与SPIHT算法采用的零树编码结构,在分辨率之间引入向下的相关性,因此所得到的嵌入式码流具有质量可伸缩性,但不具有分辨率可伸缩性.EBCOT算法采用了将子带样本分成小块,单独对每一个小块进行编码.各种相关性可能存在于一个块内,但不会存在于不同块之间.清华大学计算机系孙延奎2021EBCOT算法介绍:1.嵌入码块编码对每个码块进行编码时,编码器不用其它码块的任何信息,只用该码块本身的信息生成个嵌入码块流块本身的信息生成一个嵌入码块流。这种编码方法很容易生成图像的具有分辨率可伸缩性的编码,并可实现感兴趣区域的随机存取.清华大学计算机系孙延奎2021EBCOT图像压缩算法的两层编码结构流程图由于每个质量层必须包含辅助信息以确定每个码块对该质量层的贡献大小。当层数很多时,对给定的质量层,仅有部分码块对该层有贡献,这样在辅助信息中存在大量的冗余信息。为此EBCOT引入第二编码器(Tier2用于压缩每个质量层为引第编码器用于压缩每个质层的辅助信息。清华大学计算机系孙延奎2021第4卷第11期中国水运Vol.4No.112006年11月ChinaWaterTransportNovembdr2006收稿日期:2006-9-20作者简介:孙涛武汉理工大学土木工程与建筑学院(430070小波变换和希尔伯特—黄变换在时频分析中的应用孙涛刘晶璟孔凡万平摘要:简单介绍了时频分析的基本理论,将小波变换和希尔伯特-黄变换分别应用于几个非平稳信号的分析当中,将二者进行一个简单的比较,最终得出结论。关键词:时频分析小波变换希尔伯特-黄变换中图分类号:TN911.21文献标识码:A文章编号:1006-7973(200611-0111-03一、引言长期以来信号处理的对象局限于确定性信号或是统计量不随时间变化的平稳信号,其有效的分析工具就是Fourier分析,它是一种全局性的变换,无法表达信号的时频局部特性,但非平稳信号的广泛存在是不争的事实。由于受到信号处理理论发展的限制,对非平稳信号的分析过去人们一直是沿用平稳信号的处理方法来作近似,效果当然不够理想。随着研究的深入和科技实践的需要,针对非平稳信号的理论分析已是迫在眉睫。这就是时频分析理论产生的时代背景。时频分析实际上是将一维的时间信号映射到时频(有的是时间尺度二维,可以很好的表示出信号的频率成分随时间的化规律,而这恰恰是非平稳信号分析所需要的。二、小波变换小波变换是一种信号的时频分析方法,即在时域对信号进行离散变换,在频域进行谱分析的方法。它具有高分辨率的特点,而且在时、频两域都具有表征信号局部特征的能力。它在低频部分具有较高的频率分辨率和较低的时间分辨率,在高频部分具有较高的时间分辨率和较低的频率分辨率,很适合于探测正常信号中夹带的瞬态反常现象,所以被誉为分析信号的显微镜和望远镜。1.小波函数的定义小波(wavelet,即小区域的波,是一种特殊的长度有限、平均值为0的波形。小波函数的确切定义为:设Y(t□L2(R,若其傅立叶变换(FT:FourierTransform(ωΨ满足条件:2(Cdψψωωω+∞−∞=<∞∫(1则称(tΨ为一个基本小波或小波母函数。式(1为小波函数的可容许条件。将小波母函数(tΨ进行伸缩和平移,就可以得到函数,(atτψ:,(atτψ,;0aRaτ∈>(2式中,a为尺度因子,τ为平移因子,我们称,(atτψ为依赖于参数a、τ的小波基函数。由于尺度因子a和平移因子τ是连续变化的值,因此我们称,(atτψ