www.digitalmars.com         C & C++   DMDScript  

digitalmars.D - Re: FFT Lib?

reply jtravs <jtravs gmail.com> writes:
dsimcha Wrote:

 I'm going to need an FFT library to perform some convolutions at some point
 soon.  Two absolute, non-negotiable requirements are that it be written in
 pure D and that it be Boost or compatibly (i.e. zlib or public domain)
 licensed.  I also prefer "simple and good enough" over "has every
 micro-optimization in the book but a PITA to maintain/modify/use", as long as
 it's at least a true fft as opposed to an O(N^2) DFT.  A few questions:
 
 1.  Does anyone already have such a lib?
 
 2.  If noone has one I'll probably either write my own from scratch  or port
 some code from C if I can find code that's under a suitable license and
 written with a "simple and good enough" philosophy rather than an "every tiny
 optimization in the book" philosophy.  Could anyone recommend one to port?
 

You should have a look at kissfft. It is very small (one short c file), appropriately licensed and supports arbitrary sizes quite efficiently. I have recently ported this to the go programming language in just a few hours work. I'm willing to help with a D version too, if you want, as I might also need this soon. kissfft can be found at: http://kissfft.sourceforge.net J
Jul 28 2010
parent reply dsimcha <dsimcha yahoo.com> writes:
== Quote from jtravs (jtravs gmail.com)'s article
 dsimcha Wrote:
 I'm going to need an FFT library to perform some convolutions at some point
 soon.  Two absolute, non-negotiable requirements are that it be written in
 pure D and that it be Boost or compatibly (i.e. zlib or public domain)
 licensed.  I also prefer "simple and good enough" over "has every
 micro-optimization in the book but a PITA to maintain/modify/use", as long as
 it's at least a true fft as opposed to an O(N^2) DFT.  A few questions:

 1.  Does anyone already have such a lib?

 2.  If noone has one I'll probably either write my own from scratch  or port
 some code from C if I can find code that's under a suitable license and
 written with a "simple and good enough" philosophy rather than an "every tiny
 optimization in the book" philosophy.  Could anyone recommend one to port?


 I have recently ported this to the go programming language in just a few hours

this soon.
 kissfft can be found at: http://kissfft.sourceforge.net
 J

This looks great except that the license requires binary attribution. I've emailed the author and asked him to waive the binary attribution requirement and let me use it under the Boost license for the purpose of a D port. If he goes for it, I'll probably just port this code to D for std.numeric, and then use it where I need to in dstats. Otherwise, it's back to square 1.
Jul 28 2010
next sibling parent reply Fawzi Mohamed <fawzi gmx.ch> writes:
On 29-lug-10, at 00:01, dsimcha wrote:

 == Quote from jtravs (jtravs gmail.com)'s article
 dsimcha Wrote:
 I'm going to need an FFT library to perform some convolutions at  
 some point
 soon.  Two absolute, non-negotiable requirements are that it be  
 written in
 pure D and that it be Boost or compatibly (i.e. zlib or public  
 domain)
 licensed.  I also prefer "simple and good enough" over "has every
 micro-optimization in the book but a PITA to maintain/modify/use",  
 as long as
 it's at least a true fft as opposed to an O(N^2) DFT.  A few  
 questions:

 1.  Does anyone already have such a lib?

 2.  If noone has one I'll probably either write my own from  
 scratch  or port
 some code from C if I can find code that's under a suitable  
 license and
 written with a "simple and good enough" philosophy rather than an  
 "every tiny
 optimization in the book" philosophy.  Could anyone recommend one  
 to port?

file),

 I have recently ported this to the go programming language in just  
 a few hours

might also need this soon.
 kissfft can be found at: http://kissfft.sourceforge.net
 J

This looks great except that the license requires binary attribution. I've emailed the author and asked him to waive the binary attribution requirement and let me use it under the Boost license for the purpose of a D port. If he goes for it, I'll probably just port this code to D for std.numeric, and then use it where I need to in dstats. Otherwise, it's back to square 1.

otherwise, as I already said http://www.netlib.org/fftpack/fft.c is public domain, has reasonable performance, and supports all sizes. It might not be the always beautiful, but it is stable and works well... Fawzi
Jul 29 2010
parent reply dsimcha <dsimcha yahoo.com> writes:
== Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 otherwise, as I already said http://www.netlib.org/fftpack/fft.c is
 public domain, has reasonable performance, and supports all sizes.
 It might not be the always beautiful, but it is stable and works well...
 Fawzi

