OpenVDB 12.0.0
 
Loading...
Searching...
No Matches
BBox.h
Go to the documentation of this file.
1// Copyright Contributors to the OpenVDB Project
2// SPDX-License-Identifier: Apache-2.0
3
4#ifndef OPENVDB_MATH_BBOX_HAS_BEEN_INCLUDED
5#define OPENVDB_MATH_BBOX_HAS_BEEN_INCLUDED
6
7#include "Math.h" // for math::isApproxEqual() and math::Tolerance()
8#include "Vec3.h"
9#include <algorithm> // for std::min(), std::max()
10#include <cmath> // for std::abs()
11#include <iostream>
12#include <limits>
13#include <type_traits>
14
15
16namespace openvdb {
18namespace OPENVDB_VERSION_NAME {
19namespace math {
20
21/// @brief Axis-aligned bounding box
22template<typename Vec3T>
23class BBox
24{
25public:
26 using Vec3Type = Vec3T;
27 using ValueType = Vec3T;
28 using VectorType = Vec3T;
30
31 /// @brief The default constructor creates an invalid bounding box.
33 /// @brief Construct a bounding box that exactly encloses the given
34 /// minimum and maximum points.
35 BBox(const Vec3T& xyzMin, const Vec3T& xyzMax);
36 /// @brief Construct a bounding box that exactly encloses the given
37 /// minimum and maximum points.
38 /// @details If @a sorted is false, sort the points by their
39 /// @e x, @e y and @e z components.
40 BBox(const Vec3T& xyzMin, const Vec3T& xyzMax, bool sorted);
41 /// @brief Contruct a cubical bounding box from a minimum coordinate
42 /// and an edge length.
43 /// @note Inclusive for integral <b>ElementType</b>s
44 BBox(const Vec3T& xyzMin, const ElementType& length);
45
46 /// @brief Construct a bounding box that exactly encloses two points,
47 /// whose coordinates are given by an array of six values,
48 /// <i>x<sub>1</sub></i>, <i>y<sub>1</sub></i>, <i>z<sub>1</sub></i>,
49 /// <i>x<sub>2</sub></i>, <i>y<sub>2</sub></i> and <i>z<sub>2</sub></i>.
50 /// @details If @a sorted is false, sort the points by their
51 /// @e x, @e y and @e z components.
52 explicit BBox(const ElementType* xyz, bool sorted = true);
53
54 BBox(const BBox&) = default;
55 BBox& operator=(const BBox&) = default;
56
57 /// @brief Sort the mininum and maximum points of this bounding box
58 /// by their @e x, @e y and @e z components.
59 void sort();
60
61 /// @brief Return a const reference to the minimum point of this bounding box.
62 const Vec3T& min() const { return mMin; }
63 /// @brief Return a const reference to the maximum point of this bounding box.
64 const Vec3T& max() const { return mMax; }
65 /// @brief Return a non-const reference to the minimum point of this bounding box.
66 Vec3T& min() { return mMin; }
67 /// @brief Return a non-const reference to the maximum point of this bounding box.
68 Vec3T& max() { return mMax; }
69
70 /// @brief Return @c true if this bounding box is identical to the given bounding box.
71 bool operator==(const BBox& rhs) const;
72 /// @brief Return @c true if this bounding box differs from the given bounding box.
73 bool operator!=(const BBox& rhs) const { return !(*this == rhs); }
74
75 /// @brief Return @c true if this bounding box is empty, i.e., it has no (positive) volume.
76 bool empty() const;
77 /// @brief Return @c true if this bounding box has (positive) volume.
78 bool hasVolume() const { return !this->empty(); }
79 /// @brief Return @c true if this bounding box has (positive) volume.
80 operator bool() const { return !this->empty(); }
81
82 /// @brief Return @c true if all components of the minimum point are less than
83 /// or equal to the corresponding components of the maximum point.
84 /// @details This is equivalent to testing whether this bounding box has nonnegative volume.
85 /// @note For floating-point <b>ElementType</b>s a tolerance is used for this test.
86 bool isSorted() const;
87
88 /// @brief Return the center point of this bounding box.
90
91 /// @brief Return the extents of this bounding box, i.e., the length along each axis.
92 /// @note Inclusive for integral <b>ElementType</b>s
93 Vec3T extents() const;
94 /// @brief Return the index (0, 1 or 2) of the longest axis.
95 size_t maxExtent() const { return MaxIndex(mMax - mMin); }
96 /// @brief Return the index (0, 1 or 2) of the shortest axis.
97 size_t minExtent() const { return MinIndex(mMax - mMin); }
98
99 /// @brief Return the volume enclosed by this bounding box.
100 ElementType volume() const { Vec3T e = this->extents(); return e[0] * e[1] * e[2]; }
101
102 /// @brief Return @c true if the given point is inside this bounding box.
103 bool isInside(const Vec3T& xyz) const;
104 /// @brief Return @c true if the given bounding box is inside this bounding box.
105 bool isInside(const BBox&) const;
106 /// @brief Return @c true if the given bounding box overlaps with this bounding box.
107 bool hasOverlap(const BBox&) const;
108 /// @brief Return @c true if the given bounding box overlaps with this bounding box.
109 bool intersects(const BBox& other) const { return hasOverlap(other); }
110
111 /// @brief Pad this bounding box.
112 void expand(ElementType padding);
113 /// @brief Expand this bounding box to enclose the given point.
114 void expand(const Vec3T& xyz);
115 /// @brief Union this bounding box with the given bounding box.
116 void expand(const BBox&);
117 /// @brief Union this bounding box with the cubical bounding box with
118 /// minimum point @a xyzMin and the given edge length.
119 /// @note Inclusive for integral <b>ElementType</b>s
120 void expand(const Vec3T& xyzMin, const ElementType& length);
121
122 /// @brief Translate this bounding box by
123 /// (<i>t<sub>x</sub></i>, <i>t<sub>y</sub></i>, <i>t<sub>z</sub></i>).
124 void translate(const Vec3T& t);
125
126 /// @brief Apply a map to this bounding box.
127 template<typename MapType>
128 BBox applyMap(const MapType& map) const;
129 /// @brief Apply the inverse of a map to this bounding box
130 template<typename MapType>
131 BBox applyInverseMap(const MapType& map) const;
132
133 /// @brief Unserialize this bounding box from the given stream.
134 void read(std::istream& is) { mMin.read(is); mMax.read(is); }
135 /// @brief Serialize this bounding box to the given stream.
136 void write(std::ostream& os) const { mMin.write(os); mMax.write(os); }
137
138private:
139 Vec3T mMin, mMax;
140}; // class BBox
141
142
143////////////////////////////////////////
144
145
146template<typename Vec3T>
147inline
149 mMin( std::numeric_limits<ElementType>::max()),
150 mMax(-std::numeric_limits<ElementType>::max())
151{
152}
153
154template<typename Vec3T>
155inline
156BBox<Vec3T>::BBox(const Vec3T& xyzMin, const Vec3T& xyzMax):
157 mMin(xyzMin), mMax(xyzMax)
158{
159}
160
161template<typename Vec3T>
162inline
163BBox<Vec3T>::BBox(const Vec3T& xyzMin, const Vec3T& xyzMax, bool sorted):
164 mMin(xyzMin), mMax(xyzMax)
165{
166 if (!sorted) this->sort();
167}
168
169template<typename Vec3T>
170inline
171BBox<Vec3T>::BBox(const Vec3T& xyzMin, const ElementType& length):
172 mMin(xyzMin), mMax(xyzMin)
173{
174 // min and max are inclusive for integral ElementType
175 const ElementType size = std::is_integral<ElementType>::value ? length-1 : length;
176 mMax[0] += size;
177 mMax[1] += size;
178 mMax[2] += size;
179}
180
181template<typename Vec3T>
182inline
183BBox<Vec3T>::BBox(const ElementType* xyz, bool sorted):
184 mMin(xyz[0], xyz[1], xyz[2]),
185 mMax(xyz[3], xyz[4], xyz[5])
186{
187 if (!sorted) this->sort();
188}
189
190
191////////////////////////////////////////
192
193
194template<typename Vec3T>
195inline bool
197{
198 if (std::is_integral<ElementType>::value) {
199 // min and max are inclusive for integral ElementType
200 return (mMin[0] > mMax[0] || mMin[1] > mMax[1] || mMin[2] > mMax[2]);
201 }
202 return mMin[0] >= mMax[0] || mMin[1] >= mMax[1] || mMin[2] >= mMax[2];
203}
204
205
206template<typename Vec3T>
207inline bool
209{
210 if (std::is_integral<ElementType>::value) {
211 return mMin == rhs.min() && mMax == rhs.max();
212 } else {
213 return math::isApproxEqual(mMin, rhs.min()) && math::isApproxEqual(mMax, rhs.max());
214 }
215}
216
217
218template<typename Vec3T>
219inline void
221{
222 Vec3T tMin(mMin), tMax(mMax);
223 for (int i = 0; i < 3; ++i) {
224 mMin[i] = std::min(tMin[i], tMax[i]);
225 mMax[i] = std::max(tMin[i], tMax[i]);
226 }
227}
228
229
230template<typename Vec3T>
231inline bool
233{
234 if (std::is_integral<ElementType>::value) {
235 return (mMin[0] <= mMax[0] && mMin[1] <= mMax[1] && mMin[2] <= mMax[2]);
236 } else {
238 return (mMin[0] < (mMax[0] + t) && mMin[1] < (mMax[1] + t) && mMin[2] < (mMax[2] + t));
239 }
240}
241
242
243template<typename Vec3T>
244inline Vec3d
246{
247 return (Vec3d(mMin.asPointer()) + Vec3d(mMax.asPointer())) * 0.5;
248}
249
250
251template<typename Vec3T>
252inline Vec3T
254{
255 if (std::is_integral<ElementType>::value) {
256 return (mMax - mMin) + Vec3T(1, 1, 1);
257 } else {
258 return (mMax - mMin);
259 }
260}
261
262////////////////////////////////////////
263
264
265template<typename Vec3T>
266inline bool
267BBox<Vec3T>::isInside(const Vec3T& xyz) const
268{
269 if (std::is_integral<ElementType>::value) {
270 return xyz[0] >= mMin[0] && xyz[0] <= mMax[0] &&
271 xyz[1] >= mMin[1] && xyz[1] <= mMax[1] &&
272 xyz[2] >= mMin[2] && xyz[2] <= mMax[2];
273 } else {
275 return xyz[0] > (mMin[0]-t) && xyz[0] < (mMax[0]+t) &&
276 xyz[1] > (mMin[1]-t) && xyz[1] < (mMax[1]+t) &&
277 xyz[2] > (mMin[2]-t) && xyz[2] < (mMax[2]+t);
278 }
279}
280
281
282template<typename Vec3T>
283inline bool
285{
286 if (std::is_integral<ElementType>::value) {
287 return b.min()[0] >= mMin[0] && b.max()[0] <= mMax[0] &&
288 b.min()[1] >= mMin[1] && b.max()[1] <= mMax[1] &&
289 b.min()[2] >= mMin[2] && b.max()[2] <= mMax[2];
290 } else {
292 return (b.min()[0]-t) > mMin[0] && (b.max()[0]+t) < mMax[0] &&
293 (b.min()[1]-t) > mMin[1] && (b.max()[1]+t) < mMax[1] &&
294 (b.min()[2]-t) > mMin[2] && (b.max()[2]+t) < mMax[2];
295 }
296}
297
298
299template<typename Vec3T>
300inline bool
302{
303 if (std::is_integral<ElementType>::value) {
304 return mMax[0] >= b.min()[0] && mMin[0] <= b.max()[0] &&
305 mMax[1] >= b.min()[1] && mMin[1] <= b.max()[1] &&
306 mMax[2] >= b.min()[2] && mMin[2] <= b.max()[2];
307 } else {
309 return mMax[0] > (b.min()[0]-t) && mMin[0] < (b.max()[0]+t) &&
310 mMax[1] > (b.min()[1]-t) && mMin[1] < (b.max()[1]+t) &&
311 mMax[2] > (b.min()[2]-t) && mMin[2] < (b.max()[2]+t);
312 }
313}
314
315
316////////////////////////////////////////
317
318
319template<typename Vec3T>
320inline void
322{
323 dx = std::abs(dx);
324 for (int i = 0; i < 3; ++i) {
325 mMin[i] -= dx;
326 mMax[i] += dx;
327 }
328}
329
330
331template<typename Vec3T>
332inline void
333BBox<Vec3T>::expand(const Vec3T& xyz)
334{
335 for (int i = 0; i < 3; ++i) {
336 mMin[i] = std::min(mMin[i], xyz[i]);
337 mMax[i] = std::max(mMax[i], xyz[i]);
338 }
339}
340
341
342template<typename Vec3T>
343inline void
345{
346 for (int i = 0; i < 3; ++i) {
347 mMin[i] = std::min(mMin[i], b.min()[i]);
348 mMax[i] = std::max(mMax[i], b.max()[i]);
349 }
350}
351
352template<typename Vec3T>
353inline void
354BBox<Vec3T>::expand(const Vec3T& xyzMin, const ElementType& length)
355{
356 const ElementType size = std::is_integral<ElementType>::value ? length-1 : length;
357 for (int i = 0; i < 3; ++i) {
358 mMin[i] = std::min(mMin[i], xyzMin[i]);
359 mMax[i] = std::max(mMax[i], xyzMin[i] + size);
360 }
361}
362
363
364template<typename Vec3T>
365inline void
367{
368 mMin += dx;
369 mMax += dx;
370}
371
372template<typename Vec3T>
373template<typename MapType>
374inline BBox<Vec3T>
375BBox<Vec3T>::applyMap(const MapType& map) const
376{
377 using Vec3R = Vec3<double>;
378 BBox<Vec3T> bbox;
379 bbox.expand(map.applyMap(Vec3R(mMin[0], mMin[1], mMin[2])));
380 bbox.expand(map.applyMap(Vec3R(mMin[0], mMin[1], mMax[2])));
381 bbox.expand(map.applyMap(Vec3R(mMin[0], mMax[1], mMin[2])));
382 bbox.expand(map.applyMap(Vec3R(mMax[0], mMin[1], mMin[2])));
383 bbox.expand(map.applyMap(Vec3R(mMax[0], mMax[1], mMin[2])));
384 bbox.expand(map.applyMap(Vec3R(mMax[0], mMin[1], mMax[2])));
385 bbox.expand(map.applyMap(Vec3R(mMin[0], mMax[1], mMax[2])));
386 bbox.expand(map.applyMap(Vec3R(mMax[0], mMax[1], mMax[2])));
387 return bbox;
388}
389
390template<typename Vec3T>
391template<typename MapType>
392inline BBox<Vec3T>
393BBox<Vec3T>::applyInverseMap(const MapType& map) const
394{
395 using Vec3R = Vec3<double>;
396 BBox<Vec3T> bbox;
397 bbox.expand(map.applyInverseMap(Vec3R(mMin[0], mMin[1], mMin[2])));
398 bbox.expand(map.applyInverseMap(Vec3R(mMin[0], mMin[1], mMax[2])));
399 bbox.expand(map.applyInverseMap(Vec3R(mMin[0], mMax[1], mMin[2])));
400 bbox.expand(map.applyInverseMap(Vec3R(mMax[0], mMin[1], mMin[2])));
401 bbox.expand(map.applyInverseMap(Vec3R(mMax[0], mMax[1], mMin[2])));
402 bbox.expand(map.applyInverseMap(Vec3R(mMax[0], mMin[1], mMax[2])));
403 bbox.expand(map.applyInverseMap(Vec3R(mMin[0], mMax[1], mMax[2])));
404 bbox.expand(map.applyInverseMap(Vec3R(mMax[0], mMax[1], mMax[2])));
405 return bbox;
406}
407
408////////////////////////////////////////
409
410
411template<typename Vec3T>
412inline std::ostream&
413operator<<(std::ostream& os, const BBox<Vec3T>& b)
414{
415 os << b.min() << " -> " << b.max();
416 return os;
417}
418
419} // namespace math
420} // namespace OPENVDB_VERSION_NAME
421} // namespace openvdb
422
423#endif // OPENVDB_MATH_BBOX_HAS_BEEN_INCLUDED
General-purpose arithmetic and comparison routines, most of which accept arbitrary value types (or at...
Axis-aligned bounding box.
Definition BBox.h:24
BBox(const Vec3T &xyzMin, const ElementType &length)
Contruct a cubical bounding box from a minimum coordinate and an edge length.
Definition BBox.h:171
bool isInside(const Vec3T &xyz) const
Return true if the given point is inside this bounding box.
Definition BBox.h:267
BBox(const Vec3T &xyzMin, const Vec3T &xyzMax)
Construct a bounding box that exactly encloses the given minimum and maximum points.
Definition BBox.h:156
BBox()
The default constructor creates an invalid bounding box.
Definition BBox.h:148
size_t minExtent() const
Return the index (0, 1 or 2) of the shortest axis.
Definition BBox.h:97
const Vec3T & max() const
Return a const reference to the maximum point of this bounding box.
Definition BBox.h:64
void write(std::ostream &os) const
Serialize this bounding box to the given stream.
Definition BBox.h:136
bool hasOverlap(const BBox &) const
Return true if the given bounding box overlaps with this bounding box.
Definition BBox.h:301
void expand(const Vec3T &xyz)
Expand this bounding box to enclose the given point.
Definition BBox.h:333
void expand(const Vec3T &xyzMin, const ElementType &length)
Union this bounding box with the cubical bounding box with minimum point xyzMin and the given edge le...
Definition BBox.h:354
void sort()
Sort the mininum and maximum points of this bounding box by their x, y and z components.
Definition BBox.h:220
ElementType volume() const
Return the volume enclosed by this bounding box.
Definition BBox.h:100
bool hasVolume() const
Return true if this bounding box has (positive) volume.
Definition BBox.h:78
Vec3d getCenter() const
Return the center point of this bounding box.
Definition BBox.h:245
BBox(const Vec3T &xyzMin, const Vec3T &xyzMax, bool sorted)
Construct a bounding box that exactly encloses the given minimum and maximum points.
Definition BBox.h:163
bool isSorted() const
Return true if all components of the minimum point are less than or equal to the corresponding compon...
Definition BBox.h:232
bool empty() const
Return true if this bounding box is empty, i.e., it has no (positive) volume.
Definition BBox.h:196
void expand(ElementType padding)
Pad this bounding box.
Definition BBox.h:321
void translate(const Vec3T &t)
Translate this bounding box by (tx, ty, tz).
Definition BBox.h:366
void expand(const BBox &)
Union this bounding box with the given bounding box.
Definition BBox.h:344
size_t maxExtent() const
Return the index (0, 1 or 2) of the longest axis.
Definition BBox.h:95
bool operator==(const BBox &rhs) const
Return true if this bounding box is identical to the given bounding box.
Definition BBox.h:208
bool isInside(const BBox &) const
Return true if the given bounding box is inside this bounding box.
Definition BBox.h:284
bool operator!=(const BBox &rhs) const
Return true if this bounding box differs from the given bounding box.
Definition BBox.h:73
Vec3d ValueType
Definition BBox.h:27
Vec3T extents() const
Return the extents of this bounding box, i.e., the length along each axis.
Definition BBox.h:253
Vec3d Vec3Type
Definition BBox.h:26
const Vec3T & min() const
Return a const reference to the minimum point of this bounding box.
Definition BBox.h:62
BBox(const ElementType *xyz, bool sorted=true)
Construct a bounding box that exactly encloses two points, whose coordinates are given by an array of...
Definition BBox.h:183
Vec3T & min()
Return a non-const reference to the minimum point of this bounding box.
Definition BBox.h:66
BBox(const BBox &)=default
bool intersects(const BBox &other) const
Return true if the given bounding box overlaps with this bounding box.
Definition BBox.h:109
Vec3T & max()
Return a non-const reference to the maximum point of this bounding box.
Definition BBox.h:68
typename Vec3Type::ValueType ElementType
Definition BBox.h:29
Vec3d VectorType
Definition BBox.h:28
BBox & operator=(const BBox &)=default
BBox applyInverseMap(const MapType &map) const
Apply the inverse of a map to this bounding box.
BBox applyMap(const MapType &map) const
Apply a map to this bounding box.
void read(std::istream &is)
Unserialize this bounding box from the given stream.
Definition BBox.h:134
Definition Vec3.h:25
double ValueType
Definition Vec3.h:28
size_t MaxIndex(const Vec3T &v)
Return the index [0,1,2] of the largest value in a 3D vector.
Definition Math.h:946
std::ostream & operator<<(std::ostream &os, const BBox< Vec3T > &b)
Definition BBox.h:413
bool isApproxEqual(const Type &a, const Type &b, const Type &tolerance)
Return true if a is equal to b to within the given tolerance.
Definition Math.h:406
size_t MinIndex(const Vec3T &v)
Return the index [0,1,2] of the smallest value in a 3D vector.
Definition Math.h:930
Vec3< double > Vec3d
Definition Vec3.h:665
math::Vec3< Real > Vec3R
Definition Types.h:72
Definition Exceptions.h:13
Definition Coord.h:590
static T value()
Definition Math.h:148
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition version.h.in:121
#define OPENVDB_USE_VERSION_NAMESPACE
Definition version.h.in:218