Added 2D Fourier transformation as num.fft2 and num.fft2_inv.

[h4rm]

BTW: Scroll down for updated infos.

Hey Francesco,

I needed some 2D Fourier transformation in my work so I implemented a version. It differs a bit from the 1D calculation, namely I am not using the half complex vectors but I am working with complex numbers from the beginning. The function does 1D transformation for 1 dimension and then 1D transformation for the other dimension. The result is a complex matrix. The result should be comparably fast and probably faster than implementing it on top of fft and fft_inv. I have added the tests and the docs. Also I have cross tested the results with numpy just to make sure to have it all right. Numpy also has the 3D Fourier transformation in stock and I will also implement that for my work, following in a later commit.

Cheers, Ben

[h4rm]

One more thing:

The algorithm work for complex matrices only, so whenever real matrices are given as input I transform the matrix to a complex one with:

matrix.cnew(n1, n2, |i,j| mat:get(i,j))

Maybe there is a better solution.

[h4rm]

BTW, did some performance tests with numpy doing 1e5 2D FFTs on 16x16 matrices. Python is about 2.5 secs, the shell performs at 4.2 secs. The code is even less performant if you go to higher dimensional matrices (256x256 e.g.).

So I guess performance can still be optimized. Maybe the strides can help me to perform fast transformations.

[h4rm]

So, to be honest, after all I have read about FFT, the GSL version is not the best implementation and FFTW almost always supersedes it. It is very tempting to go for FFTW with a wrapper for https://github.com/soumith/fftw3-ffi . The FFTW libraries should be installed on every host machine due to GSL dependencies, so this should be no problem, should it? What do you think?

[soumith]

might be relevant, though not exactly what you're looking for: https://github.com/soumith/torch-signal/blob/master/init.lua

[h4rm]

Thanks soumith for the link, it helped me with the ffi.casts. So having the fftw3-ffi in place, a one dimensional fft is quite simple in gsl-shell:

fftw = require'fftw3.init'

mat = matrix.cnew(10, 1, |i,j|i)
mat2 = matrix.cnew(10,1)

input = ffi.cast("fftw_complex*",mat.data)
output = ffi.cast("fftw_complex*",mat2.data)

plan = fftw.plan_dft_1d(10, input, output, fftw.FORWARD, fftw.MEASURE)

fftw.execute(plan)

print(mat)
print(mat2)

You can also use

plan = fftw.plan_dft_1d(10, input, input, fftw.FORWARD, fftw.MEASURE)

for in-place transformation.

[h4rm]

Okay, the 2D transformation is also extremely easy in comparision to my code:

fftw = require'fftw3.init'

mat = matrix.cnew(4, 4, |i,j|i-1)
mat2 = matrix.cnew(4,4)

input = ffi.cast("fftw_complex*",mat.data)
output = ffi.cast("fftw_complex*",mat2.data)

plan = fftw.plan_dft_2d(4,4, input, output, fftw.FORWARD, fftw.MEASURE)

fftw.execute(plan)

print(mat)
print(mat2)

It produces the same results as numpys.fft2. I am going to prepare an alternative fft implementation for gsl-shell based on the ffi interface to fftw3. I will need it anyway for my projects but we can then discuss if that is a good substitute for the current implementation.

[h4rm]

This can be cancelled.