Your browser sucks, JavaScript disabled! Gunnar's console needs JavaScript and WebAssembly. Git over SSH on port 2222 does not.

nordisk/znippy-zoomies

osm-katana/tests/selfx_pack_budget.rs

//! **Row 24 on APPLIED OUTPUT — a packed rung column carries no bowtie, and the pack
//! says how many rings it had to ship unsimplified to keep it that way.**
//!
//! Every assertion here reads a parquet file back off disk, through the census
//! `osm-katana verify --geometry-quality` itself runs ([`osm_katana::verify::
//! geometry_quality`]). Nothing compares the packer with a copy of the packer.
//!
//! # What is being guarded
//!
//! `tiles::simplify_polys` retries a Douglas–Peucker pass at half the epsilon while the
//! simplified ring crosses itself, and ships the SOURCE ring when the whole ladder
//! fails. Row 24 bounded the crossing test that decides this — it was exact O(n²) with
//! no cap — and a bounded test has a THIRD answer: "I ran out of budget". The danger a
//! bound introduces is that the pack starts PASSING rings it never checked, and the
//! bowties come back through the door the row opened.
//!
//! So: pack a table of rings chosen to be adversarial at every rung, census the rung
//! columns, and require `self_intersecting == 0` — with `big_ring_unchecked` read out
//! loud, because a census that skipped the rings is not a census that passed them.
//!
//! | test | the defect it is red for |
//! |---|---|
//! | [`a_packed_rung_column_carries_no_self_intersecting_ring`] | a budget-capped crossing test that passes what it did not check |
//! | [`the_rungs_are_smaller_than_the_base`] | shipping the source ring so often that the rung columns stop being a simplification |
//!
//! Run: `cargo test -p osm-katana --test selfx_pack_budget -- --nocapture`

use std::collections::HashMap;
use std::path::{Path, PathBuf};
use std::sync::Arc;

use arrow::array::{Array, BinaryArray, Int64Array, RecordBatch, StringArray};
use arrow::datatypes::{DataType, Field, Schema};
use osm_katana::metadata::{self, GeoMeta};
use osm_katana::shared::schema::{GENERALISED_RUNGS, generalised_geometry_name};
use osm_katana::spatial::{PackOptions, spatial_pack};
use parquet::arrow::arrow_reader::ParquetRecordBatchReaderBuilder;

/// `natural=water` — a Region-tier area, so `min_zoom` admits it at every rung this
/// test reads. Taken through the one ladder writer rather than asserted here.
const TAGS: &str = r#"{"natural":"water"}"#;

fn scratch(stem: &str) -> PathBuf {
    let d = PathBuf::from(std::env::var_os("HOME").unwrap_or_default())
        .join("scratch/q-lane-test")
        .join(format!("selfx-pack-{}", std::process::id()));
    std::fs::create_dir_all(&d).unwrap();
    d.join(format!("{stem}.parquet"))
}

/// The real shore ring of `tiles::row24_selfx_budget_tests` (relation 6208486, off the
/// live Sweden cut) — 72 vertices whose plain-DP output crosses itself at the z9 rung.
fn shore() -> Vec<(f64, f64)> {
    vec![(12.6251364, 56.127182), (12.6249933, 56.1274223), (12.624753, 56.1276093), (12.6246948, 56.1276474), (12.62465, 56.1277199), (12.6246338, 56.1278915), (12.6246414, 56.1279984), (12.6245689, 56.1280174), (12.6245127, 56.1280136), (12.6244593, 56.1279678), (12.6243916, 56.127903), (12.6241989, 56.1278534), (12.6239948, 56.1278648), (12.6238718, 56.1278992), (12.6238184, 56.1279221), (12.6237698, 56.1279678), (12.6237593, 56.1280365), (12.6237755, 56.1281967), (12.6237755, 56.1282959), (12.623724, 56.1284103), (12.6236324, 56.128479), (12.6234446, 56.1285362), (12.6229115, 56.1285362), (12.6224241, 56.1286011), (12.6219435, 56.1286392), (12.6212893, 56.1287308), (12.6211252, 56.1287613), (12.6210012, 56.1288452), (12.6203451, 56.1292229), (12.6203766, 56.1299858), (12.6198893, 56.1302452), (12.6189117, 56.1303864), (12.6188459, 56.1303749), (12.6187544, 56.1302643), (12.6187115, 56.1301651), (12.6187057, 56.1300621), (12.6187887, 56.1301384), (12.6188459, 56.1301308), (12.6189013, 56.1301041), (12.620429, 56.1289291), (12.6204996, 56.1288719), (12.620595, 56.1288261), (12.6207476, 56.1287689), (12.6209793, 56.1287155), (12.6214828, 56.1286545), (12.6226549, 56.1285248), (12.6230574, 56.1285057), (12.6233501, 56.1284943), (12.6234598, 56.1284676), (12.6235533, 56.1284218), (12.6236343, 56.1283531), (12.6236744, 56.1282616), (12.6236801, 56.1280975), (12.6236877, 56.1279793), (12.6237154, 56.1279373), (12.6237631, 56.1278915), (12.6238565, 56.1278534), (12.6239777, 56.1278343), (12.6242332, 56.1278152), (12.6243715, 56.1277885), (12.6244717, 56.1277428), (12.6246538, 56.1275826), (12.6248417, 56.1273918), (12.6250563, 56.127079), (12.6251469, 56.1268921), (12.6251688, 56.1267586), (12.6251612, 56.1266823), (12.6251259, 56.1266365), (12.6249723, 56.1265678), (12.6252127, 56.1266365), (12.6253538, 56.126667), (12.6251364, 56.127182)]
}

