【正文】 %計算snr值，閾值在不同層次估計噪聲%以上的全采用的是軟閾值函數(shù)，現(xiàn)在采用硬閾值函數(shù)來處理數(shù)據(jù)title([39。用minimaxi去噪39。])xd1=wden(s1,39。minimaxi39。,39。h39。,39。one39。,lev,39。db339。)plot(xd1)ps=sum(x.^2)pn11=xxd1pn=sum(pn11.^2)c1=corrcoef(xd1,x)snr1=10*log10(ps/pn) %計算去噪后的信噪比值，閾值不做調(diào)整xd11=wden(s1,39。minimaxi39。,39。h39。,39。sln39。,lev,39。db339。)plot(xd11)ps=sum(x.^2)pn11=xxd11pn=sum(pn11.^2)c11=corrcoef(xd11,x)snr11=10*log10(ps/pn) %計算snr值，閾值根據(jù)第一層系數(shù)進行一次噪聲調(diào)整的情況xd111=wden(s1,39。minimaxi39。,39。h39。,39。mln39。,lev,39。db339。)plot(xd111)ps=sum(x.^2)pn111=xxd111pn=sum(pn111.^2)c111=corrcoef(xd11,x)snr111=10*log10(ps/pn) %計算snr值，閾值在不同層次估計噪聲%總的來看軟閾值去噪的效果要比硬閾值去噪好一點，采用的去噪的默認方法不同，則去噪的效果也不同的。這樣就可以得到：一個小波基有4種默認方法去噪，一種去噪方法有6種不同的小類型。所以一共有4*6*6=24*6=144種結(jié)果。%現(xiàn)在使用軟閾值函數(shù)去噪法來去噪了load synecg,x=synecg39。 s1=awgn(x,3,39。measured39。,39。db39。)subplot(811)plot(x)title([39。原始信號39。])subplot(812)plot(s1)title([39。含噪信號39。])lev=5[c,l]=wavedec(s1,lev,39。db339。)det1=detcoef(c,l,1)sigma=median(abs(det1))/th1=sigma*sqrt(2*log(length(x)))/log2(2)th2=sigma*sqrt(2*log(length(x)))/log2(3)th3=sigma*sqrt(2*log(length(x)))/log2(4)th4=sigma*sqrt(2*log(length(x)))/log2(5)th5=sigma*sqrt(2*log(length(x)))/log2(6)cc5=appcoef(c,l,39。db339。,5)dd5=detcoef(c,l,5)cc4=appcoef(c,l,39。db339。,4)dd4=detcoef(c,l,4)cc3=appcoef(c,l,39。db339。,3)dd3=detcoef(c,l,3)cc2=appcoef(c,l,39。db339。,2)dd2=detcoef(c,l,2)cc1=appcoef(c,l,39。db339。,1)dd1=detcoef(c,l,1)ee=[cc5 dd5 dd4 dd3 dd2 dd1]s2=waverec(ee,l,39。db339。)th=[th1,th2,th3,th4,th5][xd,cxd,lxd,perf0,perfl2]=wdencmp(39。lvd39。,s1,39。db339。,5,th,39。s39。)subplot(813)plot(xd)title([39。去后噪信號39。])ps=sum(x.^2)pn1=xxdpn=sum(pn1.^2)snr=10*log10(ps/pn)c=corrcoef(x,xd)rmse=sqrt(pn/length(x))%用硬閾值來去噪吧load synecg,x=synecg39。 s1=awgn(x,5,39。measured39。,39。db39。)subplot(811)plot(x)title([39。原始信號39。])subplot(812)plot(s1)title([39。含噪信號39。])lev=5[c,l]=wavedec(s1,lev,39。db339。)det1=detcoef(c,l,1)sigma=median(abs(det1))/th1=sigma*sqrt(2*log(length(x)))/log2(2)th2=sigma*sqrt(2*log(length(x)))/log2(3)th3=sigma*sqrt(2*log(length(x)))/log2(4)th4=sigma*sqrt(2*log(length(x)))/log2(5)th5=sigma*sqrt(2*log(length(x)))/log2(6)cc5=appcoef(c,l,39。db339。,5)dd5=detcoef(c,l,5)cc4=appcoef(c,l,39。db339。,4)dd4=detcoef(c,l,4)cc3=appcoef(c,l,39。db339。,3)dd3=detcoef(c,l,3)cc2=appcoef(c,l,39。db339。,2)dd2=detcoef(c,l,2)cc1=appcoef(c,l,39。db339。,1)dd1=detcoef(c,l,1)ee=[cc5 dd5 dd4 dd3 dd2 dd1]s2=waverec(ee,l,39。db339。)th=[th1,th2,th3,th4,th5][xd,cxd,lxd,perf0,perfl2]=wdencmp(39。lvd39。,s1,39。db339。,5,th,39。h39。)subplot(813)plot(xd)title([39。去后噪信號39。])ps=sum(x.^2)pn1=xxdpn=sum(pn1.^2)snr=10*log10(ps/pn)c=corrcoef(x,xd)rmse=sqrt(pn/length(x))%最核心的，估計各個層次的閾值的大小并且使用自己的閾值函數(shù)處理的程序模擬：load synecg,x=synecg39。 s1=awgn(x,5,39。measured39。,39。db39。)subplot(811)plot(x)title([39。原始信號39。])subplot(812)plot(s1)title([39。含噪信號39。])lev=5[c,l]=wavedec(s1,lev,39。db339。)det1=detcoef(c,l,1)sigma=median(abs(det1))/th1=sigma*sqrt(2*log(length(x)))/log2(2)th2=sigma*sqrt(2*log(length(x)))/log2(3)th3=sigma*sqrt(2*log(length(x)))/log2(4)th4=sigma*sqrt(2*log(length(x)))/log2(5)th5=sigma*sqrt(2*log(length(x)))/log2(6)cc5=appcoef(c,l,39。db339。,5)dd5=detcoef(c,l,5)cc4=appcoef(c,l,39。db339。,4)dd4=detcoef(c,l,4)cc3=appcoef(c,l,39。db339。,3)dd3=detcoef(c,l,3)cc2=appcoef(c,l,39。db339。,2)dd2=detcoef(c,l,2)cc1=appcoef(c,l,39。db339。,1)dd1=detcoef(c,l,1)ee=[cc5 dd5 dd4 dd3 dd2 dd1]s2=waverec(ee,l,39。db339。)for i=1:length(dd5)if(abs(dd5(i))th5)newdd5(i)=sign(dd5(i))*(abs(dd5(i))(2*th5)/(1+exp((abs(dd5(i))th5)/2)))。elsenewdd5(i)=0。endendsubplot(813)plot(newdd5)for i1=1:length(dd4)if(abs(dd4(i1))th4)newdd4(i1)=sign(dd4(i1))*(abs(dd4(i1))(2*th4)/(1+exp((abs(dd4(i1))th4)/2)))。elsenewdd4(i1)=0。endendsubplot(814)plot(newdd4)for i2=1:length(dd3)if(abs(dd3(i2))th3)newdd3(i2)=sign(dd3(i2))*(abs(dd3(i2))(2*th3)/(1+exp((abs(dd3(i2))th3)/2)))。elsenewdd3(i2)=0。endendsubplot(815)plot(newdd3)for j=1:length(dd2)if(abs(dd2(j))th2)newdd2(j)=sign(dd2(j))*(abs(dd2(j))(2*th2)/(1+exp((abs(dd2(j))th2)/2)))。elsenewdd2(j)=0。endendsubplot(816)plot(newdd2)for i4=1:length(dd1)if(abs(dd1(i4))th1)newdd1(i4)=sign(dd1(i4))*(abs(dd1(i4))(2*th1)/(1+exp((abs(dd1(i4))th1)/2)))。elsenewdd1(i4)=0。endendsubplot(817)plot(newdd1)newee=[cc5 newdd5 newdd4 newdd3 newdd2 newdd1]news2=waverec(newee,l,39。db339。)subplot(818)plot(news2)title([39。去后噪信號39。])ps=sum(x.^2)pn1=xnews2pn=sum(pn1.^2)snr=10*log10(ps/pn)c=corrcoef(x,news2)rmse=sqrt(pn/length(x))%得到的結(jié)果應(yīng)該是自己設(shè)計的閾值函數(shù)去噪軟閾值》硬閾值》默認去噪函數(shù)去噪。軟硬閾值折中：a=for i=1:length(dd5)if(abs(dd5(i))th5)newdd5(i)=sign(dd5(i))*(abs(dd5(i))a*th5)。elsenewdd5(i)=0。endendsubplot(813)plot(newdd5)for i1=1:length(dd4)if(abs(dd4(i1))th4)newdd4(i1)=sign(dd4(i1))*(abs(dd4(i1))a*th4)。elsenewdd4(i1)=0。endendsubplot(814)plot(newdd4)for i2=1:length(dd3)if(abs(dd3(i2))th3)newdd3(i2)=sign(dd3(i2))*(abs(dd3(i2))a*th3)。elsenewdd3(i2)=0。endendsubplot(815)plot(newdd3)for j=1:length(dd2)if(abs(dd2(j))th2)newdd2(j)=sign(dd2(j))*(abs(dd2(j))a*th2)。elsenewdd2(j)=0。endendsubplot(816)plot(newdd2)for i4=1:length(dd1)if(abs(dd1(i4))th1)newdd1(i4)=sign(dd1(i4))*(abs(dd1(i4))a*th1)。elsenewdd1(i4)=0。endendsubplot(817)plot(newdd1)newee=[cc5 newdd5 newdd4 newdd3 newdd2 newdd1]news2=waverec(newee,l,39。db339。)subplot(818)plot(news2)title([39。去后噪信號39。])ps=sum(x.^2)pn1=xnews2pn=sum(pn1.^2)snr=10*log10(ps/pn)c=corrcoef(x,news2)rmse=sqrt(pn/length(x))%得到的結(jié)果應(yīng)該是自己設(shè)計的閾值函數(shù)去噪軟閾值》硬閾值》默認去噪函數(shù)去噪。t=20:1: