geometry: nearest ray hits

crates/nanobrush-geometry/src/lib.rs

@@ -103,4 +103,5 @@ impl Budget {
     }
 }
 pub mod patch;
+pub mod ray;
 pub mod world;

crates/nanobrush-geometry/src/ray.rs

@@ -0,0 +1,318 @@
+use crate::{math::*, *};
+#[derive(Clone, Copy, Debug)]
+pub struct Ray {
+    origin: DVec3,
+    direction: DVec3,
+}
+impl Ray {
+    pub fn new(origin: DVec3, direction: DVec3) -> Result<Self> {
+        Ok(Self {
+            origin: finite(origin)?,
+            direction: unit(direction)?,
+        })
+    }
+    pub fn at(self, t: f64) -> Result<DVec3> {
+        if !t.is_finite() {
+            return Err(Error::NonFinite);
+        }
+        finite(add(self.origin, scale(self.direction, t)))
+    }
+}
+#[derive(Clone, Copy, Debug, PartialEq, Eq, PartialOrd, Ord)]
+pub enum Target {
+    Brush {
+        brush: usize,
+        plane: usize,
+    },
+    Patch {
+        patch: usize,
+        quad: usize,
+        triangle: usize,
+    },
+}
+#[derive(Clone, Copy, Debug)]
+pub struct Hit {
+    pub distance: f64,
+    pub point: DVec3,
+    pub normal: DVec3,
+    pub target: Target,
+}
+fn triangle(
+    ray: Ray,
+    a: DVec3,
+    b: DVec3,
+    c: DVec3,
+    max: f64,
+    p: Policy,
+) -> Result<Option<(f64, DVec3)>> {
+    finite(a)?;
+    finite(b)?;
+    finite(c)?;
+    let e1 = finite(sub(b, a))?;
+    let e2 = finite(sub(c, a))?;
+    let h = finite(cross(ray.direction, e2))?;
+    let det = dot(e1, h);
+    let area_scale = length(e1) * length(e2);
+    if !det.is_finite() || !area_scale.is_finite() {
+        return Err(Error::NonFinite);
+    }
+    if det.abs() <= p.angular_epsilon * area_scale {
+        return Ok(None);
+    }
+    let s = sub(ray.origin, a);
+    let u = dot(s, h) / det;
+    let q = finite(cross(s, e1))?;
+    let v = dot(ray.direction, q) / det;
+    let t = dot(e2, q) / det;
+    if !u.is_finite() || !v.is_finite() || !t.is_finite() {
+        return Err(Error::NonFinite);
+    }
+    if u < -p.angular_epsilon
+        || v < -p.angular_epsilon
+        || u + v > 1. + p.angular_epsilon
+        || t < -p.distance_epsilon
+        || t > max
+    {
+        return Ok(None);
+    }
+    Ok(Some((t.max(0.), unit(cross(e1, e2))?)))
+}
+struct Query {
+    ray: Ray,
+    max: f64,
+    policy: Policy,
+    budget: Budget,
+    hit: Option<Hit>,
+    vertices: usize,
+}
+impl Query {
+    fn vertices(&mut self, n: usize) -> Result<()> {
+        self.vertices = self.vertices.checked_add(n).ok_or(Error::Overflow)?;
+        if self.vertices > self.budget.limits.max_vertices {
+            Err(Error::Limit)
+        } else {
+            Ok(())
+        }
+    }
+    fn test(&mut self, v: [DVec3; 3], target: Target) -> Result<()> {
+        self.budget.work(1)?;
+        if let Some((distance, normal)) =
+            triangle(self.ray, v[0], v[1], v[2], self.max, self.policy)?
+        {
+            let replace = self.hit.is_none_or(|old| {
+                distance < old.distance - self.policy.distance_epsilon
+                    || ((distance - old.distance).abs() <= self.policy.distance_epsilon
+                        && target < old.target)
+            });
+            if replace {
+                self.hit = Some(Hit {
+                    distance,
+                    point: self.ray.at(distance)?,
+                    normal,
+                    target,
+                });
+            }
+        }
+        Ok(())
+    }
+}
+/// Linear double-sided nearest surface query over explicit computed geometry.
+/// No selection/visibility state exists in this model. Ray directions are unit,
+/// so returned parameters are world distances rather than caller vector scales.
+pub fn nearest(
+    world: &mb3_dat::World,
+    ray: Ray,
+    max_distance: f64,
+    policy: Policy,
+    limits: ComputeLimits,
+) -> Result<Option<Hit>> {
+    let policy = policy.validate()?;
+    if !max_distance.is_finite() || max_distance < 0. {
+        return Err(Error::InvalidDimensions);
+    }
+    if !world.computed {
+        return Err(Error::InconsistentData);
+    }
+    let mut q = Query {
+        ray,
+        max: max_distance,
+        policy,
+        budget: Budget {
+            limits,
+            stats: Stats::default(),
+        },
+        hit: None,
+        vertices: 0,
+    };
+    for (bi, b) in world.brushes.iter().enumerate() {
+        q.budget.work(1)?;
+        for (pi, p) in b.planes.iter().enumerate() {
+            q.budget.work(1)?;
+            let vertices = p.polygon.as_ref().ok_or(Error::InconsistentData)?;
+            if vertices.len() < 3 {
+                return Err(Error::InvalidBrush);
+            }
+            q.vertices(vertices.len())?;
+            for i in 1..vertices.len() - 1 {
+                q.test(
+                    [
+                        from_f32(vertices[0])?,
+                        from_f32(vertices[i])?,
+                        from_f32(vertices[i + 1])?,
+                    ],
+                    Target::Brush {
+                        brush: bi,
+                        plane: pi,
+                    },
+                )?;
+            }
+        }
+    }
+    for (pi, p) in world.patches.iter().enumerate() {
+        q.budget.work(1)?;
+        let data = p.computed.as_ref().ok_or(Error::InconsistentData)?;
+        let u = (p.anchors_u as usize / 2)
+            .checked_mul(p.quads_u as usize)
+            .and_then(|n| n.checked_add(1))
+            .ok_or(Error::Overflow)?;
+        let v = (p.anchors_v as usize / 2)
+            .checked_mul(p.quads_v as usize)
+            .and_then(|n| n.checked_add(1))
+            .ok_or(Error::Overflow)?;
+        let n = u.checked_mul(v).ok_or(Error::Overflow)?;
+        if data.points.len() != n || data.normals.len() != n {
+            return Err(Error::InconsistentData);
+        }
+        q.vertices(n)?;
+        for point in &data.points {
+            q.budget.work(1)?;
+            from_f32(*point)?;
+        }
+        for y in 0..v - 1 {
+            for x in 0..u - 1 {
+                let i = x + y * u;
+                let quad = x + y * (u - 1);
+                for (t, ids) in [[i, i + u, i + u + 1], [i, i + u + 1, i + 1]]
+                    .into_iter()
+                    .enumerate()
+                {
+                    q.test(
+                        [
+                            from_f32(data.points[ids[0]])?,
+                            from_f32(data.points[ids[1]])?,
+                            from_f32(data.points[ids[2]])?,
+                        ],
+                        Target::Patch {
+                            patch: pi,
+                            quad,
+                            triangle: t,
+                        },
+                    )?;
+                }
+            }
+        }
+    }
+    Ok(q.hit)
+}
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use mb3_dat::{Limits, World};
+    #[test]
+    fn cube_nearest_inside_miss_and_stable_ties() {
+        let source = World {
+            brushes: vec![crate::brush::cube(), crate::brush::cube()],
+            ..Default::default()
+        };
+        let (w, _) = crate::world::compute(
+            &source,
+            Policy::default(),
+            ComputeLimits::default(),
+            Limits::default(),
+        )
+        .unwrap();
+        let ray = Ray::new([0., 0., 3.], [0., 0., -2.]).unwrap();
+        let hit = nearest(&w, ray, 10., Policy::default(), ComputeLimits::default())
+            .unwrap()
+            .unwrap();
+        assert_eq!(hit.distance, 2.);
+        assert_eq!(hit.target, Target::Brush { brush: 0, plane: 5 });
+        assert_eq!(hit.normal, [0., 0., 1.]);
+        assert!(
+            nearest(&w, ray, 1., Policy::default(), ComputeLimits::default())
+                .unwrap()
+                .is_none()
+        );
+        let inside = nearest(
+            &w,
+            Ray::new([0.; 3], [1., 0., 0.]).unwrap(),
+            10.,
+            Policy::default(),
+            ComputeLimits::default(),
+        )
+        .unwrap()
+        .unwrap();
+        assert_eq!(inside.distance, 1.);
+        assert!(
+            nearest(
+                &w,
+                Ray::new([0., 0., 3.], [0., 0., 1.]).unwrap(),
+                10.,
+                Policy::default(),
+                ComputeLimits::default()
+            )
+            .unwrap()
+            .is_none()
+        );
+    }
+    #[test]
+    fn patch_hit_backfaces_and_budget() {
+        let source = World {
+            patches: vec![crate::patch::flat(3, 3)],
+            ..Default::default()
+        };
+        let (w, _) = crate::world::compute(
+            &source,
+            Policy::default(),
+            ComputeLimits::default(),
+            Limits::default(),
+        )
+        .unwrap();
+        for (origin, dir) in [
+            ([0.5, 3., 0.5], [0., -1., 0.]),
+            ([0.5, -3., 0.5], [0., 1., 0.]),
+        ] {
+            let ray = Ray::new(origin, dir).unwrap();
+            let hit = nearest(&w, ray, 10., Policy::default(), ComputeLimits::default())
+                .unwrap()
+                .unwrap();
+            assert_eq!(hit.distance, 3.);
+            assert_eq!(hit.normal, [0., 1., 0.]);
+            assert!(matches!(hit.target, Target::Patch { patch: 0, .. }));
+            assert!(
+                nearest(
+                    &w,
+                    ray,
+                    10.,
+                    Policy::default(),
+                    ComputeLimits {
+                        max_work: 1,
+                        ..Default::default()
+                    }
+                )
+                .is_err()
+            );
+        }
+        assert!(Ray::new([0.; 3], [0.; 3]).is_err());
+        assert!(
+            nearest(
+                &w,
+                Ray::new([0.; 3], [1., 0., 0.]).unwrap(),
+                f64::INFINITY,
+                Policy::default(),
+                ComputeLimits::default()
+            )
+            .is_err()
+        );
+    }
+}

