Merge pull request #101 from sfstoolbox/synthesize

Introduce synthesize() functions for soundfield superposition

doc/examples/horizontal_plane_arrays.py

@@ -29,14 +29,13 @@ def compute_and_plot_soundfield(title):
     """Compute and plot synthesized sound field."""
     print('Computing', title)
 
-    twin = tapering(a, talpha)
-    p = sfs.mono.synthesized.generic(omega, x0, n0, d * twin * a0, grid,
-                                     source=sourcetype)
+    twin = tapering(selection, talpha)
+    p = sfs.mono.synthesize(d, twin, array, secondary_source, grid=grid)
 
     plt.figure(figsize=(15, 15))
     plt.cla()
     sfs.plot.soundfield(p, grid, xnorm)
-    sfs.plot.loudspeaker_2d(x0, n0, twin)
+    sfs.plot.loudspeaker_2d(array.x, array.n, twin)
     sfs.plot.virtualsource_2d(xs)
     sfs.plot.virtualsource_2d([0, 0], npw, type='plane')
     plt.title(title)
@@ -45,110 +44,108 @@ def compute_and_plot_soundfield(title):
 
 
 # linear array, secondary point sources, virtual monopole
-x0, n0, a0 = sfs.array.linear(N, dx, center=acenter, orientation=anormal)
+array = sfs.array.linear(N, dx, center=acenter, orientation=anormal)
 
-sourcetype = sfs.mono.source.point
-a = sfs.util.source_selection_point(n0, x0, xs)
-
-d = sfs.mono.drivingfunction.wfs_3d_point(omega, x0, n0, xs)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_3d_point(
+    omega, array.x, array.n, xs)
 compute_and_plot_soundfield('linear_ps_wfs_3d_point')
 
-d = sfs.mono.drivingfunction.wfs_25d_point(omega, x0, n0, xs, xref=xnorm)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_point(
+    omega, array.x, array.n, xs, xref=xnorm)
 compute_and_plot_soundfield('linear_ps_wfs_25d_point')
 
-d = sfs.mono.drivingfunction.wfs_2d_point(omega, x0, n0, xs)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_2d_point(
+    omega, array.x, array.n, xs)
 compute_and_plot_soundfield('linear_ps_wfs_2d_point')
 
 # linear array, secondary line sources, virtual line source
-sourcetype = sfs.mono.source.line
-d = sfs.mono.drivingfunction.wfs_2d_line(omega, x0, n0, xs)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_2d_line(
+    omega, array.x, array.n, xs)
 compute_and_plot_soundfield('linear_ls_wfs_2d_line')
 
 
 # linear array, secondary point sources, virtual plane wave
-sourcetype = sfs.mono.source.point
-a = sfs.util.source_selection_plane(n0, npw)
-
-d = sfs.mono.drivingfunction.wfs_3d_plane(omega, x0, n0, npw)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_3d_plane(
+    omega, array.x, array.n, npw)
 compute_and_plot_soundfield('linear_ps_wfs_3d_plane')
 
-d = sfs.mono.drivingfunction.wfs_25d_plane(omega, x0, n0, npw, xref=xnorm)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_plane(
+    omega, array.x, array.n, npw, xref=xnorm)
 compute_and_plot_soundfield('linear_ps_wfs_25d_plane')
 
-d = sfs.mono.drivingfunction.wfs_2d_plane(omega, x0, n0, npw)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_2d_plane(
+    omega, array.x, array.n, npw)
 compute_and_plot_soundfield('linear_ps_wfs_2d_plane')
 
 
 # non-uniform linear array, secondary point sources
