Hi Charles,

Charles R Harris wrote:
> On 9/7/06, *Andrew Jaffe* <a.h.jaffe@gmail.com
> <mailto:a.h.jaffe@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.

Unless I misunderstand, I think you've got it backwards:
   - numpy.fft.rfft produces the correct complex array (with no strange
packing), with frequencies (0, 1, 2, ... n/2)

   - scipy.fftpack.rfft produces a single real array, in the correct
order, but with frequencies (0, 1, 1, 2, 2, ..., n/2) -- as given by
scipy.fftpack's rfftfreq function.

So I think you prefer numpy, not scipy.