I am calculating the numerical Fourier transform of an exponential decay exp(-|t|) and compare it to the analytically calculated result, a Lorentzian. I find that the numerically calculated spectrum contains systematically larger amplitudes than the analytical one and that this deviation increases with frequency. I'm using the python library numpy.fft.

My time array contains [-100., -99.99, ..., 0.00, ..., 100.] and the signal is correspondingly [3.72007598e-44, ..., 1.0, ..., 3.72007598e-44]. I've plotted this "original" exponential together with the back Fourier transforms below. signal

I then calculate the FFT (and adjust the phase to get a purely real result). When I compare the numerical result to the expected Lorentzian $2 \cdot \frac{1}{1^2 + (2 \pi f)^2}$, I get very good agreement at low frequencies that however gets worse towards high frequencies. The numerical result (absolute value = real part) is systematically larger than the analytical one (labelled "calc" in this plot). spectrum, all components spectrum, zoom on absolute part

This deviation is of course also visible when plotting the ratio and difference of these two curves ("calc" refers to the analytic result). The ratio of the quantities seems to be independent on the step width and the maximum of the time array. discrepancy between analytic and numerical FT

The inverse FFT of the numerical result (shown in the first plot, labelled "ift") seems to agree even a bit better with the original function than the FFT of the analytically calculated Lorentzian (labelled "icalc"). This depends however on the choice of the time array; it can also happen that the inverse FT of the Lorentzian is slightly longer close to the exponential function.

My objective is to evaluate experimental data in the frequency domain by models in the time domain. I would like to replace the usually used analytical FT (which restricts us to very idealized models) by a numerical FT. The first step is of course to show that the different FTs yield the same result... Of course, the experiment will not measure the spectrum with infinite precision, but I would very much appreciate any help to get at least 1-2 orders of magnitude better in the agreement numerical / analytical.

  • $\begingroup$ What you have discovered is that the continuous case and discrete case are not interchangeable. Intuitively, at low frequencies, the points that describe the curve look a lot like the continuous case. As you up the frequency, the resemblance weakens, as do any formulas you are applying. This is similar to calculating the area under a line. For example the formula for summing whole squares under a line (triangle numbers) is N(N+1)/2 yet the integral is N^2/2. For large N the distinction isn't as great. $\endgroup$ – Cedron Dawg Jul 17 '20 at 14:01

You have aliasing

Your continuous spectrum is not band-limited, so you are violating the sampling theorem. Sampling in time results in periodic repetition in the frequency domain so your first mirror spectra are leaking into your base band spectra.

  • $\begingroup$ Thank you. When adding the contributions of the mirror spectra to the analytical result, it resembles the numerical spectrum more and more. I just displace the spectrum by $n*2*f_{max}$ and add all contributions up, right? It takes a surprisingly large number of mirror spectra but it seems to get there eventually. Is there a way to go the other way and "correct" the numerical FT for this effect? Probably not, right, since I would have to know the shape of the spectrum outside of the range that I can calculate... $\endgroup$ – Sebastian Busch Jul 17 '20 at 18:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.