mirror of
https://github.com/storytold/vectorcraft.git
synced 2026-10-24 09:46:59 -06:00
M12: use exact crop distances and Guo–Hall stroke thinning
Compute exact squared Euclidean distances per blob with separable parabola lower envelopes. Thin with bounded Guo–Hall active sets and retain diagonal ends during staircase cleanup. Keep the whole-blob width rule and existing graph, ink, cap and fitting rules, with checked access in the moved module. Verify EDT against brute force, active sets against plain rounds, and both diagonal slopes for length and connectivity. Record the centre-lines rerun and its remaining limitations without changing feature scores. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Co-Authored-By: OpenAI <noreply@openai.com>
This commit is contained in:
andeecollard
committed by
bflatastic
bflatastic
co-authored by
Claude Opus 5.5
OpenAI
Claude Opus 5.5
OpenAI
parent
c9f0ed14f0
commit
abb7c1351b
@@ -0,0 +1,109 @@
|
||||
//! Exact squared Euclidean distance, Felzenszwalb & Huttenlocher, "Distance Transforms of
|
||||
//! Sampled Functions", §3. Separable lower envelopes of parabolas, linear time and crop storage.
|
||||
|
||||
const FAR: f64 = 1e20;
|
||||
|
||||
/// Transform one row or column. Each envelope entry is (sample, left intersection).
|
||||
fn transform(f: &[f64], out: &mut [f64], envelope: &mut Vec<(usize, f64)>) -> Option<()> {
|
||||
envelope.clear();
|
||||
for (q, &fq) in f.iter().enumerate() {
|
||||
let mut left = -FAR;
|
||||
while let Some(&(v, z)) = envelope.last() {
|
||||
let (qf, vf) = (q as f64, v as f64);
|
||||
left = ((fq + qf * qf) - (*f.get(v)? + vf * vf)) / (2.0 * (qf - vf));
|
||||
if left > z {
|
||||
break;
|
||||
}
|
||||
envelope.pop();
|
||||
}
|
||||
if envelope.is_empty() {
|
||||
left = -FAR;
|
||||
}
|
||||
envelope.push((q, left));
|
||||
}
|
||||
let mut k = 0;
|
||||
for (q, d) in out.iter_mut().enumerate() {
|
||||
while envelope.get(k + 1).is_some_and(|&(_, z)| z < q as f64) {
|
||||
k += 1;
|
||||
}
|
||||
let &(v, _) = envelope.get(k)?;
|
||||
*d = (q as f64 - v as f64).powi(2) + *f.get(v)?;
|
||||
}
|
||||
Some(())
|
||||
}
|
||||
|
||||
/// Squared centre-to-paper distances for a validated, paper-bordered component crop.
|
||||
/// Border paper also represents the space outside the raster. Returns None for invalid sizes
|
||||
/// or a non-paper border, rather than accepting an unbounded allocation or an absent seed.
|
||||
pub(super) fn squared(mask: &[bool], w: usize, h: usize) -> Option<Vec<f32>> {
|
||||
let size = super::crop_size(w, h)?;
|
||||
if mask.len() != size || mask.iter().enumerate().any(|(i, &on)| on && (i < w || i >= size - w || i % w == 0 || i % w + 1 == w)) {
|
||||
return None;
|
||||
}
|
||||
let mut grid: Vec<f64> = mask.iter().map(|&on| if on { FAR } else { 0.0 }).collect();
|
||||
let n = w.max(h);
|
||||
let (mut f, mut d, mut envelope) = (vec![0.0; n], vec![0.0; n], Vec::with_capacity(n));
|
||||
for x in 0..w {
|
||||
for (y, v) in f.get_mut(..h)?.iter_mut().enumerate() {
|
||||
*v = *grid.get(y * w + x)?;
|
||||
}
|
||||
transform(f.get(..h)?, d.get_mut(..h)?, &mut envelope)?;
|
||||
for (y, &v) in d.get(..h)?.iter().enumerate() {
|
||||
*grid.get_mut(y * w + x)? = v;
|
||||
}
|
||||
}
|
||||
for row in grid.chunks_exact_mut(w) {
|
||||
for (v, &g) in f.get_mut(..w)?.iter_mut().zip(row.iter()) {
|
||||
*v = g;
|
||||
}
|
||||
transform(f.get(..w)?, d.get_mut(..w)?, &mut envelope)?;
|
||||
for (g, &v) in row.iter_mut().zip(d.get(..w)?) {
|
||||
*g = v;
|
||||
}
|
||||
}
|
||||
Some(grid.into_iter().map(|v| v as f32).collect())
|
||||
}
|
||||
|
||||
#[cfg(test)]
|
||||
mod tests {
|
||||
use super::*;
|
||||
|
||||
#[test]
|
||||
fn edt_equals_brute_force_on_random_masks() {
|
||||
let mut seed = 19u64;
|
||||
let mut rnd = || {
|
||||
seed = seed.wrapping_mul(6364136223846793005).wrapping_add(1);
|
||||
(seed >> 33) as usize
|
||||
};
|
||||
for case in 0..200 {
|
||||
let (w, h) = (3 + rnd() % 24, 3 + rnd() % 24);
|
||||
let mask: Vec<bool> =
|
||||
(0..w * h).map(|i| i >= w && i < w * (h - 1) && i % w != 0 && i % w + 1 < w && (case % 5 == 0 || rnd() % 4 != 0)).collect();
|
||||
let d = squared(&mask, w, h).unwrap();
|
||||
for (i, &actual) in d.iter().enumerate() {
|
||||
let expected = mask
|
||||
.iter()
|
||||
.enumerate()
|
||||
.filter(|&(_, &on)| !on)
|
||||
.map(|(j, _)| {
|
||||
let dx = (i % w).abs_diff(j % w);
|
||||
let dy = (i / w).abs_diff(j / w);
|
||||
(dx * dx + dy * dy) as f32
|
||||
})
|
||||
.min_by(f32::total_cmp)
|
||||
.unwrap();
|
||||
assert_eq!(actual, expected, "case {case}, pixel {i}");
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn invalid_crops_are_refused() {
|
||||
assert!(squared(&[], 0, 0).is_none());
|
||||
assert!(squared(&[], usize::MAX, 3).is_none());
|
||||
assert!(squared(&[false; 9], 4, 3).is_none());
|
||||
assert!(squared(&[true; 9], 3, 3).is_none());
|
||||
assert!(squared(&[], crate::MAX_PIXELS as usize + 3, 3).is_none());
|
||||
assert_eq!(squared(&[false; 9], 3, 3), Some(vec![0.0; 9]));
|
||||
}
|
||||
}
|
||||
@@ -1,13 +1,12 @@
|
||||
//! Create Strokes: the parts of a colour layer no wider than the Stroke width, traced as stroked
|
||||
//! centre lines instead of filled outlines.
|
||||
//!
|
||||
//! 1. **Distance**: a chamfer (3-4) distance transform gives each pixel its distance to the nearest
|
||||
//! pixel outside the layer.
|
||||
//! 1. **Distance**: an exact Euclidean distance transform in each blob crop gives each pixel
|
||||
//! its distance to the nearest pixel outside the blob.
|
||||
//! 2. **Lines**: an 8-connected component whose widest point (twice its largest distance, less
|
||||
//! one) is within the stroke width is a line; wider ones stay filled areas.
|
||||
//! 3. **Skeleton**: the line is thinned to one pixel (Zhang–Suen, applied pixel by pixel so it
|
||||
//! keeps its topology and two-pixel-wide lines don't vanish), then staircase corners that
|
||||
//! don't connect anything are dropped.
|
||||
//! 3. **Skeleton**: Guo–Hall thins the crop to one pixel, preserving diagonal ends, then
|
||||
//! staircase corners that don't connect anything are dropped.
|
||||
//! 4. **Graph**: runs of skeleton pixels between ends and junctions (touching junction pixels
|
||||
//! merged into one point) become polylines through the pixel centres; rings without ends
|
||||
//! become closed ones. Short spurs that thinning leaves at junctions are pruned and the runs
|
||||
@@ -22,6 +21,9 @@
|
||||
//! 6. **Curves**: each polyline is fitted as outlines are ([`fit::fit_piece`]): a straight line
|
||||
//! is two anchors, corners stay corners, the rest is as few curves as Paths allows.
|
||||
|
||||
mod distance;
|
||||
mod thin;
|
||||
|
||||
use std::collections::HashSet;
|
||||
|
||||
use vectorcraft_geom::{PathData, Point};
|
||||
@@ -55,38 +57,46 @@ pub(crate) fn lines(
|
||||
ink: [f64; 3],
|
||||
style: &Style,
|
||||
) -> (Vec<bool>, Vec<Line>) {
|
||||
let mut thin = vec![false; labels.len()];
|
||||
let mut out = vec![];
|
||||
if w == 0 || h == 0 || labels.len() != w * h {
|
||||
let Some(n) = w.checked_mul(h).filter(|&n| n > 0 && n as u64 <= crate::MAX_PIXELS && n == labels.len()) else { return (vec![], out) };
|
||||
let mut thin = vec![false; n];
|
||||
if w == 0 || h == 0 || !max_width.is_finite() || max_width <= 0.0 || !ink.iter().all(|v| v.is_finite()) {
|
||||
return (thin, out);
|
||||
}
|
||||
let mask = |i: usize| labels.get(i) == Some(&label);
|
||||
let dt = distance(&mask, w, h);
|
||||
let mut seen = vec![false; w * h];
|
||||
let mut stack = vec![];
|
||||
for start in 0..w * h {
|
||||
if !mask(start) || seen[start] {
|
||||
if !mask(start) || seen.get(start).copied().unwrap_or(false) {
|
||||
continue;
|
||||
}
|
||||
let mut comp = vec![];
|
||||
seen[start] = true;
|
||||
if let Some(v) = seen.get_mut(start) {
|
||||
*v = true;
|
||||
}
|
||||
stack.push(start);
|
||||
while let Some(i) = stack.pop() {
|
||||
comp.push(i);
|
||||
for j in neighbours(i, w, h) {
|
||||
if mask(j) && !seen[j] {
|
||||
seen[j] = true;
|
||||
if mask(j) && !seen.get(j).copied().unwrap_or(false) {
|
||||
if let Some(v) = seen.get_mut(j) {
|
||||
*v = true;
|
||||
}
|
||||
stack.push(j);
|
||||
}
|
||||
}
|
||||
}
|
||||
let widest = 2.0 * f64::from(comp.iter().map(|&i| dt[i]).max().unwrap_or(0)) / 3.0 - 1.0;
|
||||
let Some(g) = Grid::from_component(&comp, w) else { continue };
|
||||
let Some(dt) = distance::squared(&g.on, g.w, g.h) else { continue };
|
||||
let widest = 2.0 * f64::from(dt.into_iter().fold(0.0f32, f32::max)).sqrt() - 1.0;
|
||||
if widest > max_width {
|
||||
continue;
|
||||
}
|
||||
let Some(line) = trace_line(&comp, w, max_width, img, ink, style) else { continue };
|
||||
let Some(line) = trace_line(g, comp.len(), max_width, img, ink, style) else { continue };
|
||||
for &i in &comp {
|
||||
thin[i] = true;
|
||||
if let Some(v) = thin.get_mut(i) {
|
||||
*v = true;
|
||||
}
|
||||
}
|
||||
out.push(line);
|
||||
}
|
||||
@@ -103,35 +113,16 @@ fn neighbours(i: usize, w: usize, h: usize) -> impl Iterator<Item = usize> {
|
||||
})
|
||||
}
|
||||
|
||||
/// Chamfer 3-4 distance (3 per pixel step, saturating) of each pixel in `mask` to the nearest
|
||||
/// pixel outside it; outside the image counts as outside.
|
||||
fn distance(mask: &impl Fn(usize) -> bool, w: usize, h: usize) -> Vec<u16> {
|
||||
let mut d = vec![0u16; w * h];
|
||||
let at = |d: &[u16], x: i64, y: i64| {
|
||||
if x < 0 || y < 0 || x >= w as i64 || y >= h as i64 { 0 } else { d.get(y as usize * w + x as usize).copied().unwrap_or(0) }
|
||||
};
|
||||
for y in 0..h {
|
||||
for x in 0..w {
|
||||
let i = y * w + x;
|
||||
if mask(i) {
|
||||
let (x, y) = (x as i64, y as i64);
|
||||
let steps = [(x - 1, y, 3), (x - 1, y - 1, 4), (x, y - 1, 3), (x + 1, y - 1, 4)];
|
||||
d[i] = steps.iter().map(|&(nx, ny, c)| at(&d, nx, ny).saturating_add(c)).min().unwrap_or(0);
|
||||
}
|
||||
}
|
||||
/// A valid crop is at most the raster limit plus a one-pixel paper border.
|
||||
fn crop_size(w: usize, h: usize) -> Option<usize> {
|
||||
if w < 3 || h < 3 {
|
||||
return None;
|
||||
}
|
||||
for y in (0..h).rev() {
|
||||
for x in (0..w).rev() {
|
||||
let i = y * w + x;
|
||||
if mask(i) {
|
||||
let (xi, yi) = (x as i64, y as i64);
|
||||
let steps = [(xi + 1, yi, 3), (xi + 1, yi + 1, 4), (xi, yi + 1, 3), (xi - 1, yi + 1, 4)];
|
||||
let v = steps.iter().map(|&(nx, ny, c)| at(&d, nx, ny).saturating_add(c)).min().unwrap_or(0);
|
||||
d[i] = d[i].min(v);
|
||||
}
|
||||
}
|
||||
let inner = (w - 2).checked_mul(h - 2)?;
|
||||
if inner as u64 > crate::MAX_PIXELS {
|
||||
return None;
|
||||
}
|
||||
d
|
||||
w.checked_mul(h).filter(|&n| n <= isize::MAX as usize)
|
||||
}
|
||||
|
||||
/// A grid around one line: its pixels with a one-pixel empty border.
|
||||
@@ -145,19 +136,49 @@ struct Grid {
|
||||
}
|
||||
|
||||
impl Grid {
|
||||
/// The 8 neighbours (N, NE, E, SE, S, SW, W, NW) of interior pixel `i`.
|
||||
fn from_component(comp: &[usize], image_width: usize) -> Option<Self> {
|
||||
if image_width == 0 {
|
||||
return None;
|
||||
}
|
||||
let xs = comp.iter().map(|&i| i % image_width);
|
||||
let ys = comp.iter().map(|&i| i / image_width);
|
||||
let (x0, x1, y0, y1) = (xs.clone().min()?, xs.max()?, ys.clone().min()?, ys.max()?);
|
||||
let (w, h) = (x1.checked_sub(x0)?.checked_add(3)?, y1.checked_sub(y0)?.checked_add(3)?);
|
||||
let n = crop_size(w, h)?;
|
||||
let mut g = Self { on: vec![false; n], w, h, x0: x0 as i64 - 1, y0: y0 as i64 - 1 };
|
||||
for &i in comp {
|
||||
*g.on.get_mut((i / image_width - y0 + 1) * w + i % image_width - x0 + 1)? = true;
|
||||
}
|
||||
Some(g)
|
||||
}
|
||||
#[cfg(test)]
|
||||
fn from_mask(mask: &[bool], w: usize, h: usize) -> Option<Self> {
|
||||
if w.checked_mul(h)? != mask.len() {
|
||||
return None;
|
||||
}
|
||||
let comp: Vec<_> = mask.iter().enumerate().filter_map(|(i, &on)| on.then_some(i)).collect();
|
||||
Self::from_component(&comp, w)
|
||||
}
|
||||
fn offsets(&self) -> [isize; 8] {
|
||||
let w = self.w as isize;
|
||||
[-w, -w + 1, 1, w + 1, w, w - 1, -1, -w - 1]
|
||||
}
|
||||
fn is_interior(&self, i: usize) -> bool {
|
||||
self.w > 0 && i < self.on.len() && !i.is_multiple_of(self.w) && i % self.w + 1 < self.w && i / self.w > 0 && i / self.w + 1 < self.h
|
||||
}
|
||||
/// N, NE, E, SE, S, SW, W, NW, with paper outside the interior.
|
||||
fn ring(&self, i: usize) -> [bool; 8] {
|
||||
let w = self.w;
|
||||
let g = |j: usize| self.on.get(j).copied().unwrap_or(false);
|
||||
[g(i - w), g(i - w + 1), g(i + 1), g(i + w + 1), g(i + w), g(i + w - 1), g(i - 1), g(i - w - 1)]
|
||||
if !self.is_interior(i) {
|
||||
return [false; 8];
|
||||
}
|
||||
self.offsets().map(|o| self.on.get(i.wrapping_add_signed(o)).copied().unwrap_or(false))
|
||||
}
|
||||
fn interior(&self) -> impl Iterator<Item = usize> + use<> {
|
||||
let (w, h) = (self.w, self.h);
|
||||
(1..h.saturating_sub(1)).flat_map(move |y| (1..w.saturating_sub(1)).map(move |x| y * w + x))
|
||||
}
|
||||
fn neighbours(&self, i: usize) -> impl Iterator<Item = usize> + '_ {
|
||||
let w = self.w;
|
||||
[i - w, i - w + 1, i + 1, i + w + 1, i + w, i + w - 1, i - 1, i - w - 1].into_iter().filter(|&j| self.on.get(j).copied().unwrap_or(false))
|
||||
self.offsets().into_iter().map(move |o| i.wrapping_add_signed(o)).filter(move |&j| self.is_interior(i) && self.on.get(j) == Some(&true))
|
||||
}
|
||||
fn degree(&self, i: usize) -> usize {
|
||||
self.neighbours(i).count()
|
||||
@@ -168,66 +189,6 @@ impl Grid {
|
||||
}
|
||||
}
|
||||
|
||||
/// Zhang–Suen thinning. Each sub-iteration picks the pixels to remove from the same state, so the
|
||||
/// skeleton stays centred; they are then removed one by one, each only if it still doesn't
|
||||
/// disconnect anything, so a two-pixel-wide line thins to one pixel instead of vanishing.
|
||||
fn thin(g: &mut Grid) {
|
||||
let removable = |p: &[bool; 8]| {
|
||||
let b = p.iter().filter(|v| **v).count();
|
||||
let a = (0..8).filter(|&k| !p[k] && p[(k + 1) % 8]).count();
|
||||
(2..=6).contains(&b) && a == 1
|
||||
};
|
||||
let mut picked = vec![];
|
||||
loop {
|
||||
let mut changed = false;
|
||||
for step in 0..2 {
|
||||
picked.clear();
|
||||
for i in g.interior() {
|
||||
if !g.on[i] {
|
||||
continue;
|
||||
}
|
||||
let p = g.ring(i);
|
||||
let (n, e, s, w) = (p[0], p[2], p[4], p[6]);
|
||||
let side = if step == 0 { !(e && s && (n || w)) } else { !(n && w && (e || s)) };
|
||||
if side && removable(&p) {
|
||||
picked.push(i);
|
||||
}
|
||||
}
|
||||
for &i in &picked {
|
||||
if removable(&g.ring(i)) {
|
||||
g.on[i] = false;
|
||||
changed = true;
|
||||
}
|
||||
}
|
||||
}
|
||||
if !changed {
|
||||
break;
|
||||
}
|
||||
}
|
||||
// Staircase corners: a pixel whose neighbours stay connected without it (Yokoi's
|
||||
// 8-connectivity number 1) and that doesn't end the line.
|
||||
loop {
|
||||
let mut changed = false;
|
||||
for i in g.interior() {
|
||||
if !g.on[i] || g.degree(i) < 2 {
|
||||
continue;
|
||||
}
|
||||
let p = g.ring(i);
|
||||
// x1 = E, x2 = NE, x3 = N, x4 = NW, x5 = W, x6 = SW, x7 = S, x8 = SE.
|
||||
let x = [p[2], p[1], p[0], p[7], p[6], p[5], p[4], p[3]];
|
||||
let off = |k: usize| !x[k % 8];
|
||||
let c8: usize = [0, 2, 4, 6].iter().map(|&k| usize::from(off(k)) - usize::from(off(k) && off(k + 1) && off(k + 2))).sum();
|
||||
if c8 == 1 {
|
||||
g.on[i] = false;
|
||||
changed = true;
|
||||
}
|
||||
}
|
||||
if !changed {
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/// A run of skeleton pixels between two nodes (`None`: none, a ring).
|
||||
#[derive(Clone)]
|
||||
struct Run {
|
||||
@@ -244,17 +205,9 @@ impl Run {
|
||||
}
|
||||
|
||||
/// The centre lines of one line component (pixel indices of the `w`-wide image).
|
||||
fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64; 3], style: &Style) -> Option<Line> {
|
||||
let xs = comp.iter().map(|&i| i % w);
|
||||
let ys = comp.iter().map(|&i| i / w);
|
||||
let (x0, x1, y0, y1) = (xs.clone().min()?, xs.max()?, ys.clone().min()?, ys.max()?);
|
||||
let (gw, gh) = (x1 - x0 + 3, y1 - y0 + 3);
|
||||
let mut g = Grid { on: vec![false; gw * gh], w: gw, h: gh, x0: x0 as i64 - 1, y0: y0 as i64 - 1 };
|
||||
for &i in comp {
|
||||
g.on[(i / w - y0 + 1) * gw + (i % w - x0 + 1)] = true;
|
||||
}
|
||||
thin(&mut g);
|
||||
let skeleton: Vec<usize> = g.interior().filter(|&i| g.on[i]).collect();
|
||||
fn trace_line(mut g: Grid, pixels: usize, max_width: f64, img: &Raster, ink: [f64; 3], style: &Style) -> Option<Line> {
|
||||
thin::thin(&mut g);
|
||||
let skeleton: Vec<usize> = g.interior().filter(|&i| g.on.get(i) == Some(&true)).collect();
|
||||
if skeleton.len() < 2 {
|
||||
return None;
|
||||
}
|
||||
@@ -263,19 +216,19 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
let mut centres: Vec<Point> = vec![];
|
||||
for &i in &skeleton {
|
||||
let d = g.degree(i);
|
||||
if d == 2 || node[i] != usize::MAX {
|
||||
if d == 2 || node.get(i) != Some(&usize::MAX) {
|
||||
continue;
|
||||
}
|
||||
let id = centres.len();
|
||||
let mut members = vec![i];
|
||||
node[i] = id;
|
||||
*node.get_mut(i)? = id;
|
||||
if d > 2 {
|
||||
let mut k = 0;
|
||||
while let Some(&p) = members.get(k) {
|
||||
k += 1;
|
||||
for q in g.neighbours(p).collect::<Vec<_>>() {
|
||||
if node[q] == usize::MAX && g.degree(q) > 2 {
|
||||
node[q] = id;
|
||||
if *node.get(q)? == usize::MAX && g.degree(q) > 2 {
|
||||
*node.get_mut(q)? = id;
|
||||
members.push(q);
|
||||
}
|
||||
}
|
||||
@@ -288,28 +241,28 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
let mut used: HashSet<(usize, usize)> = HashSet::new();
|
||||
let edge = |a: usize, b: usize| (a.min(b), a.max(b));
|
||||
let mut runs: Vec<Run> = vec![];
|
||||
for &p in skeleton.iter().filter(|&&i| node[i] != usize::MAX) {
|
||||
for &p in skeleton.iter().filter(|&&i| node.get(i) != Some(&usize::MAX)) {
|
||||
for q in g.neighbours(p).collect::<Vec<_>>() {
|
||||
if node[q] == node[p] || !used.insert(edge(p, q)) {
|
||||
if *node.get(q)? == *node.get(p)? || !used.insert(edge(p, q)) {
|
||||
continue;
|
||||
}
|
||||
let mut pts = vec![centres[node[p]], g.centre(q)];
|
||||
let mut pts = vec![*centres.get(*node.get(p)?)?, g.centre(q)];
|
||||
let (mut prev, mut cur) = (p, q);
|
||||
while node[cur] == usize::MAX {
|
||||
while *node.get(cur)? == usize::MAX {
|
||||
let Some(next) = g.neighbours(cur).find(|&n| n != prev && !used.contains(&edge(cur, n))) else { break };
|
||||
used.insert(edge(cur, next));
|
||||
pts.push(g.centre(next));
|
||||
(prev, cur) = (cur, next);
|
||||
}
|
||||
let end = (node[cur] != usize::MAX).then_some(node[cur]);
|
||||
let end = (*node.get(cur)? != usize::MAX).then_some(*node.get(cur)?);
|
||||
if let (Some(e), Some(last)) = (end, pts.last_mut()) {
|
||||
*last = centres[e];
|
||||
*last = *centres.get(e)?;
|
||||
}
|
||||
runs.push(Run { pts, ends: [Some(node[p]), end], closed: false });
|
||||
runs.push(Run { pts, ends: [Some(*node.get(p)?), end], closed: false });
|
||||
}
|
||||
}
|
||||
for &s in &skeleton {
|
||||
if node[s] != usize::MAX || g.neighbours(s).any(|n| used.contains(&edge(s, n))) {
|
||||
if *node.get(s)? != usize::MAX || g.neighbours(s).any(|n| used.contains(&edge(s, n))) {
|
||||
continue;
|
||||
}
|
||||
let mut pts = vec![g.centre(s)];
|
||||
@@ -330,7 +283,7 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
// short of each end).
|
||||
let width = |runs: &[Run]| {
|
||||
let length: f64 = runs.iter().map(Run::length).sum();
|
||||
let px = comp.len() as f64;
|
||||
let px = pixels as f64;
|
||||
((-length + (length * length + 4.0 * px).sqrt()) / 2.0).clamp(1.0, max_width.max(1.0))
|
||||
};
|
||||
let spur = width(&runs);
|
||||
@@ -341,7 +294,7 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
// How many runs end at each node: one is an end of the line.
|
||||
let mut meet = vec![0usize; centres.len()];
|
||||
for e in runs.iter().flat_map(|r| r.ends.iter().flatten()) {
|
||||
meet[*e] += 1;
|
||||
*meet.get_mut(*e)? += 1;
|
||||
}
|
||||
let (mut widths, mut soft, mut read) = (vec![], 0usize, 0usize);
|
||||
for r in &mut runs {
|
||||
@@ -356,17 +309,21 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
} else {
|
||||
None
|
||||
};
|
||||
if end.is_some_and(|e| meet[e] > 1) {
|
||||
if end.is_some_and(|e| meet.get(e).is_some_and(|&n| n > 1)) {
|
||||
continue;
|
||||
}
|
||||
let (a, b) = if r.closed { (old[(i + n - 2 % n) % n], old[(i + 2) % n]) } else { (old[i.saturating_sub(2)], old[(i + 2).min(n - 1)]) };
|
||||
let (a, b) = if r.closed {
|
||||
(*old.get((i + n - 2 % n) % n)?, *old.get((i + 2) % n)?)
|
||||
} else {
|
||||
(*old.get(i.saturating_sub(2))?, *old.get((i + 2).min(n - 1))?)
|
||||
};
|
||||
let d = b - a;
|
||||
if d.hypot() < 1e-9 {
|
||||
continue;
|
||||
}
|
||||
let upright = d.x.abs() >= d.y.abs();
|
||||
let Some(across) = ink.across(old[i], upright) else { continue };
|
||||
r.pts[i] = old[i] + if upright { kurbo::Vec2::new(0.0, across.middle) } else { kurbo::Vec2::new(across.middle, 0.0) };
|
||||
let Some(across) = ink.across(*old.get(i)?, upright) else { continue };
|
||||
*r.pts.get_mut(i)? = *old.get(i)? + if upright { kurbo::Vec2::new(0.0, across.middle) } else { kurbo::Vec2::new(across.middle, 0.0) };
|
||||
widths.push(across.width * if upright { d.x.abs() } else { d.y.abs() } / d.hypot());
|
||||
read += 1;
|
||||
soft += usize::from(across.soft);
|
||||
@@ -382,19 +339,19 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
let (mut short_by, mut told) = (0.0, 0);
|
||||
for r in runs.iter_mut().filter(|r| !r.closed && r.pts.len() >= 2) {
|
||||
for k in 0..2 {
|
||||
if r.ends[k].is_some_and(|e| meet[e] > 1) {
|
||||
if r.ends.get(k).copied().flatten().is_some_and(|e| meet.get(e).is_some_and(|&n| n > 1)) {
|
||||
continue;
|
||||
}
|
||||
let n = r.pts.len();
|
||||
let at = |i: usize| if k == 0 { i.min(n - 1) } else { n - 1 - i.min(n - 1) };
|
||||
// The last points of a wide line's skeleton turn towards a corner of its end: the end
|
||||
// lies where the line was going a width back, and they are put on that line.
|
||||
let back = (width.round() as usize + 2).max(3);
|
||||
let (tip, inner) = (r.pts[at(0)], r.pts[at(back)]);
|
||||
let back = (width.round().clamp(0.0, n as f64) as usize + 2).max(3);
|
||||
let (tip, inner) = (*r.pts.get(at(0))?, *r.pts.get(at(back))?);
|
||||
// (Only an end that has left that line: over so few pixels the line's direction is known
|
||||
// to a few degrees at best, less without anti-aliasing, and an end already on it is
|
||||
// better left where it is.)
|
||||
let before = inner - r.pts[at(2 * back)];
|
||||
let before = inner - *r.pts.get(at(2 * back))?;
|
||||
let wide = width >= 3.0 && n > 2 * back && before.cross(tip - inner).abs() > (0.25 * width).max(1.0) * before.hypot();
|
||||
let d = if wide { before } else { tip - inner };
|
||||
if d.hypot() < 1e-9 {
|
||||
@@ -404,7 +361,8 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
let on_line = |p: Point| inner + along * (p - inner).dot(along).max(0.0);
|
||||
if wide {
|
||||
for i in 0..back {
|
||||
r.pts[at(i)] = on_line(r.pts[at(i)]);
|
||||
let point = r.pts.get_mut(at(i))?;
|
||||
*point = on_line(*point);
|
||||
}
|
||||
}
|
||||
let tip = on_line(tip);
|
||||
@@ -417,7 +375,7 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
}
|
||||
let more = reach - width / 2.0;
|
||||
if more.abs() > 0.05 && more > 1.0 - (tip - inner).hypot() {
|
||||
r.pts[at(0)] = tip + along * more;
|
||||
*r.pts.get_mut(at(0))? = tip + along * more;
|
||||
}
|
||||
}
|
||||
}
|
||||
@@ -437,7 +395,7 @@ fn trace_line(comp: &[usize], w: usize, max_width: f64, img: &Raster, ink: [f64;
|
||||
if subs.is_empty() {
|
||||
return None;
|
||||
}
|
||||
Some(Line { path: PathData::new(subs), width, pixels: comp.len(), square: told > 0 && short_by / f64::from(told) < 0.07 })
|
||||
Some(Line { path: PathData::new(subs), width, pixels, square: told > 0 && short_by / f64::from(told) < 0.07 })
|
||||
}
|
||||
|
||||
/// What the image says of a line's ink.
|
||||
@@ -464,9 +422,9 @@ impl Ink<'_> {
|
||||
if x < 0 || y < 0 || x >= w || y >= h {
|
||||
return None;
|
||||
}
|
||||
let at = (y * w + x) as usize * 4;
|
||||
let p = self.img.rgba.get(at..at + 3)?;
|
||||
Some((0..3).map(|k| (f64::from(p[k]) - self.color[k]).powi(2)).sum::<f64>().sqrt())
|
||||
let at = (y as usize).checked_mul(w as usize)?.checked_add(x as usize)?.checked_mul(4)?;
|
||||
let p = self.img.rgba.get(at..at.checked_add(3)?)?;
|
||||
Some(p.iter().zip(self.color).map(|(&v, c)| (f64::from(v) - c).powi(2)).sum::<f64>().sqrt())
|
||||
}
|
||||
/// How far the paper's colour is from the line's around pixel `(x, y)`: the farthest colour
|
||||
/// within reach, along its row and its column.
|
||||
@@ -475,6 +433,9 @@ impl Ink<'_> {
|
||||
}
|
||||
/// The ink across the line through the pixel `p` lies in: along its column (`upright`) or row.
|
||||
fn across(&self, p: Point, upright: bool) -> Option<Across> {
|
||||
if !p.x.is_finite() || !p.y.is_finite() {
|
||||
return None;
|
||||
}
|
||||
let (x, y) = (p.x.floor() as i64, p.y.floor() as i64);
|
||||
let paper = self.paper(x, y);
|
||||
// Paper too like the ink says nothing.
|
||||
@@ -506,6 +467,9 @@ impl Ink<'_> {
|
||||
}
|
||||
/// How much ink (pixels) lies beyond `tip`, going `dir` (a unit vector).
|
||||
fn beyond(&self, tip: Point, dir: kurbo::Vec2) -> f64 {
|
||||
if !tip.x.is_finite() || !tip.y.is_finite() || !dir.x.is_finite() || !dir.y.is_finite() {
|
||||
return 0.0;
|
||||
}
|
||||
let paper = self.paper(tip.x.floor() as i64, tip.y.floor() as i64);
|
||||
if paper < 24.0 {
|
||||
return 0.0;
|
||||
@@ -528,19 +492,27 @@ impl Ink<'_> {
|
||||
fn prune(runs: Vec<Run>, nodes: usize, width: f64) -> Vec<Run> {
|
||||
let mut degree = vec![0usize; nodes];
|
||||
for e in runs.iter().flat_map(|r| r.ends.iter().flatten()) {
|
||||
degree[*e] += 1;
|
||||
if let Some(d) = degree.get_mut(*e) {
|
||||
*d += 1;
|
||||
}
|
||||
}
|
||||
let mut order: Vec<usize> = (0..runs.len()).collect();
|
||||
order.sort_by(|a, b| runs[*a].length().total_cmp(&runs[*b].length()));
|
||||
order.sort_by(|a, b| runs.get(*a).map_or(0.0, Run::length).total_cmp(&runs.get(*b).map_or(0.0, Run::length)));
|
||||
let mut keep = vec![true; runs.len()];
|
||||
for i in order {
|
||||
let r = &runs[i];
|
||||
let Some(r) = runs.get(i) else { continue };
|
||||
let [Some(a), Some(b)] = r.ends else { continue };
|
||||
let (tip, fork) = if degree[a] == 1 { (a, b) } else { (b, a) };
|
||||
if degree[tip] == 1 && degree[fork] >= 3 && r.length() < width.max(2.0) {
|
||||
keep[i] = false;
|
||||
degree[tip] -= 1;
|
||||
degree[fork] -= 1;
|
||||
let (tip, fork) = if degree.get(a).copied().unwrap_or(0) == 1 { (a, b) } else { (b, a) };
|
||||
if degree.get(tip).copied().unwrap_or(0) == 1 && degree.get(fork).copied().unwrap_or(0) >= 3 && r.length() < width.max(2.0) {
|
||||
if let Some(k) = keep.get_mut(i) {
|
||||
*k = false;
|
||||
}
|
||||
if let Some(d) = degree.get_mut(tip) {
|
||||
*d = d.saturating_sub(1);
|
||||
}
|
||||
if let Some(d) = degree.get_mut(fork) {
|
||||
*d = d.saturating_sub(1);
|
||||
}
|
||||
}
|
||||
}
|
||||
runs.into_iter().zip(keep).filter_map(|(r, k)| k.then_some(r)).collect()
|
||||
@@ -553,18 +525,18 @@ fn join(runs: Vec<Run>, nodes: usize) -> Vec<Run> {
|
||||
let at: Vec<(usize, usize)> = runs
|
||||
.iter()
|
||||
.enumerate()
|
||||
.flat_map(|(i, r)| r.iter().flat_map(move |r| (0..2).filter(move |&k| r.ends[k] == Some(n)).map(move |k| (i, k))))
|
||||
.flat_map(|(i, r)| r.iter().flat_map(move |r| (0..2).filter(move |&k| r.ends.get(k) == Some(&Some(n))).map(move |k| (i, k))))
|
||||
.collect();
|
||||
let [(i, ki), (j, kj)] = at[..] else { continue };
|
||||
let &[(i, ki), (j, kj)] = at.as_slice() else { continue };
|
||||
if i == j {
|
||||
if let Some(r) = runs[i].as_mut() {
|
||||
if let Some(r) = runs.get_mut(i).and_then(Option::as_mut) {
|
||||
r.pts.pop();
|
||||
r.ends = [None, None];
|
||||
r.closed = r.pts.len() >= 3;
|
||||
}
|
||||
continue;
|
||||
}
|
||||
let (Some(mut a), Some(mut b)) = (runs[i].take(), runs[j].take()) else { continue };
|
||||
let (Some(mut a), Some(mut b)) = (runs.get_mut(i).and_then(Option::take), runs.get_mut(j).and_then(Option::take)) else { continue };
|
||||
// a ends at n, b starts at n.
|
||||
if ki == 0 {
|
||||
a.pts.reverse();
|
||||
@@ -576,7 +548,9 @@ fn join(runs: Vec<Run>, nodes: usize) -> Vec<Run> {
|
||||
}
|
||||
a.pts.extend(b.pts.into_iter().skip(1));
|
||||
a.ends[1] = b.ends[1];
|
||||
runs[i] = Some(a);
|
||||
if let Some(r) = runs.get_mut(i) {
|
||||
*r = Some(a);
|
||||
}
|
||||
}
|
||||
runs.into_iter().flatten().collect()
|
||||
}
|
||||
@@ -0,0 +1,196 @@
|
||||
//! Guo–Hall A1 thinning (Guo & Hall 1989, doi:10.1145/62065.62074).
|
||||
//! Neighbourhoods are judged together in each sub-iteration; only changed ones are revisited.
|
||||
|
||||
use super::Grid;
|
||||
|
||||
/// Guo–Hall A1's test for pixel with neighbours `p` (P2..P9) in sub-iteration `step`.
|
||||
fn guo_hall_removable(p: &[bool; 8], step: usize) -> bool {
|
||||
let [p2, p3, p4, p5, p6, p7, p8, p9] = *p;
|
||||
let c = usize::from(!p2 && (p3 || p4)) + usize::from(!p4 && (p5 || p6)) + usize::from(!p6 && (p7 || p8)) + usize::from(!p8 && (p9 || p2));
|
||||
let n1 = usize::from(p9 || p2) + usize::from(p3 || p4) + usize::from(p5 || p6) + usize::from(p7 || p8);
|
||||
let n2 = usize::from(p2 || p3) + usize::from(p4 || p5) + usize::from(p6 || p7) + usize::from(p8 || p9);
|
||||
let m = if step == 0 { (p6 || p7 || !p9) && p8 } else { (p2 || p3 || !p5) && p4 };
|
||||
c == 1 && (2..=3).contains(&n1.min(n2)) && !m
|
||||
}
|
||||
|
||||
/// The Guo–Hall rounds alone (no staircase removal). Bounded: at most `w + h` rounds.
|
||||
///
|
||||
/// Only pixels whose neighbourhood changed since they were last judged in a sub-iteration are
|
||||
/// judged again (a pixel's removability depends on its 8 neighbours only): in sub-iteration 0,
|
||||
/// neighbours of every pixel removed in the last round; in sub-iteration 1, neighbours of those
|
||||
/// removed in the last round's sub-iteration 1 and this round's sub-iteration 0. The result is
|
||||
/// the plain algorithm's (tested against it below).
|
||||
fn guo_hall_rounds(g: &mut Grid) {
|
||||
let n = g.on.len();
|
||||
let ring = g.offsets();
|
||||
let mut mark = vec![0u32; n];
|
||||
let mut stamp = 0u32;
|
||||
let mut around = |g: &Grid, removed: &[&[usize]]| -> Vec<usize> {
|
||||
// Two stamps per round; the crop cap keeps this well below u32::MAX.
|
||||
stamp += 1;
|
||||
let mut out = vec![];
|
||||
for &t in removed.iter().flat_map(|r| r.iter()) {
|
||||
for o in ring {
|
||||
let i = t.wrapping_add_signed(o);
|
||||
if i < n && g.on.get(i) == Some(&true) && g.is_interior(i) && mark.get(i) != Some(&stamp) {
|
||||
if let Some(m) = mark.get_mut(i) {
|
||||
*m = stamp;
|
||||
}
|
||||
out.push(i);
|
||||
}
|
||||
}
|
||||
}
|
||||
out.sort_unstable();
|
||||
out
|
||||
};
|
||||
let (mut prev0, mut prev1): (Vec<usize>, Vec<usize>) = (vec![], vec![]);
|
||||
for round in 0..g.w + g.h {
|
||||
let cand0: Vec<usize> =
|
||||
if round == 0 { (0..n).filter(|&i| g.on.get(i) == Some(&true) && g.is_interior(i)).collect() } else { around(g, &[&prev0, &prev1]) };
|
||||
let rem0: Vec<usize> = cand0.into_iter().filter(|&i| guo_hall_removable(&g.ring(i), 0)).collect();
|
||||
for &i in &rem0 {
|
||||
if let Some(on) = g.on.get_mut(i) {
|
||||
*on = false;
|
||||
}
|
||||
}
|
||||
let cand1: Vec<usize> =
|
||||
if round == 0 { (0..n).filter(|&i| g.on.get(i) == Some(&true) && g.is_interior(i)).collect() } else { around(g, &[&prev1, &rem0]) };
|
||||
let rem1: Vec<usize> = cand1.into_iter().filter(|&i| guo_hall_removable(&g.ring(i), 1)).collect();
|
||||
for &i in &rem1 {
|
||||
if let Some(on) = g.on.get_mut(i) {
|
||||
*on = false;
|
||||
}
|
||||
}
|
||||
if rem0.is_empty() && rem1.is_empty() {
|
||||
break;
|
||||
}
|
||||
(prev0, prev1) = (rem0, rem1);
|
||||
}
|
||||
}
|
||||
|
||||
/// Thins `g` in place (Guo–Hall, then staircase removal); see the module docs.
|
||||
pub(super) fn thin(g: &mut Grid) {
|
||||
guo_hall_rounds(g);
|
||||
staircase(g);
|
||||
}
|
||||
|
||||
/// Removes staircase pixels (see the module docs). Sequential, in scan order, until stable
|
||||
/// (at most 4 passes).
|
||||
fn staircase(g: &mut Grid) {
|
||||
for _ in 0..4 {
|
||||
let mut changed = false;
|
||||
for i in 0..g.on.len() {
|
||||
// Three or more neighbours: a corner of a staircase always has them, while a pixel
|
||||
// with two touching neighbours may be the end of a line two pixels thick (removing
|
||||
// it would eat that line from its end, pixel by pixel).
|
||||
if g.on.get(i) != Some(&true) || g.degree(i) < 3 {
|
||||
continue;
|
||||
}
|
||||
let p = g.ring(i);
|
||||
// Yokoi's x1..x8 = E, NE, N, NW, W, SW, S, SE.
|
||||
let x = [p[2], p[1], p[0], p[7], p[6], p[5], p[4], p[3]];
|
||||
let off = |k: usize| !x.get(k % 8).copied().unwrap_or(false);
|
||||
let c8: usize = [0usize, 2, 4, 6].iter().map(|&k| usize::from(off(k)) - usize::from(off(k) && off(k + 1) && off(k + 2))).sum();
|
||||
if c8 == 1 {
|
||||
if let Some(on) = g.on.get_mut(i) {
|
||||
*on = false;
|
||||
}
|
||||
changed = true;
|
||||
}
|
||||
}
|
||||
if !changed {
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
#[cfg(test)]
|
||||
mod tests {
|
||||
use super::*;
|
||||
|
||||
fn skel(w: usize, h: usize, f: impl Fn(usize, usize) -> bool) -> Grid {
|
||||
let m: Vec<bool> = (0..w * h).map(|i| f(i % w, i / w)).collect();
|
||||
let mut g = Grid::from_mask(&m, w, h).unwrap();
|
||||
thin(&mut g);
|
||||
g
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn bars_thin_to_lines() {
|
||||
for thick in 1..7 {
|
||||
let g = skel(40, 12, |x, y| (3..37).contains(&x) && (3..3 + thick).contains(&y));
|
||||
let on: Vec<usize> = (0..g.on.len()).filter(|&i| g.on[i]).collect();
|
||||
assert!(on.len() >= 20, "thick {thick}: {} pixels", on.len());
|
||||
// Corner spurs may remain at the square ends (the graph prunes them).
|
||||
assert!(on.iter().filter(|&&i| g.degree(i) == 1).count() >= 2, "thick {thick}");
|
||||
}
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn diagonal_lines_two_to_five_pixels_thick_keep_their_length() {
|
||||
for thick in 2..6usize {
|
||||
for mirror in [false, true] {
|
||||
let g = skel(60, 60, |x, y| {
|
||||
let x = if mirror { 59 - x } else { x };
|
||||
(5..55).contains(&y) && x + 1 >= y && x < y + thick
|
||||
});
|
||||
let on: Vec<usize> = (0..g.on.len()).filter(|&i| g.on[i]).collect();
|
||||
let (lo, hi) = (on.iter().map(|&i| i / g.w).min().unwrap(), on.iter().map(|&i| i / g.w).max().unwrap());
|
||||
assert!(hi - lo >= 46, "thick {thick}, mirror {mirror}: shortened to {} rows", hi - lo + 1);
|
||||
assert!(on.iter().all(|&i| g.degree(i) <= 2), "thick {thick}: not one pixel wide");
|
||||
assert_eq!(on.iter().filter(|&&i| g.degree(i) == 1).count(), 2);
|
||||
let mut seen = std::collections::HashSet::new();
|
||||
let mut stack = vec![on[0]];
|
||||
while let Some(i) = stack.pop() {
|
||||
if seen.insert(i) {
|
||||
stack.extend(g.neighbours(i));
|
||||
}
|
||||
}
|
||||
assert_eq!(seen.len(), on.len(), "disconnected diagonal");
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/// Plain Guo–Hall: every pixel judged every sub-iteration.
|
||||
fn reference(g: &mut Grid) {
|
||||
loop {
|
||||
let mut any = false;
|
||||
for step in 0..2 {
|
||||
let picked: Vec<usize> =
|
||||
(0..g.on.len()).filter(|&i| g.on.get(i) == Some(&true) && g.is_interior(i) && guo_hall_removable(&g.ring(i), step)).collect();
|
||||
any |= !picked.is_empty();
|
||||
for i in picked {
|
||||
if let Some(on) = g.on.get_mut(i) {
|
||||
*on = false;
|
||||
}
|
||||
}
|
||||
}
|
||||
if !any {
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn the_active_set_gives_the_plain_algorithms_result() {
|
||||
let mut seed = 11u64;
|
||||
let mut rnd = || {
|
||||
seed = seed.wrapping_mul(6364136223846793005).wrapping_add(1442695040888963407);
|
||||
(seed >> 33) as usize
|
||||
};
|
||||
for case in 0..200 {
|
||||
let (w, h) = (3 + rnd() % 40, 3 + rnd() % 40);
|
||||
let blobs: Vec<(usize, usize, usize)> = (0..1 + rnd() % 6).map(|_| (rnd() % w, rnd() % h, 1 + rnd() % 9)).collect();
|
||||
let m: Vec<bool> = (0..w * h)
|
||||
.map(|i| {
|
||||
let (x, y) = ((i % w) as i64, (i / w) as i64);
|
||||
blobs.iter().any(|&(bx, by, r)| (x - bx as i64).pow(2) + (y - by as i64).pow(2) <= (r * r) as i64)
|
||||
})
|
||||
.collect();
|
||||
let (mut a, mut b) = (Grid::from_mask(&m, w, h).unwrap(), Grid::from_mask(&m, w, h).unwrap());
|
||||
guo_hall_rounds(&mut a);
|
||||
reference(&mut b);
|
||||
assert!(a.on == b.on, "case {case}");
|
||||
}
|
||||
}
|
||||
}
|
||||
Reference in New Issue
Block a user