← Back to branch summary

~ubuntu-branches/ubuntu/wily/qgis/wily

« back to all changes in this revision

Viewing changes to src/core/spatialindex/rtree/Index.cc

  • Committer: Bazaar Package Importer
  • Author(s): Johan Van de Wauw
  • Date: 2010-07-11 20:23:24 UTC
  • mfrom: (3.1.4 squeeze)
  • Revision ID: james.westby@ubuntu.com-20100711202324-5ktghxa7hracohmr
Tags: 1.4.0+12730-3ubuntu1
https://launchpad.net/bugs/540941
* Merge from Debian unstable (LP: #540941).
* Fix compilation issues with QT 4.7
* Add build-depends on libqt4-webkit-dev 

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/core/spatialindex/rtree/Index.cc
 
1
// Spatial Index Library
 
2
//
 
3
// Copyright (C) 2002 Navel Ltd.
 
4
//
 
5
// This library is free software; you can redistribute it and/or
 
6
// modify it under the terms of the GNU Lesser General Public
 
7
// License as published by the Free Software Foundation; either
 
8
// version 2.1 of the License, or (at your option) any later version.
 
9
//
 
10
// This library is distributed in the hope that it will be useful,
 
11
// but WITHOUT ANY WARRANTY; without even the implied warranty of
 
12
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 
13
// Lesser General Public License for more details.
 
14
//
 
15
// You should have received a copy of the GNU Lesser General Public
 
16
// License along with this library; if not, write to the Free Software
 
17
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 
18
//
 
19
//  Email:
 
20
//    mhadji@gmail.com
 
21
 
 
22
#include "../spatialindex/SpatialIndexImpl.h"
 
23
 
 
24
#include "RTree.h"
 
25
#include "Node.h"
 
26
#include "Leaf.h"
 
27
 
 
28
#include "Index.h"
 
29
 
 
30
using std::stack;
 
31
using std::vector;
 
32
using namespace SpatialIndex::RTree;
 
33
 
 
34
Index::~Index()
 
35
{
 
36
}
 
37
 
 
38
#ifdef _MSC_VER
 
39
// MSVC seems to find RTree* pTree ambiguous
 
40
Index::Index( SpatialIndex::RTree::RTree* pTree, long id, unsigned long level ) : Node( pTree, id, level, pTree->m_indexCapacity )
 
41
#else
 
42
Index::Index( RTree* pTree, long id, unsigned long level ) : Node( pTree, id, level, pTree->m_indexCapacity )
 
43
#endif//_MSC_VER
 
44
{
 
45
}
 
46
 
 
47
NodePtr Index::chooseSubtree( const Region& mbr, unsigned long insertionLevel, std::stack<long>& pathBuffer )
 
48
{
 
49
  if ( m_level == insertionLevel ) return NodePtr( this, &( m_pTree->m_indexPool ) );
 
50
 
 
51
  pathBuffer.push( m_identifier );
 
52
 
 
53
  unsigned long child = 0;
 
54
 
 
55
  switch ( m_pTree->m_treeVariant )
 
56
  {
 
57
    case RV_LINEAR:
 
58
    case RV_QUADRATIC:
 
59
      child = findLeastEnlargement( mbr );
 
60
      break;
 
61
    case RV_RSTAR:
 
62
      if ( m_level == 1 )
 
63
      {
 
64
        // if this node points to leaves...
 
65
        child = findLeastOverlap( mbr );
 
66
      }
 
67
      else
 
68
      {
 
69
        child = findLeastEnlargement( mbr );
 
70
      }
 
71
      break;
 
72
    default:
 
73
      throw Tools::NotSupportedException( "Index::chooseSubtree: Tree variant not supported." );
 
74
  }
 
75
  assert( child >= 0 );
 
76
 
 
77
  NodePtr n = m_pTree->readNode( m_pIdentifier[child] );
 
78
  NodePtr ret = n->chooseSubtree( mbr, insertionLevel, pathBuffer );
 
79
  assert( n.unique() );
 
80
  if ( ret.get() == n.get() ) n.relinquish();
 
81
 
 
82
  return ret;
 
83
}
 
84
 
 
85
NodePtr Index::findLeaf( const Region& mbr, long id, std::stack<long>& pathBuffer )
 
86
{
 
87
  pathBuffer.push( m_identifier );
 
88
 
 
89
  for ( unsigned long cChild = 0; cChild < m_children; cChild++ )
 
90
  {
 
91
    if ( m_ptrMBR[cChild]->containsRegion( mbr ) )
 
92
    {
 
93
      NodePtr n = m_pTree->readNode( m_pIdentifier[cChild] );
 
94
      NodePtr l = n->findLeaf( mbr, id, pathBuffer );
 
95
      if ( n.get() == l.get() ) n.relinquish();
 
96
      if ( l.get() != 0 ) return l;
 
97
    }
 
98
  }
 
99
 
 
100
  pathBuffer.pop();
 
101
 
 
102
  return NodePtr();
 
103
}
 
104
 
 
105
void Index::split( unsigned long dataLength, byte* pData, Region& mbr, long id, NodePtr& ptrLeft, NodePtr& ptrRight )
 
106
{
 
107
  m_pTree->m_stats.m_splits++;
 
108
 
 
109
  vector<long> g1, g2;
 
110
 
 
111
  switch ( m_pTree->m_treeVariant )
 
112
  {
 
113
    case RV_LINEAR:
 
114
    case RV_QUADRATIC:
 
115
      rtreeSplit( dataLength, pData, mbr, id, g1, g2 );
 
116
      break;
 
117
    case RV_RSTAR:
 
118
      rstarSplit( dataLength, pData, mbr, id, g1, g2 );
 
119
      break;
 
120
    default:
 
121
      throw Tools::NotSupportedException( "Index::split: Tree variant not supported." );
 
122
  }
 
123
 
 
124
  ptrLeft = m_pTree->m_indexPool.acquire();
 
125
  ptrRight = m_pTree->m_indexPool.acquire();
 
126
 
 
127
  if ( ptrLeft.get() == 0 ) ptrLeft = NodePtr( new Index( m_pTree, m_identifier, m_level ), &( m_pTree->m_indexPool ) );
 
128
  if ( ptrRight.get() == 0 ) ptrRight = NodePtr( new Index( m_pTree, -1, m_level ), &( m_pTree->m_indexPool ) );
 
129
 
 
130
  ptrLeft->m_nodeMBR = m_pTree->m_infiniteRegion;
 
131
  ptrRight->m_nodeMBR = m_pTree->m_infiniteRegion;
 
132
 
 
133
  unsigned long cIndex;
 
134
 
 
135
  for ( cIndex = 0; cIndex < g1.size(); cIndex++ )
 
136
  {
 
137
    ptrLeft->insertEntry( 0, 0, *( m_ptrMBR[g1[cIndex]] ), m_pIdentifier[g1[cIndex]] );
 
138
  }
 
139
 
 
140
  for ( cIndex = 0; cIndex < g2.size(); cIndex++ )
 
141
  {
 
142
    ptrRight->insertEntry( 0, 0, *( m_ptrMBR[g2[cIndex]] ), m_pIdentifier[g2[cIndex]] );
 
143
  }
 
144
}
 
145
 
 
146
long Index::findLeastEnlargement( const Region& r ) const
 
147
{
 
148
  double area = std::numeric_limits<double>::max();
 
149
  long best = -1;
 
150
 
 
151
  RegionPtr t = m_pTree->m_regionPool.acquire();
 
152
 
 
153
  for ( unsigned long cChild = 0; cChild < m_children; cChild++ )
 
154
  {
 
155
    m_ptrMBR[cChild]->getCombinedRegion( *t, r );
 
156
 
 
157
    double a = m_ptrMBR[cChild]->getArea();
 
158
    double enl = t->getArea() - a;
 
159
 
 
160
    if ( enl < area )
 
161
    {
 
162
      area = enl;
 
163
      best = cChild;
 
164
    }
 
165
    else if ( enl == area )
 
166
    {
 
167
      // this will rarely happen, so compute best area on the fly only
 
168
      // when necessary.
 
169
      if ( a < m_ptrMBR[best]->getArea() ) best = cChild;
 
170
    }
 
171
  }
 
172
 
 
173
  return best;
 
174
}
 
175
 
 
176
long Index::findLeastOverlap( const Region& r ) const
 
177
{
 
178
  OverlapEntry** entries = new OverlapEntry*[m_children];
 
179
 
 
180
  double leastOverlap = std::numeric_limits<double>::max();
 
181
  double me = std::numeric_limits<double>::max();
 
182
  OverlapEntry* best = 0;
 
183
 
 
184
  // find combined region and enlargement of every entry and store it.
 
185
  for ( unsigned long cChild = 0; cChild < m_children; cChild++ )
 
186
  {
 
187
    try
 
188
    {
 
189
      entries[cChild] = new OverlapEntry();
 
190
    }
 
191
    catch ( ... )
 
192
    {
 
193
      for ( unsigned long i = 0; i < cChild; i++ ) delete entries[i];
 
194
      delete[] entries;
 
195
      throw;
 
196
    }
 
197
 
 
198
    entries[cChild]->m_id = cChild;
 
199
    entries[cChild]->m_original = m_ptrMBR[cChild];
 
200
    entries[cChild]->m_combined = m_pTree->m_regionPool.acquire();
 
201
    m_ptrMBR[cChild]->getCombinedRegion( *( entries[cChild]->m_combined ), r );
 
202
    entries[cChild]->m_oa = entries[cChild]->m_original->getArea();
 
203
    entries[cChild]->m_ca = entries[cChild]->m_combined->getArea();
 
204
    entries[cChild]->m_enlargement = entries[cChild]->m_ca - entries[cChild]->m_oa;
 
205
 
 
206
    if ( entries[cChild]->m_enlargement < me )
 
207
    {
 
208
      me = entries[cChild]->m_enlargement;
 
209
      best = entries[cChild];
 
210
    }
 
211
    else if ( entries[cChild]->m_enlargement == me && entries[cChild]->m_oa < best->m_oa )
 
212
    {
 
213
      best = entries[cChild];
 
214
    }
 
215
  }
 
216
 
 
217
  if ( me < -std::numeric_limits<double>::epsilon() || me > std::numeric_limits<double>::epsilon() )
 
218
  {
 
219
    unsigned long cIterations;
 
220
 
 
221
    if ( m_children > m_pTree->m_nearMinimumOverlapFactor )
 
222
    {
 
223
      // sort entries in increasing order of enlargement.
 
224
      ::qsort( entries, m_children,
 
225
               sizeof( OverlapEntry* ),
 
226
               OverlapEntry::compareEntries );
 
227
      assert( entries[0]->m_enlargement <= entries[m_children - 1]->m_enlargement );
 
228
 
 
229
      cIterations = m_pTree->m_nearMinimumOverlapFactor;
 
230
    }
 
231
    else
 
232
    {
 
233
      cIterations = m_children;
 
234
    }
 
235
 
 
236
    // calculate overlap of most important original entries (near minimum overlap cost).
 
237
    for ( unsigned long cIndex = 0; cIndex < cIterations; cIndex++ )
 
238
    {
 
239
      double dif = 0.0;
 
240
      OverlapEntry* e = entries[cIndex];
 
241
 
 
242
      for ( unsigned long cChild = 0; cChild < m_children; cChild++ )
 
243
      {
 
244
        if ( e->m_id != cChild )
 
245
        {
 
246
          double f = e->m_combined->getIntersectingArea( *( m_ptrMBR[cChild] ) );
 
247
          if ( f != 0.0 ) dif += f - e->m_original->getIntersectingArea( *( m_ptrMBR[cChild] ) );
 
248
        }
 
249
      } // for (cChild)
 
250
 
 
251
      if ( dif < leastOverlap )
 
252
      {
 
253
        leastOverlap = dif;
 
254
        best = entries[cIndex];
 
255
      }
 
256
      else if ( dif == leastOverlap )
 
257
      {
 
258
        if ( e->m_enlargement == best->m_enlargement )
 
259
        {
 
260
          // keep the one with least area.
 
261
          if ( e->m_original->getArea() < best->m_original->getArea() ) best = entries[cIndex];
 
262
        }
 
263
        else
 
264
        {
 
265
          // keep the one with least enlargement.
 
266
          if ( e->m_enlargement < best->m_enlargement ) best = entries[cIndex];
 
267
        }
 
268
      }
 
269
    } // for (cIndex)
 
270
  }
 
271
 
 
272
  long ret = best->m_id;
 
273
 
 
274
  for ( unsigned long cChild = 0; cChild < m_children; cChild++ )
 
275
  {
 
276
    delete entries[cChild];
 
277
  }
 
278
  delete[] entries;
 
279
 
 
280
  return ret;
 
281
}
 
282
 
 
283
void Index::adjustTree( Node* n, std::stack<long>& pathBuffer )
 
284
{
 
285
  m_pTree->m_stats.m_adjustments++;
 
286
 
 
287
  // find entry pointing to old node;
 
288
  unsigned long child;
 
289
  for ( child = 0; child < m_children; child++ )
 
290
  {
 
291
    if ( m_pIdentifier[child] == n->m_identifier ) break;
 
292
  }
 
293
 
 
294
  // MBR needs recalculation if either:
 
295
  //   1. the NEW child MBR is not contained.
 
296
  //   2. the OLD child MBR is touching.
 
297
  bool bContained = m_nodeMBR.containsRegion( n->m_nodeMBR );
 
298
  bool bTouches = m_nodeMBR.touchesRegion( *( m_ptrMBR[child] ) );
 
299
  bool bRecompute = ( ! bContained || ( bTouches && m_pTree->m_bTightMBRs ) );
 
300
 
 
301
  *( m_ptrMBR[child] ) = n->m_nodeMBR;
 
302
 
 
303
  if ( bRecompute )
 
304
  {
 
305
    for ( unsigned long cDim = 0; cDim < m_nodeMBR.m_dimension; cDim++ )
 
306
    {
 
307
      m_nodeMBR.m_pLow[cDim] = std::numeric_limits<double>::max();
 
308
      m_nodeMBR.m_pHigh[cDim] = -std::numeric_limits<double>::max();
 
309
 
 
310
      for ( unsigned long cChild = 0; cChild < m_children; cChild++ )
 
311
      {
 
312
        m_nodeMBR.m_pLow[cDim] = std::min( m_nodeMBR.m_pLow[cDim], m_ptrMBR[cChild]->m_pLow[cDim] );
 
313
        m_nodeMBR.m_pHigh[cDim] = std::max( m_nodeMBR.m_pHigh[cDim], m_ptrMBR[cChild]->m_pHigh[cDim] );
 
314
      }
 
315
    }
 
316
  }
 
317
 
 
318
  m_pTree->writeNode( this );
 
319
 
 
320
  if ( bRecompute && ( ! pathBuffer.empty() ) )
 
321
  {
 
322
    long cParent = pathBuffer.top(); pathBuffer.pop();
 
323
    NodePtr ptrN = m_pTree->readNode( cParent );
 
324
    Index* p = static_cast<Index*>( ptrN.get() );
 
325
    p->adjustTree( this, pathBuffer );
 
326
  }
 
327
}
 
328
 
 
329
void Index::adjustTree( Node* n1, Node* n2, std::stack<long>& pathBuffer, byte* overflowTable )
 
330
{
 
331
  m_pTree->m_stats.m_adjustments++;
 
332
 
 
333
  // find entry pointing to old node;
 
334
  unsigned long child;
 
335
  for ( child = 0; child < m_children; child++ )
 
336
  {
 
337
    if ( m_pIdentifier[child] == n1->m_identifier ) break;
 
338
  }
 
339
 
 
340
  // MBR needs recalculation if either:
 
341
  //   1. the NEW child MBR is not contained.
 
342
  //   2. the OLD child MBR is touching.
 
343
  bool bContained = m_nodeMBR.containsRegion( n1->m_nodeMBR );
 
344
  bool bTouches = m_nodeMBR.touchesRegion( *( m_ptrMBR[child] ) );
 
345
  bool bRecompute = ( ! bContained || ( bTouches && m_pTree->m_bTightMBRs ) );
 
346
 
 
347
  *( m_ptrMBR[child] ) = n1->m_nodeMBR;
 
348
 
 
349
  if ( bRecompute )
 
350
  {
 
351
    for ( unsigned long cDim = 0; cDim < m_nodeMBR.m_dimension; cDim++ )
 
352
    {
 
353
      m_nodeMBR.m_pLow[cDim] = std::numeric_limits<double>::max();
 
354
      m_nodeMBR.m_pHigh[cDim] = -std::numeric_limits<double>::max();
 
355
 
 
356
      for ( unsigned long cChild = 0; cChild < m_children; cChild++ )
 
357
      {
 
358
        m_nodeMBR.m_pLow[cDim] = std::min( m_nodeMBR.m_pLow[cDim], m_ptrMBR[cChild]->m_pLow[cDim] );
 
359
        m_nodeMBR.m_pHigh[cDim] = std::max( m_nodeMBR.m_pHigh[cDim], m_ptrMBR[cChild]->m_pHigh[cDim] );
 
360
      }
 
361
    }
 
362
  }
 
363
 
 
364
  // No write necessary here. insertData will write the node if needed.
 
365
  //m_pTree->writeNode(this);
 
366
 
 
367
  bool bAdjusted = insertData( 0, 0, n2->m_nodeMBR, n2->m_identifier, pathBuffer, overflowTable );
 
368
 
 
369
  // if n2 is contained in the node and there was no split or reinsert,
 
370
  // we need to adjust only if recalculation took place.
 
371
  // In all other cases insertData above took care of adjustment.
 
372
  if (( ! bAdjusted ) && bRecompute && ( ! pathBuffer.empty() ) )
 
373
  {
 
374
    long cParent = pathBuffer.top(); pathBuffer.pop();
 
375
    NodePtr ptrN = m_pTree->readNode( cParent );
 
376
    Index* p = static_cast<Index*>( ptrN.get() );
 
377
    p->adjustTree( this, pathBuffer );
 
378
  }
 
379
}