matlab设计滤波器程序
- 格式:pdf
- 大小:62.09 KB
- 文档页数:2
基于MATLAB设计FIR滤波器FIR(Finite Impulse Response)滤波器是一种数字滤波器,它具有有限的冲激响应长度。
基于MATLAB设计FIR滤波器可以使用signal工具箱中的fir1函数。
fir1函数的语法如下:b = fir1(N, Wn, window)其中,N是滤波器的阶数,Wn是截止频率,window是窗函数。
要设计一个FIR低通滤波器,可以按照以下步骤进行:步骤1:确定滤波器的阶数。
阶数决定了滤波器的截止频率的陡峭程度。
一般情况下,阶数越高,滤波器的陡峭度越高,但计算复杂度也会增加。
步骤2:确定滤波器的截止频率。
截止频率是指在滤波器中将信号的频率限制在一定范围内的频率。
根据应用的需求,可以选择适当的截止频率。
步骤3:选择窗函数。
窗函数是为了在时域上窗口函数中心增加频率衰减因子而使用的函数。
常用的窗函数有Hamming、Hanning等。
窗函数可以用来控制滤波器的幅度响应特性,使得它更平滑。
步骤4:使用fir1函数设计滤波器。
根据以上步骤确定滤波器的阶数、截止频率和窗函数,可以使用fir1函数设计FIR滤波器。
具体代码如下:N=50;%设定阶数Wn=0.5;%设定截止频率window = hanning(N + 1); % 使用Hanning窗函数步骤5:使用filter函数对信号进行滤波。
设计好FIR滤波器后,可以使用filter函数对信号进行滤波。
具体代码如下:filtered_signal = filter(b, 1, input_signal);其中,input_signal是输入信号,filtered_signal是滤波后的信号。
以上,便是基于MATLAB设计FIR滤波器的简要步骤和代码示例。
根据具体需求和信号特性,可以进行相应的调整和优化。
用MATLAB信号处理工具箱进行FIR滤波器设计的三种方法摘要介绍了利用MATLAB信号处理工具箱进行FIR滤波器设计的三种方法:程序设计法、FDATool设计法和SPTool设计法,给出了详细的设计步骤,并将设计的滤波器应用到一个混和正弦波信号,以验证滤波器的性能。
关键词 MATLAB,数字滤波器,有限冲激响应,窗函数,仿真1 前言数字滤波器是一种用来过滤时间离散信号的数字系统,通过对抽样数据进行数学处理来达到频域滤波的目的。
根据其单位冲激响应函数的时域特性可分为两类:无限冲激响应(IIR)滤波器和有限冲激响应(FIR)滤波器。
与IIR滤波器相比,FIR的实现是非递归的,总是稳定的;更重要的是,FIR滤波器在满足幅频响应要求的同时,可以获得严格的线性相位特性。
因此,它在高保真的信号处理,如数字音频、图像处理、数据传输、生物医学等领域得到广泛应用。
2 FIR滤波器的窗函数设计法 FIR滤波器的设计方法有许多种,如窗函数设计法、频率采样设计法和最优化设计法等。
窗函数设计法的基本原理是用一定宽度窗函数截取无限脉冲响应序列获得有限长的脉冲响应序列,主要设计步骤为:(1) 通过傅里叶逆变换获得理想滤波器的单位脉冲响应hd(n)。
(2) 由性能指标确定窗函数W(n)和窗口长度N。
(3) 求得实际滤波器的单位脉冲响应h(n), h(n)即为所设计FIR滤波器系数向量b(n)。
(4) 检验滤波器性能。
本文将针对一个含有5Hz、15Hz和30Hz的混和正弦波信号,设计一个FIR带通滤波器,给出利用MATLAB实现的三种方法:程序设计法、 FDATool设计法和SPTool设计法。
参数要求:采样频率fs=100Hz,通带下限截止频率fc1=10 Hz,通带上限截止频率fc2=20 Hz,过渡带宽6 Hz,通阻带波动0.01,采用凯塞窗设计。
2 程序设计法MATLAB信号处理工具箱提供了各种窗函数、滤波器设计函数和滤波器实现函数。
滤波器设计MATLAB滤波器的设计在信号处理中具有重要的作用,可以用于去除噪声、增强信号等。
MATLAB是一种强大的工具,可以用于滤波器设计和分析。
本文将介绍如何使用MATLAB进行滤波器设计,并通过示例展示具体的过程。
在MATLAB中,可以使用信号处理工具箱提供的函数来设计滤波器。
常用的函数有:- `fir1`:设计FIR滤波器。
- `butter`:设计巴特沃斯滤波器。
- `cheby1`:设计切比雪夫I型滤波器。
- `cheby2`:设计切比雪夫II型滤波器。
- `ellip`:设计椭圆滤波器。
这些函数的输入参数包括滤波器类型、阶数、截止频率等。
根据具体的需求选择不同的函数来设计滤波器。
下面以设计一个低通滤波器为例,演示如何使用MATLAB进行滤波器设计。
首先,创建一个信号作为输入。
可以使用`sin`函数生成一个正弦信号作为示例。
代码如下:```matlabfs = 1000; % 采样率t = 0:1/fs:1; % 时间向量f=50;%信号频率x = sin(2*pi*f*t); % 输入信号```接下来,使用`fir1`函数设计一个低通滤波器。
该函数的输入参数`n`表示滤波器的阶数,`Wn`表示归一化的截止频率。
代码如下:```matlabn=50;%滤波器阶数Wn=0.2;%截止频率b = fir1(n, Wn);```然后,使用`filter`函数对输入信号进行滤波。
该函数的输入参数是滤波器的系数和输入信号。
代码如下:```matlaby = filter(b, 1, x);```最后,绘制原始信号和滤波后的信号的时域和频域波形。
代码如下:```matlab%时域波形subplot(2, 1, 1)plot(t, x)hold onplot(t, y)legend('原始信号', '滤波后信号') xlabel('时间 (s)')ylabel('幅值')title('时域波形')%频域波形subplot(2, 1, 2)f = linspace(-fs/2, fs/2, length(x)); X = abs(fftshift(fft(x)));Y = abs(fftshift(fft(y)));plot(f, X)hold onplot(f, Y, 'r')legend('原始信号', '滤波后信号') xlabel('频率 (Hz)')ylabel('幅值')title('频域波形')```运行以上代码,可以得到原始信号和滤波后信号的时域和频域波形图。
一、介绍频率采样法设计滤波器的背景和意义1.1 频率采样法设计滤波器的概念及其在数字信号处理中的作用 1.2 频率采样法设计滤波器与其他设计方法的比较1.3 频率采样法设计滤波器的优势和适用范围二、频率采样法设计滤波器的原理和方法2.1 频率采样法设计滤波器的基本原理2.2 频率采样法设计滤波器的设计步骤2.3 频率采样法设计滤波器的常用工具和软件三、matlab频率采样法设计滤波器的实现步骤3.1 设定滤波器的规格和要求3.2 使用matlab进行频域设计3.3 使用matlab进行时域设计3.4 验证设计的滤波器性能四、matlab频率采样法设计滤波器的案例分析4.1 案例一:低通滤波器设计4.1.1 滤波器规格要求4.1.2 频率采样法设计滤波器的实现步骤4.1.3 设计参数及性能分析4.2 案例二:带通滤波器设计4.2.1 滤波器规格要求4.2.2 频率采样法设计滤波器的实现步骤4.2.3 设计参数及性能分析五、matlab频率采样法设计滤波器的应用前景和挑战5.1 应用前景分析5.2 技术发展趋势5.3 面临的挑战和解决方案六、总结与展望6.1 频率采样法设计滤波器的优势和不足6.2 matlab工具在频率采样法设计滤波器中的应用6.3 未来发展方向和趋势在数字信号处理中,滤波器设计是一项重要的工作。
频率采样法设计滤波器是其中一种常用的设计方法,在matlab软件中进行频率采样法设计滤波器具有高效、便捷的特点。
本文将介绍频率采样法设计滤波器的原理、方法以及在matlab中的实现步骤,通过案例分析和应用前景展望来全面解析这一设计方法的优势和发展趋势。
在数字信号处理领域,滤波器设计是至关重要的一环。
而频率采样法设计滤波器作为一种常用的设计方法,在matlab软件中具有高效、便捷的特点。
接下来,我们将深入探讨频率采样法设计滤波器的原理、方法以及在matlab中的实现步骤,并通过案例分析和应用前景展望来全面解析这一设计方法的优势和发展趋势。
设计低通数字滤波器,要求在通带内频率低于0.2pirad时,允许幅度误差在1dB以内,在频率0.3pi rad~pi rad之间的阻带衰减大于15dB,用脉冲响应不变法设计数字滤波器,T=1: 切比雪夫I型模拟滤波器的设计子程序:function [b,a]=afd_chb1(Omegap,Omegar,Ar)if Omegap<=0error('通带边缘必须大于0')endif(Dt<=0)|(Ar<0)error('通带波动或阻带衰减必须大于0');endep=sqrt(10^(Dt/10)-1);A=10^(Ar/20);OmegaC=Omegap;OmegaR=Omegar/Omegap;g=sqrt(A*A-1)/ep;N=ceil(log10(g+sqrt(g*g-1))/log10(OmegaR+sqrt(OmegaR*OmegaR-1)));fprintf('\n***切比雪夫I型模拟低通滤波器阶数=%2.0f\n',N);[b,a]=u_chblap(N,Dt,OmegaC);设计非归一化切比雪夫I型模拟低通滤波器原型程序:function [b,a]=u_chblap(N,Dt,OmegaC)[z,p,k]=cheb1ap(N,Dt);a=real(poly(p));aNn=a(N+1);p=p*OmegaC;a=real(poly(p));aNu=a(N+1);k=k*aNu/aNn;b0=k;B=real(poly(z));b=k*B;直接形式转换成级联形式子程序:function [C,B,A]=sdir2cas(b,a)Na=length(a)-1;Nb=length(b)-1;b0=b(1);b=b/b0;a0=a(1);a=a/a0;C=b0/a0;p=cplxpair(roots(a));K=floor(Na/2);if K*2==NaA=zeros(K,3);for n=1:2:NaArow=p(n:1:n+1,:);Arow=poly(Arow);A((fix(n+1)/2),:)=real(Arow);elseif Na==1A=[0 real(poly(p))];elseA=zeros(K+1,3);for n=1:2:2*KArow=p(n:1:n+1,:);Arow=poly(Arow);A((fix(n+1)/2),:)=real(Arow);endA(K+1,:)=[0 real(poly(p(Na)))];endz=cplxpair(roots(b));K=floor(Nb/2);if Nb==0B=[0 0 poly(z)];elseif K*2==NbB=zeros(K,3);for n=1:2:NbBrow=z(n:1:n+1,:);Brow=poly(Brow);B((fix(n+1)/2),:)=real(Brow);endelseif Nb==1B=[0 real(poly(z))];elseB=zeros(K+1,3);for n=1:2:2*KBrow=z(n:1:n+1,:);Brow=poly(Brow);B((fix(n+1)/2),:)=real(Brow);endB(K+1,:)=[0 real(poly(z(Nb)))];End计算系统函数的幅度响应和相位响应子程序:function [db,mag,pha,w]=freqs_m(b,a,wmax)w1=0:500;w=w1*wmax/500;h=freqs(b,a,w);mag=abs(h);db=20*log10((mag+eps)/max(mag));pha=angle(h);脉冲响应不变法程序:function [b,a]=imp_invr(c,d,T)[R,p,k]=residue(c,d);p=exp(p*T);[b,a]=residuez(R,p,k);b=real(b).*T;数字滤波器响应子程序:function [db,mag,pha,grd,w]=freqz_m(b,a);[H,w]=freqz(b,a,1000,'whole');H=(H(1:501))';w=(w(1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);grd=grpdelay(b,a,w);直接转换成并联型子程序:function [C,B,A]=dir2par(b,a)M=length(b);N=length(a);[r1,p1,C]=residuez(b,a);p=cplxpair(p1,10000000*eps);x=cplxcomp(p1,p);r=r1(x);K=floor(N/2);B=zeros(K,2);A=zeros(K,3);if K*2==Nfor i=1:2:N-2br=r(i:1:i+1,:);ar=p(i:1:i+1,:);[br,ar]=residuez(br,ar,[]);B((fix(i+1)/2),:)real(br');A((fix(i+1)/2),:)real(ar');end[br,ar]=residuez(r(N-1),p(N-1),[]);B(K,:)=[real(br') 0];A(K,:)=[real(ar') 0];elsefor i=1:2:N-1br=r(i:1:i+1,:);ar=p(i:1:i+1,:);[br,ar]=residuez(br,ar,[]);B((fix(i+1)/2),:)real(br);A((fix(i+1)/2),:)real(ar);endEnd比较两个含同样标量元素但(可能)有不同下标的复数对及其相位留数向量子程序:function I=cplxcomp(p1,p2)I=[];for i=1:length(p2)for j=1:length(p1)if(abs(p1(j)-p2(i))<0.0001)I=[I,j];endendendI=I';双线性变换巴特沃斯低通滤波器设计:巴特沃思模拟滤波器的设计子程序:function [b,a]=afd_butt(wp,ws,Rp,rs)if wp<=0error('通带边缘必须大于0');endif ws<=wperror('阻带边缘必须大于通带边缘');endif(Rp<=0)|(Rs<0)error('通带波动或阻带衰减必须大于0');endN=ceil((log10((10^(Rp/10)-1)/(10^(Rs/10)-1)))/(2*log10(wp/ws))); fprintf('\n***Butterworth Filter Order=%2.0f\n',N);OmegaC=wp/((10^(Rp/10)-1)^(1/(2*N)));[b,a]=u_buttap(N,OmegaC)设计非归一化巴特沃思模拟低通滤波器原型子程序:function [b,a]=u_buttap(N,OmegaC)[z,p,k]=buttap(N);p=p*OmegaC;k=k*OmegaC^N;B=real(poly(z));b0=k;b=k*B;a=real(poly(p));直接型到级联型形式的转换:function [b0,B,A]=dir2cas(b,a)b0=b(1);b=b/b0;a0=a(1);a=a/a0;b0=b0/a0;M=length(b);N=length(a);if N>Mb=[b,zeros(1,N-M)];a=[a,zeros(1,M-N)];elseNM=0;endk=floor(N/2);B=zeros(k,3);A=zeros(k,3);if k*2==Nb=[b,0];a=[a,0];endbroots=cplxpair(roots(b));aroots=cplxpair(roots(a));for i=1:2:2*kbr=broots(i:1:i+1,:);br=real(polt(br));B((fix(i+1)/2),:)=br;ar=aroots(i:1:i+1,:);ar=real(polt(ar));A((fix(i+1)/2),:)=ar;Endfunction [db,mag,pha,grd,w]=freqz_m(b,a)[h,w]=freqz(b,a,1000,'whole');h=(h(1:501))';w=(w(1:501))';mag=abs(h);db=20*log10((mag+eps)/max(mag));pha=angle(h);grd=grdelay(b,a,w);设计一个巴特沃思高通滤波器,要求通带截止频率为0.6pi,通带内衰减不大于1dB,阻带·起始频率为0.4pi,阻带内衰减不小于15dB,T=1:>> wp=0.6*pi;ws=0.4*pi;>> Rp=1;Rs=15;T=1;>> [N,wn]=buttord(wp/pi,ws/pi,Rp,Rs) 计算巴特沃思滤波器阶数和截止频率N =4wn =>> [b,a]=butter(N,wn,'high'); 频率变换法计算巴特沃思高通滤波器>> [C,B,A]=dir2cas(b,a)C =0.0751B =1.0000 -2.0000 1.00001.0000 -2.0000 1.0000A =1.0000 0.1562 0.44881.0000 0.1124 0.0425>> [db,mag,pha,grd,w]=freqz_m(b,a);>> subplot(2,1,1);plot(w/pi,mag);>> subplot(2,1,2);plot(w/pi,db);椭圆带通滤波器的设计--ellip函数的应用:>> ws=[0.3*pi 0.75*pi]; 数字阻带边缘频率>> wp=[0.4*pi 0.6*pi]; 数字通带边缘频率>> Rp=1;Rs=40;>> Ripple=10^(-Rp/20); 通带波动>> Attn=10^(-Rs/20); 阻带衰减>> [N,wn]=ellipord(wp/pi,ws/pi,Rp,Rs) 计算椭圆滤波器参数N =4wn =0.4000 0.6000>> [b,a]=ellip(N,Rp,Rs,wn); 数字椭圆滤波器的设计>> [b0,B,A]=dir2cas(b,a) 级联形式实现b0 =0.0197B =1.0000 1.5066 1.00001.0000 0.9268 1.00001.0000 -0.9268 1.00001.0000 -1.5066 1.0000A =1.0000 0.5963 0.93991.0000 0.2774 0.79291.0000 -0.2774 0.79291.0000 -0.5963 0.9399>> figure(1);>> [db,mag,pha,grd,w]=freqz_m(b,a);>> subplot(2,2,1);plot(w/pi,mag);>> grid on;>> subplot(2,2,3);plot(w/pi,db);grid on;>> subplot(2,2,2);plot(w/pi,pha/pi);grid on;>> subplot(2,2,4);plot(w/pi,grd);设计一个巴特沃思带阻滤波器,要求通带上下截止频率为0.8pi、0.2pi,通带内衰减不大于1dB,阻带上起始频率为0.7pi、0.4pi,阻带内衰减不小于30dB:>> wp=[0.2*pi 0.8*pi];>> ws=[0.4*pi 0.7*pi];>> Rp=1;Rs=30;>> [N,wn]=buttord(wp/pi,ws/pi,Rp,Rs);>> [b,a]=butter(N,wn,'stop');>> [C,B,A]=dir2cas(b,a)C =0.0394B =1.0000 0.3559 0.99941.0000 0.3547 1.00401.0000 0.3522 0.99541.0000 0.3499 1.00461.0000 0.3475 0.99601.0000 0.3463 1.0006A =1.0000 1.3568 0.79281.0000 1.0330 0.46331.0000 0.6180 0.17751.0000 -0.2493 0.11131.0000 -0.6617 0.37551.0000 -0.9782 0.7446>> [db,mag,pha,grd,w]=freqz_m(b,a); >> subplot(2,1,1);plot(w/pi,mag);>> subplot(2,1,2);plot(w/pi);数字低通---数字带阻:function [bz,az]=zmapping(bZ,aZ,Nz,Dz) bzord=(length(bZ)-1)*(length(Nz)-1); azord=(length(aZ)-1)*(length(Dz)-1);bz=zeros(1,bzord+1);for k=0:bzordpln=[1];for i=0:k-1pln=conv(pln,Nz);endpld=[1];for i=0:bzord-k-1pld=conv(pld,Dz);endbz=bz+bZ(k+1)*conv(pln,pld); endfor k=0:azordpln=[1];for i=0:k-1pln=conv(pln,Nz);endpld=[1];for i=0:azord-k-1pld=conv(pld,Dz);endaz=az+aZ(k+1)*conv(pln,pld); endall=az(1);az=az/az1;bz=bz/az1;线性相位FIR滤波器的幅度特性:function pzkplot(num,den)hold on;axis('square');x=-1:0.01:1;y=(1-x.^2).^0.5;y1=-(1-x.^2).^0.5;plot(x,y,'b',x,y1,'b');num1=length(num);den1=length(den);if(num1>1)z=roots(num);elsez=0;endif(den1>1)p=roots(den);elsep=0;endif(num>1&den1>1)r_max_z=max(abs(real(z)));i_max_z=max(abs(imag(z)));a_max_z=max(r_max_z,i_max_z);r_max_p=max(abs(real(p)));i_max_p=max(abs(imag(p)));a_max_p=max(r_max_p,i_max_p);a_max=max(a_max_z,a_max_p);elseif (num1>1)r_max_z=max(abs(real(z)));i_max_z=max(abs(imag(z)));a_max=max(r_max_z,i_max_z);elser_max_p=max(abs(real(p)));i_max_p=max(abs(imag(p)));a_max=max(r_max_p,i_max_p);endaxis([-a_max a_max -a_max a_max]);plot([-a_max a_max],[0 0],'b');plot([0 0],[-a_max a_max],'b');plot([-a_max a_max],[a_max a_max],'b');plot([a_max a_max],[-a_max a_max],'b');Lz=length(z);for i=1:Lz;plot(real(z(i)),imag(z(i)),'bo');endLp=length(p);for j=1:Lpplot(real(p(j)),imag(p(j)),'bx');endtitle('The zeros-pole plot');xlabel('虚部');ylabel('实部');function [Hr,w,a,L]=Hr_Type1(h)M=length(h);L=(M-1)/2;a=[h(L+1) 2*h(L:-1:1)];n=[0:1:L];w=[0:1:500]'*pi/500;Hr=cos(w*n)*a';设计I型线性相位FIR滤波器:>> h=[-4 1 -1 -2 5 6 5 -2 -1 1 -4];>> M=length(h);n=0:M-1;>> [Hr,w,a,L]=Hr_Type1(h);>> amax=max(a)+1;>> amin=min(a)-1;>> subplot(2,2,1);stem(n,h);>> axis([-1 2*L+1 amin amax]);text(2*L+1.5,amin,'n'); >> xlabel('n');ylabel('h(n)');title('脉冲响应');>> subplot(2,2,3);stem(0:L,a);>> axis([-1 2*L+1 amin amax]);>> xlabel('n');ylabel('a(n)');title('a(n) 系数');>> subplot(2,2,2);plot(w/pi,Hr);>> grid on;text(1.05,-20,'频率pi');>> xlabel('频率');ylabel('Hr');title('I 型振幅响应');>> subplot(2,2,4);pzkplot(h,1);>> title('零极点分布');function [hr,w,b,L]=Hr_Type2(h)M=length(h);L=M/2;b=2*h(L:-1:1);n=[1:1:L];n=n-0.5;w=[0:1:500]'*pi/500;hr=cos(w*n)*b';II型线性相位FIR滤波器:>> h=[-4 1 -1 -2 5 6 5 -2 -1 1 -4];>> M=length(h);n=0:M-1;>> [Hr,w,b,L]=Hr_Type2(h);Warning: Integer operands are required for colon operator when used as index. > In Hr_Type2 at 2>> bmax=max(b)+1;bmin=min(b)-1;>> subplot(2,2,1);stem(n,h);axis([-1 2*L+1 bmin bmax]);text(2*L+1.5,bmin,'n');xlabel('n');ylabel('h(n)');title('脉冲响应');>> subplot(2,2,3);stem(1:L,b);axis([-1 2*L+1 bmin bmax]);xlabel('n');ylabel('b(n)');title('b(n) 系数');>> subplot(2,2,2);plot(w/pi,Hr);grid on;text(1.05,-20,'频率pi');xlabel('频率');ylabel('Hr');title('II 型振幅响应');>> subplot(2,2,4);pzkplot(h,1);title('零极点分布');function [hr,w,c,L]=Hr_Type3(h)M=length(h);L=(M-1)/2;b=2*h(L+1:-1:1);n=[1:1:L];w=[0:1:500]'*pi/500;hr=cos(w*n)*c';用MA TLAB编程绘制各种窗函数的形状。
光电图像课程设计报告书课题名称基于matlab的滤波器设计图像复原的MATLAB实现1课程设计目的〔1〕了解基于matlab的滤波器处理及其根本操作;〔2〕学习MATLAB在滤波器中的使用;〔3〕提高学习与解决问题的能力。
2课程设计根本内容2.1滤波器的根本原理设计数字滤波器的任务就是寻求一个因果稳定的线性时不变系统,并使系统函数H〔z〕具有指定的频率特性。
数字滤波器从实现的网络构造或者从单位冲激响应分类,可以分成无限单位冲激响应〔IIR〕数字滤波器和有限长单位冲激响应〔FIR〕数字滤波器。
数字滤波器频率响应的三个要素:(1)幅度平方响应(2)相位响应(3)群时延响应IIR数字滤波器:IIR数字滤波器的系统函数为有理分数,即IIR数字滤波器的逼近问题就是求解滤波器的系数和,使得在规定的物理意义上逼近所要求的特性的问题。
如果是在s平面上逼近,就得到模拟滤波器,如果是在z平面上逼近,则得到数字滤波器。
FIR数字滤波器:设FIR的单位脉冲响应h〔n〕为实数,长度为N,则其z变换和频率响应分别为按频域采样定理FIR数字滤波器的传输函数H(z)和单位脉冲响应h〔z〕可由它的N歌频域采值H(k)唯一确定。
MATLAB中提供了几个函数,分别用于实现IIR滤波器和FIR滤波器。
(1)卷积函数conv,调用格式为,c=conv〔a,b〕该格式可以计算两向量a和b的卷积,可以直接用于对有限长信号采用FIR滤波器的滤波。
(2)函数filter的调用格式为,y=filter〔b,a,*〕该格式采用数字滤波器对数据进展滤波,既可以用于IIR滤波器,也可以用于FIR滤波器。
其中向量b和a分别表示系统函数的分子,分母多项式的系数,假设a=1,此时表示FIR滤波器,否则就是IIR滤波器。
该函数就是利用给出的向量b和a,对*中的数据进展滤波,结果放入向量y。
(3)函数fftfilt的调用格式为,y=fftfilt〔b,*〕该格式是利用基于FFT的重叠相加法对数据进展滤波,这种频域滤波技术只对FIR滤波器有效。
用MATLAB 设计滤波器1 IIR 滤波器的设计freqz功能:数字滤波器的频率响应。
格式:[h ,w ]=freqz (b ,a,n )[h ,f]=freqz(b ,a ,n ,Fs)[h ,w ]=freqz(b ,a,n ,’whole')[h ,f ]=freqz(b,a ,n ,'whole ’,Fs )h=freqz (b ,a ,w)h=freqz (b,a ,f ,Fs)freqz(b ,a)说明:freqz 用于计算由矢量"和b 构成的数字滤波器H (z)=A(z)B(z)= n-1--n -1 l)z a(n ....a(2)z l l)z b(n .... b(2)z b(l)++++++++ 的复频响应H (j ω).[h ,w]=freqz (b,a ,n )可得到数字滤波器的n 点的幅频响应,这n 个点均匀地分布在上半单位圆(即0~π),并将这n 点频率记录在w 中,相应的频率响应记录在h 中。
至于n值的选择没有太多的限制,只要n 〉0的整数,但最好能选取2的幂次方,这样就可采用FFT 算法进行快速计算。
如果缺省,则n=512。
[h ,f ]二freqz(b,a,n ,Fs)允许指定采样终止频率Fs (以Hz 为单位),也即在0~Fs/2频率范围内选取n 个频率点(记录在f 中),并计算相应的频率响应h 。
[h,w]=freqz(b,a,n,’whole’)表示在0~2π之间均匀选取n个点计算频率响应.[h,f]=freqz(b,a,n,'whole',Fs)则在O~Fs之间均匀选取n个点计算频率响应.h=freqz(b,a,w)计算在矢量w中指定的频率处的频率响应,但必须注意,指定的频率必须介于0和2π之间.h=freqz(b,a,f,Fs)计算在矢量f中指定的频率处的频率响应,但指定频率必须介于0和Fs之间。
butter功能:Butterworth(比特沃思)模拟和数字滤波器设计。
matlabfir滤波器设计由于滤波器设计的复杂性,在matlab中完成滤波器设计过程可以减少大量的计算工作,下面将介绍在matlab中的FIR滤波器的设计。
一、matlab入口1. matlab工作环境:(1)安装matlab软件,在matlab中运行“dspstart”,启动数字信号处理的工作环境,包括了一些常用处理的函数和示例程序,在这种环境下可以进行滤波器设计的研究。
(2)matlab可以完成离散时间信号处理,因此,滤波器的设计需要在离散时间信号领域完成,即用离散时间信号的算法来完成滤波器设计。
2. matlab工具箱:(1)安装滤波器设计工具箱,这是一个专门用于滤波器设计的matlab库,包括了四种基本滤波器:抗谐波滤波器(notched filter),低通滤波器(low pass filter),带阻滤波器(band stop filter),以及带通滤波器(band pass filter)。
(2)安装滤波器设计高级工具箱,这是一个高级的matlab 库,可以支持滤波器的更复杂的设计,包括滤波器的参数优化,非线性滤波器,以及多种滤波器结构的数字滤波器。
3. matlab辅助性的工具:(1)安装滤波器设计助手,这是一个滤波器设计的可视化工具,可以帮助滤波器的初学者快速上手,了解滤波器设计的基本思想,并运行滤波器的例子及示例程序。
(2)安装滤波器设计指南,这是一个matlab脚本文件,可以帮助滤波器的更高级的使用者了解滤波器设计的原理,以及各种滤波器的不同之处,可以大大提高滤波器的设计效率。
二、 matlab 滤波器设计思路1. FIR滤波器的设计:(1)首先要明确滤波器的工作频率范围,记录滤波器的输入信号精度以及输出精度。
(2)设计滤波器的类型,选择适用于该频率范围内的滤波器类型,如低通滤波器,带阻滤波器等。
(3)定义滤波器的功率频谱,根据滤波器的类型,定义滤波器的功率频谱,包括最大值,超调因数以及滤波器的边带幅度。
基于MATLAB 的数字滤波器的设计1 引言数字滤波器是指完成信号滤波处理功能的,用有限精度算法实现的离散线性非时变系统,其输入是一组(由模拟信号取样和量化的)数字量,其输出是经过变换或说处理的另一组数字量。
数字滤波器具有稳定性高、精度高、灵活性大等突出优点。
这里所说的数字滤波器是指理想带通,低通等的频率选择数字滤波器。
数字滤波器设计的一个重要步骤是确定一个可实现的传输函数H(z),这个确定传输函数H(z)的过程称为数字滤波器设计。
数字滤波器的一般设计过程为:(1)按照实际需要,确定滤波器的性能要求(通常在频域内给定数字滤波的性能要求)。
(2)寻找一满足预定性能要求的离散时间线性系统。
(3)用有限精度的运算实现所设计的系统。
(4)通过模拟,验证所设计的系统是否符合给定性能要求。
2 数字滤波器的设计滤波器分为两种,分别为模拟滤波器和数字滤波器。
数字滤波器可以理解为是一个计算程序或算法,将代表输入信号的数字时间序列转化为代表输出信号的数字时间序列,并在转化的过程中,使信号按预定的形式变化。
数字滤波器有多种分类,从数字滤波器功能上分可分为低通、高通、带阻、带通滤波器,根据数字滤波器冲激响应的时域特征,可将数字滤波器分为两种,即无限长冲激响应滤波器(IIR )和有限长冲激响应滤波器(FIR )。
数字滤波器指标:一般来说,滤波器的幅频特性是分段常数的,以低通为例,在通带内逼近于1,阻带内逼近与0,实际设计的滤波器并非是锐截止的通带和阻带两个范围,两者之间总有一个过渡带。
在设计滤波器时事先给定幅频特性允许误差,在通带范围内幅度响应以误差逼近于1,在阻带内幅1σ度响应以误差逼近于0。
2σ (1)πσσ≤≤≤≤≤≤-w w e H w w e H r jwc jw ,2|)(|,1|)(|11式中wc 和wr 分别为通带边界频率和阻带边界频率,wr-wc 为过渡带。
在具体的技术指标中往往用通带波动来表示,用最小阻带衰减At 来表示,其具体的对应公式这里就不详述了。
matlab滤波器设计命令Matlab滤波器设计命令滤波器是数字信号处理中常用的工具,用于去除信号中的噪声、频率干扰或其他不需要的成分。
Matlab提供了一系列有用的滤波器设计命令,使用户能够轻松设计并应用各种类型的滤波器。
在本文中,我们将详细介绍Matlab中常用的滤波器设计命令,包括滤波器设计函数、滤波器类型和设计过程。
I. Matlab中常用的滤波器设计函数在Matlab中,有几种函数可用于设计滤波器,其中最常用的函数是`designfilt`函数和`fir1`函数。
1. designfilt函数`designfilt`函数是Matlab中最灵活和功能强大的滤波器设计函数之一,可用于设计各种类型的IIR和FIR滤波器。
它的基本语法如下:`filt = designfilt(FilterType, 'PropertyName', PropertyValue, ...)`其中,`FilterType`代表滤波器类型,包括低通滤波器(Lowpass)、高通滤波器(Highpass)、带通滤波器(Bandpass)、带阻滤波器(Bandstop)等。
`PropertyName`和`PropertyValue`是可选的参数,用于设置滤波器的各种属性,如阶数(Order)、截止频率(CutoffFrequency)、通带和阻带的最大衰减(MaximumAttenuation)等。
下面是一个使用`designfilt`函数设计低通滤波器的例子:Fs = 1000; 采样频率Fpass = 20; 通带截止频率Fstop = 30; 阻带截止频率designfilt('lowpassiir', 'FilterOrder', 4, 'PassbandFrequency', Fpass, 'StopbandFrequency', Fstop, 'SampleRate', Fs)该命令将设计一个4阶的低通IIR滤波器,其通带截止频率为20Hz,阻带截止频率为30Hz,采样频率为1000Hz。
以MATLAB 用频率采样法设计一带通数字滤波器。
低通带边缘:w1p=0.2*pi,低阻带边缘:w1s=0.35*pi,高通带边缘w2p=0.65*pi,高阻带边缘w2s=0.8*pi;设计过渡带中的频率样本值为T1=0.109021,T2=0.59417456cheb1% wp1=0.45*pi;wp2=0.65*pi;ws1=0.3*pi;ws2=0.75*pi;Rp=1;Rs=40% =============双线型变换法========================================= wp1=0.45*pi; wp2=0.65*pi;ws1=0.3*pi; ws2=0.75*pi;Rp=1; Rs=40;Wp1=tan(wp1/2); Wp2=tan(wp2/2);Ws1=tan(ws1/2); Ws2=tan(ws2/2);BW=Wp2-Wp1; W0=Wp1*Wp2; W00=sqrt(W0);WP=1; WS=WP*(W0^2-Ws1^2)/(Ws1*BW);[N,Wn]=cheb1ord(WP,WS,Rp,Rs,'s');[B,A]=cheby1(N,Rp,Wn,'s');[BT,AT]=lp2bp(B,A,W00,BW);[num,den]=bilinear(BT,AT,0.5);[h,omega]=freqz(num,den,64);subplot(2,2,1);stem(omega/pi,abs(h));xlabel('\omega/\pi');ylabel('|H(z)|');subplot(2,2,2);stem(omega/pi,20*log10(abs(h)));xlabel('\omega/\pi');ylabel('增益.dB');% =============直接法=================================wp1=0.45*pi; wp2=0.65*pi;ws1=0.3*pi; ws2=0.75*pi;Rp=1; Rs=40;Wp=[wp1/pi,wp2/pi]; Ws=[ws1/pi,ws2/pi];[N,Wn]=cheb1ord(Wp,Ws,Rp,Rs);[B,A]=cheby1(N,Rp,Wn);[h,omega]=freqz(B,A,64);subplot(2,2,3);stem(omega/pi,abs(h));xlabel('\omega/\pi');ylabel('|H(z)|');subplot(2,2,4);stem(omega/pi,20*log10(abs(h)));xlabel('\omega/\pi');ylabel('增益.dB');%cheby2%% wp1=0.45*pi;wp2=0.65*pi;ws1=0.3*pi;ws2=0.75*pi;Rp=1;Rs=40% =============双线型变换法=========================================wp1=0.45*pi; wp2=0.65*pi;ws1=0.3*pi; ws2=0.75*pi;Rp=1; Rs=40;Wp1=tan(wp1/2); Wp2=tan(wp2/2);Ws1=tan(ws1/2); Ws2=tan(ws2/2);BW=Wp2-Wp1; W0=Wp1*Wp2; W00=sqrt(W0);WP=1; WS=WP*(W0^2-Ws1^2)/(Ws1*BW);[N,Wn]=cheb2ord(WP,WS,Rp,Rs,'s');[B,A]=cheby2(N,Rs,Wn,'s');[BT,AT]=lp2bp(B,A,W00,BW);[num,den]=bilinear(BT,AT,0.5);[h,omega]=freqz(num,den,64);subplot(2,2,1);stem(omega/pi,abs(h));xlabel('\omega/\pi');ylabel('|H(z)|');subplot(2,2,2);stem(omega/pi,20*log10(abs(h)));axis([0 1 -100 0]);xlabel('\omega/\pi');ylabel('增益.dB');% =============直接法=================================wp1=0.45*pi; wp2=0.65*pi;ws1=0.3*pi; ws2=0.75*pi;Rp=1; Rs=40;Wp=[wp1/pi,wp2/pi]; Ws=[ws1/pi,ws2/pi];[N,Wn]=cheb2ord(Wp,Ws,Rp,Rs);[B,A]=cheby2(N,Rs,Wn);[h,omega]=freqz(B,A,64);subplot(2,2,3);stem(omega/pi,abs(h));xlabel('\omega/\pi');ylabel('|H(z)|');subplot(2,2,4);stem(omega/pi,20*log10(abs(h)));axis([0 1 -100 0]);xlabel('\omega/\pi');ylabel('增益.dB');%椭圆%% wp1=0.45*pi;wp2=0.65*pi;ws1=0.3*pi;ws2=0.75*pi;Rp=1;Rs=40% =============双线型变换法========================================= wp1=0.45*pi; wp2=0.65*pi;ws1=0.3*pi; ws2=0.75*pi;Rp=1; Rs=40;Wp1=tan(wp1/2); Wp2=tan(wp2/2);Ws1=tan(ws1/2); Ws2=tan(ws2/2);BW=Wp2-Wp1; W0=Wp1*Wp2; W00=sqrt(W0);WP=1; WS=WP*(W0^2-Ws1^2)/(Ws1*BW);[N,Wn]=ellipord(WP,WS,Rp,Rs,'s');[B,A]=ellip(N,Rp,Rs,Wn,'s');[BT,AT]=lp2bp(B,A,W00,BW);[num,den]=bilinear(BT,AT,0.5);[h,omega]=freqz(num,den,64);subplot(2,2,1);stem(omega/pi,abs(h));grid;xlabel('\omega/\pi');ylabel('|H(z)|');subplot(2,2,2);stem(omega/pi,20*log10(abs(h)));grid;xlabel('\omega/\pi');ylabel('增益.dB');% =============直接法================================= wp1=0.45*pi; wp2=0.65*pi;ws1=0.3*pi; ws2=0.75*pi;Rp=1; Rs=40;Wp=[wp1/pi,wp2/pi]; Ws=[ws1/pi,ws2/pi];[N,Wn]=ellipord(Wp,Ws,Rp,Rs);[B,A]=ellip(N,Rp,Rs,Wn);[h,omega]=freqz(B,A,64);subplot(2,2,3);stem(omega/pi,abs(h));grid;xlabel('\omega/\pi');ylabel('|H(z)|');subplot(2,2,4);stem(omega/pi,20*log10(abs(h)));grid;xlabel('\omega/\pi');ylabel('增益.dB');一个例程,只需要改一改参数就行了clear all;f=[0 0.19 0.2 0.3 0.31 0.59 0.6 0.8 0.81 1];% 给定频率轴分点;m=[0 0 1 1 0 0 1 1 0 0];% 给定在这些频率分点上理想的幅频响应N1=30;N2=90;% 取两种不同的滤波器长度;b1=fir2(N1,f,m);b2=fir2(N2,f,m);% 得到两个滤波器;subplot(311);stem(b1,'.');grid;subplot(312);stem(b2,'.');grid;M=128;[h1,w]=freqz(b1,1,M,1);[h2,w]=freqz(b2,1,M,1);subplot(313);plot(w,abs(h1),'b-',w,abs(h2),'g-');grid;其中,f是归依化以后的频率通过数字滤波器的采样频率算出来,根据通带和阻带算好f和m就行了看一看help,这个函数应该有窗函数的选择默认情况下是汉明窗设计一个模拟的butterworth低通滤波器,设阶数是N,OmegaC=50*pi*2;[z,p,k]=buttap(N);p=p*OmegaC;k=k*OmegaC^N;B=real(poly(z));b0=k;b=k*B;a=real(poly(p));这样一个最简单的低通滤波器就设计好了,b是传递函数的分子系数,a是分母,都是按s的降幂排列的。
在 MATLAB 中,设计滤波器通常可以使用designfilt函数或者相关的设计工具箱函数。
下面是一个简单的示例,演示如何使用designfilt函数设计一个低通滤波器:
在这个示例中,我们使用designfilt函数创建了一个低通滤波器对象,并指定了滤波器的设计参数,如滤波器类型、截止频率、设计方法和阶数。
然后,通过freqz函数绘制了滤波器的频率响应。
你可以根据需要调整滤波器的设计参数,比如滤波器类型、截止频率、设计方法和阶数等。
MATLAB 中还有其他设计工具箱函数,如fir1、butter、cheby1等,用于设计不同类型的滤波器。
请注意,滤波器设计的具体参数取决于你的应用场景和信号特性。
在选择滤波器参数时,通常需要考虑滤波器的类型、频率响应、阶数等因素。
二维高斯滤波器图像滤波MTATLAB程序设计用MATLAB设计一个3×3模板标准差为0.5,1.5,2.5的二维高斯低通滤波器(Gaussion low pass filter分别对灰度图像,真彩色图像,伪彩色图像进行滤波处理%二维高斯低通滤波器%文件为glpf.mfunction glpf(J)%处理索引图像if isind(J)[I,map]=imread(J);I=imnoise(I,'gaussian',0.005);subplot(2,2,1)imshow(I,map);title('a.原始图像');K1=filter2(fspecial('gaussian',3*3,0.5),I);K2=filter2(fspecial('gaussian',3*3,1.5),I);K3=filter2(fspecial('gaussian',3*3,2.5),I);subplot(2,2,2)imshow(K1,map);title('b.滤波器3*3,sigma=0.5')subplot(2,2,3)imshow(K2,map);title('c.滤波器3*3,sigma=1.5')subplot(2,2,4)imshow(K3,map);title('d.滤波器3*3,sigma=2.5')end%处理灰度,真彩色图像if ~isind(J)J=imread(J);%读入图像I=imnoise(J,'gaussian');%加入高斯噪声%输出原始图像subplot(2,2,1)imshow(I);title('a.原始图像')h1 = fspecial('gaussian',[3,3], 0.5);%用预定义的gaussian函数[m n p]=size(I);%[m n p] m,n代表像素点的位置,p代表空间分量%p=1,代表亮度分量,即灰度图像%p=3,代表亮度分量(Y),色差分量(Cr,Cb)if p==1%处理灰度图像I=double(I);I=conv2(I,h1,'same');%I与h的二维离散卷积endif p==3%处理真彩色I=double(I);I(:,:,1)=conv2(I(:,:,1),h1,'same');I(:,:,2)=conv2(I(:,:,2),h1,'same');I(:,:,3)=conv2(I(:,:,3),h1,'same');endI=uint8(I);%8位无符号整型%经过3*3,sigma=0.5二维高斯低通滤波器滤波后的图像subplot(2,2,2)imshow(I);title('b.滤波器3*3,sigma=0.5')h2 = fspecial('gaussian',[3,3],1.5);[m n p]=size(I);if p==1I=double(I);I=conv2(I,h2,'same');%二维离散卷积endif p==3I=double(I);I(:,:,1)=conv2(I(:,:,1),h2,'same');I(:,:,2)=conv2(I(:,:,2),h2,'same');I(:,:,3)=conv2(I(:,:,3),h2,'same');endI=uint8(I);%经过3*3,sigma=1.5二维高斯低通滤波器滤波后的图像subplot(2,2,3)imshow(I);title('c.滤波器3*3,sigma=1.5')h3 = fspecial('gaussian',[3,3], 2.5);[m n p]=size(I);if p==1I=double(I);I=conv2(I,h3,'same');endif p==3I=double(I);I(:,:,1)=conv2(I(:,:,1),h3,'same');I(:,:,2)=conv2(I(:,:,2),h3,'same');I(:,:,3)=conv2(I(:,:,3),h3,'same');endI=uint8(I);%经过3*3,sigma=2.5二维高斯低通滤波器滤波后的图像subplot(2,2,4)imshow(I);title('d.滤波器3*3,sigma=2.5')end程序调用%glpfUse.m%调用function glpf()处理灰度,真彩色,索引(伪彩色)图像并绘图glpf('填写你需要处理图像的路径')%灰度图像glpf('填写你需要处理图像的路径')%真彩色图像glpf('填写你需要处理图像的路径')%索引(伪彩色)图像本程序适用于数字图像处理多媒体课程设计。
matlabiir滤波器设计
IIR数字滤波器的设计一般采用间接法,即先设计一个模拟滤波器,再将其离散化得到所需的数字滤波器。
在MATLAB中,IIR数字滤波器的设计步骤如下:
1. 确定采样间隔Ts或采样频率fs。
2. 根据模拟频率和数字频率之间的关系,将所给出的数字滤波器的指标转化为模拟滤波器的指标(Ω=ω/Ts)。
3. 根据模拟滤波器的指标设计模拟滤波器。
4. 根据冲激响应不变法或双线性变换法,将H(s)转化为H(z)。
IIR数字滤波器的设计需要综合考虑滤波器的技术指标和设计方法,选择合适的滤波器类型和阶数,以满足实际应用的需求。
Fs=1000;
t=0.1:1/Fs:0.3;
A=1;
B=fir1(100,[0.12 0.18],'bandpass');%fir1基于窗函数的有限脉冲响应滤波器
X=3*sin(2*pi*80*t)+sqrt(5)*randn(1,length(t)); subplot(2,2,1)
plot(t,X)
title('叠加白噪声的正弦信号')
xlabel('time (seconds)')
X1=3*sin(2*pi*80*t);
subplot(2,2,2)
plot(t,X1)
title('正弦信号')
xlabel('time (seconds)')
Y=filter(B,A,X);
subplot(2,2,3)
plot(t,Y)
title('Hamming窗滤波后的信号')
xlabel('time (seconds)')
B=fir1(100,[0.12 0.18],'bandpass',kaiser(101));
Y1=filter(B,A,X);
subplot(2,2,4)
plot(t,Y1)
title('bartlett窗滤波后的信号')
xlabel('time (seconds)')
B=fir1(100,[0.12 0.18],'bandpass',barthannwin(101)); Y2=filter(B,A,X);
subplot(2,2,4)
plot(t,Y1)
title('barthannwin窗滤波后的信号')
xlabel('time (seconds)')
figure(2)
Xk=fft(X1);
Mk=abs(Xk);
subplot(2,2,1)
stem(Mk)
title('正弦信号的频谱分析')
xlabel('Frequency (Hz)')
ylabel('|X1(f)|')
Xk=fft(Y);
Mk=abs(Xk);
subplot(2,2,2)
stem(Mk)
title('Hamming窗滤波后正弦信号的频谱分析') xlabel('Frequency (Hz)')
ylabel('|Y(f)|')
Xk1=fft(Y1);
Mk1=abs(Xk1);
subplot(2,2,3)
stem(Mk1)
title('bartlett窗滤波后正弦信号的频谱分析') xlabel('Frequency (Hz)')
ylabel('|Y1(f)|')
Xk2=fft(Y2);
Mk2=abs(Xk2);
subplot(2,2,4)
stem(Mk2)
title('barthannwin窗滤波后正弦信号的频谱分析') xlabel('Frequency (Hz)')
ylabel('|Y2(f)|')。