数字信号处理实验报告王炜0518 湖南工业大学.docx
- 文档编号:25796881
- 上传时间:2023-06-14
- 格式:DOCX
- 页数:36
- 大小:564.89KB
数字信号处理实验报告王炜0518 湖南工业大学.docx
《数字信号处理实验报告王炜0518 湖南工业大学.docx》由会员分享,可在线阅读,更多相关《数字信号处理实验报告王炜0518 湖南工业大学.docx(36页珍藏版)》请在冰豆网上搜索。
数字信号处理实验报告王炜0518湖南工业大学
电气与信息工程学院
数字信号处理实验报告
老师:
舒小华
姓名:
王炜
学号:
12401720207
班级:
电子信息1202
时间:
2015年5月
实验一快速傅里叶变换及其应用
一、实验目的
(1)在理论学习的基础上,通过本实验,加强对FFT的理解,熟悉MATLAB中有关函数。
(2)应用FFT对典型信号进行频谱分析。
(3)了解应用FFT进行信号频谱分析过程中可能出现的问题,以便在实际中正确应用FFT.
(4)应用相关FFT实现序列的线性卷积和相关。
二、实验内容
实验用到的信号序列:
高斯序列
衰减正弦序列
三角波序列
反三角波序列
上机实验内容:
(1)观察高斯序列的时域和幅频特性,固定信号中参数p=8,改变q的值,使q分别等于2、4、8,观察他们的时域和幅频特性,了解当q去不同值时,对信号序列的时域和频域特性的影响;固定q=8,改变p,使p分别为等于8、13、14,观察参数p变化对信号序列的时域以及频域特性的影响,注意p等于多少时,会发生明显的泄露现象,混叠是否也随之出现?
记录实验中观察到的现象,绘制相应的时域序列和幅频特性曲线。
代码如下:
clearall;
%p=8,q=2%
[Xa1,Fa1]=gauss(8,2);
k=0:
15;
subplot(5,2,1);
plot(k,Xa1);
Xlabel('n');
Ylabel('时域特性');
text(10,0.5,'p=8,q=2');
subplot(5,2,2);
plot(k,Fa1);
Xlabel('n');
Ylabel('幅频特性');
text(8,3,'p=8,q=2');
%p=8,q=4%
[Xa1,Fa1]=gauss(8,4);
k=0:
15;
subplot(5,2,3);
plot(k,Xa1);
Xlabel('n');
Ylabel('时域特性');
text(10,0.5,'p=8,q=4');
subplot(5,2,4);
plot(k,Fa1);
Xlabel('n');
Ylabel('幅频特性');
text(8,3,'p=8,q=4');
%p=8,q=8%
[Xa1,Fa1]=gauss(8,8);
k=0:
15;
subplot(5,2,5);
plot(k,Xa1);
Xlabel('n');
Ylabel('时域特性');
text(10,0.5,'p=8,q=8');
subplot(5,2,6);
plot(k,Fa1);
Xlabel('n');
Ylabel('幅频特性');
text(8,3,'p=8,q=8');
%p=8,q=8%
[Xa1,Fa1]=gauss(13,8);
k=0:
15;
subplot(5,2,7);
plot(k,Xa1);
Xlabel('n');
Ylabel('时域特性');
text(10,0.5,'p=13,q=8');
subplot(5,2,8);
plot(k,Fa1);
Xlabel('n');
Ylabel('幅频特性');
text(8,3,'p=13,q=8');
%p=8,q=8%
[Xa1,Fa1]=gauss(14,8);
k=0:
15;
subplot(5,2,9);
plot(k,Xa1);
Xlabel('n');
Ylabel('时域特性');
text(10,0.5,'p=14,q=8');
subplot(5,2,10);
plot(k,Fa1);
Xlabel('n');
Ylabel('幅频特性');
text(8,3,'p=14,q=8');
(2)观察衰减正弦序列的时域和幅频特性,a=0.1,f=0.0625,检查峰值出现的位置是否正确,注意频谱的形状,绘出频谱特性曲线,改变发,使分别等于0.4375和0.5625,观察这两种情况下,频谱的形状和频谱出现的位置,有无混叠、泄露现象?
说明产生现象的原因。
代码如下:
clearall;
a=0.1;
f=0.0625;
i=0:
160
x=exp(-a*(i)).*sin(2*pi*f*(i));
subplot(3,2,1);
plot(i,x);
xlabel('n');
title('a=0.1,f=0.0625时域特性')
subplot(3,2,2);
G=fft(x);
plot(i,abs(G));
xlabel('k');
title('a=0.1,f=0.0625幅频特性')
f=0.4375;
x=exp(-a*(i)).*sin(2*pi*f*(i));
subplot(3,2,3);
plot(i,x);
xlabel('n');
title('a=0.1,f=0.4375时域特性')
subplot(3,2,4);
G=fft(x);
plot(i,abs(G));
xlabel('k');
title('a=0.1,f=0.4375幅频特性')
%f=0.5625;
f=0.5625;
x=exp(-a*(i)).*sin(2*pi*f*(i));
subplot(3,2,5);
plot(i,x);
xlabel('n');
title('a=0.1,f=0.5625时域特性')
subplot(3,2,6);
G=fft(x);
plot(i,abs(G));
xlabel('k');
title('a=0.1,f=0.5625幅频特性')
(3)观察三角波和反三角波序列的时域和幅频特性,用N=8点FFT分析信号序列xc(n)和xd(n)的幅频特性,观察两者的序列形状和频谱曲线有什么异同?
绘出两序列及其幅频特性曲线。
在xc(n)和xd(n)末尾补零,用N=16点FFT分析这两个信号的幅频特性,观察幅频特性发生了什么变化?
两情况的FFT频谱还有相同之处吗?
这些变化说明了什么?
N=8代码如下:
clearall;
n=0:
3;
k=0:
7;
x1(n+1)=n;
x1(n+5)=8-(n+4);
x2(n+1)=4-n;
x2(n+5)=(n+4)-4;
F1=fft(x1,8);
F2=fft(x2,8);
subplot(2,2,1);
plot(k,x1);
xlabel('n');ylabel('时域');
gridon;
subplot(2,2,2);
plot(abs(F1));
xlabel('k');ylabel('幅频');
gridon;
subplot(2,2,3);
plot(k,x2);
xlabel('n');ylabel('时域');
gridon;
subplot(2,2,4);
plot(abs(F2));
xlabel('k');ylabel('幅频');
gridon;
N=32代码如下:
clearall;
n=0:
3;
k=0:
31;
x1(n+1)=n;
x1(n+5)=8-(n+4);
x2(n+1)=4-n;
x2(n+5)=(n+4)-4;
fori=9:
32
x1(i)=0;
end
fori=9:
32
x2(i)=0;
end
F1=fft(x1,32);
F2=fft(x2,32);
subplot(2,2,1);
plot(k,x1);
xlabel('n');ylabel('时域');
gridon;
subplot(2,2,2);
plot(abs(F1));
xlabel('k');ylabel('幅频');
gridon;
subplot(2,2,3);
plot(k,x2);
xlabel('n');ylabel('时域');
gridon;
subplot(2,2,4);
plot(abs(F2));
xlabel('k');ylabel('幅频');
gridon;
(4)一个连续信号含两个频率分量,经采样得x(n)=sin(2π*0.125n)+cos(2π*(0.125+Δf)n)n=0,1……,N-1已知N=16,Δf分别为1/16和1/64,观察其频谱;当N=128时,Δf不变,其结果有何不同,为什么?
clearall;
n=0:
15;
f1=1/16;f2=1/64;
x=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f1)*n);
subplot(4,2,1);
plot(n,x);
title('f=1/16;N=16时域');
subplot(4,2,2);
y=fft(x);y=abs(y);
stem(n,y);
title('f=1/16;N=16幅频');
x1=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f2)*n);
subplot(4,2,3);
plot(n,x1);
title('f=1/64;N=16时域');
subplot(4,2,4);
y1=fft(x1);y1=abs(y1);
stem(n,y1);
title('f=1/64;N=16幅频');
n=0:
127;f1=1/16;f2=1/64;
x=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f1)*n);
subplot(4,2,5);
plot(n,x);
title('f=1/16;N=128时域');
subplot(4,2,6);
y=fft(x);y=abs(y);
stem(n,y);
title('f=1/16;N=128幅频');
x1=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f2)*n);
subplot(4,2,7);
plot(n,x1);
title('f=1/16;N=128时域');
subplot(4,2,8);
y1=fft(x1);y1=abs(y1);
stem(n,y1);
title('f=1/16;N=128幅频');
(5)用FFT分别计算(p=8,q=2)和(a=0.1,f=0.0625)的16点循环卷积和线性卷积。
相关代码如下:
clear,closeall
n=0:
15;
%p=8,q=2
xa=exp(-((n-8).^2/2));
%a=0.1,f=0.625
xb=exp((-0.1).*n).*sin(2*pi*0.0625*n);
L=max(length(xa),length(xb));
N=length(xa)+length(xb)-1;
%conv方法
y1=conv(xa,xb);
s=loopfunc(xa,xb,L);
subplot(2,2,1),stem(s);xlabel('循环卷积conv')
subplot(2,2,2),stem(y1);xlabel('线性卷积conv')
%ifft方法
ya=fft(xa);
yb=fft(xb);
y1=ya.*yb;
%ifft方法
y11=ifft(y1);
subplot(2,2,3);
stem(n,y11);xlabel('循环卷积ifft');
yaa=fft(xa,32);
ybb=fft(xb,32);
y2=yaa.*ybb;
%ifft方法
y22=ifft(y2);
subplot(2,2,4);
n=0:
31;stem(n,y22);xlabel('线性卷积ifft');
三、思考题
(1)实验中的信号序列和,在单位圆上的Z变换频谱和会相同吗?
如果不同,你能说出哪一个低频分量更多一些吗?
为什么?
答:
不同。
反三角波的低频分量更多一些。
(2)对一个有限长序列进行DFT等价于将该序列周期延拓后进行DFS展开,因为DFS也只是取其中一个周期来计算,所以FFT在一定条件下也可以用以分析周期信号序列。
如果实正弦信号sin(2πfn),f=0.1用16点FFT来做DFS运算,得到的频谱时信号本身的真实谱吗?
答:
可以把有限长非周期序列假设为一无限长周期序列的一个主直周期,即对有限长非周期序列进行周期延拓,延拓后的序列完全可以采用DFS进行处理,即采用复指数。
实验二IIR数字滤波器的设计
一、实验目的
1.掌握双线性变换法及脉冲响应不变法设计IIR数字滤波器的具体设计方法及其原理,熟悉用双线性变换法及脉冲响应不变法设计低通,高通,带通数字滤波器的急速演技编程。
2.观察双线性变换及脉冲响应不变法设计的滤波器的频域特性,了解双线性变换法及脉冲响应不变法的特点。
3.熟悉巴特沃思滤波器,切比雪夫滤波器和椭圆滤波器的频域特性。
二、实验内容
试验中有关的变量定义:
通带边界频率
阻带边界频率
通带波动
最小阻带衰减
采样频率
T采样周期
上机实验内容:
1.=0.3KHz,=0.8dB,=0.2KHz,=20dB,T=1ms;设计一切比雪夫高通滤波器,观察其通带损耗和阻带衰减是否满足要求。
具体程序如下:
%
(1)切比雪夫高通滤波器
fc=300;%fc=0.3kHz
D=0.8;%0.8dB
fr=200;%fr=0.2kHz
At=20;%20dB
fs=1000;%T=1ms
figure('NumberTitle','off','Name','
(1)切比雪夫高通滤波器');
%CHEBY1(N,R,Wp,'high','s')
%[N,Wp]=CHEB1ORD(Wp,Ws,Rp,Rs,'s')
%[NUMd,DENd]=BILINEAR(NUM,DEN,Fs),
wc=2*fs*tan(2*pi*fc/(2*fs));%omega=(2/T)*tg(w/2)w=2*pi*(f/fs)
wt=2*fs*tan(2*pi*fr/(2*fs));
[N,wn]=cheb1ord(wc,wt,D,At,'s');%N:
返回额度阶数ord:
阶wn:
通带截止频率s:
模拟
[B,A]=cheby1(N,D,wn,'high','s');%B:
分子系数A:
分母系数
[num,den]=bilinear(B,A,fs);%num:
系统函数分子系数den:
系统函数分母系数(值从workspace中得到)
[h,w]=freqz(num,den);%z:
数字
f=w/pi*(fs/2);
plot(f,20*log10(abs(h)));
axis([0,500,-80,10]);%坐标轴范围
grid;
xlabel('频率/Hz')
ylabel('幅度/dB')
=0.2KHz,=1dB,=0.3KHz,=25dB,T=1ms;分别用脉冲响应不变法及双线性变换法设计一巴特沃斯数字低通滤波器,观察所设计数字滤波器的幅频特性曲线,记录带宽和衰减量,检查是否满足要求。
比较这两种方法的优缺点。
clc;clear;
%参数
fs=1000;%T=1msfs=1/T=1kHz
fc=200;%fc=0.3kHz;
fr=300;%fr=0.2kHz
T=0.001;%T=1ms
figure('NumberTitle','off','Name','
(2)巴特沃思低通滤波器');
wp1=2*pi*fc;
wr1=2*pi*fr;
[N1,wn1]=buttord(wp1,wr1,1,25,'s');
[B1,A1]=butter(N1,wn1,'s');
[num1,den1]=impinvar(B1,A1,fs);%脉冲响应不变法
[h1,w]=freqz(num1,den1);
wp2=2*fs*tan(2*pi*fc/(2*fs));
wr2=2*fs*tan(2*pi*fr/(2*fs));
[N2,wn2]=buttord(wp2,wr2,1,25,'s');
[B2,A2]=butter(N2,wn2,'s');
[num2,den2]=bilinear(B2,A2,fs);%双线性变换法
[h2,w]=freqz(num2,den2);
f=w/(2*pi)*fs;
%前者为脉冲响应不变法曲线后者为双线性变换法曲线
plot(f,20*log10(abs(h1)),'-.',f,20*log10(abs(h2)),'-');
axis([0,500,-200,10]);grid;
xlabel('频率/hz');
ylabel('幅度/dB');
title('巴特沃思低通滤波器');
legend('脉冲响应不变法','双线性变换法');
利用双线性变换法分别设计满足下列指标的巴特沃斯型,切比雪夫型和椭圆型数字滤波器,并作图验证设计结果1.2KHz,dB,=2KHz,40dB,。
比较这三种滤波器的阶数。
具体的代码如下:
clc;clear;
fc=1200;%fc=1.2kHz
fr=2000;%fr=2kHz
rp=0.5;
rs=40;
fs=8000;%fs=8kHz
figure('NumberTitle','off','Name','(3)三种数字低通滤波器');
wc=2*pi*fc;
wr=2*pi*fr;
w1=2*fs*tan(wc/(2*fs));
w2=2*fs*tan(wr/(2*fs));
[Nb,wn]=buttord(w1,w2,rp,rs,'s');
[B,A]=butter(Nb,wn,'s');
[num1,den1]=bilinear(B,A,fs);
[h1,w]=freqz(num1,den1);
[Nc,wn]=cheb1ord(w1,w2,rp,rs,'s');%切比雪夫
[B,A]=cheby1(Nc,rp,wn,'s');
[num2,den2]=bilinear(B,A,fs);
[h2,w]=freqz(num2,den2);
[Ne,wn]=ellipord(w1,w2,rp,rs,'s');%椭圆型
[B,A]=ellip(Ne,rp,rs,wn,'low','s');
[num3,den3]=bilinear(B,A,fs);
[h3,w]=freqz(num3,den3);
f=w/pi*4000;
plot(f,20*log10(abs(h1)),'-',f,20*log10(abs(h2)),'--',f,20*log10(abs(h3)),':
');
axis([0,3000,-100,10]);grid;
%从左到右依次为:
'巴特沃思数字低通滤波器','切比雪夫数字低通滤波器','椭圆数字低通滤波器的曲线
xlabel('频率/Hz');ylabel('幅度/dB');title('三种数字低通滤波器');
legend('巴特沃思数字低通滤波器','切比雪夫数字低通滤波器','椭圆数字低通滤波器');
分别用脉冲响应不变法及双线性变换法设计一巴特沃斯性还型数字滤波器,已知=30KHz,其等效的模拟滤波器指示为3dB,2KHzf3KHz;5dB,f6KHz;20dB,f1.5KHz。
具体程序如下:
clc;clear;
figure('NumberTitle','off','Name','(4)巴特沃思数字低通滤波器');
wc=[2*pi*2000,2*pi*3000];
wr=[2*pi*1500,2*pi*6000];
rp=3;%3dB
rs=20;%20dB
fs=30000;%fs=3kHz
[N,wn]=buttord(wc,wr,rp,rs,'s');
[B,A]=butter(N,wn,'s');
[num1,den1]=impinvar(B,A,fs);%脉冲响应不变法
[h1,w]=freqz(num1,den1);
w1=2*fs*tan(2*pi*2000/(2*fs));
w2=2*fs*tan(2*pi*3000/(2*fs));
wr1=2*fs*tan(2*pi*1500/(2*fs));
wr2=2*fs*tan(2*pi*6000/(2*fs));
[N,wn]=butter(N,wn,'s');
[num2,den2]=bilinear(B,A,fs);%双线性变换法
[h2,w]=freqz(num2,den2);
f=w/pi*15000;
plot(f,20*log10(abs(h1)),'-.',f,20*log10(abs(h2)),'-');
axis([500,7000,-30,10]);grid;
xlabel('频率/Hz');ylabel('幅度/dB');title('巴特沃思数字低通滤波器');
legend('脉冲响应不变法','双线性变换法',1);
利用双线性变换法设计满足下列指标的切比雪夫型数字带阻滤波器,并作图验证设计结果:
当1KHzf2KHz时,18dB;当f500Hz以及f3KHz时,3dB,采样频率=10KHz。
具体程序如下:
clc;clear;
%1kHz<=f<=2kHzAt>=18dB
%f<=500Hzf>=3kHzd<=3dB
%fs=10kHz
d=3;
At=18;
fs=10000;
fc1=500;
fc2=3000;
fr1=1000;
fr2=2000;
wc1=2*10000*tan(2*pi*fc1/(2*fs));
wc2=2*10000*tan(2*pi*fc2/(2*fs));
wr1=2*10000*tan(2*pi*fr1/(2*fs));
wr2=2*10000*tan(2*pi*fr2/(2*fs));
figure('NumberTitle','off','Name','(5)切比雪夫数字带阻滤波器');
[N,wn]=cheb1ord([wc1,wc2],[wr1,wr2],d,At,'s');
[B,A]=cheby1(N,d,wn,'stop','s');
[num,den]=bilinear(B,A,8000);
[h,w]=freqz(num,den);
f=w/pi*fs;
plot(f,20*log10(abs(h)));
axis([0,8000,-150,10]);grid;
xlabel('频率/Hz');ylabel('幅度/dB');
%系统函数参数:
%num:
[0.190707285833214,-0.351720570024798,0.543583729019387,-0.35172057002
%4798,0.190707285833215;]
%den:
[1,-0.745790514400847,-0.153********3616,-0.247846506351433,0.46001707
%9664572;]
三、思考题
(1)双线性变换法中Ω和之间的关系是非线性的,在实验中你注意到这种非线性关系了吗?
从哪几种数字滤波器的幅频特性曲线中可以观察到这种非线性关系?
答:
在
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号处理实验报告王炜0518 湖南工业大学 数字信号 处理 实验 报告 0518 湖南 工业大学