首页笔记，计算机科学，计算机图形学

# 滤波与卷积

这里区分一下滤波和卷积，以前我也是一直认为 滤波=滤波核+卷积操作(可能是还没学数字信号处理吧，目前就学了复变，时间不够用啊，悲)

不过实际上并非如此

我们仔细梳理一下

滤波的定义我认为可以直接简单的理解为对于信号的变换，不用去考虑什么增强削弱特征之类的定义

y(n)=F{x(n)}y(n)=\mathfrak{F}\{x(n)\}

线性滤波就是满足线性性质的滤波，经常能看到线性滤波等价于卷积计算的内容，可能是工程上近似的原因吧，严谨一点还需要满足 时不变 和 有界。根本上讲就是滤波算子需要满足可与积分换序

我们来推导一下

首先是证明卷积运算是线性滤波

这很容易，从积分或求和的线性性质一步到位

然后是证明可换序的线性滤波是卷积运算

我们假设算子为 T{⋅},y(n)=T{x(n)}T\{·\},\quad y(n)=T\{x(n)\}

将x用狄拉克函数展开

y(n)=T{∫−∞∞x(t)δ(n−t)dt}y(n)=T\{\int_{-\infin}^\infin x(t) \delta(n-t)dt\}

由有界线性算子可与积分换序得到

y(n)=∫−∞∞T{x(t)δ(n−t)}dt=∫−∞∞x(t)T{δ(n−t)}dty(n)=\int_{-\infin}^\infin T\{x(t)\delta(n-t)\}dt=\int_{-\infin}^\infin x(t)T\{\delta(n-t)\}dt

由于算子是时不变的，不妨令h(t)=T{δ(t)}h(t)=T\{ \delta(t)\}，那么

y(n)=∫−∞∞x(t)h(n−t)dty(n)=\int_{-\infin}^\infin x(t) h(n-t)dt

也就是卷积，其中卷积核是 h(t)=T{δ(t)}h(t)=T\{ \delta(t)\}
离散的话就更简单了，不过思路是类似的

卷积核是对单位冲激函数应用滤波的算子

这从卷积叠加累计的角度也可以理解

现在，我们就可以用卷积的角度去理解这些滤波(有界线性滤波之后就简写了，大家应该可以理解)了，比如时域卷积等于频域乘积角度考虑等等

那么非线性滤波呢

有一些也还是可以用卷积实现的，另一些比如取中位数什么的就很难用卷积表示了

# 案例

# 高斯滤波

# 核心思想

高斯滤波是一种经典的线性平滑滤波。它根据像素与中心像素的空间几何距离来分配权重。距离越近，权重越大；距离越远，权重越小。

# 数学原理

二维高斯核的概率密度函数为：

G(x,y)=12πσ2e−x2+y22σ2G(x, y) = \frac{1}{2\pi\sigma^2} e^{-\frac{x^2 + y^2}{2\sigma^2}}

其中 (x,y)(x,y) 是邻域像素相对于中心像素的坐标，σ\sigma（标准差）决定了高斯函数的钟形曲线宽度（平滑剧烈程度）。

# 优缺点与应用

  • 优点：能有效消除高频噪声（如高斯噪声），计算速度快（具有可分离性，由于正态分布在各维度独立的情况下等价于多个一维分布的乘积，因此可拆分为一维横向和纵向滤波，O(n2)−>O(n)O(n^2)->O(n)）

  • 缺点：不保留边缘。由于它只考虑空间距离，在图像边缘处会将前景和背景的像素一起模糊掉，导致边缘变秃、细节丢失。

# 双边滤波 (Bilateral Filter)

# 核心思想

双边滤波是一种非线性滤波，旨在实现保边降噪（Edge-Preserving Smoothing）。它在分配权重时不仅考虑空间几何距离，还考虑像素色彩值相似度。

没错，这就是可以用卷积考虑的非线性滤波了

非线性性质主要由于它还考虑了亮度因素

普通卷积：

output=∑kernel×pixeloutput=∑kernel×pixel

