diff --git a/Cargo.lock b/Cargo.lock index 605f47e..381a724 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -8,6 +8,12 @@ version = "2.0.1" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "320119579fcad9c21884f5c4861d16174d0e06250625266f50fe6898340abefa" +[[package]] +name = "allocator-api2" +version = "0.2.21" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "683d7910e743518b0e34f1186f92494becacb047c7b6bf616c96772180fef923" + [[package]] name = "android_system_properties" version = "0.1.5" @@ -230,7 +236,9 @@ dependencies = [ "clap", "font8x8", "plotters", + "rayon", "reqwest", + "spade", ] [[package]] @@ -260,6 +268,31 @@ dependencies = [ "cfg-if", ] +[[package]] +name = "crossbeam-deque" +version = "0.8.6" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9dd111b7b7f7d55b72c0a6ae361660ee5853c9af73f70c3c2ef6858b950e2e51" +dependencies = [ + "crossbeam-epoch", + "crossbeam-utils", +] + +[[package]] +name = "crossbeam-epoch" +version = "0.9.18" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "5b82ac4a3c2ca9c3460964f020e1402edd5753411d7737aa39c3714ad1b5420e" +dependencies = [ + "crossbeam-utils", +] + +[[package]] +name = "crossbeam-utils" +version = "0.8.21" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d0a5c400df2834b80a4c3327b3aad3a4c4cd4de0629063962b03235697506a28" + [[package]] name = "displaydoc" version = "0.2.5" @@ -271,6 +304,18 @@ dependencies = [ "syn", ] +[[package]] +name = "either" +version = "1.15.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "48c757948c5ede0e46177b7add2e67155f70e33c07fea8284df6576da70b3719" + +[[package]] +name = "equivalent" +version = "1.0.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "877a4ace8713b0bcf2a4e7eec82529c029f1d0619886d18145fea96c3ffe5c0f" + [[package]] name = "fdeflate" version = "0.3.7" @@ -296,6 +341,12 @@ dependencies = [ "miniz_oxide", ] +[[package]] +name = "foldhash" +version = "0.1.5" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d9c4f5dac5e15c24eb999c26181a6ca40b39fe946cbe4c263c7209467bc83af2" + [[package]] name = "font8x8" version = "0.2.7" @@ -387,6 +438,17 @@ dependencies = [ "wasm-bindgen", ] +[[package]] +name = "hashbrown" +version = "0.15.5" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9229cfe53dfd69f0609a49f65461bd93001ea1ef889cd5529dd176593f5338a1" +dependencies = [ + "allocator-api2", + "equivalent", + "foldhash", +] + [[package]] name = "heck" version = "0.5.0" @@ -941,6 +1003,26 @@ dependencies = [ "getrandom 0.3.4", ] +[[package]] +name = "rayon" +version = "1.11.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "368f01d005bf8fd9b1206fb6fa653e6c4a81ceb1466406b81792d87c5677a58f" +dependencies = [ + "either", + "rayon-core", +] + +[[package]] +name = "rayon-core" +version = "1.13.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "22e18b0f0062d30d4230b2e85ff77fdfe4326feb054b9783a3460d8435c8ab91" +dependencies = [ + "crossbeam-deque", + "crossbeam-utils", +] + [[package]] name = "reqwest" version = "0.12.28" @@ -995,6 +1077,12 @@ dependencies = [ "windows-sys 0.52.0", ] +[[package]] +name = "robust" +version = "1.2.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "4e27ee8bb91ca0adcf0ecb116293afa12d393f9c2b9b9cd54d33e8078fe19839" + [[package]] name = "rustc-hash" version = "2.1.1" @@ -1137,6 +1225,18 @@ dependencies = [ "windows-sys 0.60.2", ] +[[package]] +name = "spade" +version = "2.15.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "fb313e1c8afee5b5647e00ee0fe6855e3d529eb863a0fdae1d60006c4d1e9990" +dependencies = [ + "hashbrown", + "num-traits", + "robust", + "smallvec", +] + [[package]] name = "stable_deref_trait" version = "1.2.1" diff --git a/Cargo.toml b/Cargo.toml index 3cbe801..6da5bfb 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -15,3 +15,5 @@ font8x8 = "0.2" clap = { version = "4", features = ["derive"] } chrono = { version = "0.4", default-features = false, features = ["std", "clock"] } anyhow = "1" +spade = "2" +rayon = "1" diff --git a/src/render.rs b/src/render.rs index 01c947f..19331b8 100644 --- a/src/render.rs +++ b/src/render.rs @@ -2,7 +2,8 @@ use crate::projection::OrthoProjection; use anyhow::{Context, Result}; use font8x8::{UnicodeFonts, BASIC_FONTS}; use plotters::prelude::*; -use std::collections::VecDeque; +use rayon::prelude::*; +use spade::{DelaunayTriangulation, HasPosition, Point2, Triangulation}; use std::path::Path; const IMG_WIDTH: u32 = 900; @@ -14,20 +15,8 @@ const CENTER_LON: f64 = 13.08; const HALF_W_DEG: f64 = 0.8; // horizontal half-extent in degrees of arc const HALF_H_DEG: f64 = 0.4; // vertical half-extent -/// Cloud-cover thresholds and their associated RGB colours. -/// Values below the first threshold are drawn as the background colour. -const LEVELS: [(f32, (u8, u8, u8)); 7] = [ - (1.0, (208, 232, 248)), // very thin - (2.0, (160, 200, 232)), - (5.0, (112, 168, 216)), - (10.0, ( 64, 128, 192)), - (20.0, ( 32, 80, 160)), - (50.0, ( 16, 32, 96)), - (100.0, ( 5, 10, 48)), -]; - const BACKGROUND: RGBColor = RGBColor(240, 248, 255); // pale sky-blue -const LAND_BG: RGBColor = RGBColor(232, 236, 220); // muted green for land +const LAND_BG: RGBColor = RGBColor(142, 178, 35); // strong green for land const BLACK: RGBColor = RGBColor( 0, 0, 0); const RED: RGBColor = RGBColor(220, 30, 30); @@ -42,16 +31,28 @@ const LOCATIONS: &[(&str, f64, f64)] = &[ ("Berlin (Mitte)", 13.39, 52.52), ]; +#[derive(Clone, Copy)] +struct CloudVertex { + position: Point2, + cloud: f32, +} + +impl HasPosition for CloudVertex { + type Scalar = f64; + fn position(&self) -> Point2 { + self.position + } +} + +/// Map cloud cover (0–100 %) to a colour by blending the land background toward +/// white. Returns `None` for near-zero cover so the land background shows through. fn cloud_color(cover: f32) -> Option { if cover < 1.0 { - return None; // clear sky — use background + return None; } - for &(threshold, (r, g, b)) in LEVELS.iter().rev() { - if cover >= threshold { - return Some(RGBColor(r, g, b)); - } - } - Some(RGBColor(LEVELS[0].1 .0, LEVELS[0].1 .1, LEVELS[0].1 .2)) + let t = (cover / 100.0).clamp(0.0, 1.0); + let blend = |base: u8| -> u8 { (base as f32 + t * (255.0 - base as f32)).round() as u8 }; + Some(RGBColor(blend(LAND_BG.0), blend(LAND_BG.1), blend(LAND_BG.2))) } /// Draw ASCII text using the built-in 8×8 pixel font. @@ -91,6 +92,57 @@ fn draw_pixel_text( } } +/// Separable Gaussian blur on a flat f32 grid. +/// +/// NaN pixels are treated as transparent: they contribute zero weight and the +/// result is renormalised over non-NaN neighbours only. Pixels with no non-NaN +/// neighbours remain NaN. +fn gaussian_blur(grid: &[f32], w: usize, h: usize, sigma: f32) -> Vec { + let radius = (3.0 * sigma).ceil() as usize; + let kernel: Vec = (0..=2 * radius) + .map(|i| { + let x = i as f32 - radius as f32; + (-x * x / (2.0 * sigma * sigma)).exp() + }) + .collect(); + + // Horizontal pass + let temp: Vec = (0..h * w) + .into_par_iter() + .map(|idx| { + let row = idx / w; + let col = idx % w; + let mut val = 0.0f32; + let mut wsum = 0.0f32; + for (ki, &k) in kernel.iter().enumerate() { + let c = col as i64 + ki as i64 - radius as i64; + if c < 0 || c >= w as i64 { continue; } + let v = grid[row * w + c as usize]; + if !v.is_nan() { val += k * v; wsum += k; } + } + if wsum > 0.0 { val / wsum } else { f32::NAN } + }) + .collect(); + + // Vertical pass + (0..h * w) + .into_par_iter() + .map(|idx| { + let row = idx / w; + let col = idx % w; + let mut val = 0.0f32; + let mut wsum = 0.0f32; + for (ki, &k) in kernel.iter().enumerate() { + let r = row as i64 + ki as i64 - radius as i64; + if r < 0 || r >= h as i64 { continue; } + let v = temp[r as usize * w + col]; + if !v.is_nan() { val += k * v; wsum += k; } + } + if wsum > 0.0 { val / wsum } else { f32::NAN } + }) + .collect() +} + /// Render one forecast frame as a PNG. /// /// - `lats`, `lons`, `cloud_cover` must all have the same length. @@ -119,54 +171,89 @@ pub fn render_frame( map_area.fill(&LAND_BG).context("fill map area")?; - // Nearest-neighbour fill via multi-source BFS. - // - // Each data point seeds one pixel in a flat grid. BFS then propagates each - // value outward to all 4-connected NaN neighbours, so every reachable pixel - // ends up holding the value of the closest source point (by Manhattan distance). - // Complexity: O(map_w × map_h) — ~410 k array ops, no hashing. - let grid_w = map_w as usize; - let grid_h = map_h as usize; - let mut grid = vec![f32::NAN; grid_w * grid_h]; - let mut queue: VecDeque<(i32, i32)> = VecDeque::new(); - - // Seed the grid with projected data points. + // Build a Delaunay triangulation in continuous pixel-space coordinates. + // Include all grid points that project within a small margin outside the image + // bounds so that boundary pixels interpolate cleanly without edge artefacts. + const TRI_MARGIN: f64 = 10.0; + let mut tri = DelaunayTriangulation::::new(); for i in 0..lats.len() { let cover = cloud_cover[i]; if cover.is_nan() { continue; } let Some((px, py)) = proj.project(lats[i] as f64, lons[i] as f64) else { continue }; - let Some((col, row)) = proj.to_pixel(px, py, map_w, map_h) else { continue }; - let idx = row as usize * grid_w + col as usize; - if grid[idx].is_nan() { - grid[idx] = cover; - queue.push_back((col, row)); + // Continuous pixel coordinates — no rounding. + let col_f = (px + HALF_W_DEG) / (2.0 * HALF_W_DEG) * map_w as f64; + let row_f = (-py + HALF_H_DEG) / (2.0 * HALF_H_DEG) * map_h as f64; + if col_f < -TRI_MARGIN + || col_f > map_w as f64 + TRI_MARGIN + || row_f < -TRI_MARGIN + || row_f > map_h as f64 + TRI_MARGIN + { + continue; } + let _ = tri.insert(CloudVertex { position: Point2::new(col_f, row_f), cloud: cover }); } - // BFS flood fill. - const DIRS: [(i32, i32); 4] = [(-1, 0), (1, 0), (0, -1), (0, 1)]; - while let Some((col, row)) = queue.pop_front() { - let cover = grid[row as usize * grid_w + col as usize]; - for (dc, dr) in DIRS { - let nc = col + dc; - let nr = row + dr; - if nc < 0 || nc >= map_w as i32 || nr < 0 || nr >= map_h as i32 { - continue; + // For each output pixel, locate the containing triangle and barycentric-interpolate + // the cloud cover from its three vertices. The triangulation is immutable at this + // point, so the work is embarrassingly parallel via rayon. + let grid_w = map_w as usize; + let grid_h = map_h as usize; + let cover_grid: Vec = (0..grid_h * grid_w) + .into_par_iter() + .map(|i| { + let col = (i % grid_w) as f64 + 0.5; + let row = (i / grid_w) as f64 + 0.5; + let p = Point2::new(col, row); + use spade::PositionInTriangulation::*; + match tri.locate(p) { + OnVertex(v) => tri.vertex(v).data().cloud, + OnEdge(e) => { + // Linear interpolation along the edge. + let edge = tri.directed_edge(e); + let a = edge.from(); + let b = edge.to(); + let (ax, ay) = (a.position().x, a.position().y); + let (bx, by) = (b.position().x, b.position().y); + let ab2 = (bx - ax).powi(2) + (by - ay).powi(2); + let t = if ab2 == 0.0 { + 0.0 + } else { + (((p.x - ax) * (bx - ax) + (p.y - ay) * (by - ay)) / ab2) + .clamp(0.0, 1.0) + }; + a.data().cloud * (1.0 - t) as f32 + b.data().cloud * t as f32 + } + OnFace(f) => { + let face = tri.face(f); + let [v0, v1, v2] = face.vertices(); + let (ax, ay) = (v0.position().x, v0.position().y); + let (bx, by) = (v1.position().x, v1.position().y); + let (cx, cy) = (v2.position().x, v2.position().y); + let denom = (by - cy) * (ax - cx) + (cx - bx) * (ay - cy); + if denom == 0.0 { + return v0.data().cloud; + } + let l0 = ((by - cy) * (p.x - cx) + (cx - bx) * (p.y - cy)) / denom; + let l1 = ((cy - ay) * (p.x - cx) + (ax - cx) * (p.y - cy)) / denom; + let l2 = 1.0 - l0 - l1; + v0.data().cloud * l0 as f32 + + v1.data().cloud * l1 as f32 + + v2.data().cloud * l2 as f32 + } + OutsideOfConvexHull(_) | NoTriangulation => f32::NAN, } - let nidx = nr as usize * grid_w + nc as usize; - if grid[nidx].is_nan() { - grid[nidx] = cover; - queue.push_back((nc, nr)); - } - } - } + }) + .collect(); - // Paint the filled grid. + const BLUR_SIGMA: f32 = 8.0; + let cover_grid = gaussian_blur(&cover_grid, grid_w, grid_h, BLUR_SIGMA); + + // Paint the interpolated grid. for row in 0..map_h as i32 { for col in 0..map_w as i32 { - let cover = grid[row as usize * grid_w + col as usize]; + let cover = cover_grid[row as usize * grid_w + col as usize]; if let Some(color) = cloud_color(cover) { map_area.draw_pixel((col, row), &color).ok(); } @@ -231,22 +318,17 @@ fn draw_legend( ); }; - // "< 1%" row (clear sky) - draw_entry("< 1%", LAND_BG, y as i32 + 22); - - // Remaining threshold rows — use ASCII hyphen instead of en-dash let entries = [ - (1.0f32, "1-2%"), - (2.0, "2-5%"), - (5.0, "5-10%"), - (10.0, "10-20%"), - (20.0, "20-50%"), - (50.0, "50-100%"), + (0.0f32, "< 1%"), + (25.0, "~25%"), + (50.0, "~50%"), + (75.0, "~75%"), + (100.0, "100%"), ]; - for (i, &(lo, label)) in entries.iter().enumerate() { - let row_y = y as i32 + 22 + (i as i32 + 1) * (box_h + 2); - let color = cloud_color(lo).unwrap_or(BACKGROUND); + for (i, &(cover, label)) in entries.iter().enumerate() { + let row_y = y as i32 + 22 + i as i32 * (box_h + 2); + let color = cloud_color(cover).unwrap_or(LAND_BG); draw_entry(label, color, row_y); }