Convert FFT to IFFT: Steps to Get Right Values Back

[btb4198]

If you involve scaling and Hz in the FFT debug process you will have quite a difficult time getting the phase and frequency correct, or identifying the source of bugs in the code.

I have always found the best way to debug or calibrate an FFT is to use 8 or 16 element arrays.
By generating a simple signal = Cos(2*t), then advancing gradually to things like = 1 + 2*Cos(t) – 3*Sin(2*t), you will be able to control and see what is going on in the frequency domain. Following an FFT, the inverse FFT should return your original signal.

Once you have a robust FFT core you can then apply window functions and frequency scaling external to the FFT.


Recursive FFT code is very elegant mathematically, but often spends half the CPU time pumping frames on and off the stack.
The number of characters in the source code is not really important these days.

what is a good test I should use ?

 

[btb4198]

should I just do
for (int t =0; t<16; t++)
{

list[t] =Cos(2*t);
}
?

 

[btb4198]

You don't "play" sin(2πnf/fs). You specify f in Hz (cycles/s) at this site, not 2πf in radians/s.

You synthesize the samples of a waveform with frequency f using sin(2πnf/fs).

sorry I do not understand.
I have play tone using the Naudio c# class (https://naudio.codeplex.com/)
when I make my tone I do this :
Amplitude * Math.Sin(Math.PI * Frequency * n / 44100D))
if i add the * 2D *∏
it does not sound right,
so I removed it

are you saying that is should sound good?

so I went back and add the 2D * ∏
and you were right, I am able to remove the *2 on the (bin_index/ N)*Fs

but now after the IFFT
it is just completely wrong

 

[btb4198]

olivermsun:

ok if I remove the 2D * ∏
then I have to *2 on the (bin_index/ N)*Fs
for the FFT to work...

and the IFFT works a lot better
it works but there are just some peaks that come up
like
for F = 500hz and F2 1300hz
after the IFFT
and doing thing FFT again
I get peaks at:
399.3642578125Hz
500.64697265625Hz
900.357055664063Hz
1300.06713867188Hz

 

[btb4198]

also
why does my FFT works for
http://onlinetonegenerator.com/
with the (bin_index/ N)*Fs* 2 ?
they most not be doing 2*Pi right?
is that correct ?
is that how it should be ?

 

[olivermsun]

If I play a 440 Hz tone at that website, it sounds correct.
So there may be something wrong with the way you are getting your input stream.
I noticed in an earlier thread by you that there was some confusion about the stereo samples from left and right channels are handled. Have you sorted that out?

 

[Baluncore]

Both frequency and phase are important during testing.
I have not been following details, but I would consider (bin_index/ (N-1)) rather than (bin_index/ N).

 

[btb4198]

olivermsun,

that website is different now, than what it used to be. I have a video of me testing it, how it used to look, I need to find it and I will add it on here later

 

[btb4198]

Baluncore
why
(bin_index/ (N-1)) * fs?

 

[Baluncore]

why
(bin_index/ (N-1)) * fs?

It will depend on your array base indexing.
If you get that wrong it will introduce a frequency dependent phase error.

 

[olivermsun]

I don't quite understand the last part. If he is passing in N points and receiving N FFT bins then how can the indexing be spread over N-1 frequencies?

 

[Baluncore]

Because DC is zero frequency in bin zero.

 

[olivermsun]

If you have N bins and divide by (N-1), then don't you get two bins with the same normalized frequency?

I think what you want is N bins with normalized frequency [0:N-1]/N, right?

 

[btb4198]

when you say DC do you mean Direct current?

 

[olivermsun]

Yes, that's why he said DC is zero frequency.

 

[btb4198]

ok i just did this:

1.0
-0.41614683654714241
-0.65364362086361194
0.960170286650366
-0.14550003380861354
-0.83907152907645244
0.84385395873249214
0.13673721820783361

those 8 number from
Cos(2 * n);
and that go into my fft and I get :

{(0.886399443294871, 0)}
{(0.862297845278363, 0.616189555028432)}
{(0.664289628322506, 2.35212587048179)}
{(1.42870222233886, -2.37880560416378)}
{(1.20302116482566, 0)}
{(1.42870222233886, 2.37880560416378)}
{(0.664289628322506, -2.35212587048179)}
{(0.862297845278363, -0.616189555028432)}

