请各位初学HHT入门的兄弟姐妹进来谈谈体会,请各位大牛们进来指点!

zhangnan3509

程序不会出错，肯定是自己使用的问题。没有顺利运行先看MATLAB的提示到底是什么原因，和程序的解释说明

jingrenzhi

请问AMGAUSS这个函数怎么用?

zhangnan3509

这个函数就是产生一个高斯调幅信号，对于EMD使用没有影响。我顺手发上去的。有兴趣可以作个例子！

jingrenzhi

你好,我要求取信号的包络及均值曲线的话,用哪个函数?谢谢!

jingrenzhi

在emd软件包里面,怎么没找到求取信号的包络及均值曲线的函数?各位是怎么做的?

zhangnan3509

你好好看看EMD.m程序，里面有这些程序，都做成模块了。

jingrenzhi

我对下面的程序运行了一下,结果是"mode 1 enregistre
mode 2 enregistre
mode 3 enregistre
mode 4 enregistre
mode 5 enregistre
mode 6 enregistre"
请问,问题出在哪里呢?



clc;
clear;
t=0:0.0004:1;
x=sin(2*pi*50*t)+0.3*sin(5.5*pi*50*t);
% x=sin(2*pi*50*t)+0.3*sin(6.5*pi*50*t);
tst=1;
stop = [0.01,0.1,0.01];
NBSYM = 2;
lx = length(x);
sd = stop(1);
sd2 = stop(2);
tol = stop(3);
MAXITERATIONS=2000;

sdt(lx) = 0;
sdt = sdt+sd;
sd2t(lx) = 0;
sd2t = sd2t+sd2;

% maximum number of extrema and zero-crossings in residual
ner = lx;
nzr = lx;
r = x;
imf = [];
k = 1;
% iterations counter for extraction of 1 mode
nbit=0;

% total iterations counter
NbIt=0;

while ner > 2
        
  % current mode
  m = r;
  
  % mode at previous iteration
  mp = m;
  
  sx = sd+1;
  
  % tests if enough extrema to proceed
  test = 0;
  [indmin,indmax,indzer] = extr(m);
  lm=length(indmin);
  lM=length(indmax);
  nem=lm + lM;
  nzm=length(indzer);
  
  j=1;
  
  % sifting loop
  while ( mean(sx > sd) > tol | any(sx > sd2) | (abs(nzm-nem)>1)) & (test == 0) & nbit<MAXITERATIONS
   
    if(nbit>MAXITERATIONS/5 & mod(nbit,floor(MAXITERATIONS/10))==0)
      disp(['mode ',int2str(k),' nombre d iterations : ',int2str(nbit)])
      disp(['stop parameter mean value : ',num2str(s)])
    end
   
    % boundary conditions for interpolations :
        lm=length(indmax);
        ln=length(indmin);
     if indmax(1) < indmin(1)
      if m(1) > m(indmin(1))
        lmax = fliplr(indmax(2:min(lm,NBSYM+1)));
        lmin = fliplr(indmin(1:min(ln,NBSYM)));
        lsym = indmax(1);
      else
        lmax = fliplr(indmax(1:min(lm,NBSYM)));
        lmin = [fliplr(indmin(1:min(ln,NBSYM-1))),1];
        lsym = 1;
      end
    else

      if m(1) < m(indmax(1))
        lmax = fliplr(indmax(1:min(lm,NBSYM)));
        lmin = fliplr(indmin(2:min(ln,NBSYM+1)));
        lsym = indmin(1);
      else
        lmax = [fliplr(indmax(1:min(lm,NBSYM-1))),1];
        lmin = fliplr(indmin(1:min(ln,NBSYM)));
        lsym = 1;
      end
    end
   
    if indmax(end) < indmin(end)
      if m(end) < m(indmax(end))
        rmax = fliplr(indmax(max(lm-NBSYM+1,1):end));
        rmin = fliplr(indmin(max(ln-NBSYM,1):end-1));
        rsym = indmin(end);
      else
        rmax = [lx,fliplr(indmax(max(lm-NBSYM+2,1):end))];
        rmin = fliplr(indmin(max(ln-NBSYM+1,1):end));
        rsym = lx;
      end
    else
      if m(end) > m(indmin(end))
        rmax = fliplr(indmax(max(lm-NBSYM,1):end-1));
        rmin = fliplr(indmin(max(ln-NBSYM+1,1):end));
        rsym = indmax(end);
      else
        rmax = fliplr(indmax(max(lm-NBSYM+1,1):end));
        rmin = [lx,fliplr(indmin(max(ln-NBSYM+2,1):end))];
        rsym = lx;
      end
    end
   
    tlmin = 2*t(lsym)-t(lmin);
    tlmax = 2*t(lsym)-t(lmax);
    trmin = 2*t(rsym)-t(rmin);
    trmax = 2*t(rsym)-t(rmax);
   
    % in case symmetrized parts do not extend enough
    if tlmin(1) > t(1) | tlmax(1) > t(1)
      if lsym == indmax(1)
        lmax = fliplr(indmax(1:min(lm,NBSYM)));
      else
        lmin = fliplr(indmin(1:min(ln,NBSYM)));
      end
      if lsym == 1
        error('bug')
      end
      lsym = 1;
      tlmin = 2*t(lsym)-t(lmin);
      tlmax = 2*t(lsym)-t(lmax);
    end   
   
    if trmin(end) < t(lx) | trmax(end) < t(lx)
      if rsym == indmax(end)
        rmax = fliplr(indmax(max(lm-NBSYM+1,1):end));
      else
        rmin = fliplr(indmin(max(ln-NBSYM+1,1):end));
      end
      if rsym == lx
        error('bug')
      end
      rsym = lx;
      trmin = 2*t(rsym)-t(rmin);
      trmax = 2*t(rsym)-t(rmax);
    end
         
    mlmax =m(lmax);
    mlmin =m(lmin);
    mrmax =m(rmax);
    mrmin =m(rmin);
     
    % definition of envelopes from interpolation
        
    envmax = interp1([tlmax t(indmax) trmax],[mlmax m(indmax) mrmax],t,'spline');
    envmin = interp1([tlmin t(indmin) trmin],[mlmin m(indmin) mrmin],t,'spline');

    envmoy = (envmax + envmin)/2;

    m = m - envmoy;
   
    [indmin,indmax,indzer] = extr(m);
    lm=length(indmin);
    lM=length(indmax);
    nem = lm + lM;
    nzm = length(indzer);

    % evaluation of mean zero
    sx=2*(abs(envmoy))./(abs(envmax-envmin));
    s = mean(sx);
     mp = m;
    nbit=nbit+1;
    NbIt=NbIt+1;

    if(nbit==(MAXITERATIONS-1))
      warning(['forced stop of sifting : too many iterations... mode ',int2str(k),'. stop parameter mean value : ',num2str(s)])
    end
