PSN-L Email List Message

Subject: FFT bin magnitude question
From: Randall Peters PETERS_RD@..........
Date: Wed, 11 Jul 2012 13:35:20 -0400


Randy,
      The factor of 16 difference that you mention is due to the fact that =
FFT algorithms do not 'normalize' their output.  A proper normalization req=
uires a 'sum' over the frequency domain (spectral) components to yield the =
same result as the average of the square of the time set values used to cal=
culate the Fourier transform.  Formally, this mathematical equivalence betw=
een the time domain and the frequency domain must be satisfied if Parseval'=
s theorem is to be obeyed.  In practice, I have never found a Cooley-Tukey =
(Fast Fourier) transform type for which this is done automatically.  The re=
ason involves the code-writers becoming sloppy so as to minimize memory req=
uirements placed on the computer.  In Excel's FFT calculation (Microsoft ca=
lls it 'Fourier Analysis' in their data analysis package (requiring 'add-in=
' to the 'tool-pak'), the Parseval requirement is satisfied only after  div=
iding the square of their modulus for a given row-value of FFT output [calc=
ulated using IMABS( )^2 ]  by the square of the total number of points used=
 in calculating the transform.  For the usual single-sided presentation of =
the spectrum (ignoring non-physical negative frequencies that are part of t=
he math of Fourier's gift to the world; i.e., double sided, for 'egghead' m=
athematicians) individual modulus squared values get multiplied by a factor=
 of two.   On the other hand, in Mathematica, the same modulus squared time=
s two value is treated differently.  It is instead divided by the total num=
ber of FFT points raised to the first power rather than the second power.  =
The difference between these two approaches is undoubtedly due to their dif=
ferent style of matrix manipulation.  Neither is formally correct, because =
the output, if left uncorrected (by the factor of either N^2 or N) will nev=
er satisfy Parseval's theorem.
     I will probably meet with debate over this posting, since I have gone =
'round and round' for more than a decade with the professional seismologist=
s about their need to start working with a better set of power spectral den=
sity units.  To satisfy both Parseval's theorem and to also work with somet=
hing that makes formal physics sense in terms of power-there can be only on=
e set of PSD units that are reasonable in terms of the whole world's recogn=
ition of the proper unit of power (watts).
That set of units can only be W/kg/Hz which equals m^2/s^3/Hz and not the l=
ong-standing seismology and mechanical engineering tradition of m^2/s^4/Hz =
(or g^2/Hz).
Some of them more correctly label their graphs with the words acceleration =
spectral density.  When one computes the FFT of the output from a seismomet=
er (accelerometer in general) in the units corresponding to the only thing =
to which the instrument responds; i.e., acceleration; then the ASD is simpl=
y m^2/s^4/Hz.  In other words, the ASD is very close to the single sided Fo=
urier transform spectrum discussed above, if properly normalized to satisfy=
 Parseval's theorem.  One must be very, very careful, however, to recognize=
 another subtlety that most everybody ignores.  The units as generated are =
really m^2/s^4/bin width.  The bin width is never actually one Hz, and so a=
 multiplicative factor must be considered if one is to actually use m^2/s^4=
/Hz, properly computed.
     Please, anybody, before challenging me on this subject (as I have come=
 to expect nearly every time I have brought it up elsewhere)-take a look at=
 a paper that I wrote, titled "Tutorial on power spectral density calculati=
ons for mechanical oscillators (with an exhaustive discussion of units)".  =
It can be found and read online at
http://physics.mercer.edu/hpage/psd-tutorial/psd.html
Or just type "spectral density calculations"  without the tick marks of a l=
iteral search into Google and click on the 2nd entry of the first page of t=
he sites they've catalogued.
I have spent more than a decade of intense study of this subject, turning i=
t every which way but loose; and so be prepared for a spirited debate if yo=
u think I'm wrong.
Where the attention to every one of these details is absolutely imperative =
is when one wants to calculate the total seismic power of the Earth (at lea=
st within the pass-band of an instrument, by integrating over the PSD and t=
hen (through understanding of Newton's third law) estimate the Earth's powe=
r in watts through the product of this integrated result and the mass of th=
e Earth in kg (6 x 10^24).
    Randall

Randy, &nbs= p; 

     =  The factor of 16 difference that you mention is due to the fact that = FFT algorithms do not ‘normalize’ their output.  A proper = normalization requires a ‘sum’ over the frequency domain (spect= ral) components to yield the same result as the average of the square of th= e time set values used to calculate the Fourier transform.  Formally, = this mathematical equivalence between the time domain and the frequency dom= ain must be satisfied if Parseval’s theorem is to be obeyed.  In= practice, I have never found a Cooley-Tukey (Fast Fourier) transform type = for which this is done automatically.  The reason involves the code-wr= iters becoming sloppy so as to minimize memory requirements placed on the c= omputer.  In Excel’s FFT calculation (Microsoft calls it ‘= Fourier Analysis’ in their data analysis package (requiring ‘ad= d-in’ to the ‘tool-pak’), the Parseval requirement is sat= isfied only after  dividing the square of their modulus for a given ro= w-value of FFT output [calculated using IMABS( )^2 ]  by the square of= the total number of points used in calculating the transform.  For th= e usual single-sided presentation of the spectrum (ignoring non-physical ne= gative frequencies that are part of the math of Fourier’s gift to the= world; i.e., double sided, for ‘egghead’ mathematicians) indiv= idual modulus squared values get multiplied by a factor of two.   = ;On the other hand, in Mathematica, the same modulus squared times two valu= e is treated differently.  It is instead divided by the total number o= f FFT points raised to the first power rather than the second power.  = The difference between these two approaches is undoubtedly due to their dif= ferent style of matrix manipulation.  Neither is formally correct, bec= ause the output, if left uncorrected (by the factor of either N^2 or N) wil= l never satisfy Parseval’s theorem.

     I will probably meet with debate over this posti= ng, since I have gone ‘round and round’ for more than a decade = with the professional seismologists about their need to start working with = a better set of power spectral density units.  To satisfy both Parseva= l’s theorem and to also work with something that makes formal physics= sense in terms of power—there can be only one set of PSD units that = are reasonable in terms of the whole world’s recognition of the prope= r unit of power (watts).

That set of uni= ts can only be W/kg/Hz which equals m^2/s^3/Hz and not the long-standing se= ismology and mechanical engineering tradition of m^2/s^4/Hz (or g^2/Hz).

Some of them more correctly label their gr= aphs with the words acceleration spectral density.  When one computes = the FFT of the output from a seismometer (accelerometer in general) in the = units corresponding to the only thing to which the instrument responds; i.e= .., acceleration; then the ASD is simply m^2/s^4/Hz.  In other words, t= he ASD is very close to the single sided Fourier transform spectrum discuss= ed above, if properly normalized to satisfy Parseval’s theorem. = One must be very, very careful, however, to recognize another subtlety tha= t most everybody ignores.  The units as generated are really m^2/s^4/b= in width.  The bin width is never actually one Hz, and so a multiplica= tive factor must be considered if one is to actually use m^2/s^4/Hz, proper= ly computed.    

 &n= bsp;   Please, anybody, before challenging me on this subjec= t (as I have come to expect nearly every time I have brought it up elsewher= e)—take a look at a paper that I wrote, titled “Tutorial on pow= er spectral density calculations for mechanical oscillators (with an exhaus= tive discussion of units)”.  It can be found and read online at<= o:p>

http://physics.mercer.edu/hpage/psd-tutorial/psd.= html

Or just type “spectral de= nsity calculations”  without the tick marks of a literal search = into Google and click on the 2nd entry of the first page of the = sites they’ve catalogued.

I have = spent more than a decade of intense study of this subject, turning it every= which way but loose; and so be prepared for a spirited debate if you think= I’m wrong.

Where the attention t= o every one of these details is absolutely imperative is when one wants to = calculate the total seismic power of the Earth (at least within the pass-ba= nd of an instrument, by integrating over the PSD and then (through understa= nding of Newton’s third law) estimate the Earth’s power in watt= s through the product of this integrated result and the mass of the Earth i= n kg (6 x 10^24). 

  &nb= sp; Randall

 

=

[ Top ] [ Back ] [ Home Page ]