频域滤波技术，目前主要使用的有两种类型，一种是傅立叶变换技术，还有一种是小波分析。基本逻辑就是把原始信号映射到频率空间中，使得在时域空间无法处理的信号，得以在另外一种空间体系下能够被较有效的处理。

目前，我已经完成了全部的傅立叶分析的科普性的文章，所以在这个系列里我就不再赘述这方面的基础知识，我会在这下面贴出这些文章的链接，你可以根据自己的需要去查阅对应的章节，如果想直接了解如何实现一个一维傅立叶FFT与iFFT算法的朋友，直接从6、7、8三个章节看就行了。

《漫谈傅里叶1——从无穷级数到傅里叶》
《漫谈傅里叶2——公式推导、三角函数正交性》
《漫谈傅里叶3——收敛性、非周期函数的推广应用》
《漫谈傅里叶4——全时傅里叶的缺点与短时傅里叶》
《漫谈傅里叶5——卷积与短时傅里叶的缺点》
《漫谈傅里叶6——采样与1D初步实现》
《漫谈傅里叶7——带有相位与幅值的1D实现》
《漫谈傅里叶8——傅里叶逆计算的实现》

而关于小波分析这部分，我打算近期更新一些关于量化金融方面的文章后再回头来写，你可以关注我的博客，或者新浪微博——打码的阿通。

文章目录

用OpenCV实现的图片DFT示例
- 正向 DFT 或 FFT 变换后的结果
- 逆向DFT 或从频率转回图像的结果
理想滤波器 (Ideal Filter)
Butterworth 滤波器

用OpenCV实现的图片DFT示例

尽管我知道Python有很多方便的FFT库，但是为了方便学习C/C++的朋友，我尽量用原生的OpenCV库做完了这个用例。理论的那一套在前面的一维FFT里已经讲过，所以这里不赘述，直接上完整的代码好了。