docs/GEOMETRY.md

@@ -23,3 +23,9 @@ Anchor dimensions must be odd and at least three, and quads per tile at least tw
 `world::compute` returns a new world, preserving source tables, slots, bases and anchors. It fallibly copies editable data, rebuilds every brush/patch, checks f64→f32 range, and sets computed=true only after all objects succeed. Prior derived arrays are intentionally replaced. Properties are copied unchanged, not normalized; codec-specific validation still happens at save. No source mutation or partial world escapes on a late error.
 
 Work and generated-position budgets apply across the whole world. Requested copy/geometry/output allocations are charged cumulatively (including scratch already freed, conservatively); model element/string limits also apply. CLI `convert/preflight --compute` opts into this native policy before target analysis/encoding. Defaults cap 64 planes per brush, 100,000 generated positions, 10 million work units, and the lower of codec/geometry storage limits. Empty worlds compute successfully; malformed/redundant/open brushes and degenerate patches fail instead of being dropped. Mixed cube/patch, property preservation, DAT round trip, late failure, global-budget and CLI cube-to-OBJ tests pass.
+
+## Nearest-hit queries
+
+`ray::nearest` scans explicit computed brush polygons (convex fans) and patch quads split along the first-to-third diagonal. Queries are double-sided. Ray directions are normalized, so hit parameters are world distances. Negative distances within distance epsilon clamp to zero; larger negative hits and hits beyond finite max distance are excluded. Parallel/degenerate triangles produce no hit; non-finite geometry and malformed arrays fail.
+
+Ties within distance epsilon use stable target order: brush before patch, then ascending object/plane or quad/triangle indices. Returned normals follow polygon/quad winding and do not flip toward the ray. Work/vertex caps apply to the whole linear query, including parents and patch point scans. No acceleration, visibility, selection, near-clip or callback semantics are implied; unused properties/stored normals are not revalidated. Tests cover nearest/inside/miss/max-distance, duplicate-object ties, front/back patch hits, direction normalization and invalid/budget cases.

