Fourier transform

Fourier transform is linear integral transform with the exponential kernel.

Let some complex–valued function \(A\) be defined for real values of the argument, id est,
 * \(A : \mathbb R \mapsto \mathbb C\)

Then, function \(B : \mathbb R \mapsto \mathbb C\) can be defined with
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!(1) ~ ~ ~ \displaystyle B(x)= \frac{1}{\sqrt{2\pi}} \int_{- \infty}^{\infty} \exp( - \mathrm{i} x y) A(y) ~ \mathrm{d} y \)

If the integral converges, then, function \(B\) is called Fourier transform of function \(A\).

The Fourier transform is used in various sciences and, especially, in physics. In Quantum mechanics, function \(A\) may have sense of the wave function in the coordinate representation, then \(B\) corresponds to the momentum representation. Many equations can be solved using the Fourier Transform. Especially efficient is the Fourier transform in the consideration of the paraxial approximation of the wave equations.

The Fourier-transform (FT) is often confused to the Discrete Fourier Transform (DFT). The DFT can be considered as a FT of a function that can be represented as combination of equidistant delta–functions, while the FT can be represented as a fundamental sequence of the FTs with increasing of length of the array it deals with.

Fourier operator
On the equation (1), the Fourier operator \(ふ\) can be defined in the following way:
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!(2) ~ ~ ~ \displaystyle (\hat ふ A)(x)= \frac{1}{\sqrt{2\pi}} \int_{- \infty}^{\infty} \exp( - \mathrm{i} x y) A(y) ~ \mathrm{d} y \)

Then the equation (1) can be qritten as
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!(3) ~ ~ ~ \displaystyle B=　\hat ふ A\)

without to indicate the argument.

The inverse of Fourier operator is also its Hermitian–conjugated, id est,
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!(4) ~ ~ ~ \displaystyle \hat ふ^{-1}(A)(x)=\hat ふ^{*}(A)(x)=\frac{1}{\sqrt{2\pi}} \int_{- \infty}^{\infty} \exp( \mathrm{i} x y) A(y) ~ \mathrm{d} y \)

For the integrable continuous functions \(A\) and \(B\), it is assumed that
 * if \(~A=\hat ふ(B)~\) then \(~B=\hat ふ^* (A)~\)

General properties of the Fourier operator
The Fourier operator is linear: for complex constants \(\alpha\) and \(\beta\),
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!(6) ~ ~ ~ \displaystyle \hat ふ (\alpha A+\beta B)= \alpha \hat ふ(A) + \beta \hat ふ (B)\)

On the set of continuous functions \(A\), the operatot product
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!(7) ~ ~ ~ \displaystyle \hat ふ ~\hat ふ ^*=I~\)

is interpreted as identity operator. On the language of generalized functions, this can be expressed in the following way:
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!(8) ~ ~ ~ \displaystyle \int_{-\infty}^{\infty} \exp(\mathrm {i} p x) ~\mathrm {d} x\) \(=\) \( 2\pi ~\delta(p)\)

where \(\delta\) is the Dirak's delta function. This can be seen substituting the integral representations for \(\hat ふ\) and \(\hat ふ ^*\) into the left hand side of equation (7).

Norm of the Fourier operator
The norm \(||A||\) of a function \(A\) can be defined as
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (9) ~ ~ ~ ||A||^2 = \int_{-\infty}^{\infty} |A(x)|^2 ~\mathrm{d} x\)

The norm of any function is assumed to be non–negative; \(||A||\ge 0\).

The Fourier operator is unitary, so, it preserves the norm of a function:
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (10) ~ ~ ~\) if \(~ A=\mathcal{F}(B)~\) then \(||A||=||B||\)

The iterations of the Fourier operator ho not provide a wide variety, because
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (10.1) ~ ~ ~ \mathcal{F}^4= I\)

In particular, for any absolutely integrable continuous function \(A\), the linear combination
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (10.2) ~ ~ ~ S(A)= A+ \mathcal{F} A+ \mathcal{F}^2 A+ \mathcal{F}^3 A\)

is self-Fourier function; \(~S(A)~\) is eigenfunction]] of the Fourier operator \(\mathcal{F}\) with eigenvalue unity.

Self–Fourier functions
Function \(A\) is called self-Fourier function, if
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (11) ~ ~ ~ \mathcal{F} A = A\)

All the self-Fourier functions are eigenfunction of the Fourier operator with eigenvalue unity. Any self–Fourier function can be represented as series,
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (12) ~ ~ ~ A(x) ~=~ \sum_{n=0}^{\infty} ~c_n~ \mathrm{HermiteH}_{4n}\!(x)~ \exp(-x^2/2)\)

