MATLAB Delaunay algorithm to extract the boundary of discrete points

Recently, I encountered the problem of extracting the boundary of discrete points in the project. My friend 1, who is not particularly skilled at matlab, must be confused at the beginning. Which algorithm can be used to effectively extract the contour lines of all known points? Through a large number of literature search and code experiments, I found several good contour extraction code, here to make a summary, and I hope to be able to encounter the same problem to inspire friends.

Three methods for boundary extraction of discrete points:

1.Convhull discrete point set to obtain the boundary

2.Alpha Shape algorithm detects edge points

3.Delaunay 3 Angle division algorithm

The first two methods have been summarized in a previous blog and I won't expand on them here, but I'll focus on the third algorithm.

The general idea of this algorithm is as follows:

1. Use delaunay function to perform Delaunay 3 Angle subdivision on all data points. The return value of delaunay function is a matrix of N * 3, in which N is the number of 3 Angle forms subdivided, and 3 is the ordinal number of 3 end points of each 3 Angle form.

2. Extract all the edges connected by delaunay 3 Angle slices according to triangles matrix, scan each row of triangles matrix in turn, add the edges connected by delaunay 3 Angle slices to a new matrix, and finally form a matrix of M * 2, in which M is the number of edges connected by 1 in total.

3. Obviously, the edges on the minimum convex polygon should appear only once in the above matrix. Therefore, all the edges that appear more than once in the above matrix should be removed.

4. According to the edges of the minimum convex polygon, it is easy to get the order of the nodes forming the minimum convex polygon, so as to solve the problem.

The input parameter points is a 2 * P matrix, and the number of P data points. The first line is the corresponding x coordinates of these data points, and the second line is the corresponding y coordinates. The output parameter polygon is a 2 * Q matrix. Q is the number of vertices of a convex polygon (end to end). The first row is the corresponding x coordinates of these vertices, and the second row is the corresponding y coordinates. The code implementation is as follows:

function polygon = minimal_convex_polygon(points)
 %  for  delaunay 3 Corner splitting, keeping all connected edges in the matrix  lines  In the 
 triangles = sort(delaunay(points(1, :), points(2, :)), 2);
 lines = zeros(size(triangles, 1) * 3, 2);
 for i = 1:size(triangles, 1)
 lines(3 * i - 2,:) = [triangles(i, 1), triangles(i, 2)];
 lines(3 * i - 1,:) = [triangles(i, 1), triangles(i, 3)];
 lines(3 * i,:) = [triangles(i, 2), triangles(i, 3)];
 end
 %  To get rid of  lines  More than 1 The edge of time 
 [~, IA] = unique(lines, 'rows');
 lines = setdiff(lines(IA, :), lines(setdiff(1:size(lines, 1), IA), :), 'rows');
 %  tracking  lines  The vertex number of the convex polygon is saved in  seqs  In the 
 seqs = zeros(size(lines, 1) + 1,1);
 seqs(1:2) = lines(1, :);
 lines(1, :) = [];
 for i = 3:size(seqs)
 pos = find(lines == seqs(i - 1));
 row = rem(pos - 1, size(lines, 1)) + 1;
 col = ceil(pos / size(lines, 1));
 seqs(i) = lines(row, 3 - col);
 lines(row, :) = [];
 end
 %  According to the  seqs  .   Get the vertex coordinates of the convex polygon 
 polygon = points(:, seqs);
end

The implementation function is defined and the following call is made:

plot(Pp(1,:),Pp(2,:), '*r', 'LineWidth', 4);  % Pp The first 1 behavior x Coordinates of the first 2 behavior y coordinates 
polygon = minimal_convex_polygon(Pp);
hold on;
plot(polygon(1, :), polygon(2, :), 'LineWidth', 2);

I'm not going to add the image yet, so anyone interested can try it out.