简介：拉普拉斯金字塔图像融合的MATLAB实现，面向数字图像处理初学者及需要多尺度融合算法的开发者，核心解决不同源图像在保留高频细节与边缘信息前提下的高质量合成问题。压缩包共5个文件，包含4个.m脚本和1个txt说明文档；脚本分别负责主流程控制、图像预处理、金字塔构建、逐层融合与最终重建，包体仅3KB，结构紧凑便于逐行研读。目前已有1990人学习下载。通过这套代码，读者可掌握Gaussian金字塔差分构造拉普拉斯金字塔的原理，理解逐层权重分配对融合效果的影响，并可直接借鉴其预处理、融合、重建框架。程序面向二值图像设计，使用时需注意输入格式，但整体思路同样适用于医学影像分析、遥感图像融合等广泛场景，是理论结合实践的优质参考。
1. 拉普拉斯金字塔图像融合：先搞懂它和直接叠加差在哪
做图像融合的人，十有八九第一次都试过img1 * 0.5 + img2 * 0.5。结果很明显：边缘发虚、亮度断层、拼接处能看到一条明显的缝。这不是权重没调好，而是你把像素当成了独立个体——实际上图像的“结构”出现在不同尺度上，锅里的纹理和锅边的轮廓本来就是两回事。拉普拉斯金字塔图像融合的核心思路，就是把图像拆成“低频背景 + 各级高频细节”的多层结构，每一层单独决定怎么融合，最后再倒着重构回去。
这套方法最常用的场景是医学图像融合（CT 和 MRI 互补）、多聚焦图像的景深扩展、红外与可见光叠加，以及任何两种配准好的图像需要“把各自优点合到一张图里”的任务。它不挑传感器，只要输入图像已经对齐，就能稳定工作。适合的读者是：用 OpenCV 做过基本图像处理、想摆脱简单加权、理解金字塔原理但没亲手实现过的工程师。下面从数学结构讲起，再到可跑的代码和参数调优，最后给出一个能用的验证技巧。
2. 拉普拉斯金字塔的数学结构：高斯金字塔差分与重构闭环
2.1 高斯金字塔先做“降采样 + 模糊”，拉普拉斯层是两者之差
要理解拉普拉斯金字塔，得先看它的上游——高斯金字塔。给定原图G_0，对它做一次高斯模糊再隔行隔列采样，得到G_1；对G_1重复同样操作得到G_2，以此类推。每一层面积是上一层的 1/4，像素点代表的是更大范围的局部平均。这里高斯模糊不是随便选的，5×5 高斯核（σ≈1）能有效抑制降采样产生的混叠。
拉普拉斯金字塔的每一层L_i定义为：
L_i = G_i - expand(G_{i+1})其中expand()是把低分辨率图放大回当前尺寸的操作：先插值补零（通常是双线性或最近邻），再做高斯模糊消除棋盘格效应。换句话说，拉普拉斯层保存的是“当前尺度下的细节残差”，低频部分被完全剥掉。
提示：拉普拉斯金字塔存的是“带符号的差分图”，不是自然图像，像素范围可能包含负值。后面做融合和可视化时要注意这一点。
2.2 重构就是差分的逆过程，不是简单相加
重构时，从金字塔最顶层G_N（一个极小的低频图）开始，逐层做：
G_i = L_i + expand(G_{i+1})往上一层，细节层被加回到放大的低频图上。这个闭式意味着：拉普拉斯金字塔本身不丢失任何信息——把全部分解出来的L_0, L_1, ..., L_N和G_N完整保留，理论上能无损失恢复原图（浮点误差除外）。
这个性质是融合算法的基石。假设有两张图 A 和 B，各自分解成L_A和L_B，我们做的是逐层选择或加权，然后再倒着重构。因为每一层都对应一个特定的频带，融合策略可以在不同频带独立决策——低频层决定“整体亮度跟谁走”，高频层决定“边缘细节保留谁的”。
2.3 为什么逐层选最大就能出效果，直接 max 却不行
一个常见的误区是：既然拉普拉斯金字塔保存细节，那直接对原图做np.maximum(A, B)不就行了？实际效果通常很差，因为这会引入剧烈的灰度跳变，看起来像噪声叠加。原因在于：np.maximum是像素级操作，没有区分“一个像素属于哪个尺度的结构”。
而拉普拉斯金字塔分层后，每一层都是带通滤波结果，层内的像素值本身就代表了“该尺度下的局部能量”。在同一尺度下做选大或加权，不会跨尺度打架。融合规则可以很朴素——高频层保留绝对值大的系数，低频层做加权平均——就能得到边缘清晰、背景自然的输出。
另外，金字塔的层数决定了融合的“粒度”。层数太少，细节层太粗，融合结果接近直接叠加；层数太多，最顶层只剩一个极小的 DC 分量，计算收益很低。常见选择是 4～6 层，具体和图像尺寸有关。
3. 用 OpenCV 在本地跑通拉普拉斯金字塔融合的最小实现
3.1 准备工作与输入要求
先明确输入条件：两张图像需要已经配准，即同一场景、同一尺寸、同一通道数。如果不满足，先做特征点匹配或基于互信息的配准，拉普拉斯金字塔融合本身不负责对齐。
依赖只需要 OpenCV 和 NumPy。测试时用两张多聚焦图最直观：一张前景清晰、背景模糊，另一张相反。没有现成数据的话，可以对一张清晰图分别做高斯模糊，但只模糊一半区域，再合成两张互补图。
import cv2 import numpy as np def build_gaussian_pyramid(img, levels): """构建高斯金字塔，返回从原图到最顶层的列表""" gauss = [img.copy()] for _ in range(levels - 1): img = cv2.pyrDown(img) gauss.append(img) return gauss def build_laplacian_pyramid(img, levels): """拉普拉斯金字塔：高斯层与上一层上采样之差""" gauss = build_gaussian_pyramid(img, levels) lap = [] for i in range(levels - 1): size = (gauss[i].shape[1], gauss[i].shape[0]) up = cv2.pyrUp(gauss[i + 1], dstsize=size) lap.append(cv2.subtract(gauss[i], up)) lap.append(gauss[-1]) # 最顶层直接用低频残差 return lap逻辑说明：build_gaussian_pyramid里cv2.pyrDown默认使用 5×5 高斯核，内部先做卷积再隔行采样；build_laplacian_pyramid里cv2.pyrUp先插值再卷积，dstsize必须显式指定为目标尺寸，否则 OpenCV 按默认 2 倍尺寸计算，当层数多时可能因奇偶尺寸产生偏差。
3.2 逐层融合与重构：选大策略与加权策略
def fuse_pyramids(lap_a, lap_b, mode='max'): fused = [] for la, lb in zip(lap_a, lap_b): if mode == 'max': # 高频细节保留绝对值更大的系数 mask = np.abs(la) > np.abs(lb) fused.append(np.where(mask, la, lb)) elif mode == 'avg': # 低频背景用平均，保留全局亮度 fused.append((la + lb) / 2.0) # 实际使用中可以对高层用 avg，对低层用 max，按 levels 切分 return fused def reconstruct_from_laplacian(lap): """从拉普拉斯金字塔自顶向下重构图像""" img = lap[-1] for i in range(len(lap) - 2, -1, -1): size = (lap[i].shape[1], lap[i].shape[0]) img = cv2.pyrUp(img, dstsize=size) img = cv2.add(img, lap[i]) return img # 主流程 img_a = cv2.imread('focus_left.jpg') # 左清晰 img_b = cv2.imread('focus_right.jpg') # 右清晰 levels = 5 lap_a = build_laplacian_pyramid(img_a, levels) lap_b = build_laplacian_pyramid(img_b, levels) fused_lap = [] for i in range(levels): if i < levels - 1: fused_lap.append(fuse_pyramids([lap_a[i]], [lap_b[i]], 'max')[0]) else: fused_lap.append(fuse_pyramids([lap_a[i]], [lap_b[i]], 'avg')[0]) result = reconstruct_from_laplacian(fused_lap) result = np.clip(result, 0, 255).astype(np.uint8) cv2.imwrite('fused_result.jpg', result)参数说明：levels=5对 512×512 以下的图像足够，最顶层是 32×32 的低频分量；mode='max'泛化到彩色图像时，三通道分别取绝对值的最大值即可。
3.3 细节层处理的两个隐藏难点
第一，cv2.subtract会产生截断。如果输入是uint8，G_i - expand(G_{i+1})的负值会被截断到 0，拉普拉斯系数丢失符号信息，重构时出现偏色。解决方法是使用cv2.subtract的替代——转成float32再做减法。
def build_laplacian_pyramid_float(img, levels): img = img.astype(np.float32) gauss = [] cur = img for _ in range(levels - 1): gauss.append(cur) cur = cv2.pyrDown(cur) gauss.append(cur) lap = [] for i in range(levels - 1): size = (gauss[i].shape[1], gauss[i].shape[0]) up = cv2.pyrUp(gauss[i + 1], dstsize=size) lap.append(gauss[i] - up) # 直接相减，保留负值 lap.append(gauss[-1]) return lap第二，cv2.pyrUp的插值方式在拉普拉斯重构时会产生振铃。OpenCV 的实现是双线性插值后再做高斯卷积，这不是理想的抗混叠插值。对于绝大多数应用可以接受，但如果你发现融合结果轮廓周围出现一圈亮边，可以尝试把pyrUp替换成cv2.resize配合cv2.GaussianBlur，尺寸控制更灵活。
提示：把
levels从 3 调到 8，融合结果的变化能直观体现出金字塔层数的意义——层数越深，低频分量越扁平，融合后整体色调越统一。
4. 多尺度融合实验设计：多聚焦、医学图像、红外可见光
4.1 用合成数据做回归测试：先验证代码正确性
在拿真实数据之前，先做一个闭环测试：把一张图分解再重构，计算和原图的误差。正常情况误差应小于 1e-4（float32）。如果误差很大，通常是pyrUp的dstsize没对齐，或减法截断导致。
original = cv2.imread('lena.jpg').astype(np.float32) lap = build_laplacian_pyramid_float(original, 5) recon = reconstruct_from_laplacian(lap) err = np.abs(original - recon).max() print(f"重构误差: {err:.6f}")误差来源有两类：浮点舍入误差（正常，可忽略），以及pyrDown/pyrUp组合在奇偶尺寸下无法完美恢复高频信息（异常）。后者是因为 OpenCV 的pyrDown默认核和pyrUp的核是配套的，但尺寸下取整导致最后一像素被丢弃。验证时建议先用偶数尺寸的图。
4.2 多聚焦图像：以边缘能量为准的融合效果检验
多聚焦图像是最经典的应用：把两张分别聚焦在不同深度的照片融合成全景清晰图。融合后质量如何评估？不能只看眼睛，还得看高频能量分布。一个简单指标：融合图的拉普拉斯方差（Variance of Laplacian）应该比任一输入高。
def variance_of_laplacian(img): lap = cv2.Laplacian(img, cv2.CV_64F) return lap.var() for name, img in [('A', img_a), ('B', img_b), ('Fused', result)]: print(f"{name}: VoL = {variance_of_laplacian(img):.2f}")这个值代表图像整体锐度。如果融合结果的 VoL 接近 max(VoL(A), VoL(B)) 甚至更高，说明细节保留策略有效。如果 VoL 反而下降，问题基本出在两个地方：一是levels太少，高频融合发生在过大的尺度上；二是mode='avg'被错误用在了高频层。
4.3 医学图像与红外可见光：低频策略要区分对待
医学图像融合（比如 CT 与 MRI）和红外可见光融合在低频策略上是相反的。CT 擅长骨骼等硬组织，MRI 擅长软组织，融合时低频决定整体对比度；如果低频简单平均，骨骼区域对比度会被 MRI 拉低。常见做法是低层采用“局部能量取大”：以像素为中心开 3×3 或 5×5 窗口，计算窗口内能量（像素绝对值的平方和），选能量大的那个系数。
def local_energy_weight(la, lb, win=3): """局部能量取大融合策略，win 为窗口边长""" kernel = np.ones((win, win), dtype=np.float32) energy_a = cv2.filter2D(la**2, -1, kernel) energy_b = cv2.filter2D(lb**2, -1, kernel) mask = energy_a > energy_b return np.where(mask, la, lb)参数说明：win越大，融合结果越平滑、对局部噪声越不敏感；win=3保留细节，win=7更像区域决策。这个策略适合后三层（高频），最顶层仍然用平均或能量取大。
红外可见光融合则相反，可见光图像的高频纹理丰富，红外图像低频热辐射信息更关键。常见策略是：高频层直接取可见光（或按能量加权），低频层取红外占比更高（0.6~0.7）。这就是拉普拉斯金字塔的最大优势——每一层可以独立指定不同的规则。
提示：以上三种策略对应三种实际需求：多聚焦要保边缘，医学要保结构对比度，红外要保热目标与纹理。融合代码不变，变的只是金字塔各层参数。
5. 拉普拉斯金字塔融合的必调参数：层数、核尺寸与融合规则
5.1 金字塔层数与图像尺寸的关系
层数直接决定金字塔最顶层的分辨率。最顶层不能太小，否则低频信息被过度压缩，重构时产生块状伪影；也不能太大，否则失去了多尺度分解的意义。经验公式：levels = int(log2(min(H, W))) - 3。以 512×512 为例，log2(512) = 9，所以常用 5 或 6 层。对 256×256 的医学图像，3~4 层就够。
还有一个容易被忽略的细节：第三层以上的拉普拉斯层尺寸很小（比如 64×64），这时直接做像素级选择意义不大，通常变成区域决策——两个图在这个尺度上的能量都极为稀疏，直接取平均即可。实际工程里，很多实现从levels - 2层开始全部用平均融合，可以避免顶层噪声干扰。层数对上采样步数也有影响：每多一层，pyrUp的累积插值误差就多一次迭代。非整数倍缩放时，建议层数不要超过 6。
5.2 高斯核尺寸：默认 5×5 够用的情况与必须改的情况
cv2.pyrDown和cv2.pyrUp使用内置 5×5 高斯核，这个核的 σ≈1，对大多数图像是合理默认值。但有三种情况要改用自定义核：
- 图像噪声较强：直接用默认核会让噪声在下采样时被部分保留，融合后高频层出现斑点。先用
cv2.GaussianBlur(img, (5, 5), 1.5)预处理，再进金字塔。 - 图像边缘极锐利（如屏幕截图、UI 界面）：5×5 核会造成边缘过冲，重构后出现黑白 halo。改用 3×3 核的
cv2.pyrDown变体，但 OpenCV 不支持直接传核函数给pyrDown。 - 医学图像要求各向同性：默认核近似各向同性，但如果你使用非方形像素（如某些超声图像），需要先重采样为各向同性，否则融合方向性偏差会污染结果。
另一种做法是完全抛弃pyrDown，自定义降采样函数：
def gaussian_downsample(img, kernel_size=5, sigma=1.0): blur = cv2.GaussianBlur(img, (kernel_size, kernel_size), sigma) return blur[::2, ::2] # 隔行采样对应的上采样需要先插值再卷积，插值方法对结果影响比降采样更大，推荐cv2.INTER_CUBIC而不是INTER_LINEAR，后者在边缘处会有轻微锯齿。
5.3 高频层决策的阈值化与软切换
在全部使用mode='max'时，平坦区域会出现一个视觉瑕疵——两图在该处都有很弱的噪声，但 max 选择随机切换，导致融合结果在时间序列或相邻帧之间闪烁。解决办法是引入死区（dead zone）：只有当两系数绝对值之差超过某个阈值时才选大，否则直接平均。这个思想源于小波阈值降噪，应用在图像融合上非常有效。
def threshold_fuse(la, lb, tau=5.0): """阈值化融合：差值小于 tau 时取平均，否则取绝对值大者""" diff = la - lb mask = np.abs(diff) > tau avg = (la + lb) / 2.0 return np.where(mask, np.where(np.abs(la) > np.abs(lb), la, lb), avg)tau的取值以图像灰度级为单位：输入 0-255 图像，tau=5表示灰度差小于 5 的细节系数视为噪声；tau太大（超过 30）会导致融合结果明显“偏软”，边缘过渡区被平均化。这个参数值得做一次网格搜索，是拉普拉斯金字塔融合中最需要根据数据调整的参数之一。
5.4 OpenCV 与 NumPy 的边界处理差异
最后一个埋坑点：cv2.filter2D默认边界处理是BORDER_REFLECT_101，而np.where没有边界概念；但cv2.pyrDown和cv2.pyrUp内部使用BORDER_REFLECT_101（OpenCV 的默认值）。如果你自己实现降采样时用了BORDER_CONSTANT（补零），金字塔重构后边缘一圈会发暗。
排查方法很简单：把融合空图和原图相减，查看边缘误差分布。若误差集中在边界 2~3 像素内，基本都是边界填充策略不一致导致的。统一使用cv2.BORDER_REFLECT即可。
6. 融合质量的量化评估与残差可视化验证法
6.1 从主观看效果到多光谱感知差距
人眼判断融合质量只能解决“有没有明显伪影”的问题，无法量化“结构保留度”。较实用的指标组合是：结构相似性（SSIM）、视觉信息保真度（VIF）以及梯度保留度 Qabf。SSIM 衡量亮度/对比度/结构三维相似性，VIF 更贴近人眼感知，Qabf 专门评估边缘方向的保留程度。
from skimage.metrics import structural_similarity as ssim import numpy as np def compute_ssim_with_inputs(result, img_a, img_b): ssim_a = ssim(result, img_a, multichannel=True) ssim_b = ssim(result, img_b, multichannel=True) return ssim_a, ssim_b ssim_a, ssim_b = compute_ssim_with_inputs(result, img_a, img_b) print(f"Fused vs A: {ssim_a:.4f}, Fused vs B: {ssim_b:.4f}")两个 SSIM 中较大者代表“融合结果更像哪张源图”。在多聚焦场景里，期望是两边的数值都比较高（融合图同时包含两张图的清晰区域）；如果ssim_a接近 1 而ssim_b远低于 0.8，说明融合基本没用到 B 的特征，策略失效。VIF 的实现可以用piq或sewar库，这里不做展开。
6.2 残差热力图：定位融合伪影的空间位置
SSIM 是一个标量，无法告诉你“哪块区域融合失败了”。残差图可视化是更直接的手段：
import matplotlib.pyplot as plt residual_a = np.abs(result.astype(np.float32) - img_a.astype(np.float32)) # 对残差做高斯模糊，突出区域性差异 residual_smooth = cv2.GaussianBlur(residual_a, (15, 15), 5) plt.subplot(1, 2, 1) plt.imshow(cv2.cvtColor(residual_smooth, cv2.COLOR_BGR2GRAY), cmap='jet') plt.title('Residual vs A') plt.colorbar()如果热力图中高亮区域集中在图像边缘一两像素宽的位置，说明是配准误差或金字塔重构振铃；如果高亮区域出现在大片平坦区域，说明融合策略的频率分配不合理（通常是低频层被错误使用max导致），这时应该把高亮区域对应金字塔层改为平均策略。
6.3 一个工程技巧：把金字塔各层同步可视化
调试时只看到最终融合结果太抽象，更难定位问题。建议把拉普拉斯金字塔每层归一化后拼成一张大图：
def draw_laplacian_layers(layers): """把拉普拉斯金字塔各层归一化后横向拼接""" canvas = [] for layer in layers: layer_norm = (layer - layer.min()) / (layer.max() - layer.min() + 1e-8) layer_norm = (layer_norm * 255).astype(np.uint8) canvas.append(cv2.resize(layer_norm, (256, 256))) return np.hstack(canvas)把融合前和融合后的金字塔各层各拼一张，并排对比。正常情况：高频层应该保留源图像各自最清晰的边缘，过渡区域的系数强度应逐渐变化；若发现某一层融合后能量明显低于任一输入，说明该层决策规则拿错了。利用这种分层的频域可视化，再配合 SSIM 数值辅助判断，比一次次重新生成融合图再肉眼看要高效得多。这就是拉普拉斯金字塔图像融合这套流程里最值得沉淀的调试路径。
本文还有配套的精品资源，点击获取
