Point Binning API¶

This section contains the functionality for mapping point-sampled data — satellite Level-2 swaths, station records, trajectories — onto the standard grid-doctor HEALPix representation. Point data cannot use the ESMF weight path: every granule has unique geometry (no weight reuse) and nearest source-to-destination would smear a narrow swath over the entire globe. Binning assigns each sample to its containing HEALPix cell and reduces per cell instead.

See the point data recipe for an end-to-end workflow and the shared technical decisions document for the design rationale.

Binning¶

`bin_to_healpix(ds, level, *, agg='mean', nest=True, lat_name=None, lon_name=None, source_units='auto', fill_values=None, min_count=1, with_counts=False, dense=True)` ¶

Bin point-sampled data into HEALPix cells.

Every sample is assigned to the HEALPix cell containing its coordinates (perfect-sphere geometry, consistent with all other grid-doctor output) and all samples per cell are reduced with the requested aggregation.

Parameters:

Name	Type	Description	Default
`ds`	`Dataset`	Dataset with per-sample latitude/longitude variables (swath geolocation, station coordinates, …). Variables that span the sample dimensions are binned; variables that share none of the sample dimensions are passed through unchanged; variables that overlap them only partially are skipped with a warning. Non-sample dimensions (`time`, `channel`, …) are preserved as batch dimensions.	required
`level`	`int`	Target HEALPix level. Choose it so the cell spacing (`58.6° / 2**level`) is coarser than the sample spacing — only then is the per-cell mean a faithful stand-in for conservative remapping.	required
`agg`	`BinAgg \| Mapping[str, BinAgg]`	Aggregation per variable: a single method applied to all variables or a mapping `{variable: method}`. Supported methods: `"mean"` (continuous fields; binning analogue of conservative remapping), `"mode"` (categorical fields; analogue of nearest-neighbour), `"min"`, `"max"`, and `"count"`. A per-sample sum is deliberately not offered: it scales with sample density (orbit overlap, across-track pixel count) rather than any physical integral.	`'mean'`
`nest`	`bool`	Use nested HEALPix ordering when True. Nested ordering is required for pyramid coarsening.	`True`
`lat_name`	`str \| None`	Explicit coordinate variable names. When omitted, the standard grid-doctor name lists are searched (`latitude`, `lat`, …).	`None`
`lon_name`	`str \| None`	Explicit coordinate variable names. When omitted, the standard grid-doctor name lists are searched (`latitude`, `lat`, …).	`None`
`source_units`	`SourceUnits`	Angular unit convention of the coordinates.	`'auto'`
`fill_values`	`Mapping[str, float] \| None`	Explicit per-variable fill values, e.g. `{"cloud_type": 255}`. When omitted, declared `_FillValue` / `missing_value` attributes are honoured; floating-point non-finite values are always treated as invalid. `0` is accepted as a fill value.	`None`
`min_count`	`int`	Minimum number of valid samples required for a cell to be valid. Cells with fewer samples are set to NaN (does not apply to `"count"`).	`1`
`with_counts`	`bool`	Add an `<name>_count` companion variable per binned variable holding the number of valid samples per cell. Recommended for published datasets: it makes coverage auditable and enables unbiased merging of overlapping granules downstream.	`False`
`dense`	`bool`	When True (default), return the full `12 * 4*level` cell array with standard grid-doctor coordinates, CRS variable, and attributes — directly consumable by `coarsen_healpix` and `save_pyramid`. When False*, return a compact dataset containing only the touched cells (the `cell` coordinate holds the actual HEALPix indices) — useful as a per-granule intermediate at high levels; convert with `sparse_to_dense` before coarsening or publishing.	`True`

Returns:

Type	Description
`Dataset`	Binned dataset on the HEALPix grid.

Raises:

Type	Description
`ValueError`	On unknown aggregation methods, invalid levels, or when no valid sample coordinates exist.

Examples:

hpx = gd.bin_to_healpix(
    swath,
    level=11,
    agg={"radiance": "mean", "cloud_type": "mode"},
    fill_values={"cloud_type": 255},
    with_counts=True,
)
pyramid = {11: hpx}
for lvl in range(10, -1, -1):
    pyramid[lvl] = gd.coarsen_healpix(pyramid[lvl + 1], lvl)

