Merge pull request #20 from scipp/transformation-chain-unit-fix

Fix and improve unit handling in time-dependent transformation chains

conda/meta.yaml

@@ -11,6 +11,7 @@ requirements:
     - python>=3.8
     - python-dateutil
     - scipp
+    - scipy
     - h5py
 
 test:

setup.cfg

@@ -25,6 +25,7 @@ packages = find:
 install_requires =
     python-dateutil
     scipp>=0.12
+    scipy  # we use scipp.interpolate which depends on this
     h5py
 python_requires = >=3.8
 include_package_data = True

src/scippnexus/nxobject.py

@@ -187,13 +187,13 @@ def __getitem__(self, select) -> sc.Variable:
             self._dataset.read_direct(variable.values, source_sel=index)
         else:
             variable.values = self._dataset[index]
-        if self._is_time or _is_time(variable):
+        if _is_time(variable):
             starts = []
             for name in self.attrs:
                 if (dt := _as_datetime(self.attrs[name])) is not None:
                     starts.append(dt)
             if self._is_time and len(starts) == 0:
-                starts.append(sc.epoch(unit='ns'))
+                starts.append(sc.epoch(unit=self.unit))
             if len(starts) == 1:
                 variable = convert_time_to_datetime64(
                     variable,

src/scippnexus/nxtransformations.py

@@ -90,6 +90,16 @@ def _interpolate_transform(transform, xnew):
                                    fill_value="extrapolate")(xnew=xnew)
 
 
