Signal Analysis

This page collects the theory behind the measurement chain itself: the standardized fractional-octave bands and the time-domain filter banks that implement them, the frequency weighting curves, time integration, level, event and exposure metrics, sound intensity, and the GUM uncertainty framework that underpins every measured quantity. It is part of the theory reference.

Octave Band Frequencies (ANSI S1.11 / IEC 61260)

The mid-band frequencies ($f_m$) and edges ($f_1$, $f_2$) use a base-10 ratio:

$$G=10^{0.3}$$

$G=10^{0.3}\approx 1.9953$ is the octave ratio IEC 61260-1 and ANSI S1.11 adopt in preference to an exact factor of 2, because a base-10 ratio keeps the band grid commensurate with the decade grid tabulated data uses and puts 1 kHz exactly on it. The band index $x$ is an integer, anchored so that $x=30$ lands on the 1 kHz reference, and the mid-band frequency depends on whether the fraction $b$ is odd or even:

Mid-band:

$$f_m=1000\cdot G^{(x-30)/b}(b\text{ odd}),f_m=1000\cdot G^{(2x-59)/(2b)}(b\text{ even})$$

The distinction is not cosmetic. For odd $b$ the reference sits at a band centre, so a one-third-octave bank has a band centred exactly on 1000.00 Hz; for even $b$ it sits on a band edge, so a sixth-octave bank has bands at 944.1 Hz and 1059.3 Hz and none at 1 kHz at all. A formula written only for odd $b$ cannot produce the band set the library returns for fraction=2, fraction=6 or fraction=12.

Band edges:

$$f_1=f_m\cdot G^{-1/(2b)},f_2=f_m\cdot G^{1/(2b)}$$

Exact centre against printed label. The $f_m$ above is the exact midband frequency; what an instrument prints — and what nominal_frequencies returns as its fourth item, or octave_filter(..., nominal=True) — is the ISO 266 preferred number nearest to it. The one-third-octave band whose exact centre is 1258.9 Hz is labelled 1250 Hz, and the sixth-octave band centred on 1059.3 Hz is labelled 1.06k. Compare band levels between tools by index or by exact centre, never by matching label strings: different tools round the labels differently.

Frequency Resolution vs FFT Bin Spacing

octave_filter is a time-domain fractional-octave filter bank, not an FFT or Welch spectrum estimator. Therefore, its result does not have a frequency resolution in the fs / nfft sense.

For fraction=3, the output contains one scalar level per one-third-octave band. The relevant frequency granularity is the standardized band definition: center frequency, lower edge, and upper edge. Because fractional-octave bands are logarithmically spaced, their absolute bandwidth in Hz grows with frequency while their relative bandwidth remains approximately constant.

For example, with fraction=3 and limits=[12, 20000], the exact one-third-octave band around 1 kHz is approximately:

Nominal band	Lower edge	Center	Upper edge	Bandwidth
1 kHz	891.25 Hz	1000.00 Hz	1122.02 Hz	230.77 Hz

You can inspect the exact bands with:

from phonometry import filters

fc, fl, fu, labels = filters.nominal_frequencies(fraction=3, limits=[12, 20000])
for label, center, lower, upper in zip(labels, fc, fl, fu):
    print(label, center, lower, upper, upper - lower)

If you need narrowband FFT bins for tonal inspection, run Welch/FFT on the original signal and use the phonometry band edges as masks:

import numpy as np
from scipy import signal
from phonometry import filters

fs = 100_000
# any 1D pressure signal in Pa (synthesized here so the example runs)
pressure_signal_pa = 0.02 * np.random.default_rng(0).standard_normal(fs)
x = pressure_signal_pa

# Standardized one-third-octave levels from phonometry.
levels, centers = filters.octave_filter(
    x,
    fs=fs,
    fraction=3,
    limits=[12, 20_000],
)

# Same standardized band definitions, including lower/upper edges.
fc, fl, fu, labels = filters.nominal_frequencies(fraction=3, limits=[12, 20_000])

# Narrowband Welch estimate on the original signal.
nperseg = min(2**15, len(x))
freq_bins, psd = signal.welch(
    x,
    fs=fs,
    window="hann",
    nperseg=nperseg,
    noverlap=nperseg // 2,
    scaling="density",
)

