This shows you the differences between two versions of the page.
dtmf [2016-01-24 17:32:38] |
dtmf [2016-01-24 17:32:38] (current) |
||
---|---|---|---|
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> |