shifathosn
/
effectcraft
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

M13.6: scopes, content-aware fill (PatchMatch/flow propagation/membrane), cut detection and auto-trace algorithms 

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This commit is contained in:
Brandon Thomas 2026-10-03 02:21:14 -04:00
parent afd69fb97c
commit b8a56eb4a9
6 changed files with 2137 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
crates/path/src/lib.rs
View File

@ -11,6 +11,7 @@ mod fit;
pub mod interp;
pub mod offset;
pub mod ops;
pub mod trace;
pub mod varstroke;
use effectcraft_geom::Mat3;

417
crates/path/src/trace.rs Normal file
View File

@ -0,0 +1,417 @@
//! Auto-trace: outlines of a scalar field (a layer's alpha, a colour channel or luminance) as
//! Bezier paths (Layer ▸ Auto-trace).
//!
//! 1. **Marching squares** over the samples at pixel centres, with linear interpolation of the
//! threshold crossing on each cell edge (saddles resolved by the cell's mean), on a field padded
//! with "outside" so every contour closes.
//! 2. **Douglas–Peucker** simplification of each closed loop to the Tolerance (pixels).
//! 3. **Bezier fitting**: each simplified span is fitted with one cubic (least squares with fixed
//! end tangents, Newton reparameterisation; [`crate::fit`]). The tangent at each vertex is the
//! chord direction at Corner Roundness 0 (straight segments, sharp corners) and the smooth
//! (Catmull–Rom) direction at 100 %, blended between.
//!
//! Loops smaller than Minimum Area are dropped. Each loop gets its nesting depth (0 = outermost),
//! so callers can add outer outlines and subtract holes.
use std::collections::HashMap;
use effectcraft_keyframe::ShapePath;
use kurbo::{Point, Vec2};
use crate::fit::{fit_single, unit};
/// Auto-trace options.
#[derive(Clone, Copy, Debug, PartialEq)]
pub struct TraceOpts {
/// Inside where the field is at least this (0…1).
pub threshold: f32,
/// Douglas–Peucker tolerance in pixels.
pub tolerance: f64,
/// Loops with a smaller area (square pixels) are dropped.
pub min_area: f64,
/// 0 = polygons, 1 = fully smooth.
pub corner_roundness: f64,
}
impl Default for TraceOpts {
fn default() -> Self {
TraceOpts { threshold: 0.5, tolerance: 1.0, min_area: 10.0, corner_roundness: 0.5 }
}
}
/// A traced outline.
#[derive(Clone, Debug, PartialEq)]
pub struct Traced {
pub path: ShapePath,
/// Number of other outlines that contain this one (even = filled, odd = hole).
pub depth: usize,
/// Signed area of the simplified polygon (pixels²).
pub area: f64,
}
/// Trace `field` (`w × h`, row-major, sample (x, y) at pixel centre (x + 0.5, y + 0.5)).
pub fn trace(field: &[f32], w: usize, h: usize, o: &TraceOpts) -> Vec<Traced> {
if w == 0 || h == 0 || field.len() < w * h {
return vec![];
}
let loops = contours(field, w, h, o.threshold);
let mut polys: Vec<(Vec<[f64; 2]>, Vec<[f64; 2]>, f64)> = loops
.into_iter()
.filter_map(|raw| {
let a = signed_area(&raw);
if a.abs() < o.min_area.max(0.0) || raw.len() < 3 {
return None;
}
let simp = simplify_closed(&raw, o.tolerance.max(0.01));
if simp.len() < 3 {
return None;
}
Some((raw, simp, a))
})
.collect();
// Depth: how many other loops contain a point of this loop.
let depths: Vec<usize> = (0..polys.len())
.map(|i| {
let p = polys[i].1[0];
(0..polys.len()).filter(|&j| j != i && polys[j].2.abs() > polys[i].2.abs() && point_in_poly(p, &polys[j].1)).count()
})
.collect();
let mut out: Vec<Traced> = polys
.drain(..)
.zip(depths)
.map(|((raw, simp, a), depth)| Traced { path: fit_loop(&raw, &simp, o.corner_roundness.clamp(0.0, 1.0)), depth, area: a })
.collect();
// Outer outlines first, then holes, islands…; larger first within a depth.
out.sort_by(|a, b| a.depth.cmp(&b.depth).then(b.area.abs().total_cmp(&a.area.abs())));
out
}
/// Marching squares: closed polylines (pixel coordinates).
fn contours(field: &[f32], w: usize, h: usize, thr: f32) -> Vec<Vec<[f64; 2]>> {
// Padded grid: sample (i, j) of the padded grid is pixel (i - 1, j - 1), outside = -inf.
let (gw, gh) = (w + 2, h + 2);
let val = |i: usize, j: usize| -> f32 { if i == 0 || j == 0 || i > w || j > h { f32::NEG_INFINITY } else { field[(j - 1) * w + (i - 1)] } };
let inside = |i: usize, j: usize| val(i, j) >= thr;
// Edge ids: horizontal edge (i, j)-(i+1, j) = 2 * (j * gw + i); vertical (i, j)-(i, j+1) = +1.
let pos = |i: usize, j: usize| [i as f64 - 0.5, j as f64 - 0.5];
let cross = |a: (usize, usize), b: (usize, usize)| -> [f64; 2] {
let (va, vb) = (val(a.0, a.1), val(b.0, b.1));
let t = if va.is_finite() && vb.is_finite() && (vb - va).abs() > 1e-12 {
((thr - va) / (vb - va)).clamp(0.0, 1.0) as f64
} else if va.is_finite() {
// b is padding: the crossing sits half a pixel out.
if va >= thr { 0.5 } else { 0.0 }
} else if vb >= thr {
0.5
} else {
1.0
};
let (pa, pb) = (pos(a.0, a.1), pos(b.0, b.1));
[pa[0] + (pb[0] - pa[0]) * t, pa[1] + (pb[1] - pa[1]) * t]
};
let hid = |i: usize, j: usize| 2 * (j * gw + i);
let vid = |i: usize, j: usize| 2 * (j * gw + i) + 1;
// Directed segments (from edge, to edge) keeping the inside on the left.
let mut next: HashMap<usize, usize> = HashMap::new();
let mut point: HashMap<usize, [f64; 2]> = HashMap::new();
for j in 0..gh - 1 {
for i in 0..gw - 1 {
// Corners: tl (i, j), tr (i+1, j), br (i+1, j+1), bl (i, j+1).
let (tl, tr, br, bl) = (inside(i, j), inside(i + 1, j), inside(i + 1, j + 1), inside(i, j + 1));
let code = (tl as u8) << 3 | (tr as u8) << 2 | (br as u8) << 1 | bl as u8;
if code == 0 || code == 15 {
continue;
}
let top = hid(i, j);
let bottom = hid(i, j + 1);
let left = vid(i, j);
let right = vid(i + 1, j);
let mut seg = |a: usize, b: usize| {
next.insert(a, b);
};
// Walk with the inside on the right in image coordinates (y down), which is
// counter-clockwise on screen for outer loops… orientation is fixed afterwards.
match code {
1 => seg(left, bottom),
2 => seg(bottom, right),
3 => seg(left, right),
4 => seg(right, top),
5 => {
// tr + bl inside (saddle).
let centre = (val(i, j) + val(i + 1, j) + val(i + 1, j + 1) + val(i, j + 1)) / 4.0;
if centre >= thr {
seg(left, top);
seg(right, bottom);
} else {
seg(left, bottom);
seg(right, top);
}
}
6 => seg(bottom, top),
7 => seg(left, top),
8 => seg(top, left),
9 => seg(top, bottom),
10 => {
// tl + br inside (saddle).
let centre = (val(i, j) + val(i + 1, j) + val(i + 1, j + 1) + val(i, j + 1)) / 4.0;
if centre >= thr {
seg(top, right);
seg(bottom, left);
} else {
seg(top, left);
seg(bottom, right);
}
}
11 => seg(top, right),
12 => seg(right, left),
13 => seg(right, bottom),
14 => seg(bottom, left),
_ => {}
}
for (e, a, b) in [(top, (i, j), (i + 1, j)), (bottom, (i, j + 1), (i + 1, j + 1)), (left, (i, j), (i, j + 1)), (right, (i + 1, j), (i + 1, j + 1))]
{
if inside(a.0, a.1) != inside(b.0, b.1) {
point.entry(e).or_insert_with(|| cross(a, b));
}
}
}
}
let mut out = vec![];
let mut keys: Vec<usize> = next.keys().copied().collect();
keys.sort_unstable();
for start in keys {
if !next.contains_key(&start) {
continue;
}
let mut poly = vec![];
let mut e = start;
while let Some(n) = next.remove(&e) {
if let Some(p) = point.get(&e) {
poly.push(*p);
}
e = n;
if e == start {
break;
}
}
if poly.len() >= 3 {
out.push(poly);
}
}
out
}
/// Shoelace area (positive = clockwise on screen with y down).
pub fn signed_area(p: &[[f64; 2]]) -> f64 {
let n = p.len();
let mut a = 0.0;
for i in 0..n {
let (q, r) = (p[i], p[(i + 1) % n]);
a += q[0] * r[1] - r[0] * q[1];
}
a * 0.5
}
fn point_in_poly(p: [f64; 2], poly: &[[f64; 2]]) -> bool {
let mut inside = false;
let n = poly.len();
let mut j = n - 1;
for i in 0..n {
let (a, b) = (poly[i], poly[j]);
if (a[1] > p[1]) != (b[1] > p[1]) && p[0] < (b[0] - a[0]) * (p[1] - a[1]) / (b[1] - a[1]) + a[0] {
inside = !inside;
}
j = i;
}
inside
}
fn seg_dist(p: [f64; 2], a: [f64; 2], b: [f64; 2]) -> f64 {
let (dx, dy) = (b[0] - a[0], b[1] - a[1]);
let l2 = dx * dx + dy * dy;
let t = if l2 > 0.0 { (((p[0] - a[0]) * dx + (p[1] - a[1]) * dy) / l2).clamp(0.0, 1.0) } else { 0.0 };
((p[0] - a[0] - t * dx).powi(2) + (p[1] - a[1] - t * dy).powi(2)).sqrt()
}
/// Douglas–Peucker on an open polyline: kept indices.
fn dp(p: &[[f64; 2]], lo: usize, hi: usize, tol: f64, keep: &mut Vec<usize>) {
if hi <= lo + 1 {
return;
}
let (mut best, mut bi) = (0.0, lo);
for (k, q) in p.iter().enumerate().take(hi).skip(lo + 1) {
let d = seg_dist(*q, p[lo], p[hi]);
if d > best {
best = d;
bi = k;
}
}
if best > tol {
dp(p, lo, bi, tol, keep);
keep.push(bi);
dp(p, bi, hi, tol, keep);
}
}
/// Douglas–Peucker on a closed loop: indices into `p` of the kept vertices, in order.
fn simplify_idx(p: &[[f64; 2]], tol: f64) -> Vec<usize> {
let n = p.len();
// Split at the two mutually farthest-ish points (first point and the point farthest from it).
let far = (0..n).max_by(|&a, &b| {
let d = |i: usize| (p[i][0] - p[0][0]).powi(2) + (p[i][1] - p[0][1]).powi(2);
d(a).total_cmp(&d(b))
});
let far = far.unwrap_or(0);
let mut ext: Vec<[f64; 2]> = p.to_vec();
ext.push(p[0]);
let mut keep = vec![0];
dp(&ext, 0, far, tol, &mut keep);
keep.push(far);
dp(&ext, far, n, tol, &mut keep);
keep.sort_unstable();
keep.dedup();
keep.retain(|&i| i < n);
keep
}
fn simplify_closed(p: &[[f64; 2]], tol: f64) -> Vec<[f64; 2]> {
simplify_idx(p, tol).into_iter().map(|i| p[i]).collect()
}
/// Fit cubics through the simplified vertices of a loop (`raw` holds the full contour).
fn fit_loop(raw: &[[f64; 2]], simp: &[[f64; 2]], round: f64) -> ShapePath {
let n = simp.len();
if round <= 1e-6 {
return ShapePath::polygon(simp, true);
}
// Raw index of each simplified vertex (they are raw points).
let idx: Vec<usize> = simp.iter().map(|s| raw.iter().position(|r| r == s).unwrap_or(0)).collect();
let pt = |a: [f64; 2]| Point::new(a[0], a[1]);
let v = |a: [f64; 2], b: [f64; 2]| Vec2::new(b[0] - a[0], b[1] - a[1]);
// Per vertex: (direction leaving, direction arriving), blended between chord and smooth.
let dirs: Vec<(Vec2, Vec2)> = (0..n)
.map(|i| {
let (prev, cur, nxt) = (simp[(i + n - 1) % n], simp[i], simp[(i + 1) % n]);
let smooth = unit(v(prev, nxt)).unwrap_or(Vec2::new(1.0, 0.0));
let out_chord = unit(v(cur, nxt)).unwrap_or(smooth);
let in_chord = unit(v(prev, cur)).unwrap_or(smooth);
let leave = unit(out_chord * (1.0 - round) + smooth * round).unwrap_or(out_chord);
let arrive = unit(in_chord * (1.0 - round) + smooth * round).unwrap_or(in_chord);
(leave, arrive)
})
.collect();
let mut path = ShapePath::polygon(simp, true);
for i in 0..n {
let j = (i + 1) % n;
let (a, b) = (idx[i], idx[j]);
let mut pts: Vec<Point> =
if b > a { raw[a..=b].iter().map(|q| pt(*q)).collect() } else { raw[a..].iter().chain(raw[..=b].iter()).map(|q| pt(*q)).collect() };
if pts.len() < 2 {
pts = vec![pt(simp[i]), pt(simp[j])];
}
let chord = pts[0].distance(pts[pts.len() - 1]);
let (t0, t1) = (dirs[i].0, dirs[j].1);
let (c, _, _) = if pts.len() >= 3 { fit_single(&pts, t0, t1) } else { (kurbo::CubicBez::new(pts[0], pts[0], pts[1], pts[1]), 0.0, 0) };
// Keep handles pointing along the tangents and within the chord.
let max = chord * 0.5 * round.max(0.34);
let l0 = (c.p1 - c.p0).dot(t0).clamp(0.0, max);
let l1 = (c.p3 - c.p2).dot(t1).clamp(0.0, max);
let (l0, l1) = if pts.len() < 3 { (chord / 3.0 * round, chord / 3.0 * round) } else { (l0, l1) };
path.out_tangents[i] = [t0.x * l0, t0.y * l0];
path.in_tangents[j] = [-t1.x * l1, -t1.y * l1];
}
path
}
#[cfg(test)]
mod tests {
use super::*;
use crate::{FillRule, fill_coverage, to_kurbo};
use effectcraft_geom::Mat3;
fn iou(a: &[bool], b: &[f32]) -> f64 {
let (mut i, mut u) = (0.0, 0.0);
for (x, y) in a.iter().zip(b) {
let y = *y >= 0.5;
if *x && y {
i += 1.0;
}
if *x || y {
u += 1.0;
}
}
i / u
}
/// Rasterise traced outlines with even-odd fill (as the Add / Subtract mask stack does).
fn raster(t: &[Traced], w: usize, h: usize) -> Vec<f32> {
let paths: Vec<_> = t.iter().map(|t| to_kurbo(&t.path)).collect();
fill_coverage(&paths, &Mat3::IDENTITY, w as u32, h as u32, FillRule::EvenOdd).data
}
#[test]
fn traces_a_disc_and_a_square_with_a_hole() {
let (w, h) = (120, 80);
let mut field = vec![0.0f32; w * h];
let mut truth = vec![false; w * h];
for y in 0..h {
for x in 0..w {
let (fx, fy) = (x as f64 + 0.5, y as f64 + 0.5);
let disc = (fx - 30.0).powi(2) + (fy - 40.0).powi(2) <= 20.0f64.powi(2);
let square = (60.0..110.0).contains(&fx) && (15.0..65.0).contains(&fy);
let hole = (75.0..95.0).contains(&fx) && (30.0..50.0).contains(&fy);
let on = disc || (square && !hole);
truth[y * w + x] = on;
field[y * w + x] = if on { 1.0 } else { 0.0 };
}
}
for round in [0.0, 0.5, 1.0] {
let t = trace(&field, w, h, &TraceOpts { corner_roundness: round, ..Default::default() });
assert_eq!(t.len(), 3, "disc, square, hole");
assert_eq!(t.iter().filter(|t| t.depth == 1).count(), 1, "one hole");
let r = raster(&t, w, h);
let score = iou(&truth, &r);
// Full roundness rounds the square's corners off as well.
let min = if round >= 1.0 { 0.88 } else { 0.96 };
assert!(score > min, "IoU {score} at roundness {round}");
}
// A disc stays a disc when smooth.
let disc: Vec<f32> =
(0..w * h).map(|i| if ((i % w) as f64 + 0.5 - 30.0).powi(2) + ((i / w) as f64 + 0.5 - 40.0).powi(2) <= 400.0 { 1.0 } else { 0.0 }).collect();
let t = trace(&disc, w, h, &TraceOpts { corner_roundness: 1.0, tolerance: 1.5, ..Default::default() });
let truth_disc: Vec<bool> = disc.iter().map(|v| *v > 0.5).collect();
let score = iou(&truth_disc, &raster(&t, w, h));
assert!(score > 0.97, "disc IoU {score}");
assert!(t[0].path.vertices.len() < 40, "simplified to {}", t[0].path.vertices.len());
// The square traces to 4 corners at zero roundness.
let t = trace(&field, w, h, &TraceOpts { corner_roundness: 0.0, tolerance: 1.0, ..Default::default() });
let sq = t.iter().find(|t| t.depth == 0 && t.path.vertices.iter().all(|v| v[0] > 55.0)).expect("square");
assert!(sq.path.vertices.len() <= 8, "{:?}", sq.path.vertices);
}
#[test]
fn min_area_and_threshold() {
let (w, h) = (40, 40);
let mut field = vec![0.0f32; w * h];
// A 2×2 speck and a 20×20 block at 0.6.
for y in 2..4 {
for x in 2..4 {
field[y * w + x] = 1.0;
}
}
for y in 10..30 {
for x in 10..30 {
field[y * w + x] = 0.6;
}
}
assert_eq!(trace(&field, w, h, &TraceOpts { min_area: 10.0, ..Default::default() }).len(), 1);
assert_eq!(trace(&field, w, h, &TraceOpts { min_area: 1.0, ..Default::default() }).len(), 2);
// Above 0.6 only the speck remains (and it is too small).
assert!(trace(&field, w, h, &TraceOpts { threshold: 0.7, ..Default::default() }).is_empty());
// Edge-touching shapes close along the border.
let full = vec![1.0f32; 100];
let t = trace(&full, 10, 10, &TraceOpts::default());
assert_eq!(t.len(), 1);
assert!((t[0].area.abs() - 100.0).abs() < 2.0, "{}", t[0].area);
}
}

