colr 0.2.0

A general purpose, extensible color type unifying color models and their operations at the type level.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135

//! Perceptual color difference (delta E) metrics.

use crate::Color;
use crate::illuminant::Illuminant;
use crate::model::{Lab, Oklab};

#[inline(always)]
fn euclidean(l1: f32, a1: f32, b1: f32, l2: f32, a2: f32, b2: f32) -> f32 {
    let dl = l1 - l2;
    let da = a1 - a2;
    let db = b1 - b2;
    (dl * dl + da * da + db * db).sqrt()
}

// 25^7 used in the CIEDE2000 chroma weighting factor.
const CIEDE2000_C25_7: f32 = 6103515625.0;

// Angle constants from CIE 142-2001, in radians.
const DEG_6: f32 = core::f32::consts::PI / 30.0;
const DEG_25: f32 = 25.0 * core::f32::consts::PI / 180.0;
const DEG_60: f32 = core::f32::consts::FRAC_PI_3;
const DEG_63: f32 = 63.0 * core::f32::consts::PI / 180.0;
const DEG_275: f32 = 275.0 * core::f32::consts::PI / 180.0;

#[allow(clippy::excessive_precision)]
fn ciede2000(l1: f32, a1: f32, b1: f32, l2: f32, a2: f32, b2: f32) -> f32 {
    let c1 = (a1 * a1 + b1 * b1).sqrt();
    let c2 = (a2 * a2 + b2 * b2).sqrt();
    let c_bar = (c1 + c2) * 0.5;
    let c_bar7 = c_bar * c_bar * c_bar * c_bar * c_bar * c_bar * c_bar;
    let g = 0.5 * (1.0 - (c_bar7 / (c_bar7 + CIEDE2000_C25_7))).sqrt();
    let a1p = a1 * (1.0 + g);
    let a2p = a2 * (1.0 + g);

    let c1p = (a1p * a1p + b1 * b1).sqrt();
    let c2p = (a2p * a2p + b2 * b2).sqrt();
    let h1p = (b1).atan2(a1p).rem_euclid(core::f32::consts::TAU);
    let h2p = (b2).atan2(a2p).rem_euclid(core::f32::consts::TAU);

    let dl_p = l2 - l1;
    let dc_p = c2p - c1p;
    let dh_p = if c1p * c2p == 0.0 {
        0.0
    } else {
        (h2p - h1p + core::f32::consts::PI).rem_euclid(core::f32::consts::TAU)
            - core::f32::consts::PI
    };
    let dh_big = 2.0 * (c1p * c2p).sqrt() * (dh_p * 0.5).sin();

    let l_bar = (l1 + l2) * 0.5;
    let c_barp = (c1p + c2p) * 0.5;
    let h_barp = if c1p * c2p == 0.0 {
        h1p + h2p
    } else if (h1p - h2p).abs() <= core::f32::consts::PI {
        (h1p + h2p) * 0.5
    } else if h1p + h2p < core::f32::consts::TAU {
        (h1p + h2p + core::f32::consts::TAU) * 0.5
    } else {
        (h1p + h2p - core::f32::consts::TAU) * 0.5
    };

    let t = 1.0
        - 0.17 * (h_barp - core::f32::consts::FRAC_PI_6).cos()
        + 0.24 * (2.0 * h_barp).cos()
        + 0.32 * (3.0 * h_barp + DEG_6).cos()
        - 0.20 * (4.0 * h_barp - DEG_63).cos();

    let l_bar_shifted = l_bar - 50.0;
    let sl = 1.0 + 0.015 * l_bar_shifted * l_bar_shifted
        / (20.0 + l_bar_shifted * l_bar_shifted).sqrt();
    let sc = 1.0 + 0.045 * c_barp;
    let sh = 1.0 + 0.015 * c_barp * t;

    let c_barp7 = { let c2 = c_barp * c_barp; let c4 = c2 * c2; c4 * c2 * c_barp };
    let rc = 2.0 * (c_barp7 / (c_barp7 + CIEDE2000_C25_7)).sqrt();
    let d_theta = DEG_60
        * (-((h_barp - DEG_275) / DEG_25).powi(2)).exp();
    let rt = -(2.0 * d_theta).sin() * rc;

    let term_l = dl_p / sl;
    let term_c = dc_p / sc;
    let term_h = dh_big / sh;

    (term_l * term_l + term_c * term_c + term_h * term_h + rt * term_c * term_h).sqrt()
}

macro_rules! impl_delta_e_lab {
    ($n:literal, $bound:tt) => {
        impl<$bound: Illuminant> Color<[f32; $n], Lab<$bound>> {
            /// CIE76 color difference. Euclidean distance in L*a*b*.
            ///
            /// Fast but overestimates differences in saturated regions.
            /// Prefer `delta_e2000` when perceptual accuracy matters.
            #[inline(always)]
            pub fn delta_e76(self, other: Self) -> f32 {
                let [l1, a1, b1, ..] = self.inner();
                let [l2, a2, b2, ..] = other.inner();
                euclidean(l1, a1, b1, l2, a2, b2)
            }

            /// CIEDE2000 color difference.
            ///
            /// Reference: CIE 142-2001, Luo, Cui, Rigg (2001).
            /// Parametric factors kL, kC, kH are all 1 (graphic arts default).
            pub fn delta_e2000(self, other: Self) -> f32 {
                let [l1, a1, b1, ..] = self.inner();
                let [l2, a2, b2, ..] = other.inner();
                ciede2000(l1, a1, b1, l2, a2, b2)
            }
        }
    };
}

macro_rules! impl_delta_e_oklab {
    ($n:literal) => {
        impl Color<[f32; $n], Oklab> {
            /// Oklab color difference. Euclidean distance in Oklab.
            ///
            /// Oklab is designed to be perceptually uniform, so this
            /// metric is effective for most use cases without the
            /// complexity of CIEDE2000.
            #[inline(always)]
            pub fn delta_e_ok(self, other: Self) -> f32 {
                let [l1, a1, b1, ..] = self.inner();
                let [l2, a2, b2, ..] = other.inner();
                euclidean(l1, a1, b1, l2, a2, b2)
            }
        }
    };
}

impl_delta_e_lab!(3, W);
impl_delta_e_lab!(4, W);
impl_delta_e_oklab!(3);
impl_delta_e_oklab!(4);