连续时间信号频谱分析研究及MATLAB实现

井敏英 郭佳林 李大伟 段 涛

（陕西理工学院，陕西 汉中 723001）

【摘 要】信号的频谱分析是信号与系统分析的基础。文章分析了用数值计算的方法实现确知连续时间信号的频谱分析，即采用离散傅里叶变换的快速算法实现对连续信号的频谱估计，然后在MATLAB语言工具下结合正弦信号给出了频谱分析的结果。

【关键词】频谱分析；离散傅里叶变化；MATLAB 【中图分类号】TN911.6 【文献标识码】A 【文章编号】1008-1151(2011)04-0042-02

（一）引言 机不能直接加以处理。为了实现数值计算，还需要对F(Ω)进

在信号处理过程中，频域分析方法往往比时域分析方法行离散化处理。将频率段0~1/THz划分为N个计算频率更方便和有效。对于确知连续时间信号，其频域分析可以通点，这N个离散频率点以角频率表示为： 过连续时间傅里叶变换来进行，但是，这样计算出来的结果2πkkω=ω=,k=0,1,LN−1 （2） 0仍然是连续函数，计算机不能直接加以处理。为了实现数值N

计算，还需要对其进行离散化处理，即采用离散傅里叶变换

即得到离散频率点上的近似计算式：

（DFT）进行分析。DFT的快速算法的出现，使DFT在数字通N−1

F(ω)=Tf(n)exp(−jωn) （3） 信、图像处理、功率谱估计、系统分析与仿真、雷达信号处

n=0

理、光学、医学等各个领域都得到广泛应用。本文以正弦信

对比DFT计算式，显然有 号为例，介绍用DFT的快速算法即快速傅里叶变化（FFT）实

现确知连续时间信号的频谱分析，给出了MATLAB语言工具下 （4） F(Ω)2π≈T•DFT[f(n)]=T•F(k)Ω=k的分析程序。 NT

∑

（二）信号频谱分析的原理

对于时间连续信号f(t)，其频谱分析可以通过连续时间傅里叶变换（CTFT）来进行。连续时间傅里叶变化特别适合于对时间连续信号的理论分析，但是，由于函数f(t)和其频谱函数都是连续函数，不能够直接用计算机来处理，因此在进行数值计算时必须将其离散化，然后利用离散傅里叶变换（DFT）实现近似计算。

设对连续时间信号f(t)的截取时间段长度为L，对其进行离散化的采样时间间隔为T，那么采样输出的离散时间序列f(nT)中的信号样值点数N为：⎢⎣L/T⎥⎦+1。

例如，若信号f(t)的大部分能量集中于频率段[0,fm]和时间段[t0,t1]上，那么就可以根据采样定理选取采样时间间隔为T≤1/(2fm)，并选取截取时间段长度为L≥t1−t0。这样就在离散时间点序列t=nT,n=0,1,L,N−1上对连续时间信号进行了离散化。对f(t)的连续时间傅里叶变化中的积分进行近似求和计算，即：

该式表明，利用DFT(FFT)计算连续时间傅里叶变换的频谱时，除了计算时域样点的离散傅里叶变化的频谱F(k)，还要将F(k)乘以取样时间间隔T，才能得出结果。

1965直接计算DFT的复杂度为序列长度的平方，即O(N2)。

年，库利（T.W.Cooley）和图基（J.W.Tuky）在《计算数学》杂志上发表了著名的《机器计算傅里叶级数的一种算法》论文后，桑德—图基等快速算法相继出现，又经人们进行改进，很快形成了一套高效计算方法，这就是现在的快速傅里叶算法FFT。其思想是将长序列分解为若干短序列进行DFT计算，然后通过若干旋转因子的复数乘法和复数加法合成最终的结果。FFT算法中，如果序列长度是2的幂次，可将序列长为N的DFT分割为两个长为N/2的子序列的DFT，成为基2-FFT。快算傅里叶变化算法只需要O(NlogN)的计算复杂度，FFT算法为数字信号处理技术应用于各种信号的实时处理创造了条件，大大推动了数字信号处理技术的发展，在此利用FFT计算连续时间傅里叶变化的频谱。

（三）频谱分析的若干问题讨论

通过以上原理可知，连续非周期信号频谱的数值计算必须首先对信号时域采样，得到时间离散化的信号，时域采样

t1+∞必须满足或近似满足采样定理。根据时域频域的对应关系，

F(Ω)=∫f(t)exp(−jΩt)dt≈∫f(t)exp(−jΩt)dt

−∞t0

时域采样将导致所得的抽样信号频谱周期化。然而，为了使周期化后的抽样信号频谱便于计算机处理，还必须再将其频N−1

≈Tf(nT)exp(−jΩnT) （1） 域离散化，方法是对该频谱进行频域采样。根据时域频域的

n=0

对应关系，频域离散化将对应于时域信号的周期化。因此，

但是，这样计算出来的结果F(Ω)仍然是连续函数，计算

对于连续非周期信号频谱进行数值计算时，要确定如何截取

∑

【收稿日期】2011-01-18

【作者简介】井敏英（1978－），女，陕西富平人，陕西理工学院物理系硕士，从事通信信号处理的研究。 - 42 -

