Update wfs-referencing.ipynb

- nicer code structure, format
- subplot(2,1) instead of (1,2) seems nicer
- scatter plot for ref contour
- tapering not needed
- add difference as magnitude plot

Author: fs446

doc/examples/wfs-referencing.ipynb

@@ -9,6 +9,17 @@
     "Illustrates the usage of the SFS toolbox for the simulation of different 2.5D WFS referencing schemes for contours within the listening area that exhibit amplitude correct synthesis, cf. Ch. 4.1.3 in <cite data-cite=\"Firtha2019\">[Firtha19]</cite> and <cite data-cite=\"Firtha2017\">(FFSS17)</cite>."
    ]
   },
+  {
+   "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": {},
@@ -22,19 +33,12 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "import matplotlib.pyplot as plt\n",
-    "import numpy as np\n",
-    "import sfs\n",
-    "\n",
-    "R = 1.5  # Radius of circular loudspeaker array\n",
-    "array = sfs.array.circular(N=64, R=R)\n",
+    "R = 1.5  # Radius of circular loudspeaker array, in m\n",
+    "array = sfs.array.circular(N=64, R=R)  # with N loudspeakers\n",
     "grid = sfs.util.xyz_grid([-2, 2], [-2, 2], 0, spacing=0.02)\n",
     "\n",
-    "xs = sfs.util.asarray_of_rows((-4, 0, 0))  # point source on negative x-axis\n",
-    "normalize_gain = 4 * np.pi * np.linalg.norm(xs)\n",
-    "wavelength = 1 / 4  # m\n",
-    "f = sfs.default.c / wavelength  # Hz\n",
-    "omega = 2 * np.pi * f  # rad/s"
+    "xs = -4, 0, 0  # point source on negative x-axis\n",
+    "wavelength = 1 / 4  # m"
    ]
   },
   {
@@ -43,32 +47,40 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "def sound_field(d, xref, selection,\n",
-    "                secondary_source, array, grid, tapering=True):\n",
-    "    if tapering:\n",
-    "        tapering_window = sfs.tapering.kaiser(selection, beta=1)\n",
-    "    else:\n",
-    "        tapering_window = sfs.tapering.none(selection)\n",
+    "def sound_field(\n",
+    "        d, xref, selection, secondary_source, array, grid):\n",
     "\n",
-    "    p = sfs.fd.synthesize(d, tapering_window,\n",
-    "                          array, secondary_source, grid=grid)\n",
+    "    p = sfs.fd.synthesize(\n",
+    "        d, selection, array, secondary_source, grid=grid)\n",
     "\n",
-    "    fig, axs = plt.subplots(1, 2, figsize=(10, 8))\n",
-    "    sfs.plot2d.amplitude(p, grid, vmax=2, vmin=-2, ax=axs[0])\n",
-    "    sfs.plot2d.level(p, grid, vmax=6, vmin=-6, ax=axs[1],\n",
-    "                     cmap='seismic', colorbar_kwargs={'label': 'dB'})\n",
-    "    for i in range(axs.shape[0]):\n",
-    "        sfs.plot2d.loudspeakers(array.x, array.n,\n",
-    "                                tapering_window,\n",
-    "                                size=0.125, ax=axs[i])\n",
-    "        axs[i].plot(xref[:, 0][selection],\n",
-    "                    xref[:, 1][selection],\n",
-    "                    marker='o', linestyle='', color='dimgray', ms=4)\n",
-    "        axs[i].grid(True)\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(\n",
+    "        p, grid, vmax=6, vmin=-6, ax=ax_lvl, cmap='seismic',\n",
+    "        colorbar_kwargs={'label': 'dB'})\n",
+    "    for ax in ax_amp, ax_lvl:\n",
+    "        sfs.plot2d.loudspeakers(\n",
+    "            array.x, array.n, selection, size=0.125, ax=ax)\n",
+    "        ax.scatter(*xref[selection, :2].T, s=20, c='dimgray',\n",
+    "                   zorder=3)\n",
+    "        ax.grid(True)\n",
     "    plt.tight_layout()\n",
     "    return p\n"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "xs = sfs.util.asarray_of_rows((-4, 0, 0))\n",
+    "f = sfs.default.c / wavelength  # Hz\n",
+    "omega = 2 * np.pi * f  # rad/s\n",
+    "normalize_gain = 4 * np.pi * np.linalg.norm(xs)"
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {},
@@ -89,13 +101,14 @@
     "# calc reference contour xref(x0), cf. [FFSS17, eq. (24), (31), (52)]\n",
     "# this code assumes virtual point source on x-axis\n",
     "cosbeta = (array.n @ [1, 0, 0]).reshape(-1, 1)\n",
-    "xref = array.x + (xs - array.x) * \\\n",
-    "    (xref_line + R * cosbeta) / (xs[0, 0] + R * cosbeta)\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(d * normalize_gain, xref, selection,\n",
-    "                     secondary_source, array, grid, tapering=False)"
+    "p_line = sound_field(\n",
+    "    d * normalize_gain, xref, selection,\n",
+    "    secondary_source, array, grid)"
    ]
   },
   {
@@ -123,8 +136,9 @@
     "\n",
     "d, selection, secondary_source = sfs.fd.wfs.point_25d(\n",
     "    omega, array.x, array.n, xs, xref=xref)\n",
-    "p_circ = sound_field(d * normalize_gain, xref, selection,\n",
-    "                     secondary_source, array, grid, tapering=False)"
+    "p_circ = sound_field(\n",
+    "    d * normalize_gain, xref, selection,\n",
+    "    secondary_source, array, grid)"
    ]
   },
   {
@@ -134,10 +148,26 @@
    "outputs": [],
    "source": [
     "# (complex-valued) ratio between the two soundfields tells\n",
-    "# us about the difference of both approaches\n",
-    "sfs.plot2d.level(p_line / p_circ, grid, vmax=1, vmin=-1,\n",
-    "                 cmap='seismic', colorbar_kwargs={'label': 'dB'})\n",
-    "sfs.plot2d.loudspeakers(array.x, array.n, size=0.125)\n",
+    "# us about the difference of both approaches in terms of dB\n",
+    "sfs.plot2d.level(\n",
+    "    p_line / p_circ, grid, vmax=1, vmin=-1, cmap='seismic',\n",
+    "    colorbar_kwargs={'label': 'dB'})\n",
+    "sfs.plot2d.loudspeakers(array.x, array.n, selection, size=0.125)\n",
+    "plt.grid(True)\n",
+    "plt.tight_layout()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# (complex-valued) difference between the two soundfields tells\n",
+    "# us about the difference of both approaches in terms of magnitude\n",
+    "sfs.plot2d.amplitude(\n",
+    "    p_line - p_circ, grid, vmax=1e-1, vmin=-1e-1)\n",
+    "sfs.plot2d.loudspeakers(array.x, array.n, selection, size=0.125)\n",
     "plt.grid(True)\n",
     "plt.tight_layout()"
    ]