ph-color 0.1.1

Fixed-point no_std color math for embedded targets: conversion, transfer functions, matrices, gain, and interpolation
Documentation
//! Additive `f32` color triple (`feature = "f32"`).
//!
//! Channels are `0.0..=1.0`, not UQ0.16 integers. This module does not wrap or
//! replace the fixed-point numeric core.

use core::marker::PhantomData;

use crate::color::Color;
use crate::encoding::{Encoded, Encoding, Linear};
use crate::fixed::{Q0_16, Q4_28};
use crate::space::{ColorSpace, Perceptual};

/// Maximum absolute error of the `f32` path versus frozen UQ0.16 W7 goldens,
/// in LSBs after rounding `(channel * 65535)` to nearest (half away from zero).
///
/// Host measurement on matrix, gain, LUT lookup, and lerp CSVs: max 1 LSB
/// (gain: 0). The asserted bound is 4 LSB. The fixed-point path remains
/// bit-identical when this feature is enabled.
///
/// Oklab has no `f32` entry points; that remains out of scope.
pub const F32_MAX_ERR_LSB: u16 = 4;

const U16_MAX_F: f32 = 65535.0;

/// Three `f32` channels tagged with space and encoding.
///
/// Values are linear light or encoded in `0.0..=1.0`. Enabling
/// `feature = "f32"` does not change any fixed-point result.
///
/// Oklab `f32` conversion is out of scope.
pub struct ColorF32<S: ColorSpace, E: Encoding> {
    /// Channel values in `0.0..=1.0` (`0.0` = 0, `1.0` = UQ0.16 full scale).
    pub ch: [f32; 3],
    _pd: PhantomData<fn() -> (S, E)>,
}

impl<S: ColorSpace, E: Encoding> Copy for ColorF32<S, E> {}

impl<S: ColorSpace, E: Encoding> Clone for ColorF32<S, E> {
    fn clone(&self) -> Self {
        *self
    }
}

impl<S: ColorSpace, E: Encoding> PartialEq for ColorF32<S, E> {
    fn eq(&self, other: &Self) -> bool {
        let [a0, a1, a2] = self.ch;
        let [b0, b1, b2] = other.ch;
        a0.to_bits() == b0.to_bits() && a1.to_bits() == b1.to_bits() && a2.to_bits() == b2.to_bits()
    }
}

impl<S: ColorSpace, E: Encoding> Eq for ColorF32<S, E> {}

impl<S: ColorSpace, E: Encoding> core::fmt::Debug for ColorF32<S, E> {
    fn fmt(&self, f: &mut core::fmt::Formatter<'_>) -> core::fmt::Result {
        f.debug_struct("ColorF32").field("ch", &self.ch).finish()
    }
}

impl<S: ColorSpace, E: Encoding> ColorF32<S, E> {
    /// Construct from three `0.0..=1.0` channels.
    #[must_use]
    pub const fn new(ch: [f32; 3]) -> Self {
        Self {
            ch,
            _pd: PhantomData,
        }
    }

    /// Convert to [`Q0_16`] [`Color`] by rounding to nearest and saturating.
    #[must_use]
    pub fn to_color(self) -> Color<S, E> {
        let [c0, c1, c2] = self.ch;
        Color::new([unit_to_q0_16(c0), unit_to_q0_16(c1), unit_to_q0_16(c2)])
    }
}

impl<S: ColorSpace> ColorF32<S, Linear> {
    /// Linear interpolation. `t` is `0.0..=1.0` (`0.0` returns `self`).
    #[must_use]
    pub fn lerp(self, other: Self, t: f32) -> Self {
        let [a0, a1, a2] = self.ch;
        let [b0, b1, b2] = other.ch;
        Self::new([
            lerp_unit(a0, b0, t),
            lerp_unit(a1, b1, t),
            lerp_unit(a2, b2, t),
        ])
    }
}

impl<S: ColorSpace + Perceptual> ColorF32<S, Encoded> {
    /// Interpolate encoded values only when the space is [`Perceptual`].
    #[must_use]
    pub fn lerp(self, other: Self, t: f32) -> Self {
        let [a0, a1, a2] = self.ch;
        let [b0, b1, b2] = other.ch;
        Self::new([
            lerp_unit(a0, b0, t),
            lerp_unit(a1, b1, t),
            lerp_unit(a2, b2, t),
        ])
    }
}

