基于PCNN和NSCT的图像融合MATLAB实现

基于脉冲耦合神经网络(PCNN)和非下采样轮廓波变换(NSCT)的图像融合MATLAB实现。该代码包含了NSCT分解与重构、PCNN模型实现以及融合规则设计。

function fused_image = PCNN_NSCT_Fusion(image1, image2)
    % 基于PCNN和NSCT的图像融合
    % 输入: image1, image2 - 待融合的两幅图像(灰度图像)
    % 输出: fused_image - 融合后的图像
    
    %% 参数设置
    params.levels = 4;               % NSCT分解层数
    params.directions = [4, 4, 8, 8]; % 各层方向数
    params.pcnn_beta = 0.2;           % PCNN链接强度系数
    params.pcnn_alphaL = 0.0693;      % PCNN链接输入衰减系数
    params.pcnn_alphaTheta = 0.2;     % PCNN阈值衰减系数
    params.pcnn_VL = 1.0;             % PCNN链接输入幅度
    params.pcnn_VTheta = 20;          % PCNN阈值幅度
    params.pcnn_iter = 10;            % PCNN迭代次数
    
    %% 图像预处理
    % 转换为灰度图像(如果输入是彩色图像)
    if size(image1, 3) == 3
        image1 = rgb2gray(image1);
    end
    if size(image2, 3) == 3
        image2 = rgb2gray(image2);
    end
    
    % 转换为double类型并归一化
    image1 = im2double(image1);
    image2 = im2double(image2);
    
    % 调整图像大小一致
    [rows, cols] = size(image1);
    image2 = imresize(image2, [rows, cols]);
    
    %% NSCT分解
    [low1, high1] = NSCT_Decomposition(image1, params.levels, params.directions);
    [low2, high2] = NSCT_Decomposition(image2, params.levels, params.directions);
    
    %% 低频子带融合
    low_fused = Fuse_Low_Frequency(low1, low2);
    
    %% 高频子带融合
    high_fused = cell(size(high1));
    for lev = 1:params.levels
        num_dir = length(high1{lev});
        high_fused{lev} = cell(1, num_dir);
        for dir = 1:num_dir
            band1 = high1{lev}{dir};
            band2 = high2{lev}{dir};
            high_fused{lev}{dir} = Fuse_High_Frequency(band1, band2, params);
        end
    end
    
    %% NSCT重构
    fused_image = NSCT_Reconstruction(low_fused, high_fused, params.directions);
    
    %% 后处理
    fused_image = imadjust(fused_image);  % 对比度调整
    fused_image = mat2gray(fused_image);  % 归一化到[0,1]
end

%% NSCT分解函数
function [low, high] = NSCT_Decomposition(img, levels, directions)
    % 非下采样轮廓波变换分解
    low = img;
    high = cell(levels, 1);
    
    for lev = 1:levels
        % 非下采样金字塔分解
        [lo, hi] = NSP_Decomp(low);
        low = lo;
        
        % 非下采样方向滤波器组分解
        dir_bands = cell(1, directions(lev));
        for d = 1:directions(lev)
            dir_bands{d} = NSDFB_Decomp(hi, d);
        end
        high{lev} = dir_bands;
    end
end

%% NSCT重构函数
function img = NSCT_Reconstruction(low, high, directions)
    % 非下采样轮廓波变换重构
    for lev = length(high):-1:1
        % 非下采样方向滤波器组重构
        hi = zeros(size(low));
        num_dir = length(high{lev});
        for d = 1:num_dir
            hi = hi + NSDFB_Recon(high{lev}{d}, d);
        end
        
        % 非下采样金字塔重构
        low = NSP_Recon(low, hi);
    end
    img = low;
end

%% 低频子带融合函数
function low_fused = Fuse_Low_Frequency(low1, low2)
    % 低频子带融合 - 加权平均法
    energy1 = sum(abs(low1(:)));
    energy2 = sum(abs(low2(:)));
    
    weight1 = energy1 / (energy1 + energy2 + eps);
    weight2 = energy2 / (energy1 + energy2 + eps);
    
    low_fused = weight1 * low1 + weight2 * low2;
