From f48423118a667c08b232b13d51d3ffd6072c2556 Mon Sep 17 00:00:00 2001
From: Anne Heimes <ahe@akustik.rwth-aachen.de>
Date: Fri, 25 Apr 2025 12:09:18 +0200
Subject: [PATCH 01/12] draft correlation method after mommertz

---
 docs/api_reference.rst                      |   2 +-
 docs/modules/imkar.scattering.freefield.rst |   7 +
 imkar/__init__.py                           |   7 +
 imkar/scattering/__init__.py                |   7 +
 imkar/scattering/freefield.py               | 108 +++++++++++++++
 tests/test_scattering_freefield.py          | 138 ++++++++++++++++++++
 6 files changed, 268 insertions(+), 1 deletion(-)
 create mode 100644 docs/modules/imkar.scattering.freefield.rst
 create mode 100644 imkar/scattering/__init__.py
 create mode 100644 imkar/scattering/freefield.py
 create mode 100644 tests/test_scattering_freefield.py

diff --git a/docs/api_reference.rst b/docs/api_reference.rst
index 8f0383c..49d60df 100644
--- a/docs/api_reference.rst
+++ b/docs/api_reference.rst
@@ -16,7 +16,7 @@ Modules
 .. toctree::
    :maxdepth: 1
 
-   modules/imkar
+   modules/imkar.scattering.freefield
 
 
 .. _examples gallery: https://pyfar-gallery.readthedocs.io/en/latest/examples_gallery.html
diff --git a/docs/modules/imkar.scattering.freefield.rst b/docs/modules/imkar.scattering.freefield.rst
new file mode 100644
index 0000000..e5a349e
--- /dev/null
+++ b/docs/modules/imkar.scattering.freefield.rst
@@ -0,0 +1,7 @@
+imkar.scattering.freefield
+==========================
+
+.. automodule:: imkar.scattering.freefield
+   :members:
+   :undoc-members:
+   :show-inheritance:
diff --git a/imkar/__init__.py b/imkar/__init__.py
index 0dcfbd5..ac544ff 100644
--- a/imkar/__init__.py
+++ b/imkar/__init__.py
@@ -5,3 +5,10 @@
 __author__ = """The pyfar developers"""
 __email__ = ''
 __version__ = '0.1.0'
