diff --git a/src/apps/testapps/testPolygon.c b/src/apps/testapps/testPolygon.c index 9f3015d90..df9bc39d7 100644 --- a/src/apps/testapps/testPolygon.c +++ b/src/apps/testapps/testPolygon.c @@ -48,16 +48,84 @@ SUITE(polygon) { BBox bbox; bboxFromGeoLoop(&geoloop, &bbox); + // For exact points on the polygon, we bias west and south t_assert(!pointInsideGeoLoop(&geoloop, &bbox, &sfVerts[0]), - "contains exact"); - t_assert(pointInsideGeoLoop(&geoloop, &bbox, &sfVerts[4]), - "contains exact 4"); + "does not contain exact vertex 0"); + t_assert(pointInsideGeoLoop(&geoloop, &bbox, &sfVerts[3]), + "contains exact vertex 3"); + t_assert(pointInsideGeoLoop(&geoloop, &bbox, &inside), "contains point inside"); t_assert(!pointInsideGeoLoop(&geoloop, &bbox, &somewhere), "contains somewhere else"); } + TEST(pointInsideGeoLoopCornerCases) { + LatLng verts[] = {{0, 0}, {1, 0}, {1, 1}, {0, 1}}; + GeoLoop geoloop = {.numVerts = 4, .verts = verts}; + + BBox bbox; + bboxFromGeoLoop(&geoloop, &bbox); + + LatLng point = {0, 0}; + + // Test corners. For exact points on the polygon, we bias west and + // north, so only the southeast corner is contained. + t_assert(!pointInsideGeoLoop(&geoloop, &bbox, &point), + "does not contain sw corner"); + point.lat = 1; + t_assert(!pointInsideGeoLoop(&geoloop, &bbox, &point), + "does not contain nw corner "); + point.lng = 1; + t_assert(!pointInsideGeoLoop(&geoloop, &bbox, &point), + "does not contain ne corner "); + point.lat = 0; + t_assert(pointInsideGeoLoop(&geoloop, &bbox, &point), + "contains se corner "); + } + + TEST(pointInsideGeoLoopEdgeCases) { + LatLng verts[] = {{0, 0}, {1, 0}, {1, 1}, {0, 1}}; + GeoLoop geoloop = {.numVerts = 4, .verts = verts}; + + BBox bbox; + bboxFromGeoLoop(&geoloop, &bbox); + + LatLng point; + + // Test edges. Only points on south and east edges are contained. + point.lat = 0.5; + point.lng = 0; + t_assert(!pointInsideGeoLoop(&geoloop, &bbox, &point), + "does not contain point on west edge"); + point.lat = 1; + point.lng = 0.5; + t_assert(!pointInsideGeoLoop(&geoloop, &bbox, &point), + "does not contain point on north edge"); + point.lat = 0.5; + point.lng = 1; + t_assert(pointInsideGeoLoop(&geoloop, &bbox, &point), + "contains point on east edge"); + point.lat = 0; + point.lng = 0.5; + t_assert(pointInsideGeoLoop(&geoloop, &bbox, &point), + "contains point on south edge"); + } + + TEST(pointInsideGeoLoopExtraEdgeCase) { + // This is a carefully crafted shape + point to hit an otherwise + // missed branch in coverage + LatLng verts[] = {{0, 0}, {1, 0.5}, {0, 1}}; + GeoLoop geoloop = {.numVerts = 4, .verts = verts}; + + BBox bbox; + bboxFromGeoLoop(&geoloop, &bbox); + + LatLng point = {0.5, 0.5}; + t_assert(pointInsideGeoLoop(&geoloop, &bbox, &point), + "contains inside point matching longitude of a vertex"); + } + TEST(pointInsideGeoLoopTransmeridian) { LatLng verts[] = {{0.01, -M_PI + 0.01}, {0.01, M_PI - 0.01}, diff --git a/src/apps/testapps/testPolygonToCellsReported.c b/src/apps/testapps/testPolygonToCellsReported.c index 644122f86..8e4d9a9ed 100644 --- a/src/apps/testapps/testPolygonToCellsReported.c +++ b/src/apps/testapps/testPolygonToCellsReported.c @@ -170,4 +170,37 @@ SUITE(polygonToCells_reported) { t_assert(actualNumIndexes == 4353, "got expected polygonToCells size"); free(hexagons); } + + // https://github.com/uber/h3/issues/595 + TEST(h3_136) { + H3Index center = 0x85283473fffffff; + LatLng centerLatLng; + H3_EXPORT(cellToLatLng)(center, ¢erLatLng); + + // This polygon should include the center cell. The issue here arises + // when one of the polygon vertexes is to the east of the index center, + // with exactly the same latitude + LatLng testVerts[] = {{.lat = centerLatLng.lat, -2.121207808248113}, + {0.6565301558937859, -2.1281107217935986}, + {0.6515463604919347, -2.1345342663428695}, + {0.6466583305904194, -2.1276313527973842}}; + + GeoLoop testGeoLoop = {.numVerts = 4, .verts = testVerts}; + GeoPolygon testPolygon; + testPolygon.geoloop = testGeoLoop; + testPolygon.numHoles = 0; + + int res = 5; + int64_t numHexagons; + t_assertSuccess(H3_EXPORT(maxPolygonToCellsSize)(&testPolygon, res, 0, + &numHexagons)); + H3Index *hexagons = calloc(numHexagons, sizeof(H3Index)); + + t_assertSuccess( + H3_EXPORT(polygonToCells)(&testPolygon, res, 0, hexagons)); + int64_t actualNumIndexes = countNonNullIndexes(hexagons, numHexagons); + + t_assert(actualNumIndexes == 8, "got expected polygonToCells size"); + free(hexagons); + } } diff --git a/src/h3lib/include/polygonAlgos.h b/src/h3lib/include/polygonAlgos.h index b24c99c3d..e6ce98403 100644 --- a/src/h3lib/include/polygonAlgos.h +++ b/src/h3lib/include/polygonAlgos.h @@ -92,6 +92,20 @@ bool GENERIC_LOOP_ALGO(pointInside)(const TYPE *loop, const BBox *bbox, b = tmp; } + // If the latitude matches exactly, we'll hit an edge case where + // the ray passes through the vertex twice on successive segment + // checks. To avoid this, adjust the latiude northward if needed. + // + // NOTE: This currently means that a point at the north pole cannot + // be contained in any polygon. This is acceptable in current usage, + // because the point we test in this function at present is always + // a cell center or vertex, and no cell has a center or vertex on the + // north pole. If we need to expand this algo to more generic uses we + // might need to handle this edge case. + if (lat == a.lat || lat == b.lat) { + lat += DBL_EPSILON; + } + // If we're totally above or below the latitude ranges, the test // ray cannot intersect the line segment, so let's move on if (lat < a.lat || lat > b.lat) {