数字信号处理课程设计IIR数字滤波器以实现

一、课程设计目的：

1.全面复习课程所学理论知识，巩固所学知识重点和难点，将理论与实践很

好地结合起来。

2.掌握信号分析与处理的基本方法与实现

3.提高综合运用所学知识独立分析和解决问题的能力；

4.熟练使用一种高级语言进行编程实现。

二、课程设计题目和具体设计（分题列出）；

（五）第五题

L设计题目及设计要求

设有一信号N〃）=l+cos（£n*os（空n）,设计各种HR数字滤波器以实现：

43

1）低通滤波器，滤除cos（2n）的成分，保留成分1+cos（工n）

34

2）高通滤波器，滤除l+cos（?n）的成分，保留成分cos（等n）

3）带通滤波器，滤除1+cos（二n）的成分，保留成分cos（工n）

34

4）带阻滤波器，滤除cos（工n）的成分，保留成分l+cos（生n）

43

要求：1）求出各个滤波器的阶数，设计各滤波器。画出各滤波器的幅频和相

频特性，计算滤波器的系统函数H（z）

2）画出滤波前后信号的时域、频域波形

2.设计思想及系统功能分析

本题主要考察四种基本HR数字滤波器的设需要先对时域信号进行分析，

确定各种滤波器的参数指标，然后根据matlab的一些固定调用格式对滤波器的

阶数，截止频率等进行设计最后再生成最终滤波器，然后利用生成的滤波器对信

号进行滤波。

3.关键部分的理论分析与计算

关键部分在于分析滤波器的参数指标时，要将频率都进行归一化，因为此处

设计的是数字滤波器；另外，通带衰减一般小于3dB,阻带衰减一般大于30dB；

设计带通、带阻滤波器时其截至频率分为上限截至频率和下限截至频率，书写方

式为区间形式；另外，此处根据信号形式选择通、阻带截频时，应注意在所给信

号频率的附近根据要求作适当的调整，不能刚好将所给的频率作为截止频率，因

为滤波器不可能完全理想，为了保证信号能尽量无失真的滤波出来，就要通带范

围稍微大一点。

4.程序源代码

%第五题

机1)设计数字低通滤波器

figurc(l);

T=1

wpl=O.4;wsl=0.5;apl=l;asl=30;先给出数字参数指标

[Nl,wc1]=buttord(wp1,wsl,apl,asl);%计算阶数N和3Db处的频率

N1%在命令框中显示N

[Bl,Al]=butter(Nl,wcl);%求巴特沃斯型数字滤波器

w=linspace(0,2*pi,1000);%在0-2*pi内将w分成1000个点

hl=freqz(Bl,Al,w);%调用freqz函数,求解幅频响应

maghl=abs(hl);%求出幅频特性

phalsingle(hl)外求出相频特性

subplot(321);plot(w/pi,maghl,'r');

与画图

title('低通滤波器幅频特性曲线')；xlabel('模拟频率');

subplot(322);plot(v/pi,phal,'r');

titleC低通滤波器相频特性曲线');xlabel('模拟频率');

n=0:72;%选取序列长度

w=pi*n/73;%分出与r.相对应的w

xn-l+cos(pi*n/4)+cos(2*pi*n/3);与原疗列

subplot(323);stem(n,xn,'.');%画出原序列

