@@ -1113,3 +1113,53 @@ buffer authority and reports native cost; pass `buffer_bytes=` to select a dense
 plan. `plan_concurrency` reports each candidate's strategy (or `per-tile` when
 stream subdivision is needed). Snapshot and remote session constructors do not
 currently expose this opt-in.
+
+## Parallel gathering with bounded intermediates
+
+For independently decodable Zarr shards, a local Zarrista session can opt into
+parallel region reads without allocating a dense bounding array:
+
+```python
+session = LocalReadSession(
+    "zarrista", read_strategy="parallel", gather_workers=8, buffer_bytes=8 << 20,
+)
+reader = session.open("velocity.zarr", "")
+request = reader.plan(selection)
+print(request.strategy, request.gather_workers, request.intermediate_bytes_upper_bound)
+preview = reader.read(selection)
+```
+
+The planner considers up to `gather_workers` workers (1–256), bounded further by the
+number of selected shards and the available buffer. It charges the final output,
+worker-owned selected regions, and one decoded chunk per active worker. It lowers
+the worker count before choosing ordinary native selection when fewer than two
+workers fit. Contiguous requests and requests within one region use native selection.
+Streaming, materialization and prefetch re-plan with their actual remaining buffers.
+This is an explicit policy; small requests can get slower from scheduling overhead.
+
+Selected shards are generated lazily. Each worker gathers only touched subchunks
+inside its shard; results are copied into disjoint locations in the final compact
+array. The bounded executor retains at most one region per scheduled slot, including
+the region being consumed. Every call joins its pool, including on failure, so no
+background worker outlives its read reservation. It cannot interrupt an already
+running native call. Native codec threads, shard indexes, caches and scratch remain
+outside the explicit array-buffer allowance. Outer read workers multiply inner
+worker pools; benchmark the whole concurrency configuration rather than assuming
+more threads improve throughput.
+
+`ParallelGatherBackend(base, factory, GatherWindow(region_shape, workers))` is the
+generic extension point. The factory must supply independent mutable reader state
+with the same metadata as `base`. A read-only native array handle may be shared.
+Worker results must not retain arrays larger than their selected region or one
+decoded chunk; workers must not reserve from the parent reader's shared `ReadBudget`.
+The base must guarantee a one-chunk decode window for strided reads; regions must
+align to whole chunks on every axis. This capability is not inferred for arbitrary
+engines. `LocalReadSession` currently exposes it only for Zarrista partial-shard
+reads; other engines and outer-codec shards are rejected for this strategy.
+
+`ReadStats.read_strategies` counts completed tile executions by strategy, and
+`gather_workers_peak` reports the largest admitted inner worker limit. This is
+not a measurement of active threads; a cache hit may launch none.
+`probe_reads` includes these in every trial, so the requested maximum can be
+compared with actual budget-reduced counts. Native reads report one gather worker;
+that does not count the engine's own internal threads.
@@ -1,6 +1,6 @@
 [project]
 name = "arrayroute"
-version = "0.1.11"
+version = "0.1.12"
 description = "Measured, bounded selections over scientific chunked arrays"
 requires-python = ">=3.11"
 license = "Apache-2.0"
 
@@ -1,9 +1,10 @@
 # SPDX-License-Identifier: Apache-2.0
 """Selection planning independent of a storage format, scheduler or cloud provider."""
 
-from .plan import ArraySpec, DecodeWindow, ReadPlan, Tile, plan, tiles
+from .plan import ArraySpec, DecodeWindow, GatherWindow, ReadPlan, Tile, plan, tiles
 from .reader import ArrayReader, ReadStats
 
 __all__ = [
-    "ArrayReader", "ArraySpec", "DecodeWindow", "ReadPlan", "ReadStats", "Tile", "plan", "tiles",
+    "ArrayReader", "ArraySpec", "DecodeWindow", "GatherWindow",
+    "ReadPlan", "ReadStats", "Tile", "plan", "tiles",
 ]
