Differences

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

Link to this comparison view

user:jenda:dtmf [2016/01/24 17:32] (current)
jenda created
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>​
 
Except where otherwise noted, content on this wiki is licensed under the following license: CC Attribution-Noncommercial-Share Alike 4.0 International
Recent changes RSS feed Donate Powered by PHP Valid XHTML 1.0 Valid CSS Driven by DokuWiki