% FILTFILT overrides the filtfilt function for analysis objects.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% DESCRIPTION: FILTFILT overrides the filtfilt function for analysis objects.
%              Applies the input digital IIR filter to the input analysis object
%              forwards and backwards. If the input analysis object contains a
%              time-series (tsdata) then the filter is applied using the normal
%              recursion algorithm. The output analysis object contains a tsdata
%              object.
%
%              If the input analysis object contains a frequency-series (fsdata)
%              then the response of the filter is computed and then multiplied
%              with the input frequency series. The output analysis object
%              contains a frequency series.
%
% CALL:        >> [b, filt] = filtfilt(a,pl)
%              >> b = filtfilt(a,pl)
%
% INPUTS:      pl   - a parameter list
%              a    - input analysis object
%
% OUTPUTS:     filt - a copy of the input filter object with the
%                     history values filled in.
%              b    - output analysis object containing the filtered data.
%
% <a href="matlab:utils.helper.displayMethodInfo('ao', 'filtfilt')">Parameters Description</a>
%
% VERSION:     $Id: filtfilt.m,v 1.48 2011/05/12 13:20:42 luigi Exp $
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function varargout = filtfilt(varargin)
  
  % Check if this is a call for parameters
  if utils.helper.isinfocall(varargin{:})
    varargout{1} = getInfo(varargin{3});
    return
  end
  
  import utils.const.*
  utils.helper.msg(msg.PROC3, 'running %s/%s', mfilename('class'), mfilename);
  
  % Collect input variable names
  in_names = cell(size(varargin));
  for ii = 1:nargin,in_names{ii} = inputname(ii);end
  
  % Collect all AOs and plists
  [as, ao_invars] = utils.helper.collect_objects(varargin(:), 'ao', in_names);
  [pl, pl_invars] = utils.helper.collect_objects(varargin(:), 'plist', in_names);
  [fbobjs, f_invars] = utils.helper.collect_objects(varargin(:), 'filterbank', in_names);
  [fobj, f_invars] = utils.helper.collect_objects(varargin(:), 'ltpda_filter', in_names);
  
  % Make copies or handles to inputs
  bs   = copy(as, nargout);
  
  % combine plists
  pl = parse(pl, getDefaultPlist());
  
  
  if isempty(fobj) && isempty(fbobjs)
    fobj        = find(pl, 'filter');
    if isa(fobj, 'filterbank')
    fbobjs = fobj;
    end
    f_invars{1} = class(fobj);
  end
  
  % check inputs
  if ~isa(fobj, 'miir') && ~isa(fobj, 'mfir') && ~isa(fbobjs, 'filterbank')
    error('### the filter input should be either an miir/mfir object or a filterbank object.');
  end
  
  fobj_out = [];
  fp = [];
  
  % Loop over AOs
  for jj = 1:numel(bs)
    
    % Copy filter so we can change it
    if ~isempty(fobj)
      fp = copy(fobj, 1);
    elseif ~isempty(fbobjs)
      fp = copy(fbobjs, 1);
    end
    % keep the history to suppress the history of the intermediate steps
    inhist = bs(jj).hist;
    
    if isa(bs(jj).data, 'tsdata')
      %------------------------------------------------------------------------
      %------------------------   Time-series filter   ------------------------
      %------------------------------------------------------------------------
      
      %%%%%%%%%%%%%%%%%%%%%%%%%%%     filter    %%%%%%%%%%%%%%%%%%%%%%%%%%%
      if isa(fp,'ltpda_filter')
        % redesign if needed
        fs = bs(jj).data.fs;
        if fs ~= fp.fs 
          warning('!!! Filter is designed for a different sample rate of data.');
          % Adjust/redesign if this is a standard filter
          fp = redesign(fp, fs);
        end
        % apply filter
        if isa(fp, 'miir')
          bs(jj).data.y = filtfilt(fp.a, fp.b, bs(jj).data.y);
        elseif isa(fp, 'mfir');
          bs(jj).data.y = filtfilt(fp.a, 1, bs(jj).data.y);
        else
          error('### Unknown filter object [%s]', class(fp));
        end
        % set y-units = yunits.*ounits./iunits
        bs(jj).data.setYunits(bs(jj).data.yunits.*fp.ounits./fp.iunits);
        
        %%%%%%%%%%%%%%%%%%%%%%%%%%%     filter bank     %%%%%%%%%%%%%%%%%%%%%%%%%%%
      elseif isa(fp,'filterbank')
        % redesign not implemented for filterbank
        %%%
        % utils.math routine to apply filtfilt properly
        bs(jj).data.y = utils.math.filtfilt_filterbank(bs(jj),fp);
        % not setting units yet
        %%%
      end
      
    elseif isa(bs(jj).data, 'fsdata')
      %------------------------------------------------------------------------
      %----------------------   Frequency-series filter   ---------------------
      %------------------------------------------------------------------------
      
      % apply filter
      fil_resp = resp(fp, plist('f', bs(jj)));
      bs(jj)    = bs(jj).*fil_resp.*conj(fil_resp);
    else
      error('### unknown data type.');
    end
    
    
    % name for this object
    bs(jj).name = sprintf('%s(%s)', fp.name, ao_invars{jj});
    % add history
    bs(jj).addHistory(getInfo('None'), pl, ao_invars(jj), [bs(jj).hist fp.hist]);
    % clear errors
    bs(jj).clearErrors;
    % Collect the filters for the output
    fobj_out = [fobj_out, fp];
  end
  
  % Set outputs
  if nargout == 1
    varargout{1} = bs;
  elseif nargout == 2
    varargout{1} = bs;
    varargout{2} = fobj_out;
  elseif nargout > 2
    error('### wrong number of output arguments.');
  end
end

%--------------------------------------------------------------------------
% Get Info Object
%--------------------------------------------------------------------------
function ii = getInfo(varargin)
  if nargin == 1 && strcmpi(varargin{1}, 'None')
    sets = {};
    pls  = [];
  else
    sets = {'Default'};
    pls  = getDefaultPlist;
  end
  % Build info object
  ii = minfo(mfilename, 'ao', 'ltpda', utils.const.categories.sigproc, '$Id: filtfilt.m,v 1.48 2011/05/12 13:20:42 luigi Exp $', sets, pls);
end

%--------------------------------------------------------------------------
% Get Default Plist
%--------------------------------------------------------------------------
function plout = getDefaultPlist()
  persistent pl;  
  if exist('pl', 'var')==0 || isempty(pl)
    pl = buildplist();
  end
  plout = pl;  
end

function pl = buildplist()
  pl = plist({'filter', 'The filter to apply to the data.'},  paramValue.EMPTY_STRING);
end