DSPRelated.com
Forums

FFT bin meaning

Started by hyeewang February 10, 2009
On Feb 11, 11:50�pm, kevinjmc...@netscape.net wrote:
> On Feb 11, 7:51&#4294967295;pm, robert bristow-johnson <r...@audioimagination.com> > wrote: >
...
> > > i'm not sure, but wouldn't it be the case that if you started with an > > infinitely-long continuous-time input, x(t), and multiplied it by a > > "window" that's the right sorta sinc() function > > > &#4294967295; &#4294967295; w(t) = sinc( (Fs/N)*t ) &#4294967295; &#4294967295;(sinc(u) = sin(pi*u)/(pi*u)) > > > now we have x(t)*w(t). &#4294967295;then it gets sampled and time-aliased to get a > > discrete, periodic function: > > > &#4294967295; &#4294967295; x[n] = SUM{ x((n-m*N)/Fs) * w((n-m*N)/Fs) } > > &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; m > > > and run a single period (N samples) of that x[n] into the DFT (note > > the less common scaling convention). > > > &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295;N-1 > > &#4294967295; &#4294967295; X[k] = (1/N) * SUM{ x[n] * exp(-j*2*pi*n*k/N) } > > &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295; &#4294967295;n=0 > > > wouldn't the value in bin #k, X[k] be a brick-wall filter of the sum > > of amplitudes > > > would not the bin #k contain only frequency components of x(t) that > > live strictly between &#4294967295; (k-1/2)*Fs/N &#4294967295;and &#4294967295;(k+1/2)*Fs/N ? &#4294967295;wouldn't > > that be the case? >
...
> > Let me say first to Fred and Glen: yes, you're right - the 'bin' idea > is certainly a heuristic, albeit a useful one for interpreting > things. &#4294967295;I wasn't sure of the OP's knowledge of windows/convolution so > I tried to give a simpler answer. > > As for robert - I think I undestand what you're tryng to get at, but > I've always interpreted spectral leakage as explained in the well > known Harris paper (you probably have a copy, but here's a link): > > http://web.mit.edu/xiphmont/Public/windows.pdf > > As shown, he explains spectral leakage in terms of the discontinuities > at the ends of the data record. &#4294967295;He also says that the easiest way to > diminish them is to match as many orders of the derivatives of the > weighted data to zero. > > I think that if you use a very large data record, and the right type > of window, you can come close to gettting a 'brick wall' filter. &#4294967295;But > I think that to get an actual brick wall would imply that there would > be no spectral leakage whatsoever. So if you have a single frequency > that is swept from one 'bin' to the next, then the output shows up in > the first bin when the frquency is closer to it that the next one, and > immediately disappears from the first 'bin' and instantly appears in > the next when the frequency of the input is closerr to that one than > the first. &#4294967295;I think that in order to do that, you'd have to match all > the derivatives to zero, and I've never seen a window that does that.
well, first, i have to confess that my "window" is sorta bogus. it's a sinc() function and the only thing it really has in common with a window function is that i'm multiplying r(t) by it. (i changed the name of "x(t)" to "r(t)", because i am running out of letters.) first of all, i am assuming to start that r(t) is sufficiently bandlimited. its F.T. is zero at and above Fs/2: R(f) = 0 for all |f| >= Fs/2 w(t) = sinc( (Fs/N)*t ) (sinc(u) = sin(pi*u)/(pi*u)) the F.T. of this is W(f) = (N/Fs)*rect((N/Fs)*f) where rect(u) = 1 for |u| < 1/2, 0 otherwize so now, just like with a legit window, when i multiply r(t) by w(t), i convolve the spectra R(f) and W(f) together: let v(t) = r(t) * w(t) then V(f) = R(f) (*) W(f) ( "(*)" means "convolve") this turns out to be f+Fs/(2N) V(f) = integral{ (N/Fs) * R(u) du} f-Fs/(2N) so that ends up smearing the local bunch of frequency components together. oh PISSER!! i just realized that this makes the V(f) spectrum a little bit wider so it isn't quite sample-able at rate Fs. so i will slightly modify the beginning requirement to make it a little more tight. let's bandlimit r(t) just a little more: R(f) = 0 for all |f| >= Fs/2 - Fs/(2N) so this means not only is there nothing in r(t) at or above Nyquist, but i don't want anything in there to be within 1/2 bin from Nyquist either. so now, after that convolution in the freq domain, we know that V(f) = 0 for all |f| >= Fs/2 and v(t) = r(t)*w(t) can be sampled at rate Fs. now v(t) stretches out forever and describing it with samples would require more samples than N. this reaches back to my own religious belief (that we occasionally fight over here on comp.dsp) that the DFT fundamentally transforms one periodic sequence of period N to another periodic sequence of period N. so i am gonna do something to v(t) that periodically extends it without increasing its bandwidth (before sampling): +inf x(t) = SUM{ v(t - m*N/Fs) } m=-inf so x(t) is periodic (and still bandlimited) with period N/Fs (N sample times). that ends up causing the spectrum X(f) to be sampled, having dirac delta functions at discrete frequencies of k*Fs/N up to, but not including Nyquist, +/- Fs/2. N/2-1 X(f) = SUM{ V(k*Fs/N) * delta(f - k*Fs/N) } k=-N/2+1 and V(k*Fs/N) is the average of all of R(f) (our original spectrum) from (k-1/2)*Fs/N to (k+1/2)*Fs/N . this means, going back to x(t), that N/2-1 x(t) = SUM{ V(k*Fs/N) * exp(j*2*pi*(k*Fs/N)*t) } k=-N/2+1 now i sample it: x[n] = x(n/Fs) and then we yank out one period of this periodic discrete sequence and send that to the DFT: N-1 X[k] = (1/N) * SUM{ x[n] * exp(-j*2*pi*n*k/N) } n=0 note the 1/N in the forward DFT (this is why i don't like the convention that leaves it out of the forward DFT). so the inverse DFT is N-1 x[n] = SUM{ X[k] * exp(+j*2*pi*n*k/N) } k=0 but, because of the periodicity of X[k], this is the same as (N/2)-1 x[n] = SUM{ X[k] * exp(+j*2*pi*n*k/N) } k=-(N/2) but, because of that additional bandlimiting we had to do at the beginning, we don't have anything at Nyquist: X[-N/2] = X[N/2] = 0 these non-zero components at Nyquist always fuck me up. i just don't know what to do with them. so now it's (N/2)-1 x[n] = SUM{ X[k] * exp(+j*2*pi*n*k/N) } k=-(N/2)+1 compare this to N/2-1 x(t) = SUM{ V(k*Fs/N) * exp(j*2*pi*(k*Fs/N)*t) } k=-N/2+1 or N/2-1 x(n/Fs) = SUM{ V(k*Fs/N) * exp(j*2*pi*k*n/N) } k=-N/2+1 but x(n/Fs) is x[n], so it looks like that X[k] = V(k*Fs/N) which is k*Fs/N + Fs/(2N) V(k*Fs/N) = integral{ (N/Fs) * R(f) df} k*Fs/N - Fs/(2N) or (k+1/2)*Fs/N X[k] = V(k*Fs/N) = (N/Fs)*integral{ R(f) df} (k-1/2)*Fs/N so do you see what i did? with my bogus window function and my time- aliasing (or periodic extension) of r(t)*w(t), i showed that once you sample that and DFT it, what goes into bin #k is the average of all frequency components of r(t) that are within +/- 1/2 bin width from the bin center of bin #k. it gives *some* justification to the concept of the bin collecting all of the components within its "reach". but it's pretty contrived, i will grant you. r b-j
On Feb 10, 10:39&#4294967295;pm, hyeewang <hyeew...@gmail.com> wrote:
> This means that the frequencies of all sinusoids we measure will be a > multiple of the inverse of the analysis window length - so if our > &#4294967295;nails&#4294967295; are N samples away, our STFT bins will have a spacing of > sampleRate/N Hertz. As a result, this concept imposes an artificial > frequency grid on our analysis by requiring the reference frequencies > to be an integer multiple of our signal window in period to make them > seamlessly fit into our analysis frame. > > I often see DFT bin in some DFT/stft literature,What does the term > "bin" mean exactly? > Does it refer to a specific single frequency or a frequency range? > If we do N points DFT,then have N DFT bins, where do the 0-bin and N-1 > bin refer to? > If N = 256,then where does 112.5 bin refer to? > > Thanks. > Hyee
a specific bin value will be a function of frequecy component inside and outside the bin, just more a function of the inside frequency components then the outside frequency components (depending on the window used), they all kind of get thrown into the "bin", so to speak
On Feb 12, 2:36&#4294967295;am, robert bristow-johnson <r...@audioimagination.com>
wrote:

