[数字图像处理笔记] 第四章 图像变换

1. 连续傅里叶变换

1.1 一维连续傅里叶变换

• 一维连续傅里叶正变换（\(\text{1-Dimensional Continuous Fourier Transform}\)）

  对于函数 \(f(x)\)，一维连续傅里叶变换有如下定义：

  \[\Re: \; F(u) = \int_{- \infty}^{\infty} f(x)e^{-j2 \pi u x} \text{d}x \]

  其中 \(j^2 = 1\).

• 一维连续傅里叶反变换（\(\text{1-Dimensional Continuous Fourier Inversion Transform}\)）

  上式中 \(F(u)\) 的反变换，有如下定义：

  \[\Re ^{-1}: \; f(x) = \int_{- \infty}^{\infty} F(u)e^{j2 \pi u x} \text{d}u \]

函数 \(f(x)\) 与 \(F(u)\) 被称为一个 傅里叶变换对。对于任意函数 \(f(x)\)，其傅里叶变换 \(F(u)\) 是 唯一的；反之亦然。

对于函数 \(f(x)\)，其对应的 \(F(u)\) 一般是复函数。由 欧拉公式 \(e^{ix} = \text{cos}x + i\;\text{sin}x\) 有如下定义：

• 实部

  \[R(u) = \int_{- \infty}^{\infty} f(x)\text{cos}(2\pi u x)\;\text{d}x \]

• 虚部

  \[I(u) = - \int_{- \infty}^{\infty} f(x)\text{sin}(2\pi u x)\;\text{d}x \]

• 振幅

  \[|F(u)| = [R^2(u) + I^2(u)]^{\frac{1}{2}} \]

• 能量

  \[E(u) = |F(u)|^2 = R^2(u) + I^2(u) \]

• 相位

  \[\varphi (u) = \text{arctan}\frac{I(u)}{R(u)} \]

1.2 二维连续傅里叶变换

设 \(f(x, y)\) 连续可积，且 \(F(u, v)\) 可积，则存在傅里叶变换对：

• 正变换

  \[F\{f(x, y)\} = F(u, v) = \int_{- \infty}^{\infty}f(x, y)e^{-j2\pi (ux + vy)}\; \text{d}x \text{d}y \]

• 反变换

  \[F^{-1}\{F(u, v)\} = f(x, y) = \int_{- \infty}^{\infty}F(u, v)e^{j2\pi (ux + vy)}\; \text{d}u \text{d}v \]

同样有如下定义：

• 傅里叶频谱

  \[|F(u, v)| = [R^2(u, v) + I^2(u, v)]^{\frac{1}{2}} \]

• 能量谱

  \[E(u, v) = |F(u, v)|^2 = R^2(u, v) + I^2(u, v) \]

• 相位谱

  \[\varphi (u, v) = \text{arctan}\frac{I(u, v)}{R(u, v)} \]

2. 离散傅里叶变换

2.1 一维离散傅里叶变换

\(f(x)\) 是在时域上等距采样得到的 \(N\) 点离散序列， \(x\) 是离散实变量，\(u\) 为离散频率变量。由此，离散傅里叶变换（\(\text{Discrete Fourier Transform, DFT}\)）定义如下：

• 正变换

  \[\Re: \; F(u) = \frac{1}{\sqrt{N}}\sum_{x=0}^{N-1}f(x)e^{-j2 \pi \frac{ux}{N}}, \; u = 0, 1, ... ,N - 1 \]

• 反变换

  \[\Re: \; f(x) = \frac{1}{\sqrt{N}}\sum_{x=0}^{N-1}F(u)e^{j2 \pi \frac{ux}{N}}, \; u = 0, 1, ... ,N - 1 \]

由欧拉公式，可以得到如下定义：

• 傅里叶频谱

  \[|F(u)| = [R^2(u) + I^2(u)]^{\frac{1}{2}} \]

• 能量谱

  \[E(u) = |F(u)|^2 = R^2(u) + I^2(u) \]

