Re: [Numpy-discussion] rfft different in numpy vs scipy

2006-09-11 Thread 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

-
Using Tomcat but need to do more? Need to support web services, security?
Get stuff done quickly with pre-integrated technology to make your job easier
Download IBM WebSphere Application Server v.1.0.1 based on Apache Geronimo
http://sel.as-us.falkag.net/sel?cmd=lnkkid=120709bid=263057dat=121642
___
Numpy-discussion mailing list
Numpy-discussion@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/numpy-discussion


Re: [Numpy-discussion] rfft different in numpy vs scipy

2006-09-09 Thread Steven G. Johnson
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.

(FFTW does actually support a different, packed format in its r2r
interface, but only for 1d FFTs.  The reason for this has do to with
advantages of that format for odd sizes.  We recommend that most users
employ the r2c interface, however; our r2c interface is generally
faster for even sizes.)

I hope this is helpful.

Cordially,
Steven G. Johnson


-
Using Tomcat but need to do more? Need to support web services, security?
Get stuff done quickly with pre-integrated technology to make your job easier
Download IBM WebSphere Application Server v.1.0.1 based on Apache Geronimo
http://sel.as-us.falkag.net/sel?cmd=lnkkid=120709bid=263057dat=121642
___
Numpy-discussion mailing list
Numpy-discussion@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/numpy-discussion


Re: [Numpy-discussion] rfft different in numpy vs scipy

2006-09-07 Thread Charles R Harris
On 9/7/06, Andrew Jaffe [EMAIL PROTECTED] 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 areactually real) with the frequencies equivalent to the positive part of
the fftfreq, whereas scipy returns n real numbers with the frequenciesas in rfftfreq (i.e., two real numbers at the same frequency, except forthe 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 anyunpacking after the fact, at the expense of a tiny amount of wastedspace. 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
-
Using Tomcat but need to do more? Need to support web services, security?
Get stuff done quickly with pre-integrated technology to make your job easier
Download IBM WebSphere Application Server v.1.0.1 based on Apache Geronimo
http://sel.as-us.falkag.net/sel?cmd=lnkkid=120709bid=263057dat=121642___
Numpy-discussion mailing list
Numpy-discussion@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/numpy-discussion


Re: [Numpy-discussion] rfft different in numpy vs scipy

2006-09-07 Thread Andrew Jaffe
Hi Charles,

Charles R Harris wrote:
 On 9/7/06, *Andrew Jaffe* [EMAIL PROTECTED] 
 mailto:[EMAIL PROTECTED] 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.

This is complicated by the fact that (I think) numpy.fft shows up as 
scipy.fft, so functions with the same name in scipy.fft and 
scipy.fftpack aren't actually the same!

Andrew








-
Using Tomcat but need to do more? Need to support web services, security?
Get stuff done quickly with pre-integrated technology to make your job easier
Download IBM WebSphere Application Server v.1.0.1 based on Apache Geronimo
http://sel.as-us.falkag.net/sel?cmd=lnkkid=120709bid=263057dat=121642
___
Numpy-discussion mailing list
Numpy-discussion@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/numpy-discussion


Re: [Numpy-discussion] rfft different in numpy vs scipy

2006-09-07 Thread Charles R Harris
On 9/7/06, Andrew Jaffe [EMAIL PROTECTED] wrote:
Hi Charles,Charles R Harris wrote: On 9/7/06, *Andrew Jaffe* [EMAIL PROTECTED] mailto:[EMAIL PROTECTED]
 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 strangepacking), 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 byscipy.fftpack's rfftfreq function.So I think you prefer numpy, not scipy.Ah, well then, yes. I prefer Numpy. IIRC, one way to get the Scipy ordering is to use the Hartley transform as the front to the real transform. And, now that you mention it, there was a big discussion about it in the scipy list way back when with yours truly pushing the complex form. I don't recall that anything was settled, rather the natural output of the algorithm they were using prevailed by default. Maybe you could write a front end that did the right thing?
 Chuck
-
Using Tomcat but need to do more? Need to support web services, security?
Get stuff done quickly with pre-integrated technology to make your job easier
Download IBM WebSphere Application Server v.1.0.1 based on Apache Geronimo
http://sel.as-us.falkag.net/sel?cmd=lnkkid=120709bid=263057dat=121642___
Numpy-discussion mailing list
Numpy-discussion@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/numpy-discussion