AI Flow: Reverse Engineering

//! Tented vias, read from the shape of the dome and not from brightness alone. Solder mask drapes over a via's hole
//! and rim, so a lit via is a small round bump: a highlight on the side towards the light and a shadow on the other.
//! Silkscreen, copper edges and pads are bright too, but they are long, flat, or have the wrong size.
use crate::board::RES;
use crate::copper::line_filter;
use crate::raster::{self, Plane};
use image::RgbImage;
use serde_json::{json, Value};

/// Every candidate, with what was measured on it. Pixels unless the name says otherwise.
#[derive(Clone, Debug)]
pub struct Cand {
    pub x: f64, pub y: f64,
    /// strength against the local grain
    pub z: f32,
    /// highlight above the median of the ring around it, grey levels
    pub contrast: f32,
    /// the share of the response left after opening along 1 mm straight lines: 1 for a bump, 0 on a line
    pub compact: f32,
    /// half-height blob: long and short axis (mm), area (mm2), and whether it ran out of the 0.7 mm window
    pub len: f32, pub wid: f32, pub area: f32, pub runs_off: bool,
    /// brightness dipole around the highlight, as a vector and as a share of the grain
    pub vx: f32, pub vy: f32, pub shade: f32,
    /// agreement of the dipole with the face's light (cosine), filled in once the light is known
    pub cos: f32,
    /// depth of the shadow on the far side, and the glow on the near side, against the ring (grain units)
    pub shadow: f32, pub glow: f32,
    /// distance from the highlight to the darkest point on the far side (mm)
    pub shadow_at: f32,
    /// pale pixels (min RGB over 135) in the ring around it, and in a 0.35 mm patch at it
    pub pale_ring: f32, pub pale_patch: f32,
    /// red minus blue at the highlight (gold, solder and copper are warm), and the darkest tenth of the surroundings
    pub warm: f32, pub dark: f32,
    /// ring median, grey levels
    pub ring: f32,
    /// how far the middle of the blob stands towards the light past its ends (mm): a highlight wraps round a dome, so it bends towards the light; the rim light on the far side bends away
    pub bend: f32,
    /// what stands proud once lines thinner than 0.2 mm are gone and wide areas are cancelled: a dome on a copper edge (grey levels)
    pub bump: f32,
    /// the seed's response once lines are explained, and whether it was only a seed for that
    pub thz: f32, pub th_only: bool,
    /// the blob's centre, pixels
    pub bx: f64, pub by: f64,
    /// best dark circle round the highlight: centre (pixels), radius (mm) and its strength against the grain
    pub rx: f64, pub ry: f64, pub rr: f32, pub rs: f32,
    /// how far the highlight lies towards the light from that circle's centre, as a share of the radius
    pub rside: f32,
}

pub struct FaceDomes { pub cands: Vec<Cand>, pub light: (f64, f64), pub coherence: f64 }

fn median(v: &mut [f32]) -> f32 { if v.is_empty() { return 0.0; } let k = v.len() / 2; *v.select_nth_unstable_by(k, |a, b| a.partial_cmp(b).unwrap()).1 }