206
crates/raster/src/cuts.rs Normal file
View File

@ -0,0 +1,206 @@
//! Scene edit (cut) detection (Layer ▸ Scene Edit Detection).
//!
//! Each frame is reduced to a [`Signature`]: a joint RGB colour histogram (4 × 4 × 4 bins,
//! normalised) and a small luma thumbnail. Two consecutive frames are compared with
//!
//! - the **histogram distance** (half the L1 distance, 0…1): large when the colour content
//! changes, insensitive to motion;
//! - the **motion-compensated difference**: block optical flow ([`crate::flow::block_flow`])
//! from the earlier thumbnail to the later one, the earlier one warped along it, and the mean
//! absolute luma error that remains. Inside a shot motion explains the change and the error
//! stays small; across a cut nothing matches.
//!
//! The dissimilarity is the geometric mean of the two (both must agree, which rejects flashes and
//! fast pans). A frame starts a new shot when its dissimilarity exceeds the absolute threshold
//! and is a clear local peak: several times the median of its neighbours (the adaptive threshold
//! of the classic shot-boundary literature), with a minimum shot length.
use rayon::prelude::*;
use crate::Image;
use crate::flow::block_flow;
const BINS: usize = 4;
const THUMB: u32 = 96;
/// The per-frame data compared by [`dissimilarity`].
#[derive(Clone, Debug, Default, PartialEq)]
pub struct Signature {
pub hist: Vec<f32>,
pub thumb: Image,
}
/// Compute a frame's signature.
pub fn signature(img: &Image) -> Signature {
let mut hist = vec![0f32; BINS * BINS * BINS];
let mut n = 0.0;
for p in &img.data {
let a = p[3];
let c = if a > 1e-6 { [p[0] / a * a.min(1.0), p[1] / a * a.min(1.0), p[2] / a * a.min(1.0)] } else { [0.0; 3] };
let b = |v: f32| ((v.clamp(0.0, 1.0) * BINS as f32) as usize).min(BINS - 1);
hist[b(c[0]) * BINS * BINS + b(c[1]) * BINS + b(c[2])] += 1.0;
n += 1.0;
}
if n > 0.0 {
for h in &mut hist {
*h /= n;
}
}
let long = img.width.max(img.height).max(1);
let thumb = if long > THUMB {
let s = THUMB as f64 / long as f64;
crate::resample(img, ((img.width as f64 * s).round() as u32).max(8), ((img.height as f64 * s).round() as u32).max(8))
} else {
img.clone()
};
Signature { hist, thumb }
}
fn luma(p: [f32; 4]) -> f32 {
0.2126 * p[0] + 0.7152 * p[1] + 0.0722 * p[2]
}
/// Dissimilarity of consecutive frames, 0 (same) … 1 (unrelated).
pub fn dissimilarity(a: &Signature, b: &Signature) -> f32 {
let hist: f32 = a.hist.iter().zip(&b.hist).map(|(x, y)| (x - y).abs()).sum::<f32>() * 0.5;
let (ta, tb) = (&a.thumb, &b.thumb);
let motion = if ta.width == tb.width && ta.height == tb.height && !ta.is_empty() {
let f = block_flow(tb, ta, 8, 4);
// For each pixel of B, where it came from in A.
let mut err = 0.0f32;
let mut n = 0.0f32;
for y in 0..tb.height {
for x in 0..tb.width {
let v = f.at(x as f64 + 0.5, y as f64 + 0.5);
let pa = ta.sample_bilinear_clamped(x as f64 + 0.5 + v[0] as f64, y as f64 + 0.5 + v[1] as f64);
err += (luma(pa) - luma(tb.data[tb.idx(x, y)])).abs();
n += 1.0;
}
}
// Mean absolute error of 0.25 (a quarter of full scale) counts as completely different.
(err / n.max(1.0) * 4.0).min(1.0)
} else {
1.0
};
(hist * motion).sqrt()
}
/// Cut detection options.
#[derive(Clone, Copy, Debug, PartialEq)]
pub struct CutOpts {
/// Absolute dissimilarity a cut must exceed (0…1).
pub threshold: f32,
/// A cut must be this many times the median of its neighbours.
pub peak_ratio: f32,
/// Neighbours on each side used for the local median.
pub window: usize,
/// Minimum frames between cuts.
pub min_shot: usize,
}
impl Default for CutOpts {
fn default() -> Self {
CutOpts { threshold: 0.25, peak_ratio: 3.0, window: 4, min_shot: 2 }
}
}
/// Dissimilarities of every frame to the previous one (`scores[0]` = 0).
pub fn scores(frames: &[Image]) -> Vec<f32> {
let sigs: Vec<Signature> = frames.par_iter().map(signature).collect();
let mut out = vec![0.0];
out.extend((1..sigs.len()).into_par_iter().map(|i| dissimilarity(&sigs[i - 1], &sigs[i])).collect::<Vec<_>>());
out
}
/// Frames that start a new shot (indices into `scores`; never 0).
pub fn detect(scores: &[f32], o: &CutOpts) -> Vec<usize> {
let mut cuts: Vec<usize> = vec![];
for i in 1..scores.len() {
let s = scores[i];
if s < o.threshold {
continue;
}
let lo = i.saturating_sub(o.window).max(1);
let hi = (i + o.window).min(scores.len() - 1);
let mut nb: Vec<f32> = (lo..=hi).filter(|&j| j != i).map(|j| scores[j]).collect();
nb.sort_by(f32::total_cmp);
let med = nb.get(nb.len() / 2).copied().unwrap_or(0.0);
if s < o.peak_ratio * med.max(0.02) {
continue;
}
match cuts.last() {
Some(&last) if i - last < o.min_shot.max(1) => {
// Keep the stronger of two close candidates.
if s > scores[last] {
cuts.pop();
cuts.push(i);
}
}
_ => cuts.push(i),
}
}
cuts
}
#[cfg(test)]
mod tests {
use super::*;
use crate::hash_noise;
/// A textured "shot": `kind` picks the palette and texture, panned by `dx` pixels.
fn frame(kind: u32, dx: f32, w: u32, h: u32) -> Image {
let mut img = Image::new(w, h);
for y in 0..h {
for x in 0..w {
let fx = x as f32 + dx;
let v = match kind {
0 => 0.5 + 0.4 * (fx * 0.15).sin() * (y as f32 * 0.11).cos(),
1 => {
if ((fx / 12.0).floor() as i32 + (y / 12) as i32) % 2 == 0 {
0.85
} else {
0.15
}
}
_ => hash_noise((fx.floor() as i64 / 4) as u32, y / 4, 7),
};
let c = match kind {
0 => [v, v * 0.6, 0.2],
1 => [0.1, v * 0.8, v],
_ => [v * 0.9, v * 0.9, v * 0.4 + 0.5],
};
img.set(x, y, [c[0], c[1], c[2], 1.0]);
}
}
img
}
#[test]
fn finds_synthetic_cuts_exactly() {
// Three shots with camera pans (12, 9 and 11 frames) and a brightness flicker inside the
// second.
let mut frames = vec![];
for i in 0..12 {
frames.push(frame(0, i as f32 * 2.0, 128, 72));
}
for i in 0..9 {
let mut f = frame(1, i as f32 * 3.0, 128, 72);
if i == 4 {
f.map_straight(|c| [c[0] * 1.08, c[1] * 1.08, c[2] * 1.08]);
}
frames.push(f);
}
for i in 0..11 {
frames.push(frame(2, -(i as f32) * 4.0, 128, 72));
}
let s = scores(&frames);
let cuts = detect(&s, &CutOpts::default());
assert_eq!(cuts, vec![12, 21], "scores {s:?}");
}
#[test]
fn a_static_shot_has_no_cuts() {
let frames: Vec<Image> = (0..8).map(|i| frame(2, i as f32, 64, 48)).collect();
assert!(detect(&scores(&frames), &CutOpts::default()).is_empty());
}
}