title(*xn的时域曲线')；xlabel('n');ylabelCxn);

X=abs(xn*exp(-l*j*n，*w));%原序列的序列傅里叶变换

subplot(3,2,4);stem(w,X,'.');%画出序列傅里叶变换的曲线

title(*xn的频域曲线');xlabelCwf);ylabel('X');

yl=filtcr(Bl,Al,xn);%利用filter函数恢复出滤波后的函数

subplot(325);stem(n,yl,'.nf);先画图

title('低通滤波后的时域曲线')jxlabelCn,)jylabelCyl');

Yl=abs(yl*exp(-l*j*n，*w));%恢复后的函数的序列傅里叶变换

subplot(3,2,6);stem(w,Yl,'.m')一画图

title('低通滤波后的频域曲线');xlabclCw');ylabel(,YT);

%(2)设计数字高通滤波器

figure(2);

wp2=0.5;ws2=0.4;ap2=l;as2=30;

[N2,wc2]=buttord(wp2,ws2,ap2,as2);

N2

[B2,A2]=butter(N2;wc2,'high*);

w=1inspace(O,2*pi.1()00);

h2=freqz(B2,A2,w);

magh2=abs(h2);

pha2=angle(h2);

subplot(321);plot(w/pi,magh2,'r');

titleC高通滤波器幅频特性曲线')；xlabel('模拟频率');

subplot(322);plot(w/pi,pha2,'r');

titleC高通滤波器相频特性曲线')ixlabelC模拟频率');

n-0:72;

w=pi*n/73;

xn=1+cos(pi*n/4)+cos(2*pi*n/3);

subplot(323);stem(n,xn,'.');

titleCxn的时域曲线');xlabel('n');ylabel(*xn);

X=abs(xn*exp(-1*j*n,*w));

subplot(3,2,4);stem(w,X,'.');

titleCxn的频域曲线')jxlabelCw*);ylabelCX*);

y2=fi1ter(B2,A2,xn);

subplot(325);stem(n,y2,'.m');

titleC高通滤波后的时域曲线')ixlabelCn));ylabel('y2');

Y2=abs(y2*exp(-l*j*n，*w));

subplot(3,2,6);stem(w,Y2,'.m');

titleC高通滤波后的频域曲线');xlabel(*w*);ylabel(，Y2*);

%(3)设计数字带通滤波器

figure(3);

wp3=[0.1,0.4];ws3=[0.05,0.5];ap3=l;as3=30;

[N3,wc3]=buttord(wp3,ws3,ap3,as3);

N3

[B3,A3]=butter(N3:wc3);

w=linspace(0,2*pi.1000);

h3=freqz(B3,A3,w);

magh3=abs(h3);

pha3=angle(h3);

subplot(321);plot(w/pi,magh3,'r*);

titleC带通滤波器幅频特性曲线');xlabelC模拟频率')；

subplot(322);plot(w/pi,pha3,*r*);

titleC带通滤波器相频特性曲线')ixlabel('模拟频率')；

n=0:72;

w=pi*n/73;

xn=l+cos(pi*n/4)+cos(2*pi*n/3);

subplot(323);stem(n,xn,'.');

titleCxn的时域曲线')jxlabelCn*);ylabelCxn);

X-abs(xn*exp(-1+j*n，+w));

subplot(3,2,4);stem(w,X,'.');

titleCxn的频域曲线');xlabel(,w');ylabel(?XJ);

y3=filter(B3,A3,xn);

subplot(325);stem(n,y3,'.m');

titleC带通滤波后的时域曲线');xlabeln)jylabelCy3');

Y3=abs(y3*exp(-l*j*n，*w));

subplot(3,2,6);stem(w,Y3,'.nf);

titleC带通滤波后的频域曲线')；xlabelCw');ylabelfY3,);

%(4)设计数字带阻滤波器

figure(4);

wp4=[0.05,0.5];ws4=[0,1,0.4];ap4=l;as4=30;

[N4,wc4]=buttord(wp4,ws4,ap4,as4);

N4

[B4,A4]=butter(N4.wc4,*stop*);

w=linspace(0,2*pi:1000);

h4=freqz(B4,A4,w);

magh4=abs(h4);

pha4=angle(h4);

subp1ot(321);plot(w/pi,magh4,*r*);

titleC带阻滤波器幅频特性曲线')ixlabelC模拟频率');

subplot(322);plot(v/pi,pha4,*r*);

titleC带阻滤波器相频特性曲线')；xlabel('模拟频率');

n=0:72;

w=pi*n/73;

xn=l+cos(pi*n/4)+cos(2*pi*n/3);

subplot(323);stem(n,xn,'.');

titleCxn的时域曲线');xlabel(*n');ylabelCxn);

