Digital Signal Processing Routines

Detailed Description

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

The digital filtering routines

General usage

Setup a filter using the create_filter() routine.

    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 implemented filters are :

BESSEL : Bessel IIR filter
BUTTERWORTH : Butterwordth IIR filter
CHEBYSHEV : Chebyshev IIR filter
RESONATOR : Resonators
GAUSSIAN : Non-causal Gaussian FIR filter

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 :

BESSEL : the parameter defines the ripple in dB, has to be negative !
RESONATOR : the parameter gives the Q value of the resonator, if you want to have a pure oscillator (so infinite Q), then set the parameter to a negative number or zero.
GAUSSIAN : the filter cut-off parameter, or the fraction of the gaussian convolution function below which it is set to 0. (default is 0.001)

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, wave );

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, wave, 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, wave, 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 gaussian filter is implemented as a FIR convolution with both causal and anti-causal coefficients. Note that the frequency given is treated as the -3dB level for the gaussian. There is an option to restore the definition for bandwith which was used in early ESA processing, being the gaussian sigma, use GAUSSIAN_SIGMA_BW.

The Digital Downconversion Algorithm (DDC)

The digital downconversion routine was developped to process digitised BPM waveforms and to retreive their position and amplitude. It basically implements an RF mixer in software. You need to supply it with the doublewf_t holding the waveform to mix down and the frequency for the software LO. Also you need to give a pointer to a low-pass filter in order to filter out the resulting double frequency component from the downmixing. The routine

    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

    ddc_initialise();

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

    ddc_cleanup();

Discrete (Fast) Fourier Transforms

The FFT routines in the dsp section of libbpm are based upon the General Purpose FFT Package by Takuya OOURA, 1996-2001, see http://www.kurims.kyoto-u.ac.jp/~ooura/fft.html More specifically on it's split-radix fast version (fftsg). These set of routines needs a buffer for bitswapping an a buffer to store a table with sin and cos values so they needn't be calculated for every FFT. The routine

    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

    fft_cleanup( );

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

Complex Discrete Fourier Transform

The first one is

    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

FFT_FORWARD : forward discrete Fourier transform (plus-sign)

$X[k] = \sum_{j=0}^{n-1} x[j]*\exp(2*\pi*i*j*k/n), 0<=k<n$

FFT_BACKWARD : backward discrete Fourier transform (minus-sign)

$X[k] = \sum_{j=0}^{n-1} x[j]*\exp(-2*\pi*i*j*k/n), 0<=k<n$

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

The second routine implements the real discrete Fourier transform when having FFT_FORWARD and the other way around when having FFT_BACKWARD.

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

So for FFT_FORWARD

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

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

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

Masatake MORI, Makoto NATORI, Tatuo TORII: Suchikeisan, Iwanamikouzajyouhoukagaku18, Iwanami, 1982 (Japanese)
Henri J. Nussbaumer: Fast Fourier Transform and Convolution Algorithms, Springer Verlag, 1982
C. S. Burrus, Notes on the FFT (with large FFT paper list) http://www-dsp.rice.edu/research/fft/fftnote.asc

Copyright statement for FFT routines

Copyright(C) 1996-2001 Takuya OOURA email: ooura@mmm.t.u-tokyo.ac.jp download: http://momonga.t.u-tokyo.ac.jp/~ooura/fft.html You may use, copy, modify this code for any purpose and without fee. You may distribute this ORIGINAL package.

DSP example program

There is an example program, which can be found in the examples directory under dsp. It shows how to work with the filtering and the DDC routines...

    #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_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;
    }

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 norm_phase.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, doublewf_t *w)

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

EXTERN void delete_filter (filter_t *f)

EXTERN int filter_step_response (filter_t *f, doublewf_t *w, int itrig)

EXTERN int filter_impulse_response (filter_t *f, doublewf_t *w, 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, doublewf_t *bufre, doublewf_t *bufim)

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)

EXTERN void norm_phase (double *phase)

Define Documentation

#define BESSEL

Bitmask for Bessel filter

Definition at line 384 of file bpm_dsp.h.

Referenced by create_filter(), and create_splane_representation().

#define BUTTERWORTH

Bitmask for Butterworth filter