双边：

output=∑kernel(I)×pixeloutput=∑kernel(I)×pixel

所以更精确的说应该就只是加权邻域求和，不过从计算角度考虑倒是差不多

# 数学原理

双边滤波的输出值由两个核（Kernels）的乘积决定：

Ifiltered(p)=1Wp∑q∈SGσs(∥p−q∥)⋅Gσr(∥I(p)−I(q)∥)⋅I(q)I_{\text{filtered}}(p) = \frac{1}{W_p} \sum_{q \in S} G_{\sigma_s}(\|p - q\|) \cdot G_{\sigma_r}(\|I(p) - I(q)\|) \cdot I(q)

其中 WpW_p 是归一化系数，pp 为中心像素，qq 为邻域像素。

  1. 空间邻近度因子（Space Kernel）：

    Gσs=e−∥p−q∥22σs2G_{\sigma_s} = e^{-\frac{\|p - q\|^2}{2\sigma_s^2}}

    （与高斯滤波相同，距离越近权重越高）

  2. 亮度相似度因子（Range Kernel）：

    Gσr=e−∥I(p)−I(q)∥22σr2G_{\sigma_r} = e^{-\frac{\|I(p) - I(q)\|^2}{2\sigma_r^2}}

    （关键点：如果邻域像素 qq 与中心像素 pp 的颜色差异极大，说明可能跨越了边缘，GσrG_{\sigma_r} 会趋近于 0，从而保护边缘不被模糊）

# 优缺点与应用

  • 优点：在平滑区域内降噪的同时，能够完美保留陡峭的边缘。

  • 缺点：计算复杂度高，无法拆分为一维滤波；在边缘两侧噪声较大时，Range Kernel 可能会失效，产生类似“光晕（Halos）”的伪影。

# 联合双边滤波 (Joint Bilateral Filter)

# 核心思想

考虑更多项，比如深度信息、法线等等

数学本质和上面的一样

C~p=∑q∈Ωw(p,q)⋅Cq∑q∈Ωw(p,q)\tilde{C}_p = \frac{\sum_{q \in \Omega} w(p, q) \cdot C_q}{\sum_{q \in \Omega} w(p, q)}

其中 Ω\Omega 是滤波窗口（核大小），w(p,q)w(p, q) 是综合权重。核心就在于如何计算这个权重 w(p,q)w(p, q)

为了在模糊噪点的同时不糊掉物体的边缘（Edge-Avoiding），权重 w(p,q)w(p, q) 由位置/深度、法线和亮度三者共同决定：

w(p,q)=wz(p,q)⋅wn(p,q)⋅wl(p,q)w(p, q) = w_z(p, q) \cdot w_n(p, q) \cdot w_l(p, q)

注意这里的权重不是空间权重，而是为了保留边缘考虑的相似度权重。然后考虑的权重项及其实现不唯一

# ① 几何深度项权重（Depth Weight: wzw_z）

防止在 3D 空间中距离很远、但在 2D 屏幕上相邻的两个物体（例如前景的杯子和背景的墙面）互相混色。

wz(p,q)=exp⁡(−∣z(p)−z(q)∣σz⋅∣∇z(p)⋅(p−q)∣+ϵ)w_z(p, q) = \exp \left( - \frac{|z(p) - z(q)|}{\sigma_z \cdot |\nabla z(p) \cdot (p - q)| + \epsilon} \right)

  • z(p),z(q)z(p), z(q)：像素 pp 和 qq 的线性深度值。

  • ∇z(p)\nabla z(p)：像素 pp 处的深度梯度（Depth Gradient），即深度随屏幕坐标的变化率。

  • 分母中的梯度乘项：如果表面是一个斜面，像素间的深度差自然很大。引入梯度项相当于一阶泰勒展开，相当于分母在预测如果 q 仍然位于 p 所在的这个连续表面上，那么从 p 移动到 q，深度理论上应该变化多少

  • σz\sigma_z：控制对深度差异敏感度的超参数。

  • ϵ\epsilon：防止分母为 0 的极小值。

