wolffd@0: function varargout = mirenvelope(orig,varargin)
wolffd@0: %   e = mirenvelope(x) extracts the envelope of x, showing the global shape
wolffd@0: %       of the waveform.
wolffd@0: %   mirenvelope(...,m) specifies envelope extraction method.
wolffd@0: %       Possible values:
wolffd@0: %           m = 'Filter' uses a low-pass filtering. (Default strategy)
wolffd@0: %           m = 'Spectro' uses a spectrogram.
wolffd@0: %
wolffd@0: %   Options related to the 'Filter' method:
wolffd@0: %       mirenvelope(...,'Hilbert'): performs a preliminary Hilbert
wolffd@0: %           transform.
wolffd@0: %       mirenvelope(...,'PreDecim',N) downsamples by a factor N>1, where 
wolffd@0: %           N is an integer, before the low-pass filtering (Klapuri, 1999).
wolffd@0: %           Default value: N = 1.
wolffd@0: %       mirenvelope(...,'Filtertype',f) specifies the filter type.
wolffd@0: %           Possible values are:
wolffd@0: %               f = 'IIR': filter with one autoregressive coefficient 
wolffd@0: %                   (default)
wolffd@0: %               f = 'HalfHann': half-Hanning (raised cosine) filter 
wolffd@0: %                   (Scheirer, 1998)
wolffd@0: %           Option related to the 'IIR' option:
wolffd@0: %           mirenvelope(...,'Tau',t): time constant of low-pass filter in 
wolffd@0: %               seconds.
wolffd@0: %               Default value: t = 0.02 s.
wolffd@0: %       mirenvelope(...,'PostDecim',N) downsamples  by a factor N>1, where 
wolffd@0: %           N is an integer, after the low-pass filtering.
wolffd@0: %           Default value: N = 16 if 'PreDecim' is not used, else N = 1.
wolffd@0: %       mirenvelope(...,'Trim'): trims the initial ascending phase of the
wolffd@0: %           curves related to the transitory state.
wolffd@0: %
wolffd@0: %   Options related to the 'Spectro' method:
wolffd@0: %       mirenvelope(...,b) specifies whether the frequency range is further
wolffd@0: %           decomposed into bands. Possible values:
wolffd@0: %               b = 'Freq': no band decomposition (default value)
wolffd@0: %               b = 'Mel': Mel-band decomposition 
wolffd@0: %               b = 'Bark': Bark-band decomposition 
wolffd@0: %               b = 'Cents': decompositions into cents 
wolffd@0: %       mirenvelope(...,'Frame',...) specifies the frame configuration.
wolffd@0: %           Default value: length: .1 s, hop factor: 10 %.
wolffd@0: %       mirenvelope(...,'UpSample',N) upsamples by a factor N>1, where 
wolffd@0: %           N is an integer.
wolffd@0: %           Default value if 'UpSample' called: N = 2
wolffd@0: %       mirenvelope(...,'Complex') toggles on the 'Complex' method for the
wolffd@0: %           spectral flux computation.
wolffd@0: %
wolffd@0: %   Other available for all methods:
wolffd@0: %       mirenvelope(...,'Sampling',r): resamples to rate r (in Hz).
wolffd@0: %           'Down' and 'Sampling' options cannot therefore be combined.
wolffd@0: %       mirenvelope(...,'Halfwave'): performs a half-wave rectification.
wolffd@0: %       mirenvelope(...,'Center'): centers the extracted envelope.
wolffd@0: %       mirenvelope(...,'HalfwaveCenter'): performs a half-wave
wolffd@0: %           rectification on the centered envelope.
wolffd@0: %       mirenvelope(...,'Log'): computes the common logarithm (base 10) of
wolffd@0: %           the envelope.
wolffd@0: %       mirenvelope(...,'Mu',mu): computes the logarithm of the
wolffd@0: %           envelope, before the eventual differentiation, using a mu-law
wolffd@0: %           compression (Klapuri, 2006).
wolffd@0: %               Default value for mu: 100
wolffd@0: %       mirenvelope(...,'Log'): computes the logarithm of the envelope.
wolffd@0: %       mirenvelope(...,'Power'): computes the power (square) of the
wolffd@0: %           envelope.
wolffd@0: %       mirenvelope(...,'Diff'): computes the differentation of the
wolffd@0: %           envelope, i.e., the differences between successive samples.
wolffd@0: %       mirenvelope(...,'HalfwaveDiff'): performs a half-wave
wolffd@0: %           rectification on the differentiated envelope.
wolffd@0: %       mirenvelope(...,'Normal'): normalizes the values of the envelope by
wolffd@0: %           fixing the maximum value to 1.
wolffd@0: %       mirenvelope(...,'Lambda',l): sums the half-wave rectified envelope
wolffd@0: %           with the non-differentiated envelope, using the respective
wolffd@0: %           weight 0<l<1 and (1-l). (Klapuri et al., 2006)
wolffd@0: %       mirenvelope(...,'Smooth',o): smooths the envelope using a moving
wolffd@0: %           average of order o.
wolffd@0: %           Default value when the option is toggled on: o=30
wolffd@0: %       mirenvelope(...,'Gauss',o): smooths the envelope using a gaussian
wolffd@0: %           of standard deviation o samples.
wolffd@0: %           Default value when the option is toggled on: o=30
wolffd@0: %       mirenvelope(...,'Klapuri06'): follows the model proposed in
wolffd@0: %           (Klapuri et al., 2006). 
wolffd@0: 
wolffd@0:         method.type = 'String';
wolffd@0:         method.choice = {'Filter','Spectro'};
wolffd@0:         method.default = 'Filter';
wolffd@0:     option.method = method;
wolffd@0:     
wolffd@0: %% options related to 'Filter':
wolffd@0: 
wolffd@0:         hilb.key = 'Hilbert';
wolffd@0:         hilb.type = 'Boolean';
wolffd@0:         hilb.default = 0;
wolffd@0:     option.hilb = hilb;
wolffd@0:  
wolffd@0:         decim.key = {'Decim','PreDecim'};
wolffd@0:         decim.type = 'Integer';
wolffd@0:         decim.default = 0;
wolffd@0:     option.decim = decim;
wolffd@0:     
wolffd@0:         filter.key = 'FilterType';
wolffd@0:         filter.type = 'String';
wolffd@0:         filter.choice = {'IIR','HalfHann',0};
wolffd@0:         if isamir(orig,'mirenvelope')
wolffd@0:             filter.default = 0; % no more envelope extraction, already done
wolffd@0:         else
wolffd@0:             filter.default = 'IIR';
wolffd@0:         end
wolffd@0:     option.filter = filter;
wolffd@0: 
wolffd@0:             %% options related to 'IIR': 
wolffd@0:             tau.key = 'Tau';
wolffd@0:             tau.type = 'Integer';
wolffd@0:             tau.default = .02;
wolffd@0:     option.tau = tau;
wolffd@0:     
wolffd@0:         zp.key = 'ZeroPhase'; % internal use: for manual filtfilt
wolffd@0:         zp.type = 'Boolean';
wolffd@0:         if isamir(orig,'mirenvelope')
wolffd@0:             zp.default = 0;
wolffd@0:         else
wolffd@0:             zp.default = NaN;
wolffd@0:         end
wolffd@0:     option.zp = zp;
wolffd@0: 
wolffd@0:         ds.key = {'Down','PostDecim'};
wolffd@0:         ds.type = 'Integer';
wolffd@0:         if isamir(orig,'mirenvelope')
wolffd@0:             ds.default = 1;
wolffd@0:         else
wolffd@0:             ds.default = NaN; % 0 if 'PreDecim' is used, else 16
wolffd@0:         end
wolffd@0:         ds.when = 'After';
wolffd@0:         ds.chunkcombine = 'During';
wolffd@0:     option.ds = ds;
wolffd@0: 
wolffd@0:         trim.key = 'Trim';
wolffd@0:         trim.type = 'Boolean';
wolffd@0:         trim.default = 0;
wolffd@0:         trim.when = 'After';
wolffd@0:     option.trim = trim;
wolffd@0:      
wolffd@0: %% Options related to 'Spectro':
wolffd@0: 
wolffd@0:         band.type = 'String';
wolffd@0:         band.choice = {'Freq','Mel','Bark','Cents'};
wolffd@0:         band.default = 'Freq';
wolffd@0:     option.band = band;
wolffd@0: 
wolffd@0:         up.key = {'UpSample'};
wolffd@0:         up.type = 'Integer';
wolffd@0:         up.default = 0;
wolffd@0:         up.keydefault = 2;
wolffd@0:         up.when = 'After';
wolffd@0:     option.up = up;
wolffd@0: 
wolffd@0:         complex.key = 'Complex';
wolffd@0:         complex.type = 'Boolean';
wolffd@0:         complex.default = 0;
wolffd@0:         complex.when = 'After';
wolffd@0:     option.complex = complex;    
wolffd@0:     
wolffd@0: %% Options related to all methods:
wolffd@0:     
wolffd@0:         sampling.key = 'Sampling';
wolffd@0:         sampling.type = 'Integer';
wolffd@0:         sampling.default = 0;
wolffd@0:         sampling.when = 'After';
wolffd@0:     option.sampling = sampling;
wolffd@0: 
wolffd@0:         hwr.key = 'Halfwave';
wolffd@0:         hwr.type = 'Boolean';
wolffd@0:         hwr.default = 0;
wolffd@0:         hwr.when = 'After';
wolffd@0:     option.hwr = hwr;
wolffd@0: 
wolffd@0:         c.key = 'Center';
wolffd@0:         c.type = 'Boolean';
wolffd@0:         c.default = 0;
wolffd@0:         c.when = 'After';
wolffd@0:     option.c = c;
wolffd@0:     
wolffd@0:         chwr.key = 'HalfwaveCenter';
wolffd@0:         chwr.type = 'Boolean';
wolffd@0:         chwr.default = 0;
wolffd@0:         chwr.when = 'After';
wolffd@0:     option.chwr = chwr;
wolffd@0:     
wolffd@0:         mu.key = 'Mu';
wolffd@0:         mu.type = 'Integer';
wolffd@0:         mu.default = 0;
wolffd@0:         mu.keydefault = 100;
wolffd@0:         mu.when = 'After';
wolffd@0:     option.mu = mu;
wolffd@0:     
wolffd@0:         oplog.key = 'Log';
wolffd@0:         oplog.type = 'Boolean';
wolffd@0:         oplog.default = 0;
wolffd@0:         oplog.when = 'After';
wolffd@0:     option.log = oplog;
wolffd@0: 
wolffd@0:         oppow.key = 'Power';
wolffd@0:         oppow.type = 'Integer';
wolffd@0:         oppow.default = 0;
wolffd@0:         oppow.when = 'After';
wolffd@0:     option.power = oppow;
wolffd@0: 
wolffd@0:         diff.key = 'Diff';
wolffd@0:         diff.type = 'Integer';
wolffd@0:         diff.default = 0;
wolffd@0:         diff.keydefault = 1;
wolffd@0:         diff.when = 'After';
wolffd@0:     option.diff = diff;
wolffd@0:     
wolffd@0:         diffhwr.key = 'HalfwaveDiff';
wolffd@0:         diffhwr.type = 'Integer';
wolffd@0:         diffhwr.default = 0;
wolffd@0:         diffhwr.keydefault = 1;
wolffd@0:         diffhwr.when = 'After';
wolffd@0:     option.diffhwr = diffhwr;
wolffd@0: 
wolffd@0:         lambda.key = 'Lambda';
wolffd@0:         lambda.type = 'Integer';
wolffd@0:         lambda.default = 1;
wolffd@0:         lambda.when = 'After';
wolffd@0:     option.lambda = lambda;
wolffd@0: 
wolffd@0:         aver.key = 'Smooth';
wolffd@0:         aver.type = 'Integer';
wolffd@0:         aver.default = 0;
wolffd@0:         aver.keydefault = 30;
wolffd@0:         aver.when = 'After';
wolffd@0:     option.aver = aver;
wolffd@0:         
wolffd@0:         gauss.key = 'Gauss';
wolffd@0:         gauss.type = 'Integer';
wolffd@0:         gauss.default = 0;
wolffd@0:         gauss.keydefault = 30;
wolffd@0:         gauss.when = 'After';
wolffd@0:     option.gauss = gauss;
wolffd@0: 
wolffd@0:         norm.key = 'Normal';
wolffd@0:         norm.type = 'Boolean';
wolffd@0:         norm.default = 0;
wolffd@0:         norm.when = 'After';
wolffd@0:     option.norm = norm;
wolffd@0: 
wolffd@0:         presel.type = 'String';
wolffd@0:         presel.choice = {'Klapuri06'};
wolffd@0:         presel.default = 0;
wolffd@0:     option.presel = presel;    
wolffd@0:     
wolffd@0:         frame.key = 'Frame';
wolffd@0:         frame.type = 'Integer';
wolffd@0:         frame.number = 2;
wolffd@0:         frame.default = [.1 .1];
wolffd@0:     option.frame = frame;
wolffd@0:         
wolffd@0: specif.option = option;
wolffd@0: 
wolffd@0: specif.eachchunk = 'Normal';
wolffd@0: specif.combinechunk = 'Concat';
wolffd@0: specif.extensive = 1;
wolffd@0: 
wolffd@0: varargout = mirfunction(@mirenvelope,orig,varargin,nargout,specif,@init,@main);
wolffd@0: 
wolffd@0: 
wolffd@0: function [x type] = init(x,option)
wolffd@0: type = 'mirenvelope';
wolffd@0: if isamir(x,'mirscalar')
wolffd@0:     return
wolffd@0: end
wolffd@0: if ischar(option.presel) && strcmpi(option.presel,'Klapuri06')
wolffd@0:     option.method = 'Spectro';
wolffd@0: end
wolffd@0: if not(isamir(x,'mirenvelope'))
wolffd@0:     if strcmpi(option.method,'Filter')
wolffd@0:         if isnan(option.zp)
wolffd@0:             if strcmpi(option.filter,'IIR')
wolffd@0:                 option.zp = 1;
wolffd@0:             else
wolffd@0:                 option.zp = 0;
wolffd@0:             end
wolffd@0:         end
wolffd@0:         if option.zp == 1
wolffd@0:             x = mirenvelope(x,'ZeroPhase',2,'Down',1,...
wolffd@0:                               'Tau',option.tau,'PreDecim',option.decim);
wolffd@0:         end
wolffd@0:     elseif strcmpi(option.method,'Spectro')
wolffd@0:         x = mirspectrum(x,'Frame',option.frame.length.val,...
wolffd@0:                                   option.frame.length.unit,...
wolffd@0:                                   option.frame.hop.val,...
wolffd@0:                                   option.frame.hop.unit,...
wolffd@0:                                   'Window','hanning',...
wolffd@0:                                   option.band,'Power');
wolffd@0: 
wolffd@0:     end
wolffd@0: end
wolffd@0: 
wolffd@0: 
wolffd@0: function e = main(orig,option,postoption)
wolffd@0: if iscell(orig)
wolffd@0:     orig = orig{1};
wolffd@0: end
wolffd@0: if isamir(orig,'mirscalar')
wolffd@0:     d = get(orig,'Data');
wolffd@0:     fp = get(orig,'FramePos');
wolffd@0:     for i = 1:length(d)
wolffd@0:         for j = 1:length(d{i})
wolffd@0:             d{i}{j} = reshape(d{i}{j},size(d{i}{j},2),1,size(d{i}{j},3));
wolffd@0:             p{i}{j} = mean(fp{i}{j})';
wolffd@0:         end
wolffd@0:     end
wolffd@0:     e.downsampl = 0;
wolffd@0:     e.hwr = 0;
wolffd@0:     e.diff = 0;
wolffd@0:     e.method = 'Spectro';
wolffd@0:     e.phase = {{}};
wolffd@0:     e = class(e,'mirenvelope',mirtemporal(orig));
wolffd@0:     e = set(e,'Title','Envelope','Data',d,'Pos',p,'FramePos',{{}});
wolffd@0:     postoption.trim = 0;
wolffd@0:     postoption.ds = 0;
wolffd@0:     e = post(e,postoption);
wolffd@0:     return
wolffd@0: end
wolffd@0: if isfield(option,'presel') && ischar(option.presel) && ...
wolffd@0:         strcmpi(option.presel,'Klapuri06')
wolffd@0:     option.method = 'Spectro';
wolffd@0:     postoption.up = 2;
wolffd@0:     postoption.mu = 100;
wolffd@0:     postoption.diffhwr = 1;
wolffd@0:     postoption.lambda = .8;
wolffd@0: end
wolffd@0: if isfield(postoption,'ds') && isnan(postoption.ds)
wolffd@0:     if option.decim
wolffd@0:         postoption.ds = 0;
wolffd@0:     else
wolffd@0:         postoption.ds = 16;
wolffd@0:     end
wolffd@0: end
wolffd@0: if not(isfield(option,'filter')) || not(ischar(option.filter))
wolffd@0:     e = post(orig,postoption);
wolffd@0: elseif strcmpi(option.method,'Spectro')
wolffd@0:     d = get(orig,'Data');
wolffd@0:     fp = get(orig,'FramePos');
wolffd@0:     sr = get(orig,'Sampling');
wolffd@0:     ch = get(orig,'Channels');
wolffd@0:     ph = get(orig,'Phase');
wolffd@0:     for h = 1:length(d)
wolffd@0:         sr{h} = 0;
wolffd@0:         for i = 1:length(d{h})
wolffd@0:             if size(d{h}{i},3)>1 % Already in bands (channels in 3d dim)
wolffd@0:                 d{h}{i} = permute(sum(d{h}{i}),[2 1 3]);
wolffd@0:                 ph{h}{i} = permute(ph{h}{i},[2 1 3]);
wolffd@0:             else % Simple spectrogram, frequency range sent to 3d dim
wolffd@0:                 d{h}{i} = permute(d{h}{i},[2 3 1]);
wolffd@0:                 ph{h}{i} = permute(ph{h}{i},[2 3 1]);
wolffd@0:             end
wolffd@0:             p{h}{i} = mean(fp{h}{i})';
wolffd@0:             if not(sr{h}) && size(fp{h}{i},2)>1
wolffd@0:                 sr{h} = 1/(fp{h}{i}(1,2)-fp{h}{i}(1,1));
wolffd@0:             end
wolffd@0:         end
wolffd@0:         if not(sr{h})
wolffd@0:             warning('WARNING IN MIRENVELOPE: The frame decomposition did not succeed. Either the input is of too short duration, or the chunk size is too low.');
wolffd@0:         end
wolffd@0:         ch{h} = (1:size(d{h}{1},3))';
wolffd@0:     end
wolffd@0:     e.downsampl = 0;
wolffd@0:     e.hwr = 0;
wolffd@0:     e.diff = 0;
wolffd@0:     %e.todelete = {};
wolffd@0:     e.method = 'Spectro';
wolffd@0:     e.phase = ph;
wolffd@0:     e = class(e,'mirenvelope',mirtemporal(orig));
wolffd@0:     e = set(e,'Title','Envelope','Data',d,'Pos',p,...
wolffd@0:               'Sampling',sr,'Channels',ch);
wolffd@0:     postoption.trim = 0;
wolffd@0:     postoption.ds = 0;
wolffd@0:     e = post(e,postoption);
wolffd@0: else
wolffd@0:     if isnan(option.zp)
wolffd@0:         if strcmpi(option.filter,'IIR')
wolffd@0:             option.zp = 1;
wolffd@0:         else
wolffd@0:             option.zp = 0;
wolffd@0:         end
wolffd@0:     end
wolffd@0:     if option.zp == 1
wolffd@0:         option.decim = 0;
wolffd@0:     end
wolffd@0:     e.downsampl = 1;
wolffd@0:     e.hwr = 0;
wolffd@0:     e.diff = 0;
wolffd@0:     e.method = option.filter;
wolffd@0:     e.phase = {};
wolffd@0:     e = class(e,'mirenvelope',mirtemporal(orig));
wolffd@0:     e = purgedata(e);
wolffd@0:     e = set(e,'Title','Envelope');
wolffd@0:     sig = get(e,'Data');
wolffd@0:     x = get(e,'Pos');
wolffd@0:     sr = get(e,'Sampling');
wolffd@0:     disp('Extracting envelope...')
wolffd@0:     d = cell(1,length(sig));
wolffd@0:     [state e] = gettmp(orig,e);
wolffd@0:     for k = 1:length(sig)
wolffd@0:         if option.decim
wolffd@0:             sr{k} = sr{k}/option.decim;
wolffd@0:         end
wolffd@0:         if strcmpi(option.filter,'IIR')
wolffd@0:             a2 = exp(-1/(option.tau*sr{k})); % filter coefficient 
wolffd@0:             a = [1 -a2];
wolffd@0:             b = 1-a2;
wolffd@0:         elseif strcmpi(option.filter,'HalfHann')
wolffd@0:             a = 1;
wolffd@0:             b = hann(sr{k}*.4);
wolffd@0:             b = b(ceil(length(b)/2):end);
wolffd@0:         end
wolffd@0:         d{k} = cell(1,length(sig{k}));
wolffd@0:         for i = 1:length(sig{k})
wolffd@0:             sigi = sig{k}{i};
wolffd@0:             if option.zp == 2
wolffd@0:                 sigi = flipdim(sigi,1);
wolffd@0:             end
wolffd@0:             if option.hilb
wolffd@0:                 try
wolffd@0:                     for h = 1:size(sigi,2)
wolffd@0:                         for j = 1:size(sigi,3)
wolffd@0:                             sigi(:,h,j) = hilbert(sigi(:,h,j));
wolffd@0:                         end
wolffd@0:                     end
wolffd@0:                 catch
wolffd@0:                     disp('Signal Processing Toolbox does not seem to be installed. No Hilbert transform.');
wolffd@0:                 end 
wolffd@0:             end
wolffd@0:             sigi = abs(sigi);
wolffd@0:             
wolffd@0:             if option.decim
wolffd@0:                 dsigi = zeros(ceil(size(sigi,1)/option.decim),...
wolffd@0:                               size(sigi,2),size(sigi,3));
wolffd@0:                 for f = 1:size(sigi,2)
wolffd@0:                     for c = 1:size(sigi,3)
wolffd@0:                         dsigi(:,f,c) = decimate(sigi(:,f,c),option.decim);
wolffd@0:                     end
wolffd@0:                 end
wolffd@0:                 sigi = dsigi;
wolffd@0:                 clear dsigi
wolffd@0:                 x{k}{i} = x{k}{i}(1:option.decim:end,:,:);
wolffd@0:             end
wolffd@0:             
wolffd@0:             % tmp = filtfilt(1-a,[1 -a],sigi); % zero-phase IIR filter for smoothing the envelope
wolffd@0: 
wolffd@0:             % Manual filtfilt
wolffd@0:             emptystate = isempty(state);
wolffd@0:             tmp = zeros(size(sigi));
wolffd@0:             for c = 1:size(sigi,3)
wolffd@0:                 if emptystate
wolffd@0:                     [tmp(:,:,c) state(:,c,1)] = filter(b,a,sigi(:,:,c));
wolffd@0:                 else
wolffd@0:                     [tmp(:,:,c) state(:,c,1)] = filter(b,a,sigi(:,:,c),...
wolffd@0:                                                         state(:,c,1));
wolffd@0:                 end
wolffd@0:             end
wolffd@0:             
wolffd@0:             tmp = max(tmp,0); % For security reason...
wolffd@0:             if option.zp == 2
wolffd@0:                 tmp = flipdim(tmp,1);
wolffd@0:             end
wolffd@0:             d{k}{i} = tmp;
wolffd@0:             %td{k} = round(option.tau*sr{k}*1.5); 
wolffd@0:         end
wolffd@0:     end
wolffd@0:     e = set(e,'Data',d,'Pos',x,'Sampling',sr); %,'ToDelete',td
wolffd@0:     e = settmp(e,state);
wolffd@0:     if not(option.zp == 2)
wolffd@0:         e = post(e,postoption);
wolffd@0:     end
wolffd@0: end
wolffd@0: if isfield(option,'presel') && ischar(option.presel) && ...
wolffd@0:         strcmpi(option.presel,'Klapuri06')
wolffd@0:     e = mirsum(e,'Adjacent',10);
wolffd@0: end
wolffd@0: 
wolffd@0:     
wolffd@0: function e = post(e,postoption)
wolffd@0: if isempty(postoption)
wolffd@0:     return
wolffd@0: end
wolffd@0: if isfield(postoption,'lambda') && not(postoption.lambda)
wolffd@0:     postoption.lambda = 1;
wolffd@0: end
wolffd@0: d = get(e,'Data');
wolffd@0: tp = get(e,'Time');
wolffd@0: sr = get(e,'Sampling');
wolffd@0: ds = get(e,'DownSampling');
wolffd@0: ph = get(e,'Phase');
wolffd@0: for k = 1:length(d)
wolffd@0:     if isfield(postoption,'sampling')
wolffd@0:         if postoption.sampling
wolffd@0:             newsr = postoption.sampling;
wolffd@0:         elseif isfield(postoption,'ds') && postoption.ds>1
wolffd@0:             newsr = sr{k}/postoption.ds;
wolffd@0:         else
wolffd@0:             newsr = sr{k};
wolffd@0:         end
wolffd@0:     end
wolffd@0:     if isfield(postoption,'up') && postoption.up
wolffd@0:         [z,p,gain] = butter(6,10/newsr/postoption.up*2,'low');
wolffd@0:         [sos,g] = zp2sos(z,p,gain);
wolffd@0:         Hd = dfilt.df2tsos(sos,g);
wolffd@0:     end
wolffd@0:     for i = 1:length(d{k})
wolffd@0:         if isfield(postoption,'sampling') && postoption.sampling
wolffd@0:             if and(sr{k}, not(sr{k} == postoption.sampling))
wolffd@0:                 dk = d{k}{i};
wolffd@0:                 for j = 1:size(dk,3)
wolffd@0:                     if not(sr{k} == round(sr{k}))
wolffd@0:                         mirerror('mirenvelope','The ''Sampling'' postoption cannot be used after using the ''Down'' postoption.');
wolffd@0:                     end
wolffd@0:                     rk(:,:,j) = resample(dk(:,:,j),postoption.sampling,sr{k});
wolffd@0:                 end
wolffd@0:                 d{k}{i} = rk;
wolffd@0:                 tp{k}{i} = repmat((0:size(d{k}{i},1)-1)',...
wolffd@0:                                     [1 1 size(tp{k}{i},3)])...
wolffd@0:                             /postoption.sampling + tp{k}{i}(1,:,:);
wolffd@0:                 if not(iscell(ds))
wolffd@0:                     ds = cell(length(d));
wolffd@0:                 end
wolffd@0:                 ds{k} = round(sr{k}/postoption.sampling);
wolffd@0:             end
wolffd@0:         elseif isfield(postoption,'ds') && postoption.ds>1
wolffd@0:             if not(postoption.ds == round(postoption.ds))
wolffd@0:                 mirerror('mirenvelope','The ''Down'' sampling rate should be an integer.');
wolffd@0:             end
wolffd@0:             ds = postoption.ds;
wolffd@0:             tp{k}{i} = tp{k}{i}(1:ds:end,:,:); % Downsampling...
wolffd@0:             d{k}{i} = d{k}{i}(1:ds:end,:,:);
wolffd@0:         end
wolffd@0:         if isfield(postoption,'sampling')
wolffd@0:             if not(strcmpi(e.method,'Spectro')) && postoption.trim 
wolffd@0:                 tdk = round(newsr*.1); 
wolffd@0:                 d{k}{i}(1:tdk,:,:) = repmat(d{k}{i}(tdk,:,:),[tdk,1,1]); 
wolffd@0:                 d{k}{i}(end-tdk+1:end,:,:) = repmat(d{k}{i}(end-tdk,:,:),[tdk,1,1]);
wolffd@0:             end
wolffd@0:             if postoption.log
wolffd@0:                 d{k}{i} = log10(d{k}{i});
wolffd@0:             end
wolffd@0:             if postoption.mu
wolffd@0:                 dki = max(0,d{k}{i});
wolffd@0:                 mu = postoption.mu;
wolffd@0:                 dki = log(1+mu*dki)/log(1+mu);
wolffd@0:                 dki(~isfinite(d{k}{i})) = NaN;
wolffd@0:                 d{k}{i} = dki;
wolffd@0:             end
wolffd@0:             if postoption.power
wolffd@0:                 d{k}{i} = d{k}{i}.^2;
wolffd@0:             end
wolffd@0:             if postoption.up
wolffd@0:                 dki = zeros(size(d{k}{i},1).*postoption.up,...
wolffd@0:                             size(d{k}{i},2),size(d{k}{i},3));
wolffd@0:                 dki(1:postoption.up:end,:,:) = d{k}{i};
wolffd@0:                 dki = filter(Hd,[dki;...
wolffd@0:                                 zeros(6,size(d{k}{i},2),size(d{k}{i},3))]);
wolffd@0:                 d{k}{i} = dki(1+ceil(6/2):end-floor(6/2),:,:);
wolffd@0:                 tki = zeros(size(tp{k}{i},1).*postoption.up,...
wolffd@0:                             size(tp{k}{i},2),...
wolffd@0:                             size(tp{k}{i},3));
wolffd@0:                 dt = repmat((tp{k}{i}(2)-tp{k}{i}(1))...
wolffd@0:                                 /postoption.up,...
wolffd@0:                             [size(tp{k}{i},1),1,1]);
wolffd@0:                 for j = 1:postoption.up
wolffd@0:                     tki(j:postoption.up:end,:,:) = tp{k}{i}+dt*(j-1);
wolffd@0:                 end
wolffd@0:                 tp{k}{i} = tki;
wolffd@0:                 newsr = sr{k}*postoption.up;
wolffd@0:             end
wolffd@0:             if (postoption.diffhwr || postoption.diff) && ...
wolffd@0:                     not(get(e,'Diff'))
wolffd@0:                 tp{k}{i} = tp{k}{i}(1:end-1,:,:);
wolffd@0:                 order = max(postoption.diffhwr,postoption.diff);
wolffd@0:                 if postoption.complex
wolffd@0:                     dph = diff(ph{k}{i},2);
wolffd@0:                     dph = dph/(2*pi);% - round(dph/(2*pi));
wolffd@0:                     ddki = sqrt(d{k}{i}(3:end,:,:).^2 + d{k}{i}(2:end-1,:,:).^2 ...
wolffd@0:                                               - 2.*d{k}{i}(3:end,:,:)...
wolffd@0:                                                  .*d{k}{i}(2:end-1,:,:)...
wolffd@0:                                                  .*cos(dph));
wolffd@0:                     d{k}{i} = d{k}{i}(2:end,:,:); 
wolffd@0:                     tp{k}{i} = tp{k}{i}(2:end,:,:);
wolffd@0:                 elseif order == 1
wolffd@0:                     ddki = diff(d{k}{i},1,1);
wolffd@0:                 else
wolffd@0:                     b = firls(order,[0 0.9],[0 0.9*pi],'differentiator');
wolffd@0:                     ddki = filter(b,1,...
wolffd@0:                         [repmat(d{k}{i}(1,:,:),[order,1,1]);...
wolffd@0:                          d{k}{i};...
wolffd@0:                          repmat(d{k}{i}(end,:,:),[order,1,1])]);
wolffd@0:                     ddki = ddki(order+1:end-order-1,:,:);
wolffd@0:                 end
wolffd@0:                 if postoption.diffhwr
wolffd@0:                     ddki = hwr(ddki);
wolffd@0:                 end
wolffd@0:                 d{k}{i} = (1-postoption.lambda)*d{k}{i}(1:end-1,:,:)...
wolffd@0:                             + postoption.lambda*sr{k}/10*ddki;
wolffd@0:             end
wolffd@0:             if postoption.aver
wolffd@0:                 y = filter(ones(1,postoption.aver),1,...
wolffd@0:                             [d{k}{i};zeros(postoption.aver,...
wolffd@0:                                            size(d{k}{i},2),...
wolffd@0:                                            size(d{k}{i},3))]);
wolffd@0:                 d{k}{i} = y(1+ceil(postoption.aver/2):...
wolffd@0:                              end-floor(postoption.aver/2),:,:);
wolffd@0:             end
wolffd@0:             if postoption.gauss
wolffd@0:                 sigma = postoption.gauss;
wolffd@0:                 gauss = 1/sigma/2/pi...
wolffd@0:                         *exp(- (-4*sigma:4*sigma).^2 /2/sigma^2);
wolffd@0:                 y = filter(gauss,1,[d{k}{i};zeros(4*sigma,1,size(d{k}{i},3))]);
wolffd@0:                 y = y(4*sigma:end,:,:);
wolffd@0:                 d{k}{i} = y(1:size(d{k}{i},1),:,:);
wolffd@0:             end
wolffd@0:             if postoption.chwr
wolffd@0:                 d{k}{i} = center(d{k}{i});
wolffd@0:                 d{k}{i} = hwr(d{k}{i});
wolffd@0:             end
wolffd@0:             if postoption.hwr
wolffd@0:                 d{k}{i} = hwr(d{k}{i});
wolffd@0:             end
wolffd@0:             if postoption.c
wolffd@0:                 d{k}{i} = center(d{k}{i});
wolffd@0:             end    
wolffd@0:             if postoption.norm
wolffd@0:                 d{k}{i} = d{k}{i}./repmat(max(abs(d{k}{i})),...
wolffd@0:                                          [size(d{k}{i},1),1,1]);
wolffd@0:             end        
wolffd@0:         end
wolffd@0:     end
wolffd@0:     if isfield(postoption,'sampling')
wolffd@0:         sr{k} = newsr;
wolffd@0:     end
wolffd@0: end
wolffd@0: if isfield(postoption,'ds') && postoption.ds>1
wolffd@0:     e = set(e,'DownSampling',postoption.ds,'Sampling',sr);
wolffd@0: elseif isfield(postoption,'sampling') && postoption.sampling
wolffd@0:     e = set(e,'DownSampling',ds,'Sampling',sr);
wolffd@0: elseif isfield(postoption,'up') && postoption.up
wolffd@0:     e = set(e,'Sampling',sr);
wolffd@0: end
wolffd@0: if isfield(postoption,'sampling')
wolffd@0:     if postoption.hwr
wolffd@0:         e = set(e,'Halfwave',1);
wolffd@0:     end
wolffd@0:     if postoption.diff
wolffd@0:         e = set(e,'Diff',1,'Halfwave',0,'Title','Differentiated envelope');
wolffd@0:     end
wolffd@0:     if postoption.diffhwr
wolffd@0:         e = set(e,'Diff',1,'Halfwave',1,'Centered',0);
wolffd@0:     end
wolffd@0:     if postoption.c
wolffd@0:         e = set(e,'Centered',1);
wolffd@0:     end
wolffd@0:     if postoption.chwr
wolffd@0:         e = set(e,'Halfwave',1,'Centered',1);
wolffd@0:     end
wolffd@0: end
wolffd@0: e = set(e,'Data',d,'Time',tp); 
wolffd@0: %if isfield(postoption,'frame') && isstruct(postoption.frame)
wolffd@0: %    e = mirframenow(e,postoption);
wolffd@0: %end