in src/h3lib/lib/algos.c [780:935]
int _polyfillInternal(const GeoPolygon* geoPolygon, int res, H3Index* out) {
// One of the goals of the polyfill algorithm is that two adjacent polygons
// with zero overlap have zero overlapping hexagons. That the hexagons are
// uniquely assigned. There are a few approaches to take here, such as
// deciding based on which polygon has the greatest overlapping area of the
// hexagon, or the most number of contained points on the hexagon (using the
// center point as a tiebreaker).
//
// But if the polygons are convex, both of these more complex algorithms can
// be reduced down to checking whether or not the center of the hexagon is
// contained in the polygon, and so this is the approach that this polyfill
// algorithm will follow, as it's simpler, faster, and the error for concave
// polygons is still minimal (only affecting concave shapes on the order of
// magnitude of the hexagon size or smaller, not impacting larger concave
// shapes)
//
// This first part is identical to the maxPolyfillSize above.
// Get the bounding boxes for the polygon and any holes
BBox* bboxes = H3_MEMORY(malloc)((geoPolygon->numHoles + 1) * sizeof(BBox));
assert(bboxes != NULL);
bboxesFromGeoPolygon(geoPolygon, bboxes);
// Get the estimated number of hexagons and allocate some temporary memory
// for the hexagons
int numHexagons = H3_EXPORT(maxPolyfillSize)(geoPolygon, res);
H3Index* search = H3_MEMORY(calloc)(numHexagons, sizeof(H3Index));
H3Index* found = H3_MEMORY(calloc)(numHexagons, sizeof(H3Index));
// Some metadata for tracking the state of the search and found memory
// blocks
int numSearchHexes = 0;
int numFoundHexes = 0;
// 1. Trace the hexagons along the polygon defining the outer geofence and
// add them to the search hash. The hexagon containing the geofence point
// may or may not be contained by the geofence (as the hexagon's center
// point may be outside of the boundary.)
const Geofence geofence = geoPolygon->geofence;
int failure = _getEdgeHexagons(&geofence, numHexagons, res, &numSearchHexes,
search, found);
// If this branch is reached, we have exceeded the maximum number of
// hexagons possible and need to clean up the allocated memory.
// LCOV_EXCL_START
if (failure) {
H3_MEMORY(free)(search);
H3_MEMORY(free)(found);
H3_MEMORY(free)(bboxes);
return failure;
}
// LCOV_EXCL_STOP
// 2. Iterate over all holes, trace the polygons defining the holes with
// hexagons and add to only the search hash. We're going to temporarily use
// the `found` hash to use for dedupe purposes and then re-zero it once
// we're done here, otherwise we'd have to scan the whole set on each insert
// to make sure there's no duplicates, which is very inefficient.
for (int i = 0; i < geoPolygon->numHoles; i++) {
Geofence* hole = &(geoPolygon->holes[i]);
failure = _getEdgeHexagons(hole, numHexagons, res, &numSearchHexes,
search, found);
// If this branch is reached, we have exceeded the maximum number of
// hexagons possible and need to clean up the allocated memory.
// LCOV_EXCL_START
if (failure) {
H3_MEMORY(free)(search);
H3_MEMORY(free)(found);
H3_MEMORY(free)(bboxes);
return failure;
}
// LCOV_EXCL_STOP
}
// 3. Re-zero the found hash so it can be used in the main loop below
for (int i = 0; i < numHexagons; i++) found[i] = 0;
// 4. Begin main loop. While the search hash is not empty do the following
while (numSearchHexes > 0) {
// Iterate through all hexagons in the current search hash, then loop
// through all neighbors and test Point-in-Poly, if point-in-poly
// succeeds, add to out and found hashes if not already there.
int currentSearchNum = 0;
int i = 0;
while (currentSearchNum < numSearchHexes) {
H3Index ring[MAX_ONE_RING_SIZE] = {0};
H3Index searchHex = search[i];
H3_EXPORT(kRing)(searchHex, 1, ring);
for (int j = 0; j < MAX_ONE_RING_SIZE; j++) {
if (ring[j] == H3_NULL) {
continue; // Skip if this was a pentagon and only had 5
// neighbors
}
H3Index hex = ring[j];
// A simple hash to store the hexagon, or move to another place
// if needed. This MUST be done before the point-in-poly check
// since that's far more expensive
int loc = (int)(hex % numHexagons);
int loopCount = 0;
while (out[loc] != 0) {
// If this branch is reached, we have exceeded the maximum
// number of hexagons possible and need to clean up the
// allocated memory.
// LCOV_EXCL_START
if (loopCount > numHexagons) {
H3_MEMORY(free)(search);
H3_MEMORY(free)(found);
H3_MEMORY(free)(bboxes);
return -1;
}
// LCOV_EXCL_STOP
if (out[loc] == hex) break; // Skip duplicates found
loc = (loc + 1) % numHexagons;
loopCount++;
}
if (out[loc] == hex) {
continue; // Skip this hex, already exists in the out hash
}
// Check if the hexagon is in the polygon or not
GeoCoord hexCenter;
H3_EXPORT(h3ToGeo)(hex, &hexCenter);
// If not, skip
if (!pointInsidePolygon(geoPolygon, bboxes, &hexCenter)) {
continue;
}
// Otherwise set it in the output array
out[loc] = hex;
// Set the hexagon in the found hash
found[numFoundHexes] = hex;
numFoundHexes++;
}
currentSearchNum++;
i++;
}
// Swap the search and found pointers, copy the found hex count to the
// search hex count, and zero everything related to the found memory.
H3Index* temp = search;
search = found;
found = temp;
for (int j = 0; j < numSearchHexes; j++) found[j] = 0;
numSearchHexes = numFoundHexes;
numFoundHexes = 0;
// Repeat until no new hexagons are found
}
// The out memory structure should be complete, end it here
H3_MEMORY(free)(bboxes);
H3_MEMORY(free)(search);
H3_MEMORY(free)(found);
return 0;
}