def poisson_window(window_frames: int, *, decay_frames: float | None = None) -> np.ndarray:
    """The exponentially decreasing peak-shaped window the published analysis uses.

    Parameters
    ----------
    window_frames : int
        Window length in frames [-].
    decay_frames : float, optional
        Length over which the weight falls by ``1/e`` [-].

    Returns
    -------
    np.ndarray
        Window weights, peak in the middle [-].

    Convention
    ----------
    Periodic rather than symmetric, which is what a short-time transform wants
    and what :func:`scipy.signal.spectrogram` would build from a window name.

    Defaults
    --------
    ``decay_frames = window_frames / 8`` is an assumed value.  The publications
    name the window shape but not its rate, so this is a VAFT choice, recorded
    rather than presented as published: it leaves the window edges at about 0.02
    of the peak, tapering hard enough to suppress the transform's edge artefacts
    without narrowing the effective window far below the stated 1 ms.

    Applicability
    -------------
    Machine-independent.

    Limitations
    -----------
    Because the rate is a choice, band powers from this window are not
    numerically comparable with published values; the structure they show is.

    Provenance
    ----------
    .. [1] [Jung2022]_ section 4.1: "an exponential decreasing peak-shape time
       window of 1 ms (50 images)".
    .. [2] :func:`scipy.signal.windows.exponential`.
    """
    from scipy.signal import windows

    length = int(window_frames)
    if length < 2:
        raise ValueError(f"window_frames must be at least 2; got {window_frames!r}")
    tau = float(decay_frames) if decay_frames is not None else length / 8.0
    if tau <= 0.0:
        raise ValueError(f"decay_frames must be positive; got {decay_frames!r}")
    return windows.exponential(length, tau=tau, sym=False)

def subtract_temporal_background(frames, *, window_frames: int = BACKGROUND_FRAMES_50KFPS) -> np.ndarray:
    """Remove each pixel's slowly varying emission by subtracting its local temporal mean.

    Parameters
    ----------
    frames : array_like
        Time-ordered image cube, frame axis first [any].
    window_frames : int, optional
        Width of the local mean, in frames [-].

    Returns
    -------
    np.ndarray
        ``frames`` minus its local temporal mean, same shape, float [any].

    Processing steps
    ----------------
    1. Take the centred rolling mean of every pixel over ``window_frames``.
    2. Subtract it from the frame at the centre of that window.

    Convention
    ----------
    The window is centred on the analysed frame and spans
    ``[i - w//2, i + (w-1)//2]``, clipped to the record, so the first and last
    frames average over a shrinking window rather than over zeros and the output
    keeps the input's length.  For an odd window -- both published presets are
    odd -- that span is symmetric; an even window takes one more frame from
    before than after.

    A window wider than the record is allowed and clips the same way, so every
    frame is measured against the whole record.  Refusing it would make the
    result depend on how many frames a caller happened to pass rather than on
    the window it asked for: a fourteen-frame acquisition and a fourteen-frame
    slice of a long one would come out differently.

    Frame values and their order are otherwise untouched: this subtracts, it does
    not filter, reorder or resample.

    Defaults
    --------
    ``BACKGROUND_FRAMES_50KFPS = 15`` and ``BACKGROUND_FRAMES_2500FPS = 3`` are
    literature values.  The 50 kFrames/s preset is the analysed frame plus seven
    before and seven after, which is what reproduces the published filter
    behaviour: measured at 50 kFrames/s this window passes 0.89 of a 3 kHz
    component and 1.11 of a 6 kHz one, matching the statement that it removes
    components below about 3 kHz while enhancing the several-kHz structures.  A
    15-frame window is therefore not interchangeable with a 7-frame one, which
    passes only 0.26 at 3 kHz and 0.81 at 6 kHz -- it would attenuate the MHD
    band this analysis exists to see.

    Applicability
    -------------
    Machine-independent.

    Limitations
    -----------
    The subtraction is zero-phase but the shrinking edge windows bias the first
    and last ``window_frames / 2`` frames toward the local level; a burst there is
    partly removed along with the background.  A single very bright frame pulls
    the mean toward itself and appears, inverted, in its neighbours.

    Provenance
    ----------
    .. [1] [Jung2022]_ section 3.1: "subtracts the average value of each pixel in
       a series of seven images before and after the analysed image ... removes
       slowly varying components (less than about 3 kHz)".
    .. [2] [JungThesis]_, the background-removal figure.
    """
    cube = _as_cube(frames)
    window = int(window_frames)
    if window < 2:
        raise ValueError(
            "window_frames must be at least 2; a one-frame mean is the frame "
            f"itself, so subtracting it yields exactly zero. Got {window_frames!r}"
        )
    return cube - _rolling_mean(cube, window)

