basic FFT questions

[clepsydrae]

So, trying to wrap my head around FFT.

A basic DFT takes n samples and converts to n complex result values, where the frequency represented by each slot is given by n/windowsize*srate, and with an audio input signal it only makes sense to look up to result slot n/2 as the values above would not be useful... right?

So, do I understand correctly that a series of real values as input (e.g. audio signal) still generates a series of complex values on output, but that in audio we don't generally care about the upper half of the output values?

If so, is it the case that the JS function that reaper provides -- fft() -- takes n real samples and produces n/2 pairs of result values in the same buffer, where the first is the real component and the second is the imaginary?

Thus after calling fft_permute, the 0'th pair (indices 0 and 1) is <dc_component, 0> and thereafter the i'th pair (indices 2*i and 2*i+1) where i goes from 1 to windowsize/2 - 1 is the <real, imaginary> FFT value of the signal at frequency i/windowsize*srate ? No fencepost errors there?

And magnitude would be sqrt(results[2*i]^2 + results[2*i+1]^2), and would range from 0 to 1?

I have read that phase is inverse_tan(imag/real), but in the JS I see

ang= -atan2(aa,bb);

where aa is the first results value and bb is the second... i assume the negation is just there for screen drawing purposes, but I would have expected (bb,aa) -- maybe fft() returns the values as <imaginary, real> ? Actually, it's looking increasingly like a bug in the JS. Update: yep, pretty sure; afaict -atan2(bb,aa) would be correct there. And yes, the - sign is just for drawing.

Last one: arctan only ranges from -pi/2 to pi/2, but you'd need a total of 2*pi range (+/- 1*pi) to express all possible phase shifting... I would guess that when the real part is < 0 I need to add pi to the result to get the actual phase? So that would properly be something like:

ang = atan2(bb,aa) + (aa < 0 ? 3.1415926 : 0);

Yeah?

Thanks in advance for any help!

[Tale]

I don't have time to answer all your questions right now, but I this explains atan2:

https://en.wikipedia.org/wiki/Atan2

[clepsydrae]

Originally Posted by Tale

I don't have time to answer all your questions right now, but I this explains atan2:

https://en.wikipedia.org/wiki/Atan2

Ah, got it. Unfortunately the JS documentation only says:

atan2(x,y) -- returns the Arc Tangent of x divided by y (return value is in radians).

...so I was confused. I figured it was just some kind of obscure efficiency thing. Good to know about atan2.

I will wait patiently for any other time you have for my other queries. Thanks!

[Tale]

Originally Posted by clepsydrae

If so, is it the case that the JS function that reaper provides -- fft() -- takes n real samples and produces n/2 pairs of result values in the same buffer, where the first is the real component and the second is the imaginary?

JSFX's fft(buf, n) takes n complex pairs as input, and also produces n complex pairs. If you want to input a real signal, then you need to convert it to complex first by setting the imaginary part to 0:

i = 0;
loop(n,
  complex_buf[2*i] = real_buf[i];
  complex_buf[2*i+1] = 0;
  i += 1;
);

fft(complex_buf, n);
fft_permute(complex_buf, n);

Originally Posted by clepsydrae

Thus after calling fft_permute, the 0'th pair (indices 0 and 1) is <dc_component, 0> and thereafter the i'th pair (indices 2*i and 2*i+1) where i goes from 1 to windowsize/2 - 1 is the <real, imaginary> FFT value of the signal at frequency i/windowsize*srate ? No fencepost errors there?

Yeah, that sounds about right.

Originally Posted by clepsydrae

And magnitude would be sqrt(results[2*i]^2 + results[2*i+1]^2), and would range from 0 to 1?

Yes, except I think its range is 0 to +inf.

Originally Posted by clepsydrae

maybe fft() returns the values as <imaginary, real> ?

Nope, re=buf[2*i], im=buf[2*i+1].

[clepsydrae]

Originally Posted by Tale

JSFX's fft(buf, n) takes n complex pairs as input, and also produces n complex pairs. If you want to input a real signal, then you need to convert it to complex first by setting the imaginary part to 0

D'oh! I see that now in the JSFX, missed it before, thanks.

I see this, too "Note that fft()/ifft() must NOT cross at 65,536 item boundary, so be sure to specify the offset accordingly" -- is this referring to the memory space that JSFX have access to? So I want to make sure not to have the fft workspace cross over 0[65535*k] to 0[65536*k] ?

The web docs mention rfft() and irfft() as being the "real" versions, but they are not listed in the in-program docs and are undefined functions... so, mere fantasy?

Yes, except I think its range is 0 to +inf.

Ok, so, last (?) questions regard scaling fft results magnitude to actual 0-to-1 amplitude (which I can use to get to dBFS).

Various places on the web say that where mag = sqrt(real^2 + imag^2) and n is number of samples/bins, the amplitude is then 2*mag/n, but the in-program script docs in reaper say "Your inputs or outputs will need to be scaled down by 1/size, if used", so I assume that's just mag/n.