# Example: list the Welch bins inside the one-third-octave band closest to 1 kHz.
band_index = int(np.argmin(np.abs(np.asarray(fc) - 1000.0)))
in_band = (freq_bins >= fl[band_index]) & (freq_bins <= fu[band_index])

print("Selected one-third-octave band:", labels[band_index])
print("Welch bin spacing:", freq_bins[1] - freq_bins[0], "Hz")
for f, pxx in zip(freq_bins[in_band], psd[in_band]):
    print(f, pxx)

The two objects on one axis, for a six-tone signal at 20, 100, 500, 2000, 4000 and 15 000 Hz. The grey trace is a Welch PSD ($f_s$ = 48 kHz, nperseg = 8192, so a fixed 5.86 Hz bin everywhere); the markers are the standardized one-third-octave levels of the same signal. The bin width never changes and the band width does: 4.60 Hz at the 20 Hz band, narrower than one bin, against 230.77 Hz at 1 kHz and 3657 Hz at 16 kHz. That is why the top bands each swallow hundreds of bins while the bottom ones sit inside a single one, and why the two answers cannot be converted into each other. (The PSD trace is offset vertically for legibility, so read its shape, not its level.)

This keeps the two concepts separate: phonometry gives standardized fractional-octave levels, while Welch gives narrowband FFT bins. With fs=100000 and nperseg=2**15, the Welch bin spacing is about 3.05 Hz. Window choice and overlap affect leakage and averaging variance, but they do not change the bin spacing of each FFT segment.

When sigbands=True, octave_filter can also return the time-domain waveform filtered by each band. Applying Welch/FFT to one selected filtered waveform can be useful as a diagnostic view of the content inside that filtered band, but it does not recover FFT bins from the scalar band levels.

Magnitude Responses

The magnitude response of a filter is $|H(j\omega )|$, the modulus of its transfer function evaluated on the imaginary axis — what the filter does to the amplitude of a steady sinusoid at each frequency, before any band is filtered. The library implements standard classical filter prototypes:

1. Butterworth: Maximally flat passband.

$$|H(j\omega )|=\frac{1}{\sqrt{1+(\omega /{\omega }_c{)}^{2n}}}$$

2. Chebyshev I: Equiripple in passband, steeper roll-off.

$$|H(j\omega )|=\frac{1}{\sqrt{1+{ϵ}^2{T}_n^2(\omega /{\omega }_c)}}$$

3. Chebyshev II: Inverse Chebyshev, equiripple in stopband, flat passband.

$$|H(j\omega )|=\frac{1}{\sqrt{1+\frac{1}{{ϵ}^2{T}_n^2({\omega }_{stop}/\omega )}}}$$

4. Elliptic: Equiripple in both, maximum selectivity.

$$|H(j\omega )|=\frac{1}{\sqrt{1+{ϵ}^2{R}_n^2(\omega /{\omega }_c,L)}}$$

5. Bessel: Maximally flat group delay (linear phase).

$$H(s)=\frac{{\theta }_n(0)}{{\theta }_n(s/{\omega }_0)}$$

(Where ${\theta }_n$ is the reverse Bessel polynomial)

The five are not interchangeable for a reported band level, and the table says which is which:

Prototype	Passband	Stopband	Where its edge frequency is defined	IEC 61260-1 class reachable	Where it belongs
Butterworth	Monotone, maximally flat	Monotone	−3 dB	Class 1 with the defaults	The bank default, and any reported level
Chebyshev II	Flat	Equiripple	−3 dB (mapped from the stopband edge)	Class 1 once attenuation ≥ 70 dB	A steeper skirt at the same order
Chebyshev I	Equiripple (ripple dB)	Monotone	Passband ripple edge	None	Exploratory filtering
Elliptic	Equiripple	Equiripple	Passband ripple edge	None	Exploratory filtering
Bessel	Maximally flat group delay	Slow roll-off	−3 dB	None	Preserving a burst or decay envelope