Yeah, I just looked at this code. The problem is that I wouldn't do a straight, mechanical port to D, I'd D-ify the code a little while I was at it (such as making it work on generic random access ranges if it's an out of place transform, making real-only inputs vs. complex inputs DRYer with templates, etc.). Since this code seems to have been straight-up translated from FORTRAN, it scores about a -2 the readability scale, from 1 to 10. In other words, I would probably be able to translate it and make it run, but even the slightest modification would be near impossible.
Jul 29 2010
parent reply dsimcha <dsimcha yahoo.com> writes:
== Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 On 29-lug-10, at 15:31, dsimcha wrote:
 == Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 otherwise, as I already said http://www.netlib.org/fftpack/fft.c is
 public domain, has reasonable performance, and supports all sizes.
 It might not be the always beautiful, but it is stable and works
 well...
 Fawzi

Yeah, I just looked at this code. The problem is that I wouldn't do a straight, mechanical port to D, I'd D-ify the code a little while I was at it (such as making it work on generic random access ranges if it's an out of place transform, making real-only inputs vs. complex inputs DRYer with templates, etc.). Since this code seems to have been straight-up translated from FORTRAN, it scores about a -2 the readability scale, from 1 to 10. In other words, I would probably be able to translate it and make it run, but even the slightest modification would be near impossible.

the code too bad, I classified it as a typical numerical intensive code (even rather clean where one can easily find the boundaries of the arrays, and without too many goto), and templating float, using slices as input seemed easy.

Yes, but typical numerics intensive code is basically unreadable. I think D is about the only language that has enough high-level stuff (templates, etc.) to write readable numerics code with a nice interface, but is also efficient enough to write efficient numerics code. Unfortunately, we're basically starting from scratch because C and Fortran numerics code is so crufty and unreadable that it's really, really hard to refactor it into readable D or even put a nice interface on it.
 But yes probably it is not the best basis to develop your own thing.
 Ahrg fortran ruined me ;).
 ciao
 Fawzi

I got permission from Mark Borgerding to port kissfft and have the binary attribution clause waived, but unfortunately this code is almost as difficult to read due to its extreme preprocessor abuse. (Basically he implements complex numbers as a struct with an imaginary and real component and a bunch of preprocessor macros.) Part of the preprocessor abuse is necessary for supporting fixed point FFTs. Does anyone care about these (I sure don't), or can I just get rid of them? I'm hoping I can refactor this code into a nice D interface, but it may end up being that I have to basically start from scratch. One thing I will absolutely not do is contribute code without a proper D-style (i.e. templated, works with std.complex, etc.) interface. If I can't understand kissfft well enough to do this refactoring, then maybe I'll just start from scratch using the Wikipedia FFT article and a textbook.
Jul 30 2010
next sibling parent reply Don <nospam nospam.com> writes:
dsimcha wrote:
 == Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 On 29-lug-10, at 15:31, dsimcha wrote:
 == Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 otherwise, as I already said http://www.netlib.org/fftpack/fft.c is
 public domain, has reasonable performance, and supports all sizes.
 It might not be the always beautiful, but it is stable and works
 well...
 Fawzi

a straight, mechanical port to D, I'd D-ify the code a little while I was at it (such as making it work on generic random access ranges if it's an out of place transform, making real-only inputs vs. complex inputs DRYer with templates, etc.). Since this code seems to have been straight-up translated from FORTRAN, it scores about a -2 the readability scale, from 1 to 10. In other words, I would probably be able to translate it and make it run, but even the slightest modification would be near impossible.

the code too bad, I classified it as a typical numerical intensive code (even rather clean where one can easily find the boundaries of the arrays, and without too many goto), and templating float, using slices as input seemed easy.

Yes, but typical numerics intensive code is basically unreadable. I think D is about the only language that has enough high-level stuff (templates, etc.) to write readable numerics code with a nice interface, but is also efficient enough to write efficient numerics code. Unfortunately, we're basically starting from scratch because C and Fortran numerics code is so crufty and unreadable that it's really, really hard to refactor it into readable D or even put a nice interface on it.
 But yes probably it is not the best basis to develop your own thing.
 Ahrg fortran ruined me ;).
 ciao
 Fawzi

