Implemented a very basic BVH for lists of triangles.

This isn't a good implementation by any means.  It's just to get
things started.

src/bbox.rs

@@ -5,6 +5,9 @@ use std::ops::BitOr;
 
 use math::Point;
 use lerp::{lerp, Lerp};
+use ray::Ray;
+
+const BBOX_MAXT_ADJUST: f32 = 1.00000024;
 
 /// A 3D axis-aligned bounding box.
 #[derive(Debug, Copy, Clone)]
@@ -32,6 +35,24 @@ impl BBox {
             max: max,
         }
     }
+
+    // Returns whether the given ray intersects with the bbox.
+    pub fn intersect_ray(&self, ray: &Ray) -> bool {
+        // Calculate slab intersections
+        let t1 = (self.min.co - ray.orig.co) * ray.dir_inv.co;
+        let t2 = (self.max.co - ray.orig.co) * ray.dir_inv.co;
+
+        // Find the far and near intersection
+        let hitt0 = (t1[0].min(t2[0]))
+                        .max(t1[1].min(t2[1]))
+                        .max(t1[2].min(t2[2]));
+        let hitt1 = (t1[0].max(t2[0]))
+                        .min(t1[1].max(t2[1]))
+                        .min(t1[2].max(t2[2]));
+
+        // Did we hit?
+        return hitt0.max(0.0) <= hitt1.min(ray.max_t);
+    }
 }
 
 

src/bvh.rs

@@ -0,0 +1,153 @@
+#![allow(dead_code)]
+
+use bbox::BBox;
+use math::Point;
+use ray::Ray;
+use triangle;
+use algorithm::partition;
+
+#[derive(Debug)]
+pub struct BVH {
+    nodes: Vec<BVHNode>,
+}
+
+#[derive(Debug)]
+enum BVHNode {
+    Internal {
+        bounds: BBox,
+        second_child_index: usize,
+    },
+
+    Leaf {
+        bounds: BBox,
+        triangle: (Point, Point, Point),
+    },
+}
+
+impl BVH {
+    pub fn from_triangles(triangles: &mut [(Point, Point, Point)]) -> BVH {
+        let mut bvh = BVH { nodes: Vec::new() };
+
+        bvh.recursive_build(triangles);
+
+        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 {
+        let me = self.nodes.len();
+
+        if triangles.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,
+            });
+        } else {
+            // Not a leaf node
+            self.nodes.push(BVHNode::Internal {
+                bounds: BBox::new(),
+                second_child_index: 0,
+            });
+
+            // 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);
+                }
+                bounds
+            };
+            let split_axis = {
+                let x_ext = bounds.max[0] - bounds.min[0];
+                let y_ext = bounds.max[1] - bounds.min[1];
+                let z_ext = bounds.max[2] - bounds.min[2];
+                if x_ext > y_ext && x_ext > z_ext {
+                    0
+                } else if y_ext > z_ext {
+                    1
+                } else {
+                    2
+                }
+            };
+            let split_pos = (bounds.min[split_axis] + bounds.max[split_axis]) * 0.5;
+
+            // Partition triangles based on split
+            let split_index = {
+                let mut split_i = partition(triangles, |tri| {
+                    let tb = tri_bounds(*tri);
+                    let centroid = (tb.min[split_axis] + tb.max[split_axis]) * 0.5;
+                    centroid < split_pos
+                });
+                if split_i < 1 {
+                    split_i = 1;
+                }
+
+                split_i
+            };
+
+            // Create child nodes
+            self.recursive_build(&mut triangles[..split_index]);
+            let child2_index = self.recursive_build(&mut triangles[split_index..]);
+
+            // Set node
+            self.nodes[me] = BVHNode::Internal {
+                bounds: bounds,
+                second_child_index: child2_index,
+            };
+        }
+
+        return me;
+    }
+}
+
+
+pub fn intersect_bvh(bvh: &BVH, ray: &Ray) -> bool {
+    let mut i_stack = [0; 64];
+    let mut stack_ptr = 0;
+
+    loop {
+        match bvh.nodes[i_stack[stack_ptr]] {
+            BVHNode::Internal { bounds, second_child_index } => {
+                if bounds.intersect_ray(ray) {
+                    i_stack[stack_ptr] += 1;
+                    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;
+                    }
+                    stack_ptr -= 1;
+                }
+            }
+        }
+    }
+
+    return false;
+}

