Computation of bounding sphere of triangle.

src/typed-geometry/detail/functions.hh

@@ -3,6 +3,7 @@
 #include <typed-geometry/functions/aabb.hh>
 #include <typed-geometry/functions/angle.hh>
 #include <typed-geometry/functions/area.hh>
+#include <typed-geometry/functions/bounds.hh>
 #include <typed-geometry/functions/centroid.hh>
 #include <typed-geometry/functions/closest_points.hh>
 #include <typed-geometry/functions/contains.hh>

src/typed-geometry/functions/bounds.hh

+#pragma once
+
+#include <typed-geometry/functions/interpolate.hh>
+#include <typed-geometry/tests/triangle_tests.hh>
+#include <typed-geometry/types/objects/sphere.hh>
+#include <typed-geometry/types/objects/triangle.hh>
+
+namespace tg
+{
+template <int D, class ScalarT>
+TG_NODISCARD constexpr sphere<D, ScalarT> bounding_sphere_of(triangle<D, ScalarT> const& t)
+{
+    auto e0sqr = distance_sqr(t.pos0, t.pos1);
+    auto e1sqr = distance_sqr(t.pos1, t.pos2);
+    auto e2sqr = distance_sqr(t.pos2, t.pos0);
+
+    auto const* a = &t.pos2;
+    auto const* b = &t.pos0;
+    auto const* c = &t.pos1;
+
+    // "rotate" triangle s.t. e0 is the largest edge and
+    // lies opposite of point a (equivalently, e1 <-> b, e2 <-> c)
+    if (e2sqr > e1sqr)
+    {
+        // make e1 larger than e2
+        detail::swap(e1sqr, e2sqr);
+        detail::swap(b, c);
+    }
+    if (e1sqr > e0sqr)
+    {
+        // make e0 larger than e1
+        detail::swap(e0sqr, e1sqr);
+        detail::swap(a, b);
+    }
+
+    if (e0sqr >= e1sqr + e2sqr) // triangle is obtuse or right
+    {
+        auto radius = tg::sqrt(e0sqr) / ScalarT(2);
+        auto center = mix(*b, *c, ScalarT(0.5));
+        return {center, radius};
+    }
+
+    // compute barycentric coords
+    auto alpha = e0sqr * (e1sqr + e2sqr - e0sqr);
+    auto beta = e1sqr * (e2sqr + e0sqr - e1sqr);
+    auto gamma = e2sqr * (e0sqr + e1sqr - e2sqr);
+    auto bary_sum = alpha + beta + gamma;
+
+    auto center = tg::interpolate(triangle<D, ScalarT>(*a, *b, *c), alpha, beta, gamma) / bary_sum;
+    auto radius = tg::distance(center, t.pos0);
+
+    return {center, radius};
+}
+}