简介：这份资源面向信号处理初学者与工程技术人员，提供小波分解与小波包分解的MATLAB实现范例，帮助理解非平稳信号在时间与频率上的联合分析方法。压缩包内共2个文件，均为m脚本，分别对应小波分解与小波包分解两条处理路径，包体约2KB，轻量易读，适合直接运行调试。代码覆盖数据加载、预处理、小波基选择、分解系数获取、结果可视化、阈值去噪与信号重构等环节，并对比原始信号与重构信号以评估分解效果。读者可据此掌握cwt、swt、wavedec、waverec及小波包相关函数的调用方式，通过调整小波基与分解层数观察差异，逐步建立对小波分析原理的直观认识。目前已有3852人学习下载，适合作为课程实验、故障诊断与图像分析等方向的入门参考。
1. 小波包分解到底比小波分解多算了什么
如果你只用过wavedec做小波分解，可能会觉得“分解”这件事已经到头了：一层一层往下拆低频，高频直接扔掉。但真正做过轴承故障诊断、肌电信号去噪或者电能质量分析的人都知道，小波分解只拆低频，高频部分不再细分，这就导致高频段的分辨率被浪费了。小波包分解（Wavelet Packet Decomposition）干的事，就是把每一层的高频和低频都继续往下拆，让整个频带被均匀切分。这份 MATLAB 程序包就是围绕这个逻辑写的，包含分解、重构、节点能量计算和可视化几个核心模块。它适合两类人：一类是刚接触小波包、想找一个能跑通的最小实现来改参数的新手；另一类是做信号特征提取、需要快速验证节点能量分布是否可分的老手。下面我按“先跑通、再改参数、最后避坑”的顺序拆一遍。
2. 小波包分解的数学逻辑与 MATLAB 实现路径
2.1 从小波分解到小波包：差在递归方式上
普通小波分解的递归只走低频分支。假设信号采样率是 12 kHz，做 3 层分解，低频分支的频带划分是：0–1.5 kHz、0–3 kHz、0–6 kHz，高频分支只保留 6–12 kHz、6–9 kHz 这些粗粒度结果。小波包分解则把每一层的低频和高频都当作新的输入，继续用同一组正交镜像滤波器（QMF）往下拆。3 层小波包分解会产生 2^3 = 8 个节点，每个节点占据原信号频带的 1/8，也就是 1.5 kHz 带宽。这个均匀切分特性，是它适合做频带能量特征提取的根本原因。
MATLAB 里对应的是wpdec和wprec两个函数。wpdec返回一个小波包树对象，节点命名规则是(层号, 节点序号)，比如(3,0)是第 3 层第 0 个节点，对应最低频段。节点序号从左到右按频率从低到高排列，但要注意：高频分支在每一层会翻转一次顺序，这是正交镜像滤波器组的固有特性，不是 bug。如果你发现节点能量排序和预期频带对不上，先检查这个翻转。
2.2 用 wpdec 建树并提取节点系数
下面这段代码是最小可运行版本，直接复制到 MATLAB 脚本里就能跑。我用的是一段模拟信号：两个正弦分量加白噪声，采样率 1200 Hz，做 3 层小波包分解，小波基选db4。
% 小波包分解最小示例 fs = 1200; % 采样率 t = 0:1/fs:1-1/fs; % 1 秒时间轴 f1 = 50; f2 = 180; % 两个信号分量频率 x = sin(2*pi*f1*t) + 0.6*sin(2*pi*f2*t) + 0.3*randn(size(t)); wname = 'db4'; % 小波基 level = 3; % 分解层数 T = wpdec(x, level, wname); % 小波包分解，返回树对象 % 提取第 3 层所有节点系数 numNodes = 2^level; coefs = cell(1, numNodes); for k = 0:numNodes-1 coefs{k+1} = wpcoef(T, [level, k]); end % 计算各节点能量占比 energy = cellfun(@(c) sum(c.^2), coefs); energyRatio = energy / sum(energy); disp('第3层各节点能量占比：'); disp(energyRatio);逻辑说明：wpdec的第三个参数是小波基名称，db4的消失矩是 4，对正弦类信号能提供较好的频带隔离。wpcoef按节点索引取系数，节点索引是[层号, 节点序号]的二元组。能量计算用sum(c.^2)，这是离散信号能量的标准定义。最后归一化得到占比，方便判断哪个频带占主导。
参数怎么改：采样率fs必须和实际信号一致，否则频率轴对不上。分解层数level受信号长度限制，wpdec要求信号长度能被 2^level 整除，否则会自动补零，补零会引入边缘效应。小波基选择上，dbN系列适合瞬态信号，symN系列对称性更好，coifN适合需要近似线性相位的场景。我一般先用db4试，如果重构误差大再换sym8。
2.3 节点重构与信号合成验证
分解完不验证重构，等于白做。小波包重构用wprec，它把树对象还原成时域信号。下面代码接着上面的变量继续跑：
% 全节点重构，验证分解无损性 x_rec = wprec(T); % 计算重构误差 err = norm(x - x_rec) / norm(x); fprintf('相对重构误差：%.2e\n', err); % 只重构特定节点（比如节点 3 和节点 5） nodes = [3, 5]; x_partial = wprec(T, nodes);逻辑说明：wprec(T)不带节点参数时做全重构，理论上误差应该在 1e-12 量级，如果误差大于 1e-6，说明分解层数过高或者信号长度不足导致补零过多。wprec(T, nodes)只重构指定节点，用于提取特定频带成分。注意节点序号从 0 开始，nodes = [3, 5]对应第 3 层第 3 和第 5 个节点。
这里有个容易翻车的地方：wprec的节点参数只接受同一层的节点，不能混层。如果你写[2,1; 3,4]这种跨层组合，MATLAB 会报错。正确做法是先确定你要分析的频带，统一在同一层取节点。
3. 节点能量特征提取与频带排序的实操细节
3.1 能量特征向量怎么构造才不踩频带翻转的坑
上一章算出的energyRatio是一个 8 维向量，直接拿去做分类特征行不行？行，但有个前提：你得确认节点顺序和实际频带顺序一致。前面提过，高频分支在每层会翻转一次。以 3 层分解为例，节点序号 0 到 7 对应的实际频带顺序并不是单调递增的，而是按格雷码顺序排列。具体对应关系可以用wpcoef逐节点做 FFT 来验证：
% 验证节点序号与实际频带的对应关系 for k = 0:7 c = wpcoef(T, [3, k]); C = fft(c); f = (0:length(C)-1) * fs / length(C); [~, idx] = max(abs(C(1:floor(end/2)))); fprintf('节点(3,%d) 主频约 %.1f Hz\n', k, f(idx)); end跑完你会看到节点主频不是按 0、1、2、3 顺序递增的。如果直接把energyRatio按节点序号喂给分类器，而训练集和测试集的信号长度不同导致补零位置变化，频带翻转模式可能不一致，分类精度会莫名其妙下降。稳妥做法是：先按主频排序，再重排能量向量。排序代码很简单：
% 按主频重排节点能量 nodeFreq = zeros(1, 8); for k = 0:7 c = wpcoef(T, [3, k]); C = fft(c); f = (0:length(C)-1) * fs / length(C); [~, idx] = max(abs(C(1:floor(end/2)))); nodeFreq(k+1) = f(idx); end [~, sortIdx] = sort(nodeFreq); energySorted = energyRatio(sortIdx);这样得到的energySorted才是按频带从低到高排列的特征向量。我见过不止一个人在这里翻车，特征向量顺序不一致，换一组数据就失效。
3.2 用 wpviewcf 看频带分布，比手动 FFT 快
MATLAB 小波工具箱里有个wpviewcf函数，能直接画出小波包系数的时频着色图。用法是wpviewcf(T, 1)，第二个参数是色图模式，1 表示用jet。这个图能直观看到哪个节点在哪个时间段有能量集中。对于非平稳信号，比如冲击成分，时频图比单纯的能量占比更有说服力。
不过wpviewcf有个限制：它只显示节点系数，不显示实际频率刻度。你得自己根据采样率和分解层数换算。换算公式是：第 L 层第 k 个节点的频带范围是[k*fs/(2^(L+1)), (k+1)*fs/(2^(L+1))]。注意这里分母是 2^(L+1)，不是 2^L，因为小波包分解每次滤波后做二抽取，实际频带减半。
3.3 批量处理多个信号时的参数固化策略
如果你要处理一批信号，比如 100 组轴承振动数据，每组都手动改wpdec参数不现实。常见做法是把分解参数写成一个配置结构体，循环调用：
cfg.wname = 'db4'; cfg.level = 3; cfg.fs = 12000; files = dir('*.mat'); features = zeros(length(files), 2^cfg.level); for i = 1:length(files) data = load(files(i).name); x = data.signal; T = wpdec(x, cfg.level, cfg.wname); for k = 0:2^cfg.level-1 c = wpcoef(T, [cfg.level, k]); features(i, k+1) = sum(c.^2); end end % 归一化 features = features ./ sum(features, 2);这段代码的关键是features矩阵按行归一化，保证每个样本的能量占比之和为 1。如果不归一化，信号幅值差异会直接反映到特征上，分类器会学到幅值而不是频带分布。另外dir('*.mat')的文件顺序在不同操作系统上可能不同，如果标签和文件顺序有关，建议用sort按文件名排一下。
4. 避坑与排查：小波包分解最常见的五个翻车点
4.1 重构误差大得离谱
现象：wprec重构出来的信号和原始信号波形对不上，相对误差大于 1e-3。
原因：最常见的是信号长度不是 2^level 的整数倍，wpdec自动补零，补零位置在信号末尾，重构时补零段也被还原，导致边缘失真。另一个原因是小波基的正交性不好，比如用了bior系列双正交小波，分解和重构滤波器不满足完美重构条件。
解决：先检查length(x) / 2^level是否为整数，不是的话手动截断到最近整数倍。小波基换成dbN或symN正交小波。如果必须用双正交小波，确认分解和重构用的是同一组滤波器。
4.2 节点能量占比和预期频带对不上
现象：明明信号主频在 50 Hz，但能量最大的节点序号不是最低频节点。
原因：频带翻转。前面说过，高频分支每层翻转一次，节点序号和实际频带不是单调对应。
解决：用wpcoef加 FFT 逐节点验证主频，按主频重排节点顺序。或者直接用wpviewcf看时频图，肉眼确认能量集中位置。
4.3 分解层数选多了，特征反而变差
现象：从 3 层加到 5 层，分类精度不升反降。
原因：层数越多，每个节点带宽越窄，节点系数长度越短。5 层分解时每个节点只有原信号长度的 1/32，如果原信号只有 1024 点，节点系数只有 32 点，能量估计的方差急剧增大。另外层数过多会导致频带混叠加重，相邻节点之间相关性变强。
解决：节点系数长度至少保留 64 点以上。对于 1024 点信号，3 层分解是上限。如果信号更长，比如 4096 点，可以做到 5 层。经验法则是2^level <= length(x)/64。
4.4 批量处理时内存爆掉
现象：循环处理几百个文件时 MATLAB 卡死或报Out of memory。
原因：每次循环都新建wpdec树对象，旧对象没有及时清除。小波包树对象包含所有节点的系数，内存占用是原信号的 2 倍以上。
解决：在循环末尾加clear T，或者把wpdec放在一个独立函数里，函数返回后局部变量自动释放。如果信号特别长，考虑用wpdec的'mode'参数设为'periodic'，减少边界存储。
4.5 节点序号从 0 开始还是从 1 开始
现象：wpcoef(T, [3, 0])报错索引越界。
原因：MATLAB 小波工具箱的节点序号从 0 开始，但 MATLAB 数组索引从 1 开始。写循环时容易混。
解决：统一用for k = 0:2^level-1遍历节点，取系数时用wpcoef(T, [level, k])，存数组时用coefs{k+1}。记住这个“节点从 0 起，数组从 1 起”的对应关系，能省很多调试时间。
5. 进阶技巧：用节点重构做单频带提取与包络分析
小波包分解最实用的进阶玩法，不是算能量占比，而是只重构特定节点，把某个频带单独拎出来做包络谱。这在轴承故障诊断里很常见：故障冲击成分集中在某个高频段，先用小波包把该频段重构出来，再做 Hilbert 包络，最后对包络做 FFT，看故障特征频率。
下面这段代码演示完整流程。假设采样率 12 kHz，信号长度 4096，做 4 层小波包分解，取第 4 层第 6 个节点（对应频带约 3.75–4.125 kHz）做包络分析：
fs = 12000; level = 4; wname = 'db4'; T = wpdec(x, level, wname); % 只重构节点 (4,6) node = [4, 6]; x_band = wprec(T, node); % Hilbert 包络 env = abs(hilbert(x_band)); % 包络谱 L = length(env); Env = fft(env .* hann(L)'); f = (0:L-1) * fs / L; amplitude = abs(Env(1:floor(L/2))) * 2 / L; % 画图 figure; subplot(2,1,1); plot((0:length(x_band)-1)/fs, x_band); title('节点(4,6)重构信号'); xlabel('时间 (s)'); subplot(2,1,2); plot(f(1:floor(L/2)), amplitude); xlim([0, 500]); title('包络谱'); xlabel('频率 (Hz)');逻辑说明：wprec(T, node)只重构指定节点，输出信号长度和原信号一致，但只保留该节点对应频带的成分。hilbert返回解析信号，取模得到包络。加hann窗是为了减少频谱泄漏。包络谱的频率分辨率是fs/L，12 kHz 采样 4096 点，分辨率约 2.93 Hz，足够分辨常见的轴承故障特征频率。
参数怎么调：节点选择取决于你想分析的频带。先用wpviewcf看哪个节点能量集中，再选那个节点。如果故障特征频率很低（比如 20 Hz），包络谱的横轴范围要相应缩小，xlim改成[0, 200]更清楚。小波基换成sym8有时能得到更平滑的包络，因为对称性更好。
一个我踩过的坑：wprec重构单节点时，输出信号里除了该频带成分，还会残留少量相邻频带的泄漏。如果包络谱里出现不该有的频率峰，先检查分解层数是否足够。层数越多，频带隔离越好，但节点系数越短。折中方案是层数选到节点系数长度 128 点左右。
从那以后我每次做小波包分析，都强制走一遍“全重构验证误差 → 单节点 FFT 验证频带 → 包络谱验证特征频率”这三步，少一步都可能被错误结果带偏。希望帮到你。
本文还有配套的精品资源，点击获取
