From 6932058f59d8be7efb5e6608ce667416d1c0a5f6 Mon Sep 17 00:00:00 2001
From: Christophe Riccio <christophe.riccio@unity3d.com>
Date: Tue, 22 Aug 2017 22:47:36 +0200
Subject: [PATCH] Fixed quaternion constructor from two vectors special cases
 #469

---
 glm/gtc/quaternion.inl      | 22 ++++++++++++++++++----
 readme.md                   |  1 +
 test/gtc/gtc_quaternion.cpp | 22 +++++++++++++++++-----
 3 files changed, 36 insertions(+), 9 deletions(-)

diff --git a/glm/gtc/quaternion.inl b/glm/gtc/quaternion.inl
index eed933fe..edfb6e9a 100644
--- a/glm/gtc/quaternion.inl
+++ b/glm/gtc/quaternion.inl
@@ -148,11 +148,25 @@ namespace detail
 	template<typename T, qualifier Q>
 	GLM_FUNC_QUALIFIER tquat<T, Q>::tquat(vec<3, T, Q> const& u, vec<3, T, Q> const& v)
 	{
-		vec<3, T, Q> const LocalW(cross(u, v));
-		T Dot = detail::compute_dot<vec<3, T, Q>, T, detail::is_aligned<Q>::value>::call(u, v);
-		tquat<T, Q> q(T(1) + Dot, LocalW.x, LocalW.y, LocalW.z);
+		T norm_u_norm_v = sqrt(dot(u, u) * dot(v, v));
+		T real_part = norm_u_norm_v + dot(u, v);
+		vec<3, T, Q> w;
 
-		*this = normalize(q);
+		if(real_part < static_cast<T>(1.e-6f) * norm_u_norm_v)
+		{
+			// If u and v are exactly opposite, rotate 180 degrees
+			// around an arbitrary orthogonal axis. Axis normalisation
+			// can happen later, when we normalise the quaternion.
+			real_part = static_cast<T>(0);
+			w = abs(u.x) > abs(u.z) ? vec<3, T, Q>(-u.y, u.x, static_cast<T>(0)) : vec<3, T, Q>(static_cast<T>(0), -u.z, u.y);
+		}
+		else
+		{
+			// Otherwise, build quaternion the standard way.
+			w = cross(u, v);
+		}
+
+	    *this = normalize(tquat<T, Q>(real_part, w.x, w.y, w.z));
 	}
 
 	template<typename T, qualifier Q>
diff --git a/readme.md b/readme.md
index 76e0ea5b..74043af7 100644
--- a/readme.md
+++ b/readme.md
@@ -101,6 +101,7 @@ glm::mat4 camera(float Translate, glm::vec2 const& Rotate)
 - Fixed documentation warnings
 - Fixed GLM_HAS_OPENMP when OpenMP is not enabled
 - Fixed Better follow GLSL min and max specification #372
+- Fixed quaternion constructor from two vectors special cases #469 
 
 #### Deprecation:
 - Requires Visual Studio 2013, GCC 4.7, Clang 3.4, Cuda 7, ICC 2013 or a C++11 compiler
diff --git a/test/gtc/gtc_quaternion.cpp b/test/gtc/gtc_quaternion.cpp
index f94adf5d..b006a12f 100644
--- a/test/gtc/gtc_quaternion.cpp
+++ b/test/gtc/gtc_quaternion.cpp
@@ -225,16 +225,28 @@ int test_quat_mul()
 
 int test_quat_two_axis_ctr()
 {
-	int Error(0);
+	int Error = 0;
 
-	glm::quat q1(glm::vec3(1, 0, 0), glm::vec3(0, 1, 0));
-	glm::vec3 v1 = q1 * glm::vec3(1, 0, 0);
+	glm::quat const q1(glm::vec3(1, 0, 0), glm::vec3(0, 1, 0));
+	glm::vec3 const v1 = q1 * glm::vec3(1, 0, 0);
 	Error += glm::all(glm::epsilonEqual(v1, glm::vec3(0, 1, 0), 0.0001f)) ? 0 : 1;
 
-	glm::quat q2 = q1 * q1;
-	glm::vec3 v2 = q2 * glm::vec3(1, 0, 0);
+	glm::quat const q2 = q1 * q1;
+	glm::vec3 const v2 = q2 * glm::vec3(1, 0, 0);
 	Error += glm::all(glm::epsilonEqual(v2, glm::vec3(-1, 0, 0), 0.0001f)) ? 0 : 1;
 
+	glm::quat const q3(glm::vec3(1, 0, 0), glm::vec3(-1, 0, 0));
+	glm::vec3 const v3 = q3 * glm::vec3(1, 0, 0);
+	Error += glm::all(glm::epsilonEqual(v3, glm::vec3(-1, 0, 0), 0.0001f)) ? 0 : 1;
+
+	glm::quat const q4(glm::vec3(0, 1, 0), glm::vec3(0, -1, 0));
+	glm::vec3 const v4 = q4 * glm::vec3(0, 1, 0);
+	Error += glm::all(glm::epsilonEqual(v4, glm::vec3(0, -1, 0), 0.0001f)) ? 0 : 1;
+
+	glm::quat const q5(glm::vec3(0, 0, 1), glm::vec3(0, 0, -1));
+	glm::vec3 const v5 = q5 * glm::vec3(0, 0, 1);
+	Error += glm::all(glm::epsilonEqual(v5, glm::vec3(0, 0, -1), 0.0001f)) ? 0 : 1;
+
 	return Error;
 }