X=abs(xn*cxp(T*j*n'*w));

subplot(3,2,4);stem(w,X,'.');

titleCxn的频域曲线')"label('w');ylabel('X');

y4=filter(B4,A4,xn);

subplot(325);stem(n,y4,'.m');

titleC带阻滤波后的时域曲线');xlabel(JnJ);ylabel(*y4');

Y4=abs(y4*exp(-l*j*n,*w));

subplot(3,2,6);stem(w,Y4,'.m');

titleC带阻滤波后的频域曲线');xlabelCw*);ylabelCY4J);

5.测试数据及必要的理论分析与比较

由figurel可以看出低通滤波器的幅频特性曲线和相频特性I川线，原信号的

时域曲线和频域曲线以及滤波后的信号的时域曲线和频域曲线，通过对比，可以

很清楚的看出该低通滤波器实现了题目所耍求的滤波功能；

由figure2可以看出高通滤波器的幅频特性曲线和相频特性曲线，原信号的

时域曲线和频域曲线以及滤波后的信号的时域曲线和频域曲线，通过对比，可以

很清楚的看出该高通滤波器实现了题目所要求的港波功能；

由figure3可以看出带通滤波器的幅频特性曲线和相频特性曲线，原信号的

时域曲线和频域曲线以及流波后的信号的时域曲线和频域曲线，通过对比，可以

很清楚的看出该带通滤波器实现了题目所要求的滤波功能；

由figurc4可以看出带阻滤波器的幅频特性曲线和相频特性曲线，原信号的

时域曲线和频域曲线以及滤波后的信号的时域曲线和频域曲线，通过对比，可以

很清楚的看出该带阻滤波器实现了题目所要求的滤波功能；

(六)第六题

1.设计题目及设计要求

⑴用Hanning窗设计一线性相位带通数字滤波器，要求：N=15,

叫・°,为，叼.05*。观察它的实际3dB和20dB带宽。N=45,重复这一设

计，观察幅频和相位特性的变化，注意长度N变化的影响；

(2)分别改用矩形窗和Blackman窗，设计(1)中的昔通滤波器，观察并记录窗函数

对滤波器幅频特性的影响，比较三种窗的特点；总结窗的不同长度和不同窗对

滤波器的影响

2.设计思想及系统功能分析

本题主要考察利用窗函数法设计线性相位的数字滤波器，与之前的设计方法

并无太大区别，先是确定滤波器的技术指标(归一化)，然后通过调用窗函数的

语句进行设计，只不过这里需要根据窗函数的不同选择相应的过渡带宽及选择信

号的长度和滤波器的除数；同样，带通滤波器的截止频率需要用区间形式写出；

对于不同的窗函数，过渡带宽会不同。

3.关键部分的理论分析与计算

关键部分应该有确定截止频率(注意要归一化)，还有选择不同的窗函数对

应的过渡带宽；另外就是阶数和信号长度的关系：信号长度二阶数+1,这是从资

料书上看来的，所以每个滤波器的阶数都是按照这个等式来确定的。

4.程序源代码

舟笫六题

%(1)用汉宁窗设计线性相位带通滤波器

舟当长度N1为15时

Nl=15;

wl=O.3*pi;w2=0.5*pi;

wc=[().3,0.5]$归一化后的通带截频

Bl=8*pi/Nl;与过渡带宽

Ml=N1-1;%单位脉冲响应的h(n)的长度N=M+1,此处的Ml是阶数；

nl-0:Ml;

hnl=firl(Ml,wc,hanning(Nl));为调用firl函数求解单位脉冲响应

[hl,wl]=freqz(hnl,1,512);%求其频率响应

maghl=abs(hl);%求其幅频特性

phal=angle(hl);%求其相频特性

figure(l);

subplot(221);plot(vl/pi,10*log(maghl));%画出其幅频特性

grid;

titlefN=15时,汉宁窗的幅频特性');xlabelfw*);ylabelfmagh*);

subplot(222);plot(wl/pi,phal);%画出其幅频特性

title。N=15时,汉宁窗的相频特性');xlabel(*w*);ylabel(*pha*);

grid;

舟当长度N2为45时

N2=45;

wc=[0.3,0.5]归一化后的通带截频

B2=8*pi/N2圈过渡带宽

M2=N2T;%单位脉冲响应的h(n)的长度N=M+1,此处的Ml是阶数；

n2=0:M2;

hn2=firl(M2,wc,hanning(N2));为调用firl函数求解单位脉冲响应

[h2,w2]=freqz(hn2,1,512);%求其频率响应

magh2=abs(h2);%求其幅频特性

pha2=angle(h2);%求其相频特性

subplot(223);plot(w2/pi,10*log(magh2),*g');%画出其幅频特性

grid;

title(*N=45时，汉宁窗的幅频特性');xlabelCw,);ylabel('magh*);

subplot(224);plot(w2/pi,pha2,'g*);%画出其幅频特性

titleCN=45时,汉宁窗的相频特性');xlabel(*);ylabel(Jpha');

grid：

%(2)用矩形窗设计线性相位带通滤波器

%当长度N3为15时

figure(2);

N3=15;

wc=[0.3,0.5];

B3-4*pi/N3；%过渡带宽

M3=N3-1;

n3=0:M3;

hn3=firl(M3,wc,boxcar(N3));

[h3,w3]=freqz(hn3,1,512);

magh3=abs(h3);

pha3=angle(h3);

subplot(221);plot(v3/pi,10*log(magh3));

grid;

title(,N=15时,矩形窗的幅频特性')ixlabelCw*);ylabel('magh*);

subplot(222);plot(w3/pi,pha3);

grid;

title(,N=15时,矩形窗的相频特性')ixlabeK*w*);ylabel(*pha');

%当长度N2为45时

wc=[0.3,0.5];

B4=4*pi/N4*过渡带宽

M4=N4-1;

n4=0:M4;

hn4=firl(M4,wc,boxcar(N4));

[h4,w4]=freqz(hn4,1,512);

magh4=abs(h4);

pha4=angle(h4);

subplot(223);plot(v4/pi,10*log(magh4),1g,);

grid;

title(*N=45时,矩形窗的幅频特性');xlabelCw*);ylabel(,magh,);

subplot(224);plot(w4/pi,pha4,'g');

grid：

title。N=45时,矩形窗的相频特性');xlabel(*w*);ylabel(*pha');

$(3)用Blackman窗设计线性相位带通滤波器

%当长度N3为15时

figure(3);

N5=15;

wc-[0.3,0.5];

B5=12*pi/N5;%过渡带宽

M5=N5-1;

n5=0:M5;

hn5=firl(M5,wc,blackman(N5));

[h5,w5]=freqz(hn5,1,512);

magh5=abs(h5);

pha5-angle(h5);

subplot(221);plot(v5/pi,10*log(magh5));

grid;

title(*N=15时，Blackman窗的幅频特性');xlabelw');ylabel(?magh');

subplot(222);plot(w5/pi,pha5);

grid;

titlefN=15时，Blackman窗的相频特性‘);xlabelw');ylabel(?pha');

%当长度N2为4b时

N6=45;

wc=[0.3,0.5];

B6=12*pi/N6;为过渡带宽

M6=N6-1;

n6=0:M6;

hn6=fir1(M6,wc,b1ackman(N6));

[h6,w6]=freqz(hn6,1,512);

magh6=abs(h6);

pha6=angle(h6);

subplot(223);plot(w6/pi,10*log(magh6),*g*);

grid;

ti11e(*N=45时,Blackman窗的幅频特性');xlabel[w*);ylabelCmagh*);

subplot(224);plot(w6/pi,pha6,'g');

grid;

title('N=45时，Blackman窗的相频特性’)；xlabcl['w');ylabel('pha,);

5.测试数据及必要的理论分析与比较

由figurel可以看出，N=15时的汉宁窗设计的滤波器的通带比N=45时的

通带宽，旁瓣的衰减也要小一些，所以，阶数越大，所设计出来的滤波器效果越

好,通带稳定，衰减快。同理，根据fugure2和figure3的矩形窗和Blackman

窗也可以得出相同的结论。

通过对三个图不同窗函数设计的滤波器进行比较，可以看出：矩形窗设计的

滤波器通带最稳定，衰减最快；汉宇窗次之；Blackman窗滤波效果相对最差。

所以，矩形窗更能较好的实现滤波效果。

(七)选做题7

L设计题目及设计要求

音乐信号处理：

1)获取一段音乐或语音信号，设计单回声滤波器，实现信号的单回声产生。

给出加入单回声前后的信号频谱。

2)设计多重回声滤波器，实现多重回声效果。给出加入多重回声后的信号频

