tidy3a

Tidy-style verbs for labelled multidimensional arrays: point clouds, gridded climate data, images, sensor measurements, and complex-valued arrays. Pipelines of select, filter, mutate, group_by, and summarise work on named axes and labels, in the spirit of tidy3 and dplyr. It needs only NumPy (and labelarray, its array layer, which needs only NumPy too); xarray is optional.

Status: version 0.5.0, on PyPI. What changed in each release is in CHANGELOG.md; the design is in docs/plan.md, and the rules fit on one page in docs/decisions.md.

Installing

Requires Python 3.11 or newer.

pip install tidy3a                  # NumPy and labelarray only
pip install "tidy3a[xarray]"        # with xarray and pandas: xarray inputs and outputs, lazy and file-backed data
pip install "tidy3a[xarray,fast]"   # also Numba, for compiled grouped sums and means

Results are the same with or without fast; it only makes large grouped sums and means faster.

Quick start

Read an array, name its axes, compute a value per point, and keep the points within 50 m:

import numpy as np
from tidy3a import fields, filter, formula, infer, mutate, pull, read_array

points = np.array([[3, 4, 0], [0, 0, 50], [60, 0, 0]], dtype="float32")
cloud = read_array(
    points,
    dims={"point": infer, "feature": 3},
    labels={"feature": ["x", "y", "z"]},
    field_axis="feature",
)
x, y, z, distance = fields("x", "y", "z", "distance")

nearby = (
    cloud
    >> mutate(distance = formula(r"\sqrt{x^2 + y^2 + z^2}"))
    >> filter(distance <= 50)
)
print(nearby >> pull(distance))  # [ 5. 50.]

field_axis="feature" makes the labels on that axis act as columns, so x, y, and z name fields of the array. Verbs chain with >>; pull takes one result out as a NumPy array.

Grouping

Group by a value, by an axis, or by a transform of a value, then summarise each group or compute within it:

import numpy as np
from tidy3a import cells, fields, group_by, infer, mean, mutate, n, pull, read_array, rounded, summarise

points = np.array([
    # x,    y,    z, intensity, ring
    [ 1.0,  2.0, 0.5,  10, 0],
    [ 1.5,  2.5, 0.4,  20, 0],
    [-3.0,  0.5, 0.2,  30, 1],
    [-3.2,  0.1, 0.3,  40, 1],
    [ 0.2,  4.1, 0.9,  50, 2],
], dtype="float32")
cloud = read_array(points, dims={"point": infer, "feature": 5},
                   labels={"feature": ["x", "y", "z", "intensity", "ring"]}, field_axis="feature")
x, y, z, intensity, ring = fields("x", "y", "z", "intensity", "ring")

per_ring = cloud >> group_by(ring) >> summarise(points=n(), brightness=mean(intensity))
print(per_ring >> pull("brightness"))  # [15. 35. 50.]

centred = cloud >> group_by(ring) >> mutate(centred=intensity - mean(intensity))
print(centred >> pull("centred"))  # [-5.  5. -5.  5.  0.]

per_cell = cloud >> group_by(cells(x, 2.0)) >> summarise(points=n())
print(per_cell >> pull("x_cell"))  # ['[-4, -2)' '[0, 2)']

heights = cloud >> group_by(rounded(z, 0.5)) >> summarise(points=n())
print(heights >> pull("points"))  # [1 3 1]

voxels = cloud >> group_by(cells(x, 2.0), rounded(z, 0.5), axis="voxel") >> summarise(points=n())
print(voxels >> pull("x_cell"))  # ['[-4, -2)' '[-4, -2)' '[0, 2)' '[0, 2)']
print(voxels >> pull("z_rounded"))  # [0.  0.5 0.5 1. ]
print(voxels >> pull("points"))  # [1 1 2 1]

Several value keys make one group per combination that occurs, as dplyr, pandas, and Polars do, sorted by the first key and then the next. The groups lie along one new axis, group unless axis= names it, and summarise gives each key's value per group as a variable named after the key. Empty combinations are never made, so voxels cost what is occupied, not the bounding box.