信号的时间段、如何选择时域采样率，以及在时间段上对信号进行截取的方式。截取信号的时间段长度决定了时域周期化的周期，对应于频率抽样的频率间隔，即频率分辨率；时域采样率决定了频域周期化的周期，即频谱数值计算的范围；而在某时间段上对信号进行截取的方式，即不同窗函数的应用，决定了信号频谱估计的精度和有效范围。

设要分析连续时间非周期信号f(t)在频率范围[0,fm]内的频谱，且要求分析的频谱分辨率（数值计算的频率间隔）为Δf，则首先应根据信号频率范围确定采样率，在根据所要求的频率分辨率确定截取时间长度，从而计算出所需计算FFT的序列长度（点数），最后根据信号时域波形特性选择使用不同的窗函数。

1.根据分析的信号频率范围确定采样率。要分析信号在频率范围[0,fm]内的频谱，则采样率必须满足采样定理，即

fs≥2fm，相应地，采样时间间隔T（也称为时间分辨率）满

fthamming=hamming(N).*ft./sqrt(sum(abs(hamming(N).^2

))./N);

fthann=hann(N).*ft./sqrt(sum(abs(hann(N).^2))./N); fwrectwin=T.*fft(ftrectwin,N); fwhamming=T.*fft(fthamming,N); fwhann=T.*fft(fthann,N);

subplot(3,1,1);plot(freq,10*log(abs(fwrectwin)));title('矩形窗频谱');

ylabel('幅度dB');axis([0,200,-50,0]);grid on; subplot(3,1,2);plot(freq,10*log(abs(fwhamming))); title('汉明窗频谱');ylabel('幅度dB');axis([0,200,-50,0]);grid on;

subplot(3,1,3);plot(freq,10*log(abs(fwhann)));title('汉宁窗频谱');

ylabel('幅度dB');axis([0,200,-50,0]);grid on;

程序执行后结果如图1所示：

足T=1/fs≤1/2fm。

2.根据频率分辨率要求确定分析信号f(t)的截取时间长度。要使所分析的频率分辨率达到Δf，即每隔Δf计算一个频率点，那么对信号的截取时间长度L必须满足L≥1/Δf，根据截取时间长度L和采样时间间隔T就可以计算出截取时间信号离散化之后的序列点数N，也可以由计算采样率fs和频率间隔Δf来等价计算出序列点数，即：

图1 3种加窗后的幅度频谱图

可以看出，图1中3种加窗后对信号的频谱估计与理论值都有一定的误差，这是由于在实际分析过程中，要对连续信号采样和截断。加矩形窗的频谱误差最大，那是因为事实上，加矩形窗等价于截取时不作加窗处理。从图1中3种加窗后的幅度谱估计曲线来看，应用汉明窗和汉宁窗都比应用矩形窗（等价于不加窗）的估计精度要高。

N=⎢⎥+1=⎣⎢fs/Δf⎦⎥+1 （5） ⎣L/T⎦

3.根据信号时域波形特性来应用不同的窗函数。使用窗

函数可以控制频谱主瓣宽度、旁瓣抑制度等参数，从而更好地进行波形频谱分析和滤波器参数的设计。将窗函数与信号的时域波形或频谱进行相乘的过程，就称为对信号做时域加窗和频域加窗。不同窗函数与信号时域波形相乘就是以不同窗函数对时间无限长的连续信号f(t)进行时间段截取。

（四）分析实例

试对x(t)=sin(2π50t)+sin(2π80t)进行频谱分析，要求分析的频率范围为0~100Hz，频率分辨率为1Hz。根据所要求的分析频率范围可以确定信号的采样率为fs=200Hz，采样时间间隔（时间分辨率）为T=1/fs=5ms，而根据要求的频率分辨率可以得出信号时域截断长度为L=1/Δf=1s。因此，对截断信号的采样点数为：

N=⎢⎥+1=⎣⎢fs/Δf⎦⎥+1=201 （6） ⎣L/T⎦

（五）结论

以上对带限确知连续时间信号的频谱分析的数值计算方

法进行了分析，在分析过程中要对信号进行采样和截断，因此会出现一定误差，但只要合理选择分析参数可以使误差在工程允许的范围之内。对于随机信号，由于是无限大能量的功率信号，它不满足傅里叶变换条件，而且也不存在解析表达式，因此就不能够应用确知信号的频谱计算方法去分析随机信号的频谱。但可以用以上确知信号频谱分析的方法分析随机信号的功率谱，因此同样适应于随机信号的频域研究。

【参考文献】

[1] 高西全,丁玉美.数字信号处理[M].西安:西安电子科技大学

出版社,2010.90-110.

[2] 邵玉斌.Matlab/Simulink通信系统建模与仿真实例分析[M].

北京:清华大学出版社,2008:166-180.

[3] 张登奇,杨慧银.信号的频谱分析及MATLAB实现[J].湖南

理工学院学报(自然科学版),2010,23(3):29-33.

现应用MATLAB语言工具，使用数值计算的方法对其进行频谱估计。分别用矩形窗、汉明窗和汉宁窗进行时域加窗，分析信号幅度谱曲线，程序代码如下：

fs=200;deltaf=1;T=1/fs;L=1/deltaf;N=floor(fs/deltaf)+1;t=0:T:L;

freq=0:deltaf:fs;

ft=(sin(2*pi*50*t)+sin(2*pi*80*t))';

ftrectwin=rectwin(N).*ft./sqrt(sum(abs(rectwin(N).^2))./N);

- 43 -