docs/PROGRESS_2026-10-26.md

@@ -12,5 +12,5 @@ Scope: extend native headless formats and geometry. Each stage is independently
 | 6 | P5.16 | Convex brush contours | Complete: bounded half-space intersections, recession rejection, deterministic outward contours, closed-edge and volume checks. Cube/translation/winding/invalid/redundant/budget tests pass. |
 | 7 | P5.17 | Patch9 tessellation | Complete: tiled tensor-quadratic positions, analytic derivatives, OBJ-consistent normals, shared-position seam policy and bounded allocations/work. Planar/curved/multi-tile/degenerate/limit tests pass. |
 | 8 | P5.18 | Computed world pipeline/CLI | Complete: fallible editable-data copying, aggregate work/storage/position budgets, complete-or-error derived data and explicit `--compute`. Mixed-world/property/DAT/late-failure/global-limit and CLI cube→OBJ tests pass. |
-| 9 | P5.19 | Nearest-hit ray queries | Pending |
+| 9 | P5.19 | Nearest-hit ray queries | Complete: bounded double-sided computed-surface queries, normalized rays, stable object/face/quad ties and winding normals. Cube/inside/miss/range/tie/patch/backface/budget tests pass. |
 | 10 | P5.20 | Property table edits/reindexing | Pending |

docs/TODO.md

@@ -134,6 +134,8 @@ Each task ID is stable. Split tasks into child IDs as evidence grows. Record dep
 
 - [x] P5.18 Integrate bounded complete-or-error computed worlds and explicit CLI geometry generation ([GEOMETRY.md](GEOMETRY.md)).
 
+- [x] P5.19 Implement bounded deterministic native nearest-hit ray queries over computed brush/patch surfaces ([GEOMETRY.md](GEOMETRY.md)).
+
 ## P6 — File formats and converters
 
 **Depends on:** P3 oracle where available, P5 document model. **Exit:** standalone format specifications, declared compatibility subsets, fixture corpus, and native codecs.