So I'm wondering if I need/want the 2 in there, since the reaper docs don't mention it?

Second thing is that the gfxanalyzer JSFX doesn't use n at all when determining the screen y value for drawing, even though it draws a dBFS-accurate curve.

Simplifying/reducing, it does:

floor=-120;
scale=(gfx_h-20)*20/(-floor * log(10));
dv=aa*aa+bb*bb; // this is real/imag from fft
screen_y = (-log(dv)*0.5)*scale+20;

...so it would seem that it is not using n at all in determining the screen y value so... how does that work out?

Thanks for any ideas!

[Tale]

Originally Posted by clepsydrae

I see this, too "Note that fft()/ifft() must NOT cross at 65,536 item boundary, so be sure to specify the offset accordingly" -- is this referring to the memory space that JSFX have access to? So I want to make sure not to have the fft workspace cross over 0[65535*k] to 0[65536*k] ?

Yeah, so if you have buf of n complex pairs, and if buf < 65536, then buf + 2*n must not be >= 65536 (and if buf < 2*65536, then buf + 2*n must not be >= 2*65536, etc.).

Originally Posted by clepsydrae

The web docs mention rfft() and irfft() as being the "real" versions, but they are not listed in the in-program docs and are undefined functions... so, mere fantasy?

In WDL (the open-source C++ library that presumably provides the JSFX FFT) the real FFT is commented out because it doesn't "work right". I have read that taking the real FFT can be more efficient (I guess because you don't have to convert to complex first), but the consensus seems to be that the improvement is rather minor, which is probably why is was removed from JSFX.

Originally Posted by clepsydrae

Various places on the web say that where mag = sqrt(real^2 + imag^2) and n is number of samples/bins, the amplitude is then 2*mag/n, but the in-program script docs in reaper say "Your inputs or outputs will need to be scaled down by 1/size, if used", so I assume that's just mag/n.

I think 2*mag/n is correct (i.e. so 1 == 0dB), except for DC, which is mag/n.

Originally Posted by clepsydrae

So I'm wondering if I need/want the 2 in there, since the reaper docs don't mention it?

Well, the scaling is more of a general FFT thing than a JSFX thing, i.e. scaling is a "natural" consequence of the FFT itself.

Originally Posted by clepsydrae

...so it would seem that it is not using n at all in determining the screen y value so... how does that work out?

I think gfxanalyzer prescales the window by 1/fftsize.

[clepsydrae]

Originally Posted by Tale

I think 2*mag/n is correct (i.e. so 1 == 0dB), except for DC, which is mag/n.

Finally got to testing it... turns out for reaper's fft, you don't divide by n because it's baked in the result already. So:

absvalue = 2*sqrt(real^2 + imag^2)

[clepsydrae]

Originally Posted by clepsydrae

I have read that phase is inverse_tan(imag/real), but in the JS I see

ang= -atan2(aa,bb);

where aa is the first results value and bb is the second... i assume the negation is just there for screen drawing purposes, but I would have expected (bb,aa) -- maybe fft() returns the values as <imaginary, real> ? Actually, it's looking increasingly like a bug in the JS.

Just following up that my conviction has grown that this is a bug in the gfxanalyzer JS. (Not that anyone really cares to see the phase of each bucket in an FFT result, but it was healthy for me to sort this out.) As far as I know, it should be atan2(imag, real) and that would properly yield "-atan2(bb,aa)" in the code (negated just because it's using the result to draw). (The only result of this bug in this JS is that the phase is shown shifted.)

[clepsydrae]

And an update/clarification in case anyone ever finds this thread: I missed the fact that gfxanalyzer incorporates the 1/fftsize in the windowing function, as you apparently need to do because when using a non-unity ("rectangular") window you don't scale by 1/fftsize but a different value dependent on the window (one over the sum of the window values.) You apply the scaling (or the scaling+window) on the sample values before the fft. Then when you do the fft you can calculate the magnitude (in absolute scale, 0 to 1) as shown above.

Interestingly, though, you don't need to scale back up when coming the other way. Meaning, you put your sample values in the fft workspace array (paired with 0-value complex values) at 1/fftsize scale, do the fft, do fft_permute, do whatever you want to to the data, then fft_ipermute, then ifft, and the resulting real values in the workspace do not need to be scaled up by fftsize, they can just be used as-is. This quirk led to a lot of confusion on my part.

Someone correct me if I've missed something else. :-)

[Justin]

Fixing some of the things pointed out here for the next builds, thanks!

[clepsydrae]

Originally Posted by Justin

Fixing some of the things pointed out here for the next builds, thanks!

Sweet, thanks!

I presume you're aware of this thread, too: http://forum.cockos.com/showthread.php?t=168666 ... just makin' sure. :-)

[Justin]

Originally Posted by clepsydrae

Sweet, thanks!

I presume you're aware of this thread, too: http://forum.cockos.com/showthread.php?t=168666 ... just makin' sure. :-)

Yep, thanks!