MATLAB辅助设计常用滤波器及其离散化

您所在的位置：网站首页 › fs序列 › MATLAB辅助设计常用滤波器及其离散化

MATLAB辅助设计常用滤波器及其离散化

2023-08-25 09:00| 来源: 网络整理| 查看: 265

文章目录 1. 三种滤波器的对比2. 巴特沃斯滤波器（1）低通滤波器（2）高通滤波器（3）带通滤波器（4）带阻滤波器 3. 契比雪夫I型滤波器（1）低通滤波器（2）高通滤波器（3）带通滤波器（4）带阻滤波器 4. 贝塞尔滤波器5. MATLAB中Filter type类型表示6. 滤波器的离散化

声明：本文内容整理自网络，仅作为博主学习笔记记录，各部分版权归原作者所有。

在这里插入图片描述几个概念理解：

通带与阻带交界点的频率称为截止频率通带纹波是指在滤波器的频响中通带的最大幅值和最小幅值的差，正常的滤波器一般通带纹波小于1db，不过也视情况而定，通带纹波会导致通带内的信号幅值大小有变化，对一些要求高的系统，纹波越小越好。通带纹波和滤波器的阶数有关系，阶数越大纹波越小。阻带纹波道理应该是一样的，不过好像没有人去关注阻带的纹波，主要关注的阻带参数是阻带衰减。阻带最小衰减：滤波器有部分频率是通，部分是阻。但是阻的部分，未必能够全部阻隔，而只有部分衰减，部分留下来，因此最小衰减就可以描述它阻碍该阻碍的波段的能力的高低（理想状态是100%衰减），最小衰减越大，则能力越好。同样可以理解通带最大波纹。

在这里插入图片描述 MATLAB辅助滤波器设计有多种方法：

fdatool工具箱，参考这篇文章本文使用的方法 1. 三种滤波器的对比

在这里插入图片描述

传送门：巴特沃斯、切比雪夫、贝塞尔滤波器的比较

2. 巴特沃斯滤波器

测试用原始信号设置

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % ButterWorthFilter.m % 巴特沃夫滤波器的设计 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear; close all; clc; fs = 1000; %Hz 采样频率 Ts = 1/fs; N = 1000; %序列长度 t = (0:N-1)*Ts; delta_f = 1*fs/N; f1 = 50; f2 = 100; f3 = 200; f4 = 400; x1 = 2*0.5*sin(2*pi*f1*t); x2 = 2*0.2*sin(2*pi*f2*t); x3 = 2*0.3*sin(2*pi*f3*t); x4 = 2*0.6*sin(2*pi*f4*t); x = x1 + x2 + x3 + x4; %待处理信号由四个分量组成 X = fftshift(abs(fft(x)))/N; X_angle = fftshift(angle(fft(x))); f = (-N/2:N/2-1)*delta_f; figure(1); subplot(3,1,1); plot(t,x); title('原信号'); subplot(3,1,2); plot(f,X); grid on; title('原信号频谱幅度特性'); subplot(3,1,3); plot(f,X_angle); title('原信号频谱相位特性'); grid on; （1）低通滤波器 [B,A] = butter(N,Wn,'low')

设计案例：

%设计一个巴特沃夫低通滤波器,要求把50Hz的频率分量保留,其他分量滤掉 wp = 65/(fs/2); %通带截止频率,取50~100中间的值,并对其归一化 ws = 85/(fs/2); %阻带截止频率,取50~100中间的值,并对其归一化 alpha_p = 3; %通带允许最大衰减为 db alpha_s = 20;%阻带允许最小衰减为 db %获取阶数和截止频率 [ N1 wc1 ] = buttord( wp , ws , alpha_p , alpha_s); %获得转移函数系数 [ b a ] = butter(N1,wc1,'low'); %滤波 filter_lp_s = filter(b,a,x); X_lp_s = fftshift(abs(fft(filter_lp_s)))/N; X_lp_s_angle = fftshift(angle(fft(filter_lp_s))); figure(2); freqz(b,a); %滤波器频谱特性 figure(3); subplot(3,1,1); plot(t,filter_lp_s); grid on; title('低通滤波后时域图形'); subplot(3,1,2); plot(f,X_lp_s); title('低通滤波后频域幅度特性'); subplot(3,1,3); plot(f,X_lp_s_angle); title('低通滤波后频域相位特性'); （2）高通滤波器 [B,A] = butter(N,Wn,'high')