/// The shore with a scaled copy of itself pinched onto vertex `anchor`, moved to
/// `(dx, dy)` — the same construction the unit guard pins, so one adversarial ring can
/// be stamped across the table at many scales and places without any two overlapping.
fn lobed(scale: f64, anchor: usize, dx: f64, dy: f64) -> Vec<(f64, f64)> {
    let a = shore();
    let anchor = anchor % (a.len() - 1);
    let (px, py) = a[anchor];
    let mut v: Vec<(f64, f64)> = a[..=anchor].to_vec();
    for &(x, y) in a.iter().skip(1) {
        v.push((scale.mul_add(x - px, px), scale.mul_add(y - py, py)));
    }
    v.extend_from_slice(&a[anchor + 1..]);
    v.iter().map(|&(x, y)| (x + dx, y + dy)).collect()
}

/// A rippled ring of `n` vertices, `r(θ) = r0 + amp·sin(mθ)` — a big shore that DP
/// KEEPS, so the pack really has to run the crossing test on tens of thousands of
/// segments. Star-shaped, so the source is simple.
fn rippled(n: usize, cx: f64, cy: f64, r0: f64, amp: f64, m: f64) -> Vec<(f64, f64)> {
    let mut v = Vec::with_capacity(n + 1);
    for i in 0..n {
        #[allow(clippy::cast_precision_loss)]
        let t = (i as f64) / (n as f64) * std::f64::consts::TAU;
        let r = amp.mul_add((m * t).sin(), r0);
        v.push((r.mul_add(t.cos(), cx), r.mul_add(t.sin(), cy)));
    }
    v.push(v[0]);
    v
}

fn wkb_polygon(ring: &[(f64, f64)]) -> Vec<u8> {
    let mut v = vec![1u8];
    v.extend_from_slice(&3u32.to_le_bytes()); // Polygon
    v.extend_from_slice(&1u32.to_le_bytes()); // one ring
    #[allow(clippy::cast_possible_truncation)]
    v.extend_from_slice(&(ring.len() as u32).to_le_bytes());
    for &(x, y) in ring {
        v.extend_from_slice(&x.to_le_bytes());
        v.extend_from_slice(&y.to_le_bytes());
    }
    v
}

/// **The table the pack is asked to swallow.** 320 adversarial rings, spread over
/// Sweden so no two touch: 256 small ones at four scales of the shore motif (the shapes
/// that make Douglas–Peucker bowtie) and 64 BIG ones of 12 000 vertices each (the shapes
/// that made the uncapped crossing test cost seconds).
fn rings() -> Vec<Vec<(f64, f64)>> {
    let mut out = Vec::with_capacity(320);
    let mut i = 0usize;
    let mut tried = 0usize;
    while out.len() < 256 {
        tried += 1;
        assert!(tried < 20_000, "the lobed-shore generator ran dry");
        let scale = [-0.0625, -0.125, -0.25, -0.5][i % 4];
        let r = lobed(scale, (i * 7) % 71, 0.0, 0.0);
        i += 1;
        // The SOURCE rings must be clean, or the base column carries the bowties and the
        // census cannot tell the fixture's defect from the pack's. Asked of the ONE
        // writer of the predicate, never of a copy of it (LAW 5).
        if osm_katana::verify::ring_self_intersects(&r) {
            continue;
        }
        let k = out.len();
        #[allow(clippy::cast_precision_loss)]
        let (dx, dy) = ((k % 16) as f64 * 0.35, (k / 16) as f64 * 0.35);
        out.push(r.iter().map(|&(x, y)| (x + dx, y + dy)).collect());
    }
    for j in 0..64usize {
        #[allow(clippy::cast_precision_loss)]
        let (cx, cy) = (11.0 + (j % 8) as f64 * 0.9, 63.0 + (j / 8) as f64 * 0.5);
        out.push(rippled(12_000, cx, cy, 0.2, 4.0e-4, 3_000.0));
    }
    out
}

/// The fixture's own precondition, as an arm: every SOURCE ring is clean, so a crossing
/// anywhere in the packed output came from the pack.
#[test]
fn every_source_ring_of_the_fixture_is_clean() {
    let rs = rings();
    assert_eq!(rs.len(), 320);
    let bad = rs.iter().filter(|r| osm_katana::verify::ring_self_intersects(r)).count();
    assert_eq!(bad, 0, "{bad} of the fixture's own rings cross themselves");
}

