FIR窗函数和IIR模拟、数字滤波器的MATLAB实现

导语

本文采用MATLAB R2014b对x.mat信号进行了时域频域的分析，同时创建窗函数、IIR模拟巴特沃兹、IIR数字滤波器对x.mat信号进行滤波处理。

原始信号的时频域分析

我是把x.mat信号放到了E盘目录下，其长度为100，进行时频域分析后，进行画图可以看到。当然不一定非要用x.mat，这个只是我用来测试的信号，大家自己选取就好，注意信号的长度。

load('E:\x.mat');
n=0:99;%长度
freqz(x);%求取频域
[X,w]=freqz(x);%对x信号进行频域分析
figure;
subplot(2,1,1); stem(n,x); grid;
subplot(2,1,2); plot(w/pi,abs(X));grid;%归一化处理

滤波器设计

上文我们已经分析出来他的幅频特性和相频特性，可以看出归一化后的0.2、0.5、0.9三个频点为信号的主要构成部分，因此下文进行三个滤波处理，分别滤出来三个频段信号

FIR窗函数低通滤波

这个没什么好说的，主要是对FIR窗的参数进行调用修改，可以修改的有wp1与ws1，分别对应通带与阻带频率，还有滤波器长度与矩形窗的时域响应是对应的，不同的窗有不同的特性。

%窗函数低通FIR滤波器设计  
wp1=0.3*pi;ws1=0.48*pi;   %性能指标
wdelta1=ws1-wp1;   %过渡带宽度
N1=ceil(0.92*pi/wdelta1);  %滤波器长度
Nw1=2*N1+1;     %窗口长度
wc1=(ws1+wp1)/2;  %截止频率
win1=boxcar(Nw1);%矩形窗的时域响应
b1=fir1(Nw1-1,wc1/pi,win1);%fir1是基于加窗的线性相位FIR数字滤波器设计函数
fx1=filter(b1,1,x);
figure;
subplot(2,1,1);stem(n,fx1);
fxd1=freqz(fx1);
subplot(2,1,2);plot(w/pi,abs(fxd1)); grid;
title('低通滤波器');

经过滤波后，我们再分析可以看到x信号经过低通出

如图为x信号经过LP后的时频域分析，低通滤波后对0.2频点的信号保留较好，0.5和0.9仍有残留，说明是滤波器的参数设定不够适量，但是已经满足我所需题目的要求。大家可以再去自己修改配置参数。
而滤波器本身的参数，请大家看IIR模拟的代码部分介绍。

IIR模拟带通滤波器

IIR模拟滤波器可以说，配置就相对FIR复杂很多了，我用其制作带通滤波器，设定通带、阻带频率，以及其他的频率变化法等计算需要一定的数字信号处理的知识支撑，如果不太懂的话，可以开放给大家修改的就是，四个上下限通带阻带频率与rp、rs。分别对应带内衰减，与带外衰减。

%带通滤波器   通带0.4-0.6  阻带0.3   0.7   αs=20    50->47
Fs3=1000;%给定抽样频率
wp_d=400;%通带下限截止频率
wp_s=600;%通带上限截止频率
ws_d=300;%阻带下限截止频率
ws_s=700;%阻带上限截止频率
wp3=[wp_d wp_s]/Fs3; 
ws3=[ws_d ws_s]/Fs3;%计算通带截止频率wp，阻带频率ws
rp3=3;rs3=20;%设置通带允许最大衰减设置为3dB，阻带应达到的最小衰减为20dB
[n3,wn3]=buttord(wp3,ws3,rp3,rs3);%计算巴特沃斯滤波器阶次和截止频率
%N为滤波器阶数 wn为滤波器截止频率 数字滤波器
[bz3,az3]=butter(n3,wn3,'bandpass');%频率变换法设计巴特沃斯带通滤波器
[X3,w3]=freqz(bz3,az3,256,Fs3/1500);
subplot(4,4,12); plot(w3/pi,angle(X3));
title('IIR带通滤波器相位响应');
subplot(4,4,16); plot(w3/pi,abs(X3));
title('IIR带通滤波器幅度响应');
fx3=filter(bz3,az3,x);
fxd3=freqz(fx3);
subplot(4,4,4);stem(n,fx3);
title('x信号经过IIR带通滤波器的时域分析');
subplot(4,4,8);plot(w/pi,abs(fxd3)); grid;
title('x信号经过IIR带通滤波器的频域分析');

这个地方就显示了滤波器本身的参数，在最后的六行，其先将fx计算出来，再进行freqz计算，画图。x信号经过带通出

IIR带通滤波器参数

IIR数字高通滤波器

%高通滤波器   通带0.75  阻带0.5  αs=20  
rp2=3;rs2=20;%设置通带允许最大衰减设置为3dB，阻带应达到的最小衰减为20dB
wp2=0.75*pi;%通带频率
ws2=0.5*pi;%阻带频率
FS2=1;    % let Fs=1 
wap2=tan(wp2/2); was2=tan(ws2/2);%预畸变公式
[n2,wn2]=buttord(wap2,was2,rp2,rs2,'s');%设计等效模拟滤波器
[z2,p2,k2]=buttap(n2);
[bp2,ap2]=zp2tf(z2,p2,k2);%设计模拟低通原型滤波器
[bs2,as2]=lp2hp(bp2,ap2,wap2);%将模拟低通原型滤波器转换为高通滤波器
[bz2,az2]=bilinear(bs2,as2,FS2/2);%实现双线性变换，即由模拟滤波器 得到数字滤波器
fx2=filter(bz2,az2,x);
figure;
subplot(2,1,1);stem(n,fx2);
fxd2=freqz(fx2);
subplot(2,1,2);plot(w/pi,abs(fxd2)); grid;
title('IIR数字高通滤波器');

可以看到IIR滤波器对理论知识分析要求较高，FIR只需要对不同窗函数的相关参数特性记住会修改分析就好了。x信号经过高通出

结论

看下图，在FIR与IIR的模拟与数字方面，我对IIR模拟进行了比较综合分析，大家如果想对其他两个滤波器详细分析，可以参考IIR模拟的滤波器。

代码中其他部分小的算法就不再阐述了，都有注释的
数字信号处理下的滤波器其实都是理论知识在里面的，作者做这个实验也是为了论证一下课本的一些题目以及设计综合实验而学习的，至于这些滤波器的其他LP、HP或者BP的代码，大家打开链接进行下载即可。PS：文件内的x.mat不要修改，load的位置注意看清楚。
百度网盘：滤波器代码链接
链接：https://pan.baidu.com/s/1crnYVl76MLxnQk3bJzCAug
提取码：zlet
里面可能有些注释或者什么的小问题可以在评论里询问。