Source code in grid_doctor/swath/__init__.py

def bin_to_healpix(
    ds: xr.Dataset,
    level: int,
    *,
    agg: BinAgg | Mapping[str, BinAgg] = "mean",
    nest: bool = True,
    lat_name: str | None = None,
    lon_name: str | None = None,
    source_units: SourceUnits = "auto",
    fill_values: Mapping[str, float] | None = None,
    min_count: int = 1,
    with_counts: bool = False,
    dense: bool = True,
) -> xr.Dataset:
    """Bin point-sampled data into HEALPix cells.

    Every sample is assigned to the HEALPix cell containing its
    coordinates (perfect-sphere geometry, consistent with all other
    grid-doctor output) and all samples per cell are reduced with the
    requested aggregation.

    Parameters
    ----------
    ds:
        Dataset with per-sample latitude/longitude variables (swath
        geolocation, station coordinates, …).  Variables that span the
        sample dimensions are binned; variables that share none of the
        sample dimensions are passed through unchanged; variables that
        overlap them only partially are skipped with a warning.
        Non-sample dimensions (``time``, ``channel``, …) are preserved
        as batch dimensions.
    level:
        Target HEALPix level.  Choose it so the cell spacing
        (``58.6° / 2**level``) is coarser than the sample spacing —
        only then is the per-cell mean a faithful stand-in for
        conservative remapping.
    agg:
        Aggregation per variable: a single method applied to all
        variables or a mapping ``{variable: method}``.  Supported
        methods: ``"mean"`` (continuous fields; binning analogue of
        conservative remapping), ``"mode"`` (categorical fields;
        analogue of nearest-neighbour), ``"min"``, ``"max"``, and
        ``"count"``.  A per-sample *sum* is deliberately not offered:
        it scales with sample density (orbit overlap, across-track
        pixel count) rather than any physical integral.
    nest:
        Use nested HEALPix ordering when *True*.  Nested ordering is
        required for pyramid coarsening.
    lat_name, lon_name:
        Explicit coordinate variable names.  When omitted, the standard
        grid-doctor name lists are searched (``latitude``, ``lat``, …).
    source_units:
        Angular unit convention of the coordinates.
    fill_values:
        Explicit per-variable fill values, e.g. ``{"cloud_type": 255}``.
        When omitted, declared ``_FillValue`` / ``missing_value``
        attributes are honoured; floating-point non-finite values are
        always treated as invalid.  ``0`` is accepted as a fill value.
    min_count:
        Minimum number of valid samples required for a cell to be
        valid.  Cells with fewer samples are set to NaN (does not apply
        to ``"count"``).
    with_counts:
        Add an ``<name>_count`` companion variable per binned variable
        holding the number of valid samples per cell.  Recommended for
        published datasets: it makes coverage auditable and enables
        unbiased merging of overlapping granules downstream.
    dense:
        When *True* (default), return the full ``12 * 4**level`` cell
        array with standard grid-doctor coordinates, CRS variable, and
        attributes — directly consumable by
        [`coarsen_healpix`][grid_doctor.helpers.coarsen_healpix] and
        [`save_pyramid`][grid_doctor.helpers.save_pyramid].  When
        *False*, return a compact dataset containing only the touched
        cells (the ``cell`` coordinate holds the actual HEALPix
        indices) — useful as a per-granule intermediate at high levels;
        convert with
        [`sparse_to_dense`][grid_doctor.swath.sparse_to_dense] before
        coarsening or publishing.

    Returns
    -------
    xarray.Dataset
        Binned dataset on the HEALPix grid.

    Raises
    ------
    ValueError
        On unknown aggregation methods, invalid levels, or when no
        valid sample coordinates exist.

    Examples
    --------
    ```python
    hpx = gd.bin_to_healpix(
        swath,
        level=11,
        agg={"radiance": "mean", "cloud_type": "mode"},
        fill_values={"cloud_type": 255},
        with_counts=True,
    )
    pyramid = {11: hpx}
    for lvl in range(10, -1, -1):
        pyramid[lvl] = gd.coarsen_healpix(pyramid[lvl + 1], lvl)
    ```
    """
    if not 0 <= level <= MAX_NESTED_LEVEL:
        raise ValueError(f"level must be within [0, {MAX_NESTED_LEVEL}].")
    if min_count < 1:
        raise ValueError("min_count must be at least 1.")

    lat, lon, sample_dims = resolve_point_coords(
        ds, lat_name=lat_name, lon_name=lon_name, source_units=source_units
    )
    # Fail-fast: validate every aggregation before any data is loaded.
    methods = resolve_methods(ds, sample_dims, agg)

    valid_coords = np.isfinite(lat) & np.isfinite(lon)
    if not valid_coords.any():
        raise ValueError("No samples with valid (finite) coordinates found.")

    module, kwargs = _require_healpix_geo_module(nest)
    cell_ids = np.asarray(
        module.lonlat_to_healpix(
            lon[valid_coords], lat[valid_coords], depth=level, **kwargs
        ),
        dtype=np.int64,
    )

    # Compact indexing: bin into the touched cells only, scatter later.
    unique_cells, group_idx = np.unique(cell_ids, return_inverse=True)
    group_idx = group_idx.astype(np.int64)
    n_cells = int(unique_cells.size)

    binned: dict[str, xr.DataArray] = {}
    counts: dict[str, xr.DataArray] = {}
    for name, da in ds.data_vars.items():
        var_name = str(name)
        overlap = set(sample_dims) & set(map(str, da.dims))
        if not overlap:
            binned[var_name] = da
            continue
        if overlap != set(sample_dims):
            logger.warning(
                "Skipping %r: covers only part of the sample dimensions %s.",
                var_name,
                sample_dims,
            )
            continue

        method = methods[var_name]

        batch_dims = tuple(d for d in map(str, da.dims) if d not in sample_dims)
        arranged = da.transpose(*batch_dims, *sample_dims)
        fill = (fill_values or {}).get(var_name)
        values = masked_float64(arranged, fill_value=fill)
        values = values.reshape(*arranged.shape[: len(batch_dims)], -1)
        values = values[..., valid_coords]

        valid_count = bin_count(group_idx, values, n_cells=n_cells)

        if method == "count":
            result = valid_count
        elif method == "mode":
            result = bin_mode(group_idx, values, n_cells=n_cells)
            result[valid_count < min_count] = np.nan
        else:
            result = bin_simple(
                group_idx, values, n_cells=n_cells, func=NPG_FUNC[method]
            )
            result[valid_count < min_count] = np.nan

        attrs = {
            key: value for key, value in da.attrs.items() if key not in FILL_ATTR_NAMES
        }
        attrs["grid_doctor_method"] = AGG_TO_METHOD[method]
        binned[var_name] = xr.DataArray(result, dims=(*batch_dims, "cell"), attrs=attrs)
        if with_counts and method != "count":
            counts[f"{var_name}_count"] = xr.DataArray(
                valid_count.astype(np.int32),
                dims=(*batch_dims, "cell"),
                attrs={"long_name": f"number of valid samples binned into {var_name}"},
            )

    binned.update(counts)
    result_ds = xr.Dataset(binned, attrs=ds.attrs.copy())

    # Preserve batch coordinates (time, channel, ...).
    keep_coords = {
        str(coord): ds.coords[coord]
        for coord in ds.coords
        if not set(sample_dims) & set(map(str, ds.coords[coord].dims))
    }
    result_ds = result_ds.assign_coords(keep_coords)
    result_ds.attrs["grid_doctor_method"] = _dominant_method(result_ds)

    if dense:
        return _scatter_to_dense(result_ds, unique_cells, level=level, nest=nest)
    return _attach_sparse_coords(result_ds, unique_cells, level=level, nest=nest)

Sparse intermediates¶

At high HEALPix levels a single granule touches only a tiny fraction of the global grid. bin_to_healpix(..., dense=False) returns a compact per-granule dataset that can be accumulated cheaply and expanded once before publishing.

`sparse_to_dense(ds)` ¶

Scatter a compact (sparse) binned dataset onto the full grid.