The bank default is Butterworth, order 6, and it is the only architecture that clears the IEC 61260-1:2014 class 1 mask with the default parameters; Chebyshev II joins it once the stopband attenuation is raised to 70 dB or more. Chebyshev I, elliptic and Bessel fail the mask on passband ripple or roll-off whatever the parameters, so they belong to exploratory filtering and not to a band level that will be reported. The physical reason is worth carrying: a flat passband means an in-band tone is measured with no ripple error wherever it sits, and a monotone stopband makes leakage from a loud neighbouring band predictable. An equiripple design reaches the same mask at a lower order but biases a narrowband source by up to the ripple value depending on where in the band it falls, and smears impulsive signals through group-delay distortion near the edges. Bessel does the opposite — it preserves the shape of a decay or a burst envelope, which matters when the band signal feeds a reverberation-time fit — but will not meet class 1 at any reasonable order. Raising the order tightens the mask fit and lengthens the filter transient, which the settling rule below quantifies. filter_class_compliance checks a designed bank against the masks; see Filter Class Verification.

The same 1 kHz octave band designed five ways. Look at three things: the flatness of the passband top (Butterworth and Chebyshev II flat, Chebyshev I and elliptic rippling by the design ripple), the steepness of the skirt just outside the edges, and the shape of the deep stopband. The zoom at the crossover is where the band-edge rule below becomes visible — the two equiripple designs are not at −3 dB there, because that is not where their edge is defined.

Band-edge placement

Every architecture is designed on the band edges themselves, at whichever gain that architecture defines its edge frequency: −3 dB for Butterworth, Chebyshev II and Bessel, and the passband ripple edge (ripple dB, 0.1 dB by default) for Chebyshev I and elliptic, which are equiripple in the passband and have no −3 dB point there. Two cases need special handling:

• Chebyshev II: scipy’s Wn is the stopband edge. phonometry maps the desired −3 dB edges to stopband edges analytically (the prototype transition ratio is $\cosh (\mathrm{acosh}(\sqrt{10^{A/10}-1})/N)$), applying the lowpass→bandpass transform in the pre-warped bilinear domain so the mapping stays exact for decimated bands close to Nyquist.
• Bessel: designed with norm="mag", which defines the −3 dB point exactly at Wn (the phase norm would shift the edges to roughly −10 dB).

Filter Bank Design & Numerical Stability

A single high-order band-pass at 16 Hz and 192 kHz puts its poles within about $10^{-4}$ of the unit circle, where float64 coefficients no longer resolve them and the difference equation degenerates. Two measures keep every band of the bank away from that regime:

1. Second-Order Sections (SOS): each band is a cascade of biquads rather than one high-order transfer function, so a coefficient rounding error perturbs one pole pair instead of all of them at once.
2. Multi-rate Decimation: whenever half the sample rate still clears the band’s upper edge by a factor of 1.25, the signal is downsampled (decimated) by $M=\lfloor (f_s/2)/(1.25f_{\text{upper}})\rfloor$ before filtering — which applies to most bands rather than only the low ones (29 of the 33 one-third-octave bands at 48 kHz), so each band is filtered near its own natural rate and its poles sit well inside the unit circle. The band level is computed on the decimated signal; only sigbands=True brings the band signal back to the input rate, with resample_poly(M, 1). Chebyshev II banks reserve extra decimation headroom so their stopband edges stay below the decimated Nyquist.
3. Settling, and two effects that bias a band level. A band filter needs roughly $5/f_{\text{low}}$ seconds to settle, so a 1 s record already truncates the 16 Hz band’s response and reads several tenths of a decibel low while the same record is ample above 500 Hz; the analysed record must exceed both $5/f_{\text{low}}$ and, for a decay measurement, the reverberation time of the band. mode='peak' reads the filter’s own onset overshoot — up to about 1 dB for an abruptly starting tone — so a peak band level of a gated signal should discard the first $5/f_{\text{low}}$ seconds. And zero_phase=True runs the cascade forwards and backwards: group delay disappears, which is what ISO 3382-2 clause 7.3 wants for a decay, but the effective passband narrows and the measured broadband band level falls by roughly 0.2 to 0.3 dB per band (a pure in-band tone is unaffected, sitting where both passes are ≈ 0 dB), so a level reported for conformity should come from the single forward pass.

