I am trying to calculate the first N Fourier coefficients for sampled data (harmonic + gaussian like). The Riemann sum was my first approach, however I think the results are not enough accurate. Therefore, I am trying to choose a more suitable method. The QuadGK package seems to be extremely accurate, but it only works (I think…) for analytical functions. Among the others integration packages, which would be your recommendation?

More infomation is needed.

Sampled data = experimental data? From a simulation? Noisy? Smooth function? If it is from a computer program, then can you evaluate it on demand for any desired x?

When you say “N Fourier coefficients” do you mean that the data is periodic and you want a Fourier series? If it’s not periodic, then what do you mean?

If you have equally-spaced samples of a periodic function, then a Riemann sum actually converges exponentially fast (with the number of samples) for smooth functions.

(However, adaptive methods like QuadGK can still be advantageous for functions with localized features like sharp peaks.)

And, of course, a Riemann sum of equally-spaced samples is exactly equivalent (with suitable scaling) to a discrete Fourier transform (DFT), so you can compute it very efficiently with FFT algorithms.

QuadGK and similar methods assume that you can evaluate f(x) at any requested x in the integration domain. It doesn’t require you to have an analytical formula for f(x), but you need a way to compute it on demand (e.g. via a computer program).

In addition to the above answer, if your integrant is not periodic and you have a discrete set of fixed (possibly irregular) sampling points, you can try Filon’s method or related methods. I implemented it once in Julia, but I’m not sure it is available in a package.

Edit: a quick search gives

And

Sampled data comes from a function *H(V, u)* that depends (sum) of transversal harmonic complex potential *V(x)* and nonlinear terms that depends of a localized complex function *u(x)*, where *u(x)* is a mode for a nonlinear wave equation in periodical media.

Here is the Riemann approach:

```
N # Number of modes
M # Grid size
L # Transversal length
k0 = 2*pi/L
dx = L/M
x = -L/2:dx=L/2-dx
c = zeros(N)
for n in 1:N
c[n]=dx*sum(H.*exp.(-im*k0*n*x))/L
end
```