colr 0.1.1

//! Chromatic adaptation transforms XYZ tristimulus values from one
//! illuminant to another.
//!
//! Each adaptation method defines a cone response matrix M. The adaptation
//! matrix for a given illuminant pair is computed at compile time. Callers
//! pass WHITE_POINT_XYZ values from the Illuminant trait directly to [adapt]
//! to ensure numerical consistency with specification-defined white points.
//!
//! Bradford is mandated by ICC, ACES, and CSS Color Level 4 and is the
//! correct default for all standard RGB space work.

use crate::math::Mat3;

/// A chromatic adaptation method defined by its cone response matrix.
///
/// Only M is required. The inverse is derived by [adapt] at compile time
/// via Mat3::invert, ensuring consistency and eliminating the possibility
/// of a mismatched inverse.
pub trait ChromaticAdaptation: 'static {
    /// Forward cone response matrix, row-major.
    const M: [[f32; 3]; 3];
}

/// Compute the adaptation matrix: XYZ under src_xyz to XYZ under dst_xyz.
///
/// Takes XYZ tristimulus values directly rather than xy chromaticities.
/// Use Illuminant::WHITE_POINT_XYZ to avoid the floating-point rounding
/// that occurs when deriving XYZ from xy, which can cause slight drift
/// from specification-defined adaptation matrices.
///
/// For convenience when only xy chromaticities are available, use [adapt_xy].
pub const fn adapt<A: ChromaticAdaptation>(src_xyz: [f32; 3], dst_xyz: [f32; 3]) -> Mat3 {
    let m = A::M;
    let m_inv = mat3_to_rows(Mat3::invert(&rows_to_mat3(m)));

    let cs = mul_vec(m, src_xyz);
    let cd = mul_vec(m, dst_xyz);

    let s = [cd[0] / cs[0], cd[1] / cs[1], cd[2] / cs[2]];

    let sm = [
        [s[0] * m[0][0], s[0] * m[0][1], s[0] * m[0][2]],
        [s[1] * m[1][0], s[1] * m[1][1], s[1] * m[1][2]],
        [s[2] * m[2][0], s[2] * m[2][1], s[2] * m[2][2]],
    ];

    let r = mul_mat(m_inv, sm);

    Mat3 {
        col0: [r[0][0], r[1][0], r[2][0], 0.0],
        col1: [r[0][1], r[1][1], r[2][1], 0.0],
        col2: [r[0][2], r[1][2], r[2][2], 0.0],
    }
}

/// Compute the adaptation matrix from xy chromaticity coordinates.
///
/// Converts xy to XYZ1 and delegates to [adapt]. Use this only when
/// exact XYZ values are unavailable. For standard illuminants, prefer
/// [adapt] with Illuminant::WHITE_POINT_XYZ.
pub const fn adapt_xy<A: ChromaticAdaptation>(src_xy: [f32; 2], dst_xy: [f32; 2]) -> Mat3 {
    adapt::<A>(xy_to_xyz1(src_xy), xy_to_xyz1(dst_xy))
}

const fn xy_to_xyz1(xy: [f32; 2]) -> [f32; 3] {
    let [x, y] = xy;
    [x / y, 1.0, (1.0 - x - y) / y]
}

const fn mul_vec(m: [[f32; 3]; 3], v: [f32; 3]) -> [f32; 3] {
    [
        m[0][0] * v[0] + m[0][1] * v[1] + m[0][2] * v[2],
        m[1][0] * v[0] + m[1][1] * v[1] + m[1][2] * v[2],
        m[2][0] * v[0] + m[2][1] * v[1] + m[2][2] * v[2],
    ]
}

const fn mul_mat(a: [[f32; 3]; 3], b: [[f32; 3]; 3]) -> [[f32; 3]; 3] {
    [
        [
            a[0][0] * b[0][0] + a[0][1] * b[1][0] + a[0][2] * b[2][0],
            a[0][0] * b[0][1] + a[0][1] * b[1][1] + a[0][2] * b[2][1],
            a[0][0] * b[0][2] + a[0][1] * b[1][2] + a[0][2] * b[2][2],
        ],
        [
            a[1][0] * b[0][0] + a[1][1] * b[1][0] + a[1][2] * b[2][0],
            a[1][0] * b[0][1] + a[1][1] * b[1][1] + a[1][2] * b[2][1],
            a[1][0] * b[0][2] + a[1][1] * b[1][2] + a[1][2] * b[2][2],
        ],
        [
            a[2][0] * b[0][0] + a[2][1] * b[1][0] + a[2][2] * b[2][0],
            a[2][0] * b[0][1] + a[2][1] * b[1][1] + a[2][2] * b[2][1],
            a[2][0] * b[0][2] + a[2][1] * b[1][2] + a[2][2] * b[2][2],
        ],
    ]
}

const fn rows_to_mat3(m: [[f32; 3]; 3]) -> Mat3 {
    Mat3 {
        col0: [m[0][0], m[1][0], m[2][0], 0.0],
        col1: [m[0][1], m[1][1], m[2][1], 0.0],
        col2: [m[0][2], m[1][2], m[2][2], 0.0],
    }
}

const fn mat3_to_rows(m: Mat3) -> [[f32; 3]; 3] {
    [[m.col0[0], m.col1[0], m.col2[0]], [m.col0[1], m.col1[1], m.col2[1]], [
        m.col0[2], m.col1[2], m.col2[2],
    ]]
}

/// Bradford chromatic adaptation.
///
/// A sharpened von Kries model using a non-diagonal cone response matrix.
/// More accurate than von Kries for large illuminant changes. Mandated by
/// ICC, ACES, and CSS Color Level 4.
///
/// Reference: ICC.1:2022, Annex E.
#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
pub struct Bradford;

impl ChromaticAdaptation for Bradford {
    const M: [[f32; 3]; 3] = [[0.8951, 0.2664, -0.1614], [-0.7502, 1.7135, 0.0367], [
        0.0389, -0.0685, 1.0296,
    ]];
}

/// CAT02 chromatic adaptation.
///
/// Used in CIECAM02 and ICC v4 appearance models. Provides better hue
/// constancy than Bradford in some gamut mapping contexts.
///
/// Reference: CIE 159:2004.
#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
pub struct Cat02;

impl ChromaticAdaptation for Cat02 {
    const M: [[f32; 3]; 3] = [[0.7328, 0.4296, -0.1624], [-0.7036, 1.6975, 0.0061], [
        0.0030, 0.0136, 0.9834,
    ]];
}

/// CAT16 chromatic adaptation.
///
/// Used in CAM16. Successor to CAT02 with improved accuracy for large
/// chromatic differences and reduced numerical instability.
///
/// Reference: Li et al., "Comprehensive colour appearance model (CAM16)",
/// Color Research and Application, 2017.
#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
pub struct Cat16;

impl ChromaticAdaptation for Cat16 {
    const M: [[f32; 3]; 3] = [[0.401288, 0.650173, -0.051461], [-0.250268, 1.204414, 0.045854], [
        -0.002079, 0.048952, 0.953127,
    ]];
}

/// Von Kries chromatic adaptation.
///
/// Diagonal-only scaling of cone channels. Less accurate than Bradford
/// or CAT02 for large illuminant changes. Included for completeness and
/// legacy compatibility.
///
/// Reference: Hunt, "The Reproduction of Colour", 6th edition, 2004.
#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
pub struct VonKries;

impl ChromaticAdaptation for VonKries {
    const M: [[f32; 3]; 3] = [[0.40024, 0.7076, -0.08081], [-0.2263, 1.16532, 0.0457], [
        0.0, 0.0, 0.91822,
    ]];
}