See the Filter Banks guide for usage, the Filter Gallery for the five architectures side by side, and Filter Class Verification for checking a designed bank against the IEC 61260-1 masks.

Weighting Curves (IEC 61672-1)

A frequency weighting is a fixed approximation to the ear’s frequency response at one loudness, obtained by inverting an equal-loudness contour. A follows the inverse of the 40 phon contour, C the inverse of a contour near 100 phon — hence nearly flat through the mid band, with only the band-limiting roll-offs left — and Z is a defined flat response between stated band limits, not the absence of a filter. That origin is also the domain of validity: because A is anchored at 40 phon it progressively under-reads the low-frequency content of loud sources, which is exactly why loudness models, C-weighted peaks and the G curve below exist. See Equal-loudness contours for the contours themselves.

The A-weighting transfer function:

$$R_A(f)=\frac{12194^2\cdot f^4}{(f^2+20.6^2)\sqrt{(f^2+107.7^2)(f^2+737.9^2)}(f^2+12194^2)}$$

The four corner frequencies are the pole pairs at 20.6 Hz and 12194 Hz that band- limit the response and the single poles at 107.7 Hz and 737.9 Hz that build the low-frequency slope; the $+2.00$ dB term is the normalisation that forces the response to exactly 0 dB at 1 kHz.

$$A(f)=20{\log }_{10}(R_A(f))+2.00$$

The digital filter is obtained from the analog poles/zeros via the bilinear transform. Because the bilinear transform compresses frequencies near Nyquist, the default high_accuracy mode designs and runs the filter at an internally oversampled rate (≥ 144 kHz); see Frequency Weighting.

The three IEC 61672-1 weighting curves realized by the library, with the small positive region of the A curve magnified. The special B, D and AU curves are charted in Special Weightings.

Time Integration

Implemented as a first-order IIR exponential integrator:

$$y[n]=\alpha \cdot x^2[n]+(1-\alpha )\cdot y[n-1]$$ $$\alpha =1-e^{-1/(f_s\cdot \tau )}$$

This is the digital image of the analogue RC detector of the classical sound level meter: the squared pressure charges a capacitor through a resistor, so the displayed mean square is a one-pole low pass on $p^2$ with time constant $\tau$, and $\alpha =1-e^{-1/(f_s\tau )}$ is its impulse-invariant discretisation — exact at any $f_s$, unlike the $1/(f_s\tau )$ approximation often seen. The two standard constants are $\tau =125$ ms (Fast) and $\tau =1$ s (Slow). A step input reaches within 1 dB of its final value after about $1.6\tau$ and is 0.03 dB short at $5\tau$, which is why levels read from the first few time constants of a record are biased low.

Impulse is not this filter. It cannot be: a single $\alpha$ cannot hold a peak, and the figure below shows Impulse holding one. IEC 61672-1 defines the I weighting as a fast 35 ms detector followed by a 1.5 s decay, implemented here as an asymmetric detector — it charges with $\tau =35$ ms whenever the input exceeds the current state and discharges with $\tau =1.5$ s otherwise — so a transient is captured quickly and then held for over a second. That makes I neither an r.m.s. nor an energy average, which is why it is reserved for the specific regulations that name it and must never be energy-summed across intervals.

The default initial condition is $y[-1]=0$. Use initial_state='first' to start from the first input energy, or pass a scalar/array with the previous mean-square output state. See Why phonometry for the IEC 61672-1 tone-burst verification of this implementation.

The exponential integrator at the three standard time constants: Fast follows a burst, Slow smooths it and Impulse holds its peak.

G-weighting (ISO 7196)

The G curve extends frequency weighting into the infrasound range. ISO 7196:1995 Table 1 (p. 2) defines it by four zeros at the origin and four complex-conjugate pole pairs, given as coordinates in Hz (multiplied by $2\pi$ to obtain rad/s):

$$z_{\mathrm{1..4}}=0,p=2\pi \left\{-0.707\pm j0.707,-19.27\pm j5.16,-14.11\pm j14.11,-5.16\pm j19.27\right\}\text{Hz}$$

The gain $k$ is chosen so that the response is exactly 0 dB at 10 Hz (clause 4):

