PCL
QuadTree.h
Go to the documentation of this file.
1 // ____ ______ __
2 // / __ \ / ____// /
3 // / /_/ // / / /
4 // / ____// /___ / /___ PixInsight Class Library
5 // /_/ \____//_____/ PCL 2.4.7
6 // ----------------------------------------------------------------------------
7 // pcl/QuadTree.h - Released 2020-12-17T15:46:29Z
8 // ----------------------------------------------------------------------------
9 // This file is part of the PixInsight Class Library (PCL).
10 // PCL is a multiplatform C++ framework for development of PixInsight modules.
11 //
12 // Copyright (c) 2003-2020 Pleiades Astrophoto S.L. All Rights Reserved.
13 //
14 // Redistribution and use in both source and binary forms, with or without
15 // modification, is permitted provided that the following conditions are met:
16 //
17 // 1. All redistributions of source code must retain the above copyright
18 // notice, this list of conditions and the following disclaimer.
19 //
20 // 2. All redistributions in binary form must reproduce the above copyright
21 // notice, this list of conditions and the following disclaimer in the
22 // documentation and/or other materials provided with the distribution.
23 //
24 // 3. Neither the names "PixInsight" and "Pleiades Astrophoto", nor the names
25 // of their contributors, may be used to endorse or promote products derived
26 // from this software without specific prior written permission. For written
27 // permission, please contact info@pixinsight.com.
28 //
29 // 4. All products derived from this software, in any form whatsoever, must
30 // reproduce the following acknowledgment in the end-user documentation
31 // and/or other materials provided with the product:
32 //
33 // "This product is based on software from the PixInsight project, developed
34 // by Pleiades Astrophoto and its contributors (https://pixinsight.com/)."
35 //
36 // Alternatively, if that is where third-party acknowledgments normally
37 // appear, this acknowledgment must be reproduced in the product itself.
38 //
39 // THIS SOFTWARE IS PROVIDED BY PLEIADES ASTROPHOTO AND ITS CONTRIBUTORS
40 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
41 // TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
42 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL PLEIADES ASTROPHOTO OR ITS
43 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
44 // EXEMPLARY OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, BUSINESS
45 // INTERRUPTION; PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; AND LOSS OF USE,
46 // DATA OR PROFITS) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
47 // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
48 // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
49 // POSSIBILITY OF SUCH DAMAGE.
50 // ----------------------------------------------------------------------------
51 
52 #ifndef __PCL_QuadTree_h
53 #define __PCL_QuadTree_h
54 
56 
57 #include <pcl/Defs.h>
58 
59 #include <pcl/Array.h>
60 #include <pcl/Rectangle.h>
61 #include <pcl/Vector.h>
62 
63 namespace pcl
64 {
65 
66 // ----------------------------------------------------------------------------
67 
110 template <class T>
111 class QuadTree
112 {
113 public:
114 
118  typedef T point;
119 
123  typedef typename point::component component;
124 
129 
133  typedef DRect rectangle;
134 
138  typedef rectangle::component coordinate;
139 
140  // -------------------------------------------------------------------------
141 
146  struct Node
147  {
148  rectangle rect = 0.0;
149  Node* nw = nullptr;
150  Node* ne = nullptr;
151  Node* sw = nullptr;
152  Node* se = nullptr;
153 
154 #ifdef __PCL_QUADTREE_STRUCTURAL_NODE_HAS_DATA
155  void* data = nullptr;
156 #endif
157 
161  Node() = default;
162 
166  Node( const rectangle& r )
167  : rect( r )
168  {
169  }
170 
180  bool IsLeaf() const
181  {
182  return nw == nullptr && ne == nullptr && sw == nullptr && se == nullptr;
183  }
184 
189  bool Intersects( const rectangle& r ) const
190  {
191  return rect.x1 >= r.x0 && rect.x0 <= r.x1 &&
192  rect.y1 >= r.y0 && rect.y0 <= r.y1;
193  }
194 
199  bool Includes( coordinate x, coordinate y ) const
200  {
201  return x >= rect.x0 && x <= rect.x1 &&
202  y >= rect.y0 && y <= rect.y1;
203  }
204 
209  bool Includes( const rectangle::point& p ) const
210  {
211  return Includes( p.x, p.y );
212  }
213 
217  rectangle NWRect() const
218  {
219  return rectangle( rect.TopLeft(), rect.Center() );
220  }
221 
225  rectangle NERect() const
226  {
227  return rectangle( rect.CenterTop(), rect.CenterRight() );
228  }
229 
233  rectangle SWRect() const
234  {
235  return rectangle( rect.CenterLeft(), rect.CenterBottom() );
236  }
237 
242  rectangle SERect() const
243  {
244  return rectangle( rect.Center(), rect.BottomRight() );
245  }
246  };
247 
248  // -------------------------------------------------------------------------
249 
254  struct LeafNode : public Node
255  {
263  point_list points;
264 
265 #ifndef __PCL_QUADTREE_STRUCTURAL_NODE_HAS_DATA
266 
275  void* data = nullptr;
276 
277 #endif
278 
283  LeafNode( const rectangle& r, const point_list& p )
284  : Node( r )
285  , points( p )
286  {
287  PCL_CHECK( !points.IsEmpty() )
288  }
289 
294  LeafNode( const rectangle& r, const point& p )
295  : Node( r )
296  {
297  points << p;
298  }
299 
304  int Length() const
305  {
306  return int( points.Length() );
307  }
308  };
309 
310  // -------------------------------------------------------------------------
311 
315  QuadTree() = default;
316 
329  QuadTree( const point_list& points, int bucketCapacity = 40 )
330  {
331  Build( points, bucketCapacity );
332  }
333 
349  QuadTree( const rectangle& rect, const point_list& points, int bucketCapacity = 40 )
350  {
351  Build( rect, points, bucketCapacity );
352  }
353 
359  QuadTree( const QuadTree& ) = delete;
360 
366  QuadTree& operator =( const QuadTree& ) = delete;
367 
372  : m_root( x.m_root )
373  , m_bucketCapacity( x.m_bucketCapacity )
374  , m_length( x.m_length )
375  {
376  x.m_root = nullptr;
377  x.m_length = 0;
378  }
379 
384  {
385  if ( &x != this )
386  {
387  DestroyTree( m_root );
388  m_root = x.m_root;
389  m_bucketCapacity = x.m_bucketCapacity;
390  m_length = x.m_length;
391  x.m_root = nullptr;
392  x.m_length = 0;
393  }
394  return *this;
395  }
396 
401  {
402  Clear();
403  }
404 
408  void Clear()
409  {
410  DestroyTree( m_root );
411  m_root = nullptr;
412  m_length = 0;
413  }
414 
430  void Build( const point_list& points, int bucketCapacity = 40 )
431  {
432  Clear();
433  m_bucketCapacity = Max( 1, bucketCapacity );
434 
435  if ( !points.IsEmpty() )
436  {
437  component x = points[0][0];
438  component y = points[0][1];
439  rectangle rect( x, y, x, y );
440  for ( const point& p : points )
441  {
442  x = p[0];
443  y = p[1];
444  if ( x < rect.x0 )
445  rect.x0 = x;
446  if ( y < rect.y0 )
447  rect.y0 = y;
448  if ( x > rect.x1 )
449  rect.x1 = x;
450  if ( y > rect.y1 )
451  rect.y1 = y;
452  }
453  m_root = BuildTree( rect, points );
454  }
455  }
456 
477  void Build( const rectangle& rect, const point_list& points, int bucketCapacity = 40 )
478  {
479  Clear();
480  m_bucketCapacity = Max( 1, bucketCapacity );
481  if ( !points.IsEmpty() )
482  m_root = BuildTree( rect.Ordered(), points );
483  }
484 
493  point_list Search( const rectangle& rect ) const
494  {
495  point_list found;
496  SearchTree( found, rect.Ordered(), m_root );
497  return found;
498  }
499 
516  template <class F>
517  void Search( const rectangle& rect, F callback, void* data ) const
518  {
519  SearchTree( rect.Ordered(), callback, data, m_root );
520  }
521 
525  void Insert( const point& pt )
526  {
527  if ( m_root != nullptr )
528  InsertTree( pt, m_root );
529  else
530  {
531  component x = pt[0];
532  component y = pt[1];
533  m_root = new LeafNode( rectangle( x, y, x, y ), pt );
534  }
535  }
536 
540  void Delete( const point& pt )
541  {
542  DeleteTree( pt, m_root );
543  }
544 
549  void Delete( const rectangle& rect )
550  {
551  DeleteTree( rect.Ordered(), m_root );
552  }
553 
559  int BucketCapacity() const
560  {
561  return m_bucketCapacity;
562  }
563 
568  {
569  return m_length;
570  }
571 
575  bool IsEmpty() const
576  {
577  return m_root == nullptr;
578  }
579 
587  const Node* Root() const
588  {
589  return m_root;
590  }
591 
597  {
598  return m_root;
599  }
600 
606  const LeafNode* LeafNodeAt( rectangle::point p ) const
607  {
608  return SearchLeafNode( p, m_root );
609  }
610 
616  LeafNode* LeafNodeAt( rectangle::point p )
617  {
618  return SearchLeafNode( p, m_root );
619  }
620 
630  const Node* NodeAt( rectangle::point p ) const
631  {
632  return SearchNode( p, m_root );
633  }
634 
644  Node* NodeAt( rectangle::point p )
645  {
646  return SearchNode( p, m_root );
647  }
648 
659  Node* SplitAt( rectangle::point p )
660  {
661  Node* node = SearchDeepestStructuralNodeAt( p, m_root );
662  if ( node != nullptr )
663  {
664  Node** leaf = nullptr;
665  if ( node->nw != nullptr && node->nw->Includes( p ) )
666  leaf = &node->nw;
667  else if ( node->ne != nullptr && node->ne->Includes( p ) )
668  leaf = &node->ne;
669  else if ( node->sw != nullptr && node->sw->Includes( p ) )
670  leaf = &node->sw;
671  else if ( node->se != nullptr && node->se->Includes( p ) )
672  leaf = &node->se;
673 
674  if ( leaf != nullptr ) // cannot be false!
675  {
676  Node* newNode = SplitLeafNode( *leaf );
677  if ( newNode != nullptr ) // should be true
678  {
679  delete static_cast<LeafNode*>( *leaf );
680  *leaf = newNode;
681  }
682  return newNode;
683  }
684  }
685 
686  return nullptr;
687  }
688 
708  template <class F>
709  static void Traverse( const Node* node, F f )
710  {
711  if ( node != nullptr )
712  if ( node->IsLeaf() )
713  f( node->rect,
714  static_cast<const LeafNode*>( node )->points,
715  static_cast<const LeafNode*>( node )->data );
716  else
717  {
718  Traverse( node->nw, f );
719  Traverse( node->ne, f );
720  Traverse( node->sw, f );
721  Traverse( node->se, f );
722  }
723  }
724 
744  template <class F>
745  static void Traverse( Node* node, F f )
746  {
747  if ( node != nullptr )
748  if ( node->IsLeaf() )
749  f( node->rect, static_cast<LeafNode*>( node )->points,
750  static_cast<LeafNode*>( node )->data );
751  else
752  {
753  Traverse( node->nw, f );
754  Traverse( node->ne, f );
755  Traverse( node->sw, f );
756  Traverse( node->se, f );
757  }
758  }
759 
769  template <class F>
770  void Traverse( F f ) const
771  {
772  Traverse( m_root, f );
773  }
774 
784  template <class F>
785  void Traverse( F f )
786  {
787  Traverse( m_root, f );
788  }
789 
794  static size_type NumberOfNodes( const Node* node )
795  {
796  size_type N = 0;
797  GetNumberOfNodes( N, node );
798  return N;
799  }
800 
806  {
807  return NumberOfNodes( m_root );
808  }
809 
814  static size_type NumberOfLeafNodes( const Node* node )
815  {
816  size_type N = 0;
817  GetNumberOfLeafNodes( N, node );
818  return N;
819  }
820 
825  {
826  return NumberOfLeafNodes( m_root );
827  }
828 
832  static int Height( const Node* node )
833  {
834  return TreeHeight( node, 0 );
835  }
836 
840  int Height() const
841  {
842  return Height( m_root );
843  }
844 
848  friend void Swap( QuadTree& x1, QuadTree& x2 )
849  {
850  pcl::Swap( x1.m_root, x2.m_root );
851  pcl::Swap( x1.m_bucketCapacity, x2.m_bucketCapacity );
852  pcl::Swap( x1.m_length, x2.m_length );
853  }
854 
855 private:
856 
857  Node* m_root = nullptr;
858  int m_bucketCapacity = 0;
859  size_type m_length = 0;
860 
861  Node* NewLeafNode( const rectangle& rect, const point_list& points )
862  {
863  m_length += points.Length();
864  return new LeafNode( rect, points );
865  }
866 
867  LeafNode* SearchLeafNode( const rectangle::point& p, const Node* node ) const
868  {
869  if ( node != nullptr )
870  if ( node->Includes( p ) )
871  {
872  if ( node->IsLeaf() )
873  return const_cast<LeafNode*>( static_cast<const LeafNode*>( node ) );
874 
875  LeafNode* child = SearchLeafNode( p, node->nw );
876  if ( child == nullptr )
877  {
878  child = SearchLeafNode( p, node->ne );
879  if ( child == nullptr )
880  {
881  child = SearchLeafNode( p, node->sw );
882  if ( child == nullptr )
883  child = SearchLeafNode( p, node->se );
884  }
885  }
886  return child;
887  }
888 
889  return nullptr;
890  }
891 
892  Node* SearchNode( const rectangle::point& p, const Node* node ) const
893  {
894  if ( node != nullptr )
895  if ( node->Includes( p ) )
896  {
897  if ( node->IsLeaf() )
898  return const_cast<Node*>( node );
899 
900  Node* child = SearchNode( p, node->nw );
901  if ( child == nullptr )
902  {
903  child = SearchNode( p, node->ne );
904  if ( child == nullptr )
905  {
906  child = SearchNode( p, node->sw );
907  if ( child == nullptr )
908  {
909  child = SearchNode( p, node->se );
910  if ( child == nullptr )
911  return const_cast<Node*>( node );
912  }
913  }
914  }
915  return child;
916  }
917 
918  return nullptr;
919  }
920 
921  Node* SearchDeepestStructuralNodeAt( const rectangle::point& p, const Node* node ) const
922  {
923  if ( node != nullptr )
924  if ( !node->IsLeaf() )
925  if ( node->Includes( p ) )
926  {
927  Node* child = SearchDeepestStructuralNodeAt( p, node->nw );
928  if ( child == nullptr )
929  {
930  child = SearchDeepestStructuralNodeAt( p, node->ne );
931  if ( child == nullptr )
932  {
933  child = SearchDeepestStructuralNodeAt( p, node->sw );
934  if ( child == nullptr )
935  {
936  child = SearchDeepestStructuralNodeAt( p, node->se );
937  if ( child == nullptr )
938  return const_cast<Node*>( node );
939  }
940  }
941  }
942  return child;
943  }
944 
945  return nullptr;
946  }
947 
948  Node* BuildTree( const rectangle& rect, const point_list& points )
949  {
950  if ( points.IsEmpty() )
951  return nullptr;
952 
953  if ( points.Length() <= size_type( m_bucketCapacity ) )
954  return NewLeafNode( rect, points );
955 
956  double x2 = (rect.x0 + rect.x1)/2;
957  double y2 = (rect.y0 + rect.y1)/2;
958  // Prevent geometrically degenerate subtrees. For safety, we enforce
959  // minimum region dimensions larger than twice the machine epsilon for
960  // the rectangle coordinate type.
961  if ( x2 - rect.x0 < 2*std::numeric_limits<coordinate>::epsilon() ||
962  rect.x1 - x2 < 2*std::numeric_limits<coordinate>::epsilon() ||
963  y2 - rect.y0 < 2*std::numeric_limits<coordinate>::epsilon() ||
964  rect.y1 - y2 < 2*std::numeric_limits<coordinate>::epsilon() )
965  {
966  return NewLeafNode( rect, points );
967  }
968 
969  point_list nw, ne, sw, se;
970  for ( const point& p : points )
971  {
972  component x = p[0];
973  component y = p[1];
974  if ( x <= x2 )
975  {
976  if ( y <= y2 )
977  nw << p;
978  else
979  sw << p;
980  }
981  else
982  {
983  if ( y <= y2 )
984  ne << p;
985  else
986  se << p;
987  }
988  }
989 
990  Node* node = new Node( rect );
991  node->nw = BuildTree( rectangle( rect.x0, rect.y0, x2, y2 ), nw );
992  node->ne = BuildTree( rectangle( x2, rect.y0, rect.x1, y2 ), ne );
993  node->sw = BuildTree( rectangle( rect.x0, y2, x2, rect.y1 ), sw );
994  node->se = BuildTree( rectangle( x2, y2, rect.x1, rect.y1 ), se );
995 
996  // Further degeneracies may result, e.g. if the point class is not
997  // behaving as expected. Do not allow them.
998  if ( node->IsLeaf() )
999  {
1000  delete node;
1001  return NewLeafNode( rect, points );
1002  }
1003 
1004  return node;
1005  }
1006 
1007  Node* SplitLeafNode( const Node* node )
1008  {
1009  if ( node == nullptr || !node->IsLeaf() )
1010  return nullptr;
1011 
1012  double x2 = (node->rect.x0 + node->rect.x1)/2;
1013  double y2 = (node->rect.y0 + node->rect.y1)/2;
1014  // Prevent geometrically degenerate subtrees. For safety, we enforce
1015  // minimum region dimensions larger than twice the machine epsilon for
1016  // the rectangle coordinate type.
1017  if ( x2 - node->rect.x0 < 2*std::numeric_limits<coordinate>::epsilon() ||
1018  node->rect.x1 - x2 < 2*std::numeric_limits<coordinate>::epsilon() ||
1019  y2 - node->rect.y0 < 2*std::numeric_limits<coordinate>::epsilon() ||
1020  node->rect.y1 - y2 < 2*std::numeric_limits<coordinate>::epsilon() )
1021  {
1022  return nullptr;
1023  }
1024 
1025  const LeafNode* leaf = static_cast<const LeafNode*>( node );
1026  point_list nw, ne, sw, se;
1027  for ( const point& p : leaf->points )
1028  {
1029  component x = p[0];
1030  component y = p[1];
1031  if ( x <= x2 )
1032  {
1033  if ( y <= y2 )
1034  nw << p;
1035  else
1036  sw << p;
1037  }
1038  else
1039  {
1040  if ( y <= y2 )
1041  ne << p;
1042  else
1043  se << p;
1044  }
1045  }
1046 
1047  size_type savedLength = m_length;
1048  Node* newNode = new Node( node->rect );
1049  newNode->nw = BuildTree( rectangle( node->rect.x0, node->rect.y0, x2, y2 ), nw );
1050  newNode->ne = BuildTree( rectangle( x2, node->rect.y0, node->rect.x1, y2 ), ne );
1051  newNode->sw = BuildTree( rectangle( node->rect.x0, y2, x2, node->rect.y1 ), sw );
1052  newNode->se = BuildTree( rectangle( x2, y2, node->rect.x1, node->rect.y1 ), se );
1053  m_length = savedLength;
1054 
1055  // Further degeneracies may result, e.g. if the point class is not
1056  // behaving as expected. Do not allow them.
1057  if ( newNode->IsLeaf() )
1058  {
1059  delete newNode;
1060  return nullptr;
1061  }
1062 
1063  return newNode;
1064  }
1065 
1066  void SearchTree( point_list& found, const rectangle& rect, const Node* node ) const
1067  {
1068  if ( node != nullptr )
1069  if ( node->Intersects( rect ) )
1070  if ( node->IsLeaf() )
1071  {
1072  const LeafNode* leaf = static_cast<const LeafNode*>( node );
1073  for ( const point& p : leaf->points )
1074  {
1075  component x = p[0];
1076  if ( x >= rect.x0 )
1077  if ( x <= rect.x1 )
1078  {
1079  component y = p[1];
1080  if ( y >= rect.y0 )
1081  if ( y <= rect.y1 )
1082  found << p;
1083  }
1084  }
1085  }
1086  else
1087  {
1088  SearchTree( found, rect, node->nw );
1089  SearchTree( found, rect, node->ne );
1090  SearchTree( found, rect, node->sw );
1091  SearchTree( found, rect, node->se );
1092  }
1093  }
1094 
1095  template <class F>
1096  void SearchTree( const rectangle& rect, F callback, void* data, const Node* node ) const
1097  {
1098  if ( node != nullptr )
1099  if ( node->Intersects( rect ) )
1100  if ( node->IsLeaf() )
1101  {
1102  const LeafNode* leaf = static_cast<const LeafNode*>( node );
1103  for ( const point& p : leaf->points )
1104  {
1105  component x = p[0];
1106  if ( x >= rect.x0 )
1107  if ( x <= rect.x1 )
1108  {
1109  component y = p[1];
1110  if ( y >= rect.y0 )
1111  if ( y <= rect.y1 )
1112  callback( p, data );
1113  }
1114  }
1115  }
1116  else
1117  {
1118  SearchTree( rect, callback, data, node->nw );
1119  SearchTree( rect, callback, data, node->ne );
1120  SearchTree( rect, callback, data, node->sw );
1121  SearchTree( rect, callback, data, node->se );
1122  }
1123  }
1124 
1125  void InsertTree( const point& pt, Node*& node )
1126  {
1127  PCL_CHECK( node != nullptr )
1128 
1129  component x = pt[0];
1130  component y = pt[1];
1131 
1132  if ( x < node->rect.x0 )
1133  node->rect.x0 = x;
1134  else if ( x > node->rect.x1 )
1135  node->rect.x1 = x;
1136 
1137  if ( y < node->rect.y0 )
1138  node->rect.y0 = y;
1139  else if ( y > node->rect.y1 )
1140  node->rect.y1 = y;
1141 
1142  if ( node->IsLeaf() )
1143  {
1144  LeafNode* leaf = static_cast<LeafNode*>( node );
1145  if ( leaf->Length() < m_bucketCapacity )
1146  leaf->points << pt;
1147  else
1148  {
1149  rectangle rect = leaf->rect;
1150  double x2 = (rect.x0 + rect.x1)/2;
1151  double y2 = (rect.y0 + rect.y1)/2;
1152  // Prevent geometrically degenerate subtrees. For safety, we
1153  // enforce minimum region dimensions larger than twice the
1154  // machine epsilon for the rectangle coordinate type.
1155  if ( x2 - rect.x0 < 2*std::numeric_limits<coordinate>::epsilon() ||
1156  rect.x1 - x2 < 2*std::numeric_limits<coordinate>::epsilon() ||
1157  y2 - rect.y0 < 2*std::numeric_limits<coordinate>::epsilon() ||
1158  rect.y1 - y2 < 2*std::numeric_limits<coordinate>::epsilon() )
1159  {
1160  leaf->points << pt;
1161  }
1162  else
1163  {
1164  point_list nw, ne, sw, se;
1165  for ( const point& p : leaf->points )
1166  {
1167  component x = p[0];
1168  component y = p[1];
1169  if ( x <= x2 )
1170  {
1171  if ( y <= y2 )
1172  nw << p;
1173  else
1174  sw << p;
1175  }
1176  else
1177  {
1178  if ( y <= y2 )
1179  ne << p;
1180  else
1181  se << p;
1182  }
1183  }
1184 
1185  if ( x <= x2 )
1186  {
1187  if ( y <= y2 )
1188  nw << pt;
1189  else
1190  sw << pt;
1191  }
1192  else
1193  {
1194  if ( y <= y2 )
1195  ne << pt;
1196  else
1197  se << pt;
1198  }
1199 
1200  delete leaf;
1201 
1202  node = new Node( rect );
1203 
1204  if ( !nw.IsEmpty() )
1205  node->nw = new LeafNode( rectangle( rect.x0, rect.y0, x2, y2 ), nw );
1206  if ( !ne.IsEmpty() )
1207  node->ne = new LeafNode( rectangle( x2, rect.y0, rect.x1, y2 ), ne );
1208  if ( !sw.IsEmpty() )
1209  node->sw = new LeafNode( rectangle( rect.x0, y2, x2, rect.y1 ), sw );
1210  if ( !se.IsEmpty() )
1211  node->se = new LeafNode( rectangle( x2, y2, rect.x1, rect.y1 ), se );
1212  }
1213  }
1214 
1215  ++m_length;
1216  }
1217  else
1218  {
1219  rectangle rect = node->rect;
1220  double x2 = (rect.x0 + rect.x1)/2;
1221  double y2 = (rect.y0 + rect.y1)/2;
1222  if ( pt[0] <= x2 )
1223  {
1224  if ( pt[1] <= y2 )
1225  {
1226  if ( node->nw == nullptr )
1227  node->nw = new LeafNode( rectangle( rect.x0, rect.y0, x2, y2 ), pt );
1228  else
1229  InsertTree( pt, node->nw );
1230  }
1231  else
1232  {
1233  if ( node->sw == nullptr )
1234  node->sw = new LeafNode( rectangle( rect.x0, y2, x2, rect.y1 ), pt );
1235  else
1236  InsertTree( pt, node->sw );
1237  }
1238  }
1239  else
1240  {
1241  if ( pt[1] <= y2 )
1242  {
1243  if ( node->ne == nullptr )
1244  node->ne = new LeafNode( rectangle( x2, rect.y0, rect.x1, y2 ), pt );
1245  else
1246  InsertTree( pt, node->ne );
1247  }
1248  else
1249  {
1250  if ( node->se == nullptr )
1251  node->se = new LeafNode( rectangle( x2, y2, rect.x1, rect.y1 ), pt );
1252  else
1253  InsertTree( pt, node->se );
1254  }
1255  }
1256  }
1257  }
1258 
1259  void DeleteTree( const rectangle& rect, Node*& node )
1260  {
1261  if ( node != nullptr )
1262  if ( node->Intersects( rect ) )
1263  {
1264  if ( node->IsLeaf() )
1265  {
1266  LeafNode* leaf = static_cast<LeafNode*>( node );
1267  point_list points;
1268  for ( const point& p : leaf->points )
1269  {
1270  component x = p[0];
1271  if ( x >= rect.x0 )
1272  if ( x <= rect.x1 )
1273  {
1274  component y = p[1];
1275  if ( y >= rect.y0 )
1276  if ( y <= rect.y1 )
1277  {
1278  --m_length;
1279  continue;
1280  }
1281  }
1282  points << p;
1283  }
1284 
1285  if ( points.IsEmpty() )
1286  {
1287  delete leaf;
1288  node = nullptr;
1289  }
1290  else
1291  leaf->points = points;
1292  }
1293  else
1294  {
1295  DeleteTree( rect, node->nw );
1296  DeleteTree( rect, node->ne );
1297  DeleteTree( rect, node->sw );
1298  DeleteTree( rect, node->se );
1299 
1300  if ( node->IsLeaf() )
1301  {
1302  delete node;
1303  node = nullptr;
1304  }
1305  }
1306  }
1307  }
1308 
1309  void DeleteTree( const point& pt, Node*& node )
1310  {
1311  if ( node != nullptr )
1312  {
1313  component x = pt[0];
1314  component y = pt[1];
1315  if ( node->Includes( x, y ) )
1316  if ( node->IsLeaf() )
1317  {
1318  LeafNode* leaf = static_cast<LeafNode*>( node );
1319  point_list points;
1320  for ( const point& p : leaf->points )
1321  {
1322  if ( p[0] == x )
1323  if ( p[1] == y )
1324  {
1325  --m_length;
1326  continue;
1327  }
1328  points << p;
1329  }
1330 
1331  if ( points.IsEmpty() )
1332  {
1333  delete leaf;
1334  node = nullptr;
1335  }
1336  else
1337  leaf->points = points;
1338  }
1339  else
1340  {
1341  DeleteTree( pt, node->nw );
1342  DeleteTree( pt, node->ne );
1343  DeleteTree( pt, node->sw );
1344  DeleteTree( pt, node->se );
1345 
1346  if ( node->IsLeaf() )
1347  {
1348  delete node;
1349  node = nullptr;
1350  }
1351  }
1352  }
1353  }
1354 
1355  void DestroyTree( Node* node )
1356  {
1357  if ( node != nullptr )
1358  if ( node->IsLeaf() )
1359  delete static_cast<LeafNode*>( node );
1360  else
1361  {
1362  DestroyTree( node->nw );
1363  DestroyTree( node->ne );
1364  DestroyTree( node->sw );
1365  DestroyTree( node->se );
1366  delete node;
1367  }
1368  }
1369 
1370  static void GetNumberOfNodes( size_type& N, const Node* node )
1371  {
1372  if ( node != nullptr )
1373  {
1374  ++N;
1375  GetNumberOfNodes( N, node->nw );
1376  GetNumberOfNodes( N, node->ne );
1377  GetNumberOfNodes( N, node->sw );
1378  GetNumberOfNodes( N, node->se );
1379  }
1380  }
1381 
1382  static void GetNumberOfLeafNodes( size_type& N, const Node* node )
1383  {
1384  if ( node != nullptr )
1385  {
1386  if ( node->IsLeaf() )
1387  ++N;
1388  GetNumberOfLeafNodes( N, node->nw );
1389  GetNumberOfLeafNodes( N, node->ne );
1390  GetNumberOfLeafNodes( N, node->sw );
1391  GetNumberOfLeafNodes( N, node->se );
1392  }
1393  }
1394 
1395  static int TreeHeight( const Node* node, int h )
1396  {
1397  if ( node == nullptr )
1398  return h;
1399  return h + 1 + Max( Max( Max( TreeHeight( node->nw, h ),
1400  TreeHeight( node->ne, h ) ),
1401  TreeHeight( node->sw, h ) ),
1402  TreeHeight( node->se, h ) );
1403  }
1404 };
1405 
1406 // ----------------------------------------------------------------------------
1407 
1408 } // pcl
1409 
1410 #endif // __PCL_QuadTree_h
1411 
1412 // ----------------------------------------------------------------------------
1413 // EOF pcl/QuadTree.h - Released 2020-12-17T15:46:29Z
Node * SplitAt(rectangle::point p)
Definition: QuadTree.h:659
int Length() const
Definition: QuadTree.h:304
void Search(const rectangle &rect, F callback, void *data) const
Definition: QuadTree.h:517
QuadTree()=default
rectangle SWRect() const
Definition: QuadTree.h:233
static size_type NumberOfNodes(const Node *node)
Definition: QuadTree.h:794
LeafNode * LeafNodeAt(rectangle::point p)
Definition: QuadTree.h:616
void Delete(const point &pt)
Definition: QuadTree.h:540
LeafNode(const rectangle &r, const point_list &p)
Definition: QuadTree.h:283
DRect rectangle
Definition: QuadTree.h:133
rectangle::component coordinate
Definition: QuadTree.h:138
rectangle SERect() const
Definition: QuadTree.h:242
int BucketCapacity() const
Definition: QuadTree.h:559
Node * sw
South-West child node, representing the bottom-left subregion.
Definition: QuadTree.h:151
PCL root namespace.
Definition: AbstractImage.h:76
friend void Swap(QuadTree &x1, QuadTree &x2)
Definition: QuadTree.h:848
Array< point > point_list
Definition: QuadTree.h:128
point::component component
Definition: QuadTree.h:123
void Build(const point_list &points, int bucketCapacity=40)
Definition: QuadTree.h:430
rectangle NWRect() const
Definition: QuadTree.h:217
size_type NumberOfLeafNodes() const
Definition: QuadTree.h:824
constexpr const T & Max(const T &a, const T &b) noexcept
Definition: Utility.h:119
void Swap(GenericPoint< T > &p1, GenericPoint< T > &p2) noexcept
Definition: Point.h:1404
static void Traverse(const Node *node, F f)
Definition: QuadTree.h:709
static size_type NumberOfLeafNodes(const Node *node)
Definition: QuadTree.h:814
LeafNode(const rectangle &r, const point &p)
Definition: QuadTree.h:294
rectangle NERect() const
Definition: QuadTree.h:225
void Clear()
Definition: QuadTree.h:408
rectangle rect
The rectangular region represented by this node.
Definition: QuadTree.h:148
void Traverse(F f) const
Definition: QuadTree.h:770
size_t size_type
Definition: Defs.h:605
bool IsEmpty() const noexcept
Definition: Array.h:314
QuadTree(const rectangle &rect, const point_list &points, int bucketCapacity=40)
Definition: QuadTree.h:349
64-bit floating-point rectangle in the R^2 space.
void Delete(const rectangle &rect)
Definition: QuadTree.h:549
size_type Length() const noexcept
Definition: Array.h:268
size_type Length() const
Definition: QuadTree.h:567
Node * Root()
Definition: QuadTree.h:596
Node * NodeAt(rectangle::point p)
Definition: QuadTree.h:644
Node(const rectangle &r)
Definition: QuadTree.h:166
static int Height(const Node *node)
Definition: QuadTree.h:832
bool IsEmpty() const
Definition: QuadTree.h:575
Quadtree node structure.
Definition: QuadTree.h:146
Node * nw
North-West child node, representing the top-left subregion.
Definition: QuadTree.h:149
bool Intersects(const rectangle &r) const
Definition: QuadTree.h:189
int Height() const
Definition: QuadTree.h:840
Quadtree leaf node structure.
Definition: QuadTree.h:254
static void Traverse(Node *node, F f)
Definition: QuadTree.h:745
QuadTree(QuadTree &&x)
Definition: QuadTree.h:371
const Node * NodeAt(rectangle::point p) const
Definition: QuadTree.h:630
Bucket PR quadtree for two-dimensional point data.
Definition: QuadTree.h:111
const Node * Root() const
Definition: QuadTree.h:587
size_type NumberOfNodes() const
Definition: QuadTree.h:805
QuadTree(const point_list &points, int bucketCapacity=40)
Definition: QuadTree.h:329
void Insert(const point &pt)
Definition: QuadTree.h:525
void Build(const rectangle &rect, const point_list &points, int bucketCapacity=40)
Definition: QuadTree.h:477
bool IsLeaf() const
Definition: QuadTree.h:180
Node * se
South-East child node, representing the bottom-right subregion.
Definition: QuadTree.h:152
point_list Search(const rectangle &rect) const
Definition: QuadTree.h:493
const LeafNode * LeafNodeAt(rectangle::point p) const
Definition: QuadTree.h:606
void Traverse(F f)
Definition: QuadTree.h:785
bool Includes(coordinate x, coordinate y) const
Definition: QuadTree.h:199
Node * ne
North-East child node, representing the top-right subregion.
Definition: QuadTree.h:150
QuadTree & operator=(const QuadTree &)=delete
bool Includes(const rectangle::point &p) const
Definition: QuadTree.h:209