diff --git a/packages/rust/core/src/algorithms.rs b/packages/rust/core/src/algorithms.rs index 291215a..a67e996 100644 --- a/packages/rust/core/src/algorithms.rs +++ b/packages/rust/core/src/algorithms.rs @@ -307,16 +307,37 @@ pub fn direct_map(pixels: &[u8], palette: &Palette, canonical_palette: &Palette) // ── Ordered (Bayer) dithering ───────────────────────────────────────────────── -// 4×4 Bayer matrix, normalized to [-0.5, 0.5]. Indexed as [y % 4][x % 4]. -// Values = (bayer_entry / 16.0) - 0.5 for the standard ordered Bayer matrix. -// Written as exact fractions to keep precision (all are multiples of 1/16). +// 4×4 Bayer matrix, zero-mean thresholds in (-0.5, 0.5). Indexed as [y % 4][x % 4]. +// Values = (bayer_entry + 0.5) / 16.0 - 0.5, the standard zero-mean normalization — +// the +0.5 centering makes ordered dithering brightness-neutral (mean threshold = 0). +// Written as exact fractions to keep precision (all are odd multiples of 1/32). const BAYER_4X4: [[f64; 4]; 4] = [ - [-0.5000, 0.0000, -0.3750, 0.1250], - [ 0.2500, -0.2500, 0.3750, -0.1250], - [-0.3125, 0.1875, -0.4375, 0.0625], - [ 0.4375, -0.0625, 0.3125, -0.1875], + [-0.46875, 0.03125, -0.34375, 0.15625], + [ 0.28125, -0.21875, 0.40625, -0.09375], + [-0.28125, 0.21875, -0.40625, 0.09375], + [ 0.46875, -0.03125, 0.34375, -0.15625], ]; +/// Threshold amplitude (in sRGB-fraction space) for ordered dither on this palette. +/// +/// Ordered dithering trades quantization error for a fixed threshold pattern; the +/// threshold should be scaled to the palette's quantization step. For an evenly-spaced +/// grayscale ramp of `n` levels the step is `1/(n-1)`, so a full ±0.5 threshold is +/// `2·(n-1)×` too large (e.g. 14× for 16-level gray) and swamps fine detail with noise. +/// +/// Sparse palettes (mono + color) keep the full ±0.5 spread: their transitions are +/// dominated by black↔white/ink decisions where the wide threshold is the tuned, +/// tested behavior (see the issue-#27 property test). +fn ordered_spread(palette: &Palette) -> f64 { + let is_grayscale = + palette.colors.len() >= 3 && palette.colors.iter().all(|&[r, g, b]| r == g && g == b); + if is_grayscale { + 1.0 / (palette.colors.len() - 1) as f64 + } else { + 1.0 + } +} + /// Ordered (Bayer 4×4) dither. Pixels are independent — parallelized with rayon. /// /// The Bayer threshold is added in sRGB-fraction space, not linear. Linear-space @@ -324,6 +345,9 @@ const BAYER_4X4: [[f64; 4]; 4] = [ /// (the sRGB gamma is convex, so a fixed linear ±0.5 step is huge near 0 and tiny near 1). /// sRGB-space thresholding gives uniform perceptual dot density across the tonal range, /// matching how error diffusion already accumulates error. See GitHub issue #27. +/// +/// The threshold amplitude is scaled by `ordered_spread` to the palette's quantization +/// step so dense grayscale ramps are not swamped by full-range dither noise. pub fn ordered_dither(pixels: &[u8], width: usize, palette: &Palette) -> Vec { ordered_dither_impl(pixels, width, palette, None) } @@ -344,6 +368,21 @@ fn ordered_dither_impl( canonical_palette: Option<&Palette>, ) -> Vec { let (_palette_linear, palette_lab) = build_palette_lab(palette); + let spread = ordered_spread(palette); + + // Per-threshold gamma LUT: lut[t][v] = srgb_fraction_to_linear((v/255 + threshold_t).clamp). + // Only 16 thresholds × 256 byte values exist, so this removes all three per-pixel powf + // calls. Built from the exact same inputs as the scalar path, so results are bit-identical. + let lut: Vec<[f64; 256]> = (0..16) + .map(|t| { + let threshold = BAYER_4X4[t / 4][t % 4] * spread; + let mut row = [0.0_f64; 256]; + for (v, slot) in row.iter_mut().enumerate() { + *slot = srgb_fraction_to_linear((v as f64 / 255.0 + threshold).clamp(0.0, 1.0)); + } + row + }) + .collect(); pixels .par_chunks(3) @@ -357,18 +396,12 @@ fn ordered_dither_impl( let x = i % width; let y = i / width; - - let threshold = BAYER_4X4[y % 4][x % 4]; - - // Add threshold in sRGB-fraction space, then convert to linear for OKLab match. - let r_srgb = (rgb[0] as f64 / 255.0 + threshold).clamp(0.0, 1.0); - let g_srgb = (rgb[1] as f64 / 255.0 + threshold).clamp(0.0, 1.0); - let b_srgb = (rgb[2] as f64 / 255.0 + threshold).clamp(0.0, 1.0); + let t = (y % 4) * 4 + (x % 4); let lab = rgb_to_oklab( - srgb_fraction_to_linear(r_srgb), - srgb_fraction_to_linear(g_srgb), - srgb_fraction_to_linear(b_srgb), + lut[t][rgb[0] as usize], + lut[t][rgb[1] as usize], + lut[t][rgb[2] as usize], ); match_pixel_oklab(lab, &palette_lab, WAB) as u8 }) @@ -464,7 +497,11 @@ mod tests { /// /// We measure transitions (adjacent columns with differing palette indices) per band /// in the mid-tone region and assert the max/min ratio stays modest. Empirically the - /// linear-space implementation produces ratio ≈ 13.5; sRGB-space ≈ 4.3. + /// linear-space implementation produces ratio ≈ 13.5; sRGB-space ≈ 2.2. + /// + /// The mono decision midpoint sits at OKLab L=0.5 ≈ sRGB 188, so the top two bands + /// (sRGB ≳ 192) are already in the highlight roll-off where dither activity legitimately + /// tapers; we exclude them along with the pure-black first band. #[test] fn ordered_dither_activity_is_perceptually_uniform() { const W: usize = 256; @@ -492,9 +529,10 @@ mod tests { } } - // Mid-tone bands (skip the first and last — they're near pure black/white where - // the threshold is clamped and dither activity legitimately falls off). - let mid = &transitions[1..BANDS - 1]; + // Mid-tone bands: skip the pure-black first band and the top two near-white bands + // (past the mono midpoint ≈ sRGB 188) where the threshold clamps and dither + // activity legitimately falls off. + let mid = &transitions[1..BANDS - 2]; let max = *mid.iter().max().unwrap(); let min = *mid.iter().min().unwrap(); assert!(min > 0, "every mid-tone band should have at least one transition: {transitions:?}"); @@ -507,6 +545,84 @@ mod tests { ); } + #[test] + fn bayer_matrix_is_zero_mean() { + // A biased Bayer matrix systematically darkens (or brightens) every ordered- + // dithered image. The standard (v + 0.5)/16 - 0.5 normalization sums to zero. + let sum: f64 = BAYER_4X4.iter().flatten().sum(); + assert!(sum.abs() < 1e-12, "Bayer thresholds must sum to zero, got {sum}"); + } + + #[test] + fn ordered_dither_is_brightness_neutral_at_midpoint() { + // Find the mono decision midpoint: the smallest gray value whose threshold-free + // OKLab match flips from black (0) to white (1). + let pal = ColorScheme::Mono.palette(); + let (_, pal_lab) = build_palette_lab(pal); + let midpoint = (0u16..=255) + .find(|&v| { + let lin = srgb_channel_to_linear(v as u8); + let lab = rgb_to_oklab(lin, lin, lin); + match_pixel_oklab(lab, &pal_lab, WAB) == 1 + }) + .expect("mono palette must have a black→white transition") as u8; + + // Dither a solid midpoint-gray field. With a zero-mean threshold the white + // fraction should sit near 0.5; the old -1/32 bias produced ≈ 7/16. + const N: usize = 32; + let pixels = solid_image(midpoint, midpoint, midpoint, N, N); + let out = ordered_dither(&pixels, N, pal); + let white = out.iter().filter(|&&v| v == 1).count() as f64 / out.len() as f64; + assert!( + (white - 0.5).abs() <= 1.0 / 16.0, + "ordered dither at the mono midpoint should be ~50% white, got {white}" + ); + } + + #[test] + fn ordered_dither_grayscale16_tracks_ramp_and_uses_many_levels() { + // On a horizontal sRGB ramp the full ±0.5 threshold swamps the 17/255 GS16 step + // (14× too large). After scaling by ordered_spread, block-averaged output should + // track the input ramp and exercise most of the 16 levels. + const W: usize = 256; + const H: usize = 4; + let pal = ColorScheme::Grayscale16.palette(); + + let mut pixels = Vec::with_capacity(W * H * 3); + for _ in 0..H { + for x in 0..W { + let v = x as u8; + pixels.extend_from_slice(&[v, v, v]); + } + } + let out = ordered_dither(&pixels, W, pal); + + // (a) 4-wide block average of palette level tracks the input level. + let levels: Vec = pal.colors.iter().map(|&[r, _, _]| r as f64).collect(); + let mut max_err = 0.0_f64; + for bx in 0..W / 4 { + let mut sum = 0.0; + for y in 0..H { + for dx in 0..4 { + let x = bx * 4 + dx; + sum += levels[out[y * W + x] as usize]; + } + } + let avg = sum / (H * 4) as f64; + let input = (bx * 4 + 1) as f64; // center-ish input value + max_err = max_err.max((avg - input).abs()); + } + assert!(max_err < 24.0, "block-averaged output should track the ramp, max_err={max_err}"); + + // (b) most of the 16 levels are used. + let mut seen = [false; 16]; + for &idx in &out { + seen[idx as usize] = true; + } + let used = seen.iter().filter(|&&s| s).count(); + assert!(used >= 12, "GS16 ordered dither should use ≥12 levels, used {used}"); + } + #[test] fn serpentine_first_row_matches_raster() { // Row 0 is always scanned left-to-right regardless of serpentine flag diff --git a/packages/rust/core/src/color_space.rs b/packages/rust/core/src/color_space.rs index 518e29b..f597540 100644 --- a/packages/rust/core/src/color_space.rs +++ b/packages/rust/core/src/color_space.rs @@ -15,14 +15,21 @@ pub fn srgb_fraction_to_linear(value: f64) -> f64 { } } -/// Linear [0.0–1.0] → sRGB [0–255]. Inverse of `srgb_channel_to_linear`. -pub fn linear_channel_to_srgb(linear: f64) -> u8 { - let srgb = if linear <= 0.0031308 { +/// Linear [0.0–1.0] → continuous sRGB fraction [0.0–1.0]. Inverse of +/// `srgb_fraction_to_linear`; the continuous counterpart of `linear_channel_to_srgb`, +/// useful when a gamma-encoded value is needed without quantizing to a byte (e.g. applying +/// a tone curve about the perceptual mid-gray). +pub fn linear_fraction_to_srgb(linear: f64) -> f64 { + if linear <= 0.0031308 { linear * 12.92 } else { 1.055 * linear.powf(1.0 / 2.4) - 0.055 - }; - (srgb * 255.0).round() as u8 + } +} + +/// Linear [0.0–1.0] → sRGB [0–255]. Inverse of `srgb_channel_to_linear`. +pub fn linear_channel_to_srgb(linear: f64) -> u8 { + (linear_fraction_to_srgb(linear) * 255.0).round() as u8 } diff --git a/packages/rust/core/src/color_space_lab.rs b/packages/rust/core/src/color_space_lab.rs index 8669912..fcaa5b5 100644 --- a/packages/rust/core/src/color_space_lab.rs +++ b/packages/rust/core/src/color_space_lab.rs @@ -1,17 +1,12 @@ //! OKLab color space and color matching for dithering. -// sRGB -> XYZ matrix (D65 illuminant, BruceLinbloom) -const M_RGB_XYZ: [[f64; 3]; 3] = [ - [0.4124564, 0.3575761, 0.1804375], - [0.2126729, 0.7151522, 0.0721750], - [0.0193339, 0.1191920, 0.9503041], -]; - -// XYZ -> LMS (M1, Ottosson) -const M1: [[f64; 3]; 3] = [ - [0.8189330101, 0.3618667424, -0.1288597137], - [0.0329845436, 0.9293118715, 0.0361456387], - [0.0482003018, 0.2643662691, 0.6338517070], +// Linear sRGB -> LMS (Ottosson's published combined matrix). This folds the +// sRGB→XYZ and XYZ→LMS steps into one 3×3, halving the per-pixel matrix work in the +// hot path and avoiding the drift of composing separately-sourced matrices. +const M_RGB_LMS: [[f64; 3]; 3] = [ + [0.4122214708, 0.5363325363, 0.0514459929], + [0.2119034982, 0.6806995451, 0.1073969566], + [0.0883024619, 0.2817188376, 0.6299787005], ]; // cbrt(LMS) -> OKLab (M2, Ottosson) @@ -36,15 +31,11 @@ pub struct OkLab { pub b: f64, } -/// Linear RGB → OKLab. Pipeline: RGB → XYZ → LMS → cbrt → OKLab. +/// Linear RGB → OKLab. Pipeline: RGB → LMS → cbrt → OKLab. pub fn rgb_to_oklab(r: f64, g: f64, b: f64) -> OkLab { - let x = M_RGB_XYZ[0][0] * r + M_RGB_XYZ[0][1] * g + M_RGB_XYZ[0][2] * b; - let y = M_RGB_XYZ[1][0] * r + M_RGB_XYZ[1][1] * g + M_RGB_XYZ[1][2] * b; - let z = M_RGB_XYZ[2][0] * r + M_RGB_XYZ[2][1] * g + M_RGB_XYZ[2][2] * b; - - let l = M1[0][0] * x + M1[0][1] * y + M1[0][2] * z; - let m = M1[1][0] * x + M1[1][1] * y + M1[1][2] * z; - let s = M1[2][0] * x + M1[2][1] * y + M1[2][2] * z; + let l = M_RGB_LMS[0][0] * r + M_RGB_LMS[0][1] * g + M_RGB_LMS[0][2] * b; + let m = M_RGB_LMS[1][0] * r + M_RGB_LMS[1][1] * g + M_RGB_LMS[1][2] * b; + let s = M_RGB_LMS[2][0] * r + M_RGB_LMS[2][1] * g + M_RGB_LMS[2][2] * b; let l_ = l.cbrt(); let m_ = m.cbrt(); @@ -76,18 +67,11 @@ const M2_INV: [[f64; 3]; 3] = [ [1.0, -0.0894841775, -1.2914855480], ]; -// M1_inv: LMS -> XYZ (Ottosson) -const M1_INV: [[f64; 3]; 3] = [ - [ 1.2270138511035211, -0.5577999806518222, 0.2812561489664678], - [-0.0405801784232806, 1.1122568696168302, -0.0716766786656012], - [-0.0763812845057069, -0.4214819784180127, 1.5861632204407947], -]; - -// XYZ -> linear sRGB (IEC 61966-2-1 D65) -const M_XYZ_RGB: [[f64; 3]; 3] = [ - [ 3.2404542, -1.5371385, -0.4985314], - [-0.9692660, 1.8760108, 0.0415560], - [ 0.0556434, -0.2040259, 1.0572252], +// LMS -> linear sRGB (Ottosson's published combined inverse of M_RGB_LMS). +const M_LMS_RGB: [[f64; 3]; 3] = [ + [ 4.0767416621, -3.3077115913, 0.2309699292], + [-1.2684380046, 2.6097574011, -0.3413193965], + [-0.0041960863, -0.7034186147, 1.7076147010], ]; /// OKLab → linear RGB. Output is clamped to [0, 1]. Inverse of `rgb_to_oklab`. @@ -100,13 +84,9 @@ pub fn oklab_to_rgb(lab: OkLab) -> [f64; 3] { let m = m_ * m_ * m_; let s = s_ * s_ * s_; - let x = M1_INV[0][0] * l + M1_INV[0][1] * m + M1_INV[0][2] * s; - let y = M1_INV[1][0] * l + M1_INV[1][1] * m + M1_INV[1][2] * s; - let z = M1_INV[2][0] * l + M1_INV[2][1] * m + M1_INV[2][2] * s; - - let r = M_XYZ_RGB[0][0] * x + M_XYZ_RGB[0][1] * y + M_XYZ_RGB[0][2] * z; - let g = M_XYZ_RGB[1][0] * x + M_XYZ_RGB[1][1] * y + M_XYZ_RGB[1][2] * z; - let b = M_XYZ_RGB[2][0] * x + M_XYZ_RGB[2][1] * y + M_XYZ_RGB[2][2] * z; + let r = M_LMS_RGB[0][0] * l + M_LMS_RGB[0][1] * m + M_LMS_RGB[0][2] * s; + let g = M_LMS_RGB[1][0] * l + M_LMS_RGB[1][1] * m + M_LMS_RGB[1][2] * s; + let b = M_LMS_RGB[2][0] * l + M_LMS_RGB[2][1] * m + M_LMS_RGB[2][2] * s; [r.clamp(0.0, 1.0), g.clamp(0.0, 1.0), b.clamp(0.0, 1.0)] } @@ -151,9 +131,30 @@ mod tests { #[test] fn white_l_is_one() { let lab = rgb_to_oklab(1.0, 1.0, 1.0); - assert_relative_eq!(lab.l, 1.0, epsilon = 1e-4); - assert_relative_eq!(lab.a, 0.0, epsilon = 1e-4); - assert_relative_eq!(lab.b, 0.0, epsilon = 1e-4); + // With Ottosson's combined matrix, white maps to L=1, a=b=0 to full precision. + assert_relative_eq!(lab.l, 1.0, epsilon = 1e-6); + assert_relative_eq!(lab.a, 0.0, epsilon = 1e-6); + assert_relative_eq!(lab.b, 0.0, epsilon = 1e-6); + } + + #[test] + fn oklab_roundtrip_is_identity() { + // rgb_to_oklab and oklab_to_rgb must be true inverses — catches any mis-transcribed + // coefficient in the combined forward/inverse matrices (e.g. a bad row breaks white). + for &rgb in &[ + [1.0, 1.0, 1.0_f64], + [0.0, 0.0, 0.0], + [0.5, 0.5, 0.5], + [0.8, 0.2, 0.1], + [0.1, 0.6, 0.9], + [0.2126, 0.7152, 0.0722], + ] { + let lab = rgb_to_oklab(rgb[0], rgb[1], rgb[2]); + let back = oklab_to_rgb(lab); + for c in 0..3 { + assert_relative_eq!(back[c], rgb[c], epsilon = 1e-6); + } + } } #[test] diff --git a/packages/rust/core/src/tone_map.rs b/packages/rust/core/src/tone_map.rs index 03c833c..ebfba8b 100644 --- a/packages/rust/core/src/tone_map.rs +++ b/packages/rust/core/src/tone_map.rs @@ -5,7 +5,7 @@ use rayon::prelude::*; -use crate::color_space::srgb_channel_to_linear; +use crate::color_space::{linear_fraction_to_srgb, srgb_channel_to_linear, srgb_fraction_to_linear}; use crate::color_space_lab::{oklab_to_rgb, rgb_to_oklab, OkLab}; use crate::palettes::Palette; @@ -32,6 +32,24 @@ fn palette_to_linear(palette: &Palette) -> Vec<[f64; 3]> { .collect() } +/// Scale a near-black pixel toward a target luminance while preserving its channel ratios +/// (hue/chroma). Used by both dynamic-range compressors for the `y ≈ 0` branch where the +/// `target_y / y` scale would explode. When the pixel is exactly black, all channels are +/// set to `target_y` (neutral). +fn scale_dark_pixel(pixel: &mut [f64; 3], target_y: f64) { + let max_ch = pixel[0].max(pixel[1]).max(pixel[2]); + if max_ch > 1e-12 { + let scale = target_y / max_ch; + pixel[0] = (pixel[0] * scale).clamp(0.0, 1.0); + pixel[1] = (pixel[1] * scale).clamp(0.0, 1.0); + pixel[2] = (pixel[2] * scale).clamp(0.0, 1.0); + } else { + pixel[0] = target_y; + pixel[1] = target_y; + pixel[2] = target_y; + } +} + /// Two percentiles from the same data in one sort. Returns (p_low_val, p_high_val). fn percentile_pair(values: &[f64], p_low: f64, p_high: f64) -> (f64, f64) { let mut sorted = values.to_vec(); @@ -87,24 +105,34 @@ pub fn adjust_saturation(pixels: &mut [[f64; 3]], factor: f64) { /// `shadows` controls lift on the lower half (0.0 = identity, 1.0 = strong; power = 1 − shadows). /// `highlights` controls compression on the upper half (0.0 = identity, 1.0 = strong; power = 1 + highlights). /// +/// The curve pivots about **perceptual mid-gray** (sRGB 0.5 ≈ linear 0.214), not linear 0.5 +/// (≈ sRGB 188). Pivoting in linear space would put ~74% of the perceptual tonal range in the +/// "shadows" half; gamma-space pivoting splits shadows/highlights at mid-gray as a user expects. +/// /// Each half is independent — `shadows = 0` leaves shadows alone even if `highlights > 0`. /// Hue is preserved by scaling all three channels by `Y' / Y`. pub fn apply_shadows_highlights(pixels: &mut [[f64; 3]], shadows: f64, highlights: f64) { if shadows <= 0.0 && highlights <= 0.0 { return; } + // Guard against degenerate exponents (negative → curve inversion, >1 → runaway). + let shadows = shadows.clamp(0.0, 1.0); + let highlights = highlights.clamp(0.0, 1.0); pixels.par_iter_mut().for_each(|pixel| { let y = luminance(*pixel); if y <= 1e-6 { return; } - let y_prime = if y <= 0.5 { - let t = y / 0.5; + // Apply the S-curve in gamma-encoded space so the pivot is at perceptual mid-gray. + let g = linear_fraction_to_srgb(y); + let g_prime = if g <= 0.5 { + let t = g / 0.5; 0.5 * t.powf(1.0 - shadows) } else { - let t = (y - 0.5) / 0.5; + let t = (g - 0.5) / 0.5; 0.5 + 0.5 * t.powf(1.0 + highlights) }; + let y_prime = srgb_fraction_to_linear(g_prime); let scale = (y_prime / y).clamp(0.0, 1e6); pixel[0] = (pixel[0] * scale).clamp(0.0, 1.0); pixel[1] = (pixel[1] * scale).clamp(0.0, 1.0); @@ -144,18 +172,7 @@ pub fn compress_dynamic_range(pixels: &mut [[f64; 3]], palette: &Palette, streng pixel[2] = (pixel[2] * scale).clamp(0.0, 1.0); } else { // Pixel luminance is near zero — preserve channel ratios, scale toward display black. - let blended_black = black_y * strength; - let max_ch = pixel[0].max(pixel[1]).max(pixel[2]); - if max_ch > 1e-12 { - let scale = blended_black / max_ch; - pixel[0] = (pixel[0] * scale).clamp(0.0, 1.0); - pixel[1] = (pixel[1] * scale).clamp(0.0, 1.0); - pixel[2] = (pixel[2] * scale).clamp(0.0, 1.0); - } else { - pixel[0] = blended_black; - pixel[1] = blended_black; - pixel[2] = blended_black; - } + scale_dark_pixel(pixel, black_y * strength); } }); } @@ -203,7 +220,10 @@ pub fn auto_compress_dynamic_range(pixels: &mut [[f64; 3]], palette: &Palette) { let log_range = log_max - log_min; if log_range > 1e-6 { let skew = (log_max - (l_lav + 1e-5).ln()) / log_range; - skew.powf(1.4).clamp(0.0, 1.0) + // Clamp before powf: for high-key images the log-average can slightly exceed + // log_max, making `skew` negative — `(-x).powf(1.4)` is NaN and would poison + // every output pixel. Clamping first keeps strength in [0, 1]. + skew.clamp(0.0, 1.0).powf(1.4) } else { 1.0 } @@ -211,10 +231,13 @@ pub fn auto_compress_dynamic_range(pixels: &mut [[f64; 3]], palette: &Palette) { 1.0 }; - // Remap [p_low, p_high] → [black_y, white_y] at computed strength + // Remap [p_low, p_high] → [black_y, white_y] at computed strength. + // `normalized` is clamped to [0, 1] so the darkest/brightest percentile outliers + // saturate at the display black/white points instead of extrapolating past them + // (unclamped, the bottom 2% would map below display black and crush to pure 0). pixels.par_iter_mut().for_each(|pixel| { let y = luminance(*pixel); - let normalized = (y - p_low) / image_range; + let normalized = ((y - p_low) / image_range).clamp(0.0, 1.0); let target_y_full = black_y + normalized * display_range; let target_y = y + strength * (target_y_full - y); if y > 1e-6 { @@ -223,7 +246,8 @@ pub fn auto_compress_dynamic_range(pixels: &mut [[f64; 3]], palette: &Palette) { pixel[1] = (pixel[1] * scale).clamp(0.0, 1.0); pixel[2] = (pixel[2] * scale).clamp(0.0, 1.0); } else { - *pixel = [black_y, black_y, black_y]; + // Near-black: honor strength and preserve chroma (y=0 → normalized=0 → target black_y). + scale_dark_pixel(pixel, target_y); } }); } @@ -284,6 +308,24 @@ pub fn gamut_compress(pixels: &mut [[f64; 3]], palette: &Palette, strength: f64) .collect(); let n = pal_linear.len(); + // Precompute the palette-hull edge geometry once (OKLab endpoints + linear endpoints for + // all i = (0..n) + .flat_map(|i| (i + 1..n).map(move |j| (i, j))) + .map(|(i, j)| Edge { + a_lab: [pal_lab[i].l, pal_lab[i].a, pal_lab[i].b], + b_lab: [pal_lab[j].l, pal_lab[j].a, pal_lab[j].b], + a_lin: pal_linear[i], + b_lin: pal_linear[j], + }) + .collect(); + pixels.par_iter_mut().for_each(|pixel| { let px_lab = rgb_to_oklab(pixel[0], pixel[1], pixel[2]); let px_lab_arr = [px_lab.l, px_lab.a, px_lab.b]; @@ -292,18 +334,14 @@ pub fn gamut_compress(pixels: &mut [[f64; 3]], palette: &Palette, strength: f64) let mut best_dist_sq = f64::INFINITY; let mut best_target = *pixel; - // Edges - for i in 0..n { - for j in (i + 1)..n { - let a_lab = [pal_lab[i].l, pal_lab[i].a, pal_lab[i].b]; - let b_lab = [pal_lab[j].l, pal_lab[j].a, pal_lab[j].b]; - let t = nearest_on_segment(px_lab_arr, a_lab, b_lab); - let nearest_lab = lerp3(a_lab, b_lab, t); - let d = dist_sq(px_lab_arr, nearest_lab); - if d < best_dist_sq { - best_dist_sq = d; - best_target = lerp3(pal_linear[i], pal_linear[j], t); - } + // Edges (endpoints precomputed above; vertices covered via clamped projection) + for e in &edges { + let t = nearest_on_segment(px_lab_arr, e.a_lab, e.b_lab); + let nearest_lab = lerp3(e.a_lab, e.b_lab, t); + let d = dist_sq(px_lab_arr, nearest_lab); + if d < best_dist_sq { + best_dist_sq = d; + best_target = lerp3(e.a_lin, e.b_lin, t); } } @@ -414,6 +452,55 @@ mod tests { } } + #[test] + fn auto_compress_does_not_crush_dark_outliers() { + // A high-key image (bright cluster overshooting white_y, triggering compression) with a + // couple of near-black outliers far below the 2nd percentile. Pre-fix, the unclamped + // remap sent `normalized` negative for those outliers, pulling their target luminance + // below display black and crushing them darker than they started. Post-clamp, a + // sub-percentile pixel saturates at the display black floor and is never darkened. + let pal = spectra_palette(); + + // 98 bright pixels spread 0.7..1.0 (so p2 lands ~0.7, well above the outliers), + // plus 2 near-black outliers at y ≈ 0.0004. + let mut pixels: Vec<[f64; 3]> = (0..98) + .map(|i| { + let v = 0.7 + 0.3 * (i as f64 / 97.0); + [v, v, v] + }) + .collect(); + let outlier = [0.0004_f64, 0.0004, 0.0004]; + pixels.push(outlier); + pixels.push(outlier); + + let in_y = luminance(outlier); + auto_compress_dynamic_range(&mut pixels, pal); + let out_y = luminance(pixels[98]); + assert!( + out_y >= in_y - 1e-9, + "sub-percentile dark outlier must not be crushed darker (in={in_y:.6}, out={out_y:.6})" + ); + } + + #[test] + fn auto_compress_preserves_dark_pixel_chroma() { + // A near-black but distinctly blue outlier must keep blue dominant after auto-compress, + // rather than being flattened to neutral display black. + let pal = spectra_palette(); + // Blue channel 1e-5 → luminance ≈ 7e-7 < 1e-6, so this exercises the near-zero dark + // branch, which pre-fix flattened the pixel to neutral display black. + let mut pixels = vec![[1.0_f64, 1.0, 1.0]; 96]; + pixels.extend(std::iter::repeat_n([0.0_f64, 0.0, 1e-5], 4)); + + auto_compress_dynamic_range(&mut pixels, pal); + + let [r, g, b] = pixels[96]; + assert!( + b >= r && b >= g && b > 0.0, + "blue should remain dominant after auto-compress of a dark blue pixel: [{r}, {g}, {b}]" + ); + } + #[test] fn gamut_compress_moves_out_of_gamut_pixel() { // Pure linear red [1, 0, 0] is much more vivid than the Spectra palette's desaturated diff --git a/packages/rust/core/tests/fixtures/references/cat__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/cat__burkes_spectra6_auto.bin index 6977943..10e71a4 100644 Binary files a/packages/rust/core/tests/fixtures/references/cat__burkes_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/cat__burkes_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/cat__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/cat__ordered_spectra6_auto.bin index 1712646..54ff7fd 100644 Binary files a/packages/rust/core/tests/fixtures/references/cat__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/cat__ordered_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/cat_orange__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/cat_orange__burkes_spectra6_auto.bin index 9d86b0d..812029e 100644 --- a/packages/rust/core/tests/fixtures/references/cat_orange__burkes_spectra6_auto.bin +++ b/packages/rust/core/tests/fixtures/references/cat_orange__burkes_spectra6_auto.bin @@ -1 +1 @@ - \ No newline at end of file + \ No newline at end of file diff --git a/packages/rust/core/tests/fixtures/references/cat_orange__floyd_steinberg_mono_raw.bin b/packages/rust/core/tests/fixtures/references/cat_orange__floyd_steinberg_mono_raw.bin index 506e517..68135be 100644 Binary files a/packages/rust/core/tests/fixtures/references/cat_orange__floyd_steinberg_mono_raw.bin and b/packages/rust/core/tests/fixtures/references/cat_orange__floyd_steinberg_mono_raw.bin differ diff --git a/packages/rust/core/tests/fixtures/references/cat_orange__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/cat_orange__ordered_spectra6_auto.bin index 689a080..1e1dc59 100644 Binary files a/packages/rust/core/tests/fixtures/references/cat_orange__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/cat_orange__ordered_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/frankfurt_nacht__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/frankfurt_nacht__burkes_spectra6_auto.bin index 0100f6f..2699792 100644 Binary files a/packages/rust/core/tests/fixtures/references/frankfurt_nacht__burkes_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/frankfurt_nacht__burkes_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/frankfurt_nacht__floyd_steinberg_mono_raw.bin b/packages/rust/core/tests/fixtures/references/frankfurt_nacht__floyd_steinberg_mono_raw.bin index c1f2b2c..2945549 100644 Binary files a/packages/rust/core/tests/fixtures/references/frankfurt_nacht__floyd_steinberg_mono_raw.bin and b/packages/rust/core/tests/fixtures/references/frankfurt_nacht__floyd_steinberg_mono_raw.bin differ diff --git a/packages/rust/core/tests/fixtures/references/frankfurt_nacht__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/frankfurt_nacht__ordered_spectra6_auto.bin index d14b28f..0b00851 100644 Binary files a/packages/rust/core/tests/fixtures/references/frankfurt_nacht__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/frankfurt_nacht__ordered_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/olympiapark__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/olympiapark__burkes_spectra6_auto.bin index 37ba17c..eb17991 100644 Binary files a/packages/rust/core/tests/fixtures/references/olympiapark__burkes_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/olympiapark__burkes_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/olympiapark__floyd_steinberg_mono_raw.bin b/packages/rust/core/tests/fixtures/references/olympiapark__floyd_steinberg_mono_raw.bin index 6dca27a..a2d969a 100644 Binary files a/packages/rust/core/tests/fixtures/references/olympiapark__floyd_steinberg_mono_raw.bin and b/packages/rust/core/tests/fixtures/references/olympiapark__floyd_steinberg_mono_raw.bin differ diff --git a/packages/rust/core/tests/fixtures/references/olympiapark__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/olympiapark__ordered_spectra6_auto.bin index 377b417..7e2bb81 100644 Binary files a/packages/rust/core/tests/fixtures/references/olympiapark__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/olympiapark__ordered_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/river__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/river__burkes_spectra6_auto.bin index 4661d6c..a879a16 100644 --- a/packages/rust/core/tests/fixtures/references/river__burkes_spectra6_auto.bin +++ b/packages/rust/core/tests/fixtures/references/river__burkes_spectra6_auto.bin @@ -1 +1 @@ - \ No newline at end of file + \ No newline at end of file diff --git a/packages/rust/core/tests/fixtures/references/river__floyd_steinberg_mono_raw.bin b/packages/rust/core/tests/fixtures/references/river__floyd_steinberg_mono_raw.bin index 99f3e52..1d2a64c 100644 Binary files a/packages/rust/core/tests/fixtures/references/river__floyd_steinberg_mono_raw.bin and b/packages/rust/core/tests/fixtures/references/river__floyd_steinberg_mono_raw.bin differ diff --git a/packages/rust/core/tests/fixtures/references/river__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/river__ordered_spectra6_auto.bin index c5b2070..7454af2 100644 Binary files a/packages/rust/core/tests/fixtures/references/river__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/river__ordered_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/seeed_opendisplay__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/seeed_opendisplay__burkes_spectra6_auto.bin index cc54a20..f983e9e 100644 Binary files a/packages/rust/core/tests/fixtures/references/seeed_opendisplay__burkes_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/seeed_opendisplay__burkes_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/seeed_opendisplay__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/seeed_opendisplay__ordered_spectra6_auto.bin index df0e64d..5d5b5ea 100644 Binary files a/packages/rust/core/tests/fixtures/references/seeed_opendisplay__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/seeed_opendisplay__ordered_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/ubahn_station__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/ubahn_station__burkes_spectra6_auto.bin index c9771e4..456704b 100644 Binary files a/packages/rust/core/tests/fixtures/references/ubahn_station__burkes_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/ubahn_station__burkes_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/ubahn_station__floyd_steinberg_mono_raw.bin b/packages/rust/core/tests/fixtures/references/ubahn_station__floyd_steinberg_mono_raw.bin index 33e0718..ed299a9 100644 Binary files a/packages/rust/core/tests/fixtures/references/ubahn_station__floyd_steinberg_mono_raw.bin and b/packages/rust/core/tests/fixtures/references/ubahn_station__floyd_steinberg_mono_raw.bin differ diff --git a/packages/rust/core/tests/fixtures/references/ubahn_station__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/ubahn_station__ordered_spectra6_auto.bin index 8aade53..1731b99 100644 Binary files a/packages/rust/core/tests/fixtures/references/ubahn_station__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/ubahn_station__ordered_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/unicorn__burkes_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/unicorn__burkes_spectra6_auto.bin index f6728dd..6a89750 100644 Binary files a/packages/rust/core/tests/fixtures/references/unicorn__burkes_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/unicorn__burkes_spectra6_auto.bin differ diff --git a/packages/rust/core/tests/fixtures/references/unicorn__floyd_steinberg_mono_raw.bin b/packages/rust/core/tests/fixtures/references/unicorn__floyd_steinberg_mono_raw.bin index b9b6e47..03b3ae8 100644 Binary files a/packages/rust/core/tests/fixtures/references/unicorn__floyd_steinberg_mono_raw.bin and b/packages/rust/core/tests/fixtures/references/unicorn__floyd_steinberg_mono_raw.bin differ diff --git a/packages/rust/core/tests/fixtures/references/unicorn__ordered_spectra6_auto.bin b/packages/rust/core/tests/fixtures/references/unicorn__ordered_spectra6_auto.bin index 43aadc7..d8ac04f 100644 Binary files a/packages/rust/core/tests/fixtures/references/unicorn__ordered_spectra6_auto.bin and b/packages/rust/core/tests/fixtures/references/unicorn__ordered_spectra6_auto.bin differ