$$k=\left|\frac{{\prod }_i(j{\omega }_{10}-p_i)}{{\prod }_i(j{\omega }_{10}-z_i)}\right|,{\omega }_{10}=2\pi \cdot 10\text{rad/s}$$

The four zeros against eight poles shape the characteristic response: a rise of approximately +12 dB/octave between 1 Hz and 20 Hz, with roll-offs of approximately 24 dB/octave below 1 Hz and above 20 Hz.

The shape those four zeros and four pole pairs make, against the ISO 7196 Table 2 nominal values: 0 dB at the 10 Hz anchor, the +12 dB/octave climb through the infrasound decade below it, and the two 24 dB/octave roll-offs that fence the curve off below 1 Hz and above 20 Hz.

Infrasound needs its own curve because near the hearing threshold the perceived loudness of very-low-frequency tones grows much more steeply with sound pressure level than at mid frequencies (a small dB increase above threshold produces a large loudness jump), so the A curve (anchored at 1 kHz) grossly misrepresents infrasonic annoyance.

Since G acts on 0.25 Hz – 315 Hz, far below the Nyquist frequency at audio rates, the frequency warping of the plain bilinear transform (applied without prewarping) is negligible there: about 0.014 % at 315 Hz for $f_s=48$ kHz, under 0.01 dB on the response. The internal oversampling used for the A/C designs (whose action extends to 16 kHz) is therefore not applied.

See the Special Weightings guide for usage.

Event and dose metrics

Sound exposure level (SEL; $L_{AE}$ with A-weighting, IEC 61672-1:2013) normalizes the energy of a discrete event (aircraft flyover, train pass) to a 1 s reference duration:

$$\mathrm{SEL}=L_{eq,T}+10{\log }_{10}\left(\frac{T}{T_0}\right),T_0=1\text{s}$$

What the formula does to an event: the pass-by is replaced by a one-second block of the same total energy, which is why SEL exceeds the event’s $L_{eq}$ whenever the event lasts longer than a second, and why two events of the same SEL are interchangeable in a dose even when one is loud and short and the other quiet and long.

Sound exposure $E$ (IEC 61252, 3.1) is the time integral of the squared A-weighted sound pressure, expressed in pascal-squared hours:

$$E={\int }_0^T{p}_A^2(t)dt=\overline{p_A^2}\cdot T[{\text{Pa}}^2\text{h}]$$

When the recording is a representative sample of a longer shift, $E$ scales the measured mean square by the actual exposure duration. The normalized 8 h level (IEC 61252, 3.3) converts exposure to the steady level that carries the same energy over a nominal working day:

$$L_{EX,8h}=10{\log }_{10}\left(\frac{E}{8\text{h}\cdot p_0^2}\right),p_0=20\mu \text{Pa}$$

It is identical to $L_{EP,d}$ of Directive 86/188/EEC and $L_{EX,8h}$ of ISO 1999 (IEC 61252, 3.3 NOTES 5–6). The anchor of IEC 61252 (3.3 NOTE 4): an exposure of 3.2 Pa²h corresponds to $L_{EX,8h}$ of exactly 90 dB. ISO 9612 prescribes how to sample a working day to obtain it and how to bound the result; see Occupational noise exposure.

$L_{Cpeak}$ (IEC 61672-1:2013, subclause 5.13) is the absolute maximum of the C-weighted sound pressure expressed in dB, $L_{Cpeak}=20{\log }_{10}(\max |p_C(t)|/p_0)$, the quantity behind the 135/137/140 dB(C) occupational action limits. The implementation is verified against the one-cycle and half-cycle reference responses of Table 5.

See the Levels guide for usage and the Calibration guide for absolute-scale setup.

Sound intensity (IEC 61043)

Sound intensity is the time-averaged acoustic power flux $I=\overline{pu}$. The particle velocity follows from Euler’s equation (linearized conservation of momentum):

$${\rho }_0\frac{\partial u}{\partial t}=-\frac{\partial p}{\partial r}$$

A p-p probe approximates the pressure gradient by the finite difference of two microphones a spacer distance $\Delta r$ apart (IEC 61043:1993, definition 3.2):

$$p=\frac{p_1+p_2}{2},u=-\frac{1}{{\rho }_0\Delta r}\int (p_2-p_1)dt,I=\overline{pu}$$

