Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Reduce allocations / improve performance #43

Merged
merged 5 commits into from
Sep 23, 2021
Merged

Reduce allocations / improve performance #43

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
4c3c015
add polyfill benchmark; minor cleanup based upon ReSharper feedback a…
pocketken Sep 18, 2021
2412ee0
more tweaks based on profiling
pocketken Sep 19, 2021
026c1d3
this should make a bit of a dent
pocketken Sep 19, 2021
eb7fb38
get rid of compiler warning
pocketken Sep 19, 2021
b2304fc
whee bit more
pocketken Sep 23, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
30 changes: 15 additions & 15 deletions src/H3/Algorithms/Lines.cs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
using System.Collections.Generic;
using System.Linq;

using H3.Extensions;
using H3.Model;

Expand Down Expand Up @@ -55,34 +55,34 @@ public static IEnumerable<H3Index> LineTo(this H3Index origin, H3Index destinati
startIjk = LocalCoordIJK.ToLocalIJK(origin, origin);
endIjk = LocalCoordIJK.ToLocalIJK(origin, destination);
} catch {
return Enumerable.Empty<H3Index>();
yield break;
}

// get grid distance between start/end
int distance = startIjk.GetDistanceTo(endIjk);
var distance = startIjk.GetDistanceTo(endIjk);

// Convert IJK to cube coordinates suitable for linear interpolation
startIjk.Cube();
endIjk.Cube();

double d = distance;
double iStep = distance > 0 ? (endIjk.I - startIjk.I) / d : 0.0;
double jStep = distance > 0 ? (endIjk.J - startIjk.J) / d : 0.0;
double kStep = distance > 0 ? (endIjk.K - startIjk.K) / d : 0.0;
var iStep = distance > 0 ? (endIjk.I - startIjk.I) / d : 0.0;
var jStep = distance > 0 ? (endIjk.J - startIjk.J) / d : 0.0;
var kStep = distance > 0 ? (endIjk.K - startIjk.K) / d : 0.0;

double startI = startIjk.I;
double startJ = startIjk.J;
double startK = startIjk.K;

