实验四离散信号频谱分析

一、 实验目的

1. 2.

掌握序列的傅里叶变换、离散傅里叶级数、离散傅里叶变换、快速傅里叶变换的Matlab实现；

学习用FFT对连续信号和离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT。

二、 实验原理及方法

1. 离散非周期信号的谱分析 （1） 序列的傅里叶变换

对于满足绝对可和的序列，即|x(n)|，其傅里叶变换和反变换的定义为

jjnX(e)12nx(n)e （4.1）

jx(n)X(e)ejnd （4.2）

序列x(n)是离散的，但X(ej)是以2为周期的的连续函数，为了能够在计算机上处理，需要对x(n)进行截断，对频域进行离散化，近似处理后

X(ejkn2)nn1x(n)ejkn（4.3）

k，M是对在一个周期内的采样，k的取值由读者确定，若想观察一个M周期内的频谱，k0~M1，若观察两个周期，k0~2M1，以此类推。

序列傅里叶变换的Matlab实现： n=n1:n2;

M=input(„put in the number M=‟); k=0:2*M-1; %观察两个周期

其中k2X=x*(exp(-j*2*pi/M)).^(n‟*k)；%序列的傅里叶变换 对R4(n)进行序列的傅里叶变换得到图4-1。

x(n)14|X(e)|j0arg[X(ej)]n0123202图4-1 信号及信号的幅度谱和相位谱

（2）离散傅里叶变换（DFT）

如果序列x(n)是有限长的，序列的谱分析可以采用离散傅里叶变换，其定义为：

59

N1X(k)DFT[x(n)]n0N1x(n)WNkn,0kN1 （4.4）

x(n)IDFT[X(k)]1NX(k)Wk0knN,0nN1 （4.5）

因为x(n)与X(k)都是离散的，所以可以利用计算机进行数值计算。从数学观点看，DFT表示的是对序列x(n)或X(k)的线性运算。

此处应用DFT变换近似分析采样序列的频谱。设时域序列用x表示，长度为M；x的DFT变换为X，变换区间长度为N（NM）。

M1X(k)x(n)Wn0knN k0,1,,N1

将X(k)展开，得：

XXXX(0)x(0)WNx(1)WNx(2)WNx(M1)WN(1)x(0)WNx(1)WNx(2)WNx(M1)WN2021221011120001020(M1)1(M1)2(M1)(2)x(0)WNx(1)WNx(2)WNx(M1)WN(N1)x(0)WN(N1)0

(N1)(M1)x(1)WN(N1)1x(2)WN(N1)2x(M1)WN将上式表示成矩阵的形式：

(N1)0W00　　　W10　　　WNNN(N1)10111WN　　　WN　　　WN(N1)202１2[X(0)　X(1)　X(N1)][x(0)　x(1)　x(M1)]WN　　　WN　　　WN0(M1)1(M1)(N1)(M1)WN　WN　　WN

DFT变换的Matlab实现：

n=0:M-1;

k=0:N-1;

WN=exp(-j*2*pi/N); kn=n'*k;

WNkn=WN.^kn; X=x*WNkn;

对R4(n)进行离散傅里叶变换得到图4-2。

60

14420.520001x(n)230051015-2abs(X(k))0105angle(X(k))15

图4-2 信号及信号的离散傅里叶变换

（3）快速傅里叶变换（FFT）

快速傅里叶变换并不是一种新的变换，只是离散傅里叶变换的快速算法，常用的有按时间抽取的基-2FFT算法和按频率抽取的基-2FFT算法。在Matlab中对离散信号进行FFT，可以直接调用函数。快速傅里叶变换的原理及子程序见附录。

fft(x)：利用快速算法计算x的M点DFT，其中M是x的长度。

fft(x,N)：利用快速算法计算x的N点DFT，其中N是用户指定的长度。分两种情况： ①若x的长度M>N，则将x截短为N点序列，再作N点DFT； ②若x的长度Mifft(X)：利用快速算法计算X的M点IDFT，其中M是X的长度。

ifft(X,N)：利用快速算法计算X的N点IDFT，其中N是用户指定的长度。同样分两种情况，同fft(x,N)。

对R4(n)进行16点的快速傅里叶变换得到图4-3。

10.500420020-21x(n)235abs(X(k))10150angle(X(k))51015图4-3 信号及信号的快速傅里叶变换

观察图4-1、4-2、4-3可以得到，DFT和FFT都是对序列的傅里叶变换X(ej)在[0,2]上的离散化，离散化的采样点数即为做DFT和FFT的点数。而FFT仅为DFT的快速运算，当做DFT和FFT的点数相同时，两者的结果是一样的。 2. 离散周期信号的谱分析 （1）离散傅里叶级数 定义WNej2N，周期序列的傅里叶级数（DFS）变换对为：

N1~X(k)DFS[~x(n)]n0~x(n)eN1j2NnkN12Nn0nk~x(n)WNk0,1,,N11~~x(n)IDFS[X(k)]Nk0j~X(k)enk1NN1 （4.6）

~nkX(k)WNn0,1,,N1k0式中，n和k都是离散变量。如果将n当作时间变量，k当作频率变量，则DFS[]表示时域到频域的离散傅里叶级数正变换，IDFS[]表示由频域到时域的离散傅里叶级数反变换。

61

从式（4.6）看出，只要知道周期序列一个周期的内容，其它周期的内容也都知道了。所以，这种无限长周期序列实际上只有一个周期中的N个序列值有信息。

将n和k进行周期延拓可得，

~X(k)N1n0kn~x(n)WN,k~ （4.7）

1~x(n)NN1k0~knX(k)WN，n~ （4.8）

