Repository: kmnhan/erlabpy (language: python)

Issue:
Add support for ALBA BL20 LOREA and Diamond I05
### Implement data loaders for ~~three~~ two additional beamlines:

1. **ALBA BL20 Lorea**
   `.nxs` and `.krx` files
2. **Diamond I05**
   `.nxs` files
3. ~~**NSLS-II ID21 ESM**~~
   ~~`.pxt`, `.txt`, `.ibw`, and `.zip` files (DA30-based)~~
   Implemented in 2.10.0


Discussion / hints:

Fix patch:
diff --git a/docs/source/conf.py b/docs/source/conf.py
index 4cd91aa4..521c9dca 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -171,6 +171,11 @@ def linkcode_resolve(domain, info) -> str | None:
     "ColorType": "`ColorType <matplotlib.typing.ColorType>`",
     "RGBColorType": "`RGBColorType <matplotlib.typing.RGBColorType>`",
     "RGBAColorType": "`RGBAColorType <matplotlib.typing.RGBAColorType>`",
+    "path-like": "`path-like <os.PathLike>`",
+    "NXfield": "`NXfield <nexusformat.nexus.tree.NXfield>`",
+    "NXgroup": "`NXgroup <nexusformat.nexus.tree.NXgroup>`",
+    "NXlink": "`NXlink <nexusformat.nexus.tree.NXlink>`",
+    "NXdata": "`NXdata <nexusformat.nexus.tree.NXdata>`",
 }
 napoleon_custom_sections = [("Signals", "params_style")]
 
@@ -201,6 +206,7 @@ def linkcode_resolve(domain, info) -> str | None:
     "joblib": ("https://joblib.readthedocs.io/en/latest", None),
     "panel": ("https://panel.holoviz.org", None),
     "hvplot": ("https://hvplot.holoviz.org", None),
+    "nexusformat": ("https://nexpy.github.io/nexpy", None),
 }
 
 
diff --git a/docs/source/user-guide/curve-fitting.ipynb b/docs/source/user-guide/curve-fitting.ipynb
index f3bde09e..f3416d5b 100644
--- a/docs/source/user-guide/curve-fitting.ipynb
+++ b/docs/source/user-guide/curve-fitting.ipynb
@@ -828,6 +828,30 @@
     "result_ds"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "\n",
+    "The accessor works with any `lmfit` model, including background models from\n",
+    "[lmfitxps](https://lmfitxps.readthedocs.io/). If you have\n",
+    "[lmfitxps](https://lmfitxps.readthedocs.io/) installed, you can use the `ShirleyBG`\n",
+    "model to iteratively fit a Shirley background to the data:\n",
+    "```python\n",
+    "from lmfitxps.models import ShirleyBG\n",
+    "from lmfit.models import GaussianModel\n",
+    "\n",
+    "darr.modelfit(\"alpha\", GaussianModel() + ShirleyBG())\n",
+    "```"
+   ]
+  },
   {
    "cell_type": "raw",
    "metadata": {
diff --git a/docs/source/user-guide/io.ipynb b/docs/source/user-guide/io.ipynb
index e9eb1287..60e7b066 100644
--- a/docs/source/user-guide/io.ipynb
+++ b/docs/source/user-guide/io.ipynb
@@ -860,7 +860,7 @@
     "# Generate and save cuts with different beta values\n",
     "for i, beta in enumerate(beta_coords):\n",
     "    data = make_data(beta=beta, temp=20.0, hv=50.0)\n",
-    "    filename = f\"{tmp_dir.name}/data_001_S{str(i+1).zfill(3)}.h5\"\n",
+    "    filename = f\"{tmp_dir.name}/data_001_S{str(i + 1).zfill(3)}.h5\"\n",
     "    data.to_netcdf(filename, engine=\"h5netcdf\")\n",
     "\n",
     "# Write scan coordinates to a csv file\n",
@@ -874,7 +874,7 @@
     "# Generate some cuts with different band shifts\n",
     "for i in range(4):\n",
     "    data = make_data(beta=5.0, temp=20.0, hv=50.0, bandshift=-i * 0.05)\n",
-    "    filename = f\"{tmp_dir.name}/data_{str(i+2).zfill(3)}.h5\"\n",
+    "    filename = f\"{tmp_dir.name}/data_{str(i + 2).zfill(3)}.h5\"\n",
     "    data.to_netcdf(filename, engine=\"h5netcdf\")"
    ]
   },
@@ -1148,15 +1148,10 @@
     "        # Get all valid data files in directory\n",
     "        files = {}\n",
     "        for path in erlab.io.utils.get_files(data_dir, extensions=[\".h5\"]):\n",
-    "            # Base name\n",
-    "            data_name = os.path.splitext(os.path.basename(path))[0]\n",
-    "\n",
     "            # If multiple scans, strip the _S### part\n",
-    "            name_match = re.match(r\"(.*?_\\d{3})_(?:_S\\d{3})?\", data_name)\n",
-    "            if name_match is not None:\n",
-    "                data_name = name_match.group(1)\n",
-    "\n",
-    "            files[data_name] = path\n",
+    "            name_match = re.match(r\"(.*?_\\d{3})_(?:_S\\d{3})?\", path.stem)\n",
+    "            data_name = path.stem if name_match is None else name_match.group(1)\n",
+    "            files[data_name] = str(path)\n",
     "\n",
     "        # Map dataframe column names to data attributes\n",
     "        attrs_mapping = {\n",
diff --git a/src/erlab/accessors/general.py b/src/erlab/accessors/general.py
index 3bd0ce8d..78d85721 100644
--- a/src/erlab/accessors/general.py
+++ b/src/erlab/accessors/general.py
@@ -29,9 +29,12 @@ class PlotAccessor(ERLabDataArrayAccessor):
     def __call__(self, *args, **kwargs):
         """Plot the data.
 
-        If a 2D data array is provided, it is plotted using :func:`plot_array
-        <erlab.plotting.general.plot_array>`. Otherwise, it is equivalent to calling
-        :meth:`xarray.DataArray.plot`.
+        Plots two-dimensional data using :func:`plot_array
+        <erlab.plotting.general.plot_array>`. For non-two-dimensional data, the method
+        falls back to :meth:`xarray.DataArray.plot`.
+
+        Also sets fancy labels using :func:`fancy_labels
+        <erlab.plotting.annotations.fancy_labels>`.
 
         Parameters
         ----------
@@ -99,6 +102,8 @@ def itool(self, *args, **kwargs):
     def hvplot(self, *args, **kwargs):
         """`hvplot <https://hvplot.holoviz.org/>`_-based interactive visualization.
 
+        This method is a convenience wrapper that handles importing ``hvplot``.
+
         Parameters
         ----------
         *args