# ② 几何法线项权重（Normal Weight: wnw_n）

防止朝向不同的表面（例如墙角、立方体的两个交面）互相混色，从而保留锐利的几何边缘。

wn(p,q)=exp⁡(−1−n(p)⋅n(q)σn2)=max⁡(0,n(p)⋅n(q))σnw_n(p, q) = \exp \left( - \frac{1 - n(p) \cdot n(q)}{\sigma_n^2} \right) = \max \left( 0, n(p) \cdot n(q) \right)^{\sigma_n}

  • n(p),n(q)n(p), n(q)：像素 pp 和 qq 处的单位法向量。

  • 点积 n(p)⋅n(q)n(p) \cdot n(q)：衡量两个法线的夹角。夹角越大，点积越小，权重衰减极快。

  • σn\sigma_n：控制对法线差异敏感度的超参数（通常设置得较大，使微小的法线变化不会剧烈切断模糊，从而允许平滑的曲面平滑过渡）。

# ③ 亮度项权重 (Luminance Weight)

双边滤波的项

e−∥I(p)−I(q)∥22σr2e^{-\frac{\|I(p) - I(q)\|^2}{2\sigma_r^2}}

# 卷积优化考虑

离散卷积操作的复杂度和卷积半径有关

上面我们知道高斯滤波可以通过拆分维度降低复杂度，但是双边滤波不可以，我们不禁思考还有没有什么可以优化卷积的方法，尤其是大半径的卷积核

# 渐进式增长与多级模糊 (Progressive / Multi-Pass Downsampling)

当需要的模糊半径极其巨大（例如全屏强力 Bloom 特效，核大小可能达到几百像素）时，即便是 2N2N 的一维拆分卷积依然开销巨大。这时会采用降低分辨率 + 多级渐进的策略。

# 核心步骤 (以 Kawase Blur 或 Dual Blur 为例)

  1. 降采样 (Downsample)：首先将原图分辨率降低（如降到 1/2、1/4 或 1/8）。这一步相当于变相将卷积核的感受野扩大了 2、4、8 倍，而像素总量却减少到了 1/4、1/16、1/64。

  2. 小核迭代 (Multi-Pass)：在低分辨率的 RT（Render Texture）上，频繁使用较小的核（如 3x3 或 5x5）进行多次模糊。

  3. 升采样与混合 (Upsample & Blend)：逐级拉伸放大回来，并在过程中与上一层的模糊结果进行线性插值混合。

# 数学原理

对高斯卷积的结果还是高斯，只是方差加和，而低分辨率高斯近似为高分辨率的更高边境高斯

这样一个降采样高斯卷积序列还近似等价高斯，这个用傅里叶变换很容易证明

那么升采样呢？本质上可以从离散信号重建来理解，升采样(插值)本质上是插0+低通滤波，而这里模糊本来就不用考虑高频项

# 优势

  • 成功用极小的计算代价换取了全屏级别的超大模糊感受野。

  • 避免了由于单次大核采样导致的硬件缓存不命中（Cache Miss）问题。

# 双线性采样级联优化 (Bilinear Filtering Trick)

这是一种专属于 GPU 的硬件加速黑科技，常与高斯核拆分配合使用。

# 核心原理

GPU 的纹理采样器（Texture Sampler）内置了双线性插值（Bilinear Interpolation）硬件单元。当我们采样的坐标恰好在两个像素中间时，GPU 会自动帮你把这两个像素的值按距离加权融合后返回，而这一切只消耗一次采样的性能。

也就是降低采样的次数

利用这个特性，我们可以把一维高斯核中连续的两个权重和坐标，融合成一个“虚拟坐标”和“组合权重”：

我们需要找到 W 和 t 使得

w1I1+w2I2=W[(1−t)I1+tI2]w_1I_1+w_2I_2=W[(1-t)I_1+tI_2]

易得

