initial commit

src/consts.rs

@@ -0,0 +1,43 @@
+use crate::structs::{Vec3, Y, X, Z};
+
+pub const COLOR_CHANNELS: usize = 4;
+
+pub const EPSILON: f64 = 1. / 1024.;
+pub const LIGHT_EPSILON: f64 = 1. / 512.;
+pub const MAX_DIST: f64 = 16.;
+pub const MAX_STEPS: u32 = 1024;
+
+pub const DIST_FIX_CORRECTION: f64 = 0.5;
+pub const DIST_CORRECTION: f64 = 3.;
+pub const DIST_POWER: f64 = 2.;
+
+pub const ANGLE_FIX_CORRECTION: f64 = 0.;
+pub const ANGLE_CORRECTION: f64 = 2.;
+pub const ANGLE_POWER: f64 = 2.;
+
+//pub const IMG_WIDTH: usize = 480;
+//pub const IMG_WIDTH: usize = 1280;
+pub const IMG_WIDTH: usize = 1080;
+//pub const IMG_HEIGHT: usize = 480;
+//pub const IMG_HEIGHT: usize = 720;
+pub const IMG_HEIGHT: usize = 1080;
+pub const IMG_DIM: usize = if IMG_HEIGHT > IMG_WIDTH { IMG_HEIGHT } else { IMG_WIDTH };
+pub const IMG_SIZE: usize = IMG_WIDTH * IMG_HEIGHT;
+pub const IMG_BYTE_SIZE: usize = IMG_SIZE * 3;
+
+//pub const SUPERSAMPLING: usize = 1;
+pub const SUPERSAMPLING: usize = 2;
+//pub const RAYS_PER_PIXEL: usize = 1;
+//pub const RAYS_PER_PIXEL: usize = 50;
+pub const RAYS_PER_PIXEL: usize = 500;
+//pub const MAX_BOUNCES: u32 = 1;
+//pub const MAX_BOUNCES: u32 = 4;
+//pub const MAX_BOUNCES: u32 = 8;
+pub const MAX_BOUNCES: u32 = 10;
+
+pub const THREAD_COUNT: usize = 12;
+pub const SLICES_PER_THREAD: usize = 4;
+
+pub const UP: Vec3 = Y;
+pub const RIGHT: Vec3 = X;
+pub const FORWARD: Vec3 = Z;

src/light.rs

@@ -0,0 +1,16 @@
+use crate::structs::{ColorVec, Vec3};
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Light {
+    pos: Vec3,
+    color: ColorVec,
+}
+
+impl Light {
+    pub fn new(pos: Vec3, color: ColorVec) -> Light {
+        Light { pos, color }
+    }
+
+    pub fn pos(&self) -> Vec3 { self.pos }
+    pub fn color(&self) -> ColorVec { self.color }
+}

src/main.rs

@@ -0,0 +1,84 @@
+mod material;
+mod light;
+mod object;
+mod structs;
+mod consts;
+
+use object::*;
+use crate::consts::*;
+use crate::structs::*;
+use crate::material::*;
+use image::{ColorType, ImageFormat};
+use image::codecs::png::FilterType::Sub;
+
+fn default_cam() -> Cam {
+    Cam::new_pointing(Y*3. - X*5., O, 0.5)
+}
+
+fn default_scene() -> Scene {
+    let s0 = WithMaterial::new(
+        Sphere::new_xyz(0., 0., 0., 1.),
+        RED
+    );
+    let s1 = WithMaterial::new(
+        Sphere::new_xyz(0., 1., 1., 0.5),
+        MIRROR
+    );
+    let spheres = Union::new(s0, s1);
+    let floor = WithMaterial::new(
+        Plane::new_xyz(0., 1., 0., 2.),
+        WHITE
+    );
+    let backwall = WithMaterial::new(
+        Plane::new_xyz(-1., 0., 0., 2.),
+        BLUE
+    );
+    let sidewalls = WithMaterial::new(
+        Union::new(
+            Plane::new_xyz(0., 0., -1., 2.5),
+            Plane::new_xyz(0., 0., 1., 2.5)
+        ),
+        GREEN
+    );
+    let walls = Union::new(
+        Union::new(floor, backwall),
+        sidewalls
+    );
+    let scene = Union::new(spheres, walls);
+    let light = WithMaterial::new(
+        Plane::new_xyz(0., -1., 0., 8.),
+        LIGHTSOURCE
+    );
+    let scene = Union::new(scene, light);
+    Scene::new(scene)
+}
+
+fn main() {
+    // get scene and camera
+    let scene = default_scene();
+    let cam = default_cam();
+
+    // get stats on the scene we're about to render
+    let total_rpp = (RAYS_PER_PIXEL * SUPERSAMPLING * SUPERSAMPLING) as u64;
+    let lights = scene.get_lights().len() as u32;
+    let nodes = scene.node_count();
+    let min_rays = IMG_SIZE as u64 * RAYS_PER_PIXEL as u64 * SUPERSAMPLING as u64;
+    let max_rays = min_rays * MAX_BOUNCES as u64 * (1+lights) as u64;
+    println!("Image size: {}x{} ({} pixels)", IMG_WIDTH, IMG_HEIGHT, IMG_SIZE);
+    println!("Rpp: {}, {}x supersampling ({} total rays per image pixel)", RAYS_PER_PIXEL, SUPERSAMPLING, total_rpp);
+    println!("Max bounces: {}, dist: {}, steps: {}", MAX_BOUNCES, MAX_DIST, MAX_STEPS);
+    println!("Lights: {}, Nodes: {}", lights, nodes);
+    println!("Threads: {}, {} slices per thread", THREAD_COUNT, SLICES_PER_THREAD);
+    println!("Total rays: min {} max {}", min_rays, max_rays);
+
+    // generate image and save it
+    let data = cam.render(&scene);
+    image::save_buffer_with_format(
+        "a.png",
+        &data,
+        IMG_WIDTH as u32,
+        IMG_HEIGHT as u32,
+        ColorType::Rgb8,
+        ImageFormat::Png
+    ).unwrap();
+}

src/material.rs

