/**************************************************************************
*
* ROUTINE
*               polefilt_3rd
*
* FUNCTION
*               Direct form all-pole filter
*
* SYNOPSIS
*               subroutine polefilt_3rd(a, n, z, xy, len)
*
*   formal
*
*                       data    I/O
*       name            type    type    function
*       -------------------------------------------------------------------
*       a(=3)           float   i       N+1 filter coefficients  
*       n               int     i       Filter order
*       z               float   i/o     N+1 filter delay elements
*                                       (maintained by the user)
*                                       (z[0] is a dummy delay)
*       xy              float   i/o     Input/Output data array
*       len             int     i       Number of samples to filter
*
***************************************************************************/

#include <string.h>
polefilt_3rd(a, n, z, xy, len)
int n, len;
register float a[], z[], xy[];
{
/*
  register int t, j;

  if (a[0] != 1.0)
  {
    printf("polefilt:  bad coefficients");
    exit(1);
  }
*/
  register float ar;
  size_t n_size = 16;
 

#define   STEP(K)          ar = xy[K]; \
                            ar -= a[1] * z[1];\
                            ar -= a[2] * z[2];\
                            ar -= a[3] * z[3];\
                             xy[K] = ar;\
                             z[0] = ar;\
memmove( (char *)(z+1), (char *)(z), n_size);

        STEP(0);  STEP(1);  STEP(2);  STEP(3);  STEP(4);  STEP(5);
        STEP(6);  STEP(7);  STEP(8);  STEP(9);  STEP(10); STEP(11);
        STEP(12); STEP(13); STEP(14); STEP(15); STEP(16); STEP(17);
        STEP(18); STEP(19); STEP(20); STEP(21); STEP(22); STEP(23);
        STEP(24); STEP(25); STEP(26); STEP(27); STEP(28); STEP(29);
        STEP(30); STEP(31); STEP(32); STEP(33); STEP(34); STEP(35);
        STEP(36); STEP(37); STEP(38); STEP(39); STEP(40); STEP(41);
        STEP(42); STEP(43); STEP(44); STEP(45); STEP(46); STEP(47);
        STEP(48); STEP(49); STEP(50); STEP(51); STEP(52); STEP(53);
        STEP(54); STEP(55); STEP(56); STEP(57); STEP(58); STEP(59);

/*
printf("len=%d  n=%d \n",len,n);
  for (t = 0; t < len; t++)
  {
    z[0] = xy[t];
    for (j = n; j > 0; j--) z[0] -= a[j] * z[j];
    memmove( (char *)(z+1), (char *)(z), n*sizeof(float) );
    xy[t] = z[0];
  }
*/
}