W=w1+w2,t=w2w1+w2W=w_1+w_2,\quad t=\frac{w_2}{w_1+w_2}

  • 假设原本要采样 x1x_1（权重 w1w_1）和 x2x_2（权重 w2w_2）。

  • 我们可以将采样点设在：xoffset=x1⋅w1+x2⋅w2w1+w2x_{\text{offset}} = \frac{x_1 \cdot w_1 + x_2 \cdot w_2}{w_1 + w_2}。

  • 新的加权系数为：wcombined=w1+w2w_{\text{combined}} = w_1 + w_2。

# 复杂度再减半

通过这个黑科技，一个原本需要 NN 次采样的一维卷积，可以被优化为只需要 ⌈N/2⌉\lceil N/2 \rceil 次采样。

  • 配合横纵拆分，一个 15×1515 \times 15 的二维高斯模糊，原本需要 225 次采样。

  • 拆分后需要 15+15=3015 + 15 = 30 次。

  • 引入双线性采样优化后，只需要 8+8=168 + 8 = 16 次采样！基本达到了实时渲染的极致

# 空洞卷积与暂存跳跃采样 (Dilated Convolution / Atrous)

# 核心思想

在增大卷积核半径的时候，让采样数量保持不变，而是增大采样的间距

在采样时，每隔 D−1D-1 个像素采样一次（DD 为扩张率 Dilated Rate）。

# 联级

I→D=1RT1→D=2RT2...I\overset{D=1}{\rightarrow} RT_1\overset{D=2}{\rightarrow} RT_2 ...

RT是render texture中间量

易知，采样的数量随着pass数量是线性增长的，而感受野(Receptive Field，考虑到的值)是指数增长的

# 频域变换优化 (FFT 卷积)

# 核心思想

根据卷积定理，时域/空域的卷积等于频域的乘积。

  1. 使用 快速傅里叶变换 (FFT) 将图像和超级大卷积核都变换到频域。

  2. 将两者的频域矩阵进行逐像素点乘（Hadamard Product）。

  3. 使用 逆快速傅里叶变换 (IFFT) 变回空域。

但Gpu的fft优化很差

# 复杂度对比

  • 空域大核卷积：O(W⋅H⋅N2)O(W \cdot H \cdot N^2)

  • FFT 卷积：O(W⋅Hlog⁡(W⋅H))O(W \cdot H \log(W \cdot H)) —— 其计算开销完全与卷积核的大小 NN 无关！

  • 适用场景：当卷积核极其巨大（如 N>32N > 32 或 6464），且该卷积核无法被空间拆分时（如复杂的离焦相机镜头 Bokeh 效果、复杂的环境光卷积），FFT 算法的优势才会压过它自身的变换开销。

# 考虑时间

我们考虑时间上的信号序列 xt(n)x_t(n)

可以复用时间上前面的信号来降噪

# Spatiotemporal Variance-Guided Filtering

# 方差指导的亮度项权重（Luminance Weight: wlw_l）

直接使用亮度差会导致高噪点区域停止模糊，SVGF的核心是引入历史方差来动态调节模糊力度。

wl(p,q)=exp⁡(−∣L(p)−L(q)∣σl⋅Var(L(p))+ϵ)w_l(p, q) = \exp \left( - \frac{|L(p) - L(q)|}{\sigma_l \cdot \sqrt{\text{Var}(L(p))} + \epsilon} \right)

  • L(p),L(q)L(p), L(q)：像素 pp 和 qq 的亮度（Luminance，通常由 RGB 转化为标量）。

  • Var(L(p))\text{Var}(L(p))：像素 pp 处的亮度方差（Variance），历史中亮度的方差，通过时间重投影和局部空间累积得到。

  • 方差作为分母的作用：

    • 当方差很大时（噪点极多）：分母变大，整个指数项趋近于 0，使得 wl→1w_l \to 1。这意味着算法会忽略亮度差异，强行进行大范围模糊以消除严重的噪点。

    • 当方差很小时（画面已收敛/稳定）：分母变小，亮度项变得非常敏感。只有亮度极其接近的像素才会互相模糊，从而完美保留了高频的阴影边缘和纹理细节。

