diff --git a/.gitignore b/.gitignore index 24748fc..cbec81b 100644 --- a/.gitignore +++ b/.gitignore @@ -15,3 +15,17 @@ thumbcache/ legacy/ *.pdf run-*.sh + +# demo page source assets: re-include web/img, which the broad img/ and +# *.png patterns above would otherwise swallow +!web/img/ +!web/img/* + +# demo page build artifacts, regenerated by `cargo run --release --bin +# demoassets`, `cargo run --release --bin storekeys` and +# `python3 pr_to_json.py` (needed when deploying the page) +web/flickr/ +web/pr-data.json + +# playwright-core for verify-retrieval.mjs (npm install --no-save playwright-core) +node_modules/ diff --git a/Cargo.toml b/Cargo.toml index 9e66993..c45ea46 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -52,3 +52,18 @@ required-features = ["cli"] name = "bench" path = "src/bin/bench.rs" required-features = ["cli"] + +[[bin]] +name = "demoassets" +path = "src/bin/demoassets.rs" +required-features = ["cli"] + +[[bin]] +name = "storekeys" +path = "src/bin/storekeys.rs" +required-features = ["cli"] + +[[bin]] +name = "salientlab" +path = "src/bin/salientlab.rs" +required-features = ["cli"] diff --git a/src/bin/demoassets.rs b/src/bin/demoassets.rs new file mode 100644 index 0000000..f34d94f --- /dev/null +++ b/src/bin/demoassets.rs @@ -0,0 +1,69 @@ +//! Generates the static assets for the web demo ranking pool: +//! a thumbnail per image plus hashes.json with the three hashes, +//! computed from the full-size images. +use clap::Parser; +use image_similarity::descriptors::{Descriptor, Median, DCT, PHash}; +use image_similarity::open_image; +use rayon::prelude::*; +use std::fs; +use std::io::Write; +use std::path::PathBuf; + +#[derive(Parser)] +struct Args { + /// Directory with source images + #[arg(default_value = "smolflickr")] + source: PathBuf, + /// Output directory, thumbnails go in /thumbs + #[arg(default_value = "web/flickr")] + out: PathBuf, + /// Thumbnail bounding box in pixels + #[arg(long, default_value_t = 140)] + thumb: u32, +} + +fn main() { + let args = Args::parse(); + let thumbs = args.out.join("thumbs"); + fs::create_dir_all(&thumbs).expect("creating output dir"); + + let mut names: Vec = fs::read_dir(&args.source) + .expect("reading source dir") + .filter_map(|entry| { + let entry = entry.ok()?; + entry.path().is_file().then(|| entry.file_name().to_string_lossy().into_owned()) + }) + .collect(); + names.sort(); + + let entries: Vec = names + .par_iter() + .filter_map(|name| { + let img = match open_image(args.source.join(name)) { + Ok(img) => img, + Err(e) => { + eprintln!("skipping {name}: {e}"); + return None; + } + }; + let dct = DCT::new().describe(&img); + let median = Median.describe(&img); + let phash = PHash.describe(&img); + + let thumb = img.thumbnail(args.thumb, args.thumb); + let mut buffer = std::io::Cursor::new(Vec::new()); + thumb + .to_rgb8() + .write_with_encoder(image::codecs::jpeg::JpegEncoder::new_with_quality(&mut buffer, 80)) + .expect("encoding thumbnail"); + fs::write(thumbs.join(name), buffer.get_ref()).expect("writing thumbnail"); + + Some(format!("[\"{name}\",\"{dct:016x}\",\"{median:016x}\",\"{phash:016x}\"]")) + }) + .collect(); + + let json = format!("[\n{}\n]\n", entries.join(",\n")); + let mut file = fs::File::create(args.out.join("hashes.json")).expect("creating hashes.json"); + file.write_all(json.as_bytes()).expect("writing hashes.json"); + println!("{} images -> {}", entries.len(), args.out.display()); +} diff --git a/src/bin/salientlab.rs b/src/bin/salientlab.rs new file mode 100644 index 0000000..ef89fb2 --- /dev/null +++ b/src/bin/salientlab.rs @@ -0,0 +1,853 @@ +//! Throwaway lab for the salient-point hash variants. +//! +//! Two modes: +//! - default: small-sample table of per-mutation Hamming distances plus +//! detector repeatability, quick iteration on a few hundred images +//! - --pr: the full run. Calibrates on the first --calib images, hashes +//! every image and mutant under all variants (checkpointed to +//! --hashes with resume), then computes PR curves with the same +//! tally semantics as DescriptorStore::get_stats and plots one PNG +//! per variant. The dct hash rides along as a reference method. +//! +//! Not part of the real benchmark pipeline. + +use clap::Parser; +use image_similarity::descriptors::{Descriptor, DCT}; +use image_similarity::mutators::{get_all_mutators, CropCenter, Mutator, Rotate}; +use image_similarity::open_image; +use image_similarity::salient::{ + binarize64, hash_strongest, pooled_vector, Feature, Salient, +}; +use plotters::prelude::*; +use rayon::prelude::*; +use serde::{Deserialize, Serialize}; +use std::collections::HashSet; +use std::path::{Path, PathBuf}; +use std::sync::atomic::{AtomicUsize, Ordering}; +use std::time::Instant; + +#[derive(Parser)] +#[command(about = "Salient-point hash experiments")] +struct Args { + /// Directory with base images + #[arg(long, default_value = "smolflickr")] + dir: PathBuf, + /// Number of images to sample (default: 200 for the table, all for --pr) + #[arg(long)] + count: Option, + /// Keypoints kept per image + #[arg(long, default_value_t = 100)] + features: usize, + /// Keypoints considered for the repeatability measurement + #[arg(long, default_value_t = 50)] + repeat_n: usize, + /// Full PR-curve run over the whole directory + #[arg(long, default_value_t = false)] + pr: bool, + /// Images used to calibrate codebook and thresholds in --pr mode + #[arg(long, default_value_t = 2000)] + calib: usize, + /// Checkpoint file for the hashes in --pr mode + #[arg(long, default_value = "salient-hashes.mp")] + hashes: PathBuf, + /// Skip hashing, compute PR curves from an existing checkpoint file + #[arg(long, default_value_t = false)] + stats_only: bool, +} + +const VARIANTS: [&str; 8] = [ + "pooled", "strongest", "consensus", "pooled-cal", "pooled-sim", "cons-cal", "words-any", + "words-med", +]; +/// All hashed methods: the salient variants plus the dct reference +const METHODS: usize = VARIANTS.len() + 1; + +fn method_name(v: usize) -> &'static str { + if v < VARIANTS.len() { + VARIANTS[v] + } else { + "dct-ref" + } +} + +/// Per image-version summary statistics that the hash variants read from +struct Stats { + /// Response-weighted mean descriptor, unit normalized + pooled: [f32; 64], + /// Fraction of keypoints voting 1 on each binarize64 bit position + votes: [f32; 64], + /// Keypoints assigned to each of the 64 visual words + counts: [u32; 64], +} + +struct Calibration { + /// The visual-word codebook + centroids: Vec<[f32; 64]>, + /// Per-dimension corpus median of the pooled vector + pooled_med: [f32; 64], + /// Corpus median vote fraction per bit + vote_med: [f32; 64], + /// Corpus median count per visual word + word_med: [u32; 64], + /// Random +-1 hyperplanes for the simhash variant + planes: Vec<[f32; 64]>, +} + +/// Checkpoint file contents for --pr mode. Per image: one hash per method +/// for the base image followed by every mutant, in mutator order. +#[derive(Serialize, Deserialize)] +struct HashFile { + tags: Vec, + results: Vec<(String, Vec<[u64; METHODS]>)>, +} + +fn main() { + let args = Args::parse(); + if args.pr { + run_pr(&args); + } else { + run_table(&args); + } +} + +/// The standard suite (which now carries a hard rotate45) plus a mid +/// rotation and a heavy crop, to trace the geometric falloff more finely +fn lab_mutators() -> Vec> { + let mut mutators = get_all_mutators(); + mutators.push(Box::new(Rotate { degrees: 15.0 })); + mutators.push(Box::new(CropCenter { keep_percent: 50 })); + mutators +} + +fn image_paths(dir: &Path, count: usize) -> Vec<(String, PathBuf)> { + let mut paths: Vec<(String, PathBuf)> = std::fs::read_dir(dir) + .expect("cannot read image dir") + .filter_map(|e| e.ok()) + .filter(|e| { + matches!( + e.path().extension().and_then(|x| x.to_str()), + Some("jpg") | Some("jpeg") | Some("png") + ) + }) + .filter_map(|e| e.file_name().into_string().ok().map(|n| (n, e.path()))) + .collect(); + paths.sort(); + paths.truncate(count); + assert!(!paths.is_empty(), "no images found in {:?}", dir); + paths +} + +// --------------------------------------------------------------------------- +// full PR run + +fn run_pr(args: &Args) { + let paths = image_paths(&args.dir, args.count.unwrap_or(usize::MAX)); + let salient = Salient { + params: image_similarity::salient::Params { + max_features: args.features, + ..Default::default() + }, + }; + let mutators = lab_mutators(); + let tags: Vec = mutators.iter().map(|m| m.tag()).collect(); + + let mut file = load_hashes(&args.hashes).unwrap_or(HashFile { + tags: tags.clone(), + results: Vec::new(), + }); + assert_eq!( + file.tags, tags, + "{} was made with a different mutator suite, move it away first", + args.hashes.display() + ); + + if !args.stats_only { + // calibration is deterministic, so recomputing it on resume is safe + let calib = calibrate(&paths[..args.calib.min(paths.len())], &salient); + + let done: HashSet<&String> = file.results.iter().map(|(n, _)| n).collect(); + let todo: Vec<&(String, PathBuf)> = + paths.iter().filter(|(n, _)| !done.contains(n)).collect(); + drop(done); + println!( + "{} images total, {} still to hash, {} mutations, {} methods", + paths.len(), + todo.len(), + mutators.len(), + METHODS + ); + + let t0 = Instant::now(); + let mut last_save = Instant::now(); + let mut hashed = 0usize; + for batch in todo.chunks(512) { + let out: Vec<(String, Vec<[u64; METHODS]>)> = batch + .par_iter() + .filter_map(|(name, path)| { + let img = open_image(path).ok()?; + let mut versions = Vec::with_capacity(mutators.len() + 1); + versions.push(hash_methods(&img, &salient, &calib)); + for m in &mutators { + versions.push(hash_methods(&m.mutate(&img), &salient, &calib)); + } + Some((name.clone(), versions)) + }) + .collect(); + file.results.extend(out); + hashed += batch.len(); + let rate = hashed as f32 / t0.elapsed().as_secs_f32(); + eprintln!( + " {} / {} images, {:.1} img/s, ~{:.0} min left", + hashed, + todo.len(), + rate, + (todo.len() - hashed) as f32 / rate / 60.0 + ); + if last_save.elapsed().as_secs() >= 120 { + save_hashes(&args.hashes, &file); + last_save = Instant::now(); + } + } + save_hashes(&args.hashes, &file); + println!( + "hashing done in {:.1} min, {} images in {}", + t0.elapsed().as_secs_f32() / 60.0, + file.results.len(), + args.hashes.display() + ); + } + + pr_stats(&file, &mutators); +} + +fn calibrate(paths: &[(String, PathBuf)], salient: &Salient) -> Calibration { + println!("calibrating on {} images...", paths.len()); + let t0 = Instant::now(); + let feats: Vec> = paths + .par_iter() + .filter_map(|(_, path)| Some(salient.features(&open_image(path).ok()?))) + .collect(); + let samples: Vec<[f32; 64]> = feats + .iter() + .flat_map(|f| f.iter().take(50).map(|f| f.desc)) + .collect(); + let centroids = kmeans(&samples, 64, 25); + let stats: Vec = feats + .par_iter() + .map(|f| stats_for(f, ¢roids)) + .collect(); + let calib = Calibration { + pooled_med: column_median(&stats, |s| s.pooled), + vote_med: column_median(&stats, |s| s.votes), + word_med: column_median_u32(&stats), + planes: random_planes(64), + centroids, + }; + println!( + "calibration took {:.1}s ({} descriptors)", + t0.elapsed().as_secs_f32(), + samples.len() + ); + calib +} + +fn hash_methods(img: &image::DynamicImage, salient: &Salient, calib: &Calibration) -> [u64; METHODS] { + let feats = salient.features(img); + let stats = stats_for(&feats, &calib.centroids); + let v = variant_hashes(&feats, &stats, calib); + let mut out = [0u64; METHODS]; + out[..VARIANTS.len()].copy_from_slice(&v); + out[VARIANTS.len()] = DCT::new().describe(img); + out +} + +fn load_hashes(path: &Path) -> Option { + let bytes = std::fs::read(path).ok()?; + match rmp_serde::from_slice(&bytes) { + Ok(file) => Some(file), + Err(e) => panic!("cannot parse {}: {e}", path.display()), + } +} + +fn save_hashes(path: &Path, file: &HashFile) { + let bytes = rmp_serde::to_vec(file).expect("serialize hashes"); + let tmp = path.with_extension("tmp"); + std::fs::write(&tmp, bytes).expect("write hashes"); + std::fs::rename(&tmp, path).expect("rename hashes"); +} + +/// Same tally as DescriptorStore::get_stats, done by brute force instead of +/// the BK-tree: every base compares against every other hash. A hit is a +/// true positive when it is the base's own mutant, a false positive for +/// mutator m when it is an m-mutant of another image, and a false positive +/// for every mutator when it is a different base image. +fn pr_stats(file: &HashFile, mutators: &[Box]) { + let max_threshold = 24usize; + let n = file.results.len(); + let nm = mutators.len(); + println!("\ncomputing PR curves over {n} images..."); + + for v in 0..METHODS { + let t0 = Instant::now(); + let bases: Vec = file.results.iter().map(|(_, h)| h[0][v]).collect(); + let muts: Vec> = (0..nm) + .map(|m| file.results.iter().map(|(_, h)| h[m + 1][v]).collect()) + .collect(); + + #[derive(Clone)] + struct Hists { + // 65 slots: two hashes can differ in all 64 bits + base_fp: [u64; 65], + tp: Vec<[u64; 65]>, + fp: Vec<[u64; 65]>, + } + let empty = || Hists { + base_fp: [0; 65], + tp: vec![[0; 65]; nm], + fp: vec![[0; 65]; nm], + }; + let hists = (0..n) + .into_par_iter() + .fold(empty, |mut h, i| { + let hb = bases[i]; + for (j, other) in bases.iter().enumerate() { + if j != i { + h.base_fp[(hb ^ other).count_ones() as usize] += 1; + } + } + for m in 0..nm { + let mm = &muts[m]; + h.tp[m][(hb ^ mm[i]).count_ones() as usize] += 1; + for (j, other) in mm.iter().enumerate() { + if j != i { + h.fp[m][(hb ^ other).count_ones() as usize] += 1; + } + } + } + h + }) + .reduce(empty, |mut a, b| { + for (x, y) in a.base_fp.iter_mut().zip(&b.base_fp) { + *x += y; + } + for m in 0..nm { + for (x, y) in a.tp[m].iter_mut().zip(&b.tp[m]) { + *x += y; + } + for (x, y) in a.fp[m].iter_mut().zip(&b.fp[m]) { + *x += y; + } + } + a + }); + + let cumulative = |at: &[u64; 65]| { + let mut cum = [0u64; 65]; + let mut sum = 0; + for (t, count) in at.iter().enumerate() { + sum += count; + cum[t] = sum; + } + cum + }; + let fp_base = cumulative(&hists.base_fp); + + println!( + "\nStats for salient {} ({:.1}s)", + method_name(v), + t0.elapsed().as_secs_f32() + ); + let curves: Vec<(String, Vec<(f64, f64)>)> = (0..nm) + .map(|m| { + let tp = cumulative(&hists.tp[m]); + let fp_mut = cumulative(&hists.fp[m]); + let curve: Vec<(f64, f64)> = (0..max_threshold) + .map(|t| { + let tpos = tp[t] as f64; + let fpos = (fp_mut[t] + fp_base[t]) as f64; + let recall = tpos / n as f64; + let precision = if tpos + fpos == 0.0 { + 1.0 + } else { + tpos / (tpos + fpos) + }; + (recall, precision) + }) + .collect(); + (mutators[m].info(), curve) + }) + .collect(); + for (name, curve) in &curves { + println!("{name}:"); + for (r, p) in curve { + println!("{r}, {p}"); + } + } + plot_pr(&format!("sal-{}-pr.png", method_name(v)), method_name(v), &curves); + } + println!("\nwrote sal--pr.png for {METHODS} methods"); +} + +/// Same look as the PR plots from main.rs +fn plot_pr(filename: &str, method: &str, curves: &[(String, Vec<(f64, f64)>)]) { + let root = BitMapBackend::new(filename, (1024, 768)).into_drawing_area(); + root.fill(&WHITE).unwrap(); + let mut chart = ChartBuilder::on(&root) + .caption( + format!("PR curves (salient {method})"), + ("sans-serif", 40).into_font(), + ) + .margin(25) + .set_all_label_area_size(50) + .build_cartesian_2d(0f64..1f64, 0f64..1f64) + .unwrap(); + chart + .configure_mesh() + .x_labels(10) + .y_labels(10) + .disable_mesh() + .x_desc("Recall") + .y_desc("Precision") + .x_label_formatter(&|x| format!("{x:.3}")) + .y_label_formatter(&|x| format!("{x:.3}")) + .draw() + .unwrap(); + for (i, (name, curve)) in curves.iter().enumerate() { + let color = Palette99::pick(i).to_rgba(); + chart + .draw_series(LineSeries::new(curve.clone(), color.filled()).point_size(2)) + .unwrap() + .label(name) + .legend(move |(x, y)| PathElement::new(vec![(x, y), (x + 20, y)], color)); + } + chart + .configure_series_labels() + .position(SeriesLabelPosition::MiddleLeft) + .border_style(BLACK) + .draw() + .unwrap(); +} + +// --------------------------------------------------------------------------- +// small-sample distance table + +struct MutData { + tag: String, + w: f32, + h: f32, + feats: Vec, + dct: u64, +} + +struct ImageData { + w: f32, + h: f32, + base: Vec, + dct: u64, + muts: Vec, +} + +fn run_table(args: &Args) { + let paths = image_paths(&args.dir, args.count.unwrap_or(200)); + let salient = Salient { + params: image_similarity::salient::Params { + max_features: args.features, + ..Default::default() + }, + }; + let mutators = lab_mutators(); + println!( + "{} images, {} mutations, {} keypoints per image", + paths.len(), + mutators.len(), + args.features + ); + + let t0 = Instant::now(); + let done = AtomicUsize::new(0); + let data: Vec = paths + .par_iter() + .filter_map(|(_, path)| { + let img = open_image(path).ok()?; + let dct = DCT::new(); + let base = salient.features(&img); + let muts = mutators + .iter() + .map(|m| { + let mutated = m.mutate(&img); + MutData { + tag: m.tag(), + w: mutated.width() as f32, + h: mutated.height() as f32, + feats: salient.features(&mutated), + dct: dct.describe(&mutated), + } + }) + .collect(); + let n = done.fetch_add(1, Ordering::Relaxed) + 1; + if n % 50 == 0 { + eprintln!(" {} / {} images", n, paths.len()); + } + Some(ImageData { + w: img.width() as f32, + h: img.height() as f32, + base, + dct: dct.describe(&img), + muts, + }) + }) + .collect(); + let avg_feats: f32 = + data.iter().map(|d| d.base.len() as f32).sum::() / data.len() as f32; + println!( + "extraction: {:.1}s, avg {:.0} keypoints per base image", + t0.elapsed().as_secs_f32(), + avg_feats + ); + + repeatability(&data, &mutators, args.repeat_n); + + // codebook, then per-version stats, then corpus calibration + let samples: Vec<[f32; 64]> = data + .iter() + .flat_map(|d| d.base.iter().take(50).map(|f| f.desc)) + .collect(); + println!("\ntraining 64-word codebook on {} descriptors...", samples.len()); + let centroids = kmeans(&samples, 64, 25); + + let base_stats: Vec = data + .par_iter() + .map(|d| stats_for(&d.base, ¢roids)) + .collect(); + let mut_stats: Vec> = data + .par_iter() + .map(|d| d.muts.iter().map(|m| stats_for(&m.feats, ¢roids)).collect()) + .collect(); + + let calib = Calibration { + pooled_med: column_median(&base_stats, |s| s.pooled), + vote_med: column_median(&base_stats, |s| s.votes), + word_med: column_median_u32(&base_stats), + planes: random_planes(64), + centroids, + }; + + let base_hashes: Vec<[u64; 8]> = data + .iter() + .zip(&base_stats) + .map(|(d, s)| variant_hashes(&d.base, s, &calib)) + .collect(); + let mut_hashes: Vec> = data + .iter() + .zip(&mut_stats) + .map(|(d, stats)| { + d.muts + .iter() + .zip(stats) + .map(|(m, s)| variant_hashes(&m.feats, s, &calib)) + .collect() + }) + .collect(); + + println!("\nhamming distance base vs mutant, mean/p90 (dct hash as reference):"); + print!(" {:12} {:>10}", "mutation", "dct"); + for v in VARIANTS { + print!(" {:>10}", v); + } + println!(); + for (mi, m) in mutators.iter().enumerate() { + print!(" {:12}", m.tag()); + let mut dists: Vec = data + .iter() + .map(|d| (d.dct ^ d.muts[mi].dct).count_ones()) + .collect(); + print_cell(&mut dists, false); + for v in 0..VARIANTS.len() { + let mut dists: Vec = (0..data.len()) + .map(|i| (base_hashes[i][v] ^ mut_hashes[i][mi][v]).count_ones()) + .collect(); + print_cell(&mut dists, false); + } + println!(); + } + + // impostor baseline: unrelated image pairs; p5 shows how close the + // nearest false positives come + print!(" {:12}", "impostor*"); + let mut dists: Vec = (0..data.len()) + .map(|i| (data[i].dct ^ data[(i + 1) % data.len()].dct).count_ones()) + .collect(); + print_cell(&mut dists, true); + for v in 0..VARIANTS.len() { + let mut dists: Vec = (0..data.len()) + .map(|i| { + let j = (i + 1) % data.len(); + (base_hashes[i][v] ^ base_hashes[j][v]).count_ones() + }) + .collect(); + print_cell(&mut dists, true); + } + println!("\n *impostor row is mean/p5 over unrelated pairs: match rows want"); + println!(" to sit far below it, especially below its p5"); +} + +/// One table cell: mean/p90 for match rows, mean/p5 for the impostor row +fn print_cell(dists: &mut Vec, impostor: bool) { + dists.sort(); + let mean = dists.iter().sum::() as f32 / dists.len() as f32; + let q = if impostor { + dists[dists.len() / 20] + } else { + dists[(dists.len() * 9 / 10).min(dists.len() - 1)] + }; + print!(" {:>6.1}/{:<3}", mean, q); +} + +fn repeatability(data: &[ImageData], mutators: &[Box], top: usize) { + println!("\ndetector repeatability (top {} keypoints, 4px tolerance):", top); + for (mi, m) in mutators.iter().enumerate() { + let (mut matched, mut total) = (0usize, 0usize); + for img in data { + let md = &img.muts[mi]; + let mapped: Vec<(f32, f32, f32)> = md + .feats + .iter() + .take(top) + .map(|f| map_to_base(&md.tag, img.w, img.h, md.w, md.h, f.kp.x, f.kp.y, f.kp.sigma)) + .collect(); + for f in img.base.iter().take(top) { + if !in_region(&md.tag, img.w, img.h, md.w, md.h, f.kp.x, f.kp.y) { + continue; + } + total += 1; + let hit = mapped.iter().any(|(x, y, s)| { + let (dx, dy) = (x - f.kp.x, y - f.kp.y); + let ratio = (s / f.kp.sigma).max(f.kp.sigma / s); + dx * dx + dy * dy <= 16.0 && ratio <= 1.6 + }); + if hit { + matched += 1; + } + } + } + println!( + " {:12} {:.2}", + m.tag(), + matched as f32 / total.max(1) as f32 + ); + } +} + +// --------------------------------------------------------------------------- +// shared variant machinery + +fn stats_for(feats: &[Feature], centroids: &[[f32; 64]]) -> Stats { + let pooled = pooled_vector(feats); + let mut votes = [0f32; 64]; + for f in feats { + let h = binarize64(&f.desc); + for (bit, vote) in votes.iter_mut().enumerate() { + if (h >> (63 - bit)) & 1 == 1 { + *vote += 1.0; + } + } + } + if !feats.is_empty() { + for v in votes.iter_mut() { + *v /= feats.len() as f32; + } + } + let mut counts = [0u32; 64]; + for f in feats { + counts[nearest(centroids, &f.desc)] += 1; + } + Stats { pooled, votes, counts } +} + +fn variant_hashes(feats: &[Feature], stats: &Stats, calib: &Calibration) -> [u64; 8] { + let half = [0.5f32; 64]; + let mut centered = [0f32; 64]; + for i in 0..64 { + centered[i] = stats.pooled[i] - calib.pooled_med[i]; + } + let mut simhash = 0u64; + for plane in &calib.planes { + simhash <<= 1; + let dot: f32 = plane.iter().zip(¢ered).map(|(p, v)| p * v).sum(); + if dot > 0.0 { + simhash |= 1; + } + } + let (mut any, mut med) = (0u64, 0u64); + for w in 0..64 { + any <<= 1; + med <<= 1; + if stats.counts[w] > 0 { + any |= 1; + } + if stats.counts[w] > calib.word_med[w] { + med |= 1; + } + } + [ + binarize64(&stats.pooled), + hash_strongest(feats), + threshold_bits(&stats.votes, &half), + threshold_bits(&stats.pooled, &calib.pooled_med), + simhash, + threshold_bits(&stats.votes, &calib.vote_med), + any, + med, + ] +} + +fn threshold_bits(v: &[f32; 64], t: &[f32; 64]) -> u64 { + let mut mask = 0u64; + for i in 0..64 { + mask <<= 1; + if v[i] > t[i] { + mask |= 1; + } + } + mask +} + +fn column_median [f32; 64]>(stats: &[Stats], get: F) -> [f32; 64] { + let mut out = [0f32; 64]; + for (i, o) in out.iter_mut().enumerate() { + let mut col: Vec = stats.iter().map(|s| get(s)[i]).collect(); + col.sort_by(f32::total_cmp); + *o = col[col.len() / 2]; + } + out +} + +fn column_median_u32(stats: &[Stats]) -> [u32; 64] { + let mut out = [0u32; 64]; + for (i, o) in out.iter_mut().enumerate() { + let mut col: Vec = stats.iter().map(|s| s.counts[i]).collect(); + col.sort(); + *o = col[col.len() / 2]; + } + out +} + +fn xorshift(state: &mut u64) -> u64 { + *state ^= *state << 13; + *state ^= *state >> 7; + *state ^= *state << 17; + *state +} + +/// Fixed random +-1 hyperplanes, the classic sign-random-projection simhash +fn random_planes(n: usize) -> Vec<[f32; 64]> { + let mut rng = 0x9E3779B97F4A7C15u64; + (0..n) + .map(|_| { + let mut plane = [0f32; 64]; + for p in plane.iter_mut() { + *p = if xorshift(&mut rng) & 1 == 1 { 1.0 } else { -1.0 }; + } + plane + }) + .collect() +} + +fn nearest(centroids: &[[f32; 64]], d: &[f32; 64]) -> usize { + let mut best = 0; + let mut best_dist = f32::MAX; + for (c, centroid) in centroids.iter().enumerate() { + let mut dist = 0.0; + for i in 0..64 { + let diff = centroid[i] - d[i]; + dist += diff * diff; + } + if dist < best_dist { + best_dist = dist; + best = c; + } + } + best +} + +fn kmeans(samples: &[[f32; 64]], k: usize, iters: usize) -> Vec<[f32; 64]> { + let mut rng = 0x243F6A8885A308D3u64; + let mut centroids: Vec<[f32; 64]> = (0..k) + .map(|_| samples[(xorshift(&mut rng) % samples.len() as u64) as usize]) + .collect(); + for _ in 0..iters { + let assign: Vec = samples + .par_iter() + .map(|s| nearest(¢roids, s)) + .collect(); + let mut sums = vec![[0f32; 64]; k]; + let mut counts = vec![0usize; k]; + for (s, &a) in samples.iter().zip(&assign) { + counts[a] += 1; + for i in 0..64 { + sums[a][i] += s[i]; + } + } + for c in 0..k { + if counts[c] == 0 { + centroids[c] = samples[(xorshift(&mut rng) % samples.len() as u64) as usize]; + } else { + for i in 0..64 { + centroids[c][i] = sums[c][i] / counts[c] as f32; + } + } + } + } + centroids +} + +/// Map a keypoint position in the mutant back to base image coordinates. +/// Only the geometric mutations move points; everything else is identity. +fn map_to_base( + tag: &str, + bw: f32, + bh: f32, + mw: f32, + mh: f32, + x: f32, + y: f32, + sigma: f32, +) -> (f32, f32, f32) { + if tag == ".flip." { + return (bw - 1.0 - x, y, sigma); + } + if tag.starts_with(".crop") { + // crop_imm offsets, integer division like the mutator + let ox = ((bw - mw) / 2.0).floor(); + let oy = ((bh - mh) / 2.0).floor(); + return (x + ox, y + oy, sigma); + } + if tag == ".scale50." { + return (x * bw / mw, y * bh / mh, sigma * bw / mw); + } + if let Some(deg) = rot_degrees(tag) { + // imageproc warps output through the projection inverse, so a + // mutant point o sits at base position c + R(-theta) (o - c) + let theta = deg.to_radians(); + let (s, c) = theta.sin_cos(); + let (cx, cy) = (mw / 2.0, mh / 2.0); + let (dx, dy) = (x - cx, y - cy); + return (cx + c * dx + s * dy, cy - s * dx + c * dy, sigma); + } + (x, y, sigma) +} + +fn rot_degrees(tag: &str) -> Option { + tag.strip_prefix(".rot")?.strip_suffix('.')?.parse().ok() +} + +/// Whether a base keypoint can possibly survive the mutation, for a fair +/// repeatability denominator (crops delete everything outside the window) +fn in_region(tag: &str, bw: f32, bh: f32, mw: f32, mh: f32, x: f32, y: f32) -> bool { + if tag.starts_with(".crop") { + let ox = ((bw - mw) / 2.0).floor(); + let oy = ((bh - mh) / 2.0).floor(); + return x >= ox && x < ox + mw && y >= oy && y < oy + mh; + } + true +} diff --git a/src/bin/storekeys.rs b/src/bin/storekeys.rs new file mode 100644 index 0000000..c97a2f1 --- /dev/null +++ b/src/bin/storekeys.rs @@ -0,0 +1,38 @@ +//! Dumps the distinct hashes of a store as little-endian u64s. +//! The web demo's at-scale BK-tree builds from this file. +use clap::Parser; +use image_similarity::descriptors::get_all_descriptors; +use image_similarity::store::DescriptorStore; +use std::fs; +use std::path::PathBuf; + +#[derive(Parser)] +struct Args { + /// Store to read, the file stem picks the descriptor + #[arg(default_value = "dct.store")] + store: PathBuf, + /// Output file of raw little-endian u64 keys + #[arg(default_value = "web/flickr/keys.bin")] + out: PathBuf, +} + +fn main() { + let args = Args::parse(); + let stem = args.store.file_stem().expect("store path has no file stem").to_string_lossy(); + let descriptor = get_all_descriptors() + .into_iter() + .find(|d| d.info() == stem) + .unwrap_or_else(|| panic!("no descriptor named {stem}")); + + let store = DescriptorStore::new(descriptor) + .with_file(&args.store) + .expect("loading store"); + + // sorted for a deterministic file; u64 order carries no hamming + // structure, so the resulting tree shape is as good as random + let mut keys: Vec = store.keys().collect(); + keys.sort_unstable(); + let bytes: Vec = keys.iter().flat_map(|k| k.to_le_bytes()).collect(); + fs::write(&args.out, bytes).expect("writing keys"); + println!("{} keys -> {}", keys.len(), args.out.display()); +} diff --git a/src/bk.rs b/src/bk.rs index 6530551..b9294ff 100644 --- a/src/bk.rs +++ b/src/bk.rs @@ -44,6 +44,28 @@ impl Node { fn count(&self) -> usize { 1 + self.children.iter().map(|(_, c)| c.count()).sum::() } + + /// find, but records the key of every node compared, in visit order + fn find_trace(&self, key: u64, radius: u64, found: &mut Vec<(u64, u64)>, visited: &mut Vec) { + visited.push(self.key); + let distance = hamming(self.key, key); + if distance <= radius { + found.push((self.key, distance)); + } + for (child_distance, child) in &self.children { + let child_distance = *child_distance as u64; + if child_distance + radius >= distance && child_distance <= distance + radius { + child.find_trace(key, radius, found, visited); + } + } + } + + fn walk(&self, parent: u64, distance: u8, out: &mut Vec<(u64, u64, u64)>) { + out.push((self.key, parent, distance as u64)); + for (d, child) in &self.children { + child.walk(self.key, *d, out); + } + } } #[derive(Default)] @@ -80,6 +102,27 @@ impl BkTree { } found } + + /// Like find, but also returns every key that was compared against, + /// in visit order. Keys not in the list were pruned away. + pub fn find_trace(&self, key: u64, radius: u64) -> (Vec<(u64, u64)>, Vec) { + let mut found = Vec::new(); + let mut visited = Vec::new(); + if let Some(root) = &self.root { + root.find_trace(key, radius, &mut found, &mut visited); + } + (found, visited) + } + + /// Preorder walk of the tree as (key, parent key, distance to parent) + /// triples, for visualization. The root lists itself as its parent. + pub fn walk(&self) -> Vec<(u64, u64, u64)> { + let mut out = Vec::new(); + if let Some(root) = &self.root { + root.walk(root.key, 0, &mut out); + } + out + } } #[cfg(test)] @@ -132,4 +175,38 @@ mod tests { fn empty_tree_finds_nothing() { assert!(BkTree::new().find(1, 64).is_empty()); } + + #[test] + fn trace_matches_find_and_visits_at_most_everything() { + let keys = lcg_keys(200); + let mut tree = BkTree::new(); + for &key in &keys { + tree.insert(key); + } + for radius in [0, 4, 16] { + let (found, visited) = tree.find_trace(keys[3], radius); + assert_eq!(found, tree.find(keys[3], radius)); + assert_eq!(visited[0], keys[0], "search starts at the root"); + assert!(visited.len() <= keys.len()); + // every match was necessarily compared against + for (key, _) in &found { + assert!(visited.contains(key)); + } + } + } + + #[test] + fn walk_lists_every_node_with_consistent_distances() { + let keys = lcg_keys(50); + let mut tree = BkTree::new(); + for &key in &keys { + tree.insert(key); + } + let nodes = tree.walk(); + assert_eq!(nodes.len(), tree.len()); + assert_eq!(nodes[0], (keys[0], keys[0], 0), "root is its own parent"); + for (key, parent, distance) in &nodes[1..] { + assert_eq!(hamming(*key, *parent), *distance); + } + } } diff --git a/src/lib.rs b/src/lib.rs index 137a88f..d149420 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -7,6 +7,7 @@ pub mod cache; pub mod descriptors; pub mod store; pub mod mutators; +pub mod salient; pub mod wasm; /// Opens an image by sniffing its content instead of trusting the file diff --git a/src/mutators/mod.rs b/src/mutators/mod.rs index a532f30..45fdaaa 100644 --- a/src/mutators/mod.rs +++ b/src/mutators/mod.rs @@ -281,6 +281,9 @@ pub fn get_all_mutators() -> Vec> { Box::new(Letterbox { bar_percent: 10 }), Box::new(Logo { size_percent: 10 }), Box::new(Rotate { degrees: 2.0 }), + // a hard rotation the frequency hashes cannot follow, where the + // salient-point descriptor's rotation invariance shows itself + Box::new(Rotate { degrees: 45.0 }), Box::new(Noise { stddev: 10.0 }), ] } diff --git a/src/salient.rs b/src/salient.rs new file mode 100644 index 0000000..f292585 --- /dev/null +++ b/src/salient.rs @@ -0,0 +1,644 @@ +//! Salient-point detector and descriptor, built for hashing experiments. +//! +//! The detector is SIFT-like: extrema in a Difference-of-Gaussians scale +//! space, filtered for contrast and edge responses. The descriptor is +//! SURF-like: a 4x4 grid of subregions around the keypoint, each summarized +//! by four gradient sums, 64 floats total. Every descriptor is measured in +//! the keypoint's own scale and orientation frame, which is where the scale +//! and rotation invariance comes from. +//! +//! Descriptors are canonicalized under horizontal mirroring: a flipped patch +//! produces the same 64 floats, so hashes built from them inherit flip +//! invariance the way the DCT hash does. See `mirror_descriptor` for the +//! exact symmetry. +//! +//! Not wired into the Descriptor trait or the stores. The salientlab binary +//! runs the experiments. + +use image::DynamicImage; +use std::f32::consts::PI; + +const ORI_BINS: usize = 36; + +/// Descriptor sampling grid: DESC_SAMPLES x DESC_SAMPLES points spanning +/// +-DESC_HALF_EXTENT keypoint sigmas, grouped into DESC_GRID x DESC_GRID +/// subregions. These define the exact window the descriptor reads, so the +/// browser demo draws the same box by calling Keypoint::window_corners. +const DESC_SAMPLES: usize = 20; +const DESC_HALF_EXTENT: f32 = 9.5; +pub const DESC_GRID: usize = 4; + +/// Grayscale f32 image, values in [0, 1] +pub struct GrayF32 { + w: usize, + h: usize, + data: Vec, +} + +impl GrayF32 { + pub fn from_image(img: &DynamicImage) -> GrayF32 { + let luma = img.to_luma32f(); + GrayF32 { + w: luma.width() as usize, + h: luma.height() as usize, + data: luma.into_raw(), + } + } + + fn get(&self, x: isize, y: isize) -> f32 { + let x = x.clamp(0, self.w as isize - 1) as usize; + let y = y.clamp(0, self.h as isize - 1) as usize; + self.data[y * self.w + x] + } + + fn bilinear(&self, x: f32, y: f32) -> f32 { + let x0 = x.floor(); + let y0 = y.floor(); + let fx = x - x0; + let fy = y - y0; + let (xi, yi) = (x0 as isize, y0 as isize); + let v00 = self.get(xi, yi); + let v10 = self.get(xi + 1, yi); + let v01 = self.get(xi, yi + 1); + let v11 = self.get(xi + 1, yi + 1); + v00 * (1.0 - fx) * (1.0 - fy) + + v10 * fx * (1.0 - fy) + + v01 * (1.0 - fx) * fy + + v11 * fx * fy + } + + /// Image gradient at a fractional position, central differences + fn gradient(&self, x: f32, y: f32) -> (f32, f32) { + ( + (self.bilinear(x + 1.0, y) - self.bilinear(x - 1.0, y)) * 0.5, + (self.bilinear(x, y + 1.0) - self.bilinear(x, y - 1.0)) * 0.5, + ) + } + + /// Separable gaussian blur, borders clamped + fn gauss_blur(&self, sigma: f32) -> GrayF32 { + let radius = (sigma * 3.0).ceil().max(1.0) as isize; + let mut kernel = Vec::with_capacity(2 * radius as usize + 1); + for i in -radius..=radius { + kernel.push((-(i * i) as f32 / (2.0 * sigma * sigma)).exp()); + } + let sum: f32 = kernel.iter().sum(); + for k in kernel.iter_mut() { + *k /= sum; + } + + let mut tmp = vec![0.0f32; self.w * self.h]; + for y in 0..self.h { + for x in 0..self.w { + let mut acc = 0.0; + for (ki, k) in kernel.iter().enumerate() { + acc += k * self.get(x as isize + ki as isize - radius, y as isize); + } + tmp[y * self.w + x] = acc; + } + } + let tmp = GrayF32 { w: self.w, h: self.h, data: tmp }; + let mut out = vec![0.0f32; self.w * self.h]; + for y in 0..self.h { + for x in 0..self.w { + let mut acc = 0.0; + for (ki, k) in kernel.iter().enumerate() { + acc += k * tmp.get(x as isize, y as isize + ki as isize - radius); + } + out[y * self.w + x] = acc; + } + } + GrayF32 { w: self.w, h: self.h, data: out } + } + + /// Every other pixel, for the next octave + fn half(&self) -> GrayF32 { + let w = (self.w / 2).max(1); + let h = (self.h / 2).max(1); + let mut data = Vec::with_capacity(w * h); + for y in 0..h { + for x in 0..w { + data.push(self.data[(y * 2) * self.w + x * 2]); + } + } + GrayF32 { w, h, data } + } + + fn sub(a: &GrayF32, b: &GrayF32) -> GrayF32 { + let data = a.data.iter().zip(&b.data).map(|(a, b)| a - b).collect(); + GrayF32 { w: a.w, h: a.h, data } + } +} + +#[derive(Clone, Copy, Debug)] +pub struct Keypoint { + /// Position in original image coordinates + pub x: f32, + pub y: f32, + /// Scale in original image coordinates + pub sigma: f32, + /// Dominant gradient orientation, radians + pub angle: f32, + /// Absolute DoG value at the extremum + pub response: f32, +} + +impl Keypoint { + /// Map a point from the keypoint's oriented frame (u along the + /// orientation, v perpendicular) into image coordinates. Same rotation + /// the descriptor uses to place its samples. + fn to_image(&self, u: f32, v: f32) -> (f32, f32) { + let (sin, cos) = self.angle.sin_cos(); + (self.x + u * cos - v * sin, self.y + u * sin + v * cos) + } + + /// The four corners of the oriented descriptor window in image + /// coordinates, from the (-u,-v) corner counter-clockwise. For drawing + /// the region the descriptor summarizes. + pub fn window_corners(&self) -> [(f32, f32); 4] { + let e = DESC_HALF_EXTENT * self.sigma; + [(-e, -e), (e, -e), (e, e), (-e, e)].map(|(u, v)| self.to_image(u, v)) + } + + /// Centres of the DESC_GRID x DESC_GRID subregions in image coordinates, + /// in descriptor order (row = perpendicular to orientation, column = + /// along it), so a caller can line them up with the 64-float descriptor. + pub fn subregion_centers(&self) -> [(f32, f32); DESC_GRID * DESC_GRID] { + let cell = DESC_SAMPLES / DESC_GRID; + let mut out = [(0.0, 0.0); DESC_GRID * DESC_GRID]; + for j in 0..DESC_GRID { + for i in 0..DESC_GRID { + let su = (i * cell + cell / 2) as f32 - DESC_HALF_EXTENT; + let sv = (j * cell + cell / 2) as f32 - DESC_HALF_EXTENT; + out[j * DESC_GRID + i] = self.to_image(su * self.sigma, sv * self.sigma); + } + } + out + } +} + +#[derive(Clone)] +pub struct Feature { + pub kp: Keypoint, + pub desc: [f32; 64], +} + +pub struct Params { + /// DoG scale samples per octave + pub intervals: usize, + /// Blur of the first level of each octave + pub base_sigma: f32, + /// Minimum absolute DoG value for a keypoint (image values are 0..1) + pub contrast_threshold: f32, + /// Maximum principal curvature ratio, rejects points on straight edges + pub edge_ratio: f32, + /// Keep only the strongest keypoints + pub max_features: usize, + /// Canonicalize descriptors under horizontal mirroring + pub flip_canonical: bool, +} + +impl Default for Params { + fn default() -> Self { + Params { + intervals: 3, + base_sigma: 1.6, + contrast_threshold: 0.01, + edge_ratio: 10.0, + max_features: 100, + flip_canonical: true, + } + } +} + +struct RawKp { + octave: usize, + level: usize, + x: usize, + y: usize, + bx: f32, + by: f32, + bsigma: f32, + response: f32, +} + +pub struct Salient { + pub params: Params, +} + +impl Default for Salient { + fn default() -> Self { + Self::new() + } +} + +impl Salient { + pub fn new() -> Salient { + Salient { params: Params::default() } + } + + fn k(&self) -> f32 { + 2f32.powf(1.0 / self.params.intervals as f32) + } + + /// Gaussian pyramid: per octave `intervals + 3` blur levels, so the DoG + /// stack has `intervals` usable middle layers + fn pyramid(&self, img: &DynamicImage) -> Vec> { + let s = self.params.intervals; + let k = self.k(); + let sigma0 = self.params.base_sigma; + // assume the input carries sigma 0.5, lift it to base_sigma + let gray = GrayF32::from_image(img); + let mut base = gray.gauss_blur((sigma0 * sigma0 - 0.25).max(0.01).sqrt()); + + let mut octaves = Vec::new(); + while base.w.min(base.h) >= 24 && octaves.len() < 6 { + let mut levels = vec![base]; + for i in 1..s + 3 { + let sig_prev = sigma0 * k.powi(i as i32 - 1); + let sig_total = sigma0 * k.powi(i as i32); + let sig_diff = (sig_total * sig_total - sig_prev * sig_prev).sqrt(); + levels.push(levels[i - 1].gauss_blur(sig_diff)); + } + // levels[s] carries 2 * base_sigma: halved it seeds the next octave + base = levels[s].half(); + octaves.push(levels); + } + octaves + } + + fn detect(&self, pyramid: &[Vec]) -> Vec { + let k = self.k(); + let mut raw = Vec::new(); + for (o, levels) in pyramid.iter().enumerate() { + let dogs: Vec = levels + .windows(2) + .map(|w| GrayF32::sub(&w[1], &w[0])) + .collect(); + let scale_up = (1usize << o) as f32; + for l in 1..dogs.len() - 1 { + let d = &dogs[l]; + for y in 1..d.h - 1 { + for x in 1..d.w - 1 { + let v = d.data[y * d.w + x]; + if v.abs() < self.params.contrast_threshold { + continue; + } + if !is_extremum(&dogs, l, x, y, v) { + continue; + } + if is_edge(d, x, y, self.params.edge_ratio) { + continue; + } + let osigma = self.params.base_sigma * k.powi(l as i32); + raw.push(RawKp { + octave: o, + level: l, + x, + y, + bx: x as f32 * scale_up, + by: y as f32 * scale_up, + bsigma: osigma * scale_up, + response: v.abs(), + }); + } + } + } + } + raw + } + + /// Detect, orient and describe the strongest keypoints of an image + pub fn features(&self, img: &DynamicImage) -> Vec { + let pyramid = self.pyramid(img); + let mut raw = self.detect(&pyramid); + raw.sort_by(|a, b| b.response.total_cmp(&a.response)); + + // greedy dedupe: extrema often fire in adjacent scale layers + let mut picked: Vec = Vec::new(); + for kp in raw { + if picked.len() >= self.params.max_features { + break; + } + let dup = picked.iter().any(|p| { + let (dx, dy) = (p.bx - kp.bx, p.by - kp.by); + let ratio = (p.bsigma / kp.bsigma).max(kp.bsigma / p.bsigma); + dx * dx + dy * dy < 4.0 && ratio < 1.6 + }); + if !dup { + picked.push(kp); + } + } + + picked + .iter() + .map(|kp| { + let img = &pyramid[kp.octave][kp.level]; + let osigma = self.params.base_sigma * self.k().powi(kp.level as i32); + let angle = orientation(img, kp.x, kp.y, osigma); + let mut desc = describe(img, kp.x as f32, kp.y as f32, osigma, angle); + if self.params.flip_canonical { + desc = canonicalize(desc); + } + Feature { + kp: Keypoint { + x: kp.bx, + y: kp.by, + sigma: kp.bsigma, + angle, + response: kp.response, + }, + desc, + } + }) + .collect() + } +} + +/// Strict local extremum over the 26 neighbors in space and scale +fn is_extremum(dogs: &[GrayF32], l: usize, x: usize, y: usize, v: f32) -> bool { + for layer in &dogs[l - 1..=l + 1] { + for dy in -1isize..=1 { + for dx in -1isize..=1 { + let n = layer.get(x as isize + dx, y as isize + dy); + if std::ptr::eq(layer, &dogs[l]) && dx == 0 && dy == 0 { + continue; + } + if (v > 0.0 && n >= v) || (v < 0.0 && n <= v) { + return false; + } + } + } + } + true +} + +/// Principal curvature ratio test on the 2x2 spatial Hessian of the DoG, +/// rejects points that sit on an edge and slide along it +fn is_edge(d: &GrayF32, x: usize, y: usize, r: f32) -> bool { + let (x, y) = (x as isize, y as isize); + let v = d.get(x, y); + let dxx = d.get(x + 1, y) - 2.0 * v + d.get(x - 1, y); + let dyy = d.get(x, y + 1) - 2.0 * v + d.get(x, y - 1); + let dxy = (d.get(x + 1, y + 1) - d.get(x + 1, y - 1) - d.get(x - 1, y + 1) + + d.get(x - 1, y - 1)) + / 4.0; + let tr = dxx + dyy; + let det = dxx * dyy - dxy * dxy; + det <= 0.0 || tr * tr * r >= (r + 1.0) * (r + 1.0) * det +} + +/// Dominant gradient orientation: 36-bin histogram of gradient angles in a +/// gaussian-weighted neighborhood, smoothed, peak refined by a parabola +fn orientation(img: &GrayF32, x: usize, y: usize, sigma: f32) -> f32 { + let sig_w = 1.5 * sigma; + let radius = (3.0 * sig_w).round() as isize; + let mut hist = [0f32; ORI_BINS]; + for dy in -radius..=radius { + for dx in -radius..=radius { + let (px, py) = (x as isize + dx, y as isize + dy); + let gx = (img.get(px + 1, py) - img.get(px - 1, py)) * 0.5; + let gy = (img.get(px, py + 1) - img.get(px, py - 1)) * 0.5; + let mag = (gx * gx + gy * gy).sqrt(); + if mag == 0.0 { + continue; + } + let mut phi = gy.atan2(gx); + if phi < 0.0 { + phi += 2.0 * PI; + } + let w = (-((dx * dx + dy * dy) as f32) / (2.0 * sig_w * sig_w)).exp(); + let bin = ((phi / (2.0 * PI) * ORI_BINS as f32) as usize).min(ORI_BINS - 1); + hist[bin] += w * mag; + } + } + // two passes of circular [1 2 1] smoothing + for _ in 0..2 { + let orig = hist; + for b in 0..ORI_BINS { + let l = orig[(b + ORI_BINS - 1) % ORI_BINS]; + let r = orig[(b + 1) % ORI_BINS]; + hist[b] = (l + 2.0 * orig[b] + r) / 4.0; + } + } + let peak = (0..ORI_BINS) + .max_by(|a, b| hist[*a].total_cmp(&hist[*b])) + .unwrap(); + let l = hist[(peak + ORI_BINS - 1) % ORI_BINS]; + let c = hist[peak]; + let r = hist[(peak + 1) % ORI_BINS]; + let denom = l - 2.0 * c + r; + let delta = if denom.abs() > 1e-12 { 0.5 * (l - r) / denom } else { 0.0 }; + (peak as f32 + 0.5 + delta) * 2.0 * PI / ORI_BINS as f32 +} + +/// SURF-style descriptor: 20x20 samples spaced `sigma` apart, rotated to the +/// keypoint orientation, grouped into 4x4 subregions. Each subregion gets +/// [sum du, sum dv, sum |du|, sum |dv|] where (du, dv) is the gradient +/// rotated into the keypoint frame. Unit normalized. +/// +/// Layout: index = (j * 4 + i) * 4 + c, with i along the orientation axis, +/// j perpendicular to it, c the component. +fn describe(img: &GrayF32, x: f32, y: f32, sigma: f32, angle: f32) -> [f32; 64] { + let (sin, cos) = angle.sin_cos(); + let sig_w = 3.3 * sigma; + let cell = DESC_SAMPLES / DESC_GRID; + let mut desc = [0f32; 64]; + for sj in 0..DESC_SAMPLES { + for si in 0..DESC_SAMPLES { + let u = (si as f32 - DESC_HALF_EXTENT) * sigma; + let v = (sj as f32 - DESC_HALF_EXTENT) * sigma; + let px = x + u * cos - v * sin; + let py = y + u * sin + v * cos; + let (gx, gy) = img.gradient(px, py); + let du = cos * gx + sin * gy; + let dv = -sin * gx + cos * gy; + let w = (-(u * u + v * v) / (2.0 * sig_w * sig_w)).exp(); + let base = ((sj / cell) * DESC_GRID + si / cell) * 4; + desc[base] += w * du; + desc[base + 1] += w * dv; + desc[base + 2] += w * du.abs(); + desc[base + 3] += w * dv.abs(); + } + } + let norm: f32 = desc.iter().map(|v| v * v).sum::().sqrt(); + if norm > 0.0 { + for v in desc.iter_mut() { + *v /= norm; + } + } + desc +} + +/// The descriptor a horizontally mirrored copy of the patch would produce. +/// +/// Under a mirror the dominant orientation maps to pi - theta. Working out +/// the sample grid in that frame: the u axis is preserved, the v axis +/// negates. So subregion (i, j) swaps with (i, 3 - j), the dv sum changes +/// sign, and the other three components are untouched. +pub fn mirror_descriptor(d: &[f32; 64]) -> [f32; 64] { + let mut out = [0f32; 64]; + for j in 0..4 { + for i in 0..4 { + let src = (j * 4 + i) * 4; + let dst = ((3 - j) * 4 + i) * 4; + out[dst] = d[src]; + out[dst + 1] = -d[src + 1]; + out[dst + 2] = d[src + 2]; + out[dst + 3] = d[src + 3]; + } + } + out +} + +/// Pick the mirror representative deterministically, so a patch and its +/// mirror image canonicalize to the same descriptor +fn canonicalize(d: [f32; 64]) -> [f32; 64] { + let dv_total: f32 = (0..16).map(|s| d[s * 4 + 1]).sum(); + let mirror = if dv_total.abs() > 1e-6 { + dv_total < 0.0 + } else { + // fallback: skew of the |dv| mass along the v axis + let skew: f32 = (0..16).map(|s| ((s / 4) as f32 - 1.5) * d[s * 4 + 3]).sum(); + skew < 0.0 + }; + if mirror { + mirror_descriptor(&d) + } else { + d + } +} + +/// 64 bits from one descriptor: the 32 signed sums contribute their sign, +/// the 32 magnitude sums contribute a comparison against their own median. +/// Same recipe as the DCT hash: signs plus ordinal facts, nothing absolute. +pub fn binarize64(desc: &[f32; 64]) -> u64 { + let mut mags: Vec = (0..16) + .flat_map(|s| [desc[s * 4 + 2], desc[s * 4 + 3]]) + .collect(); + mags.sort_by(f32::total_cmp); + let median = mags[16]; + let mut mask = 0u64; + for (i, &v) in desc.iter().enumerate() { + mask <<= 1; + let bit = if i % 4 < 2 { v > 0.0 } else { v > median }; + if bit { + mask |= 1; + } + } + mask +} + +/// Response-weighted mean of all descriptors, unit normalized. +/// Symmetric pooling, so it inherits every per-keypoint invariance. +pub fn pooled_vector(feats: &[Feature]) -> [f32; 64] { + let mut acc = [0f32; 64]; + for f in feats { + for (a, d) in acc.iter_mut().zip(&f.desc) { + *a += f.kp.response * d; + } + } + let norm: f32 = acc.iter().map(|v| v * v).sum::().sqrt(); + if norm > 0.0 { + for v in acc.iter_mut() { + *v /= norm; + } + } + acc +} + +/// Variant 1: binarized pooled descriptor +pub fn hash_pooled(feats: &[Feature]) -> u64 { + binarize64(&pooled_vector(feats)) +} + +/// Variant 2: binarized descriptor of the single strongest keypoint +pub fn hash_strongest(feats: &[Feature]) -> u64 { + match feats.first() { + Some(f) => binarize64(&f.desc), + None => 0, + } +} + +/// Variant 3: per-bit majority vote over all binarized descriptors +pub fn hash_consensus(feats: &[Feature]) -> u64 { + if feats.is_empty() { + return 0; + } + let hashes: Vec = feats.iter().map(|f| binarize64(&f.desc)).collect(); + let mut mask = 0u64; + for bit in (0..64).rev() { + let votes = hashes.iter().filter(|h| (*h >> bit) & 1 == 1).count(); + mask <<= 1; + if votes * 2 > hashes.len() { + mask |= 1; + } + } + mask +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::descriptors::testimg; + + fn dist(a: u64, b: u64) -> u32 { + (a ^ b).count_ones() + } + + #[test] + fn mirror_is_an_involution() { + let mut d = [0f32; 64]; + for (i, v) in d.iter_mut().enumerate() { + *v = (i as f32 * 0.37).sin(); + } + assert_eq!(mirror_descriptor(&mirror_descriptor(&d)), d); + } + + #[test] + fn canonical_ignores_mirroring() { + let mut d = [0f32; 64]; + for (i, v) in d.iter_mut().enumerate() { + *v = (i as f32 * 0.37).sin(); + } + assert_eq!(canonicalize(d), canonicalize(mirror_descriptor(&d))); + } + + #[test] + fn binarize_packs_all_bits() { + // strongest positive signal on the first and last positions + let mut d = [0f32; 64]; + d[0] = 1.0; // first sign bit, lands on bit 63 + d[63] = 1.0; // a magnitude above the (zero) median, lands on bit 0 + let h = binarize64(&d); + assert_eq!(h, (1 << 63) | 1); + } + + #[test] + fn features_deterministic() { + // synthetic gradients are too smooth to produce DoG extrema + let img = image::open("img/meowl.jpg").unwrap(); + let s = Salient::new(); + let a = s.features(&img); + let b = s.features(&img); + assert_eq!(a.len(), b.len()); + assert!(!a.is_empty()); + for (fa, fb) in a.iter().zip(&b) { + assert_eq!(fa.desc, fb.desc); + } + } + + #[test] + fn flip_invariant_on_real_photo() { + let img = image::open("img/meowl.jpg").unwrap(); + let s = Salient::new(); + let a = s.features(&img); + let b = s.features(&img.fliph()); + assert!(!a.is_empty()); + // detection mirrors exactly, orientation binning wobbles a little + let d = dist(hash_pooled(&a), hash_pooled(&b)); + assert!(d <= 6, "pooled flip distance {d}"); + let d = dist(hash_consensus(&a), hash_consensus(&b)); + assert!(d <= 6, "consensus flip distance {d}"); + } +} diff --git a/src/store.rs b/src/store.rs index 7b53314..9106369 100644 --- a/src/store.rs +++ b/src/store.rs @@ -262,6 +262,11 @@ impl DescriptorStore { self.names.len() } + /// All distinct hashes in the store + pub fn keys(&self) -> impl Iterator + '_ { + self.map.keys().copied() + } + pub fn is_empty(&self) -> bool { self.names.is_empty() } diff --git a/src/wasm/mod.rs b/src/wasm/mod.rs index f74fa52..acf7d5f 100644 --- a/src/wasm/mod.rs +++ b/src/wasm/mod.rs @@ -1,7 +1,9 @@ //! Wasm bindings for the browser demos. //! Hashes come back as BigInt, images go in and out as encoded bytes. +use crate::bk::BkTree; use crate::descriptors::{Descriptor, Median, DCT, PHash}; use crate::mutators::{self, Mutator}; +use crate::salient::{self, Params, Salient}; use wasm_bindgen::prelude::*; use std::io::Cursor; @@ -71,6 +73,227 @@ pub fn phash_lowfreq(bytes: &[u8]) -> Result, JsError> { Ok(PHash.lowfreq(&img).to_vec()) } +/// DCT hash of every mutant in the standard 16-mutator suite. +/// The image is shrunk first so the heavy mutators stay quick, which can +/// wobble the hashes a bit or two compared to a full-size run. +#[wasm_bindgen] +pub fn suite_hashes(bytes: &[u8]) -> Result, JsError> { + let img = load(bytes)?; + let img = if img.width().max(img.height()) > 512 { + img.thumbnail(512, 512) + } else { + img + }; + let dct = DCT::new(); + Ok(mutators::get_all_mutators() + .iter() + .map(|m| dct.describe(&m.mutate(&img))) + .collect()) +} + +/// Labels of the suite mutators, newline separated, same order as suite_hashes +#[wasm_bindgen] +pub fn suite_labels() -> String { + mutators::get_all_mutators() + .iter() + .map(|m| m.info()) + .collect::>() + .join("\n") +} + +/// A queryable set of hashes: the BK-tree plus a flat copy of the keys +/// as a linear scan baseline. +#[wasm_bindgen] +pub struct BkIndex { + tree: BkTree, + keys: Vec, +} + +#[wasm_bindgen] +impl BkIndex { + #[wasm_bindgen(constructor)] + pub fn new() -> BkIndex { + BkIndex { tree: BkTree::new(), keys: Vec::new() } + } + + /// Inserts one key, ignoring duplicates. Fine for demo-sized trees, + /// use insert_bytes for bulk loads. + pub fn insert(&mut self, key: u64) { + if !self.keys.contains(&key) { + self.keys.push(key); + self.tree.insert(key); + } + } + + /// Bulk load of little-endian u64 keys (the keys.bin format), + /// assumed to be distinct. + pub fn insert_bytes(&mut self, bytes: &[u8]) { + for chunk in bytes.chunks_exact(8) { + let key = u64::from_le_bytes(chunk.try_into().unwrap()); + self.keys.push(key); + self.tree.insert(key); + } + } + + pub fn len(&self) -> u32 { + self.keys.len() as u32 + } + + /// All keys within the radius, as flattened (key, distance) pairs + pub fn find(&self, key: u64, radius: u32) -> Vec { + self.tree + .find(key, radius as u64) + .into_iter() + .flat_map(|(k, d)| [k, d]) + .collect() + } + + /// Number of hashes the tree compared against during a find + pub fn find_compared(&self, key: u64, radius: u32) -> u32 { + self.tree.find_trace(key, radius as u64).1.len() as u32 + } + + /// The keys compared against during a find, in visit order + pub fn trace(&self, key: u64, radius: u32) -> Vec { + self.tree.find_trace(key, radius as u64).1 + } + + /// Linear scan baseline, same result as find + pub fn scan(&self, key: u64, radius: u32) -> Vec { + let radius = radius as u64; + self.keys + .iter() + .map(|&k| (k, ((k ^ key).count_ones()) as u64)) + .filter(|&(_, d)| d <= radius) + .flat_map(|(k, d)| [k, d]) + .collect() + } + + /// The tree as flattened (key, parent key, distance to parent) triples + /// in preorder; the root lists itself as parent. + pub fn structure(&self) -> Vec { + self.tree + .walk() + .into_iter() + .flat_map(|(k, p, d)| [k, p, d]) + .collect() + } +} + +impl Default for BkIndex { + fn default() -> Self { + Self::new() + } +} + +/// Detected salient points of an image plus their descriptors, for the +/// part-two demos. Built once, then queried for the overlay, the ranking +/// and the per-point neighbourhood description. Coordinates are in the +/// decoded image's own pixels (see width/height), scale them to the canvas. +#[wasm_bindgen] +pub struct SalientImage { + feats: Vec, + width: u32, + height: u32, +} + +#[wasm_bindgen] +impl SalientImage { + /// Detect and describe up to `max_features` points, strongest first. + #[wasm_bindgen(constructor)] + pub fn new(bytes: &[u8], max_features: u32) -> Result { + let img = load(bytes)?; + let detector = Salient { + params: Params { + max_features: max_features.max(1) as usize, + ..Default::default() + }, + }; + let feats = detector.features(&img); + Ok(SalientImage { feats, width: img.width(), height: img.height() }) + } + + /// Number of keypoints kept. + pub fn len(&self) -> u32 { + self.feats.len() as u32 + } + + pub fn is_empty(&self) -> bool { + self.feats.is_empty() + } + + /// Width of the coordinate space the keypoints live in. + pub fn width(&self) -> u32 { + self.width + } + + /// Height of the coordinate space the keypoints live in. + pub fn height(&self) -> u32 { + self.height + } + + /// Every keypoint as [x, y, sigma, angle, response], strongest first. + /// The response is the ranking score (how blob-like and high-contrast + /// the point is); it is what "importance" means here. + pub fn keypoints(&self) -> Vec { + self.feats + .iter() + .flat_map(|f| [f.kp.x, f.kp.y, f.kp.sigma, f.kp.angle, f.kp.response]) + .collect() + } + + /// The 64-float descriptor of keypoint `i`: 16 subregions, each + /// [sum du, sum dv, sum |du|, sum |dv|] in the point's oriented frame. + pub fn descriptor(&self, i: u32) -> Vec { + self.feats + .get(i as usize) + .map(|f| f.desc.to_vec()) + .unwrap_or_default() + } + + /// The four corners of keypoint `i`'s oriented sampling window in image + /// coordinates, as [x, y] x 4. Draw this to show the described region. + pub fn window(&self, i: u32) -> Vec { + self.feats + .get(i as usize) + .map(|f| { + f.kp.window_corners() + .iter() + .flat_map(|&(x, y)| [x, y]) + .collect() + }) + .unwrap_or_default() + } + + /// The 16 subregion centres of keypoint `i` in image coordinates, as + /// [x, y] x 16, in the same order as the descriptor's 16 blocks. + pub fn subregions(&self, i: u32) -> Vec { + self.feats + .get(i as usize) + .map(|f| { + f.kp.subregion_centers() + .iter() + .flat_map(|&(x, y)| [x, y]) + .collect() + }) + .unwrap_or_default() + } + + /// The 64-bit strongest-keypoint hash (the variant that survived at 25k). + pub fn hash(&self) -> u64 { + salient::hash_strongest(&self.feats) + } + + /// The 64-bit hash of an arbitrary keypoint `i`, so a demo can show a + /// specific point's fingerprint rather than only the strongest one. + pub fn descriptor_hash(&self, i: u32) -> u64 { + self.feats + .get(i as usize) + .map(|f| salient::binarize64(&f.desc)) + .unwrap_or(0) + } +} + /// Apply a single mutator to an image, returns PNG bytes. /// The meaning of `amount` depends on the kind. #[wasm_bindgen] diff --git a/web/index.html b/web/index.html index 99786f6..d93f1eb 100644 --- a/web/index.html +++ b/web/index.html @@ -1,200 +1,1308 @@ - - - - - - - - - - - -
-
-

