Orientation parameter for RealSphericalHarmonicsDirectivity module.

CHANGELOG.rst

@@ -11,6 +11,11 @@ adheres to `Semantic Versioning <http://semver.org/spec/v2.0.0.html>`_.
 `Unreleased`_
 -------------
 
+Added
+~~~~~~~
+
+- Added ``orientation`` parameter to ``RealSphericalHarmonicsDirectivity`` object.
+
 Fixed
 ~~~~~
 

pyroomacoustics/directivities/harmonics.py

@@ -52,8 +52,9 @@
 import numpy as np
 import scipy.special
 
-from ..doa import cart2spher
+from ..doa import cart2spher, spher2cart
 from .base import Directivity
+from .direction import Rotation3D
 
 
 def get_mn_in_acn_order(degree):
@@ -140,16 +141,44 @@ class RealSphericalHarmonicsDirectivity(Directivity):
         Order of the spherical harmonic.
     n: int
         Degree of the spherical harmonic.
+    orientation: Rotation3D, optional
+        A rotation to apply to the pattern. If not provided, the
+        identity rotation is used, i.e., the harmonic is evaluated in
+        its canonical frame with no reorientation.
     condon_shortley_phase: bool, optional
         If True, includes the Condon-Shortley phase factor (-1)^m. Default is False.
     """
 
-    def __init__(self, m, n, condon_shortley_phase: bool = False):
+    def __init__(self, m, n, orientation=None, condon_shortley_phase: bool = False):
         self.m = m
         self.n = n
-
         self.condon_shortley_phase = condon_shortley_phase
 
+        if orientation is None:
+            orientation = Rotation3D(
+                [0.0, 0.0, 0.0]
+            )  # identity rotation, no-op default
+        elif not isinstance(orientation, Rotation3D):
+            raise TypeError(
+                f"orientation must be a Rotation3D object, got {type(orientation).__name__}"
+            )
+        self._orientation = orientation
+
+    def set_orientation(self, orientation):
+        """
+        Set orientation of directivity pattern.
+
+        Parameters
+        ----------
+        orientation: Rotation3D
+            New direction for the directivity pattern.
+        """
+        if not isinstance(orientation, Rotation3D):
+            raise TypeError(
+                f"orientation must be a Rotation3D object, got {type(orientation).__name__}"
+            )
+        self._orientation = orientation
+
     @property
     def is_impulse_response(self):
         return True
@@ -159,18 +188,10 @@ def filter_len_ir(self):
         return 1
 
     def get_response_cartesian(self, directions, magnitude=False, frequency=None):
-        az, co, _ = cart2spher(directions.T)
-        return self.get_response(
-            az, co, magnitude=magnitude, frequency=frequency, degrees=False
-        )
+        local_dirs = self._orientation.rotate_transpose(directions.T)
 
-    def get_response(
-        self, azimuth, colatitude=None, magnitude=False, frequency=None, degrees=True
-    ):
+        azimuth, colatitude, _ = cart2spher(local_dirs)
 
-        if degrees:
-            azimuth = np.radians(azimuth)
-            colatitude = np.radians(colatitude)
         return real_sph_harm(
             self.n,
             self.m,
@@ -179,6 +200,15 @@ def get_response(
             condon_shortley_phase=self.condon_shortley_phase,
         )[:, np.newaxis]
 
+    def get_response(
+        self, azimuth, colatitude=None, magnitude=False, frequency=None, degrees=True
+    ):
+
+        # world-frame cartesian unit vectors
+        directions = spher2cart(azimuth, colatitude, degrees=degrees).T
+
+        return self.get_response_cartesian(directions, magnitude=magnitude)
+
     def sample_rays(self, n_rays, rng=None):
         """Not yet implemented."""
         raise NotImplementedError

tests/directivities/test_harmonics_directivities.py

@@ -6,6 +6,25 @@
     real_sph_harm,
 )
 
+FS = 16000
+MIC_POS = [1.0, 1.2, 1.5]
+SRC_POS = [3.5, 2.8, 1.7]
+ROOM_DIM = [6.0, 5.0, 3.0]
+
+
+def _simulate_rir(directivity, max_order=0):
+    """Run a room simulation with a single directional microphone."""
+    if max_order == 0:
+        room = pra.AnechoicRoom(fs=FS, dim=3)
+    else:
+        mat = pra.Material(energy_absorption=0.3)
+        room = pra.ShoeBox(ROOM_DIM, fs=FS, materials=mat, max_order=max_order)
+
+    room.add_source(SRC_POS)
+    room.add_microphone(MIC_POS, directivity=directivity)
+    room.compute_rir()
+    return room.rir[0][0]
+
 
 def test_harmonics_directivity():
     """
