# Nonuniform fast Fourier transforms

**URL:** <https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858>\
**Category:** Numerics\
**Tags:** nfft\
**Created:** [July 14, 2017, 4:35pm UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858 "2017-07-14T16:35:34Z")\
**Posts on this page:** 11\
**Page:** 1

<div class="post-metadata">

**Author:** ![MikaelSlevinsky](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mikaelslevinsky/32/38545_2.png) [@MikaelSlevinsky](https://discourse.julialang.org/u/MikaelSlevinsky)\
**Post date:** [July 14, 2017, 4:35pm UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/1 "2017-07-14T16:35:34Z")

</div>

[FastTransforms.jl](https://github.com/MikaelSlevinsky/FastTransforms.jl) v0.2.2 includes nonuniform fast Fourier transforms thanks to [Alex Townsend](https://github.com/ajt60gaibb):

- `nufft1` assumes uniform samples and noninteger frequencies;
- `nufft2` assumes nonuniform samples and integer frequencies;
- `nufft3 ( = nufft)` assumes nonuniform samples and noninteger frequencies;
- `inufft1` inverts an `nufft1`; and,
- `inufft2` inverts an `nufft2`.

Basic usage is like so:

```julia
julia> Pkg.add("FastTransforms")

julia> using FastTransforms

julia> n = 10^4;

julia> c = complex(rand(n));

julia> ω = collect(0:n-1) + rand(n);

julia> nufft1(c, ω, eps());

julia> p1 = plan_nufft1(ω, eps());

julia> @time p1*c;
  0.002383 seconds (6 allocations: 156.484 KiB)

```

The new algorithms are described in the following pre-print:

D. Ruiz—Antolín and A. Townsend. [A nonuniform fast Fourier transform based on low rank approximation](https://arxiv.org/abs/1701.04492), arXiv:1701.04492, 2017.

---

<div class="post-metadata">

**Author:** ![tobias.knopp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tobias.knopp/32/7551_2.png) [@tobias.knopp](https://discourse.julialang.org/u/tobias.knopp)\
**Post date:** [July 15, 2017, 7:46am UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/2 "2017-07-15T07:46:54Z")

</div>

@MikaelSlevinsky: Are you aware of [https://github.com/tknopp/NFFT.jl](https://github.com/tknopp/NFFT.jl)? It seems that nufft just allows for 1D NFFTs whereas NFFT.jl allows for d-dimensional transforms.

I am pretty interested in a performance comparison of both approaches.

---

<div class="post-metadata">

**Author:** ![tobias.knopp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tobias.knopp/32/7551_2.png) [@tobias.knopp](https://discourse.julialang.org/u/tobias.knopp)\
**Post date:** [July 15, 2017, 7:53am UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/3 "2017-07-15T07:53:27Z")

</div>

And from an interface point of view: Isn’t `nufft2` the adjoint of `nufft2` (modulo conjugation). The NFFT is quite often used as an operator when solving linear systems of equations and there you usually want that pair. In my (internal) software on MRI reconstruction we overload `A_mul_B!` and `Ac_mul_B!` for this purpose.

---

<div class="post-metadata">

**Author:** ![MikaelSlevinsky](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mikaelslevinsky/32/38545_2.png) [@MikaelSlevinsky](https://discourse.julialang.org/u/MikaelSlevinsky)\
**Post date:** [July 15, 2017, 5:19pm UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/4 "2017-07-15T17:19:24Z")

</div>

Hi Tobias, there is some support for 2D nufft’s of types I-I and II-II, but the multi-dimensional support is not fully fleshed out. You are right, `nufft1` and `nufft2` are adjoints, modulo conjugation.

I’m interested in a performance comparison, too, and the first step is actually getting a release of this approach.

There are four levels to the 1D interface:

- the functions `nufft1/2/3(args...)` plan & execute;
- the functions `plan_nufft1/2/3(args...)` construct plans;
- if `p` is an `NUFFTPlan`, then `p*c` applies the plan; and,
- `A_mul_B!(out, p, in)` applies the plan to the `in` array in-place on the `out` array.

For a fair comparison, methods should be compared at the same tolerance `ϵ`.

Note that a hierarchical approach to the Chebyshev–Legendre transform accelerated the execution phase (of any previous methods in FastTransforms.jl) by a factor of ~20. With some effort, Dutt & Rokhlin’s hierarchical nonuniform FFT method could be implemented based on the [HierarchicalMatrices.jl](https://github.com/MikaelSlevinsky/HierarchicalMatrices.jl) framework.

---

<div class="post-metadata">

**Author:** ![tobias.knopp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tobias.knopp/32/7551_2.png) [@tobias.knopp](https://discourse.julialang.org/u/tobias.knopp)\
**Post date:** [July 15, 2017, 8:19pm UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/5 "2017-07-15T20:19:57Z")

</div>

I just saw that the paper actually makes a comparison with the NFFT.jl, which is great.

I wonder if it would be good to share an interface. NFFT.jl could have different implementations that could be switchable.

---

<div class="post-metadata">

**Author:** ![Alex\_Townsend](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/alex_townsend/32/42721_2.png) [@Alex\_Townsend](https://discourse.julialang.org/u/Alex_Townsend)\
**Post date:** [July 15, 2017, 9:37pm UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/6 "2017-07-15T21:37:01Z")

</div>

Yes, in Figure 2.2 we did a quick comparison to NFFT.jl for the 1D NUFFT of type-I:

> **[1701.04492.pdf](https://arxiv.org/pdf/1701.04492.pdf)**
>
> 324.16 KB

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [July 21, 2017, 11:09am UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/7 "2017-07-21T11:09:25Z")

</div>

> [@tobias.knopp](#):
>
> I wonder if it would be good to share an interface.

Maybe a meta-package a la DifferentialEquations.jl is in order? Transforms.jl is available it looks like…there could even be an org JuliaTransforms.

---

<div class="post-metadata">

**Author:** ![tobias.knopp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tobias.knopp/32/7551_2.png) [@tobias.knopp](https://discourse.julialang.org/u/tobias.knopp)\
**Post date:** [July 21, 2017, 11:30am UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/8 "2017-07-21T11:30:15Z")

</div>

No, in this case I do not really see a need for that. It should be possible to have a single NFFT package that allows to switch different algorithms in its backend.  
FastTransforms.jl is more high level and can then integrate the low level packages (where NFFT is just one example)

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [July 21, 2017, 4:32pm UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/9 "2017-07-21T16:32:05Z")

</div>

I think that’s what I meant: a package that wraps the different algorithms for the different transforms, with a consistent naming scheme and consistent way to specify algorithms.

---

<div class="post-metadata">

**Author:** ![nantonel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nantonel/32/2889_2.png) [@nantonel](https://discourse.julialang.org/u/nantonel)\
**Post date:** [May 8, 2018, 8:46am UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/10 "2018-05-08T08:46:33Z")

</div>

Talking about a common interface maybe this package could be interesting: [AbstractOperators.jl](https://github.com/kul-forbes/AbstractOperators.jl).

This offers a common interface that allows to easily combine linear and nonlinear transformations and use them in inverse problems (through the package [StructuredOptimization](https://github.com/kul-forbes/StructuredOptimization.jl)).

In case you find this interesting feel free to contribute!

---

<div class="post-metadata">

**Author:** ![Palli](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/palli/32/3380_2.png) [@Palli](https://discourse.julialang.org/u/Palli)\
**Post date:** [July 13, 2020, 11:10am UTC](https://discourse.julialang.org/t/nonuniform-fast-fourier-transforms/4858/11 "2020-07-13T11:10:04Z")

</div>

FYI: There’s also [GitHub - ludvigak/FINUFFT.jl: Julia interface to the nonuniform FFT library FINUFFT](https://github.com/ludvigak/FINUFFT.jl) that wasn’t available during the above discussion (a wrapper for [flatironinstitute](https://github.com/flatironinstitute)’s C++ code, or I guess its C interface), nor the paper on it:

> **[A parallel non-uniform fast Fourier transform library based on an...](https://arxiv.org/abs/1808.06736)**
>
> The nonuniform fast Fourier transform (NUFFT) generalizes the FFT to off-grid data. Its many applications include image reconstruction, data analysis, and the numerical solution of differential equations. We present FINUFFT, an efficient parallel...

> We present FINUFFT, an efficient parallel library for type 1 (nonuiform to uniform), type 2 (uniform to nonuniform), or type 3 (nonuniform to nonuniform) transforms, in dimensions 1, 2, or 3. It uses minimal RAM, requires no precomputation or plan steps, and has a simple interface to several languages. […] For types 1 and 2, rigorous error bounds asymptotic in the kernel width approach the fastest known exponential rate, namely that of the Kaiser–Bessel kernel. We benchmark against several popular CPU-based libraries, showing favorable speed and memory footprint, especially in three dimensions when high accuracy and/or clustered point distributions are desired.

Also interesting: [GitHub - gwater/HexFFT.jl: Fast Fourier transform on hexagonal grids using Birdsong and Rummelt's algorithm](https://github.com/gwater/HexFFT.jl) “implements [Birdsong and Rummelt’s 2016 algorithm](http://ieeexplore.ieee.org/document/7532670/) for fast Fourier transforms on hexagonal lattices.”

> The current implementation is ca. 10 times slower than `Base.fft()` for comparable rectangular grids. However significant optimization should be possible by preallocating and reusing temporary arrays.