谱。

3)设计均衡器，使得得不同频率的混合音频信号，通过一个均衡器后，增强

或削减某些频率区域。

2.设计思想及系统功能分析

设计思想为：先从网上下载一段格式为wav的音乐信号，选择一个较长的

长度，设计一个滤波器，然后用该滤波器对音乐信号进行滤波，实质上就是将信

号拉长后，与原信号就行一次叠加，便可设计出单回声滤波器；给滤波器加入一

个阶次，让滤波器进行多次滤波然后再与原信号叠加，便可得到多重回声滤波器。

3.关键部分的理论分析与计算

设计单回声滤波器的关键在于该滤波器的系统函数的分子部分，要通过补零

的方式将滤波器进行拉长，最好还要加入一个固定常数项；设计多重回声滤波器

时，系统函数分子部分补的零根据阶数的不同会相应增加，最后的常数项也要变

为相应的阶次累，而且分母部分也应该进行补零。

4.程序源代码

%第七题

与(1)单回声滤波器

%addpathD:ProgameFi1es/MATLAB/R2011a/bin);

[x,fs,bits]=wavread(，D:\Prettyboy.wav',2cl6);%读出音频文件

%[x,fs,bits]=wavreadCD:\yue.wav',216);先读出音频文件

wavplay(x,fs);先播放音频文件

fs

hits

pause(l)；

al=0.8;

RI=5000;

Bl=[l,zeros(1,0.8];%系统函数分子

Al=[l]；

y=filter(Bl,Al,x)N滤波器函数

wavplay(y,fs);

paused)；

xl=fft(x);

figure(l)；

subplot(321);plot(abs(xl));

title('单回声前信号幅频响应')；

subplot(322);plot(angle(xl));

titleC单回声前信号相频响应')：

yl二fft(y);%快速傅里叶变换

subplot(323);plot(abs(yl),*m*);

title('单回声后信号幅频响应')；

subplot(324);plot(angle(yl),1m*)；

title(＞单回声后信号相频响应');

[hl,wl]=freqz(Bl,AD;%求滤波器的幅频响应函数

subplot(325);plot(abs(hl),*g*);

titleC单回声滤波器幅频响应');

y2=ifft(hl);

subplot(326);stem(abs(y2),*.g*);

titleC单回声源波器冲激响应')；

与多重回声滤波器

pause(2);

[x,fs,bits]=wavread(,D:\Prcttyboy.wav',2"16);%读出音频文件

%[x,fs,bits]二wavread,D:\yue.wav',216);%读出音频文件

wavplay(x,fs)；%播放音频文件

pause(0);

N=5;

a2=0.8;

K2=8U00;

B2=[1,zeros(1,N*R2-1),-0.8*N];

A2=[l,zeros(1,R2-l),-0.8];

z=filter(B2,A2,x);

wavplay(z,fs);

zl=fft(z);

xl=fft(x);

figure(2);

subplot(321);plot(abs(xl));

titleC多重回声前信号幅频响应')；

subplot(322);plot(angle(xl));

title('多重回声前信号相频响应');

subplot(323);plot(abs(zl),'m');

titleC多重回声后幅频响应');

subplot(324);plot(angle(zl),*m*)；

titleC多重回声后相频响应');

[h2,w2]=freqz(B2,A2);

subplot(325);plot(abs(h2),'g');

titleC多重回声滤波器幅频响应')；

z2二ifft(h2);

subplot(326);stem(aba(z2),*.g*);

title('多重回声滤波器冲激响应');

5.测试数据及必要的理论分析与比较

由figurcl可以看出，声音信号通过单回声流波器后整个频率域上的幅值都

发生了一些变化，变得比以前稍人一些，这是因为滤波后的信号是原信号与回声

信号叠加起来的；同理，由figure2可以看出，声音信号通过多重回声滤波器后

整个频率域上的幅值也都发生了一些变化，变得比以前更大。

通过对比figurel和figure2可以看出，单回声滤波器的幅频响应幅度比相

应的多重回声滤波器的幅度小，因为多重回声滤波器是经过N次的单回声滤波器

滤波后叠加的效果；而且单回声滤波器的冲激响应只有一个值，多重回声滤波器

的冲激响应有N个值。

（八）选做题8

1.设计题目及设计要求

对混有噪声的音乐信号进行滤波：

1）在一段音乐信号中混入两个频率的正弦型干扰信号，利用FFT计算该混合

信号的频谱并确定干扰信号的频谱；

2）利用二阶带阻滤波器设计一个能滤出干扰信号的梳状滤波器；

2.设计思想及系统功能分析

设计思想为：先从网上下载一段格式为wav的音乐信号，通过FFT对其进

行频谱分析，确定其主要的频率成分，然后通过佶号相加的形式加入两个不同频

率的正弦信号，再对混入后的信号进行FFT,分析其频谱特性，确定干扰信号的

频谱；接下来就是设计滤波器对信号进行还原，此处没有实现梳状滤波器，只实

现了二阶带阻滤波器，对信号成功的进行了还原。

3.关键部分的理论分析与计算

关键部分为利用FFT对音乐信号进行频谱分析，确定了音乐信号的主要成分

为低频，然后加入两个高频的正弦信号便能很清楚的干扰到原信号；在设计二阶

带阻滤波器的时候，要尽量将阻带做的很窄，能够准确的滤除干扰信号。

4.程序源代码

现在一段音乐信号中混入两个频率的正弦型干扰信号，利用FFT计算该混合信号的频谱并确

定干扰信号的频谱；

[x,fs,bits]=wavread(，D:\Prettyboy.wav',216);先读出音频文件

%[x,fs,bits]=wavread(>D:\S131.wav',215)；%读出音频文件

%wavplay(x,fs);先播放音频文件

a=size(x');

N=216;%截取的信号长度

n=O:l:N-l;

t=n/fs;

fsl=7000;「s2=5000;%加入的两个正弦信号的频率

xl=x*+sin(2*pi*fsl*t)+sin(2*pi*fs2*t);%混入正弦信号后的信号

y=fft(xl);%对它做fft

f1=(0:length(y)-1)**fs/length(y)将曲线的横坐标单位设为Hz

figure(l);

subplot(211);plot(fl,abs(y),'m')圈画出幅频特性

titleC叠加正弦信号的音乐信号的幅频特性');xlabelCff);ylabelC幅度');

subplot(212);plot(fl,angle(y),Jm');%画出相频特性

titleC叠加正弦信号的音乐信号的相频特性');xlabelCfl');ylabelC相位');

姐殳计二阶带阻滤波器

Nl=2;%阶数为二阶

[Bl,Al]:butter(Nl,[2*4800/fs,2*5200/fs],'stop');分带阻滤波器的系统函数系数

[B2,A2]=butter(Nl,[2*6800/fs,2*7200/fs],Jstop*)带阻滤波器的系统函数系数

Bt=conv(Bl,B2)圈带阻滤波器的传递函数系数

At=conv(A1,A2);与带阻滤波器的传递函数系数

figure(2);

[h,w]=freqz(Bt,At);%求带阻滤波器的频率响应

subplot(211);plot(w/pi.abs(h),'r');%带阻滤波器的幅频响应

titleC梳状滤波器的幅频响应')；

subplot(212);plot(w/pi:angle(h),*r*);%带阻滤波器的相频响应

titleC带阻滤波器的相频响应')；

%利用带阻滤波器滤除信号中的干扰

z=filter(Bt,At,xl);

figure(3);

subplot(311);plot(t,x);title(5原音乐信号')；%原音乐信号

subplot(312);plot(t,xl);titleC混入干扰后的音乐信号');%混入干扰后的音乐信号

subplot(313);plot(t,z);title('滤波后的音乐信号');%滤波后的音乐信号

wavplay(xl,fs);与播放处理前音乐

pause(l)；

wavplay(z,fs);先播放处理后音乐

5.测试数据及必要的理论分析与比较

由figurel可以看出声音信号为低频，而且有两个明显的高频信号(即所加

的干扰)；figure2为没计的用来滤除高频干扰的带阻滤波器的幅频特性和相频

特性；figurc3可以很清楚的看出