5
$\begingroup$

Say I have a mesh region which encloses a region. How do I get a mesh region where this region is filled? Take for example the meshregion below:

reg=MeshRegion[List[List[1.`,0.`],List[2.`,0.`],List[3.`,0.`],List[4.`,0.`],List[5.`,0.`],List[4.5`,2.5`],List[0.5`,2.5`],List[2.5`,0.5`],List[2.5`,4.5`],List[5.`,1.`],List[5.`,2.`],List[5.`,3.`],List[5.`,4.`],List[5.`,5.`],List[4.`,5.`],List[3.`,5.`],List[2.`,5.`],List[1.`,5.`],List[0.`,5.`],List[0.`,4.`],List[0.`,3.`],List[0.`,2.`],List[0.`,1.`],List[0.`,0.`],List[5.`,6.`],List[5.`,7.`],List[5.`,8.`],List[4.9`,7.`],List[4.`,8.`],List[0.5`,5.`]],List[Polygon[List[List[23,24,1],List[7,22,23],List[1,2,8],List[3,8,2],List[1,8,23],List[7,23,8],List[21,7,20],List[7,21,22],List[30,19,20],List[20,7,9],List[20,18,30],List[18,20,9],List[17,18,9],List[9,16,17],List[8,3,4],List[6,4,10],List[10,4,5],List[15,9,6],List[10,11,6],List[12,13,6],List[11,12,6],List[6,13,15],List[16,9,15],List[29,28,27],List[15,25,28],List[14,25,15],List[26,27,28],List[25,26,28],List[15,13,14],List[6,8,4]]]]]

Note that the mesh region is not concave and I want to preserve that so taking the convex hull is not solving the issue.

Related I would like to know how to get the filled mesh region of everything to a certain side of a mesh (until some cut-off).


The requested outputs would in this case be something matching: (but of course automatized)

Region@RegionUnion[reg, Rectangle[{0, 0}, {5, 5}]]

and the completion to the right (with cut-off 10)

Region@RegionUnion[reg, Rectangle[{0, 0}, {5, 5}],  Rectangle[{5, 0}, {10, 8}]]

completion to above

Region@RegionUnion[reg, Rectangle[{0, 0}, {5, 5}], 
  Rectangle[{0, 0}, {10, 10}]]

I would like to also fill holes in cases that the boundary to the outside is a point. See for example the following region:

MeshRegion[List[List[1.`,0.`],List[2.`,0.`],List[3.`,0.`],List[4.`,0.`],List[5.`,0.`],List[4.5`,2.5`],List[0.5`,2.5`],List[2.5`,0.5`],List[2.5`,4.5`],List[5.`,1.`],List[5.`,2.`],List[5.`,3.`],List[5.`,4.`],List[5.`,5.`],List[4.`,5.`],List[3.`,5.`],List[2.`,5.`],List[1.`,5.`],List[0.`,5.`],List[0.`,4.`],List[0.`,3.`],List[0.`,2.`],List[0.`,1.`],List[0.`,0.`],List[5.`,6.`],List[5.`,7.`],List[5.`,8.`],List[4.9`,7.`],List[4.`,8.`],List[0.5`,5.`]],List[Polygon[List[List[23,24,1],List[7,22,23],List[1,2,8],List[3,8,2],List[21,7,20],List[7,21,22],List[30,19,20],List[20,18,30],List[17,18,9],List[9,16,17],List[8,3,4],List[10,4,5],List[10,11,6],List[12,13,6],List[11,12,6],List[16,9,15],List[29,28,27],List[14,25,15],List[26,27,28],List[25,26,28],List[15,13,14]]]]];
$\endgroup$
6
$\begingroup$

One idea is to convert the MeshRegion to a BoundaryMeshRegion, and then to extract the bounding polygon. Your MeshRegion:

