IIR filter implementation in Direct Form 2

Question

I have implemented Direct Form 2 IIR filter. The input is the Kronecker delta function. I have written code for the response of the filter for Kronecker delta input. The code is:

clc;   
clear all;
close all;
x=[1,zeros(1,150)];%input sequence

y=[];%output sequence

b=[1,-1.5511722144889886,1];%numerator coefficients

a=[1,-1.7473520798913149,0.79762351115754226];%denominator coefficients

bs=0.099825504845559299*b; 
%...floating point.....%
w=[];

w(1)=x(1);
 
w(2)=-a(2)*w(1)+x(2);

y(1)=bs(1)*w(1);
 
y(2)=bs(1)*w(2)+bs(2)*w(1);
 for n=3:1:150
     
w(n)=x(n)-a(2)*w(n-1)-a(3)*w(n-2);
    
y(n)=bs(1)*w(n)+bs(2)*w(n-1)+bs(3)*w(n-2);
 
end

%freqz(y);

% 

% 

%.....fixed point conversion....%

%.....considering precision of 16 i.e q16 format....%

for n=1:3

    a_f(n)=int64(a(n)*pow2(16));%.....................converted into fixed......%

end

 for n=1:3

     b_f(n)=int64(bs(n)*pow2(16));

 end

for n=1:3

    a_fix2float(n)=double(a_f(n))*pow2(-16);

end

for n=1:3

    b_fix2float(n)=double(b_f(n))*pow2(-16);

end

w_f(1)=(x(1)*pow2(16));

w_f(2)=((-a_f(2)*w_f(1))*pow2(-16));

for n=3:150

    w_f(n)=((-a_f(2)*w_f(n-1))-(a_f(3)*w_f(n-2)))*pow2(-16)+x(n);

end

y_f(1)=(b_f(1)*w_f(1))*pow2(-16);

y_f(2)=((b_f(1)*w_f(2))+(b_f(2)*w_f(1)))*pow2(-16);

for n=3:150

 y_f(n)=((b_f(1)*w_f(n))+(b_f(2)*w_f(n-1))+(b_f(3)*w_f(n-2)))*pow2(-16);

 end

...converting fixed to float...%

% format long

y_fix2float(1)=double(y_f(1))*pow2(-16);

y_fix2float(2)=double(y_f(2))*pow2(-16);

for n=3:150

    y_fix2float(n)=double(y_f(n))*pow2(-16);%...convert to float..%

end

 for n=1:150

     w_fix2float(n)=double(w_f(n))*pow2(-16);

 end

   %to manually plot frequency response..%

  HH = abs(fft(y));

 yy=abs(fft(y_fix2float));

 plot(HH(1:50),'.r');

 hold on

 plot(yy(1:50),'g');

 xlabel('time');

ylabel('response');

legend('floating','fixed');

title('directform2 q16');

I have written MATLAB code in both floating and fixed point. In fixed point I considered precision of q16.16. My response is:

My question is: Is the response I got correct?

Answer 1

It seems to be wrong. The first sample of the impulse response is always equal the b0 coefficient (b(1) in Matlab, so in this case that should be 0.099825504845559299 which doesn't seem to be the case.

The correct answer can be checked with yref = filter(bs,a,x); and it's substantially different from your answer. You have shown that the fixed point version matches the floating point version, however the floating point version is wrong to start with.

Finally, Direct Form II is a horrible choice for fixed point filters. The transfer function between input and state variables is given by the pole-only transfer function. Even in this harmless example this is about 27 dB of gain. It's almost impossible to prevent clipping in Direct Form II. Use Direct Form I or Transposed Form II instead

Answer 2

Your frequency response plot is correct. Note that you labeled the x-axis as 'time', which created some confusion, because it should be 'frequency'. This is why Hilmar thought it's wrong (because he thought you plotted the impulse response).

But let me explain what I think you did, because this might be different from what you intended to do. First of all, this is no fixed-point implementation of the filter because all computations are done in floating point, and no computation results ever get quantized. The only thing you quantized are the filter coefficients. And here you didn't use 16 bits to quantize the coefficients (as you may think) but actually 18 bits, because you use 16 fractional bits, and - implicitly - as many bits before the binary point as necessary to represent the filter coefficients. And you need two non-fractional bits to represent the coefficient with the largest magnitude $a_1=-1.7474$.

So after all it is not surprising that there is no visible difference between the FFTs of the two responses. A well-behaved second order filter (i.e. with poles not close to the unit circle and not close to the real axis) with 18-bit quantized coefficients has almost the same frequency response as the corresponding system with floating point coefficients. Note that this may change drastically as soon as you implement signal quantization.