visualdynamics.core.sine_distortion

sine_distortion

The sine distortion check: what is on a channel besides the tone.

A sine controller reads its control channel through a narrow band that follows the drive (core.sine_tracking), so the level it holds is the fundamental and nothing else. Whatever else the channel carries is outside the band: harmonics from a rattle, a loose fixture bolt, a joint that claps or an amplifier at its limit, a fixture resonance excited at a multiple of the drive, noise. The controller does not see it, and the article feels it. The check reads it directly. The tones the extraction tracked are rebuilt on the record's own samples and taken away from it, and what is left is everything that is not a tone. With r the RMS of that residual and x_1 the RMS of the rebuilt tones, both over the same span of the record, the distortion is

D = r / x_1,

reported as a fraction (a percentage on screen) and in dB, 20 log10 D.

The tolerance is usually written on a ratio rather than on D: the unfiltered RMS no more than some share above the filtered one. Where the residual is orthogonal to the tones their mean squares add, so the ratio is sqrt(1 + D2) and D = sqrt(ratio2 - 1) (distortion). Five percent is an example of what test documents ask for, and it is the default here (DISTORTION_EXCESS). Five percent over is a distortion of sqrt(1.05**2 - 1) = 0.320, 32 % or -9.9 dB (distortion_limit), which is the line the check draws and the figure a line is flagged above.

The residual is read directly because D taken from two readings, the unfiltered RMS over the tone's, rests on two things that need not hold. The tones must be orthogonal to each other for their powers to add, and two tones sweeping through each other are not orthogonal where they cross. And the two readings must cover the same stretch of the record: a level read at an instant and an RMS read over a span differ wherever the level moves inside the span, as it does through a sharp resonance. The residual is read over the same samples as the tones it is measured against, and two tones that beat or cross are rebuilt together and subtracted together, so neither assumption is made.

The check is not sensitive, and it is not meant to be. A tone with a third harmonic at 10 % of it reads a ratio of 1.005, at 20 % 1.020, and only at 32 % does it reach 1.05: it catches distortion large enough to change the test. It hardly sees clipping at all. A sine clipped flat at four fifths of its peak reads 1.004, because cutting the tops of a sine puts most of the lost power back into the fundamental, so flat tops are found in the time history and its crest factor, not here.

