简介面向图像处理初学者、算法研究者及需要做小波变换实验的MATLAB开发者提供一套完整的非下采样小波变换实现。代码可将图像分解成一个低频子带近似分量和多个高频子带水平、垂直、对角细节分量在保留原图全部信息的前提下支持多层级分解与重构。资源包共7个文件包含5个可直接运行的M脚本、1份TXT说明文档和1篇PDF论文整体仅276KB。M脚本涵盖基础分解、增强处理与测试函数说明文档对参数设置和调用方式作了备注PDF则从理论角度讲解基于小波变换的图像增强算法。目前已有267人学习下载适合用于图像去噪、压缩、边缘检测与特征提取。下载后结合示例运行可掌握小波基选择、分解层数设置、系数重构等关键操作并可直接将脚本迁移到自己的图像处理任务中。1. 非下采样小波变换为什么丢掉“下采样”才是关键做图像处理的人第一次接触小波大多是从DWT离散小波变换入手的——每分解一层四个子带尺寸减半低频继续往下拆。这套流程写起来顺手但真放到图像融合、边缘增强这类任务里DWT的毛病就藏不住了平移一个像素系数能变一大截子带尺寸不对齐后处理还得做插值对齐。非下采样小波变换NSWT把“下采样”这一步直接去掉每一层子带尺寸与原始图像完全一致配合atrous空洞滤波器组实现带来了平移不变性。输出结构也很单纯一个低频子带加若干个高频子带。如果你手里正需要一个MATLAB算法做图像融合、细节增强或纹理分析这份资源比DWT系列更对口。下面我把分解原理、参数选择和坑一次性拆清楚。2. atrous分解原理与MATLAB实现从滤波器插值到完整代码2.1 去掉了下采样代价是什么DWT的标准做法是信号经过低通和高通滤波器后各自下采样2倍。下采样保证了总数据量不变但引入了两个问题——一是平移敏感性输入平移一个像素各级小波系数会发生明显变化二是子带尺寸随层数递减做融合或增强时需要对不同尺度的子带做几何对齐。NSWT的做法是每层只做滤波不抽稀。第j层使用的滤波器是在基础滤波器系数之间插入2^(j-1)-1个零等效于用尺度越来越大的滤波器对原始分辨率信号连续滤波。代价也很直观输出总数据量是原始图像的层数1倍。以三层分解为例输出四个子带每个尺寸都与原图一致整个表达是冗余的。这份资源换来的是平移不变性和子带间天然的像素级对齐——后者在图像融合里极其重要因为你可以直接在不同子带的相同坐标处做系数比较和取值不需要任何插值对齐操作。2.2 一维滤波器的atrous插值实现atrous分解第一步是构造每层的插值滤波器。假设基础低通滤波器h长度为5第j层需要把步长变为2^(j-1)也就是在原有系数之间插入2^(j-1)-1个零。MATLAB代码可以这样写function hd atrous_filter(h, j) % 生成第j层的atrous插值滤波器 % h: 基础低通滤波器行向量 % j: 当前层数从1开始 % hd: 插零后的滤波器系数 step 2^(j-1); % 新滤波器长度 (原长度-1)*步长1 len (length(h) - 1) * step 1; hd zeros(1, len); % 每隔step个位置放一个原系数 hd(1:step:end) h; end这段代码的逻辑很简单第1层step1滤波器原样使用第2层step2在相邻系数间插1个零滤波器长度从5变成9第3层step4插入3个零长度变成17。重点在于hd(1:step:end)的赋值方式——它把原滤波器系数放到插值后滤波器的第1、第1step、第12*step等位置其余位置补零。这样每一层滤波器覆盖的“感受野”越来越大但卷积之后输出的尺寸始终不变。2.3 二维可分离分解先列后行图像是二维信号二维NSWT的标准做法是通过可分离滤波实现先用一维atrous滤波器对每一列做卷积再对得到的中间结果每一行做卷积。这样做的好处是复杂度低每层只需要两次一维卷积。下面是一个完整的二维分解函数function [L, H] nswt2_decompose(img, filt, level, boundary) % 二维非下采样小波分解atrous实现 % 输入: % img - 灰度图像double类型uint8会自动转换 % filt - 一维低通滤波器如 [1 4 6 4 1]/16 % level - 分解层数正整数常见取值3~5 % boundary- 边界扩展方式: symmetric / replicate / zero % 输出: % L - 低频子带与img尺寸相同 % H - cell数组H{j}为第j层高频子带与img尺寸相同 if nargin 4 boundary symmetric; end if isa(img, uint8) img double(img); end H cell(level, 1); L img; for j 1:level % 生成当前层的atrous滤波器 step 2^(j-1); hd zeros(1, (length(filt)-1)*step 1); hd(1:step:end) filt; % 边界扩展为了避免卷积后边缘衰减 pad floor(length(hd) / 2); Lpad padarray(L, [pad pad], boundary); % 可分离二维卷积先对列卷积再对行卷积 temp conv2(Lpad, hd, same); lf conv2(temp, hd, same); % 裁剪回原始尺寸 lf lf(pad1:end-pad, pad1:end-pad); % 高频子带 当前输入 - 低频近似 H{j} L - lf; % 低频进入下一层 L lf; end end这里有两个细节需要说明。第一H{j} L - lf生成高频子带的方式对应的是atrous小波中的“提升式”实现低频近似和原始输入的差就是这一层丢失的细节信息。这样做的好处是重构时不需要设计专门的重建滤波器所有子带直接相加即可。第二padarray配合conv2的same先扩展后裁剪是为了避免conv2默认补零导致的边缘衰减。对应的重构函数更简单function rec nswt2_reconstruct(L, H) % 从NSWT系数重建图像 % 输入与nswt2_decompose的输出对应 % 输出: rec 与原始图像尺寸一致 rec L; for j 1:length(H) rec rec H{j}; end end重构时把低频子带和所有高频子带逐像素相加由于H{j}的定义就是输入减去低频近似理论上重建误差应该为零。用最大绝对误差来判断一般能到1e-10量级——如果达不到多半是数据类型或边界处理出了岔子后面避坑章会说。2.4 完整测试脚本跑通一遍把上面两个函数存成独立文件后用下面这段脚本做验证% 读取灰度图并转为double img imread(cameraman.tif); img double(img); % 基础低通滤波器5点二项式近似三次B样条 h [1 4 6 4 1] / 16; % 三层分解 level 3; [L, H] nswt2_decompose(img, h, level, symmetric); % 查看各子带尺寸是否与原图一致 fprintf(原图尺寸: %d x %d\n, size(img,1), size(img,2)); fprintf(低频尺寸: %d x %d\n, size(L,1), size(L,2)); for j 1:level fprintf(第%d层高频尺寸: %d x %d\n, j, size(H{j},1), size(H{j},2)); end % 验证重建误差 rec nswt2_reconstruct(L, H); fprintf(重建最大误差: %e\n, max(abs(rec(:) - img(:)))); % 可视化高频用绝对值自动拉伸显示 figure; subplot(2,2,1); imshow(uint8(L)); title(低频子带); for j 1:level subplot(2,2,j1); imshow(abs(H{j}), []); % 自动拉伸动态范围 title([高频子带 Level , num2str(j)]); end运行结果里三张高频子带会呈现越来越“大尺度”的边缘纹理第一层捕捉细边缘和噪声第二层是中等尺度轮廓第三层偏向大块亮暗变化。abs(H{j})的原因是高频系数有正有负直接用imshow会丢失负值信息自动拉伸后细节更清楚。这个测试脚本建议保留换滤波器或换参数时先跑一遍确认重建误差没有变大再继续下一步。3. 参数怎么选滤波器、分解层数与边界模式3.1 基础滤波器选择B样条族是默认起点NSWT对滤波器有一个隐含要求低通滤波器系数之和应为1且频响在低频端尽量平坦。最常用的基础滤波器是二项式系数归一化得到的B样条核滤波器名称系数特性3点二项式[1 2 1]/4阶数低平滑弱高频细节保留多5点二项式[1 4 6 4 1]/16折中选择最常用7点二项式[1 6 15 20 15 6 1]/64平滑强高频子带能量弱我一般从5点二项式起步。理由有三个第一系数和为1直流分量无失真通过第二频响在归一化频率0.5处有足够衰减能撑住三层以上分解而不产生明显混叠第三滤波器是对称的conv2的180度翻转问题可以忽略。如果你处理的对象是强噪声图像可以换7点二项式但要做好高频子带细节也变少的思想准备。不建议用haar或db2这种经典正交小波滤波器来做atrous。原因在于正交小波滤波器设计时依赖临界采样和下采样重建的条件直接插零使用会破坏其性质输出子带的频率分离特性并不理想。NSWT场景下B样条系天然适合atrous插值结构这是实践里的主流做法。3.2 分解层数与图像尺寸和任务目标绑定分解层数怎么定没有绝对数字但有一个硬约束滤波器等效长度不能超过图像尺寸。第j层滤波器的等效长度是(length(h)-1)*2^(j-1)1这个长度超过图像边长时边界扩展区会占据过大比例中间有效区域被压缩分解结果失真。以512×512图像和5点滤波器为例第5层等效长度为(5-1)*16165远小于512没问题但如果图像只有64×64第5层就危险了。实践里我把层数选择跟任务挂钩图像融合3~4层。高频分3层已经能覆盖从细纹理到中等轮廓的尺度范围再深层次的尺度往往和“聚焦区域”判断无关。细节增强2~3层。增强目标通常是边缘和纹理太深的层对应的是大范围亮度变化乘增益后容易破坏整体对比度。去噪4~5层。噪声集中在浅层高频深层高频主要是真实结构分层处理对噪声和信号的可区分性更好。另外提醒一点层数增加会放大边界伪影的影响范围每多一层边界扩展宽度翻倍。所以层数不是越多越好够用就行。3.3 边界模式默认symmetric特殊场景再改我在这份资源里默认的边界模式是symmetric也就是镜像对称扩展。原因很直接自然图像的边界两侧像素通常有连续性镜像扩展不会在边缘处引入阶跃生成的子带边缘伪影最少。下面是三种模式的对比模式扩展方式适用场景典型问题symmetric镜像对称绝大多数自然图像无replicate复制边缘像素纹理强烈、边缘信息关键的图像边缘区域过度平滑zero补零仅用于理论验证边缘出现黑带强烈不推荐zero模式是最容易踩的坑conv2默认就是补零卷积如果不做padarray直接卷积边缘一圈像素会被不会真实存在的零值拉低高频子带边缘会出现一条明显的黑边。我在初始版本里没有处理这个问题后来在某个边缘检测任务里发现第一层高频子带的边缘有一圈规则暗线排查后发现就是这个问题。3.4 一套推荐的参数组合如果是第一次跑NSWT我建议从这套组合开始h [1 4 6 4 1] / 16; % 5点二项式滤波器 level 3; % 三层分解 boundary symmetric; % 镜像对称边界这套组合对512×512级别的灰度图稳定子带尺寸保持512×512重建误差在1e-10量级不会出现边缘黑带或高频全黑的问题。后续再按任务需要调整图像更小就减层数噪声更强就换7点滤波器要保留更多纹理就换3点滤波器。每次改参数我都建议顺手跑一遍重构误差校验误差超过1e-6就说明参数组合里有问题别急着往下走。4. 把子带用起来去噪、融合与增强的具体套路4.1 为什么NSWT子带适合做像素级处理NSWT子带和原图尺寸一致这意味着你可以直接在子带上逐像素操作不需要上采样或插值。这一点在图像融合里尤其有价值DWT的四叉树结构要求融合策略在树的不同层间切换而NSWT只需在同一坐标下比较不同系数的绝对值大小。高频子带本身是当前输入与低频近似的差选绝对值大的系数等价于选“边缘更锐利”的那一侧——这个逻辑在各尺度上一致成立。4.2 用法一基于高频阈值去噪高斯噪声主要分布在浅层高频子带结构细节则分布在所有高频层。经典做法是对浅层高频做软阈值深层高频做保守处理或不动function denoised nswt_denoise(img, h, level, sigma) % NSWT去噪对高频子带做软阈值 % sigma: 噪声标准差估计值可用小波工具包估计 [L, H] nswt2_decompose(img, h, level, symmetric); % 各层阈值第j层阈值 sigma * sqrt(2*log(N)) / 2^(j-1) N numel(img); base_thresh sigma * sqrt(2 * log(N)); for j 1:level thr base_thresh / (2^(j-1)); % 软阈值 Hj_abs abs(H{j}); H{j} sign(H{j}) .* max(Hj_abs - thr, 0); end denoised nswt2_reconstruct(L, H); end这段代码的核心是阈值逐层递减第1层阈值最大抑制噪声最多越往后阈值越小因为深层高频中包含的真实结构占比更高。需要说明的是sigma需要外部估计——常见做法是用图像中高频区域的标准差近似或者先做一次DWT、用第一层高频的对角分量估计噪声水平。4.3 用法二多聚焦图像融合多聚焦融合的任务是两张图像一张聚焦前景、一张聚焦背景合成一幅前后都清晰的图像。NSWT的做法是逐像素比较高频系数的绝对值取大的作为融合后系数function F nswt_fuse(A, B, h, level) % NSWT融合低频取平均高频绝对值取大 % A、B: 配准好的两幅灰度图 [LA, HA] nswt2_decompose(A, h, level, symmetric); [LB, HB] nswt2_decompose(B, h, level, symmetric); % 低频子带简单平均即可 LF (LA LB) / 2; % 高频子带绝对值取大 HF cell(level, 1); for j 1:level HF{j} HA{j}; mask abs(HB{j}) abs(HA{j}); HF{j}(mask) HB{j}(mask); end F nswt2_reconstruct(LF, HF); end这个策略的原理是聚焦区域的像素梯度大高频系数绝对值大离焦区域被光学模糊抹平高频系数接近零。逐像素取大等效于在每个像素位置选择来自聚焦侧的细节不同层的判断结果自然融合。低频取平均是因为离焦和聚焦在低频上的差异不显著平均即可避免块状拼接痕迹。4.4 和DWT版本的结果差异在哪实际对比中DWT融合结果在边缘处容易出现振铃和方块伪影根源是下采样破坏了平移不变性源图像移动半个像素系数的分布就变了融合规则的选择就不稳定。NSWT版本在高频细节的连续性上明显更好代价是内存占用成倍增加和耗时变长。以512×512灰度图为例NSWT三层分解在我的机器上大约耗时0.3秒同条件下DWT不到0.05秒但融合质量提升是值得这个代价的。如果你的任务是视频流实时处理NSWT可能扛不住DWT仍是更务实的选择。5. 非下采样小波变换的五个高频坑现象、原因、解决5.1 高频子带显示一片黑现象imshow(H{1})直接输出时整张图是黑的几乎看不到任何结构。原因高频子带的系数范围远小于0~255。以8位灰度图为例像素值范围是0~255但第一层高频系数通常只有±30左右第三层甚至只有±8。imshow默认把数据最小值映射到0、最大值映射到255如果数据范围只有±30会被压缩到很窄的灰度区看起来就是一片接近黑的灰色。解决显示时用imshow(H{j}, [])让MATLAB自动把当前数据的最小值和最大值拉伸到全动态范围。或者用abs(H{j})取绝对值后再显示这样正负细节都能呈现。注意这只影响显示不影响后续系数处理。5.2 子带边缘出现规则黑带现象分解后的高频子带四周边界处有一条宽度固定的暗线内部看起来正常。原因这是conv2默认补零边界惹的祸。如果不做padarray扩展直接对原始图像卷积边缘像素周围的数据有一半来自人为补的零值。滤芯覆盖到边缘时零值把卷积结果向下拉形成一圈低值区。解决在卷积前用padarray做对称扩展卷积完成后裁剪回原尺寸。这也是我在nswt2_decompose里放pad参数的原因。还有一点容易忽略如果pad floor(length(hd)/2)计算不对裁剪位置会偏差导致子带边缘仍然异常。建议每次修改滤波器长度后都打印一下pad值核对。5.3 uint8图像直接参与卷积重建误差飙升现象把imread返回的uint8图像直接传给nswt2_decompose重建后的误差从1e-10量级变成0.5甚至更大。原因uint8做卷积时MATLAB会把参与运算的数据当作double处理但结果会先被截断为整数数值精度丢失。更隐蔽的是H{j} L - lf这一步里L是uint8、lf是double类型不匹配会导致隐式转换误差进一步放大。解决函数开头加一行类型检查发现uint8就转double。在我的实现里isa(img, uint8)判断后自动调用double(img)完成转换。如果你在外部处理也要记得分解前转double显示时再转回uint8。5.4 conv2与核方向对称滤波器不是万能缓冲现象换了一个不对称滤波器比如某些边缘提取核后分解结果的方向性出现异常水平边缘和垂直边缘的响应不对称。原因MATLAB的conv2执行的是真正的卷积操作核会被翻转180度。对对称核来说翻转等于没翻但非对称核会改变滤波器的方向响应。filter2执行的是相关操作不翻转核和conv2结果不同。解决如果确需使用非对称滤波器统一用filter2并注意方向约定或者手动将核旋转180度传给conv2。在B样条系列滤波器内部这个问题不存在但一旦你开始尝试自定义滤波器方向问题就会冒头。5.5 分解层数过多导致低频“糊成一片”现象把level设成6或7之后低频子带看起来非常模糊几乎只保留了整体亮度分布同时重建误差明显变大。原因层数过深时等效滤波器覆盖范围迅速扩大。以5点滤波器为例第6层等效长度为(5-1)*321129对于512×512图像这个尺度已经接近图像总尺寸的1/4边界扩展区域占有效区域比重过大中心区域的分解结果已经不可靠。同时深层低频包含了极低频率亮度变化这种信息对视觉感知的意义已经很低。解决先根据图像尺寸估计最大可用层数len_max floor(log2(size(img,1) / length(h))) 1一般取这个值的一半作为实际层数。以512×512图像和5点滤波器为例理论最大约7层实际取3~4层即可。6. 再进阶用NSWT高频增益做边缘增强验证重建闭环NSWT最有价值的一个特性是高低频分离后可以独立调整某层系数再重建实现尺度可控的边缘增强。这个操作在DWT里也能做但子带尺寸不同导致增强的尺度范围更难精确控制NSWT的逐像素对应让增强更直观。我常用的一种增强策略是保持低频不变对高频子带按层数加权浅层高频权重大、深层高频权重小因为浅层对应的是人眼敏感的细边缘。代码实现如下function enhanced nswt_enhance(img, h, level, gain) % NSWT细节增强 % gain: 各层高频增益向量长度应为level % 例如 gain [1.8 1.4 1.1] 表示浅层增强更强 [L, H] nswt2_decompose(img, h, level, symmetric); for j 1:level H{j} H{j} * gain(j); end enhanced nswt2_reconstruct(L, H); % 防止溢出裁剪到0~255 enhanced max(min(enhanced, 255), 0); end调用方式img double(imread(cameraman.tif)); enhanced nswt_enhance(img, [1 4 6 4 1]/16, 3, [1.8 1.4 1.1]);验证方法不能只看视觉效果还要检查重建闭环对增强结果再做一次NSWT分解观察各层高频能量的变化是否与增益设置一致。这个过程能确认增强确实发生在预期尺度上而不是某个参数错误导致的整体亮度偏移。我遇到过增益设成[2.0 1.0 1.0]时图像出现明显光晕的情况排查后发现是第三层高频接近零增益设置对它没有意义而第一层增益过强压过了原图的亮度梯度。从那以后我每次做增益增强都会先打印各层高频的RMS能量值确认每一层都还有值得增强的内容再设定增益而不是凭感觉拍一个系数。希望这个习惯对你也有用。本文还有配套的精品资源点击获取