Two spec mismatches from cellmap-flow: source_voxel_size under relabel, box under resample (0.3.1)

CHANGELOG.md

@@ -1,5 +1,20 @@
 # Changelog
 
+## 0.3.1 (unreleased)
+
+- `source_voxel_size` reports the stored level's voxel size under `on_mismatch="relabel"`,
+  as the spec says; it reported the relabelled size (cellmap-flow's report). Reads,
+  `voxel_size`, `translation` and `corner` are unchanged.
+- `box` of an array opened with a resample target is the extent of the target voxels
+  `shape` counts, as the spec's `[corner, corner + shape * voxel_size)` says and as
+  `Grid.from_meta(meta, voxel_size=t).box` already was; it was the stored data's extent, so
+  `read(box)` could return fewer voxels than `shape` (cellmap-flow's 7 voxels of 8 nm at
+  24 nm: box ended at 56 nm, now 72). No voxel moves.
+- The unit table knows kilometres (`kilometer`, `kilometre`, `km`, 1e3 m), so they
+  convert, abbreviate to `km` and are written as OME's `kilometer` (chunkmirage's ask).
+- README: the placement rule gains the corner formula, a worked example, why the centre
+  is the model, and how the read rule follows from it.
+
 ## 0.3.0 (2026-10-22)
 
 Grids and blocks, built to spec draft 3 (Behaviours 17 to 19). Nothing on the read or

README.md

@@ -106,9 +106,18 @@ s1 = wa.open_pyramid("out.zarr", mode="r+").add_level("s1", factors=2)   # place
 ## The placement rule
 
 A translation is where the centre of voxel 0 is. Index `i` sits at
-`translation + i * voxel_size`; voxel `i` covers half a voxel either side. That is how
-OME-NGFF defines it and how Neuroglancer draws it. Dialects that store the corner are
-moved half a voxel on reading:
+`translation + i * voxel_size` and covers half a voxel either side, so the array's corner
+is `corner = translation - voxel_size / 2`. With 8 nm voxels and the corner at 0, the
+translation is 4: voxel 0 covers 0 to 8 nm around its centre at 4, voxel 1 covers 8 to
+16 nm around 12. That is how OME-NGFF defines it and how Neuroglancer draws it.
+
+The centre is the model because a voxel is a sample taken at a point. Everything that
+compares or resamples grids (picking a level, reading one volume onto another's grid,
+interpolating) works on those points, and the box a voxel covers follows from its point.
+The read rule follows too: a voxel is read iff its centre lies in the box, so boxes that
+adjoin never share a voxel or drop one.
+
+Dialects that store the corner are moved half a voxel on reading:
 
 | Stored | Means | Read as `translation` |
 |---|---|---|

docs/spec.md

@@ -149,7 +149,8 @@ vol.meta           # ArrayMeta in the store's own units
 vol.grid           # ArrayMeta in the requested units, every axis
 vol.shape, vol.dtype, vol.axes, vol.spatial, vol.level_name, vol.levels, vol.remote
                    # shape is the stored shape, or with a resample target the target voxels covering the data
-vol.source_voxel_size, vol.resampled   # the stored level's voxel size; whether reads resample to a target
+vol.source_voxel_size, vol.resampled   # the stored level's voxel size, also under resample and relabel; whether reads
+                                       #   resample to a target
 vol.list_chunks()  # the chunk indices that exist in storage, by listing keys (zarr v2/v3, N5, unsharded precomputed)
 vol.path           # what the caller passed (a group, say); vol.level_path is the array actually opened
 vol.voxel_size     # tuple[float] in `units`, spatial axes, the array's order (zyx for a volume)
