2
$\begingroup$

DSP newcomer here!

I am tinkering around with a TI DSP and am trying to implement a second order IIR filter in C. Input is 16 bit 2's complement, as is the output, the accumulator is 32 bit wide. I have tried several structures: Direct form I, Transposed direct form II, Gold & Rader structure. However, I can't seem to get anything to work and I am not sure whether it is a typo in the code of all of my filter implementations or whether it is because of me misunderstanding the principle of fixed-point arithmetic. I assume - or better: hope - that it's the same mistake that is affecting all implementations.

Here's my code for the transposed df II implementation:

// state variables
int w[3] = {0};

short filter_dfii(short value) {
    w[0] = ((value * sos_filter_coeffs_b[0] + w[1]) >> 15 ) & 0xffff;
    w[1] = value * sos_filter_coeffs_b[1]
       - w[0] * sos_filter_coeffs_a[1]
       + w[2];
    w[2] = value * sos_filter_coeffs_b[2]
       - w[0] * sos_filter_coeffs_a[2];

    return (short) w[0];
}

My df I implementation looks like this:

int last_value[2] = {0};

short filter_dfi(short value) {
    w[0] =
        ((value * sos_filter_coeffs_b[0]
        + last_value[0] * sos_filter_coeffs_b[1]
        + last_value[1] * sos_filter_coeffs_b[2]
        - w[1] * sos_filter_coeffs_a[1]
        - w[2] * sos_filter_coeffs_a[2])
        >> 15) & 0xffff;
    last_value[1] = last_value[0];
    last_value[0] = value;
    w[2] = w[1];
    w[1] = w[0];
    return (short) w[0];
}

The filter coefficients are:

int sos_filter_coefficients_a = {32767, -51150, 21015};
int sos_filter_coefficients_b = {658, 1316, 658};

How did I get those coefficients? I utilized the Matlab 'butter' function and specified a 2nd order Butterworth lowpass with a cutoff frequency of 0.05*fs or 2400 kHz in my specific case. Then I multiplied each coefficient with 32767 and rounded it to an integer value. If I call Matlab freqz function with those quantized coefficients, I get a bode plot with the expected result, so I believe quantization error is not the issue here.

However, the result is not a lowpass but just garbage. Unfortunately, I lack the proper testing equipment to describe the result scientifically, but it sounds like digital oscillation that is slightly modulated by the actual input signal and is present even when there's no input at all.

So, in my opinion, this can't be an unstable filter since the poles are well within the unit circle. It can't be overflow limit cycle either because then it shouldn't be happening when there's only low input level or no input at all.. right?

If I alter the coefficients to

int sos_filter_coefficients_a = {32767, 0, 0};
int sos_filter_coefficients_b = {32767, 0, 0};

input gets returned basically unprocessed as expected, so my issues most likely aren't caused by any external factors - there must be something wrong with my coefficients, with the function itself, or both.

Any help is much appreciated!

$\endgroup$
  • $\begingroup$ why is w[current_filter][0] a two-dimensional array but w[1] and w[2] are one-dimensional? and you are returning short w[0]. is w[ ] one-dimensional or two-dimensional? i can't see how that code compiles. $\endgroup$ – robert bristow-johnson May 15 '17 at 2:14
  • $\begingroup$ Thank you for pointing that out, Robert! My original code is designed for handling multiple filters and channels, but when writing this question I left that part part out for simplicity's sake. Seems like I was sloppy with the editing. $\endgroup$ – UnbescholtenerBuerger May 15 '17 at 7:42
  • $\begingroup$ why isn't there something like w[2] = w[1]; w[1] = w[0]; in your DF2 code? $\endgroup$ – robert bristow-johnson May 15 '17 at 16:50
  • $\begingroup$ It's actually transposed df2. $\endgroup$ – UnbescholtenerBuerger May 15 '17 at 17:54
2
$\begingroup$

The " & 0xffff" operation is questionable. It works for positive numbers but for negative numbers it will mask out the sign bit as well, so you will turn everything into positive numbers.

For the DF1, it would be better to delay the shift and mask until you return the actual output and do the conversion to short. This way you can keep the state variable in full double precision.

$\endgroup$
  • $\begingroup$ Removing the & 0xffff statement indeed helped a bit. Now audio gets passed through, although the actual filter doesn't behave as specified at all, so I need to do further investigation. Regarding your DF1 suggestion - I don't see how that would be beneficial. Before multiplying the variables with the coefficients, I would need to right-shift them anyways to avoid overflow, wouldn't I? $\endgroup$ – UnbescholtenerBuerger May 15 '17 at 13:37
1
$\begingroup$

For starters, what size is "int" ? 16 bits or 32 bits? One of your coefficient is "-51150" which requires 17 bits to represent, hence the need for a 32-bit integer for coefficients.

Make sure "int" is indeed a 32-bit integer, otherwise it would not work. Or, you can maybe multiply coefficients by 16384 instead of 32768 in order to make them fit in 16 bits.

$\endgroup$
  • $\begingroup$ Yes, my int definitely holds 32 bit. Nevertheless , scaling my coefficients to 2^14 to have them fit into short sounds like an interesting idea. However, my a0 coefficient wouldn't be 1.0 (or 32767 respectively) any more. Can I just right-shift my output by 14 instead of 15 bits and everything's just fine or do I need to make any additional adjustments? $\endgroup$ – UnbescholtenerBuerger May 15 '17 at 13:31
  • $\begingroup$ Since, a[0] = 1, you scale it the way you want. You can analyze the difference in Matlab or Octave when you scale by 16384 or 32768. I'm a bit puzzled as to why you have explicitely written your a[0] coefficient. You don't use it in your code... $\endgroup$ – Ben May 15 '17 at 21:30

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.