早期视觉：滤波

从感知到处理

前几讲我们关心三件事：相机如何把三维世界投影成二维图像、光与颜色如何被传感器记录、以及如何用变换在像素之间挪动信息。从本讲开始，我们正式进入早期视觉（Early Vision）——视觉系统在「看见」之后做的第一批处理。

这一阶段的目标不是「认出物体」，而是从噪声、纹理、亮度差异中提炼出更稳定、更有用的特征：让暗部细节看得清、让噪声消失、让边缘更锐利、让相同物体在不同尺度下被同样地识别。完成这些任务的核心工具就是滤波（Filtering）。

广义的滤波涵盖三种互补的视角：

flowchart LR
    A["空间域滤波<br/>对像素邻域<br/>做数学运算"]
    B["频率域滤波<br/>修改图像的<br/>频率成分"]
    C["模板与金字塔<br/>把滤波器视作<br/>模板进行匹配"]

    A -.同一运算.- B
    B -.多尺度推广.- C

    classDef view fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    class A,B,C view

本讲沿着这三条路径走完一遍：先从最直观的空间域操作（点操作、形态学、卷积、平滑、锐化）出发；再用傅里叶变换打开频率域视角，理解卷积定理与小波分解；最后引入采样定理和图像金字塔，建立处理多尺度的统一框架。

滤波与变形

要给「图像变换」一个清晰定位，先把它和邻近的概念区分开。对图像 $f(x,y)$ 做某种数学操作得到新图像 $g(x,y)$，按改变的对象不同分为两类：

滤波改变像素的值，输出图像在每个位置上重新计算了像素值，但空间布局不变。用函数语言说，它改变的是值域（Range）：

$$g(\bm{x}) = h\{f(\bm{x})\}$$

变形（Warping）改变像素的位置，像素值本身不变但被移动到新坐标。它改变的是定义域（Domain）：

$$g(\bm{x}) = f(h\{\bm{x}\})$$

本讲专注于滤波。根据计算时是否考虑邻域，又分为点操作（仅依赖当前像素）和邻域操作（依赖周围一片区域）。点操作是接下来要看的最简单情况。

点操作

点操作（Point Processing）的定义：输出像素值仅取决于输入图像中同一位置的像素值。设原始像素值为 $x$（取值 0–255）：

操作	公式	效果
变暗	$x - c$	整体变暗
变亮	$x + c$	整体变亮
降低对比度	$x / 2$	灰度差异减小
提升对比度	$2x$	灰度差异增大
反色	$255 - x$	「底片」效果
非线性降低对比度	$255 \cdot (x/255)^{1/3}$	暗部拉伸、亮部压缩
非线性提升对比度	$255 \cdot (x/255)^{2}$	暗部压缩、亮部拉伸

非线性操作使用幂函数，亦称 Gamma 变换：指数小于 1 时图像整体变亮、对比度降低；指数大于 1 时图像整体变暗、对比度提高。这与显示器的伽马校正背后是同一个数学。

动态范围压缩

某些场景的亮度范围远超显示设备能呈现的范围（典型的如夜景中的灯光，或宇航员太空舱照片）。直接显示会让最亮处过曝、暗处全黑。对数变换能把宽动态范围压回有限区间：

$$g(x, y) = c \log\bigl(1 + f(x, y)\bigr)$$

对数曲线在小数值上斜率大、大数值上斜率小——暗部细节被拉伸放大，亮部细节被压缩在一起，恰好对应人眼对亮度对数式的感知。

直方图均衡化

点操作的更精细版本来自一个观察：图像看起来「灰扑扑」的，往往是因为像素值挤在某个狭窄的区间里。如果能把这个区间「拉开」，图像就显得有反差。灰度直方图正是定量描述像素值分布的工具。

直方图

图像的灰度直方图（Gray-level Histogram）是 1-D 离散函数：

$$h(f_k) = n_k, \quad k = 0, 1, \ldots, L-1$$

其中 $f_k$ 是第 $k$ 级灰度值，$n_k$ 是图像中该灰度的像素数，$L$ 是灰度级数（常为 256）。除以总像素数 $n$ 可得到归一化直方图（即概率密度）：

$$p(f_k) = \frac{n_k}{n}$$

直方图与图像的视觉效果有直接对应：偏暗的图像直方图集中在低灰度一侧，偏亮的集中在高灰度一侧，对比度低的集中在中间窄带，对比度高的横跨整个灰度轴。改善视觉效果可以转化为改造直方图——这就是均衡化的基本想法。

均衡化原理

直方图均衡化（Histogram Equalization）的目标是把原始直方图映射为接近均匀分布的形式，使所有灰度级被「平均」地使用，从而最大化动态范围和整体对比度。

实现方式是寻找一个单调递增的灰度映射函数 $g = T(f)$，让映射后的直方图近似均匀。理论上这个 $T$ 就是累积分布函数（Cumulative Distribution Function, CDF）：

$$T(f) = (L - 1) \cdot \sum_{k=0}^{f} p(f_k) = (L - 1) \cdot \text{CDF}(f)$$

直觉上，CDF 把「集中区间」拉伸为「均匀区间」。原直方图越密集的灰度段，CDF 在那里斜率越大，映射后被分配的输出灰度区间越宽——密集像素被「拉开」。

直方图均衡化的步骤

1. 统计原图像各灰度的像素数 $n_k$，得到归一化直方图 $p(f_k) = n_k / n$
2. 计算累积直方图 $\text{CDF}(f) = \sum_{k=0}^{f} p(f_k)$
3. 定义灰度映射 $g = \lfloor (L-1) \cdot \text{CDF}(f) \rfloor$
4. 对每个像素应用映射，得到均衡化后的图像

映射 $T$ 满足两个关键性质：(1) 单调递增，保证亮度顺序不变（亮的还是亮，暗的还是暗）；(2) 输出范围保持在 $[0, L-1]$，不溢出。

一个具体例子

考虑一张 $64 \times 64$、灰度级 $L=8$ 的图像，直方图统计如下：

$f$	0	1	2	3	4	5	6	7
$n_k$	790	1023	850	656	329	245	122	81
$p$	0.19	0.25	0.21	0.16	0.08	0.06	0.03	0.02
CDF	0.19	0.44	0.65	0.81	0.89	0.95	0.98	1.00
$T$	1	3	5	6	6	7	7	7

映射后灰度集中的低端（0,1,2）被拉到中间和高端（1,3,5），而稀疏的高端（5,6,7）被合并到 7——动态范围从 $[0,3]$ 占主导拓宽到了整个 $[1,7]$。

彩色图像的均衡化

直接对 RGB 三通道分别均衡化会破坏色调——红绿蓝的相对比例被独立改变，颜色严重失真。正确做法是：

1. 把图像转换到 HSV 或 YCbCr 等亮度-色度分离的颜色空间
2. 仅对亮度通道（V 或 Y）做均衡化
3. 保留色调和饱和度通道
4. 转换回 RGB

这样既增强对比度又保持原有色彩。

自适应均衡化与 CLAHE

全局均衡化对整张图共用一个映射，对光照不均的图像效果有限——一张照片局部曝光过度、另一处曝光不足时，单一映射难以兼顾两端。

自适应直方图均衡化（Adaptive Histogram Equalization, AHE）把图像分成若干小块（如 $8 \times 8$），每个块独立计算和应用自己的均衡化映射，再用插值消除块边界。这样每个区域都能根据局部内容调整对比度。

AHE 有两个明显的副作用：

• 噪声放大：在原本平坦的区域（直方图集中），均衡化会大幅拉伸有限的灰度差异，把不可见的微小噪声放大成显眼的颗粒
• 块边界伪影：相邻块的映射函数不同，过渡处可见接缝

对比度受限自适应均衡化（Contrast Limited Adaptive Histogram Equalization, CLAHE）是 AHE 的改进：在计算每块直方图时，把高于某个裁剪阈值（clipLimit）的部分截掉，再把截掉的「质量」均匀地分散回所有灰度。

flowchart LR
    A["原直方图<br/>峰值很高"] --> B["设置裁剪阈值<br/>clipLimit"]
    B --> C["超出阈值的部分<br/>被裁剪"]
    C --> D["把裁剪量均匀<br/>分回所有 bin"]
    D --> E["基于新直方图<br/>做均衡化"]

    classDef src fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef proc fill:#fff3e0,stroke:#ef6c00,stroke-width:2px
    classDef out fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px

    class A src
    class B,C,D proc
    class E out

裁剪限制了 CDF 的最大斜率，从而限制了局部对比度的最大增益——噪声放大被有效抑制。CLAHE 是医学影像、卫星遥感、内窥镜等低对比度成像场景的标准预处理之一。

形态学处理

直方图操作的对象是「灰度的分布」，形态学（Morphology）则关注「形状的结构」。它的基本想法是：用一个小的结构元素（Structuring Element）$B$ 在图像 $A$ 上滑动，根据 $B$ 与 $A$ 局部的相互关系生成输出。最常用于二值图像，也可推广到灰度。

腐蚀与膨胀

膨胀（Dilation）求局部最大值：对每个位置，把结构元素覆盖区域内的最大像素值赋给中心。

腐蚀（Erosion）求局部最小值：把结构元素覆盖区域内的最小像素值赋给中心。

直观地说，膨胀让亮区扩张、暗区收缩；腐蚀让暗区扩张、亮区收缩。在二值图像中：

• 膨胀：只要结构元素与目标有任何重叠（"hit"），中心置 1
• 腐蚀：只有结构元素完全包含在目标内（"fit"），中心置 1