• 相位谱

  \[\varphi (u) = \text{arctan}\frac{I(u)}{R(u)} \]

2.2 二维离散傅里叶变换

同连续函数的傅里叶变换一样，离散函数的傅里叶变换也可推广到二维的情形，其二维离散傅里叶变换定义为：

• 正变换

  \[F(u, v) = \frac{1}{\sqrt{MN}}\sum_{x=0}^{M-1} \sum_{y=0}^{N-1} f(x, y) e^{-j 2 \pi (\frac{ux}{M} + \frac{vy}{N})}, \; u = 0, 1, ... , M - 1, \; v = 0, 1, ... , N - 1 \]

• 反变换

  \[f(x, y) = \frac{1}{\sqrt{MN}}\sum_{x=0}^{M-1} \sum_{y=0}^{N-1} F(u, v) e^{j 2 \pi (\frac{ux}{M} + \frac{vy}{N})}, \; u = 0, 1, ... , M - 1, \; v = 0, 1, ... , N - 1 \]

同样有如下定义：

• 傅里叶频谱

  \[|F(u, v)| = [R^2(u, v) + I^2(u, v)]^{\frac{1}{2}} \]

• 能量谱

  \[E(u, v) = |F(u, v)|^2 = R^2(u, v) + I^2(u, v) \]

• 相位谱

  \[\varphi (u, v) = \text{arctan}\frac{I(u, v)}{R(u, v)} \]

频谱：图像的 主要成分是低频信息，它形成了图像的 基本灰度等级，对图像结构的决定作用较小；

中频信息 决定了图像的基本结构，形成了图像的 主要边缘结构；

高频信息 形成了图像的 边缘和细节，是在中频信息上对图像内容的进一步强化。

3. 快速傅里叶变换

3.1 FFT的作用

快速傅里叶变换 （\(\text{Fast Fourier Transform, FFT}\)）并不是一种新的变换，它是 离散傅里叶变换（\(\text{DFT}\)） 的一种算法。

这种方法是在分析 离散傅里叶变换（\(\text{DFT}\)） 中的多余运算的基础上，进而消除这些重复工作的思想指导下得到的，所以在运算中大大节省了工作量，达到了快速的目的。

对于一个有限长序列 \(\{f(x)\}, 0 \le x \le n - 1\)，其傅里叶变换如下：

\[F(u) = \frac{1}{\sqrt{n}} \sum_{x=0}^{n - 1}f(x) w_n^{ux} \]

其中 \(W_n = e^{-j\frac{2 \pi}{n}}\)

计算每个频率的分量，需要进行 \(n\) 次运算，计算所有频率的分量，则需要进行 \(O(n^2)\) 的时间复杂度。当 \(n\) 较大时，这显然会耗费很多时间。

而 快速傅里叶变换 (\(\text{FFT}\)) 可以做到在 \(O(n\text{log}n)\) 的时间复杂度，计算傅里叶变换对应的多项式。

3.2 单位根

以单位圆点为起点，单位圆的 \(n\) 等分点为终点，在单位圆上可以得到 \(n\) 个复数，设幅角为正且最小的复数为 \(w_n\)，称为 \(n\) 次 单位根，即:

\[w_n = \text{cos}\frac{2 \pi}{n} + i \;\text{sin} \frac{2 \pi}{n} \]

由欧拉公式可得：

\[w_n^k = \text{cos}\frac{2 k \pi}{n} + i \;\text{sin} \frac{2 k \pi}{n} \]

特别地，\(w^0_n = w^n_n = 1\)

由周期性可以得到如下性质：

\[\begin{split} w_{cn}^{ck} &= \text{cos}\frac{2 ck \pi}{cn} + i \;\text{sin} \frac{2 ck \pi}{cn} \\\\ &= w^k_n \\\\\\ w_{n}^{k + \frac{n}{2}} &= w^k_n \cdot [\text{cos}(\frac{2 \pi}{n} \cdot \frac{n}{2}) + i \;\text{sin} (\frac{2 \pi}{n} \cdot \frac{n}{2})] \\\\ &= w^k_n \cdot (\text{cos} \pi + i \; \text{sin} \pi) \\\\ & = - w^k_n \end{split} \]