滤波算法常考虑 À-Trous 算法，也就是上面说的空洞卷积

# 时间重投影与方差累积（数据源头）

在执行上述空间滤波之前，画面必须通过时间轴获取足够的数据。
首先进行的是时间重投影，利用motion vector、光流法之类的找到当前像素在上一帧的位置

# 时间累积公式（指数移动平均 - EMA）

Cˉt=α⋅Ct+(1−α)⋅Cˉt−1\bar{C}_t = \alpha \cdot C_t + (1 - \alpha) \cdot \bar{C}_{t-1}

我们实际上只需要维护关于时间的像素亮度序列的一阶矩和二阶矩

  • α\alpha（融合权重）：通常根据历史帧数 NN 动态计算：α=max⁡(1/N,0.05)\alpha = \max(1/N, 0.05)。历史帧越多，α\alpha 越小，越信任历史数据。

  • 当发生残影（历史失效）时：如果深度或法线不匹配（检测到遮挡或光源瞬变），历史帧数 NN 直接重置为 1（即 α=1\alpha = 1），此时完全断开时间累积，改由当帧的空间方差 Varspatial\text{Var}_{\text{spatial}} 接管。

因此考虑时间的降噪算法更重要的是如何决定历史有多么可信

# Adaptive Spatiotemporal Variance-Guided Filtering

对于SVGF，α\alpha小，历史权重大，噪声小但容易出现ghost；反之亦然

因此问题就在于：不同像素、不同时间的最优 α 根本不一样

比如我们考虑考虑历史长度，历史久远的像素α\alpha小，刚出现的α\alpha大

# Kalman Filter

Kalman Filter（卡尔曼滤波）不是一个空间卷积核，而是一种**随时间递归地融合“预测”和“观测”**的方法。它要解决的问题可以概括为：

真实状态看不见，只能看到带噪声的观测；如何利用过去的估计和当前观测，得到当前状态的更好估计？

例如在 Path Tracing 中，某个表面点的真实 radiance 看不见，只能得到当前帧的 Monte Carlo 采样结果。当前采样很 noisy，但上一帧的结果也可能因为物体运动、遮挡或光照变化而过时。Kalman Filter 用一个随时间变化的权重，自动决定两者各该相信多少。

# 先建立直觉：它一直维护两样东西

对每个时刻 tt，滤波器维护：

  • 状态估计 x^\hat{x}：目前认为真实状态是多少；
  • 协方差 PP：对这个估计有多不确定。

这里的 PP 不是“实际误差”，因为实际误差本来就是未知的；它是滤波器根据模型和历史观测推断出的误差不确定性。标量状态时它就是方差，向量状态时它是协方差矩阵，还能描述不同分量之间的相关性。

每一帧都做两件事：

  1. 用上一帧的估计预测这一帧的状态，同时传播不确定性；
  2. 用当前观测纠正预测，并重新计算纠正后的不确定性。

因此它不是简单地“永远取历史”或“永远取当前”，而是：预测越不确定，就越依赖观测；观测越不可靠，就越依赖预测。

# 适用条件与状态空间模型

经典 Kalman Filter 假设状态转移和观测都是线性的，并且噪声近似为零均值高斯噪声：

# 状态转移模型

xt=Atxt−1+Btut+wt,wt∼N(0,Qt)x_t = A_t x_{t-1} + B_t u_t + w_t, \qquad w_t \sim \mathcal{N}(0,Q_t)

它描述“真实状态如何从上一时刻变化到当前时刻”。

  • xtx_t：时刻 tt 的真实状态，可以是位置、速度、颜色、radiance 等；
  • AtA_t：状态转移矩阵；
  • utu_t：已知的外部输入，不一定是人为输入，例如控制信号或已知的运动量；
  • BtB_t：把输入映射到状态空间的矩阵；
  • wtw_t：过程噪声，表示状态自身的随机变化，以及模型没有描述好的部分；
  • QtQ_t：过程噪声的协方差，表示我们对状态转移模型有多不放心。

