[Breaking] Add multiscale data support

src/ndevio/nimage.py

@@ -125,7 +125,8 @@ class nImage(BioImage):
     # Class-level type hints for instance attributes
     path: str | None
     _is_remote: bool
-    _layer_data: xr.DataArray | None
+    _reference_xarray: xr.DataArray | None
+    _layer_data: list | None
 
     def __init__(
         self,
@@ -159,6 +160,7 @@ def __init__(
                 raise
 
         # Instance state
+        self._reference_xarray = None
         self._layer_data = None
         if isinstance(image, str | Path):
             import fsspec
@@ -180,17 +182,16 @@ def __init__(
             self._is_remote = False
 
     @property
-    def layer_data(self) -> xr.DataArray:
-        """
-        Image data as xarray DataArray for napari layer creation.
+    def reference_xarray(self) -> xr.DataArray:
+        """Image data as xarray DataArray for metadata determination.
 
-        Lazily loads data on first access. Uses in-memory or dask array
+        Lazily loads xarray on first access. Uses in-memory or dask array
         based on file size (determined automatically).
 
         Returns
         -------
         xr.DataArray
-            Squeezed image data.
+            Squeezed xarray for napari dimensions.
 
         Notes
         -----
@@ -199,13 +200,70 @@ def layer_data(self) -> xr.DataArray:
         (the default). No special mosaic handling needed here.
 
         """
-        if self._layer_data is None:
+        if self._reference_xarray is None:
+            # Ensure we're at the highest-res level for metadata consistency
+            current_res = self.current_resolution_level
+            self.set_resolution_level(0)
             if self._is_remote or not determine_in_memory(self.path):
-                self._layer_data = self.xarray_dask_data.squeeze()
+                self._reference_xarray = self.xarray_dask_data.squeeze()
             else:
-                self._layer_data = self.xarray_data.squeeze()
+                self._reference_xarray = self.xarray_data.squeeze()
+            self.set_resolution_level(current_res)
+        return self._reference_xarray
+
+    @property
+    def layer_data(self) -> list:
+        """Image data arrays shaped for napari, one per resolution level.
+
+        Returns a list of arrays ordered from highest to lowest resolution.
+        For single-resolution images the list has one element.
+        napari automatically treats multi-element lists as multiscale data.
+
+        Multiscale images are always dask-backed for memory efficiency.
+        Single-resolution images use numpy or dask based on file size.
+
+        Returns
+        -------
+        list[ArrayLike]
+            Squeezed image arrays (C dim retained for multichannel split
+            in :meth:`get_layer_data_tuples`).
+
+        """
+        if self._layer_data is None:
+            self._layer_data = self._build_layer_data()
         return self._layer_data
 
+    def _build_layer_data(self) -> list:
+        """Build the list of arrays for all resolution levels."""
+        current_res = self.current_resolution_level
+        levels = self.resolution_levels
+        multiscale = len(levels) > 1
+
+        # Determine which dims to keep from level 0's squeezed metadata.
+        # Using isel instead of squeeze ensures all levels have
+        # consistent ndim (lower levels may have extra singleton spatial dims
+        # that squeeze would incorrectly remove).
+        ref = self.reference_xarray
+        keep_dims = set(ref.dims)
+
+        arrays: list = []
+        for level in levels:
+            self.set_resolution_level(level)
+            if (
+                multiscale
+                or self._is_remote
+                or not determine_in_memory(self.path)
+            ):
+                xr_data = self.xarray_dask_data
+            else:
+                xr_data = self.xarray_data
+
+            indexer = {d: 0 for d in xr_data.dims if d not in keep_dims}
+            arrays.append(xr_data.isel(indexer).data)
+
+        self.set_resolution_level(current_res)
+        return arrays
+
     @property
     def path_stem(self) -> str:
         """Filename stem derived from path or URI, used in layer names.
@@ -331,7 +389,7 @@ def layer_axis_labels(self) -> tuple[str, ...]:
         ('Z', 'Y', 'X')
 
         """
-        layer_data = self.layer_data
+        layer_data = self.reference_xarray
 
         # Exclude Channel and Samples dimensions (RGB/multichannel handled separately)
         return tuple(
@@ -458,7 +516,8 @@ def get_layer_data_tuples(
         https://napari.org/dev/plugins/building_a_plugin/guides.html
 
         """
-        layer_data = self.layer_data
+        ref = self.reference_xarray
+        data = self.layer_data
         if layer_type is not None:
             channel_types = None  # Global override ignores per-channel
         names = self.layer_names
@@ -471,7 +530,7 @@ def get_layer_data_tuples(
         if 'S' in self.dims.order:
             return [
                 build_layer_tuple(
-                    layer_data.data,
+                    data,
                     layer_type='image',
                     name=names[0],
                     metadata=base_metadata,
@@ -485,7 +544,7 @@ def get_layer_data_tuples(
         channel_dim = 'C'
 
         # Single channel (no C dimension to split)
-        if channel_dim not in layer_data.dims:
+        if channel_dim not in ref.dims:
             channel_name = self.channel_names[0]
             effective_type = resolve_layer_type(
                 channel_name or '', layer_type, channel_types
@@ -497,7 +556,7 @@ def get_layer_data_tuples(
             )
             return [
                 build_layer_tuple(
-                    layer_data.data,
+                    data,
                     layer_type=effective_type,
                     name=names[0],
                     metadata=base_metadata,
@@ -510,8 +569,8 @@ def get_layer_data_tuples(
 
         # Multichannel - split into separate layers
         channel_names = self.channel_names
-        channel_axis = layer_data.dims.index(channel_dim)
-        total_channels = layer_data.shape[channel_axis]
+        channel_axis = ref.dims.index(channel_dim)
+        total_channels = ref.shape[channel_axis]
 
         tuples: list[LayerDataTuple] = []
         for i in range(total_channels):
@@ -520,10 +579,10 @@ def get_layer_data_tuples(
                 channel_name, layer_type, channel_types
             )
 
-            # Slice along channel axis
-            slices: list[slice | int] = [slice(None)] * layer_data.ndim
+            # Slice along channel axis for each resolution level
+            slices: list[slice | int] = [slice(None)] * ref.ndim
             slices[channel_axis] = i
-            channel_data = layer_data.data[tuple(slices)]
+            channel_data = [arr[tuple(slices)] for arr in data]
 
             extra_kwargs = (
                 channel_kwargs.get(channel_name) if channel_kwargs else None

src/ndevio/utils/_layer_utils.py

@@ -113,7 +113,7 @@ def determine_in_memory(
 
 
 def build_layer_tuple(
-    data: ArrayLike,
+    data: ArrayLike | list[ArrayLike],
     *,
     layer_type: str,
     name: str,
@@ -130,8 +130,9 @@ def build_layer_tuple(
 
     Parameters
     ----------
-    data : ArrayLike
-        Image data for this layer.
+    data : ArrayLike | list[ArrayLike]
+        Image data for this layer. A list of arrays is interpreted by
+        napari as multiscale data (highest to lowest resolution).
     layer_type : str
         Layer type ('image', 'labels', etc.).
     name : str

src/ndevio/widgets/_scene_widget.py

@@ -115,6 +115,7 @@ def open_scene(self) -> None:
             self.img.set_scene(scene_index)
 
             # Clear cached data so new scene is loaded
+            self.img._reference_xarray = None
             self.img._layer_data = None
 
             # Get layer tuples and add to viewer using napari's Layer.create()

tests/test_napari_reader.py

@@ -103,8 +103,9 @@ def test_reader_supported_formats(
 
     data, meta, _ = layer_data[0]
 
-    # Check layer data shape
-    assert data.shape == expected_shape
+    # Check layer data shape (data is a list of arrays, one per resolution level)
+    assert isinstance(data, list)
+    assert data[0].shape == expected_shape
 
     # Check meta has expected keys
     assert 'name' in meta
@@ -169,10 +170,13 @@ def test_for_multiscene_widget(
 
             data = viewer.layers[0].data
 
-            assert data.shape == expected_shape
+            if isinstance(data, list):
+                assert data[0].shape == expected_shape
+            else:
+                assert data.shape == expected_shape
         else:
             data, _, _ = reader(path)[0]
-            assert data.shape == expected_shape
+            assert data[0].shape == expected_shape
 
 
 def test_napari_get_reader_supported_formats_work(resources_dir: Path):

tests/test_nimage.py

@@ -29,9 +29,9 @@ def test_nImage_init(resources_dir: Path):
     # Shape is (T, C, Z, Y, X) = (1, 2, 60, 66, 85)
     assert img.data.shape == (1, 2, 60, 66, 85)
     # layer_data should not be loaded until accessed
-    assert img._layer_data is None
+    assert img._reference_xarray is None
     # Accessing the property triggers lazy loading
-    assert img.layer_data is not None
+    assert img.reference_xarray is not None
 
 
 def test_nImage_zarr(resources_dir: Path):
@@ -50,7 +50,7 @@ def test_nImage_remote_zarr():
     assert img.path == remote_zarr
     assert img._is_remote
     # original shape is (1, 2, 1, 512, 512) but layer_data is squeezed
-    assert img.layer_data.shape == (2, 512, 512)
+    assert img.reference_xarray.shape == (2, 512, 512)
 
 
 def test_nImage_ome_reader(resources_dir: Path):
@@ -115,7 +115,7 @@ def test_get_layer_data(resources_dir: Path):
     """Test loading napari layer data in memory."""
     img = nImage(resources_dir / CELLS3D2CH_OME_TIFF)
     # Access layer_data property to trigger loading
-    data = img.layer_data
+    data = img.reference_xarray
     # layer_data will be squeezed
     # Original shape (1, 2, 60, 66, 85) -> (2, 60, 66, 85)
     assert data.shape == (2, 60, 66, 85)
@@ -285,8 +285,11 @@ def test_multichannel_returns_tuple_per_channel(self, resources_dir: Path):
             # name should be a string (not a list)
             assert isinstance(meta['name'], str)
 
-            # Data shape should NOT include channel dimension
-            assert data.shape == (60, 66, 85)  # ZYX only
+            # Data should be a list of arrays (multiscale-ready)
+            assert isinstance(data, list)
+            assert len(data) == 1  # single resolution level
+            # Shape should NOT include channel dimension
+            assert data[0].shape == (60, 66, 85)  # ZYX only
 
             # Default layer type is "image" (channel names don't match label keywords)
             assert layer_type == 'image'

tests/test_sampledata/test_sampledata.py

@@ -36,9 +36,11 @@ def _validate_layer_data_tuples(
 
         data, kwargs, layer_type = layer_tuple
 
-        # Data should be array-like with shape
-        assert hasattr(data, 'shape')
-        assert len(data.shape) >= 2  # At minimum 2D
+        # Data should be a list of array-like objects (multiscale-ready)
+        assert isinstance(data, list)
+        assert len(data) >= 1
+        assert hasattr(data[0], 'shape')
+        assert len(data[0].shape) >= 2  # At minimum 2D
 
         # kwargs should be a dict
         assert isinstance(kwargs, dict)