/// Every local maximum of a sensitive small-blob response, each measured. Nothing is rejected here.
pub fn measure(img: &RgbImage) -> FaceDomes {
    let (w, h) = (img.width() as usize, img.height() as usize);
    let l = raster::luma(img);
    let fine = raster::blur(&l, 1, 3);
    let s = raster::blur(&l, 2, 3);
    let b = raster::blur(&l, 7, 3);
    let dog = Plane { w, h, d: s.d.iter().zip(&b.d).map(|(x, y)| x - y).collect() };
    // the grain around a point: the mean size of the response, clipped so a bright pad next to a dome does not count as grain
    let mut abs: Vec<f32> = dog.d.iter().step_by(7).map(|v| v.abs()).collect();
    let clip = 3.0 * median(&mut abs);
    let noise = raster::blur(&Plane { w, h, d: dog.d.iter().map(|v| v.abs().min(clip)).collect() }, (1.5 * RES) as usize, 1);
    let z = Plane { w, h, d: dog.d.iter().zip(&noise.d).map(|(d, n)| d / (n + 1.0)).collect() };
    // what straight lines explain: the brightest opening along the four directions boards are drawn in
    let lin = (1.0 * RES) as usize; let dia = (1.0 * RES / std::f64::consts::SQRT_2) as usize;
    let mut open = Plane::new(w, h, f32::MIN);
    for (dx, dy, n) in [(1i64, 0i64, lin), (0, 1, lin), (1, 1, dia), (-1, 1, dia)] {
        let o = line_filter(&line_filter(&s, dx, dy, n, false), dx, dy, n, true);
        for i in 0..w * h { open.d[i] = open.d[i].max(o.d[i]); }
    }
    // a dome on a copper edge: thin lines are opened away by a small square, then long wide areas by the same line opening
    let th = 2 * 4 + 1;
    let erode = line_filter(&line_filter(&s, 1, 0, th, false), 0, 1, th, false);
    let thick = line_filter(&line_filter(&erode, 1, 0, th, true), 0, 1, th, true);
    let mut thick_open = Plane::new(w, h, f32::MIN);
    for (dx, dy, n) in [(1i64, 0i64, lin), (0, 1, lin), (1, 1, dia), (-1, 1, dia)] {
        let o = line_filter(&line_filter(&thick, dx, dy, n, false), dx, dy, n, true);
        for i in 0..w * h { thick_open.d[i] = thick_open.d[i].max(o.d[i]); }
    }
    let bump = Plane { w, h, d: thick.d.iter().zip(&thick_open.d).map(|(a, b)| a - b).collect() };
    let r = (0.22 * RES) as usize;
    let mx = line_filter(&line_filter(&z, 1, 0, 2 * r + 1, true), 0, 1, 2 * r + 1, true);
    // what lines do not explain, against the grain: a dome joined to a copper edge peaks here beside the line's own glint
    let thz = Plane { w, h, d: s.d.iter().zip(&open.d).zip(&noise.d).map(|((a, o), n)| (a - o) / (n + 1.0)).collect() };
    let mxt = line_filter(&line_filter(&thz, 1, 0, 2 * r + 1, true), 0, 1, 2 * r + 1, true);
    let rgb = |x: usize, y: usize| img.get_pixel(x as u32, y as u32);
    let mut done = std::collections::HashSet::new();
    let mut cands = vec![]; let mut blobs: Vec<Vec<(f32, f32)>> = vec![];
    let margin = 2usize;
    for y in margin..h - margin { for x in margin..w - margin {
        let i = y * w + x;
        let seed_z = z.d[i] > 3.0 && z.d[i] >= mx.d[i];
        let seed_t = thz.d[i] > 3.0 && thz.d[i] >= mxt.d[i];
        if !seed_z && !seed_t { continue; }
        // the highlight itself
        let (mut hx, mut hy, mut best) = (x, y, f32::MIN);
        if seed_z { for dy in -6i64..=6 { for dx in -6i64..=6 { if dx * dx + dy * dy > 36 { continue; } let (u, v) = (x as i64 + dx, y as i64 + dy); if u < 0 || v < 0 || u >= w as i64 || v >= h as i64 { continue; } let q = fine.at(u as usize, v as usize); if q > best { best = q; hx = u as usize; hy = v as usize; } } } } else { best = fine.at(x, y); }
        if !done.insert((hx, hy)) { continue; }
        // the ring around it
        let (mut ring, mut pale, mut ringn, mut near) = (vec![], 0u32, 0u32, vec![]);
        let (mut patch, mut patchn) = (0u32, 0u32);
        for dy in -32i64..=32 { for dx in -32i64..=32 {
            let d2 = dx * dx + dy * dy; if d2 > 32 * 32 { continue; }
            let (u, v) = (hx as i64 + dx, hy as i64 + dy); if u < 0 || v < 0 || u >= w as i64 || v >= h as i64 { continue; }
            let (u, v) = (u as usize, v as usize); let p = rgb(u, v); let m = p[0].min(p[1]).min(p[2]);
            near.push(fine.at(u, v));
            if d2 >= 18 * 18 { ring.push(fine.at(u, v)); ringn += 1; if m > 135 { pale += 1; } }
            if d2 <= 7 * 7 { patchn += 1; if m > 135 { patch += 1; } }
        } }
        let rm = median(&mut ring); let c = best - rm;
        near.sort_by(|a, b| a.partial_cmp(b).unwrap()); let dark = near[near.len() / 10];
        // the half-height blob
        let t = rm + 0.5 * c;
        let win = 28i64; let mut seen = vec![false; ((2 * win + 1) * (2 * win + 1)) as usize];
        let idx = |dx: i64, dy: i64| ((dy + win) * (2 * win + 1) + dx + win) as usize;
        let mut st = vec![(0i64, 0i64)]; seen[idx(0, 0)] = true; let mut blob: Vec<(f32, f32)> = vec![];
        let (mut n, mut sx, mut sy, mut sxx, mut syy, mut sxy, mut off) = (0f32, 0f32, 0f32, 0f32, 0f32, 0f32, false);
        while let Some((dx, dy)) = st.pop() {
            n += 1.0; let (fx, fy) = (dx as f32, dy as f32); blob.push((fx, fy)); sx += fx; sy += fy; sxx += fx * fx; syy += fy * fy; sxy += fx * fy;
            for (ax, ay) in [(1i64, 0i64), (-1, 0), (0, 1), (0, -1)] {
                let (nx, ny) = (dx + ax, dy + ay);
                if nx.abs() > win || ny.abs() > win { off = true; continue; }
                if seen[idx(nx, ny)] { continue; }
                let (u, v) = (hx as i64 + nx, hy as i64 + ny); if u < 0 || v < 0 || u >= w as i64 || v >= h as i64 { off = true; continue; }
                if fine.at(u as usize, v as usize) >= t { seen[idx(nx, ny)] = true; st.push((nx, ny)); }
            }
        }
        let (mx0, my0) = (sx / n, sy / n); let (cxx, cyy, cxy) = (sxx / n - mx0 * mx0, syy / n - my0 * my0, sxy / n - mx0 * my0);
        let tr = cxx + cyy; let det = (cxx * cyy - cxy * cxy).max(0.0); let disc = (tr * tr / 4.0 - det).max(0.0).sqrt();
        let (l1, l2) = ((tr / 2.0 + disc).max(0.0), (tr / 2.0 - disc).max(0.0));
        // the dipole
        let (mut vx, mut vy, mut vn) = (0f32, 0f32, 0f32);
        for dy in -18i64..=18 { for dx in -18i64..=18 { let d2 = dx * dx + dy * dy; if !(36..=324).contains(&d2) { continue; } let (u, v) = (hx as i64 + dx, hy as i64 + dy); if u < 0 || v < 0 || u >= w as i64 || v >= h as i64 { continue; }
            let d = (d2 as f32).sqrt(); let q = fine.at(u as usize, v as usize) - rm; vx += q * dx as f32 / d; vy += q * dy as f32 / d; vn += 1.0; } }
        vx /= vn.max(1.0); vy /= vn.max(1.0);
        blobs.push(blob);
        let nz = noise.at(hx, hy) + 1.0;
        let hp = rgb(hx, hy);
        cands.push(Cand { x: hx as f64 + 0.5, y: hy as f64 + 0.5, z: { let mut q = f32::MIN; for dy in -6i64..=6 { for dx in -6i64..=6 { let (u, v) = (hx as i64 + dx, hy as i64 + dy); if u >= 0 && v >= 0 && u < w as i64 && v < h as i64 { q = q.max(z.at(u as usize, v as usize)); } } } q }, contrast: c, compact: ((s.at(hx, hy) - open.at(hx, hy)) / (s.at(hx, hy) - b.at(hx, hy)).max(1.0)).clamp(-1.0, 2.0),
            len: 4.0 * l1.sqrt() / RES as f32, wid: 4.0 * l2.sqrt() / RES as f32, area: n / (RES * RES) as f32, runs_off: off,
            vx, vy, shade: (vx * vx + vy * vy).sqrt() / nz, cos: 0.0, shadow: 0.0, glow: 0.0, shadow_at: 0.0,
            pale_ring: pale as f32 / ringn.max(1) as f32, pale_patch: patch as f32 / patchn.max(1) as f32, warm: hp[0] as f32 - hp[2] as f32, dark, ring: rm, bend: 0.0, thz: thz.at(hx, hy), th_only: !seed_z, bx: hx as f64 + 0.5 + mx0 as f64, by: hy as f64 + 0.5 + my0 as f64, rx: 0.0, ry: 0.0, rr: 0.0, rs: 0.0, rside: 0.0, bump: { let mut b0 = f32::MIN; for dy in -8i64..=8 { for dx in -8i64..=8 { let (u, v) = (hx as i64 + dx, hy as i64 + dy); if u >= 0 && v >= 0 && u < w as i64 && v < h as i64 { b0 = b0.max(bump.at(u as usize, v as usize)); } } } b0 } });
    } }
    // the light: where the dipoles of the plain small bumps agree
    let (mut lx, mut ly, mut m) = (0f64, 0f64, 0f64);
    for c in cands.iter().filter(|c| c.z > 5.0 && c.compact > 0.5 && !c.runs_off && c.len < 0.6) { let n = ((c.vx * c.vx + c.vy * c.vy) as f64).sqrt().max(1e-6); lx += c.vx as f64 / n; ly += c.vy as f64 / n; m += 1.0; }
    let n = (lx * lx + ly * ly).sqrt().max(1e-6);
    let coherence = if m > 0.0 { n / m } else { 0.0 };
    let (lx, ly) = (lx / n, ly / n);
    // shadow on the far side, glow on the near side
    let disk = |cx: f64, cy: f64, rd: i64| { let (mut a, mut n) = (0f32, 0f32); for dy in -rd..=rd { for dx in -rd..=rd { if dx * dx + dy * dy > rd * rd { continue; } let (u, v) = (cx as i64 + dx, cy as i64 + dy); if u < 0 || v < 0 || u >= w as i64 || v >= h as i64 { continue; } a += fine.at(u as usize, v as usize); n += 1.0; } } a / n.max(1.0) };
    for (c, blob) in cands.iter_mut().zip(&blobs) {
        c.bend = bend_towards(blob, (lx as f32, ly as f32));

        let mag = (c.vx * c.vx + c.vy * c.vy).sqrt().max(1e-6);
        c.cos = (c.vx as f64 * lx + c.vy as f64 * ly) as f32 / mag;
        let nz = noise.at(c.x as usize, c.y as usize) + 1.0;
        c.shadow = (c.ring - disk(c.x - lx * 10.0, c.y - ly * 10.0, 5)) / nz;
        c.glow = (disk(c.x + lx * 10.0, c.y + ly * 10.0, 5) - c.ring) / nz;
        let (mut bd, mut at) = (f32::MAX, 0f32);
        for k in 2..=18 { let (u, v) = (c.x - lx * k as f64, c.y - ly * k as f64); if u < 0.0 || v < 0.0 || u >= w as f64 || v >= h as f64 { break; } let q = fine.at(u as usize, v as usize); if q < bd { bd = q; at = k as f32; } }
        c.shadow_at = at / RES as f32;
        if is_dome(c) || may_be_on_track(c) {
            let (rx, ry, rr, rs) = ring_fit(&fine, c.x, c.y, nz);
            c.rx = rx; c.ry = ry; c.rr = (rr / RES) as f32; c.rs = rs;
            c.rside = (((c.x - rx) * lx + (c.y - ry) * ly) / rr.max(1.0)) as f32;
        }
    }
    FaceDomes { cands, light: (lx, ly), coherence }
}

