--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/m-toolbox/classes/@ao/lscov.m	Wed Nov 23 19:22:13 2011 +0100
@@ -0,0 +1,227 @@
+% LSCOV is a wrapper for MATLAB's lscov function.
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% DESCRIPTION: LSCOV is a wrapper for MATLAB's lscov function. It solves a
+% set of linear equations by performing a linear least-squares fit. It
+% solves the problem
+%
+%        Y = HX
+%
+%   where X are the parameters, Y the measurements, and H the linear
+%   equations relating the two.
+%
+% CALL:        X = lscov([C1 C2 ... CN], Y, pl)
+%              X = lscov(C1,C2,C3,...,CN, Y, pl)
+%
+% INPUTS:      C1...CN - AOs which represent the columns of H.
+%              Y       - AO which represents the measurement set
+%
+% Note: the length of the vectors in Ci and Y must be the same.
+% Note: the last input AO is taken as Y.
+%
+%              pl - parameter list (see below)
+%
+% OUTPUTs:     X  - A pest object with fields: 
+%                   y   - the N fitting coefficients to y_i 
+%                   dy  - the parameters' standard deviations (lscov 'STDX' vector)
+%                   cov - the parameters' covariance matrix (lscov 'COV' vector)
+%                   
+% The procinfo field of the output PEST object is filled with the following key/value
+% pairs:
+%     'MSE' - the mean-squared errors
+%
+% <a href="matlab:utils.helper.displayMethodInfo('ao', 'lscov')">Parameters Description</a>
+%
+% VERSION:     $Id: lscov.m,v 1.42 2011/04/08 08:56:17 hewitson Exp $
+%
+% EXAMPLES:
+%
+% % 1) Determine the coefficients of a linear combination of noises:
+%
+% % Make some data
+%  fs    = 10;
+%  nsecs = 10;
+%  B1 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
+%  B2 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
+%  B3 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
+%  n  = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'm'));
+%  c = [ao(1,plist('yunits','m/T')) ao(2,plist('yunits','m/T')) ao(3,plist('yunits','m T^-1'))];
+%  y = c(1)*B1 + c(2)*B2 + c(3)*B3 + n;
+%  y.simplifyYunits;
+%  % Get a fit for c
+%  p_s = lscov(B1, B2, B3, y);
+%  % do linear combination: using lincom
+%  yfit1 = lincom(B1, B2, B3, p_s);
+%  yfit1.simplifyYunits;
+%  % do linear combination: using eval
+%  yfit2 = p_s.eval(B1, B2, B3);
+% 
+%  % Plot (compare data with fit)
+%  iplot(y, yfit1, yfit2, plist('Linestyles', {'-','--'}))
+%
+% % 2) Determine the coefficients of a linear combination of noises:
+%
+% % Make some data
+%  fs    = 10;
+%  nsecs = 10;
+%  x1 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'T'));
+%  x2 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'm'));
+%  x3 = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'C'));
+%  n  = ao(plist('tsfcn', 'randn(size(t))', 'fs', fs, 'nsecs', nsecs, 'yunits', 'm'));
+%  c = [ao(1,plist('yunits','m/T')) ao(2,plist('yunits','m/m')) ao(3,plist('yunits','m C^-1'))];
+%  y = c(1)*x1 + c(2)*x2 + c(3)*x3 + n;
+%  y.simplifyYunits;
+%  % Get a fit for c
+%  p_m = lscov(x1, x2, x3, y);
+%  % do linear combination: using lincom
+%  yfit1 = lincom(x1, x2, x3, p_m);
+%  % do linear combination: using eval
+%  pl_split = plist('times', [1 5]);
+%  yfit2 = p_m.eval(plist('Xdata', {split(x1, pl_split), split(x2, pl_split), split(x3, pl_split)}));
+%  % Plot (compare data with fit)
+%  iplot(y, yfit1, yfit2, plist('Linestyles', {'-','--'}))
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+function varargout = lscov(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
+  [A, ao_invars] = utils.helper.collect_objects(varargin(:), 'ao', in_names);
+  pl             = utils.helper.collect_objects(varargin(:), 'plist', in_names);
+  
+  if numel(A) < 2
+    error('### lscov needs at least 2 inputs AOs');
+  end
+  
+  if nargout == 0
+    error('### lscov can not be used as a modifier method. Please give at least one output');
+  end
+  
+  % combine plists
+  pl = parse(pl, getDefaultPlist());
+  
+  % collect inputs names
+  argsname = A(1).name;
+  for jj = 2:numel(A)
+    argsname = [argsname ',' A(jj).name];
+  end
+  
+  % Extract parameters
+  W = find(pl, 'weights');
+  V = find(pl, 'cov');
+  
+  % Build matrices for lscov
+  C = A(1:end-1);
+  Y = A(end);
+  
+  H = C(:).y;
+  y = Y.y;
+  if isa(W,'ao'), W = W.y; end;
+  if isa(V,'ao'), V = V.y; end;
+    
+  if isempty(V)
+    [p,stdx,mse,s] = lscov(H, y, W);
+  else
+    [p,stdx,mse,s] = lscov(H, y, W, V);
+  end
+      
+  % prepare model, units, names
+  model = [];
+  for jj = 1:numel(p)
+    names{jj} = ['C' num2str(jj)];
+    units{jj} = Y.yunits / C(jj).yunits;
+    xunits{jj} = C(jj).yunits;
+    xvar{jj} = ['X' num2str(jj)];
+    if jj == 1
+      model = ['C' num2str(jj) '*X' num2str(jj) ' '];
+    else
+      model = [model ' + C' num2str(jj) '*X' num2str(jj)];
+    end
+  end
+
+  model = smodel(plist('expression', model, ...
+    'params', names, ...
+    'values', p, ...
+    'xvar', xvar, ...
+    'xunits', xunits, ...
+    'yunits', Y.yunits ...
+    ));  
+  
+  % Build the output pest object
+  X = pest;
+  X.setY(p);
+  X.setDy(stdx);
+  X.setCov(s);
+  X.setChi2(mse);
+  X.setNames(names{:});
+  X.setYunits(units{:});
+  X.setModels(model);
+  X.name = sprintf('lscov(%s)', argsname);
+  X.addHistory(getInfo('None'), pl, ao_invars, [A(:).hist]);
+  % Set procinfo object
+  X.procinfo = plist('MSE', mse);
+  % Propagate 'plotinfo'
+  plotinfo = [A(:).plotinfo];
+  if ~isempty(plotinfo)
+    X.plotinfo = combine(plotinfo);
+  end
+  
+  % Set outputs
+  varargout{1} = X;
+    
+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.op, '$Id: lscov.m,v 1.42 2011/04/08 08:56:17 hewitson Exp $', sets, pl);
+  ii.setModifier(false);
+  ii.setArgsmin(2);
+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();
+  
+  % Weights
+  p = param({'weights','An ao containing weights for the fit.'}, paramValue.EMPTY_DOUBLE);
+  pl.append(p);
+  
+  % Cov
+  p = param({'cov','An ao containing a covariance matrix for the fit.'}, paramValue.EMPTY_DOUBLE);
+  pl.append(p);
+  
+end
+% END