Definition at line 385 of file bpm_dsp.h.

Referenced by create_filter(), and create_splane_representation().

#define CHEBYSHEV

Bitmask for Chebyshev filter

Definition at line 386 of file bpm_dsp.h.

Referenced by create_filter(), and create_splane_representation().

#define RAISEDCOSINE

Bitmask for Raised Cosine filter

Definition at line 387 of file bpm_dsp.h.

#define RESONATOR

Bitmask for Resonator filter

Definition at line 388 of file bpm_dsp.h.

Referenced by create_filter(), and get_mode_response().

#define GAUSSIAN

Bitmask for Gaussian filter

Definition at line 389 of file bpm_dsp.h.

Referenced by create_filter().

#define BILINEAR_Z_TRANSFORM

Get z poles via bilinear z transform from s plane

Definition at line 391 of file bpm_dsp.h.

#define MATCHED_Z_TRANSFORM

Get z poles via matches z transform from s plane

Definition at line 392 of file bpm_dsp.h.

Referenced by zplane_transform().

#define NO_PREWARP

Don't do the prewarp correction

Definition at line 393 of file bpm_dsp.h.

Referenced by create_filter().

#define CAUSAL

Filter is purely causal (only depends on past )

Definition at line 394 of file bpm_dsp.h.

Referenced by apply_filter(), create_filter(), and print_filter().

#define ANTICAUSAL

.... purely anticausal (only depends on future)

Definition at line 395 of file bpm_dsp.h.

Referenced by apply_filter(), and print_filter().

#define NONCAUSAL

Filter is both causal and acausal

Definition at line 396 of file bpm_dsp.h.

Referenced by create_filter().

#define GAUSSIAN_SIGMA_BW

Gaussian sigma bandwidth in stead of -3 dB (def)

Definition at line 397 of file bpm_dsp.h.

Referenced by gaussian_filter_coeffs().

#define LOWPASS

Normalise filter as lowpass

Definition at line 399 of file bpm_dsp.h.

Referenced by calculate_filter_coefficients(), and normalise_filter().

#define HIGHPASS

Normalise filter as highpass

Definition at line 400 of file bpm_dsp.h.

Referenced by calculate_filter_coefficients(), and normalise_filter().

#define BANDPASS

Normalise filter as bandpass

Definition at line 401 of file bpm_dsp.h.

Referenced by calculate_filter_coefficients(), get_mode_response(), and normalise_filter().

#define BANDSTOP

Normalise filter as bandstop

Definition at line 402 of file bpm_dsp.h.

Referenced by calculate_filter_coefficients(), create_resonator_representation(), and normalise_filter().

#define NOTCH

Normalise filter as notch filter (=bandstop)

Definition at line 403 of file bpm_dsp.h.

#define ALLPASS

Normalise filter as allpass ( resonator )

Definition at line 404 of file bpm_dsp.h.

Referenced by create_resonator_representation().

#define FIR

Filter is of FIR type

Definition at line 406 of file bpm_dsp.h.

Referenced by apply_filter(), and create_filter().

#define IIR

Filter is of IIR type

Definition at line 407 of file bpm_dsp.h.

Referenced by create_filter().

#define MAXORDER

Maximum filter order

Definition at line 409 of file bpm_dsp.h.

#define MAXPZ

Maximum number of poles and zeros >2*MAXORDER

Definition at line 410 of file bpm_dsp.h.

Referenced by calculate_filter_coefficients(), create_resonator_representation(), and gaussian_filter_coeffs().

#define FILT_EPS

A small number used in bpmdsp

Definition at line 411 of file bpm_dsp.h.

Referenced by _expand_complex_polynomial(), create_resonator_representation(), and print_filter().

#define MAX_RESONATOR_ITER

Maximum iterations in resonator poles calculation

Definition at line 412 of file bpm_dsp.h.

Referenced by create_resonator_representation().

#define FFT_FORWARD

Perform FFT from time -> frequency

Definition at line 414 of file bpm_dsp.h.

Referenced by complexfft(), fft_waveform(), and realfft().

#define FFT_BACKWARD

Perform FFT from frequency -> time

Definition at line 415 of file bpm_dsp.h.

