Skip to content

visualdynamics.core.octave

octave

Proportional-band spectra: a PSD resampled onto octave bands.

A narrowband PSD answers "how much power per hertz, at this hertz". A proportional-band one answers "how much power in this band", where the bands get wider as the frequency rises — which is how a structure's response is usually specified, and how an ear hears.

The bands are the base-ten system of ANSI S1.11 / IEC 61260, which is what sdynpy uses and therefore what a result has to agree with:

  • one octave is a factor of ten to the three tenths, not a factor of two. For nth-octave the band ratio is 10 ** (3 / (10 n)) — at a sixth of an octave that is 1.122018, where a base-two reading would give 1.122462. They differ in the fourth digit, which is enough to put every band edge in a slightly different place.
  • the grid is absolute. It does not depend on the frequency range asked for, or on any specification: the range only chooses which bands of the one fixed grid are returned. 1000 Hz anchors it.
  • and which of edges or centres lands on 10 ** (3 k / (10 n)) depends on whether n is odd or even. For odd fractions — whole octaves, thirds — the band centres sit on the grid; for even ones — sixths, twelfths — the edges do, and the centres fall half a step between. Getting this backwards puts every band half a step out, which is a real disagreement and not a rounding one.
  • bands tile — each one's upper edge is the next one's lower — and a band's centre is the geometric mean of its edges.

Converting a spectrum is an integration, not a resampling. The value in a band is the mean-square content of that band divided by its width, so the area under the spectrum is unchanged and the RMS it carries is the RMS it carried before. Reading the narrowband curve at each band centre instead would throw away everything between the centres, and would not conserve anything.

Functions:

Name Description
ratio

How much wider each band is than the one below it.

edges

The edges of every band overlapping low to high.

bands

(centres, widths, edges) for every band overlapping the range.

bin_bounds

(left, right) of the bin each line stands for.

resample

values integrated onto the bands bounds describes.

Functions:

ratio

ratio(per_octave: int = PER_OCTAVE) -> float

How much wider each band is than the one below it.

Source code in src/visualdynamics/core/octave.py
def ratio(per_octave: int = PER_OCTAVE) -> float:
    """How much wider each band is than the one below it."""
    per_octave = int(per_octave)
    if per_octave < 1:
        raise ValueError('there is at least one band to an octave')
    return 10.0 ** (DECADE_FRACTION / per_octave)

edges

edges(
    low: float, high: float, per_octave: int = PER_OCTAVE
) -> ndarray

The edges of every band overlapping low to high.

One more edge than there are bands, since they tile. Snapped to the fixed grid rather than started at low: two spectra covering different ranges land on the same bands, which is the whole point of a standard grid and the reason two runs can be compared at all.

Source code in src/visualdynamics/core/octave.py
def edges(low: float, high: float,
          per_octave: int = PER_OCTAVE) -> np.ndarray:
    """The edges of every band overlapping `low` to `high`.

    One more edge than there are bands, since they tile. Snapped to the
    fixed grid rather than started at `low`: two spectra covering
    different ranges land on the same bands, which is the whole point
    of a standard grid and the reason two runs can be compared at all.
    """
    low, high = float(low), float(high)
    if not (low > 0.0 and high > low):
        raise ValueError(f'band range {low} to {high} is not a range')
    per_octave = int(per_octave)
    step = DECADE_FRACTION / per_octave
    # odd fractions put the centres on the grid, so their edges fall
    # half a step off it; even fractions put the edges on it
    offset = 0.5 if per_octave % 2 else 0.0
    # snapped before rounding: an edge that is already a grid point
    # comes back from ten-to-the-power and a logarithm a float's
    # breadth adrift, and flooring that adds a band on every re-banding
    below = np.log10(low) / step - offset
    above = np.log10(high) / step - offset
    first = int(np.floor(below + SNAP))
    last = int(np.ceil(above - SNAP))
    return 10.0 ** (step * (np.arange(first, last + 1) + offset))

bands

bands(
    low: float, high: float, per_octave: int = PER_OCTAVE
) -> tuple[ndarray, ndarray, ndarray]

(centres, widths, edges) for every band overlapping the range.

The centre is the geometric mean of the band's own edges — the arithmetic mean would sit above it, and on a log axis a band would then be drawn off-centre from the number naming it.

