Finite Impulse Filter and decimation test program
Revision 0:3425cd32678e, committed 2010-03-26
- Comitter:
- timowen
- Date:
- Fri Mar 26 12:53:06 2010 +0000
- Commit message:
Changed in this revision
| main.cpp | Show annotated file Show diff for this revision Revisions of this file |
| mbed.bld | Show annotated file Show diff for this revision Revisions of this file |
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/main.cpp Fri Mar 26 12:53:06 2010 +0000
@@ -0,0 +1,127 @@
+
+// fir_int.cpp Tm Owen 25/3/2010
+
+// FIR decimation filter implementation using integers on 3 channel data.
+// Test/timing program for implementing a simple FIR filter
+// This is a low pass filter and decimator implemented on the MBED
+// as an integer calculation in native 32 bit precision.
+// The normal filter coefficcients ( eg, see program SCOPEFIR at XXXX.com - free use for 60 days) are
+// given as floats that add up to about 1. For an integer implementation
+// the coefficients must be multiplied by 2^N to give reasonable integers.
+// Given 32 bit arithmetic, we can handle data and coefficients as follows
+// 32 = COEF + DATA + FILTER_LENGTH all expressed as bits
+// We are dealing with 13 bit data and a 24 coef long filter ( 5 bit )
+// so we can have coefficients up to 14 bits - data and coef resolution is
+// balanced. Anything bigger will overflow 32 bits somewhere.
+// SO;- multiply the coefs by 2^14 ( = 16384 ), do the sums and divide the
+// answer by 2^14 or shift it right 14 places.
+// I put in the multiplied coefs, could be done once in software.
+// DECIM is the decimation ratio - the number of input samples that will be
+// consumed for each output sample. Normally this is a power of 2, but
+// it can be any integer? You need to design the coefficients based on the
+// sampling rate and output rate - the filter cutoff frequency must be less than
+// half the output sample frequency.
+
+// I get this to run at about 8 microsecs per 3 output values - there are 72 int mults and around 125 adds.
+// in the for loop - quite fast enough for most data logging functions but a bit slow for audio.
+
+// input data is in a .csv file, data is output to a similar file.
+
+#include "mbed.h"
+Serial pc(USBTX, USBRX); // tx, rx
+DigitalOut led(LED1); // useful indicator of progress
+FILE *fp, *fp1;
+Timer t; // output write file - bin
+LocalFileSystem local("local");
+#define SHIFT 15 // coefficient multiplier as a power of 2
+#define COEF 24 // number of coefficients in fir
+#define DECIM 3 // decimation factor required
+#define BUFSIZE 400 // number of input samples to process
+
+int xarr[COEF], yarr[COEF], zarr[COEF];
+
+// The FIR filter values shown below are good for 13 bit data at an oputput cutoff of 1/8 input rate
+// and a -80db cutoff at 1/4 of input rate, e.g for 800 s/s gives 200 s/s with 100Hz -3db point.
+int coef[] = {-4,
+33,
+208,
+528,
+720,
+3121,
+-802,
+-1784,
+-1113,
+2081,
+6756,
+10245,
+10245,
+6756,
+2081,
+-1113,
+-1784,
+-802,
+312,
+720,
+528,
+208,
+33,
+-4
+ };
+int x,y,z, count =0;
+int xout=0,yout=0,zout=0;
+int once =1;
+int bufx[BUFSIZE],bufy[BUFSIZE],bufz[BUFSIZE],buf_sump =0,buf_curp;
+int obufx[BUFSIZE/DECIM],obufy[BUFSIZE/DECIM],obufz[BUFSIZE/DECIM];
+int elapse_us;
+main()
+ {
+ int ind;
+ pc.printf("\n\rFIR_INT is an integer implementation of a low pass FIR decimation filter");
+ pc.printf("\n\r...................... Tim Owen .... March 2010..........................\n\r\n\r");
+ pc.printf("\n\r Filter length = %d, Decimation factor = %d, Shift factor = 2^%d\n\r", COEF, DECIM,SHIFT);
+ pc.printf("\n\rOpen two files - one in and one out.........");
+ fp = fopen("/local/TXT2.TXT", "r") ;
+ if(!fp) pc.printf("\n\rNo go input file");
+ fp1 = fopen("/local/TXTFILT.TXT","w");
+ if(!fp1) pc.printf("\n\rNo go output file");
+ pc.printf("\n\r coefficients %d %d %d %d %d %d %d %d %d %d %d %d", coef[0],coef[1],coef[2],coef[3],coef[4], coef[5], coef[6],coef[7],coef[8],coef[9],coef[10],coef[11]);
+ t.start();
+
+ for(ind = 0; ind < BUFSIZE ; ind++)
+ { // copy data to array so we can time without reads
+ if( fscanf(fp, "%d,%d,%d",&bufx[ind],&bufy[ind],&bufz[ind]) != EOF){}; // fscanf > ptrs
+ }
+ t.reset(); // start timing
+
+ while((count < BUFSIZE/DECIM ) && (feof(fp) == 0))
+ { // go round this loop once for each OUTPUT value
+
+ obufx[count] =0; // zero the OUTPUT value
+ obufy[count] =0;
+ obufz[count] =0;
+
+ buf_sump = buf_curp; // set sumation index from current start position
+
+ for( ind = 0; ind < COEF; ind++)
+ {// do the summation
+
+ obufx[count] +=(bufx[buf_sump] * coef[ind]);
+ obufy[count] +=(bufy[buf_sump] * coef[ind]);
+ obufz[count] +=(bufz[buf_sump] * coef[ind]);
+ buf_sump++;
+ }
+ buf_curp += DECIM ; // set up start for next sumation
+ count++;
+ }
+ elapse_us = t.read_us();
+ for( ind =0; ind < BUFSIZE/DECIM ; ind++)
+ {
+ fprintf(fp1,"%d,%d,%d\r",(obufx[ind]>>SHIFT),(obufy[ind]>>SHIFT),(obufz[ind]>>SHIFT));
+ }
+ fclose(fp);
+ fclose(fp1);
+ pc.printf("\n\rSamples = %d, elapse time = %dus us/count = %d\n\r BFN..... \n\r", count,elapse_us, elapse_us/count);
+ }
+
+
+
\ No newline at end of file
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/mbed.bld Fri Mar 26 12:53:06 2010 +0000 @@ -0,0 +1,1 @@ +http://mbed.org/users/mbed_official/code/mbed/builds/49a220cc26e0