function [t,p] = csort(T,P)
% CSORT  [t,p] = csort(T,P) Sorts T, consolidates P
% Version of 4/6/97
% Modified to work with Versions 4.2 and 5.1, 5.2
% T and P matrices with the same number of elements
% The vector T(:)' is sorted:
%   * Values in T differing by no more than 10^{-13}
%     are consolidated;
%   * Corresponding values in P are added.
T = T(:)';
n = length(T);
[TS,I] = sort(T);
d  = find([1,TS(2:n) - TS(1:n-1) >1e-13]);
t  = TS(d);
m  = length(t) + 1;
P  = P(I);
F  = [0 cumsum(P(:)')];
Fd = F([d length(F)]);
p  = Fd(2:m) - Fd(1:m-1);