@@ -38,3 +57,69 @@ def test_harmonics_directivity():
     )
 
     np.allclose(rir_sum, sh, atol=1e-2)
+
+
+def test_harmonics_rotated_z():
+    """
+    Check that the `orientation` attribute of `RealSphericalHarmonicsDirectivity`
+    correctly rotates the pattern about the z-axis: an (m=1, n=1) harmonic
+    rotated by 90 degrees about z should produce the same RIR as an
+    unrotated (m=-1, n=1) harmonic at the same position.
+    """
+    max_order = 2
+
+    rot_90_z = pra.directivities.Rotation3D([90.0], "z", degrees=True)
+
+    dir_rotated = pra.directivities.RealSphericalHarmonicsDirectivity(
+        m=1, n=1, orientation=rot_90_z
+    )
+    dir_unrotated = pra.directivities.RealSphericalHarmonicsDirectivity(m=-1, n=1)
+
+    rir_rotated = _simulate_rir(dir_rotated, max_order=max_order)
+    rir_reference = _simulate_rir(dir_unrotated, max_order=max_order)
+
+    np.testing.assert_allclose(rir_rotated, rir_reference, atol=1e-9, rtol=0)
+
+
+def test_harmonics_rotated_y():
+    """
+    Check that the `orientation` attribute of `RealSphericalHarmonicsDirectivity`
+    correctly rotates the pattern about the y-axis: an (m=0, n=1) harmonic
+    rotated by 90 degrees about y should produce the same RIR as an
+    unrotated (m=1, n=1) harmonic at the same position.
+    """
+    max_order = 2
+
+    rot_90_y = pra.directivities.Rotation3D([90.0], "y", degrees=True)
+
+    dir_rotated = pra.directivities.RealSphericalHarmonicsDirectivity(
+        m=0, n=1, orientation=rot_90_y
+    )
+    dir_unrotated = pra.directivities.RealSphericalHarmonicsDirectivity(m=1, n=1)
+
+    rir_rotated = _simulate_rir(dir_rotated, max_order=max_order)
+    rir_reference = _simulate_rir(dir_unrotated, max_order=max_order)
+
+    np.testing.assert_allclose(rir_rotated, rir_reference, atol=1e-9, rtol=0)
+
+
+def test_harmonics_rotated_x():
+    """
+    Check that the `orientation` attribute of `RealSphericalHarmonicsDirectivity`
+    correctly rotates the pattern about the x-axis: an (m=0, n=1) harmonic
+    rotated by -90 degrees about x should produce the same RIR as an
+    unrotated (m=-1, n=1) harmonic at the same position.
+    """
+    max_order = 2
+
+    rot_m90_x = pra.directivities.Rotation3D([-90.0], "x", degrees=True)
+
+    dir_rotated = pra.directivities.RealSphericalHarmonicsDirectivity(
+        m=0, n=1, orientation=rot_m90_x
+    )
+    dir_unrotated = pra.directivities.RealSphericalHarmonicsDirectivity(m=-1, n=1)
+
+    rir_rotated = _simulate_rir(dir_rotated, max_order=max_order)
+    rir_reference = _simulate_rir(dir_unrotated, max_order=max_order)
+
+    np.testing.assert_allclose(rir_rotated, rir_reference, atol=1e-9, rtol=0)