flowchart TD
    subgraph 膨胀
        D1["B 与 A 有交集"] --> D2["输出 = 1"]
        D3["B 与 A 无交集"] --> D4["输出 = 0"]
    end
    subgraph 腐蚀
        E1["B 完全包含于 A"] --> E2["输出 = 1"]
        E3["B 部分超出 A"] --> E4["输出 = 0"]
    end

    classDef hit fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    classDef miss fill:#ffebee,stroke:#c62828,stroke-width:2px
    class D2,E2 hit
    class D4,E4 miss

形式化定义：

$$A \oplus B = \{ z \mid (\hat{B})_z \cap A \ne \emptyset \}, \qquad A \ominus B = \{ z \mid B_z \subseteq A \}$$

其中 $B_z$ 表示把 $B$ 平移到位置 $z$，$\hat{B}$ 表示 $B$ 关于原点的反射。

与卷积的区别

形态学的滑动操作和卷积形似，但结构元素不必是矩形——可以是线、十字、圆、菱形等任意形状，且不参与「乘法 + 求和」，只参与「集合包含」或「极值」运算。

开运算与闭运算

腐蚀和膨胀的复合运算给出更精细的形状操作。

开运算（Opening）：先腐蚀后膨胀。

$$A \circ B = (A \ominus B) \oplus B$$

效果：消除小于结构元素的细小亮点（如噪声、孤立像素）、断开窄连接、平滑较大物体的轮廓——但不显著改变其面积。

闭运算（Closing）：先膨胀后腐蚀。

$$A \bullet B = (A \oplus B) \ominus B$$

效果：填充小于结构元素的细小空洞、连接断裂的邻近像素块、平滑边界——同样不显著改变物体形状。

直观记忆：开运算像「打开」缝隙、清掉小杂物；闭运算像「闭合」孔洞、连接断点。两者都不改变物体大致形状，是形态学滤波的主力。

梯度、顶帽与黑帽

由腐蚀和膨胀可派生出一组用于结构提取的运算。

形态学梯度（Morphological Gradient）：膨胀与腐蚀之差。

$$\text{Grad}(A) = (A \oplus B) - (A \ominus B)$$

膨胀让亮区扩张、腐蚀让亮区收缩，相减后只剩下变化最剧烈的位置——也就是物体的边缘。形态学梯度是一种朴素但有效的边缘检测算子。

顶帽运算（Top-hat）：原图减开运算。

$$\text{TopHat}(A) = A - (A \circ B)$$

开运算去掉了小于结构元素的亮特征，相减后剩下的就是这些被去掉的亮特征——比邻近背景更亮的小区域被凸显。背景不均匀时，顶帽变换可以帮助恢复二值化时丢失的边缘细节。

黑帽运算（Black-hat）：闭运算减原图。

$$\text{BlackHat}(A) = (A \bullet B) - A$$

类似地凸显比邻近背景更暗的小区域，可用于在亮背景上分离暗物体。

运算	公式	用途
膨胀	$A \oplus B$	扩张亮区、连接邻近物体
腐蚀	$A \ominus B$	收缩亮区、消除毛刺
开运算	$(A \ominus B) \oplus B$	去除小亮点、平滑边界
闭运算	$(A \oplus B) \ominus B$	填充小孔洞、连接断点
梯度	$(A \oplus B) - (A \ominus B)$	提取边缘
顶帽	$A - (A \circ B)$	凸显暗背景上的亮特征
黑帽	$(A \bullet B) - A$	凸显亮背景上的暗特征

形态学操作常作为预处理或后处理出现：滤除二值图像中的噪声、提取物体轮廓、矫正光照不均、增强字符识别中的笔画等。

图像噪声

接下来转向更常被称为「滤波」的邻域操作。这之前先建立动机：图像中为什么会有噪声，又为什么需要滤掉它？

即使对静止场景拍多张照片，得到的图像也不会完全相同。噪声（Noise）来自传感器的热噪声、光子的统计波动、电路干扰、量化误差等多种物理因素，是成像过程中不可避免的副产品。

常见的噪声类型有三种：

噪声类型	数学模型	视觉表现
椒盐噪声（Salt and Pepper）	随机出现的纯黑（0）和纯白（255）像素	图像上散布着黑白斑点
脉冲噪声（Impulse）	随机出现的纯白像素	类似椒盐噪声，但只有白点
高斯噪声（Gaussian）	$g(x, y) = f(x, y) + n$，$n \sim \mathcal{N}(0, \sigma^2)$	整体颗粒感，强度均匀

去噪的一个根本观察是：真实图像中相邻像素往往相似（来自同一表面），而噪声在像素间是独立的。如果有同一场景的多张图，可以取平均消除随机噪声；只有一张图时，则可用空间邻域作为「多次观测」——这就是邻域平均的合理性。

邻域滤波

滑动平均

最简单的邻域操作是滑动平均（Moving Average）：把每个像素替换为其邻域内所有像素的均值。一维情况下：

$$g[i] = \frac{1}{2k+1} \sum_{u=-k}^{k} f[i+u]$$

效果是高频细微波动被压平，低频整体趋势保留——这正是「平滑」的数学定义。允许给邻居赋不同权重时，权重整体构成滤波器核（Filter Kernel）。

互相关

把滑动平均推广到二维并允许任意权重，得到互相关（Cross-correlation）运算：

$$g[i, j] = \sum_{u=-k}^{k} \sum_{v=-k}^{k} h[u, v] \, f[i+u, j+v]$$

记为 $g = h \otimes f$。其中 $h[u, v]$ 称为核（kernel）或掩模（mask）——它规定每个相对位置的邻居以什么权重参与加权求和。滤波的本质就是用邻域像素的线性组合替代中心像素。

当 $h$ 全部元素相等时，得到均值滤波器（Box Filter）。例如 $3 \times 3$ 均值滤波器：

$$h = \frac{1}{9} \begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{bmatrix}$$

权重之和为 1 是关键——它保证常数区域在滤波后保持不变（输出仍为该常数）。

卷积

在信号处理中，更基本的运算不是互相关而是卷积（Convolution）。一维连续卷积定义为：

$$g(x) = (f * k)(x) = \int_{-\infty}^{+\infty} f(\tau) \, k(x - \tau) \, \d \tau$$

可以把它分解为五步：

1. 变换：把输入函数自变量从 $x$ 改为 $\tau$
2. 翻转：把核 $k(\tau)$ 翻转为 $k(-\tau)$
3. 平移：得到 $k(x - \tau)$，按 $x$ 滑动
4. 相乘：在每个 $\tau$ 处计算 $f(\tau) \cdot k(x - \tau)$
5. 积分：相乘后曲线下的面积就是 $x$ 位置的卷积值

离散版本对应换积分为求和：

$$g[i] = (f * k)[i] = \sum_{l} f[l] \, k[i - l]$$

二维离散卷积也类似：

$$(f * h)[m, n] = \sum_{k, l} f[k, l] \, h[m - k, n - l]$$

直观上分三步：(1) 翻转核（上下、左右各翻一次），(2) 把翻转后的核滑动到图像每个位置，(3) 每个位置计算逐元素乘积之和。

卷积与互相关

唯一区别在于核是否翻转：

$$\text{互相关：} (f \otimes h)[m, n] = \sum_{k, l} f[m+k, n+l] \, h[k, l]$$

$$\text{卷积：} (f * h)[m, n] = \sum_{k, l} f[m-k, n-l] \, h[k, l]$$

对于左右上下都对称的核（高斯核、均值核），两者结果完全一样。深度学习中常说的「卷积」其实多数是互相关——网络要学的是核本身，翻不翻无所谓，省掉翻转这一步也无影响。

卷积的性质

翻转看似多余，但正是它赋予了卷积良好的代数性质：

性质	公式
交换律	$a * b = b * a$
结合律	$a * (b * c) = (a * b) * c$
分配律	$a * (b + c) = (a * b) + (a * c)$
数乘可提	$(ka) * b = a * (kb) = k(a * b)$
单位元	单位脉冲 $e = [\ldots, 0, 0, 1, 0, 0, \ldots]$ 满足 $a * e = a$

此外卷积还具有两个系统性质：

平移不变性（Shift Invariance）：滤波器在图像任何位置表现相同，$\text{filter}(\text{shift}(f)) = \text{shift}(\text{filter}(f))$。

线性（Linearity）：$\text{filter}(f_1 + f_2) = \text{filter}(f_1) + \text{filter}(f_2)$。

一个深刻的结论

任何线性且平移不变的算子（Linear Shift-Invariant Operator, LSI）都可以表示为卷积。这意味着卷积是所有线性空间滤波器的统一形式。

实用价值：结合律允许把一连串滤波 $((f * b_1) * b_2) * b_3$ 预编译为单核 $f * (b_1 * b_2 * b_3)$；分配律保证锐化等组合滤波可拆开计算；可分离性（下面讨论）则把二维卷积降为两次一维卷积。

一维与二维卷积的尺寸

设输入信号长度 $m$、核长度 $n$，无填充时输出长度为 $m - n + 1$——反复卷积会让信号越来越短。解决方案是零填充（Zero-padding）：在两侧各补 $p$ 个零。配合步长（Stride）$s$，输出长度为：

$$\left\lfloor \frac{m - n + 2p}{s} \right\rfloor + 1$$

二维情况下，对 $H \times W$ 输入、$k \times k$ 核、填充 $p$、步长 $s$：

$$H_{\text{out}} = \left\lfloor \frac{H - k + 2p}{s} \right\rfloor + 1, \qquad W_{\text{out}} = \left\lfloor \frac{W - k + 2p}{s} \right\rfloor + 1$$