impl<S: ColorSpace, E: Encoding> From<Color<S, E>> for ColorF32<S, E> {
    fn from(color: Color<S, E>) -> Self {
        color.to_f32()
    }
}

impl<S: ColorSpace, E: Encoding> From<ColorF32<S, E>> for Color<S, E> {
    fn from(color: ColorF32<S, E>) -> Self {
        color.to_color()
    }
}

/// Saturate to the unit interval. NaN becomes `0.0`; ±infinity clamp to the
/// nearer endpoint (`-inf` → `0.0`, `+inf` → `1.0`).
///
/// `max`/`min` rather than a comparison chain: IEEE 754 `maxNum`/`minNum`
/// return the non-NaN operand, which is exactly the NaN rule above, and a
/// target with an FPU lowers the pair to two instructions with no branch.
///
/// The trailing `+ 0.0` is load-bearing. `f32::max` may return **either**
/// operand when both are zeros of different sign, so `max(-0.0, 0.0)` is
/// `-0.0` on some targets and `+0.0` on others; the comparison chain this
/// replaced always produced `+0.0`. Adding `+0.0` normalises that one case
/// (`-0.0 + 0.0 == +0.0` under round-to-nearest) and is the exact identity
/// for every other finite input, so the result stays bit-stable across
/// targets. `sat_unit_matches_the_comparison_chain_bit_for_bit` pins it.
#[must_use]
#[allow(clippy::manual_clamp)] // `clamp` returns NaN for NaN; this must return 0.0.
pub(crate) fn sat_unit(x: f32) -> f32 {
    x.max(0.0).min(1.0) + 0.0
}

/// Unit `f32` to [`Q0_16`]: round to nearest (half away from zero), saturate.
///
/// Implemented without `libm` (`f32::round` is not in `core` on this MSRV).
#[must_use]
pub(crate) fn unit_to_q0_16(x: f32) -> Q0_16 {
    if x.is_nan() || x <= 0.0 {
        return Q0_16::ZERO;
    }
    if x >= 1.0 {
        return Q0_16::ONE;
    }
    let scaled = x * U16_MAX_F;
    if scaled >= U16_MAX_F {
        return Q0_16::ONE;
    }
    // `scaled` is in `(0, 65535)`. `as u32` truncates toward zero.
    let trunc = scaled as u32;
    let frac = scaled - (trunc as f32);
    let rounded = if frac >= 0.5 {
        trunc.saturating_add(1)
    } else {
        trunc
    };
    if rounded >= u32::from(u16::MAX) {
        Q0_16::ONE
    } else {
        // Proven bounded: `0..=65534`.
        Q0_16::from_raw(rounded as u16)
    }
}

/// Blend `a` toward `b` with unit `t`, saturating the result.
#[must_use]
pub(crate) fn lerp_unit(a: f32, b: f32, t: f32) -> f32 {
    let a = sat_unit(a);
    let b = sat_unit(b);
    let t = sat_unit(t);
    if t <= 0.0 {
        return a;
    }
    if t >= 1.0 {
        return b;
    }
    sat_unit(a + (b - a) * t)
}

