matlab手动实现IIR滤波器

matlab⼿动实现IIR滤波器
使⽤⼯具为matlab7.0
⼀，利⽤matlab fadtool⼯具获取SOS matrix 和 factor scales
1.在命令框中输⼊fdatool
2.配置滤波器信息，配置结果如下图
3.将SOS matrix 和 factor scales保存到workspace,保存的变量名为SOS和scale，并且⽣成my_filter.m⽂件，my_filter.m的⽂件内容如下function Hd = my_filter
N = 3; % Order
Fc = 150; % Cutoff Frequency
h = fdesign.lowpass('N,Fc', N, Fc, Fs);
Hd = butter(h);
⼆.将SOS matrix 和 factor scales转化为A,B系数
1.在matlab命令对话框输⼊命令
g = prod(scale)
[b,a] = sos2tf(SOS,g)
三.⼿动编写代码实现IIR滤波器
IIR差分⽅程如下：
function F = filter()
% 采样频率
fs = 1000;
x = 0:1/fs:20;
% 构造的信号
sing = 0.5sin(2pi100x)+0.8cos(2pi180x);
% 将数据进⾏fft变换
function my_plot(data)
N = length(data);
xdft = fft(data);
xdft = xdft(1:N/2+1);
% psdx = abs(xdft).^2;
freq = 0:fs/N:fs/2;
plot(freq, (2/N)abs(xdft))
end
subplot(311)
my_plot(sing)
% A,B系数
b = [0.04953,0.14860,0.14860,0.04953];
a = [1,-1.16192,0.69594,-0.13776];
len = length(sing);
% 前⼏个时刻的输⼊值
pre_x = [0, 0, 0, 0];
% 前⼏个时刻的输出值
pre_y = [0, 0, 0, 0];
out = zeros(1,len);
for i=1:len
pre_x(1) = sing(i);
% 差分⽅程
out(i) = b(1)pre_x(1)+b(2)pre_x(2)+b(3)pre_x(3)+b(4)pre_x(4)-a(2)pre_y(2)-a(3)pre_y(3)-a(4)pre_y(4);
pre_y(1) = out(i);
for j = 4 : -1 : 2
pre_x(j) = pre_x(j-1);
pre_y(j) = pre_y(j-1);
end
end
subplot(312)
my_plot(out)
% 运⾏直接保存的.m⽂件
Hd = my_filter;
output = filter(Hd, sing);
subplot(313)
my_plot(output)
end
结果如下：。