基于Blackman窗的FIR数字带通滤波器设计

课程设计课程设计名称:基于Blackman窗的FIR数字带通滤波器设计专业班级：学生姓名：学号：指导教师：课程设计时间：2010.6.21-2010.6.25数字信号处理专业课程设计任务书学生姓名专业班级学号题目基于Blackman窗的FIR数字带通滤波器设计课题性质工程技术研究课题来源自拟课题指导教师同组姓名主要内容设计一个具有如下性能指标的FIR数字带通滤波器：1p=0.45，2p=0.65，1s=0.30，2s=0.80；窗函数为Blackman窗；分析最后设计结果性能。任务要求1.温习窗函数法设计数字FIR滤波器的原理和设计方法；2.求出最后所设计出的滤波器冲激响应系数；3.绘制所设计滤波器的增益响应曲线。参考文献[1]胡广书.数字信号处理—理论、算法与实现[M].北京:清华大学出版社,1997.[2]R.Lyons.UnderstandingDigitalSignalProcessing[M].2nded.PrenticeHallPTR.,2004.[3]A.V.奥本海姆,R.W.谢弗andJ.R.巴克.离散时间信号处理[M].第二版.西安交通大学出版社,2001.[4]S.K.Mitra.DigitalSignalProcessing:AComputer-BasedApproach[M].3rded.McGraw-Hill,2005.审查意见指导教师签字：李相国教研室主任签字：年月日说明：本表由指导教师填写，由教研室主任审核后下达给选题学生，装订在设计（论文）首页1需求分析数字滤波器可以理解为是一个计算程序或算法，将代表输入信号的数字时间序列转化为代表输出信号的数字时间序列，并在转化过程中，使信号按预定的形式变化。数字滤波器有多种分类，根据数字滤波器冲激响应的时域特征，可将数字滤波器分为两种，即无限长冲激响应（iir）滤波器和有限长冲激响应（fir）滤波器。iir数字滤波器具有无限宽的冲激响应，与模拟滤波器相匹配。所以iir滤波器的设计可以采取在模拟滤波器设计的基础上进一步变换的方法。fir数字滤波器的单位脉冲响应是有限长序列。它的设计问题实质上是确定能满足所要求的转移序列或脉冲响应的常数问题，设计方法主要有窗函数法、频率采样法和等波纹最佳逼近法等。因此设计FIR滤波器的方法之一可以从时域出发，截取有限长的一段冲击响应作为H(z)的系数，冲击响应长度N就是系统函数H(z)的阶数。只要N足够长，截取的方法合理，总能满足频域的要求。一般这种时域设计、频域检验的方法要反复几个回合才能成功。要设计一个线性相位的FIR数字滤波器，首先要求理想频率响应)(jwdeH。)(jwdeH是w的周期函数，周期为2，可以展开成傅氏级数：)(jwdeH＝njwndenh)(其中)(nhd是与理想频响对应的理想单位抽样响应序列。但不能用来作为设计FIRDF用的h(n)，因为)(nhd一般都是无限长、非因果的，物理上无法实现。为了设计出频响类似于理想频响的滤波器，可以考虑用h(n)来近似)(nhd。窗函数的基本思想：先选取一个理想滤波器（它的单位抽样响应是非因果、无限长的），再截取（或加窗）它的单位抽样响应得到线性相位因果FIR滤波器。这种方法的重点是选择一个合适的窗函数和理想滤波器。设x(n)是一个长序列，w(n)是长度为N的窗函数，用w(n)截断x(n)，得到N点序列xn(n)，即xn(n)=x(n)w(n)在频域上则有ππjjjdeπ21eWeXXN由此可见，窗函数w(n)不仅仅会影响原信号x(n)在时域上的波形，而且也会影响到频域内的形状。MATLAB信号工具箱主要提供了以下几种窗函数，如表1-1所示：表1-1MATLAB窗函数窗窗函数矩形窗Boxcar巴特利特窗Barlett三角窗Triang布莱克曼窗Blackman海明窗Hamming汉宁窗Hanning凯塞窗Kaiser切比雪夫窗Chebwin加矩形窗后的频谱和理想频谱可得到以下结论：加窗使过渡带变宽，过渡带的带宽取决于窗谱的主瓣宽度。矩形窗情况下的过渡带宽是N/4。N越大，过渡带越窄、越陡；过渡带两旁产生肩峰，肩峰的两侧形成起伏振荡。肩峰幅度取决于窗谱主瓣和旁瓣面积之比。矩形窗情况下是8.95％，与N无关。工程上习惯用相对衰耗来描述滤波器，相对衰耗定义为])0(/)(lg[20])(/)(lg[20)(0HwHeHeHwAjjw这样两个肩峰点的相对衰耗分别是0.74dB和-21dB。其中（-0.0895）对应的点的值定义为阻带最小衰耗。以上的分析可见，滤波器的各种重要指标都是由窗函数决定，因此改进滤波器的关键在于改进窗函数。窗函数谱的两个最重要的指标是：主瓣宽度和旁瓣峰值衰耗。旁瓣峰值衰耗定义为：旁瓣峰值衰耗＝20lg(第一旁瓣峰值／主瓣峰值)为了改善滤波器的性能，需使窗函数谱满足：主瓣尽可能窄，以使设计出来的滤波器有较陡的过渡带；第一副瓣面积相对主瓣面积尽可能小，即能量尽可能集中在主瓣，外泄少，使设计出来的滤波器的肩峰和余振小。但上面两个条件是相互矛盾的，实际应用中，折衷处理，兼顾各项指标。本次设计要求以窗函数为Blackman窗为基础，设计一个具有如下性能指标的FIR数字带通滤波器：1p=0.45，2p=0.65，1s=0.30，2s=0.80。要求：1.温习窗函数法设计数字FIR滤波器的原理和设计方法；2．求出最后所设计出的滤波器冲激响应系数；3．绘制所设计滤波器的增益响应曲线。2概要设计如果所希望的滤波器的理想的频率响应函数为jdeH，则其对应的单位脉冲响应为deeHnhjjdd21窗函数设计法的基本原理是用有限长单位脉冲响应序列nh逼近nhd。由于nhd往往是无限长序列，而且是非因果的，所以用窗函数n将nhd截断，并进行加权处理，得到：nnhnhdnh就作为实际设计的FIR数字滤波器的单位脉冲响应序列，其频率响应函数jeH为njNnjenheH10式中，N为所选窗函数n的长度。我们知道，用窗函数法设计的滤波器性能取决于窗函数n的类型及窗口长度N的取值。设计过程中，要根据对阻带最小衰减和过渡带宽度的要求选择合适的窗函数类型和窗口长度N。主程序框图如图4.1所示。其中幅度特性要求用dB表示。布莱克曼窗函数布莱克曼窗函数的时域形式可以表示为11π4cos08.011π2cos5.042.0)(NkNkkwNk,,2,1它的频域特性为：)(W42.0RW25.01212NWNWRR1π41π404.0NWNWRR其中，)(RW为矩形窗函数的幅度频率特性函数。布莱克曼窗函数的最大旁瓣值比主瓣值低57dB，但是主瓣宽度是矩形窗函数的主瓣宽度的3倍，为12π/N。Blackman函数：生成海明窗调用方式开始读入窗口长度N计算hd(n)调用窗函数子程序求w(n)调用子程序（函数）计算H(k)=DFT[h(n)]调用绘图子程序（函数）绘制H(k)幅度相位曲线结束图1主程序框图计算h(n)=hd(n)w(n)(1)w=blackman(n)：输入参数n是窗函数的长度；输出参数w是由窗函数的值组成的n阶向量。(2)w=blackman(n,sflag)：参数sflag用来控制窗函数首尾的两个元素值；其取值为symmetric或periodic；默认值为symmetric。3运行环境软件环境：1.Windowsxp2.MATLAB6.53．其他软件根据需要选用硬件环境：计算机：LenovoY450处理器型号：Intel酷睿2双核T6600处理器参数：主频2.00GHZ主板芯片组：IntelGM45内存：2GBDDRIII硬盘：西数WDCWD3200BEVT-22ZCTOSATA显存：512MB4开发工具和编程语言开发工具:MATLAB桌面和命令窗口、历史命令窗口、编辑器和调试器、路径搜索和用于用户浏览帮助、工作空间、文件的浏览器。MATLAB一个高级的矩阵/阵列语言，它包含控制语句、函数、数据结构、输入和输出和面向对象编程特点。用户可以在命令窗口中将输入语句与执行命令同步，也可以先编写好一个较大的复杂的应用程序（M文件）后再一起运行。编程语言:汇编语言汇编语言(AssemblyLanguage)是面向机器的程序设计语言。在汇编语合中，用助记符(Memoni)代替操作码，用地址符号(Symbol)或标号(Label)代替地址码。这样用符号代替机器语言的二进制码，就把机器语言变成了汇编语言。于是汇编语言亦称为符号语言。使用汇编语言编写的程序，机器不能直接识别，要由一种程序将汇编语言翻译成机器语言，这种起翻译作用的程序叫汇编程序，汇编程序是系统软件中语言处理系统软件。汇编程序把汇编语言翻译成机器语言的过程称为汇编。汇编语言的特点:1.面向机器的低级语言，通常是为特定的计算机或系列计算机专门设计的。2.保持了机器语言的优点，具有直接和简捷的特点。3.可有效地访问、控制计算机的各种硬件设备，如磁盘、存储器、CPU、I/O端口等。4.目标代码简短，占用内存少，执行速度快，是高效的程序设计语言。5.经常与高级语言配合使用，应用十分广泛。5详细设计算法实现的源程序：（1）理想低通滤波器函数function[db,mag,pha,grd,w]=freqz_m(b,a);N=1000;[H,w]=freqz(b,a,N,'whole');H=(H(1:1:501))';w=(w(1:1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);grd=grpdelay(b,a,w);（2）定义滤波器系数向量函数functionhd=ideallp(wc,M);al=(M-1)/2;n=[0:(M-1)];m=n-al+eps;hd=sin(wc*m)./(pi*m);（3）主程序：%数字滤波器指标ws1=0.30*pi;wp1=0.45*pi;ws2=0.80*pi;wp2=0.65*pi;As=60;tr_width=min((wp1-ws1),(wp2-ws2));M=ceil(11*pi/tr_width)+1;n=[0:1:M-1];wc1=(ws1+wp1)/2;wc2=(ws2+wp2)/2;hd=ideal_lp(wc2,M)-ideal_lp(wc1,M);%生成blackman窗w_bla=(blackman(M))';h=hd.*w_bla%频域图像的绘制freqz(h,[1])figure(2);subplot(2,2,1),stem(n,hd);title('idaelimpulseresponse')axis([0M-1-0.40.5]);xlabel('n');ylabel('hd(n)')xa=0.*n;holdonplot(n,xa,'k');holdoffsubplot(2,2,2),stem(n,w_bla);title('blackmanwindow')axis([0M-101.1]);xlabel('n');ylabel('w(n)')subplot(2,2,3),stem(n,h);title('actualimpulseresponse')axis([0M-1-0.40.5]);xlabel('n');ylabel('h(n)')holdonplot(n,xa,'k');holdoffn=0:600;k=n/2000;f1=2*pi*450;f2=2*pi*600;x=sin(f1*k)+sin(f2*k);y=filter(h,[1],x);x1=sin(f1*k);x2=sin(f2*k);figure(2)subplot