1、学生实验报告开课学院及实验室： 电子楼 317 2013 年 4 月 29 日学院机械与电气工程学院年级、专业、班姓名 学号实验课程名称 数字信号处理实验 成绩实验项目名称 实验五 用 DFT（FFT ）对信号进行频谱分析 指导老师一、实验目的学习 DFT 的基本性质及对时域离散信号进行频谱分析的方法，进一步加深对频域概念和数字频率的理解，掌握 MATLAB 函数中 FFT 函数的应用。二、实验原理离散傅里叶变换（DFT）对有限长时域离散信号的频谱进行等间隔采样，频域函数被离散化了，便于信号的计算机处理。设 x(n)是一个长度为 M 的有限长序列， x(n)的 N 点傅立叶变换：其中 NjeW
2、2，它的反变换定义为： 10)()(NknkWXnx令kNz，则有： 10)()(nnkNNzx可以得到， kNWzX)(，kN是 Z 平面单位圆上幅角为kN2的点，就是将单位圆进行 N 等分以后第 K 个点。所以，X(K)是 Z 变换在单位圆上的等距采样，或者说是序列傅立叶变换的等距采样。时域采样在满足 Nyquist 定理时，就不会发生频谱混叠。DFT 是对序列傅立叶变换的等距采样，因此可以用于序列的频谱分析。如果用 FFT 对模拟信号进行谱分析，首先要把模拟信号转换成数字信号，转换时要求知道模拟信号的最高截至频率，以便选择满足采样定理的采样频率。一般选择采样频率是模拟信号中最高频率的 3
3、4 倍。另外要选择对模拟信号的观测时间，如果采样频率和观测时间确定，则采样点数也确定了。这里观测时间和对模拟信号进行谱分析的分辨率有关，最小的观测时间和分辨率成倒数关系。最小的采样点数用教材相关公式确定。要求选择的采样点数和观测时间大于它的最小值。如果要进行谱分析的模拟信号是周期信号，最好选择观测时间是信号周期的整数倍。如果不知道信号的周期，要尽量选择观测时间长一些，以减少阶段效应的影响。用 FFT 对模拟信号作谱分析是一种近似的谱分析，首先一般模拟信号（除周期信号以外）的频谱是连续谱，而用 FFT 作谱分析得到的是数字谱，因此应该取 FFT 的点数多一些，用它的包络作为模拟信号的近似谱。另外
4、，如果模拟信号不是严格的带限信号，会因为频谱混叠现象引起谱分析的误差，这种情况下可以预先将模拟信号进行滤波，或者尽量采样频率取高一些。一般频率混叠发生在折叠频率附近，分析时要注意因频率混叠引起的误差。最后要注意一般模拟信号是无限长的，分析时要截断，截断的长度和分辨率有关，但也要尽量取长一些，取得太短会21j0)()e01NkXkDFTx截断引起的误差会很大。举一个极端的例子，一个周期性正弦波，如果所取观察时间太短，例如取小于一个周期，它的波形和正弦波相差太大，肯定误差很大，但如果取得长一些，即使不是周期的倍数，这种截断效应也会小一些。如同理论课教材所讨论的，在运用 DFT 进行频谱分析的时候可
5、能有三种误差，即：（1）混叠现象当采样率不满足 Nyquist 定理，经过采样就会发生频谱混叠。这导致采样后的信号序列频谱不能真实的反映原信号的频谱。所以，在利用 DFT 分析连续信号频谱的时候，必须注意这一问题。避免混叠现象的唯一方法是保证采样的速率足够高，使频谱交叠的现象不出现。这告诉我们，在确定信号的采样频率之前，需要对频谱的性质有所了解。在一般的情况下，为了保证高于折叠频率的分量不会出现，在采样之前，先用低通模拟滤波器对信号进行滤波。（2）泄漏现象实际中的信号序列往往很长，甚至是无限长。为了方便，我们往往用截短的序列来近似它们。这样可以使用较短的 DFT 来对信号进行频谱分析。这种截短
6、等价于给原信号序列乘以一个矩形窗函数。值得一提的是，泄漏是不能和混叠完全分离开的，因为泄漏导致频谱的扩展，从而造成混叠。为了减少泄漏的影响，可以选择适当的窗函数使频谱的扩散减到最小。（3）栅栏效应因为 DFT 是对单位圆上 Z 变换的均匀采样，所以它不可能将频谱视为一个连续函数。这样就产生了栅栏效应，从某种角度看， 用 DFT 来观看频谱就好像通过一个栅栏来观看一幅景象，只能在离散点上看到真实的频谱。这样的话就会有一些频谱的峰点或谷点被“栅栏”挡住，不能被我们观察到。减小栅栏效应的一个方法是在源序列的末端补一些零值，从而变动 DFT 的点数。这种方法的实质是改变了真是频谱采样的点数和位置，相当
7、于搬动了“栅栏”的位置，从而使得原来被挡住的一些频谱的峰点或谷点显露出来。注意，这时候每根谱线所对应的频和原来的已经不相同了。从上面的分析过程可以看出，DFT 可以用于信号的频谱分析，但必须注意可能产生的误差，在应用过程中要尽可能减小和消除这些误差的影响。DFT 运算量较大，快速离散傅里叶变换算法 FFT 是解决方案。 FFT 并不是 DFT 不相同的另一种变换，而是为了减少 DFT 运算次数的一种快速算法。它是对变换式进行一次次的分解，使其成为若干小点数 DFT 的组合，从而减小运算量。常用的 FFT 是以 2 为基数的，其长度为 MN2。它的运算效率高，程序比较简单，使用也十分的方便。当需
8、要进行变换的序列的长度不是 2 的整数次方的时候，为了使用以 2 为基的 FFT，可以用末尾补零的方法，使其长度延长至 2 的整数次方。IFFT 一般也可以通过 FFT 程序来完成。三、使用仪器、材料1、硬件：计算机 2、软件：Matlab四、实验步骤(一) 离散信号给定参考实验信号如下：用 14()xnR以 8 为周期进行周期性延拓形成的周期序列。3()(1) 分别以变换区间 N8，16，32，对 14()xnR进行 DFT(FFT)，画出相应的幅频特性曲线；(2) 分别以变换区间 N4，8 ，16，对 x2(n)分别进行 DFT(FFT)，画出相应的幅频特性曲线；(3) 对 x3(n)进行
9、频谱分析，并选择变换区间，画出幅频特性曲线。（二）连续信号1. 实验信号： 1()xtRt选择 ，式中 的波形以及幅度特性如图 7.1 所示。1.5ms()Rt2()sin/8)xtft式中频率 自己选择。f3()cos816cos20xtttt2. 分别对三种模拟信号选择采样频率和采样点数。对 ，选择采样频率 ，8 ，16 ，采样点数用 . 计算。1()xtRt4sfkHzkzsf对 ，周期 ，频率 自己选择，采样频率 ，观测时间2sin/8)f1/Tff 4s， ， ，采样点数用 计算。0.5pTTps图 5.1 R(t)的波形及其幅度特性对 ，选择采用频率 ，采样点数为 16，32，64
10、。3()cos816cos20xtttt64sfHz3. 分别对它们转换成序列，按顺序用 表示。123(),()xnx4. 分别对它们进行 FFT。如果采样点数不满足 2 的整数幂，可以通过序列尾部加 0 满足。5. 计算幅度特性并进行打印。五、实验过程原始记录（数据、图表、计算等）(一) 离散信号% 14()xnRn=0:1:10;xn=ones(1,4),zeros(1,7); %输入时域序列向量 xn=R4(n)Xk8=fft(xn,8); %计算 xn 的 8 点 DFTXk16=fft(xn,16); %计算 xn 的 16 点 DFTXk32=fft(xn,32); %计算 xn
11、的 32 点 DFTk=0:7;wk=2*k/8; %产生 8 点 DFT 对应的采样点频率 (关于 归一化值)subplot(2,2,1);stem(n,xn,.);title(a) x_1 (n);xlabel(n);ylabel(x_1 (n);subplot(2,2,2);stem(wk,abs(Xk8),.); %绘制 8 点 DFT 的幅频特性图title(b) 8 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);k=0:15;wk=2*k/16; %产生 16 点 DFT 对应的采样点频率 (关于 归一化值)subplot(2,2,3);ste
12、m(wk,abs(Xk16),.); %绘制 16 点 DFT 的幅频特性图title(c) 16 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);k=0:31;wk=2*k/32; %产生 32 点 DFT 对应的采样点频率 (关于 归一化值)subplot(2,2,4);stem(wk,abs(Xk32),.); %绘制 32 点 DFT 的幅频特性图title(d) 32 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);运行结果：%x2(n)n=0:50;xn=cos(pi/4*n); %输入时域序列向量 cos(p
13、i/4*n)Xk4=fft(xn,4); %计算 xn 的 4 点 DFTXk8=fft(xn,8); %计算 xn 的 8 点 DFTXk16=fft(xn,16); %计算 xn 的 16 点 DFTk=0:3;wk=2*k/4; %产生 4 点 DFT 对应的采样点频率 (关于 归一化值)subplot(2,2,1);stem(n,xn,.); %绘制 4 点 DFT 的幅频特性图title(a) x_2 (n);xlabel(n);ylabel(x_2 (n);subplot(2,2,2);stem(wk,abs(Xk4),.); %绘制 4 点 DFT 的幅频特性图title(b)
14、4 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);axis(0,2,0,2) k=0:7;wk=2*k/8; %产生 8 点 DFT 对应的采样点频率 (关于 归一化值)subplot(2,2,3);stem(wk,abs(Xk8),.); %绘制 8 点 DFT 的幅频特性图title(c) 8 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);k=0:15;wk=2*k/16; %产生 32 点 DFT 对应的采样点频率 (关于 归一化值)subplot(2,2,4);stem(wk,abs(Xk16),.); %绘制
15、 16 点 DFT 的幅频特性图title(d) 16 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);运行结果：%x3(n)n=0:3;x=ones(1,4),zeros(1,4);xn=x(mod(n,8)+1); %输入 xn=R4(n)以 8 为周期进行周期性延拓形成的周期序列Xk8=fft(xn,8); %计算 xn 的 8 点 DFTXk16=fft(xn,16); %计算 xn 的 16 点 DFTXk32=fft(xn,32); %计算 xn 的 32 点 DFTk=0:7;wk=2*k/8; %产生 8 点 DFT 对应的采样点频率 (关于
16、 归一化值)subplot(3,1,1);stem(wk,abs(Xk8),.); %绘制 8 点 DFT 的幅频特性图title(a) 8 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);k=0:15;wk=2*k/16; %产生 16 点 DFT 对应的采样点频率 (关于 归一化值)subplot(3,1,2);stem(wk,abs(Xk16),.); %绘制 16 点 DFT 的幅频特性图title(c) 16 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);k=0:31;wk=2*k/32; %产生 32 点 D
17、FT 对应的采样点频率 (关于 归一化值)subplot(3,1,3);stem(wk,abs(Xk32),.); %绘制 32 点 DFT 的幅频特性图title(e) 32 点 DFT 的幅频特性图);xlabel(omega/pi);ylabel(幅度);运行结果：（二）连续信号%x1(n)Fs1=4000; %采样频率Fs2=8000;Fs3=16000;Tp=0.002;N1=Fs1*Tp; %计算采样点数N2=Fs2*Tp;N3=Fs3*Tp;n1=0:N1-1;n2=0:N2-1;n3=0:N3-1;n=-0.0005:0.0001:0.002;xn=(n=0); %x1(n)X
18、k1=fft(xn,8); %计算 xn 的 8 点 DFTk1=0:length(abs(Xk1)-1;wk1=2*k1/8; %产生 8 点 DFT 对应的采样点频率(关于 归一化值)Xk2=fft(xn,16); %计算 xn 的 16 点 DFTk2=0:length(abs(Xk2)-1;wk2=2*k2/16; %产生 16 点 DFT 对应的采样点频率(关于 归一化值)Xk3=fft(xn,32); %计算 xn 的 32 点 DFTk3=0:length(abs(Xk3)-1;wk3=2*k3/32; %产生 32 点 DFT 对应的采样点频率(关于 归一化值)subplot(
19、3,1,1);stem(wk1,abs(Xk1),.);xlabel(omega/pi);ylabel(幅度);title(X_1 (n) 8 点傅立叶变换(f_s =4kHz);subplot(3,1,2);stem(wk2,abs(Xk2),.);xlabel(omega/pi);ylabel(幅度);title(X_2 (n) 16 点傅立叶变换(f_s =8kHz);subplot(3,1,3);stem(wk3,abs(Xk3),.);xlabel(omega/pi);ylabel(幅度);title(X_1 (n) 32 点傅立叶变换(f_s =16kHz);运行结果：%x2(n)
20、,取 f=2kHzf=2000;Fs=4*f; %采样频率T=1/f;Tp1=0.5*T;Tp2=T;Tp3=2*T;Ts=1/Fs;N1=Fs*Tp1; %计算采样点数N2=Fs*Tp2;N3=Fs*Tp3;n1=0:N1-1;n2=0:N2-1;n3=0:N3-1;xn1=sin(2*pi*f*n1*Ts+pi/8); xn2=sin(2*pi*f*n2*Ts+pi/8);xn3=sin(2*pi*f*n3*Ts+pi/8);Xk1=fft(xn1,2);k1=0:length(abs(Xk1)-1;wk1=2*k1/2; %产生 2 点 DFT 对应的采样点频率(关于 归一化值)Xk2=
21、fft(xn2,4);k2=0:length(abs(Xk2)-1;wk2=2*k2/4; %产生 4 点 DFT 对应的采样点频率(关于 归一化值)Xk3=fft(xn3,8);k3=0:length(abs(Xk3)-1;wk3=2*k3/8; %产生 8 点 DFT 对应的采样点频率(关于 归一化值)subplot(3,1,1);stem(wk1,abs(Xk1),.);xlabel(omega/pi);ylabel(幅度);title(X_2 (n) 2 点傅立叶变换(T_p =0.5T);subplot(3,1,2);stem(wk2,abs(Xk2),.);xlabel(omega
22、/pi);ylabel(幅度);title(X_2 (n) 4 点傅立叶变换(T_p =T);subplot(3,1,3);stem(wk3,abs(Xk3),.);xlabel(omega/pi);ylabel(幅度);title(X_2 (n) 8 点傅立叶变换(T_p =2T);运行结果：%x3(n)Fs=64;T=1/Fs;N=64;n=0:N-1;x3n=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);Xk16=fft(x3n,16);k1=0:length(abs(Xk16)-1;wk3=2*k1/8; %产生 8 点 DFT 对应的采样点频
23、率(关于 归一化值)Xk32=fft(x3n,32);k2=0:length(abs(Xk32)-1;wk2=2*k2/16; %产生 16 点 DFT 对应的采样点频率(关于 归一化值)Xk64=fft(x3n,64);k3=0:length(abs(Xk64)-1;wk3=2*k3/64; %产生 64 点 DFT 对应的采样点频率(关于 归一化值)subplot(3,1,1);stem(k1,abs(Xk16),.);xlabel(omega/pi);ylabel(幅度);title(X_3 (n) 16 点傅立叶变换);subplot(3,1,2);stem(k2,abs(Xk32),
24、.);xlabel(omega/pi);ylabel(幅度);title(X_3 (n) 32 点傅立叶变换);subplot(3,1,3);stem(k3,abs(Xk64),.);xlabel(omega/pi);ylabel(幅度);title(X_3 (n) 64 点傅立叶变换);运行结果：六、实验结果及分析1.分析 DFT 的变换区间对频域分析的作用，并说明 DFT 的物理意义。答:以(一) 离散信号中的 14()xnR为例。变换区间对频域分析的作用:随着变换区间的增大,36 变换区间的频域明显比 8 区间变换的频域容易识别一个周期内部的情况。即变换区间越长，信号的基频越小。 物理意
25、义:假设 x(n)是 N 点 DFT 是 x(n)的 z 变换在单位圆上的 N 点等间隔采样，则 x(k)为 x(n)的傅里叶变换 x(e(jw)在区间0,2 上的 N 点等间隔采样。2.对周期信号的变换区间应该如何选取。如果周期信号的周期预先不知道，如何用 DFT 分析它的频谱。答：1） 变换区间的选择：对于非周期信号，假设频谱分辨率 F，而频谱分辨率直接和 FFT 的变换区间有关，因为 FFT 能够实现的频率分辨率是 2/N 因此有最小的 N2/F，根据此式可以选择 FFT 的变换区间；对于周期信号，周期信号的频谱是离散谱，只有用整数倍周期的长度作 FFT，得到的离散谱才能代表周期信号的频谱。2） 周期信号的周期预先不知道时,可先截取 M 点进行 DFT，再将截取长度扩大 1 倍截取，比较结果，如果二者的差别满足分析误差要求，则可以近似表示该信号的频谱，如果不满足误差要求就继续将截取长度加倍，重复比较，直到结果满足要求。