1）上式表明将周期序列分解成N次谐波，第k次谐波频率为k(2/N)k，

~k0,1,2,,N1，幅度为(1N)X(k)。

2）基波分量的频率是2N，幅度是(1N)X(1)。一个周期序列可以用离散傅里叶级数表示其频谱分布规律。

由（4.6）式，时域和频域都是有限长序列，可以用Matlab实现： WN=exp(-j*2*pi/N); kn=n'*k; WNkn=WN.^kn;

X=x*WNkn;

（2）快速傅里叶变换 因为周期序列以及周期序列的频谱都是有限长序列，所以可以直接调用Matlab中的fft（）和ifft（）函数进行求解。

3. 利用FFT对连续信号进行谱分析

工程实际中，经常遇到的连续信号xa(t)，其频谱函数Xa(j)也是连续函数。设时域连续信号xa(t)持续时间为Tp，最高频率为fc。xa(t)的傅里叶变换为

Xa(jf)FT[xa(t)]对xa(t)以采样频率fs2fc (TXa(if)作零阶近似(t=nT，dt=T)得

N-1~xa(t)ej2πftdt （4.9）

1fs)采样得xa(t) xa(nT)。设共采样N点，并对

X(if)Txa(nT)ej2fnT （4.10）

n0显然，X(if)仍是f的连续周期函数。对X(jf)在区间［0, fs］上等间隔采样N点，采样间隔为F，F称为频率分辨率。参数fs、Tp、N和F满足如下关系式

62

F由于NT=Tp，所以

fsN1NT （4.11）

F1Tp （4.12）

式（4.12）说明要提高频率分辨率即减小F就要增加Tp，即增加信号的有效长度。

将fkF和式（4.11）代入式（4.10）可得X(jf)的采样

N-1j2Nkn X(jkF)Txa(nT)en0， 0kN-1 （4.13）

令Xa(k)X(jkf),x(n)xa(nT)，则

N1Xa(k)Tx(n)en0j2NknTDFT[x(n)] （4.14）

式（4.14）说明连续信号的频谱特性可以通过对连续信号采样并进行DFT，再乘以T的近似

方法得到。

同理，由

xa(t)Xa(jf)ej2πftdf

可推导出

N-1x(n)xa(nT)FXa(k)en0j2NknFN[1T1NN1n0Xa(k)ej2Nkn] （4.15）

IDFT[Xa(k)]这里的DFT和IDFT都可以直接调用fft（）和ifft（）。 4.利用FFT进行谱分析的误差分析

下面分析利用FFT对信号进行傅里叶分析时可能造成的误差。 （1）频谱混叠失真

利用FFT对连续信号进行傅里叶分析时首先需要进行时域采样。按照时域采样定理，为了不产生混叠，必须满足

fs2fc （4.16） 也就是采样间隔T满足

T1fs12fc （4.17）

63

一般应取

fs(3~5f) （4.18） c如果不满足fs2fh，就会产生频谱混叠失真。  （2）栅栏效应

利用FFT计算频谱，只能给出离散点k2Nk或k2NTk上的频谱采样值，而不

可能得到连续频谱函数，这就像通过一个“栅栏”观看信号频谱，只能在离散点上看到信号频谱，称之为“栅栏效应”。这时，如果在两个离散的谱线之间有一个特别大的频谱分量，就无法检测出来了。

减小栅栏效应的一个方法就是要使频域采样更密，即增加频域采样点数N，在不改变时域数据的情况下，必然是在数据末端添加一些零值点，使一个周期内的点数增加，但并不改变原有的记录数据。频谱采样间距为

2N，N增加，必然使样点间距更近（单位圆上样点更

多），谱线更密，谱线变密后原来看不到的谱分量就有可能看到了。 

必须指出，补零以改变计算FFT的周期时，不能提高频率分辨率，这是因为数据的实际长度仍为补零前的数据长度。

（3）频谱泄漏与谱间干扰

对信号进行FFT计算，首先必须使其变成有限时宽的信号， 这就相当于信号在时域乘一个窗函数如矩形窗，窗内数据并不改变。加窗即为时域相乘，对频域的影响可用卷积公式表示

V(ej)1212X(ej)W(ejj) （4.19）

X(e)W(ej())d式（4.19）中，X(e加窗后得到的频谱V(ejj)是信号x(n)的频谱，W(ejj)是窗函数的频谱。卷积的结果使

)与原来的频谱X(e产生失真。这种失真造成频谱的“扩)不相同，

散”，这就是所谓的“频谱泄漏”。 

三、 实验内容及步骤

1. 计算序列的DTFT和DFT，观察栅栏效应 设

x(n)R4(n)，要求用MATLAB实现：

j （1）计算x(n)的傅里叶变换X(e)，并绘出其幅度谱；

64

（2）分别计算x(n)的4点DFT和8点DFT，绘出其幅度谱。并说明它们和X(ej)的关系。 （提示：DFT变换可用MATLAB提供的函数fft实现，也可以自己用C语言或matlab编写） 2.计算序列的FFT，观察频谱泄漏 已知周期为16的信号x(n)cos(1016n)cos(1216n)。

（1） 截取一个周期长度M=16点，计算其16点FFT，并绘出其幅度谱；

（2） 截取序列长度M=10点，计算其16点FFT，绘出其幅度谱，并与（1）的结果进行比

较，观察频谱泄漏现象，说明产生频谱泄漏的原因。

四、 实验报告要求

1. 2.

结合实验中所得给定典型序列幅频特性曲线，与理论结果比较，并分析说明误差产生总结实验所得主要结论。

的原因以及用FFT作谱分析时有关参数的选择方法。

65