Add function for returning the H3 indices of each endpoint of a directed edge

[nwroyer]

After looking through the code, it appears that calculating the IJK coordinates for a face is part of the process for computing a cell boundary, and functions exist internally to convert these coordinates to H3 indices. I believe exposing a function to convert a directed edge to endpoint indices directly is an easy add to the library, and I have one drafted, but want to confirm this is not already implemented in a way I'm not spotting.

The use case that sparks this need is as follows:

• I have a large dataset that's binned by H3 cells. My goal is to draw the boundaries between cells that don't share a value in the dataset - call it the color of a cell for the sake of illustration. Rendering this is trivial once I have a set of sequenced points similar to a GeoJSON LineString. For my purposes, the directed nature of the edges H3 has is fine - a boundary inside cell A, and a boundary inside cell B, where A and B are adjacent, are actually preferred.
• The easy answer to this problem would be to set up a mapping from cell colors to cell indices, and then use polygonToCells to get a boundary. This works on a small scale. Unfortunately, the dataset I have is too large for this solution.
• I'm already rendering underlying elevation data with a system that renders cells at resolution $N$ using chunks at resolution $\lfloor \sqrt{\text{max}(0, N - 1)} \rfloor$. This has proven to be performant enough for my purposes, even with elevation stored at the same resolution as cell color.
• I could, per chunk, create a mapping of cell colors to indices and use polygonToCells as mentioned above. This process works fine. This will result in borders along the edge of a chunk, though, even when the adjacent cells share the same color.
• I can, on a per chunk basis, get the boundary segments of all the directed edges that separate cells with two different colors, then assemble these into linestrings with the following algorithm:
  • Push the first segment onto the return list
  • For each subsequent segment in the raw list, check if it shares an endpoint with one in the return list, and either extend the return segment found or add the segment as a new segment in the return list.
  • Mark closed loops at the end to avoid rendering the same vertex twice (not essential.)
• The above involves a lot of vector distance comparisons, even using squared magnitudes.

Creating/exposing the function mentioned will simplify this into a series of integer equality comparisons, and I suspect there are likely other use cases for allowing users to directly get the vertices associated with edges.

Preliminary testing indicates that the below seems to do the trick. I'll put together a proper pull request if this isn't duplicating effort.

/**
 * Provides the vertices at either end of the directed edge.
 * @param edge The directed edge H3Index
 * @param vStart The starting vertex H3Index
 * @param vEnd The starting vertex H3Index
 */
