《matlab中关于数字滤波器的函数介绍.docx》由会员分享，可在线阅读，更多相关《matlab中关于数字滤波器的函数介绍.docx（10页珍藏版）》请在三一办公上搜索。
1、matlab中关于数字滤波器的函数介绍MATLAB下的数字信号处理实现示例 一 信号、系统和系统响应 1、理想采样信号序列 首先产生信号x(n),0=n=50 n=0:50; %定义序列的长度是50 A=444.128; %设置信号有关的参数 a=50*sqrt(2.0)*pi; T=0.001; %采样率 w0=50*sqrt(2.0)*pi; x=A*exp(-a*n*T).*sin(w0*n*T); %pi是MATLAB定义的，信号乘可采用“.*” close all %清除已经绘制的x(n)图形 subplot(3,1,1);stem(x); %绘制x(n)的图形 title(理想采样
2、信号序列); 绘制信号x(n)的幅度谱和相位谱 k=-25:25; W=(pi/12.5)*k; X=x*(exp(-j*pi/12.5).(n*k); magX=abs(X); %绘制x(n)的幅度谱 subplot(3,1,2);stem(magX);title(理想采样信号序列的幅度谱); angX=angle(X); %绘制x(n)的相位谱 subplot(3,1,3);stem(angX) ; title (理想采样信号序列的相位谱) 改变参数为：1,0734.2,4.0,10=TA n=0:50; %定义序列的长度是50 A=1; %设置信号有关的参数 a=0.4; T=1; %采
3、样率 w0=2.0734; x=A*exp(-a*n*T).*sin(w0*n*T); %pi是MATLAB定义的，信号乘可采用“.*” close all %清除已经绘制的x(n)图形 subplot(3,1,1);stem(x); %绘制x(n)的图形 title(理想采样信号序列); k=-25:25; W=(pi/12.5)*k; X=x*(exp(-j*pi/12.5).(n*k); magX=abs(X); %绘制x(n)的幅度谱 subplot(3,1,2);stem(magX);title(理想采样信号序列的幅度谱); angX=angle(X); %绘制x(n)的相位谱 su
4、bplot(3,1,3);stem(angX) ; title (理想采样信号序列的相位谱) 2、单位脉冲序列 在MatLab中，这一函数可以用zeros函数实现： n=1:50; %定义序列的长度是50 x=zeros(1,50); %注意：MATLAB中数组下标从1开始 x(1)=1; close all; subplot(3,1,1);stem(x);title(单位冲击信号序列); k=-25:25; X=x*(exp(-j*pi/12.5).(n*k); magX=abs(X); %绘制x(n)的幅度谱 subplot(3,1,2);stem(magX);title(单位冲击信号的幅
5、度谱); angX=angle(X); %绘制x(n)的相位谱 subplot(3,1,3);stem(angX) ; title (单位冲击信号的相位谱) 3、矩形序列 n=1:50 x=sign(sign(10-n)+1); close all; subplot(3,1,1);stem(x);title(单位冲击信号序列); k=-25:25; X=x*(exp(-j*pi/25).(n*k); magX=abs(X); %绘制x(n)的幅度谱 subplot(3,1,2);stem(magX);title(单位冲击信号的幅度谱); angX=angle(X); %绘制x(n)的相位谱 s
6、ubplot(3,1,3);stem(angX) ; title (单位冲击信号的相位谱) 4、特定冲击串：)32(5.2)1(5.2)(+=nnnnnx n=1:50; %定义序列的长度是50 x=zeros(1,50); %注意：MATLAB中数组下标从1开始 x(1)=1;x(2)=2.5;x(3)=2.5;x(4)=1; close all; subplot(3,1,1);stem(x);title(单位冲击信号序列); k=-25:25; X=x*(exp(-j*pi/12.5).(n*k); magX=abs(X); %绘制x(n)的幅度谱 subplot(3,1,2);stem(
7、magX);title(单位冲击信号的幅度谱); angX=angle(X); %绘制x(n)的相位谱 subplot(3,1,3);stem(angX) ; title (单位冲击信号的相位谱) 5、卷积计算： +=mmnhmxnhnxny)( 在MATLAB中。提供了卷积函数conv，即y=conv(x,h)，调用十分方便。例如： 系统：)32(5.2)1(5.2)(+=nnnnnhb 信号： 500),sin(0=nnTAetxnTa n=1:50; %定义序列的长度是50 hb=zeros(1,50); %注意：MATLAB中数组下标从1开始 hb(1)=1;hb(2)=2.5;hb(
8、3)=2.5;hb(4)=1; close all; subplot(3,1,1);stem(hb);title(系统hbn); m=1:50; %定义序列的长度是50 A=444.128; %设置信号有关的参数 a=50*sqrt(2.0)*pi; T=0.001; %采样率 w0=50*sqrt(2.0)*pi; x=A*exp(-a*m*T).*sin(w0*m*T); %pi是MATLAB定义的，信号乘可采用“.*” subplot(3,1,2);stem(x);title(输入信号xn); y=conv(x,hb); subplot(3,1,3);stem(y);title(输出信号
9、yn); 6、卷积定律验证 k=-25:25; X=x*(exp(-j*pi/12.5).(n*k); magX=abs(X); %绘制x(n)的幅度谱 subplot(3,2,1);stem(magX);title(输入信号的幅度谱); angX=angle(X); %绘制x(n)的相位谱 subplot(3,2,2);stem(angX) ; title (输入信号的相位谱) Hb=hb*(exp(-j*pi/12.5).(n*k); magHb=abs(Hb); %绘制hb(n)的幅度谱 subplot(3,2,3);stem(magHb);title(系统响应的幅度谱); angHb=
10、angle(Hb); %绘制hb(n)的相位谱 subplot(3,2,4);stem(angHb) ; title (系统响应的相位谱) n=1:99; k=1:99; Y=y*(exp(-j*pi/12.5).(n*k); magY=abs(Y); %绘制y(n)的幅度谱 subplot(3,2,5);stem(magY);title(输出信号的幅度谱); angY=angle(Y); %绘制y(n)的相位谱 subplot(3,2,6);stem(angY) ; title (输出信号的相位谱) %以下将验证的结果显示 XHb=X.*Hb; Subplot(2,1,1);stem(abs
11、(XHb);title(x(n)的幅度谱与hb(n)幅度谱相乘); Subplot(2,1,2);stem(abs(Y);title(y(n)的幅度谱); axis(0,60,0,8000) 二、用FFT进行信号的频谱分析 1、高斯序列： =elsenenxqpna,0150,)(2)( n=0:15; %定义序列的长度是15 p=8;q=2; x=exp(-1*(n-p).2/q); close all; subplot(3,1,1); stem(abs(fft(x) p=8;q=4; x=exp(-1*(n-p).2/q); subplot(3,1,2); stem(abs(fft(x)
12、p=8;q=8; x=exp(-1*(n-p).2/q); subplot(3,1,3); stem(abs(fft(x) 2、衰减正弦序列： =elsenfnenxnb,0150,2sin)( n=0:15; %定义序列的长度是15 a=0.1;f=0.0625; x=exp(-a*n).*sin(2*pi*f*n); close all; subplot(2,1,1); stem(x); subplot(2,1,2); stem(abs(fft(x) 3、三角波序列： +=elsennnnnxc,074,830,1)( for i=0:3 x(i)=i+1;x(i+4)=8-(i+4);
13、end for i=8:15 x(i)=0; end close all; subplot(2,1,1); stem(x); subplot(2,1,2); stem(abs(fft(x,16) 4、反三角序列： =elsennnnnxd,074,330,4)( 三、 窗函数法设计FIR滤波器 一、在MATLAB中产生窗函数十分简单： 矩形窗 调用格式：w=boxcar(n)，根据长度n产生一个矩形窗w。 三角窗 调用格式：w=triang(n) ，根据长度n产生一个三角窗w。 汉宁窗 调用格式：w=hanning(n) ，根据长度n产生一个汉宁窗w。 海明窗 调用格式：w=hamming(n
14、) ，根据长度n产生一个海明窗w。 布拉克曼窗 调用格式：w=blackman(n) ，根据长度n产生一个布拉克曼窗w。 恺撒窗 调用格式：w=kaiser(n,beta) ，根据长度n和影响窗函数旁瓣的参数产生一个恺撒窗w。 基于窗函数的FIR滤波器设计 利用MATLAB提供的函数firl 来实现 调用格式：firl (n,Wn,ftype,Window)，n为阶数、Wn是截止频率、ftype是滤波器的类型、Window是窗函数。 例 设计一个长度为8的线性相位FIR滤波器。 其理想幅频特性满足=elseeHjd,04.00,1|)(| Window=boxcar(8); b=fir1(7,
15、0.4,Window); freqz(b,1) Window=blackman(8); b=fir1(7,0.4,Window); freqz(b,1) 例 设计线性相位带通滤波器，其长度N=15，上下边带截止频率分别为W1= 0.3，w2=0.5 Window=blackman(16); b=fir1(15,0.3 0.5,Window); freqz(b,1) 设计指标为：p=0.2 Rp=0.25dB a=0.3 As=50dB 的低通数字FIR滤波器 wp=0.2*pi;ws=0.3*pi; tr_width=ws-wp; M=ceil(6.6*pi/tr_width)+1; N=0:
16、1:M-1; wc=(ws+wp)/2; hd=ideal_lp(wc,M); w_ham=(boxcar(M); h=hd.*w_ham; db,mag,pha,grd,w=freqz_m(h,1); delta_w=2*pi/1000; Rp=-(min(db(1:1:wp/delta_w+1); As=-round(max(db(ws/delta_w+1:1:501); Close all; subplot(2,2,1);stem(hd);title(理想冲击响应) axis(0 M-1 0.1 0.3);ylabel(hdn); subplot(2,2,2);stem(w_ham);t
17、itle(汉明窗); axis(0 M-1 0 1.1);ylabel(wn); subplot(2,2,3);stem(h);title(实际冲击响应); axis(0 M-1 0.1 0.3);ylabel(hn); subplot(2,2,4);plot(w/pi,db); title(衰减幅度); axis(0 1 -100 10);ylabel(Decibles); 四、IIR滤波器的实现 MATLAB中滤波器的分析和实现 1、freqs函数：模拟滤波器的频率响应 例 系统传递函数为14.013.02.0)(22+=sssssH 的模拟滤波器，在MATLAB中可以用以下程序来实现：
18、a=1 0.4 1; b=0.2 0.3 1; %w=logspace(-1,1); %产生从到之间地0个等间距点，即50个频率点110110 freqs(b,a,w) %根据输入的参数绘制幅度谱和相位谱 2、freqz函数：数字滤波器的频率响应 例 系统传递函数为224.013.02.0)(+=zzzzzH 的模拟滤波器，在MATLAB中可以用以下程序来实现： a=1 0.4 1; b=0.2 0.3 1; %根据输入的参数绘制幅度谱和相位谱，得到0到之间128个点处的频率响应 freqz(b,a,128) 3、ButterWorth模拟和数字滤波器 butterd函数：ButterWort
19、h滤波器阶数的选择。 调用格式：n,Wn=butterd(Wp,Ws,Rp,Rs)，在给定滤波器性能的情况下，计算ButterWorth滤波器的阶数n和截止频率Wn。 相同参数条件下的模拟滤波器则调用格式为：n,Wn=butterd(Wp,Ws,Rp,Rs,s) butter函数：ButterWorth滤波器设计。 调用格式：b,a=butter(n,Wn)，根据阶数n和截止频率Wn计算ButterWorth滤波器分子分母系数。 相同参数条件下的模拟滤波器则调用格式为：b,a=butter(n,Wn,s) 例 采样频率为1Hz，通带临界频率fp =0.2Hz，通带内衰减小于1dB；阻带临界频率
20、fs=0.3Hz，阻带内衰减大于25dB。设计一个数字滤波器满足以上参数。 n,Wn=buttord(0.2,0.3,1,25); b,a=butter(n,Wn); freqz(b,a,512,1); 4、Chebyshev模拟和数字滤波器 cheb1ord函数：Chebyshev型滤波器阶数计算。 调用格式：n,Wn=cheb1ord(Wp,Ws,Rp,Rs)，在给定滤波器性能的情况下，选择Chebyshev型滤波器的最小阶n和截止频率Wn。 cheby1函数：Chebyshev型滤波器设计。 调用格式：b,a=butter(n,Rp,Wn)，根据阶数n、通带内波纹Rp和截止频率Wn计算B
21、utterWorth滤波器分子分母系数。 注：Chebyshev型滤波器所用函数和型类似，分别是cheb2ord、cheby2。 例 实现上例中的滤波器 n,Wn=cheb1ord(0.2,0.3,1,25); b,a=cheby1(n,1,Wn); freqz(b,a,512,1); 脉冲响应不变法设计数字ButterWorth滤波器 调用格式：bz,az=impinvar(b,a,Fs)，再给定模拟滤波器参数b，a和取样频率Fs的前提下，计算数字滤波器的参数。两者的冲激响应不变，即模拟滤波器的冲激响应按Fs取样后等同于数字滤波器的冲激响应。 利用双线性变换法设计数字ButterWorth滤
22、波器 调用格式：bz,az=bilinearb,a,Fs，根据给定的分子b、分母系数a和取样频率Fs，根据双线性变换将模拟滤波器变换成离散滤波器，具有分子系数向量bz和分母系数向量az。 模拟域的butter函数说明与数字域的函数说明相同 b,a=butter(n,Wn,s)可以得到模拟域的Butterworth滤波器 例 采样频率为1Hz，通带临界频率fp =0.2Hz，通带内衰减小于1dB；阻带临界频率fs=0.3Hz，阻带内衰减大于25dB。设计一个数字滤波器满足以上参数。 %直接设计数字滤波器 n,Wn=buttord(0.2,0.3,1,25); b,a=butter(n,Wn); freqz(b,a,512,1); %脉冲响应不变法设计数字滤波器 n,Wn=buttord(0.2,0.3,1,25,s); b,a=butter(n,Wn,s); freqs(b,a) bz,az=impinvar(b,a,1); freqz(bz,az,512,1) %双线性变换法设计ButterWorth数字滤波器 n,Wn=buttord(0.2,0.3,1,25,s); b,a=butter(n,Wn,s); freqs(b,a) bz,az=bilinear(b,a,1); freqz(bz,az,512,1)