@@ -0,0 +1,88 @@
+use crate::structs::{ColorVec, ColorMat};
+use crate::consts::COLOR_CHANNELS;
+use crate::material::SurfaceType::{Diffuse, Reflective, Stop};
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub enum SurfaceType {
+    Diffuse,
+    Reflective,
+    Stop
+}
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Material {
+    emission: ColorVec,
+    reflection: ColorMat,
+    surface: SurfaceType
+}
+
+impl Material {
+    pub const fn new(emission: ColorVec, reflection: ColorMat, surface: SurfaceType) -> Material {
+        Material { emission, reflection, surface }
+    }
+
+    pub const fn new_from_diagonal(emission: ColorVec, reflection: ColorVec, surface: SurfaceType) -> Material {
+        Material { emission, reflection: ColorMat::new_from_diagonal(reflection), surface }
+    }
+
+    pub const fn emission(&self) -> ColorVec { self.emission }
+    pub const fn reflection(&self) -> ColorMat { self.reflection }
+    pub const fn surface(&self) -> SurfaceType { self.surface }
+}
+
+const COLOR_ZERO: ColorVec = ColorVec::new_zero();
+const COLOR_ONE: ColorVec = ColorVec::new([1.; COLOR_CHANNELS]);
+
+pub const WHITE: Material = Material::new_from_diagonal(
+    COLOR_ZERO,
+    COLOR_ONE,
+    Diffuse
+);
+pub const RED: Material = Material::new_from_diagonal(
+    COLOR_ZERO,
+    ColorVec::new([0.75, 0.25, 0.25, 0.]),
+    Diffuse
+);
+pub const GREEN: Material = Material::new_from_diagonal(
+    COLOR_ZERO,
+    ColorVec::new([0.25, 0.75, 0.25, 0.]),
+    Diffuse
+);
+pub const BLUE: Material = Material::new_from_diagonal(
+    COLOR_ZERO,
+    ColorVec::new([0.25, 0.25, 0.75, 0.]),
+    Diffuse
+);
+pub const LIGHTSOURCE: Material = Material::new_from_diagonal(
+    ColorVec::new([1., 1., 1., 0.]),
+    COLOR_ONE,
+    Diffuse
+);
+
+pub const MIRROR: Material = Material::new_from_diagonal(
+    COLOR_ZERO,
+    ColorVec::new([0.9, 0.9, 0.9, 0.9]),
+    Reflective
+);
+
+pub const STRONG_LIGHTSOURCE: Material = Material::new_from_diagonal(
+    ColorVec::new([5., 5., 5., 0.]),
+    COLOR_ZERO,
+    Stop
+);
+
+pub const UV_LIGHTSOURCE: Material = Material::new_from_diagonal(
+    ColorVec::new_one(3),
+    ColorVec::new([0.25, 0.25, 0.25, 1.]),
+    Diffuse
+);
+pub const FLUORESCENT: Material = Material::new(
+    ColorVec::new([0.25; COLOR_CHANNELS]),
+    ColorMat::new([
+        [1., 0., 0., 0.75],
+        [0., 1., 0., 0.25],
+        [0., 0., 1., 0.],
+        [0., 0., 0., 0.]
+    ]),
+    Diffuse
+);

src/object/cuboid.rs

@@ -0,0 +1,65 @@
+use crate::structs::Vec3;
+use crate::object::Obj;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Cuboid {
+    center: Vec3,
+    radius: Vec3
+}
+
+impl Cuboid {
+    pub fn new(center: Vec3, radius: Vec3) -> Cuboid {
+        Cuboid { center, radius }
+    }
+    pub fn new_xyz(cx: f64, cy: f64, cz: f64, rx: f64, ry: f64, rz: f64) -> Cuboid {
+        let center = Vec3::new(cx, cy, cz);
+        let radius = Vec3::new(rx, ry, rz);
+        Cuboid { center, radius }
+    }
+}
+
+impl Obj for Cuboid {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        let dx = point.x() - self.center.x();
+        let dy = point.y() - self.center.y();
+        let dz = point.z() - self.center.z();
+
+        if dx.abs()<self.radius.x() && dy.abs()<self.radius.y() && dz.abs()<self.radius.z() {
+            -f64::min(
+                f64::min(
+                    f64::min(
+                        f64::abs(dx-self.radius.x()),
+                        f64::abs(dx+self.radius.x())
+                    ),
+                    f64::abs(dy-self.radius.y())
+                ),
+                f64::min(
+                    f64::min(
+                        f64::abs(dy+self.radius.y()),
+                        f64::abs(dz-self.radius.z())
+                    ),
+                    f64::abs(dz+self.radius.z())
+                )
+            )
+        } else {
+            let dx0 = f64::max(0., dx.abs()-self.radius.x());
+            let dy0 = f64::max(0., dy.abs()-self.radius.y());
+            let dz0 = f64::max(0., dz.abs()-self.radius.z());
+
+            f64::sqrt(dx0*dx0 + dy0*dy0 + dz0*dz0)
+        }
+    }
+
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        let dx = (point.x()-self.center.y())/self.radius.x();
+        let dy = (point.y()-self.center.y())/self.radius.y();
+        let dz = (point.z()-self.center.z())/self.radius.z();
+        let (adx, ady, adz) = (dx.abs(), dy.abs(), dz.abs());
+
+        Vec3::new(
+            if adx>=ady && adx>adz { if dx>0. { 1. } else { -1. }} else { 0. },
+            if ady>=adz && ady>adx { if dy>0. { 1. } else { -1. }} else { 0. },
+            if adz>=adx && adz>ady { if dz>0. { 1. } else { -1. }} else { 0. },
+        )
+    }
+}

src/object/cylinder.rs

@@ -0,0 +1,23 @@
+use crate::structs::Vec3;
+use crate::object::Obj;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Cylinder {
+    center: Vec3,
+    radius: f64
+}
+
+impl Cylinder {
+    pub const fn new(center: Vec3, radius: f64) -> Cylinder {
+        Cylinder { center, radius }
+    }
+}
+
+impl Obj for Cylinder {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        Vec3::new(self.center.x(), 0., self.center.z()).distance_to(point) - self.radius
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        point - Vec3::new(self.center.x(), 0.,self.center.z())
+    }
+}

src/object/intersection.rs

@@ -0,0 +1,45 @@
+use crate::object::Obj;
+use crate::light::Light;
+use crate::structs::Vec3;
+use std::vec;
+use crate::material::Material;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Intersection<A: Obj, B: Obj> {
+    a: A,
+    b: B
+}
+
+impl<A: Obj, B: Obj> Intersection<A, B> {
+    pub fn new(a: A, b: B) -> Intersection<A, B> {
+        Intersection { a, b }
+    }
+}
+
+impl<A: Obj, B: Obj> Obj for Intersection<A, B> {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        f64::max(self.a.distance_to(point), self.b.distance_to(point))
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        if self.a.distance_to(point) > self.b.distance_to(point) {
+            self.a.normal_at(point)
+        } else {
+            self.b.normal_at(point)
+        }
+    }
+    fn material_at(&self, point: Vec3) -> Material {
+        if self.a.distance_to(point) > self.b.distance_to(point) {
+            self.a.material_at(point)
+        } else {
+            self.b.material_at(point)
+        }
+    }
+    fn get_lights(&self) -> vec::Vec<Light> {
+        let mut l = self.a.get_lights();
+        l.extend(self.b.get_lights());
+        l
+    }
+    fn node_count(&self) -> u32 {
+        self.a.node_count() + self.b.node_count() + 1
+    }
+}

src/object/mod.rs

@@ -0,0 +1,71 @@
+use crate::consts::{EPSILON, MAX_DIST};
+use crate::structs::{Vec3, O, X, Y, Z};
+use crate::material::Material;
+use crate::material;
+use crate::light::Light;
+use std::vec::Vec;
+
+mod sphere;
+mod plane;
+mod union;
+mod intersection;
+mod cuboid;
+mod cylinder;
+mod torus;
+mod waves;
+mod with_material;
+mod with_lights;
+mod transform;
+mod scene;
+
+pub trait Obj: Send + Sync {
+    fn distance_to(&self, _point: Vec3) -> f64 {
+        MAX_DIST
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        let v = self.distance_to(point);
+        let x = self.distance_to(point + X*EPSILON) - v;
+        let y = self.distance_to(point + Y*EPSILON) - v;
+        let z = self.distance_to(point + Z*EPSILON) - v;
+
+        Vec3::new(x, y, z).unit()
+    }
+    fn material_at(&self, _point: Vec3) -> Material {
+        material::WHITE
+    }
+    fn get_lights(&self) -> Vec<Light> {
+        Vec::new()
+    }
+    fn node_count(&self) -> u32 { 1 }
+}
+
+impl<T: Obj> Obj for &T {
+    fn distance_to(&self, point: Vec3) -> f64 { (*self).distance_to(point) }
+    fn normal_at(&self, point: Vec3) -> Vec3 { (*self).normal_at(point) }
+    fn material_at(&self, point: Vec3) -> Material { (*self).material_at(point) }
+    fn get_lights(&self) -> Vec<Light> { (*self).get_lights() }
+    fn node_count(&self) -> u32 { (*self).node_count() }
+}
+
+pub trait ObjClone: Obj {
+    fn clone_obj(&self) -> Box<dyn ObjClone>;
+}
+
+impl<T: 'static + Obj + Clone> ObjClone for T {
+    fn clone_obj(&self) -> Box<dyn ObjClone> {
+        Box::new(self.clone())
+    }
+}
+
+pub use sphere::Sphere;
+pub use plane::Plane;
+pub use union::Union;
+pub use intersection::Intersection;
+pub use cuboid::Cuboid;
+pub use cylinder::Cylinder;
+pub use torus::Torus;
+pub use waves::Waves;
+pub use with_material::{WithMaterial, WithDynamicMaterial};
+pub use with_lights::{WithLights, WithLight};
+pub use transform::AffineTransform;
+pub use scene::Scene;

