数字信号处理实验五

实 验 报 告

实验名称____利用DFT分析离散信号频谱 课程名称____数字信号处理________

院 系 部:电气与电子工程 专业班级：信息1002 学生姓名：王萌 学 号: 11012000219 同 组 人: 实验台号: 指导教师：范杰清 成 绩： 实验日期:

华北电力大学

一、实验目的

应用离散傅里叶变换(DFT)，分析离散信号

x[k]的频谱。深刻理解

DFT分析离散信号频谱的原理，掌握改善分析过程中产生的误差的方法。

二、实验原理

根据信号傅里叶变换建立的时域与频域之间的对应关系，可以得

到有限长序列的离散傅里叶变换(DFT)与四种确定信号傅里叶变换的之间的关系，实现由DFT分析其频谱。

Matlab中提供了fft函数，FFT是DFT的快速算法 X=fft(x)：用于计算序列x的离散傅里叶变换（DFT） X=fft(x,n)：对序列x补零或截短至n点的离散傅里叶变换。 当x的长度小于n时，在x的尾部补零使x的长度达到n点； 当x的长度大于n时，将x截短使x的长度成n点； x=ifft(X)和x=ifft(X，n)是相应的离散傅里叶反变换。 fftshift(x)将fft计算输出的零频移到输出的中心位置。

~利用DFT计算离散周期信号 x [ k ] 的频谱

分析步骤为： (1) 确定离散周期序列

~x[k] 的基本周期N；

(2) 利用fft函数求其一个周期的DFT，得到X[m];

~X[m]?X[m] (3)

。

利用DFT计算离散非周期信号x[k] 的频谱

分析步骤为：(1) 确定序列的长度M及窗函数的类型。当序列为无限长时，需要根据能量分布，进行截短。

(2) 确定作FFT的点数N；根据频域取样定理.为使时域波形 不产生

? M混叠，必须取 N 。

(3) 使用fft函数作N点FFT计算X[m]。

三、实验内容

1、利用FFT分析信号

的频谱；

(1) 确定DFT计算的参数；

(2) 进行理论值与计算值比较，讨论信号频谱分析过程中误差 原因及改善方法。 解答：

%该周期序列的周期N=32，基频W0=p/8。

N=32; k=0:N-1; x=cos(pi*3/8*k); X=fft(x,N); subplot(2,1,1);

stem(k-N/2,abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)'); subplot(2,1,2);

stem(k-N/2,angle(fftshift(X)));

ylabel('Phase'); xlabel('Frequency (rad)');

3πx[k]?cos(k), k?0,1,?318

改善方法：应该多抽取些 2. 利用的频谱；

(1) 确定DFT计算的参数；

(2) 进行理论值与计算值比较，讨论信号频谱分析过程中误差 原因及改善方法。 解答：

k=0:30; x=0.5.^k;

subplot(2,1,1);

stem(k,x); %画出序列的时域波形 subplot(2,1,2); w=k-15;

plot(w, abs(fftshift(fft(x)))); %画出序列频谱的幅度谱

FFT分析信号

kx[k]?(1)u[k] 2

改善：多抽取些 3. 有限长脉冲序列

x[k]?[2,3,3,1,0,5;k?0,1,2,3,4,5]，

利用FFT分析其频谱，并绘出其幅度谱与相位谱。 解答：

N=6; k=0:N-1;

x=[2,3,3,1,0,5]; X=fft(x,N); subplot(2,1,1);

stem(k-N/2,abs(fftshift(X))); ylabel('Magnitude'); xlabel('Frequency (rad)');

title('Magnitude—Frequency N=6 '); subplot(2,1,2);

stem(k-N/2,angle(fftshift(X))); ylabel('Phase'); xlabel('Frequency (rad)');

title('Phase—Frequency ');

x [k? cos(4. 某周期序列由3个频率组成： ] ? cos( 16 k ) 16 k ) ? cos( k ) ，27π9ππ利用FFT分析其频谱。如何选取FFT的点数N？此3个频率分别对应FFT计算结果X[m]中的哪些点？若选取的N不合适，FFT计算出的

频谱X[m]会出现什么情况？ 解答：

%该周期序列的周期N=32，N=64,基频W0=p/16。

N=32; k=0:N-1;

x=cos(pi*7/16*k)+cos(9*pi/16*k)+cos(8*pi/16*k); X=fft(x,N); subplot(4,1,1);

stem(k-N/2,abs(fftshift(X))); ylabel('Magnitude');

xlabel('Frequency (rad)');

title('Magnitude—Frequency N=32 '); subplot(4,1,2);

stem(k-N/2,angle(fftshift(X))); ylabel('Phase'); xlabel('Frequency (rad)');

title('Phase—Frequency N=32'); N=64; k=0:N-1;

x=cos(pi*7/16*k)+cos(9*pi/16*k)+cos(8*pi/16*k); X=fft(x,N); subplot(4,1,3);

stem(k-N/2,abs(fftshift(X))); ylabel('Magnitude'); xlabel('Frequency (rad)');

title('Magnitude—Frequency N=64 '); subplot(4,1,4);

stem(k-N/2,angle(fftshift(X))); ylabel('Phase'); xlabel('Frequency (rad)');

title('Phase—Frequency N=64');

