卷 IV · 相乘CH 16深度 16/20

从噪声底下把东西捞出来

前十五章处理的都是「信号在那儿，怎么看清楚」。这一章处理「信号根本看不见」——它比噪声弱一百倍，示波器上是雪花，频谱仪上什么都没有。而只要你手里有它的样板，一次相关就够了。你的手机此刻正在用这一招接收 GPS。

★★ −20 dB 也能捞7.1 σ512 ms 压成 0.375 ms
▷ 先猜一下

一段 512 毫秒的「啁啾」——频率从 200 Hz 一路扫到 3000 Hz 的那种滑音。把它埋进白噪声里，信噪比 −20 dB（信号的有效值只有噪声的百分之一）。

直接看波形：信号的峰值只有噪声标准差的 0.141 倍——完全淹没，示波器上就是一片雪花。

现在我手里有一模一样的样板（我知道那段啁啾长什么样，只是不知道它在什么时候出现）。我把样板滑过整段录音，逐点相乘再求和。

问：在信号真正出现的那个位置，相关值有多高？（用「相当于噪声本底的多少个标准差」来衡量。）

A 还是看不见。−20 dB 就是 −20 dB，相关也变不出信息 B 大约 2 σ。勉强能看出个鼓包，但达不到可靠检测的程度 C 7.1 σ，而整段噪声里最高的假峰只有 4.10 σ——清清楚楚 D 100 σ。既然信噪比差 100 倍，相关就该补回 100 倍

相关：不翻转的卷积

动作和第 13 章几乎一样，只差一个翻转：

卷积：  (a ✻ b)[k] = Σ a[i] · b[k − i]      ← b 倒着走
相关：  (a ⊛ b)[k] = Σ a[k + i] · b[i]      ← b 正着走

念出来：**把样板挪到位置 k，逐点相乘，求和。**
样板和那一段越像，这个和越大。

为什么「越像和越大」？因为这就是内积——第 3 章那个东西。内积衡量两个东西有多「同方向」，样板对齐的时候，正对正、负对负，每一项都是正的，加起来很大；错开的时候，正负乱配，互相抵消。

顺带一句，深度学习框架里的 Conv2d 实现的正是这一个（第 13 章的「术语正名」）。CNN 在做的字面意义上的事情，是拿一堆学出来的样板在图上到处找匹配。这一章讲的就是它为什么有效。

把那组数摆出来

啁啾 4096 点（512 ms @ 8 kHz），200 → 3000 Hz
埋进白噪声，输入信噪比 −20.0 dB，藏在第 12000 个采样处

  相关之前
    信号峰值 / 噪声 σ  =  0.141        ✗ 完全看不见

  相关之后
    峰的位置            =  12000        ✓ 一个样本都不差
    峰高 / 噪声本底 σ   =  7.1 σ
    整段噪声里最高的假峰 =  4.10 σ      ← 真峰高出一截，判得干脆

  实测增益  34.0 dB        理论处理增益  10·log₁₀(4096) = 36.1 dB
◆ 主线

只要你知道信号长什么样，就可以用「相关」把它从任意深的噪声里捞出来——代价是时间。

处理增益 = 10 · log₁₀(N)        N 是样板的长度（点数）

  N = 100      → 20 dB
  N = 4096     → 36 dB
  N = 1 000000 → 60 dB

道理很朴素：信号是相干的，噪声不是。把 N 个样本对齐相加，信号的幅度按 N 增长（每一项都同号），噪声的幅度只按 √N 增长（随机游走）。信噪比因此按 √N 改善，也就是每 10 倍长度换 10 dB。

这也说明了它的极限在哪：想再多 10 dB，样板得长十倍。GPS 冷启动要几十秒，就是在付这笔账。

为什么实测 34.0 dB 比理论的 36.1 dB 少了 2 dB？因为理论值假设噪声是理想白噪声且样板与信号完全对齐；实测里噪声是有限样本的（4.10 σ 那个最高假峰本身就是有限样本的产物），而且噪声本底的估计也有涨落。2 dB 的差距在这类实验里属于正常，重要的是数量级对得上。

为什么偏偏是「啁啾」

这里还藏着一个漂亮的设计。为什么雷达、声呐、LIGO 都用这种频率一路上扫的信号，而不是一个短促的脉冲？

因为它们要同时要两样互相打架的东西：

  • 时间分辨率要高（能分开两个挨得近的目标）→ 第 10 章：需要大带宽。
  • 能量要足（能打得远）→ 需要长时间发射，因为发射机的峰值功率有上限。

一个短脉冲带宽大但能量小；一个长的单频音能量大但带宽小。啁啾两样都要：它在时间上很长（能量足），同时扫过很宽的频率范围（带宽大）。而相关这一步会把它「压」回去：