src/object/plane.rs

@@ -0,0 +1,31 @@
+use crate::structs::Vec3;
+use crate::object::Obj;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Plane {
+    normal: Vec3,
+    offset: f64
+}
+
+impl Plane {
+    pub fn new(normal: Vec3, offset: f64) -> Plane {
+        let l = normal.magnitude();
+        Plane { normal: normal/l, offset: offset/l }
+    }
+    pub fn new_xyz(x: f64, y: f64, z: f64, offset: f64) -> Plane {
+        let normal = Vec3::new(x, y, z);
+        Plane::new(normal, offset)
+    }
+    pub const unsafe fn new_raw(normal: Vec3, offset: f64) -> Plane {
+        Plane { normal, offset }
+    }
+}
+
+impl Obj for Plane {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        point*self.normal + self.offset
+    }
+    fn normal_at(&self, _point: Vec3) -> Vec3 {
+        self.normal
+    }
+}

src/object/scene.rs

@@ -0,0 +1,31 @@
+use crate::object::{Obj, ObjClone};
+use crate::structs::Vec3;
+use crate::material::Material;
+use crate::light::Light;
+
+pub struct Scene {
+    scene: Box<dyn ObjClone>
+}
+
+impl Scene {
+    pub fn new<T: 'static + ObjClone>(scene: T) -> Scene {
+        Scene { scene: Box::new(scene) }
+    }
+    pub fn new_from_box(scene: Box<dyn ObjClone>) -> Scene {
+        Scene { scene }
+    }
+}
+
+impl Obj for Scene {
+    fn distance_to(&self, point: Vec3) -> f64 { self.scene.distance_to(point) }
+    fn normal_at(&self, point: Vec3) -> Vec3 { self.scene.normal_at(point) }
+    fn material_at(&self, point: Vec3) -> Material { self.scene.material_at(point) }
+    fn get_lights(&self) -> Vec<Light> { self.scene.get_lights() }
+    fn node_count(&self) -> u32 { self.scene.node_count() + 1 }
+}
+
+impl Clone for Scene {
+    fn clone(&self) -> Scene {
+        Scene { scene: self.scene.clone_obj() }
+    }
+}

src/object/sphere.rs

@@ -0,0 +1,27 @@
+use crate::structs::Vec3;
+use crate::object::Obj;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Sphere {
+    pos: Vec3,
+    radius: f64
+}
+
+impl Sphere {
+    pub const fn new(pos: Vec3, radius: f64) -> Sphere {
+        Sphere { pos, radius }
+    }
+    pub const fn new_xyz(x: f64, y: f64, z: f64, radius: f64) -> Sphere {
+        let pos = Vec3::new(x, y, z);
+        Sphere { pos, radius }
+    }
+}
+
+impl Obj for Sphere {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        self.pos.distance_to(point) - self.radius
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        (point - self.pos).unit()
+    }
+}

src/object/torus.rs

@@ -0,0 +1,34 @@
+use crate::structs::Vec3;
+use crate::object::Obj;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Torus {
+    center: Vec3,
+    radius: f64,
+    thickness: f64
+}
+
+impl Torus {
+    pub const fn new(center: Vec3, radius: f64, thickness: f64) -> Torus {
+        Torus { center, radius, thickness }
+    }
+    pub const fn new_xyz(x: f64, y: f64, z: f64, radius: f64, thickness: f64) -> Torus {
+        Torus::new(Vec3::new(x, y, z), radius, thickness)
+    }
+}
+
+impl Obj for Torus {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        let (dx, dy) = (point.x()-self.center.x(), point.y()-self.center.y());
+        let dh = f64::abs(f64::sqrt(dx*dx + dy*dy) - self.radius);
+        let dz = point.z() - self.center.z();
+        f64::sqrt(dh*dh + dz*dz) - self.thickness
+    }
+
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        let centered = point - self.center;
+        let centered_no_z= Vec3::new(centered.x(), centered.y(), 0.);
+        let closest = centered_no_z.unit() * self.radius;
+        (centered - closest).unit()
+    }
+}

src/object/transform.rs

@@ -0,0 +1,83 @@
+use crate::object::Obj;
+use crate::structs::{Mat3, Vec3, O, I3};
+use crate::material::Material;
+use crate::light::Light;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct AffineTransform<T: Obj + Clone> {
+    obj: T,
+    transform: Mat3,
+    transform_inv: Mat3,
+    translate: Vec3
+}
+
+impl<T: Obj + Clone> AffineTransform<T> {
+    pub fn new(obj: T, transform: Mat3, translate: Vec3) -> AffineTransform<T> {
+        AffineTransform { obj, transform, transform_inv: transform.invert(), translate }
+    }
+    pub fn new_linear(obj: T, transform: Mat3) -> AffineTransform<T> {
+        AffineTransform { obj, transform, transform_inv: transform.invert(), translate: O }
+    }
+    pub fn new_translate(obj: T, translate: Vec3) -> AffineTransform<T> {
+        AffineTransform { obj, transform: I3, transform_inv: I3, translate }
+    }
+
+    fn apply_rev(&self, point: Vec3) -> Vec3 {
+        self.transform_inv*point - self.translate
+    }
+    fn apply_fwd(&self, point: Vec3) -> Vec3 {
+        self.transform*point + self.translate
+    }
+}
+
+pub const SWAP_XY: Mat3 = Mat3::new(
+    0., 1., 0.,
+    1., 0., 0.,
+    0., 0., 1.,
+);
+pub const SWAP_XZ: Mat3 = Mat3::new(
+    0., 0., 1.,
+    0., 1., 0.,
+    1., 0., 0.,
+);
+pub const SWAP_YZ: Mat3 = Mat3::new(
+    1., 0., 0.,
+    0., 0., 1.,
+    0., 1., 0.,
+);
+
+pub const fn scale_xyz(x: f64, y: f64, z: f64) -> Mat3 {
+    Mat3::new(
+        x , 0., 0.,
+        0., y , 0.,
+        0., 0., z ,
+    )
+}
+pub const fn scale(k: f64) -> Mat3 { scale_xyz(k, k, k) }
+pub const fn scale_x(k: f64) -> Mat3 { scale_xyz(k, 1., 1.) }
+pub const fn scale_y(k: f64) -> Mat3 { scale_xyz(1., k, 1.) }
+pub const fn scale_z(k: f64) -> Mat3 { scale_xyz(1., 1., k) }
+
+impl<T: Obj + Clone> Obj for AffineTransform<T> {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        self.obj.distance_to(self.apply_rev(point))
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        self.apply_fwd(self.obj.normal_at(self.apply_rev(point)))
+    }
+    fn material_at(&self, point: Vec3) -> Material {
+        self.obj.material_at(self.apply_rev(point))
+    }
+    fn get_lights(&self) -> Vec<Light> {
+        self.obj
+            .get_lights()
+            .into_iter()
+            .map(|light| {
+                Light::new(self.apply_fwd(light.pos()), light.color())
+            })
+            .collect()
+    }
+    fn node_count(&self) -> u32 {
+        self.obj.node_count() + 1
+    }
+}

