Fix for bounding volume handling and calculation

- Fixed Ritter algorithm implementation
- Added notation of invalid bounding sphere (radius == -1.0)
- Handle merging of invalid bounding sphere with valid ones
- Added test cases and adjusted tests boundingsphere and
  proximityfilter
- This is necessary to ensure the correct working for viewAll and
  viewEntity

Task-number: QTBUG-78313
Change-Id: I1dc6d227cf9009f6fbd3230093c7a7a94fb05ae3
Reviewed-by: Paul Lemire <paul.lemire@kdab.com>

3 files changed, 107 insertions, 34 deletions

diff --git a/src/render/frontend/sphere.cpp b/src/render/frontend/sphere.cpp
index 4909acaef..470dbfe59 100644
--- a/src/render/frontend/sphere.cpp
+++ b/src/render/frontend/sphere.cpp
@@ -44,6 +44,7 @@
 #include <QPair>
 
 #include <math.h>
+#include <algorithm>
 
 QT_BEGIN_NAMESPACE
 
@@ -55,6 +56,9 @@ namespace {
 // returns true and intersection point q; false otherwise
 bool intersectRaySphere(const Qt3DRender::RayCasting::QRay3D &ray, const Qt3DRender::Render::Sphere &s, Vector3D *q = nullptr)
 {
+    if (s.isNull())
+        return false;
+
     const Vector3D p = ray.origin();
     const Vector3D d = ray.direction();
     const Vector3D m = p - s.center();
@@ -139,11 +143,31 @@ inline void sphereFromExtremePoints(Qt3DRender::Render::Sphere &s, const QVector
 
 inline void constructRitterSphere(Qt3DRender::Render::Sphere &s, const QVector<Vector3D> &points)
 {
-    // Calculate the sphere encompassing two axially extreme points
-    sphereFromExtremePoints(s, points);
-
-    // Now make sure the sphere bounds all points by growing if needed
-    s.expandToContain(points);
+    //def bounding_sphere(points):
+    //  dist = lambda a,b: ((a[0] - b[0])**2 + (a[1] - b[1])**2 + (a[2] - b[2])**2)**0.5
+    //  x = points[0]
+    //  y = max(points,key= lambda p: dist(p,x) )
+    //  z = max(points,key= lambda p: dist(p,y) )
+    //  bounding_sphere = (((y[0]+z[0])/2,(y[1]+z[1])/2,(y[2]+z[2])/2), dist(y,z)/2)
+    //
+    //  exterior_points = [p for p in points if dist(p,bounding_sphere[0]) > bounding_sphere[1] ]
+    //  while ( len(exterior_points) > 0 ):
+    //    pt = exterior_points.pop()
+    //    if (dist(pt, bounding_sphere[0]) > bounding_sphere[1]):
+    //      bounding_sphere = (bounding_sphere[0],dist(pt,bounding_sphere[0]))
+    //
+    //  return bounding_sphere
+
+    const Vector3D x = points[0];
+    const Vector3D y = *std::max_element(points.begin(), points.end(), [&x](const Vector3D& lhs, const Vector3D& rhs){ return (lhs - x).lengthSquared() < (rhs - x).lengthSquared(); });
+    const Vector3D z = *std::max_element(points.begin(), points.end(), [&y](const Vector3D& lhs, const Vector3D& rhs){ return (lhs - y).lengthSquared() < (rhs - y).lengthSquared(); });
+
+    const Vector3D center = (y + z) * 0.5f;
+    const Vector3D maxDistPt = *std::max_element(points.begin(), points.end(), [&center](const Vector3D& lhs, const Vector3D& rhs){ return (lhs - center).lengthSquared() < (rhs - center).lengthSquared(); });
+    const float radius = (maxDistPt - center).length();
+
+    s.setCenter(center);
+    s.setRadius(radius);
 }
 
 } // anonymous namespace
@@ -169,6 +193,12 @@ void Sphere::initializeFromPoints(const QVector<Vector3D> &points)
 
 void Sphere::expandToContain(const Vector3D &p)
 {
+    if (isNull()) {
+        m_center = p;
+        m_radius = 0.0f;
+        return;
+    }
+
     const Vector3D d = p - m_center;
     const float dist2 = d.lengthSquared();
 
@@ -184,6 +214,13 @@ void Sphere::expandToContain(const Vector3D &p)
 
 void Sphere::expandToContain(const Sphere &sphere)
 {
+    if (isNull()) {
+        *this = sphere;
+        return;
+    } else if (sphere.isNull()) {
+        return;
+    }
+
     const Vector3D d(sphere.m_center - m_center);
     const float dist2 = d.lengthSquared();
 
@@ -206,6 +243,9 @@ void Sphere::expandToContain(const Sphere &sphere)
 
 Sphere Sphere::transformed(const Matrix4x4 &mat) const
 {
+    if (isNull())
+        return *this;
+
     // Transform extremities in x, y, and z directions to find extremities
     // of the resulting ellipsoid
     Vector3D x = mat.map(m_center + Vector3D(m_radius, 0.0f, 0.0f));
diff --git a/src/render/frontend/sphere_p.h b/src/render/frontend/sphere_p.h
index 10cf92091..b7585f85a 100644
--- a/src/render/frontend/sphere_p.h
+++ b/src/render/frontend/sphere_p.h
@@ -69,7 +69,7 @@ class Q_3DRENDERSHARED_PRIVATE_EXPORT Sphere : public RayCasting::BoundingSphere
 public:
     inline Sphere(Qt3DCore::QNodeId i = Qt3DCore::QNodeId())
         : m_center()
-        , m_radius(0.0f)
+        , m_radius(-1.0f)
         , m_id(i)
     {}
 
@@ -82,7 +82,7 @@ public:
     void setCenter(const Vector3D &c);
     Vector3D center() const override;
 
-    inline bool isNull() { return m_center == Vector3D() && m_radius == 0.0f; }
+    bool isNull() const { return m_center == Vector3D() && m_radius == -1.0f; }
 
     void setRadius(float r);
     float radius() const override;
@@ -131,7 +131,9 @@ inline Vector3D Sphere::center() const
 
 inline void Sphere::setRadius(float r)
 {
-    m_radius = r;
+    Q_ASSERT(r >= 0.0f);
+    if (r >= 0.0f)
+        m_radius = r;
 }
 
 inline float Sphere::radius() const
@@ -142,11 +144,14 @@ inline float Sphere::radius() const
 inline void Sphere::clear()
 {
     m_center = Vector3D();
-    m_radius = 0.0f;
+    m_radius = -1.0f;
 }
 
 inline bool intersects(const Sphere &a, const Sphere &b)
 {
+    if (a.isNull() || b.isNull())
+        return false;
+
     // Calculate squared distance between sphere centers
     const Vector3D d = a.center() - b.center();
     const float distSq = Vector3D::dotProduct(d, d);
diff --git a/src/render/jobs/calcboundingvolumejob.cpp b/src/render/jobs/calcboundingvolumejob.cpp
index 113dc34ce..9af2f4f38 100644
--- a/src/render/jobs/calcboundingvolumejob.cpp
+++ b/src/render/jobs/calcboundingvolumejob.cpp
@@ -90,30 +90,42 @@ public:
         m_min = QVector3D(findExtremePoints.xMin, findExtremePoints.yMin, findExtremePoints.zMin);
         m_max = QVector3D(findExtremePoints.xMax, findExtremePoints.yMax, findExtremePoints.zMax);
 
-        // Calculate squared distance for the pairs of points
-        const float xDist2 = (findExtremePoints.xMaxPt - findExtremePoints.xMinPt).lengthSquared();
-        const float yDist2 = (findExtremePoints.yMaxPt - findExtremePoints.yMinPt).lengthSquared();
-        const float zDist2 = (findExtremePoints.zMaxPt - findExtremePoints.zMinPt).lengthSquared();
-
-        // Select most distant pair
-        Vector3D p = findExtremePoints.xMinPt;
-        Vector3D q = findExtremePoints.xMaxPt;
-        if (yDist2 > xDist2 && yDist2 > zDist2) {
-            p = findExtremePoints.yMinPt;
-            q = findExtremePoints.yMaxPt;
+        FindMaxDistantPoint maxDistantPointY(m_manager);
+        maxDistantPointY.setReferencePoint = true;
+        if (!maxDistantPointY.apply(positionAttribute, indexAttribute, drawVertexCount,
+                                     primitiveRestartEnabled, primitiveRestartIndex)) {
+            return false;
         }
-        if (zDist2 > xDist2 && zDist2 > yDist2) {
-            p = findExtremePoints.zMinPt;
-            q = findExtremePoints.zMaxPt;
+        if (maxDistantPointY.hasNoPoints)
+            return false;
+
+        //const Vector3D x = maxDistantPointY.referencePt;
+        const Vector3D y = maxDistantPointY.maxDistPt;
+
+        FindMaxDistantPoint maxDistantPointZ(m_manager);
+        maxDistantPointZ.setReferencePoint = false;
+        maxDistantPointZ.referencePt = y;
+        if (!maxDistantPointZ.apply(positionAttribute, indexAttribute, drawVertexCount,
+                                     primitiveRestartEnabled, primitiveRestartIndex)) {
+            return false;
+        }
+        const Vector3D z = maxDistantPointZ.maxDistPt;
+
+        const Vector3D center = (y + z) * 0.5f;
+
+        FindMaxDistantPoint maxDistantPointCenter(m_manager);
+        maxDistantPointCenter.setReferencePoint = false;
+        maxDistantPointCenter.referencePt = center;
+        if (!maxDistantPointCenter.apply(positionAttribute, indexAttribute, drawVertexCount,
+                                     primitiveRestartEnabled, primitiveRestartIndex)) {
+            return false;
         }
 
-        const Vector3D c = 0.5f * (p + q);
-        m_volume.setCenter(c);
-        m_volume.setRadius((q - c).length());
+        const float radius = (center - maxDistantPointCenter.maxDistPt).length();
 
-        ExpandSphere expandSphere(m_manager, m_volume);
-        if (!expandSphere.apply(positionAttribute, indexAttribute, drawVertexCount,
-                                primitiveRestartEnabled, primitiveRestartIndex))
+        m_volume = Qt3DRender::Render::Sphere(center, radius);
+
+        if (m_volume.isNull())
             return false;
 
         return true;
@@ -172,18 +184,34 @@ private:
         }
     };
 
-    class ExpandSphere : public Buffer3fVisitor
+    class FindMaxDistantPoint : public Buffer3fVisitor
     {
     public:
-        ExpandSphere(NodeManagers *manager, Sphere& volume)
-            : Buffer3fVisitor(manager), m_volume(volume)
+        FindMaxDistantPoint(NodeManagers *manager)
+            : Buffer3fVisitor(manager)
         { }
 
-        Sphere& m_volume;
+        float maxLengthSquared = 0.0f;
+        Vector3D maxDistPt;
+        Vector3D referencePt;
+        bool setReferencePoint = false;
+        bool hasNoPoints = true;
+
         void visit(uint ndx, float x, float y, float z) override
         {
             Q_UNUSED(ndx);
-            m_volume.expandToContain(Vector3D(x, y, z));
+            const Vector3D p = Vector3D(x, y, z);
+
+            if (hasNoPoints && setReferencePoint) {
+                maxLengthSquared = 0.0f;
+                referencePt = p;
+            }
+            const float lengthSquared = (p - referencePt).lengthSquared();
+            if ( lengthSquared >= maxLengthSquared ) {
+                maxDistPt = p;
+                maxLengthSquared = lengthSquared;
+            }
+            hasNoPoints = false;
         }
     };
 };