设计案例：

%设计一个高通滤波器,要求把400Hz的频率分量保留,其他分量滤掉 wp = 350/(fs/2); %通带截止频率,取200~400中间的值,并对其归一化 ws = 380/(fs/2); %阻带截止频率,取200~400中间的值,并对其归一化 alpha_p = 3; %通带允许最大衰减为 db alpha_s = 20;%阻带允许最小衰减为 db %获取阶数和截止频率 [ N2 wc2 ] = buttord( wp , ws , alpha_p , alpha_s); %获得转移函数系数 [ b a ] = butter(N2,wc2,'high'); %滤波 filter_hp_s = filter(b,a,x); X_hp_s = fftshift(abs(fft(filter_hp_s)))/N; X_hp_s_angle = fftshift(angle(fft(filter_hp_s))); figure(4); freqz(b,a); %滤波器频谱特性 figure(5); subplot(3,1,1); plot(t,filter_hp_s); grid on; title('高通滤波后时域图形'); subplot(3,1,2); plot(f,X_hp_s); title('高通滤波后频域幅度特性'); subplot(3,1,3); plot(f,X_hp_s_angle); title('高通滤波后频域相位特性'); （3）带通滤波器 % [B,A] = butter(N,Wn) fs = 1000 ; %信号的采样频率 wp=[8 30]*2/fs; %通带边界频率，单位为rad/s ws=[7 32]*2/fs; %阻带边界频率，单位为rad/s Rp=1; %通带最大波纹度，单位dB (不要太小) Rs=30; %表示阻带最小衰减，单位dB [N,Wn]=buttord(Wp,Ws,Rp,Rs); %巴特沃斯数字滤波器的阶数n和-3dB归一化截止频率Wn [B,A]=butter(N,Wn); %得到n阶巴特沃斯滤波的分子分母 dataOut = filter(B,A,dataIn);

设计案例：

%设计一个带通滤波器,要求把50Hz和400Hz的频率分量滤掉,其他分量保留 wp = [65 385 ] / (fs/2); %通带截止频率,50~100、200~400中间各取一个值,并对其归一化 ws = [75 375 ] / (fs/2); %阻带截止频率,50~100、200~400中间各取一个值,并对其归一化 alpha_p = 3; %通带允许最大衰减为 db alpha_s = 20;%阻带允许最小衰减为 db %获取阶数和截止频率 [ N3 wn ] = buttord( wp , ws , alpha_p , alpha_s); %获得转移函数系数 [ b a ] = butter(N3,wn,'bandpass'); %滤波 filter_bp_s = filter(b,a,x); X_bp_s = fftshift(abs(fft(filter_bp_s)))/N; X_bp_s_angle = fftshift(angle(fft(filter_bp_s))); figure(6); freqz(b,a); %滤波器频谱特性 figure(7); subplot(3,1,1); plot(t,filter_bp_s); grid on; title('带通滤波后时域图形'); subplot(3,1,2); plot(f,X_bp_s); title('带通滤波后频域幅度特性'); subplot(3,1,3); plot(f,X_bp_s_angle); title('带通滤波后频域相位特性'); （4）带阻滤波器 [B,A] = butter(N,Wn,'stop')

设计案例：