/// The dark circle a highlight sits in or on: a clearance gap round a via's pad, or the foot of its dome. Each centre and
/// radius is scored by how much darker the circle is than just inside and outside it, round most of the way (a track
/// breaks the circle where it leaves, so the extremes are trimmed).
fn ring_fit(f: &Plane, hx: f64, hy: f64, nz: f32) -> (f64, f64, f64, f32) {
    let n = 32usize;
    let dirs: Vec<(f64, f64)> = (0..n).map(|k| { let a = k as f64 / n as f64 * std::f64::consts::TAU; (a.cos(), a.sin()) }).collect();
    let at = |x: f64, y: f64| -> f32 { let (u, v) = (x as i64, y as i64); if u < 0 || v < 0 || u >= f.w as i64 || v >= f.h as i64 { 0.0 } else { f.at(u as usize, v as usize) } };
    let score = |cx: f64, cy: f64, r: f64| -> f32 {
        let mut d: Vec<f32> = dirs.iter().map(|(c, s)| (at(cx + c * (r - 3.0), cy + s * (r - 3.0)) + at(cx + c * (r + 3.0), cy + s * (r + 3.0))) / 2.0 - at(cx + c * r, cy + s * r)).collect();
        d.sort_by(|a, b| a.partial_cmp(b).unwrap());
        let t = n / 5; d[t..n - t].iter().sum::<f32>() / (n - 2 * t) as f32
    };
    let (mut best, mut bx, mut by, mut br) = (f32::MIN, hx, hy, 10.0);
    let mut dy = -20i64; while dy <= 20 { let mut dx = -20i64; while dx <= 20 { let mut r = 12i64; while r <= 20 {
        // the highlight lies inside the circle or on it
        if ((dx * dx + dy * dy) as f64).sqrt() <= r as f64 + 2.0 { let q = score(hx + dx as f64, hy + dy as f64, r as f64); if q > best { best = q; bx = hx + dx as f64; by = hy + dy as f64; br = r as f64; } }
        r += 2; } dx += 2; } dy += 2; }
    let (x0, y0, r0) = (bx, by, br);
    for dy in -2i64..=2 { for dx in -2i64..=2 { for dr in -1i64..=1 { let (x, y, r) = (x0 + dx as f64, y0 + dy as f64, (r0 + dr as f64).max(12.0)); let q = score(x, y, r); if q > best { best = q; bx = x; by = y; br = r; } } } }
    (bx, by, br, best / nz)
}

/// Fit the blob's height along the light against the square of its position across it. The result is the sagitta in mm,
/// positive when the middle stands towards the light.
fn bend_towards(blob: &[(f32, f32)], l: (f32, f32)) -> f32 {
    if blob.len() < 6 { return 0.0; }
    let pts: Vec<(f32, f32)> = blob.iter().map(|p| (p.0 * -l.1 + p.1 * l.0, p.0 * l.0 + p.1 * l.1)).collect();
    let n = pts.len() as f32;
    let (mt, ma) = (pts.iter().map(|p| p.0).sum::<f32>() / n, pts.iter().map(|p| p.1).sum::<f32>() / n);
    let mq = pts.iter().map(|p| (p.0 - mt).powi(2)).sum::<f32>() / n;
    let (mut sxy, mut sxx, mut half) = (0f32, 0f32, 0f32);
    for p in &pts { let q = (p.0 - mt).powi(2) - mq; sxy += q * (p.1 - ma); sxx += q * q; half = half.max((p.0 - mt).abs()); }
    if sxx <= 0.0 { return 0.0; }
    -(sxy / sxx) * half * half / RES as f32
}

/// Thresholds measured on the Driver photos (tuning sample of 272 labelled candidates, 2026-11-05).
/// A dome stands out from the grain, is a small bump and not part of a straight line, is not on silkscreen, and has no
/// hole in its shadow. Long thin highlights must bend towards the light.
pub fn is_dome(c: &Cand) -> bool {
    // the smallest domes are as small as grain: what tells them apart is a clear shadow on the far side from the light
    let lit = c.z >= 6.0 && c.shadow >= 2.0 && c.cos >= 0.9;
    let free = c.z >= 4.0 && c.contrast >= 55.0 && c.compact >= 0.75 && (c.len >= 0.16 || lit) && c.len >= 0.12 && c.len < 0.62 && !c.runs_off && c.area < 0.16 && c.pale_ring < 0.35 && c.shadow >= -9.0;
    // a thin straight dash (silkscreen, a piece of track) is not a highlight: those bend
    let dash = c.len >= ARC_MIN_LEN && c.wid <= 0.16 && c.bend < 0.02;
    let free = free && !dash;
    // a via with a track leaving it: the track's glint joins its highlight, but the dark circle round its pad still shows
    let on_track = may_be_on_track(c) && c.rs >= RING_ON_TRACK;
    free || on_edge(c) || on_track
}

/// A dome on a copper edge: its body stands proud of the line, wider than the line's own glint and than the rim of a ring
/// (on crops picked from a whole-board overlay of the Driver bottom, rail domes give a bump of 15 and up and a width of 0.34 and up, rims and tracks 14 and 0.27 at most).
fn on_edge(c: &Cand) -> bool { c.runs_off && c.thz >= 4.0 && c.compact >= 0.7 && c.shadow >= 1.5 && c.cos >= 0.9 && c.pale_ring < 0.3 && c.bump >= 14.5 && c.wid >= 0.3 }

fn may_be_on_track(c: &Cand) -> bool { c.z >= 4.0 && c.contrast >= 55.0 && c.compact >= 0.7 && c.pale_ring < 0.2 && c.len < 1.6 && c.wid < 0.6 && c.shadow >= -9.0 }

/// A via as the photo shows it: where its centre is and how wide it looks.
#[derive(Clone, Debug)]
pub struct Dome { pub x: f64, pub y: f64, pub d_mm: f64, pub shape: &'static str, pub z: f32, /// radius of the dark circle round it, pixels, where one was found
    pub ring_px: f64 }

/// Below this share of agreement the bumps on a face are not lit as one population, so they are not domes. On the Driver
/// photos the bottom face gives 0.6 and the flat top gives under 0.2 (the mask lies flat over the vias there).
pub const FACE_COHERENCE: f64 = 0.4;
/// A dark circle 0.3 to 0.5 mm in radius, at least this strong (grain units), is the clearance round a via's pad. On the tuning
/// sample labelled ring vias give 2.4 and up (tenth percentile), plain domes 2.0 and down (ninetieth).
const RING_STRONG: f32 = 2.2;
/// The same where a track's glint joins the highlight, so the circle is all the evidence there is.
const RING_ON_TRACK: f32 = 2.6;
/// The dark circle's centre line lies this far outside the pad's edge: half a clearance gap as the photo draws it (mm).
const RING_TO_PAD_MM: f64 = 0.05;
/// A highlight this long and this thin is the near arc of a via too wide to light in full.
const ARC_MIN_LEN: f32 = 0.38; const ARC_MAX_WID: f32 = 0.2;
/// Two domes closer than this are one (mm): the closest fields on the Driver photos have a pitch near 0.5.
const SAME_DOME_MM: f64 = 0.3;
/// A blob's half-height centre lies this far towards the light from the middle of the dome: measured by eye on 18 crops of the tuning sample.
const BLOB_TO_LIGHT_MM: f64 = 0.03;

/// Whether a face shows domes at all: the shading directions of its accepted bumps agree.
pub fn face_coherence(f: &FaceDomes) -> f64 {
    let (mut a, mut b, mut n) = (0f64, 0f64, 0f64);
    for c in f.cands.iter().filter(|c| is_dome(c) && !c.runs_off) { let m = ((c.vx * c.vx + c.vy * c.vy) as f64).sqrt().max(1e-6); a += c.vx as f64 / m; b += c.vy as f64 / m; n += 1.0; }
    if n > 0.0 { (a * a + b * b).sqrt() / n } else { 0.0 }
}

