Tile-backed virtual arrays - #80
Merged
Merged
Conversation
Port TileArray and the tile engine registry from the 0.3 line (xdas/virtual) nearly verbatim: lazy positive-step slice folding, manifest-level concatenation, streaming reductions, and per-format load engines. Two 0.2 adaptations: a lazy leading-axis expand_dims supporting the legacy concat-along-a-new-dimension path, and signature-bound streaming dispatch compatible with the DataArray reduction wrappers.
Replace the serialized-dask-graph fallback: both engines now emit tile-backed lazy DataArrays. Silixa gains time-axis push-down through TdmsReader row bounds; miniseed keeps its header semantics, with the stream method and ignore_last_sample travelling in the manifest's engine specification.
A tile-backed variable stores its manifest as a __tiles__ sibling group plus a JSON attribute holding what the arrays cannot (engine specification and dtype); open_dataarray reconstructs the TileArray. Writing dask-backed virtual arrays now emits a FutureWarning; the reader stays.
# Conflicts: # xdas/io/miniseed.py # xdas/io/silixa.py
Tile decoding becomes the load_tile half of the single per-format plugin socket xdas.io.Engine: the base class gains an abstract load_tile(path, selection, **kwargs), the silixa and miniseed tile engine classes dissolve into their io engines, and xdas.tiles.registry shrinks to a get_engine adapter resolving manifest engine names against the io registry (the 0.3 line hosts the same lookup over its own registry, keeping tilearray.py identical in both lines).
asn, febus, terra15, apsensing, prodml and the native xdas engine gain vtype="tiles": open_dataarray then backs the DataArray with a lazy TileArray instead of an HDF5 virtual source, and each engine gains the load_tile decode half (bodies shared verbatim with the 0.3 line, so saved tile views open identically there). The hdf5 vtype stays the default. Febus models a whole file as a single tile: the overlap trimming and 3-D block fusing live in load_tile, instead of one VDS mapping entry per block.
The touched blocks' trimmed windows form one rectangular hyperslab (every block keeps the same post-overlap window), so a single h5py read replaces the per-block loop and concatenation; the partial first and last blocks crop away as plain numpy slices.
The get_engine indirection existed to keep tilearray.py line-identical with 0.3; with the mirror abandoned it lost its purpose. The friendly unknown-name error moves into Engine.__class_getitem__ (where every lookup benefits) and tilearray resolves Engine[name] directly. Also reword the 0.2.9 release notes to describe deltas from 0.2.8 only, not internal churn of the 0.2.9 development.
Reads are lazy and stored views outlive the session, so a relative source path would resolve against whatever the working directory is at read time. Construction is the only moment the relative path is still trustworthy (the scan just used it), so the constructor absolutizes every entry of `paths`; manifests reloaded from disk are trusted as stored.
The version was duplicated in pyproject.toml and xdas/__init__.py, and docs/conf.py carried a third copy that had already drifted to 0.2.7. Declare it dynamic and let setuptools read xdas.__version__, which is now the only place to edit; conf.py derives its release from it too.
Set the version to 0.2.9.dev0, the canonical PEP 440 spelling: setuptools normalises devN to .devN, so the shorter 0.2.9dev0 would leave __version__ and the distribution metadata spelled differently. test_version only accepted plain digits between dots, so extend it to the pre/post/dev markers.
# Conflicts: # docs/release-notes.md # tests/test_xdas.py
The 100 000 file ceiling in `open_mfdataarray` was introduced together with the HDF5 virtual layout linking loop: building that mapping costs one libhdf5 call per source file, so both time and memory grow with the file count. The tiles vtype builds no such mapping, its manifest is a plain array write, so the ceiling does not apply to it. Resolve the engine's effective vtype up front and pick the limit from it. The message now names the real constraint, the per-file data arrays the scan holds until they are combined, and points at the way out. Document the two backends side by side in the virtual datasets guide so the trade-off is written down rather than folklore.
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## dev #80 +/- ##
===========================================
Coverage 100.00% 100.00%
===========================================
Files 44 47 +3
Lines 4772 6017 +1245
Branches 748 1016 +268
===========================================
+ Hits 4772 6017 +1245 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
A tile-backed data array stores its manifest in a `__tiles__` sibling group of its variables. `_get_depth` counted that group, so every such data array looked one level deeper than it is and the collection reader took it for a nested collection, then failed on the datasets it found where it expected a group. Any collection holding more than one tile-backed array was therefore impossible to reopen. Skip the manifest group when measuring depth.
Benchmarking the two backends over a multi-million file archive contradicted what the guide claimed. Resolving a region in the HDF5 C library does avoid all per-file Python, but its cost grows with the number of mappings the dataset holds rather than with the size of the request, so the same read gets slower as the archive grows. A tile manifest is searched, so its cost follows how many tiles the read touches and not how many the manifest contains. HDF5 is therefore the quicker reader only while the file count stays modest. Also note that the top of a virtual dataset pyramid opens fast because it defers the work, and charges the first read of each region for it. Fix the tiles API page while here: its `load_tile` entry resolved against `xdas.tiles` and so looked for `xdas.tiles.xdas.io`, which failed the whole documentation build.
A Febus file stores a stack of overlapping blocks rather than one contiguous array. The hdf5 backing has to describe every block with its own mapping, so a Febus manifest grew with the block count on top of the file count; a tile array describes the whole file as one tile and keeps the overlap trimming in the reader. Listing tiles first makes it the default for that engine. The documentation had not kept up. No engine emits dask graphs any more, yet the format table still credited Silixa and MiniSEED to dask and the guide presented it as a live backend; it is now marked deprecated, with the tables stating what each engine supports and which backing it picks by default. A test pins those defaults so the tables cannot drift.
The two module constants and the _effective_vtype helper collapse into a single MAX_OPEN_FILES dict (with a None fallback entry) and a _check_file_count helper that resolves the vtype and enforces the limit.
The open functions (open, open_dataarray, open_mfdataarray, open_mfdatatree) now declare engine, vtype and ctype explicitly, and engine accepts a configured Engine instance as well as a name. Remaining keyword arguments (**engine_kwargs) are forwarded to the engine constructor only, so format-specific parameters keep working next to the engine name but misspelled ones raise instead of being silently swallowed. Per-engine parameters (febus overlaps/offset, miniseed ignore_last_sample, prodml swapped_dims, terra15 tz, native group) move to the engine constructors, validated before any file is scanned; the callable-engine escape hatch is removed in favor of subclassing Engine. open_mfdataarray resolves the engine once up front: the file-count limit reads the resolved vtype directly and the per-file opens receive the configured instance. Along the way this fixes miniseed silently ignoring its ctype argument and RealTimeLoader defaulting to an engine name that only a doctest side effect ever registered.
The constructor now wraps a manifest dataset directly — hand-built or reopened from a stored file — while the new from_tiles classmethod builds the manifest from per-tile descriptions at scan time; the from_dataset/_setup indirection is gone. The starts kwarg goes with it: trimmed and decimated geometry is view state, arising from slicing or entering through a stored manifest.
The dtype of a tile array is not a decode target: engines decode into the element type of their sources, and the array records that type at scan time so the lazy array answers dtype without touching files. Make that role explicit: from_tiles takes (paths, sizes, dtype, engine) in the numpy order, reads verify each decoded part against the recorded dtype instead of silently casting, and the stored sidecar spec drops its dtype copy in favor of the placeholder variable's own element type. Casting stays an explicit astype step outside the tiles module.
The old repr led with the shape and the dtype, both of which the labeled array already prints one line above, and reported the engine as a quoted kwarg. Report instead what only the tiling knows: the volume the array stands for and the number of tiles, keyed by engine. Bytes render in decimal units, as xarray renders its Size: header, so the two lines agree rather than differing by the 1024/1000 base. The inline form drops the size and the dtype, which an inline row already carries, as dask's does.
Variable-length HDF5 strings pay per-object global-heap overhead on disk and decode eagerly and slowly when the manifest is reopened (xarray materializes vlen columns at open regardless of decode flags). Writing string columns as netCDF-standard fixed-width char arrays (encoding dtype S1) instead makes a 1M-tile manifest open 3x faster (629 -> 206 ms), halves its memory (289 -> 144 MB) and shrinks it 36% on disk (96 -> 62 MB), with identical values on reopen. Legacy vlen manifests still open as before.
The char-array disk encoding fixed the stored form but left the in-memory manifest an object array of str: 8-byte pointers plus a heap-allocated str per tile, python-level comparisons, and a str-to-char re-encoding step at save plus a char-to-str decode at open. Making fixed-width bytes (S kind, filesystem encoding) the canonical dtype of every string variable — normalized once in the TileArray constructor, so scan-time, hand-built and legacy stored manifests all converge — removes the save-site special case entirely (S arrays land on disk as netCDF char arrays natively) and reopens without any decode. Bytes reach str only at the engine boundary, one os.fsdecode per tile read. On a 1M-tile manifest: memory 80 -> 31 MB (-61%), save 0.44 -> 0.11 s (4.0x), open 206 -> 56 ms (3.7x, on top of the 3x the char encoding already bought), manifest fusion 46 -> 37 ms, identical 31 MB file. Scan-time build pays one os.fsencode per path (0.52 -> 0.78 s per million tiles), noise against real scan I/O.
Every string operation on the manifest went through np.frompyfunc: an object array boxing each element, one python call per tile, then a re-wrap back to fixed-width bytes. Two of the five were pure numpy work — stripping the root prefix and joining it back — and the third, rebasing paths on a new root, only looked per-tile: concat always rebases onto an ancestor directory, so where the old root sits under the new one is one constant prefix shared by every tile. Those become np.strings.slice/add ufuncs. The two that genuinely need python, os.path.abspath at scan time and encoding object arrays that may mix str and bytes, become plain comprehensions writing S arrays directly, which also drops an astype at the call site. The np.strings ufuncs size their output from the input widths, so the results are trimmed back to their longest element — otherwise the root split would stop paying for itself. Requires numpy >= 2.1 for np.strings.slice.
The scan ceiling and the streaming batch size were two numbers for one quantity: how many per-file scan products may be held at once. Fold them into MAX_OPEN_FILES and let CONSOLIDATING_VTYPES say which vtypes can shrink a drained batch, instead of naming tiles in the check. Since the batch equals the ceiling, anything that opened in one call before still takes the single-batch path, and a vtype that cannot consolidate never reaches the streaming path at all. open_mfdatacollection had the same ceiling as a bare literal; it now shares the constant.
xdas/virtual.py becomes xdas/virtual/hdf5.py (the HDF5 virtual dataset backend) and xdas/tiles.py becomes xdas/virtual/tiles.py, under one xdas.virtual package that re-exports everything both modules exposed. Tests re-mirror to tests/virtual/. No behavior change.
VirtualBackend is a marker base whose whole job is naming: backends register by passing vtype= in the class definition and are retrieved with VirtualBackend[vtype], the same registry fashion as Engine (by name) and Coordinate (by ctype). The duck-array contract stays informal: the backends share no implementation. The consolidates flag declares whether concatenation fuses scan products into one compact object.
The consolidating vtypes are no longer a hardcoded set in the routines: the multi-file scan asks the registered backend's consolidates flag, and the save-side virtual detection covers any registered backend instead of enumerating the classes.
Any non-None vtype must name a registered VirtualBackend before the per-engine support check, so a typo fails fast with the registered list — even through AutoEngine, which has no supported list of its own.
'from xdas import tiles' kept resolving through the editable install's finder to the pre-move main checkout, so the suite stayed green until the move was merged there.
Both backends gain a from_variable classmethod that exposes one stored HDF5 variable as their lazy array (a virtual source; a single tile decoded by the generic 'xdas' engine), so the native format's open path becomes VirtualBackend[vtype].from_variable(variable) instead of branching on the vtype literal.
create_variable/finalize_save form the persistence pair of the VirtualBackend contract: the first writes the variable inside the open h5netcdf handle (an HDF5 virtual dataset; a placeholder carrying the engine specification), the second appends what outlives it once the handle closes (nothing; the tile manifest as a sibling group). The native format's save path loses its per-backend branches, and the TILES_GROUP convention moves to the tiles module with the code that writes it.
VirtualBackend now states every member a backend implements — shape, dtype, __getitem__, __array__, from_variable, and the persistence pair — and hosts the properties derived from them (ndim, size, nbytes, empty), shared by both backends instead of living on the hdf5 base. TileArray exposes shape and dtype as read-only properties over the values its constructor computes.
With the whole contract declared on the base and every member overridden at class level by both backends, abstractness is now enforceable: the stubs become docstring-only abstract methods, as Coordinate already does. VirtualArray stays the (now formally) abstract hdf5 family base; its concrete forms and TileArray instantiate unchanged.
The block-crossing test built its reference with the engine default, which f08de4f flipped to tiles: the test compared the tiles path against itself, and the legacy per-block loop lost its coverage. Also flatten a double negative in the open ceiling check.
np.strings.slice (used by the tiles manifest) only exists in numpy>=2.3, which itself requires Python>=3.11. Python 3.10 reaches EOL in October 2026.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Adds
xdas.virtual.tiles, a lazy tile-backed virtual array, and wires it into the IO engines as a newtilesvtype — replacing the serialized-dask-graph fallback used by the formats that HDF5 virtual datasets cannot serve, and lifting the multi-file ceiling for the formats that adopt it.Why
HDF5 virtual datasets can only map sources that are themselves HDF5. Silixa TDMS and MiniSEED therefore fell back to a serialized dask graph: opaque to inspect, and written into files that then needed dask internals to be reopened. There was also no lazy multi-file array the engines could share — each format solved laziness its own way — and the one-mapping-per-source design put a hard 100 000-file ceiling on
open_mfdataarray.What changes
The tiles backend
TileArray. A dense rectilinear grid of file-backed tiles that behaves as one numpy-like lazy array, described by a plain manifest dataset (per-axis tile sizes, source origins and signed strides, per-tile paths and engine parameters). Manifest strings are held as fixed-width bytes and the common source directory is split off into a singlerootvalue, so manifests stay compact in memory and land on disk as plain netCDF char arrays — and a stored view relocates by editing that one root path.np.newaxisfold into the geometry and an axis map; the numpy manipulation routines whose effect is a grid rewrite dispatch lazily too — thetranspose,flip,split,stackandatleastfamilies,expand_dims,squeeze,roll,tile,delete,append/insert— and whole-array reductions (sum,mean,min,max, their nan variants,any/all) stream one tile row at a time with bounded memory. Anything the grid cannot express falls back to a bounded read of the selection.load_tile, the decode half ofio.Engine. A static method receiving a path, one source-local — possibly strided —sliceper source axis, and the manifest's engine specification as keyword arguments; it returns exactly the selected sub-box. It is deliberately state-free, so a stored manifest decodes identically everywhere.VirtualBackend, the abstract base thehdf5andtilesbackends now register on (VirtualBackend[vtype]), declares the whole backend contract — lazy array protocol, open-time wrapping, and the stored form (create_variable/finalize_save) — so the open and save paths need no per-backend branch.Engines
vtype="tiles"on every HDF5 engine — asn, terra15, apsensing, prodml and the native xdas format keephdf5as their default. Febus now defaults to tiles: a Febus file holds a stack of overlapping blocks, and the hdf5 backing needed one mapping per block where a tile array needs one per file; the touched blocks' trimmed windows are read as a single hyperslab. Silixa TDMS and MiniSEED always emit tile-backed arrays now (with time-axis push-down for Silixa).engine,vtypeandctype, andenginealso accepts a configuredxdas.io.Engineinstance. Format-specific parameters (overlaps/offsetfor febus,ignore_last_samplefor miniseed,swapped_dimsfor prodml,tzfor terra15,groupfor the native format) are engine constructor parameters, validated up front.__tiles__sibling group. Saved tile views are directly readable by the 0.3 line.Scaling to large archives
open_mfdataarrayfuses scan results every 100 000 files into compact runs instead of holding one data array per file, so memory is bounded by the batch, not the archive. Results are accumulated without coordinate simplification (lossless in any arrival order) and sorted once at the end — the same result as before, whatever the file naming; anything that opened in one call before still takes the single-batch path unchanged.xdas.sortby. Sorts a tile- or stack-backed data array along a dimension by coordinate value, lazily: blocks are permuted through the manifest, tie points gathered blockwise, no data read. This is how the streamed combine orders shuffled archives, exposed for standalone use.simplifyruns in linear time whatever the number of gaps: the reduce stage is now a one-pass greedy sleeve instead of Douglas-Peucker, which degenerated quadratically on gap-rich coordinates (a 100 000-file gappy archive simplified in minutes; now milliseconds). The deviation guarantee is unchanged — dropped tie points stay withintolerance, surviving values never move.Backward compatibility
np.stringsroutines from numpy 2.3, whose wheels requires Python 3.11+.hdf5remains the default vtype everywhere but Febus, Silixa and MiniSEED, so existing calls return what they did — Febus users can pinvtype="hdf5"to keep the old backing.enginenow raises aTypeError(subclassxdas.io.Engineinstead — see the data-formats guide, which gains a worked custom-engine example). Misspelled engine keywords now raise instead of being silently ignored, and combiningvtype/ctype/engine keywords with an already configured engine instance raises aValueError.FutureWarning; files already written that way still open.ctypeargument is now honored (it previously routed to an unused attribute and always built interpolated time coordinates).New public API is documented in
docs/api/tiles.mdand the reworkeddocs/user-guide/io/pages; the suite grows by ~700 tests (100 % coverage held).