How to find adjacent polygons of a specific polygon in a VoronoiMesh?

Question

I have developed some code for VoronoiMesh objects, and I am able to find the polygon with the highest number of adjacent polygons. How do I then find the indices of the adjacent polygons?

Please post your code. – David G. Stork Jan 29 '16 at 21:47 — David G. Stork, Jan 29 '16 at 21:47

PlatoManiac · Accepted Answer · 2016-01-31T11:00:07.130

Finding the points which correspond to the Voronoi cells with maximum number of neighboring cells is easy by sorting the VertexDegree on underlying DelanunayMesh.

Lets find the indexes of the pts with maximum vertex degree which assures that they will be surrounded by most number of Voronoi cells.

pts = RandomReal[{-1, 1}, {250, 2}];
Vmesh = VoronoiMesh[pts];
graph = With[{m = DelaunayMesh[pts]}, 
   Graph[MeshCoordinates[m], 
    MeshCells[m, 1] /. Line[{start_, end_}] -> {start, end}, 
    VertexCoordinates -> MeshCoordinates[m], 
    GraphLayout -> "PlanarEmbedding", 
    EdgeStyle -> Directive[Opacity@.5, Dashed, White], 
    VertexStyle -> White, VertexSize -> 3]];
highestdeg = Flatten@Position[VertexDegree[graph], VertexDegree[graph] // Max]

{76, 91, 247}

So we have our VoronoiMesh and the point indices highestdeg whose neighboring cells we will search with the following function. This function will return the cell index which contains the point with given index and the cell indices of it's neighboring mesh cells.

FindAdjPoly[ptIndex_, {pts_, Vmesh_}] := 
 Block[{pt, mpt, regs, lines, poly, loc, all},
  mpt = MeshCoordinates[Vmesh];
  (* Get the mesh polygons *)
  regs = MeshCells[Vmesh, 2];
  pt = pts[[ptIndex]];
  (* Select the polygon containing pt *)
  poly = SelectFirst[(regs /. Polygon[a__] :> {a, Polygon@mpt[[a]]}), 
    RegionMember[Last@#, pt] &];
  (* Find the index of the polygon in mesh *)
  loc = MeshCellIndex[Vmesh, {Polygon[First@poly]}][[1, 2]];
  (*Get the edges of the polygon *)
  lines = Partition[First@poly, 2, 1, 1];
  (* Find polygons that share the above edges *)
  all = DeleteDuplicates@
    Flatten[With[{ed = #}, 
        Position[
         regs /. Polygon[a__] :> (MemberQ[a, #] & /@ ed), {True,True}]] & /@ lines];
  {loc, Complement[all, {loc}]}
  ]

Calling the function:

polys = FindAdjPoly[#, {pts, Vmesh}] & /@ highestdeg

{{248, {11, 28, 35, 99, 149, 198, 202, 210, 234, 244}}, {249, {8, 14, 16, 59, 62, 106, 119, 169, 186, 188}}, {250, {75, 77, 80, 84, 85, 138, 147, 172, 176, 192}}}

Why not check if everything is working as expected...

Code to visualize:

Show[HighlightMesh[Vmesh, 
  Table[Prepend[#, 2] & /@ Transpose@{polys[[i, 2]]}, {i,Length@polys}]], graph, 
 Graphics[{White, Disk[#, 0.023] & /@ pts[[highestdeg]], Red,
 PointSize@.008, Point /@ pts[[highestdeg]]}]]

Hi, Thanks a lot for the code. It's a very nice code. But, when I try to visualize them, I found an error. May be because of version problem. However, I have a question. Why you use Table and transpose those values. And, It would be great for me if you describe a bit of the "FindAdjPoly" in the code. I want to learn how you do it? Thanks — Odrisso, Jan 31 '16 at 00:20
@Odrisso I have added comments to the code for FindAdjPoly. Hope this helps you. About the visualization I am working on version 10.3.1 and I see no errors. The Table and Transpose is used to process the output polys so that HighlightMesh can make sense of it. — PlatoManiac, Jan 31 '16 at 11:03

score 6 · Answer 2 · answered Jan 29 '16 at 22:57

Here are the points at the centers of the Voronoi cells:

myCenters = RandomReal[1, {100, 2}];

This creates a DelanunayMesh of the center points, i.e., the graph linking centers that is mathematically dual to the Voronoi cell representation:

myDelaunayMesh =  DelaunayMesh[myCenters];

This extracts the links that include the most surrounded point:

myLinkstomostSurroundedPoint = (Last@
    Sort@Table[
      Select[myLinks, (#[[1, 1]] == myCenters[[i]] || #[[1, 2]] == 
           myCenters[[i]])  &], {i, Length[myCenters]}]) /. 
  Line[{a__, b__}] -> {a, b};

This finds the center that appears the most often in that list, i.e., is the most-surrounded cell center:

mostSurroundedPoint = 
 Commonest@Flatten[myLinkstomostSurroundedPoint, 1][[1]];

Here are the other points in the list, i.e., the "surrounding" points:

theSurroundingPoints = 
 Complement[(Union@
    Flatten[myLinkstomostSurroundedPoint, 1]), {mostSurroundedPoint}];

Here is the position in the list of the most surrounded point:

Position[myCenters, mostSurroundedPoint[[1]]]

(*

*)

Here are the indexes of the other, surrounding points:

Flatten@(Position[myCenters, #] & /@ theSurroundingPoints)

(*

{78, 80, 43, 8, 86, 10, 84, 30, 46}

*)

Hi. Thanks. But,I still found problem. Here is my code: v = VoronoiMesh[mydata]; p1 = MeshPrimitives[v, 2] dv = Length[p1] n = Length /@ Level[pv, {2}] sv = Max[n] nv = Count[n, sv] i1 = Level[Position[n, sv], {2}] In this way, I determine the region with the highest number of adjacent regions. But, I can't determine the adjacent regions now. Please let me know about it. I run your code. But, can't understand it. — Odrisso, Jan 29 '16 at 23:50
The adjacent regions have as centers given by the last line in my code. — David G. Stork, Jan 30 '16 at 00:03
When I use the last line, it choose all the polygons. Not the surrounding ones. I think I have to Set up a For loop to find the Union of neighbouring adjacent polygons. Can you please kindly tell me about it. — Odrisso, Jan 30 '16 at 00:37
The last line gives the indexes in myCenters of the neighbors. So the two dimensional points are myCenters[[78]], myCenters[[80]], etc. Clearly not "all" the polygons. — David G. Stork, Jan 30 '16 at 00:47

score 2 · Answer 3 · edited Feb 23 '18 at 18:16

It has never been easier as with the new IGMeshCellAdjacencyGraph tool in IGraph/M, (version 0.3.93 or later), which is based on this answer. The following method works for an arbitrary MeshRegion.

Here is a mesh:

SeedRandom[1234];
pts = RandomReal[{-1, 1}, {250, 2}];
Vmesh = VoronoiMesh[pts];

This is how you find a polygon with the maxmimal number of neighbors:

Needs["IGraphM`"]
g = IGMeshCellAdjacencyGraph[Vmesh, 2];
HighlightMesh[
 Vmesh,
 Extract[
  VertexList[g],
  Position[VertexDegree[g], Max[VertexDegree[g]]]
  ]
 ]

score 2 · Answer 4 · answered Feb 01 '18 at 06:09

You can use properties of a MeshRegion object to find adjacent faces/meshes. The only problem is that MeshRegion properties return face indices, and unfortunately, face labels don't match up with face indices. So, here are two functions to convert between face labels and face indices:

fromLabel[vm_, labels_] := Ordering[Length /@ vm["Faces"]][[labels]]
toLabel[vm_, indices_] := ReplaceAll[
    indices,
    Thread[Ordering[Length /@ vm["Faces"]]->Range[Length[vm["Faces"]]]]
]

Then, the following will highlight adjacent faces:

highlightAdjacentFaces[vm_, label_] := HighlightMesh[
    vm,
    Thread[{2, toLabel[vm, vm["FaceFaceConnectivity"][[fromLabel[vm, label]]]]}]
]

Example:

SeedRandom[1];
pts = RandomReal[10, {250, 2}];

vm = VoronoiMesh[pts];
highlightAdjacentFaces[vm, 200]
highlightAdjacentFaces[vm, {100, 150, 200}]

How to find adjacent polygons of a specific polygon in a VoronoiMesh?

4 Answers4

Linked