[Numpy-discussion] rfft different in numpy vs scipy

[Andrew Jaffe]

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
> discussion.
> 
> 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
> applications.
> 
> 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?

Yours,

Andrew