QGIS API Documentation  3.8.0-Zanzibar (11aff65)
rtree.hpp
Go to the documentation of this file.
1 /*
2  * R-Tree index
3  *
4  * from http://www.superliminal.com/
5  *
6  * LICENSE : Entirely free for all uses. Enjoy!
7  * This File is in the public domain
8  *
9  * AUTORS
10  * - 1983 Original algorithm and test code by Antonin Guttman and Michael Stonebraker, UC Berkely
11  * - 1994 ANCI C ported from original test code by Melinda Green - [email protected]
12  * - 1995 Sphere volume fix for degeneracy problem submitted by Paul Brook
13  * - 2004 Templated C++ port by Greg Douglas
14  * - 2008 Portability issues fixed by Maxence Laurent
15  */
16 
17 #ifndef RTREE_H
18 #define RTREE_H
19 
20 #include "qgis.h"
21 #include <cstdio>
22 #include <cmath>
23 #include <cassert>
24 #include <QtGlobal>
25 
27 
28 #define ASSERT assert // RTree uses ASSERT( condition )
29 
30 //
31 // RTree.h
32 //
33 
34 #define RTREE_TEMPLATE template<class DATATYPE, class ELEMTYPE, int NUMDIMS, class ELEMTYPEREAL, int TMAXNODES, int TMINNODES>
35 #define RTREE_QUAL RTree<DATATYPE, ELEMTYPE, NUMDIMS, ELEMTYPEREAL, TMAXNODES, TMINNODES>
36 
37 #define RTREE_DONT_USE_MEMPOOLS // This version does not contain a fixed memory allocator, fill in lines with EXAMPLE to implement one.
38 #define RTREE_USE_SPHERICAL_VOLUME // Better split classification, may be slower on some systems
39 
40 namespace pal
41 {
42 
43 // Fwd decl
44  class RTFileStream; // File I/O helper class, look below for implementation and notes.
45 
46 
65  template < class DATATYPE, class ELEMTYPE, int NUMDIMS,
66  class ELEMTYPEREAL = ELEMTYPE, int TMAXNODES = 8, int TMINNODES = TMAXNODES / 2 >
67  class RTree
68  {
69  protected:
70 
71  struct Node; // Fwd decl. Used by other internal structs and iterator
72 
73  public:
74 
75  // These constant must be declared after Branch and before Node struct
76  // Stuck up here for MSVC 6 compiler. NSVC .NET 2003 is much happier.
77  enum
78  {
79  MAXNODES = TMAXNODES,
80  MINNODES = TMINNODES
81  };
82 
83 
84  public:
85 
86  RTree();
87  virtual ~RTree();
88 
93  void Insert( const ELEMTYPE a_min[NUMDIMS], const ELEMTYPE a_max[NUMDIMS], const DATATYPE &a_dataId );
94 
99  void Remove( const ELEMTYPE a_min[NUMDIMS], const ELEMTYPE a_max[NUMDIMS], const DATATYPE &a_dataId );
100 
108  int Search( const ELEMTYPE a_min[NUMDIMS], const ELEMTYPE a_max[NUMDIMS], bool a_resultCallback( DATATYPE a_data, void *a_context ), void *a_context );
109 
111  void RemoveAll();
112 
114  int Count();
115 
117  bool Load( const char *a_fileName );
119  bool Load( RTFileStream &a_stream );
120 
121 
123  bool Save( const char *a_fileName );
125  bool Save( RTFileStream &a_stream );
126 
131  class Iterator
132  {
133  private:
134 
135  enum { MAX_STACK = 32 }; // Max stack size. Allows almost n^32 where n is number of branches in node
136 
137  struct StackElement
138  {
139  Node *m_node = nullptr;
140  int m_branchIndex = 0;
141  };
142 
143  public:
144 
145  Iterator() { Init(); }
146 
148  bool IsNull() const { return ( m_tos <= 0 ); }
149 
151  bool IsNotNull() const { return ( m_tos > 0 ); }
152 
154  DATATYPE &operator*()
155  {
156  ASSERT( IsNotNull() );
157  StackElement &curTos = m_stack[m_tos - 1];
158  return curTos.m_node->m_branch[curTos.m_branchIndex].m_data;
159  }
160 
162  const DATATYPE &operator*() const
163  {
164  ASSERT( IsNotNull() );
165  StackElement &curTos = m_stack[m_tos - 1];
166  return curTos.m_node->m_branch[curTos.m_branchIndex].m_data;
167  }
168 
170  bool operator++() { return FindNextData(); }
171 
172  private:
173 
175  void Init() { m_tos = 0; }
176 
178  bool FindNextData()
179  {
180  for ( ;; )
181  {
182  if ( m_tos <= 0 )
183  {
184  return false;
185  }
186  StackElement curTos = Pop(); // Copy stack top cause it may change as we use it
187 
188  if ( curTos.m_node->IsLeaf() )
189  {
190  // Keep walking through data while we can
191  if ( curTos.m_branchIndex + 1 < curTos.m_node->m_count )
192  {
193  // There is more data, just point to the next one
194  Push( curTos.m_node, curTos.m_branchIndex + 1 );
195  return true;
196  }
197  // No more data, so it will fall back to previous level
198  }
199  else
200  {
201  if ( curTos.m_branchIndex + 1 < curTos.m_node->m_count )
202  {
203  // Push sibling on for future tree walk
204  // This is the 'fall back' node when we finish with the current level
205  Push( curTos.m_node, curTos.m_branchIndex + 1 );
206  }
207  // Since cur node is not a leaf, push first of next level to get deeper into the tree
208  Node *nextLevelnode = curTos.m_node->m_branch[curTos.m_branchIndex].m_child;
209  Push( nextLevelnode, 0 );
210 
211  // If we pushed on a new leaf, exit as the data is ready at TOS
212  if ( nextLevelnode->IsLeaf() )
213  {
214  return true;
215  }
216  }
217  }
218  }
219 
221  void Push( Node *a_node, int a_branchIndex )
222  {
223  m_stack[m_tos].m_node = a_node;
224  m_stack[m_tos].m_branchIndex = a_branchIndex;
225  ++m_tos;
226  ASSERT( m_tos <= MAX_STACK );
227  }
228 
230  StackElement &Pop()
231  {
232  ASSERT( m_tos > 0 );
233  --m_tos;
234  return m_stack[m_tos];
235  }
236 
237  StackElement m_stack[MAX_STACK];
238  int m_tos;
239 
240  friend class RTree; // Allow hiding of non-public functions while allowing manipulation by logical owner
241  };
242 
244  void GetFirst( Iterator &a_it )
245  {
246  a_it.Init();
247  if ( m_root && ( m_root->m_count > 0 ) )
248  {
249  a_it.Push( m_root, 0 );
250  a_it.FindNextData();
251  }
252  }
253 
255  void GetNext( Iterator &a_it ) { ++a_it; }
256 
258  bool IsNull( Iterator &a_it ) { return a_it.IsNull(); }
259 
261  DATATYPE &GetAt( Iterator &a_it ) { return *a_it; }
262 
263  protected:
264 
266  struct Rect
267  {
268  ELEMTYPE m_min[NUMDIMS];
269  ELEMTYPE m_max[NUMDIMS];
270  };
271 
275  struct Branch
276  {
277  Rect m_rect;
278  union
279  {
280  Node *m_child;
281  DATATYPE m_data;
282  };
283  };
284 
286  struct Node
287  {
288  bool IsInternalNode() { return ( m_level > 0 ); } // Not a leaf, but a internal node
289  bool IsLeaf() { return ( m_level == 0 ); } // A leaf, contains data
290 
291  int m_count;
292  int m_level;
293  Branch m_branch[MAXNODES];
294  };
295 
297  struct ListNode
298  {
299  ListNode *m_next;
300  Node *m_node;
301  };
302 
304  struct PartitionVars
305  {
306  int m_partition[MAXNODES + 1];
307  int m_total;
308  int m_minFill;
309  int m_taken[MAXNODES + 1];
310  int m_count[2];
311  Rect m_cover[2];
312  ELEMTYPEREAL m_area[2];
313 
314  Branch m_branchBuf[MAXNODES + 1];
315  int m_branchCount;
316  Rect m_coverSplit;
317  ELEMTYPEREAL m_coverSplitArea;
318  };
319 
320  Node *AllocNode();
321  void FreeNode( Node *a_node );
322  void InitNode( Node *a_node );
323  void InitRect( Rect *a_rect );
324  bool InsertRectRec( Rect *a_rect, const DATATYPE &a_id, Node *a_node, Node **a_newNode, int a_level );
325  bool InsertRect( Rect *a_rect, const DATATYPE &a_id, Node **a_root, int a_level );
326  Rect NodeCover( Node *a_node );
327  bool AddBranch( Branch *a_branch, Node *a_node, Node **a_newNode );
328  void DisconnectBranch( Node *a_node, int a_index );
329  int PickBranch( Rect *a_rect, Node *a_node );
330  Rect CombineRect( Rect *a_rectA, Rect *a_rectB );
331  void SplitNode( Node *a_node, Branch *a_branch, Node **a_newNode );
332  ELEMTYPEREAL RectSphericalVolume( Rect *a_rect );
333  ELEMTYPEREAL RectVolume( Rect *a_rect );
334  ELEMTYPEREAL CalcRectVolume( Rect *a_rect );
335  void GetBranches( Node *a_node, Branch *a_branch, PartitionVars *a_parVars );
336  void ChoosePartition( PartitionVars *a_parVars, int a_minFill );
337  void LoadNodes( Node *a_nodeA, Node *a_nodeB, PartitionVars *a_parVars );
338  void InitParVars( PartitionVars *a_parVars, int a_maxRects, int a_minFill );
339  void PickSeeds( PartitionVars *a_parVars );
340  void Classify( int a_index, int a_group, PartitionVars *a_parVars );
341  bool RemoveRect( Rect *a_rect, const DATATYPE &a_id, Node **a_root );
342  bool RemoveRectRec( Rect *a_rect, const DATATYPE &a_id, Node *a_node, ListNode **a_listNode );
343  ListNode *AllocListNode();
344  void FreeListNode( ListNode *a_listNode );
345  bool Overlap( Rect *a_rectA, Rect *a_rectB );
346  void ReInsert( Node *a_node, ListNode **a_listNode );
347  bool Search( Node *a_node, Rect *a_rect, int &a_foundCount, bool a_resultCallback( DATATYPE a_data, void *a_context ), void *a_context );
348  void RemoveAllRec( Node *a_node );
349  void Reset();
350  void CountRec( Node *a_node, int &a_count );
351 
352  bool SaveRec( Node *a_node, RTFileStream &a_stream );
353  bool LoadRec( Node *a_node, RTFileStream &a_stream );
354 
355  Node *m_root;
356  ELEMTYPEREAL m_unitSphereVolume;
357  };
358 
359 
360  // Because there is not stream support, this is a quick and dirty file I/O helper.
361  // Users will likely replace its usage with a Stream implementation from their favorite API.
362 
366  class RTFileStream
367  {
368  FILE *m_file = nullptr;
369 
370  public:
371 
372 
373  RTFileStream()
374  {
375  m_file = nullptr;
376  }
377 
378  ~RTFileStream()
379  {
380  Close();
381  }
382 
384  RTFileStream( const RTFileStream &other ) = delete;
386  RTFileStream &operator=( const RTFileStream &other ) = delete;
387 
388  bool OpenRead( const char *a_fileName )
389  {
390  if ( m_file )
391  fclose( m_file );
392 
393  m_file = fopen( a_fileName, "rb" );
394  return m_file != nullptr;
395  }
396 
397  bool OpenWrite( const char *a_fileName )
398  {
399  if ( m_file )
400  fclose( m_file );
401 
402  m_file = fopen( a_fileName, "wb" );
403  return m_file != nullptr;
404  }
405 
406  void Close()
407  {
408  if ( m_file )
409  {
410  fclose( m_file );
411  m_file = nullptr;
412  }
413  }
414 
415  template< typename TYPE >
416  size_t Write( const TYPE &a_value )
417  {
418  ASSERT( m_file );
419  return fwrite( static_cast< void * >( &a_value ), sizeof( a_value ), 1, m_file );
420  }
421 
422  template< typename TYPE >
423  size_t WriteArray( const TYPE *a_array, int a_count )
424  {
425  ASSERT( m_file );
426  return fwrite( static_cast< void * >( a_array ), sizeof( TYPE ) * a_count, 1, m_file );
427  }
428 
429  template< typename TYPE >
430  size_t Read( TYPE &a_value )
431  {
432  ASSERT( m_file );
433  return fread( static_cast< void * >( &a_value ), sizeof( a_value ), 1, m_file );
434  }
435 
436  template< typename TYPE >
437  size_t ReadArray( TYPE *a_array, int a_count )
438  {
439  ASSERT( m_file );
440  return fread( static_cast< void * >( a_array ), sizeof( TYPE ) * a_count, 1, m_file );
441  }
442 
443  };
444 
445 
446  RTREE_TEMPLATE
447  RTREE_QUAL::RTree()
448  {
449  ASSERT( MAXNODES > MINNODES );
450  ASSERT( MINNODES > 0 );
451 
452 
453  // We only support machine word size simple data type, e.g., integer index or object pointer.
454  // Since we are storing as union with non data branch
455  ASSERT( sizeof( DATATYPE ) == sizeof( void * ) || sizeof( DATATYPE ) == sizeof( int ) );
456 
457  // Precomputed volumes of the unit spheres for the first few dimensions
458  const float UNIT_SPHERE_VOLUMES[] =
459  {
460  0.000000f, 2.000000f, 3.141593f, // Dimension 0,1,2
461  4.188790f, 4.934802f, 5.263789f, // Dimension 3,4,5
462  5.167713f, 4.724766f, 4.058712f, // Dimension 6,7,8
463  3.298509f, 2.550164f, 1.884104f, // Dimension 9,10,11
464  1.335263f, 0.910629f, 0.599265f, // Dimension 12,13,14
465  0.381443f, 0.235331f, 0.140981f, // Dimension 15,16,17
466  0.082146f, 0.046622f, 0.025807f, // Dimension 18,19,20
467  };
468 
469  m_root = AllocNode();
470  m_root->m_level = 0;
471  m_unitSphereVolume = static_cast< ELEMTYPEREAL >( UNIT_SPHERE_VOLUMES[NUMDIMS] );
472  }
473 
474 
475  RTREE_TEMPLATE
476  RTREE_QUAL::~RTree()
477  {
478  Reset(); // Free, or reset node memory
479  }
480 
481 
482  RTREE_TEMPLATE
483  void RTREE_QUAL::Insert( const ELEMTYPE a_min[NUMDIMS], const ELEMTYPE a_max[NUMDIMS], const DATATYPE &a_dataId )
484  {
485 #ifdef _DEBUG
486  for ( int index = 0; index < NUMDIMS; ++index )
487  {
488  ASSERT( a_min[index] <= a_max[index] );
489  }
490 #endif //_DEBUG
491 
492  Rect rect;
493 
494  for ( int axis = 0; axis < NUMDIMS; ++axis )
495  {
496  rect.m_min[axis] = a_min[axis];
497  rect.m_max[axis] = a_max[axis];
498  }
499 
500  InsertRect( &rect, a_dataId, &m_root, 0 );
501  }
502 
503 
504  RTREE_TEMPLATE
505  void RTREE_QUAL::Remove( const ELEMTYPE a_min[NUMDIMS], const ELEMTYPE a_max[NUMDIMS], const DATATYPE &a_dataId )
506  {
507 #ifdef _DEBUG
508  for ( int index = 0; index < NUMDIMS; ++index )
509  {
510  ASSERT( a_min[index] <= a_max[index] );
511  }
512 #endif //_DEBUG
513 
514  Rect rect;
515 
516  for ( int axis = 0; axis < NUMDIMS; ++axis )
517  {
518  rect.m_min[axis] = a_min[axis];
519  rect.m_max[axis] = a_max[axis];
520  }
521 
522  RemoveRect( &rect, a_dataId, &m_root );
523  }
524 
525 
526  RTREE_TEMPLATE
527  int RTREE_QUAL::Search( const ELEMTYPE a_min[NUMDIMS], const ELEMTYPE a_max[NUMDIMS], bool a_resultCallback( DATATYPE a_data, void *a_context ), void *a_context )
528  {
529 #ifdef _DEBUG
530  for ( int index = 0; index < NUMDIMS; ++index )
531  {
532  ASSERT( a_min[index] <= a_max[index] );
533  }
534 #endif //_DEBUG
535 
536  Rect rect;
537 
538  for ( int axis = 0; axis < NUMDIMS; ++axis )
539  {
540  rect.m_min[axis] = a_min[axis];
541  rect.m_max[axis] = a_max[axis];
542  }
543 
544  // NOTE: May want to return search result another way, perhaps returning the number of found elements here.
545 
546  int foundCount = 0;
547  Search( m_root, &rect, foundCount, a_resultCallback, a_context );
548 
549  return foundCount;
550  }
551 
552 
553  RTREE_TEMPLATE
554  int RTREE_QUAL::Count()
555  {
556  int count = 0;
557  CountRec( m_root, count );
558 
559  return count;
560  }
561 
562 
563 
564  RTREE_TEMPLATE
565  void RTREE_QUAL::CountRec( Node *a_node, int &a_count )
566  {
567  if ( a_node->IsInternalNode() ) // not a leaf node
568  {
569  for ( int index = 0; index < a_node->m_count; ++index )
570  {
571  CountRec( a_node->m_branch[index].m_child, a_count );
572  }
573  }
574  else // A leaf node
575  {
576  a_count += a_node->m_count;
577  }
578  }
579 
580 
581  RTREE_TEMPLATE
582  bool RTREE_QUAL::Load( const char *a_fileName )
583  {
584  RemoveAll(); // Clear existing tree
585 
586  RTFileStream stream;
587  if ( !stream.OpenRead( a_fileName ) )
588  {
589  return false;
590  }
591 
592  bool result = Load( stream );
593 
594  stream.Close();
595 
596  return result;
597  }
598 
599 
600 
601  RTREE_TEMPLATE
602  bool RTREE_QUAL::Load( RTFileStream &a_stream )
603  {
604  // Write some kind of header
605  int _dataFileId = ( 'R' << 0 ) | ( 'T' << 8 ) | ( 'R' << 16 ) | ( 'E' << 24 );
606  int _dataSize = sizeof( DATATYPE );
607  int _dataNumDims = NUMDIMS;
608  int _dataElemSize = sizeof( ELEMTYPE );
609  int _dataElemRealSize = sizeof( ELEMTYPEREAL );
610  int _dataMaxNodes = TMAXNODES;
611  int _dataMinNodes = TMINNODES;
612 
613  int dataFileId = 0;
614  int dataSize = 0;
615  int dataNumDims = 0;
616  int dataElemSize = 0;
617  int dataElemRealSize = 0;
618  int dataMaxNodes = 0;
619  int dataMinNodes = 0;
620 
621  a_stream.Read( dataFileId );
622  a_stream.Read( dataSize );
623  a_stream.Read( dataNumDims );
624  a_stream.Read( dataElemSize );
625  a_stream.Read( dataElemRealSize );
626  a_stream.Read( dataMaxNodes );
627  a_stream.Read( dataMinNodes );
628 
629  bool result = false;
630 
631  // Test if header was valid and compatible
632  if ( ( dataFileId == _dataFileId )
633  && ( dataSize == _dataSize )
634  && ( dataNumDims == _dataNumDims )
635  && ( dataElemSize == _dataElemSize )
636  && ( dataElemRealSize == _dataElemRealSize )
637  && ( dataMaxNodes == _dataMaxNodes )
638  && ( dataMinNodes == _dataMinNodes )
639  )
640  {
641  // Recursively load tree
642  result = LoadRec( m_root, a_stream );
643  }
644 
645  return result;
646  }
647 
648 
649  RTREE_TEMPLATE
650  bool RTREE_QUAL::LoadRec( Node *a_node, RTFileStream &a_stream )
651  {
652  a_stream.Read( a_node->m_level );
653  a_stream.Read( a_node->m_count );
654 
655  if ( a_node->IsInternalNode() ) // not a leaf node
656  {
657  for ( int index = 0; index < a_node->m_count; ++index )
658  {
659  Branch *curBranch = &a_node->m_branch[index];
660 
661  a_stream.ReadArray( curBranch->m_rect.m_min, NUMDIMS );
662  a_stream.ReadArray( curBranch->m_rect.m_max, NUMDIMS );
663 
664  curBranch->m_child = AllocNode();
665  LoadRec( curBranch->m_child, a_stream );
666  }
667  }
668  else // A leaf node
669  {
670  for ( int index = 0; index < a_node->m_count; ++index )
671  {
672  Branch *curBranch = &a_node->m_branch[index];
673 
674  a_stream.ReadArray( curBranch->m_rect.m_min, NUMDIMS );
675  a_stream.ReadArray( curBranch->m_rect.m_max, NUMDIMS );
676 
677  a_stream.Read( curBranch->m_data );
678  }
679  }
680 
681  return true; // Should do more error checking on I/O operations
682  }
683 
684 
685  RTREE_TEMPLATE
686  bool RTREE_QUAL::Save( const char *a_fileName )
687  {
688  RTFileStream stream;
689  if ( !stream.OpenWrite( a_fileName ) )
690  {
691  return false;
692  }
693 
694  bool result = Save( stream );
695 
696  stream.Close();
697 
698  return result;
699  }
700 
701 
702  RTREE_TEMPLATE
703  bool RTREE_QUAL::Save( RTFileStream &a_stream )
704  {
705  // Write some kind of header
706  int dataFileId = ( 'R' << 0 ) | ( 'T' << 8 ) | ( 'R' << 16 ) | ( 'E' << 24 );
707  int dataSize = sizeof( DATATYPE );
708  int dataNumDims = NUMDIMS;
709  int dataElemSize = sizeof( ELEMTYPE );
710  int dataElemRealSize = sizeof( ELEMTYPEREAL );
711  int dataMaxNodes = TMAXNODES;
712  int dataMinNodes = TMINNODES;
713 
714  a_stream.Write( dataFileId );
715  a_stream.Write( dataSize );
716  a_stream.Write( dataNumDims );
717  a_stream.Write( dataElemSize );
718  a_stream.Write( dataElemRealSize );
719  a_stream.Write( dataMaxNodes );
720  a_stream.Write( dataMinNodes );
721 
722  // Recursively save tree
723  bool result = SaveRec( m_root, a_stream );
724 
725  return result;
726  }
727 
728 
729  RTREE_TEMPLATE
730  bool RTREE_QUAL::SaveRec( Node *a_node, RTFileStream &a_stream )
731  {
732  a_stream.Write( a_node->m_level );
733  a_stream.Write( a_node->m_count );
734 
735  if ( a_node->IsInternalNode() ) // not a leaf node
736  {
737  for ( int index = 0; index < a_node->m_count; ++index )
738  {
739  Branch *curBranch = &a_node->m_branch[index];
740 
741  a_stream.WriteArray( curBranch->m_rect.m_min, NUMDIMS );
742  a_stream.WriteArray( curBranch->m_rect.m_max, NUMDIMS );
743 
744  SaveRec( curBranch->m_child, a_stream );
745  }
746  }
747  else // A leaf node
748  {
749  for ( int index = 0; index < a_node->m_count; ++index )
750  {
751  Branch *curBranch = &a_node->m_branch[index];
752 
753  a_stream.WriteArray( curBranch->m_rect.m_min, NUMDIMS );
754  a_stream.WriteArray( curBranch->m_rect.m_max, NUMDIMS );
755 
756  a_stream.Write( curBranch->m_data );
757  }
758  }
759 
760  return true; // Should do more error checking on I/O operations
761  }
762 
763 
764  RTREE_TEMPLATE
765  void RTREE_QUAL::RemoveAll()
766  {
767  // Delete all existing nodes
768  Reset();
769 
770  m_root = AllocNode();
771  m_root->m_level = 0;
772  }
773 
774 
775  RTREE_TEMPLATE
776  void RTREE_QUAL::Reset()
777  {
778 #ifdef RTREE_DONT_USE_MEMPOOLS
779  // Delete all existing nodes
780  RemoveAllRec( m_root );
781 #else // RTREE_DONT_USE_MEMPOOLS
782  // Just reset memory pools. We are not using complex types
783  // EXAMPLE
784 #endif // RTREE_DONT_USE_MEMPOOLS
785  }
786 
787 
788  RTREE_TEMPLATE
789  void RTREE_QUAL::RemoveAllRec( Node *a_node )
790  {
791  ASSERT( a_node );
792  ASSERT( a_node->m_level >= 0 );
793 
794  if ( a_node->IsInternalNode() ) // This is an internal node in the tree
795  {
796  for ( int index = 0; index < a_node->m_count; ++index )
797  {
798  RemoveAllRec( a_node->m_branch[index].m_child );
799  }
800  }
801  FreeNode( a_node );
802  }
803 
804 
805  RTREE_TEMPLATE
806  typename RTREE_QUAL::Node *RTREE_QUAL::AllocNode()
807  {
808  Node *newNode = nullptr;
809 #ifdef RTREE_DONT_USE_MEMPOOLS
810  newNode = new Node;
811 #else // RTREE_DONT_USE_MEMPOOLS
812  // EXAMPLE
813 #endif // RTREE_DONT_USE_MEMPOOLS
814  InitNode( newNode );
815  return newNode;
816  }
817 
818 
819  RTREE_TEMPLATE
820  void RTREE_QUAL::FreeNode( Node *a_node )
821  {
822  ASSERT( a_node );
823 
824 #ifdef RTREE_DONT_USE_MEMPOOLS
825  delete a_node;
826 #else // RTREE_DONT_USE_MEMPOOLS
827  // EXAMPLE
828 #endif // RTREE_DONT_USE_MEMPOOLS
829  }
830 
831 
832 // Allocate space for a node in the list used in DeletRect to
833 // store Nodes that are too empty.
834  RTREE_TEMPLATE
835  typename RTREE_QUAL::ListNode *RTREE_QUAL::AllocListNode()
836  {
837 #ifdef RTREE_DONT_USE_MEMPOOLS
838  return new ListNode;
839 #else // RTREE_DONT_USE_MEMPOOLS
840  // EXAMPLE
841 #endif // RTREE_DONT_USE_MEMPOOLS
842  }
843 
844 
845  RTREE_TEMPLATE
846  void RTREE_QUAL::FreeListNode( ListNode *a_listNode )
847  {
848 #ifdef RTREE_DONT_USE_MEMPOOLS
849  delete a_listNode;
850 #else // RTREE_DONT_USE_MEMPOOLS
851  // EXAMPLE
852 #endif // RTREE_DONT_USE_MEMPOOLS
853  }
854 
855 
856  RTREE_TEMPLATE
857  void RTREE_QUAL::InitNode( Node *a_node )
858  {
859  a_node->m_count = 0;
860  a_node->m_level = -1;
861  }
862 
863 
864  RTREE_TEMPLATE
865  void RTREE_QUAL::InitRect( Rect *a_rect )
866  {
867  for ( int index = 0; index < NUMDIMS; ++index )
868  {
869  a_rect->m_min[index] = static_cast< ELEMTYPE >( 0 );
870  a_rect->m_max[index] = static_cast< ELEMTYPE >( 0 );
871  }
872  }
873 
874 
875 // Inserts a new data rectangle into the index structure.
876 // Recursively descends tree, propagates splits back up.
877 // Returns 0 if node was not split. Old node updated.
878 // If node was split, returns 1 and sets the pointer pointed to by
879 // new_node to point to the new node. Old node updated to become one of two.
880 // The level argument specifies the number of steps up from the leaf
881 // level to insert; e.g. a data rectangle goes in at level = 0.
882  RTREE_TEMPLATE
883  bool RTREE_QUAL::InsertRectRec( Rect *a_rect, const DATATYPE &a_id, Node *a_node, Node **a_newNode, int a_level )
884  {
885  ASSERT( a_rect && a_node && a_newNode );
886  ASSERT( a_level >= 0 && a_level <= a_node->m_level );
887 
888  int index;
889  Branch branch;
890  Node *otherNode = nullptr;
891 
892  // Still above level for insertion, go down tree recursively
893  if ( a_node->m_level > a_level )
894  {
895  index = PickBranch( a_rect, a_node );
896  if ( !InsertRectRec( a_rect, a_id, a_node->m_branch[index].m_child, &otherNode, a_level ) )
897  {
898  // Child was not split
899  a_node->m_branch[index].m_rect = CombineRect( a_rect, & ( a_node->m_branch[index].m_rect ) );
900  return false;
901  }
902  else // Child was split
903  {
904  a_node->m_branch[index].m_rect = NodeCover( a_node->m_branch[index].m_child );
905  branch.m_child = otherNode;
906  branch.m_rect = NodeCover( otherNode );
907  return AddBranch( &branch, a_node, a_newNode );
908  }
909  }
910  else if ( a_node->m_level == a_level ) // Have reached level for insertion. Add rect, split if necessary
911  {
912  branch.m_rect = *a_rect;
913  branch.m_child = reinterpret_cast< Node * >( a_id );
914  // Child field of leaves contains id of data record
915  return AddBranch( &branch, a_node, a_newNode );
916  }
917  else
918  {
919  // Should never occur
920  ASSERT( 0 );
921  return false;
922  }
923  }
924 
925 
926 // Insert a data rectangle into an index structure.
927 // InsertRect provides for splitting the root;
928 // returns 1 if root was split, 0 if it was not.
929 // The level argument specifies the number of steps up from the leaf
930 // level to insert; e.g. a data rectangle goes in at level = 0.
931 // InsertRect2 does the recursion.
932 //
933  RTREE_TEMPLATE
934  bool RTREE_QUAL::InsertRect( Rect *a_rect, const DATATYPE &a_id, Node **a_root, int a_level )
935  {
936  ASSERT( a_rect && a_root );
937  ASSERT( a_level >= 0 && a_level <= ( *a_root )->m_level );
938 #ifdef _DEBUG
939  for ( int index = 0; index < NUMDIMS; ++index )
940  {
941  ASSERT( a_rect->m_min[index] <= a_rect->m_max[index] );
942  }
943 #endif //_DEBUG
944 
945  Node *newRoot = nullptr;
946  Node *newNode = nullptr;
947  Branch branch;
948 
949  if ( InsertRectRec( a_rect, a_id, *a_root, &newNode, a_level ) ) // Root split
950  {
951  newRoot = AllocNode(); // Grow tree taller and new root
952  newRoot->m_level = ( *a_root )->m_level + 1;
953  branch.m_rect = NodeCover( *a_root );
954  branch.m_child = *a_root;
955  AddBranch( &branch, newRoot, nullptr );
956  branch.m_rect = NodeCover( newNode );
957  branch.m_child = newNode;
958  AddBranch( &branch, newRoot, nullptr );
959  *a_root = newRoot;
960  return true;
961  }
962 
963  return false;
964  }
965 
966 
967 // Find the smallest rectangle that includes all rectangles in branches of a node.
968  RTREE_TEMPLATE
969  typename RTREE_QUAL::Rect RTREE_QUAL::NodeCover( Node *a_node )
970  {
971  ASSERT( a_node );
972 
973  bool firstTime = true;
974  Rect rect;
975  InitRect( &rect );
976 
977  for ( int index = 0; index < a_node->m_count; ++index )
978  {
979  if ( firstTime )
980  {
981  rect = a_node->m_branch[index].m_rect;
982  firstTime = false;
983  }
984  else
985  {
986  rect = CombineRect( &rect, & ( a_node->m_branch[index].m_rect ) );
987  }
988  }
989 
990  return rect;
991  }
992 
993 
994 // Add a branch to a node. Split the node if necessary.
995 // Returns 0 if node not split. Old node updated.
996 // Returns 1 if node split, sets *new_node to address of new node.
997 // Old node updated, becomes one of two.
998  RTREE_TEMPLATE
999  bool RTREE_QUAL::AddBranch( Branch *a_branch, Node *a_node, Node **a_newNode )
1000  {
1001  ASSERT( a_branch );
1002  ASSERT( a_node );
1003 
1004  if ( a_node->m_count < MAXNODES ) // Split won't be necessary
1005  {
1006  a_node->m_branch[a_node->m_count] = *a_branch;
1007  ++a_node->m_count;
1008 
1009  return false;
1010  }
1011  else
1012  {
1013  ASSERT( a_newNode );
1014 
1015  SplitNode( a_node, a_branch, a_newNode );
1016  return true;
1017  }
1018  }
1019 
1020 
1021 // Disconnect a dependent node.
1022 // Caller must return (or stop using iteration index) after this as count has changed
1023  RTREE_TEMPLATE
1024  void RTREE_QUAL::DisconnectBranch( Node *a_node, int a_index )
1025  {
1026  ASSERT( a_node && ( a_index >= 0 ) && ( a_index < MAXNODES ) );
1027  ASSERT( a_node->m_count > 0 );
1028 
1029  // Remove element by swapping with the last element to prevent gaps in array
1030  a_node->m_branch[a_index] = a_node->m_branch[a_node->m_count - 1];
1031 
1032  --a_node->m_count;
1033  }
1034 
1035 
1036 // Pick a branch. Pick the one that will need the smallest increase
1037 // in area to accommodate the new rectangle. This will result in the
1038 // least total area for the covering rectangles in the current node.
1039 // In case of a tie, pick the one which was smaller before, to get
1040 // the best resolution when searching.
1041  RTREE_TEMPLATE
1042  int RTREE_QUAL::PickBranch( Rect *a_rect, Node *a_node )
1043  {
1044  ASSERT( a_rect && a_node );
1045 
1046  bool firstTime = true;
1047  ELEMTYPEREAL increase;
1048  ELEMTYPEREAL bestIncr = static_cast< ELEMTYPEREAL >( -1 );
1049  ELEMTYPEREAL area;
1050  ELEMTYPEREAL bestArea = 0;
1051  int best = 0;
1052  Rect tempRect;
1053 
1054  for ( int index = 0; index < a_node->m_count; ++index )
1055  {
1056  Rect *curRect = &a_node->m_branch[index].m_rect;
1057  area = CalcRectVolume( curRect );
1058  tempRect = CombineRect( a_rect, curRect );
1059  increase = CalcRectVolume( &tempRect ) - area;
1060  if ( ( increase < bestIncr ) || firstTime )
1061  {
1062  best = index;
1063  bestArea = area;
1064  bestIncr = increase;
1065  firstTime = false;
1066  }
1067  else if ( qgsDoubleNear( increase, bestIncr ) && ( area < bestArea ) )
1068  {
1069  best = index;
1070  bestArea = area;
1071  bestIncr = increase;
1072  }
1073  }
1074  return best;
1075  }
1076 
1077 
1078 // Combine two rectangles into larger one containing both
1079  RTREE_TEMPLATE
1080  typename RTREE_QUAL::Rect RTREE_QUAL::CombineRect( Rect *a_rectA, Rect *a_rectB )
1081  {
1082  ASSERT( a_rectA && a_rectB );
1083 
1084  Rect newRect = { {0, }, {0, } };
1085 
1086  for ( int index = 0; index < NUMDIMS; ++index )
1087  {
1088  newRect.m_min[index] = std::min( a_rectA->m_min[index], a_rectB->m_min[index] );
1089  newRect.m_max[index] = std::max( a_rectA->m_max[index], a_rectB->m_max[index] );
1090  }
1091 
1092  return newRect;
1093  }
1094 
1095 
1096 
1097 // Split a node.
1098 // Divides the nodes branches and the extra one between two nodes.
1099 // Old node is one of the new ones, and one really new one is created.
1100 // Tries more than one method for choosing a partition, uses best result.
1101  RTREE_TEMPLATE
1102  void RTREE_QUAL::SplitNode( Node *a_node, Branch *a_branch, Node **a_newNode )
1103  {
1104  ASSERT( a_node );
1105  ASSERT( a_branch );
1106 
1107  // Could just use local here, but member or external is faster since it is reused
1108  PartitionVars localVars;
1109  PartitionVars *parVars = &localVars;
1110  int level;
1111 
1112  // Load all the branches into a buffer, initialize old node
1113  level = a_node->m_level;
1114  GetBranches( a_node, a_branch, parVars );
1115 
1116  // Find partition
1117  ChoosePartition( parVars, MINNODES );
1118 
1119  // Put branches from buffer into 2 nodes according to chosen partition
1120  *a_newNode = AllocNode();
1121  ( *a_newNode )->m_level = a_node->m_level = level;
1122  LoadNodes( a_node, *a_newNode, parVars );
1123 
1124  ASSERT( ( a_node->m_count + ( *a_newNode )->m_count ) == parVars->m_total );
1125  }
1126 
1127 
1128 // Calculate the n-dimensional volume of a rectangle
1129  RTREE_TEMPLATE
1130  ELEMTYPEREAL RTREE_QUAL::RectVolume( Rect *a_rect )
1131  {
1132  ASSERT( a_rect );
1133 
1134  ELEMTYPEREAL volume = static_cast< ELEMTYPEREAL >( 1 );
1135 
1136  for ( int index = 0; index < NUMDIMS; ++index )
1137  {
1138  volume *= a_rect->m_max[index] - a_rect->m_min[index];
1139  }
1140 
1141  ASSERT( volume >= static_cast< ELEMTYPEREAL >( 0 ) );
1142 
1143  return volume;
1144  }
1145 
1146 
1147 // The exact volume of the bounding sphere for the given Rect
1148  RTREE_TEMPLATE
1149  ELEMTYPEREAL RTREE_QUAL::RectSphericalVolume( Rect *a_rect )
1150  {
1151  ASSERT( a_rect );
1152 
1153  ELEMTYPEREAL sumOfSquares = static_cast< ELEMTYPEREAL >( 0 );
1154  ELEMTYPEREAL radius;
1155 
1156  for ( int index = 0; index < NUMDIMS; ++index )
1157  {
1158  ELEMTYPEREAL halfExtent = ( static_cast< ELEMTYPEREAL >( a_rect->m_max[index] ) - static_cast< ELEMTYPEREAL >( a_rect->m_min[index] ) ) * 0.5f;
1159  sumOfSquares += halfExtent * halfExtent;
1160  }
1161 
1162  radius = static_cast< ELEMTYPEREAL >( std::sqrt( sumOfSquares ) );
1163 
1164  // Pow maybe slow, so test for common dims like 2,3 and just use x*x, x*x*x.
1165  if ( NUMDIMS == 3 )
1166  {
1167  return ( radius * radius * radius * m_unitSphereVolume );
1168  }
1169  else if ( NUMDIMS == 2 )
1170  {
1171  return ( radius * radius * m_unitSphereVolume );
1172  }
1173  else
1174  {
1175  return static_cast< ELEMTYPEREAL >( std::pow( radius, NUMDIMS ) * m_unitSphereVolume );
1176  }
1177  }
1178 
1179 
1180 // Use one of the methods to calculate retangle volume
1181  RTREE_TEMPLATE
1182  ELEMTYPEREAL RTREE_QUAL::CalcRectVolume( Rect *a_rect )
1183  {
1184 #ifdef RTREE_USE_SPHERICAL_VOLUME
1185  return RectSphericalVolume( a_rect ); // Slower but helps certain merge cases
1186 #else // RTREE_USE_SPHERICAL_VOLUME
1187  return RectVolume( a_rect ); // Faster but can cause poor merges
1188 #endif // RTREE_USE_SPHERICAL_VOLUME
1189  }
1190 
1191 
1192 // Load branch buffer with branches from full node plus the extra branch.
1193  RTREE_TEMPLATE
1194  void RTREE_QUAL::GetBranches( Node *a_node, Branch *a_branch, PartitionVars *a_parVars )
1195  {
1196  ASSERT( a_node );
1197  ASSERT( a_branch );
1198 
1199  ASSERT( a_node->m_count == MAXNODES );
1200 
1201  int index;
1202  // Load the branch buffer
1203  for ( index = 0; index < MAXNODES; ++index )
1204  {
1205  a_parVars->m_branchBuf[index] = a_node->m_branch[index];
1206  }
1207  a_parVars->m_branchBuf[MAXNODES] = *a_branch;
1208  a_parVars->m_branchCount = MAXNODES + 1;
1209 
1210  // Calculate rect containing all in the set
1211  a_parVars->m_coverSplit = a_parVars->m_branchBuf[0].m_rect;
1212  for ( index = 1; index < MAXNODES + 1; ++index )
1213  {
1214  a_parVars->m_coverSplit = CombineRect( &a_parVars->m_coverSplit, &a_parVars->m_branchBuf[index].m_rect );
1215  }
1216  a_parVars->m_coverSplitArea = CalcRectVolume( &a_parVars->m_coverSplit );
1217 
1218  InitNode( a_node );
1219  }
1220 
1221 
1222 // Method #0 for choosing a partition:
1223 // As the seeds for the two groups, pick the two rects that would waste the
1224 // most area if covered by a single rectangle, i.e. evidently the worst pair
1225 // to have in the same group.
1226 // Of the remaining, one at a time is chosen to be put in one of the two groups.
1227 // The one chosen is the one with the greatest difference in area expansion
1228 // depending on which group - the rect most strongly attracted to one group
1229 // and repelled from the other.
1230 // If one group gets too full (more would force other group to violate min
1231 // fill requirement) then other group gets the rest.
1232 // These last are the ones that can go in either group most easily.
1233  RTREE_TEMPLATE
1234  void RTREE_QUAL::ChoosePartition( PartitionVars *a_parVars, int a_minFill )
1235  {
1236  ASSERT( a_parVars );
1237 
1238  ELEMTYPEREAL biggestDiff;
1239  int group, chosen = 0, betterGroup = 0;
1240 
1241  InitParVars( a_parVars, a_parVars->m_branchCount, a_minFill );
1242  PickSeeds( a_parVars );
1243 
1244  while ( ( ( a_parVars->m_count[0] + a_parVars->m_count[1] ) < a_parVars->m_total )
1245  && ( a_parVars->m_count[0] < ( a_parVars->m_total - a_parVars->m_minFill ) )
1246  && ( a_parVars->m_count[1] < ( a_parVars->m_total - a_parVars->m_minFill ) ) )
1247  {
1248  biggestDiff = static_cast< ELEMTYPEREAL >( -1 );
1249  for ( int index = 0; index < a_parVars->m_total; ++index )
1250  {
1251  if ( !a_parVars->m_taken[index] )
1252  {
1253  Rect *curRect = &a_parVars->m_branchBuf[index].m_rect;
1254  Rect rect0 = CombineRect( curRect, &a_parVars->m_cover[0] );
1255  Rect rect1 = CombineRect( curRect, &a_parVars->m_cover[1] );
1256  ELEMTYPEREAL growth0 = CalcRectVolume( &rect0 ) - a_parVars->m_area[0];
1257  ELEMTYPEREAL growth1 = CalcRectVolume( &rect1 ) - a_parVars->m_area[1];
1258  ELEMTYPEREAL diff = growth1 - growth0;
1259  if ( diff >= 0 )
1260  {
1261  group = 0;
1262  }
1263  else
1264  {
1265  group = 1;
1266  diff = -diff;
1267  }
1268 
1269  if ( diff > biggestDiff )
1270  {
1271  biggestDiff = diff;
1272  chosen = index;
1273  betterGroup = group;
1274  }
1275  else if ( qgsDoubleNear( diff, biggestDiff ) && ( a_parVars->m_count[group] < a_parVars->m_count[betterGroup] ) )
1276  {
1277  chosen = index;
1278  betterGroup = group;
1279  }
1280  }
1281  }
1282  Classify( chosen, betterGroup, a_parVars );
1283  }
1284 
1285  // If one group too full, put remaining rects in the other
1286  if ( ( a_parVars->m_count[0] + a_parVars->m_count[1] ) < a_parVars->m_total )
1287  {
1288  if ( a_parVars->m_count[0] >= a_parVars->m_total - a_parVars->m_minFill )
1289  {
1290  group = 1;
1291  }
1292  else
1293  {
1294  group = 0;
1295  }
1296  for ( int index = 0; index < a_parVars->m_total; ++index )
1297  {
1298  if ( !a_parVars->m_taken[index] )
1299  {
1300  Classify( index, group, a_parVars );
1301  }
1302  }
1303  }
1304 
1305  ASSERT( ( a_parVars->m_count[0] + a_parVars->m_count[1] ) == a_parVars->m_total );
1306  ASSERT( ( a_parVars->m_count[0] >= a_parVars->m_minFill ) &&
1307  ( a_parVars->m_count[1] >= a_parVars->m_minFill ) );
1308  }
1309 
1310 
1311 // Copy branches from the buffer into two nodes according to the partition.
1312  RTREE_TEMPLATE
1313  void RTREE_QUAL::LoadNodes( Node *a_nodeA, Node *a_nodeB, PartitionVars *a_parVars )
1314  {
1315  ASSERT( a_nodeA );
1316  ASSERT( a_nodeB );
1317  ASSERT( a_parVars );
1318 
1319  for ( int index = 0; index < a_parVars->m_total; ++index )
1320  {
1321  ASSERT( a_parVars->m_partition[index] == 0 || a_parVars->m_partition[index] == 1 );
1322 
1323  if ( a_parVars->m_partition[index] == 0 )
1324  {
1325  AddBranch( &a_parVars->m_branchBuf[index], a_nodeA, nullptr );
1326  }
1327  else if ( a_parVars->m_partition[index] == 1 )
1328  {
1329  AddBranch( &a_parVars->m_branchBuf[index], a_nodeB, nullptr );
1330  }
1331  }
1332  }
1333 
1334 
1335 // Initialize a PartitionVars structure.
1336  RTREE_TEMPLATE
1337  void RTREE_QUAL::InitParVars( PartitionVars *a_parVars, int a_maxRects, int a_minFill )
1338  {
1339  ASSERT( a_parVars );
1340 
1341  a_parVars->m_count[0] = a_parVars->m_count[1] = 0;
1342  a_parVars->m_area[0] = a_parVars->m_area[1] = static_cast< ELEMTYPEREAL >( 0 );
1343  a_parVars->m_total = a_maxRects;
1344  a_parVars->m_minFill = a_minFill;
1345  for ( int index = 0; index < a_maxRects; ++index )
1346  {
1347  a_parVars->m_taken[index] = false;
1348  a_parVars->m_partition[index] = -1;
1349  }
1350  }
1351 
1352 
1353  RTREE_TEMPLATE
1354  void RTREE_QUAL::PickSeeds( PartitionVars *a_parVars )
1355  {
1356  int seed0 = 0, seed1 = 0;
1357  ELEMTYPEREAL worst, waste;
1358  ELEMTYPEREAL area[MAXNODES + 1];
1359 
1360  for ( int index = 0; index < a_parVars->m_total; ++index )
1361  {
1362  area[index] = CalcRectVolume( &a_parVars->m_branchBuf[index].m_rect );
1363  }
1364 
1365  worst = -a_parVars->m_coverSplitArea - 1;
1366  for ( int indexA = 0; indexA < a_parVars->m_total - 1; ++indexA )
1367  {
1368  for ( int indexB = indexA + 1; indexB < a_parVars->m_total; ++indexB )
1369  {
1370  Rect oneRect = CombineRect( &a_parVars->m_branchBuf[indexA].m_rect, &a_parVars->m_branchBuf[indexB].m_rect );
1371  waste = CalcRectVolume( &oneRect ) - area[indexA] - area[indexB];
1372  if ( waste > worst )
1373  {
1374  worst = waste;
1375  seed0 = indexA;
1376  seed1 = indexB;
1377  }
1378  }
1379  }
1380  Classify( seed0, 0, a_parVars );
1381  Classify( seed1, 1, a_parVars );
1382  }
1383 
1384 
1385 // Put a branch in one of the groups.
1386  RTREE_TEMPLATE
1387  void RTREE_QUAL::Classify( int a_index, int a_group, PartitionVars *a_parVars )
1388  {
1389  ASSERT( a_parVars );
1390  ASSERT( !a_parVars->m_taken[a_index] );
1391 
1392  a_parVars->m_partition[a_index] = a_group;
1393  a_parVars->m_taken[a_index] = true;
1394 
1395  if ( a_parVars->m_count[a_group] == 0 )
1396  {
1397  a_parVars->m_cover[a_group] = a_parVars->m_branchBuf[a_index].m_rect;
1398  }
1399  else
1400  {
1401  a_parVars->m_cover[a_group] = CombineRect( &a_parVars->m_branchBuf[a_index].m_rect, &a_parVars->m_cover[a_group] );
1402  }
1403  a_parVars->m_area[a_group] = CalcRectVolume( &a_parVars->m_cover[a_group] );
1404  ++a_parVars->m_count[a_group];
1405  }
1406 
1407 
1408 // Delete a data rectangle from an index structure.
1409 // Pass in a pointer to a Rect, the tid of the record, ptr to ptr to root node.
1410 // Returns 1 if record not found, 0 if success.
1411 // RemoveRect provides for eliminating the root.
1412  RTREE_TEMPLATE
1413  bool RTREE_QUAL::RemoveRect( Rect *a_rect, const DATATYPE &a_id, Node **a_root )
1414  {
1415  ASSERT( a_rect && a_root );
1416  ASSERT( *a_root );
1417 
1418  Node *tempNode = nullptr;
1419  ListNode *reInsertList = nullptr;
1420 
1421  if ( !RemoveRectRec( a_rect, a_id, *a_root, &reInsertList ) )
1422  {
1423  // Found and deleted a data item
1424  // Reinsert any branches from eliminated nodes
1425  while ( reInsertList )
1426  {
1427  tempNode = reInsertList->m_node;
1428 
1429  for ( int index = 0; index < tempNode->m_count; ++index )
1430  {
1431  InsertRect( & ( tempNode->m_branch[index].m_rect ),
1432  tempNode->m_branch[index].m_data,
1433  a_root,
1434  tempNode->m_level );
1435  }
1436 
1437  ListNode *remLNode = reInsertList;
1438  reInsertList = reInsertList->m_next;
1439 
1440  FreeNode( remLNode->m_node );
1441  FreeListNode( remLNode );
1442  }
1443 
1444  // Check for redundant root (not leaf, 1 child) and eliminate
1445  if ( ( *a_root )->m_count == 1 && ( *a_root )->IsInternalNode() )
1446  {
1447  tempNode = ( *a_root )->m_branch[0].m_child;
1448 
1449  ASSERT( tempNode );
1450  FreeNode( *a_root );
1451  *a_root = tempNode;
1452  }
1453  return false;
1454  }
1455  else
1456  {
1457  return true;
1458  }
1459  }
1460 
1461 
1462 // Delete a rectangle from non-root part of an index structure.
1463 // Called by RemoveRect. Descends tree recursively,
1464 // merges branches on the way back up.
1465 // Returns 1 if record not found, 0 if success.
1466  RTREE_TEMPLATE
1467  bool RTREE_QUAL::RemoveRectRec( Rect *a_rect, const DATATYPE &a_id, Node *a_node, ListNode **a_listNode )
1468  {
1469  ASSERT( a_rect && a_node && a_listNode );
1470  ASSERT( a_node->m_level >= 0 );
1471 
1472  if ( a_node->IsInternalNode() ) // not a leaf node
1473  {
1474  for ( int index = 0; index < a_node->m_count; ++index )
1475  {
1476  if ( Overlap( a_rect, & ( a_node->m_branch[index].m_rect ) ) )
1477  {
1478  if ( !RemoveRectRec( a_rect, a_id, a_node->m_branch[index].m_child, a_listNode ) )
1479  {
1480  if ( a_node->m_branch[index].m_child->m_count >= MINNODES )
1481  {
1482  // child removed, just resize parent rect
1483  a_node->m_branch[index].m_rect = NodeCover( a_node->m_branch[index].m_child );
1484  }
1485  else
1486  {
1487  // child removed, not enough entries in node, eliminate node
1488  ReInsert( a_node->m_branch[index].m_child, a_listNode );
1489  DisconnectBranch( a_node, index ); // Must return after this call as count has changed
1490  }
1491  return false;
1492  }
1493  }
1494  }
1495  return true;
1496  }
1497  else // A leaf node
1498  {
1499  for ( int index = 0; index < a_node->m_count; ++index )
1500  {
1501  if ( a_node->m_branch[index].m_child == reinterpret_cast< Node * >( a_id ) )
1502  {
1503  DisconnectBranch( a_node, index ); // Must return after this call as count has changed
1504  return false;
1505  }
1506  }
1507  return true;
1508  }
1509  }
1510 
1511 
1512 // Decide whether two rectangles overlap.
1513  RTREE_TEMPLATE
1514  bool RTREE_QUAL::Overlap( Rect *a_rectA, Rect *a_rectB )
1515  {
1516  ASSERT( a_rectA && a_rectB );
1517 
1518  for ( int index = 0; index < NUMDIMS; ++index )
1519  {
1520  if ( a_rectA->m_min[index] > a_rectB->m_max[index] ||
1521  a_rectB->m_min[index] > a_rectA->m_max[index] )
1522  {
1523  return false;
1524  }
1525  }
1526  return true;
1527  }
1528 
1529 
1530 // Add a node to the reinsertion list. All its branches will later
1531 // be reinserted into the index structure.
1532  RTREE_TEMPLATE
1533  void RTREE_QUAL::ReInsert( Node *a_node, ListNode **a_listNode )
1534  {
1535  ListNode *newListNode = nullptr;
1536 
1537  newListNode = AllocListNode();
1538  newListNode->m_node = a_node;
1539  newListNode->m_next = *a_listNode;
1540  *a_listNode = newListNode;
1541  }
1542 
1543 
1544 // Search in an index tree or subtree for all data retangles that overlap the argument rectangle.
1545  RTREE_TEMPLATE
1546  bool RTREE_QUAL::Search( Node *a_node, Rect *a_rect, int &a_foundCount, bool ( *a_resultCallback )( DATATYPE a_data, void *a_context ), void *a_context )
1547  {
1548  ASSERT( a_node );
1549  ASSERT( a_node->m_level >= 0 );
1550  ASSERT( a_rect );
1551 
1552  if ( a_node->IsInternalNode() ) // This is an internal node in the tree
1553  {
1554  for ( int index = 0; index < a_node->m_count; ++index )
1555  {
1556  if ( Overlap( a_rect, &a_node->m_branch[index].m_rect ) )
1557  {
1558  if ( !Search( a_node->m_branch[index].m_child, a_rect, a_foundCount, a_resultCallback, a_context ) )
1559  {
1560  return false; // Don't continue searching
1561  }
1562  }
1563  }
1564  }
1565  else // This is a leaf node
1566  {
1567  for ( int index = 0; index < a_node->m_count; ++index )
1568  {
1569  if ( Overlap( a_rect, &a_node->m_branch[index].m_rect ) )
1570  {
1571  DATATYPE &id = a_node->m_branch[index].m_data;
1572 
1573  // NOTE: There are different ways to return results. Here's where to modify
1574  if ( a_resultCallback )
1575  {
1576  ++a_foundCount;
1577  if ( !a_resultCallback( id, a_context ) )
1578  {
1579  return false; // Don't continue searching
1580  }
1581  }
1582  }
1583  }
1584  }
1585 
1586  return true; // Continue searching
1587  }
1588 
1589 
1590 #undef RTREE_TEMPLATE
1591 #undef RTREE_QUAL
1592 
1593 }
1594 
1596 
1597 #endif //RTREE_H
1598 
bool qgsDoubleNear(double a, double b, double epsilon=4 *std::numeric_limits< double >::epsilon())
Compare two doubles (but allow some difference)
Definition: qgis.h:265
QgsMargins operator*(const QgsMargins &margins, double factor)
Returns a QgsMargins object that is formed by multiplying each component of the given margins by fact...
Definition: qgsmargins.h:242