--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/m-toolbox/classes/@ao/bin_data.m	Wed Nov 23 19:22:13 2011 +0100
@@ -0,0 +1,298 @@
+% BIN_DATA rebins aos data, on logarithmic scale, linear scale, or arbitrarly chosen.
+% The rebinning is done taking the mean of the bins included in the range
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% DESCRIPTION: BIN_DATA rebins aos data, on logarithmic scale, linear scale, or arbitrarly chosen.
+% The rebinning is done taking the mean of the bins included in the range
+%
+% CALL:        bs = bin_data(a1,a2,a3,...,pl)
+%              bs = bin_data(as,pl)
+%              bs = as.bin_data(pl)
+%
+% INPUTS:      aN   - input analysis objects
+%              as   - input analysis objects array
+%              pl   - input parameter list
+%
+% OUTPUTS:     bs   - array of analysis objects, one for each input
+%
+% <a href="matlab:utils.helper.displayMethodInfo('ao', 'bin_data')">Parameters Description</a>
+%
+% The code is inherited from D Nicolodi, UniTN
+%
+% VERSION:    $Id: bin_data.m,v 1.20 2011/05/10 16:46:48 mauro Exp $
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+
+function varargout = bin_data(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
+  [as, ao_invars] = utils.helper.collect_objects(varargin(:), 'ao', in_names);
+  
+  % Decide on a deep copy or a modify
+  bs = copy(as, nargout);
+  
+  % Apply defaults to plist
+  usepl = applyDefaults(getDefaultPlist(), varargin{:});
+
+  x_scale     = find(usepl, 'x_scale', find(usepl, 'xscale'));
+  x_vals      = find(usepl, 'x_vals', find(usepl, 'xvals'));
+  resolution  = find(usepl, 'resolution');
+  range       = find(usepl, 'range');
+  method      = lower(find(usepl, 'method'));
+  inherit_dy  = utils.prog.yes2true(find(usepl, 'inherit-dy', find(usepl, 'inherit_dy')));
+  
+  % Loop over input AOs
+  for jj = 1:numel(bs)
+      
+    % check input data
+    if isa(bs(jj).data, 'data2D')
+      
+      w = find(usepl, 'weights');
+      
+      if isa(w, 'ao')
+        w = w.y;
+      end
+      
+      if isempty(w)
+        w =  1./(bs(jj).dy).^2;
+      end
+      
+      if isempty(x_vals)
+        if isempty(x_scale) || isempty(resolution)
+          error('### Please specify a scale and density for binning, OR the list of the values to bin around');
+        else
+          
+          switch lower(x_scale)
+            case {'lin', 'linear'}
+              % Case of linear binning
+              % number of bins in the rebinned data set
+              N = resolution; 
+
+              % maximum and minimum x
+              if ~isempty(range) && isfinite(range(1))
+                xmin = range(1);
+              else
+                xmin = min(bs(jj).x);
+              end
+              
+              if ~isempty(range) && isfinite(range(2))
+                xmax = range(2);
+              else
+                xmax = max(bs(jj).x);
+              end
+
+              dx = (xmax - xmin)/N;
+              
+              x_min = bs(jj).x(1) + dx*(0:(N-1))';
+              x_max = bs(jj).x(1) + dx*(1:N)';
+              
+            case {'log', 'logarithmic'}
+              % Case of log-based binning
+
+              % maximum and minimum x
+              if ~isempty(range) && isfinite(range(1))
+                xmin = range(1);
+              else
+                xmin = min(bs(jj).x(bs(jj).x > 0));
+              end
+              
+              if ~isempty(range) && isfinite(range(2))
+                xmax = range(2);
+              else
+                xmax = max(bs(jj).x);
+              end
+              
+              alph = 10^(1/resolution);
+                        
+              % number of bins in the rebinned data set
+              N = ceil(log10(xmax/xmin) * resolution);
+              
+              % maximum and minimum x-value for each bin
+              x_min = xmin*alph.^(0:(N-1))';
+              x_max = xmin*alph.^(1:N)';
+            otherwise
+              error(['### Unknown scaling option ' x_scale '. Please choose between ''lin'' and ''log']);
+          end
+        end
+      else
+        % number of bins in the rebinned data set
+        % If the x-scale is an AO, then take the x values
+        if isa(x_vals, 'ao')
+          if eq(x_vals.xunits, bs(jj).xunits)
+            x_vals = x_vals.x;
+          else
+            error('x_vals AO and data AO have different x-units');
+          end
+        elseif ~isnumeric(x_vals)
+          error('Unsupported x_vals object');
+        end
+        N = length(x_vals) - 1;
+        x_min = x_vals(1:N);
+        x_max = x_vals(2:N+1);
+      end
+      
+      x =  bs(jj).x;
+      y =  bs(jj).y;
+      dy = bs(jj).dy;
+      
+      % preallocate output vectors
+      xr = zeros(N, 1);
+      yr = zeros(N, size(y, 2));
+      if strcmpi(method, 'mean') || strcmpi(method, 'wmean')
+        dyr = zeros(N, size(y, 2));
+      else
+        dyr = [];
+      end
+      nr = zeros(N, 1);
+      
+      % compute the averages
+      for kk = 1:N
+        in = x >= x_min(kk) & x < x_max(kk);
+        if any(in)
+          nr(kk) = sum(in);               % number of points averaged in this bin          
+          
+          switch method
+            case {'mean', 'median', 'max', 'min', 'rms'}
+              xr(kk) = feval(method, x(in));           % rebinned x bins;
+              yr(kk) = feval(method, y(in));           % rebinned y bins;
+              if strcmpi(method, 'mean')
+                dyr(kk) = std(y(in), 0)/sqrt(nr(kk));
+                % check for zeros in the uncertainty and replace it with the individual point uncertainty
+                if dyr(kk) == 0
+                  if inherit_dy && ~isempty(dy)
+                    dyr(kk) = mean(dy(in));
+                  else
+                    dyr(kk) = Inf;
+                  end
+                end              
+              end              
+            case {'wmean'}
+              xr(kk)  = mean(x(in));                      % rebinned x bins;                                          
+              yr(kk)  = sum(y(in).*w(in))./sum(w(in));    % rebinned y bins;
+              dyr(kk) = 1./sqrt(sum(w(in)));              % rebinned dy bins;
+            otherwise
+              error(['### Unsupported method ' method]);
+          end
+        end
+      end
+      
+      % remove bins where we do not have nothing to average
+      in = nr ~= 0;
+      nr = nr(in);
+      xr = xr(in);
+      yr = yr(in,:);
+      if strcmpi(method, 'mean') || strcmpi(method, 'wmean')  
+        dyr = dyr(in,:);
+      end
+      
+      % set the new object data
+      bs(jj).setXY(xr, yr);
+      bs(jj).setDy(dyr);
+      
+      % nr goes into the procinfo
+      bs(jj).procinfo = plist('navs', nr);
+      
+      % set name
+      bs(jj).name = sprintf('bin_data(%s)', ao_invars{jj});
+      % Add history
+      bs(jj).addHistory(getInfo('None'), usepl, ao_invars(jj), bs(jj).hist);
+    else
+      warning('### Ignoring input AO number %d (%s); it is not a 2D data object.', jj, bs(jj).name)
+    end
+  end % loop over analysis objects
+  
+  % Set output
+  varargout = utils.helper.setoutputs(nargout, bs);
+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: bin_data.m,v 1.20 2011/05/10 16:46:48 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();
+  
+  % method
+  p = param({'method',['method for binning. Choose from:<ul>', ...
+    '<li>mean</li>', ...
+    '<li>median</li>', ...
+    '<li>max</li>', ...
+    '<li>min</li>', ...
+    '<li>rms</li>', ...
+    '<li>weighted mean (weights can be input or are taken from data dy)</li></ul>']}, ...
+    {1, {'MEAN', 'MEDIAN', 'MAX', 'MIN', 'RMS', 'WMEAN'}, paramValue.SINGLE});
+  pl.append(p);
+  
+  % x-scale
+  p = param({'xscale',['scaling of binning. Choose from:<ul>', ...
+    '<li>log - logaritmic</li>', ...
+    '<li>lin - linear</li></ul>']}, {1, {'LOG', 'LIN'}, paramValue.SINGLE});
+  pl.append(p);
+  
+  % resolution
+  p = param({'resolution',['When setting logaritmic x scale, it sets the number of points per decade.<br>' ...
+    'When setting linear x scale, it sets the number of points.']}, paramValue.DOUBLE_VALUE(10));
+  pl.append(p);
+  
+  % x_vals
+  p = param({'xvals',['List of x values to evaluate the binning between.<br>', ...
+    'It may be a vector or an ao, in which case it will take the x field']}, paramValue.DOUBLE_VALUE([]));
+  pl.append(p);
+  
+  % weights
+  p = param({'weights', ['List of weights for the case of weighted mean.<br>', ...
+    'If empty, weights will be taken from object(s) dy field as w = 1/dy^2']}, paramValue.DOUBLE_VALUE([]));
+  pl.append(p);
+
+  % range
+  p = param({'range', ['Range of x where to operate.<br>', ...
+    'If empty, the whole data set will be used']}, paramValue.DOUBLE_VALUE([]));
+  pl.append(p);
+  
+  % inherit_dy
+  p = param({'inherit_dy', ['Choose what to do in the case of mean, and bins with only one point. Choose from:<ul>', ...
+    '<li>''yes''  - take the uncertainty from the original data, if defined</li>', ...
+    '<li>''no''   - set it to Inf so it weighs 0 in averaged means</li></ul>' ...
+    ]}, paramValue.YES_NO);
+  pl.append(p);
+end
+% END