Skip to content

Animation of time-domain WFS #193

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open

Animation of time-domain WFS #193

Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
115 changes: 115 additions & 0 deletions doc/examples/wfs-animation.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,115 @@
{
"cells": [
{
"cell_type": "code",
"execution_count": null,
"id": "77fc3672",
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"import sfs\n",
"from matplotlib import animation\n",
"from IPython.display import HTML\n",
"from functools import partial\n",
"from scipy.signal import unit_impulse\n",
"\n",
"\n",
"# Point source\n",
"xs = 0, 2, 0\n",
"rs = np.linalg.norm(xs) # distance from origin\n",
"ts = rs / sfs.default.c # time-of-arrival at origin\n",
"\n",
"# Impulsive excitation\n",
"fs = 8000 # Adjust this to change the shape (width) of the impulse\n",
"signal = unit_impulse(512), fs # Band-limited pulse (e.g. sinc) can be used instead\n",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is unit_impulse really necessary here? What about using this?

signal = [1], fs

"\n",
"# Circular loudspeaker array\n",
"N = 32 # number of loudspeakers\n",
"R = 1.5 # radius\n",
"array = sfs.array.circular(N, R)\n",
"\n",
"grid = sfs.util.xyz_grid([-2, 2], [-2, 2], 0, spacing=0.02)\n",
"\n",
"delays, weights, selection, secondary_source = \\\n",
" sfs.td.wfs.point_25d(array.x, array.n, xs)\n",
"d = sfs.td.wfs.driving_signals(delays, weights, signal)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "e3560c30",
"metadata": {
"scrolled": false
},
"outputs": [],
"source": [
"# Animation\n",
"def plot(d, selection, secondary_source, t=0, ax=None, **kw):\n",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In general, I'm a fan of reusable plot functions, especially if users can copy it and use it in their own code. However, this doesn't work here, because grid and array are not passed as function arguments.

I'm also not sure if it helps the understanding here, in the context of the animation, what do you think?

I think you should either add more arguments or completely dissolve the function. I could imagine a few notebook cells with easily digestible steps like this:

Create loudspeaker driving signals.

delays, weights, selection, secondary_source = sfs.td.wfs.point_25d(
    array.x, array.n, xs)
d = sfs.td.wfs.driving_signals(delays, weights, signal)

Create one frame of the initial sound field for the initial plot.
The observation_time doesn't matter,
because this frame will not be used in the final animation.

p = sfs.td.synthesize(
    d, selection, array, secondary_source, grid=grid, observation_time=0)

Create the initial plot (but don't show it yet).

fig, ax = plt.subplots(figsize=(5, 5))
im = sfs.plot2d.amplitude(p, grid, ax=ax, vmin=-0.01, vmax=0.01)
sfs.plot2d.loudspeakers(array.x, array.n, selection * array.a, size=0.15)
plt.close()

Define a function which will update im for each frame of the animation.

def update_frame_pressure(t):
    p = sfs.td.synthesize(
        d, selection, array, secondary_source, grid=grid, observation_time=t)
    im.set_array(p)
    return [im]

Put everything together to create an animation.

ani = animation.FuncAnimation(
    fig, update_frame_pressure, frames=time_stamps,
    interval=interval, blit=True)

" p = sfs.td.synthesize(d, selection, array, secondary_source, grid=grid,\n",
" observation_time=t)\n",
" im = sfs.plot2d.amplitude(p, grid, ax=ax, **kw)\n",
" sfs.plot2d.loudspeakers(array.x, array.n, selection * array.a, size=0.15)\n",
" return im\n",
"\n",
"def update_frame_pressure(i, time_stamps):\n",
" t_i = time_stamps[i]\n",
" p = sfs.td.synthesize(d, selection, array, secondary_source, grid=grid,\n",
" observation_time=t_i)\n",
" im.set_array(p)\n",
" return [im]\n",
"\n",
"\n",
"time_stamps = np.linspace(0.5/343, 5/343, 100) # Time sampling is different from fs defined above\n",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The number 343 is suspicious here. If you want to use the speed of sound, you should use sfs.default.c, if not, it would probably be better to use different values.

"frames = 100\n",
"interval = 150\n",
"\n",
"\n",
"fig, ax = plt.subplots(figsize=(5, 5))\n",
"p = sfs.td.synthesize(d, selection, array, secondary_source, grid=grid,\n",
" observation_time=0)\n",
Comment on lines +71 to +72
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think that this does anything, right?

Suggested change
"p = sfs.td.synthesize(d, selection, array, secondary_source, grid=grid,\n",
" observation_time=0)\n",

It just creates p which is then never used?

"im = plot(d, selection, secondary_source, t=ts, ax=ax, vmin=-0.01, vmax=0.01)\n",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This frame is not used in the final animation, so the value of t is not really relevant, and defining ts is not necessary.

"\n",
"ani = animation.FuncAnimation(\n",
" fig, partial(update_frame_pressure, time_stamps=time_stamps),\n",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is a bit strange that you use partial function application here, when you have full control over the function. You could simply move the time_stamps into the function.

But even better, you could use the time_stamps as frames:

def update_frame_pressure(t):
    p = sfs.td.synthesize(
        d, selection, array, secondary_source, grid=grid, observation_time=t)
    im.set_array(p)
    return [im]

ani = animation.FuncAnimation(
    fig, update_frame_pressure, frames=time_stamps,
    interval=interval, blit=True)

" frames=frames, interval=interval, blit=True)\n",
"plt.close()\n",
"HTML(ani.to_jshtml())"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "9f77a060",
"metadata": {},
"outputs": [],
"source": [
"# Save as gif file - This might take a few minutes.\n",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Whenever I write that something takes a long time, I also try to show how long it actually takes. You could use the %%time magic here:

%%time
ani.save("wfs-25d-td.gif", writer='imagemagick', fps=10, dpi=200)

"ani.save(\"wfs-25d-td.gif\", writer='imagemagick',fps=10, dpi=200)"
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

After creating the GIF file, you could also include it in a Markdown cell like this:

![](wfs-25d-td.gif)

BTW, why are you specifying fps here and not above when using to_jshtml()?

]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"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.9.5"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Loading