construct quaternion from angle/axis

include/cif++/Point.hpp

@@ -378,7 +378,9 @@ Point Nudge(Point p, float offset);
 
 Quaternion Normalize(Quaternion q);
 
+Quaternion ConstructFromAngleAxis(float angle, Point axis);
 std::tuple<double, Point> QuaternionToAngleAxis(Quaternion q);
+
 Point Centroid(const std::vector<Point> &Points);
 Point CenterPoints(std::vector<Point> &Points);
 

src/Point.cpp

@@ -295,6 +295,19 @@ Quaternion Normalize(Quaternion q)
 
 // --------------------------------------------------------------------
 
+Quaternion ConstructFromAngleAxis(float angle, Point axis)
+{
+	auto q = std::cos((angle * mmcif::kPI / 180) / 2);
+	auto s = std::sqrt(1 - q * q);
+
+	axis.normalize();
+
+	return Normalize(Quaternion{q,
+		static_cast<float>(s * axis.mX),
+		static_cast<float>(s * axis.mY),
+		static_cast<float>(s * axis.mZ)});
+}
+
 std::tuple<double,Point> QuaternionToAngleAxis(Quaternion q)
 {
 	if (q.R_component_1() > 1)

test/unit-test.cpp

@@ -1745,6 +1745,47 @@ BOOST_AUTO_TEST_CASE(t1)
 	// std::cout << "rmsd: " << RMSd(p1, p2) << std::endl;
 }
 
+BOOST_AUTO_TEST_CASE(t2)
+{
+	Point p[] = {
+		{ 1, 1, 0 },
+		{ 2, 1, 0 },
+		{ 1, 2, 0 }
+	};
+
+	Point xp = mmcif::CrossProduct(p[1] - p[0], p[2] - p[0]);
+
+	Quaternion q = ConstructFromAngleAxis(45, xp); //mmcif::Normalize(Quaternion{45 * mmcif::kPI / 180, xp.mX, xp.mY, xp.mZ});
+
+	auto &&[angle, axis] = mmcif::QuaternionToAngleAxis(q);
+
+	BOOST_TEST(angle == 45, tt::tolerance(0.01));
+}
+
+BOOST_AUTO_TEST_CASE(t3)
+{
+	Point p[] = {
+		{ 1, 1, 0 },
+		{ 2, 1, 0 },
+		{ 1, 2, 0 }
+	};
+
+	Point xp = mmcif::CrossProduct(p[1] - p[0], p[2] - p[0]);
+
+	Quaternion q = ConstructFromAngleAxis(45, xp); //mmcif::Normalize(Quaternion{45 * mmcif::kPI / 180, xp.mX, xp.mY, xp.mZ});
+
+	Point v = p[1];
+	v -= p[0];
+	v.rotate(q);
+	v += p[0];
+
+	std::cout << v << std::endl;
+
+	double a = mmcif::Angle(v, p[0], p[1]);
+
+	BOOST_TEST(a == 45, tt::tolerance(0.01));
+}
+
 BOOST_AUTO_TEST_CASE(parser_test_1)
 {
 	auto data1 = R"(