% emd_principe.m
%
% P. Flandrin

N = 1000;
t = 1:N;

trend = sin(2*pi*0.0002*t).*hanning(N)';

rap = 10;
b = randn(1,N/rap);
noise = interp(b,rap);

signal = trend+noise/40;

[emin,emax,ezer] = extr(signal);

envsup = spline(emax,signal(emax),t);
envinf = spline(emin,signal(emin),t);
envmoy = (envsup+envinf)/2;
ysup = envsup(501:570);
maxenv = max(ysup);
yinf = envinf(501:570);
minenv = min(yinf);
ymoy = (ysup+yinf)/2;
ychop = signal(501:570);
BB = [1 70 minenv-(maxenv-minenv)/10 maxenv+(maxenv-minenv)/10];
BBB = [1 70 mean(ychop-ymoy)-3*(maxenv-minenv)/5 mean(ychop-ymoy)+3*(maxenv-minenv)/5];

figure(1)
clf

subplot(131)

plot(ychop,'k','LineWidth',2)
title('signal')
axis(BB)
set(gca,'YTick',[]);set(gca,'XTick',[])
axis('square')


pause

hold on

plot(ymoy)
title('signal =')

hold off

subplot(132)

plot(ymoy,'LineWidth',2)
title('slow oscillation ...')

axis(BB)
set(gca,'YTick',[]);set(gca,'XTick',[])

axis('square')

pause

subplot(133)

plot(ychop-ymoy,'r','LineWidth',2)

axis(BBB)
set(gca,'YTick',[]);set(gca,'XTick',[])

axis('square')

title('+ fast oscillation')