若N选取不合适，则无法区分频率

5. 某离散序列 析其频谱。

(1) 对x[k]做64点FFT，绘出信号频谱，能分辨出其中的两个频率吗？

(2) 对x[k]补零到256点后计算FFT，能分辨出其中的两个频率吗？ (3) 选用非矩形窗计算FFT，能够分辨出其中的两个频率吗？ (4) 若不能够很好地分辨出其中的两个频谱，应采取哪些措施？ 解答：

（1）、该周期序列的周期N=64，基频W0=p/15 N=64;

k=0:N-1;

2π2.3πx[k]?cos(k)?0.75cos( 15 15 k ),

? 63 0 ? k ，利用FFT分

x=cos(pi*2/15*k)+0.75*cos(2.3*pi/15*k); X=fft(x,N); subplot(2,1,1);

stem(k-N/2,abs(fftshift(X))); ylabel('Magnitude'); xlabel('Frequency (rad)');

title('Magnitude—Frequency N=64 '); subplot(2,1,2);

stem(k-N/2,angle(fftshift(X))); ylabel('Phase'); xlabel('Frequency (rad)');

title('Phase—Frequency N=64');

（2）对x[k]补零到256点后计算FFT

N=256;

k=0:N-1;

x=cos(pi*2/15*k)+0.75*cos(2.3*pi/15*k); X=fft(x,N); subplot(2,1,1);

stem(k-N/2,abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)');

title('Magnitude—Frequency N=256'); subplot(2,1,2);

stem(k-N/2,angle(fftshift(X))); ylabel('Phase'); xlabel('Frequency (rad)');

title('Phase—Frequency N=256');

（3）

Hanning:

N=64; k=0:N-1;

w=1/2*(1-cos(2*pi*k/(N-1)));

x=cos(2*pi/15*k).*w+0.75*cos(2.3*pi/15*k).*w; X=fft(x,N); subplot(2,1,1);

stem(k,abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)'); subplot(2,1,2);

stem(k,angle(fftshift(X)));

ylabel('Phase'); xlabel('Frequency (rad)');

Blackman:

N=64; k=0:N-1;

w=0.42-0.5*cos(2*pi*k/(N-1))+0.08*cos(4*pi*k/(N-1)); x=cos(2*pi/15*k).*w+0.75*cos(2.3*pi/15*k).*w; X=fft(x,N); subplot(2,1,1);

stem(k,abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)'); subplot(2,1,2);

stem(k,angle(fftshift(X)));

ylabel('Phase'); xlabel('Frequency (rad)');

??(0.1k)?2,6. 已知序列x[k]??e??0,2k?50 others利用FFT分析下列信号的幅频特性，频率范围为 Ω?[?π,π) ，N=500点。

(1) y[k]?x[2k]； (2) g[k]?x[4k]；

(3) 若将上述x[k]乘以cos(pk/2) ，重做(1)和(2)。 解答：(1)%y[k]?x[2k]

N=500;

K=-50:2:50;

w=-(0.1*K).^2/2;x=exp(w); X=fft(x,N); subplot(2,1,1);

stem(abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)'); subplot(2,1,2);

stem(angle(fftshift(X)));

ylabel('Phase'); xlabel('Frequency (rad)

%g[k]?x[4k]

N=500; K=-50:4:50; w=-(0.1*K).^2/2; x=exp(w); X=fft(x,N); subplot(2,1,1);

stem(abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)'); subplot(2,1,2);

stem(angle(fftshift(X)));

ylabel('Phase'); xlabel('Frequency (rad)');

(2) 将上述x[k]乘以cos(pk/2) %y[k]?x[2k]

N=500; K=-50:2:50; w=-(0.1*K).^2/2;

x=exp(w).*cos(pi*K/2); X=fft(x,N); subplot(2,1,1);

stem(abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)'); subplot(2,1,2);

stem(angle(fftshift(X)));

ylabel('Phase'); xlabel('Frequency (rad)');

%g[k]?x[4k]

N=500; K=-50:4:50; w=-(0.1*K).^2/2;

x=exp(w).*cos(pi*K/2); X=fft(x,N); subplot(2,1,1);

stem(abs(fftshift(X)));

ylabel('Magnitude'); xlabel('Frequency (rad)'); subplot(2,1,2);

stem(angle(fftshift(X)));

ylabel('Phase'); xlabel('Frequency (rad)');

四、实验思考题

1. 既然可直接由DTFT定义计算序列DTFT，为何利用DFT分析序列的频谱？

答：DFT用于分析有限长序列更方便

2. 若序列持续时间无限长，且无解析表达式，如何利用DFT分析其频谱？ 答：利用窗函数

3. 在利用DFT分析离散信号频谱时，会出现哪些误差？如何克服或改善？

答：频率分不出，应该多抽取些

4. 在利用DFT分析离散信号频谱时，如何选择窗函数？

答：根据窗函数的不同特性

5. 序列补零和增加序列长度都可以提高频谱分辨率吗？ 两者有何本质区别？

答：可以。补零是DFT中经常采用的方法，补零只能使序列的频谱变得细致,但不能提高序列的频率分辨率,只有采集更多的有效数据,才能得到序列的高分辨率频谱。

本文来源：https://www.bwwdw.com/article/3gf.html