%默認閾值
%使用wdencmp進行降噪
load data;
x=data;
figure(1);
subplot(211),plot(x),title('噪聲信號');
%尋找默認值
[thr,sorh,keepapp]=ddencmp('den','wv',x)
%使用全局閾值進行降噪
xd=wdencmp('gbl',x,'db3',5,thr,sorh,keepapp);
a=xd;
subplot(212),plot(xd),title('默認閾值降噪信號')

 data=a;
%hectorload data.txt;
%x=load(hectorload);
%data=data(:,2);
[n,p]=size(data);
t=1:n;
%plot(t,data)
 %hold on;
mu=1/2*mean(data)
avg=double(zeros(n,1));
avg1=double(zeros(n,1));
signal=double(zeros(n,1));
noise=double(zeros(n,1));
ii=1;
t=0;
k=0;
sum1=0.0;%信號加權平均
sum2=0.0;%噪聲加權平均
sum3=0.0;%信號平均
sum4=0.0;%噪聲平均
aa=0.5;
bb=2;
while (ii<=n)
   while (ii<=n & data(ii)<=mu )    %求信號加權平均或者直接平均
       m=1;
      while ( ii<=n & data(ii)<=mu )
         
         signal(m)=data(ii);
           m=m+1;
           ii=ii+1;
      end
        sum3=sum3+sum(signal);
        m=m-1; 
        avg(ii-m:ii-1)=sum(signal)/m;
        
        %將每段信號再分成3段,設定權重
        r=fix(m/3);
        
        signal(1:r)=aa*signal(1:r);   %加權
        signal(r+1:2*r)=bb*signal(r+1:2*r);%加權
        signal(2*r+1:m)=aa*signal(2*r+1:m);%加權
        avg1(ii-m:ii-1)=sum(signal(1:m))/m ;
        sum1=sum1+sum(signal(1:m)) ;
        
       
   end
   while (ii<=n & data(ii)>mu )   %求噪聲
              b=1;
             
       while(ii<=n & data(ii)>mu)
          noise(b)=data(ii);
        
            b=b+1;
           ii=ii+1;
       end
        sum4=sum4+sum(noise);
        b=b-1;
        avg(ii-b:ii)=sum(noise)/b; 
        r=fix(b/3);
        noise(1:r)=aa* noise(1:r);
        noise(r+1:2*r)=bb* noise(r+1:2*r);
        noise(2*r+1:b)=aa* noise(2*r+1:b);
        avg1(ii-b:ii-1)=sum(noise)/b ;
        sum2=sum2+sum(noise);
   end
end
for ii=1:n                     %計算信號、噪聲點數
    if (data(ii)<mu)
        k=k+1;
    else t=t+1;
    end 
end
signal_weighing=sum1/k  %信號加權平均
signal=sum3/k    %信號平均
noise_weighing=sum2/t    %噪聲加權平均
noise=sum4/t    %噪聲平均

s=1:n;                 %畫圖
avg2=avg(1:n);
figure(2)
plot(s,avg2,'LineWidth',1,'color','k')  %直接平均圖像
hold on
plot(s,avg1,'LineWidth',1,'color','r')  %加權平均圖像
hold on