emd_blx.m

B样条曲线生成程序,HHT插件等。rilling版本

t=0:1:2000;
s=sin(2*pi*0.08*t)+sin(2*pi*0.11*t)+sin(2*pi*0.15*t)
%%%%%
s1=s;
sd=1;
for n=1:29
if(sd>0.3)
   jd=find(diff(sign(diff(s1)))==-2)+1;  %找h1极大值的位置
   jx=find(diff(sign(diff(s1)))==2)+1; %找h1极小值的位置
   ah=length(jd);
   al=length(jx);
   jdz(ah)=0;
   for i=1:ah
       bh=jd(i);
       jdz(i)=s1(bh);
   end                     %找出h1极大值对应函数值
   jxz(al)=0;
   for i=1:al
       bl=jx(i);
       jxz(i)=s1(bl);
   end                     %找出h1极小值对应函数值
   jsbl=spline(jd,jdz,t);  %极大值拟和的上包络
   jxbl=spline(jx,jxz,t);  %极小值拟和的下包络
   m1=(jsbl+jxbl)/2;       %上下包络均值
   h1=s1-m1; 
   sd=sum(((s1-h1)./h1).^2);
end
end
%%%%%%%
[indmin,indmax,indzer]=extr(s,t)
tt1=indmax
tt2=indmin
max=s(indmax)
min=s(indmin)
tt=indmax(1):indmax(end)
TT=indmin(1):indmin(end)
max1=spline(tt1,max,tt)
min1=spline(tt2,min,TT)
%%%%%%%
[up,down] = envelope_s(t,s1,'linear')
%%%%%%
figure(1)
plot(s1,'c')
hold on
plot(jsbl,'r')
plot(jxbl,'r')
plot(jd,jdz,'o')
plot(tt,max1,'m')
plot(TT,min1,'m')
plot(t,up,'b-.');
plot(t,down,'b-.');
hold off;