How is the polyarea function input order to receive correct results?

Hi everyone,

In the help of the polyarea function it says that the value returned can depend on the order of the input values?

Does anybody know more about that? When does it depend any when not? how do i hve to order my input values so that the result will be correct?

v1=[0 0 1 1];
v2=[0 1 0 1];
close all
figure
patch(v1,v2,rand(1,3))
A = polyarea(v1,v2) ;

This will lead to 2 triangles, somehow one with a negativ and one with positiv area. At least that is whats happening i think, because the result is 0, even though that these 4 point form a square.

v1=[0 0 1 1];
v2=[0 1 1 0];
close all
figure
patch(v1,v2,rand(1,3))
A = polyarea(v1,v2) ;

This instead will lead to a square with the area 1.

Im trying to use voronoin function, to get n articicial grin structure and i later have to dd more polyeders inside the existing "mesh" and the areas of those hould be calculated correctly, so i need to know how to order the input values.

Many thanks in advance

Best regards

Answer 1

[John D'Errico]

Yes. You do need to think about what you are doing.

I'll create a quasi-triangular region as an example.

px = [0, 1, .5, .49];
py = [0, 0, 1, 1-0.1];

plot(px([1 2 3 4 1]),py([1 2 3 4 1]),'r-')

So in the sequence I listed them, then connecting back to the beginning, the area is:

polyarea(px([1 2 3 4 1]),py([1 2 3 4 1]))
ans =
      0.48

This is exactly as I would expect, just a bit less than 1/2, which would be the area of the triangle that forms the convex hull. However, in a different sequence, the included area is seen as:

polyarea(px([1 4 2 3 1]),py([1 4 2 3 1]))
ans =
                      0.05
                      
plot(px([1 4 2 3 1]),py([1 4 2 3 1]),'g-')

So the same points, but in a different sequence, we got a completely different result, and that is expected.

Can you know which sequence was intended? Of course not. That is impossible, nor should a tool like polyarea be expected to know which you intended.

In the case of a convex polygon, IF the region is KNOWN to be convex, then you can just use a sort on polar angle.

t = rand(1,50)*2*pi;
px = cos(t);
py = sin(t);

polyarea(px,py)
ans =
         0.968571643648727

As it turns out, I would have expected an area that was reasonably close to pi.

mux = mean(px);
muy = mean(py);
[th,r] = cart2pol(px,py); % I could have used atan2 here also
[~,ind] = sort(th);

polyarea(px(ind),py(ind))
ans =
          3.11074834389343

So the area is now seen as quite close to pi. Which is correct? As I said, IF WE KNEW the region was a convex polygon, then the latter is easy. But given only the points from some general possibly non-convex domain in unsorted order, the problem is essentially impossible to solve without knowing the desired sequence of the points.

plot(px,py,'r-',px([ind,ind(1)]),py([ind,ind(1)]),'b-')

Is the area that of the blue curve in the last plot, or the area "inside" the red curve? See that when you have a Catholic polygon, that is a polygon that crosses itself multiple time,s then some parts of the polygon have positive area, some have a negative area. This effectively has some of those areas cancel each other out. (Think right hand rule, for those who know what a cross product is and what it tells you.)

[Marc Laub]

Ok, thanks.

so since voronoi leads to only convex polygons the polar angle sort is propaply done here, since in the ouput where the vertices per polygon are listed they are always in the correct order.

[John D'Errico]

My guess without thinking, is I would use a cross product to determine the orientation of a triangle. However, since voronoi uses delaunay to trianglute the set, I'm pretty sure that delaunay has the property that all triangles are oriented in a consistent direction.

px = [0 1 1 0];
py = [0 0 1 1];
tri = delaunay(px,py)
tri =
     3     4     2
     4     1     2

That is both triangles were created to have a cross product that points along the positive z-axis. We can think of it as if they were sorted in a counter-clockwise sequence, so increasing in polar angle.