数字图像处理笔记（五）

傅里叶变换

为什么要在频域上研究图像增强

• 可以利用频率成分和图像外表之间的对应关系。一些在空间域表述困难的增强任务，在频率域中变得非常普通
• 滤波在频率域更为直观，他可以解释空间域滤波的某些性质
• 可以在频率域指定滤波器，做反变换，然后在空间域使用结果滤波器作为空间域滤波器的指导
• 一旦通过频率域试验选择了空间滤波，通常实施都在空间域进行

  频域滤波与空间滤波的关系

  傅里叶变换可以将图像从空间域变换到频率域，而傅里叶反变换则可以将图像的频谱逆变换维空域图像，即人可以直接识别的图像。这样一来，我们可以利用空域图像与频谱之间的对应关系，尝试将空间域滤波变为频域滤波，然后再将频域滤波处理后的图像反变换回空间域，从而达到图像增强的目的。这样做的一个最主要的吸引力在于频域滤波的直观性特点。

  傅里叶变换及其反变换

1. 一维连续傅里叶变换及其反变换

• 单变量连续函数f(x)的傅里叶变换$F(\mu)$定义为 $$F(\mu)=\int_{-∞}^{∞}f(x)e^{-j2\pi\mu x}dx$$ 其中，$j=\sqrt{-1}$
• 给定$F(\mu)$,通过傅里叶反变换可以得到f(x) $$f(x)=\int_{-∞}^{∞}F(\mu)e^{j2\pi\mu x}d\mu$$

1. 二维连续傅里叶变换及其反变换

• 二维连续函数f(x,y)的傅里叶变换$F(\mu,v)$定义为 $$F(\mu,v)=\int_{-∞}^{∞}\int_{-∞}^{∞}f(x,y)e^{-j2\pi(\mu x+vy)}dxdy$$
• 给定$F(\mu,v)$,通过傅里叶反变换可以得到f(x,y) $$f(x,y)=\int_{-∞}^{∞}\int_{-∞}^{∞}F(\mu,v)e^{j2\pi(\mu x+vy)}d\mu dv$$

1. 一维离散傅里叶变换(DFT)及反变换

• 单变量离散函数f(x)(x=0,1,2,…,M-1)的傅里叶变换$F(\mu)$定义为 $$F(\mu)=\frac{1}{M}\sum_{x=0}^{M-1}f(x)e^{-j2\pi\mu x/M}$$ $\mu=0,1,2,…,M-1$
• 给定$F(\mu)$,通过傅里叶反变换可以得到f(x) $$f(x)=\sum_{\mu=0}^{M-1}F(\mu)e^{j2\pi\mu x/M}$$ x=0,1,2,…,M-1

从欧拉公式$e^{j\theta}=cos\theta+jsin\theta$

$$F(\mu)=\frac{1}{M}\sum_{x=0}^{M-1}f(x)e^{j(-2\pi\mu x)/M}$$

傅里叶变换的极坐标表示：

$$F(\mu)=|F(\mu)|e^{-j\phi(\mu)}$$

• 幅度或频率谱为 $$|F(\mu)|=[R(\mu)^2+I(\mu)^2]^{(\frac{1}{2})}$$ $R(\mu)和I(\mu)$分别是$F(\mu)$的实部和虚部
• 相角或相位谱为 $$\phi(\mu)=arctan[\frac{I(\mu)}{R(\mu)}]$$

1. 二维DFT的极坐标表示

• 功率谱 $$P(\mu,v)=|F(\mu,v)|^2=R(\mu,v)^2+I(\mu,v)^2$$

F(\mu,v)的原点变换

用$(-1)^{x+y}乘以f(x,y),将F(\mu,v)原点变换到频率坐标下的(M/2,N/2),它是M×N区域的中心$
$\mu=0,1,2,…,M-1 $
$v=0,1,2,….,N-1$

• F(0,0)表示 $$F(0,0)=\frac{1}{MN}\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x,y)$$ 这说明：假设f(x,y)是一幅图像，在原点的傅里叶变换等价于图像的平均灰度级。
• 如果f(x,y)是实函数，它的傅里叶变换是对称的，即 $$F(\mu,v)=F(-\mu,-v)$$
• 傅里叶变换的频率谱是对称的 $$|F(\mu,v)|=|F(-\mu,-v)|$$

  傅里叶变换的性质

