Skip to content

Smoothing

spectrakit.smooth.smooth_savgol

smooth_savgol(
    intensities: ndarray,
    window_length: int = DEFAULT_WINDOW_LENGTH,
    polyorder: int = DEFAULT_POLYORDER,
) -> np.ndarray

Apply Savitzky-Golay smoothing filter.

Fits successive sub-sets of adjacent data points with a low-degree polynomial by the method of linear least squares.

Parameters:

Name Type Description Default
intensities ndarray

Spectral intensities, shape (W,) or (N, W).

required
window_length int

Length of the filter window (must be odd and greater than polyorder).

DEFAULT_WINDOW_LENGTH
polyorder int

Order of the polynomial used to fit the samples.

DEFAULT_POLYORDER

Returns:

Type Description
ndarray

Smoothed intensities, same shape as input.

Raises:

Type Description
SpectrumShapeError

If input is not 1-D or 2-D.

EmptySpectrumError

If input has zero elements.

ValueError

If window_length or polyorder are invalid.

Source code in src/spectrakit/smooth/savgol.py
def smooth_savgol(
    intensities: np.ndarray,
    window_length: int = DEFAULT_WINDOW_LENGTH,
    polyorder: int = DEFAULT_POLYORDER,
) -> np.ndarray:
    """Apply Savitzky-Golay smoothing filter.

    Fits successive sub-sets of adjacent data points with a low-degree
    polynomial by the method of linear least squares.

    Args:
        intensities: Spectral intensities, shape ``(W,)`` or ``(N, W)``.
        window_length: Length of the filter window (must be odd and
            greater than ``polyorder``).
        polyorder: Order of the polynomial used to fit the samples.

    Returns:
        Smoothed intensities, same shape as input.

    Raises:
        SpectrumShapeError: If input is not 1-D or 2-D.
        EmptySpectrumError: If input has zero elements.
        ValueError: If ``window_length`` or ``polyorder`` are invalid.
    """
    if window_length < 1:
        raise ValueError(f"window_length must be >= 1, got {window_length}")
    if window_length % 2 == 0:
        raise ValueError(f"window_length must be odd, got {window_length}")
    if polyorder < 0:
        raise ValueError(f"polyorder must be >= 0, got {polyorder}")
    if polyorder >= window_length:
        raise ValueError(
            f"polyorder ({polyorder}) must be less than window_length ({window_length})"
        )

    intensities = ensure_float64(intensities)
    validate_1d_or_2d(intensities)
    warn_if_not_finite(intensities)

    return apply_along_spectra(
        _smooth_savgol_1d,
        intensities,
        window_length=window_length,
        polyorder=polyorder,
    )

spectrakit.smooth.smooth_whittaker

smooth_whittaker(
    intensities: ndarray,
    lam: float = DEFAULT_LAMBDA,
    differences: int = DEFAULT_DIFFERENCES,
    wavenumbers: ndarray | None = None,
) -> np.ndarray

Apply Whittaker smoother (penalized least squares).

Minimizes the sum of squared residuals plus a penalty on the roughness (measured by finite differences) of the fitted curve.

When wavenumbers are provided, the penalty matrix accounts for non-uniform spacing between spectral points. This is important for spectra measured on instruments with variable point spacing, as the standard uniform-spacing penalty would over- or under-smooth regions with different point densities.

Parameters:

Name Type Description Default
intensities ndarray

Spectral intensities, shape (W,) or (N, W).

required
lam float

Smoothness parameter (lambda). Larger values produce smoother results. Typical range: 1e2 to 1e8.

DEFAULT_LAMBDA
differences int

Order of the difference penalty. 2 penalizes curvature (default), 1 penalizes slope.

DEFAULT_DIFFERENCES
wavenumbers ndarray | None

Wavenumber axis, shape (W,). If provided, builds a non-uniform finite-difference penalty matrix that accounts for the actual point spacing.

None

Returns:

Type Description
ndarray

Smoothed intensities, same shape as input.

Raises:

Type Description
SpectrumShapeError

If input is not 1-D or 2-D.

EmptySpectrumError

If input has zero elements.

ValueError

If lam or differences are invalid, or if wavenumbers length does not match spectral width.

Source code in src/spectrakit/smooth/whittaker.py
def smooth_whittaker(
    intensities: np.ndarray,
    lam: float = DEFAULT_LAMBDA,
    differences: int = DEFAULT_DIFFERENCES,
    wavenumbers: np.ndarray | None = None,
) -> np.ndarray:
    """Apply Whittaker smoother (penalized least squares).

    Minimizes the sum of squared residuals plus a penalty on the
    roughness (measured by finite differences) of the fitted curve.

    When *wavenumbers* are provided, the penalty matrix accounts for
    non-uniform spacing between spectral points. This is important for
    spectra measured on instruments with variable point spacing, as the
    standard uniform-spacing penalty would over- or under-smooth
    regions with different point densities.

    Args:
        intensities: Spectral intensities, shape ``(W,)`` or ``(N, W)``.
        lam: Smoothness parameter (lambda). Larger values produce
            smoother results. Typical range: 1e2 to 1e8.
        differences: Order of the difference penalty. 2 penalizes
            curvature (default), 1 penalizes slope.
        wavenumbers: Wavenumber axis, shape ``(W,)``. If provided,
            builds a non-uniform finite-difference penalty matrix that
            accounts for the actual point spacing.

    Returns:
        Smoothed intensities, same shape as input.

    Raises:
        SpectrumShapeError: If input is not 1-D or 2-D.
        EmptySpectrumError: If input has zero elements.
        ValueError: If ``lam`` or ``differences`` are invalid, or if
            ``wavenumbers`` length does not match spectral width.
    """
    if lam <= 0:
        raise ValueError(f"lam (smoothness) must be positive, got {lam}")
    if differences < 1:
        raise ValueError(f"differences must be >= 1, got {differences}")

    intensities = ensure_float64(intensities)
    validate_1d_or_2d(intensities)
    warn_if_not_finite(intensities)

    if wavenumbers is not None:
        wavenumbers = ensure_float64(wavenumbers)
        expected_w = intensities.shape[-1]
        if wavenumbers.shape[0] != expected_w:
            raise ValueError(
                f"wavenumbers length {wavenumbers.shape[0]} does not match "
                f"intensities spectral width {expected_w}"
            )

    # Pre-compute the penalty matrix once for the entire batch.
    n = intensities.shape[-1]

    if wavenumbers is not None and differences <= 2:
        if differences == 1:
            D = _build_nonuniform_diff1(wavenumbers)
        else:
            D = _build_nonuniform_diff2(wavenumbers)
    else:
        D = sparse.eye(n, format="csc")
        for _ in range(differences):
            D = D[1:] - D[:-1]

    penalty_z = sparse.eye(n, format="csc") + lam * D.T @ D

    return apply_along_spectra(
        _smooth_whittaker_1d,
        intensities,
        penalty_z=penalty_z,
    )