Loading...
Searching...
No Matches
AABBTree.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015 OpenFOAM Foundation
9 Copyright (C) 2016-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
29#include "AABBTree.H"
30#include "bitSet.H"
31
32// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
33
34template<class Type>
36(
37 const bool leavesOnly,
38 const bool writeLinesOnly,
39 const treeBoundBox& bb,
40 const label nodeI,
41 const List<Pair<treeBoundBox>>& bbs,
42 const List<Pair<label>>& nodes,
43 label& vertI,
44 Ostream& os
45) const
46{
47 if (!leavesOnly || nodeI < 0)
48 {
49 AABBTreeBase::writeOBJ(os, bb, vertI, writeLinesOnly);
50 }
51
52 // recurse to find leaves
53 if (nodeI >= 0)
54 {
56 (
57 leavesOnly,
58 writeLinesOnly,
59 bbs[nodeI].first(),
60 nodes[nodeI].first(),
61 bbs,
62 nodes,
63 vertI,
64 os
65 );
67 (
68 leavesOnly,
69 writeLinesOnly,
70 bbs[nodeI].second(),
71 nodes[nodeI].second(),
72 bbs,
73 nodes,
74 vertI,
76 );
77 }
78}
79
80
81template<class Type>
83(
84 const bool equalBinSize,
85 const label level,
86 const UList<Type>& objects,
87 const pointField& points,
88 const labelUList& objectIDs,
89 const treeBoundBox& bb,
90 const label nodeI,
91
92 DynamicList<Pair<treeBoundBox>>& bbs,
93 DynamicList<labelPair>& nodes,
94 DynamicList<labelList>& addressing
95) const
96{
97 // Determine which direction to divide the box
98
99 const direction maxDir = bb.maxDir();
100 const scalar maxSpan = bb.span()[maxDir];
101
102 scalar pivotValue;
103
104 if (equalBinSize)
105 {
106 // Pick up points used by this set of objects
107
108 bitSet isUsedPoint(points.size());
109 DynamicList<scalar> component(points.size());
110
111 for (const label objI : objectIDs)
112 {
113 const Type& obj = objects[objI];
114
115 for (const label pointI : obj)
116 {
117 if (isUsedPoint.set(pointI))
118 {
119 component.push_back(points[pointI][maxDir]);
120 }
121 }
122 }
123
124 // Determine the median
125
126 Foam::sort(component);
127
128 pivotValue = component[component.size()/2];
129 }
130 else
131 {
132 // Geometric middle
133 pivotValue = bb.min()[maxDir] + 0.5*maxSpan;
134 }
135
136
137 const scalar divMin = pivotValue + tolerance_*maxSpan;
138 const scalar divMax = pivotValue - tolerance_*maxSpan;
139
140
141 // Assign the objects to min or max bin
142
143 DynamicList<label> minBinObjectIDs(objectIDs.size());
144 DynamicList<label> maxBinObjectIDs(objectIDs.size());
145
146 treeBoundBox minBb;
147 treeBoundBox maxBb;
148
149 for (const label objI : objectIDs)
150 {
151 const Type& obj = objects[objI];
152
153 bool addMin = false;
154 bool addMax = false;
155
156 for (const label pointi : obj)
157 {
158 const scalar& cmptValue = points[pointi][maxDir];
159
160 addMin = addMin || (cmptValue < divMin);
161 addMax = addMax || (cmptValue > divMax);
162
163 if (addMin && addMax) break;
164 }
165
166 // Note: object is inserted into both min/max child boxes (duplicated)
167 // if it crosses the bin boundaries
168 if (addMin)
169 {
170 minBinObjectIDs.push_back(objI);
171 minBb.add(points, obj);
172 }
173 if (addMax)
174 {
175 maxBinObjectIDs.push_back(objI);
176 maxBb.add(points, obj);
177 }
178 }
179
180 // Inflate box in case geometry reduces to 2-D
181 if (minBinObjectIDs.size())
182 {
183 minBb.inflate(0.01);
184 }
185 if (maxBinObjectIDs.size())
186 {
187 maxBb.inflate(0.01);
188 }
189
190 minBinObjectIDs.shrink();
191 maxBinObjectIDs.shrink();
192
193 bool addMin = (minBinObjectIDs.size() > minLeafSize_ && level < maxLevel_);
194 bool addMax = (maxBinObjectIDs.size() > minLeafSize_ && level < maxLevel_);
195
196 // Since bounding boxes overlap, verify that splitting was effective
197
198 if
199 (
200 objectIDs.size() <= (minBinObjectIDs.size() + minLeafSize_/2)
201 || objectIDs.size() <= (maxBinObjectIDs.size() + minLeafSize_/2)
202 )
203 {
204 addMin = addMax = false;
205 }
206
207 label minI;
208 if (addMin)
209 {
210 // New leaf
211 minI = nodes.size();
212 nodes.emplace_back(-1, -1);
213 }
214 else
215 {
216 // Update existing leaf
217 minI = -addressing.size() - 1;
218 addressing.push_back(minBinObjectIDs);
219 }
220
221 label maxI;
222 if (addMax)
223 {
224 // New leaf
225 maxI = nodes.size();
226 nodes.emplace_back(-1, -1);
227 }
228 else
229 {
230 // Update existing leaf
231 maxI = -addressing.size() - 1;
232 addressing.push_back(maxBinObjectIDs);
233 }
234
235 nodes(nodeI) = labelPair(minI, maxI);
236 bbs(nodeI) = Pair<treeBoundBox>(minBb, maxBb);
237
238 // Recurse
239 if (minI >= 0)
240 {
241 createBoxes
242 (
243 equalBinSize,
244 level + 1,
245 objects,
246 points,
247 minBinObjectIDs,
248 minBb,
249 minI,
250 bbs,
251 nodes,
252 addressing
253 );
254 }
255 if (maxI >= 0)
256 {
257 createBoxes
258 (
259 equalBinSize,
260 level + 1,
261 objects,
262 points,
263 maxBinObjectIDs,
264 maxBb,
265 maxI,
266 bbs,
267 nodes,
268 addressing
269 );
271}
272
273
274// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
275
276template<class Type>
278:
279 maxLevel_(0),
281 boundBoxes_(),
283{}
284
285
286template<class Type>
288(
289 const UList<Type>& objects,
290 const pointField& points,
291 const bool equalBinSize,
292 label maxLevel,
293 label minLeafSize
294)
295:
296 maxLevel_(maxLevel),
297 minLeafSize_(minLeafSize),
298 boundBoxes_(),
299 addressing_()
300{
301 if (objects.empty())
302 {
303 return;
304 }
305
306
307 DynamicList<Pair<treeBoundBox>> bbs(maxLevel);
308 DynamicList<labelPair> nodes(maxLevel);
309 DynamicList<labelList> addr(maxLevel);
310
311 nodes.emplace_back(-1, -1);
312 treeBoundBox topBb(points);
313 topBb.inflate(0.01);
314
315 labelList objectIDs(identity(objects.size()));
316
318 (
319 equalBinSize,
320 0, // starting at top level
321 objects,
322 points,
323 objectIDs,
324 topBb,
325 0, // starting node
326
327 bbs,
328 nodes,
329 addr
330 );
331
332
333 //{
334 // OFstream os("tree.obj");
335 // label vertI = 0;
336 // writeOBJ
337 // (
338 // true, // leavesOnly
339 // false, // writeLinesOnly
340 //
341 // topBb,
342 // 0,
343 // bbs,
344 // nodes,
345 // vertI,
346 // os
347 // );
348 //}
349
350
351 // transfer flattened tree to persistent storage
353 DynamicList<labelList> addressing(2*addr.size());
354
355 forAll(nodes, nodeI)
356 {
357 if (nodes[nodeI].first() < 0)
358 {
359 boundBoxes.push_back(bbs[nodeI].first());
360 addressing.push_back(addr[-(nodes[nodeI].first() + 1)]);
361 }
362 if (nodes[nodeI].second() < 0)
363 {
364 boundBoxes.push_back(bbs[nodeI].second());
365 addressing.push_back(addr[-(nodes[nodeI].second() + 1)]);
366 }
367 }
368
369 boundBoxes_.transfer(boundBoxes);
370 addressing_.transfer(addressing);
371
372
373 if (0)
374 {
375 bitSet checked(objects.size());
376 for (const auto& box : addressing_)
377 {
378 for (const auto& id : box)
379 {
380 checked.set(id);
381 }
382 }
383
384 const label unsetSize = checked.count(false);
385
386 if (unsetSize)
387 {
388 Info<< "*** Problem: IDs not set: " << unsetSize << endl;
389 }
391}
392
393
394// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
395
396template<class Type>
397void Foam::AABBTree<Type>::writeOBJ(Ostream& os) const
398{
399 label vertIndex(0);
400
401 for (const treeBoundBox& bb : boundBoxes_)
402 {
403 // writeLinesOnly=false
404 AABBTreeBase::writeOBJ(os, bb, vertIndex, false);
405 }
406}
407
408
409template<class Type>
410bool Foam::AABBTree<Type>::pointInside(const point& pt) const
411{
412 for (const treeBoundBox& bb : boundBoxes_)
413 {
414 if (bb.contains(pt))
415 {
416 return true;
417 }
419
420 return false;
421}
422
423
424template<class Type>
425bool Foam::AABBTree<Type>::overlaps(const boundBox& bbIn) const
426{
427 for (const treeBoundBox& bb : boundBoxes_)
428 {
429 if (bb.overlaps(bbIn))
430 {
431 return true;
432 }
433 }
434
435 return false;
436}
437
438
439// * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * //
440
441template<class Type>
442Foam::Ostream& Foam::operator<<(Ostream& os, const AABBTree<Type>& tree)
443{
444 os << tree.boundBoxes_ << tree.addressing_;
445
447 return os;
448}
449
450
451template<class Type>
452Foam::Istream& Foam::operator>>(Istream& is, AABBTree<Type>& tree)
453{
454 is >> tree.boundBoxes_ >> tree.addressing_;
455
457 return is;
458}
459
460
461// ************************************************************************* //
static scalar tolerance_
Relative tolerance.
Definition AABBTree.H:80
static void writeOBJ(Ostream &os, const treeBoundBox &bb, label &vertIndex, const bool writeLinesOnly=false)
Write treeBoundBox in OBJ format.
Templated tree of axis-aligned bounding boxes (AABB).
Definition AABBTree.H:116
const List< labelList > & addressing() const noexcept
Return the contents addressing.
Definition AABBTree.H:217
void createBoxes(const bool equalBinSize, const label level, const UList< Type > &objects, const pointField &points, const labelUList &objectIDs, const treeBoundBox &bb, const label nodeI, DynamicList< Pair< treeBoundBox > > &bbs, DynamicList< labelPair > &nodes, DynamicList< labelList > &addressing) const
Create the bounding boxes by interrogating points.
Definition AABBTree.C:76
bool overlaps(const boundBox &bbIn) const
Determine whether a bounding box overlaps the tree bounding boxes.
Definition AABBTree.C:418
label minLeafSize_
Minimum points per leaf.
Definition AABBTree.H:129
AABBTree()
Default construct.
Definition AABBTree.C:270
const List< treeBoundBox > & boundBoxes() const noexcept
Return the bounding boxes making up the tree.
Definition AABBTree.H:209
void writeOBJ(const bool leavesOnly, const bool writeLinesOnly, const treeBoundBox &bb, const label nodeI, const List< Pair< treeBoundBox > > &bbs, const List< Pair< label > > &nodes, label &vertI, Ostream &os) const
Write OBJ for all bounding boxes.
Definition AABBTree.C:29
label maxLevel_
Maximum tree level.
Definition AABBTree.H:124
bool pointInside(const point &pt) const
Determine whether a point is inside the bounding boxes.
Definition AABBTree.C:403
List< labelList > addressing_
Leaf addressing.
Definition AABBTree.H:139
List< treeBoundBox > boundBoxes_
Bounding boxes making up the tree.
Definition AABBTree.H:134
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition DynamicList.H:68
T & emplace_back(Args &&... args)
Construct an element at the end of the list, return reference to the new list element.
void push_back(const T &val)
Copy append an element to the end of this list.
DynamicList< T, SizeMin > & shrink()
Calls shrink_to_fit() and returns a reference to the DynamicList.
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition IOstream.C:45
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition Istream.H:60
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition List.H:72
void push_back(const T &val)
Append an element at the end of the list.
Definition ListI.H:221
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition Ostream.H:59
An ordered pair of two objects of type <T> with first() and second() elements.
Definition Pair.H:66
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition UList.H:89
bool empty() const noexcept
True if List is empty (ie, size() is zero).
Definition UList.H:701
void size(const label n)
Older name for setAddressableSize.
Definition UList.H:118
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition bitSet.H:61
unsigned int count(const bool on=true) const
Count number of bits set.
Definition bitSetI.H:420
void set(const bitSet &bitset)
Set specified bits from another bitset.
Definition bitSetI.H:502
A bounding box defined in terms of min/max extrema points.
Definition boundBox.H:71
const point & min() const noexcept
Minimum describing the bounding box.
Definition boundBoxI.H:162
void add(const boundBox &bb)
Extend to include the second box.
Definition boundBoxI.H:323
direction maxDir() const
Direction (X/Y/Z) of the largest span (for empty box: 0).
Definition boundBoxI.H:254
void inflate(const scalar factor)
Expand box by factor*mag(span) in all dimensions.
Definition boundBoxI.H:381
vector span() const
The bounding box span (from minimum to maximum).
Definition boundBoxI.H:192
Standard boundBox with extra functionality for use in octree.
OBJstream os(runTime.globalPath()/outputName)
const pointField & points
#define FUNCTION_NAME
Pair< label > labelPair
A pair of labels.
Definition Pair.H:54
List< label > labelList
A List of labels.
Definition List.H:62
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
messageStream Info
Information stream (stdout output on master, null elsewhere).
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces).
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition Ostream.H:519
Istream & operator>>(Istream &, directionInfo &)
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
void sort(UList< T > &list)
Sort the list.
Definition UList.C:283
uint8_t direction
Definition direction.H:49
vector point
Point is a vector.
Definition point.H:37
vectorField pointField
pointField is a vectorField.
UList< label > labelUList
A UList of labels.
Definition UList.H:75
Tree tree(triangles.begin(), triangles.end())
#define forAll(list, i)
Loop across all elements in list.
Definition stdFoam.H:299