Utilities

Various utility functions.

sfs.util.rotation_matrix(n1, n2)[source]

Compute rotation matrix for rotation from n1 to n2.

Parameters:n1, n2 ((3,) array_like) – Two vectors. They don’t have to be normalized.
Returns:(3, 3) numpy.ndarray – Rotation matrix.
sfs.util.wavenumber(omega, c=None)[source]

Compute the wavenumber for a given radial frequency.

sfs.util.direction_vector(alpha, beta=1.5707963267948966)[source]

Compute normal vector from azimuth, colatitude.

sfs.util.sph2cart(alpha, beta, r)[source]

Spherical to cartesian coordinate transform.

\[\begin{split}x = r \cos \alpha \sin \beta \\ y = r \sin \alpha \sin \beta \\ z = r \cos \beta\end{split}\]

with \(\alpha \in [0, 2\pi), \beta \in [0, \pi], r \geq 0\)

Parameters:
  • alpha (float or array_like) – Azimuth angle in radiants
  • beta (float or array_like) – Elevation angle in radiants (with 0 denoting North pole)
  • r (float or array_like) – Radius
Returns:

  • x (float or array_like) – x-component of Cartesian coordinates
  • y (float or array_like) – y-component of Cartesian coordinates
  • z (float or array_like) – z-component of Cartesian coordinates

sfs.util.cart2sph(x, y, z)[source]

Cartesian to spherical coordinate transform.

\[\begin{split}\alpha = \arctan \left( \frac{y}{x} \right) \\ \beta = \arccos \left( \frac{z}{r} \right) \\ r = \sqrt{x^2 + y^2 + z^2}\end{split}\]

with \(\alpha \in [0, 2\pi), \beta \in [0, \pi], r \geq 0\)

Parameters:
  • x (float or array_like) – x-component of Cartesian coordinates
  • y (float or array_like) – y-component of Cartesian coordinates
  • z (float or array_like) – z-component of Cartesian coordinates
Returns:

  • alpha (float or array_like) – Azimuth angle in radiants
  • beta (float or array_like) – Elevation angle in radiants (with 0 denoting North pole)
  • r (float or array_like) – Radius

sfs.util.asarray_1d(a, **kwargs)[source]

Squeeze the input and check if the result is one-dimensional.

Returns a converted to a numpy.ndarray and stripped of all singleton dimensions. Scalars are “upgraded” to 1D arrays. The result must have exactly one dimension. If not, an error is raised.

sfs.util.asarray_of_rows(a, **kwargs)[source]

Convert to 2D array, turn column vector into row vector.

Returns a converted to a numpy.ndarray and stripped of all singleton dimensions. If the result has exactly one dimension, it is re-shaped into a 2D row vector.

sfs.util.as_xyz_components(components, **kwargs)[source]

Convert components to XyzComponents of numpy.ndarrays.

The components are first converted to NumPy arrays (using numpy.asarray()) which are then assembled into an XyzComponents object.

Parameters:
  • components (triple or pair of array_like) – The values to be used as X, Y and Z arrays. Z is optional.
  • **kwargs – All further arguments are forwarded to numpy.asarray(), which is applied to the elements of components.
sfs.util.as_delayed_signal(arg, **kwargs)[source]

Make sure that the given argument can be used as a signal.

Parameters:
  • arg (sequence of 1 array_like followed by 1 or 2 scalars) – The first element is converted to a NumPy array, the second element is used as the sampling rate (in Hertz) and the optional third element is used as the starting time of the signal (in seconds). Default starting time is 0.
  • **kwargs – All keyword arguments are forwarded to numpy.asarray().
Returns:

DelayedSignal – A named tuple consisting of a numpy.ndarray containing the audio data, followed by the sampling rate (in Hertz) and the starting time (in seconds) of the signal.

Examples

Typically, this is used together with tuple unpacking to assign the audio data, the sampling rate and the starting time to separate variables:

>>> import sfs
>>> sig = [1], 44100
>>> data, fs, signal_offset = sfs.util.as_delayed_signal(sig)
>>> data
array([1])
>>> fs
44100
>>> signal_offset
0
sfs.util.strict_arange(start, stop, step=1, endpoint=False, dtype=None, **kwargs)[source]

Like numpy.arange(), but compensating numeric errors.

Unlike numpy.arange(), but similar to numpy.linspace(), providing endpoint=True includes both endpoints.

Parameters:
  • start, stop, step, dtype – See numpy.arange().

  • endpoint – See numpy.linspace().

    Note

    With endpoint=True, the difference between start and end value must be an integer multiple of the corresponding spacing value!

  • **kwargs – All further arguments are forwarded to numpy.isclose().

Returns:

numpy.ndarray – Array of evenly spaced values. See numpy.arange().