/// From accepted candidates to vias. A strong circle gives the centre and width outright, and takes in every highlight inside
/// it. A long thin arc lies on the near side of a via: its centre is measured from the strong circles on this face. A round
/// blob is the dome itself.
pub fn finish(f: &FaceDomes) -> Vec<Dome> {
    let (lx, ly) = f.light;
    let good: Vec<&Cand> = f.cands.iter().filter(|c| is_dome(c)).collect();
    let is_ring = |c: &Cand| c.rs >= if c.runs_off || c.len >= 0.62 { RING_ON_TRACK } else { RING_STRONG };
    let is_arc = |c: &Cand| c.len >= ARC_MIN_LEN && c.wid <= ARC_MAX_WID && !c.runs_off;
    // what the strong rings say about the arcs on this face
    let mut along: Vec<f32> = vec![]; let mut dia: Vec<f32> = vec![];
    for c in good.iter().filter(|c| is_ring(c) && is_arc(c) && c.len < 0.62) { along.push(((c.bx - c.rx) * lx + (c.by - c.ry) * ly) as f32); dia.push(2.0 * (c.rr - RING_TO_PAD_MM as f32)); }
    let (arc_off, arc_d) = if along.len() >= 5 { (median(&mut along) as f64, median(&mut dia) as f64) } else { (0.0, 0.0) };
    let mut blob_d: Vec<f32> = good.iter().filter(|c| !c.runs_off && !is_ring(c) && !is_arc(c)).map(|c| (c.len + c.wid) / 2.0).collect();
    let edge_d = median(&mut blob_d) as f64;
    let mut rings: Vec<&Cand> = good.iter().copied().filter(|c| is_ring(c) && c.rside >= -0.25).collect();
    rings.sort_by(|a, b| b.rs.partial_cmp(&a.rs).unwrap());
    let mut out: Vec<Dome> = vec![];
    for c in rings {
        if out.iter().any(|o| crate::geom::dist((o.x, o.y), (c.rx, c.ry)) < SAME_DOME_MM * RES) { continue; }
        out.push(Dome { x: c.rx + 0.5, y: c.ry + 0.5, d_mm: 2.0 * (c.rr as f64 - RING_TO_PAD_MM), shape: "ring", z: c.z, ring_px: c.rr as f64 * RES });
    }
    let n_rings = out.len();
    let mut rest: Vec<&Cand> = good.iter().copied().filter(|c| !is_ring(c) || c.rside < -0.25).collect();
    rest.sort_by(|a, b| b.z.partial_cmp(&a.z).unwrap());
    for c in rest {
        // inside a via already found: its rim light or the highlight on its hole
        if out[..n_rings].iter().any(|o| crate::geom::dist((o.x, o.y), (c.x, c.y)) < o.ring_px + 2.0) { continue; }
        // a highlight that runs into a line stands for a via only as a dome on a copper edge
        if (c.runs_off || c.len >= 0.62) && !on_edge(c) { continue; }
        // the rim on the far side of a circle is not a via of its own
        if c.rs >= RING_STRONG && c.rside < -0.25 && is_arc(c) { continue; }
        let (x, y, d, shape) = if c.runs_off { (c.x, c.y, (c.wid as f64).clamp(0.2, edge_d.max(0.45)), "edge") }
            else if is_arc(c) && arc_d > 0.0 { (c.bx - lx * arc_off, c.by - ly * arc_off, arc_d, "arc") }
            else { (c.bx - lx * BLOB_TO_LIGHT_MM * RES, c.by - ly * BLOB_TO_LIGHT_MM * RES, ((c.len + c.wid) / 2.0) as f64, "blob") };
        if out.iter().any(|o| crate::geom::dist((o.x, o.y), (x, y)) < SAME_DOME_MM * RES) { continue; }
        out.push(Dome { x, y, d_mm: d, shape, z: c.z, ring_px: 0.0 });
    }
    out
}

pub fn cand_json(c: &Cand, side: &str) -> Value {
    let (mx, my) = crate::board::to_mm(c.x, c.y);
    let r = |v: f32| (v as f64 * 1000.0).round() / 1000.0;
    let rnd = |x: f64, y: f64| -> Value { let m = crate::board::to_mm(x, y); json!([(m.0 * 1000.0).round() / 1000.0, (m.1 * 1000.0).round() / 1000.0]) };
    json!({"side": side, "pos_mm": [(mx * 1000.0).round() / 1000.0, (my * 1000.0).round() / 1000.0], "z": r(c.z), "contrast": r(c.contrast), "compact": r(c.compact), "len": r(c.len), "wid": r(c.wid), "area": r(c.area), "runsOff": c.runs_off,
        "shade": r(c.shade), "cos": r(c.cos), "shadow": r(c.shadow), "glow": r(c.glow), "shadowAt": r(c.shadow_at), "paleRing": r(c.pale_ring), "palePatch": r(c.pale_patch), "warm": r(c.warm), "dark": r(c.dark), "ring": r(c.ring), "bump": r(c.bump), "bend": r(c.bend), "dome": is_dome(c), "blob_mm": rnd(c.bx, c.by), "ring_mm": rnd(c.rx, c.ry), "thz": r(c.thz), "thOnly": c.th_only, "rr": r(c.rr), "rs": r(c.rs), "rside": r(c.rside)})
}

/// Does a face show domes at all? Average its picture over many known via positions and compare the middle with the
/// surroundings. Grain cancels, so anything a via does to the picture stands out. The Driver bottom scores +56 standard
/// errors and the top -3.4, inside the limit and not a bump. This tests for domes, not for clearance rings: most vias sit in
/// pours and have no ring. It is a conservative gate on calling a via inward-only. Returns the score and the count used.
pub fn face_signature(img: &RgbImage, at_mm: &[(f64, f64)]) -> (f64, usize) {
    let l = raster::luma(img); let (w, h) = (l.w as i64, l.h as i64); let k = crate::prep::mask_tone(img).exposure();
    let mut diffs: Vec<f64> = vec![];
    for (x, y) in at_mm {
        let (u, v) = crate::board::to_px(*x, *y); let (u, v) = (u as i64, v as i64);
        if u < 24 || v < 24 || u >= w - 24 || v >= h - 24 { continue; }
        let (mut a, mut an, mut b, mut bn, mut all) = (0f64, 0f64, 0f64, 0f64, vec![]);
        for dy in -24i64..=24 { for dx in -24i64..=24 { let d = dx * dx + dy * dy; if d > 576 { continue; } let q = l.at((u + dx) as usize, (v + dy) as usize); all.push(q); if d <= 64 { a += q as f64; an += 1.0; } else if d >= 144 { b += q as f64; bn += 1.0; } } }
        // silkscreen or a light part over the spot says nothing about the via under it
        if median(&mut all) < VIEW_MAX * k { diffs.push(a / an - b / bn); }
    }
    if diffs.len() < 20 { return (0.0, diffs.len()); }
    let n = diffs.len() as f64; let m = diffs.iter().sum::<f64>() / n;
    let sd = (diffs.iter().map(|d| (d - m).powi(2)).sum::<f64>() / (n - 1.0)).sqrt().max(1e-6);
    (m / (sd / n.sqrt()), diffs.len())
}
/// A 1.2 mm patch whose median is brighter than this is not bare mask (Driver photos, 8 bit luma).
const VIEW_MAX: f32 = 110.0;
/// Standard errors above which a face shows domes. A judgement, set well above the Driver top and well below its bottom.
pub const SIGNATURE_MIN: f64 = 5.0;