Fucking around with perceptual hashes

-

Introduction

-

TODO: why this project exists, the master's thesis that wasn't, discovering pHash after the fact.

-

TODO: the rules of the game: 64 bits per image, hamming distance, nothing else.

-
+ + + + + + Fucking around with perceptual hashes + + + + + + + +
+
+
+

Fucking around with perceptual hashes

+

+ Quite some time ago I started a research project: Near-copy + detection of images using 64-bit image descriptors. For + various reasons, some personal and some technical, this + research project never saw the light of day. One of the + technical reasons is that I discovered halfway in that my + best (and only) original idea mostly already existed under + the name pHash + [5] +

+

+ I finally decided to get all of this out of my system and + write it up as a little blog post. So here it is: Perceptual + hashes, what is it, how they work, and what I fucked around + with. +

+ -
-

Perceptual hashing

-

TODO: cryptographic hash vs perceptual hash. One flipped pixel: md5 avalanches, a perceptual hash shrugs.

-

Terminology

-

TODO: descriptor / fingerprint / perceptual hash / signature all mean roughly the same thing depending on which corner of the literature you are standing in. Copy detection vs near-duplicate detection vs LSH.

-

Subjective image similarity

-

TODO: what "the same image" even means. Same pixels? Same scene? Same vibe?