return Enumerable.Range(0, distance + 1)
.Select(n => {
var currentIjk = CoordIJK.CubeRound(
startI + iStep * n,
startJ + jStep * n,
startK + kStep * n
).Uncube();
return LocalCoordIJK.ToH3Index(origin, currentIjk);
});
for (var n = 0; n < distance + 1; n += 1) {
CoordIJK.CubeRound(
startI + iStep * n,
startJ + jStep * n,
startK + kStep * n,
endIjk
).Uncube();
yield return LocalCoordIJK.ToH3Index(origin, endIjk);
}
}

}
Expand Down
126 changes: 60 additions & 66 deletions src/H3/Algorithms/Polyfill.cs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,7 +1,8 @@
using System;
using System.Collections.Generic;
using System.Linq;
using H3.Extensions;
using H3.Model;
using static H3.Constants;
using NetTopologySuite.Algorithm.Locate;
using NetTopologySuite.Geometries;
using NetTopologySuite.LinearReferencing;
Expand All @@ -11,6 +12,7 @@
namespace H3.Algorithms {

internal sealed class PositiveLonFilter : ICoordinateSequenceFilter {

public bool Done => false;

public bool GeometryChanged => true;
Expand All @@ -19,9 +21,11 @@ public void Filter(CoordinateSequence seq, int i) {
double x = seq.GetX(i);
seq.SetOrdinate(i, Ordinate.X, x < 0 ? x + 360.0 : x);
}

}

internal sealed class NegativeLonFilter : ICoordinateSequenceFilter {

public bool Done => false;

public bool GeometryChanged => true;
Expand All @@ -30,42 +34,53 @@ public void Filter(CoordinateSequence seq, int i) {
double x = seq.GetX(i);
seq.SetOrdinate(i, Ordinate.X, x > 0 ? x - 360.0 : x);
}

}

/// <summary>
/// Polyfill algorithms for H3Index.
/// </summary>
public static class Polyfill {

private static readonly ICoordinateSequenceFilter _negativeLonFilter = new NegativeLonFilter();

private static readonly ICoordinateSequenceFilter _positiveLonFilter = new PositiveLonFilter();

/// <summary>
/// Returns all of the H3 indexes that are contained within the provided
/// Polygon at the specified resolution. Supports Polygons with holes.
/// </summary>
/// <param name="polygon">Containment polygon</param>
/// <param name="resolution">H3 resolution</param>
/// <returns>Indicies where center point is contained within polygon</returns>
public static IEnumerable<H3Index> Fill(this Polygon polygon, int resolution) {
public static IEnumerable<H3Index> Fill(this Geometry polygon, int resolution) {
bool isTransMeridian = polygon.IsTransMeridian();
var testPoly = isTransMeridian ? SplitPolygon(polygon) : polygon;
var testPoly = isTransMeridian ? SplitGeometry(polygon) : polygon;

HashSet<H3Index> searched = new();
HashSet<ulong> searched = new();

Stack<H3Index> toSearch = new(GetIndicies(testPoly.Coordinates, resolution));
Stack<H3Index> toSearch = new(TraceCoordinates(testPoly.Coordinates, resolution));
if (toSearch.Count == 0 && !testPoly.IsEmpty) {
toSearch.Push(H3Index.FromPoint(testPoly.InteriorPoint, resolution));
toSearch.Push(testPoly.InteriorPoint.Coordinate.ToH3Index(resolution));
}

IndexedPointInAreaLocator locator = new(testPoly);
var coordinate = new Coordinate();
var faceIjk = new FaceIJK();

while (toSearch.Count != 0) {
var index = toSearch.Pop();

if (index != H3Index.Invalid) {
foreach (var neighbour in GetKRingInPolygon(index, locator, searched)) {
if (neighbour == H3Index.Invalid) continue;
yield return neighbour;
toSearch.Push(neighbour);
}
foreach (var neighbour in index.GetNeighbours()) {
if (searched.Contains(neighbour)) continue;
searched.Add(neighbour);

var location = locator.Locate(neighbour.ToCoordinate(coordinate, faceIjk));
if (location != Location.Interior)
continue;

yield return neighbour;
toSearch.Push(neighbour);
}
}
}
Expand All @@ -78,49 +93,55 @@ public static IEnumerable<H3Index> Fill(this Polygon polygon, int resolution) {
/// <param name="resolution"></param>
/// <returns></returns>
public static IEnumerable<H3Index> Fill(this LineString polyline, int resolution) =>
polyline.Coordinates.GetIndicies(resolution);
polyline.Coordinates.TraceCoordinates(resolution);

/// <summary>
/// Gets all of the H3 indices that define the provided set of Coordinates.
/// Gets all of the H3 indices that define the provided set of <see cref="Coordinate"/>s.
/// </summary>
/// <param name="coordinates"></param>
/// <param name="resolution"></param>
/// <returns></returns>
public static IEnumerable<H3Index> GetIndicies(this Coordinate[] coordinates, int resolution) {
public static IEnumerable<H3Index> TraceCoordinates(this Coordinate[] coordinates, int resolution) {
HashSet<H3Index> indicies = new();

// trace the coordinates
int coordLen = coordinates.Length - 1;
for (int c = 0; c < coordLen; c += 1) {
var coordLen = coordinates.Length - 1;
FaceIJK faceIjk = new();
GeoCoord v1 = new();
GeoCoord v2 = new();
Vec3d v3d = new();
for (var c = 0; c < coordLen; c += 1) {
// from this coordinate to next/first
var vertA = coordinates[c];
var vertB = coordinates[c + 1];
var vA = coordinates[c];
var vB = coordinates[c + 1];
v1.Longitude = vA.X * M_PI_180;
v1.Latitude = vA.Y * M_PI_180;
v2.Longitude = vB.X * M_PI_180;
v2.Latitude = vB.Y * M_PI_180;

// estimate number of indicies between points, use that as a
// number of segments to chop the line into
var count = GeoCoord.FromCoordinate(vertA)
.LineHexEstimate(GeoCoord.FromCoordinate(vertB), resolution);
var count = v1.LineHexEstimate(v2, resolution);

for (int j = 1; j < count; j += 1) {
// interpolate line
var interpolated = LinearLocation.PointAlongSegmentByFraction(vertA, vertB, j / count);
var index = H3Index.FromCoordinate(interpolated, resolution);
if (!indicies.Contains(index)) indicies.Add(index);
var interpolated = LinearLocation.PointAlongSegmentByFraction(vA, vB, (double)j / count);
indicies.Add(interpolated.ToH3Index(resolution, faceIjk, v3d));
}
}

return indicies;
}

/// <summary>
/// Determines whether or not the polygon is flagged as transmeridian;
/// Determines whether or not the geometry is flagged as transmeridian;
/// that is, has an arc > 180 deg lon.
/// </summary>
/// <param name="polygon"></param>
/// <param name="geometry"></param>
/// <returns></returns>
public static bool IsTransMeridian(this Polygon polygon) {
if (polygon.IsEmpty) return false;
var coords = polygon.Envelope.Coordinates;
public static bool IsTransMeridian(this Geometry geometry) {
if (geometry.IsEmpty) return false;
var coords = geometry.Envelope.Coordinates;
return Math.Abs(coords[0].X - coords[2].X) > 180.0;
}

Expand All @@ -129,45 +150,18 @@ public static bool IsTransMeridian(this Polygon polygon) {
/// a multipolygon by clipping coordinates on either side of it and
/// then unioning them back together again.
/// </summary>
/// <param name="originalPolygon"></param>
/// <param name="originalGeometry"></param>
/// <returns></returns>
private static Geometry SplitPolygon(Polygon originalPolygon) {
var left = originalPolygon.Copy();
left.Apply(new NegativeLonFilter());
var right = originalPolygon.Copy();
right.Apply(new PositiveLonFilter());

var polygon = left.Union(right);
return polygon.IsEmpty ? originalPolygon : polygon;
private static Geometry SplitGeometry(Geometry originalGeometry) {
var left = originalGeometry.Copy();
left.Apply(_negativeLonFilter);
var right = originalGeometry.Copy();
right.Apply(_positiveLonFilter);

var geometry = left.Union(right);
return geometry.IsEmpty ? originalGeometry : geometry;
}

/// <summary>
/// Executes a k = 1 neighbour search for the provided H3 index, returning
/// any neighbours that have center points contained within the provided
/// polygon and that are not already present within the provided search
/// history hashset.
/// </summary>
/// <param name="index">H3 index to get neighbours for</param>
/// <param name="locator">IndexedPointInAreaLocator to use for point-in-poly
/// checks</param>
/// <param name="searched">Hashset of previously searched indicies; will
/// be updated to include any newly discovered neighbours automatically.
/// </param>
/// <returns>Neighbouring H3 indicies who's center points are contained
/// within the provided polygon</returns>
private static IEnumerable<H3Index> GetKRingInPolygon(H3Index index, IndexedPointInAreaLocator locator, HashSet<H3Index> searched) =>
index.GetKRing(1)
.Where(cell => {
if (searched.Contains(cell.Index)) {
return false;
}
searched.Add(cell.Index);
var coord = cell.Index.ToPoint().Coordinate;
var location = locator.Locate(coord);
return location == Location.Interior;
})
.Select(cell => cell.Index);

}

}
}
18 changes: 9 additions & 9 deletions src/H3/Algorithms/Rings.cs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -129,23 +129,23 @@ public static IEnumerable<RingCell> GetKRingSlow(this H3Index origin, int k) {

// since k >= 0, start with origin
Queue<RingCell> queue = new();
Dictionary<H3Index, int> searched = new();
Dictionary<ulong, int> searched = new();
queue.Enqueue(new RingCell { Index = origin, Distance = 0 });

while (queue.Count != 0) {
var cell = queue.Dequeue();
yield return cell;

var nextK = cell.Distance + 1;
if (nextK <= k) {
for (int i = 0; i < 6; i += 1) {
var neighbour = cell.Index.GetDirectNeighbour(LookupTables.CounterClockwiseDirections[i]).Item1;
if (neighbour == origin || neighbour == H3Index.Invalid || (searched.TryGetValue(neighbour, out int previousK) && previousK <= nextK)) {
continue;
}
searched[neighbour] = nextK;
queue.Enqueue(new RingCell { Index = neighbour, Distance = nextK });
if (nextK > k)
continue;

foreach (var neighbour in cell.Index.GetNeighbours()) {
if (neighbour == origin || searched.TryGetValue(neighbour, out int previousK) && previousK <= nextK) {
continue;
}
searched[neighbour] = nextK;
queue.Enqueue(new RingCell { Index = neighbour, Distance = nextK });
}
}
}
Expand Down