[SciPy-User] fft not giving Hermitian output for real input

Aronne Merrelli aronne.merrelli@gmail....
Sat Feb 25 11:24:11 CST 2012


On Thu, Feb 23, 2012 at 8:15 PM, Jeff Alstott <jeffalstott@gmail.com> wrote:

> I have a particular data set, d0, which is real:
>
> In [72]: d0
> Out[72]:
> array([ 0.00907105,  0.0372916 ,  0.01402867, ..., -0.04779497,
>        -0.07054817, -0.0436582 ])
>
> It doesn't produce Hermitian output from numpy.fft.fft.
>
>
>    In [87]: y0 = fft.fft(d0)
>    In [88]: y0
>    Out[88]:
>    array([  7.77156117e-14 +0.00000000e+00j,
>             6.89226454e-13 -1.56319402e-13j,
>             1.72140080e-13 +0.00000000e+00j, ...,
>            -7.95807864e-13 -1.13686838e-13j,
>            -3.41060513e-13 +1.13686838e-13j,   1.25055521e-12
>    -3.41060513e-13j])
>    In [89]: y0[1].conj()==y0[-1]
>    Out[89]: False
>    In [90]: y0[1].conj()==y0[-2]
>    Out[90]: False
>
>
I might be missing something, but this seems like floating point rounding
error to me. I don't know what your d0 really looks like, but with some
random numbers,

In [1]: d0 = np.random.normal(0,1,1024)
In [2]: y0 = np.fft.fft(d0)
In [3]: np.allclose( y0[1:].conj(), y0[-1:0:-1] )
Out[3]: True
In [4]: np.max( (y0[1:].conj() - y0[-1:0:-1]).imag )
Out[4]: 4.3520742565306136e-13
In [5]: np.max( (y0[1:].conj() - y0[-1:0:-1]).real )
Out[5]: 2.957634137601417e-13

I tried this in MATLAB and the two halves are exactly equal. I would assume
this means that the MATLAB FFT implementation does something extra to
eliminate the rounding error for the special case of real input, whereas
the NumPy FFT does not take that step. I'm pretty sure the NumPy FFT is
just calling FFTPACK, so you'd need to check the implementation details
there.

Cheers,
Aronne
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20120225/40ffee20/attachment.html 


More information about the SciPy-User mailing list