[Numpy-discussion] rfft different in numpy vs scipy
Charles R Harris
charlesr.harris at gmail.com
Thu Sep 7 18:03:59 CDT 2006
On 9/7/06, Andrew Jaffe <a.h.jaffe at gmail.com> wrote:
>
> Hi all,
>
> It seems that scipy and numpy define rfft differently.
>
> numpy returns n/2+1 complex numbers (so the first and last numbers are
> actually real) with the frequencies equivalent to the positive part of
> the fftfreq, whereas scipy returns n real numbers with the frequencies
> as in rfftfreq (i.e., two real numbers at the same frequency, except for
> the highest and lowest) [All of the above for even n; but the difference
> between numpy and scipy remains for odd n.]
>
> I think the numpy behavior makes more sense, as it doesn't require any
> unpacking after the fact, at the expense of a tiny amount of wasted
> space. But would this in fact require scipy doing extra work from
> whatever the 'native' real_fft (fftw, I assume) produces?
>
> Anyone else have an opinion?
Yes, I prefer the scipy version because the result is actually a complex
array and can immediately be use as the coefficients of the fft for
frequencies <= Nyquist. I suspect, without checking, that what you get in
numpy is a real array with f[0] == zero frequency, f[1] + 1j* f[2] as the
coefficient of the second frequency, etc. This makes it difficult to
multiply by a complex transfer function or phase shift the result to rotate
the original points by some fractional amount.
As to unpacking, for some algorithms the two real coefficients are packed
into the real and complex parts of the zero frequency so all that is needed
is an extra complex slot at the end. Other algorithms produce what you
describe. I just think it is more convenient to think of the real fft as an
efficient complex fft that only computes the coefficients <= Nyquist because
Hermitean symmetry determines the rest.
Chuck
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://projects.scipy.org/pipermail/numpy-discussion/attachments/20060907/08e5a4f5/attachment.html
More information about the Numpy-discussion
mailing list