啁啾长度        512.0 ms
相关峰的包络半高宽   0.375 ms          ← 压缩了 1365 倍
理论宽度 1/B     0.357 ms （B = 2800 Hz）

时间带宽积 BT = 2800 × 0.512 = 1434     ← 压缩比 ≈ BT

这叫脉冲压缩：发射时把能量摊在半秒里，接收时用相关把它挤回三分之一毫秒的一根针。你既拿到了长信号的能量，又拿到了短脉冲的分辨率——第 10 章那条不等式没有被违反，因为那根针的带宽确实是 2800 Hz。

这也顺带回答了「为什么峰这么窄」：相关峰的宽度由带宽决定（1/B），和信号本身有多长完全无关。

✎ 术语正名

匹配滤波（matched filter）和这一章的「相关」是同一件事的两种说法。

把相关写成卷积：a ⊛ b 等于 a ✻ b̃，其中 b̃ 是 b 的时间反转（复信号还要取共轭）。所以「拿样板去相关」就是「用一个冲激响应等于样板时间反转的滤波器去滤」——那个滤波器就叫匹配滤波器。

它有一条很强的最优性定理：在白噪声下，所有线性滤波器里，匹配滤波器让输出峰值信噪比最大，没有之一。注意那三个前提：线性、白噪声、只关心峰值信噪比。噪声不白的时候要先做「白化」（这正是 LIGO 数据处理的第一步）；关心别的指标（比如同时要压低旁瓣）时，最优的就不是它了。

还有一个实现上的关键：相关也可以走频域。按第 14 章的定理，相关 = 频域里共轭相乘再逆变换。GPS 接收机搜索卫星要在几千个「码相位 × 多普勒频移」组合里找峰，全靠这条才算得动。

⌨ 自己跑一遍

把整个实验重跑一遍，二十行：

import numpy as np

fs, N, L, at = 8000, 4096, 32768, 12000
t = np.arange(N) / fs
T = N / fs
s = np.cos(2*np.pi*(200*t + 2800*t*t/(2*T)))          # 200 → 3000 Hz 啁啾

rng = np.random.default_rng(424242)
sigma = s.std() * 10                                   # SNR = −20 dB
rec = rng.standard_normal(L) * sigma
rec[at:at+N] += s

print('相关前  信号峰 / 噪声σ = %.3f' % (s.max()/sigma))     # 0.141

c = np.correlate(rec, s, 'valid')                      # 样板滑过去
pk = int(np.argmax(np.abs(c)))
bg = c[np.abs(np.arange(len(c)) - at) > 64].std()      # 排除峰附近估本底
print('峰在', pk, ' 峰高 %.1f σ' % (abs(c[pk])/bg))     # 12000  7.1 σ
print('噪声里最高假峰 %.2f σ' % (np.abs(c[np.abs(np.arange(len(c))-at)>64]).max()/bg))

把 sigma 的系数从 10 改成 30（−30 dB），峰会掉到 2σ 左右，检测失败；这时把 N 从 4096 改成 40960（样板长十倍），它又回来了。这一步能让你亲手把「处理增益 = 10·log₁₀(N)」摸一遍。

走频域的版本（大 N 时快几个数量级）：

M = 1
while M < L + N: M <<= 1
c2 = np.real(np.fft.ifft(np.fft.fft(rec, M) * np.conj(np.fft.fft(s, M))))
print('频域路线的峰位置：', int(np.argmax(c2)))          # 12000

python3 -c "import numpy as np;print('%.1f dB' % (10*np.log10(4096)))"

在线跑：Google Colab。想看真实数据：LIGO 开放数据中心（gwosc.org）提供 GW150914 的原始应变数据和官方教程，那个教程做的就是这一章的事，只是多了一步噪声白化。

▸ 在现实里
  • GPS。卫星在两万公里外，功率二十几瓦，到你手机时信号强度大约比热噪声底低 20–30 dB——用任何频谱仪去看，那个频段上什么都没有，就是噪声。手机能收到，全靠这一章：每颗卫星有一段 1023 位的伪随机码，接收机拿本地生成的同一段码去相关，处理增益 30 dB 把它拉出水面。「用一段谁都知道的码去换 30 dB」是 GPS 能存在的全部原因。
  • LIGO 与引力波。2015 年 9 月 14 日那次并合，探测器的臂长变化只有质子直径的千分之一。做法是准备了几十万个不同质量组合的模板（波形都是广义相对论算出来的啁啾），拿每一个去和数据相关。这一章的方法直接换来了一次诺贝尔奖。顺带说，那个波形是啁啾不是巧合：双星越靠近转得越快，频率自然一路上扫。
  • 雷达与超声。脉冲压缩是现代雷达的标配。医用超声也一样——探头发出的是编码的长脉冲，接收后压缩，这样既保证穿透深度又保证分辨率。
  • 听歌识曲。Shazam 不做全波形相关（太贵），但它的哈希匹配在精神上是同一件事：拿一组已知的特征去和数据库对齐，靠大量弱证据的一致对齐换取判决可信度。
  • 生物信息学。DNA 序列比对里的 Smith–Waterman 与这一章同源；而基因组序列的快速比对确实有基于 FFT 的实现——把碱基编码成复数，相关一遍，找峰。
