FIR Filters (Part 3)

Part 123

Note: As usual, the code samples are extracted from a standardized library (actually PlainFIR). Please check this page if you are interested in the code..

From the previous readings, we understand that the most complex part of the whole process is the filter design. For sake of simplicity, I decided to program the window method which is simple and efficient. The pseudo code looks like:

  • Save memory space for the vector which will contain the filter coefficients
  • Compute coefficients depending upon the filter type
  • Apply windowing to the coefficients

Here is the function which does the job. It has a low level of optimization and may be improved depending the type of application and the available memory space (Still this RAM limitation…)

void PlainFIR::SetFilter(uint8_t filterType, uint16_t order, uint16_t samplingFrequency, uint8_t windowType, uint16_t transition1, uint16_t transition2) 
/* Order shall be an even number in order to simplify the code */
{
	_order  = order;
	uint16_t taps = order + 1;
	_vFilter = (double *)malloc(taps * sizeof(double)); /* allocate memory for n taps buffer */
	double normTransFreq1 = transition1 / samplingFrequency;
	double normTransFreq2 = transition2 / samplingFrequency;
	/* Compute half + 1 weighing factors, because the filter is symetric */
	for (uint16_t n = 0; n < (order >> 1) + 1; n++) {
		double a = M_PI * (n - (order >> 1));
		double weigthingFactor;
		/* Compute weighing factor */
		switch(filterType){
		case FIR_FIL_TYP_LOW_PASS:
			if (n != (order >> 1)) {
				weigthingFactor = sin(2.0 * normTransFreq1 * a) / a;
			}
			else {
				weigthingFactor = 2.0 * normTransFreq1;
			}
			break;		
		case FIR_FIL_TYP_HIG_PASS:
			if (n != (order >> 1)) {
				weigthingFactor = - sin(2.0 * normTransFreq1 * a) / a;
			}
			else {
				weigthingFactor = 1.0 - (2.0 * normTransFreq1);
			}
			break;		
		case FIR_FIL_TYP_BAN_PASS:
			if (n != (order >> 1)) {
				weigthingFactor = (sin(2.0 * normTransFreq2 * a) - sin(2.0 * normTransFreq1 * a)) / a;
			}
			else {
				weigthingFactor = 2.0 * (normTransFreq2 - normTransFreq1);
			}
			break;		
		case FIR_FIL_TYP_BAN_STOP:
			if (n != (order >> 1)) {
				weigthingFactor = (sin(2.0 * normTransFreq1 * a) - sin(2.0 * normTransFreq2 * a)) / a;
			}
			else {
				weigthingFactor = 1.0 - (2.0 * (normTransFreq2 - normTransFreq1));
			}
			break;			
		};
		/* Apply windowing */
		switch(windowType){
		case FIR_WIN_TYP_BARLETT:
			weigthingFactor *= 1.0 - ((2.0 * abs(n - (order >> 1))) / order);
			break;
		case FIR_WIN_TYP_BLACKMAN:
			weigthingFactor *= 0.42 - (0.5 * cos((2.0 * M_PI * n) / order)) + (0.08 * cos((4.0 * M_PI * n) / order));
			break;
		case FIR_WIN_TYP_RECTANGLE:
			weigthingFactor *= 1.0;
			break;
		case FIR_WIN_TYP_HAMMING:
			weigthingFactor *= 0.5 - (0.46 * cos((2.0 * M_PI * n) / order));
			break;
		case FIR_WIN_TYP_HANN:
			weigthingFactor *= 0.5 - (0.5 * cos((2.0 * M_PI * n) / order));
			break;
		};
		/* Record weighing factors in filter vector */
		_vFilter[n] = weigthingFactor;
		_vFilter[taps - (n +1)] = weigthingFactor;
	}
};

Notes:

  • For sake of code simplicity and standardization, I decided to use odd number of coefficients (filter taps), while this property is mandatory for high-pass filters.
  • M_PI is a constant from the math.h AVR library

Links:
You want to exercise the filter design on your PC. I found that WindowsFIR is a quite convenient (freeware) tool, as long as you understand FIR and that you find the tricky way to export coefficients.

One Comment

  1. […] Next post on same subject Tags: DSP, Programming Comment (RSS)  |  Trackback […]

Leave a Reply

You must be logged in to post a comment.