Implement local fixed-frame Trajectory

src/scene_synthesis/__init__.py

@@ -6,3 +6,4 @@
 from .microphones import Microphone as Microphone
 from .scene import Scene as Scene
 from .sources import Source as Source
+from .trajectory import Trajectory as Trajectory

src/scene_synthesis/sources.py

@@ -1,10 +1,11 @@
 """Acoustic source definition and properties."""
 
 import numpy as np
-from acoular import SignalGenerator, Trajectory
+from acoular import SignalGenerator
 from traits.api import CArray, HasStrictTraits, Instance
 
 from scene_synthesis.directivities import Directivity
+from scene_synthesis.trajectory import Trajectory
 
 
 class Source(HasStrictTraits):
@@ -14,10 +15,10 @@ class Source(HasStrictTraits):
     --------
     Instantiate a simple source with a sine signal and a default trajectory:
 
-    >>> from acoular import SineGenerator, Trajectory
-    >>> from scene_synthesis.sources import Source
+    >>> from acoular import SineGenerator
+    >>> from scene_synthesis import Source, Trajectory
     >>> signal = SineGenerator(freq=1000, sample_freq=44100, num_samples=44100)
-    >>> trajectory = Trajectory()
+    >>> trajectory = Trajectory(points={0.0: (0.0, 0.0, 0.0), 1.0: (1.0, 0.0, 0.0)})
     >>> source = Source(signal=signal, trajectory=trajectory)
     """
 

src/scene_synthesis/trajectory.py

@@ -0,0 +1,125 @@
+"""Trajectory definition with a fixed frame of reference."""
+
+import numpy as np
+from scipy.interpolate import splev, splprep
+from traits.api import Dict, Float, HasStrictTraits, Property, Tuple, cached_property, property_depends_on
+
+
+class Trajectory(HasStrictTraits):
+    """Represent a source trajectory in a fixed frame of reference.
+
+    The trajectory is specified by a mapping from time instants to sampled
+    ``(x, y, z)`` positions. A spline is fit through those samples and can
+    then be evaluated at arbitrary times to obtain positions or time
+    derivatives such as velocity.
+
+    Notes
+    -----
+    - The spline order is chosen automatically based on the number of
+      available points, up to cubic interpolation.
+    - The frame of reference is fixed, i.e. the sampled positions are
+      interpreted directly as global coordinates.
+
+    Examples
+    --------
+    >>> import scene_synthesis as ss
+    >>> trajectory = ss.Trajectory(points={0.0: (0.0, 0.0, 0.0), 1.0: (1.0, 0.0, 0.0)})
+    >>> trajectory.location(0.5)
+    [array(0.5), array(0.), array(0.)]
+    """
+
+    #: Dictionary mapping time instants to sampled ``(x, y, z)`` positions.
+    points = Dict(
+        key_trait=Float,
+        value_trait=Tuple(Float, Float, Float),
+    )
+
+    #: Start and end time of the trajectory as ``(t_min, t_max)``.
+    interval = Property()
+
+    #: Internal spline representation returned by :func:`scipy.interpolate.splprep`.
+    tck = Property()
+
+    @property_depends_on(['points[]'])
+    def _get_interval(self):
+        return np.sort(list(self.points.keys()))[np.r_[0, -1]]
+
+    @cached_property
+    @property_depends_on(['points[]'])
+    def _get_tck(self):
+        t = np.sort(list(self.points.keys()))
+        xp = np.array([self.points[i] for i in t]).T
+        k = min(3, len(self.points) - 1)
+        tcku = splprep(xp, u=t, s=0, k=k)
+        return tcku[0]
+
+    def location(self, t, der=0):
+        """Evaluate the trajectory or one of its derivatives.
+
+        Parameters
+        ----------
+        t : float or array-like of float
+            Time instant or time instants at which to evaluate the trajectory.
+        der : int, optional
+            Derivative order. Use ``0`` for position, ``1`` for velocity,
+            ``2`` for acceleration, and so on. Defaults to ``0``.
+
+        Returns
+        -------
+        list[numpy.ndarray]
+            Three arrays representing the ``x``, ``y``, and ``z`` components
+            at the requested times.
+
+        Examples
+        --------
+        >>> import scene_synthesis as ss
+        >>> trajectory = ss.Trajectory(points={0.0: (0.0, 0.0, 0.0), 1.0: (1.0, 0.0, 0.0)})
+        >>> trajectory.location(0.5)
+        [array(0.5), array(0.), array(0.)]
+        """
+        return splev(t, self.tck, der)
+
+    def traj(self, t_start, t_end=None, delta_t=None, der=0):
+        """Iterate through trajectory samples over a time range.
+
+        Parameters
+        ----------
+        t_start : float
+            Start time of the iteration. If ``delta_t`` is omitted, this value
+            is interpreted as the step size and the full trajectory interval is
+            used.
+        t_end : float, optional
+            End time of the iteration. Defaults to the end of
+            :attr:`interval`.
+        delta_t : float, optional
+            Time step between yielded samples. If omitted, ``t_start`` is used
+            as the step size for traversing the full trajectory interval.
+        der : int, optional
+            Derivative order to evaluate. Defaults to ``0``.
+
+        Yields
+        ------
+        tuple[numpy.float64, numpy.float64, numpy.float64]
+            The interpolated ``(x, y, z)`` values at each sampled time.
+
+        Examples
+        --------
+        >>> import scene_synthesis as ss
+        >>> trajectory = ss.Trajectory(points={0.0: (0.0, 0.0, 0.0), 1.0: (1.0, 0.0, 0.0)})
+        >>> samples = list(trajectory.traj(0.5))
+        >>> samples[0]
+        (np.float64(0.0), np.float64(0.0), np.float64(0.0))
+        >>> samples[1]
+        (np.float64(0.5), np.float64(0.0), np.float64(0.0))
+        >>> samples = list(trajectory.traj(0.0, 1.0, 0.5))
+        >>> samples[0]
+        (np.float64(0.0), np.float64(0.0), np.float64(0.0))
+        >>> samples[1]
+        (np.float64(0.5), np.float64(0.0), np.float64(0.0))
+        """
+        if not delta_t:
+            delta_t = t_start
+            t_start, t_end = self.interval
+        if not t_end:
+            t_end = self.interval[1]
+        yield from zip(*self.location(np.arange(t_start, t_end, delta_t), der), strict=True)

tests/cases_trajectory.py

@@ -1,11 +1,11 @@
 """Test cases for source trajectories."""
 
-import acoular as ac
 import numpy as np
+import scene_synthesis as ss
 
 
 class Trajectories:
-    """Test cases for all :class:`acoular.trajectories.Trajectory`-derived classes.
+    """Test cases for all :class:`scene_synthesis.trajectory.Trajectory` objects.
 
     New trajectories should be added here.
     """