end

%% 高频子带融合函数 (PCNN)
function band_fused = Fuse_High_Frequency(band1, band2, params)
    % 高频子带融合 - PCNN方法
    [rows, cols] = size(band1);
    
    % 初始化PCNN参数
    F = cat(3, abs(band1), abs(band2));  % 双通道输入
    Y = zeros(rows, cols, 2);             % 神经元输出
    U = zeros(rows, cols, 2);             % 内部状态
    L = zeros(rows, cols, 2);             % 链接输入
    Theta = ones(rows, cols, 2);          % 动态阈值
    
    % PCNN迭代
    for iter = 1:params.pcnn_iter
        for ch = 1:2
            % 链接输入更新
            L(:,:,ch) = params.pcnn_VL * conv2(Y(:,:,mod(ch,2)+1), ones(3)/9, 'same');
            
            % 内部状态更新
            U(:,:,ch) = F(:,:,ch) .* (1 + params.pcnn_beta * L(:,:,ch));
            
            % 脉冲输出更新
            Y(:,:,ch) = U(:,:,ch) > Theta(:,:,ch);
            
            % 阈值更新
            Theta(:,:,ch) = exp(-params.pcnn_alphaTheta) * Theta(:,:,ch) + ...
                            params.pcnn_VTheta * Y(:,:,ch);
        end
    end
    
    % 计算点火次数
    ignition_count = sum(Y, 3);
    
    % 融合规则:选择点火次数多的系数
    mask = ignition_count(:,:,1) > ignition_count(:,:,2);
    band_fused = band1;
    band_fused(~mask) = band2(~mask);
end

%% 非下采样金字塔分解 (简化实现)
function [lo, hi] = NSP_Decomp(img)
    % 使用高斯滤波器进行分解
    h = fspecial('gaussian', [5 5], 1);
    lo = imfilter(img, h, 'replicate');
    hi = img - lo;
end

%% 非下采样金字塔重构 (简化实现)
function img = NSP_Recon(lo, hi)
    % 金字塔重构
    img = lo + hi;
end

%% 非下采样方向滤波器组分解 (简化实现)
function band = NSDFB_Decomp(hi, direction)
    % 使用方向滤波器进行分解
    switch direction
        case 1 % 0度方向
            filt = [1 0 -1; 2 0 -2; 1 0 -1]/4;
        case 2 % 45度方向
            filt = [0 1 2; -1 0 1; -2 -1 0]/4;
        case 3 % 90度方向
            filt = [-1 2 -1; 0 0 0; 1 -2 1]/4;
        case 4 % 135度方向
            filt = [-2 -1 0; -1 0 1; 0 1 2]/4;
        otherwise % 其他方向使用随机滤波器
            filt = randn(3,3)/9;
    end
    band = imfilter(hi, filt, 'replicate');
end

%% 非下采样方向滤波器组重构 (简化实现)
function hi = NSDFB_Recon(band, direction)
    % 方向滤波器组重构
    hi = band; % 简化处理
end

%% 主函数 - 图像融合演示
function Demo_PCNN_NSCT_Fusion()
    % 读取示例图像
    img1 = imread('cameraman.tif');
    img2 = imread('rice.png');
    
    % 执行融合
    fused_img = PCNN_NSCT_Fusion(img1, img2);
    
    % 显示结果
    figure;
    subplot(2,2,1); imshow(img1); title('源图像1');
    subplot(2,2,2); imshow(img2); title('源图像2');
    subplot(2,2,3); imshow(fused_img); title('融合图像');
    
    % 计算性能指标
    metrics = Evaluate_Fusion(fused_img, img1, img2);
    subplot(2,2,4); text(0.1, 0.5, sprintf(...
        '信息熵: %.4f\n互信息: %.4f\nSSIM: %.4f\nQ^{AB/F}: %.4f', ...
        metrics.EN, metrics.MI, metrics.SSIM, metrics.Qabf));
    axis off; title('性能指标');