@@ -0,0 +1,72 @@
+# SPDX-License-Identifier: Apache-2.0
+"""Bounded region gathering over independent, single-chunk decoding readers."""
+from __future__ import annotations
+
+from collections.abc import Callable
+from dataclasses import replace
+
+import numpy as np
+
+from .backends import _grid_regions
+from .execution import ordered_map
+from .plan import GatherWindow, ReadPlan, Selection
+from .reader import Backend
+
+
+class ParallelGatherBackend:
+    """Read disjoint chunk-aligned regions with bounded, worker-owned intermediates.
+
+    Each factory result must have the same metadata and single-chunk decode guarantee
+    as base. It must own its mutable state; read-only native array handles may be shared.
+    Results must not retain arrays larger than their selected region or one decoded
+    chunk. Workers must not reserve from the parent reader's shared ReadBudget.
+    Only planned strided calls parallelize. No thread pool survives a call, including
+    failures. Native reads are joined before the output reservation is released.
+    """
+
+    def __init__(self, base: Backend, factory: Callable[[], Backend], window: GatherWindow):
+        if base.spec.gather_window is not None:
+            raise ValueError("Nested gather backends are not supported")
+        self.base, self.factory, self.window = base, factory, window
+        self.dtype = base.dtype
+        self.spec = replace(base.spec, gather_window=window)
+
+    def read(self, selection: Selection) -> np.ndarray:
+        return self.base.read(selection)
+
+    def read_plan(self, request: ReadPlan) -> np.ndarray:
+        if request.strategy == 'native':
+            return self.base.read(request.selection)
+        if (request.strategy != 'parallel-gather'
+                or not 2 <= request.gather_workers <= self.window.workers):
+            raise ValueError('Parallel gather requires its planned worker count')
+        output = np.empty(request.output_shape, dtype=self.dtype)
+
+        def make_worker():
+            backend = self.factory()
+            if backend.spec != self.base.spec or backend.dtype != self.dtype:
+                raise ValueError('Gather factory changed array metadata')
+
+            def read(region):
+                return region, backend.read(region)
+
+            return read
+
+        regions = _grid_regions(request.selection, self.window.region_shape)
+        with ordered_map(make_worker, regions, workers=request.gather_workers,
+                         max_pending=request.gather_workers) as rows:
+            for region, value in rows:
+                targets = []
+                shape = []
+                for part, original in zip(region, request.selection, strict=True):
+                    if isinstance(original, slice):
+                        assert isinstance(part, slice)
+                        count = len(range(part.start, part.stop, part.step))
+                        begin = (part.start - original.start) // original.step
+                        targets.append(slice(begin, begin + count))
+                        shape.append(count)
+                if value.shape != tuple(shape) or value.dtype != self.dtype:
+                    raise ValueError('Gather worker returned an unexpected shape or dtype')
+                output[tuple(targets)] = value
+                del value
+        return output
@@ -1,5 +1,5 @@
 # SPDX-License-Identifier: Apache-2.0
-"""O(rank) selection costs; lazy subdivision without a global chunk inventory."""
+"""Bounded selection costs; lazy subdivision without a global chunk inventory."""
 
 from __future__ import annotations
 
@@ -34,13 +34,30 @@ def __post_init__(self) -> None:
             raise TypeError("strided_only must be a boolean")
 
 
+@dataclass(frozen=True)
+class GatherWindow:
+    """Independent grid regions with one live decoded chunk per worker."""
+
+    region_shape: tuple[int, ...]
+    workers: int
+
+    def __post_init__(self) -> None:
+        if (isinstance(self.workers, bool) or not isinstance(self.workers, int)
+                or not 1 <= self.workers <= 256):
+            raise ValueError("Gather workers must be an integer from 1 to 256")
+        if (not self.region_shape or any(isinstance(v, bool) or not isinstance(v, int)
+                                        or v < 1 for v in self.region_shape)):
+            raise ValueError("Gather regions require positive integer dimensions")
+
+
 @dataclass(frozen=True)
 class ArraySpec:
     shape: tuple[int, ...]
     chunks: tuple[int, ...]
     itemsize: int
     decode_window: DecodeWindow | None = None
     dense_gather: bool = False
+    gather_window: GatherWindow | None = None
 
     def __post_init__(self) -> None:
         if not isinstance(self.dense_gather, bool):
@@ -53,6 +70,17 @@ def __post_init__(self) -> None:
             if any(isinstance(v, bool) or not isinstance(v, int) or v < minimum for v in values):
                 raise ValueError("shape, chunks and itemsize require valid integer dimensions")
 
+        gather = self.gather_window
+        if gather is not None:
+            if not isinstance(gather, GatherWindow):
+                raise TypeError("gather_window must be a GatherWindow or None")
+            if (self.dense_gather or self.decode_window is None
+                    or self.decode_window.max_chunks != 1):
+                raise ValueError("Parallel gathering requires a single-chunk native decode window")
+            if (len(gather.region_shape) != len(self.chunks)
+                    or any(r % c for r, c in zip(gather.region_shape, self.chunks, strict=True))):
+                raise ValueError("Gather regions must be a whole number of chunks on each axis")
+
     @property
     def logical_bytes(self) -> int:
         return math.prod(self.shape) * self.itemsize