/// Whether a drilled pad's hole is open on this photo: its middle is dark, where a soldered pin's tip is bright.
pub fn hole_is_open(img: &RgbImage, x: f64, y: f64, drill: f64) -> bool {
    let (u, v) = crate::board::to_px(x, y); let r = (drill / 2.0 * 0.6 * RES) as i64; let mut q = vec![];
    for dy in -r..=r { for dx in -r..=r { let (a, b) = (u as i64 + dx, v as i64 + dy); if a >= 0 && b >= 0 && a < img.width() as i64 && b < img.height() as i64 && dx * dx + dy * dy <= r * r { let p = img.get_pixel(a as u32, b as u32); q.push(0.299 * p[0] as f32 + 0.587 * p[1] as f32 + 0.114 * p[2] as f32); } } }
    median(&mut q) < OPEN_HOLE_MAX
}
/// The middle of an open hole on the Driver photos is deep brown, under this 8 bit luma; a pin tip is well over it.
const OPEN_HOLE_MAX: f32 = 70.0;

/// A via as a pit: what a tented via looks like under flat, frontal light, where nothing is lit from one side and no dome highlight forms. The mask sinks into the hole, so the via
/// reads as a small dark round spot, darker than the mask at any level of the copper, with the camera's sharpening often drawing a light ring round it (CM5 photos, both faces). Everything measured
/// on every candidate is kept, so a rejection can be checked.
#[derive(Clone, Debug)]
pub struct Pit { pub x: f64, pub y: f64, /// the half-depth blob's long and short axis, mm
    pub len: f32, pub wid: f32, /// the ring's median brightness (green) minus the darkest point, and that as a share of the ring
    pub depth: f32, pub share: f32, /// against the local grain
    pub z: f32, /// share of the ring that is plain mask, and the darkest point's G - R over the mask's own
    pub ring_mask: f32, pub green: f32, /// share of the ring that is pale: print or a part's end
    pub ring_pale: f32, pub runs_off: bool, pub keep: bool }

/// What a pit is, against nothing but this photo. The blob must close inside the window and be round (no gap line or track end), its surroundings must be mask (no part or print), and it must stay green at
/// its darkest (a gap between parts is neutral). The bars on roundness, ring, depth, strength and colour are read from the CM5 everypin photos, tile set A of the 5 mm checkerboard only, by how many pits on one face have a pit within
/// 0.1 mm on the other against the same count with one face mirrored top to bottom, which shows chance: 114 of 167 top pits against 1 at these values. The counter of a printed letter and the shadow beside a part's bright end are dark too, and the
/// pale share of the ring is what tells them, at some cost: it gives up 8 of 122 confirmed pits to lose 18 of 185. The size comes from the face's own population.
const PIT_ROUND: f32 = 2.0; const PIT_RING_MASK: f32 = 0.3; const PIT_MIN_SHARE: f32 = 0.12; const PIT_MIN_Z: f32 = 2.5; const PIT_GREEN: f32 = 0.05; const PIT_MAX_PALE: f32 = 0.4;
/// A population is a via field only if it is many of one size: this many at least, with the middle half of their widths within this share of the median.
const PIT_MIN_COUNT: usize = 40; const PIT_SPREAD: f32 = 0.5;
/// The size band round the population's median, as shares of it.
const PIT_BAND: (f32, f32) = (0.6, 1.6);

pub struct PitFace { pub pits: Vec<Pit>, pub median_mm: f64, pub spread: f64, pub shows: bool }

/// Read the pits on a face. `shows`: they form a population of one size, so this face shows its vias this way.
pub fn read_pits(img: &RgbImage) -> PitFace {
    let (w, h) = (img.width() as usize, img.height() as usize);
    let tone = crate::prep::mask_tone(img);
    let g = Plane { w, h, d: img.pixels().map(|p| p[1] as f32).collect() };
    let fine = raster::blur(&g, 1, 1);
    let (s, b) = (raster::blur(&g, 2, 2), raster::blur(&g, 8, 2));
    let dog = Plane { w, h, d: b.d.iter().zip(&s.d).map(|(x, y)| x - y).collect() };
    let mut abs: Vec<f32> = dog.d.iter().step_by(7).map(|v| v.abs()).collect();
    let clip = 3.0 * median(&mut abs);
    let noise = raster::blur(&Plane { w, h, d: dog.d.iter().map(|v| v.abs().min(clip)).collect() }, (1.5 * RES) as usize, 1);
    let r = 6usize;
    let mx = line_filter(&line_filter(&dog, 1, 0, 2 * r + 1, true), 0, 1, 2 * r + 1, true);
    let plain = |x: usize, y: usize| { let p = img.get_pixel(x as u32, y as u32); (p[1] as f32 - p[0] as f32) > 0.3 * tone.gr50 && p[1] > p[2] };
    let mut pits = vec![]; let win = 18i64;
    for y in 20..h.saturating_sub(20) { for x in 20..w.saturating_sub(20) {
        let i = y * w + x;
        if dog.d[i] < mx.d[i] || dog.d[i] <= 0.0 || dog.d[i] / (noise.d[i] + 1.0) < 2.5 { continue; }
        // the darkest point near the peak, and the ring round it
        let (mut cx, mut cy, mut low) = (x, y, f32::MAX);
        for dy in -4i64..=4 { for dx in -4i64..=4 { let (u, v) = ((x as i64 + dx) as usize, (y as i64 + dy) as usize); if fine.at(u, v) < low { low = fine.at(u, v); cx = u; cy = v; } } }
        let (mut ring, mut on, mut tot, mut pale) = (vec![], 0u32, 0u32, 0u32);
        for dy in -17i64..=17 { for dx in -17i64..=17 { let d2 = dx * dx + dy * dy; if !(144..=289).contains(&d2) { continue; } let (u, v) = (cx as i64 + dx, cy as i64 + dy); if u < 0 || v < 0 || u >= w as i64 || v >= h as i64 { continue; } ring.push(fine.at(u as usize, v as usize)); tot += 1; if plain(u as usize, v as usize) { on += 1; } let q = img.get_pixel(u as u32, v as u32); if q[0].min(q[1]).min(q[2]) > 135 { pale += 1; } } }
        let rm = median(&mut ring); let depth = rm - low;
        if depth <= 0.0 { continue; }
        let t = low + 0.5 * depth;
        let mut seen = vec![false; ((2 * win + 1) * (2 * win + 1)) as usize]; let idx = |dx: i64, dy: i64| ((dy + win) * (2 * win + 1) + dx + win) as usize;
        let mut st = vec![(0i64, 0i64)]; seen[idx(0, 0)] = true;
        let (mut n, mut sx, mut sy, mut sxx, mut syy, mut sxy, mut off) = (0f32, 0f32, 0f32, 0f32, 0f32, 0f32, false);
        while let Some((dx, dy)) = st.pop() {
            n += 1.0; let (fx, fy) = (dx as f32, dy as f32); sx += fx; sy += fy; sxx += fx * fx; syy += fy * fy; sxy += fx * fy;
            for (ax, ay) in [(1i64, 0i64), (-1, 0), (0, 1), (0, -1)] { let (nx, ny) = (dx + ax, dy + ay);
                if nx.abs() > win || ny.abs() > win { off = true; continue; }
                if seen[idx(nx, ny)] { continue; }
                let (u, v) = (cx as i64 + nx, cy as i64 + ny); if u < 0 || v < 0 || u >= w as i64 || v >= h as i64 { off = true; continue; }
                if fine.at(u as usize, v as usize) <= t { seen[idx(nx, ny)] = true; st.push((nx, ny)); } }
        }
        let (m0, m1) = (sx / n, sy / n); let (cxx, cyy, cxy) = (sxx / n - m0 * m0, syy / n - m1 * m1, sxy / n - m0 * m1);
        let tr = cxx + cyy; let det = (cxx * cyy - cxy * cxy).max(0.0); let disc = (tr * tr / 4.0 - det).max(0.0).sqrt();
        let (l1, l2) = ((tr / 2.0 + disc).max(0.0), (tr / 2.0 - disc).max(0.0));
        let centre = img.get_pixel(cx as u32, cy as u32);
        pits.push(Pit { x: cx as f64 + 0.5 + m0 as f64, y: cy as f64 + 0.5 + m1 as f64, len: 4.0 * l1.sqrt() / RES as f32, wid: 4.0 * l2.sqrt() / RES as f32, depth, share: depth / rm.max(1.0), z: dog.d[i] / (noise.d[i] + 1.0),
            ring_mask: on as f32 / tot.max(1) as f32, green: (centre[1] as f32 - centre[0] as f32) / tone.gr50.max(1.0), ring_pale: pale as f32 / tot.max(1) as f32, runs_off: off, keep: false });
    } }
    // close peaks of one spot: the deepest stands
    pits.sort_by(|a, b| b.depth.partial_cmp(&a.depth).unwrap());
    let mut uniq: Vec<Pit> = vec![];
    for p in pits { if !uniq.iter().any(|q| crate::geom::dist((q.x, q.y), (p.x, p.y)) < 6.0) { uniq.push(p); } }
    let passes = |p: &Pit| !p.runs_off && p.z >= PIT_MIN_Z && p.share >= PIT_MIN_SHARE && p.ring_mask >= PIT_RING_MASK && p.ring_pale <= PIT_MAX_PALE && p.green >= PIT_GREEN && p.wid > 0.0 && p.len / p.wid <= PIT_ROUND;
    let mut d: Vec<f32> = uniq.iter().filter(|p| passes(p)).map(|p| (p.len + p.wid) / 2.0).collect();
    d.sort_by(|a, b| a.partial_cmp(b).unwrap());
    let (med, q1, q3) = if d.is_empty() { (0.0, 0.0, 0.0) } else { (d[d.len() / 2], d[d.len() / 4], d[d.len() * 3 / 4]) };
    for p in uniq.iter_mut() { let dm = (p.len + p.wid) / 2.0; p.keep = !p.runs_off && p.z >= PIT_MIN_Z && p.share >= PIT_MIN_SHARE && p.ring_mask >= PIT_RING_MASK && p.ring_pale <= PIT_MAX_PALE && p.green >= PIT_GREEN && p.wid > 0.0 && p.len / p.wid <= PIT_ROUND && dm >= PIT_BAND.0 * med && dm <= PIT_BAND.1 * med; }
    let spread = if med > 0.0 { (q3 - q1) / med } else { 9.0 };
    // the same via seen in a pour's clearance fails the plain-mask ring above, and most vias under a BGA sit in one
    for r in round_pits(img) {
        match uniq.iter_mut().find(|q| crate::geom::dist((q.x, q.y), (r.x, r.y)) < 4.0) { Some(q) => q.keep = true, None => uniq.push(r) }
    }
    let n = uniq.iter().filter(|p| p.keep).count();
    PitFace { pits: uniq, median_mm: med as f64, spread: spread as f64, shows: n >= PIT_MIN_COUNT && spread <= PIT_SPREAD }
}

