convfilter.cpp

Adaptive digital Filters in C++

// convfilter.c++	a general purpose convolution, filter

//                      (c) Copyright 1995, Everett F. Carter Jr.
//                      Permission is granted by the author to use
//			this software for any application provided this
//			copyright notice is preserved.

        /* data is organized as follows:                           *
         *                                                         *
	 *  x[2]     --->  x[n]         latest value               *
         *  x[1]     --->  x[n-1]       previous value             *
         *  x[0]     --->  x[n-2]       value before x[n-1]        *
         *                                                         *
	 *  y[] is organized in the same way                       *
         *                                                         *
	 */


static const char rcsid[] = "@(#)convfilter.c++	1.6 16:36:05 5/16/95   EFC";

#include <convfilter.hpp>

static const complex zero = 0.0;
static const complex  one = 1.0;

void ConvolutionFilter::free_x_space()
{
	if ( nx > 0 )
	{
		delete []x;
		delete []alpha;
		err--;
	}
}

void ConvolutionFilter::free_y_space()
{

	if ( ny > 0 )
	{
		delete []y;
		delete []beta;
		err--;
	}
}

void ConvolutionFilter::allocate_x_space(const int wx)
{

	if ( nx != wx )
	{
		free_x_space();

		nx = wx;

		x     = new float[nx];
		alpha = new float[nx];
	}

	if (x == NULL || alpha == NULL )
				err--;
	else
				err++;
	nm1 = nx - 1;

	invals = 0;

	
}

void ConvolutionFilter::allocate_y_space(const int wy)
{

	if ( ny != wy )
	{
		free_y_space();

		if ( wy <= 0 )
			ny = nx - 1;
		else
			ny = wy;

		if ( ny < 1 )
			return;
	
		y     = new float[ny];
		beta  = new float[ny + 1];

	}

	if (y == NULL || beta == NULL )
				err--;
	else
				err++;

	mm1 = ny - 1;
	invals = 0;

	
}

void ConvolutionFilter::shift(float *vals, const int last)
{
          /* shuffle the data down one */
          for (int k = 0; k < last; k++)
			 vals[k] = vals[k + 1];

}

float ConvolutionFilter::convolve()
{
	 int k;

	 for (k = 0, out_val = 0.0; k < nx; k++)
			 out_val += x[nm1 - k] * alpha[k];

	 for (k = 0; k < ny; k++)
			 out_val -= y[mm1 - k] * beta[k+1];

	return out_val;

}

complex ConvolutionFilter::h(const float z,const int inv)
{
	int i;
	complex num, denom, transfer;

	for (i = 0, num = zero; i < nx; i++)
	{
		num = num + complex( alpha[i] * cos( i * z ) ,
				 inv * alpha[i] * sin( i * z ) );
	}

	if ( ny < 1 )
		denom = one;
	else
	for (i = 0, denom = zero; i <= ny; i++)
	{
		denom = denom + complex( beta[i] * cos( i * z ) ,
				   inv * beta[i] * sin( i * z ) );
	}

	transfer = num / denom;

	return ( transfer );

}

ostream& ConvolutionFilter::coefs(ostream& os) const
{
	int k;

	if ( alpha )
	{
		for (k = 0; k < nx; k++)
			os << alpha[k] << "  ";
		os << '\n';
	}

	if ( beta )
	{
		for (k = 0; k < ny; k++)
			os << beta[k+1] << "  ";
		os << '\n';
	}


	return os;
}