@@ -17,25 +17,25 @@ def case_none(self):
     def case_static(self):
         """Static trajectory test case."""
         points = {0.0: (1.0, 0.0, 0.0), 1.0: (1.0, 0.0, 0.0)}
-        return [ac.Trajectory(points=points)]
+        return [ss.Trajectory(points=points)]
 
     def case_linear_pass(self):
         """Linear pass trajectory (fly by) test case."""
         points = {0.0: (-1.0, -1.0, 1.0), 1.0: (1.0, 1.0, 1.0)}
-        return [ac.Trajectory(points=points)]
+        return [ss.Trajectory(points=points)]
 
     def case_linear_approach(self):
         """Linear approach trajectory (fly at) test case."""
         points = {0.0: (5.0, 0.0, 0.0), 1.0: (0.5, 0.0, 0.0)}
-        return [ac.Trajectory(points=points)]
+        return [ss.Trajectory(points=points)]
 
     def case_circular(self):
         """Circular trajectory (fly around) test case."""
         n = 3600
         points = {i / n: (1.0 * np.cos(2 * np.pi * i / n), 1.0 * np.sin(2 * np.pi * i / n), 0.0) for i in range(n + 1)}
-        return [ac.Trajectory(points=points)]
+        return [ss.Trajectory(points=points)]
 
     def case_static_array(self):
         """Static array trajectory test case."""
         points_array = [{0.0: (x, y, 1.0), 1.0: (x, y, 1.0)} for x, y in [(1.0, 1.0), (1.0, 0.0), (0.0, 0.0)]]
-        return [ac.Trajectory(points=points) for points in points_array]
+        return [ss.Trajectory(points=points) for points in points_array]

tests/test_trajectory.py

@@ -0,0 +1,31 @@
+"""Unit tests for the scene-synthesis trajectory class."""
+
+import numpy as np
+import scene_synthesis as ss
+
+
+def test_trajectory_interval_uses_point_bounds():
+    """``Trajectory.interval`` should expose the first and last point times."""
+    trajectory = ss.Trajectory(points={0.5: (0.0, 0.0, 0.0), 2.0: (1.0, 0.0, 0.0), 1.0: (0.5, 0.0, 0.0)})
+
+    np.testing.assert_allclose(trajectory.interval, np.array([0.5, 2.0]))
+
+
+def test_trajectory_location_matches_sampled_points():
+    """``Trajectory.location`` should pass through the sampled points."""
+    trajectory = ss.Trajectory(points={0.0: (0.0, 0.0, 0.0), 1.0: (1.0, 2.0, 0.0), 2.0: (2.0, 4.0, 0.0)})
+
+    location = np.array(trajectory.location(1.0))
+
+    np.testing.assert_allclose(location, np.array([1.0, 2.0, 0.0]))
+
+
+def test_trajectory_traj_iterates_over_requested_range():
+    """``Trajectory.traj`` should iterate over positions with the requested step size."""
+    trajectory = ss.Trajectory(points={0.0: (0.0, 0.0, 0.0), 1.0: (1.0, 0.0, 0.0)})
+
+    samples = list(trajectory.traj(0.0, 1.0, 0.25))
+
+    assert len(samples) == 4
+    np.testing.assert_allclose(samples[0], (0.0, 0.0, 0.0))
+    np.testing.assert_allclose(samples[-1], (0.75, 0.0, 0.0))