From 1c660dda13ab9fc5a05e3a6d82ac43c9e62935c5 Mon Sep 17 00:00:00 2001 From: Nathan Vegdahl Date: Thu, 31 Dec 2015 00:32:38 -0800 Subject: [PATCH] Generalized the BVH to... pretty much anything. The BVH is now generic over any kind of data. The building function takes in a closure that can bound the given data type in 3d space, and the rest just works. --- src/bvh.rs | 119 +++++++++++++++++++++++++++++------------------- src/main.rs | 67 +++++++++++++++++---------- src/math/mod.rs | 14 ++++-- 3 files changed, 125 insertions(+), 75 deletions(-) diff --git a/src/bvh.rs b/src/bvh.rs index 7943932..512372b 100644 --- a/src/bvh.rs +++ b/src/bvh.rs @@ -7,8 +7,10 @@ use triangle; use algorithm::partition; #[derive(Debug)] -pub struct BVH { +pub struct BVH<'a, T: 'a> { nodes: Vec, + objects: &'a [T], + depth: usize, } #[derive(Debug)] @@ -20,36 +22,51 @@ enum BVHNode { Leaf { bounds: BBox, - triangle: (Point, Point, Point), + object_index: usize, }, } -impl BVH { - pub fn from_triangles(triangles: &mut [(Point, Point, Point)]) -> BVH { - let mut bvh = BVH { nodes: Vec::new() }; +impl<'a, T> BVH<'a, T> { + pub fn from_objects(objects: &'a mut [T], bounder: F) -> BVH<'a, T> + where F: Fn(&T) -> BBox + { + let mut bvh = BVH { + nodes: Vec::new(), + objects: &[], + depth: 0, + }; - bvh.recursive_build(triangles); + bvh.recursive_build(0, 0, objects, &bounder); + bvh.objects = objects; + + println!("BVH Depth: {}", bvh.depth); bvh } - // Recursively builds the BVH starting at the given node with the given - // first and last primitive indices (in bag). - fn recursive_build(&mut self, triangles: &mut [(Point, Point, Point)]) -> usize { + + fn recursive_build(&mut self, + offset: usize, + depth: usize, + objects: &mut [T], + bounder: &F) + -> usize + where F: Fn(&T) -> BBox + { let me = self.nodes.len(); - if triangles.len() == 1 { + if objects.len() == 0 { + return 0; + } else if objects.len() == 1 { // Leaf node - let tri = triangles[0]; - self.nodes.push(BVHNode::Leaf { - bounds: { - let minimum = tri.0.min(tri.1.min(tri.2)); - let maximum = tri.0.max(tri.1.max(tri.2)); - BBox::from_points(minimum, maximum) - }, - triangle: tri, + bounds: bounder(&objects[0]), + object_index: offset, }); + + if self.depth < depth { + self.depth = depth; + } } else { // Not a leaf node self.nodes.push(BVHNode::Internal { @@ -58,18 +75,10 @@ impl BVH { }); // Determine which axis to split on - fn tri_bounds(tri: (Point, Point, Point)) -> BBox { - let minimum = tri.0.min(tri.1.min(tri.2)); - let maximum = tri.0.max(tri.1.max(tri.2)); - BBox { - min: minimum, - max: maximum, - } - } let bounds = { let mut bounds = BBox::new(); - for tri in &triangles[..] { - bounds = bounds | tri_bounds(*tri); + for obj in objects.iter() { + bounds = bounds | bounder(obj); } bounds }; @@ -87,10 +96,10 @@ impl BVH { }; let split_pos = (bounds.min[split_axis] + bounds.max[split_axis]) * 0.5; - // Partition triangles based on split + // Partition objects based on split let split_index = { - let mut split_i = partition(triangles, |tri| { - let tb = tri_bounds(*tri); + let mut split_i = partition(objects, |obj| { + let tb = bounder(obj); let centroid = (tb.min[split_axis] + tb.max[split_axis]) * 0.5; centroid < split_pos }); @@ -102,8 +111,11 @@ impl BVH { }; // Create child nodes - self.recursive_build(&mut triangles[..split_index]); - let child2_index = self.recursive_build(&mut triangles[split_index..]); + self.recursive_build(offset, depth + 1, &mut objects[..split_index], bounder); + let child2_index = self.recursive_build(offset + split_index, + depth + 1, + &mut objects[split_index..], + bounder); // Set node self.nodes[me] = BVHNode::Internal { @@ -117,11 +129,18 @@ impl BVH { } -pub fn intersect_bvh(bvh: &BVH, ray: &Ray) -> bool { - let mut i_stack = [0; 64]; - let mut stack_ptr = 0; +pub fn intersect_bvh(bvh: &BVH<(Point, Point, Point)>, ray: &mut Ray) -> Option<(f32, f32, f32)> { + if bvh.nodes.len() == 0 { + return None; + } - loop { + let mut i_stack = [0; 65]; + let mut stack_ptr: usize = 1; + let mut hit = false; + let mut u = 0.0; + let mut v = 0.0; + + while stack_ptr > 0 { match bvh.nodes[i_stack[stack_ptr]] { BVHNode::Internal { bounds, second_child_index } => { if bounds.intersect_ray(ray) { @@ -129,25 +148,29 @@ pub fn intersect_bvh(bvh: &BVH, ray: &Ray) -> bool { i_stack[stack_ptr + 1] = second_child_index; stack_ptr += 1; } else { - if stack_ptr == 0 { - break; - } stack_ptr -= 1; } } - BVHNode::Leaf{bounds: _, triangle: tri} => { - if let Some(_) = triangle::intersect_ray(ray, tri) { - return true; - } else { - if stack_ptr == 0 { - break; + BVHNode::Leaf { bounds: _, object_index } => { + if let Some((t, tri_u, tri_v)) = + triangle::intersect_ray(ray, bvh.objects[object_index]) { + if t < ray.max_t { + hit = true; + ray.max_t = t; + u = tri_u; + v = tri_v; } - stack_ptr -= 1; } + + stack_ptr -= 1; } } } - return false; + if hit { + return Some((ray.max_t, u, v)); + } else { + return None; + } } diff --git a/src/main.rs b/src/main.rs index e52d20e..fc4ee07 100644 --- a/src/main.rs +++ b/src/main.rs @@ -20,6 +20,7 @@ use docopt::Docopt; use image::Image; use math::{Point, Vector, fast_logit}; use ray::Ray; +use bbox::BBox; // ---------------------------------------------------------------- @@ -74,18 +75,32 @@ fn main() { // Generate a scene of triangles let mut triangles = { let mut triangles = Vec::new(); - for x in 0..10 { - for y in 0..10 { - let cx = x as f32 * 32.0; - let cy = y as f32 * 32.0; + let xres = 512; + let yres = 512; + let xinc = 512.0 / (xres as f32); + let yinc = 512.0 / (yres as f32); + for x in 0..xres { + for y in 0..yres { + let cx = x as f32 * xinc; + let cy = y as f32 * yinc; triangles.push((Point::new(cx, cy, 1.0), - Point::new(cx + 32.0, cy, 1.0), - Point::new(cx, cy + 32.0, 1.0))); + Point::new(cx + xinc, cy, 1.1), + Point::new(cx, cy + yinc, 1.2))); + triangles.push((Point::new(cx + xinc, cy + yinc, 1.0), + Point::new(cx, cy + yinc, 1.1), + Point::new(cx + xinc, cy, 1.2))); } } triangles }; - let scene = bvh::BVH::from_triangles(&mut triangles[..]); + let scene = bvh::BVH::from_objects(&mut triangles[..], |tri| { + let minimum = tri.0.min(tri.1.min(tri.2)); + let maximum = tri.0.max(tri.1.max(tri.2)); + BBox { + min: minimum, + max: maximum, + } + }); println!("Scene built."); // Write output image of ray-traced triangle @@ -96,23 +111,29 @@ fn main() { let mut g = 0.0; let mut b = 0.0; let offset = hash_u32(((x as u32) << 16) ^ (y as u32), 0); - const SAMPLES: usize = 64; + const SAMPLES: usize = 16; for si in 0..SAMPLES { - let ray = Ray::new(Point::new(x as f32 + - fast_logit(halton::sample(0, offset + si as u32), - 1.5), - y as f32 + - fast_logit(halton::sample(3, offset + si as u32), - 1.5), - 0.0), - Vector::new(0.0, 0.0, 1.0)); - if bvh::intersect_bvh(&scene, &ray) { - r += 1.0; - g += 1.0; - b += 1.0; - // r += u; - // g += v; - // b += (1.0 - u - v).max(0.0); + let mut ray = Ray::new(Point::new(x as f32 + + fast_logit(halton::sample(0, + offset + si as u32), + 1.5), + y as f32 + + fast_logit(halton::sample(3, + offset + si as u32), + 1.5), + 0.0), + Vector::new(0.0, 0.0, 1.0)); + if let Some((_, u, v)) = bvh::intersect_bvh(&scene, &mut ray) { + r += u; + g += v; + b += (1.0 - u - v).max(0.0); + // r += 1.0; + // g += 1.0; + // b += 1.0; + } else { + r += 0.1; + g += 0.1; + b += 0.1; } } r *= 255.0 / SAMPLES as f32; diff --git a/src/math/mod.rs b/src/math/mod.rs index d4e81de..40671a3 100644 --- a/src/math/mod.rs +++ b/src/math/mod.rs @@ -40,14 +40,20 @@ pub fn fast_ln(x: f32) -> f32 { return y - 87.989971088; } -pub fn fast_logit(p: f32, width: f32) -> f32 { - let n = 0.001 + (p * 0.998); - fast_ln((n / (1.0 - n))) * width * (0.6266 / 4.0) -} +/// The logit function, scaled to approximate the probit function. +/// +/// We use this as a close approximation to the gaussian inverse CDF, +/// since the gaussian inverse CDF (probit) has no analytic formula. pub fn logit(p: f32, width: f32) -> f32 { let n = 0.001 + (p * 0.998); (n / (1.0 - n)).ln() * width * (0.6266 / 4.0) } + +pub fn fast_logit(p: f32, width: f32) -> f32 { + let n = 0.001 + (p * 0.998); + + fast_ln((n / (1.0 - n))) * width * (0.6266 / 4.0) +}