From 9af7218f239813851ae37d1770c8d200a96346dd Mon Sep 17 00:00:00 2001 From: Nick Rabinowitz Date: Mon, 27 Nov 2023 11:28:47 -0800 Subject: [PATCH] Add support for full containment and overlapping modes in polygonToCellsExperimental (#796) Adds support for full containment mode and overlapping modes in polygonToCellsExperimental. * Defines an enum for the containment modes. My assumption here is that the lowest 4 bits of the 32-bit flag are reserved for containment modes (not sure that we'd need more than 7, but wanted to be on the safe side), and that we'd use them as an integer, not as bitwise flags, since the modes are mutually exclusive. * Implements the check for containment in iterStepPolygonCompact. This was very straightforward as we already have cellBoundaryInsidePolygon polygon available. * Implements the check for overlapping cells in iterStepPolygonCompact. This checks center point inclusion first, then shares code with the full containment check, substituting a new cellBoundaryCrossesPolygon check for cellBoundaryInsidePolygon. --- scripts/make_countries.js | 17 +- .../benchmarkPolygonToCellsExperimental.c | 79 +++++- src/apps/testapps/testH3Memory.c | 13 +- src/apps/testapps/testPolyfillInternal.c | 77 +++++- src/apps/testapps/testPolygonInternal.c | 24 +- src/apps/testapps/testPolygonToCells.c | 10 +- .../testapps/testPolygonToCellsExperimental.c | 231 ++++++++++++++---- .../testPolygonToCellsReportedExperimental.c | 44 ++-- src/h3lib/include/polygon.h | 24 +- src/h3lib/lib/algos.c | 10 +- src/h3lib/lib/polyfill.c | 84 ++++++- src/h3lib/lib/polygon.c | 47 +++- 12 files changed, 537 insertions(+), 123 deletions(-) diff --git a/scripts/make_countries.js b/scripts/make_countries.js index 4e8a495d4..8e810dff3 100644 --- a/scripts/make_countries.js +++ b/scripts/make_countries.js @@ -131,6 +131,7 @@ async function makeCountries(sourceUrl, targetPath) { #include "h3api.h" #include "mathExtensions.h" #include "polyfill.h" +#include "polygon.h" const GeoPolygon COUNTRIES[${polygons.length}] = {${ polygons.map((poly, i) => formatGeoPolygon(poly, names[i])).join(',') @@ -149,13 +150,25 @@ for (int res = 0; res < MAX_RES + 1; res++) { BENCHMARK(polygonToCells_AllCountries1, 5, { for (int index = 0; index < ${polygons.length}; index++) { - H3_EXPORT(polygonToCells)(&COUNTRIES[index], res, 0, hexagons); + H3_EXPORT(polygonToCells)(&COUNTRIES[index], res, CONTAINMENT_CENTER, hexagons); } }); BENCHMARK(polygonToCells_AllCountries2, 5, { for (int index = 0; index < ${polygons.length}; index++) { - H3_EXPORT(polygonToCellsExperimental)(&COUNTRIES[index], res, 0, hexagons); + H3_EXPORT(polygonToCellsExperimental)(&COUNTRIES[index], res, CONTAINMENT_CENTER, hexagons); + } + }); + + BENCHMARK(polygonToCells_AllCountries3, 5, { + for (int index = 0; index < ${polygons.length}; index++) { + H3_EXPORT(polygonToCellsExperimental)(&COUNTRIES[index], res, CONTAINMENT_FULL, hexagons); + } + }); + + BENCHMARK(polygonToCells_AllCountries4, 5, { + for (int index = 0; index < ${polygons.length}; index++) { + H3_EXPORT(polygonToCellsExperimental)(&COUNTRIES[index], res, CONTAINMENT_OVERLAPPING, hexagons); } }); diff --git a/src/apps/benchmarks/benchmarkPolygonToCellsExperimental.c b/src/apps/benchmarks/benchmarkPolygonToCellsExperimental.c index b0837c7fc..1034f3372 100644 --- a/src/apps/benchmarks/benchmarkPolygonToCellsExperimental.c +++ b/src/apps/benchmarks/benchmarkPolygonToCellsExperimental.c @@ -17,6 +17,7 @@ #include "benchmark.h" #include "h3api.h" #include "polyfill.h" +#include "polygon.h" // Fixtures LatLng sfVerts[] = { @@ -122,24 +123,84 @@ southernGeoPolygon.geoloop = southernGeoLoop; int64_t numHexagons; H3Index *hexagons; -BENCHMARK(polygonToCellsSF, 500, { - H3_EXPORT(maxPolygonToCellsSize)(&sfGeoPolygon, 9, 0, &numHexagons); +BENCHMARK(polygonToCellsSF_Center, 500, { + H3_EXPORT(maxPolygonToCellsSize) + (&sfGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons); hexagons = calloc(numHexagons, sizeof(H3Index)); - H3_EXPORT(polygonToCellsExperimental)(&sfGeoPolygon, 9, 0, hexagons); + H3_EXPORT(polygonToCellsExperimental) + (&sfGeoPolygon, 9, CONTAINMENT_CENTER, hexagons); free(hexagons); }); -BENCHMARK(polygonToCellsAlameda, 500, { - H3_EXPORT(maxPolygonToCellsSize)(&alamedaGeoPolygon, 9, 0, &numHexagons); +BENCHMARK(polygonToCellsSF_Full, 500, { + H3_EXPORT(maxPolygonToCellsSize) + (&sfGeoPolygon, 9, CONTAINMENT_FULL, &numHexagons); hexagons = calloc(numHexagons, sizeof(H3Index)); - H3_EXPORT(polygonToCellsExperimental)(&alamedaGeoPolygon, 9, 0, hexagons); + H3_EXPORT(polygonToCellsExperimental) + (&sfGeoPolygon, 9, CONTAINMENT_FULL, hexagons); free(hexagons); }); -BENCHMARK(polygonToCellsSouthernExpansion, 10, { - H3_EXPORT(maxPolygonToCellsSize)(&southernGeoPolygon, 9, 0, &numHexagons); +BENCHMARK(polygonToCellsSF_Overlapping, 500, { + H3_EXPORT(maxPolygonToCellsSize) + (&sfGeoPolygon, 9, CONTAINMENT_OVERLAPPING, &numHexagons); hexagons = calloc(numHexagons, sizeof(H3Index)); - H3_EXPORT(polygonToCellsExperimental)(&southernGeoPolygon, 9, 0, hexagons); + H3_EXPORT(polygonToCellsExperimental) + (&sfGeoPolygon, 9, CONTAINMENT_OVERLAPPING, hexagons); + free(hexagons); +}); + +BENCHMARK(polygonToCellsAlameda_Center, 500, { + H3_EXPORT(maxPolygonToCellsSize) + (&alamedaGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons); + hexagons = calloc(numHexagons, sizeof(H3Index)); + H3_EXPORT(polygonToCellsExperimental) + (&alamedaGeoPolygon, 9, CONTAINMENT_CENTER, hexagons); + free(hexagons); +}); + +BENCHMARK(polygonToCellsAlameda_Full, 500, { + H3_EXPORT(maxPolygonToCellsSize) + (&alamedaGeoPolygon, 9, CONTAINMENT_FULL, &numHexagons); + hexagons = calloc(numHexagons, sizeof(H3Index)); + H3_EXPORT(polygonToCellsExperimental) + (&alamedaGeoPolygon, 9, CONTAINMENT_FULL, hexagons); + free(hexagons); +}); + +BENCHMARK(polygonToCellsAlameda_Overlapping, 500, { + H3_EXPORT(maxPolygonToCellsSize) + (&alamedaGeoPolygon, 9, CONTAINMENT_OVERLAPPING, &numHexagons); + hexagons = calloc(numHexagons, sizeof(H3Index)); + H3_EXPORT(polygonToCellsExperimental) + (&alamedaGeoPolygon, 9, CONTAINMENT_OVERLAPPING, hexagons); + free(hexagons); +}); + +BENCHMARK(polygonToCellsSouthernExpansion_Center, 10, { + H3_EXPORT(maxPolygonToCellsSize) + (&southernGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons); + hexagons = calloc(numHexagons, sizeof(H3Index)); + H3_EXPORT(polygonToCellsExperimental) + (&southernGeoPolygon, 9, CONTAINMENT_CENTER, hexagons); + free(hexagons); +}); + +BENCHMARK(polygonToCellsSouthernExpansion_Full, 10, { + H3_EXPORT(maxPolygonToCellsSize) + (&southernGeoPolygon, 9, CONTAINMENT_FULL, &numHexagons); + hexagons = calloc(numHexagons, sizeof(H3Index)); + H3_EXPORT(polygonToCellsExperimental) + (&southernGeoPolygon, 9, CONTAINMENT_FULL, hexagons); + free(hexagons); +}); + +BENCHMARK(polygonToCellsSouthernExpansion_Overlapping, 10, { + H3_EXPORT(maxPolygonToCellsSize) + (&southernGeoPolygon, 9, CONTAINMENT_OVERLAPPING, &numHexagons); + hexagons = calloc(numHexagons, sizeof(H3Index)); + H3_EXPORT(polygonToCellsExperimental) + (&southernGeoPolygon, 9, CONTAINMENT_OVERLAPPING, hexagons); free(hexagons); }); diff --git a/src/apps/testapps/testH3Memory.c b/src/apps/testapps/testH3Memory.c index 96627b543..54a7c688e 100644 --- a/src/apps/testapps/testH3Memory.c +++ b/src/apps/testapps/testH3Memory.c @@ -26,6 +26,7 @@ #include "h3api.h" #include "latLng.h" #include "polyfill.h" +#include "polygon.h" #include "test.h" #include "utility.h" @@ -230,22 +231,22 @@ SUITE(h3Memory) { sfGeoPolygon.numHoles = 0; int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&sfGeoPolygon, 9, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); resetMemoryCounters(0); failAlloc = true; - H3Error err = H3_EXPORT(polygonToCellsExperimental)(&sfGeoPolygon, 9, 0, - hexagons); + H3Error err = H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 9, CONTAINMENT_CENTER, hexagons); t_assert(err == E_MEMORY_ALLOC, "polygonToCellsExperimental failed (1)"); t_assert(actualAllocCalls == 1, "alloc called once"); t_assert(actualFreeCalls == 0, "free not called"); resetMemoryCounters(1); - err = H3_EXPORT(polygonToCellsExperimental)(&sfGeoPolygon, 9, 0, - hexagons); + err = H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 9, CONTAINMENT_CENTER, hexagons); t_assert(err == E_SUCCESS, "polygonToCellsExperimental succeeded (1)"); t_assert(actualAllocCalls == 1, "alloc called one time"); t_assert(actualFreeCalls == 1, "free called one time"); diff --git a/src/apps/testapps/testPolyfillInternal.c b/src/apps/testapps/testPolyfillInternal.c index db393cc30..4e44255ad 100644 --- a/src/apps/testapps/testPolyfillInternal.c +++ b/src/apps/testapps/testPolyfillInternal.c @@ -14,11 +14,14 @@ * limitations under the License. */ +#include + #include "bbox.h" #include "h3Index.h" #include "h3api.h" #include "latLng.h" #include "polyfill.h" +#include "polygon.h" #include "test.h" #include "utility.h" @@ -33,16 +36,24 @@ static GeoPolygon sfGeoPolygon = { {0.6599990002976, -2.1376771158464}}}, .numHoles = 0}; +static GeoPolygon invalidGeoPolygon = { + .geoloop = {.numVerts = 4, + .verts = (LatLng[]){{NAN, -2.1364398519396}, + {0.6595011102219, NAN}, + {NAN, -2.1354884206045}, + {0.6581220034068, NAN}}}, + .numHoles = 0}; + SUITE(polyfillInternal) { TEST(iterInitPolygonCompact_errors) { IterCellsPolygonCompact iter; - iter = iterInitPolygonCompact(&sfGeoPolygon, -1, 0); + iter = iterInitPolygonCompact(&sfGeoPolygon, -1, CONTAINMENT_CENTER); t_assert(iter.error == E_RES_DOMAIN, "Got expected error for invalid res"); t_assert(iter.cell == H3_NULL, "Got null output for invalid res"); - iter = iterInitPolygonCompact(&sfGeoPolygon, 16, 0); + iter = iterInitPolygonCompact(&sfGeoPolygon, 16, CONTAINMENT_CENTER); t_assert(iter.error == E_RES_DOMAIN, "Got expected error for invalid res"); t_assert(iter.cell == H3_NULL, "Got null output for invalid res"); @@ -53,11 +64,12 @@ SUITE(polyfillInternal) { t_assert(iter.cell == H3_NULL, "Got null output for invalid flags"); } - TEST(iterStepPolygonCompact_errors) { + TEST(iterStepPolygonCompact_invalidCellErrors) { IterCellsPolygonCompact iter; H3Index cell; - iter = iterInitPolygonCompact(&sfGeoPolygon, 9, 0); + iter = iterInitPolygonCompact(&sfGeoPolygon, 9, CONTAINMENT_CENTER); + t_assertSuccess(iter.error); // Give the iterator a cell with a bad base cell cell = 0x85283473fffffff; @@ -69,7 +81,37 @@ SUITE(polyfillInternal) { "Got expected error for invalid cell"); t_assert(iter.cell == H3_NULL, "Got null output for invalid cell"); - iter = iterInitPolygonCompact(&sfGeoPolygon, 9, 0); + iter = iterInitPolygonCompact(&sfGeoPolygon, 9, CONTAINMENT_CENTER); + t_assertSuccess(iter.error); + + // Give the iterator a cell with a bad base cell, at the target res + cell = 0x89283470003ffff; + H3_SET_BASE_CELL(cell, 123); + iter.cell = cell; + + iterStepPolygonCompact(&iter); + t_assert(iter.error == E_CELL_INVALID, + "Got expected error for invalid cell"); + t_assert(iter.cell == H3_NULL, + "Got null output for invalid cell at res"); + + iter = iterInitPolygonCompact(&sfGeoPolygon, 9, CONTAINMENT_FULL); + t_assertSuccess(iter.error); + + // Give the iterator a cell with a bad base cell, at the target res + // (full containment) + cell = 0x89283470003ffff; + H3_SET_BASE_CELL(cell, 123); + iter.cell = cell; + + iterStepPolygonCompact(&iter); + t_assert(iter.error == E_CELL_INVALID, + "Got expected error for invalid cell"); + t_assert(iter.cell == H3_NULL, + "Got null output for invalid cell at res"); + + iter = iterInitPolygonCompact(&sfGeoPolygon, 9, CONTAINMENT_CENTER); + t_assertSuccess(iter.error); // Give the iterator a cell that's too fine for a child check, // and a target resolution that allows this to run. This cell has @@ -84,9 +126,28 @@ SUITE(polyfillInternal) { t_assert(iter.cell == H3_NULL, "Got null output for invalid cell"); } + TEST(iterStepPolygonCompact_invalidPolygonErrors) { + IterCellsPolygonCompact iter; + + // Start with a good polygon, otherwise we error out early + iter = + iterInitPolygonCompact(&sfGeoPolygon, 5, CONTAINMENT_OVERLAPPING); + t_assertSuccess(iter.error); + + // Give the iterator a bad polygon and a cell at target res + iter._polygon = &invalidGeoPolygon; + iter.cell = 0x85283473fffffff; + + iterStepPolygonCompact(&iter); + t_assert(iter.error == E_LATLNG_DOMAIN, + "Got expected error for invalid polygon"); + t_assert(iter.cell == H3_NULL, "Got null output for invalid cell"); + } + TEST(iterDestroyPolygonCompact) { IterCellsPolygonCompact iter = - iterInitPolygonCompact(&sfGeoPolygon, 9, 0); + iterInitPolygonCompact(&sfGeoPolygon, 9, CONTAINMENT_CENTER); + t_assertSuccess(iter.error); iterDestroyPolygonCompact(&iter); t_assert(iter.error == E_SUCCESS, "Got success for destroyed iterator"); @@ -101,7 +162,9 @@ SUITE(polyfillInternal) { } TEST(iterDestroyPolygon) { - IterCellsPolygon iter = iterInitPolygon(&sfGeoPolygon, 9, 0); + IterCellsPolygon iter = + iterInitPolygon(&sfGeoPolygon, 9, CONTAINMENT_CENTER); + t_assertSuccess(iter.error); iterDestroyPolygon(&iter); t_assert(iter.error == E_SUCCESS, "Got success for destroyed iterator"); diff --git a/src/apps/testapps/testPolygonInternal.c b/src/apps/testapps/testPolygonInternal.c index 4b0f80286..42d19d14a 100644 --- a/src/apps/testapps/testPolygonInternal.c +++ b/src/apps/testapps/testPolygonInternal.c @@ -626,55 +626,55 @@ SUITE(polygonInternal) { H3_EXPORT(destroyLinkedMultiPolygon)(&polygon); } - TEST(lineIntersectsLine) { + TEST(lineCrossesLine) { LatLng lines1[4] = {{0, 0}, {1, 1}, {0, 1}, {1, 0}}; t_assert( - lineIntersectsLine(&lines1[0], &lines1[1], &lines1[2], &lines1[3]), + lineCrossesLine(&lines1[0], &lines1[1], &lines1[2], &lines1[3]), "diagonal intersection"); LatLng lines2[4] = {{1, 1}, {0, 0}, {1, 0}, {0, 1}}; t_assert( - lineIntersectsLine(&lines2[0], &lines2[1], &lines2[2], &lines2[3]), + lineCrossesLine(&lines2[0], &lines2[1], &lines2[2], &lines2[3]), "diagonal intersection, reverse vertexes"); LatLng lines3[4] = {{0.5, 0}, {0.5, 1}, {0, 0.5}, {1, 0.5}}; t_assert( - lineIntersectsLine(&lines3[0], &lines3[1], &lines3[2], &lines3[3]), + lineCrossesLine(&lines3[0], &lines3[1], &lines3[2], &lines3[3]), "horizontal/vertical intersection"); LatLng lines4[4] = {{0.5, 1}, {0.5, 0}, {1, 0.5}, {0, 0.5}}; t_assert( - lineIntersectsLine(&lines4[0], &lines4[1], &lines4[2], &lines4[3]), + lineCrossesLine(&lines4[0], &lines4[1], &lines4[2], &lines4[3]), "horizontal/vertical intersection, reverse vertexes"); LatLng lines5[4] = {{0, 0}, {0.4, 0.4}, {0, 1}, {1, 0}}; t_assert( - !lineIntersectsLine(&lines5[0], &lines5[1], &lines5[2], &lines5[3]), + !lineCrossesLine(&lines5[0], &lines5[1], &lines5[2], &lines5[3]), "diagonal non-intersection, below"); LatLng lines6[4] = {{0.6, 0.6}, {1, 1}, {0, 1}, {1, 0}}; t_assert( - !lineIntersectsLine(&lines6[0], &lines6[1], &lines6[2], &lines6[3]), + !lineCrossesLine(&lines6[0], &lines6[1], &lines6[2], &lines6[3]), "diagonal non-intersection, above"); LatLng lines7[4] = {{0.5, 0}, {0.5, 1}, {0, 0.5}, {0.4, 0.5}}; t_assert( - !lineIntersectsLine(&lines7[0], &lines7[1], &lines7[2], &lines7[3]), + !lineCrossesLine(&lines7[0], &lines7[1], &lines7[2], &lines7[3]), "horizontal/vertical non-intersection, below"); LatLng lines8[4] = {{0.5, 0}, {0.5, 1}, {0.6, 0.5}, {1, 0.5}}; t_assert( - !lineIntersectsLine(&lines8[0], &lines8[1], &lines8[2], &lines8[3]), + !lineCrossesLine(&lines8[0], &lines8[1], &lines8[2], &lines8[3]), "horizontal/vertical non-intersection, above"); LatLng lines9[4] = {{0.5, 0}, {0.5, 0.4}, {0, 0.5}, {1, 0.5}}; t_assert( - !lineIntersectsLine(&lines9[0], &lines9[1], &lines9[2], &lines9[3]), + !lineCrossesLine(&lines9[0], &lines9[1], &lines9[2], &lines9[3]), "horizontal/vertical non-intersection, left"); LatLng lines10[4] = {{0.5, 0.6}, {0.5, 1}, {0, 0.5}, {1, 0.5}}; - t_assert(!lineIntersectsLine(&lines10[0], &lines10[1], &lines10[2], - &lines10[3]), + t_assert(!lineCrossesLine(&lines10[0], &lines10[1], &lines10[2], + &lines10[3]), "horizontal/vertical non-intersection, right"); } diff --git a/src/apps/testapps/testPolygonToCells.c b/src/apps/testapps/testPolygonToCells.c index 22a0cdecf..b35dedd58 100644 --- a/src/apps/testapps/testPolygonToCells.c +++ b/src/apps/testapps/testPolygonToCells.c @@ -461,19 +461,21 @@ SUITE(polygonToCells) { TEST(invalidFlags) { int64_t numHexagons; - for (uint32_t flags = 1; flags <= 32; flags++) { + for (uint32_t flags = 3; flags <= 32; flags++) { t_assert( H3_EXPORT(maxPolygonToCellsSize)( &sfGeoPolygon, 9, flags, &numHexagons) == E_OPTION_INVALID, - "Flags other than 0 are invalid for maxPolygonToCellsSize"); + "Flags other than polyfill modes are invalid for " + "maxPolygonToCellsSize"); } t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&sfGeoPolygon, 9, 0, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - for (uint32_t flags = 1; flags <= 32; flags++) { + for (uint32_t flags = 3; flags <= 32; flags++) { t_assert(H3_EXPORT(polygonToCells)(&sfGeoPolygon, 9, flags, hexagons) == E_OPTION_INVALID, - "Flags other than 0 are invalid for polygonToCells"); + "Flags other than polyfill modes are invalid for " + "polygonToCells"); } free(hexagons); } diff --git a/src/apps/testapps/testPolygonToCellsExperimental.c b/src/apps/testapps/testPolygonToCellsExperimental.c index f389456ab..24464597d 100644 --- a/src/apps/testapps/testPolygonToCellsExperimental.c +++ b/src/apps/testapps/testPolygonToCellsExperimental.c @@ -22,6 +22,7 @@ #include "h3Index.h" #include "latLng.h" #include "polyfill.h" +#include "polygon.h" #include "test.h" #include "utility.h" @@ -99,7 +100,7 @@ static void fillIndex_assertions(H3Index h) { H3Index *polygonToCellsOut = calloc(polygonToCellsSize, sizeof(H3Index)); t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( - &polygon, nextRes, 0, polygonToCellsOut)); + &polygon, nextRes, CONTAINMENT_CENTER, polygonToCellsOut)); int64_t polygonToCellsCount = countNonNullIndexes(polygonToCellsOut, polygonToCellsSize); @@ -160,26 +161,56 @@ SUITE(polygonToCells) { TEST(polygonToCells) { int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&sfGeoPolygon, 9, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&sfGeoPolygon, 9, - 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 9, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == 1253, "got expected polygonToCells size"); free(hexagons); } + TEST(polygonToCells_FullContainment) { + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 9, CONTAINMENT_FULL, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 9, CONTAINMENT_FULL, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 1175, + "got expected polygonToCells size (full containment mode)"); + free(hexagons); + } + + TEST(polygonToCells_Overlapping) { + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 9, CONTAINMENT_OVERLAPPING, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 9, CONTAINMENT_OVERLAPPING, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 1334, + "got expected polygonToCells size (overlapping mode)"); + free(hexagons); + } + TEST(polygonToCellsHole) { int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&holeGeoPolygon, 9, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &holeGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&holeGeoPolygon, - 9, 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &holeGeoPolygon, 9, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == 1214, @@ -187,14 +218,45 @@ SUITE(polygonToCells) { free(hexagons); } + TEST(polygonToCellsHoleFullContainment) { + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &holeGeoPolygon, 9, CONTAINMENT_FULL, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &holeGeoPolygon, 9, CONTAINMENT_FULL, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert( + actualNumIndexes == 1118, + "got expected polygonToCells size (hole, full containment mode)"); + free(hexagons); + } + + TEST(polygonToCellsHoleOverlapping) { + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &holeGeoPolygon, 9, CONTAINMENT_OVERLAPPING, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &holeGeoPolygon, 9, CONTAINMENT_OVERLAPPING, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 1311, + "got expected polygonToCells size (hole, overlapping mode)"); + free(hexagons); + } + TEST(polygonToCellsEmpty) { int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&emptyGeoPolygon, 9, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &emptyGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&emptyGeoPolygon, - 9, 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &emptyGeoPolygon, 9, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == 0, @@ -202,6 +264,77 @@ SUITE(polygonToCells) { free(hexagons); } + TEST(polygonToCellsContainsPolygon) { + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 4, CONTAINMENT_CENTER, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 4, CONTAINMENT_CENTER, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 0, "got expected polygonToCells size"); + free(hexagons); + } + + TEST(polygonToCellsContainsPolygon_CenterContained) { + // Contains the center point of a res 4 cell + LatLng centerVerts[] = {{0.6595645, -2.1353315}, + {0.6595645, -2.1353314}, + {0.6595644, -2.1353314}, + {0.6595644, -2.1353314265}}; + GeoLoop centerGeoLoop = {.numVerts = 4, .verts = centerVerts}; + GeoPolygon centerGeoPolygon; + centerGeoPolygon.geoloop = centerGeoLoop; + centerGeoPolygon.numHoles = 0; + + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + ¢erGeoPolygon, 4, CONTAINMENT_CENTER, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + ¢erGeoPolygon, 4, CONTAINMENT_CENTER, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 1, "got expected polygonToCells size"); + t_assert(hexagons[0] == 0x8428309ffffffff, "got expected hexagon"); + + free(hexagons); + } + + TEST(polygonToCellsContainsPolygon_FullContainment) { + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 4, CONTAINMENT_FULL, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 4, CONTAINMENT_FULL, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 0, + "got expected polygonToCells size (full containment mode)"); + free(hexagons); + } + + TEST(polygonToCellsContainsPolygon_Overlapping) { + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 4, CONTAINMENT_OVERLAPPING, &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &sfGeoPolygon, 4, CONTAINMENT_OVERLAPPING, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 1, + "got expected polygonToCells size (overlapping mode)"); + t_assert(hexagons[0] == 0x8428309ffffffff, "got expected hexagon"); + free(hexagons); + } + TEST(polygonToCellsExact) { LatLng somewhere = {1, 2}; H3Index origin; @@ -223,15 +356,27 @@ SUITE(polygonToCells) { someHexagon.numHoles = 0; int64_t numHexagons; - t_assertSuccess( - H3_EXPORT(maxPolygonToCellsSize)(&someHexagon, 9, 0, &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &someHexagon, 9, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&someHexagon, 9, - 0, hexagons)); - int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + int64_t actualNumIndexes; + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &someHexagon, 9, CONTAINMENT_CENTER, hexagons)); + actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + t_assert(actualNumIndexes == 1, + "got expected polygonToCells size for center containment (1)"); + + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &someHexagon, 9, CONTAINMENT_FULL, hexagons)); + actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + t_assert(actualNumIndexes == 1, + "got expected polygonToCells size for full containment (1)"); + + // TODO: CONTAINMENT_OVERLAPPING yields 7 cells, presumably due to FPE + // in the various cell boundaries - t_assert(actualNumIndexes == 1, "got expected polygonToCells size (1)"); free(hexagons); free(verts); } @@ -272,11 +417,11 @@ SUITE(polygonToCells) { expectedSize = 4228; int64_t numHexagons; t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( - &primeMeridianGeoPolygon, 7, 0, &numHexagons)); + &primeMeridianGeoPolygon, 7, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( - &primeMeridianGeoPolygon, 7, 0, hexagons)); + &primeMeridianGeoPolygon, 7, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == expectedSize, @@ -287,11 +432,11 @@ SUITE(polygonToCells) { // differences in hex size and grid offset between the two cases expectedSize = 4238; t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( - &transMeridianGeoPolygon, 7, 0, &numHexagons)); + &transMeridianGeoPolygon, 7, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagonsTM = calloc(numHexagons, sizeof(H3Index)); t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( - &transMeridianGeoPolygon, 7, 0, hexagonsTM)); + &transMeridianGeoPolygon, 7, CONTAINMENT_CENTER, hexagonsTM)); actualNumIndexes = countNonNullIndexes(hexagonsTM, numHexagons); t_assert(actualNumIndexes == expectedSize, @@ -300,21 +445,23 @@ SUITE(polygonToCells) { // Transmeridian filled hole case -- only needed for calculating hole // size t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( - &transMeridianFilledHoleGeoPolygon, 7, 0, &numHexagons)); + &transMeridianFilledHoleGeoPolygon, 7, CONTAINMENT_CENTER, + &numHexagons)); H3Index *hexagonsTMFH = calloc(numHexagons, sizeof(H3Index)); t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( - &transMeridianFilledHoleGeoPolygon, 7, 0, hexagonsTMFH)); + &transMeridianFilledHoleGeoPolygon, 7, CONTAINMENT_CENTER, + hexagonsTMFH)); int64_t actualNumHoleIndexes = countNonNullIndexes(hexagonsTMFH, numHexagons); // Transmeridian hole case t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( - &transMeridianHoleGeoPolygon, 7, 0, &numHexagons)); + &transMeridianHoleGeoPolygon, 7, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagonsTMH = calloc(numHexagons, sizeof(H3Index)); t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( - &transMeridianHoleGeoPolygon, 7, 0, hexagonsTMH)); + &transMeridianHoleGeoPolygon, 7, CONTAINMENT_CENTER, hexagonsTMH)); actualNumIndexes = countNonNullIndexes(hexagonsTMH, numHexagons); t_assert(actualNumIndexes == expectedSize - actualNumHoleIndexes, @@ -337,12 +484,12 @@ SUITE(polygonToCells) { GeoPolygon polygon = {.geoloop = geoloop, .numHoles = 0}; int64_t numHexagons; - t_assertSuccess( - H3_EXPORT(maxPolygonToCellsSize)(&polygon, 4, 0, &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &polygon, 4, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess( - H3_EXPORT(polygonToCellsExperimental)(&polygon, 4, 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &polygon, 4, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); @@ -389,12 +536,12 @@ SUITE(polygonToCells) { polygon.numHoles = 0; int64_t numHexagons; - t_assertSuccess( - H3_EXPORT(maxPolygonToCellsSize)(&polygon, 9, 0, &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &polygon, 9, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess( - H3_EXPORT(polygonToCellsExperimental)(&polygon, 9, 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &polygon, 9, CONTAINMENT_CENTER, hexagons)); int found = 0; int numPentagons = 0; @@ -413,19 +560,21 @@ SUITE(polygonToCells) { TEST(invalidFlags) { int64_t numHexagons; - for (uint32_t flags = 1; flags <= 32; flags++) { + for (uint32_t flags = 3; flags <= 32; flags++) { t_assert( H3_EXPORT(maxPolygonToCellsSize)( &sfGeoPolygon, 9, flags, &numHexagons) == E_OPTION_INVALID, - "Flags other than 0 are invalid for maxPolygonToCellsSize"); + "Flags other than polyfill modes are invalid for " + "maxPolygonToCellsSize"); } - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&sfGeoPolygon, 9, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &sfGeoPolygon, 9, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - for (uint32_t flags = 1; flags <= 32; flags++) { + for (uint32_t flags = 3; flags <= 32; flags++) { t_assert(H3_EXPORT(polygonToCellsExperimental)( &sfGeoPolygon, 9, flags, hexagons) == E_OPTION_INVALID, - "Flags other than 0 are invalid for polygonToCells"); + "Flags other than polyfill modes are invalid for " + "polygonToCells"); } free(hexagons); } diff --git a/src/apps/testapps/testPolygonToCellsReportedExperimental.c b/src/apps/testapps/testPolygonToCellsReportedExperimental.c index 58e688504..b3fe8e5b8 100644 --- a/src/apps/testapps/testPolygonToCellsReportedExperimental.c +++ b/src/apps/testapps/testPolygonToCellsReportedExperimental.c @@ -21,6 +21,7 @@ #include "h3Index.h" #include "latLng.h" #include "polyfill.h" +#include "polygon.h" #include "test.h" #include "utility.h" @@ -42,23 +43,26 @@ SUITE(polygonToCells_reported) { for (int res = 0; res < 3; res++) { int64_t polygonToCellsSize; t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( - &worldGeoPolygon, res, 0, &polygonToCellsSize)); + &worldGeoPolygon, res, CONTAINMENT_CENTER, + &polygonToCellsSize)); H3Index *polygonToCellsOut = calloc(polygonToCellsSize, sizeof(H3Index)); t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( - &worldGeoPolygon, res, 0, polygonToCellsOut)); + &worldGeoPolygon, res, CONTAINMENT_CENTER, polygonToCellsOut)); int64_t actualNumIndexes = countNonNullIndexes(polygonToCellsOut, polygonToCellsSize); int64_t polygonToCellsSize2; t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( - &worldGeoPolygon2, res, 0, &polygonToCellsSize2)); + &worldGeoPolygon2, res, CONTAINMENT_CENTER, + &polygonToCellsSize2)); H3Index *polygonToCellsOut2 = calloc(polygonToCellsSize2, sizeof(H3Index)); t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( - &worldGeoPolygon2, res, 0, polygonToCellsOut2)); + &worldGeoPolygon2, res, CONTAINMENT_CENTER, + polygonToCellsOut2)); int64_t actualNumIndexes2 = countNonNullIndexes(polygonToCellsOut2, polygonToCellsSize2); @@ -105,12 +109,12 @@ SUITE(polygonToCells_reported) { int res = 7; int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&testPolygon, res, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &testPolygon, res, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&testPolygon, res, - 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &testPolygon, res, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == 4499, @@ -134,12 +138,12 @@ SUITE(polygonToCells_reported) { int res = 7; int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&testPolygon, res, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &testPolygon, res, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&testPolygon, res, - 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &testPolygon, res, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == 4609, @@ -160,12 +164,12 @@ SUITE(polygonToCells_reported) { int res = 13; int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&testPolygon, res, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &testPolygon, res, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&testPolygon, res, - 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &testPolygon, res, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == 4353, "got expected polygonToCells size"); @@ -193,12 +197,12 @@ SUITE(polygonToCells_reported) { int res = 5; int64_t numHexagons; - t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&testPolygon, res, 0, - &numHexagons)); + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)( + &testPolygon, res, CONTAINMENT_CENTER, &numHexagons)); H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); - t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)(&testPolygon, res, - 0, hexagons)); + t_assertSuccess(H3_EXPORT(polygonToCellsExperimental)( + &testPolygon, res, CONTAINMENT_CENTER, hexagons)); int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); t_assert(actualNumIndexes == 8, "got expected polygonToCells size"); diff --git a/src/h3lib/include/polygon.h b/src/h3lib/include/polygon.h index 387ed034a..d7d24ffbc 100644 --- a/src/h3lib/include/polygon.h +++ b/src/h3lib/include/polygon.h @@ -40,18 +40,38 @@ /** Macro: Whether a GeoLoop is empty */ #define IS_EMPTY_GEOFENCE(geoloop) geoloop->numVerts == 0 +/** + * Values representing polyfill containment modes, to be used in + * the `flags` bit field. + */ +typedef enum { + CONTAINMENT_CENTER = 0, ///< Cell center is contained in the shape + CONTAINMENT_FULL = 1, ///< Cell is fully contained in the shape + CONTAINMENT_OVERLAPPING = 2, ///< Cell overlaps the shape at any point + CONTAINMENT_INVALID = 3 ///< This mode is invalid and should not be used +} ContainmentMode; + +// 1s in the 4 bits defining the polyfill containment mode, 0s elsewhere +#define FLAG_CONTAINMENT_MODE_MASK ((uint32_t)(15)) +#define FLAG_GET_CONTAINMENT_MODE(flags) (flags & FLAG_CONTAINMENT_MODE_MASK) + // Defined directly in polygon.c: +H3Error validatePolygonFlags(uint32_t flags); void bboxesFromGeoPolygon(const GeoPolygon *polygon, BBox *bboxes); bool pointInsidePolygon(const GeoPolygon *geoPolygon, const BBox *bboxes, const LatLng *coord); bool cellBoundaryInsidePolygon(const GeoPolygon *geoPolygon, const BBox *bboxes, const CellBoundary *boundary, const BBox *boundaryBBox); +bool cellBoundaryCrossesPolygon(const GeoPolygon *geoPolygon, + const BBox *bboxes, + const CellBoundary *boundary, + const BBox *boundaryBBox); bool cellBoundaryCrossesGeoLoop(const GeoLoop *geoloop, const BBox *loopBBox, const CellBoundary *boundary, const BBox *boundaryBBox); -bool lineIntersectsLine(const LatLng *a1, const LatLng *a2, const LatLng *b1, - const LatLng *b2); +bool lineCrossesLine(const LatLng *a1, const LatLng *a2, const LatLng *b1, + const LatLng *b2); // The following functions are created via macro in polygonAlgos.h, // so their signatures are documented here: diff --git a/src/h3lib/lib/algos.c b/src/h3lib/lib/algos.c index 79666d843..949988a08 100644 --- a/src/h3lib/lib/algos.c +++ b/src/h3lib/lib/algos.c @@ -774,8 +774,9 @@ H3Error H3_EXPORT(gridRingUnsafe)(H3Index origin, int k, H3Index *out) { */ H3Error H3_EXPORT(maxPolygonToCellsSize)(const GeoPolygon *geoPolygon, int res, uint32_t flags, int64_t *out) { - if (flags != 0) { - return E_OPTION_INVALID; + H3Error flagErr = validatePolygonFlags(flags); + if (flagErr) { + return flagErr; } // Get the bounding box for the GeoJSON-like struct BBox bbox; @@ -890,8 +891,9 @@ H3Error _getEdgeHexagons(const GeoLoop *geoloop, int64_t numHexagons, int res, */ H3Error H3_EXPORT(polygonToCells)(const GeoPolygon *geoPolygon, int res, uint32_t flags, H3Index *out) { - if (flags != 0) { - return E_OPTION_INVALID; + H3Error flagErr = validatePolygonFlags(flags); + if (flagErr) { + return flagErr; } // One of the goals of the polygonToCells algorithm is that two adjacent // polygons with zero overlap have zero overlapping hexagons. That the diff --git a/src/h3lib/lib/polyfill.c b/src/h3lib/lib/polyfill.c index 7ba5a28b5..bde5eedbc 100644 --- a/src/h3lib/lib/polyfill.c +++ b/src/h3lib/lib/polyfill.c @@ -354,8 +354,9 @@ IterCellsPolygonCompact iterInitPolygonCompact(const GeoPolygon *polygon, return iter; } - if (flags != 0) { - iterErrorPolygonCompact(&iter, E_OPTION_INVALID); + H3Error flagErr = validatePolygonFlags(flags); + if (flagErr) { + iterErrorPolygonCompact(&iter, flagErr); return iter; } @@ -406,23 +407,80 @@ void iterStepPolygonCompact(IterCellsPolygonCompact *iter) { iter->_started = true; } + ContainmentMode mode = FLAG_GET_CONTAINMENT_MODE(iter->_flags); + while (cell) { int cellRes = H3_GET_RESOLUTION(cell); // Target res: Do a fine-grained check if (cellRes == iter->_res) { - // Check if the cell is in the polygon - // TODO: Handle other polyfill modes here - LatLng center; - H3Error centerErr = H3_EXPORT(cellToLatLng)(cell, ¢er); - if (NEVER(centerErr != E_SUCCESS)) { - iterErrorPolygonCompact(iter, centerErr); - return; + if (mode == CONTAINMENT_CENTER || mode == CONTAINMENT_OVERLAPPING) { + // Check if the cell center is inside the polygon + LatLng center; + H3Error centerErr = H3_EXPORT(cellToLatLng)(cell, ¢er); + if (centerErr != E_SUCCESS) { + iterErrorPolygonCompact(iter, centerErr); + return; + } + if (pointInsidePolygon(iter->_polygon, iter->_bboxes, + ¢er)) { + // Set to next output + iter->cell = cell; + return; + } } - if (pointInsidePolygon(iter->_polygon, iter->_bboxes, ¢er)) { - // Set to next output - iter->cell = cell; - return; + if (mode == CONTAINMENT_OVERLAPPING) { + // For overlapping, we need to do a quick check to determine + // whether the polygon is wholly contained by the cell. We check + // the first polygon vertex, which if it is contained could also + // mean we simply intersect. + H3Index polygonCell; + H3Error polygonCellErr = H3_EXPORT(latLngToCell)( + &(iter->_polygon->geoloop.verts[0]), cellRes, &polygonCell); + if (polygonCellErr != E_SUCCESS) { + iterErrorPolygonCompact(iter, polygonCellErr); + return; + } + if (polygonCell == cell) { + // Set to next output + iter->cell = cell; + return; + } + } + if (mode == CONTAINMENT_FULL || mode == CONTAINMENT_OVERLAPPING) { + CellBoundary boundary; + H3Error boundaryErr = + H3_EXPORT(cellToBoundary)(cell, &boundary); + if (boundaryErr != E_SUCCESS) { + iterErrorPolygonCompact(iter, boundaryErr); + return; + } + BBox bbox; + H3Error bboxErr = cellToBBox(cell, &bbox, false); + if (NEVER(bboxErr != E_SUCCESS)) { + // Should be unreachable - invalid cells would be caught in + // the previous boundaryErr + iterErrorPolygonCompact(iter, bboxErr); + return; + } + // Check if the cell is fully contained by the polygon + if (mode == CONTAINMENT_FULL && + cellBoundaryInsidePolygon(iter->_polygon, iter->_bboxes, + &boundary, &bbox)) { + // Set to next output + iter->cell = cell; + return; + } + // For overlap, we've already checked for center point inclusion + // above; if that failed, we only need to check for line + // intersection + else if (mode == CONTAINMENT_OVERLAPPING && + cellBoundaryCrossesPolygon( + iter->_polygon, iter->_bboxes, &boundary, &bbox)) { + // Set to next output + iter->cell = cell; + return; + } } } diff --git a/src/h3lib/lib/polygon.c b/src/h3lib/lib/polygon.c index a92fecfba..1944e5bda 100644 --- a/src/h3lib/lib/polygon.c +++ b/src/h3lib/lib/polygon.c @@ -42,6 +42,20 @@ #undef ITERATE #undef IS_EMPTY +/** + * Whether the flags for the polyfill operation are valid + * TODO: Move to polyfill.c when the old algo is removed + * @param flags Flags to validate + * @return Whether the flags are valid + */ +H3Error validatePolygonFlags(uint32_t flags) { + if (flags & (~FLAG_CONTAINMENT_MODE_MASK) || + FLAG_GET_CONTAINMENT_MODE(flags) >= CONTAINMENT_INVALID) { + return E_OPTION_INVALID; + } + return E_SUCCESS; +} + /** * Create a bounding box from a GeoPolygon * @param polygon Input GeoPolygon @@ -115,6 +129,33 @@ bool cellBoundaryInsidePolygon(const GeoPolygon *geoPolygon, const BBox *bboxes, return true; } +/** + * Whether any part of a cell boundary crosses a polygon. Crossing in this case + * means whether any line segments intersect; it does not include containment. + * @param geoPolygon The polygon to test + * @param bboxes The bboxes for the main geoloop and each of its holes + * @param boundary The cell boundary to test + * @return Whether the cell boundary is contained + */ +bool cellBoundaryCrossesPolygon(const GeoPolygon *geoPolygon, + const BBox *bboxes, + const CellBoundary *boundary, + const BBox *boundaryBBox) { + // Check for line intersections with outer loop + if (cellBoundaryCrossesGeoLoop(&(geoPolygon->geoloop), &bboxes[0], boundary, + boundaryBBox)) { + return true; + } + // Check for line intersections with any hole + for (int i = 0; i < geoPolygon->numHoles; i++) { + if (cellBoundaryCrossesGeoLoop(&(geoPolygon->holes[i]), &bboxes[i + 1], + boundary, boundaryBBox)) { + return true; + } + } + return false; +} + /** * Whether a cell boundary crosses a geo loop. Crossing in this case means * whether any line segments intersect; it does not include containment. @@ -167,7 +208,7 @@ bool cellBoundaryCrossesGeoLoop(const GeoLoop *geoloop, const BBox *loopBBox, } for (int j = 0; j < normalBoundary.numVerts; j++) { - if (lineIntersectsLine( + if (lineCrossesLine( &loop1, &loop2, &normalBoundary.verts[j], &normalBoundary.verts[(j + 1) % normalBoundary.numVerts])) { return true; @@ -187,8 +228,8 @@ bool cellBoundaryCrossesGeoLoop(const GeoLoop *geoloop, const BBox *loopBBox, * @param b2 End of line B * @return Whether the lines intersect */ -bool lineIntersectsLine(const LatLng *a1, const LatLng *a2, const LatLng *b1, - const LatLng *b2) { +bool lineCrossesLine(const LatLng *a1, const LatLng *a2, const LatLng *b1, + const LatLng *b2) { double denom = ((b2->lng - b1->lng) * (a2->lat - a1->lat) - (b2->lat - b1->lat) * (a2->lng - a1->lng)); if (!denom) return false;