西安交通大学实验报告
第 页 共 页
课 程 数字信号处理
系 别 实 验 日 期 年 月 日
专业班级 组别 交报告 日 期 年 月 日
姓 名 学号 报 告 退 发 (订正、重做)
同 组 者 教师审批签字
实验名称 用DFT计算连续信号的谱
一、实验目的:
1.学习用Matlab完成信号的DFT计算。
2.了解用DFT对连续周期信号进行谱分析可能产生的误差。
二、实验内容及要求:
1. 有一单频信号x(t)=sin(2πft)
若 ①f=f1=15Hz ②f=f2=2Hz
分别用DFT求x(t)的谱。
要求:①抽样频率fs以不发生混叠为宜。
②记录长度tp应取整数周期。
③用stem语句绘出幅度谱,横坐标为模拟频率f。用plot画出时域波形。
2. 有一复频信号x(t)=sin(2πf1t)+ sin(2πf2t)
其中 f1=15Hz ,f2=2Hz
当 T=0.01s时,求:
① N=100时x(t)的幅度谱并图示;(用stem语句绘出幅度谱)
② N=50时x(t)的幅度谱并图示。(用stem语句绘出幅度谱)
请分别预测①与②是否逼近真实的谱?如有误差,请分析原因。
三、实验内容程序及结果
1、运行程序如下:
n=0:5;
n1=0:0.01:6;
n1=n1./45;
fs1=45;T1=1/fs1;
x=sin(2*pi*15*n.*T1);
x1=sin(2*pi*15*n1);
subplot(221);plot(n1,x1)
xlabel('时间');ylabel('信号');
y=abs(fft(x));
y=y(1:4);
f=(0:3)*15/2;
subplot(222);stem(f,y)
xlabel('频率(Hz)');ylabel('幅值');
n=0:5;
n1=(0:0.01:6)/6;
fs2=6;T2=1/fs2;
x=sin(2*pi*2*n.*T2);
x1=sin(2*pi*2*n1);
subplot(223);plot(n1,x1)
xlabel('时间');ylabel('信号');
y=abs(fft(x));
y=y(1:4);
f=(0:3)*2/2;
subplot(224);stem(f,y)
xlabel('频率(Hz)');ylabel('幅值');
运行结果如下图所示:
实验结果分析:
①要使得在抽样过程中不发生混叠现象,那么根据抽样定理,抽样频率fs应大于2*f,故取fs=3*f,即fs1=45 Hz ,fs2=6 Hz,这样抽样后就不会发生混叠的问题。
② 记录长度tp取了2个整周期,因取抽样频率为原始频率的3倍,因此每个周期应该含有3个点,2个周期含有6个点,故设置n=0:5。如果tp不为整数个周期,则会出现频谱泄漏现象。
③ 要求横坐标为模拟频率f,因为由DFT运算出来的结果为时域、频域都是离散周期的序列,所以要求把频域横坐标化为模拟频率。设计算出来的模拟频率里每个样点间的间隔为f0,取的样点个数为N, 则有:2*pi*f/fs=2*pi*k/N,
可得:f= fs* k/N,其中N=6,fs=3f0,所以f=k* f0/2
2、运行程序如下:
N1=0:99;
N2=0:49;
T=0.01;
f1=15;
f2=2;
x=sin(2*pi*f1*T*N1)+sin(2*pi*f2*T*N1);
y=fft(x);
y=abs(y);
subplot(121);stem(N1,y)
xlabel('频率');ylabel('幅值');
x=sin(2*pi*f1*T*N2)+sin(2*pi*f2*T*N2);
y=fft(x);
y=abs(y);
subplot(122);stem(N2,y)
xlabel('频率');ylabel('幅值')
运行结果如下图所示:
问题:请分别预测①与②是否逼近真实的谱?如有误差,请分析原因。
答:由已知的时域t轴上的最小间隔T和采样点数N,即可确定抽样频率和记录长度。由实验结果可得,①与②之间有较大的区别,其中①更符合真真实的谱子。因为当N=100时,取的为复频信号的整数个周期,而当N=50时不为整数周期,因此会产生频谱泄漏现象,所以①的结果逼近真实的谱,而②则不逼近真实情况。
实验结果分析:这道题目很直观的表现出频谱泄漏的现象,让我明白了这种现象出现的原因,和应该避免的措施。
四、实验心得体会
通过本次实验,我熟悉了进行DFT运算中的常用指令,掌握了运用MATLAB语言快速实现信号的DFT运算的命令程序,也更加深入的领会了频谱泄漏现象的表现。处理信号转换中繁琐的过程通过Matlab的简便语句就能得出完美的结果,从其中可以很明显地感觉到Matlab的强大之处,也更加促使我进一步掌握Matlab的兴趣和信心。
第二篇:数字信号处理实验报告
数字信号处理实验报告
(基于MATLAB的DFT频谱分析)
学 院:通信与信息工程学院
班 级:通信工程1104班
姓 名:杨 文 刚
学 号:1107020420
20##年11月20日
一、 实验目的
1.掌握DFT函数的用法。
2. 利用DFT进行信号检测及谱分析。
3.了解信号截取长度对谱分析的影响。
二、 实验内容
1.利用DFT计算信号功率谱。
实验程序:
t=0:0.001:0.6;
x=sin(2*pi*50*t)+sin(2*pi*120*t)+randn(1,length(t));
Y=dft(x,512);
P=Y.*conj(Y)/512;
f=1000*(0:255)/512;
plot(f,P(1:256))
运行截图
2. 进行信号检测,分析信号频谱所对应频率轴的数字频率和频率之间的关系。
模拟信号,以 进行取样,求N点DFT的幅值谱。
实验程序:
subplot(2,2,1)
N=45;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
y=dft(x,N);
plot(q,abs(y));title('DFT N=45')
subplot(2,2,2)
N=50;n=0:N-1;t=0.01*n; q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
y=dft(x,N);
plot(q,abs(y));title('DFT N=50')
subplot(2,2,3)
N=55;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
y=dft(x,N);
plot(q,abs(y));title('DFT N=55')
subplot(2,2,4)
N=60;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
y=dft(x,N);
plot(q,abs(y));title('DFT N=60')
运行截图
3. 对2,进一步增加截取长度和DFT点数,如N加大到256,观察信号频谱的变化,分析产生这一变化的原因。在截取长度不变的条件下改变采样频率,观察信号频谱的变化,分析产生这一变化的原因。
N加大到256时的程序:
N=256;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
y=dft(x,N);
plot(q,abs(y));title('DFT N=256')
运行截图
分析原因:在T=0.01s的情况下,第一个序列的周期是100,第二个序列的周期是50,所以当取样点数小于100时,频率分辨率不够,不能够区分出两个信号。当采样点数足够多(256)时,频率分辨率增加,能够区分出两个频率的信号。
将采样间隔变为T=0.1s时,N仍为45的程序:
N=45;n=0:N-1;t=0.1*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
y=dft(x,N);
plot(q,abs(y));title('DFT N=45')
运行图如右图
分析原因:在T=0. 1s的情况下,第一个序列的周期是10,第二个序列的周期是5,所以当取样点数为45时,能够区分出两个信号。
参数同上,N取64,并在信号中加入噪声w(t)。
figure(2)
subplot(2,1,1)
N=64;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t); y=dft(x,N);
plot(q,abs(y));title('DFT N=64')
subplot(2,1,2)
N=64;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t)+0.8*randn(1,N); y=dft(x,N);
plot(q,abs(y));title('DFT N=64(with noise)')
由图可以看出这种噪音不影响信号检测。
4. 对3,加大噪声到2*randn(1,N)和8*randn(1,N),画出并比较不同噪声下时域波形和频谱。
subplot(2,1,1)
N=64;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t)+2*randn(1,N); y=dft(x,N);
plot(q,abs(y));title('DFT N=64(with noise2)')
subplot(2,1,2)
N=64;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t)+8*randn(1,N); y=dft(x,N);
plot(q,abs(y));title('DFT N=64(with noise8)')
subplot(3,2,1)
N=64;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t)+0.8*randn(1,N);
plot(x);title('噪声为0.8*w的信号')
y=dft(x,N);
subplot(3,2,2)
plot(q,abs(y));title('噪声为0.8*w时的频谱')
subplot(3,2,3)
N=64;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t)+2*randn(1,N);
plot(x);title('噪声为2*w时的信号')
y=dft(x,N);
subplot(3,2,4)
plot(q,abs(y));title('噪声为2*w时的频谱')
subplot(3,2,5)
N=64;n=0:N-1;t=0.01*n;
q=n*2*pi/N;
x=2*sin(4*pi*t)+5*cos(8*pi*t)+8*randn(1,N);
plot(x);title('噪声为8*w时的信号')
y=dft(x,N);
subplot(3,2,6)
plot(q,abs(y));title('噪声为8*w时的频谱')
实验分析:当噪声较小时,不影响信号的检测,但当噪声较大时,就看不出原信号的频率成分了,可以继续加大噪声,可看到其频谱杂乱无章了。
5. 用一个N点DFT计算两个长度为N的实序列N点离散傅里叶变换,并将结果和直接使用两个N点DFT得到的结果进行比较。
x=[1 2 3 4 5 6];
y=[6 5 4 3 2 1];
[a,b]=dft_2(x,y)
a =
Columns 1 through 3
21.0000 -3.0000 + 5.1962i -3.0000 + 1.7321i
Columns 4 through 6
-3.0000 -3.0000 - 1.7321i -3.0000 - 5.1962i
b =
Columns 1 through 3
21.0000 3.0000 - 5.1962i 3.0000 - 1.7321i
Columns 4 through 6
3.0000 3.0000 + 1.7321i 3.0000 + 5.1962i
函数文件如下:
function [y1,y2]=dft_2(a,b)
N=length(a);
x=zeros(1,N);
x=a+j*b;
X=dft(x,N);
X0=conj(fliplr(X));
X0=[X0(N) X0(1:N-1)];
y1=(X+X0)./2;
y2=(X-X0)./2./j;
直接运行计算 :
dft(x)
ans =
Columns 1 through 3
21.0000 -3.0000 + 5.1962i -3.0000 + 1.7321i
Columns 4 through 6
-3.0000 -3.0000 - 1.7321i -3.0000 - 5.1962i
dft(y)
ans =
Columns 1 through 3
21.0000 3.0000 - 5.1962i 3.0000 - 1.7321i
Columns 4 through 6
3.0000 3.0000 + 1.7321i 3.0000 + 5.1962i
6.比较DFT和DFT的运算时间。(计时函数 tic, toc)
N分别取256,512,1024,2048,4096,…
程序如下:
N=256;
N=4096;
x=randn(1,N);
tic
y=dft(x,N);
toc
tic
z=dft(x);
toc
N=256:Elapsed time is 0.172000 seconds.
Elapsed time is 0.015000 seconds.
N=512:Elapsed time is 0.687000 seconds.
Elapsed time is 0.000000 seconds.
N=1024:Elapsed time is 3.031000 seconds.
Elapsed time is 0.047000 seconds.
N=2048:Elapsed time is 13.375000 seconds.
Elapsed time is 0.063000 seconds.
N=4096:Elapsed time is 59.250000 seconds.
Elapsed time is 0.125000 seconds.
7.对给定语音信号进行谱分析,写出采样频率,画出语音信号的波形及频谱,并分析语音信号的频率分布特点。
(1)画时域波形并对整个语音序列做DFT
[x,fs]=wavread('C:\ai1.wav');
subplot(2,1,1)
N=length(x);
n=0:N-1;
plot(n,x);
xlabel('n');
ylabel('x');
title('时域波形');
subplot(2,1,2);
N=length(x);
n=0:N-1;
t=0.01*n;
q=n*2*pi/N;
y=dft(x,N);
plot(q,abs(y));
xlabel('n');
ylabel('ai1');
title('DFT');
(2)并分别求出k=300,3500所对应的信号频率(Hz)
[x,fs]=wavread('C:\ai1.wav')
N=length(x);
n=0:N-1;
t=n*(1/fs);
q=n*2*pi/N;
n1=300;
q1=n1*2*pi/N;
f1=q1*fs/(2*pi);
n2=3500;
q2=n2*2*pi/N;
f2=q2*fs/(2*pi);
fs =
16000
(3)从语音中截取一段语音(256点)做DFT,得出频谱,画出时域波形及频谱。
程序如下:
[x,fs]=wavread('C:\ai1.wav');
subplot(2,1,1);
N=256;
n=0:N-1;
x=x(1:256);
plot(n,x);
xlabel('n');
ylabel('x');
title('256点时域波形');
subplot(2,1,2);
N=256;
n=0:N-1;
t=0.01*n;
q=n*2*pi/N;
x=x(1:256);
y=dft(x,N);
plot(q,abs(y));
xlabel('n');
ylabel('ai1');
title('DFT-256');
三、 实验心得
本次数字信号处理实验,我使用matlab工具做了对于一些DFT信号的频谱分析,由于我们并没有学习过matlab的使用方法,对这个工具软件的使用我并不熟悉。我在网上查找了用matlab做DFT频谱分析的一些相关资料,去图书馆借阅了相关的书籍,作为参考,做了这次实验。通过这次实验,首先我对matlab这个软件有了初步的了解,同时也感受到了这个软件的强大,是一个很好的数学工具。同时,我也对离散傅里叶变换有了更进一步的认识。由于这学期我们的课程是英文版的,对我们平时的学习带来了一定的困难,本次实验我认为是一个很好的机会去实际动手利用工具来操作和实践,认识和理解数字信号处理的一些方法。