Added ctis.inverters.MartInverter, an implementation of the Multiplicative Algebraic Reconstruction Technique (MART).

.github/workflows/tests.yml

@@ -32,6 +32,8 @@ jobs:
           pip install setuptools wheel
           pip install -e .[test]
       - name: Test with pytest
+        env:
+          MPLBACKEND: "agg"
         run: |
           pip install pytest pytest-cov
           pytest --cov=. --cov-report=xml

ctis/__init__.py

@@ -5,8 +5,10 @@
 
 from . import scenes
 from . import instruments
+from . import inverters
 
 __all__ = [
     "scenes",
     "instruments",
+    "inverters",
 ]

ctis/instruments/_instruments.py

@@ -114,12 +114,19 @@ def uncertainty(self) -> Callable[[na.ScalarArray], na.ScalarArray]:
         for a given number of photons.
         """
 
+    @property
+    @abc.abstractmethod
+    def channel(self):
+        """
+        Human-readable name of each independent CTIS channel.
+        """
+
     @property
     @abc.abstractmethod
     def axis_channel(self) -> str | tuple[str, ...]:
         """
         The logical axis or axes of this instrument corresponding to
-        the different dispersion magnitudes and angles.
+        the different CTIS channels.
         """
 
     @property
@@ -391,6 +398,11 @@ class IdealInstrument(
     A grid of wavelength and position coordinates on the sensor plane.
     """
 
