Sick ass images

rhai/cow.rhai

@@ -1,35 +0,0 @@
-
-let scene = Scene();
-
-let distance = 3.0;
-let falloff = V(0.0,0.0,0.1);
-let camera = Camera( P(0.0,0.0,distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+Z Cam", camera);
-let camera = Camera( P(0.0,distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+Y Cam", camera);
-let camera = Camera( P(distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+X Cam", camera);
-let camera = Camera( P(0.0,0.0,-distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-Z Cam", camera);
-let camera = Camera( P(0.0,-distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-Y Cam", camera);
-let camera = Camera( P(-distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-X Cam", camera);
-
-let height = 4.0;
-let spacing = 13.0;
-let light = Light(P(0.0,height,spacing), V(0.0,0.3,0.3), falloff);
-scene.addLight("blue", light);
-let light = Light(P(0.0,height,0.0), V(0.0,0.6,0.0), falloff);
-scene.addLight("green", light);
-let light = Light(P(0.0,height,-spacing), V(0.3,0.0,0.0), falloff);
-scene.addLight("red", light);
-
-let material = Material(V(0.2,0.2,0.2), V(0.2, 0.8, 0.8), 10.0);
-scene.addMaterial("bluegreen", material);
-
-let mesh = Mesh("obj/cow.obj", material);
-let mesh_node = Node(mesh);
-scene.addNode("mesh", mesh_node);
-
-scene

rhai/mesh.rhai

@@ -1,66 +0,0 @@
-let scene = Scene();
-
-let material = Material(V(0.2,0.2,0.2), V(0.2, 0.2, 0.2), 10.0);
-scene.addMaterial("material", material);
-let material2 = Material(V(0.2,0.7,0.2), V(0.2, 0.2, 0.2), 10.0);
-scene.addMaterial("mat2", material2);
-
-let camera = Camera(P(0.0,0.0,2.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("Cam", camera);
-
-let falloff = V(0.0,0.0,0.0);
-let light = Light(P(6.0,6.0,6.0), V(0.4,0.4,0.4), falloff);
-light.active(false);
-scene.addLight("white", light);
-let light = Light(P(2.0,0.0,0.0), V(0.0,1.0,0.0), V(0.1, 0.01, 0.001));
-light.active(false);
-scene.addLight("green", light);
-let light = Light(P(-2.0,0.0,0.0), V(1.0,0.0,0.0), V(0.1, 0.01, 0.001));
-light.active(false);
-scene.addLight("red", light);
-
-let light = Ambient(V(0.3,0.3,0.3));
-scene.addLight("ambient", light);
-
-let tri = TriangleUnit();
-let tri_node = Node(tri, material);
-tri_node.active(false);
-scene.addNode("tri", tri_node);
-
-let circle = CircleUnit();
-let circle_node = Node(circle, material);
-circle_node.active(false);
-scene.addNode("circle", circle_node);
-
-let cone = ConeUnit();
-let cone_node = Node(cone, material);
-cone_node.active(false);
-scene.addNode("cone", cone_node);
-
-let torus = Torus(0.5, 1.5);
-let torus_node = Node(torus, material);
-torus_node.active(true);
-scene.addNode("torus", torus_node);
-
-
-let sphere = SphereUnit();
-let sphere_node = Node(sphere, material);
-sphere_node.translate(0.0,0.0,0.0);
-sphere_node.active(false);
-scene.addNode("sphere", sphere_node);
-
-let ground = SphereUnit();
-let ground_node = Node(ground, material2);
-let scale = 2.0;
-ground_node.translate(0.0,-scale*2.0,0.0);
-ground_node.scale(scale,scale,scale);
-ground_node.active(false);
-scene.addNode("ground", ground_node);
-
-let mesh = Mesh("obj/cat.obj");
-let mesh_node = Node(mesh, material);
-mesh_node.active(false);
-scene.addNode("mesh", mesh_node);
-
-
-scene

rhai/scene.rhai

@@ -1,66 +1,174 @@
 let scene = Scene();
 
-let distance = 10.0;
+let distance = 0.99;
 let camera = Camera( P(0.0,0.0,distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
 scene.addCamera("+Z Cam", camera);
-let camera = Camera( P(0.0,distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
+// let camera = Camera( P(0.0,distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+Y Cam", camera);
+// scene.addCamera("+Y Cam", camera);
-let camera = Camera( P(distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
+// let camera = Camera( P(distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+X Cam", camera);
+// scene.addCamera("+X Cam", camera);
-let camera = Camera( P(0.0,0.0,-distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
+// let camera = Camera( P(0.0,0.0,-distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-Z Cam", camera);
+// scene.addCamera("-Z Cam", camera);
-let camera = Camera( P(0.0,-distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
+// let camera = Camera( P(0.0,-distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-Y Cam", camera);
+// scene.addCamera("-Y Cam", camera);
-let camera = Camera( P(-distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
+// let camera = Camera( P(-distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-X Cam", camera);
+// scene.addCamera("-X Cam", camera);
 
-let material = Material(V(0.2,0.2,0.2), V(0.2, 0.8, 0.8), 10.0);
+let falloff = V(0.1, 0.1, 0.15);
-scene.addMaterial("bluegreen", material);
 
-let light = Light(P(0.0,7.0,0.0), V(0.0,0.0,1.0), V(0.1, 0.01, 0.001));
+let colour = V(0.0,0.5,0.5);
-light.active(false);
+let pos = P(-0.5,0.9,0.5);
+let light = Light(pos, colour, falloff);
+light.active(true);
 scene.addLight("blue", light);
 
-let light = Light( P(2.0,7.0,0.0), V(0.0,1.0,0.0), V(0.1, 0.01, 0.001));
+let colour = V(0.0,1.0,0.0);
+let pos = P(-0.5,0.9,-0.5);
+let light = Light(pos, colour, falloff);
 light.active(false);
 scene.addLight("green", light);
 
-let light = Light( P(2.0,7.0,2.0), V(1.0,0.0,0.0), V(0.1, 0.01, 0.001));
+let colour = V(1.0,0.0,0.0);
+let light = Light(pos, colour, falloff);
+light.active(false);
 scene.addLight("red", light);
 
+let colour = V(0.7,0.7,0.7);
+let pos = P(0.0,0.9,0.0);
+let light = Light(pos, colour, falloff);
+light.active(true);
+scene.addLight("white", light);
+
 let light = Ambient(V(0.1,0.1,0.1));
+light.active(true);
 scene.addLight("ambient", light);
 
-
-let sphere = Sphere(P(0.0,0.0,0.0), 1.0 );
-let sphere_node = Node(sphere, material);
- scene.addNode("sphere", sphere_node);
-
 //let mesh = Mesh("obj/cow.obj" );
 //let mesh_node = Node(mesh);
 //scene.addNode("mesh", mesh_node);
 
-for i in 0..6 {
-    let sphere = Sphere(P(0.0,0.0,0.0), 2.0 );
-     let sphere_node = Node(sphere, material);
-     sphere_node.translate(4.0*cos(i.to_float()), -4.0, 4.0*sin(i.to_float()));
-     scene.addNode(i.to_string(), sphere_node);
- }
 // let child = sphere_node.child(sphere);
 // child.translate(V(1.0,1.0,1.0));
 //scene.addNode(child);
 
+//let sphere2= SphereUnit();
+//let sphere2_node = Node( sphere2, material2);
+// sphere2_node.rotate(0.1,0.1,0.0);
+// sphere2_node.translate(0.0,1.0,0.0);
+//scene.addNode("sphere2",  sphere2_node);
+
+let kd = V(1.0, 1.0, 1.0); // Diffuse color (white)
+let ks = V(0.0,0.0,0.0); // Specular color (no specular reflection)
+let kr = V(0.0,0.0,0.0); // Reflection color (no reflection)
+let white_wall = Material(kd, ks, kr, 10.0);
+scene.addMaterial("white_wall", white_wall);
+
+let kd = V(1.0, 0.0, 0.0); // Diffuse color (white)
+let ks = V(0.0,0.0,0.0); // Specular color (no specular reflection)
+let kr = V(0.0,0.0,0.0); // Reflection color (no reflection)
+let red_wall = Material(kd, ks, kr, 10.0);
+scene.addMaterial("red_wall", red_wall);
+
+let kd = V(0.0, 1.0, 0.0); // Diffuse color (white)
+let ks = V(0.0,0.0,0.0); // Specular color (no specular reflection)
+let kr = V(0.0,0.0,0.0); // Reflection color (no reflection)
+let green_wall = Material(kd, ks, kr, 10.0);
+scene.addMaterial("green_wall", green_wall);
+
+let kd = V(0.0, 0.0, 1.0); // Diffuse color (white)
+let ks = V(0.0,0.0,0.0); // Specular color (no specular reflection)
+let kr = V(0.0,0.0,0.0); // Reflection color (no reflection)
+let blue_wall = Material(kd, ks, kr, 10.0);
+scene.addMaterial("blue_wall", blue_wall);
+
+
+//Rear wall
+let rectangle1 = RectangleUnit();
+let rectangle_node1 = Node(rectangle1, white_wall);
+rectangle_node1.rotate(0.0, 0.0, 0.0);
+rectangle_node1.translate(0.0, 0.0, -1.0);
+rectangle_node1.active(true);
+scene.addNode("rectangle1", rectangle_node1);
+
+//Behind wall
+// let rectangle6 = RectangleUnit();
+// let rectangle_node6 = Node(rectangle6, white_wall);
+// rectangle_node6.rotate(0.0, 180.0, 0.0);
+// rectangle_node6.translate(0.0, 0.0, 1.0);
+// rectangle_node6.active(true);
+// scene.addNode("rectangle6", rectangle_node6);
+
+//Right wall
+let rectangle2 = RectangleUnit();
+let rectangle_node2 = Node(rectangle2, green_wall);
+rectangle_node2.rotate(0.0, -90.0, 0.0);
+rectangle_node2.translate(1.0, 0.0, 0.0);
+rectangle_node2.active(true);
+scene.addNode("rectangle2", rectangle_node2);
+
+//Floor
+let rectangle3 = RectangleUnit();
+let rectangle_node3 = Node(rectangle3, red_wall);
+rectangle_node3.rotate(0.0, 90.0, 0.0);
+rectangle_node3.translate(-1.0, 0.0, 0.0);
+rectangle_node3.active(true);
+scene.addNode("rectangle3", rectangle_node3);
+
+//Left wall
+let rectangle4 = RectangleUnit();
+let rectangle_node4 = Node(rectangle4, white_wall);
+rectangle_node4.rotate(90.0, 0.0, 0.0);
+rectangle_node4.translate(0.0, 1.0, 0.0);
+rectangle_node4.active(true);
+scene.addNode("rectangle4", rectangle_node4);
+
+//Ceiling
+let rectangle5 = RectangleUnit();
+let rectangle_node5 = Node(rectangle5, white_wall);
+rectangle_node5.rotate(-90.0, 0.0, 0.0);
+
+rectangle_node5.translate(0.0, -1.0, 0.0);
+rectangle_node5.active(true);
+scene.addNode("rectangle5", rectangle_node5);
+
+
+let kd = V(0.0, 0.0, 0.0); // Diffuse color (white)
+let ks = V(0.0,0.0,0.0); // Specular color (no specular reflection)
+let kr = V(1.0,1.0,1.0); // Reflection color (no reflection)
+let reflective = Material(kd, ks, kr, 10.0);
+scene.addMaterial("reflective", reflective);
+
+let sphere = Sphere(P(0.0,0.0,0.0), 0.4 );
+let sphere_node = Node( sphere, reflective);
+sphere_node.translate(0.4, -0.6, 0.0);
+scene.addNode("sphere",sphere_node);
+
+let kd = V(0.3, 0.3, 0.3); // Diffuse color (white)
+let ks = V(0.3,0.3,0.0); // Specular color (no specular reflection)
+let kr = V(0.0,0.0,1.0); // Reflection color (no reflection)
+let shiny = Material(kd, ks, kr, 2.0);
+scene.addMaterial("shiny", shiny);
+
 let cube = CubeUnit();
-let cube_node = Node(cube, material);
+let cube_node = Node( cube, shiny);
-scene.addNode("cube", cube_node);
+cube_node.translate(-0.5,-0.6,0.0);
+cube_node.scale(0.3,0.2,0.2);
+cube_node.rotate(0.0,45.0,30.0);
+scene.addNode("cube",cube_node);
 
-//let gnonom = Gnonom();
+let gnonom = Gnonom();
-//let gnonom_node = Node(gnonom);
+let gnonom_node = Node(gnonom, shiny);
-//scene.addNode("gnonom", gnonom_node);
+gnonom_node.scale(0.2,0.2,0.2);
+gnonom_node.translate(0.0, 0.-0.7, 0.8);
+gnonom_node.rotate(0.0, 45.0, 0.0);
+gnonom_node.active(false);
+scene.addNode("gnonom", gnonom_node);
 
-//let cylinder = Cylinder(2.0,1.0 );
+// let cylinder = Cylinder(2.0, 1.0);
-//let cylinder_node = Node(cylinder);
+// let cylinder_node = Node(cylinder, material);
-//cylinder_node.scale(1.0,1.0,1.0);
+// cylinder_node.scale(1.0, 1.0, 1.0);
-//scene.addNode("cylinder",cylinder_node);
+// scene.addNode("cylinder", cylinder_node);
+
+//let cone 
 
 scene

rhai/shapes.rhai

@@ -2,7 +2,6 @@
 let scene = Scene();
 
 let distance = 3.0;
-let falloff = V(0.0,0.0,0.1);
 let camera = Camera( P(0.0,0.0,distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
 scene.addCamera("+Z Cam", camera);
 let camera = Camera( P(0.0,distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
@@ -16,8 +15,29 @@ scene.addCamera("-Y Cam", camera);
 let camera = Camera( P(-distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
 scene.addCamera("-X Cam", camera);
 
+let kd = V(0.6, 0.3, 0.9); // Diffuse color (white)
+let ks = V(0.0,0.0,0.0); // Specular color (no specular reflection)
+let kr = V(1.0,1.0,1.0); // Reflection color (no reflection)
+let material = Material(kd, ks, kr, 10.0);
+scene.addMaterial("mattee", material);
+
+let kd = V(0.3, 0.3, 0.3); // Diffuse color (white)
+let ks = V(0.7,0.7,0.7); // Specular color (no specular reflection)
+let kr = V(1.0,1.0,1.0); // Reflection color (no reflection)
+let material1 = Material(kd, ks, kr, 10.0);
+scene.addMaterial("mattee", material1);
+
+let kd = V(0.4, 0.0, 0.8); // Diffuse color (white)
+let ks = V(0.0,0.7,0.7); // Specular color (no specular reflection)
+let kr = V(1.0,1.0,1.0); // Reflection color (no reflection)
+let material2 = Material(kd, ks, kr, 10.0);
+scene.addMaterial("reflect", material2);
+
 let height = 4.0;
 let spacing = 4.0;
+
+let falloff = V(0.0,0.0,0.01);
+
 let blue = V(0.0,0.0,0.6);
 let light = Light(P(0.0,height,spacing), blue, falloff);
 scene.addLight("blue", light);
@@ -28,36 +48,43 @@ let red = V(0.6,0.0,0.0);
 let light = Light(P(0.0,height,-spacing), red, falloff);
 scene.addLight("red", light);
 
-let material = Material(V(0.2,0.2,0.2), V(0.2, 0.8, 0.8), 10.0);
-scene.addMaterial("bluegreen", material);
-
 let steiner = Steiner();
-let steiner_node = Node(steiner, material);
+let steiner_node = Node(steiner, material2);
+steiner_node.rotate(90.0,0.0,0.0);
+steiner_node.translate(0.0,0.0,1.0);
  scene.addNode("steiner", steiner_node);
 
 let steiner2 = Steiner2();
 let steiner2_node = Node(steiner2, material2);
+steiner2_node.active(false);
  scene.addNode("steiner2", steiner2_node);
 
 let crosscap = CrossCap();
 let crosscap_node = Node(crosscap, material);
+crosscap_node.active(false);
  scene.addNode("crosscap", crosscap_node);
 
-let p = 1.0;
+let p = 0.9;
-let q = 1.0;
+let q = 0.1;
 let crosscap2 = CrossCap2(p, q);
 let crosscap2_node = Node(crosscap2, material);
+crosscap2_node.active(true);
+crosscap2_node.translate(0.0,0.0,-1.5);
+crosscap2_node.rotate(140.0,0.0,90.0);
  scene.addNode("crosscap2", crosscap2_node);
 
-let k = 0.5;
+let k = 2.0;
 let roman = Roman(k );
 let roman_node = Node(roman, material);
+roman_node.active(false);
  scene.addNode("roman", roman_node);
 
 let inner_rad = 1.0;
-let outer_rad = 0.5;
+let outer_rad = 1.2;
 let torus = Torus(inner_rad, outer_rad );
 let torus_node = Node(torus, material);
+torus_node.scale(0.2,0.2,0.2);
+torus_node.rotate(0.0,70.0,0.0);
  scene.addNode("torus", torus_node);
 
 scene

rhai/space.rhai

@@ -0,0 +1,34 @@
+let scene = Scene();
+
+let falloff = V(0.0,0.0,0.0);
+
+//CAMERAS
+// let camera = Camera( P(100.0,100.0,100.0), P(500.0,500.0,500.0), V(0.0,1.0,0.0));
+// scene.addCamera("Main Camera", camera);
+
+// //Light for the sun
+// let light = Light(P(800.0, 800.0, 250.0), V(1.0, 1.0, 0.929), falloff);
+// scene.addLight("Sun", light);
+
+// //Ball for the sun 
+// let material = Material(V(1.0, 1.0, 0.9), V(0.9, 0.9, 0.9), 10.0);
+// scene.addMaterial("material_sun", material);
+// let sphere = Sphere(P(800.0, 800.0, 200.0), 50.0);
+// let sphere_node = Node(sphere, material);
+//  scene.addNode("sphere", sphere_node);
+
+// //Ball for the planet 
+// let material = Material(V(0.2,0.8,0.2), V(0.2, 0.8, 0.8), 10.0);
+// scene.addMaterial("material_planet", material);
+// let sphere = Sphere(P(500.0, 500.0, 500.0), 50.0);
+// let sphere_node = Node(sphere, material);
+//  scene.addNode("sphere", sphere_node);
+
+let material = 
+
+let steiner = Steiner();
+let steiner_node = Node(steiner, material);
+ scene.addNode("sphere", steiner);
+
+
+scene

rhai/sphere.rhai

@@ -1,34 +0,0 @@
-let scene = Scene();
-
-let distance = 3.0;
-let falloff = V(0.0,0.0,0.1);
-let camera = Camera( P(0.0,0.0,distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+Z Cam", camera);
-let camera = Camera( P(0.0,distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+Y Cam", camera);
-let camera = Camera( P(distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("+X Cam", camera);
-let camera = Camera( P(0.0,0.0,-distance), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-Z Cam", camera);
-let camera = Camera( P(0.0,-distance,0.1), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-Y Cam", camera);
-let camera = Camera( P(-distance,0.0,0.0), P(0.0,0.0,0.0), V(0.0,1.0,0.0));
-scene.addCamera("-X Cam", camera);
-
-let height = 4.0;
-let spacing = 13.0;
-let light = Light(P(0.0,height,spacing), V(0.0,0.3,0.3), falloff);
-scene.addLight("blue", light);
-let light = Light(P(0.0,height,0.0), V(0.0,0.6,0.0), falloff);
-scene.addLight("green", light);
-let light = Light(P(0.0,height,-spacing), V(0.3,0.0,0.0), falloff);
-scene.addLight("red", light);
-
-let material = Material(V(0.2,0.2,0.2), V(0.2, 0.8, 0.8), 10.0);
-scene.addMaterial("bluegreen", material);
-
-let sphere = Steiner( material);
-let sphere_node = Node(sphere);
- scene.addNode("sphere", sphere_node);
-
-scene

src/bvh.rs

@@ -1,11 +1,21 @@
-use crate::{ray::*, EPSILON};
+use crate::{node::Node, ray::*, EPSILON};
-use nalgebra::{Point3, Vector3};
+use nalgebra::{distance, Matrix4, Point3, Vector3};
+use std::collections::HashMap;
+use std::fmt;
+
+// Debuging statics
+static mut STATIC0: i32 = 0;
+static mut STATIC1: i32 = 0;
+static mut STATIC2: i32 = 0;
+static mut STATIC3: i32 = 0;
+static mut STATIC4: i32 = 0;
 
 // BOUNDING BOX -----------------------------------------------------------------
 #[derive(Clone)]
 pub struct AABB {
     pub bln: Point3<f64>,
     pub trf: Point3<f64>,
+    pub centroid: Point3<f64>,
 }
 
 impl AABB {
@@ -13,10 +23,28 @@ impl AABB {
     pub fn new(bln: Point3<f64>, trf: Point3<f64>) -> AABB {
         let bln = bln + Vector3::new(EPSILON, EPSILON, EPSILON);
         let trf = trf - Vector3::new(EPSILON, EPSILON, EPSILON);
-        AABB { bln, trf }
+        let centroid = bln + (trf - bln) / 2.0;
+        AABB { bln, trf, centroid }
+    }
+    //Empty box
+    pub fn empty() -> AABB {
+        AABB {
+            bln: Point3::new(f64::MAX, f64::MAX, f64::MAX),
+            trf: Point3::new(f64::MIN, f64::MIN, f64::MIN),
+            centroid: Point3::new(0.0, 0.0, 0.0),
+        }
+    }
+    //Apply a matrix transformation to a box
+    pub fn transform_mut(&mut self, mat: &Matrix4<f64>) {
+        let bln = &mut self.bln;
+        let trf = &mut self.trf;
+        let centroid = &mut self.centroid;
+        self.bln = mat.transform_point(bln);
+        self.trf = mat.transform_point(trf);
+        self.centroid = mat.transform_point(centroid);
     }
     // Intersect bounding box exactly
-    pub fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    pub fn intersect_ray(&self, ray: &Ray) -> bool {
         let bln = &self.bln;
         let trf = &self.trf;
         let t1 = (bln - ray.a).component_div(&ray.b);
@@ -42,7 +70,7 @@ impl AABB {
         false
     }
     // Intersect way with some epsilon term
-    pub fn intersect_bounding_box_aprox(&self, ray: &Ray) -> bool {
+    pub fn intersect_ray_aprox(&self, ray: &Ray) -> bool {
         let bln = &self.bln;
         let trf = &self.trf;
         let t1 = (bln - ray.a).component_div(&ray.b);
@@ -69,7 +97,7 @@ impl AABB {
     }
     // Get the center of this bounding box
     fn get_centroid(&self) -> Point3<f64> {
-        self.bln + (self.trf - self.bln) / 2.0
+        self.centroid
     }
     // Make a new AABB that contains both
     pub fn join(&self, other: &AABB) -> AABB {
@@ -86,6 +114,19 @@ impl AABB {
             ),
         )
     }
+    //Join mutably
+    pub fn join_mut(&mut self, other: &AABB) {
+        self.bln = Point3::new(
+            self.bln.x.min(other.bln.x),
+            self.bln.y.min(other.bln.y),
+            self.bln.z.min(other.bln.z),
+        );
+        self.trf = Point3::new(
+            self.trf.x.max(other.trf.x),
+            self.trf.y.max(other.trf.y),
+            self.trf.z.max(other.trf.z),
+        );
+    }
     //Grow the AABB to contain the cover the point
     pub fn grow(&self, other: &Point3<f64>) -> AABB {
         AABB::new(
@@ -101,34 +142,322 @@ impl AABB {
             ),
         )
     }
+    //Grow mutably
+    pub fn grow_mut(&mut self, other: &Point3<f64>) {
+        self.bln = Point3::new(
+            self.bln.x.min(other.x),
+            self.bln.y.min(other.y),
+            self.bln.z.min(other.z),
+        );
+        self.trf = Point3::new(
+            self.trf.x.max(other.x),
+            self.trf.y.max(other.y),
+            self.trf.z.max(other.z),
+        );
+    }
     // Size of AABB
     pub fn size(&self) -> Vector3<f64> {
         self.trf - self.bln
-    }
+    } //Surface area of AABB
-    //Surface area of AABB
     pub fn surface_area(&self) -> f64 {
         let size = self.size();
         2.0 * (size.x * size.y + size.x * size.z + size.y * size.z)
     }
+    pub fn area(&self) -> f64 {
+        let extent = self.trf - self.bln;
+        return extent.x * extent.y + extent.y * extent.z + extent.z * extent.x;
+    }
     // Volume of the AABB
     pub fn volume(&self) -> f64 {
         let size = self.size();
         size.x * size.y * size.z
     }
 }
-
+impl fmt::Display for AABB {
-pub enum BVHNode<'a> {
+    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
-    Leaf {
+        if self.bln[0] == f64::MAX || self.trf[0] == f64::MIN {
-        parent: &'a BVHNode<'a>,
+            writeln!(f, "Empty aabb")
-        bounding_box: AABB,
+        } else {
-        depth: u32,
+            writeln!(f, "bln: {}\ntrf: {}", self.bln, self.trf)
-    },
+        }
-    Node {
+    }
-        parent: Option<&'a BVHNode<'a>>,
+}
-        child_l: &'a BVHNode<'a>,
+#[derive(Clone)]
-        child_r: &'a BVHNode<'a>,
+pub struct BVHNode {
-        depth: u32,
+    aabb: AABB,        //The nodes bounding box
-    },
+    l_idx: usize,      //Child node l, the right node is alway l_idx + 1
+    first_prim: usize, //First primitive that the node encapsulates
+    prim_count: usize, //Number of primitives the node encapsulates
 }
 
-impl<'a> BVHNode<'a> {}
+impl BVHNode {
+    pub fn default() -> BVHNode {
+        BVHNode {
+            aabb: AABB::empty(),
+            l_idx: 0,
+            first_prim: 0,
+            prim_count: 0,
+        }
+    }
+}
+
+impl fmt::Display for BVHNode {
+    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
+        writeln!(f, "l_idx: {}", self.l_idx)?;
+        writeln!(f, "First Prim: {}", self.first_prim)?;
+        writeln!(f, "Prim Count: {}", self.prim_count)?;
+        writeln!(f, "aabb: {}", self.aabb)
+    }
+}
+
+pub struct BVH {
+    bvh_nodes: Vec<BVHNode>, //BVH nodes with AABBs
+    nodes: Vec<Node>,        //Nodes with primitives
+    nodes_used: usize,
+}
+
+impl BVH {
+    //Build a bvh by subdividing recursively
+    pub fn build(in_nodes: &HashMap<String, Node>) -> BVH {
+        /*
+        Make our own vec of nodes so that we can refer to it by index
+        This might be expensive so another method is preferred
+        */
+        let mut nodes = vec![];
+        for (_, node) in in_nodes {
+            nodes.push(node.clone());
+        }
+
+        //A BVH tree will be maximum size of 2*n + 1
+        //Initialise an empty BVHNode with empty AABB
+        let n = nodes.len();
+        let bvh_nodes: Vec<BVHNode> = vec![BVHNode::default(); 2 * n + 1];
+
+        //Begin constructing our BVH tree
+        //One node used to begin with (The root node)
+        let nodes_used = 1;
+        let mut tree = BVH {
+            nodes,
+            bvh_nodes,
+            nodes_used,
+        };
+        // Get the root node at index 0
+        let root = &mut tree.bvh_nodes[0];
+        root.l_idx = 0; //Root node has no left or right child to begin
+        (root.first_prim, root.prim_count) = (0, n); //Make root include all n nodes
+        tree.update_bvh_node_aabb(0); //Create the root nodes AABB on the n primitives
+        tree.subdivide(0); //Sub divide the root node
+        tree
+    }
+    // Will update the node's AABB at bvh_nodes[index]
+    fn update_bvh_node_aabb(&mut self, index: usize) {
+        // We will make his node bound all its primitives
+        let bvh_node = &mut self.bvh_nodes[index]; // Current BVHNode
+        let bvh_node_aabb = &mut bvh_node.aabb; //Current node AABB
+
+        let first_prim = bvh_node.first_prim; //Start index of prim
+        let prim_count = bvh_node.prim_count; //Number of primitives within the nodes aabb
+
+        for i in 0..prim_count {
+            let node = &self.nodes[first_prim + i]; //Get the node from the Vec<Node>
+            bvh_node_aabb.join_mut(&node.aabb); //Join it with the BVH node's AABB
+        }
+
+        // unsafe {
+        //     println!("UPDATE TO AABB ---- {STATIC0}");
+        //     STATIC0 += 1;
+        //     let bvh_node = &mut self.bvh_nodes[index]; //Get the BVHNode we are working on
+        //     println!("{bvh_node}");
+        // }
+    }
+    // Subdivision, will subdivide a split
+    fn subdivide(&mut self, index: usize) {
+        //Get the bvh_node we will be altering
+        // Determine the axis and position of the split plane
+        // Split the group of primitives in two halves using the split plane
+        // Create child nodes for each half
+        // Recurse into each of the child nodes.
+
+        //Leaf node case, we cannot sub-divide any more
+        if self.bvh_nodes[index].prim_count == 1 {
+            return;
+        };
+
+        /* ------------ SUBDIVIDE BY LONGEST AXIS ------------ */
+        //Get information about the node we want to subdivide
+        let (bln, trf) = (
+            self.bvh_nodes[index].aabb.bln,
+            self.bvh_nodes[index].aabb.trf,
+        );
+        let extent = trf - bln;
+        let mut axis = 0; // Assume that x is longest
+        if extent.y > extent.x {
+            axis = 1; // Split y if longest
+        };
+        if extent.z > extent[axis] {
+            axis = 2; // Split z if longest
+        };
+        let split_pos = bln[axis] + extent[axis] * 0.5; // Final split down the middle of AABB
+
+        /* --------- SUBDIVIDE BY Surface Area Heuristic ---------*/
+        // let mut best_axis: Option<usize> = None;
+        // let mut best_pos = 0.0;
+        // let mut best_cost = 1e30;
+        // let first_prim_idx = self.bvh_nodes[index].first_prim;
+        // for axis in 0..2 {
+        //     for i in 0..self.bvh_nodes[index].prim_count {
+        //         let node = &self.nodes[first_prim_idx + i];
+        //         //Get the centroid of the bounding box
+        //         let centroid = node.aabb.get_centroid();
+        //         //Get the candidate position
+        //         let candidate_pos = world_centroid[axis];
+        //         let cost = self.evaluate_sah(&self.bvh_nodes[index], axis, candidate_pos);
+        //         if cost < best_cost {
+        //             best_pos = candidate_pos;
+        //             best_axis = Some(axis);
+        //             best_cost = cost;
+        //         }
+        //     }
+        // }
+        // let axis = match best_axis {
+        //     Some(axis) => axis,
+        //     None => 0,
+        // };
+        // let split_pos = best_pos;
+
+        let left_count;
+        let right_count;
+        let mut i;
+        let mut j;
+        {
+            let bvh_node = &mut self.bvh_nodes[index];
+            i = bvh_node.first_prim; //Start of array
+            j = i + bvh_node.prim_count - 1; //End of array
+            while i <= j {
+                //Perform a quicksort dependent on location
+                let node = &self.nodes[i]; // Node we would like to sort
+                let centroid = node.aabb.get_centroid(); //Centroid of node we would like to sort
+                if centroid[axis] < split_pos {
+                    i += 1; // On Left-Hand-Side
+                } else {
+                    self.nodes.swap(i, j);
+                    j -= 1; // On Right-Hand-Side
+                }
+            }
+            //Now we have two splits
+            //The lhs of the array is in the left split  0..left_count
+            //The rhs of the array is on the right split left_count + 1..n
+            left_count = i - bvh_node.first_prim; //Number of prims on lhs
+            right_count = bvh_node.prim_count - left_count;
+            //println!("SPLIT INTO: {left_count} {right_count}");
+            if left_count == 0 || left_count == bvh_node.prim_count {
+                //Split did nothing
+                return;
+            }
+        }
+        // unsafe {
+        //     println!("SUBDIVIDE: {STATIC1}");
+        //     println!("SPLIT INTO: {left_count} ");
+        //     STATIC1 += 1;
+        // }
+
+        let l_idx = self.nodes_used; //Left child
+        self.bvh_nodes[index].l_idx = l_idx;
+        self.nodes_used = self.nodes_used + 2;
+
+        //Set left node information
+        self.bvh_nodes[l_idx].first_prim = self.bvh_nodes[index].first_prim; //Left split begins at parent split
+        self.bvh_nodes[l_idx].prim_count = left_count; // Left prims
+
+        //Set right node information
+        self.bvh_nodes[l_idx + 1].first_prim = i; // Right split start index
+        self.bvh_nodes[l_idx + 1].prim_count = right_count;
+
+        //Current node is not a leaf node
+        self.bvh_nodes[index].prim_count = 0;
+
+        self.update_bvh_node_aabb(l_idx); //Update AABB for left of split
+        self.update_bvh_node_aabb(l_idx + 1); //Update AABB for right of split
+
+        //Recurse
+        self.subdivide(l_idx); // Subdivide left index
+        self.subdivide(l_idx + 1); // SUbdivide right index
+    }
+    // Traverse the BVH, 0 will be needed to start at root node
+    pub fn traverse(&self, ray: &Ray, idx: usize) -> Option<(&Node, Intersection)> {
+        let bvh_node = &self.bvh_nodes[idx];
+        if !bvh_node.aabb.intersect_ray(&ray) {
+            // No intersection with BVH in world coordinates
+            return None;
+        }
+        if bvh_node.prim_count != 0 {
+            // Leaf node intersection
+            let node_idx = bvh_node.first_prim;
+            let node = &self.nodes[node_idx];
+            if !node.active {
+                return None;
+            }
+            if let Some(intersect) = node.intersect_ray(&ray) {
+                if intersect.distance < EPSILON {
+                    return None;
+                } else {
+                    return Some((node, intersect));
+                }
+            }
+            return None;
+        } else {
+            //Recurse down the BVH
+            //Recurse down the BVH right node
+            let intersect_l = self.traverse(ray, bvh_node.l_idx);
+            let intersect_r = self.traverse(ray, bvh_node.l_idx + 1);
+
+            match (intersect_l, intersect_r) {
+                (None, None) => return None,
+                (Some(intersect), None) => return Some(intersect),
+                (None, Some(intersect)) => return Some(intersect),
+                (Some((node_l, inter_l)), Some((node_r, inter_r))) => {
+                    //Compare intersect distance
+                    let dist_l = distance(&ray.a, &inter_l.point);
+                    let dist_r = distance(&ray.a, &inter_r.point);
+                    if dist_l < dist_r {
+                        return Some((node_l, inter_l));
+                    } else {
+                        return Some((node_r, inter_r));
+                    }
+                }
+            }
+        }
+    }
+    fn evaluate_sah(&self, node: &BVHNode, axis: usize, pos: f64) -> f64 {
+        // determine triangle counts and bounds for this split candidate
+        let mut l_aabb = AABB::empty();
+        let mut r_aabb = AABB::empty();
+        let mut l_count = 0;
+        let mut r_count = 0;
+        for i in 0..node.prim_count {
+            let aabb = self.nodes[node.first_prim + i].get_world_aabb();
+            if aabb.trf[axis] < pos {
+                l_count += 1;
+                l_aabb.grow_mut(&aabb.trf);
+            } else {
+                r_count += 1;
+                r_aabb.grow_mut(&aabb.bln);
+            }
+        }
+        let cost = l_count as f64 * l_aabb.area() + r_count as f64 * r_aabb.area();
+        match cost > 0.0 {
+            true => 0.0,
+            false => 1e30,
+        }
+    }
+}
+
+impl fmt::Display for BVH {
+    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
+        for (i, node) in self.bvh_nodes.iter().enumerate() {
+            writeln!(f, "Node: {i}")?;
+            writeln!(f, "{node}")?;
+        }
+        write!(f, "")
+    }
+}

src/gui.rs

@@ -5,7 +5,7 @@ use crate::{
     node::*,
     primitive::*,
     scene::*,
-    state::{INIT_FILE, SAVE_FILE},
+    state::{RaytracingOption, INIT_FILE, SAVE_FILE},
 };
 use imgui::*;
 use nalgebra::{Point3, Vector3};
@@ -14,14 +14,28 @@ use rhai::Engine;
 use std::time::Instant;
 
 //BUFFER CONSTANTS
-const BUFFER_PROPORTION_INIT: f32 = 0.2;
 const BUFFER_PROPORTION_MIN: f32 = 0.1;
 const BUFFER_PROPORTION_MAX: f32 = 1.0;
 
 //RAY CONSTANTS
-const RAYS_INIT: i32 = 100;
+const MIN_THREADS: u32 = 1;
-const RAYS_MIN: i32 = 100;
+const MAX_THREADS: u32 = 12;
-const RAYS_MAX: i32 = 10000;
+const RAYS_MIN: u32 = 100;
+const RAYS_MAX: u32 = 10000;
+const MIN_DEPTH: u8 = 1;
+const MAX_DEPTH: u8 = 10;
+const MIN_SAMPLES: u32 = 1;
+const MAX_SAMPLES: u32 = 10;
+const MIN_RANDOM: f64 = 100.0;
+const MAX_RANDOM: f64 = 1000.0;
+const MIN_EPSILON: f64 = 1e-11;
+const MAX_EPSILON: f64 = 1.0;
+
+//DIFFUSE CONSTANTS
+const MIN_DIFFUSE_RAYS: u8 = 1;
+const MAX_DIFFUSE_RAYS: u8 = 10;
+const MIN_DIFFUSE_COEFFICIENT: f32 = 0.0;
+const MAX_DIFFUSE_COEFFICIENT: f32 = 1.0;
 
 //MATERIAL CONSTANTS
 const MIN_D: f32 = 0.0;
@@ -47,14 +61,14 @@ const MAX_ROTATION: f64 = 180.0;
 const MAX_TRANSLATE: f64 = 10.0;
 
 // CAMERA CONSTANTS
-const MIN_FOV: f32 = 10.0;
+const MIN_FOV: f64 = 10.0;
-const MAX_FOV: f32 = 160.0;
+const MAX_FOV: f64 = 160.0;
 //const CAMERA_INIT: f32 = 5.0;
 
 /// Manages all state required for rendering Dear ImGui over `Pixels`test.
 pub enum GuiEvent {
-    BufferResize(f32, f32),
+    RaytracerOption(RaytracingOption),
-    CameraUpdate(Camera, f32),
+    CameraUpdate(Camera),
     SceneLoad(Scene),
     SaveImage(String),
 }
@@ -72,12 +86,9 @@ pub struct Gui {
     engine: Engine,
     scene: Scene,
 
-    pub ray_num: i32,
+    raytracing_option: RaytracingOption,
-
-    buffer_proportion: f32,
 
     camera: Camera,
-    camera_fov: f32,
 
     image_filename: String,
 }
@@ -122,7 +133,7 @@ impl Gui {
         let renderer = imgui_wgpu::Renderer::new(&mut imgui, device, queue, config);
 
         // Return GUI context
-        Self {
+        let mut gui = Self {
             imgui,
             platform,
             renderer,
@@ -135,14 +146,29 @@ impl Gui {
             engine: init_engine(),
             scene: Scene::empty(),
 
-            ray_num: RAYS_INIT,
+            raytracing_option: RaytracingOption::default(),
-            buffer_proportion: BUFFER_PROPORTION_INIT,
 
             camera: Camera::unit(),
-            camera_fov: 110.0,
 
             image_filename: String::from(SAVE_FILE),
+        };
+
+        // ------------ TESTING CODE (LOAD SCENE ON START) -----------------
+        match std::fs::read_to_string(&mut gui.script_filename) {
+            Ok(script) => {
+                gui.script = script;
             }
+            Err(e) => println!("{}", e),
+        }
+        match gui.engine.eval(&gui.script) {
+            Ok(scene) => {
+                gui.scene = scene;
+                gui.event = Some(GuiEvent::SceneLoad(gui.scene.clone()));
+            }
+            Err(e) => println!("{e}"),
+        }
+        // ------------ TESTING CODE (LOAD SCENE ON START) -----------------
+        gui
     }
 
     /// Prepare Dear ImGui.
@@ -183,23 +209,83 @@ impl Gui {
 
         //Raytracing options -------------------------------------------
         if CollapsingHeader::new("Raytracer").build(ui) {
-            // Numbers of rays to render
+            ui.slider(
-            ui.slider("# Rays: ", RAYS_MIN, RAYS_MAX, &mut self.ray_num);
+                "Threads",
+                MIN_THREADS,
+                MAX_THREADS,
+                &mut self.raytracing_option.threads,
+            );
+            // Numbers of rays to render per pass
+            ui.slider(
+                "Rays Per Pass",
+                RAYS_MIN,
+                RAYS_MAX,
+                &mut self.raytracing_option.pixels_per_thread,
+            );
             // Proportion of the window the buffer occupies
             ui.slider(
                 "% Buffer: ",
                 BUFFER_PROPORTION_MIN,
                 BUFFER_PROPORTION_MAX,
-                &mut self.buffer_proportion,
+                &mut self.raytracing_option.buffer_proportion,
+            );
+            //Clear colour for scene
+            ui.slider_config("Clear Colour", 0, 255)
+                .build_array(&mut self.raytracing_option.clear_color);
+            //Clear colour if no intersect
+            ui.slider_config("Pixel Clear Colour", 0, 255)
+                .build_array(&mut self.raytracing_option.pixel_clear);
+            //Ray depth slider
+            ui.slider(
+                "Ray Depth",
+                MIN_DEPTH,
+                MAX_DEPTH,
+                &mut self.raytracing_option.ray_depth,
+            );
+            //Ray samples slider
+            ui.slider(
+                "Ray Samples",
+                MIN_SAMPLES,
+                MAX_SAMPLES,
+                &mut self.raytracing_option.ray_samples,
+            );
+            //Ray randomness
+            ui.slider(
+                "Ray Randomness",
+                MIN_RANDOM,
+                MAX_RANDOM,
+                &mut self.raytracing_option.ray_randomness,
+            );
+            //Number of diffuse rays
+            ui.slider(
+                "Diffuse Rays",
+                MIN_DIFFUSE_RAYS,
+                MAX_DIFFUSE_RAYS,
+                &mut self.raytracing_option.diffuse_rays,
+            );
+            //Diffuse Coefficient
+            ui.slider(
+                "Diffuse Coefficient",
+                MIN_DIFFUSE_COEFFICIENT,
+                MAX_DIFFUSE_COEFFICIENT,
+                &mut self.raytracing_option.diffuse_coefficient,
             );
             // Fov of the buffer
-            ui.slider("fov", MIN_FOV, MAX_FOV, &mut self.camera_fov);
+            ui.slider(
+                "fov",
+                MIN_FOV,
+                MAX_FOV,
+                &mut self.raytracing_option.buffer_fov,
+            );
+            // Enable BVH
+            ui.checkbox("Enable BVH", &mut self.raytracing_option.bvh_active);
+            ui.checkbox("Enable Shadows", &mut self.raytracing_option.shadows);
+            ui.checkbox("Enable Reflections", &mut self.raytracing_option.reflect);
+            ui.checkbox("Enable Specular", &mut self.raytracing_option.specular);
+            ui.checkbox("Enable Diffuse", &mut self.raytracing_option.diffuse);
             // Apply stored changes
             if ui.button("Apply") {
-                self.event = Some(GuiEvent::BufferResize(
+                self.event = Some(GuiEvent::RaytracerOption(self.raytracing_option.clone()));
-                    self.buffer_proportion,
-                    self.camera_fov,
-                ));
             };
         }
         // CAMERA OPTIONS ----------------------------------------
@@ -214,7 +300,7 @@ impl Gui {
                 .build_array(self.camera.up.as_mut_slice());
             if ui.button("Apply Camera") {
                 println!("Camera changed");
-                self.event = Some(GuiEvent::CameraUpdate(self.camera.clone(), self.camera_fov));
+                self.event = Some(GuiEvent::CameraUpdate(self.camera.clone()));
             }
         }
         // SCRIPTING --------------------------------------------
@@ -266,9 +352,7 @@ impl Gui {
         // SCENE --------------------------------------------
         if CollapsingHeader::new("Scene").build(ui) {
             if ui.button("Update Scene") {
-                for (_, node) in &mut self.scene.nodes {
+                self.scene.compute();
-                    node.compute();
-                }
                 self.event = Some(GuiEvent::SceneLoad(self.scene.clone()));
             }
             // Edit transformation of nodes
@@ -318,7 +402,7 @@ impl Gui {
                 for (label, camera) in &self.scene.cameras {
                     if ui.button(label) {
                         self.camera = camera.clone();
-                        self.event = Some(GuiEvent::CameraUpdate(camera.clone(), self.camera_fov));
+                        self.event = Some(GuiEvent::CameraUpdate(camera.clone()));
                     }
                 }
             }
@@ -450,8 +534,8 @@ pub fn init_engine() -> Engine {
         .register_type::<Mesh>()
         .register_fn("Mesh", Mesh::from_file);
     engine
-        .register_type::<Rectangle>()
+        .register_type::<RectangleXY>()
-        .register_fn("Rectange", Rectangle::new)
+        .register_fn("Rectange", RectangleXY::new)
-        .register_fn("RectangleUnit", Rectangle::unit);
+        .register_fn("RectangleUnit", RectangleXY::unit);
     engine
 }

src/main.rs

@@ -1,10 +1,11 @@
 use crate::state::run;
 use error_iter::ErrorIter;
 
-const EPSILON: f64 = 1e-8;
+const EPSILON: f64 = 1e-7;
 const INFINITY: f64 = 1e10;
 
 use log::error;
+//use nalgebra::{Matrix4, RowVector4, Vector3, Vector4};
 use std::env;
 use std::error::Error;
 
@@ -21,13 +22,159 @@ mod state;
 
 fn main() {
     env_logger::init();
-
     env::set_var("RUST_BACKTRACE", "1");
+    //let args: Vec<String> = env::args().collect();
+
+    // let vec = Vector3::new(1.0, 1.0, 1.0);
+    // let translation = Vector3::new(1.0, 1.0, 1.0);
+    // let translation_matrix = Matrix4::new_translation(&translation);
+    // println!(
+    //     "{}, {}",
+    //     translation_matrix,
+    //     translation_matrix.transform_vector(&vec)
+    // );
+    // let mut translation_matrix = translation_matrix.transpose();
+    // translation_matrix.set_row(3, &RowVector4::new(0.0, 0.0, 0.0, 0.0));
+    // println!(
+    //     "{}, {}", //     translation_matrix, //     translation_matrix.transform_vector(&vec)
+    // );
     if let Err(e) = run() {
         println!("Error at runtime: {}", e);
     };
+
+    // if args.len() == 6 {
+    //     let width: usize = args[1].parse().unwrap();
+    //     let height: usize = args[2].parse().unwrap();
+    //     let fovy = args[3].parse::<f64>().unwrap();
+    //     let filename = &args[4];
+    //     let savefile = &args[5];
+    //     headless(
+    //         width,
+    //         height,
+    //         fovy,
+    //         filename.to_string(),
+    //         savefile.to_string(),
+    //     );
+    // } else {
+    //}
 }
 
+// fn headless(width: usize, height: usize, fovy: f64, filename: String, savefile: String) {
+//     let options = Arc::new(RaytracingOption {
+//         threads: 12,
+//         ray_samples: 1,
+//         ray_randomness: 100.0,
+//         clear_color: [0x22, 0x00, 0x11, 0x55],
+//         pixel_clear: [0x55, 0x00, 0x22, 0x55],
+//         pixels_per_thread: 200,
+//         buffer_proportion: 1.0,
+//         buffer_fov: 110.0,
+//         ray_depth: 5,
+//         diffuse_rays: 3,
+//         diffuse_coefficient: 0.8,
+//         bvh_active: false,
+//     });
+//     //Read script from file
+//     let script = match std::fs::read_to_string(&filename) {
+//         Ok(in_script) => in_script,
+//         Err(e) => {
+//             println!("{}", e);
+//             return;
+//         }
+//     };
+//     //Evaluate scene in file
+//     let engine = init_engine();
+//     let scene: Arc<Scene> = match engine.eval(&script) {
+//         Ok(in_scene) => Arc::new(in_scene),
+//         Err(e) => {
+//             println!("{e}");
+//             return;
+//         }
+//     };
+//     //Set the camera
+//     let mut camera = Camera::unit();
+//     for (_, in_camera) in &scene.cameras {
+//         camera = in_camera.clone();
+//     }
+//     //Cast the rays
+//     let rays = Arc::new(Ray::cast_rays(
+//         &camera.eye,
+//         &camera.target,
+//         &camera.up,
+//         fovy,
+//         width as u32,
+//         height as u32,
+//     ));
+//     //Enable bounding volume heirarchy
+//     let bvh;
+//     match options.bvh_active {
+//         true => bvh = Arc::new(Some(BVH::build(&scene.nodes))),
+//         false => bvh = Arc::new(None),
+//     }
+//     //Create our frame and indexer
+//     let size = width * height;
+//     let frame_mutex = Arc::new(Mutex::new(vec![0; size * 4]));
+//     //Multithreading
+//     let mut handles = vec![];
+
+//     for index in 0..size {
+//         for _ in 0..options.threads {
+//             //Get random index from queue
+//             //Create a nre thread for this pixel
+//             let handle = thread::spawn({
+//                 let rays = rays.clone();
+//                 let scene = scene.clone();
+//                 let options = options.clone();
+//                 let bvh = bvh.clone();
+//                 let rays = rays.clone();
+//                 let frame_mutex = frame_mutex.clone();
+//                 move || {
+//                     //Shade colour for selected ray
+//                     let mut colour: Vector3<f32> = Vector3::zeros();
+//                     //Get the ray we want to make
+//                     let shot_ray = &rays[index];
+//                     //Send out ray_samples rays
+//                     for _ in 0..options.ray_samples {
+//                         let point = shot_ray.a;
+//                         let dir = shot_ray.b;
+//                         //Generate a random ray
+//                         let rx = (random::<f64>() - 0.5) / options.ray_randomness;
+//                         let ry = (random::<f64>() - 0.5) / options.ray_randomness;
+//                         let rz = (random::<f64>() - 0.5) / options.ray_randomness;
+//                         let nx = dir.x + rx;
+//                         let ny = dir.y + ry;
+//                         let nz = dir.z + rz;
+//                         let rand_ray = Ray::new(point, Vector3::new(nx, ny, nz));
+
+//                         if let Some(ray_colour) = rand_ray.shade_ray(&scene, 0, &options, &bvh) {
+//                             colour += ray_colour;
+//                         }
+//                     }
+//                     colour = (colour / options.ray_samples as f32) * 255.0;
+//                     let rgba = [colour.x as u8, colour.y as u8, colour.z as u8, 0xff];
+//                     {
+//                         let frame = &mut frame_mutex.lock().unwrap();
+//                         frame[index * 4..(index + 1) * 4].copy_from_slice(&rgba);
+//                     }
+//                 }
+//             });
+//             handles.push(handle);
+//         }
+//         for handle in handles.drain(..) {
+//             handle.join().unwrap();
+//         }
+//     }
+//     use std::path::Path;
+//     image::save_buffer(
+//         Path::new(&savefile),
+//         &frame_mutex.lock().unwrap(),
+//         width as u32,
+//         height as u32,
+//         image::ColorType::Rgba8,
+//     )
+//     .unwrap();
+// }
+
 fn log_error<E: Error + 'static>(method_name: &str, err: E) {
     error!("{method_name}() failed: {err}");
     for source in err.sources().skip(1) {

src/material.rs

@@ -10,10 +10,10 @@ pub struct Material {
 }
 
 impl Material {
-    pub fn new(kd: Vector3<f64>, ks: Vector3<f64>, shininess: f64) -> Material {
+    pub fn new(kd: Vector3<f64>, ks: Vector3<f64>, kr: Vector3<f64>, shininess: f64) -> Material {
         let kd = kd.cast();
         let ks = ks.cast();
-        let kr = ks.cast();
+        let kr = kr.cast();
         let shininess = shininess as f32;
         Material {
             kd,

src/node.rs

@@ -1,12 +1,19 @@
-use crate::{material::Material, primitive::*};
+use crate::{
-use nalgebra::{Matrix4, Vector3};
+    bvh::AABB,
-use std::rc::Rc;
+    material::Material,
+    primitive::*,
+    ray::{Intersection, Ray},
+    EPSILON,
+};
+use nalgebra::{distance, Matrix3, Matrix4, Vector3};
+use std::sync::Arc;
 
 #[derive(Clone)]
 pub struct Node {
     //Primitive
-    pub primitive: Rc<dyn Primitive>,
+    pub primitive: Arc<dyn Primitive>,
     pub material: Material,
+    pub aabb: AABB,
     //Transformations
     pub rotation: [f64; 3],
     pub scale: [f64; 3],
@@ -14,27 +21,30 @@ pub struct Node {
     //Model matricies
     pub model: Matrix4<f64>,
     pub inv_model: Matrix4<f64>,
-
+    pub inv_transpose_model: Matrix3<f64>,
+    //If the node is active
     pub active: bool,
 }
 
 impl Node {
     //New node with no transformations
-    pub fn new(primitive: Rc<dyn Primitive>, material: Material) -> Node {
+    pub fn new(primitive: Arc<dyn Primitive>, material: Material) -> Node {
+        let aabb = primitive.get_aabb();
         Node {
             primitive,
             material,
+            aabb,
             rotation: [0.0, 0.0, 0.0],
             scale: [1.0, 1.0, 1.0],
             translation: [0.0, 0.0, 0.0],
             model: Matrix4::identity(),
             inv_model: Matrix4::identity(),
-
+            inv_transpose_model: Matrix3::identity(),
             active: true,
         }
     }
     //New node with parent transformations
-    pub fn child(self, primitive: Rc<dyn Primitive>) -> Node {
+    pub fn child(self, primitive: Arc<dyn Primitive>) -> Node {
         let mut child = self.clone();
         child.primitive = primitive;
         child
@@ -46,12 +56,6 @@ impl Node {
 
     //Rotate a mesh by adding to its rotation
     pub fn rotate(&mut self, roll: f64, pitch: f64, yaw: f64) {
-        //Convert to radians
-        let roll = roll.to_radians();
-        // Convert pitch and yaw to radians
-        let pitch = pitch.to_radians();
-        let yaw = yaw.to_radians();
-
         // Add the roll, pitch, and yaw to the current rotation
         self.rotation[0] += roll;
         self.rotation[1] += pitch;
@@ -71,9 +75,9 @@ impl Node {
     }
     // Scale a mesh by adding to its current scale
     pub fn scale(&mut self, x: f64, y: f64, z: f64) {
-        self.scale[0] += x;
+        self.scale[0] = x;
-        self.scale[1] += y;
+        self.scale[1] = y;
-        self.scale[2] += z;
+        self.scale[2] = z;
 
         // Recompute the model and inverse model matrices
         self.compute();
@@ -88,10 +92,30 @@ impl Node {
         let scale_matrix = Matrix4::new_nonuniform_scaling(&scale);
         // Rotation matrix
         let (roll, pitch, yaw) = (self.rotation[0], self.rotation[1], self.rotation[2]);
-        let rotation_matrix = Matrix4::from_euler_angles(roll, pitch, yaw);
+        let rotation_matrix =
+            Matrix4::from_euler_angles(roll.to_radians(), pitch.to_radians(), yaw.to_radians());
         // Compute the model matrix by combining the translation, rotation, and scale matrices
         self.model = (translation_matrix * rotation_matrix * scale_matrix).cast();
         // Compute the inverse model matrix by inverting the model matrix
         self.inv_model = self.model.try_inverse().unwrap();
+        self.inv_transpose_model = self.inv_model.transpose().remove_row(3).remove_column(3);
+        self.aabb.transform_mut(&self.model);
+    }
+    // Intersection of a ray, will convert to model coords and check
+    pub fn intersect_ray(&self, ray: &Ray) -> Option<Intersection> {
+        let ray = ray.transform(&self.inv_model); //Transform from world coordinates
+        if let Some(mut intersect) = self.primitive.intersect_ray(&ray) {
+            if intersect.distance < EPSILON {
+                return None;
+            }
+            intersect.transform_mut(&self.model, &self.inv_transpose_model); //Transform to world coords
+            intersect.distance = distance(&intersect.point, &ray.a);
+            return Some(intersect);
+        }
+        return None;
+    }
+    //Gets the bounding box in world coords
+    pub fn get_world_aabb(&self) -> AABB {
+        return self.aabb.clone();
     }
 }

src/primitive.rs

@@ -9,11 +9,11 @@ use nalgebra::{distance, Point3, Vector3};
 use roots::{find_roots_quadratic, find_roots_quartic, Roots};
 use std::fs::File;
 use std::io::{BufRead, BufReader};
-use std::rc::Rc;
+use std::sync::Arc;
 // PRIMITIVE TRAIT -----------------------------------------------------------------
-pub trait Primitive {
+pub trait Primitive: Send + Sync {
     fn intersect_ray(&self, ray: &Ray) -> Option<Intersection>;
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool;
+    fn get_aabb(&self) -> AABB;
 }
 
 // SPHERE -----------------------------------------------------------------
@@ -21,23 +21,14 @@ pub trait Primitive {
 pub struct Sphere {
     position: Point3<f64>,
     radius: f64,
-    bounding_box: AABB,
 }
 
 impl Sphere {
-    pub fn new(position: Point3<f64>, radius: f64) -> Rc<dyn Primitive> {
+    pub fn new(position: Point3<f64>, radius: f64) -> Arc<dyn Primitive> {
-        let radius_vec = Vector3::new(radius, radius, radius);
+        Arc::new(Sphere { position, radius })
-        let bln = position - radius_vec;
-        let trf = position + radius_vec;
-        let bounding_box = AABB::new(bln, trf);
-        Rc::new(Sphere {
-            position,
-            radius,
-            bounding_box,
-        })
     }
 
-    pub fn unit() -> Rc<dyn Primitive> {
+    pub fn unit() -> Arc<dyn Primitive> {
         Sphere::new(Point3::new(0.0, 0.0, 0.0), 1.0)
     }
 }
@@ -78,8 +69,12 @@ impl Primitive for Sphere {
         })
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        return self.bounding_box.intersect_bounding_box(ray);
+        let radius = self.radius;
+        let radius_vec = Vector3::new(radius, radius, radius);
+        let bln = self.position - radius_vec;
+        let trf = self.position + radius_vec;
+        AABB::new(bln, trf)
     }
 }
 
@@ -90,28 +85,32 @@ pub struct Circle {
     radius: f64,
     normal: Vector3<f64>,
     constant: f64,
-    bounding_box: AABB,
 }
 
 impl Circle {
-    pub fn new(position: Point3<f64>, radius: f64, normal: Vector3<f64>) -> Rc<dyn Primitive> {
+    pub fn new(position: Point3<f64>, radius: f64, normal: Vector3<f64>) -> Arc<dyn Primitive> {
-        let radius_vec = Vector3::new(radius, radius, radius);
-        let bln = position - radius_vec;
-        let trf = position + radius_vec;
-        let bounding_box = AABB::new(bln, trf);
-
         let normal = normal.normalize();
         let constant = normal.dot(&position.coords);
-        Rc::new(Circle {
+        Arc::new(Circle {
             position,
             radius,
             normal,
             constant,
-            bounding_box,
         })
     }
 
-    pub fn unit() -> Rc<dyn Primitive> {
+    pub fn new_unboxed(position: Point3<f64>, radius: f64, normal: Vector3<f64>) -> Circle {
+        let normal = normal.normalize();
+        let constant = normal.dot(&position.coords);
+        Circle {
+            position,
+            radius,
+            normal,
+            constant,
+        }
+    }
+
+    pub fn unit() -> Arc<dyn Primitive> {
         let position = Point3::new(0.0, 0.0, 0.0);
         let normal = Vector3::new(0.0, 0.0, -1.0);
         let radius = 1.0;
@@ -144,8 +143,13 @@ impl Primitive for Circle {
         }
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let radius = self.radius;
+        let position = self.position;
+        let radius_vec = Vector3::new(radius, radius, radius);
+        let bln = position - radius_vec;
+        let trf = position + radius_vec;
+        AABB::new(bln, trf)
     }
 }
 
@@ -154,31 +158,27 @@ impl Primitive for Circle {
 pub struct Cylinder {
     radius: f64,
     height: f64,
-    base_circle: Rc<dyn Primitive>,
+    base_circle: Circle,
-    top_circle: Rc<dyn Primitive>,
+    top_circle: Circle,
-    bounding_box: AABB,
 }
 
 impl Cylinder {
-    pub fn new(radius: f64, height: f64) -> Rc<dyn Primitive> {
+    pub fn new(radius: f64, height: f64) -> Arc<dyn Primitive> {
-        let base_circle = Circle::new(
+        let base_circle = Circle::new_unboxed(
             Point3::new(0.0, 0.0, 0.0),
             radius,
             Vector3::new(0.0, -1.0, 0.0),
         );
-        let top_circle = Circle::new(
+        let top_circle = Circle::new_unboxed(
             Point3::new(0.0, height, 0.0),
             radius,
             Vector3::new(0.0, 1.0, 0.0),
         );
-        let bln = Point3::new(-radius, 0.0, -radius);
+        Arc::new(Cylinder {
-        let trf = Point3::new(radius, height, radius);
-        Rc::new(Cylinder {
             radius,
             height,
             base_circle,
             top_circle,
-            bounding_box: AABB { bln, trf },
         })
     }
 }
@@ -261,8 +261,12 @@ impl Primitive for Cylinder {
         }
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let radius = self.radius;
+        let height = self.height;
+        let bln = Point3::new(-radius, 0.0, -radius);
+        let trf = Point3::new(radius, height, radius);
+        AABB::new(bln, trf)
     }
 }
 
@@ -271,28 +275,24 @@ impl Primitive for Cylinder {
 pub struct Cone {
     height: f64,
     constant: f64,
-    circle: Rc<dyn Primitive>,
+    circle: Circle,
-    bounding_box: AABB,
 }
 
 impl Cone {
-    pub fn new(radius: f64, height: f64) -> Rc<dyn Primitive> {
+    pub fn new(radius: f64, height: f64) -> Arc<dyn Primitive> {
-        let circle = Circle::new(
+        let circle = Circle::new_unboxed(
             Point3::new(0.0, 0.0, 0.0),
             radius,
             Vector3::new(0.0, -1.0, 0.0),
         );
-        let bln = Point3::new(-radius, 0.0, -radius);
-        let trf = Point3::new(radius, height, radius);
         let constant = radius * radius / (height * height);
-        Rc::new(Cone {
+        Arc::new(Cone {
             height,
             constant,
             circle,
-            bounding_box: AABB { bln, trf },
         })
     }
-    pub fn unit() -> Rc<dyn Primitive> {
+    pub fn unit() -> Arc<dyn Primitive> {
         Cone::new(0.5, 1.0)
     }
 
@@ -363,86 +363,59 @@ impl Primitive for Cone {
         }
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let height = self.height;
+        let radius = (self.constant * height * height).sqrt();
+        let bln = Point3::new(-radius, 0.0, -radius);
+        let trf = Point3::new(radius, height, radius);
+        AABB::new(bln, trf)
     }
 }
 
 // RECTANGLE -----------------------------------------------------------------
+// Normal is (0.0, 0.0, 1.0) always facing towards camera at positive z axis
 #[derive(Clone)]
-pub struct Rectangle {
+pub struct RectangleXY {
-    position: Point3<f64>,
+    bl: Point3<f64>,
-    normal: Vector3<f64>,
+    tr: Point3<f64>,
-    width_direction: Vector3<f64>,
-    width: f64,
-    height: f64,
-    bounding_box: AABB,
 }
 
-impl Rectangle {
+impl RectangleXY {
-    pub fn new(
+    pub fn new(bl: Point3<f64>, tr: Point3<f64>) -> Arc<dyn Primitive> {
-        position: Point3<f64>,
+        Arc::new(RectangleXY { bl, tr })
-        normal: Vector3<f64>,
-        width_direction: Vector3<f64>,
-        width: f64,
-        height: f64,
-    ) -> Rc<dyn Primitive> {
-        let normal = normal.normalize();
-        let width_direction = width_direction.normalize();
-        let height_direction = width_direction.cross(&normal);
-        let bln = position - width / 2.0 * width_direction - height / 2.0 * height_direction;
-        let trf = position + width / 2.0 * width_direction + height / 2.0 * height_direction;
-        Rc::new(Rectangle {
-            position,
-            normal: normal.normalize(),
-            width_direction: width_direction.normalize(),
-            width,
-            height,
-            bounding_box: AABB { bln, trf },
-        })
     }
-    pub fn unit() -> Rc<dyn Primitive> {
+    pub fn unit() -> Arc<dyn Primitive> {
-        Rectangle::new(
+        RectangleXY::new(Point3::new(-1.0, -1.0, 0.0), Point3::new(1.0, 1.0, 0.0))
-            Point3::new(0.0, 0.0, 0.0),
-            Vector3::new(0.0, 1.0, 0.0),
-            Vector3::new(1.0, 0.0, 0.0),
-            2.0,
-            2.0,
-        )
     }
 }
 
-impl Primitive for Rectangle {
+impl Primitive for RectangleXY {
     fn intersect_ray(&self, ray: &Ray) -> Option<Intersection> {
-        let constant = self.position.coords.dot(&self.normal);
+        let z = self.bl.z;
-        let denominator = ray.b.dot(&self.normal);
+        let az = ray.a.z;
-        let t = (constant - ray.a.coords.dot(&self.normal)) / denominator;
+        let bz = ray.b.z;
-
+        let t = (z - az) / bz;
         if t > INFINITY {
             return None;
         }
-
         let intersect = ray.at_t(t);
-        let height_direction = self.width_direction.cross(&self.normal);
+        let (ix, iy) = (intersect.x, intersect.y);
-        let (w2, h2) = (self.width / 2.0, self.height / 2.0);
-        let r1 = w2 * self.width_direction;
-        let r2 = h2 * height_direction;
-        let pi = intersect - self.position;
-        let pi_dot_r1 = pi.dot(&r1);
-        let pi_dot_r2 = pi.dot(&r2);
 
-        if pi_dot_r1 >= -w2 && pi_dot_r1 <= w2 && pi_dot_r2 >= -h2 && pi_dot_r2 <= h2 {
+        if (ix < self.bl.x) || (ix > self.tr.x) || (iy < self.bl.y) || (iy > self.tr.y) {
-            return Some(Intersection {
+            return None;
+        }
+
+        Some(Intersection {
             point: intersect,
-                normal: self.normal,
+            normal: Vector3::new(0.0, 0.0, 1.0),
             distance: t,
-            });
+        })
-        }
-        None
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let bl = self.bl + Vector3::new(0.0, 0.0, -0.1);
+        let tr = self.tr + Vector3::new(0.0, 0.0, 0.1);
+        AABB::new(bl, tr)
     }
 }
 
@@ -451,19 +424,18 @@ impl Primitive for Rectangle {
 pub struct Cube {
     bln: Point3<f64>,
     trf: Point3<f64>,
-    bounding_box: AABB,
 }
 
 impl Cube {
-    pub fn new(bln: Point3<f64>, trf: Point3<f64>) -> Rc<dyn Primitive> {
+    pub fn new(bln: Point3<f64>, trf: Point3<f64>) -> Arc<dyn Primitive> {
-        Rc::new(Cube {
+        Arc::new(Cube { bln, trf })
-            bln,
-            trf,
-            bounding_box: AABB { bln, trf },
-        })
     }
 
-    pub fn unit() -> Rc<dyn Primitive> {
+    pub fn new_unboxed(bln: Point3<f64>, trf: Point3<f64>) -> Cube {
+        Cube { bln, trf }
+    }
+
+    pub fn unit() -> Arc<dyn Primitive> {
         let bln = Point3::new(-1.0, -1.0, -1.0);
         let trf = Point3::new(1.0, 1.0, 1.0);
         Cube::new(bln, trf)
@@ -488,12 +460,12 @@ impl Primitive for Cube {
             let intersect = ray.at_t(tmin);
 
             // Check if the intersection is outside the box
-            if intersect.x < bln.x
+            if intersect.x < bln.x - EPSILON
-                || intersect.x > trf.x
+                || intersect.x > trf.x + EPSILON
-                || intersect.y < bln.y
+                || intersect.y < bln.y - EPSILON
-                || intersect.y > trf.y
+                || intersect.y > trf.y + EPSILON
-                || intersect.z < bln.z
+                || intersect.z < bln.z - EPSILON
-                || intersect.z > trf.z
+                || intersect.z > trf.z + EPSILON
             {
                 return None; // Intersection is outside the box
             }
@@ -516,7 +488,7 @@ impl Primitive for Cube {
 
             Some(Intersection {
                 point: intersect,
-                normal: normal,
+                normal,
                 distance: tmin,
             })
         } else {
@@ -524,8 +496,8 @@ impl Primitive for Cube {
         }
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        AABB::new(self.bln, self.trf)
     }
 }
 
@@ -539,27 +511,17 @@ pub struct Triangle {
     v: Point3<f64>,
     w: Point3<f64>,
     normal: Vector3<f64>,
-    bounding_box: AABB,
 }
 
 impl Triangle {
-    pub fn new(u: Point3<f64>, v: Point3<f64>, w: Point3<f64>) -> Rc<dyn Primitive> {
+    pub fn new(u: Point3<f64>, v: Point3<f64>, w: Point3<f64>) -> Arc<dyn Primitive> {
         let uv = v - u;
         let uw = w - u;
         let normal = uw.cross(&uv).normalize();
-        let bln = u.inf(&v).inf(&w);
+        Arc::new(Triangle { u, v, w, normal })
-        let trf = u.sup(&v).sup(&w);
-        let bounding_box = AABB { bln, trf };
-        Rc::new(Triangle {
-            u,
-            v,
-            w,
-            normal,
-            bounding_box,
-        })
     }
     #[allow(dead_code)]
-    pub fn unit() -> Rc<dyn Primitive> {
+    pub fn unit() -> Arc<dyn Primitive> {
         let u = Point3::new(-1.0, -1.0, 0.0);
         let v = Point3::new(0.0, 1.0, 0.0);
         let w = Point3::new(1.0, -1.0, 0.0);
@@ -608,8 +570,13 @@ impl Primitive for Triangle {
         None
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let u = self.u;
+        let v = self.v;
+        let w = self.w;
+        let bln = u.inf(&v).inf(&w);
+        let trf = u.sup(&v).sup(&w);
+        AABB::new(bln, trf)
     }
 }
 
@@ -617,17 +584,13 @@ impl Primitive for Triangle {
 #[derive(Clone)]
 pub struct Mesh {
     triangles: Vec<Triangle>,
-    bounding_box: AABB,
 }
 
 impl Mesh {
-    pub fn new(triangles: Vec<Triangle>) -> Rc<dyn Primitive> {
+    pub fn new(triangles: Vec<Triangle>) -> Arc<dyn Primitive> {
         // Calculate the bounding box for the entire mesh based on the bounding boxes of individual triangles
-        let bounding_box = Mesh::compute_bounding_box(&triangles);
+        let _bounding_box = Mesh::compute_bounding_box(&triangles);
-        Rc::new(Mesh {
+        Arc::new(Mesh { triangles })
-            triangles,
-            bounding_box,
-        })
     }
 
     fn compute_bounding_box(triangles: &Vec<Triangle>) -> AABB {
@@ -641,10 +604,10 @@ impl Mesh {
             trf = trf.sup(&triangle.v);
             trf = trf.sup(&triangle.w);
         }
-        AABB { bln, trf }
+        AABB::new(bln, trf)
     }
 
-    pub fn from_file(filename: &str) -> Rc<dyn Primitive> {
+    pub fn from_file(filename: &str) -> Arc<dyn Primitive> {
         let mut triangles: Vec<Triangle> = Vec::new();
         let mut vertices: Vec<Point3<f64>> = Vec::new();
 
@@ -685,16 +648,7 @@ impl Mesh {
                                 let uv = u - v;
                                 let uw = w - v;
                                 let normal = uv.cross(&uw).normalize();
-                                let bln = u.inf(&v).inf(&w);
+                                triangles.push(Triangle { u, v, w, normal });
-                                let trf = u.sup(&v).sup(&w);
-                                let bounding_box = AABB { bln, trf };
-                                triangles.push(Triangle {
-                                    u,
-                                    v,
-                                    w,
-                                    normal,
-                                    bounding_box,
-                                });
                             }
                         }
                         _ => {}
@@ -727,8 +681,8 @@ impl Primitive for Mesh {
         closest_intersect
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        Mesh::compute_bounding_box(&self.triangles)
     }
 }
 
@@ -737,18 +691,14 @@ impl Primitive for Mesh {
 pub struct Torus {
     inner_rad: f64,
     outer_rad: f64,
-    bounding_box: AABB,
 }
 
 impl Torus {
-    pub fn new(inner_rad: f64, outer_rad: f64) -> Rc<dyn Primitive> {
+    pub fn new(inner_rad: f64, outer_rad: f64) -> Arc<dyn Primitive> {
         // I need to find the bounding box for this shape
-        let trf = Point3::new(1.0, 1.0, 1.0);
+        Arc::new(Torus {
-        let bln = Point3::new(-1.0, -1.0, -1.0);
-        Rc::new(Torus {
             inner_rad,
             outer_rad,
-            bounding_box: AABB { bln, trf },
         })
     }
 }
@@ -854,53 +804,42 @@ impl Primitive for Torus {
         })
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        //TODO!
+        let trf = Point3::new(1.0, 1.0, 1.0);
+        let bln = Point3::new(-1.0, -1.0, -1.0);
+        AABB::new(bln, trf)
     }
 }
 
 // GNOMON -----------------------------------------------------------------
 #[derive(Clone)]
 pub struct Gnonom {
-    x_cube: Rc<dyn Primitive>,
+    x_cube: Cube,
-    y_cube: Rc<dyn Primitive>,
+    y_cube: Cube,
-    z_cube: Rc<dyn Primitive>,
+    z_cube: Cube,
-    bounding_box: AABB,
 }
 
 impl Gnonom {
     const GNONOM_WIDTH: f64 = 0.1;
     const GNONOM_LENGTH: f64 = 2.0;
-    pub fn new() -> Rc<dyn Primitive> {
+    pub fn new() -> Arc<dyn Primitive> {
-        let x_cube = Cube::new(
+        let x_cube = Cube::new_unboxed(
             Point3::new(0.0, -Self::GNONOM_WIDTH, -Self::GNONOM_WIDTH),
             Point3::new(Self::GNONOM_LENGTH, Self::GNONOM_WIDTH, Self::GNONOM_WIDTH),
         );
-        let y_cube = Cube::new(
+        let y_cube = Cube::new_unboxed(
             Point3::new(-Self::GNONOM_WIDTH, 0.0, -Self::GNONOM_WIDTH),
             Point3::new(Self::GNONOM_WIDTH, Self::GNONOM_LENGTH, Self::GNONOM_WIDTH),
         );
-        let z_cube = Cube::new(
+        let z_cube = Cube::new_unboxed(
             Point3::new(-Self::GNONOM_WIDTH, -Self::GNONOM_WIDTH, 0.0),
             Point3::new(Self::GNONOM_WIDTH, Self::GNONOM_WIDTH, Self::GNONOM_LENGTH),
         );
-        let bounding_box = AABB {
+        Arc::new(Gnonom {
-            bln: Point3::new(
-                -Self::GNONOM_WIDTH,
-                -Self::GNONOM_WIDTH,
-                -Self::GNONOM_WIDTH,
-            ),
-            trf: Point3::new(
-                Self::GNONOM_LENGTH,
-                Self::GNONOM_LENGTH,
-                Self::GNONOM_LENGTH,
-            ),
-        };
-        Rc::new(Gnonom {
             x_cube,
             y_cube,
             z_cube,
-            bounding_box,
         })
     }
 }
@@ -922,25 +861,30 @@ impl Primitive for Gnonom {
         None
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        AABB::new(
+            Point3::new(
+                -Self::GNONOM_WIDTH,
+                -Self::GNONOM_WIDTH,
+                -Self::GNONOM_WIDTH,
+            ),
+            Point3::new(
+                Self::GNONOM_LENGTH,
+                Self::GNONOM_LENGTH,
+                Self::GNONOM_LENGTH,
+            ),
+        )
     }
 }
 
 // CROSS CAP ---------
 #[derive(Clone)]
-pub struct CrossCap {
+pub struct CrossCap {}
-    bounding_box: AABB,
-}
 
 impl CrossCap {
-    pub fn new() -> Rc<dyn Primitive> {
+    pub fn new() -> Arc<dyn Primitive> {
         // I need to find the bounding box for this shape
-        let trf = Point3::new(1.0, 1.0, 1.0);
+        Arc::new(CrossCap {})
-        let bln = Point3::new(-1.0, -1.0, -1.0);
-        Rc::new(CrossCap {
-            bounding_box: AABB { bln, trf },
-        })
     }
 }
 
@@ -1014,8 +958,10 @@ impl Primitive for CrossCap {
         })
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let trf = Point3::new(1.0, 1.0, 1.0);
+        let bln = Point3::new(-1.0, -1.0, -1.0);
+        AABB::new(bln, trf)
     }
 }
 
@@ -1024,19 +970,12 @@ impl Primitive for CrossCap {
 pub struct CrossCap2 {
     p: f64,
     q: f64,
-    bounding_box: AABB,
 }
 
 impl CrossCap2 {
-    pub fn new(p: f64, q: f64) -> Rc<dyn Primitive> {
+    pub fn new(p: f64, q: f64) -> Arc<dyn Primitive> {
         // I need to find the bounding box for this shape
-        let trf = Point3::new(1.0, 1.0, 1.0);
+        Arc::new(CrossCap2 { p, q })
-        let bln = Point3::new(-1.0, -1.0, -1.0);
-        Rc::new(CrossCap2 {
-            p,
-            q,
-            bounding_box: AABB { bln, trf },
-        })
     }
 }
 
@@ -1135,25 +1074,21 @@ impl Primitive for CrossCap2 {
         })
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let trf = Point3::new(1.0, 1.0, 1.0);
+        let bln = Point3::new(-1.0, -1.0, -1.0);
+        AABB::new(bln, trf)
     }
 }
 
 //  Steiner  ---------
 #[derive(Clone)]
-pub struct Steiner {
+pub struct Steiner {}
-    bounding_box: AABB,
-}
 
 impl Steiner {
-    pub fn new() -> Rc<dyn Primitive> {
+    pub fn new() -> Arc<dyn Primitive> {
         // I need to find the bounding box for this shape
-        let trf = Point3::new(1.0, 1.0, 1.0);
+        Arc::new(Steiner {})
-        let bln = Point3::new(-1.0, -1.0, -1.0);
-        Rc::new(Steiner {
-            bounding_box: AABB { bln, trf },
-        })
     }
 }
 
@@ -1216,25 +1151,21 @@ impl Primitive for Steiner {
         })
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let trf = Point3::new(1.0, 1.0, 1.0);
+        let bln = Point3::new(-1.0, -1.0, -1.0);
+        AABB::new(bln, trf)
     }
 }
 
 //  Steiner 2 ---------
 #[derive(Clone)]
-pub struct Steiner2 {
+pub struct Steiner2 {}
-    bounding_box: AABB,
-}
 
 impl Steiner2 {
-    pub fn new() -> Rc<dyn Primitive> {
+    pub fn new() -> Arc<dyn Primitive> {
         // I need to find the bounding box for this shape
-        let trf = Point3::new(1.0, 1.0, 1.0);
+        Arc::new(Steiner2 {})
-        let bln = Point3::new(-1.0, -1.0, -1.0);
-        Rc::new(Steiner2 {
-            bounding_box: AABB { bln, trf },
-        })
     }
 }
 
@@ -1308,8 +1239,10 @@ impl Primitive for Steiner2 {
         })
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let trf = Point3::new(1.0, 1.0, 1.0);
+        let bln = Point3::new(-1.0, -1.0, -1.0);
+        AABB::new(bln, trf)
     }
 }
 
@@ -1317,18 +1250,12 @@ impl Primitive for Steiner2 {
 #[derive(Clone)]
 pub struct Roman {
     k: f64,
-    bounding_box: AABB,
 }
 
 impl Roman {
-    pub fn new(k: f64) -> Rc<dyn Primitive> {
+    pub fn new(k: f64) -> Arc<dyn Primitive> {
         // I need to find the bounding box for this shape
-        let trf = Point3::new(1.0, 1.0, 1.0);
+        Arc::new(Roman { k })
-        let bln = Point3::new(-1.0, -1.0, -1.0);
-        Rc::new(Roman {
-            k,
-            bounding_box: AABB { bln, trf },
-        })
     }
 }
 
@@ -1419,8 +1346,10 @@ impl Primitive for Roman {
         })
     }
 
-    fn intersect_bounding_box(&self, ray: &Ray) -> bool {
+    fn get_aabb(&self) -> AABB {
-        self.bounding_box.intersect_bounding_box(ray)
+        let trf = Point3::new(1.0, 1.0, 1.0);
+        let bln = Point3::new(-1.0, -1.0, -1.0);
+        AABB::new(bln, trf)
     }
 }
 

src/ray.rs

@@ -1,11 +1,7 @@
-use crate::{node::Node, scene::Scene, EPSILON};
+use crate::{bvh::BVH, light::Light, node::Node, scene::Scene, state::RaytracingOption, EPSILON};
-use nalgebra::{distance, Matrix4, Point3, Vector3};
+use nalgebra::{distance, Matrix3, Matrix4, Point3, Vector3};
 use rand;
 
-const MAX_DEPTH: u8 = 5;
-const DIFFUSE_RAYS: i8 = 5;
-const DIFFUSE_COEFFICIENT: f32 = 0.5;
-
 fn random_vec() -> Vector3<f64> {
     Vector3::new(rand::random(), rand::random(), rand::random())
 }
@@ -22,15 +18,17 @@ pub struct Intersection {
 }
 //Intersection point including point and normal
 impl Intersection {
-    pub fn transform(&self, trans: &Matrix4<f64>, inv_trans: &Matrix4<f64>) -> Intersection {
+    pub fn transform(&mut self, trans: &Matrix4<f64>, inv_trans: &Matrix4<f64>) -> Intersection {
-        let point = trans.transform_point(&self.point);
-        let normal = inv_trans.transpose().transform_vector(&self.normal);
         Intersection {
-            point,
+            point: trans.transform_point(&self.point),
-            normal,
+            normal: inv_trans.transpose().transform_vector(&self.normal),
             distance: self.distance,
         }
     }
+    pub fn transform_mut(&mut self, trans: &Matrix4<f64>, inv_transpose: &Matrix3<f64>) {
+        self.point = trans.transform_point(&self.point);
+        self.normal = inv_transpose * self.normal;
+    }
 }
 
 // Ray struct represents a ray in 3D space with a starting point 'a' and a direction 'b'
@@ -66,45 +64,45 @@ impl Ray {
             b: trans.transform_vector(&self.b),
         }
     }
+    //Transform mutably
+    pub fn transform_mut(&mut self, trans: &Matrix4<f64>) {
+        self.a = trans.transform_point(&self.a);
+        self.b = trans.transform_vector(&self.b);
+    }
     //This function will determine if the ray hits an object in the scene
-    pub fn hit_scene(&self, scene: &Scene) -> bool {
+    //This is not optimised as it does not include bounding boxes
+    pub fn hit_scene(ray: &Ray, scene: &Scene) -> bool {
         for (_, node) in &scene.nodes {
             if !node.active {
                 continue;
             }
             // Transform ray into local model cordinates
-            let ray = self.transform(&node.inv_model);
+            if node.intersect_ray(&ray).is_some() {
-            // Check bounding box intersection
-            if node.primitive.intersect_bounding_box(&ray) {
-                // Check primitive intersection
-                if node.primitive.intersect_ray(&ray).is_some() {
                 return true;
             }
         }
-        }
         false
     }
     //This function find the closest intersection point of a ray with an object in the scene
-    pub fn closest_intersect<'a>(&'a self, scene: &'a Scene) -> Option<(&Node, Intersection)> {
+    //Also not optimised, as it does not include bounding boxes
+    pub fn closest_intersect<'a>(
+        ray: &'a Ray,
+        scene: &'a Scene,
+    ) -> Option<(&'a Node, Intersection)> {
         let mut closest_distance = f64::MAX;
         let mut closest_intersect: Option<(&Node, Intersection)> = None;
+        let ray_a = ray.a;
         for (_, node) in &scene.nodes {
+            //position of ray in world coords
             if !node.active {
                 continue;
             }
-            // Transform ray into local model cordinates
+
-            let ray = self.transform(&node.inv_model);
+            if node.aabb.intersect_ray(&ray) {
-            // Check bounding box intersection
+                //Check node intersection
-            if node.primitive.intersect_bounding_box(&ray) {
+                if let Some(intersect) = node.intersect_ray(&ray) {
-                // Check primitive intersection
-                if let Some(intersect) = node.primitive.intersect_ray(&ray) {
-                    // Dont intersect with itself
-                    if intersect.distance < EPSILON {
-                        continue;
-                    }
                     // Check for closest distance by converting to world coords
-                    let intersect = intersect.transform(&node.model, &node.inv_model);
+                    let distance = distance(&ray_a, &intersect.point);
-                    let distance = distance(&ray.a, &intersect.point);
                     if distance < closest_distance {
                         closest_distance = distance;
                         closest_intersect = Some((node, intersect));
@@ -115,19 +113,42 @@ impl Ray {
         closest_intersect
     }
     // This function takes a scene and returns the color of the point where the ray intersects the scene
-    pub fn shade_ray(&self, scene: &Scene, depth: u8) -> Option<Vector3<f32>> {
+    pub fn shade_ray(
-        if depth == MAX_DEPTH {
+        &self,
+        scene: &Scene,
+        depth: u8,
+        options: &RaytracingOption,
+        sbvh: &Option<BVH>,
+    ) -> Option<Vector3<f32>> {
+        //If we have exceeded depth then return
+        if depth == options.ray_depth {
             return None;
         }
-        match self.closest_intersect(scene) {
+        match sbvh {
+            //We have a bvh so use bvh traversal
+            Some(bvh) => {
+                //Intersect the scene with the bvh
+                if let Some((node, intersect)) = bvh.traverse(self, 0) {
+                    return Some(Ray::phong_shade_point(
+                        &scene, &self, &node, &intersect, depth, options, sbvh,
+                    ));
+                }
+                return None;
+            }
+            //We dont have a bvh so use generic algorithm
+            None => {
+                //No BVH given so intersect normally
+                match Ray::closest_intersect(self, scene) {
                     Some((node, intersect)) => {
                         Some(Ray::phong_shade_point(
-                    &scene, &self, &node, &intersect, depth,
+                            &scene, &self, &node, &intersect, depth, options, sbvh,
                         )) // If there is an intersection, shade it
                     }
                     None => None, // If there is no intersection, return None
                 }
             }
+        }
+    }
 
     // Function to shade a point in the scene using Phong shading model
     pub fn phong_shade_point(
@@ -136,11 +157,12 @@ impl Ray {
         node: &Node,
         intersect: &Intersection,
         depth: u8,
+        options: &RaytracingOption,
+        bvh: &Option<BVH>,
     ) -> Vector3<f32> {
         let normal = &intersect.normal;
-        let point = intersect.point;
+        let point = &intersect.point;
         let incidence = &ray.b;
-
         let material = &node.material;
 
         // Compute the ambient light component and set it as base colour
@@ -161,67 +183,103 @@ impl Ray {
             let to_light = to_light.normalize();
 
             //Niave Shadows
-            let to_light_ray = Ray::new(point, to_light);
+            if options.shadows {
-            if to_light_ray.light_blocked(scene, node) {
+                let to_light_ray = Ray::new(*point, to_light);
+                if to_light_ray.light_blocked(scene, light, bvh) {
                     continue;
                 }
+            }
 
             let n_dot_l = normal.dot(&to_light).max(0.0) as f32;
 
             //Reflected component
             let mut reflect = Vector3::zeros();
+            if options.reflect {
                 let reflect_dir = incidence - 2.0 * incidence.dot(&normal) * normal;
-            let reflect_ray = Ray::new(point, reflect_dir);
+                let reflect_ray = Ray::new(*point, reflect_dir);
-            if let Some(col) = reflect_ray.shade_ray(scene, depth + 1) {
+                if let Some(col) = reflect_ray.shade_ray(scene, depth + 1, options, bvh) {
                     reflect += col.component_mul(&material.kr)
                 }
+            }
 
             //Diffuse component (Lambertian)
             let mut diffuse = Vector3::zeros();
+            if options.diffuse {
                 diffuse += material.kd * n_dot_l;
-            for _ in 0..DIFFUSE_RAYS {
+                for _ in 0..options.diffuse_rays {
                     let diffuse_dir = random_unit_vec();
                     let diffuse_ray = Ray::new(point.clone(), diffuse_dir + normal);
-                if let Some(col) = diffuse_ray.shade_ray(scene, depth + 1) {
+                    if let Some(col) = diffuse_ray.shade_ray(scene, depth + 1, options, bvh) {
-                    diffuse += col * DIFFUSE_COEFFICIENT;
+                        diffuse += col * options.diffuse_coefficient;
+                    }
                 }
             }
 
             //Specular component
             let mut specular = Vector3::zeros();
-            if n_dot_l < 0.0 {
+            if options.specular {
+                if n_dot_l > 0.0 {
                     let h = (to_light - incidence).normalize();
                     let n_dot_h = normal.dot(&h).max(0.0) as f32;
                     specular = material.ks * n_dot_h.powf(material.shininess);
                 }
+            }
 
             //Falloff
-            // let falloff = 1.0
+            let mut falloff = 1.0;
-            //     / (1.0
+            if options.falloff {
-            //         + light.falloff[0]
+                falloff = 1.0
-            //         + light.falloff[1] * light_distance
+                    / ((1.0 + light.falloff[0])
-            //         + light.falloff[2] * light_distance * light_distance);
+                        + light.falloff[1] * light_distance
+                        + light.falloff[2] * light_distance * light_distance);
+            }
 
-            let intensity = light.colour.component_mul(&(diffuse + reflect + specular));
+            let intensity = light.colour.component_mul(&(diffuse + reflect + specular)) * falloff;
             colour += &intensity;
         }
 
         colour
     }
 
-    pub fn light_blocked(&self, scene: &Scene, _node: &Node) -> bool {
+    pub fn light_blocked(&self, scene: &Scene, light: &Light, bvh: &Option<BVH>) -> bool {
+        let light_distance = distance(&self.a, &light.position);
+        match bvh {
+            Some(bvh) => {
+                //We have a bvh so use bvh traversal
                 for (_, node) in &scene.nodes {
                     if !node.active {
                         continue;
                     }
-            let ray = self.transform(&node.inv_model);
+                    match bvh.traverse(self, 0) {
-            if node.primitive.intersect_bounding_box(&ray) {
+                        Some((_, intersect)) => {
-                if node.primitive.intersect_ray(&ray).is_some() {
+                            if intersect.distance < light_distance {
                                 return true;
                             }
                         }
+                        None => continue,
                     }
-        false
+                }
+                return false;
+            }
+            None => {
+                for (_, node) in &scene.nodes {
+                    if !node.active {
+                        continue;
+                    }
+                    if node.aabb.intersect_ray(self) {
+                        match node.intersect_ray(self) {
+                            Some(intersect) => {
+                                if intersect.distance < light_distance {
+                                    return true;
+                                }
+                            }
+                            None => continue,
+                        }
+                    }
+                }
+            }
+        }
+        return false;
     }
     //Cast a set of rays
     pub fn cast_rays(

src/scene.rs

@@ -35,4 +35,10 @@ impl Scene {
     pub fn add_camera(&mut self, label: String, camera: Camera) {
         self.cameras.insert(label, camera);
     }
+    // Compute all matricies for nodes
+    pub fn compute(&mut self) {
+        for (_, node) in &mut self.nodes {
+            node.compute();
+        }
+    }
 }

src/state.rs

@@ -1,16 +1,19 @@
 //Use linear algebra module
 
+use crate::bvh::BVH;
 use crate::camera::Camera;
 use crate::ray::Ray;
 use crate::{gui::Gui, scene::Scene};
 use crate::{gui::GuiEvent, log_error};
 use std::path::Path;
+use std::thread;
 
 use nalgebra::Vector3;
 use rand::seq::SliceRandom;
 use rand::{random, thread_rng};
 
 use std::error::Error;
+use std::sync::Arc;
 
 use anyhow::Result;
 use pixels::{Pixels, SurfaceTexture};
@@ -21,16 +24,57 @@ use winit::window::{Window, WindowBuilder};
 
 const START_WIDTH: i32 = 1200;
 const START_HEIGHT: i32 = 700;
-const RAY_SAMPLES: i8 = 5;
-const RAY_RANDOMNESS: f64 = 100.0;
-const COLOUR_CLEAR: [u8; 4] = [0x22, 0x00, 0x11, 0x55];
-const PIXEL_CLEAR: [u8; 4] = [0x55, 0x00, 0x22, 0x55];
 
 pub const INIT_FILE: &str = "rhai/scene.rhai";
 pub const SAVE_FILE: &str = "img.png";
 
+#[derive(Clone)]
+pub struct RaytracingOption {
+    pub threads: u32,
+    pub ray_samples: u32,
+    pub ray_randomness: f64,
+    pub clear_color: [u8; 4],
+    pub pixel_clear: [u8; 4],
+    pub pixels_per_thread: u32,
+    pub buffer_proportion: f32,
+    pub buffer_fov: f64,
+    pub ray_depth: u8,
+    pub diffuse_rays: u8,
+    pub diffuse_coefficient: f32,
+    pub bvh_active: bool,
+    pub shadows: bool,
+    pub diffuse: bool,
+    pub reflect: bool,
+    pub specular: bool,
+    pub falloff: bool,
+}
+impl RaytracingOption {
+    pub fn default() -> RaytracingOption {
+        RaytracingOption {
+            threads: 12,
+            ray_samples: 1,
+            ray_randomness: 700.0,
+            clear_color: [0x22, 0x00, 0x11, 0x55],
+            pixel_clear: [0x11, 0x00, 0x22, 0x55],
+            pixels_per_thread: 100,
+            buffer_proportion: 1.0,
+            buffer_fov: 70.0,
+            ray_depth: 1,
+            diffuse_rays: 3,
+            diffuse_coefficient: 0.1,
+            bvh_active: false,
+            shadows: true,
+            diffuse: true,
+            reflect: true,
+            specular: true,
+            falloff: true,
+        }
+    }
+}
+
 pub struct State {
-    scene: Scene,
+    scene: Arc<Scene>,
+    bvh: Arc<Option<BVH>>,
     camera: Camera,
     window: Window,
 
@@ -40,52 +84,61 @@ pub struct State {
     pixels: Pixels,
     gui: Gui,
 
-    rays: Vec<Ray>,
+    rays: Arc<Vec<Ray>>,
     ray_queue: Vec<usize>,
+    raytracing_options: Arc<RaytracingOption>,
 }
 
 impl State {
     pub fn new(window: Window, pixels: Pixels, gui: Gui) -> Self {
-        let scene = Scene::empty();
+        let scene = Arc::new(Scene::empty());
         let window_size = window.inner_size();
+        let pixels = pixels;
         let camera = Camera::unit();
-        let rays = Vec::new();
+        let rays = Arc::new(Vec::new());
 
         Self {
             scene,
+            bvh: Arc::new(None),
             camera,
             window,
             buffer_width: window_size.width as u32,
             buffer_height: window_size.height as u32,
-            pixels: pixels,
+            pixels,
             gui,
             rays,
             ray_queue: Vec::new(),
+            raytracing_options: Arc::new(RaytracingOption::default()),
         }
     }
 
     fn update(&mut self) -> Result<(), Box<dyn Error>> {
         if let Some(event) = self.gui.event.take() {
             match event {
-                GuiEvent::BufferResize(proportion, fov) => {
+                GuiEvent::RaytracerOption(options) => {
-                    self.resize_buffer(proportion, fov as f64)?
+                    self.raytracing_options = Arc::new(options);
+                    match self.raytracing_options.bvh_active {
+                        true => self.bvh = Arc::new(Some(BVH::build(&self.scene.nodes))),
+                        false => self.bvh = Arc::new(None),
                     }
-                GuiEvent::CameraUpdate(camera, fovy) => {
+                    self.resize_buffer()?
-                    self.rays = Ray::cast_rays(
+                }
+                GuiEvent::CameraUpdate(camera) => {
+                    self.rays = Arc::new(Ray::cast_rays(
                         &camera.eye,
                         &camera.target,
                         &camera.up,
-                        fovy as f64,
+                        self.raytracing_options.buffer_fov,
                         self.buffer_width,
                         self.buffer_height,
-                    );
+                    ));
                     self.camera = camera;
-                    self.clear()?;
+                    self.clear_buffer()?;
                     self.reset_queue();
                 }
                 GuiEvent::SceneLoad(scene) => {
-                    self.scene = scene;
+                    self.scene = Arc::new(scene);
-                    self.clear()?;
+                    self.clear_buffer()?;
                     self.reset_queue();
                 }
                 GuiEvent::SaveImage(filename) => {
@@ -103,30 +156,32 @@ impl State {
         Ok(())
     }
 
-    fn resize_buffer(&mut self, proportion: f32, fovy: f64) -> Result<(), Box<dyn Error>> {
+    fn resize_buffer(&mut self) -> Result<(), Box<dyn Error>> {
         // Calculate new buffer dimensions based on proportion
         let size = self.window.inner_size();
+        let proportion = &self.raytracing_options.buffer_proportion;
+        let fovy = self.raytracing_options.buffer_fov;
         self.buffer_width = (size.width as f32 * proportion) as u32;
         self.buffer_height = (size.height as f32 * proportion) as u32;
 
         // Clear the buffer and reset the ray queue
-        self.clear()?;
+        self.clear_buffer()?;
         self.reset_queue();
 
         // Recalculate rays with new buffer dimensions
-        self.rays = Ray::cast_rays(
+        self.rays = Arc::new(Ray::cast_rays(
             &self.camera.eye,
             &self.camera.target,
             &self.camera.up,
             fovy,
             self.buffer_width,
             self.buffer_height,
-        );
+        ));
 
         // Resize buffer and surface
-        let pixels = &mut self.pixels;
+        self.pixels.resize_surface(size.width, size.height)?;
-        pixels.resize_surface(size.width, size.height)?;
+        self.pixels
-        pixels.resize_buffer(self.buffer_width, self.buffer_height)?;
+            .resize_buffer(self.buffer_width, self.buffer_height)?;
 
         Ok(())
     }
@@ -148,48 +203,102 @@ impl State {
 
     fn draw(&mut self) -> Result<(), Box<dyn Error>> {
         //Draw ray_num in a block
-        let frame = self.pixels.frame_mut();
+        let randomness = self.raytracing_options.ray_randomness;
-        for _ in 0..self.gui.ray_num {
+        let samples = self.raytracing_options.ray_samples;
-            //Get random index from queue
+        let samples_f32 = samples as f32;
-            let index = match self.ray_queue.pop() {
+
-                Some(index) => index,
+        let num_threads = self.raytracing_options.threads;
+        let pixels_per_thread = self.raytracing_options.pixels_per_thread;
+
+        let mut handles = vec![];
+
+        for _ in 0..num_threads {
+            //Get necessary variables to render
+            let rays = self.rays.clone();
+            let scene = self.scene.clone();
+            let options = self.raytracing_options.clone();
+            let bvh = self.bvh.clone();
+
+            //Get the workload for a thread
+            let mut load = vec![];
+            for _ in 0..pixels_per_thread {
+                match self.ray_queue.pop() {
+                    Some(index) => load.push(index),
                     None => break,
-            };
+                }
-            //Shade colour for selected ray
+            }
-            let mut colour = Vector3::zeros();
+            //The finished queue of the thread
-            for _ in 0..RAY_SAMPLES {
+            let mut finished = vec![];
-                let ray = &self.rays[index];
+
+            //Create a new thread for these pixels
+            let handle = thread::spawn({
+                move || {
+                    for index in &load {
+                        //Shade colour for selected index
+                        let mut colour: Vector3<f32> = Vector3::zeros();
+                        let ray = &rays[*index];
+                        for _ in 0..samples {
+                            //Generate a ray in a random direction
                             let point = ray.a;
                             let dir = ray.b;
-                let rx = (random::<f64>() - 0.5) / RAY_RANDOMNESS;
+                            let rx = (random::<f64>() - 0.5) / randomness;
-                let ry = (random::<f64>() - 0.5) / RAY_RANDOMNESS;
+                            let ry = (random::<f64>() - 0.5) / randomness;
-                let rz = (random::<f64>() - 0.5) / RAY_RANDOMNESS;
+                            let rz = (random::<f64>() - 0.5) / randomness;
                             let nx = dir.x + rx;
                             let ny = dir.y + ry;
                             let nz = dir.z + rz;
 
                             let rand_ray = Ray::new(point, Vector3::new(nx, ny, nz));
 
-                if let Some(ray_colour) = rand_ray.shade_ray(&self.scene, 0) {
+                            if let Some(ray_colour) = rand_ray.shade_ray(&scene, 0, &options, &bvh)
+                            {
                                 colour += ray_colour;
-                };
                             }
-            colour = (colour / RAY_SAMPLES as f32) * 255.0;
+                        }
+                        colour = (colour / samples_f32) * 255.0;
                         let rgba = [colour.x as u8, colour.y as u8, colour.z as u8, 0xff];
+                        finished.push(rgba);
+                    }
+                    return (load, finished);
+                }
+            });
+            handles.push(handle);
+        }
+
+        let mut all_results = vec![];
+
+        for handle in handles.drain(..) {
+            let (load, finished) = handle
+                .join()
+                .map_err(|e| format!("Thread panicked: {:?}", e))?;
+            let thread_results: Vec<_> = load.into_iter().zip(finished.into_iter()).collect();
+            all_results.extend(thread_results);
+        }
+
+        //Now we have two vectors will all the indicies and rgba values, we can upload them to the bufer
+
+        let frame = self.pixels.frame_mut();
+        for result in all_results {
+            let index = result.0;
+            let rgba = result.1;
             frame[index * 4..(index + 1) * 4].copy_from_slice(&rgba);
         }
         Ok(())
     }
 
-    fn clear(&mut self) -> Result<(), Box<dyn Error>> {
+    fn clear_buffer(&mut self) -> Result<(), Box<dyn Error>> {
         let frame = self.pixels.frame_mut();
         for pixel in frame.chunks_exact_mut(4) {
-            pixel.copy_from_slice(&COLOUR_CLEAR);
+            pixel.copy_from_slice(&self.raytracing_options.pixel_clear);
         }
         Ok(())
     }
 
     fn reset_queue(&mut self) {
+        match self.raytracing_options.bvh_active {
+            true => self.bvh = Arc::new(Some(BVH::build(&self.scene.nodes))),
+            false => self.bvh = Arc::new(None),
+        }
         let size = self.buffer_height as usize * self.buffer_width as usize;
         let mut ray_queue: Vec<usize> = (0..size).collect();
         ray_queue.shuffle(&mut thread_rng());
@@ -200,14 +309,12 @@ impl State {
         // Update state
         self.update()?;
         // Draw rays if we have remaining rays in queue
-        if !self.ray_queue.is_empty() {
         match self.draw() {
             Err(e) => {
                 println!("ERROR: {}", e);
             }
             _ => {}
         }
-        }
         // Render Gui
         self.gui
             .prepare(&self.window)
@@ -234,7 +341,7 @@ pub fn run() -> Result<(), Box<dyn Error>> {
     let gui = Gui::new(&window, &pixels);
 
     let mut state = State::new(window, pixels, gui);
-    state.resize_buffer(1.0, 90.0)?;
+    state.resize_buffer()?;
 
     event_loop.run(move |event, _, control_flow| {
         state.gui.handle_event(&state.window, &event);