indexedOctree.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::indexedOctree
26 
27 Description
28  Non-pointer based hierarchical recursive searching
29 
30 SourceFiles
31  indexedOctree.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef indexedOctree_H
36 #define indexedOctree_H
37 
38 #include "treeBoundBox.H"
39 #include "pointIndexHit.H"
40 #include "FixedList.H"
41 #include "Ostream.H"
42 #include "HashSet.H"
43 #include "labelBits.H"
44 #include "PackedList.H"
45 #include "volumeType.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 // Forward declaration of classes
53 template<class Type> class indexedOctree;
54 template<class Type> Ostream& operator<<(Ostream&, const indexedOctree<Type>&);
55 class Istream;
56 
57 
58 /*---------------------------------------------------------------------------*\
59  Class indexedOctreeName Declaration
60 \*---------------------------------------------------------------------------*/
61 
62 TemplateName(indexedOctree);
63 
64 
65 /*---------------------------------------------------------------------------*\
66  Class indexedOctree Declaration
67 \*---------------------------------------------------------------------------*/
68 
69 template<class Type>
70 class indexedOctree
71 :
72  public indexedOctreeName
73 {
74 public:
75 
76  // Data types
77 
78  //- Tree node. Has up pointer and down pointers.
79  class node
80  {
81  public:
82 
83  //- Bounding box of this node
85 
86  //- Parent node (index into nodes_ of tree)
87  label parent_;
88 
89  //- IDs of the 8 nodes on all sides of the mid point
91 
92  friend Ostream& operator<< (Ostream& os, const node& n)
93  {
94  return os << n.bb_ << token::SPACE
95  << n.parent_ << token::SPACE << n.subNodes_;
96  }
97 
98  friend Istream& operator>> (Istream& is, node& n)
99  {
100  return is >> n.bb_ >> n.parent_ >> n.subNodes_;
101  }
103  friend bool operator==(const node& a, const node& b)
104  {
105  return
106  a.bb_ == b.bb_
107  && a.parent_ == b.parent_
108  && a.subNodes_ == b.subNodes_;
109  }
111  friend bool operator!=(const node& a, const node& b)
112  {
113  return !(a == b);
114  }
115  };
116 
117 
118 private:
119 
120  // Static data
121 
122  //- Relative perturbation tolerance. Determines when point is
123  // considered to be close to face/edge of bb of node.
124  // The tolerance is relative to the bounding box of the smallest
125  // node.
126  static scalar perturbTol_;
127 
128 
129  // Private Data
130 
131  //- Underlying shapes for geometric queries.
132  const Type shapes_;
133 
134  //- List of all nodes
135  List<node> nodes_;
136 
137  //- List of all contents (referenced by those nodes that are contents)
138  labelListList contents_;
139 
140  //- Per node per octant whether is fully inside/outside/mixed.
141  mutable PackedList<2> nodeTypes_;
142 
143  // Private Member Functions
144 
145  //- Helper: does bb intersect a sphere around sample? Or is any
146  // corner point of bb closer than nearestDistSqr to sample.
147  // (bb is implicitly provided as parent bb + octant)
148  static bool overlaps
149  (
150  const treeBoundBox& parentBb,
151  const direction octant,
152  const scalar nearestDistSqr,
153  const point& sample
154  );
155 
156  // Construction
157 
158  //- Split list of indices into 8 bins
159  // according to where they are in relation to mid.
160  void divide
161  (
162  const labelList& indices,
163  const treeBoundBox& bb,
164  labelListList& result
165  ) const;
166 
167  //- Subdivide the contents node at position contentI.
168  // Appends to contents.
169  node divide
170  (
171  const treeBoundBox& bb,
173  const label contentI
174  ) const;
175 
176  //- Split any contents node with more than minSize elements.
177  void splitNodes
178  (
179  const label minSize,
181  DynamicList<labelList>& contents
182  ) const;
183 
184  //- Reorder contents to be in same order as nodes.
185  // Returns number of nodes on the compactLevel.
186  static label compactContents
187  (
188  DynamicList<node>& nodes,
189  DynamicList<labelList>& contents,
190  const label compactLevel,
191  const label nodeI,
192  const label level,
193  List<labelList>& compactedContents,
194  label& compactI
195  );
196 
197  //- Determine inside/outside per node (mixed if cannot be
198  // determined). Only valid for closed shapes.
199  volumeType calcVolumeType(const label nodeI) const;
200 
201  //- Search cached volume type.
202  volumeType getVolumeType(const label nodeI, const point&) const;
203 
204 
205  // Query
206 
207  //- Find nearest point to line.
208  template<class FindNearestOp>
209  void findNearest
210  (
211  const label nodeI,
212  const linePointRef& ln,
213 
214  treeBoundBox& tightest,
215  label& nearestShapeI,
216  point& linePoint,
217  point& nearestPoint,
218 
219  const FindNearestOp& fnOp
220  ) const;
221 
222  //- Return bbox of octant
223  treeBoundBox subBbox
224  (
225  const label parentNodeI,
226  const direction octant
227  ) const;
228 
229  //- Helper: take a point on/close to face of bb and push it
230  // inside or outside of bb.
231  static point pushPoint
232  (
233  const treeBoundBox&,
234  const point&,
235  const bool pushInside
236  );
237 
238  //- Helper: take a point on face of bb and push it
239  // inside or outside of bb.
240  static point pushPoint
241  (
242  const treeBoundBox&,
243  const direction,
244  const point&,
245  const bool pushInside
246  );
247 
248  //- Helper: take point on face(s) of bb and push it away from
249  // edges of face.
250  // Guarantees that if pt is on a face it gets perturbed
251  // so it is away from the face edges.
252  // If pt is not on a face does nothing.
253  static point pushPointIntoFace
254  (
255  const treeBoundBox& bb,
256  const vector& dir, // end-start
257  const point& pt
258  );
259 
260  //- Walk to parent of node+octant.
261  bool walkToParent
262  (
263  const label nodeI,
264  const direction octant,
265 
266  label& parentNodeI,
267  label& parentOctant
268  ) const;
269 
270  //- Walk tree to neighbouring node. Return false if edge of tree
271  // hit.
272  bool walkToNeighbour
273  (
274  const point& facePoint,
275  const direction faceID, // direction to walk in
276  label& nodeI,
277  direction& octant
278  ) const;
279 
280  //- Debug: return verbose the bounding box faces
281  static word faceString(const direction faceID);
282 
283  //- Traverse a node. If intersects a triangle return first
284  // intersection point.
285  // findAny=true : return any intersection
286  // findAny=false: return nearest (to start) intersection
287  template<class FindIntersectOp>
288  void traverseNode
289  (
290  const bool findAny,
291  const point& treeStart,
292  const vector& treeVec,
293 
294  const point& start,
295  const point& end,
296  const label nodeI,
297  const direction octantI,
298 
299  pointIndexHit& hitInfo,
300  direction& faceID,
301 
302  const FindIntersectOp& fiOp
303  ) const;
304 
305  //- Find any or nearest intersection
306  template<class FindIntersectOp>
307  pointIndexHit findLine
308  (
309  const bool findAny,
310  const point& treeStart,
311  const point& treeEnd,
312  const label startNodeI,
313  const direction startOctantI,
314  const FindIntersectOp& fiOp,
315  const bool verbose = false
316  ) const;
317 
318  //- Find any or nearest intersection of line between start and end.
319  template<class FindIntersectOp>
320  pointIndexHit findLine
321  (
322  const bool findAny,
323  const point& start,
324  const point& end,
325  const FindIntersectOp& fiOp
326  ) const;
327 
328  //- Find all elements intersecting box.
329  void findBox
330  (
331  const label nodeI,
332  const treeBoundBox& searchBox,
333  labelHashSet& elements
334  ) const;
335 
336 
337  //- Find all elements intersecting sphere.
338  void findSphere
339  (
340  const label nodeI,
341  const point& centre,
342  const scalar radiusSqr,
343  labelHashSet& elements
344  ) const;
345 
346 
347  template<class CompareOp>
348  static void findNear
349  (
350  const scalar nearDist,
351  const bool okOrder,
352  const indexedOctree<Type>& tree1,
353  const labelBits index1,
354  const treeBoundBox& bb1,
355  const indexedOctree<Type>& tree2,
356  const labelBits index2,
357  const treeBoundBox& bb2,
358  CompareOp& cop
359  );
360 
361 
362  // Other
363 
364  //- Count number of elements on this and sublevels
365  label countElements(const labelBits index) const;
366 
367  //- Dump node+octant to an obj file
368  void writeOBJ(const label nodeI, const direction octant) const;
369 
370  //- From index into contents_ to subNodes_ entry
371  static labelBits contentPlusOctant
372  (
373  const label i,
374  const direction octant
375  )
376  {
377  return labelBits(-i - 1, octant);
378  }
379 
380  //- From index into nodes_ to subNodes_ entry
381  static labelBits nodePlusOctant
382  (
383  const label i,
384  const direction octant
385  )
386  {
387  return labelBits(i + 1, octant);
388  }
389 
390  //- From empty to subNodes_ entry
391  static labelBits emptyPlusOctant
392  (
393  const direction octant
394  )
395  {
396  return labelBits(0, octant);
397  }
398 
399 public:
400 
401  //- Get the perturbation tolerance
402  static scalar& perturbTol();
403 
404 
405  // Constructors
406 
407  //- Construct null
408  indexedOctree(const Type& shapes);
409 
410  //- Construct from components
412  (
413  const Type& shapes,
414  const List<node>& nodes,
415  const labelListList& contents
416  );
417 
418  //- Construct from shapes
420  (
421  const Type& shapes,
422  const treeBoundBox& bb,
423  const label maxLevels, // maximum number of levels
424  const scalar maxLeafRatio, // how many elements per leaf
425  const scalar maxDuplicity // in how many leaves is a shape on
426  // average
427  );
428 
429  //- Construct from Istream
430  indexedOctree(const Type& shapes, Istream& is);
431 
432  //- Clone
434  {
436  (
437  new indexedOctree<Type>(*this)
438  );
439  }
440 
441  // Member Functions
442 
443  // Access
444 
445  //- Reference to shape
446  const Type& shapes() const
447  {
448  return shapes_;
449  }
450 
451  //- List of all nodes
452  const List<node>& nodes() const
453  {
454  return nodes_;
455  }
456 
457  //- List of all contents (referenced by those nodes that are
458  // contents)
459  const labelListList& contents() const
460  {
461  return contents_;
462  }
463 
464  //- Top bounding box
465  const treeBoundBox& bb() const
466  {
467  if (nodes_.empty())
468  {
470  << "Tree is empty" << abort(FatalError);
471  }
472  return nodes_[0].bb_;
473  }
474 
475 
476  // Conversions for entries in subNodes_.
478  static bool isContent(const labelBits i)
479  {
480  return i.val() < 0;
481  }
483  static bool isEmpty(const labelBits i)
484  {
485  return i.val() == 0;
486  }
488  static bool isNode(const labelBits i)
489  {
490  return i.val() > 0;
491  }
493  static label getContent(const labelBits i)
494  {
495  if (!isContent(i))
496  {
498  << abort(FatalError);
499  }
500  return -i.val()-1;
501  }
503  static label getNode(const labelBits i)
504  {
505  if (!isNode(i))
506  {
508  << abort(FatalError);
509  }
510  return i.val() - 1;
511  }
513  static direction getOctant(const labelBits i)
514  {
515  return i.bits();
516  }
517 
518 
519  // Queries
520 
521  pointIndexHit findNearest
522  (
523  const point& sample,
524  const scalar nearestDistSqr
525  ) const;
526 
527  //- Calculate nearest point on nearest shape.
528  // Returns
529  // - bool : any point found nearer than nearestDistSqr
530  // - label: index in shapes
531  // - point: actual nearest point found
532  template<class FindNearestOp>
533  pointIndexHit findNearest
534  (
535  const point& sample,
536  const scalar nearestDistSqr,
537 
538  const FindNearestOp& fnOp
539  ) const;
540 
541  //- Low level: calculate nearest starting from subnode.
542  template<class FindNearestOp>
543  void findNearest
544  (
545  const label nodeI,
546  const point&,
547 
548  scalar& nearestDistSqr,
549  label& nearestShapeI,
550  point& nearestPoint,
551 
552  const FindNearestOp& fnOp
553  ) const;
554 
555  //- Find nearest to line.
556  // Returns
557  // - bool : any point found?
558  // - label: index in shapes
559  // - point: actual nearest point found
560  // sets:
561  // - linePoint : corresponding nearest point on line
562  pointIndexHit findNearest
563  (
564  const linePointRef& ln,
565  treeBoundBox& tightest,
566  point& linePoint
567  ) const;
568 
569  template<class FindNearestOp>
570  pointIndexHit findNearest
571  (
572  const linePointRef& ln,
573  treeBoundBox& tightest,
574  point& linePoint,
575 
576  const FindNearestOp& fnOp
577  ) const;
578 
579  //- Find nearest intersection of line between start and end.
580  pointIndexHit findLine
581  (
582  const point& start,
583  const point& end
584  ) const;
585 
586  //- Find any intersection of line between start and end.
588  (
589  const point& start,
590  const point& end
591  ) const;
592 
593  //- Find nearest intersection of line between start and end.
594  template<class FindIntersectOp>
595  pointIndexHit findLine
596  (
597  const point& start,
598  const point& end,
599  const FindIntersectOp& fiOp
600  ) const;
601 
602  //- Find any intersection of line between start and end.
603  template<class FindIntersectOp>
605  (
606  const point& start,
607  const point& end,
608  const FindIntersectOp& fiOp
609  ) const;
610 
611  //- Find (in no particular order) indices of all shapes inside or
612  // overlapping bounding box (i.e. all shapes not outside box)
613  labelList findBox(const treeBoundBox& bb) const;
614 
615  //- Find (in no particular order) indices of all shapes inside or
616  // overlapping a bounding sphere (i.e. all shapes not outside
617  // sphere)
618  labelList findSphere
619  (
620  const point& centre,
621  const scalar radiusSqr
622  ) const;
623 
624  //- Find deepest node (as parent+octant) containing point. Starts
625  // off from starting index in nodes_ (use 0 to start from top)
626  // Use getNode and getOctant to extract info, or call findIndices.
627  labelBits findNode(const label nodeI, const point&) const;
628 
629  //- Find shape containing point. Only implemented for certain
630  // shapes.
631  label findInside(const point&) const;
632 
633  //- Find the shape indices that occupy the result of findNode
634  const labelList& findIndices(const point&) const;
635 
636  //- Determine type (inside/outside/mixed) for point. unknown if
637  // cannot be determined (e.g. non-manifold surface)
638  volumeType getVolumeType(const point&) const;
639 
640  //- Helper function to return the side. Returns outside if
641  // outsideNormal&vec >= 0, inside otherwise
642  static volumeType getSide
643  (
644  const vector& outsideNormal,
645  const vector& vec
646  );
647 
648  //- Helper: does bb intersect a sphere around sample? Or is any
649  // corner point of bb closer than nearestDistSqr to sample.
650  static bool overlaps
651  (
652  const point& bbMin,
653  const point& bbMax,
654  const scalar nearestDistSqr,
655  const point& sample
656  );
657 
658  //- Find near pairs and apply CompareOp to them.
659  // tree2 can be *this or different tree.
660  template<class CompareOp>
661  void findNear
662  (
663  const scalar nearDist,
664  const indexedOctree<Type>& tree2,
665  CompareOp& cop
666  ) const;
667 
668 
669  // Write
670 
671  //- Print tree. Either print all indices (printContent = true) or
672  // just size of contents nodes.
673  void print
674  (
676  const bool printContents,
677  const label
678  ) const;
679 
680  bool write(Ostream& os) const;
681 
682 
683  // IOstream Operators
684 
685  friend Ostream& operator<< <Type>(Ostream&, const indexedOctree<Type>&);
686 };
687 
688 
689 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
690 
691 } // End namespace Foam
692 
693 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
694 
695 #ifdef NoRepository
696  #include "indexedOctree.C"
697 #endif
698 
699 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
700 
701 #endif
702 
703 // ************************************************************************* //
void divide(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
TemplateName(blendedSchemeBase)
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:325
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A line primitive.
Definition: line.H:56
treeBoundBox bb_
Bounding box of this node.
Definition: indexedOctree.H:83
static direction getOctant(const labelBits i)
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:54
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
pointIndexHit findLineAny(const point &start, const point &end) const
Find any intersection of line between start and end.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
static bool isEmpty(const labelBits i)
uint8_t direction
Definition: direction.H:45
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
const labelListList & contents() const
List of all contents (referenced by those nodes that are.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
static scalar & perturbTol()
Get the perturbation tolerance.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
label parent_
Parent node (index into nodes_ of tree)
Definition: indexedOctree.H:86
static volumeType getSide(const vector &outsideNormal, const vector &vec)
Helper function to return the side. Returns outside if.
const List< node > & nodes() const
List of all nodes.
static bool isNode(const labelBits i)
labelBits findNode(const label nodeI, const point &) const
Find deepest node (as parent+octant) containing point. Starts.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
Version of OSstream which prints a prefix on each line.
A class for handling words, derived from string.
Definition: word.H:59
direction bits() const
Definition: labelBits.H:102
const Type & shapes() const
Reference to shape.
const treeBoundBox & bb() const
Top bounding box.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Tree node. Has up pointer and down pointers.
Definition: indexedOctree.H:78
bool ln(const fileName &src, const fileName &dst)
Create a softlink. dst should not exist. Returns true if successful.
Definition: POSIX.C:908
indexedOctree(const Type &shapes)
Construct null.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
friend bool operator!=(const node &a, const node &b)
void print(prefixOSstream &, const bool printContents, const label) const
Print tree. Either print all indices (printContent = true) or.
static bool isContent(const labelBits i)
label findInside(const point &) const
Find shape containing point. Only implemented for certain.
label facePoint(const int facei, const block &block, const label i, const label j)
autoPtr< indexedOctree< Type > > clone() const
Clone.
static label getContent(const labelBits i)
friend bool operator==(const node &a, const node &b)
const labelList & findIndices(const point &) const
Find the shape indices that occupy the result of findNode.
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
static label getNode(const labelBits i)
FixedList< labelBits, 8 > subNodes_
IDs of the 8 nodes on all sides of the mid point.
Definition: indexedOctree.H:89
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
A 29bits label and 3bits direction packed into single label.
Definition: labelBits.H:51
label n
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
friend Istream & operator>>(Istream &is, node &n)
Definition: indexedOctree.H:97
friend Ostream & operator<<(Ostream &os, const node &n)
Definition: indexedOctree.H:91
bool write(Ostream &os) const
label val() const
Definition: labelBits.H:97
Namespace for OpenFOAM.