6 template <
class Po
int,
class Base>
9 memset(_children, 0,
sizeof(_children));
12 template <
class Po
int,
class Base>
14 _axis(cell._axis), _splitValue(cell._splitValue)
16 if (cell._children[0])
24 if (cell._children[1])
26 _children[1] =
new AAKdCell<Point, Base>(*cell._children[1]);
34 template <
class Po
int,
class Base>
37 for (
unsigned int i = 0; i < 2; ++i)
44 template <
class Po
int,
class Base>
48 return _children[
index];
51 template <
class Po
int,
class Base>
55 return _children[
index];
58 template <
class Po
int,
class Base>
62 _children[
index] = child;
65 template <
class Po
int,
class Base>
72 template <
class Po
int,
class Base>
79 template <
class Po
int,
class Base>
86 template <
class Po
int,
class Base>
95 template <
class Strategies>
105 std::list<CellType*> stack;
106 stack.push_back(root);
111 if (ShouldSubdivide(
c))
118 for (
unsigned int i = 0; i < CellType::NChildren; ++i)
120 stack.push_back((*
c)[i]);
126 template <
class Strategies>
130 std::vector<size_t>* points)
const
137 PointsInSphere(Root(),
center, radius, points);
140 template <
class Strategies>
143 std::vector<size_t>* points)
const
152 PointsInAACube(Root(),
center, cube.
Width() / 2, points);
155 template <
class Strategies>
159 std::vector<size_t>* points)
const
166 PointsInAACube(Root(),
center, width / 2, points);
169 template <
class Strategies>
181 std::list<CellType*> stack;
189 if (!ShouldSubdivide(
c))
191 for (
unsigned int i = 0; i < CellType::NChildren; ++i)
200 for (
unsigned int i = 0; i < CellType::NChildren; ++i)
202 stack.push_back((*
c)[i]);
206 else if (ShouldSubdivide(
c))
209 for (
unsigned int i = 0; i < CellType::NChildren; ++i)
211 stack.push_back((*
c)[i]);
215 StrategyBaseType::RefreshWithNewTreeData(bc);
218 template <
class Strategies>
224 *dist = std::numeric_limits<ScalarType>::infinity();
232 if (*dist < std::numeric_limits<ScalarType>::infinity())
238 template <
class Strategies>
242 std::vector<NN>* neighbors,
245 *dist = std::numeric_limits<ScalarType>::infinity();
254 neighbors->reserve(k + 1);
255 KNearestNeighbors(*Root(), rootBox, p, k, neighbors, &worstIdx, dist);
256 if (*dist < std::numeric_limits<ScalarType>::infinity())
262 template <
class Strategies>
268 ScalarType* dist2)
const
273 for (
size_t i = 0; i < cell.Size(); ++i)
275 size_t idx = Translate(cell.Points() + i);
276 PointType diff = PointTranslated(idx) - p;
277 ScalarType d2 = diff * diff;
287 unsigned int nearerChild, fartherChild;
288 if (p[cell.Axis()] <= cell.Split())
298 AABox<PointType> subBoxes[2];
299 box.
Split(cell.Axis(), cell.Split(), &subBoxes[0], &subBoxes[1]);
300 NearestNeighbor(*cell[nearerChild], subBoxes[nearerChild], p, neighbor, dist2);
303 if (Intersect(subBoxes[fartherChild], p, *dist2))
304 NearestNeighbor(*cell[fartherChild], subBoxes[fartherChild], p, neighbor, dist2);
307 template <
class Strategies>
310 const AABox<PointType>& box,
313 std::vector<NN>* neighbors,
315 ScalarType* dist2)
const
320 for (
size_t i = 0; i < cell.Size(); ++i)
322 size_t idx = Translate(cell.Points() + i);
323 PointType diff = PointTranslated(idx) - p;
324 ScalarType d2 = diff * diff;
325 if (neighbors->size() < k)
327 neighbors->push_back(NN(idx, d2));
330 else if (d2 < *dist2)
332 (*neighbors)[*worstIdx] = NN(idx, d2);
339 unsigned int nearerChild, fartherChild;
340 if (p[cell.Axis()] <= cell.Split())
350 AABox<PointType> subBoxes[2];
351 box.Split(cell.Axis(), cell.Split(), &subBoxes[0], &subBoxes[1]);
354 *cell[nearerChild], subBoxes[nearerChild], p, k, neighbors, worstIdx, dist2);
357 if (neighbors->size() < k || Intersect(subBoxes[fartherChild], p, *dist2))
359 *cell[fartherChild], subBoxes[fartherChild], p, k, neighbors, worstIdx, dist2);
362 template <
class Strategies>
364 AAKdTree<Strategies>::ReadjustData(CellType* cell)
371 pmax = pmin = PointTranslated(Translate(cell->Points()));
372 for (
size_t i = 1; i < cell->Size(); ++i)
374 const PointType& p = PointTranslated(Translate(cell->Points() + i));
375 for (
unsigned int j = 0; j < Dim; ++j)
381 else if (pmax[j] < p[j])
389 unsigned int axis = 0;
390 ScalarType length = pdiff[0];
391 for (
unsigned int j = 1; j < Dim; ++j)
393 if (pdiff[j] > length)
400 ScalarType
split = (pmax[axis] + pmin[axis]) / 2;
402 CellType *
left = (*cell)[0], *
right = (*cell)[1];
408 template <
class Strategies>
413 std::vector<size_t>* points)
const
417 ScalarType radius2 = radius * radius;
418 for (
size_t i = 0; i < cell->Size(); ++i)
420 size_t ind = Translate(cell->Points() + i);
421 const PointType& p = PointTranslated(ind);
423 if ((d * d) <= radius2)
425 points->push_back(ind);
431 ScalarType d =
center[cell->Axis()] - cell->Split();
434 PointsInSphere((*cell)[1],
center, radius, points);
438 PointsInSphere((*cell)[0],
center, radius, points);
443 template <
class Strategies>
448 std::vector<HandleType>* points)
const
452 for (
size_t i = 0; i < cell->Size(); ++i)
454 size_t ind = Translate(cell->Points() + i);
455 const PointType& p = PointTranslated(ind);
457 bool isInside =
true;
458 for (
unsigned int j = 0; j < PointType::Dim; ++j)
468 points->push_back(ind);
474 ScalarType d =
center[cell->Axis()] - cell->Split();
477 PointsInAACube((*cell)[1],
center, radius, points);
481 PointsInAACube((*cell)[0],
center, radius, points);
486 template <
class Strategies>
488 AAKdTree<Strategies>::ShouldSubdivide(
const CellType* cell)
const
490 if (cell->Size() > 10)
497 template <
class Strategies>
499 AAKdTree<Strategies>::Subdivide(CellType* cell)
502 pmax = pmin = PointTranslated(Translate(cell->Points()));
503 for (
size_t i = 1; i < cell->Size(); ++i)
505 const PointType& p = PointTranslated(Translate(cell->Points() + i));
506 for (
unsigned int j = 0; j < Dim; ++j)
512 else if (pmax[j] < p[j])
520 unsigned int axis = 0;
521 ScalarType length = pdiff[0];
522 for (
unsigned int j = 1; j < Dim; ++j)
524 if (pdiff[j] > length)
531 ScalarType
split = (pmax[axis] + pmin[axis]) / 2;
535 right =
new CellType;
538 if (
left->Size() == cell->Size() ||
right->Size() == cell->Size())
546 cell->Child(0,
left);
547 cell->Child(1,
right);
550 template <
class Strategies>
552 AAKdTree<Strategies>::Absolute(ScalarType f)
const