-
+

Isn't that a bit small?

+

+ Putting an image into a 64-bit string is quite the + challenge. A typical image can vary from a couple hundred + KiB to several MiB. Putting all that in a 64-bit value is + like compressing an image to ~0.0005% of its original size. +

+

+ Even disregarding the challenges of trying to capture that + much information in that few bits, you very quickly run into + the + pigeonhole principle. We are trying to put an infinite amount of pigeons into + 264 pigeonholes, which is going to involve + creative bookkeeping, a lot of trimmed feathers, and + possibly the violation of some animal welfare laws. +

+

+ Our saving grace is that we don't really need to reconstruct + the original image from our 64-bit string representation. + What we want is some representation of the image that is + going to be resistant to small changes in the + input. We want to be able to say that two images are the + same (or very similar) when the Hamming distance + d between their representations is + small: +

+ + + d(a,b) + = + popcount(ab) + + +

+ I.e. how many bits are different between a and b? The + advantage of this simple distance function is that one + compare is a one xor and one + popcnt on x86. Depending on your CPU that means + you can easily compare upwards of a million of these pairs + per second in a single thread. +

+ +

+ Why limit ourselves when currently there are datacenters + floating around the earth that are processing terabytes of + data per second to generate a picture of a kitten falling + over? +

+

+ Part of it is that competing with + all those big scary algorithms + is too intimidating for me, but also working around such a + constraint is fun for me, as well as that I like my software + to be efficient. More seriously, every machine language, + database, and OS has their own corresponding primitive + 64-bit value (usually some type of integer). This makes a + 64-bit hash a natural fit for almost any system or algorithm + that might want to use it. +