src/main.rs

@@ -10,6 +10,7 @@ mod bbox;
 mod data_tree;
 mod image;
 mod triangle;
+mod bvh;
 mod halton;
 
 use std::path::Path;
@@ -17,7 +18,6 @@ use std::path::Path;
 use docopt::Docopt;
 
 use image::Image;
-use data_tree::DataTree;
 use math::{Point, Vector, fast_logit};
 use ray::Ray;
 
@@ -71,10 +71,24 @@ fn main() {
         return;
     }
 
+    // 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;
+                triangles.push((Point::new(cx, cy, 1.0),
+                                Point::new(cx + 32.0, cy, 1.0),
+                                Point::new(cx, cy + 32.0, 1.0)));
+            }
+        }
+        triangles
+    };
+    let scene = bvh::BVH::from_triangles(&mut triangles[..]);
+    println!("Scene built.");
+
     // Write output image of ray-traced triangle
-    let p1 = Point::new(10.0, 80.0, 1.0);
-    let p2 = Point::new(420.0, 40.0, 1.0);
-    let p3 = Point::new(235.0, 490.0, 1.0);
     let mut img = Image::new(512, 512);
     for y in 0..img.height() {
         for x in 0..img.width() {
@@ -82,7 +96,7 @@ 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 = 16;
+            const SAMPLES: usize = 64;
             for si in 0..SAMPLES {
                 let ray = Ray::new(Point::new(x as f32 +
                                               fast_logit(halton::sample(0, offset + si as u32),
@@ -92,10 +106,13 @@ fn main() {
                                                          1.5),
                                               0.0),
                                    Vector::new(0.0, 0.0, 1.0));
-                if let Some((_, u, v)) = triangle::intersect_ray(&ray, (p1, p2, p3)) {
+                if bvh::intersect_bvh(&scene, &ray) {
-                    r += u;
+                    r += 1.0;
-                    g += v;
+                    g += 1.0;
-                    b += (1.0 - u - v).max(0.0);
+                    b += 1.0;
+                    // r += u;
+                    // g += v;
+                    // b += (1.0 - u - v).max(0.0);
                 }
             }
             r *= 255.0 / SAMPLES as f32;
@@ -106,19 +123,4 @@ fn main() {
         }
     }
     let _ = img.write_binary_ppm(Path::new(&args.arg_imgpath));
-
-    let test_string = r##"
-        Thing $yar { # A comment
-            Obj [Things and stuff\]]
-        }
-
-        Thing { # A comment
-            Obj [23]
-            Obj [42]
-            Obj ["The meaning of life!"]
-        }
-    "##;
-    let tree = DataTree::from_str(test_string);
-
-    println!("{:#?}", tree);
 }

src/ray.rs

@@ -2,12 +2,14 @@
 
 use std;
 
+use float4::Float4;
 use math::{Vector, Point, Matrix4x4};
 
 #[derive(Debug, Copy, Clone)]
 pub struct Ray {
     pub orig: Point,
     pub dir: Vector,
+    pub dir_inv: Vector,
     pub max_t: f32,
     pub time: f32,
 }
@@ -17,6 +19,7 @@ impl Ray {
         Ray {
             orig: orig,
             dir: dir,
+            dir_inv: Vector { co: Float4::new(1.0, 1.0, 1.0, 1.0) / dir.co },
             max_t: std::f32::INFINITY,
             time: 0.0,
         }

todo.txt

@@ -1,2 +1 @@
-- Implement a basic camera.
+- Implement a basic camera.
-- Implement a BVH for lists of triangles.