Image Trace: a palette of the picture's own colours

Color and Grayscale chose their colours from every pixel, so the blends along anti-aliased edges took palette
entries and traced as slivers, and with more colours asked for than the art has, k-means split flat colours
along their noise and the halves traced as specks all over each other.

- The palette is made from the picture's own colours: a pixel that lies between the pixels on opposite sides
  of it is a blend and counts for nothing (lines one pixel wide and noise are no blends).
- Colours that noise made of one colour merge: close colours whose pixels are mixed through each other. The
  bands of a gradient keep to themselves and stay apart.
- A blend pixel whose nearest palette colour is a third colour (halfway from navy to yellow lies an olive that
  is nearer to a green elsewhere) goes to the nearer of the two colours it lies between.
- Boxes are cut where their two sides differ most (Otsu's criterion) instead of at the median, which split the
  largest colour of flat art along its noise and left two smaller colours in one box.
- Colours are compared in Oklab when pixels are assigned and in k-means.

On each preset's own inputs, all sizes and variants (docs/image-trace-quality.md), before -> after:
- 16 Colors: 69.9 -> 4.2 paths, 1743 -> 66 anchors, outlines 5.08 -> 0.13 px off on average (largest 85.9 ->
  1.09), kinks 1376 -> 4.4, paths in a colour the art doesn't have 5.9 -> 0
- 6 Colors: 21.2 -> 4.0 paths, 4.44 -> 0.14 px off; 3 Colors: colour error 4.86 -> 0.43, 0.33 -> 0.09 px off
- High Fidelity Photo: 3636 -> 639 paths, 81,597 -> 25,604 anchors, colour error 2.50 -> 2.37; a 1 megapixel
  photo-like image 1.98 -> 0.42 s (51,389 -> 4,298 paths), a 4 megapixel one traced (27,583 paths) where it was
  refused
- Low Fidelity Photo: 128 -> 55 paths, colour error 2.59 -> 2.67; Shades of Gray: 91 -> 73 paths, 9.29 -> 9.33
This commit is contained in:
bflatastic
2026-10-13 03:48:30 +07:00
parent a891530ef2
commit 4204ce36db
4 changed files with 556 additions and 223 deletions
+322 -78
View File
@@ -1,4 +1,5 @@
//! Colour reduction and noise removal.
//! Colour reduction and noise removal. Colours are compared in Oklab, where equal distances look
//! equally different.
use crate::{Mode, Palette, Raster, TraceParams};
@@ -37,12 +38,18 @@ pub fn quantize(img: &Raster, params: &TraceParams) -> Quantized {
Quantized { labels, palette: vec![[255, 255, 255], [0, 0, 0]] }
}
Mode::Grayscale => {
let mut hist = [0u64; 256];
for p in img.rgba.as_chunks::<4>().0 {
if p[3] >= 128 {
hist[luma([p[0], p[1], p[2], p[3]]) as usize] += 1;
let own = own_colours(img);
let histogram = |only_own: bool| {
let mut hist = [0u64; 256];
for (p, &f) in img.rgba.as_chunks::<4>().0.iter().zip(&own) {
if p[3] >= 128 && (f || !only_own) {
hist[luma([p[0], p[1], p[2], p[3]]) as usize] += 1;
}
}
}
hist
};
let (every, calm) = (histogram(false), histogram(true));
let hist = if palette_from_own(calm.iter().sum::<u64>() as f64, every.iter().sum::<u64>() as f64) { calm } else { every };
let levels = kmeans_1d(&hist, params.colors.clamp(2, 256) as usize);
let mut lut = [0u16; 256];
for (v, slot) in lut.iter_mut().enumerate() {
@@ -53,36 +60,251 @@ pub fn quantize(img: &Raster, params: &TraceParams) -> Quantized {
let g = l.round().clamp(0.0, 255.0) as u8;
[g, g, g]
});
compact(labels, palette.collect())
settled(merge_mixed(compact(labels, palette.collect()), img.width as usize, img.height as usize), img, &own)
}
// Logo mode traces without quantising; asked directly, it quantises like Color.
Mode::Color | Mode::Logo => {
// 15-bit colour histogram with per-bin sums (exact bin means), fitted to the image's
// range per channel so that close tones (a dark, low-contrast photo) stay apart.
let hk = HistKey::new(img);
let mut bins: Vec<Bin> = vec![Bin::default(); 1 << 15];
for p in img.rgba.as_chunks::<4>().0 {
let own = own_colours(img);
let (mut every, mut calm) = (vec![Bin::default(); 1 << 15], vec![Bin::default(); 1 << 15]);
for (p, &f) in img.rgba.as_chunks::<4>().0.iter().zip(&own) {
if p[3] >= 128 {
let b = &mut bins[hk.key(p[0], p[1], p[2])];
b.n += 1;
b.sum[0] += p[0] as u64;
b.sum[1] += p[1] as u64;
b.sum[2] += p[2] as u64;
let k = hk.key(p[0], p[1], p[2]);
every[k].add(p);
if f {
calm[k].add(p);
}
}
}
let used: Vec<Used> = bins.iter().enumerate().filter(|(_, b)| b.n > 0).map(|(i, b)| (i, b.mean(), b.n as f64)).collect();
let centers = color_palette(&used, params);
let mut lut = vec![0u16; 1 << 15];
for (i, c, _) in &used {
lut[*i] = nearest_rgb(&centers, *c) as u16;
}
let labels = px.map(|p| if p[3] < 128 { TRANSPARENT } else { lut[hk.key(p[0], p[1], p[2])] }).collect();
let palette = centers.iter().map(|c| [c[0].round() as u8, c[1].round() as u8, c[2].round() as u8]).collect();
compact(labels, palette)
let used =
|bins: &[Bin]| -> Vec<Used> { bins.iter().enumerate().filter(|(_, b)| b.n > 0).map(|(i, b)| (i, b.mean(), b.n as f64)).collect() };
let (every, calm) = (used(&every), used(&calm));
let count = |u: &[Used]| u.iter().map(|b| b.2).sum::<f64>();
let from = if palette_from_own(count(&calm), count(&every)) { &calm } else { &every };
// The picture with `colors` colours asked for; fewer come back when it has fewer, or
// when noise made several of one (see `merge_mixed`).
let reduce = |colors: u32| {
let centers = color_palette(from, &TraceParams { colors, ..params.clone() });
let mut lut = vec![0u16; 1 << 15];
let center_labs = labs(&centers);
for (i, c, _) in &every {
lut[*i] = nearest_lab(&center_labs, oklab(*c)) as u16;
}
let labels =
img.rgba.as_chunks::<4>().0.iter().map(|p| if p[3] < 128 { TRANSPARENT } else { lut[hk.key(p[0], p[1], p[2])] }).collect();
let palette = centers.iter().map(|c| [c[0].round() as u8, c[1].round() as u8, c[2].round() as u8]).collect();
let q = compact(labels, palette);
// A swatch library's colours are used as they are.
if params.palette == Palette::DocumentLibrary && !params.swatches.is_empty() {
q
} else {
merge_mixed(q, img.width as usize, img.height as usize)
}
};
let q = reduce(params.colors);
settled(q, img, &own)
}
}
}
/// The opaque pixels whose colour is the picture's own, not a blend: a pixel is a blend when it
/// lies between the two pixels on opposite sides of it (left and right, above and below, or
/// across a diagonal) that differ by more than 20 levels, which is what anti-aliasing and blur
/// leave along an edge. A line one pixel wide is no blend (the pixels on its two sides are alike),
/// nor is noise.
fn own_colours(img: &Raster) -> Vec<bool> {
let (w, h) = (img.width as isize, img.height as isize);
let px = img.rgba.as_chunks::<4>().0;
let at = |x: isize, y: isize| {
let p = (x >= 0 && y >= 0 && x < w && y < h).then(|| px.get((y * w + x) as usize)).flatten()?;
(p[3] >= 128).then(|| [f64::from(p[0]), f64::from(p[1]), f64::from(p[2])])
};
let mut out = Vec::with_capacity(px.len());
for y in 0..h {
for x in 0..w {
let Some(p) = at(x, y) else {
out.push(false);
continue;
};
let blend = [(1, 0), (0, 1), (1, 1), (1, -1)].into_iter().any(|(dx, dy)| {
let (Some(a), Some(b)) = (at(x - dx, y - dy), at(x + dx, y + dy)) else { return false };
let d = [0, 1, 2].map(|k| b[k] - a[k]);
let len2: f64 = d.iter().map(|v| v * v).sum();
if len2 <= 400.0 {
return false;
}
let t = (0..3).map(|k| (p[k] - a[k]) * d[k]).sum::<f64>() / len2;
let off2: f64 = (0..3).map(|k| (p[k] - a[k] - t * d[k]).powi(2)).sum();
t > 0.06 && t < 0.94 && off2 < len2 / 16.0
});
out.push(!blend);
}
}
out
}
/// Whether the palette is made from the picture's own colours alone (`own` of `all` opaque
/// pixels), so that the blends along edges take none of its colours. A picture with next to none
/// is all soft edges: there every pixel counts.
fn palette_from_own(own: f64, all: f64) -> bool {
own >= 0.05 * all && own >= 16.0
}
/// sRGB (0–255) to Oklab × 100 (Björn Ottosson, "A perceptual color space for image processing",
/// 2020): equal distances look equally different, and 2 or 3 is about the least one sees.
fn oklab(c: [f64; 3]) -> [f64; 3] {
let lin = c.map(|v| {
let v = (v / 255.0).clamp(0.0, 1.0);
if v <= 0.04045 { v / 12.92 } else { ((v + 0.055) / 1.055).powf(2.4) }
});
let l = (0.412_221_470_8 * lin[0] + 0.536_332_536_3 * lin[1] + 0.051_445_992_9 * lin[2]).cbrt();
let m = (0.211_903_498_2 * lin[0] + 0.680_699_545_1 * lin[1] + 0.107_396_956_6 * lin[2]).cbrt();
let s = (0.088_302_461_9 * lin[0] + 0.281_718_837_6 * lin[1] + 0.629_978_700_5 * lin[2]).cbrt();
[
100.0 * (0.210_454_255_3 * l + 0.793_617_785 * m - 0.004_072_046_8 * s),
100.0 * (1.977_998_495_1 * l - 2.428_592_205 * m + 0.450_593_709_9 * s),
100.0 * (0.025_904_037_1 * l + 0.782_771_766_2 * m - 0.808_675_766 * s),
]
}
/// `q` with the colours that noise made of one colour merged, largest mix first.
///
/// Asked for more colours than a picture has, k-means splits a flat colour along its noise, and
/// the two halves would trace as specks of one all over the other. Such halves are close (under 30
/// levels apart: noise is about as strong in a dark colour as in a light one, so the distance is
/// taken in RGB) and thoroughly mixed: the pixels of the lesser one have, on average, more than
/// half a neighbour of the other each; or they are closer still (under 12 levels) and mixed in
/// patches (compression blocks). The bands of a gradient are close too, but each keeps to itself
/// and they stay apart.
fn merge_mixed(q: Quantized, w: usize, h: usize) -> Quantized {
let k = q.palette.len();
if k < 2 || q.labels.len() != w * h {
return q;
}
// How many pixels each colour has, and how many neighbouring pairs each pair of colours.
let mut count = vec![0.0f64; k];
let mut touch = vec![0.0f64; k * k];
let label = |i: usize| q.labels.get(i).map(|&l| l as usize).filter(|&l| l < k);
for i in 0..w * h {
let Some(a) = label(i) else { continue };
count[a] += 1.0;
for j in [(i % w + 1 < w).then_some(i + 1), (i + w < w * h).then_some(i + w)].into_iter().flatten() {
if let Some(b) = label(j).filter(|&b| b != a) {
touch[a * k + b] += 1.0;
touch[b * k + a] += 1.0;
}
}
}
let mut color: Vec<[f64; 3]> = q.palette.iter().map(|c| c.map(f64::from)).collect();
// Which colour each has merged into (itself when none).
let mut into: Vec<usize> = (0..k).collect();
loop {
let mut best: Option<(f64, usize, usize)> = None;
for a in 0..k {
for b in a + 1..k {
if into[a] != a || into[b] != b || count[a] <= 0.0 || count[b] <= 0.0 {
continue;
}
let mixed = touch[a * k + b] / count[a].min(count[b]);
let apart = (0..3).map(|c| (color[a][c] - color[b][c]).powi(2)).sum::<f64>().sqrt();
if (mixed > 0.5 && apart < 30.0 || mixed > 0.15 && apart < 12.0) && best.is_none_or(|x| mixed > x.0) {
best = Some((mixed, a, b));
}
}
}
let Some((_, a, b)) = best else { break };
color[a] = [0, 1, 2].map(|c| (color[a][c] * count[a] + color[b][c] * count[b]) / (count[a] + count[b]));
count[a] += count[b];
count[b] = 0.0;
into[b] = a;
for c in 0..k {
if c != a && c != b {
touch[a * k + c] += touch[b * k + c];
touch[c * k + a] = touch[a * k + c];
}
}
}
if into.iter().enumerate().all(|(i, &t)| i == t) {
return q;
}
// A colour may have merged into one that merged again: follow each to where it ended up.
for i in 0..k {
let mut to = into[i];
while into[to] != to {
to = into[to];
}
into[i] = to;
}
let labels = q.labels.iter().map(|&l| into.get(l as usize).map_or(l, |&t| t as u16)).collect();
compact(labels, color.iter().map(|c| [c[0].round() as u8, c[1].round() as u8, c[2].round() as u8]).collect())
}
/// `q` with its blend pixels settled, and the colours that leaves unused dropped.
///
/// A pixel on an edge is a blend of the colours on its two sides, and the palette colour nearest
/// to a blend can be a third colour altogether (halfway from navy to yellow lies an olive that is
/// nearer to a green elsewhere in the picture), which would trace as a sliver of that colour along
/// the edge. A blend pixel (see [`own_colours`]) whose colour is none of those of the picture's
/// own pixels within two pixels of it takes the side it is nearer to of the pair of those colours
/// it is best a blend of, when that explains it twice as well as its own colour does. With one
/// colour around it, it keeps its own.
fn settled(mut q: Quantized, img: &Raster, own: &[bool]) -> Quantized {
let (w, h) = (img.width as usize, img.height as usize);
if q.labels.len() != w * h || own.len() != w * h {
return q;
}
let before = q.labels.clone();
let color = |l: u16| q.palette.get(l as usize).map(|c| c.map(f64::from));
let dist = |a: [f64; 3], b: [f64; 3]| (0..3).map(|k| (a[k] - b[k]).powi(2)).sum::<f64>().sqrt();
let mut changed = false;
for i in 0..w * h {
let mine_l = before[i];
if own[i] || mine_l == TRANSPARENT {
continue;
}
let (x, y) = (i % w, i / w);
// The colours of the picture's own pixels within two pixels (seldom more than three).
let mut near: Vec<u16> = Vec::new();
for yy in y.saturating_sub(2)..=(y + 2).min(h - 1) {
for xx in x.saturating_sub(2)..=(x + 2).min(w - 1) {
let j = yy * w + xx;
if own[j] && before[j] != TRANSPARENT && !near.contains(&before[j]) {
near.push(before[j]);
}
}
}
let (Some(p), Some(mine)) = (img.rgba.get(i * 4..i * 4 + 3).map(|p| [f64::from(p[0]), f64::from(p[1]), f64::from(p[2])]), color(mine_l))
else {
continue;
};
if near.len() < 2 || near.contains(&mine_l) {
continue;
}
let mut best: Option<(f64, u16)> = None;
for (n, &a) in near.iter().enumerate() {
for &b in &near[n + 1..] {
let (Some(ca), Some(cb)) = (color(a), color(b)) else { continue };
let d = [0, 1, 2].map(|k| cb[k] - ca[k]);
let len2: f64 = d.iter().map(|v| v * v).sum();
let t = if len2 > 0.0 { ((0..3).map(|k| (p[k] - ca[k]) * d[k]).sum::<f64>() / len2).clamp(0.0, 1.0) } else { 0.0 };
let off = dist(p, [0, 1, 2].map(|k| ca[k] + t * d[k]));
if best.is_none_or(|x| off < x.0) {
best = Some((off, if t < 0.5 { a } else { b }));
}
}
}
if let Some((off, to)) = best
&& off < 0.5 * dist(p, mine)
{
q.labels[i] = to;
changed = true;
}
}
if changed { compact(q.labels, q.palette) } else { q }
}
/// Drop unused palette entries and merge duplicates.
fn compact(mut labels: Vec<u16>, palette: Vec<[u8; 3]>) -> Quantized {
let mut used = vec![false; palette.len()];
@@ -120,6 +342,12 @@ struct Bin {
}
impl Bin {
fn add(&mut self, p: &[u8; 4]) {
self.n += 1;
for (sum, v) in self.sum.iter_mut().zip(p) {
*sum += u64::from(*v);
}
}
fn mean(&self) -> [f64; 3] {
let n = self.n.max(1) as f64;
[self.sum[0] as f64 / n, self.sum[1] as f64 / n, self.sum[2] as f64 / n]
@@ -208,16 +436,11 @@ fn kmeans_1d(hist: &[u64; 256], k: usize) -> Vec<f64> {
levels
}
fn dist2(a: [f64; 3], b: [f64; 3]) -> f64 {
// Slightly perceptual weights.
let (dr, dg, db) = (a[0] - b[0], a[1] - b[1], a[2] - b[2]);
2.0 * dr * dr + 4.0 * dg * dg + 3.0 * db * db
}
fn nearest_rgb(centers: &[[f64; 3]], c: [f64; 3]) -> usize {
/// The index of the colour of `labs` (Oklab) nearest to `c` (Oklab): the one that looks most like it.
fn nearest_lab(labs: &[[f64; 3]], c: [f64; 3]) -> usize {
let mut best = (f64::INFINITY, 0);
for (i, &m) in centers.iter().enumerate() {
let d = dist2(m, c);
for (i, m) in labs.iter().enumerate() {
let d = (m[0] - c[0]).powi(2) + (m[1] - c[1]).powi(2) + (m[2] - c[2]).powi(2);
if d < best.0 {
best = (d, i);
}
@@ -225,6 +448,10 @@ fn nearest_rgb(centers: &[[f64; 3]], c: [f64; 3]) -> usize {
best.1
}
fn labs(colors: &[[f64; 3]]) -> Vec<[f64; 3]> {
colors.iter().map(|c| oklab(*c)).collect()
}
/// The most colours Full Tone (and Automatic on continuous-tone images) makes.
pub const FULL_TONE_MAX: usize = 256;
@@ -295,8 +522,9 @@ fn library_colors(used: &[Used], swatches: &[[u8; 3]], k: usize) -> Vec<[f64; 3]
}
}
let mut weight = vec![0.0; all.len()];
let all_labs = labs(&all);
for (_, c, n) in used {
if let Some(w) = weight.get_mut(nearest_rgb(&all, *c)) {
if let Some(w) = weight.get_mut(nearest_lab(&all_labs, oklab(*c))) {
*w += n;
}
}
@@ -321,60 +549,73 @@ fn kmeans_rgb(used: &[Used], k: usize) -> Vec<[f64; 3]> {
refine(used, centers(&median_cut(used, k, 0.0)))
}
/// Median cut: repeatedly split the box with the largest weighted extent, while there are fewer
/// than `k` boxes and some box spans more than `spread` levels of a channel.
/// Divisive clustering in the manner of median cut, but each box is cut where that leaves its two
/// sides most alike within themselves (the cut, along one channel, with the largest between-side
/// variance: Otsu's criterion), and the box cut next is the one whose cut gains most. A cut at the
/// median, as in classic median cut, splits the largest colour of flat art in two along its noise
/// and leaves two smaller colours in one box. Boxes are cut while there are fewer than `k` and
/// some box spans more than `spread` levels of a channel.
fn median_cut(used: &[Used], k: usize, spread: f64) -> Vec<ColorBox> {
if used.is_empty() {
return vec![ColorBox { mean: [255.0; 3], weight: 0.0, peak: 0.0 }];
}
let mut boxes: Vec<Vec<usize>> = vec![(0..used.len()).collect()];
while boxes.len() < k {
let mut best: Option<(f64, usize, usize)> = None;
for (bi, b) in boxes.iter().enumerate() {
if b.len() < 2 {
continue;
}
let (mut lo, mut hi) = ([f64::INFINITY; 3], [f64::NEG_INFINITY; 3]);
let mut wsum = 0.0;
for &i in b {
for ch in 0..3 {
lo[ch] = lo[ch].min(used[i].1[ch]);
hi[ch] = hi[ch].max(used[i].1[ch]);
/// A box's bins and its best cut: (gain, bins sorted along the cut's channel, bins on the first side).
struct Cut {
bins: Vec<usize>,
best: Option<(f64, Vec<usize>, usize)>,
}
let plan = |bins: Vec<usize>| {
let mut best: Option<(f64, Vec<usize>, usize)> = None;
let widest = (0..3)
.map(|ch| {
let v = bins.iter().map(|&i| used[i].1[ch]);
v.clone().fold(f64::NEG_INFINITY, f64::max) - v.fold(f64::INFINITY, f64::min)
})
.fold(0.0, f64::max);
if bins.len() >= 2 && widest > spread.max(0.0) {
for ch in 0..3 {
let mut sorted = bins.clone();
sorted.sort_by(|x, y| used[*x].1[ch].total_cmp(&used[*y].1[ch]));
let (total, sum): (f64, f64) = sorted.iter().fold((0.0, 0.0), |a, &i| (a.0 + used[i].2, a.1 + used[i].1[ch] * used[i].2));
let (mut w, mut s) = (0.0, 0.0);
let mut cut: Option<(f64, usize)> = None;
for (j, &i) in sorted.iter().enumerate().take(sorted.len() - 1) {
w += used[i].2;
s += used[i].1[ch] * used[i].2;
let rest = total - w;
if w > 0.0 && rest > 0.0 {
let gain = w * rest / total * (s / w - (sum - s) / rest).powi(2);
if cut.is_none_or(|c| gain > c.0) {
cut = Some((gain, j + 1));
}
}
}
if let Some((gain, at)) = cut
&& best.as_ref().is_none_or(|b| gain > b.0)
{
best = Some((gain, sorted, at));
}
wsum += used[i].2;
}
let (ch, range) = (0..3).map(|c| (c, hi[c] - lo[c])).fold((0, -1.0), |a, b| if b.1 > a.1 { b } else { a });
if range <= spread.max(0.0) {
continue;
}
let score = range * wsum.sqrt();
if best.is_none_or(|x| score > x.0) {
best = Some((score, bi, ch));
}
}
let Some((_, bi, ch)) = best else { break };
let mut b = boxes.swap_remove(bi);
b.sort_by(|x, y| used[*x].1[ch].total_cmp(&used[*y].1[ch]));
let total: f64 = b.iter().map(|&i| used[i].2).sum();
let mut acc = 0.0;
let mut cut = 1;
for (j, &i) in b.iter().enumerate() {
acc += used[i].2;
if acc >= total / 2.0 {
cut = (j + 1).clamp(1, b.len() - 1);
break;
}
}
let tail = b.split_off(cut);
boxes.push(b);
boxes.push(tail);
Cut { bins, best }
};
let mut boxes = vec![plan((0..used.len()).collect())];
while boxes.len() < k {
let gain = |c: &Cut| c.best.as_ref().map_or(f64::NEG_INFINITY, |b| b.0);
let Some(bi) = (0..boxes.len()).filter(|&i| boxes[i].best.is_some()).max_by(|&a, &b| gain(&boxes[a]).total_cmp(&gain(&boxes[b]))) else {
break;
};
let Some((_, mut sorted, at)) = boxes.swap_remove(bi).best else { break };
let tail = sorted.split_off(at);
boxes.push(plan(sorted));
boxes.push(plan(tail));
}
boxes
.iter()
.map(|b| {
let mut s = [0.0; 3];
let (mut w, mut peak) = (0.0, 0.0f64);
for &i in b {
for &i in &b.bins {
for (acc, v) in s.iter_mut().zip(used[i].1) {
*acc += v * used[i].2;
}
@@ -386,13 +627,16 @@ fn median_cut(used: &[Used], k: usize, spread: f64) -> Vec<ColorBox> {
.collect()
}
/// Lloyd refinement of `centers` over the bins `used`.
/// Lloyd refinement of `centers` over the bins `used`: each bin goes to the centre that looks most
/// like it (nearest in Oklab), and each centre becomes the mean colour of its bins.
fn refine(used: &[Used], mut centers: Vec<[f64; 3]>) -> Vec<[f64; 3]> {
let bins = labs(&used.iter().map(|u| u.1).collect::<Vec<_>>());
for _ in 0..12 {
let center_labs = labs(&centers);
let mut sum = vec![[0.0; 3]; centers.len()];
let mut cnt = vec![0.0; centers.len()];
for (_, c, n) in used {
let j = nearest_rgb(&centers, *c);
for ((_, c, n), lab) in used.iter().zip(&bins) {
let j = nearest_lab(&center_labs, *lab);
for ch in 0..3 {
sum[j][ch] += c[ch] * n;
}
@@ -402,7 +646,7 @@ fn refine(used: &[Used], mut centers: Vec<[f64; 3]>) -> Vec<[f64; 3]> {
for j in 0..centers.len() {
if cnt[j] > 0.0 {
let m = [sum[j][0] / cnt[j], sum[j][1] / cnt[j], sum[j][2] / cnt[j]];
moved = moved.max(dist2(m, centers[j]));
moved = moved.max((0..3).map(|ch| (m[ch] - centers[j][ch]).powi(2)).sum());
centers[j] = m;
}
}