26 #ifndef __TERRALIB_SAM_KDTREE_INTERNAL_INDEX_H
27 #define __TERRALIB_SAM_KDTREE_INTERNAL_INDEX_H
30 #include "../../geometry/Envelope.h"
67 template<
class KdTreeNode>
class Index
72 typedef typename KdTreeNode::kdKey
kdKey;
75 typedef typename KdTreeNode::kdData
kdData;
109 const std::size_t&
size()
const
138 const std::size_t last = dataSet.size() - 1;
155 erase(node->getLeft());
158 erase(node->getRight());
173 node->getData().push_back(data);
177 inline void search(
const te::gm::Envelope& e, KdTreeNode* node,
const char& level, std::vector<KdTreeNode*>& report)
const;
180 KdTreeNode*
buildOptimized(std::vector<std::pair<kdKey, kdData> >& dataSet,
const std::size_t& first,
const std::size_t& last)
182 const std::size_t kth = (last - first + 1) / 2;
184 KdTreeNode* newNode =
new KdTreeNode(dataSet[first + kth].first);
186 newNode->setData(dataSet[first + kth].second);
190 if(first + kth > first)
191 newNode->setLeft(
buildOptimized(dataSet, first, first + kth - 1));
193 if(first + kth < last)
194 newNode->setRight(
buildOptimized(dataSet, first + kth + 1, last));
222 template<
class KdTreeNode>
227 m_root =
new KdTreeNode(key);
236 KdTreeNode* x = m_root;
245 if(key.getX() > x->getKey().getX())
250 else if(key.getX() < x->getKey().getX())
255 else if(key.getY() == x->getKey().getY())
270 if(key.getY() > x->getKey().getY())
275 else if(key.getY() < x->getKey().getY())
280 else if(key.getX() == x->getKey().getX())
295 KdTreeNode* newNode =
new KdTreeNode(key);
302 y->setRight(newNode);
308 template<
class KdTreeNode>
311 if((node->getKey().getX() >= e.
m_llx) && (node->getKey().getX() <= e.
m_urx) &&
312 (node->getKey().getY() >= e.
m_lly) && (node->getKey().getY() <= e.
m_ury))
314 report.push_back(node);
320 if(node->getKey().getX() >= e.
m_llx)
321 search(e, node->getLeft(),
'y', report);
324 if(node->getKey().getX() <= e.
m_urx)
325 search(e, node->getRight(),
'y', report);
330 if(node->getKey().getY() >= e.
m_lly)
331 search(e, node->getLeft(),
'x', report);
334 if(node->getKey().getY() <= e.
m_ury)
335 search(e, node->getRight(),
'x', report);
343 if (data.getX() > e.
m_urx)
347 if (e.
m_llx > data.getX())
351 if (data.getY() > e.
m_ury)
355 if (e.
m_lly > data.getY())
379 typedef typename KdTreeNode::kdKey
kdKey;
382 typedef typename KdTreeNode::kdData
kdData;
417 const std::size_t&
size()
const
453 void build(std::vector<std::pair<kdKey, kdDataItem> >& dataSet)
469 for(std::size_t i = 0; i < k; ++i)
470 sqrDists.push_back(std::numeric_limits<double>::max());
472 te::gm::Envelope e(-std::numeric_limits<double>::max(), -std::numeric_limits<double>::max(),
473 +std::numeric_limits<double>::max(), +std::numeric_limits<double>::max());
495 erase(node->getLeft());
498 erase(node->getRight());
504 inline KdTreeNode*
build(std::vector<std::pair<kdKey, kdDataItem> >& dataSet,
double averageValue,
const te::gm::Envelope& mbr);
513 inline void update(KdTreeNode* node,
const kdKey& key, std::vector<kdDataItem>& report, std::vector<double>& sqrDists,
te::gm::Envelope& e)
const;
516 double average(std::vector<std::pair<kdKey, kdDataItem> >& dataSet,
const char& discriminator)
const
518 const std::size_t
size = dataSet.size();
520 double medianValue = 0.0;
522 if(discriminator ==
'x')
524 for(
unsigned int i = 0; i <
size; ++i)
525 medianValue += dataSet[i].first.getX();
527 return medianValue /
size;
531 for(
unsigned int i = 0; i <
size; ++i)
532 medianValue += dataSet[i].first.getY();
534 return medianValue /
size;
562 template<
class KdTreeNode>
565 std::vector<KdTreeNode*> reportNodes;
567 search(e, reportNodes);
569 std::size_t nNodes = reportNodes.size();
571 for(std::size_t i = 0; i < nNodes; ++i)
573 std::size_t nElements = reportNodes[i]->getData().size();
575 for(std::size_t j = 0; j < nElements; ++j)
577 if(
Intersects((reportNodes[i])->getData()[j], e))
578 report.push_back((reportNodes[i])->getData()[j]);
583 template<
class KdTreeNode>
588 if(dataSet.size() <= m_bucketSize)
590 KdTreeNode* node =
new KdTreeNode(averageValue);
592 node->setDiscriminator(
'l');
594 std::size_t size = dataSet.size();
596 for(std::size_t i = 0; i < size; ++i)
597 node->getData().push_back(dataSet[i].second);
605 char discriminator =
'x';
607 std::vector<std::pair<kdKey, kdDataItem> > leftDataSet;
608 std::vector<std::pair<kdKey, kdDataItem> > rightDataSet;
614 averageValue = average(dataSet,
'x');
617 newMbr1.
m_urx = averageValue;
618 newMbr2.
m_llx = averageValue;
620 std::size_t size = dataSet.size();
622 for(std::size_t i = 0; i < size; ++ i)
624 if(dataSet[i].first.getX() <= averageValue)
625 leftDataSet.push_back(dataSet[i]);
627 rightDataSet.push_back(dataSet[i]);
635 averageValue = average(dataSet,
'y');
638 newMbr1.
m_ury = averageValue;
639 newMbr2.
m_lly = averageValue;
641 std::size_t size = dataSet.size();
643 for(std::size_t i = 0; i < size; ++ i)
645 if(dataSet[i].first.getY() <= averageValue)
646 leftDataSet.push_back(dataSet[i]);
648 rightDataSet.push_back(dataSet[i]);
654 KdTreeNode* node =
new KdTreeNode(averageValue);
656 if(rightDataSet.size() == 0)
658 node->setDiscriminator(
'l');
660 std::size_t size = leftDataSet.size();
662 for(std::size_t i = 0; i < size; ++i)
663 node->getData().push_back(leftDataSet[i].second);
665 else if(leftDataSet.size() == 0)
667 node->setDiscriminator(
'l');
669 std::size_t size = rightDataSet.size();
671 for(std::size_t i = 0; i < size; ++i)
672 node->getData().push_back(rightDataSet[i].second);
676 node->setDiscriminator(discriminator);
677 node->setLeft(build(leftDataSet, averageValue, newMbr1));
678 node->setRight(build(rightDataSet, averageValue, newMbr2));
684 template<
class KdTreeNode>
687 if(node->getDiscriminator() ==
'l')
689 update(node, key, report, sqrDists, e);
691 else if(node->getDiscriminator() ==
'x')
693 if(key.getX() <= node->getKey())
695 nearestNeighborSearch(node->getLeft(), key, report, sqrDists, e);
697 if((e.
m_llx < node->getKey()) && (node->getKey() < e.
m_urx))
698 nearestNeighborSearch(node->getRight(), key, report, sqrDists, e);
702 nearestNeighborSearch(node->getRight(), key, report, sqrDists, e);
704 if((e.
m_llx < node->getKey()) &&(node->getKey() < e.
m_urx))
705 nearestNeighborSearch(node->getLeft(), key, report, sqrDists, e);
708 else if(node->getDiscriminator() ==
'y')
710 if(key.getY() <= node->getKey())
712 nearestNeighborSearch(node->getLeft(), key, report, sqrDists, e);
714 if((e.
m_lly < node->getKey()) &&(node->getKey() < e.
m_ury))
715 nearestNeighborSearch(node->getRight(), key, report, sqrDists, e);
719 nearestNeighborSearch(node->getRight(), key, report, sqrDists, e);
721 if((e.
m_lly < node->getKey()) &&(node->getKey() < e.
m_ury))
722 nearestNeighborSearch(node->getLeft(), key, report, sqrDists, e);
727 template<
class KdTreeNode>
730 if(node->getDiscriminator() ==
'x')
733 if(e.
m_llx <= node->getKey())
734 search(e, node->getLeft(), report);
737 if(e.
m_urx >= node->getKey())
738 search(e, node->getRight(), report);
740 else if(node->getDiscriminator() ==
'y')
743 if(e.
m_lly <= node->getKey())
744 search(e, node->getLeft(), report);
747 if(e.
m_ury >= node->getKey())
748 search(e, node->getRight(), report);
751 report.push_back(node);
754 template<
class KdTreeNode>
757 const std::size_t size = node->getData().size();
759 const std::size_t nNeighbors = report.size();
762 for(std::size_t i = 0; i < size; ++i)
764 double dx = (key.getX() - node->getData()[i].getX());
765 double dy = (key.getY() - node->getData()[i].getY());
767 double dkp = (dx * dx) + (dy * dy);
770 if(dkp < sqrDists[nNeighbors - 1])
774 for(std::size_t j = 0; j < nNeighbors; ++j)
776 if(dkp < sqrDists[j])
779 for(std::size_t k = nNeighbors - 1; k > j; --k)
781 report[k] = report[k - 1];
782 sqrDists[k] = sqrDists[k - 1];
786 report[j] = node->getData()[i];
796 double maxDist = sqrDists[nNeighbors - 1];
798 if(maxDist != std::numeric_limits<double>::max())
799 maxDist = sqrt(maxDist);
801 e.
m_llx = key.getX() - maxDist;
802 e.
m_lly = key.getY() - maxDist;
803 e.
m_urx = key.getX() + maxDist;
804 e.
m_ury = key.getY() + maxDist;
811 #endif // __TERRALIB_SAM_KDTREE_INTERNAL_INDEX_H