Numerical computation of continuous Fourier transforms

Are there any packages available for computation of Fourier transforms of the form F(x_k) = \int_{-\infty}^\infty f(t)\mathrm{e}^{-\mathrm{i}tx_k}\,\mathrm{d}t? e.g. the method developed by DBailey and Swarztrauber (The Society for Industrial and Applied Mathematics). The function f(t) is rapidly decaying.

You could always use a generic quadrature routine, though that will have suboptimal efficiency if you have many frequencies x_k at which you want to compute F(x_k).

Technically, that paper describes a Bluestein’s (1968) algorithm for the “Zoom FFT”, generalized to the chirp-z transform by Rabiner and Schafer (1969), which is not the same as the continuous Fourier transform (note: Bailey and Swarztrauber’s terminology was not widely adopted, because a fractional Fourier transform already meant something else). (It computes a sum, not an integral.)

You can use this sum to approximate a continuous Fourier transform of a rapidly decaying function if you sample finely enough and normalize correctly, of course.

I don’t know of a Julia package for this offhand, but it’s pretty easily to implement given an fft function (e.g. from the FFTW.jl package), since it’s just a convolution that can be written in a few lines of code.

More generally, if your frequencies are not equally spaced, you could use a nonuniform FFT, e.g. FINUFFT.jl or FastTransforms.jl.

4 Likes

For non-equidistant NFFTs please have a look at GitHub - JuliaMath/NFFT.jl: Julia implementation of the Non-equidistant Fast Fourier Transform (NFFT), which is nowadays a package family around NFFTs. One can use AbstractNFFTs as the base package and then we have different implementations (Julia based, Julia/Cuda based, NFFT3, FINUFFT). Performance-wise the Julia implementation is one of the fastest: Performance · NFFT

3 Likes

You may find

https://github.com/JuliaApproximation/OscillatoryIntegrals.jl

useful. It gives a method whose computational cost is independent of the frequency.

2 Likes