How to find a figure's centroid inside a polyshape?

I have the following problem:

I am plotting the boundaries of several municipalities and I need to list each one.

As there are many municipalities, I created a script that does this numbering automatically.

I used the centroid function to find a better position automatically, but some municipalities, like the one shown below, have a concave shape and the centroi is outside the limits of the municipalities.

My question is: there is a automatic way to find a position as central as possible within the city limits?

This is the result I get when I use the centroid of each polyshape to put the number. But for some municipalities with irregular shapes the numbers may be outside the limits (not well centralized):

Ps.: The boundaries of the municipalities are polyshapes and I have all the information regarding the location of the borders (I attached the polyshape for this example).

[Matt J]

Where would you consider the "center" of the above shape to be? Could you mark it on your posted figure?

[David Franco]

In fact Matt, anywhere within the limits is already viable for me.

Answer 1

[Matt J]

In fact Matt, anywhere within the limits is already viable for me.

If so, you could use nearestvertex command to project the centroid to the nearest point within city limits,

https://www.mathworks.com/help/matlab/ref/polyshape.nearestvertex.html

[David Franco]

Using this option I get a point (in black below) that is in the border.

I need the number to be all inside the limits.

[Matt J]

Shrink the boundaries of the polyshape using polybuffer (see example below) by whatever buffer distance from the city limits that you need. Then project to the nearest vertex like before.

>> pgonBig=polyshape([0 0; 0 1; 1 0]);
>> pgonSmall=polybuffer(pgonBig, -0.1);
>> plot([pgonBig,pgonSmall])

[David Franco]

Very clever option. It's working.

Any ideas on how to set the reduction parameter as a percentage of the polyshape?

[Matt J]

Perhaps as follows,

>> city=polyshape([0 0; 0 1; 1 0]);
>> A=area(city);
>> p = .8; %percentage

>> [d,res]=fzero(@(d)area(polybuffer(city,d))-p*A, 0)

        d =
        
           -0.0309
        
        
        res =
        
           8.0499e-09

>> innercity=polybuffer(city,d);

[David Franco]

Perfect! Thank you!

Answer 2

[Chad Greene]

Also answered here, but reposting my answer in case anyone stumbles upon this thread.

I just ran into this problem when trying to place a text label in the middle of a crescent-shaped ice shelf. The centroid or the mean or median of the coordinates of the ice shelf polygon are all outside the bounds of the ice shelf. Here's the best solution I could come up with:

% Convert the outline to a polyshape:
P = polyshape(x,y); 
% And get the delaunay triangulation of the polygon:  
T = triangulation(P);
% Now find the center points of all the triangles:
[C,r] = circumcenter(T);
% Get the index of the centerpoint that has the largest radius from the boundary:
[~,ind] = max(r); 
% These center coordinates are in the center of the fattest part of the polygon:  
xc = C(ind,1); 
yc = C(ind,2); 

[David Franco]

Thank you for your answer Chad! But I got two points completely outside of the map...

Here is the code:

% Find all centroids

load mapapr.mat

centroids = zeros(399,2);

warning off

tic

for ii = 1:399

A = Snew(ii).X;

B = Snew(ii).Y;

% Convert the outline to a polyshape:

P = polyshape(A,B);

% And get the delaunay triangulation of the polygon:

T = triangulation(P);

% Now find the center points of all the triangles:

[C,r] = circumcenter(T);

% Get the index of the centerpoint that has the largest radius from the boundary:

[~,ind] = max(r);

% These center coordinates are in the center of the fattest part of the polygon:

centroids(ii,1) = C(ind,1);

centroids(ii,2) = C(ind,2);

end

toc

warning on

mapshow(Snew,'FaceColor','c');

hold on

plot(centroids(:,1),centroids(:,2),'k*')

hold off

Attached is my data.

Thank you!

[Chad Greene]

Oh, that's crazy. I'm able to replicate the problem (thanks for providing the data, minimal working example, and the figure!!) but I can't say I understand why this is an issue for those polygons. Nonetheless, here's a fix that works. Just below the circumcenter line, add this:

% If it's not inside the polygon, make sure it doesn't have the maximum radius: 
r(~isinterior(P,C(:,1),C(:,2))) = -1;

[Chad Greene]

And here's what it looks like with text labels:

text(centroids(:,1),centroids(:,2),{Snew.MUNICIPIO},...
    'horiz','center','vert','middle','fontsize',7,'clipping','on')

[Chad Greene]

Alrighty, @David Franco, you inspired me to add a new polycenter function to the Climate Data Toolbox for Matlab.

load mapapr.mat

[xc,yc] = polycenter(Snew); 

mapshow(Snew,'FaceColor','c');
hold on
text(xc,yc,num2str((1:399)'),'fontsize',7,'horiz','center','vert','middle')

[David Franco]

Nice @Chad Greene!!! Thank you!!!