clear
sqrt_snr=4;init=2055415866;%設(shè)置信噪比和隨機(jī)種子
[x,y]=wnoise(1,11,sqrt_snr,init);%產(chǎn)生原始信號x和被高斯白噪聲污染的信號xnoised
figure(1)
subplot(211);plot(x);title('參考信號','fontsize',8)
subplot(212);plot(y);title('染噪信號','fontsize',8)

%用sym6小波對原始信號作4層分解并提取系數(shù)

lev=4;wname='sym6';
[c,l]=wavedec(y,lev,wname);
a4=appcoef(c,l,wname,lev);
d4=detcoef(c,l,4);
d3=detcoef(c,l,3);
d2=detcoef(c,l,2);
d1=detcoef(c,l,1);

%1.對信號進(jìn)行強(qiáng)制性消噪處理
dd4=zeros(1,length(d4));
dd3=zeros(1,length(d3));
dd2=zeros(1,length(d2));
dd1=zeros(1,length(d1));
c1=[a4 dd4 dd3 dd2 dd1];
y1=waverec(c1,l,wname);

figure(2)
subplot(311);plot(y1)
title('強(qiáng)制消噪','fontsize',8)
grid on
%求原始信號與估計信號之間的均方根誤差(RMSE)和經(jīng)小波去噪后的估計信號的信噪比(SNR)
sum=0;sum1=0;
n=length(x);
for i=1:n
    sum=sum+(x(i)-y1(i)).^2;
    sum1=sum1+x(i).^2;
end
RMSE1=sqrt(sum/n);
SNR1=10*log10(sum1/sum);


%2.用默認(rèn)閾值對信號進(jìn)行消噪處理
%用ddencmp函數(shù)獲得信號的默認(rèn)閾值

[thr,sorh,keepapp]=ddencmp('den','wv',y);
y2=wdencmp('gbl',y,wname,lev,thr,sorh,keepapp);

subplot(312);plot(y2)
title('默認(rèn)閾值消噪','fontsize',8)
grid on
%求原始信號與估計信號之間的均方根誤差(RMSE)和經(jīng)小波去噪后的估計信號的信噪比(SNR)
sum=0;sum1=0;
n=length(x);
for i=1:n
    sum=sum+(x(i)-y2(i)).^2;
    sum1=sum1+x(i).^2;
end
RMSE2=sqrt(sum/n);
SNR2=10*log10(sum1/sum);


%3.用給定閾值進(jìn)行消噪處理(這里給定軟閾值)

lev=4;wname='sym6';
[c,l]=wavedec(y,lev,wname);
a4=appcoef(c,l,wname,lev);
d4=detcoef(c,l,4);
d3=detcoef(c,l,3);
d2=detcoef(c,l,2);
d1=detcoef(c,l,1);
softd1=wthresh(d1,'h',4.39);
softd2=wthresh(d2,'h',4.35);
softd3=wthresh(d3,'h',4.31);
softd4=wthresh(d4,'h',4.36);
c2=[a4 softd4 softd3 softd2 softd1];
y3=waverec(c2,l,wname);

subplot(313);plot(y3)
title('給定軟閾值消噪','fontsize',8)
grid on
%求原始信號與估計信號之間的均方根誤差(RMSE)和經(jīng)小波去噪后的估計信號的信噪比(SNR)
sum=0;sum1=0;
n=length(x);
for i=1:n
    sum=sum+(x(i)-y3(i)).^2;
    sum1=sum1+x(i).^2;
end
RMSE3=sqrt(sum/n);
SNR3=10*log10(sum1/sum);