图像去噪序列——BM3D图像去噪模型实现
生活随笔
收集整理的這篇文章主要介紹了
图像去噪序列——BM3D图像去噪模型实现
小編覺得挺不錯的,現在分享給大家,幫大家做個參考.
1. BM3D模型簡介
BM3D模型是一個兩階段圖像去噪方法，主要包含兩個步驟:
(1) 在噪聲圖像上，利用局部區域搜索相似塊，并進行堆疊，在變換域(DCT域、FFT域）利用硬閾值去噪方法對堆疊的圖像塊進行去噪，獲得堆疊相似塊的估計值，最后，根據均值權重進行聚合；
(2) 通過步驟(1) 獲取初步估計的圖像，在初步估計的圖像上進行相似塊的聚合； 然后，利用維納協同濾波進行圖像去噪，從而，獲取最后的去噪結果
2. 模型實現(代碼參考網絡實現)：
% BM3D_Color_Demo
% BM3D 在彩色圖像上去噪
% Author: HSW
% Date: 2018-05-06
%
clc;
close all;
clear all;
img_org = imread('timg.png');
figure(1);
imshow(img_org);
title('原圖像');
% 加噪聲
sigma = 25;
img_noise = double(img_org)+sigma * randn(size(img_org));
figure;
imshow(img_noise / 255, []);
title('噪聲圖像');
img_denoise = BM3D_Color(img_noise, 0, sigma, 0, 1);
figure;
imshow(img_denoise / 255, []);
title('去噪圖像'); % BM3D_Gray_Demo
% BM3D 在灰度圖像上去噪
% Author: HSW
% Date: 2018-05-06
%
clc;
close all;
clear all;
img_org = imread('timg.png');
img_gray = rgb2gray(img_org);
figure(1);
imshow(img_gray);
title('原圖像');
% 加噪聲
sigma = 25;
img_noise = double(img_gray)+sigma * randn(size(img_gray));
figure;
imshow(img_noise / 255, []);
title('噪聲圖像');
img_denoise = BM3D_Gray(img_noise, 0, sigma, 1);
figure;
imshow(img_denoise / 255, []);
title('去噪圖像'); function img_denoise = BM3D_Color(img_noise, tran_mode, sigma, color_mode, isDisplay)
% BM3D實現去噪
% Inputs:
% img_noise: 噪聲圖像
% tran_mode: 變換方法： 默認值為0, tran_mode: = 0, fft; = 1, dct; = 2, dwt, = 3, db1
% sigma: 噪聲水平，默認值為10
% color_mode: 彩色圖像去噪時采用的顏色空間, 默認值為0, color_mode: = 0, YUV; = 1, YCbCr; = 2, OPP
% Ouputs:
% img_out: 去噪圖像
% 參考文獻：An Analysis and Implementation of the BM3D Image Denoising Method
% Inputs:
% img_in: 噪聲圖像，必須為矩形方陣
% tran_mode: = 0, FFT; = 1, DCT; = 2, DWT, = 3, db1
% Outputs:
% img_denoise: 去噪圖像
%
%
if ~exist('isDisplay', 'var')
isDisplay = 0;
end
if ~exist('color_mode', 'var')
color_mode = 0;
end
if ~exist('sigma', 'var')
sigma = 10;
end
if ~exist('tran_mode', 'var')
tran_mode = 0;
end
[row, col, dims] = size(img_noise);
img_trans = rgb2other(img_noise, color_mode);
% First Step 參數
kHard = 8; % 塊大小
pHard = 4; % 塊移動間隔
lambda_distHard = 0; % 求相似的距離時，變換后，收縮的閾值
nHard = 40; % 搜索窗口大小
NHard = 28; % 最多相似塊個數
tauHard = 5000; % 最大的相似距離for fft
% kaiser窗口的參數，實際上并沒有特別大的影響
beta=2;
Wwin2D = kaiser(kHard, beta) * kaiser(kHard, beta)';
% Second Step參數
kWien = kHard;
pWien = pHard;
lambda_distWien = lambda_distHard;
nWien = nHard;
NWien = NHard;
tauWien = tauHard;
sigma2 = sigma*sigma;
if tran_mode == 0
% FFT
lambda2d=400;
lambda1d=500;
lambda2d_wie=50;
lambda1d_wie=500;
elseif tran_mode == 1
% DCT
lambda2d=50;
lambda1d=80;
lambda2d_wie=20;
lambda1d_wie=60;
elseif tran_mode == 2
% DWT
lambda2d=50;
lambda1d=80;
lambda2d_wie=20;
lambda1d_wie=60;
end
fprintf('BM3D: First Stage Start...\n');
%block為原始圖像塊， tran_block為FFT變換且硬閾值截斷后的頻域系數(頻域， 計算距離的時候采用的是變換塊)
[block_ch1, tran_block_ch1, block2row_idx_ch1, block2col_idx_ch1] = im2block(img_trans(:,:,1), kHard, pHard, lambda_distHard, 0);
[block_ch2, tran_block_ch2, block2row_idx_ch2, block2col_idx_ch2] = im2block(img_trans(:,:,2), kHard, pHard, lambda_distHard, 0);
[block_ch3, tran_block_ch3, block2row_idx_ch3, block2col_idx_ch3] = im2block(img_trans(:,:,3), kHard, pHard, lambda_distHard, 0);
%bn_r和bn_c為行和列上的圖像塊個數
bn_r = floor((row - kHard) / pHard) + 1;
bn_c = floor((col - kHard) / pHard) + 1;
%基礎估計的圖像
img_basic_sum = zeros(row, col, 3);
img_basic_weight = zeros(row, col, 3);
%對每個塊遍歷
for i=1:bn_r
for j=1:bn_c
% 利用亮度通道進行相似塊搜索
[sim_blk_ch1, sim_num, sim_blk_idx] = search_similar_block(i, j, block_ch1, tran_block_ch1, floor(nHard/pHard), bn_r, bn_c, tauHard, NHard);
% 進行亮度通道處理
% 協同濾波： 公式(2)
tran3d_blk_shrink_ch1 = transform_3d(sim_blk_ch1, tran_mode, lambda2d, lambda1d);
tran3d_blk_shrink_ch2 = transform_3d(block_ch2(:,:,sim_blk_idx), tran_mode, lambda2d, lambda1d);
tran3d_blk_shrink_ch3 = transform_3d(block_ch3(:,:,sim_blk_idx), tran_mode, lambda2d, lambda1d);
% 聚合： 公式(3)中的說明
NHard_P_ch1 = nnz(tran3d_blk_shrink_ch1);
NHard_P_ch2 = nnz(tran3d_blk_shrink_ch2);
NHard_P_ch3 = nnz(tran3d_blk_shrink_ch3);
if NHard_P_ch1 > 1
wHard_P_ch1 = 1 / NHard_P_ch1;
else
wHard_P_ch1 = 1;
end
if NHard_P_ch2 > 1
wHard_P_ch2 = 1 / NHard_P_ch2;
else
wHard_P_ch2 = 1;
end
if NHard_P_ch3 > 1
wHard_P_ch3 = 1 / NHard_P_ch3;
else
wHard_P_ch3 = 1;
end
blk_est_ch1 = inv_transform_3d(tran3d_blk_shrink_ch1,tran_mode);
blk_est_ch1 = real(blk_est_ch1);
blk_est_ch2 = inv_transform_3d(tran3d_blk_shrink_ch2, tran_mode);
blk_est_ch2 = real(blk_est_ch2);
blk_est_ch3 = inv_transform_3d(tran3d_blk_shrink_ch3, tran_mode);
blk_est_ch3 = real(blk_est_ch3);
% 公式(3): 對亮度通道，即第1個通道
for k=1:sim_num
idx = sim_blk_idx(k);
ir = block2row_idx_ch1(idx);
jr = block2col_idx_ch1(idx);
img_basic_sum(ir:ir+kHard-1, jr:jr+kHard-1, 1) = img_basic_sum(ir:ir+kHard-1, jr:jr+kHard-1, 1) + wHard_P_ch1 * blk_est_ch1(:, :, k);
img_basic_weight(ir:ir+kHard-1, jr:jr+kHard-1, 1) = img_basic_weight(ir:ir+kHard-1, jr:jr+kHard-1, 1) + wHard_P_ch1;
img_basic_sum(ir:ir+kHard-1, jr:jr+kHard-1, 2) = img_basic_sum(ir:ir+kHard-1, jr:jr+kHard-1, 2) + wHard_P_ch2 * blk_est_ch2(:, :, k);
img_basic_weight(ir:ir+kHard-1, jr:jr+kHard-1, 2) = img_basic_weight(ir:ir+kHard-1, jr:jr+kHard-1, 2) + wHard_P_ch2;
img_basic_sum(ir:ir+kHard-1, jr:jr+kHard-1, 3) = img_basic_sum(ir:ir+kHard-1, jr:jr+kHard-1, 3) + wHard_P_ch3 * blk_est_ch3(:, :, k);
img_basic_weight(ir:ir+kHard-1, jr:jr+kHard-1, 3) = img_basic_weight(ir:ir+kHard-1, jr:jr+kHard-1, 3) + wHard_P_ch3;
end
end
end
img_basic = img_basic_sum ./ img_basic_weight;
if isDisplay
figure;
img_rgb = other2rgb(img_basic, color_mode);
imshow(img_rgb / 255.0 ,[]);
title('BM3D:Fist Stage Result');
end
fprintf('BM3D: First Stage End...\n');
fprintf('BM3D: Second Stage Start...\n');
[block_basic_ch1,tran_block_basic_ch1,block2row_idx_basic_ch1,block2col_idx_basic_ch1] = im2block(img_basic(:, :, 1), kWien, pWien, lambda_distWien, 0);
[block_basic_ch2,tran_block_basic_ch2,block2row_idx_basic_ch3,block2col_idx_basic_ch2] = im2block(img_basic(:, :, 2), kWien, pWien, lambda_distWien, 0);
[block_basic_ch3,tran_block_basic_ch3,block2row_idx_basic_ch3,block2col_idx_basic_ch3] = im2block(img_basic(:, :, 3), kWien, pWien, lambda_distWien, 0);
bn_r = floor((row - kWien) / pWien) + 1;
bn_c = floor((col - kWien) / pWien) + 1;
img_wien_sum = zeros(row, col, 3);
img_wien_weight = zeros(row, col, 3);
for i=1:1:bn_r
for j=1:1:bn_c
% 公式(5), 利用亮度進行相似性搜索
[sim_blk_basic_ch1, sim_num, sim_blk_basic_idx] = search_similar_block(i, j, block_basic_ch1, tran_block_basic_ch1, floor(nWien/pWien), bn_r, bn_c, tauWien, NWien);
% 公式(6)
tran3d_blk_basic_ch1 = transform_3d(sim_blk_basic_ch1, tran_mode, lambda2d_wie, lambda1d_wie);
tran3d_blk_basic_ch2 = transform_3d(block_basic_ch2(:, :, sim_blk_basic_idx), tran_mode, lambda2d_wie, lambda1d_wie);
tran3d_blk_basic_ch3 = transform_3d(block_basic_ch3(:, :, sim_blk_basic_idx), tran_mode, lambda2d_wie, lambda1d_wie);
omega_P_ch1 = (tran3d_blk_basic_ch1.^2) ./ ((tran3d_blk_basic_ch1.^2) + sigma2);
omega_P_ch2 = (tran3d_blk_basic_ch2.^2) ./ ((tran3d_blk_basic_ch2.^2) + sigma2);
omega_P_ch3 = (tran3d_blk_basic_ch3.^2) ./ ((tran3d_blk_basic_ch3.^2) + sigma2);
% 公式(7)
tran3d_blk_ch1 = transform_3d(block_ch1(:, :, sim_blk_basic_idx), tran_mode, lambda2d_wie, lambda1d_wie);
tran3d_blk_ch2 = transform_3d(block_ch2(:, :, sim_blk_basic_idx), tran_mode, lambda2d_wie, lambda1d_wie);
tran3d_blk_ch3 = transform_3d(block_ch3(:, :, sim_blk_basic_idx), tran_mode, lambda2d_wie, lambda1d_wie);
blk_est_ch1 = inv_transform_3d(omega_P_ch1 .* tran3d_blk_ch1, tran_mode);
blk_est_ch2 = inv_transform_3d(omega_P_ch2 .* tran3d_blk_ch2, tran_mode);
blk_est_ch3 = inv_transform_3d(omega_P_ch3 .* tran3d_blk_ch3, tran_mode);
blk_est_ch1 = real(blk_est_ch1);
blk_est_ch2 = real(blk_est_ch2);
blk_est_ch3 = real(blk_est_ch3);
NWien_P_ch1 = nnz(omega_P_ch1);
NWien_P_ch2 = nnz(omega_P_ch2);
NWien_P_ch3 = nnz(omega_P_ch3);
if NWien_P_ch1 > 1
wWien_P_ch1 = 1 / (NWien_P_ch1);
else
wWien_P_ch1 = 1;
end
if NWien_P_ch2 > 1
wWien_P_ch2 = 1/(NWien_P_ch2);
else
wWien_P_ch2 = 1;
end
if NWien_P_ch3 > 1
wWien_P_ch3 = 1 / (NWien_P_ch3);
else
wWien_P_ch3 = 1;
end
% 公式(8)
for k=1:sim_num
idx=sim_blk_basic_idx(k);
ir=block2row_idx_basic_ch1(idx);
jr=block2col_idx_basic_ch1(idx);
img_wien_sum(ir:ir+kWien-1, jr:jr+kWien-1, 1) = img_wien_sum(ir:ir+kWien-1, jr:jr+kWien-1, 1) + wWien_P_ch1 * blk_est_ch1(:, :, k);
img_wien_weight(ir:ir+kWien-1, jr:jr+kWien-1, 1) = img_wien_weight(ir:ir+kWien-1, jr:jr+kWien-1, 1) + wWien_P_ch1;
img_wien_sum(ir:ir+kWien-1, jr:jr+kWien-1, 2) = img_wien_sum(ir:ir+kWien-1, jr:jr+kWien-1, 2) + wWien_P_ch2 * blk_est_ch2(:, :, k);
img_wien_weight(ir:ir+kWien-1, jr:jr+kWien-1, 2) = img_wien_weight(ir:ir+kWien-1, jr:jr+kWien-1, 2) + wWien_P_ch2;
img_wien_sum(ir:ir+kWien-1, jr:jr+kWien-1, 3) = img_wien_sum(ir:ir+kWien-1, jr:jr+kWien-1, 3) + wWien_P_ch3 * blk_est_ch3(:, :, k);
img_wien_weight(ir:ir+kWien-1, jr:jr+kWien-1, 3) = img_wien_weight(ir:ir+kWien-1, jr:jr+kWien-1, 3) + wWien_P_ch3;
end
end
end
img_other = img_wien_sum ./ img_wien_weight;
img_denoise = other2rgb(img_other, color_mode);
fprintf('BM3D: Second Stage End\n');
function img_denoise = BM3D_Gray(img_noise, tran_mode, sigma, isDisplay)
% 參考文獻：An Analysis and Implementation of the BM3D Image Denoising Method
% Inputs:
% img_noise: 灰度噪聲圖像，必須為矩形方陣
% tran_mode: = 0, fft; = 1, dct; = 2, dwt, = 3, db1
% Outputs:
% img_denoise: 去噪圖像
%
if ~exist('tran_mode', 'var')
tran_mode = 0;
end
if ~exist('sigma', 'var')
sigma = 10;
end
if ~exist('isDisplay', 'var')
isDisplay = 0;
end
[row,col] = size(img_noise);
% First Step 參數
kHard = 8; % 塊大小
pHard = 4; % 塊移動間隔
lambda_distHard = 0; % 求相似的距離時，變換后，收縮的閾值
nHard = 40; % 搜索窗口大小
NHard = 28; % 最多相似塊個數
tauHard = 5000; % 最大的相似距離for fft
% kaiser窗口的參數，實際上并沒有特別大的影響
beta=2;
Wwin2D = kaiser(kHard, beta) * kaiser(kHard, beta)';
% Second Step參數
kWien = kHard;
pWien = pHard;
lambda_distWien = lambda_distHard;
nWien = nHard;
NWien = NHard;
tauWien = tauHard;
sigma2 = sigma*sigma;
if(tran_mode==0) %fft
lambda2d=400;
lambda1d=500;
lambda2d_wie=50;
lambda1d_wie=500;
elseif(tran_mode == 1) %dct
lambda2d=50;
lambda1d=80;
lambda2d_wie=20;
lambda1d_wie=60;
elseif(tran_mode == 2) %dwt
lambda2d=50;
lambda1d=80;
lambda2d_wie=20;
lambda1d_wie=60;
end
%block為原始圖像塊， tran_block為FFT變換且硬閾值截斷后的頻域系數(頻域， 計算距離的時候采用的是變換塊)
[block,tran_block,block2row_idx,block2col_idx]=im2block(img_noise,kHard,pHard,lambda_distHard,0);
%bn_r和bn_c為行和列上的圖像塊個數
bn_r=floor((row-kHard)/pHard)+1;
bn_c=floor((col-kHard)/pHard)+1;
%基礎估計的圖像
img_basic_sum=zeros(row,col);
img_basic_weight=zeros(row,col);
%basic處理
fprintf('BM3D: First Stage Start...\n');
%對每個塊遍歷
for i=1:bn_r
for j=1:bn_c
[sim_blk,sim_num,sim_blk_idx]=search_similar_block(i,j,block,tran_block,floor(nHard/pHard),bn_r,bn_c,tauHard,NHard);
% 協同濾波： 公式(2)
tran3d_blk_shrink=transform_3d(sim_blk,tran_mode,lambda2d,lambda1d);
% 聚合： 公式(3)中的說明
NHard_P=nnz(tran3d_blk_shrink);
if(NHard_P >1)
wHard_P=1/NHard_P;
else
wHard_P=1;
end
blk_est =inv_transform_3d(tran3d_blk_shrink,tran_mode);
blk_est=real(blk_est);
% 公式(3)
for k=1:sim_num
idx=sim_blk_idx(k);
ir=block2row_idx(idx);
jr=block2col_idx(idx);
img_basic_sum(ir:ir+kHard-1,jr:jr+kHard-1) = img_basic_sum(ir:ir+kHard-1,jr:jr+kHard-1) + wHard_P*blk_est(:,:,k);
img_basic_weight(ir:ir+kHard-1,jr:jr+kHard-1) = img_basic_weight(ir:ir+kHard-1,jr:jr+kHard-1) + wHard_P;
end
end
end
fprintf('BM3D: First Stage End...\n');
img_basic=img_basic_sum./img_basic_weight;
if isDisplay
figure;
imshow(img_basic,[]);
title('BM3D:Fist Stage Result');
end
[block_basic,tran_block_basic,block2row_idx_basic,block2col_idx_basic] = im2block(img_basic,kWien,pWien,lambda_distWien,0);
bn_r=floor((row-kWien)/pWien)+1;
bn_c=floor((col-kWien)/pWien)+1;
img_wien_sum=zeros(row,col);
img_wien_weight=zeros(row,col);
fprintf('BM3D: Second Stage Start...\n');
for i=1:1:bn_r
for j=1:1:bn_c
% 公式(5)
[sim_blk_basic,sim_num,sim_blk_basic_idx] = search_similar_block(i,j,block_basic,tran_block_basic,floor(nWien/pWien),bn_r,bn_c,tauWien,NWien);
% 公式(6)
tran3d_blk_basic = transform_3d(sim_blk_basic,tran_mode,lambda2d_wie,lambda1d_wie);
omega_P=(tran3d_blk_basic.^2)./((tran3d_blk_basic.^2)+sigma2);
% 公式(7)
tran3d_blk = transform_3d(block(:,:,sim_blk_basic_idx),tran_mode,lambda2d_wie,lambda1d_wie);
blk_est=inv_transform_3d(omega_P.*tran3d_blk,tran_mode);
blk_est=real(blk_est);
NWien_P=nnz(omega_P);
if(NWien_P >1)
wWien_P=1/(NWien_P);
else
wWien_P=1;
end
% 公式(8)
for k=1:sim_num
idx=sim_blk_basic_idx(k);
ir=block2row_idx_basic(idx);
jr=block2col_idx_basic(idx);
img_wien_sum(ir:ir+kWien-1,jr:jr+kWien-1) = img_wien_sum(ir:ir+kWien-1,jr:jr+kWien-1) + wWien_P*blk_est(:,:,k);
img_wien_weight(ir:ir+kWien-1,jr:jr+kWien-1) = img_wien_weight(ir:ir+kWien-1,jr:jr+kWien-1) + wWien_P;
end
end
end
fprintf('BM3D: Second Stage End\n');
img_denoise = img_wien_sum./img_wien_weight;
function [block,transform_block,block2row_idx,block2col_idx] =im2block(img,k,p,lambda2D,delta)
% 實現圖像分塊
% Inputs：
% k: 塊大小
% p: 塊移動步長
% lambda_2D: 收縮閾值
% delta: 收縮閾值
% Outputs:
% block: 返回的塊
% transform_block: 變換后的塊
% block2row_idx: 塊索引與圖像塊的左上角行坐標對應關系
% block2col_idx: 塊索引與圖像塊的左上角列坐標對應關系
%
[row,col] = size(img);
% 頻域去噪中的硬閾值，實際上原文中，對于噪聲方差小于40時thres = 0, 具體見公式(1)的說明第2點(即距離計算)
thres = lambda2D*delta*sqrt(2*log(row*col));
% r_num 和 c_num分別表示行和列上可以采集的塊的數目
r_num = floor((row-k)/p)+1;
c_num = floor((col-k)/p)+1;
block = zeros(k,k,r_num*c_num);
block2row_idx = [];
block2col_idx = [];
cnt = 1;
for i = 0:r_num-1
rs = 1+i*p;
for j = 0:c_num-1
cs = 1+j*p;
block(:,:,cnt) = img(rs:rs+k-1,cs:cs+k-1);
block2row_idx(cnt) = rs;
block2col_idx(cnt) = cs;
tr_b = fft2(block(:,:,cnt));
idx = find(abs(tr_b)<thres);
tr_b(idx) = 0;
transform_block(:,:,cnt) = tr_b;
cnt = cnt+1;
end
end
endfunction [blk_est]=inv_transform_3d(blk_tran3d,tran_mode)
% 3D 逆變換
% Inputs:
% blk_tran3d: 在頻域中，硬閾值濾波的圖像塊
% tran_mode: 變換方法
% Outputs:
% blk_est:
%
global blk_tran1d_s;
global blk_2d_s;
[m,n,blk_num]=size(blk_tran3d);
blk_invtran1d=zeros(m,n,blk_num);
blk_est=zeros(m,n,blk_num);
if(tran_mode==0) %fft
for i=1:1:m
for j=1:1:n
blk_invtran1d(i,j,:)=ifft(blk_tran3d(i,j,:));
end
end
for i=1:1:blk_num
blk_est(:,:,i)=ifft2(blk_invtran1d(:,:,i));
end
elseif(tran_mode==1) %dct
for i=1:1:m
for j=1:1:n
blk_invtran1d(i,j,:)=idct(blk_tran3d(i,j,:));
end
end
for i=1:1:blk_num
blk_est(:,:,i)=idct2(blk_invtran1d(:,:,i));
end
elseif(tran_mode==2) %dwt
blk_num=length(blk_2d_s);
blk_c=waverec2(blk_tran3d,blk_tran1d_s,'haar');
blk_est=[];
for i=1:1:blk_num
blk_est(:,:,i)=waverec2(blk_c(:,i),blk_2d_s{i},'Bior1.5');
end
else
error('tran_mode error');
end
endfunction img_trans = other2rgb(img_in, color_mode)
% 將RGB顏色空間轉為其他顏色空間
% Inputs:
% img_in: RGB顏色空間圖像
% color_mode: 彩色圖像去噪時采用的顏色空間, 默認值為0, color_mode: = 0, YUV; = 1, YCbCr; = 2, OPP
% Outputs:
% img_trans: 其他顏色空間
%
% Author: HSW
% Date: 2018-05-06
img_trans = zeros(size(img_in));
[row, col, dims] = size(img_in);
if color_mode == 0
color_tran = [0.30, 0.59, 0.11; -0.15, -0.29, 0.44; 0.61, -0.51, -0.10];
color_tran_inv = inv(color_tran);
for i = 1:row
for j = 1:col
other = [img_in(i, j, 1); img_in(i, j, 2); img_in(i, j, 3)];
img_trans(i, j, :) = color_tran_inv * other;
end
end
elseif color_mode == 1
color_tran = [0.30, 0.59, 0.11; -0.17, -0.33, 0.50; 0.50, -0.42, -0.08];
color_tran_inv = inv(color_tran);
for i = 1:row
for j = 1:col
other = [img_in(i, j, 1); img_in(i, j, 2); img_in(i, j, 3)];
img_trans(i, j, :) = color_tran_inv * other;
end
end
elseif color_mode == 2
color_tran = [1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0; 1.0 / 2.0, 0, -1.0 / 2.0; 1.0 / 4.0, -1.0 / 2.0, 1.0 / 4.0];
color_tran_inv = inv(color_tran);
for i = 1:row
for j = 1:col
other = [img_in(i, j, 1); img_in(i, j, 2); img_in(i, j, 3)];
img_trans(i, j, :) = color_tran_inv * other;
end
end
end
endfunction img_trans = rgb2other(img_in, color_mode)
% 將RGB顏色空間轉為其他顏色空間
% Inputs:
% img_in: RGB顏色空間圖像
% color_mode: 彩色圖像去噪時采用的顏色空間, 默認值為0, color_mode: = 0, YUV; = 1, YCbCr; = 2, OPP
% Outputs:
% img_trans: 其他顏色空間
%
% Author: HSW
% Date: 2018-05-06
img_trans = zeros(size(img_in));
[row, col, dims] = size(img_in);
if color_mode == 0
color_tran = [0.30, 0.59, 0.11; -0.15, -0.29, 0.44; 0.61, -0.51, -0.10];
for i = 1:row
for j = 1:col
rgb = [img_in(i, j, 1); img_in(i, j, 2); img_in(i, j, 3)];
img_trans(i, j, :) = (color_tran * rgb)';
end
end
elseif color_mode == 1
color_tran = [0.30, 0.59, 0.11; -0.17, -0.33, 0.50; 0.50, -0.42, -0.08];
for i = 1:row
for j = 1:col
rgb = [img_in(i, j, 1); img_in(i, j, 2); img_in(i, j, 3)];
img_trans(i, j, :) = (color_tran * rgb)';
end
end
elseif color_mode == 2
color_tran = [1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0; 1.0 / 2.0, 0, -1.0 / 2.0; 1.0 / 4.0, -1.0 / 2.0, 1.0 / 4.0];
for i = 1:row
for j = 1:col
rgb = [img_in(i, j, 1); img_in(i, j, 2); img_in(i, j, 3)];
img_trans(i, j, :) = (color_tran * rgb)';
end
end
end
endfunction [sim_blk,sim_num,sim_blk_idx]=search_similar_block(ik,jk,block,tran_block,np,bn_r,bn_c,tau,max_sim_num)
% 搜索相似塊
% Inputs:
% ik, jk： 待搜索相似塊的索引
% block: 圖像塊集合
% tran_block： 圖像塊FFT硬閾值過濾后的FFT系數
% k: 圖像塊大小
% np: floor(nHard / pHard)， 其中nHard表示圖像的搜索區域大小， pHard表示塊的移動步長
% bn_r, bn_c: 圖像總的行/列可以采集圖像塊的數目
% tau: 圖像塊相似性判斷閾值，見公式(1)
% max_sim_num: 最多保留相似塊的數目
% Ouputs:
% sim_blk:
% sim_num:
% sim_blk_idx:
%
% 搜索窗口的左上角，右下角的塊索引
in_s = max(ik-floor(np/2),1);
jn_s = max(jk-floor(np/2),1);
in_e = min(ik+floor(np/2),bn_r);
jn_e = min(jk+floor(np/2),bn_c);
% 當前參考塊
ref_blk = tran_block(:,:,((ik-1)*bn_c+jk));
ii = in_s:1:in_e;
jj = jn_s:1:jn_e;
[II,JJ] = meshgrid(ii,jj);
IDX = (II-1)*bn_c+JJ;
blk_idx=IDX(:);
% 收縮范圍內的全部圖像塊
cur_blk=tran_block(:,:,blk_idx);
cnt=size(cur_blk,3);
ref_blk_mat=repmat(ref_blk,[1,1,cnt]);
delta_blk=cur_blk-ref_blk_mat;
dist=sum(sum(delta_blk.*delta_blk,1),2);
[dist_sort,dist_idx]=sort(dist);
% 最大相似塊是真實相似塊和目標參數相似塊的最小值
max_num=min(cnt,max_sim_num);
if(dist_sort(max_num)<tau)
sim_num=max_num;
else
sim_num=sum(dist_sort(1:max_num)<tau);
end
cnt_idx=dist_idx(1:sim_num);
sim_blk_idx=blk_idx(cnt_idx);
sim_blk=block(:,:,sim_blk_idx);
endfunction [val]=thres_shrink(data,thres)
% 進行閾值截斷： 即 data(i) < thres ? data(i) = 0 : data(i) = data(i)
% Inputs:
% data: 閾值截斷前的數據
% thres: 閾值
% Ouputs:
% val: 閾值截斷后的數據
%
val=data;
idx=find(abs(data)<thres);
val(idx)=0;
endfunction blk_tran3d = transform_3d(blk_3d,tran_mode,lambda2d,lambda1d)
% 進行3D變換，即Collaborative Filtering: 在圖像塊內進行2D變換，在圖像塊間進行1D變換
% 公式(2)
% Inputs:
% blk_3d:
% tran_mode:
% Ouputs:
%
global blk_tran1d_s;
global blk_2d_s;
[m,n,blk_num]=size(blk_3d);
%變換不同時，可能需要修改？？
blk_2d_shrink=zeros(m,n,blk_num);
blk_1d_shrink=zeros(m,n,blk_num);
if(tran_mode==0) %fft
for i=1:1:blk_num
blk_tran2d = fft2(blk_3d(:,:,i));
blk_2d_shrink(:,:,i) = thres_shrink(blk_tran2d,lambda2d);
end
for i=1:1:m
for j=1:1:n
blk_tran1d = fft(blk_2d_shrink(i,j,:));
blk_1d_shrink(i,j,:) = thres_shrink(blk_tran1d,lambda1d);
end
end
blk_tran3d=blk_1d_shrink;
elseif(tran_mode==1) %dct
for i=1:1:blk_num
blk_tran2d=dct2(blk_3d(:,:,i));
blk_2d_shrink(:,:,i)=thres_shrink(blk_tran2d,lambda2d);
end
for i=1:1:m
for j=1:1:n
blk_tran1d=dct(blk_2d_shrink(i,j,:));
blk_1d_shrink(i,j,:)=thres_shrink(blk_tran1d,lambda1d);
end
end
blk_tran3d=blk_1d_shrink;
elseif(tran_mode==2) %dwt
blk_2d_s={};
blk_2d_shrink=[];%zeros()
for i=1:1:blk_num
[blk_tran2d_c,blk_tran2d_s]=wavedec2(blk_3d(:,:,i),2,'Bior1.5');
blk_2d_shrink(:,i)=thres_shrink(blk_tran2d_c,lambda2d);
blk_2d_s{i}=blk_tran2d_s;
end
%這里應該用 wavedec.因為是對1維？？
[blk_tran1d_c,blk_tran1d_s]=wavedec2(blk_2d_shrink,1,'haar');
blk_tran3d=thres_shrink(blk_tran1d_c,lambda1d);
% elseif(strcmp(tran_mode,'db1')) %還未實現
% blk_2d_s={};
% blk_2d_shrink=[];%zeros()
% for i=1:1:blk_num
% [blk_tran2d_cA,blk_tran2d_cH,blk_tran2d_cV,blk_tran2d_cD]=...
% dwt2(blk_3d(:,:,i),'db1');
% blk_2d_shrink(:,i)=thres_shrink(blk_tran2d_c,lambda2d);
% blk_2d_s{i}=blk_tran2d_s;
% end
% [blk_tran1d_c,blk_tran1d_s]=wavedec2(blk_2d_shrink,1,'haar');
% blk_tran3d=thres_shrink(blk_tran1d_c,lambda1d);
else
error('tran_mode error');
end
end3. 模型效果：
3.1 灰度圖像
3.2 彩色圖像
總結
以上是生活随笔為你收集整理的图像去噪序列——BM3D图像去噪模型实现的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： MIUI12极简模式怎么打开 极简模式打
- 下一篇： 斗地主的春天是什么意思？斗地主的春天含义