libbpm: Digital Signal Processing Routines

This module contains the definitions for the digital signal processing routines for libbpm.

    filter_t *filter = create_filter( "the_filter", RESONATOR | , 0,
                                       nsamples, 40.*kHz, 8.*kHz, 0., 200. );

The arguments the filter expects is a name for the filter (just for esthetic purposes when printing the filter), the filter options, which are explained below, the order of the filter, where it is meaning full (e.g. Butterworth, Bessel, Chebyshev). Then it needs the number of samples in the waveforms which will be filtered by this filter, the sampling frequency and one (optionally two) frequency parameter. For lowpass/highpass filters and the resonater, only the first frequency defines respectively the -3dB frequency level for the low/high pass and the resonance frequency for the resonator (the witdh is defined by the Q value in this case). For bandpass/stop filters the two frequencies are required and define the -3dB level which defines the bandwidth of the filter, with f1 being the lower end frequency and f2 the higher end.

The IIR Bessel, Butterworth and Chebyshev filters can be normalised as lowpass (option LOWPASS) which is the default, highpass (option HIGHPASS), bandstop (option BANDSTOP) or bandpass (option BANDPASS) filters. They are designed with poles and zeros in the s plane that are transformed to the z plane either by bilinear z transform (option BILINEAR_Z_TRANSFORM) or matched z transform (option MATCHED_Z_TRANSFORM). Just "OR" the options together to setup the filter, e.g. :

    filter_t *filter = create_filter( "lp", BESSEL | HIGHPASS | MATCHED_Z_TRANSFORM, 0,
                                       ns, 40.*kHz, 8.*kHz, 0., 200. );

The resonators are designed directly with their 2 poles and 2 zeros in the z plane and can be normalised either as BANDPASS (default), BANDSTOP (or NOTCH) or ALLPASS resonators.

The last argument to the create_filter() routine is a parameter which can optionally be given to the filter. It depends on the filter chosen, currently the parameter has meaning for the following filters :

The filter coefficients for the difference equation are calculated and checked for consistency, upon which they are stored in the filter structure. Once this is done and the filter is setup, application to various waveforms is fairly straightforward. Note that you only have to define your filter once during initialisation. Once setup, it can be used to filter any number of waveforms of the same type.

    apply_filter( filter, wf );

To get an impulse response from the filter into the secified waveform, where the impulse is given at sample 1000, the following routine is implemented.

    filter_impulse_response( filter, wf, 1000 );