bins(field, edges) groups by explicit edges, or by a duration for times (bins(col("time"), "1 day")). Cell and bin edges fall where the printed settings put them: with a size of 0.1, the value 0.3 is in [0.3, 0.4).

A missing key (NaN, or NaT for times) is one group, left out by default with a note saying how many positions it held; group_by(..., dropna=False) keeps it as the last group:

import numpy as np
from tidy3a import fields, group_by, infer, n, pull, read_array, summarise

values = np.array([[10, 0], [20, 0], [30, 1], [40, np.nan]], dtype="float32")
data = read_array(values, dims={"point": infer, "feature": 2},
                  labels={"feature": ["intensity", "ring"]}, field_axis="feature")
intensity, ring = fields("intensity", "ring")

print(data >> group_by(ring) >> summarise(points=n()) >> pull("points"))  # [2 1]
print(data >> group_by(ring, dropna=False) >> summarise(points=n()) >> pull("ring"))  # [ 0.  1. nan]

Working with xarray

tidy() wraps an xarray Dataset or DataArray, and collect() gives the result back as a Dataset:

import numpy as np
import xarray as xr
from tidy3a import bins, col, group_by, mean, pull, summarise, tidy

times = np.datetime64("2026-01-01T00", "ns") + np.arange(8) * np.timedelta64(6, "h")
air = xr.Dataset(
    {"air": (("time", "lat"), np.arange(16, dtype="float64").reshape(8, 2) + 270)},
    coords={"time": times, "lat": [10.0, 20.0]},
)

per_lat = tidy(air) >> group_by("lat") >> summarise(mean_air=mean(col("air")))
print(per_lat >> pull("mean_air"))  # [277. 278.]

daily = tidy(air) >> group_by(bins(col("time"), "1 day"), "lat") >> summarise(mean_air=mean(col("air")))
result = daily.collect()
print(dict(result.sizes))  # {'time_bin': 2, 'lat': 2}

Labels on a shared axis must match exactly: the same labels in another order are matched by label (and the change is noted), and different labels need a join (left_join, inner_join, full_join). The integer label 2**53 + 1 never matches the float 2**53.

Explaining and reusing pipelines

Each verb records what it did as it runs. explain() lists the source and every step with its effect, from those records, computing nothing; verbs also compose before any data into a pipeline to run, preview, or reuse:

import numpy as np
from tidy3a import fields, filter, group_by, infer, mean, mutate, n, pull, read_array, summarise

points = np.array([[3, 4, 0, 10, 1], [0, 0, 50, 20, 2], [60, 0, 0, 30, 1], [1, 1, 1, 40, 2]], dtype="float32")
cloud = read_array(points, dims={"point": infer, "feature": 5},
                   labels={"feature": ["x", "y", "z", "intensity", "ring"]}, field_axis="feature")
x, y, z, intensity, ring, d = fields("x", "y", "z", "intensity", "ring", "d")

steps = mutate(d=x + y + z) >> filter(d <= 50) >> group_by(ring) >> summarise(points=n(), m=mean(intensity))
per_ring = cloud >> steps
print(per_ring >> pull("m"))  # [10. 30.]
print(len(per_ring.plan()))  # 4

print(per_ring.explain()) shows:

<plan source: NumPy array>
  values (point: 4, feature: 5) float32
1. mutate(d = (x + y) + z)
     + d (point: 4) float32
2. filter(d <= 50)
     point: 4 → 3
3. group_by(ring)
     groups: ring (2)
4. summarise(points = n(), m = mean(intensity))
     - axis point
     - axis feature
     + axis ring: 2
     - values
     - d
     + points (ring: 2) int64
     + m (ring: 2) float32
     field axis: feature → none
     groups: none

steps.explain(cloud, head=1000) runs a pipeline on the first 1,000 positions of each observation axis and explains that, to check a long pipeline on a large array quickly. per_ring.plan() gives the verbs that made a result as a pipeline, so other_cloud >> per_ring.plan() repeats them on other data.

Tables: Polars and tidy3

