diff --git a/crates/render/src/camera.rs b/crates/render/src/camera.rs index 8b0117f..3e33bc8 100644 --- a/crates/render/src/camera.rs +++ b/crates/render/src/camera.rs @@ -5,6 +5,7 @@ use serde::{Deserialize, Serialize}; +use crate::pick::Projector; use crate::{Aabb3, Vec3}; /// Column-major 4×4 matrix. @@ -206,13 +207,9 @@ impl Camera { self.proj(aspect, n, f).mul(&self.view()) } /// World point → pixel coordinates (x right, y down) and depth (0 near … 1 far). + /// For many points, a [`Projector`] computes the projection once. pub fn to_screen(&self, p: Vec3, w: f64, h: f64, radius: f64) -> Option<(f64, f64, f64)> { - let c = self.view_proj(w / h.max(1.0), radius).apply(p); - if c[3].abs() < 1e-12 { - return None; - } - let (x, y, z) = (c[0] / c[3], c[1] / c[3], c[2] / c[3]); - (c[3] > 0.0).then_some(((x + 1.0) * 0.5 * w, (1.0 - y) * 0.5 * h, z)) + Projector::new(self, w, h, radius).screen(p).map(|[x, y, z]| (x, y, z)) } /// Picking ray through pixel (x, y): origin and unit direction. pub fn ray(&self, x: f64, y: f64, w: f64, h: f64) -> (Vec3, Vec3) { diff --git a/crates/render/src/display.rs b/crates/render/src/display.rs index c71130f..2f406ae 100644 --- a/crates/render/src/display.rs +++ b/crates/render/src/display.rs @@ -14,6 +14,7 @@ use fieldcraft_model::{ElementKind, Model}; +use crate::pick::{Hit, first_hit}; use crate::{Aabb3, Vec3}; /// One drawn triangle of an element face. @@ -25,6 +26,18 @@ pub struct DisplayTri { /// The face's unit normal, pointing out of its element (a shell's own orientation). pub normal: Vec3, pub kind: ElementKind, + /// The element it is a face of, and the face's number (from 1, as + /// [`ElementKind::faces`] numbers them; a shell's is 1, its top side). + pub element: u32, + pub face: u8, +} + +impl DisplayTri { + /// The element face seen looking along `dir`: a shell's top side (1) from in front and its + /// bottom side (2) from behind; a solid's face as it is. + pub fn seen_face(&self, dir: Vec3) -> u8 { + if self.kind.dimension() == 2 && self.normal.dot(dir) > 0.0 { 2 } else { self.face } + } } /// An element edge. @@ -41,6 +54,7 @@ pub struct DisplayLine { pub b: Vec3, pub nodes: [u32; 2], pub kind: ElementKind, + pub element: u32, } /// A node. @@ -57,6 +71,8 @@ pub struct SurfaceTri { pub normal: Vec3, /// Which surface (body) it belongs to: the order [`DisplayMesh::add_surface`] got them in. pub body: u32, + /// The face of that surface it lies on (its `tri_face` tag; `u32::MAX` when it has none). + pub face: u32, } /// Display data for a model and its geometry. @@ -80,7 +96,9 @@ pub struct DisplayMesh { } /// The faces of an element as corner indices (0-based into its node list), each wound -/// counter-clockwise seen from outside for a positively oriented element. +/// counter-clockwise seen from outside for a positively oriented element. A solid's faces come in +/// the order [`ElementKind::faces`] numbers them (with the same corners); a shell draws its top +/// side, face 1. fn faces(kind: ElementKind) -> &'static [&'static [usize]] { match kind { ElementKind::Rod2 | ElementKind::Rod3 | ElementKind::Beam2 => &[], @@ -142,7 +160,7 @@ impl DisplayMesh { inner.dedup(); drop(keys); let mut edges: Vec<(u32, u32)> = Vec::new(); - for e in model.elements().values() { + for (&label, e) in model.elements() { // Corners: one past the highest face index (line elements: their two ends). let corners = faces(e.kind).iter().flat_map(|f| f.iter()).max().map_or(2, |m| m + 1).min(MAX_CORNERS); let mut at = [Vec3::ZERO; MAX_CORNERS]; @@ -159,13 +177,14 @@ impl DisplayMesh { let at = at.get(..corners).unwrap_or(&[]); if e.kind.dimension() == 1 { if let (Some(&a), Some(&b), Some(&na), Some(&nb)) = (at.first(), at.get(1), e.nodes.first(), e.nodes.get(1)) { - out.lines.push(DisplayLine { a, b, nodes: [na, nb], kind: e.kind }); + out.lines.push(DisplayLine { a, b, nodes: [na, nb], kind: e.kind, element: label }); } continue; } let solid = e.kind.dimension() == 3; let element_centre = centre(at); - for f in faces(e.kind) { + for (number, f) in faces(e.kind).iter().enumerate() { + let face = u8::try_from(number + 1).unwrap_or(u8::MAX); if solid && face_key(&e.nodes, f).is_some_and(|k| inner.binary_search(&k).is_ok()) { continue; } @@ -179,7 +198,7 @@ impl DisplayMesh { if let (Some(&a), Some(&b), Some(&c), Some(&la), Some(&lb), Some(&lc)) = (p.first(), p.get(k), p.get(k + 1), labels.first(), labels.get(k), labels.get(k + 1)) { - out.tris.push(DisplayTri { p: [a, b, c], nodes: [la, lb, lc], normal, kind: e.kind }); + out.tris.push(DisplayTri { p: [a, b, c], nodes: [la, lb, lc], normal, kind: e.kind, element: label, face }); } } for (j, &i) in f.iter().enumerate() { @@ -229,8 +248,8 @@ impl DisplayMesh { for p in [a, b, c] { self.bounds.add(p); } - self.surface_tris.push(SurfaceTri { p: [a, b, c], normal, body: self.bodies }); let face = tri_face.get(t).copied().unwrap_or(u32::MAX); + self.surface_tris.push(SurfaceTri { p: [a, b, c], normal, body: self.bodies, face }); let [v0, v1, v2] = [vertex(i0), vertex(i1), vertex(i2)]; for (u, v) in [(v0, v1), (v1, v2), (v2, v0)] { sides.push((u.min(v), u.max(v), face)); @@ -263,30 +282,24 @@ impl DisplayMesh { self.nodes.binary_search_by_key(&label, |n| n.label).ok().and_then(|i| self.nodes.get(i)).map(|n| n.p) } + /// The element face triangle ([`DisplayMesh::tris`]) a ray (`origin`, `dir`) hits first + /// ([`pick::first_hit`]). + pub fn hit_element_face(&self, origin: Vec3, dir: Vec3) -> Option { + first_hit(origin, dir, self.tris.iter().map(|t| t.p)) + } + + /// The geometry surface triangle ([`DisplayMesh::surface_tris`]) a ray hits first. + pub fn hit_surface(&self, origin: Vec3, dir: Vec3) -> Option { + first_hit(origin, dir, self.surface_tris.iter().map(|t| t.p)) + } + /// The node nearest to where a ray (`origin`, `dir`) first hits an element face ahead of /// its origin: the hit face's closest corner, and the hit point. Line elements and geometry /// are not hit. pub fn pick(&self, origin: Vec3, dir: Vec3) -> Option<(u32, Vec3)> { - let mut best: Option<(f64, &DisplayTri)> = None; - for t in &self.tris { - // Möller–Trumbore, both sides. - let (e1, e2) = (t.p[1] - t.p[0], t.p[2] - t.p[0]); - let h = dir.cross(e2); - let det = e1.dot(h); - if !(det.abs() > 1e-300) { - continue; - } - let s = origin - t.p[0]; - let u = s.dot(h) / det; - let q = s.cross(e1); - let v = dir.dot(q) / det; - let d = e2.dot(q) / det; - if u >= 0.0 && v >= 0.0 && u + v <= 1.0 && d >= 0.0 && d.is_finite() && best.is_none_or(|(b, _)| d < b) { - best = Some((d, t)); - } - } - let (d, t) = best?; - let hit = origin + dir * d; + let h = self.hit_element_face(origin, dir)?; + let t = self.tris.get(h.index)?; + let hit = origin + dir * h.t; let (i, _) = t.p.iter().enumerate().min_by(|a, b| a.1.dist(hit).total_cmp(&b.1.dist(hit)))?; Some((*t.nodes.get(i)?, hit)) } @@ -461,6 +474,33 @@ mod tests { assert_eq!(d.pick(Vec3::new(f64::NAN, 0.0, 0.0), Vec3::new(0.0, 0.0, -1.0)), None); } + /// A solid's drawn faces are numbered as the model numbers them, so a picked triangle names + /// the element face a set or load refers to. + #[test] + fn drawn_faces_carry_the_models_face_numbers() { + for kind in [ElementKind::Tet4, ElementKind::Tet10, ElementKind::Wedge6, ElementKind::Hex8, ElementKind::Hex20] { + let drawn = faces(kind); + assert_eq!(drawn.len(), kind.faces().len(), "{kind:?}"); + for (d, m) in drawn.iter().zip(kind.faces()) { + let corners = d.len(); + let (mut a, mut b) = (d.to_vec(), m.get(..corners).unwrap().to_vec()); + a.sort_unstable(); + b.sort_unstable(); + assert_eq!(a, b, "{kind:?}"); + } + } + let m = block(2, 1, 1); + let d = DisplayMesh::from_model(&m); + for t in &d.tris { + let face = m.face_nodes(fieldcraft_model::ElementFace { element: t.element, face: t.face }).unwrap(); + assert!(t.nodes.iter().all(|n| face.contains(n)), "{t:?} {face:?}"); + } + // The ray from +x hits element 2's face 4 (corners 1, 2, 6, 5: x = 2). + let h = d.hit_element_face(Vec3::new(10.0, 0.3, 0.6), Vec3::new(-1.0, 0.0, 0.0)).unwrap(); + let t = d.tris[h.index]; + assert_eq!((t.element, t.face, h.t), (2, 4, 8.0)); + } + #[test] fn empty_and_deterministic() { assert_eq!(DisplayMesh::from_model(&Model::default()), DisplayMesh::default()); diff --git a/crates/render/src/lib.rs b/crates/render/src/lib.rs index ed56dce..3d900ef 100644 --- a/crates/render/src/lib.rs +++ b/crates/render/src/lib.rs @@ -1,7 +1,8 @@ //! Display data for FieldCraft, with no GPU dependency: what a finite element model looks like //! on screen ([`DisplayMesh`]: the visible faces, element edges, line elements and nodes), an -//! orbit camera with perspective or orthographic projection, and a CPU rasterizer for headless -//! snapshots. The desktop viewport draws the same [`Scene`] on the GPU. +//! orbit camera with perspective or orthographic projection, picking ([`pick`]: what a ray through +//! a pixel hits), and a CPU rasterizer for headless snapshots. The desktop viewport draws the +//! same [`Scene`] on the GPU. #![deny(clippy::unwrap_used, clippy::expect_used, clippy::panic, clippy::unimplemented, clippy::todo, clippy::unreachable)] #![forbid(unsafe_code)] @@ -10,6 +11,7 @@ pub mod contour; mod display; mod font; pub mod glyphs; +pub mod pick; mod raster; mod vec; @@ -17,6 +19,6 @@ pub use camera::{Camera, Mat4, Orientation, StandardView}; pub use display::{DisplayEdge, DisplayLine, DisplayMesh, DisplayNode, DisplayTri, SurfaceTri}; pub use raster::{ Canvas, LABEL_OFFSET, MAX_DIM, MAX_PIXELS, MAX_PRIMITIVES, RenderError, Rgb, Scene, SceneLabel, SceneLine, ScenePoint, SceneTri, fit_size, png, - render, render_png, shade, + render, render_layers, render_layers_png, render_png, shade, }; pub use vec::{Aabb3, Vec3}; diff --git a/crates/render/src/pick.rs b/crates/render/src/pick.rs new file mode 100644 index 0000000..ff791ed --- /dev/null +++ b/crates/render/src/pick.rs @@ -0,0 +1,295 @@ +//! Picking: what a ray through a pixel hits first, and what lies near a pixel or inside a box on +//! screen. +//! +//! Rays are tested against triangles with the watertight ray–triangle intersection of Woop, +//! Benthin and Wald ("Watertight Ray/Triangle Intersection", Journal of Computer Graphics +//! Techniques 2(1), 2013): the ray is sheared onto the +z axis and the triangle's edges are tested +//! with 2D edge functions. An edge shared by two triangles is evaluated from the same operands in +//! both (once each way round, which IEEE arithmetic negates exactly), so a ray through a shared +//! edge or vertex hits at least one of them and never slips between. Hits are two-sided. The +//! nearest hit wins; equal distances go to the lower index, so the answer is the same on every +//! run and platform. +//! +//! Screen-space picks (nodes, edges, box selection) go through a [`Projector`] (the camera's +//! view-projection, computed once) and, for box selection, a [`DepthMap`] of the drawn faces that +//! tells which points are in front. + +use crate::raster::raster_tri; +use crate::{Camera, Mat4, Vec3}; + +/// The distance along a ray (`origin`, `dir`; `dir` need not be unit) at which it hits triangle +/// `p`, from either side; `None` when it misses, runs in the triangle's plane, or hits behind the +/// origin. Non-finite input misses. +pub fn ray_triangle(origin: Vec3, dir: Vec3, p: &[Vec3; 3]) -> Option { + let d = [dir.x, dir.y, dir.z]; + // The dominant axis of the direction becomes z; x and y keep the winding. + let kz = if d[0].abs() >= d[1].abs() && d[0].abs() >= d[2].abs() { + 0 + } else if d[1].abs() >= d[2].abs() { + 1 + } else { + 2 + }; + let (mut kx, mut ky) = ((kz + 1) % 3, (kz + 2) % 3); + let dz = *d.get(kz)?; + if !(dz.abs() > 0.0) || !dz.is_finite() { + return None; + } + if dz < 0.0 { + std::mem::swap(&mut kx, &mut ky); + } + let (sx, sy, sz) = (d.get(kx)? / dz, d.get(ky)? / dz, 1.0 / dz); + let rel = |q: Vec3| { + let r = q - origin; + [r.x, r.y, r.z] + }; + let (a, b, c) = (rel(p[0]), rel(p[1]), rel(p[2])); + let shear = |v: &[f64; 3]| -> Option<(f64, f64, f64)> { + let z = *v.get(kz)?; + Some((v.get(kx)? - sx * z, v.get(ky)? - sy * z, sz * z)) + }; + let ((ax, ay, az), (bx, by, bz), (cx, cy, cz)) = (shear(&a)?, shear(&b)?, shear(&c)?); + // Edge functions: e(P, Q) = Px·Qy − Py·Qx for the edges (C, B), (A, C), (B, A). + let u = cx * by - cy * bx; + let v = ax * cy - ay * cx; + let w = bx * ay - by * ax; + if !(u.is_finite() && v.is_finite() && w.is_finite()) { + return None; + } + if (u < 0.0 || v < 0.0 || w < 0.0) && (u > 0.0 || v > 0.0 || w > 0.0) { + return None; + } + let det = u + v + w; + if det == 0.0 { + return None; + } + let t = (u * az + v * bz + w * cz) / det; + (t >= 0.0 && t.is_finite()).then_some(t) +} + +/// A ray's first hit among triangles: the triangle's index and the distance along the ray. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct Hit { + pub index: usize, + pub t: f64, +} + +/// The first triangle `tris` (in order) that the ray hits: the nearest, the lowest index among +/// equally near ones. +pub fn first_hit(origin: Vec3, dir: Vec3, tris: impl IntoIterator) -> Option { + let mut best: Option = None; + for (index, p) in tris.into_iter().enumerate() { + if let Some(t) = ray_triangle(origin, dir, &p) + && best.is_none_or(|b| t < b.t) + { + best = Some(Hit { index, t }); + } + } + best +} + +/// The distance on screen from `p` to the segment `a`–`b` (pixels), and where along it (0 at +/// `a`, 1 at `b`) the nearest point lies. +pub fn segment_distance(p: [f64; 2], a: [f64; 2], b: [f64; 2]) -> (f64, f64) { + let (dx, dy) = (b[0] - a[0], b[1] - a[1]); + let len2 = dx * dx + dy * dy; + let k = if len2 > 0.0 { (((p[0] - a[0]) * dx + (p[1] - a[1]) * dy) / len2).clamp(0.0, 1.0) } else { 0.0 }; + let (qx, qy) = (a[0] + dx * k, a[1] + dy * k); + (((p[0] - qx).powi(2) + (p[1] - qy).powi(2)).sqrt(), k) +} + +/// A camera's projection for one viewport size, computed once: world points to pixels. +#[derive(Clone, Copy, Debug)] +pub struct Projector { + cam: Camera, + vp: Mat4, + w: f64, + h: f64, +} + +impl Projector { + /// The projection of `cam` onto a `w × h` viewport whose scene has radius `radius` (the + /// clipping planes', as drawn). + pub fn new(cam: &Camera, w: f64, h: f64, radius: f64) -> Projector { + let (w, h) = (if w.is_finite() { w.max(1.0) } else { 1.0 }, if h.is_finite() { h.max(1.0) } else { 1.0 }); + Projector { cam: *cam, vp: cam.view_proj(w / h, radius.max(1e-3)), w, h } + } + + /// The pixel (x right, y down) and depth (0 near … 1 far) of a world point in front of the + /// eye ([`Camera::to_screen`]). + pub fn screen(&self, p: Vec3) -> Option<[f64; 3]> { + let c = self.vp.apply(p); + if !(c[3] > 1e-12) { + return None; + } + let (x, y, z) = (c[0] / c[3], c[1] / c[3], c[2] / c[3]); + let s = [(x + 1.0) * 0.5 * self.w, (1.0 - y) * 0.5 * self.h, z]; + s.iter().all(|v| v.is_finite()).then_some(s) + } + + /// The picking ray through pixel (x, y): origin and unit direction ([`Camera::ray`]). + pub fn ray(&self, x: f64, y: f64) -> (Vec3, Vec3) { + self.cam.ray(x, y, self.w, self.h) + } + + pub fn size(&self) -> [f64; 2] { + [self.w, self.h] + } +} + +/// Most pixels a [`DepthMap`] holds; a larger viewport is sampled more coarsely. +pub const MAX_DEPTH_PIXELS: usize = 2048 * 2048; +/// Depth (of the 0…1 range) a point may lie behind the drawn surface and still count as on it. +const DEPTH_ALLOWANCE: f64 = 2e-3; + +/// The depth of the nearest drawn face at every pixel: which points are in front. +#[derive(Clone, Debug)] +pub struct DepthMap { + proj: Projector, + /// Map pixels per viewport pixel (1, or less for a large viewport). + scale: f64, + w: usize, + h: usize, + depth: Vec, +} + +impl DepthMap { + /// Rasterise `tris` with `proj`. `None` when the buffer cannot be allocated. + pub fn new(proj: &Projector, tris: impl IntoIterator) -> Option { + let [vw, vh] = proj.size(); + let scale = ((MAX_DEPTH_PIXELS as f64) / (vw * vh)).sqrt().min(1.0); + let (w, h) = (((vw * scale).ceil() as usize).max(1), ((vh * scale).ceil() as usize).max(1)); + let n = w.checked_mul(h)?; + let mut depth = Vec::new(); + depth.try_reserve_exact(n).ok()?; + depth.resize(n, f32::INFINITY); + let mut map = DepthMap { proj: *proj, scale, w, h, depth }; + for p in tris { + map.add(&p); + } + Some(map) + } + + fn add(&mut self, p: &[Vec3; 3]) { + let s = |q: Vec3| self.proj.screen(q).map(|v| [v[0] * self.scale, v[1] * self.scale, v[2]]); + let (Some(a), Some(b), Some(c)) = (s(p[0]), s(p[1]), s(p[2])) else { return }; + let (w, depth) = (self.w, &mut self.depth); + raster_tri(a, b, c, w, self.h, |x, y, z| { + if let Some(d) = depth.get_mut(y * w + x) { + *d = d.min(z); + } + }); + } + + /// Whether a world point is on screen and not behind a drawn face (a point on a face is in + /// front of it). A point whose pixel no face covers, on a face's outline, is compared with + /// the nearest face around it, so what lies behind an outline stays hidden. + pub fn visible(&self, p: Vec3) -> bool { + let Some([x, y, z]) = self.proj.screen(p) else { return false }; + let [vw, vh] = self.proj.size(); + if !(0.0..=vw).contains(&x) || !(0.0..=vh).contains(&y) || !(0.0..=1.0).contains(&z) { + return false; + } + let (mx, my) = (((x * self.scale) as usize).min(self.w - 1), ((y * self.scale) as usize).min(self.h - 1)); + let at = |x: usize, y: usize| self.depth.get(y * self.w + x).copied().unwrap_or(f32::INFINITY); + let mut d = at(mx, my); + if d == f32::INFINITY { + for y in my.saturating_sub(1)..=(my + 1).min(self.h - 1) { + for x in mx.saturating_sub(1)..=(mx + 1).min(self.w - 1) { + d = d.min(at(x, y)); + } + } + } + z <= f64::from(d) + DEPTH_ALLOWANCE + } +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::StandardView; + + const SQUARE: [[Vec3; 3]; 2] = [ + [Vec3::new(0.0, 0.0, 0.0), Vec3::new(1.0, 0.0, 0.0), Vec3::new(1.0, 1.0, 0.0)], + [Vec3::new(0.0, 0.0, 0.0), Vec3::new(1.0, 1.0, 0.0), Vec3::new(0.0, 1.0, 0.0)], + ]; + + #[test] + fn rays_hit_from_either_side_and_report_the_distance() { + let down = Vec3::new(0.0, 0.0, -1.0); + assert_eq!(ray_triangle(Vec3::new(0.7, 0.2, 5.0), down, &SQUARE[0]), Some(5.0)); + assert_eq!(ray_triangle(Vec3::new(0.7, 0.2, -2.0), Vec3::Z, &SQUARE[0]), Some(2.0)); + // Behind the origin, beside the triangle, in its plane. + assert_eq!(ray_triangle(Vec3::new(0.7, 0.2, -2.0), down, &SQUARE[0]), None); + assert_eq!(ray_triangle(Vec3::new(0.2, 0.7, 5.0), down, &SQUARE[0]), None); + assert_eq!(ray_triangle(Vec3::new(-1.0, 0.5, 0.0), Vec3::X, &SQUARE[0]), None); + // A slanted ray, a zero or non-finite one. + let t = ray_triangle(Vec3::new(0.0, 0.0, 1.0), Vec3::new(0.5, 0.25, -1.0), &SQUARE[0]).unwrap(); + assert!((t - 1.0).abs() < 1e-15); + assert_eq!(ray_triangle(Vec3::new(0.5, 0.2, 1.0), Vec3::ZERO, &SQUARE[0]), None); + assert_eq!(ray_triangle(Vec3::new(f64::NAN, 0.2, 1.0), down, &SQUARE[0]), None); + assert_eq!(ray_triangle(Vec3::new(0.5, 0.2, 1.0), Vec3::new(0.0, 0.0, f64::INFINITY), &SQUARE[0]), None); + let huge = [Vec3::new(1e300, 0.0, 0.0), Vec3::new(-1e300, 1e300, 0.0), Vec3::new(-1e300, -1e300, 0.0)]; + assert_eq!(ray_triangle(Vec3::new(0.0, 0.0, 1.0), down, &huge), None, "overflowing edge functions miss"); + } + + /// Rays through the shared diagonal (inside the square) and corners of two triangles, at + /// every angle, always hit one of them. + #[test] + fn no_ray_slips_through_a_shared_edge_or_vertex() { + for k in 1..400 { + let s = k as f64 / 400.0; + for dir in [Vec3::new(0.0, 0.0, -1.0), Vec3::new(0.3, -0.1, -1.0), Vec3::new(-1e-9, 2e-9, -1.0)] { + let at = Vec3::new(s, s, 0.0) - dir * 3.0; + assert!(first_hit(at, dir, SQUARE).is_some(), "{s} {dir:?}"); + } + } + for corner in [Vec3::ZERO, Vec3::new(1.0, 1.0, 0.0)] { + assert!(first_hit(corner + Vec3::Z, Vec3::new(0.0, 0.0, -1.0), SQUARE).is_some()); + } + } + + #[test] + fn the_nearest_hit_wins_and_ties_go_to_the_lower_index() { + let lift = |t: [Vec3; 3], z: f64| t.map(|p| p + Vec3::new(0.0, 0.0, z)); + let tris = [lift(SQUARE[0], 0.0), lift(SQUARE[0], 2.0), lift(SQUARE[0], 1.0)]; + let hit = first_hit(Vec3::new(0.7, 0.2, 10.0), Vec3::new(0.0, 0.0, -1.0), tris).unwrap(); + assert_eq!(hit, Hit { index: 1, t: 8.0 }); + let same = [SQUARE[0], SQUARE[0]]; + assert_eq!(first_hit(Vec3::new(0.7, 0.2, 1.0), Vec3::new(0.0, 0.0, -1.0), same).map(|h| h.index), Some(0)); + assert_eq!(first_hit(Vec3::new(0.7, 0.2, 1.0), Vec3::new(0.0, 0.0, -1.0), []), None); + } + + #[test] + fn screen_distances_to_segments() { + assert_eq!(segment_distance([1.0, 1.0], [0.0, 0.0], [2.0, 0.0]), (1.0, 0.5)); + assert_eq!(segment_distance([3.0, 0.0], [0.0, 0.0], [2.0, 0.0]), (1.0, 1.0)); + assert_eq!(segment_distance([0.0, 2.0], [0.0, 0.0], [0.0, 0.0]), (2.0, 0.0)); + } + + #[test] + fn the_projector_agrees_with_the_camera_and_the_depth_map_hides_what_is_behind() { + let mut cam = Camera { distance: 10.0, ..Camera::default() }; + cam.set_view(StandardView::Top); + cam.target = Vec3::new(0.5, 0.5, 0.0); + let p = Projector::new(&cam, 400.0, 300.0, 2.0); + let q = Vec3::new(0.25, 0.75, 0.0); + let (x, y, z) = cam.to_screen(q, 400.0, 300.0, 2.0).unwrap(); + assert_eq!(p.screen(q), Some([x, y, z])); + // The ray through a point's pixel passes through the point. + let (o, d) = p.ray(x, y); + assert!((q - o).cross(d).len() < 1e-9); + let map = DepthMap::new(&p, SQUARE).unwrap(); + assert!(map.visible(q), "a point on the face"); + assert!(map.visible(q + Vec3::Z * 0.5), "a point above it"); + assert!(!map.visible(q - Vec3::Z * 0.5), "a point under it"); + assert!(map.visible(Vec3::new(1.5, 0.5, -0.5)), "beside it, nothing covers the point"); + assert!(!map.visible(Vec3::new(500.0, 0.5, 0.0)), "off screen"); + assert!(!map.visible(Vec3::new(f64::NAN, 0.5, 0.0))); + // A huge viewport is sampled coarsely. + let big = DepthMap::new(&Projector::new(&cam, 8192.0, 8192.0, 2.0), SQUARE).unwrap(); + assert!(big.depth.len() <= MAX_DEPTH_PIXELS + 2 * 2048 + 1); + let hostile = Projector::new(&cam, f64::NAN, -1.0, f64::INFINITY); + assert!(DepthMap::new(&hostile, SQUARE).is_some()); + } +} diff --git a/crates/render/src/raster.rs b/crates/render/src/raster.rs index 8dc6125..b3d50d1 100644 --- a/crates/render/src/raster.rs +++ b/crates/render/src/raster.rs @@ -5,6 +5,7 @@ //! Adapted from SolidCraft's `crates/render/src/raster.rs` (MIT OR Apache-2.0), see //! ATTRIBUTION.md. +use crate::pick::Projector; use crate::{Camera, Vec3, font}; #[derive(Clone, Copy, Debug, PartialEq, Eq)] @@ -179,10 +180,22 @@ impl Canvas { /// 1..=[`MAX_DIM`] per side or over [`MAX_PIXELS`], and scenes over [`MAX_PRIMITIVES`], are /// refused ([`fit_size`] gives an accepted size). pub fn render(scene: &Scene, cam: &Camera, w: usize, h: usize) -> Result { + render_layers(scene, None, cam, w, h) +} + +/// Depth (of the 0…1 range) an overlay face may lie behind the scene's and still be drawn: it +/// covers the faces it repeats. +const OVERLAY_DEPTH: f32 = 2e-4; + +/// [`render`] with an overlay drawn over the scene: its faces over the scene's faces they +/// coincide with (a highlight repeats the faces it marks), then its lines and points. The +/// overlay's background and radius are not used. +pub fn render_layers(scene: &Scene, overlay: Option<&Scene>, cam: &Camera, w: usize, h: usize) -> Result { if !size_ok(w, h) { return Err(RenderError::Size { w, h }); } - let primitives = scene.tris.len().saturating_add(scene.lines.len()).saturating_add(scene.points.len()).saturating_add(scene.labels.len()); + let count = |s: &Scene| s.tris.len().saturating_add(s.lines.len()).saturating_add(s.points.len()).saturating_add(s.labels.len()); + let primitives = count(scene).saturating_add(overlay.map_or(0, count)); if primitives > MAX_PRIMITIVES { return Err(RenderError::SceneTooLarge { primitives }); } @@ -197,48 +210,23 @@ pub fn render(scene: &Scene, cam: &Camera, w: usize, h: usize) -> Result Option<[f64; 3]> { - let q = vp.apply(p); - if q[3] <= 1e-9 { - return None; - } - Some([(q[0] / q[3] + 1.0) * 0.5 * sw as f64, (1.0 - q[1] / q[3]) * 0.5 * sh as f64, q[2] / q[3]]) - }; + let proj = Projector::new(cam, sw as f64, sh as f64, scene.radius); + let project = |p: Vec3| proj.screen(p); let view_dir = cam.back(); - for t in &scene.tris { + for (t, allowance) in scene.tris.iter().map(|t| (t, 0.0)).chain(overlay.into_iter().flat_map(|o| &o.tris).map(|t| (t, OVERLAY_DEPTH))) { let (Some(a), Some(b), Some(cc)) = (project(t.p[0]), project(t.p[1]), project(t.p[2])) else { continue }; let col = shade(t.normal, view_dir, t.color.f()); - let (x0, x1) = (a[0].min(b[0]).min(cc[0]).floor().max(0.0), a[0].max(b[0]).max(cc[0]).ceil().min(sw as f64 - 1.0)); - let (y0, y1) = (a[1].min(b[1]).min(cc[1]).floor().max(0.0), a[1].max(b[1]).max(cc[1]).ceil().min(sh as f64 - 1.0)); - if !(x0 <= x1 && y0 <= y1) { - continue; - } - let area = (b[0] - a[0]) * (cc[1] - a[1]) - (b[1] - a[1]) * (cc[0] - a[0]); - if !(area.abs() >= 1e-12) { - continue; - } - for y in y0 as usize..=y1 as usize { - for x in x0 as usize..=x1 as usize { - let (px, py) = (x as f64 + 0.5, y as f64 + 0.5); - let w0 = ((b[0] - px) * (cc[1] - py) - (b[1] - py) * (cc[0] - px)) / area; - let w1 = ((cc[0] - px) * (a[1] - py) - (cc[1] - py) * (a[0] - px)) / area; - let w2 = 1.0 - w0 - w1; - if w0 < -1e-9 || w1 < -1e-9 || w2 < -1e-9 { - continue; - } - let z = (w0 * a[2] + w1 * b[2] + w2 * cc[2]) as f32; - let i = y * sw + x; - if let Some(d) = c.depth.get_mut(i) - && z < *d - { - *d = z; - c.put(x, y, col, 1.0); - } + raster_tri(a, b, cc, sw, sh, |x, y, z| { + if let Some(d) = c.depth.get_mut(y * sw + x) + && z < *d + allowance + { + *d = d.min(z); + c.put(x, y, col, 1.0); } - } + }); } - for l in &scene.lines { + let layers = || std::iter::once(scene).chain(overlay); + for l in layers().flat_map(|s| &s.lines) { let (Some(a), Some(b)) = (project(l.a), project(l.b)) else { continue }; let len = ((b[0] - a[0]).powi(2) + (b[1] - a[1]).powi(2)).sqrt(); let steps = if len.is_finite() { (len.ceil() as usize).clamp(1, MAX_LINE_STEPS) } else { continue }; @@ -249,7 +237,7 @@ pub fn render(scene: &Scene, cam: &Camera, w: usize, h: usize) -> Result Result= 1e-12) { + return; + } + for y in y0 as usize..=y1 as usize { + for x in x0 as usize..=x1 as usize { + let (px, py) = (x as f64 + 0.5, y as f64 + 0.5); + let w0 = ((b[0] - px) * (c[1] - py) - (b[1] - py) * (c[0] - px)) / area; + let w1 = ((c[0] - px) * (a[1] - py) - (c[1] - py) * (a[0] - px)) / area; + let w2 = 1.0 - w0 - w1; + if w0 < -1e-9 || w1 < -1e-9 || w2 < -1e-9 { + continue; + } + pixel(x, y, (w0 * a[2] + w1 * b[2] + w2 * c[2]) as f32); + } + } +} + /// Text above and right of screen point `p`, over everything, with a halo of the opposite /// lightness so it reads on any background. fn label(c: &mut Canvas, p: [f64; 3], text: &str, col: [f64; 3]) { @@ -350,7 +366,12 @@ fn downsample(c: &Canvas, w: usize, h: usize) -> Result { /// Render to PNG bytes. pub fn render_png(scene: &Scene, cam: &Camera, w: usize, h: usize) -> Result, RenderError> { - let c = render(scene, cam, w, h)?; + render_layers_png(scene, None, cam, w, h) +} + +/// [`render_layers`] as PNG bytes. +pub fn render_layers_png(scene: &Scene, overlay: Option<&Scene>, cam: &Camera, w: usize, h: usize) -> Result, RenderError> { + let c = render_layers(scene, overlay, cam, w, h)?; png(c.w, c.h, &c.rgba).ok_or(RenderError::Encode) } @@ -407,6 +428,25 @@ mod tests { assert!(p[2] > 200 && p[0] < 60, "{p:?}"); } + /// An overlay face repeating a scene face covers it; one behind another face does not show. + #[test] + fn an_overlay_covers_the_faces_it_repeats() { + let scene = red_triangle(); + let mut cam = Camera { distance: 50.0, ..Default::default() }; + cam.set_view(StandardView::Front); + let px = |c: &Canvas| c.rgba[(32 * c.w + 32) * 4..(32 * c.w + 32) * 4 + 3].to_vec(); + let green = |s: &Scene| Scene { tris: s.tris.iter().map(|t| SceneTri { color: Rgb(0, 255, 0), ..*t }).collect(), ..Scene::default() }; + let over = render_layers(&scene, Some(&green(&scene)), &cam, 64, 64).unwrap(); + assert!(px(&over)[1] > px(&over)[0] + 50, "{:?}", px(&over)); + // The same face moved behind the red one (front view looks along +y) stays hidden. + let mut behind = green(&scene); + for t in &mut behind.tris { + t.p = t.p.map(|p| p + Vec3::Y * 5.0); + } + let hidden = render_layers(&scene, Some(&behind), &cam, 64, 64).unwrap(); + assert_eq!(px(&hidden), px(&render(&scene, &cam, 64, 64).unwrap())); + } + #[test] fn png_bytes_and_hostile_sizes() { let scene = red_triangle();