1// Copyright (C) 2014 Klaralvdalens Datakonsult AB (KDAB).
2// SPDX-License-Identifier: LicenseRef-Qt-Commercial OR LGPL-3.0-only OR GPL-2.0-only OR GPL-3.0-only
3
4#include "sphere_p.h"
5
6#include <Qt3DRender/private/qray3d_p.h>
7
8#include <QPair>
9
10#include <math.h>
11#include <algorithm>
12
13QT_BEGIN_NAMESPACE
14
15namespace {
16
17// Algorithms taken from Real-time collision detection, p178-179
18
19// Intersects ray r = p + td, |d| = 1, with sphere s and, if intersecting,
20// returns true and intersection point q; false otherwise
21bool intersectRaySphere(const Qt3DRender::RayCasting::QRay3D &ray, const Qt3DRender::Render::Sphere &s, Vector3D *q = nullptr)
22{
23 if (s.isNull())
24 return false;
25
26 const Vector3D p = ray.origin();
27 const Vector3D d = ray.direction();
28 const Vector3D m = p - s.center();
29 const float c = Vector3D::dotProduct(a: m, b: m) - s.radius() * s.radius();
30
31 // If there is definitely at least one real root, there must be an intersection
32 if (q == nullptr && c <= 0.0f)
33 return true;
34
35 const float b = Vector3D::dotProduct(a: m, b: d);
36 // Exit if r’s origin outside s (c > 0) and r pointing away from s (b > 0)
37 if (c > 0.0f && b > 0.0f)
38 return false;
39
40 const float discr = b*b - c;
41 // A negative discriminant corresponds to ray missing sphere
42 if (discr < 0.0f)
43 return false;
44
45 // If we don't need the intersection point, return early
46 if (q == nullptr)
47 return true;
48
49 // Ray now found to intersect sphere, compute smallest t value of intersection
50 float t = -b - sqrt(x: discr);
51
52 // If t is negative, ray started inside sphere so clamp t to zero
53 if (t < 0.0f)
54 t = 0.0f;
55
56 *q = p + t * d;
57 return true;
58}
59
60inline void constructRitterSphere(Qt3DRender::Render::Sphere &s, const QList<Vector3D> &points)
61{
62 //def bounding_sphere(points):
63 // dist = lambda a,b: ((a[0] - b[0])**2 + (a[1] - b[1])**2 + (a[2] - b[2])**2)**0.5
64 // x = points[0]
65 // y = max(points,key= lambda p: dist(p,x) )
66 // z = max(points,key= lambda p: dist(p,y) )
67 // bounding_sphere = (((y[0]+z[0])/2,(y[1]+z[1])/2,(y[2]+z[2])/2), dist(y,z)/2)
68 //
69 // exterior_points = [p for p in points if dist(p,bounding_sphere[0]) > bounding_sphere[1] ]
70 // while ( len(exterior_points) > 0 ):
71 // pt = exterior_points.pop()
72 // if (dist(pt, bounding_sphere[0]) > bounding_sphere[1]):
73 // bounding_sphere = (bounding_sphere[0],dist(pt,bounding_sphere[0]))
74 //
75 // return bounding_sphere
76
77 const Vector3D x = points[0];
78 const Vector3D y = *std::max_element(first: points.begin(), last: points.end(), comp: [&x](const Vector3D& lhs, const Vector3D& rhs){ return (lhs - x).lengthSquared() < (rhs - x).lengthSquared(); });
79 const Vector3D z = *std::max_element(first: points.begin(), last: points.end(), comp: [&y](const Vector3D& lhs, const Vector3D& rhs){ return (lhs - y).lengthSquared() < (rhs - y).lengthSquared(); });
80
81 const Vector3D center = (y + z) * 0.5f;
82 const Vector3D maxDistPt = *std::max_element(first: points.begin(), last: points.end(), comp: [&center](const Vector3D& lhs, const Vector3D& rhs){ return (lhs - center).lengthSquared() < (rhs - center).lengthSquared(); });
83 const float radius = (maxDistPt - center).length();
84
85 s.setCenter(center);
86 s.setRadius(radius);
87}
88
89} // anonymous namespace
90
91namespace Qt3DRender {
92
93namespace Render {
94
95const float Sphere::ms_epsilon = 1.0e-7f;
96
97Sphere Sphere::fromPoints(const QList<Vector3D> &points)
98{
99 Sphere s;
100 s.initializeFromPoints(points);
101 return s;
102}
103
104void Sphere::initializeFromPoints(const QList<Vector3D> &points)
105{
106 if (!points.isEmpty())
107 constructRitterSphere(s&: *this, points);
108}
109
110void Sphere::expandToContain(const Vector3D &p)
111{
112 if (isNull()) {
113 m_center = p;
114 m_radius = 0.0f;
115 return;
116 }
117
118 const Vector3D d = p - m_center;
119 const float dist2 = d.lengthSquared();
120
121 if (dist2 > m_radius * m_radius) {
122 // Expand radius so sphere also contains p
123 const float dist = sqrt(x: dist2);
124 const float newRadius = 0.5f * (m_radius + dist);
125 const float k = (newRadius - m_radius) / dist;
126 m_radius = newRadius;
127 m_center += k * d;
128 }
129}
130
131void Sphere::expandToContain(const Sphere &sphere)
132{
133 if (isNull()) {
134 *this = sphere;
135 return;
136 } else if (sphere.isNull()) {
137 return;
138 }
139
140 const Vector3D d(sphere.m_center - m_center);
141 const float dist2 = d.lengthSquared();
142
143 const float dr = sphere.m_radius - m_radius;
144 if (dr * dr >= dist2) {
145 // Larger sphere encloses the smaller. Set our size to the larger
146 if (m_radius > sphere.m_radius)
147 return;
148 else
149 *this = sphere;
150 } else {
151 // The spheres are overlapping or disjoint
152 const float dist = sqrt(x: dist2);
153 const float newRadius = 0.5f * (dist + m_radius + sphere.m_radius);
154 if (dist > ms_epsilon)
155 m_center += d * (newRadius - m_radius) / dist;
156 m_radius = newRadius;
157 }
158}
159
160Sphere Sphere::transformed(const Matrix4x4 &mat) const
161{
162 if (isNull())
163 return *this;
164
165 // Transform extremities in x, y, and z directions to find extremities
166 // of the resulting ellipsoid
167 Vector3D x = mat.map(point: m_center + Vector3D(m_radius, 0.0f, 0.0f));
168 Vector3D y = mat.map(point: m_center + Vector3D(0.0f, m_radius, 0.0f));
169 Vector3D z = mat.map(point: m_center + Vector3D(0.0f, 0.0f, m_radius));
170
171 // Transform center and find maximum radius of ellipsoid
172 Vector3D c = mat.map(point: m_center);
173 float rSquared = qMax(a: qMax(a: (x - c).lengthSquared(), b: (y - c).lengthSquared()), b: (z - c).lengthSquared());
174 return Sphere(c, sqrt(x: rSquared), id());
175}
176
177Qt3DCore::QNodeId Sphere::id() const
178{
179 return m_id;
180}
181
182bool Sphere::intersects(const RayCasting::QRay3D &ray, Vector3D *q, Vector3D *uvw) const
183{
184 Q_UNUSED(uvw);
185 return intersectRaySphere(ray, s: *this, q);
186}
187
188Sphere::Type Sphere::type() const
189{
190 return RayCasting::QBoundingVolume::Sphere;
191}
192
193#ifndef QT_NO_DEBUG_STREAM
194
195QDebug operator<<(QDebug dbg, const Sphere &sphere)
196{
197 QDebugStateSaver saver(dbg);
198 dbg.nospace() << "Sphere(center("
199 << sphere.center().x() << ", " << sphere.center().y() << ", "
200 << sphere.center().z() << ") - radius(" << sphere.radius() << "))";
201 return dbg;
202}
203
204#endif
205
206} // Render
207
208} // Qt3DRender
209
210QT_END_NAMESPACE
211

source code of qt3d/src/render/frontend/sphere.cpp