这是后续 CNN 计算输出尺寸的同一公式。

边缘效应

当滤波器窗口移到图像边缘，部分窗口落在图像外。处理方式有几种：

方法	策略	效果
忽略边缘	不计算边缘像素	输出尺寸 < 输入
零填充（Zero）	外部视为 0	边缘出现暗边和梯度
边缘填充（Edge）	复制最近的边缘像素	自然但可能引入伪影
反射填充（Reflect）	沿边缘对称翻转	通常视觉效果最佳
环绕填充（Wrap）	视为周期信号	适合周期性图像或频域处理

输出尺寸有三种约定：

• full：输出包含所有部分重叠区域，尺寸为 $(H + k - 1) \times (W + k - 1)$
• same：输出与输入同尺寸（最常用，需配合填充）
• valid：仅保留核完全在图像内的位置，尺寸为 $(H - k + 1) \times (W - k + 1)$

用核做几件具体的事

不同核做不同的事，从直觉上认识几个例子：

核	作用
$\begin{bmatrix} 0 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 0 \end{bmatrix}$	单位脉冲——原样输出
$\begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 0 \end{bmatrix}$	把图像向左平移 1 像素
$\frac{1}{9}\begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{bmatrix}$	均值滤波——平滑模糊
$\begin{bmatrix} 0 & 0 & 0 \\ 0 & 2 & 0 \\ 0 & 0 & 0 \end{bmatrix} - \frac{1}{9}\begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{bmatrix}$	锐化滤波——强调与局部均值之差

最后一行的锐化滤波利用了卷积的分配律：原图（用单位脉冲乘以 2 得到的「双倍原图」）减去其平滑版本，等于原图加上「细节」。后面会专门讨论。

高斯滤波

为什么不用均值滤波

均值滤波给所有邻居相同权重，这在直觉上不合理——距离中心越近的像素，应该对输出贡献越大。均值滤波在原图细节较密时还会产生网格状伪影（边缘附近出现明暗条纹），称为振铃效应（Ringing Artifact）。

高斯核

高斯滤波器（Gaussian Filter）的核由二维高斯函数定义：

$$G(u, v) = \frac{1}{2\pi\sigma^2} \, \e^{-\frac{u^2 + v^2}{2\sigma^2}}$$

其中 $\sigma$ 是标准差，控制函数宽度。前面的常数因子保证连续高斯积分为 1；离散化后通常忽略它，先按公式计算权重再归一化为和为 1。一个常用的 $5 \times 5$、$\sigma = 1$ 的离散近似是：

$$G \approx \begin{bmatrix} 0.003 & 0.013 & 0.022 & 0.013 & 0.003 \\ 0.013 & 0.059 & 0.097 & 0.059 & 0.013 \\ 0.022 & 0.097 & 0.159 & 0.097 & 0.022 \\ 0.013 & 0.059 & 0.097 & 0.059 & 0.013 \\ 0.003 & 0.013 & 0.022 & 0.013 & 0.003 \end{bmatrix}$$

中心权重最大、向外迅速衰减，符合「近邻影响大」的直觉。

关键参数

标准差 $\sigma$ 控制平滑程度：$\sigma$ 越大，高斯越宽越平坦，纳入更多远处像素，平滑更强。$\sigma$ 是高斯滤波器最重要的参数，常被称为「尺度」（scale）。

核宽度 由经验法则决定：将核半宽设为约 $3\sigma$，覆盖 99.7% 的高斯面积——更窄会因截断而失真，更宽则徒增计算量。

import numpy as np
from scipy.ndimage import gaussian_filter

img_smoothed = gaussian_filter(img, sigma=1.5)

可分离性

如果二维核能写成一个列向量与一个行向量的外积，则称该核是可分离的（Separable）。此时二维卷积可以拆为两次一维卷积——先沿行用一维核滤波，再沿列用一维核滤波。

均值滤波器可分离：

$$\frac{1}{9} \begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{bmatrix} = \frac{1}{3} \begin{bmatrix} 1 \\ 1 \\ 1 \end{bmatrix} \cdot \frac{1}{3} \begin{bmatrix} 1 & 1 & 1 \end{bmatrix}$$

二维高斯天然可分离：

$$\frac{1}{2\pi\sigma^2} \e^{-\frac{x^2 + y^2}{2\sigma^2}} = \left[\frac{1}{\sqrt{2\pi}\sigma} \e^{-\frac{x^2}{2\sigma^2}}\right] \cdot \left[\frac{1}{\sqrt{2\pi}\sigma} \e^{-\frac{y^2}{2\sigma^2}}\right]$$

复杂度的差异显著：对 $n \times n$ 图像、$m \times m$ 核，直接卷积是 $O(n^2 m^2)$，分离后是 $O(n^2 m)$——核越大，加速比越明显。一个 $11 \times 11$ 的高斯核，可分离方案能把每像素 121 次乘加降为 22 次。

平滑滤波器的共性

均值滤波和高斯滤波都属于平滑滤波器（Smoothing Filter），共同性质：

• 核中所有权重为正
• 权重之和为 1（保证常数区域不变）
• 平滑程度与核大小成正比
• 抑制高频、保留低频——属于低通滤波器（Low-pass Filter）

后面在频率域将会看到，「保留低频、抑制高频」是平滑的等价描述。

中值滤波

线性滤波器在处理高斯噪声时表现良好，但面对椒盐噪声会失效——极端的 0/255 像素被「稀释」却向四周扩散，整张图变得污秽模糊。这时需要非线性滤波器。

中值滤波器（Median Filter）的规则极简：把邻域内所有像素值排序，取中间值作为输出。

中值滤波的计算

对于 $3 \times 3$ 邻域

$$\begin{bmatrix} 10 & 15 & 20 \\ 23 & 90 & 27 \\ 33 & 31 & 30 \end{bmatrix}$$

排序得 $10, 15, 20, 23, \boxed{27}, 30, 31, 33, 90$，中值为 27。异常值 90 被完全排除。

中值滤波有三个独特优势：

不引入新值。输出始终是邻域中已有的某个像素值，不会产生线性滤波器特有的「介于两值之间」的模糊过渡。

有效去除脉冲噪声。极端值排序后位于两端，中值不受影响——这正是线性滤波器难以做到的。

保持边缘锐利。在边缘处邻域大部分像素属于边缘某一侧，中值落在多数像素那一侧；不像均值那样把两侧值混合成渐变。

中值滤波是非线性的

中值不满足叠加性。反例：

$$A = \begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 2 \\ 2 & 2 & 2 \end{bmatrix}, \quad B = \begin{bmatrix} 0 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 0 \end{bmatrix}$$

$\text{med}(A) + \text{med}(B) = 1 + 0 = 1$，但 $\text{med}(A + B) = 2 \ne 1$。

因此中值滤波不能表示为卷积，也不能用 FFT 加速。

双边滤波

高斯滤波只考虑像素间的空间距离，平滑纹理的同时也会模糊边缘。双边滤波器（Bilateral Filter）的设计目标是平滑纹理而保留边缘——同时考虑空间邻近度和像素值相似度。

直觉来自一维信号：在阶跃边缘附近，普通高斯模糊把两侧像素一视同仁地平均，让锐利跳变变成斜坡。双边滤波加上一条规则——若邻居的像素值与中心差异很大（说明跨越了边缘），则降低其权重。

普通高斯模糊的公式：

$$I_{\bm{p}}^{\text{b}} = \sum_{\bm{q} \in \mathcal{S}} G_{\sigma_s}(\|\bm{p} - \bm{q}\|) \, I_{\bm{q}}$$

权重仅取决于像素 $\bm{p}, \bm{q}$ 之间的空间距离。双边滤波在此基础上增加值域权重：

$$I_{\bm{p}}^{\text{bf}} = \frac{1}{W_{\bm{p}}^{\text{bf}}} \sum_{\bm{q} \in \mathcal{S}} \underbrace{G_{\sigma_s}(\|\bm{p} - \bm{q}\|)}_{\text{空间权重}} \,\underbrace{G_{\sigma_r}(|I_{\bm{p}} - I_{\bm{q}}|)}_{\text{值域权重}} \, I_{\bm{q}}$$

其中 $W_{\bm{p}}^{\text{bf}}$ 是归一化因子。

• $G_{\sigma_s}(\|\bm{p} - \bm{q}\|)$：空间上越近权重越大，由 $\sigma_s$ 控制
• $G_{\sigma_r}(|I_{\bm{p}} - I_{\bm{q}}|)$：像素值越相似权重越大，由 $\sigma_r$ 控制

flowchart LR
    subgraph 高斯模糊
        A1["仅空间距离"] --> A2["边缘也被平滑<br/>→ 边缘模糊"]
    end
    subgraph 双边滤波
        B1["空间 + 值域"] --> B2["跨边缘像素<br/>值域权重 ≈ 0<br/>→ 边缘保留"]
    end

    classDef gauss fill:#ffebee,stroke:#c62828,stroke-width:2px
    classDef bilat fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    class A1,A2 gauss
    class B1,B2 bilat

双边滤波的核形状是空间变化的（spatially varying）——平坦区域中值域权重几乎不起作用，核近似标准高斯；边缘处值域权重大幅削弱另一侧像素的贡献，核变得不对称，避免跨边缘平滑。

双边滤波不擅长椒盐噪声

椒盐噪声的像素值与邻居差异极大，值域权重会把这些噪声点视为「边缘」而保留。对椒盐噪声仍应优先用中值滤波。