+

+ Comparing two hashes costs a couple of instructions. A + million images index into 8 MB of RAM, and the whole thing + can run on your phone, an old netbook or anything else you + can find. +

+

+ The hash is even small enough to use as a filename. 64 + characters in a binary string converts to 16 characters hex + or 11 characters in + base58, + while keeping compatibility with most filesystems. This + works because + 264=1616 + and + 5811>264>5810. +

+

+ That covers comparing two images. Searching a whole + collection with one query is its own problem, with its own + clever data structure, and gets its own chapter near the + bottom, once we have built a hash worth searching with. +

+
-
-

Demo

-

Select an image to run this demo with. Don't worry, nothing will be sent to any server! All calculations are done in the browser.

-

Unanimated images with no transparancy work best, but feel free to experiment.

- - or pick a sample: - -
-
+
+

The demo image

+

+ Every demo on this page runs on one image, live in your + browser. A sample is preselected. Swap in another sample or + one of your own whenever you like. Selecting a new image + will recalculate all examples in this demo with the new + image. +

+

+ This entire demo runs in your browser. It is the original + code compiled to WebAssembly, no server involved, and I + never get your image. +

+

+

+ +

+

+ Most samples come from the MIRFLICKR-25000 collection[2], the same set used for the experiments at the bottom. The + first (and default) image is one I took myself of my cat + Spook, who was the best boy ever. +

+
+
+
+
+ +
+

Why a normal hash gets you nowhere

+

+ The naive approach is to use any standard hash function that + you would normally use to validate file equality. One of the + most well-known examples is MD5, which is a + cryptographic hash, which is a fancy way of saying + it is easy to calculate but hard to reverse*. +

+

+ Cryptographic hashes find byte-identical copies of files and + nothing else. Hashes built for integrity checking are + engineered for the avalanche effect: flip one input + bit and every output bit flips with probability one half. + Avalanche is exactly what you want when verifying a download + or storing a secret, and exactly what you do not want when + looking for pictures. +

+

+ To demonstrate this, the above image was altered slightly: + we did +1 to the red channel of the single pixel at the + center. A change so small that on a typical image without a + magnifying glass you cannot see it. MD5 outputs a 128 + "digest". Here are the first 64 bits of the original and + altered images' MD5 digest (as 16 characters of hex): +

+
+
+

MD5, original

+
+ +
+
+

MD5, one pixel changed

+
+ +
+
+

Changed bits are marked red.

+

+

+ For every change, the changed output bit is an + independent coin flip, so we expect around 32 + changed bits from the original. If you got a result + that significantly differs from 32 changed bits, you + might have gotten (un)lucky! +

+
+
+

+ Any operation on the original image, such as a + recompression, adjustment of the metadata, or the changing + of a single pixel, will completely change the hash value. +

+ +
+ +
+

Locality sensitive hashing

+

+ Hashes that don't display this avalanche behavior are called + locality-sensitive hashes (LSH). It is still a hash + function, i.e. it maps some arbitrary input domain to a + fixed-size output domain, but by some measure of similarity + the function clusters similar input. +

+

+ Unlike a cryptographic hash, LSH's are not designed to be + irreversible (though they may be). +

+

+ We will introduce the perceptual hashing algorithm + that I developed in the next steps, but as a preview + example, here is that perceptual hash applied to three + images. The first image is the original, the second is the + same image "hue rotated", and the third is an unrelated + photo. +

+
+
+

Your image

+
+ the selected image +
+
+
+
+

Hue shifted 90°

+
+ the selected image with shifted colors +
+
+

+ distance +

+
+
+

Unrelated photo

+
+ an unrelated photo +
+
+

+ distance +

+
+
+

+ Note that the second image has + every pixel different from the original, something + the md5 would definitely see as a change. But our perceptual + hash treats it as almost the same image. Depending on your + input image you should expect an exact match or at least a + very small distance. +

+

+ For the unrelated image we can see the same type of + difference as we were getting before, depending on the + visual similarity to the input image we would expect at + least upwards of 10 bits difference, but more likely about + half of the total 64 bits again. +

+

+ The thing to note here is that distance between hashes + suddenly carries meaning. +

+

+ So how do you build a perceptual hash? Every bit is a yes/no + question about the image. The art is picking questions whose + answers survive "mutation". We will discuss mutation in a + later chapter, but it can be any operation on the image that + preserves the perceptual similarity: recompression, + rotation, mirroring, resizing, color grading, contrast + adjustment, and more. +

+
+ +
+

Step 1: throw almost everything away

+

Scale down to 8x8. Fuck aspect ratio.

+

+ As it turns out, how you convert to grayscale matters. Luma + weights: +

+ + + Y= 0.2126R+ + 0.7152G+ 0.0722B + + +

+ 64 pixels left, which is also our bit-budget. How can that + now? What a wild coinkydink. +

+

+ Your original image, compared to what the hash gets to work + with: +

+
+
+ +
+

Step 2: Median threshold

+

+ Take the median of the 64 pixels and emit one bit per pixel: +

+ + + bi + = + [ + pi + > + median(p) + ] + + +

+ 1 when the condition holds, 0 otherwise. Median from Thomee + et al.[1] + + Does a couple of things well. Brightness, gamma, and color + shifts. +

+
+
+

8×8 input

+
+
+
+

Bits

+
+
+
+

Hash

+ +
+
+

+ The bits are tied to pixel positions, so a mirrored copy + scrambles them completely. See the shitty flip recall below. + Also unrelated images might share this course layout. + Imagine that many landscape fotographs might have the same + general silhouette. +

+
+ +
+

Step 3: DCT

+

+ The discrete cosine transform rewrites the 8×8 + thumbnail as a weighted sum of 64 fixed cosine patterns. The + weights are the coefficients: +

+ + + C(u,v) + = + + + α(u)α(v) + + 4 + + x=07 + y=07 + (pxy128) + cos((2x+1)uπ16) + cos((2y+1)vπ16) + = + + +

+ TODO: Dit is waarschijnlijk beter als het meer simpel is. + Hoop van die termen doen er niet echt toe hier. +

+ +
+
+

The 64 patterns

+
+ +
+
+
+

8×8 input

+
+
+
+

Your coefficients

+
+
+ +
+

+ +

+
+
+

+ Top-left is average brightness. Horizontal frequency + increases to the right, vertical frequency increases + downward. +

+

+ One of the key observations of JPEG image compression is + that typically speaking most images have most of their + "energy" concentrated in the lower-frequency patterns. That + is why the zigzag ordering is so effective: it effectively + gives a MSB ordering of the coefficients. After that, JPEG + uses a quantization table to reduce the higher-frequency + components to less bits. For our case, we can just discard + the LSBs. +

+

Putting it back together

+

+ The proof is that the weighted patterns sum back to the + picture. Add them one at a time, in zigzag order: +

+
+
+
+
the 8×8
+
+
+ +
+ first 8 of 64 patterns +
+
+
+
+
+ +
+
+ coefficients added, in zigzag order +
+
+
+

+ + + the hash stops reading at 36 +

+

+ Note that the general shape of the picture comes quite + quickly with the lower frequency components. The higher + frequency components simply refine existing shapes and + edges. That is why cutting the zigzag at 36 loses so little: + the animation pauses there, and the second half barely + changes the picture. +

+
+ +
+

Step 4: 64 bits

+

+ Oké maar dat is dus nog steeds veel te veel data. 64 floats. + ja doei. Dus we doen signs: +

+ + + sk= [c~k<0] + + +

+ En de rest dan op volgorde: Is het volgende patroontje + sterker dan de huidige? +

+ + + ok= + [ + |ck+1| + > + |ck| + ] + + +

+ 36 + 28 = 64, qed. Dan nu het slimme. + + c~ + als normalisatie stap. We gaan altijd uit van dat het eerste + patroon niet geinverteerd is. Als ie dat wel is draaien we + hem alsnog om. Dus dan maakt het niet uit of de patroontjes + omgedraaid zijn of niet. +

+

Snippet:

+
let flip = dct[1].signum();              // sign of C(1,0)
+for &i in ZIGZAG.iter().take(36) {
+    let mut sign = dct[i].signum();
+    if i % 2 != 0 { sign *= flip }       // odd horizontal frequency
+    sign_mask = sign_mask << 1 | (sign < 0.0) as u64;
+}
+for pair in ZIGZAG[..=28].windows(2) {
+    let bigger = dct[pair[1]].abs() > dct[pair[0]].abs();
+    ordinal_mask = ordinal_mask << 1 | bigger as u64;
+}
+let hash = sign_mask << 28 | ordinal_mask;
+
+
+

Sign mask (36)

+
+
+
+

Ordinal mask (28)

+
+
+
+

Hash

+ +
+
+

+ Kim (2003)[3] + + is de Cosine baseline in Thomee et al.[1] +

+

+ Dit is zo mogelijk het enige echt unieke aan mijn oplossing. + whoop. +

+

What we remember

+

+ We can reconstruct the original 8x8 image from these values. + We store no values, only signs and cardinality. Still, we + can make a rough approximation of the 8x8 starting point: +

+
+
+
+
the 8×8
+
+
+ +
rebuilt from the 64 bits
+
+
+

+ A ghost, but a recognizable one, out of 8 bytes. One quirk: + because of the flip normalization the bits genuinely cannot + tell left from right, so the ghost sometimes comes out + mirrored. +

+

+ +

+
+ +
+

Step 5: how pHash does it

+

+ +

+
+
+

32×32

+
+
+
+

Low-frequency block

+
+
+
+

Bits

+
+
+
+

Hash

+ +
+
+ main difference is the mid frequencies that get saved. same + idea, different execution. makes it better on rotations and logo + insertions. improvements to dct possible. + +

+ pHash is the work of the + pHash.org project, + which also hosts an + online demo in the + same spirit as this page. The implementation here follows + the widely used + imagehash + recipe. +

+
+ +
+

Mutations

+

+ + + +

+
+ + + + + + + + + + + + + +
dctmedianphash
Hamming distance
+

+ Distances of 4 or less are green: [motivate this with PR + curves?] +

+
+ +
+

The flip trick

+ + + C(u,v) + + (1)u + C(u,v) + + +

+ Magnitudes never change, so ordinal bits are flip-invariant + for free. Only the signs of odd horizontal frequencies flip, + and they all flip together. So normalize them + against one of their own: +

+ + + c~(u,v) + = + C(u,v) + · + + sign(C(1,0)) + u + + + +

+ +

+
+ the selected image, flippable +
+ + + + + + + + + + + + + + + + + + + + + + + + + +
originalflippeddistance
dct +
+
+
+
median +
+
+
+
phash +
+
+
+
+

Bits that changed under the flip are red.

+
+ +
+

Find the copy

+

+ The pool below holds 1,111 images from MIRFLICKR-25000[2], hashed offline from the full-size originals. Click any + image to make it the query, or upload your own. Each method + then returns its ten nearest neighbours by Hamming distance. +

+

+ Mutate the query: + + or query your own image: + +

+
+
+
+

Query

+
+
+
+

dct

+
+
+
+

median

+
+
+
+

phash

+
+
+
+

+ no mutation should always yield d 0, as the algorithms are + deterministic. with mutation you can see that some + algorithms still find the original for some mutations. try + to figure out which algorithms are invariant to which + changes. +

+