+    channel: str | na.AbstractScalar = dataclasses.MISSING
+    """
+    Human-readable name of each independent CTIS channel.
+    """
+
     axis_channel: str | tuple[str, ...] = dataclasses.MISSING
     """
     The logical axis or axes of this instrument corresponding to
@@ -492,10 +504,6 @@ def _coordinates_output(self) -> na.AbstractSpectralPositionalVectorArray:
 
         coordinates_output = coordinates_output.cell_centers(self.axis_wavelength)
 
-        p = coordinates_output.position
-        coordinates_output.position.x = na.random.normal(p.x, 1e-3 * u.pix)
-        coordinates_output.position.y = na.random.normal(p.y, 1e-3 * u.pix)
-
         return coordinates_output
 
     @functools.cached_property
@@ -518,12 +526,12 @@ def weights_transpose(self):
         coordinates_input = self._coordinates_input
         coordinates_output = self._coordinates_output
 
-        return na.regridding.weights(
-            coordinates_input=coordinates_output.position,
-            coordinates_output=coordinates_input.position,
-            axis_input=self.axis_sensor_xy,
-            axis_output=self.axis_scene_xy,
-            method="conservative",
+        return na.regridding.transpose_weights_conservative(
+            weights=self.weights,
+            coordinates_input=coordinates_input.position,
+            coordinates_output=coordinates_output.position,
+            axis_input=self.axis_scene_xy,
+            axis_output=self.axis_sensor_xy,
         )
 
     def image(

ctis/instruments/_instruments_test.py

@@ -7,6 +7,8 @@
 
 velocity = na.linspace(-500, 500, axis="wavelength", num=21) * u.km / u.s
 
+wavelength_rest = 171 * u.AA
+
 position_scene = na.Cartesian2dVectorLinearSpace(
     start=-20 * u.arcsec,
     stop=20 * u.arcsec,
@@ -19,15 +21,18 @@
     y=na.arange(0, 64, axis="sensor_y") * u.pix,
 )
 
-coordinates_scene = na.SpectralPositionalVectorArray(velocity, position_scene)
-coordinates_sensor = na.SpectralPositionalVectorArray(velocity, position_sensor)
-
-gaussians = ctis.scenes.gaussians(
-    inputs=coordinates_scene,
-    width=na.SpectralPositionalVectorArray(30 * u.km / u.s, 1 * u.arcsec),
+coordinates_scene = na.DopplerPositionalVectorArray.from_velocity(
+    velocity=velocity,
+    wavelength_rest=wavelength_rest,
+    position=position_scene,
+)
+coordinates_sensor = na.DopplerPositionalVectorArray.from_velocity(
+    velocity=velocity,
+    wavelength_rest=wavelength_rest,
+    position=position_sensor,
 )
 
-wavelength_rest = 171 * u.AA
+gaussians = ctis.scenes.gaussians(coordinates_scene)
 
 AA = dict(
     unit=u.AA,
@@ -40,16 +45,19 @@
 dispersion = 200 * u.km / u.s
 dispersion = (dispersion.to(**AA) - wavelength_rest) / u.pix
 
+angle = na.linspace(0, 360, axis="channel", num=3, endpoint=False)
+
 instrument_ideal = ctis.instruments.IdealInstrument(
     area_effective=1 * u.cm**2,
     timedelta_exposure=10 * u.s,
     plate_scale=2 * u.arcsec / u.pix,
     dispersion=dispersion,
-    angle=na.linspace(0, 360, axis="channel", num=3, endpoint=False),
+    angle=angle,
     wavelength_ref=wavelength_rest,
     position_ref=32 * u.pix,
     coordinates_scene=coordinates_scene,
     coordinates_sensor=coordinates_sensor,
+    channel=angle,
     axis_channel="channel",
     axis_wavelength="wavelength",
     axis_scene_xy=("scene_x", "scene_y"),

ctis/inverters/__init__.py

@@ -0,0 +1,20 @@
+"""Inversion algorithms which can reconstruct scenes from observed images."""
+
+from . import merit
+from ._results import AbstractInversionResult, InversionResult
+from ._inverters import AbstractInverter
+from ._iterative import (
+    AbstractIterativeInverter,
+    MartInverter,
+    IterativeInversionResult,
+)
+
+__all__ = [
+    "merit",
+    "AbstractInverter",
+    "AbstractIterativeInverter",
+    "MartInverter",
+    "AbstractInversionResult",
+    "InversionResult",
+    "IterativeInversionResult",
+]

ctis/inverters/_inverters.py

@@ -0,0 +1,48 @@
+import abc
+import dataclasses
+import named_arrays as na
+import ctis
+from ._results import InversionResult
+
+__all__ = [
+    "AbstractInverter",
+]
+
+
+@dataclasses.dataclass
+class AbstractInverter(
+    abc.ABC,
+):
+    """
+    An interface describing an algorithm which can invert CTIS observations
+    to yield a reconstruction of the observed scene.
+    """
+
+    @property
+    @abc.abstractmethod
+    def instrument(self) -> ctis.instruments.AbstractInstrument:
+        """
+        A model of a CTIS instrument which transforms the radiance of an observed
+        scene to photons measured by the sensors.
+        """
+
+    @abc.abstractmethod
+    def __call__(
+        self,
+        images: na.FunctionArray[na.SpectralPositionalVectorArray, na.ScalarArray],
+        **kwargs,
+    ) -> InversionResult:
+        """
+        Reconstruct a scene using the observed images.
+
+        Parameters
+        ----------
+        images
+            The observed images used to calculate the reconstruction.
+            Must be evaluated on the same coordinates as
+            :attr:`~ctis.instruments.AbstractInstrument.coordinates_sensor`
+            attribute of :attr:`instrument`.
+        kwargs
+            Additional keyword arguments which can be used by subclass
+            implementations.
+        """

ctis/inverters/_inverters_test.py

@@ -0,0 +1,31 @@
+import abc
+import numpy as np
+import named_arrays as na
+import ctis
+
+
+class AbstractTestAbstractInverter(
+    abc.ABC,
+):
+
+    def test_instrument(self, a: ctis.inverters.AbstractInverter):
+        result = a.instrument
+        assert isinstance(result, ctis.instruments.AbstractInstrument)
+
+    def test__call__(
+        self,
+        a: ctis.inverters.AbstractInverter,
+        images: na.FunctionArray[na.SpectralPositionalVectorArray, na.ScalarArray],
+        **kwargs,
+    ) -> ctis.inverters.AbstractInversionResult:
+        result = a(images, **kwargs)
+
+        assert isinstance(result, ctis.inverters.AbstractInversionResult)
+
+        assert result.solution.sum() > 0
+        assert isinstance(result.success, bool)
+        assert isinstance(result.message, str)
+        assert np.all(result.images == images)
+        assert result.inverter == a
+
+        return result

ctis/inverters/_iterative/__init__.py

@@ -0,0 +1,8 @@
+from ._iterative import AbstractIterativeInverter, IterativeInversionResult
+from ._mart import MartInverter
+
+__all__ = [
+    "AbstractIterativeInverter",
+    "IterativeInversionResult",
+    "MartInverter",
+]