QGIS API Documentation  2.7.0-Master
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
qgsspatialindex.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgsspatialindex.cpp - wrapper class for spatial index library
3  ----------------------
4  begin : December 2006
5  copyright : (C) 2006 by Martin Dobias
6  email : wonder.sk at gmail dot com
7  ***************************************************************************
8  * *
9  * This program is free software; you can redistribute it and/or modify *
10  * it under the terms of the GNU General Public License as published by *
11  * the Free Software Foundation; either version 2 of the License, or *
12  * (at your option) any later version. *
13  * *
14  ***************************************************************************/
15 
16 #include "qgsspatialindex.h"
17 
18 #include "qgsgeometry.h"
19 #include "qgsfeature.h"
20 #include "qgsfeatureiterator.h"
21 #include "qgsrectangle.h"
22 #include "qgslogger.h"
23 
24 #include "SpatialIndex.h"
25 
26 using namespace SpatialIndex;
27 
28 
29 
30 // custom visitor that adds found features to list
31 class QgisVisitor : public SpatialIndex::IVisitor
32 {
33  public:
34  QgisVisitor( QList<QgsFeatureId> & list )
35  : mList( list ) {}
36 
37  void visitNode( const INode& n )
38  { Q_UNUSED( n ); }
39 
40  void visitData( const IData& d )
41  {
42  mList.append( d.getIdentifier() );
43  }
44 
45  void visitData( std::vector<const IData*>& v )
46  { Q_UNUSED( v ); }
47 
48  private:
49  QList<QgsFeatureId>& mList;
50 };
51 
52 class QgsSpatialIndexCopyVisitor : public SpatialIndex::IVisitor
53 {
54  public:
55  QgsSpatialIndexCopyVisitor( SpatialIndex::ISpatialIndex* newIndex )
56  : mNewIndex( newIndex ) {}
57 
58  void visitNode( const INode& n )
59  { Q_UNUSED( n ); }
60 
61  void visitData( const IData& d )
62  {
63  SpatialIndex::IShape* shape;
64  d.getShape( &shape );
65  mNewIndex->insertData( 0, 0, *shape, d.getIdentifier() );
66  delete shape;
67  }
68 
69  void visitData( std::vector<const IData*>& v )
70  { Q_UNUSED( v ); }
71 
72  private:
73  SpatialIndex::ISpatialIndex* mNewIndex;
74 };
75 
76 
78 class QgsFeatureIteratorDataStream : public IDataStream
79 {
80  public:
82  QgsFeatureIteratorDataStream( const QgsFeatureIterator& fi ) : mFi( fi ), mNextData( 0 )
83  {
84  readNextEntry();
85  }
86 
88  {
89  delete mNextData;
90  }
91 
93  virtual IData* getNext()
94  {
95  RTree::Data* ret = mNextData;
96  mNextData = 0;
97  readNextEntry();
98  return ret;
99  }
100 
102  virtual bool hasNext() { return mNextData != 0; }
103 
105  virtual uint32_t size() { Q_ASSERT( 0 && "not available" ); return 0; }
106 
108  virtual void rewind() { Q_ASSERT( 0 && "not available" ); }
109 
110  protected:
112  {
113  QgsFeature f;
114  SpatialIndex::Region r;
115  QgsFeatureId id;
116  while ( mFi.nextFeature( f ) )
117  {
118  if ( QgsSpatialIndex::featureInfo( f, r, id ) )
119  {
120  mNextData = new RTree::Data( 0, 0, r, id );
121  return;
122  }
123  }
124  }
125 
126  private:
127  QgsFeatureIterator mFi;
128  RTree::Data* mNextData;
129 };
130 
131 
133 class QgsSpatialIndexData : public QSharedData
134 {
135  public:
137  {
138  initTree();
139  }
140 
142  {
143  QgsFeatureIteratorDataStream fids( fi );
144  initTree( &fids );
145  }
146 
148  : QSharedData( other )
149  {
150  initTree();
151 
152  // copy R-tree data one by one (is there a faster way??)
153  double low[] = { DBL_MIN, DBL_MIN };
154  double high[] = { DBL_MAX, DBL_MAX };
155  SpatialIndex::Region query( low, high, 2 );
156  QgsSpatialIndexCopyVisitor visitor( mRTree );
157  other.mRTree->intersectsWithQuery( query, visitor );
158  }
159 
161  {
162  delete mRTree;
163  delete mStorage;
164  }
165 
166  void initTree( IDataStream* inputStream = 0 )
167  {
168  // for now only memory manager
169  mStorage = StorageManager::createNewMemoryStorageManager();
170 
171  // R-Tree parameters
172  double fillFactor = 0.7;
173  unsigned long indexCapacity = 10;
174  unsigned long leafCapacity = 10;
175  unsigned long dimension = 2;
176  RTree::RTreeVariant variant = RTree::RV_RSTAR;
177 
178  // create R-tree
179  SpatialIndex::id_type indexId;
180 
181  if ( inputStream )
182  mRTree = RTree::createAndBulkLoadNewRTree( RTree::BLM_STR, *inputStream, *mStorage, fillFactor, indexCapacity,
183  leafCapacity, dimension, variant, indexId );
184  else
185  mRTree = RTree::createNewRTree( *mStorage, fillFactor, indexCapacity,
186  leafCapacity, dimension, variant, indexId );
187  }
188 
190  SpatialIndex::IStorageManager* mStorage;
191 
193  SpatialIndex::ISpatialIndex* mRTree;
194 };
195 
196 // -------------------------------------------------------------------------
197 
198 
200 {
201  d = new QgsSpatialIndexData;
202 }
203 
205 {
206  d = new QgsSpatialIndexData( fi );
207 }
208 
210  : d( other.d )
211 {
212 }
213 
215 {
216 }
217 
219 {
220  if ( this != &other )
221  d = other.d;
222  return *this;
223 }
224 
226 {
227  double pt1[2], pt2[2];
228  pt1[0] = rect.xMinimum();
229  pt1[1] = rect.yMinimum();
230  pt2[0] = rect.xMaximum();
231  pt2[1] = rect.yMaximum();
232  return Region( pt1, pt2, 2 );
233 }
234 
235 bool QgsSpatialIndex::featureInfo( const QgsFeature& f, SpatialIndex::Region& r, QgsFeatureId &id )
236 {
237  QgsGeometry *g = f.geometry();
238  if ( !g )
239  return false;
240 
241  id = f.id();
242  r = rectToRegion( g->boundingBox() );
243  return true;
244 }
245 
246 
248 {
249  Region r;
250  QgsFeatureId id;
251  if ( !featureInfo( f, r, id ) )
252  return false;
253 
254  // TODO: handle possible exceptions correctly
255  try
256  {
257  d->mRTree->insertData( 0, 0, r, FID_TO_NUMBER( id ) );
258  return true;
259  }
260  catch ( Tools::Exception &e )
261  {
262  Q_UNUSED( e );
263  QgsDebugMsg( QString( "Tools::Exception caught: " ).arg( e.what().c_str() ) );
264  }
265  catch ( const std::exception &e )
266  {
267  Q_UNUSED( e );
268  QgsDebugMsg( QString( "std::exception caught: " ).arg( e.what() ) );
269  }
270  catch ( ... )
271  {
272  QgsDebugMsg( "unknown spatial index exception caught" );
273  }
274 
275  return false;
276 }
277 
279 {
280  Region r;
281  QgsFeatureId id;
282  if ( !featureInfo( f, r, id ) )
283  return false;
284 
285  // TODO: handle exceptions
286  return d->mRTree->deleteData( r, FID_TO_NUMBER( id ) );
287 }
288 
289 QList<QgsFeatureId> QgsSpatialIndex::intersects( QgsRectangle rect ) const
290 {
291  QList<QgsFeatureId> list;
292  QgisVisitor visitor( list );
293 
294  Region r = rectToRegion( rect );
295 
296  d->mRTree->intersectsWithQuery( r, visitor );
297 
298  return list;
299 }
300 
301 QList<QgsFeatureId> QgsSpatialIndex::nearestNeighbor( QgsPoint point, int neighbors ) const
302 {
303  QList<QgsFeatureId> list;
304  QgisVisitor visitor( list );
305 
306  double pt[2];
307  pt[0] = point.x();
308  pt[1] = point.y();
309  Point p( pt, 2 );
310 
311  d->mRTree->nearestNeighborQuery( neighbors, p, visitor );
312 
313  return list;
314 }
315 
316 QAtomicInt QgsSpatialIndex::refs() const
317 {
318  return d->ref;
319 }
QgsFeatureId id() const
Get the feature id for this feature.
Definition: qgsfeature.cpp:100
Wrapper for iterator of features from vector data provider or vector layer.
A rectangle specified with double values.
Definition: qgsrectangle.h:35
SpatialIndex::IStorageManager * mStorage
storage manager
QgsSpatialIndex & operator=(const QgsSpatialIndex &other)
implement assignment operator
static SpatialIndex::Region rectToRegion(QgsRectangle rect)
SpatialIndex::ISpatialIndex * mRTree
R-tree containing spatial index.
void initTree(IDataStream *inputStream=0)
bool deleteFeature(const QgsFeature &f)
remove feature from index
double yMaximum() const
Get the y maximum value (top side of rectangle)
Definition: qgsrectangle.h:188
#define QgsDebugMsg(str)
Definition: qgslogger.h:33
QgsGeometry * geometry() const
Get the geometry object associated with this feature.
Definition: qgsfeature.cpp:112
QgsSpatialIndexData(const QgsFeatureIterator &fi)
void visitNode(const INode &n)
The feature class encapsulates a single feature including its id, geometry and a list of field/values...
Definition: qgsfeature.h:113
QgsSpatialIndexCopyVisitor(SpatialIndex::ISpatialIndex *newIndex)
virtual void rewind()
sets the stream pointer to the first entry, if possible.
QList< QgsFeatureId > intersects(QgsRectangle rect) const
returns features that intersect the specified rectangle
double x() const
Definition: qgspoint.h:126
virtual bool hasNext()
returns true if there are more items in the stream.
QgsSpatialIndexData(const QgsSpatialIndexData &other)
QgisVisitor(QList< QgsFeatureId > &list)
QgsFeatureIteratorDataStream(const QgsFeatureIterator &fi)
constructor - needs to load all data to a vector for later access when bulk loading ...
double yMinimum() const
Get the y minimum value (bottom side of rectangle)
Definition: qgsrectangle.h:193
double xMaximum() const
Get the x maximum value (right side of rectangle)
Definition: qgsrectangle.h:178
#define FID_TO_NUMBER(fid)
Definition: qgsfeature.h:82
Data of spatial index that may be implicitly shared.
QList< QgsFeatureId > nearestNeighbor(QgsPoint point, int neighbors) const
returns nearest neighbors (their count is specified by second parameter)
QgsSpatialIndex()
constructor - creates R-tree
Utility class for bulk loading of R-trees.
void visitData(std::vector< const IData * > &v)
A class to represent a point.
Definition: qgspoint.h:63
void visitData(const IData &d)
void visitData(const IData &d)
QAtomicInt refs() const
get reference count - just for debugging!
QgsRectangle boundingBox()
Returns the bounding box of this feature.
bool insertFeature(const QgsFeature &f)
add feature to index
~QgsSpatialIndex()
destructor finalizes work with spatial index
virtual uint32_t size()
returns the total number of entries available in the stream.
void visitData(std::vector< const IData * > &v)
virtual IData * getNext()
returns a pointer to the next entry in the stream or 0 at the end of the stream.
qint64 QgsFeatureId
Definition: qgsfeature.h:30
double y() const
Definition: qgspoint.h:134
void visitNode(const INode &n)
static bool featureInfo(const QgsFeature &f, SpatialIndex::Region &r, QgsFeatureId &id)
double xMinimum() const
Get the x minimum value (left side of rectangle)
Definition: qgsrectangle.h:183