120 files changed, 2013 insertions, 1404 deletions

diff --git a/libdimension/Makefile.am b/libdimension/Makefile.am
index c488016..5c14363 100644
--- a/libdimension/Makefile.am
+++ b/libdimension/Makefile.am
@@ -17,129 +17,144 @@
 ## along with this program.  If not, see <http://www.gnu.org/licenses/>. ##
 ###########################################################################
 
-SUBDIRS = .                                                                    \
-          bench                                                                \
-          tests
+SUBDIRS = . \
+  bench \
+  tests
 
 AM_CFLAGS = $(PTHREAD_CFLAGS)
 # Make dmnsn_error() backtraces useful
 AM_LDFLAGS = $(PTHREAD_LIBS) -rdynamic
 CC = $(PTHREAD_CC)
 
-nobase_include_HEADERS = dimension.h                                           \
-                         dimension/array.h                                     \
-                         dimension/camera.h                                    \
-                         dimension/cameras.h                                   \
-                         dimension/canvas.h                                    \
-                         dimension/color.h                                     \
-                         dimension/compiler.h                                  \
-                         dimension/csg.h                                       \
-                         dimension/dictionary.h                                \
-                         dimension/error.h                                     \
-                         dimension/finish.h                                    \
-                         dimension/finishes.h                                  \
-                         dimension/future.h                                    \
-                         dimension/geometry.h                                  \
-                         dimension/gl.h                                        \
-                         dimension/interior.h                                  \
-                         dimension/light.h                                     \
-                         dimension/lights.h                                    \
-                         dimension/malloc.h                                    \
-                         dimension/map.h                                       \
-                         dimension/object.h                                    \
-                         dimension/objects.h                                   \
-                         dimension/pattern.h                                   \
-                         dimension/patterns.h                                  \
-                         dimension/pigment.h                                   \
-                         dimension/pigments.h                                  \
-                         dimension/png.h                                       \
-                         dimension/polynomial.h                                \
-                         dimension/pool.h                                      \
-                         dimension/ray_trace.h                                 \
-                         dimension/scene.h                                     \
-                         dimension/tcolor.h                                    \
-                         dimension/texture.h                                   \
-                         dimension/timer.h
+nobase_include_HEADERS = dimension.h \
+  dimension/base.h \
+  dimension/base/compiler.h \
+  dimension/base/common.h \
+  dimension/base/error.h \
+  dimension/base/malloc.h \
+  dimension/base/pool.h \
+  dimension/base/array.h \
+  dimension/base/dictionary.h \
+  dimension/concurrency.h \
+  dimension/concurrency/future.h \
+  dimension/platform.h \
+  dimension/platform/timer.h \
+  dimension/math.h \
+  dimension/math/aabb.h \
+  dimension/math/matrix.h \
+  dimension/math/ray.h \
+  dimension/math/scalar.h \
+  dimension/math/vector.h \
+  dimension/color.h \
+  dimension/color/color.h \
+  dimension/color/tcolor.h \
+  dimension/canvas.h \
+  dimension/canvas/canvas.h \
+  dimension/canvas/gl.h \
+  dimension/canvas/png.h \
+  dimension/pattern.h \
+  dimension/pattern/pattern.h \
+  dimension/pattern/patterns.h \
+  dimension/pattern/map.h \
+  dimension/model.h \
+  dimension/model/camera.h \
+  dimension/model/cameras.h \
+  dimension/model/csg.h \
+  dimension/model/finishes.h \
+  dimension/model/finish.h \
+  dimension/model/interior.h \
+  dimension/model/light.h \
+  dimension/model/lights.h \
+  dimension/model/object.h \
+  dimension/model/objects.h \
+  dimension/model/pigment.h \
+  dimension/model/pigments.h \
+  dimension/model/scene.h \
+  dimension/model/texture.h \
+  dimension/render.h \
+  dimension/render/render.h
 
 lib_LTLIBRARIES = libdimension.la
 
-libdimension_la_SOURCES = $(nobase_include_HEADERS)                            \
-                          array.c                                              \
-                          bvh.c                                                \
-                          bvh.h                                                \
-                          camera.c                                             \
-                          canvas.c                                             \
-                          canvas_pigment.c                                     \
-                          checker.c                                            \
-                          compiler-internal.h                                  \
-                          cone.c                                               \
-                          cube.c                                               \
-                          csg.c                                                \
-                          dictionary.c                                         \
-                          dimension-internal.h                                 \
-                          error.c                                              \
-                          finish.c                                             \
-                          future.c                                             \
-                          future-internal.h                                    \
-                          geometry.c                                           \
-                          gradient.c                                           \
-                          inline.c                                             \
-                          interior.c                                           \
-                          lambertian.c                                         \
-                          leopard.c                                            \
-                          light.c                                              \
-                          malloc.c                                             \
-                          map.c                                                \
-                          object.c                                             \
-                          pattern.c                                            \
-                          perspective.c                                        \
-                          phong.c                                              \
-                          pigment.c                                            \
-                          pigment_map.c                                        \
-                          plane.c                                              \
-                          platform.c                                           \
-                          platform.h                                           \
-                          profile.h                                            \
-                          point_light.c                                        \
-                          polynomial.c                                         \
-                          pool.c                                               \
-                          prtree.c                                             \
-                          prtree.h                                             \
-                          ray_trace.c                                          \
-                          reflection.c                                         \
-                          rgba.c                                               \
-                          rgba.h                                               \
-                          scene.c                                              \
-                          solid_pigment.c                                      \
-                          sphere.c                                             \
-                          texture.c                                            \
-                          threads.c                                            \
-                          threads.h                                            \
-                          timer.c                                              \
-                          torus.c                                              \
-                          triangle.c                                           \
-                          triangle_fan.c
+libdimension_la_SOURCES = $(nobase_include_HEADERS) \
+  base/array.c \
+  base/dictionary.c \
+  base/error.c \
+  base/inline.c \
+  base/malloc.c \
+  base/pool.c \
+  bvh/bvh.c \
+  bvh/prtree.c \
+  canvas/canvas.c \
+  canvas/rgba.c \
+  concurrency/future.c \
+  concurrency/threads.c \
+  dimension.h \
+  internal.h \
+  internal/bvh.h \
+  internal/compiler.h \
+  internal/future.h \
+  internal/platform.h \
+  internal/polynomial.h \
+  internal/profile.h \
+  internal/prtree.h \
+  internal/rgba.h \
+  internal/threads.h \
+  math/matrix.c \
+  math/polynomial.c \
+  model/camera.c \
+  model/cameras/perspective.c \
+  model/finish.c \
+  model/finishes/lambertian.c \
+  model/finishes/phong.c \
+  model/finishes/reflection.c \
+  model/interior.c \
+  model/light.c \
+  model/lights/point_light.c \
+  model/object.c \
+  model/objects/cone.c \
+  model/objects/csg.c \
+  model/objects/cube.c \
+  model/objects/plane.c \
+  model/objects/sphere.c \
+  model/objects/torus.c \
+  model/objects/triangle.c \
+  model/objects/triangle_fan.c \
+  model/pigment.c \
+  model/pigments/canvas_pigment.c \
+  model/pigments/pigment_map.c \
+  model/pigments/solid_pigment.c \
+  model/scene.c \
+  model/texture.c \
+  pattern/checker.c \
+  pattern/gradient.c \
+  pattern/leopard.c \
+  pattern/map.c \
+  pattern/pattern.c \
+  platform/platform.c \
+  platform/timer.c \
+  render/render.c
 libdimension_la_CFLAGS  = $(AM_CFLAGS)
 libdimension_la_LDFLAGS = -version-info 0:0:0 -no-undefined $(AM_LDFLAGS)
 libdimension_la_LIBADD  =
 
 if PNG
-libdimension_la_SOURCES += png.c
+libdimension_la_SOURCES += canvas/png.c
 libdimension_la_CFLAGS  += $(libpng_CFLAGS)
 libdimension_la_LIBADD  += $(libpng_LIBS)
 else
-libdimension_la_SOURCES += png-stubs.c
+libdimension_la_SOURCES += canvas/png-stubs.c
 endif
 
 if GL
-libdimension_la_SOURCES += gl.c
+libdimension_la_SOURCES += canvas/gl.c
 libdimension_la_LIBADD  += -lGL
 else
-libdimension_la_SOURCES += gl-stubs.c
+libdimension_la_SOURCES += canvas/gl-stubs.c
 endif
 
 if PROFILE
-libdimension_la_SOURCES += profile.c
+libdimension_la_SOURCES += base/profile.c
 endif
 
 pkgconfigdir = $(libdir)/pkgconfig
diff --git a/libdimension/array.c b/libdimension/base/array.c
index deea035..d8faed7 100644
--- a/libdimension/array.c
+++ b/libdimension/base/array.c
@@ -23,7 +23,7 @@
  * Non-inline array functions.
  */
 
-#include "dimension.h"
+#include "dimension/base.h"
 
 void
 dmnsn_array_cleanup(void *ptr)
diff --git a/libdimension/dictionary.c b/libdimension/base/dictionary.c
index 84ebedf..6e99abd 100644
--- a/libdimension/dictionary.c
+++ b/libdimension/base/dictionary.c
@@ -23,7 +23,7 @@
  * Associative arrays, implemented with PATRICIA tries.
  */
 
-#include "dimension.h"
+#include "dimension/base.h"
 
 struct dmnsn_dictionary {
   size_t obj_size;       ///< The size of the objects in the trie.
diff --git a/libdimension/error.c b/libdimension/base/error.c
index e9e7cf1..6b8d18e 100644
--- a/libdimension/error.c
+++ b/libdimension/base/error.c
@@ -23,7 +23,8 @@
  * Error handling.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/platform.h"
 #include <errno.h>
 #include <pthread.h>
 #include <stdatomic.h>
diff --git a/libdimension/inline.c b/libdimension/base/inline.c
index b0622fe..b0622fe 100644
--- a/libdimension/inline.c
+++ b/libdimension/base/inline.c
diff --git a/libdimension/malloc.c b/libdimension/base/malloc.c
index 2e4969d..6aaf68c 100644
--- a/libdimension/malloc.c
+++ b/libdimension/base/malloc.c
@@ -23,7 +23,7 @@
  * Dynamic memory.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
 #include <stdlib.h>
 #include <string.h>
 #include <stdatomic.h>
diff --git a/libdimension/pool.c b/libdimension/base/pool.c
index a7d6b4f..4bfdd7d 100644
--- a/libdimension/pool.c
+++ b/libdimension/base/pool.c
@@ -23,7 +23,8 @@
  * Memory pool implementation.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/concurrency.h"
 #include <stdatomic.h>
 
 /// Number of pointers per block, we want a block to fit in a single page.
diff --git a/libdimension/profile.c b/libdimension/base/profile.c
index 164e002..87f27c1 100644
--- a/libdimension/profile.c
+++ b/libdimension/base/profile.c
@@ -23,10 +23,11 @@
  * Branch profile accounting.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/concurrency.h"
+#include <inttypes.h>
 #include <pthread.h>
 #include <stdio.h>
-#include <inttypes.h>
 
 /// Information on one predicted branch.
 typedef struct {
diff --git a/libdimension/bench/future.c b/libdimension/bench/future.c
index fd26363..c8fb529 100644
--- a/libdimension/bench/future.c
+++ b/libdimension/bench/future.c
@@ -17,9 +17,9 @@
  * along with this program.  If not, see <http://www.gnu.org/licenses/>. *
  *************************************************************************/
 
-#include "../platform.c"
-#include "../future.c"
-#include "../threads.c"
+#include "../platform/platform.c"
+#include "../concurrency/future.c"
+#include "../concurrency/threads.c"
 #include <sandglass.h>
 #include <stdlib.h>
 
diff --git a/libdimension/bench/geometry.c b/libdimension/bench/geometry.c
index 68c9885..59d27e3 100644
--- a/libdimension/bench/geometry.c
+++ b/libdimension/bench/geometry.c
@@ -17,17 +17,18 @@
  * along with this program.  If not, see <http://www.gnu.org/licenses/>. *
  *************************************************************************/
 
-#include "dimension.h"
+#include "dimension/math.h"
 #include <sandglass.h>
 #include <stdlib.h>
+#include <stdio.h>
 
 int
 main(void)
 {
   dmnsn_vector vector, vector2;
   dmnsn_matrix matrix, matrix2;
-  dmnsn_line line;
-  dmnsn_bounding_box box;
+  dmnsn_ray ray;
+  dmnsn_aabb box;
   double result;
 
   sandglass_t sandglass;
@@ -74,19 +75,19 @@ main(void)
   });
   printf("dmnsn_rotation_matrix(): %ld\n", sandglass.grains);
 
-  // dmnsn_new_line()
+  // dmnsn_new_ray()
   vector2 = dmnsn_new_vector(3.0, 2.0, 1.0);
   sandglass_bench_fine(&sandglass, {
-    line = dmnsn_new_line(vector, vector2);
+    ray = dmnsn_new_ray(vector, vector2);
   });
-  printf("dmnsn_new_line(): %ld\n", sandglass.grains);
+  printf("dmnsn_new_ray(): %ld\n", sandglass.grains);
 
-  // dmnsn_new_bounding_box()
+  // dmnsn_new_aabb()
   vector2 = dmnsn_new_vector(3.0, 4.0, 5.0);
   sandglass_bench_fine(&sandglass, {
-    box = dmnsn_new_bounding_box(vector, vector2);
+    box = dmnsn_new_aabb(vector, vector2);
   });
-  printf("dmnsn_new_bounding_box(): %ld\n", sandglass.grains);
+  printf("dmnsn_new_aabb(): %ld\n", sandglass.grains);
 
   // dmnsn_vector_add()
   sandglass_bench_fine(&sandglass, {
@@ -175,23 +176,23 @@ main(void)
   });
   printf("dmnsn_transform_normal(): %ld\n", sandglass.grains);
 
-  // dmnsn_transform_line()
+  // dmnsn_transform_ray()
   sandglass_bench_fine(&sandglass, {
-    line = dmnsn_transform_line(matrix, line);
+    ray = dmnsn_transform_ray(matrix, ray);
   });
-  printf("dmnsn_transform_line(): %ld\n", sandglass.grains);
+  printf("dmnsn_transform_ray(): %ld\n", sandglass.grains);
 
-  // dmnsn_transform_bounding_box()
+  // dmnsn_transform_aabb()
   sandglass_bench_fine(&sandglass, {
-    box = dmnsn_transform_bounding_box(matrix, box);
+    box = dmnsn_transform_aabb(matrix, box);
   });
-  printf("dmnsn_transform_bounding_box(): %ld\n", sandglass.grains);
+  printf("dmnsn_transform_aabb(): %ld\n", sandglass.grains);
 
-  // dmnsn_line_point()
+  // dmnsn_ray_point()
   sandglass_bench_fine(&sandglass, {
-    vector = dmnsn_line_point(line, result);
+    vector = dmnsn_ray_point(ray, result);
   });
-  printf("dmnsn_line_point(): %ld\n", sandglass.grains);
+  printf("dmnsn_ray_point(): %ld\n", sandglass.grains);
 
   return EXIT_SUCCESS;
 }
diff --git a/libdimension/bench/polynomial.c b/libdimension/bench/polynomial.c
index bf3d2b0..33fd390 100644
--- a/libdimension/bench/polynomial.c
+++ b/libdimension/bench/polynomial.c
@@ -17,7 +17,8 @@
  * along with this program.  If not, see <http://www.gnu.org/licenses/>. *
  *************************************************************************/
 
-#include "dimension.h"
+#define DMNSN_INLINE extern inline
+#include "../math/polynomial.c"
 #include <sandglass.h>
 #include <stdlib.h>
 
diff --git a/libdimension/bench/prtree.c b/libdimension/bench/prtree.c
index 6c5d3ef..6b5e7c1 100644
--- a/libdimension/bench/prtree.c
+++ b/libdimension/bench/prtree.c
@@ -17,19 +17,19 @@
  * along with this program.  If not, see <http://www.gnu.org/licenses/>. *
  *************************************************************************/
 
-#include "../platform.c"
-#include "../threads.c"
-#include "../future.c"
-#include "../bvh.c"
-#include "../prtree.c"
+#include "../platform/platform.c"
+#include "../concurrency/threads.c"
+#include "../concurrency/future.c"
+#include "../bvh/bvh.c"
+#include "../bvh/prtree.c"
 #include <sandglass.h>
 #include <stdlib.h>
 
 static bool
-dmnsn_fake_intersection_fn(const dmnsn_object *object, dmnsn_line line,
+dmnsn_fake_intersection_fn(const dmnsn_object *object, dmnsn_ray ray,
                            dmnsn_intersection *intersection)
 {
-  intersection->t = (object->bounding_box.min.z - line.x0.z)/line.n.z;
+  intersection->t = (object->aabb.min.z - ray.x0.z)/ray.n.z;
   intersection->normal = dmnsn_x;
   return true;
 }
@@ -40,7 +40,7 @@ dmnsn_fake_inside_fn(const dmnsn_object *object, dmnsn_vector point)
   return true;
 }
 
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_fake_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   dmnsn_vector a, b;
@@ -53,7 +53,7 @@ dmnsn_fake_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
   b.y = 2.0*((double)rand())/RAND_MAX - 1.0;
   b.z = 2.0*((double)rand())/RAND_MAX - 1.0;
 