src/object/union.rs

@@ -0,0 +1,45 @@
+use crate::object::Obj;
+use crate::light::Light;
+use crate::structs::Vec3;
+use std::vec;
+use crate::material::Material;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Union<A: Obj, B: Obj> {
+    a: A,
+    b: B
+}
+
+impl<A: Obj, B: Obj> Union<A, B> {
+    pub fn new(a: A, b: B) -> Union<A, B> {
+        Union { a, b }
+    }
+}
+
+impl<A: Obj, B: Obj> Obj for Union<A, B> {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        f64::min(self.a.distance_to(point), self.b.distance_to(point))
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        if self.a.distance_to(point) < self.b.distance_to(point) {
+            self.a.normal_at(point)
+        } else {
+            self.b.normal_at(point)
+        }
+    }
+    fn material_at(&self, point: Vec3) -> Material {
+        if self.a.distance_to(point) < self.b.distance_to(point) {
+            self.a.material_at(point)
+        } else {
+            self.b.material_at(point)
+        }
+    }
+    fn get_lights(&self) -> vec::Vec<Light> {
+        let mut l = self.a.get_lights();
+        l.extend(self.b.get_lights());
+        l
+    }
+    fn node_count(&self) -> u32 {
+        self.a.node_count() + self.b.node_count() + 1
+    }
+}

src/object/waves.rs

@@ -0,0 +1,22 @@
+use crate::object::Obj;
+use crate::structs::Vec3;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Waves {}
+
+impl Waves {
+    pub const fn new() -> Waves {
+        Waves {}
+    }
+}
+
+impl Obj for Waves {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        let dist = point.x().sin() + point.y().sin() + point.z();
+        dist / f64::sqrt(3.)
+    }
+
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        Vec3::new(-point.x().cos(), -point.y().cos(), 1.).unit()
+    }
+}

src/object/with_lights.rs

@@ -0,0 +1,44 @@
+use crate::light::Light;
+use crate::object::Obj;
+use crate::structs::Vec3;
+use crate::material::Material;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct WithLights<T: Obj, const N: usize> {
+    obj: T,
+    lights: [Light; N]
+}
+
+impl<T: Obj, const N: usize> WithLights<T, N> {
+    pub fn new(obj: T, lights: [Light; N]) -> WithLights<T, N> {
+        WithLights { obj, lights }
+    }
+}
+
+impl<T: Obj, const N: usize> Obj for WithLights<T, N> {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        self.obj.distance_to(point)
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        self.obj.normal_at(point)
+    }
+    fn material_at(&self, point: Vec3) -> Material {
+        self.obj.material_at(point)
+    }
+    fn get_lights(&self) -> Vec<Light> {
+        let mut l = self.obj.get_lights();
+        l.extend(&self.lights);
+        l
+    }
+    fn node_count(&self) -> u32 {
+        self.obj.node_count() + 1
+    }
+}
+
+pub type WithLight<T> = WithLights<T, 1>;
+
+impl<T: Obj> WithLight<T> {
+    pub fn new_one(obj: T, light: Light) -> WithLight<T> {
+        WithLight { obj, lights: [light; 1] }
+    }
+}

src/object/with_material.rs

@@ -0,0 +1,64 @@
+use crate::object::Obj;
+use crate::material::Material;
+use crate::structs::Vec3;
+use crate::light::Light;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct WithMaterial<T: Obj> {
+    obj: T,
+    material: Material
+}
+
+impl<T: Obj> WithMaterial<T> {
+    pub fn new(obj: T, material: Material) -> WithMaterial<T> {
+        WithMaterial { obj, material }
+    }
+}
+
+impl<T: Obj> Obj for WithMaterial<T> {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        self.obj.distance_to(point)
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        self.obj.normal_at(point)
+    }
+    fn material_at(&self, point: Vec3) -> Material {
+        self.material
+    }
+    fn get_lights(&self) -> Vec<Light> {
+        self.obj.get_lights()
+    }
+    fn node_count(&self) -> u32 {
+        self.obj.node_count() + 1
+    }
+}
+
+#[derive(Debug, Clone, PartialEq)]
+pub struct WithDynamicMaterial<T: Obj, F: Clone + Send + Sync + Fn(Vec3) -> Material> {
+    obj: T,
+    fun: F
+}
+
+impl<T: Obj, F: Clone + Send + Sync + Fn(Vec3) -> Material> WithDynamicMaterial<T, F> {
+    pub fn new(obj: T, fun: F) -> WithDynamicMaterial<T, F> {
+        WithDynamicMaterial { obj, fun }
+    }
+}
+
+impl<T: Obj, F: Clone + Send + Sync + Fn(Vec3) -> Material> Obj for WithDynamicMaterial<T, F> {
+    fn distance_to(&self, point: Vec3) -> f64 {
+        self.obj.distance_to(point)
+    }
+    fn normal_at(&self, point: Vec3) -> Vec3 {
+        self.obj.normal_at(point)
+    }
+    fn material_at(&self, point: Vec3) -> Material {
+        (self.fun)(point)
+    }
+    fn get_lights(&self) -> Vec<Light> {
+        self.obj.get_lights()
+    }
+    fn node_count(&self) -> u32 {
+        self.obj.node_count() + 1
+    }
+}

src/structs/cam.rs