958
crates/raster/src/inpaint.rs Normal file
View File

@ -0,0 +1,958 @@
//! Content-Aware Fill: filling the transparent area of a layer across a range of frames.
//!
//! Written from the published methods:
//!
//! - **Membrane (harmonic) fill**: the hole takes the smoothest values that match its boundary
//! (Laplace's equation with Dirichlet boundary values), solved coarse to fine with
//! known-weighted pyramids and Jacobi relaxation. Edge Blend uses it directly; it also
//! initialises the other methods and completes flow fields.
//! - **PatchMatch** (Barnes, Shechtman, Finkelstein, Goldman 2009) nearest-neighbour fields
//! inside the **space-time completion** EM loop of Wexler, Shechtman and Irani (2007), applied
//! per frame: every patch overlapping the hole is matched to a fully known patch (random
//! initialisation, propagation, random search); hole pixels become the weighted vote of the
//! overlapping matches; repeated over a multi-scale pyramid, coarse to fine.
//! - **Flow-guided propagation** (the propagation stage of flow-guided video inpainting, e.g. Xu
//! et al. 2019 / Huang et al. 2016): optical flow between neighbouring frames is computed on
//! membrane-filled frames, *completed* inside the hole from the surrounding motion, and each
//! missing pixel follows the flow forwards and backwards until it lands on a frame where that
//! point of the scene is visible; the two candidates are blended by temporal distance.
//! **Lighting Correction** adds the membrane interpolation of the mismatch along the hole's
//! border (Poisson-style seamless blending), scaled by the chosen strength.
//!
//! Methods (After Effects' Content-Aware Fill panel):
//! - **Object**: propagation, then PatchMatch synthesis for what no frame shows; each frame's
//! synthesis starts from the previous frame's result warped along the flow, so the fill stays
//! temporally stable.
//! - **Surface**: propagation, then one synthesised fill carried through the range by the flow
//! (flat, rigid surfaces).
//! - **Edge Blend**: the membrane fill of each frame from its own border (fast; textureless
//! surfaces).
//!
//! Reference frames are complete images of given frames (for example a clean plate): they count as
//! frames where every pixel is visible, so propagation pulls their content into the hole.
use rayon::prelude::*;
use crate::Image;
use crate::flow::{Flow, block_flow};
/// Fill method.
#[derive(Clone, Copy, Debug, Default, PartialEq, Eq)]
pub enum FillMethod {
#[default]
Object,
Surface,
EdgeBlend,
}
impl FillMethod {
pub fn from_name(s: &str) -> Option<FillMethod> {
match s.to_ascii_lowercase().replace([' ', '_', '-'], "").as_str() {
"object" => Some(FillMethod::Object),
"surface" => Some(FillMethod::Surface),
"edgeblend" | "edge" => Some(FillMethod::EdgeBlend),
_ => None,
}
}
pub fn name(self) -> &'static str {
match self {
FillMethod::Object => "object",
FillMethod::Surface => "surface",
FillMethod::EdgeBlend => "edgeBlend",
}
}
}
/// Fill options.
#[derive(Clone, Copy, Debug, PartialEq)]
pub struct FillOpts {
pub method: FillMethod,
/// Lighting Correction strength 0 (None) … 1 (Strong).
pub lighting: f32,
/// Patch size (odd).
pub patch: usize,
/// Seed of the PatchMatch random search.
pub seed: u32,
}
impl Default for FillOpts {
fn default() -> Self {
FillOpts { method: FillMethod::Object, lighting: 0.0, patch: 7, seed: 1 }
}
}
// ------------------------------------------------------------------ helpers
#[inline]
fn rnd(seed: u32, a: u32, b: u32) -> u32 {
let mut h = seed.wrapping_mul(0x9E37_79B9) ^ a.wrapping_mul(0x85EB_CA6B) ^ b.wrapping_mul(0xC2B2_AE35);
h ^= h >> 16;
h = h.wrapping_mul(0x7FEB_352D);
h ^= h >> 15;
h = h.wrapping_mul(0x846C_A68B);
h ^ (h >> 16)
}
/// Grow a hole mask by `r` pixels (square structuring element, separable).
pub fn dilate(hole: &[bool], w: usize, h: usize, r: usize) -> Vec<bool> {
if r == 0 {
return hole.to_vec();
}
let mut tmp = vec![false; w * h];
for y in 0..h {
for x in 0..w {
let lo = x.saturating_sub(r);
let hi = (x + r).min(w - 1);
tmp[y * w + x] = (lo..=hi).any(|k| hole[y * w + k]);
}
}
let mut out = vec![false; w * h];
for y in 0..h {
for x in 0..w {
let lo = y.saturating_sub(r);
let hi = (y + r).min(h - 1);
out[y * w + x] = (lo..=hi).any(|k| tmp[k * w + x]);
}
}
out
}
/// Peak signal-to-noise ratio (dB) of the RGB of `a` against `b` over `region` (all pixels when
/// `None`), for values in 0…1.
pub fn psnr(a: &Image, b: &Image, region: Option<&[bool]>) -> f64 {
let (mut se, mut n) = (0.0f64, 0.0f64);
for i in 0..a.data.len().min(b.data.len()) {
if region.is_some_and(|r| !r[i]) {
continue;
}
for c in 0..3 {
let d = (a.data[i][c] - b.data[i][c]) as f64;
se += d * d;
}
n += 3.0;
}
if se <= 0.0 || n == 0.0 {
return 99.0;
}
10.0 * (1.0 / (se / n)).log10()
}
// ------------------------------------------------------------------ membrane
/// Harmonic fill of `data` (`w × h`, 4 channels) where `known` is false.
fn membrane_core(data: &mut [[f32; 4]], known: &[bool], w: usize, h: usize) {
if !known.iter().any(|k| !*k) {
return;
}
if !known.iter().any(|k| *k) {
return;
}
// Coarse level: known-weighted 2×2 averages.
if w > 8 && h > 8 {
let (cw, ch) = (w.div_ceil(2), h.div_ceil(2));
let mut cd = vec![[0.0f32; 4]; cw * ch];
let mut ck = vec![false; cw * ch];
for y in 0..ch {
for x in 0..cw {
let mut s = [0.0f32; 4];
let mut n = 0.0;
for (dx, dy) in [(0, 0), (1, 0), (0, 1), (1, 1)] {
let (fx, fy) = (2 * x + dx, 2 * y + dy);
if fx < w && fy < h && known[fy * w + fx] {
let p = data[fy * w + fx];
for c in 0..4 {
s[c] += p[c];
}
n += 1.0;
}
}
if n > 0.0 {
cd[y * cw + x] = s.map(|v| v / n);
ck[y * cw + x] = true;
}
}
}
membrane_core(&mut cd, &ck, cw, ch);
for y in 0..h {
for x in 0..w {
if !known[y * w + x] {
data[y * w + x] = cd[(y / 2).min(ch - 1) * cw + (x / 2).min(cw - 1)];
}
}
}
}
// Jacobi relaxation of the unknown pixels (in place, Gauss–Seidel order).
let iters = if w * h <= 4096 { 200 } else { 40 };
let unknown: Vec<usize> = (0..w * h).filter(|&i| !known[i]).collect();
for _ in 0..iters {
let mut delta = 0.0f32;
for &i in &unknown {
let (x, y) = (i % w, i / w);
let mut s = [0.0f32; 4];
let mut n = 0.0;
for (nx, ny) in [(x.wrapping_sub(1), y), (x + 1, y), (x, y.wrapping_sub(1)), (x, y + 1)] {
if nx < w && ny < h {
let p = data[ny * w + nx];
for c in 0..4 {
s[c] += p[c];
}
n += 1.0;
}
}
let v = s.map(|v| v / n);
delta = delta.max((v[0] - data[i][0]).abs());
data[i] = v;
}
if delta < 1e-5 {
break;
}
}
}
/// Membrane fill: `img` with the `hole` pixels replaced by the harmonic interpolation of their
/// border.
pub fn membrane_fill(img: &Image, hole: &[bool]) -> Image {
let (w, h) = (img.width as usize, img.height as usize);
let mut out = img.clone();
let known: Vec<bool> = hole.iter().map(|h| !*h).collect();
membrane_core(&mut out.data, &known, w, h);
out
}
// ------------------------------------------------------------------ PatchMatch
struct Level {
w: usize,
h: usize,
img: Vec<[f32; 4]>,
hole: Vec<bool>,
}
fn downsample(l: &Level) -> Level {
let (w, h) = (l.w.div_ceil(2), l.h.div_ceil(2));
let mut img = vec![[0.0f32; 4]; w * h];
let mut hole = vec![false; w * h];
for y in 0..h {
for x in 0..w {
let mut s = [0.0f32; 4];
let mut n = 0.0;
let mut any_hole = false;
for (dx, dy) in [(0, 0), (1, 0), (0, 1), (1, 1)] {
let (fx, fy) = ((2 * x + dx).min(l.w - 1), (2 * y + dy).min(l.h - 1));
let i = fy * l.w + fx;
if l.hole[i] {
any_hole = true;
} else {
for c in 0..4 {
s[c] += l.img[i][c];
}
n += 1.0;
}
}
hole[y * w + x] = any_hole;
if n > 0.0 {
img[y * w + x] = s.map(|v| v / n);
}
}
}
Level { w, h, img, hole }
}
/// Summed-area table of hole pixels, for "is this patch fully known" queries.
struct HoleSat {
w: usize,
s: Vec<u32>,
}
impl HoleSat {
fn new(hole: &[bool], w: usize, h: usize) -> HoleSat {
let sw = w + 1;
let mut s = vec![0u32; sw * (h + 1)];
for y in 0..h {
for x in 0..w {
s[(y + 1) * sw + x + 1] = hole[y * w + x] as u32 + s[y * sw + x + 1] + s[(y + 1) * sw + x] - s[y * sw + x];
}
}
HoleSat { w: sw, s }
}
/// Hole pixels in [x0, x1) × [y0, y1).
fn count(&self, x0: usize, y0: usize, x1: usize, y1: usize) -> u32 {
self.s[y1 * self.w + x1] + self.s[y0 * self.w + x0] - self.s[y0 * self.w + x1] - self.s[y1 * self.w + x0]
}
}
struct Pm<'a> {
w: usize,
h: usize,
r: usize,
img: &'a [[f32; 4]],
sat: HoleSat,
}
impl Pm<'_> {
fn valid_source(&self, x: i64, y: i64) -> bool {
let r = self.r as i64;
x >= r
&& y >= r
&& x + r < self.w as i64
&& y + r < self.h as i64
&& self.sat.count((x - r) as usize, (y - r) as usize, (x + r + 1) as usize, (y + r + 1) as usize) == 0
}
/// SSD between the patch at target (tx, ty) (clipped to the image) and source (sx, sy).
fn dist(&self, tx: i64, ty: i64, sx: i64, sy: i64, cap: f32) -> f32 {
let r = self.r as i64;
let mut d = 0.0f32;
let mut n = 0.0f32;
for dy in -r..=r {
let (ty2, sy2) = (ty + dy, sy + dy);
if ty2 < 0 || ty2 >= self.h as i64 {
continue;
}
for dx in -r..=r {
let tx2 = tx + dx;
if tx2 < 0 || tx2 >= self.w as i64 {
continue;
}
let a = self.img[ty2 as usize * self.w + tx2 as usize];
let b = self.img[sy2 as usize * self.w + (sx + dx) as usize];
for c in 0..4 {
let e = a[c] - b[c];
d += e * e;
}
n += 1.0;
}
if d > cap * n.max(1.0) {
return f32::INFINITY;
}
}
d / n.max(1.0)
}
}
/// One Wexler/PatchMatch level: refine `img` (hole pixels) with `em` EM rounds, starting from the
/// nearest-neighbour field `nnf` (per pixel, source centre; `None` = random).
fn pm_level(l: &mut Level, nnf: &mut Vec<Option<[i32; 2]>>, em: usize, pm_iters: usize, patch: usize, seed: u32) {
let (w, h) = (l.w, l.h);
let r = (patch / 2).min((w.min(h).saturating_sub(1)) / 2).max(1);
// Target centres: pixels whose patch touches the hole.
let near = dilate(&l.hole, w, h, r);
let targets: Vec<usize> = (0..w * h).filter(|&i| near[i]).collect();
if targets.is_empty() {
return;
}
let sat = HoleSat::new(&l.hole, w, h);
// Sources: centres of fully known patches.
let sources: Vec<[i32; 2]> = {
let pm = Pm { w, h, r, img: &l.img, sat: HoleSat::new(&l.hole, w, h) };
(0..w * h).filter(|&i| pm.valid_source((i % w) as i64, (i / w) as i64)).map(|i| [(i % w) as i32, (i / w) as i32]).collect()
};
if sources.is_empty() {
return;
}
nnf.resize(w * h, None);
for e in 0..em {
let pm = Pm { w, h, r, img: &l.img, sat: HoleSat { w: sat.w, s: sat.s.clone() } };
// Initialise / validate.
let mut field: Vec<([i32; 2], f32)> = vec![([0, 0], f32::INFINITY); w * h];
targets
.par_iter()
.map(|&i| {
let (x, y) = ((i % w) as i64, (i / w) as i64);
let s = nnf[i]
.filter(|s| pm.valid_source(s[0] as i64, s[1] as i64))
.unwrap_or_else(|| sources[rnd(seed ^ e as u32, i as u32, 0) as usize % sources.len()]);
(i, (s, pm.dist(x, y, s[0] as i64, s[1] as i64, f32::INFINITY)))
})
.collect::<Vec<_>>()
.into_iter()
.for_each(|(i, v)| field[i] = v);
let is_target = &near;
for it in 0..pm_iters {
let fwd = it % 2 == 0;
let snap = field.clone();
let rows: Vec<Vec<(usize, ([i32; 2], f32))>> = (0..h)
.into_par_iter()
.map(|y| {
let mut row: Vec<([i32; 2], f32)> = snap[y * w..(y + 1) * w].to_vec();
let mut out = vec![];
let xs: Box<dyn Iterator<Item = usize>> = if fwd { Box::new(0..w) } else { Box::new((0..w).rev()) };
for x in xs {
let i = y * w + x;
if !is_target[i] {
continue;
}
let (mut best, mut bd) = row[x];
let (xi, yi) = (x as i64, y as i64);
let try_c = |c: [i32; 2], best: &mut [i32; 2], bd: &mut f32| {
if c != *best && pm.valid_source(c[0] as i64, c[1] as i64) {
let d = pm.dist(xi, yi, c[0] as i64, c[1] as i64, *bd);
if d < *bd {
*bd = d;
*best = c;
}
}
};
// Propagation from the scan-order neighbours.
let step: i64 = if fwd { -1 } else { 1 };
let nx = xi + step;
if nx >= 0 && nx < w as i64 && is_target[y * w + nx as usize] {
let s = row[nx as usize].0;
try_c([s[0] - step as i32, s[1]], &mut best, &mut bd);
}
let ny = yi + step;
if ny >= 0 && ny < h as i64 && is_target[ny as usize * w + x] {
let s = snap[ny as usize * w + x].0;
try_c([s[0], s[1] - step as i32], &mut best, &mut bd);
}
// Random search around the best match, halving the radius.
let mut rad = w.max(h) as i64;
let mut k = 0u32;
while rad >= 1 {
let a = rnd(seed.wrapping_add(it as u32 * 977 + e as u32 * 31), i as u32, k);
let ox = (a % (2 * rad as u32 + 1)) as i64 - rad;
let oy = ((a >> 16) % (2 * rad as u32 + 1)) as i64 - rad;
try_c([(best[0] as i64 + ox) as i32, (best[1] as i64 + oy) as i32], &mut best, &mut bd);
rad /= 2;
k += 1;
}
// And a random source anywhere (helps escape bad starts).
let any = sources[rnd(seed ^ 0xABCD, i as u32, it as u32 + e as u32 * 64) as usize % sources.len()];
try_c(any, &mut best, &mut bd);
row[x] = (best, bd);
out.push((i, (best, bd)));
}
out
})
.collect();
for row in rows {
for (i, v) in row {
field[i] = v;
}
}
}
// Vote.
let mut ds: Vec<f32> = targets.iter().map(|&i| field[i].1).filter(|d| d.is_finite()).collect();
ds.sort_by(f32::total_cmp);
let sigma2 = ds.get(ds.len() * 3 / 4).copied().unwrap_or(1.0).max(1e-5);
let ri = r as i64;
let new: Vec<(usize, [f32; 4])> = (0..w * h)
.into_par_iter()
.filter(|&i| l.hole[i])
.map(|i| {
let (x, y) = ((i % w) as i64, (i / w) as i64);
let mut acc = [0.0f32; 4];
let mut wsum = 0.0f32;
for dy in -ri..=ri {
for dx in -ri..=ri {
let (px, py) = (x + dx, y + dy);
if px < 0 || py < 0 || px >= w as i64 || py >= h as i64 {
continue;
}
let (s, d) = field[py as usize * w + px as usize];
if !d.is_finite() {
continue;
}
let (sx, sy) = (s[0] as i64 - dx, s[1] as i64 - dy);
let wt = (-d / (2.0 * sigma2)).exp();
let v = l.img[sy as usize * w + sx as usize];
for c in 0..4 {
acc[c] += v[c] * wt;
}
wsum += wt;
}
}
(i, if wsum > 0.0 { acc.map(|v| v / wsum) } else { l.img[i] })
})
.collect();
for (i, v) in new {
l.img[i] = v;
}
for &i in &targets {
nnf[i] = Some(field[i].0);
}
}
}
/// Fill the `hole` of `img` with PatchMatch synthesis (multi-scale Wexler EM). `init`: a starting
/// guess for the hole (e.g. the previous frame's fill warped here); with it only the finest
/// scale is refined.
pub fn patch_fill(img: &Image, hole: &[bool], init: Option<&Image>, patch: usize, seed: u32) -> Image {
let (w, h) = (img.width as usize, img.height as usize);
if !hole.iter().any(|b| *b) || w == 0 || h == 0 {
return img.clone();
}
let patch = patch.max(3) | 1;
let fine = Level { w, h, img: img.data.clone(), hole: hole.to_vec() };
if let Some(init) = init {
let mut l = fine;
for (i, v) in l.img.iter_mut().enumerate() {
if l.hole[i] {
*v = init.data[i];
}
}
let mut nnf = vec![];
pm_level(&mut l, &mut nnf, 2, 4, patch, seed);
return Image { width: img.width, height: img.height, data: l.img };
}
// Pyramid: stop while the image is still a few patches across.
let mut levels = vec![fine];
while levels.len() < 6 {
let last = levels.last().expect("level");
if last.w.min(last.h) / 2 < patch * 4 || !last.hole.iter().any(|b| *b) {
break;
}
let d = downsample(last);
// Keep at least one fully known patch.
if !d.hole.iter().any(|b| !*b) {
break;
}
levels.push(d);
}
// Coarsest: membrane start.
let top = levels.len() - 1;
{
let l = &mut levels[top];
let known: Vec<bool> = l.hole.iter().map(|b| !*b).collect();
membrane_core(&mut l.img, &known, l.w, l.h);
}
let mut nnf: Vec<Option<[i32; 2]>> = vec![];
for k in (0..=top).rev() {
if k < top {
// Upsample the coarser result into this level's hole and the NNF (×2).
let (coarse_img, cw, ch, coarse_nnf) = {
let c = &levels[k + 1];
(c.img.clone(), c.w, c.h, nnf.clone())
};
let l = &mut levels[k];
for y in 0..l.h {
for x in 0..l.w {
let i = y * l.w + x;
if l.hole[i] {
// Bilinear sample of the coarse image.
let fx = ((x as f32 + 0.5) / 2.0 - 0.5).clamp(0.0, (cw - 1) as f32);
let fy = ((y as f32 + 0.5) / 2.0 - 0.5).clamp(0.0, (ch - 1) as f32);
let (x0, y0) = (fx.floor() as usize, fy.floor() as usize);
let (x1, y1) = ((x0 + 1).min(cw - 1), (y0 + 1).min(ch - 1));
let (tx, ty) = (fx - x0 as f32, fy - y0 as f32);
let g = |xx: usize, yy: usize| coarse_img[yy * cw + xx];
let mut v = [0.0f32; 4];
for c in 0..4 {
let top = g(x0, y0)[c] + (g(x1, y0)[c] - g(x0, y0)[c]) * tx;
let bot = g(x0, y1)[c] + (g(x1, y1)[c] - g(x0, y1)[c]) * tx;
v[c] = top + (bot - top) * ty;
}
l.img[i] = v;
}
}
}
nnf = (0..l.w * l.h)
.map(|i| {
let (x, y) = (i % l.w, i / l.w);
let ci = (y / 2).min(ch - 1) * cw + (x / 2).min(cw - 1);
coarse_nnf.get(ci).copied().flatten().map(|s| [s[0] * 2 + (x % 2) as i32, s[1] * 2 + (y % 2) as i32])
})
.collect();
}
let em = if k == top { 8 } else { 4 };
pm_level(&mut levels[k], &mut nnf, em, 5, patch, seed.wrapping_add(k as u32 * 101));
}
let l = levels.swap_remove(0);
Image { width: img.width, height: img.height, data: l.img }
}
// ------------------------------------------------------------------ flow-guided propagation
/// Block flow from `a` to `b` with the tiles that touch either hole replaced by the membrane
/// interpolation of the surrounding motion (flow completion).
fn completed_flow(a: &Image, b: &Image, hole_a: &[bool], hole_b: &[bool]) -> Flow {
let mut f = block_flow(a, b, 8, 6);
if f.v.is_empty() {
return f;
}
let (w, bs) = (a.width as usize, f.block as usize);
let mut known = vec![true; f.v.len()];
for (i, k) in known.iter_mut().enumerate() {
let (c, r) = (i % f.cols, i / f.cols);
// A tile is unreliable when any pixel within it (and a block's margin) is a hole.
let (x0, y0) = ((c * bs).saturating_sub(bs / 2), (r * bs).saturating_sub(bs / 2));
let (x1, y1) = (((c + 1) * bs + bs / 2).min(w), ((r + 1) * bs + bs / 2).min(a.height as usize));
'tile: for y in y0..y1 {
for x in x0..x1 {
if hole_a[y * w + x] || hole_b[y * w + x] {
*k = false;
break 'tile;
}
}
}
}
if known.iter().all(|k| !*k) {
f.v.iter_mut().for_each(|v| *v = [0.0; 2]);
return f;
}
let mut d: Vec<[f32; 4]> = f.v.iter().map(|v| [v[0], v[1], 0.0, 0.0]).collect();
membrane_core(&mut d, &known, f.cols, f.rows);
for (v, x) in f.v.iter_mut().zip(d) {
*v = [x[0], x[1]];
}
f
}
/// The input of a sequence fill.
pub struct FillInput {
/// Frames (premultiplied RGBA, all the same size).
pub frames: Vec<Image>,
/// Per frame, `true` where the pixel is missing.
pub holes: Vec<Vec<bool>>,
/// Complete images for some frames: (frame index, image).
pub references: Vec<(usize, Image)>,
}
/// Trace pixel (x, y) of frame `t` through `flows` (`step` = +1 forward with `fwd`, −1 backward
/// with `bwd`) to the nearest frame where it is visible: (sampled colour, frames travelled).
fn trace(
t: usize,
x: f64,
y: f64,
step: i64,
fwd: &[Flow],
bwd: &[Flow],
known: &dyn Fn(usize, i64, i64) -> bool,
frames: &[&Image],
max: usize,
) -> Option<([f32; 4], usize)> {
let n = frames.len();
let (w, h) = (frames[0].width as f64, frames[0].height as f64);
let (mut px, mut py) = (x, y);
let mut k = t as i64;
for d in 1..=max {
let v = if step > 0 {
if k as usize + 1 >= n {
return None;
}
fwd[k as usize].at(px, py)
} else {
if k == 0 {
return None;
}
bwd[k as usize - 1].at(px, py)
};
px += v[0] as f64;
py += v[1] as f64;
k += step;
if px < 0.0 || py < 0.0 || px >= w || py >= h {
return None;
}
// Visible when the bilinear support is known.
let (x0, y0) = ((px - 0.5).floor() as i64, (py - 0.5).floor() as i64);
let ok = [(x0, y0), (x0 + 1, y0), (x0, y0 + 1), (x0 + 1, y0 + 1)]
.iter()
.all(|&(xx, yy)| known(k as usize, xx.clamp(0, w as i64 - 1), yy.clamp(0, h as i64 - 1)));
if ok {
return Some((frames[k as usize].sample_bilinear_clamped(px, py), d));
}
}
None
}
/// Fill a range of frames. `progress(done, total)` returns `false` to cancel (→ `None`).
pub fn fill_sequence(input: &FillInput, o: &FillOpts, progress: &mut dyn FnMut(usize, usize) -> bool) -> Option<Vec<Image>> {
let n = input.frames.len();
if n == 0 {
return Some(vec![]);
}
let (w, h) = (input.frames[0].width as usize, input.frames[0].height as usize);
// Steps: flows (1) + propagation (1) + per-frame synthesis (n).
let total = n + 2;
if o.method == FillMethod::EdgeBlend {
let mut out = Vec::with_capacity(n);
for (i, (f, hl)) in input.frames.iter().zip(&input.holes).enumerate() {
out.push(opaque_fill(membrane_fill(f, hl), hl));
if !progress(i + 1, n) {
return None;
}
}
return Some(out);
}
// Reference frames make their frame fully known.
let mut frames: Vec<Image> = input.frames.clone();
let mut holes: Vec<Vec<bool>> = input.holes.clone();
let mut is_ref = vec![false; n];
for (k, img) in &input.references {
if *k < n && img.width as usize == w && img.height as usize == h {
frames[*k] = img.clone();
holes[*k] = vec![false; w * h];
is_ref[*k] = true;
}
}
// 1. Flows between membrane-filled neighbours, completed inside the holes.
let pre: Vec<Image> = frames.par_iter().zip(&holes).map(|(f, hl)| membrane_fill(f, hl)).collect();
let pairs: Vec<(Flow, Flow)> = (0..n.saturating_sub(1))
.into_par_iter()
.map(|i| (completed_flow(&pre[i], &pre[i + 1], &holes[i], &holes[i + 1]), completed_flow(&pre[i + 1], &pre[i], &holes[i + 1], &holes[i])))
.collect();
let (fwd, bwd): (Vec<Flow>, Vec<Flow>) = pairs.into_iter().unzip();
if !progress(1, total) {
return None;
}
// 2. Propagation (and the lighting-correction border ring).
let refs: Vec<&Image> = frames.iter().collect();
let known = |k: usize, x: i64, y: i64| !holes[k][y as usize * w + x as usize];
let ring_r = 2usize;
let propagated: Vec<(Image, Vec<bool>)> = (0..n)
.into_par_iter()
.map(|t| {
let mut img = frames[t].clone();
let mut left = holes[t].clone();
if is_ref[t] || !holes[t].iter().any(|b| *b) {
return (img, left);
}
let sample = |i: usize| -> Option<[f32; 4]> {
let (x, y) = ((i % w) as f64 + 0.5, (i / w) as f64 + 0.5);
let f = trace(t, x, y, 1, &fwd, &bwd, &known, &refs, n);
let b = trace(t, x, y, -1, &fwd, &bwd, &known, &refs, n);
match (f, b) {
(Some((a, da)), Some((c, dc))) => {
let (wa, wc) = (1.0 / da as f32, 1.0 / dc as f32);
Some(std::array::from_fn(|k| (a[k] * wa + c[k] * wc) / (wa + wc)))
}
(Some((a, _)), None) | (None, Some((a, _))) => Some(a),
_ => None,
}
};
let filled: Vec<(usize, [f32; 4])> = (0..w * h).filter(|&i| holes[t][i]).filter_map(|i| sample(i).map(|v| (i, v))).collect();
for (i, v) in &filled {
img.data[*i] = *v;
left[*i] = false;
}
// Lighting correction: membrane of the border mismatch, added to propagated pixels.
if o.lighting > 0.0 && !filled.is_empty() {
let grown = dilate(&holes[t], w, h, ring_r);
let mut diff = vec![[0.0f32; 4]; w * h];
let mut dk = vec![true; w * h];
for i in 0..w * h {
if holes[t][i] {
dk[i] = false;
} else if grown[i]
&& let Some(v) = sample(i)
{
let a = frames[t].data[i];
diff[i] = std::array::from_fn(|c| if c < 3 { a[c] - v[c] } else { 0.0 });
}
}
membrane_core(&mut diff, &dk, w, h);
for (i, _) in &filled {
for c in 0..3 {
img.data[*i][c] += o.lighting * diff[*i][c];
}
}
}
(img, left)
})
.collect();
if !progress(2, total) {
return None;
}
// 3. Synthesis of what no frame shows.
let mut out: Vec<Option<Image>> = vec![None; n];
let first = (0..n).find(|&t| propagated[t].1.iter().any(|b| *b));
let Some(t0) = first else {
return Some(propagated.into_iter().zip(&holes).map(|((img, _), hl)| opaque_fill(img, hl)).collect());
};
// Frames before t0 are complete already.
for (t, slot) in out.iter_mut().enumerate().take(t0) {
*slot = Some(propagated[t].0.clone());
}
let mut done = 2;
out[t0] = Some(patch_fill(&propagated[t0].0, &propagated[t0].1, None, o.patch, o.seed));
done += 1;
if !progress(done, total) {
return None;
}
for t in t0 + 1..n {
let (img, left) = &propagated[t];
if !left.iter().any(|b| *b) {
out[t] = Some(img.clone());
} else {
// Warp the previous result here along the backward flow (t → t − 1).
let prev = out[t - 1].as_ref().expect("previous frame");
let mut warped = img.clone();
for i in 0..w * h {
if left[i] {
let (x, y) = ((i % w) as f64 + 0.5, (i / w) as f64 + 0.5);
let v = bwd[t - 1].at(x, y);
warped.data[i] = prev.sample_bilinear_clamped(x + v[0] as f64, y + v[1] as f64);
}
}
out[t] = Some(match o.method {
FillMethod::Surface => warped,
_ => patch_fill(img, left, Some(&warped), o.patch, o.seed),
});
}
done += 1;
if !progress(done, total) {
return None;
}
}
Some(out.into_iter().zip(&holes).map(|(img, hl)| opaque_fill(img.unwrap_or_default(), hl)).collect())
}
/// Filled pixels are opaque.
fn opaque_fill(mut img: Image, hole: &[bool]) -> Image {
for (p, hl) in img.data.iter_mut().zip(hole) {
if *hl && p[3] > 1e-6 {
let a = p[3];
*p = [p[0] / a, p[1] / a, p[2] / a, 1.0];
}
}
img
}
#[cfg(test)]
mod tests {
use super::*;
use crate::hash_noise;
/// A textured background (periodic stripes + blocky noise), shifted by `dx`.
fn background(w: u32, h: u32, dx: f32) -> Image {
let mut img = Image::new(w, h);
for y in 0..h {
for x in 0..w {
let fx = x as f32 + dx;
let s = 0.5 + 0.3 * (fx * std::f32::consts::TAU / 10.0).sin();
let n = hash_noise(((fx.floor() as i64).rem_euclid(20) / 4) as u32, (y % 20) / 4, 3);
img.set(x, y, [s, 0.3 + 0.4 * n, 0.6 - 0.3 * s * n, 1.0]);
}
}
img
}
fn cut(img: &Image, hole: &[bool]) -> Image {
let mut o = img.clone();
for (p, h) in o.data.iter_mut().zip(hole) {
if *h {
*p = [0.0; 4];
}
}
o
}
fn rect_hole(w: u32, h: u32, x0: u32, y0: u32, s: u32) -> Vec<bool> {
(0..w * h).map(|i| (x0..x0 + s).contains(&(i % w)) && (y0..y0 + s).contains(&(i / w))).collect()
}
#[test]
fn membrane_matches_a_linear_ramp() {
let (w, h) = (40u32, 30u32);
let mut img = Image::new(w, h);
for y in 0..h {
for x in 0..w {
let v = x as f32 / w as f32;
img.set(x, y, [v, v, v, 1.0]);
}
}
let hole = rect_hole(w, h, 10, 8, 14);
let f = membrane_fill(&cut(&img, &hole), &hole);
assert!(psnr(&f, &img, Some(&hole)) > 40.0, "{}", psnr(&f, &img, Some(&hole)));
}
#[test]
fn patchmatch_beats_the_membrane_on_texture() {
let (w, h) = (96u32, 72u32);
let gt = background(w, h, 0.0);
let hole = rect_hole(w, h, 40, 28, 16);
let src = cut(&gt, &hole);
let mem = membrane_fill(&src, &hole);
let pm = patch_fill(&src, &hole, None, 7, 1);
let (pm_db, mem_db) = (psnr(&pm, &gt, Some(&hole)), psnr(&mem, &gt, Some(&hole)));
assert!(pm_db > mem_db + 3.0 && pm_db > 20.0, "PatchMatch {pm_db:.1} dB vs membrane {mem_db:.1} dB");
// Known pixels are untouched.
for i in 0..hole.len() {
if !hole[i] {
assert_eq!(pm.data[i], gt.data[i]);
}
}
}
/// A moving object removed from a (static or panning) textured background.
fn sequence(n: usize, pan: f32) -> (FillInput, Vec<Image>) {
let (w, h) = (96u32, 64u32);
let mut frames = vec![];
let mut holes = vec![];
let mut gts = vec![];
for t in 0..n {
let gt = background(w, h, t as f32 * pan);
let hole = rect_hole(w, h, 10 + 8 * t as u32, 22, 18);
frames.push(cut(&gt, &hole));
holes.push(hole);
gts.push(gt);
}
(FillInput { frames, holes, references: vec![] }, gts)
}
#[test]
fn object_fill_recovers_background_and_is_stable() {
for pan in [0.0, 1.0] {
let (input, gts) = sequence(8, pan);
let out = fill_sequence(&input, &FillOpts::default(), &mut |_, _| true).expect("fill");
for t in 0..out.len() {
let db = psnr(&out[t], &gts[t], Some(&input.holes[t]));
assert!(db > 26.0, "frame {t} pan {pan}: {db:.1} dB");
assert_eq!(out[t].data[0][3], 1.0);
}
if pan == 0.0 {
// Temporal consistency: on a static background the filled pixels agree between
// frames.
for t in 1..out.len() {
let both: Vec<bool> = input.holes[t].iter().zip(&input.holes[t - 1]).map(|(a, b)| *a || *b).collect();
let db = psnr(&out[t], &out[t - 1], Some(&both));
assert!(db > 30.0, "frames {}–{t} differ: {db:.1} dB", t - 1);
}
}
}
}
#[test]
fn surface_and_edge_blend_and_references() {
// A hole that never moves: nothing to propagate, so synthesis must carry it.
let (w, h) = (80u32, 60u32);
let gt = background(w, h, 0.0);
let hole = rect_hole(w, h, 30, 20, 14);
let input = FillInput { frames: vec![cut(&gt, &hole); 4], holes: vec![hole.clone(); 4], references: vec![] };
let surf = fill_sequence(&input, &FillOpts { method: FillMethod::Surface, ..Default::default() }, &mut |_, _| true).expect("fill");
for t in 1..4 {
// Surface carries one fill through: identical frames.
assert!(psnr(&surf[t], &surf[0], Some(&hole)) > 40.0);
}
let edge = fill_sequence(&input, &FillOpts { method: FillMethod::EdgeBlend, ..Default::default() }, &mut |_, _| true).expect("fill");
assert!(psnr(&edge[0], &membrane_fill(&input.frames[0], &hole), Some(&hole)) > 60.0);
// A clean reference frame recovers the background exactly.
let mut with_ref = FillInput { frames: input.frames.clone(), holes: input.holes.clone(), references: vec![(2, gt.clone())] };
let out = fill_sequence(&with_ref, &FillOpts::default(), &mut |_, _| true).expect("fill");
for (t, f) in out.iter().enumerate() {
assert!(psnr(f, &gt, Some(&hole)) > 35.0, "frame {t}: {:.1}", psnr(f, &gt, Some(&hole)));
}
// Cancel.
with_ref.references.clear();
assert!(fill_sequence(&with_ref, &FillOpts::default(), &mut |d, _| d < 2).is_none());
}
#[test]
fn lighting_correction_absorbs_a_brightness_change() {
// The background brightens over time; propagated pixels come from darker frames.
let (mut input, mut gts) = sequence(6, 0.0);
for t in 0..6 {
let k = 1.0 + 0.06 * t as f32;
for img in [&mut input.frames[t], &mut gts[t]] {
for p in &mut img.data {
for c in 0..3 {
p[c] *= k;
}
}
}
}
let plain = fill_sequence(&input, &FillOpts::default(), &mut |_, _| true).expect("fill");
let lit = fill_sequence(&input, &FillOpts { lighting: 1.0, ..Default::default() }, &mut |_, _| true).expect("fill");
let score = |o: &[Image]| (0..6).map(|t| psnr(&o[t], &gts[t], Some(&input.holes[t]))).sum::<f64>();
assert!(score(&lit) > score(&plain), "lit {} vs plain {}", score(&lit), score(&plain));
}
}

