geometry: Faster AABB transformations.

3 files changed, 55 insertions, 42 deletions

diff --git a/libdimension/bench/geometry.c b/libdimension/bench/geometry.c
index 56c8b58..68c9885 100644
--- a/libdimension/bench/geometry.c
+++ b/libdimension/bench/geometry.c
@@ -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 Benchmark Suite.                   *
  *                                                                       *
@@ -27,6 +27,7 @@ main(void)
   dmnsn_vector vector, vector2;
   dmnsn_matrix matrix, matrix2;
   dmnsn_line line;
+  dmnsn_bounding_box box;
   double result;
 
   sandglass_t sandglass;
@@ -41,7 +42,7 @@ main(void)
   });
   printf("dmnsn_new_vector(): %ld\n", sandglass.grains);
 
-  // dmnsn_matrix_construct()
+  // dmnsn_new_matrix()
   sandglass_bench_fine(&sandglass, {
     matrix = dmnsn_new_matrix(1.0, 1.0, 0.0, 0.0,
                               1.0, 1.0, 1.0, 0.0,
@@ -80,6 +81,13 @@ main(void)
   });
   printf("dmnsn_new_line(): %ld\n", sandglass.grains);
 
+  // dmnsn_new_bounding_box()
+  vector2 = dmnsn_new_vector(3.0, 4.0, 5.0);
+  sandglass_bench_fine(&sandglass, {
+    box = dmnsn_new_bounding_box(vector, vector2);
+  });
+  printf("dmnsn_new_bounding_box(): %ld\n", sandglass.grains);
+
   // dmnsn_vector_add()
   sandglass_bench_fine(&sandglass, {
     vector = dmnsn_vector_add(vector, vector2);
@@ -173,6 +181,12 @@ main(void)
   });
   printf("dmnsn_transform_line(): %ld\n", sandglass.grains);
 
+  // dmnsn_transform_bounding_box()
+  sandglass_bench_fine(&sandglass, {
+    box = dmnsn_transform_bounding_box(matrix, box);
+  });
+  printf("dmnsn_transform_bounding_box(): %ld\n", sandglass.grains);
+
   // dmnsn_line_point()
   sandglass_bench_fine(&sandglass, {
     vector = dmnsn_line_point(line, result);
diff --git a/libdimension/dimension/geometry.h b/libdimension/dimension/geometry.h
index 26dd608..b2ba551 100644
--- a/libdimension/dimension/geometry.h
+++ b/libdimension/dimension/geometry.h
@@ -170,6 +170,13 @@ dmnsn_new_matrix4(dmnsn_vector a, dmnsn_vector b, dmnsn_vector c,
   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);
 
diff --git a/libdimension/geometry.c b/libdimension/geometry.c
index 1ec2801..c93dd11 100644
--- a/libdimension/geometry.c
+++ b/libdimension/geometry.c
@@ -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.                           *
  *                                                                       *
@@ -358,48 +358,40 @@ dmnsn_bounding_box
 dmnsn_transform_bounding_box(dmnsn_matrix trans, dmnsn_bounding_box box)
 {
   // Infinite/zero bounding box support
-  if (isinf(box.min.x))
+  if (isinf(box.min.x)) {
     return box;
+  }
+
+  // Taking the "absolute value" of the matrix saves some min/max calculations
+  for (int i = 0; i < 3; ++i) {
+    for (int j = 0; j < 3; ++j) {
+      trans.n[i][j] = fabs(trans.n[i][j]);
+    }
+  }
+
+  dmnsn_vector Mx = dmnsn_matrix_column(trans, 0);
+  dmnsn_vector xmin = dmnsn_vector_mul(box.min.x, Mx);
+  dmnsn_vector xmax = dmnsn_vector_mul(box.max.x, Mx);
+
+  dmnsn_vector My = dmnsn_matrix_column(trans, 1);
+  dmnsn_vector ymin = dmnsn_vector_mul(box.min.y, My);
+  dmnsn_vector ymax = dmnsn_vector_mul(box.max.y, My);
+
+  dmnsn_vector Mz = dmnsn_matrix_column(trans, 2);
+  dmnsn_vector zmin = dmnsn_vector_mul(box.min.z, Mz);
+  dmnsn_vector zmax = dmnsn_vector_mul(box.max.z, Mz);
+
+  dmnsn_vector Mt = dmnsn_matrix_column(trans, 3);
 
-  dmnsn_vector corner;
   dmnsn_bounding_box ret;
-  ret.min = dmnsn_transform_point(trans, box.min);
-  ret.max = ret.min;
-
-  corner  = dmnsn_new_vector(box.min.x, box.min.y, box.max.z);
-  corner  = dmnsn_transform_point(trans, corner);
-  ret.min = dmnsn_vector_min(ret.min, corner);
-  ret.max = dmnsn_vector_max(ret.max, corner);
-
-  corner  = dmnsn_new_vector(box.min.x, box.max.y, box.min.z);
-  corner  = dmnsn_transform_point(trans, corner);
-  ret.min = dmnsn_vector_min(ret.min, corner);
-  ret.max = dmnsn_vector_max(ret.max, corner);
-
-  corner  = dmnsn_new_vector(box.min.x, box.max.y, box.max.z);
-  corner  = dmnsn_transform_point(trans, corner);
-  ret.min = dmnsn_vector_min(ret.min, corner);
-  ret.max = dmnsn_vector_max(ret.max, corner);
-
-  corner  = dmnsn_new_vector(box.max.x, box.min.y, box.min.z);
-  corner  = dmnsn_transform_point(trans, corner);
-  ret.min = dmnsn_vector_min(ret.min, corner);
-  ret.max = dmnsn_vector_max(ret.max, corner);
-
-  corner  = dmnsn_new_vector(box.max.x, box.min.y, box.max.z);
-  corner  = dmnsn_transform_point(trans, corner);
-  ret.min = dmnsn_vector_min(ret.min, corner);
-  ret.max = dmnsn_vector_max(ret.max, corner);
-
-  corner  = dmnsn_new_vector(box.max.x, box.max.y, box.min.z);
-  corner  = dmnsn_transform_point(trans, corner);
-  ret.min = dmnsn_vector_min(ret.min, corner);
-  ret.max = dmnsn_vector_max(ret.max, corner);
-
-  corner  = dmnsn_new_vector(box.max.x, box.max.y, box.max.z);
-  corner  = dmnsn_transform_point(trans, corner);
-  ret.min = dmnsn_vector_min(ret.min, corner);
-  ret.max = dmnsn_vector_max(ret.max, corner);
+
+  ret.min = dmnsn_vector_add(xmin, ymin);
+  ret.min = dmnsn_vector_add(ret.min, zmin);
+  ret.min = dmnsn_vector_add(ret.min, Mt);
+
+  ret.max = dmnsn_vector_add(xmax, ymax);
+  ret.max = dmnsn_vector_add(ret.max, zmax);
+  ret.max = dmnsn_vector_add(ret.max, Mt);
 
   return ret;
 }