[Numpy-discussion] rfft different in numpy vs scipy
a.h.jaffe at bakerjaffe.plus.com
Mon Sep 11 10:02:08 CDT 2006
Steven G. Johnson wrote:
> Andrew Jaffe wrote:
>> 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?
> As an author of FFTW, let me interject a couple of points into this
> First, if you are using FFTW, then its real-input r2c routines
> "natively" produce output in the "unpacked" numpy format as described
> above: an array of n/2+1 complex numbers. Any "packed" format would
> require some data permutations. Other FFT implementations use a
> variety of formats.
> Second, the *reason* why FFTW's r2c routines produce unpacked output is
> largely because "packed" formats do not generalize well to
> multi-dimensional FFTs, while the "unpacked" format does. (Packed
> formats are *possible* for multidimensional transforms, but become
> increasingly intricate as you add more dimensions.) Additionally, I
> personally find the unpacked format more convenient in most
> I hope this is helpful.
OK -- so it appears that all (three) of the votes so far are in support
of the numpy convention -- a complex result.
If this is something that can be done in pure python, I'm willing to
give it a stab, but I'm probably not capable of handling any python/C
issues. Does anyone out there understand the interaction between fftpack
(C & python?), fftw (C), scipy and numpy in this context well enough to
give some advice?
More information about the Numpy-discussion