The rhp_wrappers Module

Fig. 21 The H3-style wrappers over New Zealand: polyfill (left) returns the cells whose centroids fall inside a polygon; linetrace (right) returns the cells a line passes through, in order. The color change partway down the South Island is the boundary between the equatorial R cells (quads) and the south polar S cells (skew quads).

Fig. 22 linetrace across the north polar cap, viewed from above the pole (Arctic coastlines in grey: northern Greenland at lower left, Svalbard below, Franz Josef Land at lower right; the resolution 4 sub-grid is drawn as a faint outline). A linestring’s segments are straight in longitude-latitude coordinates, so at these latitudes each leg curves around the pole rather than hopping across it, passing through the circular cap cell N444 and the dart and skew-quad cells surrounding it – the cell shapes for which no longitude-latitude rectangle is a valid approximation.

Fig. 23 cell_ring: away from the cube’s corners (left, around Q44) the ring at distance k has the familiar 8k cells. Around a cell at a cube corner (right, around N00, corner marked) the rings are smaller – the k=1 ring has 7 cells – because only three cells meet at the corner, so the ring spans three cube faces without a fourth, diagonally-opposite cell.

Fig. 24 The wrap_antimeridian flag on linetrace. Coordinates are treated as plate carree, so the segment (179, 10) to (-179, 10) literally crosses the whole plane through longitude 0 (top, the default). With wrap_antimeridian=True (bottom) the segment is interpreted as taking the short way across the antimeridian. Users who split their geometries at the antimeridian beforehand never need the flag.

Line semantics

Two points do not by themselves say which line joins them, and the three answers give different cells. line names the one wanted: "geodesic" for the shortest path on the ellipsoid, "plate_carree" for a line straight in longitude and latitude, and "plane" for one straight in the rHEALPix plane. The older plane flag still selects between the last two.

A linestring digitised in a geographic CRS usually means geodesic segments between its vertices, which is what line="geodesic" traces, exactly and without densifying. wrap_antimeridian does not apply to a geodesic, which already takes the short way, and an antipodal pair raises ValueError because no single geodesic between them is shortest.

Fig. 25 Each line is straight in exactly one of the three spaces, and the panels are those spaces. The plate carree line stays flat near 47 degrees north while the geodesic and the planar line arc over the pole. All three trace sixteen cells here, but not the same sixteen.

Fig. 26 Once a line crosses faces the readings separate further. The geodesic and the plate carree line stay close; the straight planar line does not, and part of it leaves the grid’s cross-shaped planar image altogether, so the cells it returns are not a connected path. cells_from_line warns when that happens.

Attributing length to cells

line_crossings() takes the same arguments and returns, per cell in the order cells_from_line gives them, the fraction of the line’s length on the ellipsoid lying in that cell, so a feature’s length can be shared out among the cells it crosses. A cell visited twice gets an entry per visit, and the fractions sum to 1.

cell_area(rhpindex, unit, ] =, plane, dggs)	Returns the area of a cell in the requested unit (or None if rhpindex is invalid).
cell_ring(rhpindex, k, dggs)	Returns the ring of cell indices around rhpindex at distance k (in graph distance -- shortest number of edge/corner hops -- not straight-line distance), or None if rhpindex is invalid.
geo_to_rhp(lat, lng, resolution, plane, dggs)	Turn a latitute and longitude (in degrees) into an rHEALPix cell address at the requested resolution.
k_ring(rhpindex, k, dggs)	Returns the k-ring (the centre cell together with every ring out to distance k, in graph distance -- see cell_ring()) of cell indices around rhpindex, or None if rhpindex is invalid.
linetrace(geometry, res, plane, verbose, ...)	Returns the list of cell indices touched by a shapely linestring or multilinestring at the requested resolution.
polyfill(geometry, res, plane, compress, ...)	Turns the area contained in a shapely polygon or multipolygon into a set of cell indices at the requested resolution.
rhp_get_base_cell(rhpindex, dggs)	Returns the resolution 0 cell id of a given cell index (or None if invalid).
rhp_get_resolution(rhpindex, dggs)	Returns the resolution of a given cell index (or None if invalid).
rhp_is_valid(rhpindex, dggs)	Checks if the given cell address is valid within the DGGS, i.e. whether dggs.parse_index accepts it.
rhp_to_center_child(rhpindex, res, verbose, dggs)	Returns central child of rhpindex at resolution res (immediate central child if res == None).
rhp_to_geo(rhpindex, geo_json, plane, dggs)	Look up the centroid (in degrees) of the cell identified by rhpindex: the cell's direct position in the sense of OGC Topic 21 (see Cell.centroid).
rhp_to_geo_boundary(rhpindex, geo_json, ...)	Extract the corner coordinates of a cell at a given cell ID and returns them as a tuple of coordinate pairs (in degrees).
rhp_to_parent(rhpindex, res, verbose, dggs)	Returns parent of rhpindex at resolution res (immediate parent if res == None).