双边滤波的代价是计算量大：每像素核都不同，无法用 FFT 加速；暴力实现复杂度 $O(n^2 m^2)$，大图上可达数分钟。Durand 等人提出的快速近似算法能显著加速但伴随精度损失。

滤波器	处理高斯噪声	处理椒盐噪声	是否保边	是否线性
均值	良	差	否	是
高斯	优	差	否	是
中值	中	优	是	否
双边	优	差	是	否

锐化滤波

平滑滤波去掉高频信息；锐化的目标看起来像是反过来——但「只保留高频」会得到一张几乎全黑的图（高频本身能量小）。锐化的真正含义是增强高频信息，而非舍弃低频。

反锐化掩模

观察一个等式：原图 $f$ 减去其平滑版本 $f * g$，得到的就是「细节」（高频）：

$$\text{detail} = f - f * g$$

把细节加回原图，得到锐化结果：

$$f_{\text{sharp}} = f + \alpha \cdot (f - f * g) = (1 + \alpha) f - \alpha (f * g)$$

利用卷积的分配律和单位元（$f = f * e$），整理为单一卷积：

$$f_{\text{sharp}} = f * \bigl((1 + \alpha) e - \alpha g\bigr)$$

这种方法称为反锐化掩模（Unsharp Masking）——名字来源于摄影暗房，先做一张模糊（unsharp）的副本作为掩模，再用它增强原图。$\alpha$ 控制锐化强度。

构造一个 $3 \times 3$ 的近似锐化核（用单位脉冲乘以 2 减去均值核）：

$$\begin{bmatrix} 0 & 0 & 0 \\ 0 & 2 & 0 \\ 0 & 0 & 0 \end{bmatrix} - \frac{1}{9}\begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{bmatrix} = \begin{bmatrix} -1/9 & -1/9 & -1/9 \\ -1/9 & 17/9 & -1/9 \\ -1/9 & -1/9 & -1/9 \end{bmatrix}$$

核中心为正、四周为负，权重之和仍为 1（不改变常数区域亮度）——这是「拉普拉斯型」核的典型形态。

高斯的拉普拉斯算子

更精确的锐化使用高斯的拉普拉斯算子（Laplacian of Gaussian, LoG）作为核：先用高斯平滑，再求二阶导。LoG 在边缘处呈典型的「墨西哥帽」形态——中心负、环绕正、远处归零，正是高频成分的检测器。锐化等价于把图像与「单位脉冲减去 LoG 形状的核」做卷积。

LoG 在拉普拉斯金字塔和后续边缘检测中将多次出现。

混合图像

锐化的双重视角——把高频从原图中「分离出来」并叠加——直接催生了一个有趣的应用：混合图像（Hybrid Images）。

将图像 $A$ 的低频与图像 $B$ 的高频叠加：

$$I_{\text{hybrid}} = \text{LowPass}(A) + \text{HighPass}(B)$$

低通用高斯模糊实现，高通通过原图减去其低通版本得到。在频域中观察可以清晰看到：低频部分来自 $A$，高频部分来自 $B$。

这种图的奇妙之处在于：近看是 $B$，远看是 $A$。原因来自人眼的对比敏感度函数（Contrast Sensitivity Function, CSF）——人眼对中等空间频率最敏感。靠近图像时，物理高频对应较低的「每度视角周期数」，落在敏感区，于是 $B$ 的细节主导感知；远离时高频成分超出感知能力，仅低频可见，$A$ 浮现。

flowchart LR
    A["图像 A"] --> LP["高斯低通<br/>(LowPass)"]
    B["图像 B"] --> HP["原图减低通<br/>(HighPass)"]
    LP --> SUM["叠加"]
    HP --> SUM
    SUM --> H["混合图像<br/>近看 B，远看 A"]

    classDef src fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef proc fill:#fff3e0,stroke:#ef6c00,stroke-width:2px
    classDef out fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px

    class A,B src
    class LP,HP,SUM proc
    class H out

经典示例是 Oliva 等人 SIGGRAPH 2006 论文中爱因斯坦与玛丽莲·梦露的混合——近看是爱因斯坦，远看（或眯起眼睛）变成梦露。这种「同一张图、两种观感」的特性优雅地证明了空间频率与视觉感知之间的紧密联系。

频率域

到这里我们一直在空间域操作像素和邻域；接下来切换视角——把图像看作不同频率分量的叠加，在频率域（Frequency Domain）操作这些分量。这一视角能回答空间域分析难以解答的问题：为什么高斯滤波比均值滤波平滑更自然？降采样为什么会出现摩尔纹？混合图像为什么近看远看不同？

时域与频域

时域（Time Domain）以时间或空间作为自变量描述信号——能直观看到信号的起伏、峰值、持续时间。声音的波形图、心电图都是时域表示。

频域（Frequency Domain）以频率作为自变量——揭示信号的内部结构。一首音乐在时域是连续的声波起伏，在频域则是不同频率成分（低音、中音、高音）的强度分布。

Fourier 的核心思想

1807 年，Fourier 提出：任何周期函数都可以表示为不同频率的正弦和余弦函数的加权和。这一想法在当时备受争议——Lagrange、Laplace、Poisson 都曾质疑——直到 1878 年才被翻译为英文。但它（在满足一定技术条件下）基本正确，被称为傅里叶级数（Fourier Series）。

直觉上，正弦波 $A \sin(\omega x + \phi)$ 是基本构建块：$A$ 控制振幅、$\omega$ 控制频率、$\phi$ 控制相位。叠加足够多的不同正弦波，就能逼近任意信号。

方波的傅里叶级数

方波是奇数次谐波的叠加：

$$f(t) = \sum_{k=1, 3, 5, \ldots}^{\infty} \frac{1}{k} \sin(kt)$$

只用 $k=1$ 时是简单正弦曲线；叠加 $k=3$ 后波形顶部开始变平；继续叠加更高谐波，逐步逼近方波——「尖角」的凝聚需要无穷多个分量。

正余弦波在频域的几何意义是圆周运动在直线上的投影——频域的基本单元可以理解为以特定频率始终旋转的圆。

复指数

把正余弦统一在一个紧凑框架中的工具是欧拉公式（Euler's Formula）：

$$\e^{\mathrm{i}\theta} = \cos\theta + \mathrm{i}\sin\theta$$

复指数 $\e^{\mathrm{i}\theta}$ 在复平面上对应单位圆上一个旋转点——实部是余弦、虚部是正弦。离散时间下，长度 $N$ 的复指数信号为 $\e^{\mathrm{i} 2\pi k n / N}$，其中 $k$ 是频率（每 $N$ 个采样完成 $k$ 个周期）。复指数同时编码了振幅与相位，比单独用正弦更简洁——这也是它取代正弦成为频谱分析基本信号的原因。

推广到二维：

$$\e^{\mathrm{i} 2\pi (un / N + vm / M)} = \e^{\mathrm{i} 2\pi un/N} \cdot \e^{\mathrm{i} 2\pi vm/M}$$

二维复指数可分离——这是后续 2D DFT 可分离实现的基础。

傅里叶级数与傅里叶变换

傅里叶级数把周期信号分解为直流分量加一系列复指数信号之和：

$$f(t) = \sum_{k=-\infty}^{+\infty} c_k \e^{\mathrm{i} 2\pi k t / T}, \quad c_k = \frac{1}{T} \int_0^T f(t) \e^{-\mathrm{i} 2\pi k t / T} \d t$$

傅里叶变换把非周期的连续信号转换为频域的连续信号——可视作把周期推到无穷大的极限：

$$F(\omega) = \int_{-\infty}^{+\infty} f(t) \e^{-\mathrm{i} \omega t} \d t \quad \text{（正变换）}$$

$$f(t) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} F(\omega) \e^{+\mathrm{i} \omega t} \d \omega \quad \text{（逆变换）}$$

正逆变换的形式几乎对称——只差指数符号和归一化系数。

离散傅里叶变换

1D DFT

数字图像处理的对象是离散信号，对应离散傅里叶变换（Discrete Fourier Transform, DFT）。对长度 $N$ 的离散信号 $f[n]$：

DFT 与逆 DFT

正变换（空间域 → 频率域）：

$$F[u] = \sum_{n=0}^{N-1} f[n] \, \e^{-\mathrm{i} 2\pi u n / N}, \quad u = 0, 1, \ldots, N-1$$

逆变换（频率域 → 空间域）：

$$f[n] = \frac{1}{N} \sum_{u=0}^{N-1} F[u] \, \e^{+\mathrm{i} 2\pi u n / N}, \quad n = 0, 1, \ldots, N-1$$

$u$ 是频率索引，$F[u]$ 是复数，编码了频率分量 $u$ 的振幅和相位。直觉上 $F[u]$ 是信号 $f$ 与频率 $u$ 的复指数波的内积——内积越大，该频率分量越强。

一个简单的例子

取 $N = 4$，$f = [1, 2, 3, 4]$。$F[0] = 1+2+3+4 = 10$，除以 $N$ 得平均值 $2.5$——这就是直流分量（DC component）。$F[1]$ 衡量「每 4 采样点完成 1 周期」的频率成分有多强。

DFT 可以写成矩阵形式 $\hat{\bm{f}} = \bm{W} \bm{f}$，其中 $\bm{W}$ 的 $(u, n)$ 元为 $\e^{-\mathrm{i} 2\pi u n / N}$。矩阵每一行是一个频率的基函数：第 0 行是常数（直流），后续各行振荡频率递增。傅里叶变换本质上是一种基变换——把信号从「像素基」转换到「频率基」。

