function [sorted,indices] = sort(x,dim,mode)

% Copyright 2004 The MathWorks, Inc.

eml_assert(nargin  > 0, 'error', 'Not enough input arguments.');

% Validate data
if isempty(x)
    sorted = x;
    if nargout == 2,
        indices = zeros(size(x));
    end
    return
end

% Validate DIM input
if nargin == 1
    % Sort the row if it's a row vector, otherwise sort each column
    if size(x,1)==1
        dim = 2;
    else
        dim = 1;
    end
else 
    eml_assert(~isempty(dim) && isscalar(dim) && (dim==1 || dim==2),...
      'DIM must be 1 or 2.'); 
end
if dim == 1,
    dim0 = 2;
else
    dim0 = 1;
end

% Validate the sort MODE input
if nargin > 2
    % Flip the indexing of the sorted results when using a descending sort
    if strcmp(mode,'ascend')
        direction = 1;
    elseif strcmp(mode,'descend')
        direction = -1;
    else
        eml_assert(false, 'error', 'sorting direction must be ''ascend'' or ''descend''.');
    end
else
    direction = 1;
end

if isreal(x),
    sorted = zeros(size(x),class(x));
else
    sorted = zeros(size(x),class(x)) + 0i;
end
if nargout > 1,
    indices = zeros(size(x));
end
n = size(x,dim0);
for i=1:n,
    if dim == 2,
        v = x(i,:);
    else
        v = x(:,i);
    end
    if nargout > 1,
        [sorted_v,sorted_indices] = sort0(v,direction);
    else
        sorted_v = sort0(v,direction);
    end
    if dim == 2,
        sorted(i,:) = sorted_v;
        if nargout > 1,
            indices(i,:) = sorted_indices;
        end
    else
        sorted(:,i) = sorted_v;
        if nargout > 1,
            indices(:,i) = sorted_indices;
        end
    end
end
        
function [sorted_v,sorted_indices] = sort0(v,direction)
    n = length(v);
    if isreal(v),
        v0 = zeros(size(v),class(v));
    else
        v0 = zeros(size(v),class(v)) + 0i;
    end
    if nargout > 1,
        indices0 = zeros(size(v));
        indices = indices0;
        for i=1:n
            indices(i) = i;
        end
    end
    i = 1;
    while i <= n,
        i2 = 2*i;
        j = 1;
        % i % disp
        while j <= n,
            p = j;
            p_end = j + i;
            if p_end > n + 1,
                p_end = n + 1;
            end
            q = p_end;
            q_end = j + i2;
            if q_end > n + 1,
                q_end = n + 1;
            end
            k = 1;
            k_end = q_end - j;
            % j % disp
            while k <= k_end,
                b = q == q_end;
                if ~b && p < p_end,
                    if isreal(v(p))
                        if direction == 1,
                            b = v(p) <= v(q) || isnan(v(q));
                        else
                            b = v(q) < v(p) || isnan(v(p));
                        end
                    else
                        d = abs(v(p)) - abs(v(q));
                        if d == 0,
                            d = angle(v(p)) - angle(v(q));
                        end
                        if direction == 1,
                            b = d <= 0;
                        else
                            b = d > 0;
                        end
                    end
                end
                        
                if b,
                    v0(k) = v(p);
                    if nargout > 1,
                        indices0(k) = indices(p);
                    end
                    p = p + 1;
                else
                    v0(k) = v(q);
                    if nargout > 1,
                        indices0(k) = indices(q);
                    end
                    q = q + 1;
                end
                k = k + 1;
            end;
            % v0 % disp
            k = 1;
            p = j;
            while k <= k_end,
                v(p) = v0(k);
                if nargout > 1,
                    indices(p) = indices0(k);
                end
                k = k + 1;
                p = p + 1;
            end
            % v % disp
            j = q_end;
        end
        i = i2;
    end
    sorted_v = v;
    if nargout > 1,
        sorted_indices = indices;
    end