@@ -0,0 +1,323 @@
+use crate::structs::{Vec3, ColorVec, COLOR_ZERO, Ray, Y, X};
+use crate::object::Obj;
+use crate::material::{Material, SurfaceType};
+use crate::consts::*;
+use crate::light::Light;
+use std::f64::consts::PI;
+use rand::prelude::*;
+use crossbeam_channel::unbounded;
+use std::ops::{Mul, Add, Sub};
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Cam {
+    pos: Vec3,
+    dir: Vec3,
+    up: Vec3,
+    right: Vec3
+}
+
+impl Cam {
+    pub fn new(pos: Vec3, dir: Vec3, fov: f64) -> Cam {
+        let right = (dir ^ UP).unit() * fov;
+        let up = (right ^ dir).unit() * fov;
+        let dir = dir.unit();
+        Cam { pos, dir, up, right }
+    }
+
+    pub fn new_pointing(pos: Vec3, pointing: Vec3, fov: f64) -> Cam {
+        Self::new(pos, pointing - pos, fov)
+    }
+}
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+struct RayData {
+    pub pos: Vec3,
+    pub dir: Vec3,
+    pub normal: Vec3,
+    pub material: Material,
+    pub dist: f64,
+    pub steps: u32
+}
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+struct LightRayData<'a> {
+    pub light: &'a Light,
+    pub dir: Vec3,
+    pub dist: f64
+}
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+struct Slice {
+    pub x: usize,
+    pub y: usize,
+    pub w: usize,
+    pub h: usize
+}
+
+type CastData = [Option<RayData>; MAX_BOUNCES as usize];
+const NONE_CAST_DATA: CastData = [None; MAX_BOUNCES as usize];
+pub type Image = [u8; IMG_BYTE_SIZE];
+
+fn l2s(l: f64) -> u8 {
+    (f64::powf(f64::clamp(l, 0., 1.), 1./2.2) * 255. + 0.5) as u8
+}
+
+fn randab(a: f64, b: f64) -> f64 {
+    random::<f64>() * (b-a) + a
+}
+
+fn shoot_ray<T: Obj>(scene: &T, pos: Vec3, dir: Vec3) -> Option<RayData> {
+    let dir = dir.unit();
+    let mut steps = 0;
+    let mut dist = 0.;
+    let ray = Ray::new(pos, dir);
+    while steps < MAX_STEPS && dist < MAX_DIST {
+        let pos = ray.point(dist);
+        let d = scene.distance_to(pos);
+        if d < EPSILON {
+            return Some(RayData {
+                pos,
+                dir,
+                normal: scene.normal_at(pos),
+                material: scene.material_at(pos),
+                dist,
+                steps
+            })
+        }
+        dist += d;
+        steps += 1;
+    }
+    None
+}
+
+fn shoot_ray_at<T: Obj>(scene: &T, pos: Vec3, dest: Vec3) -> bool {
+    let dir = (dest - pos).unit();
+    let mut steps = 0;
+    let mut dist = EPSILON; //XXX 0.
+    let ray = Ray::new(pos, dir);
+    while steps < MAX_STEPS && dist < MAX_DIST {
+        let point = ray.point(dist);
+        let sd = scene.distance_to(point);
+        let dd = dest.distance_to(point);
+        if sd < 0. {
+            return false;
+        }
+        if dd < sd {
+            return true;
+        }
+        dist += if sd < EPSILON { EPSILON } else { sd };
+        steps += 1;
+    }
+    false
+}
+
+impl Cam {
+    pub fn render<T: Obj + Clone>(&self, scene: &T) -> Image {
+        if THREAD_COUNT > 1 {
+            self.render_multithreaded(scene)
+        } else {
+            self.render_singlethreaded(scene)
+        }
+    }
+
+    pub fn render_singlethreaded<T: Obj>(&self, scene: &T) -> Image {
+        let lights = scene.get_lights();
+        let mut pixels = [0; IMG_BYTE_SIZE];
+        for y in 0..IMG_HEIGHT {
+            for x in 0..IMG_WIDTH {
+                let field_pos = (x + y*IMG_WIDTH) * 3;
+                self.render_single(scene, &lights, x, y, &mut pixels[field_pos..(field_pos+3)]);
+            }
+        }
+        pixels
+    }
+
+    pub fn render_multithreaded<T: Obj + Clone>(&self, scene: &T) -> Image {
+        let mut pixels = [0; IMG_BYTE_SIZE];
+
+        crossbeam::scope(|scope| {
+            // initialize everything for render
+            let lights = scene.get_lights();
+            let slice_width = IMG_WIDTH / THREAD_COUNT;
+            let slice_height = IMG_HEIGHT / SLICES_PER_THREAD;
+            let slice_byte_size = slice_width * slice_height * 3;
+
+            let (slice_tx, slice_rx) = unbounded::<Slice>();
+            let (data_tx, data_rx) = unbounded();
+
+            // spawn a bunch of threads
+            for _i in 0..THREAD_COUNT {
+                let cam = self.clone();
+                let scene = scene.clone();
+                let lights = lights.clone();
+                let slice_rx = slice_rx.clone();
+                let data_tx = data_tx.clone();
+
+                scope.spawn(move |_| {
+                    for slice in slice_rx.iter() {
+                        let mut pixels = Vec::with_capacity(slice_byte_size);
+                        pixels.resize(slice_byte_size, 0 as u8);
+                        let data = pixels.as_mut_slice();
+
+                        for x in 0..slice.w {
+                            for y in 0..slice.h {
+                                let i = (x + y*slice_width) * 3;
+                                cam.render_single(&scene, &lights, slice.x+x, slice.y+y, &mut data[i..(i+3)]);
+                            }
+                        }
+
+                        data_tx.send((slice, pixels)).unwrap();
+                    }
+                });
+            }
+
+            // don't forget to close the ends of the channel we don't use
+            drop(slice_rx);
+            drop(data_tx);
+
+            // send the slice data
+            let mut y = 0;
+            while y < IMG_HEIGHT {
+                let mut x = 0;
+                while x < IMG_WIDTH {
+                    let w = usize::min(slice_width, IMG_WIDTH - x);
+                    let h = usize::min(slice_height, IMG_HEIGHT - y);
+                    slice_tx.send(Slice { x, y, w, h }).unwrap();
+                    x += slice_width;
+                }
+                y += slice_height;
+            }
+            drop(slice_tx);
+
+            // merge stuff as we get it
+            for (slice, data) in data_rx {
+                let data = data.as_slice();
+                for sy in 0..slice.h {
+                    let py = sy + slice.y;
+                    for sx in 0..slice.w {
+                        let px = sx + slice.x;
+
+                        let pi = (px + py * IMG_WIDTH) * 3;
+                        let si = (sx + sy * slice_width) * 3;
+                        pixels[pi + 0] = data[si + 0];
+                        pixels[pi + 1] = data[si + 1];
+                        pixels[pi + 2] = data[si + 2];
+                    }
+                }
+            }
+        }).unwrap();
+
+        pixels
+    }
+
+    fn render_single<T: Obj>(&self, scene: &T, lights: &Vec<Light>, x: usize, y: usize, field: &mut [u8]) -> () {
+        let mut color = COLOR_ZERO;
+        for xs in 0..SUPERSAMPLING {
+            for ys in 0..SUPERSAMPLING {
+                for _i in 0..RAYS_PER_PIXEL {
+                    let cast = self.cast_single(scene, x * SUPERSAMPLING + xs, y * SUPERSAMPLING + ys);
+                    color = color + self.color_single(scene, &lights, cast);
+                }
+            }
+        }
+        color = color / (SUPERSAMPLING * SUPERSAMPLING * RAYS_PER_PIXEL) as f64;
+
+        field[0] = l2s(color.v(0));
+        field[1] = l2s(color.v(1));
+        field[2] = l2s(color.v(2));
+    }
+
+    fn color_single<'a, T: Obj>(&self, scene: &'a T, lights: &'a Vec<Light>, cast: CastData) -> ColorVec {
+        let mut color = COLOR_ZERO;
+        let mut dist: f64 = 1.;
+
+        for i in (0..MAX_BOUNCES).rev() {
+            if let Some(ray) = cast[i as usize] {
+                match ray.material.surface() {
+                    SurfaceType::Diffuse => {
+                        color = color / dist.powf(DIST_POWER);
+
+                        color = ray.material.reflection() * color + ray.material.emission();
+                        for lightray in Self::get_lights(scene, lights, ray.pos + (scene.normal_at(ray.pos)- ray.dir) * LIGHT_EPSILON) {
+                            let angle = f64::max(0., ray.dir.angle_to(lightray.dir).abs() / PI / ANGLE_CORRECTION + 1. - ANGLE_FIX_CORRECTION);
+                            let dist = f64::max(1., lightray.dist / DIST_CORRECTION - DIST_FIX_CORRECTION);
+                            color = color + ray.material.reflection() * lightray.light.color() / angle.powf(ANGLE_POWER) / dist.powf(DIST_POWER);
+                        }
+                    },
+                    SurfaceType::Reflective => {
+                        color = ray.material.reflection() * color + ray.material.emission();
+                    },
+                    SurfaceType::Stop => {
+                        color = ray.material.emission();
+                    }
+                }
+                dist = f64::max(1., ray.dist / DIST_CORRECTION - DIST_FIX_CORRECTION);
+            }
+        }
+
+        color / dist.powf(DIST_POWER)
+    }
+
+    fn cast_single<T: Obj>(&self, scene: &T, x: usize, y: usize) -> CastData {
+        let mut cast = NONE_CAST_DATA;
+        let mut pos = self.pos;
+        let mut dir = self.dir_for(x, y);
+        let mut i = 0;
+        while i < MAX_BOUNCES {
+            let hit = shoot_ray(scene, pos, dir);
+            cast[i as usize] = hit;
+            match hit {
+                Some(hit) => {
+                    i += 1;
+                    pos = hit.pos - dir * EPSILON;
+                    match hit.material.surface() {
+                        SurfaceType::Stop => break,
+                        SurfaceType::Diffuse => {
+                            let r1 = randab(0., 2.*PI);
+                            let r2 = randab(0., 1.);
+                            let w = hit.normal;
+                            let u = ((if w.x().abs() > 0.5 { Y } else { X }) ^ w).unit();
+                            let v = w ^ u;
+                            let x = r1.cos()*r2.sqrt();
+                            let y = r1.sin()*r2.sqrt();
+                            let z = (1.-r2).sqrt();
+                            if w.x().abs() < 0.5 {
+                                dir = (u*x + v*y + w*z).unit();
+                            } else {
+                                dir = (u*x + v*y + w*z).unit();
+                            }
+                        },
+                        SurfaceType::Reflective => {
+                            dir = dir - hit.normal * (hit.normal * dir) * 2.;
+                        }
+                    }
+                },
+                None => break,
+            }
+        }
+        cast
+    }
+
+    fn get_lights<'a, T: Obj>(scene: &'a T, lights: &'a Vec<Light>, pos: Vec3) -> Vec<LightRayData<'a>> {
+        lights
+            .into_iter()
+            .filter_map(|light| {
+                if shoot_ray_at(scene, pos, light.pos()) {
+                    Some(LightRayData {
+                        light,
+                        dir: (pos - light.pos()),
+                        dist: pos.distance_to(light.pos())
+                    })
+                } else {
+                    None
+                }
+            })
+            .collect()
+    }
+
+    fn dir_for(&self, x: usize, y: usize) -> Vec3 {
+        let x = (((x + (IMG_DIM - IMG_WIDTH) * SUPERSAMPLING / 2) as f64) / ((IMG_DIM * SUPERSAMPLING) as f64) - 0.5) * 2.;
+        let y = (0.5 - ((y + (IMG_DIM - IMG_HEIGHT) * SUPERSAMPLING / 2) as f64) / ((IMG_DIM * SUPERSAMPLING) as f64)) * 2.;
+        (self.dir + self.right * x + self.up * y).unit()
+    }
+}