rhealpixdggs.rhp_wrappers.cell_area(rhpindex: str, unit: Literal['km^2', 'm^2']='km^2', plane: bool = True, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → float | None

Returns the area of a cell in the requested unit (or None if rhpindex is invalid).

TODO: investigate use case where unit is ‘rads^2’

EXAMPLES:

>>> cell_area('S001450634', unit='m^2', plane=False)
219428.08748562282
>>> cell_area('S001450635', unit='m^2', plane=False)
219428.08748562282
>>> cell_area('S00145063', unit='km^2', plane=False)
1.9748527873706059
>>> cell_area('INVALID', unit='km^2', plane=False)

rhealpixdggs.rhp_wrappers.cell_ring(rhpindex: str, k: int = 1, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → list[str] | None

Returns the ring of cell indices around rhpindex at distance k (in graph distance – shortest number of edge/corner hops – not straight-line distance), or None if rhpindex is invalid.

Also returns None if k < 0.

Returns [ rhpindex ] if k == 0 (by convention).

Built by growing the ring outward one shell at a time from neighbor() (edge-adjacent) and diagonal_neighbor() (corner-adjacent) steps, so this is correct near a cube corner too: exactly 3 cells meet at each of the cube’s 8 corners rather than 4, so a ring passing near one has fewer cells there than the usual 8*k for an interior ring, rather than an incorrect or duplicate cell.

EXAMPLES:

>>> cell_ring('S001450634')
['S001450631', 'S001450632', 'S001450635', 'S001450638', 'S001450637', 'S001450636', 'S001450633', 'S001450630']
>>> cell_ring('S001450634', k=2)
['S001450607', 'S001450608', 'S001450606', 'S001450616', 'S001450640', 'S001450643', 'S001450646', 'S001450670', 'S001450662', 'S001450661', 'S001450660', 'S001442882', 'S001442858', 'S001442855', 'S001442852', 'S001442828']
>>> cell_ring('S001450634', k=0)
['S001450634']
>>> cell_ring('S', k=1)
['O', 'P', 'Q', 'R']
>>> cell_ring('S', k=-1)

rhealpixdggs.rhp_wrappers.geo_to_rhp(lat: float, lng: float, resolution: int, plane: bool = True, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → str | None

Turn a latitute and longitude (in degrees) into an rHEALPix cell address at the requested resolution.

Mostly passes through the parameters to the function turning coordinate points into cells, but converts the address tuple from the resulting cell into a string.

Returns None if no cell matching the coordinates is found.

EXAMPLES::

>>> geo_to_rhp(-43.738058, -176.260506, 9, plane=False)
'S001450634'
>>> geo_to_rhp(-43.738379,-176.258807, 9, plane=False)
'S001450635'

rhealpixdggs.rhp_wrappers.k_ring(rhpindex: str, k: int = 1, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → list[str] | None

Returns the k-ring (the centre cell together with every ring out to distance k, in graph distance – see cell_ring()) of cell indices around rhpindex, or None if rhpindex is invalid.

Also returns None if k < 0.

EXAMPLES:

>>> k_ring('S001450634')
['S001450634', 'S001450631', 'S001450632', 'S001450635', 'S001450638', 'S001450637', 'S001450636', 'S001450633', 'S001450630']
>>> k_ring('S001450634', k=2)
['S001450634', 'S001450631', 'S001450632', 'S001450635', 'S001450638', 'S001450637', 'S001450636', 'S001450633', 'S001450630', 'S001450607', 'S001450608', 'S001450606', 'S001450616', 'S001450640', 'S001450643', 'S001450646', 'S001450670', 'S001450662', 'S001450661', 'S001450660', 'S001442882', 'S001442858', 'S001442855', 'S001442852', 'S001442828']
>>> k_ring('S001450634', k=0)
['S001450634']
>>> k_ring('S001450634', k=-1)
>>> k_ring('INVALID')

rhealpixdggs.rhp_wrappers.linetrace(geometry: LineString | MultiLineString, res: int, plane: bool = True, verbose: bool = False, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>, wrap_antimeridian: bool = False, line: Literal['plane', 'plate_carree', 'geodesic'] | None=None) → list[str] | None

Returns the list of cell indices touched by a shapely linestring or multilinestring at the requested resolution. Removes internal sequences of duplicate cells before returning result.

Returns None if the geom_type field in the input geometry is anything other than ‘LineString’ or ‘MultiLineString’.

Returns None if the geometry is empty, or if it has no length.

Returns None if no cells match the geometry for some reason.

Returns None if the geometry is invalid in other ways, e.g. if a linestring contains self intersecting segments.

A linestring’s segments are straight in the input coordinate space (shapely geometries are planar): planar rHEALPix coordinates if plane = True, longitude-latitude coordinates otherwise – so longitude-latitude segments are plate carree lines, not geodesics, and by default do not wrap around the antimeridian (a segment from longitude 179 to -179 runs the long way around, through longitude 0 – the literal planar reading). Splitting inputs at the antimeridian, as GeoJSON’s RFC 7946 prescribes for data producers, avoids the ambiguity entirely; alternatively pass wrap_antimeridian = True to trace segments spanning more than half a turn of longitude the short way, across the antimeridian, which traces the same cells as splitting them there. To trace a geodesic, densify it into short segments first (e.g. with pyproj.Geod). See RHEALPixDGGS.cells_from_line, which traces each segment exactly, for the details.

EXAMPLES:

>>> from shapely import LineString
>>> line = LineString([[-176.260506, -43.738058], [-176.258807, -43.738379]])
>>> linetrace(line, res=9, plane=False)
['S001450634', 'S001450635']

rhealpixdggs.rhp_wrappers.polyfill(geometry: Polygon | MultiPolygon, res: int, plane: bool = True, compress: bool = False, verbose: bool = False, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>, containment: Literal['center', 'full', 'overlapping']='center', min_res: int = 0) → set[str] | None

Turns the area contained in a shapely polygon or multipolygon into a set of cell indices at the requested resolution. Which cells count is chosen by containment, following H3’s modes:

• "center" (default): a cell is included if its centroid is inside the geometry defined by the boundaries and holes.

• "full": a cell is included if the whole cell lies within the geometry (its boundary may touch the geometry’s).

• "overlapping": a cell is included if any part of it meets the geometry, so the result covers the geometry.

full and overlapping test the cell’s boundary polygon. In the plane that polygon is exact. On the ellipsoid, equatorial cells are longitude-latitude aligned and exact too; polar cells have curved edges, which are represented by 6 points per edge (Cell.boundary(n=6)), so a decision there is made against that piecewise-linear approximation. A cap cell, which spans every longitude, is fully inside only a geometry that covers the whole latitude band from the cap’s boundary to the pole, and overlaps any geometry reaching that latitude.

The candidate cells are those of the geometry’s bounding box (RHEALPixDGGS.cells_in_box), and every test runs on all candidates at once: centroids from RHEALPixDGGS.centroids, boundaries from RHEALPixDGGS.boundary_array, predicates from shapely’s vectorised functions.

Returns an empty set if no cell centroids fall within the input geometry.

Returns None if the geom_type field in the input geometry is anything other than ‘Polygon’ or ‘MultiPolygon’.

Returns None if the geometry is empty, or if it has no area.

Returns None if no cells match the geometry for some reason.

Returns None if the geometry is invalid in other ways, e.g. if a point on a hole boundary is outside the exterior boundary of its polygon, or if two polygons in a multipolygon overlap.

Polygon edges are straight in the input coordinate space (shapely geometries are planar) and do not wrap around the antimeridian; a region straddling it must be split into a MultiPolygon along the antimeridian by the caller first, as is standard for planar GIS geometry.

The fill descends the hierarchy from resolution min_res: a cell wholly inside the geometry is kept whole, a cell missing it is dropped, and only the cells the boundary passes through are split, down to resolution res, where the remaining cells are decided by containment. So the work grows with the length of the boundary, not the area. With compress = False every kept coarse cell is expanded into its resolution-res descendants, giving exactly the cells a cell-by-cell test at res would give. With compress = True the coarse cells are returned as they are, a mixed-resolution cover with no cell coarser than min_res and no complete group of siblings left unmerged; so every cell of the result has exactly one ancestor at resolution min_res, which suits partitioning by that resolution. “Wholly inside” is judged by the cell’s boundary polygon for full and overlapping where that is exact (the plane, equatorial cells), and otherwise by the cell’s longitude-latitude bounding box, which contains the cell: for polar cells, whose 6-point polygon is only approximate, and for center, where the box also contains every descendant’s centroid, so the centroid rule is kept exactly.

polyfill_array returns the same cells as a sorted numpy string array, which costs less than half the memory per cell for large fills.

EXAMPLES::

>>> from shapely import Polygon
>>> coords = ((0., 0.), (0., 1.), (1., 1.), (1., 0.), (0., 0.))
>>> polygon = Polygon(coords)
>>> polyfill(polygon, res=1, plane=False)
set()
>>> polyfill(polygon, res=2, plane=False)
set()
>>> polyfill(polygon, res=3, plane=False)
set()
>>> polyfill(polygon, res=4, plane=False)
{'Q3330'}
>>> sorted(polyfill(polygon, res=5, plane=False))
['Q33303', 'Q33304', 'Q33305', 'Q33306', 'Q33307', 'Q33308', 'Q33330', 'Q33331', 'Q33332']
>>> sorted(polyfill(polygon, res=6, plane=False))
['Q333033', 'Q333034', 'Q333035', 'Q333036', 'Q333037', 'Q333038', 'Q333043', 'Q333044', 'Q333045', 'Q333046', 'Q333047', 'Q333048', 'Q333053', 'Q333054', 'Q333056', 'Q333057', 'Q333060', 'Q333061', 'Q333062', 'Q333063', 'Q333064', 'Q333065', 'Q333066', 'Q333067', 'Q333068', 'Q333070', 'Q333071', 'Q333072', 'Q333073', 'Q333074', 'Q333075', 'Q333076', 'Q333077', 'Q333078', 'Q333080', 'Q333081', 'Q333083', 'Q333084', 'Q333086', 'Q333087', 'Q333300', 'Q333301', 'Q333302', 'Q333303', 'Q333304', 'Q333305', 'Q333306', 'Q333307', 'Q333308', 'Q333310', 'Q333311', 'Q333312', 'Q333313', 'Q333314', 'Q333315', 'Q333316', 'Q333317', 'Q333318', 'Q333320', 'Q333321', 'Q333323', 'Q333324', 'Q333326', 'Q333327', 'Q333330', 'Q333331', 'Q333332', 'Q333333', 'Q333334', 'Q333335', 'Q333340', 'Q333341', 'Q333342', 'Q333343', 'Q333344', 'Q333345', 'Q333350', 'Q333351', 'Q333353', 'Q333354']

Four of those res-6 groups have all 9 siblings present (Q33306x, Q33307x, Q33330x, Q33331x), so compress=True collapses them one level to res-5 parents. The remaining partial groups stay at res-6, giving a mixed-resolution result:

>>> sorted(c for c in polyfill(polygon, res=6, plane=False, compress=True) if len(c) < 7)
['Q33306', 'Q33307', 'Q33330', 'Q33331']

At res=7 the same four cells are again fully covered (all 81 of each cell’s res-7 grandchildren fall inside the polygon), so compaction runs two levels: res-7 siblings collapse to res-6 parents, those in turn collapse to res-5:

>>> result7 = polyfill(polygon, res=7, plane=False, compress=True)
>>> 'Q33306' in result7  # res-5 cell present despite res=7 request
True
>>> 'Q333060' not in result7  # intermediate res-6 cell absorbed
True
>>> 'Q3330600' not in result7  # original res-7 cell absorbed
True

rhealpixdggs.rhp_wrappers.polyfill_array(geometry: Polygon | MultiPolygon, res: int, plane: bool = True, compress: bool = False, verbose: bool = False, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>, containment: Literal['center', 'full', 'overlapping']='center', min_res: int = 0) → ndarray | None

polyfill returning the cell indices as a numpy string array instead of a set: the same cells for the same arguments, sorted (by base cell, then digit by digit, the order Cell objects sort in; with compress a coarse cell comes where its descendants would) and without duplicates, or None where polyfill returns None. An empty array means no cell qualified.

The array costs about 4 * (res + 2) bytes per cell against roughly a hundred for a set entry, and feeds array-based callers such as RHEALPixDGGS.boundary_array, RHEALPixDGGS.centroids or a data frame column directly.

EXAMPLES::

>>> from shapely import Polygon
>>> polygon = Polygon(((0., 0.), (0., 1.), (1., 1.), (1., 0.), (0., 0.)))
>>> polyfill_array(polygon, res=5, plane=False).tolist()
['Q33303', 'Q33304', 'Q33305', 'Q33306', 'Q33307', 'Q33308', 'Q33330', 'Q33331', 'Q33332']
>>> polyfill_array(polygon, res=2, plane=False).size
0
>>> polyfill_array(Polygon(), res=2, plane=False) is None
True

rhealpixdggs.rhp_wrappers.rhp_get_base_cell(rhpindex: str, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → str | None

Returns the resolution 0 cell id of a given cell index (or None if invalid).

EXAMPLES::

>>> rhp_get_base_cell('S001450634')
'S'
>>> rhp_get_base_cell('INVALID')

rhealpixdggs.rhp_wrappers.rhp_get_resolution(rhpindex: str, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → int | None

Returns the resolution of a given cell index (or None if invalid).

EXAMPLES:

>>> rhp_get_resolution('S001450634')
9
>>> rhp_get_resolution('INVALID')

rhealpixdggs.rhp_wrappers.rhp_is_valid(rhpindex: str, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → bool

Checks if the given cell address is valid within the DGGS, i.e. whether dggs.parse_index accepts it. Like every function taking an index string, raises ValueError if dggs has no index strings (N_side other than 2 or 3).

EXAMPLES::

>>> rhp_is_valid('S001450634')
True
>>> rhp_is_valid('N')
True
>>> rhp_is_valid('INVALID')
False
>>> rhp_is_valid(None)
False
>>> rhp_is_valid(123456)
Traceback (most recent call last):
...
TypeError: object of type 'int' has no len()

rhealpixdggs.rhp_wrappers.rhp_to_center_child(rhpindex: str, res: int | None = None, verbose: bool = True, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → str | None

Returns central child of rhpindex at resolution res (immediate central child if res == None).

Returns None if the cell index is invalid.

Returns None if the DGGS has an even number of cells on a side.

EXAMPLES:

>>> rhp_to_center_child('S001450634')
'S0014506344'
>>> rhp_to_center_child('S001450634', res=13)
'S0014506344444'
>>> rhp_to_center_child('S001450634', res=9)
'S001450634'
>>> rhp_to_center_child('S001450634', res=10)
'S0014506344'
>>> rhp_to_center_child('INVALID')

rhealpixdggs.rhp_wrappers.rhp_to_geo(rhpindex: str, geo_json: bool = True, plane: bool = True, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → tuple[float, float] | None

Look up the centroid (in degrees) of the cell identified by rhpindex: the cell’s direct position in the sense of OGC Topic 21 (see Cell.centroid). With plane = True, the default, that is the centre of the planar square, which is also the cell’s nucleus; with plane = False it is the area-weighted centroid on the ellipsoid.

Returns None if the cell index is invalid.

If geojson is requested as the output format:

• Will return a (longitude, latitude) coordinate pair.

if geojson is NOT requested as the output format:

• Will return a (latitude, longitude) coordinate pair in order to be consistent with h3 coordinate ordering.

EXAMPLES:

>>> rhp_to_geo('S001450634', True, False)
(-176.2606635452476, -43.7365450535837)
>>> rhp_to_geo('S001450635', True, False)
(-176.25592420875034, -43.7365450535837)
>>> rhp_to_geo('NotACellId', True, False)

rhealpixdggs.rhp_wrappers.rhp_to_geo_boundary(rhpindex: str, geo_json: bool = True, plane: bool = True, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → tuple[tuple[float, float], ...] | None

Extract the corner coordinates of a cell at a given cell ID and returns them as a tuple of coordinate pairs (in degrees).

Returns None if the cell index is invalid.

If geojson is requested as the output format:

• Will return (longitude, latitude) coordinate pairs.

• Will repeat the first vertex and insert it at the end if geojson is requested as the output format.

If geojson is NOT requested as the output format:

• Will return (latitude, longitude) coordinate pairs in order to be consistent with rHEALPix coordinate ordering.

EXAMPLES:

>>> rhp_to_geo_boundary('S001450634', False, False)
((-43.73395872598705, -176.26086040756147), (-43.73395872598705, -176.25612132062557), (-43.73913136381169, -176.26046658591815), (-43.73913136381169, -176.2652061719943))
>>> rhp_to_geo_boundary('S001450634', True, False)
((-176.26086040756147, -43.73395872598705), (-176.25612132062557, -43.73395872598705), (-176.26046658591815, -43.73913136381169), (-176.2652061719943, -43.73913136381169), (-176.26086040756147, -43.73395872598705))
>>> rhp_to_geo_boundary('INVALID')

rhealpixdggs.rhp_wrappers.rhp_to_parent(rhpindex: str, res: int | None = None, verbose: bool = True, dggs: RHEALPixDGGS = <rhealpixdggs.dggs.RHEALPixDGGS object>) → str | None

Returns parent of rhpindex at resolution res (immediate parent if res == None).

Returns None if the cell index is invalid.

EXAMPLES:

>>> rhp_to_parent('S001450634', res=0)
'S'
>>> rhp_to_parent('S001450634')
'S00145063'
>>> rhp_to_parent('INVALID')