没有控制输入时通常写成：

xt=Atxt−1+wtx_t = A_t x_{t-1} + w_t

# 观测模型

zt=Htxt+vt,vt∼N(0,Rt)z_t = H_t x_t + v_t, \qquad v_t \sim \mathcal{N}(0,R_t)

它描述“真实状态经过观测系统后，变成了什么样的带噪观测”。

  • ztz_t：当前实际拿到的观测；
  • HtH_t：把状态空间映射到观测空间的矩阵；
  • vtv_t：观测噪声；
  • RtR_t：观测噪声的协方差，表示当前观测有多不可靠。

当噪声确实是高斯噪声时，Kalman Filter 给出的后验均值是均方误差意义下的最优估计，同时也是后验分布的均值和 MAP 估计。如果噪声不是高斯分布，公式仍然可以作为线性最小均方估计，但不一定是完整贝叶斯意义下的最优解。

# 在图形学中的对应关系

以时间降噪为例，必须先通过 motion vector、深度和法线等信息，把当前像素重投影到上一帧中对应的表面点。重投影成功后，常见的简化是：

xt=xt−1+wt,zt=xt+vtx_t = x_{t-1} + w_t, \qquad z_t = x_t + v_t

也就是 At=IA_t=I、Ht=IH_t=I。

  • xtx_t：同一个表面点在当前时刻的真实 radiance；
  • ztz_t：当前帧的路径追踪采样结果；
  • RtR_t：当前采样的噪声方差，可以根据采样数、样本方差等估计。如果 ztz_t 是 NN 个样本的平均值，通常应使用类似 s2/Ns^2/N 的估计，而不是直接使用单个样本的方差 s2s^2；
  • QtQ_t：radiance 发生变化、重投影不完全准确、模型遗漏等带来的不确定性。

注意：同一个屏幕坐标不一定对应同一个世界表面点。如果不先重投影，就把当前像素和上一帧的同位置像素直接融合，很容易产生拖影和 ghosting。

# 两个阶段：先预测，再更新

记号 x^t∣t−1\hat{x}_{t\mid t-1} 表示“只使用 t−1t-1 及以前的信息，对 tt 时刻的状态做出的估计”；x^t∣t\hat{x}_{t\mid t} 表示“已经使用当前观测 ztz_t 后的估计”。竖线左边是要估计的时刻，右边是已经使用的信息截止时刻。

# 1. 预测阶段（Predict / Time Update）

使用上一时刻的后验估计预测当前状态：

x^t∣t−1=Atx^t−1∣t−1+Btut\hat{x}_{t\mid t-1} = A_t\hat{x}_{t-1\mid t-1} + B_tu_t

同时传播不确定性：

Pt∣t−1=AtPt−1∣t−1AtT+QtP_{t\mid t-1} = A_tP_{t-1\mid t-1}A_t^T + Q_t

这里的 Pt∣t−1P_{t\mid t-1} 叫先验协方差。即使上一帧很确定，只要 QtQ_t 不为零，经过状态转移后当前预测也会变得更不确定；这就是公式中要加上 QtQ_t 的原因。

# 2. 更新阶段（Update / Measurement Update）

现在把当前观测 ztz_t 拿来修正先验预测。

首先计算innovation，也就是实际观测和预测观测之间的差：

rt=zt−Htx^t∣t−1r_t = z_t - H_t\hat{x}_{t\mid t-1}

innovation的协方差为：

St=HtPt∣t−1HtT+RtS_t = H_tP_{t\mid t-1}H_t^T + R_t

它同时包含两部分不确定性：预测本身的不确定性，以及观测噪声。然后计算 Kalman 增益：

Kt=Pt∣t−1HtTSt−1K_t = P_{t\mid t-1}H_t^TS_t^{-1}

Kalman 增益就是“应该沿着创新修正多少”的权重。用它更新状态：

x^t∣t=x^t∣t−1+Ktrt\hat{x}_{t\mid t} = \hat{x}_{t\mid t-1} + K_tr_t

展开就是：