1. 平移性质 $$f(x,y)e^{j2\pi(\mu_0x/M+v_0y/N)}\Leftrightarrow F(\mu-\mu_0,v-v_0)$$ $$f(x-x_0,y-y_0)\Leftrightarrow F(\mu,v)e^{-j2\pi(\mu x_0/M+vy_0/N)}$$ 公式（1）表明将f(x,y)与一个指数项相乘就相当于把其变换后的频域中心移动到新的位置
   公式(2)表明将$F(\mu,v)$与一个指数项相乘就相当于把其变换后的空域中心移动到新的位置
   公式2表明f(x,y)的平移不影响其傅里叶变换的幅值

当$\mu_0=M/2且v_0=N/2,$

$$e^{j2\pi(\mu_0x/M+v_0y/N)}=e^{j\pi(x+y)}=(-1)^{x+y}$$

带入(1)和(2),得到

$$f(x,y)(-1)^{x+y}\Leftrightarrow F(\mu-M/2,v-N/2)$$ $$f(x-M/2,y-N/2)\Leftrightarrow F(\mu,v)(-1)^{\mu+v}$$

1. 分配律
   傅里叶变换只对加法满足分配律，但对乘法则不满足
2. 尺度变换(缩放)
   给定2个标量a和b，可以证明对傅里叶变换下列2个公式成立 $$af(x,y)\Leftrightarrow aF(\mu,v)$$ $$f(ax,by)\Leftrightarrow \frac{F(\mu/a,v/b)}{|ab|}$$

3. 旋转性
   引入极坐标$x=rcos\theta$,$y=rsin\theta$,$\mu=wcos\phi,v=wsin\phi$
   将f(x,y)和$F(\mu,v)$转换为$f(r,\theta)$和$F(w,\phi)$.将其带入傅里叶变换得到
   f(r,\theta+\theta_0)\Leftrightarrow F(w,\phi+\theta_0)
   f(x,y)旋转角度$\theta_0$,F(\mu,v)也将转过相同的角度
   $F(\mu,v)$旋转角度$\theta_0$,f(x,y)也将转过相同的角度

4. 周期性和共轭对称性

   $$F(\mu,v)=F(\mu+M,v)=F(\mu,v+N)=F(\mu+M,v+N)$$ $$f(x,y)=f(x+M,y)=f(x,y+N)=f(x+M,y+N)$$

   尽管$F(\mu,v)$对无穷多个$\mu$和v的值重复出现，但只需根据在任一个周期里的N个值就可以从$F(\mu,v)$得到f(x,y)

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

$$F(\mu,v)=F^{* }(-\mu,-v)$$ $$|F(\mu,v)|=|F(-\mu,-v)|$$

其中$F^{ }(-\mu,-v)$为$F(\mu,v)$的复共轭。
*当两个复数实部相等，虚部互为相反数时，这两个复数叫做互为共轭复数
周期性和共轭对称性举例

1. 平均值
   由二维傅里叶变换的定义 $$F(\mu,v)=\frac{1}{MN}\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x,y)e^{-j2\pi(\mu x/M+vy/N)}$$ 所以$F(0,0)=\frac{1}{MN}\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x,y)$

而$\overline{f}(x,y)=\frac{1}{MN}\sum_{x=0}^{M-1}\sum_{y=0}^{N-1}f(x,y)$

1. 可分性 $$F(\mu,v)=\frac{1}{M}\sum_{x=0}^{M-1}e^{-j2\pi\mu x/M}\frac{1}{N}\sum_{y=0}^{N-1}f(x,y)e^{-j2\pi vy/N}=\frac{1}{M}\sum_{x=0}^{M-1}e^{-j2\pi\mu x/M}F(x,v)$$

F(x,v)是沿着f(x,y)的一行所进行的傅里叶变换。当x=0,1,…,M-1,沿着f(x,y)的所有行计算傅里叶变换

分离性——二维傅里叶变换的全过程

1. 卷机性
   大小为M×N的两个函数f(x,y)和h(x,y)的离散卷积 $$f(x,y)* h(x,y)=\frac{1}{MN}\sum_{m=0}^{M-1}\sum_{n=0}^{N-1}f(m,n)h(x-m,y-n)$$

卷积定理

$$f(x,y)* h(x,y)\Leftrightarrow F(\mu,v)H(\mu,v)$$ $$f(x,y)h(x,y)\Leftrightarrow F(\mu,v)* H(\mu,v)$$

1. 相关性
   大小为M×N的两个函数f(x,y)和h(x,y)的相关定义为 $$f(x,y)。h(x,y)=\frac{1}{MN}\sum_{m=0}^{M-1}\sum_{n=0}^{N-1}f^{* }(m,n)h(x+m,y+n)$$

