angleSort3d

ANGLESORT3D Sort 3D coplanar points according to their angles in plane.

   PTS2 = angleSort3d(PTS);
   Considers all points are located on the same plane, and sort them
   according to the angle on plane. PTS is a [Nx2] array. Note that the
   result depends on the plane orientation: points can be in reverse order
   compared to expected. The reference plane is computed based on the
   first three points.

   PTS2 = angleSort3d(PTS, PTS0);
   Computes angles between each point of PTS and PT0. By default, uses
   centroid of points.

   PTS2 = angleSort3d(PTS, PTS0, PTS1);
   Specifies the point which will be used as a start.

   [PTS2, I] = angleSort3d(...);
   Also return in I the indices of PTS, such that PTS2 = PTS(I, :);

   See also
   points3d, angles3d, angleSort

0001 function varargout = angleSort3d(pts, varargin)
0002 %ANGLESORT3D Sort 3D coplanar points according to their angles in plane.
0003 %
0004 %   PTS2 = angleSort3d(PTS);
0005 %   Considers all points are located on the same plane, and sort them
0006 %   according to the angle on plane. PTS is a [Nx2] array. Note that the
0007 %   result depends on the plane orientation: points can be in reverse order
0008 %   compared to expected. The reference plane is computed based on the
0009 %   first three points.
0010 %
0011 %   PTS2 = angleSort3d(PTS, PTS0);
0012 %   Computes angles between each point of PTS and PT0. By default, uses
0013 %   centroid of points.
0014 %
0015 %   PTS2 = angleSort3d(PTS, PTS0, PTS1);
0016 %   Specifies the point which will be used as a start.
0017 %
0018 %   [PTS2, I] = angleSort3d(...);
0019 %   Also return in I the indices of PTS, such that PTS2 = PTS(I, :);
0020 %
0021 %   See also
0022 %   points3d, angles3d, angleSort
0023 
0024 % ------
0025 % Author: David Legland
0026 % E-mail: david.legland@inrae.fr
0027 % Created: 2005-11-24
0028 % Copyright 2005-2024 INRA - CEPIA Nantes - MIAJ (Jouy-en-Josas)
0029 
0030 % default values
0031 pt0     = mean(pts, 1);
0032 pt1     = pts(1,:);
0033 
0034 if isscalar(varargin)
0035     pt0 = varargin{1};
0036 elseif length(varargin)==2
0037     pt0 = varargin{1};
0038     pt1 = varargin{2};
0039 end
0040 
0041 % create support plane
0042 plane   = createPlane(pts(1:3, :));
0043 
0044 % project points onto the plane
0045 pts2d   = planePosition(pts, plane);
0046 pt0     = planePosition(pt0, plane);
0047 pt1     = planePosition(pt1, plane);
0048 
0049 % compute origin angle
0050 theta0  = atan2(pt1(2)-pt0(2), pt1(1)-pt0(1));
0051 theta0  = mod(theta0 + 2*pi, 2*pi);
0052 
0053 % translate to reference point
0054 n       = size(pts, 1);
0055 pts2d   = pts2d - repmat(pt0, [n 1]);
0056 
0057 % compute angles
0058 angle   = atan2(pts2d(:,2), pts2d(:,1));
0059 angle   = mod(angle - theta0 + 4*pi, 2*pi);
0060 
0061 % sort points according to angles
0062 [angle, I] = sort(angle); %#ok<ASGLU>
0063 
0064 
0065 % format output
0066 if nargout<2
0067     varargout{1} = pts(I, :);
0068 elseif nargout==2
0069     varargout{1} = pts(I, :);
0070     varargout{2} = I;
0071 end
0072