-  return dmnsn_new_bounding_box(dmnsn_vector_min(a, b), dmnsn_vector_max(a, b));
+  return dmnsn_new_aabb(dmnsn_vector_min(a, b), dmnsn_vector_max(a, b));
 }
 
 static dmnsn_object_vtable dmnsn_fake_vtable = {
@@ -99,7 +99,7 @@ main(void)
   printf("dmnsn_new_bvh(DMNSN_BVH_PRTREE): %ld\n", sandglass.grains);
 
   // dmnsn_bvh_intersection()
-  dmnsn_line ray = dmnsn_new_line(
+  dmnsn_ray ray = dmnsn_new_ray(
     dmnsn_new_vector( 1.0,  1.0, -2.0),
     dmnsn_new_vector(-0.5, -0.5,  1.0)
   );
diff --git a/libdimension/bench/triangle.c b/libdimension/bench/triangle.c
index b2b91d0..4a5a456 100644
--- a/libdimension/bench/triangle.c
+++ b/libdimension/bench/triangle.c
@@ -43,21 +43,21 @@ main(void)
   dmnsn_object_precompute(triangle);
 
   dmnsn_intersection intersection;
-  dmnsn_line line;
+  dmnsn_ray ray;
   bool intersected;
 
   // Intersecting case
-  line = dmnsn_new_line(dmnsn_new_vector(2.0, 1.0, -1.0), dmnsn_z);
+  ray = dmnsn_new_ray(dmnsn_new_vector(2.0, 1.0, -1.0), dmnsn_z);
   sandglass_bench_fine(&sandglass, {
-    intersected = dmnsn_object_intersection(triangle, line, &intersection);
+    intersected = dmnsn_object_intersection(triangle, ray, &intersection);
   });
   dmnsn_assert(intersected, "Didn't intersect");
   printf("dmnsn_triangle_intersection(true): %ld\n", sandglass.grains);
 
   // Non-intersecting case
-  line = dmnsn_new_line(dmnsn_new_vector(3.0, 3.0, -1.0), dmnsn_z);
+  ray = dmnsn_new_ray(dmnsn_new_vector(3.0, 3.0, -1.0), dmnsn_z);
   sandglass_bench_fine(&sandglass, {
-    intersected = dmnsn_object_intersection(triangle, line, &intersection);
+    intersected = dmnsn_object_intersection(triangle, ray, &intersection);
   });
   dmnsn_assert(!intersected, "Intersected");
   printf("dmnsn_triangle_intersection(false): %ld\n", sandglass.grains);
diff --git a/libdimension/bvh.c b/libdimension/bvh/bvh.c
index d1bd0fc..eab2c28 100644
--- a/libdimension/bvh.c
+++ b/libdimension/bvh/bvh.c
@@ -23,7 +23,10 @@
  * BVH implementation.  These are the hottest code paths in libdimension.
  */
 
-#include "dimension-internal.h"
+#include "internal/bvh.h"
+#include "internal/concurrency.h"
+#include "internal/prtree.h"
+#include <pthread.h>
 
 /// Implementation for DMNSN_BVH_NONE: just stick all objects in one node.
 static dmnsn_bvh_node *
@@ -48,9 +51,9 @@ struct dmnsn_bvh {
 
 /// A flat BVH node for storing in an array for fast pre-order traversal.
 typedef struct dmnsn_flat_bvh_node {
-  dmnsn_bounding_box bounding_box; // The bounding box of this node.
-  dmnsn_object *object;            // The referenced object, for leaf nodes.
-  ptrdiff_t skip;                  // Displacement to the next sibling.
+  dmnsn_aabb aabb;      ///< The bounding box of this node.
+  dmnsn_object *object; ///< The referenced object, for leaf nodes.
+  ptrdiff_t skip;       ///< Displacement to the next sibling.
 } dmnsn_flat_bvh_node;
 
 /// Add an object or its children, if any, to an array.
@@ -86,7 +89,7 @@ dmnsn_split_unbounded(dmnsn_array *objects)
   dmnsn_object **array = dmnsn_array_first(objects);
   size_t i, skip;
   for (i = 0, skip = 0; i < dmnsn_array_size(objects); ++i) {
-    if (dmnsn_bounding_box_is_infinite(array[i]->bounding_box)) {
+    if (dmnsn_aabb_is_infinite(array[i]->aabb)) {
       dmnsn_array_push(unbounded, &array[i]);
       ++skip;
     } else {
@@ -106,8 +109,8 @@ dmnsn_flatten_bvh_recursive(dmnsn_bvh_node *node, dmnsn_array *flat)
   dmnsn_array_resize(flat, currenti + 1);
   dmnsn_flat_bvh_node *flatnode = dmnsn_array_at(flat, currenti);
 
-  flatnode->bounding_box = node->bounding_box;
-  flatnode->object       = node->object;
+  flatnode->aabb = node->aabb;
+  flatnode->object = node->object;
 
   for (size_t i = 0; i < node->nchildren && node->children[i]; ++i) {
     dmnsn_flatten_bvh_recursive(node->children[i], flat);
@@ -171,49 +174,48 @@ dmnsn_delete_bvh(dmnsn_bvh *bvh)
   }
 }
 
-/// A line with pre-calculated reciprocals to avoid divisions.
-typedef struct dmnsn_optimized_line {
-  dmnsn_vector x0;    ///< The origin of the line.
-  dmnsn_vector n_inv; ///< The inverse of each component of the line's slope
-} dmnsn_optimized_line;
+/// A ray with pre-calculated reciprocals to avoid divisions.
+typedef struct dmnsn_optimized_ray {
+  dmnsn_vector x0;    ///< The origin of the ray.
+  dmnsn_vector n_inv; ///< The inverse of each component of the ray's slope
+} dmnsn_optimized_ray;
 
 /// Precompute inverses for faster ray-box intersection tests.
-static inline dmnsn_optimized_line
-dmnsn_optimize_line(dmnsn_line line)
+static inline dmnsn_optimized_ray
+dmnsn_optimize_ray(dmnsn_ray ray)
 {
-  dmnsn_optimized_line optline = {
-    .x0    = line.x0,
-    .n_inv = dmnsn_new_vector(1.0/line.n.x, 1.0/line.n.y, 1.0/line.n.z)
+  dmnsn_optimized_ray optray = {
+    .x0    = ray.x0,
+    .n_inv = dmnsn_new_vector(1.0/ray.n.x, 1.0/ray.n.y, 1.0/ray.n.z)
   };
-  return optline;
+  return optray;
 }
 
 /// Ray-AABB intersection test, by the slab method.  Highly optimized.
 static inline bool
-dmnsn_ray_box_intersection(dmnsn_optimized_line optline,
-                           dmnsn_bounding_box box, double t)
+dmnsn_ray_box_intersection(dmnsn_optimized_ray optray, dmnsn_aabb box, double t)
 {
   // This is actually correct, even though it appears not to handle edge cases
-  // (line.n.{x,y,z} == 0).  It works because the infinities that result from
-  // dividing by zero will still behave correctly in the comparisons.  Lines
+  // (ray.n.{x,y,z} == 0).  It works because the infinities that result from
+  // dividing by zero will still behave correctly in the comparisons.  Rays
   // which are parallel to an axis and outside the box will have tmin == inf
-  // or tmax == -inf, while lines inside the box will have tmin and tmax
+  // or tmax == -inf, while rays inside the box will have tmin and tmax
   // unchanged.
 
-  double tx1 = (box.min.x - optline.x0.x)*optline.n_inv.x;
-  double tx2 = (box.max.x - optline.x0.x)*optline.n_inv.x;
+  double tx1 = (box.min.x - optray.x0.x)*optray.n_inv.x;
+  double tx2 = (box.max.x - optray.x0.x)*optray.n_inv.x;
 
   double tmin = dmnsn_min(tx1, tx2);
   double tmax = dmnsn_max(tx1, tx2);
 
-  double ty1 = (box.min.y - optline.x0.y)*optline.n_inv.y;
-  double ty2 = (box.max.y - optline.x0.y)*optline.n_inv.y;
+  double ty1 = (box.min.y - optray.x0.y)*optray.n_inv.y;
+  double ty2 = (box.max.y - optray.x0.y)*optray.n_inv.y;
 
   tmin = dmnsn_max(tmin, dmnsn_min(ty1, ty2));
   tmax = dmnsn_min(tmax, dmnsn_max(ty1, ty2));
 
-  double tz1 = (box.min.z - optline.x0.z)*optline.n_inv.z;
-  double tz2 = (box.max.z - optline.x0.z)*optline.n_inv.z;
+  double tz1 = (box.min.z - optray.x0.z)*optray.n_inv.z;
+  double tz2 = (box.max.z - optray.x0.z)*optray.n_inv.z;
 
   tmin = dmnsn_max(tmin, dmnsn_min(tz1, tz2));
   tmax = dmnsn_min(tmax, dmnsn_max(tz1, tz2));
@@ -250,8 +252,7 @@ dmnsn_get_intersection_cache(const dmnsn_bvh *bvh)
 
 /// Test for a closer object intersection than we've found so far.
 static inline bool
-dmnsn_closer_intersection(dmnsn_object *object, dmnsn_line ray,
-                          dmnsn_intersection *intersection, double *t)
+dmnsn_closer_intersection(dmnsn_object *object, dmnsn_ray ray, dmnsn_intersection *intersection, double *t)
 {
   dmnsn_intersection local_intersection;
   if (dmnsn_object_intersection(object, ray, &local_intersection)) {
@@ -265,8 +266,7 @@ dmnsn_closer_intersection(dmnsn_object *object, dmnsn_line ray,
 }
 
 DMNSN_HOT bool
-dmnsn_bvh_intersection(const dmnsn_bvh *bvh, dmnsn_line ray,
-                       dmnsn_intersection *intersection, bool reset)
+dmnsn_bvh_intersection(const dmnsn_bvh *bvh, dmnsn_ray ray, dmnsn_intersection *intersection, bool reset)
 {
   double t = INFINITY;
 
@@ -276,7 +276,7 @@ dmnsn_bvh_intersection(const dmnsn_bvh *bvh, dmnsn_line ray,
   }
 
   // Precalculate 1.0/ray.n.{x,y,z} to save time in intersection tests
-  dmnsn_optimized_line optline = dmnsn_optimize_line(ray);
+  dmnsn_optimized_ray optray = dmnsn_optimize_ray(ray);
 
   // Search the intersection cache
   dmnsn_intersection_cache *cache = dmnsn_get_intersection_cache(bvh);
@@ -287,7 +287,7 @@ dmnsn_bvh_intersection(const dmnsn_bvh *bvh, dmnsn_line ray,
   if (dmnsn_likely(cache->i < DMNSN_INTERSECTION_CACHE_SIZE)) {
     cached = cache->objects[cache->i];
   }
-  if (cached && dmnsn_ray_box_intersection(optline, cached->bounding_box, t)) {
+  if (cached && dmnsn_ray_box_intersection(optray, cached->aabb, t)) {
     if (dmnsn_closer_intersection(cached, ray, intersection, &t)) {
       found = cached;
     }
@@ -297,7 +297,7 @@ dmnsn_bvh_intersection(const dmnsn_bvh *bvh, dmnsn_line ray,
   dmnsn_flat_bvh_node *node = dmnsn_array_first(bvh->bounded);
   dmnsn_flat_bvh_node *last = dmnsn_array_last(bvh->bounded);
   while (node <= last) {
-    if (dmnsn_ray_box_intersection(optline, node->bounding_box, t)) {
+    if (dmnsn_ray_box_intersection(optray, node->aabb, t)) {
       if (node->object && node->object != cached) {
         if (dmnsn_closer_intersection(node->object, ray, intersection, &t)) {
           found = node->object;
@@ -331,7 +331,7 @@ dmnsn_bvh_inside(const dmnsn_bvh *bvh, dmnsn_vector point)
   dmnsn_flat_bvh_node *node = dmnsn_array_first(bvh->bounded);
   dmnsn_flat_bvh_node *last = dmnsn_array_last(bvh->bounded);
   while (node <= last) {
-    if (dmnsn_bounding_box_contains(node->bounding_box, point)) {
+    if (dmnsn_aabb_contains(node->aabb, point)) {
       if (node->object && dmnsn_object_inside(node->object, point)) {
         return true;
       }
@@ -344,16 +344,16 @@ dmnsn_bvh_inside(const dmnsn_bvh *bvh, dmnsn_vector point)
   return false;
 }
 
-dmnsn_bounding_box
-dmnsn_bvh_bounding_box(const dmnsn_bvh *bvh)
+dmnsn_aabb
+dmnsn_bvh_aabb(const dmnsn_bvh *bvh)
 {
   if (dmnsn_array_size(bvh->unbounded) > 0) {
-    return dmnsn_infinite_bounding_box();
+    return dmnsn_infinite_aabb();
   } else if (dmnsn_array_size(bvh->bounded) > 0) {
     dmnsn_flat_bvh_node *root = dmnsn_array_first(bvh->bounded);
-    return root->bounding_box;
+    return root->aabb;
   } else {
-    return dmnsn_zero_bounding_box();
+    return dmnsn_zero_aabb();
   }
 }
 
@@ -361,7 +361,7 @@ dmnsn_bvh_node *
 dmnsn_new_bvh_node(unsigned int max_children)
 {
   dmnsn_bvh_node *node = dmnsn_malloc(sizeof(dmnsn_bvh_node) + max_children*sizeof(dmnsn_bvh_node *));
-  node->bounding_box = dmnsn_zero_bounding_box();
+  node->aabb = dmnsn_zero_aabb();
   node->object = NULL;
   node->nchildren = 0;
   node->max_children = max_children;
@@ -372,7 +372,7 @@ dmnsn_bvh_node *
 dmnsn_new_bvh_leaf_node(dmnsn_object *object)
 {
   dmnsn_bvh_node *node = DMNSN_MALLOC(dmnsn_bvh_node);
-  node->bounding_box = object->bounding_box;
+  node->aabb = object->aabb;
   node->object = object;
   node->nchildren = 0;
   node->max_children = 0;
@@ -396,9 +396,7 @@ dmnsn_bvh_node_add(dmnsn_bvh_node *parent, dmnsn_bvh_node *child)
   dmnsn_assert(parent->nchildren < parent->max_children,
                "Too many BVH children inserted.");
 
-  parent->bounding_box.min = dmnsn_vector_min(parent->bounding_box.min,
-                                              child->bounding_box.min);
-  parent->bounding_box.max = dmnsn_vector_max(parent->bounding_box.max,
-                                              child->bounding_box.max);
+  parent->aabb.min = dmnsn_vector_min(parent->aabb.min, child->aabb.min);
+  parent->aabb.max = dmnsn_vector_max(parent->aabb.max, child->aabb.max);
   parent->children[parent->nchildren++] = child;
 }
diff --git a/libdimension/prtree.c b/libdimension/bvh/prtree.c
index bb7beec..c8e4e54 100644
--- a/libdimension/prtree.c
+++ b/libdimension/bvh/prtree.c
@@ -23,7 +23,9 @@
  * Priority R-tree implementation.
  */
 
-#include "dimension-internal.h"
+#include "internal/platform.h"
+#include "internal/prtree.h"
+#include "internal/concurrency.h"
 #include <stdlib.h>
 
 /// Number of children per PR-node.
@@ -72,48 +74,48 @@ enum {
 static int
 dmnsn_xmin_comp(const void *l, const void *r)
 {
-  double lval = (*(const dmnsn_colored_prnode **)l)->node->bounding_box.min.x;
-  double rval = (*(const dmnsn_colored_prnode **)r)->node->bounding_box.min.x;
+  double lval = (*(const dmnsn_colored_prnode **)l)->node->aabb.min.x;
+  double rval = (*(const dmnsn_colored_prnode **)r)->node->aabb.min.x;
   return (lval > rval) - (lval < rval);
 }
 
 static int
 dmnsn_ymin_comp(const void *l, const void *r)
 {
-  double lval = (*(const dmnsn_colored_prnode **)l)->node->bounding_box.min.y;
-  double rval = (*(const dmnsn_colored_prnode **)r)->node->bounding_box.min.y;
+  double lval = (*(const dmnsn_colored_prnode **)l)->node->aabb.min.y;
+  double rval = (*(const dmnsn_colored_prnode **)r)->node->aabb.min.y;
   return (lval > rval) - (lval < rval);
 }
 
 static int
 dmnsn_zmin_comp(const void *l, const void *r)
 {
-  double lval = (*(const dmnsn_colored_prnode **)l)->node->bounding_box.min.z;
-  double rval = (*(const dmnsn_colored_prnode **)r)->node->bounding_box.min.z;
+  double lval = (*(const dmnsn_colored_prnode **)l)->node->aabb.min.z;
+  double rval = (*(const dmnsn_colored_prnode **)r)->node->aabb.min.z;
   return (lval > rval) - (lval < rval);
 }
 
 static int
 dmnsn_xmax_comp(const void *l, const void *r)
 {
-  double lval = (*(const dmnsn_colored_prnode **)l)->node->bounding_box.max.x;
-  double rval = (*(const dmnsn_colored_prnode **)r)->node->bounding_box.max.x;
+  double lval = (*(const dmnsn_colored_prnode **)l)->node->aabb.max.x;
+  double rval = (*(const dmnsn_colored_prnode **)r)->node->aabb.max.x;
   return (lval < rval) - (lval > rval);
 }
 
 static int
 dmnsn_ymax_comp(const void *l, const void *r)
 {
-  double lval = (*(const dmnsn_colored_prnode **)l)->node->bounding_box.max.y;
-  double rval = (*(const dmnsn_colored_prnode **)r)->node->bounding_box.max.y;
+  double lval = (*(const dmnsn_colored_prnode **)l)->node->aabb.max.y;
+  double rval = (*(const dmnsn_colored_prnode **)r)->node->aabb.max.y;
   return (lval < rval) - (lval > rval);
 }
 
 static int
 dmnsn_zmax_comp(const void *l, const void *r)
 {
-  double lval = (*(const dmnsn_colored_prnode **)l)->node->bounding_box.max.z;
-  double rval = (*(const dmnsn_colored_prnode **)r)->node->bounding_box.max.z;
+  double lval = (*(const dmnsn_colored_prnode **)l)->node->aabb.max.z;
+  double rval = (*(const dmnsn_colored_prnode **)r)->node->aabb.max.z;
   return (lval < rval) - (lval > rval);
 }
 
diff --git a/libdimension/canvas.c b/libdimension/canvas/canvas.c
index 2aade47..350cd60 100644
--- a/libdimension/canvas.c
+++ b/libdimension/canvas/canvas.c
@@ -23,7 +23,7 @@
  * Canveses.
  */
 
-#include "dimension-internal.h"
+#include "dimension/canvas.h"
 
 dmnsn_canvas *
 dmnsn_new_canvas(dmnsn_pool *pool, size_t width, size_t height)
diff --git a/libdimension/gl-stubs.c b/libdimension/canvas/gl-stubs.c
index df31308..df31308 100644
--- a/libdimension/gl-stubs.c
+++ b/libdimension/canvas/gl-stubs.c
diff --git a/libdimension/gl.c b/libdimension/canvas/gl.c
index 77099d2..900f53c 100644
--- a/libdimension/gl.c
+++ b/libdimension/canvas/gl.c
@@ -23,8 +23,10 @@
  * OpenGL import/export.
  */
 
-#include "dimension-internal.h"
+#include "internal/rgba.h"
+#include "dimension/canvas.h"
 #include <GL/gl.h>
+#include <math.h>
 #include <stdlib.h>
 #include <stdint.h>
 
diff --git a/libdimension/png-stubs.c b/libdimension/canvas/png-stubs.c
index 9c752a5..9c752a5 100644
--- a/libdimension/png-stubs.c
+++ b/libdimension/canvas/png-stubs.c
diff --git a/libdimension/png.c b/libdimension/canvas/png.c
index 6b73738..b1d29a3 100644
--- a/libdimension/png.c
+++ b/libdimension/canvas/png.c
@@ -23,7 +23,10 @@
  * PNG import/export.
  */
 
-#include "dimension-internal.h"
+#include "internal/concurrency.h"
+#include "internal/platform.h"
+#include "internal/rgba.h"
+#include "dimension/canvas.h"
 #include <png.h>
 #include <errno.h>
 #include <setjmp.h>
diff --git a/libdimension/rgba.c b/libdimension/canvas/rgba.c
index 9889189..4abf457 100644
--- a/libdimension/rgba.c
+++ b/libdimension/canvas/rgba.c
@@ -23,7 +23,7 @@
  * 16-bit RGBA canvas optimizer.
  */
 
-#include "dimension-internal.h"
+#include "internal/rgba.h"
 #include <stdint.h>
 
 void
diff --git a/libdimension/future.c b/libdimension/concurrency/future.c
index 8133812..90ffa24 100644
--- a/libdimension/future.c
+++ b/libdimension/concurrency/future.c
@@ -23,7 +23,9 @@
  * Future objects.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/concurrency.h"
+#include "internal/future.h"
 #include <pthread.h>
 
 /**
diff --git a/libdimension/threads.c b/libdimension/concurrency/threads.c
index 1fee01d..93d2ea9 100644
--- a/libdimension/threads.c
+++ b/libdimension/concurrency/threads.c
@@ -23,7 +23,9 @@
  * Background threading.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/concurrency.h"
+#include "internal/future.h"
 #include <pthread.h>
 
 /// The payload to pass to the pthread callback.
diff --git a/libdimension/dimension.h b/libdimension/dimension.h
index fc3a681..db67a48 100644
--- a/libdimension/dimension.h
+++ b/libdimension/dimension.h
@@ -22,6 +22,9 @@
  * @file
  * The main \#include file for the Dimension library.  This file declares all
  * of its public functions and types.
+ *
+ * To import only a subset of the libdimension's functionality, you can include
+ * the headers for submodules directly.
  */
 
 /**
@@ -31,55 +34,18 @@
  * all rendering-related tasks for Dimension.
  */
 
-#ifndef DIMENSION_H
-#define DIMENSION_H
-
-/* Include compiler.h first for DMNSN_CXX */
-#include <dimension/compiler.h>
+#ifndef DMNSN_H
+#define DMNSN_H
 
-#if DMNSN_CXX
-/* We've been included from a C++ file; mark everything here as extern "C" */
-extern "C" {
-#endif
-
-/* Include all the libdimension headers */
-#include <dimension/common.h>
-#include <dimension/error.h>
-#include <dimension/malloc.h>
-#include <dimension/pool.h>
-#include <dimension/array.h>
-#include <dimension/dictionary.h>
-#include <dimension/future.h>
-#include <dimension/timer.h>
+/* Include all modules. */
+#include <dimension/base.h>
+#include <dimension/platform.h>
+#include <dimension/concurrency.h>
 #include <dimension/math.h>
-#include <dimension/geometry.h>
-#include <dimension/polynomial.h>
 #include <dimension/color.h>
-#include <dimension/tcolor.h>
 #include <dimension/canvas.h>
-#include <dimension/gl.h>
-#include <dimension/png.h>
 #include <dimension/pattern.h>
-#include <dimension/patterns.h>
-#include <dimension/map.h>
-#include <dimension/pigment.h>
-#include <dimension/pigments.h>
-#include <dimension/finish.h>
-#include <dimension/finishes.h>
-#include <dimension/texture.h>
-#include <dimension/interior.h>
-#include <dimension/object.h>
-#include <dimension/objects.h>
-#include <dimension/csg.h>
-#include <dimension/light.h>
-#include <dimension/lights.h>
-#include <dimension/camera.h>
-#include <dimension/cameras.h>
-#include <dimension/scene.h>
-#include <dimension/ray_trace.h>
-
-#if DMNSN_CXX
-}
-#endif
+#include <dimension/model.h>
+#include <dimension/render.h>
 
-#endif /* DIMENSION_H */
+#endif /* DMNSN_H */
diff --git a/libdimension/dimension/base.h b/libdimension/dimension/base.h
new file mode 100644
index 0000000..d9f063f
--- /dev/null
+++ b/libdimension/dimension/base.h
@@ -0,0 +1,46 @@
+/*************************************************************************
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Basic functionality: compiler abstractions, error reporting, memory
+ * management, data structures, etc.
+ */
+
+#ifndef DMNSN_BASE_H
+#define DMNSN_BASE_H
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#include <dimension/base/compiler.h>
+#include <dimension/base/common.h>
+#include <dimension/base/error.h>
+#include <dimension/base/malloc.h>
+#include <dimension/base/pool.h>
+#include <dimension/base/array.h>
+#include <dimension/base/dictionary.h>
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* DMNSN_BASE_H */
diff --git a/libdimension/dimension/array.h b/libdimension/dimension/base/array.h
index 9261a0e..a54d0e5 100644
--- a/libdimension/dimension/array.h
+++ b/libdimension/dimension/base/array.h
@@ -23,6 +23,10 @@
  * Simple dynamic arrays.
  */
 
+#ifndef DMNSN_BASE_H
+#error "Please include <dimension/base.h> instead of this header directly."
+#endif
+
 #include <stddef.h> /* For size_t */
 #include <stdlib.h> /* For qsort() */
 #include <string.h> /* For memcpy() */
diff --git a/libdimension/dimension/common.h b/libdimension/dimension/base/common.h
index 15bafd8..b947b4a 100644
--- a/libdimension/dimension/common.h
+++ b/libdimension/dimension/base/common.h
@@ -23,6 +23,10 @@
  * Common types and utilities.
  */
 
+#ifndef DMNSN_BASE_H
+#error "Please include <dimension/base.h> instead of this header directly."
+#endif
+
 /**
  * Generic callback type.
  * @param[in,out] ptr  A pointer to an object to act on.
diff --git a/libdimension/dimension/compiler.h b/libdimension/dimension/base/compiler.h
index af3f4c8..a83f1b9 100644
--- a/libdimension/dimension/compiler.h
+++ b/libdimension/dimension/base/compiler.h
@@ -23,6 +23,10 @@
  * Compiler abstractions.
  */
 
+#ifndef DMNSN_BASE_H
+#error "Please include <dimension/base.h> instead of this header directly."
+#endif
+
 /**
  * @internal
  * @def DMNSN_C_VERSION
diff --git a/libdimension/dimension/dictionary.h b/libdimension/dimension/base/dictionary.h
index 887b171..7360b45 100644
--- a/libdimension/dimension/dictionary.h
+++ b/libdimension/dimension/base/dictionary.h
@@ -23,6 +23,10 @@
  * Simple associative arrays.
  */
 
+#ifndef DMNSN_BASE_H
+#error "Please include <dimension/base.h> instead of this header directly."
+#endif
+
 /** A string-object associative array. */
 typedef struct dmnsn_dictionary dmnsn_dictionary;
 
diff --git a/libdimension/dimension/error.h b/libdimension/dimension/base/error.h
index 0561b8a..d039081 100644
--- a/libdimension/dimension/error.h
+++ b/libdimension/dimension/base/error.h
@@ -20,12 +20,13 @@
 
 /**
  * @file
- * Error reporting interface.  Errors are reported at a given severity by the
- * dmnsn_error() macro at a given severity, which prints a warning if it is
- * below the set resilience, or prints an error and exits if it's at or above
- * the set resilience.
+ * Error reporting.
  */
 
+#ifndef DMNSN_BASE_H
+#error "Please include <dimension/base.h> instead of this header directly."
+#endif
+
 #include <stdbool.h>
 
 /**
diff --git a/libdimension/dimension/malloc.h b/libdimension/dimension/base/malloc.h
index 742e3a2..fcc492e 100644
--- a/libdimension/dimension/malloc.h
+++ b/libdimension/dimension/base/malloc.h
@@ -25,6 +25,10 @@
  * fails, they instead call dmnsn_error().
  */
 
+#ifndef DMNSN_BASE_H
+#error "Please include <dimension/base.h> instead of this header directly."
+#endif
+
 #include <stddef.h> /* For size_t */
 
 /**
diff --git a/libdimension/dimension/pool.h b/libdimension/dimension/base/pool.h
index 164bbbc..9983f1d 100644
--- a/libdimension/dimension/pool.h
+++ b/libdimension/dimension/base/pool.h
@@ -25,6 +25,10 @@
  * once once a scene is rendered (for example).
  */
 
+#ifndef DMNSN_BASE_H
+#error "Please include <dimension/base.h> instead of this header directly."
+#endif
+
 #include <stddef.h> /* For size_t */
 
 /* Forward-declare dmnsn_pool. */
diff --git a/libdimension/dimension/canvas.h b/libdimension/dimension/canvas.h
index e325364..53bbda3 100644
--- a/libdimension/dimension/canvas.h
+++ b/libdimension/dimension/canvas.h
@@ -1,5 +1,5 @@
 /*************************************************************************
- * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
  *                                                                       *
  * This file is part of The Dimension Library.                           *
  *                                                                       *
@@ -20,108 +20,26 @@
 
 /**
  * @file
- * A canvas which is rendered to.
+ * Using, importing, and exporting canvases.
  */
 
-#include <stddef.h>
+#ifndef DMNSN_CANVAS_H
+#define DMNSN_CANVAS_H
 
-/** A canvas, or image. */
-typedef struct dmnsn_canvas {
-  size_t width;  /**< Canvas width. */
-  size_t height; /**< Canvas height. */
+#ifdef __cplusplus
+extern "C" {
+#endif
 
-  /** An array of <tt>dmnsn_canvas_optimizer</tt>s. */
-  dmnsn_array *optimizers;
+#include <dimension/base.h>
+#include <dimension/concurrency.h>
+#include <dimension/color.h>
 
-  /**
-   * @internal
-   * Stored in first-quadrant representation (origin is bottom-left).  The pixel
-   * at (a,b) is accessible as pixels[b*width + a].
-   */
-  dmnsn_tcolor *pixels;
-} dmnsn_canvas;
+#include <dimension/canvas/canvas.h>
+#include <dimension/canvas/png.h>
+#include <dimension/canvas/gl.h>
 
-/* Forward-declare dmnsn_canvas_optimizer. */
-typedef struct dmnsn_canvas_optimizer dmnsn_canvas_optimizer;
-
-/**
- * Canvas optimizer callback type.
- * @param[in] optimizer  The canvas optimizer.
- * @param[in] canvas  The canvas that was just updated.
- * @param[in] x  The x-coordinate that was just updated.
- * @param[in] y  The y-coordinate that was just updated.
- */
-typedef void dmnsn_canvas_optimizer_fn(dmnsn_canvas_optimizer *optimizer, const dmnsn_canvas *canvas, size_t x, size_t y);
-
-/** Canvas optimizer. */
-struct dmnsn_canvas_optimizer {
-  dmnsn_canvas_optimizer_fn *optimizer_fn; /**< Optimizer callback. */
-};
-
-/**
- * Allocate a new canvas.
- * @param[in] pool  The memory pool to allocate from.
- * @param[in] width  The width of the canvas to allocate (in pixels).
- * @param[in] height  The height of the canvas to allocate (in pixels).
- * @return The allocated canvas.
- */
-dmnsn_canvas *dmnsn_new_canvas(dmnsn_pool *pool, size_t width, size_t height);
-
-/**
- * Initialize a dmnsn_canvas_optimizer field
- * @param[in] optimizer  The optimizer to initialize.
- */
-void dmnsn_init_canvas_optimizer(dmnsn_canvas_optimizer *optimizer);
-
-/**
- * Set a canvas optimizer
- * @param[in,out] canvas    The canvas to optimize.
- * @param[in]     optimizer The optimizer to use.
- */
-void dmnsn_canvas_optimize(dmnsn_canvas *canvas,
-                           const dmnsn_canvas_optimizer *optimizer);
-
-/**
- * Find a canvas optimizer by its callback.
- * @param[in] canvas        The canvas to check.
- * @param[in] optimizer_fn  The callback to search for for.
- * @return A pointer to the canvas optimizer with the callback \p optimizer_fn,
- *         or NULL if none is found.
- */
-dmnsn_canvas_optimizer *
-dmnsn_canvas_find_optimizer(const dmnsn_canvas *canvas,
-                            dmnsn_canvas_optimizer_fn *optimizer_fn);
-
-/* Pixel accessors */
-
-/**
- * Get the color of a pixel.
- * @param[in] canvas  The canvas to access.
- * @param[in] x       The x coordinate.
- * @param[in] y       The y coordinate.
- * @return The color of the pixel at (\p x, \p y).
- */
-DMNSN_INLINE dmnsn_tcolor
-dmnsn_canvas_get_pixel(const dmnsn_canvas *canvas, size_t x, size_t y)
-{
-  dmnsn_assert(x < canvas->width && y < canvas->height,
-               "Canvas access out of bounds.");
-  return canvas->pixels[y*canvas->width + x];
+#ifdef __cplusplus
 }
+#endif
 
-/**
- * Set the value of a pixel.
- * @param[in,out] canvas  The canvas to modify.
- * @param[in]     x       The x coordinate of the pixel.
- * @param[in]     y       The y coordinate of the pixel.
- * @param[in]     tcolor  The value to set the pixel at (\p x, \p y) to.
- */
-void dmnsn_canvas_set_pixel(dmnsn_canvas *canvas, size_t x, size_t y,
-                            dmnsn_tcolor tcolor);
-
-/**
- * Clear a canvas uniformly with a given color.
- * @param[in,out] canvas  The canvas to erase.
- * @param[in]     tcolor  The color to paint it with.
- */
-void dmnsn_canvas_clear(dmnsn_canvas *canvas, dmnsn_tcolor tcolor);
+#endif /* DMNSN_CANVAS_H */
diff --git a/libdimension/dimension/canvas/canvas.h b/libdimension/dimension/canvas/canvas.h
new file mode 100644
index 0000000..07d4fde
--- /dev/null
+++ b/libdimension/dimension/canvas/canvas.h
@@ -0,0 +1,131 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * A canvas which is rendered to.
+ */
+
+#ifndef DMNSN_CANVAS_H
+#error "Please include <dimension/canvas.h> instead of this header directly."
+#endif
+
+#include <stddef.h>
+
+/** A canvas, or image. */
+typedef struct dmnsn_canvas {
+  size_t width;  /**< Canvas width. */
+  size_t height; /**< Canvas height. */
+
+  /** An array of <tt>dmnsn_canvas_optimizer</tt>s. */
+  dmnsn_array *optimizers;
+
+  /**
+   * @internal
+   * Stored in first-quadrant representation (origin is bottom-left).  The pixel
+   * at (a,b) is accessible as pixels[b*width + a].
+   */
+  dmnsn_tcolor *pixels;
+} dmnsn_canvas;
+
+/* Forward-declare dmnsn_canvas_optimizer. */
+typedef struct dmnsn_canvas_optimizer dmnsn_canvas_optimizer;
+
+/**
+ * Canvas optimizer callback type.
+ * @param[in] optimizer  The canvas optimizer.
+ * @param[in] canvas  The canvas that was just updated.
+ * @param[in] x  The x-coordinate that was just updated.
+ * @param[in] y  The y-coordinate that was just updated.
+ */
+typedef void dmnsn_canvas_optimizer_fn(dmnsn_canvas_optimizer *optimizer, const dmnsn_canvas *canvas, size_t x, size_t y);
+
+/** Canvas optimizer. */
+struct dmnsn_canvas_optimizer {
+  dmnsn_canvas_optimizer_fn *optimizer_fn; /**< Optimizer callback. */
+};
+
+/**
+ * Allocate a new canvas.
+ * @param[in] pool  The memory pool to allocate from.
+ * @param[in] width  The width of the canvas to allocate (in pixels).
+ * @param[in] height  The height of the canvas to allocate (in pixels).
+ * @return The allocated canvas.
+ */
+dmnsn_canvas *dmnsn_new_canvas(dmnsn_pool *pool, size_t width, size_t height);
+
+/**
+ * Initialize a dmnsn_canvas_optimizer field
+ * @param[in] optimizer  The optimizer to initialize.
+ */
+void dmnsn_init_canvas_optimizer(dmnsn_canvas_optimizer *optimizer);
+
+/**
+ * Set a canvas optimizer
+ * @param[in,out] canvas    The canvas to optimize.
+ * @param[in]     optimizer The optimizer to use.
+ */
+void dmnsn_canvas_optimize(dmnsn_canvas *canvas,
+                           const dmnsn_canvas_optimizer *optimizer);
+
+/**
+ * Find a canvas optimizer by its callback.
+ * @param[in] canvas        The canvas to check.
+ * @param[in] optimizer_fn  The callback to search for for.
+ * @return A pointer to the canvas optimizer with the callback \p optimizer_fn,
+ *         or NULL if none is found.
+ */
+dmnsn_canvas_optimizer *
+dmnsn_canvas_find_optimizer(const dmnsn_canvas *canvas,
+                            dmnsn_canvas_optimizer_fn *optimizer_fn);
+
+/* Pixel accessors */
+
+/**
+ * Get the color of a pixel.
+ * @param[in] canvas  The canvas to access.
+ * @param[in] x       The x coordinate.
+ * @param[in] y       The y coordinate.
+ * @return The color of the pixel at (\p x, \p y).
+ */
+DMNSN_INLINE dmnsn_tcolor
+dmnsn_canvas_get_pixel(const dmnsn_canvas *canvas, size_t x, size_t y)
+{
+  dmnsn_assert(x < canvas->width && y < canvas->height,
+               "Canvas access out of bounds.");
+  return canvas->pixels[y*canvas->width + x];
+}
+
+/**
+ * Set the value of a pixel.
+ * @param[in,out] canvas  The canvas to modify.
+ * @param[in]     x       The x coordinate of the pixel.
+ * @param[in]     y       The y coordinate of the pixel.
+ * @param[in]     tcolor  The value to set the pixel at (\p x, \p y) to.
+ */
+void dmnsn_canvas_set_pixel(dmnsn_canvas *canvas, size_t x, size_t y,
+                            dmnsn_tcolor tcolor);
+
+/**
+ * Clear a canvas uniformly with a given color.
+ * @param[in,out] canvas  The canvas to erase.
+ * @param[in]     tcolor  The color to paint it with.
+ */
+void dmnsn_canvas_clear(dmnsn_canvas *canvas, dmnsn_tcolor tcolor);
diff --git a/libdimension/dimension/gl.h b/libdimension/dimension/canvas/gl.h
index 284a248..c2d77a1 100644
--- a/libdimension/dimension/gl.h
+++ b/libdimension/dimension/canvas/gl.h
@@ -23,6 +23,10 @@
  * OpenGL export/import of canvases.
  */
 
+#ifndef DMNSN_CANVAS_H
+#error "Please include <dimension/canvas.h> instead of this header directly."
+#endif
+
 /**
  * Optimize a canvas for GL drawing
  * @param[in] pool  The memory pool to allocate from.
diff --git a/libdimension/dimension/png.h b/libdimension/dimension/canvas/png.h
index f4d2c1e..57a68c2 100644
--- a/libdimension/dimension/png.h
+++ b/libdimension/dimension/canvas/png.h
@@ -23,6 +23,10 @@
  * PNG import/export of canvases.
  */
 
+#ifndef DMNSN_CANVAS_H
+#error "Please include <dimension/canvas.h> instead of this header directly."
+#endif
+
 #include <stdio.h>
 
 /**
diff --git a/libdimension/dimension/color.h b/libdimension/dimension/color.h
index adbd396..74bed59 100644
--- a/libdimension/dimension/color.h
+++ b/libdimension/dimension/color.h
@@ -1,5 +1,5 @@
 /*************************************************************************
- * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
  *                                                                       *
  * This file is part of The Dimension Library.                           *
  *                                                                       *
@@ -20,174 +20,23 @@
 
 /**
  * @file
- * Colors.
+ * Color handling.
  */
 
-#include <stdbool.h>
+#ifndef DMNSN_COLOR_H
+#define DMNSN_COLOR_H
 
-/** A color value. */
-typedef struct {
-  double R; /**< Red component. */
-  double G; /**< Green component. */
-  double B; /**< Blue component. */
-} dmnsn_color;
+#ifdef __cplusplus
+extern "C" {
+#endif
 
-/** A standard format string for colors. */
-#define DMNSN_COLOR_FORMAT "Color<%g, %g, %g>"
-/** The appropriate arguements to printf() a color. */
-#define DMNSN_COLOR_PRINTF(c) (c).R, (c).G, (c).B
+#include <dimension/math.h>
 
-/** Construct a new color. */
-DMNSN_INLINE dmnsn_color
-dmnsn_new_color(double R, double G, double B)
-{
-  dmnsn_color ret = { R, G, B };
-  return ret;
-}
-
-/** Apply sRGB gamma */
-DMNSN_INLINE double
-dmnsn_sRGB_gamma(double Clinear)
-{
-  /*
-   * If C represents R, G, and B, then the sRGB values are now found as follows:
-   *
-   *         { 12.92*Clinear,                Clinear <= 0.0031308
-   * Csrgb = {                1/2.4
-   *         { (1.055)*Clinear      - 0.055, Clinear >  0.0031308
-   */
-  if (Clinear == 1.0) {
-    return 1.0; /* Map 1.0 to 1.0 instead of 0.9999999999999999 */
-  } else if (Clinear > 0.0031308) {
-    return 1.055*pow(Clinear, 1.0/2.4) - 0.055;
-  } else {
-    return 12.92*Clinear;
-  }
-}
-
-/** Convert to sRGB space. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_to_sRGB(dmnsn_color color)
-{
-  return dmnsn_new_color(
-    dmnsn_sRGB_gamma(color.R),
-    dmnsn_sRGB_gamma(color.G),
-    dmnsn_sRGB_gamma(color.B)
-  );
-}
-
-/** Remove sRGB gamma */
-DMNSN_INLINE double
-dmnsn_sRGB_inverse_gamma(double CsRGB)
-{
-  /*
-   * If C represents R, G, and B, then the Clinear values are now found as
-   * follows:
-   *
-   *           { Csrgb/12.92,                Csrgb <= 0.04045
-   * Clinear = {                        2.4
-   *           { ((Csrgb + 0.055)/1.055)   , Csrgb >  0.04045
-   */
-  if (CsRGB == 1.0) {
-    return 1.0; /* Map 1.0 to 1.0 instead of 0.9999999999999999 */
-  } else if (CsRGB <= 0.040449936) {
-    return CsRGB/12.92;
-  } else {
-    return pow((CsRGB + 0.055)/1.055, 2.4);
-  }
-}
-
-/** Convert from sRGB space. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_from_sRGB(dmnsn_color color)
-{
-  return dmnsn_new_color(
-    dmnsn_sRGB_inverse_gamma(color.R),
-    dmnsn_sRGB_inverse_gamma(color.G),
-    dmnsn_sRGB_inverse_gamma(color.B)
-  );
-}
-
-/** Greyscale color intensity. */
-DMNSN_INLINE double
-dmnsn_color_intensity(dmnsn_color color)
-{
-  return 0.2126*color.R + 0.7152*color.G + 0.0722*color.B;
-}
-
-/** Add two colors together. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_add(dmnsn_color lhs, dmnsn_color rhs)
-{
-  return dmnsn_new_color(lhs.R + rhs.R, lhs.G + rhs.G, lhs.B + rhs.B);
-}
-
-/** Subtract two colors. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_sub(dmnsn_color lhs, dmnsn_color rhs)
-{
-  return dmnsn_new_color(lhs.R - rhs.R, lhs.G - rhs.G, lhs.B - rhs.B);
-}
-
-/** Scale a color's intensity. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_mul(double n, dmnsn_color color)
-{
-  return dmnsn_new_color(n*color.R, n*color.G, n*color.B);
-}
-
-/** Return the color at \p n on a gradient from \p c1 at 0 to \p c2 at 1. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_gradient(dmnsn_color c1, dmnsn_color c2, double n)
-{
-  return dmnsn_new_color(
-    n*(c2.R - c1.R) + c1.R,
-    n*(c2.G - c1.G) + c1.G,
-    n*(c2.B - c1.B) + c1.B
-  );
-}
-
-/** Illuminate \p color with \p light. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_illuminate(dmnsn_color light, dmnsn_color color)
-{
-  return dmnsn_new_color(light.R*color.R, light.G*color.G, light.B*color.B);
-}
+#include <dimension/color/color.h>
+#include <dimension/color/tcolor.h>
 
-/** Saturate the color components to [0.0, 1.0]. */
-DMNSN_INLINE dmnsn_color
-dmnsn_color_clamp(dmnsn_color color)
-{
-  color.R = dmnsn_clamp(color.R, 0.0, 1.0);
-  color.G = dmnsn_clamp(color.G, 0.0, 1.0);
-  color.B = dmnsn_clamp(color.B, 0.0, 1.0);
-  return color;
+#ifdef __cplusplus
 }
+#endif
 
-/** Return whether a color contains any NaN components. */
-DMNSN_INLINE bool
-dmnsn_color_isnan(dmnsn_color color)
-{
-  return dmnsn_isnan(color.R) || dmnsn_isnan(color.G) || dmnsn_isnan(color.B);
-}
-
-/* Standard colors */
-
-/** Black. */
-#define dmnsn_black   dmnsn_new_color(0.0, 0.0, 0.0)
-/** White. */
-#define dmnsn_white   dmnsn_new_color(1.0, 1.0, 1.0)
-/** Red. */
-#define dmnsn_red     dmnsn_new_color(1.0, 0.0, 0.0)
-/** Green. */
-#define dmnsn_green   dmnsn_new_color(0.0, 1.0, 0.0)
-/** Blue. */
-#define dmnsn_blue    dmnsn_new_color(0.0, 0.0, 1.0)
-/** Magenta. */
-#define dmnsn_magenta dmnsn_new_color(1.0, 0.0, 1.0)
-/** Orange. */
-#define dmnsn_orange  dmnsn_new_color(1.0, 0.21404114048223255, 0.0)
-/** Yellow. */
-#define dmnsn_yellow  dmnsn_new_color(1.0, 1.0, 0.0)
-/** Cyan. */
-#define dmnsn_cyan    dmnsn_new_color(0.0, 1.0, 1.0)
+#endif /* DMNSN_COLOR_H */
diff --git a/libdimension/dimension/color/color.h b/libdimension/dimension/color/color.h
new file mode 100644
index 0000000..84e66ea
--- /dev/null
+++ b/libdimension/dimension/color/color.h
@@ -0,0 +1,198 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Colors.
+ */
+
+#ifndef DMNSN_COLOR_H
+#error "Please include <dimension/color.h> instead of this header directly."
+#endif
+
+#include <math.h>
+#include <stdbool.h>
+
+/** A color value. */
+typedef struct {
+  double R; /**< Red component. */
+  double G; /**< Green component. */
+  double B; /**< Blue component. */
+} dmnsn_color;
+
+/** A standard format string for colors. */
+#define DMNSN_COLOR_FORMAT "Color<%g, %g, %g>"
+/** The appropriate arguements to printf() a color. */
+#define DMNSN_COLOR_PRINTF(c) (c).R, (c).G, (c).B
+
+/** Construct a new color. */
+DMNSN_INLINE dmnsn_color
+dmnsn_new_color(double R, double G, double B)
+{
+  dmnsn_color ret = { R, G, B };
+  return ret;
+}
+
+/** Apply sRGB gamma */
+DMNSN_INLINE double
+dmnsn_sRGB_gamma(double Clinear)
+{
+  /*
+   * If C represents R, G, and B, then the sRGB values are now found as follows:
+   *
+   *         { 12.92*Clinear,                Clinear <= 0.0031308
+   * Csrgb = {                1/2.4
+   *         { (1.055)*Clinear      - 0.055, Clinear >  0.0031308
+   */
+  if (Clinear == 1.0) {
+    return 1.0; /* Map 1.0 to 1.0 instead of 0.9999999999999999 */
+  } else if (Clinear > 0.0031308) {
+    return 1.055*pow(Clinear, 1.0/2.4) - 0.055;
+  } else {
+    return 12.92*Clinear;
+  }
+}
+
+/** Convert to sRGB space. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_to_sRGB(dmnsn_color color)
+{
+  return dmnsn_new_color(
+    dmnsn_sRGB_gamma(color.R),
+    dmnsn_sRGB_gamma(color.G),
+    dmnsn_sRGB_gamma(color.B)
+  );
+}
+
+/** Remove sRGB gamma */
+DMNSN_INLINE double
+dmnsn_sRGB_inverse_gamma(double CsRGB)
+{
+  /*
+   * If C represents R, G, and B, then the Clinear values are now found as
+   * follows:
+   *
+   *           { Csrgb/12.92,                Csrgb <= 0.04045
+   * Clinear = {                        2.4
+   *           { ((Csrgb + 0.055)/1.055)   , Csrgb >  0.04045
+   */
+  if (CsRGB == 1.0) {
+    return 1.0; /* Map 1.0 to 1.0 instead of 0.9999999999999999 */
+  } else if (CsRGB <= 0.040449936) {
+    return CsRGB/12.92;
+  } else {
+    return pow((CsRGB + 0.055)/1.055, 2.4);
+  }
+}
+
+/** Convert from sRGB space. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_from_sRGB(dmnsn_color color)
+{
+  return dmnsn_new_color(
+    dmnsn_sRGB_inverse_gamma(color.R),
+    dmnsn_sRGB_inverse_gamma(color.G),
+    dmnsn_sRGB_inverse_gamma(color.B)
+  );
+}
+
+/** Greyscale color intensity. */
+DMNSN_INLINE double
+dmnsn_color_intensity(dmnsn_color color)
+{
+  return 0.2126*color.R + 0.7152*color.G + 0.0722*color.B;
+}
+
+/** Add two colors together. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_add(dmnsn_color lhs, dmnsn_color rhs)
+{
+  return dmnsn_new_color(lhs.R + rhs.R, lhs.G + rhs.G, lhs.B + rhs.B);
+}
+
+/** Subtract two colors. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_sub(dmnsn_color lhs, dmnsn_color rhs)
+{
+  return dmnsn_new_color(lhs.R - rhs.R, lhs.G - rhs.G, lhs.B - rhs.B);
+}
+
+/** Scale a color's intensity. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_mul(double n, dmnsn_color color)
+{
+  return dmnsn_new_color(n*color.R, n*color.G, n*color.B);
+}
+
+/** Return the color at \p n on a gradient from \p c1 at 0 to \p c2 at 1. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_gradient(dmnsn_color c1, dmnsn_color c2, double n)
+{
+  return dmnsn_new_color(
+    n*(c2.R - c1.R) + c1.R,
+    n*(c2.G - c1.G) + c1.G,
+    n*(c2.B - c1.B) + c1.B
+  );
+}
+
+/** Illuminate \p color with \p light. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_illuminate(dmnsn_color light, dmnsn_color color)
+{
+  return dmnsn_new_color(light.R*color.R, light.G*color.G, light.B*color.B);
+}
+
+/** Saturate the color components to [0.0, 1.0]. */
+DMNSN_INLINE dmnsn_color
+dmnsn_color_clamp(dmnsn_color color)
+{
+  color.R = dmnsn_clamp(color.R, 0.0, 1.0);
+  color.G = dmnsn_clamp(color.G, 0.0, 1.0);
+  color.B = dmnsn_clamp(color.B, 0.0, 1.0);
+  return color;
+}
+
+/** Return whether a color contains any NaN components. */
+DMNSN_INLINE bool
+dmnsn_color_isnan(dmnsn_color color)
+{
+  return dmnsn_isnan(color.R) || dmnsn_isnan(color.G) || dmnsn_isnan(color.B);
+}
+
+/* Standard colors */
+
+/** Black. */
+#define dmnsn_black   dmnsn_new_color(0.0, 0.0, 0.0)
+/** White. */
+#define dmnsn_white   dmnsn_new_color(1.0, 1.0, 1.0)
+/** Red. */
+#define dmnsn_red     dmnsn_new_color(1.0, 0.0, 0.0)
+/** Green. */
+#define dmnsn_green   dmnsn_new_color(0.0, 1.0, 0.0)
+/** Blue. */
+#define dmnsn_blue    dmnsn_new_color(0.0, 0.0, 1.0)
+/** Magenta. */
+#define dmnsn_magenta dmnsn_new_color(1.0, 0.0, 1.0)
+/** Orange. */
+#define dmnsn_orange  dmnsn_new_color(1.0, 0.21404114048223255, 0.0)
+/** Yellow. */
+#define dmnsn_yellow  dmnsn_new_color(1.0, 1.0, 0.0)
+/** Cyan. */
+#define dmnsn_cyan    dmnsn_new_color(0.0, 1.0, 1.0)
diff --git a/libdimension/dimension/tcolor.h b/libdimension/dimension/color/tcolor.h
index 8c48800..b4b4167 100644
--- a/libdimension/dimension/tcolor.h
+++ b/libdimension/dimension/color/tcolor.h
@@ -23,6 +23,10 @@
  * Colors with transparency information.
  */
 
+#ifndef DMNSN_COLOR_H
+#error "Please include <dimension/color.h> instead of this header directly."
+#endif
+
 /** A transparent color. */
 typedef struct dmnsn_tcolor {
   dmnsn_color c; /**< Color. */
diff --git a/libdimension/dimension/concurrency.h b/libdimension/dimension/concurrency.h
new file mode 100644
index 0000000..6c8c29c
--- /dev/null
+++ b/libdimension/dimension/concurrency.h
@@ -0,0 +1,39 @@
+/*************************************************************************
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Concurrency concerns.
+ */
+
+#ifndef DMNSN_CONCURRENCY_H
+#define DMNSN_CONCURRENCY_H
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#include <dimension/concurrency/future.h>
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* DMNSN_CONCURRENCY_H */
diff --git a/libdimension/dimension/future.h b/libdimension/dimension/concurrency/future.h
index 9ba28b1..24d5ee2 100644
--- a/libdimension/dimension/future.h
+++ b/libdimension/dimension/concurrency/future.h
@@ -26,6 +26,12 @@
  * is returned as an int from dmnsn_finish_progress().
  */
 
+#ifndef DMNSN_CONCURRENCY_H
+#error "Please include <dimension/concurrency.h> instead of this header directly."
+#endif
+
+#include <stdbool.h>
+
 /** A future object. */
 typedef struct dmnsn_future dmnsn_future;
 
diff --git a/libdimension/dimension/geometry.h b/libdimension/dimension/geometry.h
deleted file mode 100644
index 2ea10ca..0000000
--- a/libdimension/dimension/geometry.h
+++ /dev/null
@@ -1,500 +0,0 @@
-/*************************************************************************
- * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
- *                                                                       *
- * This file is part of The Dimension Library.                           *
- *                                                                       *
- * The Dimension Library is free software; you can redistribute it and/  *
- * or modify it under the terms of the GNU Lesser General Public License *
- * as published by the Free Software Foundation; either version 3 of the *
- * License, or (at your option) any later version.                       *
- *                                                                       *
- * The Dimension Library is distributed in the hope that it will be      *
- * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
- * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
- * Lesser General Public License for more details.                       *
- *                                                                       *
- * You should have received a copy of the GNU Lesser General Public      *
- * License along with this program.  If not, see                         *
- * <http://www.gnu.org/licenses/>.                                       *
- *************************************************************************/
-
-/**
- * @file
- * Core geometric types like vectors, matricies, and rays.
- */
-
-#include <math.h>
-#include <stdbool.h>
-
-/** A vector in 3 dimensions. */
-typedef struct dmnsn_vector {
-  double x; /**< The x component. */
-  double y; /**< The y component. */
-  double z; /**< The z component. */
-} dmnsn_vector;
-
-/** A standard format string for vectors. */
-#define DMNSN_VECTOR_FORMAT "<%g, %g, %g>"
-/** The appropriate arguements to printf() a vector. */
-#define DMNSN_VECTOR_PRINTF(v) (v).x, (v).y, (v).z
-
-/** A 4x4 affine transformation matrix, with implied [0 0 0 1] bottom row. */
-typedef struct dmnsn_matrix {
-  double n[3][4]; /**< The matrix elements in row-major order. */
-} dmnsn_matrix;
-
-/** A standard format string for matricies. */
-#define DMNSN_MATRIX_FORMAT                     \
-  "[%g\t%g\t%g\t%g]\n"                          \
-  "[%g\t%g\t%g\t%g]\n"                          \
-  "[%g\t%g\t%g\t%g]\n"                          \
-  "[%g\t%g\t%g\t%g]"
-/** The appropriate arguements to printf() a matrix. */
-#define DMNSN_MATRIX_PRINTF(m)                          \
-  (m).n[0][0], (m).n[0][1], (m).n[0][2], (m).n[0][3],   \
-  (m).n[1][0], (m).n[1][1], (m).n[1][2], (m).n[1][3],   \
-  (m).n[2][0], (m).n[2][1], (m).n[2][2], (m).n[2][3],   \
-  0.0, 0.0, 0.0, 1.0
-
-/** A line, or ray. */
-typedef struct dmnsn_line {
-  dmnsn_vector x0; /**< A point on the line. */
-  dmnsn_vector n;  /**< A normal vector; the direction of the line. */
-} dmnsn_line;
-
-/** A standard format string for lines. */
-#define DMNSN_LINE_FORMAT "(<%g, %g, %g> + t*<%g, %g, %g>)"
-/** The appropriate arguements to printf() a line. */
-#define DMNSN_LINE_PRINTF(l)                                    \
-  DMNSN_VECTOR_PRINTF((l).x0), DMNSN_VECTOR_PRINTF((l).n)
-
-/** An axis-aligned bounding box (AABB). */
-typedef struct dmnsn_bounding_box {
-  dmnsn_vector min; /**< The coordinate-wise minimum extent of the box. */
-  dmnsn_vector max; /**< The coordinate-wise maximum extent of the box. */
-} dmnsn_bounding_box;
-
-/** A standard format string for bounding boxes. */
-#define DMNSN_BOUNDING_BOX_FORMAT "(<%g, %g, %g> ==> <%g, %g, %g>)"
-/** The appropriate arguements to printf() a bounding box. */
-#define DMNSN_BOUNDING_BOX_PRINTF(box)                                  \
-  DMNSN_VECTOR_PRINTF((box).min), DMNSN_VECTOR_PRINTF((box).max)
-
-/* Constants */
-
-/** The zero vector. */
-static const dmnsn_vector dmnsn_zero = { 0.0, 0.0, 0.0 };
-/** The x vector. */
-static const dmnsn_vector dmnsn_x = { 1.0, 0.0, 0.0 };
-/** The y vector. */
-static const dmnsn_vector dmnsn_y = { 0.0, 1.0, 0.0 };
-/** The z vector. */
-static const dmnsn_vector dmnsn_z = { 0.0, 0.0, 1.0 };
-
-/* Shorthand for vector/matrix construction */
-
-/** Construct a new vector. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_new_vector(double x, double y, double z)
-{
-  dmnsn_vector v = { x, y, z };
-  return v;
-}
-
-/** Construct a new transformation matrix. */
-DMNSN_INLINE dmnsn_matrix
-dmnsn_new_matrix(double a0, double a1, double a2, double a3,
-                 double b0, double b1, double b2, double b3,
-                 double c0, double c1, double c2, double c3)
-{
-  dmnsn_matrix m = { { { a0, a1, a2, a3 },
-                       { b0, b1, b2, b3 },
-                       { c0, c1, c2, c3 } } };
-  return m;
-}
-
-/** Construct a new transformation matrix from column vectors. */
-DMNSN_INLINE dmnsn_matrix
-dmnsn_new_matrix4(dmnsn_vector a, dmnsn_vector b, dmnsn_vector c,
-                  dmnsn_vector d)
-{
-  dmnsn_matrix m = { { { a.x, b.x, c.x, d.x },
-                       { a.y, b.y, c.y, d.y },
-                       { a.z, b.z, c.z, d.z } } };
-  return m;
-}
-
-/** Extract column vectors from a matrix. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_matrix_column(dmnsn_matrix M, unsigned int i)
-{
-  return dmnsn_new_vector(M.n[0][i], M.n[1][i], M.n[2][i]);
-}
-
-/** Return the identity matrix. */
-dmnsn_matrix dmnsn_identity_matrix(void);
-
-/**
- * A scale transformation.
- * @param[in] s  A vector with components representing the scaling factor in
- *               each axis.
- * @return The transformation matrix.
- */
-dmnsn_matrix dmnsn_scale_matrix(dmnsn_vector s);
-/**
- * A translation.
- * @param[in] d  The vector to translate by.
- * @return The transformation matrix.
- */
-dmnsn_matrix dmnsn_translation_matrix(dmnsn_vector d);
-/**
- * A left-handed rotation.
- * @param[in] theta  A vector representing an axis and angle.
- *                   @f$ axis = \vec{\theta}/|\vec{\theta}| @f$,
- *                   @f$ angle = |\vec{\theta}| @f$
- * @return The transformation matrix.
- */
-dmnsn_matrix dmnsn_rotation_matrix(dmnsn_vector theta);
-/**
- * An alignment matrix.
- * @param[in] from   The initial vector.
- * @param[in] to     The desired direction.
- * @param[in] axis1  The first axis about which to rotate.
- * @param[in] axis2  The second axis about which to rotate.
- * @return A transformation matrix that will rotate \p from to \p to.
- */
-dmnsn_matrix dmnsn_alignment_matrix(dmnsn_vector from, dmnsn_vector to,
-                                    dmnsn_vector axis1, dmnsn_vector axis2);
-
-/**
- * Construct a new line.
- * @param[in] x0  A point on the line.
- * @param[in] n   The direction of the line.
- * @return The new line.
- */
-DMNSN_INLINE dmnsn_line
-dmnsn_new_line(dmnsn_vector x0, dmnsn_vector n)
-{
-  dmnsn_line l = { x0, n };
-  return l;
-}
-
-/**
- * Construct a new bounding box.
- * @param[in] min  The minimal extent of the bounding box.
- * @param[in] max  The maximal extent of the bounding box.
- * @return The new bounding box.
- */
-DMNSN_INLINE dmnsn_bounding_box
-dmnsn_new_bounding_box(dmnsn_vector min, dmnsn_vector max)
-{
-  dmnsn_bounding_box box = { min, max };
-  return box;
-}
-
-/** Return the bounding box which contains nothing. */
-DMNSN_INLINE dmnsn_bounding_box
-dmnsn_zero_bounding_box(void)
-{
-  dmnsn_bounding_box box = {
-    {  DMNSN_INFINITY,  DMNSN_INFINITY,  DMNSN_INFINITY },
-    { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY }
-  };
-  return box;
-}
-
-/** Return the bounding box which contains everything. */
-DMNSN_INLINE dmnsn_bounding_box
-dmnsn_infinite_bounding_box(void)
-{
-  dmnsn_bounding_box box = {
-    { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY },
-    {  DMNSN_INFINITY,  DMNSN_INFINITY,  DMNSN_INFINITY }
-  };
-  return box;
-}
-
-/* Vector and matrix arithmetic */
-
-/** Negate a vector. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_negate(dmnsn_vector rhs)
-{
-  /* 3 negations */
-  dmnsn_vector v = { -rhs.x, -rhs.y, -rhs.z };
-  return v;
-}
-
-/** Add two vectors. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_add(dmnsn_vector lhs, dmnsn_vector rhs)
-{
-  /* 3 additions */
-  dmnsn_vector v = { lhs.x + rhs.x, lhs.y + rhs.y, lhs.z + rhs.z };
-  return v;
-}
-
-/** Subtract two vectors. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_sub(dmnsn_vector lhs, dmnsn_vector rhs)
-{
-  /* 3 additions */
-  dmnsn_vector v = { lhs.x - rhs.x, lhs.y - rhs.y, lhs.z - rhs.z };
-  return v;
-}
-
-/** Multiply a vector by a scalar. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_mul(double lhs, dmnsn_vector rhs)
-{
-  /* 3 multiplications */
-  dmnsn_vector v = { lhs*rhs.x, lhs*rhs.y, lhs*rhs.z };
-  return v;
-}
-
-/** Divide a vector by a scalar. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_div(dmnsn_vector lhs, double rhs)
-{
-  /* 3 divisions */
-  dmnsn_vector v = { lhs.x/rhs, lhs.y/rhs, lhs.z/rhs };
-  return v;
-}
-
-/** Return the dot product of two vectors. */
-DMNSN_INLINE double
-dmnsn_vector_dot(dmnsn_vector lhs, dmnsn_vector rhs)
-{
-  /* 3 multiplications, 2 additions */
-  return lhs.x*rhs.x + lhs.y*rhs.y + lhs.z*rhs.z;
-}
-
-/** Return the cross product of two vectors. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_cross(dmnsn_vector lhs, dmnsn_vector rhs)
-{
-  /* 6 multiplications, 3 additions */
-  dmnsn_vector v = { lhs.y*rhs.z - lhs.z*rhs.y,
-                     lhs.z*rhs.x - lhs.x*rhs.z,
-                     lhs.x*rhs.y - lhs.y*rhs.x };
-  return v;
-}
-
-/** Return the projection of \p u onto \p d. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_proj(dmnsn_vector u, dmnsn_vector d)
-{
-  /* 1 division, 9 multiplications, 4 additions */
-  return dmnsn_vector_mul(dmnsn_vector_dot(u, d)/dmnsn_vector_dot(d, d), d);
-}
-
-/** Return the magnitude of a vector. */
-DMNSN_INLINE double
-dmnsn_vector_norm(dmnsn_vector n)
-{
-  /* 1 sqrt, 3 multiplications, 2 additions */
-  return sqrt(dmnsn_vector_dot(n, n));
-}
-
-/** Return the direction of a vector. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_normalized(dmnsn_vector n)
-{
-  /* 1 sqrt, 3 divisions, 3 multiplications, 2 additions */
-  return dmnsn_vector_div(n, dmnsn_vector_norm(n));
-}
-
-/** Return the component-wise minimum of two vectors. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_min(dmnsn_vector a, dmnsn_vector b)
-{
-  return dmnsn_new_vector(
-    dmnsn_min(a.x, b.x),
-    dmnsn_min(a.y, b.y),
-    dmnsn_min(a.z, b.z)
-  );
-}
-
-/** Return the component-wise maximum of two vectors. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_vector_max(dmnsn_vector a, dmnsn_vector b)
-{
-  return dmnsn_new_vector(
-    dmnsn_max(a.x, b.x),
-    dmnsn_max(a.y, b.y),
-    dmnsn_max(a.z, b.z)
-  );
-}
-
-/** Invert a matrix. */
-dmnsn_matrix dmnsn_matrix_inverse(dmnsn_matrix A);
-
-/** Multiply two matricies. */
-dmnsn_matrix dmnsn_matrix_mul(dmnsn_matrix lhs, dmnsn_matrix rhs);
-
-/** Transform a point by a matrix. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_transform_point(dmnsn_matrix T, dmnsn_vector v)
-{
-  /* 9 multiplications, 9 additions */
-  dmnsn_vector r;
-  r.x = T.n[0][0]*v.x + T.n[0][1]*v.y + T.n[0][2]*v.z + T.n[0][3];
-  r.y = T.n[1][0]*v.x + T.n[1][1]*v.y + T.n[1][2]*v.z + T.n[1][3];
-  r.z = T.n[2][0]*v.x + T.n[2][1]*v.y + T.n[2][2]*v.z + T.n[2][3];
-  return r;
-}
-
-/** Transform a direction by a matrix. */
-DMNSN_INLINE dmnsn_vector
-dmnsn_transform_direction(dmnsn_matrix T, dmnsn_vector v)
-{
-  /* 9 multiplications, 6 additions */
-  dmnsn_vector r;
-  r.x = T.n[0][0]*v.x + T.n[0][1]*v.y + T.n[0][2]*v.z;
-  r.y = T.n[1][0]*v.x + T.n[1][1]*v.y + T.n[1][2]*v.z;
-  r.z = T.n[2][0]*v.x + T.n[2][1]*v.y + T.n[2][2]*v.z;
-  return r;
-}
-
-/**
- * Transform a pseudovector by a matrix.
- * @param[in] Tinv  The inverse of the transformation matrix.
- * @param[in] v     The pseudovector to transform
- * @return The transformed pseudovector.
- */
-DMNSN_INLINE dmnsn_vector
-dmnsn_transform_normal(dmnsn_matrix Tinv, dmnsn_vector v)
-{
-  /* Multiply by the transpose of the inverse
-     (9 multiplications, 6 additions) */
-  dmnsn_vector r;
-  r.x = Tinv.n[0][0]*v.x + Tinv.n[1][0]*v.y + Tinv.n[2][0]*v.z;
-  r.y = Tinv.n[0][1]*v.x + Tinv.n[1][1]*v.y + Tinv.n[2][1]*v.z;
-  r.z = Tinv.n[0][2]*v.x + Tinv.n[1][2]*v.y + Tinv.n[2][2]*v.z;
-  return r;
-}
-
-/** Transform a bounding box by a matrix. */
-dmnsn_bounding_box dmnsn_transform_bounding_box(dmnsn_matrix T,
-                                                dmnsn_bounding_box box);
-
-/**
- * Transform a line by a matrix.
- * \f$ n' = T(l.\vec{x_0} + l.\vec{n}) - T(l.\vec{x_0}) \f$,
- * \f$ \vec{x_0}' = T(l.\vec{x_0}) \f$
- */
-DMNSN_INLINE dmnsn_line
-dmnsn_transform_line(dmnsn_matrix T, dmnsn_line l)
-{
-  /* 18 multiplications, 15 additions */
-  dmnsn_line ret;
-  ret.x0 = dmnsn_transform_point(T, l.x0);
-  ret.n  = dmnsn_transform_direction(T, l.n);
-  return ret;
-}
-
-/**
- * Return the point at \p t on a line.
- * The point is defined by \f$ l.\vec{x_0} + t \cdot l.\vec{n} \f$
- */
-DMNSN_INLINE dmnsn_vector
-dmnsn_line_point(dmnsn_line l, double t)
-{
-  return dmnsn_vector_add(l.x0, dmnsn_vector_mul(t, l.n));
-}
-
-/** Add epsilon*l.n to l.x0, to avoid self-intersections. */
-DMNSN_INLINE dmnsn_line
-dmnsn_line_add_epsilon(dmnsn_line l)
-{
-  return dmnsn_new_line(
-    dmnsn_vector_add(
-      l.x0,
-      dmnsn_vector_mul(1.0e3*dmnsn_epsilon, l.n)
-    ),
-    l.n
-  );
-}
-
-/**
- * Construct a new symmetric bounding box.
- * @param[in] r  The extent of the bounding box from the origin.
- * @return The new bounding box.
- */
-DMNSN_INLINE dmnsn_bounding_box
-dmnsn_symmetric_bounding_box(dmnsn_vector r)
-{
-  dmnsn_vector minus_r = dmnsn_vector_negate(r);
-  dmnsn_bounding_box box = {
-    dmnsn_vector_min(r, minus_r),
-    dmnsn_vector_max(r, minus_r)
-  };
-  return box;
-}
-
-/** Return whether \p p is within the axis-aligned bounding box. */
-DMNSN_INLINE bool
-dmnsn_bounding_box_contains(dmnsn_bounding_box box, dmnsn_vector p)
-{
-  return (p.x >= box.min.x && p.y >= box.min.y && p.z >= box.min.z)
-      && (p.x <= box.max.x && p.y <= box.max.y && p.z <= box.max.z);
-}
-
-/** Return whether a bounding box is infinite. */
-DMNSN_INLINE bool
-dmnsn_bounding_box_is_infinite(dmnsn_bounding_box box)
-{
-  return box.min.x == -DMNSN_INFINITY;
-}
-
-/**
- * Expand a bounding box to contain a point
- * @param[in] box  The bounding box to expand.
- * @param[in] point  The point to swallow.
- * @return The expanded bounding box.
- */
-DMNSN_INLINE dmnsn_bounding_box
-dmnsn_bounding_box_swallow(dmnsn_bounding_box box, dmnsn_vector point)
-{
-  dmnsn_bounding_box ret = {
-    dmnsn_vector_min(box.min, point),
-    dmnsn_vector_max(box.max, point)
-  };
-  return ret;
-}
-
-/** Return whether a vector contains any NaN components. */
-DMNSN_INLINE bool
-dmnsn_vector_isnan(dmnsn_vector v)
-{
-  return dmnsn_isnan(v.x) || dmnsn_isnan(v.y) || dmnsn_isnan(v.z);
-}
-
-/** Return whether a matrix contains any NaN components. */
-DMNSN_INLINE bool
-dmnsn_matrix_isnan(dmnsn_matrix m)
-{
-  size_t i, j;
-  for (i = 0; i < 3; ++i) {
-    for (j = 0; j < 4; ++j) {
-      if (dmnsn_isnan(m.n[i][j])) {
-        return true;
-      }
-    }
-  }
-  return false;
-}
-
-/** Return whether a line contains any NaN entries. */
-DMNSN_INLINE bool
-dmnsn_line_isnan(dmnsn_line l)
-{
-  return dmnsn_vector_isnan(l.x0) || dmnsn_vector_isnan(l.n);
-}
-
-/** Return whether a bounding box has any NaN components. */
-DMNSN_INLINE bool
-dmnsn_bounding_box_isnan(dmnsn_bounding_box box)
-{
-  return dmnsn_vector_isnan(box.min) || dmnsn_vector_isnan(box.max);
-}
diff --git a/libdimension/dimension/math.h b/libdimension/dimension/math.h
index 597be36..603373f 100644
--- a/libdimension/dimension/math.h
+++ b/libdimension/dimension/math.h
@@ -20,80 +20,26 @@
 
 /**
  * @file
- * Useful math functions.
+ * Mathematical functions and types.
  */
 
-#include <math.h>
-#include <stdbool.h>
+#ifndef DMNSN_MATH_H
+#define DMNSN_MATH_H
 
-/** The smallest value considered non-zero by some numerical algorithms. */
-#define dmnsn_epsilon 1.0e-10
-
-/**
- * @def DMNSN_INFINITY
- * Expands to floating-point infinity.
- */
-#if defined(INFINITY) || DMNSN_C99
-  #define DMNSN_INFINITY INFINITY
-#else
-  #define DMNSN_INFINITY HUGE_VAL
+#ifdef __cplusplus
+extern "C" {
 #endif
 
-/** Find the minimum of two values. */
-DMNSN_INLINE double
-dmnsn_min(double a, double b)
-{
-  return a < b ? a : b;
-}
-
-/** Find the maximum of two values. */
-DMNSN_INLINE double
-dmnsn_max(double a, double b)
-{
-  return a > b ? a : b;
-}
-
-/** Clamp a value to an interval. */
-DMNSN_INLINE double
-dmnsn_clamp(double n, double min, double max)
-{
-  return dmnsn_min(dmnsn_max(n, min), max);
-}
-
-/** Convert degrees to radians. */
-DMNSN_INLINE double
-dmnsn_radians(double degrees)
-{
-  return degrees*(atan(1.0)/45.0);
-}
+#include <dimension/base.h>
 
-/** Convert radians to degrees. */
-DMNSN_INLINE double
-dmnsn_degrees(double radians)
-{
-  return radians*(45.0/atan(1.0));
-}
+#include <dimension/math/scalar.h>
+#include <dimension/math/vector.h>
+#include <dimension/math/ray.h>
+#include <dimension/math/aabb.h>
+#include <dimension/math/matrix.h>
 
-/** Signum function: return the sign of a value. */
-DMNSN_INLINE int
-dmnsn_sgn(double n)
-{
-  if (n > 0.0) {
-    return 1;
-  } else if (n < 0.0) {
-    return -1;
-  } else {
-    return 0;
-  }
+#ifdef __cplusplus
 }
-
-/** Return whether a value is NaN. */
-DMNSN_INLINE bool
-dmnsn_isnan(double n)
-{
-#if DMNSN_C99
-  return isnan(n);
-#else
-  return n != n;
 #endif
-}
+
+#endif /* DMNSN_MATH_H */
diff --git a/libdimension/dimension/math/aabb.h b/libdimension/dimension/math/aabb.h
new file mode 100644
index 0000000..14cc575
--- /dev/null
+++ b/libdimension/dimension/math/aabb.h
@@ -0,0 +1,129 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Axis-aligned bounding boxes.
+ */
+
+#ifndef DMNSN_MATH_H
+#error "Please include <dimension/math.h> instead of this header directly."
+#endif
+
+/** An axis-aligned bounding box. */
+typedef struct dmnsn_aabb {
+  dmnsn_vector min; /**< The coordinate-wise minimum extent of the box. */
+  dmnsn_vector max; /**< The coordinate-wise maximum extent of the box. */
+} dmnsn_aabb;
+
+/** A standard format string for bounding boxes. */
+#define DMNSN_AABB_FORMAT "(<%g, %g, %g> ==> <%g, %g, %g>)"
+/** The appropriate arguements to printf() a bounding box. */
+#define DMNSN_AABB_PRINTF(box)                                          \
+  DMNSN_VECTOR_PRINTF((box).min), DMNSN_VECTOR_PRINTF((box).max)
+
+/**
+ * Construct a new bounding box.
+ * @param[in] min  The minimal extent of the bounding box.
+ * @param[in] max  The maximal extent of the bounding box.
+ * @return The new bounding box.
+ */
+DMNSN_INLINE dmnsn_aabb
+dmnsn_new_aabb(dmnsn_vector min, dmnsn_vector max)
+{
+  dmnsn_aabb box = { min, max };
+  return box;
+}
+
+/** Return the bounding box which contains nothing. */
+DMNSN_INLINE dmnsn_aabb
+dmnsn_zero_aabb(void)
+{
+  dmnsn_aabb box = {
+    {  DMNSN_INFINITY,  DMNSN_INFINITY,  DMNSN_INFINITY },
+    { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY }
+  };
+  return box;
+}
+
+/** Return the bounding box which contains everything. */
+DMNSN_INLINE dmnsn_aabb
+dmnsn_infinite_aabb(void)
+{
+  dmnsn_aabb box = {
+    { -DMNSN_INFINITY, -DMNSN_INFINITY, -DMNSN_INFINITY },
+    {  DMNSN_INFINITY,  DMNSN_INFINITY,  DMNSN_INFINITY }
+  };
+  return box;
+}
+
+/**
+ * Construct a new symmetric bounding box.
+ * @param[in] r  The extent of the bounding box from the origin.
+ * @return The new bounding box.
+ */
+DMNSN_INLINE dmnsn_aabb
+dmnsn_symmetric_aabb(dmnsn_vector r)
+{
+  dmnsn_vector minus_r = dmnsn_vector_negate(r);
+  dmnsn_aabb box = {
+    dmnsn_vector_min(r, minus_r),
+    dmnsn_vector_max(r, minus_r)
+  };
+  return box;
+}
+
+/** Return whether \p p is within the axis-aligned bounding box. */
+DMNSN_INLINE bool
+dmnsn_aabb_contains(dmnsn_aabb box, dmnsn_vector p)
+{
+  return (p.x >= box.min.x && p.y >= box.min.y && p.z >= box.min.z)
+      && (p.x <= box.max.x && p.y <= box.max.y && p.z <= box.max.z);
+}
+
+/** Return whether a bounding box is infinite. */
+DMNSN_INLINE bool
+dmnsn_aabb_is_infinite(dmnsn_aabb box)
+{
+  return box.min.x == -DMNSN_INFINITY;
+}
+
+/**
+ * Expand a bounding box to contain a point
+ * @param[in] box  The bounding box to expand.
+ * @param[in] point  The point to swallow.
+ * @return The expanded bounding box.
+ */
+DMNSN_INLINE dmnsn_aabb
+dmnsn_aabb_swallow(dmnsn_aabb box, dmnsn_vector point)
+{
+  dmnsn_aabb ret = {
+    dmnsn_vector_min(box.min, point),
+    dmnsn_vector_max(box.max, point)
+  };
+  return ret;
+}
+
+/** Return whether a bounding box has any NaN components. */
+DMNSN_INLINE bool
+dmnsn_aabb_isnan(dmnsn_aabb box)
+{
+  return dmnsn_vector_isnan(box.min) || dmnsn_vector_isnan(box.max);
+}
diff --git a/libdimension/dimension/math/matrix.h b/libdimension/dimension/math/matrix.h
new file mode 100644
index 0000000..7471bf5
--- /dev/null
+++ b/libdimension/dimension/math/matrix.h
@@ -0,0 +1,193 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Affine transformation matrices.
+ */
+
+#ifndef DMNSN_MATH_H
+#error "Please include <dimension/math.h> instead of this header directly."
+#endif
+
+/** A 4x4 affine transformation matrix, with implied [0 0 0 1] bottom row. */
+typedef struct dmnsn_matrix {
+  double n[3][4]; /**< The matrix elements in row-major order. */
+} dmnsn_matrix;
+
+/** A standard format string for matricies. */
+#define DMNSN_MATRIX_FORMAT                     \
+  "[%g\t%g\t%g\t%g]\n"                          \
+  "[%g\t%g\t%g\t%g]\n"                          \
+  "[%g\t%g\t%g\t%g]\n"                          \
+  "[%g\t%g\t%g\t%g]"
+/** The appropriate arguements to printf() a matrix. */
+#define DMNSN_MATRIX_PRINTF(m)                          \
+  (m).n[0][0], (m).n[0][1], (m).n[0][2], (m).n[0][3],   \
+  (m).n[1][0], (m).n[1][1], (m).n[1][2], (m).n[1][3],   \
+  (m).n[2][0], (m).n[2][1], (m).n[2][2], (m).n[2][3],   \
+  0.0, 0.0, 0.0, 1.0
+
+/** Construct a new transformation matrix. */
+DMNSN_INLINE dmnsn_matrix
+dmnsn_new_matrix(double a0, double a1, double a2, double a3,
+                 double b0, double b1, double b2, double b3,
+                 double c0, double c1, double c2, double c3)
+{
+  dmnsn_matrix m = { { { a0, a1, a2, a3 },
+                       { b0, b1, b2, b3 },
+                       { c0, c1, c2, c3 } } };
+  return m;
+}
+
+/** Construct a new transformation matrix from column vectors. */
+DMNSN_INLINE dmnsn_matrix
+dmnsn_new_matrix4(dmnsn_vector a, dmnsn_vector b, dmnsn_vector c,
+                  dmnsn_vector d)
+{
+  dmnsn_matrix m = { { { a.x, b.x, c.x, d.x },
+                       { a.y, b.y, c.y, d.y },
+                       { a.z, b.z, c.z, d.z } } };
+  return m;
+}
+
+/** Extract column vectors from a matrix. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_matrix_column(dmnsn_matrix M, unsigned int i)
+{
+  return dmnsn_new_vector(M.n[0][i], M.n[1][i], M.n[2][i]);
+}
+
+/** Return the identity matrix. */
+dmnsn_matrix dmnsn_identity_matrix(void);
+
+/**
+ * A scale transformation.
+ * @param[in] s  A vector with components representing the scaling factor in
+ *               each axis.
+ * @return The transformation matrix.
+ */
+dmnsn_matrix dmnsn_scale_matrix(dmnsn_vector s);
+/**
+ * A translation.
+ * @param[in] d  The vector to translate by.
+ * @return The transformation matrix.
+ */
+dmnsn_matrix dmnsn_translation_matrix(dmnsn_vector d);
+/**
+ * A left-handed rotation.
+ * @param[in] theta  A vector representing an axis and angle.
+ *                   @f$ axis = \vec{\theta}/|\vec{\theta}| @f$,
+ *                   @f$ angle = |\vec{\theta}| @f$
+ * @return The transformation matrix.
+ */
+dmnsn_matrix dmnsn_rotation_matrix(dmnsn_vector theta);
+
+/**
+ * An alignment matrix.
+ * @param[in] from   The initial vector.
+ * @param[in] to     The desired direction.
+ * @param[in] axis1  The first axis about which to rotate.
+ * @param[in] axis2  The second axis about which to rotate.
+ * @return A transformation matrix that will rotate \p from to \p to.
+ */
+dmnsn_matrix dmnsn_alignment_matrix(dmnsn_vector from, dmnsn_vector to,
+                                    dmnsn_vector axis1, dmnsn_vector axis2);
+
+/** Invert a matrix. */
+dmnsn_matrix dmnsn_matrix_inverse(dmnsn_matrix A);
+
+/** Multiply two matricies. */
+dmnsn_matrix dmnsn_matrix_mul(dmnsn_matrix lhs, dmnsn_matrix rhs);
+
+/** Transform a point by a matrix. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_transform_point(dmnsn_matrix T, dmnsn_vector v)
+{
+  /* 9 multiplications, 9 additions */
+  dmnsn_vector r;
+  r.x = T.n[0][0]*v.x + T.n[0][1]*v.y + T.n[0][2]*v.z + T.n[0][3];
+  r.y = T.n[1][0]*v.x + T.n[1][1]*v.y + T.n[1][2]*v.z + T.n[1][3];
+  r.z = T.n[2][0]*v.x + T.n[2][1]*v.y + T.n[2][2]*v.z + T.n[2][3];
+  return r;
+}
+
+/** Transform a direction by a matrix. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_transform_direction(dmnsn_matrix T, dmnsn_vector v)
+{
+  /* 9 multiplications, 6 additions */
+  dmnsn_vector r;
+  r.x = T.n[0][0]*v.x + T.n[0][1]*v.y + T.n[0][2]*v.z;
+  r.y = T.n[1][0]*v.x + T.n[1][1]*v.y + T.n[1][2]*v.z;
+  r.z = T.n[2][0]*v.x + T.n[2][1]*v.y + T.n[2][2]*v.z;
+  return r;
+}
+
+/**
+ * Transform a pseudovector by a matrix.
+ * @param[in] Tinv  The inverse of the transformation matrix.
+ * @param[in] v     The pseudovector to transform
+ * @return The transformed pseudovector.
+ */
+DMNSN_INLINE dmnsn_vector
+dmnsn_transform_normal(dmnsn_matrix Tinv, dmnsn_vector v)
+{
+  /* Multiply by the transpose of the inverse
+     (9 multiplications, 6 additions) */
+  dmnsn_vector r;
+  r.x = Tinv.n[0][0]*v.x + Tinv.n[1][0]*v.y + Tinv.n[2][0]*v.z;
+  r.y = Tinv.n[0][1]*v.x + Tinv.n[1][1]*v.y + Tinv.n[2][1]*v.z;
+  r.z = Tinv.n[0][2]*v.x + Tinv.n[1][2]*v.y + Tinv.n[2][2]*v.z;
+  return r;
+}
+
+/**
+ * Transform a ray by a matrix.
+ * \f$ n' = T(l.\vec{x_0} + l.\vec{n}) - T(l.\vec{x_0}) \f$,
+ * \f$ \vec{x_0}' = T(l.\vec{x_0}) \f$
+ */
+DMNSN_INLINE dmnsn_ray
+dmnsn_transform_ray(dmnsn_matrix T, dmnsn_ray l)
+{
+  /* 18 multiplications, 15 additions */
+  dmnsn_ray ret;
+  ret.x0 = dmnsn_transform_point(T, l.x0);
+  ret.n  = dmnsn_transform_direction(T, l.n);
+  return ret;
+}
+
+/** Transform a bounding box by a matrix. */
+dmnsn_aabb dmnsn_transform_aabb(dmnsn_matrix T, dmnsn_aabb box);
+
+/** Return whether a matrix contains any NaN components. */
+DMNSN_INLINE bool
+dmnsn_matrix_isnan(dmnsn_matrix m)
+{
+  size_t i, j;
+  for (i = 0; i < 3; ++i) {
+    for (j = 0; j < 4; ++j) {
+      if (dmnsn_isnan(m.n[i][j])) {
+        return true;
+      }
+    }
+  }
+  return false;
+}
diff --git a/libdimension/dimension/math/ray.h b/libdimension/dimension/math/ray.h
new file mode 100644
index 0000000..1dd98f7
--- /dev/null
+++ b/libdimension/dimension/math/ray.h
@@ -0,0 +1,83 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Lines in 3-D space.
+ */
+
+#ifndef DMNSN_MATH_H
+#error "Please include <dimension/math.h> instead of this header directly."
+#endif
+
+/** A line, or ray. */
+typedef struct dmnsn_ray {
+  dmnsn_vector x0; /**< The origin of the ray. */
+  dmnsn_vector n;  /**< The direction of the ray. */
+} dmnsn_ray;
+
+/** A standard format string for rays. */
+#define DMNSN_RAY_FORMAT "(<%g, %g, %g> + t*<%g, %g, %g>)"
+/** The appropriate arguements to printf() a ray. */
+#define DMNSN_RAY_PRINTF(l)                                     \
+  DMNSN_VECTOR_PRINTF((l).x0), DMNSN_VECTOR_PRINTF((l).n)
+
+/**
+ * Construct a new ray.
+ * @param[in] x0  The origin of the ray.
+ * @param[in] n   The direction of the ray.
+ * @return The new ray.
+ */
+DMNSN_INLINE dmnsn_ray
+dmnsn_new_ray(dmnsn_vector x0, dmnsn_vector n)
+{
+  dmnsn_ray l = { x0, n };
+  return l;
+}
+
+/**
+ * Return the point at \p t on a ray.
+ * The point is defined by \f$ l.\vec{x_0} + t \cdot l.\vec{n} \f$
+ */
+DMNSN_INLINE dmnsn_vector
+dmnsn_ray_point(dmnsn_ray l, double t)
+{
+  return dmnsn_vector_add(l.x0, dmnsn_vector_mul(t, l.n));
+}
+
+/** Add epsilon*l.n to l.x0, to avoid self-intersections. */
+DMNSN_INLINE dmnsn_ray
+dmnsn_ray_add_epsilon(dmnsn_ray l)
+{
+  return dmnsn_new_ray(
+    dmnsn_vector_add(
+      l.x0,
+      dmnsn_vector_mul(1.0e3*dmnsn_epsilon, l.n)
+    ),
+    l.n
+  );
+}
+
+/** Return whether a ray contains any NaN entries. */
+DMNSN_INLINE bool
+dmnsn_ray_isnan(dmnsn_ray l)
+{
+  return dmnsn_vector_isnan(l.x0) || dmnsn_vector_isnan(l.n);
+}
diff --git a/libdimension/dimension/math/scalar.h b/libdimension/dimension/math/scalar.h
new file mode 100644
index 0000000..3887c14
--- /dev/null
+++ b/libdimension/dimension/math/scalar.h
@@ -0,0 +1,103 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Mathematical functions of one variable.
+ */
+
+#ifndef DMNSN_MATH_H
+#error "Please include <dimension/math.h> instead of this header directly."
+#endif
+
+#include <math.h>
+#include <stdbool.h>
+
+/** The smallest value considered non-zero by some numerical algorithms. */
+#define dmnsn_epsilon 1.0e-10
+
+/**
+ * @def DMNSN_INFINITY
+ * Expands to floating-point infinity.
+ */
+#if defined(INFINITY) || DMNSN_C99
+  #define DMNSN_INFINITY INFINITY
+#else
+  #define DMNSN_INFINITY HUGE_VAL
+#endif
+
+/** Find the minimum of two values. */
+DMNSN_INLINE double
+dmnsn_min(double a, double b)
+{
+  return a < b ? a : b;
+}
+
+/** Find the maximum of two values. */
+DMNSN_INLINE double
+dmnsn_max(double a, double b)
+{
+  return a > b ? a : b;
+}
+
+/** Clamp a value to an interval. */
+DMNSN_INLINE double
+dmnsn_clamp(double n, double min, double max)
+{
+  return dmnsn_min(dmnsn_max(n, min), max);
+}
+
+/** Convert degrees to radians. */
+DMNSN_INLINE double
+dmnsn_radians(double degrees)
+{
+  return degrees*(atan(1.0)/45.0);
+}
+
+/** Convert radians to degrees. */
+DMNSN_INLINE double
+dmnsn_degrees(double radians)
+{
+  return radians*(45.0/atan(1.0));
+}
+
+/** Signum function: return the sign of a value. */
+DMNSN_INLINE int
+dmnsn_sgn(double n)
+{
+  if (n > 0.0) {
+    return 1;
+  } else if (n < 0.0) {
+    return -1;
+  } else {
+    return 0;
+  }
+}
+
+/** Return whether a value is NaN. */
+DMNSN_INLINE bool
+dmnsn_isnan(double n)
+{
+#if DMNSN_C99
+  return isnan(n);
+#else
+  return n != n;
+#endif
+}
diff --git a/libdimension/dimension/math/vector.h b/libdimension/dimension/math/vector.h
new file mode 100644
index 0000000..8eacee9
--- /dev/null
+++ b/libdimension/dimension/math/vector.h
@@ -0,0 +1,183 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Vectors in 3-D space.
+ */
+
+#ifndef DMNSN_MATH_H
+#error "Please include <dimension/math.h> instead of this header directly."
+#endif
+
+#include <math.h>
+#include <stdbool.h>
+
+/** A vector in 3 dimensions. */
+typedef struct dmnsn_vector {
+  double x; /**< The x component. */
+  double y; /**< The y component. */
+  double z; /**< The z component. */
+} dmnsn_vector;
+
+/** A standard format string for vectors. */
+#define DMNSN_VECTOR_FORMAT "<%g, %g, %g>"
+/** The appropriate arguements to printf() a vector. */
+#define DMNSN_VECTOR_PRINTF(v) (v).x, (v).y, (v).z
+
+/* Constants */
+
+/** The zero vector. */
+static const dmnsn_vector dmnsn_zero = { 0.0, 0.0, 0.0 };
+/** The x vector. */
+static const dmnsn_vector dmnsn_x = { 1.0, 0.0, 0.0 };
+/** The y vector. */
+static const dmnsn_vector dmnsn_y = { 0.0, 1.0, 0.0 };
+/** The z vector. */
+static const dmnsn_vector dmnsn_z = { 0.0, 0.0, 1.0 };
+
+/* Shorthand for vector construction */
+
+/** Construct a new vector. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_new_vector(double x, double y, double z)
+{
+  dmnsn_vector v = { x, y, z };
+  return v;
+}
+
+/* Vector arithmetic */
+
+/** Negate a vector. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_negate(dmnsn_vector rhs)
+{
+  /* 3 negations */
+  dmnsn_vector v = { -rhs.x, -rhs.y, -rhs.z };
+  return v;
+}
+
+/** Add two vectors. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_add(dmnsn_vector lhs, dmnsn_vector rhs)
+{
+  /* 3 additions */
+  dmnsn_vector v = { lhs.x + rhs.x, lhs.y + rhs.y, lhs.z + rhs.z };
+  return v;
+}
+
+/** Subtract two vectors. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_sub(dmnsn_vector lhs, dmnsn_vector rhs)
+{
+  /* 3 additions */
+  dmnsn_vector v = { lhs.x - rhs.x, lhs.y - rhs.y, lhs.z - rhs.z };
+  return v;
+}
+
+/** Multiply a vector by a scalar. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_mul(double lhs, dmnsn_vector rhs)
+{
+  /* 3 multiplications */
+  dmnsn_vector v = { lhs*rhs.x, lhs*rhs.y, lhs*rhs.z };
+  return v;
+}
+
+/** Divide a vector by a scalar. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_div(dmnsn_vector lhs, double rhs)
+{
+  /* 3 divisions */
+  dmnsn_vector v = { lhs.x/rhs, lhs.y/rhs, lhs.z/rhs };
+  return v;
+}
+
+/** Return the dot product of two vectors. */
+DMNSN_INLINE double
+dmnsn_vector_dot(dmnsn_vector lhs, dmnsn_vector rhs)
+{
+  /* 3 multiplications, 2 additions */
+  return lhs.x*rhs.x + lhs.y*rhs.y + lhs.z*rhs.z;
+}
+
+/** Return the cross product of two vectors. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_cross(dmnsn_vector lhs, dmnsn_vector rhs)
+{
+  /* 6 multiplications, 3 additions */
+  dmnsn_vector v = { lhs.y*rhs.z - lhs.z*rhs.y,
+                     lhs.z*rhs.x - lhs.x*rhs.z,
+                     lhs.x*rhs.y - lhs.y*rhs.x };
+  return v;
+}
+
+/** Return the projection of \p u onto \p d. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_proj(dmnsn_vector u, dmnsn_vector d)
+{
+  /* 1 division, 9 multiplications, 4 additions */
+  return dmnsn_vector_mul(dmnsn_vector_dot(u, d)/dmnsn_vector_dot(d, d), d);
+}
+
+/** Return the magnitude of a vector. */
+DMNSN_INLINE double
+dmnsn_vector_norm(dmnsn_vector n)
+{
+  /* 1 sqrt, 3 multiplications, 2 additions */
+  return sqrt(dmnsn_vector_dot(n, n));
+}
+
+/** Return the direction of a vector. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_normalized(dmnsn_vector n)
+{
+  /* 1 sqrt, 3 divisions, 3 multiplications, 2 additions */
+  return dmnsn_vector_div(n, dmnsn_vector_norm(n));
+}
+
+/** Return the component-wise minimum of two vectors. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_min(dmnsn_vector a, dmnsn_vector b)
+{
+  return dmnsn_new_vector(
+    dmnsn_min(a.x, b.x),
+    dmnsn_min(a.y, b.y),
+    dmnsn_min(a.z, b.z)
+  );
+}
+
+/** Return the component-wise maximum of two vectors. */
+DMNSN_INLINE dmnsn_vector
+dmnsn_vector_max(dmnsn_vector a, dmnsn_vector b)
+{
+  return dmnsn_new_vector(
+    dmnsn_max(a.x, b.x),
+    dmnsn_max(a.y, b.y),
+    dmnsn_max(a.z, b.z)
+  );
+}
+
+/** Return whether a vector contains any NaN components. */
+DMNSN_INLINE bool
+dmnsn_vector_isnan(dmnsn_vector v)
+{
+  return dmnsn_isnan(v.x) || dmnsn_isnan(v.y) || dmnsn_isnan(v.z);
+}
diff --git a/libdimension/dimension/model.h b/libdimension/dimension/model.h
new file mode 100644
index 0000000..5367729
--- /dev/null
+++ b/libdimension/dimension/model.h
@@ -0,0 +1,59 @@
+/*************************************************************************
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Scene modeling.
+ */
+
+#ifndef DMNSN_MODEL_H
+#define DMNSN_MODEL_H
+
+#include "dimension/base.h"
+#include "dimension/color.h"
+#include "dimension/canvas.h"
+#include "dimension/pattern.h"
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#include <dimension/platform.h>
+
+#include <dimension/model/pigment.h>
+#include <dimension/model/pigments.h>
+#include <dimension/model/finish.h>
+#include <dimension/model/finishes.h>
+#include <dimension/model/texture.h>
+#include <dimension/model/interior.h>
+#include <dimension/model/object.h>
+#include <dimension/model/objects.h>
+#include <dimension/model/csg.h>
+#include <dimension/model/light.h>
+#include <dimension/model/lights.h>
+#include <dimension/model/camera.h>
+#include <dimension/model/cameras.h>
+#include <dimension/model/scene.h>
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* DMNSN_MODEL_H */
diff --git a/libdimension/dimension/camera.h b/libdimension/dimension/model/camera.h
index 3c6494c..37f80b9 100644
--- a/libdimension/dimension/camera.h
+++ b/libdimension/dimension/model/camera.h
@@ -33,8 +33,7 @@ typedef struct dmnsn_camera dmnsn_camera;
  * @param[in] y       The y coordinate of the pixel (in [0, 1]).
  * @return The ray through (\p x, \p y).
  */
-typedef dmnsn_line dmnsn_camera_ray_fn(const dmnsn_camera *camera,
-                                       double x, double y);
+typedef dmnsn_ray dmnsn_camera_ray_fn(const dmnsn_camera *camera, double x, double y);
 
 /** A camera. */
 struct dmnsn_camera {
@@ -64,4 +63,4 @@ void dmnsn_init_camera(dmnsn_camera *camera);
  * @param[in] y       The y coordinate of the pixel (in [0, 1]).
  * @return The ray through (\p x, \p y).
  */
-dmnsn_line dmnsn_camera_ray(const dmnsn_camera *camera, double x, double y);
+dmnsn_ray dmnsn_camera_ray(const dmnsn_camera *camera, double x, double y);
diff --git a/libdimension/dimension/cameras.h b/libdimension/dimension/model/cameras.h
index 9ef2646..9ef2646 100644
--- a/libdimension/dimension/cameras.h
+++ b/libdimension/dimension/model/cameras.h
diff --git a/libdimension/dimension/csg.h b/libdimension/dimension/model/csg.h
index b2ce83f..b2ce83f 100644
--- a/libdimension/dimension/csg.h
+++ b/libdimension/dimension/model/csg.h
diff --git a/libdimension/dimension/finish.h b/libdimension/dimension/model/finish.h
index d975877..d975877 100644
--- a/libdimension/dimension/finish.h
+++ b/libdimension/dimension/model/finish.h
diff --git a/libdimension/dimension/finishes.h b/libdimension/dimension/model/finishes.h
index e1f7b44..e1f7b44 100644
--- a/libdimension/dimension/finishes.h
+++ b/libdimension/dimension/model/finishes.h
diff --git a/libdimension/dimension/interior.h b/libdimension/dimension/model/interior.h
index 0ff697d..0ff697d 100644
--- a/libdimension/dimension/interior.h
+++ b/libdimension/dimension/model/interior.h
diff --git a/libdimension/dimension/light.h b/libdimension/dimension/model/light.h
index 218611d..218611d 100644
--- a/libdimension/dimension/light.h
+++ b/libdimension/dimension/model/light.h
diff --git a/libdimension/dimension/lights.h b/libdimension/dimension/model/lights.h
index e7de4cc..e7de4cc 100644
--- a/libdimension/dimension/lights.h
+++ b/libdimension/dimension/model/lights.h
diff --git a/libdimension/dimension/object.h b/libdimension/dimension/model/object.h
index c8e9d86..287d838 100644
--- a/libdimension/dimension/object.h
+++ b/libdimension/dimension/model/object.h
@@ -30,8 +30,8 @@ typedef struct dmnsn_object dmnsn_object;
 
 /** A type to represent a ray-object intersection. */
 typedef struct dmnsn_intersection {
-  dmnsn_line ray; /**< The ray that intersected. */
-  double t;       /**< The line index that intersected. */
+  dmnsn_ray ray; /**< The ray that intersected. */
+  double t;      /**< The ray index that intersected. */
 
   /** The surface normal at the intersection point. */
   dmnsn_vector normal;
@@ -43,12 +43,12 @@ typedef struct dmnsn_intersection {
 /**
  * Ray-object intersection callback.
  * @param[in]  object        The object to test.
- * @param[in]  line          The line to test.
+ * @param[in]  ray           The ray to test.
  * @param[out] intersection  Where to store the intersection details of the
  *                           closest (if any) intersection.
- * @return Whether \p line intersected \p object.
+ * @return Whether \p ray intersected \p object.
  */
-typedef bool dmnsn_object_intersection_fn(const dmnsn_object *object, dmnsn_line line, dmnsn_intersection *intersection);
+typedef bool dmnsn_object_intersection_fn(const dmnsn_object *object, dmnsn_ray ray, dmnsn_intersection *intersection);
 
 /**
  * Object inside callback.
@@ -62,8 +62,9 @@ typedef bool dmnsn_object_inside_fn(const dmnsn_object *object, dmnsn_vector poi
  * Object bounding callback.
  * @param[in,out] object  The object to bound.
  * @param[in,out] trans  The effective transformation for the object.
+ * @return A world-coordinate bounding box computed for the object.
  */
-typedef dmnsn_bounding_box dmnsn_object_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans);
+typedef dmnsn_aabb dmnsn_object_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans);
 
 /**
  * Object precomputation callback.
@@ -96,7 +97,7 @@ struct dmnsn_object {
   bool precomputed; /**< @internal Whether the object is precomputed yet. */
   dmnsn_matrix trans_inv; /**< Inverse of the transformation matrix. */
   dmnsn_matrix pigment_trans; /**< Inverse transformation for the texture. */
-  dmnsn_bounding_box bounding_box; /**< Bounding box in world coordinates. */
+  dmnsn_aabb aabb; /**< Bounding box in world coordinates. */
 };
 
 /**
@@ -121,19 +122,19 @@ void dmnsn_object_precompute(dmnsn_object *object);
 /**
  * Appropriately transform a ray, then test for an intersection.
  * @param[in]  object        The object to test.
- * @param[in]  line          The ray to test.
+ * @param[in]  ray           The ray to test.
  * @param[out] intersection  Where to store the intersection details.
  * @return Whether there was an intersection.
  */
 DMNSN_INLINE bool
-dmnsn_object_intersection(const dmnsn_object *object, dmnsn_line line,
+dmnsn_object_intersection(const dmnsn_object *object, dmnsn_ray ray,
                           dmnsn_intersection *intersection)
 {
-  dmnsn_line line_trans = dmnsn_transform_line(object->trans_inv, line);
+  dmnsn_ray ray_trans = dmnsn_transform_ray(object->trans_inv, ray);
   intersection->object = NULL;
-  if (object->vtable->intersection_fn(object, line_trans, intersection)) {
+  if (object->vtable->intersection_fn(object, ray_trans, intersection)) {
     /* Get us back into world coordinates */
-    intersection->ray = line;
+    intersection->ray = ray;
     intersection->normal = dmnsn_vector_normalized(
       dmnsn_transform_normal(object->trans_inv, intersection->normal)
     );
diff --git a/libdimension/dimension/objects.h b/libdimension/dimension/model/objects.h
index 2865d82..2865d82 100644
--- a/libdimension/dimension/objects.h
+++ b/libdimension/dimension/model/objects.h
diff --git a/libdimension/dimension/pigment.h b/libdimension/dimension/model/pigment.h
index 14d8bae..14d8bae 100644
--- a/libdimension/dimension/pigment.h
+++ b/libdimension/dimension/model/pigment.h
diff --git a/libdimension/dimension/pigments.h b/libdimension/dimension/model/pigments.h
index 100016d..100016d 100644
--- a/libdimension/dimension/pigments.h
+++ b/libdimension/dimension/model/pigments.h
diff --git a/libdimension/dimension/scene.h b/libdimension/dimension/model/scene.h
index dd0c1ba..dd0c1ba 100644
--- a/libdimension/dimension/scene.h
+++ b/libdimension/dimension/model/scene.h
diff --git a/libdimension/dimension/texture.h b/libdimension/dimension/model/texture.h
index df08a4a..df08a4a 100644
--- a/libdimension/dimension/texture.h
+++ b/libdimension/dimension/model/texture.h
diff --git a/libdimension/dimension/pattern.h b/libdimension/dimension/pattern.h
index 2caed07..17a43a9 100644
--- a/libdimension/dimension/pattern.h
+++ b/libdimension/dimension/pattern.h
@@ -1,5 +1,5 @@
 /*************************************************************************
- * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
  *                                                                       *
  * This file is part of The Dimension Library.                           *
  *                                                                       *
@@ -20,43 +20,25 @@
 
 /**
  * @file
- * Patterns.  Patterns are functions which map vectors to scalars, which are
- * used for pigments and normals.
+ * Patterns.
  */
 
-/* Forward-declare dmnsn_pattern */
-typedef struct dmnsn_pattern dmnsn_pattern;
+#ifndef DMNSN_PATTERN_H
+#define DMNSN_PATTERN_H
 
-/**
- * Pattern callback.
- * @param[in] pattern  The pattern itself.
- * @param[in] v        The point at which to evaluate the pattern.
- * @return The value of the pattern at \p v.
- */
-typedef double dmnsn_pattern_fn(const dmnsn_pattern *pattern, dmnsn_vector v);
+#ifdef __cplusplus
+extern "C" {
+#endif
 
-/** A pattern. */
-struct dmnsn_pattern {
-  dmnsn_pattern_fn *pattern_fn; /**< The pattern callback. */
-};
+#include <dimension/base.h>
+#include <dimension/math.h>
 
-/**
- * Allocate a dummy pattern.
- * @param[in] pool  The memory pool to allocate from.
- * @return A pattern with no callbacks set.
- */
-dmnsn_pattern *dmnsn_new_pattern(dmnsn_pool *pool);
+#include <dimension/pattern/pattern.h>
+#include <dimension/pattern/patterns.h>
+#include <dimension/pattern/map.h>
 
-/**
- * Initialize a dmnsn_pattern field.
- * @param[out] pattern  The pattern to initialize.
- */
-void dmnsn_init_pattern(dmnsn_pattern *pattern);
+#ifdef __cplusplus
+}
+#endif
 
-/**
- * Invoke the pattern callback.
- * @param[in] pattern  The pattern to evaluate.
- * @param[in] v        The point to get the pattern value for.
- * @return The value of the pattern at \p v.
- */
-double dmnsn_pattern_value(const dmnsn_pattern *pattern, dmnsn_vector v);
+#endif /* DMNSN_PATTERN_H */
diff --git a/libdimension/dimension/map.h b/libdimension/dimension/pattern/map.h
index 7229a24..7229a24 100644
--- a/libdimension/dimension/map.h
+++ b/libdimension/dimension/pattern/map.h
diff --git a/libdimension/dimension/pattern/pattern.h b/libdimension/dimension/pattern/pattern.h
new file mode 100644
index 0000000..2caed07
--- /dev/null
+++ b/libdimension/dimension/pattern/pattern.h
@@ -0,0 +1,62 @@
+/*************************************************************************
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Patterns.  Patterns are functions which map vectors to scalars, which are
+ * used for pigments and normals.
+ */
+
+/* Forward-declare dmnsn_pattern */
+typedef struct dmnsn_pattern dmnsn_pattern;
+
+/**
+ * Pattern callback.
+ * @param[in] pattern  The pattern itself.
+ * @param[in] v        The point at which to evaluate the pattern.
+ * @return The value of the pattern at \p v.
+ */
+typedef double dmnsn_pattern_fn(const dmnsn_pattern *pattern, dmnsn_vector v);
+
+/** A pattern. */
+struct dmnsn_pattern {
+  dmnsn_pattern_fn *pattern_fn; /**< The pattern callback. */
+};
+
+/**
+ * Allocate a dummy pattern.
+ * @param[in] pool  The memory pool to allocate from.
+ * @return A pattern with no callbacks set.
+ */
+dmnsn_pattern *dmnsn_new_pattern(dmnsn_pool *pool);
+
+/**
+ * Initialize a dmnsn_pattern field.
+ * @param[out] pattern  The pattern to initialize.
+ */
+void dmnsn_init_pattern(dmnsn_pattern *pattern);
+
+/**
+ * Invoke the pattern callback.
+ * @param[in] pattern  The pattern to evaluate.
+ * @param[in] v        The point to get the pattern value for.
+ * @return The value of the pattern at \p v.
+ */
+double dmnsn_pattern_value(const dmnsn_pattern *pattern, dmnsn_vector v);
diff --git a/libdimension/dimension/patterns.h b/libdimension/dimension/pattern/patterns.h
index 59b7fec..59b7fec 100644
--- a/libdimension/dimension/patterns.h
+++ b/libdimension/dimension/pattern/patterns.h
diff --git a/libdimension/dimension/platform.h b/libdimension/dimension/platform.h
new file mode 100644
index 0000000..4c01709
--- /dev/null
+++ b/libdimension/dimension/platform.h
@@ -0,0 +1,39 @@
+/*************************************************************************
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Platform abstractions.
+ */
+
+#ifndef DMNSN_PLATFORM_H
+#define DMNSN_PLATFORM_H
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#include <dimension/platform/timer.h>
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* DMNSN_PLATFORM_H */
diff --git a/libdimension/dimension/timer.h b/libdimension/dimension/platform/timer.h
index 55d66ff..bcf4a9d 100644
--- a/libdimension/dimension/timer.h
+++ b/libdimension/dimension/platform/timer.h
@@ -23,6 +23,10 @@
  * A platform-agnostic timer abstraction.
  */
 
+#ifndef DMNSN_PLATFORM_H
+#error "Please include <dimension/platform.h> instead of this header directly."
+#endif
+
 /** A platform-agnotic timer. */
 typedef struct dmnsn_timer {
   double real;   /**< Wall-clock time. */
diff --git a/libdimension/dimension/render.h b/libdimension/dimension/render.h
new file mode 100644
index 0000000..296d7a3
--- /dev/null
+++ b/libdimension/dimension/render.h
@@ -0,0 +1,41 @@
+/*************************************************************************
+ * Copyright (C) 2014 Tavian Barnes <tavianator@tavianator.com>          *
+ *                                                                       *
+ * This file is part of The Dimension Library.                           *
+ *                                                                       *
+ * The Dimension Library is free software; you can redistribute it and/  *
+ * or modify it under the terms of the GNU Lesser General Public License *
+ * as published by the Free Software Foundation; either version 3 of the *
+ * License, or (at your option) any later version.                       *
+ *                                                                       *
+ * The Dimension Library is distributed in the hope that it will be      *
+ * useful, but WITHOUT ANY WARRANTY; without even the implied warranty   *
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU  *
+ * Lesser General Public License for more details.                       *
+ *                                                                       *
+ * You should have received a copy of the GNU Lesser General Public      *
+ * License along with this program.  If not, see                         *
+ * <http://www.gnu.org/licenses/>.                                       *
+ *************************************************************************/
+
+/**
+ * @file
+ * Rendering.
+ */
+
+#ifndef DMNSN_RENDER_H
+#define DMNSN_RENDER_H
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#include <dimension/model.h>
+
+#include <dimension/render/render.h>
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* DMNSN_RENDER_H */
diff --git a/libdimension/dimension/ray_trace.h b/libdimension/dimension/render/render.h
index 82fb759..69f52c2 100644
--- a/libdimension/dimension/ray_trace.h
+++ b/libdimension/dimension/render/render.h
@@ -1,5 +1,5 @@
 /*************************************************************************
- * Copyright (C) 2009-2011 Tavian Barnes <tavianator@tavianator.com>     *
+ * Copyright (C) 2009-2014 Tavian Barnes <tavianator@tavianator.com>     *
  *                                                                       *
  * This file is part of The Dimension Library.                           *
  *                                                                       *
@@ -20,18 +20,22 @@
 
 /**
  * @file
- * Ray-trace a scene.
+ * Rendering.
  */
 
+#ifndef DMNSN_RENDER_H
+#error "Please include <dimension/render.h> instead of this header directly."
+#endif
+
 /**
- * Render a scene by ray tracing.
+ * Render a scene.
  * @param[in,out] scene  The scene to render.
  */
-void dmnsn_ray_trace(dmnsn_scene *scene);
+void dmnsn_render(dmnsn_scene *scene);
 
 /**
  * Render a scene in the background.
  * @param[in,out] scene  The scene to render.
  * @return A \p dmnsn_future object.
  */
-dmnsn_future *dmnsn_ray_trace_async(dmnsn_scene *scene);
+dmnsn_future *dmnsn_render_async(dmnsn_scene *scene);
diff --git a/libdimension/dimension-internal.h b/libdimension/internal.h
index 53fa24c..3db2612 100644
--- a/libdimension/dimension-internal.h
+++ b/libdimension/internal.h
@@ -20,21 +20,16 @@
 
 /**
  * @file
- * The internal libdimension API.  These functions and types are used to
+ * The basic internal libdimension API.  These functions and types are used to
  * implement libdimension, but are not part of its public API.
  */
 
-#ifndef DIMENSION_INTERNAL_H
-#define DIMENSION_INTERNAL_H
+#ifndef DMNSN_INTERNAL_H
+#define DMNSN_INTERNAL_H
 
-#include "dimension.h"
-#include "compiler-internal.h"
-#include "profile.h"
-#include "platform.h"
-#include "future-internal.h"
-#include "threads.h"
-#include "bvh.h"
-#include "prtree.h"
-#include "rgba.h"
+#include "dimension/base.h"
 
-#endif // DIMENSION_INTERNAL_H
+#include "internal/compiler.h"
+#include "internal/profile.h"
+
+#endif // DMNSN_INTERNAL_H
diff --git a/libdimension/bvh.h b/libdimension/internal/bvh.h
index f965632..c97e2ba 100644
--- a/libdimension/bvh.h
+++ b/libdimension/internal/bvh.h
@@ -25,6 +25,12 @@
  * traversal algorithm.
  */
 
+#ifndef DMNSN_INTERNAL_BVH_H
+#define DMNSN_INTERNAL_BVH_H
+
+#include "internal.h"
+#include "dimension/math.h"
+#include "dimension/model.h"
 #include <stdbool.h>
 
 /// A bounding volume hierarchy.
@@ -43,18 +49,15 @@ DMNSN_INTERNAL dmnsn_bvh *dmnsn_new_bvh(const dmnsn_array *objects,
 DMNSN_INTERNAL void dmnsn_delete_bvh(dmnsn_bvh *bvh);
 
 /// Find the closest ray-object intersection in the tree.
-DMNSN_INTERNAL bool dmnsn_bvh_intersection(const dmnsn_bvh *bvh,
-                                           dmnsn_line ray,
-                                           dmnsn_intersection *intersection,
-                                           bool reset);
+DMNSN_INTERNAL bool dmnsn_bvh_intersection(const dmnsn_bvh *bvh, dmnsn_ray ray, dmnsn_intersection *intersection, bool reset);
 /// Determine whether a point is inside any object in the tree.
 DMNSN_INTERNAL bool dmnsn_bvh_inside(const dmnsn_bvh *bvh, dmnsn_vector point);
 /// Return the bounding box of the whole hierarchy.
-DMNSN_INTERNAL dmnsn_bounding_box dmnsn_bvh_bounding_box(const dmnsn_bvh *bvh);
+DMNSN_INTERNAL dmnsn_aabb dmnsn_bvh_aabb(const dmnsn_bvh *bvh);
 
 /// A non-flat BVH representation, used by BVH implementations.
 typedef struct dmnsn_bvh_node {
-  dmnsn_bounding_box bounding_box;   ///< The bounding box of this node.
+  dmnsn_aabb aabb;                   ///< The bounding box of this node.
   dmnsn_object *object;              ///< The object, for leaf nodes.
   unsigned int nchildren;            ///< How many children this node has.
   unsigned int max_children;         ///< Maximum number of children.
@@ -72,3 +75,5 @@ DMNSN_INTERNAL void dmnsn_delete_bvh_node(dmnsn_bvh_node *node);
 
 /// Add a child to a BVH node.
 DMNSN_INTERNAL void dmnsn_bvh_node_add(dmnsn_bvh_node *parent, dmnsn_bvh_node *child);
+
+#endif // DMNSN_INTERNAL_BVH_H
diff --git a/libdimension/compiler-internal.h b/libdimension/internal/compiler.h
index 40e749d..e3555a4 100644
--- a/libdimension/compiler-internal.h
+++ b/libdimension/internal/compiler.h
@@ -23,33 +23,12 @@
  * Internally-used compiler abstractions.
  */
 
-#include <stdalign.h>
-
-/**
- * @def dmnsn_likely
- * Indicate that a test is likely to succeed.
- * @param test  The test to perform.
- * @return The truth value of \p test.
- */
-/**
- * @def dmnsn_unlikely
- * Indicate that a test is unlikely to succeed.
- * @param test  The test to perform.
- * @return The truth value of \p test.
- */
-#ifdef DMNSN_PROFILE
-  #define dmnsn_likely(test)                                      \
-    dmnsn_expect(!!(test), true, DMNSN_FUNC, __FILE__, __LINE__)
-  #define dmnsn_unlikely(test)                                    \
-    dmnsn_expect(!!(test), false, DMNSN_FUNC, __FILE__, __LINE__)
-#elif DMNSN_GNUC
-  #define dmnsn_likely(test)   __builtin_expect(!!(test), true)
-  #define dmnsn_unlikely(test) __builtin_expect(!!(test), false)
-#else
-  #define dmnsn_likely(test)   (!!(test))
-  #define dmnsn_unlikely(test) (!!(test))
+#ifndef DMNSN_INTERNAL_H
+#error "Please include \"internal.h\" instead of this header directly."
 #endif
 
+#include <stdalign.h>
+
 /**
  * @def DMNSN_HOT
  * Mark a function as a hot path.
@@ -82,7 +61,7 @@
 /// Synonym for _Atomic that stdatomic.h doesn't define for some reason
 #define atomic _Atomic
 
-/// C11-compliant alloca variant
+/// Nearly-compliant alloca variant
 #define DMNSN_ALLOCA(var, size) DMNSN_ALLOCA_IMPL(var, size, __LINE__)
 
 #define DMNSN_ALLOCA_IMPL(var, size, ctr) DMNSN_ALLOCA_IMPL2(var, size, ctr)
diff --git a/libdimension/threads.h b/libdimension/internal/concurrency.h
index 839299e..e9c024a 100644
--- a/libdimension/threads.h
+++ b/libdimension/internal/concurrency.h
@@ -20,11 +20,30 @@
 
 /**
  * @file
- * Background threading interface.
+ * Internal threading abstraction.
  */
 
+#ifndef DMNSN_INTERNAL_CONCURRENCY_H
+#define DMNSN_INTERNAL_CONCURRENCY_H
+
+#include "internal.h"
+#include "dimension/concurrency.h"
 #include <pthread.h>
 
+/// Allocate a new future object.
+DMNSN_INTERNAL dmnsn_future *dmnsn_new_future(void);
+
+/// Set the total number of loop iterations.
+DMNSN_INTERNAL void dmnsn_future_set_total(dmnsn_future *future, size_t total);
+/// Increment the progress of a background task.
+DMNSN_INTERNAL void dmnsn_future_increment(dmnsn_future *future);
+/// Instantly complete the background teask.
+DMNSN_INTERNAL void dmnsn_future_finish(dmnsn_future *future);
+/// Set the number of worker threads.
+DMNSN_INTERNAL void dmnsn_future_set_nthreads(dmnsn_future *future, unsigned int nthreads);
+/// Notify completion of a worker thread.
+DMNSN_INTERNAL void dmnsn_future_finish_thread(dmnsn_future *future);
+
 /**
  * Thread callback type.
  * @param[in,out] ptr  An arbitrary pointer.
@@ -193,8 +212,7 @@ DMNSN_INTERNAL void dmnsn_once(pthread_once_t *once, dmnsn_once_fn *once_fn);
  * @param[out] key         The key to initialize.
  * @param[in]  destructor  An optional destructor callback.
  */
-DMNSN_INTERNAL void dmnsn_key_create(pthread_key_t *key,
-                                     dmnsn_callback_fn *destructor);
+DMNSN_INTERNAL void dmnsn_key_create(pthread_key_t *key, dmnsn_callback_fn *destructor);
 
 /**
  * Set a thread-specific pointer, bailing out on failure.
@@ -214,3 +232,5 @@ DMNSN_INTERNAL void dmnsn_key_delete(pthread_key_t key);
  * @param[in,out] thread  The thread to join.
  */
 DMNSN_INTERNAL void dmnsn_join_thread(pthread_t thread, void **retval);
+
+#endif // DMNSN_INTERNAL_CONCURRENCY_H
diff --git a/libdimension/future-internal.h b/libdimension/internal/future.h
index 7ca318a..047523e 100644
--- a/libdimension/future-internal.h
+++ b/libdimension/internal/future.h
@@ -23,22 +23,10 @@
  * Future object implementation.
  */
 
-#include <pthread.h>
-
-/// Allocate a new future object.
-DMNSN_INTERNAL dmnsn_future *dmnsn_new_future(void);
+#ifndef DMNSN_INTERNAL_FUTURE_H
+#define DMNSN_INTERNAL_FUTURE_H
 
-/// Set the total number of loop iterations.
-DMNSN_INTERNAL void dmnsn_future_set_total(dmnsn_future *future, size_t total);
-/// Increment the progress of a background task.
-DMNSN_INTERNAL void dmnsn_future_increment(dmnsn_future *future);
-/// Instantly complete the background teask.
-DMNSN_INTERNAL void dmnsn_future_finish(dmnsn_future *future);
-/// Set the number of worker threads.
-DMNSN_INTERNAL void dmnsn_future_set_nthreads(dmnsn_future *future,
-                                              unsigned int nthreads);
-/// Notify completion of a worker thread.
-DMNSN_INTERNAL void dmnsn_future_finish_thread(dmnsn_future *future);
+#include <pthread.h>
 
 struct dmnsn_future {
   size_t progress; ///< Completed loop iterations.
@@ -69,3 +57,5 @@ struct dmnsn_future {
   /// Condition variable for waiting for npaused == 0.
   pthread_cond_t resume_cond;
 };
+
+#endif // DMNSN_INTERNAL_FUTURE_H
diff --git a/libdimension/platform.h b/libdimension/internal/platform.h
index c03ae52..e037d21 100644
--- a/libdimension/platform.h
+++ b/libdimension/internal/platform.h
@@ -20,9 +20,14 @@
 
 /**
  * @file
- * Platform abstractions.
+ * Internal platform abstractions.
  */
 
+#ifndef DMNSN_INTERNAL_PLATFORM_H
+#define DMNSN_INTERNAL_PLATFORM_H
+
+#include "internal.h"
+#include "dimension/platform.h"
 #include <stdbool.h>
 #include <stddef.h>
 #include <stdio.h>
@@ -58,3 +63,5 @@ DMNSN_INTERNAL size_t dmnsn_ncpus(void);
  *                    times.
  */
 DMNSN_INTERNAL void dmnsn_get_times(dmnsn_timer *timer);
+
+#endif // DMNSN_INTERNAL_PLATFORM_H
diff --git a/libdimension/dimension/polynomial.h b/libdimension/internal/polynomial.h
index ab95d50..7af3e5a 100644
--- a/libdimension/dimension/polynomial.h
+++ b/libdimension/internal/polynomial.h
@@ -25,6 +25,7 @@
  * coefficient on x^i.  In general, we are only interested in positive roots.
  */
 
+#include "internal.h"
 #include <stddef.h>
 #include <stdio.h>
 
@@ -36,7 +37,7 @@
  * @param[in] degree  The degree of the polynomial.
  * @param[in] x       The value of the variable at which to evaluate.
  */
-DMNSN_INLINE double
+DMNSN_INTERNAL DMNSN_INLINE double
 dmnsn_polynomial_evaluate(const double poly[], size_t degree, double x)
 {
   double ret = poly[degree];
@@ -53,9 +54,8 @@ dmnsn_polynomial_evaluate(const double poly[], size_t degree, double x)
  * @param[in] degree  The degree of the polynomial.
  * @param[in] x       The value of the variable at which to evaluate.
  */
-DMNSN_INLINE double
-dmnsn_polynomial_evaluate_derivative(const double poly[], size_t degree,
-                                     double x)
+DMNSN_INTERNAL DMNSN_INLINE double
+dmnsn_polynomial_evaluate_derivative(const double poly[], size_t degree, double x)
 {
   double ret = poly[degree]*degree;
   size_t i;
@@ -73,7 +73,7 @@ dmnsn_polynomial_evaluate_derivative(const double poly[], size_t degree,
  *                     dimension \p degree.
  * @return The number of positive roots stored in \c x[].
  */
-size_t dmnsn_polynomial_solve(const double poly[], size_t degree, double x[]);
+DMNSN_INTERNAL size_t dmnsn_polynomial_solve(const double poly[], size_t degree, double x[]);
 
 /**
  * Output a polynomial.  The polynomial is printed as a function of x suitable
@@ -82,4 +82,4 @@ size_t dmnsn_polynomial_solve(const double poly[], size_t degree, double x[]);
  * @param[in]     poly    The coefficients of the polynomial to print.
  * @param[in]     degree  The degree of the polynomial.
  */
-void dmnsn_polynomial_print(FILE *file, const double poly[], size_t degree);
+DMNSN_INTERNAL void dmnsn_polynomial_print(FILE *file, const double poly[], size_t degree);
diff --git a/libdimension/profile.h b/libdimension/internal/profile.h
index 4817482..0f59b8f 100644
--- a/libdimension/profile.h
+++ b/libdimension/internal/profile.h
@@ -23,9 +23,38 @@
  * Built-in branch profiler.
  */
 
+#ifndef DMNSN_INTERNAL_H
+#error "Please include \"internal.h\" instead of this header directly."
+#endif
+
 #include <stdbool.h>
 
 /**
+ * @def dmnsn_likely
+ * Indicate that a test is likely to succeed.
+ * @param test  The test to perform.
+ * @return The truth value of \p test.
+ */
+/**
+ * @def dmnsn_unlikely
+ * Indicate that a test is unlikely to succeed.
+ * @param test  The test to perform.
+ * @return The truth value of \p test.
+ */
+#ifdef DMNSN_PROFILE
+  #define dmnsn_likely(test)                                      \
+    dmnsn_expect(!!(test), true, DMNSN_FUNC, __FILE__, __LINE__)
+  #define dmnsn_unlikely(test)                                    \
+    dmnsn_expect(!!(test), false, DMNSN_FUNC, __FILE__, __LINE__)
+#elif DMNSN_GNUC
+  #define dmnsn_likely(test)   __builtin_expect(!!(test), true)
+  #define dmnsn_unlikely(test) __builtin_expect(!!(test), false)
+#else
+  #define dmnsn_likely(test)   (!!(test))
+  #define dmnsn_unlikely(test) (!!(test))
+#endif
+
+/**
  * Record a test and its expected result.  Called by dmnsn_[un]likely();
  * don't call directly.
  * @param[in] result    The result of the test.
diff --git a/libdimension/prtree.h b/libdimension/internal/prtree.h
index 76a4051..7d85f44 100644
--- a/libdimension/prtree.h
+++ b/libdimension/internal/prtree.h
@@ -26,5 +26,13 @@
  * from B-trees.
  */
 
+#ifndef DMNSN_INTERNAL_PRTREE_H
+#define DMNSN_INTERNAL_PRTREE_H
+
+#include "internal.h"
+#include "internal/bvh.h"
+
 /// Create a PR-tree.
 DMNSN_INTERNAL dmnsn_bvh_node *dmnsn_new_prtree(const dmnsn_array *objects);
+
+#endif // DMNSN_INTERNAL_PRTREE_H
diff --git a/libdimension/rgba.h b/libdimension/internal/rgba.h
index b43ca3f..964d622 100644
--- a/libdimension/rgba.h
+++ b/libdimension/internal/rgba.h
@@ -23,6 +23,11 @@
  * RGBA canvas optimizer interface, used by image optimizers.
  */
 
+#ifndef DMNSN_INTERNAL_RGBA_H
+#define DMNSN_INTERNAL_RGBA_H
+
+#include "internal.h"
+#include "dimension/canvas.h"
 #include <stdint.h>
 
 /// RGBA8 optimizer type.
@@ -46,3 +51,5 @@ DMNSN_INTERNAL void dmnsn_rgba16_optimize_canvas(dmnsn_pool *pool, dmnsn_canvas
 DMNSN_INTERNAL void dmnsn_rgba8_optimizer_fn(dmnsn_canvas_optimizer *optimizer, const dmnsn_canvas *canvas, size_t x, size_t y);
 /// RGBA16 optimizer callback.
 DMNSN_INTERNAL void dmnsn_rgba16_optimizer_fn(dmnsn_canvas_optimizer *optimizer, const dmnsn_canvas *canvas, size_t x, size_t y);
+
+#endif // DMNSN_INTERNAL_RGBA_H
diff --git a/libdimension/geometry.c b/libdimension/math/matrix.c
index 4a16a3c..25590d8 100644
--- a/libdimension/geometry.c
+++ b/libdimension/math/matrix.c
@@ -20,10 +20,11 @@
 
 /**
  * @file
- * Geometrical function implementations.
+ * Matrix function implementations.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/math.h"
 #include <math.h>
 
 // Identity matrix
@@ -353,8 +354,8 @@ dmnsn_matrix_mul(dmnsn_matrix lhs, dmnsn_matrix rhs)
 }
 
 // Give an axis-aligned box that contains the given box transformed by `lhs'
-dmnsn_bounding_box
-dmnsn_transform_bounding_box(dmnsn_matrix trans, dmnsn_bounding_box box)
+dmnsn_aabb
+dmnsn_transform_aabb(dmnsn_matrix trans, dmnsn_aabb box)
 {
   // Infinite/zero bounding box support
   if (isinf(box.min.x)) {
@@ -369,7 +370,7 @@ dmnsn_transform_bounding_box(dmnsn_matrix trans, dmnsn_bounding_box box)
   }
 
   dmnsn_vector Mt = dmnsn_matrix_column(trans, 3);
-  dmnsn_bounding_box ret = { Mt, Mt };
+  dmnsn_aabb ret = { Mt, Mt };
 
   dmnsn_vector Mz = dmnsn_matrix_column(trans, 2);
   ret.min = dmnsn_vector_add(ret.min, dmnsn_vector_mul(box.min.z, Mz));
diff --git a/libdimension/polynomial.c b/libdimension/math/polynomial.c
index 609364a..09e9603 100644
--- a/libdimension/polynomial.c
+++ b/libdimension/math/polynomial.c
@@ -24,7 +24,9 @@
  * Akritas, Johnson, Krandick, and Rouillier and Zimmerman.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/polynomial.h"
+#include "dimension/math.h"
 #include <math.h>
 
 /// Get the real degree of a polynomial, ignoring leading zeros.
diff --git a/libdimension/camera.c b/libdimension/model/camera.c
index 550b41c..8d05c28 100644
--- a/libdimension/camera.c
+++ b/libdimension/model/camera.c
@@ -23,7 +23,8 @@
  * Cameras.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 // Allocate a new dummy camera
@@ -43,9 +44,9 @@ dmnsn_init_camera(dmnsn_camera *camera)
 }
 
 // Invoke the camera ray function
-dmnsn_line
+dmnsn_ray
 dmnsn_camera_ray(const dmnsn_camera *camera, double x, double y)
 {
-  dmnsn_line ray = camera->ray_fn(camera, x, y);
-  return dmnsn_transform_line(camera->trans, ray);
+  dmnsn_ray ray = camera->ray_fn(camera, x, y);
+  return dmnsn_transform_ray(camera->trans, ray);
 }
diff --git a/libdimension/perspective.c b/libdimension/model/cameras/perspective.c
index 518e52e..6f1a9fc 100644
--- a/libdimension/perspective.c
+++ b/libdimension/model/cameras/perspective.c
@@ -23,14 +23,14 @@
  * Perspective cameras.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 /// Perspective camera ray callback.
-static dmnsn_line
+static dmnsn_ray
 dmnsn_perspective_camera_ray_fn(const dmnsn_camera *camera, double x, double y)
 {
-  dmnsn_line l = dmnsn_new_line(
+  dmnsn_ray l = dmnsn_new_ray(
     dmnsn_zero,
     dmnsn_new_vector(x - 0.5, y - 0.5, 1.0)
   );
diff --git a/libdimension/finish.c b/libdimension/model/finish.c
index 13adefb..a5dbaf9 100644
--- a/libdimension/finish.c
+++ b/libdimension/model/finish.c
@@ -23,7 +23,8 @@
  * Finishes.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/model.h"
 
 dmnsn_ambient *
 dmnsn_new_ambient(dmnsn_pool *pool, dmnsn_color ambient_light)
diff --git a/libdimension/lambertian.c b/libdimension/model/finishes/lambertian.c
index cced4a7..9e041ed 100644
--- a/libdimension/lambertian.c
+++ b/libdimension/model/finishes/lambertian.c
@@ -23,7 +23,7 @@
  * Diffuse finish.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 #include <math.h>
 #include <stdlib.h>
 
diff --git a/libdimension/phong.c b/libdimension/model/finishes/phong.c
index 10e5b1a..6b1b1e1 100644
--- a/libdimension/phong.c
+++ b/libdimension/model/finishes/phong.c
@@ -23,7 +23,7 @@
  * Phong highlights.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 /// Phone specular type.
diff --git a/libdimension/reflection.c b/libdimension/model/finishes/reflection.c
index 6765bef..afaae8f 100644
--- a/libdimension/reflection.c
+++ b/libdimension/model/finishes/reflection.c
@@ -23,7 +23,7 @@
  * Reflective finish.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 #include <math.h>
 #include <stdlib.h>
 
diff --git a/libdimension/interior.c b/libdimension/model/interior.c
index d725b6f..82db3ea 100644
--- a/libdimension/interior.c
+++ b/libdimension/model/interior.c
@@ -23,7 +23,8 @@
  * Interiors.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 // Allocate an interior
diff --git a/libdimension/light.c b/libdimension/model/light.c
index 14917de..6143b3c 100644
--- a/libdimension/light.c
+++ b/libdimension/model/light.c
@@ -23,7 +23,8 @@
  * Lights.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 // Allocate a new dummy light
diff --git a/libdimension/point_light.c b/libdimension/model/lights/point_light.c
index b8c2cfb..08aae67 100644
--- a/libdimension/point_light.c
+++ b/libdimension/model/lights/point_light.c
@@ -23,7 +23,7 @@
  * Point lights.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 /// Point light type.
diff --git a/libdimension/object.c b/libdimension/model/object.c
index a67b8cb..0473f54 100644
--- a/libdimension/object.c
+++ b/libdimension/model/object.c
@@ -23,7 +23,8 @@
  * Objects.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 // Allocate a dummy object
@@ -93,7 +94,7 @@ dmnsn_object_precompute_recursive(dmnsn_object *object, dmnsn_matrix pigment_tra
   object->pigment_trans = pigment_trans;
   object->trans_inv = dmnsn_matrix_inverse(total_trans);
   if (vtable->bounding_fn) {
-    object->bounding_box = vtable->bounding_fn(object, total_trans);
+    object->aabb = vtable->bounding_fn(object, total_trans);
   }
   if (vtable->precompute_fn) {
     vtable->precompute_fn(object);
diff --git a/libdimension/cone.c b/libdimension/model/objects/cone.c
index 1e95b0a..26e59ca 100644
--- a/libdimension/cone.c
+++ b/libdimension/model/objects/cone.c
@@ -23,7 +23,8 @@
  * Cones/cylinders.
  */
 
-#include "dimension.h"
+#include "internal/polynomial.h"
+#include "dimension/model.h"
 #include <math.h>
 
 /// Cone type.
@@ -34,7 +35,7 @@ typedef struct dmnsn_cone {
 
 /// Intersection callback for a cone.
 static bool
-dmnsn_cone_intersection_fn(const dmnsn_object *object, dmnsn_line l,
+dmnsn_cone_intersection_fn(const dmnsn_object *object, dmnsn_ray l,
                            dmnsn_intersection *intersection)
 {
   const dmnsn_cone *cone = (const dmnsn_cone *)object;
@@ -55,14 +56,14 @@ dmnsn_cone_intersection_fn(const dmnsn_object *object, dmnsn_line l,
     dmnsn_vector p;
     if (n == 2) {
       t = dmnsn_min(t, x[1]);
-      p = dmnsn_line_point(l, t);
+      p = dmnsn_ray_point(l, t);
 
       if (p.y <= -1.0 || p.y >= 1.0) {
         t = dmnsn_max(x[0], x[1]);
-        p = dmnsn_line_point(l, t);
+        p = dmnsn_ray_point(l, t);
       }
     } else {
-      p = dmnsn_line_point(l, t);
+      p = dmnsn_ray_point(l, t);
     }
 
     if (t >= 0.0 && p.y >= -1.0 && p.y <= 1.0) {
@@ -89,14 +90,14 @@ dmnsn_cone_inside_fn(const dmnsn_object *object, dmnsn_vector point)
 }
 
 /// Cone bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_cone_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   const dmnsn_cone *cone = (const dmnsn_cone *)object;
 
   double rmax = dmnsn_max(cone->r1, cone->r2);
-  dmnsn_bounding_box box = dmnsn_symmetric_bounding_box(dmnsn_new_vector(rmax, 1.0, rmax));
-  return dmnsn_transform_bounding_box(trans, box);
+  dmnsn_aabb box = dmnsn_symmetric_aabb(dmnsn_new_vector(rmax, 1.0, rmax));
+  return dmnsn_transform_aabb(trans, box);
 }
 
 /// Cone vtable.
@@ -114,14 +115,14 @@ typedef struct dmnsn_cone_cap {
 
 /// Cone cap intersection function.
 static bool
-dmnsn_cone_cap_intersection_fn(const dmnsn_object *object, dmnsn_line l,
+dmnsn_cone_cap_intersection_fn(const dmnsn_object *object, dmnsn_ray l,
                                dmnsn_intersection *intersection)
 {
   if (l.n.y != 0.0) {
     const dmnsn_cone_cap *cap = (const dmnsn_cone_cap *)object;
     double r = cap->r;
     double t = -l.x0.y/l.n.y;
-    dmnsn_vector p = dmnsn_line_point(l, t);
+    dmnsn_vector p = dmnsn_ray_point(l, t);
     if (t >= 0.0 && p.x*p.x + p.z*p.z <= r*r) {
       intersection->t      = t;
       intersection->normal = dmnsn_new_vector(0.0, -1.0, 0.0);
@@ -140,12 +141,12 @@ dmnsn_cone_cap_inside_fn(const dmnsn_object *object, dmnsn_vector point)
 }
 
 /// Cone cap bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_cone_cap_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   const dmnsn_cone_cap *cap = (const dmnsn_cone_cap *)object;
-  dmnsn_bounding_box box = dmnsn_symmetric_bounding_box(dmnsn_new_vector(cap->r, 0.0, cap->r));
-  return dmnsn_transform_bounding_box(trans, box);
+  dmnsn_aabb box = dmnsn_symmetric_aabb(dmnsn_new_vector(cap->r, 0.0, cap->r));
+  return dmnsn_transform_aabb(trans, box);
 }
 
 /// Cone cap vtable.
diff --git a/libdimension/csg.c b/libdimension/model/objects/csg.c
index e31cb0e..15008c0 100644
--- a/libdimension/csg.c
+++ b/libdimension/model/objects/csg.c
@@ -23,7 +23,9 @@
  * Constructive solid geometry.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/bvh.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 ////////////
@@ -38,11 +40,11 @@ typedef struct {
 /// CSG union intersection callback.
 static bool
 dmnsn_csg_union_intersection_fn(const dmnsn_object *object,
-                                dmnsn_line line,
+                                dmnsn_ray ray,
                                 dmnsn_intersection *intersection)
 {
   const dmnsn_csg_union *csg = (const dmnsn_csg_union *)object;
-  return dmnsn_bvh_intersection(csg->bvh, line, intersection, true);
+  return dmnsn_bvh_intersection(csg->bvh, ray, intersection, true);
 }
 
 /// CSG union inside callback.
@@ -62,7 +64,7 @@ dmnsn_csg_union_precompute_fn(dmnsn_object *object)
 
   dmnsn_bvh *bvh = dmnsn_new_bvh(csg->object.children, DMNSN_BVH_PRTREE);
   csg->bvh = bvh;
-  csg->object.bounding_box = dmnsn_bvh_bounding_box(bvh);
+  csg->object.aabb = dmnsn_bvh_aabb(bvh);
 }
 
 /// CSG union vtable.
@@ -100,16 +102,16 @@ dmnsn_new_csg_union(dmnsn_pool *pool, dmnsn_array *objects)
 /**
  * Generic CSG intersection callback.
  * @param[in]  csg           The CSG object.
- * @param[in]  line          The intersection ray.
+ * @param[in]  ray           The intersection ray.
  * @param[out] intersection  The intersection data.
  * @param[in]  inside1       Whether the first object is allowed inside the
  *                           second object.
  * @param[in]  inside2       Whether the second object is allowed inside the
  *                           first object.
- * @return Whether \p line intersected \p csg.
+ * @return Whether \p ray intersected \p csg.
  */
 static bool
-dmnsn_csg_intersection_fn(const dmnsn_object *csg, dmnsn_line line,
+dmnsn_csg_intersection_fn(const dmnsn_object *csg, dmnsn_ray ray,
                           dmnsn_intersection *intersection,
                           bool inside1, bool inside2)
 {
@@ -117,21 +119,21 @@ dmnsn_csg_intersection_fn(const dmnsn_object *csg, dmnsn_line line,
   const dmnsn_object *B = *(dmnsn_object **)dmnsn_array_last(csg->children);
 
   dmnsn_intersection i1, i2;
-  bool is_i1 = dmnsn_object_intersection(A, line, &i1);
-  bool is_i2 = dmnsn_object_intersection(B, line, &i2);
+  bool is_i1 = dmnsn_object_intersection(A, ray, &i1);
+  bool is_i2 = dmnsn_object_intersection(B, ray, &i2);
 
   double oldt = 0.0;
   while (is_i1) {
-    i1.ray = line;
+    i1.ray = ray;
     i1.t += oldt;
     oldt = i1.t + dmnsn_epsilon;
 
-    dmnsn_vector point = dmnsn_line_point(i1.ray, i1.t);
+    dmnsn_vector point = dmnsn_ray_point(i1.ray, i1.t);
     if (inside2 ^ dmnsn_object_inside(B, point)) {
-      dmnsn_line newline = line;
-      newline.x0 = dmnsn_line_point(line, i1.t);
-      newline    = dmnsn_line_add_epsilon(newline);
-      is_i1 = dmnsn_object_intersection(A, newline, &i1);
+      dmnsn_ray newray = ray;
+      newray.x0 = dmnsn_ray_point(ray, i1.t);
+      newray    = dmnsn_ray_add_epsilon(newray);
+      is_i1 = dmnsn_object_intersection(A, newray, &i1);
     } else {
       break;
     }
@@ -139,16 +141,16 @@ dmnsn_csg_intersection_fn(const dmnsn_object *csg, dmnsn_line line,
 
   oldt = 0.0;
   while (is_i2) {
-    i2.ray = line;
+    i2.ray = ray;
     i2.t += oldt;
     oldt = i2.t + dmnsn_epsilon;
 
-    dmnsn_vector point = dmnsn_line_point(i2.ray, i2.t);
+    dmnsn_vector point = dmnsn_ray_point(i2.ray, i2.t);
     if (inside1 ^ dmnsn_object_inside(A, point)) {
-      dmnsn_line newline = line;
-      newline.x0 = dmnsn_line_point(line, i2.t);
-      newline    = dmnsn_line_add_epsilon(newline);
-      is_i2 = dmnsn_object_intersection(B, newline, &i2);
+      dmnsn_ray newray = ray;
+      newray.x0 = dmnsn_ray_point(ray, i2.t);
+      newray    = dmnsn_ray_add_epsilon(newray);
+      is_i2 = dmnsn_object_intersection(B, newray, &i2);
     } else {
       break;
     }
@@ -178,10 +180,10 @@ dmnsn_csg_intersection_fn(const dmnsn_object *csg, dmnsn_line line,
 /// CSG intersection intersection callback.
 static bool
 dmnsn_csg_intersection_intersection_fn(const dmnsn_object *csg,
-                                       dmnsn_line line,
+                                       dmnsn_ray ray,
                                        dmnsn_intersection *intersection)
 {
-  return dmnsn_csg_intersection_fn(csg, line, intersection, true, true);
+  return dmnsn_csg_intersection_fn(csg, ray, intersection, true, true);
 }
 
 /// CSG intersection inside callback.
@@ -201,8 +203,8 @@ dmnsn_csg_intersection_precompute_fn(dmnsn_object *csg)
   dmnsn_object *B = *(dmnsn_object **)dmnsn_array_last(csg->children);
 
   csg->trans_inv = dmnsn_identity_matrix();
-  csg->bounding_box.min = dmnsn_vector_max(A->bounding_box.min, B->bounding_box.min);
-  csg->bounding_box.max = dmnsn_vector_min(A->bounding_box.max, B->bounding_box.max);
+  csg->aabb.min = dmnsn_vector_max(A->aabb.min, B->aabb.min);
+  csg->aabb.max = dmnsn_vector_min(A->aabb.max, B->aabb.max);
 }
 
 /// CSG intersection vtable.
@@ -232,10 +234,10 @@ dmnsn_new_csg_intersection(dmnsn_pool *pool, dmnsn_object *A, dmnsn_object *B)
 /// CSG difference intersection callback.
 static bool
 dmnsn_csg_difference_intersection_fn(const dmnsn_object *csg,
-                                     dmnsn_line line,
+                                     dmnsn_ray ray,
                                      dmnsn_intersection *intersection)
 {
-  return dmnsn_csg_intersection_fn(csg, line, intersection, true, false);
+  return dmnsn_csg_intersection_fn(csg, ray, intersection, true, false);
 }
 
 /// CSG difference inside callback.
@@ -254,7 +256,7 @@ dmnsn_csg_difference_precompute_fn(dmnsn_object *csg)
   dmnsn_object *A = *(dmnsn_object **)dmnsn_array_first(csg->children);
 
   csg->trans_inv = dmnsn_identity_matrix();
-  csg->bounding_box = A->bounding_box;
+  csg->aabb = A->aabb;
 }
 
 /// CSG difference vtable.
@@ -284,10 +286,10 @@ dmnsn_new_csg_difference(dmnsn_pool *pool, dmnsn_object *A, dmnsn_object *B)
 /// CSG merge intersection callback.
 static bool
 dmnsn_csg_merge_intersection_fn(const dmnsn_object *csg,
-                                dmnsn_line line,
+                                dmnsn_ray ray,
                                 dmnsn_intersection *intersection)
 {
-  return dmnsn_csg_intersection_fn(csg, line, intersection, false, false);
+  return dmnsn_csg_intersection_fn(csg, ray, intersection, false, false);
 }
 
 /// CSG merge inside callback.
@@ -307,8 +309,8 @@ dmnsn_csg_merge_precompute_fn(dmnsn_object *csg)
   dmnsn_object *B = *(dmnsn_object **)dmnsn_array_last(csg->children);
 
   csg->trans_inv = dmnsn_identity_matrix();
-  csg->bounding_box.min = dmnsn_vector_min(A->bounding_box.min, B->bounding_box.min);
-  csg->bounding_box.max = dmnsn_vector_max(A->bounding_box.max, B->bounding_box.max);
+  csg->aabb.min = dmnsn_vector_min(A->aabb.min, B->aabb.min);
+  csg->aabb.max = dmnsn_vector_max(A->aabb.max, B->aabb.max);
 }
 
 /// CSG merge vtable.
diff --git a/libdimension/cube.c b/libdimension/model/objects/cube.c
index c455467..7d6fe0f 100644
--- a/libdimension/cube.c
+++ b/libdimension/model/objects/cube.c
@@ -23,21 +23,21 @@
  * Cubes.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 #include <math.h>
 
 /// Intersection callback for a cube.
 static bool
-dmnsn_cube_intersection_fn(const dmnsn_object *cube, dmnsn_line line,
+dmnsn_cube_intersection_fn(const dmnsn_object *cube, dmnsn_ray ray,
                            dmnsn_intersection *intersection)
 {
-  // Clip the given line against the X, Y, and Z slabs
+  // Clip the given ray against the X, Y, and Z slabs
 
   dmnsn_vector nmin, nmax;
   double tmin, tmax;
 
-  double tx1 = (-1.0 - line.x0.x)/line.n.x;
-  double tx2 = (+1.0 - line.x0.x)/line.n.x;
+  double tx1 = (-1.0 - ray.x0.x)/ray.n.x;
+  double tx2 = (+1.0 - ray.x0.x)/ray.n.x;
 
   if (tx1 < tx2) {
     tmin = tx1;
@@ -54,8 +54,8 @@ dmnsn_cube_intersection_fn(const dmnsn_object *cube, dmnsn_line line,
   if (tmin > tmax)
     return false;
 
-  double ty1 = (-1.0 - line.x0.y)/line.n.y;
-  double ty2 = (+1.0 - line.x0.y)/line.n.y;
+  double ty1 = (-1.0 - ray.x0.y)/ray.n.y;
+  double ty2 = (+1.0 - ray.x0.y)/ray.n.y;
 
   if (ty1 < ty2) {
     if (ty1 > tmin) {
@@ -80,8 +80,8 @@ dmnsn_cube_intersection_fn(const dmnsn_object *cube, dmnsn_line line,
   if (tmin > tmax)
     return false;
 
-  double tz1 = (-1.0 - line.x0.z)/line.n.z;
-  double tz2 = (+1.0 - line.x0.z)/line.n.z;
+  double tz1 = (-1.0 - ray.x0.z)/ray.n.z;
+  double tz2 = (+1.0 - ray.x0.z)/ray.n.z;
 
   if (tz1 < tz2) {
     if (tz1 > tmin) {
@@ -130,11 +130,11 @@ dmnsn_cube_inside_fn(const dmnsn_object *cube, dmnsn_vector point)
 }
 
 /// Boundary callback for a cube.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_cube_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
-  dmnsn_bounding_box box = dmnsn_symmetric_bounding_box(dmnsn_new_vector(1.0, 1.0, 1.0));
-  return dmnsn_transform_bounding_box(trans, box);
+  dmnsn_aabb box = dmnsn_symmetric_aabb(dmnsn_new_vector(1.0, 1.0, 1.0));
+  return dmnsn_transform_aabb(trans, box);
 }
 
 /// Cube vtable.
diff --git a/libdimension/plane.c b/libdimension/model/objects/plane.c
index 5428485..b34d8aa 100644
--- a/libdimension/plane.c
+++ b/libdimension/model/objects/plane.c
@@ -23,7 +23,7 @@
  * Planes.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 #include <math.h>
 #include <stdlib.h>
 
@@ -33,17 +33,17 @@ typedef struct {
   dmnsn_vector normal;
 } dmnsn_plane;
 
-/// Returns the closest intersection of `line' with `plane'.
+/// Returns the closest intersection of `ray' with `plane'.
 static bool
-dmnsn_plane_intersection_fn(const dmnsn_object *object, dmnsn_line line,
+dmnsn_plane_intersection_fn(const dmnsn_object *object, dmnsn_ray ray,
                             dmnsn_intersection *intersection)
 {
   const dmnsn_plane *plane = (const dmnsn_plane *)object;
   dmnsn_vector normal = plane->normal;
 
-  double den = dmnsn_vector_dot(line.n, normal);
+  double den = dmnsn_vector_dot(ray.n, normal);
   if (den != 0.0) {
-    double t = -dmnsn_vector_dot(line.x0, normal)/den;
+    double t = -dmnsn_vector_dot(ray.x0, normal)/den;
     if (t >= 0.0) {
       intersection->t      = t;
       intersection->normal = normal;
@@ -62,10 +62,10 @@ dmnsn_plane_inside_fn(const dmnsn_object *object, dmnsn_vector point)
 }
 
 /// Plane bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_plane_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
-  return dmnsn_infinite_bounding_box();
+  return dmnsn_infinite_aabb();
 }
 
 /// Plane vtable.
diff --git a/libdimension/sphere.c b/libdimension/model/objects/sphere.c
index d8fae05..e1ca784 100644
--- a/libdimension/sphere.c
+++ b/libdimension/model/objects/sphere.c
@@ -23,11 +23,13 @@
  * Spheres.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/polynomial.h"
+#include "dimension/model.h"
 
 /// Sphere intersection callback.
 static bool
-dmnsn_sphere_intersection_fn(const dmnsn_object *sphere, dmnsn_line l,
+dmnsn_sphere_intersection_fn(const dmnsn_object *sphere, dmnsn_ray l,
                              dmnsn_intersection *intersection)
 {
   // Solve (x0 + nx*t)^2 + (y0 + ny*t)^2 + (z0 + nz*t)^2 == 1
@@ -48,7 +50,7 @@ dmnsn_sphere_intersection_fn(const dmnsn_object *sphere, dmnsn_line l,
   }
 
   intersection->t = t;
-  intersection->normal = dmnsn_line_point(l, t);
+  intersection->normal = dmnsn_ray_point(l, t);
   return true;
 }
 
@@ -71,14 +73,14 @@ dmnsn_implicit_dot(double row[4])
 }
 
 /// Sphere bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_sphere_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   // Get a tight bound using the quadric representation of a sphere.  For
   // details, see
   // http://tavianator.com/2014/06/exact-bounding-boxes-for-spheres-ellipsoids
 
-  dmnsn_bounding_box box;
+  dmnsn_aabb box;
 
   double cx = trans.n[0][3];
   double dx = sqrt(dmnsn_implicit_dot(trans.n[0]));
diff --git a/libdimension/torus.c b/libdimension/model/objects/torus.c
index 30acbd8..b4baebd 100644
--- a/libdimension/torus.c
+++ b/libdimension/model/objects/torus.c
@@ -23,7 +23,8 @@
  * Torii.  A special case of a quartic.
  */
 
-#include "dimension.h"
+#include "internal/polynomial.h"
+#include "dimension/model.h"
 
 /// Torus type.
 typedef struct {
@@ -33,7 +34,7 @@ typedef struct {
 
 /// Bound the torus in a cylindrical shell.
 static inline bool
-dmnsn_torus_bound_intersection(const dmnsn_torus *torus, dmnsn_line l)
+dmnsn_torus_bound_intersection(const dmnsn_torus *torus, dmnsn_ray l)
 {
   double R = torus->major, r = torus->minor;
   double rmax = R + r, rmin = R - r;
@@ -42,8 +43,8 @@ dmnsn_torus_bound_intersection(const dmnsn_torus *torus, dmnsn_line l)
   // Try the caps first
   double tlower = (-r - l.x0.y)/l.n.y;
   double tupper = (+r - l.x0.y)/l.n.y;
-  dmnsn_vector lower = dmnsn_line_point(l, tlower);
-  dmnsn_vector upper = dmnsn_line_point(l, tupper);
+  dmnsn_vector lower = dmnsn_ray_point(l, tlower);
+  dmnsn_vector upper = dmnsn_ray_point(l, tupper);
   double ldist2 = lower.x*lower.x + lower.z*lower.z;
   double udist2 = upper.x*upper.x + upper.z*upper.z;
   if ((ldist2 < rmin2 || ldist2 > rmax2) && (udist2 < rmin2 || udist2 > rmax2)) {
@@ -61,7 +62,7 @@ dmnsn_torus_bound_intersection(const dmnsn_torus *torus, dmnsn_line l)
 
     size_t i;
     for (i = 0; i < n; ++i) {
-      dmnsn_vector p = dmnsn_line_point(l, x[i]);
+      dmnsn_vector p = dmnsn_ray_point(l, x[i]);
       if (p.y >= -r && p.y <= r)
         break;
     }
@@ -77,7 +78,7 @@ dmnsn_torus_bound_intersection(const dmnsn_torus *torus, dmnsn_line l)
 
 /// Torus intersection callback.
 static bool
-dmnsn_torus_intersection_fn(const dmnsn_object *object, dmnsn_line l,
+dmnsn_torus_intersection_fn(const dmnsn_object *object, dmnsn_ray l,
                             dmnsn_intersection *intersection)
 {
   const dmnsn_torus *torus = (const dmnsn_torus *)object;
@@ -119,7 +120,7 @@ dmnsn_torus_intersection_fn(const dmnsn_object *object, dmnsn_line l,
     return false;
   }
 
-  dmnsn_vector p = dmnsn_line_point(l, t);
+  dmnsn_vector p = dmnsn_ray_point(l, t);
   dmnsn_vector center = dmnsn_vector_mul(
     R,
     dmnsn_vector_normalized(dmnsn_new_vector(p.x, 0.0, p.z))
@@ -141,14 +142,14 @@ dmnsn_torus_inside_fn(const dmnsn_object *object, dmnsn_vector point)
 }
 
 /// Torus bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_torus_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   const dmnsn_torus *torus = (const dmnsn_torus *)object;
 
   double extent = torus->major + torus->minor;
-  dmnsn_bounding_box box = dmnsn_symmetric_bounding_box(dmnsn_new_vector(extent, torus->minor, extent));
-  return dmnsn_transform_bounding_box(trans, box);
+  dmnsn_aabb box = dmnsn_symmetric_aabb(dmnsn_new_vector(extent, torus->minor, extent));
+  return dmnsn_transform_aabb(trans, box);
 }
 
 /// Torus vtable.
diff --git a/libdimension/triangle.c b/libdimension/model/objects/triangle.c
index fdd2b96..5af3301 100644
--- a/libdimension/triangle.c
+++ b/libdimension/model/objects/triangle.c
@@ -25,11 +25,12 @@
  * for a description of the intersection algorithm.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/model.h"
 
 /// Optimized ray/triangle intersection test.
 static inline bool
-dmnsn_ray_triangle_intersection(dmnsn_line l, double *t, double *u, double *v)
+dmnsn_ray_triangle_intersection(dmnsn_ray l, double *t, double *u, double *v)
 {
   // See the change of basis in dmnsn_triangle_basis()
   *t = -l.x0.z/l.n.z;
@@ -40,7 +41,7 @@ dmnsn_ray_triangle_intersection(dmnsn_line l, double *t, double *u, double *v)
 
 /// Triangle intersection callback.
 DMNSN_HOT static bool
-dmnsn_triangle_intersection_fn(const dmnsn_object *object, dmnsn_line l,
+dmnsn_triangle_intersection_fn(const dmnsn_object *object, dmnsn_ray l,
                                dmnsn_intersection *intersection)
 {
   double t, u, v;
@@ -61,16 +62,16 @@ dmnsn_triangle_inside_fn(const dmnsn_object *object, dmnsn_vector point)
 }
 
 /// Triangle bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_triangle_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   dmnsn_vector a = dmnsn_transform_point(trans, dmnsn_zero);
   dmnsn_vector b = dmnsn_transform_point(trans, dmnsn_x);
   dmnsn_vector c = dmnsn_transform_point(trans, dmnsn_y);
 
-  dmnsn_bounding_box box = dmnsn_new_bounding_box(a, a);
-  box = dmnsn_bounding_box_swallow(box, b);
-  box = dmnsn_bounding_box_swallow(box, c);
+  dmnsn_aabb box = dmnsn_new_aabb(a, a);
+  box = dmnsn_aabb_swallow(box, b);
+  box = dmnsn_aabb_swallow(box, c);
   return box;
 }
 
@@ -89,7 +90,7 @@ typedef struct {
 
 /// Smooth triangle intersection callback.
 DMNSN_HOT static bool
-dmnsn_smooth_triangle_intersection_fn(const dmnsn_object *object, dmnsn_line l,
+dmnsn_smooth_triangle_intersection_fn(const dmnsn_object *object, dmnsn_ray l,
                                       dmnsn_intersection *intersection)
 {
   const dmnsn_smooth_triangle *triangle = (const dmnsn_smooth_triangle *)object;
diff --git a/libdimension/triangle_fan.c b/libdimension/model/objects/triangle_fan.c
index 9940614..93768a9 100644
--- a/libdimension/triangle_fan.c
+++ b/libdimension/model/objects/triangle_fan.c
@@ -25,7 +25,8 @@
  * for a description of the intersection algorithm.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/model.h"
 
 /// Triangle fan type.
 typedef struct {
@@ -56,11 +57,11 @@ dmnsn_change_normal_basis(const double coeffs[6], dmnsn_vector n)
   );
 }
 
-/// Change basis from one triangle to the next for a line
-static inline dmnsn_line
-dmnsn_change_line_basis(const double coeffs[6], dmnsn_line l)
+/// Change basis from one triangle to the next for a ray
+static inline dmnsn_ray
+dmnsn_change_ray_basis(const double coeffs[6], dmnsn_ray l)
 {
-  return dmnsn_new_line(dmnsn_change_basis(coeffs, l.x0), dmnsn_change_basis(coeffs, l.n));
+  return dmnsn_new_ray(dmnsn_change_basis(coeffs, l.x0), dmnsn_change_basis(coeffs, l.n));
 }
 
 /// Store the compressed incremental matrix.
@@ -97,7 +98,7 @@ dmnsn_triangle_basis(dmnsn_vector a, dmnsn_vector ab, dmnsn_vector ac)
 
 /// Optimized ray/triangle intersection test.
 static inline bool
-dmnsn_ray_triangle_intersection(dmnsn_line l, double *t, double *u, double *v)
+dmnsn_ray_triangle_intersection(dmnsn_ray l, double *t, double *u, double *v)
 {
   *t = -l.x0.z/l.n.z;
   *u = l.x0.x + (*t)*l.n.x;
@@ -107,7 +108,7 @@ dmnsn_ray_triangle_intersection(dmnsn_line l, double *t, double *u, double *v)
 
 /// Triangle fan intersection callback.
 DMNSN_HOT static bool
-dmnsn_triangle_fan_intersection_fn(const dmnsn_object *object, dmnsn_line l, dmnsn_intersection *intersection)
+dmnsn_triangle_fan_intersection_fn(const dmnsn_object *object, dmnsn_ray l, dmnsn_intersection *intersection)
 {
   const dmnsn_triangle_fan *fan = (const dmnsn_triangle_fan *)object;
 
@@ -123,7 +124,7 @@ dmnsn_triangle_fan_intersection_fn(const dmnsn_object *object, dmnsn_line l, dmn
 
   for (size_t i = 0; i < fan->ncoeffs; ++i) {
     const double *coeffs = fan->coeffs[i];
-    l = dmnsn_change_line_basis(coeffs, l);
+    l = dmnsn_change_ray_basis(coeffs, l);
     normal = dmnsn_change_normal_basis(coeffs, normal);
 
     if (dmnsn_ray_triangle_intersection(l, &t, &u, &v) && t < best_t) {
@@ -149,33 +150,33 @@ dmnsn_triangle_fan_inside_fn(const dmnsn_object *object, dmnsn_vector point)
 }
 
 /// Computes the bounding box for the first triangle
-static inline dmnsn_bounding_box
+static inline dmnsn_aabb
 dmnsn_bound_first_triangle(dmnsn_matrix trans)
 {
   dmnsn_vector a = dmnsn_transform_point(trans, dmnsn_zero);
   dmnsn_vector b = dmnsn_transform_point(trans, dmnsn_x);
   dmnsn_vector c = dmnsn_transform_point(trans, dmnsn_y);
 
-  dmnsn_bounding_box box = dmnsn_new_bounding_box(a, a);
-  box = dmnsn_bounding_box_swallow(box, b);
-  box = dmnsn_bounding_box_swallow(box, c);
+  dmnsn_aabb box = dmnsn_new_aabb(a, a);
+  box = dmnsn_aabb_swallow(box, b);
+  box = dmnsn_aabb_swallow(box, c);
 
   return box;
 }
 
 /// Triangle fan bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_triangle_fan_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   const dmnsn_triangle_fan *fan = (const dmnsn_triangle_fan *)object;
 
-  dmnsn_bounding_box box = dmnsn_bound_first_triangle(trans);
+  dmnsn_aabb box = dmnsn_bound_first_triangle(trans);
 
   for (size_t i = 0; i < fan->ncoeffs; ++i) {
     dmnsn_matrix incremental = dmnsn_decompress_coeffs(fan->coeffs[i]);
     trans = dmnsn_matrix_mul(trans, dmnsn_matrix_inverse(incremental));
     dmnsn_vector vertex = dmnsn_transform_point(trans, dmnsn_y);
-    box = dmnsn_bounding_box_swallow(box, vertex);
+    box = dmnsn_aabb_swallow(box, vertex);
   }
 
   return box;
@@ -198,7 +199,7 @@ typedef struct dmnsn_smooth_triangle_fan {
 
 /// Smooth triangle fan intersection callback.
 DMNSN_HOT static bool
-dmnsn_smooth_triangle_fan_intersection_fn(const dmnsn_object *object, dmnsn_line l, dmnsn_intersection *intersection)
+dmnsn_smooth_triangle_fan_intersection_fn(const dmnsn_object *object, dmnsn_ray l, dmnsn_intersection *intersection)
 {
   const dmnsn_smooth_triangle_fan *fan = (const dmnsn_smooth_triangle_fan *)object;
 
@@ -216,7 +217,7 @@ dmnsn_smooth_triangle_fan_intersection_fn(const dmnsn_object *object, dmnsn_line
 
   for (size_t i = 0; i < fan->ncoeffs; ++i) {
     const double *coeffs = fan->coeffs[i];
-    l = dmnsn_change_line_basis(coeffs, l);
+    l = dmnsn_change_ray_basis(coeffs, l);
     nab = nac;
     nac = dmnsn_new_vector(coeffs[6], coeffs[7], coeffs[8]);
 
@@ -236,18 +237,18 @@ dmnsn_smooth_triangle_fan_intersection_fn(const dmnsn_object *object, dmnsn_line
 }
 
 /// Smooth triangle fan bounding callback.
-static dmnsn_bounding_box
+static dmnsn_aabb
 dmnsn_smooth_triangle_fan_bounding_fn(const dmnsn_object *object, dmnsn_matrix trans)
 {
   const dmnsn_smooth_triangle_fan *fan = (const dmnsn_smooth_triangle_fan *)object;
 
-  dmnsn_bounding_box box = dmnsn_bound_first_triangle(trans);
+  dmnsn_aabb box = dmnsn_bound_first_triangle(trans);
 
   for (size_t i = 0; i < fan->ncoeffs; ++i) {
     dmnsn_matrix incremental = dmnsn_decompress_coeffs(fan->coeffs[i]);
     trans = dmnsn_matrix_mul(trans, dmnsn_matrix_inverse(incremental));
     dmnsn_vector vertex = dmnsn_transform_point(trans, dmnsn_y);
-    box = dmnsn_bounding_box_swallow(box, vertex);
+    box = dmnsn_aabb_swallow(box, vertex);
   }
 
   return box;
diff --git a/libdimension/pigment.c b/libdimension/model/pigment.c
index 5fae4f6..a8c7c0b 100644
--- a/libdimension/pigment.c
+++ b/libdimension/model/pigment.c
@@ -23,7 +23,7 @@
  * Pigments.
  */
 
-#include "dimension-internal.h"
+#include "dimension/model.h"
 
 // Allocate a dummy pigment
 dmnsn_pigment *
diff --git a/libdimension/canvas_pigment.c b/libdimension/model/pigments/canvas_pigment.c
index 2b1e2b6..bb83b0a 100644
--- a/libdimension/canvas_pigment.c
+++ b/libdimension/model/pigments/canvas_pigment.c
@@ -23,7 +23,7 @@
  * Image maps.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 
 /// Canvas pigment type.
 typedef struct dmnsn_canvas_pigment {
diff --git a/libdimension/pigment_map.c b/libdimension/model/pigments/pigment_map.c
index 71de4e0..a74a66a 100644
--- a/libdimension/pigment_map.c
+++ b/libdimension/model/pigments/pigment_map.c
@@ -23,7 +23,7 @@
  * Pigment-mapped pigment patterns.
  */
 
-#include "dimension.h"
+#include "dimension/model.h"
 
 /// Initialize a pigment in a pigment map.
 static void
diff --git a/libdimension/solid_pigment.c b/libdimension/model/pigments/solid_pigment.c
index 8c94e1c..bdbd724 100644
--- a/libdimension/solid_pigment.c
+++ b/libdimension/model/pigments/solid_pigment.c
@@ -23,8 +23,7 @@
  * Solid color pigments.
  */
 
-#include "dimension.h"
-#include <stdlib.h>
+#include "dimension/model.h"
 
 // Create a solid color
 dmnsn_pigment *
diff --git a/libdimension/scene.c b/libdimension/model/scene.c
index 6c9c495..875c376 100644
--- a/libdimension/scene.c
+++ b/libdimension/model/scene.c
@@ -23,7 +23,9 @@
  * Scenes.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/platform.h"
+#include "dimension/model.h"
 #include <stdlib.h>
 
 // Allocate an empty scene
diff --git a/libdimension/texture.c b/libdimension/model/texture.c
index 515e260..b7eb7ef 100644
--- a/libdimension/texture.c
+++ b/libdimension/model/texture.c
@@ -23,7 +23,7 @@
  * Textures.
  */
 
-#include "dimension-internal.h"
+#include "dimension/model.h"
 
 dmnsn_texture *
 dmnsn_new_texture(dmnsn_pool *pool)
diff --git a/libdimension/checker.c b/libdimension/pattern/checker.c
index 740abe8..cce9623 100644
--- a/libdimension/checker.c
+++ b/libdimension/pattern/checker.c
@@ -23,7 +23,7 @@
  * Checker pattern.
  */
 
-#include "dimension-internal.h"
+#include "dimension/pattern.h"
 
 /// Checker pattern callback.
 static double
diff --git a/libdimension/gradient.c b/libdimension/pattern/gradient.c
index 74e8045..f7a2b97 100644
--- a/libdimension/gradient.c
+++ b/libdimension/pattern/gradient.c
@@ -23,7 +23,7 @@
  * Gradient pattern.
  */
 
-#include "dimension.h"
+#include "dimension/pattern.h"
 
 /// Gradient pattern type.
 typedef struct dmnns_gradient {
diff --git a/libdimension/leopard.c b/libdimension/pattern/leopard.c
index a4000d2..1a7bce0 100644
--- a/libdimension/leopard.c
+++ b/libdimension/pattern/leopard.c
@@ -23,7 +23,7 @@
  * Leopard pattern.
  */
 
-#include "dimension-internal.h"
+#include "dimension/pattern.h"
 #include <math.h>
 
 /// Leopard pattern callback.
diff --git a/libdimension/map.c b/libdimension/pattern/map.c
index f5454b1..ac07960 100644
--- a/libdimension/map.c
+++ b/libdimension/pattern/map.c
@@ -23,7 +23,8 @@
  * Generic maps.
  */
 
-#include "dimension-internal.h"
+#include "internal.h"
+#include "dimension/pattern.h"
 
 /// dmnsn_map definition.
 struct dmnsn_map {
diff --git a/libdimension/pattern.c b/libdimension/pattern/pattern.c
index 22f5c13..60e1ae7 100644
--- a/libdimension/pattern.c
+++ b/libdimension/pattern/pattern.c
@@ -23,7 +23,7 @@
  * Patterns.
  */
 
-#include "dimension-internal.h"
+#include "dimension/pattern.h"
 
 dmnsn_pattern *
 dmnsn_new_pattern(dmnsn_pool *pool)
diff --git a/libdimension/platform.c b/libdimension/platform/platform.c
index b5d23cc..84d9bf0 100644
--- a/libdimension/platform.c
+++ b/libdimension/platform/platform.c
@@ -24,7 +24,9 @@
  */
 
 #define _GNU_SOURCE
-#include "dimension-internal.h"
+#include "internal.h"
+#include "internal/platform.h"
+#include "dimension/platform.h"
 
 #if HAVE_UNISTD_H
   #include <unistd.h>
diff --git a/libdimension/timer.c b/libdimension/platform/timer.c
index 35bf05d..2d91469 100644
--- a/libdimension/timer.c
+++ b/libdimension/platform/timer.c
@@ -23,7 +23,8 @@
  * Performance counter.
  */
 
-#include "dimension-internal.h"
+#include "internal/platform.h"
+#include "dimension/platform.h"
 
 void
 dmnsn_timer_start(dmnsn_timer *timer)
diff --git a/libdimension/ray_trace.c b/libdimension/render/render.c
index 6f4ce33..842b41e 100644
--- a/libdimension/ray_trace.c
+++ b/libdimension/render/render.c
@@ -23,7 +23,9 @@
  * The ray-tracing algorithm.
  */
 
-#include "dimension-internal.h"
+#include "internal/bvh.h"
+#include "internal/concurrency.h"
+#include "dimension/render.h"
 #include <stdlib.h>
 
 ////////////////////////////////////
@@ -35,45 +37,45 @@ typedef struct {
   dmnsn_future *future;
   dmnsn_scene *scene;
   dmnsn_bvh *bvh;
-} dmnsn_ray_trace_payload;
+} dmnsn_render_payload;
 
 // Ray-trace a scene
 void
-dmnsn_ray_trace(dmnsn_scene *scene)
+dmnsn_render(dmnsn_scene *scene)
 {
-  dmnsn_future *future = dmnsn_ray_trace_async(scene);
+  dmnsn_future *future = dmnsn_render_async(scene);
   if (dmnsn_future_join(future) != 0) {
     dmnsn_error("Error occured while ray-tracing.");
   }
 }
 
 /// Background thread callback.
-static int dmnsn_ray_trace_scene_thread(void *ptr);
+static int dmnsn_render_scene_thread(void *ptr);
 
 // Ray-trace a scene in the background
 dmnsn_future *
-dmnsn_ray_trace_async(dmnsn_scene *scene)
+dmnsn_render_async(dmnsn_scene *scene)
 {
   dmnsn_future *future = dmnsn_new_future();
 
-  dmnsn_ray_trace_payload *payload = DMNSN_MALLOC(dmnsn_ray_trace_payload);
+  dmnsn_render_payload *payload = DMNSN_MALLOC(dmnsn_render_payload);
   payload->future = future;
   payload->scene  = scene;
 
-  dmnsn_new_thread(future, dmnsn_ray_trace_scene_thread, payload);
+  dmnsn_new_thread(future, dmnsn_render_scene_thread, payload);
 
   return future;
 }
 
 /// Worker thread callback.
-static int dmnsn_ray_trace_scene_concurrent(void *ptr, unsigned int thread,
+static int dmnsn_render_scene_concurrent(void *ptr, unsigned int thread,
                                             unsigned int nthreads);
 
 // Thread callback -- set up the multithreaded engine
 static int
-dmnsn_ray_trace_scene_thread(void *ptr)
+dmnsn_render_scene_thread(void *ptr)
 {
-  dmnsn_ray_trace_payload *payload = ptr;
+  dmnsn_render_payload *payload = ptr;
 
   // Pre-calculate bounding box transformations, etc.
   dmnsn_scene_initialize(payload->scene);
@@ -89,7 +91,7 @@ dmnsn_ray_trace_scene_thread(void *ptr)
   // Time the render itself
   dmnsn_timer_start(&payload->scene->render_timer);
     int ret = dmnsn_execute_concurrently(payload->future,
-                                         dmnsn_ray_trace_scene_concurrent,
+                                         dmnsn_render_scene_concurrent,
                                          payload, payload->scene->nthreads);
   dmnsn_timer_stop(&payload->scene->render_timer);
 
@@ -135,15 +137,14 @@ typedef struct dmnsn_rtstate {
 static inline void
 dmnsn_rtstate_initialize(dmnsn_rtstate *state,
                          const dmnsn_intersection *intersection);
-/// Main helper for dmnsn_ray_trace_scene_concurrent - shoot a ray.
-static dmnsn_tcolor dmnsn_ray_shoot(dmnsn_rtstate *state, dmnsn_line ray);
+/// Main helper for dmnsn_render_scene_concurrent - shoot a ray.
+static dmnsn_tcolor dmnsn_ray_shoot(dmnsn_rtstate *state, dmnsn_ray ray);
 
 // Actually ray-trace a scene
 static int
-dmnsn_ray_trace_scene_concurrent(void *ptr, unsigned int thread,
-                                 unsigned int nthreads)
+dmnsn_render_scene_concurrent(void *ptr, unsigned int thread, unsigned int nthreads)
 {
-  const dmnsn_ray_trace_payload *payload = ptr;
+  const dmnsn_render_payload *payload = ptr;
   dmnsn_future *future = payload->future;
   dmnsn_scene *scene = payload->scene;
   dmnsn_bvh *bvh = payload->bvh;
@@ -158,7 +159,7 @@ dmnsn_ray_trace_scene_concurrent(void *ptr, unsigned int thread,
   for (size_t y = thread; y < scene->canvas->height; y += nthreads) {
     for (size_t x = 0; x < scene->canvas->width; ++x) {
       // Get the ray corresponding to the (x,y)'th pixel
-      dmnsn_line ray = dmnsn_camera_ray(
+      dmnsn_ray ray = dmnsn_camera_ray(
         scene->camera,
         ((double)(x + scene->region_x))/(scene->outer_width - 1),
         ((double)(y + scene->region_y))/(scene->outer_height - 1)
@@ -187,7 +188,7 @@ dmnsn_rtstate_initialize(dmnsn_rtstate *state,
   state->texture      = intersection->object->texture;
   state->interior     = intersection->object->interior;
 
-  state->r = dmnsn_line_point(intersection->ray, intersection->t);
+  state->r = dmnsn_ray_point(intersection->ray, intersection->t);
   state->pigment_r = dmnsn_transform_point(
     intersection->object->pigment_trans,
     state->r
@@ -206,7 +207,7 @@ dmnsn_rtstate_initialize(dmnsn_rtstate *state,
 }
 
 /// Calculate the background color.
-static void dmnsn_trace_background(dmnsn_rtstate *state, dmnsn_line ray);
+static void dmnsn_trace_background(dmnsn_rtstate *state, dmnsn_ray ray);
 /// Calculate the base pigment at the intersection.
 static void dmnsn_trace_pigment(dmnsn_rtstate *state);
 /// Handle light, shadow, and shading.
@@ -218,7 +219,7 @@ static void dmnsn_trace_transparency(dmnsn_rtstate *state);
 
 // Shoot a ray, and calculate the color
 static dmnsn_tcolor
-dmnsn_ray_shoot(dmnsn_rtstate *state, dmnsn_line ray)
+dmnsn_ray_shoot(dmnsn_rtstate *state, dmnsn_ray ray)
 {
   if (state->reclevel == 0
       || dmnsn_color_intensity(state->adc_value) < state->scene->adc_bailout)
@@ -254,7 +255,7 @@ dmnsn_ray_shoot(dmnsn_rtstate *state, dmnsn_line ray)
 }
 
 static void
-dmnsn_trace_background(dmnsn_rtstate *state, dmnsn_line ray)
+dmnsn_trace_background(dmnsn_rtstate *state, dmnsn_ray ray)
 {
   dmnsn_pigment *background = state->scene->background;
   if (state->scene->quality & DMNSN_RENDER_PIGMENT) {
@@ -340,12 +341,12 @@ dmnsn_evaluate_diffuse(const dmnsn_rtstate *state)
 static bool
 dmnsn_trace_light_ray(dmnsn_rtstate *state, const dmnsn_light *light)
 {
-  dmnsn_line shadow_ray = dmnsn_new_line(
+  dmnsn_ray shadow_ray = dmnsn_new_ray(
     state->r,
     light->direction_fn(light, state->r)
   );
   // Add epsilon to avoid hitting ourselves with the shadow ray
-  shadow_ray = dmnsn_line_add_epsilon(shadow_ray);
+  shadow_ray = dmnsn_ray_add_epsilon(shadow_ray);
 
   // Check if we're casting a shadow on ourself
   if ((dmnsn_vector_dot(shadow_ray.n, state->intersection->normal)
@@ -455,8 +456,8 @@ dmnsn_trace_reflection(dmnsn_rtstate *state)
 {
   const dmnsn_reflection *reflection = state->texture->finish.reflection;
   if (reflection) {
-    dmnsn_line refl_ray = dmnsn_new_line(state->r, state->reflected);
-    refl_ray = dmnsn_line_add_epsilon(refl_ray);
+    dmnsn_ray refl_ray = dmnsn_new_ray(state->r, state->reflected);
+    refl_ray = dmnsn_ray_add_epsilon(refl_ray);
 
     dmnsn_rtstate recursive_state = *state;
 
@@ -481,8 +482,8 @@ dmnsn_trace_transparency(dmnsn_rtstate *state)
   if (state->pigment.T >= dmnsn_epsilon) {
     const dmnsn_interior *interior = state->interior;
 
-    dmnsn_line trans_ray = dmnsn_new_line(state->r, state->intersection->ray.n);
-    trans_ray = dmnsn_line_add_epsilon(trans_ray);
+    dmnsn_ray trans_ray = dmnsn_new_ray(state->r, state->intersection->ray.n);
+    trans_ray = dmnsn_ray_add_epsilon(trans_ray);
 
     dmnsn_vector r = dmnsn_vector_normalized(trans_ray.n);
     dmnsn_vector n = state->intersection->normal;
diff --git a/libdimension/tests/future.c b/libdimension/tests/future.c
index ee2511a..a223a44 100644
--- a/libdimension/tests/future.c
+++ b/libdimension/tests/future.c
@@ -22,9 +22,9 @@
  * Tests for dmnsn_future.
  */
 
-#include "../platform.c"
-#include "../future.c"
-#include "../threads.c"
+#include "../platform/platform.c"
+#include "../concurrency/future.c"
+#include "../concurrency/threads.c"
 #include "tests.h"
 #include <stdatomic.h>
 
diff --git a/libdimension/tests/polynomial.c b/libdimension/tests/polynomial.c
index 492c22c..a983bd4 100644
--- a/libdimension/tests/polynomial.c
+++ b/libdimension/tests/polynomial.c
@@ -22,8 +22,8 @@
  * Basic tests of the polynomial root-finder.
  */
 
+#include "../math/polynomial.c"
 #include "tests.h"
-#include "../polynomial.c"
 #include <stdarg.h>
 
 #define DMNSN_CLOSE_ENOUGH 1.0e-6
diff --git a/libdimension/tests/pool.c b/libdimension/tests/pool.c
index 1722d66..6090c36 100644
--- a/libdimension/tests/pool.c
+++ b/libdimension/tests/pool.c
@@ -22,10 +22,10 @@
  * Tests for memory pools.
  */
 
-#include "../dimension-internal.h"
+#include "../concurrency/future.c"
+#include "../concurrency/threads.c"
+#include "../internal.h"
 #include "tests.h"
-#include "../future.c"
-#include "../threads.c"
 
 static dmnsn_pool *pool;
 
diff --git a/libdimension/tests/prtree.c b/libdimension/tests/prtree.c
index acdd426..8824e93 100644
--- a/libdimension/tests/prtree.c
+++ b/libdimension/tests/prtree.c
@@ -22,28 +22,28 @@
  * Basic tests of PR-trees.
  */
 
-#include "../platform.c"
-#include "../threads.c"
-#include "../future.c"
-#include "../bvh.c"
-#include "../prtree.c"
+#include "../platform/platform.c"
+#include "../concurrency/threads.c"
+#include "../concurrency/future.c"
+#include "../bvh/bvh.c"
+#include "../bvh/prtree.c"
 #include <stdio.h>
 #include <stdlib.h>
 
 unsigned int calls = 0;
 
 static bool
-dmnsn_fake_intersection_fn(const dmnsn_object *object, dmnsn_line line,
+dmnsn_fake_intersection_fn(const dmnsn_object *object, dmnsn_ray ray,
                            dmnsn_intersection *intersection)
 {
-  intersection->t = (object->bounding_box.min.z - line.x0.z)/line.n.z;
+  intersection->t = (object->aabb.min.z - ray.x0.z)/ray.n.z;
   intersection->normal = dmnsn_x;
   ++calls;
   return true;
 }
 
 static void
-dmnsn_randomize_bounding_box(dmnsn_object *object)
+dmnsn_randomize_aabb(dmnsn_object *object)
 {
   dmnsn_vector a, b;
 
@@ -55,8 +55,8 @@ dmnsn_randomize_bounding_box(dmnsn_object *object)
   b.y = 2.0*((double)rand())/RAND_MAX - 1.0;
   b.z = 2.0*((double)rand())/RAND_MAX - 1.0;
 
-  object->bounding_box.min = dmnsn_vector_min(a, b);
-  object->bounding_box.max = dmnsn_vector_max(a, b);
+  object->aabb.min = dmnsn_vector_min(a, b);
+  object->aabb.max = dmnsn_vector_max(a, b);
 }
 
 static const dmnsn_object_vtable dmnsn_fake_vtable = {
@@ -67,7 +67,7 @@ static dmnsn_object *
 dmnsn_new_fake_object(dmnsn_pool *pool)
 {
   dmnsn_object *object = dmnsn_new_object(pool);
-  dmnsn_randomize_bounding_box(object);
+  dmnsn_randomize_aabb(object);
   object->vtable = &dmnsn_fake_vtable;
   object->trans_inv = dmnsn_identity_matrix();
   return object;
@@ -92,7 +92,7 @@ main(void)
   dmnsn_bvh *bvh = dmnsn_new_bvh(objects, DMNSN_BVH_PRTREE);
 
   dmnsn_intersection intersection;
-  dmnsn_line ray = dmnsn_new_line(
+  dmnsn_ray ray = dmnsn_new_ray(
     dmnsn_new_vector(0.0, 0.0, -2.0),
     dmnsn_new_vector(0.0, 0.0, 1.0)
   );
diff --git a/libdimension/tests/render.c b/libdimension/tests/render.c
index 56a80d3..fd1e96d 100644
--- a/libdimension/tests/render.c
+++ b/libdimension/tests/render.c
@@ -350,7 +350,7 @@ main(void)
   // Render the scene
 
   printf("Rendering scene\n");
-  dmnsn_future *future = dmnsn_ray_trace_async(scene);
+  dmnsn_future *future = dmnsn_render_async(scene);
 
   // Display the scene as it's rendered
   if (display) {