@@ -313,13 +318,13 @@ def __call__(
         Parameters
         ----------
         indexers
-            Dictionary specifying the dimensions and their values or slices.
-            Position along a dimension can be specified in three ways:
+            Dictionary specifying the dimensions and their values or slices. Position
+            along a dimension can be specified in three ways:
 
             - As a scalar value: `alpha=-1.2`
 
               If no width is specified, the data is selected along the nearest value. It
-              is equivalent to `xarray.DataArray.sel` with `method='nearest'`.
+              is equivalent to :meth:`xarray.DataArray.sel` with `method='nearest'`.
 
             - As a value and width: `alpha=5, alpha_width=0.5`
 
@@ -342,10 +347,23 @@ def __call__(
         DataArray
             The selected and averaged data.
 
-        Raises
-        ------
-        ValueError
-            If a specified dimension is not present in the data.
+
+        Note
+        ----
+        Unlike :meth:`xarray.DataArray.sel`, this method treats all dimensions without
+        coordinates as equivalent to having coordinates assigned from 0 to ``n-1``,
+        where ``n`` is the size of the dimension. For example:
+
+        .. code-block:: python
+
+            da = xr.DataArray(np.random.rand(10), dims=("x",))
+
+            da.sel(x=slice(2, 3))  # This works
+
+            da.sel(x=slice(2.0, 3.0))  # This raises a TypeError
+
+            da.qsel(x=slice(2.0, 3.0))  # This works
+
         """
         indexers = either_dict_or_kwargs(indexers, indexers_kwargs, "qsel")
 
@@ -359,9 +377,16 @@ def __call__(
                     raise ValueError(
                         f"Slice not allowed for width of dimension `{dim}`"
                     )
+
                 bin_widths[dim] = float(width)
                 if dim not in self._obj.dims:
-                    raise ValueError(f"Dimension `{dim}` not found in data.")
+                    raise ValueError(f"Dimension `{dim}` not found in data")
+            else:
+                target_dim = str(dim).removesuffix("_width")
+                if target_dim not in indexers:
+                    raise ValueError(
+                        f"`{target_dim}_width` was specified without `{target_dim}`"
+                    )
 
         scalars: dict[Hashable, float] = {}
         slices: dict[Hashable, slice] = {}
@@ -381,35 +406,53 @@ def __call__(
                         f"Slice not allowed for value of dimension `{dim}` "
                         "with width specified"
                     )
-                slices[dim] = slice(value - width / 2, value + width / 2)
+                slices[dim] = slice(value - width / 2.0, value + width / 2.0)
                 avg_dims.append(dim)
 
+        unindexed_dims: list[Hashable] = [
+            k for k in slices | scalars if k not in self._obj.indexes
+        ]
+
+        if len(unindexed_dims) >= 1:
+            out = self._obj.assign_coords(
+                {k: np.arange(self._obj.sizes[k]) for k in unindexed_dims}
+            )
+        else:
+            out = self._obj
+
         if len(scalars) >= 1:
             for k, v in scalars.items():
-                if v < self._obj[k].min() or v > self._obj[k].max():
+                if v < out[k].min() or v > out[k].max():
                     warnings.warn(
                         f"Selected value {v} for `{k}` is outside coordinate bounds",
                         stacklevel=2,
                     )
-            out = self._obj.sel(
-                {str(k): v for k, v in scalars.items()}, method="nearest"
-            )
-        else:
-            out = self._obj
+            out = out.sel({str(k): v for k, v in scalars.items()}, method="nearest")
 
         if len(slices) >= 1:
             out = out.sel(slices)
 
-            lost_coords = {k: out[k].mean() for k in avg_dims}
+            lost_coords = {
+                k: out[k].mean() for k in avg_dims if k not in unindexed_dims
+            }
             out = out.mean(dim=avg_dims, keep_attrs=True)
             out = out.assign_coords(lost_coords)
 
         if verbose:
-            print(
-                f"Selected data with {scalars} and {slices}, averaging over {avg_dims}"
-            )
+            out_str = "Selected data with "
+            if len(scalars) >= 1:
+                out_str = out_str + f"{scalars}"
+            if len(slices) >= 1:
+                if len(scalars) >= 1:
+                    out_str = out_str + " and "
+                out_str = out_str + f"{slices}"
 
-        return out
+            if len(avg_dims) >= 1:
+                out_str = out_str + f", averaging over {avg_dims}"
+
+            print(out_str)
+
+        return out.drop_vars(unindexed_dims, errors="ignore")
 
     def around(
         self, radius: float | dict[Hashable, float], *, average: bool = True, **sel_kw
diff --git a/src/erlab/analysis/fit/functions/dynamic.py b/src/erlab/analysis/fit/functions/dynamic.py
index 6042c616..cc729fba 100644
--- a/src/erlab/analysis/fit/functions/dynamic.py
+++ b/src/erlab/analysis/fit/functions/dynamic.py
@@ -154,12 +154,16 @@ class MultiPeakFunction(DynamicFunction):
         `offset`, each corresponding to the Fermi level, temperature in K, and constant
         background.
     background
-        The type of background to include in the model. The options are: ``'constant'``,
-        ``'linear'``, ``'polynomial'``, or ``'none'``. If ``'constant'``, adds a
-        ``const_bkg`` parameter. If ``'linear'``, adds a ``lin_bkg`` parameter and a
-        ``const_bkg`` parameter. If ``'polynomial'``, adds  ``c0``, ``c1``, ...
-        corresponding to the polynomial coefficients. The polynomial degree can be
-        specified with `degree`. If ``'none'``, no background is added.
+        The type of background to include in the model. Possible values are:
+
+        ============  ============================================================
+        Value         Additional parameters
+        ============  ============================================================
+        'none'        None
+        'constant'    ``const_bkg``
+        'linear'      ``lin_bkg``, ``const_bkg``
+        'polynomial'  ``c0``, ``c1``, ... depending on `degree`
+        ============  ============================================================
     degree
         The degree of the polynomial background. Only used if `background` is
         ``'polynomial'``. Default is 2.
diff --git a/src/erlab/analysis/fit/models.py b/src/erlab/analysis/fit/models.py
index e7c68b57..7b876d98 100644
--- a/src/erlab/analysis/fit/models.py
+++ b/src/erlab/analysis/fit/models.py
@@ -267,7 +267,7 @@ def __init__(
         npeaks: int = 1,
         peak_shapes: list[str] | str | None = None,
         fd: bool = True,
-        background: Literal["constant", "linear", "polynomial", "none"] = "linear",
+        background: Literal["none", "constant", "linear", "polynomial"] = "linear",
         degree: int = 2,
         convolve: bool = True,
         **kwargs,
diff --git a/src/erlab/analysis/image.py b/src/erlab/analysis/image.py
index 3843eeba..28870471 100644
--- a/src/erlab/analysis/image.py
+++ b/src/erlab/analysis/image.py
@@ -99,28 +99,29 @@ def gaussian_filter(
 
     Parameters
     ----------
-    darr
+    darr : DataArray
         The input DataArray.
-    sigma
+    sigma : float or Sequence of floats or dict
         The standard deviation(s) of the Gaussian filter in data dimensions. If a float,
         the same value is used for all dimensions, each scaled by the data step. If a
         dict, the value can be specified for each dimension using dimension names as
         keys. The filter is only applied to the dimensions specified in the dict. If a
         sequence, the values are used in the same order as the dimensions of the
         DataArray.
-    order
+    order : int or Sequence of ints or dict
         The order of the filter along each dimension. If an int, the same order is used
         for all dimensions. See Notes below for other options. Defaults to 0.
-    mode
+    mode : str or Sequence of str or dict
         The boundary mode used for the filter. If a str, the same mode is used for all
         dimensions. See Notes below for other options. Defaults to 'nearest'.
     cval
         Value to fill past edges of input if mode is 'constant'. Defaults to 0.0.
     truncate
         The truncation value used for the Gaussian filter. Defaults to 4.0.
-    radius
-        The radius of the Gaussian filter in data units. See Notes below. Defaults to
-        None.
+    radius : float or Sequence of floats or dict, optional
+        The radius of the Gaussian filter in data units. See Notes below. If specified,
+        the size of the kernel along each axis will be ``2*radius + 1``, and `truncate`
+        is ignored.
 
     Returns
     -------
@@ -367,22 +368,22 @@ def ndsavgol(
 
     Parameters
     ----------
-    arr
+    arr : array-like
         The input N-dimensional array to be filtered. The array will be cast to float64
         before filtering.
-    window_shape
+    window_shape : int or tuple of ints
         The shape of the window used for filtering. If an integer, the same size will be
         used across all axes.
-    polyorder
+    polyorder : int
         The order of the polynomial used to fit the samples. `polyorder` must be less
         than the minimum of `window_shape`.
-    deriv
+    deriv : int or tuple of ints
         The order of the derivative to compute given as a single integer or a tuple of
         integers. If an integer, the derivative of that order is computed along all
         axes. If a tuple of integers, the derivative of each order is computed along the
         corresponding dimension. The default is 0, which means to filter the data
         without differentiating.
-    delta
+    delta : float or tuple of floats
         The spacing of the samples to which the filter will be applied. If a float, the
         same value is used for all axes. If a tuple, the values are used in the same
         order as in `deriv`. The default is 1.0.
@@ -390,7 +391,7 @@ def ndsavgol(
         Must be 'mirror', 'constant', 'nearest', or 'wrap'. This determines the type of
         extension to use for the padded signal to which the filter is applied.  When
         `mode` is 'constant', the padding value is given by `cval`.
-    cval
+    cval : float
         Value to fill past the edges of the input if `mode` is 'constant'. Default is
         0.0.
     method
@@ -524,8 +525,8 @@ def _calc_savgol(values_ptr, len_values, result, data) -> int:
 
 def gradient_magnitude(
     arr: npt.NDArray[np.float64],
-    dx: np.float64,
-    dy: np.float64,
+    dx: float,
+    dy: float,
     mode: str = "nearest",
     cval: float = 0.0,
 ) -> npt.NDArray[np.float64]:
diff --git a/src/erlab/interactive/imagetool/__init__.py b/src/erlab/interactive/imagetool/__init__.py
index f646dbfd..6a11dac4 100644
--- a/src/erlab/interactive/imagetool/__init__.py
+++ b/src/erlab/interactive/imagetool/__init__.py
@@ -40,7 +40,7 @@
 from erlab.interactive.utils import DictMenuBar, copy_to_clipboard, generate_code
 
 if TYPE_CHECKING:
-    from collections.abc import Callable, Collection
+    from collections.abc import Callable, Collection, Hashable
 
     from erlab.interactive.imagetool.core import ImageSlicerState
     from erlab.interactive.imagetool.slicer import ArraySlicer
@@ -518,6 +518,11 @@ def _generate_menu_kwargs(self) -> dict:
                         "toggled": self._set_colormap_options,
                         "sep_after": True,
                     },
+                    "Normalize": {"triggered": self._normalize},
+                    "Reset": {
+                        "triggered": self._reset_filters,
+                        "sep_after": True,
+                    },
                 },
             },
             "editMenu": {
@@ -671,6 +676,15 @@ def _rotate(self) -> None:
         dialog = RotationDialog(self.slicer_area)
         dialog.exec()
 
+    @QtCore.Slot()
+    def _normalize(self) -> None:
+        dialog = NormalizeDialog(self.slicer_area)
+        dialog.exec()
+
+    @QtCore.Slot()
+    def _reset_filters(self) -> None:
+        self.slicer_area.apply_func(lambda x: x)
+
     def _set_colormap_options(self) -> None:
         self.slicer_area.set_colormap(
             reverse=self.colorAct[0].isChecked(),
@@ -784,10 +798,8 @@ def _to_hdf5(darr: xr.DataArray, file: str, **kwargs) -> None:
             fn(self.slicer_area._data, files[0], **kargs)
 
 
-class _DataEditDialog(QtWidgets.QDialog):
+class _DataManipulationDialog(QtWidgets.QDialog):
     title: str | None = None
-    prefix: str = ""
-    suffix: str = ""
     show_copy_button: bool = False
 
     def __init__(self, slicer_area: ImageSlicerArea) -> None:
@@ -802,10 +814,6 @@ def __init__(self, slicer_area: ImageSlicerArea) -> None:
 
         self.setup_widgets()
 
-        self.new_window_check = QtWidgets.QCheckBox("Open in New Window")
-        self.new_window_check.setChecked(True)
-        self.layout_.addRow(self.new_window_check)
-
         self.buttonBox = QtWidgets.QDialogButtonBox(
             QtWidgets.QDialogButtonBox.StandardButton.Ok
             | QtWidgets.QDialogButtonBox.StandardButton.Cancel
@@ -831,6 +839,35 @@ def layout_(self) -> QtWidgets.QFormLayout:
     def array_slicer(self) -> ArraySlicer:
         return self.slicer_area.array_slicer
 
+    def setup_widgets(self) -> None:
+        # Overridden by subclasses
+        pass
+
+    def process_data(self) -> xr.DataArray:
+        # Overridden by subclasses
+        return self.slicer_area.data
+
+    def make_code(self) -> str:
+        # Overridden by subclasses
+        return ""
+
+
+class _DataTransformDialog(_DataManipulationDialog):
+    """
+    Parent class for implementing data changes that affect both shape and values.
+
+    These changes are destructive and cannot be undone.
+    """
+
+    prefix: str = ""
+    suffix: str = ""
+
+    def __init__(self, slicer_area: ImageSlicerArea) -> None:
+        super().__init__(slicer_area)
+        self.new_window_check = QtWidgets.QCheckBox("Open in New Window")
+        self.new_window_check.setChecked(True)
+        self.layout_.insertRow(-1, self.new_window_check)
+
     @QtCore.Slot()
     def accept(self) -> None:
         if self.slicer_area.data.name is not None:
@@ -844,20 +881,45 @@ def accept(self) -> None:
             self.slicer_area.set_data(self.process_data().rename(new_name))
         super().accept()
 
-    def setup_widgets(self) -> None:
-        # Overridden by subclasses
-        pass
 
-    def process_data(self) -> xr.DataArray:
-        # Overridden by subclasses
-        return self.slicer_area.data
+class _DataFilterDialog(_DataManipulationDialog):
+    """Parent class for implementing data changes that affect only the values."""
 
-    def make_code(self) -> str:
-        # Overridden by subclasses
-        return ""
+    enable_preview: bool = True
+
+    def __init__(self, slicer_area: ImageSlicerArea) -> None:
+        super().__init__(slicer_area)
+        self._previewed: bool = False
+
+        if self.enable_preview:
+            self.preview_button = QtWidgets.QPushButton("Preview")
+            self.preview_button.clicked.connect(self._preview)
+            self.buttonBox.addButton(
+                self.preview_button, QtWidgets.QDialogButtonBox.ButtonRole.ActionRole
+            )
+
+    @QtCore.Slot()
+    def _preview(self):
+        self._previewed = True
+        self.slicer_area.apply_func(self.func)
+
+    @QtCore.Slot()
+    def reject(self) -> None:
+        if self._previewed:
+            self.slicer_area.apply_func(None)
+        super().reject()
+
+    @QtCore.Slot()
+    def accept(self) -> None:
+        self.slicer_area.apply_func(self.func)
+        super().accept()
+
+    def func(self, data: xr.DataArray) -> xr.DataArray:
+        # Implement this method in subclasses
+        raise NotImplementedError
 
 
-class RotationDialog(_DataEditDialog):
+class RotationDialog(_DataTransformDialog):
     suffix = " Rotated"
     show_copy_button = True
 
@@ -933,3 +995,65 @@ def make_code(self) -> str:
         return generate_code(
             rotate, [f"|{placeholder}|"], self._rotate_params, module="era.transform"
         )
+
+
+class NormalizeDialog(_DataFilterDialog):
+    title = "Normalize"
+    show_copy_button = False
+
+    def setup_widgets(self) -> None:
+        dim_group = QtWidgets.QGroupBox("Dimensions")
+        dim_layout = QtWidgets.QVBoxLayout()
+        dim_group.setLayout(dim_layout)
+
+        self.dim_checks: dict[Hashable, QtWidgets.QCheckBox] = {}
+
+        for d in self.slicer_area.data.dims:
+            self.dim_checks[d] = QtWidgets.QCheckBox(str(d))
+            dim_layout.addWidget(self.dim_checks[d])
+
+        option_group = QtWidgets.QGroupBox("Options")
+        option_layout = QtWidgets.QVBoxLayout()
+        option_group.setLayout(option_layout)
+
+        self.opts: list[QtWidgets.QRadioButton] = []
+        self.opts.append(QtWidgets.QRadioButton("Data/Area"))
+        self.opts.append(QtWidgets.QRadioButton("(Data−m)/(M−m)"))
+        self.opts.append(QtWidgets.QRadioButton("Data−m"))
+        self.opts.append(QtWidgets.QRadioButton("(Data−m)/Area"))
+
+        self.opts[0].setChecked(True)
+        for opt in self.opts:
+            option_layout.addWidget(opt)
+
+        self.layout_.addRow(dim_group)
+        self.layout_.addRow(option_group)
+
+    def func(self, data: xr.DataArray) -> xr.DataArray:
+        norm_dims = tuple(k for k, v in self.dim_checks.items() if v.isChecked())
+        if len(norm_dims) == 0:
+            return data
+
+        calc_area: bool = self.opts[0].isChecked() or self.opts[3].isChecked()
+        calc_minimum: bool = not self.opts[0].isChecked()
+        calc_maximum: bool = self.opts[1].isChecked()
+
+        if calc_area:
+            area = data.mean(norm_dims)
+
+        if calc_minimum:
+            minimum = data.min(norm_dims)
+
+        if calc_maximum:
+            maximum = data.max(norm_dims)
+
+        if self.opts[0].isChecked():
+            return data / area
+
+        if self.opts[1].isChecked():
+            return (data - minimum) / (maximum - minimum)
+
+        if self.opts[2].isChecked():
+            return data - minimum
+
+        return (data - minimum) / area
diff --git a/src/erlab/interactive/imagetool/core.py b/src/erlab/interactive/imagetool/core.py
index b3e824a7..3c4c1382 100644
--- a/src/erlab/interactive/imagetool/core.py
+++ b/src/erlab/interactive/imagetool/core.py
@@ -31,8 +31,6 @@
 if TYPE_CHECKING:
     from collections.abc import Callable, Iterable, Sequence
 
-    from pyqtgraph.graphicsItems.ViewBox import ViewBoxMenu
-
     from erlab.interactive.imagetool.slicer import ArraySlicerState
 
 
@@ -123,9 +121,6 @@ def __init__(
     def getPlotItem(self) -> ItoolPlotItem:
         return self.plotItem
 
-    def getPlotItemViewBox(self) -> pg.ViewBox:
-        return self.getPlotItem().vb
-
 
 def suppress_history(method: Callable | None = None):
     """Ignore history changes made within the decorated method."""
@@ -605,10 +600,6 @@ def linked_slicers(self) -> set[ImageSlicerArea]:
     def colormap(self) -> str | pg.ColorMap:
         return self.colormap_properties["cmap"]
 
-    @colormap.setter
-    def colormap(self, cmap: str | pg.ColorMap) -> None:
-        self.set_colormap(cmap)
-
     @property
     def levels_locked(self) -> bool:
         return self.colormap_properties["levels_locked"]
@@ -709,8 +700,10 @@ def history_suppressed(self):
             self._write_history = original
 
     def on_close(self) -> None:
-        self.array_slicer.clear_cache()
-        self.data.close()
+        if hasattr(self, "array_slicer"):
+            self.array_slicer.clear_cache()
+        if hasattr(self, "data"):
+            self.data.close()
         if hasattr(self, "_data") and self._data is not None:
             self._data.close()
             del self._data
@@ -975,6 +968,35 @@ def update_values(
             self.array_slicer.clear_val_cache(include_vals=True)
             self.refresh_all(only_plots=True)
 
+            # This will update colorbar limits if visible
+            self.lock_levels(self.levels_locked)
+
+    def apply_func(self, func: Callable[[xr.DataArray], xr.DataArray] | None) -> None:
+        """Apply a function to the data.
+
+        The function must accept the data as the first argument and return a new
+        DataArray. The returned DataArray must have the same dimensions and coordinates
+        as the original data.
+
+        This action is not recorded in the history, and the data is not affected. Only
+        one function can be applied at a time.
+
+        Parameters
+        ----------
+        func
+            The function to apply to the data. if None, the data is restored.
+
+        """
+        # self._data is original data passed to `set_data`
+        # self.data is the current data transformed by ArraySlicer
+        if self._data is None:
+            return
+
+        if func is None:
+            self.update_values(self._data)
+        else:
+            self.update_values(func(self._data))
+
     @QtCore.Slot(int, int)
     @link_slicer
     @record_history
@@ -1539,12 +1561,23 @@ def __init__(
         for act in ["Transforms", "Downsample", "Average", "Alpha", "Points"]:
             self.setContextMenuActionVisible(act, False)
 
+        self.vb.setCursor(QtGui.QCursor(QtCore.Qt.CursorShape.CrossCursor))
+
+        for i in (0, 1):
+            # Hide unnecessary menu items
+            self.vb.menu.ctrl[i].linkCombo.setVisible(False)
+            self.vb.menu.ctrl[i].label.setVisible(False)
+
+        self.vb.menu.addSeparator()
+
         save_action = self.vb.menu.addAction("Save data as HDF5")
         save_action.triggered.connect(self.save_current_data)
 
         copy_code_action = self.vb.menu.addAction("Copy selection code")
         copy_code_action.triggered.connect(self.copy_selection_code)
 
+        self.vb.menu.addSeparator()
+
         if image:
             goldtool_action = self.vb.menu.addAction("Open in goldtool")
             self._goldtool: None | QtWidgets.QWidget = None
@@ -1554,10 +1587,7 @@ def __init__(
             self._dtool: None | QtWidgets.QWidget = None
             dtool_action.triggered.connect(self.open_in_dtool)
 
-        for i in (0, 1):
-            self.getViewBoxMenu().ctrl[i].linkCombo.setVisible(False)
-            self.getViewBoxMenu().ctrl[i].label.setVisible(False)
-        self.getViewBox().setCursor(QtGui.QCursor(QtCore.Qt.CursorShape.CrossCursor))
+            self.vb.menu.addSeparator()
 
         self.slicer_area = slicer_area
         self.display_axis = display_axis
@@ -1718,9 +1748,6 @@ def getMenu(self) -> QtWidgets.QMenu:
     def getViewBox(self) -> pg.ViewBox:
         return self.vb
 
-    def getViewBoxMenu(self) -> ViewBoxMenu:
-        return self.getViewBox().menu
-
     def mouseDragEvent(
         self, ev: mouseEvents.MouseDragEvent | mouseEvents.MouseClickEvent
     ) -> None:
@@ -1958,9 +1985,11 @@ def _print_angle():
             self._guideline_angle = line.angle_effective
             self._guideline_offset = [line_pos.x(), line_pos.y()]
             for i in range(2):
-                self._guideline_offset[i] = np.round(
-                    self._guideline_offset[i],
-                    self.array_slicer.get_significant(self.display_axis[i]),
+                self._guideline_offset[i] = float(
+                    np.round(
+                        self._guideline_offset[i],
+                        self.array_slicer.get_significant(self.display_axis[i]),
+                    )
                 )
 
             self.setTitle(
@@ -2176,4 +2205,4 @@ def setVisible(self, visible: bool) -> None:
         self.cb._span.blockSignals(not visible)
 
         if visible:
-            self.cb._span.setRegion(self.cb.limits)
+            self.cb.setSpanRegion(self.cb.limits)
diff --git a/src/erlab/interactive/imagetool/manager.py b/src/erlab/interactive/imagetool/manager.py
index 341eef8d..25d19786 100644
--- a/src/erlab/interactive/imagetool/manager.py
+++ b/src/erlab/interactive/imagetool/manager.py
@@ -19,7 +19,6 @@
     "main",
     "show_in_manager",
 ]
-
 import contextlib
 import gc
 import os
@@ -30,6 +29,7 @@
 import tempfile
 import threading
 import time
+import traceback
 import uuid
 from multiprocessing import shared_memory
 from typing import TYPE_CHECKING, Any, cast
@@ -149,8 +149,11 @@ def run(self) -> None:
                 EOFError,
                 ImportError,
                 IndexError,
-            ) as e:
-                print("Failed to unpickle data:", e)
+            ):
+                print(
+                    f"Failed to unpickle data due to the following error:\n"
+                    f"{traceback.format_exc()}"
+                )
 
             conn.close()
 
@@ -738,6 +741,7 @@ def is_running() -> bool:
     except FileExistsError:
         return True
     else:
+        shm.close()
         shm.unlink()
         return False
 
@@ -750,6 +754,7 @@ def main() -> None:
     qapp = QtWidgets.QApplication(sys.argv)
     qapp.setStyle("Fusion")
     qapp.setWindowIcon(QtGui.QIcon(os.path.join(os.path.dirname(__file__), "icon.png")))
+    qapp.setApplicationDisplayName("ImageTool Manager")
 
     while is_running():
         dialog = _InitDialog()
diff --git a/src/erlab/interactive/imagetool/slicer.py b/src/erlab/interactive/imagetool/slicer.py
index 14ca3b62..ce7c4942 100644
--- a/src/erlab/interactive/imagetool/slicer.py
+++ b/src/erlab/interactive/imagetool/slicer.py
@@ -324,7 +324,7 @@ def validate_array(data: xr.DataArray) -> xr.DataArray:
             The converted data.
 
         """
-        data = data.squeeze()
+        data = data.copy().squeeze()
 
         if data.ndim < 2:
             raise ValueError("Data must have at least two dimensions.")
diff --git a/src/erlab/io/__init__.py b/src/erlab/io/__init__.py
index 9aafa423..b132bea8 100644
--- a/src/erlab/io/__init__.py
+++ b/src/erlab/io/__init__.py
@@ -15,6 +15,7 @@
    dataloader
    utils
    igor
+   nexusutils
    exampledata
    characterization
 
diff --git a/src/erlab/io/dataloader.py b/src/erlab/io/dataloader.py
index af532cd4..49047a58 100644
--- a/src/erlab/io/dataloader.py
+++ b/src/erlab/io/dataloader.py
@@ -354,6 +354,7 @@ def load(
         single: bool = False,
         combine: bool = True,
         parallel: bool | None = None,
+        load_kwargs: dict[str, Any] | None = None,
         **kwargs,
     ) -> (
         xr.DataArray
@@ -404,6 +405,9 @@ def load(
         parallel
             Whether to load multiple files in parallel. If not specified, files are
             loaded in parallel only when there are more than 15 files to load.
+        load_kwargs
+            Additional keyword arguments to be passed to :meth:`load_single
+            <erlab.io.dataloader.LoaderBase.load_single>`.
         **kwargs
             Additional keyword arguments are passed to `identify`.
 
@@ -455,6 +459,8 @@ def load(
         """
         if self.always_single:
             single = True
+        if load_kwargs is None:
+            load_kwargs = {}
 
         if isinstance(identifier, int):
             if data_dir is None:
@@ -470,18 +476,20 @@ def load(
             if len(file_paths) == 1:
                 # Single file resolved
                 data: xr.DataArray | xr.Dataset | DataTree = self.load_single(
-                    file_paths[0]
+                    file_paths[0], **load_kwargs
                 )
             else:
                 # Multiple files resolved
                 if combine:
                     data = self.combine_multiple(
-                        self.load_multiple_parallel(file_paths, parallel=parallel),
+                        self.load_multiple_parallel(
+                            file_paths, parallel=parallel, **load_kwargs
+                        ),
                         coord_dict,
                     )
                 else:
                     return self.load_multiple_parallel(
-                        file_paths, parallel=parallel, post_process=True
+                        file_paths, parallel=parallel, post_process=True, **load_kwargs
                     )
 
         else:
@@ -507,12 +515,18 @@ def load(
 
                     new_kwargs = kwargs | additional_kwargs
                     return self.load(
-                        new_identifier, new_dir, single=single, **new_kwargs
+                        new_identifier,
+                        new_dir,
+                        single=single,
+                        combine=combine,
+                        parallel=parallel,
+                        **new_kwargs,
                     )
+
                 # On failure, assume single file
                 single = True
 
-            data = self.load_single(identifier)
+            data = self.load_single(identifier, **load_kwargs)
 
         data = self.post_process_general(data)
 
@@ -932,7 +946,7 @@ def load_single(
 
     def identify(
         self, num: int, data_dir: str | os.PathLike
-    ) -> tuple[list[str], dict[str, Sequence]] | None:
+    ) -> tuple[list[str | os.PathLike], dict[str, Sequence]] | None:
         """Identify the files and coordinates for a given scan number.
 
         This method takes a scan index and transforms it into a list of file paths and
@@ -961,7 +975,7 @@ def identify(
 
         Returns
         -------
-        files : list[str]
+        files : list of str or path-like
             A list of file paths.
         coord_dict : dict[str, Sequence]
             A dictionary mapping scan axes names to scan coordinates. For scans spread
@@ -1354,6 +1368,7 @@ def load_multiple_parallel(
         file_paths: list[str],
         parallel: bool | None = None,
         post_process: bool = False,
+        **kwargs,
     ) -> list[xr.DataArray] | list[xr.Dataset] | list[DataTree]:
         """Load multiple files in parallel.
 
@@ -1365,6 +1380,9 @@ def load_multiple_parallel(
             If `True`, data loading will be performed in parallel using `dask.delayed`.
         post_process
             Whether to post-process each data object after loading.
+        **kwargs
+            Additional keyword arguments to be passed to :meth:`load_single
+            <erlab.io.dataloader.LoaderBase.load_single>`.
 
         Returns
         -------
@@ -1376,10 +1394,12 @@ def load_multiple_parallel(
         if post_process:
 
             def _load_func(filename):
-                return self.load(filename, single=True)
+                return self.load(filename, single=True, load_kwargs=kwargs)
 
         else:
-            _load_func = self.load_single
+
+            def _load_func(filename):
+                return self.load_single(filename, **kwargs)
 
         if parallel:
             import joblib
diff --git a/src/erlab/io/nexusutils.py b/src/erlab/io/nexusutils.py
new file mode 100644
index 00000000..ae982f8c
--- /dev/null
+++ b/src/erlab/io/nexusutils.py
@@ -0,0 +1,444 @@
+"""Utilities for reading NeXus files into xarray objects."""
+
+import importlib
+
+if not importlib.util.find_spec("nexusformat"):
+    raise ImportError("The `nexusformat` package is required to read NeXus files")
+
+import os
+from collections.abc import Callable, Hashable, Iterable, Mapping, Sequence
+from typing import Any
+
+import h5py
+import numpy as np
+import xarray as xr
+from nexusformat.nexus import nxload
+from nexusformat.nexus.tree import NXattr, NXentry, NXfield, NXgroup, NXlink
+
+
+def _parse_value(value):
+    """Convert numpy scalars or bytes to native Python types and regular strings."""
+    if isinstance(value, bytes | bytearray):
+        return value.decode("utf-8")
+    if isinstance(value, NXattr):
+        return _parse_value(value.nxdata)
+    if isinstance(value, np.ndarray) and value.size == 1:
+        return _parse_value(np.atleast_1d(value)[0])
+    if isinstance(value, np.number):
+        # Convert to native Python type
+        return value.item()
+    return value
+
+
+def _parse_field(field):
+    """Convert a NeXus field to a native Python type."""
+    if not isinstance(field, NXfield):
+        return field
+
+    if field.size == 1:
+        return _parse_value(field.nxdata)
+
+    return field.nxdata
+
+
+def _remove_axis_attrs(attrs: Mapping[str, Any]) -> dict[str, Any]:
+    out: dict[str, Any] = {}
+    for k, v in attrs.items():
+        if k in ("axis", "primary", "target"):
+            continue
+        out[k] = _parse_value(v)
+
+    return out
+
+
+def _parse_h5py(obj: h5py.Group | h5py.Dataset, out: dict) -> None:
+    if isinstance(obj, h5py.Group):
+        for v in obj.values():
+            _parse_h5py(v, out)
+    elif isinstance(obj, h5py.Dataset):
+        out[obj.name] = _parse_value(obj[()])
+
+
+def _parse_group(
+    group: NXgroup, out: dict[str, Any], exclude: Sequence[str], parse: bool
+) -> None:
+    """Recursively parse a NeXus group and its items into a nested dictionary.
+
+    The keys are the absolute nxpath of the items, and the values are the corresponding
+    NXfield objects. The function is called recursively for each sub-group. NXlink
+    objects are skipped.
+
+    Parameters
+    ----------
+    group
+        The NeXus group to be parsed.
+    out
+        The dictionary to store the parsed items.
+    exclude
+        List of full paths to exclude from the output.
+    parse
+        Wheter to parse the values of NXfields to native Python types.
+
+    Note
+    ----
+    Sometimes the items in the group may be invalid NXobjects. In that case, the
+    function will attempt to get the values from h5py as a fallback.
+
+    """
+    for item in group.values():
+        if item.nxpath in exclude:
+            continue
+
+        if isinstance(item, NXgroup):
+            _parse_group(item, out, exclude, parse)
+
+        elif isinstance(item, NXlink):
+            # Skip links
+            continue
+
+        elif isinstance(item, NXfield):
+            if parse:
+                out[item.nxpath] = _parse_field(item)
+            else:
+                out[item.nxpath] = item
+
+        else:
+            # Sometimes item may be an invalid NXobject
+            # In that case, get the values from h5py as a fallback
+            _parse_h5py(item.nxfile.get(item.nxpath), out)
+
+
+def _get_primary_coords(group: NXgroup, out: list[NXfield]):
+    for item in group.values():
+        if isinstance(item, NXgroup):
+            _get_primary_coords(item, out)
+        elif isinstance(item, NXlink):
+            # Skip links
+            continue
+        elif (
+            isinstance(item, NXfield)
+            and "primary" in item.attrs
+            and int(item.primary) == 1
+        ):
+            out.append(item)
+
+
+def _get_non_primary_coords(group: NXgroup, out: list[NXfield]):
+    for item in group.values():
+        if isinstance(item, NXgroup):
+            _get_non_primary_coords(item, out)
+        elif isinstance(item, NXlink):
+            # Skip links
+            continue
+        elif (
+            isinstance(item, NXfield)
+            and "axis" in item.attrs
+            and ("primary" not in item.attrs or int(item.primary) != 1)
+        ):
+            out.append(item)
+
+
+def get_primary_coords(group: NXgroup) -> list[NXfield]:
+    """Get all primary coordinates in a group.
+
+    Retrieves all fields with the attribute `primary=1` in the group and its subgroups
+    recursively. The output list is sorted by the `axis` attribute of the fields.
+
+    Parameters
+    ----------
+    group
+        The group to search for the coordinates.
+
+    Returns
+    -------
+    list of NXfield
+
+    """
+    out: list[NXfield] = []
+    _get_primary_coords(group, out)
+    return sorted(out, key=lambda field: int(field.axis))
+
+
+def get_non_primary_coords(group: NXgroup) -> list[NXfield]:
+    """Get all non-primary coordinates in a group.
+
+    Retrieves all fields with the attribute `axis` in the group and its subgroups
+    recursively. The output list in the order of traversal.
+
+    Parameters
+    ----------
+    group
+        The group to search for the coordinates.
+
+    Returns
+    -------
+    list of NXfield
+
+    """
+    out: list[NXfield] = []
+    _get_non_primary_coords(group, out)
+    return out
+
+
+def get_primary_coord_dict(
+    fields: list[NXfield],
+) -> tuple[tuple[str, ...], dict[str, xr.DataArray | tuple]]:
+    """Generate a dictionary of primary coordinates from a list of NXfields.
+
+    The output dictionary contains the absolute NXpath of the field as the key and the
+    value is either a `xarray.DataArray` or a tuple of the associated dimensions, data,
+    and attributes.
+
+    The output arguments can be directly used as arguments for the `xarray.DataArray`
+    constructor.
+
+    If there are multiple fields with the same axis, one of the fields will depend on
+    the other field.
+
+    Parameters
+    ----------
+    fields : list of NXfield
+        The list of primary coordinates to be converted. The list can be obtained from a
+        NXgroup using :func:`get_primary_coords`.
+
+    Returns
+    -------
+    dims : tuple of str
+        The dimension names of the primary coordinates.
+    coords : dict of str to DataArray or tuple
+        The dictionary of primary coordinates in the group.
+    """
+    coords: dict[str, xr.DataArray | tuple] = {}
+
+    # Dict to store processed axes
+    # Ensure that if there are multiple fields with the same axis, only one is used
+    processed_axes: dict[int, NXfield] = {}
+
+    for field in fields:
+        if field.ndim == 1:
+            if int(field.axis) in processed_axes:
+                # Axis already added, make this coordinate depend on that dimension
+                coords[field.nxpath] = (
+                    processed_axes[int(field.axis)].nxpath,
+                    field.nxdata,
+                    _remove_axis_attrs(field.attrs),
+                )
+                continue
+
+            coords[field.nxpath] = nxfield_to_xarray(field)
+
+        else:
+            # Multidimensional coordinate (like kinetic energy in hv-dependent scans)
+            assoc_dims: list[str] = []
+            for num in field.shape:
+                for field in fields:
+                    if field.ndim == 1 and field.shape == (num,):
+                        assoc_dims.append(field.nxpath)
+                        break
+                else:
+                    # Shape not in list, must be new dimension
+                    assoc_dims.append(field.nxpath)
+
+            coords[field.nxpath] = (
+                tuple(assoc_dims),
+                field.nxdata,
+                _remove_axis_attrs(field.attrs),
+            )
+
+        processed_axes[int(field.axis)] = field
+
+    dims = tuple(f.nxpath for f in dict(sorted(processed_axes.items())).values())
+
+    return dims, coords
+
+
+def get_coord_dict(
+    group: NXgroup,
+) -> tuple[tuple[str, ...], dict[str, xr.DataArray | tuple]]:
+    """Generate a dictionary of coordinates from a NeXus group.
+
+    Parameters
+    ----------
+    group : NXgroup
+        The NeXus group from which to extract coordinates.
+
+    Returns
+    -------
+    dims : tuple of str
+        The dimension names of the primary coordinates.
+    coords : dict of str to DataArray or tuple
+        The dictionary of all coordinates in the group.
+    """
+    fields_primary: list[NXfield] = get_primary_coords(group)
+    fields_non_primary: list[NXfield] = get_non_primary_coords(group)
+
+    dims, coords = get_primary_coord_dict(fields_primary)
+
+    for field in fields_non_primary:
+        if field.size == 1:
+            # Scalar coordinate
+            coords[field.nxpath] = nxfield_to_xarray(field, no_dims=True)
+
+        else:
+            # Coord depends on some other primary coordinate
+            associated_primary: NXfield = group.nxroot[dims[int(field.axis) - 1]]
+            coords[field.nxpath] = (
+                associated_primary.nxpath,
+                field.nxdata,
+                _remove_axis_attrs(field.attrs),
+            )
+
+    return dims, coords
+
+
+def nexus_group_to_dict(
+    group: NXgroup,
+    exclude: Sequence[str] | None,
+    relative: bool = True,
+    replace_slash: bool = True,
+    parse: bool = False,
+) -> dict[str, Any]:
+    """Convert a NeXus group to a dictionary.
+
+    This function takes a NeXus group and converts it into a dictionary where the keys
+    are the paths of the group's items relative to the group's root path, and the values
+    are the corresponding items.
+
+    Parameters
+    ----------
+    group
+        The NeXus group to be converted.
+    exclude
+        List of paths to exclude from the output.
+    relative
+        Wheter to use the relative or absolute paths of the items. If `True`, the keys
+        are the paths of the items relative to the path of the group. If `False`, the
+        keys are the absolute paths of the items relative to the root of the NeXus file.
+    replace_slash
+        Wheter to replace the slashes in the paths with dots.
+    parse
+        Wheter to coerce the values of NXfields to native Python types.
+
+    """
+    if exclude is None:
+        exclude = []
+
+    prefix: str = group.nxpath
+    out: dict[str, Any] = {}
+    _parse_group(group, out, exclude, parse)
+
+    if relative:
+        out = {k.removeprefix(prefix).lstrip("/"): v for k, v in out.items()}
+
+    if replace_slash:
+        out = {k.replace("/", "."): v for k, v in out.items()}
+
+    return out
+
+
+def nxfield_to_xarray(field: NXfield, no_dims: bool = False) -> xr.DataArray:
+    """Convert a coord-like 1D NeXus field to a single `xarray.DataArray`.
+
+    Parameters
+    ----------
+    field
+        The NeXus field to be converted.
+
+    """
+    attrs = _remove_axis_attrs(field.attrs)
+
+    if no_dims:
+        if not field.size == 1:
+            raise ValueError(
+                "The field must have a single value "
+                "to be converted to a non-dimensional array"
+            )
+        return xr.DataArray(field.nxdata.item(), attrs=attrs)
+
+    if field.ndim > 1:
+        raise ValueError("The field must have a single dimension to be converted")
+
+    return xr.DataArray(field, dims=(field.nxpath,), attrs=attrs)
+
+
+def nxgroup_to_xarray(
+    group: NXgroup, data: str | Callable[[NXgroup], NXfield]
+) -> xr.DataArray:
+    """Convert a NeXus group to an xarray DataArray.
+
+    Parameters
+    ----------
+    group : NXgroup
+        The NeXus group to be converted.
+    data : str or callable
+        The location of the data values which can be specified in two ways:
+
+        - If a string, it must be the relative path from ``group`` to the `NXfield
+          <nexusformat.nexus.tree.NXfield>` containing the data values.
+
+        - If a callable, it must be a function that takes ``group`` as an argument and
+          returns the `NXfield <nexusformat.nexus.tree.NXfield>` containing the data
+          values.
+
+    Returns
+    -------
+    xarray.DataArray
+        The DataArray containing the data. Dimension and coordinate names are the
+        relative paths of the corresponding NXfields, with the slashes replaced by dots.
+
+    """
+    if callable(data):
+        values: NXfield = data(group)
+    else:
+        values = group[data]
+
+    if isinstance(values, NXlink):
+        values = values.nxlink
+
+    dims, coords = get_coord_dict(group)
+    excluded_paths = [*list(coords.keys()), values.nxpath]
+
+    attrs = nexus_group_to_dict(group, exclude=excluded_paths, parse=True)
+
+    # Strip path prefix from dimensions and coordinates
+    prefix: str = group.nxpath
+
+    def _make_relative(s: Hashable) -> str:
+        return str(s).removeprefix(prefix).lstrip("/").replace("/", ".")
+
+    def _make_coord_relative(t: xr.DataArray | tuple) -> xr.DataArray | tuple:
+        if isinstance(t, xr.DataArray):
+            if len(t.dims) == 1:
+                return t.rename({t.dims[0]: _make_relative(t.dims[0])})
+
+            return t
+
+        _tmp_list = list(t)
+        if isinstance(t[0], str):
+            _tmp_list[0] = _make_relative(t[0])
+        elif isinstance(t[0], Iterable):
+            _tmp_list[0] = tuple(_make_relative(s) for s in t[0])
+        return tuple(_tmp_list)
+
+    dims = tuple(_make_relative(d) for d in dims)
+    coords = {_make_relative(k): _make_coord_relative(v) for k, v in coords.items()}
+    return xr.DataArray(values, dims=dims, coords=coords, attrs=attrs)
+
+
+def get_entry(filename: str | os.PathLike, entry: str | None = None) -> NXentry:
+    """Get an NXentry object from a NeXus file.
+
+    Parameters
+    ----------
+    filename
+        The path to the NeXus file.
+    entry
+        The path of the entry to get. If `None`, the first entry in the file is
+        returned.
+
+    """
+    root = nxload(filename)
+    if entry is None:
+        return next(iter(root.entries.values()))
+
+    return root[entry]
diff --git a/src/erlab/io/plugins/esm.py b/src/erlab/io/plugins/esm.py
index 3a2f6483..0ccbdc91 100644
--- a/src/erlab/io/plugins/esm.py
+++ b/src/erlab/io/plugins/esm.py
@@ -2,7 +2,6 @@
 
 import os
 import re
-from collections.abc import Sequence
 from typing import ClassVar
 
 import erlab.io.utils
@@ -18,24 +17,24 @@ class ESMLoader(DA30Loader):
 
     additional_attrs: ClassVar[dict] = {"configuration": 3}
 
-    def identify(
-        self, num: int, data_dir: str | os.PathLike
-    ) -> tuple[list[str], dict[str, Sequence]] | None:
+    def identify(self, num: int, data_dir: str | os.PathLike):
         for file in erlab.io.utils.get_files(
             data_dir, extensions=(".ibw", ".pxt", ".zip")
         ):
-            if file.endswith(".zip"):
-                match = re.match(r"(.*?)" + str(num).zfill(4) + r".zip", file)
+            match file.suffix:
+                case ".zip":
+                    m = re.match(r"(.*?)" + str(num).zfill(4) + r".zip", file.name)
 
-            elif file.endswith(".pxt"):
-                match = re.match(r"(.*?)" + str(num).zfill(4) + r".pxt", file)
+                case ".pxt":
+                    m = re.match(r"(.*?)" + str(num).zfill(4) + r".pxt", file.name)
 
-            elif file.endswith(".ibw"):
-                match = re.match(
-                    r"(.*?)" + str(num).zfill(4) + str(num).zfill(3) + r".ibw", file
-                )
+                case ".ibw":
+                    m = re.match(
+                        r"(.*?)" + str(num).zfill(4) + str(num).zfill(3) + r".ibw",
+                        file.name,
+                    )
 
-            if match is not None:
+            if m is not None:
                 return [file], {}
 
         return None
diff --git a/src/erlab/io/plugins/i05.py b/src/erlab/io/plugins/i05.py
new file mode 100644
index 00000000..ac700b96
--- /dev/null
+++ b/src/erlab/io/plugins/i05.py
@@ -0,0 +1,75 @@
+"""Data loader for beamline I05 at Diamond."""
+
+from typing import ClassVar
+
+import xarray as xr
+
+from erlab.io.dataloader import LoaderBase
+from erlab.io.nexusutils import get_entry, nxgroup_to_xarray
+
+
+class I05Loader(LoaderBase):
+    name = "i05"
+
+    aliases = ("diamond_i05",)
+
+    name_map: ClassVar[dict] = {
+        "eV": "instrument.analyser.energies",
+        "alpha": "instrument.analyser.angles",
+        "beta": [
+            "instrument.deflector_x.deflector_x",
+            "instrument.analyser.deflector_x",
+        ],
+        "delta": "instrument.manipulator.saazimuth",  # azi
+        "chi": "instrument.manipulator.sapolar",  # polar
+        "xi": "instrument.manipulator.satilt",  # tilt
+        "x": "instrument.manipulator.sax",
+        "y": "instrument.manipulator.say",
+        "z": "instrument.manipulator.saz",
+        "hv": "instrument.monochromator.energy",
+        "temp_sample": "sample.temperature",
+    }
+
+    coordinate_attrs = ("beta", "delta", "chi", "xi", "hv", "x", "y", "z")
+    additional_attrs: ClassVar[dict] = {"configuration": 3}
+
+    skip_validate: bool = True
+    always_single: bool = True
+
+    @property
+    def file_dialog_methods(self):
+        return {"Diamond I05 Raw Data (*.nxs)": (self.load, {})}
+
+    def load_single(self, file_path) -> xr.DataArray:
+        out = nxgroup_to_xarray(get_entry(file_path), "analyser/data").squeeze()
+
+        if (
+            "instrument.centre_energy.centre_energy" in out.dims
+            and "instrument.monochromator.energy" in out.coords
+        ):
+            # swap_dims not working for 2D coords, workaround
+            out = (
+                out.assign_coords(
+                    {
+                        "instrument.monochromator.energy": out[
+                            "instrument.centre_energy.centre_energy"
+                        ],
+                        "instrument.centre_energy.centre_energy": out[
+                            "instrument.monochromator.energy"
+                        ],
+                    }
+                )
+                .rename(
+                    {
+                        "instrument.monochromator.energy": "energy_temp",
+                        "instrument.centre_energy.centre_energy": "hv_temp",
+                    }
+                )
+                .rename(
+                    {
+                        "hv_temp": "instrument.monochromator.energy",
+                        "energy_temp": "instrument.centre_energy.centre_energy",
+                    }
+                )
+            )
+        return out
diff --git a/src/erlab/io/plugins/kriss.py b/src/erlab/io/plugins/kriss.py
index 4d6af047..a4f65112 100644
--- a/src/erlab/io/plugins/kriss.py
+++ b/src/erlab/io/plugins/kriss.py
@@ -2,7 +2,6 @@
 
 import os
 import re
-from collections.abc import Sequence
 from typing import ClassVar
 
 import erlab.io.utils
@@ -22,15 +21,14 @@ class KRISSLoader(DA30Loader):
     def name_map(self):
         return super().name_map | {"chi": "ThetaY", "xi": "ThetaX"}
 
-    def identify(
-        self, num: int, data_dir: str | os.PathLike
-    ) -> tuple[list[str], dict[str, Sequence]] | None:
+    def identify(self, num: int, data_dir: str | os.PathLike):
         for file in erlab.io.utils.get_files(data_dir, extensions=(".ibw", ".zip")):
-            if file.endswith(".zip"):
-                match = re.match(r"(.*?)" + str(num).zfill(4) + r".zip", file)
+            if file.suffix == ".zip":
+                match = re.match(r"(.*?)" + str(num).zfill(4) + r".zip", str(file))
             else:
                 match = re.match(
-                    r"(.*?)" + str(num).zfill(4) + str(num).zfill(3) + r".ibw", file
+                    r"(.*?)" + str(num).zfill(4) + str(num).zfill(3) + r".ibw",
+                    str(file),
                 )
 
             if match is not None:
diff --git a/src/erlab/io/plugins/lorea.py b/src/erlab/io/plugins/lorea.py
new file mode 100644
index 00000000..d29100b3
--- /dev/null
+++ b/src/erlab/io/plugins/lorea.py
@@ -0,0 +1,179 @@
+"""Data loader for beamline 20 LOREA at ALBA."""
+
+import pathlib
+import re
+import struct
+from typing import ClassVar
+
+import numpy as np
+import xarray as xr
+
+import erlab.io
+from erlab.io.dataloader import LoaderBase
+from erlab.io.nexusutils import get_entry, nxgroup_to_xarray
+
+
+def _get_data(group):
+    return group.get_default()[group.get_default().signal]
+
+
+class LOREALoader(LoaderBase):
+    name = "lorea"
+
+    aliases = ("alba_bl20",)
+
+    name_map: ClassVar[dict] = {
+        "eV": ["instrument.analyser.energies", "energies"],
+        "alpha": ["instrument.analyser.angles", "angles"],
+        "beta": ["instrument.analyser.defl_angles", "defl_angles"],
+        "delta": "instrument.manipulator.saazimuth",  # azi
+        "chi": "instrument.manipulator.sapolar",  # polar
+        "xi": "instrument.manipulator.satilt",  # tilt
+        "x": "instrument.manipulator.sax",
+        "y": "instrument.manipulator.say",
+        "z": "instrument.manipulator.saz",
+        "hv": "instrument.monochromator.energy",
+        "temp_sample": "sample.temperature",
+    }
+
+    coordinate_attrs = ("beta", "delta", "chi", "xi", "hv", "x", "y", "z")
+    additional_attrs: ClassVar[dict] = {"configuration": 3}
+
+    skip_validate: bool = True
+    always_single: bool = True
+
+    @property
+    def file_dialog_methods(self):
+        return {"ALBA BL20 LOREA Raw Data (*.nxs, *.krx)": (self.load, {})}
+
+    def load_single(self, file_path) -> xr.DataArray:
+        if pathlib.Path(file_path).suffix == ".krx":
+            return self._load_krx(file_path)
+
+        return nxgroup_to_xarray(get_entry(file_path), _get_data)
+
+    def identify(self, num, data_dir, krax=False):
+        file = None
+        if krax:
+            for f in erlab.io.utils.get_files(data_dir, ".krx"):
+                if re.match(rf".+-\d-{str(num).zfill(5)}_\d.krx", f.name):
+                    file = f
+                    break
+        else:
+            for f in erlab.io.utils.get_files(data_dir, ".nxs"):
+                if re.match(rf"{str(num).zfill(3)}_.+.nxs", f.name):
+                    file = f
+                    break
+
+        if file is None:
+            return None
+
+        return [file], {}
+
+    def _load_krx(self, file_path):
+        # Adapted from `load_krax_FS.ipf` Igor procedure by Felix Baumberger
+
+        def read_header(file, pos, size):
+            file.seek(pos)
+            return file.read(size).decode("utf-8")
+
+        def header_to_dict(header_str):
+            out = {}
+            for row in header_str.strip().split("\r\n"):
+                if row.startswith("DATA:"):
+                    break
+                k, v = row.split("\t")
+                try:
+                    v = float(v)
+                except ValueError:
+                    pass
+                else:
+                    if v.is_integer():
+                        v = int(v)
+                out[k] = v
+            return out
+
+        with open(file_path, "rb") as file:
+            file.seek(4)
+            v0 = struct.unpack("<I", file.read(4))[0]
+            is_64bit = v0 == 0
+
+            if is_64bit:
+                fmt = "<Q"
+                num_bytes = 8
+            else:
+                fmt = "<I"
+                num_bytes = 4
+
+            file.seek(0)
+            v1 = struct.unpack(fmt, file.read(num_bytes))[0]
+            n_images = v1 // 3
+
+            # File-position of first image
+            file.seek(num_bytes)
+            image_pos = struct.unpack(fmt, file.read(num_bytes))[0]
+
+            # Parallel detection angle
+            file.seek(2 * num_bytes)
+            image_sizeY = struct.unpack(fmt, file.read(num_bytes))[0]
+
+            # Energy coordinate
+            file.seek(3 * num_bytes)
+            image_sizeX = struct.unpack(fmt, file.read(num_bytes))[0]
+
+            # Autodetect header format and get wave scaling from first header
+            header_pos = (image_pos + image_sizeX * image_sizeY + 1) * 4
+            header = read_header(file, header_pos, 1200)
+            header = header[: header.find("DATA:")]
+            header_dict = header_to_dict(header)
+
+            if header.startswith("Lines\t"):
+                e0 = header_dict.get("Start K.E.", None)
+                e1 = header_dict.get("End K.E.", None)
+                x0 = header_dict.get("ScaleMin", None)
+                x1 = header_dict.get("ScaleMax", None)
+                y0 = header_dict.get("MapStartX", None)
+                y1 = header_dict.get("MapEndX", None)
+            else:
+                # Old header format
+                e0 = header_dict.get("Start K.E.", None)
+                e1 = header_dict.get("End K.E.", None)
+                x0 = header_dict.get("XScaleMin", None)
+                x1 = header_dict.get("XScaleMax", None)
+                y0 = header_dict.get("YScaleMin", None)
+                y1 = header_dict.get("YScaleMax", None)
+
+            data = []
+            for ii in range(n_images):
+                file.seek((ii * 3 + 1) * num_bytes)
+                image_pos = struct.unpack(fmt, file.read(num_bytes))[0]
+
+                file.seek(image_pos * 4)
+                buffer = np.frombuffer(
+                    file.read(image_sizeX * image_sizeY * 4), dtype=np.int32
+                ).reshape((image_sizeY, image_sizeX))
+
+                header_pos = (image_pos + image_sizeX * image_sizeY + 1) * 4
+                file.seek(header_pos)
+                header = file.read(1200).decode("utf-8")
+                header = header[: header.find("DATA:")]
+                header_dict = header_to_dict(header)
+
+                data.append(
+                    xr.DataArray(
+                        buffer,
+                        dims=["angles", "energies"],
+                        coords={
+                            "angles": np.linspace(x0, x1, image_sizeY),
+                            "energies": np.linspace(e0, e1, image_sizeX),
+                        },
+                        attrs=header_dict,
+                    )
+                )
+
+            if n_images == 1:
+                return data[0]
+
+            return xr.concat(data, dim="defl_angles").assign_coords(
+                defl_angles=np.linspace(y0, y1, n_images)
+            )
diff --git a/src/erlab/io/plugins/merlin.py b/src/erlab/io/plugins/merlin.py
index 2329e32b..bf988333 100644
--- a/src/erlab/io/plugins/merlin.py
+++ b/src/erlab/io/plugins/merlin.py
@@ -67,7 +67,7 @@ class MERLINLoader(LoaderBase):
 
     @property
     def file_dialog_methods(self):
-        return {"ALS BL4.0.3 Raw Data (*.pxt, *.ibw)": (self.load, {})}
+        return {"ALS BL4.0.3 Raw Data (*.pxt *.ibw)": (self.load, {})}
 
     def load_single(self, file_path: str | os.PathLike) -> xr.DataArray:
         if os.path.splitext(file_path)[1] == ".ibw":
@@ -173,19 +173,18 @@ def generate_summary(
     ) -> pd.DataFrame:
         files: dict[str, str] = {}
 
-        for path in erlab.io.utils.get_files(data_dir, extensions=(".pxt",)):
-            data_name = os.path.splitext(os.path.basename(path))[0]
+        for pth in erlab.io.utils.get_files(data_dir, extensions=(".pxt",)):
+            data_name = pth.stem
             name_match = re.match(r"(.*?_\d{3})_(?:_S\d{3})?", data_name)
             if name_match is not None:
                 data_name = name_match.group(1)
-            files[data_name] = path
+            files[data_name] = str(pth)
 
         if not exclude_live:
             for file in os.listdir(data_dir):
                 if file.endswith(".ibw"):
                     data_name = os.path.splitext(file)[0]
-                    path = os.path.join(data_dir, file)
-                    files[data_name] = path
+                    files[data_name] = os.path.join(data_dir, file)
 
         summary_attrs: dict[str, str] = {
             "Lens Mode": "Lens Mode",
diff --git a/src/erlab/io/plugins/ssrl52.py b/src/erlab/io/plugins/ssrl52.py
index 43a88603..fe503935 100644
--- a/src/erlab/io/plugins/ssrl52.py
+++ b/src/erlab/io/plugins/ssrl52.py
@@ -248,10 +248,10 @@ def identify(
                 data_dir, extensions=(".h5",), notcontains="zap"
             )
 
-        for file in target_files:
-            match = re.match(r"(.*?)_" + str(num).zfill(4) + r".h5", file)
+        for path in target_files:
+            match = re.match(r"(.*?)_" + str(num).zfill(4) + r".h5", path.name)
             if match is not None:
-                return [file], {}
+                return [path], {}
 
         return None
 
@@ -270,9 +270,8 @@ def generate_summary(
         else:
             target_files = erlab.io.utils.get_files(data_dir, extensions=(".h5",))
 
-        for pth in target_files:
-            base_name = os.path.splitext(os.path.basename(pth))[0]
-            files[base_name] = pth
+        for path in target_files:
+            files[path.stem] = str(path)
 
         summary_attrs: dict[str, str] = {
             "Type": "Description",
diff --git a/src/erlab/io/utils.py b/src/erlab/io/utils.py
index e72e8417..317c7c8f 100644
--- a/src/erlab/io/utils.py
+++ b/src/erlab/io/utils.py
@@ -41,11 +41,11 @@ def showfitsinfo(path: str | os.PathLike) -> None:
 
 
 def get_files(
-    directory,
+    directory: str | os.PathLike,
     extensions: Sequence[str] | str | None = None,
     contains: str | None = None,
     notcontains: str | None = None,
-) -> set[str]:
+) -> set[pathlib.Path]:
     """Return file names in a directory with the given extension(s).
 
     Directories are ignored.
@@ -63,11 +63,11 @@ def get_files(
 
     Returns
     -------
-    files : set of str
-        Set of file names in the directory.
+    files : set of pathlib.Path
+        Set of file objects in the directory.
 
     """
-    files: set[str] = set()
+    files: set[pathlib.Path] = set()
 
     if isinstance(extensions, str):
         extensions = [extensions]
@@ -84,7 +84,7 @@ def get_files(
         ):
             continue
 
-        files.add(str(f))
+        files.add(f)
 
     return files
 
diff --git a/src/erlab/plotting/atoms.py b/src/erlab/plotting/atoms.py
index d1c144f4..dce6b823 100644
--- a/src/erlab/plotting/atoms.py
+++ b/src/erlab/plotting/atoms.py
@@ -10,7 +10,7 @@
 import contextlib
 import functools
 import itertools
-from collections.abc import Callable, Iterable, Mapping, Sequence
+from collections.abc import Callable, Iterable, Mapping
 from typing import Literal, cast
 
 import matplotlib.collections
@@ -42,7 +42,7 @@ def sunflower_sphere(n: int = 100) -> npt.NDArray[np.float64]:
 
 
 def projected_length(
-    ax: mpl_toolkits.mplot3d.Axes3D, length: np.float64 | Sequence[np.float64]
+    ax: mpl_toolkits.mplot3d.Axes3D, length: float | npt.NDArray[np.float64]
 ):
     if np.iterable(length):
         return np.asarray([projected_length(ax, d) for d in np.asarray(length).flat])
@@ -218,7 +218,9 @@ def __init__(self, segments, *, scale_linewidths: bool = True, **kwargs) -> None
 
     def set_linewidth(self, lw) -> None:
         super().set_linewidth(lw)
-        self.linewidths_orig = np.asarray(lw) if np.iterable(lw) else lw
+        self.linewidths_orig: npt.NDArray[np.float64] | float = (
+            np.asarray(lw, dtype=np.float64) if np.iterable(lw) else float(lw)
+        )
 
     def draw(self, renderer) -> None:
         if self._scale_linewidths:  # and np.isfinite(self.axes._focal_length):


Test patch:
diff --git a/tests/accessors/test_general.py b/tests/accessors/test_general.py
index 5956b3cb..7027a4e1 100644
--- a/tests/accessors/test_general.py
+++ b/tests/accessors/test_general.py
@@ -56,3 +56,86 @@ def test_ds_qshow_fit(plot_components: bool):
     assert isinstance(
         result_ds.qshow(plot_components=plot_components), panel.layout.base.Column
     )
+
+
+@pytest.mark.parametrize(
+    ("indexers", "expected_shape"),
+    [
+        ({"x": 2.0}, (5,)),
+        ({"x": slice(1.0, 3.0)}, (3, 5)),
+        ({"x": 2.0, "x_width": 1.0}, (5,)),
+        ({"x": slice(1.0, 3.0), "y": 2.0}, (3,)),
+    ],
+)
+def test_qsel_shape(indexers, expected_shape):
+    dat = xr.DataArray(np.arange(25).reshape(5, 5), dims=("x", "y"))
+    result = dat.qsel(indexers)
+    assert result.shape == expected_shape
+
+
+def test_qsel_verbose(capfd):
+    dat = xr.DataArray(np.arange(25).reshape(5, 5), dims=("x", "y"))
+    dat.qsel({"x": 2.0, "x_width": 1.0}, verbose=True)
+    captured = capfd.readouterr()
+    assert (
+        "Selected data with {'x': slice(1.5, 2.5, None)}, averaging over ['x']"
+        in captured.out
+    )
+
+
+def test_qsel_invalid_dimension():
+    dat = xr.DataArray(np.arange(25).reshape(5, 5), dims=("x", "y"))
+    with pytest.raises(ValueError, match="Dimension `z` not found in data"):
+        dat.qsel({"z": 2.0})
+
+
+def test_qsel_slice_with_width():
+    dat = xr.DataArray(np.arange(25).reshape(5, 5), dims=("x", "y"))
+    with pytest.raises(
+        ValueError,
+        match="Slice not allowed for value of dimension `x` with width specified",
+    ):
+        dat.qsel({"x": slice(1.0, 3.0), "x_width": 1.0})
+
+
+def test_qsel_value_outside_bounds():
+    dat = xr.DataArray(
+        np.arange(25).reshape(5, 5),
+        dims=("x", "y"),
+        coords={"x": np.arange(5), "y": np.arange(5)},
+    )
+    with pytest.warns(
+        UserWarning, match="Selected value 10.0 for `x` is outside coordinate bounds"
+    ):
+        dat.qsel({"x": 10.0})
+
+
+def test_qsel_drop_unindexed_dims():
+    dat = xr.DataArray(np.arange(25).reshape(5, 5), dims=("x", "y"))
+    result = dat.qsel({"x": 2.0})
+
+    xr.testing.assert_equal(result, dat.isel(x=2))
+
+
+def test_qsel_around():
+    dat = xr.DataArray(np.arange(25).reshape(5, 5), dims=("x", "y"))
+    result = dat.qsel.around(radius=2.0, x=2.0, y=2.0, average=False)
+    xr.testing.assert_identical(
+        result,
+        xr.DataArray(
+            np.array(
+                [
+                    [np.nan, np.nan, 2.0, np.nan, np.nan],
+                    [np.nan, 6.0, 7.0, 8.0, np.nan],
+                    [10.0, 11.0, 12.0, 13.0, 14.0],
+                    [np.nan, 16.0, 17.0, 18.0, np.nan],
+                    [np.nan, np.nan, 22.0, np.nan, np.nan],
+                ]
+            ),
+            dims=("x", "y"),
+        ),
+    )
+
+    xr.testing.assert_identical(
+        dat.qsel.around(radius=2.0, x=2.0, y=2.0, average=True), xr.DataArray(12.0)
+    )
diff --git a/tests/conftest.py b/tests/conftest.py
index a60058e1..5ef724de 100644
--- a/tests/conftest.py
+++ b/tests/conftest.py
@@ -9,10 +9,10 @@
 
 from erlab.io.exampledata import generate_data_angles, generate_gold_edge
 
-DATA_COMMIT_HASH = "4503bd2a789d5182c2aee434b53bf12ae3925974"
+DATA_COMMIT_HASH = "985457233a02eace67da27a0d5450fc0301e55c3"
 """The commit hash of the commit to retrieve from `kmnhan/erlabpy-data`."""
 
-DATA_KNOWN_HASH = "7426d2654fe9d00ea00d9f11dc346180ebc5e14353d7adbc756b6c2e73713c32"
+DATA_KNOWN_HASH = "691d69af64feaa92f928365050392f08c640d3f11a04243c9bdf323c1bc3f0b0"
 """The hash of the `.tar.gz` file."""
 
 
diff --git a/tests/interactive/test_imagetool.py b/tests/interactive/test_imagetool.py
index a09ec8d0..d8cc222d 100644
--- a/tests/interactive/test_imagetool.py
+++ b/tests/interactive/test_imagetool.py
@@ -11,12 +11,15 @@
 from qtpy import QtCore, QtWidgets
 
 import erlab.analysis.transform
-from erlab.interactive.imagetool import itool
+from erlab.interactive.imagetool import ImageTool, _parse_input, itool
 from erlab.interactive.imagetool.manager import ImageToolManager
 
 
 def accept_dialog(
-    dialog_trigger: Callable, time_out: int = 5, pre_call: Callable | None = None
+    dialog_trigger: Callable,
+    time_out: int = 5,
+    pre_call: Callable | None = None,
+    accept_call: Callable | None = None,
 ) -> QtWidgets.QDialog:
     """Accept a dialog during testing.
 
@@ -33,6 +36,8 @@ def accept_dialog(
     pre_call
         Callable that takes the dialog as a single argument. If provided, it is executed
         before calling ``.accept()`` on the dialog.
+    accept_call
+        If provided, it is called instead of ``.accept()`` on the dialog.
     """
     dialog = None
     start_time = time.time()
@@ -48,7 +53,10 @@ def dialog_creation():
         if isinstance(dialog, QtWidgets.QDialog):
             if pre_call is not None:
                 pre_call(dialog)
-            if (
+
+            if accept_call is not None:
+                accept_call(dialog)
+            elif (
                 isinstance(dialog, QtWidgets.QMessageBox)
                 and dialog.defaultButton() is not None
             ):
@@ -112,6 +120,7 @@ def test_itool_dtypes(qtbot, val_dtype, coord_dtype):
     with qtbot.waitExposed(win):
         win.show()
         win.activateWindow()
+    time.sleep(0.5)
 
     move_and_compare_values(qtbot, win, [12.0, 7.0, 6.0, 11.0])
     win.close()
@@ -176,7 +185,7 @@ def _go_to_file(dialog: QtWidgets.QFileDialog):
 
 def test_itool(qtbot):
     data = xr.DataArray(np.arange(25).reshape((5, 5)), dims=["x", "y"])
-    win = itool(data, execute=False)
+    win = itool(data, execute=False, cmap="terrain_r")
     qtbot.addWidget(win)
 
     with qtbot.waitExposed(win):
@@ -303,21 +312,22 @@ def test_itool_tools(qtbot, gold):
     del win
 
 
-def test_itool_ds(qtbot):
+def test_parse_input():
     # If no 2D to 4D data is present in given Dataset, ValueError is raised
     with pytest.raises(
         ValueError, match="No valid data for ImageTool found in the Dataset"
     ):
-        itool(
+        _parse_input(
             xr.Dataset(
                 {
                     "data1d": xr.DataArray(np.arange(5), dims=["x"]),
                     "data0d": 1,
                 }
-            ),
-            execute=False,
+            )
         )
 
+
+def test_itool_ds(qtbot):
     data = xr.Dataset(
         {
             "data1d": xr.DataArray(np.arange(5), dims=["x"]),
@@ -342,6 +352,7 @@ def test_itool_ds(qtbot):
 
     # Check if properly linked
     assert wins[0].slicer_area._linking_proxy == wins[1].slicer_area._linking_proxy
+    assert wins[0].slicer_area.linked_slicers == {wins[1].slicer_area}
 
     wins[0].slicer_area.unlink()
     wins[1].slicer_area.unlink()
@@ -351,15 +362,49 @@ def test_itool_ds(qtbot):
     del wins
 
 
+def test_itool_multidimensional(qtbot):
+    win = itool(
+        xr.DataArray(np.arange(25).reshape((5, 5)), dims=["x", "y"]), execute=False
+    )
+    qtbot.addWidget(win)
+
+    with qtbot.waitExposed(win):
+        win.show()
+        win.activateWindow()
+
+    win.slicer_area.set_data(
+        xr.DataArray(np.arange(125).reshape((5, 5, 5)), dims=["x", "y", "z"])
+    )
+    move_and_compare_values(qtbot, win, [62.0, 37.0, 32.0, 57.0])
+
+    win.slicer_area.set_data(
+        xr.DataArray(np.arange(625).reshape((5, 5, 5, 5)), dims=["x", "y", "z", "t"])
+    )
+    move_and_compare_values(qtbot, win, [312.0, 187.0, 162.0, 287.0])
+
+    win.close()
+
+
 def test_value_update(qtbot):
     win = itool(
         xr.DataArray(np.arange(25).reshape((5, 5)), dims=["x", "y"]), execute=False
     )
+    qtbot.addWidget(win)
+
+    with qtbot.waitExposed(win):
+        win.show()
+        win.activateWindow()
 
     new_vals = -np.arange(25).reshape((5, 5)).astype(float)
     win.slicer_area.update_values(new_vals)
     assert_almost_equal(win.array_slicer.point_value(0), -12.0)
 
+    win.close()
+
+
+def test_value_update_errors():
+    win = ImageTool(xr.DataArray(np.arange(25).reshape((5, 5)), dims=["x", "y"]))
+
     with pytest.raises(ValueError, match="DataArray dimensions do not match"):
         win.slicer_area.update_values(
             xr.DataArray(np.arange(24).reshape((2, 2, 6)), dims=["x", "y", "z"])
@@ -375,8 +420,6 @@ def test_value_update(qtbot):
     with pytest.raises(ValueError, match="^Data shape does not match.*"):
         win.slicer_area.update_values(np.arange(24).reshape((4, 6)))
 
-    win.close()
-
 
 def test_sync(qtbot):
     manager = ImageToolManager()
@@ -542,3 +585,59 @@ def _set_dialog_params(dialog):
     assert not win.slicer_area.main_image.is_guidelines_visible
 
     win.close()
+
+
+def normalize(data, norm_dims, option):
+    area = data.mean(norm_dims)
+    minimum = data.min(norm_dims)
+    maximum = data.max(norm_dims)
+
+    match option:
+        case 0:
+            return data / area
+        case 1:
+            return (data - minimum) / (maximum - minimum)
+        case 2:
+            return data - minimum
+        case _:
+            return (data - minimum) / area
+
+
+@pytest.mark.parametrize("option", [0, 1, 2, 3])
+def test_itool_normalize(qtbot, option):
+    data = xr.DataArray(
+        np.arange(25).reshape((5, 5)).astype(float),
+        dims=["x", "y"],
+        coords={"x": np.arange(5), "y": np.arange(5)},
+    )
+    win = itool(data, execute=False)
+    qtbot.addWidget(win)
+
+    # Test dialog
+    def _set_dialog_params(dialog):
+        dialog.dim_checks["x"].setChecked(True)
+        dialog.opts[option].setChecked(True)
+
+        # Preview
+        dialog.preview_button.click()
+
+    accept_dialog(win.mnb._normalize, pre_call=_set_dialog_params)
+
+    # Check if the data is normalized
+    xarray.testing.assert_identical(
+        win.slicer_area.data, normalize(data, ("x",), option)
+    )
+
+    # Reset normalization
+    win.slicer_area.apply_func(None)
+    xarray.testing.assert_identical(win.slicer_area.data, data)
+
+    # Check if canceling the dialog does not change the data
+    accept_dialog(
+        win.mnb._normalize,
+        pre_call=_set_dialog_params,
+        accept_call=lambda d: d.reject(),
+    )
+    xarray.testing.assert_identical(win.slicer_area.data, data)
+
+    win.close()
diff --git a/tests/io/plugins/test_i05.py b/tests/io/plugins/test_i05.py
new file mode 100644
index 00000000..5218af05
--- /dev/null
+++ b/tests/io/plugins/test_i05.py
@@ -0,0 +1,24 @@
+import pytest
+import xarray as xr
+
+import erlab.io
+
+
+@pytest.fixture(scope="module")
+def data_dir(test_data_dir):
+    erlab.io.set_loader("i05")
+    erlab.io.set_data_dir(test_data_dir / "i05")
+    return test_data_dir / "i05"
+
+
+@pytest.fixture(scope="module")
+def expected_dir(data_dir):
+    return data_dir / "expected"
+
+
+@pytest.mark.parametrize("fname", ["core", "cut", "fs", "kz"])
+def test_load(expected_dir, fname):
+    xr.testing.assert_identical(
+        erlab.io.load(f"{fname}.nxs"),
+        xr.load_dataarray(expected_dir / f"{fname}.h5", engine="h5netcdf"),
+    )
diff --git a/tests/io/plugins/test_lorea.py b/tests/io/plugins/test_lorea.py
new file mode 100644
index 00000000..2f636749
--- /dev/null
+++ b/tests/io/plugins/test_lorea.py
@@ -0,0 +1,35 @@
+import pytest
+import xarray as xr
+
+import erlab.io
+
+
+@pytest.fixture(scope="module")
+def data_dir(test_data_dir):
+    erlab.io.set_loader("lorea")
+    erlab.io.set_data_dir(test_data_dir / "lorea")
+    return test_data_dir / "lorea"
+
+
+@pytest.fixture(scope="module")
+def expected_dir(data_dir):
+    return data_dir / "expected"
+
+
+@pytest.mark.parametrize(
+    ("args", "expected"),
+    [
+        ({"identifier": 3, "krax": True}, "sample-1-00003_0.h5"),
+        ({"identifier": 1}, "001_sample.h5"),
+        ({"identifier": 2}, "002_sample.h5"),
+        ("001_sample.nxs", "001_sample.h5"),
+        ("002_sample.nxs", "002_sample.h5"),
+        ("sample-1-00003_0.krx", "sample-1-00003_0.h5"),
+    ],
+)
+def test_load(expected_dir, args, expected):
+    loaded = erlab.io.load(**args) if isinstance(args, dict) else erlab.io.load(args)
+
+    xr.testing.assert_identical(
+        loaded, xr.load_dataarray(expected_dir / expected, engine="h5netcdf")
+    )
diff --git a/tests/io/test_dataloader.py b/tests/io/test_dataloader.py
index 398a9da4..77e63d80 100644
--- a/tests/io/test_dataloader.py
+++ b/tests/io/test_dataloader.py
@@ -184,15 +184,10 @@ def generate_summary(self, data_dir):
             # Get all valid data files in directory
             files = {}
             for path in erlab.io.utils.get_files(data_dir, extensions=[".h5"]):
-                # Base name
-                data_name = os.path.splitext(os.path.basename(path))[0]
-
                 # If multiple scans, strip the _S### part
-                name_match = re.match(r"(.*?_\d{3})_(?:_S\d{3})?", data_name)
-                if name_match is not None:
-                    data_name = name_match.group(1)
-
-                files[data_name] = path
+                name_match = re.match(r"(.*?_\d{3})_(?:_S\d{3})?", path.stem)
+                data_name = path.stem if name_match is None else name_match.group(1)
+                files[data_name] = str(path)
 
             # Map dataframe column names to data attributes
             attrs_mapping = {