where \(\mathrm{HermiteH}\) are the Hermit polynomials , and \(c\) are arbitrary complex coefficients; the sequence \(c\) should decay quickly enough to provide the convergence of the series .

In particular, the first three self-Fourier functions \(A_0, A_4,A_8\) can be defined with
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (13) ~ ~ ~ A_0(x)=\exp(-x^2/2)\)
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (14) ~ ~ ~ A_4(x)=\exp(-x^2/2) (x^4-3x^2)\)
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (15) ~ ~ ~ A_8(x)=\exp(-x^2/2) (x^8-14x^6+35x^4)\)

The self-Fourier functions are good for testing of the numerical implementations of the Fourier operator \(\mathcal{F}\).

Numerical implementation of the Fourier operator


In order to approximate the Fourier operator, the Discrete Fourier operator is used. In this section, the numerical implementation is discussed. The C++ implementation fafo.cin of the Fourier operator \(\mathcal {F}\) is used to plot figures. Up to year 2011, it seems to be simplest and the most portable routine. This implementation is described below. The codes used to plot figures are supplied in the descriptions of the figures. (In TORI, the figures are clickable.)

There is no need to describe functions \(f(x)\!=\!\sin(2x)\), \(f(x)\!=\!\sin(kx)\) , \(f(x)\!=\!\sin(2 \pi x/T)\) and so on separately, it is sufficient to have function \(\sin\). In the similar way, it is sufficient to describe the discrete implementation of the Fourier operator at some grid. Then, other implementations can be built-up with minimal modifications of the algorithm, shifting and scaling the argument of the function and applying the corresponding modifications to the Fourier transform.

For the nimber \(N=2^m\) grid points, the step \(d\) of the grid is chosen as follows:
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (16) ~ ~ ~ d = \sqrt{2 \pi /N}\)

The grid points \(x_j\), \(j=0..N\) are chosen in the following way:
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (17) ~ ~ ~ x_j = (j-N/2) d ~\) for natural number  \(~j~\); \(~0\le j < N\)

At such a discretization, the node with number \(N/2\) always corresponds to the argument zero. The 0th node is assumed at the left hand side of the sampled range. However, it can be placed also at the right hand side as well, and interpreted as \(x_{N}\). The algorithm makes no difference between these points.

The same grid is used both for the function and for its Fourier–transform. Functions \(A\) and \(B=F(A)\) are represented with the arrays;
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (18) ~ ~ ~ A(x_j)=A_j~\) ; \(~B(x_j)=B_j~\)

The approcimation of \(B\) corresponds to the replacement of the integral to the discrete sum:
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (19) ~ ~ ~

B_j=\frac{d}{\sqrt{2 \pi}} \sum_{n=0}^{N-1} A_n \exp( -\mathrm{i} ~ x_j ~  x_n) \)

It is convenient to chose \(N=2^m\) for some natural number \(m\); then the fast implementation of such a summation is especially simple.

For the testing of the algorithm, the transform of discrete analogy of the Gaussian exponential \(A(x)\!=\!\exp(-x^2/2)\)  is suggested at the figure above. Such a Gaussian is self-Fourier function, and it coincides with its Fourier-transform. However, at small number of the grid points, \(N\!=\!4\), the deviation of the blue segmented line (that represents the array \(A\)) deviates from the red line, that represents the array \(B\) as the numeric transform of \(A\). At larger number \(N\) of grid points, the strong zooming-in is necessary to see the deviation.

A little bit more complicated example with self-Fourier function
 * \(\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (20) ~ ~ ~ A(x)\!=\!\exp(-x^2/2)(−3x^2+x^4)\)

is shown in figure at right in the same notations, as in previous figure, but \(N\!=\!16\). In this case, the segmented lines for arrays \(A\) and \(B\) overlap so perfectly, that the deviation is not seen. This deviation \(A\!-\!B\) scaled with factor 100 is shown with the saw-like line along the abscissa axis. The deviation is of order of \(10^{-3}\), and only in the 0th node of the grid it is a little bit larger, than \(10^{-3}\).

The examples show, that the numerical approximation of the Fourier operator is much better, than the 1st order approximation of the function and its Fourier transform with the arrays \(A\) and \(B\), the deviation of the segmented line from the smooth Gaussian is much larger than the deviation of the initial array from its discrete transform.

The figures are clickable, and the generators are copypasted in the descriptions. Everyone may load the code and play with the discrete implementation of the Fourier operator for various input arrays, approximating various functions.

Keywords
CosFT, FFT, SinFT