+ Images: the + MIRFLICKR-25000 + collection (Huiskes and Lew, MIR '08)[2], Creative Commons photography collected from Flickr. +

+
+ +
+

PR CURVES

+
+
+ +
+
+
+
dct
+
+
+
+
phash
+
+
+
+
median
+
+
+
+ +

+ Rendered live from the 25,000-image run of 2026-07-04; line + style marks the mutation category. The dots mark the + slider's threshold. Data: pr-data.json, extracted from the + run log by pr_to_json.py. +

+
+ +
+

Retrieval

+

+ Comparing two is nice, but ranking is the real deal. Two + questions: When do we declare two hashes to be the same + image? And how do we collect that subset without comparing + the query against every hash we have? +

+ +

When are two images the same?

+

+ pick a threshold t and call everything + at distance t or less the same image. + your image against two groups: its own 16 mutated copies + from the experiment suite, and the 1,111 unrelated pool + images from the ranking demo. +

+
+
+
+ Hamming distance from your image to its 16 mutated + copies (green, hover for the mutation) and to the 1,111 + pool images (gray bars, square-root count scale). +
+
+

+ + +

+

+ The two groups keep a comfortable distance from each other. + Copies near zero, strangers pile up just under half of the 64 + bits*. + In between the kingdom of mutated hits and the people's republic of unrelated images is the eehm.. the.. Federated Islands of false positives. + When a mutation goes too far, the hash distance increases and the fingerprint gets banished from the kingdom and might end up in the republic of strangers. + + Set your threshold too high and you start catching some of these strangers. + + The tradeoff then is how serious you want to be about it: do you accept some false positives or would you rather get false negatives? + It all depends on what you might use these copy-detectors for. + If it is the first step in reducing an image set in order for the Big Guns to take over, you might be more willing to accept false positives. + If you want to quickly see if a given image is already in your folder of holiday pictures, you might not mind a false negative or two, but would rather not have a false positive. + + The threshold is what makes this tradeoff, and the PR-curve is where it shows. + +

+ + +

Skipping most of the work

+

+ The obvious retrieval algorithm is a linear scan: compare + the query against all N stored hashes, + keep everything within t. At a million comparisons per second per thread that is + genuinely fine for a while. But the cost grows with every + stored image and is paid again on every query, and "compare + against everything" should offend you a little when the + answer is almost always "no". +

+

+ The + Burkhard-Keller tree[6] fixes + this for any metric distance, and Hamming distance is one. + Pick any stored hash as the root. Every other hash goes into + a subtree based on its distance to the root: all hashes at + distance 7 from the root share subtree 7, and inside each + subtree the same rule repeats. To search, compare the query + to the root, giving some distance d. A match can only hide in a subtree whose label lies + between + dt + and + d+t, that is the triangle inequality doing its thing. Recurse + into the surviving subtrees, ignore the rest forever. +

+

+ + + +

+
+
+
+ A BK-tree over the 16 mutant hashes of your image plus + 16 pool strangers. Mutants with identical hashes + collapse into one node, so the tree is usually smaller + than 32. Edge labels are distances to the parent. Query: + your image's hash. Green: within the radius. Outlined: + compared, too far. Dimmed: pruned, never even looked at. +
+
+

+ At small t the search drops straight + into the branch where the copies cluster and skips most + stranger branches without computing a single distance in + them. Raise t and the band + d±t + widens, fewer branches get pruned, and the search slowly + degrades back into visiting everyone. +

+ +

Does it scale?

+

+ The toy tree has 32 hashes, the experiment above produced + 425,000 (25,000 images, 16 mutants each). Those collapse to + 239,689 distinct hashes, identical images simply share one + entry. Below, all of them sit in a BK-tree in your browser's + memory, and your image queries it live. +

+
+
+

In the tree

+

+

distinct hashes

+
+
+

Found

+

+

+
+
+

Compared

+

+

+
+
+

Compared vs radius

+
+
+
+

+ +

+

+ Query with one of the MIRFLICKR samples from the strip at + the top and its mutants come right back. Query with Spook + and nothing comes back, because my cat is not among the + 25,000. +

+

+ At + t=4 + the tree answers after touching a few percent of the hashes, + roughly a 30× saving in comparisons. Now look at the + stopwatch: the dumb scan is still competitive, and at wide + radii it wins outright. Marching sequentially through memory + doing one xor and one popcnt per + hash is about the kindest thing you can do to a CPU, while + the tree spends its savings on hopping through pointers. + The comparisons saved only turn into time saved when + the collection outgrows this demo by an order of + magnitude or two, since the scan grows linearly and the + visited slice of the tree does not. + Also note how quickly the pruning decays in the curve: + perceptual hashes cluster, so a wide radius keeps almost + every branch alive. A BK-tree only earns its keep at small + radii, which is conveniently the only place our threshold + wants to be. +

+ +

Do we expect it to be perfect?

+

+ No, and it does not have to be. The histogram showed copies + that drift out of reach and the occasional stranger inside + the radius, the precision-recall curves put numbers on both. + If a wrong answer is expensive, treat the whole thing as a + preselection filter: the hash plus BK-tree reduces 425,000 + candidates to a handful in microseconds, and whatever + heavyweight comparison you actually trust (full-resolution + diffing, feature matching, a neural embedding, a human) only + runs on that handful. A cheap filter in front of an + expensive judge is a classic setup, and a 64-bit hash is + about the cheapest filter there is. +

+
+ + + +
+

References

+
    +
  1. + B. Thomee, M. Huiskes, E. Bakker, M. Lew. + Large scale image copy detection evaluation. + MIR '08. The mutation taxonomy and the Median and Cosine + baselines come from here. +
  2. +
  3. + M. Huiskes, M. Lew. + The MIR Flickr retrieval evaluation. MIR '08. + press.liacs.nl/mirflickr +
  4. +
  5. + C. Kim. Content-based image copy detection. + Signal Processing: Image Communication, 2003. Ordinal + measures over DCT coefficients. +
  6. +
  7. + M. Charikar. + Similarity estimation techniques from rounding + algorithms. + STOC '02. The random hyperplane bound. +
  8. +
  9. + pHash.org and the + Python + imagehash + library. +
  10. +
  11. + W. Burkhard, R. Keller. + Some approaches to best-match file searching. + Communications of the ACM, 1973. The BK-tree. +
  12. +
+

+ Everything on this page is Rust compiled to WebAssembly plus + plain JS and CSS. Source: link to repo. +

+
-
- -
-

Step 1: Resize

-

- The first step is to size the image down, and remove all color information. - Converting the image to grayscale is done by simply averaging the pixels. - Since we are aiming for a 64-bit descriptor, the logical target size is 8 by 8 since that will give us 64 pixels to work with. -

-

- It is obvious, but worth noting, that this action is destructive. - We are throwing away a lot of information here, especially regarding the finer details. - For now, this is a good thing. -

-

- Here we use a nearest neighbour downscaling algorithm that does not preserve aspect ratio. - Press the button to resize your image:
- -

-
-
-
-
-
- -
-

Step 2: The simplest hash that could possibly work

-

TODO: median hash: one bit per pixel, brighter than the median or not. This is the "Median" baseline from Thomee et al.

-
-
-

Bits

-
-
-
-

Hash

- -
-
-

TODO: why this breaks: global median shifts, flips scramble everything. Foreshadow the moon/plate anecdote.

-
- -
-

Step 3: Frequency space

-

TODO: what the DCT does. JPEG uses the same trick. Low frequencies = global shape, high frequencies = detail we already threw away.

-
-
-

DCT coefficients

-
-
-
-

TODO: reading the coefficient grid: top-left is the average, first row is horizontal waves, first column vertical waves.

-

From coefficients to bits

-

TODO: sign bits of the first 36 zigzag coefficients + 28 ordinal comparisons between neighbours in zigzag order. The odd-column signs get multiplied by the sign of coefficient (1,0), which buys horizontal flip invariance for one bit.

-
-
-

Sign mask (36)

-
-
-
-

Ordinal mask (28)

-
-
-
-

Hash

- -
-
-
- -
-

Step 4: How pHash does it

-

TODO: same idea, different route: resize to 32x32 instead of 8x8, DCT, keep the top-left 8x8 block of low frequencies, threshold against the median coefficient.

-
-
-

32×32

-
-
-
-

Low-frequency block

-
-
-
-

Bits

-
-
-
-

Hash

- -
-
-

TODO: what's genuinely different between our dct hash and phash (bit extraction, flip invariance) and what isn't (everything else).

-
- -
-

Mutations

-

TODO: a copy is rarely byte-identical. Recoding, resampling, content processing, framing, inserted logos (the Thomee et al. taxonomy). Sweep the slider and watch which hash survives what.

-

- - - -

-
- original - mutated -
-

- -

- - - - - - - - -
dctmedianphash
hamming distance
-

TODO: what counts as "the same" now? Thresholds. pHash uses 22 of 64 bits, we will measure our own.

-
- -
-

Find the copy

-

TODO: the fun part: a tiny search engine. Add a pile of images, click one, get the nearest neighbours per algorithm. Mention the white plate that matched the moon.

-

- - -

-
-

Results

-
-
- -
-

At scale

-

TODO: the browser demo is anecdote, this section is data. 25k images, 16 mutations each, precision-recall over the hamming threshold.

-
- PR curves for the dct hash -
TODO: dct-pr.png from the experiment run
-
-
- PR curves for median and phash -
TODO: median-pr.png and phash-pr.png
-
-

TODO: how the PR curves are computed, what a false positive means here, where the thresholds land per algorithm.

-
- -
-

Loose ends

-

TODO: vertical flips and 180 rotations (same trick, one more bit). Coefficient stability near zero. Video: keyframes vs temporally averaged frames vs 3D-DCT. The seen-images daemon idea. Sorting a folder by visual similarity.

-
-
- - - + + diff --git a/web/script.js b/web/script.js index c836061..4da9950 100644 --- a/web/script.js +++ b/web/script.js @@ -24,8 +24,21 @@ worker.onmessage = (e) => { } }; +if (window.hljs) hljs.highlightAll(); + +// Background fade-in, same trick as the landing page +{ + const bg = document.querySelector('#bg img'); + if (bg.complete) bg.classList.add('loaded'); + else bg.addEventListener('load', () => bg.classList.add('loaded')); +} + // Helpers +function cssColor(name) { + return getComputedStyle(document.documentElement).getPropertyValue(name).trim(); +} + function pngUrl(bytes) { return URL.createObjectURL(new Blob([bytes], { type: 'image/png' })); } @@ -44,54 +57,242 @@ function hex64(hash) { return hash.toString(16).padStart(16, '0'); } -// Renders the top `bits` bits of a hash, msb first -function bitGrid(container, hash, bits = 64) { +// Renders the top `bits` bits of a hash, msb first. +// Bits set in `diff` get marked as changed. +function bitGrid(container, hash, bits = 64, diff = 0n) { container.replaceChildren(); for (let i = bits - 1; i >= 0; i--) { const cell = document.createElement('div'); cell.className = (hash >> BigInt(i)) & 1n ? 'bit on' : 'bit off'; + if ((diff >> BigInt(i)) & 1n) cell.classList.add('diff'); container.appendChild(cell); } } -// 8x8 grid of coefficient magnitudes, log scale, red negative / blue positive +// 8x8 grid of coefficient magnitudes, log scale. +// Green positive, red negative (sign is the semantic part), brightness is magnitude. function heatGrid(container, values) { container.replaceChildren(); const max = Math.max(...values.map(Math.abs), 1e-9); for (const value of values) { const cell = document.createElement('div'); const strength = Math.log1p(Math.abs(value)) / Math.log1p(max); - const hue = value < 0 ? 4 : 215; - cell.style.background = `hsl(${hue} 70% ${15 + strength * 45}%)`; + const hue = value < 0 ? 4 : 140; + cell.style.background = `hsl(${hue} 55% ${10 + strength * 48}%)`; cell.title = value.toFixed(1); container.appendChild(cell); } } -function show(section) { +function pixelatedImg(url, size = 176) { + const img = document.createElement('img'); + img.src = url; + img.className = 'pixelated'; + img.width = size; + img.height = size; + return img; +} + +function show() { document.querySelectorAll('.needs-image').forEach((el) => el.classList.add('visible')); - if (section) document.getElementById(section).scrollIntoView({ behavior: 'smooth' }); +} + +// MD5, for the avalanche demo. crypto.subtle dropped it long ago, so a plain +// implementation (verified against md5sum). Returns 16 bytes. +function md5(bytes) { + const K = new Uint32Array(64); + for (let i = 0; i < 64; i++) K[i] = Math.floor(Math.abs(Math.sin(i + 1)) * 2 ** 32); + const S = [7, 12, 17, 22, 5, 9, 14, 20, 4, 11, 16, 23, 6, 10, 15, 21]; + const n = bytes.length; + const total = (((n + 8) >> 6) + 1) << 6; + const buf = new Uint8Array(total); + buf.set(bytes); + buf[n] = 0x80; + const dv = new DataView(buf.buffer); + dv.setUint32(total - 8, (n * 8) >>> 0, true); + dv.setUint32(total - 4, Math.floor(n / 536870912), true); + let a0 = 0x67452301, b0 = 0xefcdab89, c0 = 0x98badcfe, d0 = 0x10325476; + const M = new Uint32Array(16); + for (let off = 0; off < total; off += 64) { + for (let j = 0; j < 16; j++) M[j] = dv.getUint32(off + j * 4, true); + let A = a0, B = b0, C = c0, D = d0; + for (let i = 0; i < 64; i++) { + let F, g; + if (i < 16) { F = (B & C) | (~B & D); g = i; } + else if (i < 32) { F = (D & B) | (~D & C); g = (5 * i + 1) % 16; } + else if (i < 48) { F = B ^ C ^ D; g = (3 * i + 5) % 16; } + else { F = C ^ (B | ~D); g = (7 * i) % 16; } + F = (F + A + K[i] + M[g]) | 0; + A = D; D = C; C = B; + const s = S[((i >> 4) << 2) | (i & 3)]; + B = (B + ((F << s) | (F >>> (32 - s)))) | 0; + } + a0 = (a0 + A) | 0; b0 = (b0 + B) | 0; c0 = (c0 + C) | 0; d0 = (d0 + D) | 0; + } + const out = new Uint8Array(16); + const odv = new DataView(out.buffer); + odv.setUint32(0, a0 >>> 0, true); + odv.setUint32(4, b0 >>> 0, true); + odv.setUint32(8, c0 >>> 0, true); + odv.setUint32(12, d0 >>> 0, true); + return out; +} + +// Comparison slider: two stacked images, a draggable vertical divider. + +function makeCompare(container, { pixelated = false } = {}) { + const a = document.createElement('img'); + const b = document.createElement('img'); + a.className = 'compare-a'; + b.className = 'compare-b'; + if (pixelated) b.classList.add('pixelated'); + const labelA = document.createElement('span'); + const labelB = document.createElement('span'); + labelA.className = 'compare-label a'; + labelB.className = 'compare-label b'; + labelA.innerText = container.dataset.labelA || ''; + labelB.innerText = container.dataset.labelB || ''; + const divider = document.createElement('div'); + divider.className = 'compare-divider'; + const handle = document.createElement('div'); + handle.className = 'compare-handle'; + handle.innerText = '↔'; + divider.appendChild(handle); + container.append(a, b, labelA, labelB, divider); + + let pos = 50; + let sweeping = null; + function apply() { + b.style.clipPath = `inset(0 0 0 ${pos}%)`; + divider.style.left = `${pos}%`; + } + apply(); + + function positionFromEvent(e) { + const rect = container.getBoundingClientRect(); + return Math.max(0, Math.min(100, ((e.clientX - rect.left) / rect.width) * 100)); + } + container.addEventListener('pointerdown', (e) => { + cancelAnimationFrame(sweeping); + container.setPointerCapture(e.pointerId); + pos = positionFromEvent(e); + apply(); + }); + container.addEventListener('pointermove', (e) => { + if (!container.hasPointerCapture(e.pointerId)) return; + pos = positionFromEvent(e); + apply(); + }); + + // One slow reveal sweep so it is obvious the thing moves + function sweep() { + cancelAnimationFrame(sweeping); + const start = performance.now(); + const duration = 1200; + function tick(now) { + const t = Math.min(1, (now - start) / duration); + const ease = 1 - Math.pow(1 - t, 3); + pos = 2 + ease * 48; + apply(); + if (t < 1) sweeping = requestAnimationFrame(tick); + } + sweeping = requestAnimationFrame(tick); + } + + // The reveal sweep only plays when the base image changes, not on + // every update of the overlaid result. + let lastA = null; + return { + set(aUrl, bUrl) { + a.src = aUrl; + b.src = bUrl; + if (aUrl !== lastA) sweep(); + lastA = aUrl; + }, + }; +} + +// DCT helpers shared by the reconstruction demos + +// Full 64-entry zigzag, matching the order in dct.rs (which stops at 36) +const FULL_ZIGZAG = (() => { + const order = []; + for (let s = 0; s < 15; s++) { + const diagonal = []; + for (let u = 0; u <= s; u++) { + const v = s - u; + if (u < 8 && v < 8) diagonal.push(v * 8 + u); + } + if (s % 2 === 1) diagonal.reverse(); + order.push(...diagonal); + } + return order; +})(); +const ZIGZAG = FULL_ZIGZAG.slice(0, 36); + +// Inverse of the transform in dct.rs: orthonormal scaling, +128 recentering +function idct(values) { + const out = new Array(64); + for (let k = 0; k < 64; k++) { + const x = k % 8; + const y = Math.floor(k / 8); + let sum = 0; + for (let u = 0; u < 8; u++) { + for (let v = 0; v < 8; v++) { + let alpha = 1; + if (u === 0) alpha /= Math.SQRT2; + if (v === 0) alpha /= Math.SQRT2; + sum += alpha * values[v * 8 + u] * + Math.cos(((2 * x + 1) * u * Math.PI) / 16) * + Math.cos(((2 * y + 1) * v * Math.PI) / 16); + } + } + out[k] = Math.max(0, Math.min(255, Math.round(128 + 0.25 * sum))); + } + return out; +} + +function drawPixels(canvas, pixels) { + const ctx = canvas.getContext('2d'); + const data = ctx.createImageData(8, 8); + for (let i = 0; i < 64; i++) { + data.data[i * 4] = pixels[i]; + data.data[i * 4 + 1] = pixels[i]; + data.data[i * 4 + 2] = pixels[i]; + data.data[i * 4 + 3] = 255; + } + ctx.putImageData(data, 0, 0); } // Demo state + let baseBuffer = null; // ArrayBuffer of the selected image let baseUrl = null; -let baseHashes = null; // BigUint64Array [dct, median, phash] +let baseHashes = null; // [dct, median, phash] as BigInt +let baseCoefficients = null; // 64 DCT coefficients of the 8x8 +let loadToken = 0; -// Pipeline demo +// Image picker and pipeline const imageContainers = document.getElementsByClassName('image-original'); const formImage = document.getElementById('dctimage'); +const resizeCompare = makeCompare(document.getElementById('resize-compare'), { pixelated: true }); +const mutateCompare = makeCompare(document.getElementById('mutate-compare')); function reset() { for (const container of imageContainers) { container.replaceChildren(); } - document.getElementById('image-resize').replaceChildren(); - document.getElementById('resize').classList.remove('resize'); + for (const name of ['dct', 'median', 'phash']) { + document.getElementById(`flip-b-${name}`).replaceChildren(); + document.getElementById(`flip-d-${name}`).innerText = ''; + } + document.getElementById('flip-image').classList.remove('flipped'); + flipHashes = null; } async function loadImage(buffer, url) { + const token = ++loadToken; reset(); baseBuffer = buffer; baseUrl = url; @@ -100,19 +301,21 @@ async function loadImage(buffer, url) { img.src = url; container.appendChild(img); } + document.getElementById('flip-image').src = url; const data = await rpc('pipeline', { buffer }); + if (token !== loadToken) return; baseHashes = data.hashes; + baseCoefficients = data.coefficients; const [dct, median, phash] = data.hashes; // resize step - const resized = document.createElement('img'); - resized.src = pngUrl(data.resize8); - document.getElementById('image-resize').replaceChildren(resized); - const original = document.getElementById('resize-original').children[0]; - if (original && original.width) { - original.width = original.width; - original.height = original.height; + const resize8Url = pngUrl(data.resize8); + resizeCompare.set(url, resize8Url); + + // small 8x8 reference thumbs sprinkled through the sections + for (const id of ['median-input', 'dct-input', 'recon-original', 'ghost-original']) { + document.getElementById(id).replaceChildren(pixelatedImg(resize8Url)); } // median @@ -124,6 +327,9 @@ async function loadImage(buffer, url) { bitGrid(document.getElementById('dct-sign-bits'), data.masks[0], 36); bitGrid(document.getElementById('dct-ordinal-bits'), data.masks[1], 28); document.getElementById('dct-hash').innerText = hex64(dct); + heatGrid(document.getElementById('recon-grid'), data.coefficients); + renderReconstruction(); + renderGhost(data.masks); // phash const small = document.createElement('img'); @@ -134,7 +340,14 @@ async function loadImage(buffer, url) { bitGrid(document.getElementById('phash-bits'), phash); document.getElementById('phash-hash').innerText = hex64(phash); + // flip demo, original column + const names = ['dct', 'median', 'phash']; + names.forEach((name, i) => bitGrid(document.getElementById(`flip-a-${name}`), data.hashes[i])); + show(); + runMd5Demo(buffer, token); + runLshDemo(buffer, token); + runRetrieval(buffer, token); resetMutation(); runMutation(); } @@ -145,27 +358,333 @@ formImage.addEventListener('change', () => { file.arrayBuffer().then((buffer) => loadImage(buffer, URL.createObjectURL(file))); }); -document.getElementById('resize-button').addEventListener('click', () => { - document.getElementById('resize').classList.add('resize'); -}); - -// Sample images, also used to seed the ranking pool -const samples = ['img/moon1.jpg', 'img/moon2.jpg', 'img/sunflower1.jpg', 'img/sunflower2.jpg']; +// Sample images. The first one doubles as the preselected default. +const samples = [ + 'img/spook.png', + 'img/im21050.jpg', + 'img/im21644.jpg', + 'img/im21068.jpg', + 'img/im21123.jpg', + 'img/im21694.jpg', + 'img/moon1.jpg', + 'img/moon2.jpg', + 'img/sunflower1.jpg', + 'img/sunflower2.jpg', +]; const sampleContainer = document.getElementById('sample-images'); for (const src of samples) { + // buttons are created up front so the strip keeps a deterministic order + const button = document.createElement('img'); + button.className = 'sample'; + button.title = src.split('/').pop(); + sampleContainer.appendChild(button); fetch(src) .then((response) => (response.ok ? response.arrayBuffer() : Promise.reject(response.status))) .then((buffer) => { - const button = document.createElement('img'); button.src = src; - button.className = 'sample'; button.addEventListener('click', () => loadImage(buffer, src)); - sampleContainer.appendChild(button); - addToPool(src.split('/').pop(), buffer, src); + if (src === samples[0] && baseBuffer === null) loadImage(buffer, src); }) - .catch(() => {}); + .catch(() => button.remove()); } +// MD5 avalanche demo: the same PNG twice, one pixel one unit apart. +// Also quietly hashes both files perceptually for the payoff in the mutation section. + +async function runMd5Demo(buffer, token) { + const bitmap = await createImageBitmap(new Blob([buffer])); + const scale = Math.min(1, 1024 / Math.max(bitmap.width, bitmap.height)); + const width = Math.max(1, Math.round(bitmap.width * scale)); + const height = Math.max(1, Math.round(bitmap.height * scale)); + const canvas = document.createElement('canvas'); + canvas.width = width; + canvas.height = height; + const ctx = canvas.getContext('2d', { willReadFrequently: true }); + ctx.drawImage(bitmap, 0, 0, width, height); + bitmap.close(); + + const toPng = () => + new Promise((resolve) => canvas.toBlob((blob) => resolve(blob.arrayBuffer()), 'image/png')); + + const pngA = await toPng(); + const data = ctx.getImageData(0, 0, width, height); + const center = (Math.floor(height / 2) * width + Math.floor(width / 2)) * 4; + data.data[center] = data.data[center] === 255 ? 254 : data.data[center] + 1; + ctx.putImageData(data, 0, 0); + const pngB = await toPng(); + + const digest = (bytes) => new DataView(md5(new Uint8Array(bytes)).buffer).getBigUint64(0); + const md5A = digest(pngA); + const md5B = digest(pngB); + const [hashesA, hashesB] = await Promise.all([ + rpc('hashes', { buffer: pngA }), + rpc('hashes', { buffer: pngB }), + ]); + if (token !== loadToken) return; + + const diff = md5A ^ md5B; + bitGrid(document.getElementById('md5-bits-a'), md5A); + bitGrid(document.getElementById('md5-bits-b'), md5B, 64, diff); + document.getElementById('md5-hex-a').innerText = hex64(md5A) + '…'; + document.getElementById('md5-hex-b').innerText = hex64(md5B) + '…'; + document.getElementById('md5-diff').innerText = `${hamming(md5A, md5B)} of 64 bits changed`; + + // the payoff, delivered where the perceptual hashes have been introduced. + // The target span lives in prose that is sometimes commented out. + const payoff = document.getElementById('onepixel'); + if (payoff) { + const distances = [0, 1, 2].map((i) => hamming(hashesA.hashes[i], hashesB.hashes[i])); + payoff.innerText = + `The three hashes built on this page put them at distance ` + + `${distances[0]}, ${distances[1]} and ${distances[2]}.`; + } +} + +// LSH demo: the hash we are about to build, as a black box. +// Near copy agrees on almost all bits, a stranger on about half. + +const strangerCache = new Map(); +async function runLshDemo(buffer, token) { + // pick a stranger that cannot be the current sample + const stranger = baseUrl === 'img/im21694.jpg' ? 'img/sunflower1.jpg' : 'img/im21694.jpg'; + if (!strangerCache.has(stranger)) { + const response = await fetch(stranger); + strangerCache.set(stranger, await response.arrayBuffer()); + } + const [hue, other] = await Promise.all([ + rpc('mutate', { buffer, kind: 'hue', amount: 90 }), + rpc('hashes', { buffer: strangerCache.get(stranger) }), + ]); + if (token !== loadToken) return; + + document.getElementById('lsh-img-a').src = baseUrl; + document.getElementById('lsh-img-b').src = pngUrl(hue.png); + document.getElementById('lsh-img-c').src = stranger; + const base = baseHashes[0]; + bitGrid(document.getElementById('lsh-bits-a'), base); + bitGrid(document.getElementById('lsh-bits-b'), hue.hashes[0], 64, base ^ hue.hashes[0]); + bitGrid(document.getElementById('lsh-bits-c'), other.hashes[0], 64, base ^ other.hashes[0]); + for (const [id, hash] of [['lsh-d-b', hue.hashes[0]], ['lsh-d-c', other.hashes[0]]]) { + const distance = hamming(base, hash); + const cell = document.getElementById(id); + cell.innerText = distance; + cell.className = distance <= 4 ? 'match' : 'nomatch'; + } +} + +// The 64 DCT basis patterns, rendered once + +(function basisGrid() { + const canvas = document.getElementById('dct-basis'); + const ctx = canvas.getContext('2d'); + const block = 27; + const cell = 3; + + function channel(css, fallback) { + const value = css.match(/#([0-9a-f]{6})/i); + if (!value) return fallback; + const n = parseInt(value[1], 16); + return [(n >> 16) & 255, (n >> 8) & 255, n & 255]; + } + const bg = channel(cssColor('--background'), [16, 16, 24]); + const fg = channel(cssColor('--foreground'), [225, 225, 218]); + + for (let u = 0; u < 8; u++) { + for (let v = 0; v < 8; v++) { + for (let x = 0; x < 8; x++) { + for (let y = 0; y < 8; y++) { + const value = + Math.cos(((2 * x + 1) * u * Math.PI) / 16) * + Math.cos(((2 * y + 1) * v * Math.PI) / 16); + const t = (value + 1) / 2; + const rgb = bg.map((b, i) => Math.round(b + (fg[i] - b) * t)); + ctx.fillStyle = `rgb(${rgb[0]} ${rgb[1]} ${rgb[2]})`; + ctx.fillRect(u * block + x * cell, v * block + y * cell, cell, cell); + } + } + } + } +})(); + +// Hovering a basis pattern fills its (u,v) into the transform formula, +// together with the coefficient it produces for the current image. + +(function basisHover() { + const wrap = document.getElementById('dct-basis-wrap'); + const canvas = document.getElementById('dct-basis'); + const formula = document.getElementById('dct-formula'); + const us = formula.querySelectorAll('.var-u'); + const vs = formula.querySelectorAll('.var-v'); + const value = formula.querySelector('.formula-value'); + const box = document.createElement('div'); + box.className = 'basis-hover'; + wrap.appendChild(box); + let markedCell = null; + + canvas.addEventListener('mousemove', (e) => { + const rect = canvas.getBoundingClientRect(); + const u = Math.min(7, Math.floor(((e.clientX - rect.left) / rect.width) * 8)); + const v = Math.min(7, Math.floor(((e.clientY - rect.top) / rect.height) * 8)); + const block = rect.width / 8; + box.style.display = 'block'; + box.style.left = `${u * block}px`; + box.style.top = `${v * block}px`; + box.style.width = `${block}px`; + box.style.height = `${block}px`; + us.forEach((n) => (n.textContent = u)); + vs.forEach((n) => (n.textContent = v)); + if (baseCoefficients) value.textContent = baseCoefficients[v * 8 + u].toFixed(1); + formula.classList.add('live'); + // mark the matching cell in the coefficient grid + markedCell?.classList.remove('hover-cell'); + markedCell = document.getElementById('image-dct').children[v * 8 + u] || null; + markedCell?.classList.add('hover-cell'); + }); + canvas.addEventListener('mouseleave', () => { + box.style.display = 'none'; + us.forEach((n) => (n.textContent = 'u')); + vs.forEach((n) => (n.textContent = 'v')); + formula.classList.remove('live'); + markedCell?.classList.remove('hover-cell'); + markedCell = null; + }); +})(); + +// Zigzag trace over the coefficient grid + +(function zigzagTrace() { + const svg = document.getElementById('dct-zigzag'); + const line = document.createElementNS('http://www.w3.org/2000/svg', 'polyline'); + const points = ZIGZAG.map((i) => `${(i % 8) * 32 + 16},${Math.floor(i / 8) * 32 + 16}`); + line.setAttribute('points', points.join(' ')); + line.setAttribute('fill', 'none'); + svg.appendChild(line); + const length = line.getTotalLength() || 1200; + line.style.strokeDasharray = length; + line.style.strokeDashoffset = length; + + function play() { + line.classList.remove('play'); + line.style.strokeDashoffset = length; + void line.getBoundingClientRect(); + line.classList.add('play'); + line.style.strokeDashoffset = 0; + } + + document.getElementById('zigzag-replay').addEventListener('click', play); + const observer = new IntersectionObserver((entries) => { + if (entries.some((entry) => entry.isIntersecting)) { + play(); + observer.disconnect(); + } + }); + observer.observe(svg); +})(); + +// Progressive reconstruction: sum the patterns back up in zigzag order. +// A copy of the coefficient grid traces the zigzag along, so image space +// and frequency space fill up together. + +const reconSlider = document.getElementById('recon-k'); +const reconCanvas = document.getElementById('recon-canvas'); + +// Two polylines over the recon grid: the 36 coefficients the hash reads, +// and the tail it ignores. +const reconZigzag = (() => { + const svg = document.getElementById('recon-zigzag'); + const cell = 176 / 8; + const points = FULL_ZIGZAG.map((i) => [(i % 8) * cell + cell / 2, Math.floor(i / 8) * cell + cell / 2]); + const cumulative = [0]; + for (let i = 1; i < 64; i++) { + cumulative.push(cumulative[i - 1] + Math.hypot(points[i][0] - points[i - 1][0], points[i][1] - points[i - 1][1])); + } + function polyline(pts, className) { + const line = document.createElementNS('http://www.w3.org/2000/svg', 'polyline'); + line.setAttribute('points', pts.map(([x, y]) => `${x},${y}`).join(' ')); + line.setAttribute('fill', 'none'); + line.setAttribute('class', className); + svg.appendChild(line); + return line; + } + const head = polyline(points.slice(0, 36), 'zig-head'); + const tail = polyline(points.slice(35), 'zig-tail'); + const headLength = cumulative[35]; + const tailLength = cumulative[63] - cumulative[35]; + head.style.strokeDasharray = headLength; + tail.style.strokeDasharray = tailLength; + + return function update(k) { + const drawn = cumulative[k - 1]; + head.style.strokeDashoffset = headLength - Math.min(drawn, headLength); + tail.style.strokeDashoffset = tailLength - Math.max(0, drawn - headLength); + }; +})(); + +function renderReconstruction() { + if (!baseCoefficients) return; + const k = Number(reconSlider.value); + document.getElementById('recon-count').innerText = k; + const partial = new Array(64).fill(0); + for (let i = 0; i < k; i++) { + partial[FULL_ZIGZAG[i]] = baseCoefficients[FULL_ZIGZAG[i]]; + } + drawPixels(reconCanvas, idct(partial)); + const grid = document.getElementById('recon-grid'); + for (let i = 0; i < 64; i++) { + const cell = grid.children[FULL_ZIGZAG[i]]; + if (cell) cell.style.opacity = i < k ? 1 : 0.15; + } + reconZigzag(k); +} + +let reconPlayToken = 0; + +reconSlider.addEventListener('input', () => { + reconPlayToken++; + renderReconstruction(); +}); + +// Slow pass over the 36 coefficients the hash uses, a beat at the cut, +// then the leftovers: they barely change the picture. +document.getElementById('recon-play').addEventListener('click', async () => { + const token = ++reconPlayToken; + const sleep = (ms) => new Promise((resolve) => setTimeout(resolve, ms)); + for (let k = 1; k <= 64; k++) { + if (token !== reconPlayToken) return; + reconSlider.value = k; + renderReconstruction(); + if (k === 36) await sleep(900); + else await sleep(k < 36 ? 150 : 90); + } +}); + +// The ghost: reconstruction from the 64 hash bits alone. +// Signs come from the sign mask; magnitudes are guessed by walking the +// ordinal chain (each bit says bigger/smaller than the previous, so the +// walk gives every coefficient a relative rank). + +function renderGhost(masks) { + const [signMask, ordinalMask] = masks; + const ranks = [0]; + for (let k = 0; k < 28; k++) { + const bigger = (ordinalMask >> BigInt(27 - k)) & 1n; + ranks.push(ranks[k] + (bigger ? 1 : -1)); + } + const floor = Math.min(...ranks) - 1; + const coefficients = new Array(64).fill(0); + for (let k = 0; k < 36; k++) { + const negative = (signMask >> BigInt(35 - k)) & 1n; + const rank = k < ranks.length ? ranks[k] : floor; + const magnitude = Math.min(400, 24 * Math.pow(1.35, rank)); + coefficients[ZIGZAG[k]] = (negative ? -1 : 1) * magnitude; + } + drawPixels(document.getElementById('ghost-canvas'), idct(coefficients)); +} + +document.getElementById('ghost-flip').addEventListener('click', () => { + document.getElementById('ghost-canvas').classList.toggle('flipped'); +}); + // Mutation demo const mutators = { @@ -211,19 +730,27 @@ function resetMutation() { } let mutationTimer = null; +let mutationToken = 0; async function runMutation() { if (!baseBuffer) return; + const token = ++mutationToken; const kind = kindSelect.value; const amount = amountSlider.disabled ? 0 : Number(amountSlider.value); amountValue.innerText = amountSlider.disabled ? '' : amount; + const busy = [document.getElementById('mutate-compare'), document.querySelector('.distances')]; + busy.forEach((el) => el.classList.add('busy')); const data = await rpc('mutate', { buffer: baseBuffer, kind, amount }); + if (token !== mutationToken || !baseHashes) return; + busy.forEach((el) => el.classList.remove('busy')); - document.getElementById('compare-a').src = baseUrl; - document.getElementById('compare-b').src = pngUrl(data.png); + mutateCompare.set(baseUrl, pngUrl(data.png)); const names = ['dct', 'median', 'phash']; names.forEach((name, i) => { - document.getElementById(`dist-${name}`).innerText = hamming(baseHashes[i], data.hashes[i]); + const distance = hamming(baseHashes[i], data.hashes[i]); + const cell = document.getElementById(`dist-${name}`); + cell.innerText = distance; + cell.className = distance <= 4 ? 'match' : 'nomatch'; }); } @@ -239,62 +766,637 @@ kindSelect.addEventListener('change', () => { amountSlider.addEventListener('input', scheduleMutation); resetMutation(); -document.getElementById('compare-slider').addEventListener('input', (e) => { - document.getElementById('compare-b').style.clipPath = `inset(0 0 0 ${e.target.value}%)`; +// Flip invariance demo + +let flipHashes = null; +document.getElementById('flip-button').addEventListener('click', async () => { + if (!baseBuffer) return; + document.getElementById('flip-image').classList.toggle('flipped'); + if (flipHashes) return; + const token = loadToken; + const data = await rpc('mutate', { buffer: baseBuffer, kind: 'flip', amount: 0 }); + if (token !== loadToken) return; + flipHashes = data.hashes; + const names = ['dct', 'median', 'phash']; + names.forEach((name, i) => { + const diff = baseHashes[i] ^ flipHashes[i]; + bitGrid(document.getElementById(`flip-b-${name}`), flipHashes[i], 64, diff); + const distance = hamming(baseHashes[i], flipHashes[i]); + const cell = document.getElementById(`flip-d-${name}`); + cell.innerText = distance; + cell.className = 'flip-distance ' + (distance <= 4 ? 'match' : 'nomatch'); + }); }); -// Ranking demo +// Ranking demo over the precomputed MIRFLICKR pool -const pool = []; // {name, url, hashes} +const rankingMutations = [ + { key: 'none', label: 'no mutation' }, + { key: 'flip', label: 'horizontal flip', amount: 0 }, + { key: 'jpeg', label: 'JPEG quality 30', amount: 30 }, + { key: 'blur', label: 'Gaussian blur 3', amount: 3 }, + { key: 'hue', label: 'hue shift 60', amount: 60 }, + { key: 'scale', label: 'rescale to 25%', amount: 25 }, + { key: 'crop', label: 'crop to 80%', amount: 80 }, + { key: 'logo', label: 'logo insert 15%', amount: 15 }, + { key: 'noise', label: 'Gaussian noise 20', amount: 20 }, +]; + +const pool = []; // {name, url, hashes: [dct, median, phash]} const poolContainer = document.getElementById('ranking-pool'); -const resultsContainer = document.getElementById('ranking-results'); -const algoSelect = document.getElementById('ranking-algo'); -let queryIndex = null; +const rankingMutationSelect = document.getElementById('ranking-mutation'); +const queryContainer = document.getElementById('ranking-query'); +const bufferCache = new Map(); +let currentQuery = null; // {name, url, hashes, buffer?} +let rankToken = 0; -async function addToPool(name, buffer, url) { +for (const { key, label } of rankingMutations) { + const option = document.createElement('option'); + option.value = key; + option.innerText = label; + rankingMutationSelect.appendChild(option); +} + +const poolReady = fetch('flickr/hashes.json') + .then((response) => response.json()) + .then((entries) => { + const fragment = document.createDocumentFragment(); + for (const [name, dct, median, phash] of entries) { + const url = `flickr/thumbs/${name}`; + const index = pool.length; + pool.push({ + name, + url, + hashes: [BigInt('0x' + dct), BigInt('0x' + median), BigInt('0x' + phash)], + }); + const img = document.createElement('img'); + img.src = url; + img.loading = 'lazy'; + img.decoding = 'async'; + img.title = name; + img.dataset.index = index; + fragment.appendChild(img); + } + poolContainer.appendChild(fragment); + return pool; + }) + .catch(() => null); + +poolContainer.addEventListener('click', (e) => { + const index = e.target.dataset?.index; + if (index === undefined) return; + poolContainer.querySelector('.query')?.classList.remove('query'); + e.target.classList.add('query'); + currentQuery = pool[Number(index)]; + runRanking(); +}); + +document.getElementById('ranking-file').addEventListener('change', async (e) => { + const file = e.target.files[0]; + if (!file) return; + const buffer = await file.arrayBuffer(); const { hashes } = await rpc('hashes', { buffer }); - const index = pool.length; - pool.push({ name, url, hashes }); + poolContainer.querySelector('.query')?.classList.remove('query'); + currentQuery = { name: file.name, url: URL.createObjectURL(file), hashes, buffer }; + runRanking(); +}); - const thumb = document.createElement('figure'); +rankingMutationSelect.addEventListener('change', runRanking); + +async function queryBuffer() { + if (currentQuery.buffer) return currentQuery.buffer; + if (!bufferCache.has(currentQuery.url)) { + const response = await fetch(currentQuery.url); + bufferCache.set(currentQuery.url, await response.arrayBuffer()); + } + return bufferCache.get(currentQuery.url); +} + +function thumbFigure(url, caption, hit = false) { + const figure = document.createElement('figure'); const img = document.createElement('img'); img.src = url; - const caption = document.createElement('figcaption'); - caption.innerText = name; - thumb.append(img, caption); - thumb.addEventListener('click', () => { - queryIndex = index; - rank(); + img.loading = 'lazy'; + const figcaption = document.createElement('figcaption'); + figcaption.innerText = caption; + if (hit) figure.className = 'hit'; + figure.append(img, figcaption); + return figure; +} + +async function runRanking() { + if (!currentQuery) return; + const token = ++rankToken; + const mutation = rankingMutations.find((m) => m.key === rankingMutationSelect.value); + + let queryHashes = currentQuery.hashes; + let mutatedUrl = null; + if (mutation.key !== 'none') { + const buffer = await queryBuffer(); + const data = await rpc('mutate', { buffer, kind: mutation.key, amount: mutation.amount }); + if (token !== rankToken) return; + queryHashes = data.hashes; + mutatedUrl = pngUrl(data.png); + } + + queryContainer.replaceChildren(thumbFigure(currentQuery.url, currentQuery.name)); + if (mutatedUrl) queryContainer.appendChild(thumbFigure(mutatedUrl, mutation.label)); + + ['dct', 'median', 'phash'].forEach((name, algo) => { + const ranked = pool + .map((entry) => ({ entry, distance: hamming(queryHashes[algo], entry.hashes[algo]) })) + .sort((a, b) => a.distance - b.distance) + .slice(0, 10); + const container = document.getElementById(`results-${name}`); + container.replaceChildren(); + for (const { entry, distance } of ranked) { + const caption = `${distance} · ${entry.name.replace('.jpg', '')}`; + container.appendChild(thumbFigure(entry.url, caption, entry.name === currentQuery.name)); + } }); - poolContainer.appendChild(thumb); } -function rank() { - if (queryIndex === null) return; - const algo = Number(algoSelect.value); - const query = pool[queryIndex]; - const ranked = pool - .map((entry) => ({ entry, distance: hamming(query.hashes[algo], entry.hashes[algo]) })) - .sort((a, b) => a.distance - b.distance) - .slice(0, 10); +// Precision-recall charts, rendered from the run data. +// +// 16 series cannot get 16 distinguishable colors, least of all from a pywal +// palette that changes with the wallpaper. So color encodes the mutation +// *category* (fixed assignment) and a per-category dash pattern is the +// secondary encoding; individual mutators are identified by the legend, +// hover isolation and the tooltip, never by color alone. - resultsContainer.replaceChildren(); - for (const { entry, distance } of ranked) { - const thumb = document.createElement('figure'); - const img = document.createElement('img'); - img.src = entry.url; - const caption = document.createElement('figcaption'); - caption.innerText = `${distance} ${entry.name}`; - if (entry === query) thumb.className = 'query'; - thumb.append(img, caption); - resultsContainer.appendChild(thumb); +(async function prCharts() { + let data; + try { + data = await (await fetch('pr-data.json')).json(); + } catch { + return; + } + + const CATEGORIES = [ + { name: 'recoding & resampling', color: '--color2', dash: '', muts: ['JPEG quality 90', 'JPEG quality 50', 'JPEG quality 20', 'Rescale 50%'] }, + { name: 'content processing', color: '--color4', dash: '7 4', muts: ['Hue shift 30', 'Gaussian blur 1.5', 'Unsharp mask', 'Contrast +25', 'Brightness +30', 'Gaussian noise 10'] }, + { name: 'framing', color: '--color3', dash: '2 4', muts: ['Crop to 90%', 'Crop to 70%', 'Letterbox 10%'] }, + { name: 'insertion', color: '--color1', dash: '11 4 2 4', muts: ['Logo insert 10%'] }, + { name: 'flip & rotation', color: '--foreground', dash: '15 5', muts: ['Horizontal flip', 'Rotate 2'] }, + ]; + const styleOf = {}; + for (const category of CATEGORIES) { + for (const mut of category.muts) styleOf[mut] = category; + } + + const ALGOS = ['dct', 'phash', 'median']; + const SIZE = 330; + const MARGIN = { left: 42, right: 12, top: 10, bottom: 38 }; + const px = (r) => MARGIN.left + r * (SIZE - MARGIN.left - MARGIN.right); + const py = (p) => SIZE - MARGIN.bottom - p * (SIZE - MARGIN.top - MARGIN.bottom); + const NS = 'http://www.w3.org/2000/svg'; + const tooltip = document.getElementById('pr-tooltip'); + const tSlider = document.getElementById('pr-t'); + const series = []; // {algo, mut, line, dot, points} + + function el(name, attrs, parent) { + const node = document.createElementNS(NS, name); + for (const [key, value] of Object.entries(attrs)) node.setAttribute(key, value); + if (parent) parent.appendChild(node); + return node; + } + + for (const algo of ALGOS) { + const svg = el('svg', { viewBox: `0 0 ${SIZE} ${SIZE}`, class: 'pr-svg' }); + document.getElementById(`pr-chart-${algo}`).appendChild(svg); + + // recessive grid and axes, text in text tokens + for (const tick of [0, 0.25, 0.5, 0.75, 1]) { + el('line', { x1: px(tick), y1: py(0), x2: px(tick), y2: py(1), class: 'pr-grid' }, svg); + el('line', { x1: px(0), y1: py(tick), x2: px(1), y2: py(tick), class: 'pr-grid' }, svg); + const xLabel = el('text', { x: px(tick), y: py(0) + 14, class: 'pr-tick', 'text-anchor': 'middle' }, svg); + xLabel.textContent = tick; + const yLabel = el('text', { x: px(0) - 5, y: py(tick) + 3, class: 'pr-tick', 'text-anchor': 'end' }, svg); + yLabel.textContent = tick; + } + const xTitle = el('text', { x: px(0.5), y: SIZE - 8, class: 'pr-axis', 'text-anchor': 'middle' }, svg); + xTitle.textContent = 'recall'; + const yTitle = el('text', { + x: 12, y: py(0.5), class: 'pr-axis', 'text-anchor': 'middle', + transform: `rotate(-90 12 ${py(0.5)})`, + }, svg); + yTitle.textContent = 'precision'; + + for (const [mut, points] of Object.entries(data[algo])) { + const category = styleOf[mut]; + if (!category) continue; + const color = cssColor(category.color); + const line = el('polyline', { + points: points.map(([r, p]) => `${px(r)},${py(p)}`).join(' '), + class: 'pr-line', + stroke: color, + 'stroke-dasharray': category.dash, + }, svg); + const dot = el('circle', { r: 4, class: 'pr-dot', fill: color }, svg); + series.push({ algo, mut, line, dot, points }); + } + + svg.addEventListener('mousemove', (event) => { + const rect = svg.getBoundingClientRect(); + const x = ((event.clientX - rect.left) / rect.width) * SIZE; + const y = ((event.clientY - rect.top) / rect.height) * SIZE; + let best = null; + for (const s of series) { + if (s.algo !== algo || s.line.classList.contains('faded')) continue; + s.points.forEach(([r, p], t) => { + const distance = Math.hypot(px(r) - x, py(p) - y); + if (distance < 14 && (!best || distance < best.distance)) { + best = { distance, mut: s.mut, t, r, p }; + } + }); + } + if (best) { + tooltip.hidden = false; + tooltip.innerText = + `${best.mut} · t=${best.t} · ` + + `recall ${(best.r * 100).toFixed(1)}% · precision ${(best.p * 100).toFixed(1)}%`; + tooltip.style.left = `${event.clientX + 14}px`; + tooltip.style.top = `${event.clientY + 14}px`; + } else { + tooltip.hidden = true; + } + }); + svg.addEventListener('mouseleave', () => { + tooltip.hidden = true; + }); + } + + // legend: grouped by category, hover isolates, click pins + const legend = document.getElementById('pr-legend'); + const pinned = new Set(); + let hovered = null; + const chips = new Map(); + + function applyHighlight() { + const active = pinned.size ? pinned : hovered ? new Set([hovered]) : null; + for (const s of series) { + s.line.classList.toggle('faded', Boolean(active && !active.has(s.mut))); + s.dot.classList.toggle('faded', Boolean(active && !active.has(s.mut))); + } + for (const [mut, chip] of chips) { + chip.classList.toggle('dim', Boolean(active && !active.has(mut))); + chip.classList.toggle('pinned', pinned.has(mut)); + } + } + + for (const category of CATEGORIES) { + const group = document.createElement('span'); + group.className = 'pr-group'; + for (const mut of category.muts) { + const chip = document.createElement('button'); + chip.className = 'pr-chip'; + chip.type = 'button'; + const swatch = el('svg', { width: 26, height: 10, class: 'pr-swatch' }); + el('line', { + x1: 1, y1: 5, x2: 25, y2: 5, + stroke: cssColor(category.color), + 'stroke-width': 2, + 'stroke-dasharray': category.dash, + }, swatch); + const label = document.createElement('span'); + label.innerText = mut; + chip.append(swatch, label); + chip.addEventListener('mouseenter', () => { + hovered = mut; + applyHighlight(); + }); + chip.addEventListener('mouseleave', () => { + hovered = null; + applyHighlight(); + }); + chip.addEventListener('click', () => { + if (pinned.has(mut)) pinned.delete(mut); + else pinned.add(mut); + applyHighlight(); + }); + chips.set(mut, chip); + group.appendChild(chip); + } + legend.appendChild(group); + } + + function placeDots() { + const t = Number(tSlider.value); + document.getElementById('pr-t-value').innerText = t; + for (const s of series) { + const [r, p] = s.points[Math.min(t, s.points.length - 1)]; + s.dot.setAttribute('cx', px(r)); + s.dot.setAttribute('cy', py(p)); + } + } + + tSlider.addEventListener('input', placeDots); + placeDots(); +})(); + +// Retrieval demos: threshold histogram, toy BK-tree, and the full index. +// The trees themselves live in wasm memory inside the worker; this side +// only draws structures and traces it gets back. + +const SVG_NS = 'http://www.w3.org/2000/svg'; + +function svgNode(name, attrs, parent) { + const node = document.createElementNS(SVG_NS, name); + for (const [key, value] of Object.entries(attrs)) node.setAttribute(key, value); + if (parent) parent.appendChild(node); + return node; +} + +function fmtTime(ms) { + return ms < 1 ? `${Math.round(ms * 1000)} µs` : `${ms.toFixed(1)} ms`; +} + +// The 425k-hash index is built once, off the critical path +const bigReady = fetch('flickr/keys.bin') + .then((response) => (response.ok ? response.arrayBuffer() : Promise.reject(response.status))) + .then((buffer) => rpc('big_build', { buffer })) + .then(({ size }) => size) + .catch(() => null); + +const retrieval = { + suite: null, // [{hash, label}] of the current image's 16 mutants + labelOf: new Map(), // hex hash -> human name, for toy node tooltips + toy: null, // {byKey, size} of the visualized tree + profile: null, // compared/found counts per radius on the big tree + bigSize: null, +}; + +const distSlider = document.getElementById('dist-t'); +const toySlider = document.getElementById('toy-t'); +const bigSlider = document.getElementById('big-t'); + +async function runRetrieval(buffer, token) { + const data = await rpc('suite', { buffer }); + if (token !== loadToken) return; + retrieval.suite = Array.from(data.hashes).map((hash, i) => ({ hash, label: data.labels[i] })); + + const poolData = await poolReady; + if (token !== loadToken) return; + if (poolData) renderHistogram(); + + // Toy tree: 16 spread-out pool strangers go in first (so the root sits + // in stranger country), then the mutants, which cluster in one branch + const strangers = poolData ? poolData.filter((_, i) => i % 70 === 0).slice(0, 16) : []; + retrieval.labelOf = new Map(); + for (const entry of strangers) retrieval.labelOf.set(entry.hashes[0].toString(16), entry.name); + for (const { hash, label } of retrieval.suite) { + const hex = hash.toString(16); + const existing = retrieval.labelOf.get(hex); + retrieval.labelOf.set(hex, existing ? `${existing} / ${label}` : label); + } + const keys = [...strangers.map((e) => e.hashes[0]), ...retrieval.suite.map((m) => m.hash)]; + const built = await rpc('toy_build', { keys }); + if (token !== loadToken) return; + renderToyTree(built.structure, built.size); + runToyQuery(false); + + const size = await bigReady; + if (token !== loadToken || size === null) return; + retrieval.bigSize = size; + document.getElementById('big-size').innerText = size.toLocaleString('en'); + retrieval.profile = await rpc('big_profile', { key: baseHashes[0], max: Number(bigSlider.max) }); + if (token !== loadToken) return; + renderBigCurve(); + runBigQuery(); +} + +// Distance histogram: strangers as bars, copies as dots, threshold as a line + +function renderHistogram() { + if (!retrieval.suite || !baseHashes || !pool.length) return; + const t = Number(distSlider.value); + document.getElementById('dist-t-value').innerText = t; + const base = baseHashes[0]; + + const strangerCounts = new Array(65).fill(0); + for (const entry of pool) strangerCounts[hamming(base, entry.hashes[0])]++; + const copies = retrieval.suite.map((m) => ({ d: hamming(base, m.hash), label: m.label })); + + const W = 640, H = 210; + const M = { left: 34, right: 10, top: 16, bottom: 30 }; + const step = (W - M.left - M.right) / 65; + const x = (d) => M.left + (d + 0.5) * step; + const maxCount = Math.max(...strangerCounts, 1); + // square root scale, so the lone stranger bars in the overlap zone stay visible + const y = (c) => H - M.bottom - Math.sqrt(c / maxCount) * (H - M.top - M.bottom); + + const svg = svgNode('svg', { viewBox: `0 0 ${W} ${H}`, class: 'pr-svg dist-svg' }); + svgNode('rect', { + x: M.left, y: M.top, + width: Math.max(0, x(t) + step / 2 - M.left), + height: H - M.top - M.bottom, + class: 'dist-region', + }, svg); + for (let d = 0; d <= 64; d += 8) { + const label = svgNode('text', { x: x(d), y: H - M.bottom + 14, class: 'pr-tick', 'text-anchor': 'middle' }, svg); + label.textContent = d; + } + const axis = svgNode('text', { x: x(32), y: H - 2, class: 'pr-axis', 'text-anchor': 'middle' }, svg); + axis.textContent = 'hamming distance'; + svgNode('line', { x1: M.left, y1: H - M.bottom, x2: W - M.right, y2: H - M.bottom, class: 'pr-grid' }, svg); + + for (let d = 0; d <= 64; d++) { + if (!strangerCounts[d]) continue; + const bar = svgNode('rect', { + x: x(d) - step * 0.4, y: y(strangerCounts[d]), + width: step * 0.8, height: H - M.bottom - y(strangerCounts[d]), + class: 'dist-bar', + }, svg); + svgNode('title', {}, bar).textContent = `${strangerCounts[d]} pool images at distance ${d}`; + } + const stacked = new Map(); + for (const copy of copies) { + const k = stacked.get(copy.d) || 0; + stacked.set(copy.d, k + 1); + const dot = svgNode('circle', { + cx: x(copy.d), cy: H - M.bottom - 6 - k * 11, r: 4.5, class: 'dist-copy', + }, svg); + svgNode('title', {}, dot).textContent = `${copy.label}, distance ${copy.d}`; + } + svgNode('line', { + x1: x(t) + step / 2, y1: M.top - 4, x2: x(t) + step / 2, y2: H - M.bottom, + class: 'dist-thresh', + }, svg); + + document.getElementById('dist-chart').replaceChildren(svg); + const copiesIn = copies.filter((c) => c.d <= t).length; + const strangersIn = strangerCounts.slice(0, t + 1).reduce((a, b) => a + b, 0); + document.getElementById('dist-readout').innerText = + `catches ${copiesIn} of ${copies.length} copies, lets in ${strangersIn} of ${pool.length} strangers`; +} + +distSlider.addEventListener('input', renderHistogram); + +// Toy BK-tree: drawn from the wasm tree's own structure dump, +// colored by the visit trace of a real query + +function renderToyTree(structure, size) { + const byKey = new Map(); + let root = null; + for (let i = 0; i < structure.length; i += 3) { + const [key, parent, distance] = [structure[i], structure[i + 1], structure[i + 2]]; + const node = { key, distance: Number(distance), children: [], parent: null, el: {} }; + if (root === null) { + root = node; + } else { + node.parent = byKey.get(parent.toString(16)); + node.parent.children.push(node); + } + byKey.set(key.toString(16), node); + } + retrieval.toy = { byKey, size }; + if (!root) return; + + for (const node of byKey.values()) node.children.sort((a, b) => a.distance - b.distance); + let nextLeaf = 0; + let maxDepth = 0; + (function layout(node, depth) { + node.depth = depth; + maxDepth = Math.max(maxDepth, depth); + if (!node.children.length) { + node.slot = nextLeaf++; + return; + } + node.children.forEach((child) => layout(child, depth + 1)); + node.slot = (node.children[0].slot + node.children[node.children.length - 1].slot) / 2; + })(root, 0); + + const W = 960; + const rowH = 58; + const H = maxDepth * rowH + 60; + const px = (node) => 30 + (node.slot * (W - 60)) / Math.max(nextLeaf - 1, 1); + const py = (node) => 30 + node.depth * rowH; + const svg = svgNode('svg', { viewBox: `0 0 ${W} ${H}`, class: 'pr-svg toy-svg' }); + + // edges under nodes + for (const node of byKey.values()) { + if (!node.parent) continue; + node.el.edge = svgNode('line', { + x1: px(node.parent), y1: py(node.parent), x2: px(node), y2: py(node), class: 'toy-edge', + }, svg); + node.el.edgeLabel = svgNode('text', { + x: (px(node.parent) + px(node)) / 2, + y: (py(node.parent) + py(node)) / 2 - 3, + class: 'toy-edge-label', 'text-anchor': 'middle', + }, svg); + node.el.edgeLabel.textContent = node.distance; + } + for (const node of byKey.values()) { + node.el.circle = svgNode('circle', { cx: px(node), cy: py(node), r: 8, class: 'toy-node' }, svg); + const hex = node.key.toString(16); + svgNode('title', {}, node.el.circle).textContent = + `${retrieval.labelOf.get(hex) || 'pool image'} · ${hex64(node.key)}`; + } + document.getElementById('toy-tree').replaceChildren(svg); +} + +function setToyState(node, state) { + node.el.circle?.setAttribute('class', `toy-node ${state}`); + const edgeState = state === 'pruned' ? ' pruned' : ''; + node.el.edge?.setAttribute('class', `toy-edge${edgeState}`); + node.el.edgeLabel?.setAttribute('class', `toy-edge-label${edgeState}`); +} + +let toyPaintToken = 0; + +async function runToyQuery(animate) { + if (!retrieval.toy || !baseHashes) return; + const token = loadToken; + const paintToken = ++toyPaintToken; + const t = Number(toySlider.value); + document.getElementById('toy-t-value').innerText = t; + const { found, visited } = await rpc('toy_query', { key: baseHashes[0], radius: t }); + if (token !== loadToken || paintToken !== toyPaintToken) return; + + const foundSet = new Set(); + for (let i = 0; i < found.length; i += 2) foundSet.add(found[i].toString(16)); + document.getElementById('toy-stats').innerText = + `compared ${visited.length} of ${retrieval.toy.size} hashes, found ${found.length / 2}`; + + for (const node of retrieval.toy.byKey.values()) setToyState(node, 'pruned'); + const sleep = (ms) => new Promise((resolve) => setTimeout(resolve, ms)); + for (const key of visited) { + if (animate) await sleep(140); + if (paintToken !== toyPaintToken) return; + const hex = key.toString(16); + setToyState(retrieval.toy.byKey.get(hex), foundSet.has(hex) ? 'match' : 'seen'); } } -algoSelect.addEventListener('change', rank); - -document.getElementById('ranking-files').addEventListener('change', (e) => { - for (const file of e.target.files) { - file.arrayBuffer().then((buffer) => addToPool(file.name, buffer, URL.createObjectURL(file))); - } +let toyTimer = null; +toySlider.addEventListener('input', () => { + clearTimeout(toyTimer); + toyTimer = setTimeout(() => runToyQuery(false), 120); +}); +document.getElementById('toy-run').addEventListener('click', () => runToyQuery(true)); + +// The full index: instant numbers from the precomputed profile, +// timings measured in the worker on demand + +function renderBigCurve() { + const { compared } = retrieval.profile; + const size = retrieval.bigSize; + const W = 260, H = 150; + const M = { left: 40, right: 10, top: 8, bottom: 26 }; + const maxT = compared.length - 1; + const x = (t) => M.left + (t * (W - M.left - M.right)) / maxT; + const y = (frac) => H - M.bottom - frac * (H - M.top - M.bottom); + + const svg = svgNode('svg', { viewBox: `0 0 ${W} ${H}`, class: 'pr-svg big-curve-svg' }); + for (const frac of [0, 0.5, 1]) { + svgNode('line', { x1: M.left, y1: y(frac), x2: W - M.right, y2: y(frac), class: 'pr-grid' }, svg); + const label = svgNode('text', { x: M.left - 4, y: y(frac) + 3, class: 'pr-tick', 'text-anchor': 'end' }, svg); + label.textContent = `${frac * 100}%`; + } + for (const t of [0, 4, 8, 12, 16]) { + if (t > maxT) continue; + const label = svgNode('text', { x: x(t), y: H - M.bottom + 13, class: 'pr-tick', 'text-anchor': 'middle' }, svg); + label.textContent = t; + } + const axis = svgNode('text', { x: x(maxT / 2), y: H - 2, class: 'pr-axis', 'text-anchor': 'middle' }, svg); + axis.textContent = 'radius t'; + svgNode('polyline', { + points: compared.map((c, t) => `${x(t)},${y(c / size)}`).join(' '), + class: 'pr-line big-curve-line', + }, svg); + retrieval.curveDot = svgNode('circle', { r: 4, class: 'pr-dot big-curve-dot' }, svg); + retrieval.curveX = x; + retrieval.curveY = y; + document.getElementById('big-curve').replaceChildren(svg); +} + +let bigTimer = null; +let bigQueryToken = 0; + +async function runBigQuery() { + if (!retrieval.profile || !baseHashes) return; + const t = Number(bigSlider.value); + document.getElementById('big-t-value').innerText = t; + const { compared, found } = retrieval.profile; + const size = retrieval.bigSize; + document.getElementById('big-found').innerText = found[t].toLocaleString('en'); + document.getElementById('big-found-note').innerText = `hashes within distance ${t}`; + document.getElementById('big-compared').innerText = compared[t].toLocaleString('en'); + const pct = (100 * compared[t]) / size; + document.getElementById('big-times').innerText = + `${pct < 10 ? pct.toFixed(1) : Math.round(pct)}% of the tree`; + retrieval.curveDot?.setAttribute('cx', retrieval.curveX(t)); + retrieval.curveDot?.setAttribute('cy', retrieval.curveY(compared[t] / size)); + + const token = ++bigQueryToken; + const outer = loadToken; + const data = await rpc('big_query', { key: baseHashes[0], radius: t }); + if (token !== bigQueryToken || outer !== loadToken) return; + document.getElementById('big-times').innerText = + `${pct < 10 ? pct.toFixed(1) : Math.round(pct)}% of the tree · tree ${fmtTime(data.treeMs)}, scan ${fmtTime(data.scanMs)}`; +} + +bigSlider.addEventListener('input', () => { + clearTimeout(bigTimer); + bigTimer = setTimeout(runBigQuery, 120); }); diff --git a/web/style.css b/web/style.css index d437cb4..b410f9c 100644 --- a/web/style.css +++ b/web/style.css @@ -2,19 +2,33 @@ body { margin: 0; padding: 0; color: var(--foreground); - background: var(--background) url('bg.jpeg'); - background-size: cover; + background: var(--background); font-family: sans-serif; - height: 100%; - overflow: hidden; font-size: 14px; } -#wrapper { - overflow-x: hidden; - overflow-y: auto; +/* Fixed background layer, same construction as the landing page */ + +#bg { + position: fixed; + width: 100vw; + height: 100vh; + z-index: -1; + background: var(--background); +} + +#bg img { + display: block; width: 100%; height: 100%; + object-fit: cover; + object-position: 50% 50%; + opacity: 0; + transition: opacity 3s linear; +} + +#bg img.loaded { + opacity: 1; } h1, h2, h3, h4, h5, h6 { @@ -25,7 +39,7 @@ h1, h2, h3, h4, h5, h6 { h1 { margin-top: 0; - font-size: 24px; + font-size: 28px; } h2 { @@ -34,45 +48,98 @@ h2 { p { margin-top: 2px; + line-height: 1.5; } a { color: var(--color2); + text-decoration: none; + font-weight: bold; } -strong { - color: var(--color4); +a:hover { + color: var(--color3); } em { color: var(--color3); } +sup.fn { + font-size: 10px; +} + +sup.fn a { + font-weight: normal; +} + pre { - margin: 0; + margin: 8px 0; + max-width: 640px; } pre > code { - border: 1px solid; + border: 1px solid color-mix(in srgb, var(--foreground) 25%, transparent); + display: block; + padding: 8px 12px; + overflow-x: auto; + font-size: 13px; } article { max-width: 1024px; - min-height: 128px; margin: 25px auto; padding: 20px; background: color-mix(in srgb, var(--background) 90%, transparent); box-shadow: 5px 5px 10px rgba(0, 0, 0, 0.9); } -code.block { - display: block; +@media (max-width: 1064px) { + article { + margin: 15px 12px; + } } img { max-width: 100%; } +math[display="block"] { + margin: 12px 0 12px 24px; + font-size: 17px; + color: var(--color4); +} + +aside.explainer { + border-left: 3px solid var(--color4); + background: color-mix(in srgb, var(--foreground) 5%, transparent); + padding: 8px 12px; + margin: 12px 0; + font-size: 13px; +} + +.rewrite-target { + text-decoration: underline dashed var(--color3); + text-underline-offset: 3px; +} + +.rewrite-target::after { + content: ' ✎'; + color: var(--color3); + font-size: 11px; +} + +.attribution { + font-size: 12px; + opacity: 0.7; +} + +.big-number { + font-size: 26px; + color: var(--color1); + margin: 8px 0 4px 0; +} + /* Demo */ .needs-image { @@ -88,135 +155,107 @@ img { min-width: 256px; } -/* Resize demo */ +.sample { + height: 48px; + margin: 0 4px 4px 0; + cursor: pointer; + vertical-align: middle; +} -#resize { +.sample:hover { + outline: 2px solid var(--color1); +} + +/* Comparison slider */ + +.compare { position: relative; - height: 512px; + max-width: 512px; + overflow: hidden; + touch-action: none; + cursor: ew-resize; + user-select: none; +} + +.compare img { + display: block; width: 100%; - transition: height 2s linear; + pointer-events: none; } -#resize.resize { - height: 256px; -} - -#resize .image-original { +/* The top image stretches to the base image's box, so differently + sized results (8x8 thumbs, crops) stay comparable under the divider */ +.compare .compare-b { position: absolute; - width: 100%; - height: 512px; + top: 0; + left: 0; + height: 100%; + object-fit: fill; + clip-path: inset(0 0 0 50%); } -#resize .image-original img { - max-height: 50vh; - min-height: 257px; - max-width: 100%; - min-width: 257px; - /* height: 512px; */ -} - -#resize-original { - position: relative; - height: 512px; -} - -#resize-original img { +.compare .compare-divider { position: absolute; - max-width: 90vw; - max-height: 512px; - min-width: 257px; - min-height: 257px; - opacity: 1; - transition: - width 2s linear, - height 2s linear, - filter 3s linear 2s, - opacity 3s ease-in 5s; + top: 0; + bottom: 0; + left: 50%; + width: 2px; + margin-left: -1px; + background: var(--color1); } -#resize.resize #resize-original img { - width: 256px !important; - height: 256px !important; - filter: grayscale(1) blur(10px); - opacity: 0; -} - -#resize-original img.grayscale { - filter: grayscale(1) blur(10px); -} - -#resize-original img.fade { - opacity: 0; -} - -#image-resize img { - width: 256px; - height: 256px; - image-rendering: pixelated; - image-rendering: -moz-crisp-edges; -} - -#image-dct { - margin-left: 32px; - display:grid; - grid-template-columns: repeat(8, 32px); -} - -#image-dct > div { - height: 32px; - line-height:32px; - font-size: 11px; +.compare .compare-handle { + position: absolute; + top: 50%; + left: 50%; + transform: translate(-50%, -50%); + width: 28px; + height: 28px; + line-height: 28px; text-align: center; - border: 1px solid var(--color1); - box-sizing: border-box; + border-radius: 50%; + background: var(--color1); + color: var(--background); + font-size: 15px; } -.similarity-example { - display: flex; - flex-flow: row nowrap; - align-items: flex-end; +.compare .compare-label { + position: absolute; + top: 6px; + padding: 2px 6px; + font-size: 11px; + background: color-mix(in srgb, var(--background) 75%, transparent); + color: var(--foreground); } -figure { - margin: 10px; - width: 50%; +.compare .compare-label.a { + left: 6px; } -/* Loader */ - -#loader.show { - border-radius: 50%; - width: 256px; - height: 256px; - font-size: 10px; - position: fixed; - z-index: 99; - top: calc(50vh - 128px); - left: calc(50vw - 128px); - border: 16px solid transparent; - border-top: 16px solid darkred; - transform: translateZ(0); - animation: loader 2s infinite ease-in-out; +.compare .compare-label.b { + right: 6px; } -@keyframes loader { - 0% { - transform: rotate(0deg); - } - 50% { - transform: rotate(360deg); - } - 100% { - transform: rotate(360deg); - } -} -/* Skeleton */ +/* Busy indicator while the worker renders a mutation */ -.todo { - opacity: 0.55; - font-style: italic; - border-left: 3px solid var(--color3); - padding-left: 8px; +.compare.busy::after { + content: 'mutating…'; + position: absolute; + right: 8px; + bottom: 8px; + padding: 3px 9px; + background: color-mix(in srgb, var(--background) 80%, transparent); + color: var(--color1); + font-size: 12px; + animation: busy-pulse 1s ease-in-out infinite; +} + +.distances.busy td { + opacity: 0.35; +} + +@keyframes busy-pulse { + 50% { opacity: 0.35; } } /* Hash visualisation */ @@ -233,6 +272,14 @@ figure { letter-spacing: 2px; } +.hash-value.small { + display: block; + margin-top: 4px; + font-size: 12px; + letter-spacing: 1px; + opacity: 0.8; +} + .bit-grid { display: grid; grid-template-columns: repeat(8, 20px); @@ -248,6 +295,14 @@ figure { border: 1px solid var(--background); } +.bit-grid.small { + grid-template-columns: repeat(8, 13px); +} + +.bit-grid.small .bit { + height: 13px; +} + .bit-grid .bit.on { background: var(--foreground); } @@ -256,6 +311,20 @@ figure { background: color-mix(in srgb, var(--foreground) 12%, transparent); } +/* Changed bits keep their on/off brightness and get a red tint + border, + so light-changed and dark-changed cells stay distinguishable */ +.bit-grid .bit.diff { + border-color: hsl(4 70% 55%); +} + +.bit-grid .bit.on.diff { + background: color-mix(in srgb, var(--foreground) 55%, hsl(4 70% 50%)); +} + +.bit-grid .bit.off.diff { + background: hsl(4 60% 50% / 0.25); +} + .heat-grid { display: grid; grid-template-columns: repeat(8, 32px); @@ -266,7 +335,15 @@ figure { box-sizing: border-box; } -.pixelated, #image-resize img, #phash-resize img { +.heat-grid.mini { + grid-template-columns: repeat(8, 22px); +} + +.heat-grid.mini > div { + height: 22px; +} + +.pixelated, #phash-resize img, .input-thumb img { image-rendering: pixelated; image-rendering: -moz-crisp-edges; } @@ -276,71 +353,520 @@ figure { height: 128px; } -.sample { - height: 48px; - margin: 0 4px; - cursor: pointer; - vertical-align: middle; +.input-thumb img { + width: 176px; + height: 176px; } -/* Mutation compare */ +/* DCT extras */ -.compare { +.grid-overlay { position: relative; - max-width: 640px; + width: 256px; } -.compare img { - display: block; - width: 100%; +.grid-overlay.mini { + width: 176px; } -.compare #compare-b { +#dct-basis-wrap { + position: relative; + display: inline-block; +} + +.basis-hover { + display: none; + position: absolute; + outline: 2px solid var(--color1); + outline-offset: -1px; + pointer-events: none; +} + +#dct-basis { + cursor: crosshair; +} + +#dct-formula .formula-result { + display: none; +} + +#dct-formula.live .formula-result { + display: inline; +} + +#dct-formula.live .var-u, +#dct-formula.live .var-v, +#dct-formula.live .formula-value { + color: var(--color1); +} + +#image-dct .hover-cell { + outline: 2px solid var(--color1); + outline-offset: -2px; +} + +#dct-zigzag { position: absolute; top: 0; left: 0; - clip-path: inset(0 0 0 50%); + pointer-events: none; } -#compare-slider, #mutator-amount { - width: 320px; +#dct-zigzag polyline { + stroke: var(--color1); + stroke-width: 3; + stroke-linejoin: round; +} + +#dct-zigzag polyline.play { + transition: stroke-dashoffset 2.5s ease-in-out; +} + +#dct-basis { + background: var(--background); +} + +/* Reconstruction demos */ + +.recon-row { + display: flex; + flex-flow: row wrap; + gap: 24px; + align-items: flex-start; + margin: 10px 0; +} + +.recon-row figure { + margin: 0; + text-align: center; +} + +.recon-row figcaption { + font-size: 12px; + opacity: 0.8; + margin-top: 4px; +} + +.recon-row canvas { + width: 176px; + height: 176px; + image-rendering: pixelated; + image-rendering: -moz-crisp-edges; +} + +#recon-k { + width: 256px; + vertical-align: middle; +} + +.recon-note { + font-size: 12px; + opacity: 0.7; + margin-left: 8px; +} + +#recon-zigzag { + position: absolute; + top: 0; + left: 0; + pointer-events: none; +} + +#recon-zigzag polyline { + stroke-width: 2.5; + stroke-linejoin: round; + transition: stroke-dashoffset 0.12s linear; +} + +#recon-zigzag .zig-head { + stroke: var(--color1); +} + +#recon-zigzag .zig-tail { + stroke: color-mix(in srgb, var(--foreground) 45%, transparent); +} + +/* Tables */ + +table.facts { + border-collapse: collapse; + margin: 10px 0; +} + +table.facts th, table.facts td { + text-align: left; + padding: 4px 14px 4px 0; + border-bottom: 1px solid color-mix(in srgb, var(--foreground) 20%, transparent); + font-size: 13px; +} + +table.facts th { + color: var(--color1); + font-weight: normal; } .distances td { text-align: center; font-size: 18px; min-width: 64px; +} + +.distances td.match { + color: hsl(135 50% 45%); +} + +.distances td.nomatch { + color: hsl(4 70% 55%); +} + +/* LSH demo */ + +.lsh-thumb img { + height: 110px; + margin-bottom: 6px; + display: block; +} + +.lsh-distance { + font-size: 16px; +} + +.lsh-distance .match { + font-size: 22px; + color: hsl(135 50% 45%); +} + +.lsh-distance .nomatch { + font-size: 22px; + color: hsl(4 70% 55%); +} + +/* Flip demo */ + +#flip-stage img { + max-width: 256px; + max-height: 256px; + transition: transform 0.7s ease-in-out; +} + +#flip-stage img.flipped { + transform: scaleX(-1); +} + +#ghost-canvas { + transition: transform 0.4s ease-in-out; +} + +#ghost-canvas.flipped { + transform: scaleX(-1); +} + +.flip-table th { color: var(--color1); + font-weight: normal; + text-align: left; + padding-right: 12px; +} + +.flip-table td { + vertical-align: top; + padding: 6px 12px 6px 0; +} + +.flip-distance { + font-size: 24px; +} + +.flip-distance.match { + color: hsl(135 50% 45%); +} + +.flip-distance.nomatch { + color: hsl(4 70% 55%); } /* Ranking */ -.thumb-grid { +.pool { + display: grid; + grid-template-columns: repeat(auto-fill, minmax(72px, 1fr)); + gap: 4px; + max-height: 330px; + overflow-y: auto; + margin: 10px 0; +} + +.pool img { + width: 100%; + aspect-ratio: 1; + object-fit: cover; + cursor: pointer; + display: block; +} + +.pool img:hover { + outline: 2px solid var(--color1); +} + +.pool img.query { + outline: 3px solid var(--color1); +} + +.ranking-results h3 { + margin-top: 14px; +} + +.thumb-row { display: flex; flex-flow: row wrap; gap: 8px; } -.thumb-grid figure { - width: 96px; +.thumb-row figure { + width: 84px; margin: 0; - cursor: pointer; text-align: center; } -.thumb-grid figure.query { - outline: 2px solid var(--color1); +.thumb-row figure.hit { + outline: 2px solid hsl(135 50% 45%); } -.thumb-grid img { - width: 96px; - height: 96px; +.thumb-row img { + width: 84px; + height: 84px; object-fit: cover; + display: block; } -.thumb-grid figcaption { +.thumb-row figcaption { font-size: 11px; overflow: hidden; text-overflow: ellipsis; white-space: nowrap; } + +/* Precision-recall charts */ + +.pr-controls { + margin: 10px 0; +} + +#pr-legend { + display: flex; + flex-flow: row wrap; + gap: 4px 14px; + margin-bottom: 8px; +} + +.pr-group { + display: flex; + flex-flow: row wrap; + gap: 4px; +} + +.pr-chip { + display: inline-flex; + align-items: center; + gap: 5px; + border: 1px solid color-mix(in srgb, var(--foreground) 20%, transparent); + background: none; + color: var(--foreground); + font-size: 11px; + padding: 2px 7px; + cursor: pointer; + transition: opacity 0.15s; +} + +.pr-chip.dim { + opacity: 0.3; +} + +.pr-chip.pinned { + border-color: var(--color1); +} + +.pr-slider { + display: inline-block; + margin: 4px 0; + color: var(--color4); +} + +.pr-slider input { + width: 220px; + vertical-align: middle; + margin-left: 8px; +} + +.pr-charts { + display: flex; + flex-flow: row wrap; + gap: 12px; +} + +.pr-charts figure { + margin: 0; + flex: 1 1 300px; + max-width: 340px; +} + +.pr-charts figcaption { + color: var(--color1); + font-size: 14px; + margin-bottom: 2px; +} + +.pr-svg { + width: 100%; + height: auto; + background: color-mix(in srgb, var(--foreground) 3%, transparent); +} + +.pr-grid { + stroke: color-mix(in srgb, var(--foreground) 12%, transparent); + stroke-width: 1; +} + +.pr-tick { + fill: color-mix(in srgb, var(--foreground) 55%, transparent); + font-size: 9px; +} + +.pr-axis { + fill: color-mix(in srgb, var(--foreground) 70%, transparent); + font-size: 11px; +} + +.pr-line { + fill: none; + stroke-width: 1.8; + transition: opacity 0.15s; +} + +.pr-dot { + stroke: var(--background); + stroke-width: 2; + transition: opacity 0.15s; +} + +.pr-line.faded, .pr-dot.faded { + opacity: 0.12; +} + +#pr-tooltip { + position: fixed; + z-index: 10; + background: color-mix(in srgb, var(--background) 92%, transparent); + border: 1px solid var(--color1); + color: var(--foreground); + font-size: 12px; + padding: 4px 8px; + pointer-events: none; + white-space: nowrap; +} + +/* Retrieval demos */ + +.dist-figure, .toy-figure { + margin: 10px 0; +} + +.dist-figure figcaption, .toy-figure figcaption { + font-size: 12px; + opacity: 0.8; + margin-top: 4px; +} + +.dist-svg { + max-width: 720px; +} + +.dist-bar { + fill: color-mix(in srgb, var(--foreground) 45%, transparent); +} + +.dist-copy { + fill: hsl(135 50% 45%); + stroke: var(--background); + stroke-width: 1; +} + +.dist-region { + fill: color-mix(in srgb, var(--color1) 10%, transparent); +} + +.dist-thresh { + stroke: var(--color1); + stroke-width: 2; +} + +.dist-readout, .toy-stats { + margin-left: 12px; + color: var(--color4); +} + +.toy-node { + fill: color-mix(in srgb, var(--foreground) 30%, transparent); + transition: opacity 0.3s, fill 0.3s; +} + +.toy-node.seen { + fill: color-mix(in srgb, var(--foreground) 10%, var(--background)); + stroke: var(--foreground); + stroke-width: 2; +} + +.toy-node.match { + fill: hsl(135 50% 45%); +} + +.toy-node.pruned { + opacity: 0.15; +} + +.toy-edge { + stroke: color-mix(in srgb, var(--foreground) 45%, transparent); + stroke-width: 1.5; + transition: opacity 0.3s; +} + +.toy-edge-label { + fill: color-mix(in srgb, var(--foreground) 70%, transparent); + font-size: 9px; + transition: opacity 0.3s; +} + +.toy-edge.pruned, .toy-edge-label.pruned { + opacity: 0.15; +} + +.big-stats .stat-note { + font-size: 12px; + opacity: 0.7; + margin: 0; + max-width: 220px; +} + +.big-curve-svg { + width: 250px; +} + +.big-curve-line { + stroke: var(--color2); +} + +.big-curve-dot { + fill: var(--color1); +} + +/* Lists */ + +ul.loose-ends li, ol.references li { + margin-bottom: 6px; + line-height: 1.5; +} + +ol.references li:target { + outline: 1px dashed var(--color1); + outline-offset: 4px; +} diff --git a/web/worker.js b/web/worker.js index 178aceb..153e891 100644 --- a/web/worker.js +++ b/web/worker.js @@ -1,15 +1,32 @@ -import init, { hash_all, resize_preview, dct_coefficients, dct_masks, phash_lowfreq, mutate } from './pkg/image_similarity.js'; +import init, { hash_all, resize_preview, dct_coefficients, dct_masks, phash_lowfreq, mutate, suite_hashes, suite_labels, BkIndex } from './pkg/image_similarity.js'; const ready = init(); +// The retrieval demos keep their trees in wasm memory between calls: +// a small one that gets visualized and the full experiment index. +let toyIndex = null; +let bigIndex = null; + +// Average runtime of fn, sampled until enough wall time passed to +// outrun the clamped performance.now resolution +function time(fn) { + const start = performance.now(); + let reps = 0; + do { + fn(); + reps++; + } while (performance.now() - start < 20); + return (performance.now() - start) / reps; +} + onmessage = async (e) => { const { id, op, args } = e.data; await ready; try { - const bytes = new Uint8Array(args.buffer); let result; switch (op) { - case 'pipeline': + case 'pipeline': { + const bytes = new Uint8Array(args.buffer); result = { resize8: resize_preview(bytes, 8), resize32: resize_preview(bytes, 32), @@ -19,14 +36,63 @@ onmessage = async (e) => { hashes: hash_all(bytes), }; break; + } case 'hashes': - result = { hashes: hash_all(bytes) }; + result = { hashes: hash_all(new Uint8Array(args.buffer)) }; break; case 'mutate': { - const png = mutate(bytes, args.kind, args.amount); + const png = mutate(new Uint8Array(args.buffer), args.kind, args.amount); result = { png, hashes: hash_all(png) }; break; } + case 'suite': + result = { + hashes: suite_hashes(new Uint8Array(args.buffer)), + labels: suite_labels().split('\n'), + }; + break; + case 'toy_build': { + toyIndex?.free(); + toyIndex = new BkIndex(); + for (const key of args.keys) toyIndex.insert(key); + result = { structure: toyIndex.structure(), size: toyIndex.len() }; + break; + } + case 'toy_query': + result = { + found: toyIndex.find(args.key, args.radius), + visited: toyIndex.trace(args.key, args.radius), + }; + break; + case 'big_build': { + bigIndex?.free(); + bigIndex = new BkIndex(); + bigIndex.insert_bytes(new Uint8Array(args.buffer)); + result = { size: bigIndex.len() }; + break; + } + case 'big_query': { + const found = bigIndex.find(args.key, args.radius); + result = { + found, + compared: bigIndex.find_compared(args.key, args.radius), + treeMs: time(() => bigIndex.find(args.key, args.radius)), + scanMs: time(() => bigIndex.scan(args.key, args.radius)), + }; + break; + } + // compared counts and found counts for every radius up to max, + // so the radius slider works without a round trip per step + case 'big_profile': { + const compared = []; + const found = []; + for (let t = 0; t <= args.max; t++) { + compared.push(bigIndex.find_compared(args.key, t)); + found.push(bigIndex.find(args.key, t).length / 2); + } + result = { compared, found }; + break; + } default: throw new Error(`unknown op ${op}`); }