src/structs/mat3.rs

@@ -0,0 +1,217 @@
+use crate::structs::vec3::Vec3;
+use std::ops::{Add, Sub, Neg, Mul, Div};
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Mat3 {
+    a: f64, b: f64, c: f64,
+    d: f64, e: f64, f: f64,
+    g: f64, h: f64, i: f64,
+}
+
+pub const I3: Mat3 = Mat3 {
+    a: 1., b: 0., c: 0.,
+    d: 0., e: 1., f: 0.,
+    g: 0., h: 0., i: 1.,
+};
+
+pub const O3: Mat3 = Mat3 {
+    a: 0., b: 0., c: 0.,
+    d: 0., e: 0., f: 0.,
+    g: 0., h: 0., i: 0.,
+};
+
+pub const COFACTORS: Mat3 = Mat3 {
+    a:  1., b: -1., c:  1.,
+    d: -1., e:  1., f: -1.,
+    g:  1., h: -1., i:  1.,
+};
+
+impl Mat3 {
+    pub const fn new(a: f64, b: f64, c: f64, d: f64, e: f64, f: f64, g: f64, h: f64, i: f64) -> Mat3 {
+        Mat3 {
+            a, b, c,
+            d, e, f,
+            g, h, i,
+        }
+    }
+
+    pub const fn new_cols(a: Vec3, b: Vec3, c: Vec3) -> Mat3 {
+        Mat3 {
+            a: a.x(), b: b.x(), c: c.x(),
+            d: a.y(), e: b.y(), f: c.y(),
+            g: a.z(), h: b.z(), i: c.z(),
+        }
+    }
+
+    pub const fn new_rows(a: Vec3, d: Vec3, g: Vec3) -> Mat3 {
+        Mat3 {
+            a: a.x(), b: a.y(), c: a.z(),
+            d: d.x(), e: d.y(), f: d.z(),
+            g: g.x(), h: g.y(), i: g.z(),
+        }
+    }
+
+    pub const fn a(&self) -> f64 { self.a }
+    pub const fn b(&self) -> f64 { self.b }
+    pub const fn c(&self) -> f64 { self.c }
+    pub const fn d(&self) -> f64 { self.d }
+    pub const fn e(&self) -> f64 { self.e }
+    pub const fn f(&self) -> f64 { self.f }
+    pub const fn g(&self) -> f64 { self.g }
+    pub const fn h(&self) -> f64 { self.h }
+    pub const fn i(&self) -> f64 { self.i }
+
+    pub fn det(&self) -> f64 {
+        let (a, b, c) = (self.a(), self.b(), self.c());
+        let (d, e, f) = (self.d(), self.e(), self.f());
+        let (g, h, i) = (self.g(), self.h(), self.i());
+
+        a*(e*i - f*h) - b*(d*i - f*g) + c*(d*h - e*g)
+    }
+
+    pub const fn trans(&self) -> Mat3 {
+        Mat3 {
+            a: self.a, b: self.d, c: self.g,
+            d: self.b, e: self.e, f: self.h,
+            g: self.c, h: self.f, i: self.i,
+        }
+    }
+
+    pub fn minor(&self) -> Mat3 {
+        let (a, b, c) = (self.a(), self.b(), self.c());
+        let (d, e, f) = (self.d(), self.e(), self.f());
+        let (g, h, i) = (self.g(), self.h(), self.i());
+
+        Mat3 {
+            a: e*i - f*h, b: d*i - f*g, c: d*h - e*g,
+            d: b*i - c*h, e: a*i - c*g, f: a*h - b*g,
+            g: b*f - c*e, h: a*f - c*d, i: a*e - b*d,
+        }
+    }
+
+    pub fn invert(&self) -> Mat3 {
+        (self.minor().member_mul(COFACTORS)).trans() / self.det()
+    }
+
+    pub fn member_mul(&self, other: Mat3) -> Mat3 {
+        Mat3 {
+            a: self.a*other.a, b: self.b*other.b, c: self.c*other.c,
+            d: self.d*other.d, e: self.e*other.e, f: self.f*other.f,
+            g: self.g*other.g, h: self.h*other.h, i: self.i*other.i,
+        }
+    }
+}
+
+impl Add for Mat3 {
+    type Output = Mat3;
+
+    fn add(self, other: Mat3) -> Mat3 {
+        Mat3 {
+            a: self.a+other.a, b: self.b+other.b, c: self.c+other.c,
+            d: self.d+other.d, e: self.e+other.e, f: self.f+other.f,
+            g: self.g+other.g, h: self.h+other.h, i: self.i+other.i,
+        }
+    }
+}
+
+impl Sub for Mat3 {
+    type Output = Mat3;
+
+    fn sub(self, other: Mat3) -> Mat3 {
+        Mat3 {
+            a: self.a-other.a, b: self.b-other.b, c: self.c-other.c,
+            d: self.d-other.d, e: self.e-other.e, f: self.f-other.f,
+            g: self.g-other.g, h: self.h-other.h, i: self.i-other.i,
+        }
+    }
+}
+
+impl Mul<Mat3> for Mat3 {
+    type Output = Mat3;
+
+    fn mul(self, other: Mat3) -> Mat3 {
+        Mat3 {
+            a: self.a*other.a + self.b*other.d + self.c*other.g,
+            b: self.a*other.b + self.b*other.e + self.c*other.h,
+            c: self.a*other.c + self.b*other.f + self.c*other.i,
+
+            d: self.d*other.a + self.e*other.d + self.f*other.g,
+            e: self.d*other.b + self.e*other.e + self.f*other.h,
+            f: self.d*other.c + self.e*other.f + self.e*other.i,
+
+            g: self.g*other.a + self.h*other.d + self.i*other.g,
+            h: self.g*other.b + self.h*other.e + self.i*other.h,
+            i: self.g*other.c + self.h*other.f + self.i*other.i,
+        }
+    }
+}
+
+impl Neg for Mat3 {
+    type Output = Mat3;
+
+    fn neg(self) -> Mat3 {
+        Mat3 {
+            a: -self.a, b: -self.b, c: -self.c,
+            d: -self.d, e: -self.e, f: -self.f,
+            g: -self.h, h: -self.h, i: -self.i,
+        }
+    }
+}
+
+impl Mul<f64> for Mat3 {
+    type Output = Mat3;
+
+    fn mul(self, k: f64) -> Mat3 {
+        Mat3 {
+            a: self.a*k, b: self.b*k, c: self.c*k,
+            d: self.d*k, e: self.e*k, f: self.f*k,
+            g: self.g*k, h: self.g*k, i: self.i*k,
+        }
+    }
+}
+
+impl Mul<Vec3> for Mat3 {
+    type Output = Vec3;
+
+    fn mul(self, other: Vec3) -> Vec3 {
+        let (x, y, z) = (other.x(), other.y(), other.z());
+        let (a, b, c) = (self.a(), self.b(), self.c());
+        let (d, e, f) = (self.d(), self.e(), self.f());
+        let (g, h, i) = (self.g(), self.h(), self.i());
+
+        Vec3::new(
+            a*x + b*y + c*z,
+            d*x + e*y + f*z,
+            g*x + h*y + i*z,
+        )
+    }
+}
+
+impl Mul<Mat3> for Vec3 {
+    type Output = Vec3;
+
+    fn mul(self, other: Mat3) -> Vec3 {
+        let (x, y, z) = (self.x(), self.y(), self.z());
+        let (a, b, c) = (other.a(), other.b(), other.c());
+        let (d, e, f) = (other.d(), other.e(), other.f());
+        let (g, h, i) = (other.g(), other.h(), other.i());
+
+        Vec3::new(
+            x*a + y*d + z*g,
+            x*b + y*e + z*h,
+            x*c + y*f + z*i,
+        )
+    }
+}
+
+impl Div<f64> for Mat3 {
+    type Output = Mat3;
+
+    fn div(self, k: f64) -> Mat3 {
+        Mat3 {
+            a: self.a/k, b: self.b/k, c: self.c/k,
+            d: self.d/k, e: self.e/k, f: self.f/k,
+            g: self.g/k, h: self.h/k, i: self.i/k,
+        }
+    }
+}