%设计一个带阻滤波器,要求把50Hz和400Hz的频率分量保留,其他分量滤掉 wp = [65 385 ] / (fs/2); %通带截止频率?,50~100、200~400中间各取一个值,并对其归一化 ws = [75 375 ] / (fs/2); %阻带截止频率?,50~100、200~400中间各取一个值,并对其归一化 alpha_p = 3; %通带允许最大衰减为 db alpha_s = 20;%阻带允许最小衰减为 db %获取阶数和截止频率 [ N4 wn ] = buttord( wp , ws , alpha_p , alpha_s); %获得转移函数系数 [ b a ] = butter(N4,wn,'stop'); %滤波 filter_bs_s = filter(b,a,x); X_bs_s = fftshift(abs(fft(filter_bs_s)))/N; X_bs_s_angle = fftshift(angle(fft(filter_bs_s))); figure(8); freqz(b,a); %滤波器频谱特性 figure(9); subplot(3,1,1); plot(t,filter_bs_s); grid on; title('带阻滤波后时域图形'); subplot(3,1,2); plot(f,X_bs_s); title('带阻滤波后频域幅度特性'); subplot(3,1,3); plot(f,X_bs_s_angle); title('带阻滤波后频域相位特性'); 3. 契比雪夫I型滤波器

在这里插入图片描述

