mvtd.m

这是改进的经验模式分解方法

function m = mvtd(s,inds,MODE_COMPLEX)
% 用中值定理法求极值处的局部均值(MVTD)
% e(te) = 1/(t(i+1)-t(i)+1)*sum(s(t(l)) l=i..i+1
% s(t(l))--信号 t(i),t(i+1) 两个极值点的时刻
% te= (t(i)+t(i+1))/2
% m(t(i+1)) = h(t(i))*e(te)+h(t(i+1))*e(te+1)
% h(t(i)) = (t(i+2))-t(i+1))/(t(i+2)-t(i))
% h(t(i+1)) = (t(i+1)-t(i))/(t(i+2)-t(i))

if size(s,1) > 1
  s = s.';     % 向量需为[1,n]维
end;

if size(inds,1) > 1
  ins = ins';  % 向量需为[1,n]维
end

len = length(inds);
e = zeros(1,len);  % [1,n]

if MODE_COMPLEX		% 复信号情况
  switch MODE_COMPLEX
    case 1,2
    m = s;
  end;
else                % 实信号情况
   for i = 1:len-1
      t = floor((inds(i)+inds(i+1))/2);%% te= (t(i)+t(i+1))/2 {floor不超过x 的最大整数.(高斯取整) floor( [3.12 -3.12])ans = 3    -4 }
      e(t) = 1/(inds(i+1)-inds(i)+1)*(s(i)+s(i+1));% e(te) = 1/(t(i+1)-t(i)+1)*sum(s(t(l)) l=i..i+1
   end;
   for i = 1:len-2
      t = floor((inds(i)+inds(i+1))/2);
      h1(inds(i)) = (inds(i+2)-inds(i+1))/(inds(i+2)-inds(i));
      h2(inds(i)) = (inds(i+1)-inds(i))/(inds(i+2)-inds(i));
      m(inds(i+1)) = h1(i)*e(t)+h2(i)*e(t+1);
   end
end;