then I reverse and send it back to the fft and dived by N
and I get this :

0.99999999999999989
0.41614683654714235
0.65364362086361194
0.96017028665036586
0.14550003380861354
0.83907152907645233
0.84385395873249214
0.13673721820783363

I lost my - value , will they all are + now... but it is because I am taking the
Magnitude
I am doing

double[] returntemp = new double[R];
           for (int v = 0; v < R; v++)
           {
               F[v] = F[v].Magnitude / N;
               returntemp[v] = F[v].Magnitude;
           }
           return returntemp;

I did not know if I should use the real or imaginary numbers, so I just used the
Magnitude
is that ok ?

 

[olivermsun]

If you want to verify that your FFT and IFFT are working correctly, you need to output both the real and imaginary parts. There will probably be a small imaginary part due to numerical errors but it should be negligible.

Your magnitudes look good.

 

[btb4198]

ok this is my real:

0.99999999999999989
-0.41614683654714235
-0.65364362086361194
0.96017028665036586
-0.14550003380861354
-0.83907152907645233
0.84385395873249214
0.13673721820783363

so I should use my real and not the Magnitude?

 

[btb4198]

my img is all like 0 and 0.00000000000000017478323165292
which is 0

 

[btb4198]

here are some video i make




see my fft does work for that site with the
(bin index/N) * fs * 2

also I did test it again with the new version of this site and it still work

 

[olivermsun]

Then that's good. Your original signal was all real, so it makes sense that any imaginary part that you get after your IFFT is just numerical error. The real part is the part you want to keep.

 

[olivermsun]

here are some video i make




see my fft does work for that site with the
(bin index/N) * fs * 2

also I did test it again with the new version of this site and it still work

Did you ever figure out if your samples contain interleaved L+R stereo channels?

 

[btb4198]

yes it is Stereo Channels.
I am read in Stereo data from the sound card

so my FFT and IFFT it working now...
thanks for all the help

but question about the *2

as you can see in the youtube video my fft is working the * 2

you think this has something to do with the fact I am using Stereo data ?

 

[olivermsun]

I don't know how your library is extracting the stereo data, but sure it could explain the 2* factor. Interleaving the L+R channels would mean you have twice as many samples per second as you might think. If the signals look more or less the same, the interleaved data could look similar to a (noisy) sine wave of half the frequency.

 

[btb4198]

yes I think that is what it is doing and that is how it plays that will... that is why I can add 2*PI to play it.

ok I make a low pass filter:

       public void lowpass_filter( int frequency)
        {
 for (int i = 0; i < R; i++)
 {      
     if((i / N) * Fs >= frequency)
     {
         F[i] = 0D;
     }

 }  }

and it is working good from what I have been testing ...
but my high pass is not working?

