/* *  -*- Mode ANSI C -*- *  $Id: filter.c,v 1.4 1996/08/16 17:25:47 fernande Exp $ *  $Source: /sgi.acct/sweldens/cvs/liftpack/Lifting/filter.c,v $ *  Author: Gabriel Fernandez * *  This module finds the filter coefficients used in the wavelet *  transform to generate the wavelet coefficients (gamma) using *  the contribution of the lambda values. The process uses the *  Neville's algorithm (neville.c) for the polynomial interpolation *  and a given value, N, to know how many filter coefficients are *  needed. The results are returned in a Matrix of N/2+1 rows by *  N columns in a linear representation. *//* do not edit anyhting above this line *//* System header files */#include <assert.h>/* FLWT header files */#include <flwtdef.h>   /* Matrix type and error codes definitions */#include <filter.h>    /* external declaration of GetFilter() */#include <mem.h>#include <neville.h>   /* Neville() function */#include <util.h>      /* ODD() macro *//* code *//* * GetFilter function: finds the filter coefficients used in the *                     Gamma coefficients prediction routine. The *                     Neville polynomial interpolation algorithm *                     is used to find all the possible values for *                     the filter coefficients. Thanks to this *                     process, the boundaries are correctly treated *                     without including artifacts. *                     Results are ordered in matrix as follows: *                         0 in the left and N   in the right *                         1 in the left and N-1 in the right *                         2 in the left and N-2 in the right *                                        . *                                        . *                                        . *                         N/2 in the left and N/2 in the right *                     For symmetry, the cases from *                         N/2+1 in the left and N/2-1 *                     to *                         N in the left and 0 in the right *                     are the same as the ones shown above, but with *                     switched sign. */VectorGetFilter ( const int N ){    static Flt xa[40], ya[40], x;    static int row, col, cc, Nrows;    Vector filter, ptr;    assert ( N < sizeof xa / sizeof *xa );   /* N is greater than dim(xa) */    /* Check for validity of N */    if ( ODD (N) || (N == 0) )        return (Vector)NULL;   /* N is not even and greater than zero */    /* Number of cases for filter calculations */    Nrows = (N>>1) + 1;    /* N/2 + 1 */    /* Allocate memory for filter matrix */    filter = vector( 0, (long)(Nrows*N - 1) );    ptr = filter;    /* Generate values of xa */    xa[0] = (Flt)0.5*(Flt)(1 - N);    /* -N/2 + 0.5 */    for (col = 1 ; col < N ; col++)        xa[col] = xa[col-1] + (Flt)1;    /* Find filter coefficient values */    filter += ( (Nrows*N) - 1 );   /* go to last position in filter matrix */    for (row = 0 ; row < Nrows ; row++) {        x = (Flt)row;        for (col = 0 ; col < N ; col++) {            for (cc = 0 ; cc < N ; cc++)                ya[cc] = (Flt)0;            ya[col] = (Flt)1;            *(filter--) = Neville (xa, ya, N, x);        }    }    return ptr;}