treeBoundBox.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2015 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::treeBoundBox
26 
27 Description
28  Standard boundBox + extra functionality for use in octree.
29 
30  Numbering of corner points is according to octant numbering.
31 
32  On the back plane (z=0):
33 
34  \verbatim
35  Y
36  ^
37  |
38  +--------+
39  |2 3|
40  | |
41  | |
42  | |
43  |0 1|
44  +--------+->X
45  \endverbatim
46 
47  For the front plane add 4 to the point labels.
48 
49 
50 SourceFiles
51  treeBoundBoxI.H
52  treeBoundBox.C
53  treeBoundBoxTemplates.C
54 
55 \*---------------------------------------------------------------------------*/
56 
57 #ifndef treeBoundBox_H
58 #define treeBoundBox_H
59 
60 #include "boundBox.H"
61 #include "direction.H"
62 #include "pointField.H"
63 #include "faceList.H"
64 
65 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
66 
67 namespace Foam
68 {
69 
70 class Random;
71 
72 /*---------------------------------------------------------------------------*\
73  Class treeBoundBox Declaration
74 \*---------------------------------------------------------------------------*/
75 
76 class treeBoundBox
77 :
78  public boundBox
79 {
80 
81 private:
82 
83  //- To initialise edges.
84  static edgeList calcEdges(const label[12][2]);
85 
86  //- To initialise faceNormals.
87  static FixedList<vector, 6> calcFaceNormals();
88 
89 public:
90 
91  // Static data members
92 
93  //- The great value used for greatBox and invertedBox
94  static const scalar great;
95 
96  //- As per boundBox::greatBox, but with GREAT instead of VGREAT
97  static const treeBoundBox greatBox;
98 
99  //- As per boundBox::invertedBox, but with GREAT instead of VGREAT
100  static const treeBoundBox invertedBox;
101 
102  //- Bits used for octant/point coding.
103  // Every octant/corner point is the combination of three faces.
104  enum octantBit
105  {
106  RIGHTHALF = 0x1 << 0,
107  TOPHALF = 0x1 << 1,
108  FRONTHALF = 0x1 << 2
109  };
110 
111  //- Face codes
112  enum faceId
113  {
114  LEFT = 0,
115  RIGHT = 1,
116  BOTTOM = 2,
117  TOP = 3,
118  BACK = 4,
119  FRONT = 5
120  };
121 
122  //- Bits used for face coding
123  enum faceBit
124  {
125  NOFACE = 0,
126  LEFTBIT = 0x1 << LEFT, //1
127  RIGHTBIT = 0x1 << RIGHT, //2
128  BOTTOMBIT = 0x1 << BOTTOM, //4
129  TOPBIT = 0x1 << TOP, //8
130  BACKBIT = 0x1 << BACK, //16
131  FRONTBIT = 0x1 << FRONT, //32
132  };
133 
134  //- Edges codes.
135  // E01 = edge between 0 and 1.
136  enum edgeId
137  {
138  E01 = 0,
139  E13 = 1,
140  E23 = 2,
141  E02 = 3,
143  E45 = 4,
144  E57 = 5,
145  E67 = 6,
146  E46 = 7,
148  E04 = 8,
149  E15 = 9,
150  E37 = 10,
151  E26 = 11
152  };
153 
154  //- Face to point addressing
155  static const faceList faces;
156 
157  //- Edge to point addressing
158  static const edgeList edges;
159 
160  //- Per face the unit normal
161  static const FixedList<vector, 6> faceNormals;
162 
163 
164  // Constructors
165 
166  //- Construct null setting points to zero
167  inline treeBoundBox();
168 
169  //- Construct from a boundBox
170  explicit inline treeBoundBox(const boundBox& bb);
171 
172  //- Construct from components
173  inline treeBoundBox(const point& min, const point& max);
174 
175  //- Construct as the bounding box of the given pointField.
176  // Local processor domain only (no reduce as in boundBox)
177  explicit treeBoundBox(const UList<point>&);
178 
179  //- Construct as subset of points
180  // Local processor domain only (no reduce as in boundBox)
181  treeBoundBox(const UList<point>&, const labelUList& indices);
182 
183  //- Construct as subset of points
184  // The indices could be from edge/triFace etc.
185  // Local processor domain only (no reduce as in boundBox)
186  template<unsigned Size>
188  (
189  const UList<point>&,
190  const FixedList<label, Size>& indices
191  );
192 
193 
194  //- Construct from Istream
195  inline treeBoundBox(Istream&);
196 
197 
198  // Member functions
199 
200  // Access
201 
202  //- Typical dimension length,height,width
203  inline scalar typDim() const;
204 
205  //- Vertex coordinates. In octant coding.
206  tmp<pointField> points() const;
207 
208 
209  // Check
210 
211  //- Corner point given octant
212  inline point corner(const direction) const;
213 
214  //- Sub box given by octant number. Midpoint calculated.
215  treeBoundBox subBbox(const direction) const;
216 
217  //- Sub box given by octant number. Midpoint provided.
218  treeBoundBox subBbox(const point& mid, const direction) const;
219 
220  //- Returns octant number given point and the calculated midpoint.
221  inline direction subOctant
222  (
223  const point& pt
224  ) const;
225 
226  //- Returns octant number given point and midpoint.
227  static inline direction subOctant
228  (
229  const point& mid,
230  const point& pt
231  );
232 
233  //- Returns octant number given point and the calculated midpoint.
234  // onEdge set if the point is on edge of subOctant
235  inline direction subOctant
236  (
237  const point& pt,
238  bool& onEdge
239  ) const;
240 
241  //- Returns octant number given point and midpoint.
242  // onEdge set if the point is on edge of subOctant
243  static inline direction subOctant
244  (
245  const point& mid,
246  const point& pt,
247  bool& onEdge
248  );
249 
250  //- Returns octant number given intersection and midpoint.
251  // onEdge set if the point is on edge of subOctant
252  // If onEdge, the direction vector determines which octant to use
253  // (acc. to which octant the point would be if it were moved
254  // along dir)
255  static inline direction subOctant
256  (
257  const point& mid,
258  const vector& dir,
259  const point& pt,
260  bool& onEdge
261  );
262 
263  //- Calculates optimal order to look for nearest to point.
264  // First will be the octant containing the point,
265  // second the octant with boundary nearest to the point etc.
266  inline void searchOrder
267  (
268  const point& pt,
269  FixedList<direction, 8>& octantOrder
270  ) const;
271 
272  //- Overlaps other bounding box?
273  using boundBox::overlaps;
274 
275  //- Intersects segment; set point to intersection position and face,
276  // return true if intersection found.
277  // (pt argument used during calculation even if not intersecting).
278  // Calculates intersections from outside supplied vector
279  // (overallStart, overallVec). This is so when
280  // e.g. tracking through lots of consecutive boxes
281  // (typical octree) we're not accumulating truncation errors. Set
282  // to start, (end-start) if not used.
283  bool intersects
284  (
285  const point& overallStart,
286  const vector& overallVec,
287  const point& start,
288  const point& end,
289  point& pt,
290  direction& ptBits
291  ) const;
292 
293  //- Like above but does not return faces point is on
294  bool intersects
295  (
296  const point& start,
297  const point& end,
298  point& pt
299  ) const;
300 
301  //- Contains point or other bounding box?
302  using boundBox::contains;
303 
304  //- Contains point (inside or on edge) and moving in direction
305  // dir would cause it to go inside.
306  bool contains(const vector& dir, const point&) const;
307 
308  //- Code position of point on bounding box faces
309  direction faceBits(const point&) const;
310 
311  //- Position of point relative to bounding box
312  direction posBits(const point&) const;
313 
314  //- Calculate nearest and furthest (to point) vertex coords of
315  // bounding box
316  void calcExtremities
317  (
318  const point& pt,
319  point& nearest,
320  point& furthest
321  ) const;
322 
323  //- Returns distance point to furthest away corner.
324  scalar maxDist(const point&) const;
325 
326  //- Compare distance to point with other bounding box
327  // return:
328  // -1 : all vertices of my bounding box are nearer than any of
329  // other
330  // +1 : all vertices of my bounding box are further away than
331  // any of other
332  // 0 : none of the above.
333  label distanceCmp(const point&, const treeBoundBox& other) const;
334 
335  //- Return slightly wider bounding box
336  // Extends all dimensions with s*span*Random::scalar01()
337  // and guarantees in any direction s*mag(span) minimum width
338  inline treeBoundBox extend(Random&, const scalar s) const;
339 
340  // Friend Operators
341 
342  friend bool operator==(const treeBoundBox&, const treeBoundBox&);
343  friend bool operator!=(const treeBoundBox&, const treeBoundBox&);
344 
345  // IOstream operator
346 
347  friend Istream& operator>>(Istream& is, treeBoundBox&);
348  friend Ostream& operator<<(Ostream& os, const treeBoundBox&);
349 
350 };
351 
352 
353 //- Data associated with treeBoundBox type are contiguous
354 template<>
355 inline bool contiguous<treeBoundBox>() {return contiguous<boundBox>();}
356 
357 
358 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
359 
360 } // End namespace Foam
361 
362 #include "treeBoundBoxI.H"
363 
364 #ifdef NoRepository
365 # include "treeBoundBoxTemplates.C"
366 #endif
367 
368 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
369 
370 #endif
371 
372 // ************************************************************************* //
void calcExtremities(const point &pt, point &nearest, point &furthest) const
Calculate nearest and furthest (to point) vertex coords of.
Definition: treeBoundBox.C:513
bool intersects(const point &overallStart, const vector &overallVec, const point &start, const point &end, point &pt, direction &ptBits) const
Intersects segment; set point to intersection position and face,.
Definition: treeBoundBox.C:260
Simple random number generator.
Definition: Random.H:49
unsigned char direction
Definition: direction.H:43
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
edgeId
Edges codes.
Definition: treeBoundBox.H:135
faceId
Face codes.
Definition: treeBoundBox.H:111
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
point corner(const direction) const
Corner point given octant.
Definition: treeBoundBoxI.H:63
bool contains(const vector &dir, const point &) const
Contains point (inside or on edge) and moving in direction.
Definition: treeBoundBox.C:402
friend bool operator!=(const treeBoundBox &, const treeBoundBox &)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
direction subOctant(const point &pt) const
Returns octant number given point and the calculated midpoint.
Definition: treeBoundBoxI.H:75
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
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:75
point nearest(const point &) const
Return the nearest point on the boundBox to the supplied point.
Definition: boundBox.C:303
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:53
tmp< pointField > points() const
Vertex coordinates. In octant coding.
Definition: treeBoundBox.C:168
Namespace for OpenFOAM.
static const edgeList edges
Edge to point addressing.
Definition: treeBoundBox.H:157
static faceList faces()
Return faces with correct point order.
Definition: boundBox.C:167
static const scalar great
The great value used for greatBox and invertedBox.
Definition: treeBoundBox.H:93
treeBoundBox()
Construct null setting points to zero.
Definition: treeBoundBoxI.H:31
octantBit
Bits used for octant/point coding.
Definition: treeBoundBox.H:103
treeBoundBox extend(Random &, const scalar s) const
Return slightly wider bounding box.
void searchOrder(const point &pt, FixedList< direction, 8 > &octantOrder) const
Calculates optimal order to look for nearest to point.
friend bool operator==(const treeBoundBox &, const treeBoundBox &)
faceBit
Bits used for face coding.
Definition: treeBoundBox.H:122
static const FixedList< vector, 6 > faceNormals
Per face the unit normal.
Definition: treeBoundBox.H:160
direction faceBits(const point &) const
Code position of point on bounding box faces.
Definition: treeBoundBox.C:442
const point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
friend Ostream & operator<<(Ostream &os, const treeBoundBox &)
direction posBits(const point &) const
Position of point relative to bounding box.
Definition: treeBoundBox.C:476
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:55
bool contains(const point &) const
Contains point? (inside or on edge)
Definition: boundBoxI.H:170
scalar maxDist(const point &) const
Returns distance point to furthest away corner.
Definition: treeBoundBox.C:560
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
Direction is an integer type used to represent the Cartesian directions etc. Currently it is a typede...
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const treeBoundBox greatBox
As per boundBox::greatBox, but with GREAT instead of VGREAT.
Definition: treeBoundBox.H:96
treeBoundBox subBbox(const direction) const
Sub box given by octant number. Midpoint calculated.
Definition: treeBoundBox.C:183
bool contiguous< boundBox >()
Data associated with boundBox type are contiguous.
Definition: boundBox.H:248
scalar typDim() const
Typical dimension length,height,width.
Definition: treeBoundBoxI.H:57
static const treeBoundBox invertedBox
As per boundBox::invertedBox, but with GREAT instead of VGREAT.
Definition: treeBoundBox.H:99
label distanceCmp(const point &, const treeBoundBox &other) const
Compare distance to point with other bounding box.
Definition: treeBoundBox.C:573
friend Istream & operator>>(Istream &is, treeBoundBox &)
bool overlaps(const boundBox &) const
Overlaps/touches boundingBox?
Definition: boundBoxI.H:120
A class for managing temporary objects.
Definition: PtrList.H:118
bool contiguous< treeBoundBox >()
Data associated with treeBoundBox type are contiguous.
Definition: treeBoundBox.H:354