变换矩阵 $\bm{W}$ 是酉矩阵（Unitary Matrix）：$\bm{W}^{-1} = \bm{W}^*$（逆等于共轭转置）。酉变换保持向量长度不变——能量不增不减，变换完全可逆，不丢失也不创造信息。这是后面会用到的Parseval 定理的根源。

振幅与相位

DFT 输出 $F[u]$ 是复数，由实部 $R[u]$、虚部 $I[u]$ 组成，可分解为：

$$|F[u]| = \sqrt{R[u]^2 + I[u]^2} \quad \text{（振幅）}$$

$$\phi[u] = \arctan\frac{I[u]}{R[u]} \quad \text{（相位）}$$

振幅衡量频率 $u$ 的分量有多强，相位编码该分量在空间上的偏移。$F[u] = |F[u]| \e^{\mathrm{i} \phi[u]}$ 完整描述了该频率。

常见变换对

掌握几对常用变换有助于快速判断滤波器的频率特性：

空间域	频率域	说明
脉冲 $\delta$	常数 1	脉冲包含所有频率且振幅相等
常数 1	脉冲 $\delta$	常数信号仅有零频率分量
方波（Box）	sinc 函数	sinc 有无穷振荡的旁瓣
高斯函数	高斯函数	高斯在变换下保持形状
三角波（Tent）	$\text{sinc}^2$	两方波卷积；卷积对应频域乘法
梳状函数（Shah）	梳状函数	采样函数的变换仍是采样函数

高斯函数在傅里叶变换下仍为高斯——空间域宽高斯对应频率域窄高斯，反之亦然。这种「无旁瓣」的特性是高斯滤波器在图像处理中广受青睐的关键原因。

二维傅里叶变换

2D DFT 的定义

将 1D DFT 推广到 2D。对 $N \times M$ 图像 $f[n, m]$：

$$F[u, v] = \sum_{n=0}^{N-1} \sum_{m=0}^{M-1} f[n, m] \, \e^{-\mathrm{i} 2\pi (un/N + vm/M)}$$

逆变换：

$$f[n, m] = \frac{1}{NM} \sum_{u=0}^{N-1} \sum_{v=0}^{M-1} F[u, v] \, \e^{+\mathrm{i} 2\pi (un/N + vm/M)}$$

$(u, v)$ 是二维频率坐标——$u$ 对应垂直方向频率，$v$ 对应水平方向。由于复指数可分离，2D DFT 等价于先对每行做 1D DFT，再对每列做 1D DFT。

二维波

二维复指数 $\e^{\mathrm{i} 2\pi (un/N + vm/M)}$ 在空间中表现为特定方向的平面波：

• $(u, v) = (2, 0)$：纯垂直方向的波，垂直方向有 2 个完整周期
• $(u, v) = (3, 1)$：斜方向的波
• $(u, v) = (7, -5)$：更高频率的斜向波

模 $\sqrt{u^2 + v^2}$ 决定总频率高低，$\arctan(v/u)$ 决定波的传播方向。

振幅谱与相位谱

2D DFT 的结果通常可视化为振幅谱和相位谱两张图。常用约定是把零频率（直流分量）置于中心、频率向外递增（fftshift 操作）。直流分量通常远大于其他频率，一般取对数 $\log(1 + |F[u, v]|)$ 增强可视性。

import numpy as np
import matplotlib.pyplot as plt

img = plt.imread('image.png')[:, :, 0]      # 灰度化
F = np.fft.fft2(img)                         # 二维 DFT
F_shift = np.fft.fftshift(F)                 # 零频率移到中心
magnitude = np.log(1 + np.abs(F_shift))      # 对数振幅谱
phase = np.angle(F_shift)                    # 相位谱

观察振幅谱可以得到几个直觉：

平移不改变振幅谱。移动图像只改变各频率分量的相位，不改变振幅——图像的频率组成与空间位置无关。

空间尺度与频率互逆。图像中尺度越小的细节（精细纹理、锐利边缘），对应频率域中越高的频率。

方向性。含强水平边缘的图像，振幅谱在垂直轴上能量集中——水平边缘意味着亮度沿垂直方向快速变化。

一个有趣的实验是在频率域中选择性地把某些方向的频率分量置零再逆变换：把所有水平频率分量置零，逆变换后图像中的垂直边缘消失；反之亦然。这直观展示了不同方向的频率分量各自编码图像中对应方向的结构信息。

相位比振幅更重要

把图像 A 的振幅谱与图像 B 的相位谱组合并逆变换，得到的图像视觉上几乎完全像图像 B。这说明相位编码了物体的空间结构和轮廓——人类视觉最依赖的信息。振幅更多影响整体纹理和对比度，对「像什么」的贡献较小。

卷积定理

频域视角的最大威力来自一个核心定理：

卷积定理

空间域的卷积等价于频率域的逐元素乘法：

$$\mathcal{F}[g * h] = \mathcal{F}[g] \cdot \mathcal{F}[h]$$

反之，频率域的乘法等价于空间域的卷积：

$$\mathcal{F}^{-1}[G \cdot H] = \mathcal{F}^{-1}[G] * \mathcal{F}^{-1}[H]$$

flowchart LR
    subgraph 空间域
        A["图像 f"] -- "卷积 *" --> C["输出"]
        B["滤波器 h"] --> C
    end
    subgraph 频率域
        D["F = FFT(f)"] -- "逐元素相乘" --> F["G = F · H"]
        E["H = FFT(h)"] --> F
        F -- "IFFT" --> G["输出"]
    end

    A -.-> D
    B -.-> E

    classDef spatial fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef freq fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    classDef result fill:#fff3e0,stroke:#ef6c00,stroke-width:2px
    class A,B spatial
    class D,E freq
    class C,G result
    class F freq

实际意义

计算加速。直接卷积复杂度 $O(N^2 M^2)$（$N$ 图像边长、$M$ 核边长）。快速傅里叶变换（Fast Fourier Transform, FFT）能在 $O(N^2 \log N)$ 内做完正反变换。当核较大（$M^2 > \log N$）时，FFT 路线更快。

频率域设计滤波器。在空间域应用滤波器等价于在频率域将图像频谱乘以滤波器的频率响应（Frequency Response）。这给出了设计滤波器的新思路——直接指定要保留或抑制的频率范围（理想低通、理想高通、带通……）。

频率域视角下的滤波器

卷积定理为前面的观察提供了根本解释：为什么高斯滤波平滑自然，而均值滤波会产生方块状伪影？

高斯滤波器的傅里叶变换仍是高斯——一条平滑的钟形曲线，从低频到高频逐渐衰减、无任何振荡。高频被平滑地抑制。

均值滤波器的傅里叶变换是 sinc 函数——主瓣外有无穷个振荡的旁瓣。这些旁瓣意味着某些高频分量未被抑制甚至可能被放大，导致振铃效应——边缘附近出现波纹状伪影。

高斯函数的独特之处是它在空间域和频率域都没有旁瓣。此外空间域 $\sigma$ 与频率域宽度成反比：$\sigma$ 越大，频率域高斯越窄，低通效果越强——与「$\sigma$ 越大平滑越强」的直觉一致。

频率域滤波示例

在频域显式构造一个 Sobel 类的边缘检测核：

$$h = \begin{bmatrix} 1 & 0 & -1 \\ 2 & 0 & -2 \\ 1 & 0 & -1 \end{bmatrix}$$

它对水平方向上的亮度差异敏感（垂直边缘）。其 FFT 表现为水平方向上为低、两侧为高的带通形态——把图像的垂直边缘（水平方向的高频）提取出来。

相位的重要性再现

若只对图像振幅谱做某种调整（如裁剪低频或高频）后保留相位再逆变换，仍能识别原图；但若打乱相位、只保留振幅谱，逆变换得到的几乎是噪声。在频域设计滤波器时，振幅响应决定了「保留哪些频率」，但相位响应必须保持「线性相位」（否则边缘错位）——这是设计实用滤波器的隐性约束。

小波变换

傅里叶变换强大但有一个根本缺陷：它告诉你信号包含哪些频率，却不告诉你这些频率出现在哪里。整个时间或空间维度被「积分掉」了。对平稳信号（频率组成不随时间变化）这无所谓，但图像是充满突变的——同一张图像中，平坦区域只含低频、边缘区域含高频，频率组成随空间位置而变。傅里叶变换无法表达这种局部性。

flowchart LR
    A["纯傅里叶变换<br/>只有频率分辨率"] --> B["不知道频率<br/>出现在哪里"]
    A --> C["对突变和<br/>非平稳信号效果差"]

    classDef weak fill:#ffebee,stroke:#c62828,stroke-width:2px
    class A,B,C weak

短时傅里叶变换的折中

一个直观的解决方案是短时傅里叶变换（Short-Time Fourier Transform, STFT）：用一个窗函数把信号切成小段，对每段做傅里叶变换。结果是一个二维表示——横轴时间、纵轴频率——能告诉你「在某段时间内信号包含哪些频率」。

但 STFT 的窗口大小固定，存在 Heisenberg 式的取舍：窗口窄则时间分辨率高、但频率分辨率低（无法分辨低频）；窗口宽则反之。一个窗口尺寸适合的频率范围未必适合另一段。

小波变换（Wavelet Transform）的关键想法是：用可缩放、可平移的「小波形」（wavelet）作为基函数——分析高频时用窄窗（小尺度小波）保证时间分辨率，分析低频时用宽窗（大尺度小波）保证频率分辨率。不同尺度上「时间-频率」的精度自适应分配。

连续小波变换

