数字信号处理IIR滤波器设计.docx
《数字信号处理IIR滤波器设计.docx》由会员分享，可在线阅读，更多相关《数字信号处理IIR滤波器设计.docx（45页珍藏版）》请在冰豆网上搜索。
数字信号处理IIR滤波器设计
课程设计报告
课程名称数字信号处理
课题名称IIR滤波器设计
专业通信工程
班级通信工程
学号
姓名
指导教师胡瑛张细政
2013年9月8日
湖南工程学院
课程设计任务书
课程名称数字信号处理
课题IIR滤波器设计
专业班级通信工程
学生姓名
学号
指导老师胡瑛张细政
审批
任务书下达日期2013年9月1日
任务完成日期2013年9月8日
《数字信号处理》课程设计任务书
一、设计目的
综合运用数字信号处理的理论知识进行频谱分析和滤波器设计，通过理论推导得出相应结论，再利用MATLAB作为编程工具进行计算机实现，从而加深对所学知识的理解，建立概念。
二、设计要求
1、MATLAB的使用，掌握MATLAB的程序设计方法。
2、Windows环境下语音信号采集的方法。
3、数字信号处理的基本概念、基本理论和基本方法。
4、TLAB设计FIR和nR数字滤波器的方法。
5、用MATLAB对信号进行分析和处理
6、计报告4000以上，含程序设计说明，用户使用说明，源程序清单及程序框图。
7、机演示。
8、有详细的文档。
文档中包括设计思路、设计仿真程序、仿真结果及相应的分析与结论。
三、进度安排
第一周星期一：
课题讲解，查阅资料
星期二:
总体设计，详细设计
星期三：
编程，上机调试、修改程序
星期四：
上机调试、完善程序
星期五：
答辩
星期六-星期天：
撰写课程设计报告
附：
课程设计报告装订顺序：
封面、任务书、目录、正文、评分、附件（A4大小的图纸及程序清单）。
正文的格式:
一级标题用3号黑体,二级标题用四号宋体加粗,正文用小四号宋体;行距为22。
正文的内容:
一、课题的主要功能；二、课题的功能模块的划分；三、主要功能的实现；四、程序调试；五、总结；六、附件（所有程序的原代码，要求对程序写出必要的注释）；七、评分表。
1、课题的主要功能
1.1问题描述
采用IIR滤波器的完全设计函数来设计各型IIR滤波器（Butterworth型与切比雪夫I型、切比雪夫Ⅱ型和椭圆型的低通、高通、带通、带阻中的至少2种类型）,绘出滤波器的频域响应；并用这些数字滤波器对含噪语音信号分别进行滤波处理，比较不同方法下设计出来的数字滤波器的滤波效果，并从理论上进行分析（或解释）。
1.2功能要求
1、语音信号采集（每个同学必须录制本人的语音信号）
录制一段课程设计学生的语音信号并保存为文件，要求长度不小于10秒，并对录制的信号进行采样；录制时可以使用Windows自带的录音机，或者使用其它专业的录音软件，录制时需要配备录音硬件（如麦克风），为便于比较，需要在安静、干扰小的环境下录音。
2、语音信号分析
使用MATLAB绘出采样后的语音信号的时域波形和频谱图。
根据频谱图求出其带宽，并说明语音信号的采样频率不能低于多少赫兹。
3、含噪语音信号合成
在MATLAB软件平台下，给原始的语音信号叠加上噪声，噪声类型分为如下几种：
（1）白噪声；
（2）单频噪色（正弦干扰）；（3）多频噪声（多正弦干扰）；（4）其它干扰，可设置为低频、高频、带限噪声，或Chirp干扰、冲激干扰。
绘出叠加噪声后的语音信号时域和频谱图，在视觉上与原始语音信号图形对比，也可通过Windows播放软件从听觉上进行对比，分析并体会含噪语音信号频谱和时域波形的改变。
4、数字滤波器设计及滤波
利用不同的原型低通滤波器（Butterworth型与切比雪夫I型、切比雪夫Ⅱ型和椭圆型）来设计各型IIR滤波器（低通、高通、带通、带阻中的至少3种类型），绘出滤波器的频域响应；并用这些数字滤波器对含噪语音信号分别进行滤波处理，比较不同方法下设计出来的数字滤波器的滤波效果，并从理论上进行分析（或解释）。
5、回放语音信号
对滤波后的语音信号进行回放，感觉滤波前后语音信号的变化。
2.、课题的功能模块的划分
2.1语音信号采集（Butterworth型低通滤波器为例）
[x1,fs,bits]=wavread（'d:
/baojing2'）;
sound（x1,fs,bits）;
y1=fft（x1,32768）;
2.2语音信号分析
原始信号
[x1,fs,bits]=wavread（'d:
/baojing2'）;
sound（x1,fs,bits）;
y1=fft（x1,32768）;
figure
（1）
subplot（2,1,1）
plot（x1）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y1））;
title（'原始语音信号频谱'）
2.3含噪语音信号合成
加噪后的语音信号
N=length（x1）-1;
t=0:
1/fs:
N/fs;
d=[0.8*cos（2*pi*5000*t）]';
x3=x1+d;
%sound（x3,fs,bits）;
figure
（2）
subplot（2,1,1）
plot（x3）
title（'加单频余弦信号的语音信号时域波形'）
y3=fft（x3,32768）;
subplot（2,1,2）
plot（abs（y3））;
title（'加单频余弦信号的语音信号频谱'）
2.4数字滤波器设计及滤波
加噪后的语音信号进行滤波处理
Fs=fs;%采样频率
wp=4000*2/Fs;ws=4300*2/Fs;%根据采样频率将滤波器边界频率进行转换
Rp=1;Rs=20;%通带波纹和阻带衰减
Nn=128;%显示滤波器频率特性的数据长度
[N,Wn]=buttord（wp,ws,Rp,Rs）;%求得数字滤波器的最小阶数和截止频率（归一化频率）
[b,a]=butter（N,Wn）;%设计Butterworth高通数字滤波器
figure（3）
[H,f]=freqz（b,a,Nn,Fs）;%用Nn点绘出频率特性
subplot（2,1,1）,plot（f,20*log10（abs（H）））;
xlabel（'频率/Hz'）;ylabel（'振幅/dB'）;gridon;
subplot（2,1,2）,plot（f,180/pi*unwrap（angle（H）））
xlabel（'频率/Hz'）;ylabel（'相位/^o'）;gridon;
figure（4）
y=filter（b,a,x3）;%对输入信号进行滤波
2.5回放语音信号
subplot（2,1,1）
plot（y）;
title（'滤波语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y4））;
title（'滤波语音信号频谱'）%
sound（y,fs,bits）;
3、滤波器实现的原理
IIR数字滤波器设计原理：
IIR数字滤波器是通过因果稳定的Ha（s）映射成因果稳定的H（z），即s平面的左半平面必须映射到z平面单位圆的内部。
H（z）的频率响应能模仿Ha（s）的频率响应，即s平面的虚轴必须映射到z平面的单位圆上。
变换前后的滤波器在时域或频域的主要特征（频率响应或单位冲激响应等）应尽可能相同或接近。
将传输函数Ha（s）从s平面转换到z平面的方法有多种，主要有冲激不变法和双线性变换法。
如果不考虑混叠现象，利用冲击不变法实现的数字滤波器会很好地重现原模拟滤波器的频率特性。
而且数字滤波器的单位脉冲响应完全模仿模拟滤波器的单位冲激响应，时域特性逼近好。
在本设计中，用冲激函数作为系统激励信号，用各种数字滤波器作为测试系统。
冲激函数具有无限宽广的频谱，用冲激函数做激励信号相当于对测试系统输入所有频率的信号，系统必然有对应的输出。
用Transfer函数计算出系统输出与输入的傅立叶变换之比，从而得到系统的频率响应函数。
IIR数字滤波器的系统函数为
的有理分式：
设计IIR滤波器的系统函数，就是要确定
的阶数N及分子分母多项式的系数
和
，使其
满足指定的频率特性。
由于模拟滤波器的设计有许多简单而严谨的设计公式和大量的图表可以利用，因此IIR滤波器设计的方法之一是：
先设计一个合适的模拟滤波器，然后将模拟滤波器通过适当的变换转换成满足给定指标的数字滤波器。
1、Butterworth模拟低通滤波器
幅度平方函数：
其中，N为滤波器的阶数，
为通带截止频率。
2．Chebyshev模拟低通滤波器
IIR数字滤波器经典设计法的一般步骤为：
根据给定的性能指标和方法不同，首先对设计性能指标中的频率指标，如数字边界频率进行变换，转换后的模拟频率指标作为模拟滤波器原型设计的性能指标。
估计模拟滤波器最小阶数和截止频率，利用MATLAB工具函数buttord、cheb1ord、cheb2ord、ellipord等。
设计模拟低通滤波器原型。
利用MATLAB工具函数buttap、cheb1ap、cheb2ap、ellipap等。
由模拟原型低通滤波器经频率变换获得模拟滤波器（低通、高通、带通、带阻等），利用MATLAB工具函数lp2lp、lp2hp、lp2bp、lp2bs。
将模拟滤波器离散化获得IIR数字滤波器，利用MATLAB工具函数bilinear或impinvar。
后面的几个步骤在前面的章节中已经论述，这里主要介绍第一步，关于设计性能指标的转换。
设计IIR滤波器时，给出的性能指标通常分数字指标和模拟指标两种。
数字性能指标给出通带截止频率
，阻带起始频率
，通带波纹Rp，阻带衰减Rs等。
数字频率
和
的取值范围为0~
，单位弧度。
而MATLAB工具函数常采用归一化频率，
和
的取值范围为0~1，对应于0～
，此时需进行转换。
模拟性能指标给出通带截止频率
，阻带起始频率
，通带波纹Rp，阻带衰减Rs等。
模拟频率
和
单位为弧度/秒（rad/s）。
MATLAB信号处理工具箱中，设计性能指标的转换应根据不同设计方法进行不同处理。
4、程序调试
4.1Butterworth型低通滤波器
图4-1-1图4-1-2
图4-1-3图4-1-4
4.2Butterworth型高通滤波器
图4-2-1图4-2-2
图4-2-3图4-2-4
4.3切比雪夫I型低通滤波器
图4-3-1图4-3-2
图4-3-3图4-3-4
4.4切比雪夫I型带通滤波器
图4-4-1图4-4-2
图4-4-3图4-4-4
4.5切比雪夫Ⅱ型带通滤波器
图4-5-1图4-5-2
图4-5-3图4-5-4
4.6切比雪夫Ⅱ型低通滤波器
图4-6-1图4-6-2
图4-6-3图4-6-4
4.7椭圆型的低通
图4-7-1图4-7-2
图4-7-3图4-7-4
4.8椭圆型的高通
图4-8-1图4-8-2
图4-8-3图4-8-4
5、总结
通过这次课程设计，我对matlab软件更加了解，更加熟练，用得更加顺手，对数字信号中的滤波器也有了新的认识，本次我是采用IIR滤波器的完全设计函数来设计各型IIR滤波器（Butterworth型与切比雪夫I型、切比雪夫Ⅱ型和椭圆型的低通、高通、带通、带阻中的至少2种类型）,并要求绘出滤波器的频域响应；并用这些数字滤波器对含噪语音信号分别进行滤波处理，它的设计方法：
（1）根据给定的性能指标和方法不同，首先对设计性能指标中的频率指标，如数字边界频率进行变换，转换后的模拟频率指标作为模拟滤波器原型设计的性能指标。
（2）估计模拟滤波器最小阶数和截止频率，利用MATLAB工具函数buttord、cheb1ord、cheb2ord、ellipord等。
（3）设计模拟低通滤波器原型。
利用MATLAB工具函数buttap、cheb1ap、cheb2ap、ellipap等。
（4）由模拟原型低通滤波器经频率变换获得模拟滤波器（低通、高通、带通、带阻等），利用MATLAB工具函数lp2lp、lp2hp、lp2bp、lp2bs。
（5）将模拟滤波器离散化获得IIR数字滤波器，利用MATLAB工具函数bilinear或impinvar。
数字滤波器可以通过编程实现各种不同系统，满足不同的需要，也可以随时改动系数，调整滤波器参数，选择最佳方案。
使用LabVIEW软件平台开发电气参数测量仪等虚拟仪器，实现了更高的效率，节省了更多的硬件开销，方便了系统的维护和减轻了仪器更新的负担。
使用虚拟仪器逐步代替传统仪器已经成为测试领域发展的趋势。
但是在实际应用中，仍要根据具体情况进行程序的优化和软硬件的结合，使虚拟仪器发挥更高的性能。
6、附件
Butterworth型低通滤波器
[x1,fs,bits]=wavread（'d:
/baojing2'）;
%sound（x1,fs,bits）;
y1=fft（x1,32768）;
figure
（1）
subplot（2,1,1）
plot（x1）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y1））;
title（'原始语音信号频谱'）
N=length（x1）-1;
t=0:
1/fs:
N/fs;
d=[0.8*cos（2*pi*5000*t）]';
x3=x1+d;
%sound（x3,fs,bits）;
figure
（2）
subplot（2,1,1）
plot（x3）
title（'加单频余弦信号的语音信号时域波形'）
y3=fft（x3,32768）;
subplot（2,1,2）
plot（abs（y3））;
title（'加单频余弦信号的语音信号频谱'）
Fs=fs;%采样频率
wp=4000*2/Fs;ws=4300*2/Fs;%根据采样频率将滤波器边界频率进行转换
Rp=1;Rs=20;%通带波纹和阻带衰减
Nn=128;%显示滤波器频率特性的数据长度
[N,Wn]=buttord（wp,ws,Rp,Rs）;%求得数字滤波器的最小阶数和截止频率（归一化频率）
[b,a]=butter（N,Wn）;%设计Butterworth高通数字滤波器
figure（3）
[H,f]=freqz（b,a,Nn,Fs）;%用Nn点绘出频率特性
subplot（2,1,1）,plot（f,20*log10（abs（H）））;
xlabel（'频率/Hz'）;ylabel（'振幅/dB'）;gridon;
subplot（2,1,2）,plot（f,180/pi*unwrap（angle（H）））
xlabel（'频率/Hz'）;ylabel（'相位/^o'）;gridon;
figure（4）
y=filter（b,a,x3）;%对输入信号进行滤波
y4=fft（y,32768）;
subplot（2,1,1）
plot（y）;
title（'滤波语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y4））;
title（'滤波语音信号频谱'）%
sound（y,fs,bits）;
Butterworth型高通滤波器
[x1,fs,bits]=wavread（'d:
/baojing2'）;
sound（x1,fs,bits）;
y1=fft（x1,32768）;
figure
（1）
subplot（2,1,1）
plot（x1）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y1））;
title（'原始语音信号频谱'）
N=length（x1）-1;
t=0:
1/fs:
N/fs;
d=[0.8*cos（2*pi*40*t）]';
x3=x1+d;
%sound（x3,fs,bits）;
figure
（2）
subplot（2,1,1）
plot（x3）
title（'加单频余弦信号的语音信号时域波形'）
y3=fft（x3,32768）;
subplot（2,1,2）
plot（abs（y3））;
title（'加单频余弦信号的语音信号频谱'）
Fs=1000;%采样频率
wp=50*2/Fs;ws=40*2/Fs;%根据采样频率将滤波器边界频率进行转换
Rp=1;Rs=20;%通带波纹和阻带衰减
Nn=128;%显示滤波器频率特性的数据长度
[N,Wn]=buttord（wp,ws,Rp,Rs）;%求得数字滤波器的最小阶数和截止频率（归一化频率）
[b,a]=butter（N,Wn,'high'）;%设计Butterworth高通数字滤波器
figure（3）
[H,f]=freqz（b,a,Nn,Fs）;%用Nn点绘出频率特性
subplot（2,1,1）,plot（f,20*log10（abs（H）））;
xlabel（'频率/Hz'）;ylabel（'振幅/dB'）;gridon;
subplot（2,1,2）,plot（f,180/pi*unwrap（angle（H）））
xlabel（'频率/Hz'）;ylabel（'相位/^o'）;gridon;
figure（4）
y=filter（b,a,x3）;%对输入信号进行滤波
y4=fft（y,32768）;
subplot（2,1,1）
plot（y）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y4））;
title（'原始语音信号频谱'）%
sound（y,fs,bits）;
切比雪夫I型带通滤波器
[x1,fs,bits]=wavread（'d:
/baojing2'）;
sound（x1,fs,bits）;
y1=fft（x1,32768）;
figure
（1）
subplot（2,1,1）
plot（x1）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y1））;
title（'原始语音信号频谱'）
N=length（x1）-1;
t=0:
1/fs:
N/fs;
d=[0.8*cos（2*pi*5000*t）]';
x3=x1+d;
%sound（x3,fs,bits）;
figure
（2）
subplot（2,1,1）
plot（x3）
title（'加单频余弦信号的语音信号时域波形'）
y3=fft（x3,32768）;
subplot（2,1,2）
plot（abs（y3））;
title（'加单频余弦信号的语音信号频谱'）
Fs=fs;%采样频率
wp=[504000]*2/Fs;%通带边界频率（归一化频率）（6-20式）
ws=[305000]*2/Fs;%阻带边界频率（归一化频率）（6-20式）
Rp=1;Rs=30;Nn=128;%通带波纹和阻带衰减以及绘制频率特性的数据点数
[N,Wn]=cheb1ord（wp,ws,Rp,Rs）;%求得数字滤波器的最小阶数和归一化截止频率
[b,a]=cheby1（N,Rp,Wn）;%按最小阶数、通带波纹和截止频率设计数字滤波器
figure（3）
[H,f]=freqz（b,a,Nn,Fs）;%求得滤波器的频率特性
subplot（2,1,1）,plot（f,20*log10（abs（H）））;
xlabel（'频率/Hz'）;ylabel（'振幅/dB'）;gridon;
subplot（2,1,2）,plot（f,180/pi*unwrap（angle（H）））
xlabel（'频率/Hz'）;ylabel（'相位/^o'）;gridon;
figure（4）
y=filtfilt（b,a,x3）;%对输入信号进行滤波
y4=fft（y,32768）;
subplot（2,1,1）
plot（y）;
title（'滤波语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y4））;
title（'滤波语音信号频谱'）%
sound（y,fs,bits）;
切比雪夫I型低通滤波器
[x1,fs,bits]=wavread（'d:
/baojing2'）;
sound（x1,fs,bits）;
y1=fft（x1,32768）;
figure
（1）
subplot（2,1,1）
plot（x1）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y1））;
title（'原始语音信号频谱'）
N=length（x1）-1;
t=0:
1/fs:
N/fs;
d=[0.8*cos（2*pi*5000*t）]';
x3=x1+d;
%sound（x3,fs,bits）;
figure
（2）
subplot（2,1,1）
plot（x3）
title（'加单频余弦信号的语音信号时域波形'）
y3=fft（x3,32768）;
subplot（2,1,2）
plot（abs（y3））;
title（'加单频余弦信号的语音信号频谱'）
Fs=fs;%采样频率
wp=4000*2/Fs;%通带边界频率（归一化频率）（6-20式）
ws=4500*2/Fs;%阻带边界频率（归一化频率）（6-20式）
Rp=1;Rs=30;Nn=128;%通带波纹和阻带衰减以及绘制频率特性的数据点数
[N,Wn]=cheb1ord（wp,ws,Rp,Rs）;%求得数字滤波器的最小阶数和归一化截止频率
[b,a]=cheby1（N,Rp,Wn）;%按最小阶数、通带波纹和截止频率设计数字滤波器
figure（3）
[H,f]=freqz（b,a,Nn,Fs）;%求得滤波器的频率特性
subplot（2,1,1）,plot（f,20*log10（abs（H）））;
xlabel（'频率/Hz'）;ylabel（'振幅/dB'）;gridon;
subplot（2,1,2）,plot（f,180/pi*unwrap（angle（H）））
xlabel（'频率/Hz'）;ylabel（'相位/^o'）;gridon;
figure（4）
y=filtfilt（b,a,x3）;%对输入信号进行滤波
y4=fft（y,32768）;
subplot（2,1,1）
plot（y）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y4））;
title（'原始语音信号频谱'）%
sound（y,fs,bits）;
切比雪夫Ⅱ型低通滤波器
[x1,fs,bits]=wavread（'d:
/baojing2'）;
sound（x1,fs,bits）;
y1=fft（x1,32768）;
figure
（1）
subplot（2,1,1）
plot（x1）;
title（'原始语音信号时域波形'）;
subplot（2,1,2）
plot（abs（y1））;
ti
升级会员 