/// One matrix row: decoded [`Q4_28`] × unit channels, saturate to `0.0..=1.0`.
#[must_use]
pub(crate) fn apply_row_f32(coefs: [Q4_28; 3], ch: [f32; 3]) -> f32 {
    let [k0, k1, k2] = coefs;
    let [c0, c1, c2] = ch;
    let c0 = sat_unit(c0);
    let c1 = sat_unit(c1);
    let c2 = sat_unit(c2);
    sat_unit(k0.to_f32() * c0 + k1.to_f32() * c1 + k2.to_f32() * c2)
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::color::Color;
    use crate::encoding::{Encoded, Linear};
    use crate::space::Srgb;

    #[test]
    fn sat_unit_matches_the_comparison_chain_bit_for_bit() {
        // The branchless form replaced `is_nan() || x <= 0.0 ...`; pin the
        // cases where IEEE `maxNum`/`minNum` could plausibly differ.
        fn reference(x: f32) -> f32 {
            if x.is_nan() || x <= 0.0 {
                0.0
            } else if x >= 1.0 {
                1.0
            } else {
                x
            }
        }
        for x in [
            f32::NAN,
            -f32::NAN,
            f32::INFINITY,
            f32::NEG_INFINITY,
            -0.0,
            0.0,
            f32::MIN_POSITIVE,
            -f32::MIN_POSITIVE,
            1.0,
            1.0 - f32::EPSILON,
            1.0 + f32::EPSILON,
            0.5,
            -1.0,
            1e30,
            -1e30,
        ] {
            assert_eq!(
                sat_unit(x).to_bits(),
                reference(x).to_bits(),
                "x={x:?} bits={:08x}",
                x.to_bits()
            );
        }
    }

    #[test]
    fn positive_infinity_saturates_high() {
        assert_eq!(sat_unit(f32::INFINITY), 1.0);
        assert_eq!(sat_unit(f32::NEG_INFINITY), 0.0);
        assert_eq!(sat_unit(f32::NAN), 0.0);
        assert_eq!(unit_to_q0_16(f32::INFINITY), Q0_16::ONE);
        assert_eq!(unit_to_q0_16(f32::NEG_INFINITY), Q0_16::ZERO);
        assert_eq!(unit_to_q0_16(f32::NAN), Q0_16::ZERO);
        let lut = crate::interp::InterpLut::<17>::from_knots(
            Q0_16::array_from_raw([
                0, 4096, 8192, 12288, 16384, 20480, 24576, 28672, 32768, 36864, 40960, 45056,
                49152, 53248, 57344, 61440, 65535,
            ]),
            0,
        );
        assert_eq!(lut.lookup_f32(f32::INFINITY), 1.0);
        assert_eq!(lut.lookup_f32(f32::NEG_INFINITY), 0.0);
        let one = Q4_28::ONE;
        let zero = Q4_28::ZERO;
        let m = crate::matrix::Matrix3::<Srgb, Srgb>::from_q428([
            [one, zero, zero],
            [zero, one, zero],
            [zero, zero, one],
        ]);
        let mixed = ColorF32::<Srgb, Linear>::new([0.5, f32::NAN, 0.5]);
        assert_eq!(m.apply_f32(mixed).ch, [0.5, 0.0, 0.5]);
    }

    #[test]
    fn endpoints_round_trip_to_u16() {
        let z = Color::<Srgb, Linear>::new(Q0_16::array_from_raw([0, 0, 0]));
        let f = Color::<Srgb, Linear>::new(Q0_16::array_from_raw([65535, 65535, 65535]));
        assert_eq!(ColorF32::from(z).to_color(), z);
        assert_eq!(ColorF32::from(f).to_color(), f);
    }

    #[test]
    fn lerp_endpoints() {
        let a = ColorF32::<Srgb, Linear>::new([0.0, 0.0, 0.0]);
        let b = ColorF32::<Srgb, Linear>::new([1.0, 0.5, 0.0]);
        assert_eq!(a.lerp(b, 0.0), a);
        assert_eq!(a.lerp(b, 1.0), b);
        let hi = ColorF32::<Srgb, Linear>::new([1.0, 1.0, 1.0]);
        let lo = ColorF32::<Srgb, Linear>::new([1e-8, 1e-8, 1e-8]);
        assert_eq!(hi.lerp(lo, 1.0), lo);
        let nan = ColorF32::<Srgb, Linear>::new([f32::NAN, f32::NAN, f32::NAN]);
        let one = ColorF32::<Srgb, Linear>::new([1.0, 1.0, 1.0]);
        assert_eq!(nan.lerp(one, 0.5).ch, [0.5, 0.5, 0.5]);
    }

    fn assert_send_sync<T: Send + Sync>() {}

    #[test]
    fn color_f32_is_send_sync() {
        assert_send_sync::<ColorF32<Srgb, Linear>>();
        assert_send_sync::<ColorF32<Srgb, Encoded>>();
    }
}