Last active
February 27, 2024 16:18
-
-
Save dfaivre/acfef42cdbf411555956e9eba65dd30d to your computer and use it in GitHub Desktop.
PolyLabel (pole of inaccessibility of a polygon) - C#/NTS
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
using System; | |
using System.Collections.Generic; | |
using System.Linq; | |
using GeoAPI.Geometries; | |
using NetTopologySuite.Geometries; | |
namespace Foo | |
{ | |
/// <summary> | |
/// Finds the largest area inside a polygon (ie, for a label) | |
/// | |
/// See: https://github.com/mapbox/polylabel/blob/master/polylabel.js | |
/// </summary> | |
public static class PolyLabel | |
{ | |
public static SearchCell FindMax( | |
IPolygon polygon, double precision = 1) | |
{ | |
var boundingBox = polygon.EnvelopeInternal; | |
var cellSize = System.Math.Min(boundingBox.Width, boundingBox.Height); | |
var searchSize = cellSize / 2; | |
// I'm not sure why the source has a search size of zero for the centroid cell. Oh well? | |
// DMF 2018-05-16 | |
var bestCell = new SearchCell(polygon.Centroid, searchSize: 0, polygon); | |
var cellCandidates = new SortedList<SearchCell, SearchCell>(); | |
for (var x = boundingBox.MinX; x < boundingBox.MaxX; x += cellSize) | |
{ | |
for (var y = boundingBox.MinY; y < boundingBox.MaxY; y += cellSize) | |
{ | |
var cell = new SearchCell(new Point(x + searchSize, y + searchSize), searchSize, polygon); | |
cellCandidates.Add(cell, cell); | |
} | |
} | |
while (cellCandidates.Any()) | |
{ | |
// pick the most promising cell from the queue | |
var cell = cellCandidates.Keys[0]; | |
cellCandidates.RemoveAt(0); | |
// update the best cell if we found a better one | |
if (cell.CentroidDistance > bestCell.CentroidDistance) | |
{ | |
bestCell = cell; | |
} | |
// do not drill down further if there's no chance of a better solution | |
if (cell.MaxPotentialDistance - bestCell.CentroidDistance <= precision) continue; | |
// else: split the cell and serach | |
var splitCells = cell.Split(); | |
foreach (var splitCell in splitCelrrrls) | |
{ | |
cellCandidates.Add(splitCell, splitCell); | |
} | |
} | |
return bestCell; | |
} | |
public class SearchCell : IComparable<SearchCell> | |
{ | |
private readonly IPolygon _searchPolygon; | |
private static readonly double Sqrt2 = System.Math.Sqrt(2); | |
public SearchCell(IPoint point, double searchSize, double centroidDistance, double maxPotentialDistance) | |
{ | |
Point = point; | |
SearchSize = searchSize; | |
CentroidDistance = centroidDistance; | |
MaxPotentialDistance = maxPotentialDistance; | |
} | |
public SearchCell(IPoint point, double searchSize, IPolygon searchPolygon) | |
{ | |
_searchPolygon = searchPolygon; | |
Point = point; | |
SearchSize = searchSize; | |
CentroidDistance = PointToPolygonDistanceSigned(point, searchPolygon); | |
MaxPotentialDistance = CentroidDistance + searchSize * Sqrt2; | |
} | |
public IPoint Point { get; } | |
public double SearchSize { get; } | |
public double CentroidDistance { get; } | |
public double MaxPotentialDistance { get; } | |
public SearchCell[] Split() | |
{ | |
var splitSearchSize = SearchSize / 2; | |
SearchCell CreateSplitCell(double x, double y) | |
{ | |
return new SearchCell( | |
new Point(x, y), splitSearchSize, _searchPolygon); | |
} | |
return new[] | |
{ | |
CreateSplitCell(Point.X - splitSearchSize, Point.Y - splitSearchSize), | |
CreateSplitCell(Point.X + splitSearchSize, Point.Y - splitSearchSize), | |
CreateSplitCell(Point.X + splitSearchSize, Point.Y + splitSearchSize), | |
CreateSplitCell(Point.X - splitSearchSize, Point.Y + splitSearchSize), | |
}; | |
} | |
public int CompareTo(SearchCell other) => MaxPotentialDistance > other.MaxPotentialDistance ? -1 : 1; | |
} | |
private static double PointToPolygonDistanceSigned( | |
IPoint point, IPolygon polygon) | |
{ | |
var inside = polygon.Intersects(point); | |
var minDistance = | |
polygon.InteriorRings | |
.Select(r => r.Distance(point)) | |
.Concat(new[] {polygon.ExteriorRing.Distance(point)}) | |
.Min(); | |
return (inside ? 1 : -1) * minDistance; | |
} | |
} | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Hi @dfaivre
Thanks for sharing this. Any idea why this is not working for the following polygon?
I am showing the marker at the center point of the Search Cell. I would like the marker to be in the visual center of the polygon.