Skip to content

Bidirectional gridPathCells #1111

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

Open
Zaczero wants to merge 1 commit into
base: master
Choose a base branch
from
+144 −37
Open

Bidirectional gridPathCells #1111

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
48 changes: 44 additions & 4 deletions src/apps/testapps/testGridPathCells.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@
* usage: `testGridPathCells`
*/

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
Expand All @@ -29,6 +30,20 @@
#include "test.h"
#include "utility.h"

static void assertPathValid(H3Index start, H3Index end, const H3Index *path,
int64_t size) {
t_assert(size > 0, "path size is positive");
t_assert(path[0] == start, "path starts with start index");
t_assert(path[size - 1] == end, "path ends with end index");

for (int64_t i = 1; i < size; i++) {
int isNeighbor;
t_assertSuccess(
H3_EXPORT(areNeighborCells)(path[i], path[i - 1], &isNeighbor));
t_assert(isNeighbor, "path is contiguous");
}
}

SUITE(gridPathCells) {
TEST(gridPathCells_acrossMultipleFaces) {
H3Index start = 0x85285aa7fffffff;
Expand All @@ -40,15 +55,40 @@ SUITE(gridPathCells) {
"Line not computable across multiple icosa faces");
}

TEST(gridPathCells_pentagon) {
TEST(gridPathCells_pentagonReverseInterpolation) {
H3Index start = 0x820807fffffffff;
H3Index end = 0x8208e7fffffffff;

int64_t size;
t_assertSuccess(H3_EXPORT(gridPathCellsSize)(start, end, &size));
H3Index *path = calloc(sizeof(H3Index), size);
t_assert(H3_EXPORT(gridPathCells)(start, end, path) == E_PENTAGON,
"Line not computable due to pentagon distortion");
H3Index *path = calloc(size, sizeof(H3Index));

t_assertSuccess(H3_EXPORT(gridPathCells)(start, end, path));
assertPathValid(start, end, path, size);

free(path);
}

TEST(gridPathCells_knownFailureNotCoveredByReverseInterpolation) {
// Known limitation case:
//
// There are still pairs where `gridDistance` succeeds but interpolation
// fails in both origin-anchored local IJK charts (anchored at `start`
// and anchored at `end`). Since `gridPathCells` only performs these two
// interpolation attempts, it returns an error.
//
// This test keeps a pinned input pair to demonstrate that the current
// approach is not complete.
H3Index start = 0x8411b61ffffffff;
H3Index end = 0x84016d3ffffffff;

int64_t size;
t_assertSuccess(H3_EXPORT(gridPathCellsSize)(start, end, &size));
H3Index *path = calloc(size, sizeof(H3Index));

H3Error err = H3_EXPORT(gridPathCells)(start, end, path);
t_assert(err != E_SUCCESS, "Expected gridPathCells to fail");

free(path);
}
}
133 changes: 100 additions & 33 deletions src/h3lib/lib/localij.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -22,9 +22,12 @@
#include <assert.h>
#include <inttypes.h>
#include <math.h>
#include <stdbool.h>
#include <stdlib.h>
#include <string.h>

#include "algos.h"
#include "alloc.h"
#include "baseCells.h"
#include "faceijk.h"
#include "h3Assert.h"
Expand Down Expand Up @@ -398,9 +401,10 @@ H3Error localIjkToCell(H3Index origin, const CoordIJK *ijk, H3Index *out) {
for (int i = 0; i < pentagonRotations; i++) {
dir = _rotate60ccw(dir);
}
// The pentagon rotations are being chosen so that dir is not the
// deleted direction. If it still happens, it means we're moving
// into a deleted subsequence, so there is no index here.
// The pentagon rotations are chosen to avoid the deleted direction
// (the missing neighbor direction around a pentagon). If we still
// land on it, the coordinate would cross pentagon distortion and
// cannot be represented.
if (dir == K_AXES_DIGIT) {
return E_PENTAGON;
}
Expand Down Expand Up @@ -603,13 +607,15 @@ H3Error H3_EXPORT(gridDistance)(H3Index origin, H3Index index, int64_t *out) {

/**
* Number of indexes in a line from the start index to the end index,
* to be used for allocating memory. Returns a negative number if the
* line cannot be computed.
* to be used for allocating memory.
*
* On success, sets `*size` to `gridDistance(start, end) + 1`
* (including both endpoints).
*
* @param start Start index of the line
* @param end End index of the line
* @param size Size of the line
* @returns 0 on success, or another value on error
* @param size Output size of the line
* @returns 0 on success, otherwise the error from `gridDistance`
*/
H3Error H3_EXPORT(gridPathCellsSize)(H3Index start, H3Index end,
int64_t *size) {
Expand Down Expand Up @@ -654,34 +660,32 @@ static void cubeRound(double i, double j, double k, CoordIJK *ijk) {
}

/**
* Given two H3 indexes, return the line of indexes between them (inclusive).
*
* This function may fail to find the line between two indexes, for
* example if they are very far apart. It may also fail when finding
* distances for indexes on opposite sides of a pentagon.
* Attempts to generate a shortest-length path by interpolating through an
* origin-anchored local IJK coordinate space.
*
* Notes:
* This helper implements the interpolation-based path construction used by
* `gridPathCells`. It can fail if interpolation lands on intermediate IJK
* coordinates that cannot be mapped back to valid H3 cells. This can occur
* because the origin-anchored local IJ(K) space is not globally continuous, and
* some intermediate coordinates do not have an inverse mapping back to a cell
* in the chosen chart (for example, due to discontinuities or warping near
* pentagons).
*
* - The specific output of this function should not be considered stable
* across library versions. The only guarantees the library provides are
* that the line length will be `gridDistance(start, end) + 1` and that
* every index in the line will be a neighbor of the preceding index.
* - Lines are drawn in grid space, and may not correspond exactly to either
* Cartesian lines or great arcs.
* The output is written to `out[outOffset + outStep * n]`, allowing callers to
* fill the path in either direction without an intermediate buffer.
*
* @param start Start index of the line
* @param end End index of the line
* @param out Output array, which must be of size gridPathCellsSize(start, end)
* @return 0 on success, or another value on failure.
* @param start Origin cell for local IJK conversions
* @param end Target cell
* @param distance Expected edge distance between `start` and `end`
* @param out Output buffer
* @param outOffset Output index for the first element
* @param outStep Output stride (+1 for forward fill, -1 for reverse fill)
* @return E_SUCCESS if all intermediate steps convert successfully, otherwise
* the first encountered conversion error
*/
H3Error H3_EXPORT(gridPathCells)(H3Index start, H3Index end, H3Index *out) {
int64_t distance;
H3Error distanceError = H3_EXPORT(gridDistance)(start, end, &distance);
// Early exit if we can't calculate the line
if (distanceError) {
return distanceError;
}

static H3Error gridPathCellsInterpolate(H3Index start, H3Index end,
int64_t distance, H3Index *out,
int64_t outOffset, int64_t outStep) {
// Get IJK coords for the start and end. We've already confirmed
// that these can be calculated with the distance check above.
CoordIJK startIjk = {0};
Expand All @@ -703,7 +707,7 @@ H3Error H3_EXPORT(gridPathCells)(H3Index start, H3Index end, H3Index *out) {
ijkToCube(&startIjk);
ijkToCube(&endIjk);

double invDistance = distance ? 1.0 / (double)distance : 0;
double invDistance = 1.0 / (double)distance;
double iStep = (double)(endIjk.i - startIjk.i) * invDistance;
double jStep = (double)(endIjk.j - startIjk.j) * invDistance;
double kStep = (double)(endIjk.k - startIjk.k) * invDistance;
Expand All @@ -715,7 +719,8 @@ H3Error H3_EXPORT(gridPathCells)(H3Index start, H3Index end, H3Index *out) {
(double)startIjk.k + kStep * n, &currentIjk);
// Convert cube -> ijk -> h3 index
cubeToIjk(&currentIjk);
H3Error currentError = localIjkToCell(start, &currentIjk, &out[n]);
const int64_t idx = outOffset + outStep * n;
H3Error currentError = localIjkToCell(start, &currentIjk, &out[idx]);
if (currentError) {
// The cells between `start` and `end` may fall in pentagon
// distortion.
Expand All @@ -725,3 +730,65 @@ H3Error H3_EXPORT(gridPathCells)(H3Index start, H3Index end, H3Index *out) {

return E_SUCCESS;
}

/**
* Given two H3 indexes, return the line of indexes between them (inclusive).
*
* This function relies on `gridDistance(start, end)` to determine the expected
* path length, and will return the same error if `gridDistance` fails.
*
* Path construction is performed by straight-line interpolation in the
* origin-anchored local IJK coordinate space:
*
* - First, interpolate using `start` as the local IJK origin.
* - If interpolation fails, retry using `end` as the local IJK origin and
* reverse the resulting sequence into `out`.
*
* If both interpolation attempts fail, this function returns the error from the
* first attempt.
*
* Notes:
*
* - The specific output of this function should not be considered stable
* across library versions. The only guarantees the library provides are
* that the line length will be `gridDistance(start, end) + 1` and that
* every index in the line will be a neighbor of the preceding index.
* - Lines are drawn in grid space, and may not correspond exactly to either
* Cartesian lines or great arcs.
*
* @param start Start index of the line
* @param end End index of the line
* @param out Output array, which must be of size gridPathCellsSize(start, end)
* @return 0 on success, or another value on failure.
*/
H3Error H3_EXPORT(gridPathCells)(H3Index start, H3Index end, H3Index *out) {
int64_t distance;
H3Error distanceError = H3_EXPORT(gridDistance)(start, end, &distance);
// Early exit if we can't calculate the line
if (distanceError) {
return distanceError;
}

if (distance == 0) {
out[0] = start;
return E_SUCCESS;
}

// Straight-line interpolation in local IJK space anchored at `start`.
H3Error interpolateErr =
gridPathCellsInterpolate(start, end, distance, out, 0, 1);
if (!interpolateErr) {
return E_SUCCESS;
}

// Retry interpolation anchored at `end` and reverse the output.
// This can resolve cases where the local IJK chart is discontinuous
// relative to one origin but not the other.
H3Error reverseErr =
gridPathCellsInterpolate(end, start, distance, out, distance, -1);
if (!reverseErr) {
return E_SUCCESS;
}

return interpolateErr;
}