• 大小: 2KB
    文件类型: .zip
    金币: 1
    下载: 0 次
    发布日期: 2021-06-19
  • 语言: 其他
  • 标签: ceemd  改进版  

资源简介

改进版程序ceemd,加入噪声,抑制噪声,有不错的效果。

资源截图

代码片段和文件信息

function [imftresidtspecwaddmeancountimfnresidn] = ceemdnew(xthrshsthreshNalpha)
% Empiricial Mode Decomposition (Hilbert-Huang Transform)
% [imfresidspecwaphsmeancount] = emd(xthrshsthreshNalpha)
%   Inputs:
%         thrsh=thershold variance [0.2 to 0.3 usually]
 %        sthrh=threshold for first iteration [between 4/8]
%         N=number of trials for the ensemble emd
%         alpha=multiplier for white noise (defaults to .1*std(x) as  suggested by Huang)
%   Outputs:
%         imf= intrinsic mode functions (imf(:end))
%         resid=residue; spec is the spectrum (A*e^jwt) w is the
%         frequency a is the ampitlude phs is the instphase
x   = transpose(x(:));
imf = [];
k=0;
xt=x;
if nargin<5
%     alpha=.1*std(xt);
alpha=0.5*std(xt);% 按照需要修改,huang的建议是0.1,但是0.1会造成端点噪声大。
end
cx=0;
s=0;
i=1;
lim=1;
imfprev=[];
while i    if i~=1
        Noise=alpha.*randn(size(xt));
    else
        Noise=0;
    end
    for t=1:lim
        count=0;
        x=xt+(-1)^t*Noise;
        while ~ismonotonic(x)
            count=count+1;
            if i==1
                k=k+1;
                imf{count}=zeros(size(x));
                meancount(k)=0;
            end
           x1 = x;
           VA = Inf;
           cn=0;
           while (VA > thrsh) | ~isimf(x1)
               cn=cn+1;
              s1 = getspline(x1);
              s2 = -getspline(-x1);
              x2 = x1-(s1+s2)/2;
              VA = var((s1+s2)./(2*x1));
              x1 = x2;
              if i==1 cx=cx+1; elseif cn>cx*2 break; end
              if i==1 && cx>numel(xt)/2
                  s=s+1;
                  if s>sthresh
                      x1=xt+alpha.*randn(size(xt));
                      s=0;
                      cx=0;
                      disp(‘Restarting Sifting‘);
                  end
              end
           end
%            if cn>cx*2 && i~=1
%                disp(‘Too long for ensemble member breaking count‘)
%                flag=1;
%                break;
%            else
%                flag=0;
%            end
           
           if count<=k 
               if i==1 
                   imf{count} = x1+imf{count}; 
               else
                   imf{tcount}=x1+imf{tcount}; 
               end;
           else
               break;
           end
           x=x-x1;
        end
        
        if i==1
            resid=repmat(x21);
            imf=repmat(imf21); 
%         elseif flag==0
        else
            resid(t:) = x+resid(t:);
            if t==lim
                meancount=meancount+1;
            end
%             residprev=x;
%         else
%             if t==lim
%                 resid(1:)=resid(1:)-residprev;
%                 for l=1:numel(imf(1:))
%                     imf{1l}=imf{1l}-imfprev{l};
%                 end
%             end
%             i=i-1;
        end


 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----
     文件        4623  2019-05-13 11:39  ceemd (1).m

评论

共有 条评论