raytracers/rtiow/renderer/src/triangles.rs

use std::f32::EPSILON;

use stl::STL;

use crate::{
    aabb::AABB,
    hitable::{Hit, HitRecord},
    material::Material,
    ray::Ray,
    vec3::{cross, dot, Vec3},
};

#[derive(Debug)]
pub struct Triangle {
    normal: Vec3,
    verts: [Vec3; 3],
    // Precomputed data
}

#[derive(Debug)]
pub struct Triangles<M>
where
    M: Material,
{
    pub triangles: Vec<Triangle>,
    pub bbox: AABB,
    material: M,
}

impl<M> Triangles<M>
where
    M: Material,
{
    pub fn new(stl: &STL, material: M, scale_factor: f32) -> Triangles<M> {
        let triangles: Vec<_> = stl
            .triangles
            .iter()
            .map(|t| {
                let v0 = t.verts[0] * scale_factor;
                let v1 = t.verts[1] * scale_factor;
                let v2 = t.verts[2] * scale_factor;
                Triangle {
                    normal: t.normal,
                    verts: [v0, v1, v2],
                }
            })
            .collect();
        let (min, max) = triangles.iter().fold(
            (Vec3::from(f32::MAX), Vec3::from(f32::MIN)),
            |(min, max), t| {
                let t_min_x = t.verts[0].x.min(t.verts[1].x.min(t.verts[2].x));
                let t_min_y = t.verts[0].y.min(t.verts[1].y.min(t.verts[2].y));
                let t_min_z = t.verts[0].z.min(t.verts[1].z.min(t.verts[2].z));
                let t_max_x = t.verts[0].x.max(t.verts[1].x.max(t.verts[2].x));
                let t_max_y = t.verts[0].y.max(t.verts[1].y.max(t.verts[2].y));
                let t_max_z = t.verts[0].z.max(t.verts[1].z.max(t.verts[2].z));

                (
                    Vec3::from([min.x.min(t_min_x), min.y.min(t_min_y), min.z.min(t_min_z)]),
                    Vec3::from([max.x.max(t_max_x), max.y.max(t_max_y), max.z.max(t_max_z)]),
                )
            },
        );
        let bbox = AABB::new(min, max);
        Triangles {
            triangles,
            bbox,
            material,
        }
    }
}

/// Based on https://www.scratchapixel.com/lessons/3d-basic-rendering/ray-tracing-rendering-a-triangle/moller-trumbore-ray-triangle-intersection.html
fn ray_triangle_intersect_moller_trumbore(r: Ray, tri: &Triangle) -> Option<RayTriangleResult> {
    // #ifdef MOLLER_TRUMBORE
    //     Vec3f v0v1 = v1 - v0;
    //     Vec3f v0v2 = v2 - v0;
    //     Vec3f pvec = dir.crossProduct(v0v2);
    //     float det = v0v1.dotProduct(pvec);
    // #ifdef CULLING
    //     // if the determinant is negative, the triangle is 'back facing'
    //     // if the determinant is close to 0, the ray misses the triangle
    //     if (det < kEpsilon) return false;
    // #else
    //     // ray and triangle are parallel if det is close to 0
    //     if (fabs(det) < kEpsilon) return false;
    // #endif
    //     float invDet = 1 / det;
    //
    //     Vec3f tvec = orig - v0;
    //     u = tvec.dotProduct(pvec) * invDet;
    //     if (u < 0 || u > 1) return false;
    //
    //     Vec3f qvec = tvec.crossProduct(v0v1);
    //     v = dir.dotProduct(qvec) * invDet;
    //     if (v < 0 || u + v > 1) return false;
    //
    //     t = v0v2.dotProduct(qvec) * invDet;
    //
    let v0 = tri.verts[0];
    let v1 = tri.verts[1];
    let v2 = tri.verts[2];

    let v0v1 = v1 - v0;
    let v0v2 = v2 - v0;
    let p = cross(r.direction, v0v2);
    let det = dot(v0v1, p);
    if det < EPSILON {
        return None;
    }

    let inv_det = 1. / det;

    let t = r.origin - v0;
    let u = dot(t, p) * inv_det;
    if u < 0. || u > 1. {
        return None;
    }

    let q = cross(t, v0v1);
    let v = dot(r.direction, q) * inv_det;
    if v < 0. || u + v > 1. {
        return None;
    }
    let t = dot(v0v2, q) * inv_det;

    if t > EPSILON {
        return Some(RayTriangleResult {
            t,
            p: r.origin + r.direction * t,
        });
    }
    None
}

impl<M> Hit for Triangles<M>
where
    M: Material,
{
    fn hit(&self, r: Ray, _t_min: f32, _t_max: f32) -> Option<HitRecord> {
        // TODO(wathiede): add an acceleration structure to more cheaply skip some triangles.
        for tri in &self.triangles {
            if let Some(RayTriangleResult { t, p }) = ray_triangle_intersect_moller_trumbore(r, tri)
            {
                //if let Some(RayTriangleResult { t, p }) = ray_triangle_intersect_geometric(r, tri) {
                // We don't support UV (yet?).
                let uv = (0.5, 0.5);
                return Some(HitRecord {
                    t,
                    uv,
                    p,
                    normal: tri.normal,
                    material: &self.material,
                });
            }
        }
        None
    }

    fn bounding_box(&self, _t_min: f32, _t_max: f32) -> Option<AABB> {
        Some(self.bbox)
    }
}

struct RayTriangleResult {
    t: f32,
    p: Vec3,
}