diff --git a/crates/cms/README.md b/crates/cms/README.md --- a/crates/cms/README.md +++ b/crates/cms/README.md @@ -30,6 +30,11 @@ Buffers are checked once per call, before any pixel is written. * `Transform::proof(src, proof, display, intent, bpc, simulate_paper)`: soft-proof chain. * `GamutCheck::new(src, proof, ΔE)`: out-of-gamut test (gamut warning). + Uses ΔE76 and the supplied threshold; `gamut::DEFAULT_THRESHOLD` is a constant of 4.0. + Opt into CIEDE2000 with + `GamutCheck::with_metric(src, proof, gamut::Metric::De2000, threshold)`; + `delta_e` and the threshold use the selected metric's units. Both metrics use the same + grid sizes and interpolation, so values between grid nodes are approximate. * `Lut3d::from_transform(&t, 33)?.to_rgba16f_bytes()`: display transform as a 3D texture. * `transform::cached(src, dst, opts)`: process-wide transform cache. * `lutfile::parse(name, bytes)` / `lutfile::write_cube(&lut)`: `.cube`, `.3dl` and `.look` LUT diff --git a/crates/cms/src/gamut.rs b/crates/cms/src/gamut.rs --- a/crates/cms/src/gamut.rs +++ b/crates/cms/src/gamut.rs @@ -4,26 +4,47 @@ use crate::clut::Clut; use crate::pipeline::{Pipeline, Stage}; use crate::profile::Profile; use crate::transform::link_stages; -use crate::{CmsError, Intent, check_len}; +use crate::{CmsError, Intent, check_len, math}; /// Default ΔE76 threshold above which a colour counts as out of gamut. pub const DEFAULT_THRESHOLD: f32 = 4.0; +/// Colour-difference metric used for a proof round trip. +#[derive(Clone, Copy, Debug, PartialEq, Eq)] +pub enum Metric { + /// CIE 1976 ΔE*ab (the default for [`GamutCheck::new`]). + De76, + /// CIEDE2000 with kL = kC = kH = 1. + De2000, +} + /// Precomputed out-of-gamut test for colours of `src` against `proof`. /// /// A colour is out of gamut when sending it to the proof device (relative colorimetric, -/// clipping) and back to Lab moves it by more than the threshold (ΔE76). The ΔE is tabulated -/// over the source device space (33³ grid for 3 channels, 17⁴ for 4, exact otherwise). +/// clipping) and back to Lab moves it by more than the threshold in the selected metric. +/// The ΔE is tabulated over the source device space (256 points for 1 channel, 33³ for 3, +/// 17⁴ for 4, exact otherwise). Between grid nodes the result is interpolated. #[derive(Debug)] pub struct GamutCheck { inputs: usize, grid: Option, exact: (Pipeline, Pipeline), + metric: Metric, + /// Warning threshold in the selected metric's units. pub threshold: f32, } impl GamutCheck { + /// Builds a gamut check using ΔE76 and the supplied threshold. pub fn new(src: &Profile, proof: &Profile, threshold: f32) -> Result { + Self::with_metric(src, proof, Metric::De76, threshold) + } + + /// Builds a gamut check with `threshold` in the chosen metric's units. + /// + /// The metric is fixed when the lookup table is built. Thresholds are not interchangeable + /// between metrics; [`DEFAULT_THRESHOLD`] remains the ΔE76 default. + pub fn with_metric(src: &Profile, proof: &Profile, metric: Metric, threshold: f32) -> Result { let lab = crate::Builtin::LabD50.profile(); let direct = Pipeline::new(src.channels(), link_stages(src, lab, Intent::RelativeColorimetric, false)?)?; let mut via = link_stages(src, proof, Intent::RelativeColorimetric, false)?; @@ -34,8 +55,13 @@ impl GamutCheck { let (mut a, mut b) = ([0.0f32; 16], [0.0f32; 16]); direct.run(p, &mut a); via.run(p, &mut b); - let d = |i: usize, s: f32| (a[i] - b[i]) * s; - o[0] = (d(0, 100.0).powi(2) + d(1, 255.0).powi(2) + d(2, 255.0).powi(2)).sqrt(); + o[0] = match metric { + Metric::De76 => { + let d = |i: usize, s: f32| (a[i] - b[i]) * s; + (d(0, 100.0).powi(2) + d(1, 255.0).powi(2) + d(2, 255.0).powi(2)).sqrt() + } + Metric::De2000 => delta_e2000(&a, &b), + }; }; let n = src.channels(); let grid = match n { @@ -44,14 +70,15 @@ impl GamutCheck { 4 => Some(Clut::sample(vec![17; 4], 1, de)?), _ => None, }; - Ok(GamutCheck { inputs: n, grid, exact: (direct, via), threshold }) + Ok(GamutCheck { inputs: n, grid, exact: (direct, via), metric, threshold }) } pub fn inputs(&self) -> usize { self.inputs } - /// ΔE76 between the colour (at least `inputs` normalised values) and its proof round trip. + /// ΔE in the selected metric between the colour (at least `inputs` normalised values) + /// and its proof round trip. pub fn delta_e(&self, px: &[f32]) -> Result { check_len("gamut check input", px.len(), self.inputs)?; let mut o = [0.0f32; 1]; @@ -67,7 +94,10 @@ impl GamutCheck { let (mut a, mut b) = ([0.0f32; 16], [0.0f32; 16]); self.exact.0.run(px, &mut a); self.exact.1.run(px, &mut b); - o[0] = ((a[0] - b[0]) * 100.0).hypot((a[1] - b[1]) * 255.0).hypot((a[2] - b[2]) * 255.0); + o[0] = match self.metric { + Metric::De76 => ((a[0] - b[0]) * 100.0).hypot((a[1] - b[1]) * 255.0).hypot((a[2] - b[2]) * 255.0), + Metric::De2000 => delta_e2000(&a, &b), + }; } } Ok(o[0]) @@ -77,3 +107,101 @@ impl GamutCheck { Ok(self.delta_e(px)? > self.threshold) } } + +fn delta_e2000(a: &[f32], b: &[f32]) -> f32 { + let lab = |v: &[f32]| [v[0] as f64 * 100.0, v[1] as f64 * 255.0 - 128.0, v[2] as f64 * 255.0 - 128.0]; + math::delta_e2000(lab(a), lab(b)) as f32 +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::profile::{ColorSpace, Lut, LutKind, ProfileClass}; + use crate::{Builtin, math}; + + fn constant_lab_profile(inputs: usize, lab: [f64; 3]) -> Profile { + let mut p = Builtin::LabD50.profile().clone(); + p.class = ProfileClass::Input; + p.color_space = if inputs == 1 { ColorSpace::Gray } else { ColorSpace::Color(inputs as u8) }; + p.a2b[1] = Some(Lut { + kind: LutKind::Ab, + inputs, + outputs: 3, + stages: vec![Stage::Matrix { + rows: 3, + cols: inputs, + m: vec![0.0; 3 * inputs], + offset: vec![lab[0] / 100.0, (lab[1] + 128.0) / 255.0, (lab[2] + 128.0) / 255.0], + }], + }); + p + } + + #[test] + fn selected_metric_reaches_sampled_and_exact_paths() -> Result<(), CmsError> { + // Sharma reference pair 1: the same known difference through each channel-count route. + let a = [50.0, 2.6772, -79.7751]; + let b = [50.0, 0.0, -82.7485]; + let mut proof = constant_lab_profile(3, b); + proof.color_space = ColorSpace::Lab; + for inputs in [1, 2, 3, 4] { + let src = constant_lab_profile(inputs, a); + let legacy = GamutCheck::new(&src, &proof, DEFAULT_THRESHOLD)?; + let de76 = GamutCheck::with_metric(&src, &proof, Metric::De76, DEFAULT_THRESHOLD)?; + let mut de2000 = GamutCheck::with_metric(&src, &proof, Metric::De2000, DEFAULT_THRESHOLD)?; + let grid = match inputs { + 1 => Some(vec![256]), + 3 => Some(vec![33; 3]), + 4 => Some(vec![17; 4]), + _ => None, + }; + assert_eq!(de2000.grid.as_ref().map(|g| g.grid().to_vec()), grid); + assert_eq!(de2000.inputs(), inputs); + for value in [0.0, 0.37, 1.0] { + let px = vec![value; inputs]; + assert_eq!(legacy.delta_e(&px)?, de76.delta_e(&px)?); + assert!((de76.delta_e(&px)? as f64 - math::delta_e76(a, b)).abs() < 1e-4); + assert!((de2000.delta_e(&px)? - 2.0425).abs() < 1e-4); + assert!(legacy.out_of_gamut(&px)?); + assert!(!de2000.out_of_gamut(&px)?); + } + let px = vec![0.5; inputs]; + de2000.threshold = de2000.delta_e(&px)?; + assert!(!de2000.out_of_gamut(&px)?); + de2000.threshold -= 0.01; + assert!(de2000.out_of_gamut(&px)?); + assert!(matches!(de2000.delta_e(&px[..inputs - 1]), Err(CmsError::InvalidArgument(_)))); + } + Ok(()) + } + + #[test] + fn de2000_grid_tracks_exact_display_p3_to_srgb_roundtrips() -> Result<(), CmsError> { + let g = GamutCheck::with_metric(Builtin::DisplayP3.profile(), Builtin::Srgb.profile(), Metric::De2000, 2.0)?; + let exact = |px: &[f32]| { + let (mut a, mut b) = ([0.0; 3], [0.0; 3]); + g.exact.0.run(px, &mut a); + g.exact.1.run(px, &mut b); + let decode = |l: [f32; 3]| [l[0] as f64 * 100.0, l[1] as f64 * 255.0 - 128.0, l[2] as f64 * 255.0 - 128.0]; + math::delta_e2000(decode(a), decode(b)) as f32 + }; + for r in [0.0, 0.5, 1.0] { + for g0 in [0.0, 0.5, 1.0] { + for b in [0.0, 0.5, 1.0] { + let px = [r, g0, b]; + assert!((g.delta_e(&px)? - exact(&px)).abs() < 1e-4, "node {px:?}"); + } + } + } + for r in [0.03, 0.17, 0.37, 0.61, 0.79, 0.97] { + for g0 in [0.03, 0.17, 0.37, 0.61, 0.79, 0.97] { + for b in [0.03, 0.17, 0.37, 0.61, 0.79, 0.97] { + let px = [r, g0, b]; + let error = (g.delta_e(&px)? - exact(&px)).abs(); + assert!(error < 0.5, "interpolated {px:?}: error {error}"); + } + } + } + Ok(()) + } +} diff --git a/crates/cms/src/math.rs b/crates/cms/src/math.rs --- a/crates/cms/src/math.rs +++ b/crates/cms/src/math.rs @@ -234,10 +234,19 @@ mod tests { #[test] fn delta_e2000_matches_sharma_reference_pairs() { - // Sharma, Wu & Dalal (2005), test data pairs 1, 7, 17, 25 and 34. + // Sharma, Wu & Dalal (2005), pairs 1, 7, 9–17, 25 and 29, including hue-wrap boundaries. + // https://hajim.rochester.edu/ece/sites/gsharma/ciede2000/dataNprograms/ciede2000testdata.txt let pairs = [ ([50.0, 2.6772, -79.7751], [50.0, 0.0, -82.7485], 2.0425), ([50.0, 0.0, 0.0], [50.0, -1.0, 2.0], 2.3669), + ([50.0, 2.49, -0.001], [50.0, -2.49, 0.0009], 7.1792), + ([50.0, 2.49, -0.001], [50.0, -2.49, 0.001], 7.1792), + ([50.0, 2.49, -0.001], [50.0, -2.49, 0.0011], 7.2195), + ([50.0, 2.49, -0.001], [50.0, -2.49, 0.0012], 7.2195), + ([50.0, -0.001, 2.49], [50.0, 0.0009, -2.49], 4.8045), + ([50.0, -0.001, 2.49], [50.0, 0.001, -2.49], 4.8045), + ([50.0, -0.001, 2.49], [50.0, 0.0011, -2.49], 4.7461), + ([50.0, 2.5, 0.0], [50.0, 0.0, -2.5], 4.3065), ([50.0, 2.5, 0.0], [73.0, 25.0, -18.0], 27.1492), ([60.2574, -34.0099, 36.2677], [60.4626, -34.1751, 39.4387], 1.2644), ([22.7233, 20.0904, -46.6940], [23.0331, 14.9730, -42.5619], 2.0373), @@ -248,6 +257,26 @@ mod tests { } } + #[test] + fn delta_e2000_is_finite_symmetric_and_zero_for_identical_lab() { + let mut colours = Vec::new(); + for l in [0.0, 50.0, 100.0] { + for a in [-128.0, 0.0, 127.0] { + for b in [-128.0, 0.0, 127.0] { + colours.push([l, a, b]); + } + } + } + for &a in &colours { + assert_eq!(delta_e2000(a, a), 0.0); + for &b in &colours { + let d = delta_e2000(a, b); + assert!(d.is_finite() && d >= 0.0, "{a:?} {b:?}: {d}"); + assert!((d - delta_e2000(b, a)).abs() < 1e-12); + } + } + } + #[test] fn bradford_maps_white() { let d65 = xy_to_xyz(D65_XY);