Skip to content

Spectral Transforms

spectrakit.transform.transform_kubelka_munk

transform_kubelka_munk(reflectance: ndarray) -> np.ndarray

Convert diffuse reflectance to Kubelka-Munk units.

Applies the transformation: K/S = (1 - R)^2 / (2R)

where R is the diffuse reflectance (0 to 1 scale). This linearizes the relationship between concentration and spectral response for diffuse reflectance measurements.

Parameters:

Name Type Description Default
reflectance ndarray

Diffuse reflectance values in [0, 1], shape (W,) or (N, W).

required

Returns:

Type Description
ndarray

Kubelka-Munk values (K/S), same shape as input.

Raises:

Type Description
SpectrumShapeError

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

EmptySpectrumError

If input has zero elements.

Source code in src/spectrakit/transform/kubelka_munk.py
def transform_kubelka_munk(reflectance: np.ndarray) -> np.ndarray:
    """Convert diffuse reflectance to Kubelka-Munk units.

    Applies the transformation: K/S = (1 - R)^2 / (2R)

    where R is the diffuse reflectance (0 to 1 scale). This
    linearizes the relationship between concentration and spectral
    response for diffuse reflectance measurements.

    Args:
        reflectance: Diffuse reflectance values in [0, 1], shape
            ``(W,)`` or ``(N, W)``.

    Returns:
        Kubelka-Munk values (K/S), same shape as input.

    Raises:
        SpectrumShapeError: If input is not 1-D or 2-D.
        EmptySpectrumError: If input has zero elements.
    """
    reflectance = ensure_float64(reflectance)
    validate_1d_or_2d(reflectance)
    warn_if_not_finite(reflectance)

    if np.any(reflectance < 0) or np.any(reflectance > 1):
        warnings.warn(
            "Reflectance values outside [0, 1] detected. "
            "Ensure you are passing reflectance, not absorbance.",
            stacklevel=2,
        )

    # Clamp reflectance to avoid division by zero
    r = np.clip(reflectance, EPSILON, 1.0 - EPSILON)

    return (1.0 - r) ** 2 / (2.0 * r)

spectrakit.transform.transform_atr_correction

transform_atr_correction(
    intensities: ndarray,
    wavenumbers: ndarray,
    n_crystal: float = DEFAULT_N_CRYSTAL,
    n_sample: float = DEFAULT_N_SAMPLE,
    angle: float = DEFAULT_ANGLE,
) -> np.ndarray

Apply ATR path-length correction to infrared spectra.

Corrects for the wavenumber-dependent penetration depth in ATR measurements. The effective path length in ATR varies with wavenumber, making peaks at lower wavenumbers appear stronger. This correction normalizes for that effect.

Parameters:

Name Type Description Default
intensities ndarray

ATR spectral intensities, shape (W,) or (N, W).

required
wavenumbers ndarray

Wavenumber axis in cm^-1, shape (W,).

required
n_crystal float

Refractive index of the ATR crystal. Common values: diamond = 2.4, ZnSe = 2.4, Ge = 4.0.

DEFAULT_N_CRYSTAL
n_sample float

Refractive index of the sample. Typical organic samples: 1.4-1.6.

DEFAULT_N_SAMPLE
angle float

Angle of incidence in degrees.

DEFAULT_ANGLE

Returns:

Type Description
ndarray

ATR-corrected 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 physics parameters are invalid (non-positive refractive indices, angle outside (0, 90), or angle below the critical angle for the given crystal/sample pair).

Source code in src/spectrakit/transform/atr_correction.py
def transform_atr_correction(
    intensities: np.ndarray,
    wavenumbers: np.ndarray,
    n_crystal: float = DEFAULT_N_CRYSTAL,
    n_sample: float = DEFAULT_N_SAMPLE,
    angle: float = DEFAULT_ANGLE,
) -> np.ndarray:
    """Apply ATR path-length correction to infrared spectra.

    Corrects for the wavenumber-dependent penetration depth in ATR
    measurements. The effective path length in ATR varies with
    wavenumber, making peaks at lower wavenumbers appear stronger.
    This correction normalizes for that effect.

    Args:
        intensities: ATR spectral intensities, shape ``(W,)`` or ``(N, W)``.
        wavenumbers: Wavenumber axis in cm^-1, shape ``(W,)``.
        n_crystal: Refractive index of the ATR crystal. Common values:
            diamond = 2.4, ZnSe = 2.4, Ge = 4.0.
        n_sample: Refractive index of the sample. Typical organic
            samples: 1.4-1.6.
        angle: Angle of incidence in degrees.

    Returns:
        ATR-corrected intensities, same shape as input.

    Raises:
        SpectrumShapeError: If input is not 1-D or 2-D.
        EmptySpectrumError: If input has zero elements.
        ValueError: If physics parameters are invalid (non-positive
            refractive indices, angle outside (0, 90), or angle below
            the critical angle for the given crystal/sample pair).
    """
    intensities = ensure_float64(intensities)
    wavenumbers = ensure_float64(wavenumbers)
    validate_1d_or_2d(intensities)
    warn_if_not_finite(intensities)

    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}"
        )

    # Validate physics parameters
    if n_crystal <= 0:
        raise ValueError(f"n_crystal must be positive, got {n_crystal}")
    if n_sample <= 0:
        raise ValueError(f"n_sample must be positive, got {n_sample}")
    if not 0 < angle < 90:
        raise ValueError(f"angle must be in (0, 90) degrees, got {angle}")

    n_ratio = n_sample / n_crystal
    if n_ratio >= 1.0:
        raise ValueError(
            f"n_sample ({n_sample}) must be less than n_crystal ({n_crystal}) "
            "for total internal reflection"
        )

    # Compute penetration depth factor: dp ∝ 1 / (ν * sqrt(sin²θ - (n2/n1)²))
    theta = np.radians(angle)
    sin2_theta = np.sin(theta) ** 2

    discriminant = sin2_theta - n_ratio**2
    if discriminant <= 0:
        raise ValueError(
            f"Angle {angle}° is below the critical angle for n_crystal={n_crystal}, "
            f"n_sample={n_sample}. No total internal reflection occurs."
        )

    # Penetration depth: dp = λ / (2π * n1 * sqrt(sin²θ - (n2/n1)²))
    # Since λ = 1/ν (in cm⁻¹), dp ∝ 1 / (ν * sqrt(discriminant))
    # Correction: multiply absorbance by ν * sqrt(discriminant) to remove
    # the path-length dependence, then normalize so the correction factor
    # at the maximum wavenumber equals 1.
    dp_inv = wavenumbers * np.sqrt(discriminant)
    correction = dp_inv / np.max(dp_inv)

    if intensities.ndim == 1:
        return intensities * correction  # type: ignore[no-any-return]
    return intensities * correction[np.newaxis, :]  # type: ignore[no-any-return]