Room Acoustics: Image sources#
Contact: m.berzborn@tue.nl
This assignment is part of the course “Architectural Acoustics” at the Eindhoven University of Technology. The goal of the assignment is to predict the required order of image sources for a given room. Sabine’s equation is used to predict the reverberation time of the room, which is then used to calculate the required order of image sources.
Duration: 30 Minutes
Requirements: Basic knowledge of geometrical acoustics and the image source model. Familiarity with Python and the used libraries is helpful but not required.
References
Kuttruff, Room acoustics, Taylor & Francis.
Vorländer, Auralization, Springer-Verlag Berlin Heidelberg, 2008.
Dependencies
If you are running this notebook locally or in Google Colab, make sure to install the required dependencies. You can do this by executing the following command in an arbitrary code cell:
!pip install numpy matplotlib pyfar ipympl pyroomacoustics
[ ]:
import pyfar as pf
import numpy as np
import pyroomacoustics as pra
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle
%matplotlib ipympl
The following cell contains code required to create the visualization of the image sources up to a specified order. Please run the cell to make sure the visualization works later in the notebook. You do not need to understand the code in this cell to complete the assignment.
[ ]:
def mirror_plot_content_n_times(ax, room_size, points, n_orders, draw_room=True):
"""Mirror 2D points across a rectangular room boundary up to order n_orders.
Uses a diamond-shaped neighborhood based on Manhattan distance, i.e.
abs(i) + abs(j) <= n_orders.
"""
l_x, l_y = room_size
for i in range(-n_orders, n_orders + 1):
for j in range(-n_orders, n_orders + 1):
if abs(i) + abs(j) > n_orders:
continue
room_origin = (i * l_x, j * l_y)
if draw_room:
ax.add_patch(Rectangle(room_origin, l_x, l_y, fill=False, alpha=1))
for point in points:
is_original = i == 0 and j == 0
# Skip mirrored positions if point should only appear in original room
if point.get("only_original", False) and not is_original:
continue
x, y = point["xy"]
x_mirror = i * l_x + (x if i % 2 == 0 else l_x - x)
y_mirror = j * l_y + (y if j % 2 == 0 else l_y - y)
marker = point.get("marker_original", "o") if is_original else point.get("marker_image", "o")
ax.scatter(
x_mirror,
y_mirror,
label=point.get("label") if is_original else None,
marker=marker,
s=90 if is_original else 40,
facecolors=point.get("color", "black") if is_original else "none",
edgecolors=point.get("color", "black"),
linewidths=1.5,
)
def draw_reflection_paths(ax, source_xy, receiver_xy, room_size, color="0.25"):
"""Draw first- and second-order reflection paths for opposite walls x=L_x and y=L_y."""
source_xy = np.asarray(source_xy, dtype=float)
receiver_xy = np.asarray(receiver_xy, dtype=float)
l_x, l_y = room_size
label_used = {"construction": False, "first": False, "second": False}
def _line_intersection_with_vertical(p0, p1, x_const):
direction = p1 - p0
if np.isclose(direction[0], 0):
return None
t = (x_const - p0[0]) / direction[0]
return p0 + t * direction
def _line_intersection_with_horizontal(p0, p1, y_const):
direction = p1 - p0
if np.isclose(direction[1], 0):
return None
t = (y_const - p0[1]) / direction[1]
return p0 + t * direction
def _within(value, lower, upper, eps=1e-9):
return (lower - eps) <= value <= (upper + eps)
def _draw_first_order(source_image, wall_type, wall_value):
if wall_type == "x":
reflection_point = _line_intersection_with_vertical(source_image, receiver_xy, wall_value)
is_valid = reflection_point is not None and _within(reflection_point[1], 0, l_y)
else:
reflection_point = _line_intersection_with_horizontal(source_image, receiver_xy, wall_value)
is_valid = reflection_point is not None and _within(reflection_point[0], 0, l_x)
if not is_valid:
return
ax.plot(
[receiver_xy[0], source_image[0]],
[receiver_xy[1], source_image[1]],
linestyle=":",
color=color,
linewidth=1.2,
label="Construction lines" if not label_used["construction"] else None,
)
label_used["construction"] = True
ax.plot(
[source_xy[0], reflection_point[0]],
[source_xy[1], reflection_point[1]],
linestyle="--",
color=color,
linewidth=1.8,
label="1st-order reflection path" if not label_used["first"] else None,
)
label_used["first"] = True
ax.plot(
[reflection_point[0], receiver_xy[0]],
[reflection_point[1], receiver_xy[1]],
linestyle="--",
color=color,
linewidth=1.8,
)
source_image_x = np.array([2 * l_x - source_xy[0], source_xy[1]])
source_image_y = np.array([source_xy[0], 2 * l_y - source_xy[1]])
_draw_first_order(source_image_x, wall_type="x", wall_value=l_x)
_draw_first_order(source_image_y, wall_type="y", wall_value=l_y)
source_image_xy = np.array([2 * l_x - source_xy[0], 2 * l_y - source_xy[1]])
p_y = _line_intersection_with_horizontal(source_image_xy, receiver_xy, l_y)
if p_y is not None and _within(p_y[0], 0, l_x):
p_x = _line_intersection_with_vertical(source_image_x, p_y, l_x)
if p_x is not None and _within(p_x[1], 0, l_y):
ax.plot(
[receiver_xy[0], source_image_xy[0]],
[receiver_xy[1], source_image_xy[1]],
linestyle=":",
color="0.45",
linewidth=1.2,
label=None,
)
ax.plot(
[source_xy[0], p_x[0], p_y[0], receiver_xy[0]],
[source_xy[1], p_x[1], p_y[1], receiver_xy[1]],
linestyle="-.",
color="0.1",
linewidth=2.0,
label="2nd-order reflection path" if not label_used["second"] else None,
)
label_used["second"] = True
Room description#
Assume a rectangular room with the dimensions \([L_x, L_y, L_z] = [6, 3.5, 2.5]\) m. The speed of sound is \(c = 343\) m/s.
[ ]:
L_x = 6
L_y = 3.5
L_z = 2.5
L = [L_x, L_y, L_z] # m
c = 343 # m/s
Image source visualization#
By executing the following cell, you can visualize the image sources up to a specified order. This is chosen as 2 by default. A sound source is located at \([x_s, y_s, z_s] = [1, 1, 1.8]\) m and a microphone is initially located at \([x_m, y_m, z_m] = [3, 2.5, 1.25]\) m.
Task: Move the microphone to a different position and observe how the reflection paths change. Finally, identify a position in the room where the second order reflection is no longer audible. Explain why this is the case. Hint: When the reflection is no longer audible, the reflection path disappears from the figure.
[ ]:
source_position = [1, 1, 1.8]
microphone_position = [3.5, 2.75, 1.25]
plot_image_source_order = 2
fig, ax = plt.subplots(figsize=(10, 6))
mirror_plot_content_n_times(
ax=ax,
room_size=(L_x, L_y),
points=[
{
"xy": source_position[:2],
"color": "C0",
"label": "Source",
},
{
"xy": microphone_position[:2],
"color": "C3",
"label": "Microphone",
"only_original": True,
},
],
n_orders=plot_image_source_order,
)
draw_reflection_paths(ax, source_position[:2], microphone_position[:2], room_size=(L_x, L_y))
ax.set_aspect("equal")
ax.set_xlim(-(plot_image_source_order + 0.25) * L_x, (plot_image_source_order + 1.25) * L_x)
ax.set_ylim(-(plot_image_source_order + 0.25) * L_y, (plot_image_source_order + 1.25) * L_y)
ax.legend(loc="upper left", bbox_to_anchor=(1.02, 1.0), borderaxespad=0)
ax.set_xticks(np.arange((-plot_image_source_order)*L_x, (plot_image_source_order+2)*L_x, L_x))
ax.set_xticklabels([f"{n} $L_x$" for n in np.arange((-plot_image_source_order), (plot_image_source_order+2))])
ax.set_yticks(np.arange((-plot_image_source_order)*L_y, (plot_image_source_order+2)*L_y, L_y))
ax.set_yticklabels([f"{n} $L_y$" for n in np.arange((-plot_image_source_order), (plot_image_source_order+2))])
fig.tight_layout()
plt.show()
Reverberation time prediction#
To ensure that a sufficiently high order of image sources is considered, a reverberation time prediction is required. Sabine’s equation can be used to predict the reverberation time of the room, which is then used to calculate the required order of image sources. Make sure to save the predicted reverberation time in a variable named reverberation_time. Task: Implement Sabine’s equation to predict the reverberation time of the room.
[ ]:
average_absorption_coefficient = 0.15
# YOUR CODE HERE
raise NotImplementedError()
print(f"Predicted reverberation time: {reverberation_time:.2f} s")
Task: Now estimate the required order of image sources using the predicted reverberation time, the speed of sound and the room dimensions. Make sure to store the solution in a variable named maximum_image_source_order. Hint: The required order of image sources is not the same as the total number of image sources. The figure above shows all image sources up to second order if you did not change the respective code cell.
[ ]:
# YOUR CODE HERE
raise NotImplementedError()
Room impulse response simulation#
Finally, you can compute the room impulse response using the image source model. This is done in the following cell which you can run after implementing the previous tasks. No changes are required to run the following cells.
[ ]:
sampling_rate = 8000
room = pra.ShoeBox(
L, fs=sampling_rate, max_order=maximum_image_source_order,
materials=pra.Material(average_absorption_coefficient, 0.),
air_absorption=True,
ray_tracing=False,
)
room.add_source(source_position)
room.add_microphone(microphone_position)
room.compute_rir()
room_impulse_response = pf.dsp.normalize(pf.Signal(room.rir[0][0], sampling_rate))
To visualize the simulated room impulse response, run the following cell. You can also use the visualization to make sure the chosen image source order is sufficient.
[ ]:
plt.figure(figsize=(8, 4))
ax = pf.plot.time(room_impulse_response, dB=True)
ax.set_xlim(0, 3)
ax.set_ylim(-65, 5)
plt.show()
License notice#
This notebook © 2026 by Marco Berzborn is licensed under CC BY 4.0
Watermark#
The following watermark might help others to install specific package versions that might be required to run the notebook. Please give at least the versions of Python, IPython, numpy , and scipy, major third party packagers (e.g., pytorch), and all used pyfar packages.
[ ]:
%load_ext watermark
%watermark -v -m -p numpy,scipy,pyfar,sofar,watermark