Referenced by complexfft(), and realfft().

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 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 get_mode_response().

EXTERN int apply_filter	(	filter_t *	f,
		doublewf_t *	w
	)

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::nyc, filter_t::options, doublewf_t::wf, filter_t::wfbuffer, filter_t::xc, filter_t::xc_ac, filter_t::xv, filter_t::xv_ac, filter_t::yc, filter_t::yv, and filter_t::yv_ac.

Referenced by ddc(), filter_impulse_response(), filter_step_response(), generate_diodesignal(), 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 ANTICAUSAL, bpm_error(), 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::nxc_ac, filter_t::nyc, filter_t::options, print_filter_representation(), filter_t::xc, filter_t::xc_ac, and filter_t::yc.

EXTERN void delete_filter ( filter_t * f )

Clears the memory that was allocated on the heap for the filter f.

Parameters:

f

a pointer to the filter

Returns:: void

Definition at line 7 of file delete_filter.c.

References filter_t::cplane, and filter_t::wfbuffer.

Referenced by get_mode_response().

EXTERN int filter_step_response	(	filter_t *	f,
		doublewf_t *	w,
		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(), filter_t::ns, and doublewf_t::wf.

EXTERN int filter_impulse_response	(	filter_t *	f,
		doublewf_t *	w,
		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(), filter_t::ns, and doublewf_t::wf.

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:

f

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 BESSEL, bpm_error(), BUTTERWORTH, filter_t::cheb_ripple, CHEBYSHEV, filterrep_t::npoles, filter_t::options, filter_t::order, and filterrep_t::pole.

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:

f

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(), ALLPASS, filter_t::alpha1, BANDSTOP, bpm_error(), 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(), 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 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 BANDPASS, BANDSTOP, bpm_error(), 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:

f

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, 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:

f

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(), 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.

Referenced by calculate_filter_coefficients(), and create_resonator_representation().

EXTERN int ddc_initialise	(	int	ns,
		double	fs
	)

Initialises and allocates memory for the DDC buffers with the correct number of samples and sampling frequency

Parameters:

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

Returns:: BPM_SUCCESS upon success, BPM_FAILURE upon failure

Definition at line 50 of file ddc.c.

References bpm_error(), and doublewf().

EXTERN void ddc_cleanup ( void )

Clears up and frees the buffer memory for the ddc routines

Definition at line 70 of file ddc.c.

References doublewf_delete().

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

Do a digital downconversion on the waveform f. The routine returns a complex DC waveform "wdc". If the buffer arguments are NULL pointers, the DDC routine will use an internal buffer. This is a good option when all the BPMs in the system have the same sampling frequency and number of samples.

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
	bufre	The real ddc buffer
	bufim	The imaginary ddc buffer

Returns:: BPM_SUCCESS upon success, BPM_FAILURE upon failure

Definition at line 78 of file ddc.c.

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

Referenced by ddc_waveform().

EXTERN int fft_gen_tables ( void )

Regenerates the sin/cos tables that are needed for the fast DFT algorithm.

Definition at line 116 of file discrete_fourier_transforms.c.

References bpm_error().

Referenced by fft_initialise().

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:

ns

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 void fft_cleanup ( void )

This routine frees up the memory used by the FFT buffers

Definition at line 163 of file discrete_fourier_transforms.c.

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 bpm_error(), bpm_warning(), 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 230 of file discrete_fourier_transforms.c.

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

Referenced by fft_waveform().

EXTERN void norm_phase ( double * phase )

Normalises the phase, to the interval [0,2pi[

Parameters:

phase

Pointer to the phase value to normalise

Definition at line 8 of file norm_phase.c.

Referenced by complexwf_getphase(), complexwf_getphase_new(), postprocess_waveform(), process_caltone(), and process_waveform().


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	norm_phase.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, doublewf_t w)
EXTERN void	print_filter (FILE of, filter_t f)
EXTERN void	delete_filter (filter_t *f)
EXTERN int	filter_step_response (filter_t f, doublewf_t w, int itrig)
EXTERN int	filter_impulse_response (filter_t f, doublewf_t w, 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, doublewf_t bufre, doublewf_t *bufim)
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)
EXTERN void	norm_phase (double *phase)