Saturday 17 January 2015

Let's Do Some More Digital Signal Processing - Infinite Impulse Response (IIR) Filters


Following on from my earlier blog : Let's Do Some Digital Signal Processing (Part 2) - One Pole Filter I thought it was about time that I took a look at Infinite Impulse Response (IIR) filters.

IIR filters give much higher performance than Finite Impulse Response filters because they use feedback as well as feedforward in order to implement the filter response.

This code shows how to implement an IIR biquad filter structure.

#include <stdio.h>

#define COUNT 8  /* Number of samples to process */

#define FILTER_STAGES   1           /* Single 2nd ordrer biquad */

// Fs = 44.1 kHz, Fc = 1 kHz LPF
double LPFCoefficientArray [] =
{
3.43448050764769180000e-003, 6.86896101529538370000e-003, 3.43448050764769180000e-003,
-1.89857508921973460000e+000, 9.12313011250325380000e-001
};

float impulseSrc[] = {
 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0
 };
float stepSrc[] = {
 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0
 };
float Dst[] = {
 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0
 };

double FilterState [2 * FILTER_STAGES];

void iir_init (double *pState,
  const int NumberOfStages);

double iir_filter (const double Source,
  double *pState,
  const double *pCoeffs,
  const int NumberOfStages);

void main (void)

{
  int i;

  iir_init (FilterState, FILTER_STAGES);      /* Initialize the filter */
  for (i = 0; i < COUNT; i++)
    Dst[i] = iir_filter (impulseSrc[i], FilterState, LPFCoefficientArray, FILTER_STAGES);
  printf ("One Pole Impulse = %f, %f, %f, %f, %f, %f, %f, %f\n", Dst[0], Dst[1], Dst[2], Dst[3], Dst[4], Dst[5], Dst[6], Dst[7]);

  iir_init (FilterState, FILTER_STAGES);      /* Initialize the filter */
  for (i = 0; i < COUNT; i++)
    Dst[i] = iir_filter (stepSrc[i], FilterState, LPFCoefficientArray, FILTER_STAGES);
  printf ("One Pole Step = %f, %f, %f, %f, %f, %f, %f, %f\n", Dst[0], Dst[1], Dst[2], Dst[3], Dst[4], Dst[5], Dst[6], Dst[7]);

}


void iir_init (double *pState,
  const int NumberOfStages)

{
  int  i;

      /* Initialise the filter state array to 0 */
  for (i = 0; i < (NumberOfStages * 2); i++)
  {
    *pState++ = 0.0;
  }

}   /* End of iir_init () */


double iir_filter (const double Source,
  double *pState,
  const double *pCoeffs,
  const int NumberOfStages)

{
  register double       FeedbackSumOfProducts;
  register int i;
  register double       TempInputData;

  TempInputData = Source;

  for (i = 0; i < NumberOfStages; i++)
  {
    FeedbackSumOfProducts = TempInputData - (*(pCoeffs+3) * *pState) -
        (*(pCoeffs+4) * *(pState+1));       /* Feedback */
    TempInputData = (*pCoeffs * FeedbackSumOfProducts) +
        (*(pCoeffs + 1) * *pState) +
        (*(pCoeffs + 2) * *(pState+1));     /* Feedforward and save result for next time round */
    
    *(pState+1) = *pState;                  /* Move delayed samples */
    *pState = FeedbackSumOfProducts;

    pState += 2;        /* Increment array pointers */
    pCoeffs += 5;
  }

  return (TempInputData);                     /* Save output */

}   /* End of iir_filter () */


The result of running this code is :
One Pole Impulse = 0.003434, 0.013390, 0.025722, 0.036620, 0.046059, 0.054038, 0.060575, 0.065706
One Pole Step = 0.003434, 0.016824, 0.042546, 0.079167, 0.125226, 0.179264, 0.239839, 0.305545

The filter implemented in this example is a low-pass filter with 1 kHz cut-off frequency and a sampling rate of 44.1 kHz (CD frequency). This can be observed by modifying the above code to generate a longer impulse response (e.g. 1024) and then passing this data to a Fast Fourier Transform (FFT).

If you have found this solution useful then please do hit the Google (+1) button so that others may be able to find it as well.
Numerix-DSP Libraries : http://www.numerix-dsp.com/eval/

No comments:

Post a Comment