reg = MeshRegion[
    {
    {1.,0.},{2.,0.},{3.,0.},{4.,0.},{5.,0.},{4.5,2.5},{0.5,2.5},{2.5,0.5},
    {2.5,4.5},{5.,1.},{5.,2.},{5.,3.},{5.,4.},{5.,5.},{4.,5.},{3.,5.},
    {2.,5.},{1.,5.},{0.,5.},{0.,4.},{0.,3.},{0.,2.},{0.,1.},{0.,0.},
    {5.,6.},{5.,7.},{5.,8.},{4.9,7.},{4.,8.},{0.5,5.}
    },
    {Polygon[{
        {23,24,1},{7,22,23},{1,2,8},{3,8,2},{1,8,23},{7,23,8},{21,7,20},
        {7,21,22},{30,19,20},{20,7,9},{20,18,30},{18,20,9},{17,18,9},
        {9,16,17},{8,3,4},{6,4,10},{10,4,5},{15,9,6},{10,11,6},{12,13,6},
        {11,12,6},{6,13,15},{16,9,15},{29,28,27},{15,25,28},{14,25,15},
        {26,27,28},{25,26,28},{15,13,14},{6,8,4}
    }]
    }
];

The equivalent BoundaryMeshRegion:

boundary = BoundaryMesh[reg]

enter image description here

Extract the bounding polygon:

p = First @ boundary["BoundaryPolygons"];

Visualization:

Region @ p

enter image description here

| improve this answer | |
$\endgroup$
  • $\begingroup$ +1. It seems to work for the more difficult case as well when points are on the boundary too: RegionUnion @@ (RegionUnion @@ MeshPrimitives[reg, 2])["BoundaryPolygons"]. It's a shame I can't find this documented anywhere though. $\endgroup$ – flinty Jul 30 at 18:11
  • $\begingroup$ I can't say for sure, but I think we can use the "BoundaryNesting" property to extract the outer polygons programatically. e.g. BoundaryDiscretizeRegion[ RegionUnion @@ Extract[boundary["BoundaryPolygons"], Position[boundary["BoundaryNesting"], {_, _, 0, _}]]] $\endgroup$ – Chip Hurst Jul 31 at 0:09
4
$\begingroup$

First, let's write the data like this:

coords = {{1.`, 0.`}, {2.`, 0.`}, {3.`, 0.`}, {4.`, 0.`}, {5.`, 
    0.`}, {4.5`, 2.5`}, {0.5`, 2.5`}, {2.5`, 0.5`}, {2.5`, 
    4.5`}, {5.`, 1.`}, {5.`, 2.`}, {5.`, 3.`}, {5.`, 4.`}, {5.`, 
    5.`}, {4.`, 5.`}, {3.`, 5.`}, {2.`, 5.`}, {1.`, 5.`}, {0.`, 
    5.`}, {0.`, 4.`}, {0.`, 3.`}, {0.`, 2.`}, {0.`, 1.`}, {0.`, 
    0.`}, {5.`, 6.`}, {5.`, 7.`}, {5.`, 8.`}, {4.9`, 7.`}, {4.`, 
    8.`}, {0.5`, 5.`}};
poly = Polygon[{{23, 24, 1}, {7, 22, 23}, {1, 2, 8}, {3, 8, 2}, {1, 8,
      23}, {7, 23, 8}, {21, 7, 20}, {7, 21, 22}, {30, 19, 20}, {20, 7,
      9}, {20, 18, 30}, {18, 20, 9}, {17, 18, 9}, {9, 16, 17}, {8, 3, 
     4}, {6, 4, 10}, {10, 4, 5}, {15, 9, 6}, {10, 11, 6}, {12, 13, 
     6}, {11, 12, 6}, {6, 13, 15}, {16, 9, 15}, {29, 28, 27}, {15, 25,
      28}, {14, 25, 15}, {26, 27, 28}, {25, 26, 28}, {15, 13, 14}, {6,
      8, 4}}];

We need to get the coordinates into 3D to use RepairMesh so we can fill the hole:

