BRL-CAD
BRNode.cpp
Go to the documentation of this file.
1 /* B R N O D E . C P P
2  * BRL-CAD
3  *
4  * Copyright (c) 2014 United States Government as represented by
5  * the U.S. Army Research Laboratory.
6  *
7  * This library is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public License
9  * version 2.1 as published by the Free Software Foundation.
10  *
11  * This library is distributed in the hope that it will be useful, but
12  * WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14  * Lesser General Public License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with this file; see the file named COPYING for more
18  * information.
19  */
20 #include "common.h"
21 #include <algorithm> // for std::max
22 #define NOMINMAX
23 extern "C" {
24 #include "bu/log.h"
25 }
26 #include "brep.h"
27 
28 namespace brlcad {
29 BRNode::~BRNode()
30 {
31  /* delete the children */
32  for (size_t i = 0; i < m_children.size(); i++) {
33  delete m_children[i];
34  }
35 }
36 
37 
38 int
39 BRNode::depth()
40 {
41  int d = 0;
42  for (size_t i = 0; i < m_children.size(); i++) {
43  d = 1 + std::max(d, m_children[i]->depth());
44  }
45  return d;
46 }
47 
48 void
49 BRNode::getLeaves(std::list<BRNode *> &out_leaves)
50 {
51  if (m_children.size() > 0) {
52  for (size_t i = 0; i < m_children.size(); i++) {
53  m_children[i]->getLeaves(out_leaves);
54  }
55  } else {
56  out_leaves.push_back(this);
57  }
58 }
59 
60 BRNode *
61 BRNode::closer(const ON_3dPoint &pt, BRNode *left, BRNode *right)
62 {
63  double ldist = pt.DistanceTo(left->m_estimate);
64  double rdist = pt.DistanceTo(right->m_estimate);
65  TRACE("\t" << ldist << " < " << rdist);
66  if (ldist < rdist) {
67  return left;
68  } else {
69  return right;
70  }
71 }
72 
73 int
74 BRNode::isTrimmed(const ON_2dPoint &uv, double &trimdist) const
75 {
76  point_t bmin, bmax;
77  BRNode::GetBBox(bmin, bmax);
78  if ((bmin[X] <= uv[X]) && (uv[X] <= bmax[X])) /* if check trim and in BBox */
79  {
80  fastf_t v = getCurveEstimateOfV(uv[X], 0.0000001);
81  trimdist = v - uv[Y];
82  if (uv[Y] <= v) {
83  if (m_XIncreasing) {
84  return 1;
85  } else {
86  return 0;
87  }
88  } else if (uv[Y] > v) {
89  if (!m_XIncreasing) {
90  return 1;
91  } else {
92  return 0;
93  }
94  } else {
95  return 1;
96  }
97  } else {
98  trimdist = -1.0;
99  if (m_trimmed) {
100  return 1;
101  } else {
102  return 0;
103  }
104  }
105 }
106 
107 ON_2dPoint
108 BRNode::getClosestPointEstimate(const ON_3dPoint &pt)
109 {
110  ON_Interval u, v;
111  return getClosestPointEstimate(pt, u, v);
112 }
113 
114 ON_2dPoint
115 BRNode::getClosestPointEstimate(const ON_3dPoint &pt, ON_Interval &u, ON_Interval &v)
116 {
117  if (isLeaf()) {
118  double uvs[5][2] = {{m_u.Min(), m_v.Min()}, /* include the corners for an easy refinement */
119  {m_u.Max(), m_v.Min()},
120  {m_u.Max(), m_v.Max()},
121  {m_u.Min(), m_v.Max()},
122  {m_u.Mid(), m_v.Mid()}
123  }; /* include the estimate */
124  ON_3dPoint corners[5];
125  const ON_Surface *surf = m_face->SurfaceOf();
126 
127  u = m_u;
128  v = m_v;
129 
130  /* ??? should we pass these in from SurfaceTree::curveBBox()
131  * to avoid this recalculation?
132  */
133  if (!surf->EvPoint(uvs[0][0], uvs[0][1], corners[0]) ||
134  !surf->EvPoint(uvs[1][0], uvs[1][1], corners[1]) ||
135  !surf->EvPoint(uvs[2][0], uvs[2][1], corners[2]) ||
136  !surf->EvPoint(uvs[3][0], uvs[3][1], corners[3]))
137  {
138  throw new std::exception(); /* FIXME */
139  }
140  corners[4] = BRNode::m_estimate;
141 
142  /* find the point on the curve closest to pt */
143  size_t mini = 0;
144  double mindist = pt.DistanceTo(corners[mini]);
145  double tmpdist;
146  for (size_t i = 1; i < 5; i++) {
147  tmpdist = pt.DistanceTo(corners[i]);
148  TRACE("\t" << mindist << " < " << tmpdist);
149  if (tmpdist < mindist) {
150  mini = i;
151  mindist = tmpdist;
152  }
153  }
154  TRACE("Closest: " << mindist << "; " << PT2(uvs[mini]));
155  return ON_2dPoint(uvs[mini][0], uvs[mini][1]);
156  } else {
157  if (m_children.size() > 0) {
158  BRNode *closestNode = m_children[0];
159  for (size_t i = 1; i < m_children.size(); i++) {
160  closestNode = closer(pt, closestNode, m_children[i]);
161  }
162  return closestNode->getClosestPointEstimate(pt, u, v);
163  } else {
164  throw new std::exception();
165  }
166  }
167 }
168 
169 fastf_t
170 BRNode::getLinearEstimateOfV(fastf_t u)
171 {
172  fastf_t v = m_start[Y] + m_slope * (u - m_start[X]);
173  return v;
174 }
175 
176 fastf_t
177 BRNode::getCurveEstimateOfV(fastf_t u, fastf_t tol) const
178 {
179  ON_3dVector tangent;
180  point_t A, B;
181  double Ta, Tb;
182 
183  if (m_start[X] < m_end[X]) {
184  VMOVE(A, m_start);
185  VMOVE(B, m_end);
186  Ta = m_t[0];
187  Tb = m_t[1];
188  } else {
189  VMOVE(A, m_end);
190  VMOVE(B, m_start);
191  Ta = m_t[1];
192  Tb = m_t[0];
193  }
194 
195  fastf_t dU = B[X] - A[X];
196  if (NEAR_ZERO(dU, tol)) /* vertical */
197  {
198  return A[Y];
199  }
200 
201  ON_3dVector Tan_start = m_trim->TangentAt(Ta);
202  ON_3dVector Tan_end = m_trim->TangentAt(Tb);
203 
204  fastf_t dT = Tb - Ta;
205  fastf_t guess;
206  ON_3dPoint p;
207 
208  /* Use quick binary subdivision until derivatives at end points in 'u' are within 5 percent */
209  while (!NEAR_ZERO(dU, tol) && !NEAR_ZERO(dT, tol)) {
210  guess = Ta + dT / 2;
211  p = m_trim->PointAt(guess);
212 
213  if (UNLIKELY(NEAR_EQUAL(p[X], u, SMALL_FASTF))) {
214  /* hit 'u' exactly, done deal */
215  return p[Y];
216  }
217 
218  if (p[X] > u) {
219  /* v is behind us, back up the end */
220  Tb = guess;
221  VMOVE(B, p);
222  Tan_end = m_trim->TangentAt(Tb);
223  } else {
224  /* v is in front, move start forward */
225  Ta = guess;
226  VMOVE(A, p);
227  Tan_start = m_trim->TangentAt(Ta);
228  }
229  dT = Tb - Ta;
230  dU = B[X] - A[X];
231  }
232 
233  dU = B[X] - A[X];
234  if (NEAR_ZERO(dU, tol)) /* vertical */
235  {
236  return A[Y];
237  }
238 
239  guess = Ta + (u - A[X]) * dT / dU;
240  p = m_trim->PointAt(guess);
241 
242  int cnt = 0;
243  while ((cnt < 1000) && (!NEAR_EQUAL(p[X], u, tol))) {
244  if (p[X] < u) {
245  Ta = guess;
246  VMOVE(A, p);
247  } else {
248  Tb = guess;
249  VMOVE(B, p);
250  }
251  dU = B[X] - A[X];
252  if (NEAR_ZERO(dU, tol)) { /* vertical */
253  return A[Y];
254  }
255 
256  dT = Tb - Ta;
257  guess = Ta + (u - A[X]) * dT / dU;
258  p = m_trim->PointAt(guess);
259  cnt++;
260  }
261  if (cnt > 999) {
262  bu_log("getCurveEstimateOfV(): estimate of 'v' given a trim curve and "
263  "'u' did not converge within iteration bound(%d).\n", cnt);
264  }
265  return p[Y];
266 }
267 
268 fastf_t
269 BRNode::getCurveEstimateOfU(fastf_t v, fastf_t tol) const
270 {
271  ON_3dVector tangent;
272  point_t A, B;
273  double Ta, Tb;
274 
275  if (m_start[Y] < m_end[Y]) {
276  VMOVE(A, m_start);
277  VMOVE(B, m_end);
278  Ta = m_t[0];
279  Tb = m_t[1];
280  } else {
281  VMOVE(A, m_end);
282  VMOVE(B, m_start);
283  Ta = m_t[1];
284  Tb = m_t[0];
285  }
286 
287  fastf_t dV = B[Y] - A[Y];
288  if (NEAR_ZERO(dV, tol)) /* horizontal */
289  {
290  return A[X];
291  }
292 
293  ON_3dVector Tan_start = m_trim->TangentAt(Ta);
294  ON_3dVector Tan_end = m_trim->TangentAt(Tb);
295 
296  fastf_t dT = Tb - Ta;
297  fastf_t guess;
298  ON_3dPoint p;
299 
300  /* Use quick binary subdivision until derivatives at end points in 'u' are within 5 percent */
301  while (!NEAR_ZERO(dV, tol) && !NEAR_ZERO(dT, tol)) {
302  guess = Ta + dT / 2;
303  p = m_trim->PointAt(guess);
304  if (p[Y] < v) {
305  Ta = guess;
306  VMOVE(A, p);
307  Tan_start = m_trim->TangentAt(Ta);
308  } else {
309  Tb = guess;
310  VMOVE(B, p);
311  Tan_end = m_trim->TangentAt(Tb);
312  }
313  dT = Tb - Ta;
314  dV = B[Y] - A[Y];
315  }
316 
317  dV = B[Y] - A[Y];
318  if (NEAR_ZERO(dV, tol)) /* horizontal */
319  {
320  return A[X];
321  }
322 
323  guess = Ta + (v - A[Y]) * dT / dV;
324  p = m_trim->PointAt(guess);
325 
326  int cnt = 0;
327  while ((cnt < 1000) && (!NEAR_EQUAL(p[Y], v, tol))) {
328  if (p[Y] < v) {
329  Ta = guess;
330  VMOVE(A, p);
331  } else {
332  Tb = guess;
333  VMOVE(B, p);
334  }
335  dV = B[Y] - A[Y];
336  if (NEAR_ZERO(dV, tol)) { /* horizontal */
337  return A[X];
338  }
339  dT = Tb - Ta;
340  guess = Ta + (v - A[Y]) * dT / dV;
341  p = m_trim->PointAt(guess);
342  cnt++;
343  }
344  if (cnt > 999) {
345  bu_log("getCurveEstimateOfV(): estimate of 'u' given a trim curve and "
346  "'v' did not converge within iteration bound(%d).\n", cnt);
347  }
348  return p[X];
349 }
350 }
351 
352 // Local Variables:
353 // tab-width: 8
354 // mode: C++
355 // c-basic-offset: 4
356 // indent-tabs-mode: t
357 // c-file-style: "stroustrup"
358 // End:
359 // ex: shiftwidth=4 tabstop=8
void bu_log(const char *,...) _BU_ATTR_PRINTF12
Definition: log.c:176
#define SMALL_FASTF
Definition: defines.h:342
Header file for the BRL-CAD common definitions.
Definition: color.c:49
#define NEAR_ZERO(val, epsilon)
Definition: color.c:55
#define A
Definition: msr.c:51
double fastf_t
Definition: defines.h:300
Definition: color.c:50
#define UNLIKELY(expression)
Definition: common.h:282