@@ -212,7 +213,8 @@ window and a spatial box at once); anything with `.begin` and `.end`, or `.offse
 `shift`, `snap(voxel_size, translation, mode="grow" | "shrink" | "closest")` to a grid's
 voxel corners, `to_index` (the centre rule) and `Box.from_index`; `vol.box` is the array's
 own extent on the spatial axes, `[corner, corner + shape * voxel_size)`, what funlib calls
-the array's roi. So `vol.read(vol.box.intersect(roi).snap(vol.voxel_size, vol.translation))`
+the array's roi; with a resample target, the target voxels `shape` counts, so `read(vol.box)`
+returns `shape` and equals `Grid.from_meta(meta, voxel_size=t).box`. So `vol.read(vol.box.intersect(roi).snap(vol.voxel_size, vol.translation))`
 is the migration of a clipped, snapped Roi read, without funlib.
 
 There is no `offset` attribute anywhere in this API, because that word means the corner in
@@ -393,7 +395,7 @@ this document and that one differ, fix one of them.
 3. **Units are kept as stored and converted only on the axes a caller asks to convert**
    (`ArrayMeta.in_units(space=, time=)`, `in_nm()`, the `units` argument), through one table:
    for length, American and British spellings, `nm`, `um`, `µm` (U+00B5), `μm` (U+03BC),
-   `mm`, `cm`, `m`, `angstrom`/`Å`, `pm`; for time, `s`/`second`, `ms`, `us`/`µs`, `ns`,
+   `mm`, `cm`, `m`, `km`, `angstrom`/`Å`, `pm`; for time, `s`/`second`, `ms`, `us`/`µs`, `ns`,
    `minute`, `hour`, `day`. Any unit string on any axis is accepted and never rejected until
    a conversion is asked for (`"s"`, `"m"`, `"degrees_north"`, `""`). A spatial axis with no
    unit is taken as nanometres, with a **warning once per array** (every repo assumed this
@@ -549,8 +551,8 @@ this document and that one differ, fix one of them.
    reads that level voxel for voxel and reports `voxel_size` as the requested size with
    **voxel `i` staying voxel `i`**: every world quantity, the translation and the corner
    included, scales by the ratio of the sizes (cellmap-flow's rule, which chunkmirage's
-   server relies on), so a box of `m` requested voxels returns `m` stored voxels of the
-   level; `"keep"` reads the chosen level at its own voxel size (a level pick, nothing
+   server relies on), so a box of `m` requested voxels returns `m` stored voxels of the level,
+   and `source_voxel_size` still reports the level's own size; `"keep"` reads the chosen level at its own voxel size (a level pick, nothing
    more). `relabel` and `nearest` are different things and both exist: one is a grid
    policy, the other a resampling kernel. `level` also takes a position, `0` or `"first"`
    for the first listed level. The pyramid an array belongs to is found through

pyproject.toml

@@ -1,6 +1,6 @@
 [project]
 name = "worldarray"
-version = "0.3.0"
+version = "0.3.1"
 description = "Read zarr, N5 and Neuroglancer precomputed volumes in world coordinates through tensorstore: correct voxel placement, pyramid levels, padded region reads, resampling. numpy + tensorstore only."
 readme = "README.md"
 license = { file = "LICENSE" }

src/worldarray/__init__.py

@@ -28,7 +28,7 @@
     parse_chunk_key,
 )
 
-__version__ = "0.3.0"
+__version__ = "0.3.1"
 
 __all__ = [
     "ArrayMeta",

src/worldarray/array.py

@@ -139,6 +139,8 @@ def __init__(
             self._meta = self._rescaled_to_base(base_voxel_size, base_units)
         space, time = _unit_targets(units)
         self._grid = self._meta.in_units(space=space, time=time, strict=strict_units)
+        # the stored level's spatial voxel size, before any relabel renames the grid
+        self._stored_voxel_size = tuple(self._grid.voxel_size[i] for i in self._grid.spatial)
 
         self._target: tuple[float, ...] | None = None
         if voxel_size is not None:
@@ -299,8 +301,9 @@ def voxel_size(self) -> tuple[float, ...]:
 
     @property
     def source_voxel_size(self) -> tuple[float, ...]:
-        """The stored level's voxel size on the spatial axes, in ``units``."""
-        return self._source_voxel_size()
+        """The stored level's voxel size on the spatial axes, in ``units``: what the level
+        is, whether reads resample it or ``on_mismatch="relabel"`` renames its grid."""
+        return self._stored_voxel_size
 
     @property
     def resampled(self) -> bool:
@@ -337,8 +340,22 @@ def corner(self) -> tuple[float, ...]:
     def box(self) -> Box:
         """The array's extent on the spatial axes as a world box, ``[corner, corner + shape *
         voxel_size)``: what funlib calls the array's roi. ``Box.from_index`` gives the same
-        over every axis."""
+        over every axis. With a resample target it is the extent of the target voxels that
+        ``shape`` counts, which may reach past the stored data by less than one target voxel,
+        so ``read(box)`` returns ``shape``."""
         sp = self._grid.spatial
+        if self._target is not None:
+            begin, end = [], []
+            for k, axis in enumerate(sp):
+                first, last = _resample.target_range(
+                    self._grid.shape[axis],
+                    self._grid.voxel_size[axis],
+                    self._target[k],
+                    self._origins()[k],
+                )
+                begin.append(self.anchor[k] + first * self._target[k])
+                end.append(self.anchor[k] + last * self._target[k])
+            return Box(tuple(begin), tuple(end))
         return Box.from_index(
             (0,) * len(sp),
             tuple(self._grid.shape[i] for i in sp),

src/worldarray/resample.py

@@ -153,14 +153,22 @@ def apply_axis(block: np.ndarray, axis: int, plan: AxisPlan) -> np.ndarray:
     raise ValueError(plan.kind)
 
 
+def target_range(n: int, s: float, t: float, origin: float = 0.0) -> tuple[int, int]:
+    """``[first, last)``: the target voxels of size ``t`` that overlap ``n`` source voxels of
+    size ``s``, indexed on a target grid whose voxel 0 corner sits ``origin`` source voxels
+    from the source corner."""
+    ratio = t / s
+    first = int(np.floor(snap_integral(-origin / ratio)))
+    last = int(np.ceil(snap_integral((n - origin) / ratio)))
+    return first, max(last, first)
+
+
 def target_shape(n: int, s: float, t: float, origin: float = 0.0) -> int:
     """How many target voxels of size ``t`` cover ``n`` source voxels of size ``s``, counting
     from a target grid whose voxel 0 corner sits ``origin`` source voxels from the source
     corner: the first and last target voxels that overlap the data, inclusive."""
-    ratio = t / s
-    first = np.floor(snap_integral(-origin / ratio))
-    last = np.ceil(snap_integral((n - origin) / ratio))
-    return int(max(last - first, 0))
+    first, last = target_range(n, s, t, origin)
+    return last - first
 
 
 def resample(

src/worldarray/units.py

@@ -19,6 +19,7 @@
     "millimeter": 1e-3, "millimetre": 1e-3, "mm": 1e-3,
     "centimeter": 1e-2, "centimetre": 1e-2, "cm": 1e-2,
     "meter": 1.0, "metre": 1.0, "m": 1.0,
+    "kilometer": 1e3, "kilometre": 1e3, "km": 1e3,
 }  # fmt: skip
 TIME_S = {
     "nanosecond": 1e-9, "ns": 1e-9,
@@ -81,7 +82,7 @@ def factor(
 # entry through the table above, so "nanometre", "nm" and "Nanometers" agree.
 _SHORT = {
     "picometer": "pm", "angstrom": "Å", "nanometer": "nm", "micrometer": "µm",
-    "millimeter": "mm", "centimeter": "cm", "meter": "m",
+    "millimeter": "mm", "centimeter": "cm", "meter": "m", "kilometer": "km",
     "nanosecond": "ns", "microsecond": "µs", "millisecond": "ms", "second": "s",
     "minute": "min", "hour": "h", "day": "d",
 }  # fmt: skip
@@ -103,6 +104,7 @@ def abbreviation(unit: str) -> str:
 # The American spelling is OME's.
 _OME_NAMES = (
     "picometer", "angstrom", "nanometer", "micrometer", "millimeter", "centimeter", "meter",
+    "kilometer",
     "nanosecond", "microsecond", "millisecond", "second", "minute", "hour", "day",
 )  # fmt: skip
 _OME_BY_SIZE = {_lookup(name): name for name in _OME_NAMES}

tests/test_grid.py

@@ -299,6 +299,7 @@ def test_level_factors_stop_when_a_level_fits():
     "mm": ("millimeter", "millimetre", "mm"),
     "cm": ("centimeter", "centimetre", "cm"),
     "m": ("meter", "metre", "m"),
+    "km": ("kilometer", "kilometre", "km"),
     "ns": ("nanosecond", "ns"),
     "µs": ("microsecond", "us", "µs", "μs"),
     "ms": ("millisecond", "ms"),

tests/test_read.py

@@ -550,6 +550,28 @@ def test_relabel_rescales_everything_by_the_ratio(tmp_path):
         and vol.translation == (20.0, 20.0, 20.0)
         and vol.corner == (16.0, 16.0, 16.0)
     )
+    assert vol.source_voxel_size == (4.0, 4.0, 4.0)  # what the level is (cellmap-flow's report)
+    resampled = wa.open(str(group / "s0"), voxel_size=(8, 8, 8), on_mismatch="resample")
+    assert resampled.source_voxel_size == (4.0, 4.0, 4.0)
+    assert resampled.voxel_size == (8.0, 8.0, 8.0)
+
+
+def test_a_resampled_box_is_the_target_voxels_shape_counts(tmp_path):
+    """cellmap-flow's 7 voxels of 8 nm read at 24 nm: 3 target voxels, box 0 to 72, not 56."""
+    group = tmp_path / "rs.zarr"
+    write_zarr2_group(group, ome_multiscales("zyx", [("s0", (8, 8, 8), (4, 4, 4))]))
+    write_zarr2(group / "s0", pattern((7, 7, 7)))
+    vol = wa.open(str(group / "s0"), voxel_size=(24, 24, 24), on_mismatch="resample")
+    assert vol.shape == (3, 3, 3) and vol.box == wa.Box((0.0,) * 3, (72.0,) * 3)
+    assert vol.read(vol.box).shape == vol.shape
+    g = wa.Grid.from_meta(wa.open(str(group / "s0")).grid, voxel_size=(24, 24, 24))
+    assert g.box == vol.box  # the planned grid and the opened array agree
+    # target corners at 4 + 24 k: the voxels from -20 cover the data from 0
+    anchored = wa.open(
+        str(group / "s0"), voxel_size=(24, 24, 24), on_mismatch="resample", anchor=(4, 4, 4)
+    )
+    assert anchored.shape == (4, 4, 4) and anchored.box == wa.Box((-20.0,) * 3, (76.0,) * 3)
+    assert anchored.read(anchored.box).shape == anchored.shape
 
 
 def test_post_transform_runs_after_resampling_and_before_padding(tmp_path):

tests/test_write.py

@@ -57,6 +57,7 @@ def test_units_are_written_in_ome_names_and_unknown_ones_as_given():
         assert wa.units.ome_name(spelling) == "nanometer"
     assert wa.units.ome_name("um") == wa.units.ome_name("\u00b5m") == "micrometer"
     assert wa.units.ome_name("s") == "second"
+    assert wa.units.ome_name("km") == wa.units.ome_name("Kilometres") == "kilometer"
     with pytest.raises(ValueError, match="unknown unit"):
         wa.units.ome_name("pixel")
     doc = wa.ome_multiscales(("z", "y", "x"), "pixel", [("s0", (1, 1, 1), (0.5,) * 3)])