小波去噪主程序资料

1、精选优质文档-倾情为你奉上小波去噪主程序%xiaoboquzao.m%对彩色图像进行去噪I = imread(菊花.png,png); % 读入图像X = im2double(I); % 转换成双精度类型x_noise = imnoise(X, gaussian, 0.01); % 加入高斯噪声%提取三个通道信息xr = x_noise(:, :, 1); % R通道xg = x_noise(:, :, 2); % G通道xb = x_noise(:, :, 3); % B通道%估计三个通道的阈值Cr, Sr = wavedec2(xr, 2, sym4);Cg, Sg = wavedec2(

2、xg, 2, sym4);Cb, Sb = wavedec2(xb, 2, sym4);thr_r = Donoho(xr); % R通道全局阈值thr_g = Donoho(xg); % G通道全局阈值thr_b = Donoho(xb); % B通道全局阈值%对三个通道分别进行去噪% Donoho全局阈值 软阈值公式x_soft_r = wdenoise(xr, gbl, s, thr_r, sym4, 2);x_soft_g = wdenoise(xg, gbl, s, thr_g, sym4, 2);x_soft_b = wdenoise(xb, gbl, s, thr_b, sym4

3、, 2);% Donoho全局阈值 硬阈值公式-x_hard_r = wdenoise(xr, gbl, h, thr_r, sym4, 2);x_hard_g = wdenoise(xg, gbl, h, thr_g, sym4, 2);x_hard_b = wdenoise(xb, gbl, h, thr_b, sym4, 2);% Birge-Massart策略 软阈值公式-thr_lvd_r=thr_rx_soft_lvd_r = wdenoise(xr, lvd, s, thr_lvd_r, sym4, 1);thr_lvd_g=thr_gx_soft_lvd_g = wdenois

4、e(xg, lvd, s, thr_lvd_g, sym4, 1);thr_lvd_b=thr_bx_soft_lvd_b = wdenoise(xb, lvd, s, thr_lvd_b, sym4, 1);%半软阈值-x1_r = den1(xr, sym4, 2, thr_r);x1_g = den1(xg, sym4, 2, thr_g);x1_b = den1(xb, sym4, 2, thr_b);%半软阈值 + 均值滤波-x1_5_r = den1_5_1(xr, sym4, 2, thr_r, 0.5*thr_r);x1_5_g = den1_5_1(xg, sym4, 2,

5、thr_g, 0.5*thr_g);x1_5_b = den1_5_1(xb, sym4, 2, thr_b, 0.5*thr_b);%自适应阈值-x4_r = den4(xr, sym4, 2);x4_g = den4(xg, sym4, 2);x4_b = den4(xb, sym4, 2);% 恢复去噪后的图像x_soft = cat(3, x_soft_r, x_soft_g, x_soft_b); % Donoho 软阈值x_hard = cat(3, x_hard_r, x_hard_g, x_hard_b); % Donoho 硬阈值x1 = cat(3, x1_r, x1_g,

6、 x1_b); % 半软阈值x4 = cat(3, x4_r, x4_g, x4_b); % 自适应阈值x1_5= cat(3,x1_5_r,x1_5_g,x1_5_b); % 半软阈值 + 均值滤波%计算去噪图像与原图像峰值信噪比psnr_soft = PSNR_color(x_soft, X)psnr_hard = PSNR_color(x_hard, X)psnr1 = PSNR_color(x1, X)psnr1_5 = PSNR_color(x1_5, X)psnr4 = PSNR_color(x4, X)%显示去噪后的图像figure; imshow(X); title(原图像);

7、figure; imshow(x_noise); title(带噪声图像);figure; imshow(x_soft); title(Donoho 软阈值)figure; imshow(x_hard); title(Donoho 硬阈值);figure; imshow(x1); title(半软阈值);figure; imshow(x4); title(自适应阈值);figure; imshow(x1_5); title(半软阈值加均值滤波);半软阈值去噪方法%den1.mfunction X = den1(x, wname, n, thr)thr1 = 0.5 * thr;C, S = w

