做数据压缩作业时，我把“完全重建QAM滤波器组”这组概念翻来覆去折腾了两周。刚开始顶着“滤波器组不是通信课的东西吗，跟压缩有什么关系”的疑问硬啃，等到代码跑通、频谱图一张张摆出来，才真正明白子带编码的底气全部系在“重建质量”四个字上。这篇文章就是完整的过程复盘：从子带编码为什么要拆信号，到QAM滤波器组的调制结构，再到完全重建的数学条件、工程实现与踩坑排查，最后串到压缩系统里的比特分配。内容按一条可复现的路径来写，适合正在学多率信号处理、研究子带编码，或者想搞懂音频和图像压缩底层原理的读者。代码以Python为基础，思路不依赖特定平台，换个语言也能平移过去。
1. 压缩为什么要先拆信号：子带编码的基本逻辑
1.1 不是所有频带都值得用同样的码率去编码
直接对采样点做PCM编码，等于默认每个样本的重要性完全一样，每个时刻的频率成分都要花同样多的比特去表示。但真实信号几乎从不长这样：语音的主要能量集中在低频，图像的细节集中在中高频但占比很小，一段音乐里可能某个频段噪声很大而另一个频段几乎空白。把信号在频域切开，对每个子带单独处理，就能把比特花在真正需要的地方。
用一个最简单的公式说明收益：若整段信号用R比特/样本编码，误差功率约为σ²·2^(-2R)。如果把信号分成M个子带，给第k个子带分配R_k比特、并保持ΣR_k = MR，那么整体误差是各子带误差之和。由于子带方差σ_k²差异往往很大，让σ_k²大的子带多拿比特、小的少拿甚至不拿比特，总误差会明显小于直接PCM。这个“能量集中带来压缩增益”的现象，就是子带编码的第一推动力。
1.2 拆得开还得合得上：分解环节不能提前丢信息
滤波器组是把信号拆开再拼起来的“切割机”。问题在于：如果切割本身就已经让信息丢失，后面量化做得再精细也救不回来。很多新手做子带编码实验时，把所有注意力放在量化和熵编码上，结果先被分析-综合级联本身的失真吓了一跳——明明没有量化，输入和输出就已经对不上了。
这里要区分两个概念：混叠失真和量化噪声。混叠失真来源于下采样时的频谱搬移，如果综合端没有做对应的对消，高频会被折叠到低频，造成类似“镜像鬼影”的听觉或视觉伪影；量化噪声则是压缩时比特不足导致的误差。完全重建滤波器组解决的是前者：让“拆-合”这个过程本身不引入任何失真，相当于把干净信号交给量化环节。如果这一步不满足，后续所有压缩优化都是在错误的地基上盖楼。
1.3 用什么指标衡量重建质量
工程上最常用的指标是重建信噪比（SNR）：
[ SNR = 10\log_{10}\frac{\sum x[n]^2}{\sum (x[n-d]-\hat{x}[n])^2} ]
d是系统固有的群延迟。注意这里的延迟必须对齐，否则信噪比会被人为拉低。除了SNR，还要看误差频谱是否均匀：完全重建系统在无量化时误差频谱应该是平坦的或接近数值零；一旦出现某个频段的尖峰，基本可以断定是混叠对消不彻底。对于音频，人耳对某些频段的误差更敏感，有人会用感知加权或掩蔽阈值来修正指标，但那属于心理声学范畴，先把基础SNR做扎实再谈。
2. QAM滤波器组在数学上到底做了什么
2.1 一个低通原型，复制出M路带通滤波器
QAM滤波器组的核心思想非常朴素：设计一个低通原型滤波器p[n]，然后把它调制到M个不同的中心频率上，就得到M个并行的带通滤波器。调制公式是：
[ h_k[n] = p[n]\cdot e^{j\frac{2\pi}{M}kn}, \quad k=0,1,\dots,M-1 ]
这里e的指数项本质上是一个复正弦序列，它把p[n]的频谱从零频搬到2πk/M处。p[n]的带宽决定每个子带占用多宽：如果p[n]是理想低通、截止频率为π/M，M个通道拼接起来正好覆盖[-π,π]的完整数字频率范围。实际中p[n]只能做到近似，过渡带会带来子带间的重叠，这也是后面所有混叠问题的来源。
这种做法的巧妙之处在于：M个子带滤波器不需要分别设计，只要设计好一个p[n]，其余的全部由调制生成，设计自由度一下子压缩到一个滤波器上。这也解释了为什么“原型滤波器的质量决定整个滤波器组的质量”——后续的完全重建条件，本质上都是对p[n]提要求。
2.2 “QAM”这个名字是怎么回事
第一次看到“QAM滤波器组”的人容易想起通信里的正交幅度调制，两者确实有血缘关系。QAM调制把符号分成同相分量和虚部正交分量，两个分量相位差90度叠加到载波上，从而在一个频带内传输两路信息。QAM滤波器组则把这种思想用在滤波器构造上：用复指数调制生成M个复系数滤波器后，每个子带的输出都是复信号，可以拆成实部I和虚部Q两路处理。
为什么要绕这么一圈用复信号而不是直接做M个实的带通滤波器？因为一个实带通滤波器对应正负两个频率区间，只覆盖半边频带，要铺满[-π,π]需要2M个实滤波器；而复调制滤波器只占据单边频带，M个通道就够了。在数字实现时，复子带被称为I/Q两路信号，虽然每路仍是实数，但两路共享同一组多相计算和DFT变换，结构和效率都更规整。所以“QAM滤波器组”这个叫法的核心含义是：用类似QAM调制的方式构造滤波器组，换来的是频谱利用率和实现上的对称性。
2.3 和余弦调制、MDCT、DFT滤波器组的关系
如果把复指数换成cos(2πkn/M)，就得到余弦调制滤波器组（CMFB）。它是实数滤波器组，不产生复子带，临界采样时可以用更严格的数学条件实现完全重建。我们熟悉的MDCT（改进离散余弦变换）本质上是一种余弦调制滤波器组，配合时间域混叠消除（TDAC）就能做到完整重建，所以它长期霸占音频编码的舞台。DFT滤波器组则是复指数调制的直接版本，也就是最朴素的QAM滤波器组。
把它们放在一起看就清楚了：MDCT/CMFB是实数分支，DFT/QAM滤波器组是复数分支。理论研究时经常先设计一个余弦调制完全重建系统，再通过解析延拓或并联I/Q两个实数滤波器组的方式构造复数版本，所以QAM滤波器组可以理解为余弦调制滤波器组的复扩展。理解了这层关系，看文献时就不会被一堆名称绕晕：它们都是原型低通滤波器加正弦/余弦/指数调制这个大框架下的不同实例。
2.4 关键参数一览
| 参数 | 含义 | 典型取值 |
|---|---|---|
| M | 子带数量 | 8、16、32、64 |
| p[n] | 原型低通滤波器 | FIR，长度N |
| L = N/M | 每个多相分量的长度 | 4~8 |
| ω_k | 第k个子带的中心频率 | 2πk/M |
| 子带带宽 | 约等于p[n]的带宽 | 接近2π/M |
| 输出子带采样率 | 降低M倍 | f_s/M |
N必须是M的整数倍，否则多相分解无法均匀切分，完全重建条件会变得非常难满足。这是设计时最容易忽略的约束，我在第一次实现时随意选了个N=50，结果多相分支长短不一，后面对齐一直出问题。后来老老实实N=ML，一切顺利很多。
3. 完全重建的硬约束：从输入输出方程到三个条件
3.1 分析-综合级联的长相
把分析滤波器组和综合滤波器组串起来，输入x[n]先经过分析、下采样M倍、再上采样M倍、经过综合滤波器，最后叠加。在z域写出来是：
[ \hat{X}(z) = \sum_{l=0}^{M-1} A_l(z) X(zW^l), \quad W=e^{-j2\pi/M} ]
其中
[ A_l(z) = \frac{1}{M}\sum_{k=0}^{M-1} H_k(zW^l) G_k(z) ]
这里X(zW^l)不是新信号，而是x[n]乘上一个复指数后的z变换，代表频谱搬移。l=0时是原始频谱本身，l≠0时就是下采样造成的镜像频谱。整个公式直观含义是：输出信号等于原始信号加M-1个搬移副本，各自乘以一个系数A_l(z)。完全重建的目标就是让搬移副本全部消失，只留下l=0项。
3.2 三个条件：混叠对消、幅度全通、相位线性
把目标拆开，完全重建要求三条同时成立：
第一，对所有l≠0，A_l(z)=0。这叫混叠对消条件。它保证任何频谱搬移产生的镜像分量在综合端被清除。
第二，A_0(z)的幅频响应是全通的，即|A_0(e^{jω})|是常数。这保证信号各个频率成分不会被滤波器组“夹生”——有的频段被放大，有的被衰减。
第三，A_0(z)的相位是线性的，也就是群延迟恒定。否则不同频率的时延不一致，波形会散开，主观听感甚至比轻微幅频失真更差。
实际系统到不了理想状态，设计目标就变成：在通带内幅度尽量平、相位尽量线性，阻带部分尽量把镜像压下去。教科书里所谓“近似完全重建”，指的就是这三条在工程误差允许范围内成立。
3.3 怎么用匹配设计把混叠项干掉
一个很自然的做法是让综合滤波器与分析滤波器“互相对着干”：分析滤波器把频谱搬到某个位置，综合滤波器再从那个位置搬回来，同时做时间反转来补偿相位。具体地，令
[ G_k(z) = z^{-(N-1)} H_k^*(z^{-1}) ]
H_k^*(z^{-1})表示把H_k的系数时间反转并取共轭。此时H_k(z)G_k(z)变成z^{-(N-1)}|H_k(e^{jω})|²，相位部分变得干净，只剩下幅度平方。把所有通道相加后，失真传递函数为：
[ T(e^{j\omega}) = \frac{1}{M}e^{-j\omega(N-1)}\sum_{k=0}^{M-1}|H_k(e^{j\omega})|^2 ]
注意|H_k(e^{jω})|²等于|p(e^{j(ω-2πk/M)})|²，所以括号里是对原型滤波器的幅度谱做M个等间隔搬移后求和。要让这个和为常数，原型滤波器必须满足“M通道功率互补”条件。这个条件直观上就是：原型滤波器每个频点上的能量，经过M次搬移后恰好无缝覆盖、不重不漏。这一点可以通过p[n]的多相序列正交性来保证：
[ \sum_m p[n+mM],p[l+mM] = \frac{1}{M}\delta[n-l] ]
换句话说，M个多相分量两两正交且能量归一。很多设计方法（窗函数法、特征滤波法、迭代优化法）最终都是为了逼近这个关系。
3.4 有限长滤波器的残差混叠有多大
理想条件下无限长原型才可能阻带为零，实际FIR滤波器总有一段过渡带和阻带泄漏。因此临界采样的DFT/QAM滤波器组几乎不可能做到严格完全重建，残差混叠的幅度大致与原型滤波器的阻带衰减同量级。举个例子：M=8、原型滤波器阻带衰减60dB时，无量化的重建SNR通常就是45到50dB左右，瓶颈不是计算精度，而是混叠没有完全消除。
这个结论提醒我们：如果压缩系统的量化误差目标在40dB以下，阻带衰减40~50dB的原型可能勉强能用；但如果想追求无损或高保真重建，就必须用过渡带更陡、阻带更深的设计，或者改用余弦调制滤波器组或MDCT这类可以严格PR的路线。认清这一点，后面选型才有依据。
4. 实现与调通：一套可复现的QAM滤波器组
4.1 原型滤波器的设计选择
工程上最省事的原型设计是用窗函数法。M=8，取L=6，原型长度N=48，Kaiser窗参数β取8左右，截止频率设为1/M（对应π/M）：
import numpy as np from scipy.signal import firwin M, L = 8, 6 N = M * L # 窗函数法设计原型低通，截止归一化频率1.0/M p = firwin(N, 1.0/M, window=('kaiser', 8.0)) # 粗略归一化，让功率互补条件尽量满足 p = p * np.sqrt(M)这里有个细节：firwin设计的是线性相位滤波器，中心对称，这对相位条件有好处。截止频率取1/M意味着每个子带宽度接近2π/M，相邻子带边缘刚好相接。但窗函数法得到的滤波器功率互补性一般，残差混叠会偏高；想做得更好就要用专门设计平方根奈奎斯特型原型。作为作业验证，上面的代码足够跑通全流程，把注意力放在滤波器组本身而不是滤波器设计的细节上。
4.2 生成分析和综合滤波器组
有了原型之后，生成M个分析滤波器再做时间反转得到M个综合滤波器：
n = np.arange(N) # 分析滤波器：复指数调制 H = np.array([p * np.exp(2j * np.pi * k * n / M) for k in range(M)]) # 综合滤波器：时间反转 + 复共轭调制 G = np.array([H[k][::-1].conj() for k in range(M)]) # 延迟补偿：由于时间反转，系统总延迟是 N-1 delay = N - 1综合滤波器选择H[k]的时间反转共轭，对应前面说的匹配设计。这样每对分析-综合滤波器的乘积在频域是模平方，相位被对齐成线性延迟。实际运算时复信号的处理要小心：分析得到的每个子带都是复数，后面量化时要分别处理实部和虚部，或者转换成幅度和相位再量化。
4.3 脉冲响应验证：先别急着压数据
任何滤波器组实现完成后，第一件事是用单位脉冲做测试。输入x[n]只在n=0为1，其余为零，观察输出应该在n=N-1处出现一个脉冲，其他位置应当接近零。下面这段是信号经过完整分析-综合级联的验证：
N_signal = 2048 x = np.zeros(N_signal) x[0] = 1.0 # 分析 subbands = np.array([np.convolve(x, hk)[:N_signal] for hk in H]) subbands = subbands[:, ::M] # 下采样M倍 # 综合 up = np.zeros((M, subbands.shape[1] * M)) for k in range(M): up[k, ::M] = subbands[k] y = np.sum([np.convolve(up[k], G[k])[:N_signal] for k in range(M)], axis=0) # 对比 align = y[delay] error_power = np.sum(np.abs(y - np.roll(x, delay))**2)如果error_power相对能量在-40dB以下，说明级联基本正常。如果看到明显的大误差，先不要怀疑量化，优先检查原型滤波器长度是否M整数倍、综合端是否真的做了时间反转、延迟补偿是不是算错了。这三个问题几乎覆盖了初学阶段的现场犯罪。
4.4 多相结构：把计算量降到可以接受
直接按卷积实现分析滤波器组，每个输出样本要算M×L次乘加，M=8时还好，M=64时就是好几倍的浪费。标准做法是多相分解加DFT：输入按M相抽帧，每相经过一个长度为L的子滤波器，再做M点DFT得到子带信号。复杂度从O(M²L)降到O(ML+M log M)。综合端镜像操作即可。
多相结构还有一个额外好处：它能天然保证N=ML的整数倍约束，代码层面不容易写出“多相分支长度不一致”的低级错误。建议封装成两个函数，一个analysis、一个synthesis，内部用多相切分和DFT完成，而不是反复调np.convolve，否则长信号处理时内存和时间都会很难看。我第一版就是用笨方法卷积跑的，处理一分钟音频等了半天；改成多相后快了一个数量级。
5. 实测中的失真来源与排查顺序
5.1 症状一：单频测试出现镜像鬼影
用单一频率正弦作为输入，假设是1kHz，采样率16kHz，M=8。正常输出应该只有1kHz这一根谱线加一点数值噪声。如果输出频谱在镜像位置（如3kHz、5kHz、7kHz）出现明显谱线，说明混叠对消失效。排查方向按顺序来：
- 原型滤波器阻带衰减够不够？窗函数法Kaiser窗β=8大约给到60dB以上的阻带衰减，若发现镜像谱线在-40dB以上，先加长N或加大β。
- 原型带宽是不是太宽？截止频率如果大于1/M，相邻子带重叠严重，功率互补条件被破坏，混叠自然抬头。
- 综合端有没有做时间反转？漏掉[::-1]会让相位匹配失效，镜像谱线会非常整齐地出现在所有l≠0位置。
单频测试是排查混叠最快的手段，比看白噪声误差功率直观得多。建议把它写进回归测试，参数一改就重跑。
5.2 症状二：波形形状没问题但整体错位
如果误差是“输出很像输入，但在时间轴上偏移了几个点”，这是延迟补偿没对齐。分析-综合级联中时间反转带来的延迟是N-1，但多相结构可能还会引入额外延迟，具体和实现方式有关。碰到这种情况，不要拍脑袋改代码，做一个交叉相关：
corr = np.correlate(y, x, mode='full') peak_pos = np.argmax(np.abs(corr)) estimated_delay = len(x) - 1 - peak_pos拿计算出来的实际延迟去对齐，再用对齐后的信号算SNR。另一个容易出问题的地方是：如果你在子带端做了重采样、滤波或帧处理，每一步都可能带来额外延迟，最终系统延迟不再是N-1。复杂系统里建议每加一个模块就重新测一次延迟，比最后阶段一次性找问题省力得多。
5.3 症状三：浮点没问题，一改定点底噪飙升
纯Python/Numpy浮点环境里，重建SNR能做到非常接近理论值。但实际压缩系统尤其是嵌入式平台会用定点运算，滤波器系数一旦量化成16位整数，底噪可能直接抬到-60dB以上，甚至压过弱信号的混叠水平。测试下来，系数用24位比较稳妥，子带谱值建议保留32位累加精度。常规经验：每减少一位精度，底噪大约恶化6dB。做定点模型时先按浮点测出基线SNR，再逐级降低精度，观察误差曲线何时开始偏离，以这个拐点作为最小位宽。
5.4 症状四：压缩后拼起来的信号咔哒作响
滤波器组按块处理时，边界处容易出现点击声。根因是独立处理每个块时，滤波器状态没有跨块延续，块边界处的信号被生硬截断，综合端又产生了瞬态。解决思路有两种：一是带重叠的处理，块之间留一部分样本重叠，综合时用窗函数渐入渐出；二是保持滤波器的状态向量在块间传递，让分析-综合过程连续运行。第一种实现简单但会带来少量码率开销，第二种保持临界采样效率但状态管理更复杂。音频编码里普遍用重叠变换就是出于这个原因，MDCT的50%重叠既保证了完全重建，又天然解决了块边界问题。
| 症状 | 主要原因 | 排查方法 |
|---|---|---|
| 镜像频谱出现 | 混叠未对消 | 单频输入，检查原型阻带与时间反转 |
| 波形错位 | 延迟未补偿 | 交叉相关估计实际延迟 |
| 底噪偏高 | 系数位宽不足 | 浮点基线对比，逐级降位宽 |
| 帧边界咔哒 | 块处理不连续 | 改用重叠块或传递滤波器状态 |
| 高频细节模糊 | 子带带宽过宽/过窄 | 画各子带幅频响应检查覆盖 |
这张表基本覆盖了我在这个项目里踩过的坑，照着从前往后查，大多数问题都能定位。
6. 把完全重建的价值变现：压缩系统里的量化和比特分配
6.1 从分解到编码的完整链路
一套完整的滤波器组压缩系统长这样：
- 分析滤波器组：把信号拆成M个子带，每个子带降采样M倍；
- 对每个子带计算能量/方差，决定比特分配；
- 按分配的比特数对子带信号做量化（常用标量量化或带死区的均匀量化）；
- 可选的熵编码：对量化索引做Huffman或算术编码；
- 解码端量化值还原，子带上采样、综合滤波器组叠加；
- 延迟对齐，得到重建信号。
完全重建保证的是第1步和第5步本身不引入失真，让整个系统的误差尽量只来自第3步的量化。这样误差分析就被简化成量化噪声问题，理论干净、实现也好预测。
6.2 最优比特分配：钱要花在方差大的子带
给定总比特预算B，每个子带分多少比特？在均方误差准则下，最优解是对数量化器背景下各子带独立编码的经典结果：
[ R_k = \frac{B}{M} + \frac{1}{2}\log_2 \frac{\sigma_k^2}{\left(\prod_{j=0}^{M-1}\sigma_j^2\right)^{1/M}} ]
直观理解：方差大的子带多分比特，方差小的少分；所有子带平均下来仍等于总预算。给一个具体数字例子：M=8，子带方差从高到低设为[100, 50, 20, 10, 5, 2, 1, 0.5]，总预算每样本4比特，那么高方差子带可能分到6~7比特，低方差子带只分到1~2比特，甚至直接置零。对比全带等比特分配，这种方式通常能多拿到3~6dB的SNR增益。计算时记得所有权重的单位要保持一致，实际实现一般用拉格朗日乘子法逐带搜索。
6.3 QAM滤波器组的参数权衡：M取大还是取小
子带数量M是设计系统的核心旋钮。M大，频域分辨率高，能把能量集中到更窄的频带里，比特分配更精细；但代价是原型滤波器要做得更陡，否则相邻子带混叠严重，需要更长的N，计算量和延迟一起涨。M小，实现简单、延迟低，但频带切得粗，一个子带里可能同时包含信号和噪声，压缩效率上不去。
音频编码常见的经验值是32到64个子带，正好对应听觉临界频带的粗略划分。图像编码则更多用树形结构或小波包，本质上是一种非均匀的滤波器组划分。QAM复滤波器组的好处是子带信号带相位信息，方便做能量感知或空间线索估计，坏处是复数运算量和I/Q两路量化开销比实数滤波器组大。如果你的系统不太需要相位信息，单纯追求高保真压缩，余弦调制或MDCT路线往往更划算。
6.4 滤波器组和DCT、小波不是竞争对手，是一个框架
JPEG用的DCT是8×8分块变换，本质上是一个固定基底的块变换；JPEG 2000用小波变换，可以看作可递归分解的双通道滤波器组；音频里的MDCT又是带时间混叠消除的余弦调制滤波器组。把这些技术放在“完全重建滤波器组”这个框架下看，它们都是不同约束下对同一个问题的解。QAM滤波器组是这个家族里偏数学完备、偏通信血统的一支，它的价值在于帮你理解所有变换编码的共同底层逻辑：先找一组能完全重建的基底或通道，把信号拆到不同的频域/时域位置，再按感知和能量重要性分配比特。
做压缩选型时不必过度纠结“谁更高级”。DCT简单、块效应明显；小波多分辨率好、适合图像；MDCT重叠变换适合音频；QAM复滤波器组适合需要同时做信道化分析和压缩的场景。理解了完全重建条件之后，这些选项在你眼里就不是黑盒了，而是同一道数学题的不同边界条件。
最后再说两句实在的
做完这轮作业，最大的体会是：调试滤波器组一定要用脉冲、单频、白噪声三种测试信号分别过一遍。脉冲测延迟和时域对齐，单频测混叠对消程度，白噪声测整体SNR和误差频谱平坦度。这套流程帮我定位过至少四个问题，以后设计其他滤波器组也沿用，非常省事。
另一个心得是：如果只是想给压缩系统配个前端，不建议一上来就用临界采样的DFT/QAM滤波器组，那个路线的严格完全重建条件确实难伺候。先跑通余弦调制或MDCT这类成熟方案，理解清楚“时间混叠消除”怎么运作，再回头看复指数版本，很多困惑会自然解开。滤波器组理论里的坑往往不是单个概念看不懂，而是概念之间的关系没理顺；顺着“原型滤波器-调制-多相-完全重建”这条线走，比东看一眼西看一眼高效得多。
