数字信号处理课程
设计
领导形象设计圆作业设计ao工艺污水处理厂设计附属工程施工组织设计清扫机器人结构设计
一、课程设计目的
1、掌握MATLAB用于数字滤波器的设计技术和
方法
快递客服问题件处理详细方法山木方法pdf计算方法pdf华与华方法下载八字理论方法下载
；
2、掌握FIR滤波器的设计方法和步骤；
3、掌握IIR滤波器的设计方法和步骤。
二、课程设计任务及步骤
题目1、双线性变换法设计数字滤波器
（1） 采用双线性变换法设计一个巴特沃斯数字低通滤波器，要求：wp=0.25
，Rp=1dB，ws=0.4
，As=15dB，滤波器的采样频率为Fs=100Hz。
设计步骤：
1、用MATLAB信号处理工具箱函数buttord，butter是巴特沃斯滤波器设计函数，该程序中调用以下格式：
[N,wc]=buttord(wp,ws,Rp,As)
该格式用于计算巴特沃斯数字滤波器的阶数N和3dB截止频率wc，调用参数wp和ws分别为巴特沃斯滤波器的通带边界频率和阻带边界频率的归一化值，0<
表
关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf
示为设计高通滤波器。
开始
2、滤波器设计流程图
读入要求设计指标的fp,fs,As,Rp,Fs
用双极性变换法转换wp,ws为wp1,ws1
调用工具箱函数ellipord, ellip求N,wpo,B,A
调用bilinear将模拟转为数字滤波器参数sh
调用子程序(函数)计算H，w，pha
调用子程序(函数)绘出H，pha，损耗曲线
结束
3、程序实现代码如下：
%双极性变换法设计IIR椭圆数字高通滤波器
close all;clear all;
%============================================
Rp=1; As=20; fp=250; fs=150; Fs=1000;T=1/Fs;
%============================================
wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; %设置指标参数
wp1=(2/T)*tan(wp/2); ws1=(2/T)*tan(ws/2);%双极性变换法边界频率转换指标
H1=10^(-Rp/20); %通带波纹幅度 Rp=-20lg(H1)
H2=10^(-As/20); %阻带波纹幅度 As=-20lg(H2)
[N,wpo]=ellipord(wp1,ws1,Rp,As,'s');%计算椭圆高通滤波器最低阶数和通带边界频率
[B,A] =ellip(N,Rp,As,wpo,'high','s');%计算滤波器系统函数的系数
[Bd,Ad]=bilinear(B,A,Fs);%用双极性变换法转换成数字滤波器
[H,w] =freqz(Bd,Ad);%计算H(z)的幅频响应
subplot(2,2,1);plot(w/(2*pi)*Fs,abs(H)); xlabel('频率(Hz)'); ylabel('幅度(dB)')
set(gca,'XTick',[0,fs,fp,Fs/2]); set(gca,'YTick',[0,H2,H1,1]);
text(fp,H1,'<--通带最大衰减','FontSize',8);%文本注释
text(fs,H2,'<--阻带最小衰减','FontSize',8);%文本注释
title('椭圆高通滤波器');axis([0,Fs/2,0,1.1]);grid on;
subplot(2,2,2);pha=angle(H)*180/pi;plot(w/2/pi*Fs,pha);grid on
set(gca,'XTick',[0,fs,fp,Fs/2]); set(gca,'YTick',[-180,0,180]);title('相位响应曲线')
subplot(2,2,3);myplot(Bd,Ad,As);grid onset(gca,'XTick',[0,0.3,0.5,1]);
4、运行结果及相应波形
N =
3
wpo =
2.0000e+003
（3）采用双极性变换法设计一个切比雪夫I型数字高通滤波器，要求：通带wp1=0.3
，wp2=0.7
，Rp=1dB，阻带ws1=0.2
，ws2=0.8
，As=20dB，滤波器采样频率Fs=1000Hz。
设计步骤：
1、用MATLAB信号处理工具箱函数cheb1ord，cheby1是切比雪夫I型滤波器设计函数，该程序中调用以下格式：
[N,wc]=cheb1ord(wp,ws,Rp,As,’s’)
该格式用于计算切比雪夫I型数字滤波器的最低阶数N和通带边界频率wc，
[B,A]=cheby1(N,Rp,wc,’s’)
该格式用于计算切比雪夫I型数字滤波器系统函数分子和分母多项式的字数向量B,A。调用参数N,和wc分别为切比雪夫I型数字滤波器的阶数和通带边界频率wc。
开始
2、滤波器设计流程图
读入要求设计指标的wp1,ws1,wp2,ws2,As,Rp,Fs
用双极性变换法转换wp1,ws1,wp2,ws2为Wp1,Ws1,Wp2,Ws2
调用工具箱函数cheb1ord, cheby1求N,wc,B,A
调用bilinear将模拟转换为数字滤波器参数
调用子程序(函数)计算H，w，pha
调用子程序(函数)绘出H，pha，损耗曲线
结束
3、程序实现代码如下：
%双极性变换法设计IIR切比雪夫I型数字带通滤波器
close all;clear all;
%=============要求设计指标===================
wp1=0.3*pi; wp2=0.7*pi; ws1=0.2*pi; ws2=0.8*pi; Rp=1;
As=20; Fs=1000;T=1/Fs;
%============================================
Wp1=(2/T)*tan(wp1/2); Wp2=(2/T)*tan(wp2/2);
Ws1=(2/T)*tan(ws1/2); Ws2=(2/T)*tan(ws2/2);
ws=[Ws1,Ws2];wp=[Wp1,Wp2];%通带宽度
H1=10^(-Rp/20);%通带波纹幅度 Rp=-20lg(H1)
H2=10^(-As/20);%阻带波纹幅度 As=-20lg(H2)
[N,wc]=cheb1ord(wp,ws,Rp,As,'s') %计算滤波器阶数和3dB截止频率wc
[B,A] =cheby1(N,Rp,wc,'s'); %计算滤波器系统函数的系数
[Bz,Az]=bilinear(B,A,Fs); %用双极性变换法转换成数字滤波器
[H,w]=freqz(Bz,Az); %计算H(z)的幅频响应
subplot(2,2,1);plot(w/2/pi*Fs,abs(H),'k');
xlabel('频率(Hz)'); ylabel('幅度(dB)') ;axis([0,500,0,1.1]);
set(gca,'XTick',[0,ws1*Fs/(2*pi),wp1*Fs/(2*pi),wp2*Fs/(2*pi),ws2*Fs/(2*pi),Fs/2]);
set(gca,'YTick',[0,H2,H1,1]); title('切比雪夫I型数字带通滤波器');grid on
subplot(2,2,2),plot(w/2/pi*Fs,angle(H)/pi*180,'k');
ylabel('\phi');title('相位函数曲线');axis([0,Fs/2,-180,180]);
set(gca,'XTick',[0,ws1*Fs/(2*pi),wp1*Fs/(2*pi),wp2*Fs/(2*pi),ws2*Fs/(2*pi),Fs/2]);
set(gca,'YTick',[-180,0,180]);grid on
subplot(2,2,3);myplot(Bz,Az,As);grid on;set(gca,'XTick',[0,ws1/pi,wp1/pi,wp2/pi,ws2/pi]);
4、运行结果及相应波形
N =
3
wc =
1.0e+003 *
1.0191 3.9252
（4） 采用双极性变换法设计切比雪夫I型数字带阻滤波器，要求：下通带wp1=0.2
，上通带wp2=0.8
，Rp=1dB，阻带下限ws1=0.3
，阻带上限ws2=0.7
，As=20dB，滤波器的采样频率Fs=1000Hz。
设计步骤：
1、用MATLAB信号处理工具箱函数cheb1ord，cheby1是切比雪夫I型滤波器设计函数，该程序中调用以下格式：
[N,wc]=cheb1ord(wp,ws,Rp,As,’s’)
该格式用于计算切比雪夫I型数字滤波器的最低阶数N和通带边界频率wc，
[B,A]=cheby1(N,Rp,wc,’stop’,’s’)
该格式用于计算切比雪夫I型数字滤波器系统函数分子和分母多项式的字数向量B,A。调用参数N,和wc分别为切比雪夫I型数字滤波器的阶数和通带边界频率wc。stop表示为设计带阻滤波器。
2、程序实现代码如下：
%双极性变换法设计IIR切比雪夫I型数字带阻滤波器
close all;clear all;
%=============要求设计指标===================
wp1=0.2*pi; wp2=0.8*pi; ws1=0.3*pi; ws2=0.7*pi;
Rp=1; As=20; T=0.001;Fs=1/T;
%============================================
Wp1=(2/T)*tan(wp1/2); Wp2=(2/T)*tan(wp2/2);
Ws1=(2/T)*tan(ws1/2); Ws2=(2/T)*tan(ws2/2);
ws=[Ws1,Ws2];wp=[Wp1,Wp2];%通带宽度
H1=10^(-Rp/20);%通带波纹幅度 Rp=-20lg(H1)
H2=10^(-As/20); %阻带波纹幅度 As=-20lg(H2)
[N,wc]=cheb1ord(wp,ws,Rp,As,'s');%计算滤波器阶数和3dB截止频率wc
[B,A] =cheby1(N,Rp,wc,'stop','s');%计算滤波器系统函数的系数
[Bz,Az]=bilinear(B,A,Fs);%用双极性变换法转换成数字滤波器
[H,w]=freqz(Bz,Az);%计算H(z)的幅频响应
subplot(2,2,1),plot(w/2/pi*Fs,abs(H),'k');
xlabel('频率(Hz)'); ylabel('幅度(dB)') ;axis([0,Fs/2,0,1.1]);
set(gca,'XTick',[0,100,150,350,400,500]); set(gca,'YTick',[0,H2,H1,1]);
title('切比雪夫I型数字带通滤波器');grid on
subplot(2,2,2);plot(w/2/pi*Fs,angle(H)/pi*180,'k');
ylabel('\phi');title('相位函数曲线');axis([0,Fs/2,-180,180]);
set(gca,'XTick',[0,100,150,350,400,500]);%获取当前X轴句柄值
set(gca,'YTick',[-180,0,180]);grid on
subplot(2,2,3);myplot(Bz,Az,As);grid on;set(gca,'XTick',[0,0.2,0.3,0.7,0.8]);
3、运行结果及相应波形
N =
3
wc =
1.0e+003 *
0.6498 6.1554
题目2、窗函数法设计数字FIR滤波器
（1） 用矩形窗设计一个FIR数字低通滤波器，要求：N=64,截止频率为wp=0.45
，描绘理想和实际滤波器的脉冲相应，窗函数及滤波器的幅频特性响应曲线。
设计步骤：
1、用MATLAB信号处理工具箱函数fir1是用于设计线性相位FIR数字滤波器的工具箱函数，实现线性相位FIR数字滤波器的标准窗函数法设计，调用格式如下：
hn=fir1(M,wc,’low’,boxcar(N))
用于计算6dB截止频率为wc的M阶单位脉冲响应hn。
理想的单位脉冲响应函数：hd=sin(wc*(n-a)/(pi*(n-a)));
b=boxcar(N);%返回长度为N的矩形窗函数b
开始
2、滤波器设计流程图
读入窗口的长度
计算hd(n)
调用窗函数子程序求b
计算h(n)=hd(n)*b’
调用子程序(函数)计算理想和实际的H，w，pha
调用子程序(函数)绘出理想和实际的H，pha，损耗曲线
结束
3、程序实现代码如下：
%用矩形窗设计FIR数字低通滤波器
close all;clear all;
%============要求设计指标=========
N=64;wc=0.4*pi;
%=================================
n=0:N-1;a=(N-1)/2;hd=sin(wc*(n-a)/(pi*(n-a)));
b=boxcar(N);%返回长度为N的矩形窗函数b
hn=hd.*(b)';%计算理想滤波器的单位脉冲响应
[H,w]=freqz(hn);%计算理想H(z)的幅频响应
pha=angle(H);%相频特性
subplot(3,2,1);stem(n,hn,'.');xlabel('n');ylabel('h(n)');
title('理想低通滤波器脉冲响应');
subplot(3,2,2);hn1=fir1(N-1,wc/pi,'low',boxcar(N));%计算hn
stem(n,hn1,'.');axis([0,70,-0.1,0.4]);grid on
title('实际低通滤波器脉冲响应');xlabel('n');ylabel('hn1');
subplot(3,2,3);plot(w,abs(H));grid on ;set(gca,'XTick',[0,0.4]);
title('理想幅频响应曲线');xlabel('\omega/\pi');ylabel('幅值/dB');
subplot(3,2,4);[H1,w1]=freqz(hn1);%计算实际H(z)的幅频响应
plot(w/pi,abs(H1));axis([0,1,0,1.1]);grid on;title('实际幅频响应曲线');
xlabel('\omega/\pi');ylabel('幅值/dB');set(gca,'XTick',[0,0.4]);
subplot(3,2,5);stem(n,b,'.');axis([0,70,0,1.2]);
title('矩形窗函数曲线');xlabel('n');ylabel('幅度');subplot(3,2,6)
A=1;myplot(hn1,A);set(gca,'XTick',[0,0.4]);
4、运行结果及相应波形
（2） 用合适的窗函数设计一个FIR数字高通滤波器，要求：通带截止频率为wp=0.45
，Rp=0.5dB，阻带截止频率为ws=0.3
，As=20dB，描绘该滤波器的脉冲响应，窗函数机滤波器的幅频响应曲线和相频响应曲线。
设计步骤：
1、用MATLAB信号处理工具箱函数fir1是用于设计线性相位FIR数字滤波器的工具箱函数，实现线性相位FIR数字滤波器的标准窗函数法设计，调用格式如下：
hn=fir1(M,wc,’high’, bartlett (N))
用于计算6dB截止频率为wc的M阶单位脉冲响应hn。
b= bartlett (N);%返回长度为N的三角形窗函数b
N1=ceil(6.1*pi/Bt);%计算三角形窗hn所需长度
2、程序实现代码如下：
%用三角形窗设计FIR数字高通滤波器
close all;clear all;
%==========要求设计指标================
wp=0.45*pi;ws=0.3*pi;As=20;
%======================================
%由As=20dB，查表得知使用三角形窗函数符合设计要求
wc=(wp+ws)/2; %计算理想高通滤波器的通带截止频率
Bt=wp-ws; %过渡带宽度
N1=ceil(6.1*pi/Bt);%计算三角形窗hn所需长度
N=N1+mod(N1+1,2)%如果N为偶数加1，保证N为奇数
hn=fir1(N-1,wc/pi,'high',bartlett(N));%计算hn
subplot(3,2,1);yn='hn';tstem(hn,yn);%绘制频率响应hn波形
title('三角形高通滤波器频率响应函数');grid on;
subplot(3,2,2);A=1;myplot(hn,A)%绘制损耗函数波形
set(gca,'XTick',[0,ws/pi,wp/pi]); %获取当前X轴句柄值
subplot(3,2,3);[H,w]=freqz(hn);plot(w/pi,abs(H));
set(gca,'XTick',[0,ws/pi,wp/pi]);grid on; ylabel('幅值/dB');
axis([0.25,0.8,0,1.2]);title('幅频响应曲线');xlabel('\omega/\pi');
subplot(3,2,4);plot(w/pi,angle(H));
set(gca,'XTick',[0,ws/pi,wp/pi]); axis([0,1,-3,3]);title('相频响应曲线');
grid on;xlabel('\omega/\pi');ylabel('相位/rad');
subplot(3,2,5);wn=bartlett(N);%返回长度为N的三角形窗函数wn
n=0:N-1;stem(n,wn,'.');axis([0,45,0,1.2]);title('三角形窗函数曲线');
xlabel('n');ylabel('幅度');grid on;
3、运行结果及相应波形
N =
41
（3） 选择合适的窗函数设计一个FIR数字带通滤波器，要求：下通带截止频率ws1=0.2
，As=65dB，通带低端截止频率wp1=0.3
，Rp=0.05dB；通带高端截止频率wp2=0.8
，As=65dB，描绘实际滤波器的脉冲响应，窗函数机滤波器的幅频响应曲线机相频响应曲线。
设计步骤：
1、用MATLAB信号处理工具箱函数fir1是用于设计线性相位FIR数字滤波器的工具箱函数，实现线性相位FIR数字滤波器的标准窗函数法设计，调用格式如下：
hn=fir1(M,wc,blackman(N)); 用于计算6dB截止频率为wc的M阶单位脉冲响应hn。
b=blackman(N);%返回长度为N的布莱克曼窗函数b
N1=ceil(11*pi/Bt);%计算布莱克曼窗hn所需长度
2、程序实现代码如下：
%用布莱克曼窗设计FIR数字带通滤波器
close all;clear all;
%====================要求设计指标============
wp1=0.3*pi;wp2=0.7*pi;ws1=0.2*pi;ws2=0.8*pi;As=65;
%============================================
%由As=65dB，查表得知使用布莱克曼窗函数符合设计要求
wc=[(wp1+ws1)/2,(wp2+ws2)/2];
Bt=wp1-ws1;%过渡带宽度
N=ceil(11*pi/Bt)%计算布莱克曼窗hn所需长度
hn=fir1(N-1,wc/pi,blackman(N));%计算hn
subplot(3,2,1);yn='hn';tstem(hn,yn);%绘制hn波形
title('布莱克曼带通滤波器频率响应函数');grid on;
subplot(3,2,2);A=1;myplot(hn,A);
set(gca,'XTick',[0,ws1/pi,wp1/pi,wp2/pi,ws2/pi]);
subplot(3,2,3);[H,w]=freqz(hn);
plot(w/pi,abs(H));axis([0,1,0,1.1]);grid on
set(gca,'XTick',[0,ws1/pi,wp1/pi,wp2/pi,ws2/pi]);
title('幅频响应曲线');xlabel('w');ylabel('幅值/dB');
subplot(3,2,4);plot(w/pi,angle(H));axis([0,1,-3,3]);grid on;
set(gca,'XTick',[0,ws1/pi,wp1/pi,wp2/pi,ws2/pi]);
xlabel('w');ylabel('相位/rad');title('相频响应曲线');
subplot(3,2,5);wn=blackman(N);%返回长度为N的布莱克曼窗函数wn
n=0:N-1;stem(n,wn,'.');axis([0,120,0,1.2]);
title('布莱克曼窗函数曲线');xlabel('n');ylabel('幅度');grid on;
3、运行结果及相应波形
N =
110
（4） 选择合适的窗函数设计一个FIR数字带阻滤波器，要求：下通带截止频率wp1=0.2
，
Rp=0.1dB；阻带低端截止频率ws1=0.3π，As=40dB，阻带高端截止频率ws2=0.7π，As=40dB，上通带截止频率wp2=0.8π，Rp=0.1dB；描绘实际滤波器的脉冲响应，窗函数机滤波器的幅频响应曲线及相频响应曲线。
1、用MATLAB信号处理工具箱函数fir1是用于设计线性相位FIR数字滤波器的工具箱函数，实现线性相位FIR数字滤波器的标准窗函数法设计，调用格式如下：
hn=fir1(M,wc, 'stop',hamming(N)); 用于计算6dB截止频率为wc的M阶单位脉冲响应hn。
b=hamming(N);%返回长度为N的哈明窗函数b
N1=ceil(6.6*pi/Bt);%计算哈明窗hn所需长度
2、程序实现代码如下：
%用哈明窗设计FIR数字带阻滤波器
close all;clear all;
%====================要求设计指标================
wp1=0.2*pi;wp2=0.8*pi;ws1=0.3*pi;ws2=0.7*pi;As=65;
%================================================
%由As=40dB，查表得知使用哈明窗函数符合设计要求
Bt=ws1-wp1;%计算理想高通滤波器的通带截止频率
N=ceil(6.6*pi/Bt)-1%计算哈明窗hn所需长度
wp=[(wp1+ws1)/2,(wp2+ws2)/2];
hn=fir1(N-1,wp/pi,'stop',hamming(N)');%计算hn
subplot(3,2,1);yn='hn';tstem(hn,yn);%绘制hn波形
title('哈明带阻滤波器频率响应曲线');grid on;
subplot(3,2,2);A=1;myplot(hn,A)%绘制损耗函数波形
set(gca,'XTick',[0,wp1/pi,ws1/pi,ws2/pi,wp2/pi]);
subplot(3,2,3);[H,w]=freqz(hn);plot(w/pi,abs(H));axis([0,1,0,1.1]);grid on
set(gca,'XTick',[0,wp1/pi,ws1/pi,ws2/pi,wp2/pi]);
title('幅频响应曲线');xlabel('w');ylabel('幅值/dB');
subplot(3,2,4);plot(w/pi,angle(H));
set(gca,'XTick',[0,wp1/pi,ws1/pi,ws2/pi,wp2/pi]); axis([0,1,-3,3]);grid on;
title('相频响应曲线');xlabel('w');ylabel('相位/rad');
subplot(3,2,5);wn=hamming(N);%返回长度为N的哈明窗函数wn
n=0:N-1;stem(n,wn,'.');;axis([0,70,0,1.2]);title('哈明窗函数曲线');
xlabel('n');ylabel('幅度');grid on;
3、运行结果及相应波形
N =
65
浙公网安备 33021202002483