digital audio mod11 / Mbed 2 deprecated equalizerseperate

Important changes to repositories hosted on mbed.com

Mbed hosted mercurial repositories are deprecated and are due to be permanently deleted in July 2026.

To keep a copy of this software download the repository Zip archive or clone locally using Mercurial.

It is also possible to export all your personal repositories from the account settings page.

Dependencies:   mbed

FFT.cpp@1:099f1a4c5fc8, 2017-03-29 (annotated)

Committer:
rik
Date:
Wed Mar 29 08:42:31 2017 +0000
Revision:
1:099f1a4c5fc8
Parent:
0:8c128e047ec9 
first compile (full of bugs?)

Who changed what in which revision?

UserRevisionLine numberNew contents of line
rik 1:099f1a4c5fc8 1 #include "FFT.h"
rik 1:099f1a4c5fc8 2
rik 1:099f1a4c5fc8 3
rik 1:099f1a4c5fc8 4 // Calculates and subtracts the signal mean
rik 1:099f1a4c5fc8 5 void subMean(float* data, int datalength){
rik 1:099f1a4c5fc8 6 float mean = 0;
rik 1:099f1a4c5fc8 7 float total = 0;
rik 1:099f1a4c5fc8 8 for (int i = 0; i < datalength; i++)
rik 1:099f1a4c5fc8 9 total += data[i]/datalength;
rik 1:099f1a4c5fc8 10 for (int i = 0; i < datalength; i++)
rik 1:099f1a4c5fc8 11 data[i] = data[i] - mean;
rik 1:099f1a4c5fc8 12 }
rik 1:099f1a4c5fc8 13
rik 1:099f1a4c5fc8 14 void reverseCooleyTukeyRadix2(float* data, complex_num* dataout, int datalength){
rik 1:099f1a4c5fc8 15
rik 1:099f1a4c5fc8 16 complex_num fftbuffer[256];
rik 1:099f1a4c5fc8 17 complex_num twiddlebuffer;
rik 1:099f1a4c5fc8 18
rik 1:099f1a4c5fc8 19 complex_num minusone;
rik 1:099f1a4c5fc8 20 minusone.real = -1;
rik 1:099f1a4c5fc8 21 minusone.imaginary = 0;
rik 1:099f1a4c5fc8 22
rik 1:099f1a4c5fc8 23 // copy iput to output first
rik 1:099f1a4c5fc8 24 for (int i = 0; i < datalength; i++){
rik 1:099f1a4c5fc8 25 dataout[i].real = data[i];
rik 1:099f1a4c5fc8 26 dataout[i].imaginary = 0;
rik 1:099f1a4c5fc8 27 }
rik 1:099f1a4c5fc8 28
rik 1:099f1a4c5fc8 29 for (int stage = 8; stage > -1; stage--){
rik 1:099f1a4c5fc8 30
rik 1:099f1a4c5fc8 31 unsigned int step = 1<<(stage+1);
rik 1:099f1a4c5fc8 32
rik 1:099f1a4c5fc8 33 for (unsigned int i = 0; i < 512; i = i+step){ // Crack it up in parts, i.e. for stage0 -> 512, stage1 -> 256 , 2>128, 3>64, 4>32, 5>16, 6>8, 7>4, 8>2, 9>1 --> 1024/2^stage OR IN STEPS: stage0 -> 2, stage1 -> 4, 2>8, 3>16, 4,32, 5>64, ......... 9>512 --> 1024 - 1024/2^stage+1
rik 1:099f1a4c5fc8 34 // Buffer the first 'half' of the thing, which is the first half step
rik 1:099f1a4c5fc8 35 for (unsigned int k = 0; k < step / 2; k++){
rik 1:099f1a4c5fc8 36 fftbuffer[k].real = dataout[i+k].real;
rik 1:099f1a4c5fc8 37 fftbuffer[k].imaginary = dataout[i+k].imaginary;
rik 1:099f1a4c5fc8 38 }
rik 1:099f1a4c5fc8 39
rik 1:099f1a4c5fc8 40 // Apply twiddle factors to the seconds half (half step) depending on stage and numbering
rik 1:099f1a4c5fc8 41 for (unsigned int k = 0; k < step / 2; k++){
rik 1:099f1a4c5fc8 42 // Read in the correct twiddle factor
rik 1:099f1a4c5fc8 43 // Twiddle factor is recalculated with k*(1024 / 2^stage+1) which can be recalculated to 1<<
rik 1:099f1a4c5fc8 44 twiddlebuffer.imaginary = twiddleLookupTable_Imaginary[k * (1<<(9-stage))];
rik 1:099f1a4c5fc8 45 twiddlebuffer.real = twiddleLookupTable_Real[k * (1 << (9 - stage))];
rik 1:099f1a4c5fc8 46
rik 1:099f1a4c5fc8 47 dataout[i + k + step / 2] = complex_multiply(dataout[i + k + step / 2], twiddlebuffer);
rik 1:099f1a4c5fc8 48 }
rik 1:099f1a4c5fc8 49
rik 1:099f1a4c5fc8 50 // Recalculate the first half by addition of the first and recalculated half
rik 1:099f1a4c5fc8 51 for (unsigned int k = 0; k < step / 2; k++)
rik 1:099f1a4c5fc8 52 dataout[i + k] = complex_add(dataout[i+k], dataout[i+k+step/2]);
rik 1:099f1a4c5fc8 53
rik 1:099f1a4c5fc8 54
rik 1:099f1a4c5fc8 55 // Multiply the second half with minus 1
rik 1:099f1a4c5fc8 56 for (unsigned int k = 0; k < step / 2; k++)
rik 1:099f1a4c5fc8 57 dataout[i + step / 2 + k] = complex_multiply(dataout[i + step/2+k], minusone);
rik 1:099f1a4c5fc8 58
rik 1:099f1a4c5fc8 59 // Add the buffered old values, with the new second half values
rik 1:099f1a4c5fc8 60 for (unsigned int k = 0; k < step / 2; k++)
rik 1:099f1a4c5fc8 61 dataout[i + step / 2 + k] = complex_add(dataout[i + step/2 + k], fftbuffer[k]);
rik 1:099f1a4c5fc8 62 }
rik 1:099f1a4c5fc8 63 }
rik 1:099f1a4c5fc8 64 }
rik 1:099f1a4c5fc8 65
rik 1:099f1a4c5fc8 66 void cooleyTukeyBack (float* dataout, complex_num* datain, int datalength){
rik 1:099f1a4c5fc8 67 complex_num fftbuffer[256];
rik 1:099f1a4c5fc8 68 complex_num twiddlebuffer;
rik 1:099f1a4c5fc8 69
rik 1:099f1a4c5fc8 70 complex_num minusone;
rik 1:099f1a4c5fc8 71 minusone.real = -1;
rik 1:099f1a4c5fc8 72 minusone.imaginary = 0;
rik 1:099f1a4c5fc8 73
rik 1:099f1a4c5fc8 74 for (int stage = 0; stage < 9; stage++){
rik 1:099f1a4c5fc8 75
rik 1:099f1a4c5fc8 76 unsigned int step = 1<<(stage+1);
rik 1:099f1a4c5fc8 77
rik 1:099f1a4c5fc8 78 for (unsigned int i = 0; i < 512; i = i+step){ // Crack it up in parts, i.e. for stage0 -> 512, stage1 -> 256 , 2>128, 3>64, 4>32, 5>16, 6>8, 7>4, 8>2, 9>1 --> 1024/2^stage OR IN STEPS: stage0 -> 2, stage1 -> 4, 2>8, 3>16, 4,32, 5>64, ......... 9>512 --> 1024 - 1024/2^stage+1
rik 1:099f1a4c5fc8 79 // Buffer the first 'half' of the thing, which is the first half step
rik 1:099f1a4c5fc8 80 for (unsigned int k = 0; k < step / 2; k++){
rik 1:099f1a4c5fc8 81 fftbuffer[k].real = datain[i+k].real;
rik 1:099f1a4c5fc8 82 fftbuffer[k].imaginary = datain[i+k].imaginary;
rik 1:099f1a4c5fc8 83 }
rik 1:099f1a4c5fc8 84
rik 1:099f1a4c5fc8 85 // Apply twiddle factors to the seconds half (half step) depending on stage and numbering
rik 1:099f1a4c5fc8 86 for (unsigned int k = 0; k < step / 2; k++){
rik 1:099f1a4c5fc8 87 // Read in the correct twiddle factor
rik 1:099f1a4c5fc8 88 // Twiddle factor is recalculated with N-k*(1024 / 2^stage+1) which can be recalculated to 1<<
rik 1:099f1a4c5fc8 89 twiddlebuffer.imaginary = twiddleLookupTable_Imaginary[512-(k * (1<<(9-stage)))];
rik 1:099f1a4c5fc8 90 twiddlebuffer.real = twiddleLookupTable_Real[512-(k * (1 << (9 - stage)))];
rik 1:099f1a4c5fc8 91
rik 1:099f1a4c5fc8 92 datain[i + k + step / 2] = complex_multiply(datain[i + k + step / 2], twiddlebuffer);
rik 1:099f1a4c5fc8 93 }
rik 1:099f1a4c5fc8 94
rik 1:099f1a4c5fc8 95 // Recalculate the first half by addition of the first and recalculated half
rik 1:099f1a4c5fc8 96 for (unsigned int k = 0; k < step / 2; k++)
rik 1:099f1a4c5fc8 97 datain[i + k] = complex_add(datain[i+k], datain[i+k+step/2]);
rik 1:099f1a4c5fc8 98
rik 1:099f1a4c5fc8 99
rik 1:099f1a4c5fc8 100 // Multiply the second half with minus 1
rik 1:099f1a4c5fc8 101 for (unsigned int k = 0; k < step / 2; k++)
rik 1:099f1a4c5fc8 102 datain[i + step / 2 + k] = complex_multiply(datain[i + step/2+k], minusone);
rik 1:099f1a4c5fc8 103
rik 1:099f1a4c5fc8 104 // Add the buffered old values, with the new second half values
rik 1:099f1a4c5fc8 105 for (unsigned int k = 0; k < step / 2; k++)
rik 1:099f1a4c5fc8 106 datain[i + step / 2 + k] = complex_add(datain[i + step/2 + k], fftbuffer[k]);
rik 1:099f1a4c5fc8 107 }
rik 1:099f1a4c5fc8 108 }
rik 1:099f1a4c5fc8 109 //first convert complex data back to magnitude data then divide by 1/N
rik 1:099f1a4c5fc8 110 for (int i = 0; i < 512; i++)
rik 1:099f1a4c5fc8 111 dataout[i] = (complex_magnitude(datain[i])/512);
rik 1:099f1a4c5fc8 112 }
rik 1:099f1a4c5fc8 113
rik 1:099f1a4c5fc8 114 void equalizer(complex_num* data){
rik 1:099f1a4c5fc8 115
rik 1:099f1a4c5fc8 116 }
rik 1:099f1a4c5fc8 117
rik 1:099f1a4c5fc8 118 void performEqualizer(float* datain, float* dataout, complex_num* freqdata, int datalength){
rik 1:099f1a4c5fc8 119 if (datalength == 512){
rik 1:099f1a4c5fc8 120 // - Calculating and substracting the mean -
rik 1:099f1a4c5fc8 121 subMean(datain, datalength);
rik 1:099f1a4c5fc8 122 // - Applying a hanning window - N multiplications
rik 1:099f1a4c5fc8 123 applyHanningTable(datain, datalength);
rik 1:099f1a4c5fc8 124 // It then applies a radix-2 Cooley Tukey FFT to create a fourier transform
rik 1:099f1a4c5fc8 125 reverseCooleyTukeyRadix2(datain, freqdata, datalength);
rik 1:099f1a4c5fc8 126 // It gets equalized
rik 1:099f1a4c5fc8 127 equalizer(freqdata);
rik 1:099f1a4c5fc8 128 // back to time domain
rik 1:099f1a4c5fc8 129 cooleyTukeyBack(dataout, freqdata, datalength);
rik 1:099f1a4c5fc8 130 }
rik 1:099f1a4c5fc8 131
rik 1:099f1a4c5fc8 132 }