Old 12-06-2015, 07:17 PM   #1
clepsydrae
Human being with feelings
 
clepsydrae's Avatar
 
Join Date: Nov 2011
Posts: 1,759
Default basic FFT questions

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!

Last edited by clepsydrae; 12-20-2015 at 07:24 PM.
clepsydrae is offline   Reply With Quote
Old 12-07-2015, 12:09 AM   #2
Tale
Human being with feelings
 
Tale's Avatar
 
Join Date: Jul 2008
Location: Holland
Posts: 2,761
Default

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

https://en.wikipedia.org/wiki/Atan2
Tale is offline   Reply With Quote
Old 12-07-2015, 02:14 PM   #3
clepsydrae
Human being with feelings
 
clepsydrae's Avatar
 
Join Date: Nov 2011
Posts: 1,759
Default

Quote:
Originally Posted by Tale View Post
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!
clepsydrae is offline   Reply With Quote
Old 12-07-2015, 02:46 PM   #4
Tale
Human being with feelings
 
Tale's Avatar
 
Join Date: Jul 2008
Location: Holland
Posts: 2,761
Default

Quote:
Originally Posted by clepsydrae View Post
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:

Code:
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);
Quote:
Originally Posted by clepsydrae View Post
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.

Quote:
Originally Posted by clepsydrae View Post
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.

Quote:
Originally Posted by clepsydrae View Post
maybe fft() returns the values as <imaginary, real> ?
Nope, re=buf[2*i], im=buf[2*i+1].
Tale is offline   Reply With Quote
Old 12-07-2015, 10:55 PM   #5
clepsydrae
Human being with feelings
 
clepsydrae's Avatar
 
Join Date: Nov 2011
Posts: 1,759
Default

Quote:
Originally Posted by Tale View Post
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?

Quote:
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!
clepsydrae is offline   Reply With Quote
Old 12-08-2015, 02:09 AM   #6
Tale
Human being with feelings
 
Tale's Avatar
 
Join Date: Jul 2008
Location: Holland
Posts: 2,761
Default

Quote:
Originally Posted by clepsydrae View Post
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.).

Quote:
Originally Posted by clepsydrae View Post
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.

Quote:
Originally Posted by clepsydrae View Post
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.

Quote:
Originally Posted by clepsydrae View Post
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.

Quote:
Originally Posted by clepsydrae View Post
...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.
Tale is offline   Reply With Quote
Old 12-19-2015, 01:50 AM   #7
clepsydrae
Human being with feelings
 
clepsydrae's Avatar
 
Join Date: Nov 2011
Posts: 1,759
Default

Quote:
Originally Posted by Tale View Post
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 is offline   Reply With Quote
Old 12-20-2015, 07:22 PM   #8
clepsydrae
Human being with feelings
 
clepsydrae's Avatar
 
Join Date: Nov 2011
Posts: 1,759
Default

Quote:
Originally Posted by clepsydrae View Post
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 is offline   Reply With Quote
Old 12-28-2015, 11:55 PM   #9
clepsydrae
Human being with feelings
 
clepsydrae's Avatar
 
Join Date: Nov 2011
Posts: 1,759
Default

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. :-)
clepsydrae is offline   Reply With Quote
Old 12-29-2015, 07:26 PM   #10
Justin
Administrator
 
Justin's Avatar
 
Join Date: Jan 2005
Location: NYC
Posts: 10,051
Default

Fixing some of the things pointed out here for the next builds, thanks!
Justin is online now   Reply With Quote
Old 12-29-2015, 07:44 PM   #11
clepsydrae
Human being with feelings
 
clepsydrae's Avatar
 
Join Date: Nov 2011
Posts: 1,759
Default

Quote:
Originally Posted by Justin View Post
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. :-)
clepsydrae is offline   Reply With Quote
Old 12-29-2015, 08:24 PM   #12
Justin
Administrator
 
Justin's Avatar
 
Join Date: Jan 2005
Location: NYC
Posts: 10,051
Default

Quote:
Originally Posted by clepsydrae View Post
Sweet, thanks!

I presume you're aware of this thread, too: http://forum.cockos.com/showthread.php?t=168666 ... just makin' sure. :-)
Yep, thanks!
Justin is online now   Reply With Quote
Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off

Forum Jump


All times are GMT -7. The time now is 07:44 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2017, vBulletin Solutions, Inc.