文章编号:1673-0291(2011)06-0049-08零相位数字滤波器在非平稳信号处理中的应用常 广,鄢素云,王 毅(北京交通大学电气工程学院,北京100044)摘 要:研究零相位数字滤波器在进行非平稳信号滤波时的特点.选用一种典型带通零相位数字滤波器,以非平稳调幅信号作为滤波器输入,进行仿真分析.将零相位数字滤波器与小波包分解重构和经验模态分解方法的滤波能力进行了比较.探讨了零相位数字滤波器在处理非平稳调幅信号时存在的过渡过程,及对滤波结果幅值和频率的影响.论述了滤波误差与滤波器参数、输入信号特性和信噪比等因素的关系.提出了分段零相位滤波器方法,改善了滤波器性能.最后,以一个实测的振荡信号对上述分析进行了验证.为在非平稳信号处理中,正确使用零相位数字滤波器提供了参考.关键词:数字滤波器;零相位;调幅信号;非平稳信号;分段零相位数字滤波中图分类号:TM 930 文献标志码:AApplication of zero -phase digital filter onnon -stationary signal processingCHAN G Guang ,YAN Suyun ,WAN G Yi(School of Electrical Eng ineering,Beijing Jiaotong U niversity,Beijing 100044,China)Abstract:T he characteristics of zero -phase dig ital filter w hen being utilized to process the non -station -ary signals are studied.Ty pical band -pass zero -phase digital filters are simulated.And non -stationary amplitude -modulation signals are selected to be input of the simulation.Wavelet packet decomposition and reconstruction,empirical mode decomposition and the zero -phase dig ital filter are applied to com -pare their band -pass filter capabilities.The simulation demonstrates the transition process in non -sta -tionary signal filtering.And it clarifies amplitude characteristics,and frequency characteristics existing in filtering the non -stationary am plitude -modulation signal in detail.This article also discusses the rela -tionship between error and filter parameters,characteristics of input sig nal and signal to noise ratio ofinput sig nal.A segment zero -phase dig ital filter m ethod is proposed in this paper to enhance the perfor -mance of the normal zero -phase dig ital filter.The segment zero -phase digital filter is em ployed in ex -tracting the main component from a real oscillation signal to verify the validity of the new zero -phase filtering method.The study prov ides support for proper usage of zero -phase digital filter applied on non -stationary signal processing.Key words:dig ital filter;zero -phase;amplitude -modulation signal;non -stationary sig nal;segment ze -ro -phase dig ital filter收稿日期:2011-05-10基金项目:国家自然科学基金资助项目(60674013); 十一五 国家科技支撑计划(2009G09-1-5)作者简介:常广(1978 ),男,湖南汨罗人,博士生,主要研究方向为智能电器、机电系统状态检修.email:guang -chang@.王毅(1958 ),男,辽宁沈阳人,教授,博士,博士生导师.email:yw ang5@.数字滤波是数字信号处理的常用手段.普通的数字滤波器在滤波时存在一定的相移.为解决该问题,零相位数字滤波器被引入到数字信号处理领域中.依据正向序列和翻转序列所处位置的不同,主要第35卷第6期2011年12月北 京 交 通 大 学 学 报JOU RN AL O F BEIJIN G JIAOT O NG U N IV ERSI T YV ol.35N o.6Dec.2011有FRR(Forw ard-filter,Reverse-filter,and Reverse-output)和RRF(Reverse-filter,Reverse-filter,and Forw ard-output)两种常用的基本零相位数字滤波器算法[1-2].前人对此开展了大量的研究工作.通过选择合适的滤波器初始状态[1,3]、信号延拓[3]、FRR与RRF相结合[4]等方法,较好地改善了零相位数字滤波器的边界问题,在平稳信号滤波中获得了很好的效果.在对非平稳信号进行滤波时,与小波分解重构和经验模态分解等无相移的非平稳信号处理方法相比,零相位数字滤波器具有截止频率明确、计算量小的优点.尤其在计算能力较弱的嵌入式应用系统中,零相位数字滤波器的这一优点更加突出.调幅信号、调频信号是两类典型的非平稳信号,很多工程实际信号可以借助这两种信号进行模拟.例如,电力系统的次同步振荡、电压暂态振荡、开关设备动作时的振动、发电机组轴系的扭振等就可用调幅信号模拟;旋转机械调速过程中的转速则多为调频信号.文献[5]利用线性调频信号对零相位数字滤波器进行了分析.因此,本文作者将主要针对调幅信号来讨论零相位数字滤波器在非平稳信号处理中的特点.利用调幅信号对一种典型的带通零相位数字滤波器进行了仿真分析;对小波包分解重构、经验模态分解及零相位数字滤波器方法的滤波效果进行了比较;探讨了在调幅信号处理中存在的过渡过程,及其对零相位数字滤波的幅值特性以及频率特性的影响;论述了影响滤波误差的因素;提出了利用分段零相位滤波的方法改善零相位数字滤波器过渡过程的方法;用实测的次同步振荡信号对两种方法进行了验证.1 基本原理对于一维时间序列,零相位数字滤波器的基本原理是,利用正向时间序列和翻转时间序列通过滤波器时的相移相互抵消,从而实现滤波结果的零相移.利用傅里叶变换的时域位移性质和时间翻转性质,即可从理论上论证其正确性[1-2].以FRR方法为例,将长度为N+1的单边时间序列x[n]通过冲激响应为h[n]的滤波器,得到一次滤波序列y1[n],再将一次滤波序列进行时序反转、平移得到数据序列y2[n],y2[n]通过滤波器h[n]得到二次滤波序列y3[n],最后将二次滤波序列进行时序反转、平移得到最终滤波输出y[n].从而,有y1[n]=x[n]*h[n](1)y2[n]=y1[N-n](2)y3[n]=y2[n]*h[n](3)y[n]=y3[N-n](4)设x[n] X( )、y1[n] Y1( )、y2[n] Y2( )、y3[n] Y3( )、y[n] Y( )、h[n] H( )构成傅里叶变换对,且考虑到对实序列h[n]有H(- )=H*( ).依式(1~4),FRR滤波输出的频域表达如下式推导Y( )=e-j N Y3(- )=e-j N Y2(- )H(- )=e-j N e j N Y1( )H(- )=Y1( )H(- )=X( )H( )H(- )=X( )|H( )|2(5)因而,FRR方法实现了零相移滤波,其传递函数频域表示如下H FRR( )=H( )H(- )(6)同理可得,RRF方法实现零相移滤波时,其传递函数频域表示如下H RRF( )=H(- )H( )(7) 选定冲激响应序列为h[n]的普通数字滤波器,实际使用时可以通过两个分立的滤波过程来实现零相位滤波,也可以直接构造零相位数字滤波器,其冲激响应序列分别为h FRR[n]和h RRF[n].h FRR[n]=h[n]*h[-n](8)h RRF[n]=h[-n]*h[n](9)其中h[-n]为h[n]的时间翻转序列.h[n]的选择可有FIR滤波器和I IR滤波器两种类型.在零相位滤波器的应用中,IIR滤波器在相频特性上的劣势可以忽略,且在相同阶数下IIR滤波器比FIR滤波器具有更好的截止特性,因而h[n]一般选用IIR数字滤波器.2 3种实现零相位处理方法的比较为比较小波包分解重构、经验模态分解和零相位数字滤波器3种方法的差异.选取含有3个衰减振荡分量的输入信号f(t),信号采样率为1kH z.设定3个分量的振荡频率分别为10H z、23H z、39Hz,其计算公式如下f1(t)=e-3(t-1)u(t-1)sin(20 (t-1))(10)f2(t)=e-2(t-1)u(t-1)sin(46 (t-1))(11)f3(t)=e-(t-1)u(t-1)sin(78 (t-1))(12)输入信号f(t)的计算公式如下f(t)=f1(t)+f2(t)+f3(t)(13)50北 京 交 通 大 学 学 报 第35卷1)小波包分解重构.利用小波包可将信号分解到不同节点,各节点对应了不同频带的子波,对目标节点系数进行重构即可得到所需的滤波输出.考虑到重构信号的不失真和光滑,采用了bior 、sym 和db 小波基对多分量信号进行处理,最终选定db8小波基对f (t)进行5层分解重构,由于信号采样率为1kHz,小波包5层分解的各子波频带宽度为15 625H z,使10H z 、23Hz 、39H z 3个频率成分位于[5,0]、[5,1]、[5,3]节点对应的频带内[6],滤波结果如图1所示(图中虚线为滤波输出,实线为输入分量).小波包分解重构的滤波输出,突变点附近存在过渡过程,其余部分信号不光滑,在幅值和相位上都存在抖动,其中23H z 和39Hz 分量抖动比较严重.图1 小波包分解重构F ig.1 Wavelet packet decomposit ion and reconstructio n2)经验模态分解.经验模态分解方法是自适应的将信号分解为不同时间尺度的本征模态函数(IMF)分量和趋势项,本征模态函数是一类零均值的调幅调频信号.对含随机宽带噪声的信号,经验模态分解是类小波的二进滤波器组,第一个本征模态函数的频带最高[7].由于m (t)由3个调幅正弦信号组成,经验模态分解得到的前3个本征模态函数分量即对应其滤波输出.其中,第一个本征模态函数对应39Hz 分量、第二个本征模态函数对应23H z 分量、第三个本征模态函数对应10Hz 分量,滤波结果如图2所示(图中虚线为滤波输出,实线为输入分量).图2 经验模态分解F ig.2 Empir ical mode decomposit ion经验模态分解方法的滤波输出,突变点附近存在过渡过程,但其余部分信号很光滑,幅值和相位的稳定性强于小波包分解重构方法,但是仍然存在抖动.我们又采用标准差为0 005、聚合次数为50次的聚合经验模态分解对f (t)进行滤波,发现能改善过渡过程.但是,同时造成了39Hz 分量混叠到23H z 分量中,为减少混叠效应需要增加聚合次数,极大的增加了计算量.3)零相位数字滤波器.数字滤波器在对有限长信号进行滤波时,相当于给信号加了一个矩形窗,在滤波输出信号上造成衰减振荡.零相位数字滤波器进行了两次滤波,在输出信号的首尾两端形成失真[3].FRR 方法第二次滤波的首端是输入信号的尾端,滤波输出的首端失真较小;RRF 方法第二次滤波的首端是输入信号的首端,滤波输出的尾端失真较小[4].对于本文考虑的调幅正弦信号,通过在输入信号首尾两端进行信号延拓的方法,来消除端点效应.依据对端点效应误差的要求,延拓的数据长度应比拟于过渡过程时间,本文取1s 长度的数据点.消除端点效应的FRR 方法和RRF,在使用中的差异可忽略不计.本文的后续讨论中,零相位数字滤波器均采用FRR 方法实现,且采用了首尾端的端点延拓,不再特别说明.综合考虑滤波器阶跃响应的过渡过程时间和振荡次数等因素,选定用4阶巴特沃斯带通滤波器(直接型系数见表1)构造零相位数字滤波器,中心频率分别为10H z 、23Hz 、39H z,带宽为3Hz,记为h bp10、h bp23、h bp39.滤波结果见图3(图中虚线为滤波输出,实线为输入分量).表1 带通滤波器系数T ab.1 Band -pass filters coefficients10Hz23Hz39Hz反向系数前向系数反向系数前向系数反向系数前向系数1 00000-3 96551 5 90483-3 913000 97369 0 00009 0 00000-0 00018 0 000000 00009 1 00000-3 93192 5 83851-3 879860 97369 0 00009 0 00000-0 00018 0 000000 00009 1 00000-3 85464 5 68806-3 803600 97369 0 00009 0 00000-0 00018 0 000000 00009滤波输出很光滑,除了突变点附近,幅值和相位与原始信号重合,效果好于小波包分解重构和经验模态分解方法.不考虑信号的突变,则零相位数字滤波器实现了较严格的零相移滤波.3种方法都存在过渡过程.以10H z 分量的误差为例来考察过渡过程,将3种方法所得结果前3s的误差示于图4中.从上至下依次为小波包分解重构、经验模态分解、零相位滤波器误差,实线为原始51第6期 常 广等:零相位数字滤波器在非平稳信号处理中的应用图3 零相位数字滤波器Fig.3 Zero -phase digital filter误差,虚线为幅值放大100倍后的误差.其中,经验模态分解误差极值最大,零相位滤波器与小波包分解重构的误差极值较小;零相位滤波与经验模态分解误差稳定时间较短,小波包分解重构的误差稳定时间较长;但零相位滤波器反时序的特性造成突变点前有信号的提前输出.图4 10Hz 分量的误差F ig.4 Error of 10Hz component下面我们将对零相位滤波器的过渡过程详细探讨,并研究其幅值和频率特性,从而对其过渡过程性能进行改善.3 过渡过程3 1 阶跃信号零相位数字滤波器解决了滤波中存在的相移问题,但是在对非平稳信号进行滤波时,仍需注意滤波中可能存在的过渡过程问题.尤其是在信号发生突变的位置影响较大,这是因为突变点的存在类似于将信号在该点截断,引入了新的端点.考察一个在第2秒时刻幅值由0跃变为1在第8秒再跳变为0的10Hz 正弦信号x (t),其计算公式如下x (t )=sin (20 (t -2))[u(t -2)-u (t -8)](14)其中u (t)为单位阶跃信号.将x (t )通过h bp10进行滤波,滤波结果如图5所示,存在过渡过程.为了更清楚地表达过渡过程,将其幅值响应曲线示于图6.输入信号的突变在滤波输出的对应位置上存在过渡过程,稳定到 2%的时间约0 48s,最大超调量为3 4%;进入平稳阶段后,滤波输出与输入信号能很好的吻合,误差均值在10-5量级,而误差的方差在10-7量级,平稳阶段时间越长,误差越小.对这类幅值发生阶跃的正弦信号,可以采用滤波输出信号最大幅值的中值时刻来近似定位突变点,从而去除突变点前的虚假信号.采用该方法估计x (t )的突变点的估计误差只有0 005s,实际应用中受波形形态及噪声的影响,误差变大.图5 信号突变点处的过渡过程F ig.5 T ransition pro cess at abrupt disco ntinuouspoints图6 信号突变点处的幅值响应F ig.6 Amplitude r esponses at abrupt discontinuous points事实上,不同特点的信号通过同一零相位数字滤波器,滤波信号幅值和频率表现出的特性也会不同.为进行具体的仿真分析,选取衰减振荡、发散振荡、发散振荡经采取措施后又衰减振荡这3种情况进行研究,利用幅值衰减、幅值发散、幅值先发散后衰减这3种不同调幅正弦信号进行模拟.3 2 3类调幅信号先构造调幅正弦信号的幅值包络a 1(t )、a 2(t )、a 3(t),其中a 1(t)为衰减幅值、a 2(t)为发散幅值、a 3(t )为先发散后衰减幅值,其计算公式如下a 1(t)=e -(t-1)u (t -1)(15)a 2(t)=e (t-9)[1-u(t -9)](16)a 3(t)=e (t-5)[1-u(t -5)]+e -(t -5)u (t -5)(17)对应的调幅正弦信号分别为m 1(t)、m 2(t )、m 3(t ),其计算公式如下m 1(t)=a 1(t )sin (20 t)(18)m 2(t)=a 2(t )sin (20 t)(19)52北 京 交 通 大 学 学 报 第35卷m 3(t )=a 3(t)sin (20 t )(20)采用一个高斯白噪声信号w n (t )作为噪声干扰信号,w n (t )的统计分布为N (0,0 01).对应的含噪声的输入信号分别为s 1(t)、s 2(t)、s 3(t),其表达式分别为下列公式所示s 1(t)=m 1(t )+w n (t)(21)s 2(t)=m 2(t )+w n (t)(22)s 3(t)=m 3(t )+w n (t)(23)s 1(t)、s 2(t)、s 3(t )的波形如图7所示,经过零相位带通数字滤波器h bp10滤波后,滤波输出信号如图8所示.可明显看出,这3种信号的滤波输出的最大幅值小于的x (t)的滤波输出的最大幅值;过渡过程时间也小于x (t)的过渡过程时间;其中s 3(t)的滤波结果最好,误差均值和方差都在10-5量级.可见输入信号的形态特点对滤波结果有影响.图7 含噪声的输入信号Fig.7 Noisy inputsignals图8 滤波输出Fig.8 Filtered signals4 幅值与频率特性4 1 幅值与频率的获取滤波输出的幅值和频率可利用希尔伯特变换求得.其基本原理如下.对信号 (t),其希尔伯特变换由下式定义(t )=1+ - ( )t -d (24)对应的幅值包络a (t)可按下式计算a(t)= 2(t)+ 2(t )(25) 对应的相位 (t)可按下式计算(t)=arctan (t )(t )(26)对应的瞬时频率定义为f (t)=12 d (t )d (t )(27)下面利用希尔伯特变换获取滤波结果的幅值包络和频率,对零相位数字滤波器的幅值响应和频率响应特点进行详细分析.4 2 幅值特性s 1(t)、s 2(t)、s 3(t)经过零相位带通数字滤波器h bp10滤波后,理想的幅值应为a 1(t )、a 2(t )、a 3(t ).实际按式(25)计算出的滤波输出信号幅值分别如图9~图11所示.可明显看出,s 1(t)、s 2(t)、s 3(t)的滤波输出的最大幅值分别只有0 81、0 80和0 93,明显小于的x (t)的滤波输出的最大幅值1 034;s 3(t)的最大值没有明显延迟,s 1(t )、s 2(t)过渡过程时间为0 22s 、0 22s,也明显小于对x (t)滤波时的过渡过程时间0 48s;信号幅值较小的部分噪声影响比较大.由于不像x (t)有稳定的信号阶段可以使滤波输出完成过渡达到稳定,上述3类信号衰减系数越大幅值损失越大.将不同衰减系数的10H z 分量滤波输出的幅值示于图12中,其结果印证了该结论.图9 含噪声的衰减信号幅值Fig.9 Amplitudes of the noisy decaying signal图10含噪声的发散信号幅值Fig.10 Amplitudes of the noisy grow ing signal图11 含噪声的先发散后衰减信号幅值Fig.11 Amplitudes of the noisy grow ing -decay ing signal53第6期 常 广等:零相位数字滤波器在非平稳信号处理中的应用图12 不同衰减系数下的幅值Fig.12 Comparison of different decaying expo nents对于零相位带通数字滤波器而言,与普通带通数字滤波器相似,其幅值响应具有如下特点:1)滤波器类型一定时,较宽的频带可以获得较大的幅值和较小的时延;提高滤波器阶数可加速阻频带衰减,减少临近频带的干扰,但也会带来通带幅值的振荡.2)输入信号幅值包络变化时,滤波输出的时延从幅值包络变化率大的部分偏向变化率小的部分;幅值包络变化越平缓,滤波后的时延也越小;滤波器阶数大于等于4时,信号幅值包络越接近阶跃,超调量越大3)因为白噪声无法完全滤除,幅值较小的部分噪声影响比较明显.4 3 频率特性按式(27)式定义对m 1(t )、m 2(t)、m 3(t )的滤波输出信号进行计算,得到的瞬时频率分别示于图13~图15中.其中,首尾两端的振荡是由于希尔伯特变换的边界效应造成的.另外需要注意的地方是,在信号突变点附近(约为一个过渡过程时间),滤波输出存在较大的频率突变,说明滤波输出信号在该位置存在畸变,在图13~图14中分别位于0 755s 、9 246s 处.取中间频率稳定的3~8s 时间段进行分析,可得频率均值为10 000H z,方差小于10-5量级.其中,m 3(t)在中间突变点处的频率处于9 99Hz 与10 01H z 之间,明显好于m 1(t )、m 2(t)在突变点处的滤波结果.图13 衰减信号的频率F ig.13 Frequency o f the filtered decay ing signal为考察噪声对频率响应的影响,取s 1(t )信号为例.图16给出了s 1(t)经滤波后的瞬时频率曲线,可见滤波结果变差,信号幅值较小的部位上噪声造成相位的波动,使波形相对发生了较大的畸变.仍取3~8s 时间段进行分析,可得频率均值为9 962Hz,方差在10-1量级.图14 发散信号的频率Fig.14 F requency of the filt ered g rowingsignal图15 先发散后衰减信号的频率Fig.15 F requency of the filter ed grow ing -decayingsignal图16 含噪声的衰减信号的频率F ig.16 Frequency of the noisy filtered decaying sig nal零相位带通数字滤波器的频率响应具有如下特点:1)中心频率和频带宽度变化对稳定的频率段影响很小,主要对信号突变点附近的频率有影响 2)当输入信号的幅值包络变化时,滤波输出的频率畸变点位置与幅值时延的方向相反,位于突变点偏幅值包络变化率较大的一侧3)噪声的存在,造成信号的畸变,滤波输出频率存在波动,幅值较小的部分噪声影响比较明显.5 分段零相位滤波由上述仿真结果可知,m 3(t )信号的滤波结果明显好于m 1(t)与m 2(t).m 3(t)与m 1(t)、m 2(t)的不同之处在于,其突变点两边的幅值是对称的.基于此,作者提出分段零相位滤波方法,用以改善零相位滤波器在处理存在突变的调幅信号时的过渡过程性能,其具体步骤如下:1)对信号进行零相位滤波54北 京 交 通 大 学 学 报 第35卷2)在滤波幅值包络的中值点附近,搜索滤波输出中绝对值最小的点,进行突变点的近似定位 3)将输入信号按定位结果分段,对分段后的信号进行端点延拓,再分别进行零相位滤波4)将分段滤波结果组合,形成滤波输出.零相位滤波器采用h bp10,利用零相位滤波和分段零相位滤波方法分别对s 1(t )进行滤波,滤波结果示于图17中.从图中可以明显观察到,分段零相位滤波的结果不但比原始的零相位滤波结果过渡过程缩短,而且也解决了反时序提前输出的问题.可见,分段零相位滤波能缩短过渡过程时间,消除反时序提前输出,极大改善零相位滤波器的滤波效果.其缺点是,增加了滤波次数和端点延拓次数,增加了计算量.图17 分段零相位数字滤波Fig.17 Comparisons of ZPDF and Segment ZPDF6 次同步振荡分量的提取为从次同步振荡中分离出主要的振荡分量,将该方法应用到次同步振荡信号的滤波中.在某次试验中测得振荡信号,采样率为1kHz,数据长度为95s,如图18所示.图18 次同步振荡信号Fig.18 Subsy nchr onous oscillation signal通过时频分析知该振荡信号中主要包含有13Hz 、50H z 、25Hz 等几个较大的频率成分,其中在次同步频域内13Hz 振幅最大.以13Hz 分量为例,仍采用4阶巴特沃斯带通滤波器(直接型的中心频率为13H z,带宽3Hz),构造零相位数字滤波器、小波包分解重构和经验模态分解方法分别从振荡信号中分离出该分量.如图19曲线所示,从上到下依次为零相位数字滤波器、小波包分解重构和经验模态分解的分离结果.实际被激起的次同步振荡信号为多段衰减的正弦信号,其幅值包络应比较光滑,零相位数字滤波器结果最好,小波包分解重构所得主分量的幅值剧烈抖动,而经验模态分解结果则几乎完全失真.为进一步比较零相位滤波器和分段零相位滤波器的滤波效果,取其中15 8~16 8s 时间段作局部放大示于图20中.可见主分量很好地滤除了振荡信号中的其他频率成分,且相位相对于原振荡信号几乎没有移动.分段零相位数字滤波器滤波结果的过渡过程优于普通的零相位数字滤波器,且减小了反时序的提前输出.但是,由于没有精确定位到原信号的突变点上,在突变点附近存在较小的畸变;如计算能力允许,应采用小波分解的模极大值等方法进行突变点的精确定位.图19 不同方法获得的主分量Fig.19 M ain component obtained v ia differ ent methods图20 实际信号的分段零相位数字滤波F ig.20 Comparisons of ZP DF and seg ment ZPDFapplied on subsynchronous oscillat ion signal7 结论零相位数字滤波器能有效解决滤波中存在的相移问题,但是在处理非平稳信号时可能存在过渡过程.针对调幅信号通过仿真分析和实例验证,阐明了零相位数字滤波的误差与滤波器参数、输入信号的幅值形态和信噪比等因素的关系:1)滤波误差主要受过渡过程和信噪比的影响,过渡过程与滤波器参数、信号幅值形态有关.2)对一定类型的零相位带通滤波器,滤波器阶数和频带宽度的设置会改变过渡过程响应,从而影响滤波误差.3)输入信号的幅值包络不同,滤波误差存在差55第6期 常 广等:零相位数字滤波器在非平稳信号处理中的应用异;幅值包络变化越剧烈,误差越大.4)提高信噪比能减小滤波误差.为解决零相位数字滤波器在处理存在突变的调幅振荡信号时的过渡过程问题,本文提出了分段零相位滤波的方法.通过与小波包分解重构、经验模态分解方法对实际信号的处理进行对比,验证了该方法的有效性.实际使用中可根据信号特点选取合适的滤波方法.参考文献(References):[1]Gustafsson F.Determining the initial states in forw ar d-backward filter ing[J].I EEE T rans on Signal Processing, 1996,44(4):988-992.[2]纪跃波,秦树人,汤宝平.零相位数字滤波器[J].重庆大学学报:自然科学版,2000,23(6):4-7.JI Yuebo,Q IN Shuren,T A NG Baoping.Digital filtering wit h zero phase erro r[J].Jour nal of Cho ngqing U nivers-i ty:N atural Science Edit ion,2000,23(6):4-7.(in Ch-i nese)[3]陈淑珍,杨涛.零相移滤波器的改进及实现方法[J].武汉大学学报:理学版,2001,47(3):373-376.CHEN Shuzhen,YAN G T ao.Improvement and realiza-tion of the zer o-phase filter[J].Journal of Wuhan U niv er-sity:Natural Science Editio n,2001,47(3):373-376.(inChinese)[4]尚秀辉,郭爱煌,李广宇.基于零相位数字滤波器的边界问题的分析[J].电子测量技术,2010,33(4):25-27.SHA NG X iuhui,GU O Aihuang,LI Guang yu.Analysis of the boundar y problems based on zer o-phase digital filter [J].Electro nic Measurement T echnology,2010,33(4): 25-27.(in Chinese)[5]李之雄,郭瑜,郑华文.零相位滤波器在非平稳信号分析中的应用研究[J].昆明理工大学学报:理工版,2007,32(5):18-22.L I Zhix iong,GU O Y u,ZHENG Huawen.Study on the application o f zero-phase filter to non-stationary signal ana-l ysis[J].Journal of Kunming U niversity of Science and T echnology:Science and T echnolo gy,2007,32(5):18-22.(in Chinese)[6]曾宪伟,赵卫明,盛菊琴.小波包分解树结点与信号子空间频带的对应关系及其应用[J].地震学报,2008,30(1):90-96.ZENG Xianw ei,ZHAO Weiming,SH EN G Juqin.Corre-spo nding r elationships between nodes of decomposition tree of w avelet packet and fr equency bands o f signal subspace [J].Acta Seismolog ica Sinica,2008,30(1):90-96.(in Chinese)[7]Flandr in P,Rilling G,Goncalv es P.Empirical mode de-composition as a filter bank[J].I EEE Sig nal P rocessing L ett,2004,11(2):112-114.56北 京 交 通 大 学 学 报 第35卷。