基于窗函数设计法线性相位型FIR数字低通滤波器设计.doc

实 验 报 告课程名称 数字信号处理 实验项目 数字信号处理课程设计- 7.基于窗函数设计法线性相位型FIR数字低通滤波器设计系 别 信息与通信工程 专业/班级 姓名/学号 实验日期 成 绩 _ 指导教师 _ 数字信号处理课程设计任务书题目7基于窗函数设计法线性相位型FIR数字低通滤波器设计主要内容1、根据ALPF指标，采用三种窗函数，设计线性相位型FIR滤波器；2、生成一个时间序列，对之滤波，并分析三种滤波器的滤波效果。设计要求1、根据ALPF指标，构建物理可实现的线性相位FIR滤波器的冲击响应函数；2、采用三种窗函数，设计线性相位型FIR滤波器；3、生成一个时间序列，并利用已设计的线性相位型FIR滤波器对之滤波；4、对比分析三种窗函数法设计的数字滤波器性能。主要仪器设备1、计算机1台，安装MATLAB软件主要参考文献美维纳.K.恩格尔，约翰.G.普罗科斯著，刘树棠译.数字信号处理使用MATLABM.西安：西安交通大学出版社，2002. 飞思科技产品研发中心编著.MATLAB7辅助信号处理技术与应用M.北京：电子工业出版社，2005.课程设计进度安排（起止时间、工作内容）课程设计共设8个设计题目，每班3至4人为1组，1人1套设备，每组选作不同的题目，4个班1批。完整课程设计共20学时，为期1周，具体进度如下：5学时 学习题目相关知识，掌握实现原理；5学时 用MATLAB语言实现题目要求；5学时 进一步完善功能，现场检查、答辩；5学时 完成并提交课程设计报告。课程设计开始日期2010.12.13课程设计完成日期2010.12.17课程设计实验室名称健翔桥校区计算中心地 点计算中心资料下载地址各班公共邮箱摘要本次实验主要利用窗函数设计法线性相位型FIR数字低通滤波器。根据ALPF指标，采用三种窗函数，设计线性相位型FIR滤波器，并比较它们在相同条件下的滤波效果。最后生成一个时间序列，对之滤波，并分析两种滤波器的滤波效果。关键词：Kaiser窗 Blackman窗 Hamming窗 FIR 基于窗函数设计法线性相位型FIR数字低通滤波器设计目的： 熟练掌握MATLAB工具软件在工程设计中的使用； 熟练掌握窗函数法设计线性相位型FIR滤波器的方法。要求： 根据给定ALPF幅频特性指标，生成与之相对应的模拟滤波器h(t)； 根据h(t)取得与之相对应的数字滤波器冲击响应函数h(n)； 选择适当的窗函数，构建线性相位型FIR的冲击响应函数hd(n)； 设计与之相对应的DLPF，对比分析DLPF幅频特性是否符合要求； 产生一个有干扰频率的时域序列（借助FFT分析说明其有干扰），使之通过所设计的DLPF，对滤波输出结果作出分析，说明输出结果。步骤： 根据给定设计指标（通带截止频率fb1 000 Hz，通带最大衰减p=1dB，阻带截止频率fc1 200 Hz，阻带最小衰减s=100dB，令T=2）生成与之相对应的模拟滤波器h(t)； 根据h(t)取得与之相对应的数字滤波器冲击响应函数h(n)； 选择适当的窗函数，构建线性相位型FIR的冲击响应函数hd(n)； 设计与之相对应的DLPF，给出窗函数及所设计滤波器的幅度特性，并对比理想冲激响应和实际冲击响应，对比分析DLPF幅频特性是否符合要求； 试说明不同窗函数对所设计数字滤波器性能的影响； 产生一个有干扰频率的时域序列（借助FFT分析说明其有干扰），使之通过所设计的DLPF，对滤波输出结果作出分析，说明输出结果。目录一、前言5二、 FIR滤波器的设计原理6三、FIR滤波器的窗函数设计法7（一）Hamming窗7（二）Blackman窗9（三）Kaiser-Basel窗11四、实验分析13五、实验总结13一、前言 数字信号处理（Digital Signal Processing，简称DSP）是一门涉及许多学科而又广泛应用于许多领域的新兴学科。20世纪60年代以来，随着计算机和信息技术的飞速发展，数字信号处理技术应运而生并得到迅速的发展。在过去的二十多年时间里，数字信号处理已经在通信等领域得到极为广泛的应用。数字信号处理是利用计算机或专用处理设备，以数字形式对信号进行采集、变换、滤波、估值、增强、压缩、识别等处理，以得到符合人们需要的信号形式。数字信号处理是围绕着数字信号处理的理论、实现和应用等几个方面发展起来的。数字信号处理在理论上的发展推动了数字信号处理应用的发展。反过来，数字信号处理的应用又促进了数字信号处理理论的提高。而数字信号处理的实现则是理论和应用之间的桥梁。虽然数字信号处理的理论发展迅速，但在20世纪80年代以前，由于实现方法的限制，数字信号处理的理论还得不到广泛的应用。直到20世纪70年代末80年代初世界上第一片单片可编程DSP芯片的诞生，才将理论研究结果广泛应用到低成本的实际系统中，并且推动了新的理论和应用领域的发展。可以毫不夸张地说，DSP芯片的诞生及发展对近20年来通信、计算机、控制等领域的技术发展起到十分重要的作用。二、 FIR滤波器的设计原理FIR滤波器通常采用窗函数方法来设计。窗设计的基本思想是，首先选择一个适当的理想选频滤波器（它总是具有一个非因果，无限持续时间脉冲响应），然后街区（加窗）它的脉冲响应得到线性相位和因果FIR滤波器。我们用Hd（ejw）表示理想的选频滤波器，它在通带上具有单位增益和线性相位，在阻带上具有零响应。一个带宽wcpi的低通滤波器由下式给定：为了从hd(n)得到一个FIR滤波器，必须同时在两边截取hd（n）。而要得到一个因果的线性相位滤波器，它的h(n)长度为N，必须有：这种操作叫做加窗，h(n)可以看做是hd(n)与窗函数w（n）的乘积： h(n)=hd(n)w(n)其中 根据w（n）的不同定义，可以得到不同的窗结构。 在频域中，因果FIR滤波器响应H(ejw)由Hd(ejw)和窗响应W(ejw)的周期卷积得到，即 常用的窗函数有矩形窗、巴特利特（BARTLETT）窗、汉宁（HANNING）窗、海明（HAMMING）窗、布莱克曼（BLACKMAN）窗、凯泽（KAISER）窗等。三、FIR滤波器的窗函数设计法 FIR滤波器的设计方法有许多种，如窗函数设计法、频率采样设计法和最优化设计法等。窗函数设计法的基本原理是用一定宽度窗函数截取无限脉冲响应序列获得有限长的脉冲响应序列，主要设计步骤为： (1)通过傅里叶逆变换获得理想滤波器的单位脉冲响应hd(n)。 (2)由性能指标确定窗函数W(n)和窗口长度N。 (3)求得实际滤波器的单位脉冲响应h(n)， h(n)即为所设计FIR滤波器系数向量b(n)。 （一）Hamming窗hamming窗程序fp=1000;fs=1200;Fs=10000;rs=100;%给定设计指标wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;rp=1;B=ws-wp;%过渡带M=ceil(rs-8)/2.285/Bt);%求Nwc=(wp+ws)/2/pi;x=-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0;%一个有干扰频率的时域序列subplot(3,2,1);plot(x);%画出干扰序列时域图title(干扰序列时域图);y=fft(x,1024);%让x做快速傅里叶变换subplot(3,2,2);plot(abs(y);%画出x频域图像title(干扰序列频域图);hn1=fir1(M,wc,hamming(M+1); %求hamming窗hn1 n=0:N;subplot(3,2,3);stem(n,hn1,.);title(海明窗);Hw,w=freqz(hn1);subplot(3,2,4);plot(w*Fs/20*pi,20*log10(abs(Hw);grid on;%画hn1频域的图形title(海明窗频域图);hn2=ifft(fft(hn1,1024).*fft(x,1024),1024);%让干扰序列通过滤波器subplot(3,2,5);plot(hn2);title(通过滤波器后时域图);grid on%画hn2图形hn3=fft(hn2);%让hn2做fftsubplot(3,2,6);plot(abs(hn3);title(通过滤波器后频域图);grid on;%画出hn2频域的图像db=20*log10(abs(hn3)+eps)/max(abs(hn3); %求FIR滤波器频响的dB值delta_w=2*pi/1000; %将2pi等分1000份Rp1=-min(db(1:1:(wp/delta_w+1) %求Passband RippleAs1=-max(db(ws/delta_w+1:1:501) Rp1 = 4.2440As1 = 66.5416（二）Blackman窗Blackman窗程序fp=1000;fs=1200;Fs=10000;rs=100;%给定设计指标wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;rp=1;B=ws-wp;%过渡带N=ceil(rs-8)/2.285/Bt);%求Nwc=(wp+ws)/2/pi;x=-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0;%一个有干扰频率的时域序列subplot(3,2,1);plot(x);%画出干扰序列时域图title(干扰序列时域图);y=fft(x,1024);%让x做快速傅里叶变换subplot(3,2,2);plot(abs(y);%画出x频域图像title(干扰序列频域图);hn1=fir1(N,wc,blackman(N+1); %求Blackman窗hn1 n=0:N;subplot(3,2,3);stem(n,hn1,.);title(Blackman窗);Hw,w=freqz(hn1);subplot(3,2,4);plot(w*Fs/20*pi,20*log10(abs(Hw);title(Blackman窗频域图);grid on;%画hn1频域的图形hn2=ifft(fft(hn1,1024).*fft(x,1024),1024);%让干扰序列通过滤波器subplot(3,2,5);plot(hn2);title(通过滤波器后时域图);grid on;%画hn2图形hn3=fft(hn2);%让hn2做fftsubplot(3,2,6);plot(abs(hn3);title(通过滤波器后频域图);grid on;%画出hn2频域的图像db=20*log10(abs(hn3)+eps)/max(abs(hn3); %求FIR滤波器频响的dB值delta_w=2*pi/1000; %将2pi等分1000份Rp1=-min(db(1:1:(wp/delta_w+1) %求Passband RippleAs1=-max(db(ws/delta_w+1:1:501) Rp1 = 4.2329Warning: Integer operands are required for colon operator when used as index As1 = 58.8327（三）Kaiser-Basel窗Kaiser-Basel窗程序fp=1000;fs=1200;Fs=10000;rs=100;%给定设计指标wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;rp=1;Bt=ws-wp;%计算过渡带alph=0.112*(rs-8.7);%计算Kaiser窗参数M=ceil(rs-8)/2.285/Bt);%计算Mwc=(wp+ws)/2/pi;x=-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0;%一个有干扰频率的时域序列subplot(3,2,1);plot(x);%画出干扰序列的时域图像title(干扰序列时域图);y=fft(x,1024);%做1024点快速傅里叶变换subplot(3,2,2);plot(abs(y);%画出fft后的频谱图title(干扰序列频域图);hn=fir1(M,wc,kaiser(M+1,alph); %求Kaiser窗hn n=0:M;subplot(3,2,3);stem(n,hn,.);%画出hn图形title(Kaiser窗);Hw,w=freqz(hn);subplot(3,2,4);plot(w*Fs/20*pi,20*log10(abs(Hw);title(Kaiser窗频域图);grid on;%hn频域图想hn2=ifft(fft(hn,1024).*fft(x,1024),1024);%让干扰序列通过低通滤波器subplot(3,2,5);plot(hn2);title(通过滤波器后时域图);grid on;%画滤波之后时域hn2hn3=fft(hn2);%hn2做f