For stationary signals the same estimator has an exact frequency-domain form through the imaginary part of the one-sided cross spectrum $G_{12}$ of the two pressures; the implementation estimates it with Welch-averaged, Hann-windowed segments:

$$I(f)=-\frac{\mathrm{Im}\{G_{12}(f)\}}{2\pi f{\rho }_0\Delta r}$$

The finite difference underestimates the true plane-wave intensity by the factor

$$\frac{\sin (k\Delta r)}{k\Delta r},k=\frac{2\pi f}{c}$$

IEC 61043 clause 7.3 specifies the probe intensity response with exactly this argument and Table 3 tabulates it (e.g. −10.5 dB at 6.3 kHz for a 25 mm spacer). Below $f=0.1c/\Delta r$ (i.e. $k\Delta r$ under 0.63) the bias stays within about 0.3 dB; bias_correction provides the reciprocal factor per band and max_valid_frequency the bound.

The usable band is bounded below as well, and the library exposes no bound for it. Halving the spacer doubles the top of the band and simultaneously halves the true pressure difference between the two microphones, because that difference shrinks as $k\Delta r$. Any residual inter-channel phase mismatch is therefore amplified in the same proportion, and at low frequency it eventually swamps the true intensity — which is why a probe ships with a set of spacers rather than one. As orders of magnitude, a 50 mm spacer covers roughly the 50 Hz to 1.25 kHz thirds and a 12 mm spacer roughly 200 Hz to 5 kHz, so a full-range survey is two passes stitched at the overlap. The consequence for validity is that a band fails not because it lies above max_valid_frequency alone, but whenever the field reactivity climbs to within the bias-error factor of the instrument’s residual index — which happens first at the low-frequency end.

That residual index is measured, not looked up. A residual intensity testing device presents the same sound pressure to both microphones at once, both ports opening onto one small cavity, so whatever intensity the probe then indicates is pure phase-mismatch residual and ${\delta }_{pI0}=L_p-L_I$ in that state, per band. It is part of the checks made before a survey, with the microphones in the same ports and orientation they will be used in, and any phase compensation applied there stays applied for the measurement.

The pressure-intensity index ${\delta }_{pI}=L_p-L_I$ measures how reactive the field is: in a free plane progressive wave it equals $10{\log }_{10}({\rho }_{0}c/400)=0.14$ dB, while large values flag reactive or noisy fields in which the inter-channel phase error dominates. ISO 9614-1:1993 Annex A generalizes it over a measurement surface as the indicator F2 (with F3 for negative partial power and F4 for field non-uniformity), and the instrument’s dynamic capability $L_d={\delta }_{pI0}-K$ (pressure-residual intensity index minus the bias error factor: 10 dB for grades 1/2, 7 dB for grade 3) must exceed F2 for the measurement to be valid (criterion 1).

See the Sound Intensity guide for usage.

The p-p estimator in the two limiting fields: the gap between $L_p$ and $L_I$ is the pressure-intensity index that flags reactive fields.

Measurement uncertainty (ISO/IEC Guide 98-3: GUM and Supplement 1)

Domain budgets like ISO 12999-1 and ISO 9612 Annex C are instances of the general framework of the GUM (ISO/IEC Guide 98-3:2008). Given a measurement model $y=f(x_1,\dots ,x_N)$, the law of propagation of uncertainty (clause 5) combines the input standard uncertainties through sensitivity coefficients:

$$u_c^2(y)=\sum _{i=1}^N{\left(\frac{\partial f}{\partial x_i}\right)}^2{u}^2(x_i),$$

generalized to $(c\odot u{)}^{\top }r(c\odot u)$ for correlated inputs. The sensitivities are obtained by central differences on the user’s model callable (step scaled to $10^{-3}$ of each input uncertainty), so no hand-derived partials are needed. Type B inputs enter through the clause 4.3 half-width rules: rectangular $a/\sqrt{3}$ (4.3.7), triangular $a/\sqrt{6}$ (4.3.9), U-shaped $a/\sqrt{2}$. The expanded uncertainty $U=k{u}_c$ takes $k$ from the t-distribution at the Welch–Satterthwaite effective degrees of freedom (Annex G.4):

