Spectra and resonances¶
Extract resonances from a ring-down by harmonic inversion, and complete a truncated spectrum when a run is stopped before the fields have fully decayed.
Resonances¶
ResonanceAnalysis(
freq_window: Tuple[float, float],
init_num_freqs: int = 200,
rcond: float = 1e-08,
)
Extract resonances (f, decay, Q, amplitude, phase) from a time signal.
Parameters¶
freq_window : (float, float)
(f_min, f_max) in Hz. The solver is seeded with an even grid of
trial frequencies in this window; it then iteratively optimizes and
prunes them. A narrow window around a few resonances is more accurate
than a broad one. Resonances slightly outside the window may still be
returned.
init_num_freqs : int, default 200
Number of trial frequencies. Larger finds more (and fainter)
resonances at higher cost; the number returned is typically smaller.
rcond : float, default 1e-8
Relative cutoff for the generalized-eigenproblem conditioning: basis
directions whose overlap-matrix eigenvalue magnitude is below
rcond * max are dropped. A real probe signal carries each mode's
mirror pole at -f, outside the seeded window; the kept directions
must be able to model it, or it leaks into the in-window pole. At the
former default of 1e-4 that leak moved a high-Q fit by up to 36 % in
Q, erratically with the record length, and dropped a mode 60 dB
below a neighbour; at 1e-8 the high-Q fit holds to about 1e-6 and
the weak mode is returned to about 1e-3 on float32 probe signals.
Larger values return fewer poles: on a record dominated by noise
(white noise of 1e-2 of the peak or more) a larger rcond returns
fewer spurious in-window poles and, on longer records, less scatter
in Q.
Notes¶
The input signal must be uniformly sampled and should begin after the
sources have turned off -- the source-injection transient is not part of
the ring-down model and degrades the fit. run and
run_time_series gate it with t_start= or simulation= and
warn when the first sample fitted precedes the end of the last source.
Find resonances in a 1-D time series sampled every dt seconds.
Parameters¶
signal : array-like, shape (n_samples,) Real or complex time series (an FDTD field probe is real). Must be uniformly sampled and start after the sources are off. dt : float Sample spacing in seconds.
Returns¶
xarray.Dataset
Indexed by freq (Hz), with data variables decay (1/s),
Q, amplitude, phase (rad) and error (a rough
per-mode reliability estimate). Sorted by freq. The raw set may
include faint or spurious modes -- use select_resonances to
keep the physical ones.
run_time_series(
dataarray: DataArray,
*,
t_start: Optional[float] = None,
simulation=None,
) -> Dataset
Find resonances in a TimeMonitor xarray.DataArray.
dt is derived from the t coordinate, which must be uniformly
spaced -- except for the one short final interval the engine's monitor
cadence produces when n_steps is not a multiple of
interval_steps; that extra final
sample is dropped with a warning.
If the array has a component dimension, the electric components are
summed, falling back to the magnetic ones if no E field is present (the
standard convention -- never E and H mixed); pass a single-component
array to control this exactly.
t_start and simulation gate the fit as in run.
run(
sim_data: Mapping[str, DataArray],
monitors: Union[str, Sequence[str]],
fields: Optional[Sequence[str]] = None,
*,
t_start: Optional[float] = None,
simulation=None,
) -> Dataset
Find resonances in one or more TimeMonitor outputs.
Parameters¶
sim_data : RunResult or mapping name -> DataArray
The run output (anything indexable by monitor name).
monitors : str or sequence of str
Monitor name(s) to read. Multiple monitors are summed sample by
sample into one signal; they must share one time axis. Summing
can also cancel a mode (one antisymmetric between two probes
vanishes from their sum), so check a mode you miss on each probe
alone.
fields : sequence of str, optional
Field components to use (e.g. ["Ez"]). Default: sum the electric
components present, falling back to the magnetic ones -- the standard
convention.
t_start : float, optional
Fit only the samples at t >= t_start (seconds). The ring-down
model has no term for the source pulse, so the fit should start
after the last source ends.
simulation : Simulation, optional
The run's simulation. Without t_start, the fit starts at the
end of its last source (max source_end_time_s, the instant
auto-shutoff measures from). A fit whose first sample
precedes that instant warns; the check also uses the simulation a
RunResult carries when none is passed.
Returns¶
xarray.Dataset
Same schema as run_raw_signal.
select_resonances(
resonances: Dataset,
*,
freq_window: Optional[Tuple[float, float]] = None,
min_amplitude: Optional[float] = None,
min_amplitude_rel: Optional[float] = None,
min_q: Optional[float] = None,
max_error: Optional[float] = None,
require_decay: bool = True,
sort_by: Optional[str] = "amplitude",
) -> Dataset
Filter and rank a raw resonance xarray.Dataset.
The filter-diagonalization output deliberately includes faint and spurious
modes (low amplitude, high error, or non-physical negative decay). This
keeps the physical ones.
Parameters¶
resonances : xarray.Dataset
Output of ResonanceAnalysis.run_raw_signal / run.
freq_window : (float, float), optional
Keep only resonances with f_min <= freq <= f_max.
min_amplitude, min_q, max_error : float, optional
Lower bound on amplitude / Q, upper bound on error.
min_amplitude is in RAW FIELD UNITS, a unit PointDipole
ring-down probes at ~1e-6 V/m scale, so an absolute threshold that
looks small (1e-3) can silently reject every physical mode. Unless the
signal was pre-normalized, prefer min_amplitude_rel.
min_amplitude_rel : float, optional
Lower amplitude bound as a FRACTION of the strongest mode in the
(freq_window-restricted) set, scale-free, so it works on raw probe
signals without knowing the field magnitude. E.g. 1e-3 keeps modes
within 60 dB of the dominant one. Combines with min_amplitude
(both must pass).
require_decay : bool, default True
Drop modes with decay <= 0 (growing in time). NOTE: in a CLOSED
lossless scene (PEC box, undamped cavity) the physical mode is
undamped, and its fitted decay can land at a tiny NEGATIVE value, the default then drops the very mode you are after. Pass
require_decay=False there and treat |decay| * T_window << 1
as "undamped within the window's resolution".
sort_by : str or None, default "amplitude"
Data variable to sort by, descending (e.g. "amplitude" or
"Q"). None leaves the freq ordering. The strongest mode
comes first by default: a spurious pole fitted to noise can carry a
higher Q than any physical mode, so ranking by Q over a broad window
can put it first.
Returns¶
xarray.Dataset The filtered (and optionally ranked) subset.
Spectral completion¶
SpectrumCompleter(
freq_window: Tuple[float, float],
*,
fit_start: Optional[float] = None,
holdout_fraction: float = 0.25,
residual_tol: float = 0.02,
min_amplitude_rel: float = 0.001,
min_decay_resolved: float = 0.1,
strict: bool = False,
resonance_analysis: Optional[ResonanceAnalysis] = None,
)
Complete a truncated resonator spectrum from its fitted ringdown.
Parameters¶
freq_window : (float, float)
Frequency band (Hz) searched for resonant poles, passed to
ResonanceAnalysis. Should generously cover the band of
freqs you will complete.
fit_start : float, optional
Absolute time (seconds, on the record's own time axis) where the
pure-ringdown model becomes valid, i.e. after every source has
switched off. Defaults to the midpoint of the record; setting it
explicitly to just after your source's end is always better.
holdout_fraction : float, default 0.25
Fraction of the fit window [fit_start, T] held out for the
acceptance gate. The model is fitted on the LAST
1 - holdout_fraction of the window and extrapolated backward
over the first holdout_fraction, which the fit never saw.
residual_tol : float, default 0.02
Acceptance threshold on the held-out window: relative RMS residual
between the backward-extrapolated model and the recorded signal.
min_amplitude_rel : float, default 1e-3
Poles whose FDM amplitude is below this fraction of the strongest
pole's are dropped before the amplitude re-fit (spurious-mode guard).
min_decay_resolved : float, default 0.1
Reject unless every kept pole decays by at least this many decay
constants across the fit segment (alpha * span >= this). The
tail is an extrapolation of the fitted decay RATE; a record over
which a mode's amplitude changes by under ~10% has not measured that
rate, it has assumed it, the limiting case being a pure sinusoid,
whose fitted alpha is numerical noise and whose "tail" would be
a near-divergent spike the holdout gate cannot catch (the model fits
the held-out window perfectly).
strict : bool, default False
If True, a rejected gate raises CompletionRejected instead
of warning and returning the truncated spectrum.
resonance_analysis : ResonanceAnalysis, optional
Pre-configured pole finder; by default one is built from
freq_window with its default basis size and rcond=1e-4. The
completer keeps that former default: at ResonanceAnalysis's
current 1e-8, heavily damped spurious poles clear
min_amplitude_rel on noisy records, and extrapolated back over
the holdout they reject records the former cutoff accepted.
complete(
sim_data: Mapping[str, DataArray],
monitor: str,
freqs: Sequence[float],
fields: Optional[Sequence[str]] = None,
) -> Dataset
Complete the spectrum of a TimeMonitor output.
Parameters mirror ResonanceAnalysis.run; the monitor's own time
coordinate supplies both the sample times and the spacing, so a
decimated probe (interval_steps > 1) is handled consistently.
complete_raw_signal(
signal: ndarray,
dt: float,
freqs: Sequence[float],
*,
t0: float = 0.0,
) -> Dataset
Complete the spectrum of a raw, uniformly sampled REAL signal.
Parameters¶
signal : array
Real samples u(t_n) at t_n = t0 + n*dt.
dt : float
Sample spacing (seconds). For an engine probe this is the
monitor's recording interval, not necessarily the solver step.
freqs : sequence of float
Frequencies (Hz) at which to evaluate the spectrum.
t0 : float, default 0.0
Absolute time of the first sample.
Returns¶
xarray.Dataset
spectrum (complex, completed), truncated (complex, the
recorded-samples-only DFT), tail (complex, the analytic
remainder actually added, zero when rejected), each over
freq; per-mode diagnostics mode_freq / mode_decay /
mode_q over mode; and attrs accepted (1/0),
holdout_residual, residual_tol, fit_start,
holdout_end, record_end, reject_reason ("" when
accepted).
Bases: RuntimeError
The held-out-window gate refused the ringdown model (strict=True).