def summed_region_signal(frames, *, region: Sequence[int] | None = None) -> np.ndarray:
    """Sum the intensity inside one image region, giving a single time series.

    Parameters
    ----------
    frames : array_like
        Time-ordered image cube, frame axis first [any].
    region : sequence of int, optional
        Pixel box as ``(row_start, row_stop, column_start, column_stop)``, the
        stops exclusive; ``None`` sums the whole frame [-].

    Returns
    -------
    np.ndarray
        One value per frame [any].

    Convention
    ----------
    Bounds are half-open and in row-then-column order, the order an image array
    is indexed with -- which is the opposite of the column-then-row order
    :func:`vaft.process.camera_geometry.project_points` returns pixels in.  A
    region is stated explicitly rather than inferred, because which part of the
    image faces the low-field side depends on the camera pose.

    Applicability
    -------------
    Machine-independent.

    Limitations
    -----------
    A plain sum, not a mean: two regions of different size are not comparable in
    magnitude. Nothing checks that the region holds plasma rather than wall.

    Provenance
    ----------
    .. [1] [Jung2022]_ figure 3, the summed camera pixel intensity viewing the
       low-field side whose spectrogram is compared with the magnetic probe's.
    """
    cube = _as_cube(frames)
    if region is None:
        return cube.reshape(cube.shape[0], -1).sum(axis=1)
    if len(region) != 4:
        raise ValueError(
            "region must be (row_start, row_stop, column_start, column_stop); "
            f"got {region!r}"
        )
    r0, r1, c0, c1 = (int(value) for value in region)
    if r0 >= r1 or c0 >= c1:
        raise ValueError(f"region bounds must be increasing; got {region!r}")
    window = cube[:, r0:r1, c0:c1]
    if window.size == 0:
        raise ValueError(
            f"region {region!r} selects no pixels of a {cube.shape[1:]} frame"
        )
    return window.reshape(window.shape[0], -1).sum(axis=1)

def track_reference_frequency(
    spectrogram,
    *,
    search_range: Sequence[float] = REFERENCE_SEARCH_RANGE_HZ,
) -> np.ndarray:
    """Follow the dominant MHD frequency a magnetic probe reports, window by window.

    Parameters
    ----------
    spectrogram : Any
        Any object carrying ``time``, ``frequency`` and a ``(frequency, time)``
        ``magnitude``, such as a magnetics spectrogram result [-].
    search_range : sequence of float, optional
        Lowest and highest frequency the peak may be found at [Hz].

    Returns
    -------
    np.ndarray
        One frequency per window, ``nan`` where no bin lies in range [Hz].

    Convention
    ----------
    The frequency is tracked, not assumed.  The published filtering is "centred
    at the dominant MHD mode frequency measured by the magnetic probe", which
    moves during a shot; a fixed 6 kHz is the value of one example, not the
    method.

    Defaults
    --------
    ``REFERENCE_SEARCH_RANGE_HZ = (3000, 15000)`` is a legacy compatibility
    value.  The legacy code searched the first fifteen bins of a magnetics
    transform whose resolution was 1 kHz and rejected anything below the fourth,
    which is this range expressed in hertz rather than in bins, so it means the
    same thing at any window length.  The lower edge matters: without it the
    peak sits on the equilibrium field's own slow variation instead of the mode.

    Applicability
    -------------
    Machine-independent.

    Limitations
    -----------
    A single argmax per window, with no continuity requirement between windows,
    so it can jump between two comparable modes from one window to the next. It
    reports the strongest component in range whether or not one is really there;
    a window with no mode still returns its loudest bin.
    :func:`vaft.process.fluctuation.track_dominant_frequency`, which this calls,
    offers a floor and a continuity limit for callers that want neither.

    Provenance
    ----------
    .. [1] [CameraFFT]_: ``[~, idx_max] = max(P1_mir(1:15, iFT))`` with windows
       whose ``idx_max`` falls below four discarded.
    .. [2] [Jung2022]_ section 4.1, which specifies the magnetic probe as the
       source of the centre frequency.
    """
    from vaft.process.fluctuation import track_dominant_frequency

    # The shared ridge tracker with no floor and no continuity requirement is
    # exactly the legacy per-window argmax; only an empty window reports NaN.
    return track_dominant_frequency(
        spectrogram, search_range=search_range, floor_ratio=0.0
    ).frequency