$${\nu }_{\mathrm{eff}}=\frac{u_c^4}{{\sum }_i{u}_i^4/{\nu }_i}.$$

Supplement 1 (ISO/IEC Guide 98-3-1:2008) propagates the full distributions instead: $10^6$ Monte Carlo draws (clause 6.4) through the same model give $u(y)$ and the probabilistically symmetric coverage interval from the $\frac{1}{2}(1\mp p)$ fractiles (clause 7.7): the route when the model is non-linear or the output visibly non-Gaussian. The Guides’ own examples are reproduced: the four-term additive model gives $u_c=2.0$ and the Monte Carlo 95 % interval $\pm 3.88$ of Supplement 1 clause 9.2/Table 3 (four rectangular inputs; the output is nearly trapezoidal, not Gaussian, so the interval is narrower than $\pm 1.96u$), and the GUM Annex H.1 end-gauge example gives $k=t_{0.99}({\nu }_{\mathrm{eff}}=16)=2.92$ and $U_{99}=93$ nm.

The two routes on one problem — an A-weighted level, not the Supplement 1 four-term example quoted above. Left is the law of propagation as a budget: one bar per input, so the term worth reducing is visible. Right is the Supplement 1 route: the Monte Carlo output distribution with the GUM Gaussian drawn over it and the 95 % coverage interval shaded. Here the two agree, which is what clause 8 calls validation; where the model is non-linear or the output visibly non-Gaussian the histogram departs from the curve and the interval is read off the fractiles instead.

See the GUM Uncertainty guide for usage.

References

• Fahy, F. J. (1995). Sound intensity (2nd ed.). E&FN Spon. https://doi.org/10.4324/9780203475386ISBN 978-0-419-19810-9. The physics of the p-p estimator: active intensity, the finite-difference bias and the phase-mismatch error budget.
• International Electrotechnical Commission. (1993). Electroacoustics — Instruments for the measurement of sound intensity — Measurements with pairs of pressure sensing microphones (IEC 61043:1993). The p-p instrument standard, adopted in Europe as EN 61043:1994: the cross-spectral estimator and the pressure-residual intensity index behind the dynamic capability.
• International Electrotechnical Commission. (2013). Electroacoustics — Sound level meters — Part 1: Specifications (IEC 61672-1:2013). The A/C/Z weighting curves, the exponential time integration and the SEL and LCpeak definitions of the event-metric section.
• International Electrotechnical Commission. (2014). Electroacoustics — Octave-band and fractional-octave-band filters — Part 1: Specifications (IEC 61260-1:2014). The base-10 mid-band and band-edge definitions of the octave-band section and the class masks the banks are verified against.
• International Organization for Standardization. (1993). Acoustics — Determination of sound power levels of noise sources using sound intensity — Part 1: Measurement at discrete points (ISO 9614-1:1993). The surface indicators F2–F4 that generalize the pressure-intensity index over a measurement surface.
• Joint Committee for Guides in Metrology. (2008). Evaluation of measurement data — Guide to the expression of uncertainty in measurement (JCGM 100:2008, the GUM). BIPM. https://doi.org/10.59161/JCGM100-2008EThe law of propagation of uncertainty and the Welch–Satterthwaite expanded uncertainty of the GUM section. The linked PDF is the free download.
• Joint Committee for Guides in Metrology. (2008). Evaluation of measurement data — Supplement 1 to the "Guide to the expression of uncertainty in measurement" — Propagation of distributions using a Monte Carlo method (JCGM 101:2008). BIPM. https://doi.org/10.59161/JCGM101-2008The Monte Carlo propagation of distributions and its coverage-interval construction. The linked PDF is the free download.
• Oppenheim, A. V., & Schafer, R. W. (2010). Discrete-time signal processing (3rd ed.). Pearson. ISBN 978-0-13-198842-2. The digital-filter theory behind the SOS cascades, the bilinear transform and the multirate decimation of the filter-bank design section.
• Smith, J. O. (n.d.). Introduction to digital filters with audio applications. Center for Computer Research in Music and Acoustics (CCRMA), Stanford University. A free online-book companion treatment of the classical filter prototypes and their magnitude responses.