3
crates/raster/src/lib.rs
View File

@ -5,7 +5,10 @@
pub mod blur;
pub mod channels3d;
pub mod cuts;
pub mod flow;
pub mod inpaint;
pub mod scopes;
pub mod warp;
use effectcraft_color::{BlendMode, blend_pixel};

552
crates/raster/src/scopes.rs Normal file
View File

@ -0,0 +1,552 @@
//! Video scopes (the Lumetri Scopes panel): waveform (RGB, Luma, YC), vectorscope (YUV, HLS),
//! histogram and parade (RGB, YUV), computed from a frame as density plots.
//!
//! The maths is the textbook broadcast definition: luma `Y = Kr·R + Kg·G + Kb·B` with the
//! ITU-R BT.601 / BT.709 / BT.2020 coefficients, colour differences `Cb = (B − Y) / (2(1 − Kb))`
//! and `Cr = (R − Y) / (2(1 − Kr))` (each in −0.5…0.5). Values are read straight (un-premultiplied)
//! from the frame. Every scope is a [`Scope`]: a `width × height` grid of three density planes
//! (one per trace colour) normalised so the busiest cell is 1, ready to be tinted and drawn.
//!
//! Scales: **8-bit** maps 0…1 to the plot (0…255 graticule); **float** maps −0.25…1.25 so
//! super-white and sub-black values stay visible. **Clamp** limits the signal to 0…1 first.
use rayon::prelude::*;
use crate::Image;
/// Luma / colour-difference coefficients.
#[derive(Clone, Copy, Debug, Default, PartialEq, Eq)]
pub enum ColorStandard {
Rec601,
#[default]
Rec709,
Rec2020,
}
impl ColorStandard {
pub const ALL: [ColorStandard; 3] = [ColorStandard::Rec601, ColorStandard::Rec709, ColorStandard::Rec2020];
/// (Kr, Kb).
pub fn coefficients(self) -> (f32, f32) {
match self {
ColorStandard::Rec601 => (0.299, 0.114),
ColorStandard::Rec709 => (0.2126, 0.0722),
ColorStandard::Rec2020 => (0.2627, 0.0593),
}
}
pub fn name(self) -> &'static str {
match self {
ColorStandard::Rec601 => "rec601",
ColorStandard::Rec709 => "rec709",
ColorStandard::Rec2020 => "rec2020",
}
}
pub fn from_name(s: &str) -> Option<ColorStandard> {
match s.to_ascii_lowercase().replace([' ', '.', '-', '_'], "").as_str() {
"rec601" | "601" | "bt601" => Some(ColorStandard::Rec601),
"rec709" | "709" | "bt709" => Some(ColorStandard::Rec709),
"rec2020" | "2020" | "bt2020" => Some(ColorStandard::Rec2020),
_ => None,
}
}
/// (Y, Cb, Cr) of a straight RGB colour.
pub fn ycbcr(self, rgb: [f32; 3]) -> [f32; 3] {
let (kr, kb) = self.coefficients();
let y = kr * rgb[0] + (1.0 - kr - kb) * rgb[1] + kb * rgb[2];
[y, (rgb[2] - y) / (2.0 * (1.0 - kb)), (rgb[0] - y) / (2.0 * (1.0 - kr))]
}
}
/// Which scope.
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
pub enum ScopeKind {
WaveformRgb,
WaveformLuma,
/// Luma (plane 0) and chroma amplitude `√(Cb² + Cr²)` (plane 2).
WaveformYc,
/// Cb across, Cr up.
VectorscopeYuv,
/// Hue around the circle, saturation outwards (HLS).
VectorscopeHls,
/// R, G, B histograms (bins across, counts up).
Histogram,
ParadeRgb,
ParadeYuv,
}
impl ScopeKind {
pub const ALL: [ScopeKind; 8] = [
ScopeKind::WaveformRgb,
ScopeKind::WaveformLuma,
ScopeKind::WaveformYc,
ScopeKind::VectorscopeYuv,
ScopeKind::VectorscopeHls,
ScopeKind::Histogram,
ScopeKind::ParadeRgb,
ScopeKind::ParadeYuv,
];
pub fn name(self) -> &'static str {
match self {
ScopeKind::WaveformRgb => "waveformRgb",
ScopeKind::WaveformLuma => "waveformLuma",
ScopeKind::WaveformYc => "waveformYc",
ScopeKind::VectorscopeYuv => "vectorscopeYuv",
ScopeKind::VectorscopeHls => "vectorscopeHls",
ScopeKind::Histogram => "histogram",
ScopeKind::ParadeRgb => "paradeRgb",
ScopeKind::ParadeYuv => "paradeYuv",
}
}
pub fn label(self) -> &'static str {
match self {
ScopeKind::WaveformRgb => "Waveform (RGB)",
ScopeKind::WaveformLuma => "Waveform (Luma)",
ScopeKind::WaveformYc => "Waveform (YC)",
ScopeKind::VectorscopeYuv => "Vectorscope YUV",
ScopeKind::VectorscopeHls => "Vectorscope HLS",
ScopeKind::Histogram => "Histogram",
ScopeKind::ParadeRgb => "Parade (RGB)",
ScopeKind::ParadeYuv => "Parade (YUV)",
}
}
pub fn from_name(s: &str) -> Option<ScopeKind> {
let n = s.to_ascii_lowercase().replace([' ', '(', ')', '_', '-'], "");
Self::ALL.iter().copied().find(|k| k.name().to_ascii_lowercase() == n || k.label().to_ascii_lowercase().replace([' ', '(', ')'], "") == n)
}
/// Vectorscopes are square.
pub fn is_vectorscope(self) -> bool {
matches!(self, ScopeKind::VectorscopeYuv | ScopeKind::VectorscopeHls)
}
}
/// Scope options.
#[derive(Clone, Copy, Debug, PartialEq)]
pub struct ScopeOpts {
pub standard: ColorStandard,
/// Float scale (−0.25…1.25) instead of 8-bit (0…1).
pub float: bool,
/// Clamp the signal to 0…1 before plotting.
pub clamp: bool,
/// Plot size.
pub width: u32,
pub height: u32,
}
impl Default for ScopeOpts {
fn default() -> Self {
ScopeOpts { standard: ColorStandard::Rec709, float: false, clamp: true, width: 256, height: 256 }
}
}
impl ScopeOpts {
/// Signal range shown on the value axis.
pub fn range(&self) -> (f32, f32) {
if self.float { (-0.25, 1.25) } else { (0.0, 1.0) }
}
/// Plot row (0 = top) of a value, or `None` when off-scale.
fn row(&self, v: f32) -> Option<usize> {
let (lo, hi) = self.range();
let f = (v - lo) / (hi - lo);
if !(0.0..=1.0).contains(&f) {
return None;
}
let r = ((1.0 - f) * (self.height - 1) as f32).round() as usize;
Some(r.min(self.height as usize - 1))
}
}
/// A computed scope: three density planes (`planes[c][y * width + x]`, 0…1).
#[derive(Clone, Debug, Default, PartialEq)]
pub struct Scope {
pub width: u32,
pub height: u32,
pub planes: [Vec<f32>; 3],
/// Raw counts before normalisation (same layout), for statistics and tests.
pub counts: [Vec<u32>; 3],
}
impl Scope {
fn new(w: u32, h: u32) -> Scope {
let n = w as usize * h as usize;
Scope { width: w, height: h, planes: [vec![0.0; n], vec![0.0; n], vec![0.0; n]], counts: [vec![0; n], vec![0; n], vec![0; n]] }
}
fn merge(mut self, o: Scope) -> Scope {
for c in 0..3 {
for (a, b) in self.counts[c].iter_mut().zip(&o.counts[c]) {
*a += *b;
}
}
self
}
fn normalise(&mut self) {
// Square-root density so faint traces stay visible next to dense ones.
let max = self.counts.iter().flat_map(|p| p.iter()).copied().max().unwrap_or(0).max(1) as f32;
for c in 0..3 {
self.planes[c] = self.counts[c].iter().map(|&n| (n as f32 / max).sqrt()).collect();
}
}
/// Total count of plane `c`.
pub fn total(&self, c: usize) -> u64 {
self.counts[c].iter().map(|&n| n as u64).sum()
}
/// The cell of plane `c` with the most samples, as (x, y).
pub fn peak(&self, c: usize) -> (u32, u32) {
let (i, _) = self.counts[c].iter().enumerate().max_by_key(|(_, n)| **n).unwrap_or((0, &0));
(i as u32 % self.width.max(1), i as u32 / self.width.max(1))
}
/// Straight RGBA8 picture: plane 0 red, 1 green, 2 blue (`mono`: all planes white).
pub fn to_rgba8(&self, mono: bool) -> Vec<u8> {
let n = self.width as usize * self.height as usize;
let mut out = vec![0u8; n * 4];
for i in 0..n {
let (r, g, b) = (self.planes[0][i], self.planes[1][i], self.planes[2][i]);
let (r, g, b) = if mono {
let m = r.max(g).max(b);
(m, m, m)
} else {
(r, g, b)
};
let a = r.max(g).max(b);
out[i * 4] = (r.min(1.0) * 255.0) as u8;
out[i * 4 + 1] = (g.min(1.0) * 255.0) as u8;
out[i * 4 + 2] = (b.min(1.0) * 255.0) as u8;
out[i * 4 + 3] = (a.min(1.0) * 255.0) as u8;
}
out
}
}
/// Straight RGB of pixel `i` (transparent pixels read as black), clamped when asked.
#[inline]
fn rgb(img: &Image, i: usize, clamp: bool) -> [f32; 3] {
let p = img.data[i];
let mut c = if p[3] > 1e-6 { [p[0] / p[3], p[1] / p[3], p[2] / p[3]] } else { [0.0; 3] };
// Composite over black, as a monitor shows the frame.
let a = p[3].clamp(0.0, 1.0);
for v in &mut c {
*v *= a;
if clamp {
*v = v.clamp(0.0, 1.0);
}
}
c
}
/// Compute a scope of `img`.
pub fn compute(img: &Image, kind: ScopeKind, o: &ScopeOpts) -> Scope {
let (w, h) = (o.width.max(8), o.height.max(8));
let o = ScopeOpts { width: w, height: h, ..*o };
if img.is_empty() {
return Scope::new(w, h);
}
let iw = img.width as usize;
let rows: Vec<usize> = (0..img.height as usize).collect();
let mut s = rows
.par_chunks(16)
.map(|chunk| {
let mut s = Scope::new(w, h);
for &y in chunk {
for x in 0..iw {
plot(&mut s, kind, &o, rgb(img, y * iw + x, o.clamp), x as f32 / iw as f32);
}
}
s
})
.reduce(|| Scope::new(w, h), Scope::merge);
s.normalise();
s
}
#[inline]
fn bump(s: &mut Scope, c: usize, x: usize, y: usize) {
let i = y * s.width as usize + x.min(s.width as usize - 1);
s.counts[c][i] += 1;
}
/// Add one pixel (straight RGB, at horizontal position `fx` 0…1 in the frame).
fn plot(s: &mut Scope, kind: ScopeKind, o: &ScopeOpts, c: [f32; 3], fx: f32) {
let w = s.width as usize;
let col = |fx: f32, x0: usize, span: usize| x0 + ((fx * span as f32) as usize).min(span - 1);
match kind {
ScopeKind::WaveformRgb => {
let x = col(fx, 0, w);
for (ch, v) in c.iter().enumerate() {
if let Some(y) = o.row(*v) {
bump(s, ch, x, y);
}
}
}
ScopeKind::WaveformLuma => {
let x = col(fx, 0, w);
let yv = o.standard.ycbcr(c)[0];
if let Some(y) = o.row(yv) {
for ch in 0..3 {
bump(s, ch, x, y);
}
}
}
ScopeKind::WaveformYc => {
let x = col(fx, 0, w);
let [yv, cb, cr] = o.standard.ycbcr(c);
if let Some(y) = o.row(yv) {
bump(s, 0, x, y);
bump(s, 1, x, y);
}
if let Some(y) = o.row((cb * cb + cr * cr).sqrt()) {
bump(s, 2, x, y);
}
}
ScopeKind::ParadeRgb | ScopeKind::ParadeYuv => {
let third = (w / 3).max(1);
let vals = if kind == ScopeKind::ParadeRgb {
c
} else {
let [y, cb, cr] = o.standard.ycbcr(c);
// Colour differences centred on mid-scale.
[y, cb + 0.5, cr + 0.5]
};
for (ch, v) in vals.iter().enumerate() {
let x = col(fx, ch * third, third);
if let Some(y) = o.row(*v) {
// RGB parade: tinted per channel; YUV: Y white, Cb blue, Cr red.
match (kind, ch) {
(ScopeKind::ParadeYuv, 0) => {
for p in 0..3 {
bump(s, p, x, y);
}
}
(ScopeKind::ParadeYuv, 1) => bump(s, 2, x, y),
(ScopeKind::ParadeYuv, _) => bump(s, 0, x, y),
_ => bump(s, ch, x, y),
}
}
}
}
ScopeKind::VectorscopeYuv => {
let [_, cb, cr] = o.standard.ycbcr(c);
if let Some((x, y)) = vec_cell(s, cb, cr) {
for ch in 0..3 {
bump(s, ch, x, y);
}
}
}
ScopeKind::VectorscopeHls => {
let (hue, _, sat) = hls(c);
let a = hue.to_radians();
if let Some((x, y)) = vec_cell(s, sat * 0.5 * a.cos(), sat * 0.5 * a.sin()) {
for ch in 0..3 {
bump(s, ch, x, y);
}
}
}
ScopeKind::Histogram => {
let (lo, hi) = o.range();
for (ch, v) in c.iter().enumerate() {
let f = (v - lo) / (hi - lo);
if (0.0..=1.0).contains(&f) {
let x = ((f * (w - 1) as f32).round() as usize).min(w - 1);
// Histograms keep their counts in row 0; the bars are drawn by `histogram_bars`.
bump(s, ch, x, 0);
}
}
}
}
}
/// Vectorscope cell of (u, v) in −0.5…0.5 (the disc fills the plot; +v up).
fn vec_cell(s: &Scope, u: f32, v: f32) -> Option<(usize, usize)> {
let (w, h) = (s.width as f32, s.height as f32);
let x = ((u + 0.5) * (w - 1.0)).round();
let y = ((0.5 - v) * (h - 1.0)).round();
if x < 0.0 || y < 0.0 || x >= w || y >= h {
return None;
}
Some((x as usize, y as usize))
}
/// (hue degrees, lightness, saturation) of straight RGB (HLS).
pub fn hls(c: [f32; 3]) -> (f32, f32, f32) {
let (r, g, b) = (c[0], c[1], c[2]);
let max = r.max(g).max(b);
let min = r.min(g).min(b);
let l = (max + min) * 0.5;
let d = max - min;
if d < 1e-6 {
return (0.0, l, 0.0);
}
let s = if l > 0.5 { d / (2.0 - max - min).max(1e-6) } else { d / (max + min).max(1e-6) };
let h = if max == r {
((g - b) / d).rem_euclid(6.0)
} else if max == g {
(b - r) / d + 2.0
} else {
(r - g) / d + 4.0
};
(h * 60.0, l, s.min(1.0))
}
/// Histogram bars: per channel, `bins` counts across the value range (counts in a
/// [`ScopeKind::Histogram`] scope's first row).
pub fn histogram_bars(s: &Scope) -> [Vec<u32>; 3] {
let w = s.width as usize;
[s.counts[0][..w].to_vec(), s.counts[1][..w].to_vec(), s.counts[2][..w].to_vec()]
}
/// Where a vectorscope target (75 % colour bars' primaries and secondaries) sits, as (u, v) in
/// −0.5…0.5, for the graticule.
pub fn vector_targets(standard: ColorStandard, hls_mode: bool) -> Vec<(&'static str, [f32; 2])> {
let cols: [(&str, [f32; 3]); 6] = [
("R", [0.75, 0.0, 0.0]),
("Mg", [0.75, 0.0, 0.75]),
("B", [0.0, 0.0, 0.75]),
("Cy", [0.0, 0.75, 0.75]),
("G", [0.0, 0.75, 0.0]),
("Yl", [0.75, 0.75, 0.0]),
];
cols.iter()
.map(|(n, c)| {
if hls_mode {
let (h, _, s) = hls(*c);
let a = h.to_radians();
(*n, [s * 0.5 * a.cos(), s * 0.5 * a.sin()])
} else {
let [_, cb, cr] = standard.ycbcr(*c);
(*n, [cb, cr])
}
})
.collect()
}
#[cfg(test)]
mod tests {
use super::*;
fn solid(w: u32, h: u32, c: [f32; 3]) -> Image {
Image::filled(w, h, [c[0], c[1], c[2], 1.0])
}
#[test]
fn waveform_of_flat_grey_is_one_line() {
let img = solid(64, 32, [0.5, 0.5, 0.5]);
let o = ScopeOpts { width: 64, height: 101, ..Default::default() };
let s = compute(&img, ScopeKind::WaveformLuma, &o);
// 0.5 sits on row 50 of 0..=100; every sample lands there.
assert_eq!(s.total(0), 64 * 32);
for x in 0..64 {
assert_eq!(s.counts[0][50 * 64 + x], 32);
}
}
#[test]
fn waveform_rgb_ramps_follow_columns() {
// Horizontal ramp: column x has value x / (w - 1).
let w = 101;
let mut img = Image::new(w, 4);
for y in 0..4 {
for x in 0..w {
let v = x as f32 / (w - 1) as f32;
img.set(x, y, [v, 1.0 - v, 0.0, 1.0]);
}
}
let o = ScopeOpts { width: w, height: 101, ..Default::default() };
let s = compute(&img, ScopeKind::WaveformRgb, &o);
for x in 0..w as usize {
// Red: row 100 - x; green: row x.
assert_eq!(s.counts[0][(100 - x) * w as usize + x], 4, "red column {x}");
assert_eq!(s.counts[1][x * w as usize + x], 4, "green column {x}");
}
}
#[test]
fn luma_coefficients_per_standard() {
let red = [1.0, 0.0, 0.0];
assert!((ColorStandard::Rec601.ycbcr(red)[0] - 0.299).abs() < 1e-6);
assert!((ColorStandard::Rec709.ycbcr(red)[0] - 0.2126).abs() < 1e-6);
assert!((ColorStandard::Rec2020.ycbcr(red)[0] - 0.2627).abs() < 1e-6);
// Pure red has Cr = +0.5 in every standard, Cb negative.
for st in ColorStandard::ALL {
let [_, cb, cr] = st.ycbcr(red);
assert!((cr - 0.5).abs() < 1e-6);
assert!(cb < 0.0);
// White and grey have no chroma.
let [_, cb, cr] = st.ycbcr([0.6, 0.6, 0.6]);
assert!(cb.abs() < 1e-6 && cr.abs() < 1e-6);
}
}
#[test]
fn vectorscope_places_red_on_its_target() {
let img = solid(16, 16, [1.0, 0.0, 0.0]);
let o = ScopeOpts { width: 101, height: 101, ..Default::default() };
let s = compute(&img, ScopeKind::VectorscopeYuv, &o);
let [_, cb, cr] = ColorStandard::Rec709.ycbcr([1.0, 0.0, 0.0]);
let (x, y) = s.peak(0);
assert_eq!(x, ((cb + 0.5) * 100.0).round() as u32);
assert_eq!(y, ((0.5 - cr) * 100.0).round() as u32);
// Grey sits in the centre.
let g = compute(&solid(8, 8, [0.3, 0.3, 0.3]), ScopeKind::VectorscopeYuv, &o);
assert_eq!(g.peak(0), (50, 50));
// HLS: fully saturated green at 120° on the rim.
let h = compute(&solid(8, 8, [0.0, 1.0, 0.0]), ScopeKind::VectorscopeHls, &o);
let a = 120f32.to_radians();
assert_eq!(h.peak(0), (((0.5 * a.cos() + 0.5) * 100.0).round() as u32, ((0.5 - 0.5 * a.sin()) * 100.0).round() as u32));
}
#[test]
fn histogram_counts_and_clamp() {
// Half black, half white.
let mut img = Image::new(10, 10);
for y in 0..10 {
for x in 0..10 {
let v = if x < 5 { 0.0 } else { 1.0 };
img.set(x, y, [v, v, v, 1.0]);
}
}
let o = ScopeOpts { width: 256, height: 64, ..Default::default() };
let s = compute(&img, ScopeKind::Histogram, &o);
let bars = histogram_bars(&s);
for b in &bars {
assert_eq!(b[0], 50);
assert_eq!(b[255], 50);
assert_eq!(b.iter().sum::<u32>(), 100);
}
// Super-white: off scale in 8-bit unclamped, on the float scale, pinned to 1 clamped.
let hot = solid(4, 4, [1.2, 1.2, 1.2]);
let unclamped = ScopeOpts { clamp: false, ..o };
assert_eq!(compute(&hot, ScopeKind::Histogram, &unclamped).total(0), 0);
let fl = ScopeOpts { clamp: false, float: true, ..o };
let f = histogram_bars(&compute(&hot, ScopeKind::Histogram, &fl));
let expect = (((1.2 + 0.25) / 1.5) * 255.0f32).round() as usize;
assert_eq!(f[0][expect], 16);
let c = histogram_bars(&compute(&hot, ScopeKind::Histogram, &o));
assert_eq!(c[0][255], 16);
}
#[test]
fn parade_splits_channels_into_thirds() {
let img = solid(30, 3, [1.0, 0.5, 0.0]);
let o = ScopeOpts { width: 90, height: 101, ..Default::default() };
let s = compute(&img, ScopeKind::ParadeRgb, &o);
// Red trace in the first third at the top, green mid-height in the second, blue at the
// bottom of the third.
let (rx, ry) = s.peak(0);
let (gx, gy) = s.peak(1);
let (bx, by) = s.peak(2);
assert!(rx < 30 && ry == 0);
assert!((30..60).contains(&gx) && gy == 50);
assert!((60..90).contains(&bx) && by == 100);
let yuv = compute(&solid(30, 3, [0.5, 0.5, 0.5]), ScopeKind::ParadeYuv, &o);
// Grey: Y = 0.5 and both colour differences at mid-scale.
assert_eq!(yuv.peak(1).1, 50);
assert_eq!(yuv.peak(2).1, 50);
}
#[test]
fn yc_waveform_has_chroma_trace() {
let img = solid(8, 8, [1.0, 0.0, 0.0]);
let o = ScopeOpts { width: 8, height: 101, ..Default::default() };
let s = compute(&img, ScopeKind::WaveformYc, &o);
let [y, cb, cr] = ColorStandard::Rec709.ycbcr([1.0, 0.0, 0.0]);
assert_eq!(s.peak(0).1, ((1.0 - y) * 100.0).round() as u32);
assert_eq!(s.peak(2).1, ((1.0 - (cb * cb + cr * cr).sqrt()) * 100.0).round() as u32);
}
}