x^t∣t=x^t∣t−1+Kt(zt−Htx^t∣t−1)\hat{x}_{t\mid t} = \hat{x}_{t\mid t-1} + K_t\left(z_t-H_t\hat{x}_{t\mid t-1}\right)

最后更新后验协方差：

Pt∣t=(I−KtHt)Pt∣t−1P_{t\mid t} = (I-K_tH_t)P_{t\mid t-1}

这个简洁形式在数学推导中很常见。实际浮点实现中，为了更好地保持协方差矩阵的对称性和半正定性，通常使用等价的 Joseph 形式：

Pt∣t=(I−KtHt)Pt∣t−1(I−KtHt)T+KtRtKtTP_{t\mid t} = (I-K_tH_t)P_{t\mid t-1}(I-K_tH_t)^T + K_tR_tK_t^T

整个流程可以压缩成：

上一帧后验→A,Q当前先验→z,H,R当前后验\text{上一帧后验}\xrightarrow{A,Q}\text{当前先验}\xrightarrow{z,H,R}\text{当前后验}

# 最重要的特例：标量、直接观测

在图像降噪中，先看最容易理解的 A=1A=1、H=1H=1 情况。此时：

Pt−=Pt−1+QP_t^- = P_{t-1}+Q

Kt=Pt−Pt−+RK_t = \frac{P_t^-}{P_t^-+R}

x^t=x^t−+Kt(zt−x^t−)=(1−Kt)x^t−+Ktzt\hat{x}_t = \hat{x}_t^- + K_t(z_t-\hat{x}_t^-) = (1-K_t)\hat{x}_t^-+K_tz_t

Pt=(1−Kt)Pt−P_t = (1-K_t)P_t^-

这几行公式已经揭示了 Kalman Filter 的核心：它就是一个权重会根据不确定性自动变化的加权平均。

  • 若 Pt−≫RP_t^-\gg R，预测很不可靠，Kt≈1K_t\approx 1，更多相信当前观测；
  • 若 R≫Pt−R\gg P_t^-，当前观测很 noisy，Kt≈0K_t\approx 0，更多相信历史预测；
  • 增大 QQ 会增大先验不确定性，从而增大 KtK_t，让滤波器更快响应变化；
  • 增大 RR 会减小 KtK_t，让滤波器更平滑，但也更容易滞后或产生拖影。

为什么是“按方差的倒数”加权？在 H=1H=1 的情况下，后验估计也可以看作最小化：

(x−x^t−)2Pt−+(zt−x)2R\frac{(x-\hat{x}_t^-)^2}{P_t^-}+\frac{(z_t-x)^2}{R}

不确定性越小，分母越小，对应的惩罚越大，因此越应该相信它。矩阵情况下同样的思想写成：

min⁡x(x−x^t−)T(Pt−)−1(x−x^t−)+(zt−Htx)TRt−1(zt−Htx)\min_x\;(x-\hat{x}_t^-)^T(P_t^-)^{-1}(x-\hat{x}_t^-) +(z_t-H_tx)^TR_t^{-1}(z_t-H_tx)

# 一个数值例子

假设预测值为 x^t−=10\hat{x}_t^-=10，先验方差为 Pt−=4P_t^-=4；当前观测为 zt=14z_t=14，观测噪声方差为 R=9R=9。则：

Kt=44+9=413≈0.308K_t=\frac{4}{4+9}=\frac{4}{13}\approx0.308

所以：

x^t=10+0.308(14−10)≈11.23\hat{x}_t=10+0.308(14-10)\approx11.23

结果不是 1010 和 1414 的普通平均值 1212，因为当前观测的方差 99 大于预测的方差 44，滤波器认为观测不够可靠。更新后：

Pt=(1−0.308)×4≈2.77P_t=(1-0.308)\times4\approx2.77

融合观测后，不确定性从 44 降到了约 2.772.77。

# 初始化与历史失效

第一帧没有历史信息，因此不能凭空得到 x^0∣−1\hat{x}_{0\mid-1}。一种简单的初始化方式是直接使用当前观测：

