mortie.convert

Address-space conversions between geographic coordinates, packed morton words, UNIQ cell numbers and HEALPix NESTED ids — plus mort2bbox / mort2polygon, which turn a word into a bounding box or a ring. Split out of mortie.tools by domain (issue #159) so the Python surface mirrors the Rust tree (geo2mort.rs, morton.rs, cell_geom.rs); the names stay flat on the package (mortie.geo2mort, mortie.mort2polygon).

Every geographic entry point here takes a keyword-only latitude= argument. Its default, "authalic", maps WGS84 geodetic latitude to authalic latitude on the way into the spherical kernel and back on the way out, so cells are equal-area on the ellipsoid; latitude="geodetic-spherical" is the pre-0.10 escape. The two conventions are non-corresponding partitions — see specification.md §9. The geodetic_to_authalic / authalic_to_geodetic pair below exposes that latitude→latitude mapping on its own.

Address-space conversions between geographic, morton, UNIQ and HEALPix.

The X2Y family: :func:geo2mort / :func:mort2geo and :func:geo2uniq / :func:uniq2geo across the geographic boundary, :func:norm2mort / :func:mort2norm and :func:norm2uniq / :func:unique2parent across the normalized-address boundary, and :func:mort2healpix out to NESTED cell ids. :func:mort2bbox and :func:mort2polygon belong here too: from the caller's side they turn a word into a bounding box or a ring, which is a conversion -- even though their kernels live in src_rust/src/cell_geom.rs rather than in geo2mort.rs / morton.rs with the rest of this module's twins.

Split out of mortie.tools (issue #159) so the Python surface mirrors the Rust tree's own decomposition. The names stay flat on the package (mortie.geo2mort, mortie.mort2polygon): this module is where they live, not how they are spelled.

`geo2mort(lats, lons, order=None, points=None, *, latitude='authalic')`

Compute morton indices from geographic coordinates.

The entire pipeline runs in Rust via the healpix crate — no Python HEALPix backend is needed.

lat/lon inputs are treated as points by default (indeterminate resolution, encoded at max precision), so a bare geo2mort(lats, lons) returns order-29 Kind::Point words. Passing an explicit order asks for an area cell at that resolution instead (points inferred False). The two flags resolve as:

order=None, points=None (bare call) -> order-29 point words;
an explicit order with points unset -> area cell at order;
points=True -> order-29 point words (order-29-only; an explicit order != 29 raises ValueError, matching :meth:MortonIndexArray.from_latlon);
points=False -> area cell at order (order=None -> 29).

Non-finite lat/lon encode to the reserved empty word 0 (base cell 0 is the null sentinel) on both the area and point routes.

Parameters:

Name	Type	Description	Default
`lats`	`array - like`	Latitude(s) in degrees.	required
`lons`	`array - like`	Longitude(s) in degrees.	required
`order`	`int`	HEALPix order (0-29). Defaults to 29. An explicit value implies an area cell unless `points=True` is also given.	`None`
`points`	`bool`	Encode `Kind::Point` (order-29) vs `Kind::Area` words. Defaults to `True` for a bare call and `False` when an `order` is given.	`None`
`latitude`	`str`	Latitude convention of the input (issue #186): `"authalic"` (default; geodetic latitudes are converted so cells are equal-area on the WGS84 ellipsoid) or `"geodetic-spherical"` (legacy: geodetic latitude fed to the spherical kernel as-is). Cell ids under the two conventions are non-corresponding partitions — never mix them.	`'authalic'`

Returns:

Type	Description
`ndarray`	Packed `uint64` morton word(s), same shape family as the input (scalar in -> length-1 ndarray).

Raises:

Type	Description
`ValueError`	If `points=True` is combined with an explicit `order != 29`, or latitude is not a valid convention.

Source code in mortie/convert.py

def geo2mort(lats, lons, order=None, points=None, *, latitude="authalic"):
    """Compute morton indices from geographic coordinates.

    The entire pipeline runs in Rust via the ``healpix`` crate — no
    Python HEALPix backend is needed.

    lat/lon inputs are treated as **points** by default (indeterminate
    resolution, encoded at max precision), so a bare ``geo2mort(lats, lons)``
    returns order-29 ``Kind::Point`` words. Passing an explicit ``order`` asks
    for an **area** cell at that resolution instead (``points`` inferred
    ``False``). The two flags resolve as:

    * ``order=None, points=None`` (bare call) -> order-29 **point** words;
    * an explicit ``order`` with ``points`` unset -> **area** cell at ``order``;
    * ``points=True`` -> order-29 point words (order-29-only; an explicit
      ``order != 29`` raises ``ValueError``, matching
      :meth:`MortonIndexArray.from_latlon`);
    * ``points=False`` -> area cell at ``order`` (``order=None`` -> 29).

    Non-finite ``lat``/``lon`` encode to the reserved empty word ``0`` (base
    cell 0 is the null sentinel) on both the area and point routes.

    Parameters
    ----------
    lats : array-like
        Latitude(s) in degrees.
    lons : array-like
        Longitude(s) in degrees.
    order : int, optional
        HEALPix order (0-29). Defaults to 29. An explicit value implies an area
        cell unless ``points=True`` is also given.
    points : bool, optional
        Encode ``Kind::Point`` (order-29) vs ``Kind::Area`` words. Defaults to
        ``True`` for a bare call and ``False`` when an ``order`` is given.
    latitude : str, optional
        Latitude convention of the input (issue #186): ``"authalic"``
        (default; geodetic latitudes are converted so cells are equal-area on
        the WGS84 ellipsoid) or ``"geodetic-spherical"`` (legacy: geodetic
        latitude fed to the spherical kernel as-is).  Cell ids under the two
        conventions are non-corresponding partitions — never mix them.

    Returns
    -------
    ndarray
        Packed ``uint64`` morton word(s), same shape family as the input
        (scalar in -> length-1 ndarray).

    Raises
    ------
    ValueError
        If ``points=True`` is combined with an explicit ``order != 29``, or
        *latitude* is not a valid convention.
    """
    # Resolve the point/area mode: a bare call encodes points; an explicit order
    # implies an area cell at that resolution unless the caller forces points.
    if points is None:
        points = order is None
    if order is None:
        order = MAX_ORDER
    if points and int(order) != MAX_ORDER:
        raise ValueError(
            "points=True encodes an order-29 point; pass order=29 "
            "(the default) or omit it"
        )
    # Ensure contiguous arrays for Rust FFI
    if not np.isscalar(lats):
        lats = np.ascontiguousarray(lats, dtype=np.float64)
        lons = np.ascontiguousarray(lons, dtype=np.float64)
    result = _rust_geo2mort(lats, lons, int(order), points, latitude)
    # Always return a contiguous uint64 ndarray. The scalar Rust path hands back
    # a Python int (which np would otherwise infer as int64), so coerce to keep
    # the dtype uint64 regardless of scalar-vs-array input or hemisphere.
    return np.ascontiguousarray(np.atleast_1d(result), dtype=np.uint64)

`mort2geo(morton, *, latitude='authalic')`

Convert morton index to lat/lon of pixel center.

This is the inverse of geo2mort, returning the center coordinates of the HEALPix cell identified by the morton index.

Mixed-order arrays are supported (issue #116): elements are grouped by order (:func:orders_of), each group runs the uniform kernel, and the results scatter back to input positions. Point words (spec §4) are order 29 by definition and group with order 29 — a point's location is exactly what mort2geo returns.

Parameters:

Name	Type	Description	Default
`morton`	`int or array - like`	Morton index (mixed orders allowed).	required
`latitude`	`str`	Latitude convention of the returned coordinates (issue #186): `"authalic"` (default) converts the kernel-frame latitude back to WGS84 geodetic; `"geodetic-spherical"` returns the legacy spherical latitude as-is. Pass the same convention the words were encoded under.	`'authalic'`

Returns:

Name	Type	Description
`lat`	`float or array`	Latitude in degrees
`lon`	`float or array`	Longitude in degrees

Source code in mortie/convert.py

def mort2geo(morton, *, latitude="authalic"):
    """Convert morton index to lat/lon of pixel center.

    This is the inverse of geo2mort, returning the center coordinates
    of the HEALPix cell identified by the morton index.

    Mixed-order arrays are supported (issue #116): elements are grouped by
    order (:func:`orders_of`), each group runs the uniform kernel, and the
    results scatter back to input positions. Point words (spec §4) are order
    29 by definition and group with order 29 — a point's location is exactly
    what mort2geo returns.

    Parameters
    ----------
    morton : int or array-like
        Morton index (mixed orders allowed).
    latitude : str, optional
        Latitude convention of the **returned** coordinates (issue #186):
        ``"authalic"`` (default) converts the kernel-frame latitude back to
        WGS84 geodetic; ``"geodetic-spherical"`` returns the legacy spherical
        latitude as-is.  Pass the same convention the words were encoded
        under.

    Returns
    -------
    lat : float or array
        Latitude in degrees
    lon : float or array
        Longitude in degrees
    """
    _check_latitude(latitude)
    # Handle scalar vs array input to match geo2mort behavior
    input_is_scalar = np.isscalar(morton)

    # Group-by-order dispatch for mixed-order input (issue #116).
    if not input_is_scalar:
        words = np.atleast_1d(np.asarray(morton, dtype=np.uint64))
        orders = orders_of(words)
        unique_orders = np.unique(orders)
        if unique_orders.size > 1:
            lat = np.empty(words.size, dtype=np.float64)
            lon = np.empty(words.size, dtype=np.float64)
            for order in unique_orders:
                mask = orders == order
                lat[mask], lon[mask] = mort2geo(words[mask], latitude=latitude)
            return lat, lon

    # Decode morton to normalized address and parent
    normed, parent, order = mort2norm(morton)

    # Convert to UNIQ
    uniq = norm2uniq(normed, parent, order)

    # Convert to lat/lon (uniq2geo decodes the order from the UNIQ value and
    # applies the egress latitude conversion)
    lat, lon = uniq2geo(uniq, latitude=latitude)

    # Return array to match geo2mort behavior
    if input_is_scalar:
        return np.array([lat]), np.array([lon])
    return lat, lon

`mort2bbox(morton, *, latitude='authalic')`

Convert morton index to bounding box of the pixel.

For pixels touching the antimeridian, vertex longitudes at ±180° are normalized to use consistent representation based on hemisphere voting, preventing bbox misinterpretation as spanning the entire globe.

Mixed-order arrays are supported (issue #116): elements are grouped by order (:func:orders_of), each group runs the uniform kernel, and the results scatter back to input positions. Point words (spec §4) are order 29 by definition and group with order 29 — a point yields the bounding box of its containing order-29 cell (the cell that contains the point), which is exactly the bbox of the order-29 area word at the same location. A group of points therefore covers a well-defined area, element by element.

Parameters:

Name	Type	Description	Default
`morton`	`int or array - like`	Morton index (mixed orders allowed).	required
`latitude`	`str`	Latitude convention of the returned box (issue #186): `"authalic"` (default) converts vertex latitudes back to WGS84 geodetic; `"geodetic-spherical"` returns legacy spherical latitudes. Pass the convention the words were encoded under.	`'authalic'`

Returns:

Name	Type	Description
`bbox`	`dict or list of dicts`	Bounding box in format suitable for STAC/CMR: {"west": min_lon, "south": min_lat, "east": max_lon, "north": max_lat}

Source code in mortie/convert.py

def mort2bbox(morton, *, latitude="authalic"):
    """Convert morton index to bounding box of the pixel.

    For pixels touching the antimeridian, vertex longitudes at ±180° are
    normalized to use consistent representation based on hemisphere voting,
    preventing bbox misinterpretation as spanning the entire globe.

    Mixed-order arrays are supported (issue #116): elements are grouped by
    order (:func:`orders_of`), each group runs the uniform kernel, and the
    results scatter back to input positions. Point words (spec §4) are order
    29 by definition and group with order 29 — a point yields the bounding box
    of its containing order-29 cell (the cell that contains the point), which
    is exactly the bbox of the order-29 **area** word at the same location. A
    group of points therefore covers a well-defined area, element by element.

    Parameters
    ----------
    morton : int or array-like
        Morton index (mixed orders allowed).
    latitude : str, optional
        Latitude convention of the **returned** box (issue #186):
        ``"authalic"`` (default) converts vertex latitudes back to WGS84
        geodetic; ``"geodetic-spherical"`` returns legacy spherical
        latitudes.  Pass the convention the words were encoded under.

    Returns
    -------
    bbox : dict or list of dicts
        Bounding box in format suitable for STAC/CMR:
        {"west": min_lon, "south": min_lat, "east": max_lon, "north": max_lat}
    """
    _check_latitude(latitude)
    morton = np.atleast_1d(morton)
    is_scalar = len(morton) == 1

    words = np.asarray(morton, dtype=np.uint64)
    # Group-by-order dispatch for mixed-order input (issue #116).
    orders = orders_of(words)
    unique_orders = np.unique(orders)
    if unique_orders.size > 1:
        bboxes = [None] * words.size
        for order in unique_orders:
            (idx,) = np.nonzero(orders == order)
            group = mort2bbox(words[idx], latitude=latitude)
            if idx.size == 1:
                bboxes[idx[0]] = group  # length-1 call returns the bare dict
            else:
                for i, bbox in zip(idx, group):
                    bboxes[i] = bbox
        return bboxes

    # First get the pixel center
    normed, parent, order = mort2norm(morton)
    uniq = norm2uniq(normed, parent, order)

    nside = 2**order
    nest = uniq - 4 * (nside**2)

    # Get pixel boundaries: (N, 3, 4) — cell in axis 0, xyz in axis 1, the 4
    # corners in axis 2.  A single cell comes back 2-D (3, 4); promote it.
    boundaries = hp.boundaries(order, nest)
    if boundaries.ndim == 2:
        boundaries = boundaries[np.newaxis, ...]
    n = len(morton)

    # One batched vec2ang over every cell's corners (one Rust round-trip instead
    # of one per cell), then reshape to (N, 4).
    verts = np.transpose(boundaries, (0, 2, 1)).reshape(-1, 3)
    theta, phi = hp.vec2ang(verts)
    lats_all = (90 - np.degrees(theta)).reshape(n, 4)
    if latitude == "authalic":  # egress: kernel frame -> geodetic (issue #186)
        lats_all = authalic_to_geodetic(lats_all.ravel()).reshape(n, 4)
    lons_all = np.degrees(phi)
    lons_all = np.where(lons_all > 180, lons_all - 360, lons_all).reshape(n, 4)

    bboxes = []
    for i in range(n):
        lats = lats_all[i]
        lons = lons_all[i]

        # Normalize antimeridian representation
        # Check if bbox touches antimeridian with mixed ±180°
        ANTIMERIDIAN_TOLERANCE = 1e-6
        on_antimeridian = np.abs(np.abs(lons) - 180.0) < ANTIMERIDIAN_TOLERANCE

        if np.any(on_antimeridian) and (np.max(lons) - np.min(lons)) > 180:
            # Count vertices in each hemisphere (excluding those on antimeridian)
            non_antimeridian = ~on_antimeridian
            if np.any(non_antimeridian):
                western_count = np.sum(lons[non_antimeridian] < -0.1)
                eastern_count = np.sum(lons[non_antimeridian] > 0.1)

                # Determine target longitude for antimeridian vertices
                if western_count > eastern_count:
                    target_lon = -180.0
                elif eastern_count > western_count:
                    target_lon = 180.0
                else:
                    # Use median of non-antimeridian lons
                    median_lon = np.median(lons[non_antimeridian])
                    target_lon = -180.0 if median_lon < 0 else 180.0

                # Normalize antimeridian vertices
                lons = lons.copy()
                lons[on_antimeridian] = target_lon

        # Create bounding box
        bbox = {
            "west": float(np.min(lons)),
            "south": float(np.min(lats)),
            "east": float(np.max(lons)),
            "north": float(np.max(lats))
        }
        bboxes.append(bbox)

    if is_scalar:
        return bboxes[0]
    return bboxes

`mort2polygon(morton, step=1, *, latitude='authalic')`

Convert morton index to polygon representation.

Parameters:

Name	Type	Description	Default
`morton`	`int or array - like`	Morton index.	required
`step`	`int`	Points per side for the cell boundary (default 1 = 4 corners). Use step=32 for 128 boundary points that accurately trace curved cell edges, important for polar cells where 4-corner polygons poorly approximate the true HEALPix boundary.	`1`
`latitude`	`str`	Latitude convention of the returned ring (issue #186): `"authalic"` (default) converts vertex latitudes back to WGS84 geodetic; `"geodetic-spherical"` returns legacy spherical latitudes. Pass the convention the words were encoded under.	`'authalic'`

Returns:

Name	Type	Description
`polygon`	`list or list of lists`	Polygon coordinates as [[lat, lon], ...] in standard geographic order. The polygon is closed (first point repeated at end). Note: Returns [lat, lon] pairs, NOT [lon, lat]. This is the standard geographic coordinate order used by most spatial analysis libraries.

Notes

Polygons that touch the antimeridian (±180° longitude) are automatically normalized to use consistent longitude representation (-180 or +180) based on which hemisphere contains the majority of vertices. This prevents spatial libraries from misinterpreting touching polygons as crossing polygons.

Mixed-order arrays are supported (issue #116): elements are grouped by order (:func:orders_of), each group runs the uniform kernel, and the results scatter back to input positions (rings are 4step+1 vertices at every order, so mixed orders do not change the output shape). Point words (spec §4) are order 29 by definition and group with order 29 — a point yields the polygon ring of its containing order-29 cell, exactly the ring of the order-29 area* word at the same location.

Source code in mortie/convert.py

def mort2polygon(morton, step=1, *, latitude="authalic"):
    """Convert morton index to polygon representation.

    Parameters
    ----------
    morton : int or array-like
        Morton index.
    step : int, optional
        Points per side for the cell boundary (default 1 = 4 corners).
        Use step=32 for 128 boundary points that accurately trace
        curved cell edges, important for polar cells where 4-corner
        polygons poorly approximate the true HEALPix boundary.
    latitude : str, optional
        Latitude convention of the **returned** ring (issue #186):
        ``"authalic"`` (default) converts vertex latitudes back to WGS84
        geodetic; ``"geodetic-spherical"`` returns legacy spherical
        latitudes.  Pass the convention the words were encoded under.

    Returns
    -------
    polygon : list or list of lists
        Polygon coordinates as [[lat, lon], ...] in standard geographic order.
        The polygon is closed (first point repeated at end).

        **Note**: Returns [lat, lon] pairs, NOT [lon, lat]. This is the standard
        geographic coordinate order used by most spatial analysis libraries.

    Notes
    -----
    Polygons that touch the antimeridian (±180° longitude) are automatically
    normalized to use consistent longitude representation (-180 or +180) based
    on which hemisphere contains the majority of vertices. This prevents spatial
    libraries from misinterpreting touching polygons as crossing polygons.

    Mixed-order arrays are supported (issue #116): elements are grouped by
    order (:func:`orders_of`), each group runs the uniform kernel, and the
    results scatter back to input positions (rings are 4*step+1 vertices at
    every order, so mixed orders do not change the output shape). Point words
    (spec §4) are order 29 by definition and group with order 29 — a point
    yields the polygon ring of its containing order-29 cell, exactly the ring
    of the order-29 **area** word at the same location.
    """
    _check_latitude(latitude)
    morton = np.atleast_1d(morton)
    is_scalar = len(morton) == 1

    words = np.asarray(morton, dtype=np.uint64)
    # Group-by-order dispatch for mixed-order input (issue #116).
    orders = orders_of(words)
    unique_orders = np.unique(orders)
    if unique_orders.size > 1:
        polygons = [None] * words.size
        for order in unique_orders:
            (idx,) = np.nonzero(orders == order)
            group = mort2polygon(words[idx], step=step, latitude=latitude)
            if idx.size == 1:
                polygons[idx[0]] = group  # length-1 call returns the bare ring
            else:
                for i, polygon in zip(idx, group):
                    polygons[i] = polygon
        return polygons

    # Get pixel information
    normed, parent, order = mort2norm(morton)
    uniq = norm2uniq(normed, parent, order)

    nside = 2**order
    nest = uniq - 4 * (nside**2)

    # Get pixel boundaries: (N, 3, 4*step) — cell in axis 0, xyz in axis 1, the
    # boundary points in axis 2.  A single cell comes back 2-D (3, ncols);
    # promote it.
    boundaries = hp.boundaries(order, nest, step=step)
    if boundaries.ndim == 2:
        boundaries = boundaries[np.newaxis, ...]
    n = len(morton)
    ncols = 4 * step

    # One batched vec2ang over every cell's boundary points (one Rust round-trip
    # instead of one per cell), then reshape to (N, ncols).
    verts = np.transpose(boundaries, (0, 2, 1)).reshape(-1, 3)
    theta, phi = hp.vec2ang(verts)
    lats_all = (90 - np.degrees(theta)).reshape(n, ncols)
    if latitude == "authalic":  # egress: kernel frame -> geodetic (issue #186)
        lats_all = authalic_to_geodetic(lats_all.ravel()).reshape(n, ncols)
    lons_all = np.degrees(phi)
    lons_all = np.where(lons_all > 180, lons_all - 360, lons_all).reshape(n, ncols)

    polygons = []
    for i in range(n):
        lats = lats_all[i]
        lons = lons_all[i]

        # Create polygon as list of [lat, lon] pairs (standard geographic order)
        # Close the polygon by repeating first point
        polygon = [[float(lats[j]), float(lons[j])] for j in range(len(lons))]
        polygon.append(polygon[0])  # Close the polygon

        # Normalize antimeridian representation to prevent misinterpretation
        polygon = _normalize_antimeridian_polygon(polygon)

        polygons.append(polygon)

    if is_scalar:
        return polygons[0]
    return polygons

`mort2healpix(morton)`

Convert morton index to HEALPix cell ID and order.

Parameters:

Name	Type	Description	Default
`morton`	`int or array - like`	Morton index.	required

Returns:

Name	Type	Description
`cell_ids`	`int or ndarray`	HEALPix cell ID(s) in NESTED scheme
`order`	`int`	HEALPix order (resolution level)

Raises:

Type	Description
`ValueError`	If the words are at mixed orders (propagated from :func:`mort2norm`, which enforces the same-order precondition below).

Notes

The function converts morton indices to HEALPix NESTED scheme cell IDs. All input morton indices must be at the same order.

Examples:

>>> import mortie
>>> m = mortie.geo2mort(-80.0, 120.0, order=6)[0]
>>> cell_id, order = mortie.mort2healpix(m)
>>> print(f"HEALPix cell {cell_id} at order {order}")
HEALPix cell 37010 at order 6

Source code in mortie/convert.py

def mort2healpix(morton):
    """Convert morton index to HEALPix cell ID and order.

    Parameters
    ----------
    morton : int or array-like
        Morton index.

    Returns
    -------
    cell_ids : int or ndarray
        HEALPix cell ID(s) in NESTED scheme
    order : int
        HEALPix order (resolution level)

    Raises
    ------
    ValueError
        If the words are at mixed orders (propagated from :func:`mort2norm`,
        which enforces the same-order precondition below).

    Notes
    -----
    The function converts morton indices to HEALPix NESTED scheme cell IDs.
    All input morton indices must be at the same order.

    Examples
    --------
    >>> import mortie
    >>> m = mortie.geo2mort(-80.0, 120.0, order=6)[0]
    >>> cell_id, order = mortie.mort2healpix(m)
    >>> print(f"HEALPix cell {cell_id} at order {order}")
    HEALPix cell 37010 at order 6
    """
    # Check if input is scalar before converting to array
    is_scalar = np.isscalar(morton)
    morton = np.atleast_1d(morton)

    # Get normalized morton and order
    normed, parent, order = mort2norm(morton)

    # Convert to UNIQ indexing
    uniq = norm2uniq(normed, parent, order)

    # Convert UNIQ to HEALPix NESTED cell ID
    # UNIQ = 4 * nside^2 + nest_index
    nside = 2**order
    cell_ids = uniq - 4 * (nside**2)

    # Ensure arrays for consistent handling
    cell_ids = np.atleast_1d(cell_ids).astype(np.int64)
    order = np.atleast_1d(order)

    if is_scalar:
        return int(cell_ids[0]), int(order[0])

    # For array input, return single order if all are the same
    order_val = int(order[0]) if len(np.unique(order)) == 1 else order
    return cell_ids, order_val

`mort2norm(morton)`

Convert morton index back to normalized address and parent cell.

Parameters:

Name	Type	Description	Default
`morton`	`int or array - like`	Packed morton word(s) (`uint64`; base cells 7-11 set bit 63).	required

Returns:

Name	Type	Description
`normed`	`int or array`	Normalized HEALPix address
`parent`	`int or array`	Parent base cell (0-11)
`order`	`int or array`	HEALPix order inferred from morton index

Raises:

Type	Description
`ValueError`	If the words are at mixed orders — the return contract carries a single scalar order, so use :func:`orders_of` for per-element orders.

Notes

Empty input returns two empty int64 arrays and order == 0.

Source code in mortie/convert.py

def mort2norm(morton):
    """Convert morton index back to normalized address and parent cell.

    Parameters
    ----------
    morton : int or array-like
        Packed morton word(s) (``uint64``; base cells 7-11 set bit 63).

    Returns
    -------
    normed : int or array
        Normalized HEALPix address
    parent : int or array
        Parent base cell (0-11)
    order : int or array
        HEALPix order inferred from morton index

    Raises
    ------
    ValueError
        If the words are at mixed orders — the return contract carries a single
        scalar order, so use :func:`orders_of` for per-element orders.

    Notes
    -----
    Empty input returns two empty ``int64`` arrays and ``order == 0``.
    """
    morton = np.atleast_1d(np.asarray(morton, dtype=np.uint64))
    is_scalar = len(morton) == 1

    # Empty input: nothing to decode. Return empty int64 arrays (matching the
    # array-path dtype) and order 0.
    if morton.size == 0:
        empty = np.empty(0, dtype=np.int64)
        return empty, empty.copy(), 0

    # The packed-u64 kernel decodes each word to (nested, depth); the depth is
    # the HEALPix order (no decimal-digit scan). Reject mixed orders: the
    # return contract is a single scalar order (the geo kernels above this
    # dispatch group-by-order and never hit this — issue #116).
    nested, depths = _rust_mort2nested(np.ascontiguousarray(morton))
    if np.any(depths != depths[0]):
        raise ValueError(
            f"Mixed orders in morton array: {sorted(set(int(d) for d in depths))}; "
            "use orders_of for per-element orders"
        )

    order = int(depths[0])
    # nested ids are HEALPix cell ids (<< 2^58 for order <= 29), so int64 is safe
    # arithmetic here and keeps normed/parent signed for downstream callers.
    nested = nested.astype(np.int64)
    nside_sq = np.int64(1) << np.int64(2 * order)
    parent = nested // nside_sq
    normed = nested % nside_sq

    if is_scalar:
        return normed[0], parent[0], order
    return normed, parent, order

`norm2mort(normed, parent, order)`

Convert a normalized HEALPix address + base cell to a packed morton word.

The exact inverse of :func:mort2norm: mort2norm(norm2mort(n, p, o)) returns (n, p, o). Born order-29-native (issue #48) — there is no order cap beyond the kernel's MAX_ORDER of 29. The returned uint64 is the packed decimal_morton word (issue #58; the prefix is base+1, so bit 63 is set — a large unsigned value — for base cells 7-11), not the retired decimal encoding.

Parameters:

Name	Type	Description	Default
`normed`	`int or array`	Normalized HEALPix address (the in-base z-order, `0 <= normed < 4**order`).	required
`parent`	`int or array`	Parent base cell (0-11).	required
`order`	`int`	HEALPix order (0-29).	required

Returns:

Name	Type	Description
`morton`	`uint64 or ndarray`	Packed morton word(s).

Source code in mortie/convert.py

def norm2mort(normed, parent, order):
    """Convert a normalized HEALPix address + base cell to a packed morton word.

    The exact inverse of :func:`mort2norm`: ``mort2norm(norm2mort(n, p, o))``
    returns ``(n, p, o)``. Born order-29-native (issue #48) — there is no order
    cap beyond the kernel's ``MAX_ORDER`` of 29. The returned ``uint64`` is the
    packed ``decimal_morton`` word (issue #58; the prefix is ``base+1``, so bit 63
    is set — a large unsigned value — for base cells 7-11), not the retired
    decimal encoding.

    Parameters
    ----------
    normed : int or array
        Normalized HEALPix address (the in-base z-order, ``0 <= normed < 4**order``).
    parent : int or array
        Parent base cell (0-11).
    order : int
        HEALPix order (0-29).

    Returns
    -------
    morton : uint64 or ndarray
        Packed morton word(s).
    """
    normed = np.atleast_1d(np.asarray(normed, dtype=np.int64))
    parent = np.atleast_1d(np.asarray(parent, dtype=np.int64))
    is_scalar = normed.size == 1 and parent.size == 1
    # nested = parent * nside^2 + normed; pack via the kernel bridge.
    nested = (parent.astype(np.uint64) << np.uint64(2 * order)) | normed.astype(
        np.uint64
    )
    n = max(normed.size, parent.size)
    nested = np.ascontiguousarray(np.broadcast_to(nested, (n,)))
    depths = np.full(nested.size, order, dtype=np.uint8)
    morton = _rust_nested2mort(nested, depths)
    if is_scalar:
        return np.uint64(morton[0])
    return morton

`geodetic_to_authalic(lats)`

Convert WGS84 geodetic latitude(s) to authalic latitude (degrees).

The forward half of the issue #186 convention change: authalic latitude substituted into the spherical HEALPix formulas makes mortie's cells equal-area on the WGS84 ellipsoid by construction. The conversion is a 5-harmonic trigonometric series with coefficients derived from the pinned WGS84 constants (a = 6378137, 1/f = 298.257223563); it is exact to <= 1e-13 rad (~0.6 um on the ground). The equator and poles are fixed points; the divergence peaks in the +/-45-degree band, where the authalic latitude is ~0.12830 degrees (~14.26 km of meridian arc) closer to the equator. Longitude is unaffected by the convention, so there is no lons argument.

Every mortie entry point applies this conversion internally under its default latitude="authalic"; this function is the standalone spelling for callers who need the raw latitude mapping (e.g. to reproduce a binning decision or to label an external dataset).

Parameters:

Name	Type	Description	Default
`lats`	`float or array - like`	Geodetic latitude(s) in degrees.	required

Returns:

Type	Description
`float or ndarray`	Authalic latitude(s) in degrees (scalar in -> scalar out).

`authalic_to_geodetic(lats)`

Convert authalic latitude(s) back to WGS84 geodetic latitude (degrees).

The inverse of :func:geodetic_to_authalic, exact to the same <= 1e-13 rad series bound — see there for the convention background (issue #186).

Parameters:

Name	Type	Description	Default
`lats`	`float or array - like`	Authalic latitude(s) in degrees.	required

Returns:

Type	Description
`float or ndarray`	Geodetic latitude(s) in degrees (scalar in -> scalar out).