sfs.util.xyz_grid(x, y, z, spacing, endpoint=True, **kwargs)[source]

Create a grid with given range and spacing.

Parameters:
  • x, y, z (float or pair of float) – Inclusive range of the respective coordinate or a single value if only a slice along this dimension is needed.
  • spacing (float or triple of float) – Grid spacing. If a single value is specified, it is used for all dimensions, if multiple values are given, one value is used per dimension. If a dimension (x, y or z) has only a single value, the corresponding spacing is ignored.
  • endpoint (bool, optional) – If True (the default), the endpoint of each range is included in the grid. Use False to get a result similar to numpy.arange(). See strict_arange().
  • **kwargs – All further arguments are forwarded to strict_arange().
Returns:

XyzComponents – A grid that can be used for sound field calculations.

sfs.util.normalize(p, grid, xnorm)[source]

Normalize sound field wrt position xnorm.

sfs.util.probe(p, grid, x)[source]

Determine the value at position x in the sound field p.

sfs.util.broadcast_zip(*args)[source]

Broadcast arguments to the same shape and then use zip().

sfs.util.normalize_vector(x)[source]

Normalize a 1D vector.

sfs.util.displacement(v, omega)[source]

Particle displacement

\[d(x, t) = \int_0^t v(x, t) dt\]
sfs.util.db(x, power=False)[source]

Convert x to decibel.

Parameters:
  • x (array_like) – Input data. Values of 0 lead to negative infinity.
  • power (bool, optional) – If power=False (the default), x is squared before conversion.
class sfs.util.XyzComponents(components)[source]

Triple (or pair) of components: x, y, and optionally z.

Instances of this class can be used to store coordinate grids (either regular grids like in xyz_grid() or arbitrary point clouds) or vector fields (e.g. particle velocity).

This class is a subclass of numpy.ndarray. It is one-dimensional (like a plain list) and has a length of 3 (or 2, if no z-component is available). It uses dtype=object in order to be able to store other numpy.ndarrays of arbitrary shapes but also scalars, if needed. Because it is a NumPy array subclass, it can be used in operations with scalars and “normal” NumPy arrays, as long as they have a compatible shape. Like any NumPy array, instances of this class are iterable and can be used, e.g., in for-loops and tuple unpacking. If slicing or broadcasting leads to an incompatible shape, a plain numpy.ndarray with dtype=object is returned.

To make sure the components are NumPy arrays themselves, use as_xyz_components().

Parameters:components ((3,) or (2,) array_like) – The values to be used as X, Y and Z data. Z is optional.
x

x-component.

y

y-component.

z

z-component (optional).

apply(func, *args, **kwargs)[source]

Apply a function to each component.

The function func will be called once for each component, passing the current component as first argument. All further arguments are passed after that. The results are returned as a new XyzComponents object.

class sfs.util.DelayedSignal

A tuple of audio data, sampling rate and start time.

This class (a collections.namedtuple) is not meant to be instantiated by users.

To pass a signal to a function, just use a simple tuple or list containing the audio data and the sampling rate (in Hertz), with an optional starting time (in seconds) as a third item. If you want to ensure that a given variable contains a valid signal, use sfs.util.as_delayed_signal().

data

Alias for field number 0

samplerate

Alias for field number 1

time

Alias for field number 2

sfs.util.image_sources_for_box(x, L, N, prune=True)[source]

Image source method for a cuboid room.

The classical method by Allen and Berkley [AB79].

Parameters:
  • x ((D,) array_like) – Original source location within box. Values between 0 and corresponding side length.
  • L ((D,) array_like) – side lengths of room.
  • N (int) – Maximum number of reflections per image source, see below.
  • prune (bool, optional) – selection of image sources:
    • If True (default): Returns all images reflected up to N times. This is the usual interpretation of N as “maximum order”.
    • If False: Returns reflected up to N times between individual wall pairs, a total number of \(M := (2N+1)^D\). This larger set is useful e.g. to select image sources based on distance to listener, as suggested by [Bor84].
Returns:

  • xs ((M, D) array_like) – original & image source locations.
  • wall_count ((M, 2D) array_like) – number of reflections at individual walls for each source.

sfs.util.spherical_hn2(n, z)[source]

Spherical Hankel function of 2nd kind.

Defined as http://dlmf.nist.gov/10.47.E6,

\[\hankel{2}{n}{z} = \sqrt{\frac{\pi}{2z}} \Hankel{2}{n + \frac{1}{2}}{z},\]

where \(\Hankel{2}{n}{\cdot}\) is the Hankel function of the second kind and n-th order, and \(z\) its complex argument.

Parameters:
  • n (array_like) – Order of the spherical Hankel function (n >= 0).
  • z (array_like) – Argument of the spherical Hankel function.