Convert Set of (x,y) Coordinates Into Polygon

I converted the attached coordinates (x,y) into an alphashape. See image.

shp = alphaShape(coordinates(:,1),coordinates(:,2),'HoleThreshold',50);

What I need are the coordinates at the polygon vertices (shown in red), and importantly, in the proper order shown with the numbers. I would like to be able to use the polyshape function next...

I have toyed with boundaryfacet, delauney but with no luck.

Any suggestions?

Thanks.

Answer 1

[Matt J]

This uses tspsearch from the File Exchange,

https://www.mathworks.com/matlabcentral/fileexchange/71226-tspsearch

load coordinates

shp = alphaShape(coords(:,1),coords(:,2),'HoleThreshold',50);

[bf,P]=boundaryFacets(shp);

p=tspsearch(P,5); P=P(p,:);

pgon1=polyshape(P); pgon2=polyshape(P,Simplify=true);

Warning: Polyshape has duplicate vertices, intersections, or other inconsistencies that may produce inaccurate or unexpected results. Input data has been modified to create a well-defined polyshape.

Warning: Polyshape has duplicate vertices, intersections, or other inconsistencies that may produce inaccurate or unexpected results. Input data has been modified to create a well-defined polyshape.

idx= ismember(pgon1.Vertices, pgon2.Vertices,'rows');

V=flipud(pgon1.Vertices(idx,:));

[~,s]=min(sum(V,2)); V=circshift(V,1-s);

V(4,1)=V(3,1); V(5,1)=V(6,1); V([3,6],:)=[]; %Fix bad corners

plot(pgon1,FaceColor='none') ;hold on
scatlabel( scatter(V(:,1),V(:,2)) );hold off

●

[Paul Safier]

Hey @Star Strider and @Matt J , these are both really nice solutions! Thank you. Matt, nice use of the TSP! I don't know which one to accept as they're both very appreciated.

[Matt J]

Thanks @Paul Safier. You can Accept whichever solution you end up using, and you can upvote the other solutions.

Answer 2

[Star Strider]

Try this —

LD = load('coordinates.mat')

LD = struct with fields:
    coords: [17650x2 double]

coordinates = LD.coords;

shp = alphaShape(coordinates(:,1),coordinates(:,2),'HoleThreshold',50);

shpx = shp.Points(:,1);
shpy = shp.Points(:,2);

[minx,maxx] = bounds(shpx);
[miny,maxy] = bounds(shpy);

minxv = shpx == minx;
% nnz(minxv)

maxxv = shpx == maxx;
% nnz(maxxv)

minyv = shpy == miny;
% nnz(minyv)
minyidx = find(minyv);
endsidx = find(diff(minyidx) > 151);
corner = minyidx(endsidx:endsidx+1);
midpt = round(mean(corner));
[minmidpty,maxmidpty] = bounds(shpy(midpt+(0:151)));

VertexPoints = table([1; 8; 6; 7; 2; 5; 3; 4], [minx; minx; maxx; maxx; shpx(corner(1)); shpx(corner(2)); shpx(corner(1)); shpx(corner(2))], [miny; maxy; miny; maxy;  shpy(corner(1)); shpy(corner(2)); minmidpty; minmidpty], VariableNames=["Corner","X","Y"] );
VertexPoints = sortrows(VertexPoints,1)

VertexPoints = 8x3 table
    Corner      X         Y   
    ______    ______    ______

      1       124.41    25.084
      2       173.58    25.084
      3       173.58    126.42
      4       225.75    126.42
      5       225.75    25.084
      6       274.92    25.084
      7       274.92    175.59
      8       124.41    175.59

figure
plot(shp)
hold on
plot(minx, miny, 'rs', MarkerFaceColor='r')
plot(minx, maxy, 'rs', MarkerFaceColor='r')
plot(maxx, miny, 'rs', MarkerFaceColor='r')
plot(maxx, maxy, 'rs', MarkerFaceColor='r')
plot(shpx(corner), shpy(corner), 'rs', MarkerFaceColor='r')
plot(shpx(corner(1)), minmidpty, 'rs', MarkerFaceColor='r')
plot(shpx(corner(2)), minmidpty, 'rs', MarkerFaceColor='r')
hold off
axis([100  300    20  180])
text(VertexPoints{:,2}, VertexPoints{:,3}, compose('%2d',VertexPoints{:,1}), Color='r', FontWeight='bold', Vert='middle', Horiz='left', FontSize=12)

●

EDIT — Added text call.

.

Answer 3

[Walter Roberson]

load coordinates
coordinates = coords;
k = boundary(coordinates);
plot(coordinates(k,1), coordinates(k,2))

[Paul Safier]

@Walter Roberson curious that using boundary alone results in those two chopped corners. Any idea why considering the inital data points clearly had a corners there?

[Matt J]

It just depends what shrink factor you use. It is not the case that boundary() alone produces chopped corners. Your original alphashape has them as well, but they are smaller and less noticeable because the alpha radius you used happened to be a more optimal choice. We could tune boundary()'s shrink factor as well for better results:

load coordinates
coordinates = coords;
k = boundary(coordinates,0.99);
plot(coordinates(k,1), coordinates(k,2))

●

The main incompleteness in Walter's solution, however, is that it doesn't actually find the vertices. Just the boundary points:

numel(k)

ans = 801

[Paul Safier]

@Matt J agreed, I'm wondering if bypassing the use of alphashape in place of boundary might add some robustness for the case of other shapes. If boundary is more forgiving than alphashape with a less-than-optimal factor, then it might be preferable.

[Paul]

load coordinates
coordinates = coords;
shp = alphaShape(coordinates(:,1),coordinates(:,2),'HoleThreshold',50);
[bf,P] = boundaryFacets(shp);
pgon = polyshape(P);  % KeepCollinearPoints = false by default

Warning: Polyshape has duplicate vertices, intersections, or other inconsistencies that may produce inaccurate or unexpected results. Input data has been modified to create a well-defined polyshape.

figure
plot(pgon)
hold on
plot(pgon.Vertices(:,1),pgon.Vertices(:,2),'o')  % only 10 points

●

Not clear to me what the expectations are for handling the "double vertices" on those inner corners, which I think are also there using the boundary() solution shown above.

boundary doesn't have an option to remove colinear points AFAICT, so I suppose one could use boundary() with whatever shrink factor and then covert that result to a polyshape to get rid of colinear points.

Answer 4

[Catalytic]

load coordinates

shp = alphaShape(coords(:,1),coords(:,2),'HoleThreshold',50);

[~,P]=boundaryFacets(shp);

chull=convhull( polyshape(coords));

Warning: Polyshape has duplicate vertices, intersections, or other inconsistencies that may produce inaccurate or unexpected results. Input data has been modified to create a well-defined polyshape.

concavity =subtract(chull   , polyshape(P));

Warning: Polyshape has duplicate vertices, intersections, or other inconsistencies that may produce inaccurate or unexpected results. Input data has been modified to create a well-defined polyshape.

[xlim,ylim]=boundingbox(concavity);

innerV=table2array(combinations(xlim,ylim));

finalPolyshape=subtract( chull , polyshape(innerV([1,2,4,3],:))   );

V=flipud(finalPolyshape.Vertices);

plot(finalPolyshape); 
hold on; 
scatter(V(:,1),V(:,2),'ro','filled'); 
scatter(coords(1:2:end,1),coords(1:2:end,2),'.k','SizeData',1)
hold off

●