diff --git a/sasdata/abscissa.py b/sasdata/abscissa.py new file mode 100644 index 00000000..a0c16a1e --- /dev/null +++ b/sasdata/abscissa.py @@ -0,0 +1,123 @@ +from abc import ABC, abstractmethod + +import numpy as np +from exceptions import InterpretationError +from numpy._typing import ArrayLike +from quantities.quantity import Quantity +from util import is_increasing + + +class Abscissa(ABC): + + def __init__(self, axes: list[Quantity]): + self._axes = axes + self._dimensionality = len(axes) + + @property + def dimensionality(self) -> int: + """ Dimensionality of this data """ + return self._dimensionality + + @property + @abstractmethod + def is_grid(self) -> bool: + """ Are these coordinates using a grid representation + (as opposed to a general list representation) + + is_grid = True: implies that the corresponding ordinate is n-dimensional tensor + is_grid = False: implies that the corresponding ordinate is a 1D list + + If the data is one dimensional, is_grid=True + """ + + @property + def axes(self) -> list[Quantity]: + """ Axes of the data: + + If it's an (n1-by-n2-by-n3...) grid (is_grid=True): give the values for each axis, returning a list like + [Quantity(length n1), Quantity(length n2), Quantity(length n3) ... ] + + If it is not grid data (is_grid=False), but n points on a general mesh, give one array for each dimension + [Quantity(length n), Quantity(length n), Quantity(length n) ... ] + """ + + return self._axes + + @staticmethod + def _determine_error_message(axis_arrays: list[np.ndarray], ordinate_shape: tuple): + """ Error message for the `.determine` function""" + + shape_string = ", ".join([str(axis.shape) for axis in axis_arrays]) + + return f"Cannot interpret array shapes axis: [{shape_string}], ordinate: {ordinate_shape}" + + @staticmethod + def determine(axis_data: list[Quantity[ArrayLike]], ordinate_data: Quantity[ArrayLike]) -> "Abscissa": + """ Get an Abscissa object that fits the combination of axes and data""" + + # Different posibilites: + # 1: axes_data[i].shape == axes_data[j].shape == ordinate_data.shape + # 1a: axis_data[i] is 1D => + # 1a-i: len(axes_data) == 1 => Grid type or Scatter type depending on sortedness + # 1a-ii: len(axes_data) != 1 => Scatter type + # 1b: axis_data[i] dimensionality matches len(axis_data) => Meshgrid type + # 1c: other => Error + # 2: (len(axes_data[0]), len(axes_data[1])... ) == ordinate_data.shape => Grid type + # 3: None of the above => Error + + ordinate_shape = np.array(ordinate_data.value).shape + axis_arrays = [np.array(axis.value) for axis in axis_data] + + # 1: + if all([axis.shape == ordinate_shape for axis in axis_arrays]): + # 1a: + if all([len(axis.shape)== 1 for axis in axis_arrays]): + # 1a-i: + if len(axis_arrays) == 1: + # Is it sorted + if is_increasing(axis_arrays[0]): + return GridAbscissa(axis_data) + else: + return ScatterAbscissa(axis_data) + # 1a-ii + else: + return ScatterAbscissa(axis_data) + # 1b + elif all([len(axis.shape) == len(axis_arrays) for axis in axis_arrays]): + + return MeshgridAbscissa(axis_data) + + else: + raise InterpretationError(Abscissa._determine_error_message(axis_arrays, ordinate_shape)) + + elif all([len(axis.shape)== 1 for axis in axis_arrays]) and \ + tuple([axis.shape[0] for axis in axis_arrays]) == ordinate_shape: + + # Require that they are sorted + if all([is_increasing(axis) for axis in axis_arrays]): + + return GridAbscissa(axis_data) + + else: + raise InterpretationError("Grid axes are not sorted") + + else: + raise InterpretationError(Abscissa._determine_error_message(axis_arrays, ordinate_shape)) + +class GridAbscissa(Abscissa): + + @property + def is_grid(self): + return True + +class MeshgridAbscissa(Abscissa): + + @property + def is_grid(self): + return True + +class ScatterAbscissa(Abscissa): + + @property + def is_grid(self): + return False diff --git a/sasdata/data.py b/sasdata/data.py index 73e64e61..afb280d8 100644 --- a/sasdata/data.py +++ b/sasdata/data.py @@ -11,6 +11,58 @@ class SasData: + """ General object containing data in the SasView ecosystem""" + + def __init__(self, + name: str, + ordinate: Quantity, + mask: Quantity, + abscissae: list[Quantity], + dependents: list["SasData"]): + + self._ordinate = ordinate + self._abscissae = abscissae + self._mask = mask + + @property + def ordinate(self) -> Quantity: + return self._ordinate + + @property + def abscissae(self) -> list[Quantity]: + return self._abscissae + + @property + def mask(self) -> Quantity: + return self._mask + + def scatter_data(self): + """ Return data in the coordinate/value form [(x1, x2, x3, y)...]""" + + +class SasDerivedMeasurement(SasData): + """ General object sas measurement that has not come directly from a file, + for example, the difference between two datasets""" + + + def __init__(self, + name: str, + ordinate: Quantity, + abscissae: list[Quantity], + dependents: list["SasData"], + metadata: DerivedMetadata): + + super().__init__( + name=name, + ordinate=ordinate, + abscissae=abscissae, + dependents=dependents) + + self.metadata = metadata + + +class SasMeasurement(SasData): + def __init__( self, name: str, @@ -119,7 +171,6 @@ def _save_h5(self, sasentry: HDF5Group): for idx, (key, sasdata) in enumerate(self._data_contents.items()): sasdata.as_h5(group, key) - @staticmethod def save_h5(data: dict[str, typing.Self], path: str | typing.BinaryIO): with h5py.File(path, "w") as f: @@ -130,7 +181,6 @@ def save_h5(data: dict[str, typing.Self], path: str | typing.BinaryIO): data._save_h5(sasentry) - class SasDataEncoder(MetadataEncoder): def default(self, obj): match obj: diff --git a/sasdata/exceptions.py b/sasdata/exceptions.py new file mode 100644 index 00000000..6fc3e545 --- /dev/null +++ b/sasdata/exceptions.py @@ -0,0 +1,2 @@ +class InterpretationError(Exception): + """ Error interpreting data """ diff --git a/sasdata/util.py b/sasdata/util.py index 8decc68b..5dd0d404 100644 --- a/sasdata/util.py +++ b/sasdata/util.py @@ -1,6 +1,8 @@ from collections.abc import Callable from typing import TypeVar +import numpy as np + T = TypeVar("T") def cache[T](fun: Callable[[], T]): @@ -16,3 +18,15 @@ def wrapper() -> T: return cache_state[1] return wrapper + +def is_increasing(data: np.ndarray): + """ Check if a 1D array is sorted in strictly increasing order""" + return np.all(data[1:] > data[:-1]) + +def is_decreasing(data: np.ndarray): + """ Check if a 1D array is sorted in strictly decreasing order""" + return np.all(data[1:] < data[:-1]) + +def is_sorted(data: np.ndarray): + """ Check if a 1D array is strictly sorted """ + return is_increasing(data) or is_decreasing(data) diff --git a/test/utest_asbscissa.py b/test/utest_asbscissa.py new file mode 100644 index 00000000..49debb3d --- /dev/null +++ b/test/utest_asbscissa.py @@ -0,0 +1,71 @@ +import numpy as np +import pytest +from abscissa import Abscissa, GridAbscissa, MeshgridAbscissa, ScatterAbscissa +from exceptions import InterpretationError +from quantities.quantity import Quantity +from quantities.units import none + + +def test_deterimine_1d_grid(): + """ Test that 1D ordered data is grid type""" + q = Quantity(np.arange(10), units=none) + + determined = Abscissa.determine([q], q) + + assert isinstance(determined, GridAbscissa) + +def test_deterimine_1d_scatter(): + """ Test that 1D unordered data is scatter type """ + a = Quantity(np.array([1, 2, 3, 4, 5, 0, 9, 8, 7, 6]), units=none) + d = Quantity(np.arange(10), units=none) + + determined = Abscissa.determine([a], d) + + assert isinstance(determined, ScatterAbscissa) + +def test_2D_scatter(): + """ Test the nD scatter case with 2D data """ + q = Quantity(np.arange(10), units=none) + + determined = Abscissa.determine([q, q], q) + + assert isinstance(determined, ScatterAbscissa) + +def test_2D_meshgrid(): + """ Test the meshgrid case with 2x5""" + q = Quantity(np.arange(10).reshape(2, 5), units=none) + + determined = Abscissa.determine([q, q], q) + + assert isinstance(determined, MeshgridAbscissa) + + +def test_2D_grid(): + """ Test the nD grid case with 2x5 """ + a1 = Quantity(np.arange(2), units=none) + a2 = Quantity(np.arange(5), units=none) + + d = Quantity(np.arange(10).reshape(2, 5), units=none) + + determined = Abscissa.determine([a1, a2], d) + + assert isinstance(determined, GridAbscissa) + +def test_2D_grid_axis_error(): + """ Test the nD grid case with bad axes """ + + a1 = Quantity(np.arange(2), units=none) + a2 = Quantity(np.arange(5), units=none) + + d = Quantity(np.arange(10).reshape(5, 2, 1), units=none) + + with pytest.raises(InterpretationError): + Abscissa.determine([a1, a2], d) + +def test_2D_meshgrid_error_mismatched_dimensionality(): + """ Test the nD meshgrid case with bad axes """ + + q = Quantity(np.arange(10).reshape(5, 2), units=none) + + with pytest.raises(InterpretationError): + Abscissa.determine([q, q, q], q) # three axes, each 2D