/// Bars for `round_pits`, in 8-bit green at 40 px/mm. Set on tile set A of the CM5 close-up set (2026-11-07) from where pairs across the faces stop rising while chance falls: 103 pairs against about 1 by chance.
/// A via's rim stands this far above its centre.
const ROUND_DEPTH: f32 = 22.0;
/// Spread round the circle as a share of that depth.
const ROUND_SPREAD: f32 = 0.9;
/// The centre's radius at half depth, px: a via is not a speck and not a part's shadow.
const ROUND_RADIUS: (usize, usize) = (2, 7);
/// Share of what surrounds it that is pale (print, a part's end) or black (a part).
const ROUND_PALE: f32 = 0.2; const ROUND_BLACK: f32 = 0.2;
/// Flooded out from the centre at half depth a via's core stays this small and this round; a gap between tracks runs away.
const ROUND_CORE_R: f32 = 8.0; const ROUND_CORE_ASPECT: f32 = 1.6;

/// Vias by what a via looks like whatever lies round it: a dark centre with a round rim. Unlike the test above this needs no plain mask nearby.
pub fn round_pits(img: &RgbImage) -> Vec<Pit> {
    let (w, h) = (img.width() as usize, img.height() as usize);
    let g = Plane { w, h, d: img.pixels().map(|p| p[1] as f32).collect() };
    let (coarse, fine) = (raster::blur(&g, 1, 2), raster::blur(&g, 1, 1));
    let neg = Plane { w, h, d: coarse.d.iter().map(|v| -*v).collect() };
    let lowest = line_filter(&line_filter(&neg, 1, 0, 9, true), 0, 1, 9, true);
    let rays: Vec<(f32, f32)> = (0..24).map(|k| { let a = k as f32 / 24.0 * std::f32::consts::TAU; (a.cos(), a.sin()) }).collect();
    let mut out = vec![];
    for y in 16..h.saturating_sub(16) { for x in 16..w.saturating_sub(16) {
        if neg.at(x, y) < lowest.at(x, y) { continue; }
        let (mut m, mut sd) = ([0f32; 15], [0f32; 15]);
        let (mut pale, mut black, mut seen) = (0u32, 0u32, 0u32);
        for r in 0..15usize {
            let (mut sum, mut sq) = (0f32, 0f32);
            for (c, si) in &rays {
                let (u, v) = ((x as f32 + r as f32 * c) as u32, (y as f32 + r as f32 * si) as u32);
                let p = img.get_pixel(u, v); let gv = p[1] as f32;
                sum += gv; sq += gv * gv;
                if (3..12).contains(&r) { seen += 1; if p[0] as f32 > 0.98 * gv || p[2] as f32 > 0.98 * gv { pale += 1; } if gv < 25.0 { black += 1; } }
            }
            let n = rays.len() as f32; m[r] = sum / n; sd[r] = (sq / n - m[r] * m[r]).max(0.0).sqrt();
        }
        let centre = (m[0] + m[1]) / 2.0;
        let depth = m[3..7].iter().cloned().fold(f32::MIN, f32::max) - centre;
        if depth <= ROUND_DEPTH { continue; }
        if sd[2..9].iter().sum::<f32>() / 7.0 / depth.max(1.0) >= ROUND_SPREAD { continue; }
        let half = centre + depth / 2.0;
        let rad = m.iter().position(|v| *v >= half).unwrap_or(0);
        if rad < ROUND_RADIUS.0 || rad > ROUND_RADIUS.1 { continue; }
        if pale as f32 / seen as f32 >= ROUND_PALE || black as f32 / seen as f32 >= ROUND_BLACK { continue; }
        let c = img.get_pixel(x as u32, y as u32);
        if c[1] as f32 - c[0] as f32 <= -5.0 { continue; }
        // flood the core at half depth
        if fine.at(x, y) > half { continue; }
        let win = 14i64; let idx = |dx: i64, dy: i64| ((dy + win) * (2 * win + 1) + dx + win) as usize;
        let mut mark = vec![false; ((2 * win + 1) * (2 * win + 1)) as usize]; mark[idx(0, 0)] = true;
        let mut st = vec![(0i64, 0i64)];
        let (mut rmax, mut x0, mut x1, mut y0, mut y1) = (0f32, 0i64, 0i64, 0i64, 0i64);
        while let Some((dx, dy)) = st.pop() {
            rmax = rmax.max(((dx * dx + dy * dy) as f32).sqrt()); x0 = x0.min(dx); x1 = x1.max(dx); y0 = y0.min(dy); y1 = y1.max(dy);
            for (ax, ay) in [(1i64, 0i64), (-1, 0), (0, 1), (0, -1)] { let (nx, ny) = (dx + ax, dy + ay);
                if nx.abs() > win || ny.abs() > win || mark[idx(nx, ny)] { continue; }
                if fine.at((x as i64 + nx) as usize, (y as i64 + ny) as usize) <= half { mark[idx(nx, ny)] = true; st.push((nx, ny)); } }
        }
        let (bw, bh) = ((x1 - x0 + 1) as f32, (y1 - y0 + 1) as f32);
        if rmax > ROUND_CORE_R || bw.max(bh) / bw.min(bh).max(1.0) > ROUND_CORE_ASPECT { continue; }
        let d = 2.0 * rad as f32 / RES as f32;
        out.push(Pit { x: x as f64 + 0.5, y: y as f64 + 0.5, len: d, wid: d, depth, share: depth / (centre + depth).max(1.0), z: 0.0, ring_mask: 0.0, green: 0.0, ring_pale: pale as f32 / seen as f32, runs_off: false, keep: true });
    } }
    out
}