@@ -67,8 +95,9 @@ class ReadPlan:
     decoded_bytes_upper_bound: int
     peak_decoded_bytes_upper_bound: int
     intermediate_bytes_upper_bound: int = 0
-    strategy: Literal["native", "dense-gather"] = "native"
+    strategy: Literal["native", "dense-gather", "parallel-gather"] = "native"
     read_selection: Selection | None = None
+    gather_workers: int = 1
 
     @property
     def decode_amplification(self) -> float:
@@ -129,9 +158,29 @@ def plan(spec: ArraySpec, selection: Selection, *, buffer_bytes: int | None = No
         touched * math.prod(spec.chunks) * spec.itemsize,
         live_chunks * math.prod(spec.chunks) * spec.itemsize,
     )
-    if (not spec.dense_gather or buffer_bytes is None or not request.output_bytes
+    if (buffer_bytes is None or not request.output_bytes
             or not any(isinstance(s, slice) and s.step > 1 for s in normalized)):
         return request
+    gather = spec.gather_window
+    if gather is not None:
+        regions = plan(ArraySpec(spec.shape, gather.region_shape, spec.itemsize),
+                       request.selection).chunks_touched
+        region_output = spec.itemsize * math.prod(
+            min(_length(s), (r + s.step - 1) // s.step)
+            for s, r in zip(request.selection, gather.region_shape, strict=True)
+            if isinstance(s, slice)
+        )
+        for workers in range(min(gather.workers, regions), 1, -1):
+            candidate = replace(
+                request, strategy="parallel-gather", gather_workers=workers,
+                peak_decoded_bytes_upper_bound=workers * request.peak_decoded_bytes_upper_bound,
+                intermediate_bytes_upper_bound=min(request.output_bytes, workers * region_output),
+            )
+            if candidate.buffer_bytes_upper_bound <= buffer_bytes:
+                return candidate
+        return request
+    if not spec.dense_gather:
+        return request
     bounds = tuple(
         s if isinstance(s, int) else slice(s.start, s.start + (_length(s) - 1) * s.step + 1, 1)
         for s in normalized
 
@@ -8,6 +8,7 @@
 import random
 import statistics
 import time
+from collections import Counter
 from collections.abc import Callable, Sequence
 from typing import cast
 
@@ -131,7 +132,9 @@ def read(index: int) -> dict:
                                        buffersize=max(1, 32768 // dtype.itemsize)):
                     digest.update(cast(np.ndarray, block).tobytes())
                 return {'sha256': digest.hexdigest(), 'bytes': value.nbytes,
-                        'requests': stats.requests, 'read_seconds': stats.seconds}
+                        'requests': stats.requests, 'read_seconds': stats.seconds,
+                        'read_strategies': stats.read_strategies,
+                        'gather_workers_peak': stats.gather_workers_peak}
 
             return read
 
@@ -152,11 +155,16 @@ def read(index: int) -> dict:
             if any((r['sha256'], r['bytes']) != (e['sha256'], e['bytes'])
                    for r, e in zip(rows, reference, strict=True)):
                 raise ValueError("Pilot values changed from the sequential reference")
+            strategies = Counter()
+            for row in rows:
+                strategies.update(row['read_strategies'])
             report['trials'].append({
                 'repetition': repetition, 'workers': count, 'elapsed_seconds': elapsed,
                 'output_bytes': sum(r['bytes'] for r in rows),
                 'read_seconds_sum': sum(r['read_seconds'] for r in rows),
                 'requests': sum(r['requests'] for r in rows), 'verified': True,
+                'read_strategies': dict(strategies),
+                'gather_workers_peak': max(r['gather_workers_peak'] for r in rows),
             })
         medians = {count: statistics.median(t['elapsed_seconds'] for t in report['trials']
                                            if t['workers'] == count) for count in allowed}
 
@@ -7,7 +7,7 @@
 from collections.abc import Callable, Generator, Iterator
 from concurrent.futures import ThreadPoolExecutor
 from contextlib import contextmanager, nullcontext
-from dataclasses import dataclass
+from dataclasses import dataclass, field
 from typing import Protocol
 
 import numpy as np
@@ -25,10 +25,17 @@ def read(self, selection: Selection) -> np.ndarray: ...
 
 @dataclass
 class ReadStats:
+    """Completed tile plans; worker limits are admitted maxima, not observed active threads.
+
+    A decoded cache hit can satisfy a parallel plan without starting any workers.
+    """
+
     requests: int = 0
     output_bytes: int = 0
     decoded_bytes_upper_bound: int = 0
     seconds: float = 0.0
+    read_strategies: dict[str, int] = field(default_factory=dict)
+    gather_workers_peak: int = 0
 
 
 class ArrayReader:
@@ -42,8 +49,9 @@ def __init__(
         self.buffer_bytes = (min(buffer_bytes, shared_budget.max_bytes)
                              if shared_budget is not None else buffer_bytes)
         self.shared_budget = shared_budget
-        if backend.spec.dense_gather and not callable(getattr(backend, "read_plan", None)):
-            raise TypeError("Dense gather metadata requires a backend with read_plan")
+        if ((backend.spec.dense_gather or backend.spec.gather_window is not None)
+                and not callable(getattr(backend, "read_plan", None))):
+            raise TypeError("Gather metadata requires a backend with read_plan")
         self._execute: Callable[[ReadPlan], np.ndarray] = getattr(
             backend, "read_plan", lambda request: backend.read(request.selection))
         if not callable(self._execute):
@@ -73,6 +81,10 @@ def stream(
                     raise ValueError("Backend returned an unexpected shape or dtype")
                 if stats is not None:
                     stats.requests += 1
+                    strategy = tile.plan.strategy
+                    stats.read_strategies[strategy] = stats.read_strategies.get(strategy, 0) + 1
+                    stats.gather_workers_peak = max(
+                        stats.gather_workers_peak, tile.plan.gather_workers)
                     stats.seconds += time.perf_counter() - start
                     stats.output_bytes += value.nbytes
                     stats.decoded_bytes_upper_bound += tile.plan.decoded_bytes_upper_bound
 
@@ -11,7 +11,8 @@
 
 from .backends import DenseGatherBackend, TensorStoreBackend, ZarrBackend, ZarristaBackend
 from .budget import ReadBudget
-from .plan import Selection
+from .gather import ParallelGatherBackend
+from .plan import GatherWindow, Selection
 from .ranges import RangePolicy
 from .reader import ArrayReader, Backend, ReadStats
 
@@ -49,6 +50,7 @@ def __init__(
         range_policy: RangePolicy | None = None,
         read_budget: ReadBudget | None = None,
         read_strategy: str = "native",
+        gather_workers: int = 4,
     ):
         if not isinstance(immutable, bool):
             raise TypeError("immutable must be a boolean")
@@ -58,8 +60,16 @@ def __init__(
         )
         if engine not in {"zarr-python", "tensorstore", "zarrista"}:
             raise ValueError(f"Unknown array engine {engine!r}")
-        if read_strategy not in {"native", "dense"}:
-            raise ValueError("read_strategy must be native or dense")
+        if read_strategy not in {"native", "dense", "parallel"}:
+            raise ValueError("read_strategy must be native, dense or parallel")
+        if (isinstance(gather_workers, bool) or not isinstance(gather_workers, int)
+                or not 1 <= gather_workers <= 256):
+            raise ValueError("Gather workers must be an integer from 1 to 256")
+        if read_strategy == "parallel" and engine != "zarrista":
+            raise ValueError("Parallel gathering currently requires the Zarrista engine")
+        if read_strategy != "parallel" and gather_workers != 4:
+            raise ValueError("gather_workers requires the parallel read strategy")
+        self.gather_workers = gather_workers
         self.read_strategy = read_strategy
         if (
             min(buffer_bytes, max_open_arrays, tensorstore_threads) < 1
@@ -141,9 +151,16 @@ def open(self, store: str | Path, array_path: str) -> ArrayReader:
         else:
             import zarrista
 
-            backend = ZarristaBackend(
-                zarrista.Array.open(zarrista.store.FilesystemStore(key[0]), path="/" + key[1])
-            )
+            array = zarrista.Array.open(zarrista.store.FilesystemStore(key[0]), path="/" + key[1])
+            native = ZarristaBackend(array)
+            backend = native
+            if self.read_strategy == "parallel":
+                if not native._partial_shards:
+                    raise ValueError("Parallel gathering requires independently decodable shards")
+                backend = ParallelGatherBackend(
+                    native, lambda: ZarristaBackend(array),
+                    GatherWindow(tuple(array.chunk_shape([0] * array.ndim)), self.gather_workers),
+                )
         if self.read_strategy == "dense":
             backend = DenseGatherBackend(backend)
         reader = ArrayReader(