乔建华+张雪英
摘 要: 谱分析是数字信号处理中的一个重要问题,初学者普遍对连续信号谱分析理解不深,尤其是在误差分析时缺乏统一示例,更容易产生困惑。介绍了用离散傅里叶变换(DFT)对连续信号进行谱分析的过程,并详细说明了误差产生的原因和减小误差的方法。而且通过对模拟信号谱分析的实例全面说明了各项误差的影响及解决方案,并应用Matlab直观地进行了分析和对比验证。
关键词: 谱分析; 误差分析; DFT; Matlab
中图分类号: TN911.7?34 文献标识码: A 文章编号: 1004?373X(2014)13?0053?04
Error of continuous signal spectrum analysis based on DFT
QIAO Jian?hua1,2, ZHANG Xue?ying2
(1. College of Electronic & Information Engineering, Taiyuan University of Science & Technology, Taiyuan 030024, China;
2. College of Information Engineering, Taiyuan University of Technology, Taiyuan 030024, China)
Abstract: Spectrum analysis is one of the important missions in digital signal processing. The understanding of most beginners is not deep to continuous signal spectrum analysis, especially lacking of a uniform example in the error analysis, so it is easier to be confused. The process of continuous time signal spectrum analysis using the discrete Fourier transform (DFT) is introduced in this paper. The reasons of generating the errors and the methods of minimizing the errors are illustrated in detail. The influence and solutions of various errors are elaborated through a concrete example of analog signal spectrum analysis. It is verified intuitively through analysis and contrast with Matlab software.
Keywords: spectrum analysis; analysis error; DFT; Matlab
0 引 言
连续信号的谱分析是数字信号处理的一个重要内容。对信号进行谱分析,关键是得到信号的傅里叶变换。由于离散傅里叶变换(DFT)有快速算法FFT,因此经常用DFT(FFT)对信号进行谱分析。但是根据傅里叶变换理论,若信号持续时间有限长,则其频谱无限宽;若信号的频谱有限宽,则其持续时间无限长[1]。因此,用有限长序列的有限点的DFT对连续信号进行谱分析必然是近似的,要明确误差产生的原因,并能够采取合理措施来降低误差。文献[2?5]也以不同的侧重点介绍了DFT在谱分析中的仿真实现,物理实现,具体应用以及针对各种不同信号的谱分析。而对于实际中广泛存在的连续信号,初学者普遍对其应用DFT的谱分析的原理和过程理解不深,尤其对误差产生原因不明确,对减小误差方案不清晰,在具体的应用中不知如何整合这些方案,这都是对整个谱分析过程研究不透彻的缘故。本文从具体的模拟信号谱分析的实例入手,来充分说明这一点,并通过Matlab演示来加深对整个谱分析过程中误差的认识。
1 用DFT对连续信号谱分析
假设模拟信号为[xa(t)],根据Nyquist采样定理,使用模拟信号带宽2倍以上的采样率才能保证频谱不混叠,但实际上理想的带限信号时域上总是无限长的,因此在工程应用时对于带宽外能量占总能量比例很小的信号,可以近似看成是带限信号,认为抽样信号是不失真的。
例如设[xa(t)=e-1 000|t|,]其傅里叶变换为:
[Xa(jΩ)=FT[xa(t)]=-∞∞xa(t)e-jΩtdt=-∞∞e-1 000|t|e-jΩtdt=2 0001 0002+Ω2] (1)
可以看到这个信号不是一个理想带限的信号。但是随着[Ω]增加,[Xa(jΩ)]会逐渐减小,可以规定一个带宽[Ωc,]使得[Xa(jΩc)?1,]认为当采样频率[fs>Ωcπ]时可以无失真地恢复原模拟信号,如本例中取[fc=Ωc2π=2]kHz,此即造成了误差的产生。[xa(t)]及其幅频特性曲线如图1(a),(b)所示。
对[xa(t)]以采样频率[fs=5 ]kHz采样,即采样间隔[T=][0.2 ]ms,得采样信号[xa(nT)=xa(t)t=nT=x(n),]在式(1)中,由于[t→nT,dt→T,-∞∞dt→n=-∞∞T,]因此得采样信号的频谱:
[Xa(jΩ)≈T?n=-∞∞xa(nT)e-jΩnT] (2)
假设将采样序列[xa(nT)]截取成长度为[M]的有限长序列,此截断即为误差产生之二。若频率用[f]表示,则式(2)写为:
[Xa(jf)≈T?n=0M-1xa(nT)e-j2πfnT] (3)
显然[Xa(jf)]仍是[f]的连续周期函数,[xa(nT)]和[Xa(jf)]如图1(c),(d)所示。
为了进行数值运算,频域上也要离散化,因此在频域的一个周期[0, fs]内等间隔采样[N]点,如图1(f)所示。由于间隔采样有可能漏掉重要的频域信息,成为误差产生之三。设频域采样间隔为[F,]则[F=fsN=1(NT),]将[f=kF]代入式(3)中可得对采样信号频谱[Xa(j f)]的采样:
[Xa(jkF)≈T?n=0M-1xa(nT)e-j2πNkn] (4)
根据频域采样定理[6],如果序列[x(n)]的长度为[M,]只有当频域采样点数[N≥M]时,才可由频域采样[X(k)]恢复原序列[x(n),]否则会产生时域混叠。因此为简便计算,将序列也以[N]点长度来截取,则式(4) 的求和上限为[N-1,]并令:[Xa(k)=Xa(jkF),x(n)=xa(nT)],则:
[Xa(k)≈T?n=0N-1x(n)e-j2πNkn=T?DFT[x(n)]N, k=0,1,…,N-1] (5)
可见,对近似持续时间有限的带限连续信号,其频谱特性可以通过对连续信号采样[N]点并进行DFT再乘以[T]的近似方法得到。因此用DFT对模拟信号谱分析,产生误差是必然的。其具体分析过程是首先确定模拟信号的最高频率[fc]和频率分辨率[F,]再选取参数,包括采样频率[fs,]采样点数[N,]记录时间[Tp]等。其取值分别为:
[fs≥2fcN=fsF≥2fcFTp=NT=Nfs≥1F] (6)
然后根据计算的参数,通过采样得到信号序列[x(n)],即可根据式(5)进行运算分析。
图1 DFT对连续信号谱分析过程
下面具体分析误差的影响及其减小误差的措施。
2 误差分析
2.1 频谱混叠
在对模拟信号[xa(t)]进行采样时,必须满足采样定理,即采样频率[fs≥2fc,]而通过上面分析时域有限的信号不可能是锐截止的,并且信号中不可避免地有一些高频杂散信号,因此在采样之前,一般都要对模拟信号进行滤波,滤除高频杂散信号。如本例中,[fc]取2 kHz,当[fs]=5 kHz时,采样信号的频谱混叠很小,见图1(d);但当选[fs]=1 kHz时,其频谱如图1(h)所示,混叠严重。以折叠频率[fs2]处的频谱幅度来比较[T=]0.2 ms和[T=]1 ms时的混叠程度,此处[k=256,]键入Matlab语句,得到:
abs(X1(k))/max(abs(X1))
ans=0.009 4
abs(X2(k))/max(abs(X2))
ans=0.212 6
而在[fc]处,模拟信号的幅度比值为:
abs(Xa(2k))/max(abs(Xa))
ans=0.006 6
可以看出,随着采样频率的减小,混叠现象加大。
因此,要减小混叠,必须满足乃奎斯特采样定理,并且在采样前进行预滤波,滤除高于折叠频率[fs2]的频率成分,一般取采样频率[6][fs≥(3~5)fc。]
2.2 截断效应
对无限长的模拟信号,用DFT进行谱分析时,必须先进行截断,通过采样才能得到有限点的序列,这样必然产生误差。截断可以理解为加窗,即:
[y(n)=x(n)w(n)] (7)
式中:[x(n)]为模拟信号经采样得到的时域离散信号;[w(n)]为窗函数序列。根据频域卷积定理,加窗后的信号频谱为:
[Y(ejω)=12πX(ejω)*W(ejω)] (8)
显然与原序列的频谱是不同的。
文献[7]对常用的五种窗函数的时频域波形有详细的对比,是对窗函数本身的性能做了说明。下面仅以矩形窗和hamming窗为例来从加窗对信号频谱的影响方面来进行分析。图2所示为矩形窗和hamming窗的时域和频域波形。由时域波形可见,矩形窗是对原信号的原样截取,hamming窗是在截取的同时对信号做了加权。由频域波形可见,矩形窗的主瓣窄,其宽度[8]为[4πN,]但旁瓣幅度也高;hamming窗是以主瓣展宽为代价换取了旁瓣幅度的减小,其主瓣宽度[8]为[8πN,]而主瓣宽度影响的就是频率分辨率。可见,通过加窗对信号截断使得对模拟信号的谱分析产生误差,表现为原本离散的谱线展宽(高频泄漏)和谱间干扰等现象,即所谓截断效应。从而降低了谱的分辨率,使频率相近的两个信号不易分清。
图2 矩形窗和hamming窗的时频曲线
因此,为了减小截断效应的影响,一种方法是增大截取长度,通过改变[N,]使窗函数的主瓣变窄,提高频率分辨率。另一种方法是改变截断窗函数的形状。
下面通过一个例子,来更清楚地展示泄漏和谱间干扰的影响,以及通过截取长度和加不同的窗函数对频谱的改善作用[9],设:
[x(n)=cos(2πf1n)+sin(2πf2n)+cos(2πf3n)]
其中[f1]=25 Hz,[f2]=50 Hz,[f3]=100 Hz,可见其最高频率为[f3,]频率分辨率至少为25 Hz,则根据式(6),截取长度[Tp≥1F=]0.04 s,取采样频率[fs]=400 Hz,则最小采样点[N=][Tp?][fs]=16点,图3所示为其加矩形窗和hamming窗以[N,][4N,][8N]点分别截断的频谱效果。由图可以清晰地看到,当截取长度太短的时候,泄漏和谱间干扰都非常严重,25 Hz和50 Hz的两条谱线已无法分辨,只有增加截取长度[N,]才使得主瓣变窄,提高了频率分辨率;而采用hamming窗又降低了旁瓣的影响,减小了谱间干扰。而且可以看出,当[N]一定时,旁瓣越小的窗函数,其主瓣就越宽。增加[N]使主瓣变窄,但旁瓣的相对幅度并不减小。矩形窗的主瓣较窄,但是旁瓣幅度也较高。可见,当截取长度一定时,频率分辨率和谱间干扰是相互抵消的,只能以一种优势来换取另一种性能的降低。
图3 矩形窗和hamming窗以不同长度截断效果
2.3 栅栏效应
栅栏效应是指用DFT对连续信号进行谱分析时,由于DFT的离散特性,使离散点之间的频谱无法得到,相当于透过栅栏观察频谱,只看到等间隔的离散点的频线,其他的频谱都被栅栏挡住了,故称之为栅栏效应。因此用DFT得到的离散谱线的包络只能是近似谱。
为了减小栅栏效应,可以多增加些栅栏的缝隙,即增大DFT变换点数,这一方面可以通过在原序列尾部补零来实现[10]。例如对下式信号进行分析:
[x(n)=sin(2πf1n)+sin(2πf2n)+sin(2πf3n)]
其中[f1]=20 Hz,[f2]=20.5 Hz,[f3]=40 Hz,由其最高频率为[f3],频率分辨率至少要达到0.5 Hz,根据式(6)选择参数,若采样频率[fs]取80 Hz,则采样点数[N≥][800.5=]160。若[N]取160点,再进行DFT,可得到图4波形[x1(n)]及幅度谱[X1(k),]由[X1(k)]曲线可见,无法分辨20 Hz和20.5 Hz的两条谱线,这就是由于栅栏效应造成的。如果在序列后补零加长到512点,如图4中[x2(n),]再看其幅度谱[X2(k),]频谱被细化,挨近的两条谱线已能有所区分。
图4 补零抑制栅栏效应
但是,如果对原序列采样128点,即使补零加长到512点,如图5所示,其频谱[X2(k)]也无法分辨出这两条谱线。这是因为截断效应已经使得频谱变模糊,即使再增加栅栏的缝隙,看到的也是模糊的频谱,补零并没有真正提高频率分辨率。这种情况下,只能通过增加截取长度,首先提高频谱分辨率。如对信号直接采样512点,观察图5的[X3(k)]幅频特性曲线,此时,这两条谱线可以清晰地分辨出来。
图5 增加截取长度抑制栅栏效应
总之,对模拟信号的DFT谱分析主要集中于两个方面:首先是对信号的采样,其中涉及的问题有频谱混叠和截断效应,因此要选择合适的采样频率和窗函数,足够的记录时间;然后是对得到的时域离散信号[x(n)]进行DFT运算,这时要注意栅栏效应,要计算足够点数的DFT,避免漏掉邻近的谱线。因此在选择参数的时候,在依据式(6)的基础上,实际应用中要适当加大取值。
3 结 语
本文针对数字信号处理中的一个普遍应用——采用DFT对连续信号谱分析中的误差问题,从基本原理上说明了利用DFT对连续信号谱分析的近似性,对造成误差的三种主要因素(频谱混叠、截断效应和栅栏效应)的产生做了深入分析,提出减小误差的相应措施,以及参数选择的基本原则,而且以详实细致的图例对相应问题作了说明。需要提出的是如果利用近代谱估计的方法,可以得到一些更好的效果,下一阶段可以在这方面继续进行研究。
参考文献
[1] OPPENHEIM A V, SCHAFER R W.数字信号处理[M].董士嘉,杨耀增,译.北京:科学出版社,1980.
[2] 谢海霞,孙志雄.连续非周期信号频谱分析及Matlab实现[J].现代电子技术,2013,36(11):53?56.
[3] 郭连平,田书林,王志刚.数字信号处理过程中信号截位误差抑制方法研究[J].信号处理,2013(5):544?546.
[4] 温万惠.在生理信号频谱分析中掌握离散傅里叶变换法[J].西南师范大学学报:自然科学版,2010(6):227?230.
[5] 林爱英,滕红丽,袁超,等.DFT在信号谱分析中的应用[J].安徽工业大学学报:自然科学版,2011(2):192?196.
[6] 高西全,丁玉美.数字信号处理[M].3版.西安:西安电子科技大学出版社,2008.
[7] 汪伟,谢皓臣,梁光明,等.加窗离散傅里叶变换性能分析和比对[J].现代电子技术,2012,35(3):115?118.
[8] ORFANIDIS S J.信号处理导论(影印版)[M].北京:清华大学出版社,1999.
[9] 陈怀琛.Matlab及在电子信息课程中的应用[M].3版.北京:电子工业出版社,2006.
[10] 薛年喜.Matlab在数字信号处理中的应用[M].2版.北京:清华大学出版社,2008.
图3 矩形窗和hamming窗以不同长度截断效果
2.3 栅栏效应
栅栏效应是指用DFT对连续信号进行谱分析时,由于DFT的离散特性,使离散点之间的频谱无法得到,相当于透过栅栏观察频谱,只看到等间隔的离散点的频线,其他的频谱都被栅栏挡住了,故称之为栅栏效应。因此用DFT得到的离散谱线的包络只能是近似谱。
为了减小栅栏效应,可以多增加些栅栏的缝隙,即增大DFT变换点数,这一方面可以通过在原序列尾部补零来实现[10]。例如对下式信号进行分析:
[x(n)=sin(2πf1n)+sin(2πf2n)+sin(2πf3n)]
其中[f1]=20 Hz,[f2]=20.5 Hz,[f3]=40 Hz,由其最高频率为[f3],频率分辨率至少要达到0.5 Hz,根据式(6)选择参数,若采样频率[fs]取80 Hz,则采样点数[N≥][800.5=]160。若[N]取160点,再进行DFT,可得到图4波形[x1(n)]及幅度谱[X1(k),]由[X1(k)]曲线可见,无法分辨20 Hz和20.5 Hz的两条谱线,这就是由于栅栏效应造成的。如果在序列后补零加长到512点,如图4中[x2(n),]再看其幅度谱[X2(k),]频谱被细化,挨近的两条谱线已能有所区分。
图4 补零抑制栅栏效应
但是,如果对原序列采样128点,即使补零加长到512点,如图5所示,其频谱[X2(k)]也无法分辨出这两条谱线。这是因为截断效应已经使得频谱变模糊,即使再增加栅栏的缝隙,看到的也是模糊的频谱,补零并没有真正提高频率分辨率。这种情况下,只能通过增加截取长度,首先提高频谱分辨率。如对信号直接采样512点,观察图5的[X3(k)]幅频特性曲线,此时,这两条谱线可以清晰地分辨出来。
图5 增加截取长度抑制栅栏效应
总之,对模拟信号的DFT谱分析主要集中于两个方面:首先是对信号的采样,其中涉及的问题有频谱混叠和截断效应,因此要选择合适的采样频率和窗函数,足够的记录时间;然后是对得到的时域离散信号[x(n)]进行DFT运算,这时要注意栅栏效应,要计算足够点数的DFT,避免漏掉邻近的谱线。因此在选择参数的时候,在依据式(6)的基础上,实际应用中要适当加大取值。
3 结 语
本文针对数字信号处理中的一个普遍应用——采用DFT对连续信号谱分析中的误差问题,从基本原理上说明了利用DFT对连续信号谱分析的近似性,对造成误差的三种主要因素(频谱混叠、截断效应和栅栏效应)的产生做了深入分析,提出减小误差的相应措施,以及参数选择的基本原则,而且以详实细致的图例对相应问题作了说明。需要提出的是如果利用近代谱估计的方法,可以得到一些更好的效果,下一阶段可以在这方面继续进行研究。
参考文献
[1] OPPENHEIM A V, SCHAFER R W.数字信号处理[M].董士嘉,杨耀增,译.北京:科学出版社,1980.
[2] 谢海霞,孙志雄.连续非周期信号频谱分析及Matlab实现[J].现代电子技术,2013,36(11):53?56.
[3] 郭连平,田书林,王志刚.数字信号处理过程中信号截位误差抑制方法研究[J].信号处理,2013(5):544?546.
[4] 温万惠.在生理信号频谱分析中掌握离散傅里叶变换法[J].西南师范大学学报:自然科学版,2010(6):227?230.
[5] 林爱英,滕红丽,袁超,等.DFT在信号谱分析中的应用[J].安徽工业大学学报:自然科学版,2011(2):192?196.
[6] 高西全,丁玉美.数字信号处理[M].3版.西安:西安电子科技大学出版社,2008.
[7] 汪伟,谢皓臣,梁光明,等.加窗离散傅里叶变换性能分析和比对[J].现代电子技术,2012,35(3):115?118.
[8] ORFANIDIS S J.信号处理导论(影印版)[M].北京:清华大学出版社,1999.
[9] 陈怀琛.Matlab及在电子信息课程中的应用[M].3版.北京:电子工业出版社,2006.
[10] 薛年喜.Matlab在数字信号处理中的应用[M].2版.北京:清华大学出版社,2008.
图3 矩形窗和hamming窗以不同长度截断效果
2.3 栅栏效应
栅栏效应是指用DFT对连续信号进行谱分析时,由于DFT的离散特性,使离散点之间的频谱无法得到,相当于透过栅栏观察频谱,只看到等间隔的离散点的频线,其他的频谱都被栅栏挡住了,故称之为栅栏效应。因此用DFT得到的离散谱线的包络只能是近似谱。
为了减小栅栏效应,可以多增加些栅栏的缝隙,即增大DFT变换点数,这一方面可以通过在原序列尾部补零来实现[10]。例如对下式信号进行分析:
[x(n)=sin(2πf1n)+sin(2πf2n)+sin(2πf3n)]
其中[f1]=20 Hz,[f2]=20.5 Hz,[f3]=40 Hz,由其最高频率为[f3],频率分辨率至少要达到0.5 Hz,根据式(6)选择参数,若采样频率[fs]取80 Hz,则采样点数[N≥][800.5=]160。若[N]取160点,再进行DFT,可得到图4波形[x1(n)]及幅度谱[X1(k),]由[X1(k)]曲线可见,无法分辨20 Hz和20.5 Hz的两条谱线,这就是由于栅栏效应造成的。如果在序列后补零加长到512点,如图4中[x2(n),]再看其幅度谱[X2(k),]频谱被细化,挨近的两条谱线已能有所区分。
图4 补零抑制栅栏效应
但是,如果对原序列采样128点,即使补零加长到512点,如图5所示,其频谱[X2(k)]也无法分辨出这两条谱线。这是因为截断效应已经使得频谱变模糊,即使再增加栅栏的缝隙,看到的也是模糊的频谱,补零并没有真正提高频率分辨率。这种情况下,只能通过增加截取长度,首先提高频谱分辨率。如对信号直接采样512点,观察图5的[X3(k)]幅频特性曲线,此时,这两条谱线可以清晰地分辨出来。
图5 增加截取长度抑制栅栏效应
总之,对模拟信号的DFT谱分析主要集中于两个方面:首先是对信号的采样,其中涉及的问题有频谱混叠和截断效应,因此要选择合适的采样频率和窗函数,足够的记录时间;然后是对得到的时域离散信号[x(n)]进行DFT运算,这时要注意栅栏效应,要计算足够点数的DFT,避免漏掉邻近的谱线。因此在选择参数的时候,在依据式(6)的基础上,实际应用中要适当加大取值。
3 结 语
本文针对数字信号处理中的一个普遍应用——采用DFT对连续信号谱分析中的误差问题,从基本原理上说明了利用DFT对连续信号谱分析的近似性,对造成误差的三种主要因素(频谱混叠、截断效应和栅栏效应)的产生做了深入分析,提出减小误差的相应措施,以及参数选择的基本原则,而且以详实细致的图例对相应问题作了说明。需要提出的是如果利用近代谱估计的方法,可以得到一些更好的效果,下一阶段可以在这方面继续进行研究。
参考文献
[1] OPPENHEIM A V, SCHAFER R W.数字信号处理[M].董士嘉,杨耀增,译.北京:科学出版社,1980.
[2] 谢海霞,孙志雄.连续非周期信号频谱分析及Matlab实现[J].现代电子技术,2013,36(11):53?56.
[3] 郭连平,田书林,王志刚.数字信号处理过程中信号截位误差抑制方法研究[J].信号处理,2013(5):544?546.
[4] 温万惠.在生理信号频谱分析中掌握离散傅里叶变换法[J].西南师范大学学报:自然科学版,2010(6):227?230.
[5] 林爱英,滕红丽,袁超,等.DFT在信号谱分析中的应用[J].安徽工业大学学报:自然科学版,2011(2):192?196.
[6] 高西全,丁玉美.数字信号处理[M].3版.西安:西安电子科技大学出版社,2008.
[7] 汪伟,谢皓臣,梁光明,等.加窗离散傅里叶变换性能分析和比对[J].现代电子技术,2012,35(3):115?118.
[8] ORFANIDIS S J.信号处理导论(影印版)[M].北京:清华大学出版社,1999.
[9] 陈怀琛.Matlab及在电子信息课程中的应用[M].3版.北京:电子工业出版社,2006.
[10] 薛年喜.Matlab在数字信号处理中的应用[M].2版.北京:清华大学出版社,2008.