+
+
+from . import scattering
+
+__all__ = [
+    "scattering",
+]
diff --git a/imkar/scattering/__init__.py b/imkar/scattering/__init__.py
new file mode 100644
index 0000000..ce0af71
--- /dev/null
+++ b/imkar/scattering/__init__.py
@@ -0,0 +1,7 @@
+"""Imkar scattering module."""
+
+from . import freefield
+
+__all__ = [
+    "freefield",
+]
diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
new file mode 100644
index 0000000..e88224c
--- /dev/null
+++ b/imkar/scattering/freefield.py
@@ -0,0 +1,108 @@
+"""Scattering calculation functions based on free-field data."""
+import numpy as np
+import pyfar as pf
+
+
+def correlation_method(
+        sample_pressure, reference_pressure, microphone_weights):
+    r"""
+    Calculate the incident-dependent free-field scattering coefficient.
+
+    This function uses the Mommertz correlation method [#]_ to compute the
+    scattering coefficient of the input data:
+
+    .. math::
+        s = 1 -
+            \frac{|\sum_w \underline{p}_{\text{sample}}(\vartheta,\varphi)
+            \cdot \underline{p}_{\text{reference}}^*(\vartheta,\varphi)
+            \cdot w(\vartheta,\varphi)|^2}
+            {\sum_w |\underline{p}_{\text{sample}}(\vartheta,\varphi)|^2
+            \cdot w(\vartheta,\varphi) \cdot \sum_w
+            |\underline{p}_{\text{reference}}(\vartheta,\varphi)|^2
+            \cdot w(\vartheta,\varphi) }
+
+    where:
+        - :math:`\underline{p}_{\text{sample}}` is the reflected sound
+          pressure of the sample under investigation.
+        - :math:`\underline{p}_{\text{reference}}` is the reflected sound
+          pressure from the reference sample.
+        - :math:`w` represents the area weights of the sampling, and
+          :math:`\vartheta` and :math:`\varphi` are the ``colatitude``
+          and ``azimuth`` angles from the
+          :py:class:`~pyfar.classes.coordinates.Coordinates` object.
+
+    The test sample is assumed to lie in the x-y-plane.
+
+    Parameters
+    ----------
+    sample_pressure : :py:class:`~pyfar.classes.audio.FrequencyData`
+        Reflected sound pressure or directivity of the test sample. Its cshape
+        must be (..., microphone_weights.size) and broadcastable to the
+        cshape of ``reference_pressure``.
+    reference_pressure : :py:class:`~pyfar.classes.audio.FrequencyData`
+        Reflected sound pressure or directivity of the reference sample. Its
+        cshape must be (..., microphone_weights.size) and broadcastable to the
+        cshape of ``sample_pressure``.
+    microphone_weights : array_like
+        1D array containing the area weights for the microphone positions.
+        No normalization is required. Its shape must match the last dimension
+        in the cshape of ``sample_pressure`` and ``reference_pressure``.
+
+    Returns
+    -------
+    scattering_coefficients : :py:class:`~pyfar.classes.audio.FrequencyData`
+        The scattering coefficient for each incident direction as a function
+        of frequency.
+
+    References
+    ----------
+    .. [#]  E. Mommertz, "Determination of scattering coefficients from the
+            reflection directivity of architectural surfaces," Applied
+            Acoustics, vol. 60, no. 2, pp. 201-203, June 2000,
+            doi: 10.1016/S0003-682X(99)00057-2.
+
+    """
+    # check input types
+    if not isinstance(sample_pressure, pf.FrequencyData):
+        raise TypeError("sample_pressure must be of type pyfar.FrequencyData")
+    if not isinstance(reference_pressure, pf.FrequencyData):
+        raise TypeError(
+            "reference_pressure must be of type pyfar.FrequencyData")
+    microphone_weights = np.atleast_1d(
+        np.asarray(microphone_weights, dtype=float))
+
+    # check input dimensions
+    if sample_pressure.cshape[-1] != microphone_weights.size:
+        raise ValueError(
+            "The last dimension of sample_pressure must match the size of "
+            "microphone_weights")
+    if reference_pressure.cshape[-1] != microphone_weights.size:
+        raise ValueError(
+            "The last dimension of reference_pressure must match the size of "
+            "microphone_weights")
+
+    if sample_pressure.cshape[:-1] != reference_pressure.cshape[:-1]:
+        raise ValueError(
+            "The cshape of sample_pressure and reference_pressure must be "
+            "broadcastable except for the last dimension")
+
+    # prepare data
+    microphone_weights = microphone_weights[:, np.newaxis]
+    p_sample = sample_pressure.freq
+    p_reference = reference_pressure.freq
+
+    # calculate according to mommertz correlation method Equation (5)
+    p_sample_sum = np.sum(microphone_weights * np.abs(p_sample)**2, axis=-2)
+    p_ref_sum = np.sum(microphone_weights * np.abs(p_reference)**2, axis=-2)
+    p_cross_sum = np.sum(
+            p_sample * np.conj(p_reference) * microphone_weights, axis=-2)
+
+    data_scattering_coefficient \
+        = 1 - ((np.abs(p_cross_sum)**2)/(p_sample_sum*p_ref_sum))
+
+    # create pyfar.FrequencyData object
+    scattering_coefficients = pf.FrequencyData(
+        data_scattering_coefficient,
+        sample_pressure.frequencies)
+
+    return scattering_coefficients
diff --git a/tests/test_scattering_freefield.py b/tests/test_scattering_freefield.py
new file mode 100644
index 0000000..960a7a3
--- /dev/null
+++ b/tests/test_scattering_freefield.py
@@ -0,0 +1,138 @@
+import pyfar as pf
+import numpy as np
+import numpy.testing as npt
+import pytest
+from imkar.scattering import freefield as sff
+
+
+def plane_wave(amplitude, direction, sampling):
+    f = 5000
+    c = 343
+    x = sampling
+    direction.cartesian = direction.cartesian/direction.radius
+    dot_product = direction.x*x.x+direction.y*x.y+direction.z*x.z
+    dot_product = dot_product[..., np.newaxis]
+    f = np.atleast_1d(f)
+    return pf.FrequencyData(
+        amplitude*np.exp(-1j*2*np.pi*f/c*dot_product),
+        frequencies=f,
+    )
+
+
+def test_correlation_method_0():
+    sampling = pf.samplings.sph_equal_area(5000)
+    sampling.weights = pf.samplings.calculate_sph_voronoi_weights(sampling)
+    sampling = sampling[sampling.z>0]
+    sample_pressure = plane_wave(1, pf.Coordinates(0, 0, 1), sampling)
+    reference_pressure = plane_wave(1, pf.Coordinates(0, 0, 1), sampling)
+    s = sff.correlation_method(
+        sample_pressure, reference_pressure, sampling.weights,
+    )
+    npt.assert_almost_equal(s.freq, 0)
+
+
+@pytest.mark.parametrize("s_scatter", [0.1, 0.5, 0.9])
+@pytest.mark.parametrize("Phi_scatter_deg", [30, 60, 90, 120, 150, 42])
+def test_correlation_method_with_plane_waves(s_scatter, Phi_scatter_deg):
+    s_spec = 1-s_scatter
+    Phi_spec = 45/180*np.pi
+    Phi_scatter = Phi_scatter_deg/180*np.pi
+    R_spec = np.sqrt(s_spec)
+    R_scatter = np.sqrt(np.abs(s_scatter*np.sin(Phi_spec)/np.sin(Phi_scatter)))
+    sampling = pf.samplings.sph_equal_area(10000)
+    sampling.weights = pf.samplings.calculate_sph_voronoi_weights(sampling)
+    sampling = sampling[sampling.z>0]
+    sample_pressure = plane_wave(
+        R_spec,
+        pf.Coordinates.from_spherical_front(np.pi/2, Phi_spec, 1), sampling)
+    sample_pressure += plane_wave(
+        R_scatter,
+        pf.Coordinates.from_spherical_front(np.pi/2, Phi_scatter, 1), sampling)
+    reference_pressure = plane_wave(
+        1, pf.Coordinates.from_spherical_front(np.pi/2, Phi_spec, 1), sampling)
+    sd_spec = 1-sff.correlation_method(
+        sample_pressure, reference_pressure, sampling.weights,
+    )
+    reference_pressure = plane_wave(
+        1, pf.Coordinates.from_spherical_front(
+            np.pi/2, Phi_scatter, 1), sampling)
+    sd_scatter = 1-sff.correlation_method(
+        sample_pressure, reference_pressure, sampling.weights,
+    )
+    npt.assert_almost_equal(sd_spec.freq, s_spec, 1)
+    npt.assert_almost_equal(sd_scatter.freq, s_scatter, 1)
+
+    reference_pressure = plane_wave(
+        1, pf.Coordinates.from_spherical_front(
+            np.pi/2, Phi_spec+5/180*np.pi, 1), sampling)
+    sd_scatter_0 = 1-sff.correlation_method(
+        sample_pressure, reference_pressure, sampling.weights,
+    )
+    npt.assert_almost_equal(sd_scatter_0.freq, 0, 1)
+
+
+def test_correlation_method():
+    sampling = pf.samplings.sph_equal_area(5000)
+    sampling.weights = pf.samplings.calculate_sph_voronoi_weights(sampling)
+    sampling = sampling[sampling.z>0]
+    sample_pressure = plane_wave(1, pf.Coordinates(0, 0, 1), sampling)
+    reference_pressure = plane_wave(1, pf.Coordinates(0, 0, 1), sampling)
+    s = sff.correlation_method(
+        sample_pressure, reference_pressure, sampling.weights,
+    )
+    npt.assert_almost_equal(s.freq, 0)
+
+
+def test_correlation_method_invalid_sample_pressure_type():
+    reference_pressure = pf.FrequencyData(np.array([1, 2, 3]), [100, 200, 300])
+    microphone_weights = np.array([0.5, 0.5, 0.5])
+    with pytest.raises(
+            TypeError, match="sample_pressure must be of type "
+            "pyfar.FrequencyData"):
+        sff.correlation_method(
+            "invalid_type", reference_pressure, microphone_weights)
+
+
+def test_correlation_method_invalid_reference_pressure_type():
+    sample_pressure = pf.FrequencyData(np.array([1, 2, 3]), [100, 200, 300])
+    microphone_weights = np.array([0.5, 0.5, 0.5])
+    with pytest.raises(
+            TypeError, match="reference_pressure must be of type "
+            "pyfar.FrequencyData"):
+        sff.correlation_method(
+            sample_pressure, "invalid_type", microphone_weights)
+
+
+def test_correlation_method_mismatched_sample_pressure_weights():
+    sample_pressure = pf.FrequencyData(np.array([[1, 2, 3]]), [100, 200, 300])
+    reference_pressure = pf.FrequencyData(
+        np.array([[1, 2, 3]]), [100, 200, 300])
+    microphone_weights = np.array([0.5, 0.5])
+    with pytest.raises(
+            ValueError, match="The last dimension of sample_pressure must "
+            "match the size of microphone_weights"):
+        sff.correlation_method(
+            sample_pressure, reference_pressure, microphone_weights)
+
+
+def test_correlation_method_mismatched_reference_pressure_weights():
+    sample_pressure = pf.FrequencyData(np.array([[1, 2, 3]]), [100, 200, 300])
+    reference_pressure = pf.FrequencyData(
+        np.array([[1, 2, 3]]), [100, 200, 300])
+    microphone_weights = np.array([0.5, 0.5])
+    with pytest.raises(
+            ValueError, match="The last dimension of sample_pressure must "
+            "match the size of microphone_weights"):
+        sff.correlation_method(
+            sample_pressure, reference_pressure, microphone_weights)
+
+
+def test_correlation_method_mismatched_sample_reference_shapes():
+    sample_pressure = pf.FrequencyData(np.array([[1, 2, 3]]), [100, 200, 300])
+    reference_pressure = pf.FrequencyData(np.array([[1, 2]]), [100, 200])
+    microphone_weights = np.array([0.5, 0.5, 0.5])
+    with pytest.raises(
+            ValueError, match="The last dimension of sample_pressure must "
+            "match the size of microphone_weights"):
+        sff.correlation_method(
+            sample_pressure, reference_pressure, microphone_weights)