H3Error H3_EXPORT(directedEdgeToVertices)(H3Index edge, H3Index* vStart, H3Index* vEnd) {
    // Get the origin and neighbor direction from the edge
    Direction direction = H3_GET_RESERVED_BITS(edge);
    H3Index origin;
    H3Error originResult = H3_EXPORT(getDirectedEdgeOrigin)(edge, &origin);
    if (originResult) {
        return originResult;
    }

    // Get the start vertex for the edge
    int startVertex = vertexNumForDirection(origin, direction);
    if (startVertex == INVALID_VERTEX_NUM) {
        // This is not actually an edge (i.e. no valid direction),
        // so return no vertices.
        return E_DIR_EDGE_INVALID;
    }

    H3Error err = H3_EXPORT(cellToVertex)(origin, startVertex, vStart);
    if (err) {
        return err;
    }
    err = H3_EXPORT(cellToVertex)(origin, startVertex + 1, vEnd);
    if (err) {
        return err;
    }
    return E_SUCCESS;

[nwroyer]

Further experimentation shows that the above doesn't account for the differing order of vertices and edges. To solve this, I believe it's possible to just get the starting vertex of the next edge clockwise around the cell.

The below code seems to work fine at multiple resolutions from 0 to 7 in my application; I'll spin up proper tests and a previous-edge function as well.

static const Direction nextEdgeClockwiseDirectionsHex[8] = {
    // Direction 0, center
    CENTER_DIGIT,

    // Direction 1, k-axes direction
    JK_AXES_DIGIT, // K_AXES_DIGIT -> JK_AXES_DIGIT,
    // Direction 2, j-axes direction
    IJ_AXES_DIGIT, // J_AXES_DIGIT -> IJ_AXES_DIGIT
    // Direction 3, j == k direction
    J_AXES_DIGIT, // JK_AXES_DIGIT -> J_AXES_DIGIT
    // Direction 4, i-axes direction
    IK_AXES_DIGIT, // I_AXES_DIGIT -> IK_AXES_DIGIT
    // Direction 5, i == k direction
    K_AXES_DIGIT, // IK_AXES_DIGIT -> K_AXES_DIGIT
    // Direction 6, i == j direction
    I_AXES_DIGIT, // IJ_AXES_DIGIT -> I_AXES_DIGIT

    // Invalid digit (7)
    INVALID_DIGIT // INVALID_DIGIT -> INVALID_DIGIT
};

static const Direction nextEdgeClockwiseDirectionsPent[8] = {
    // Direction 0, center
    CENTER_DIGIT,

    // Direction 1, k-axes direction
    IK_AXES_DIGIT, // K_AXES_DIGIT -> IK_AXES_DIGIT
    // Direction 2, j-axes direction
    K_AXES_DIGIT, // J_AXES_DIGIT -> K_AXES_DIGIT
    // Direction 3, j == k direction
    I_AXES_DIGIT, // JK_AXES_DIGIT -> I_AXES_DIGIT
    // Direction 4, i-axes direction
    J_AXES_DIGIT, // I_AXES_DIGIT -> J_AXES_DIGIT
    // Direction 5, i == k direction
    JK_AXES_DIGIT, // IK_AXES_DIGIT -> JK_AXES_DIGIT
    // Direction 6, i == j direction
    INVALID_DIGIT, // IJ_AXES_DIGIT (invalid for pentagons) -> INVALID_DIGIT

    // Invalid digit (7)
    INVALID_DIGIT // INVALID_DIGIT -> INVALID_DIGIT
};

// ...

/**
 * Returns the next edge in the clockwise order around the origin cell
 * @param edge The directed edge H3Index
 * @param nextClockwiseEdge The next clockwise edge output
 */
H3Error H3_EXPORT(directedEdgeGetNextClockwise)(H3Index edge, H3Index* nextClockwiseEdge) {
    Direction direction = H3_GET_RESERVED_BITS(edge);
    if (direction == CENTER_DIGIT || direction == INVALID_DIGIT) {
        return E_DIR_EDGE_INVALID;
    }

    H3Index origin;
    H3Error originResult = H3_EXPORT(getDirectedEdgeOrigin)(edge, &origin);
    if (originResult) {
        return originResult;
    }

    *nextClockwiseEdge = edge;
    int isPent = H3_EXPORT(isPentagon)(origin);
    if (isPent) {
        Direction nextDirection = nextEdgeClockwiseDirectionsPent[direction];
        H3_SET_RESERVED_BITS(*nextClockwiseEdge, nextDirection);
    } else {
        Direction nextDirection = nextEdgeClockwiseDirectionsHex[direction];
        H3_SET_RESERVED_BITS(*nextClockwiseEdge, nextDirection);
    }

    return E_SUCCESS;
}

/**
 * Provides the vertices at either end of the directed edge.
 * @param edge The directed edge H3Index
 * @param vStart The starting vertex H3Index
 * @param vEnd The starting vertex H3Index
 */
H3Error H3_EXPORT(directedEdgeToVertices)(H3Index edge, H3Index* vStart, H3Index* vEnd) {
    // Get the origin and neighbor direction from the edge
    Direction direction = H3_GET_RESERVED_BITS(edge);
    H3Index origin;
    H3Error originResult = H3_EXPORT(getDirectedEdgeOrigin)(edge, &origin);
    if (originResult) {
        return originResult;
    }

    // Get the start vertex for the edge
    int startVertex = vertexNumForDirection(origin, direction);
    if (startVertex == INVALID_VERTEX_NUM) {
        // This is not actually an edge (i.e. no valid direction),
        // so return no vertices.
        return E_DIR_EDGE_INVALID;
    }

    H3Error firstVertexErr = H3_EXPORT(cellToVertex)(origin, startVertex, vStart);
    if (firstVertexErr) {
        return firstVertexErr;
    }

    // Handle the end vertex for the edge
    int isPent = H3_EXPORT(isPentagon)(origin);
    int nextDirection;
    if (isPent) {
        nextDirection = nextEdgeClockwiseDirectionsPent[direction];
    } else {
        nextDirection = nextEdgeClockwiseDirectionsHex[direction];
    }

    int endVertex = vertexNumForDirection(origin, nextDirection);
    H3Error secondVertexErr = H3_EXPORT(cellToVertex)(origin, endVertex, vEnd);
    if (secondVertexErr) {
        return secondVertexErr;
    }

    return E_SUCCESS;
}

[nrabinowitz]

If I'm understanding correctly, this is already implemented as directedEdgeToBoundary. Is there something that function doesn't provide that you need for your use case?

[nrabinowitz]

My goal is to draw the boundaries between cells that don't share a value in the dataset

Another option here is to use cellsToMultiPolygon, which will draw a single polygon for multiple cells. You'd need to cluster your cells by shared value, then send each cluster to cellsToMultiPolygon to get the appropriate outline. This is likely easier than dealing with edges directly.

[nwroyer]

Neither of the existing functions you mention really get at what I need, or would require significant additional processing to achieve the result I'm looking for.

Both functions mentioned give me, ultimately, the boundary around a group of cells, where the only information on what category adjacent cells belong to is "is this is in the provided set or not?" What it doesn't give me is the set of boundary edges or vertices between cells with different values in the provided set, where I don't care about the edges between cells in the provided set and outside.

Consider the following diagram:

Here I'm only interested in retrieving the purple edges. This is already doable by iterating over the provided set of cells and retaining only those edges that are both:

• Between cells that are both in the set
• Between cells of different values

Turning this into a rendering solution is where the problem comes in. The size of the dataset means that I have an interest in rendering these borders with as few vertices in the boundary as possible. directedEdgeToBoundary or cellsToMultipolgyon means that to filter out the edges not in purple, I'll need to potentially iterate over a large number of points in the result of the function that will then need to be discarded.

Furthermore, filtering out coincident points on the purple edges will mean two floating-point comparisons per vertex, instead of a single equality check if I'm able to get vertices associated with each edge directly, impacting performance.

The second problem comes with large datasets. I can't fit in memory an array with all the cells in a particular set - think of binning the world by nations at resolution 7+ or so. What I can do is split the entire world up into chunks that I know do fit in memory (say, children of cells at resolution 2) and then compute the interior (purple) borders per chunk.

Testing this approach with country data has shown that it's viable - when rendering a slice of the world at resolution 7, I'm able to pull only those cells in resolution 2 within view and then use this chunk-wise method to compute borders.

I've actually achieved what I set out to do using the code from my last post, and am still preparing the tests necessary for a proper pull request. I believe there's probably value in users being able to associate vertices with edge data anyway, even beyond my admittedly niche use case. We can associated edges with cells, and cells with vertices, but unless I'm mistaken there's no way to associate edges with vertices directly.

[nrabinowitz]

Ok, I think I understand now. I apologize, I skimmed over this at first and didn't realize you were dealing directly with vertex-mode vertexes instead of lat/lng vertexes. I'm still not 100% sure I understand the use case, which I still suspect could be solved as efficiently without the new method, but I can definitely see that directedEdgeToVertexes is something we'd want to add (I think it was originally planned but we never added it in after adding vertex mode), and we'd be happy to accept a PR. I have a few notes on the code you've posted above, but it would be easier to add them in code review.