+def _smaller_unit(a, b):
+    if a.unit == b.unit:
+        return a.unit
+    ratio = sc.scalar(1.0, unit=a.unit).to(unit=b.unit)
+    if ratio.value < 1.0:
+        return a.unit
+    else:
+        return b.unit
+
+
 def get_full_transformation(depends_on: Field) -> Union[None, sc.DataArray]:
     """
     Get the 4x4 transformation matrix for a component, resulting
@@ -105,12 +115,15 @@ def get_full_transformation(depends_on: Field) -> Union[None, sc.DataArray]:
     for transform in transformations:
         if isinstance(total_transform, sc.DataArray) and isinstance(
                 transform, sc.DataArray):
-            time = sc.concat([
-                total_transform.coords["time"].to(unit='ns', copy=False),
-                transform.coords["time"].to(unit='ns', copy=False)
-            ],
+            unit = _smaller_unit(transform.coords['time'],
+                                 total_transform.coords['time'])
+            total_transform.coords['time'] = total_transform.coords['time'].to(
+                unit=unit, copy=False)
+            transform.coords['time'] = transform.coords['time'].to(unit=unit,
+                                                                   copy=False)
+            time = sc.concat([total_transform.coords["time"], transform.coords["time"]],
                              dim="time")
-            time = sc.datetimes(values=np.unique(time.values), dims=["time"], unit='ns')
+            time = sc.datetimes(values=np.unique(time.values), dims=["time"], unit=unit)
             total_transform = _interpolate_transform(transform, time) \
                 * _interpolate_transform(total_transform, time)
         else:

tests/nxtransformations_test.py

@@ -52,6 +52,33 @@ def test_Transformation_with_single_value(nxroot):
     assert sc.identical(t[()], expected)
 
 
+def test_chain_with_single_values_and_different_unit(nxroot):
+    detector = create_detector(nxroot)
+    detector.create_field('depends_on', sc.scalar('/detector_0/transformations/t1'))
+    transformations = detector.create_class('transformations',
+                                            NX_class.NXtransformations)
+    value = sc.scalar(6.5, unit='mm')
+    offset = sc.spatial.translation(value=[1, 2, 3], unit='mm')
+    vector = sc.vector(value=[0, 0, 1])
+    t = value.to(unit='m') * vector
+    value1 = transformations.create_field('t1', value)
+    value1.attrs['depends_on'] = 't2'
+    value1.attrs['transformation_type'] = 'translation'
+    value1.attrs['offset'] = offset.values
+    value1.attrs['offset_units'] = str(offset.unit)
+    value1.attrs['vector'] = vector.value
+    value2 = transformations.create_field('t2', value.to(unit='cm'))
+    value2.attrs['depends_on'] = '.'
+    value2.attrs['transformation_type'] = 'translation'
+    value2.attrs['vector'] = vector.value
+
+    expected = sc.spatial.affine_transform(value=np.identity(4), unit=t.unit)
+    expected = expected * sc.spatial.translations(
+        dims=t.dims, values=2 * t.values, unit=t.unit)
+    expected = expected * sc.spatial.translation(value=[0.001, 0.002, 0.003], unit='m')
+    assert sc.identical(detector[...].coords['depends_on'], expected)
+
+
 def test_Transformation_with_multiple_values(nxroot):
     detector = create_detector(nxroot)
     detector.create_field('depends_on', sc.scalar('/detector_0/transformations/t1'))
@@ -81,3 +108,71 @@ def test_Transformation_with_multiple_values(nxroot):
     assert sc.identical(t.offset, offset)
     assert sc.identical(t.vector, vector)
     assert sc.identical(t[()], expected)
+
+
+def test_chain_with_multiple_values(nxroot):
+    detector = create_detector(nxroot)
+    detector.create_field('depends_on', sc.scalar('/detector_0/transformations/t1'))
+    transformations = detector.create_class('transformations',
+                                            NX_class.NXtransformations)
+    log = sc.DataArray(
+        sc.array(dims=['time'], values=[1.1, 2.2], unit='m'),
+        coords={'time': sc.array(dims=['time'], values=[11, 22], unit='s')})
+    log.coords['time'] = sc.epoch(unit='ns') + log.coords['time'].to(unit='ns')
+    offset = sc.spatial.translation(value=[1, 2, 3], unit='m')
+    vector = sc.vector(value=[0, 0, 1])
+    t = log * vector
+    t.data = sc.spatial.translations(dims=t.dims, values=t.values, unit=t.unit)
+    value1 = transformations.create_class('t1', NX_class.NXlog)
+    value1['time'] = log.coords['time'] - sc.epoch(unit='ns')
+    value1['value'] = log.data
+    value1.attrs['depends_on'] = 't2'
+    value1.attrs['transformation_type'] = 'translation'
+    value1.attrs['offset'] = offset.values
+    value1.attrs['offset_units'] = str(offset.unit)
+    value1.attrs['vector'] = vector.value
+    value2 = transformations.create_class('t2', NX_class.NXlog)
+    value2['time'] = log.coords['time'] - sc.epoch(unit='ns')
+    value2['value'] = log.data
+    value2.attrs['depends_on'] = '.'
+    value2.attrs['transformation_type'] = 'translation'
+    value2.attrs['vector'] = vector.value
+
+    expected = sc.spatial.affine_transform(value=np.identity(4), unit=t.unit)
+    expected = t * (t * (offset * expected))
+    assert sc.identical(detector[...].coords['depends_on'].value, expected)
+
+
+def test_chain_with_multiple_values_and_different_time_unit(nxroot):
+    detector = create_detector(nxroot)
+    detector.create_field('depends_on', sc.scalar('/detector_0/transformations/t1'))
+    transformations = detector.create_class('transformations',
+                                            NX_class.NXtransformations)
+    # Making sure to not use nanoseconds since that is used internally and may thus
+    # mask bugs.
+    log = sc.DataArray(
+        sc.array(dims=['time'], values=[1.1, 2.2], unit='m'),
+        coords={'time': sc.array(dims=['time'], values=[11, 22], unit='s')})
+    log.coords['time'] = sc.epoch(unit='us') + log.coords['time'].to(unit='us')
+    offset = sc.spatial.translation(value=[1, 2, 3], unit='m')
+    vector = sc.vector(value=[0, 0, 1])
+    t = log * vector
+    t.data = sc.spatial.translations(dims=t.dims, values=t.values, unit=t.unit)
+    value1 = transformations.create_class('t1', NX_class.NXlog)
+    value1['time'] = log.coords['time'] - sc.epoch(unit='us')
+    value1['value'] = log.data
+    value1.attrs['depends_on'] = 't2'
+    value1.attrs['transformation_type'] = 'translation'
+    value1.attrs['offset'] = offset.values
+    value1.attrs['offset_units'] = str(offset.unit)
+    value1.attrs['vector'] = vector.value
+    value2 = transformations.create_class('t2', NX_class.NXlog)
+    value2['time'] = log.coords['time'].to(unit='ms') - sc.epoch(unit='ms')
+    value2['value'] = log.data
+    value2.attrs['depends_on'] = '.'
+    value2.attrs['transformation_type'] = 'translation'
+    value2.attrs['vector'] = vector.value
+
+    expected = sc.spatial.affine_transform(value=np.identity(4), unit=t.unit)
+    expected = t * (t * (offset * expected))
+    assert sc.identical(detector[...].coords['depends_on'].value, expected)