定义一个母小波（Mother Wavelet）$\psi(t)$——一个均值为零、能量集中在有限区间的函数。通过缩放和平移得到一族小波：

$$\psi_{a, b}(t) = \frac{1}{\sqrt{|a|}} \, \psi\!\left(\frac{t - b}{a}\right)$$

其中 $a$ 是尺度参数（scale，控制小波的「宽度」）、$b$ 是位置参数（position，控制平移）。连续小波变换（Continuous Wavelet Transform, CWT）定义为：

$$C(a, b) = \int_{-\infty}^{+\infty} f(t) \, \psi_{a, b}(t) \, \d t$$

每个 $(a, b)$ 对应一个小波系数——衡量信号在位置 $b$、尺度 $a$ 上与小波的相似度。

离散小波变换

实际计算中使用离散小波变换（Discrete Wavelet Transform, DWT），$a, b$ 取离散值（典型为 2 的幂次：$a = 2^j$、$b = k \cdot 2^j$），构成正交基。Mallat 1989 年提出的多分辨率分解（Multiresolution Decomposition）框架是 DWT 的核心：通过一对镜像滤波器（低通 $h$ 和高通 $g$）+ 二倍下采样，把信号迭代分解为「概貌」（approximation）和「细节」（detail）两部分；下一层只对概貌继续分解。

flowchart LR
    S["信号 s"] --> H1["低通 h"]
    S --> G1["高通 g"]
    H1 --> D1["↓2"]
    G1 --> E1["↓2"]
    D1 --> A1["概貌 a₁"]
    E1 --> X1["细节 d₁"]
    A1 --> H2["低通 h"]
    A1 --> G2["高通 g"]
    H2 --> D2["↓2"]
    G2 --> E2["↓2"]
    D2 --> A2["概貌 a₂"]
    E2 --> X2["细节 d₂"]

    classDef sig fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef filt fill:#fff3e0,stroke:#ef6c00,stroke-width:2px
    classDef out fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    class S sig
    class H1,G1,H2,G2,D1,E1,D2,E2 filt
    class A1,X1,A2,X2 out

每层下采样使尺寸减半——经过 $J$ 层得到的系数总数等于原信号长度（无冗余的正交分解）。重建时反向操作：上采样后用对偶滤波器卷积叠加。

Haar 小波

最简单也最有教学价值的小波是 Haar 小波。1909 年提出，比 Fourier 变换更早地具备「时间局部性」。

定义两个基本函数：

尺度函数（scaling function，又称父小波）$\phi(x)$：

$$\phi(x) = \begin{cases} 1 & 0 \le x < 1 \\ 0 & \text{其他} \end{cases}$$

小波函数（wavelet function，又称母小波）$\psi(x)$：

$$\psi(x) = \begin{cases} 1 & 0 \le x < 1/2 \\ -1 & 1/2 \le x < 1 \\ 0 & \text{其他} \end{cases}$$

通过缩放（$2^j$ 倍压缩）和平移（$i$ 个单位）生成一族基函数：

$$\phi^j_i(x) = \phi(2^j x - i), \quad \psi^j_i(x) = \psi(2^j x - i), \quad i = 0, 1, \ldots, 2^j - 1$$

这族函数构成 $L^2[0, 1]$ 的一组正交基。

Haar 变换的过程

Haar 变换有一个直观的算法：对信号每两个相邻值取「均值 + 半差」。

离散 Haar 变换的步骤

考虑信号 $[9, 7, 3, 5]$（分辨率 4）。

分辨率	平均值	细节系数
4	$[9, 7, 3, 5]$	—
2	$[8, 4]$	$[1, -1]$
1	$[6]$	$[2]$

每一步两两组合：

• 第 1 步：$(9+7)/2 = 8$，$(9-7)/2 = 1$；$(3+5)/2 = 4$，$(3-5)/2 = -1$
• 第 2 步：$(8+4)/2 = 6$，$(8-4)/2 = 2$

最终的 Haar 系数为 $[6, 2, 1, -1]$——一个直流系数 + 三个细节系数。

每一层的「平均」对应概貌系数，「半差」对应细节系数。展开为基函数表示：

$$I(x) = c_0^0 \phi_0^0(x) + d_0^0 \psi_0^0(x) + d_0^1 \psi_0^1(x) + d_1^1 \psi_1^1(x)$$

最初用 4 个 $\phi^2_i$ 表示，逐层「升级」：先把两个 $\phi^2_i$ 合并成一个 $\phi^1_i$ 加一个 $\psi^1_i$，再合并到 $\phi^0_0$ 加 $\psi^0_0$。每次合并，「细节」从概貌中被分离出来存到对应尺度的小波系数。

变换是完美可逆的——给定所有概貌和细节系数，可以从最粗尺度逐层恢复出原始信号。

二维 Haar 小波变换

二维 Haar 变换有两种分解策略。

标准分解（Standard Decomposition）：先对每一行做完整的一维变换，再对每一列做完整的一维变换。

procedure StandardDecomposition(C: array[1..h, 1..w] of reals)
    for row = 1 to h do
        Decomposition(C[row, 1..w])         // 行方向变换
    for col = 1 to w do
        Decomposition(C[1..h, col])         // 列方向变换
end procedure

非标准分解（Non-standard Decomposition）：每次行列各做一步变换，然后把图像分成 4 块（LL/LH/HL/HH），只对 LL 子块继续递归。

procedure NonstandardDecomposition(C: array[1..h, 1..h] of reals)
    C ← C / h
    while h > 1 do
        for row = 1 to h do DecompositionStep(C[row, 1..h])
        for col = 1 to h do DecompositionStep(C[1..h, col])
        h ← h / 2
end procedure

非标准分解的结果在视觉上更直观——把图像分解为四个子带：

flowchart TD
    subgraph 一层分解
        LL["LL<br/>低-低<br/>水平+垂直<br/>都低频<br/>（缩略图）"]
        LH["LH<br/>低-高<br/>水平低频<br/>垂直高频<br/>（水平边缘）"]
        HL["HL<br/>高-低<br/>水平高频<br/>垂直低频<br/>（垂直边缘）"]
        HH["HH<br/>高-高<br/>水平+垂直<br/>都高频<br/>（对角细节）"]
    end

    classDef ll fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef lh fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    classDef hl fill:#fff3e0,stroke:#ef6c00,stroke-width:2px
    classDef hh fill:#f3e5f5,stroke:#4a148c,stroke-width:2px
    class LL ll
    class LH lh
    class HL hl
    class HH hh

LL 子块继续递归分解，给出多分辨率金字塔——这正是 JPEG 2000 等图像压缩标准使用小波而非 DCT 的原因之一。

import pywt
import numpy as np

img = np.random.rand(256, 256)
coeffs = pywt.wavedec2(img, 'haar', level=3)
# coeffs[0] 是最粗的近似系数 LL3
# coeffs[1..] 是各层的细节系数 (LH, HL, HH)

reconstructed = pywt.waverec2(coeffs, 'haar')

Gabor 小波

Haar 小波在边缘附近能精确定位但在频率上不锐利（方波在频域的 sinc 旁瓣很大）。Gabor 小波（Gabor Wavelet）是另一类——用高斯包络调制的复指数：

$$\psi(x, y) = \e^{-(x^2 + y^2) / (2\sigma^2)} \cdot \e^{\mathrm{i} 2\pi u_0 x}$$

取实部和虚部分别得到余弦（偶对称）和正弦（奇对称）两种滤波器：

$$\psi_e(x, y) = \e^{-(x^2 + y^2)/(2\sigma^2)} \cos(2\pi u_0 x)$$

$$\psi_s(x, y) = \e^{-(x^2 + y^2)/(2\sigma^2)} \sin(2\pi u_0 x)$$

参数解读：$\sigma$ 控制空间范围（高斯包络宽度），$u_0$ 控制空间频率（正弦波密度）。通过旋转坐标 $x' = x\cos\alpha + y\sin\alpha$ 可得到任意方向的 Gabor 滤波器；不同 $u_0$ 对应不同空间频率。

Gabor 滤波器的独特价值在于它在空间-频率联合不确定性原理下达到了最小不确定性——既不能再精确定位空间也不能再精确定位频率，是两者乘积的下界。这让它成为图像分析中「同时关心位置和频率」的最佳工具之一。

Gabor 与人脑

实验表明，人类初级视觉皮层 V1 简单细胞的感受野与二维 Gabor 函数的拟合几乎完美——拟合残差与噪声无法区分。这说明人脑在视觉处理早期阶段本质上就是在做多尺度、多方向的带通滤波。Gabor 滤波器组是纹理分析、人脸识别（如 Gabor jet）等任务的经典选择。

采样与混叠

频域视角还能解决另一个根本问题：怎么正确地把图像变小？

降采样的问题

降采样（Downsampling）是减小图像分辨率的操作。最朴素的 2 倍降采样是丢弃每隔一行一列的像素，得到原图 $1/4$ 大小的图像。但直接降采样可能产生严重伪影——混叠（Aliasing），表现为锯齿、摩尔纹、电影中车轮倒转等。

频域看采样

从频域理解为什么会出现混叠：

1. 原始信号的频谱占据一定带宽——假设最高频率为 $f_{\max}$，频谱在 $[-f_{\max}, f_{\max}]$ 范围内
2. 采样操作等价于在空间域将信号乘以一个梳状函数（Shah / Dirac comb，间隔 $T_s$ 的脉冲序列）。由卷积定理，空间域乘法对应频率域卷积；梳状函数的傅里叶变换仍是梳状函数，因此采样后的频谱变成原始频谱沿频率轴的周期性复制——以采样频率 $f_s = 1/T_s$ 为间隔无限重复
3. 采样率足够（$f_s > 2 f_{\max}$）时，频谱副本之间有间隔不会重叠；可以用一个理想低通滤波器「切出」中心那一份频谱，完美重建原始信号
4. 采样率不足（$f_s < 2 f_{\max}$）时，频谱副本发生重叠——高频分量「伪装」成低频混入原始频谱，无法再被分离。这就是混叠，一种不可逆的信息损坏