From a865d2caf035256216eb3b66c29878994c7fdfb2 Mon Sep 17 00:00:00 2001
From: Anne Heimes <ahe@akustik.rwth-aachen.de>
Date: Fri, 25 Apr 2025 12:21:35 +0200
Subject: [PATCH 02/12] remove unused files

---
 docs/modules/imkar.rst | 7 -------
 1 file changed, 7 deletions(-)
 delete mode 100644 docs/modules/imkar.rst

diff --git a/docs/modules/imkar.rst b/docs/modules/imkar.rst
deleted file mode 100644
index 0ea5287..0000000
--- a/docs/modules/imkar.rst
+++ /dev/null
@@ -1,7 +0,0 @@
-imkar
-=====
-
-.. automodule:: imkar
-   :members:
-   :undoc-members:
-   :show-inheritance:

From ae1a7c398424ac397f1b6b03e3eea17ddad26e69 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Tue, 27 May 2025 01:07:51 +0200
Subject: [PATCH 03/12] Apply suggestions from code review

---
 imkar/scattering/freefield.py | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index e88224c..5b6d2e2 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -85,7 +85,9 @@ def correlation_method(
         raise ValueError(
             "The cshape of sample_pressure and reference_pressure must be "
             "broadcastable except for the last dimension")
-
+    # Test whether the objects are able to perform arithmetic operations.
+    # e.g. does the frequency vectors match
+    _ = sample_pressure + p_reference
     # prepare data
     microphone_weights = microphone_weights[:, np.newaxis]
     p_sample = sample_pressure.freq

From 026f8049a90c9680ff70ddc399e708f7a1b4a207 Mon Sep 17 00:00:00 2001
From: Anne Heimes <ahe@akustik.rwth-aachen.de>
Date: Tue, 27 May 2025 01:10:46 +0200
Subject: [PATCH 04/12] fix tests

---
 imkar/scattering/freefield.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index 5b6d2e2..bc9880c 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -87,7 +87,8 @@ def correlation_method(
             "broadcastable except for the last dimension")
     # Test whether the objects are able to perform arithmetic operations.
     # e.g. does the frequency vectors match
-    _ = sample_pressure + p_reference
+    _ = sample_pressure + reference_pressure
+
     # prepare data
     microphone_weights = microphone_weights[:, np.newaxis]
     p_sample = sample_pressure.freq

From fba1ada994aa7a59c589709362acecb7bad5a031 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Tue, 29 Jul 2025 10:30:44 +0200
Subject: [PATCH 05/12] apply review, enhance tests and tests doc

---
 imkar/scattering/freefield.py      | 59 +++++++++++++++---------------
 tests/test_scattering_freefield.py | 52 +++++++++++++++++++++++---
 2 files changed, 75 insertions(+), 36 deletions(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index bc9880c..ff5b5de 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -6,43 +6,41 @@
 def correlation_method(
         sample_pressure, reference_pressure, microphone_weights):
     r"""
-    Calculate the incident-dependent free-field scattering coefficient.
+    Calculate the scattering coefficient from free-field reflection pattern.
 
     This function uses the Mommertz correlation method [#]_ to compute the
-    scattering coefficient of the input data:
+    scattering coefficient by averaging over the microphone positions:
 
     .. math::
         s = 1 -
-            \frac{|\sum_w \underline{p}_{\text{sample}}(\vartheta,\varphi)
-            \cdot \underline{p}_{\text{reference}}^*(\vartheta,\varphi)
-            \cdot w(\vartheta,\varphi)|^2}
-            {\sum_w |\underline{p}_{\text{sample}}(\vartheta,\varphi)|^2
-            \cdot w(\vartheta,\varphi) \cdot \sum_w
-            |\underline{p}_{\text{reference}}(\vartheta,\varphi)|^2
-            \cdot w(\vartheta,\varphi) }
+            \frac{|\sum_{\Omega_r} \underline{p}_{\text{sample}}(\Omega_r)
+            \cdot \underline{p}_{\text{reference}}^*(\Omega_r)
+            \cdot w(\Omega_r)|^2}
+            {\sum_{\Omega_r} |\underline{p}_{\text{sample}}(\Omega_r)|^2
+            \cdot w(\Omega_r) \cdot \sum_{\Omega_r}
+            |\underline{p}_{\text{reference}}(\Omega_r)|^2
+            \cdot w(\Omega_r) }
 
     where:
         - :math:`\underline{p}_{\text{sample}}` is the reflected sound
           pressure of the sample under investigation.
         - :math:`\underline{p}_{\text{reference}}` is the reflected sound
           pressure from the reference sample.
-        - :math:`w` represents the area weights of the sampling, and
-          :math:`\vartheta` and :math:`\varphi` are the ``colatitude``
-          and ``azimuth`` angles from the
-          :py:class:`~pyfar.classes.coordinates.Coordinates` object.
-
-    The test sample is assumed to lie in the x-y-plane.
+        - :math:`\Omega_r` represents the solid angle of the microphone
+          positions and :math:`w(\Omega_r)` represents its area weights.
 
     Parameters
     ----------
-    sample_pressure : :py:class:`~pyfar.classes.audio.FrequencyData`
+    sample_pressure : :py:class:`~pyfar.FrequencyData`, :py:class:`~pyfar.Signal`
         Reflected sound pressure or directivity of the test sample. Its cshape
         must be (..., microphone_weights.size) and broadcastable to the
-        cshape of ``reference_pressure``.
-    reference_pressure : :py:class:`~pyfar.classes.audio.FrequencyData`
+        cshape of ``reference_pressure``. The frequency vectors of both
+        ``sample_pressure`` and ``reference_pressure`` must match.
+    reference_pressure : :py:class:`~pyfar.FrequencyData`, :py:class:`~pyfar.Signal`
         Reflected sound pressure or directivity of the reference sample. Its
         cshape must be (..., microphone_weights.size) and broadcastable to the
-        cshape of ``sample_pressure``.
+        cshape of ``sample_pressure``. The frequency vectors of both
+        ``sample_pressure`` and ``reference_pressure`` must match.
     microphone_weights : array_like
         1D array containing the area weights for the microphone positions.
         No normalization is required. Its shape must match the last dimension
@@ -51,7 +49,9 @@ def correlation_method(
     Returns
     -------
     scattering_coefficients : :py:class:`~pyfar.classes.audio.FrequencyData`
-        The scattering coefficient for each incident direction as a function
+        The scattering coefficient of the broadcasted shape of
+        ``sample_pressure`` and ``reference_pressure`` excluding the
+        last dimension.
         of frequency.
 
     References
@@ -63,11 +63,14 @@ def correlation_method(
 
     """
     # check input types
-    if not isinstance(sample_pressure, pf.FrequencyData):
-        raise TypeError("sample_pressure must be of type pyfar.FrequencyData")
-    if not isinstance(reference_pressure, pf.FrequencyData):
+    if not isinstance(sample_pressure, (pf.FrequencyData, pf.Signal)):
+        raise TypeError(
+            "sample_pressure must be of type pyfar.FrequencyData or "
+            "pyfar.Signal")
+    if not isinstance(reference_pressure, (pf.FrequencyData, pf.Signal)):
         raise TypeError(
-            "reference_pressure must be of type pyfar.FrequencyData")
+            "reference_pressure must be of type pyfar.FrequencyData or "
+            "pyfar.Signal")
     microphone_weights = np.atleast_1d(
         np.asarray(microphone_weights, dtype=float))
 
@@ -81,10 +84,6 @@ def correlation_method(
             "The last dimension of reference_pressure must match the size of "
             "microphone_weights")
 
-    if sample_pressure.cshape[:-1] != reference_pressure.cshape[:-1]:
-        raise ValueError(
-            "The cshape of sample_pressure and reference_pressure must be "
-            "broadcastable except for the last dimension")
     # Test whether the objects are able to perform arithmetic operations.
     # e.g. does the frequency vectors match
     _ = sample_pressure + reference_pressure
@@ -100,8 +99,8 @@ def correlation_method(
     p_cross_sum = np.sum(
             p_sample * np.conj(p_reference) * microphone_weights, axis=-2)
 
-    data_scattering_coefficient \
-        = 1 - ((np.abs(p_cross_sum)**2)/(p_sample_sum*p_ref_sum))
+    data_scattering_coefficient = 1 - np.abs(p_cross_sum)**2/(
+        p_sample_sum*p_ref_sum)
 
     # create pyfar.FrequencyData object
     scattering_coefficients = pf.FrequencyData(
diff --git a/tests/test_scattering_freefield.py b/tests/test_scattering_freefield.py
index 960a7a3..ddc73c1 100644
--- a/tests/test_scattering_freefield.py
+++ b/tests/test_scattering_freefield.py
@@ -6,6 +6,26 @@
 
 
 def plane_wave(amplitude, direction, sampling):
+    """Generate a plane wave for a given direction and sampling.
+
+    This function is just used for testing purposes, so the frequency is
+    fixed to 5000 Hz.
+
+    Parameters
+    ----------
+    amplitude : float
+        The amplitude of the plane wave.
+    direction : pf.Coordinates
+        The direction of the plane wave.
+    sampling : pf.Sampling
+        The sampling grid for the plane wave.
+
+    Returns
+    -------
+    pf.FrequencyData
+        The generated plane wave in the frequency domain.
+    """
+
     f = 5000
     c = 343
     x = sampling
@@ -19,7 +39,7 @@ def plane_wave(amplitude, direction, sampling):
     )
 
 
-def test_correlation_method_0():
+def test_correlation_zero_scattering():
     sampling = pf.samplings.sph_equal_area(5000)
     sampling.weights = pf.samplings.calculate_sph_voronoi_weights(sampling)
     sampling = sampling[sampling.z>0]
@@ -33,12 +53,23 @@ def test_correlation_method_0():
 
 @pytest.mark.parametrize("s_scatter", [0.1, 0.5, 0.9])
 @pytest.mark.parametrize("Phi_scatter_deg", [30, 60, 90, 120, 150, 42])
-def test_correlation_method_with_plane_waves(s_scatter, Phi_scatter_deg):
+def test_correlation_fractional_scattering(s_scatter, Phi_scatter_deg):
+    """
+    Test analytic scattering coefficient calculation for non 0 or 1.
+
+    for the sample pressure, two plane waves are used. The scattered and the
+    specular wave. The reflection factor of both waves is calculated based on
+    the scattering coefficient and the scattering angle.
+    """
+    # calculate the specular and scattered reflection factors
     s_spec = 1-s_scatter
     Phi_spec = 45/180*np.pi
     Phi_scatter = Phi_scatter_deg/180*np.pi
     R_spec = np.sqrt(s_spec)
     R_scatter = np.sqrt(np.abs(s_scatter*np.sin(Phi_spec)/np.sin(Phi_scatter)))
+
+    # define the sampling grid and generate the plane waves
+    # for the sample pressure and the reference pressure
     sampling = pf.samplings.sph_equal_area(10000)
     sampling.weights = pf.samplings.calculate_sph_voronoi_weights(sampling)
     sampling = sampling[sampling.z>0]
@@ -48,20 +79,29 @@ def test_correlation_method_with_plane_waves(s_scatter, Phi_scatter_deg):
     sample_pressure += plane_wave(
         R_scatter,
         pf.Coordinates.from_spherical_front(np.pi/2, Phi_scatter, 1), sampling)
+
+    # Calculate the scattering coefficient with for the specular wave
+    # in other words, the reference pressure is the specular wave (R=1)
     reference_pressure = plane_wave(
         1, pf.Coordinates.from_spherical_front(np.pi/2, Phi_spec, 1), sampling)
     sd_spec = 1-sff.correlation_method(
         sample_pressure, reference_pressure, sampling.weights,
     )
+    npt.assert_almost_equal(sd_spec.freq, s_spec, 1)
+
+    # Calculate the scattering coefficient with for the scattered wave
+    # in other words, the reference pressure is the scattered wave (R=1)
+    # to the result will be 1-s_scatter
     reference_pressure = plane_wave(
         1, pf.Coordinates.from_spherical_front(
             np.pi/2, Phi_scatter, 1), sampling)
     sd_scatter = 1-sff.correlation_method(
         sample_pressure, reference_pressure, sampling.weights,
     )
-    npt.assert_almost_equal(sd_spec.freq, s_spec, 1)
     npt.assert_almost_equal(sd_scatter.freq, s_scatter, 1)
 
+    # Calculate the scattering coefficient with for 5° more then the specular
+    # wave, than the scattering coefficient should be 0
     reference_pressure = plane_wave(
         1, pf.Coordinates.from_spherical_front(
             np.pi/2, Phi_spec+5/180*np.pi, 1), sampling)
@@ -71,16 +111,16 @@ def test_correlation_method_with_plane_waves(s_scatter, Phi_scatter_deg):
     npt.assert_almost_equal(sd_scatter_0.freq, 0, 1)
 
 
-def test_correlation_method():
+def test_correlation_one_scattering():
     sampling = pf.samplings.sph_equal_area(5000)
     sampling.weights = pf.samplings.calculate_sph_voronoi_weights(sampling)
     sampling = sampling[sampling.z>0]
-    sample_pressure = plane_wave(1, pf.Coordinates(0, 0, 1), sampling)
+    sample_pressure = plane_wave(1, pf.Coordinates(0, 1, 0), sampling)
     reference_pressure = plane_wave(1, pf.Coordinates(0, 0, 1), sampling)
     s = sff.correlation_method(
         sample_pressure, reference_pressure, sampling.weights,
     )
-    npt.assert_almost_equal(s.freq, 0)
+    npt.assert_almost_equal(s.freq, 1, decimal=3)
 
 
 def test_correlation_method_invalid_sample_pressure_type():

From 67de44575a60b02943481c510d5be3f7210e62e2 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Tue, 29 Jul 2025 10:32:50 +0200
Subject: [PATCH 06/12] fix ruff

---
 imkar/scattering/freefield.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index ff5b5de..3091dfc 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -61,7 +61,7 @@ def correlation_method(
             Acoustics, vol. 60, no. 2, pp. 201-203, June 2000,
             doi: 10.1016/S0003-682X(99)00057-2.
 
-    """
+    """  # noqa: E501
     # check input types
     if not isinstance(sample_pressure, (pf.FrequencyData, pf.Signal)):
         raise TypeError(

From 051fb263ab144ac76a761b036e09594e7a5df318 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Tue, 29 Jul 2025 10:35:41 +0200
Subject: [PATCH 07/12] improve docstring

---
 imkar/scattering/freefield.py | 5 ++---
 1 file changed, 2 insertions(+), 3 deletions(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index 3091dfc..d96fd02 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -49,10 +49,9 @@ def correlation_method(
     Returns
     -------
     scattering_coefficients : :py:class:`~pyfar.classes.audio.FrequencyData`
-        The scattering coefficient of the broadcasted shape of
-        ``sample_pressure`` and ``reference_pressure`` excluding the
+        The scattering coefficient of the broadcasted cshape of
+        ``sample_pressure`` and ``reference_pressure``, excluding the
         last dimension.
-        of frequency.
 
     References
     ----------

From 205f4818d435e12db4a5152e08dff7300dde55e1 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Tue, 29 Jul 2025 10:42:42 +0200
Subject: [PATCH 08/12] fix tests

---
 tests/test_scattering_freefield.py | 21 +++++----------------
 1 file changed, 5 insertions(+), 16 deletions(-)

diff --git a/tests/test_scattering_freefield.py b/tests/test_scattering_freefield.py
index ddc73c1..67671c9 100644
--- a/tests/test_scattering_freefield.py
+++ b/tests/test_scattering_freefield.py
@@ -144,9 +144,9 @@ def test_correlation_method_invalid_reference_pressure_type():
 
 
 def test_correlation_method_mismatched_sample_pressure_weights():
-    sample_pressure = pf.FrequencyData(np.array([[1, 2, 3]]), [100, 200, 300])
+    sample_pressure = pf.FrequencyData(np.array([[1, 2, 3]]).T, [100])
     reference_pressure = pf.FrequencyData(
-        np.array([[1, 2, 3]]), [100, 200, 300])
+        np.array([[1, 2]]).T, [100])
     microphone_weights = np.array([0.5, 0.5])
     with pytest.raises(
             ValueError, match="The last dimension of sample_pressure must "
@@ -156,23 +156,12 @@ def test_correlation_method_mismatched_sample_pressure_weights():
 
 
 def test_correlation_method_mismatched_reference_pressure_weights():
-    sample_pressure = pf.FrequencyData(np.array([[1, 2, 3]]), [100, 200, 300])
+    sample_pressure = pf.FrequencyData(np.array([[1, 2]]).T, [100])
     reference_pressure = pf.FrequencyData(
-        np.array([[1, 2, 3]]), [100, 200, 300])
+        np.array([[1, 2, 3]]).T, [100])
     microphone_weights = np.array([0.5, 0.5])
     with pytest.raises(
-            ValueError, match="The last dimension of sample_pressure must "
-            "match the size of microphone_weights"):
-        sff.correlation_method(
-            sample_pressure, reference_pressure, microphone_weights)
-
-
-def test_correlation_method_mismatched_sample_reference_shapes():
-    sample_pressure = pf.FrequencyData(np.array([[1, 2, 3]]), [100, 200, 300])
-    reference_pressure = pf.FrequencyData(np.array([[1, 2]]), [100, 200])
-    microphone_weights = np.array([0.5, 0.5, 0.5])
-    with pytest.raises(
-            ValueError, match="The last dimension of sample_pressure must "
+            ValueError, match="The last dimension of reference_pressure must "
             "match the size of microphone_weights"):
         sff.correlation_method(
             sample_pressure, reference_pressure, microphone_weights)

From 0e7714a4efbf8cdd8628ccd3e2fbd79aabfaa5b5 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Tue, 29 Jul 2025 12:07:24 +0200
Subject: [PATCH 09/12] make pyfar.FrequencyData link to pyfar doc

---
 imkar/scattering/freefield.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index d96fd02..6955f3f 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -31,12 +31,12 @@ def correlation_method(
 
     Parameters
     ----------
-    sample_pressure : :py:class:`~pyfar.FrequencyData`, :py:class:`~pyfar.Signal`
+    sample_pressure : pyfar.FrequencyData, pyfar.Signal
         Reflected sound pressure or directivity of the test sample. Its cshape
         must be (..., microphone_weights.size) and broadcastable to the
         cshape of ``reference_pressure``. The frequency vectors of both
         ``sample_pressure`` and ``reference_pressure`` must match.
-    reference_pressure : :py:class:`~pyfar.FrequencyData`, :py:class:`~pyfar.Signal`
+    reference_pressure : pyfar.FrequencyData, pyfar.Signal
         Reflected sound pressure or directivity of the reference sample. Its
         cshape must be (..., microphone_weights.size) and broadcastable to the
         cshape of ``sample_pressure``. The frequency vectors of both
@@ -48,7 +48,7 @@ def correlation_method(
 
     Returns
     -------
-    scattering_coefficients : :py:class:`~pyfar.classes.audio.FrequencyData`
+    scattering_coefficients : pyfar.FrequencyData
         The scattering coefficient of the broadcasted cshape of
         ``sample_pressure`` and ``reference_pressure``, excluding the
         last dimension.
@@ -60,7 +60,7 @@ def correlation_method(
             Acoustics, vol. 60, no. 2, pp. 201-203, June 2000,
             doi: 10.1016/S0003-682X(99)00057-2.
 
-    """  # noqa: E501
+    """
     # check input types
     if not isinstance(sample_pressure, (pf.FrequencyData, pf.Signal)):
         raise TypeError(

From 06581bafd3891372e406968012677f6f07e73170 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Tue, 29 Jul 2025 12:12:48 +0200
Subject: [PATCH 10/12] Apply suggestions from code review

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>
---
 tests/test_scattering_freefield.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/tests/test_scattering_freefield.py b/tests/test_scattering_freefield.py
index 67671c9..2923a63 100644
--- a/tests/test_scattering_freefield.py
+++ b/tests/test_scattering_freefield.py
@@ -57,7 +57,7 @@ def test_correlation_fractional_scattering(s_scatter, Phi_scatter_deg):
     """
     Test analytic scattering coefficient calculation for non 0 or 1.
 
-    for the sample pressure, two plane waves are used. The scattered and the
+    For the sample pressure, two plane waves are used. The scattered and the
     specular wave. The reflection factor of both waves is calculated based on
     the scattering coefficient and the scattering angle.
     """
@@ -91,7 +91,7 @@ def test_correlation_fractional_scattering(s_scatter, Phi_scatter_deg):
 
     # Calculate the scattering coefficient with for the scattered wave
     # in other words, the reference pressure is the scattered wave (R=1)
-    # to the result will be 1-s_scatter
+    # so the result will be 1-s_scatter
     reference_pressure = plane_wave(
         1, pf.Coordinates.from_spherical_front(
             np.pi/2, Phi_scatter, 1), sampling)
@@ -100,8 +100,8 @@ def test_correlation_fractional_scattering(s_scatter, Phi_scatter_deg):
     )
     npt.assert_almost_equal(sd_scatter.freq, s_scatter, 1)
 
-    # Calculate the scattering coefficient with for 5° more then the specular
-    # wave, than the scattering coefficient should be 0
+    # Calculate the scattering coefficient for 5° more than the specular
+    # wave, then the scattering coefficient should be 0
     reference_pressure = plane_wave(
         1, pf.Coordinates.from_spherical_front(
             np.pi/2, Phi_spec+5/180*np.pi, 1), sampling)

From 51e459b008aa3cf19337c2c58e4ee1a71d2a3872 Mon Sep 17 00:00:00 2001
From: Anne Heimes <64446926+ahms5@users.noreply.github.com>
Date: Fri, 22 Aug 2025 10:39:48 +0200
Subject: [PATCH 11/12] Apply suggestions from code review

Co-authored-by: sikersten <70263411+sikersten@users.noreply.github.com>
---
 imkar/scattering/freefield.py      | 16 ++++++++--------
 tests/test_scattering_freefield.py |  2 +-
 2 files changed, 9 insertions(+), 9 deletions(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index 6955f3f..6fb0477 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -27,29 +27,29 @@ def correlation_method(
         - :math:`\underline{p}_{\text{reference}}` is the reflected sound
           pressure from the reference sample.
         - :math:`\Omega_r` represents the solid angle of the microphone
-          positions and :math:`w(\Omega_r)` represents its area weights.
+          positions and :math:`w(\Omega_r)` represents its area weight.
 
     Parameters
     ----------
     sample_pressure : pyfar.FrequencyData, pyfar.Signal
-        Reflected sound pressure or directivity of the test sample. Its cshape
-        must be (..., microphone_weights.size) and broadcastable to the
-        cshape of ``reference_pressure``. The frequency vectors of both
+        Reflected sound pressure or directivity of the test sample. Its `cshape`
+        must be ``(..., microphone_weights.size)`` and broadcastable to the
+        `cshape` of ``reference_pressure``. The frequency vectors of both
         ``sample_pressure`` and ``reference_pressure`` must match.
     reference_pressure : pyfar.FrequencyData, pyfar.Signal
         Reflected sound pressure or directivity of the reference sample. Its
-        cshape must be (..., microphone_weights.size) and broadcastable to the
-        cshape of ``sample_pressure``. The frequency vectors of both
+        `cshape` must be (..., microphone_weights.size) and broadcastable to the
+        `cshape` of ``sample_pressure``. The frequency vectors of both
         ``sample_pressure`` and ``reference_pressure`` must match.
     microphone_weights : array_like
         1D array containing the area weights for the microphone positions.
         No normalization is required. Its shape must match the last dimension
-        in the cshape of ``sample_pressure`` and ``reference_pressure``.
+        in the `cshape` of ``sample_pressure`` and ``reference_pressure``.
 
     Returns
     -------
     scattering_coefficients : pyfar.FrequencyData
-        The scattering coefficient of the broadcasted cshape of
+        The scattering coefficient of the broadcasted `cshape` of
         ``sample_pressure`` and ``reference_pressure``, excluding the
         last dimension.
 
diff --git a/tests/test_scattering_freefield.py b/tests/test_scattering_freefield.py
index 2923a63..e630cbc 100644
--- a/tests/test_scattering_freefield.py
+++ b/tests/test_scattering_freefield.py
@@ -89,7 +89,7 @@ def test_correlation_fractional_scattering(s_scatter, Phi_scatter_deg):
     )
     npt.assert_almost_equal(sd_spec.freq, s_spec, 1)
 
-    # Calculate the scattering coefficient with for the scattered wave
+    # Calculate the scattering coefficient for the scattered wave
     # in other words, the reference pressure is the scattered wave (R=1)
     # so the result will be 1-s_scatter
     reference_pressure = plane_wave(

From 232bc973e72f47d9543d6009b927ab810ecffc33 Mon Sep 17 00:00:00 2001
From: Anne Heimes <ahe@akustik.rwth-aachen.de>
Date: Fri, 22 Aug 2025 11:36:24 +0200
Subject: [PATCH 12/12] apply review

---
 imkar/scattering/freefield.py      | 21 +++++++++++----------
 tests/test_scattering_freefield.py |  1 +
 2 files changed, 12 insertions(+), 10 deletions(-)

diff --git a/imkar/scattering/freefield.py b/imkar/scattering/freefield.py
index 6fb0477..907b32b 100644
--- a/imkar/scattering/freefield.py
+++ b/imkar/scattering/freefield.py
@@ -32,14 +32,14 @@ def correlation_method(
     Parameters
     ----------
     sample_pressure : pyfar.FrequencyData, pyfar.Signal
-        Reflected sound pressure or directivity of the test sample. Its `cshape`
-        must be ``(..., microphone_weights.size)`` and broadcastable to the
-        `cshape` of ``reference_pressure``. The frequency vectors of both
-        ``sample_pressure`` and ``reference_pressure`` must match.
+        Reflected sound pressure or directivity of the test sample. Its
+        `cshape` must be ``(..., microphone_weights.size)`` and broadcastable
+        to the `cshape` of ``reference_pressure``. The frequency vectors
+        of both ``sample_pressure`` and ``reference_pressure`` must match.
     reference_pressure : pyfar.FrequencyData, pyfar.Signal
         Reflected sound pressure or directivity of the reference sample. Its
-        `cshape` must be (..., microphone_weights.size) and broadcastable to the
-        `cshape` of ``sample_pressure``. The frequency vectors of both
+        `cshape` must be (..., microphone_weights.size) and broadcastable
+        to the `cshape` of ``sample_pressure``. The frequency vectors of both
         ``sample_pressure`` and ``reference_pressure`` must match.
     microphone_weights : array_like
         1D array containing the area weights for the microphone positions.
@@ -82,10 +82,11 @@ def correlation_method(
         raise ValueError(
             "The last dimension of reference_pressure must match the size of "
             "microphone_weights")
-
-    # Test whether the objects are able to perform arithmetic operations.
-    # e.g. does the frequency vectors match
-    _ = sample_pressure + reference_pressure
+    # check frequency vectors
+    if not np.allclose(
+            reference_pressure.frequencies, sample_pressure.frequencies,
+            atol=1e-15):
+        raise ValueError("The frequencies do not match.")
 
     # prepare data
     microphone_weights = microphone_weights[:, np.newaxis]
diff --git a/tests/test_scattering_freefield.py b/tests/test_scattering_freefield.py
index e630cbc..1cc240f 100644
--- a/tests/test_scattering_freefield.py
+++ b/tests/test_scattering_freefield.py
@@ -112,6 +112,7 @@ def test_correlation_fractional_scattering(s_scatter, Phi_scatter_deg):
 
 
 def test_correlation_one_scattering():
+    """Two plane waves in the different directions should result in s=1."""
     sampling = pf.samplings.sph_equal_area(5000)
     sampling.weights = pf.samplings.calculate_sph_voronoi_weights(sampling)
     sampling = sampling[sampling.z>0]