Skip to content

signals.correlation

Correlation analysis and time-delay estimation.

Auto- and cross-correlation estimators with the three standard normalizations, generalized cross-correlation (GCC) time-delay estimation with the Knapp & Carter weightings, and sub-sample delay location for impulse responses, following Bendat & Piersol, Random Data: Analysis and Measurement Procedures (4th ed., 2010) and Knapp & Carter (1976):

  • correlation estimates computed via FFT with zero padding so the circular product never wraps (B&P Section 11.4.2, Eq. 11.95), with the biased, unbiased (Eq. 11.96) and coefficient (Eq. 5.16) normalizations, plus the large-T normalized random error of the estimate for bandwidth-limited Gaussian data, (Eqs. 8.109/8.112), exposed as correlation_random_error;
  • time-delay estimation: the peak of the cross-correlation locates the delay of a common signal between two sensors (B&P Section 5.1.4, Eq. 5.21). time_delay implements the direct correlator, the weighted-phase-slope estimator of the cross-spectrum (Eq. 5.101b) and the generalized cross-correlation of Knapp & Carter (1976): the averaged cross-spectrum is weighted by before the inverse transform, with the Table I processors 'roth' (), 'scot' (), 'phat' () and the maximum-likelihood 'ml' (Hannan-Thomson) weighting that attains the Cramér-Rao bound;
  • sub-sample peak location by three-point parabolic interpolation, optionally after band-limited local upsampling of the correlation around its peak, and the peak-location uncertainty of the delay estimate, (Eq. 8.129), with the 95 % interval (Eq. 8.130);
  • impulse-response utilities: the sub-sample arrival time of a single IR (its peak is the cross-correlation with an ideal impulse) and the alignment of an IR pair by the estimated delay, applied as an exact fractional shift in the frequency domain.

The GCC estimators run on the same Welch core (segmentation, tapering, overlap policy) as phonometry.signals.spectra, so a GCC and a cross-spectral density computed with the same segment length are mutually consistent bin by bin.

Auto-generated from the source docstrings by scripts/generate_api_docs.py (make api-docs). Do not edit by hand.

align_impulse_responses(
ir: NDArray[np.float64] | list[float],
reference: NDArray[np.float64] | list[float],
fs: float,
*,
interpolation: _Interpolation = 'parabolic',
upsample: int = 8,
) -> AlignedImpulseResponseResult

Align an impulse response onto a reference by its estimated delay.

Estimates the sub-sample delay of ir relative to reference (impulse_response_delay) and removes it with an exact band-limited fractional shift (frequency-domain phase ramp over a zero-padded record). Use it to average IR ensembles or to compare measurements taken at slightly different distances.

Parameters

NameDescription
irImpulse response to align, 1-D.
referenceReference impulse response, same length.
fsSample rate, in Hz.
interpolation'parabolic' (default) or 'none'.
upsampleInteger local-upsampling factor (default 8).

Returns: An AlignedImpulseResponseResult.

Raises

ExceptionWhen
ValueErrorIf the inputs or parameters are invalid.
AlignedImpulseResponseResult(
aligned: NDArray[np.float64],
reference: NDArray[np.float64],
delay: float,
delay_samples: float,
fs: float,
)

An impulse response aligned onto a reference.

Attributes

NameDescription
alignedThe input IR advanced by the estimated delay (exact fractional shift applied as a frequency-domain phase ramp over a zero-padded record, so nothing wraps around).
referenceThe reference IR.
delayEstimated delay removed from the IR, in seconds.
delay_samplesThe same delay in (fractional) samples.
fsSample rate, in Hz.
AlignedImpulseResponseResult.plot(
ax: Axes | None = None,
*,
language: str = 'en',
**kwargs: Any,
) -> Axes

Plot the reference and the aligned impulse response.

Parameters

NameDescription
languageLabel language, "en" (default) or "es".
correlation(
x: NDArray[np.float64] | list[float],
y: NDArray[np.float64] | list[float] | None = None,
fs: float = 1.0,
*,
normalization: _Normalization = 'unbiased',
max_lag: float | None = None,
) -> CorrelationResult

Auto- or cross-correlation estimate with a chosen normalization.

Computed via zero-padded FFT so the circular product never wraps (B&P Section 11.4.2). y=None gives the autocorrelation of x. The sign convention follows B&P Eq. 5.19-5.20: with the estimate peaks at .

Normalizations:

  • 'biased' - the raw lag sums divided by N; tapers toward the record ends and stays bounded by ;
  • 'unbiased' - divided by (Eq. 11.96), an unbiased estimate of whose variance grows toward the ends;
  • 'coefficient' - the correlation coefficient function over the mean-removed records (Eq. 5.16).

