clear all;        
load woman;          %载入带噪声的图像
init=2055615866;         %生成含噪图像并显示
randn('seed',init);
XX=X+2*randn(size(X));
[c,l]=wavedec2(XX,2,'db2');     %对图像进行消噪处理，用db2小波函数对x进行两层分解a2=wrcoef2(‘a’,c,l, ’db2’,2);     %重构第二层图像的近似系数
a2=wrcoef2('a',c,l,'db2',2);
n=[1,2];          %设置尺度向量
p=[10.28,24.08];             %设置阈值向量
nc = wthcoef2('t',c,l,n,p,'s');      %对高频小波系数进行阈值处理
mc = wthcoef2('t',nc,l,n,p,'s');    %再次对高频小波系数进行阈值处理
X2=waverec2(mc,l,'db2');     %图像的二维小波重构
colormap(map);
subplot(131),image(XX),axis square;
title('含噪图像');
subplot(132),image(a2),axis square;
title('小波分解去噪');
subplot(133),image(X2),axis square;
title('小波阈值去噪');
ps = sum(sum((X-mean(mean(X))).^2));  %计算信噪比
Pn = sum(sum((a2-X).^2));
disp('利用小波2层分解去噪的信噪比');
snr1 = 10*log10(Ps/Pn)
disp('利用小波阈值去噪的信噪比');
Pn1 = sum(sum((X2-X).^2));
snr2 = 10*log10(Ps/Pn1);