src/structs/mat_n.rs

@@ -0,0 +1,133 @@
+use crate::structs::vec_n::Vec;
+use std::ops::{Add, Sub, Neg, Mul, Div};
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Mat<const N: usize> {
+    v: [[f64; N]; N],
+}
+
+impl<const N: usize> Mat<N> {
+    pub const fn new(v: [[f64; N]; N]) -> Mat<N> {
+        Mat { v }
+    }
+
+    pub const fn new_zero() -> Mat<N> {
+        Mat { v: [[0.; N]; N] }
+    }
+
+    pub const fn new_unit() -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        let mut i = 0;
+        while i<N {
+            v[i][i] = 1.;
+            i += 1;
+        }
+        Mat { v }
+    }
+
+    pub const fn new_one(i: usize, j: usize) -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        v[i][j] = 1.;
+        Mat { v }
+    }
+
+    pub const fn new_from_diagonal(d: Vec<N>) -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        let mut i = 0;
+        while i<N {
+            v[i][i] = d.v(i);
+            i += 1;
+        }
+        Mat { v }
+    }
+
+    pub const fn v(&self, i: usize, j: usize) -> f64 {
+        self.v[i][j]
+    }
+}
+
+impl<const N: usize> Add for Mat<N> {
+    type Output = Mat<N>;
+
+    fn add(self, other: Mat<N>) -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        for i in 0..N {
+            for j in 0..N {
+                v[i][j] = self.v[i][j] + other.v[i][j];
+            }
+        }
+        Mat { v }
+    }
+}
+
+impl<const N: usize> Sub for Mat<N> {
+    type Output = Mat<N>;
+
+    fn sub(self, other: Mat<N>) -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        for i in 0..N {
+            for j in 0..N {
+                v[i][j] = self.v[i][j] - other.v[i][j];
+            }
+        }
+        Mat { v }
+    }
+}
+
+impl<const N: usize> Neg for Mat<N> {
+    type Output = Mat<N>;
+
+    fn neg(self) -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        for i in 0..N {
+            for j in 0..N {
+                v[i][j] = -self.v[i][j];
+            }
+        }
+        Mat { v }
+    }
+}
+
+impl<const N: usize> Mul<f64> for Mat<N> {
+    type Output = Mat<N>;
+
+    fn mul(self, k: f64) -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        for i in 0..N {
+            for j in 0..N {
+                v[i][j] = self.v[i][j] * k;
+            }
+        }
+        Mat { v }
+    }
+}
+
+impl<const N: usize> Mul<Vec<N>> for Mat<N> {
+    type Output = Vec<N>;
+
+    fn mul(self, other: Vec<N>) -> Vec<N> {
+        let mut v = [0.; N];
+        for i in 0..N {
+            let mut x = 0.;
+            for j in 0..N {
+                x += self.v[i][j] * other.v(i);
+            }
+            v[i] = x;
+        }
+        Vec::new(v)
+    }
+}
+
+impl<const N: usize> Div<f64> for Mat<N> {
+    type Output = Mat<N>;
+
+    fn div(self, k: f64) -> Mat<N> {
+        let mut v = [[0.; N]; N];
+        for i in 0..N {
+            for j in 0..N {
+                v[i][j] = self.v[i][j] / k;
+            }
+        }
+        Mat { v }
+    }
+}

src/structs/mod.rs

@@ -0,0 +1,20 @@
+use crate::consts::COLOR_CHANNELS;
+
+mod vec3;
+mod mat3;
+mod vec_n;
+mod mat_n;
+
+mod ray;
+mod cam;
+
+pub use vec3::{Vec3, X, Y, Z, O};
+pub use mat3::{Mat3, I3, O3};
+pub use vec_n::Vec;
+pub use mat_n::Mat;
+pub use ray::Ray;
+pub use cam::{Cam, Image};
+
+pub type ColorVec = Vec<COLOR_CHANNELS>;
+pub const COLOR_ZERO: ColorVec = ColorVec::new_zero();
+pub type ColorMat = Mat<COLOR_CHANNELS>;

