# Why does the filter function in matlab get different results from the IIR filters using direct form?

I have 6 first order filters, each one has 2 coefficients (B(x1,x1), B(x2,x2) for the numerator and A(y1,y1),A(y2,y2) for the denominator).

B = [0.6 0 ; 1 0.3 ; 1 -1; 1 -1 ; 1 -1; 1 -1]; % b(1,1) ajustado para 0.6
A = [1 -0.2025 ; 1 -0.2025 ; 1 -0.9860 ; 1 -0.9097 ; 1 -0.9973 ; 1 -0.9973];

[H1,f1] = freqz([B(1,1) B(1,2)] ,[A(1,1) A(1,2)],512,48000);
[H2,f2] = freqz([B(2,1) B(2,2)] ,[A(2,1) A(2,2)],512,48000);
[H3,f3] = freqz([B(3,1) B(3,2)] ,[A(3,1) A(3,2)],512,48000);
[H4,f4] = freqz([B(4,1) B(4,2)] ,[A(4,1) A(4,2)],512,48000);
[H5,f5] = freqz([B(5,1) B(5,2)] ,[A(5,1) A(5,2)],512,48000);
[H6,f6] = freqz([B(6,1) B(6,2)] ,[A(6,1) A(6,2)],512,48000);

mag = abs( H1.*H2.*H3.*H4.*H5.*H6); % pus um ajuste de 0.5

magdb = 20*log10(mag);
semilogx(f2,magdb)

ylim([-70 40])
title('Digital implementation of the A-weighting filter (fs = 48 kHz)')
xlabel('Frequency (Hz)')
ylabel('Gain(dB)')
grid on


This is the transfer function i obtain (ignore the title):

And here are all the poles (each pole corresponds to one first order filter):

P =

0.202500000000000

P =

0.202500000000000

P =

0.986000000000000

P =

0.909700000000000

P =

0.997300000000000

P =

0.997300000000000


Every pole is inside the unit circle.

This is the code and frequency response of the third first order filter:

%%% H3
figure()
mag = abs( H3);
magdb = 20*log(mag);
semilogx(f1,magdb);

ylim([-70 40])
title('Digital implementation of the A-weighting filter (fs = 48 kHz)')
xlabel('Frequency (Hz)')
ylabel('Gain(dB)')


Using the filter function i obtain the following for a 10000Hz sin wave:

t = linspace(0,1,fs);
x = sin(2*pi*10000*t);

y = filter([B(3,1) B(3,2)],[A(3,1) A(3,2)],x);

figure()
plot(t,x)
hold on
plot(t,y)
legend('Input Data','Filtered Data')


This is exactly what i'm expecting.

But now if i use the IIR Filter coeficcent using the past values i get something different:

x=[0 0];
y=[0 0];

t = linspace(0,1,fs);

yy = zeros(1,fs);

for c= 1:fs
x(1) = sin(2*pi*10000*t(c));

y(1) = (1/A(3,1))*(B(3,1)*x(1) + B(3,2)*x(2) + A(3,2)*y(2) );

yy(c) = y(1);
% update x and y data vectors
for i = 2:-1:1
x(i+1) = x(i); % store xi
y(i+1) = y(i); % store yi
end

end

figure()
plot(t,yy)


Why does this happen? How can i get the same response as the filter function?

I think the y(1) equation should be

y(1) = (1/A(3,1))*(B(3,2)*x(1) + B(3,1)*x(2) + A(3,1)*y(2) );


then I get

Note also, the update should just be

x(2) = x(1)
y(2) = y(1)


otherwise the x and y matrices grow to be $$3 \times 1$$ after the first time through the inner for loop.

• I agree with the inner loop part, but if i change the coefficients the way u did it i always get amplitude 1 even for low frequencies, that should not happen Commented Dec 24, 2022 at 21:04

Both your formula and the other answer formula for y(1) is wrong. Here is the correct formula:

y(1) = (1/A(3,1))*(B(3,1)*x(1) + B(3,2)*x(2) - A(3,2)*y(2));


Relative to the formula in your question, notice the subtraction rather than addition of the final term.

An alternative way to implement it that matches Matlab's filter function can be found in the algorithm section of the corresponding documentation.