Fft() function does not return the correct analytic result

Specific Domains Signal and Image Processing
,
2 
import FFTW
using PyPlot

Nx = 1024
x = range(-10, 10, length=Nx)
dx = x[2] - x[1]
k = 2 * pi * FFTW.fftfreq(Nx, 1/dx)
k = FFTW.ifftshift(k)

a = 1.0
f = [abs(i) <= a / 2 ? 1 : 0 for i in x]

sa = @. a * sinc(k * a / 2 / pi)   # sinc(x) = sin(pi * x) / (pi * x)

s = FFTW.fftshift(f)   # compensate osciallations
s = FFTW.fft(s) * dx
s = FFTW.ifftshift(s)

plot(k, sa)
plot(k, real.(s))
show()
  1. Note that sinc function in Julia is defined as normalized sinc function: sinc(x) = sin(pi * x) / (pi * x).
  2. Spectrum has extra oscillations due to the way how fft algorithm stores arrays. You should apply extra fftshift before calling fft in order to avoid these oscillations (shift in space = frequency modulation in spectrum).
  3. By default FFTW does not normalize the integrals by 2pi.
2 Likes