$f^ $表示f的复共轭。对于实函数，$f^{ }=f$
相关定理：

• 卷积和相关性理论总结
  卷积是空间域过滤和频率域过滤之间的纽带
  相关的重要应用在于匹配:确定是否有感兴趣的物体区域

1. f(x,y)是原始图像
2. h(x,y)作为感兴趣的物体或区域(模板)
3. 如果匹配,两个函数的相关值会在h找到f中相应点的位置上达到最大

快速傅里叶变换

只考虑一维的情况，根据傅里叶变换的分离性可知，二维傅里叶变换可由连续2次一维傅里叶变换得到。

为什么需要快速傅里叶变换(FFT)

$$F(\mu)=\frac{1}{M}\sum_{x=0}^{M-1}f(x)e^{-j2\pi\mu x/M} \ \ \ \ \ \mu=0,1,2,...,M-1$$

对$\mu$的M个值中的每一个都需进行M次复数乘法(将f(x)与$e^{-j2\pi\mu x/M}$相乘)和M-1次加法，即时间复杂度为$O(M^2)$.
而快速傅里叶变换FFT时间复杂度为$O(Mlog_2M)$

FFT算法基本思想

FFT算法基于一个叫做逐次加倍的方法。通过推导将原始傅里叶变换成两个递推公式

$$F(\mu)=\frac{1}{M}\sum_{x=0}^{M-1}f(x)e^{-j2\pi\mu x/M} \ \ \ \ \ \mu=0,1,2,...,M-1$$ $$F(\mu)=\frac{1}{2}[F_{even}(\mu)+F_{odd}(\mu)W_{2k}^{\mu}]$$ $$F(\mu+K)=\frac{1}{2}[F_{even}(\mu)-F_{odd}(\mu)W_{2k}^{\mu}]$$

其中M=2K
$F_{even}(\mu)、F_{odd}(\mu)$是K个点的傅里叶值

假设M的形式是$M=2^n$
n为正整数。因此，M可以表示为M=2K
将M=2K带入上面式子有

$$F(\mu)=\frac{1}{2K}\sum_{x=0}^{2K-1}f(x)W_{2K}^{\mu x}=\frac{1}{2}[\frac{1}{K}\sum_{x=0}^{K-1}f(2x)W_{2K}^{\mu(2x)}+\frac{1}{K}\sum_{x=0}^{K-1}f(2x+1)W_{2K}^{\mu(2x+1)}]$$

推导：
因为 $W_M=e^{-j2\pi/M}$
所以

$$W_{2K}^{2\mu x}=e^{-j2\pi(2\mu x)/2K}=e^{-j2\pi(\mu x)/K}=W_{K}^{\mu x}$$

代入上式有

$$F(\mu)=\frac{1}{2}[\frac{1}{K}\sum_{x=0}^{K-1}f(2x)W_K^{\mu x}+\frac{1}{K}\sum_{x=0}^{K-1}f(2x+1)W_K^{\mu x}W_{2K}^{\mu}]$$

定义两个符号

$$F_{even}(\mu)=\frac{1}{K}\sum_{x=0}^{K-1}f(2x)W_K^{\mu x} \ \ \ \ \ \mu=0,1,2,...,K-1$$ $$F_{odd}(\mu)=\frac{1}{K}\sum_{x=0}^{K-1}f(2x+1)W_K^{\mu x}$$

得到FFT的第一个公式

$$F(\mu)=\frac{1}{2}[F_{even}(\mu)+F_{odd}(\mu)W_{2K}^{\mu}]$$

推导:

$$W_K^{\mu+K}=e^{-j2\pi(\mu+K)/K=e^{-j2\pi\mu/K}e^{-j2\pi}}=W_K^{\mu}$$ $$W_{2K}^{\mu+K}=e^{-j2\pi(\mu+K)/2K}=e^{-j2\pi\mu/2K}e^{-j\pi}=-W_{2K}^{\mu}$$

由此可得

$$F(\mu+K)=\frac{1}{2K}\sum_{x=0}^{2K-1}f(x)W_{2K}^{(\mu+K)x}=\frac{1}{2}[F_{even}(\mu)-F_{odd}(\mu)W_{2K}^{\mu}]$$

由此可以得到FFT的第二个公式

$$F(\mu+K)=\frac{1}{2}[F_{even}(\mu)-F_{odd}(\mu)W_{2K}^{\mu}]$$