利用matlab设计巴特沃斯低通滤波器

三峡大学课程设计报告专业班级 20091421课程数字信号处理课程设计学号 2009142116 学生姓名姜祥奔指导教师王露2012年 5 月平时成绩（20%）报告成绩（40%）答辩成绩（40%）总成绩数字信号处理课程设计实验一：用双线性变换法和脉冲响应不变法设计IIR 数字滤波器采用双线性变换法设计一个巴特沃斯数字低通滤波器，要求：通带截止频率100Hz ，阻带截止频率200Hz ，通带衰减指标Rp 小于2dB ，阻带衰减Rs 大于15dB ，滤波器采样频率Fs=500Hz 。

绘制频率响应曲线。

理论部分：（原理及设计过程）第一步：论ω和f 的关系及数字域性能的公式表示。

模拟频率与数字频率之间为线性关系，T f T fT T s,101,24-===Ω=πω为抽样周期，故 ππω4.0500110021=⨯⨯==c c kHz f 对应于 修正后)2t a n (21c T ω=Ω ππω8.0500120021=⨯⨯==st st kHz f 对应于 修正后)2tan(21st T ω=Ω 按衰减的定义2)()(log 204.0010≤πj j e H e H 15)()(log 208.0010≥πj j e H e H 设0=ω处频率响应幅度归一化为1,即1)(0=j e H ，则上两式变成2)(log 204.010-≥πj e H (1)15)(log 208.010-≤πj e H (2)这就是数字滤波器的性能指标的表达式。

2 下面把数字低通滤波器的性能要求转变为“样本”模拟低通滤波器的性能要求。

由T Ω=ω，按修正式)()(T jH e H a j ωω≈，设没有混叠效应（即混叠效应设计完成后再进行校验）则有πωωω≤Ω==),()()(j H T j H e H a a j (3)利用（3）式，由（1）、（2）式可写出模拟低通滤波器的指标为2)102(log 20)4.0(log 2021010-≥⨯=ππj H Tj H a a …….（4） 15)104(log 20)8.0(log 2021010-≤⨯=ππj H T jH a a …（5） 3 计算“样本”模拟低通滤波器所需的阶数N 及3dB 截止频率C Ω。

巴特沃思低通滤波器的幅度平方函数是N C a j H 22)(11)(ΩΩ+=Ω以分贝形式表示上式，即⎥⎦⎤⎢⎣⎡ΩΩ+-=ΩN C a j H 21010)(1log 10)(log 20 ……………（6） 把求出的性能指标关系（4）式、(5)式代入(6)式得2)(1log 102110-≥⎥⎦⎤⎢⎣⎡ΩΩ+-N C 15)(1log 102210-≤⎥⎦⎤⎢⎣⎡ΩΩ+-N C 先用等号来满足指标，可得N C21)(1ΩΩ+＝2.010 N C22)(1ΩΩ+＝5.110 解此两方程，得N=1.3709, N 是滤波器阶次，必须取整数，为了满足或超过给定指标，故应选取比求出的N 大一点的整数，故取N ＝2,代入通带条件（6-118a ）式，得C Ω＝1.038310⨯,即N ＝2, C Ω＝1.038310⨯ (7)显然，利用（7）式的N 和C Ω，阻带指标正好满足，而通带指标则低于2dB ，这正好对减少冲激响应不变法造成的频谱混叠效应是有利的。

程序部分：wp=100*2*pi; %数字滤波器的通带截止频率ws=200*2*pi; %数字滤波器的阻带截止频率Rp=2;As=15; %输入滤波器的通阻带衰减指标%转换为模拟滤波器指标Fs=500;T=1/Fs;wp1=wp*T;ws1=ws*T;Omgp=(2/T)*tan(wp1/2); %原型通带频率预修正Omgs=(2/T)*tan(ws1/2); %原型通阻频率预修正%模拟原型滤波器计算[n, Omgc]=buttord(Omgp, Omgs,Rp,As,'s'); %计算阶数n 和截止频率[z0,p0,k0]=buttap(n); %归一化原型设计[ba,aa]=zp2tf(z0,p0,k0); %把滤波器零极点转化为传递函数模型[ba1,aa1]=lp2lp(ba,aa, Omgc); %变换为模拟低通滤波器系数b,a%或者[ba1,aa1]=butter(n, Omgc,’s ’);%用双线性变换法计算数字滤波器系数[bd,ad]=bilinear(ba1,aa1,Fs) ; %用双线性变换法求数字滤波器系数b,a%求数字系统的频率特性[H,w]=freqz(bd,ad);plot(w*Fs/(2*pi),20*log10(abs(H)));ylabel('增益/dB');xlabel('频率（Hz ）');title('数字滤波器幅度响应');grid on;数字滤波器的频率响应图：实验二：用窗函数法设计FIR 数字滤波器选择合适的窗函数设计一个FIR 数字低通滤波器，要求：带通截止频率为Wp=0.2π，阻带截止频率Ws=0.3π，阻带衰减不小于40dB ，通带衰减不大于3dB 。

描绘滤波器的幅频响应曲线。

理论部分：（原理及设计过程）1．数字频率。

通带截止频率：π2.0=p w阻带截止频率：π3.0=s w阻带衰减：dB 402=δ2．求)(n h d 。

设)(jw d e H 为理想线性相位滤波器⎪⎩⎪⎨⎧≤-其它,0,)(c jw jwd w we e H τ 首先由所需低通滤波器的过渡带求理想低通滤波器的截止频率c Ω（由于c Ω为两个肩峰值处的频率中点，而由p Ω到s Ω之间的过渡带宽并非两个肩峰间的频率差，因而以下求出的c Ω有一定的近似），并将其转化为对应的数字频率为π25.0)(21=+=s p c w w w 由此可得[]⎰--⎪⎪⎩⎪⎪⎨⎧=≠--==c c w w c c n jw d n w n n w n dw e n h ππτττππτ,,)(sin )(121)()( 其中τ为线性相位所必需的移位，且21-=N τ。

3．求窗函数。

由阻带衰减dB 402=δ，查课本中的表7-3可选汉宁窗，其阻带最小衰减44dB 满足要求。

所要求的过度带宽（数字频域）π1.0=-=∆p s w w w 由于汉宁窗过渡带宽满足N w π2.6=∆ 所以622.6=∆=wN π 3121=-=N τ 4．求)(n h 。

由汉宁窗表达式)(n w 确定FIR 滤波器的)(n h 。

汉宁窗 )()]12cos(1[21)(n R N n n w N --=π )21()]21(sin[)(----=N n N n w n h c d π 所以 )(*)()(n w n h n h d =5．由)(n h 求得)(jw e H ，并检验各项指标是否满足要求。

如不满足要求，则要改变N ，或改变窗形状（或两者都改变），然后重新计算。

\程序部分：子函数1：function hd=ideal_lp(wc,M);alpha = (M-1)/2;n = [0:1:(M-1)];m = n - alpha +eps; % add smallest number to avoi divided by zerohd = sin(wc*m)./(pi*m);子函数2：function[db,mag,pha,w]=freqz_m4(b,a)[H,w]=freqz(b,a,1000,'whole');H=(H(1:1:501));w=(w(1:1:501));mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);主程序：clear allWp=0.2*pi;Ws=0.3*pi;tr_width=Ws-Wp; %过渡带宽度N=ceil(6.2*pi/tr_width)+1; %滤波器长度n=0:1:N-1;Wc=(Ws+Wp)/2; %理想低通滤波器的截止频率hd=ideal_lp(Wc,N); %理想低通滤波器的单位冲激响应w_ham=(hann(N))'; %汉宁窗h=hd.*w_ham; %截取得到实际的单位脉冲响应[db,mag,pha,w]=freqz_m4(h,[1]); %计算实际滤波器的幅度响应delta_w=2*pi/1000;Ap=-(min(db(1:1:Wp/delta_w+1))); %实际通带纹波As=-round(max(db(Ws/delta_w+1:1:501))); %实际阻带纹波subplot(221)stem(n,hd) %火柴棍图title('理想单位脉冲响应hd(n)')subplot(222)stem(n,w_ham)title('汉宁窗w(n)')subplot(223)stem(n,h)title('实际单位脉冲响应h(n)')subplot(224)plot(w/pi,db)title('幅度响应(dB)')axis([0,1,-100,10])FIR低通滤波器特性图：心得体会：不得不承认，这一次数字信号处理的课程设计过程是艰苦的，但同时又是收获丰富的。

因为我们是上学期修的数字信号处理，而课程设计却安排在下学期，经过半年的时间，很多学过的东西都变得模糊的。

刚选定课题时，自己甚至连题目的意思都看不懂，“什么是双线性变换法？”“什么是巴特沃思数字低通滤波器？”，一系列的问号在头脑中打转。

Matlab程序设计经验的缺乏也曾一度令自己失去持续做程题的勇气。

现在，虽然不能拍着胸脯说自己弄通了如何用双线性变换法设计巴特沃思数字滤波器，也不敢说完全通悟利用窗函数设计一个FIR 数字低通滤波器的方法，但庆幸的是自己在课程设计这一环节中学会了独立学习知识的方法，更体验到了从“迷惘——坚定——求索——明悟”的乐趣，也发现了以前学习数字信号处理时一直被自己忽视的许多问题。

如模拟频率f，模拟角频率Ω及数字频率ω三者的概念区别及相互换算问题，在这一次课程设计中就是必须弄清楚的问题，不然一切寸步难行。