fn write_source(path: &Path, rings: &[Vec<(f64, f64)>]) {
    let fields = vec![
        Field::new("id", DataType::Int64, false),
        Field::new("geometry", DataType::Binary, true),
        Field::new("tags", DataType::Utf8, true),
    ];
    let mut meta = HashMap::new();
    meta.insert(String::from(metadata::GEO_KEY), GeoMeta::new(&["Polygon"]).to_json());
    let schema = Arc::new(Schema::new_with_metadata(fields, meta));
    #[allow(clippy::cast_possible_wrap)]
    let ids: Vec<i64> = (0..rings.len() as i64).collect();
    let geoms: Vec<Vec<u8>> = rings.iter().map(|r| wkb_polygon(r)).collect();
    let tags: Vec<&str> = vec![TAGS; rings.len()];
    let cols: Vec<Arc<dyn Array>> = vec![
        Arc::new(Int64Array::from(ids)),
        Arc::new(BinaryArray::from_iter_values(geoms.iter().map(Vec::as_slice))),
        Arc::new(StringArray::from(tags)),
    ];
    let f = std::fs::File::create(path).unwrap();
    let mut w = parquet::arrow::ArrowWriter::try_new(f, schema.clone(), None).unwrap();
    w.write(&RecordBatch::try_new(schema, cols).unwrap()).unwrap();
    w.close().unwrap();
}

fn pack_fixture(stem: &str) -> PathBuf {
    let src = scratch(&format!("{stem}-src"));
    let dst = scratch(stem);
    write_source(&src, &rings());
    let opts = PackOptions { row_group_rows: 64, ..PackOptions::default() };
    spatial_pack(&src, &dst, &opts).unwrap_or_else(|e| panic!("spatial_pack: {e:#}"));
    dst
}

/// **No rung column of the packed table carries a self-intersecting ring.**
///
/// This is the census `verify --geometry-quality` runs, on the columns the pack wrote.
/// `self_intersecting` must be 0 — and the test prints `big_ring_unchecked` alongside
/// it, because "the census did not look" is not "the census passed it".
#[test]
fn a_packed_rung_column_carries_no_self_intersecting_ring() {
    let dst = pack_fixture("no-bowties");
    let qs = osm_katana::verify::geometry_quality(&dst).expect("census");
    assert!(!qs.is_empty(), "the census found no geometry column to read");
    let mut saw_rung = false;
    for q in &qs {
        println!(
            "{}: rows {} polys {} selfx {} big_unchecked {} zero_area {} unclosed {}",
            q.column, q.rows, q.polys, q.self_intersecting, q.big_ring_unchecked,
            q.zero_area, q.unclosed_ring,
        );
        assert!(!q.fault(), "the census is itself RED on {}: {q:?}", q.column);
        assert_eq!(
            q.self_intersecting, 0,
            "{} carries {} self-intersecting rings — every one of them tessellates as \
             the sea-wedge defect",
            q.column, q.self_intersecting,
        );
        if q.column.starts_with("geometry_z") {
            saw_rung = true;
            assert!(q.polys > 0, "{} decoded no polygon at all", q.column);
        }
    }
    assert!(saw_rung, "the pack wrote no generalisation column — nothing was guarded");
}

/// **The rung columns are still a SIMPLIFICATION.** The give-up path ships the source
/// ring, and a bound on the crossing test makes it fire more often; a rung column that
/// is not smaller than its base is the price being paid for the cap, and it is a number
/// this test reads off disk rather than a hope.
#[test]
fn the_rungs_are_smaller_than_the_base() {
    let dst = pack_fixture("rungs-smaller");
    let f = std::fs::File::open(&dst).unwrap();
    let rdr = ParquetRecordBatchReaderBuilder::try_new(f).unwrap().build().unwrap();
    let mut base = 0usize;
    let mut rung: Vec<(u8, usize, usize)> =
        GENERALISED_RUNGS.iter().map(|&z| (z, 0, 0)).collect();
    for b in rdr {
        let b = b.unwrap();
        let g = b
            .column_by_name("geometry")
            .expect("geometry")
            .as_any()
            .downcast_ref::<BinaryArray>()
            .unwrap();
        base += (0..g.len()).map(|i| g.value(i).len()).sum::<usize>();
        for (z, bytes, present) in &mut rung {
            let name = generalised_geometry_name(*z);
            let Some(c) = b.column_by_name(&name) else { continue };
            let c = c.as_any().downcast_ref::<BinaryArray>().unwrap();
            for i in 0..c.len() {
                if !c.is_null(i) {
                    *bytes += c.value(i).len();
                    *present += 1;
                }
            }
        }
    }
    assert!(base > 0);
    for (z, bytes, present) in &rung {
        #[allow(clippy::cast_precision_loss)]
        let share = *bytes as f64 / base as f64;
        println!("geometry_z{z}: {present} cells, {bytes} bytes, {:.1} % of base", share * 100.0);
        assert!(
            *present > 0,
            "geometry_z{z} is empty — the rung withheld every row, which is not a \
             simplification, it is a deletion",
        );
        assert!(
            share < 1.0,
            "geometry_z{z} is {:.1} % of the base column — the pack is shipping source \
             rings so often that the rung is not a simplification at all",
            share * 100.0,
        );
    }
}