（1）低通滤波器 fs = 1000; % Hz采样频率 wp = 55/(fs/2); %通带截止频率,取50~100中间的值,并对其归一化 ws = 90/(fs/2); %阻带截止频率,取50~100中间的值,并对其归一化 Rp = 3; %通带允许最大衰减为 db Rs = 40; %阻带允许最小衰减为 db [n,Wn]=cheb1ord(Wp,Ws,Rp,Rs); % 获取阶数和截止频率 [b,a]=cheby1(n,Rp,Wn, 'low'); %获得转移函数系数 dataOut = filter(b,a,dataIn); %信号滤波运算 function y=lowp(x,f1,f3,rp,rs,Fs) %低通滤波 %使用注意事项：通带或阻带的截止频率的选取范围是不能超过采样率的一半 %即，f1,f3的值都要小于 Fs/2 %x:需要带通滤波的序列 % f 1：通带截止频率 % f 3：阻带截止频率 %rp：边带区衰减DB数设置 %rs：截止区衰减DB数设置 %FS：序列x的采样频率 % rp=0.1;rs=30;%通带边衰减DB值和阻带边衰减DB值 % Fs=2000;%采样率 % wp=2*pi*f1/Fs; ws=2*pi*f3/Fs; % 设计切比雪夫滤波器； [n,wn]=cheb1ord(wp/pi,ws/pi,rp,rs); [bz1,az1]=cheby1(n,rp,wp/pi); %查看设计滤波器的曲线 [h,w]=freqz(bz1,az1,256,Fs); h=20*log10(abs(h)); figure;plot(w,h);title('所设计滤波器的通带曲线');grid on; % y=filter(bz1,az1,x);%对序列x滤波后得到的序列y end （2）高通滤波器 function y=highp(x,f1,f3,rp,rs,Fs) %高通滤波 %使用注意事项：通带或阻带的截止频率的选取范围是不能超过采样率的一半 %即，f1,f3的值都要小于 Fs/2 %x:需要带通滤波的序列 % f 1：通带截止频率 % f 2：阻带截止频率 %rp：边带区衰减DB数设置 %rs：截止区衰减DB数设置 %FS：序列x的采样频率 % rp=0.1;rs=30;%通带边衰减DB值和阻带边衰减DB值 % Fs=2000;%采样率 % wp=2*pi*f1/Fs; ws=2*pi*f3/Fs; % 设计切比雪夫滤波器； [n,wn]=cheb1ord(wp/pi,ws/pi,rp,rs); [bz1,az1]=cheby1(n,rp,wp/pi,'high'); %查看设计滤波器的曲线 [h,w]=freqz(bz1,az1,256,Fs); h=20*log10(abs(h)); figure;plot(w,h);title('所设计滤波器的通带曲线');grid on; y=filter(bz1,az1,x); end （3）带通滤波器 function y=bandp(x,f1,f3,fsl,fsh,rp,rs,Fs) %带通滤波 %使用注意事项：通带或阻带的截止频率与采样率的选取范围是不能超过采样率的一半 %即，f1,f3,fs1,fsh,的值小于 Fs/2 %x:需要带通滤波的序列 % f 1：通带左边界 % f 3：通带右边界 % fs1：衰减截止左边界 % fsh：衰变截止右边界 %rp：边带区衰减DB数设置 %rs：截止区衰减DB数设置 %FS：序列x的采样频率 % f1=300;f3=500;%通带截止频率上下限 % fsl=200;fsh=600;%阻带截止频率上下限 % rp=0.1;rs=30;%通带边衰减DB值和阻带边衰减DB值 % Fs=2000;%采样率 % wp1=2*pi*f1/Fs; wp3=2*pi*f3/Fs; wsl=2*pi*fsl/Fs; wsh=2*pi*fsh/Fs; wp=[wp1 wp3]; ws=[wsl wsh]; % % 设计切比雪夫滤波器； [n,wn]=cheb1ord(ws/pi,wp/pi,rp,rs); [bz1,az1]=cheby1(n,rp,wp/pi); %查看设计滤波器的曲线 [h,w]=freqz(bz1,az1,256,Fs); h=20*log10(abs(h)); figure;plot(w,h);title('所设计滤波器的通带曲线');grid on; y=filter(bz1,az1,x); end （4）带阻滤波器 function y=bands(x,f1,f3,fsl,fsh,rp,rs,Fs) %带阻滤波 %使用注意事项：通带或阻带的截止频率与采样率的选取范围是不能超过采样率的一半 %即，f1,f3,fs1,fsh,的值小于 Fs/2 %x:需要带通滤波的序列 % f 1：通带左边界 % f 3：通带右边界 % fs1：衰减截止左边界 % fsh：衰变截止右边界 %rp：边带区衰减DB数设置 %rs：截止区衰减DB数设置 %FS：序列x的采样频率 % f1=300;f3=500;%通带截止频率上下限 % fsl=200;fsh=600;%阻带截止频率上下限 % rp=0.1;rs=30;%通带边衰减DB值和阻带边衰减DB值 % Fs=2000;%采样率 % wp1=2*pi*f1/Fs; wp3=2*pi*f3/Fs; wsl=2*pi*fsl/Fs; wsh=2*pi*fsh/Fs; wp=[wp1 wp3]; ws=[wsl wsh]; % % 设计切比雪夫滤波器； [n,wn]=cheb1ord(ws/pi,wp/pi,rp,rs); [bz1,az1]=cheby1(n,rp,wp/pi,'stop'); %查看设计滤波器的曲线 [h,w]=freqz(bz1,az1,256,Fs); h=20*log10(abs(h)); figure;plot(w,h);title('所设计滤波器的通带曲线');grid on; y=filter(bz1,az1,x); end 4. 贝塞尔滤波器 % [b,a] = besself(n,Wn) %低通滤波器（Wn为二元向量时，为带通） % Wn是一个实际频率值，不是相对量 [b,a] = besself(n,Wn,’ftype’); % ftype类型low，hight，stop dataOut = filter(b,a,dataIn); %信号滤波运算 5. MATLAB中Filter type类型表示 'low''bandpass''high''stop' 6. 滤波器的离散化

传送门：点这里~ 在这里插入图片描述

参考文献：

https://blog.csdn.net/xiaokun19870825/article/details/78777325https://blog.csdn.net/cqfdcw/article/details/84939698https://www.mathworks.com/help/signal/ref/cheby1.htmlhttps://www.mathworks.com/help/signal/ref/cheby1.html#bucqk89-ftypehttps://www.mathworks.com/help/matlab/ref/filter.htmlhttp://blog.sina.com.cn/s/blog_49c02a8c0100yszh.htmlhttps://blog.csdn.net/zhoufan900428/article/details/8969470

【本文地址】

MATLAB辅助设计常用滤波器及其离散化

MATLAB辅助设计常用滤波器及其离散化

今日新闻

推荐新闻