(完整word版)用MATLAB设计滤波器
使用MATLAB设计FIR滤波器

使⽤MATLAB设计FIR滤波器1. 采⽤fir1函数设计,fir1函数可以设计低通、带通、⾼通、带阻等多种类型的具有严格线性相位特性的FIR滤波器。
语法形式:b = fir1(n, wn)b = fir1(n, wn, ‘ftype’)b = fir1(n, wn, ‘ftype’, window)b = fir1(n, wn, ‘ftype’, window, ‘noscale’)参数的意义及作⽤:b:返回的FIR滤波器单位脉冲响应,脉冲响应为偶对称,长度为n+1;n:滤波器的介数;wn:滤波器的截⽌频率,取值范围为0<wn<1,1对应信号采样频率⼀半。
如果wn是单个数值,且ftype参数为low,则表⽰设计截⽌频率为wn的低通滤波器,如果ftype参数为high,则表⽰设计截⽌频率为wn的⾼通滤波器;如果wn是有两个数组成的向量[wn1wn2],ftype为stop,则表⽰设计带阻滤波器,ftype为bandpass,则表⽰设计带通滤波器;如果wn是由多个数组成的向量,则根据ftype的值设计多个通带或阻带范围的滤波器,ftype为DC-1,表⽰设计的第⼀个频带为通带,ftype为DC-0,表⽰设计的第⼀个频带为阻带;window:指定使⽤的窗函数,默认为海明窗;noscale:指定是否归⼀化滤波器的幅度。
⽰例:N=41; %滤波器长度fs=2000; %采样频率%各种滤波器的特征频率fc_lpf=200;fc_hpf=200;fp_bandpass=[200 400];fc_stop=[200 400];%以采样频率的⼀半,对频率进⾏归⼀化处理wn_lpf=fc_lpf*2/fs;wn_hpf=fc_hpf*2/fs;wn_bandpass=fp_bandpass*2/fs;wn_stop=fc_stop*2/fs;%采⽤fir1函数设计FIR滤波器b_lpf=fir1(N-1,wn_lpf);b_hpf=fir1(N-1,wn_hpf,'high');b_bandpass=fir1(N-1,wn_bandpass,'bandpass');b_stop=fir1(N-1,wn_stop,'stop');%求滤波器的幅频响应m_lpf=20*log(abs(fft(b_lpf)))/log(10);m_hpf=20*log(abs(fft(b_hpf)))/log(10);m_bandpass=20*log(abs(fft(b_bandpass)))/log(10);m_stop=20*log(abs(fft(b_stop)))/log(10);%设置幅频响应的横坐标单位为Hzx_f=0:(fs/length(m_lpf)):fs/2;%绘制单位脉冲响应%绘制单位脉冲响应subplot(421);stem(b_lpf);xlabel('n');ylabel('h(n)');subplot(423);stem(b_hpf);xlabel('n');ylabel('h(n)');subplot(425);stem(b_bandpass);xlabel('n');ylabel('h(n)');subplot(427);stem(b_stop);xlabel('n');ylabel('h(n)');%绘制幅频响应曲线subplot(422);plot(x_f,m_lpf(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);subplot(424);plot(x_f,m_hpf(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);subplot(426);plot(x_f,m_bandpass(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);subplot(428);plot(x_f,m_stop(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);2. 采⽤fir2函数设计,函数算法是:⾸先根据要求的幅频响应向量形式进⾏插值,然后进⾏傅⾥叶变换得到理想滤波器的单位脉冲响应,最后利⽤窗函数对理想滤波器的单位脉冲响应激进型截断处理,由此得到FIR滤波器系数。
matlab滤波器设计(源代码)

某合成信号,表达式如下:f=10cos(2pi*30t)+cos(2pi*150t)+5cos(2pi*600t),请设计三个滤波器,分别提取出信号中各频率分量,并分别绘制出通过这三个滤波器后信号的时域波形和频谱这个信号的频率分量分别为30、150和600Hz,因此可分别设计一个低通、带通和高通的滤波器来提取。
以FIR滤波器为例,程序如下:clear;fs=2000;t=(1:1000)/fs;x=10*cos(2*pi*30*t)+cos(2*pi*150*t)+5*cos(2*pi*600*t);L=length(x);N=2^(nextpow2(L));Hw=fft(x,N);figure(1);subplot(2,1,1);plot(t,x);grid on;title('滤波前信号x');xlabel('时间/s');% 原始信号subplot(2,1,2);plot((0:N-1)*fs/L,abs(Hw));% 查看信号频谱grid on;title('滤波前信号频谱图');xlabel('频率/Hz');ylabel('振幅|H(e^jw)|');%% x_1=10*cos(2*pi*30*t)Ap=1;As=60;% 定义通带及阻带衰减dev=[(10^(Ap/20)-1)/(10^(Ap/20)+1),10^(-As/20)];% 计算偏移量mags=[1,0];% 低通fcuts=[60,100];% 边界频率[N,Wn,beta,ftype]=kaiserord(fcuts,mags,dev,fs);% 估算FIR滤波器阶数hh1=fir1(N,Wn,ftype,kaiser(N+1,beta));% FIR滤波器设计x_1=filter(hh1,1,x);% 滤波x_1(1:ceil(N/2))=[];% 群延时N/2,删除无用信号部分L=length(x_1);N=2^(nextpow2(L));Hw_1=fft(x_1,N);figure(2);subplot(2,1,1);plot(t(1:L),x_1);grid on;title('x_1=10*cos(2*pi*30*t)');xlabel('时间/s');subplot(2,1,2);plot((0:N-1)*fs/L,abs(Hw_1));% 查看信号频谱grid on;title('滤波后信号x_1频谱图');xlabel('频率/Hz');ylabel('振幅|H(e^jw)|');%% x_2=cos(2*pi*150*t)Ap=1;As=60;% 定义通带及阻带衰减dev=[10^(-As/20),(10^(Ap/20)-1)/(10^(Ap/20)+1),10^(-As/20)];% 计算偏移量mags=[0,1,0];% 带通fcuts=[80,120,180,220];% 边界频率[N,Wn,beta,ftype]=kaiserord(fcuts,mags,dev,fs);% 估算FIR滤波器阶数hh2=fir1(N,Wn,ftype,kaiser(N+1,beta));% FIR滤波器设计x_2=filter(hh2,1,x);% 滤波x_2(1:ceil(N/2))=[];% 群延时N/2,删除无用信号部分L=length(x_2);N=2^(nextpow2(L));Hw_2=fft(x_2,N);figure(3);subplot(2,1,1);plot(t(1:L),x_2);grid on;title('x_2=cos(2*pi*150*t)');xlabel('时间/s');subplot(2,1,2);plot((0:N-1)*fs/L,abs(Hw_2));% 查看信号频谱grid on;title('滤波后信号x_2频谱图');xlabel('频率/Hz');ylabel('振幅|H(e^jw)|');%% x_3=5*cos(2*pi*600*t)Ap=1;As=60;% 定义通带及阻带衰减dev=[10^(-As/20),(10^(Ap/20)-1)/(10^(Ap/20)+1)];% 计算偏移量mags=[0,1];% 高通fcuts=[500,550];% 边界频率[N,Wn,beta,ftype]=kaiserord(fcuts,mags,dev,fs);% 估算FIR滤波器阶数hh2=fir1(N,Wn,ftype,kaiser(N+1,beta));% FIR滤波器设计x_3=filter(hh2,1,x);% 滤波x_3(1:ceil(N/2))=[];% 群延时N/2,删除无用信号部分L=length(x_3);N=2^(nextpow2(L));Hw_3=fft(x_3,N);figure(4);subplot(2,1,1);plot(t(1:L),x_3);grid on;title('x_3=5*cos(2*pi*600*t)');xlabel('时间/s');subplot(2,1,2);plot((0:N-1)*fs/L,abs(Hw_3));% 查看信号频谱grid on;title('滤波后信号x_3频谱图');xlabel('频率/Hz');ylabel('振幅|H(e^jw)|');。
matlab滤波器设计

实验四 IIR 滤波器的实现(IIR 滤波器是无穷脉冲响应数字滤波器的简称)一、 讲课目的熟悉用双线性变换法设计IIR 数字滤波器和方式。
二、讲课内容(一)大体概念1、 数字滤波器是指输入、输出均为数字信号,通过必然运算关系,改变输入信号所含频率成份的相对照例或滤除某些频率成份的器件。
数字滤波器与模拟滤波器概念相同,只是信号的形式和实现滤波的方式不同。
一样数字滤波器从现实的网络结构或从单位脉冲响应分类,能够分成无穷脉冲脉冲响应(IIR)滤波器和有限脉冲响应(FIR)滤波器。
2、 IIR 滤波器设计的要紧方式是先设计低通模拟滤波器,进行频率变换,将其转换为相应的高通、带通模拟滤波器,再将模拟滤波器转换为相应的数字滤波器。
3、 N 阶模拟滤波器系统函数的一样形式为 11101110...()()...()M M M M N N N N b s b s b s b B s H s a s a s a s a A s ----++++==++++ 4、 N 阶数字模拟滤波器系统函数的一样形式为12(1)012112(1)121...()()1...()M M M M N N N N b b z b z b z b z B z H z a z a z a z a z A z ------------+++++==+++++(二)求滤波器频率响应的函数一、freqs 函数:求模拟滤波器的频率响应freqs(b,a,w)计算由向量w (rad/s )指定的频率点上(或W 直接为采样点数)频率响应,其中b 和a 别离为模拟滤波器系统函数H(s)的分子与分母。
freqs 函数将自动绘出幅频和相频曲线。
实验4-1:系统传输函数为220.20.31()0.41s s H s s s ++=++的模拟滤波器,绘出其幅频和相频谱。
a=[1,,1];b=[,,1];w=logspace(-1,1); %产生从10-1到101之间地50个等间距点,即50个频率点运行得幅频特性和相频特性如下:二、freqz函数:求数字滤波器的频率响应freqz(b,a,w)计算由向量w指定的频率点上模拟滤波器H(z)的频率响应,其中b和a别离为数字滤波器系统函数H(z)的分子与分母。
(完整word版)巴特沃斯带阻数字滤波器设计matlab程序及仿真图 - 副本

fs=15000;T= 1/fs;rp=1;rs=40;wp1=0.11*pi;wp2=0.81*pi;ws1=0.31*pi;ws2=0.61*pi;%数字带阻滤波器技术指标wc1=(2/T)*tan(wp1/2);%频率预畸变wc2=(2/T)*tan(wp2/2);wr1=(2/T)*tan(ws1/2);wr2=(2/T)*tan(ws2/2);w0=sqrt(wc1*wc2);B=wc2-wc1;wp=1;%归一化通带截止频率ws=wp*(wr1*B) / (w0^2-wr1^2) ; %归一化阻带截止频率[N,wc]=buttord(wp,ws,rp,rs,'s')%求滤波器阶数和3dB截止频率[Z,P,K]=buttap(N)%设计模拟低通滤波器[Md,Nd]=zp2tf(Z,P,K)%将零极点形式转换为传输函数形式[M,N]=lp2bs(Md,Nd,w0,B)%对低通滤波器进行频率变换,转换为带阻滤波器[h,w]=freqs(M,N);%模拟带阻滤波器的幅频响应plot(w/(2*pi),abs(h));grid;xlabel('频率Hz');ylabel('幅度');title('模拟带阻滤波器');[b,a]=bilinear(M,N,15000)%对模拟滤波器双线性变换figure(1);freqz(b,a);[H,W]=freqz(b,a); %绘出频率响应;axis([0,1,-100,20]);figure(2);plot(W*fs/(2*pi),abs(H));grid on;xlabel('频率/Hz');ylabel('幅值');n=0:199;t=n/fs;x=sin(2*pi*400*t)+3*sin(2*pi*3000*t)+2*sin(2*pi*5000*t);figure(3);subplot(311);plot(t,x);axis([0,0.01,-5,5]);title('输入信号');grid on;y=filter(b,a,x);subplot(312);stem(y,'.');title('输出序列');grid on;ya=y*sinc(fs*(ones(length(n),1)*t-(n/fs)'*ones(1,length(t))));subplot(313);plot(t,ya);axis([0,0.01,-3,3]);title('输出波形');grid on;t=(0:100)/fs;figure(4)fs=1.5*10000;n=(0:100)/fs;f=sin(2*pi*400*t)+3*sin(2*pi*3000*t)+2*sin(2*pi*5000*t);y=fftfilt(b,x);[H1,f1]=freqz(f,[1]);[H2,f2]=freqz(y,[1]);f1=f1/pi*fs/2;f2=f2/pi*fs/2;subplot(2,1,1);plot(f1,abs(H1));title('输入信号的频谱');subplot(2,1,2);plot(f2,abs(H2));title('输出信号的频谱');基于Matlab 的带阻滤波器设计.10.20.30.40.50.60.70.80.91-800-600-400-2000N o r m a l i z e d Fre q u⨯π r a d /s a m p l e Ph a se(d e g r e e s )00.10.20.30.40.50.60.70.80.91-100-50N o r m a l i z e d Fr e q u⨯π r a d /s a m p l e M a g n i tu d e1000200030004000500060007000800000.20.40.60.811.21.4频率/Hz幅值00.0010.0020.0030.0040.0050.0060.0070.0080.0090.01-505输入信号020406080100120140160180200-22输出序列0.0010.0020.0030.0040.0050.0060.0070.0080.0090.01-202输出波形01000200030004000500060007000800050100150200输入信号的频谱010002000300040005000600070008000102030输出信号的频谱N =4wc =1.7947b =0.0186 -0.0410 0.1082 -0.1355 0.1810 -0.1355 0.1082 -0.0410 0.0186a =1.0000 -0.6707 -1.3750 0.5678 1.1964 -0.2996 -0.4631 0.0496 0.0762>。
matlab带通滤波器 (2)

MATLAB带通滤波器1. 简介带通滤波器是一种数字信号处理中常用的滤波器。
它可以选择特定的频率范围内的信号并传递,同时抑制其他频率范围的信号。
在MATLAB中,可以使用信号处理工具箱中的函数来设计和实现带通滤波器。
本文档将介绍如何使用MATLAB设计和使用带通滤波器,包括滤波器的设计方法和常见的应用场景。
2. 带通滤波器的设计带通滤波器的设计过程可以分为以下几个步骤:2.1 滤波器类型选择MATLAB中提供了多种带通滤波器类型的设计方法,包括巴特沃斯滤波器、切比雪夫滤波器和椭圆滤波器等。
根据需求选择合适的滤波器类型。
2.2 滤波器规格确定确定滤波器的通带范围、阻带范围和过渡带宽等规格参数。
2.3 滤波器设计根据滤波器类型和规格参数,使用相应的MATLAB函数进行滤波器设计。
常用的函数包括butter、cheby1和ellip等。
2.4 滤波器特性分析设计完成的滤波器可以通过频率响应、相位响应和零极点分布等特性进行分析。
MATLAB提供了函数来绘制和分析滤波器的特性曲线。
3. MATLAB中的带通滤波器函数MATLAB提供了多个函数用于设计和实现带通滤波器,下面介绍其中几个常用的函数:3.1 butter函数butter函数可用于设计巴特沃斯滤波器。
它的语法为:[b, a] = butter(n, Wn, 'bandpass')其中,n表示滤波器的阶数,Wn为通带范围,可以是一个长度为2的向量表示最低频率和最高频率的范围。
'bandpass'表示带通滤波器。
3.2 cheby1函数cheby1函数可用于设计切比雪夫滤波器。
它的语法为:[b, a] = cheby1(n, Rp, Wn, 'bandpass')其中,n表示滤波器的阶数,Rp为通带中允许的最大衰减量,Wn为通带范围,可以是一个长度为2的向量表示最低频率和最高频率的范围。
'bandpass'表示带通滤波器。
滤波器matlab

>> n=31;%定义滤波器阶数32fs=12.8*10^3;fc1=49;fc2=51;w1=2*pi*fc1/fs;w2=2*pi*fc2/fs;%参数转换,将模拟滤波器的技术指标转换为数字滤波器的技术指标window=hanning(n+1);%使用hanning窗函数q=fir1(n,[w1/pi w2/pi],hanning(n+1));%滤波器时域函数,使用标准响应的加窗设计函数fir1w=linspace(0,pi,512);h1=freqz(q,1,512);%进行512个点的傅里叶变换figure(2);plot(w/pi,20*log10(abs(h1)));title('滤波器频谱图');xlabel('频率');ylabel('幅度');grid ;设计FIR低通滤波器,系统频率为50MHz,通带截止频率Fpass为1MHz,阻带截止频率Fstop为4MHz,通带最大衰减Apass为1dB,阻带最小衰减Astop 为30dB。
程序和必要的程序注释谢谢最佳答案只要用一个公式就行。
library IEEE;use IEEE.STD_LOGIC_1164.ALL;use IEEE.STD_LOGIC_ARITH.ALL;use IEEE.STD_LOGIC_UNSIGNED.ALL;entity fir isPort (clk: in std_logic;reset: in std_logic;inpx: in std_logic_vector(11 downto 0);outy: out std_logic_vector(11 downto 0));end fir;architecture beh of fir issignal x0,x1,x2,x3: std_logic_vector(11 downto 0);constant c0:integer :=-1234*32768/1000;constant c1:integer :=2345*32768/10000;constant c2:integer :=5*32768;constant c3:integer :=-3*32768/10000;signal p0,p1,p2,p3:integer;signal sum: integer;beginsample_delay_line:process(clk)beginif rising_edge(clk) thenif reset='1' thenx3 <=(others=>'0');x2 <=(others=>'0');x1 <=(others=>'0');x0 <=(others=>'0');elsex3 <=x2;x2 <=x1;x1 <=x0;x0 <=inpx;end if;end if;end process;p0 <= conv_integer(x0)*c0;p1 <= conv_integer(x1)*c1;p2 <= conv_integer(x2)*c2;p3 <= conv_integer(x3)*c3;sum <=p0+p1+p2+p3;outy <=conv_std_logic_vector(sum/32768,12);end beh;4.1 数字滤波器简介数字滤波器是一种用来过滤时间离散信号的数字系统,通过对抽样数据进行数学处理来达到频域滤波的目的。
滤波器设计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('频域波形')```运行以上代码,可以得到原始信号和滤波后信号的时域和频域波形图。
完整版使用MATLAB设计ISE中FIR滤波器系数方法

使用 MATLAB 设计 ISE中 FIR 滤波器系数的方法
1、翻开 MATLAB,在命令行窗口输入“fdatool”,翻开“Filter Designer& Analysis Tool”工具。
以下列图所示:
2、因为 FPGA中滤波器的系数需要为整数,所以需要在此处将系数设置为“Fixed
-point ”种类。
点击上图中红色方框内的按钮,在新出现的页面中将“ Filter arithmetic ”设置为“ Fixed -point ”。
设置达成后以下列图所示:
3、点击上图中红色方框内的按钮,进入滤波器参数设置页面,在此中设置采样
频次( Fs)、通带频次( Fpass)、阻带频次( Fstop )以及阻带衰减( Astop )等参数,并按最下边的“ Design Filter”按钮生成滤波器系数。
以下列图所示,采样频次为 62MHz,通带频次为 2MHz,阻带频次为 4MHz,阻带衰减为 -80dB。
4、而后导出 coe 文件,点击下列图方框中的按钮即可导出coe 文件:
5、在 ISE 中新建一个 FIR 滤波器 IP 核,在第一页设置中将“ Select Source ”改为“ COE File ”,而后在下边选择上一步生成的 coe 文件即可。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
用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(比特沃思)模拟和数字滤波器设计。
格式:有六种[b,a]=butter(n,Wn)[b,a]=butter(n,Wn,’ftype’)[b,a]=butter(n,Wn,’s’)[b,a]=butter(n,Wn,'ftyPe',’s')[z,p,k]=butter(…)[A,B,C,D]=butter(…)说明:butter函数可设计低通、带通、高通和带阻的数字和模拟滤波器,其特性为使通带内的幅度响应最大限度地平坦,这会损失截止频率处的下降斜度.在期望通带平滑的情况下,可使用butter函数,但在期望下降斜度大的场合,应使用椭圆和Chebyshev(切比雪夫)滤波器。
butter函数可设计出数字域和模拟域的Butterworth滤波器。
(1)数字域[b,a]=butter(n,Wn)可设计出截止频率为Wn的n阶低通Butterworth滤波器,其滤波器为H(z)= A(z)B(z)=n -1--n -1l)z a(n ....a(2)z l l)z b(n .... b(2)z b(l)++++++++ 截止频率是滤波器幅度响应下降至1/2处的频率,Wn ∈[0,1],其中1相应于0。
5f s (取样频率,即Nyquist 频率)。
当Wn=[Wl W2](Wl 〈W2)时,butter 函数产生一2n 阶的数字带通滤波器,其通带为Wl<ω<W2。
[b ,a]=butter (n ,Wn ,'ftype ’)可设计出高通或带阻滤波器:·当ftype=high 时,可设计出截止频率为Wn 的高通滤波器;·当ftype=stop 时,可设计出带阻滤波器,这时Wn=[Wl W2],且阻带为Wl 〈ω〈W2利用输出变量个数的不同,可得到滤波器的另外两种表示:零极点增益和状态方程。
·[z ,p,k]=butter (n ,Wn)或[z,p ,k]=butter (n ,Wn,'ftype')可得到滤波器的零极点增益表示;·[A,B ,C ,D ]=butter (n,Wn )或[A ,B,C ,D]=butter(n ,Wn ,'ftype ’)可得到滤波器的状态空间表示。
(2)模拟域[b ,a]=butter (n,Wn ,’s ’)可设计截止频率为Wn 的n 阶低通模拟Butterworth 滤波器为H (s)= A(s)B(s)=l)a(n ....a(2)s s l)b(n .... b(2)s b(l)s 1-n n 1-n n ++++++++其中截止频率Wn 〉0。
模拟域的butter 函数说明与数字域的完全相同,它也有六种形式。
例1:设数据采样频率为900Hz ,现要设计一9阶的高通Butterworth 滤波器,截止频率为300Hz 。
%ex104.m[b,a]=butter(9,300/500,'high’);freqz(b,a,128,1000)例2:设计一10阶的带通Butterworth滤波器,带通为100--200 Hz,并画出滤波器的冲击响应.%ex105.mn=5;Wn=[100 200]/500;[b,a]=butter(n,Wn);figure(1);freqz(b,a,128,1000)figure(2);impz(b,a,101)例3利用双线性变换设计数字Butterworth滤波器πω2.0=pdBRP1=πω3.0=SdBAS15=源程序:Wp=0。
2*pi;Ws=0。
3*pi;Rp=1;As=15;T=1;Fs=1/T;OmegaP=(2/T)*tan(Wp/2); OmegaS=(2/T)*tan(Ws/2);ep=sqrt(10^(Rp/10)—1);Ripple=sqrt(1/(1+ep*ep));Attn=1/10^(As/20);[cs,ds]=adf_butt(OmegaP,OmegaS,Rp,As);[b,a]=bilinear(cs,ds,T);[C,B,A]=sdir2cas(b,a)[cs,ds]=s521(OmegaP,OmegaS,Rp,As);%调用1[b,a]=bilinear(cs,ds,T);[b0,B,A]=s522(b,a);%调用2%%plotfigure(1);[db,mag,pha,grd,w]=freqz_m(b,a);subplot(2,2,1);plot(w/pi,mag);title('Magnitude Response');xlabel(’frequency in pi units');ylabel(’|H|');axis([0,1,0,1]);set(gca,’XTickMode’,'manual’,'XTick’,[0,0。
2,0。
3,1]);set(gca,’YTickMode’,’manual','YTick’,[0,Attn,Ripple,1]);grid;subplot(2,2,3);plot(w/pi,db);title(’Magnitude in dB’);xlabel(’frequency in pi units’);ylabel('decibels');axis([0,1,-40,5]);set(gca,’XTickMode','manual’,'XTick',[0,0.2,0。
3,1]);set(gca,'YTickMode',’manual','YTick’,[-50,-15,—1,0]);grid ;%set(gca,'YTickLabelMode','manual’,'YTickLabels’,['50’;'15’;'1’;'0’]);subplot(2,2,2);plot(w/pi,pha/pi);title(’Phase Response’)xlabel(’frequency in pi units’);ylabel(’pi units');axis([0,1,-1,1]);set(gca,’XTickMode’,'manual’,’XTick',[0,0.2,0。
3,1]);set(gca,’YTickMode’,’manual’,'YTick',[—1,0,1]);gridsubplot(2,2,4);plot(w/pi,grd);title('Group Delay');xlabel(’frequency in pi units');ylabel('samples’);axis([0,1,0,10]);set(gca,’XTickMode','manual','XTick’,[0,0。
2,0.3,1]);set(gca,’YTickMode’,’manual',’YTick',[0:2:10]);grid子程序1: freqz的修正%freqz_m.mfunction[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);子程序2变滤波器直接形式为级联形式%sdir2cas。
mfunction [C,B,A]=sdir2cas(b,a)%DIRECT—form to CASCADE-form conversin in s-plane%-——-————-———--———---———--——-—————-————--—%[C,B,A]=sdir2cas(b,a)%C=gain coefficient%B=K by 3 matrix of real coefficients containing bk's%A=K by 3 matrix of real coefficients containing ak’s%b=numberator polynomial coefficients of DIRECT form%a=denominator polynomial coefficients of DIRECT form%Na=length(a)-1;Nb=length(b)-1;%compute gain coefficient Cb0=b(1);b=b/b0;a0=a(1);a=a/a0;C=b0/a0;%%Debnominator second—order sections:p=cplxpair(roots(a));K=floor(Na/2);if K*2==Na %Computation when Na is evenA=zeros(K,3);for n=1:2:NaArow=p(n:1:n+1,:);Arow=poly(Arow);A(fix((n+1)/2),:)=real(Arow);endelseif Na==1 %Computation when Na =1A=[0 real(poly(p))];else %Computation when Na is odd and>1 A=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)))];end%numerator second—order sections:z=cplxpair(roots(b));K=floor(Nb/2);if Nb==0 %computation when Nb=0B=[0 0 poly(z)];elseif K*2==Nb %computation when Nb is even B=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==1 %computation when Nb=1B=[0 real(poly(z))];else %Computation when Nb is odd and>1B=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)))];endchebyl功能:Chebyshev(切比雪夫)I型滤波器设计(通带等波纹)。