logo资料库

滤波器设计实例及Matlab代码.pdf

发布时间：2022-06-13 发布人：admin 分类：说明书 资料大小：0.99M 资料格式：pdf 举报 版权申诉
第1页 / 共29页
第2页 / 共29页
第3页 / 共29页
第4页 / 共29页
第5页 / 共29页
第6页 / 共29页
第7页 / 共29页
第8页 / 共29页
资料共29页，剩余部分请下载后查看
    数字滤波器 1.从功能分：低通、带通、高通、带阻（每一种又可分：Digital 和 Analog filter） 2.从实现方法分：FIR、IIR 3.从设计方法分：Chebyshev(切比雪夫)、Butterworth(巴特沃斯) 、椭圆滤波 4.从处理信号分：经典滤波器、现代滤波器 （模拟滤波器用于处理模拟信号，即连续时间信号； 数字滤波器用于处理离散时间信号。） 经典数字滤波器按照单位取样响应 h(n)的时域特性可分为无限冲激响应(IIR， Infinite Impulse Response) 系 统 和 有 限 冲 激 响 应 (FIR ， Finite Impulse Response)系统。如果信号和噪声的频谱相互重叠，那么经典滤波器无能为力。 现代滤波器把信号和噪声都视为随机信号，利用他们的统计特征(如自相关函数、 功率谱等)导出一套最佳估值算法，然后用硬件或软件予以实现。如维纳滤波器。 2.1 FIR 滤波器 FIR 滤波器的优点在于有线性的相位特征，差分方程为： 有三种设计方法：窗函数法、频率采样法、切比雪夫等波纹逼近法 介绍窗函数法滤波器的 Matlab 程序： 如果所希望的滤波器的理想频率响应函数为 ，则其对应的单位脉冲响应为： 窗函数设计法的基本原理是用有限长单位脉冲响应序列 h(n)逼近 hd(n)。由于 hd(n)往往是无 限长序列，且是非因果的，做法是：只取其中的某些项，即只截取 hd(n)中的一部分，比如 n=0，1，2，3，...，N-1，N 为整数。这种处理相当于将 hd(n)与函数 w(n)相乘，w(n)为： w(n)相当与一个矩形，称之为矩形窗。即我们可用矩形窗函数 w(n)将无限脉冲响应 hd(n)截 取一段 h(n)来近似 hd(n)，这种截取在数学上表示为乘积形式（变换为频率域为褶积形式）： 这里加窗函数不是可有可无的，而是变为物理可实现所必须的。h(n)就作为实际设计的 FIR 数字滤波器的单位脉冲响应序列，其频率响应函数 为： （离散时间傅立叶取有限项）
    式中，N 为所选窗函数 w(n)的长度。 2.1.1 窗函数 这些窗函数的基本参数如表 1 所示
    应用类型： 幂窗--如矩形、三角形、梯形，应用于幂函数 三角函数窗--如汉宁窗、海明窗，应用于三角函数 指数窗--如高斯窗，应用于指数函数 2.1.2 fir1 函数 设计标准响应 FIR 滤波器可使用 fir1 函数，以经典方法实现加窗线性相位 FIR 滤波器设计， 可设计出：低通、带通、高通、带阻滤波器。形式为： b-----滤波器系数；长度为 n+1 n-----滤波器阶数； Wc-----截止频率；（归一化） ftype-----’high’ ‘stop’ ‘low’ ‘bandpass’; Window-----窗函数，长度为 FIR 滤波器系数个数，即 n+1 说明：Window 先进行 Z 变换，再根据 Wc 截至频率，确定频率域的定义域，然后和理想滤 波器相乘，就得到 b(z)，取其系数就得到 b。 2.1.3 freqz 函数 freqz 函数计算线性系统的频率响应,包括幅频响应和相频响应,基本输入为线性系统的 AMMA 模型系数向量,一个典型的 AMMA 模型为 % jw -jw -jmw % jw B(e) b(1) + b(2)e + .... + b(m+1)e % H(e) = ---- = ------------------------------------ % jw -jw -jnw % A(e) a(1) + a(2)e + .... + a(n+1)e 其中,向量 A 为 IIR 系统分母向量,对 FIR 系统为1;向量 B 为分子向量函数使用离散傅立叶变 换(fft)计算系统响应,计算范围为0~pi,即单位圆的上半部分,因为默认情况下对于实数系 数系统响应关于 DC 对称.而 fft 是在2Pi 单位圆。 ，w----等分频率的点；h----对应的频率响应；n----默 调用格式： 认值是512。 2.1.4 各个窗函数响应曲线 程序如下： 取窗函数宽度 N=150 程序代码如下： clear all; wls=0.2*pi;wlp=0.35*pi; whp=0.65*pi;whs=0.8*pi;%参数设置 delta_w=min((wlp-wls),(whs-whp));%求两个过渡带中的小者 wc1=(wls+wlp)/2;wc2=(whp+whs)/2;%截至频率取通带阻带边界频率的平均值
    N=150; %矩形窗 w1=boxcar(N); hn1=fir1(N-1,[wc1,wc2]/pi,boxcar(N));%检验设计出的滤波器单位脉冲响应，即：离散傅立叶变换 [h1,w]=freqz(hn1,1);%验证设计出的滤波器的频率响应，默认点数为512，在0-pi范围 %Hamming窗 w2=hamming(N); hn2=fir1(N-1,[wc1,wc2]/pi,hamming(N)); [h2,w]=freqz(hn2,1); %Blackman窗 w3=blackman(N); hn3=fir1(N-1,[wc1,wc2]/pi,blackman(N)); [h3,w]=freqz(hn3,1); %Kaiser窗(beta=7.865) w4=kaiser(N); hn4=fir1(N-1,[wc1,wc2]/pi,kaiser(N,7.865)); [h4,w]=freqz(hn4,1); figure(3) n=0:N-1; plot(n,w1,':b',n,w2,'-.k',n,w3,'--r',n,w4,'-k'); axis([-10 160 0 1.1]); xlabel('n');ylabel('w(n)'); title('各种窗函数(N=150)'); legend('矩形窗带通滤波器','Hamming窗带通滤波器','Blackman窗带通滤波器','Kaiser窗带通滤波 器'); figure(4) plot(w/pi,20*log10(abs(h1)),':b',w/pi,20*log10(abs(h2)),'-.k',w/pi,20*log10(abs(h3)),'--r',w/pi,20*log1 0(abs(h4)),'-k'); axis([0 1 -150 5]); xlabel('归一化角频率');ylabel('幅度（单位：分贝）');grid on; title('各种窗函数幅频特征比较(N=150)') legend('矩形窗带通滤波器','Hamming窗带通滤波器','Blackman窗带通滤波器','Kaiser窗带通滤波 器'); 图形如下：
    2.1.5 设计低通、带通、高通、带阻滤波器
    ○1 采用 Kaiser 窗设计一个低通 FIR 滤波器 要求：采样频率为 8kHz 通带：0Hz—1kHz，带内波动小于 5% 阻带：1.5kHz，带内最小衰减：Rs=40dB 思路分析：通带截至频率为 wp=0.25pi，阻带截止频率为 ws=0.375pi，低通滤波器的时域表 ， 其 中 α 关 于 (N-1)/2 对 称 ， 这 样 滤 波 器 就 得 到 为 达 式 为 ) ( ) α − ) π α w n ( ( n − sin n h ( ) d = h ( ) n = d n h ( ) ( ) nω ，最后利用函数 freqz 得到加窗后的滤波器的幅频响应和相频响应。 低通 FIR 滤波器 程序如下： %采样频率为8kHz %通带：0Hz-1kHz,带内波动小于5%；wp=0.25pi %阻带：1.5Hz，带内最小衰减：Rs=40dB；ws=0.375pi % clc; clear all; Rs=40;%带内最小衰减 wp=0.25*pi; ws=0.375*pi; delta_w=ws-wp; N=ceil(10*pi/delta_w); beta=0.5842*(Rs-21)^0.4+0.07886*(Rs-21); % wc=(wp+ws)/2; alpha=(N-1)/2; n=[0:N-1]; m1=n-alpha+eps;%加一个小数以避免零作除数 hd=sin(wc*m1)./(pi*m1);%滤波器在时域系统的冲击响应 B=kaiser(N,beta);%凯泽窗 h=hd.*(B)';%加窗后 %以上相当于 h=fir1(N-1,wc/pi,’low’,kaiser(N,beta)); % [H,m]=freqz(h,[1],1024,'whole');%获取频率响应,’whole’表示 m在0-2pi范围内被分成1024段 mag=abs(H);%幅值 db=20*log10((mag+eps)/max(mag));%分贝数 pha=angle(H);%相位 %绘图 w=m/pi; figure(1); subplot(2,2,1);stem(n,hd); xlabel('n');ylabel('hd');
    title('滤波器时域'); subplot(2,2,2); plot(w,mag); %加窗后时间域的响应 xlabel('w');ylabel('h'); title('加窗后幅度响应');％画错了，应该w的范围是：0-1 subplot(2,2,3); plot(w,db); %实际低通滤波器频率域响应 xlabel('归一化w');ylabel('db'); title('分贝数'); axis([0 1 -100 0]); subplot(2,2,4); plot(w,pha);%实际低通滤波器相位响应 xlabel('归一化w');ylabel('相位'); title('相频响应'); axis([0 1 -4 4]); ○2 采用 Hamming 窗设计一个高通 FIR 滤波器 要求：通带截止频率 wp=0.6pi 阻带截止频率 ws=0.4pi 0.25 通带最大衰减 α = dB p 阻带最小衰减 分析： w p α p = = / f F 2 π s p 20log (1 10 w 2 π = δ α = + − p s s ) f F / s s 20log 10 δ s α = s 50 dB
    通带偏差 pδ = 0.0292 阻带偏差 sδ = 0.0032 通带边沿频率 fp=1000Hz 阻带边沿频率 fs=600Hz 选择窗函数 W(n),计算窗函数长度 N，由已知条件知：阻带最小衰减 α = s 40 dB 汉宁窗和哈明窗都满足要求，这里选择汉宁窗，过渡带宽 B w t ≤ p w −= s 0.2 π 汉宁窗的精确带宽 tB Nπ= 6.2 / ， 故要求 tB =≤ 6.2 / Nπ 0.2 π 。 解得：N>= 31 对于高通滤波器和带阻滤波器，N 必须取奇数。 程序代码： % clc; clear all; wp=0.6*pi; ws=0.4*pi; tr_width=wp-ws; N=ceil(6.2*pi/tr_width); if (mod(N,2)==0) N=N+1; end %高通滤波器取奇数 n=0:1:N-1; wc=(ws+wp)/2; % alpha=(N-1)/2; n=[0:N-1]; m1=n-alpha+eps;%加一个小数以避免零作除数 hd_0=sin(pi*m1)./(pi*m1); hd_h=sin(wc*m1)./(pi*m1); hd=hd_0-hd_h;%高通滤波器的时域的冲击响应 w_han=(hanning(N))'; h=hd.*w_han; %equal to h=fir1(N-1,wc/pi,'high',hanning(N)); % [H,w]=freqz(h,1,1000,'whole');%注：'whole'表示此时w在2pi范围内被分成1000份,离散傅立叶变换 H=(H(1:1:501))'; w=(w(1:1:501))'; mag=abs(H); db=20*log10((mag+eps)/max(mag)); pha=angle(H); % delta_w=2*pi/1000; Ap=-(min(db(wp/delta_w+1:1:501)));
    分享到：
    赞 收藏

    相关推荐

    资料库

    课程资源

    共收录17145份资料，累计13个分类，关注成员有19位，主要包括：PHP,网络管理,网页制作,Java,.Net,数据库,3G/移动开发,C/C++,游戏开发,嵌入式,讲义,软件测试,专业指导

    热门标签

    PHP 网络管理 网页制作 Java .Net 数据库 3G/移动开发 C/C++ 游戏开发 嵌入式 讲义 软件测试 专业指导

    最新资料