The tones are rebuilt by the estimator the extraction uses (core.sine.vold_kalman): every tone the levels hold that plays over the checked tone's stretch, solved jointly on the record's own samples with the extraction's smoothing, each tone's complex envelope on its own sweep argument. The envelope follows the level as it moves, a resonance included, so the rebuilt tone is the tone the controller's band would pass, sample by sample. Where the extraction corrected the sweep clock (SineLevel.drift_hz not zero) the rebuild fits the same kind of correction to its own envelope's phase and solves again, so a drifted sweep is not rebuilt off its tone. The residual and the tones are then read over cycles cycles of the tone's phase centered on each line's instant (the extraction's smoothing by default), moved inward where a span would run off either end of the tone's stretch.

The record's mean over the tone's stretch is taken out first: a sensor's offset is not motion the article feels. Anything slower than the span that is not constant is counted, as it should be, since it is not the tone either. Content up to the Nyquist frequency is counted in full; above it the record holds nothing, so a sweep whose third harmonic passes the Nyquist frequency says so in the notes.

A line whose span reaches a tone of the specification the levels do not hold cannot be judged, and neither can a line the extraction found under its noise floor (SineLevel.below_floor), nor a line with no tone on the channel to judge against (a dead channel), nor any line of a channel whose stretch holds a sample that is not a finite number. All are kept, counted in the notes and left out of the verdict (DistortionCheck.judged), and a channel with no judged line does not pass.

References

The definition of the distortion and the estimator of the tone it is read against are set out in these.

1. Department of Defense (2019). MIL-STD-810H, Environmental Engineering Considerations and Laboratory Tests. Method 514.8, Table 514.8-III, note 3: the distortion check on a sinusoidal test, the unfiltered RMS within 5 percent of the filtered RMS, stated as about 32 percent distortion.
2. Brandt, A. (2011). Noise and Vibration Analysis: Signal Analysis and Experimental Procedures. Wiley. Order tracking, the Vold-Kalman filter among it: a component on a known frequency history read as a slowly varying amplitude by least squares.
3. Bendat, J. S., & Piersol, A. G. (2010). Random Data: Analysis and Measurement Procedures, 4th ed. Wiley. The mean square of a sum of orthogonal components as the sum of their mean squares.

Classes:

Name	Description
DistortionCheck	One tone's distortion along its sweep, on each channel read: what

Functions:

Name	Description
distortion	The distortion from two RMS readings, where what is not the tone
distortion_limit	The distortion a tolerance on the ratio allows: sqrt((1 +
decibels	A distortion in dB, 20 log10 D; minus infinity where it is zero.
check_distortion	Check one tone's distortion along its sweep: the RMS of the record

Classes

DistortionCheck dataclass

DistortionCheck(*, tone: str, dofs: list[str], dimensions: list[str], units: list[str | None], frequency: ndarray, seconds: ndarray, unfiltered: ndarray, fundamental: ndarray, residual: ndarray, judged: ndarray, excess: float, cycles: float, notes: tuple[str, ...] = ())

One tone's distortion along its sweep, on each channel read: what check_distortion returns.

frequency and seconds are the extraction's lines, ascending in frequency, each with the second it was read, from the recording's first sample. unfiltered, fundamental and residual are channels x lines of RMS in SI over cycles cycles about each line: the channel itself, the tones rebuilt on it, and the channel less the tones. judged is channels x lines, False where a line cannot be judged (the tone under its noise floor, a tone the levels do not hold live beside it, no tone on the channel, or a sample that is not a finite number). excess is the tolerance on the ratio. notes says what was not judged and why, empty when every line was. Two of these compare equal only if they are the same object, since the arrays have no single truth value.

Methods:

Name	Description
row	The row of channel dof; the first when None. Raises
worst	The largest judged distortion and the frequency it is at, on
describe	The verdict in one line, such as 'Up at 101Z+: worst 12.3 %

Attributes:

Name	Type	Description
distortion	ndarray	Channels x lines: D as a fraction, the residual's RMS over the
ratio	ndarray	Channels x lines: the ratio the distortion stands for,
decibels	ndarray	Channels x lines: D in dB (decibels).
limit	float	The distortion the tolerance allows, as a fraction
over	ndarray	Channels x lines, True where a judged line's distortion is
passed	bool	Whether every channel read has a judged line and no judged

Attributes

distortion property

distortion: ndarray

Channels x lines: D as a fraction, the residual's RMS over the rebuilt tones'; not a number where there is no tone.

ratio property

ratio: ndarray

Channels x lines: the ratio the distortion stands for, sqrt(1 + D**2), the unfiltered RMS over the tone's where the rest is orthogonal to the tones; the number a tolerance like 'within 5 %' is written on.

decibels property

decibels: ndarray

Channels x lines: D in dB (decibels).

limit property

limit: float

The distortion the tolerance allows, as a fraction (distortion_limit).

over property

over: ndarray

Channels x lines, True where a judged line's distortion is above the limit.

passed property

passed: bool

Whether every channel read has a judged line and no judged line on any channel is over the limit. A channel with nothing judged (dead, not a number, or every line unjudged) has not passed: there is no reading to pass.

Methods:

row

row(dof: str | None = None) -> int

The row of channel dof; the first when None. Raises ValueError when the check did not read it.

Source code in src/visualdynamics/core/sine_distortion.py

def row(self, dof: str | None = None) -> int:
    """The row of channel `dof`; the first when None. Raises
    ValueError when the check did not read it."""
    if dof is None:
        return 0
    if str(dof) not in self.dofs:
        raise ValueError(f'{dof} was not read; the check holds '
                         f'{self.dofs}')
    return self.dofs.index(str(dof))

worst

worst(dof: str | None = None) -> tuple[float, float] | None

The largest judged distortion and the frequency it is at, on channel dof or, when None, on any channel: (Hz, D as a fraction), or None when no line was judged.

Parameters:

Name	Type	Description	Default
dof	str	The channel; every channel read when omitted.	None

Returns:

Type	Description
tuple of float, or None	(frequency in Hz, distortion as a fraction).

Source code in src/visualdynamics/core/sine_distortion.py

def worst(self, dof: str | None = None) -> tuple[float, float] | None:
    """The largest judged distortion and the frequency it is at, on
    channel `dof` or, when None, on any channel: (Hz, D as a
    fraction), or None when no line was judged.

    Parameters
    ----------
    dof : str, optional
        The channel; every channel read when omitted.

    Returns
    -------
    tuple of float, or None
        (frequency in Hz, distortion as a fraction).
    """
    rows = (range(len(self.dofs)) if dof is None else [self.row(dof)])
    best = None
    values = self.distortion
    for r in rows:
        seen = self.judged[r] & np.isfinite(values[r])
        if not seen.any():
            continue
        k = int(np.flatnonzero(seen)[np.argmax(values[r][seen])])
        if best is None or values[r, k] > best[1]:
            best = (float(self.frequency[k]), float(values[r, k]))
    return best

describe

describe(dof: str | None = None) -> str

The verdict in one line, such as 'Up at 101Z+: worst 12.3 % (-18.2 dB) at 412 Hz, within the 32.0 % limit'. Every channel's worst when dof is None, naming any channel with no line judged.

Source code in src/visualdynamics/core/sine_distortion.py

def describe(self, dof: str | None = None) -> str:
    """The verdict in one line, such as 'Up at 101Z+: worst 12.3 %
    (-18.2 dB) at 412 Hz, within the 32.0 % limit'. Every channel's
    worst when `dof` is None, naming any channel with no line
    judged."""
    where = self.tone + ('' if dof is None else f' at {dof}')
    worst = self.worst(dof)
    if worst is None:
        return f'{where}: no line could be judged'
    rows = (range(len(self.dofs)) if dof is None else [self.row(dof)])
    over = int(sum(self.over[r].sum() for r in rows))
    judged = int(sum(self.judged[r].sum() for r in rows))
    frequency, value = worst
    words = (f'{where}: worst {100.0 * value:.1f} % '
             f'({float(decibels(value)):.1f} dB) at {frequency:.4g} Hz, ')
    if over:
        words += (f'{over} of {judged} lines over the '
                  f'{100.0 * self.limit:.1f} % limit')
    else:
        words += f'within the {100.0 * self.limit:.1f} % limit'
    silent = [self.dofs[r] for r in rows if not self.judged[r].any()]
    if silent:
        words += f'; no line judged on {", ".join(silent)}'
    return words

Functions:

distortion

distortion(unfiltered: Any, fundamental: Any) -> ndarray

The distortion from two RMS readings, where what is not the tone is orthogonal to it: sqrt((unfiltered / fundamental)**2 - 1), the RMS of everything that is not the tone over the RMS of the tone.

This is the conversion between a tolerance written on the ratio and a distortion. The check itself reads the residual directly (check_distortion), since two tones that cross are not orthogonal and two readings over different stretches are not comparable.

A ratio under one, which only the scatter of the two readings gives (the unfiltered record holds the tone and more), reads as no distortion rather than as the root of a negative number. A fundamental of zero reads as not a number.

Parameters:

Name	Type	Description	Default
unfiltered	float or array - like	The channel's RMS, unfiltered.	required
fundamental	float or array - like	The tone's RMS, in the same units.	required

Returns:

Type	Description
ndarray	D as a fraction (0.32 is 32 %), the shape the two broadcast to.

Source code in src/visualdynamics/core/sine_distortion.py

def distortion(unfiltered: Any, fundamental: Any) -> np.ndarray:
    """The distortion from two RMS readings, where what is not the tone
    is orthogonal to it: sqrt((unfiltered / fundamental)**2 - 1), the RMS
    of everything that is not the tone over the RMS of the tone.

    This is the conversion between a tolerance written on the ratio and
    a distortion. The check itself reads the residual directly
    (`check_distortion`), since two tones that cross are not orthogonal
    and two readings over different stretches are not comparable.

    A ratio under one, which only the scatter of the two readings gives
    (the unfiltered record holds the tone and more), reads as no
    distortion rather than as the root of a negative number. A
    fundamental of zero reads as not a number.

    Parameters
    ----------
    unfiltered : float or array-like
        The channel's RMS, unfiltered.
    fundamental : float or array-like
        The tone's RMS, in the same units.

    Returns
    -------
    ndarray
        D as a fraction (0.32 is 32 %), the shape the two broadcast to.
    """
    unfiltered = np.asarray(unfiltered, dtype=float)
    fundamental = np.asarray(fundamental, dtype=float)
    with np.errstate(divide='ignore', invalid='ignore'):
        ratio = np.where(fundamental > 0.0, unfiltered / fundamental, np.nan)
    return np.sqrt(np.maximum(ratio * ratio - 1.0, 0.0))

distortion_limit

distortion_limit(excess: float = DISTORTION_EXCESS) -> float

The distortion a tolerance on the ratio allows: sqrt((1 + excess)**2 - 1), 0.320 for 5 %.

Raises ValueError unless excess is more than zero.

Parameters:

Name	Type	Description	Default
excess	float	The share the unfiltered RMS may stand above the fundamental's.	`DISTORTION_EXCESS`

Returns:

Type	Description
float	D as a fraction.

Source code in src/visualdynamics/core/sine_distortion.py

def distortion_limit(excess: float = DISTORTION_EXCESS) -> float:
    """The distortion a tolerance on the ratio allows: sqrt((1 +
    excess)**2 - 1), 0.320 for 5 %.

    Raises ValueError unless `excess` is more than zero.

    Parameters
    ----------
    excess : float, default `DISTORTION_EXCESS`
        The share the unfiltered RMS may stand above the fundamental's.

    Returns
    -------
    float
        D as a fraction.
    """
    excess = float(excess)
    if not excess > 0.0:
        raise ValueError(f'a tolerance of {excess:g} on the ratio must be '
                         'more than zero')
    return float(np.sqrt((1.0 + excess) ** 2 - 1.0))

decibels

decibels(value: Any) -> ndarray

A distortion in dB, 20 log10 D; minus infinity where it is zero.

Parameters:

Name	Type	Description	Default
value	float or array - like	D as a fraction.	required

Returns:

Type	Description
ndarray	dB, the shape of value.

Source code in src/visualdynamics/core/sine_distortion.py

def decibels(value: Any) -> np.ndarray:
    """A distortion in dB, 20 log10 D; minus infinity where it is zero.

    Parameters
    ----------
    value : float or array-like
        D as a fraction.

    Returns
    -------
    ndarray
        dB, the shape of `value`.
    """
    with np.errstate(divide='ignore'):
        return 20.0 * np.log10(np.asarray(value, dtype=float))

check_distortion

check_distortion(history: Any, specification: SineSweepSpecification, levels: SineLevelSet | None = None, *, tone: str | None = None, channels: Sequence[str] | None = None, excess: float = DISTORTION_EXCESS, cycles: float | None = None) -> DistortionCheck

Check one tone's distortion along its sweep: the RMS of the record less the tracked tones against the RMS of the tones, at every line of the extraction.

The tones come from levels, a SineLevelSet that extract_sine (Extract Sine Levels) read from this recording against this specification; with none, the extraction is run here with the recording's own sine_extraction setting (or its suggestion), as the project verb runs it. Over the checked tone's stretch, every tone the levels hold is rebuilt on each channel's samples by the joint Vold-Kalman solve with the extraction's smoothing, the record's mean taken out, and the rebuilt tones subtracted. At each line the residual and the rebuilt tones are read over cycles cycles of the tone centered on the line's instant, and the distortion is the one RMS over the other. A line is flagged where it is above the limit excess sets (distortion_limit).

Raises ValueError when excess or cycles is not more than zero, when the recording is not evenly sampled, when the levels hold no such tone or channel, when the recording carries no record of a channel in the quantity its levels were read in, when the levels do not say when each line was read or were read at times the recording does not hold, or when the recording holds less of the tone than one span.

Parameters:

Name	Type	Description	Default
history	TimeHistory	The recording, evenly sampled.	required
specification	SineSweepSpecification	The sweep the drive followed.	required
levels	SineLevelSet	The extracted levels; extracted here when omitted.	None
tone	str	Which tone, by name; the first the levels hold by default.	None
channels	sequence of str	Which control channels; every channel the levels hold by default.	None
excess	float	The share the unfiltered RMS may stand above the fundamental's; 5 %, a distortion of 32 %, by default.	`DISTORTION_EXCESS`
cycles	float	The span the RMS readings are taken over, in cycles of the tone; the extraction's smoothing by default.	None

Returns:

Type	Description
DistortionCheck	The three RMS readings at every line on every channel read, which lines could be judged, and the tolerance.

Source code in src/visualdynamics/core/sine_distortion.py

def check_distortion(history: Any, specification: SineSweepSpecification,
                     levels: SineLevelSet | None = None, *,
                     tone: str | None = None,
                     channels: Sequence[str] | None = None,
                     excess: float = DISTORTION_EXCESS,
                     cycles: float | None = None) -> DistortionCheck:
    """Check one tone's distortion along its sweep: the RMS of the record
    less the tracked tones against the RMS of the tones, at every line of
    the extraction.

    The tones come from `levels`, a `SineLevelSet` that `extract_sine`
    (Extract Sine Levels) read from this recording against this
    specification; with none, the extraction is run here with the
    recording's own `sine_extraction` setting (or its suggestion), as the
    project verb runs it. Over the checked tone's stretch, every tone the
    levels hold is rebuilt on each channel's samples by the joint
    Vold-Kalman solve with the extraction's smoothing, the record's mean
    taken out, and the rebuilt tones subtracted. At each line the residual
    and the rebuilt tones are read over `cycles` cycles of the tone
    centered on the line's instant, and the distortion is the one RMS
    over the other. A line is flagged where it is above the limit
    `excess` sets (`distortion_limit`).

    Raises ValueError when `excess` or `cycles` is not more than zero, when
    the recording is not evenly sampled, when the levels hold no such tone
    or channel, when the recording carries no record of a channel in the
    quantity its levels were read in, when the levels do not say when
    each line was read or were read at times the recording does not
    hold, or when the recording holds less of the
    tone than one span.

    Parameters
    ----------
    history : TimeHistory
        The recording, evenly sampled.
    specification : SineSweepSpecification
        The sweep the drive followed.
    levels : SineLevelSet, optional
        The extracted levels; extracted here when omitted.
    tone : str, optional
        Which tone, by name; the first the levels hold by default.
    channels : sequence of str, optional
        Which control channels; every channel the levels hold by default.
    excess : float, default `DISTORTION_EXCESS`
        The share the unfiltered RMS may stand above the fundamental's;
        5 %, a distortion of 32 %, by default.
    cycles : float, optional
        The span the RMS readings are taken over, in cycles of the tone;
        the extraction's smoothing by default.

    Returns
    -------
    DistortionCheck
        The three RMS readings at every line on every channel read, which
        lines could be judged, and the tolerance.
    """
    distortion_limit(excess)                  # refuses a bad tolerance
    lines = _tone_lines(history, specification, levels, tone=tone,
                        channels=channels, cycles=cycles,
                        doing='Checking sine distortion', reading='the RMS')
    level, chosen, cycles = lines.level, lines.chosen, lines.cycles
    wanted, picks, rows = lines.wanted, lines.picks, lines.rows
    first, last, phase = lines.first, lines.last, lines.phase
    ends, starts, seconds = lines.ends, lines.starts, lines.seconds
    shared, laid = lines.shared, lines.laid
    dt = lines.dt

    count = len(seconds)
    unfiltered = np.full((len(rows), count), np.nan)
    fundamental = np.full((len(rows), count), np.nan)
    residual = np.full((len(rows), count), np.nan)
    notes = []
    for c, row in enumerate(rows):
        x = np.real(np.asarray(history.ordinate[row, first:last],
                               dtype=float))
        bad = int((~np.isfinite(x)).sum())
        if bad:
            notes.append(f'{wanted[c]}: {bad} sample{"s" * (bad != 1)} over '
                         "the tone's stretch "
                         f'{"are" if bad != 1 else "is"} not a finite '
                         'number, so no line is judged')
            continue
        x = x - np.mean(x)
        arguments, offsets, envelopes = _rebuilt(
            x, laid, dt, lines.smoothing,
            [drifted for *_rest, drifted in laid])
        tones = np.zeros(len(x))
        for argument, offset, envelope in zip(arguments, offsets,
                                              envelopes, strict=True):
            tones[offset:offset + len(argument)] += np.real(
                envelope * np.exp(1j * argument))
        unfiltered[c] = np.sqrt(_mean_squares(x, phase, ends, starts))
        fundamental[c] = np.sqrt(_mean_squares(tones, phase, ends, starts))
        residual[c] = np.sqrt(np.maximum(
            _mean_squares(x - tones, phase, ends, starts), 0.0))

    floor = (np.zeros((len(rows), count), dtype=bool)
             if level.below_floor is None
             else np.asarray(level.below_floor, dtype=bool)[picks])
    toneless = np.isfinite(fundamental) & ~(fundamental > 0.0)
    judged = (~floor & ~shared[None, :] & np.isfinite(residual)
              & (fundamental > 0.0))
    if floor.any():
        notes.append(f'{int(floor.sum())} line{"s" * int(floor.sum() != 1)} '
                     'under the noise floor, not judged')
    if shared.any():
        notes.append(f'{int(shared.sum())} line{"s" * int(shared.sum() != 1)} '
                     'shared with a tone the levels do not hold, not judged')
    for c in range(len(rows)):
        none = int(toneless[c].sum())
        if none:
            notes.append(f'{wanted[c]}: {none} line{"s" * (none != 1)} '
                         'with no tone on the channel to judge against, '
                         'not judged')
    reach = 0.5 / dt / HARMONIC
    if float(np.max(level.abscissa)) > reach:
        notes.append(f'above {reach:.4g} Hz the {HARMONIC_WORD} harmonic is '
                     'past the Nyquist frequency, so the record cannot hold '
                     'it and the check cannot see it')
    return DistortionCheck(
        tone=chosen.name, dofs=wanted,
        dimensions=[history.ordinate_dim[row] for row in rows],
        units=[history.ordinate_unit[row] for row in rows],
        frequency=np.asarray(level.abscissa, dtype=float), seconds=seconds,
        unfiltered=unfiltered, fundamental=fundamental, residual=residual,
        judged=judged, excess=float(excess), cycles=cycles,
        notes=tuple(notes))