数字信号处理实验报告.docx
《数字信号处理实验报告.docx》由会员分享,可在线阅读,更多相关《数字信号处理实验报告.docx(34页珍藏版)》请在冰豆网上搜索。
数字信号处理实验报告
南京邮电大学
实验报告
实验名称熟悉MATLAB环境
快速傅里叶变换(FFT)及其应用
IIR数字滤波器的设计
FIR数字滤波器的设计
课程名称数字信号处理
班级学号99(B130904)
姓名许永耀
开课时间2014/2015学年,第二学期
实验一熟悉MATLAB环境
一、实验目的
(1)熟悉MATLAB的主要操作命令。
(2)学会简单的矩阵输入和数据读写。
(3)掌握简单的绘图命令。
(4)用MATLAB编程并学会创建函数。
(5)观察离散系统的频率响应。
二、实验内容
(1)数组的加、减、乘、除和乘方运算。
输入A=[1234],B=[3,4,5,6],求
C=A+B,D=A-B,E=A.*B,F=A./B,G=A.^B。
并用stem语句画出A、B、C、D、
E、F、G。
n=0:
1:
3;
A=[1234];
subplot(7,1,1)
stem(n,A)
xlabel('n')
ylabel('A')
B=[3,4,5,6];
subplot(7,1,2)
stem(n,B)
xlabel('n')
ylabel('B')
C=A+B;
subplot(7,1,3)
stem(n,C)
xlabel('n')
ylabel('C')
figure
C=A+B
D=A-B
E=A.*B
F=A./B
G=A.^B
C=
46810
D=
-2-2-2-2
E=
381524
F=
0.33330.50000.60000.6667
G=
1162434096
n=0:
1:
3;
A=[1234];
B=[3,4,5,6];
D=A-B;
subplot(4,1,1)
stem(n,D)
xlabel('n')
ylabel('D')
E=A.*B;
subplot(4,1,2)
stem(n,E)
xlabel('n')
ylabel('E')
F=A./B;
subplot(4,13)
stem(n,F)
xlabel('n')
ylabel('F')
G=A.^B;
subplot(4,1,4)
stem(n,G)
xlabel('n')
ylabel('G')
(2)用MATLAB实现下列序列:
a)
n=0:
1:
15;
x1=0.8.^n;
stem(n,x1)
xlabel('n')
ylabel('x(n)')
title('2(a)')
b)
n=0:
1:
15;
i=sqrt(-1);
a=0.2+3*i;
x2=exp(a*n);
figure
subplot(1,2,1)
stem(n,real(x2))
xlabel('n')
ylabel('x(n)实部')
subplot(1,2,2)
stem(n,imag(x2))
xlabel('n')
ylabel('x(n)虚部')
c)
n=0:
1:
15;
x3=3*cos(0.125*pi*n+0.2*pi)+2*sin(0.25*pi*n+0.1*pi);
stem(n,x3)
xlabel('n')
ylabel('x(n)')
(4)绘出下列时间函数的图形,对x轴、y轴以及图形上方均须加上适当的标注:
a)
t=0:
0.001:
10;
x=sin(2*pi*t);
plot(t,x,'r-')
xlabel('t'),ylabel('x(t)'),title('sin(2\pit)')
b)
t=0:
0.001:
4;
x=cos(100*pi*t).*sin(pi*t);
plot(t,x,'r-')
xlabel('t'),ylabel('x(t)'),title('cos(100*pi*t).*sin(pi*t)')
(6)给定一因果系统
,求出并绘制H(z)的幅频响应和相频响应。
b=[1,sqrt
(2),1];a=[1,-0.67,0.9];
[H,w]=freqz(b,a,512,'whole');magH=abs(H);phaH=angle(H);
subplot(2,1,1),plot(w/pi,magH);grid
xlabel('');ylabel('幅度');title('幅频响应')
subplot(2,1,2),plot(w/pi,phaH/pi);grid
xlabel('频率(单位:
pi)');ylabel('相位(单位:
pi)');title('相频响应')
(7)计算序列{8-2-123}和序列{23-1-3}的离散卷积,并作图表示卷积结果。
n=[0:
8-1];
x=[8,-2,-1,2,3];h=[2,3,-1,-3];
y=conv(x,h);stem(n,y);
xlabel('n')
ylabel('y');
figure
y=conv(x,h)
y=
1620-16-211910-9-9
(8)求以下差分方程所描述系统的单位脉冲响应h(n),
a=[1,0.1,-0.06];
b=[1,-2];
N=50;
x=[1,zeros(1,N-1)];
h=filter(b,a,x);
n=0:
1:
50;
stem(h,'.');
每一小题均给出实验过程与结果(含实验程序、运行的数据结果和图形);
实验二快速傅里叶变换(FFT)及其应用
一、实验目的
(1)在理论学习的基础上,通过本实验,加深对FFT的理解,熟悉MATLAB中的有关函数。
(2)应用FFT对典型信号进行频谱分析。
(3)了解应用FFT进行信号频谱分析过程中可能出现的问题,以便在实际中正确应用FFT。
(4)应用FFT实现序列的线性卷积和相关。
二、实验内容
实验中用到的信号序列
a) 高斯序列
b) 衰减正弦序列
c) 三角波序列
d) 反三角波序列
(1)观察高斯序列的时域和幅频特性,固定信号
中参数p=8,改变q的值,使q分别等于2,4,8,观察它们的时域和幅频特性,了解当q取不同值时,对信号序列的时域幅频特性的影响;固定q=8,改变p,使p分别等于8,13,14,观察参数p变化对信号序列的时域及幅频特性的影响,观察p等于多少时,会发生明显的泄漏现象,混叠是否也随之出现?
记录实验中观察到的现象,绘出相应的时域序列和幅频特性曲线。
n=0:
15;
p=8;q=2;
xa=exp(-(n-p).^2/q);
subplot(3,2,1);
stem(n,xa,'.');
title('q=2高斯序列时域特性');
subplot(3,2,2);stem(fft(xa),'.');
title('q=2高斯序列幅频特性');
n=0:
15;
p=8;q=4;
xa=exp(-(n-p).^2/q);
subplot(3,2,3);stem(n,xa,'.');
title('q=4高斯序列时域特性');
subplot(3,2,4);stem(fft(xa),'.');
title('q=4高斯序列幅频特性');
n=0:
15;p=8;q=8;
xa=exp(-(n-p).^2/q);
subplot(3,2,5);stem(n,xa,'.');
title('q=8高斯序列时域特性');
subplot(3,2,6);
stem(fft(xa),'.');
title('q=8高斯序列幅频特性');
n=0:
15;
p=8;q=8;
xa=exp(-(n-p).^2/q);
subplot(3,2,1);
stem(n,xa,'.');
title('p=8高斯序列时域特性');
subplot(3,2,2);stem(fft(xa),'.');
title('p=8高斯序列幅频特性');
n=0:
15;
p=13;q=8;
xa=exp(-(n-p).^2/q);
subplot(3,2,3);
stem(n,xa,'.');
title('p=13高斯序列时域特性');
subplot(3,2,4);
stem(fft(xa),'.');
title('p=13高斯序列幅频特性');
n=0:
15;
p=14;q=8;
xa=exp(-(n-p).^2/q);
subplot(3,2,5);
stem(n,xa,'.');
title('p=14高斯序列时域特性');
subplot(3,2,6);stem(fft(xa),'.');
title('p=14高斯序列幅频特性');
n=0:
15;
p=8;q=8;
xa=exp(-(n-p).^2/q);
subplot(3,2,1);
stem(n,xa,'.');
title('p=8高斯序列时域特性');
subplot(3,2,2);
plot(abs(fft(xa)));
title('p=8高斯序列幅频特性');
n=0:
15;
p=13;q=8;
xa=exp(-(n-p).^2/q);
subplot(3,2,3);
stem(n,xa,'.');
title('p=13高斯序列时域特性');
subplot(3,2,4);
plot(abs(fft(xa)));
title('p=13高斯序列幅频特性');
n=0:
15;
p=14;q=8;
xa=exp(-(n-p).^2/q);
subplot(3,2,5);
stem(n,xa,'.');
title('p=14高斯序列时域特性');
subplot(3,2,6);
plot(abs(fft(xa)));
title('p=14高斯序列幅频特性');
(3)观察三角波和反三角波序列的时域和幅频特性,用N=8点FFT分析信号序列
和
的幅频特性,观察两者的序列形状和频谱曲线有什么异同?
绘出两序列及其幅频特性曲线。
forn=1:
4
xc(n)=n-1;
end;
forn=5:
8;
xc(n)=8-(n-1);
end;
%xc=[0,1,2,3,4,3,2,1];
n=0:
1:
7;
subplot(2,1,1)
stem(n,xc)
N=8;%fft点数
%N=32;%%改变fft点数
Xc=fft(xc,N);
k=0:
1:
N-1;
subplot(2,1,2)
plot(k,abs(Xc))
forn=1:
4
xd(n)=5-n;
end;
forn=5:
8;
xd(n)=n-5;
end;
n=0:
1:
7;
subplot(2,1,1)
stem(n,xd)
N=8;
Xd=fft(xd,N);
k=0:
1:
N-1;
subplot(2,1,2)
plot(k,abs(Xd))
在
和
末尾补零,用N=32点FFT分析这两个信号的幅频特性,观察幅频特性发生了什么变化?
两种情况的FFT频谱还有相同之处吗?
这些变化说明了什么?
forn=1:
4
xc(n)=n-1;
end;
forn=5:
8;
xc(n)=8-(n-1);
end;
forn=9:
32
xc(n)=0;
end;
n=0:
1:
31;
subplot(2,1,1)
stem(n,xc)
N=32;
Xc=fft(xc,N);
k=0:
1:
N-1;
subplot(2,1,2)
plot(k,abs(Xc))
forn=1:
4
xd(n)=5-n;
end;
forn=5:
8;
xd(n)=n-5;
end;
forn=9:
32
xd(n)=0;
end;
n=0:
1:
31;
subplot(2,1,1)
stem(n,xd)
N=32;
Xd=fft(xd,N);
k=0:
1:
N-1;
subplot(2,1,2)
plot(k,abs(Xd))
(5)用FFT分别实现
(p=8,q=2)和
(a=0.1,f=0.0625)的16点循环卷积和线性卷积。
n1=0:
1:
15;p1=8;q1=2;
x=exp(-(n1-p1).^2/q1);
n2=0:
1:
15;a=0.1;f2=0.0625;
y=(exp(-a*n2)).*sin(2*pi*f2*n2);
N=length(x);n=0:
N-1;n3=0:
30;
X=fft(x);Y=fft(y);x32=[xzeros(1,16)];y32=[yzeros(1,16)];
X32=fft(x32);Y32=fft(y32);
z16=ifft(X.*Y);
z32=ifft(X32.*Y32);
subplot(2,2,1);plot(n,z16,'-*');xlabel('n');ylabel('z(n)');title('循环卷积结果');subplot(2,2,2);plot(n3,z32(1:
2*N-1),'-o');xlabel('n');
ylabel('z(n)');title('线性卷积结果');
rm16=real(ifft(conj(X).*Y));
rm32_0=real(ifft(conj(X32).*Y32));rm32=[rm32_0(N+2:
2*N)rm32_0(1:
N)];
m=n;subplot(2,2,3);plot(m,rm16,'--');xlabel('m');ylabel('rm');title('循环相关结果');m=-(N-1):
N-1;subplot(2,2,4);plot(m,rm32,'--o');xlabel('m');ylabel('rm');title('线性相关结果');
每一小题均给出实验过程与结果(含实验程序、运行的数据结果和图形);
实验三IIR数字滤波器的设计
一、实验目的
(1)掌握双线性变换法及脉冲响应不变法设计IIR数字滤波器的具体设计方法及其原理,熟悉用双线性变换法及脉冲响应不变法设计低通、高通和带通IIR数字滤波器的计算机编程。
(2)观察双线性变换及脉冲响应不变法设计的滤波器的频域特性,了解双线性变换法及脉冲响应不变法的特点。
(3)熟悉巴特沃思滤波器、切比雪夫滤波器和椭圆滤波器的频率特性。
二、实验内容
(1)P162例4.4
设采样周期T=250
(采样频率
=4kHz),分别用脉冲响应不变法和双线性变换法设计一个三阶巴特沃思低通滤波器,其3dB边界频率为
=1kHz。
fc=1000;
fs=4000;
[B,A]=butter(3,2*pi*fc,'s');
[num1,den1]=impinvar(B,A,fs);
[h1,w]=freqz(num1,den1);
f=w/pi*fs/2;
%双线性变换法
fc=1000;
fs=4000;
[B,A]=butter(3,2*fs*tan(pi*fc/fs),'s');
[num2,den2]=bilinear(B,A,fs);
[h2,w]=freqz(num2,den2);
f=w/pi*fs/2;
plot(f,abs(h1),'-.',f,abs(h2),'-');
gridon
legend('脉冲响应不变法','双线性变换法')
xlabel('频率/Hz')
ylabel('幅度')
(2)
=0.2kHz,
=1dB,
=0.3kHz,At=25dB,T=1ms;分别用脉冲响应不变法及双线性变换法设计一巴特沃思数字低通滤波器,观察所设计数字滤波器的幅频特性曲线,记录带宽和衰减量,检查是否满足要求。
比较这两种方法的优缺点。
fc=200;
Rp=1;
fr=300;
Rs=25;
T=0.001;
fs=1/T;
%脉冲响应不变法
Wp=2*pi*fc;
Ws=2*pi*fr;
[N,Wn]=buttord(Wp,Ws,Rp,Rs,'s')
[B,A]=butter(N,Wn,'s');
[num1,den1]=impinvar(B,A,fs);
[h1,w]=freqz(num1,den1);
%双线性变换法
Wp=2*fs*tan(pi*fc/fs);
Ws=2*fs*tan(pi*fr/fs);
[N,Wn]=buttord(Wp,Ws,Rp,Rs,'s')
[B,A]=butter(N,Wn,'s');
[num2,den2]=bilinear(B,A,fs);
[h2,w]=freqz(num2,den2);
f=w/pi*fs/2;
plot(f,abs(h1),'-.',f,abs(h2),'-');
gridon
legend('脉冲响应不变法','双线性变换法')
xlabel('频率/Hz')
ylabel('幅度')
N=
9
Wn=
1.3693e+003
N=
6
Wn=
1.7043e+003
(3)利用双线性变换法设计满足下列指标的巴特沃思数字低通滤波器,并作图验证设计结果:
=1.2kHz,
=2kHz,At
=8kHz。
fc=1200;
Rp=0.5;
fr=2000;
Rs=40;
fs=8000;
T=1/fs;
Wp=2*fs*tan(pi*fc/fs);
Ws=2*fs*tan(pi*fr/fs);
[N,Wn]=buttord(Wp,Ws,Rp,Rs,'s')
[B,A]=butter(N,Wn,'s');
[num2,den2]=bilinear(B,A,fs);
[h2,w]=freqz(num2,den2);
f=w/pi*fs/2;
plot(f,abs(h2),'-');gridon
legend('双线性变换法')
xlabel('频率/Hz')
ylabel('幅度')
每一小题均给出实验过程与结果(含实验程序、运行的数据结果和图形);
实验四FIR数字滤波器的设计
一、实验目的
(1)掌握用窗函数法,频率采样法及优化设计法设计FIR滤波器的原理及方法,熟悉相应的计算机编程;
(2)熟悉线性相位FIR滤波器的幅频特性和相频特性;
(3)了解各种不同窗函数对滤波器性能的影响。
二、实验内容
(1)N=45,计算并画出矩形窗、汉明窗、布莱克曼窗的归一化的幅度谱。
N=45;
n=0:
1:
N-1;
win=boxcar(N);%矩形窗
[H,w]=freqz(win,1);
H=H/H
(1);
stem(n,win)
figure
plot(w,abs(H));
N=45;
n=0:
1:
N-1;
win=hamming(N);%汉明窗
[H,w]=freqz(win,1);
H=H/H
(1);
stem(n,win)
figure
plot(w,abs(H));
N=45;
n=0:
1:
N-1;
win=blackman(N);%布莱克曼窗
[H,w]=freqz(win,1);
H=H/H
(1);
stem(n,win)
figure
plot(w,abs(H));
(2)用矩形窗设计一个21阶的线性相位低通FIR数字滤波器,截止频率Wc=0.25π,求出滤波器系数,并绘出滤波器的幅频特性。
修改程序,分别得到阶次为N=41,61的滤波器,并显示其各自的幅频曲线。
a)在上面所得的几幅图中,在截止频率两边可以观察到幅频响应的摆动行为。
请问波纹的数量与滤波器脉冲响应的长度之间有什么关系?
Wc=0.25*pi;
N=21;
M=(N-1)/2;%位移量
forn=0:
(N-1)
if(n==fix(M))%中间的点单独算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));
end;
end;
win=boxcar(N);
%%%不同窗函数
h=hd.*win';
[H,w]=freqz(h,1);
n=0:
1:
N-1;
subplot(3,1,1);
stem(n,h)
xlabel('n'),ylabel('h(n)')
subplot(3,1,2);
plot(w,abs(H));
xlabel('\omega'),ylabel('幅度谱')
subplot(3,1,3);
plot(w,angle(H));
xlabel('\omega'),ylabel('相位谱')
Wc=0.25*pi;
N=41;
M=(N-1)/2;%位移量
forn=0:
(N-1)
if(n==fix(M))%中间的点单独算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));
end;
end;
win=boxcar(N);
%%%不同窗函数
h=hd.*win';
[H,w]=freqz(h,1);
n=0:
1:
N-1;
subplot(3,1,1);
stem(n,h)
xlabel('n'),ylabel('h(n)')
subplot(3,1,2);
plot(w,abs(H));
xlabel('\omega'),ylabel('幅度谱')
subplot(3,1,3);
plot(w,angle(H));
xlabel('\omega'),ylabel('相位谱')
Wc=0.25*pi;
N=61;
M=(N-1)/2;%位移量
forn=0:
(N-1)
if(n==fix(M))%中间的点单独算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));
end;
end;
win=boxcar(N);
%%%不同窗函数
h=hd.*win';
[H,w]=freqz(h,1);
n=0:
1:
N-1;
subplot(3,1,1);
stem(n,h)
xlabel('n'),ylabel('h(n)')
subplot(3,1,2);
plot(w,abs(H));
xlabel('\omega'),ylabel('幅度谱')
subplot(3,1,3);
plo