奈奎斯特-香农采样定理

采样定理

以离散间隔对连续信号采样时，采样频率必须大于信号最高频率的两倍：

$$f_s \ge 2 f_{\max}$$

满足此条件即可从采样信号完美重建原始信号。临界频率 $2 f_{\max}$ 称为奈奎斯特频率（Nyquist Frequency）。

直觉上，捕获频率为 $f$ 的正弦波至少需要每周期 2 个采样点（一峰一谷），否则该正弦波与一个更低频率的正弦波不可区分——这正是混叠的本质。

混叠在日常生活中很常见：

• 电影中车轮看似向后转——时间采样率（帧率）不足以捕捉旋转频率
• 渲染中棋盘格远处崩解为乱花——空间采样率不足以表示高频纹理
• 旧电视上条纹衬衫出现奇怪花纹——衬衫的空间频率超出显示器采样能力

抗混叠

混叠源于信号中存在超出采样能力的高频分量。两种解决方案：

提高采样率：采集更多样本，但受硬件限制，往往不可行。

预滤波（Pre-filtering）：在降采样前用低通滤波去除超过新采样频率一半的频率分量。有控制地丢失高频远好于让高频混叠到低频产生伪影——这就是抗混叠（Anti-aliasing）。

2 倍降采样的标准算法：

from scipy.ndimage import gaussian_filter

# 1. 高斯低通预滤波（抗混叠）
blurred = gaussian_filter(img, sigma=1.0)
# 2. 每隔一行一列取样
downsampled = blurred[::2, ::2]

不做预滤波时，连续 1/2、1/4、1/8 降采样会让图像出现严重锯齿和摩尔纹；做了高斯预滤波后图像虽然更模糊但视觉自然、无伪影。

模板匹配

回到本讲开头提到的滤波第三种视角：把滤波器看作模板（Template），将图像中的每个局部区域与模板比较，找到最相似的位置——这就是计算机视觉中「找东西」最朴素的方法。

核心问题是：如何衡量两个图像块的相似度？

互相关匹配

最直接的方案：把模板 $g$ 当作滤波器，对图像 $f$ 做互相关：

$$h[m, n] = \sum_{k, l} g[k, l] \, f[m+k, n+l]$$

$h[m, n]$ 给出每个位置上模板与对应图像块的「匹配得分」——得分越高越相似。

但直接互相关有个致命问题：对亮度敏感。一个明亮的区域即使与模板毫不相似，也可能因像素值本身大而得到很高得分——大范围天空区域往往比真正的目标得分更高。

零均值互相关

改进：先减去模板的均值，让匹配关注形状而非亮度：

$$h[m, n] = \sum_{k, l} (g[k, l] - \bar{g}) \, f[m+k, n+l]$$

减去均值后，模板的正值区域表示「比平均亮」、负值区域表示「比平均暗」，互相关变成了形状匹配。仍有误检——它没考虑图像局部区域自身的亮度和对比度。

平方差之和

平方差之和（Sum of Squared Differences, SSD）从另一角度衡量相似度——直接计算逐像素差异：

$$h[m, n] = \sum_{k, l} \bigl(g[k, l] - f[m+k, n+l]\bigr)^2$$

SSD 越小越相似。展开为：

$$h[m, n] = \sum g^2 - 2 \sum g f + \sum f^2$$

第一项常数、第二项可用线性滤波计算、第三项依赖图像局部能量需额外计算。SSD 的弱点是对整体亮度变化敏感——区域整体偏亮或偏暗时，即使形状完全匹配 SSD 也可能很大。

归一化互相关

归一化互相关（Normalized Cross-Correlation, NCC）是模板匹配的「黄金标准」——同时消除亮度和对比度影响：

$$h[m, n] = \frac{\displaystyle \sum_{k, l} (g[k, l] - \bar{g})\bigl(f[m+k, n+l] - \bar{f}_{m, n}\bigr)}{\left[\displaystyle \sum_{k, l} (g[k, l] - \bar{g})^2 \cdot \sum_{k, l} \bigl(f[m+k, n+l] - \bar{f}_{m, n}\bigr)^2\right]^{0.5}}$$

其中 $\bar{f}_{m, n}$ 是以 $(m, n)$ 为中心的图像块局部均值。NCC 取值在 $[-1, 1]$，1 为完美匹配。分子减均值消除亮度偏移，分母归一化消除对比度差异。

方法	速度	抗亮度变化	抗对比度变化
互相关	最快	差	差
零均值互相关	快	中	差
SSD	较快	差	差
NCC	最慢	好	好

选择哪种方法

没有一种方法在所有场景下都最优。光照稳定时，零均值互相关或 SSD 的速度优势更重要；光照变化大时，NCC 是最可靠的选择。MATLAB 的 normxcorr2 提供了 NCC 的开箱即用实现。

模板匹配假设模板与目标尺度相同，对尺度变化无鲁棒性——这正是下一节图像金字塔要解决的问题。

高斯金字塔

为什么需要多尺度

模板匹配的根本问题是它假设模板与图像中的目标尺度相同。但现实中同一物体可能出现在不同距离上，在图像中占据不同大小。仅靠平移不变性不够，还需要尺度不变性（Scale Invariance）。

朴素方案是准备多个不同大小的模板，但模板尺度连续变化无法穷举。更聪明的做法是固定模板大小，把图像本身缩放到不同分辨率，在每个分辨率上做匹配——这就是图像金字塔（Image Pyramid）的思想。

构建高斯金字塔

高斯金字塔（Gaussian Pyramid）的构建非常简单：对原图反复执行「高斯模糊 + 2 倍降采样」。

import cv2
import numpy as np

def build_gaussian_pyramid(img, levels):
    pyramid = [img]
    for _ in range(levels - 1):
        blurred = cv2.GaussianBlur(pyramid[-1], (5, 5), 1)
        downsampled = blurred[::2, ::2]
        pyramid.append(downsampled)
    return pyramid

用 $g_k$ 表示第 $k$ 层、$G_k$ 表示「模糊+降采样」算子，则：

$$g_1 = G_0 \, g_0, \quad g_2 = G_1 \, g_1, \quad g_3 = G_2 \, g_2, \quad \ldots$$

经典的滤波核采用 $[1, 4, 6, 4, 1] / 16$——一个近似高斯的可分离核。每一层的分辨率是上一层的一半，图像越来越小，故称「金字塔」。

金字塔的应用

由粗到精配准（Coarse-to-Fine Registration）是高斯金字塔最经典的应用：先在最粗（最小）层做全局搜索——图像很小、搜索空间小、速度快；然后逐层向上传递，仅在局部范围内微调。这种策略不仅快得多，还能避免陷入局部最优。

flowchart LR
    L4["最粗层 g₄<br/>32×32"] --> M4["全局搜索<br/>粗略对齐"]
    M4 --> L3["g₃<br/>64×64"]
    L3 --> M3["局部微调"]
    M3 --> L2["g₂<br/>128×128"]
    L2 --> M2["进一步微调"]
    M2 --> L1["原图<br/>512×512"]
    L1 --> M1["最终精对齐"]

    classDef coarse fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef fine fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    class L4,L3,L2 coarse
    class L1 fine
    class M4,M3,M2,M1 coarse

多尺度模板匹配：

1. 在当前尺度上匹配模板
2. 对图像降采样（实践中缩放步长取 1.1 ~ 1.2 而非 2，以避免漏掉中间尺度）
3. 重复 1-2 直到图像很小
4. 对所有响应做阈值筛选和非极大值抑制

从图像金字塔到特征金字塔

经典金字塔在每个尺度上独立提取特征效率低。深度学习时代提出了更巧妙的方案。

CNN 本身天然产生层次化特征——浅层分辨率高但语义弱、深层分辨率低但语义强。特征金字塔网络（Feature Pyramid Network, FPN）利用这一点：通过自顶向下的路径将深层强语义信息传递回高分辨率层，让每一层都兼具高分辨率和强语义。

flowchart LR
    subgraph 自底向上
        C2["C₂<br/>高分辨率<br/>弱语义"] --> C3["C₃"] --> C4["C₄"] --> C5["C₅<br/>低分辨率<br/>强语义"]
    end
    subgraph 自顶向下
        P5["P₅"] --> P4["P₄"] --> P3["P₃"] --> P2["P₂<br/>高分辨率<br/>强语义"]
    end

    C5 -.-> P5
    C4 -.-> P4
    C3 -.-> P3
    C2 -.-> P2

    classDef bottom fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef top fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    class C2,C3,C4,C5 bottom
    class P2,P3,P4,P5 top

FPN 的关键是把深层特征 2 倍上采样后加上经过 $1 \times 1$ 卷积对齐通道数的浅层特征——这个简单的 add 操作就完成了语义信息传递。金字塔结构后来也被引入 Vision Transformer：Pyramid Vision Transformer（PVT）逐阶段降低分辨率构建层次化特征，Swin Transformer 通过 patch merging 实现类似效果——这些架构证明了多尺度表示不仅适用于 CNN，也是 Transformer 处理密集预测任务的基本需求。

拉普拉斯金字塔

从低通到带通

