1. 这篇笔记要解决的问题写这篇笔记前我先把前面几篇的思路捋一捋我们已经搞定了离散傅里叶变换的原理、快速傅里叶变换的用法、频率谱的解读低频中间、高频外围、亮度代表能量以及频域滤波的基本流程——把图像从空间域变到频率域乘以一个滤波器再变回去。这个系列到了第 5 篇核心问题就一个滤波器本身长什么样也就是说在 H(u, v) 这一步上我们到底有哪些选择。不同滤波器带来的效果差异非常大有的平滑得像磨皮有的锐利得像拉满对比度有的会把图像弄出一圈一圈的水波纹有的却干净利落。这篇笔记想解决的就是把常见的频域滤波器全部过一遍讲清楚每个的原理、数学表达式、参数怎么调、适合用在什么地方以及哪些坑是你写代码时一定会踩的。适合看这篇笔记的人有两类。一类是正在学数字图像处理课程的学生课本上对这部分讲得比较零散公式堆在那但不知道怎么落到代码里另一类是工作中真正要做图像去噪、增强、纹理分析的开发者需要在传统方法里挑一个能用的方案。这篇笔记既照顾理论推导也把完整的 Python 代码贴出来大家可以直接抄作业。我的目标只有一个看完这篇你能对着任意一张图像快速决定该用哪个滤波器、参数设多少、代码怎么写而不是翻着书发愣。2. 频域滤波的框架回顾你还缺哪块拼图2.1 频域滤波的本质乘法就是卷积先花半分钟回顾一下频域滤波的核心公式。空间域的卷积对应频率域的乘法这是卷积定理告诉我们的g(x, y) IFFT[ F(u, v) × H(u, v) ]其中 F(u, v) 是图像经过傅里叶变换后的频域表示H(u, v) 就是我们设计的滤波器传递函数g 是滤波结果。整件事的本质就一句话你想保留哪些频率、滤掉哪些频率就把对应的频点乘以 1 或乘以 0或者中间值。但光知道这个框架是不够的。实际动手时你会发现真正的技术含量完全在 H(u, v) 的设计上。公式写起来都简单参数选不好结果天差地别。这也是为什么这个系列需要单独用一整篇笔记来讲滤波器本身。以我做图像处理的经验为例很多人一开始以为频域滤波就是把高频砍掉就模糊了把低频砍掉就锐化了。这个理解大方向没错但真的上手后要么图像出现诡异的振铃要么整张图变暗要么边缘出现条状伪影。这些现象大多不是流程写错了而是滤波器选型或者参数设置的问题。2.2 一张图说清滤波器分类常见的频域滤波器按功能可以分成四个基本类型低通、高通、带通、带阻。每个类型又有多种实现形式常见的就是理想型、巴特沃斯型、高斯型这三种。低通保留低频抑制高频用于平滑、去噪高通保留高频抑制低频用于锐化、边缘提取带通保留某个频带内的成分用于提取特定尺度纹理带阻抑制某个频带内的成分用于去除周期性噪声而在每种类型内部区别主要在于过渡带的处理方式。理想型在截止频率处直接一刀切巴特沃斯型有一个平滑但可调的过渡带高斯型则是完全平滑没有振铃。这三种之间的细微差别就是我在实际工作中最看重的部分。3. 理想滤波器理论上的标准答案实践中的踩坑之王3.1 理想低通与理想高通的数学定义理想低通滤波器ILPF的传递函数是H(u, v) 1, 如果 D(u, v) D0 H(u, v) 0, 如果 D(u, v) D0其中 D(u, v) 是频率点 (u, v) 到频率中心也就是零频的欧几里得距离D0 是截止半径。理想高通IHPF刚好反过来H(u, v) 0, 如果 D(u, v) D0 H(u, v) 1, 如果 D(u, v) D0这组定义好理解也好看但它是一个典型的数学上完美、工程上翻车的方案。原因在于理想的频率响应意味着时域空间域的脉冲响应是无限长的 sinc 函数而对我们这种有限尺寸的图像截断必然产生振荡这就是振铃效应ringing的根源。我最早学这块的时候第一次写出理想低通滤波的代码兴冲冲把一张人像图拿去跑结果出来的图像边缘一圈一圈的水波纹当时还以为是自己代码写错了。查了半天才发现这就是理想滤波器的天性不是 bug是 feature 里带的问题。3.2 为什么理想滤波器几乎不会在实际中使用用吉布斯现象最容易解释这个问题。理想低通在截止频率处从 1 直接跳到 0相当于在频域乘了一个矩形窗对应在空间域就是与一个无限长的 sinc 核做卷积。sinc 函数的两侧有正负振荡的旁瓣卷积结果在图像边缘和细节附近就会出现过冲和振铃。凡是说明书里写着可获得最锐利的截止特性这类话的滤波器你都要警惕。实际图像处理任务中几乎没有场景能容忍振铃伪影。医学影像里出现这种伪影会被怀疑是病变工业视觉里这直接导致误判。所以我的建议很明确理想滤波器只作为理论教学的标尺不做实际选型。如果你非要用理想型做实验唯一需要注意的实操点就是 D0 的选择。这个值一定要配合频谱图的能量分布来看先画出 log 频谱看看能量集中区域大概多大再定 D0不要拍脑袋给一个数。即便如此结果大概率还是会有振铃。4. 巴特沃斯滤波器平滑过渡的实用主义者4.1 巴特沃斯低通的公式与阶数含义巴特沃斯低通滤波器BLPF的传递函数为H(u, v) 1 / (1 [D(u, v) / D0]^(2n))这里 n 是滤波器的阶数正整数D0 是截止频率。注意在 D D0 这个点上无论 n 取多少分母都是 112所以 H 0.5也就是说巴特沃斯滤波器的截止频率定义在增益衰减 3dB 的位置这和理想滤波器的定义方式完全不同看曲线的同学千万别搞混。阶数 n 控制了过渡带的陡峭程度。n 越大过渡带越窄频率响应越接近理想滤波器n 越小过渡带越宽振铃越弱但过滤效果也没那么干净。n1 时基本没有振铃n2 时轻微振铃但可接受n3 时振铃逐渐明显。这是一个典型的敏感度-平滑度博弈。4.2 从低通到高通只需要把比值倒过来巴特沃斯高通滤波器BHPF的公式是H(u, v) 1 / (1 [D0 / D(u, v)]^(2n))其实就是把 D 和 D0 的比值倒过来。这个变换非常优雅完全不需要重新设计滤波器。当 D 远大于 D0 时H 趋近于 1高频通过当 D 趋近于 0 时H 趋近于 0低频被抑制。我在做高通滤波时个人的默认选择就是巴特沃斯 n2兼顾平滑和保留细节。相比于后面要讲的高斯高通巴特沃斯高通在同样 D0 下对中间频段的保留更慷慨一点更适合需要同时保留部分低频信息的场景。当然这只是经验值具体还要看你的图。4.3 巴特沃斯滤波器的实战参数建议直接说几个我实测过的参数范围D0 在频率谱归一化半径的 0.1~0.3 之间适合做平滑去噪n1 适合对图像质量要求高、不能容忍任何伪影的场景n2 是通用默认值多数情况下效果和效率均衡n4 就没什么必要了振铃会逐步回归你还不如用理想滤波器省事有一点必须提醒随着 n 增大巴特沃斯滤波器越来越接近理想滤波器但永远达不到理想那种陡峭所以振铃也只是逼近通常不会像理想型那样夸张。这个特性让它成为很多教科书里折中方案的首选推荐。5. 高斯滤波器实测中最稳的选择5.1 高斯低通与高通为什么没有振铃高斯低通滤波器GLPF的传递函数是H(u, v) exp(-D²(u, v) / (2 * D0²))高斯高通滤波器GHPF则是H(u, v) 1 - exp(-D²(u, v) / (2 * D0²))高斯函数有一个绝妙的自相似属性它的傅里叶变换仍然是高斯函数。这意味着频域的高斯滤波对应空间域的卷积核也是高斯形状而高斯核是一个非负、单峰的平滑核不存在负旁瓣。旁瓣是振铃的根源没有旁瓣自然就没有振铃。我私心最喜欢高斯滤波器原因就是它干净。不管 D0 设得多极端输出图像不会出现一圈圈的伪影最多就是过度模糊或者过度锐化。对于实际工程来说宁可效果差一点也不能让输出结果出现奇怪的假信息在这个原则下高斯型几乎是必选。5.2 高斯带通与高斯带阻的构造思路高斯带通和高斯带阻不像低通高通那样有标准的单公式更常见的做法是通过低通和高通组合构造。带通就是两个低通的差或者两个高通的组合带阻则是反过来用 1 减去带通。以带阻为例如果想滤除中心频率为 D0、带宽为 W 的频带一个常见的高斯带阻形式是H(u, v) 1 - exp(-[(D²(u, v) - D0²) / (D(u, v) * W)]²)这里 D0 是待抑制频带的中心W 是带宽参数。这个公式的含义是当 D 接近 D0 时指数项接近 exp(-0)1H 接近 0该频带被抑制偏离 D0 越远H 越接近 1。带阻滤波器在去周期性噪声上特别好用。比如扫描图片上有一层周期性的网纹或者照片上有灯光造成的摩尔纹它们反映在频谱上就是几个对称的亮斑点用带阻滤波器或者说陷波滤波器把这些斑点盖住反变换回来网纹就消失了。这是频域滤波相对空间域滤波最有优势的场景之一空间域做这种事极其痛苦频域只要在频谱图上圈几个点就行。5.3 同态滤波高斯高通最有价值的应用同态滤波是我实际工作中用得最多的频域增强手段它本质上是把图像先取对数再做高通滤波再做指数还原。数学依据是照射-反射模型f(x, y) i(x, y) * r(x, y)照射分量 i 变化缓慢集中在低频反射分量 r 变化剧烈集中在高频。取对数后乘法变成加法于是可以用高通滤波把二者分开压缩照射分量的动态范围增强反射分量的细节。一个标准的同态滤波传递函数长这样H(u, v) (γH - γL) * (1 - exp(-c * (D²(u, v) / D0²))) γL其中 γL 控制低频增益通常 1γH 控制高频增益通常 1c 控制过渡带的陡峭程度D0 是截止频率。我常用的初值是 γL0.5、γH1.5、c2然后根据结果微调。同态滤波对背光照片、雾天图像、光照不均的图像有奇效在我看来是每个做图像增强的人必须掌握的工具。下表是三种滤波器在几个关键维度上的对比方便大家选型滤波器类型过渡带特性振铃风险可调参数主要用途理想型极陡突变严重D0理论教学、效果对比巴特沃斯型可调n控制中低D0、n通用平滑/锐化高斯型平滑自然无D0去噪、增强、同态滤波6. 频域滤波与空间域滤波的对应关系6.1 卷积定理视角下的等价性前面提到卷积定理这里需要正面回答一个很多人纠结的问题既然频域滤波本质是空间卷积那为什么不直接在空间域做卷积算了答案是要看滤波器尺寸。如果你设计的等效空间核很小比如 3x3 或 5x5那空间域卷积速度快频域反而要来回 FFT浪费时间。但如果等效核很大比如几十上百像素半径空间域卷积的时间复杂度是 O(N²×M²)而频域只需要两次 FFT 加一次逐点相乘O(N²logN)优势就非常明显了。以高斯滤波为例一个标准差为 20 像素的高斯核对应空间域核尺寸可能要 100x100 以上这时候空间域卷积慢得离谱频域实现几乎是瞬时的。很多项目里我都建议用频域实现大核高斯滤波这是实打实的性能优化。6.2 空间核的尺寸与边界问题频域滤波隐含了一个周期性假设DFT 认为图像在水平和垂直方向上是无限周期重复的。这带来两个后果。第一如果图像左右边界本来就不连续频域会出现十字形的亮线伪影第二的结果在图像边界处可能产生折叠效应。缓解办法主要有两个一是先对图像做边界填充用镜像填充或复制填充而不是零填充把边界的不连续性尽量消除二是用窗函数对原始图像做预处理让图像边缘渐变为零。这两种方式都能减少频谱上的十字亮线我在实操中优先推镜像填充因为简单有效且不改变图像整体亮度分布。6.3 什么时候用频域什么时候用空间域给你一个判断标准如果处理的是大尺寸核、全局性操作去噪、去周期纹理、同态增强用频域如果是局部小窗口操作边缘检测、局部对比度、锐化细节用空间域。空间域滤波的优点是局部性好、计算开销低、适合实时处理频域滤波的优点是全局性强、效果直观直接在频谱上画形状、适合分析和定制。两者各有各的舞台不冲突甚至经常配合使用。7. Python 实操从零实现四种滤波器7.1 基础工具函数中心化与距离矩阵进入代码环节。我们的基础环境是 Python NumPy OpenCV可选 Matplotlib。首先准备两个基础工具一个是频谱中心化把零频移到画面中心另一个是构建距离矩阵。import numpy as np import matplotlib.pyplot as plt def fft2_shift(img): # 图像转 float32执行 FFT 并中心化 f np.fft.fft2(img.astype(np.float32)) f_shift np.fft.fftshift(f) return f, f_shift def get_distance_matrix(rows, cols): # 返回每个频点到中心的欧氏距离矩阵 cx, cy rows // 2, cols // 2 u np.arange(rows) - cx v np.arange(cols) - cy U, V np.meshgrid(u, v) return np.sqrt(U**2 V**2) def ifft2_shift(f_shift): f np.fft.ifftshift(f_shift) img np.fft.ifft2(f) img np.real(img) return np.clip(img, 0, 255)这里有个新手特别容易忽略的细节fftshift和ifftshift不是同一个函数的互逆吗是但有个微妙的区别。fftshift把零频从左上角搬到中心ifftshift把中心搬回左上角。对于偶数尺寸的图像两者效果一致但如果是奇数尺寸两者结果不完全一样所以安全起见正变换用 fftshift逆变换用 ifftshift不要混用。这个细节我曾经在某个深夜调了半小时最后发现就是这里错了。7.2 四种滤波器的完整实现下面直接给出四种滤波器的实现代码。我刻意把输入输出设计成统一的接口输入原始图像和参数直接返回滤波后的图。def ideal_lowpass_filter(img, D0): rows, cols img.shape D get_distance_matrix(rows, cols) H (D D0).astype(np.float32) _, f_shift fft2_shift(img) return ifft2_shift(f_shift * H) def ideal_highpass_filter(img, D0): rows, cols img.shape D get_distance_matrix(rows, cols) H (D D0).astype(np.float32) _, f_shift fft2_shift(img) return ifft2_shift(f_shift * H) def butterworth_lowpass_filter(img, D0, n2): rows, cols img.shape D get_distance_matrix(rows, cols) H 1.0 / (1.0 np.power(D / D0, 2 * n)) _, f_shift fft2_shift(img) return ifft2_shift(f_shift * H) def butterworth_highpass_filter(img, D0, n2): rows, cols img.shape D get_distance_matrix(rows, cols) D np.maximum(D, 1e-6) # 避免除零 H 1.0 / (1.0 np.power(D0 / D, 2 * n)) _, f_shift fft2_shift(img) return ifft2_shift(f_shift * H) def gaussian_lowpass_filter(img, D0): rows, cols img.shape D get_distance_matrix(rows, cols) H np.exp(-(D**2) / (2 * D0**2)) _, f_shift fft2_shift(img) return ifft2_shift(f_shift * H) def gaussian_highpass_filter(img, D0): rows, cols img.shape D get_distance_matrix(rows, cols) H 1.0 - np.exp(-(D**2) / (2 * D0**2)) _, f_shift fft2_shift(img) return ifft2_shift(f_shift * H)注意巴特沃斯高通里我加了np.maximum(D, 1e-6)因为距离矩阵在中心点是 0直接代入 D0/D 会得到无穷大导致结果全部变成 NaN。这个坑不多说踩过的人都懂。7.3 实操对比同一张图的滤波效果分析下面我们加载一张灰度图分别用六种滤波器跑一遍并对比结果。from skimage import data, img_as_float img img_as_float(data.camera()) fig, axes plt.subplots(2, 4, figsize(12, 8)) axes[0][0].imshow(img, cmapgray) axes[0][0].set_title(Original) filters [ (Ideal LP D030, ideal_lowpass_filter(img, 30)), (Butterworth LP n2 D030, butterworth_lowpass_filter(img, 30, 2)), (Gaussian LP D030, gaussian_lowpass_filter(img, 30)), (Ideal HP D030, ideal_highpass_filter(img, 30)), (Butterworth HP n2 D030, butterworth_highpass_filter(img, 30, 2)), (Gaussian HP D030, gaussian_highpass_filter(img, 30)), ] for i, (title, result) in enumerate(filters): ax axes[(i 1) // 4][(i 1) % 4] ax.imshow(result, cmapgray) ax.set_title(title) plt.tight_layout() plt.show()实测效果通常是这样高斯低通最平滑柔和巴特沃斯低通边缘保留了稍多细节但略有振铃理想低通振铃最明显高通的结果在图像边缘和纹理处响应强烈理想高通会呈现明显的黑白振铃圈高斯高通最自然。在实际项目中几乎每次我都会优先看高斯型的结果再根据需求决定要不要换成巴特沃斯增加过渡带的锐利度。7.4 D0 与阶数 n 的调参经验关于参数我直接给一组经过多次试验的参考区间D0 的选择先看频谱图的能量分布。把 log(1|F|) 的频谱用 imshow 显示出来观察亮区半径。如果亮区半径大约是图像短边的 1/10那 D0 取这个值做低通就是保留主要能量如果要做高通D0 取更小一点比如 10~20不然会把整张图变成黑底白边。n 的选择n1 最稳n2 通用n3 以上慎用。做对比实验时建议把 n1、2、3 的结果放一起看振铃的出现是非常直观的信号一有振铃就说明你调过头了。要特别强调D0 的单位是像素距离也就是频域坐标里的半径不是你图像里的像素尺寸。很多人第一次上手会把 D0 设成和图像边长一样大结果滤波后图像几乎没变化因为截止频率太高了所有频率都被放行了。我自己的经验是 D0 先从小到大多扫几组10、30、50、80肉眼观察输出质量的拐点这样最靠谱。8. 常见问题与避坑指南8.1 中心化操作fshift 与 ifftshift 的完全指南新手最容易出的问题就是中心化。FFT 的结果默认零频在左上角直接画频谱图会看到四个角各有一块亮斑而不是中心一块。这时候要先用fftshift把零频挪到中间可视化才正常。滤波时滤波器 H 通常也是中心对称的所以你必须在移位后的频谱上乘 H然后ifftshift再逆变换。流程顺序错了结果完全不对。给一个判断方法如果逆变换后的图像出现周期性错位也就是图像内容看起来被切成了四块并发生了位移那基本就是fftshift和ifftshift用混了。还有一个常见表现是低频滤波器输出的图像整体发灰、亮度明显下降这大概率是零频DC 分量没有被正确处理——要么你乘的滤波器在中心点的值是 0直接把图像的均值干掉了。8.2 数据类型与显示范围FFT 的输入最好是 float32不要直接给 uint8。uint8 的像素范围是 0~255FFT 之后频域的幅度可能非常大中间值可能超过 1e8直接存成 uint8 会截断成一片白或者一片黑。所以代码里所有频域操作都应该在 float32 或 float64 下进行最后输出时再 clip 到 0~255。频谱可视化的经典操作是对幅度取对数np.log(1 np.abs(f_shift))。直接看原始幅度谱除了中心那个大亮斑其他细节几乎都是黑的因为动态范围太大了。取对数可以把范围压缩让高频细节显现出来。这也是为什么教材里的频谱图看起来总是灰蒙蒙的那不是原图是 log 变换后的图。还有一个容易忽略的问题逆变换后取实部np.real而不是幅值np.abs。如果取幅值会把负的像素值变成正的导致图像对比度失真。我见过有人用np.abs反变换后图像看起来发白就是这个原因。8.3 振铃与边界伪影怎么识别振铃的识别特征很清晰图像边缘附近出现一圈一圈交替的明暗条纹像水波纹扩散。如果你看到这种图案先检查是不是用了理想滤波器如果确实用了理想型换成巴特沃斯或高斯型即可解决。如果用的是巴特沃斯但 n 太大同样要降低 n。边界伪影则是另一回事图像边缘出现横纵交叉的亮线尤其在做高通滤波后会放大。这是因为 DFT 隐含周期延拓图像左右和上下边缘在周期延拓后不连续产生了沿 x 轴和 y 轴的十字亮线。缓解方法是先对图像做边界填充我用过的最有效的两种镜像填充np.pad(img, pad_width, modereflect)有效消除边缘跳变不影响频谱主体结构窗函数在滤波前把图像与汉宁窗相乘使边缘渐变为 0伪影消除彻底但会略微改变图像局部统计信息具体选哪个取决于图像内容。风景图像的整体灰度变化缓慢镜像填充更好显微图像或工业图样需要保留细小特征窗函数更稳。两种都试一遍用肉眼判断最简单直接。8.4 非线性噪声求你别用频域滤波去处理这里我要泼一盆冷水频域滤波只擅长处理加性平稳噪声比如高斯噪声、周期噪声。对于椒盐噪声这种脉冲型噪声频域滤波的效果非常差。原因好理解椒盐噪声的频域特性是宽频谱它和图像的频谱混叠在一起你不可能在不伤及图像细节的情况下把它去掉。所以如果你需要对含椒盐噪声的图像去噪请去用中值滤波或者现在更流行的非局部均值、双边滤波。频域方法在这个场景下是杀鸡用牛刀而且杀不干净。判断噪声类型的方法很简单先看频谱图。如果频谱上出现规律排布的亮星点说明存在周期性噪声用带阻/陷波滤波如果频谱整体均匀分布、没有明显星点那就是随机噪声用低通或小波去噪。观察频谱再做决定这个习惯能让你的项目少走很多弯路。9. 最后的个人心得写到这回想起我做图像处理项目从入门到现在的过程最有价值的一条经验是不要一上来就套滤波器先花几分钟看一下频谱长什么样。频谱是图像在频率域的体检报告该滤哪个频带、不该动哪个频带它都会告诉你。很多调参难题在看了频谱之后都会变得特别简单因为答案就写在那里。另外关于这两年被反复讨论的频域滤波是不是过时了——我的看法是深度学习虽然在图像增强去噪上打败了很多传统方法但频域的思想从来没有退场甚至连深度学习模型里都在用频域做注意力机制。传统频域滤波依然是理解图像全局结构的捷径也是很多工程场景里成本最低、可解释性最强的一招。如果你要继续往下深入我建议接下来研究两个方向。一是陷波滤波器在去周期性噪声上的精细应用多频带噪声怎么逐一处理二是频域滤波与金字塔结合的多尺度增强它在 HDR 和医学图像里非常好用。这个系列我也会继续往后写下一篇应该会专门把同态滤波和陷波滤波展开这两个应用场景最多、也最容易出效果。