-x0, n0, a0 = sfs.array.linear_diff(N//3 * [dx] + N//3 * [dx/2] + N//3 * [dx],
+array = sfs.array.linear_diff(N//3 * [dx] + N//3 * [dx/2] + N//3 * [dx],
                                    center=acenter, orientation=anormal)
 
-d = sfs.mono.drivingfunction.wfs_25d_point(omega, x0, n0, xs, xref=xnorm)
-a = sfs.util.source_selection_point(n0, x0, xs)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_point(
+    omega, array.x, array.n, xs, xref=xnorm)
 compute_and_plot_soundfield('linear_nested_ps_wfs_25d_point')
 
-d = sfs.mono.drivingfunction.wfs_25d_plane(omega, x0, n0, npw, xref=xnorm)
-a = sfs.util.source_selection_plane(n0, npw)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_plane(
+    omega, array.x, array.n, npw, xref=xnorm)
 compute_and_plot_soundfield('linear_nested_ps_wfs_25d_plane')
 
 
 # random sampled linear array, secondary point sources
-x0, n0, a0 = sfs.array.linear_random(N, dx/2, 1.5*dx, center=acenter,
+array = sfs.array.linear_random(N, dx/2, 1.5*dx, center=acenter,
                                      orientation=anormal)
 
-d = sfs.mono.drivingfunction.wfs_25d_point(omega, x0, n0, xs, xref=xnorm)
-a = sfs.util.source_selection_point(n0, x0, xs)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_point(
+    omega, array.x, array.n, xs, xref=xnorm)
 compute_and_plot_soundfield('linear_random_ps_wfs_25d_point')
 
-d = sfs.mono.drivingfunction.wfs_25d_plane(omega, x0, n0, npw, xref=xnorm)
-a = sfs.util.source_selection_plane(n0, npw)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_plane(
+    omega, array.x, array.n, npw, xref=xnorm)
 compute_and_plot_soundfield('linear_random_ps_wfs_25d_plane')
 
 
 # rectangular array, secondary point sources
-x0, n0, a0 = sfs.array.rectangular((N, N//2), dx, center=acenter, orientation=anormal)
-d = sfs.mono.drivingfunction.wfs_25d_point(omega, x0, n0, xs, xref=xnorm)
-a = sfs.util.source_selection_point(n0, x0, xs)
+array = sfs.array.rectangular((N, N//2), dx, center=acenter, orientation=anormal)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_point(
+    omega, array.x, array.n, xs, xref=xnorm)
 compute_and_plot_soundfield('rectangular_ps_wfs_25d_point')
 
-d = sfs.mono.drivingfunction.wfs_25d_plane(omega, x0, n0, npw, xref=xnorm)
-a = sfs.util.source_selection_plane(n0, npw)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_plane(
+    omega, array.x, array.n, npw, xref=xnorm)
 compute_and_plot_soundfield('rectangular_ps_wfs_25d_plane')
 
 
 # circular array, secondary point sources
 N = 60
-x0, n0, a0 = sfs.array.circular(N, 1, center=acenter)
-d = sfs.mono.drivingfunction.wfs_25d_point(omega, x0, n0, xs, xref=xnorm)
-a = sfs.util.source_selection_point(n0, x0, xs)
+array = sfs.array.circular(N, 1, center=acenter)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_point(
+    omega, array.x, array.n, xs, xref=xnorm)
 compute_and_plot_soundfield('circular_ps_wfs_25d_point')
 
-d = sfs.mono.drivingfunction.wfs_25d_plane(omega, x0, n0, npw, xref=xnorm)
-a = sfs.util.source_selection_plane(n0, npw)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_plane(
+    omega, array.x, array.n, npw, xref=xnorm)
 compute_and_plot_soundfield('circular_ps_wfs_25d_plane')
 
 
 # circular array, secondary line sources, NFC-HOA
-x0, n0, a0 = sfs.array.circular(N, 1)
+array = sfs.array.circular(N, 1)
 xnorm = [0, 0, 0]
 talpha = 0  # switches off tapering
-sourcetype = sfs.mono.source.line
 
-d = sfs.mono.drivingfunction.nfchoa_2d_plane(omega, x0, 1, npw)
-a = sfs.util.source_selection_all(N)
+d, selection, secondary_source = sfs.mono.drivingfunction.nfchoa_2d_plane(
+    omega, array.x, 1, npw)
 compute_and_plot_soundfield('circular_ls_nfchoa_2d_plane')
 
 
 # circular array, secondary point sources, NFC-HOA
-x0, n0, a0 = sfs.array.circular(N, 1)
+array = sfs.array.circular(N, 1)
 xnorm = [0, 0, 0]
 talpha = 0  # switches off tapering
-sourcetype = sfs.mono.source.point
 
-d = sfs.mono.drivingfunction.nfchoa_25d_point(omega, x0, 1, xs)
-a = sfs.util.source_selection_all(N)
+d, selection, secondary_source = sfs.mono.drivingfunction.nfchoa_25d_point(
+    omega, array.x, 1, xs)
 compute_and_plot_soundfield('circular_ps_nfchoa_25d_point')
 
-d = sfs.mono.drivingfunction.nfchoa_25d_plane(omega, x0, 1, npw)
-a = sfs.util.source_selection_all(N)
+d, selection, secondary_source = sfs.mono.drivingfunction.nfchoa_25d_plane(
+    omega, array.x, 1, npw)
 compute_and_plot_soundfield('circular_ps_nfchoa_25d_plane')

doc/examples/mirror-image-source-model.ipynb

@@ -93,7 +93,7 @@
    "outputs": [],
    "source": [
     "grid = sfs.util.xyz_grid([0, L[0]], [0, L[1]], 1.5, spacing=0.02)\n",
-    "P = sfs.mono.source.point_image_sources(omega, x0, [1, 0, 0], grid, L,\n",
+    "P = sfs.mono.source.point_image_sources(omega, x0, grid, L,\n",
     "                                        max_order, coeffs=coeffs)"
    ]
   },

doc/examples/modal-room-acoustics.ipynb

@@ -36,7 +36,6 @@
     "x0 = 1, 3, 1.80  # source position\n",
     "L = 6, 6, 3  # dimensions of room\n",
     "deltan = 0.01  # absorption factor of walls\n",
-    "n0 = 1, 0, 0  # normal vector of source (only for compatibility)\n",
     "N = 20  # maximum order of modes"
    ]
   },
@@ -88,7 +87,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "p = sfs.mono.source.point_modal(omega, x0, n0, grid, L, N=N, deltan=deltan)"
+    "p = sfs.mono.source.point_modal(omega, x0, grid, L, N=N, deltan=deltan)"
    ]
   },
   {
@@ -136,7 +135,7 @@
     "\n",
     "receiver = 1, 1, 1.8\n",
     "\n",
-    "p = [sfs.mono.source.point_modal(om, x0, n0, receiver, L, N=N, deltan=deltan)\n",
+    "p = [sfs.mono.source.point_modal(om, x0, receiver, L, N=N, deltan=deltan)\n",
     "     for om in omega]\n",
     "     \n",
     "plt.plot(f, sfs.util.db(p))\n",

doc/examples/plot_particle_density.py

@@ -31,12 +31,12 @@ def plot_particle_displacement(title):
 
 
 # point source
-v = sfs.mono.source.point_velocity(omega, xs, npw, grid)
+v = sfs.mono.source.point_velocity(omega, xs, grid)
 amplitude = 1.5e6
 plot_particle_displacement('particle_displacement_point_source')
 
 # line source
-v = sfs.mono.source.line_velocity(omega, xs, npw, grid)
+v = sfs.mono.source.line_velocity(omega, xs, grid)
 amplitude = 1.3e6
 plot_particle_displacement('particle_displacement_line_source')
 

doc/examples/sound_field_synthesis.py

@@ -25,63 +25,57 @@
 # angular frequency
 omega = 2 * np.pi * f
 # normal vector of plane wave
-npw = sfs.util.direction_vector(np.radians(pw_angle), np.radians(90))
+npw = sfs.util.direction_vector(np.radians(pw_angle))
 
 
 # === get secondary source positions ===
-#x0, n0, a0 = sfs.array.linear(N, dx, center=[-1, 0, 0])
-#x0, n0, a0 = sfs.array.linear_random(N, 0.2*dx, 5*dx)
-#x0, n0, a0 = sfs.array.rectangular(N, dx, orientation=sfs.util.direction_vector(0*np.pi/4, np.pi/2))
-#x0, n0, a0 = sfs.array.circular(N, R)
-x0, n0, a0 = sfs.array.load('../../data/arrays/university_rostock.csv')
+#array = sfs.array.linear(N, dx, center=[-1, 0, 0])
+#array = sfs.array.linear_random(N, 0.2*dx, 5*dx)
+array = sfs.array.rectangular(N, dx, orientation=sfs.util.direction_vector(0*np.pi/4))
+#array = sfs.array.circular(N, R)
+#array = sfs.array.load('../../data/arrays/university_rostock.csv')
 
-#x0, n0, a0 = sfs.array.planar(N, dx, orientation=sfs.util.direction_vector(np.radians(0),np.radians(180)))
-#x0, n0, a0 = sfs.array.cube(N, dx, orientation=sfs.util.direction_vector(0, np.pi/2))
+#array = sfs.array.planar(N, dx, orientation=sfs.util.direction_vector(np.radians(0), np.radians(180)))
+#array = sfs.array.cube(N, dx, orientation=sfs.util.direction_vector(0, np.pi/2))
 
-#x0, n0, a0 = sfs.array.sphere_load('/Users/spors/Documents/src/SFS/data/spherical_grids/equally_spaced_points/006561points.mat', 1, center=[.5,0,0])
+#array = sfs.array.sphere_load('/Users/spors/Documents/src/SFS/data/spherical_grids/equally_spaced_points/006561points.mat', 1, center=[.5,0,0])
 
 
-# === compute driving function ===
-#d = sfs.mono.drivingfunction.delay_3d_plane(omega, x0, n0, npw)
+# === compute driving function and determine active secondary sources ===
+#d, selection, secondary_source = sfs.mono.drivingfunction.delay_3d_plane(omega, array.x, array.n, npw)
 
-#d = sfs.mono.drivingfunction.wfs_2d_line(omega, x0, n0, xs)
+#d, selection, secondary_source = sfs.mono.drivingfunction.wfs_2d_line(omega, array.x, array.n, xs)
 
-#d = sfs.mono.drivingfunction.wfs_2d_plane(omega, x0, n0, npw)
-d = sfs.mono.drivingfunction.wfs_25d_plane(omega, x0, n0, npw, xref)
-#d = sfs.mono.drivingfunction.wfs_3d_plane(omega, x0, n0, npw)
+#d, selection, secondary_source = sfs.mono.drivingfunction.wfs_2d_plane(omega, array.x, array.n, npw)
+d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_plane(omega, array.x, array.n, npw, xref)
+#d, selection, secondary_source = sfs.mono.drivingfunction.wfs_3d_plane(omega, array.x, array.n, npw)
 
-#d = sfs.mono.drivingfunction.wfs_2d_point(omega, x0, n0, xs)
-#d = sfs.mono.drivingfunction.wfs_25d_point(omega, x0, n0, xs)
-#d = sfs.mono.drivingfunction.wfs_3d_point(omega, x0, n0, xs)
+#d, selection, secondary_source = sfs.mono.drivingfunction.wfs_2d_point(omega, array.x, array.n, xs)
+#d, selection, secondary_source = sfs.mono.drivingfunction.wfs_25d_point(omega, array.x, array.n, xs)
+#d, selection, secondary_source = sfs.mono.drivingfunction.wfs_3d_point(omega, array.x, array.n, xs)
 
-#d = sfs.mono.drivingfunction.nfchoa_2d_plane(omega, x0, R, npw)
+#d, selection, secondary_source = sfs.mono.drivingfunction.nfchoa_2d_plane(omega, array.x, R, npw)
 
-#d = sfs.mono.drivingfunction.nfchoa_25d_point(omega, x0, R, xs)
-#d = sfs.mono.drivingfunction.nfchoa_25d_plane(omega, x0, R, npw)
-
-# === determine active secondary sources ===
-a = sfs.util.source_selection_plane(n0, npw)
-#a = sfs.util.source_selection_point(n0, x0, xs)
-#a = sfs.util.source_selection_all(len(x0))
+#d, selection, secondary_source = sfs.mono.drivingfunction.nfchoa_25d_point(omega, array.x, R, xs)
+#d, selection, secondary_source = sfs.mono.drivingfunction.nfchoa_25d_plane(omega, array.x, R, npw)
 
 
 # === compute tapering window ===
-#twin = sfs.tapering.none(a)
-#twin = sfs.tapering.kaiser(a, 8.6)
-twin = sfs.tapering.tukey(a,.3)
+#twin = sfs.tapering.none(selection)
+#twin = sfs.tapering.kaiser(selection, 8.6)
+twin = sfs.tapering.tukey(selection, 0.3)
 
 # === compute synthesized sound field ===
-p = sfs.mono.synthesized.generic(omega, x0, n0, d * twin * a0 , grid,
-                                 source=sfs.mono.source.point)
+p = sfs.mono.synthesize(d, twin, array, secondary_source, grid=grid)
 
 
 # === plot synthesized sound field ===
 plt.figure(figsize=(10, 10))
 sfs.plot.soundfield(p, grid, [0, 0, 0])
-sfs.plot.loudspeaker_2d(x0, n0, twin)
+sfs.plot.loudspeaker_2d(array.x, array.n, twin)
 plt.grid()
 plt.savefig('soundfield.png')
 
 
-#sfs.plot.loudspeaker_3d(x0, n0, twin)
+#sfs.plot.loudspeaker_3d(array.x, array.n, twin)
 #plt.savefig('loudspeakers.png')

doc/examples/soundfigures.py

@@ -24,16 +24,16 @@
 grid = sfs.util.xyz_grid([-3, 3], [-3, 3], 0, spacing=0.02)
 
 # get secondary source positions
-x0, n0, a0 = sfs.array.cube(N, dx)
+array = sfs.array.cube(N, dx)
 
 # driving function for sound figure
 figure = np.array(Image.open('figures/tree.png'))  # read image from file
 figure = np.rot90(figure)  # turn 0deg to the top
-d = sfs.mono.soundfigure.wfs_3d_pw(omega, x0, n0, figure, npw=npw)
+d, selection, secondary_source = sfs.mono.soundfigure.wfs_3d_pw(
+    omega, array.x, array.n, figure, npw=npw)
 
 # compute synthesized sound field
-p = sfs.mono.synthesized.generic(omega, x0, n0, d * a0, grid,
-                                 source=sfs.mono.source.point)
+p = sfs.mono.synthesize(d, selection, array, secondary_source, grid=grid)
 
 # plot and save synthesized sound field
 plt.figure(figsize=(10, 10))