This routine creates an impulse function (zero everywhere, except at the sample you enter, where it's value is 1) and puts it through the filter. The FFT of this impulse response gives you the filter characteristic in frequency domain. Also you can check the filter's response to a step function, it's so-called step response :

    filter_step_response( filter, wf, 1000 )

The step response is defined as the response of the filter to an input function which is zero at the beginning and 1 for samples >= the sample you specify.

The Bessel, Butterworth and Chebyshev filters

The Resonator filter

The gaussian filter

The Digital Downconversion Algorithm (DDC)

    int ddc( doublewf_t *w, double f, filter_t *filter, complexwf_t *dcw );

returns then the complex DC waveform (dcw), where it's amplitude and phase can then be used in further calculations for beam position and slope in the BPM. We recommend the usage of a GAUSSIAN low-pass filter for the double frequency filtering as this shows the best phase behaviour combined with linearity (see create_filter()).

For fast execution, the DDC routine comes with a buffer which it only allocates once by doing

This buffer is used in the filtering routine, you can clean up after the execution of the buffer by having

Discrete (Fast) Fourier Transforms

    fft_initialise( int ns )

intialises the buffers for waveforms of a certain sample length ns. Note that ns has to be a power of 2. You can clear the FFT buffers by issuing

Then two wrapper routines are implemented which take doublewf_t and complexwf_t data.

Complex Discrete Fourier Transform

    int complexfft( complexwf_t *z, int fft_mode );

which takes a complex waveform and performs an FFT in place. The fft_mode argument can be either

Note the backward and forward FFT's have a factor of n inbetween them, so to get the orginal wf back after applying both the backward and the forward FFT, you need to divdide by the number of samples z->n.

Real Discrete Fourier Transform

    int realfft( doublewf_t *y, int fft_mode, complexwf_t *z );

and FFT_BACKWARD takes the input frmo the first half (n/2) of the complexwf_t and FFTs it, expanding to a doublewf_t of length n.

$X[k] = \frac{(Re(x[0]) + Re(x[n/2])*\cos(\pi*k) )}{2} + \sum_{j=1}^{n/2-1} Re(x[j])*\cos(2*\pi*j*k/n) + \sum_{j=1}^{n/2-1} Im(x[j])*\sin(2*\pi*j*k/n), 0<=k<n$

Reference for FFT routines

DSP example program

    #include <stdio.h>
    #include <stdlib.h>
    #include <string.h>
    #include <math.h>

    #include <iostream>
    
    #include <TROOT.h>
    #include <TFile.h>
    #include <TTree.h>
    
    #include <bpm/bpm_process.h>
    #include <bpm/bpm_units.h>
    #include <bpm/bpm_simulation.h>
    #include <bpm/bpm_nr.h>
    #include <bpm/bpm_rf.h>
    #include <bpm/bpm_alloc.h>
    #include <bpm/bpm_dsp.h>
    #include <bpm/bpm_wf.h>
    
    using namespace std;

    int main( int argc, char **argv ) {
    
      cout << "Welcome to the libbpm DSP sandbox" << endl;
    
      int ns    = 256;
      double fs = 119.*MHz;
      
      doublewf_t *w = doublewf( ns, fs );
      doublewf_t *s = doublewf_sample_series( ns, fs );
      
      doublewf_t *ddc_amp   = doublewf( ns, fs );
      doublewf_t *ddc_phase = doublewf( ns, fs );
      
      // setup the root trees...
      TFile *rootfile = new TFile( "dsp.root", "recreate" );
      TTree *roottree = new TTree( "dsp", "libbpm dsp tests" );
      
      int evt;
      double amp, phase;
      double gen_amp, gen_phase;
      
      // setup the branches in the tree
      roottree->Branch( "evt",       &evt,          "evt/I"            );
      roottree->Branch( "wf",        w->wf,         "wf[256]/D"        );
      roottree->Branch( "s",         s->wf,         "s[256]/D"         );
      roottree->Branch( "gen_amp",   &gen_amp,      "gen_amp/D"        );
      roottree->Branch( "gen_phase", &gen_phase,    "gen_phase/D"      );
      roottree->Branch( "ddc_amp",   ddc_amp->wf,   "ddc_amp[256]/D"   );
      roottree->Branch( "ddc_phase", ddc_phase->wf, "ddc_phase[256]/D" );
      
      complexwf_t *ddcwf = complexwf( ns, fs ); 

      filter_t *gauss  = create_filter( "gauss", GAUSSIAN,0,ns,fs,6.*MHz,0.,0.001);
      filter_t *butter = create_filter( "butter", BUTTERWORTH | LOWPASS,4,ns,fs,6.*MHz,0.,0.);
      filter_t *bessel = create_filter( "bessel", BESSEL | LOWPASS,4,ns, fs,6.*MHz,0., 0.);
      filter_t *cheby  = create_filter( "cheby",  CHEBYSHEV | LOWPASS,4,ns,fs,6.*MHz,0.,-10.);
      
      // init the DDC
      ddc_initialise( ns, fs );
      
      for ( evt = 1; evt<=1000; evt++ ) {
      
        // Make the waveform
        gen_amp   = (double) evt * 10.;
        gen_phase = PI / (double) evt;
        
        // reset the w to 0... quite important :D
        doublewf_reset( w );
        
        doublewf_add_dcywave( w, gen_amp, gen_phase, 21.4*MHz, 0.15*usec, 0.2*usec, 0. );
        
        // do the DDC :)
        if ( ddc( w, 21.4*MHz, gauss, ddcwf ) ) return 1;
        
        // want to try differen filters ?
        //if ( ddc( w, 21.4*MHz, butter, ddcwf ) ) return 1;
        //if ( ddc( w, 21.4*MHz, bessel, ddcwf ) ) return 1;
        //if ( ddc( w, 21.4*MHz, cheby, ddcwf ) ) return 1;
        
        
        // get amplitude and phase from complex wf
        complexwf_getamp( ddc_amp, ddcwf );
        complexwf_getphase( ddc_phase, ddcwf );
        
        // fill the tree...
        roottree->Fill();
        
        if ( evt % 100 == 0 ) cout << "Simulated " << evt << " events." << endl;
      }    

      // clear the DDC memory buffers
      ddc_cleanup();
      
      rootfile->Write();
      rootfile->Close();
      
      delete_filter( gauss );
      delete_filter( butter );
      delete_filter( bessel);
      delete_filter( cheby );
      
      complexwf_delete( ddcwf );
      
      doublewf_delete( w );
      doublewf_delete( s );
      doublewf_delete( ddc_amp );
      doublewf_delete( ddc_phase );
      
      return 0;
    }

Define Documentation

Function Documentation

EXTERN filter_t* create_filter	(	char	name[],
		unsigned int	options,
		int	order,
		int	ns,
		double	fs,
		double	f1,
		double	f2,
		double	par
	)

Creates the filter.

Parameters:

	name	a name for the filter
	options	filter specification and options bitword
	order	filter order
	ns	number of samples of the waveforms
	fs	sampling frequency
	f1	first frequency
	f2	optional second frequency ( bandpass/bandstop )
	par	optional parameter for chebyshev : ripple in dB for resonator : Q factor

Returns:: A pointer to the created filter structure, memory is allocated on the heap inside this routine, the user has to take of deleting it using delete_filter().

Definition at line 10 of file create_filter.c.

References alloc_simple_wave_double(), filter_t::alpha1, filter_t::alpha2, BESSEL, bpm_error(), bpm_warning(), BUTTERWORTH, calculate_filter_coefficients(), CAUSAL, filter_t::cheb_ripple, CHEBYSHEV, filter_t::cplane, create_resonator_representation(), create_splane_representation(), filter_t::f1, filter_t::f2, FIR, filter_t::fs, filter_t::gauss_cutoff, GAUSSIAN, gaussian_filter_coeffs(), IIR, filter_t::name, NO_PREWARP, NONCAUSAL, normalise_filter(), filterrep_t::npoles, filter_t::ns, filter_t::options, filter_t::order, filter_t::Q, RESONATOR, filter_t::w_alpha1, filter_t::w_alpha2, filter_t::wfbuffer, filter_t::yc, and zplane_transform().

Referenced by add_mode_response(), and get_mode_response().

EXTERN int apply_filter	(	filter_t *	f,
		double *	wf
	)

Apply the filter to the given waveform. Note that the filter is applied in place, the user has to make a copy of the waveform if he/she wants to keep the original before applying the filter. The number of samples in the waveform has to be set in advance when creating the filter, it is stored in the filter structure (f->ns).

Parameters:

	f	pointer to a filter that was created using create_filter
	wf	an array containing the waveform to be filtered

Returns:: BPM_SUCCESS upon success and BPM_FAILURE upon failure

Definition at line 19 of file apply_filter.c.

References ANTICAUSAL, bpm_error(), CAUSAL, FIR, filter_t::gain, filter_t::ns, filter_t::nxc, filter_t::nxc_ac, filter_t::options, filter_t::wfbuffer, filter_t::xc, filter_t::xc_ac, filter_t::xv, filter_t::xv_ac, filter_t::yv, and filter_t::yv_ac.

Referenced by add_mode_response(), ddc(), filter_impulse_response(), filter_step_response(), and get_mode_response().

EXTERN void print_filter	(	FILE *	of,
		filter_t *	f
	)

Prints the filter to the given file pointer.

Parameters:

	of	the filepointer, use "stdout" to print to the terminal
	f	the filter to be printed

Returns:: void

Definition at line 8 of file print_filter.c.

References bpm_error(), c_abs(), c_arg(), CAUSAL, filter_t::cplane, filter_t::dc_gain, filter_t::fc_gain, FILT_EPS, filter_t::gain, filter_t::hf_gain, filter_t::name, filter_t::nxc, filter_t::options, print_filter_representation(), and filter_t::xc.

EXTERN int filter_step_response	(	filter_t *	f,
		double *	wf,
		int	itrig
	)

This routine fills the given wf with the step response of the filter. The step response is defined as wf[i] = 0. for i < itrig and wf[i] = 1. for i >= itrig.

Parameters:

	f	a pointer to the filter to use
	wf	pointer to a waveform which will be overwritten with the step response
	itrig	the sample number in the waveform which will have the step

Returns:: BPM_SUCCESS upon succes and BPM_FAILURE upon failure

Produces a stepresponse for the filter, step is defined by the trigger sample number the starting level and the endlevel

Definition at line 8 of file filter_step_response.c.

References apply_filter(), bpm_error(), and filter_t::ns.

EXTERN int filter_impulse_response	(	filter_t *	f,
		double *	wf,
		int	itrig
	)

This routine fills the given wf with the impulse response of the filter. The impulse response is defined as wf[i] = 1. for i == itrig and wf[i] = 0. elsewhere.

Parameters:

	f	a pointer to the filter to use
	wf	pointer to a waveform which will be overwritten with the impulse response
	itrig	the sample number in the waveform which will have the impulse

Returns:: BPM_SUCCESS upon succes and BPM_FAILURE upon failure

Produces an impulse response for the filter, step is defined by the trigger sample number the starting level and the endlevel

Definition at line 7 of file filter_impulse_response.c.

References apply_filter(), bpm_error(), and filter_t::ns.

EXTERN filterrep_t* create_splane_representation ( filter_t * f )

This routine returns a pointer to a filter representation filterrep_t in the s plane for Butterworth, Chebyshev and Bessel filters. It need an initialised filter structure which has the filter type and the order set. Memory is allocated for this routine on the heap, so the user is responsible to delete this memory using free().

Parameters:

the initialised filter with the correct options in f->options

Returns:: the filter representation in the s plane

Definition at line 32 of file create_splane_representation.c.

References _add_splane_pole(), BESSEL, bpm_error(), BUTTERWORTH, c_conj(), c_exp(), filter_t::cheb_ripple, CHEBYSHEV, complex(), filterrep_t::npoles, filter_t::options, and filter_t::order.

Referenced by create_filter().

EXTERN filterrep_t* create_resonator_representation ( filter_t * f )

This routine returns a pointer to a filter representation filterrep_t in the z plane for resonance filters. It needs an initialised filter structure which has the filter type and the Q factor set. Memory is allocated for this routine on the heap, so the user is responsible to delete this memory using free().

Parameters:

the initialised filter with the correct options in f->options

Returns:: the filter representation in the z plane

Definition at line 15 of file create_resonator_representation.c.

References _eval_complex_polynomial(), _expand_complex_polynomial(), _reflect(), ALLPASS, filter_t::alpha1, BANDSTOP, bpm_error(), c_conj(), c_div(), c_exp(), complex(), FILT_EPS, complex_t::im, MAX_RESONATOR_ITER, MAXPZ, filterrep_t::npoles, filterrep_t::nzeros, filter_t::options, filterrep_t::pole, filter_t::Q, complex_t::re, and filterrep_t::zero.

Referenced by create_filter().

EXTERN filterrep_t* zplane_transform	(	filter_t *	f,
		filterrep_t *	s
	)

This routine transforms the poles and zeros for Bessel, Chebyshev and Butterworth filters to the z plane either via matched z transform or bilinear z transform. This is set in f->options. Memory is allocated for this routine on the heap, so the user is responsible to delete this memory using free().

Parameters:

	f	the filter, needs the options from it to check how to transform
	s	filter s plane poles and zeros

Returns:: a pointer to the z plane representation

Definition at line 8 of file zplane_transform.c.

References bpm_error(), c_div(), c_exp(), c_scale(), c_sum(), complex(), MATCHED_Z_TRANSFORM, filterrep_t::npoles, filterrep_t::nzeros, filter_t::options, filterrep_t::pole, and filterrep_t::zero.

Referenced by create_filter().

EXTERN void print_filter_representation	(	FILE *	of,
		filterrep_t *	r
	)

Prints the filter representation in terms of poles and zeros to the filepointer.

Parameters:

	of	the filepointer, use "stdout" to print to the terminal
	r	the filter representation to be printed

Returns:: void

Display filter representation

Definition at line 8 of file print_filter_representation.c.

References c_imag(), c_real(), filterrep_t::npoles, filterrep_t::nzeros, filterrep_t::pole, and filterrep_t::zero.

Referenced by print_filter().

EXTERN int normalise_filter	(	filter_t *	f,
		filterrep_t *	s
	)

Normalises the Butterworth, Chebyshev or Bessel filters to be Bandpass/stop or Low/Highpass

Parameters:

	f	the filter
	s	the filter's representation in the s plane

Returns:: BPM_SUCCESS upon success or BPM_FAILURE upon failure.

Definition at line 7 of file normalise_filter.c.

References bpm_error(), c_div(), c_scale(), complex(), HIGHPASS, LOWPASS, filterrep_t::npoles, filterrep_t::nzeros, filter_t::options, filterrep_t::pole, filter_t::w_alpha1, filter_t::w_alpha2, and filterrep_t::zero.

Referenced by create_filter().

EXTERN int calculate_filter_coefficients ( filter_t * f )

Calculates the filter coefficients from the z plane representation for Butterworth, Chebyshev, Bessel and Resonators. Before this routine is called, one has to make sure that the member cplane, which holds a pointer to the filter's representation in the complex plane is set. This routine than calculates the filter coefficients and stores them in f->xc ( coefficients of x[n], x[n-1], x[n-2]...) and f->yc ( coefficients of y[n-1], y[n-2], y[n-3], ... in case of IIR filters ).

Parameters:

the filter, having it's f->cplane member set to the z plan representation

Returns:: BPM_SUCCESS upon success or BPM_FAILURE upon failure.

Calculates the filter coefficients from the poles and zeros in the cplane representation... Also calculates the filter gains...

Definition at line 56 of file calculate_filter_coefficients.c.

References _eval_complex_polynomial(), _expand_complex_polynomial(), filter_t::alpha1, filter_t::alpha2, BANDPASS, BANDSTOP, c_abs(), c_div(), c_mult(), c_real(), c_sqrt(), complex(), filter_t::cplane, filter_t::dc_gain, filter_t::fc_gain, filter_t::gain, filter_t::hf_gain, HIGHPASS, LOWPASS, MAXPZ, filterrep_t::npoles, filter_t::nxc, filter_t::nyc, filterrep_t::nzeros, filter_t::options, filterrep_t::pole, filter_t::xc, filter_t::yc, and filterrep_t::zero.

Referenced by create_filter().

EXTERN int gaussian_filter_coeffs ( filter_t * f )

Calculates the gaussian filter coefficients from the original gaussian filter implementation in the digital downconversion algortithm in Yury's code. Note that this filter is implemented as a FIR non-causal filter.

Parameters:

the filter structure with the coefficients to fill

Returns:: BPM_SUCCESS upon success or BPM_FAILURE upon failure.

Definition at line 8 of file gaussian_filter_coeffs.c.

References bpm_error(), dround(), filter_t::f1, filter_t::fs, filter_t::gain, filter_t::gauss_cutoff, GAUSSIAN_SIGMA_BW, MAXPZ, filter_t::ns, filter_t::nxc, filter_t::nxc_ac, filter_t::options, filter_t::xc, and filter_t::xc_ac.

Referenced by create_filter().

EXTERN int _expand_complex_polynomial	(	complex_t *	w,
		int	n,
		complex_t *	a
	)

Helper routine to expand a complex polynomial from a set of zeros.

Parameters:

	w	array of complex zeros for the polynomial
	n	nunber of zeros
	a	array of coeffiecients for the polynomial that is returned

Returns:: BPM_SUCCESS upon success or BPM_FAILURE upon failure.

Calculate the polynomial coefficients in a0 + a1 * z + a2 * z^2 + a3 * z^3 + ... = (z-w1)(z-w2)(z-w3)... from the n polynomial's zero's "w" returns the results in a, the array of coefficients...

Definition at line 8 of file calculate_filter_coefficients.c.

References bpm_error(), c_imag(), c_mult(), c_neg(), c_sum(), complex(), and FILT_EPS.

Referenced by calculate_filter_coefficients(), and create_resonator_representation().

EXTERN complex_t _eval_complex_polynomial	(	complex_t *	a,
		int	n,
		complex_t	z
	)

Helper routine to evaluate a complex polynomial for value z

Parameters:

	a	array of coeffiecients for the polynomial that is returned
	n	number of zeros
	z	the value for which to evalute the polynomial

Returns:: the value of the polynomial for z ( complex_t )

Definition at line 44 of file calculate_filter_coefficients.c.

References c_mult(), c_sum(), and complex().

Referenced by calculate_filter_coefficients(), and create_resonator_representation().

int ddc	(	doublewf_t *	w,
		double	f,
		filter_t *	filter,
		complexwf_t *	dcw
	)

Do a digital downconversion on the waveform f. The routine returns a complex DC waveform "wdc".

Parameters:

	w	The waveform of doubles to process
	f	The frequency of the digital local oscillator
	filter	The lowpass filter to get rid of the 2omega component
	dcw	The complex DC waveform

Returns:: BPM_SUCCESS upon success, BPM_FAILURE upon failure

Definition at line 78 of file ddc.c.

References _check_ddc_buffers(), apply_filter(), complexwf_setimag(), complexwf_setreal(), complexwf_t::fs, doublewf_t::fs, complexwf_t::ns, doublewf_t::ns, and doublewf_t::wf.

Referenced by ddc_sample_waveform(), and ddc_waveform().

EXTERN int fft_initialise ( int ns )

This one initialised the FFT buffers, checks whether they are large enough for the given number of samples and frees and re-allocates memory where necessary

Parameters:

The number of samples in the waveforms to be transformed

Returns:: BPM_SUCCESS upon succes, BPM_FAILURE upon failure

Definition at line 130 of file discrete_fourier_transforms.c.

References bpm_error(), and fft_gen_tables().

EXTERN int complexfft	(	complexwf_t *	z,
		int	fft_mode
	)

Executes a complex fast fourier transform in line. See the reference guide for details.

Parameters:

	z	The complex waveform to transform (original waveform is destroyed) Note that the number of samples need to be a power of 2.
	fft_mode	Specifies whether to do the forward or backward transform

Returns:: BPM_SUCCESS upon succes, BPM_FAILURE upon failure

Definition at line 178 of file discrete_fourier_transforms.c.

References _check_fft_buffers(), _is_pow2(), bpm_error(), bpm_warning(), cdft(), FFT_BACKWARD, FFT_FORWARD, complex_t::im, complexwf_t::ns, complex_t::re, and complexwf_t::wf.

EXTERN int realfft	(	doublewf_t *	y,
		int	fft_mode,
		complexwf_t *	z
	)

Executes a real fast fourier transform, between the real waveform y and the complex waveform z. See documentation for further explanation.

Parameters:

	y	Pointer to the real wavefrom
	fft_mode	Specifies whether to do the forward or backward transform
	z	Pointer to the complex waveform

Returns:: BPM_SUCCESS upon succes, BPM_FAILURE upon failure

Definition at line 225 of file discrete_fourier_transforms.c.

References _check_fft_buffers(), _is_pow2(), bpm_error(), bpm_warning(), FFT_BACKWARD, FFT_FORWARD, complex_t::im, complexwf_t::ns, rdft(), complex_t::re, doublewf_t::wf, and complexwf_t::wf.

Digital Signal Processing Routines

Detailed Description

The digital filtering routines

General usage

The Bessel, Butterworth and Chebyshev filters

The Resonator filter

The gaussian filter

The Digital Downconversion Algorithm (DDC)

Discrete (Fast) Fourier Transforms

Complex Discrete Fourier Transform

Real Discrete Fourier Transform

Reference for FFT routines

Copyright statement for FFT routines

DSP example program

Files

Data Structures

Defines

Functions

Define Documentation

Function Documentation


Files
file	bpm_dsp.h
	libbpm digital signal processing routines
file	calculate_filter_coefficients.c
file	create_filter.c
file	create_resonator_representation.c
file	create_splane_representation.c
file	ddc.c
file	delete_filter.c
file	discrete_fourier_transforms.c
file	filter_impulse_response.c
file	filter_step_response.c
file	gaussian_filter_coeffs.c
file	normalise_filter.c
file	print_filter.c
file	print_filter_representation.c
file	zplane_transform.c
Data Structures
struct	filterrep_t
struct	filter_t
Defines
#define	BESSEL
#define	BUTTERWORTH
#define	CHEBYSHEV
#define	RAISEDCOSINE
#define	RESONATOR
#define	GAUSSIAN
#define	BILINEAR_Z_TRANSFORM
#define	MATCHED_Z_TRANSFORM
#define	NO_PREWARP
#define	CAUSAL
#define	ANTICAUSAL
#define	NONCAUSAL
#define	GAUSSIAN_SIGMA_BW
#define	LOWPASS
#define	HIGHPASS
#define	BANDPASS
#define	BANDSTOP
#define	NOTCH
#define	ALLPASS
#define	FIR
#define	IIR
#define	MAXORDER
#define	MAXPZ
#define	FILT_EPS
#define	MAX_RESONATOR_ITER
#define	FFT_FORWARD
#define	FFT_BACKWARD
Functions
EXTERN filter_t *	create_filter (char name[], unsigned int options, int order, int ns, double fs, double f1, double f2, double par)
EXTERN int	apply_filter (filter_t f, double wf)
EXTERN void	print_filter (FILE of, filter_t f)
EXTERN void	delete_filter (filter_t *f)
EXTERN int	filter_step_response (filter_t f, double wf, int itrig)
EXTERN int	filter_impulse_response (filter_t f, double wf, int itrig)
EXTERN filterrep_t *	create_splane_representation (filter_t *f)
EXTERN filterrep_t *	create_resonator_representation (filter_t *f)
EXTERN filterrep_t *	zplane_transform (filter_t f, filterrep_t s)
EXTERN void	print_filter_representation (FILE of, filterrep_t r)
EXTERN int	normalise_filter (filter_t f, filterrep_t s)
EXTERN int	calculate_filter_coefficients (filter_t *f)
EXTERN int	gaussian_filter_coeffs (filter_t *f)
EXTERN int	_expand_complex_polynomial (complex_t w, int n, complex_t a)
EXTERN complex_t	_eval_complex_polynomial (complex_t *a, int n, complex_t z)
EXTERN int	ddc_initialise (int ns, double fs)
EXTERN void	ddc_cleanup (void)
int	ddc (doublewf_t w, double f, filter_t filter, complexwf_t *dcw)
EXTERN int	fft_gen_tables (void)
EXTERN int	fft_initialise (int ns)
EXTERN void	fft_cleanup (void)
EXTERN int	complexfft (complexwf_t *z, int fft_mode)
EXTERN int	realfft (doublewf_t y, int fft_mode, complexwf_t z)

	ns	Nuber of samples in waveforms to be processed
	fs	The sampling frequency of the waveforms