Source code in src/visualdynamics/core/octave.py
def bands(low: float, high: float, per_octave: int = PER_OCTAVE
          ) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    """(centres, widths, edges) for every band overlapping the range.

    The centre is the geometric mean of the band's own edges — the
    arithmetic mean would sit above it, and on a log axis a band would
    then be drawn off-centre from the number naming it.
    """
    bounds = edges(low, high, per_octave)
    lower, upper = bounds[:-1], bounds[1:]
    return np.sqrt(lower * upper), upper - lower, bounds

bin_bounds

bin_bounds(
    centres: ArrayLike, widths: ArrayLike | None = None
) -> tuple[ndarray, ndarray]

(left, right) of the bin each line stands for.

With no widths the bins are the midpoints between neighbours, which is exact for evenly spaced FFT lines. With them the edges follow from c = sqrt(l u) and w = u - l:

u = (w + sqrt(w^2 + 4 c^2)) / 2,   l = u - w

the positive root, which is where a proportional band's edges actually are. Taking c +/- w/2 instead would be assuming the centre is the arithmetic mean of the edges, and it is not.

Source code in src/visualdynamics/core/octave.py
def bin_bounds(centres: ArrayLike, widths: ArrayLike | None = None
               ) -> tuple[np.ndarray, np.ndarray]:
    """(left, right) of the bin each line stands for.

    With no widths the bins are the midpoints between neighbours, which
    is exact for evenly spaced FFT lines. With them the edges follow
    from `c = sqrt(l u)` and `w = u - l`:

        u = (w + sqrt(w^2 + 4 c^2)) / 2,   l = u - w

    the positive root, which is where a proportional band's edges
    actually are. Taking `c +/- w/2` instead would be assuming the
    centre is the arithmetic mean of the edges, and it is not.
    """
    centres = np.asarray(centres, dtype=float)
    if widths is not None:
        widths = np.asarray(widths, dtype=float)
        upper = (widths + np.sqrt(widths * widths
                                  + 4.0 * centres * centres)) / 2.0
        return upper - widths, upper
    spacing = np.gradient(centres)
    return centres - spacing / 2.0, centres + spacing / 2.0

resample

resample(
    frequencies: ArrayLike,
    values: ArrayLike,
    bounds: ArrayLike,
    widths: ArrayLike | None = None,
    source: ArrayLike | None = None,
) -> ndarray

values integrated onto the bands bounds describes.

Each line stands for its own bin, flat across it, which is what a discrete density is. A band's content is the part of every bin that falls inside it, and the band's value is that content over the band's width — so the area under the spectrum survives and a band straddling the end of the data is credited only for the part of it that has data underneath.

Complex values pass through as complex: the cross terms of a CPSD average over a band the same way the diagonal does.

Source code in src/visualdynamics/core/octave.py
def resample(frequencies: ArrayLike, values: ArrayLike, bounds: ArrayLike,
             widths: ArrayLike | None = None,
             source: ArrayLike | None = None) -> np.ndarray:
    """`values` integrated onto the bands `bounds` describes.

    Each line stands for its own bin, flat across it, which is what a
    discrete density is. A band's content is the part of every bin that
    falls inside it, and the band's value is that content over the
    band's width — so the area under the spectrum survives and a band
    straddling the end of the data is credited only for the part of it
    that has data underneath.

    Complex values pass through as complex: the cross terms of a CPSD
    average over a band the same way the diagonal does.
    """
    frequencies = np.asarray(frequencies, dtype=float)
    values = np.atleast_2d(np.asarray(values))
    if frequencies.size < 2:
        raise ValueError('a spectrum needs at least two lines to be banded')

    # the bins the source actually has. Inferring them from the centres
    # is exact for evenly spaced lines and wrong for a spectrum that is
    # already banded — the guess runs a fraction of a percent wide
    # through the middle and six percent wide at the first band, which
    # smears content into its neighbours every time one is re-banded.
    left, right = source if source is not None else bin_bounds(frequencies)
    lower, upper = bounds[:-1], bounds[1:]

    # (bands, lines): how much of each line's bin lies in each band
    overlap = np.clip(
        np.minimum(right[None, :], upper[:, None])
        - np.maximum(left[None, :], lower[:, None]), 0.0, None)
    held = np.nan_to_num(values, nan=0.0)
    content = held @ overlap.T                      # (records, bands)
    widths = (upper - lower) if widths is None else np.asarray(widths)
    with np.errstate(divide='ignore', invalid='ignore'):
        banded = np.where(widths > 0.0, content / widths, 0.0)
    # a band with no data under it says nothing rather than nothing-much
    covered = overlap.sum(axis=1)
    return np.where(covered > 0.0, banded, np.nan)