reg = MeshRegion[Append[#, 0] & /@ coords, poly];

Then we fill the hole. Notice how the mesh has a quite poor triangulation, even though it filled the hole properly:

filled = RepairMesh[reg, "HoleEdges"]

3d mesh filled bad triangulation

So we'll rectify that by getting it back into 2D and re-discretizing it. We'll get the polygons and drop the z coordinates, then convert to a Graphics and finally call DiscretizeGraphics:

gr = Graphics[Polygon[#[[1, All, 1 ;; 2]]] & /@ MeshPrimitives[filled, 2]];
(* re-descretize to get new clean mesh *)
DiscretizeGraphics@gr

better mesh

This cleans up some of the triangulation issues, but there's still a problem. The faces that filled the hole are actually on top of the 2D mesh and not well connected. Any region operations on this mesh could produce spurious lines and connectivity issues. Unfortunately, Mathematica doesn't provide a way to set a tolerance in RegionUnion, otherwise I would have just union'd all the polygons to begin with.

To fix this I can rasterize the graphics first at a very high resolution and then use ImageMesh:

gr = Rasterize[
   Graphics[{White, 
     Polygon[#[[1, All, 1 ;; 2]]] & /@ MeshPrimitives[filled, 2]}, 
    Background -> Black], ImageSize -> {2048, 2048}];
(* re-descretize to get cleaner mesh *)
GraphicsRow[{ImageMesh[gr], TriangulateMesh@ImageMesh[gr]}]

best mesh

Note 1: Using the rasterize approach will cause the scale to change uniformly. If you need the original scale, you will need to use FindGeometricTransform on some select boundary points to find the scale/translate matrix that returns the region to the original size.

Note 2: The geometry produced by RepairMesh is not very good and introduces some extra polygons that shouldn't be there as shown below.

repair mesh y u no wrk

You could also accomplish the filling using the raster method using this much simpler one-liner:

reg = MeshRegion[coords, poly];
ImageMesh@
 FillingTransform[
  Graphics[{White, reg, ImageSize -> {2048, 2048}}, 
   Background -> Black]]

It's also possible to create a Graph of the RegionBoundary and find connected component subgraphs, then use FindShortestTour to get their polygons. This approach seems better to me because there's no scaling problems and it also gives you both the hole and the filled outer polygon:

reg = MeshRegion[coords, poly];
gr = Graph[
   MeshPrimitives[RegionBoundary[reg], 1] /. 
    Line[x_] :> UndirectedEdge @@ x];

With[{cgc = ConnectedGraphComponents[gr]},
 Graphics[{Thick,
   Riffle[
    RandomColor[
     Length[cgc]], (EdgeList[#] /. 
        UndirectedEdge[x_, y_] :> Line[{x, y}]) & /@ cgc]}]
 ]
polys = Polygon[Last[FindShortestTour[#]]] & /@ 
  ConnectedGraphComponents[gr]
MeshRegion[polys[[1]]]

good mesh graph

polygons from graph

| improve this answer | |
$\endgroup$
  • $\begingroup$ Thanks! I am going through it now. How did you extract the poly = Polygon[...] from the MeshRegion? (For the coordinates I understand to use MeshCoordinates). (Okay I managed: Polygon@Flatten[((MeshCells@reg)[[3]] /. Polygon -> List), {1, 2}]) $\endgroup$ – Kvothe Jul 30 at 16:59
  • $\begingroup$ I just took your data and rewrote it manually. From a mesh region you can use MeshPrimitives[reg,2] to get the polygons, if your version of Mathematica supports it. Regions in Mathematica are really broken and there are all sorts of weird hacks to get around inadequacies in things like RegionIntersection and RegionUnion. If I come up with another way to solve it I'll post another answer or extend this one. $\endgroup$ – flinty Jul 30 at 17:03
  • $\begingroup$ Great! Thank you! Both work great. I hope they scale well with bigger initial meshes. (I will try.) One problem I encounter is that if a hole is separated only by a point somewhere the hole is not filled. Is there anything simple that can be done about that? (I added an example region of that form to the question.) $\endgroup$ – Kvothe Jul 30 at 17:44
  • $\begingroup$ @Kvothe the RepairMesh method is not that great - thanks Wolfram! - I've just noticed that in the top right where it's narrow RepairMesh introduces an extra triangle and makes it a bit thicker. I think the graph method is the the best. As for the single point on the boundary - the rasterize approach is probably the easiest, but the graph approach could be extended to use find graph partitions. This is getting very difficult though. $\endgroup$ – flinty Jul 30 at 17:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.