treeDataEdge.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 | Copyright (C) 2011-2013 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 \*---------------------------------------------------------------------------*/
25 
26 #include "treeDataEdge.H"
27 #include "indexedOctree.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 defineTypeNameAndDebug(treeDataEdge, 0);
34 }
35 
36 
37 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
38 
39 Foam::treeBoundBox Foam::treeDataEdge::calcBb(const label edgeI) const
40 {
41  const edge& e = edges_[edgeI];
42  const point& p0 = points_[e[0]];
43  const point& p1 = points_[e[1]];
44 
45  return treeBoundBox(min(p0, p1), max(p0, p1));
46 }
47 
48 
49 void Foam::treeDataEdge::update()
50 {
51  if (cacheBb_)
52  {
53  bbs_.setSize(edgeLabels_.size());
54 
55  forAll(edgeLabels_, i)
56  {
57  bbs_[i] = calcBb(edgeLabels_[i]);
58  }
59  }
60 }
61 
62 
63 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
64 
66 (
67  const bool cacheBb,
68  const edgeList& edges,
69  const pointField& points,
70  const labelUList& edgeLabels
71 )
72 :
73  edges_(edges),
74  points_(points),
75  edgeLabels_(edgeLabels),
76  cacheBb_(cacheBb)
77 {
78  update();
79 }
80 
81 
83 (
84  const bool cacheBb,
85  const edgeList& edges,
86  const pointField& points,
87  const Xfer<labelList>& edgeLabels
88 )
89 :
90  edges_(edges),
91  points_(points),
92  edgeLabels_(edgeLabels),
93  cacheBb_(cacheBb)
94 {
95  update();
96 }
97 
98 
100 (
101  const indexedOctree<treeDataEdge>& tree
102 )
103 :
104  tree_(tree)
105 {}
106 
107 
109 (
110  const indexedOctree<treeDataEdge>& tree
111 )
112 {}
113 
114 
115 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
116 
118 {
119  pointField eMids(edgeLabels_.size());
120 
121  forAll(edgeLabels_, i)
122  {
123  const edge& e = edges_[edgeLabels_[i]];
124 
125  eMids[i] = e.centre(points_);
126  }
127  return eMids;
128 }
129 
130 
131 //- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
132 // Only makes sense for closed surfaces.
134 (
135  const indexedOctree<treeDataEdge>& oc,
136  const point& sample
137 ) const
138 {
139  return volumeType::UNKNOWN;
140 }
141 
142 
143 // Check if any point on shape is inside cubeBb.
145 (
146  const label index,
147  const treeBoundBox& cubeBb
148 ) const
149 {
150  const edge& e = edges_[edgeLabels_[index]];
151 
152  const point& start = points_[e.start()];
153  const point& end = points_[e.end()];
154 
155  point intersect;
156 
157  return cubeBb.intersects(start, end, intersect);
158 }
159 
160 
161 // Check if any point on shape is inside sphere.
163 (
164  const label index,
165  const point& centre,
166  const scalar radiusSqr
167 ) const
168 {
169  const edge& e = edges_[edgeLabels_[index]];
170 
171  const pointHit nearHit = e.line(points_).nearestDist(centre);
172 
173  const scalar distSqr = sqr(nearHit.distance());
174 
175  if (distSqr <= radiusSqr)
176  {
177  return true;
178  }
179 
180  return false;
181 }
182 
183 
184 void Foam::treeDataEdge::findNearestOp::operator()
185 (
186  const labelUList& indices,
187  const point& sample,
188 
189  scalar& nearestDistSqr,
190  label& minIndex,
191  point& nearestPoint
192 ) const
193 {
194  const treeDataEdge& shape = tree_.shapes();
195 
196  forAll(indices, i)
197  {
198  const label index = indices[i];
199 
200  const edge& e = shape.edges()[shape.edgeLabels()[index]];
201 
202  pointHit nearHit = e.line(shape.points()).nearestDist(sample);
203 
204  scalar distSqr = sqr(nearHit.distance());
205 
206  if (distSqr < nearestDistSqr)
207  {
208  nearestDistSqr = distSqr;
209  minIndex = index;
210  nearestPoint = nearHit.rawPoint();
211  }
212  }
213 }
214 
215 
216 void Foam::treeDataEdge::findNearestOp::operator()
217 (
218  const labelUList& indices,
219  const linePointRef& ln,
220 
221  treeBoundBox& tightest,
222  label& minIndex,
223  point& linePoint,
224  point& nearestPoint
225 ) const
226 {
227  const treeDataEdge& shape = tree_.shapes();
228 
229  // Best so far
230  scalar nearestDistSqr = magSqr(linePoint - nearestPoint);
231 
232  forAll(indices, i)
233  {
234  const label index = indices[i];
235 
236  const edge& e = shape.edges()[shape.edgeLabels()[index]];
237 
238  // Note: could do bb test ? Worthwhile?
239 
240  // Nearest point on line
241  point ePoint, lnPt;
242  scalar dist = e.line(shape.points()).nearestDist(ln, ePoint, lnPt);
243  scalar distSqr = sqr(dist);
244 
245  if (distSqr < nearestDistSqr)
246  {
247  nearestDistSqr = distSqr;
248  minIndex = index;
249  linePoint = lnPt;
250  nearestPoint = ePoint;
251 
252  {
253  point& minPt = tightest.min();
254  minPt = min(ln.start(), ln.end());
255  minPt.x() -= dist;
256  minPt.y() -= dist;
257  minPt.z() -= dist;
258  }
259  {
260  point& maxPt = tightest.max();
261  maxPt = max(ln.start(), ln.end());
262  maxPt.x() += dist;
263  maxPt.y() += dist;
264  maxPt.z() += dist;
265  }
266  }
267  }
268 }
269 
270 
271 bool Foam::treeDataEdge::findIntersectOp::operator()
272 (
273  const label index,
274  const point& start,
275  const point& end,
276  point& result
277 ) const
278 {
280  (
281  "treeDataEdge::intersects(const label, const point&,"
282  "const point&, point&)"
283  );
284  return false;
285 }
286 
287 
288 // ************************************************************************* //
PointHit< Point > nearestDist(const Point &p) const
Return nearest distance to line from a given point.
Definition: lineI.H:95
const Point & rawPoint() const
Return point with no checking.
Definition: PointHit.H:158
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
bool overlaps(const label index, const treeBoundBox &sampleBb) const
Does (bb of) shape at index overlap bb.
Definition: treeDataEdge.C:145
volumeType getVolumeType(const indexedOctree< treeDataEdge > &, const point &) const
Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
Definition: treeDataEdge.C:134
vector point
Point is a vector.
Definition: point.H:41
findNearestOp(const indexedOctree< treeDataEdge > &tree)
Definition: treeDataEdge.C:100
dimensioned< scalar > magSqr(const dimensioned< Type > &)
A simple container for copying or transferring objects of type <T>.
Definition: Xfer.H:85
findIntersectOp(const indexedOctree< treeDataEdge > &tree)
Definition: treeDataEdge.C:109
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
const edgeList & edges() const
Definition: treeDataEdge.H:169
static const Vector max
Definition: Vector.H:82
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
Namespace for OpenFOAM.
A line primitive.
Definition: line.H:56
const Cmpt & y() const
Definition: VectorI.H:71
bool ln(const fileName &src, const fileName &dst)
Create a softlink. dst should not exist. Returns true if successful.
Definition: POSIX.C:855
linePointRef line(const pointField &) const
Return edge line.
Definition: edgeI.H:169
const pointField & points() const
Definition: treeDataEdge.H:174
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define forAll(list, i)
Definition: UList.H:421
const Cmpt & x() const
Definition: VectorI.H:65
const Cmpt & z() const
Definition: VectorI.H:77
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
Holds data for octree to work on an edges subset.
Definition: treeDataEdge.H:53
static const Vector min
Definition: Vector.H:83
label end() const
Return end vertex label.
Definition: edgeI.H:92
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
label start() const
Return start vertex label.
Definition: edgeI.H:81
const labelList & edgeLabels() const
Definition: treeDataEdge.H:179
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
#define notImplemented(functionName)
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:356
treeDataEdge(const bool cacheBb, const edgeList &edges, const pointField &points, const labelUList &edgeLabels)
Construct from selected edges. !Holds references to edges and points.
Definition: treeDataEdge.C:66
pointField shapePoints() const
Get representative point cloud for all shapes inside.
Definition: treeDataEdge.C:117
defineTypeNameAndDebug(combustionModel, 0)
point centre(const pointField &) const
Return centre (centroid)
Definition: edgeI.H:151