3.3 快速傅里叶变换

对于傅里叶变换计算表达式

\[F(k) = \frac{1}{\sqrt{n}} \sum_{x=0}^{n - 1}f(x) w_n^{kx} \]

\[w_n^{kx} = e^{-j\frac{2 \pi}{n}kx} \]

我们令 \(t = w_n^k, a_x = f(x)\)

由此可化简为：

\[F(k) = \frac{1}{\sqrt{n}} \sum_{x=0}^{n - 1}a_x t^k \]

展开傅里叶变换表达式：

\[\begin{split} F(t) &= (a_0 + a_2 t^2 + ... + a_{n-2}t^{n-2}) \\ &+ (a_1t + a_3t^{3} + ... + a_{n-1}t^{(n-1)}) \end{split} \]

令

\[\begin{split} G(t) = a_0 + a_2t + a_4t^{2} + ... + a_{n-2}t^{\frac{n}{2} - 1} \\\\ H(t) = a_1 + a_3t + a_5t^{2} + ... + a_{n-1}t^{\frac{n}{2} - 1} \end{split} \]

由此可得：

\[F(t) = G(t^2) + t\; H(t^2) \]

将 \(w^k_n = t\) 回代到等式中，由 单元根 性质可得：

\[\begin{split} F(w^k_n) &= G(w_n^{2k}) + w^k_n \cdot H(w_n^{2k}) \\\\ &= G(w^k_{\frac{n}{2}}) + w^k_n \cdot H(w^k_{\frac{n}{2}}) \\\\\\ F(w^{k + \frac{n}{2}}_n) &= G(w^{2k + n}_n) + w^{k + \frac{n}{2}}_n \cdot H(w_n^{2k + n}) \\\\ &= G(w^{2k}_n \cdot w^n_n) - w_n^k \cdot H(w^{2k}_n \cdot w^n_n) \\\\ &= G(w^{2k}_n) - w^k_n \cdot H(w^{2k}_n) \\\\ &= G(w_{\frac{n}{2}}^k) - w^k_n \cdot H(w_{\frac{n}{2}}^k) \end{split} \]

可以发现，只要求出 \(G(w_{\frac{n}{2}}^k)\) 和 \(H(w_{\frac{n}{2}}^k)\) 就可以进而求出 \(F(w^k_n)\) 和 \(F(w_n^{k + \frac{n}{2}})\)。

对于 \(G(w_{\frac{n}{2}}^k)\) 和 \(H(w_{\frac{n}{2}}^k)\) 的求解，只需要转换求解为 \(G(w^k_{\frac{n}{4}})\) 和 \(H(w^k_{\frac{n}{4}})\) 的问题即可。

直到最后，转换为求解 \(G(w^k_1) = 1\) 和 \(H(w_1^k) = 1\) 的问题。

所以，只需要不断对 \(G, H\) 进行 递归 求解即可。 求解每个分量对应的 \(F(w_n^k)\) 只需要 \(O(\text{log} n)\) 的时间复杂度。求解所有分量则需要 \(O(n \; \text{log} n)\) 的时间复杂度。

涉猎过算法竞赛的都知道，\(O(n \; \text{log} n)\) 比 \(O(n ^ 2)\) 在 \(n\) 较大时会快很多。

4. 傅里叶变换的性质

4.1 可分离性

二维离散傅里叶变换对，可以写成如下形式：

