Add dielectric material

src/material.c

@@ -1,5 +1,24 @@
 #include "material.h"
 
+#include <math.h>
+#include <stdio.h>
+
+static vec3_t reflect(vec3_t v, vec3_t n)
+{
+	const vec3_t dn = vec3_scale(n, -2 * vec3_dot(n, v));
+	return vec3_unit(vec3_add(v, dn));
+}
+
+static vec3_t refract(vec3_t v, vec3_t n, double eta_ratio)
+{
+	const double cos_theta = fmin(-vec3_dot(v, n), 1.0);
+	const vec3_t perp
+	    = vec3_scale(vec3_add(v, vec3_scale(n, cos_theta)), eta_ratio);
+	const vec3_t para
+	    = vec3_scale(n, -sqrt(fabs(1.0 - vec3_len_squared(perp))));
+	return vec3_unit(vec3_add(perp, para));
+}
+
 bool scatter_lambertian(
     material_params_t params, hit_t hit, rng_t *rng, ray_t *ray,
     vec3_t *atten_out)
@@ -17,10 +36,30 @@ bool scatter_reflective(
 {
 	(void)rng;
 
-	const vec3_t dn
-	    = vec3_scale(hit.normal, -2 * vec3_dot(hit.normal, ray->dir));
 	ray->orig = hit.point;
-	ray->dir = vec3_unit(vec3_add(ray->dir, dn));
+	ray->dir = reflect(ray->dir, hit.normal);
 	*atten_out = params.reflective.tint;
 	return true;
 }
+
+bool scatter_dielectric(
+    material_params_t params, hit_t hit, rng_t *rng, ray_t *ray,
+    vec3_t *atten_out)
+{
+	const double eta_ratio
+	    = hit.front ? 1.0 / params.dielectric.eta : params.dielectric.eta;
+	const vec3_t n = hit.front ? hit.normal : vec3_scale(hit.normal, -1);
+	const double cos_theta = fmin(-vec3_dot(ray->dir, n), 1.0);
+	const double sin_theta = sqrt(1.0 - cos_theta * cos_theta);
+	const double r = pow((1 - eta_ratio) / (1 + eta_ratio), 2);
+	const double reflectance = r + (1 - r) * pow(1 - cos_theta, 5);
+
+	ray->orig = hit.point;
+	if (sin_theta * eta_ratio > 1.0 || rng_canon(rng) < reflectance)
+		ray->dir = reflect(ray->dir, n);
+	else
+		ray->dir = refract(ray->dir, n, eta_ratio);
+	*atten_out = (vec3_t) { 1.0, 1.0, 1.0 };
+
+	return true;
+}