function [imf,ort,nbit] = emd_online(x,t,stop,nbpresift,tst,tst2)%EMD_ONLINE  (On Line Empirical Mode Decomposition) computes on-line EMD%% IMPORTANT: EMD_ONLINE does not truly apply EMD on-line but it does AS IF.% It is rather a demonstration that EMD can be applied on-line.%% stopping criterion for sifting : %   at each point : mean amplitude < threshold*envelope amplitude %                   if mean amplitude > max(envelope amplitude)/tolerance%   &%   at each point : mean amplitude < threshold2*enveloppe amplitude %                   if mean amplitude > max(envelope amplitude)/tolerance2%%[imf,ort,nbits] = EMD_ONLINE(x,t,stop,nbpresift,tst,tst2)% inputs: %         - x : analyzed signal %         - t (optional) : sampling times (default : 1:length(x))%         - stop (optional) : threshold, and threshold2 (optional)%                             tolerance, and tolerance2 (both optional)%                             for sifting stopping criterion %                             default : [0.05,0.5,20,100]%         - nbpresift (optional) : number of sifting by blocks iterations (default 4)%         - tst (optional) : if equals to 1 shows sifting steps%         - tst2 (optional) : if equals to 1 shows sifting by blocks steps                   %% outputs:%         - imf : intrinsic mode functions (last line = residual)%         - ort : index of orthogonality%         - nbits : number of iterations for each mode%% calls:%         - extr (finds extrema and zero-crossings)%         - io : computes the index of orthogonality% % G. Rilling, July 2002% gabriel.rilling@ens-lyon.frDEFSTOP = [0.05,0.5,20,100];% default parameters for sifting stopNBPRESIFT = 4;%number of sifting iterations per blockif(nargin==1)  t = 1:length(x);  stop = DEFSTOP;  tst = 0;  tst2 = 0;endif(nargin==2)  stop = DEFSTOP;  tst = 0;  tst2 = 0;endif (nargin==3)  tst=0;  tst2 = 0;endif (nargin==4)  tst=0;  tst2 = 0;endif (nargin==5)  tst2 = 0;endif nargin > 3    NBPRESIFT = nbpresift;endS = size(x);if ((S(1) > 1) & (S(2) > 1)) | (length(S) > 2)  error('x must have only one row or one column')endif S(1) > 1  x = x';endS = size(t);if ((S(1) > 1) & (S(2) > 1)) | (length(S) > 2)  error('t must have only one row or one column')endif S(1) > 1  t = t';endif (length(t)~=length(x))  error('x and t must have the same length')endS = size(stop);if ((S(1) > 1) & (S(2) > 1)) | (S(1) > 4) | (S(2) > 4) | (length(S) > 2)  error('stop must have only one row or one column of max four elements')endif S(1) > 1  stop = stop';  S = size(stop);endif S(2) < 4    stop(4) = DEFSTOP(4);endif S(2) < 3    stop(3) = DEFSTOP(3);endif S(2) < 2    stop(2) = DEFSTOP(2);endif S(2) == 1  stop=[stop, DEFSTOP(2)];endsd = stop(1);sd2 = stop(2);tol = stop(3);tol2 = stop(4);if tst  figure  figures(1) = gcf;  figure  figures(3) = gcf;endif tst2  figure  figures(2) = gcf;end  MAXITERATIONS=10000;LARGMIN = 5;NBSYM = 2;% maximum number of symmetrized points for interpolationsLARGTRANS = 10;LARGTRANSPS = 5;PAS = 20;STEP = 5;% maximal number of iterations on a modeLX = length(x);% for displaysdt(LX) = 0;sdt = sdt+sd;sd2t(LX) = 0;sd2t = sd2t+sd2;% number of minima and maxima on the considered zonelm = 0;lM = 0;% number of minima and maxima right of the considered zone, % after "stop" or "stopps"lmr = 0;lMr = 0;% same, but left before "start"lml = 0;lMl = 0;% total number of extrema, left and rightnem = 0;nemr = 0;neml = 0;k = 1;nbit = 0;% number of modes, and number of modes on which block siftings are completednbmodes = 1;nbmodes_psdone = 0;start = 1;% end of the constant part of the windowstop = min(PAS+1,LX);% start and end of the considered zonestopr = 1;startl = 1;% end of available data on the considered zonefin = 1;% start of the considered zone for block siftinglimpsl(1,1:NBPRESIFT) = 1;% start and end of segment to which sifting is appliedstartps = 1;stopps(1,NBPRESIFT) = 0;stopps = stopps + 1;% end of available data for block siftingfinps(1,1:NBPRESIFT) = 1;finps(1,1) = 10*PAS;% tests if all data are available for an iteration of block siftinglafinps(1,NBPRESIFT) = 0;% allows to interrupt a mode extraction to process to the next mode% interrupts also if not enough data availablesuspps(1,NBPRESIFT) = 0;% tests for the termination of one iteration of block sifting stoptestps(1,NBPRESIFT) = 0;indmin = [];indmax = [];% tests if all data are available for on-line siftinglafin = 0;% tests if a mode is entirely extractedstoptest = 0;% allows to interrupt a mode extraction to process to the next mode% interrupts also if not enough data availablesusp = 0;% tells if the considered zone has to be moved forward for having enough extremaneedextr = 1;% idem for block siftingneedextrps(1,1:NBPRESIFT) = 1;% tells how many iterations of block sifting have been initiatednbstartedpresift = 1;% modes concerned by block siftingmps = x;% mode concerned by on-line siftingm(LX) = 0;trig = 0;if tst | tst2  disp('appuyer sur une touche pour commencer')  pauseend  while sum(stoptest) < nbmodes % global loop  for k = 1:nbmodes	  nsteps = 0;    waittest = 0;    if k == 1 & trig        suspps(1,1) = 0;        trig = 0;    end            while  sum(stoptestps(k,:)) < NBPRESIFT & ~waittest & sum(suspps(k,:)) < nbstartedpresift(k) % boucle de presifting      for i = 1:nbstartedpresift(k)          	if needextrps(k,i) == 1	  [indmintmp,indmaxtmp] = extr(mps(k,max([(limpsl(k,i)-1),1]):finps(k,i),i));	  nb = sum(indmintmp > stopps(k,i))+sum(indmaxtmp > stopps(k,i));	  stoprps(k,i) = finps(k,i);      	  if nb < 8*LARGTRANSPS & finps(k,i) < LX          suspps(k,i) = 1;          if k == 1 & i == 1              finps(1,1) = min(LX,finps(1,1) + 10*PAS);              [indmintmp,indmaxtmp] = extr(mps(k,max([(limpsl(k,i)-1),1]):finps(k,i),i));              nb = sum(indmintmp > stopps(k,i))+sum(indmaxtmp > stopps(k,i));              stoprps(k,i) = finps(k,i);              trig = 1;          end      else	    lmt = length(indmintmp);	    lMt = length(indmaxtmp);	    if lmt > 0	      indminps(k,1:lmt,i) = indmintmp + max([(limpsl(k,i)-1),1])-1;	    end	    if lMt > 0	      indmaxps(k,1:lMt,i) = indmaxtmp + max([(limpsl(k,i)-1),1])-1;	    end	    if lmt < size(indminps,2)	      indminps(k,length(indmintmp)+1:end,i) = 0;	    end	    if lMt < size(indmaxps,2)	      indmaxps(k,length(indmaxtmp)+1:end,i) = 0;	    end	    needextrps(k,i) = 0;	    	  end	  	  if stoprps(k,i) >= LX	    lafinps(k,i) = 1;	    needextrps(k,i) = 0;	  end	end		if ~suspps(k,i)	  curindminps = indminps(k,find(indminps(k,:,i) >= limpsl(k,i)),i);	  curindmaxps = indmaxps(k,find(indmaxps(k,:,i) >= limpsl(k,i)),i);	  nemps = length(curindminps) + length(curindmaxps);	  	end      	% loop of block (pre)sifting	while (~needextrps(k,i) | lafinps(k,i)) & ~stoptestps(k,i) & ~waittest & ~suspps(k,i)	  if nemps < 3 & lafinps(k,i)	    stoptestps(k,:) = 1;	    stoptest(k) = 1;	    m(k,:) = mps(k,:,i);	    if i > 1	      m(k+1,:) = x - sum(m(1:k,:));	    end	    break	  end	  	  if limpsl(k,i) == 1	    startps(k,i) = 1;	  else	    startps(k,i) = stopps(k,i);	  end	  if lafinps(k,i)	    stopps(k,i) = LX;	    stoptestps(k,i) = 1;	  else	    stopps(k,i) = min(curindminps(max([1,end - LARGTRANSPS+1])),curindmaxps(max([1,end - LARGTRANSPS+1]))); % si ~lafinps(k,i)	  end	        if startps(k,i) == stopps(k,i)          pause             needextrps(k,i) = 1;          break      end      	  lmr = sum(curindminps > stopps(k,i));	  lMr = sum(curindmaxps > stopps(k,i));	  nemrps(k,i) = lmr + lMr;	  	  if nemrps(k,i) < 8*LARGTRANSPS	    needextrps(k,i) = 1;	  end	  	  	  if limpsl(k,i) == 1	    margeml = 0;	    margeMl = 0;	             if curindmaxps(1) < curindminps(1)		if mps(k,1,i) > mps(k,curindminps(1),i)	      lmax = fliplr(curindmaxps(2:min(end,NBSYM+1)));	      lmin = fliplr(curindminps(1:min(end,NBSYM)));	      lsym = curindmaxps(1);		else	      lmax = fliplr(curindmaxps(1:min(end,NBSYM)));	      lmin = [fliplr(curindminps(1:min(end,NBSYM-1))),1];	      lsym = 1;		end	      else		if mps(k,1,i) < mps(k,curindmaxps(1),i)	      lmax = fliplr(curindmaxps(1:min(end,NBSYM)));	      lmin = fliplr(curindminps(2:min(end,NBSYM+1)));	      lsym = curindminps(1);		else	      lmax = [fliplr(curindmaxps(1:min(end,NBSYM-1))),1];	      lmin = fliplr(curindminps(1:min(end,NBSYM)));	      lsym = 1;		end	      end	      	    tlmin = 2*t(lsym)-t(lmin);	    tlmax = 2*t(lsym)-t(lmax);	        % in case symmetrized parts do not extend enough    if tlmin(1) > t(1) | tlmax(1) > t(1)      if lsym == curindmaxps(1)	lmax = fliplr(curindmaxps(1:min(end,NBSYM)));      else	lmin = fliplr(curindminps(1:min(end,NBSYM)));      end      if lsym == 1	error('bug')      end      lsym = 1;      tlmin = 2*t(lsym)-t(lmin);      tlmax = 2*t(lsym)-t(lmax);    end	    	    	  else % limpsl ~= 1	    	    lmin = curindminps(find(curindminps <= startps(k,i)));	    lmax = curindmaxps(find(curindmaxps <= startps(k,i)));	    if length(lmin) < 5 |length(lmax) < 5	      error('souci')	    end	    tlmin = t(lmin);	    tlmax = t(lmax);	    margeml = length(lmin);	    margeMl = length(lmax);	    	  end % if limpsl...	    	  if lafinps(k,i)	    	    margemr = 0;	    margeMr = 0;	      	      if curindmaxps(end) < curindminps(end)		if mps(k,LX,i) > mps(k,curindmaxps(end),i)	rmax = [LX,fliplr(curindmaxps(max(end-NBSYM+2,1):end))];	rmin = fliplr(curindminps(max(end-NBSYM+1,1):end));	rsym = LX;		else	rmax = fliplr(curindmaxps(max(end-NBSYM+1,1):end));	rmin = fliplr(curindminps(max(end-NBSYM,1):end-1));	rsym = curindminps(end);		end