KissFFT output of kiss_fftr

Question

I'm receiving PCM data trough socket connection in packets containing 320 samples. Sample rate of sound is 8000 samples per second. I am doing with it something like this:

int size = 160 * 2;//160;
int isinverse = 1;
kiss_fft_scalar zero;
memset(&zero,0,sizeof(zero));
kiss_fft_cpx fft_in[size];
kiss_fft_cpx fft_out[size];
kiss_fft_cpx fft_reconstructed[size];

kiss_fftr_cfg fft = kiss_fftr_alloc(size*2 ,0 ,0,0);
kiss_fftr_cfg ifft = kiss_fftr_alloc(size*2,isinverse,0,0);

for (int i = 0; i < size; i++) {
    fft_in[i].r = zero;
    fft_in[i].i = zero;
    fft_out[i].r = zero;
    fft_out[i].i = zero;
    fft_reconstructed[i].r = zero;
    fft_reconstructed[i].i = zero;
}

// got my data through socket connection

for (int i = 0; i < size; i++) {
     // samples are type of short
     fft_in[i].r = samples[i];
     fft_in[i].i = zero;
     fft_out[i].r = zero;
     fft_out[i].i = zero;
 }

 kiss_fftr(fft, (kiss_fft_scalar*) fft_in, fft_out);
 kiss_fftri(ifft, fft_out, (kiss_fft_scalar*)fft_reconstructed);

 // lets normalize samples
 for (int i = 0; i < size; i++) {
     short* samples = (short*) bufTmp1;
     samples[i] = rint(fft_reconstructed[i].r/(size*2));
 }

After that I fill OpenAL buffers and play them. Everything works just fine but I would like to do some filtering of audio between kiss_fftr and kiss_fftri. Starting point as I think for this is to convert sound from time domain to frequency domain, but I don't really understand what kind of data I'm receiving from kiss_fftr function. What information is stored in each of those complex number, what its real and imaginary part can tell me about frequency. And I don't know which frequencies are covered (what frequency span) in fft_out - which indexes corresponds to which frequencies.

I am total newbie in signal processing and Fourier transform topics.

Any help?

Answer 1

What you might want to investigate is FFT fast convolution using overlap add or overlap save algorithms. You will need to expand the length of each FFT by the length of the impulse of your desired filter. This is because (1) FFT/IFFT convolution is circular, and (2) each index in the FFT array result corresponds to almost all frequencies (a Sinc shaped response), not just one (even if mostly near one), so any single bin modification will leak throughout the entire frequency response (except certain exact periodic frequencies).

Answer 2

Before you jump in with both feet into a C implementation, get familiar with digital filters, esp FIR filters.

You can design the FIR filter using something like GNU Octave's signal toolbox. Look at the command fir1(the simplest), firls, or remez. Alternately, you might be able to design a FIR filter through a web page. A quick web search for "online fir filter design" found this (I have not used it, but it appears to use the equiripple design used in the remez or firpm command )

Try implementing your filter first with a direct convolution (without FFTs) and see if the speed is acceptable -- that is an easier path. If you need an FFT-based approach, there is a sample implementation of overlap-save in the kissfft/tools/kiss_fastfir.c file.

Answer 3

I will try to answer your questions directly.

// a) the real and imaginary components of the output need to be combined to calculate the amplitude at each frequency. 

float ar,ai,scaling; 

scaling=1.0/(float)size;

// then for each output [i] from the FFT...

ar = fft_out[i].r; 
ai = fft_out[i].i; 
amplitude[i] = 2.0 * sqrtf( ar*ar + ai*ai ) * scaling ;

// b) which index refers to which frequency? This can be calculated as follows. Only the first half of the FFT results are needed (assuming your 8KHz sampling rate) 

for(i=1;i<(size/2);i++) freq = (float)i / (1/8000) / (float)size ; 

// c)  phase (range +/- PI) for each frequency is calculated like this: 

phase[i] = phase = atan2(fft_out[i].i / fft_out[i].r);