x^0∣0=z0,P0∣0=R0\hat{x}_{0\mid0}=z_0,\qquad P_{0\mid0}=R_0

也可以根据实际情况给 P0∣0P_{0\mid0} 一个更保守的较大值。之后滤波器会通过连续的预测和更新逐渐修正初始值。

如果重投影失败，或者深度、法线、材质等检查表明历史已经不再对应当前表面点，就应该丢弃旧的 x^\hat{x} 和 PP，重新用当前观测初始化。不要只把旧颜色清掉而继续沿用旧协方差，因为那会让滤波器错误地认为“当前历史仍然很可靠”。

# 它和 EMA、SVGF 的关系

在 A=H=1A=H=1 时，Kalman 更新式为：

x^t=(1−Kt)x^t−1+Ktzt\hat{x}_t=(1-K_t)\hat{x}_{t-1}+K_tz_t

这和指数移动平均（EMA）完全同形：

xˉt=(1−α)xˉt−1+αzt\bar{x}_t=(1-\alpha)\bar{x}_{t-1}+\alpha z_t

所以可以把 Kalman 增益 KtK_t 看成自适应的 EMA 系数 α\alpha。区别在于：

  • EMA 由我们手动指定固定或启发式的 α\alpha；
  • Kalman Filter 根据 PP、QQ、RR 推导出当前应该使用的 KtK_t；
  • QQ 和 RR 变化时，KtK_t 也会随像素、随时间变化。

放到 Path Tracing 的时域降噪中，可以这样理解：

  • 当前帧样本方差大，增大 RtR_t，降低 KtK_t，更多使用历史；
  • 预期 radiance 变化大，增大 QtQ_t，提高 KtK_t，更快接受当前帧；
  • 深度、法线或运动矢量检查失败时，历史不再代表同一个表面点，应当丢弃或重置历史，而不是继续套用旧的 PP；
  • Kalman Filter 主要负责时间上的融合，SVGF 或双边滤波仍可以负责空间上的边缘感知平滑。

这里要区分两种方差：当前样本的统计方差通常用来估计观测噪声 RtR_t；PtP_t 则是滤波器对“融合后的状态估计”仍有多不确定的判断，它们不是同一个东西。

# 工程实现时的伪代码

// x_prev, P_prev: 上一帧后验估计和协方差
// z: 当前观测；A, H, Q, R: 当前时刻的模型和噪声参数

x_pred = A * x_prev + B * u
P_pred = A * P_prev * transpose(A) + Q

innovation = z - H * x_pred
S = H * P_pred * transpose(H) + R
K = P_pred * transpose(H) * inverse(S)

x = x_pred + K * innovation
P = (I - K * H) * P_pred

实际代码中不建议显式计算 inverse(S)，而应使用线性方程求解器；如果数值稳定性比较重要，协方差更新使用 Joseph 形式。对于 RGB 状态，可以把 xx 设为三维向量，PP 可以是每个通道独立的对角矩阵，也可以保留通道之间的相关性。

# 常见误区

  1. 把 PP 当成真实误差。 PP 是模型对误差的预测，模型参数错了时它也会错。
  2. 把 QQ 和 RR 混为一谈。 QQ 描述状态模型的不确定性，RR 描述当前观测的不确定性；前者发生在状态空间，后者发生在观测空间。
  3. 令 Q=0Q=0 却希望滤波器快速跟随变化。 Q=0Q=0 等价于认为状态永远严格按照模型变化，连续融合后 PP 会越来越小，KK 也会变小，最终几乎不再相信新观测。
  4. 不做重投影就复用历史像素。 屏幕位置相同不代表世界中的表面点相同，这是产生 ghosting 的常见原因。
  5. 认为 Kalman Filter 能解决所有噪声。 它对线性、高斯、近似独立的噪声最合适；遇到明显非线性可考虑 EKF/UKF，遇到离群值则需要创新值裁剪、统计门限或鲁棒滤波。

待续...

慢慢补信号的知识了...