Parameters

NameDescription
xFirst signal, 1-D.
ySecond signal, same length, or None for autocorrelation.
fsSample rate, in Hz.
normalizationSee above (default 'unbiased').
max_lagLargest lag magnitude to keep, in seconds (default: the full N-1 samples).

Returns: A CorrelationResult.

Raises

ExceptionWhen
ValueErrorIf the inputs or parameters are invalid.
correlation_random_error(
coefficient: float,
signal_bandwidth: float,
duration: float,
) -> float

Normalized random error of a correlation estimate.

B&P Eqs. 8.109/8.112: for bandwidth-limited Gaussian data of bandwidth B observed for T seconds, with the correlation coefficient at the lag of interest. At the zero lag of an autocorrelation () this is (Eq. 8.111). Valid for and (Section 8.4.1).

For the two-detector time-delay problem , (Section 8.4.2) the peak coefficient is , which reproduces the book’s Example 8.5: Hz, s, give .

Parameters

NameDescription
coefficientCorrelation coefficient at the lag.
signal_bandwidthSignal bandwidth B, in Hz.
durationRecord length T, in seconds.

Returns: Normalized random error (dimensionless; inf at ).

Raises

ExceptionWhen
ValueErrorIf the bandwidth or duration is not positive, or the coefficient is outside [-1, 1].
CorrelationResult(
lags: NDArray[np.float64],
values: NDArray[np.float64],
coefficient: NDArray[np.float64],
normalization: str,
kind: str,
fs: float,
n_samples: int,
duration: float,
)

Auto- or cross-correlation estimate (B&P Sections 5.1, 8.4, 11.4).

Attributes

NameDescription
lagsLag axis , in seconds, symmetric about zero. Positive lag means the second signal is delayed relative to the first (, B&P Eq. 5.19-5.20 convention).
valuesCorrelation estimate on lags with the requested normalization.
coefficientCorrelation coefficient function on the same lags (Eq. 5.16; equals values when normalization='coefficient').
normalization'biased' (), 'unbiased' (, Eq. 11.96) or 'coefficient'.
kind'autocorrelation' or 'cross-correlation'.
fsSample rate, in Hz.
n_samplesRecord length N, in samples.
durationRecord length , in seconds.
CorrelationResult.plot(
ax: Axes | None = None,
*,
language: str = 'en',
**kwargs: Any,
) -> Axes

Plot the correlation estimate against the lag in seconds.

Parameters

NameDescription
languageLabel language, "en" (default) or "es".
CorrelationResult.random_error(
signal_bandwidth: float,
) -> NDArray[np.float64]

Normalized random error of the estimate at each lag.

For bandwidth-limited Gaussian data of bandwidth B and record length T (B&P Eqs. 8.109 and 8.112, identical in form for the cross- and autocorrelation, with the measured coefficient in place of the true value):

so (Eq. 8.111). The large-T approximation behind it assumes and (Section 8.4.1). Lags where the measured coefficient is zero return inf.

Parameters

NameDescription
signal_bandwidthSignal bandwidth B, in Hz.

Returns: Normalized random error per lag (dimensionless).

Raises

ExceptionWhen
ValueErrorIf the bandwidth is not positive.
impulse_response_delay(
ir: NDArray[np.float64] | list[float],
fs: float,
*,
reference: NDArray[np.float64] | list[float] | None = None,
interpolation: _Interpolation = 'parabolic',
upsample: int = 8,
) -> float

Sub-sample delay of an impulse response, in seconds.

Without a reference, the arrival time of the IR itself: the cross-correlation of an IR with an ideal unit impulse is the IR, so its peak magnitude location - refined to sub-sample resolution by band-limited local upsampling plus parabolic interpolation - is the delay relative to . With a reference IR, the delay of ir relative to it, from the peak of their full-record cross-correlation with the same refinement (one-shot transients are not stationary records, so the direct correlator is used rather than the Welch-averaged GCC).

Sub-sample accuracy presumes the IR is band-limited below Nyquist; the synthetic fractional-delay tests pin the achievable accuracy (about 1e-3 samples for a band-limited pulse at the default upsample=8).

Parameters

NameDescription
irImpulse response, 1-D.
fsSample rate, in Hz.
referenceOptional reference IR (same length) the delay is measured against.
interpolation'parabolic' (default) or 'none'.
upsampleInteger local-upsampling factor (default 8).

Returns: Delay in seconds (relative to or to reference).

Raises

