[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[Help-gsl] Requesting feedback on the FFT interface

From: Dimitris Papavasiliou
Subject: [Help-gsl] Requesting feedback on the FFT interface
Date: Sun, 14 May 2006 19:13:24 +0300
User-agent: Mail/News 1.5 (X11/20060228)

Hello again,

Upon trying to implement a dct transform with the GSL fft routines I noticed a couple of things regarding the fft interface:

1) There are no routines accepting gsl_vector arguments instead of (data, stride, n) triplets. Unless I'm missing something this seems almost a bug to me.

2) There's a weird assymetry in the interface. Complex routines have *_forward, *_transform, *_inverse and *_backward flavors where *_transform accept a direction argument. Real routines on the other hand only *_transform and *_inverse flavors and *_transform simply does an forward transform.

3) Well I suppose there's a good reason for this but the halfcomplex arrays produced by radix2 and mixed-radix routines are different. This tends to make postprocessing of the results of FFTs harder if you need to support both radix2 and mixed-radix routines. For example I thought I'd implement DCTs for both routines and it was enough of a headache figuring out how these halfcomplex sequences where stored, how to multiply the real and imaginary parts with certain sines and cosines (in-place) and where to put the results to end up with the real result of the DCT without having to do it twice.

I'd therefore like to propose the following changes, either as a patch or as an extension, whichever you prefer:

1) Adding gsl_fft_*_*_vector routines.  I have code for that already.

2) Changing both real and complex transforms to a single scheme, either *_forward, *_transform, *_inverse, or *_transform, *_inverse.

3) In order to avoid rewriting the final modulation step of the DCT for radix2 and mixed transforms I have written a small function that permutes the mixed-radix version of the halfcomplex array (the FFTPACK version) to the radix2 version (and back). The radix2 version is actually (at least in my case) better since you can process the array a pair at a time, possibly overwriting each pair without destroying data from other coefficients (well this probably doesn't make sense but I can't put it much better). It is technically O(n^2) but it only performs element swaps and it probably could be optimized. So I think it might be better to use such a function to permute the output of the real mixed-radix functions so that the end-user won't have to a) learn and get accustomed to two halfcomplex formats and b) potentially write code to do the same thing with both formats. We can use my function or a more optimized one but in any case, I doubt the overhead will be a problem.

This one is not a problem with the current interface but:

4) We can add a dct transform pair. It (hopefully) works for arbitrary length sequences.

If the library maintainers and users agree to any of these then I can handle 1) and 4) I suppose but I'd need some help from the original implementor of the FFT routines for 2) and 3) as I don't have a lot of spare time right now. I'd also need some help testing this stuff.

Sorry for the endless blabbering,
Dimitris P.

reply via email to

[Prev in Thread] Current Thread [Next in Thread]