6

The following code gets all vertices of all polygons (mesh cells) of VoronoiMesh[pts]:

SeedRandom[3]; 
pts = RandomReal[{-1, 1}, {25, 2}];
mesh = VoronoiMesh[pts];
vertices = MeshCoordinates[mesh];
Show[mesh, Graphics[{Black, Point[pts], Red, Point[vertices]}]]

This outputs:

Voronoi

My question

How can I get a list of vertices for each polygon and compute the area of each polygon using the Shoelace formula?

The output should be similar to:

Output

So, by clicking on the polygon number, it should show its vertices and its size.

I found this tool-tip image in Finding the perimeter, area and number of sides of a Voronoi cell

J. M.'s missing motivation
  • 124,525
  • 11
  • 401
  • 574
Eman
  • 199
  • 7

3 Answers3

4

Use MeshPrimitives like this:

Show[Graphics[{FaceForm@RGBColor[
    0.666, 0.776, 0.952], 
   Table[Tooltip[p, 
     Grid@{{"Perimeter", Perimeter@p}, {"Area", Area@p}, {"Edges", 
        Length @@ p}}], {p, MeshPrimitives[mesh, 2]}]}], 
 Graphics[{Black, Point[pts], Red, Point[vertices]}]]

enter image description here

M.R.
  • 31,425
  • 8
  • 90
  • 281
  • Thanks so much for your help. But, if I want the vertices of each polygon to be shown also with area,edges and Perimeter. How to do that?? – Eman Sep 08 '18 at 18:41
  • 1
    How are you ordering them? – M.R. Sep 08 '18 at 18:47
  • Thanks so much for your help and your reply. What did you mean by them ? Did you mean the vertices?? If you mean the vertices, I don't order them. the code get all vertices of all voronoi polygons. I want to get the vertices of each polygon, separately. So, by clicking on each polygon; I can get its vertices. – Eman Sep 08 '18 at 19:12
  • 2
    Note that PropertyValue[{mesh, 2}, MeshCellMeasure] is a faster way to get all of the areas. However I don't think the other properties can be computed in this way. – Greg Hurst Sep 08 '18 at 19:52
3

Here's an efficient way to implement the shoelace formula, assuming no self intersections:

ShoelaceArea[Polygon[pts_?MatrixQ]] := 
  0.5 * #1.(RotateLeft[#2] - RotateRight[#2])& @@ Transpose[pts]

A comparison:

shoeareas = ShoelaceArea /@ MeshPrimitives[mesh, 2]; // AbsoluteTiming

 {0.000233, Null}

areas = PropertyValue[{mesh, 2}, MeshCellMeasure]; // AbsoluteTiming

 {0.000013, Null}

Max[Abs[shoeareas - areas]]

3.33067*10^-16
Greg Hurst
  • 35,921
  • 1
  • 90
  • 136
2
polygons = Join @@ MeshCells[mesh, 2, "Multicells" -> True][[All, 1]];
polygondata = With[{x = MeshCoordinates[mesh]}, Map[
    p \[Function] Partition[x[[p]], 2, 1, 1],
    polygons
    ]];
areas = 0.5 Total[Map[Det, polygondata, {2}], {2}];
circumferences = Total[Map[Norm, Differences /@ polygondata, {2}], {2}];

For the tooltipping, you can also use the option MeshCellLabel of MeshRegion, but that's are a bit unwieldy:

MeshRegion[mesh, MeshCellLabel -> Map[
   i \[Function] ({2, i} -> Tooltip[
       i,
       Grid[{
         {"Vertices", polygons[[i]]},
         {"Vertex Coordinates", polygondata[[i, All, 1]]},
         {"Area", areas[[i]]},
         {"Perimeter", circumferences[[i]]}
         },
        Alignment -> {Left, Top}
        ]
       ]
     ),
   Range[MeshCellCount[mesh, 2]]
   ]
 ]
Henrik Schumacher
  • 106,770
  • 7
  • 179
  • 309