To both glen and robert: great posts - you really got me thinking
about this.  First, let me correct an error when I wrote "a single
frequency that is swept from one 'bin' to the next..."  This would
obviously produce smearing, so I should have written "step the input
frequency by some small increment and recompute using the incremented
frequency."

Glen, I'll have a lot to say about your post (Gaussian windows,
minimum time-bandwidth product).  There's a great technique for
estimating tonal inputs based on the use of a Gaussian window, and it
gives a result that's almost as if you used brick wall filters.  I'll
refer to this a little later.

As for robert, your 'contrived' result may not be so contrived after
all.  Think of it this way.  Suppose you want 64 frequency 'bins' for
your input, and you did so by creating 64 band pass FIR filters to
process it. The transition zones of those filters would have to be
incredibly sharp, so you'd need an awful lot of taps.  Now it seems
that if you have an input that was periodic or nearly periodic over N
(as you said: sufficiently band limited), then the above would work.
But what really caught my eye in your post was "now v(t) stretches out
forever and describing it with samples would require more samples than
N."  That would seem to be similar to the use of large FIR filters,
where the impulse response would be much larger than N.

Doing the above with 64 large band pass FIR filters would require an
awful lot of calculation, but perhaps there's some saving to be had by
making the number of taps some multiple of N.  And all the band pass
filters would be uniformly separated in frequency.  I can't help but
think that someone somewhere must have looked at this before.  It's an
interesting thing to think about.