I got permission from Mark Borgerding to port kissfft and have the binary attribution clause waived, but unfortunately this code is almost as difficult to read due to its extreme preprocessor abuse. (Basically he implements complex numbers as a struct with an imaginary and real component and a bunch of preprocessor macros.) Part of the preprocessor abuse is necessary for supporting fixed point FFTs. Does anyone care about these (I sure don't), or can I just get rid of them? I'm hoping I can refactor this code into a nice D interface, but it may end up being that I have to basically start from scratch. One thing I will absolutely not do is contribute code without a proper D-style (i.e. templated, works with std.complex, etc.) interface. If I can't understand kissfft well enough to do this refactoring, then maybe I'll just start from scratch using the Wikipedia FFT article and a textbook.

This is a useful read on FFT performance. http://cnx.org/content/m16336/latest/
Jul 31 2010
parent bearophile <bearophileHUGS lycos.com> writes:
Don:
 This is a useful read on FFT performance.
 http://cnx.org/content/m16336/latest/

That looks good, but it's also a large amount of information. Register allocation (spilling strategy) is critical in not uber-optimized FFT routines. And GPUs can help :-) Bye, bearophile
Jul 31 2010
prev sibling parent "Mike James" <foo bar.com> writes:
"dsimcha" <dsimcha yahoo.com> wrote in message 
news:i2umov$2esa$1 digitalmars.com...
 == Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 On 29-lug-10, at 15:31, dsimcha wrote:
 == Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 otherwise, as I already said http://www.netlib.org/fftpack/fft.c is
 public domain, has reasonable performance, and supports all sizes.
 It might not be the always beautiful, but it is stable and works
 well...
 Fawzi

Yeah, I just looked at this code. The problem is that I wouldn't do a straight, mechanical port to D, I'd D-ify the code a little while I was at it (such as making it work on generic random access ranges if it's an out of place transform, making real-only inputs vs. complex inputs DRYer with templates, etc.). Since this code seems to have been straight-up translated from FORTRAN, it scores about a -2 the readability scale, from 1 to 10. In other words, I would probably be able to translate it and make it run, but even the slightest modification would be near impossible.

the code too bad, I classified it as a typical numerical intensive code (even rather clean where one can easily find the boundaries of the arrays, and without too many goto), and templating float, using slices as input seemed easy.

Yes, but typical numerics intensive code is basically unreadable. I think D is about the only language that has enough high-level stuff (templates, etc.) to write readable numerics code with a nice interface, but is also efficient enough to write efficient numerics code. Unfortunately, we're basically starting from scratch because C and Fortran numerics code is so crufty and unreadable that it's really, really hard to refactor it into readable D or even put a nice interface on it.
 But yes probably it is not the best basis to develop your own thing.
 Ahrg fortran ruined me ;).
 ciao
 Fawzi

I got permission from Mark Borgerding to port kissfft and have the binary attribution clause waived, but unfortunately this code is almost as difficult to read due to its extreme preprocessor abuse. (Basically he implements complex numbers as a struct with an imaginary and real component and a bunch of preprocessor macros.) Part of the preprocessor abuse is necessary for supporting fixed point FFTs. Does anyone care about these (I sure don't), or can I just get rid of them?

Fixed-point FFTs could be important if your number-crunching data that could have been collected by, say, a microcontroller-based system with 16-bit ADCs. Integer-only FFTs are useful.
 I'm hoping I can refactor this code into a nice D interface, but it may 
 end up
 being that I have to basically start from scratch.  One thing I will 
 absolutely
 not do is contribute code without a proper D-style (i.e. templated, works 
 with
 std.complex, etc.) interface.  If I can't understand kissfft well enough 
 to do
 this refactoring, then maybe I'll just start from scratch using the 
 Wikipedia FFT
 article and a textbook. 

Aug 02 2010
prev sibling next sibling parent Fawzi Mohamed <fawzi gmx.ch> writes:
On 29-lug-10, at 15:31, dsimcha wrote:

 == Quote from Fawzi Mohamed (fawzi gmx.ch)'s article
 otherwise, as I already said http://www.netlib.org/fftpack/fft.c is
 public domain, has reasonable performance, and supports all sizes.
 It might not be the always beautiful, but it is stable and works  
 well...
 Fawzi

Yeah, I just looked at this code. The problem is that I wouldn't do a straight, mechanical port to D, I'd D-ify the code a little while I was at it (such as making it work on generic random access ranges if it's an out of place transform, making real-only inputs vs. complex inputs DRYer with templates, etc.). Since this code seems to have been straight-up translated from FORTRAN, it scores about a -2 the readability scale, from 1 to 10. In other words, I would probably be able to translate it and make it run, but even the slightest modification would be near impossible.

eheh this must meant that too much fortran inured me, I did not find the code too bad, I classified it as a typical numerical intensive code (even rather clean where one can easily find the boundaries of the arrays, and without too many goto), and templating float, using slices as input seemed easy. But yes probably it is not the best basis to develop your own thing. Ahrg fortran ruined me ;). ciao Fawzi
Jul 30 2010
prev sibling parent Fawzi Mohamed <fawzi gmx.ch> writes:
On 30-lug-10, at 16:14, dsimcha wrote:

 [...]
 Part of the preprocessor abuse is necessary for supporting
 fixed point FFTs.  Does anyone care about these (I sure don't), or  
 can I just get
 rid of them?

I don't care about them
 I'm hoping I can refactor this code into a nice D interface, but it  
 may end up
 being that I have to basically start from scratch.  One thing I will  
 absolutely
 not do is contribute code without a proper D-style (i.e. templated,  
 works with
 std.complex, etc.) interface.  If I can't understand kissfft well  
 enough to do
 this refactoring, then maybe I'll just start from scratch using the  
 Wikipedia FFT
 article and a textbook.

The main thing that you should probably do (apart documenting the factors and order that you assume/generate) is splitting fft in initialization of a given size (creating a struct/ class environment) and executing an fft of that size (which might be a method of the environment (probably better) or a free function taking the env as parameter). Both fftpack (http://www.netlib.org/fftpack/doc) and fftw (http://fftw.org/doc/ ) have that structure. This because you can prepare/precompute things to then do the actual fft in the optimal way. Then you can still provide easy to use wrappers that call the two things in sequence. ciao Fawzi
Jul 30 2010