\[\begin{split} & F(u, v) = \frac{1}{\sqrt{MN}}\sum_{x=0}^{M-1}e^{-j 2\pi \frac{ux}{M}}f(x, y)\sum_{y=0}^{N-1}e^{-j 2\pi \frac{vy}{N}} \\\\ & f(x, y) = \frac{1}{\sqrt{MN}}\sum_{u=0}^{M-1}e^{j2 \pi \frac{ux}{M}} \sum_{v=0}^{N-1} F(u, v) e^{j 2\pi \frac{vy}{N}} \\\\ & u = 0, 1, ... , M-1; \; v = 0, 1, ... , N - 1 \end{split} \]

一个 \(\text{2-D}\) 傅里叶变换可由连续两次运用 \(1-D\) 傅里叶变换来实现：

\[\begin{split} F(u, y) &= \frac{1}{\sqrt{M}}\sum_{x=0}^{M-1}f(x, y) e^{-j 2 \pi \frac{ux}{M}} \\\\ F(u, v) &= \frac{1}{\sqrt{N}}\sum_{y=0}^{N-1}F(u, y) e^{-j 2 \pi \frac{vy}{N}} \end{split} \]

当 \(M = N\) 时，变换过程如图所示：

4.2 平移性质

如果 \(F(u,v)\) 的频率变量 \(u\)，\(v\) 各移动了 \(u_0\)，\(v_0\) 距离，\(f(x, y)\) 的实变量 \(x\)，\(y\) 各移动了 \(x_0\)，\(y_0\)距离，则：

\[\begin{split} & f(x, y)e^{j 2 \pi (\frac{u_0x}{M} + \frac{v_0y}{N})} \Leftrightarrow F(u - u_0, v - v_0) \\\\ & f(x - x_0, y - y_0) \Leftrightarrow F(u, v)e^{-j 2 \pi (\frac{ux_0}{M} + \frac{vy_0}{N})} \end{split} \]

特别地，当 \(u_0 = \frac{M}{2}, v_0 = \frac{N}{2}\) 时，\(e^{j 2 \pi (\frac{u_0x}{M} + \frac{v_0y}{N})} = e^{j \pi (x + y)} = (-1)^{x + y}\)

则：

\[\begin{split} & f(x, y)(-1)^{x+y} \Leftrightarrow F(u - \frac{M}{2}, v - \frac{N}{2}) \\\\ & f(x - \frac{M}{2}, y - \frac{N}{2}) \Leftrightarrow F(u, v)(-1)^{u + v} \end{split} \]

这条性质被称为 移中性 。

如图所示，图像平移 不会改变图像的频谱形状（fftshift后），但是会 改变相位：

相关代码如下：

I = imread('3.png');
I = rgb2gray(I);
I = im2double(I);
F = fftshift(fft2(I));

subplot(2,3,1), imshow(I,[]), title('Original Image');
subplot(2,3,2), imshow(log(1+abs(F)),[]), title('Fourier Transform - Magnitude');
subplot(2,3,3), imshow(angle(F),[-pi pi]), title('Fourier Transform - Phase');

translation = [400, 400];
I_translate = circshift(I, translation);
F_translate = fftshift(fft2(I_translate));

subplot(2,3,4), imshow(I_translate,[]), title('Translated Image');
subplot(2,3,5), imshow(log(1+abs(F_translate)),[]), title('Fourier Transform of Translated Image - Magnitude');
subplot(2,3,6), imshow(angle(F_translate),[-pi pi]), title('Fourier Transform of Translated Image - Phase');

4.3 周期性和共轭对称性

傅里叶变换和反变换的周期性：

\[\begin{split} F(u, v) &= F(u + M, v) = F(u, v + N) = F(u + M, v + M) \\\\ f(x, y) &= f(x + M, y) = f(x, y + N) = f(x + M, y + N) \end{split} \]

如果 \(f(x, y)\) 是实函数，则它的傅里叶变换具有共轭对称性：

\[\begin{split} F(u, v) &= F^*(-u, -v) \\\\ |F(u, v)| &= |F(-u, -v)| \end{split} \]

其中 \(F^*(u, v)\) 为 \(F(u, v)\) 的 复共轭 。

4.4 旋转性质

