new turn style

[blank.git] / src / model / geometry.cpp

diff --git a/src/model/geometry.cpp b/src/model/geometry.cpp
index c2816b44a0f585f67781ef7e6bbe4e811dcd4fe9..a4ac5006d55a019544f122a035af5d98d21cf308 100644 (file)
--- a/src/model/geometry.cpp
+++ b/src/model/geometry.cpp
@@ -1,10 +1,41 @@
 #include "geometry.hpp"
 
 #include <limits>
+#include <glm/gtx/matrix_cross_product.hpp>
+#include <glm/gtx/optimum_pow.hpp>
+#include <glm/gtx/transform.hpp>
 
 
 namespace blank {
 
+glm::mat3 find_rotation(const glm::vec3 &a, const glm::vec3 &b) noexcept {
+       glm::vec3 v(cross(a, b));
+       if (iszero(v)) {
+               // a and b are parallel
+               if (iszero(a - b)) {
+                       // a and b are identical
+                       return glm::mat3(1.0f);
+               } else {
+                       // a and b are opposite
+                       // create arbitrary unit vector perpendicular to a and
+                       // rotate 180° around it
+                       glm::vec3 arb(a);
+                       if (std::abs(a.x - 1.0f) > std::numeric_limits<float>::epsilon()) {
+                               arb.x += 1.0f;
+                       } else {
+                               arb.y += 1.0f;
+                       }
+                       glm::vec3 axis(normalize(cross(a, arb)));
+                       return glm::mat3(glm::rotate(PI, axis));
+               }
+       }
+       float mv = length_squared(v);
+       float c = dot(a, b);
+       float f = (1 - c) / mv;
+       glm::mat3 vx(matrixCross3(v));
+       return glm::mat3(1.0f) + vx + (pow2(vx) * f);
+}
+
 bool Intersection(
        const Ray &ray,
        const AABB &aabb,
@@ -47,20 +78,17 @@ bool Intersection(
                *dist = t_min;
        }
        if (normal) {
-               glm::vec4 norm(0.0f);
                if (min_all.x > min_all.y) {
                        if (min_all.x > min_all.z) {
-                               norm.x = t2.x < t1.x ? 1 : -1;
+                               normal->x = t2.x < t1.x ? 1 : -1;
                        } else {
-                               norm.z = t2.z < t1.z ? 1 : -1;
+                               normal->z = t2.z < t1.z ? 1 : -1;
                        }
                } else if (min_all.y > min_all.z) {
-                       norm.y = t2.y < t1.y ? 1 : -1;
+                       normal->y = t2.y < t1.y ? 1 : -1;
                } else {
-                       norm.z = t2.z < t1.z ? 1 : -1;
+                       normal->z = t2.z < t1.z ? 1 : -1;
                }
-               norm = M * norm;
-               *normal = glm::vec3(norm);
        }
        return true;
 }