ExceptionWhen
ValueErrorIf the inputs or parameters are invalid.
time_delay(
x: NDArray[np.float64] | list[float],
y: NDArray[np.float64] | list[float],
fs: float,
*,
method: _Method = 'gcc',
weighting: _Weighting = 'phat',
window: str = 'hann',
nperseg: int | None = None,
overlap: float = 0.5,
max_delay: float | None = None,
interpolation: _Interpolation = 'parabolic',
upsample: int = 1,
signal_bandwidth: float | None = None,
) -> TimeDelayResult

Time delay of y relative to x (TDE).

Three estimators of the delay in the two-sensor model (B&P Section 5.1.4):

  • 'direct' - the peak of the full-record correlation coefficient function (Eq. 5.21);
  • 'gcc' - the peak of the generalized cross-correlation of Knapp & Carter (1976): the Welch-averaged cross-spectrum (shared core with cross_spectral_density) is weighted by before the inverse transform. Weightings (Table I): 'none' (plain correlator), 'roth' (, suppresses bands where the first sensor is noisy), 'scot' (, prewhitens both channels), 'phat' (: for uncorrelated noises the ideal GCC is a delta at the delay, but errors are accentuated wherever signal power is small), and 'ml' (Hannan-Thomson, ): the maximum-likelihood processor that weights the phase by its coherence-derived reliability and attains the Cramér-Rao bound - provided the coherence estimate is averaged over at least two segments, the discrete form of Knapp & Carter’s existence condition (enforced here). The delay must fit within half a segment; raise nperseg for long delays.
  • 'phase' - the -weighted least-squares slope of the cross-spectrum phase (Eq. 5.101b); accurate for clean, moderate delays where the unwrapped phase is unambiguous, and independent of any peak interpolation (interpolation/upsample do not apply).

The correlation-peak methods refine the sample peak to sub-sample resolution by three-point parabolic interpolation, optionally after band-limited local upsampling (upsample > 1); this presumes the signals are band-limited below Nyquist so the peak is oversampled. With signal_bandwidth given, the result carries the B&P peak-location uncertainty (Eq. 8.129) and its interval (Eq. 8.130).

Parameters

NameDescription
xReference record, 1-D.
yDelayed record, 1-D, same length.
fsSample rate, in Hz.
method'gcc' (default), 'direct' or 'phase'.
weightingGCC weighting (default 'phat'; ignored otherwise).
windowWelch taper for 'gcc'/'phase' (default Hann).
npersegWelch segment length for 'gcc'/'phase' (None picks the shared default).
overlapWelch overlap fraction for 'gcc'/'phase'.
max_delayLargest delay magnitude searched, in seconds.
interpolation'parabolic' (default) or 'none'.
upsampleInteger local-upsampling factor (default 1: off).
signal_bandwidthSignal bandwidth B in Hz for the Eq. 8.129 delay uncertainty (None: no error reported).

Returns: A TimeDelayResult.

Raises

ExceptionWhen
ValueErrorIf the inputs or parameters are invalid.
TimeDelayResult(
delay: float,
delay_samples: float,
method: str,
weighting: str | None,
lags: NDArray[np.float64],
correlation: NDArray[np.float64],
peak_correlation: float,
delay_std: float | None,
delay_interval: tuple[float, float] | None,
signal_bandwidth: float | None,
fs: float,
)

Time-delay estimate between two records.

Attributes

NameDescription
delayEstimated delay of the second record relative to the first, in seconds (positive: y lags x).
delay_samplesThe same delay in (fractional) samples.
method'direct', 'gcc' or 'phase'.
weightingGCC weighting name (None unless method='gcc').
lagsLag axis of correlation, in seconds.
correlationThe correlation function whose peak was located: the correlation coefficient for 'direct', the weighted GCC (normalized to unit peak magnitude) for 'gcc', and the unweighted equivalent for 'phase' (whose estimate comes from Eq. 5.101b, not from this curve).
peak_correlationPlain correlation coefficient at the estimated delay (rounded to the nearest sample) - the quantity entering the B&P error formulas, whatever the method.
delay_stdStandard deviation of the peak-location estimate, (Eq. 8.129), in seconds; None unless signal_bandwidth was given.
delay_intervalApproximate 95 % confidence interval (Eq. 8.130), in seconds; None without a bandwidth.
signal_bandwidthThe bandwidth B used for the error, Hz.
fsSample rate, in Hz.
TimeDelayResult.plot(
ax: Axes | None = None,
*,
language: str = 'en',
**kwargs: Any,
) -> Axes

Plot the correlation function with the estimated delay marked.

Parameters

NameDescription
languageLabel language, "en" (default) or "es".