[bz,az]=impinvar(bs,as,1/T); [H,w]=freqz(bz,az,w); plot(w/pi,abs(H));
10.90.80.70.60.50.40.30.20.1000.10.20.30.40.50.60.70.80.91
②双线性z变换法
编程实现：设计巴特沃斯低通滤波器，fp=100Hz，fs=300Hz，ap=3dB，as=20dB，抽样频率Fs=1000Hz，画出设计后的数字低通滤波器的幅频特性曲线。
提示：采用[N, OmgC]=buttord(wap,was,Rp,As,‘s’)、[bs,as]=butter(N, OmgC,‘s’)，说明：若Ω和w在映射的过程中乘了系数2/T，调用[Bz,Az]=bilinear(bs,as,Fs)，否则调用 [Bz,Az]=bilinear(bs,as,Fs/2)和[H,w]=freqz(Bz,Az,w)函数。 fs=1000; fp=100; fst=300; Rp=3; As=20; w=0:0.005*pi:pi;
wap=tan(fp*2*pi/fs/2)*fs*2; was=tan(fst*2*pi/fs/2)*fs*2; [N,OmgC]=buttord(wap,was,Rp,As,'s');
[bs,as]=butter(N,OmgC,'s'); [bz,az]=bilinear(bs,as,fs); [H,w]=freqz(bz,az,w); plot(w,abs(H));
10.90.80.70.60.50.40.30.20.1000.511.522.533.5
四．思考题
1.冲激响应不变法设计滤波器的步骤？优缺点？ 2.双线性z变换法设计滤波器的步骤？优缺点？
实验六 数字滤波器组的设计——FIR数字滤波器的设计
一、实验目的
1、掌握用窗函数法和频率采样法设计FIR数字滤波器的原理和方法； 2、熟悉线性相位FIR滤波器的幅频特性和相频特性； 3、了解不同窗函数对滤波器性能的影响。 4、了解过渡带采样点对滤波器性能的影响。 二、实验原理
?nH(z)?h(n)z?FIR数字滤波器的系统函数为，该滤波器总是因果稳定的系统，且
n?0N?1当满足h(n)?h(N?1?n)或h(n)??h(N?1?n)的对称条件时，FIR滤波器具有线性相位。其在数据通信、图像处理、语音信号处理等实际领域得到广泛应用。FIR数字滤波器的设计方法主要有窗函数法和频率采样法。
（1）窗函数法
滤波器设计是以理想滤波器的特性为基准，使待设计的实际滤波器的特性逼近理想滤波器的特性。以低通滤波器为例，一个截止频率为?c、相位响应为零的理想低通滤波器的频率响应为：
??1Hd(e)????0j??c????
???c求离散时间傅里叶反变换可得理想冲激响应为：
sin?cn1?1?cj?nj?j?nhd(n)??Hd(e)ed???ed?????n???
-?-?c2?2??n理想冲激响应hd(n)是时间上的无限长序列，信号的绝大部分能量集中在零时刻附近，随着时间的增加，能量越来越小，最后趋于零。因此，可以用一个窗函数对hd(n)进行截短，原则上尽量保留hd(n)的最高能量部分，以达到用时间上的有限长序列逼近理想冲激响应的目的，这就是窗函数法设计FIR数字滤波器的基本思想。
若窗函数用w(n)表示，则有：
hw(n)?hd(n)w(n)
从物理可实现的角度考虑，此时得到的有限长序列hw(n)并非因果序列，需要将其进行移位，最终得到FIR数字低通滤波器的冲击响应序列h(n)，即：
h(n)?hw(n?N)
若相位响应不为零，而是令?(?)????，即?(?)具有线性相位，则：
-j???e???c?N?1Hd(ej?)????，其中，N为滤波器阶数。 0?????2?c?1hd(n)?2???e-c?c?j??sin[?c(n??)]ed???(n??)
j?n窗函数法的设计思想是从时域出发，截取有限长的一段冲激响应作为H(z)系数，冲激响应长度N就是系统函数H(z)的阶数。只要N足够长，并且截取方法合理，总能够满足频域的要求。窗函数法设计FIR数字滤波器的基本步骤为：
Step1：确定数字滤波器理想特性Hd(e)；
j?Step2：由Hd(e)求出单位冲激响应hd(n)；
。
Step3：由过渡带宽及阻带最小衰减的要求选择适当的窗函数，并根据线性相位条件确定窗函数的长度N；
如h(n)偶对称，N为奇数，或h(n)奇对称，N为偶数。
2??3.3????p??s，对于海明窗，N0???，若选h(n)偶对称，N为奇数，在MATLAB中，可用下面语句实现：
deltaw=ws-wp;%求过渡带宽度
N0=ceil(6.6*pi/deltaw);%/窗宽由过渡带决定，课本342页表7-3， N=N+mod(N0+1,2);
wdham=(hamming(N))';%求窗函数
在MATLAB中，可由w=boxcar(N)（矩形窗，返回值w是一个长度为N的矩形窗序列）、w=hanning(N)（汉宁窗）、w=hamming(N)（汉明窗）、w=Blackman(N)（布莱克曼窗）、w=Kaiser(N,beta)（凯塞窗，最有用的窗结构之一，对于给定的阻带衰减，它能够提供最陡峭的过渡带，beta是窗函数的形状参数，通过改变beta可以对主瓣宽度）等函数来实现窗函数设计法中所需的窗函数。
Step4：求出单位冲激响应h(n)：h(n)=hd(n)w(n) 0≤n≤N-1
Step5：对h(n)作离散时间傅立叶变换，得到H(ej?)，分析幅频特性是否满足要求。
j?（2）频率采样法
频率采样法是先对理想频响
Hd(ej?)进行采样，得到采样值
H(k)，即
H(k)?Hdk(?)Hdej?(2???kN)k?,N0?,1，,再对频域序列,1H(k)做IDFT得到时域的滤波器
序列：h(n)?IDFT[H(k)]，这种设计过程是在频域进行抽样得到离散的频域点，称之为频率抽样法。
频率抽样法的设计思路：
Hd(ej?)频率取样H(k)?Hd(2k?/N)IDFTh(n)
FIR数字滤波器的频率抽样法设计思路图
频率采样法设计FIR数字滤波器的基本步骤为：
Step1：根据设计要求选择滤波器的种类，确定合适的滤波器阶数N[0,2?)； Step2：然后对其在上进行N点等间隔采样得到H(k)；
Step3：将H(k)带入内插公式得到所设计滤波器的频率响应H(ej?)；
在 MATLAB中可由函数h=real(ifft(H,N))和实现。频率采样法得到滤波器，在采样点上滤波器的实际频率响应是严格地和理想频率响应数值相等的。但是在采样点之间的频响则是由各采样点的加权内插函数的延伸叠加而成的，因而有一定的逼近误差，误差大小取决于理想频率响应曲线形状。理想频率响应特性变化越平缓，则内插值越接近理想值，逼近误差越小。为解决上述问题，可在频率响应的过渡带内插入一个、两个或三个采样点，增加过渡带，减小频带边缘的突变，减小通带和阻带的波动，因而可增大阻带最小衰减。 三、实验内容
1、数字滤波器的技术指标如下：?P?0.4?、?P?0.5dB、?s?0.6?、?S?50dB，采用窗函数法来设计一个FIR数字滤波器。画出各种滤波器的幅频和相频特性曲线。
（1）选择一个合适的窗函数，取 N＝15，观察所设计滤波器的幅频特性，分析是否满足设计要求；
（2）取N＝35，重复上述设计，观察幅频和相频特性的变化，分析长度N变化的影响；
clear all;
N=35; wp=0.4*pi; ws=0.6*pi; wc=1/2*(wp+ws); n=[0:(N-1)]; tao=(N-1)/2; m=n-tao+eps; hd=sin(wc*m)./(pi*m); deltaw=ws-wp; N0=ceil(6.6*pi/deltaw); N=N+mod(N0+1,2); wdham=(hamming(N))'; h=hd.*wdham; [H,w]=freqz(h,[1],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); subplot(221);stem(n,hd); subplot(222);stem(n,wdham); subplot(223);stem(n,h); subplot(224);plot(w,db);
0.60.40.20-0.201020304010.80.60.40.200102030400.60.40.20-0.2010203040500-50-100-15001234
（3）保持 N＝35不变，改变窗函数为blackman窗，观察并记录窗函数对滤波器幅频特性的影响，比较两种窗的特点。
提示：m=n-tao+eps; %eps为一个很小的数，避免下一步中分母为0； hd=sin(wc*m)./(pi*m); %求理想脉冲响应。 clear all;
N=35; wp=0.4*pi; ws=0.6*pi; wc=1/2*(wp+ws); n=[0:(N-1)]; tao=(N-1)/2; m=n-tao+eps;