Merge remote-tracking branch 'origin/main' into removals

.azure-pipelines/ci-conda-env.txt

@@ -1,18 +1,18 @@
 conda-forge::boost
 conda-forge::boost-cpp
 conda-forge::bzip2
-conda-forge::c-compiler
+conda-forge::c-compiler<1.5
 conda-forge::cctbx-base==2022.2
 conda-forge::conda
-conda-forge::cxx-compiler
+conda-forge::cxx-compiler<1.5
 conda-forge::python-dateutil
 conda-forge::dials-data
 conda-forge::docutils
 conda-forge::eigen
 conda-forge::future
 conda-forge::h5py>=3.1
 conda-forge::hdf5plugin
-conda-forge::hdf5=1.12
+conda-forge::hdf5>=1.12,<1.13
 conda-forge::jpeg
 conda-forge::matplotlib-base
 conda-forge::mrcfile

CHANGELOG.rst

@@ -1,3 +1,12 @@
+DIALS 3.11.1 (2022-09-02)
+=========================
+
+Bugfixes
+--------
+
+- ``dxtbx.dlsnxs2cbf``: Fix bug introduced by #572. (`#545 </~https://github.com/cctbx/dxtbx/issues/545>`_)
+
+
 dxtbx 3.11.0 (2022-08-24)
 =========================
 

newsfragments/538.feature

@@ -0,0 +1 @@
+Add spectrum support to FormatNXmx

newsfragments/541.misc

@@ -0,0 +1 @@
+Add dataset for testing the SwissFEL JF16M detector using NeXus

newsfragments/546.bugfix

@@ -0,0 +1 @@
+``NXmx``: Cope more gracefully with scalar NXtransformations values.

newsfragments/547.misc

@@ -0,0 +1 @@
+Fix ``test_get_dependency_chain``: Minor change to output as a result of #546

newsfragments/548.bugfix

@@ -0,0 +1 @@
+``dxtbx.dlsnxs2cbf``: Fix distance and pixel size bugs.

src/dxtbx/format/FormatNXmx.py

@@ -9,8 +9,6 @@
 class FormatNXmx(FormatNexus):
     """Read NXmx-flavour NeXus-format HDF5 data."""
 
-    _cached_file_handle = None
-
     @staticmethod
     def understand(image_file):
         with h5py.File(image_file) as handle:
@@ -30,49 +28,47 @@ def __init__(self, image_file, **kwargs):
     def _start(self):
         self._static_mask = None
 
-        with h5py.File(self._image_file, swmr=True) as fh:
-            nxmx = self._get_nxmx(fh)
-            nxentry = nxmx.entries[0]
-            nxsample = nxentry.samples[0]
-            nxinstrument = nxentry.instruments[0]
-            nxdetector = nxinstrument.detectors[0]
-            nxbeam = nxinstrument.beams[0]
-            self._goniometer_model = dxtbx.nexus.get_dxtbx_goniometer(nxsample)
-            self._beam_factory = dxtbx.nexus.CachedWavelengthBeamFactory(nxbeam)
-            wavelength = self._beam_factory.make_beam(index=0).get_wavelength()
-            self._detector_model = dxtbx.nexus.get_dxtbx_detector(
-                nxdetector, wavelength
-            )
-            self._scan_model = dxtbx.nexus.get_dxtbx_scan(nxsample, nxdetector)
-            self._static_mask = dxtbx.nexus.get_static_mask(nxdetector)
-            self._bit_depth_readout = nxdetector.bit_depth_readout
-
-            if self._scan_model:
-                self._num_images = len(self._scan_model)
+        self._cached_file_handle = h5py.File(self._image_file, swmr=True)
+        nxmx = self._get_nxmx(self._cached_file_handle)
+        nxentry = nxmx.entries[0]
+        nxsample = nxentry.samples[0]
+        nxinstrument = nxentry.instruments[0]
+        nxdetector = nxinstrument.detectors[0]
+        nxbeam = nxinstrument.beams[0]
+        self._goniometer_model = dxtbx.nexus.get_dxtbx_goniometer(nxsample)
+        self._beam_factory = dxtbx.nexus.CachedWavelengthBeamFactory(nxbeam)
+        wavelength = self._beam_factory.make_beam(index=0).get_wavelength()
+        self._detector_model = dxtbx.nexus.get_dxtbx_detector(nxdetector, wavelength)
+        self._scan_model = dxtbx.nexus.get_dxtbx_scan(nxsample, nxdetector)
+        self._static_mask = dxtbx.nexus.get_static_mask(nxdetector)
+        self._bit_depth_readout = nxdetector.bit_depth_readout
+
+        if self._scan_model:
+            self._num_images = len(self._scan_model)
+        else:
+            nxdata = nxmx.entries[0].data[0]
+            if nxdata.signal:
+                data = nxdata[nxdata.signal]
             else:
-                nxdata = nxmx.entries[0].data[0]
-                if nxdata.signal:
-                    data = nxdata[nxdata.signal]
-                else:
-                    data = list(nxdata.values())[0]
-                self._num_images, *_ = data.shape
+                data = list(nxdata.values())[0]
+            self._num_images, *_ = data.shape
 
     def _get_nxmx(self, fh: h5py.File):
         return dxtbx.nexus.nxmx.NXmx(fh)
 
     def _beam(self, index: int | None = None) -> dxtbx.model.Beam:
         return self._beam_factory.make_beam(index=index or 0)
 
+    def get_spectrum(self, index=None):
+        return self._beam_factory.make_spectrum(index=index or 0)
+
     def get_num_images(self) -> int:
         return self._num_images
 
     def get_static_mask(self, index=None, goniometer=None):
         return self._static_mask
 
     def get_raw_data(self, index):
-        if self._cached_file_handle is None:
-            self._cached_file_handle = h5py.File(self._image_file, swmr=True)
-
         nxmx = self._get_nxmx(self._cached_file_handle)
         nxdata = nxmx.entries[0].data[0]
         nxdetector = nxmx.entries[0].instruments[0].detectors[0]

src/dxtbx/format/nexus.py

@@ -566,7 +566,10 @@ def load_model(self, index=None):
         primary_key = "incident_wavelength"
         wavelength = self.obj.handle[primary_key]
         spectrum_wavelengths = wavelength
-        spectrum_weights = self.obj.handle.get(primary_key + "_weight")
+        spectrum_weights = self.obj.handle.get(primary_key + "_weights")
+        if spectrum_weights is None:
+            # Handle deprecation: /~https://github.com/nexusformat/definitions/issues/837
+            spectrum_weights = self.obj.handle.get(primary_key + "_weight")
 
         # If the wavelength array does not represent spectra, look for spectra
         # in the variant chain
@@ -576,7 +579,10 @@ def load_model(self, index=None):
             if "variant" in variant_test.attrs:
                 variant_key = variant_test.attrs["variant"]
                 variant_wavelengths = self.obj.handle[variant_key]
-                variant_weights = self.obj.handle.get(variant_key + "_weight")
+                variant_weights = self.obj.handle.get(variant_key + "_weights")
+                if variant_weights is None:
+                    # Handle deprecation: /~https://github.com/nexusformat/definitions/issues/837
+                    variant_weights = self.obj.handle.get(variant_key + "_weight")
                 if variant_weights is None:
                     variant_test = variant_wavelengths  # Keep looking
                 else:
@@ -635,7 +641,7 @@ def get_wavelength(wavelength):
         return self.model
 
 
-def get_change_of_basis(transformation):
+def get_change_of_basis(transformation, include_setting=True):
     """
     Get the 4x4 homogenous coordinate matrix for a given NXtransformation.
     """
@@ -645,7 +651,7 @@ def get_change_of_basis(transformation):
     axis_type = h5str(transformation.attrs["transformation_type"])
 
     vector = n2i_cob * col(transformation.attrs["vector"]).normalize()
-    setting = transformation[0]
+    setting = transformation[0] if include_setting else 0
     units = h5str(transformation.attrs["units"])
 
     if "offset" in transformation.attrs:
@@ -924,7 +930,11 @@ def set_frame(pg, transformation):
                 fast = matrix.col([-fast[0], fast[1], -fast[2]])
                 slow = slow_pixel_direction_handle.attrs["vector"]
                 slow = matrix.col([-slow[0], slow[1], -slow[2]])