/// Through vias show on both faces, so real pits line up with what the other face shows and grain does not. How many of `a` have something of `b` within 0.1 mm, how many would by chance (the same count with `b` mirrored
/// across the frame each way, the higher one), and whether that is decisive: at least a fifth of them and five times chance. CM5 everypin held-out tiles: about two thirds against under one percent. Driver top, whose
/// mask is flat with no via to see: 3 of 684 against 9, so its dark specks are grain. `wh` is the frame in pixels; `a` and `b` are in pixels.
pub fn pits_confirmed(a: &[(f64, f64)], b: &[(f64, f64)], wh: (f64, f64)) -> (usize, usize, bool) {
    let near = 0.1 * RES;
    let count = |other: &[(f64, f64)]| a.iter().filter(|p| other.iter().any(|q| crate::geom::dist(**p, *q) < near)).count();
    let hit = count(b);
    let mx: Vec<(f64, f64)> = b.iter().map(|q| (wh.0 - q.0, q.1)).collect();
    let my: Vec<(f64, f64)> = b.iter().map(|q| (q.0, wh.1 - q.1)).collect();
    let chance = count(&mx).max(count(&my));
    (hit, chance, hit * 5 >= a.len() && hit >= 5 * (chance + 1))
}

pub fn pit_json(p: &Pit, side: &str) -> Value {
    let (mx, my) = crate::board::to_mm(p.x, p.y); let r = |v: f32| (v as f64 * 1000.0).round() / 1000.0;
    json!({"side": side, "pos_mm": [(mx * 1000.0).round() / 1000.0, (my * 1000.0).round() / 1000.0], "len": r(p.len), "wid": r(p.wid), "depth": r(p.depth), "share": r(p.share), "z": r(p.z), "ringMask": r(p.ring_mask), "ringPale": r(p.ring_pale), "green": r(p.green), "runsOff": p.runs_off, "pit": p.keep})
}

#[cfg(test)]
mod tests {
    use super::*;
    #[test]
    fn pits_count_only_when_the_other_face_agrees() {
        let a: Vec<(f64, f64)> = (0..100).map(|k| (37.0 + (k * 211 % 1500) as f64, 23.0 + (k * 97 % 1100) as f64)).collect();
        let b: Vec<(f64, f64)> = a.iter().take(60).map(|p| (p.0 + 1.0, p.1 - 1.0)).collect();
        let (hit, _, sure) = pits_confirmed(&a, &b, (1600.0, 1200.0));
        assert_eq!(hit, 60); assert!(sure);
        let far: Vec<(f64, f64)> = a.iter().map(|p| (p.0 + 11.0, p.1 + 13.0)).collect();
        assert!(!pits_confirmed(&a, &far, (1600.0, 1200.0)).2, "nothing lines up");
        assert!(!pits_confirmed(&a, &[], (1600.0, 1200.0)).2, "one face cannot confirm itself");
    }

    #[test]
    fn a_via_in_a_clearance_is_read_and_a_gap_between_tracks_is_not() {
        // a dark centre, a pale pad and a dark halo; and beside it a long dark line
        let img = RgbImage::from_fn(300, 200, |x, y| {
            let d = (((x as f32 - 80.0).powi(2) + (y as f32 - 100.0).powi(2)) as f32).sqrt();
            if d < 3.5 { image::Rgb([20, 58, 30]) } else if d < 9.0 { image::Rgb([70, 140, 80]) } else if d < 14.0 { image::Rgb([24, 80, 50]) }
            else if (195..199).contains(&x) { image::Rgb([24, 80, 50]) } else { image::Rgb([70, 140, 80]) } });
        let p = round_pits(&img);
        assert_eq!(p.len(), 1, "{:?}", p.iter().map(|q| (q.x, q.y)).collect::<Vec<_>>());
        assert!((p[0].x - 80.5).abs() < 2.0 && (p[0].y - 100.5).abs() < 2.0);
    }

    #[test]
    fn pits_are_read_and_lines_are_not() {
        // a bright mask with 100 dark round spots, a dark gap line, and a grey part: the spots are vias, the line and part are not
        let (w, h) = (1600u32, 1200u32);
        let mut img = RgbImage::from_pixel(w, h, image::Rgb([60, 190, 140]));
        let mut at = vec![];
        for k in 0..100i64 { at.push((60 + (k % 10) * 90, 60 + (k / 10) * 90)); }
        for (cx, cy) in &at { for dy in -6i64..=6 { for dx in -6i64..=6 { let d = ((dx * dx + dy * dy) as f64).sqrt(); if d <= 5.0 { img.put_pixel((*cx + dx) as u32, (*cy + dy) as u32, image::Rgb([10, 70, 50])); } } } }
        for x in 1000..1500u32 { for y in 300..305u32 { img.put_pixel(x, y, image::Rgb([20, 100, 70])); } }
        for x in 1100..1300u32 { for y in 600..800u32 { img.put_pixel(x, y, image::Rgb([50, 50, 50])); } }
        let f = read_pits(&img);
        assert!(f.shows, "{} {}", f.median_mm, f.spread);
        let kept: Vec<&Pit> = f.pits.iter().filter(|p| p.keep).collect();
        assert_eq!(kept.len(), 100, "{}", kept.len());
        assert!(kept.iter().all(|p| p.x < 1000.0), "only the spots");
    }

    use image::Rgb;
    /// A mask-coloured photo with grain, for what is drawn on it.
    fn board(w: u32, h: u32, f: impl Fn(f64, f64) -> Option<f32>) -> RgbImage {
        RgbImage::from_fn(w, h, |x, y| {
            let g = (((x * 7919 + y * 104729) % 23) as f32) - 11.0;
            let v = f(x as f64, y as f64).unwrap_or(70.0) + g;
            Rgb([(v * 0.6).clamp(0.0, 255.0) as u8, v.clamp(0.0, 255.0) as u8, (v * 0.8).clamp(0.0, 255.0) as u8])
        })
    }
    fn dome_at(cx: f64, cy: f64, x: f64, y: f64) -> Option<f32> {
        let (dx, dy) = (x - cx, y - cy); let d = (dx * dx + dy * dy).sqrt();
        if d < 5.0 { Some(210.0) } else if d < 10.0 && dy > 0.0 { Some(35.0) } else { None }
    }
    #[test]
    fn domes_are_kept_and_found_in_the_right_place() {
        let spots = [(60.0, 60.0), (120.0, 60.0), (60.0, 120.0), (120.0, 120.0), (180.0, 60.0), (180.0, 120.0)];
        let img = board(260, 200, |x, y| spots.iter().find_map(|s| dome_at(s.0, s.1, x, y)));
        let d = finish(&measure(&img));
        assert_eq!(d.len(), spots.len(), "{d:?}");
        for s in spots { assert!(d.iter().any(|v| crate::geom::dist((v.x, v.y), s) < 3.0), "{s:?} {d:?}"); }
    }
    #[test]
    fn a_straight_line_and_a_wide_patch_are_not_domes() {
        let img = board(320, 200, |x, y| if (y - 60.0).abs() < 3.0 && x > 20.0 && x < 300.0 { Some(220.0) } else if x > 100.0 && x < 200.0 && y > 110.0 && y < 170.0 { Some(200.0) } else { None });
        let d = finish(&measure(&img));
        assert!(d.is_empty(), "{d:?}");
    }
    #[test]
    fn a_face_without_vias_shows_no_signature() {
        let at: Vec<(f64, f64)> = (0..30).map(|k| (102.0 + k as f64 * 0.2, 102.0)).collect();
        let flat = board(400, 200, |_, _| None);
        assert!(face_signature(&flat, &at).0.abs() < SIGNATURE_MIN);
        let dots: Vec<(f64, f64)> = at.iter().map(|p| crate::board::to_px(p.0, p.1)).collect();
        let lit = board(400, 200, |x, y| if dots.iter().any(|d| ((x - d.0).powi(2) + (y - d.1).powi(2)).sqrt() < 3.0) { Some(140.0) } else { None });
        assert!(face_signature(&lit, &at).0 > SIGNATURE_MIN);
    }
}