As for Glen's post. I've often used a Gaussian window.  There's a
great technique for estimating the frequency and amplitude of tonal
inputs using a Gaussian window:

http://www.edn.com/archives/1994/030394/05df1.htm

I've used it many times, and I've been impressed by the results.  As
noted in the above, when you pass a sinusoid of amplitude &#4294967295;A&#4294967295; and
frequency &#4294967295;F&#4294967295; through a Gaussian filter, you get an attenuated copy of
the signal that is a Gaussian function of frequency.  When we apply a
Gaussian window when doing an FFT, we are in fact creating the
equivalent of a series of Gaussian band pass filters in the frequency
domain, where the center of each filter is located at each one of the
FFT&#4294967295;s output points.  Now if we consider an input tone that lies
somewhere between 2 consecutive filters, the amplitude output of each
filter (call them a(F) and b(F)) are:

a(F) = Ae-(F - ndf)2/cdf2
b(F) = Ae-(F - (n + 1)df)2/cdf2

Where c is a constant, n and n + 1 are the nth and (n +1)th filter
(i.e.: windowed FFT output point), and df is the frequency bin spacing
of the FFT outputs. Now if we take the natural logarithm of the ratio a
(F)/b(F) and manipulate results, we get:

F = ndf + df/2 + (df/2)c( ln(b(F)) - ln(a(F)) )

[Note: the author in the above link has a negative last term because
he uses ln(a(F)) - ln(b(F)) in his equation, but ln(b(F)) - ln(a(F))
in his code].

So the frequency of a tone that lies between 2 output points is a
simple function of the difference between the natural logarithm of the
difference between adjacent FFT outputs.  The author of the link then
goes on to provide a C example using 3 tones: one of amplitude 1 at
10.5 cps, one of amplitude .01 at 16 cps, and one of amplitude .1 at
25.21 cps (the 1978 HARRIS paper used just the first two). The
estimates produced by the technique are quite good.

There are several caveats about using it, though.  As noted by the
author, the tone(s) in the input can't be rapidly changing, and the
input can't have too much noise, and there shouldn't be any tones very
close in frequency ...   He doesn't quantify the noise limit - you get
processing gain from the FFT, but I haven't done any experiments using
various noise levels to test where the technique starts to fail.

You also have to be careful with some of the numbers used in the
above.  For instance, the author used a sample rate and N such that
his 'df' was = 1.  And his equations and code don't quite match (the
code works fine, but he does some things at different places in the
code than when he shows the equations).  It should also be pointed out
that one needs to compute the normalized power spectrum and examine it
for peaks (no sense in computing an amplitude and frequency between
output bins if the only thing in it is noise).

The overall process is: window the input data, compute the FFT,
compute the normalized power spectrum, examine the power spectrum for
peaks (a peak indicates the presence of a tone), take the natural log
of the two FFT outputs that correspond to the location of the peak
(i.e.: if the power spectrum shows a peak between the 5th and 6th
output, then compute the frequency and amplitude of the tone using the
natural logs of the 5th and 6th outputs).

It works almost as if you have brick wall band pass filters in the
frequency domain, and it's useful for tonal inputs.  I'm not quite
sure if it'd be useful if you have some broadband component in the
input as well.

(p.s. re: steve's comment: yep, a 'bin' will actually contain
something from a large frequency range - the tails of the filter due
to the windowing function extend much futher on either side of the
center of the bin).

Kevin