高斯金字塔的每一层都是原图的低通版本——包含该尺度及以下的所有频率信息。但若想知道每一层新增了哪些信息——只属于该尺度的频率分量——就需要做减法。

拉普拉斯金字塔（Laplacian Pyramid）的第 $k$ 层定义为高斯金字塔相邻两层之差：

$$l_k = g_k - \text{Upsample}(g_{k+1})$$

由于 $g_{k+1}$ 分辨率只有 $g_k$ 的一半，需要先上采样到相同大小再做差。上采样的方法是先插入零值（使尺寸翻倍），再用平滑核插值——这一步本身又是一次卷积。

拉普拉斯金字塔的本质

每一层 $l_k$ 捕获的是 $g_k$ 相比 $g_{k+1}$ 多出的细节信息——高斯金字塔逐步丢弃高频，拉普拉斯金字塔则把这些被丢弃的高频「存」起来。在频率域中，每一层对应一个环形带通滤波器（annular band-pass filter），不同层的频段互不重叠。

用矩阵形式更清楚。令 $G_k$ 为模糊+降采样算子、$F_k$ 为上采样+模糊算子，则：

$$l_k = (I - F_k G_k) \, g_k$$

$I - F_k G_k$ 就是带通滤波矩阵——对 $g_k$ 先降采样再上采样会丢失细节，原图减去这个「丢失细节后的版本」就得到了细节本身。

LoG 与 DoG

「拉普拉斯」金字塔的名字来源于高斯的拉普拉斯算子（Laplacian of Gaussian, LoG）：先用高斯平滑再求二阶拉普拉斯导数。

$$\text{LoG}\{I; \sigma\} = \nabla^2 (G_\sigma * I) = (\nabla^2 G_\sigma) * I, \quad \nabla^2 = \frac{\pd^2}{\pd x^2} + \frac{\pd^2}{\pd y^2}$$

LoG 在频域是带通的（既抑制低频也抑制高频），在空间域形似「墨西哥帽」。实际计算中常用高斯差分（Difference of Gaussians, DoG）来近似 LoG：

$$\text{DoG}\{I; \sigma_1, \sigma_2\} = G_{\sigma_1} * I - G_{\sigma_2} * I = (G_{\sigma_1} - G_{\sigma_2}) * I$$

两个宽度不同的低通滤波器相减就是带通滤波——这正是相邻金字塔层相减的几何含义。DoG 计算更简单，是半倍频 LoG 的良好近似（SIFT 等特征检测器的核心算子）。

完美重建

拉普拉斯金字塔的关键性质是可完美重建（Perfect Reconstruction）：给定所有拉普拉斯层 $l_0, l_1, \ldots, l_{N-1}$ 和最粗层的高斯残差 $g_N$，可以精确恢复原图。

重建过程自底向上：

$$g_k = l_k + \text{Upsample}(g_{k+1})$$

从最粗层 $g_N$ 开始，逐层加上对应的拉普拉斯细节，最终恢复原图 $g_0$。拉普拉斯金字塔是一种无损的多尺度分解——把图像拆成不同频段的分量，任何时候都能完整拼回去。

flowchart LR
    subgraph 分解 Analysis
        I["原图 g₀"] --> B0["模糊+降采样"]
        B0 --> G1["g₁"]
        G1 --> U1["上采样"]
        U1 --> S1["I - U(g₁)"]
        I --> S1
        S1 --> L0["拉普拉斯 l₀"]
    end
    subgraph 重建 Synthesis
        G1B["g₁"] --> U1B["上采样"]
        L0B["l₀"] --> A1["相加"]
        U1B --> A1
        A1 --> R0["恢复 g₀"]
    end

    classDef ana fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef syn fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    class I,B0,G1,U1,S1,L0 ana
    class G1B,U1B,L0B,A1,R0 syn

图像融合

拉普拉斯金字塔最经典的应用之一是图像融合（Image Blending）。

直接拼接两张图像，接缝处会出现明显边界。最朴素的方案是用渐变 $\alpha$ 遮罩：$I = m \cdot I^A + (1 - m) \cdot I^B$。但渐变窗口大小很难选——窗口太宽，两张图相互「透」过来产生鬼影（ghosting）；窗口太窄，接缝仍可见。

深刻洞察：不同频率的内容需要不同宽度的混合窗口。高频细节（如边缘）需要窄窗口才不会模糊，低频内容（如整体亮度）需要宽窗口才能平滑过渡。拉普拉斯金字塔恰好把图像按频段拆分——每个频段可以用各自合适的窗口分别混合。

拉普拉斯金字塔融合算法

1. 分别构建 $I^A$、$I^B$ 的拉普拉斯金字塔 $L^A$、$L^B$
2. 构建遮罩 $m$ 的高斯金字塔 $M$（遮罩在每层都被平滑了不同程度）
3. 在每层用对应遮罩混合：$l_k = l_k^A \cdot m_k + l_k^B \cdot (1 - m_k)$
4. 从混合后的拉普拉斯金字塔重建出最终图像

这就是 Burt 和 Adelson 1983 年的经典「苹果 + 橙子 = orapple」实验：苹果和橙子各自分解为金字塔，每个频段用不同宽度的遮罩混合，重建后得到自然无缝过渡——肉眼难以察觉的「水果嫁接」。

混合图像的金字塔视角

前面讨论的混合图像也可以用拉普拉斯金字塔重新理解：把混合图像分解为拉普拉斯金字塔后，高频层显示高频源图像（如猫的轮廓），低频层显示低频源图像（如骷髅的整体形态）——不同频段来自不同原图。两张图像的「拼接」是在频域完成的。

可控金字塔

方向的缺失

高斯和拉普拉斯金字塔只分解了尺度信息——每一层对所有方向一视同仁。但许多任务需要同时感知尺度和方向：检测特定朝向的纹理、分析边缘走向、估计运动方向等。拉普拉斯金字塔的带通滤波器在频域中是环形的（各向同性），要分辨方向就需要把环切成若干「扇形」。

高斯导数与方向

二维高斯函数 $g(x, y; \sigma)$ 对 $x$ 求偏导：

$$g_x(x, y; \sigma) = \frac{\pd g}{\pd x} = \frac{-x}{2\pi \sigma^4} \, \e^{-(x^2 + y^2) / (2\sigma^2)}$$

$g_x$ 对垂直方向的边缘有强响应；类似地，$g_y$ 对水平方向的边缘有强响应。其他方向呢？

可控滤波器

关键的数学事实是：沿任意方向 $\alpha$ 的平滑方向梯度，是 $g_x$ 和 $g_y$ 响应的线性组合。

$$\bm{u}^{\intercal} \nabla g \otimes I = \cos\alpha \cdot (g_x \otimes I) + \sin\alpha \cdot (g_y \otimes I)$$

这意味着只需计算两次滤波（$g_x \otimes I$ 和 $g_y \otimes I$），就能合成出任意方向的滤波结果。这种性质称为可控性（Steerability）。

合成 30° 方向的滤波器

$$G_{30°}^1 = \cos(30°) \cdot G_0^1 + \sin(30°) \cdot G_{90°}^1$$

只需预计算 0° 和 90° 两个基滤波器的响应图，然后对每像素按上式线性组合，即可得到 30° 方向的响应——不需要重新做完整滤波。

可控金字塔

可控金字塔（Steerable Pyramid）将可控滤波器与多尺度分解结合：在每个尺度上应用一组不同方向的带通滤波器（$B_0, B_1, \ldots, B_n$），得到该尺度下各方向的分量；然后对残余低频信号降采样，进入下一层继续分解。

可控金字塔有几个重要性质：

保持所有信息——是一种过完备表示（Overcomplete Representation），系数总数多于原始像素数。但好处是可从系数完美重建原图（自逆性），且各频段/方向的子带之间没有混叠。

独立通道——每个子带表示特定尺度和方向上的信息，可独立操作。比如增强某方向的纹理而不影响其他方向。

方向分析——通过比较各方向子带的响应强度，可以得到每个位置的主方向（dominant orientation），生成方向场。

flowchart TD
    G["高斯金字塔<br/>逐步模糊+降采样<br/>→ 尺度不变性"]
    L["拉普拉斯金字塔<br/>相邻层之差<br/>→ 带通分解 + 融合"]
    S["可控金字塔<br/>尺度 × 方向<br/>→ 纹理/特征分析"]

    G --> L --> S

    classDef gauss fill:#e3f2fd,stroke:#1565c0,stroke-width:2px
    classDef lap fill:#e8f5e8,stroke:#2e7d32,stroke-width:2px
    classDef steer fill:#f3e5f5,stroke:#4a148c,stroke-width:2px
    class G gauss
    class L lap
    class S steer

可控金字塔的应用包括纹理合成、噪声去除、运动分析、运动放大（Eulerian Video Magnification）等。

三种视角的统一

回顾整讲，可以看到滤波的三种视角并非独立——它们是同一回事的不同表述：

表示	空间信息	频率信息	适用场景
像素	精确	无	直接显示、像素级操作
傅里叶变换	无	精确	全局频率分析、滤波器设计
金字塔 / 小波	局部	频段	多尺度分析、纹理、特征提取

像素表示和傅里叶变换走向了两个极端——一个只有空间信息、一个只有频率信息。金字塔和小波在两者之间找到平衡，既保留空间又保留频率信息（虽然都不精确）——这种「既知道在哪里，又知道是什么频率」的能力，正是它们在计算机视觉中如此重要的原因。

下一讲我们将把这些工具用于一个具体的视觉任务：边缘检测——从图像中提取最朴素也最基本的几何结构。