/// A via's pad: at each radius the median brightness round the via, so a track leaving it covers only a few sectors. Where the surroundings are laminate the median falls off the pad's edge, and the half
/// height of that fall is the pad's radius. Returns the diameter in mm, or None where the surroundings are copper and nothing can be read.
fn pad_diameter(g: &crate::raster::Plane, x: f64, y: f64, min_step: f32) -> Option<f64> {
    const SECTORS: usize = 48;
    let at = |r: f64, k: usize| { let a = k as f64 / SECTORS as f64 * std::f64::consts::TAU; g.sample((x + r * a.cos()) as f32, (y + r * a.sin()) as f32) };
    let med: Vec<f32> = (0..40).map(|r| { let mut v: Vec<f32> = (0..SECTORS).map(|k| at(r as f64, k)).collect(); v.sort_by(|a, b| a.partial_cmp(b).unwrap()); v[SECTORS / 2] }).collect();
    let mid = |a: usize, b: usize| { let mut v = med[a..b].to_vec(); v.sort_by(|p, q| p.partial_cmp(q).unwrap()); v[v.len() / 2] };
    let (inner, outer) = (mid(6, 11), mid(26, 40));
    if inner - outer < min_step { return None; }
    let half = (inner + outer) / 2.0;
    let k = (6..40).find(|i| med[*i] < half)?;
    let f = (med[k - 1] - half) / (med[k - 1] - med[k] + 1e-6);
    Some(2.0 * (k as f64 - 1.0 + f as f64) / crate::board::RES)
}

/// Via pads on the face where domes show. A dome's width is the dimple over its hole and a dark ring is the clearance round the pad, so neither is the pad. A via whose surroundings are laminate shows the pad
/// itself. Those are the evidence: on the Driver bottom they cluster round one size whatever the via was read as, so every via takes their median, to the nearest 0.05 mm (the photo cannot tell 0.48 from 0.50).
/// The drill stays assumed: `drill` from the rules, cut down to leave a 0.1 mm ring.
pub fn size_vias(img: &image::RgbImage, second: Option<&image::RgbImage>, vias: &mut serde_json::Value, drill: f64) -> (serde_json::Value, String) {
    let (w, h) = (img.width() as usize, img.height() as usize);
    let mean = crate::raster::Plane { w, h, d: (0..w * h).map(|i| { let (x, y) = ((i % w) as u32, (i / w) as u32); let a = img.get_pixel(x, y)[1] as f32; let b = second.map(|s| s.get_pixel(x, y)[1] as f32).unwrap_or(a); ((a + b) / 2.0).max(2.0) }).collect() };
    let lg = { let b = crate::raster::blur(&mean, 1, 2); crate::raster::Plane { w, h, d: b.d.iter().map(|v| v.ln()).collect() } };
    // the fall must be a fifth of a log unit, about half the Driver bottom's copper-to-laminate separation; a bright mask has its own separation (CM5: 0.16), and half of that is the same bar
    let tone = crate::prep::mask_tone(img);
    let min_step = if tone.bright() { let plain: Vec<f32> = (0..w * h).step_by(3).filter(|i| { let (x, y) = ((*i % w) as u32, (*i / w) as u32); let p = img.get_pixel(x, y); (p[1] as f32 - p[0] as f32) > 0.3 * tone.gr50 && p[1] > p[2] }).map(|i| lg.d[i]).collect(); crate::trace::two_modes(&plain).map(|m| 0.5 * (m.1 - m.0)).unwrap_or(0.2) } else { 0.2 };
    let mut read = vec![];
    for v in vias.as_array().into_iter().flatten() { let (Some(x), Some(y)) = (v["pos_mm"][0].as_f64(), v["pos_mm"][1].as_f64()) else { continue }; let (u, t) = crate::board::to_px(x, y); if u < 45.0 || t < 45.0 || u > w as f64 - 45.0 || t > h as f64 - 45.0 { continue; }
        if let Some(d) = pad_diameter(&lg, u, t, min_step) { read.push(serde_json::json!({"id": v["id"], "shape": v["shape"], "padMm": (d * 1000.0).round() / 1000.0})); } }
    let mut ds: Vec<f64> = read.iter().filter_map(|r| r["padMm"].as_f64()).collect(); ds.sort_by(|a, b| a.partial_cmp(b).unwrap());
    if ds.len() < 10 { return (serde_json::json!({"read": read, "note": "too few vias show their pad: sizes left as read"}), format!("{} vias show their pad: too few, sizes left as read", ds.len())); }
    let q = |p: f64| ds[((ds.len() - 1) as f64 * p) as usize];
    // one size is the premise: on the Driver bottom the tenth and ninetieth percentiles differ by 1.55 times. Where they differ by more than double, what fell off was not one pad size (CM5 bottom: 0.32 to 1.17 mm, tracks and pour edges among them), so nothing is written
    if q(0.9) > 2.0 * q(0.1) { return (serde_json::json!({"read": read, "p10Mm": q(0.1), "p90Mm": q(0.9), "note": "the falls read are not one size: sizes left as read"}), format!("{} vias show a fall, from {:.2} to {:.2} mm: not one pad size, sizes left as read", ds.len(), q(0.1), q(0.9))); }
    let med = q(0.5); let pad = (med * 20.0).round() / 20.0; let dr = drill.min(pad - 0.2).max(0.1);
    if let Some(a) = vias.as_array_mut() { for v in a.iter_mut() { v["diameter_mm"] = serde_json::json!(pad); v["drill_mm"] = serde_json::json!((dr * 1000.0).round() / 1000.0); v["sizeTags"] = serde_json::json!({"diameter": "PADS_ON_LAMINATE_MEDIAN", "drill": "ASSUMED"}); } }
    (serde_json::json!({"padMm": pad, "medianMm": med, "p10Mm": q(0.1), "p90Mm": q(0.9), "n": ds.len(), "read": read, "method": "half height of the fall of median brightness round each via on the bottom photos, at vias whose surroundings are laminate; see domes::pad_diameter"}), format!("via pads: {} vias show theirs, median {:.3} mm (10th to 90th percentile {:.3} to {:.3}), written as {pad:.2} mm", ds.len(), med, q(0.1), q(0.9)))
}

#[cfg(test)]
mod pad_size_tests {
    use super::*;
    #[test]
    fn a_pad_on_laminate_reads_its_diameter_and_one_in_copper_reads_nothing() {
        let (w, h) = (200usize, 200usize);
        let paint = |lit: f32, ground: f32| crate::raster::Plane { w, h, d: (0..w * h).map(|i| { let (x, y) = ((i % w) as f32 - 100.0, (i / w) as f32 - 100.0); if (x * x + y * y).sqrt() <= 10.0 { lit } else { ground } }).collect() };
        let d = pad_diameter(&paint(4.3, 3.7), 100.0, 100.0, 0.2).expect("a fall to read");
        assert!((d - 0.5).abs() < 0.05, "{d}");
        assert!(pad_diameter(&paint(4.3, 4.3), 100.0, 100.0, 0.2).is_none());
    }
}