Evaluating the continuous Fourier transform of a constant, and matching it up with the FFT result

Question

I am following my optics textbook (Optics, by Eugene Hecht), throughout which are given various exact analytical results for the diffraction patterns that result from light passing through differently shaped apertures. These diffraction patterns are often calculated by taking continuous Fourier transforms. Ultimately, I am trying to get to the stage where I can use a FFT to numerically calculate these patterns when an analytical solution is not possible. Therefore, as a first sanity-check step I would like to check that my FFT results match up with the analytical Fourier transforms in some simple cases including the correct scaling factors.

Unfortunately, I have failed at the first hurdle, in the simplest possible case I can think of!

The book gives the following definition for the continuous Fourier transform pair:

\begin{align} f(x) &= \frac{1}{2\pi} \int_{-\infty}^{+\infty} F(k)\; e^{-ikx}dk\hspace{1cm} &\text{Inverse Transform} \\ F(k) &= \int_{-\infty}^{+\infty} f(x)\; e^{ikx}dx\hspace{1cm} &\text{Forward Transform} \end{align}

The book then goes on to show some examples of functions along with their transforms - in particular, the two simple cases shown below. In Figure (a) it is shown that the spectrum of a constant $f(x) = A$ has a peak located at $k=0$. The peak is shown as an arrow, with a label of $2\pi A$:

$\hskip1in$

I tried to arrive at this by proceeding as follows:

\begin{align} F(k) &= \int_{-\infty}^{+\infty} f(x)\; e^{ikx}dx \\ &= \int_{-\infty}^{+\infty} A\; e^{ikx}dx \\ &= A\int_{-\infty}^{+\infty}\cos(kx)\;dx\; +\; iA\int_{-\infty}^{+\infty}\sin(kx)\;dx \end{align}

Now, because the cosine and sine waves extend from $\pm\infty$, I know that their areas will average out and the integrals will evaluate to zero EXCEPT for the case when $k=0$. In this special case, the sine term will still evaluate to zero (because $\sin(0)=0$), however the cosine term will evaluate to $\infty$, because the area under $\cos(0)=1$ will be infinite. This behaviour is captured by the Dirac delta function: $$\delta(k) = \left\{ \begin{array}{1 1} +\infty & \quad k = 0\\ 0 & \quad k\neq 0 \end{array} \right.$$ Therefore, I would say that the transform should evaluate to: $$F(k) = A\;\delta(k)$$ and when $k=0$ this should take a value of $+\infty$. I don't understand (a) why the value is finite, and (b) why there is a factor of $2\pi$ in the book's graph.

I also then tried to calculate numerically the value of the spectrum coefficient using Matlab's FFT function, which implements the DFT, using a value of $A=3$. According to the textbook, I guess I should see a value of $2\pi\cdot 3=18.8$.

num_points = 10;
dx = 0.5;                     % Sampling interval
x = 0:dx:(dx*num_points-dx);  % Generate position vector
A = 3;                        % Value of constant function
f = A*ones(num_points,1);     % Input function f(x) = A
F = fft(f)*dx;                % Scale FFT

disp(['DC Component = ',num2str(F(1))])

I understand there are some disagreements between people for right way to scale the FFT output (see here for example). Since the DFT is defined as $$Y_k = \sum^{N-1}_{n=0}y_n\; e^{-2\pi i k n/N}$$

it seems to make sense to me that this is a summing approximation to the continuous integral version of the FT above, and so should be scaled by $dx$ in order to give results that can be compared directly - so that is what I have done in the Matlab snippet. Unfortunately the script shows a value for the DC peak of 15 (and depends on the number of samples, which must be wrong).

Can anyone help me resolve these conflicts? As I say, I am looking to be able to calculate numerical FFTs and have them match the continuous analytical versions. I would like to move on to try a more complex function that can still be done analytically (such as a pure cosine, or a Gaussian), but if I can't do the simple constant value I don't think I should move on just yet!

---

Similar SE questions for reference:

Fourier Transform and Delta Function

About Fourier transform of periodic signal

Analytical Fourier transform vs FFT of functions in Matlab

how fft points should be rescaled to get the same results as the analytical solution?

Matlab FFT for gaussian function

Answer 1

Welcome! Regarding your questions:

a) Its value is not finite. $A$ (or $2\pi A$, which is the correct answer) is not its value, it's the delta function's "area". $\delta(k)$ is not a proper function, it's a distribution or a generalized function. When we sketch such functions, we denote their area as a scaling factor. Indeed the actual value of the delta function is $+\infty$ at $k=0$.

b) You recognize that the forward Transform is somewhat problematic. Try the inverse one to see that the FT pair holds: $$f(x) = \frac{1}{2\pi}\int_{-\infty}^{+\infty}F(k)e^{-jkx}dk = \frac{1}{2\pi}2\pi A\int_{-\infty}^{+\infty}\delta(k)e^{-jkx}dk= Ae^{-jkx}\Big|_{k=0} = A$$

So when you use MATLAB to approximate this pair, you can never get the $A=3$ number. I would use another piece of code to approximate the FT of the constant, which gives better control on scaling etc.

 dx = 0.01;              % x-step
 x = -10:dx:10;          % x-axis
 A = 3;                  % constant value
 f = A*ones(size(x));    % just an approximation of the constant
 dk = 0.01;              % k-step
 k = -10:dk:10;          % k-axis
 F = dx*f*exp(j*k'*x).'; % Riemann approximation of CTFT
 plot(k, F);             % plot

As you can see from the plot, you can only approximate constant $A$ by some time-limited rectangular pulse. In the $k-$domain you get something that approaches the behavior of the delta function more and more as you increase the samples of your x-axis and the duration of your rectangular pulse.

You'll never get a perfect dirac delta because you would need infinitely many samples of your $k$ and $x$ axes. What you really get in the transform is the so-called sinc function that approaches the delta function under the notion of limit.