-                parent, cob = get_cumulative_change_of_basis(depends_on)
+                cob = get_change_of_basis(
+                    fast_pixel_direction_handle, include_setting=False
+                ) * get_change_of_basis(
+                    slow_pixel_direction_handle, include_setting=False
+                )
                 origin = matrix.col((cob * matrix.col((0, 0, 0, 1)))[0:3])
 
                 p.set_local_frame(fast.elems, slow.elems, origin.elems)

src/dxtbx/nexus/__init__.py

@@ -6,11 +6,13 @@
 import h5py
 import numpy as np
 
+import cctbx
 from cctbx import eltbx
 from scitbx.array_family import flex
 
 import dxtbx.model
 from dxtbx import flumpy
+from dxtbx.nexus.nxmx import units
 
 from . import nxmx
 
@@ -65,18 +67,105 @@ class CachedWavelengthBeamFactory:
     """Defer Beam generation whilst caching the wavelength value"""
 
     def __init__(self, nxbeam: nxmx.NXbeam):
-        self.incident_wavelength = nxbeam.incident_wavelength.to("angstrom").magnitude
+        self.handle = nxbeam._handle
+        self.index = None
+        self.model = None
+        self.spectrum = None
 
     def make_beam(self, index: int = 0) -> dxtbx.model.Beam:
-        if np.isscalar(self.incident_wavelength):
-            wavelength = self.incident_wavelength
+        self.read_models(index)
+        return self.model
+
+    def make_spectrum(self, index: int = 0) -> dxtbx.model.Spectrum:
+        self.read_models(index)
+        return self.spectrum
+
+    def read_models(self, index: int = 0):
+        # Cached model
+        if self.model is not None and index == self.index:
+            return
+
+        # Get the items from the NXbeam class
+        # Note it would be better if there were a general way to read weights
+        # and variants from the nxmx classes
+        # See /~https://github.com/cctbx/dxtbx/issues/549
+        primary_key = "incident_wavelength"
+        wavelength = self.handle[primary_key]
+        spectrum_wavelengths = wavelength
+        spectrum_weights = self.handle.get(primary_key + "_weights")
+        if spectrum_weights is None:
+            # Handle deprecation: /~https://github.com/nexusformat/definitions/issues/837
+            spectrum_weights = self.handle.get(primary_key + "_weight")
+
+        # If the wavelength array does not represent spectra, look for spectra
+        # in the variant chain
+        variant_test = wavelength
+        has_variant_spectra = False
+        while spectrum_weights is None:
+            if "variant" in variant_test.attrs:
+                variant_key = variant_test.attrs["variant"]
+                variant_wavelengths = self.handle[variant_key]
+                variant_weights = self.handle.get(variant_key + "_weights")
+                if variant_weights is None:
+                    # Handle deprecation: /~https://github.com/nexusformat/definitions/issues/837
+                    variant_weights = self.handle.get(variant_key + "_weight")
+                if variant_weights is None:
+                    variant_test = variant_wavelengths  # Keep looking
+                else:
+                    # Found spectra
+                    spectrum_wavelengths = variant_wavelengths
+                    spectrum_weights = variant_weights  # cause while loop to end
+                    has_variant_spectra = True
+            else:
+                break
+
+        if index is None:
+            index = 0
+        self.index = index
+
+        def get_wavelength(wavelength):
+            if wavelength.shape in ((), (1,)):
+                wavelength_value = wavelength[()]
+            else:
+                wavelength_value = wavelength[index]
+            wavelength_units = units(wavelength)
+            wavelength_value = float(
+                (wavelength_value * wavelength_units).to("angstrom").magnitude
+            )
+            return wavelength_value
+
+        if spectrum_weights is None:
+            # Construct the beam model
+            wavelength_value = get_wavelength(wavelength)
+            self.model = dxtbx.model.Beam(
+                direction=(0, 0, 1), wavelength=wavelength_value
+            )
         else:
