# Help needed with FFT and NFFT
**URL:** https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449
**Category:** Specific Domains
**Tags:** question, nfft, fft
**Created:** [January 24, 2023, 10:18am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449 "2023-01-24T10:18:18Z")
**Posts on this page:** 15
**Page:** 1
### Author:  [@Karajan](https://discourse.julialang.org/u/Karajan)
#### Post date: [January 24, 2023, 10:18am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/1 "2023-01-24T10:18:18Z")
Hi, I think I could use some help understanding nonuniform discrete Fourier transforms (NFFT). Either I’m not using the tools correctly or I’m misunderstanding how they work in general.
I read that when my signal has nonuniform time steps I should use a NFFT (type II, I think) instead of a FFT. Since my results looked all kind of weird, I tried with a small example, shown below. This is similar to what I want to use it on: multiple periods of sine waves with uneven time steps. The goal is to detect additional harmonics (of course there are none in this example).
```julia
x = range(0, 360; length=200)
x´ = x[1:end-1] # skip 360° so it's periodic
c = cosd.(x´)
using FFTW
rfft(c) .|> abs |> plot
using FastTransforms
nufft2(complex.(c), x´, eps()) .|> abs |> plot
using NFFT
k = x´ ./ x[end] ./ 2 # nfft wants k to be between -0.5 and 0.5
nfft(k, complex.(c)) .|> abs |> plot
```
I don’t understand what’s going on with `nfft`, while the results so `nufft2` make somewhat sense to me. However, if I use `length=201` `nufft2` is just zeros.
When I add the line
```julia
x´ = x´ .+ randn(length(x´))
```
before running `cosd` to simulate the uneven time steps, `rfft` still looks great, `nfft` looks better-ish than before, while `nufft2` is all over the place.
I could even use `10randn` and the `rfft` results still look good – shouldn’t I just always use the normal FFT? Breaking the FFT isn’t that hard (quadratically spaced time steps will do it) so it would be good to know how to work with NFFTs nonetheless.
If helpful I can also write a bit more about my actual application.
---
### Author:  [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)
#### Post date: [January 24, 2023, 12:58pm UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/2 "2023-01-24T12:58:26Z")
> [@Karajan](#):
>
> ```julia
> using NFFT
> k = x´ ./ x[end] ./ 2
>
> ```
This maps `k` over [0, 0.5). Mapping it over [-0.5, 0.5) produces good results.
---
### Author:  [@Karajan](https://discourse.julialang.org/u/Karajan)
#### Post date: [January 25, 2023, 6:08pm UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/3 "2023-01-25T18:08:07Z")
Ah yes, that must have been the remnants of some experiment I tried.
So the results I’m seeing here look okay (I need to get better at reading FFT from complex values…). However, if I now go to `length=201` the results look off again, which doesn’t happen with `rfft`. Is this expected?
---
### Author:  [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)
#### Post date: [January 25, 2023, 7:02pm UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/4 "2023-01-25T19:02:15Z")
> [@Karajan](#):
>
> if I now go to `length=201` the results look off again
In the [NFFT.jl](https://github.com/JuliaMath/NFFT.jl) package’s [reference paper linked](https://arxiv.org/abs/2208.00049), it is specified that the size N, of the NDFT in the equidistant sampling domain, is restricted to even numbers.
---
### Author:  [@tknopp](https://discourse.julialang.org/u/tknopp)
#### Post date: [January 25, 2023, 10:17pm UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/5 "2023-01-25T22:17:03Z")
That issue is resolved. We now also support odd N. I need to update the arXiv paper…
---
### Author:  [@Karajan](https://discourse.julialang.org/u/Karajan)
#### Post date: [January 26, 2023, 7:45am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/6 "2023-01-26T07:45:35Z")
I did a quick comparison, with less points so it clearer. (`length` refers here to `c` and `x´`, not `x`).
And one where the `x´` are shifted by`randn(length(x´))`.
(Here I calculated `k = (x´ .- minimum(x´)) ./ (maximum(x´) - minimum(x´)) .- 0.5`)
Are these the correct results?
To my eyes in all cases `fft` looks better. What would be the propper use case for NFFT then?
---
### Author:  [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)
#### Post date: [January 26, 2023, 8:05am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/7 "2023-01-26T08:05:11Z")
> [@Karajan](#):
>
> `k = (x´ .- minimum(x´)) ./ (maximum(x´) - minimum(x´)) .- 0.5`)
I think here you may be mapping onto [-0.5, 0.5].
As indicated, if you map over [-0.5, 0.5) you get better NFFT results, which overlay the FFT perfectly:
---
### Author:  [@Karajan](https://discourse.julialang.org/u/Karajan)
#### Post date: [January 26, 2023, 8:26am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/8 "2023-01-26T08:26:59Z")
Ah, thanks for spelling it out for me. How do you make that happen though, could you share your code?
Having `k = -0.5:0.03447285262855501:0.4997127262280953` for example didn’t make it any better.
---
### Author:  [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)
#### Post date: [January 26, 2023, 8:37am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/9 "2023-01-26T08:37:17Z")
The code for the plot above is provided here below:
```julia
using Plots; gr()
n = 30
x = 360/n * (iseven(n) ? (-n÷2:(n÷2-1)) : (-(n-1)÷2:(n-1)÷2))
c = cosd.(x)
using FFTW
using Plots; gr()
fft(c) |> fftshift .|> abs |> plot
using NFFT
k = x / 360 # k over [-0.5, 0.5) if even, (-0.5, 0.5) if odd
nfft(k, complex.(c)) .|> abs |> x -> plot!(x; ls=:dash)
```
_ **EDIT:** _
_The symmetry of support x above has been improved for the n-odd case_
---
### Author:  [@Karajan](https://discourse.julialang.org/u/Karajan)
#### Post date: [January 26, 2023, 1:05pm UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/10 "2023-01-26T13:05:17Z")
That looks much better, thank you.
Is this the general rule to form `k`? If so, would it make sense to have `nfft` calculate this automatically? If not, how would I do this with my uneven timesteps? I don’t really understand why this is correct and the other `k` aren’t.
The documentation says
and in the examples uses
```julia
k = range(-0.4, stop=0.4, length=J)
k = rand(2, J) .- 0.5
k = Float32.(LinRange(-0.5,0.5,64))
```
in different places. From there I didn’t know any better than “it doesn’t really matter” which clearly isn’t true. If all of your answers are basic stuff I’m not aware of, where could I read up on this?
---
### Author:  [@tknopp](https://discourse.julialang.org/u/tknopp)
#### Post date: [January 28, 2023, 9:02am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/11 "2023-01-28T09:02:06Z")
> [@Karajan](#):
>
> Is this the general rule to form `k`?
No, there is no general rule for forming `k` and it makes no sense to let the NFFT automatically chose `k` in an equidistant way.
The point here really is that the NFFT and the FFT are two different things. Only in case that the nodes are taken equidistantly, they are the same. Otherwise not.
In the examples you show we just use some exemplary set of nodes, which you need to adapt to the actual nodes you have in your application.
> [@Karajan](#):
>
> From there I didn’t know any better than “it doesn’t really matter” which clearly isn’t true.
yes, shifting nodes in frequency space has a huge impact.
> [@Karajan](#):
>
> If all of your answers are basic stuff I’m not aware of, where could I read up on this?
Well, one needs some deeper experience with Fourier mathematics in order to use the NFFT correctly. Since you aim at calculating the Fourier coefficients from a time signal you likely need to investigate the inverse NFFT, which is different from the adjoint NFFT. Probably this documentation [Tools · NFFT](https://juliamath.github.io/NFFT.jl/dev/tools/)
or this paper
[https://downloads.hindawi.com/journals/ijbi/2007/024727.pdf](https://downloads.hindawi.com/journals/ijbi/2007/024727.pdf)
are helpful. In MRI the NFFT is regularly used, usually under the name “gridding”.
---
### Author:  [@Karajan](https://discourse.julialang.org/u/Karajan)
#### Post date: [January 30, 2023, 5:52pm UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/12 "2023-01-30T17:52:46Z")
> [@tknopp](#):
>
> No, there is no general rule for forming `k` and it makes no sense to let the NFFT automatically chose `k` in an equidistant way.
Just so there is no misunderstanding, I meant the process of scaling of the time steps I have so they are `[-0.5, 0.5)`. It’s clear that I can’t “invent” an equidistant grid if I want to take advantage of NFFT.
> [@tknopp](#):
>
> Well, one needs some deeper experience with Fourier mathematics in order to use the NFFT correctly.
Well, that is unfortunate. I took a look at the NFFT paper ([https://arxiv.org/pdf/2208.00049.pdf](https://arxiv.org/pdf/2208.00049.pdf)) and got immediately lost, but I’ll take a look at the other resources, too. I understand the basic priciple behind the normal fourier transform and can use a FFT, but that’s about where it ends.
I was hoping for a solution where could give a function the time steps in addition to the signal and get a more correct result than using a normal FFT on the signal (or at least equally straight-forward). If that isn’t possible and it takes getting deep in the mathematical weeds I probably can’t justify spending the time, sadly.
---
### Author:  [@tknopp](https://discourse.julialang.org/u/tknopp)
#### Post date: [February 2, 2023, 10:05am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/13 "2023-02-02T10:05:43Z")
> [@Karajan](#):
>
> Just so there is no misunderstanding, I meant the process of scaling of the time steps I have so they are `[-0.5, 0.5)`. It’s clear that I can’t “invent” an equidistant grid if I want to take advantage of NFFT.
Ok, this scaling is no magic. Within the FFT ([Fast Fourier transform - Wikipedia](https://en.wikipedia.org/wiki/Fast_Fourier_transform)) one has the term `kn/N` in the exponential. The `k/N` leads to sampling points being located between [0, 1). If we wrap this at 0.5 (open interval) we get the interval [-1/2, 1/2).
---
### Author:  [@marius311](https://discourse.julialang.org/u/marius311)
#### Post date: [February 2, 2023, 10:46am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/14 "2023-02-02T10:46:07Z")
Not totally your question but this really helped me understand the NFFT, so figured I may as well share somewhere. With an FFT, you’re rewriting your function sampled at N grid points into a series sum of N sines and cosines. If you were to inverse FFT, you’re evaluating those sines and cosines on the original grid (getting exactly your function back). But with the NFFT you can evaluate that series at any arbitrary points. This gets you something like interpolation (since its basically assuming theres no power beyond the nyquist of the original grid):
```julia
using NFFT, Distributions, Random, PyPlot
N_coarse = 32
N_fine = 256
grid_coarse = range(-0.5, 0.5, length=N_coarse+1)[1:end-1]
grid_fine = sort(rand(Uniform(-0.5, 0.5), N_fine))
f = randn(N_coarse)
f_resampled = real.(nfft(grid_fine, nfft_adjoint(grid_coarse, N_coarse, f)) / N_coarse)
plot(grid_coarse, f, "o-", label="original")
plot(grid_fine, f_resampled, ".-", label="resampled")
legend()
ylim(-2, 2)
```
---
### Author:  [@Karajan](https://discourse.julialang.org/u/Karajan)
#### Post date: [February 12, 2023, 11:36am UTC](https://discourse.julialang.org/t/help-needed-with-fft-and-nfft/93449/15 "2023-02-12T11:36:10Z")
Sorry for the delay.
This helped very much indeed. I wasn’t really aware that an NFFT is more like an inverse FFT rather than a FFT (that thought didn’t even cross my mind, and it didn’t seem terribly obvious from the stuff I read).
That makes this quote make much more sense to me as well
> [@tknopp](#):
>
> Since you aim at calculating the Fourier coefficients from a time signal you likely need to investigate the inverse NFFT, which is different from the adjoint NFFT. Probably this documentation [Tools · NFFT](https://juliamath.github.io/NFFT.jl/dev/tools/)
– I started digging into it, hopefully I can make some progress.