end

%% 融合性能评估函数
function metrics = Evaluate_Fusion(fused, img1, img2)
    % 计算融合图像的性能指标
    metrics = struct();
    
    % 信息熵 (EN)
    metrics.EN = entropy(fused);
    
    % 互信息 (MI)
    metrics.MI = immutualinfo(fused, img1) + immutualinfo(fused, img2);
    
    % 结构相似性 (SSIM)
    metrics.SSIM = ssim(fused, img1) + ssim(fused, img2);
    metrics.SSIM = metrics.SSIM / 2; % 平均值
    
    % 边缘保持度 (Q^{AB/F})
    metrics.Qabf = 0;
    % 这里简化处理,实际实现需要更复杂的计算
end

%% 辅助函数 - 显示图像
function Show_Images(img1, img2, fused)
    figure;
    subplot(1,3,1); imshow(img1); title('源图像1');
    subplot(1,3,2); imshow(img2); title('源图像2');
    subplot(1,3,3); imshow(fused); title('融合图像');
end

使用示例

% 读取图像
image1 = imread('medical1.jpg');
image2 = imread('medical2.jpg');

% 执行融合
fused_image = PCNN_NSCT_Fusion(image1, image2);

% 显示结果
figure;
subplot(1,3,1); imshow(image1); title('源图像1');
subplot(1,3,2); imshow(image2); title('源图像2');
subplot(1,3,3); imshow(fused_image); title('融合图像');

% 保存结果
imwrite(fused_image, 'fused_result.jpg');

算法说明

1. NSCT分解与重构

  • 非下采样金字塔(NSP):实现多尺度分解,捕获不同尺度的图像特征

  • 非下采样方向滤波器组(NSDFB):实现多方向分解,捕获不同方向的边缘信息

  • 优势:平移不变性、良好的频率选择性、多尺度多方向分析能力

2. PCNN模型

  • 脉冲耦合神经网络:模拟视觉皮层神经元的同步脉冲发放特性

  • 核心方程:

    F_ij[n] = S_ij                   // 反馈输入
    L_ij[n] = e^{-α_L} L_ij[n-1] + V_L Σ M_kl Y_kl[n-1]  // 链接输入
    U_ij[n] = F_ij[n](1 + βL_ij[n])   // 内部状态
    Y_ij[n] = {1 if U_ij[n] > θ_ij[n], else 0}  // 输出脉冲
    θ_ij[n] = e^{-α_θ} θ_ij[n-1] + V_θ Y_ij[n-1]  // 动态阈值
    

3. 融合规则

  • 低频子带:加权平均法(能量自适应)

    L_F(i,j) = w_A(i,j)L_A(i,j) + w_B(i,j)L_B(i,j)
    w_A = E_A/(E_A+E_B), w_B = E_B/(E_A+E_B)
    E = ∑|L(i,j)| (子带能量)
    
  • 高频子带:PCNN点火次数决定

    • 将两幅源图像的高频系数作为PCNN输入

    • 计算每个系数的点火次数

    • 选择点火次数多的系数作为融合结果

参考代码 利用pcnn和nsct实现图像融合 www.youwenfan.com/contentcns/52378.html

应用场景

  1. 医学图像融合

    • 多模态医学图像(CT/MRI/PET)融合

    • 手术导航中的实时影像融合

    • 病理分析与诊断

  2. 遥感图像处理

    • 多光谱与全色图像融合

    • 多时相遥感图像变化检测

    • 灾害评估与监测

  3. 安防监控系统

    • 红外与可见光图像融合

    • 多摄像头视角融合

    • 低照度环境增强

  4. 计算机视觉

  • 多聚焦图像融合

  • 多曝光图像融合

  • 三维重建中的图像融合

posted @ 2026-03-10 16:33  风一直那个吹  阅读(2)  评论(0)    收藏  举报