line source
+ − % FILTSUBTRACT subtracts a frequency dependent noise contribution from an input ao.
+ − %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ − %
+ − % DESCRIPTION: FILTSUBTRACT subtracts a frequency dependent noise contribution from an input ao.
+ − % The method computes the transfer function between both input AOs and
+ − % fits a miir model to it. The frequency band is applied is set by a
+ − % threshold in the coherence that the user defines as an input
+ − % parameter.
+ − %
+ − % CALL: c = filtSubtract(a,b pl)
+ − %
+ − % INPUTS: a - AO from where subtract linear contributions
+ − % b - AOs with noise contributions
+ − % pl - parameter list (see below)
+ − %
+ − % OUTPUTs: c - output AO with contributions subtracted (tsdata)
+ − %
+ − % <a href="matlab:utils.helper.displayMethodInfo('ao', 'filtSubtract')">Parameters Description</a>
+ − %
+ − %
+ − % VERSION: $Id: filtSubtract.m,v 1.19 2011/05/22 22:26:56 mauro Exp $
+ − %
+ − % TODO: handling errors
+ − % split by coherence function
+ − %
+ − %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ −
+ − function varargout = filtSubtract(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);
+ −
+ − % Method can not be used as a modifier
+ − if nargout == 0
+ − error('### filtSubtract cannot be used as a modifier. Please give an output variable.');
+ − end
+ −
+ − % 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 = utils.helper.collect_objects(varargin(:), 'plist', in_names);
+ −
+ − % Only two inputs ao's accepted
+ − if numel(as) ~= 2
+ − error('### filtSubtract only accepts two inputs AOs.');
+ − end
+ −
+ − % Decide on a deep copy or a modify
+ − bs = copy(as, nargout);
+ −
+ − % combine plists
+ − pl = combine(pl, getDefaultPlist());
+ −
+ − % get parameters
+ − times = find(pl,'times');
+ − filttimes = find(pl,'times postfilter');
+ − fs = find(pl,'fs');
+ − filt = find(pl,'filt');
+ − fspl = find(pl,'frequencies');
+ −
+ − % resample and consolidate
+ − if isempty(fs)
+ − bs1r = bs(1);
+ − pl = pset(pl,'fs', bs(1).fs);
+ − else
+ − if fs ~= bs(1).fs
+ − bs1r = resample(bs(1), plist('fsout',fs));
+ − else
+ − bs1r = bs(1);
+ − end
+ − end
+ −
+ − % consolidate and split
+ − c = consolidate(bs1r, bs(2), pl);
+ − if ~isempty(times)
+ − cs = split(c(1), plist('times', times));
+ − cn = split(c(2), plist('times', times));
+ − else
+ − cs = c(1);
+ − cn = c(2);
+ − end
+ −
+ − if isempty(filt)
+ − % Transfer functions
+ − tf = ltfe(cn, cs, pl);
+ − % split transfer function to relevant frequencies
+ − if ~isempty(fspl)
+ − tf_spl = split(tf, plist('frequencies', [fspl(1) fspl(2)]));
+ − else
+ − tf_spl = tf;
+ − end
+ −
+ − % get filter from transfer function
+ − fp = zDomainFit(tf_spl, pl);
+ − fp.filters.setIunits(cn.yunits);
+ − fp.filters.setOunits(cs.yunits);
+ − else
+ − fp = filt;
+ − end
+ −
+ − % get noise contribution
+ − cs_cont = filter(cn,fp);
+ − cs_cont.simplifyYunits();
+ −
+ − % remove filter transient
+ − if ~isempty(filttimes)
+ − cs_cont = split(cs_cont, plist('times', filttimes));
+ − cs = split(cs, plist('times', filttimes));
+ − end
+ −
+ − % subtraction
+ − cs_subt = detrend(cs) - detrend(cs_cont);
+ −
+ − % new tsdata
+ − fsd = tsdata(cs_subt.x, cs_subt.y, cs_subt.fs);
+ − % make output analysis object
+ − cs = ao(fsd);
+ − % set name
+ − cs.name = sprintf('filtSubtract(%s)', ao_invars{1});
+ − % set units
+ − cs.setYunits(cs_subt.yunits);
+ − % t0
+ − if ~isempty(times) && ~isempty(filttimes)
+ − cs.setT0(bs(1).t0 + times(1) + filttimes(1));
+ − elseif isempty(times) && ~isempty(filttimes)
+ − cs.setT0(bs(1).t0 + filttimes(1));
+ − elseif ~isempty(times) && isempty(filttimes)
+ − cs.setT0(bs(1).t0 + times(1));
+ − else
+ − cs.setT0(bs(1).t0);
+ − end
+ − % Add procinfo
+ − cs.procinfo = plist('filter', fp);
+ − % Add history
+ − cs.addHistory(getInfo('None'), pl, ao_invars, [bs(:).hist]);
+ − % Set output
+ − varargout{1} = cs;
+ −
+ − end
+ −
+ −
+ − %--------------------------------------------------------------------------
+ − % Get Info Object
+ − %--------------------------------------------------------------------------
+ − function ii = getInfo(varargin)
+ − if nargin == 1 && strcmpi(varargin{1}, 'None')
+ − sets = {};
+ − pl = [];
+ − else
+ − sets = {'Default'};
+ − pl = getDefaultPlist();
+ − end
+ − % Build info object
+ − ii = minfo(mfilename, 'ao', 'ltpda', utils.const.categories.sigproc, '$Id: filtSubtract.m,v 1.19 2011/05/22 22:26:56 mauro Exp $', sets, pl);
+ − end
+ −
+ −
+ − %--------------------------------------------------------------------------
+ − % Get Default Plist
+ − %--------------------------------------------------------------------------
+ − function plout = getDefaultPlist()
+ − persistent pl;
+ − if ~exist('pl', 'var') || isempty(pl)
+ − pl = buildplist();
+ − end
+ − plout = pl;
+ − end
+ −
+ − function pl = buildplist()
+ −
+ − pl = plist();
+ −
+ − % fs
+ − p = param({'fs','target sampling frequency to resample the data.'}, paramValue.EMPTY_DOUBLE);
+ − pl.append(p);
+ −
+ − % times
+ − p = param({'times', 'selects the interval where the subtraction is applied.'},...
+ − paramValue.EMPTY_STRING);
+ − pl.append(p);
+ −
+ − % times
+ − p = param({'times postfilter', 'selects the filter transient intervals to be removed.'},...
+ − paramValue.EMPTY_STRING);
+ − pl.append(p);
+ −
+ − % filtspl
+ − p = param({'frequencies', 'selects the frequency band where the transfer<br>'...
+ − 'function is fitted.'}, paramValue.EMPTY_STRING);
+ − pl.append(p);
+ −
+ − % filt
+ − p = param({'filt', ['a miir/mfir object which will be used as a<br>'...
+ − 'transfer function. If this option is selected<br>'...
+ − ' the fit is avoided.']}, paramValue.EMPTY_STRING);
+ − pl.append(p);
+ −
+ − end
+ −