src/structs/ray.rs

@@ -0,0 +1,28 @@
+use crate::structs::Vec3;
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Ray {
+    from: Vec3,
+    dir: Vec3,
+}
+
+impl Ray {
+    pub fn new(from: Vec3, dir: Vec3) -> Ray {
+        Ray { from, dir }
+    }
+
+    pub fn new_from_to(from: Vec3, to: Vec3) -> Ray {
+        Ray { from, dir: to - from }
+    }
+
+    pub fn unit(&self) -> Ray {
+        Ray { from: self.from, dir: self.dir.unit() }
+    }
+
+    pub fn from(&self) -> Vec3 { self.from }
+    pub fn dir(&self) -> Vec3 { self.dir }
+
+    pub fn point(&self, steps: f64) -> Vec3 {
+        self.from + self.dir * steps
+    }
+}

src/structs/vec3.rs

@@ -0,0 +1,127 @@
+use std::ops::{Add, Sub, Neg, Mul, Div, BitXor};
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Vec3 {
+    x: f64,
+    y: f64,
+    z: f64,
+}
+
+pub const X: Vec3 = Vec3 { x: 1., y: 0., z: 0. };
+pub const Y: Vec3 = Vec3 { x: 0., y: 1., z: 0. };
+pub const Z: Vec3 = Vec3 { x: 0., y: 0., z: 1. };
+pub const O: Vec3 = Vec3 { x: 0., y: 0., z: 0. };
+
+impl Vec3 {
+    pub const fn new(x: f64, y: f64, z: f64) -> Vec3 {
+        Vec3 { x, y, z }
+    }
+
+    pub const fn x(&self) -> f64 { self.x }
+    pub const fn y(&self) -> f64 { self.y }
+    pub const fn z(&self) -> f64 { self.z }
+
+    pub fn magnitude(self) -> f64 {
+        f64::sqrt(self.x*self.x + self.y*self.y + self.z*self.z)
+    }
+
+    pub fn unit(self) -> Vec3 {
+        self / self.magnitude()
+    }
+
+    pub fn distance_to(self, other: Vec3) -> f64 {
+        (self - other).magnitude()
+    }
+
+    pub fn angle_to(self, other: Vec3) -> f64 {
+        ((self * other) / (self.magnitude() * other.magnitude())).acos()
+    }
+
+    pub fn member_mul(self, other: Vec3) -> Vec3 {
+        Vec3 {
+            x: self.x*other.x,
+            y: self.y*other.y,
+            z: self.z*other.z,
+        }
+    }
+}
+
+impl Add for Vec3 {
+    type Output = Vec3;
+
+    fn add(self, other: Vec3) -> Vec3 {
+        Vec3 {
+            x: self.x + other.x,
+            y: self.y + other.y,
+            z: self.z + other.z,
+        }
+    }
+}
+
+impl Sub for Vec3 {
+    type Output = Vec3;
+
+    fn sub(self, other: Vec3) -> Vec3 {
+        Vec3 {
+            x: self.x - other.x,
+            y: self.y - other.y,
+            z: self.z - other.z,
+        }
+    }
+}
+
+impl Neg for Vec3 {
+    type Output = Vec3;
+
+    fn neg(self) -> Vec3 {
+        Vec3 {
+            x: -self.x,
+            y: -self.y,
+            z: -self.z,
+        }
+    }
+}
+
+impl BitXor for Vec3 {
+    type Output = Vec3;
+
+    fn bitxor(self, other: Vec3) -> Vec3 {
+        Vec3 {
+            x: self.y * other.z - self.z * other.y,
+            y: self.z * other.x - self.x * other.z,
+            z: self.x * other.y - self.y * other.x,
+        }
+    }
+}
+
+impl Mul<Vec3> for Vec3 {
+    type Output = f64;
+
+    fn mul(self, other: Vec3) -> f64 {
+        self.x*other.x + self.y*other.y + self.z*other.z
+    }
+}
+
+impl Mul<f64> for Vec3 {
+    type Output = Vec3;
+
+    fn mul(self, k: f64) -> Vec3 {
+        Vec3 {
+            x: self.x * k,
+            y: self.y * k,
+            z: self.z * k,
+        }
+    }
+}
+
+impl Div<f64> for Vec3 {
+    type Output = Vec3;
+
+    fn div(self, k: f64) -> Vec3 {
+        Vec3 {
+            x: self.x / k,
+            y: self.y / k,
+            z: self.z / k,
+        }
+    }
+}

src/structs/vec_n.rs

@@ -0,0 +1,111 @@
+use std::ops::{Add, Sub, Neg, Mul, Div};
+
+#[derive(Debug, Copy, Clone, PartialEq)]
+pub struct Vec<const N: usize> {
+    v: [f64; N],
+}
+
+impl<const N: usize> Vec<N> {
+    pub const fn new(v: [f64; N]) -> Vec<N> {
+        Vec { v }
+    }
+
+    pub const fn new_zero() -> Vec<N> {
+        Vec { v: [0.; N] }
+    }
+
+    pub const fn new_one(i: usize) -> Vec<N> {
+        let mut v = [0.; N];
+        v[i] = 1.;
+        Vec { v }
+    }
+
+    pub const fn v(&self, i: usize) -> f64 {
+        self.v[i]
+    }
+
+    pub fn magnitude(self) -> f64 {
+        let mut x = 0.;
+        for i in 0..N {
+            let v = self.v[i];
+            x += v*v;
+        }
+        f64::sqrt(x)
+    }
+
+    pub fn unit(self) -> Vec<N> {
+        self / self.magnitude()
+    }
+
+    pub fn distance_to(self, other: Vec<N>) -> f64 {
+        (self - other).magnitude()
+    }
+
+    pub fn member_mul(self, other: Vec<N>) -> Vec<N> {
+        let mut v = [0.; N];
+        for i in 0..N {
+            v[i] = self.v[i] * other.v[i];
+        }
+        Vec { v }
+    }
+}
+
+impl<const N: usize> Add for Vec<N> {
+    type Output = Vec<N>;
+
+    fn add(self, other: Vec<N>) -> Vec<N> {
+        let mut v = [0.; N];
+        for i in 0..N {
+            v[i] = self.v[i] + other.v[i];
+        }
+        Vec { v }
+    }
+}
+
+impl<const N: usize> Sub for Vec<N> {
+    type Output = Vec<N>;
+
+    fn sub(self, other: Vec<N>) -> Vec<N> {
+        let mut v = [0.; N];
+        for i in 0..N {
+            v[i] = self.v[i] - other.v[i];
+        }
+        Vec { v }
+    }
+}
+
+impl<const N: usize> Neg for Vec<N> {
+    type Output = Vec<N>;
+
+    fn neg(self) -> Vec<N> {
+        let mut v = [0.; N];
+        for i in 0..N {
+            v[i] = -self.v[i];
+        }
+        Vec { v }
+    }
+}
+
+impl<const N: usize> Mul<f64> for Vec<N> {
+    type Output = Vec<N>;
+
+    fn mul(self, k: f64) -> Vec<N> {
+        let mut v = [0.; N];
+        for i in 0..N {
+            v[i] = self.v[i] * k;
+        }
+        Vec { v }
+    }
+}
+
+impl<const N: usize> Div<f64> for Vec<N> {
+    type Output = Vec<N>;
+
+    fn div(self, k: f64) -> Vec<N> {
+        let mut v = [0.; N];
+        for i in 0..N {
+            v[i] = self.v[i] / k;
+        }
+        Vec { v }
+    }
+}