✗ 这个直觉是错的

「信噪比小于 1 就没救了。信号比噪声还弱，说明它的信息已经被淹没了。」

「信噪比」这个数没有说清楚是在什么带宽、什么时长上测的，而这一章的全部内容就是：换一个时长去测，它就变了。

上面那个信号在单个采样点上的信噪比是 −20 dB；在整段 4096 点上的信噪比是 −20 + 36.1 = +16 dB。同一个信号，同一份数据，两个数差了 36 dB——差别只在于你把多长的一段当成「一次观测」。

这里被淹没的从来不是信息，是某一种看法。信息一直在，只是分散在四千个样本里，每个样本里只有一点点。相关这个动作把它们相干地加了起来。

当然有前提，而且很硬：你必须知道信号长什么样。这一章的一切建立在「手里有样板」上。没有样板的话（比如你要找的是一个未知形状的信号），处理增益就没了——这时只能去用别的东西：能量检测、周期图、或者干脆去找一族可能的模板挨个试（LIGO 就是这么干的，代价是几十万次相关）。判据是：先问「我知道我在找什么吗」——知道，就还有 10·log₁₀(N) 可以拿。

◇ 揭晓

正确答案是 C：7.1 σ，而整段噪声里最高的假峰只有 4.10 σ。峰的位置一个样本都不差。

A 「相关也变不出信息」——相关确实变不出信息，它做的是把已有的信息重新聚拢。信号的能量原本摊在 4096 个样本里，每个样本里都被噪声压着；相关把这 4096 份微弱的证据相干地加了起来，而噪声只能非相干地加。这个「相干 vs 非相干」的差别就是全部收益。 B 「大约 2 σ」——量级差了一档。可以粗算一下：输入信噪比 −20 dB，处理增益 36.1 dB，输出约 +16 dB，也就是幅度上大约 6.3 倍。实测 7.1 σ，对得上。这类问题养成先算一遍 10·log₁₀(N) 的习惯，比凭感觉猜可靠得多。 D 「补回 100 倍」——混淆了「增益」和「信噪比」。处理增益是 36 dB（约 63 倍幅度），不是恰好抵消那 100 倍，也不是无限。而且要注意信噪比不是「峰高几个 σ」：后者还要看你比较的对象（这里是相关输出的本底涨落）。这两个数在实践中经常被混着报，看到一个 dB 数先问它的分母是什么。
⟳ 换个域看
时域里，相关是「拿样板滑过去，每停一次算一次内积」——一个直观但昂贵的动作，L 个位置 × N 点样板。 频域里，它是逐点共轭相乘：C[k] = X[k] · conj(S[k])。而这句话还有一个漂亮的读法——共轭就是把相位取反，也就是「把样板在每个频率上的相位抵消掉」。

于是匹配滤波在频域里的意思变得极其清楚：它把信号里所有频率成分的相位对齐到同一个点上，让它们在那一刻同时到达波峰、叠出一根针；而噪声的相位是随机的，对齐不了，只能继续摊平。

这也一口气解释了脉冲压缩：啁啾的各个频率成分本来是先后到达的（低频先、高频后），匹配滤波给每个频率补上一个相反的延迟，把它们全部拉到同一时刻。512 毫秒压成 0.375 毫秒，用的就是第 4 章那半张脸——相位。

这一章的一句话

只要你知道要找的东西长什么样，就能用相关把它从比它强一百倍的噪声里捞出来，代价是 10·log₁₀(N) 这条明码标价的兑换率；而做到这件事的机制，是把它所有频率分量的相位重新对齐到同一刻。

卷 IV 到此结束。你现在有了一条定理（卷积 = 相乘）和它的两个大用场（滤波、检测）。

但这一切有个前提：你得算得动。按定义式算一次 N 点的变换要 N² 次运算，本机实测 2²⁰ 个点要跑 6.2 小时——那样的话，前面十六章全是纸上谈兵，没有一样能进你的耳机。下一卷从那个把 6.2 小时变成 29 毫秒的算法开始。