diff --git a/doc/examples.rst b/doc/examples.rst index c04994c4..80d60fe0 100644 --- a/doc/examples.rst +++ b/doc/examples.rst @@ -13,6 +13,7 @@ Examples .. nbgallery:: examples/sound-field-synthesis + examples/wfs-referencing examples/modal-room-acoustics examples/mirror-image-source-model examples/animations-pulsating-sphere diff --git a/doc/examples/wfs-referencing.ipynb b/doc/examples/wfs-referencing.ipynb new file mode 100644 index 00000000..ea406a40 --- /dev/null +++ b/doc/examples/wfs-referencing.ipynb @@ -0,0 +1,252 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# 2.5D WFS Referencing Schemes\n", + "\n", + "This notebook illustrates the usage of the SFS toolbox for the simulation of different 2.5D WFS referencing schemes.\n", + "A dedicated referencing scheme allows correct amplitude alongside a reference contour within the listening area.\n", + "For the theory please check\n", + "Ch 3.1-3.3 in (Start1997),\n", + "Ch. 4.1.3 in (Firtha2019) and\n", + "(Firtha2017)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "import sfs" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Circular loudspeaker arrays" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "R = 1.5 # radius [m] of circular loudspeaker array\n", + "N = 64 # loudspeakers\n", + "array = sfs.array.circular(N=N, R=R)\n", + "grid = sfs.util.xyz_grid([-2, 2], [-2, 2], 0, spacing=0.02)\n", + "\n", + "xs = -4, 0, 0 # virtual point source on negative x-axis\n", + "wavelength = 1 / 4 # m" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def sound_field(d, selection, array, secondary_source, grid, xref):\n", + " p = sfs.fd.synthesize(d, selection, array, secondary_source, grid=grid)\n", + " fig, [ax_amp, ax_lvl] = plt.subplots(2, 1, sharex=True)\n", + " fig.set_figheight(fig.get_figwidth() * 3/2)\n", + " sfs.plot2d.amplitude(p, grid, vmax=2, vmin=-2, ax=ax_amp)\n", + " sfs.plot2d.level(p, grid, vmax=12, vmin=-12, ax=ax_lvl)\n", + " sfs.plot2d.level_contour(p, grid, levels=[0], colors='w', ax=ax_lvl)\n", + " xref = np.broadcast_to(xref, array.x.shape)\n", + " for ax in ax_amp, ax_lvl:\n", + " sfs.plot2d.loudspeakers(array.x, array.n, selection, size=0.125, ax=ax)\n", + " ax_lvl.scatter(*xref[selection, :2].T, marker='o', s=20, c='lightsalmon',\n", + " zorder=3)\n", + " plt.tight_layout()\n", + " return p" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "xs = sfs.util.asarray_of_rows(xs)\n", + "frequency = sfs.default.c / wavelength # Hz\n", + "omega = 2 * np.pi * frequency # rad/s\n", + "normalize_gain = 4 * np.pi * np.linalg.norm(xs)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Line as reference contour\n", + "\n", + "The reference contour is calculated according to eqs. (24), (31), (52) in (Firtha2017). \n", + "The code assumes a virtual point source on x-axis.\n", + "The reference contour is a straight line on y-axis." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "xref_line = 0\n", + "cosbeta = (array.n @ [1, 0, 0]).reshape(-1, 1)\n", + "xref = array.x + \\\n", + " (xs - array.x) * (xref_line + R * cosbeta) / (xs[0, 0] + R * cosbeta)\n", + "\n", + "d, selection, secondary_source = sfs.fd.wfs.point_25d(\n", + " omega, array.x, array.n, xs, xref=xref)\n", + "p_line = sound_field(\n", + " d * normalize_gain, selection, array, secondary_source, grid, xref)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The level plot includes a white 0 dB isobar curve.\n", + "The orange-like dots represent the stationary phase points at which amplitude correct synthesis is to be expected.\n", + "These dots shape the line reference contour.\n", + "Note that the isobar curve is not perfectly aligned along line reference contour due to diffraction artifacts." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Circle as reference contour\n", + "\n", + "This reference contour is a circle with its origin at xs and a radius |xs|. This contour is obtained with more straightforward vector calculus than the previous example." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# reference contour is a circle with origin xs and radius |xs|\n", + "xref_dist = np.linalg.norm(xs)\n", + "# calc reference contour xref(x0), cf. [Firtha19, eq. (24), (31)]\n", + "xref = xs + xref_dist * sfs.util.normalize_rows(array.x - xs)\n", + "d, selection, secondary_source = sfs.fd.wfs.point_25d(\n", + " omega, array.x, array.n, xs, xref=xref)\n", + "p_circ = sound_field(\n", + " d * normalize_gain, selection, array, secondary_source, grid, xref)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Reference point" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The default handling in\n", + "`point_25d(omega, x0, n0, xs, xref=[0, 0, 0], c=None, omalias=None)`\n", + "uses just a reference point xref, and more specifically this default point is the origin of the coordinate system.\n", + "This single point xref, the virtual source position xs and the loudspeaker array geometry together determine the reference contour without further user access to it.\n", + "This handling is chosen due to convenience and practical relevance when working with circular loudspeaker arrays.\n", + "\n", + "The example below shows the resulting reference contour for the default case.\n", + "In the example it looks similar to the line reference contour, but is in general not exactly the same.\n", + "For example, please try a virtual point source that is far away from the array." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "d, selection, secondary_source = sfs.fd.wfs.point_25d(\n", + " omega, array.x, array.n, xs)\n", + "p_point = sound_field(\n", + " d * normalize_gain, selection, array, secondary_source,\n", + " grid, [0, 0, 0])" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Points with amplitude correct synthesis need to be stationary phase points, theoretically.\n", + "Within the listening area, these points are found on rays that start at the virtual point source and intersect with active loudspeakers.\n", + "The chosen points together shall shape a smooth contour, i.e. the reference contour.\n", + "\n", + "The example below shows a reference point xref that does not meet any ray (the gray lines in the level plot) alongside the stationary phase holds with its corresponding loudspeaker.\n", + "\n", + "The single point referencing scheme results in 0 dB isobar curve that closely passes the chosen xref point.\n", + "In practice this typically works with sufficient precision once the position of xref is appropriately chosen (i.e. not too close, not too far, not to off-center from the active loudspeakers etc.)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "xref = 0, 0.1175, 0 # intentionally no stationary phase point\n", + "# we don't forget to normalize the point source's amplitude\n", + "# to this new reference point:\n", + "normalize_gain = 4 * np.pi * np.linalg.norm(xs - xref)\n", + "d, selection, secondary_source = sfs.fd.wfs.point_25d(\n", + " omega, array.x, array.n, xs, xref=xref)\n", + "p_point = sound_field(\n", + " d * normalize_gain, selection, array, secondary_source,\n", + " grid, xref)\n", + "\n", + "# plot stationary phase rays\n", + "# one ray connects the virtual source with one activate loudspeaker\n", + "spa = array.x + 3*R * sfs.util.normalize_rows(array.x - xs)\n", + "plt.plot(\n", + " np.vstack((array.x[selection, 0], spa[selection, 0])),\n", + " np.vstack((array.x[selection, 1], spa[selection, 1])),\n", + " color='gray')\n", + "plt.xlim(-2, 2)\n", + "plt.ylim(-2, 2);" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "A plane wave like sound field, e.g. by setting `xs = -100, 0, 0`, for all above examples reveals some further interesting implications of the different referencing schemes." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "sfs", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.13.7" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/doc/references.bib b/doc/references.bib index 8b149bac..09a92a15 100644 --- a/doc/references.bib +++ b/doc/references.bib @@ -96,3 +96,10 @@ @phdthesis{Schultz2016 year = {2016}, doi = {10.18453/rosdok_id00001765} } +@phdthesis{Firtha2019, + author = {Firtha, G.}, + title = {{A Generalized Wave Field Synthesis Framework with Application + for Moving Virtual Sources}}, + school = {Budapest University of Technology and Economics}, + year = {2019} +} diff --git a/sfs/fd/wfs.py b/sfs/fd/wfs.py index 840b59e9..c6132d22 100644 --- a/sfs/fd/wfs.py +++ b/sfs/fd/wfs.py @@ -205,7 +205,7 @@ def point_25d(omega, x0, n0, xs, xref=[0, 0, 0], c=None, omalias=None): is implemented. The theoretical link of `point_25d()` and `point_25d_legacy()` was - introduced as *unified WFS framework* in :cite:`Firtha2017`. + introduced as *unified WFS framework* in :cite:`Firtha2017`, :cite:`Firtha2019`. Examples -------- @@ -217,6 +217,12 @@ def point_25d(omega, x0, n0, xs, xref=[0, 0, 0], c=None, omalias=None): normalize_gain = 4 * np.pi * np.linalg.norm(xs) plot(normalize_gain * d, selection, secondary_source) + .. nblinkgallery:: + :caption: Further Example + :name: wfs-referencing-link-gallery + + examples/wfs-referencing + """ x0 = _util.asarray_of_rows(x0) n0 = _util.asarray_of_rows(n0) @@ -292,7 +298,7 @@ def point_25d_legacy(omega, x0, n0, xs, xref=[0, 0, 0], c=None, omalias=None): \e{-\i\wc |\x_0-\x_\text{s}|} The theoretical link of `point_25d()` and `point_25d_legacy()` was - introduced as *unified WFS framework* in :cite:`Firtha2017`. + introduced as *unified WFS framework* in :cite:`Firtha2017`, :cite:`Firtha2019`. Also cf. Eq. (2.145)-(2.147) :cite:`Schultz2016`. Examples