public void hightpass_filter(int frequency)
        {
            for (int i = 0; i < R; i++)
            {
                if ((i / N) * Fs <= frequency)
                {
                    F[i] = 0D;
                }

            }

I still get both frequencies back
high and low
frequencies back


the lowpass filter is working so good
why it the high pass filter not?

 

[olivermsun]

The problem is what I mentioned to you earlier — the frequencies in an FFT output appear in pairs. So instead of thinking of the frequencies as [0, 1, 2, …, N-1]/N, you need to think of them as [0, 1, …, N/2-1, -N/2, -N/2+1, …, -2, -1]/N. You always need to filter positive and negative frequency bins in pairs.

The lowpass kind of works because you are zeroing out all the bins higher than a certain number. However, you won't get back a pure real signal because all you have left are unpaired (low numbered) bins. The high pass doesn't work at all because you are removing a few low numbered bins, but the corresponding negative frequency bins are untouched.

Also, if you are dealing with stereo channels, what you really need to do is de-interleave the stream and then FFT/filter the two channels separately.

 

[btb4198]

ok how do I find the negative frequency bins?
I tried to google it, but I do not see anything ...
can you show me an example ?
I get that there are two part of a signal way, one - and +. but what I do not get is how with the high pass filter I can till getting the low Frequency.
I mean you use bin Index to get the frequencies right and and every bin is different so no two bins should get you the same Frequency right?

 

[olivermsun]

look carefully at the [0, 1, 2, …, N-1]/N I posted above.

Let's say you have 8 samples at sampling rate fs = 8 Hz. The "normal" way you'd label your frequency bins is
[0, 1, 2, …, N-1]/N * fs = [0, 1, 2, …, 7]/8 * 8 Hz = [0, 1, 2, …, 7] Hz.

Since the frequencies wrap around the 8 Hz sampling frequency, these are exactly the same as
[0, 1, 2, …, N/2-1, -N/2, …,-2, -1,]/N * fs = [0, 1, 2, 3, -4, -3, -2, -1]/8 * 8 Hz = [0, 1, 2, 3, -4, -3, -2, -1] Hz.

If you do a quick test, you'll notice that the Fourier coefficients contained in the 1 and -1 Hz bins are just complex conjugates of each other. Similarly for 3, -3 Hz and every other pair of bins. This is because you have a real-valued signal. If you want to so something to your signal but keep it real-valued, you have to maintain this symmetry.

For example, to high pass keeping >= 2 Hz signals, you can zero the 0 Hz and +1, -1 Hz bins and then inverse transform.

 

[btb4198]

ok based off your last post I make this function :

 public void hightpass_filter(int frequency)
        {
          int i = 0;
          int t = 0;
            for (; i < R/2; i++)
            {
                if ((i / N) * Fs <= frequency)
                {
                    F[i] = 0D;
                }

            }
            for (; i < R; i++)
            {
                if ((t / N) * Fs <= frequency)
                {
                    F[i] = 0D;
                }

                t++;
            }
      
      }

but it is this not working...
I still see both Frequencies
so what did I do wrong ?
from 0 to N/2 is positive
and from N/2 to (N-1) is negative
right?
or did I not understand that right ?

 

[olivermsun]

Try running the example from my post through this function. Just use double F[] = {0, 1, 2, 3, -4, -3, -2, -1} or something easily identifiable to fill the bins. I think you'll discover the issue.

 

[btb4198]

what Fs should i use?
because I am using 44100 and that is too high for only 8 points
after 1/N * 44100 is too high...

(never mind)

 

[olivermsun]

Use any Fs for your test. You're just testing with synthetic data that you fully understand so that you can verify your code. I think my example above used 8 points at 8 Hz.

 

[btb4198]

ok I did what you said, but I am not see why it is now working:
ok I am using:
double[] list = { 0, 1, 2, 3, -4, -3, -2, -1 };
and
and I set the high pass filter to 2
and after doing it

I get
{(0, 0)}
{(0, 0)}
{(0, 0)}
{(4, -1.65685424949238)}
{(0, 0)}
{(0, 0)}
{(0, 0)}
{(4, 9.65685424949238)}

oh I am using fs = 8;

this is what I get before the filter:
1 10.452503719011
2 5.65685424949238
3 4.32956880116958

and this is what I get after the filter:
1 5.22625185950551
3 2.16478440058479

the 2 is done, but the 1 is this there
how 1 hz is bin 1 and bin 6 right?
are am I wrong about this ?

I also tried with Fs = 44100
I tried to filter out anything lower then 11000
i got:
{(0, 0)}
{(0, 0)}
{(-4, 4)}
{(4, -1.65685424949238)}
{(0, 0)}
{(0, 0)}
{(-4, -4)}
{(4, 9.65685424949238)}
this is before the filter:

5512.5 10.452503719011
11025 5.65685424949238
16537.5 4.32956880116958

and after the filter:
5512.5 5.22625185950551
11025 5.65685424949238
16537.5 2.16478440058479

5512.5 is divide by half
...
what i am missing here ?

 

[olivermsun]

Look carefully. If your FFT output is { 0, 1, 2, 3, -4, -3, -2, -1 } and your frequencies are [ 0, 1, 2, 3, -4, -3, -2, -1 ] Hz, then a high pass that keeps only frequencies greater than 2 Hz should leave you with { 0, 0, 0, 3, -4, -3, 0, 0 }, which you then pass to your IFFT.

 

[btb4198]

hey,
i think it is working!
thanks,

 

[olivermsun]

Great!
You're welcome.