-            assert len(self.incident_wavelength) > index, len(self.incident_wavelength)
-            wavelength = self.incident_wavelength[index]
-        return dxtbx.model.BeamFactory.make_beam(
-            sample_to_source=(0, 0, 1),
-            wavelength=wavelength,
-        )
+            self.model = dxtbx.model.Beam()
+            self.model.set_direction((0, 0, 1))
+
+            wavelength_units = units(spectrum_wavelengths)
+
+            if len(spectrum_wavelengths.shape) > 1:
+                spectrum_wavelengths = spectrum_wavelengths[index]
+            else:
+                spectrum_wavelengths = spectrum_wavelengths[()]
+            if len(spectrum_weights.shape) > 1:
+                spectrum_weights = spectrum_weights[index]
+            else:
+                spectrum_weights = spectrum_weights[()]
+
+            spectrum_wavelengths = (
+                (spectrum_wavelengths * wavelength_units).to("angstrom").magnitude
+            )
+            spectrum_energies = cctbx.factor_ev_angstrom / spectrum_wavelengths
+            self.spectrum = dxtbx.model.Spectrum(spectrum_energies, spectrum_weights)
+
+            if has_variant_spectra:
+                wavelength_value = get_wavelength(wavelength)
+                self.model.set_wavelength(wavelength_value)
+            else:
+                self.model.set_wavelength(self.spectrum.get_weighted_wavelength())
 
 
 def get_dxtbx_scan(
@@ -216,7 +305,14 @@ def get_dxtbx_detector(
             # Hierarchical detector model
             fast_axis = MCSTAS_TO_IMGCIF @ module.fast_pixel_direction.vector
             slow_axis = MCSTAS_TO_IMGCIF @ module.slow_pixel_direction.vector
-            origin = np.array((0.0, 0.0, 0.0))
+            origin = MCSTAS_TO_IMGCIF @ (
+                module.fast_pixel_direction.offset.to("mm").magnitude
+                if module.fast_pixel_direction.offset is not None
+                else np.array([0.0, 0.0, 0.0])
+                + module.slow_pixel_direction.offset.to("mm").magnitude
+                if module.slow_pixel_direction.offset is not None
+                else np.array([0.0, 0.0, 0.0])
+            )
         else:
             # Flat detector model
 
@@ -338,7 +434,10 @@ def get_static_mask(nxdetector: nxmx.NXdetector) -> tuple[flex.bool, ...] | None
     if pixel_mask is None or not pixel_mask.size or pixel_mask.ndim != 2:
         return None
     all_slices = get_detector_module_slices(nxdetector)
-    return tuple(flumpy.from_numpy(pixel_mask[slices]) == 0 for slices in all_slices)
+    return tuple(
+        flumpy.from_numpy(np.ascontiguousarray(pixel_mask[slices])) == 0
+        for slices in all_slices
+    )
 
 
 def _dataset_as_flex(

src/dxtbx/nexus/nxmx.py

@@ -6,7 +6,7 @@
 import operator
 from collections import abc, namedtuple
 from functools import cached_property, reduce
-from typing import Iterable, Iterator, Sequence, Union, cast, overload
+from typing import Iterable, Iterator, Sequence, Union, overload
 
 import dateutil.parser
 import h5py
@@ -433,7 +433,7 @@ def depends_on(self) -> NXtransformationsAxis | None:
         return None
 
     def __getitem__(self, key) -> pint.Quantity:
-        return self._handle[key] * self.units
+        return np.atleast_1d(self._handle)[key] * self.units
 
     @cached_property
     def end(self) -> NXNumber | None:
@@ -451,8 +451,7 @@ def increment_set(self) -> pint.Quantity | None:
 
     @cached_property
     def matrix(self) -> np.ndarray:
-
-        values = cast(pint.Quantity, np.atleast_1d(self[()]))
+        values = self[()]
         if np.any(values):
             values = (
                 values.to("mm").magnitude
@@ -1166,8 +1165,7 @@ def get_rotation_axes(dependency_chain: DependencyChain) -> Axes:
     for transformation in dependency_chain:
         if transformation.transformation_type != "rotation":
             continue
-        values = cast(pint.Quantity, np.atleast_1d(transformation[()]))
-        values = values.to("degrees").magnitude
+        values = transformation[()].to("degrees").magnitude
         is_scan = len(values) > 1 and not np.all(values == values[0])
         axes.append(transformation.vector)
         angles.append(values[0])