import cv2from Common import *def load_img_grayscale(file: str):image = cv2.imread(file, cv2.IMREAD_GRAYSCALE)return imagedef shift_spectrum(gamma_matrix):# 和C++语言的OpenCV库中获取矩阵大小方式类似# 可以用如下方法获得# int rows = gamma_matrix.rows# int cols = gamma_matrix.colsrows, cols = gamma_matrix.shapecr = int(rows / 2)cc = int(cols / 2)# 计算四个象限子矩阵# 可以使用C/C++中通常使用 Rect 来获取原矩阵的四个象限范围# (0,0)----------------------------->cols# Q2             | Q1                 |# (0r,0c)        |                    |  #                |(cr,cc)             |# ---------------|------------------- |# Q3             | Q4                 |#                |                    |#                |                    |#                                    (rows,cols)  quadrant_1 = gamma_matrix[0:cr, cc:cols]quadrant_2 = gamma_matrix[0:cr, 0:cc]quadrant_3 = gamma_matrix[cr:rows, 0:cc]quadrant_4 = gamma_matrix[cr:rows, cc:cols]temp = quadrant_1.copy()quadrant_1 = quadrant_3quadrant_3 = temptemp = quadrant_2.copy()quadrant_2 = quadrant_4quadrant_4 = temp# 将修改好的数值写入到output矩阵中output = np.zeros(gamma_matrix.shape)output[0:cr, cc:cols] = quadrant_1output[0:cr, 0:cc] = quadrant_2output[cr:rows, 0:cc] = quadrant_3output[cr:rows, cc:cols] = quadrant_4return outputdef gen_spectrum(real_matrix, imag_matrix):# 创建一个gamma变化后的矩阵gamma_matrix = np.zeros(real_matrix.shape, dtype=np.float32)cv2.magnitude(real_matrix, imag_matrix, gamma_matrix)cv2.log(gamma_matrix + 1, gamma_matrix)# 归一化：normalizeoutput = np.zeros(real_matrix.shape, dtype=np.float32)cv2.normalize(gamma_matrix, output, 0, 1, cv2.NORM_MINMAX)return outputdef do_fft():# 获取原始图像的灰度图img = load_img_grayscale("Data/Illustrations/1.jpg")# 获取图像的长宽# C++代码里，需要用 img.rows, img.cols来获取对应的数据rows, cols = img.shape# 对长宽进行优化，DFT算法需要输入的数组长度为2^nm = cv2.getOptimalDFTSize(rows)n = cv2.getOptimalDFTSize(cols)# 按照要求，对数据进行填充，不足的部分填充0# C++中使用# cv2.copyMakeBorder(src, dst, 0, m - rows, 0, n - cols, cv2.BORDER_CONSTANT)fft2 = np.zeros((m, n, 2), np.float32)fft2[:rows, :cols, 0] = img# 为DFT创建实数部和虚数部# C++对应的方法：# Mat planes[] = {Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F)};# Mat complex;# cv2::merge(planes, 2, complex):# dft(complexImg, complexImg, cv2.DFT_COMPLEX_OUTPUT)# 使用DFT算法，进行快速傅立叶变换cv2.dft(fft2, fft2, cv2.DFT_COMPLEX_OUTPUT)# 分离复数矩阵# C/C++中，用的是 split(fft2, real_mat, imag_mat) 进行分离real_mat, imag_mat = cv2.split(fft2)# 生成频率数据，由于一部分数据属于不可见数据，所以需要进行gamma变换# 并执行归一化操作gamma = gen_spectrum(real_mat, imag_mat)# 为gamma图执行中央转换# 方法是1，3象限对调，2，4象限对调shift = shift_spectrum(gamma_spectrum)# 把原始图像、FFT频率图、转化后的频率图全部显示出来plt = PltImageCache()plt.add(img, "origin")plt.add(gamma, "dft computed")plt.add(shift, "dft shifted")plt.plots(1, 3)# 返回复数矩阵供下一步操作return fft2def do_ifft(fft2):# 获取原始图像的灰度图img = load_img_grayscale("Data/Illustrations/1.jpg")ifft2 = np.zeros(fft2.shape[:2], np.float32)cv2.dft(fft2, ifft2, cv2.DFT_REAL_OUTPUT | cv2.DFT_INVERSE | cv2.DFT_SCALE)# 把原始图像、FFT频率图、转化后的频率图全部显示出来plt = PltImageCache()plt.add(img, "original")plt.add(ifft2, "idft")plt.plots(1, 2)if __name__ == "__main__":dft_matrix = do_fft()do_ifft(dft_matrix)

输出的结果就是这样的了：

正向 DFT 或 FFT 变换后的结果

逆向DFT 或从频率转回图像的结果

这里你会发现转换回来的照片多了一圈黑色的边框，那个就是为了进行DFT运算，而调整的新图片大小，长宽大小都是 2n2^n2n。

通常，在做完FFT运算后，它输出的结果是对称的，对于2维的也是这样，不做任何处理的频率图，频率范围从四边角到中央逐渐递增，所以为了方便观察，会把四个象限的子矩阵做一下调整。

所以最终的图片也就成了从中央一个点向四周呈现放射状的形式。因此频率范围的调整就变成了计算圆半径R的问题。

通常，我们为了过滤图片中的一些信息，就会使用一些滤波器，下面介绍一些比较常用的频域滤波器。

理想滤波器 (Ideal Filter)

类似于二极管的作用，只放走设定域值内的频率，而域值外的频率一律设置为0