8、avedec2(x, n, wname); %对图像进行小波分解dcoef = C( prod(S(1, :) + 1 : end); %高频部分系数ind = find( abs(dcoef) = thr1 & abs(dcoef) = thr ) + prod(S(1, :); %大于thr的系数% C(ind) = sign(C(ind) .* ( abs(C(ind) - alpha * thr );X = waverec2(C, S, wname); %重构图像改进的半软阈值去噪方法%den1_5_1.mfunction X = den1_5_1(x, wname, n, thr,

9、thr1)%对一层的重构图像进行均值值滤波% thr1 = 0.6 * thr;C, S = wavedec2(x, n, wname); %对图像进行小波分解dcoef = C( prod(S(1, :) + 1 : end); %高频部分系数ind = find( abs(dcoef) = thr1 & abs(dcoef) = thr ) + prod(S(1, :); %大于thr的系数% C(ind) = sign(C(ind) .* ( abs(C(ind) - alpha * thr );%重构至第1层A1 = wrcoef2(a, C, S, wname, 1); H1 = w

10、rcoef2(h, C, S, wname, 1);V1 = wrcoef2(v, C, S, wname, 1);D1 = wrcoef2(d, C, S, wname, 1);%对三个子图像进行均值滤波h_fir = 1 1 1 1 1 / 5; %水平方向滤波器v_fir = 1 1 1 1 1 / 5; %垂直方向滤波器d_fir = 0 0 1 0 0; 0 0 1 0 0; 1 1 1 1 1; 0 0 1 0 0; 0 0 1 0 0 / 9; %对角线方向滤波器H1 = filter2(h_fir, H1);V1 = filter2(v_fir, V1);D1 = filter

11、2(d_fir, D1);%中值-% H1 = medfilt2(H1, 3 3);% V1 = medfilt2(V1, 3 3);% D1 = medfilt2(D1, 3 3);%重构图像X = A1 + H1 + V1 + D1;自适应阈值去噪方法%den4.mfunction X = den4(x, wname, n)% Feature Adaptive Wavelet Shrinkage for Image Denoising%初始化参数值R = 5; %窗口大小alpha= 0.1; %控制小波系数缩减的程度beta= 0.3;delta= DELTA(x); %噪方差lambd

12、a2 = 4 * delta2 * log(R); %局部阈值 2 C, S = wavedec2(x, n, wname); %对图像进行小波分解%提取每层系数并进行处理for i = n : -1 : 1 cH = detcoef2(h, C, S, i); %水平细节系数 cV = detcoef2(v, C, S, i); %垂直细节系数 cD = detcoef2(d, C, S, i); %对角线细节系数 dim = size(cH); %分别处理三个方向的系数 for j = 1 : dim(1) for k = 1 : dim(2) S_jk2 = energy(cH, j,

13、k, R); cH(j, k) = shrink(cH(j, k), S_jk2, alpha, beta, lambda2); S_jk2 = energy(cV, j, k, R); cV(j, k) = shrink(cV(j, k), S_jk2, alpha, beta, lambda2); S_jk2 = energy(cD, j, k, R); cD(j, k) = shrink(cD(j, k), S_jk2, alpha, beta, lambda2); end end %再把系数放回去 k = size(S,1) - i;first = S(1,1)*S(1,2) + 3

14、* sum(S(2:k-1, 1).*S(2:k-1, 2) + 1; %起始位置add = S(k,1)*S(k,2); %系数长度 C(first : first + add - 1) = reshape(cH, 1, add); C(first + add : first + 2*add - 1) = reshape(cV, 1, add); C(first + 2*add : first + 3*add - 1) = reshape(cD, 1, add);endX = waverec2(C, S, wname); %重构图像% deltafunction delta = DELTA(

15、x)%估计噪声方差 C, S = wavedec2(x, 1, db1); %小波分解d = C( prod( S(1,:) ) + 2 * prod( S(2,:) ) + 1 : end); % HH子带系数delta = median( abs(d) ) / 0.6745; %计算delta%energyfunction S_jk2 = energy(cM, j, k, R)%计算小波系数附近的能量dim = size(cM);%边界判断row_min = (j-1 = fix(R/2) * fix(-R/2);row_max = (dim(1)-j = fix(R/2) * fix(R

16、/2);col_min = (k-1 = fix(R/2) * fix(-R/2);col_max = (dim(2)-k = fix(R/2) * fix(R/2);s = 0;for m = row_min : row_max for n = col_min : col_max s = cM(j + m, k + n)2 + s; endendS_jk2 = s / R2;% shrink%function d_jk = shrink(d, S_jk2, alpha, beta, lambda2)%处理小波系数if S_jk2 = beta * lambda2 d_jk = d * (1

17、- alpha * lambda2 / S_jk2);else d_jk = 0;endDonoho阈值去噪算法%Donoho.mfunction thr = Donoho(x)%用Donoho通用阈值公式计算阈值 x为要进行处理的图像% thr = delta * sqrt( 2 * log(n)% n为信号的长度或尺寸% delta = MAD / 0.6745 %经验公式，其中MAD为小波分解后高子带系数的中值n = prod( size(x) ); %图像尺寸%计算deltaC, S = wavedec2(x, 1, db1); %小波分解d = C( prod( S(1,:) ) +

18、 2 * prod( S(2,:) ) + 1 : end); %HH子带系数delta = median( abs(d) ) / 0.6745;%计算阈值thr = delta * sqrt(2*log(n);计算彩色图像(含RGB三个通道)的均方差%MSE_color.mfunction mse_color = MSE_color(I, J)% 计算彩色图像(含RGB三个通道)的均方差 m, n, o = size(I);sigma = 0;for i = 1:m for j = 1:n for k = 1:o sigma = sigma + ( I(i,j,k) - J(i,j,k) )

19、2; end endendmse_color = sigma / (m*n) / 3;峰值信噪比计算函数%PSNR_color.mfunction psnr = PSNR_color(I, J)%计算消噪前后图像的峰值信噪比psnr = 10 * log10( 2552 / MSE_color(I, J) );%wdenoise.mfunction X = wdenoise(x, measure, sorh, thr, wname, n)% 阈值去噪函数% x为带噪声图像 % measure表示全局或局部% sorh表示软硬阈值方法%thr为阈值% wname为小波函数名% n为分解层次 C, S = wavedec2(x, n, wname); %对图像进行小波分解switch measure case gbl %全局阈值方法 dcoef = C( prod(S(1, :) + 1 : end); %提