end
  imf(k,:) = m;
  if tst
    disp(['mode ',int2str(k),' enregistre'])
  end
  nbits(k) = nbit;
  k = k+1;
  r = r - m;
  [indmin,indmax,indzer] = extr(r);
  ner = length(indmin) + length(indmax);
  nzr = length(indzer);
  nbit=1;

  if (max(r) - min(r)) < (1e-10)*(max(x) - min(x))
    if ner > 2
      warning('forced stop of EMD : too small amplitude')
    else

      disp('forced stop of EMD : too small amplitude')
    end
    break
  end
  
end

imf(k,:) = r;

if tst
  close
end

jingrenzhi

为啥不是图?

eight

原帖由 jingrenzhi 于 2007-5-17 21:02 发表
为啥不是图?

你贴的程序结果是正确的，至于如何画图，建议先看看 emd_fmsin.m

[ 本帖最后由 eight 于 2007-5-17 21:09 编辑 ]

jingrenzhi

我的软件包里面,没发现有这样一个函数呀?

eight

原帖由 jingrenzhi 于 2007-5-17 21:18 发表
我的软件包里面,没发现有这样一个函数呀?

package_emd\examples\NSIP2003

实在不明白，为何你找不到，不就几个文件夹吗？不想自己找的话用 windows 的搜索命令该不会找不到吧

jingrenzhi

哦!不好意思!我是真的没找到!现在找到了!谢谢:loveliness:

jingrenzhi

周末了，版面显的冷清多了!

jingrenzhi

您好!我看了你说的 emd_fmsin.m程序,还是不知道怎么将这个结果画成图?请指点!谢谢!

eight

原帖由 jingrenzhi 于 2007-5-19 15:41 发表
您好!我看了你说的 emd_fmsin.m程序,还是不知道怎么将这个结果画成图?请指点!谢谢!

请学会自己调试一下。其实就是调用了一个函数，而且肯定是在调用了 emd 之后的，不可能找不到