function [unique_positions, value_sums] = unique_totals(positions, values)
%UNIQUE_TOTALS return sums of values found at each unique position.
%
% [unique_positions, value_sums] = unique_totals(positions, values)
%
%
% Inputs:
%            positions NxD - may not provide DxN unless D==1
%               values Nx1 - can provide 1xN
%
% Outputs:
%     unique_positions MxD - M depends on input
%           value_sums Mx1 

% Iain Murray, November 2006.
% Updated to use accumarray in February 2008 after reading:
% http://blogs.mathworks.com/steve/2008/01/31/counting-occurrences-of-image-colors/

% Sanitize input
if isempty(positions) && isempty(values)
    unique_positions = [];
    value_sums = [];
    return
end
if min(size(values)) ~= 1
    error('Values must be a vector.');
end
values = values(:);
[N, D] = size(positions);
if (N == 1) && (D == length(values))
    positions = positions';
    [N, D] = size(positions);
end
if length(values) ~= N
    error('imurray:unique_totals:invalid_input', ...
        ['You must have the same number of positions and values.\n',...
        'Is the NxD parity correct on the positions matrix?']);
end

if D > 1
    [unique_positions, dummy, position_indexes] = unique(positions, 'rows');
    value_sums = accumarray(position_indexes, values);
else
    % A more complex, but faster implementation for D == 1. Matlab's unique
    % function seems to be a slower than it should be, and it's easier to avoid
    % when D == 1. Also, when I last tried, Octave didn't work with D == 1 and
    % the implementation for D > 1.

    % Sorting by position makes the values that need to be summed adjacent
    [sorted_positions, idx] = sort(positions);
    sorted_values = values(idx);
    % Find the blocks and do the sums
    position_block_ends = [find(diff(sorted_positions)); N];
    if D == 1
        unique_positions = sorted_positions(position_block_ends);
    end
    cum_values = cumsum(sorted_values);
    value_sums = diff([0; cum_values(position_block_ends)]);
end