当前位置：首页 > news >正文

MATLAB实现图像PCA降噪

news 2025/11/13 10:26:15

一、算法实现

function denoised_img = pca_denoise(img_path, block_size, n_components)% 读取图像并灰度化img = imread(img_path);if size(img,3)==3img = rgb2gray(img);endimg = im2double(img);[h,w] = size(img);% 分块参数pad = block_size//2;padded = padarray(img, [pad pad], 'replicate');% 初始化输出denoised = zeros(h,w);% 滑动窗口处理for i = 1:block_size:h-padfor j = 1:block_size:w-pad% 提取当前块block = padded(i:i+block_size-1, j:j+block_size-1);% 数据标准化block_centered = block - mean(block(:));% PCA降维[U,S,V] = svd(block_centered, 'econ');components = U(:,1:n_components);% 重构图像块recon_block = components * S(1:n_components,1:n_components) * components';recon_block = recon_block + mean(block(:));% 边界处理denoised(i:i+block_size-1, j:j+block_size-1) = recon_block;endend% 峰值信噪比计算psnr_value = psnr(denoised, img);denoised_img = uint8(denoised);fprintf('PSNR: %.2f dB\n', psnr_value);
end

二、参数优化

1. 分块大小选择

块尺寸	适用场景	计算效率
8x8	快速处理（实时性要求高）	高
16x16	平衡质量与效率	中
32x32	高分辨率图像细节保留	低

2. 主成分数量确定

% 自动选择保留95%方差的成分数
function k = auto_select_components(block)cov_mat = cov(block(:)');[~,S] = eig(cov_mat);S = diag(S);k = find(cumsum(S)/sum(S) >= 0.95, 1);
end

3. 预处理增强

% 添加高斯噪声模拟
noisy_img = imnoise(img, 'gaussian', 0, 0.05);% 中值滤波预处理
preprocessed = medfilt2(noisy_img, [3 3]);

三、完整处理流程示例

%% 加载测试图像
img = imread('lena.png');
img = rgb2gray(img);
img = im2double(img);%% 参数设置
block_size = 16;    % 分块尺寸
n_components = 20;  % 主成分数量%% 执行PCA降噪
denoised_img = pca_denoise(img, block_size, n_components);%% 结果可视化
figure;
subplot(1,2,1);
imshow(img);
title('原始图像');
subplot(1,2,2);
imshow(denoised_img);
title(sprintf('PCA降噪 (PSNR=%.2f dB)', psnr(denoised_img,img)));%% 保存结果
imwrite(denoised_img, 'denoised_lena.png');

四、优化改进

重叠分块处理：

function denoised = overlapped_pca(img, block_size, overlap)pad = overlap//2;padded = padarray(img, [pad pad], 'replicate');[h,w] = size(img);denoised = zeros(h,w);for i = 1:block_size:h-padfor j = 1:block_size:w-padblock = padded(i:i+block_size-1, j:j+block_size-1);% PCA处理...% 加权平均重叠区域weight = hamming(block_size);denoised(i-overlap:i-overlap+block_size, j-overlap:j-overlap+block_size) = ...denoised(i-overlap:i-overlap+block_size, j-overlap:j-overlap+block_size) + ...recon_block .* weight;endend
end

多尺度PCA：

% 先进行小波分解再进行PCA
coeffs = wavedec2(img, 2, 'db1');
denoised_coeffs = cell(size(coeffs));
for k = 1:length(coeffs)denoised_coeffs{k} = pca_denoise_block(coeffs{k});
end
denoised_img = waverec2(denoised_coeffs, 'db1');

残差增强：

residual = img - denoised_img;
enhanced = denoised_img + wiener2(residual, [5 5]);

五、效果评估脚本

function metrics = evaluate_denoising(original, denoised, noisy)metrics.PSNR = psnr(denoised, original);metrics.SSIM = ssim(denoised, original);metrics.MSE = mean((original(:) - denoised(:)).^2);metrics.NR = mean((denoised(:) - original(:)).^2);
end% 使用示例
original = imread('lena_clean.png');
denoised = imread('denoised_lena.png');
noisy = imread('noisy_lena.png');metrics = evaluate_denoising(original, denoised, noisy);
disp('=== 性能评估 ===');
for name = fieldnames(metrics)'fprintf('%s: %.4f\n', name{1}, metrics.(name{1}));
end