数字信号处理实验报告.docx
《数字信号处理实验报告.docx》由会员分享，可在线阅读，更多相关《数字信号处理实验报告.docx（27页珍藏版）》请在冰豆网上搜索。
数字信号处理实验报告
中南大学
课程名称：
数字信号处理
指导老师：
陈科文
专业班级:
电子信息工程
姓名:
学号:
实验时间:
2014·10·8-------2014·10·17
目录
一、课程设计题目；
二、课程设计要求；
三、课程设计思想介绍；
四、课程设计流图；
五、课程设计结果及部分代码；
六、课程设计心得与体会；
七、参考文献；
八、附录（源代码）；
一、课程设计题目
1信号发生器
用户根据测试需要，可任选以下两种方式之一生成测试信号：
、直接输入（或从文件读取）测试序列；
②、输入由多个不同频率正弦信号叠加组合而成的模拟信号公式（如式1-1所示）、采样频率（Hz）、采样点数，动态生成该信号的采样序列，作为测试信号。
2频谱分析
使用FFT对产生的测试信号进行频谱分析并展示其幅频特性与相频特性，指定需要滤除的频带，通过选择滤波器类型（IIR/FIR），确定对应的滤波器（低通、高通）技术指标。
3滤波器设计
根据以上技术指标（通带截止频率、通带最大衰减、阻带截止频率、阻带最小衰减），设计数字滤波器，生成相应的滤波器系数，并画出对应的滤波器幅频特性与相频特性。
①IIRDF设计：
可选择滤波器基型（巴特沃斯或切比雪夫型）；
②FIRDF设计：
使用窗口法（可选择窗口类型，并比较分析基于不同窗口、不同阶数所设计数字滤波器的特点）。
4数字滤波
根据设计的滤波器系数，对测试信号进行数字滤波，展示滤波后信号的幅频特性与相频特性，分析是否满足滤波要求（对同一滤波要求，对比分析各类滤波器的差异）。
①IIRDF：
要求通过差分方程迭代实现滤波（未知初值置零处理）；
②FIRDF：
要求通过快速卷积实现滤波（对于长序列，可以选择使用重叠相加或重叠保留法进行卷积运算）。
5选做内容
将一段语音作为测试信号，通过频谱展示和语音播放，对比分析滤波前后语音信号的变化，进一步加深对数字信号处理的理解。
二、课程设计要求；
（1）、使用MATLAB（或其它开发工具）编程实现上述内容，写出课程设计报告。
（2）、课程设计报告的内容包括：
①课程设计题目和题目设计要求；
②设计思想和系统功能结构及功能说明；
③设计中关键部分的详细描述和介绍，采用流程图描述关键模块的设计思路；
④总结，包括设计过程中遇到的问题和解决方法，心得体会等；
⑤参考文献；
6程序源代码清单。
三、课程设计介绍；
用双线性变换法设计无限脉冲响应数字滤波器（IIRDF）时，先把数字滤波器指标转换成模拟滤波器的指标，然后根据模拟滤波器的指标设计模拟滤波器，再经过线性变换把模拟滤波器转换成数字滤波器。
该部分要能够设计巴特沃斯型低通滤波器，并能够输出该滤波器的幅频特性与滤波后结果。
该部分的关键是滤波器的设计部分，按照双线性变换法设计滤波器的步骤进行设计即可。
设计线性相位有限脉冲数字滤波器（FIRDF）时，先根据阻带衰减及过渡带的指标要求，选择窗函数类型，计算其窗口长度M；再构造希望逼近的频响函数Hdg（ejw），再计算hd（n）;最后进行加窗。
该部分要求能够设计出FIR滤波器，并输出窗函数、hd（n）函数的图形及损耗曲线与输出图形。
该部分的关键是设计窗函数及绘制损耗曲线，可通过查找资料进行完成。
四、课程设计流图；
IIR滤波器参数设置
IIR滤波器幅频输出
IIR滤波器设计
信号产生
IIR滤波
信号输入
信号输出
采样输出
幅频输出
采样输出
幅频输出
相频输出
相频输出
IIR滤波器流程图
h（n）函数输出
FIR滤波器参数设置
损耗曲线输出
FIR滤波器幅频输出
FIR滤波器设计
信号产生
FIR滤波
信号输入
信号输出
采样输出
幅频输出
采样输出
幅频输出
相频输出
相频输出
FIR滤波器流程图
五、课程设计结果及部分代码；
IIR滤波器运行结果及滤波器部分代码：
%滤波器设计%截止频率600Hz;
T=1;
Fr=1;
fs=20000;
wpz=0.05;wsz=0.06;
wp=2*tan（wpz*pi/2）;
ws=2*tan（wsz*pi/2）;
rp=1;rs=15;%频率低于500Hz内，衰减小于1db；频率大于600Hz时，衰减大于15db
[N,wc]=buttord（wp,ws,rp,rs,'s'）;
[B,A]=butter（N,wc,'s'）;
[Bz,Az]=bilinear（B,A,Fr）;
[h,w]=freqz（Bz,Az,256,fs*1000）;
[Nd,wdc]=buttord（wpz,wsz,rp,rs）;
[Bdz,Adz]=butter（Nd,wdc）;
figure;
plot（w,abs（h））;
xlabel（'频率'）;ylabel（'幅度'）;
title（'滤波器幅频相应'）;
FIR滤波器结果及滤波器部分代码：
（1）凯塞窗：
%滤波器设计
fp=1500;fs=2500;rs=40;%截止频率2500；
Fs=20000;
wp=2*pi*fp/Fs;
ws=2*pi*fs/Fs;
Bt=ws-wp;
alph=0.5842*（rs-21）^0.4+0.07886*（rs-21）;
M=ceil（（rs-8）/2.285/Bt）;
wc=（wp+ws）/2/pi;
hn=fir1（M,wc,kaiser（M+1,alph））;
figure;
%窗函数w（n）波形；
n=0:
M;
subplot（4,1,1）;
stem（n,kaiser（M+1,alph））;gridon;
xlabel（'n'）;ylabel（'w（n）'）;
title（'窗函数w（n）波形'）;
%h（n）波形；
subplot（4,1,2）;
stem（n,hn,'.'）;gridon;
xlabel（'n'）;ylabel（'h（n）'）;
title（'h（n）波形'）;
%损耗函数曲线；
[h1,w1]=freqz（hn,1）;
subplot（4,1,3）;
plot（w1/pi,20*log10（abs（h1）））;%通过求对数得到以db为单位的幅度特性
gridon;
axis（[0,1,-100,10]）;
xlabel（'归一化频率/pi'）;ylabel（'幅度/dB'）;
title（'损耗函数曲线'）;
（2）汉宁窗：
%滤波器设计
fp=1500;fs=2500;rs=40;%截止频率2500；
Fs=20000;
wp=2*pi*fp/Fs;
ws=2*pi*fs/Fs;
Bt=ws-wp;
alph=0.5842*（rs-21）^0.4+0.07886*（rs-21）;
M=ceil（（rs-8）/2.285/Bt）;
wc=（wp+ws）/2/pi;
hn=fir1（M,wc,kaiser（M+1,alph））;
figure;
%窗函数w（n）波形；
n=0:
M;
subplot（4,1,1）;
stem（n,kaiser（M+1,alph））;gridon;
xlabel（'n'）;ylabel（'w（n）'）;
title（'窗函数w（n）波形'）;
%h（n）波形；
subplot（4,1,2）;
stem（n,hn,'.'）;gridon;
xlabel（'n'）;ylabel（'h（n）'）;
title（'h（n）波形'）;
%损耗函数曲线；
[h1,w1]=freqz（hn,1）;
subplot（4,1,3）;
plot（w1/pi,20*log10（abs（h1）））;%通过求对数得到以db为单位的幅度特性
gridon;
axis（[0,1,-100,10]）;
xlabel（'归一化频率/pi'）;ylabel（'幅度/dB'）;
title（'损耗函数曲线'）;
六、课程设计心得与体会；
通过本次课程设计，终于让我在课堂上所有的收获得到了很好的表达，我在课程设计期间不断地去回忆课堂老师所讲的内容和自己在课后所了解的知识，在写代码的过程中也去翻看matlab书籍，这对我本次课设有了极大地帮助。
在实验中，我也遇到了很多种困难，像画幅频特性时就出现了横纵坐标不相对的问题，我改变采样点数及采样频率并检查频率公式经多次的调试才终于得到了正确的图像；在滤波器的设计中，也出现参数设置的问题，我自己设计的参数我也是经过多次计算才最终确定其截止频率与衰减。
在本次课设中，我不仅将理论与实际相结合，而且又在这个过程中学到了很多新的东西，而且对matlab编程的能力也有了很大的提高，真是个一举多得的事情。
我想这对我的实践动手能力及经验都是一个很好的挑战，对我以后的发展及工作都会是一个跳板，我很希望有机会接触这样有挑战性的实验，并努力提高自己各方面的综合素质。
七、参考文献：
《数字信号处理（第三版）》，高息全丁玉美编著西安电子科技大学出版社；
《数字信号处理教程》，陈怀琛编著，电子工业出版社；
《信号与线性系统分析（第四版）》，吴大正主编，高等教育出版社；
《数字信号处理实验指导书（matlab版）》，SanjitK.Mitra编著，电子工业出版社；
《MATLAB程序设计与应用（第二版）》，刘卫国编著，高等教育出版社；
《MATLAB7辅助信号处理技术与应用》，飞科技术产品研发中心编著，电子工业出版社；
八、附录（源代码）；
IIR巴特沃斯滤波器代码：
%滤波器设计%截止频率600Hz;
T=1;
Fr=1;
fs=20000;
wpz=0.05;wsz=0.06;
wp=2*tan（wpz*pi/2）;
ws=2*tan（wsz*pi/2）;
rp=1;rs=15;%频率低于500Hz内，衰减小于1db；频率大于600Hz时，衰减大于15db
[N,wc]=buttord（wp,ws,rp,rs,'s'）;
[B,A]=butter（N,wc,'s'）;
[Bz,Az]=bilinear（B,A,Fr）;
[h,w]=freqz（Bz,Az,256,fs*1000）;
[Nd,wdc]=buttord（wpz,wsz,rp,rs）;
[Bdz,Adz]=butter（Nd,wdc）;
figure;
plot（w,abs（h））;
xlabel（'频率'）;ylabel（'幅度'）;
title（'滤波器幅频相应'）;
%信号
x=0;
fs=20000;
N=1024;
n=0:
1023;
t=n/fs;
f0=10;
whilef0~=0
f0=input（'请输入频率值：
'）;
a=100.*sin（2*pi*f0*t）;
x=x+a;
end
figure;
N=1024;
fs=20000;
n=0:
N-1;
t=n/fs;
%画出信号波形
subplot（4,2,1）;plot（t,x）;gridon;
xlabel（'t'）;ylabel（'f（t）'）;
title（'f（t）波形'）;
%画出信号采样后波形
subplot（4,2,3）;stem（t,x,'.'）;gridon;
xlabel（'n'）;ylabel（'f（n）'）;
title（'f（n）波形'）;
%画出幅频特性
%T=1;
%w=linspace（0,2*pi,N）./pi;
%dt=t
（2）-t
（1）;
%f=1/dt;
%ft=fft（x,N）;
%F=ft（1:
N/2+1）;
%f=f*（0:
N/2）/N;
%subplot（4,2,5）;plot（f,abs（F）,'--'）;gridon;
%axis（[01200010000]）;
ft=fft（x,N）;
w=linspace（0,2*pi,1024）./pi;
f=[0:
N-1]*fs/N;
hag=abs（ft）;
subplot（4,2,5）;
plot（f（1:
ceil（N/2））,hag（1:
ceil（N/2）））;gridon;
axis（[01200010000]）;
xlabel（'f'）;ylabel（'|F（K）|'）;
title（'幅频特性'）;
%画出相频特性
ang=angle（ft）;
subplot（4,2,7）;
plot（w（1:
1024）,ang（1:
1024））;gridon;
axis（[02-44]）;
xlabel（'w'）;ylabel（'Φ（w）'）;
title（'相频特性'）;
%滤波
y=filter（Bdz,Adz,x）;
%画出滤波后信号
subplot（4,2,2）;
plot（t,y）;gridon;
xlabel（'t'）,ylabel（'f''（t）'）;
title（'滤波后f（t）波形'）;
%画出滤波采样后波形
subplot（4,2,4）;stem（t,y,'.'）;gridon;
xlabel（'n'）;ylabel（'f''（n）'）;
title（'滤波后f（n）波形'）;
%画出滤波后幅频特性
fw=fft（y,N）;
ha=abs（fw）;
subplot（4,2,6）;
plot（f（1:
ceil（N/2））,ha（1:
ceil（N/2）））;gridon;
axis（[01200010000]）;
xlabel（'f'）,ylabel（'|F（K）|'）;
title（'滤波后抽样信号幅频特性'）;
%画出滤波后相频特性
an=angle（fw）;
subplot（4,2,8）;
plot（w（1:
1024）,an（1:
1024））;gridon;
axis（[02-44]）;
xlabel（'w'）,ylabel（'Φ''（w）'）;
title（'滤波后抽样信号相频特性'）;
FIR凯塞窗滤波器代码：
%信号
x=0;
f0=10;
whilef0~=0
f0=input（'请输入频率值：
'）;
a=100.*sin（2*pi*f0*t）;
x=x+a;
end
%滤波器设计
fp=1500;fs=2500;rs=40;%截止频率2500；
Fs=20000;
wp=2*pi*fp/Fs;
ws=2*pi*fs/Fs;
Bt=ws-wp;
alph=0.5842*（rs-21）^0.4+0.07886*（rs-21）;
M=ceil（（rs-8）/2.285/Bt）;
wc=（wp+ws）/2/pi;
hn=fir1（M,wc,kaiser（M+1,alph））;
figure;
%窗函数w（n）波形；
n=0:
M;
subplot（4,1,1）;
stem（n,kaiser（M+1,alph））;gridon;
xlabel（'n'）;ylabel（'w（n）'）;
title（'窗函数w（n）波形'）;
%h（n）波形；
subplot（4,1,2）;
stem（n,hn,'.'）;gridon;
xlabel（'n'）;ylabel（'h（n）'）;
title（'h（n）波形'）;
%损耗函数曲线；
[h1,w1]=freqz（hn,1）;
subplot（4,1,3）;
plot（w1/pi,20*log10（abs（h1）））;%通过求对数得到以db为单位的幅度特性
gridon;
axis（[0,1,-100,10]）;
xlabel（'归一化频率/pi'）;ylabel（'幅度/dB'）;
title（'损耗函数曲线'）;
%对信号进行滤波
Y=filter（hn,2,x）;
%滤波后信号波形；
subplot（4,1,4）;
plot（t,Y）;gridon;
axis（[00.01-400400]）;
xlabel（'t'）;ylabel（'f''（t）'）;
title（'滤波后的信号'）;
figure;
%定义变量
fs=20000;
N=1024;
n=0:
1023;
t=n/fs;
%画出信号波形
subplot（4,2,1）;plot（t,x）;gridon;
axis（[00.01-400400]）;
xlabel（'t'）;ylabel（'f（t）'）;
title（'f（t）波形'）;
%画出信号采样后波形
subplot（4,2,3）;stem（t,x,'.'）;gridon;
axis（[00.01-400400]）;
xlabel（'n'）;ylabel（'f（n）'）;
title（'f（n）波形'）;
%画出幅频特性
%T=1;
%w=linspace（0,2*pi,N）./pi;
%dt=t
（2）-t
（1）;
%f=1/dt;
%ft=fft（x,N）;
%F=ft（1:
N/2+1）;
%f=f*（0:
N/2）/N;
%subplot（2,5,3）;plot（f,abs（F）,'--'）;gridon;
%axis（[01200010000]）;
ft=fft（x,N）;
w=linspace（0,2*pi,1024）./pi;
f=[0:
N-1]*fs/N;
hag=abs（ft）;
subplot（4,2,5）;
plot（f（1:
ceil（N/2））,hag（1:
ceil（N/2）））;gridon;
axis（[010000010000]）;
xlabel（'f'）;ylabel（'|X（e（jw））|'）;
title（'幅频特性'）;
%画出相频特性
ang=angle（ft）;
subplot（4,2,7）;
stem（w（1:
1024）,ang（1:
1024）,'.'）;gridon;
axis（[01-44]）;
xlabel（'w'）;ylabel（'Φ（w）'）;
title（'相频特性'）;
%画出滤波后信号波形
subplot（4,2,2）;plot（t,Y）;gridon;
axis（[00.01-400400]）;
xlabel（'t'）;ylabel（'f''（t）'）;
title（'滤波后f（t）波形'）;
%画出滤波后信号采样波形
subplot（4,2,4）;stem（t,Y,'.'）;gridon;
axis（[00.01-400400]）;
xlabel（'n'）;ylabel（'f''（n）'）;
title（'滤波后f（n）波形'）;
%画出滤波后信号幅频特性
fw=fft（Y,N）;
w=linspace（0,2*pi,1024）./pi;
f=[0:
N-1]*fs/N;
ha=abs（fw）;
subplot（4,2,6）;
plot（f（1:
ceil（N/2））,ha（1:
ceil（N/2）））;gridon;
axis（[010000010000]）;
xlabel（'f'）;ylabel（'|X''（e（jw））|'）;
title（'滤波后幅频特性'）;
%画出滤波后相频特性
an=angle（fw）;
subplot（4,2,8）;
stem（w（1:
1024）,an（1:
1024）,'.'）;gridon;
axis（[01-44]）;
xlabel（'w'）;ylabel（'Φ''（w）'）;
title（'滤波后相频特性'）;
FIR汉宁窗滤波器代码：
%信号
x=0;
f0=10;
whilef0~=0
f0=input（'请输入频率值：
'）;
a=100.*sin（2*pi*f0*t）;
x=x+a;
end
%滤波器设计
fp=1500;fs=2500;rs=40;%截止频率2500；
Fs=20000;
wp=2*pi*fp/Fs;
ws=2*pi*fs/Fs;
Bt=ws-wp;
alph=0.5842*（rs-21）^0.4+0.07886*（rs-21）;
M=ceil（（rs-8）/2.285/Bt）;
wc=（wp+ws）/2/pi;
hn=fir1（M,wc,kaiser（M+1,alph））;
figure;
%窗函数w（n）波形；
n=0:
M;
subplot（4,1,1）;
stem（n,kaiser（M+1,alph））;gridon;
xlabel（'n'）;ylabel（'w（n）'）;
title（'窗函数w（n）波形'）;
%h（n）波形；
subplot（4,1,2）;
stem（n,hn,'.'）;gridon;
xlabel（'n'）;ylabel（'h（n）'）;
title（'h（n）波形'）;
%损耗函数曲线；
[h1,w1]=freqz（hn,1）;
subplot（4,1,3）;
plot（w1/pi,20*log10（abs（h1）））;%通过求对数得到以db为单位的幅度特性
gridon;
axis（[0,1,-100,10]）;
xlabel（'归一化频率/pi'）;ylabel（'幅度/dB'）;
title（'损耗函数曲线'）;
%对信号进行滤波
Y=filter（hn,2,x）;
%滤波后信号波形；
subplot（4,1,4）;
plot（t,Y）;gridon;
axis（[00.01-400400]）;
xlabel（'t'）;ylabel（'f''（t）'）;
title（'滤波后的信号'）;
figure;
%定义变量
fs=20000;
N=1024;
n=0:
1023;
t=n/fs;
%画出信号波形
%subplot（4,2,1）;plot（t,x）;gridon;
subplot（4,2,1）;plot（t,x）;gridon;
axis（[00.01-400400]）;
xlabel（'t'）;ylabel（'f（t）'）;
title（'f（t）波形'）;
%画出信号采样后波形
subplot（4,2,3）;stem（t,x,'.'）;gridon;
axis（[00.01-400400]）;
xlabel（'n'）;ylabel（'f（n）'）;
title（'f（n）波形'）;
%画出幅频特性
%T=1;
%w=linspace（0,2*pi,N）./pi;
%dt=t
（2）-t
（1）;
%f=1/dt;
%ft=fft（x,N）;
%F=ft（1:
N/2+1）;
%f=f*（0:
N/2）/N;
%subplot（2,5,3）;plot（f,abs（F）,'--'）;gridon;
%axis（[01200010000]）;
升级会员 