f(x)={g(x)a<x≤b0elsef(x) = \left\{\begin{matrix} g(x) & a < x \leq b \\ 0 & else \\ \end{matrix}\right. f(x)={g(x)0a<x≤belse

而在设代码结构上与椒盐噪音的方式是相似的，其实现代码大概是这样的：

def generate_ideal_mask(rows, cols, alpha, beta=0):""":param rows: rows of spectrum:param cols: columns of spectrum:param alpha: the frequency upper limit:param beta: the frequency lower limit:return:"""import numpy as npoutput = np.zeros((rows, cols), np.float32)# center coordinatecr = rows / 2cc = cols / 2for r in range(rows):for c in range(cols):distance = np.sqrt((r - cr) ** 2 + (c - cc) ** 2)if alpha > distance >= beta:output[r, c] = 1.else:output[r, c] = 0.return output

而输出结果如下：

Butterworth 滤波器

其数学解析式如下：

f(x)=11+[D(u,v)Do]2nf(x) = \frac{1}{1 +[\frac{D(u,v)}{D_o}]^{2n}} f(x)=1+[DoD(u,v)]2n1

代码的实现方式大体如下：

def generate_butterworth_mask(rows, cols, n, d, flip=False):""":param rows: rows of spectrum:param cols: columns of spectrum:param n: the filter adjustment factor:param d: the d0:return:"""import numpy as npoutput = np.zeros((rows, cols), np.float32)# center coordinatecr = rows / 2cc = cols / 2for r in range(rows):for c in range(cols):distance = np.sqrt((r - cr) ** 2 + (c - cc) ** 2)frac = (distance / d) ** (2 * n)output[r, c] = 1 / (1 + frac)if flip:output = 1 - outputreturn output

得到的效果大概是这样的：

数字图像学笔记——10. 频域与傅里叶分析方法相关推荐

数字图像学笔记 —— 16. 图像退化与复原（自适应滤波之「最小均方差滤波」）
文章目录图像恢复的一般运算过程什么是「最小均方差滤波」实现步骤实现代码最后的结果图像恢复的一般运算过程我们从前几章的基本理论出发,退化信号恢复成原始信号的步骤,可以概括成两步基本公式.对 ...
数字图像学笔记——13. 图像退化与复原（退化函数的评估方法：观察法、实验法、数学建模法与湍流导致的退化）
在对受到多种原因影响的图像进行复原时,我们经常需要先行评估对图像质量产生影响的退化函数,有时甚至需要尝试建模.通过这些手段,能够最大程度上恢复图像上的噪音,并重建高清的图像细节. 文章目录线性位置不 ...
数字图像学笔记——14. 图像退化与复原（线性退化）
文章目录运动导致的退化(线性退化) 水平运动导致的退化垂直运动导致的退化运动导致的退化(线性退化) 在上一章 <数字图像学笔记--13. 图像退化与复原(退化函数的评估方法:观察法.实验法 ...
数字图像学笔记——7. 噪音生成（泊松噪音生成方法）
文章目录泊松噪音 Knuth算法散列生成算法生成泊松噪音的图像泊松噪音 Knuth算法首先,回顾泊松分布的函数: P(x=k)=e−λλkk!P(x=k) = \frac{e^{- \lam ...
数字图像学笔记——3.彩色转黑白
文章目录一些说明关于示例代码关于依赖环境关于教材灰度图.亮度图(Gray Image) 彩色图转灰度图一般亮度转换(luminosity method) 亮度优先转换(luminosity ...
数字图像学笔记——4. 直方图计算、线性变换、对数变换、Gamma变换
文章目录灰度直方图(Gray Histogram) 直方图的计算方法简单的图像转换方法线性变换 / 图像翻转(Image Nagatives) 对数变换(Log Transformation) ...
Linux学习笔记10
Linux学习笔记10 Linux学习笔记10 正则表达式源码包约定目录 Shell脚本约定目录 Shell脚本的创建与执行 date命令同步时间 Shell脚本预设变量与用户交互数学计算 S ...
台大李宏毅Machine Learning 2017Fall学习笔记 (10)Tips for Deep Learning
台大李宏毅Machine Learning 2017Fall学习笔记 (10)Tips for Deep Learning 注:本博客主要参照 http://blog.csdn.net/xzy_thu ...
Autoware感知瞎学笔记（一）lidar_kf_contour_track
Autoware感知瞎学笔记(一)lidar_kf_contour_track 目录代码分析: 一.雷达目标Kalman滤波器 1. lidar_kf_contour_track.cpp 2. li ...

数字图像学笔记——10. 频域与傅里叶分析方法