The mathematics behind Fourier Transform for Image Processing I am following http://homepages.inf.ed.ac.uk/rbf/HIPR2/fourier.htm . I understand the application of Fourier Transform behind Image Processing, but right now, I am curious about the mathematics behind it, and it is giving me a bit of a hard time.
For example:

In this formula, where do all these equations come from? Could somebody please elaborate the mathematics behind the scene in layman's term?
 A: In my opinion it is quite okay to understand Fourier transform as orthogonal basis matrices to evaluate the certain frequencies for a given image. However I have the following link which will be helpful for you to further understand:
http://sharp.bu.edu/~slehar/fourier/fourier.html#harmonics (dead link)
The working one
Good luck.
A: The discrete Fourier transform (which is a function from $\mathbb C^N$ to $\mathbb C^N$) simply changes basis, to a special basis known as the discrete Fourier basis.
And what's so special about this basis?  It's a basis of eigenvectors for the shift operator $S$, which maps $\begin{bmatrix} x_0 & x_1 & x_2 &\ldots & x_{N-1} \end{bmatrix}^T$ to $\begin{bmatrix} x_{N-1} & x_0 & x_1& \ldots & x_{N-2} \end{bmatrix}^T$.
Each basis vector is constructed by taking an $N$th root of unity $\omega$ and forming the vector $v = \begin{bmatrix} 1 & \omega & \omega^2 & \ldots & \omega^{N-1} \end{bmatrix}^T$.  It's easy to check that $v$ is an eigenvector of $S$:
\begin{align*}
S(v) &= \begin{bmatrix} \omega^{N-1} & 1 & \omega & \ldots & \omega^{N-2} \end{bmatrix}^T \\
&= \omega^{-1} v.
\end{align*}
We have one discrete Fourier basis vector for each $N$th root of unity.
Because $S$ is unitary (it clearly preserves norms), we have that $S$ is normal, so there is an orthonormal basis of eigenvectors for $S$.  This explains why the discrete Fourier basis (once normalized) is orthonormal and why the discrete Fourier transform preserves norms.
Any convolution operator on $\mathbb C^N$ can be expressed as a linear combination of powers of $S$.  This explains why the discrete Fourier basis diagonalizes any convolution operator on $\mathbb C^N$.
The 2D discrete Fourier transform is analogous.  We have two shift operators on $\mathbb C^{M \times N}$, $S_1$ (which shifts each row to the right) and $S_2$ (which shifts each column down).  $S_1$ and $S_2$ are normal operators and they commute, so there is an orthonormal basis of eigenvectors that simultaneously diagonalizes $S_1$ and $S_2$.  This gives us the 2D discrete Fourier basis.
A: I think the easiest way of understanding the math, Is understanding what it means. It's also easier to start in 1D and only then move on to higher dimensions.
To see what the sum you mentioned is, try inserting a wave with constant frequency. You see that the summation only gives a peak at the "correct" frequency. If you have a super-position of waves, the sum will give you a number of peaks each corresponding to a different part of the super-position.
Thus, the sum gives you the spectrum of frequencies in the wave or image.
A: In layman's terms:
First let's start with the Fourier series, a method that Fourier wrote for the first time in a paper about heat diffusion modelling. The idea is that any continuous function can be approximated by adding up lots of sine an cosine functions. The more terms you use, the more accurate will be the approximation:
$f(x) = a_0\cos\frac{\pi y}{2}+a_1\cos 3\frac{\pi y}{2}+a_2\cos5\frac{\pi y}{2}+\cdots + b_0\sin\frac{\pi y}{2}+b_1\sin 3\frac{\pi y}{2}+b_2\sin5\frac{\pi y}{2}+\cdots$
In order to find the coefficients, the following trick was used for cosine functions:
$$
a_n = \displaystyle\frac{1}{\pi}\int_{-\pi}^\pi f(x) \cos(nx)\, dx,
$$
and similarly for sine functions:
$$
b_n = \displaystyle\frac{1}{\pi}\int_{-\pi}^\pi f(x) \sin(nx)\, dx.
$$
These formulae are known as the Fourier sine and cosine transforms. 
Euler's formula $e^{i\theta} = \cos(\theta) + i \sin(\theta)$ shows a relationship between exponential and cosine and sine functions. The Fourier sine and cosine transforms can thus be combined into a single transform:
$$
\hat{f}(\xi) = \int_{-\infty}^{\infty} f(x)\ e^{- 2\pi i x \xi} \, dx,
$$
and this explains why you see exponential functions in the Fourier transform instead of cosine and sine functions.
Now, in image processing we are typically working with two-dimensional images. So the Fourier transform has been done twice:
$$
\displaystyle F(u,v)=\int_{-\infty}^{\infty} \int_{-\infty}^{\infty} f(x,y)\ e^{-j2\pi(ux+vy)} \, dx \, dy. 
$$
Since images actually come in discrete pixels rather than continuous the integrals are replaced by summations. 
A: I have given a similar explanation of Fourier series here and here. 
The Fourier series of a function $f : \mathbb R / \mathbb Z \to \mathbb C$ is a sum $\sum_{k = -\infty}^\infty \hat{f}(e^{2 \pi i k x}) e^{2 \pi i k x}$. Here, $e^{2 \pi i k x} = \chi_k (x)$ denotes the $k$-th character of the topological group $\mathbb R / \mathbb Z $ and $\hat{f} : \widehat{\mathbb R / \mathbb Z} \to \mathbb C$ denotes the Fourier transform of $f$. 
The setting is: You have a topological group $G$, then you consider a space of functions $G \to \mathbb C$ on it, say for example $L^2(G)$. This is a Hilbert space and comes with an inner product $\langle \cdot , \cdot \rangle$ which lets you define orthogonality of functions in the space, and characters in particular. Those characters form an orthonormal basis for your functions which means that you can approximate every function in $L^2(G)$ as $\| f - \sum_{k = -N}^N \hat{f}(\chi_k) \chi_k (x) \|_2 \xrightarrow{N \to \infty} 0$.