借助极坐标 \(x = r \text{cos} \theta, y = r\text{sin}\theta, u = w\text{cos}\phi, v = w\text{sin}\phi\)：

\[\begin{split} & f(x, y) \Leftrightarrow F(u, v) \\\\ & f(r\text{cos}\theta, r \text{sin}\theta) \Leftrightarrow F(w \text{cos}\phi, w \text{sin}\phi) \\\\ & f(r, \theta + \theta_0) \Leftrightarrow F(w, \phi + \theta_0) \end{split} \]

如图所示，旋转之后，傅里叶频谱也随之旋转，相位也发生了改变。会 改变图像的频谱，也会 改变相位：

相关代码如下：

I = imread('3.png');
I = rgb2gray(I);
I = im2double(I);
F = fftshift(fft2(I));

subplot(2,3,1), imshow(I,[]), title('Original Image');
subplot(2,3,2), imshow(log(1+abs(F)),[]), title('Fourier Transform - Magnitude');
subplot(2,3,3), imshow(angle(F),[-pi pi]), title('Fourier Transform - Phase');

I_translate = imrotate(I, -45, 'bilinear', 'crop');
F_translate = fftshift(fft2(I_translate));

subplot(2,3,4), imshow(I_translate,[]), title('Translated Image');
subplot(2,3,5), imshow(log(1+abs(F_translate)),[]), title('Fourier Transform of Translated Image - Magnitude');
subplot(2,3,6), imshow(angle(F_translate),[-pi pi]), title('Fourier Transform of Translated Image - Phase');

4.5 分配律

傅里叶变换和反变换对 加法满足分配律，根据傅里叶变换对的定义可得到：

\[F\{f_1(x, y) + f_2(x, y)\} = F\{f_1(x, y)\} + F\{f_2(x, y)\} \]

乘法则不满足 ：

\[F\{f_1(x, y) \cdot f_2(x, y)\} \not = F\{f_1(x, y)\} \cdot F\{f_2(x, y)\} \]

4.6 尺度变换

尺度变换描述了函数自变量的尺度变化对其傅里叶变换的作用，有如下表达式：

\[af(x, y) \Leftrightarrow aF(u, v) \]

可证明得到：

\[f(ax, by) \Leftrightarrow \frac{1}{|ab|}F(\frac{u}{a}, \frac{v}{b}) \]

如图所示：

相关代码：

img = imread('img.jpg'); 
img_norm = im2double(img);

figure;
subplot(2, 2, 1);
imshow(img_norm);
title('initial image');

img_scaled = img_norm * 0.1;
subplot(2, 2, 2);
imshow(img_scaled);
title('scaled image');

fft_img = fft2(img_norm);
fft_img = fftshift(fft_img);
subplot(2, 2, 3);
imshow(log(1 + abs(fft_img)), []);
title('initial image fft');

fft_scaled_img = fft2(img_scaled);
fft_scaled_img = fftshift(fft_scaled_img);
subplot(2, 2, 4);
imshow(log(1 + abs(fft_scaled_img)), []);
title('scaled image fft');

4.7 平均值

对于一个二维离散函数，其平均值为：

\[\bar f(x, y) = \frac{1}{MN}\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x, y) \]

傅里叶变换在原点 \((0, 0)\) 的频谱分量：

\[\begin{split} F(0, 0) &= \frac{1}{\sqrt{MN}}\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x, y)e^{-j2 \pi 0} \\\\ &= \sqrt{MN}[\frac{1}{MN}\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x, y)] \\\\ &= \sqrt{MN} \bar f(x, y) \end{split} \]

4.8 卷积定理

• 一维傅里叶变换

  \[f(x) * g(x) = \int_{- \infty}^{\infty}f(z)g(x - z)\text{d}z \Leftrightarrow F(u)G(u) \]

• 二维傅里叶变换

  \[f(x, y) * g(x, y) \Leftrightarrow F(u, v)G(u, v) \]