Re: [Matplotlib-users] Possible error in psd, specgram and company

[Eric F. <ef...@ha...>]

On 2013/11/21 8:45 AM, Eric Larson wrote:
> Hey Luis,
>
> I think in general the application of normalizations in DFT / FFT (or
> even the continuous transforms, for that matter) can vary from
> implementation to another, in terms of whether they're applied during
> the forward transform, inverse transform, or half on each. That sounds
> like what you're running into here. MATLAB and other packages (FFTW, and
> apparently whatever matplotlib uses internally for FFTs) tend to
> normalize during the inverse transform, so you would need to put the
> normalization in yourself to see Parseval's theorem work out, which is
> essentially what you're showing.
>

It is not a matter of how the FFT is normalized, but rather of whether 
one is calculating a power spectrum, or power spectral density (psd).

For example, see:

http://pubman.mpdl.mpg.de/pubman/item/escidoc:152164:1/component/escidoc:152163/395068.pdf

> To my mind, this isn't a bug, so much as something users should be made
> aware of so they understand the issues at stake. It would be good to add
> a note to the PSD computation about these normalization issues, assuming
> one doesn't already exist. I have not contributed to matplotlib before,
> but I'm happy to take a stab at this if people think it would be
> helpful. My guess is that matplotlib won't be changed to include the
> suggested normalization by default (but maybe it could be added as an
> option?) because this would break backward-compatibility for users that
> have taken this factor into account in their own code.

I think the mpl psd is consistent with other psd implementations and 
with the definition of the psd.

Eric

>
> Eric
>
>
>
> On Thu, Nov 21, 2013 at 10:30 AM, Luis Miguel García-Cuevas González
> <lui...@gm... <mailto:lui...@gm...>> wrote:
>
>     Hello!
>
>     I have a doubt about the way matplotlib computes (and plots) a psd.
>       I have
>     the following example code (I've imported pylab for convenience here):
>
>        from pylab import *
>
>        A = 10. ** (10. / 20.) / sqrt(2.) * 2. # This amplitude should
>     give 10 dB
>        frequency = 5.
>        sampling_frequency = 1000.
>        t_max = 12.
>        t = arange(0., t_max, 1. / sampling_frequency)
>        x = A * cos(2. * pi * frequency * t)
>
>        Pxx, freqs = mlab.psd(x, Fs=sampling_frequency, NFFT=len(x),
>            window=mlab.window_none)
>
>        dB_actual = 20. * log10(rms_flat(x)) # dB_actual = 10., as expected
>        dB_psd = 10. * log10(Pxx.max()) # dB_psd = 20.791812460476248
>
>
>     As I understand, this should not happen this way.  dB_psd should be
>     equal to
>     10, as dB_actual.  I can get it to compute the correct value with this:
>
>        Pxx, freqs = mlab.psd(x / sqrt(t_max), Fs=sampling_frequency,
>     NFFT=len(x),
>            window=mlab.window_none)
>
>     (i.e., dividing the variable by the square root of the time). Is
>     this really
>     the expected behavior or this is a bug?
>
>     If it does matter, I'm using using matplotlib 1.3.1.
>
>     Thanks for your help!
>
>     --
>     Luis Miguel García-Cuevas González
>
>     ------------------------------------------------------------------------------
>     Shape the Mobile Experience: Free Subscription
>     Software experts and developers: Be at the forefront of tech innovation.
>     Intel(R) Software Adrenaline delivers strategic insight and
>     game-changing
>     conversations that shape the rapidly evolving mobile landscape. Sign
>     up now.
>     http://pubads.g.doubleclick.net/gampad/clk?id=63431311&iu=/4140/ostg.clktrk
>     _______________________________________________
>     Matplotlib-users mailing list
>     Mat...@li...
>     <mailto:Mat...@li...>
>     https://lists.sourceforge.net/lists/listinfo/matplotlib-users
>
>
>
>
> ------------------------------------------------------------------------------
> Shape the Mobile Experience: Free Subscription
> Software experts and developers: Be at the forefront of tech innovation.
> Intel(R) Software Adrenaline delivers strategic insight and game-changing
> conversations that shape the rapidly evolving mobile landscape. Sign up now.
> http://pubads.g.doubleclick.net/gampad/clk?id=63431311&iu=/4140/ostg.clktrk
>
>
>
> _______________________________________________
> Matplotlib-users mailing list
> Mat...@li...
> https://lists.sourceforge.net/lists/listinfo/matplotlib-users
>