nordisk/facett

facett-map3d-flight/benches/heavy_terrain.rs

//! **HEAVY terrain-DEM benchmark** — the real-elevation terrain pipeline driven at
//! full cores + max memory over a large tile (the Bencher LAW: bench the PARTS + the
//! end-to-end build, with element/triangle throughput).
//!
//! Steps benched, each with a throughput:
//! - `hgt_load`          — parse a 1-arc-second SRTM `.hgt` tile (3601² big-endian i16)
//! - `dem_sample`        — bilinear `elevation_at` throughput over a random point cloud
//! - `procedural_sample` — the synthetic-fallback heightfield sampling throughput
//! - `build_mesh_dem`    — build the 3D terrain mesh (verts+normals+indices) from a DEM
//! - `build_mesh_synth`  — the same over the procedural fallback (no DEM in memory)
//!
//! Grid side is `MAP3D_FLIGHT_BENCH_GRID` (default 1024 → 1024² ≈ 1.05 M vertices,
//! ~2.1 M triangles per build). Run:
//!   CARGO_TARGET_DIR=/mnt/work4t/cargo-target/facett-map3d-flight \
//!   cargo bench -p facett-map3d-flight

use std::hint::black_box;

use criterion::{Criterion, Throughput, criterion_group, criterion_main};

use facett_map3d_flight::origin::FloatingOrigin;
use facett_map3d_flight::terrain::{Dem, Elevation, GeoRegion, ProceduralTerrain, build_terrain_mesh};
use facett_map3d_flight::WorldPos;

/// Terrain vertex-grid side (`n×n` vertices). `MAP3D_FLIGHT_BENCH_GRID` overrides it.
fn grid_side() -> usize {
    std::env::var("MAP3D_FLIGHT_BENCH_GRID").ok().and_then(|v| v.parse().ok()).unwrap_or(1024)
}

/// Synthesise a 1-arc-second SRTM `.hgt` byte buffer (`n×n` big-endian i16) with a
/// deterministic ridged heightfield — the real on-disk format the loader parses.
fn hgt_bytes(n: usize) -> Vec<u8> {
    let mut bytes = Vec::with_capacity(n * n * 2);
    for r in 0..n {
        for c in 0..n {
            let fx = c as f64 / n as f64;
            let fy = r as f64 / n as f64;
            let h = (900.0 + 800.0 * (fx * 21.0).sin() * (fy * 17.0).cos()) as i16;
            bytes.extend_from_slice(&h.to_be_bytes());
        }
    }
    bytes
}

/// A large real DEM (synthetic heights, real grid) over a 1° Sweden cell.
fn heavy_dem(n: usize) -> Dem {
    Dem::synthetic(17.0, 60.0, 1.0, 1.0, n, n, |lon, lat| {
        600.0 + 500.0 * (lon * 9.0).sin() as f32 * (lat * 7.0).cos() as f32 + 120.0 * (lon * 31.0).cos() as f32
    })
}

fn bench_terrain(c: &mut Criterion) {
    let side = grid_side();
    let region = GeoRegion::new(17.1, 59.1, 17.9, 59.9);
    let origin = FloatingOrigin::snapped(WorldPos::from_lonlat(17.5, 59.5, 4000.0), 512.0);

    // ── hgt_load: parse a real 1-arc-second tile (3601²) ──────────────────────────
    {
        let n = 3601usize; // SRTM 1-arc-second
        let bytes = hgt_bytes(n);
        let mut g = c.benchmark_group("hgt_load");
        g.throughput(Throughput::Elements((n * n) as u64));
        g.bench_function("srtm_1arcsec_3601", |b| {
            b.iter(|| {
                let dem = Dem::from_hgt_bytes(black_box(&bytes), 17.0, 59.0).unwrap();
                black_box(dem.void_count)
            });
        });
        g.finish();
    }

    // ── dem_sample / procedural_sample: elevation query throughput ────────────────
    {
        let dem = heavy_dem(side);
        let proc = ProceduralTerrain::default();
        // A deterministic point cloud inside the region.
        let pts: Vec<(f64, f64)> = (0..1_000_000)
            .map(|i| {
                let fx = ((i * 2654435761u64 as usize) % 100_000) as f64 / 100_000.0;
                let fy = ((i * 40503) % 100_000) as f64 / 100_000.0;
                (17.0 + fx, 59.0 + fy)
            })
            .collect();

        let mut g = c.benchmark_group("dem_sample");
        g.throughput(Throughput::Elements(pts.len() as u64));
        g.bench_function("bilinear_1M", |b| {
            b.iter(|| {
                let mut acc = 0.0f32;
                for &(lon, lat) in black_box(&pts) {
                    acc += dem.elevation_at(lon, lat);
                }
                black_box(acc)
            });
        });
        g.finish();

        let mut g = c.benchmark_group("procedural_sample");
        g.throughput(Throughput::Elements(pts.len() as u64));
        g.bench_function("fbm_1M", |b| {
            b.iter(|| {
                let mut acc = 0.0f32;
                for &(lon, lat) in black_box(&pts) {
                    acc += proc.elevation_at(lon, lat);
                }
                black_box(acc)
            });
        });
        g.finish();
    }

    // ── build_mesh_dem / build_mesh_synth: the full terrain mesh build ────────────
    {
        let dem = heavy_dem(side);
        let proc = ProceduralTerrain::default();
        let verts = (side * side) as u64;
        let tris = ((side - 1) * (side - 1) * 2) as u64;

        let mut g = c.benchmark_group("build_mesh_dem");
        g.throughput(Throughput::Elements(tris));
        g.bench_function(format!("dem_{side}x{side}"), |b| {
            b.iter(|| {
                let m = build_terrain_mesh(black_box(&dem), region, side, side, origin, 1.0);
                black_box((m.vertex_count(), m.triangle_count()))
            });
        });
        g.finish();

        let mut g = c.benchmark_group("build_mesh_synth");
        g.throughput(Throughput::Elements(tris));
        g.bench_function(format!("synth_{side}x{side}"), |b| {
            b.iter(|| {
                let m = build_terrain_mesh(black_box(&proc), region, side, side, origin, 1.0);
                black_box(m.vertex_count())
            });
        });
        g.finish();

        eprintln!("[heavy_terrain] grid {side}×{side} → {verts} vertices / {tris} triangles per build");
    }
}

criterion_group!(benches, bench_terrain);
criterion_main!(benches);