Differences

This shows you the differences between two versions of the page.

@@ Line 1: / Line 1: @@
+====== FFTW demo - DTMF decoder ======
+Hanging on IRC, we decided to try how fast can I write a DTMF decoder with library FFT. I estimated 30 minutes. Finally, it took me exactly one hour, and the code quality is terrible.
+Yes, I know that computing these 8 bins from definition would be faster and would not require fftw.
+Now the whole world sees how I code for competitions where only right output, not the code quality, is evaluated.
+<code c>
+#include <stdlib.h>
+#include <stdio.h>
+#include <inttypes.h>
+#include <string.h>
+#include <unistd.h>
+#include <math.h>
+#include <fftw3.h>
+#include <fcntl.h>
+#include <sys/mman.h>
+#include <sys/types.h>
+#include <unistd.h>
+#include <sys/stat.h>
+#include <sys/ioctl.h>
+#include <string.h>
+#include <linux/fs.h>
+double hamming(int i, int length) {
+  double a, b, w, N1;
+  a = 25.0/46.0;
+  b = 21.0/46.0;
+  N1 = (double)(length-1);
+  w = a - b*cos(2*i*M_PI/N1);
+  return w;
+}
+int main(int argc, char **argv) {
+  int i;
+  float* acc;
+  float* fftbuf;
+  int fftsize = 512;
+  acc = (float*) calloc(fftsize * 2, sizeof(float));
+  fftbuf = (float*) calloc(fftsize * 2, sizeof(float));
+  /* Setup FFTW */
+  int N = fftsize;
+  fftwf_complex *fftw_in, *fftw_out;
+  fftwf_plan p;
+  fftw_in = (fftwf_complex*) fftwf_malloc(sizeof(fftwf_complex) * N);
+  fftw_out = (fftwf_complex*) fftwf_malloc(sizeof(fftwf_complex) * N);
+  p = fftwf_plan_dft_1d(N, fftw_in, fftw_out, FFTW_FORWARD, FFTW_ESTIMATE);
+  float * win = malloc(fftsize * sizeof(float));
+  for(i=0; i<fftsize; i++) {
+    win[i] = hamming(i, fftsize);
+  }
+  int16_t * data;
+  int fd = open("/tmp/nahravka.wav", O_RDONLY);
+  size_t dsize = 168204;
+  data = (int16_t*)mmap(NULL, dsize, PROT_READ, MAP_PRIVATE, fd, 0);
+  int samples = dsize/2;
+  /* compute smooth func */
+  int16_t * mx = malloc(sizeof(int16_t) * dsize);
+  for(int i = 0; i<samples-50; i++) {
+    int max = 0;
+    for(int j = i; j<i+30; j++) {
+      if(abs(data[j]) > max) {
+        max = abs(data[j]);
+      }
+    }
+    mx[i] = max;
+    /*if(i%20 == 1) {
+      printf("max %i\n", max);
+    }*/
+  }
+#define sigth1 15000
+#define sigth2 10000
+#define sigstart 32000
+//printf("HERE\n");
+  int pos = sigstart;
+  while(1) {
+    /* find start of burst */
+    while(pos++) {
+//      printf("pos %i mx %i\n", pos, mx[pos]);
+      if(mx[pos] > sigth1) {
+        break;
+      }
+    }
+    /* find end of burst */
+    int endpos = pos+10;
+    while(endpos++) {
+//      printf("endpos %i mx %i\n", endpos, mx[endpos]);
+      if(mx[endpos] < sigth2) {
+        break;
+      }
+    }
+    /* align */
+    int center = (pos+endpos)/2;
+    int dur = endpos-pos;
+    int start = center-fftsize/2;
+    //printf("len = %i (%i - %i), start %i\n", endpos-pos, pos, endpos, start);
+    /* Read data */
+    for(i=0; i<fftsize; i++) {
+      float sample = (float) ( data[i+start] );
+      fftbuf[i*2] = sample * win[i];
+      fftbuf[i*2+1] = 0.0;
+    }
+    // Transform
+    memcpy(fftw_in, fftbuf, 2*fftsize * sizeof(float));
+    fftwf_execute(p);
+    float* fout = (float*)fftw_out;
+    // Read output
+    float* rout = (float*)malloc(sizeof(float)*fftsize);
+    for(int i = 0; i<fftsize; i++) {
+      float orig = hypot(fftbuf[i*2],fftbuf[i*2+1]);
+      float ret = hypot(fout[i*2], fout[i*2+1]);
+      rout[i] = ret;
+      //printf("%f %f\n", orig, ret);
+    }
+    /* find two maximums */
+    int max1 = 8;
+    int max2 = 8;
+    rout[511] = 0.0;
+    float dtmf[] = { 1209,  1336,  1477,  1633, 697, 770, 852, 941 };
+    int bins[9];
+    for(int i = 0; i<8; i++) {
+      bins[i] = (float)fftsize * dtmf[i]/8000.0;
+      //printf("bin %i = %f\n", bins[i], rout[bins[i]]);
+    }
+    bins[8] = 511;
+    /* horiz */
+    for(int i = 0; i<4; i++) {
+      if(rout[bins[i]] > rout[bins[max1]]) {
+        max1 = i;
+      }
+    }
+    /* vert */
+    for(int i = 4; i<8; i++) {
+      //printf("%i cmp %f %f\n", i, rout[bins[i]], rout[bins[max2]]);
+      if(rout[bins[i]] > rout[bins[max2]]) {
+        max2 = i;
+      }
+    }
+    max2 -= 4;
+    //printf("%i %i\n", max1, max2);
+    int res[4][4] = {{1, 2, 3, 'A'},
+                   {4, 5, 6, 'B'},
+                   {7, 8, 9, 'C'},
+                   {'*', 0, '#', 'D'} };
+    int rt = res[max2][max1];
+    if(rt>10) {
+      printf(" *** %c ***\n", rt);
+    } else {
+      printf(" *** %i ***\n", rt);
+    }
+    pos = endpos+10;
+  }
+}
+</code>

wiki

User Tools

Site Tools

Differences

Page Tools