With Polars installed (pip install "tidy3a[polars]"), as_table() gives data with one observation axis as a Polars DataFrame, one row per position, and tidy(table, axis=...) reads a table back:

import numpy as np
from tidy3a import as_table, col, filter, infer, pull, read_array, tidy

points = np.array([[3, 4, 0], [0, 0, 50], [60, 0, 0]], dtype="float32")
cloud = read_array(points, dims={"point": infer, "feature": 3},
                   labels={"feature": ["x", "y", "z"]}, field_axis="feature")
table = cloud >> filter(col("x") < 10) >> as_table()
print(table.columns)  # ['point', 'x', 'y', 'z']
print(table["point"].to_list())  # [1, 2]

back = tidy(table, axis="point", field_axis="feature", fields=["x", "y", "z"])
print(back >> pull("z"))  # [ 0. 50.]

The axis's labels become a column named after it (original positions after a filter, as pull gives them), and field-axis labels and the other variables along the axis become columns. Data with more than one observation axis needs as_table(long=True), one row per cell. Reading a table, each column becomes a variable along axis=, a column of that name becomes its labels, and fields= packs columns into one array along field_axis=. A LazyFrame is collected, with a note; Polars nulls become NaN, promoting integer and boolean columns to floats, also with a note.

With tidy3 installed, tidy3.tidy(cloud) gives a tidy3 frame (tidy3a defines the __tidy__ hook tidy3 calls), and tidy(frame, axis="point") reads one back.

Backends

New data uses the NumPy backend wherever it holds the data without loss, and the xarray backend for the rest: lazy, chunked, and memory-mapped values, coordinates other than an axis's labels, multi-level labels, time zones, and data with xarray storage encodings, so a netCDF file read and written back keeps them. Results are the same on both. set_backend("xarray") chooses xarray for all new data, as before 0.3.0; docs/backends.md lists the few differences.

Numbers

Float32 sums and means add in float64 and round once. Every grouped one, and every ungrouped one over 8 or more values, is correctly rounded: the float32 nearest the exact result, the same for any memory layout, on either backend, with or without fast, and on lazy, chunked (Dask) data whatever the chunking: each chunk is summed in its own task, the chunk sums combine in any order, and the few results an error bound cannot settle get a second pass for their exact sums. Float64 data follows NumPy. Integer arithmetic follows NumPy too, including unsigned subtraction wrapping around: convert first (nir * 1.0 - red) when a result can be negative.

The API at a glance

An array from the calling code, such as a rotation matrix, can be an xarray DataArray or a one-variable read_array(...) result.

Laziness: read_array(..., lazy=True) memory-maps a file and tidy() keeps a lazy xarray input lazy until collect(), but the verbs themselves run eagerly: arithmetic allocates results and filter evaluates its mask. read_binary's options are in the reader contract.

Notebooks

Install IPython separately to use the notebook extension. In a notebook cell:

%load_ext tidy3a.jupyter
from tidy3a import *

The extension enables bare field names inside verbs, so the examples above can omit fields(...). Multi-line pipes need parentheses. Sharing a notebook with tidy3 requires its ownership-aware transformer; the prerequisite is recorded in the plan.

Examples and development

The point-cloud example runs on real nuScenes Mini data and checks itself against plain NumPy:

python examples/nuscenes_sweep.py /path/to/nuscenes-mini

Three public sample datasets (gridded air temperature, ERA-Interim wind, and stars' Landsat scene) run against plain xarray; see examples/sample_data.py for the files, optional reader dependencies, and checks.

From a checkout, install the development dependencies and run the tests (the suite runs every test on both backends):

python -m pip install -e ".[dev]"
python -m pytest -q

Before pushing, run the CI checks in fresh environments with only declared dependencies (requires uv):

scripts/check.sh

scripts/bench.sh times both backends on the reference workloads in a pinned environment (--fast adds the compiled kernels). The release procedure describes the GitHub Release workflow. tidy3a's arrays, exact labels, grouped reductions, and float32 sums come from labelarray, a library of their own that needs only NumPy and has its own tests; the name tidy3a.labelarray still imports it.

License

MIT